不完全CI下的Wolbachia動力學行為

邱俊雄, 郭文亮, 鄭 波

(廣州大學 數學與信息科學學院, 廣東 廣州 510006)

不完全CI下的Wolbachia動力學行為

邱俊雄, 郭文亮, 鄭 波*

(廣州大學 數學與信息科學學院, 廣東 廣州 510006)

登革熱已經成為世界上最常見的媒介傳播病毒性疾病之一, 而利用感染Wolbachia蚊子去控制登革熱已經在實踐中證明是一項有效的策略. 文章討論了不完全CI下的Wolbachia的動力學行為，系統地分析了平衡點的存在條件和穩定性，精確地給出了不完全CI給成功傳播Wolbachia的閾值帶來的影響，所得的結論對釋放策略的制定具有一定的指導作用.

登革熱; Wolbachia; 不完全CI; 吸引域; 閾值

登革熱是登革病毒經伊蚊傳播引起的急性傳染病,主要包括登革出血熱和登革休克綜合癥.埃及伊蚊和白紋伊蚊是其主要載體.目前估計,在全球范圍內的病例數從1997年的1億增加到5億,使得它成為一個日益嚴重的公共健康問題. 但不幸的是,醫學上目前仍缺乏有效的藥物和疫苗來控制登革熱.而且因為抗藥性、環境等因素的影響,傳統的化學策略未能實現消滅傳播登革熱病毒的蚊子,所以現在我們把重點放在具有巨大生殖能力和基因靈活性的蚊子上. 2005年，密歇根州立大學奚志勇首次實現通過胚胎顯微把Wolbachia注射到埃及伊蚊體內,建立了穩定鏈,并且在實驗室已證實Wolbachia可以以較低的引入率擴散到整個蚊群[1].隨后,他們也成功地將Wolbachia注射到傳播瘧疾的斯氏按蚊中[2],從而使得運用Wolbachia控制蟲媒傳播疾病具有了更堅實的應用基礎和更廣闊的前景.實驗發現:①Wolbachia是一種垂直傳播的共生菌,攜帶Wolbachia的雌蚊的后代一般都攜帶Wolbachia; ②細胞質不親和性(CI)在攜帶Wolbachia的雄蚊與不攜帶Wolbachia的雌蚊交配時產生,造成的直接后果是這時候的受精卵不能完全甚至完全不能孵化.而攜帶Wolbachia的雌蚊不管與哪種雄蚊交配,受精卵都能正常孵化.這就給攜帶Wolbachia的雌蚊帶來了一定的生殖優勢.于是,在利用Wolbachia控制蟲媒傳播疾病時,出現了2種策略:①種群替換:同時釋放攜帶Wolbachia的雌蚊與雄蚊,使得整個野外蚊群被攜帶Wolbachia的蚊群完全替換;②種群壓制:只釋放攜帶Wolbachia的雄蚊,使得野外雌蚊與釋放的雄蚊交配后不能正常孵化后代,從而達到“殺死”或者壓制野外蚊群的目的. 而不管是哪一種釋放策略,數學模型的建立與分析都對實際釋放工作有著不可或缺的指導作用.

早在1959年,CASPARI等首次從數學的角度分析了細胞質不親和性給蚊子帶來的進化優勢[3].隨后,許多學者,如FINE, HOFFMANN, TURELLI和TURELLI等建立了離散世代不重疊的數學模型[4-7]. KEELING[8]與FARKAS[9]建立了常微分方程組來反映Wolbachia在蚊群中的傳播動力學行為,分析了平衡點及其穩定性,利用數值模擬對吸引域做了估計,但是仍舊缺少對“如何保證Wolbachia成功傳播”的數學分析與刻畫.另外,在制定釋放策略時,釋放公母蚊子的比例是有技巧的.因為只有母蚊子會咬人,所以希望放的母蚊子數量是盡可能少的.由此,將蚊子按照雌雄與攜帶Wolbachia分類,并考慮到不同性別蚊子的成熟期因素,ZHENG等[10]建立時滯微分方程模型來研究Wolbachia在種群中傳播的動力學行為,同時數值模擬顯示該模型能與2011年發表在Nature上的連續2篇論文的實驗數據吻合得非常好[11-12].在此基礎上,本文考慮不完全細胞質不親和(不完全CI)對“保證Wolbachia成功傳播的閾值”的數學刻畫.與完全CI不同,不完全CI導致未感染雌蚊與感染雄蚊交配后,仍有少部分未感染后代出現.CI的不完全會提高種群替換的閾值,降低種群替換的速度,從而了解其對Wolbachia的動力學行為的影響至關重要.

