Skip to main content

probability or statistics - How to remove outliers from data


I have a data with noise which some times includes significant outliers. The position of the outliers are random.


For example:


data1 = Table[PDF[NormalDistribution[3.5, .8], i], {i, -5, 15, .01}] +
RandomReal[{100, 500}];
noise = RandomReal /@ RandomReal[{-0.2, .2}, Length[data1]];
data2 = data1 + noise;
n = RandomInteger[{1, Length[data2]}, RandomInteger[{2, 10}]];
data2[[n]] = data2[[n]]*1.01;
ListPlot[{data2}, PlotRange -> All]


enter image description here


One solution is to use the average of the data but because the position of the outliers are random, the non-noisy data is hard to extract. The whole level of the data is random which means I can not use fixed reference to check and remove outliers.


Any idea how to remove these points using Mathematica?


Thanks.



Answer



I will give you two similar methods. But, I will rewrite one of the comments above just to make sure it is read.




You've been given some fine answers, but be absolutely sure that removing the outliers is Doing The Right Thing™. You might want to consider "robust" methods that can deal with the presence of outliers. – Guess who it is.





Simple Gaussian Threshold


The simplest way is to remove the moving mean of the data, then compute its standard deviation ($\sigma$), then pick a level at which you want to reject the data, say at 1%, so you can remove any points that vary more than $ 3\times \sigma$ . If you know how the data is distributed about its mean values, then you can pick a different method. You can also remove the median since that would be less sensitive to the distribution.


SeedRandom[1245]; 
data1 =
Table[PDF[NormalDistribution[3.5, .8], i], {i, -5, 15, .01}] +
RandomReal[{100, 500}];
noise = RandomReal /@ RandomReal[{-0.2, .2}, Length[data1]];
data2 = data1 + noise;

n = RandomInteger[{1, Length[data2]}, RandomInteger[{2, 10}]];
data2[[n]] = data2[[n]]*1.01;
ListPlot[{data2}, PlotRange -> All]

enter image description here


We have about 8 outliers. We compute the moving average,


movingAvg=ArrayPad[MovingAverage[data2, 5],{5-1,0},"Fixed"]

Here we subtract the moving mean,


subtractedmean = (Subtract @@@Transpose[{data2, movingAvg}]);


enter image description here


Now find the locations of the outliers:


 outpos=Position[subtractedmean, x_ /; x>StandardDeviation[subtractedmean]*3];

Length[outpos]


8




looks like we got the right number of outliers. Removing them.


 newdata=Delete[data2,outpos]
ListPlot[newdata, PlotRange -> All]

enter image description here


To give you an idea of "Threshold" line in this case,


 dathreshold = 
ConstantArray[StandardDeviation[subtractedmean]*3, Length[data2]] + movingAvg;

Here is the "Threshold" line drawn along with the points removed,



 Show[ListPlot[data2, PlotRange -> All, AspectRatio -> 1], 
ListPlot[dathreshold, Joined -> True, PlotStyle -> {Thick, Purple}],
Graphics[{Red,Circle[#, {100, 0.5}] & /@
Thread[{First /@ outpos, data2[[First /@ outpos]]}]}]]

enter image description here


By derivatives


A second way to remove outliers, is by looking at the Derivatives, then threshold on them. Differences in the data are more likely to behave gaussian then the actual distributions.


 diff=Abs@Differences[data2,2];


ListPlot[diff, PlotRange -> All, Joined -> True]

enter image description here


Now you do the same threshold, (based on the standard deviation) on these peaks. Note that the outliers are now really well separated from the actual data. You can find the peak positions that are above the threshold you set, in our case we will keep using $3 \times \sigma$.


You can probably use the peak finding function from V10 (not sure if there is a way to threshold the peaks), but since I stuck in V9 I do the poor's man way.


 newpos=Flatten[Position[Partition[diff, 3, 1], 
x_ /; ((x[[1]] < x[[2]] > x[[3]]) && (x[[2]] > 3*sddif)), {1}]] + 2

newpos===First /@ outpos



True



Comments

Popular posts from this blog

functions - Get leading series expansion term?

Given a function f[x] , I would like to have a function leadingSeries that returns just the leading term in the series around x=0 . For example: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x)] x and leadingSeries[(1/x + 2 + (1 - 1/x^3)/4)/(4 + x)] -(1/(16 x^3)) Is there such a function in Mathematica? Or maybe one can implement it efficiently? EDIT I finally went with the following implementation, based on Carl Woll 's answer: lds[ex_,x_]:=( (ex/.x->(x+O[x]^2))/.SeriesData[U_,Z_,L_List,Mi_,Ma_,De_]:>SeriesData[U,Z,{L[[1]]},Mi,Mi+1,De]//Quiet//Normal) The advantage is, that this one also properly works with functions whose leading term is a constant: lds[Exp[x],x] 1 Answer Update 1 Updated to eliminate SeriesData and to not return additional terms Perhaps you could use: leadingSeries[expr_, x_] := Normal[expr /. x->(x+O[x]^2) /. a_List :> Take[a, 1]] Then for your examples: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x), x] leadingSeries[Exp[x], x] leadingSeries[(1/x + 2 + (1 - 1/x...

How to thread a list

I have data in format data = {{a1, a2}, {b1, b2}, {c1, c2}, {d1, d2}} Tableform: I want to thread it to : tdata = {{{a1, b1}, {a2, b2}}, {{a1, c1}, {a2, c2}}, {{a1, d1}, {a2, d2}}} Tableform: And I would like to do better then pseudofunction[n_] := Transpose[{data2[[1]], data2[[n]]}]; SetAttributes[pseudofunction, Listable]; Range[2, 4] // pseudofunction Here is my benchmark data, where data3 is normal sample of real data. data3 = Drop[ExcelWorkBook[[Column1 ;; Column4]], None, 1]; data2 = {a #, b #, c #, d #} & /@ Range[1, 10^5]; data = RandomReal[{0, 1}, {10^6, 4}]; Here is my benchmark code kptnw[list_] := Transpose[{Table[First@#, {Length@# - 1}], Rest@#}, {3, 1, 2}] &@list kptnw2[list_] := Transpose[{ConstantArray[First@#, Length@# - 1], Rest@#}, {3, 1, 2}] &@list OleksandrR[list_] := Flatten[Outer[List, List@First[list], Rest[list], 1], {{2}, {1, 4}}] paradox2[list_] := Partition[Riffle[list[[1]], #], 2] & /@ Drop[list, 1] RM[list_] := FoldList[Transpose[{First@li...

front end - keyboard shortcut to invoke Insert new matrix

I frequently need to type in some matrices, and the menu command Insert > Table/Matrix > New... allows matrices with lines drawn between columns and rows, which is very helpful. I would like to make a keyboard shortcut for it, but cannot find the relevant frontend token command (4209405) for it. Since the FullForm[] and InputForm[] of matrices with lines drawn between rows and columns is the same as those without lines, it's hard to do this via 3rd party system-wide text expanders (e.g. autohotkey or atext on mac). How does one assign a keyboard shortcut for the menu item Insert > Table/Matrix > New... , preferably using only mathematica? Thanks! Answer In the MenuSetup.tr (for linux located in the $InstallationDirectory/SystemFiles/FrontEnd/TextResources/X/ directory), I changed the line MenuItem["&New...", "CreateGridBoxDialog"] to read MenuItem["&New...", "CreateGridBoxDialog", MenuKey["m", Modifiers-...