http:∥xuebao.gxnu.edu.cn

數(shù)是使用 Tanimoto 系數(shù)及其化學(xué)結(jié)構(gòu)的乘積[45]計算得到的);

⑥ 蛋白質(zhì)序列相似性網(wǎng)絡(luò)(基于基因組序列使用成對的史密斯-沃特曼得分[46]獲得蛋白質(zhì)相似性網(wǎng)

絡(luò),用于訓(xùn)練卷積網(wǎng)絡(luò)以更新藥物靶標(biāo)特征)。

2.2 模型性能

現(xiàn)有的藥物靶標(biāo)預(yù)測主要將藥物靶標(biāo)的已知相互作用作為陽性實例,未知相互作用作為陰性實例。

實驗中采用 10 倍交叉驗證,并隨機(jī)選擇 10%的數(shù)據(jù)作為測試集,其余 90%的數(shù)據(jù)作為訓(xùn)練集。 實驗將

NGDTI 與 NeoDTI

[31]

、DTINet

[22]

、BLMNI

[12]

、NetLapRLS

[47]和 HNM

[48] 等 5 種方法進(jìn)行比較。 在實驗中使

用 AUPR(精確召回曲線下的面積)來衡量 NGDTI 的預(yù)測效果。 從圖 3( a)可知,NGDTI 比其他方法效果

更好,其中 AUPR 比最佳方法高 0.01。 與 DTINet 相比,盡管使用的都是降維之后的網(wǎng)絡(luò)擴(kuò)散特征,但

NGDTI 方法在此基礎(chǔ)上進(jìn)一步使用圖卷積模型優(yōu)化特征,從而獲得更好的結(jié)果。 在正負(fù)樣本比例設(shè)置為

1 ∶ 1 的情況下,實驗結(jié)果顯示 NeoDTI 比 DTINet 更糟糕。 為了驗證 NGDTI 在稀疏陽性樣本下的性能,實

驗對樣本數(shù)進(jìn)行修改,并為陽性和陰性樣本指定了 1 ∶ 10 的比例。 從圖 3(b)可知,每種算法的性能都下

降了。 相比之下,NGDTI 仍然取得了最佳的預(yù)測性能,且比對比方法中最好的 AUPR 值多 0.015。 這表

明,即使在標(biāo)記稀疏的情況下,對比方法的預(yù)測性能仍然不如 NGDTI 方法。

由于數(shù)據(jù)可能是冗余的,例如,數(shù)據(jù)集中存在一種蛋白質(zhì)的多種同源蛋白質(zhì)或一種藥物的多種高度相

似的藥物,這可能會對預(yù)測性能產(chǎn)生負(fù)面影響。 因此,本文采用與 Luo 等[22]相同的策略,通過在藥物靶標(biāo)

矩陣中刪除那些具有相似藥物或靶標(biāo)的藥物-靶標(biāo)關(guān)聯(lián)來減少數(shù)據(jù)冗余的影響。 實驗還在消除藥物-靶標(biāo)

關(guān)聯(lián)情況下進(jìn)行測試,其中關(guān)聯(lián)網(wǎng)絡(luò)中的 Jaccard 相似度大于 0.6,藥物化學(xué)相似性網(wǎng)絡(luò)中的結(jié)構(gòu)相似性得

分超過 0.6,蛋白質(zhì)-蛋白質(zhì)序列相似性網(wǎng)絡(luò)中的同一性得分超過 0.4。 在這些實驗中,陰性和陽性樣本的

比例保持為 1 ∶ 1。 從圖 3(c)(d)(e)(f)的實驗結(jié)果來看,雖然冗余的藥物靶標(biāo)關(guān)聯(lián)刪除后,NGDTI 性能

下降,但其仍優(yōu)于其他預(yù)測方法。

圖 3 NGDTI 方法與其他方法的比較

Fig. 3 Comparison of NGDTI with other baseline methods

97

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

2.3 模型參數(shù)的影響

本文提出的多網(wǎng)絡(luò)集成的藥物靶標(biāo)預(yù)測算法 NGDTI 的核心是使用變分圖自編碼器(VGAE)來更新

藥物和靶點特征。 與文獻(xiàn)[49]中使用圖卷積神經(jīng)網(wǎng)絡(luò)(GCN)的目的不同,NGDTI 是通過 GCN 學(xué)習(xí)藥物

和蛋白質(zhì)的網(wǎng)絡(luò)特征,使得藥物和蛋白質(zhì)的特征編碼含更豐富的生物信息,而不是通過 GCN 來學(xué)習(xí)每種

藥物的分子結(jié)構(gòu)。 VGAE 能夠聚合鄰域特征來進(jìn)一步提高特征的可用性,它使用基于譜的圖卷積網(wǎng)絡(luò)

(GCN)方法從圖信號處理的角度引入濾波器來定義圖卷積,其中的圖卷積操作可以被認(rèn)為是從圖信號中

去除噪聲。 為了驗證 VGAE 部件的有效性,實驗中實現(xiàn)了一個不包含 VGAE 部件的多網(wǎng)絡(luò)集成框架來預(yù)

測藥物-靶標(biāo)相互作用。 實驗首先對 NGDTI 中是否有 VGAE 部件、藥物特征維數(shù)、蛋白質(zhì)特征維數(shù)進(jìn)行不

同設(shè)置,并比較在不同條件下的預(yù)測效果,實驗結(jié)果如表 1 所示。 從表 1 可以看出,在藥物和蛋白質(zhì)特征

維數(shù)相同的情況下,有 VGAE 部件的 NGDTI 方法預(yù)測效果更好,而且在藥物特征維數(shù)為 100、蛋白質(zhì)特征

維數(shù)為 400 時預(yù)測效果最好。

表 1 NGDTI 在不同設(shè)置下的預(yù)測性能(正負(fù)樣本比例為 1 ∶ 1)

Tab. 1 Prediction performance of NGDTI under different settings (ratio of positive and negative samples is 1 ∶ 1)

是否有 VGAE 部件 藥物特征維數(shù) 蛋白質(zhì)特征維數(shù) AUPR AUROC

否 100 200 0.889 0.863

是 100 200 0.901 0.880

否 200 200 0.894 0.875

是 200 200 0.914 0.895

否 100 400 0.924 0.904

是 100 400 0.943 0.901

否 200 400 0.921 0.900

是 200 400 0.928 0.910

在之后的實驗中,主要評估參數(shù)的影響和 NGDTI 的魯棒性。 這些實驗通過更改與藥物或靶標(biāo)相關(guān)的

網(wǎng)絡(luò)數(shù)量以及 NGDTI 的超參數(shù)來測試 NGDTI 的魯棒性。 所有實驗結(jié)果均通過多次實驗取平均獲得。

首先,通過實驗驗證聚合多個異構(gòu)網(wǎng)絡(luò)對預(yù)測結(jié)果的影響。 在實驗中只使用部分網(wǎng)絡(luò)的情況下進(jìn)行

性能評估,同時將預(yù)測結(jié)果與使用所有網(wǎng)絡(luò)的預(yù)測結(jié)果進(jìn)行對比,結(jié)果如圖 4(a)所示。 通過結(jié)果可以觀

察到隨著多數(shù)據(jù)源數(shù)據(jù)的整合,預(yù)測性能顯著提高。 在添加了疾病和副作用相關(guān)信息的網(wǎng)絡(luò)之后,模型的

預(yù)測效果也得到改進(jìn),這也表明整合多數(shù)據(jù)源數(shù)據(jù)的有效性(NGDTI 可以整合多種蛋白質(zhì)或藥物相關(guān)數(shù)

據(jù)來改善預(yù)測性能)。

圖 4 整合更多與藥物或靶標(biāo)相關(guān)的信息的效果和重啟隨機(jī)游走概率 p 的影響

Fig. 4 Effect of integrating more information related to the drug or target and the effect of restarting the random walk probability

98

http:∥xuebao.gxnu.edu.cn

此外,本文還探索模型超參數(shù)對實驗性能的影響。 在這里,主要研究隨機(jī)游走重啟概率 p 對實驗結(jié)果

的影響。 如圖 4(b)所示,在測試中改變重新啟動概率值為 0.4 到 0.7,以觀察不同概率下的性能穩(wěn)定性。

在圖 4(b)中可以看出,當(dāng)重啟概率為 0.4 至 0.7 時,NGDTI 實現(xiàn)了穩(wěn)定的性能。 從以上實驗可以得出,模

型的參數(shù)對實驗性能的影響較小。

2.4 NGDTI 預(yù)測的藥物靶標(biāo)相互作用

最終的預(yù)測結(jié)果選取可信度排名前 10 位的藥物-靶標(biāo)相互作用,其中有 4 個藥物-靶標(biāo)相互作用有相

關(guān)文獻(xiàn)研究的支持。 例如,nifedipine 是一種被批準(zhǔn)用于輔助治療高原肺水腫的藥物,而 NR3C1 的多態(tài)性

與高原肺水腫有著重要的關(guān)聯(lián)。 這一預(yù)測可以被先前的一項研究支持,該研究表明 NR3C1 多態(tài)性與高原

肺水腫的易感性有關(guān)[50]

。 此前有研究表明硝苯地平是一種可以抑制自發(fā)性心律失常的藥物,而 SCN5A

在心律失常中起重要作用[51]

,這一關(guān)系也被 NGDTI 所預(yù)測。 此外,sorafenib 被批準(zhǔn)用于治療晚期腎細(xì)胞

癌,NGDTI 預(yù)測 sorafenib 與集落刺激因子受體(CSF1R)存在相互作用,已有研究也證實了 CSF1R 確實在

乳腺的發(fā)展和乳腺癌變中起著重要作用[52-53]

。 最后,NGDTI 預(yù)測 rivastigmine 和 CES1 的相互作用也在文

獻(xiàn)[54]中得到支持。 簡而言之,NGDTI 預(yù)測的藥物-靶標(biāo)相互關(guān)系有一部分是有文獻(xiàn)研究支持的,這進(jìn)一

步表明了 NGDTI 優(yōu)秀的預(yù)測能力。

3 結(jié)語

本文提出一種名為 NGDTI 的模型,用于集成來自不同異構(gòu)網(wǎng)絡(luò)的信息來預(yù)測新的藥物-靶標(biāo)相互作

用。 NGDTI 可以通過網(wǎng)絡(luò)擴(kuò)散過程從異構(gòu)網(wǎng)絡(luò)中提取低維隱藏特征信息,為了對節(jié)點特征進(jìn)行平滑和降

噪,該模型添加了圖卷積編碼來獲得更加有效的節(jié)點特征。 從實驗結(jié)果看,NGDTI 獲得了比其他基準(zhǔn)方

法更好的預(yù)測性能,而且 NGDTI 具有很強的魯棒性。 此外,NGDTI 是可擴(kuò)展的框架,其他有關(guān)藥物和靶標(biāo)

的更多信息也可以輕松地納入到當(dāng)前框架中。 因此,NGDTI 可以為加強藥物開發(fā)和藥物靶標(biāo)預(yù)測提供有

用的工具。 后續(xù)將進(jìn)一步優(yōu)化 NGDTI 模型,整合更多異構(gòu)信息,并改善模型的預(yù)測結(jié)果。 在本研究中,

NGDTI 模型僅用于預(yù)測未知的藥物-靶標(biāo)相互作用,但 NGDTI 模型也可以擴(kuò)展應(yīng)用于其他研究領(lǐng)域。

參 考 文 獻(xiàn)

[1] 徐國保, 陳媛曉, 王驥. 基于圖卷積網(wǎng)絡(luò)的藥物靶標(biāo)關(guān)聯(lián)預(yù)測算法[J]. 計算機(jī)應(yīng)用, 2021, 41(5): 1522-1526. DOI:

10.11772 / j.issn.1001-9081.2020081186.

[2] 余冬華, 郭茂祖, 劉曉燕, 等. 藥物靶標(biāo)作用關(guān)系預(yù)測結(jié)果評價及查詢驗證[J]. 計算機(jī)研究與發(fā)展, 2019, 56(9):

1881-1888. DOI: 10.7544 / issn1000-1239.2019.20180830.

[3] 宮莉, 康經(jīng)武. 基于化學(xué)蛋白質(zhì)組學(xué)的藥物靶標(biāo)鑒定[J]. 色譜, 2020, 38(8): 877-879.

[4] 李梢. 網(wǎng)絡(luò)靶標(biāo): 中藥方劑網(wǎng)絡(luò)藥理學(xué)研究的一個切入點[J]. 中國中藥雜志, 2011, 36(15): 2017-2020.

[5] PAUL S M, MYTELKA D S, DUNWIDDIE C T, et al. How to improve R&D productivity: the pharmaceutical industry’ s

grand challenge[J]. Nature Reviews Drug Discovery, 2010, 9(3): 203-214. DOI: 10.1038 / nrd3078.

[6] 彭利紅, 劉海燕, 任日麗, 等. 基于多標(biāo)記學(xué)習(xí)預(yù)測藥物-靶標(biāo)相互作用[ J]. 計算機(jī)工程與應(yīng)用, 2017, 53(15):

260-265.

[7] MORRIS G M, HUEY R, LINDSTROM W, et al. AutoDock4 and AutoDockTools4: automated docking with selective

receptor flexibility[J]. Journal of Computational Chemistry, 2009, 30(16): 2785-2791. DOI: 10.1002 / jcc.21256.

[8] KEISER M J, ROTH B L, ARMBRUSTER B N, et al. Relating protein pharmacology by ligand chemistry[ J]. Nature

Biotechnology, 2007, 25(2): 197-206. DOI: 10.1038 / nbt1284.

99

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

[9] KLIPP E, WADE R C, KUMMER U. Biochemical network-based drug-target prediction[J]. Current opinion in

biotechnology, 2010, 21(4): 511-516. DOI: 10.1016 / j.copbio.2010.05.004.

[10] YAMANISHI Y, ARAKI M, GUTTERIDGE A, et al. Prediction of drug-target interaction networks from the integration of

chemical and genomic spaces[J]. Bioinformatics, 2008, 24(13): i232-i240. DOI: 10.1093 / bioinformatics/ btn162.

[11] BLEAKLEY K, YAMANISHI Y. Supervised prediction of drug-target interactions using bipartite local models [ J ].

Bioinformatics, 2009, 25(18): 2397-2403. DOI: 10.1093 / bioinformatics/ btp433.

[12] MEI J P, KWOH C K, YANG P, et al. Drug-target interaction prediction by learning from local information and neighbors

[J]. Bioinformatics, 2013, 29(2): 238-245. DOI: 10.1093 / bioinformatics/ bts670.

[13] VAN LAARHOVEN T, NABUURS S B, MARCHIORI E. Gaussian interaction profile kernels for predicting drug-target

interaction[J]. Bioinformatics, 2011, 27(21): 3036-3043. DOI: 10.1093 / bioinformatics/ btr500.

[14] PAHIKKALA T, AIROLA A, PIETIL? S, et al. Toward more realistic drug-target interaction predictions[ J]. Briefings in

Bioinformatics, 2015, 16(2): 325-337. DOI: 10.1093 / bib / bbu010.

[15] HAO M, WANG Y L, BRYANT S H. Improved prediction of drug-target interactions using regularized least squares

integrating with kernel fusion technique[J]. Analytica Chimica Acta, 2016, 909: 41-50. DOI: 10.1016 / j.aca.2016.01.014.

[16] NASCIMENTO A C A, PRUDêNCIO R B C, COSTA I G. A multiple kernel learning algorithm for drug-target interaction

prediction[J]. BMC Bioinformatics, 2016, 17: 46. DOI: 10.1186 / s12859-016-0890-3.

[17] LIU Y, WU M, MIAO C Y, et al. Neighborhood regularized logistic matrix factorization for drug-target interaction prediction

[J]. PLOS Computational Biology, 2016, 12(2): e1004760. DOI: 10.1371 / journal.pcbi.1004760.

[18] G?NEN M. Predicting drug-target interactions from chemical and genomic kernels using Bayesian matrix factorization[ J].

Bioinformatics, 2012, 28(18): 2304-2310. DOI: 10.1093 / bioinformatics/ bts360.

[19] ZHENG X D, DING H, MAMITSUKA H, et al. Collaborative matrix factorization with multiple similarities for predicting

drug-target interactions[C]∥ Proceedings of the 19th ACM SIGKDD International Conference on Knowledge Discovery and

Data Mining. New York, NY: ACM Press, 2013: 1025-1033. DOI: 10.1145 / 2487575.2487670.

[20] HAO M, BRYANT S H, WANG Y L. Predicting drug-target interactions by dual-network integrated logistic matrix

factorization[J]. Scientific Reports, 2017, 7: 40376. DOI: 10.1038 / srep40376.

[21] MIZUTANI S, PAUWELS E, STOVEN V, et al. Relating drug-protein interaction network with drug side effects [ J].

Bioinformatics, 2012, 28(18): i522-i528. DOI: 10.1093 / bioinformatics/ bts383.

[22] LUO Y N, ZHAO X B, ZHOU J T, et al. A network integration approach for drug-target interaction prediction and

computational drug repositioning from heterogeneous information [ J]. Nature Communications, 2017, 8: 537. DOI: 10.

1038 / s41467-017-00680-8.

[23] NATARAJAN N, DHILLON I S. Inductive matrix completion for predicting gene-disease associations[ J]. Bioinformatics,

2014, 30(12): i60-i68. DOI: 10.1093 / bioinformatics/ btu269.

[24] CHEN Q F, LAI D H, LAN W, et al. ILDMSF: inferring associations between long non-coding RNA and disease based on

multi-similarity fusion[J]. IEEE/ ACM Transactions on Computational Biology and Bioinformatics, 2021, 18( 3): 1106-

1112. DOI: 10.1109 / TCBB.2019.2936476.

[25] LAN W, LAI D H, CHEN Q F, et al. LDICDL: LncRNA-disease association identification based on collaborative deep

learning[J/ OL]. IEEE/ ACM Transactions on Computational Biology and Bioinformatics, 2020 [ 2021-07-23]. https:∥

ieeexplore.ieee.org / document / 9246263. DOI: 10.1109 / TCBB.2020.3034910.

[26] 李琳, 梁永全, 劉廣明. 基于重啟隨機(jī)游走的圖自編碼器[ J]. 計算機(jī)應(yīng)用研究, 2021, 38(10): 3009-3013. DOI:

10.19734 / j.issn.1001-3695.2021.03.0083.

[27] 翟正利, 梁振明, 周煒, 等. 變分自編碼器模型綜述[ J]. 計算機(jī)工程與應(yīng)用, 2019, 55(3): 1-9. DOI: 10.3778 / j.

issn.1002-8331.1810-0284.

[28] 郭景峰, 董慧, 張庭瑋, 等. 主題關(guān)注網(wǎng)絡(luò)的表示學(xué)習(xí)[J]. 計算機(jī)應(yīng)用, 2020, 40(2): 441-447. DOI: 10.11772 / j.

issn.1001-9081.2019081529.

[29] 王佩琪, 高原, 劉振宇, 等. 深度卷積神經(jīng)網(wǎng)絡(luò)的數(shù)據(jù)表示方法分析與實踐[J]. 計算機(jī)研究與發(fā)展, 2017, 54(6):

100

http:∥xuebao.gxnu.edu.cn

1348-1356. DOI: 10.7544 / issn1000-1239.2017.20170098.

[30] 劉思, 劉海, 陳啟買, 等. 基于網(wǎng)絡(luò)表示學(xué)習(xí)與隨機(jī)游走的鏈路預(yù)測算法[ J]. 計算機(jī)應(yīng)用, 2017, 37(8): 2234-

2239. DOI: 10.11772 / j.issn.1001-9081.2017.08.2234.

[31] WAN F P, HONG L X, XIAO A, et al. NeoDTI: neural integration of neighbor information from a heterogeneous network for

discovering new drug-target interactions[J]. Bioinformatics, 2019, 35(1): 104-111. DOI: 10.1093/ bioinformatics/ bty543.

[32] GILMER J, SCHOENHOLZ S S, RILEY P F, et al. Neural message passing for quantum chemistry[C]∥ Proceedings of

the 34th International Conference on Machine Learning: Volume 70. Sydney: PMLR, 2017: 1263-1272.

[33] ROSSI R A, ZHOU R, AHMED N K. Deep inductive graph representation learning[ J]. IEEE Transactions on Knowledge

and Data Engineering, 2020, 32(3): 438-452. DOI: 10.1109 / TKDE.2018.2878247.

[34] WU Z H, PAN S R, CHEN F W, et al. A comprehensive survey on graph neural networks[ J]. IEEE Transactions on

Neural Networks and Learning Systems, 2021, 32(1): 4-24. DOI: 10.1109 / TNNLS.2020.2978386.

[35] 焦李成, 楊淑媛, 劉芳, 等. 神經(jīng)網(wǎng)絡(luò)七十年: 回顧與展望[J]. 計算機(jī)學(xué)報, 2016, 39(8): 1697-1716.

[36] SUN C, CAO Y K, WEI J M, et al. Autoencoder-based drug-target interaction prediction by preserving the consistency of

chemical properties and functions of drugs[J]. Bioinformatics, 2021, 37(20): 3618-3625. DOI: 10.1093 / bioinformatics/

btab384.

[37] XUAN P, ZHANG Y, CUI H, et al. Integrating multi-scale neighbouring topologies and cross-modal similarities for drugprotein interaction prediction[J]. Briefings in Bioinformatics, 2021, 22(5): bbab119. DOI: 10.1093 / bib / bbab119.

[38] GAO K Y, FOKOUE A, LUO H, et al. Interpretable drug target prediction using deep neural representation [ C]∥

Proceedings of the Twenty-Seventh International Joint Conference on Artificial Intelligence. Marina del Rey, CA: IJCAI,

2018: 3371-3377. DOI: 10.24963 / ijcai.2018 / 468.

[39] WANG S, CHO H H, ZHAI C X, et al. Exploiting ontology graph for predicting sparsely annotated gene function[ J].

Bioinformatics, 2015, 31(12): i357-i364. DOI: 10.1093 / bioinformatics/ btv260.

[40] KIPF T N, WELLING M. Variational graph auto-encoders[EB/ OL]. (2016-11-21)[2021-07-23]. https:∥arxiv.org / abs/

1611.07308.

[41] KNOX C, LAW V, JEWISON T, et al. DrugBank 3.0: a comprehensive resource for ‘Omics’ research on drugs[ J].

Nucleic Acids Research, 2011, 39(suppl_1): D1035-D1041. DOI: 10.1093 / nar/ gkq1126.

[42] PRASAD T S K, GOEL R, KANDASAMY K, et al. Human protein reference database: 2009 update[ J]. Nucleic Acids

Research, 2009, 37(suppl_1): D767-D772. DOI: 10.1093 / nar/ gkn892.

[43] DAVISA P, MURPHY C G, JOHNSON R, et al. The comparative toxicogenomics database: update 2013[ J]. Nucleic

Acids Research, 2013, 41(D1): D1104-D1114. DOI: 10.1093 / nar/ gks994.

[44] KUHN M, CAMPILLOS M, LETUNIC I, et al. A side effect resource to capture phenotypic effects of drugs[J]. Molecular

Systems Biology, 2010, 6(1): 343. DOI: 10.1038 / msb.2009.98.

[45] ROGERS D, HAHN M. Extended-connectivity fingerprints[ J]. Journal ofChemical Information and Modeling, 2010, 50

(5): 742-754. DOI: 10.1021 / ci100050t.

[46] SMITH T F, WATERMAN M S. Identification of common molecular subsequences[J]. Journal of Molecular Biology, 1981,

147(1): 195-197. DOI: 10.1016 / 0022-2836(81)90087-5.

[47] XIA Z, WU L Y, ZHOU X B, et al. Semi-supervised drug-protein interaction prediction from heterogeneous biological

spaces[J]. BMC Systems Biology, 2010, 4(Suppl 2): S6. DOI: 10.1186 / 1752-0509-4-S2-S6.

[48] WANG W H, YANG S, ZHANG X, et al. Drug repositioning by integrating target information through a heterogeneous

