鼠疫（plague）是一种由鼠疫耶尔森菌（Yersinia pestis，鼠疫菌）引发的自然疫源性疾病，因疫源地宿主、媒介、病原的差异而呈现丰富的多态性^[1-4]。目前，在中国的鼠疫自然疫源地中，其中4处山地疫源地位于新疆维吾尔自治区（新疆）^[5]。鼠疫在自然疫源地中的流行基于鼠疫菌及其宿主和媒介构成的生态系统，而蚤类被认为是疫源地鼠疫持续流行的关键^[6]。以往研究表明，新疆各山地鼠疫疫源地蚤类群落的种类、数量、分布均存在差异^[7]，而生态环境对传染性疾病流行的影响一直备受关注，如生态因素对汉坦病毒（Hantavirus）、克氏锥虫（Trypanosoma curzi）、苍白密螺旋体（Treponema pallidum）流行的影响^[8-9]。为深入了解鼠疫流行规律，本文从新疆各山地鼠疫疫源地蚤类群落结构入手，选择鼠疫主要宿主动物血清抗体阳性率水平代表鼠疫流行强度，分析疫源地蚤的数量、密度、蚤群落物种的多样性、丰富度、均匀度和优势度在4类疫源地分布的差异及不同蚤生态结构对鼠疫流行的影响作用，为开展鼠疫防治、监测提供更多的科学依据。

数据来源于新疆疾病预防控制中心2001－2013年的鼠疫监测结果，并已在《中国地方病防治杂志》《新疆疾病预防通报》等杂志公开发表。本研究数据选择新疆4类山地鼠疫疫源地鼠疫监测结果，4类疫源地分别为北天山灰旱獭（Marmota baibacina）和长尾黄鼠（Spermophilus undulatus）、南天山灰旱獭、帕米尔高原阿赖山红旱獭（M. caudata）和昆仑山喜马拉雅旱獭（M. himalayana）鼠疫自然疫源地，研究对象为各疫源地主要鼠疫宿主动物及其体蚤。

疫源地动物鼠疫监测按照《全国鼠疫监测方案》执行。旱獭采用铁丝套法捕获，长尾黄鼠采用弓形夹法捕获，采集心脏血液用于血清学检测，选择间接血凝试验检测鼠疫菌F1抗体，并计算鼠疫菌血清抗体阳性率；单只啮齿动物梳检媒介蚤，并进行形态学分类鉴定，按同一时间、地点、种类进行分组送检，计算蚤指数和染蚤率，计算4类疫源地蚤群落多样性、丰富度、均匀度和优势度水平。蚤指数、染蚤率等宿主体蚤群落指标按以下公式计算：

多样性指数采用Shannon⁃Wienner（1948）多样性指数，均匀度采用Pielon（1966）公式，多样性指数（H'）、均匀度（E）、和优势度（D）按下列公式计算：

式中，P_i为物种i个体数占群落中总个体数的比例，S为群落中物种种类数（丰富度），N_i为物种i的个体数，N为群落中所有物种的总个体数。

采用SPSS 25.0软件进行数据整理和统计分析。山地鼠疫疫源地宿主体蚤群落蚤指数、染蚤率、多样性、优势度、均匀度和丰富度两组间差异性分析选择Wilcoxon秩和检验；多组间差异性分析选择Kruskal-Wallis H秩和检验。应用χ²检验各疫源地宿主抗体阳性率差异性，选择Spearman秩相关检验上述疫源地宿主抗体阳性率与蚤类群落指标间的相关性。P < 0.05为差异有统计学意义。

应用SPSS 25.0软件对相关性分析存在统计学意义的指标进行主成分分析（principal component analysis，PCA）。主成分分析旋转方法选择最大方差法^[10]，主成分数量的确定采用均数法。

首先，将主成分作为自变量，利用R 4.3.0软件将宿主抗体阳性率作为应变量构建广义相加模型（generalized additive model，GAM）。模型中，总体样本70%划分为训练集，30%划分为测试集^[11]。模型公式如下：

式中，g（）是连接函数，E_Y表示目标变量的期望值，α为截距，f₁（X₁）到f_p（X_p）是平滑函数，表示主成分的非线性关系，ε表示残差。评估模型拟合效果的指标采用决定系数（R²）、均方根误差（root mean square error，RMSE）和平均绝对误差（mean absolute error，MAE）。

