数学模型与python基础——稳定状态数学模型

数学模型与python基础（4）

之前讲了如何利用微分方程建立的数学模型来描述一些量动态变化的过程，微分方程能够具体描绘出每一时刻的各个量的值，但是在一些情况下，我们并不需要知道一些量每时每刻的变化，只想知道其稳定状态的情况以及趋势性。这种情况下可以利用微分方程稳定性理论，而不需要对微分方程进行具体求解。

微分方程稳定性理论

自治

如果一个常微分方程或方程组不含时间变量，即
$$\frac{dx}{dt}=F(x,t)=\left[\begin{array}{c}{f}_1(x,t)\\ ⋮\\ f_N(t)\end{array}\right]=F(x)$$
则称该常微分方程（组）是自治的，该系统为自治系统。自治表示该系统的变化速度与该系统所处的状况有关，与时间无关，假如与时间有关，则称非自治的。

相空间

假定有一个系统$\frac{dx}{dt}=F(x)$，则定义以点$(x_1,x_2,x_3,...,x_n)$为坐标的空间$R^n$，当n=2时，称相空间为相平面。满足微分方程（组）的点组成的点集称为系统的轨线，所有轨线在相空间的分布图称为相图。

奇点（平衡点）与稳定性

相空间中满足$F(x_0)=0$的点$x_0$称为系统的奇点或平衡点。奇点可以是孤立的，也可以是连续的点集。

设$F(x)$为实解析函数（处处可微），且$x_0$为系统$\frac{dx}{dt}=F(x)$的奇点。若$F(x)$在点$x_0$处的雅各比矩阵$J\left(x_0\right)=\left[\frac{\partial f_i}{\partial x_j}\right]$是非奇异的，则$x_0$是该系统的孤立奇点。

如果$x_0$是$\frac{dx}{dt}=F(x)$的奇点，如果对于任意一个给定的$ϵ>0$，存在一个$\delta >0$，使得如果$|x(0)-x_0|<\delta$，则$|x(t)-x_0|<ϵ$对所有的t都成立，则称奇点$x_0$是稳定的。即任何初始条件在平衡态附近的轨迹均能维持在平衡态附近，则称奇点是稳定的。否则称奇点是不稳定的。

如果一个稳定的奇点满足$\underset{t\to \infty }{\lim }|x(t)-x_0|=0$，则称该奇点是渐进稳定的。

线性常系数系统的系数矩阵A对应的行列式不等于零时，即$detA\ne 0$时，(0,0)是系统唯一的平衡点，其稳定性由特征方程 $det(A-\lambda I)=0$ 的特征根$\lambda$决定

一阶方程的平衡点与稳定性

设有方程
$$\dot{x}(t)=f(x)$$
f(x)无自变量t，因此该方程为自治方程，该方程对应的系统为自治系统

方程$f(x)=0$的解为系统奇点（平衡点）

判断该系统在奇点处是否稳定的方法是，将$f(x)$在平衡点$x_e$处做泰勒展开，只取一阶项，方程近似为
$$\dot{x}(t)=f^{\prime }(x_e)(x-x_e)$$
该方程称为原方程的近似线性方程。

若$f^{\prime }(x_e)<0$则$x_e$对于原方程和近似线性方程都是稳定的

二阶方程的平衡点与稳定性