network model[J]. Bioinformatics, 2014, 30(20): 2923-2930. DOI: 10.1093 / bioinformatics/ btu403.

[49] DUVENAUD D, MACLAURIN D, AGUILERA-IPARRAGUIRRE J, et al. Convolutional networks on graphs for learning

molecular fingerprints[EB/ OL]. (2015-11-03)[2021-07-23]. https:∥arxiv.org / abs/ 1509.09292v2.

[50] YANG Y Z, DU H, LI Y H, et al. NR3C1 gene polymorphisms are associated with high-altitude pulmonary edema in Han

Chinese[J]. Journal of Physiological Anthropology, 2019, 38: 4. DOI: 10.1186 / s40101-019-0194-1.

[51] LIU G X, REMME C A, BOUKENS B J, et al. Overexpression of SCN5A in mouse heart mimics human syndrome of

101

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

enhanced atrioventricular nodal conduction[ J]. Heart Rhythm, 2015, 12( 5): 1036-1045. DOI: 10.1016 / j. hrthm.2015.

01.029.

[52] THOMAS G, GURUNG I S, KILLEEN M J, et al. Effects of L-type Ca

2+

channel antagonism on ventricular

arrhythmogenesis in murine hearts containing a modification in the Scn5a gene modelling human long QT syndrome 3[ J].

The Journal of Physiology, 2007, 578(1): 85-97. DOI: 10.1113 / jphysiol.2006.121921.

[53] ULLRICH K, WURSTER K D, LAMPRECHT B, et al. BAY 43-9006 / Sorafenib blocks CSF1R activity and induces

apoptosis in various classical Hodgkin lymphoma cell lines[ J]. British Journal of Haematology, 2011, 155(3): 398-402.

DOI: 10.1111 / j.1365-2141.2011.08685.x.

[54] TSURKAN L G, HATFIELD M J, EDWARDS C C, et al. Inhibition of human carboxylesterases hCE1 and hiCE by

cholinesterase inhibitors[J]. Chemico-biological interactions, 2013, 203(1): 226-230. DOI: 10.1016 / j.cbi.2012.10.018.

Study on Multi-information Integration for Drug Target Prediction

TAN Kai

1

, LI Yongjie

1

, PAN Haiming

1

, HUANG Kexin

2

, QIU Jie

2

, CHEN Qingfeng

1?

(1. School of Computer, Electronics and Information, Guangxi University, Nanning Guangxi 530004, China;

2. Guangxi Medical University, Nanning Guangxi 530021, China;

3. School of Computer Science and Engineering, Yulin Normal University, Yulin Guangxi 537000, China)

Abstract: Accurate determination of drug-target interactions is crucial in drug discovery process and

repositioning. Traditional methods for DTI prediction are either time-consuming ( simulation-based methods) or

heavily dependent on domain expertise ( similarity-based and feature-based methods). Existing computationbased methods using single data information or sparse data, always suffer from high false positive rates. Although

integrating multiple heterogeneous networks has been prevalent for drug target prediction, how to retain as much

structural information as possible is still a big challenge. This paper proposes a novel framework NGDTI, which

extracts relevant biological properties and association information from the network while maintaining the topology

information. Further, the graph neural network is applied to update the extracted feature information. The learned

topology-preserving representations of drugs and targets promote DTI prediction. Compared with the state-of-theart methods, NGDTI increases the AUPR value by nearly 0.01. The results demonstrate that NGDTI is promising

for drug development and repositioning.

Keywords: drug target association prediction; network embedding; network integration; matrix decomposition;

graph neural network

(責(zé)任編輯 黃 勇)

102

第 40 卷 第 2 期

2022 年 3 月

廣西師范大學(xué)學(xué)報(自然科學(xué)版)

Journal of Guangxi Normal University (Natural Science Edition)

Vol. 40 No. 2

Mar. 2022

DOI: 10.16088 / j.issn.1001-6600.2021052801 http: xuebao.gxnu.edu.cn

徐王軍, 曹進(jìn)德, 伍代勇, 等.一類具有遷移和 Allee 效應(yīng)的食餌-捕食者系統(tǒng)穩(wěn)定性[J]. 廣西師范大學(xué)學(xué)報(自然科學(xué)版), 2022, 40(2):

103-115. XU W J, CAO J D, WU D Y, et al. Stability of a prey-predator model with migration and Allee effects [ J]. Journal of Guangxi Normal

University(Natural Science Edition), 2022, 40(2): 103-115.

一類具有遷移和 Allee 效應(yīng)的食餌-捕食者

系統(tǒng)穩(wěn)定性

徐王軍1

, 曹進(jìn)德2

, 伍代勇1

, 申傳勝1?

(1. 安慶師范大學(xué) 數(shù)理學(xué)院, 安徽 安慶 246133; 2. 東南大學(xué) 數(shù)學(xué)學(xué)院, 江蘇 南京 211189)

摘 要: 研究一類食餌具有 Allee 效應(yīng)且捕食者具有人工控制遷移的食餌-捕食者系統(tǒng), 該系統(tǒng)具有平方根項的功能性反

應(yīng)函數(shù)。 首先通過定性分析, 證明解的有界性, 分析平衡點的存在性, 得到系統(tǒng)平衡點的局部穩(wěn)定性的充分條件。 接著

討論平衡點的 Hopf 分岔存在性, 并通過計算第一李雅普諾夫系數(shù), 研究平衡點 Hopf 分岔的穩(wěn)定性和方向。 最后通過數(shù)

值模擬驗證所得結(jié)論的正確性, 結(jié)果表明 Allee 效應(yīng)和人工控制遷移率對食餌種群和捕食者種群的生存與滅絕具有重要

意義。

關(guān)鍵詞: 食餌-捕食者模型; Allee 效應(yīng); 遷移率; Hopf 分岔; 穩(wěn)定性

中圖分類號: Q141; O175 文獻(xiàn)標(biāo)志碼: A 文章編號: 1001-6600(2022)02-0103-13

在生態(tài)系統(tǒng)中,任何種群都不是獨立存在的,而是與生物群落中其他種群相互聯(lián)系、互相共存。 同時,

種群之間一般存在著種內(nèi)競爭、種間競爭、捕食與被捕食及互惠等關(guān)系。 其中,關(guān)于捕食與被捕食系統(tǒng),人

們研究了多種功能性反應(yīng)函數(shù),如 Holling 型功能性反應(yīng)函數(shù)[1-6]

、Beddington-DeAngelis 型功能性反應(yīng)函

數(shù)[7]

,以及帶有平方根項的功能性反應(yīng)函數(shù)[8]

。 這里平方根項是基于食餌集體防御行為,即食餌種群常

常抱團(tuán)形成一個圓形防御區(qū),處于圓形區(qū)域最外圍的個體數(shù) R(x)將成為被捕食對象。 注意到 R(x)與食

餌種群圓形防御區(qū)域的周長 l 成比例,即 R(x)≈l。 假設(shè)該圓形防御區(qū)域的面積為 S, 則有 l = 2 πS 。一

般地,食餌種群密度 x, 即單位表面積內(nèi)食餌的個體數(shù)與其所占據(jù)的區(qū)域面積 S 成正比, 即 x∝S。 因此,

有 R(x)∝ x 。 然而目前關(guān)于平方根項功能性反應(yīng)函數(shù)的研究鮮有報道。 英國人口學(xué)家 Malthus(1766—

1834)根據(jù)人口統(tǒng)計資料,于 1798 年提出人口指數(shù)增長模型[9]

,該模型假設(shè)單位時間內(nèi)人口的增長量與

當(dāng)年的人口數(shù)量成正比。 后來生態(tài)學(xué)家們又提出了密度制約的 Logistic 增長模型[10]

。 實際上,當(dāng)種群數(shù)

量低于某一個閾值時,由于環(huán)境的作用,種群密度往往呈現(xiàn)負(fù)增長,例如自然界中一些密度很小的珍稀動

物(大熊貓、金絲猴、華南虎),如果人類不予保護(hù),順其自然,即使不加捕殺也會滅絕,這種現(xiàn)象被稱為

Allee 效應(yīng)。 近年來,Allee 效應(yīng)[11-16]受到廣泛關(guān)注,例如 Petrovskii 等[17] 研究了具有 Allee 效應(yīng)的生物入

侵機(jī)制的食餌-捕食者系統(tǒng),Rao 等[18]討論了具有 Allee 效應(yīng)的食餌-捕食者擴(kuò)散模型。

近年來,農(nóng)作物病蟲害綠色防控方法與措施受到高度關(guān)注,尤其是生物防治技術(shù)受到生態(tài)學(xué)家們的普

遍歡迎,如利用釋放禽類等害蟲天敵來防控病蟲害。 Barclay

[19] 提出利用捕食者釋放、棲息地管理和農(nóng)藥

釋放組合防治害蟲的模型。 Tang 等[20]研究了蟲害綜合治理的最佳時機(jī):天敵釋放和殺蟲劑施用率的模

擬。 成定平[21]建立鼠類與天敵系統(tǒng)漸近穩(wěn)定性的數(shù)學(xué)分析。 Chen 等[22]討論了具有捕食者遷移的時滯捕

食者-食餌模型,該模型中捕食者的遷移率依賴于有效食餌,當(dāng)捕食者捕獲率較低時,不能有效控制食餌的

數(shù)量,從而導(dǎo)致食餌急劇泛濫,故通過人為控制遷移捕食者達(dá)到控制捕食者的種群密度。 最近,文獻(xiàn)[23]

收稿日期: 2021-05-28 修回日期: 2021-06-28

基金項目: 國家自然科學(xué)基金(11975025); 國家自然科學(xué)基金委員會與英國皇家學(xué)會合作交流項目(12011530158)

通信作者: 申傳勝(1975—), 男, 安徽六安人, 安慶師范大學(xué)教授, 博士。 E-mail: csshen@mail.ustc.edu.cn

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

在三種群的食餌-捕食者模型中引入食餌的時變遷移速率,發(fā)現(xiàn)遷移誘導(dǎo)出現(xiàn) 2 個 Hopf 分叉和 2 個極限

環(huán),有趣的是,第 2 個 Hopf 分叉后的滅絕點出現(xiàn)類似于不連續(xù)的一級相變。 關(guān)于其他密度依賴的遷移對

食餌-捕食者模型穩(wěn)定性的影響研究參見文獻(xiàn)[24-27]。

文獻(xiàn)[22]建立的模型沒有考慮 Allee 效應(yīng)以及群防御行為。 實際上,當(dāng)食餌的數(shù)量較大時,由于被捕

食者群防御行為,捕食者捕獲到的食餌往往不能呈線性增長。 其次,捕獲率還會受到其他因素影響,如時

間、環(huán)境因素等。 為了更貼近實際生態(tài)意義,本文基于文獻(xiàn)[22]建立一類同時具有 Allee 效應(yīng)和人工控制

遷移率的 2 種群食餌-捕食者模型,分析了模型解的有界性、平衡點的存在性、局部穩(wěn)定性以及 Hopf 分岔。

此外,通過數(shù)值模擬驗證了以上動力學(xué)性質(zhì)的正確性。

1 模型的建立

經(jīng)典兩物種食餌-捕食者模型為

dx

dt

= h0(x)-R(x)y,

dy

dt

= -uy+eR(x)y。

ì

?

í

?

?

?

?

(1)

式中:x 和 y 分別表示食餌種群的密度和捕食者種群的密度;h0(x)是遵循 Logistic 增長的函數(shù),即 h0(x)=

rx(1-x / K),這里 r 和 K 分別表示食餌的內(nèi)稟自然增長率和環(huán)境容納量;u 和 e 分別表示捕食者的死亡率

和捕食者捕獲食餌后的轉(zhuǎn)化率;R(x)= α x ,α 表示捕食者的捕獲率,R(x)描述了在食餌種群中觀察到的

群體防御行為[8]

。

對模型(1)中第 1 個方程引入 Allee 效應(yīng)項, 得到

h(x)= rx 1-

x

K

( ) 1-

A+c

x+c

( ) ,

式中 A 和 c 分別表示 Allee 效應(yīng)閾值和輔助參數(shù),這里 0<A<K。 顯然,當(dāng)食餌種群密度小于 Allee 閾值時,

食餌呈現(xiàn)負(fù)增長,可見 Allee 效應(yīng)對食餌種群的增長具有很大的影響。

對模型(1)中第 2 個方程引入人工控制遷移函數(shù)

P(x,y)= ε(x-ρy),

式中:ε 表示捕食者的遷移率,ρ 表示單位時間內(nèi)捕食者捕獲食餌的消耗率。 因此,得到如下模型:

dx

dt

= h(x)-α x y,

dy

dt

= -uy+eα x y+P(x,y)。

ì

?

í

?

?

?

?

(2)

為了簡便起見,對模型(2)引入無量綱變換:

?t = rt,x?=

x

K

,y?=

y

K

,α?=

α K

r

,c?=

c

K

,?e = e,u?=

u

r

,ε?=

ε

r

,ρ?= ρ,θ =

A

K

。

式中 0<θ<1。 于是,去掉“短橫”,得到系統(tǒng)

dx

dt

=

x

x+c

(1-x)(x-θ)-α x y,

dy

dt

= -uy+eα x y+ε(x-ρy)。

ì

?

í

?

?

?

?

(3)

考慮到生物學(xué)的實際意義,系統(tǒng)(3)的可行域為 E ={(x,y) x≥0,y≥0} 。 系統(tǒng)中的所有參數(shù)都是正常數(shù)。

104

http:∥xuebao.gxnu.edu.cn

2 系統(tǒng)解的有界性

定理 1 系統(tǒng)(3)初值滿足(x0 ,y0 )∈E 的解非負(fù)且有界。

證明 對于任何初始值(x0 ,y0 )∈E,并且對任意 t>0,由系統(tǒng)(3)可知:

x(t) = x(0)exp ∫

t

0

(1 - x(s))(x(s) - θ)

x(s) + c

-

αy(s)

x(s)

é

?

ê

ê

ù

?

ú

ú

ds ( ) ≥ 0,

y(t) = y(0)exp ∫

t

0

- u + eα x(s) +

εx(s)

y(s)

- ερ

é

?

ê

ê

ù

?

ú

ú ( ds) ≥ 0。

ì

?

í

?

?

?

?

因此,對任意 t>0,都有 x(t)≥0,y(t)≥0。 系統(tǒng)(3)的解非負(fù)。

當(dāng) x∈(0,θ)∪(1,+∞ )時,恒有

dx

dt

<0,即有l(wèi)im

t→∞

x(t)= 0,而當(dāng) x∈(θ,1)時,有

dx

dt

≤

x(1-x)(x-θ)

x+c

≤(1-x)(x-θ)。 (4)

考慮方程

dz

dt

= (1-z)(z-θ),

由方程解得lim sup

t→∞

z(t)≤1。 通過比較定理,得lim sup

t→∞

x(t)≤1。

定義一個關(guān)于解的和函數(shù) L(t)= x(t)+

y(t)

e

,其中 e>0,那么

dL

dt

=

x

x+c

(1-x)(x-θ)-

u

e

y+

ε

e

(x-ρy)≤(1-x)(x-θ)-

u+ερ

e

( ) y+

ε

e

x =

-x

2+(θ+1)x-θ+

ε

e

x-(u+ερ)(L-x)= -x

2+(θ+

ε

e

+u+ερ+1)x-θ-(u+ερ)L。 (5)

記

g(x)= -x

2+(θ+

ε

e

+u+ερ+1)x-θ,

則有

dg

dx

= -2x+θ+

ε

e

+u+ερ+1。 (6)

令

dg

dx

= 0,解得

xa

=

(θ+u+ρε+1)e+ε

2e

,

且

d

2

g

dx

2

= -2<0。

(ⅰ)當(dāng) xa<1 時,g(x)在 x = xa 取得最大值。

g(xa )=

(θe+ue+ερe+e+ε)

2

4e

2

- θ,

即有

dL

dt

≤g(xa )-(u+ερ)L,

從而

105

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

lim sup

t→∞

L(t)≤

g(xa )

u+ερ

。 (7)

(ⅱ)當(dāng) xa≥1 時,g(x)在 x = 1 取得最大值

g(1)=

ε+ue+ρεe

e

,

即有

dL

dt

≤g(1)-(u+ερ)L,

從而

lim sup

t→∞

L(t)≤

g(1)

u+ερ

。 (8)

綜合式(7)和式(8),可得

lim sup

t→∞

x(t)+

y(t)

e

( ) ≤

g(xa )

u+ερ

, xa<1,

g(1)

u+ερ

, xa≥1。

ì

?

í

?

?

?

?

(9)

于是結(jié)論成立。 證畢。

3 平衡點分析

3.1 平衡點的存在性

確定系統(tǒng)(3)的非負(fù)平衡點,令(x

?

,y

?

)為其平衡點,則:

x

?

(1-x

?

)(x

? -θ)-α(x

? +c) x

?

y

? = 0,

-uy

? +eα x

?

y

? +ε(x

? -ρy

?

)= 0。 { (10)

1)系統(tǒng)(3)恒存在滅絕平衡點 E0(0,0)。

2)當(dāng) ε = 0,系統(tǒng)(3)有 2 個軸向平衡點 E11(1,0)和 E12( θ,0)以及正平衡點 E1( x

?

1 ,y

?

1 ),記 x

?

1

= η,

y

?

1

= ξ。 當(dāng) η∈(θ,1),由方程組(10)可得

η =

u

eα

( )

2

,ξ =

u(1-η)(η-θ)

eα

2

(η+c)

。 (11)

3)當(dāng) ε≠0,系統(tǒng)(3)有正平衡點 E2(x

?

,y

?

),即:

y

? =

x

?

(1-x

?

)(x

? -θ)

α(x

? +c)

,

(eα x

? -u-ρε)y

? +εx

? = 0。

ì

?

í

?

?

?

?

(12)

從而

x

?

(1-x

?

)(x

? -θ)

α(x

? +c)

=

-εx

?

eα x

? -u-ρε

。 (13)

記 δ = x

?

,由式(13)得

eαδ

5-(u+ρε)δ

4-α(e(θ+1)+ε)δ

3+(u+ρε)(θ+1)δ

2-α(εc-eθ)δ-(u+ρε)θ = 0。 (14)

再令

M(δ)= a0

δ

5+a1

δ

4+a2

δ

3+a3

δ

2+a4

δ+a5 , (15)

式中:

a0

= eα>0,a1

= -(u+ρε)<0,a2

= -α(e(θ+1)+ε)<0,a3

= (u+ρε)(θ+1)>0,

106

http:∥xuebao.gxnu.edu.cn

a4

= -α(εc-eθ),a5

= -(u+ρε)θ<0。

若 M(δ)= 0 的 5 個根 δ 是實數(shù),由笛卡爾符號法則可知,則 M(δ)= 0 正根的個數(shù)等于它的系數(shù)序列的變

號數(shù),即 M(δ)= 0 有 5 個正實數(shù)根,但是它的系數(shù)序列的變號數(shù)為 3,這與結(jié)論相矛盾。 因此,M(δ)= 0 至

多存在 3 個正實數(shù)根滿足 δ = x

?

>0 且 x

? ∈( θ,1),即系統(tǒng)(3)至多存在 3 個正平衡點 E2( x

?

,y

?

),同

時,因為 a0>0 且 a5<0,所以 M(δ)的系數(shù)序列的變號數(shù)至少存在一個變號數(shù),即系統(tǒng)(3)至少存在一個正

平衡點 E2(x

?

,y

?

)。

3.2 穩(wěn)定性分析

對系統(tǒng)(3)的所有非負(fù)平衡點進(jìn)行局部穩(wěn)定性分析,還要注意,對任意 x∈(0,θ],恒有 dx / dt<0,隨著

時間的不斷演化,食餌最終會趨于滅絕(x = 0),由于食餌滅絕,捕食者的增長率呈現(xiàn)負(fù)增長,即有dy / dt<0,

捕食者也趨于滅絕(y = 0),因此,滅絕平衡點 E0 是穩(wěn)定的。 至于系統(tǒng)(3)的其他平衡點的局部穩(wěn)定性,可

以通過計算其 Jacobian 矩陣的特征值,并運用 Lyapunov 第一種方法和 Routh-Hurwitz 準(zhǔn)則確定各平衡點局

部漸近穩(wěn)定的充分條件。

令

dx

dt

=

x

x+c

(1-x)(x-θ)-α x y:=Q(x,y),

dy

dt

= -uy+eα x y+ε(x-ρy):=W(x,y)。

ì

?

í

?

?

?

?

(16)

于是,有

?Q

?x

=

(2θ-3c+1)x

2-2x

3+2c(θ+1)x-cθ

(x+c)

2

-

αy

2 x

,

?Q

?y

= -α x ,

?W

?x

=

eαy

2 x

+ε,

?W

?y

= -u+eα x -ερ。

ì

?

í

?

??

?

?

(17)

(ⅰ)系統(tǒng)(3)在軸向平衡點 E11處的 Jacobian 矩陣為

J11

=

θ-1

c+1

-α

0 eα-u

( ) 。 (18)

因此,軸向平衡點 E11對應(yīng)的特征方程為

λθ-1

c+1

( ) (λ-eα+u)= 0, (19)

它的特征值為 λ1

=

θ-1

c+1

,λ2

= eα-u。 顯然,λ1<0。 若 u>eα,即 λ2<0,則 E11局部漸近穩(wěn)定,此時 E11為穩(wěn)定結(jié)

點;若 u<eα,即 λ2>0,則 E11為鞍點。

(ⅱ)系統(tǒng)(3)在軸向平衡點 E12處的 Jacobian 矩陣為

J12

=

θ(1-θ)

θ+c

-α θ

0 eα θ -u

( )

。 (20)

因此,軸向平衡點 E12對應(yīng)的特征方程為

λθ(1-θ)

θ+c

( ) (λ-eα θ +u)= 0, (21)

它的特征值為 λ1

=

θ(1-θ)

c+1

,λ2

= eα θ -u。 顯然,λ1>0。 若 θ<

u

eα

( )

2

,即 λ2<0,則 E12是鞍點;若 θ>

u

eα

( )

2

,即

λ2>0,則 E12是不穩(wěn)定結(jié)點。

(ⅲ)系統(tǒng)(3)在正平衡點 E1 處的 Jacobian 矩陣為

107

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

J1

=

a11 a12

a21 0 ( ) 。 (22)

式中:

a11

=

(θ-3c+1)η

2-2η

3+2c(θ+1)η-cθ

(η+c)

2

-

αξ

2 η

,a12

= -

u

e

,a21

=

e

2

α

2

ξ

2u

。 (23)

因此,E1 對應(yīng)的特征方程為

λ

2-a11λ-a12 a21

= 0, (24)

設(shè)其特征值為 λi(i = 1,2),若

θ<

3η

3+(5c-1)η

2-3cη

η

2+(3c+1)η-c

, (25)

則有

λ1

+λ2

= a11<0,

λ1λ2