再将主成分作为自变量，构建随机森林（random forest，RF）模型，用于衡量每个主成分对模型拟合的贡献程度^[12]，检验水准α=0.05。

北天山疫源地主要宿主动物为长尾黄鼠和灰旱獭，南天山疫源地为灰旱獭，帕米尔高原疫源地为阿赖山红旱獭，昆仑山疫源地为喜马拉雅旱獭。2001－2013年，4类疫源地记录捕获鼠疫宿主83 017只，其中54 023只染蚤，平均染蚤率为65.07%，其中1 836份鼠血清鼠疫菌抗体检测呈阳性，平均阳性率为2.21%。发现蚤9种，分别为谢氏山蚤（Oropsylla silantiewi）、斧形盖蚤（Callopsylla dolabris）、矩凹黄鼠蚤（Citellophilus lebedewi）、方形黄鼠蚤（Ci. tesquorum）、似升额蚤（Frontopsylla elatoides）、宽新蚤（Neopsylla mana）、腹窦纤蚤（Rhadinopsylla li）、人蚤（Pulex irritans）和阿州山蚤（O. alaskensis），共计338 008匹，平均蚤指数为4.07。其中北天山主要蚤类为谢氏山蚤和方形黄鼠蚤，分别占该疫源地捕获蚤类样本的57.50%和11.86%；南天山主要蚤类为谢氏山蚤，占该疫源地捕获蚤类样本的88.92%；帕米尔高原疫源地主要蚤类为矩凹黄鼠蚤，占该疫源地捕获蚤类样本的68.41%；昆仑山疫源地主要蚤类为谢氏山蚤，占该疫源地捕获蚤类样本总数的75.43%。根据蚤群落中蚤类构成情况，计算各疫源地宿主体蚤群落结构指标。各疫源地宿主血清抗体阳性率及各蚤群落生态学指标见表 1。

表 1 新疆维吾尔自治区4类山地鼠疫疫源地宿主体蚤群落生态学指标[M（P₂₅，P₇₅）] Table 1 Ecological indicators of host flea community in 4 types of mountain plague foci in Xinjiang Uygur Autonomous Region

表选项

蚤指数、染蚤率、蚤类多样性、丰富度、蚤群落均匀度和优势度在4类疫源地分布的差异均有统计学意义（H=30.146，P < 0.001；H=31.602，P < 0.001；H=33.443，P < 0.001；H=34.532，P < 0.001；H=17.900，P < 0.001；H=26.686，P < 0.001）。北天山鼠疫疫源地宿主动物染蚤率、蚤指数、蚤类多样性和丰富度水平最高，中位数分别为73.55%、4.67、1.40和7.00（表 1），各疫源地间差异分析结果显示，北天山疫源地宿主体蚤指数、染蚤率、蚤类多样性和丰富度与南天山和昆仑山疫源地差异存在统计学意义，北天山疫源地与帕米尔疫源地宿主染蚤率和蚤类丰富度差异存在统计学意义。南天山疫源地与昆仑山疫源地宿主染蚤率差异存在统计学意义。南天山疫源地蚤群落均匀度最低，与其他疫源地差异均存在统计学意义，南天山蚤群落优势度最高，与北天山和帕米尔疫源地差异存在统计学意义（表 2）。北天山疫源地宿主抗体阳性率最高，中位数为2.67%，与南天山、帕米尔和昆仑山疫源地差异存在统计学意义（χ²=138.711，P < 0.001；χ²=103.990，P < 0.001；χ²=75.808，P < 0.001）。

表 2 新疆维吾尔自治区4类山地鼠疫宿主体蚤群落生态学指标差异性分析 Table 2 Analysis on the difference of ecological indexes of the four types of plague flea communities in Xinjiang Uygur Autonomous Region

表选项

蚤类群落结构与宿主抗体阳性率呈正相关的指标有染蚤率（r_s=0.682，P < 0.001）、蚤指数（r_s=0.686，P < 0.001）、丰富度（r_s=0.529，P < 0.001）和Shannon多样性指数（r_s=0.339，P=0.023），优势度和均匀度指标与宿主抗体阳性率相关性无统计学意义（r_s=-0.201，P=0.186；r_s=-0.130，P=0.393）。

