摘" 要： 主要研究一類病原體與中性粒細胞（PMNs）相互作用的免疫應(yīng)答模型，詳細分析了其正平衡點的存在條件和類型，并通過MATLAB軟件包對所得到的結(jié)論進行了數(shù)值模擬。

關(guān)鍵詞： 免疫應(yīng)答模型；正平衡點；尖點

中圖分類號： O175

文獻標(biāo)識碼： A" 文章編號： 2096-3998（2024）02-0079-07

收稿日期：2023-04-25" 修回日期：2023-06-27

基金項目：國家自然科學(xué)基金項目（61763024）

作者簡介：孫得寧（1998—），男，甘肅張掖人，碩士研究生，主要研究方向為微分方程與動力系統(tǒng)；張存華（1972—），女，甘肅武山人，博士，教授，主要研究方向為微分方程與動力系統(tǒng)。

引用格式：孫得寧，張存華.一類免疫應(yīng)答模型的分支分析.陜西理工大學(xué)學(xué)報（自然科學(xué)版），2024，40（2）：79-85.

人體內(nèi)主要通過先天免疫系統(tǒng)對病原體作出反應(yīng)。病原體入侵肺部后，首先巨噬細胞殺死病原體，當(dāng)巨噬細胞檢測到自身無法遏制的威脅時，中性粒細胞（PMNs）會對病原體作出反應(yīng)。中性粒細胞與淋巴細胞的比率（NLR）已被視為預(yù)測肺炎感染的可靠因子，Wang Dawei等發(fā)現(xiàn)新冠患者在危重期PMNs計數(shù)增加且淋巴細胞計數(shù)減少（即NLR增加），這表明患者病情可能處于危重狀態(tài)。

Pugliese等提出了一個關(guān)于病原體與免疫細胞相互作用的系統(tǒng)，但相互作用項非常具體；D’onofrio對該系統(tǒng)進行了優(yōu)化。從文獻來看，PMNs的活性在正常水平下較低，隨著病原體入侵肺部后迅速增加，因此，Young等構(gòu)建了如下系統(tǒng)：

=αx-mx1+bx-βxz，

=δx-γxz+η-μz，（1）

其中，x（t）表示t時刻病原體數(shù)量的密度，z（t）表示t時刻PMNs肺部中載量的密度，α表示病原體的復(fù)制速率，mx1+bx表示通過先天免疫系統(tǒng)清除病原體的數(shù)量，-βxz表示PMNs吞噬病原體的數(shù)量，-γxz表示單個獨立的免疫細胞抑制病原體的數(shù)目有限，η表示正常水平下PMNs的流入量，-μz表示PMNs的自然死亡率，參數(shù)α、m、b、β、δ、γ、η、 μ均為正。

Young等發(fā)現(xiàn)當(dāng)少量病原體入侵人體后，能夠迅速被PMNs清除，即

α＜m，（2）

而他們認為巨噬細胞清除病原體的能力有限，則分析了系統(tǒng)（1）在條件αβ＞δγ下平衡點的類型，且得到邊界平衡點是穩(wěn)定結(jié)點，正平衡點是鞍點。然而，參數(shù)的實際值因人而異，較大病原體可能不會克服先天免疫系統(tǒng)，所以又考慮了相反的情況

αβ＜δγ，（3）

并分析了系統(tǒng)（1）在條件（3）下平衡點的個數(shù)，但沒有提供具體的理論分析。Shi Shujing等分析了系統(tǒng)（1）在條件（3）下平衡點的類型及其穩(wěn)定性，并得到隨著參數(shù)值的變化，系統(tǒng)（1）展現(xiàn)了一系列分支現(xiàn)象：包括余維一鞍結(jié)點分支和Hopf分支。然而，隨著年齡的增長人體的免疫系統(tǒng)開始衰老和減弱，PMNs的自然死亡率隨之升高，故提出如下系統(tǒng)：

=αx-mx1+bx-βxz，

=δx-γxz+η-μz2。（4）

接下來主要討論系統(tǒng)（4）的正平衡點及其類型。

1" 正平衡點及其類型

首先做如下尺度變換：

x=α2βδx1，" z=αβz1，" t=1αt1，

用x、z、t表示x1、z1、t1，則系統(tǒng)（4）化為

=x-p1x1+p2x-xz，

=x-p3xz+p4-p5z2，（5）

其中，

p1=mα，

p2=α2bβδ，

p3=αγβδ，

p4=βηα2，

p5=μβ。

由式（2）和式（3）知p1gt;1，p3lt;1，而p1、p2、p3、p4、p5均為正參數(shù)，因此，僅在下列條件下研究系統(tǒng)（5）：

0lt;p3lt;1

p1gt;1

p2，p4，p5gt;0。（6）

