乙肝免疫球蛋白阻斷乙肝病毒在母嬰間傳播的動力學模型

李冬梅, 劉偉華, 汪 琪, 郭美靜

(1- 哈爾濱理工大學理學院，哈爾濱 150080; 2- 黑龍江省醫院兒科，哈爾濱 150036)

1 引言

乙型病毒(Hepatitis B Virus, HBV)性肝炎，簡稱乙肝，是由HBV 感染引起的傳染病.母嬰傳播是HBV 傳播的三種方式之一.實施母嬰傳播阻斷策略能夠顯著降低HBV 感染，且已經成為現階段我國乙肝預防工作的重點[1].應用乙肝免疫球蛋白(Hepatitis B Immunoglobulin, HBIG)可有效的阻斷HBV 通過胎盤對胎兒的感染，進一步提高出生人口的健康質量[2].

眾多數學工作者通過研究HBV 的發病機制和治療機理建立了數學模型，并希望通過定量研究HBV 變化規律對防治乙肝疾病帶來一定的幫助[3–5].Lewin 和Ribeiro[6]通過數學模型進一步研究了抗病毒藥物治療對乙肝疾病的控制問題.Stanca 等人建立了具有抗體免疫反應的數學模型，這使得人們更加精確地刻畫了宿主體內的HBV 變化規律[7].實際上，在臨床醫學中是通過觀察母體內的體征情況來預測胎兒的患病情況.周煜等人[8]運用Meta 分析的方法對臨床檢驗的數據分析來評價藥物阻斷HBV 母嬰宮內免疫策略的可行性，發現HBV 攜帶孕婦注射HBIG 可降低新生兒的感染率，但目前對于母嬰之間HBV 傳播的定量模型研究結果很少.

本文綜合考慮HBV 進入母體后的傳播機理，以及藥物對病毒的抑制作用，建立了母體懷孕期間的藥物阻斷HBV 在母嬰間傳播的動力學模型，研究了母體接受HBIG 來控制體內HBV 的數量，從而減少HBV 對胎兒的感染機率.再根據臨床用藥量數據，運用Matlab 模擬胎兒體內的病毒數量變化曲線，從理論上優化給藥方案，為指導臨床治療提供了可參考的依據.

2 乙肝免疫球蛋白阻斷乙肝病毒在母嬰間傳播的動力學模型的建立

研究發現，當個體感染了HBV 后，其作為抗原刺激宿主細胞的免疫系統產生抗體來阻滯感染HBV.若免疫系統不強無法抵抗HBV 復制，就會帶來慢性持續性感染，從而演變成乙肝患者.乙肝免疫球蛋白(HBIG)是預防HBV 入侵復制的被動免疫制劑，當人體接受這種外源性抗體后，可以在短期內中和并清除血清中游離的HBV.臨床上通常采用HBIG 阻斷HBV 在母嬰間的傳播，有關方面的研究已取得了可觀的結果[9,10].這為建立數學模型提供了可依據的理論基礎.在妊娠初期，由于胎盤的屏障作用，母體與胎兒血液分開，互不干擾，除了可以進行選擇性的物質交換外，胎盤還可以阻止母體內的一部分病毒進入胚胎、保障胎兒正常生長發育.但在孕婦妊娠20 周后，胎盤開始傳遞母體抗體給胎兒的功能，但胎兒體內的抗體不足以抵抗母體傳遞的病毒，為此要控制母體內的病毒數量.通過應用HBIG 來降低母體內HBV 數量，可以有效阻斷HBV 感染胎兒[7].由文獻[6]從數學和醫學角度做出如下假設：

1) 胎兒體內的HBV 主要是來源于母體內HBV 輸入；

2) 在胎兒的生長發育過程中，其免疫系統尚不完善，不足以抵御外來病毒的入侵.為了降低胎兒的感染可能性，藥物可以減少母體內HBV，從而阻斷HBV 在母嬰間的垂直傳播.因而不考慮胎兒免疫系統功能；

3) 每個個體乙肝傳播機理相同；

4) 認為HBV 攜帶者與感染者之間存在著臨界值.若個體的HBV 超過臨界值，可診斷為HBV 感染者.

基于上述原因與假設，將母體與胎兒看成兩個個體，胎盤連接兩個個體，胎兒是否感染主要是受母體HBV 量的影響.下面分別討論母體及胎兒體內的易感肝細胞(Ti)、感染肝細胞(Ii)、HBV(Vi)和免疫抗體(C1)的情況，其中i= 1,2，下角標為1 代表母體、下角標為2 代表胎兒.考慮到HBIG 對HBV 的作用機理，建立了如下的動力學模型

模型(1)中相關參數意義，如表1 所示.

表1 模型(1)相關參數表

3 主要研究結果