= -a12 a21>0。 { (26)

從而它的特征值 λi<0(i = 1,2),于是可得 E1 是局部漸近穩(wěn)定的。

(ⅳ)系統(tǒng)(3)在正平衡點 E2 的 Jacobian 矩陣為

J2

=

b11

b12

b21

b22

( ) , (27)

式中:

b11

=

(2θ-3c+1)x

? 2-2x

? 3+2c(θ+1)x

? -cθ

(x

? +c)

2

-

αy

?

2 x

?

,

b12

= -α x

?

,

b21

=

e1αy

?

2 x

?

+ε,

b22

= eα x

? -u-ερ。

ì

?

í

?

?

?

?

?

?

?

?

??

(28)

于是,它所對應(yīng)的特征方程為

λ

2+?1λ+?2

= 0, (29)

式中:

?1

= -(b11

+b22 ),

?2

= b11

b22

-b12

b21 。 { (30)

根據(jù) Routh-Hurwitz 判據(jù)準(zhǔn)則,有 Δ1

=?1>0,Δ2

=

?1 1

0 ?2

=?1?2>0。 如果 ?1>0 和 ?2>0 成立,那么 E2 是

局部漸近穩(wěn)定的。

3.3 Hopf 分岔

以 Allee 效應(yīng)閾值 θ 為分岔參數(shù),討論在正平衡點 E1 處發(fā)生 Hopf 分岔的條件。 設(shè) λ( θ) = λR(θ) +

iλI(θ)為特征方程(24)的特征值,代入特征方程(24),然后將其實部和虛部分離,得到:

λ

2

R

-λ

2

I

-a11λR

-a12 a21

= 0,

2λRλI

-a11λI

= 0。 { (31)

當(dāng)特征值出現(xiàn)實部為零的復(fù)共軛時,正平衡點 E1 發(fā)生 Hopf 分岔失去其穩(wěn)定性。 此時,Hopf 分岔點定義為

θH ,則有 λR(θH )= 0,將其代入式(31),得到:

-λ

2

I

-a12 a21

= 0,

a11λI

= 0(λI≠0)。 { (32)

108

http:∥xuebao.gxnu.edu.cn

從方程(32)解得 a11(θH )= 0,而 λI(θH )= -a12(θH )a21(θH ) >0,即有

det(J1 ) θ=θH

= -a12(θH )a21(θH )>0。

聯(lián)合式(23),得到

θH

=

αξ η (η+c)

2+4η

4-4cη

2+2(3c-1)η

3

2η(η

2+2cη-c)

。 (33)

接著對方程(31)兩邊關(guān)于 θ 求導(dǎo),注意 λR(θ)= 0,于是有

-a11

dλR

dθ

-2λI

dλI

dθ

=

d(a12 a21 )

dθ

,

2λI

dλR

dθ

-a11

dλI

dθ

=λI

d(a11 )

dθ

。

ì

?

í

?

??

?

?

(34)

解得

dλR

dθ θ=θH

=

2λ

2

I da11

dθ

-

a11 d(a12 a21 )

dθ

a

2

11

+4λ

2

I θ=θH

≠0 ?

2λ

2

I da11

dθ θ=θH

≠

a11 d(a12 a21 )

dθ θ=θH

。 (35)

定理 2 如果det(J1 ) θ=θH

= -a12(θH )a21(θH ) >0 且

2λ

2

I da11

dθ θ=θH

≠

a11 d(a12 a21 )

dθ θ=θH

成立,那么當(dāng) θ<θH

時,系統(tǒng)(3)的正平衡點 E1 局部漸近穩(wěn)定;并且 E1 在 θ = θH 處發(fā)生 Hopf 分岔。

定理 3 定義 l

1 為

l

1 := 3sxxx

s

2

y

-(sxy

+2sxxy

syη)?sxx

sy

+3sxxx

s

2

yη- -syhx

sy

sxxy

ξ」, (36)

式中:

sy

= a12

= -α η ,sxx

=

-2η

3-6cη

2-6c

2

η+2c

2

θ+2c

2+2cθ

(η+c)

3

+

αξ

4η η

,sxy

=

-α

2 η

,

sxxx

=

6c

2

η-6c

2

θ-6c

2-6cθ

(η+c)

4

-

3αξ

8η

2

η

,sxxy

=

α

4η η

,hx

= a21

=

e

2

α

2

ξ

2u

。

ì

?

í

?

??

?

??

(37)

若 l

1>0,則正平衡點 E1 的 Hopf 分岔是亞臨界的;若 l

1<0,則正平衡點 E1 的 Hopf 分岔是超臨界的。

證明 為了研究 Hopf 分岔的穩(wěn)定性和方向,下面計算第一李雅普諾夫系數(shù),令 u = x-η,v = y-ξ,那么

系統(tǒng)(3)轉(zhuǎn)變?yōu)?/p>

du

dt

=

u+η

u+η+c

(1-u-η)(u+η-θ)-α u+η (ξ+v):= s(u,v),

dv

dt

= -u(ξ+v)+eα u+η (ξ+v)+ε(u+η-ρ(ξ+v)):= h(u,v)。

ì

?

í

?

?

?

?

(38)

現(xiàn)在對方程(38)在(u,v)= (0,0)考慮其三階泰勒展開式,得到:

du

dt

= a11 u+a12

v+s1(u,v),

dv

dt

= a21 u+a22

v+h1(u,v)。

ì

?

í

?

?

?

?

(39)

這里 s1(u,v)和 h1(u,v)是關(guān)于 u 和 v 的高階項,則有

s1(u,v)= suu u

2+suvuv+svv

v

2+suuu u

3+suuvu

2

v+suvvuv

2+svvv

v

3

,

h1(u,v)= huu u

2+huvuv+hvv

v

2+huuu u

3+huuvu

2

v+huvvuv

2+hvvv

v

3

。 { (40)

式中:

109

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

su

= a11

=

(θ-3c+1)η

2-2η

3+2c(θ+1)η-cθ

(η+c)

2

-

αξ

2 η

,sv

= a12

= -α η ,

suu

=

-2η

3-6cη

2-6c

2

η+2c

2

θ+2c

2+2cθ

(η+c)

3

+

αξ

4η η

,suv

=

-α

2 η

,

svv

= 0,suuu

=

6c

2

η-6c

2

θ-6c

2-6cθ

(η+c)

4

-

3αξ

8η

2

η

,suuv

=

α

4η η

,suvv

= 0,svvv

= 0,

hu

= a21

=

e

2

α

2

ξ

2u

,hv

= a22

= 0,huu

= -

eαξ

4η η

,huv

=

eα

2 η

,

hvv

= 0,huuu

=

3eαξ

8η

2

η

,huuv

=

-eα

4η η

,huvv

= 0,hvvv

= 0。

ì

?

í

?

?

?

?

?

?

??

?

?

?

?

?

?

?

(41)

方程(40)中所有偏導(dǎo)數(shù)都是在分岔點計算的,即(u,v)= (0,0)。 因此系統(tǒng)(39)可以寫成

P? = JE1

P+D(P), (42)

式中:

P = (u,v)

T

,D = (s1(u,v),h1(u,v))

T

,s1(u,v)= suu u

2+suvuv+suuu u

3+suuvu

2

v,

h1(u,v)= huu u

2+huvuv+huuu u

3+huuvu

2

v。

當(dāng) su

= a11

= 0 時,正平衡點 E1 在 Hopf 分岔點的特征值是純虛數(shù),即 λ = i -a12 a21

= i -svhu 。 再令該

特征值對應(yīng)的特征向量為 u?= ( u1 ,u2 )

T

。 記正平衡點 E1 對應(yīng)的 Jacobian 矩陣為 A =

0 sv

hu 0 ( ) ,則有 Au?=

λu?,解得 u?= (sv,i -svhu )

T

。

定義

F= (Re(u?)-Im(u?))=

sv 0

0 - -svhu

é

?

ê

ê

ù

?

ú

ú

,

令 P =FZ,其中 Z = (z1 ,z2 )

T

,得

u

v

( ) =

sv 0

0 - -svhu

é

?

ê

ê

ù

?

ú

ú

z1

z2

é

?

ê

ê

ù

?

ú

ú

,

從而

u = sv

z1 ,

v = - -svhu

z2 。 {

因此,通過變換得到如下系統(tǒng)

Z? = (F

-1

J1F)Z+F

-1D(FZ),

即

z?

1

z?

2

é

?

ê

ê

ù

?

ú

ú

=

0 - -svhu

-svhu 0

é

?

ê

ê

ê

ù

?

ú

ú

ú

z1

z2

é

?

ê

ê

ù

?

ú

ú

+

D1(z1 ,z2 )

D2(z1 ,z2 )

é

?

ê

ê

ù

?

ú

ú

, (43)

式中:

D1(z1 ,z2 )=

1

sv

[suu (sv

z1 )

2- -svhu

suv

sv

z1

z2

+suuu(sv

z1 )

3-suuv(sv

z1 )

2 -svhu

z2 ],

D2(z1 ,z2 )= -

1

-svhu

huu(sv

z1 )

2- -svhu huv

sv

z1

z2

+huuu(sv

z1 )

3- -svhu huuv(sv

z1 )

2

z2

[ ] 。

ì

?

í

?

??

?

?

(44)

Hopf 分岔的方向由第一李雅普諾夫系數(shù)的符號決定,于是有

l

1 :=

1

16

?

3D1

?z

3

1

+

?

3D2

?z

2

1 ?z2

é

?

ê

ê

ù

?

ú

ú

+

1

16 -svhu

?

2D1

?z1 ?z2

?

2D1

?z

2

1

-

?

2D2

?z1 ?z2

?

2D2

?z

2

1

-

?

2D1

?z

2

1

?

2D2

?z

2

1

é

?

ê

ê

ù

?

ú

ú

, (45)

110

http:∥xuebao.gxnu.edu.cn

式中:

?

3D1

?z

3

1

= 6suuu

s

2

v ,

?

3D1

?z1 ?z

2

2

= 0,

?

3D2

?z

2

1 ?z2

= 2huuv

s

2

v ,

?

3D2

?z

3

2

= 0,

?

2D1

?z1 ?z2

= - -svhu

suv

-2 -svhu

suuvu,

?

2D2

?z1 ?z2

= svhuv

+2huuv

svu,

?

2D1

?z

2

1

= 2suu

sv

+6suuu

svu+2vsv

suuv,

?

2D1

?z

2

2

= 0,

?

2D2

?z

2

1

= -

2huu

s

2

v

-svhu

-

6huuu

s

2

v u

-svhu

-

2huuv

s

2

v

v

-svhu

,

?

2D2

?z

2

2

= 0。

ì

?

í

?

?

?

?

?

?

?

?

?

?

?

?

(46)

通過簡化得到 l

1 的表達(dá)式為

l

1 := 3suuu

s

2

v

-(suv

+2suuv

svη)[suu

sv

+3suuu

s

2

v η- -svhu

sv

suuv

ξ]。 (47)

若 l

1>0,則正平衡點 E1 的 Hopf 分岔是亞臨界的;若 l

1<0,則正平衡點 E1 的 Hopf 分岔是超臨界的。

4 數(shù)值模擬

下面,通過選擇不同的 Allee 效應(yīng)閾值 θ 和人工控制遷移參數(shù) ε,采用數(shù)值模擬來驗證各個平衡點存

在條件及其穩(wěn)定性。

在圖 1 中,固定參數(shù) r = 0.9,c = 1.5。 紅色曲線和藍(lán)色曲線指食餌種群的環(huán)境容納量為 15 和 25。 當(dāng)食

餌種群不受 Allee 效應(yīng)的影響(θ = 0)且環(huán)境容納量較大時,h0(x)的最大值也較大(見圖 1(a));當(dāng)食餌種

群受到 Allee 效應(yīng)的影響時,食餌種群在 0<x<5 呈現(xiàn)負(fù)增長,在 5<x<15 呈現(xiàn)正增長。 與沒有 Allee 效應(yīng)相

比,當(dāng)食餌種群的環(huán)境容納量相同時,食餌種群受到 Allee 效應(yīng)影響后的增長率會減小(見圖 1( b))。 因

此,Allee 效應(yīng)對種群的存活和滅絕具有非常重要的作用。

(a)沒有 Allee 效應(yīng)時 (b)具有 Allee 效應(yīng)時

圖 1 參數(shù) r = 0.9,c = 1.5 時食餌種群的增長率 h0(x)和 h(x)

Fig. 1 Growth rate h0(x) and h(x) of the prey population when r = 0.9,c = 1.5

在圖 2 中,固定參數(shù) α= 0.24,u = 0.05,e = 0.32,c = 0.01,根據(jù)方程(11)可得 η = 0.423 85。 當(dāng) θ = 0.04 滿

足條件(25)時,食餌種群和捕食者種群的相圖曲線趨于一個穩(wěn)定點 E1(0.423 85,1.382 7),此時,E1 是局

部漸近穩(wěn)定的(如圖 2(a)); 當(dāng) θ = 0.075不滿足條件(25)時,食餌種群和捕食者種群的相圖曲線趨于一個

穩(wěn)定的極限環(huán)(如圖 2(b))。 因此,從圖 2(a)和圖 2(b)可以看出,當(dāng) Allee 閾值 θ 較小時,則系統(tǒng)趨于正

平衡點的可能性較大。

111

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

圖 2 系統(tǒng)(3)在 ε= 0 時的相圖

Fig. 2 Phase diagram of system (3) at ε= 0

在圖 3 中,假設(shè)系統(tǒng)(3)沒有人為控制遷移。 固定參數(shù) α = 0.24,u = 0.05,e = 0.32,c = 0.01, 根據(jù)方程

(11)解得 η = 0.423 85,再代入方程(33),得到 θH≈0.06,由方程(36)解得 l

1

= -0.031 8<0。 圖 3(a)繪制了

參數(shù)為 θ 的食餌種群分岔圖,其中插圖是臨界值(θc

= 0.075 6)附近的放大。 由圖 3( a)可知,隨著 θ 的增

大,系統(tǒng)經(jīng)歷從穩(wěn)定狀態(tài)到 Hopf 分岔(HB),以及極限環(huán)的不穩(wěn)定狀態(tài)和急速下降,直到食餌種群的密度

低于 θc 時,食餌進(jìn)入滅絕狀態(tài)。 可見,模擬得到的分岔點與理論計算閾值(θH≈0.06)一致。 若 θ>θc,則滅

絕平衡點 E0 是局部漸近穩(wěn)定的。 系統(tǒng)的正平衡點 E1 在 θ = θH 發(fā)生 Hopf 分岔,當(dāng) 0<θ<θH 時,正平衡點 E1

是局部漸近穩(wěn)定的。 圖 3(b)繪制了參數(shù)為 θ 的捕食者種群分岔圖。 由圖 3(b)可知,捕食者種群經(jīng)歷的

狀態(tài)和食餌種群經(jīng)歷的狀態(tài)相同,但是捕食者種群在穩(wěn)定狀態(tài)(θ∈(0,0.06))的密度是單調(diào)遞減的,而食

餌種群的密度一直保持穩(wěn)定不變,即 x = 0.423 85。 因此,食餌種群受到 Allee 效應(yīng)的影響越小,對保護(hù)一

些密度較小的珍稀動物的種群多樣性是有益的。

圖 3 系統(tǒng)(3)中參數(shù)為 θ 且 ε= 0 的食餌和捕食者分岔圖

Fig. 3 Bifurcation diagram of prey population and predator population in system (3) with parameter θ and ε= 0

在圖 4 中,假設(shè)系統(tǒng)(3)存在人為控制遷移。 固定參數(shù) α= 0.24,u = 0.05,e = 0.32,c = 0.01。 圖 4(a)繪

制了參數(shù)為 θ 的食餌種群分岔圖,其中圖 4(a)插圖分別是分岔點(eθh

= 0.125)附近的放大和臨界值(θ1

=

0.135)附近的放大。 由圖 4(a)可知,隨著 θ 的增大,系統(tǒng)經(jīng)歷從穩(wěn)定狀態(tài)到 Hopf 分岔(HB),以及極限環(huán)

的不穩(wěn)定狀態(tài)和急速下降,直到食餌種群的密度低于 θ1 時,食餌進(jìn)入滅絕狀態(tài)。 若 θ>θ1 ,則滅絕平衡點

E0 是局部漸近穩(wěn)定的。 系統(tǒng)的正平衡點 E2 在 θ = θh 發(fā)生 Hopf 分岔,當(dāng) 0<θ<θh 時,隨著 θ 的增大,食餌種

群的密度在不斷遞減,此時 E2 在穩(wěn)定狀態(tài)是局部漸近穩(wěn)定的。 圖 4(b)繪制了參數(shù)為 θ 的捕食者種群分

岔圖, 其中圖 4(b)插圖是臨界值(θ2

= 0.156)附近的放大。 由圖 4(b)可知,系統(tǒng)經(jīng)歷穩(wěn)定狀態(tài)和急劇下

降,直到捕食者種群滅絕,并且捕食者種群在穩(wěn)定狀態(tài)的密度也在逐漸減小。 因此, 生物控制遷移對維護(hù)

生態(tài)系統(tǒng)的動態(tài)平衡有著非常重要的意義。

112

http:∥xuebao.gxnu.edu.cn

圖 4 系統(tǒng)(3)中參數(shù)為 θ 且 ε= 0.2 的食餌和捕食者分岔

Fig. 4 Bifurcation diagram of prey population and predator population in system (3) with parameter θ and ε= 0.2

在圖 5 中,固定參數(shù) α= 0.3,u = 0.2,e = 1,c = 0.001,ρ = 0.5,并解出 θ = 0.103 不滿足條件(25)。 當(dāng)系統(tǒng)

(3)沒有人為控制遷移時,食餌種群和捕食者種群的相圖曲線趨于一個穩(wěn)定的極限環(huán)(如圖 5( a));當(dāng)系

統(tǒng)(3)存在人為控制遷移時,食餌種群和捕食者種群的相圖曲線趨于一個穩(wěn)定點 E2

= (0.46,0.95)(如圖 5

(b))。 因此,人為控制遷移有利于生態(tài)系統(tǒng)的穩(wěn)定性。

圖 5 系統(tǒng)(3)在 θ= 0.103 時的相圖

Fig.5 Phase diagram of system (3) at θ= 0.103

在圖 6 中,固定參數(shù) α= 0.24,u = 0.05,c = 0.01,ρ = 0.5,θ = 0.1。 繪制參數(shù)為 e 的捕食者種群分岔圖,其

中圖 6(a)插圖分別是第 1 個分岔點(eH1

= 0.16)附近的放大、第 2 個分岔點(eH2

= 0.306)附近的放大和臨

界值(ec

= 0.312)附近的放大,圖 6(b)插圖是臨界值(e0

= 0.443)附近的放大。 當(dāng)系統(tǒng)(3)沒有人為控制遷

移時,由圖 6(a)可知,e 從 0 增加到 0.4,系統(tǒng)(3)經(jīng)歷從滅絕狀態(tài)到 Hopf 分岔(HB),以及極限環(huán)的不穩(wěn)定

狀態(tài)和急劇上升,再到經(jīng)歷微小的 Hopf 分岔(HB)和急劇下降,直到捕食者種群趨于滅絕,即 y = 0。 若捕

食者轉(zhuǎn)化率高于臨界值 ec,則滅絕平衡點 E0 是局部漸近穩(wěn)定的;若捕食者轉(zhuǎn)化率低于分岔點 eH1

,則軸向

平衡點 E11是局部漸近穩(wěn)定的;若捕食者轉(zhuǎn)化率 e∈(0.24,0.306),則正平衡點 E1 是局部漸近穩(wěn)定的,如圖

6(a)所示。 其次,如果系統(tǒng)(3)受到人為控制遷移時,系統(tǒng)經(jīng)歷穩(wěn)定狀態(tài)和急劇下降,直到捕食者種群滅

絕。 若捕食者轉(zhuǎn)化率 e∈(0,0.443),則捕食者種群的密度隨著 e 的增大而逐漸增加,并且正平衡點 E2 是

局部漸近穩(wěn)定的,如圖 6(b)所示。 因此,適當(dāng)?shù)娜藶榭刂七w移不僅提高了捕食者的轉(zhuǎn)化率,而且有利于持

久維護(hù)物種的多樣性。

113

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

圖 6 系統(tǒng)(3)中參數(shù)為 e 的捕食者種群分岔

Fig. 6 Bifurcation diagram of predator population with parameterein system (3)

5 結(jié)論

數(shù)值模擬結(jié)果和分析表明,食餌-捕食者系統(tǒng)在 Allee 效應(yīng)和人工控制遷移作用下,出現(xiàn)了穩(wěn)定的極限

環(huán)、穩(wěn)定點以及 Hopf 分岔等復(fù)雜的動力學(xué)行為。 食餌種群受到 Allee 效應(yīng)的影響越小越有利于維持其物

種多樣性,從而對保護(hù)一些密度較小的珍稀動物種群是有益的,例如北極熊、金絲猴等;生態(tài)系統(tǒng)中有時需

要適當(dāng)?shù)娜藶榭刂七w移來維護(hù)系統(tǒng)的穩(wěn)定性,如自然界中一些鳥類在捕獲害蟲時, 因為鳥類自身捕獲率

較低,不能抑制害蟲數(shù)量的增長,所以需要人為控制遷移一些鳥類種群來捕獲害蟲,從而使得鳥類和害蟲

構(gòu)成的生態(tài)系統(tǒng)保持動態(tài)平衡。

參 考 文 獻(xiàn)

[1] MA Z H, WANG S F, WANG T T, et al. Stability analysis of prey-predator system with Holling type functional response

and prey refuge[J]. Advances in Difference Equations, 2017, 2017: 243.

[2] LIU B, ZHANG Y J, CHEN L S. Dynamic complexities of a Holling I predator-prey model concerning periodic biological

and chemical control[J]. Chaos, Solitons & Fractals, 2004, 22(1): 123-134.

[3] ZHANG Y J, XU Z L, LIU B, et al. Dynamic analysis of a Holling I predator-prey system with mutual interference

concerning pest control[J]. Journal of Biological Systems, 2005, 13(1): 45-58.

[4] ZU J, MIMURA M. The impact of Allee effect on a predator-prey system with Holling type II functional response [ J].

Applied Mathematics and Computation, 2010, 217(7): 3542-3556.

[5] SUN G, SARWARDI S, PAL P J, et al. The spatial patterns through diffusion driven instability in modified Leslie-Gower

and Holling-type II predator-prey model[J]. Journal of Biological Systems, 2010, 18(3): 593-603.

[6] YI F Q, WEI J J, SHI J P. Bifurcation and spatiotemporal patterns in a homogeneous diffusive predator-prey system[ J].

Journal of Differential Equations, 2009, 246(5): 1944-1977.

[7] FAN M, KUANG Y. Dynamics of a nonautonomous predator-prey system with the Beddington-DeAngelis functional response

[J]. Journal of Mathematical Analysis and Applications, 2004, 295(1): 15-39.

[8] AJRALDI V, PITTAVINO M, VENTURINO E. Modeling herd behavior in population systems[J]. Nonlinear Analysis: Real

World Applications, 2011, 12(4): 2319-2338.

[9] KOSOBUD R F, O'NEILL W D. On the dependence of population growth on income: New results in a ricardian-Malthus

model[J]. De Economist, 1981, 129(2): 206-223.

[10] TSOULARIS A, WALLACE J. Analysis of logistic growth models[J]. Mathematical Biosciences, 2002, 179(1): 21-55.

[11] BOUKAL D S, SABELIS M W, BEREC L. How predator functional responses and Allee effects in prey affect the paradox of

enrichment and population collapses[J]. Theoretical Population Biology, 2007, 72(1): 136-147.

114

http:∥xuebao.gxnu.edu.cn

[12] KRAMER A M, DENNIS B, LIEBHOLD A M, et al. The evidence for Allee effects[ J]. Population Ecology, 2009, 51

(3): 341-354.

[13] BEREC L, ANGULO E, COURCHAMP F. Multiple Allee effects and population management[ J]. Trends in Ecology &

Evolution, 2007, 22(4): 185-191.

[14] BISWAS S, SAIFUDDIN M, SASMAL S K, et al. A delayed prey-predator system with prey subject to the strong Allee

effect and disease[J]. Nonlinear Dynamics, 2016, 84(3): 1569-1594.

[15] COURCHAMP F, BEREC L, GASCOIGNE J. Allee effects in ecology and conservation[ J]. Environmental Conservation,

2008, 36(1): 80-85.

[16] WU D Y, ZHAO H Y. Complex dynamics of a discrete predator-prey model with the prey subject to the Allee effect[ J].

Journal of Difference Equations Applications, 2017, 23(11): 1765-1806.

[17] PETROVSKII S, MOROZOV A, LI B L. Regimes of biological invasion in a predator-prey system with the Allee effect[J].

Bulletin of Mathematical Biology, 2005, 67(3): 637-661.

[18] RAO F, KANG Y. The complex dynamics of a diffusive prey-predator model with an Allee effect in prey[ J]. Ecological

Complexity, 2016, 28: 123-144.

[19] BARCLAY H J. Models for pest control using predator release, habitat management and pesticide release in combination

[J]. Journal of Applied Ecology, 1982, 19(2): 337-348.

[20] TANG S Y, TANG G Y, CHEKE R A. Optimum timing for integrated pest management: Modelling rates of pesticide

application and natural enemy releases[J]. Journal of Theoretical Biology, 2010, 264(2): 623-638.

[21] 成定平. 鼠類-天敵系統(tǒng)漸近穩(wěn)定性的數(shù)學(xué)分析[J]. 生物數(shù)學(xué)學(xué)報, 2003, 18(3): 283-286.

[22] CHEN Y M, ZHANG F Q. Dynamics of a delayed predator-prey model with predator migration[ J]. Applied Mathematical

Modelling, 2013, 37(3): 1400-1412.

[23] JIN M, XU F, SHEN C S, et al. Nontrivial effect of time varying migration on the three species prey-predator system[ J].

Communications in Theoretical Physics, 2019, 71(1): 127-131.

[24] ABDLLAOUI A E, AUGER P, KOOI B W, et al. Effects of density-dependent migrations on stability of a two-patch

predator-prey model[J]. Mathematical Biosciences, 2007, 210(1): 335-354.

[25] HUANG Y, DIEKMANN O. Predator migration in response to prey density: What are the consequences? [ J]. Journal of

Mathematical Biology, 2001, 43(6): 561-581.

[26] SUN G Q, JIN Z, LIU Q X, et al. Dynamical complexity of a spatial predator-prey model with migration[ J]. Ecological

Modelling, 2008, 219(1 / 2): 248-255.

[27] PETROVSKII S, LI B L. An exactly solvable model of population dynamics with density-dependent migrations and the Allee

effect[J]. Mathematical Biosciences, 2003, 186(1): 79-91.

Stability of a Prey-predator Model with Migration and Allee Effects

XU Wangjun

1

, CAO Jinde

2

, WU Daiyong

1

, SHEN Chuansheng

1?

(1. School of Mathematics and Physics, Anqing Normal University, Anqing Anhui 246133, China;

2. School of Mathematics, Southeast University, Nanjing Jiangsu 211189, China)

Abstract: A kind of prey-predator system with Allee effect and artificially controlled migration of predators is

studied. The system has a square root functional response function. Firstly, by qualitative analysis of the model,

the boundedness of the solution is proved, and the existence of the equilibrium point is analyzed. Sufficient

conditions for the local stability of the equilibrium point of the system are obtained. Then, the existence of the

Hopf-bifurcation of the equilibrium point is discussed, and the stability and direction of the equilibrium Hopfbifurcation are studied by calculating the first Lyapunov coefficient. Finally, the correctness of the conclusion is

verified by numerical simulation. The results indicate that the Allee effect and artificially controlled migration rate

are important for the survival and extinction of prey and predator populations.

Keywords: prey-predator model; Allee effect; migration rate; Hopf bifurcation; stability

(責(zé)任編輯 蘇凱敏)

115

第 40 卷 第 2 期

2022 年 3 月

廣西師范大學(xué)學(xué)報(自然科學(xué)版)

Journal of Guangxi Normal University (Natural Science Edition)

Vol. 40 No. 2

Mar. 2022

DOI: 10.16088 / j.issn.1001-6600.2021052401 http: xuebao.gxnu.edu.cn

蔣群群, 王林峰.一類非線性 p-Laplace 方程的 Liouville 定理[ J]. 廣西師范大學(xué)學(xué)報(自然科學(xué)版), 2022, 40(2): 116-124. JIANG Q Q,

WANG L F. Liouville theorems for a nonlinear p-Laplace equation[J]. Journal of Guangxi Normal University (Natural Science Edition), 2022, 40

(2): 116-124.

一類非線性 p-Laplace 方程的 Liouville 定理

蔣群群, 王林峰?

(南通大學(xué) 理學(xué)院, 江蘇 南通 226019)

摘 要: 在帶有適當(dāng)曲率條件的完備流形上研究非線性 p-Laplace 方程 Δpu+au

p-1

ln u+λu

p-1 = 0, 式中 a、 λ 和 p>1 為給定

常數(shù)。 通過考慮幾何量沿 p-Laplace 方程的演化, 在 Ricci 曲率有下界的緊致流形上建立上述方程的微分不等式。 借助截

斷函數(shù)及 Hessian 比較定理, 在截面曲率有下界的非緊流形上也建立類似不等式。 作為應(yīng)用得到了 Liouville 定理。

關(guān)鍵詞: 微分不等式; 非線性; p-Laplace 方程; Liouville 定理; Ricci 孤立子

中圖分類號: O186.1 文獻(xiàn)標(biāo)志碼: A 文章編號: 1001-6600(2022)02-0116-09

假設(shè) M 是一個帶度量 g 的 n 維完備流形。 本文研究如下非線性 p-Laplace 方程

Δpu+au

p-1

ln u+λu

p-1 = 0, (1)

式中 p>1,且 Δpu = div( | u |

p-2

u)。

如果 p = 2, 式(1)為

Δu+auln u+λu = 0。 (2)

文獻(xiàn)[1]得到 Ricci 曲率為負(fù)下界的完備流形上正調(diào)和函數(shù)的一個最優(yōu)微分不等式。 λ 為正時方程 Δu+

λu = 0 的微分不等式的證明方法類似于文獻(xiàn)[1],也可參看文獻(xiàn)[2]。 文獻(xiàn)[1]不等式中的等號能夠取到,

因此這個微分不等式是最優(yōu)的,這也意味著著名的 Yau 的 Liouville 定理[3]

, 即在具有非負(fù) Ricci 曲率的完

備流形上不存在非常數(shù)正的調(diào)和函數(shù)。 光滑度量測量空間上類似的最優(yōu)微分不等式在文獻(xiàn)[4] 中被

建立。

當(dāng) a<0 時,Ma

[5]得到式(2)的正解的一個微分不等式,還通過觀察式(2)與膨脹梯度 Ricci 孤立子之

間的關(guān)系,說明該微分不等式最優(yōu)。 膨脹梯度 Ricci 孤立子由下式定義[6]

Ric+Hess f =

a

2

g, (3)

式中: f 為勢函數(shù);a<0。 Ma 觀察到如果 f 是式(3)的勢函數(shù),則對某個常數(shù) λ,u = e

-f滿足式 (2)。

注意式(2)也與 Perelman 的 W 函數(shù)有關(guān)。 定義

λ = inf W(u) u > 0,u ∈ C

∞

(M),∫

M

u

2 { dx = 1} ,

式中

W(u) = ∫

M

(4 u

2 - u

2

ln u

2

)dx。

文獻(xiàn)[7]證明在緊流形上 λ 可以由滿足

4Δu+2uln u+λu = 0

的正光滑函數(shù) u 取到。 文獻(xiàn)[6]利用這個最小化結(jié)果來排除緊致流形上非平凡收縮的梯度 Ricci 孤立子

的存在性。 在其他情況下 W 泛函的研究可以參看文獻(xiàn)[8-9]。 由文獻(xiàn)[7]可知,式(2)也與著名的 Gross

收稿日期: 2021-05-24 修回日期: 2021-07-08

基金項目: 國家自然科學(xué)基金(11771223)

通信作者: 王林峰(1973—), 男, 江蘇南通人, 南通大學(xué)教授, 博士。 E-mail: wlf711178@ntu.edu.cn

http:∥xuebao.gxnu.edu.cn

對數(shù) Sobolev 不等式密切相關(guān)。

對 p>1,Kotschwar 等[10]在截面曲率有下界的假設(shè)下建立 p-調(diào)和函數(shù)的局部微分不等式。 由這個微分

不等式能推出 Liouville 定理,該定理表明在截面曲率非負(fù)的完備流形上不存在非常數(shù)正的 p-調(diào)和函數(shù)。

Wang 等[11]在截面曲率有下界的非緊流形上證明方程 Δpu+λu

p-1 = 0 正解的一個最優(yōu)微分不等式。 與 pLaplace 相關(guān)的 Liouville 定理的研究還可以參看文獻(xiàn)[12-13]及其中的參考文獻(xiàn)。

對于 p>1 和某個給定常數(shù) a,定義

λa,p

= inf Wa,p(u) u > 0,u ∈ C

∞

(M),∫

M

u

p { dx = 1} ,

式中

Wa,p(u) = ∫

M

u

p -

a

p

u

p

lnu

p

( ) dx。

在緊致流形上能使 λa,p取到的函數(shù) u 滿足 Euler-Lagrange 方程(1),式中 λ =λa,p。 關(guān)于 p-Laplace 方程的微

分不等式已有不少研究[14-15]

。 本文將建立非線性 p-Laplace 方程(1)正解的微分不等式。

1 本文主要結(jié)果及預(yù)備知識

為方便起見,設(shè)

? =

2

n+2

, a>0,

6(p-1)

6(p-1)+(4(p-1)+1)n

, a<0,

ì

?

í

?

?

?

?

(4)

以及

A =

(p-1)

1

4

(p-1)+

1

2

( ) n,a>0,

4(p-1)

2+8(p-1)+1

12

n, a<0。

ì

?

í

?

?

?

?

(5)

下面結(jié)果是緊流形上的微分不等式。

定理 1 設(shè) M 為 Ricci 曲率以-(n-1)K 為下界的 n 維緊致流形,式中 K≥0。 假設(shè) u 是式(1)的正解,

令 u

? = v,并且 M= sup

x∈M

{u(x)} ,則

① 對 a>0,如果 p≥2,則

v ≤ p

1

p a

1

p (2A)

-

1

p ?

p-2

p +A

-

1

2 ((n-1)K)

1

[ 2 ] M

?

, (6)

如果 1<p≤2,又有

v ≤ p

1

2 (p-1)

-

1

2 a

1

p (2A)

-

1

p ?

p-2

p +(p-1)

-

1

2 A

-

1

2 ((n-1)K)

1

[ 2 ] M

?

。 (7)

② 對 a<0,如果 p≥2,則

v ≤max pa(2A)

-1

?

p-2+A

-

p

2 ((n-1)K)

p

{ 2 ,0}

1

p M

?

, (8)

如果 1<p≤2,又有

v ≤A

-

1

2 ((n-1)K)

1

2 M

?

。 (9)

定理 1 能推出如下 Liouville 定理。

定理 2 設(shè) M 為 Ricci 曲率以-(n-1)K 為下界的 n 維緊致流形,式中 K≥0。 假設(shè) u 是 a<0 時式(1)

的正解。 對 p≥2,如果

pa(2A)

-1

?

p-2+A

-

p

2 ((n-1)K)

p

2 ≤0, (10)

或者當(dāng) 1<p<2 時 K = 0,則 u≡e

-

λ

a 是常數(shù)。

117

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

流形非緊時也可以建立式(1)的微分不等式。 引入截斷函數(shù)時計算中包含距離函數(shù)的 Hessian,所以

假設(shè)截面曲率有下界[10]

。

定理 3 設(shè) M 為 n 維完備非緊流形,截斷曲率以-K 為下界,式中 K≥0。 假設(shè) 1<p≤2,u 是式(1)的正

解,設(shè) M= sup

x∈M

{u(x)} ,如果 a>0,那么 v = u

? 滿足

v ≤ p

1

2 (p-1)

-

1

2 a

1

p (2A)

-

1

p ?

p-2

p +(p-1)

-

1

2 A

-

1

2 ((n-1)K)

1

[ 2 ] M

?

, (11)

如果 a<0,則

v ≤A

-

1

2 ((n-1)K)

1

2 M

?

。 (12)

如第 3 章推論 2 所述,由定理 3 能得到關(guān)于截面曲率為負(fù)的非緊流形的 Liouville 定理。

假設(shè) u>0 是式(1)的連續(xù)弱解[10,15]

(以下簡稱解),并且 M= sup

x∈M

{u(x)} 。 對任意 δ>0,存在只依賴于 p

和 δ 的常數(shù) C(p,δ),使得

u

p-1

ln u ≤C(p,δ)max 1,u

p-1+δ

{ }≤C(p,δ)max 1,M

p-1+δ

{ } ,

則文獻(xiàn)[16]中的梯度估計意味對某些 α>0, u

2∈C

α

,且對某些 β>1 有 u

2∈W

1,β

loc 。 實際上,u 在

{ u≠0}處光滑。 由 v = u

?

,得 v = ?u

?-1

u,ln v = ?ln u,且式(1)可改寫為

Δp

v =

(p-1)(?-1)

?

v

p

v

-a?

p-2

v

p-1

ln v-λ?

p-1

v

p-1

。 (13)

對給定的滿足式(13)的函數(shù) v,定義算子 L 為

L(w)= div( v

p-2 w)+(p-2)div( v

p-4

( v· w) v)。

下面引理可以看作是 Bochner 公式的推廣(見文獻(xiàn)[17]式(2.5))。

引理 1 在任意 v≠0 點處,有

L( v

p

)= p v

2p-4 ‖Hess v‖2

g?

( +Ric( v, v) ) +p v

p-2 Δp

v· v, (14)

式中

g?= g+(p-2)

dv?dv

v

2

,

‖Hess v‖2

?g

= Hess v

2

g

+

p-2

2

v

2 2

v

2

+

(p-2)

2

4

v· v

2 2

v

4

。

基于引理 1 可以推導(dǎo)出以下估計。

引理 2 設(shè) M 為 Ricci 曲率以-( n-1)K 為下界的 n 維完備非緊流形,式中 K≥0。 假設(shè) v = u

? 是式

(13)的解,則

L( v

p

)+(n-1)Kp v

2p-2≥

pA

v

2p

v

2

-ap?

p-2

v

p-2

v

p+

p(p-1)(?-1)

?

υ

p-2

v

p· v

v

(15)

在所有 v≠0 的點處成立,這里 ? 和 A 分別見式(4)、(5)。

證明 由式(13)得

Δp

v· v =

(p-1)(?-1)

?

v

p· v

v

-

v

p-2

v

( 2 ) -

(p-1)a?

p-2

ln v+a?

p-2+(p-1)λ?

p-1

[ ] v

p-2

v

2

。 (16)

易知,

‖Hess v‖2

g?≥

1

n

(Trg?(Hess v))

2 =

1

n

Δv+

p-2

2

v· v

2

v

( 2 )

2

=

1

n

( v

2-pΔp

v)

2 =

1

n

v

4-2p

(Δp

v)

2 =

1

n

v

4-2p (p-1)(?-1)

?

v

p

v

- a?

p-2

v

p-1

ln v-λ?

p-1

v

é p-1

?

ê

ê

ù

?

ú

ú

2

。 (17)

把式(16)、(17)代入式(14),得

L( v

p

)≥

p(p-1)

2

(?-1)

2

n?

2

-

p(p-1)(?-1)

? ( )

v

2p

v

2

-

118

http:∥xuebao.gxnu.edu.cn

p(p-1) 1+

2(?-1)

n?

( )

v

p

v

(a?

p-2

v

p-1

ln v+λ?

p-1

v

p-1

)+

p

n

(a?

p-2

v

p-1

ln v+λ?

p-1

v

p-1

)

2-(n-1)Kp v

2p-2-

ap?

p-2

v

p-2

v

p+

p(p-1)(?-1)

?

v

p-2

v

p· v

v

。

? 的定義意味著

a?

p-2

1+

2(?-1)

n?

( ) ≥0。 (18)

假設(shè) δ 滿足

a

n

?

p-2-

p-1

δ

1+

2(?-1)

n?

( ) = 0, (19)

在點 x∈M,如果 a?

p-2

v

p-1

ln v+λ?

p-1

v

p-1≤

a?

p-2

δ v

p

v

,則

L( v

p

)+(n-1)Kp v

2p-2≥

p(p-1)

(p-1)(?-1)

2

n?

2

-

?-1

?

- a?

p-2

δ 1+

2(?-1)

n?

( )

é

?

ê

ê

ù

?

ú

ú

v

2p

v

2

-

ap?

p-2

v

p-2

v

p+

p(p-1)(?-1)

?

v

p-2

v

p· v

v

。 (20)

在點 x 處,如果 a?

p-2

v

p-1

ln v+λ?

p-1

v

p-1≥

a?

p-2

δ v

p

v

,則

L( v

p

)+(n-1)Kp v

2p-2≥

p(p-1)

2

(?-1)

2

n?

2

-

p(p-1)(?-1)

? ( )

v

2p

v

2

+

p(p-1)(?-1)

?

v

p-2

v

p· v

v

+

p

n

-

p(p-1)

a?

p-2

δ

1+

2(?-1)

n?

( )

é

?

ê

ê

ù

?

ú

ú

(a?

p-2

v

p-1

ln v+λ?

p-1

ν

p-1

)

2-ap?

p-2

v

p-2

v

p =

p(p-1)

2

(?-1)

2

n?

2

-

p(p-1)(?-1)

? ( )

v

2p

v

2

- ap?

p-2

v

p-2

v

p+

p(p-1)(?-1)

?

v

p-2

v

p· v

v

,

這里最后一個等式來自式(19)。 因此式(20)在所有點都成立。 由式(19),

L( v

p

)+(n-1)Kp v

2p-2≥

p(p-1) -

3(p-1)

n

?-1

?

( )

2

-(4(p-1)+1)

?-1

?

- n(p-1)

é

?

ê

ê

ù

?

ú

ú

v

2p

v

2

-

ap?

p-2

v

p-2

v

p+

p(p-1)(?-1)

?

v

p-2

v

p· v

v

=

pA

v

2p

v

2

- ap?

p-2

v

p-2

v

p+

p(p-1)(?-1)

?

v

p-2

v

p· v

v

,

從而得到式(15)。 證畢。

注 1 式(4)中選擇的 ? 是最優(yōu)的,因為

-

3(p-1)

n

?-1

?

( )

2

-(4(p-1)+1)

?-1

?

- n(p-1)

在式(18)的條件下達(dá)到最大值。

119

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

2 緊致流形

本章證明定理 1 和定理 2。

定理 1 的證明 假設(shè) v

p 在 x0∈M 取到最大值,并且 v

p

(x0 )>0。 因此,在 x0 處,有 v

p = 0,

L( v

p

)≤0。 由式(15)得

pA

v

2p

v

2 ≤(n-1)Kp v

2p-2+ap?

p-2

v

p-2

v

p

。 (21)

首先考慮 a>0 的情況。 如果 p≥2,容易看出

(n-1)K v

p-2≤

(p-2)A

pv

2

v

p+2((n-1)K)

p

2 p

-

p

2 (p-2)

p-2

2

(p-2)A

pv

( 2 )

-

p

2

+1

。 (22)

把式(22)代入式(21),得

v

p

(x0 )≤ pa(2A)

-1

?

p-2+A

-

p

2 ((n-1)K)

p

[ 2 ] v

p

(x0 )≤

pa(2A)

-1

?

p-2+A

-

p

2 ((n-1)K)

p

[ 2 ] M

?p

。

因此,對任意 x∈M,

v

p

(x)≤ pa(2A)

-1

?

p-2+A

-

p

2 ((n-1)K)

p

[ 2 ] M

?p

。

從而得式(6)。

如果 1<p<2,容易看出

a?

p-2

v

p-2

v

2-p≤(2-p)Av

-2

v

2+p2

-

2

p (a?

p-2

)

2

p A

1-

2

p 。 (23)

將式(23)代入式(21),并與 p≥2 的情況類似進(jìn)行討論,得到式(7)。

對 a<0 的情況,當(dāng) p≥2 時,將式(22)代入式(21),得

v

p

(x0 )≤ pa(2A)

-1

?

p-2+A

-

p

2 ((n-1)K)

p

[ 2 ] v

p

。

當(dāng) 1<p<2 時,由式(21) 直接得 pA

v

2p

v

2 ≤( n-1) Kp v

2p-2在 x0 處成立,于是可以得到式(8)、(9)。

證畢。

現(xiàn)在證明定理 2。

定理 2 的證明 當(dāng) a<0 時,若對 p≥2,K 滿足式(10),或者對 1<p<2,K = 0,從定理 1 可得 v≡0,因

此,u 是常數(shù)。 由式(1)知 u≡e

-

λ

a 。 證畢。

在定理 2 中令 p = 2 得到對式(2)的 Liouville 定理。

推論 1 設(shè) M 為 Ricci 曲率以-(n-1)K 為下界的 n 維緊致流形,式中 K≥0。 假設(shè) a<0 時 u 是式(2)

的正解。 如果 a+(n-1)K≤0,那么 u≡e

-

λ

a 是常數(shù)。

3 非緊流形

本章將利用截斷函數(shù)給出非緊致流形上式(1) 正解的微分不等式。 由于計算涉及到距離函數(shù)的

Hessian,需要假設(shè)截斷曲率有下界。 首先證明定理 3。

定理 3 的證明 設(shè) O∈M 為定點,r( x) = dist(O,x) 為由 O 決定的距離函數(shù)。 對足夠大的 R>0,用

B(O,R)表示中心為 O 半徑為 R 的測地球。

120

http:∥xuebao.gxnu.edu.cn

考慮如下函數(shù)[15]

θ(t):[0,+∞ )→[0,1],θ(t)=

1,0≤t≤1,

0,t≥2, {

使得

-10θ

1

2 ≤θ′≤0,θ″>-10θ。 (24)

定義截斷函數(shù) φ:M→R 為 φ(x,t)= θ

r(x)

R

( ) ,計算表明:

L(φ v

p

)= div( v

p-2

(φ v

p

))+(p-2)div( v

p-4

( v· (φ v

p

)) v)=

φL( v

p

)+2(p-1) v

2p-4

v

2· φ+(p-2) v

p-2

v· φΔp

v +

v

2p-2 Δφ+

(p-2)Hess φ( v, v)

v

( 2 ) +

(p-1)(p-2) v

2p-6

( v· v

2

)( v· φ)。 (25)

設(shè) φ v

p =G,假設(shè) G 在 x0∈M 時取到最大值,不妨假設(shè) G(x0 )>0,那么在 x0 點,φ>0, v≠0, G= 0,即

v

2 = -

2

pφ

v

2

φ。 (26)

在 x0 點,還有 L(G)≤0。 對任意 σ∈(0,1),有

‖Hess v‖2

g?≥

1

n

( v

2-pΔp

v)

2 =

σ

n

v

4-2p

(Δp

v)

2+

1-σ

n

v

4-2p

(Δp

v)

2 =

σ

n

v

4-2p (p-1)(?-1)

?

v

p

v

- a?

p-2

v

p-1

ln v-λ?

p-1

v

é p-1

?

ê

ê

ù

?

ú

ú

2

+

1-σ

n

v

4-2p

(Δp

v)

2

。 (27)

令

?σ

=

2σ

n+2σ

, a≥0,

6(p-1)σ

6(p-1)σ+(4(p-1)+1)n

, a<0,

ì

?

í

?

?

?

?

以及

Aσ

=

(p-1)

1

4

(p-1)+

1

2

( )

n

σ

, a≥0,

4(p-1)

2+8(p-1)+1

12

n

σ

, a<0。

ì

?

í

?

?

?

?

假設(shè) δσ 滿足

aσ

n

?

p-2

σ

-

p-1

δσ

1+

2σ(?σ

-1)

n?σ

( ) = 0。

與式(15)類似,有

L( v

p

)+(n-1)Kp v

2p-2-

p(1-σ)

n

(Δp

v)

2≥

pAσ

v

2p

v

2

- ap?

p-2

σ

v

p-2

v

p+

p(p-1)(?σ

-1)

?σ

v

p-2

v

p· v

v

。

由式(25)、(26)、 (27),以及在 x0 處 L(G)≤0,得在 x0 處有

0≥pAσ φ

v

2p

v

2

- apφ?

p-2

σ

v

p-2

v

p-(n-1)Kpφ v

2p-2+

121

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

p(1-σ)φ

n

(Δp

v)

2-

p(p-1)(?σ

-1)

?σ

v

2p-2

v· φ

v

-

4(p-1)

p

v

2p-2 φ

2

φ

+(p-2) v

p-2

v· φΔp

v+

v

2p-2 Δφ+

(p-2)Hessφ( v, v)

v

( 2 ) -

2(p-1)(p-2)

p

v

2p-4

v· φ

2

φ

。

易見

p(1-σ)φ

n

(Δp

v)

2+(p-2) v

p-2

v· φΔp

v≥-

(p-2)

2

n

4p(1-σ)

v

2p-4

v· φ

2

φ

,

以及 v· φ ≤ v φ ,則在 x0 處有

0≥pAσ φ

v

2p

v

2

- apφ?

p-2

σ

v

p-2

v

p- (n-1)Kpφ v

2p-2+

p(1-σ)φ

n

(Δp

v)

2-

p(p-1)(?σ

-1)

?σ

v

2p-2

v· φ

v

-

4(p-1)

p

v

2p-2 φ

2

φ

+(p-2) v

p-2

v· φΔp

v+

v

2p-2 Δφ+

(p-2)Hess φ( v, v)

v

( 2 ) -

2(p-1)(p-2)

p

v

2p-4

v· φ

2

φ

,

及

Δφ+(p-2) v

-2Hess φ( v, v)=

θ″

R

2

g+

p-2

v

2

dv?dv ( ) ( r, r)+

θ′

R

<g+

p-2

v

2

dv?dv,Hess r>。

由式(24)、 Hessian 比較定理[18]

,Hess r≤

1+Kr

r

g,以及 g?= g+

p-2

v

2

dv?dv 正定,得

Δφ+

(p-2)Hess φ( v, v)

v

2 ≥-

10φ

R

2

-

10(n+p-2) φ

R

1+KR

R

。

易見 φ =

θ′

R

≤

10 φ

R

。 那么在 x0 處有

0≥pAσ

G

2

v

2

φ

- ap?

p-2

σ

v

p-2G-(n-1)Kpφ

2

p

- 1G

2 -

2

p +

10p(p-1)(?σ

-1)

?σ

φ

1

p

-

3

2

Rv

G

2 -

1

p -

400(p-1)

pR

2

φ

2

p

-2G

2 -

2

p -

10φ

R

2

+

10(n+p-2) φ

R

1+ K R

R

é

?

ê

ê

ù

?

ú

ú

φ

2

p

- 2G

2 -

2

p -

2(p-1) p-2

p

+

(p-2)

2

n

4p(1-σ)

é

?

ê

ê

ù

?

ú

ú

100

R

2

φ

2

p

-2G

2 -

2

p 。

由 0≤φ≤1,得

pAσ

G

2

v

2 ≤ap?

p-2

σ

v

p-2G+(n-1)KpG

2 -

2

p +

C1φ

1

p

-

1

2

Rv

G

2 -

1

p +

C2

+C3R

R

2

φ

2

p

-1G

2 -

2

p , (28)

式中 Ci,i = 1,2,3,…是依賴于 n、p、σ、K、?σ 的常數(shù)。 由 1<p<2,0≤φ≤1,得

pAσ

G

2

v

2 ≤ap?

p-2

σ

v

p-2G+(n-1)KpG

2-

2

p +

C1

Rv

G

2-

1

p +

C2

+C3R

R

2

G

2-

2

p 。 (29)

首先考慮 a>0 的情形。 易知

a?

p-2

σ

v

p-2G≤(2-p)Aσ

v

-2G

2+2

-

2

p p(a?

p-2

σ )

2

p A

1-

2

p

σ G

2-

2

p 。 (30)

122

http:∥xuebao.gxnu.edu.cn

容易看出,對充分小的 δ?>0,

C1

Rv

G

2-

1

p ≤δ?Aσ

G

2

v

2

+

C

2

1

4δ?Aσ R

2

G

2-

2

p 。 (31)

將式(30)、(31)代入式(29),得

(p(p-1)-δ?)Aσ G

2

p (x0 )≤ p

2

2

-

2

p (a?

p-2

σ )

2

p A

1-

2

p

σ

+(n-1)Kp+

C

2

1

4δ?Aσ R

2

+

C2

+C3R

R

2

é

?

ê

ê

ù

?

ú

ú

M

2

。

因此,對所有 x∈B(O,R),

(p(p-1)-δ?)Aσ

v

2

(x)= (p(p-1)-δ?)Aσ G

2

p (x)≤(p(p-1)-δ?)Aσ G

2

p (x0 )≤

p

2

2

-

2

p (a?

p-2

σ )

2

p A

1-

2

p

σ

+(n-1)Kp+

C

2

1

4δ?Aσ R

2

+

C2

+C3R

R

2

é

?

ê

ê

ù

?

ú

ú

M

2

。

令 R→∞ ,對所有 x∈M,

(p(p-1)-δ?)Aσ

v

2

(x)≤ p

2

2

-

2

p (a?

p-2

σ )

2

p A

1-

2

p

σ

[ +(n-1)Kp] M

2

。

令 δ?↘0,得

(p-1)Aσ

v

2

(x)≤ p2

-

2

p (a?

p-2

σ )

2

p A

1-

2

p

σ

[ +(n-1)K] M

2

。 (32)

注意當(dāng) σ↗1 時,?σ→?,Aσ→A,在式(32)中令 σ↗1,得

(p-1)A v

2

(x)≤ p2

-

2

p (a?

p-2

)

2

p A

1-

2

[ p +(n-1)K] M

2

。

這個不等式意味著式(11),通過類似討論可以得到式(12)。 證畢。

推論 2 在具有非負(fù)截面曲率的 n 維非緊流形上,當(dāng) a<0,且 1<p≤2 時,式(1)無非平凡的正有界解。

參 考 文 獻(xiàn)

[1] LI P, WANG J P. Complete manifolds with positive spectrum, II[J]. Journal of Differential Geometry, 2002, 62: 143-162.

[2] LI P. Geometric analysis[M]. Cambridge: Cambridge University Press, 2012.

[3] YAU S T. Harmonic functions on complete Riemannian manifolds[J]. Communications on Pure and Applied Mathematics,

1975, 28(2): 201-228.

[4] WANG L F. The upper bound of the L

2

μ

spectrum[J]. Annals of Global Analysis and Geometry, 2010, 37(4) : 393-402.

[5] MA L. Gradient estimates for a simple elliptic equation on complete non-compact Riemannian manifolds[ J]. Journal of

Functional Analysis, 2006, 241(1): 374-382.

[6] PERELMAN G. The entropy formula for the Ricci flow and its geometric applications[EB/ OL]. (2002-11-11) [2021-05-

24]. https:∥arxiv.org / abs/ math / 0211159.

[7] ROTHAUS S. Logarithmic Sobolev inequalities and the spectrum of Schr?dinger operators [ J ]. Journal of Functional

Analysis, 1981, 42(1): 110-120.

[8] NI L. The entropy formula for linear heat equation[J]. Journal of Geometric Analysis, 2004, 14(1): 87-100.

[9] WANG L F. A new entropy formula for the linear heat equation[J]. Archiv der Mathematik, 2011, 96(5): 473-481.

[10] KOTSCHWAR B, NI L. Local gradient estimates of p-harmonic functions, 1 / H-flow, and an entropy formula[ J]. Annales

Scientifiques de l′école Normale Supérieure, 2009, 42(1): 1-36.

[11] WANG L F, ZHU Y P. A sharp gradient estimate for the weighted p-Laplacian[ J]. Applied Mathematics-A Journal of

Chinese Universities, 2012, 27(4): 462-474.

[12] 周俞潔, 張澤宇, 王林峰. 黎曼流形上 p-Laplace 算子的 Liouville 定理[ J]. 西南大學(xué)學(xué)報(自然科學(xué)版), 2017,

123

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

39(10): 62-68.

[13] 魏公明. 具奇系數(shù)發(fā)展型 p-Laplace 不等方程整體解的不存在性[J]. 數(shù)學(xué)年刊 A 輯, 2007, 28(3): 387-394.

[14] 汪悅. Riemann 流形上改進(jìn)的 p-Laplace 方程的梯度估計[J]. 中國科學(xué)(數(shù)學(xué)), 2014, 44(3): 288-294.

[15] 王林峰. p-Laplace Schr?dinger 熱方程的橢圓型梯度估計[J]. 數(shù)學(xué)學(xué)報, 2010, 53(4): 643-654.

[16] TOLKSDORF P. Regularity for a more general class of quasilinear elliptic equations[ J]. Journal of Differential Equations,

1984, 51(1): 126-150.

[17] WANG L F. Gradient estimates on the weighted p-Laplace heat equation [ J]. Journal of Differential Equations, 2018,

264(1): 506-524.

[18] SCHEON R, YAU S T. Lectures on differential geometry[M]. Cambridge, MA: International Press, 1994.

Liouville Theorems for a Nonlinear p-Laplace Equation

JIANG Qunqun, WANG Linfeng

?

(School of Sciences, Nantong University, Nantong Jiangsu 226019, China)

Abstract: In this paper the nonlinear p-Laplace equation Δpu + au

p-1

lnu + λu

p-1 = 0 is studied on complete

manifolds with some suitable curvature condition, where a, λ and p > 1 are some given constants. Differential

inequalities for the p-Laplace equation on compact manifolds with Ricci curvature bounded from below are

established, based on the evolution of the geometric quantity along the p-Laplace equation. Similar inequalities

can also be established on a noncompact manifold whose sectional curvature is bounded from below, based on the

skills of cut off function and the Hessian comparison theorem. As an application, Liouville theorems are obtained.

Keywords: differential inequality; nonlinear; p-Laplace equation; Liouville theorem; Ricci soliton

(責(zé)任編輯 吳佃華)

124

第 40 卷 第 2 期

2022 年 3 月

廣西師范大學(xué)學(xué)報(自然科學(xué)版)

Journal of Guangxi Normal University (Natural Science Edition)

Vol. 40 No. 2

Mar. 2022

DOI: 10.16088 / j.issn.1001-6600.2021061201 http: xuebao.gxnu.edu.cn

王涵, 張映輝. 模擬趨化現(xiàn)象的三維雙曲-拋物系統(tǒng)的最優(yōu)衰減率[ J]. 廣西師范大學(xué)學(xué)報(自然科學(xué)版), 2022, 40( 2): 125-131. WANG

H, ZHANG Y H. Optimal time-decay rates of the hyperbolic-parabolic system modeling chemotaxis in R

3

[ J]. Journal of Guangxi Normal University

(Natural Science Edition), 2022, 40(2): 125-131.

模擬趨化現(xiàn)象的三維雙曲-拋物系統(tǒng)的最優(yōu)衰減率

王 涵, 張映輝?

(廣西師范大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院, 廣西 桂林 541006)

摘 要: 本文研究一個模擬趨化現(xiàn)象的三維雙曲-拋物系統(tǒng)的 Cauchy 問題解的大時間行為, 得到其解及其各階空間導(dǎo)數(shù)

的最優(yōu)時間衰減率。 跟已有結(jié)果相比, 本文主要創(chuàng)新在于給出解的最高階空間導(dǎo)數(shù)的最優(yōu)衰減率, 且該衰減率與熱方程

的衰減率一樣。 研究方法主要基于高頻-低頻分解和精細(xì)的能量估計。

關(guān)鍵詞: 雙曲-拋物系統(tǒng); 最優(yōu)衰減速率; 高頻-低頻分解; 大時間行為; 趨化現(xiàn)象

中圖分類號: O29 文獻(xiàn)標(biāo)志碼: A 文章編號: 1001-6600(2022)02-0125-07

1 本文主要成果

本文研究一個模擬趨化現(xiàn)象的三維雙曲-拋物系統(tǒng)解的時間衰減率,該系統(tǒng)可表示為如下形式:

vt

- v = 0,x∈R

3

, t>0,

ut

- ·(uv)= DΔu,x∈R

3

, t>0。 { (1)

初值為

(v, u)(x,0)= (v0 ,u0 )(x)→(0,u?), | x |→∞ , (2)

式中常數(shù) u?>0。 式(1) 與以下系統(tǒng)密切相關(guān),

?p

?t

=D · p ln

p

Φ(w)

( ( ) ) ,

?w

?t

= βpw。

ì

?

í

?

?

?

?

(3)

該模型由 Othmer 等[1]以及 Levine 等[2]通過大量生物學(xué)估計和數(shù)值計算得到。 式中:p( x,t)表示粒子密

度;w(x,t)表示化學(xué)物質(zhì)濃度;D 表示粒子擴(kuò)散速率且 D>0;Φ 表示化學(xué)勢能。

類似文獻(xiàn)[3-9],假設(shè) Φ(w)= w

-α

,α 是大于零的常數(shù),系統(tǒng)(3)可以改寫為

pt

=DΔp+Dα · p

w

w

( ) ,

wt

= βpw,

ì

?

í

??

?

(4)

進(jìn)一步,假設(shè)

q = (ln w)=

w

w

,

系統(tǒng)(4)可改寫為

pt

=DΔp+Dα ·(pq),

qt

= β p。 { (5)

收稿日期: 2021-06-12 修回日期: 2021-07-02

基金項目: 國家 自 然 科 學(xué) 基 金 ( 11771150, 11571280 ); 廣 西 自 然 科 學(xué) 基 金 ( 2019JJG110003, 2019AC20214,

2019JJA110071)

通信作者: 張映輝(1981—), 男, 湖南祁陽人, 廣西師范大學(xué)教授, 博士。 E-mail: yinghuizhang@gxnu.edu.cn

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

最后,正常數(shù) A、B 和 c1 由以下式子給出。 令 τ = At,ξ =Bx,u = p,v = c1

q, 則將系統(tǒng)(5)改寫為

vτ

=

βBc1

A

ξu,

uτ

=

DB

2

A

Δξu+

DαB

Ac1

ξ·(uv)。

ì

?

í

?

?

?

?

(6)

若在式(6)中取

βBc1

A

= 1,

B

2

A

= 1,

DαB

Ac1

= 1,

ì

?

í

?

?

?

?

?

?

?

?

即,

A =Dαβ>0,B = Dαβ >0,c1

=

Dα

β

>0,

那么易知 u 和 v 滿足

vτ

- ξu = 0,

uτ

- ξ·(uv)= DΔξu。 { (7)

用變量(x,t)替換(τ,ξ),則式(7)恰好是式(1)。

對系統(tǒng)(1)的一維情形, 文獻(xiàn)[5,10]和[3]分別考慮了初邊值問題和 Cauchy 問題光滑解的存在性和

漸近性;文獻(xiàn)[5]考慮了系統(tǒng)(1)的初邊值問題, 當(dāng) u0

-1

2

H2+ v0

2

H2足夠小時,證明系統(tǒng)(1)光滑解的全局

存在性;文獻(xiàn)[3]證明在初始值充分大的條件下,系統(tǒng)(1)Cauchy 問題光滑解的全局存在性;文獻(xiàn)[11]考

慮帶流量限制的趨化模型的初邊值問題, 證明經(jīng)典解整體存在且關(guān)于時間一致有界;文獻(xiàn)[12]研究非線

性擴(kuò)散-趨化聚集模型,證明生物趨化模型弱解的全局有界性及唯一性;文獻(xiàn)[13] 得到具有不同分?jǐn)?shù)階擴(kuò)

散趨化模型的衰減估計;文獻(xiàn)[14] 研究一個耦合雙曲-拋物系統(tǒng)的全局光滑解;對于高維情形,文獻(xiàn)

[9,15-17]分別研究系統(tǒng)(1)初邊值問題和 Cauchy 問題的光滑小解全局存在性。 但是, 對高維趨化模型

全局解的大時間行為方面的研究結(jié)果很少。 文獻(xiàn)[4]首次詳細(xì)證明在定常狀態(tài)下, 初始擾動的 H

2 范數(shù)足

夠小,并且 L

1 范數(shù)有界時,三維問題強解的全局存在性和最優(yōu)時間衰減速率

(v,u)(t) L

2≤C0(1+t)

-

3

4 和 (v,u)(t) H1≤C0(1+t)

-

5

4 。 (8)

最近, 文獻(xiàn)[14]得到當(dāng)(v0 ,u0 )∈H

N∩L

p

(p∈[1,2]), 且 (v0 ,u0 ) H2≤δ0 時,L

p -L

2 型時間衰減速率,其

主要結(jié)果為

?(v, u)(t) HN-?≤C0(1+t)

-

3

2

1

p

-

1

2 ( ) +

?

2 ,? = 0,…,N-1。 (9)

一方面,容易發(fā)現(xiàn)式(8)中解(v,u)的 2 階(即最高階)空間導(dǎo)數(shù)的 L

2 衰減速率是(1+t)

-

5

4 ,與它的 1

階空間導(dǎo)數(shù)具有相同的衰減速率,并且比熱方程的 L

2 衰減速率(1+t)

-

7

4 慢。 另一方面,如果在式(9)中取

p = 1,? =N-1,那么(v,u)的 N 階(即最高階)空間導(dǎo)數(shù)的 L

2 衰減速率為(1+ t)

-

3+2(N-1)

4 ,與它的 N-1 階空間

導(dǎo)數(shù)的衰減速率相同,并且比熱方程的 L

2 衰減速率(1+t)

-

3+2N

4 慢。 因此,在這種意義下,式(8)、(9)中解的

最高階空間導(dǎo)數(shù)的衰減率都不是最優(yōu)。

本文旨在對上述問題給出一個明確答案,更精確地說,本文得到 Cauchy 問題(1)、(2)的解及其各階

空間導(dǎo)數(shù)(1 階到最高階 N 階)的最優(yōu)時間衰減速率。 值得一提的是,本文得到解的任一階空間導(dǎo)數(shù)的衰

減速率都與熱方程的衰減速率一樣,證明方法是基于高頻-低頻分解和精細(xì)的能量估計。

類似文獻(xiàn)[4],重新改寫 Cauchy 問題(1)、(2)。 假設(shè)

126

http:∥xuebao.gxnu.edu.cn

λ =

1

u?

,λ1

= u?。

作變量替換(v,u)→(λv,u+u?),在(0,u?)附近將系統(tǒng)(1)線性化,Cauchy 問題(1)、(2)改寫為

vt

-λ1 u = 0,

ut

-λ1 ·v-DΔu =λ ·(uv),

(v,u)(x,0)= (v0 ,u0 )(x)→(0,0), | x |→∞ 。

ì

?

í

?

?

??

(10)

在本處以及后續(xù)工作中,為了簡單起見,仍用(v,u)表示重新構(gòu)造的變量。 為方便起見,本文采用以下記

號:H

s

(R

3

)(s∈R)表示通常 Sobolev 空間,具有范數(shù) · Hs;L

p

(R

3

) (1≤p≤∞ )表示通常的 L

p 空間,且范

數(shù)為 · L

p;

? 表示任意 ? 階空間導(dǎo)數(shù),式中 ?≥0;當(dāng) ?<0 或 ? 不是正整數(shù)時,

? 表示 Λ

? 并且 Λ

?

f: =

F

-1

( | ξ |

iF f),式中 F 表示通常的 Fourier 變換算子,F

-1表示其逆算子。

本文假設(shè) C 為大于零的常數(shù),符號 a?b 表示 a≤Cb,式中 C 只依賴于問題參數(shù),且 C>0。 設(shè)存在徑向

函數(shù) φ∈C

∞

0 (R

3

ξ),當(dāng)| ξ |≤1 時滿足 φ(ξ)= 1,當(dāng) | ξ | ≥2 時滿足 φ(ξ)= 0。 定義 f 的低頻部分和高頻部分

分別為

f

l =F

-1

[φ(ξ)

^

f], f

h=F

-1

[(1-φ(ξ))

^

f]。

為了簡單起見, 記 (A,B) X := A X

+ B X 。

本文主要結(jié)果如下。

定理 1 假設(shè) ×v0

= 0,(v0 ,u0 )∈H

N∩ L

1

,整數(shù) N≥2,v0 、u0∈L

1

,存在常數(shù) δ0 ,使得

(v0 ,u0 ) H2≤δ0 , (11)

則 Cauchy 問題(10)存在一個唯一全局解(v,u),使得對所有 t≥0,以下衰減估計成立:

k(v,u)(t) HN-k?(1+t)

-

3+2k

4 ,k = 0,…,N。 (12)

注 1 與文獻(xiàn)[18]的結(jié)果(式(9))相比, 本文主要結(jié)果(式(11))中關(guān)于解(v,u)的最高階(即 N)空

間導(dǎo)數(shù)的 L

2 衰減速率的結(jié)果是全新的。 特別地,本文結(jié)果蘊含解(v,u)的最高階(N)空間導(dǎo)數(shù)的 L

2 衰減

速率為(1+t)

-

3+2N

4 ,比文獻(xiàn)[18]的結(jié)果(式(9))中的(1+t)

-

3+2(N-1)

4 要快, 并且該衰減速率與熱方程的衰減速

率相同。 因此,本文的衰減速率是最優(yōu)的。

現(xiàn)在簡單給出本文分析要點。 基于文獻(xiàn)[18]結(jié)果,只需要證明式(12)中關(guān)于解( v,u)的最高階(N)

空間導(dǎo)數(shù)的 L

2 衰減速率。 證明主要包括如下 3 個步驟。

首先得到解的能量估計:

d

dt

N(v,u)

2

L

2+C

N+1 u

2

L

2?δ(

N

v

2

L

2+ N+1 u

2

L

2 )。 (13)

其次,注意到上述能量估計中不包含 v 的耗散。 因此,為了揭示 v 的耗散,可以利用文獻(xiàn)[14]方法,

即通過構(gòu)造 v 和 u 之間的相互作用能量泛函,得

-

d

dt

∫

R3

N-1 u

N

vdx + C

N

v

2

L

2 ? Nu

2

H1 。 (14)

注意到式(14)中時間導(dǎo)數(shù)的積分項包含 N-1 u,所以,通過式(13)、(14)似乎不可能封閉解的 N 階能量估

計。 為此將充分利用高頻-低頻分解好的性質(zhì)來克服這一困難。 具體地說,代替構(gòu)造 v 和 u 之間的相互作

用能量泛函,本文構(gòu)造 v

h 和 u

h 之間新的相互作用能量泛函,得

-

d

dt

∫

R3

N-1 u

h N

v

h

dx + C

N

v

h 2

L

2 ? Nu

h 2

L

2 + N+1 u

h 2

L

2 + (1 + t)

-

3+2N

2 。 (15)

于是,取足夠大的 T0 和 D0 ,定義能量泛函

E(t) = D0

N(v,u)

2

L

2 - ∫

R3

N-1 u

h N

v

h

dx,

當(dāng) t≥T0 時,由于 D0 充分大,則 E(t)等價于 N(v,u)

2

L

2 。 因此,結(jié)合式(13)、(15),可得

d

dt

E(t)+CE(t)?(1+t)

-

3+2N

2 + N(v

l

,u

l

)

2

L

2 。 (16)

127

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

最后, 利用 Plancherel 定理、Hausdorff-Young 不等式和 Duhamel 原理可得到解的最高階空間導(dǎo)數(shù)的低

頻衰減速率。

N(v

l

,u

l

)

2

L

2?(1+t)

-

3+2N

2 。 (17)

結(jié)合式(16)、(17),并利用 Gronwall 不等式可得到解的最高階空間導(dǎo)數(shù)的最優(yōu)衰減速率,從而完成定理 1

的證明。

2 非線性能量估計

本節(jié)將推導(dǎo)式(10)的先驗非線性能量估計。 假設(shè)對足夠小的 δ>0,有先驗估計

(v,u)(t) H2≤δ。

本節(jié)將廣泛使用 Gagliardo-Nirenberg Sobolev 插值不等式。

以下引理是文獻(xiàn)[19]第 125 頁定理的特例。

引理 1 假設(shè) 0≤i,j≤k,則

i

f

L

p? Δ

j

f

1-θ

L

q

k

f

θ

L

r,

式中 θ 滿足

i

3

-

1

p

=

j

3

-

1

q

( ) (1-θ)+

k

3

-

1

r

( ) θ。

引理 2 假設(shè)整數(shù) k≥1,則

k(fg) L

p? f

L

p

1

kg L

p

2

+ k

f

L

p

3

g L

p

4

,

式中 p、p2 、p3∈[1,+∞ ],且

1

p

=

1

p1

+

1

p2

=

1

p3

+

1

p4

。

證明 對 p = p2

= p3

= 2,可以利用引理 1 證明,詳情可參看文獻(xiàn)[20]第 98 頁。 對于一般情形,參看文

獻(xiàn)[21]中引理 3.1。

引理 3 若對任意 2≤p≤∞ , f∈L

p

(R

3

),則

f

l

L

p+ f

h

L

p? f

L

p。

證明 對 2≤p≤∞ , 根據(jù)卷積的 Young 不等式,對于低頻可得

f

l

L

p? F

-1

φ L

1 f

L

p? f

L

p,

因此,

f

h

L

p? f

L

p+ f

l

L

p? f

L

p。

證畢。

接下來給出一類包含 u 耗散估計的能量估計。

引理 4 假設(shè) 0≤k≤N-1, 則

d

dt

k+1(v,u)

2

L

2+2D k+2 u

2

L

2?δ(

k+1

v

2

L

2+ k+2 u

2

L

2 )。 (18)

證明 對 0≤k≤N-1,將 k+1作用于式(10)的第 1、2 式子,且將所得等式分別乘以 k+1

v 和 k+1 u,再將

其結(jié)果求和,然后在 R

3 上進(jìn)行分部積分,可得

1

2

d

dt

∫

R3

|

k+1(v,u) |

2

dx + D k+2 u

2

L

2 = ∫

R3

k+1(λ ·(uv))

k+1 udx。 (19)

下面估計等式(19)的右邊項,由分部積分、H?lder 不等式、Cauchy-Schwarz 不等式和引理 1、2,可得

∫

R3

k+1(λ ·(uv))

k+1 udx = - λ∫

R3

k+1(uv)

k+2 udx ? k+1(vu) L

2

k+2 u L

2 ? (

k+1 u L

6 v L

3 +

128

http:∥xuebao.gxnu.edu.cn

k+1

v L

2 u L∞ )

k+2 u L

2 ? ( v H1

k+2 u L

2 + u H2

k+1

v L

2 )

k+2 u L

2 ? δ(

k+1

v

2

L

2 + k+2 u

2

L

2 )。

(20)

結(jié)合式(19)、(20)即可證得式(18)成立。 證畢。

最后,考慮以下線性系統(tǒng):

vt

-λ1 u = 0,

ut

-λ1 ·v-DΔu = 0,

(v,u)(x,0)= (v0 ,u0 )(x)→(0,0), | x |→∞ 。

ì

?

í

?

?

??

(21)

給出線性系統(tǒng)(21)的 L

2 時間衰減速率。

以下引理 5 見文獻(xiàn)[4]。

引理 5 假設(shè)(v?,u?)是系統(tǒng)(21)的解,初值(v?0 ,u?0 )∈H

N∩L

1

,則對 0≤k≤N,有

k(v?

l

,u?

l

)(t) L

2?(1+t)

-

3

4

-

k

2 ( (v?0 ,u?0 ) L

1+ k(v?0 ,u?0 ) L

2 )。

引理 6 在定理 1 的假設(shè)條件下,非線性系統(tǒng)(10)的解(v,u)滿足以下衰減估計

k(v

l

,u

l

)(t) L

2?(1+t)

-

3

4

-

k

2 ,k = 0,…,N。 (22)

證明 定義 S = ( S

1

,S

2

)

t = (0,λ ·( uv))

t

,根據(jù)引理 5、引理 2、引理 3、Plancherel 定理、HausdorffYoung 不等式、Duhamel 原理和 H?lder 不等式,可得

k(v

l

,u

l

)(t) L

2 ? (1 + t)

-

3

4

-

k

2 (v,u)(0) L

1 + ∫

1

2

0

(1 + t - τ)

-

3

4

-

k

2 S(τ) L

1 dτ +

∫

t

t

2

(1 + t - τ)

-

5

4

| ξ |

k-1

S^

l

(τ) L∞ dτ ? (1 + t)

-

3

4

-

k

2 + ∫

t

t

2

(1 + t - τ)

-

5

4

| ξ |

k-1

S^

l

(τ) L∞ dτ。 (23)

另一方面,估計 | ξ |

k-1

S^

l

(τ) L∞ 。

| ξ |

k-1

S^

l

(τ) L∞ ? k-2(λ ·(uv)(τ) L

1 =λ

k-1(uv)(τ) L

1?λ( u(τ) L

2

k-1

v(τ) L

2+

v(τ) L

2

k-1 u(τ) L

2 )?λ( (v,u)(τ) L

2

k-1(v,u)(τ) L

2 )?(1+τ)

-

3

4 (1+τ)

-

3

4

-

k-1

2 ?(1+τ)

-1-

k

2 。 (24)

將式(24)代入式(23)即可證式(22)。 證畢。

3 定理 1 的證明

本章將利用低頻和高頻分解完成定理 1 的證明。 定理 1 中 l≤N-1 的情形在文獻(xiàn)[14]已證明,只需

證明當(dāng) l = N 時

N(v,u)(t) L

2?(1+t)

-

3+2N

4 (25)

成立,即可完成定理 1 的證明。

步驟 1

N

v

h 的耗散。 將 N-1F

-1

(1-φ(ξ))作用到式(10)的第 2 個式子,所得結(jié)果與 N

v

h 相乘,并

在 R

3 上積分,得

∫

R3

N-1 u

h

t

N

v

h

dx - λ1 ∫

R3

N

v

h N

v

h

dx - D∫

R3

N+1 u

h N

v

h

dx = λ∫

R3

N(uv)

h N

v

h

dx。

又因為

∫

R3

N-1 u

h

t

N

v

h

dx =

d

dt

∫

R3

N-1 u

h N

v

h

dx - ∫

R3

N-1 u

h N

v

h

t dx =

d

dt

∫

R3

N-1 u

h N

v

h

dx - λ1 ∫

R3

N-1 u

h N+1 u

h

dx。

于是,

-

d

dt

∫

R3

N-1 u

h N

v

h

dx + λ1

N

v

h 2

L

2 = - λ1 ∫

R3

N-1 u

h N+1 u

h

dx - D∫

R3

N+1 u

h N

v

h

dx -

λ∫

R3

N(uv)

h N

v

h

dx: = I1

+ I2

+ I3 。 (26)

下面對等式(25) 右邊進(jìn)行估計。 首先估計 I1 。 由分部積分得

129

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

| I1

| = λ1

| ∫

R3

N-1 u

h N+1 u

h

dx | ? Nu

h 2

L

2 。

其次估計 I2 。 根據(jù) H?lder 不等式和 Young 不等式,得

| I2

| = D | ∫

R3

N+1 u

h N

v

h

dx | ? N+1 u

h

L

2

N

v

h

L

2 ? N+1 u

h 2

L

2 +

λ1

4

N

v

h 2

L

2 。

最后估計 I3 。 根據(jù) H?lder 不等式、引理 3、引理 2、引理 1、式(8)、(9)和 Young 不等式,得

|I3

| =λ | 〈

N(uv)

h

,

N

v

h

〉 |? N(uv)

h

L

2

N

v

h

L

2? (v,u) L∞

N(v,u) L

2

N

v

h

L

2?

(v,u)

1

2

L

2

2(v,u)

1

2

L

2

N(v,u) L

2

N

v

h

L

2?(1+t)

-

5

4

×

1

2 (1+t)

-

7

4

×

1

2 (1+t)

-

3+2(N-1)

4 N

v

h

L

2?

(1+t)

-

3+2N

4 (1+t)

-1 N

v

h

L

2?(1+t)

-

3+2N

2 +(1+t)

-2 N

v

h 2

L

2 。

因此,

-

d

dt

∫

R3

N-1u

h N

v

h

dx + λ1

N

v

h 2

L

2 ?

λ1

4

+ (1 + t)

-2

( )

N

v

h 2

L

2 + Nu

h 2

L

2 + N+1u

h 2

L

2 + (1 + t)

-

3+2N

2 。 (27)

步驟 2 式(12)的證明。 在式(18)中取 k =N-1,得

d

dt

N(v,u)

2

L2

+2D N+1 u

2

L

2?δ(

N

v

2

L

2+ N+1 u

2

L

2 。 (28)

取足夠大的 T1 和 D1 ,定義能量泛函

E(t) = D1

N(v,u)

2

L

2 - ∫

R3

N-1 u

h N

v

h

dx。

當(dāng) t≥T1 時,由于 D1 充分大,則 E(t)等價于 N(v,u)

2

L

2 。 將式(28)乘以 D1 ,然后將所得結(jié)果與式(27)相

加,注意到當(dāng) t≥T1 時,由于 δ 充分小,再結(jié)合引理 3 的高頻性質(zhì)和引理 6,整理可得

d

dt

E(t)+C1E(t)?(1+t)

-

3+2N

2 + N(v

l

,u

l

)

2

L

2 。

結(jié)合式(22)與 Gronwall 不等式,得

E(t)?(1+t)

-

3+2N

2 。

故式(25)成立。

綜上所述,定理 1 得證。

參 考 文 獻(xiàn)

[1] OTHMER H G, STEVENS A. Aggregation, blowup, and collapse: the ABCs of taxis in reinforced random walks[J]. SIAM

Journal on Applied Mathematics, 1997, 57(4): 1044-1081.

[2] LEVINE H A, SLEEMAN B D. A system of reaction diffusion equations arising in the theory of reinforced random walks[J].

SIAM Journal on Applied Mathematics, 1997, 57(3): 683-730.

[3] GUO J, XIAO J X, ZHAO H J, et al. Global solutions to a hyperbolic-parabolic coupled system with large initial data[J].

Acta Mathematica Scientia, 2009, 29(3): 629-641.

[4] XIE W J, ZHANG Y H, XIAO Y D, et al. Global existence and convergence rates for thestrong solutions in H

2

to the 3D

chemotaxis model[J]. Journal of Applied Mathematics, 2013, 2013(1): 1-9.

[5] ZHANG M, ZHU C J. Global existence of solutions to a hyperbolicparabolic system[ J]. Proceedings of the American

Mathematical Society, 2007, 135(4): 1017-1027.

[6] 張映輝, 譚忠, 賴柏順.一個模擬趨化現(xiàn)象的廣義雙曲-拋物系統(tǒng)的光滑解的全局分析[ J]. 數(shù)學(xué)年刊 A 輯, 2012,

33(1): 27-38.

[7] ZHANG Y H, TAN Z, SUN M B. Global existence and asymptotic behavior of smooth solutions to a coupled hyperbolicparabolic system[J]. Nonlinear Analysis: Real World Applications, 2013, 14 (1): 465-482.

[8] ZHANG Y H, DENG H Y, SUN M B. Global analysis of smooth solutions to a hyperbolic-parabolic coupled system[ J].

130

http:∥xuebao.gxnu.edu.cn

Frontiers of Mathematics in China, 2013, 8(6): 1437-1460.

[9] ZHANG Y H, XIE W J. Global existence and exponential stability for the strong solutions in H

2

to the 3-D chemotaxis model

[J]. Boundary Value Problems, 2015, 116: 1-13.

[10] LI D, PAN R H, ZHAO K. Quantitative decay of a one-dimensional hybrid chemotaxis model with large data [ J].

Nonlinearity, 2015, 28(7): 2181-2210.

[11] 涂馨予, 穆春來, 鄭攀. 帶有流暈限制的擬線性趨化模型解的整體有界性[ J]. 中國科學(xué): 數(shù)學(xué), 2021, 51(6):

1003-1012.

[12] 張婕燕, 辛巧, 穆春來. 具有非線性擴(kuò)散的趨化模型弱解的有界性[ J]. 西南大學(xué)學(xué)報(自然科學(xué)版), 2021,

43(5): 88-95.

[13] 王曦, 劉祖漢, 周玲. 具不同分?jǐn)?shù)階擴(kuò)散趨化模型的衰減估計[J]. 數(shù)學(xué)年刊 A 輯, 2020, 41(2): 175-200.

[14] 張映輝, 譚忠, 孫明保.一個耦合雙曲-拋物系統(tǒng)的全局光滑解[J]. 數(shù)學(xué)年刊 A 輯, 2013, 34(1): 29-46.

[15] LI D, PAN R H, ZHAO K. On a hybrid type chemotaxis model on bounded domains with large data[J]. SIAM Journal on

Applied Mathematics, 2012, 72(7): 417-443.

[16] HAO C C. Global well-posedness for a multidimensional chemotaxis model in critical Besov spaces [ J]. Zeitschrift für

Angewandte Mathematik und Physik, 2012, 63: 825-834.

[17] LI D, LI T, ZHAO K. On a hyperbolic-parabolic system modeling chemotaxis[ J]. Mathematical Models and Methods in

Applied Sciences, 2011, 21(8): 1631-1650.

[18] ZHANG Y H, LI C, XIE W J. Decay of a 3-D hyperbolic parabolic system modeling chemotaxis[J]. Journal of Information

and Optimization Sciences, 2018, 39(7): 1505-1525.

[19] NIRENBERG L. On elliptic partial differential equations [ J]. Annali della Scuola Normale Superiore di Pisa-Classe di

Scienze, 1959, 13(2): 115-162.

[20] MAJDA A J, BERTOZZI A L. Vorticity and incompressible flow[M]. Cambridge: Cambridge University Press, 2002.

[21] JU N. Existence and uniqueness of the solution to the dissipative 2D quasi-geostrophic equations in the Sobolev space[ J].

Communications in Mathematical Physics, 2004, 251: 365-376.

Optimal Time-decay Rates of the Hyperbolic-parabolic

System Modeling Chemotaxis in R

3

WANG Han, ZHANG Yinghui

?

(School of Mathematics and Statistics, Guangxi Normal University, Guilin Guangxi 541006, China)

Abstract: The large-time behavior of solutions to the Cauchy problem of a 3D hyperbolic-parabolic system

modeling chemotaxis is investigated. The optimal time decay rates of the higher-order spatial derivatives of the

solutions are obtained. Compared with previous results, the main innovation of this paper is to give the highest

order spatial derivative of the solutions which is the same as that of the heat equation. The proof is mainly based

on low-frequency and high-frequency decomposition and delicate energy estimates.

Keywords: hyperbolic-parabolic system; optimal decay rates; low-frequency and high-frequency decomposition;

large-time behavior; chemotaxis

(責(zé)任編輯 吳佃華)

131

第 40 卷 第 2 期

2022 年 3 月

廣西師范大學(xué)學(xué)報(自然科學(xué)版)

Journal of Guangxi Normal University (Natural Science Edition)

Vol. 40 No. 2

Mar. 2022

DOI: 10.16088 / j.issn.1001-6600. 2021070902 http: xuebao.gxnu.edu.cn

張琬婧, 林支桂. 增長區(qū)域上一類寄生蟲-宿主模型的 Turing 不穩(wěn)定[J]. 廣西師范大學(xué)學(xué)報(自然科學(xué)版), 2022, 40(2): 132-139. ZHANG

W J, LIN Z G. Turing instability of a parasite-host model on growing domains[J]. Journal of Guangxi Normal University (Natural Science Edition),

2022, 40(2): 132-139.

增長區(qū)域上一類寄生蟲-宿主模型的 Turing 不穩(wěn)定

張琬婧, 林支桂?

(揚州大學(xué) 數(shù)學(xué)科學(xué)學(xué)院, 江蘇 揚州 225002)

摘 要: 為理解區(qū)域演化對寄生蟲-宿主 Turing 不穩(wěn)定的影響, 本文以寄生蟲-宿主傳染病模型為主體, 研究增長區(qū)域上

的反應(yīng)擴(kuò)散問題, 通過線性化和譜分析給出模型產(chǎn)生 Turing 不穩(wěn)定的條件, 再利用數(shù)值模擬驗證理論結(jié)果。 結(jié)果表明擴(kuò)

散系數(shù)的增加有利于 Turing 斑圖的形成, 而區(qū)域增長對 Turing 斑圖形成起破壞作用。

關(guān)鍵詞: 寄生蟲-宿主模型; 增長區(qū)域; Turing 不穩(wěn)定; 數(shù)值模擬; 傳染病

中圖分類號: R181; O175 文獻(xiàn)標(biāo)志碼: A 文章編號: 1001-6600(2022)02-0132-08

近期 H1N1

[1]

、H7N9

[2]

、登革熱[3]等傳染病反復(fù)發(fā)作,新型冠狀病毒肺炎[4] 不斷流行,傳染病的傳播

引起流行病學(xué)專家以及社會廣泛關(guān)注。 傳染病受季節(jié)、溫度等因素影響,其傳播范圍隨時間變化,意味著

受傳染病影響的區(qū)域與時間有關(guān)。 考慮到種群在空間上的遷徙,在模型中研究種群密度的時空演化能更

好地描述傳染病的變化規(guī)律,故而在常微分方程的基礎(chǔ)上引入擴(kuò)散作用,從而得到反應(yīng)擴(kuò)散方程。 此類方

程描繪了生物種群隨時空變化而產(chǎn)生的一系列規(guī)律[5]

。 但是引入擴(kuò)散作用后,原本模型的一系列穩(wěn)定性

質(zhì)就會被打亂,其中 Turing 不穩(wěn)定就是一個例子:由均一的平衡狀態(tài)變?yōu)椴痪粻顟B(tài),這一工作得到了人

們的極大關(guān)注[6-9]

。 由此,為理解在增長的棲息地區(qū)域中,擴(kuò)散對物種生存的影響,本文在增長區(qū)域下討

論一類寄生蟲-宿主傳染病模型的 Turing 不穩(wěn)定。

1 擴(kuò)散的寄生蟲-宿主傳染病模型

2004 年,在 Kermack 等[10]建立的經(jīng)典倉室模型基礎(chǔ)上,Berezovsky 等[11]提出 1 個包含可變?nèi)丝凇⒔佑|

傳播、因病死亡等因素的寄生蟲-宿主模型。

dS

dt

= rN 1-

N

K

( ) -β

SI

N

-(μ+m)S,

dI

dt

= β

SI

N

-(μ+d)I,

ì

?

í

?

?

?

?

(1)

式中:N 表示物種總數(shù),它包含易感人群 S 和染病人群 I,且滿足 N= S+I;r 表示易感人群 S 的內(nèi)稟增長率;

K 表示易感人群 S 的承載能力;β 表示染病人群 I 的接觸傳播率;μ 表示自然死亡率;d 表示病死率;m 表示

非染病者的人均移民率。

Wang 等[12]將模型(1)進(jìn)行無量綱化,并引入擴(kuò)散

收稿日期: 2021-07-09 修回日期: 2021-08-16

基金項目: 國家自然科學(xué)基金(11771381, 11911540464)

通信作者: 林支桂(1965—), 男, 江蘇興化人, 揚州大學(xué)教授, 博士。 E-mail: zglin@yzu.edu.cn

http:∥xuebao.gxnu.edu.cn

?S

?t

-d1ΔS = νRd(S+I)(1-S-I)-R0

SI

S+I

-νS, x∈Ω,t>0,

?I

?t

-d2ΔI =R0

SI

S+I

-I, x∈Ω,t>0,

?S

?υ

= 0,

?I

?υ

= 0, x∈?Ω,t>0,

S(x,0)= S0(x)≥0,I(x,0)= I0(x)≥0, x∈Ω,

ì

?

í

?

?

?

??

?

?

?

?

(2)

式中:R0

=

β

μ+d

為基本再生數(shù);Rd

=

r

μ+m

為人口再生數(shù);ν =

μ+m

μ+d

表示易感染者與感染者平均壽命的比值;Ω

是具有光滑邊界的有界開區(qū)域;υ 表示邊界上的單位外法向量;正常數(shù) d1 、d2 為擴(kuò)散系數(shù);S0(x)、I0(x)為

非負(fù)光滑且不恒為零的函數(shù)。

由偏微分方程理論,易知系統(tǒng)(2)的解存在唯一,且當(dāng) Rd >1 時,系統(tǒng)(2)存在一個無病平衡點 E0

=

(S

Δ

,I

Δ

)= 1-

1

Rd

( ,0) ;當(dāng) Rd >R

?

d : =

R0

+v-1

vR0

,且 R0 >1 時,系統(tǒng)(2) 存在一個染病平衡點 E

? = ( S

?

,I

?

),

式中:

S

? =

vR0Rd

-R0

+1-v

vR

2

0Rd

; I

? = (R0

-1)

vR0Rd

-R0

+1-v

vR

2

0Rd

= (R0

-1)S

?

。

其對應(yīng)的常微分方程的正平衡解所對應(yīng)的 Jacobi 矩陣為

J =

-

R

2

0

+νR0Rd

+νR0

-4R0

-2ν+3

R0

2R0

+2ν-νR0Rd

-3

R0

(R0

-1)

2

R0

1

R0

-1

?

è

?

?

?

?

?

?

?

÷

÷

÷

÷

÷

△=

a11 a12

a21 a22

( ) 。

下面討論系統(tǒng)(2)解的 Turing 失穩(wěn)條件。

定理 1 當(dāng) 1<R0<2 時,如果

Rd>

1

R0(2-R0 )

,

R0

-1

R0Rd

-1

<ν<

(R0

-1)(3-R0 )

R0

-2+R0Rd

,

(d1 a22

+d2 a11 )

2-4d1 d2(a11 a22

-a12 a21 )>0,

且存在某個 δi>0,使得不等式 0<k1<δi<k2 成立,那么系統(tǒng)(2)在 E

? 處 Turing 不穩(wěn)定,式中:δi 是 Ω 中帶有

齊次 Neumann 邊界條件的-Δ 算子的特征值;

k1

=

d1 a22

+d2 a11

- (d1 a22

+d2 a11 )

2-4d1 d2 det(J)

2d1 d2

;

k2

=

d1 a22

+d2 a11

+ (d1 a22

+d2 a11 )

2-4d1 d2 det(J)

2d1 d2

。

證明 易知,當(dāng) R0>1,Rd>

R0

-1+ν

vR0

時,E

?對系統(tǒng)(2)所對應(yīng)的常微分系統(tǒng)穩(wěn)定。

下面設(shè) 0 = δ0<δ1<δ2<…是 Ω 中帶有齊次 Neumann 邊界條件的-Δ 算子的特征值,在 C

1

(?Ω)中 δi 對應(yīng)

的特征空間為 E(δi),設(shè)

X = {U= (φ,ψ)∈C

1

(?Ω)×C

1

(?Ω) | ?ηφ= ?ηψ= 0,x∈?Ω},

{φij,j = 1,2,…,dim E(δi)}是 E(δi)的一組標(biāo)準(zhǔn)正交基,

Xij

= c·φij

| c∈R

2

{ } ,Xi

=?dimE(δI

)

j = 1 Xij,X=?+∞

i = 1Xi。

μ 是 Xi 上的特征值當(dāng)且僅當(dāng) μ 是矩陣

Ai

=

-d1

δi

+a11 a12

a21

-d2

δi

+a22

( )

133

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

的特征值。 此時特征方程為

μ

2+[(d1

+d2 )δi

-a11

-a22 ]μ+[d1 d2

δ

2

i

-(d1 a22

+d2 a11 )δi

+a11 a22

-a12 a21 ] = 0。

顯然,要證明 E

?對于系統(tǒng)(2)不穩(wěn)定,只需要 det(Ai)<0,也就是特征方程至少有 1 個正解。 因為 δi >0,所

以不穩(wěn)定第一個必要條件就是 d1 a22

+d2 a11>0,這等價于 a11>-d1 a22

/ d2>0。

要使 a11

= -

R

2

0

+vR0Rd

+vR0

-4R0

-2v+3

R0

>0,即要 ν<

(R0

-1)(3-R0 )

R0

-2+R0Rd

,R0 <2,Rd >

1

R0(2-R0 )

。 下面證明要

使特征方程存在 1 個正解,則必須滿足

(d1 a22

+d2 a11 )

2-4d1 d2(a11 a22

-a12 a21 )>0。

綜上所述,要使系統(tǒng)(2)在平衡解 E

?處 Turing 不穩(wěn)定,那么 Rd 、R0 、ν 必須滿足

1<R0<2,Rd>

1

R0(2-R0 )

,

R0

-1

R0Rd

-1

<ν<

(R0

-1)(3-R0 )

R0

-2+R0Rd

,

(d1 a22

+d2 a11 )

2-4d1 d2(a11 a22

-a12 a21 )>0。

事實上,此情形還不能保證一定存在 δi>0,使得

h(δi)△=d1 d2

δ

2

i

-(d1 a22

+d2 a11 )δi

+a11 a22

-a12 a21<0。 (3)

所以直接計算

h(δi)= δi(d1

δi

-a11 ) d2

-

M-d1 a22

δi

δi(a11

-d1

δi)

( ) ,

式中 M= det J = a11 a22

-a12 a21 。

進(jìn)一步,假設(shè) a11>0,且 d2 >Di: =

M-d1 a22

δi

δi(a11

-d1

δi)

成立,h(δi) <0。 于是特征方程至少有 1 個正根,故平衡

解 E

?對于系統(tǒng)(2)不穩(wěn)定。 如果對所有 δi>0 都有 di<Di,那么平衡解 E

?對于系統(tǒng)(2)穩(wěn)定。

在滿足上述條件的情況下,方程 h(δi)= 0 存在 2 個正解:

k1

=

d1 a22

+d2 a11

- (d1 a22

+d2 a11 )

2-4d1 d2(a11 a22

-a12 a21 )

2d1 d2

,

k2

=

d1 a22

+d2 a11

+ (d1 a22

+d2 a11 )

2-4d1 d2(a11 a22

-a12 a21 )

2d1 d2

。

從而有結(jié)論:如果存在 δi>0,使不等式 0<k1<δi <k2 成立,那么系統(tǒng)(2)在染病平衡點 E

? 處 Turing 不穩(wěn)定。

證畢。

注 1 假設(shè) Ω = (0,lπ),容易計算出( 0,lπ) 上具有齊次 Neumann 邊界條件的-Δ 算子特征值是

i-1

l

( )

2

,特征函數(shù)是 cos

i-1

l

x,i = 1,2,3,…。 由此得 k1<

i-1

l

( )

2

<k2 ,從而,

l k1 <i-1<l k2 ,

最后只需要滿足(k1

+k2 )-2 k1

k2 >1 / l

2

。

定理 2 假設(shè) E

?是系統(tǒng)(2)的染病平衡點,如果 Ω= (0,lπ),滿足

Rd>

1

R0(2-R0 )

,

R0

-1

R0Rd

-1

<ν<

(R0

-1)(3-R0 )

R0

-2+R0Rd

,

d1(1-R0 )-d2(R

2

0

+vR0Rd

+vR0

-4R0

-2v+3)-2 d1 d2(vR0Rd

-R0

+1-v)(R0

-1) >

R0 d1 d2

l

2

,

那么系統(tǒng)(2)在平衡解 E

?處 Turing 不穩(wěn)定。

故而當(dāng) i→∞ 時,δi→∞ ,并且特征值 δi

{ }

+∞

i = 1會隨區(qū)域 Ω 的大小變化而相應(yīng)連續(xù)發(fā)生變化。 當(dāng)區(qū)域 Ω

134

http:∥xuebao.gxnu.edu.cn

足夠大時,不等式 0<k1<δi<k2 自然滿足,從而定理 1 自然成立;當(dāng)區(qū)域 Ω 充分小時,不等式 0<k1 <δi <k2 不

能滿足,即定理 1 不成立。 因此,空間大小會影響平衡解的穩(wěn)定性。

2 增長區(qū)域上的寄生蟲-宿主傳染病模型

第 1 章通過線性化方法給出了固定區(qū)域下寄生蟲-宿主傳染病模型的 Turing 不穩(wěn)定的條件,結(jié)果表

明,擴(kuò)散會導(dǎo)致此模型的 Turing 不穩(wěn)定。 本章將考慮增長區(qū)域上具 Neumann 邊界條件的寄生蟲-宿主模

型。 記

f(S,I)= νRd(S+I)(1-S-I)-R0

SI

S+I

-νS,

g(S,I)= R0

SI

S+I

-I。

ì

?

í

?

?

?

?

運用 Plaza 等[13]建立的框架,將區(qū)域的增長率和曲率融入傳染病動力學(xué)模型,即設(shè)二維平面 X 嵌入三維空

間 R

3 中,表示為

X(ξ,η,t)= (x(ξ,η,t),y(ξ,η,t),z(ξ,η,t)),

并定義

h1

= Xξ , h2

= Xη , (4)

式中 · 是 Euclidean 范數(shù),下標(biāo)表示偏導(dǎo)數(shù)。 由此,系統(tǒng)(2)轉(zhuǎn)化為

?S

?t

=

d1

h1 h2

h2

h1

( Sξ )

ξ

+

h1

h2

( Sη )

η

é

?

ê

ê

ù

?

ú

ú

-?t(ln(h1 h2 ))S+f(S,I), X(t)∈Ω(t), t>0,

?I

?t

=

d2

h1 h2

h2

h1

I ( ξ )

ξ

+

h1

h2

I ( η )

η

é

?

ê

ê

ù

?

ú

ú

-?t(ln(h1 h2 ))I+g(S,I), X(t)∈Ω(t), t>0,

?S

?υ

= 0,

?I

?υ

= 0, X(t)∈?Ω(t), t>0,

S(X(0),t)= S0(X(0))≥0,I(X(0),t)= I0(X(0))≥0, X(t)∈Ω(0)。

ì

?

í

?

?

?

?

?

?

?

?

?

?

(5)

式中:S0(X(0))和 I0(X(0))是正有界函數(shù);Ω(0)為初始區(qū)域。 由此得到增長區(qū)域下具 Neumann 邊界條

件下的寄生蟲-宿主傳染病模型。

為刻畫區(qū)域演化率的影響,考慮棲息地的演化是各向同性的。 首先探索限制在平面區(qū)域上增長時的

情形,數(shù)學(xué)上表達(dá)式為

X(ξ,η,t)= ρ(t)Y(ξ,η)= ρ(t)[ξ,η,0]

T

,

式中 ρ(t)是增長函數(shù),滿足 ρ(t)在[0,+∞ )上連續(xù)可微,

ρ(0)= 1,

?ρ(t)≥0,lim

t→∞

ρ(t)= ρ∞ >1,lim

t→∞

?ρ(t)= 0。

經(jīng)計算,得

h

2

1

= Xξ

2 = ρ

2

(t),h

2

2

= Xη

2 = ρ

2

(t),h1

/ h2

= 1。

由此,式(5)轉(zhuǎn)化為平面增長區(qū)域下的寄生蟲-宿主傳染病擴(kuò)散模型

?S

?t

=

d1

ρ

2

(t)

(Sξξ

+Sηη )+f(S,I)-

2

?ρ(t)

ρ(t)

S, (ξ,η)∈Ω(0),t>0,

?I

?t

=

d2

ρ

2

(t)

(Iξξ

+Iηη )+g(S,I)-

2

?ρ(t)

ρ(t)

I, (ξ,η)∈Ω(0),t>0,

?S

?υ

= 0,

?I

?υ

= 0, (ξ,η)∈?Ω(0),t>0,

S(ξ,η,0)= S0(ξ,η)≥0,I(ξ,η,0)= I0(ξ,η)≥0, (ξ,η)∈Ω(0)。

ì

?

í

?

?

?

?

?

?

?

?

??

(6)

135

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

3 增長區(qū)域上的 Turing 不穩(wěn)定

本章將通過線性化和譜分析得出正常數(shù)平衡點的穩(wěn)定性,從而給出此情形下所需要的 Turing 不穩(wěn)定

條件。 為具體展現(xiàn)區(qū)域演化對 Turing 不穩(wěn)定的影響,取增長函數(shù) ρ( t) = e

kt

,則式(6) 中的稀釋項即為

?ρ(t) / ρ(t)= k。 為簡化起見,記

?f(S,I)= f(S,I)-2kS, g?(S,I)= g(S,I)-2kI, (7)

從而式(6)變?yōu)?/p>

?S

?t

=

d1

ρ

2

(t)

(Sξξ

+Sηη )+?f(S,I), (ξ,η)∈Ω(0),t>0,

?I

?t

=

d2

ρ

2

(t)

(Iξξ

+Iηη )+g?(S,I), (ξ,η)∈Ω(0),t>0,

?S

?υ

= 0,

?I

?υ

= 0, (ξ,η)∈?Ω(0),t>0,

S(ξ,η,0)= S0(ξ,η)≥0,I(ξ,η,0)= I0(ξ,η)≥0, (ξ,η)∈Ω(0)。

ì

?

í

?

?

?

?

?

?

?

?

??

(8)

式(8)對應(yīng)的動力系統(tǒng)為

dS

dt

=?f(S,I),

dI

dt

= g?(S,I),

S(0)= S0≥0, I(0)= I0≥0。

ì

?

í

??

?

(9)

定理 3 當(dāng) R0>1+2k,Rd >

(R0

+ν-1)(1+2k)

νR0

時,式(9) 對應(yīng)的正平衡點為 E

Δ = ( S

Δ

,I

Δ

),式中:S

Δ =

(1+2k) vR0Rd

-(R0

[ -1+v)(1+2k) ]

vR

2

0Rd

;I

Δ =

(R0

-1-2k)S

Δ

1+2k

,并且正平衡點 E

Δ 局部漸近穩(wěn)定。

證明 為了得到具有每種分支類型(Hopf、Turing 和 Turing-Hopf)下的參數(shù)值的條件,根據(jù)文獻(xiàn)[14],

先做下列計算。 式(7)對應(yīng)的 ODE 系統(tǒng)的零斜率線為

I1(S)= -

νRdR

2

0

(1+2k)

3

S

2+

νRdR0

(1+2k)

2

-

v+2k

1+2k ( ) S 和 I2(S)=

R0 S

1+2k

-S。 (10)

I1 和 I2 相交所滿足的等式為

νRdR

2

0

(1+2k)

3

S

2-

νRdR0

(1+2k)

2

-

ν+2k+R0

1+2k

+1 ( ) S = 0。 (11)

因此,給出條件

νRdR0

1+2k

-v-R0

+1>0 和 R0>1+2k, (12)

從而可以保證式(11)中 S

Δ

>0,I

Δ

>0。 故 E

Δ = (S

Δ

,I

Δ

)是正平衡解。

其次,記(S

Δ

,I

Δ

)處的 Jacobi 矩陣為 J(E

Δ

),這里

J(E

Δ

)=

?f

S?f

I

g?S g?I

( ) , (13)

式中:

?f

S

= -

R

2

0

+νR0Rd

+νR0

-(4+6k)R0

-2(1+2k)v+(1+2k)(2k+3)

R0

;

?f

I

=

-νR0Rd

+(2R0

+2ν-3-2k)(1+2k)

R0

;

g?S

=

(R0

-1-2k)

2

R0

, g?I

=

(1+2k)

R0

2

-1-2k<0。

136

http:∥xuebao.gxnu.edu.cn

考慮其特征方程

μE-J(E

Δ

) = μ

2-(?f

S

+g?I)μ+?f

S g?I

-?f

Ig?S

= 0, (14)

經(jīng)計算,

?f

S

+g?I

= -

(R0

-1-2k)(R0

-1+ν)+νR0Rd

-(R0

-1+ν)(1+2k)

R0

,

?f

S g?I

-?f

Ig?S

= detJ(E

Δ

)=

(R0

-1-2k) (1+2k)(-R0

-ν+1)+vR0Rd

[ ]

R0

。

由條件 R0>1+2k,Rd>

(R0

+v-1)(1+2k)

vR0

,得 ?f

S

+g?I<0 且?f

S g?I

-?f

Ig?S >0,故特征值 μ 具有負(fù)實部。 根據(jù)常微分

方程穩(wěn)定性理論,平衡解(S

Δ

,I

Δ

)對式(9)局部漸近穩(wěn)定。 證畢。

類似定理 1 的證明,通過線性化得特征矩陣,討論特征方程,給出其至少存在 1 個正根的必要條件,即

得系統(tǒng)產(chǎn)生 Turing 不穩(wěn)定的必要條件。

定理 4 式(8)產(chǎn)生 Turing 不穩(wěn)定的必要條件為

f

S

+gI

-4k<0,

f

S gI

-f

IgS

-2k(f

S

+gI)+4k

2

>0,

(d2

f

S

+d1 gI)-2k(d1

+d2 )>0,

- (d2

f

S

+d1 gI)-2k(d1

[ +d2 ) ]

2+4d1 d2

f

S gI

-f

IgS

-2k(f

S

+gI)+4k

2

[ ] <0。

把 f、g 的表達(dá)式代入計算可得

1+2k<R0<2+4k, (15)

Rd>

(1+2k)

2

R0(2-R0

+2k)

, (16)

(R0

-1)(1+2k)

R0Rd

-1-2k

<v<

(R0

-1-2k)(1+2k)+(-R0

+4+4k)(R0

-1)

R0

-2-4k+R0Rd

, (17)

- -

d2 [R

2

0

+νR0Rd

+νR0

-(4+4k)R0

-2(1+2k)ν+(1+2k)(2k+3)]

R0

+d1

(1+2k)

2

R0

( -1-4k)

é

?

ê

ê

ù

?

ú

ú

2

+

4d1 d2

(R0

-1-2k) [R0((1+2k)(4k+3-ν-3R0 )+νR0Rd )+(1+2k)8kν]

R

2

0

é

?

ê

ê

ù

?

ú

ú

<0。 (18)

由此,通過分析發(fā)生 Turing 不穩(wěn)定的必要條件,可以發(fā)現(xiàn):擴(kuò)散系數(shù) d2 越大,模型更容易出現(xiàn)斑圖模

式;而區(qū)域增長快,即 k 大時,上述必要條件不容易滿足,所以對 Turing 斑圖起破壞作用。

4 數(shù)值模擬及解釋

本章將通過數(shù)值模擬來研究固定區(qū)域上 Turing 不穩(wěn)定在此類寄生蟲-宿主模型中的動力學(xué)行為。 在

實際流行病學(xué)中,空間上的擴(kuò)散影響物種初始活動空間范圍[15]

。 為研究模型非平凡解的動態(tài)分布,選取

參數(shù) R0

= 1.3,Rd

= 1.7,ν = 0.3,d1

= 0.1。

對時間和空間都采用有限差分方法,且取時間步長 Δt = 0.01,空間步長 Δx = 1.5,最終時間 t = 1 000,初

始條件是種群初值的隨機(jī)分布。

由圖 1 可以得出,當(dāng)擴(kuò)散系數(shù) d2 越小,即(d1 a22

+d2 a11 )

2-4d1 d2(a11 a22

-a12 a21 )<0 時,種群最后會收斂

到一個平衡解,不發(fā)生 Turing 不穩(wěn)定;當(dāng)擴(kuò)散系數(shù) d2 較大時,Turing 不穩(wěn)定的充分條件滿足,因此種群在

不同位置的值不同,于是常數(shù)解不穩(wěn)定。 圖 2 則是產(chǎn)生 Turing 不穩(wěn)定時種群 S 和 I 所產(chǎn)生的斑圖。

137

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

圖 1 不同擴(kuò)散系數(shù)下種群 S 在位置(25,25)、(50,50)、(75,75)時的密度

Fig. 1 Densities of population S at locations (25, 25), (50, 50) and (75, 75)

with different diffusion coefficients

圖 2 d2

= 5 時種群 S 和 I 所形成的斑圖

Fig. 2 Spatial patterns of population S and I when d2

= 5

參 考 文 獻(xiàn)

[1] NISHIURA H, WILSON N, BAKER M. Estimating the reproduction number of the novel influenza A virus(H1N1) in a

Southern Hemisphere setting: preliminary estimate in New Zealand [ J ]. The New Zealand Medical Journal, 2009,

122(1299): 73-77.

[2] ZHOU J F, WANG D Y, GAO R B, et al. Biological features of novel avian influenza A (H7N9) virus[J]. Nature, 2013,

499(7459): 500-503.

[3] GUZMAN M G, HALSTEAD S B, ARTSOB H, et al. Dengue: a continuing global threat[J]. Nature Reviews Microbiology,

2010, 8(12 Suppl): S7-S16.

[4] 周濤, 劉權(quán)輝, 楊紫陌, 等. 新型冠狀病毒肺炎基本再生數(shù)的初步預(yù)測[ J]. 中國循證醫(yī)學(xué)雜志, 2020, 20( 3):

359-364.

[5] KANG H, RUAN S G. Nonlinear age-structured population models with nonlocal diffusion and nonlocal boundary conditions

[J]. Journal of Differential Equations, 2021, 278: 430-462.

[6] 丁亞君, 張存華. Lengyel-Epstein 反應(yīng)擴(kuò)散系統(tǒng)的穩(wěn)定性和 Turing 不穩(wěn)定性[ J]. 高校應(yīng)用數(shù)學(xué)學(xué)報 A 輯, 2018,

33(3): 272-280.

[7] 林支桂. 數(shù)學(xué)生態(tài)學(xué)導(dǎo)引[M]. 北京: 科學(xué)出版社, 2013: 94-107.

[8] ZU?IGA-GALINDO W A. Reaction-diffusion equations on complex networks and Turing patterns, via p-adic analysis[ J].

138

http:∥xuebao.gxnu.edu.cn

Journal of Mathematical Analysis and Applications, 2020, 491(1): 124239.

[9] 張道祥, 趙李鮮, 胡偉.一類三種群食物鏈模型中交錯擴(kuò)散引起的 Turing 不穩(wěn)定[ J]. 山東大學(xué)學(xué)報(理學(xué)版),

2017, 52(1): 88-97.

[10] KERMACK W O, McKendrick A G. A contribution to the mathematical theory of epidemics[J]. Proceedings of the Royal

Society A: Mathematical, Physical and Engineering Sciences, 1927, 115(772): 700-721.

[11] BEREZOVSKY F, KAREV G, SONG B, et al. A simple epidemic model with surprising dynamics [ J]. Mathematical

Biosciences and Engineering, 2005, 2(1): 133-152.

[12] WANG W M, CAI Y L, WU M J, et al. Complex dynamics of a reaction-diffusion epidemic model[J]. Nonlinear Analysis:

Real World Applications, 2012, 13(5): 2240-2258.

[13] PLAZA R G, SáNCHEZ-GARDU?O F, PADILLA P, et al. The effect of growth and curvature on pattern formation [ J].

Journal of Dynamics and Differential Equations, 2004, 16(4): 1093-1121.

[14] SANCHEZ-GARDUNO F, KRAUSE A L, CASTILLO J A, et al. Turing-Hopf patterns on growing domains: The torus and

the sphere[J]. Journal of Theoretical Biology, 2019, 481: 136-150.

[15] 王瑋明, 蔡永麗. 生物數(shù)學(xué)模型斑圖動力學(xué)[M]. 北京: 科學(xué)出版社, 2020: 8-11.

Turing Instability of a Parasite-host Model on Growing Domains

ZHANG Wanjing, LIN Zhigui

?

(School of Mathematical Science, Yangzhou University, Yangzhou Jiangsu 225002, China)

Abstract: In order to study the influence of growing domain on Turing instability of the parasite-host model, this

paper takes the parasite-host infectious disease model as the main body and analyzes its Turing instability

conditions through linearization and spectral analysis. Numerical simulations are used to verify the theoretical

results. The results show that the increase of the diffusion coefficient is benefitial for formation of Turing pattern,

but regional growth has a destructive effect on Turing pattern formation.

Keywords: parasite-host model; growing domain; Turing instability; numerical simulation; infectious diseases

(責(zé)任編輯 吳佃華)

139

第 40 卷 第 2 期

2022 年 3 月

廣西師范大學(xué)學(xué)報(自然科學(xué)版)

Journal of Guangxi Normal University (Natural Science Edition)

Vol. 40 No. 2

Mar. 2022

DOI: 10.16088 / j.issn.1001-6600.2021041901 http: xuebao.gxnu.edu.cn

劉奇文,李丹,黃小芳,等. 納米金催化甲酸還原磷鉬酸耦合電置換反應(yīng)-共振瑞利散射測定痕量汞[ J]. 廣西師范大學(xué)學(xué)報(自然科學(xué)版),

2022, 40(2): 140-148. LIU Q W, LI D, HUANG X F, et al. A new strategy for the determination of trace mercury by resonance Rayleigh scattering

method based on nano-gold catalytic amplification and galvanic replacement reaction-phosphomolybdic acid[J]. Journal of Guangxi Normal University

(Natural Science Edition), 2022, 40(2): 140-148.

納米金催化甲酸還原磷鉬酸耦合電置換反應(yīng)共振瑞利散射測定痕量汞

劉奇文1,2

, 李 丹1,2

, 黃小芳1,2

, 梁愛惠1,2?

, 蔣治良1,2?

(1. 珍稀瀕危動植物生態(tài)與環(huán)境保護(hù)教育部重點實驗室(廣西師范大學(xué)), 廣西 桂林 541006;

2. 廣西環(huán)境污染控制理論與技術(shù)重點實驗室(廣西師范大學(xué)), 廣西 桂林 541006)

摘 要: 在 pH= 3.1 的 HCOOH-HCOONa 緩沖液中, 磷鉬酸粒子在 450 nm 處產(chǎn)生一個共振瑞利散射(RRS)峰。 金納米粒

子(AuNPs)可催化磷鉬酸-甲酸反應(yīng)生成磷鉬藍(lán), 使得 450 nm 處磷鉬酸的 RRS 強度線性降低。 Hg

2+可與 AuNPs 發(fā)生電置

換反應(yīng), 從而抑制 AuNPs 的催化作用, RRS 峰增強。 在 2.5×10

-4

~ 3.5 μmol / L, 隨著 Hg

2+濃度的增加, AuNPs 的催化作

用逐漸減弱, 反應(yīng)液的顏色逐漸從藍(lán)色變?yōu)闊o色, 體系在 450 nm 處的 RRS 峰值(ΔI) 線性增高, 其線性方程為 ΔI =

0.32C+46.1, 檢出限為 0.18 nmol / L。 該法用于廢水中 Hg

2+的檢測, 結(jié)果令人滿意。

關(guān)鍵詞: 汞; 納米金催化; 靜電取代; 磷鉬藍(lán); 共振瑞利散射

中圖分類號: O657.3 文獻(xiàn)標(biāo)志碼: A 文章編號: 1001-6600(2022)02-0140-09

電置換反應(yīng)(GR)指一種較活潑金屬被一種較不活潑的金屬侵蝕,反應(yīng)較快,常在幾分鐘內(nèi)發(fā)生[1]

。

由于納米表面的高催化活性,納米粒子表面的 GR 反應(yīng)更快。 GR 反應(yīng)由于具有高度的可調(diào)性而引起人們

的極大興趣,可用于研究金屬納米結(jié)構(gòu)表面錯綜復(fù)雜的合金化和去合金化反應(yīng)[2]

;GR 反應(yīng)還提供一種非

常簡易、多功能路徑制備可控空心和多孔壁納米結(jié)構(gòu)材料方法[3]

。 Sun 等[4] 采用 Ag 作為模板金屬,和

Au

3+發(fā)生氧化還原反應(yīng),制備了 Au 中空結(jié)構(gòu),并對其反應(yīng)機(jī)理進(jìn)行了解釋。 Liu 等[5] 利用超小 ssDNA-模

板化銀納米簇與 Cu(II)發(fā)生拮抗電置換反應(yīng)制備了 Ag-Cu 合金,并用光散射技術(shù)檢測反應(yīng)的發(fā)生。 Bi

等[6]通過 Ag / AgCl 核殼納米線和 H2PtCl

6發(fā)生電置換反應(yīng)得到較高電催化活性的納米管。 近年來,GR 反

應(yīng)在納米分析檢測技術(shù)中的應(yīng)用也較活躍[7-10]

。 Netzer 等[7] 將 GR 反應(yīng)與銀納米線 Langmuir-Blodgett 膜

技術(shù)有機(jī)結(jié)合成功地制備了高活性 SERS 基底,基于 1 394 cm

-1處的靈敏拉曼峰可用于檢測 8 nmol / L 4-氨

基硫酚。 共振瑞利散射(RRS)不僅是一種簡便靈敏的分子光譜分析技術(shù),也是一種研究納米微粒反應(yīng)的

靈敏光譜技術(shù)[11-13]

,納米催化放大是提高分析方法靈敏度的重要途徑之一[14]

。 RRS 與納米催化技術(shù)結(jié)合

構(gòu)建了重金屬環(huán)境污染物 RRS 分析新方法。 Wang 等[15]以摻鈀共價有機(jī)骨架(TpPaPd)為催化劑,催化次

磷酸鈉還原 Ni(II)形成 Ni-P 合金,該產(chǎn)物具有 RRS 信號,將納米催化反應(yīng)與基于 Pb

2+的 DNA 酶反應(yīng)相

結(jié)合,建立了檢測 0.001~0.100 nmol / L Pb

2+的 RRS 分析方法,其檢出限為 0.4 pmol / L。 Zhang 等[16]利用摻

金碳點(CDAu )對 AgNO3

-葡萄糖反應(yīng)具有強催化,其產(chǎn)物銀納米粒子(AgNPs)具有 SPR 效應(yīng)。 結(jié)合 Apt

特異性反應(yīng)與納米催化放大信號,建立了 RRS 法檢測水中痕量 As

3+

,線性范圍為 0.025~0.750 μg / L ,檢出

限為 0.01 μg / L。 但尚未見基于納米金催化甲酸還原磷鉬酸耦合金納米粒子(AuNPs)-汞離子電置換反

應(yīng),RRS 測定痕量汞的研究報道。

收稿日期: 2021-04-19 修回日期: 2021-04-29

基金項目: 國家自然科學(xué)基金(21767004)

通信作者: 梁愛惠(1965—), 女, 廣西岑溪人, 廣西師范大學(xué)研究員。 E-mail: ahliang2008@163.com

蔣治良(1965—), 男, 廣西全州人, 廣西師范大學(xué)教授, 博導(dǎo)。 E-mail: zljiang@gxnu.edu.cn

http:∥xuebao.gxnu.edu.cn

汞離子是重金屬中最重要的污染物之一,其不可生物降解,容易通過食物鏈在魚類和貝類中積累[17]

。

此外,汞離子在自然環(huán)境中可以轉(zhuǎn)化為高毒性甲基汞,比 Hg

2+更容易穿過細(xì)胞膜,對腦組織、腎臟、神經(jīng)系

統(tǒng)、泌尿系統(tǒng)和內(nèi)分泌系統(tǒng)造成嚴(yán)重?fù)p害[18]