系統(tǒng)的正不變區(qū)域為

Ω={（x，z）|x≥0，z≥0}，

為了討論系統(tǒng)（5）的正平衡點，令

x-p1x1+p2x-xz=0，

x-p3xz+p4-p5z2=0，（7）

由此得

z=1-p11+p2x，

（p22-p22p3）x3+（2p2-2p2p3+p1p2p3+p22p4-p22p5）x2+

（1-p3+p1p3+2p2p4-2p2p5+2p1p2p5）x+p4-p5+2p1p5-p21p5=0，

令

f（x）=（p22-p22p3）x3+（2p2-2p2p3+p1p2p3+p22p4-p22p5）x2+

（1-p3+p1p3+2p2p4-2p2p5+2p1p2p5）x+p4-p5+2p1p5-p21p5，

則

f′（x）=3（p22-p22p3）x2+（4p2-4p2p3+2p1p2p3+2p22p4-2p22p5）x+

1-p3+p1p3+2p2p4-2p2p5+2p1p2p5，

令

a=3p22-3p22p3，

b=4p2-4p2p3+2p1p2p3+2p22p4-2p22p5，

c=1-p3+p1p3+2p2p4-2p2p5+2p1p2p5，

Δ=b2-4ac，

則當(dāng)Δgt;0且blt;0時， f′（x）有兩個正根，即f（x）存在極值點μ1、 μ2；當(dāng)Δlt;0時， f′（x）不存在正根，即f（x）是單調(diào)遞增的。

由式（6）得

p22-p22p3gt;0。

那么得到f（x）根的5種情況（見圖1）。

（a） 三個不同的正根：E1，E2，E3""" （b） 兩個正根：重根E*，單根E3""" （c） 兩個正根：單根E1，重根E*

（d） 一個正根：重根E*""""""""""""" （e） 不存在正根

由圖1可知，首先，當(dāng)f（0）lt;0，blt;0，Δgt;0， f（μ1）gt;0且f（μ2）lt;0（μ2gt;μ1gt;0）時， f（x）有三個正根（見圖1（a））。

其次，當(dāng)f（0）lt;0，blt;0，Δgt;0， f（μ1）=0或f（μ2）=0時， f（x）有一個正的單根和一個正的重根（見圖1（b）（c））。

再次，當(dāng)f（0）gt;0，blt;0，Δgt;0， f（μ1）=0時， f（x）有一個正的重根（見圖1（d））。

最后，當(dāng)f（0）gt;0，blt;0，Δgt;0， f（μ2）gt;0或f（0）gt;0，Δlt;0時， f（x）不存在正根（見圖1（e））。

接下來研究系統(tǒng)（5）正平衡點的類型，系統(tǒng)（5）在E00，p4p5處的Jacobian矩陣為

J（E0）=

1-p1-p4p50

1-p3p4p5-2p4p5

〗，

該矩陣有兩個特征值λ1=1-p1-p4p5和λ2=-2p4p5，再由（6）知p1gt;1，則λ1lt;0，因此E0是穩(wěn)定雙曲結(jié)點。

系統(tǒng)（5）在正平衡點E（x，z）處的Jacobian矩陣為

J（E）=

p1p2x（1+p2x）2-x

1-p3+p1p31+p2x-p3x-2p5+2p1p51+p2x

〗，

由Jacobian矩陣得到J（E）的行列式為

Det（J（E））=

2p21p2p5x-（p1p2p3x2+2p1p2p5x）（1+p2x）

（1+p2x）3

，

J（E）的跡為

Tr（J（E））=

p1p2x（1+p2x）2+

2p1p51+p2x-

p3x-2p5。

根據(jù)Det（J（E））與f′（x）之間的關(guān)系得到

Det（J（E））=

xf′（x）（1+p2x）2。（8）

令

p*1=（1+p2x）p2x

則當(dāng)p1lt;p*1（p1gt;p*1）時，Tr（J（E））gt;0（Tr（J（E））lt;0），由以上結(jié)論得到以下定理。

定理1" 當(dāng)（6）中的條件滿足時，系統(tǒng)（5）總是存在一個無病平衡點E00，p4p5，它是穩(wěn)定雙曲結(jié)點。此外，

（Ⅰ）若f（0）gt;0、blt;0、Δgt;0、 f（μ2）gt;0，或f（0）gt;0、Δlt;0，則系統(tǒng)（5）不存在正平衡點；

（Ⅱ）若f（0）gt;0，blt;0，Δgt;0， f（μ2）=0，則系統(tǒng)（5）存在唯一正平衡點E*（x*，z*），它是一個退化平衡點；