在臨床實踐中，對母體應用HBIG 阻斷HBV 在母嬰之間傳播，藥物可以有效控制母親體內病毒的數量，使母體內的病毒不傳播給胎兒.因此，先考慮母體在藥物作用下病毒數量變化趨勢，再分析胎兒體內HBV感染可能性.可將模型(1)分解成如下的母體模型、胎兒模型

將母體模型(2)的第一、二、四方程相加，令m=min{dT1,dI1,dC1}，則由

可知，存在t1，當t ≥t1時，有

再由模型(2)的第三個方程可求得，存在t2>t1，當t ≥t2時，有

且Γ 是模型(2)的不變集.

3.1 母體模型(2)的穩定性

模型(2)在平衡點處滿足如下方程

由式(4)求得

令

將式(5)代入式(4)的第三個方程，得到

其中

由式(5)、式(6)易求得無病毒平衡點

由二次函數性質知，式(6)確定的函數為f(V1) = 0.當R0> 1 時，有唯一正根V ?1> 0，使得f(V ?1) = 0.再將V ?1代入式(5)，可計算出模型(2)有病毒平衡點Q?1=(T1?,I1?,V1?,C1?).

定理1當R0< 1 時，模型(2)存在局部漸近穩定的無病毒平衡點Q?0；當R0>1 時，模型(2)無病毒平衡點Q?0是不穩定.

證明 模型(2)的Jacobi 矩陣為

模型(2)在Q?0處線性近似方程對應的特征方程為

當R0< 1 時，由Routh-Hurwitz 定理[11]知，式(8)的所有特征根均為負特征根，無病毒平衡點是局部漸近穩定點；當R0>1 時，無病毒平衡點是不穩定點.

考慮系統

其中X1∈1, X2∈, X=(X1,X2), X?=(,0)為平衡點.假設下列條件：

H2:·X1=A1(X)(X1?)系統(9)的子系統存在唯一的全局漸近穩定點；

H3: 對任意的X ∈Γ，A2(X)是Metzler 矩陣且不可約；

H4: 定義在Γ 上的A2(X)存在一個上界矩陣，即存在∈Γ，使得

H5: 矩陣的譜半徑α()<0.

引理1[12]如果系統(9)滿足條件H1–H5，那么平衡點X?是系統(9)的全局漸近穩定點.

為了研究模型(2)無病毒平衡點的全局穩定性，將模型(2)寫成具有形式(9)的矩陣方程形式.令

定理2當R0<1 時，模型(2)的無病毒平衡點是全局漸近穩定點.

證明 下面驗證引理1 的條件.

H1: 由于Γ 是模型(2)的不變集，因而模型(2)在Γ 內是耗散的.

H2: 對于子系統

·X1=A1(X)(X1?X?1)，即有

顯然

是子系統(10)的全局漸近穩定點.

H3: 對于矩陣A2(X)，非主對角線上元素均為正，故A2(X)是Metzler 矩陣且不可約[13].

H4: 由A2(X)的定義及式(10)的穩定性可知，A2(X)在Γ 內存在上界矩陣，取=(X?1,0)，可得

H5: ˉA2是常數矩陣，要證α()<0，等價于α(D ?CA?1B)<0[11].

由引理1 知，當R0<1 時，模型(2)的無病毒平衡點是全局漸近穩定點.

定理3當R0>1 時，模型(2)在Γ 內是一致持久的.

證明 由定理1 知，當R0>1 時，在區域Γ 內模型(2)存在著兩個平衡點Q?0和Q?1，且Q?0是不穩定點.根據式(8)計算出模型(2)的特征值知，只有?Γ 上初值(T10,0,0,C10)出發的軌線最終趨于唯一的平衡點Q?0，其它軌線都將遠離Q?0.因此Q?0是?Γ 上孤立的最大緊的不變集.因而，存在正的常數a，對于從在Γ 的內部出發的模型(2)解有

由此模型(2)的解不會跑出?Γ，從而，有Ws({})∩?Γ=?，其中

綜上，當R0>1 時，模型(2)在空間區域Γ 內是一致持久的[13].

現考慮系統

其中F為Rn的連續函數，D ?Rn是連續的開集，X ∈Rn.令X(t,X0)是系統(11)的解，初始條件為X(0)=X0，假設：

1) 系統(11)存在吸引子集K ?D；

2) 系統(11)在D內有唯一的平衡點X?.

引理2[14]如果假設(1)、(2)成立，且存在Lyapunov 函數L(X)、有界可導函數G(X)、正常數a1、g和a2滿足：

(i)a1|X|≤L ≤a2|X|;

(ii)L′(X)≤(G′(X)?g)L(X);

則系統(11)的平衡點X?是全局漸近穩定點.

定理4當R0>1 時，模型(2)的有病毒平衡點Q?1是全局漸近穩定點.

