- 非线性函数不转化成线性函数
- 直接使用 s c i p y . o p t i m i z e scipy.optimize scipy.optimize 下面的 c u r v e curve curve_ f i t ( ) fit() fit()拟合
scipy.optimize.curve_fit(f, xdata, ydata, p0=None, sigma=None)
输入参数：
f
f
f: 拟合函数
x
d
a
t
a
xdata
xdata,
y
d
a
t
a
ydata
ydata: 被拟合数据点
p
0
p0
p0: 拟合参数初始值
s
i
g
m
a
sigma
sigma: 数据点
y
d
a
t
a
ydata
ydata的误差
返回值 :
p
o
p
t
popt
popt: 拟合参数
p
c
o
v
pcov
pcov: 拟合参数的协方差
矩阵:
参数标准偏差 :
perr = np.diag(pcov)
实例：
发射源的发射强度具有指数形式
I
=
I
0
e
−
α
t
I = I_0e^{-alpha t}
I=I0e−αt.现有观测数据如下，试确定
I
0
I_0
I0和
α
alpha
α.
t
i
0.2
0.3
0.4
0.5
0.6
0.7
0.8
I
i
3.16
2.38
1.75
1.34
1
0.74
0.56
begin{array}{|c|l|l|l|l|l|l|l|} hline boldsymbol{t}_{boldsymbol{i}} & 0.2 & 0.3 & 0.4 & 0.5 & 0.6 & 0.7 & 0.8 \ hline I_{i} & 3.16 & 2.38 & 1.75 & 1.34 & 1 & 0.74 & 0.56 \ hline end{array}
tiIi0.23.160.32.380.41.750.51.340.610.70.740.80.56
Python代码实现：
import numpy as np
from scipy.optimize import curve_fit
x = np.arange(0.2, 0.9, 0.1)
y = np.array([3.16,2.38,1.75,1.34,1,0.74,
0.56])
def func(x, a, b): #定义拟合模型
return a * np.exp(b * x)
pars, pcov = curve_fit(func, x, y)
print(pars) #拟合参数
print(np.sqrt(np.diag(pcov))) #拟合误差
运行结果：
[ 5.63606121 -2.89059333] [0.03875286 0.01987778]置信区间
拟合参数的置信区间：
β
0
^
=
β
0
±
t
α
2
(
n
−
p
)
σ
β
0
β
1
^
=
β
1
±
t
α
2
(
n
−
p
)
σ
β
1
begin{array}{l} widehat{beta_{0}}=beta_{0} pm t_{frac{alpha}{2}}(n-p) sigma_{beta_{0}} \ widehat{beta_{1}}=beta_{1} pm t_{frac{alpha}{2}}(n-p) sigma_{beta_{1}} end{array}
β0
=β0±t2α(n−p)σβ0β1
=β1±t2α(n−p)σβ1
푛: 数据点的个数
푝: 拟合参数的个数
σ
β
0
sigma_{beta_0}
σβ0,
σ
β
1
sigma_{beta_1}
σβ1, 拟合参数的标准偏差
拟合参数的置信区间（Python代码）
import numpy as np
from scipy.stats.distributions import t
from scipy.optimize import curve_fit
x = np.arange(0.2, 0.9, 0.1)
y = np.array([3.16,2.38,1.75,1.34,1,0.74,
0.56])
def func(x, a, b): #定义拟合模型
return a * np.exp(b * x)
pars, pcov = curve_fit(func, x, y)
print(pars) #拟合参数
print(np.sqrt(np.diag(pcov))) #拟合误差
n = len(y) # number of data points
p = len(pars) # number of parameters
alpha = 0.05 # 95% confidence interval = 100*(1-alpha)
dof = n - p #自由度
# student-t value for the dof and confidence level
tval = t.ppf(1.0-alpha/2., dof) #t_0.975(n-p)
sigma = np.sqrt(np.diag(pcov)) #标准偏差
for i in range(p):
print('p{0}: {1:5.3f} [{2:5.3f} {3:5.3f}]'
.format(i,pars[i], #参数
pars[i]-sigma[i]*tval, #置信下限
pars[i]+sigma[i]*tval)) #置信上限
运行结果：
[ 5.63606121 -2.89059333] [0.03875286 0.01987778] p0: 5.636 [5.536 5.736] p1: -2.891 [-2.942 -2.839]