（Ⅲ）若f（0）lt;0，blt;0，Δgt;0， f（μ1）=0或f（μ2）=0，則系統(tǒng)（5）存在兩個不等的正平衡點：一個退化平衡點E*（x*，z*）和一個雙曲平衡點Ei（i=1，3）。當(dāng)p1lt;p*1時，Ei（i=1，3）為穩(wěn)定的雙曲結(jié)點或焦點；當(dāng)p1gt;p*1時，Ei（i=1，3）為不穩(wěn)定的雙曲結(jié)點或焦點；當(dāng)p1=p*1時，Ei（i=1，3）為弱焦點或中心；

（Ⅳ）若f（0）lt;0，blt;0，Δgt;0， f（μ1）gt;0且f（μ2）lt;0，則系統(tǒng)（5）存在三個不等的正平衡點：E2為雙曲鞍點，Ei（i=1，3）為雙曲結(jié)點或焦點；當(dāng)p1lt;p*1時，Ei為穩(wěn)定的雙曲結(jié)點或焦點；當(dāng)p1gt;p*1時，Ei為不穩(wěn)定的雙曲結(jié)點或焦點；當(dāng)p1=p*1時，Ei為弱焦點或中心。

證明" 根據(jù)等式（8）與圖1可知Det（J（Ei））gt;0（i=1，3），Det（J（E2））lt;0，Det（J（E*））=0，則E1、E3為雙曲結(jié)點或焦點，E2為雙曲鞍點，E*為退化平衡點。

定理2" 當(dāng)f（0）gt;0，blt;0，Δgt;0， f（μ2）=0與（6）中的條件滿足時，系統(tǒng)（5）存在唯一正平衡點E*（x*，z*）。此外，

（Ⅰ）當(dāng)p1lt;p*1（或p1gt;p*1）時，E*（x*，z*），是包含一個穩(wěn)定拋物扇形（或不穩(wěn)定拋物扇形）的鞍結(jié)點；

（Ⅱ）當(dāng)p1=p*1時，E*（x*，z*）是一個尖點，此外，

（?。┤鬎≠0，則E*（x*，z*）是一個余維二的尖點；

（ⅱ）若F=0，則E*（x*，z*）是一個至少余維三的尖點。

證明" 令Tr（J（E*））=0得到

p1=p*1=（1+p2x）p2x，

那么當(dāng)p1lt;p*1（p1gt;p*1）時，Tr（J（E*））lt;0（Tr（J（E*））gt;0）。

下面來討論E*（x*，z*）的具體類型，令X=x-x*，Z=z-z*，將平衡點移到原點并將該式進行泰勒展開，則系統(tǒng)（5）化為

dXdt=a11X+a12Z+a13X2-XZ，

dZdt=a21X+a22Z-p3XZ，

（9）

其中，

a11=p1p2x*（1+p2x*）2，

a12=-x*，

a13=p1p2（1+p2x*）4，

a21=p1p31+p2x*-p3+1，

a22=2p1p51+p2x*-p3x*-2p5。

令

X=-x*μ+p1p2x*（1+p2x*）2v，

Y=-p1p2x*（1+p2x*）2μ+1-p3+p1p31+p2x*v，

τ=p1p2x*（1+p2x*）2+2p1p51+p2x*-p3x*-2p5t。

用t表示τ，則系統(tǒng)（9）化為

=b11μ2+b12μv+b13v2，

v·=v+b21μ2+b22μv+b23v2，（10）

其中，

b11=a12a13a21a22+a11a21a22-p3a211a22（a11a21+a21a22）（a11+a22），

b12=（2a11a12a13-a12a21+a211+p3a11a12-p3a211a12）（a11+a22）a12a22（a11+a22），

b13=a211a13a22-a11a21a22+p3a211a22（a12a22+a11a12）（a11+a22），

b21=p3a11a12a11a21+a21a22+a11a12a13a21a22+a211a21a22-p3a311a22（a11a221+a221a22）（a11+a22），

b22=p3a211-p3a12a21a21（a11+a22）+2a211a12a13-a11a12a21+a311+p3a211a12-p3a311a12a12a21a22，

b23=a311a13a22-a211a21a22+p3a311a22（a12a21a22+a11a12a21）（a11+a22）-p3a11a11+a22，

經(jīng)計算判斷b11≠0，依據(jù)文獻中11.3節(jié)的中心流形理論，得到限制在中心流形上的方程為

dμdt=b11μ2+o（|μ|3），（11）

依據(jù)文獻中的定理7.1，結(jié)論（Ⅰ）成立。

令X=x-x*，Z=z-z*，p1=p*1，Det（J（E*））=0，則系統(tǒng)（5）化為

dXdt=c11X+c12Z+c13X2-XZ，

dZdt=c21X+c22Z-p3XZ，（12）

其中，

c11=p1p2x*（1+p2x*）2，

