1
$\begingroup$

I have a huge data file which I can't ListPlot.

This code generates similar kind of data:

datatest =RandomSample[Join[RandomReal[{0.5, 15}, 20], RandomReal[.1, 10000]]];
datatest2 = 5 + Riffle[datatest, -datatest];

I want to filter (delete) the part of the data that is not necessary as follows:

peaks = FindPeaks[datatest2, 0, 0, 5.2];
ListPlot[datatest2, PlotRange -> All, Joined -> True, 
 Epilog -> {Red, PointSize[0.01], Point[peaks]}]

enter image description here Currently I am using some kind of long way to do the task.

Is there any signal processing functionality in MMA that can do this easily?

Thank you

$\endgroup$
  • 2
    $\begingroup$ I don't quite understand. Isn't the output of FindPeaks exactly what you're looking for? $\endgroup$ – yohbs Feb 1 '15 at 22:32
  • $\begingroup$ I want to keep the peaks and some data around. In another word I want to delete the marked parts in the plot. $\endgroup$ – Algohi Feb 1 '15 at 22:34
  • $\begingroup$ I have more than 1 Billion data point which I can't visualize. most of the data is similar to the parts marked in the plot. I want to delete these parts so that I can visualize the data $\endgroup$ – Algohi Feb 1 '15 at 22:36
  • $\begingroup$ I think FindPeaks[] is the best way to do this. If it is taking too long, split datatest2 into $k$ equal portions (e.g., $k=8$ or $16$) and then Parallelize[] your FindPeaks[] operation. $\endgroup$ – David G. Stork Feb 1 '15 at 22:42
  • 1
    $\begingroup$ The question is a bit unclear (squiggly marks on a graph is not a substitute for a precise specification) but this previous question seems related: mathematica.stackexchange.com/q/22528/484 $\endgroup$ – Rahul Feb 2 '15 at 1:50
3
$\begingroup$

I am not an expert (understatemnt of the year!) in signal processing but you can use the band function and create a sparse matrix that has as many 1s as you want around the positions of the peaks. I am not sure the following is the best way to do this but it works:

With[{width = 200},
  spArray = SparseArray[
    Thread[
     (Band[# - width, # + width] & /@ peaks[[All, 1]]) -> 1], 
    Length[datatest2]]
  ];

so this means that I have made a sparse array with 1s in 400 samples ($\pm 200$) around each of your peaks and if I multiply this with your original array dataset2 these are gonna be the only elements that will survive:

ListPlot[spArray datatest2, PlotRange -> All, Joined -> True, 
 Epilog -> {Red, PointSize[0.01], Point[peaks]}]

enter image description here

$\endgroup$
  • $\begingroup$ This is a good approach however the number of points still the same which is the main problem. you can image how impossible to plot data contains 1 Billion point. I will think of way to remove those zeros without affecting the zeros around the peaks. thanks $\endgroup$ – Algohi Feb 2 '15 at 21:31
  • $\begingroup$ This will do the job: pos = Flatten[Position[Normal[spArray], 1]]; datatest2[[pos]]; $\endgroup$ – Algohi Feb 2 '15 at 21:46
  • $\begingroup$ @Algohi I don't quite understand what the connection is between the filter and the plotting bit. If you need to plot, why not use the method I posted and subsample the resulting array. Is it billions of data with only a few peaks or with millions of peaks? I am sure we (SE) think of something if you re edit to clarify the question a bit $\endgroup$ – gpap Feb 2 '15 at 22:18
  • $\begingroup$ It is billion with hundreds of peaks. when I said filter I did not know what the correct word is. that is why I add delete after. I hope this clarify the issue and sorry for any ambiguity:) $\endgroup$ – Algohi Feb 2 '15 at 22:23
  • $\begingroup$ @Algohi, I didn't mean it in a frustrated way :) - and by no means don't feel pressured to accept if it's not working - but there's a few things I don't understand still. Does FindPeaks work at all with your data? I would doubt that it does with a billion points. So do you know where the peaks are? Do you know how many there are if not? If you don't that's a (difficult) question in itself: "how to find hundreds of peaks among billions of points?". $\endgroup$ – gpap Feb 3 '15 at 10:21
1
$\begingroup$

I'm not sure I understand what you need, but here's my try: If you want to keep only 2d+1 data points around each peak, you can use

toKeep = Map[# + Range[-d, d] &, peaks[[All, 1]]];
choppedData = Map[Part[datatest2, #] &, toKeep];

choppedData is a list of lists, the i-th list contains the 2d+1 values around the i-th peak.

To speed up you can use ParallelMap instead of Map.

$\endgroup$
  • $\begingroup$ Thanks but note here that, peaks[[;;,1]] does not necessary to be Integer number. you may add Floor or Ceiling. Second, there will be intersection if d is larger than the minimum distance between peaks. This will create additional work to do. Any way thanks for the answer. $\endgroup$ – Algohi Feb 1 '15 at 22:56
1
$\begingroup$

A common way to remove outliers is with the Median filter. What you want to do is the opposite: to keep the outliers and remove the inliers. Subtracting the data from the median, then clipping the result and selecting all those larger than a threshold is one way to proceed.

short = Select[Chop[datatest2 - MedianFilter[datatest2, 5], 0.5], Abs[#] > 0.1 &]; 
ListPlot[short, PlotRange -> All]

enter image description here

One downside is that this removes the mean value (about five in your data). You may wish to add this back in, and also to fiddle with the parameters to get the width you are looking for.

$\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.