Нелинейная регрессия e^(-x) с использованием scipy, python, numpy

Приведенный ниже код дает мне плоскую линию для линии наилучшего соответствия, а не красивую кривую вдоль модели e ^ (-x), которая соответствовала бы данным. Может ли кто-нибудь показать мне, как исправить приведенный ниже код, чтобы он соответствовал моим данным?

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

def _eNegX_(p,x):
    x0,y0,c,k=p  
    y = (c * np.exp(-k*(x-x0))) + y0
    return y

def _eNegX_residuals(p,x,y):
    return y - _eNegX_(p,x)

def Get_eNegX_Coefficients(x,y):
    print 'x is:  ',x  
    print 'y is:  ',y 

    # Calculate p_guess for the vectors x,y.  Note that p_guess is the
    # starting estimate for the minimization.
    p_guess=(np.median(x),np.min(y),np.max(y),.01)

    # Calls the leastsq() function, which calls the residuals function with an initial 
    # guess for the parameters and with the x and y vectors.  Note that the residuals
    # function also calls the _eNegX_ function.  This will return the parameters p that
    # minimize the least squares error of the _eNegX_ function with respect to the original
    # x and y coordinate vectors that are sent to it.
    p, cov, infodict, mesg, ier = scipy.optimize.leastsq(  
        _eNegX_residuals,p_guess,args=(x,y),full_output=1,warning=True)

    # Define the optimal values for each element of p that were returned by the leastsq() function. 
    x0,y0,c,k=p  
    print('''Reference data:\  
    x0 = {x0}
    y0 = {y0}
    c = {c}
    k = {k}
    '''.format(x0=x0,y0=y0,c=c,k=k))  

    print 'x.min() is:  ',x.min()
    print 'x.max() is:  ',x.max()
    # Create a numpy array of x-values
    numPoints = np.floor((x.max()-x.min())*100)
    xp = np.linspace(x.min(), x.max(), numPoints)
    print 'numPoints is:  ',numPoints
    print 'xp is:  ',xp
    print 'p is:  ',p
    pxp=_eNegX_(p,xp)
    print 'pxp is:  ',pxp

    # Plot the results  
    plt.plot(x, y, '>', xp, pxp, 'g-')
    plt.xlabel('BPM%Rest') 
    plt.ylabel('LVET/BPM',rotation='vertical')
    plt.xlim(0,3)
    plt.ylim(0,4)
    plt.grid(True) 
    plt.show()

    return p

# Declare raw data for use in creating regression equation 
x = np.array([1,1.425,1.736,2.178,2.518],dtype='float')  
y = np.array([3.489,2.256,1.640,1.043,0.853],dtype='float')  

p=Get_eNegX_Coefficients(x,y)

person MedicalMath    schedule 20.12.2010    source источник


Ответы (1)


Похоже, это проблема с вашими первоначальными предположениями; что-то вроде (1, 1, 1, 1) работает нормально:график, который выглядит хорошо
У вас есть

p_guess=(np.median(x),np.min(y),np.max(y),.01)

для функции

def _eNegX_(p,x):
    x0,y0,c,k=p  
    y = (c * np.exp(-k*(x-x0))) + y0
    return y

Итак, test_data_maxe^(-.01(x - test_data_median)) + test_data_min

Я не очень разбираюсь в искусстве выбора хороших стартовых параметров, но могу кое-что сказать. leastsq здесь находит локальный минимум - ключ к выбору этих значений - найти правильную гору, на которую нужно взобраться, а не пытаться сократить работу, которую должен выполнять алгоритм минимизации. Ваше первоначальное предположение выглядит так (green): (1.736, 0.85299999999999998, 3.4889999999999999, 0.01) alt text

что приводит к вашей плоской линии (синяя): (-59.20295956, 1.8562 , 1.03477144, 0.69483784)

При регулировке высоты линии был получен больший выигрыш, чем при увеличении значения k. Если вы знаете, что подходите для таких данных, используйте большее значение k. Если вы не знаете, я думаю, вы могли бы попытаться найти приличное значение k путем выборки ваших данных или обратного отсчета от наклона между средним значением первой и второй половины, но я не знаю, как это сделать. об этом.

Редактировать: вы также можете начать с нескольких предположений, запустить минимизацию несколько раз и выбрать линию с наименьшими остатками.

person Thomas    schedule 21.12.2010
comment
Другой простой способ угадать k — просто использовать y.ptp() / x.ptp() (т. е. наклон точек). В любом случае, он должен попасть на стадион. Определенно лучше исходного фиксированного начального предположения 0,01, если значение для k обычно не близко к 0,01! В любом случае, +1 от меня... Это определенно проблема со стартовым предположением. - person Joe Kington; 21.12.2010
comment
Это действительно помогло мне, спасибо. Я считаю этот вопрос решенным. - person MedicalMath; 21.12.2010
comment
@MedicalMath, тогда почему бы не щелкнуть зеленую галочку и не наградить его за ответ. - person ptomato; 21.12.2010