反应动力学参数拟合与停留时间分布函数——基于Python实现

1 非线性回归Python实现方法

在很多领域中往往需要对实验数据（或者调查数据）通过某种方程进行拟合，得到方程中的关键参数。在统计、物理、化学、生物等领域都是非常常见的方法。
非线性方程是指方程中存在非线性项，例如$e^x$，$\ln x$，$x^3$，$\sin x$等项时，求解过程通过人工计算很难得到准确的解。

1.1 多项式拟合 polyfit()

Python的numpy库中的函数polyfit()能够拟合一元多项式
$$y=a_1{x}^m+a_2{x}^{m-1}+...+a_{m}x+a_{m+1}$$numpy.polyfit()专门用于拟合这种形式的多项式，看下面的例子

import numpy as np
import matplotlib.pyplot as plt

# 生成数据
x = np.linspace(0, 5, 100)
y = np.exp(x)
# 方差初始容器，辅助内容
i_range = [i for i in range(3,10)]
rsds = []
# 拟合并绘制不同项数多项式
for i in i_range:
    p = np.polyfit(x, y, i, full=True)
    y_hat = np.polyval(p[0], x)
    rsds.append(p[1])
    plt.plot(x, y_hat, lw=2, label=str(i))
# 绘制原始数据
plt.plot(x[::10], y[::10], 'r*', ms=20, label='exp(x)')
plt.legend()
# 字体大小设置
font = {
   'size':22}
plt.xlabel('x', font)
plt.ylabel('y', font)
plt.tick_params(labelsize=19)
plt.show()
# 绘制方差图
plt.plot(i_range, np.log10(rsds), '-o')
plt.xlabel('n', font)
plt.ylabel(r'$\lg \sigma^2$', font)
plt.tick_params(labelsize=19)

图1-1多项式拟合

图1-2 多项式拟合方差随多项式次数的变化规律
简单举个例子，$y=e^x$这个函数也可以拟合成多项式的形式，随着多项式次数的增加，拟合结果方差逐渐减小。其原理是连续函数的泰勒公式，高数里证明了连续可导的函数总能表示成多项式的形式。
多项式拟合能在一定程度上获得变化规律，但是也存在着无法忽略的问题，超出原始数据范围的预测效果差，而且随着n的增加，有可能出现过拟合的情况（噪音信号被带入拟合结果）。

---

1.2 非线性拟合 curve_fit()

curve_fit()是scipy.optimize的函数，用于拟合已知形式的函数，接下来还是看个例子。仍然以$y=\exp (x)$为例说明，通过随机数函数增加生成一些噪音，使噪音范围为[-0.2, 0.2），然后分别用curve_fit()和polyfit()拟合并对比。用下面的函数拟合，a和n为待拟合参数。
$$y=a\exp (nx)$$
代码如下：

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

# 生成数据
x = np.linspace(0, 5, 100)
y = np.exp(x)
# 随机函数生成噪音
eps = 0.4*np.array([random() for i in range(len(x))]) - 0.2
y = y * (1 + eps)

# 拟合函数形式
def f(x, a, n):
    return a*np.exp(n*x)

# 非线性拟合
popt,pc