。 鑒于 Hg

2+在環(huán)境和生物系統(tǒng)中的毒性,因此,開發(fā)一種快

速、高效、靈敏檢測痕量 Hg

2+的方法至關(guān)重要。 最近,已有學(xué)者給出了一些用于檢測 Hg

2+的傳感平臺,包

括比色法[19]

、電化學(xué)法[20]

、熒光法[21]和 RRS 法[22]

。 在這些方法中,RRS 法因其高靈敏度、高選擇性、操

作簡單和響應(yīng)速度快等獨特優(yōu)勢而備受關(guān)注。 Zhang 等[23] 合成了二維納米復(fù)合物 rGO/ PEI/ Pd,借助其

過氧化物酶活性,開發(fā)了一種快速、高度選擇性和超痕量肉眼比色法檢測水溶液中 Hg

2+的通用策略。 主

要基于在汞離子存在下,rGO/ PEI/ Pd 可以促進(jìn) 3,3’,5,5’-四甲基聯(lián)苯胺(TMB)的有效氧化,使其顏色變

為肉眼和吸收光譜法可檢測到的深藍(lán)色,最低檢測濃度為 0.39 nmol / L Hg

2+

。 Tan 等[24]制備并表征了電化

學(xué)衍生的還原型氧化石墨烯化學(xué)阻抗傳感器,該傳感器對 Hg

