2014-11-14 2 views
8

필자는 데이터에 다중 가우스 피팅을 수행하는 방법을 찾고있었습니다. 지금까지 발견 한 예제 중 대부분은 난수 생성을 위해 정규 분포를 사용합니다. 그러나 나는 데이터의 플롯을보고 1-3 피크가 있는지 확인하는 데 관심이 있습니다.파이썬 데이터로드 및 멀티 가우스 적합성

하나의 피크에 대해이 작업을 수행 할 수 있지만 더 많은 작업을 수행하는 방법을 모르겠습니다.

예를 들어,이 데이터를 가지고 : 나는 lmfit를 사용하여 시도 http://www.filedropper.com/data_11

을, 그리고 물론 scipy의,하지만 좋은 결과.

도움 주셔서 감사합니다.

+1

귀하의 질문은 명확하지 않습니다. 귀하의 (다소 시끄러운) 데이터에 가우스를 맞추고 싶습니까? 맥시마의 위치를 ​​찾고 싶습니까? 데이터가 1 ​​~ 3 가우스의 합계입니까? 각각의 평균 및 표준 편차를 구하고 싶습니까? –

+1

안녕하세요. 회신 주셔서 감사합니다 :) 각 피크에 대해 하나의 가우스를 맞추고 싶습니다. – astromath

+0

"데이터가 1 ​​~ 3 가우스의 합계이고 각각의 평균 및 표준 편차를 얻고 싶습니까?" 정확하게! – astromath

답변

11

단일 가우시안 합계의 매개 변수화 된 모델 함수를 간단히 작성하십시오. 초기 추측에 대한 좋은 가치를 선택하십시오 (이것은 매우 중요한 단계입니다). scipy.optimize은 그 숫자를 약간 조정합니다. 여기

는 당신이 그것을 할 수있는 방법은 다음과 같습니다

import numpy as np 
import matplotlib.pyplot as plt 
from scipy import optimize 

data = np.genfromtxt('data.txt') 
def gaussian(x, height, center, width, offset): 
    return height*np.exp(-(x - center)**2/(2*width**2)) + offset 
def three_gaussians(x, h1, c1, w1, h2, c2, w2, h3, c3, w3, offset): 
    return (gaussian(x, h1, c1, w1, offset=0) + 
     gaussian(x, h2, c2, w2, offset=0) + 
     gaussian(x, h3, c3, w3, offset=0) + offset) 

def two_gaussians(x, h1, c1, w1, h2, c2, w2, offset): 
    return three_gaussians(x, h1, c1, w1, h2, c2, w2, 0,0,1, offset) 

errfunc3 = lambda p, x, y: (three_gaussians(x, *p) - y)**2 
errfunc2 = lambda p, x, y: (two_gaussians(x, *p) - y)**2 

guess3 = [0.49, 0.55, 0.01, 0.6, 0.61, 0.01, 1, 0.64, 0.01, 0] # I guess there are 3 peaks, 2 are clear, but between them there seems to be another one, based on the change in slope smoothness there 
guess2 = [0.49, 0.55, 0.01, 1, 0.64, 0.01, 0] # I removed the peak I'm not too sure about 
optim3, success = optimize.leastsq(errfunc3, guess3[:], args=(data[:,0], data[:,1])) 
optim2, success = optimize.leastsq(errfunc2, guess2[:], args=(data[:,0], data[:,1])) 
optim3 

plt.plot(data[:,0], data[:,1], lw=5, c='g', label='measurement') 
plt.plot(data[:,0], three_gaussians(data[:,0], *optim3), 
    lw=3, c='b', label='fit of 3 Gaussians') 
plt.plot(data[:,0], two_gaussians(data[:,0], *optim2), 
    lw=1, c='r', ls='--', label='fit of 2 Gaussians') 
plt.legend(loc='best') 
plt.savefig('result.png') 

result of fitting

당신이 볼 수 있듯이,이 두 맞는 (시각) 사이에 거의 차이가 없다. 3,

err3 = np.sqrt(errfunc3(optim3, data[:,0], data[:,1])).sum() 
err2 = np.sqrt(errfunc2(optim2, data[:,0], data[:,1])).sum() 
print('Residual error when fitting 3 Gaussians: {}\n' 
    'Residual error when fitting 2 Gaussians: {}'.format(err3, err2)) 
# Residual error when fitting 3 Gaussians: 3.52000910965 
# Residual error when fitting 2 Gaussians: 3.82054499044 

을이 경우 : 소스의 현재 또는 단지 2. 그러나 3 가우시안이 있다면 당신이 추측을해야한다면 그래서 당신은 작은 잔류를 확인 후, 확실히 알 수 없다 Gaussians은 더 나은 결과를 제공하지만, 초기 추정도 상당히 정확합니다.

+0

안녕하십니까. 답변 해 주셔서 대단히 감사합니다. 나는 두 개의 분리 된 Gaussians을 얻은 ​​다음 결합하는 방법을 시도했지만 해결책을 본 후에 잘못된 생각 이었음을 이해합니다. "센터"와 " 오프셋 "매개 변수는 무엇입니까? – astromath

+0

이것들은 Gaussian의 평균과 그것에 주어진 수직 오프셋이 될 것입니다. 확인하려면'Gaussian'의 정의를 확인하십시오 ;-) Stackoverflow에 오신 것을 환영합니다. 그것이 도움이된다면 [upvote or accept the answer] (https://stackoverflow.com/help/someone-answers)를 잊지 마십시오. –

+0

Oliver 안녕하세요. 다시 한 번 감사드립니다.) 스크립트를 이해할 수는 있지만 어떻게 의미를 추측 했습니까? 그 점에 대해서는 numpy를 사용 하겠지만 더 간단하고 빠른 것이 있다고 생각합니다. 다시 한번 감사드립니다 :) – astromath