如何获得多个峰下的面积值
How to get value of area under multiple peaks
我有一些来自生物分析仪的数据,这些数据给出了时间(x 轴)和吸光度值(y 轴)。时间是每 0.05 秒一次,从 32 秒到 138 秒,所以你可以想象我有多少数据点。我已经使用 plotly 和 matplotlib 创建了一个图形,这样我就有更多的库可以使用来找到解决方案,所以任何一个库中的解决方案都可以!我想要做的是让我的脚本找到每个峰值下方的区域和 return 我的值。
def create_plot(sheet_name):
sample = book.sheet_by_name(sheet_name)
data = [[sample.cell_value(r, c) for r in range(sample.nrows)] for c in range(sample.ncols)]
y = data[2][18:len(data[2]) - 2]
x = np.arange(32, 138.05, 0.05)
indices = peakutils.indexes(y, thres=0.35, min_dist=0.1)
peaks = [y[i] for i in indices]
此代码段获取我的 Y 值、X 值和峰索引。现在有没有办法得到每条曲线下的面积?假设有 15 个索引。
图表如下所示:
自动回答
给定一组 x
和 y
值以及一组 peaks
(峰的 x
坐标),以下是您可以如何自动找到每个峰下的面积。我假设 x
、y
和 peaks
都是 Numpy 数组:
import numpy as np
# find the minima between each peak
ixpeak = x.searchsorted(peaks)
ixmin = np.array([np.argmin(i) for i in np.split(y, ixpeak)])
ixmin[1:] += ixpeak
mins = x[ixmin]
# split up the x and y values based on those minima
xsplit = np.split(x, ixmin[1:-1])
ysplit = np.split(y, ixmin[1:-1])
# find the areas under each peak
areas = [np.trapz(ys, xs) for xs,ys in zip(xsplit, ysplit)]
输出:
示例数据已设置,因此每个峰下的面积(或多或少)保证为 1.0
,因此底部图中的结果是正确的。绿色 X 标记是每两个峰之间的最小值的位置。 "belonging" 到每个峰的曲线部分被确定为与每个峰相邻的最小值之间的曲线部分。
完整代码
这是我用来生成示例数据的完整代码:
import scipy as sp
import scipy.stats
prec = 1e5
n = 10
N = 150
r = np.arange(0, N+1, N//n)
# generate some reasonable fake data
peaks = np.array([np.random.uniform(s, e) for s,e in zip(r[:-1], r[1:])])
x = np.linspace(0, N + n, num=int(prec))
y = np.max([sp.stats.norm.pdf(x, loc=p, scale=.4) for p in peaks], axis=0)
以及我用来制作图表的代码:
import matplotlib.pyplot as plt
# plotting stuff
plt.figure(figsize=(5,7))
plt.subplots_adjust(hspace=.33)
plt.subplot(211)
plt.plot(x, y, label='trace 0')
plt.plot(peaks, y[ixpeak], '+', c='red', ms=10, label='peaks')
plt.plot(mins, y[ixmin], 'x', c='green', ms=10, label='mins')
plt.xlabel('dep')
plt.ylabel('indep')
plt.title('Example data')
plt.ylim(-.1, 1.6)
plt.legend()
plt.subplot(212)
plt.bar(np.arange(len(areas)), areas)
plt.xlabel('Peak number')
plt.ylabel('Area under peak')
plt.title('Area under the peaks of trace 0')
plt.show()
我有一些来自生物分析仪的数据,这些数据给出了时间(x 轴)和吸光度值(y 轴)。时间是每 0.05 秒一次,从 32 秒到 138 秒,所以你可以想象我有多少数据点。我已经使用 plotly 和 matplotlib 创建了一个图形,这样我就有更多的库可以使用来找到解决方案,所以任何一个库中的解决方案都可以!我想要做的是让我的脚本找到每个峰值下方的区域和 return 我的值。
def create_plot(sheet_name):
sample = book.sheet_by_name(sheet_name)
data = [[sample.cell_value(r, c) for r in range(sample.nrows)] for c in range(sample.ncols)]
y = data[2][18:len(data[2]) - 2]
x = np.arange(32, 138.05, 0.05)
indices = peakutils.indexes(y, thres=0.35, min_dist=0.1)
peaks = [y[i] for i in indices]
此代码段获取我的 Y 值、X 值和峰索引。现在有没有办法得到每条曲线下的面积?假设有 15 个索引。
图表如下所示:
自动回答
给定一组 x
和 y
值以及一组 peaks
(峰的 x
坐标),以下是您可以如何自动找到每个峰下的面积。我假设 x
、y
和 peaks
都是 Numpy 数组:
import numpy as np
# find the minima between each peak
ixpeak = x.searchsorted(peaks)
ixmin = np.array([np.argmin(i) for i in np.split(y, ixpeak)])
ixmin[1:] += ixpeak
mins = x[ixmin]
# split up the x and y values based on those minima
xsplit = np.split(x, ixmin[1:-1])
ysplit = np.split(y, ixmin[1:-1])
# find the areas under each peak
areas = [np.trapz(ys, xs) for xs,ys in zip(xsplit, ysplit)]
输出:
示例数据已设置,因此每个峰下的面积(或多或少)保证为 1.0
,因此底部图中的结果是正确的。绿色 X 标记是每两个峰之间的最小值的位置。 "belonging" 到每个峰的曲线部分被确定为与每个峰相邻的最小值之间的曲线部分。
完整代码
这是我用来生成示例数据的完整代码:
import scipy as sp
import scipy.stats
prec = 1e5
n = 10
N = 150
r = np.arange(0, N+1, N//n)
# generate some reasonable fake data
peaks = np.array([np.random.uniform(s, e) for s,e in zip(r[:-1], r[1:])])
x = np.linspace(0, N + n, num=int(prec))
y = np.max([sp.stats.norm.pdf(x, loc=p, scale=.4) for p in peaks], axis=0)
以及我用来制作图表的代码:
import matplotlib.pyplot as plt
# plotting stuff
plt.figure(figsize=(5,7))
plt.subplots_adjust(hspace=.33)
plt.subplot(211)
plt.plot(x, y, label='trace 0')
plt.plot(peaks, y[ixpeak], '+', c='red', ms=10, label='peaks')
plt.plot(mins, y[ixmin], 'x', c='green', ms=10, label='mins')
plt.xlabel('dep')
plt.ylabel('indep')
plt.title('Example data')
plt.ylim(-.1, 1.6)
plt.legend()
plt.subplot(212)
plt.bar(np.arange(len(areas)), areas)
plt.xlabel('Peak number')
plt.ylabel('Area under peak')
plt.title('Area under the peaks of trace 0')
plt.show()