c12=-x*，

c13=p1p2（1+p2x*）4，

c21=p1p31+p2x*-p3+1，

c22=2p1p51+p2x*-p3x*-2p5。

令

X=-x*μ-v，

Z=-p1p2x*（1+p2x*）2μ，

τ=（1+p2x*）2p1p2x*（1+p2x*）t。

用t表示τ，則系統(tǒng)（12）化為

=v+d11μ2+d12μv，

v·=d21μ2+d22μv+d23v2，（13）

其中，

d11=p3c11c12c21，

d12=p3c11c21，

d21=p3c11c212-c11c212c13-c211c12c21，

d22=c211-p3c11+2c11c12c13c21，

d23=-c11c13c21。

根據(jù)文獻中的引理2.4，用x、z表示μ、v，得到系統(tǒng)（13）在原點的小鄰域內(nèi)等價于系統(tǒng)（14）：

=z，

=Dx2+Fxz+o（|x，z|3），（14）

其中，D=d21，F(xiàn)=d22+2d11。經(jīng)計算得D=d21≠0，則當(dāng)F≠0時，E*（x*，z*）是一個余維二的尖點；而當(dāng)F=0時，E*（x*，z*）是一個至少余維三的尖點。

2" 數(shù)值模擬

通過MATLAB軟件包對所獲得的結(jié)論進行數(shù)值驗證。取p1=3 462/2 500+0.1，p2=23，p3=7/10，p4=194/12 000，p5=13/25，滿足條件p1lt;p*1，得到E*（x*，z*）是包含一個穩(wěn)定拋物扇形的鞍結(jié)點（見圖2（a））；取p1=57 598/237 862+0.1，p2=58/91，p3=78/62，p4=86/91-0.046，p5=65/93，滿足條件p1gt;p*1，得到E*（x*，z*）是包含一個不穩(wěn)定拋物扇形的鞍結(jié)點（見圖2（b））。

（a） p1lt;p*1時，E*（x*，z*）是包含一個穩(wěn)定拋物扇形的鞍結(jié)點

（b） p1gt;p*1時，E*（x*，z*）是包含一個不穩(wěn)定拋物扇形的鞍結(jié)點

3" 結(jié)語

本文主要研究了一類描述病原體與宿主免疫反應(yīng)之間相互作用的免疫應(yīng)答模型，并討論了它的正平衡點及其類型，得到該系統(tǒng)至多存在三個正平衡點，當(dāng)其中兩個平衡點合為二重平衡點時，該二重平衡點可能為余維二尖點或余維三尖點。最后通過MATLAB軟件包對所獲得的結(jié)論進行了數(shù)值驗證。

［" 參" 考" 文" 獻" ］

WANG Dawei，HU Bo，HU Chang，et al.Clinical characteristics of 138 hospitalized patients with 2019 novel coronavirus-infected pneumonia in Wuhan，China.Jama，2020，323（11）：1061-1069.

PUGLIESE A，GANDOLFI A.A simple model of pathogen-immune dynamics including specific and nonspecific immuni-ty.Mathematical Biosciences，2008，214（1/2）：73-80.

D’ONOFRIO A.On the interaction between the immune system and an exponentially replicating pathogen.Mathematical Biosciences And Engineering，2010（7）：579-602.

YOUNG T R，BUCKALEW R，MAY A K，et al.A low dimensional dynamical model ofthe initial pulmonary innate reponse to infection.Mathematical Biosciences，2012，235（2）：189-200.

SHI Shujing，HUANG Jicai，WEN Jing.Bifurcation Analysis of a Dynamical Model for the Innate Immune Response to Initial Pulmonary Infections.International Journal of Bifurcation and Chaos，2020，30（16）：1-22.

張錦炎，馮貝葉.常微分方程幾何理論與分支問題.北京：北京大學(xué)出版社，1997：329-339.

張芷芬，丁同仁，黃文灶，等.微分方程定性理論.北京：科學(xué)出版社，1997：125-140.

閆夢月.腫瘤細胞與效應(yīng)細胞相互作用模型的動力學(xué)分析.武漢：華中師范大學(xué)，2022.

PERKO L.Differential Equations and Dynamical Systems.New York： Springer，2001：147-163.

［責(zé)任編輯：魏 強］

Bifurcation analysis of a class of immune response models

SUN Dening，" ZHANG Cunhua

School of Mathematics and Physics， Lanzhou Jiaotong University， Lanzhou 730070， China

Abstract：" We mainly study the immune response model of the interaction between a class of pathogens and neutrophils（PMNs）， and analyze the conditions and type of their positive equilibrium point in detail， and carry out numerical simulation of the conclusions obtained by MATLAB software package.

Key words：" immune response model; positive equilibrium point; sharp points