檢測峯寬的魯棒算法

Question

我問how to programmatically judge spectrum bands和@detly使用FWHM（半峯全寬），以確定峯值的寬度建議。我四處搜尋，發現FWHM可以用於擬合模型（我幾乎是一個普通人），特別是高斯模型。具體而言，對於高斯模型，2.354 * sigma是the width。

我在找一個高斯適合因爲存在高峯。這張圖中有4個格局良好的山峯。然後是「雙」峯（儘管兩者可能都很重要）和兩個展開峯。他們將證明對一個天真的FWHM不可能的挑戰。

由於它的x座標是大概位置，您可以幫助生成高斯在Scip/Numpy中的高斯擬合（用於FWHM計算）嗎？如果Guassian是一個不好的選擇，那麼一些其他的計劃。

Answer 1

擬合高斯是一個很好的方法。如果你對初始參數值有好的猜測，你可以嘗試一下猜測它們。一個很大的問題是噪聲，實際上你可能想要將每個峯值單獨進行擬合（即只考慮一個給定峯值的範圍），或者在數據上得到基線噪聲曲線並將其減去。

這裏是一些代碼，以適應多個高斯。我已經爲我認爲是最突出的8個峯值輸入了一些非常鬆散的參數，再加上一個非常寬泛的高斯帶來嘗試捕捉背景噪音。在處理之前，我會稍微清理已發佈的圖像，以幫助獲取數據（刪除鼠標光標和軸邊緣，並倒轉圖像）。

代碼：

import Image 
from scipy import * 
from scipy.optimize import leastsq 
import numpy as np 

im = Image.open("proc.png") 
im = im.convert('L') 
h, w = im.size 
#Extract data from the processed image: 
im = np.asarray(im) 
y_vals, junk = np.mgrid[w:0:-1, h:0:-1] 
y_vals[im < 255] = 0 
y_vals = np.amax(y_vals,axis=0) 

def gaussian(x, A, x0, sig): 
    return A*exp(-(x-x0)**2/(2.0*sig**2)) 

def fit(p,x): 
    return np.sum([gaussian(x, p[i*3],p[i*3+1],p[i*3+2]) 
        for i in xrange(len(p)/3)],axis=0) 

err = lambda p, x, y: fit(p,x)-y 

#params are our intitial guesses for fitting gaussians, 
#(Amplitude, x value, sigma): 
params = [[50,40,5],[50,110,5],[100,160,5],[100,220,5], 
      [50,250,5],[100,260,5],[100,320,5], [100,400,5], 
      [30,300,150]] # this last one is our noise estimate 
params = np.asarray(params).flatten() 

x = xrange(len(y_vals)) 
results, value = leastsq(err, params, args=(x, y_vals)) 

for res in results.reshape(-1,3): 
    print "amplitude, position, sigma", res 

import pylab 
pylab.subplot(211, title='original data') 
pylab.plot(y_vals) 
pylab.subplot(212, title='guassians fit') 
y = fit(results, x) 
pylab.plot(x, y) 
pylab.savefig('fig2.png') 
pylab.show() 

這些擬合輸出高斯參數：

#Output: 
amplitude, position, sigma [ 23.47418352 41.27086097 5.91012897] 
amplitude, position, sigma [ 16.26370263 106.39664543 3.67827824] 
amplitude, position, sigma [ 59.74500239 163.11210316 2.89866545] 
amplitude, position, sigma [ 57.91752456 220.24444939 2.91145375] 
amplitude, position, sigma [ 39.78579841 251.25955921 2.74646334] 
amplitude, position, sigma [ 86.50061756 260.62004831 3.08367483] 
amplitude, position, sigma [ 79.26849867 319.64343319 3.57224402] 
amplitude, position, sigma [ 127.97009966 399.27833126 3.14331212] 
amplitude, position, sigma [ 20.21224956 379.41066063 195.47122312]