数学模型课程期末复习提纲（下）

MathModeling(Ⅱ)

主编： 😄 审阅：🔥

第五章：预测算法

1. 回归预测

   拟合最常用的方法为最小二乘法，其原理是：因变量的实际值与拟合值之差称为残差，将所有残差平方后相加，即得残差平方和，”最好“拟合效果就是使残差平方和最小，于是运用极值原理，将”泛函问题“转换为”求残差平方和最小”问题（具体思想可参考数值分析相关内容）

   • 多元线性回归
     $$y={\beta }_0+{\beta }_1{x}_1+\cdots +{\beta }_p{x}_p+\epsilon$$
     其中 $x_i,i=1,\cdots ,p$ 称为回归变量（自变量），$y$ 称为回归因变量（ $n$ 维列向量），${\beta }_i,i=1,\cdots ,p$ 称为回归系数，$\epsilon$ 应大致服从均值为 $0$ 的正态分布。特别地，当 $p=1$ 时，退化为一元线性回归。

     • $\text{Matlab}$ 实现

       b=regress(y,X)
       [b,bint]=regress(y,X)
       [b,bint,r]=regress(y,X)
       [b,bint,r,rint]=regress(y,X)
       [b,bint,r,rint,stats]=regress(y,X)
       []=regress(y,X,alpha)
       

       输入：

       y	$n\times 1$ 数值向量
       X	$n\times p$ 数值矩阵
       $\alpha$	置信水平（缺省值 $0.05$ ）

       输出：

       b	系数估计，$p\times 1$ 数值向量
       bint	系数估计的置信区间
       r	残差，$p\times 1$ 数值向量
       rint	残差的置信区间
       stats	包括 $R^2$ （决定系数，越大越好），$F$ （方差齐性检验的统计量，$F>F_{\alpha }(p+1,n)$ ），$p$ （ $F$ 检验的 $p$ 值，$p<\alpha$ ），$s^2$ （误差方差的估计，越小越好）

       $F_{\alpha }(p+1,n)$ ：利用 $R$ 实现

       qf(1-alpha,p+1,n) #求F分布的1-alpha分位数
       

   • 多项式拟合
     $$y=a_0+a_{1}x+a_2{x}^2+\cdots +a_m{x}^m$$

     • $\text{Matlab}$ 实现——拟合曲线相关参数

       p=polyfit(x,y,n);
       [p,s]=polyfit(x,y,n);
       [p,s,mu]=polyfit(x,y,n)
       

       输入：

       x	观测点，数值向量
       y	观测值，数值向量
       n	多项式拟合的程度，即利用n次多项式进行拟合

       输出：

       p	系数估计，数值向量，元素按照次数从小到大排列
       s	残差估计的结构，$\text{S.R}$ （三角形因子），$\text{S.df}$（自由度），$\text{S.normr}$ （残差范数）
       mu	两元素向量 ，$mu(1)$ ：$x$ 的平均值，$mu(2)$ ：$x$ 的标准差

     • $\text{Matlab}$实现——可视化

       y=polyval(p,x)
       [y,delta]=ppolyval(p,x,s)
       y=polyval(p,x,[],mu)
       [y,delta]=polyval(p,x,s,mu)
       

       输入：

       x	观测点，数值向量
       p	多项式系数估计值
       s	残差估计的结构
       mu	两元素向量

       输出：

       y	拟合函数值
       delta	残差的标准偏差

       说明：

       • 多项式计算时用到 $mu$ ，即表明利用 $\hat{x}=\frac{x-mu(1)}{mu(2)}$ 代替 $x$ 。
       • 多项式拟合的次数不宜过高，以免出现 Runge 现象。

     • 常用的数据变换方式

       曲线类型	变换	直线方程
       幂函数 $y=a{x}^b$	$Y=\log (y),X=\log (x)$	$Y=\log (a)+bX$
       指数函数 $y=a{e}^{bx}$	$Y=\log (y),X=x$	$Y=\log (a)+bX$
       双曲函数 $\frac{1}{y}=a+\frac{b}{x}$	$Y=\frac{1}{y},X=\frac{1}{x}$	$Y=a+bX$
       $S$ 型函数 $y=\frac{1}{a+b{e}^{-x}}$	$Y=\frac{1}{y},X=e^{-x}$	$Y=a+bX$

     • 参考线

       直线 $y=ax+b$ ：$\text{refline(a,b)}$

       多项式曲线 $y=p_n{x}^n+\cdots +p_{1}x+p_0$ ：$\text{refcurve}([p_n,\cdots ,p_1,p_0])$

   • 一般曲线拟合

     曲线形式 $y=f(a,x)$ 已知，其中 $a=[a_1,a_2,\cdots ]$ 为未知的待定系数。

     • $\text{lsqcurvefit}$