2+表現(xiàn)出選擇性響應(yīng),將其用于水樣中 Hg

2+

的檢測,最低檢測濃度為 0.5 nmol / L。 Yu 等[25]以金華佛手柑為碳源,水熱法制備了具有較高光致發(fā)光性

水溶性熒光碳點,基于 Hg

2+會導(dǎo)致碳點發(fā)生熒光猝滅,建立一個熒光檢測 0.01~100 μmol / L Hg

2+的分析方

法,其檢出限為 5.5 nmol / L。 Ngernpimai 等[26]利用 Hg

2+可誘導(dǎo)單鏈 DNA 的構(gòu)象變化,進(jìn)一步導(dǎo)致金納米

棒(GNR)聚集,引起 RRS 強度的增強,建立共振瑞利散射法檢測 Hg

2+的策略,其檢出限為 0.23 nmol / L。

但上述方法有的靈敏度欠佳,有的分析過程復(fù)雜。

本文利用 AuNPs 催化磷鉬酸-甲酸反應(yīng),結(jié)合 AuNPs 與 Hg

2+的 GR 反應(yīng),以磷鉬酸顆粒為 RRS 指示

劑,建立一種簡便、靈敏的檢測 Hg

2+的 RRS 分析新方法,檢出限為 0.18 nmol / L。

