一、 寻找变化率具线性关系的等式
基于数学方法进行分类的模型中,微分方程模型是非常重要的一类(其他的还有几何模型、图论模型、概率模型、最优控制模型、规划论模型、排队论模型和马氏链模型等)。但是怎样建立,却并不那么自明。回忆一下初中解应用题的时候,我们建立方程都是找变量间的一个等式。根据这个思路,微分方程的等式可以这样寻找。
因变量的变化率 = 自变量的线性关系
用形式化语言可以抽象表达为:
方程左边使用的是增量符号Δ,这是便于建模时好理解(一般情况下都是如此,对这个方法熟悉后可以直接使用“率”来理解),毕竟有
X 可以是一个常量,也可以是一个变量,还可以是一个向量,花体字 ℒ 表示一个线性关系。
这就是我们要建立微分方程的一个基本思路和出发点。以下通过两个🌰来说明怎样寻找等式关系。
二、 速度和人口的两个栗子
(一) 匀速直线运动
假设位移为 s ,时间为 t ,匀速直线运动就是相同的时间内位移相同。也就是取极限后得到
这是一个非常简单的可分离变量的微分方程,其通解为 s = kt + C 。C是任意常数,由于 t = 0 时,s = 0 ,这是初始条件,带入上式得到特解 s = kt ,如果把k换成v,就是我们熟悉的 s = vt
你可能会说,这么简单的公式我早就会了,还需要你来说?且慢,一个可以从另外一个角度去看看我们习以为常的一些事情有别的理解。另外,也是为理解以下更复杂模型打一个基础。
对于没有阻力的自由落体运动,该怎样求其s和t之间的关系呢?
提示:根据实验,相同时间内,速度的变化量相同。于是,方程可以建立为:
解出上述公式后,将 s = vt 带入,就能得到s 和 t 的平方之间的线性关系模型。
(二) 人口模型
1. 人口增长率为常数
基本假设 : 人口自然增长率(出生率-死亡率) r 是常数, x(t) ~时刻t的人口,于是可以得到下列等式。
(t 时刻 + Δt时间)后人口变化数量 = t 时刻人口数量人口增长率Δt
即:x(t+Δt)-x(t) = x(t)rΔt, 移项, [x(t+Δt)-x(t)]/ Δt = x(t)r ,
X(0) = x_0 ,方程求解后得到
②就是马尔萨斯在1798年提出的著名的指数增长模型,其实,只要是可以表示为等比数列的关系,都可以通过这个方式建立微分方程。由于资源有限,环境承载能力有上限,故以上模型在短期、人口压力不大的时候可以使用,但是长期一定有一个增长的极限,资源、环境等外部因素对人口增长有阻滞作用人口,反应到增长率是一个随时间变化的量。
2. 人口增长率为变量
当到了某个临界点后,人口数量越多,各种资源越来越难以满足人的需求,出生率一定会下降,死亡率假定不变情况下,人口自然增长率r 是x 的减函数,为了计算简单,我们还是假定 r 和 x 存在线性关系的减函数。
假设: r(x) = r – sx …③ (r,s >0) r~固有增长率(x很小时)
x_M ~人口容量(资源、环境能容纳的最大人口数量),此时,人口不再增加,即 r(x_M) = 0, 于是 s = r/ x_M. , 带入③式得到
将④带入①得
求解⑤得到
图形如下
图 1 阻滞模型
⑥式为著名的logistic 模型,广泛应用于生物、科技、社会等各个领域,这个模型在动力学系统中也叫寻的模型。
三、 一种流行病模型简单理解
按照单一群体和复合群体两种基于群体进行的测算,使用的都是如下几类基本模型,具备S(易感者)、E(潜伏者)、I(感染者)和R(抵抗者)的全部或者部分对象人群,以这四类人群为基础,可以选择全部或者部分构成不同的模型类型,如适合普通流感的SI模型、适合急性传染病的SIR模型等。
(一) SEIR模型的基本构成
单一群体模型采用宏观视角建模,关注整个人群状态的变化,一般采用微分方程描述。最流行的单一群体模型是仓室模型,所有处于相同状态的人构成一个仓室,随着状态的变化,人员在仓室之间移动,仓室模型的基本假设是每个人都在相同人群中均匀混合,接触是瞬时的,并与历史无关。每个仓室的人口数量足够大,在传染病流行过程中,感染率和恢复率是常数。具体的抽象模型见下图:
图 2 传染病一般发病形式对应的专业路线
表 1 符号说明
(一) SEIR模型理解
SEIR模型由⑦、⑧、⑨、⑩四个方程构成的一个方程组和⑪这样一个约束条件组成。方程组中,左边表示不同人群的变化率。∂是求偏导数符号,表示变化量的符号∆趋于0时,就使用∂来表示。如:抵抗者增加的数量/天可以表示为 ΔR/Δt 。当 ΔR→0 与 Δt→0 时,就得到了 ∂R/∂t, 也就是抵抗者增加率。下面分别对四个方程进行说明。
⑦式:易感者S增加的数量/天= -感染概率β(易感者接触感染者后转变为潜伏者的概率)×当天易感者S数量×感染者占比(感染者I/总人数N)。SI/N 表示易感者S接触的人中恰好是感染者I(假设接触的人中只有感染者I具有传染性,其他三类S、E、R都不具备传染性)概率。负号表示易感者每天减少数量。
⑧式:潜伏者E增加的数量/天=当天由易感者转化为潜伏者的数量-转化率σ×当天潜伏者E数量。
⑨式:感染者I增加的数量/天=转化率σ×当天潜伏者数量-治愈率γ×当天感染者数量
⑩式:抵抗者R增加的数量/天=治愈率γ×当天感染者数量I
以上四个方程结合图2传染病一般发病形式对应的专业路线就更清晰了。
s.t N=S+E+I+R……⑪
只有极少部分(偏)微分方程有解析解,SEIR模型也不能得到解析解,通过数值算法求解。具体为用差商代替导数,将微分方程离散化后成为差分方程,利用4阶龙格-库塔(R-K)方法插值求解差分方程,这个过程就是 RK 算法。在R中编程计算还是很方便的,绘图效果也不错。
四、 微分方程模型常见应用场景
(一) 静态优化模型
存储模型,价格模型,消费者均衡模型、万有引力定律的发现
(二) 动态模型
人口模型、传染病模型、经济增长模型、药物在体内的分布与排除、烟雾的扩散与消失
描述对象特征随时间(空间)的演变过程
分析对象特征的变化规律
预报对象特征的未来性质
研究控制对象特征的手段
(三) 稳态模型
种群的竞争模型,军备竞赛模型、正规战与游击战、香烟过滤嘴的作用