1 模型的建立及其平衡點的穩定性

Wolbachia的注射會改變蚊群的生殖能力與壽命[1-2], 為此, 令b1,b2分別為感染蚊子與未感染蚊子的出生率;δ1與δ2為感染蚊子與未感染蚊子的死亡率;令q(0

(1)

(2)

引理1系統(1)～(2)的從第一象限出發的任意解都是有界的.而且,若令

κ1=b1/δ1,κ2=b2/δ2

(3)

對于點(x,y)滿足x≥0,y≥0,其ω極限集

ω(x,y)?D:={(x,y):0≤x≤κ1,0≤y≤κ2}.

證明記x(t)和y(t)是系統(1)～(2)在第一象限出發的任意解.有

x′(t)=b1x-δ1x(x+y)

同時,

當x(t)>κ1,則x′(t)<0;若y(t)>κ2,則y′(t)<0.從而κ1(κ2)可視為(未)感染蚊群的環境容納量,且(x,y)的ω極限集包含在D中.引理證畢!

注1當κ1≥κ2,環境有利于感染蚊群,此時感染了Wolbachia的蚊群擁有適合度優勢(fitnessbenefit).當κ1<κ2,環境有利于未感染蚊群,此時感染了Wolbachia的蚊群擁有適合度劣勢(fitnesscost).

引理2系統(1)～(2)在x軸上存在平衡點E0=(0,0)和E1=(κ1,0),在y軸上存在平衡點E2=(0,κ2).其中,E0是一個源; 當κ2(1-q)<κ1,E1是一個局部漸近穩定點,當κ2(1-q)>κ1,E1是一個鞍點;當κ1<κ2,E2是一個漸近穩定點,當κ1>κ2,E2是一個鞍點.

證明易得系統存在下面的平衡點:E0=(0,0),E1=(κ1,0),E2=(0,κ2). 注意到在x軸和y軸上,系統(1)～(2)退化為

x′=δ1x(κ1-x),y′=δ2y(κ2-y)

(4)

由式(4)可知,從x軸上與y軸上出發的解(原點附近出發)將遠離E0=(0,0),因此，由常微分方程解對初值的連續依賴性可知E0=(0,0)是一個源.

系統(1)～(2)的Jacobian矩陣:

J(x,y)=

從應用的角度考慮,以下只考慮κ2(1-q)<κ1<κ2,即保證E1和E2都是漸近穩定的平衡點,對應于競爭系統感染蚊群替換掉未感染蚊群的成功與否.

引理3若κ2(1-q)<κ1<κ2,系統(1)～(2)存在唯一的內部鞍點

證明令式(1)的左邊等于0,除去x=0,有

x+y=κ1

(5)

把式(5)代入式(2),可得

(6)

由此可得唯一的內部平衡點為

E*=(x*,y*)=

在E*的Jacobian矩陣是:

進而可以得到保證E*是鞍點的必要條件:

detJ(E*)=

-δ1x*δ2(κ1-κ2(1-q))<0.

2 保證Wolbachia成功傳播的閾值分析——E1的吸引域

定理1(1)存在E*的穩定曲線H:y=h(x),x≥0滿足

曲線在0x*部分定義了E*的穩定曲線.

(2)定義

則可進一步估計y=h(x):當δ1=δ2時,有

h0(x)=h(x)=h1(x).

當δ1>δ2時,有:

h0(x)

h1(x)x*.

當δ1<δ2時,有:

h1(x)

h0(x)x*.

注2定理1說明,若記E1和E2的吸引域分別為B(E1)和B(E2)，則有

B(E1)={(x,y):x≥0,0≤y

B(E2)={(x,y):x≥0,y>h(x)}.

證明(1)x-零漸近線x+y=κ1和y-零漸近線b2(1-qx/(x+y))=δ2(x+y)將第一象限分成了4個部分,見圖1、圖2.

圖1 向量場示意圖Fig.1 Vector fields schematic

圖2 估計穩定曲線示意圖Fig.2 Estimation of stable curve

通過分析向量場的方向,能夠知道區域D4有x′(t)<0,y′(t)>0.當t→∞,從區域D4出發的任意解(x(t),y(t))將會趨向E2,且D4?Β(E2).同理,D2?Β(E1).由D4?Β(E2)和D2?Β(E1),E*的穩定曲線必須位于區域D1和D3中.在區域D1和D3中,發現dy/dx>0,因此，y隨x的增大而增大,從而穩定的曲線能夠表示成一條光滑且嚴格單調遞增的函數y=h(x).當0

(2)只證明δ1=δ2與δ1>δ2的情況.δ1<δ2的情形與δ1>δ2類似,略去.

為此,記連接E0和E*的直線:

L:l(x,y)=(κ1-κ2(1-q))x-(κ2-κ1)y=0.

則系統(1)～(2)沿著上述直線的導數:

(κ2-κ1)y(b1-δ1(x+y))-

κ2q(δ2-δ1)(x-x*)y

(7)

當δ1=δ2，表明L是系統(1)～(2)的解.同時也意味著函數h0(x)，h(x)和h1(x)是相等的.

接下來證明δ1>δ2的情況.將證明:當0x*時,連接鞍點E*與原點的穩定曲線的下界是y=h1(x),上界是y=h0(x).由式(7)可以得到:

(8)

(9)

進而由向量場分析易知(圖2)

h(x)>h0(x) 若0

h(x)x*

(10)

另一方面,令

由于δ1>δ2,從而通過鞍點E*的曲線z(t)=c與直線y=h0(x)的位置關系由圖2確定,則系統(1)～(2)的解沿著z(t)的導數滿足

另一方面,直線l(x,y)=0的斜率為

所以,

(11)

若0h0(x),從而必有z′(t)<0,這意味著z(t)沿著異宿軌y=h(x)移動的時候是單調遞減的,即

log((x*)1/δ1(y*)-1/δ2).

這等價于:

從而,若δ1>δ2,且0

h(x)

(12)

而若x>x*,由式(10)已知h(x)0(圖2),意味著z(t)沿著異宿軌y=h(x)移動的時候是遞增的,即

log((x*)1/δ1(y*)-1/δ2).

這等價于:

(13)

結合式(10),式(12)和式(13)，定理證畢！

通過釋放Wolbachia雄蚊，利用CI的特性去壓制野外蚊群時，一個顯然的結論是CI強度直接決定了種群壓制的速度.而在種群替換時，CI強度不僅會影響替換速度，更重要的是，它也直接影響了使替換成功的最小釋放數目.為此，采用文獻[10]中的2種Wolbachia株(wMelPop, wMel)的參數，分析在不同CI的參數條件下吸引域的變化情況.

例1實驗發現[1-2]，wMel對于雌蚊的生殖力幾乎沒有影響，相比于未感染的雌蚊，感染了wMel的雌蚊的壽命減少了大概10%，采用文獻[10]的參數，令

b1=b2=0.198 8,δ1=9.448 2×10-6,

δ2=8.503 4×10-6.

而wMelPop不光大大降低了雌蚊的生殖力，而且感染了wMelPop的雌蚊的壽命減少了大約40%，此時，文獻[10]中的參數如下：

b1=0.107 7,b2=0.198 8,

δ1=1.417 2×10-5,δ2=8.503 4×10-6.

圖3顯示在不同CI強度下，通過數值模擬得到的曲線是E1和E2的吸引域的分界線.可知：隨著q值(CI強度)的增大，E1的吸引域在變大，而E2的吸引域在變小.這也意味著在未感染蚊子數量是一定的時候，CI強度越大，要達到種群替換需要釋放感染的蚊子的數量越少.

圖3 CI強度對Wolbachia感染的吸引域的影響Fig.3 Effect of CI intensity on the attraction basin of Wolbachia infection

[1] XI Z, KHOO C, DOBSON S L. Wolbachia establishment and invasion in an Aedes aegypti laboratory population[J]. Science, 2005, 310: 326-328.

[2] BIAN G, JOSHI D, DONG Y, et al. Wolbachia invades Anopheles stephensi populations and induces refractoriness to Plasmodium infection[J]. Science, 2013,340: 748-751.

[3] CASPARI E, WATSON G S. On the evolutionary importance of cytoplasmic sterility in mosquitoes[J]. Evolution, 1959,13(4): 568-570.

[4] FINE P E. On the dynamics of symbiote-dependent cytoplasmic incompatibility in culicine mosquitoes[J]. J Invert Pathol, 1978,31(1): 10-18.

[5] TURELLI M, HOFFMANN A A. Rapid spread of an inherited incompatibility factor in California Drosophila[J]. Nature,1991, 353:440-442.

[6] TURELLI M, HOFFMANN A A. Microbe-induced cytoplasmic incompatibility as a mechanism for introducing transgenes into arthropod populations[J]. Insect Mol Biol, 1999,8(2): 243-255.

[7] TURELLI M. Cytoplasmic incompatibility in populations with overlapping generations[J]. Evolution, 2010,64(1): 232-241.

[8] KEELING M J, JIGGINS F M, READ J M. The invasion and coexistence of competing Wolbachia strains[J]. Heredity, 2003,91(4):382-388.

[9] FARKAS J Z, HINOW P. Structured and unstructured continuous models for Wolbachia infections[J]. Bull Math Biol, 2010,72(8):2067-2088.

[10]ZHENG B, TANG M, YU J. Modeling Wolbachia spread in mosquitoes through delay differential equations[J]. SIAM J Appl Math, 2014,74(3):743-770.

[11]HOFFMANN A A, MONTGOMERY B L, POPOVICI J, et al. Successful establishment of Wolbachia in Aedes populations to suppress dengue transmission[J]. Nature, 2011, 476: 454-457.

[12]WALKER T, JOHNSON P H, MOREIRA L A, et al. The wMel Wolbachia strain blocks dengue and invades caged Aedes aegypti populations[J]. Nature, 2011, 476: 450-453.

【責任編輯： 周 全】

Dynamics of Wolbachia under incomplete CI

QIU Jun-xiong, GUO Wen-liang, ZHENG Bo

(School of Mathematics and Information Sciences, Guangzhou University, Guangzhou 510006, China)

Dengue fever has become one of the world’s most common vector-borne viral diseases. A novel method to control dengue is to release Wolbachia mosquitoes which has been proved to be effective. In this paper, we discuss the dynamics of Wolbachia with incomplete CI. After analyzing the stability and instability of equilibrium points, we sharply describe the effect of incomplete CI on the threshold value for successful Wolbachia invasion which could offer an intuitive guide for mosquito releasing to ensure the success.

dengue fever; Wolbachia; incomplete CI; attraction basin; threshold value

2016-04-09;

2016-05-04

國家自然科學基金資助項目(11301103)；廣東省高等學校優秀青年教師培養計劃資助項目(Yq2013131)

邱俊雄(1990-),男,碩士研究生.E-mail：405259872@qq.com

*通信作者.E-mail：zhengbo611@outlook.com

1671- 4229(2016)06-0034-05

O 175.1