1 實驗部分

1.1 主要儀器和試劑

日立 F-7000 熒光分光光度計(日立高新技術(shù)公司);TU-1901 型雙光束紫外可見分光光度計(北京普

析通用儀器有限責(zé)任公司);SYZ-550 型石英亞沸蒸餾水器(江蘇晶玻儀器廠);KQ3200DB 型數(shù)控超聲波

清洗器(昆山市超聲儀器有限公司,功率 150 W,工作頻率 40 kHz);pH 計(梅特勒-托利多儀器上海有限公

司);Nano-23 型納米粒度與電位分析儀(英國 Malvern 公司);FEI Quanta 200-場發(fā)射掃描電子顯微鏡(FEI

公司)。

5.0 mmol / L 磷鉬酸(中國菱湖化工試劑廠,分析純);0.1 mol / L 甲酸、質(zhì)量分?jǐn)?shù) 0.05% NaBH4(汕頭市

西隴化工廠);0.1 mol / L NaOH(天津市百世化工有限公司);質(zhì)量分?jǐn)?shù) 1%檸檬酸三鈉(汕頭市西隴化工

廠);0.01 mol / L HAuCl

4(國藥集團(tuán)化學(xué)試劑公司);0.01 mol / L HgSO4(汕頭市西隴化工廠)儲備液;pH3.1

