3
$\begingroup$

I am trying to detect outliers/noise as indicated on the diagram below from sensor data. Can anyone advice how to go about it? I can only do this in python, so are there libraries in python that I can leverage? Is there an example that can be given. Sample data am using has timestamps and the value. Also how can I also get the frequency of the noise/outlier occurence? enter image description here

$\endgroup$
2
$\begingroup$

There are several ways you can do this, however this is highly specific to the data you are analyzing. I will describe some of very general (not domain/data specific) approaches below:

One possible approach I can think of that'd deal with this sudden, HF noise is to run a window with desired length through the signal. In every iteration, you'd take several samples for this window (i.e. for t=0 to t=10, then t=1 to t=11 and so on) and measure the rate of change within the window. If the rate of change is below some threshold you defined, you'd draw a line from the last point where the rate of change is below the threshold to the beginning on the next window where the rate of change is below the threshold again. This might or might not work, however I described the more general and popular ways below:

  1. One of the simplest and most naive ways would be to apply a high order low pass filter that would remove higher frequencies (sudden changes are high frequency peaks in Fourier spectrum) - this would simply smooth out your signal but it would also alter it, so it could be no longer useful. Obviously you'd have to set up the cutoff frequency properly to get a satisfactory result. If the HF peaks have some specific, known frequency, you can also apply a band stop filter (or a few of them, pretty much like an audio equalizer!). It's pretty trivial to implement these filters as IIR and a bit more complex if you need FIR. You can find a lot of information about how to implement them in the Internet and there are most likely ready scripts for Python that readily implement them.

  2. Another - and perhaps a better approach - would be to use Kalman filter, but again, it's not obvious to tune and you'd need some more understanding of your signal and the Kalman filter to implement it succesfuly.

  3. Yet another (more sophisticated) approach is to use sample-based noise removal techniques, however the signal you presented doesn't seem to have just a single feature that could be processed with such the filter, so you'd have to apply it several times.

If you need help understanding/implementing any of the solutions described below, feel free to ask!

Regarding the question about how to measure the noise frequency - you can manually select the time range in which the noise occurs and draw a FFT spectrum / Histogram of it and find the peak value. That's the most prevalent frequency in your noise. As far as I can see, your signal is very smooth and thus low frequency, so it won't have much influence on this kind of measurement. If you don't feel like doing this in Python, you can always export it to a text format acceptable by Matlab or Scilab (free) and these programs have ready tools for doing this kind of analysis.

$\endgroup$
  • $\begingroup$ Thank you so much for your detailed response and I would like to understand more and try to implement all and see the difference in all. I have never worked with IIR filters before, do you mind giving more on how to implement 1. I would like to get started with the data set that I have and see what I get. $\endgroup$ – user5559 Mar 12 '16 at 7:09
3
$\begingroup$

Regarding your next question, I will try to give you some long description explainig a general way of transforming filters defined by their transfer function into discrete domain. If you are not interested in the details, you can skip to the bottom, where I wrote a C code (tested).

The transfer function of a first-order low pass filter (with unit gain) is defined as follows: $$K(s) = \frac{1}{1+sT}$$ And its angular cutoff frequency defined as $$\omega_c = \frac{1}{T}$$ Therefore the cutoff frequency expressed in Hertz becomes $$f_c=\frac{1}{2 \pi T}$$ $$T=\frac{1}{2 \pi f_c}$$ Whatever goes above this frequency, will be cut off with gain equal to -3dB/decade for the first order low pass filter.

Either using bilinear transform or by expressing the gain in form of differential equations and then discretizing them with Euler scheme you can obtain a discrete domain equation that expresses the gain. I prefer bilinear transform as it is more straightforward. In order to apply bilinear transform we perform the following substitution: $$ s = \frac{1}{T_s} \frac{1-z^{-1}}{1+z^{-1}}$$ Where Ts (not to be mistaken with T) is sampling period, so in your case it is the number of seconds between each sample. I.e. for 44100Hz sampling frequency, Ts becomes 1/44100 seconds.

Now, after performing this substitution and simplifying terms, the discrete time transfer function (Z transform) of your filter becomes: $$\frac{1+z^{-1}}{(1+2T_s/T) + (1-2T_s/T)z^{-1}}$$ Now lets recall the definition of Z transform. In general form it is expressed as the following: $$ H(z) = \frac{b_0 + b_1z^{-1}+ b_2z^{-2} + ... + b_Mz^{-M}}{1 + a_1z^{-1}+ a_2z^{-2} + ... + a_Nz^{-N}} = \frac{Y(z)}{X(z)}$$ Where bare numerator coefficients from 1 to M (highest power) and a are denominator coefficients.

Now, a difference equation can be derived in the following way: $$y(n) = b_0x(n) + b_1x(n-1) + b_2x(n-2) + ... + b_Mx(n-M) - a_1y(n-1) - a_2y(n-2) - ... - a_Ny(n-N) $$ Where n means n-th sample, y is your output and x is your input. After applying the above, the simplified form of difference equation for your filter becomes: $$y(n) = (1-\alpha)y(n-1) + \alpha x(n)$$ Where $$ \alpha = \frac{T_s}{T + T_s}$$

This means that n-th sample (the current one) in filtered signal, called y(n) here, equals one minus alpha times the previous sample in filtered signal plus alpha times n-th sample in your source (noised) signal. In C language:

double lowpass (double x, double y_prev, double alpha)
{
    return (1.0-alpha)*y_prev + alpha*x;
}

int main()
{
    double x[10] = {1,2,3,4,5,6,7,8,9,10}; //signal with noise
    double y[10]; //filtered signal (will be generated)
    double fc = 500.0; //500Hz cutoff frequency
    double Ts = 1.0/44100.0; //for example
    double T  = 1.0/(2.0*3.1415*fc);
    double alpha = Ts/(T+Ts);

    y[0] = x[0]; //initialize first y sample as equal to first x sample

    printf("Input\tOutput\n");
    int i;
    for(i=1; i<10; ++i)
    {
        y[i] = lowpass(x[i], y[i-1], alpha);
        printf("%.4f\t%.4f\n", x[i], y[i]);
    }
}

Without complicating this code, for higher order filtering, you can swap y with x once it is filtered and repeat what's in the loop. It's slow and naive, but otherwise, you'd have to take a power of the transfer function (i.e. K(s)^2 for second order filter) and derive the discrete equation again. Thankfully with some practice, you will see a pattern in this and (even better!) you can find ready derivations of higher order low pass filters in the Internet :)

Sources:
https://en.wikipedia.org/wiki/Low-pass_filter
https://en.wikipedia.org/wiki/Bilinear_transform https://ccrma.stanford.edu/~jos/fp/Z_Transform_Difference_Equations.html

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.