一个二阶方程可以用两个一阶方程组成的方程组表示为
$$\{\begin{array}{l}{\dot{x}}_1=f(x_1,x_2)\\ \\ {\dot{x}}_2=g(x_1,x_2)\end{array}$$
$f(x_1,x_2),g(x_1,x_2)$均不含t，因此都是自治方程

如下方程的解为系统的奇点（平衡点）
$$\{\begin{array}{l}{\dot{x}}_1=f(x_1,x_2)=0\\ \\ {\dot{x}}_2=g(x_1,x_2)=0\end{array}$$
对于线性常系数方程（非齐次方程通过平移换元来齐次）
$$\{\begin{array}{l}{\dot{x}}_1=a_1{x}_1+a_2{x}_2\\ \\ {\dot{x}}_2=b_1{x}_1+b_2{x}_2\end{array}$$
二阶线性常系数方程显然只有一个平衡点(0,0)

系数矩阵为
$$A=\left[\begin{array}{cc}{a}_1& a_2\\ b_1& b_2\end{array}\right]$$
若A的行列式不为0

则可根据微分方程组的特征方程特征根决定
$$det(A-\lambda I)=0$$
即
$$\{\begin{array}{l}{\lambda }^2+p\lambda +q=0\\ \\ p=-(a_1+b_2)\\ \\ q=detA\end{array}$$
特征根${\lambda }_{1,2}=\frac{1}{2}(-p\pm \sqrt{p^2-4q})$

当${\lambda }_1\ne {\lambda }_2$时，方程一般解为$c_1{e}^{{\lambda }_{1}t}+c_2{e}^{{\lambda }_{2}t}$

当${\lambda }_1={\lambda }_2$时，方程一般解为$c_1{e}^{{\lambda }_{1}t}+c_{2}t{e}^{{\lambda }_{1}t}$

当${\lambda }_{1,2}$均有负实部的时候，(0,0)为稳定平衡点，任意一个有正实部的时候，(0,0)为不稳定平衡点

因此$p>0\text{且}q>0$时，平衡点稳定，$p<0\text{或}q<0$时，平衡点不稳定

并且对于特征根的分布和奇点类型的关系有如下结论（具体不同类型的奇点周边的相轨迹图形，可以参见我的自动控制的非线性系统分析博客）

（1）${\lambda }_1<{\lambda }_2<0$时，原点是稳定节点

（2）${\lambda }_1={\lambda }_2<0$时，原点是稳定退化节点

（3）${\lambda }_1>{\lambda }_2>0$时，原点是不稳定节点

（4）${\lambda }_1={\lambda }_2>0$时，原点是不稳定退化节点

（5）${\lambda }_1<0<{\lambda }_2$时，原点是不稳定鞍点

（6）${\lambda }_{1,2}=\alpha \pm \beta i,\alpha <0$时，原点是稳定焦点

（7）${\lambda }_{1,2}=\alpha \pm \beta i,\alpha <0$时，原点是稳定焦点

（6）${\lambda }_{1,2}=\alpha \pm \beta i,\alpha =0$时，原点是不稳定焦点

（8）${\lambda }_{1,2}=\alpha \pm \beta i,\alpha >0$时，原点是不稳定中心点

其余系统的稳定性判别方法

稳定性的判断方法：

如果常系数齐次线性系统的系数矩阵A的所有特征根都有负实部，则其渐进稳定，如果有正实部则不稳定。如果实部非正，但是有零实部，则可能稳定，但不可能渐进稳定。

非线性系统通过在奇点附近近似线性化，根据这个近似线性系统的奇点稳定性来判断非线性系统的奇点稳定性。

对系统$\frac{dx}{dt}=F(x)$，如果F在$x_0$处的雅各比矩阵是非奇异的，$detJ(x_0)\ne 0$，且在$x=0$的某邻域连续，有直到2阶的连续偏导数，则根据泰勒展开，系统的线性近似可以表示为
$$\frac{dx}{dt}=Ax$$
其中
$$A=\left[\begin{array}{ccc}\frac{\partial f_1}{\partial x_1}& \cdots & \frac{\partial f_1}{\partial x_n}\\ ⋮& ⋮& ⋮\\ \frac{\partial f_n}{\partial x_1}& \cdots & \frac{\partial f_n}{\partial x_n}\end{array}\right]$$

如果线性近似系统的零解是渐进稳定的或者不稳定的，则原系统也是渐进稳定的或者不稳定的。

如果线性近似系统的零解稳定，则无法推知原系统的稳定性

设非线性系统
$$\{\begin{array}{l}\frac{dx}{dt}=ax+by+\phi (x,y)\\ \\ \frac{dy}{dt}=ax+by+\psi (x,y)\end{array}$$
中的$\phi (x,y)$和$\psi (x,y)$满足条件：

（1）在点O的某邻域内存在连续的一阶偏导数

（2）存在常数$\delta >0$，使得
$$\underset{r\to 0}{\lim }\frac{\phi (x,y)}{r^{1+\delta }}=\underset{r\to 0}{\lim }\frac{\psi (x,y)}{r^{1+\delta }}=0r=\sqrt{x^2+y^2}$$
（3）该非线性系统的一次近似系统的特征方程的根没有零实部，则其与线性近似系统的奇点类型、稳定性相同

种群增长模型

一种动物种群或者其他具有生长性质的事物可能满足以下假设。

（1）由于总量足够大，因此可以将离散的数量变化看做连续可微的变化

（2）增长率之和环境性质有关，和时间无关，是自治的

（3）由于环境总资源有限，因此种群含有的个体数量会对增长率有一个成正比的抑制作用

满足以上假设，得到Logistic模型，微分方程如下
$$\{\begin{array}{l}\dot{x}(t)=rx(1-\frac{x}{N})\\ \\ x(0)=N_0\end{array}$$
x(t)表示t时刻的个体数量，r为固有增长率，N是环境容许的最大数量

和上一篇一样，我们可以通过python求解这个微分方程

import sympy as sy
N, r, N0, t, C1 = sy.symbols('N, r, N0, t, C1') 
x = sy.Function('x') 
ode = sy.diff(x(t), t, 1)-r*x(t)*(1-x(t)/N)
ode_sol = sy.simplify(sy.dsolve(ode, x(t)))  
sy.pprint(ode_sol)
ics = {x(0): N0} 
C_sol = sy.solve(ode_sol.subs(t, 0).subs(ics),C1)
sy.pprint(sy.simplify(ode_sol.subs({C1:C_sol[0]}))) 

输出如下

           C₁⋅N + r⋅t 
        N⋅ℯ
x(t) = ───────────────
        C₁⋅N + r⋅t
       ℯ           - 1
                r⋅t    
          N⋅N₀⋅ℯ
x(t) = ────────────────
               r⋅t
       N + N₀⋅ℯ    - N₀

可知其通解为
$$x(t)=\frac{N\cdot e^{C_1\cdot N+r\cdot t}}{e^{C_1\cdot N+r\cdot t}-1}$$
可知其特解为
$$x(t)=\frac{N\cdot N_0\cdot e^{rt}}{N_0\cdot e^{rt}+N-N_0}=\frac{N}{1+\frac{N-N_0}{N_0}{e}^{-rt}}$$
求出这个微分方程的特解之后，可以通过$t\to \infty$的时候$x(t)=N$来判断最后种群中的个体数会稳定在什么值

不过非线性系统可能有多个稳定位置，不同的初始条件对应的稳定值是不同的，因此如果通过这个方法，在初始条件不确定的情况下很难确定系统最后的情况。

通过稳定性理论可以在不求解微分方程的情况下获知，步骤如下

根据微分方程可以得知，$\dot{x}=0$时，$x=0\text{或}N$，这两点为系统的平衡点

对该系统
$$f(x)=rx(1-\frac{x}{N})=-\frac{r}{N}\cdot x^2+rx$$

$$f^{\prime }(x)=-\frac{2r}{N}\cdot x+r$$

在平衡点附近得到近似线性方程
$$\dot{x}(t)=f^{\prime }(x_e)(x-x_e)$$
对于平衡点$x_{e1}=0$，$f^{\prime }(x_e)=r>0$，近似线性系统和原系统均不稳定

对于平衡点$x_{e2}=N$，$f^{\prime }(x_e)=-r<0$，近似线性系统和原系统均稳定

因此系统最后会稳定至$x=N$

种群开发模型

假如我们要对一个种群进行采集、捕捉和开发，开发过程不能竭泽而渔，因此，开发的本质是在保证种群个体量稳定的情况下获取最大的经济效益，也是一个稳定状态数学模型。

我们认为单位时间能够完成的开发量，比如捕猎量，采摘量等，和环境中的种群个体量$x(t)$成正比。为了能够方便我们进行控制和规划，因此我们令开发系数为开发强度E和单位开发强度的开发率q的乘积，于是，在种群个体量为x的情况下，单位时间的开发量为
$$h(x)=Eqx(t)$$
结合上述种群增长模型，得到Scheafer模型
$$\dot{x}(t)=rx(1-\frac{x}{N})-Eqx$$
根据微分方程可以得知，$\dot{x}=0$时，
$$x(\frac{rx}{N}-r+Eq)=0$$
得$x=0\text{或}\frac{(r-Eq)N}{r}$，这两点为系统的平衡点

对该系统
$$f(x)=rx(1-\frac{x}{N})-Eqx=-\frac{r}{N}\cdot x^2+(r-Eq)x$$

$$f^{\prime }(x)=-\frac{2r}{N}\cdot x+r-Eq$$

在平衡点附近得到近似线性方程
$$\dot{x}(t)=f^{\prime }(x_e)(x-x_e)$$
对于平衡点$x_{e1}=0$，$f^{\prime }(x_e)=r-Eq$

对于平衡点$x_{e2}=\frac{(r-Eq)N}{r}$，$f^{\prime }(x_e)=-r+Eq$

因此

当固有增长率r大于开发系数Eq时，$r-Eq>0$，此时平衡点$x_{e1}=0$不稳定，平衡点$x_{e2}=\frac{(r-Eq)N}{r}$稳定，系统稳定状态时，$x=\frac{(r-Eq)N}{r}$

当固有增长率r小于开发系数Eq时，$r-Eq<0$，此时平衡点$x_{e1}=0$稳定，平衡点$x_{e2}=\frac{(r-Eq)N}{r}$不稳定，系统稳定状态时，$x=0$

因此，为了能够持续开采，我们首先要保证Eq<r，在这种情况下稳定的单位时间开发量为
$$h(\frac{(r-Eq)N}{r})=Eq\frac{(r-Eq)N}{r}$$
要求h达到最大值，即求$Eq\frac{(r-Eq)N}{r}=qN(-\frac{q}{r}{E}^2+E)$在$E$取何值时，最大

通过简单的二次函数求极值，可以求得$E=\frac{r}{2q}$时，$h_{max}=\frac{rN}{4}$

依照这个开发强度，可以获得最大的开发量

不过种群捕捉开发的目的往往是获得最大的经济效益而不是获得最大的开发量，由于开发强度同时也意味着开发成本，因此对于经济效益的规划问题就要建立收支模型。令单位开发强度会对应成本c，单位种群个体的收益为p，

则单位时间的收入$T$为
$$T=ph(x)=pEqx$$
单位时间的支出$S$为
$$S=cE$$
单位时间的利润$R$为
$$R=T-S=pEqx-cE$$
在稳定状态时，
$$R=pEq\frac{(r-Eq)N}{r}-cE$$
同样可以通过二次函数求极值的方法求解，在此不具体论述了

种群相互竞争的模型

当两种种群分别在有限环境中增长时，各自满足Logistic模型
$$\dot{x}(t)=rx(1-\frac{x}{N})$$
但是当有限环境中有两个类似的种群各自增长时，因为其对环境的消耗也会抑制另外的种群增长，所以需要在其增长模型里减去由另一种种群增长对环境的消耗带来的抑制作用

引入指标${\sigma }_1,{\sigma }_2$，表示抑制作用，假设这两个指标相互独立，增长模型变为
$$\{\begin{array}{l}{\dot{x}}_1=r_1{x}_1(1-\frac{x_1}{N_1}-{\sigma }_1\frac{x_2}{N_2})\\ \\ {\dot{x}}_2=r_2{x}_2(1-{\sigma }_2\frac{x_1}{N_1}-\frac{x_2}{N_2})\end{array}$$
分析其增长结果，求平衡点
$$\{\begin{array}{l}{\dot{x}}_1=r_1{x}_1(1-\frac{x_1}{N_1}-{\sigma }_1\frac{x_2}{N_2})=0\\ \\ {\dot{x}}_2=r_2{x}_2(1-{\sigma }_2\frac{x_1}{N_1}-\frac{x_2}{N_2})=0\end{array}$$
得到平衡点有$x_{e1}(0,0),x_{e2}(0,N_2),x_{e3}(N_1,0),x_{e4}(\frac{N_1(1-{\sigma }_1)}{1-{\sigma }_1{\sigma }_2},\frac{N_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2})$

利用上文提到的一阶线性化方法，得到系数矩阵A
$$A=\left[\begin{array}{cc}\frac{\partial f_1}{\partial x_1}& \frac{\partial f_1}{\partial x_2}\\ & \\ \frac{\partial f_2}{\partial x_1}& \frac{\partial f_2}{\partial x_2}\end{array}\right]=\left[\begin{array}{cc}{r}_1-\frac{2r_1{x}_1}{N_1}-\frac{r_1{\sigma }_1{x}_2}{N_2}& -\frac{r_1{\sigma }_1{x}_1}{N_2}\\ & \\ -\frac{r_2{\sigma }_2{x}_2}{N_1}& r_2-\frac{2r_2{x}_2}{N_2}-\frac{r_2{\sigma }_2{x}_1}{N_1}\end{array}\right]$$
对$x_{e1}(0,0)$得
$$A=\left[\begin{array}{cc}{r}_1& 0\\ 0& r_2\end{array}\right]$$
根据上文二阶系统稳定性分析的方法，
$$p=-(r_1+r_2)$$

$$q=detA=r_1{r}_2$$

$${\lambda }_{1,2}=\frac{1}{2}(r_1+r_2\pm |r_1-r_2|)$$

因为$p<0$，所以该奇点不稳定

对$x_{e2}(0,N_2)$得
$$A=\left[\begin{array}{cc}{r}_1-r_1{\sigma }_1& 0\\ -r_2{\sigma }_2& -r_2\end{array}\right]$$
根据上文二阶系统稳定性分析的方法，
$$p=-(r_1-r_1{\sigma }_1-r_2)=-r_1(1-{\sigma }_1)+r_2$$

$$q=detA=-r_1{r}_2+r_1{r}_2{\sigma }_1=-r_1{r}_2(1-{\sigma }_1)$$

$${\lambda }_{1,2}=\frac{1}{2}(r_1(1-{\sigma }_1)-r_2\pm |r_1(1-{\sigma }_1)+r_2|)$$

假如系统稳定，即$p>0\text{且}q>0$，
$$\{\begin{array}{l}-r_1(1-{\sigma }_1)+r_2>0\\ \\ -r_1{r}_2(1-{\sigma }_1)>0\end{array}$$
得${\sigma }_1>1$

对$x_{e3}(N_1,0)$得
$$A=\left[\begin{array}{cc}-r_1& -r_1{\sigma }_1\\ & \\ 0& (1-{\sigma }_2)r_2\end{array}\right]$$
根据上文二阶系统稳定性分析的方法，
$$p=-[-r_1+(1-{\sigma }_2)r_2]=r_1-(1-{\sigma }_2)r_2$$

$$q=detA=-r_1{r}_2(1-{\sigma }_2)$$

假如系统稳定，即$p>0\text{且}q>0$，
$$\{\begin{array}{l}{r}_1-(1-{\sigma }_2)r_2>0\\ \\ -r_1{r}_2(1-{\sigma }_2)>0\end{array}$$
得${\sigma }_2>1$

对$x_{e4}(\frac{N_1(1-{\sigma }_1)}{1-{\sigma }_1{\sigma }_2},\frac{N_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2})$得
$$A=\left[\begin{array}{cc}{r}_1(-\frac{1-{\sigma }_1}{1-{\sigma }_1{\sigma }_2})& -\frac{r_1{\sigma }_1}{N_2}\frac{N_1(1-{\sigma }_1)}{1-{\sigma }_1{\sigma }_2}\\ & \\ -\frac{r_2{\sigma }_2}{N_1}\frac{N_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2}& r_2(-\frac{1-{\sigma }_2}{1-{\sigma }_1{\sigma }_2})\end{array}\right]$$
根据上文二阶系统稳定性分析的方法，
$$\begin{array}{l}p=r_1(\frac{1-{\sigma }_1}{1-{\sigma }_1{\sigma }_2})+r_2(\frac{1-{\sigma }_2}{1-{\sigma }_1{\sigma }_2})\\ \\ =\frac{(r_1+r_2)-r_1{\sigma }_1-r_2{\sigma }_2}{1-{\sigma }_1{\sigma }_2}\\ \\ =\frac{r_1(1-{\sigma }_1)+r_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2}\end{array}$$

$$\begin{array}{l}q=detA\\ =r_1(\frac{1-{\sigma }_1}{1-{\sigma }_1{\sigma }_2})r_2(\frac{1-{\sigma }_2}{1-{\sigma }_1{\sigma }_2})-\frac{r_1{\sigma }_1}{N_2}\frac{N_1(1-{\sigma }_1)}{1-{\sigma }_1{\sigma }_2}\frac{r_2{\sigma }_2}{N_1}\frac{N_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2}\\ \\ =\frac{r_1{r}_2(1-{\sigma }_1)(1-{\sigma }_2)}{(1-{\sigma }_1{\sigma }_2{)}^2}-\frac{r_1{r}_2{\sigma }_1{\sigma }_2(1-{\sigma }_1)(1-{\sigma }_2)}{(1-{\sigma }_1{\sigma }_2{)}^2}\\ \\ =\frac{r_1{r}_2(1-{\sigma }_1)(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2}\end{array}$$

假如系统稳定，即$p>0\text{且}q>0$，
$$\{\begin{array}{l}\frac{r_1(1-{\sigma }_1)+r_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2}>0\\ \\ \frac{r_1{r}_2(1-{\sigma }_1)(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2}>0\end{array}$$

得${\sigma }_1<1,{\sigma }_2<1$

由此可以分情况讨论

（1）${\sigma }_1<1,{\sigma }_2<1$此时，$x_{e1},x_{e2},x_{e3}$不稳定， $x_{e4}$稳定，因此最后两种种群共存，保持在$x_{e4}(\frac{N_1(1-{\sigma }_1)}{1-{\sigma }_1{\sigma }_2},\frac{N_2(1-{\sigma }_2)}{1-{\sigma }_1{\sigma }_2})$

（2）${\sigma }_1>1,{\sigma }_2<1$此时，$x_{e1},x_{e3},x_{e4}$不稳定， $x_{e2}$稳定，此时，第二种种群增长对第一种种群增长的抑制作用过大，导致第一种种群最后灭绝，只剩第二种种群，稳定在$x_{e2}(0,N_2)$

（3）${\sigma }_1<1,{\sigma }_2>1$此时，$x_{e1},x_{e2},x_{e4}$不稳定， $x_{e3}$稳定，此时，第一种种群增长对第二种种群增长的抑制作用过大，导致第二种种群最后灭绝，只剩第一种种群，稳定在$x_{e3}(N_1,0)$

（4）${\sigma }_1>1,{\sigma }_2>1$此时，$x_{e1},x_{e4}$不稳定， $x_{e2},x_{e3}$稳定，此时，稳定在这两个平衡点都有可能，根据初始时刻两种种群的量而定。

不过在一个只有两种种群的情况下，往往${\sigma }_1<1$代表第二种种群增长对第一种种群增长的抑制能力弱，由于两种种群在争夺同样的环境要素，所以，${\sigma }_1<1$表示第一种种群抢夺环境的能力比第二种强，一般情况下${\sigma }_2\approx \frac{1}{{\sigma }_1}>1$。所以稳定在$x_{e4}$以及同时可能稳定在两个平衡点的情况并不会出现。

种群捕食与被捕食的模型

为了简化模型，我们暂且认为，在会被捕食的情况下，种群个体不会泛滥，种群增长带来的对增长率的抑制作用可以暂且忽略不计，同时环境中只有两种种群，且其只存在捕食与被捕食的关系。

由于会被捕食，所以捕食者的数量会对被捕食种群的增长产生抑制作用，我们定义被捕食者个体数为$x_1$，抑制系数为${\lambda }_1$，没有捕食者时的自然增长率为$r_1$，则
$${\dot{x}}_1(t)=x_1(r_1-{\lambda }_1{x}_2)$$
由于捕食者并不能从环境中直接汲取资源使种群增长，因此当没有被捕食者存在时，捕食者数量会自然下降。而被捕食者的数量会捕食者种群的增长产生促进作用，定义捕食者个体数为$x_2$，促进系数为${\lambda }_2$，在没有被捕食者时的自然衰减率为$r_2$，则
$${\dot{x}}_2(t)=x_2(-r_2+{\lambda }_2{x}_1)$$
于是该系统的微分方程组为
$$\{\begin{array}{l}{\dot{x}}_1(t)=x_1(r_1-{\lambda }_1{x}_2)\\ \\ {\dot{x}}_2(t)=x_2(-r_2+{\lambda }_2{x}_1)\end{array}$$
分析稳定状态时的种群情况，求其平衡点
$$\{\begin{array}{l}{\dot{x}}_1(t)=x_1(r_1-{\lambda }_1{x}_2)=0\\ \\ {\dot{x}}_2(t)=x_2(-r_2+{\lambda }_2{x}_1)=0\end{array}$$
得到平衡点$x_{e1}(0,0)$和$x_{e2}(\frac{r_2}{{\lambda }_2},\frac{r_1}{{\lambda }_1})$

利用上文提到的一阶线性化方法，得到系数矩阵A
$$A=\left[\begin{array}{cc}\frac{\partial f_1}{\partial x_1}& \frac{\partial f_1}{\partial x_2}\\ & \\ \frac{\partial f_2}{\partial x_1}& \frac{\partial f_2}{\partial x_2}\end{array}\right]=\left[\begin{array}{cc}{r}_1-{\lambda }_1{x}_2& -x_1{\lambda }_1\\ & \\ x_2{\lambda }_2& -r_2+{\lambda }_2{x}_1\end{array}\right]$$
对$x_{e1}(0,0)$得
$$A=\left[\begin{array}{cc}{r}_1& 0\\ 0& -r_2\end{array}\right]$$
根据上文二阶系统稳定性分析的方法，
$$p=-(r_1-r_2)$$

$$q=detA=-r_1{r}_2$$

因为$q<0$，所以该奇点不稳定

对$x_{e2}(\frac{r_2}{{\lambda }_2},\frac{r_1}{{\lambda }_1})$得
$$A=\left[\begin{array}{cc}0& -\frac{r_2{\lambda }_1}{{\lambda }_2}\\ & \\ \frac{r_1{\lambda }_2}{{\lambda }_1}& 0\end{array}\right]$$
根据上文二阶系统稳定性分析的方法，
$$p=0$$

$$q=detA=r_1{r}_2$$

特征根为
$$\pm \frac{1}{2}\sqrt{-4r_1{r}_2}$$
因此该奇点为不稳定中心点

由此我们可以得出结论，捕食与被捕食的关系会导致两个种群数量形成一个生态循环，此消彼长地往复，而不会趋于一个固定的值。