指数曲线拟合的置信区间 [英] Confidence interval for exponential curve fit

问题描述

我正在尝试对某些x,y数据进行指数拟合以获取置信区间(在此处 ).这是我必须找到最适合数据的指数的MWE:

I'm trying to obtain a confidence interval on an exponential fit to some x,y data (available here). Here's the MWE I have to find the best exponential fit to the data:

from pylab import *
from scipy.optimize import curve_fit

# Read data.
x, y = np.loadtxt('exponential_data.dat', unpack=True)

def func(x, a, b, c):
    '''Exponential 3-param function.'''
    return a * np.exp(b * x) + c

# Find best fit.
popt, pcov = curve_fit(func, x, y)
print popt

# Plot data and best fit curve.
scatter(x, y)
x = linspace(11, 23, 100)
plot(x, func(x, *popt), c='r')
show()

产生:

我该如何使用纯python，numpy或scipy(这是我已经安装的软件包)获得此拟合值的95％(或其他值)的置信区间?

How can I obtain the 95% (or some other value) confidence interval on this fit preferably using either pure python, numpy or scipy (which are the packages I already have installed)?

推荐答案

Gabriel的答案不正确.此处以红色表示由GraphPad Prism计算得出的他的数据的95％置信带:

Gabriel's answer is incorrect. Here in red the 95% confidence band for his data as calculated by GraphPad Prism:

背景:拟合曲线的置信区间"通常称为置信带.对于95％的置信带，可以95％的置信度包含真实曲线. (这与上面以灰色显示的预测带不同.预测带与将来的数据点有关.有关更多详细信息，请参见例如此

Background: the "confidence interval of a fitted curve" is typically called confidence band. For a 95% confidence band, one can be 95% confident that it contains the true curve. (This is different from prediction bands, shown above in gray. Prediction bands are about future data points. For more details, see, e.g., this page of the GraphPad Curve Fitting Guide.)

在Python中， kmpfit 可以计算非-线性最小二乘法.以下是加百列的例子:

In Python, kmpfit can calculate the confidence band for non-linear least squares. Here for Gabriel's example:

from pylab import *
from kapteyn import kmpfit

x, y = np.loadtxt('_exp_fit.txt', unpack=True)

def model(p, x):
  a, b, c = p
  return a*np.exp(b*x)+c

f = kmpfit.simplefit(model, [.1, .1, .1], x, y)
print f.params

# confidence band
a, b, c = f.params
dfdp = [np.exp(b*x), a*x*np.exp(b*x), 1]
yhat, upper, lower = f.confidence_band(x, dfdp, 0.95, model)

scatter(x, y, marker='.', s=10, color='#0000ba')
ix = np.argsort(x)
for i, l in enumerate((upper, lower, yhat)):
  plot(x[ix], l[ix], c='g' if i == 2 else 'r', lw=2)
show()

dfdp是关于每个参数p(即a，b和c)的模型f = a * e ^(b * x)+ c的偏导数∂f/∂p.有关背景，请参见 kmpfit教程或此页面 GraphPad曲线拟合指南. (与我的示例代码不同，kmpfit教程不使用 confidence_band() 来自库，但其实现略有不同.)

The dfdp are the partial derivatives ∂f/∂p of the model f = a*e^(b*x) + c with respect to each parameter p (i.e., a, b, and c). For background, see the kmpfit Tutorial or this page of the GraphPad Curve Fitting Guide. (Unlike my sample code, the kmpfit Tutorial does not use confidence_band() from the library but its own, slightly different, implementation.)

最后，Python图与棱镜1相匹配:

Finally, the Python plot matches the Prism one:

这篇关于指数曲线拟合的置信区间的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！