證明 由式(7)，可得模型(2)在Q?1處線性近似方程對應的特征方程為

則由式(12)知，模型(2)的所有特征根均為負特征根.

由Routh-Hurwitz 定理[11]知，當R0>1 時，有病毒平衡點Q?1是局部漸近穩定點.

再由式(7)可計算出模型(2)的第三加性復合矩陣為[15]

其中

a11=?b1V1?dT1?dI1?εC1?dV1, a22=?b1V1?dT1?dI1?kV1?dC1

a33=?b1V1?dT1?εC1?dV1?kV1?dC1, a44=?εC1?dI1?dV1?kV1?dC1.

式(13)對應的線性復合系統為

構造Lyapunov 函數

L(t,X,Y,W,Z)=max{C1|X|, V1|Y|, I1|W+Z|, I1|W|, I1|Z|}.

由定理3，在Γ 內存在a1>0, a2>0，使得

下面計算L的導數.

1) 若取L=C1|X|時，有V1|Y|≤C1|X|，可計算

其中G1=lnC1, g1=dT1+dI1+dV1.

2) 若取L=V1|Y|時，有C1|X|≤V1|Y|，則可計算

其中G2=lnI1+lnV1, g2=dT1+dC1.

3) 若取L=I1|W+Z|時，有V1|Y|≤I1|W+Z|，則可計算

其中g3=dT1+dI1+dC1.

4) 若取L=I1|W|時，有V1|Y|≤I1|W|，則可計算

其中g4=dT1+dC1.

5) 若取L=I1|Z|時，則可計算

其中G3=lnI1, g5=dI1+dV1+dC1.

由式(16)–(20)，可有

D+L ≤[G′?g]L,

其中G′=max{G1,G2}, g=min{g1,g2,g3,g4,g5}.

由引理2 知，模型(2)的有病毒平衡點Q?2是全局漸近穩定點.

3.2 胎兒模型(3)的穩定性

在假設(1)的約定下，通過藥物治療后，母體內的HBV 最終趨于V ?1，可以認為感染的母體輸入給胎兒的HBV 量為ˉη=ηV ?1，由式(3)可得如下胎兒體內HBV 模型

同母體可行域Γ 的分析，可以求得方程(21)可行域為

其中

方程(21)的平衡點滿足下列方程

整理得

將式(23)代入式(22)中第二個方程，得

由二次函數性質知，f(I2)=0 存在唯一正根，使得f()=0.再將代入式(23)中，可得到方程(21)唯一的有病毒平衡點.

定理5當R0>1 時，模型(21)存在唯一的有病毒平衡點是全局漸近穩定的.

證明 模型(21)的Jacobi 矩陣為

計算

顯然有a1> 0, a3> 0，且a1a2?a3> 0.由Routh-Hurwitz 定理[11]知，有病毒平衡點是局部漸近穩定的.

模型(21)對應以ω為周期的任意周期解P(t) = (T2(t),T2(t),V2(t))的二階復合系統為[15]

做Lyapunove 函數

由區域D內只含有穩定平衡點Q?2，故存在c1，使得L ≥c1sup{|X|,|Y|,|Z|}.

將對(24)式求右導數，有

其中α=min{dT2,dI2}.

由式(25)得

其中

再由式(21)可知

由式(26)可得

取

由式(27)知

由此可知，方程(24)零解是穩定的.

另一方面，取H= diag(?1,1,?1)可證明的非對角線元素非正，模型(21)是競爭系統，即模型(21)滿足Poincare-Bendixson 性質，由文獻[11]知，模型(21)有病毒平衡點是全局漸近穩定.

4 模型應用

目前，用HBIG 阻斷攜帶HBV 的母體感染胎兒的治療方案沒有明確統一的標準，但已取得了階段性的進展.臨床實踐中，通過檢測治療過程中母體的HBV-DNA 數據(其與HBV 存在正相關關系)來確定治療效果.并認為HBV 攜帶者與感染者之間存在著臨界值，若檢測結果沒有超過臨界值，可以視為HBV 轉陰，無傳染性可能.因此，用藥原則遵循著母體內HBV 不應超過臨界[7].根據檢測到母體HBV-DNA 數據，可以定性推測胎兒感染HBV 可能性.用數學模型可以定量分析母嬰體內HBV 數量變化趨勢，預測胎兒體內HBV 數量的變化情況.并且通過調整模型的參數，進一步優選用藥量使得母體的HBV 控制在臨界值之下.其中HBV-DNA 的臨界值為277，HBV 的臨界值為1000[7].模型(1)中參數見表2[16].

表2 模型相關參數數值表

下面根據模型(1)結果和表2 中乙肝患者的相關參數選取，針對于HBIG 治療HBV 孕婦案例，分析模型結果的合理性，數值模擬HBIG 劑量對母體內HBV 的數量影響，預測胎兒感染HBV 數量的變化情況.

