微分方程参数估计：如何基于观测数据估算a与b？

嘿，这个问题其实可以通过把非线性微分方程转化为线性回归问题来解决，非常直观，我来一步步给你拆解：

第一步：把微分方程转化为线性形式 #

首先看原方程：
$$\frac{dy}{dt} = (a + bx)y$$
我们可以对两边做个简单变形——同时除以y（这里要注意y必须大于0，要是观测数据里有y≤0的情况，得先检查数据合理性：原方程里y=0是平衡点，后续y都应该保持0，大概率是异常值）：
$$\frac{1}{y}\frac{dy}{dt} = a + bx$$
左边其实就是$\frac{d(\ln y)}{dt}$，所以我们令$z(t) = \ln y(t)$，方程就变成了：
$$\frac{dz}{dt} = a + bx(t)$$
这就变成了一个关于z的导数的线性方程，接下来我们用离散观测数据来近似这个导数。

第二步：离散化处理观测数据 #

假设我们有n组观测数据$(t_1,y_1,x_1), (t_2,y_2,x_2), ..., (t_n,y_n,x_n)$，我们用有限差分近似导数$\frac{dz}{dt}$：
$$\frac{z_i - z_{i-1}}{t_i - t_{i-1}} \approx a + bx_i$$
这里$z_i = \ln y_i$，$t_i - t_{i-1}$是相邻两个时刻的时间差，记为$\Delta t_i$。把式子整理一下：
$$z_i - z_{i-1} = a \cdot \Delta t_i + b \cdot (x_i \cdot \Delta t_i)$$
我们再定义几个新变量简化表达：

• $\Delta z_i = z_i - z_{i-1}$（相邻z值的差）
• $w_i = x_i \cdot \Delta t_i$（x_i和时间差的乘积）

现在式子就变成了标准的线性回归模型：
$$\Delta z_i = a \cdot \Delta t_i + b \cdot w_i + \epsilon_i$$
其中$\epsilon_i$是误差项，代表离散近似带来的误差或者观测噪声。

第三步：用最小二乘法求解参数a和b #

现在问题转化为求解线性回归的系数，我们可以用普通最小二乘法（OLS）——最小化所有误差的平方和：
$$S = \sum_{i=2}^n (\Delta z_i - a \Delta t_i - b w_i)^2$$
对a和b分别求偏导并令其等于0，就能得到正规方程组：
$$
\begin{cases}
a \sum_{i=2}^n (\Delta t_i)^2 + b \sum_{i=2}^n (\Delta t_i w_i) = \sum_{i=2}^n (\Delta t_i \Delta z_i) \
a \sum_{i=2}^n (\Delta t_i w_i) + b \sum_{i=2}^n (w_i)^2 = \sum_{i=2}^n (w_i \Delta z_i)
\end{cases}
$$
解这个二元一次方程组就能得到a和b的估计值。

当然你也不用手动解方程，用工具就能直接算：

• 在Python里，可以用scikit-learn的LinearRegression，把$\Delta t_i$和$w_i$组成特征矩阵，$\Delta z_i$作为目标变量直接拟合；
• 在R里，用lm()函数，把公式写成Delta_z ~ Delta_t + w - 1（因为没有截距项）。

注意事项 #

• 如果时间间隔是固定的（比如$\Delta t_i = \Delta t$为常数），式子会更简单：$\frac{\Delta z_i}{\Delta t} = a + bx_i$，直接对$\frac{\Delta z_i}{\Delta t}$和$x_i$做线性回归，截距就是a，斜率就是b；
• 如果有y≤0的情况，不能直接取对数，这时候可以先检查数据是否有异常，或者考虑引入小常数偏移后再取对数（比如$z_i = \ln(y_i + c)$，不过会引入新参数，复杂度更高）；
• 如果误差项不满足OLS的假设（比如异方差），可以考虑用加权最小二乘法来优化结果。

内容的提问来源于stack exchange，提问作者L.bronze