姬文娟 李斌 崔明 張燁
【摘要】我們考慮HIV1病毒動力學(xué)模型的特殊情況,一類三維競爭系統(tǒng),通過常微分方程中競爭系統(tǒng)和一致持續(xù)系統(tǒng)的性質(zhì),證明了高維自治系統(tǒng)周期解的存在性和軌道漸近穩(wěn)定性.若R0<1,體內(nèi)病毒可以被清除,疾病最終消失;若R0>1,病毒在宿主體內(nèi)持續(xù)存在,模型解會趨于地方病平衡點或者一個周期解,后者只在一定條件下排除存在的可能性.同時,我們選擇合適的參數(shù)和閾值,針對幾類模型的理論分析進行數(shù)值模擬,得到了較為合理的模擬結(jié)果并且直觀反映了理論分析的合理性.
【關(guān)鍵詞】HIV;細(xì)胞間傳播;穩(wěn)定性;周期解;競爭系統(tǒng)
本文我們主要證明一類三維競爭系統(tǒng)周期解的存在性和穩(wěn)定性,同時在模型中不考慮宿主體內(nèi)的細(xì)胞直接接觸感染k2項,是為了利用三維競爭系統(tǒng)的性質(zhì)來得到結(jié)論.從而得到下面這個模型:
dTdt=s+pT(1-TTmax)-dTT-(1-ηRT)k1VTdT*dt=(1-ηRT)k1VT-δT*dVdt=NδT*-cV(1)
s表示宿主體內(nèi)產(chǎn)生健康T細(xì)胞的比例,當(dāng)然已有的T細(xì)胞也能夠產(chǎn)生健康T細(xì)胞,p表示最長生命周期率,Tmax表示T細(xì)胞在其生命周期內(nèi)產(chǎn)生的T細(xì)胞種群密度.假設(shè)健康細(xì)胞T的自然死亡率dT和被感染細(xì)胞T*的死亡率δ是恒定的,c表示自由病毒粒子V的死亡率.k1VT表示T細(xì)胞接觸病毒粒子V以后被感染的數(shù)量比,k1是恒定的感染率.N表示被感染細(xì)胞生命周期內(nèi)釋放的病毒粒子的總數(shù).記f(T)=s+pT(1-TTmax)-dTT.
1.周期解的存在性和軌道漸近穩(wěn)定性
定理5 假如系統(tǒng)(1)存在周期解,那么解軌道是漸近穩(wěn)定的.
證明.設(shè)這個周期解是p(t)=(p1(t),p2(t),p3(t))T,并且假定其最小正周期是ω>0.根據(jù)解的有界性,0≤p1(t)≤T0,對t∈[0,ω].
證明周期解的穩(wěn)定性需要利用第二次復(fù)合方程的方法,構(gòu)造下面這個系統(tǒng):
z·=f[2]x(p(t))z,z=(z1(t),z2(t),z3(t))T.(4)
f[2]x=a11+a22a23-a13a32a11+a33a12-a31a21a22+a33
=f′(T)-(1-ηRT)k1V-δ(1-ηRT)k1T(1-ηRT)k1TNδf′(T)-(1-ηRT)k1V-c00(1-ηRT)k1V-δ-c
aij表示系統(tǒng)Jacobian矩陣中第i行第j列的值.(4)是原系統(tǒng)的線性變分方程,也被稱為第二次復(fù)合方程,f[2]x稱為f的Jacobian矩陣fx的第二次復(fù)合矩陣.假如系統(tǒng)(4)漸近穩(wěn)定,那么原系統(tǒng)的周期解p(t)也是漸近穩(wěn)定的.這是Poincaré關(guān)于平面自治系統(tǒng)存在周期解時解的軌道漸近穩(wěn)定性準(zhǔn)則的推廣.建立系統(tǒng)的一個Lyapunov函數(shù):
V(z1,z2,z3;p(t)):=sup{|z1|,p2(t)p3(t)(|z2|+|z3|)},顯然上面建立的V函數(shù)是恒正的,但是并非處處可微.用V的右導(dǎo)數(shù)修正,表示為D+V(t):
D+(|z1(t)|)≤(f′(p1)-(1-ηRT)k1p3(t)-δ)|z1(t)|+(1-ηRT)k1p1(t)p3(t)p2(t)p2(t)p3(t)(|z2(t)|+|z3(t)|)
D+(p2(t)p3(t)(|z2(t)|+|z3(t)|))=(p·2(t)p2(t)-p·3(t)p3(t))*p2(t)p3(t)(|z2(t)|+|z3(t)|)+p2(t)p3(t)D+(|z2(t)|+|z3(t)|)≤p2(t)p3(t)Nδ|z1(t)|+p2(t)p3(t)(|z2(t)|+|z3(t)|)(p·2(t)p2(t)-p·3(t)p3(t)-c)-p2(t)p3(t)(-f′(p1(t))|z2(t)|+δ|z3(t)|) ≤p2(t)p3(t)Nδ|z1(t)|+p2(t)p3(t)(|z2(t)|+|z3(t)|)(p·2(t)p2(t)-p·3(t)p3(t)-c-min(α*,δ))
0<α*=-maxT∈[0,T0]f′(T).定義這樣一個函數(shù):
g1(t)=f′(p1)-(1-ηRT)k1p3(t)-δ+(1-ηRT)k1p1(t)p3(t)/p2(t)=f′(p1)-(1-ηRT)k1p3(t)+p·2(t)/p2(t).
這里
g2(t)=p2(t)p3(t)Nδ+p·2(t)p2(t)-p·3(t)p3(t)-c-min(α*,δ)
=p·2(t)p2(t)-min(α*,δ).
根據(jù)上面兩個等式g1(t)和g2(t),以及p(t)滿足模型,我們易得:D+V(t)≤sup(g1(t),g2(t))V(t).因為g1(t)≤-α*+p2(t)/p2(t)·,g1(t)≤g2(t),所以有D+V(t)≤g2(t)V(t).為了證明D+V(t)<0,只須要有∫ω0g2(t)dt<0.
∫ω0g2(t)dt=∫ω0(p·2(t)p2(t)-min(α*,δ))dt<0,∫ω0p·2(t)p2(t)dt=∫ω0dp2(t)p2(t)=0
最終得到D+V(t)<0.因此我們證明了模型具有周期解時,解軌道的漸近穩(wěn)定性
定理6 當(dāng)R0>1,但模型(1)的地方病平衡點E1=(T1,T1*,V)不穩(wěn)定時,模型(1)在D內(nèi)部存在軌道漸近穩(wěn)定的周期解.
證明.利用RouthHurwitz準(zhǔn)則,容易證得模型(1)在R0>1和f′(T1)<0的條件下地方病平衡點是局部漸近穩(wěn)定的(證明過程見定理4).但如果E1不穩(wěn)定,對模型(1)作一個變換,令w=(X,Y,Z)=(-T,T*,-V),則(1)轉(zhuǎn)換為:
dXdt=-s+pX1+XTmax-dTX+(1-ηRT)k1XZdYdt=(1-ηRT)k1XZ-δYdVdt=-NδY-cZ(5)
模型(5)在w點處的Jacobian矩陣記為t>T(w0)>0(T(0),T*(0),V(0))=(800,500,1000)N=100k1=0.0003w(t,w0)∈R.
記集合W:={(X,Y,Z):X,Z≤0,Y≥0},在集合W中J(w)除對角線以外的元素均為非正,所以模型(5)在W中是三維競爭系統(tǒng).令w*=(X*,Y*,Z*),根據(jù)題設(shè)中的條件易得,w*是不穩(wěn)定的,并且det(J(w*))<0.另外,對于w0∈intW,W內(nèi)部都存在緊集R,使得當(dāng)t>T(w0)時,對w(t,w0)∈R,都有T(w0)>0.引用文獻[10]中的定理1.1和1.2,因為w*是不穩(wěn)定的,所以模型(7)至少存在一個軌道漸近穩(wěn)定的周期解.同時由于模型(5)是(1)的線性變換,所以模型(1)至少存在一個軌道漸近穩(wěn)定的周期解.
2.數(shù)值模擬
下圖反映了初值?。═(0),T*(0),V(0))=(800,500,1000)時,模型(1)存在一個軌道漸近穩(wěn)定的周期解.
T、T*和V隨時間的變化
圖中的參數(shù)取值為:s=5,dT=0.008,δ=0.3,Tmax=1500,N=100,c=2.5,ηRT=0.1,k1=0.0003.
本文小結(jié)
利用RouthHurwitz準(zhǔn)則和構(gòu)造Lyapunov函數(shù)的方法可證明了模型(1)中無病平衡點和地方病平衡點的穩(wěn)定性.模型(1)是一個三維競爭系統(tǒng),因此我們可以利用競爭系統(tǒng)的特性和論證周期解的軌道漸近穩(wěn)定性證明地方病平衡點的性質(zhì),并且利用一些小技巧找到了滿足周期解存在的參數(shù)條件,模擬出了地方病平衡點不穩(wěn)定時的周期解軌道.