案例：孕婦懷孕第20 周，孕期檢查時發現孕婦體內攜帶HBV，檢測HBV-DNA 沒有超過臨界值為277，認為還沒有發病.為避免將病毒傳染給胎兒，在家屬知情的情況下，對孕婦進行注射HBIG 的阻斷治療.治療期間檢測到母體內HBV-DNA 含量如下表3.采取的治療方案是標準用藥方案，即每月注射一針HBIG，每次劑量200 IU，直至分娩[17].

表3 孕婦血清中HBV-DNA 含量

由表3 中數據可以繪制出治療期間孕婦血清中HBV-DNA 含量的曲線如下圖1.

由圖1 可以看出，在懷孕20 周時開始用藥，母體內HBV-DNA 含量先增加，到了第22 周左右時達到峰值，隨著時間的推移，母體內HBV-DNA 含量逐漸減少，直至40 周后趨于平穩.

圖1 母體內HBV-DNA 含量

4.1 模型的合理性分析

結合案例數據，用模型(1)數值模擬母體HBV 數量變化趨勢，分析模型結果的合理性.選取母體內的HBV 數量的初值V1(20)<，將表2 中數據代入模型(1)，計算出R0= 1.05> 1，由定理4 知病毒平衡點是全局漸近穩定點.應用Matlab 軟件對模型(1)分別數值模擬出母體、胎兒體內的HBV 數量隨時間變化的曲線圖像如下圖2 和圖3.

圖2 當μ2 =200IU 時，母體內V1(t)隨時間的變化曲線

圖3 當μ2 =200IU 時，胎兒體內V2(t)隨時間的變化曲線

從圖2 可以看出，在20 周時接受治療，HBV 數量增加到第22 周時達到峰值，隨后逐漸降低.圖2 的HBV 數量曲線變化趨勢與圖1 反映的HBV-DNA 含量曲線變化趨勢基本相同.在治療期間內模型給出HBV 數量變化趨勢與臨床實際監測HBV-DNA 含量變化趨勢一致，說明模型是合理的.同時，由圖2 的數據可知，母體內HBV 達到的峰值數量900，低于在治療期間(20 周–45 周)平均HBV 數量為483，HBV 數量在可控的范圍內，降低了胎兒的感染HBV 的風險.由圖3 數據可計算出孕期胎兒的HBV 數量是362.

4.2 藥物劑量對攜帶HBV 孕婦病毒數量的影響

當孕婦攜帶HBV 的數量沒有超過臨界值數量時，胎兒感染HBV 風險很大.臨床上應用HBIG 阻斷HBV 在母嬰間傳播，4.1 中給出了標準用藥方案下孕婦、胎兒體內HBV 數量，能夠有效控制HBV 數量不超出臨界數值.在不超出最小中毒劑量下，當增加用藥劑量至時400IU[17]，分析數值模擬HBV數量的變化趨勢.選取母體內的HBV數量的初值V1(20) = 600<，其它參數仍取表2，計算出R0= 0.64<1.由定理3 知，無病毒平衡點是全局漸近穩定點.由于孕周期時間是有限的，在此只考慮45 周之內孕婦體內HBV 數量變化情況，用Matlab 數值模擬出母體、胎兒體內HBV 數量隨時間變化的曲線圖像如下圖4 和圖5.

由圖4 數值可知，當改變藥物劑量為時400IU，母體內HBV 數量峰值795 低于1000，平均病毒數量367.對比圖2、圖4，發現在孕期內HBV 數量均下降，說明調整劑量后的用藥方案能更好的控制體內的病毒數量.由圖5 數值可以計算出胎兒體內的HBV 的數量是340.由此說明增加藥劑量可以快速降低母體內HBV 數量，胎兒的感染風險減小，達到更好的治療效果.

圖4 當μ2 =400IU 時，母體內V1(t)隨時間的變化曲線

圖5 當μ2 =400IU 時，胎兒體內V2(t)隨時間的變化曲線

4.3 藥物劑量對感染HBV 孕婦病毒數量的影響

圖6 不同劑量下，母體內V1(t)隨時間的變化曲線

圖7 不同劑量下，胎兒體內V2(t)隨時間的變化曲線

從圖6 可以看出，隨著藥物劑量的增大，母體內HBV 的數量下降的速度越快，且在第12 周時HBV 數量分別為1050、830、810，說明給藥劑量為200IU 時，不能嚴格的控制母體HBV 傳染給胎兒，不能達到阻斷的效果，通過劑量為400IU 與500IU 對比知，500IU 的母體內的HBV 數量下降得更快，則在臨床實踐中，在患者經濟允許情況下，可適當地加大藥物的治療劑量，做出針對性的治療方案.