Пики подгонки с Scipy curve_fit, оптимальные параметры ошибки не найдены - PullRequest
0 голосов
/ 18 октября 2018

Я недавно начал с Python, потому что у меня есть огромное количество данных, где я хочу автоматически подгонять гауссиан к пикам в спектрах.Ниже приведен пример трех пиков, которые я хочу совместить с тремя отдельными пиками.

Я нашел вопрос, где кто-то ищет что-то очень похожее, Как я могу подогнать несколько гауссовых кривых к масс-спектрометрии?данные в Python? , и приняли его к моему сценарию.

Я добавил свой код внизу, и когда я запускаю последний раздел, я получаю сообщение об ошибке «RuntimeError: Оптимальные параметры не найдены: число вызовов функции достигло maxfev = 800».Чего мне не хватает?

Данные можно скачать по адресу https://www.dropbox.com/s/zowawljcjco70yh/data_so.h5?dl=0

enter image description here

#%%
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import sparse
from scipy.sparse.linalg import spsolve
from scipy.optimize import curve_fit
#%% Read data
path = 'D:/Python/data_so.h5'
f = pd.read_hdf(path, mode = 'r')
t = f.loc[:, 'Time stamp']
d = f.drop(['Time stamp', 'Name spectrum'], axis = 1)
#%% Extract desired wavenumber range
wn_st=2000
wn_ed=2500
ix_st=np.argmin(abs(d.columns.values-wn_st))
ix_ed=np.argmin(abs(d.columns.values-wn_ed))
d = d.iloc[:, ix_st:ix_ed+1]
#%% AsLS baseline correction
spectrum = 230
y = d.iloc[spectrum]

niter = 10
lam = 200000
p = 0.005

L = len(y)
D = sparse.diags([1,-2,1],[0,-1,-2], shape=(L,L-2))
w = np.ones(L)
for i in range(niter):
    W = sparse.spdiags(w, 0, L, L)
    Z = W + lam * D.dot(D.transpose())
    z = spsolve(Z, w*y)
    w = p * (y > z) + (1-p) * (y < z)

corr = d.iloc[spectrum,:] - z
#%% Plot spectrum, baseline and corrected spectrum
plt.clf()
plt.plot(d.columns, d.iloc[spectrum,:])
plt.plot(d.columns, z)
plt.plot(d.columns, corr)
plt.gca().invert_xaxis()
plt.show()
#%%
x = d.columns.values
def gauss(x, a, mu, sig):
    return a*np.exp(-(x.astype(float)-mu)**2/(2*sig**2))

fitx = x[(x>2232)*(x<2252)]
fity = y[(x>2232)*(x<2252)]
mu=np.sum(fitx*fity)/np.sum(fity)
sig=np.sqrt(np.sum(fity*(fitx-mu)**2)/np.sum(fity))

popt, pcov = curve_fit(gauss, fitx, fity, p0=[max(fity),mu, sig])
plt.plot(x, gauss(x, popt[0],popt[1],popt[2]), 'r-', label='fit')
...