1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
|
PeakUtils tutorial
==================
This tutorial shows the basic usage of PeakUtils to detect the peaks of
1D data.
Importing the libraries
-----------------------
.. code:: python
import numpy
import peakutils
from peakutils.plot import plot as pplot
from matplotlib import pyplot
%matplotlib inline
Preparing the data
------------------
Lets generate some noisy data from two Gaussians:
.. code:: python
centers = (30.5, 72.3)
x = numpy.linspace(0, 120, 121)
y = (peakutils.gaussian(x, 5, centers[0], 3) +
peakutils.gaussian(x, 7, centers[1], 10) +
numpy.random.rand(x.size))
pyplot.figure(figsize=(10,6))
pyplot.plot(x, y)
pyplot.title("Data with noise")
.. image:: _static/tut_a_1.png
Getting a first estimate of the peaks
-------------------------------------
By using peakutils.indexes, we can get the indexes of the peaks from the
data. Due to the noise, it will be just a rough approximation.
.. code:: python
indexes = peakutils.indexes(y, thres=0.5, min_dist=30)
print(indexes)
print(x[indexes], y[indexes])
pyplot.figure(figsize=(10,6))
pplot(x, y, indexes)
pyplot.title('First estimate')
.. parsed-literal::
[31 74]
[ 31. 74.] [ 5.67608909 7.79403394]
.. image:: _static/tut_a_2.png
Enhancing the resolution by interpolation
-----------------------------------------
We can enhance the resolution by using interpolation. We will try to fit
a Gaussian near each previously detected peak.
.. code:: python
peaks_x = peakutils.interpolate(x, y, ind=indexes)
print(peaks_x)
.. parsed-literal::
[ 30.58270223 72.34348214]
Estimating and removing the baseline
------------------------------------
It is common for data to have an undesired baseline.
*PeakUtils* implements a function for estimating the baseline by using
an iterative polynomial regression algorithm.
.. code:: python
y2 = y + numpy.polyval([0.002,-0.08,5], x)
pyplot.figure(figsize=(10,6))
pyplot.plot(x, y2)
pyplot.title("Data with baseline")
.. image:: _static/tut_a_3.png
.. code:: python
base = peakutils.baseline(y2, 2)
pyplot.figure(figsize=(10,6))
pyplot.plot(x, y2-base)
pyplot.title("Data with baseline removed")
.. image:: _static/tut_a_4.png
Related functionality in SciPy
------------------------------
SciPy also implements functions that can be used for peak detection.
Some examples:
- `scipy.signal.find\_peaks\_cwt <http://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.find_peaks_cwt.html>`__
- `scipy.signal.savgol\_filter <http://docs.scipy.org/doc/scipy/reference/generated/scipy.signal.savgol_filter.html>`__
|