- $\text{lsqnonlin}$

- $\text{nlinfit}$

 - 曲线拟合工具箱

2. Markov 预测

   Markov 预测，是一种预测事件发生的概率的方法。它基于 Markov 链，根据事件的目前状况预测其将来各个时刻（或时期）变动状况的一种预测方法。其基本要求时状态转移概率矩阵必须具有一定的稳定性，且必须建立在大量的统计数据基础上。

   • 相关概念

     状态：某一事件在某个时刻（或时期）出现的某种结果，记为 $E_i,i=1,2,\cdots ,n$ 。

     状态转移：事件的发展从一种状态转变为另一状态，记为 $E_i\to E_j$ 。

     Markov 过程：在事件的发展过程中，若每次状态的转移都仅与前一时刻的状态有关，与过去的状态无关，或者说状态转移是无后效性的。

     状态转移概率：在事件的发展变化过程中，从有一种状态出发，下一时刻转移到其他状态的可能性，记为 P { E i ∣ E j } = ⋅ p i j \mathbb{P}\{E_i|E_j\}\overset{\cdot}= p_{ij} P{Ei​∣Ej​}=⋅pij​ 。

     状态转移概率矩阵： P = { p i j } n × n P=\{p_{ij}\}_{n\times n} P={pij​}n×n​ ，且满足
     $$\begin{gathered} 0\le p_{ij}\le 1,i,j=1,2,\cdots ,n \\ \sum _{i=1}^n{p}_{ij}=1 \end{gathered}$$

   • 状态概率的计算

     用 ${\pi }_j(k)$ 表示事件在第 $k$ 时刻处于状态 $E_j$ 的概率。显然，
     $$\sum _{j=1}^m{\pi }_j(k)=1$$
     根据 Markov 过程的无后效性及 Bayes 条件概率公式，有
     $${\pi }_j(k)=\sum _{i=1}^n{\pi }_i(k-1)p_{ij},i,j=1,\cdots ,n$$
     记 $\pi (k)=[{\pi }_1(k),\cdots ,{\pi }_n(k)]$ 为第 $k$ 时刻的状态概率向量，则
     $$\pi (k)=\pi (k-1)P=\cdots =\pi (0)P^k$$
     其中 $\pi (0)$ 为初始状态概率向量。

   • 终极状态概率预测
     $$\underset{k\to \infty }{\lim }\pi (k)=[\underset{k\to \infty }{\lim }{\pi }_1(k),\cdots ,\underset{k\to \infty }{\lim }{\pi }_n(k)]=[{\pi }_1,\cdots ,{\pi }_n]=\pi$$
     计算方法：解下列方程组即可得到终极状态概率
     $$\{\begin{array}{l}\pi =\pi P\\ \sum _{i=1}^n{\pi }_i=1\end{array}$$

   • 拓展内容

     • **正则链：**从任意状态出发经过有限次转移都能达到另外的任意状态

       • 定义：

         一个有 $k$ 个状态的马氏链如果存在正整数 $N$ ，使从任意状态 $i$ 经 $N$ 次转移都以大于零的概率到达状态 $j(i,j=1,\cdots ,k)$ ，则称为正则链。

       • 检验定理：

         若马氏链的转移矩阵为 $P$ ，则它是正则链的充要条件是，存在正整数 $N$ ，使 $P^N>0$ （指 $P^N$ 的每一元素大于零）

       • 极限定理：

         正则链存在唯一的极限状态概率 $w=(w_1,\cdots ,w_k)$ ，使得当 $n\to \infty$ 时状态概率 $a(n)\to w$ ，$w$ 与初始状态概率 $a(0)$ 无关，$w$ 又称稳态概率，满足
         $$\begin{gathered} wP=w \\ \sum _{i=1}^k{w}_i=1 \end{gathered}$$
         当转移矩阵 $P$ 给定后，我们可以通过联立上述两式构成方程组，解出极限状态概率 $w$ 的具体值。

         通过上面的陈述，我们可以看出 $\underset{n\to \infty }{\lim }{P}^n$ 存在，记作 $P^{\infty }$ ，并且 $P^{\infty }$ 的每一行都是稳态概率 $w$ 。如果记 P ∞ = { p i j ( ∞ ) } P^\infty=\{p_{ij}^{(\infty)}\} P∞={pij(∞)​} ，那么有 $p_{ii}^{(\infty )}=w_i$ 。

     • 首达概率与平均转移次数

       从状态 $i$ 出发经过 $n$ 次转移，第一次到达状态 $j$ 的概率称为 $i$ 到 $j$ 的首达概率，记作 $f_{ij}(n)$ 。于是
       $${\mu }_{ij}=\sum _{i=1}^{\infty }n{f}_{ij}(n)$$
       为由状态 $i$ 第一次到达状态 $j$ 的平均转移次数。特别的，${\mu }_{ii}$ 是状态 $i$ 首次返回的平均转移次数。

       • 定理：（${\mu }_{ii}$ 与 $w$ 的相互转换）

         对于正则链，
         $${\mu }_{ii}=\frac{1}{w_i}$$
         

         解：

       ​ （1）由题意，状态转移矩阵 $P=\left(\begin{array}{ccc}0.6& 0.2& 0.2\\ 0.1& 0.8& 0.1\\ 0.1& 0.4& 0.5\end{array}\right)$，初始状态为$\pi (0)=[0.40.40.2]$

       ​ $\therefore \pi (1)=\pi (0)P=[0.30.480.22]$

       ​ $\therefore$ 2018年，预测B厂家的市场占有率最高

       ​ （2）由题意，设终极状态概率向量为$\pi =\pi P,{\sum }_{i=1}^3{\pi }_i=1$可得，
       $$\{\begin{array}{l}{\pi }_1=0.6{\pi }_1+0.1{\pi }_2+0.1{\pi }_3\\ {\pi }_2=0.2{\pi }_1+0.8{\pi }_2+0.4{\pi }_3\\ {\pi }_3=0.2{\pi }_1+0.1{\pi }_2+0.5{\pi }_3\\ {\pi }_1+{\pi }_2+{\pi }_3=1\end{array}解得\{\begin{array}{l}{\pi }_1=0.2\\ {\pi }_2=0.6\\ {\pi }_3=0.2\end{array}$$
       ​ $\therefore$ 很长时间后，B厂家的市场占有率最高

​ 3.哎，我就是玩（备考请跳过 😘）

• 数据标准化

  • 原因：

    • 属性不同：有的越大越好，有的越小越好，有的在某个区间最佳
    • 量纲不同：计量单位不同，不同列的数值大小差异很大
    • 极差不同：最大值最小值之差不同，不具有可比性

  • 数据标准化方法：

    vecStd.m		向量归一化，转换到区间[0,1]
    linearStd.m		线性比例变换法，分正向型和反向型分别转换到区间[0,1]
    minmaxStd.m		极差归一化，分正向型和反向型分别转换到区间[0,1]
    optValue.m		最优值为给定数值的标准化，转换到区间[0,1]
    optInterval.m	最优值为给定区间的标准化，转换到区间[0,1]
    zscoreStd.m		Z标准变换法，化为标准正态分布数值
    

第六章：图论

最短路径的 Dijkstra 算法，Floyd 算法；最小生成树的 Prime 算法只考其一。因此下面着重介绍这三种算法及其原理。

1. 最短路径

   • 基本概念：

     顶点集：顶点组成的集合，记为 $V$

     边集：顶点之间的连接组成的集合，记为 $E$

     图：顶点集和边集组成，记为 $G=(V,E)$

     网络：被赋予具体的含义和权的图，记为 $G=(V,E,W)$

     权值：任意两顶点 $u,v$ 之间的权重，记为 $⟨u,v⟩$

     有（无）向图：边的连接时有（无）方向的图

   • Dijkstra 算法（只适用于正权边的图）

     • 算法原理：

       将图中顶点集分为两组：

       已求出最短路径的顶点集合 $S$ ：初始时 $S$ 只有一个起始点（或称为源点 $v$ ），以后每求得一条最短路径，就将其加入集合 $S$ ，直到全部顶点都加入 $S$ 中，算法结束。

       未确定最短路径的顶点集合 $U$ ：按最短路径长度的递增次序依次把 $U$ 中顶点加入 $S$ 中，每次保持从源点 $v$ 到 $S$ 中各顶点的最短路径长度不大于从 $v$ 到 $U$ 中任何顶点的最短路径长度。

       $S$ 中的顶点的距离：从 $v$ 到此顶点的最短路径长度

       $U$ 中的顶点的距离：从 $v$ 到此顶点只包含 $S$ 中的顶点为中间顶点的当前最短路径长度

     • 算法步骤：

       • 初始时， S = { v } S=\{v\} S={v} ，$v$ 的距离为 $0$ ； U = { 其余顶点 } U=\{\text{其余顶点}\} U={其余顶点} ，若 $v$ 与 $U$ 中顶点 $u$ 有边，则距离为其权值 $⟨u,v⟩$ ，否则设为 $\infty$ ；
       • 从 $U$ 中选取一个距离 $v$ 最小的顶点 $k$ ，并加入 $S$ 中；
       • 以 $k$ 为新考虑的中间点，修改 $U$ 中各顶点的距离，若从起始点 $v$ 到顶点 $u$（经过顶点 $k$ ）比原来距离（不经过顶点 $k$ ）短，则修改顶点 $u$ 的距离值为顶点 $k$ 的距离加上 $⟨k,u⟩$ 的长度；
       • 重复第二步和第三步，直到所有顶点都包含在 $S$ 中。

     例题：用Dijkstra 法求出以A 为起始点到其余各点的最短路径及距离。
     
     

• Floyd 算法（适用于所有图）

  • 算法原理：

    假设 $\text{Dis}(i,j)$ 为任意两个节点 $i,j$ 之间的最短路径的距离，对于每一个节点 $k$ ，检查
    $$\text{Dis}(i,k)+\text{Dis}(k,j)<\text{Dis}(i,j)$$
    是否成立，若成立，则设置 $\text{Dis}(i,j)=\text{Dis}(i,k)+\text{Dis}(k,j)$ ，依次类推，当遍历完所有节点 $k$ ，$\text{Dis}(i,j)$ 中记录的就是节点 $i$ 到节点 $j$ 的最短路径距离。

  • 算法步骤：

    • 从任意一条单边路径开始，所有两点之间的距离是边的权值，若没有相连则权值为无穷大；
    • 对于每一对顶点 $(u,v)$ ，检查是否存在顶点 $w$ 使得 $u\to w\to v$ 的距离更短，若是则更新它。

2. 最小生成树

   • 基本概念：

     树：连通且不含圈的无向图，记为 $T$

     树枝：树中的边

     树叶：树中度数为 $1$ 的顶点

     生成树：若 $T$ 是包含图 $G$ 的全部顶点的子图，且又是树

     权：$T$ 中全部边上的权数之和，记为 $w(T)=\sum _{e\in E_1}w(e)$

     最小生成树：所有生成树的权中最小者

   • Prim 算法

     • 算法原理：

       给定连通赋权图 $G=(V,E,W)$ ，其中 $W$ 为邻接矩阵，设两个集合 $P$ 和 $Q$ ，其中 $P$ 用于存放 $G$ 的最小生成树中的节点，集合 $Q$ 存放 $G$ 的最小生成树中的边。令集合 $P$ 的初值为 P = { v 1 } P=\{v_1\} P={v1​} ，集合 $Q$ 的初值为 Q = { } Q=\{\} Q={} 。

     • 算法步骤：

       • 假设集 { v 1 } \{v_1\} {v1​} 是起始点，初始化 P = { v 1 } , Q = { } P=\{v_1\},Q=\{\} P={v1​},Q={}

       • while $P\sim =V$

         找到最小边 $⟨p,v⟩$ ，其中 $p\in P,v\in V\setminus P$

         P = P + { v } P=P+\{v\} P=P+{v}

         Q = Q + { ⟨ p , v ⟩ } Q=Q+\{\langle p,v \rangle\} Q=Q+{⟨p,v⟩}

         end

   • Kruskal 算法

     • 算法原理：

       设无向连通赋值图$G=(V,E)$，令$G$的最小生成树为$T=(U,TE)$，其初始状态为$U=V$和 T E = { } TE=\{\} TE={}，这样$T$中各顶点各自构成一个连通分量。按照权值从小到大的顺序$E_{(1)}\le E_{(2)}\le \cdots \le E_{(m)}$,有选择性地考察，最终使得$T$的连通分量个数为1，即得到$G$的最小生成树。

     • 算法步骤：

       • 初始化 U = V , T E = { } U=V,TE=\{\} U=V,TE={}

       • 将$E$中的元素按权值从小到大重新排列，得$E_{(1)}\le E_{(2)}\le \cdots \le E_{(m)}$，且对应的顶点为

         $V_{(1)},V_{(2)},\cdots ,V_{(m)}$

       • for i= 1 : m

         ​ while $f(U)\sim =1$（$f(X)$表示顶点集$X$中的连通分量个数）

         ​ if g ( T E + { E ( i ) } ) g(TE+\{E_{(i)}\}) g(TE+{E(i)​})（$g(X)$判断边集$X$构成的树中是否有圈）

         ​ continue;

         ​ else

         ​ T E = T E + { E ( i ) } ; TE=TE+\{E_{(i)}\}; TE=TE+{E(i)​};

         ​ end

         ​ end

         end

         ​

3. 指派问题

   某公司准备分派 $n$ 个工人 $X_1,\cdots ,X_n$ 做 $n$ 件工作 $Y_1,\cdots ,Y_n$ ，但由于任务性质和个人专长不同，因此各人完成每件工作的效率也就不同，于是产生了一个问题：应指派哪个人去完成哪一件工作，使工作效率最高？

   用效率矩阵表示：
   $$C=\left(\begin{array}{ccc}{c}_{11}& \cdots & c_{1n}\\ ⋮& \ddots & ⋮\\ c_{n1}& \cdots & c_{nn}\end{array}\right)$$
   其中 $c_{ij}$ 表示指派第 $i$ 个人去完成第 $j$ 项任务时的效率（或时间、成本等），$i,j=1,\cdots ,n$ 。令
   $$x_{ij}=\{\begin{array}{l}1,\text{分配第 i 个人去完成第 j 项任务}\\ 0,\text{不分配第 i 个人去完成第 j 项任务}\end{array}$$
   指派问题的一般数学模型为：
   $$\begin{gathered} \min Z=\sum _{i=1}^n\sum _{j=1}^n{c}_{ij}{x}_{ij} \\ s.t.\sum _{j=1}^n{x}_{ij}=1(一项工作只能由一个人完成) \\ \sum _{i=1}^n{x}_{ij}=1(一个人只能完成一项工作) \\ {x}_{ij}=0\text{或}1i,j=1,\cdots ,n \end{gathered}$$

   • Hungary 算法原理

     • 定理 $1$ ：

       若从指派问题效率矩阵 { c i j } n × n \{c_{ij}\}_{n\times n} {cij​}n×n​ 的每一行元素中分别减去（或加上）一个常数 $u_i$（称为该行的位势），从每一列元素中分别减去（或加上）一个常数 $v_j$（称为该列的列势），得到一个新的效率矩阵 { b i j } n × n \{b_{ij}\}_{n\times n} {bij​}n×n​ ，则 { b i j } n × n \{b_{ij}\}_{n\times n} {bij​}n×n​ 对应的指派问题与原分派问题有相同的最优解，只是最优值相差一个常数。

     • 定理 $2$ ：

       指派问题效率矩阵 { c i j } n × n \{c_{ij}\}_{n\times n} {cij​}n×n​ 中独立零元素的最多个数等于能覆盖所有零元素的最少线数。

       只有一个零元素的行或列将该零元素标注为独立零元素。

     • **例：**利用 Hungary 算法分析如下效率矩阵的指派问题
       $$C=\left(\begin{array}{ccccc}12& 7& 9& 7& 9\\ 8& 9& 6& 6& 6\\ 7& 17& 12& 14& 9\\ 15& 14& 6& 6& 10\\ 4& 10& 7& 10& 9\end{array}\right)$$
       第一步：对矩阵 $C$ 做行规约（每行元素减去该行的最小元素）、列规约（每列元素减去该列的最小元素），使得每一行和每一列都含有零：
       $$B=\left(\begin{array}{ccccc}5& 0& 2& 0& 2\\ 2& 3& 0& 0& 0\\ 0& 10& 5& 7& 2\\ 9& 8& 0& 0& 4\\ 0& 6& 3& 6& 5\end{array}\right)$$
       第二步：从最少“ $0$ ”数的行或列开始，将“ $0$ ”圈起来，并划掉它所在行和列的其他 “ $0$ ”；反复这样做，直到所有“ $0$ ”圈起来或被划掉为止：

       

       圈起的“ 0 ”的数量$m=4<5=n$，故不是最优解，也就是说如果在这一步$m=n$，后面就不用算了

       第三步：在没有圈起“ $0$ ”所在行打“√”；在打“√”行中所有“ $\varnothing$”所在的列打“√”，在打“√”列中含有圈起“ $0$ ”的行上打“√”，反复执行这两步，知道不能打“√”为止，然后用直线划去打“√”的列和不打“√”的行，这样没有划去的行构成调整行，划去的列构成调整列：

       

       第四步：在所有调整行中寻找没被划线的最小的元素作为调整量，将调整行各元素减去调整量，对调整列中各元素加上调整量，然后执行圈“ $0$ ”和划“ $0$ ”的操作，并循环以上步骤，直到圈起的“ $0$ ”数等于 $n$ 为止：
       2. - -

       

       可以看到调整行中最小元素为2.

       注意: 要明确调整列和调整行

       

       在红色区域里的元素加调整量, 在蓝色区域里的元素减调整量, 行列交叉处的0不变

       第五步：留下圈0的元素为1, 剩下所有元素为零, 即得到最优矩阵

       

       得到一个最优解：
       $$x_{12}=x_{23}=x_{35}=x_{43}=x_{51}=1$$
       此时效率成本为：
       $$7+6+9+6+4=32$$

       即: 第一个人做第2个任务; 第二个人做第3个任务; 第三个人做第5个任务; 第四个人做第4个任务;
       第五个人做第1个任务.

第七章：动态模型

1. 微分方程模型

   • 适用范围

     • 物体的运动、震动、受力形变（物体的运动轨迹、运动目标的跟踪与拦截、 建筑的防震/防强风设计、物体受力形变与控制等）；
     • 动植物/微生物的数量或密度的变化（植物生长、多种生物间相互依存、传 染病的传播与控制等）；
     • 物质、能量扩散与传递（粉尘、烟雾、化学物质在空气、水、土壤中的扩散 与沉积，化学反应过程的描述，热传导问题等）；
     • 消费品在市场上的销售过程；
     • 信息的扩散与传播。

   • 基本概念

     解析解（精确解）：适用于线性和少量非线性系统；

     数值解（近似解）：适用于多数线性和非线性系统，但不能对系统的行为提供一个定性解释；

     定性解：用定性理论（适用于二维或三维系统）和稳定性理论（适用于高维系统）分析系统在局部和全局的动态行为。

   • 建模方法

     根据规律建模；用微元法建模；用模拟近似法建模

   • 建模原理

     转化；准确性和总体特征；约束条件

     注意事项：

     • 所建立的方程或方程组应满足守恒定律；
     • 尽量简化方程；
     • 量纲一致；
     • 比较理想化的建模方法，适用于定性讨论或精度要求不高的情形。

   • 微分方程的解析解

     S=dsolve(eqn)
     S=dsolve(eqn,cond)
     S=dsolve(eqn,cond,Name,Value)
     [y1,……,yN]=dsolve(problem)
     

     输入：

     eqn	微分方程或方程组
     cond	初值或边界条件
     Name，Value	名称-值对参数

     输出：

     S	微分方程的解
     y1，……，yN	多个变量的微分方程的解

     

   • 微分方程的数值解

     

     [t,y]=ode45(odefun,tspan,y0)
     [t,y]=ode45(odefun,tspan,y0,options)
     [t,y,te,ye,ie]=ode45(odefun,tspan,y0,options)
     sol=ode(problem)
     

     输入：

     odefun	待解方程
     tspan	积分区间
     y0	初始值
     option	选项结构体

     输出：

     t	时间演化点
     y	数值解
     te	事件的时间
     ye	事件时间处的解
     ie	指定发生了哪个事件
     sol	用于计算的结构体

   • 传染病模型

     为了简单起见，假设该地区总人口不变，记为 $N$ 。

     • 模型目标：描述传染病的传播过程，分析受感染人数的变化规律，预报传染病高峰到来的时刻，得到控制和根除预防传染病的方法。

     • 基本假设：传染病是由病人通过“接触”健康人进行传播的。疾病流行区域内的人分为三类：

       易感人群（Susceptible individuals）— S 类

       病人（Infected individuals）— I 类

       移出者（Recovered individuals）— R 类

     • 传染病模型 $1$ ：SI 模型

       • 模型假设：

         假设 $1$ ：只考虑传染病（病人 $I$ ）和不传染病（易感 $S$ ）两类人群，且其数量只与时间有关。记 $S(t)$ 为 $t$ 时刻健康人占总人口的比例，$I(t)$ 为 $t$ 时刻病人占总人口的比例；

         假设 $2$ ：人群数量足够大，只考虑传播过程中的平均效应，即函数 $S(t)$ 和 $I(t)$ 可以视为连续且可微的；

         假设 $3$ ：每个 $I$ 类的人每天“有效接触”的人数（包括宾根与健康人）为常数 $\lambda$ ，即传染率，反映该地区的卫生水平；

         假设 $4$ ：不考虑出生与死亡，以及人群的迁入迁出因素。

       • 模型构建：

         考虑到 $t$ 到 $t+\Delta t$ 时间内病人人数的变化，根据假设 $1$ ，可得在 $\Delta t$ 时间内受感染的人数为
         $$NI(t+\Delta t)-NI(t)=NI(t)\cdot \lambda S(t)\Delta t$$
         令 $\Delta t\to 0$ ，可得
         $$\{\begin{array}{l}\frac{dI(t)}{dt}=\lambda I(t)S(t)=\lambda I(t)(1-I(t))\\ I(0)=I_0\end{array}$$
         注： SI 模型的缺点很明显，当传播时间足够长时，该地区的所有人都会感染。

     • 传染病模型 $2$ ：SIS 模型

       在 SI 模型的四条假设的基础上，增加一条假设：

       假设 $5$ ：病人以固定的比率痊愈，再次成为易感人群。每天被治愈的病人数占病人总数的比例为 $\mu$ 。

       • 模型构建

         考虑到 $t$ 到 $t+\Delta t$ 时间内病人人数的变化，根据假设 $1$ ，可得在 $\Delta t$ 时间内受感染的人数为
         $$NI(t+\Delta t)-NI(t)=NI(t)\cdot \lambda S(t)\Delta t-N\cdot \mu I(t)\Delta t$$
         令 $\Delta t\to 0$ ，可得
         $$\{\begin{array}{l}\frac{dI(t)}{dt}=\lambda I(t)(1-I(t))-\mu I(t)\\ I(0)=I_0\end{array}$$
         若令 $\sigma =\lambda /\mu$ ，即一个病人在整个患病期间有效接触的平均人数，称为接触数。此时，得到 SIS 模型的另一个等价形式：
         $$\{\begin{array}{l}\frac{dI(t)}{dt}=\lambda I(t)[-I(t)+(1-\frac{1}{\sigma })]\\ I(0)=I_0\end{array}$$
         注：
         $$\underset{t\to \infty }{\lim }I(t)=\{\begin{array}{l}1-1/\sigma ,\sigma >1\\ 0,\sigma <1\end{array}$$

     • 传染病模型 $3$ ：SIR 模型

       引入 $R$ 类人群（免疫者或者治愈后免疫），分别记 $S(t),I(t),R(t)$ 为病人、易感人群、移出者在总人口中所占的比例，即
       $$S(t)+I(t)+R(t)=1$$

       • 模型构建

         考虑到 $t$ 到 $t+\Delta t$ 时间内病人和移出者人数人数的变化，可得在 $\Delta t$ 时间内受感染的人数为
         KaTeX parse error: No such environment: eqnarray at position 8: \begin{̲e̲q̲n̲a̲r̲r̲a̲y̲}̲ NI(t+\Delta t)…
         令 $\Delta t\to 0$ ，可得
         $$\{\begin{array}{l}\frac{dI(t)}{dt}=\lambda I(t)S(t)-\mu I(t)\\ \frac{dR(t)}{dt}=\mu I(t)\end{array}$$
         这也表明
         $$\frac{dS(t)}{dt}=-\lambda I(t)S(t)$$
         注：

         • $\underset{t\to \infty }{\lim }I(t)=0$ ，表明传染病最终将消失；

         • 若 $S_0>1/\sigma$ ，则 $I(t)$ 先升后降至 $0$ ，表示传染病蔓延；

         • 若 $S_0<1/\sigma$ ，则 $I(t)$ 单调下降至 $0$ ，表示传染病不会蔓延；

       

   • 生物种群模型

     • 生物种群模型 $1$ :种群竞争模型

       考虑甲、乙两个种群在争夺同一食物来源和生存空间而进行生存竞争时的变化 规律。分别用 $N1(t)$ 和 $N2(t)$ 表示在 $t$ 时刻两个种群的数量。 假设种群甲单独存在时，收到食物和生活空间的限制，遵循种群变化的 Logistic 规律，即种群平均增长率为
       $$r=r_1\left(1-\frac{N_1}{K_1}\right)$$
       现增加一个竞争对手，并假设所导致的平均增长率的降低量正比于竞争对手的 数量，则种群甲的平均增长率变为：
       $$r=b_1-a_{11}{N}_1-a_{12}{N}_2$$
       其中，$b_1$ 表示种群甲单独存在时的平均增长率，$a_{12}$ 表示种群乙对种群甲的竞争力。于是，可得两种群竞争模型
       $$\{\begin{array}{l}\frac{d{N}_1}{dt}=N_1(b_1-a_{11}{N}_1-a12N_2)\\ \frac{d{N}_2}{dt}=N_2(b_2-a_{21}{N}_1-a22N_2)\end{array}$$

     • 生物种群模型 $2$ : 食饵-捕食者模型

       设种群甲依靠丰富的自然资源生长，种群乙依靠捕食种群甲为生，考虑两个种 群的变化规律。由于资源丰富，故种群甲单独存在时的平均增长率为 $r=b_1$ 。 存在捕食者是，假设由于被捕食者造成种群甲的平均增长率的降低量正比于捕 食者的数量，此时种群甲的平均增长率为：
       $$r=b_1-a_{12}{N}_2$$
       其中 $a_{12}$ 表示单位捕食者造成食饵增长率的减少量。

       由于种群乙依靠捕食种群甲为生，故种群乙单独存在时的增长率为 $r=-b_2$ 。 由于捕食，假设种群乙的平均增长率正比于食饵数量，则种群乙的增长率为：
       $$r=-b_2+a_{21}{N}_1$$
       于是，可得两种群的食饵-捕食者模型
       $$\{\begin{array}{l}\frac{d{N}_1}{dt}=N_1(b_1-a_{12}{N}_2)\\ \frac{d{N}_2}{dt}=N_2(-b_2+a_{21}{N}_1)\end{array}$$
       注：

       • 存在两个不稳定的平衡态 $(0,0),(b_2/a_{21},b_1/a_{12})$
       • 初值微小扰动，轨线发生变化

     • 生物种群模型 $3$ : 修正的食饵-捕食者模型

       假设在捕食者与食饵没有相遇时，种群甲乙都依靠有限的自然资源生长，即根据 Logistic 规律，种群甲乙单独存在时的平均增长率分别为：
       $$b_1=r_1\left(1-\frac{N_1}{K_1}\right),b_2=r_2\left(1-\frac{N_2}{K_2}\right)$$
       假设由于捕食造成种群甲的平均增长率降低 $a1N2$ ，种群乙的平均增长率增 加 $a2N1$ ，则得到修正的食饵- 捕食者模型：
       $$\{\begin{array}{l}\frac{d{N}_1}{dt}=r_1\left(1-\frac{N_1}{K_1}\right)-a_1{N}_1{N}_2\\ \frac{d{N}_2}{dt}=r_2\left(1-\frac{N_2}{K_2}\right)+a_2{N}_1{N}_2\end{array}$$
       该模型有四个平衡态 $(x^{\ast },y^{\ast })$ ，若正平衡态 $(x^{\ast }>0,y^{\ast }>0)$ 是渐近稳定的，则该模型的解都将随时间趋于平衡态，这表明捕食者与食饵可以长期共存，符合实际。

       • 二阶自治微分方程的平衡点与稳定性

         

         注意到，
         $$p=-tr(A){|}_{P^{\ast }(x_1^{\ast },x_2^{\ast })},q=\det (A)=(f_{x_1}{g}_{x_2}-f_{x_2}{g}_{x_1}){|}_{P^{\ast }(x_1^{\ast },x_2^{\ast })}$$

       考虑修正的食饵-捕食者模型
       $$\{\begin{array}{l}\frac{d{N}_1}{dt}=r_1\left(1-\frac{N_1}{K_1}\right)-a_1{N}_1{N}_2\stackrel{\Delta }{=}f(N_1,N_2)\\ \frac{d{N}_2}{dt}=r_2\left(1-\frac{N_2}{K_2}\right)+a_2{N}_1{N}_2\stackrel{\Delta }{=}g(N_1,N_2)\end{array}$$

       • 首先，计算平衡点：
         $$\{\begin{array}{l}f(N_1,N_2)=0\\ g(N_1,N_2)=0\end{array}⟺\{\begin{array}{l}{N}_1=0或r_1\left(1-\frac{N_1}{K_1}\right)-a_1{N}_1{N}_2=0\\ N_2=0或r_2\left(1-\frac{N_2}{K_2}\right)+a_2{N}_1{N}_2=0\end{array}$$
         解得
         $$\{\begin{array}{l}{N}_1=0\\ N_2=0\end{array}或\{\begin{array}{l}{N}_1=0\\ N_2=K_2\end{array}或\{\begin{array}{l}{N}_1=K_1\\ N_2=0\end{array}或\{\begin{array}{l}{N}_1=\frac{r_2{K}_1(r_1-K_2{a}_1)}{r_1{r}_2+K_1{K}_2{a}_1{a}_2}\\ N_1=\frac{r_1{K}_2(r_2+K_1{a}_2)}{r_1{r}_2+K_1{K}_2{a}_1{a}_2}\end{array}$$
         因此有平衡点：
         $$P_1(0,0),P_2(0,K_2),P_3(K_1,0),P_4(\frac{r_2{K}_1(r_1-K_2{a}_1)}{r_1{r}_2+K_1{K}_2{a}_1{a}_2},\frac{r_1{K}_2(r_2+K_1{a}_2)}{r_1{r}_2+K_1{K}_2{a}_1{a}_2})$$
         其中，当 $r_1\ge K_2{a}_1$ 时，$P_4$ 才有意义，且等号成立时退化为 $P_2$ 。

       • 其次，计算对应的 Jacobi 矩阵：
         $$A=\left(\begin{array}{cc}{f}_{N_1}& f_{N_2}\\ g_{N_1}& g_{N_2}\end{array}\right)=\left(\begin{array}{cc}\frac{-2r_1{N}_1+K_1(r_1-N_2{a}_1)}{K_1}& -a_1{N}_1\\ a_2{N}_2& \frac{-2r_2{N}_2+K_2(r_2+N_1{a}_2)}{K_2}\end{array}\right)$$

       • 列表，判断稳定条件：

       

       综上，当 $r_1<K_2{a}_1$ 时，因捕食者能力过高，$P_2$ 稳定，即食饵灭绝，捕食者达到最大容量；当 $r_1>K_2{a}_1$ ，$P_4$ 稳定，即二者共存，分别趋向非零的有限值。