与抗体阳性率存在相关性的4个蚤群落结构指标——染蚤率、蚤指数、多样性和丰富度纳入主成分分析进行主成分提取。保留特征根值> 1所对应的前2个主成分，见图 1。其中，主成分1（PC1）特征值为2.04，贡献率为50.95%；主成分2（PC2）特征值为1.74，贡献率为43.38%。2个主成分累计贡献率达到94.33%，满足均数法确定主成分个数的原则。经KMO（Kaiser-Meyer-Olkin）及Bartlett检验，KMO值等于0.672，Bartlett值等于152.420，P < 0.001，模型拟合效果良好。

图 1 新疆维吾尔自治区4类山地鼠疫疫源地宿主体蚤群落生态学指标主成分碎石图 Figure 1 The main component lithotripsy map of the ecological index of the host flea community in 4 types of mountain plague foci in Xinjiang Uygur Autonomous Region

图选项

根据主成分系数矩阵的特征向量，2个主成分可以表示为：

(1)

(2)

式（1）和式（2）中x₁、x₂、x₃、x₄分别表示染蚤率、蚤指数、多样性和丰富度。

从主成分系数可以看出，PC1由染蚤率和蚤指数决定，反映宿主体蚤群落丰度；PC2由蚤类多样性和丰富度决定，反映宿主体蚤群落的多样性。见表 3。

表 3 新疆维吾尔自治区4类山地鼠疫疫源地鼠疫菌抗体阳性率相关因子主成分分析相关系数矩阵特征向量值 Table 3 Eigenvector values of correlation coefficient matrix of principal component analysis

表选项

将PC1和PC2纳入GAM，结果显示，宿主抗体阳性率随宿主体蚤指数和染蚤率升高而不断增高，PC1与抗体血清阳性率之间存在线性正相关关系（图 2A）；宿主抗体阳性率随多样性和丰富度升高呈先增强后减弱的波动形态，PC2与宿主抗体阳性率之间存在非线性关系（图 2B）。

图 2 基于广义相加模型的主成分1（A）和主成分2（B）与宿主抗体阳性率间关系 Figure 2 The relationship between principal components 1(A) and 2(B) and the positive rate of host antibodies by the generalized additive model

图选项

回归方程为：PV=α+s（Y₁）+s（Y₂）+ε

式中，PV为宿主抗体阳性率，α为截距，s（）为平滑函数，表示主成分的非线性关系，Y₁和Y₂分别为主成分1和主成分2，ε表示残差。

模型模拟效果检验结果显示，R²=0.743、RMSE=1.026、MAE=0.876，模型拟合效果良好。将2个重要主成分纳入随机森林模型，PC1重要性得分为28.84分，PC2为4.52分，显示宿主抗体阳性率对宿主体蚤群落的丰度（蚤指数与染蚤率）变化更为敏感。

鼠疫是一种自然疫源性疾病，鼠疫的流行与其所处的生态系统密切相关，而鼠疫自然疫源地的类型因其生态地理景观、宿主、媒介的不同，乃至病原体遗传差异性而存在丰富的多样性^[13]，不同生态系统常存在动物鼠疫流行强度上的差异^[14-15]。相关研究表明，鼠疫生态系统中的宿主、媒介影响着鼠疫自然疫源地中动物鼠疫的自然流行^[16-18]，鼠疫菌检出情况、主要宿主动物或鼠疫流行指示动物（如牧犬）的鼠疫菌血清抗体阳性率等鼠疫监测指标可以反映鼠疫流行强度，但鼠疫菌的检出受多种人为因素影响，犬类血清阳性率会受到饲养方式等因素的干扰，并且可获得的样本有限，鉴于这2种指标反映鼠疫流行的能力有限，所以本文选择主要鼠疫宿主动物血清抗体阳性率反映鼠疫流行强度，分析其与蚤类群落结构之间的关系，这对深入研究鼠疫的发生和发展有着重要意义。