HCOOH-HCOONa 緩沖液:10 mL 離心管中加入 2.0 mL 0.1 mol / L NaOH 溶液和 8.0 mL 0.1 mol / L HCOOH

溶液,充分混合,得濃度為 5.0 mmol / L 的 HCOOH-HCOONa 緩沖溶液(濃度以醋酸根計)。

金納米粒子(AuNPs):常溫下,量取 40.0 mL 二次水于錐形瓶,在攪拌下加入 0.5 mL 質(zhì)量分?jǐn)?shù) 1%

HAuCl

4 溶液和 3.5 mL 質(zhì)量分?jǐn)?shù) 1%檸檬酸鈉水溶液,充分混合后,在攪拌下緩慢滴加 4.0 mL 質(zhì)量分?jǐn)?shù)

0.05%的 NaBH4 溶液,加完后繼續(xù)攪拌 10 min,最后定容至 50.0 mL,其濃度為 58 mg / L AuNPs。

1.2 實驗方法

在 5.0 mL 的刻度試管中,依次移取 170 μL 5 mmol / L 磷鉬酸、500 μL pH = 3.1 的 HCOOH-HCOONa 緩

沖液、50 μL 58 mg / L 的 AuNPs 及一定濃度的 Hg

2+

,定容至 2.0 mL,混勻。 80 ℃ 水浴 35 min 后,冰水終止