第一，鼠疫自然疫源地宿主体蚤群落结构与鼠疫流行强度相关。本研究发现，新疆各山地鼠疫疫源地宿主体蚤群落存在丰富的多样性，且疫源地鼠疫的流行强度与蚤类多样性相关性存在统计学意义。各鼠疫宿主体蚤群落结构指标中，多样性、丰富度、蚤指数和染蚤率在北天山灰旱獭-长尾黄鼠鼠疫自然疫源地呈现最高水平，与其他疫源地相比差异有统计学意义，同时，宿主抗体阳性率亦显著高于其他疫源地。宿主体蚤群落指标相关因素分析发现，有些指标对疫源地鼠疫流行强度影响是共性的，比如蚤指数、染蚤率、多样性和丰富度，在4类山地疫源地中均与鼠疫流行强度呈正相关关系，有些则因疫源地的不同而存在差异，如北天山灰旱獭-长尾黄鼠宿主抗体阳性率与优势度之间存在正相关关系，而在昆仑山喜马拉雅旱獭鼠疫自然疫源地呈现负相关关系。其他学者亦有类似的发现，如陈永明等^[19]和张璐等^[20]通过对内蒙古自治区和河北省长爪沙鼠鼠疫自然疫源地调查研究发现，携带并传播鼠疫风险更大的长爪沙鼠（Meriones unguiculatus）比疫源地中其他啮齿动物有着更高的体蚤多样性与丰富度。张渝疆等^[15]对准噶尔盆地鼠疫研究亦发现，宿主体蚤的丰富度差异与该地区鼠疫流行强度之间存在相关性。

第二，宿主体蚤群落生态结构对鼠疫流行产生复合型组合和动态的影响作用。通过各群落因子的主成分分析及构建的GAM和RF模型结果，发现疫源地动物鼠疫流行强度与PC1因子，即宿主体蚤指数与染蚤率呈线性正相关关系，鼠疫流行强度随蚤指数、染蚤率升高而不断增强，RF模型的拟合结果显示主成分1重要性得分（28.84分）远高于主成分2（4.52分），这表明鼠疫流行强度的改变对蚤指数与染蚤率的变化更为灵敏。这一结果可解释Enscore等^[21]和白雪薇等^[22]的发现，即发生疫情的地区或者年份所捕获的啮齿动物染蚤率与蚤指数明显高于非疫区或无疫情发生年份^[2]，也在理论上解释了新疆各山地鼠疫疫源地中北天山鼠疫疫源地疫情流行强度最高的原因^[13]。Amatre等^[23]与Eisen等^[2]在研究乌干达森林疫源地时发现，疫源地环境中存在3种受鼠疫侵害严重的啮齿动物，它们之间共享多种蚤类，且有较高水平体蚤多样性与丰富度，证明宿主体蚤多样性和丰富度会影响鼠疫的流行^{[16, 23]}。而本研究结果显示，鼠疫流行强度随PC2因子，即宿主体蚤多样性与丰富度的增加呈现先升高后降低的动态变化过程，说明宿主体蚤群落多样性和丰富度对动物鼠疫的波动平衡具有调节作用。宿主体蚤多样性和丰富度水平较低时，蚤类多样性和丰富度升高会加剧鼠疫流行。在多样性和丰富度较高的情况下蚤类丰富度和多样性继续增加则会抑制鼠疫流行。此种宿主体蚤多样性和丰富度与鼠疫流行强度间的非线性关系，很好地反映出鼠疫在自然界中的流行不受物种多样性简单稀释或放大调节作用，而是一个复杂的动态的过程^[24-25]。在理论上，宿主体蚤群落多样性和丰富度增高到一定程度可能会降低宿主动物体蚤群落中主要传播媒介的占比，从而抑制动物鼠疫在宿主间的传播。当然这一假设需要更多的实验和调查研究予以证实，不过在鼠疫自然疫源地宿主体蚤群落的研究中亦发现存在此种宿主体蚤群落动态变化的现象。如，新疆北天山灰旱獭-长尾黄鼠鼠疫自然疫源地中，长尾黄鼠体蚤群落会随着季节的不同而发生变化，从春末到秋初，长尾黄鼠主要媒介蚤——方形黄鼠蚤占比逐渐减少，随之而来的是体蚤群落多样性和丰富度增加，以及动物鼠疫流行强度的下降^[26]。

综上所述，通过对新疆山地鼠疫自然疫源地宿主动物体蚤群落结构的研究，发现新疆4类山地鼠疫自然疫源地宿主体蚤群落结构有着非常大的差异，并且这种差异直接影响着动物鼠疫的流行。但需要注意的是，这种宿主体蚤群落差异与疫源地鼠疫流行强度的相关性，只是鼠疫与环境相互作用的复杂生态系统的一部分，并未包括宿主动物、微生物群落、植被等生物因素和气象、土壤等非生物因素与鼠疫流行的关系。因此，要更全面深入探究鼠疫流行规律依然任重而道远。

利益冲突  无