反應(yīng)。 取溶液于石英皿內(nèi),在 400 V,激發(fā)和發(fā)射狹縫 10 nm,用熒光分光光度計掃描,獲得共振瑞利散射

光譜;不加 Hg

2+做空白,測定溶液 450 nm 處的 RRS 峰值,計算 ΔI = I450 nm

-(I450 nm )0 。

2 結(jié)果與討論

2.1 分析原理

在實驗條件下,磷鉬酸具有共振瑞利散射信號,可作為該分析反應(yīng)的指示組分。 AuNPs 可以催化磷

鉬酸與甲酸反應(yīng),生成磷鉬藍(lán)而吸收瑞利散射光,使得體系的 RRS 信號降低。 Hg

2+與 AuNP 發(fā)生電置換反

141

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

應(yīng)生成穩(wěn)定的金汞齊(AuHg)而覆蓋在其表面,抑制了 AuNPs 的催化活性,使得生成的磷鉬藍(lán)降低,導(dǎo)致

體系的 RRS 信號增強,據(jù)此可以建立測定 Hg

2+的共振瑞利散射新方法(圖 1)。

圖 1 納米金催化甲酸還原磷鉬酸耦合電置換反應(yīng)-RRS 測定汞原理

Fig. 1 Analytical principle of the RRS measurement of Hg

2+

by coupling AuNP-phosphomolybdic

acid-formic acid catalytic amplification with GR reaction of AuNP-Hg(Ⅱ)

2.2 RRS 光譜

在 pH = 3.1 的 HCOOH-HCOONa 緩沖溶液及常溫下,磷鉬酸與甲酸不反應(yīng),即使在沸水加熱條件下反

應(yīng)也極難進(jìn)行。 AuNPs 可催化磷鉬酸-甲酸反應(yīng)生成磷鉬藍(lán),隨著 AuNPs 濃度的增加,反應(yīng)液的顏色逐漸

從無色變?yōu)樗{(lán)色,體系在 450 nm 處的共振瑞利散射峰逐漸降低(圖 2A),系生成藍(lán)色的磷鉬藍(lán)對 RRS 光

吸收所致。 Hg

2+可與 Au 發(fā)生電置換反應(yīng),生成的金汞齊覆蓋在 AuNPs 表面抑制 AuNPs 的催化作用,在

2.5×10

-4

~3.5 μmol / L 隨著 Hg

2+濃度的增加,AuNPs 的催化作用逐漸減弱,反應(yīng)液的顏色逐漸從藍(lán)色變?yōu)?/p>

無色,體系在 450 nm 處的共振瑞利散射峰線性升高(圖 2B)。

a: 0.425 mmol / L 磷鉬酸+5.0 mmol / L HCOOHHCOONa; b: a+0.29 mg / L AuNPs; c: a+0.58

mg / L AuNPs; d: a+ 1.16 mg / L AuNPs; e: a+

1.45 mg / L AuNPs; f: a+2.03 mg / L AuNPs; g:

a+2.61 mg / L AuNPs

a: 0. 425 mmol / L 磷鉬酸 + 5. 0 mmol / L HCOOHHCOONa+ 1. 45 mg / L AuNPs; b: a + 0. 25 nmol / L

Hg

2+

; c: a + 5. 0 nmol / L Hg

2+

; d: a + 1. 0 μmol / L

Hg

2+

; e: a + 1. 5 μmol / L Hg

2+

; f: a + 2. 5 μmol / L

Hg

2+

; g: a+3.5 μmol / L Hg

2+

圖 2 磷鉬酸-HCOOH-HCOONa-AuNPs-Hg

2+體系 RRS 光譜

Fig. 2 RRS spectra of phosphomolybdic acid-HCOOH-HCOONa-AuNPs-Hg

2+

system

142

http:∥xuebao.gxnu.edu.cn

2.3 紫外吸收光譜

由圖 3A 可知,AuNPs 可催化磷鉬酸-甲酸反應(yīng),其吸光度隨著 AuNPs 濃度的增加線性增強。 隨著

Hg

2+濃度的增加,AuNPs 對磷鉬酸-甲酸反應(yīng)的催化作用逐漸減弱,生成的磷鉬藍(lán)逐漸減少,在 1.0×10

-2

~

3.5 μmol / L 范圍內(nèi),隨著 Hg

2+濃度的增大,體系在 700 nm 處的 Abs 信號線性降低(圖 3B),此結(jié)論與 RRS

光譜結(jié)果一致。 雖然紫外吸收光譜也可用于汞離子測定,但不及 RRS 靈敏。

2.4 掃描電鏡(SEM)

圖 4A 是磷鉬酸的掃描電鏡圖,由圖可見其顆粒均勻分散,這是由于磷鉬酸分子具有較強的疏水性并

聚集形成粒徑約為 8 μm 的顆粒。 隨著 Hg

2+的加入,AuNP-磷鉬酸-甲酸鈉催化反應(yīng)逐漸被抑制,但磷鉬酸

顆粒的形狀變化不大(圖 4B)。

a: 0. 425 mmol / L 磷鉬酸+ 5. 0 mmol / L HCOOHHCOONa; b: a+0.29 mg / L AuNPs; c: a+0.58 mg /

L AuNPs; d: a+0.87 mg / L AuNPs; e: a+1.16 mg /

L AuNPs; f: a+1.45 mg / L AuNPs; g: a+2.03 mg /

L AuNPs; h: a+2.32 mg / L AuNPs; i: a+2.61 mg /

L AuNPs

a: 0.425 mmol / L 磷鉬酸+ 5. 0 mmol / L HCOOHHCOONa+ 1. 45 mg / L AuNPs; b: a + 10 nmol / L

Hg

2+

; c: a+ 0.2 μmol / L Hg

2+

; d: a + 1. 0 mmol / L

Hg

2+

; e: a + 1. 5 mmol / L Hg

2+

; f: a + 2. 5 mmol / L

Hg

2+

; g: a+3.5 mmol / L Hg

2+

圖 3 磷鉬酸-AuNPs-HCOOH-HCOONa-Hg

2+體系吸收光譜

Fig. 3 Absorption spectra of phosphomolybdic acid-AuNPs-HCOOH-HCOONa-Hg

2+

system

2.5 激光散射

按實驗方法制得反應(yīng)液,用納米粒度與 Zeta 電位分析儀記錄其粒度。 從實驗結(jié)果可以看出,無 Hg

2+

存在時,磷鉬酸-甲酸鈉反應(yīng)在 AuNPs 的催化作用下反應(yīng)充分,體系中的磷鉬酸少量聚集,顆粒粒度較小

(圖 5a),其粒徑分布在 100~500 nm;隨著 Hg

2+的加入,磷鉬酸-甲酸鈉反應(yīng)逐漸被抑制,體系的顆粒分布

在 140~400 nm(圖 5b)。 此激光散射粒徑與 SEM 粒徑不一致的主要原因是電鏡制備樣品需要脫水,此過

程會導(dǎo)致顆粒聚集。

2.6 分析條件優(yōu)化

考察磷鉬酸濃度對體系 RRS 信號的影響,加入 5.0 mmol / L HCOOH-HCOONa 緩沖溶液,1.45 mg / L

AuNPs,1.0 μmol / L HgSO4 和不同濃度的磷鉬酸溶液,其余按照 1. 2 節(jié)的實驗步驟操作,當(dāng)加入 0. 425

143

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

圖 4 體系掃描電鏡

Fig. 4 Scanning electron micrograph of the system

a: 0.425 mmol / L 磷鉬酸+ 5.0 mmol / L HCOOH-HCOONa+1.45 mg / L AuNPs; b: a+0.35 μmol / L Hg

2+

圖 5 粒度分析

Fig. 5 Particle size analysis

mmol / L 磷鉬酸時 ΔI 達(dá)到最大,故選用 0.425 mmol / L 磷鉬酸(圖 6A)。 隨著 AuNPs 濃度的增加,其催化作

用逐漸增強,當(dāng)加入 1. 45 mg / L AuNPs 時 ΔI 達(dá)到最大,故選用 1. 45 mg / L AuNPs ( 圖 6B)。 考察了

HCOOH-HCOONa 緩沖溶液 pH 對體系 RRS 信號的影響,當(dāng)加入 pH = 3.1 的 HCOOH-HCOONa 緩沖溶液時

ΔI 達(dá)到最大,故選用 pH = 3.1 的 HCOOH-HCOONa 緩沖溶液(圖 6C)。 HCOOH-HCOONa 緩沖溶液的濃度

同樣會對體系 RRS 信號產(chǎn)生影響,當(dāng)加入 5.0 mmol / L HCOOH-HCOONa 緩沖溶液時 ΔI 達(dá)到最大,故選用

5.0 mmol / L 的 HCOOH-HCOONa 緩沖溶液(圖 6D)。 考察反應(yīng)溫度對體系 RRS 信號的影響,當(dāng)反應(yīng)溫度

為 80 ℃時 ΔI 達(dá)到最大,故選用 80 ℃作為反應(yīng)溫度(圖 6E)。 考察反應(yīng)時間對體系 RRS 信號的影響,當(dāng)

反應(yīng)時間為 35 min 時 ΔI 達(dá)到最大,故選用 35 min 作為反應(yīng)時間(圖 6 F)。

2.7 工作曲線

在最佳實驗條件下,AuNPs 表現(xiàn)出對于磷鉬酸-甲酸鈉 RRS 體系的催化作用,在 0.29 ~ 2. 61 mg / L

AuNPs,450 nm 處的 RRS 強度變化 ΔI 與 AuNPs 濃度呈線性關(guān)系, AuNPs 的催化線性方程為 y =

288.39CAuNPs

+26.9,線性相關(guān)系數(shù) R

2 為 0.994(圖 2A)。 基于磷鉬酸-AuNPs-甲酸鈉體系表現(xiàn)出對 Hg

2+濃

度的響應(yīng),建立一個測定 Hg

2+的共振瑞利散射方法,在 2.5×10

-4

~3.5 μmol / L Hg

2+

,450 nm 處的 RRS 強度

變化 ΔI 與 Hg

2+濃度呈線性關(guān)系,線性方程為 y = 0.32CHg

2+ +46.1,線性相關(guān)系數(shù) R

2 為 0.994 6(圖 2B)。 對

144

http:∥xuebao.gxnu.edu.cn

圖 6 磷鉬酸-AuNPs-HCOOH-HCOONa-Hg

2+體系反應(yīng)條件優(yōu)化

Fig. 6 Optimization of the conditions of phosphomolybdic acid-AuNPs-HCOOH-HCOONa-Hg

2+

system

于磷鉬酸-甲酸鈉 UV 體系,在 0.29~2.61 mg / L AuNP,700 nm 處的 UV 強度變化 ΔA 與 AuNPs 濃度呈線性

關(guān)系,線性方程為 y = 0.241CAuNPs

+0.003,線性相關(guān)系數(shù) R

2 為 0.994 8(圖 3A)。 已知磷鉬酸-AuNPs-甲酸鈉

體系可用于檢測 Hg

2+

,在 1.0×10

-2

~3.5 μmol / L Hg

2+

,700 nm 處的 UV 強度變化 ΔA 與 Hg

2+濃度呈線性關(guān)

系,線性方程為 y = 0.000 1CHg

2+ +0.028,線性相關(guān)系數(shù) R

2 為 0.986 9(圖 3B)。 表 1 是本方法與部分文獻(xiàn)報

道方法的對比。

表 1 本法與已報道方法分析特性比較

Tab. 1 Comparison of analysis characteristics between this method and the reported method

方法 線性范圍/ (nmol·L

-1

) 檢出限/ (nmol·L

-1

) 注解 參考文獻(xiàn)

熒光法 50~ 20 000 13.2 高靈敏但底物合成繁瑣 [27]

伏安法 2.4~ 220 0.8 選擇性好但電極修飾繁雜 [28]

比色法 10

4

~ 10

5

3.6 穩(wěn)定性好但檢測范圍不寬 [29]

比色法 5 000~ 35 000 2 370 快速方便但檢測靈敏度低 [30]

本法 0.25~ 3 500 0.18 靈敏、簡便、快速 本文

2.8 干擾試驗

按實驗方法考察了共存離子對磷鉬酸-AuNPs-HCOOH-HCOONa-Hg

2+ RRS 測定 0.5 μmol / L Hg

2+的干

擾情況。 由表 2 表明,5 μmol / L Mg

2+

、SO

2-

4 、Fe

3+

、Pb

2+

;25 μmol / L Ca

2+

、Br

-

、Cr

3+

、Al

3+

、HSA、血紅蛋白、血

漿蛋白;50 μmol / L IO

-

3 、Na

+

、Zn

2+

、Ba

2+

、NH

+

4 、NO

-

3 、Co

2+

、Cu

2+

、Mn

2+

、BSA、鯡魚精 DNA、牛白蛋白、核糖核酸

酶不干擾測定,表明本法具有較好的選擇性。

2.9 樣品分析

從污水廠取得 3 份工業(yè)廢水,過濾處理后分別取樣品 200 μL,按 1.2 節(jié)試驗方法進(jìn)行測定,檢測結(jié)果

見表 3。 水樣中 Hg

2+含量為 1.9 ~ 5.1 μmol / L,相對標(biāo)準(zhǔn)偏差(RSD) 為 2.6% ~ 7.0%,回收率為 93.7% ~

98.4%,具有較好的回收率和重現(xiàn)性。

145

廣西師范大學(xué)學(xué)報(自然科學(xué)版),2022,40(2)

表 2 干擾離子對體系的影響

Tab. 2 Influence of interfering ions on the system

干擾離子 相對倍數(shù) 相對誤差/ % 干擾離子 相對倍數(shù) 相對誤差/ %

Mg

2+

10 5.4 SO

2-

4 10 -7.0

Fe

3+

10 1.6 Pb

2+

10 3.5

Ca

2+

50 9.5 Br

-

50 4.5

Cr

3+

50 -8.8 Al

3+

50 2.8

HSA 50 9.4 血紅蛋白 50 -5.2

血漿蛋白 50 -3.4 IO

-

3 100 8.8

Na

+

100 5.6 Zn

2+

100 4.6

Ba

2+

100 1.9 NH

+

4 100 9.1

NO

-

3 100 2.9 Co

2+

100 4.4

Cu

2+

100 -8.8 Mn

2+

100 -8.7

BSA 100 -3.1 鯡魚精 DNA 100 3.3

牛白蛋白 100 -0.2 核糖核酸酶 100 3.1

表 3 樣品的 RRS 測定結(jié)果

Tab. 3 RRS measurement results of the samples

樣品

平均值

/ (nmol·L

-1

)(n = 5)

加入 Hg

2+

/ (nmol·L

-1

)

測得值/

(nmol·L

-1

)

回收率/ % RSD/ %

Hg

2+含量/

(μmol·L

-1

)

水樣 1 192.6 200 386.2 96.8 7.0 1.9

水樣 2 261.7 200 449.1 93.7 4.9 2.6

水樣 3 506.9 200 703.7 98.4 2.6 5.1

3 結(jié)論

本文以磷鉬酸顆粒具有 RRS 效應(yīng)作為分析反應(yīng)的指示劑,借助 AuNPs 的催化放大作用以及 AuNP 與

Hg

2+的電置換反應(yīng),建立了一個簡單、靈敏檢測痕量 Hg

2+ 的共振瑞利散射分析方法,線性范圍為 0.25 ~

3 500 nmol / L,檢測限為 0.18 nmol / L。

參 考 文 獻(xiàn)

[1] ZHANG C L, LUO L, LUO J, et al. A process-analysis microsystem based on density gradient centrifugation and its application in

the study of the galvanic replacement mechanism of Ag nanoplates with HAuCl

4[ J]. Chemical Communications, 2012, 48

(58): 7241-7243.

[2] JAIN P K, HUANG X H, EL-SAYED I H, et al. Noble metals on the nanoscale: optical and photothermal properties and

some applications in imaging, sensing, biology, and medicine[J]. Accounts of Chemical Research, 2008, 41: 1578-1586.

[3] LU X M, CHEN J Y, SKRABALAK S E, et al. Galvanic replacement reaction: a simple and powerful route to hollow and

porous metal nanostructures[J]. Proceedings of the Institution of Mechanical Engineers, Part N: Journal of Nanomaterials,

Nanoengineering and Nanosystems, 2007, 221(1): 1-16.

146