基于阶段结构蚊子种群的登革热传播动力学模型构建与分析

基于阶段结构蚊子种群的登革热传播动力学模型构建与分析一、引言1.1研究背景与意义登革热作为全球传播最为广泛的蚊媒传染病之一，正给人类健康和社会经济发展带来沉重负担。世界卫生组织（WHO）数据显示，全球约有一半人口面临登革热感染风险，每年有1亿至4亿人感染，其中320万例为有症状感染，约50-100万人需入院治疗，严重者可发展为登革出血热或登革休克综合征，甚至危及生命。在过去50年里，登革热发病率急剧上升了30倍，其传播范围不断扩大，不仅在热带和亚热带地区广泛流行，如东南亚、太平洋岛屿和加勒比海地区，还在一些原本无登革热流行的地区出现疫情，成为全球性的严重公共卫生问题。登革热的传播主要依赖于埃及伊蚊和白纹伊蚊等媒介蚊虫。蚊子的生命周期包含卵、幼虫、蛹和成虫四个阶段，不同阶段的蚊子在生存、繁殖和传播病毒的能力上存在显著差异。在温带气候下，伊蚊冬季仍可存活，但主要活跃于热带和亚热带地区。在这些适宜的环境中，蚊子繁殖速度极快，使得登革热的传播速度大大加快，短时间内即可引发大规模疫情。例如，在2014年，广东省暴发了史上最严重的登革热疫情，截止到当年10月21日，患病者登记病例就已达36889例。又如2024年，全球多个国家和地区遭受登革热侵袭，巴西卫生部5月4日公布的数据显示，巴西2024年登革热疑似及确诊病例上升至4270624例，死亡病例2197例，另有2276名死亡病例待核查；巴拉圭公共卫生和社会福利部当地时间5月3日发布的最新报告显示，从去年9月至当时，该国累计确诊登革热病例42775例，死亡病例数增至98例。鉴于蚊子种群在登革热传播中的关键作用，深入研究基于阶段结构蚊子种群的登革热传播动力学模型具有极其重要的意义。通过构建和分析这样的模型，我们能够更加准确地揭示登革热的传播机制，明确不同阶段蚊子种群对病毒传播的影响程度。这不仅有助于我们从理论层面深入理解登革热的传播规律，还能为制定科学有效的防控策略提供坚实的理论依据，从而减少登革热的传播范围和发病率，保护人类健康，减轻社会经济负担。1.2国内外研究现状在登革热传播动力学模型的研究领域，国内外学者已取得了丰硕的成果。早期的研究主要集中在构建简单的传染病模型来描述登革热在人群中的传播，如经典的SIR（易感者-感染者-康复者）模型及其衍生模型。这些模型为理解登革热的传播过程提供了基础框架，能够初步分析疾病的传播趋势和基本再生数等关键指标。随着研究的深入，考虑到蚊子种群在登革热传播中的重要作用，越来越多的学者开始将蚊子种群纳入模型中。一些研究构建了包含蚊子种群的SIR型模型，将蚊子分为易感蚊、感染蚊和带毒蚊等不同类别，分析蚊子种群动态与登革热传播之间的相互关系。还有研究进一步考虑了蚊子的繁殖、死亡等因素对登革热传播的影响，使得模型更加贴近实际情况。近年来，随着对蚊子生物学特性的深入了解，考虑蚊子种群阶段结构的登革热传播动力学模型逐渐成为研究热点。这类模型将蚊子的生命周期划分为卵、幼虫、蛹和成虫等不同阶段，分别考虑各阶段的生长、发育、死亡以及它们与登革热传播的关联。例如，部分研究通过实验和数据分析，确定了不同阶段蚊子的发育时间、存活率以及成虫的繁殖能力等参数，并将这些参数纳入模型中，以更准确地模拟登革热的传播过程。然而，目前的研究仍存在一些不足之处。一方面，虽然考虑了蚊子种群的阶段结构，但对于不同阶段蚊子之间的相互作用以及它们对环境因素的响应机制，研究还不够深入。例如，环境温度、湿度等因素如何影响蚊子不同阶段的发育速率和存活率，以及这些影响如何进一步作用于登革热的传播，尚未得到充分的探讨。另一方面，在模型的应用方面，虽然已有研究尝试利用模型来评估防控策略的效果，但由于模型的复杂性和实际情况的多样性，模型的预测结果与实际防控效果之间仍存在一定的差距。因此，如何进一步完善模型，提高其对实际情况的预测能力和指导作用，是未来研究需要解决的重要问题。1.3研究内容与方法本论文旨在深入研究基于阶段结构蚊子种群的登革热传播动力学模型，具体研究内容如下：模型构建：综合考虑登革热病毒在人群和蚊子种群中的传播过程，以及蚊子种群的阶段结构，构建更加符合实际情况的登革热传播动力学模型。在模型中，详细描述人群中易感者、感染者和康复者之间的转换关系，以及蚊子种群中卵、幼虫、蛹和成虫各阶段的生长、发育、死亡和感染过程。平衡点分析：对构建的模型进行平衡点分析，确定模型的无病平衡点和地方病平衡点，并研究这些平衡点的稳定性。通过分析平衡点的稳定性，判断登革热在人群中是否会爆发以及在何种条件下会保持稳定传播。数值模拟：利用数值计算方法对模型进行模拟，分析不同参数对登革热传播的影响。通过改变蚊子种群各阶段的发育时间、繁殖率、死亡率以及人群的感染率、康复率等参数，观察登革热的传播趋势和发病率的变化，从而为防控策略的制定提供依据。本研究拟采用以下方法：数学建模方法：基于传染病动力学原理，结合蚊子种群的生物学特性和登革热的传播特点，建立数学模型来描述登革热的传播过程。在建模过程中，充分考虑各种因素之间的相互关系，确保模型的合理性和准确性。理论分析方法：运用常微分方程理论、稳定性理论等数学工具，对模型的平衡点进行分析，研究其稳定性条件。通过理论分析，揭示登革热传播的内在规律，为数值模拟和防控策略的制定提供理论基础。数值计算方法：采用数值计算软件，如MATLAB等，对模型进行数值求解和模拟。通过数值模拟，直观地展示登革热的传播过程和不同参数对传播的影响，为研究结果的分析和讨论提供数据支持。二、蚊子种群阶段结构与登革热传播概述2.1蚊子种群的阶段划分及特征蚊子作为完全变态昆虫，其生命周期历经卵、幼虫、蛹和成虫四个显著不同的阶段，各阶段在形态、生活习性和存活时间上展现出独特的特征。在卵阶段，蚊子卵的形态因蚊种而异，如按蚊卵呈船形，伊蚊卵为圆锥形，它们常成簇地附着在水面或湿润边缘。卵的孵化对环境条件要求较为苛刻，需适宜的温度（20-30℃）与湿度。以埃及伊蚊为例，其卵具备特殊的“防水涂层”，这使其能够在干燥环境中休眠长达数月之久，一旦积水形成，便迅速孵化，这种特性极大地确保了种群的延续。蚊子幼虫，俗称孑孓，在水中以螺旋状方式游动。不同蚊种的幼虫呼吸方式有所不同，伊蚊幼虫通过尾部的呼吸吸管获取氧气，而按蚊幼虫则借助体侧呼吸孔进行呼吸。它们主要以藻类、细菌或有机碎屑为食，在生长过程中会经历3-4次蜕皮，进而发育为蛹。幼虫对水体类型的选择具有高度敏感性，库蚊偏好污浊积水，按蚊倾向于清洁流水，白纹伊蚊作为登革热的主要传播媒介之一，常利用花盆托盘、轮胎积水等城市微环境进行繁殖。当水体受到震动时，幼虫会迅速下沉躲避天敌，这种“逃生技能”有助于它们在复杂的生存环境中存活。蛹期的蚊子静止悬浮于水面，通过体表呼吸孔维持基本的代谢活动，但不再摄食。此时的蚊子对环境变化极为敏感，若水体出现扰动或干涸的情况，蛹可能会提前羽化。在非洲季节性河流中，蛹能够在短暂的干涸期存活，待雨季来临，便会大量羽化为成虫，这体现了蛹期作为抵御环境变化缓冲阶段的重要作用。在羽化前，蛹的背部会裂开，成虫像“破茧蝴蝶”般钻出，整个过程仅需短短几秒。成虫阶段的蚊子，雌雄个体在习性上存在明显分化。雄蚊仅吸食植物汁液，而雌蚊为了获取蛋白质以发育卵巢，需要吸食人类或动物的血液。多数蚊种在羽化后1-2天便完成交配，雌蚊单次产卵量可达100-300枚，部分种类还可多次产卵。不同蚊种的活动节律差异显著，按蚊多在黄昏或清晨吸血，埃及伊蚊则白天活跃，这种活动习性使其能够更好地适应城市人类的作息，增加了与人类接触的机会，从而提高了登革热病毒的传播风险。在适宜条件下，雄蚊通常可活1-3周，雌蚊在夏天一般能活1个月左右，而越冬的雌蚊寿命可达数月之久。2.2登革热传播机制登革热病毒主要通过埃及伊蚊和白纹伊蚊等蚊子的叮咬在人和蚊子之间传播，整个传播过程涉及多个复杂的环节。当蚊子叮咬登革热患者时，带有登革病毒的血液被吸入蚊子体内。病毒首先在蚊子的中肠上皮细胞中感染和复制，中肠上皮细胞含有接受登革病毒的膜受体，如寨卡病毒受体和Laminin受体，这些受体促进了病毒的进入和复制。蚊子的中肠屏障在调节病毒感染中起着关键作用，它由一群紧密连接的上皮细胞组成，有助于防止病原体的入侵，并可能限制病毒在蚊子体内的传播。然而，对于易感染登革病毒的蚊子种类，病毒能够突破中肠屏障，进入蚊子的血淋巴，并进一步感染其他组织，如唾液腺。在蚊子的唾液腺中，病毒大量复制并积累。当感染病毒的蚊子再次叮咬健康人时，蚊子唾液中的病毒会随着唾液注入人体皮肤。病毒随后通过皮肤上的伤口或毛孔进入人体血液，进而感染人体细胞。登革病毒具有广泛的宿主范围，在人体内，它主要感染单核巨噬细胞、内皮细胞和淋巴细胞等。病毒进入细胞后，利用宿主细胞的物质和能量进行复制和组装，产生大量的子代病毒。这些子代病毒释放后，又可以继续感染周围的细胞，从而导致病毒在体内的扩散。登革热的传播效率受到多种因素的影响。环境因素方面，温度和湿度对蚊子的繁殖、发育和活动能力有着重要影响。适宜的温度（25-32℃）和较高的湿度（60%-80%）有利于蚊子的生长和繁殖，同时也能提高病毒在蚊子体内的复制效率。在高温高湿的环境下，蚊子的繁殖速度加快，寿命延长，叮咬人类的频率增加，从而使得登革热病毒的传播风险大大提高。此外，降雨也会影响蚊子的滋生地，雨水形成的积水为蚊子提供了更多的繁殖场所。宿主免疫力也是影响登革热传播的关键因素。人体的免疫系统在感染登革病毒后会启动一系列的免疫反应，包括细胞免疫和体液免疫。初次感染登革病毒时，免疫系统能够识别病毒抗原并产生相应的抗体，这些抗体可以中和病毒，阻止病毒的进一步感染和扩散。然而，登革病毒有4种血清型，当人体再次感染不同血清型的登革病毒时，可能会发生抗体依赖的增强作用（ADE）。ADE效应会导致免疫系统对病毒的清除能力下降，反而促进病毒在体内的复制和感染，使得病情加重，增加了登革热传播的复杂性和严重性。2.3蚊子种群阶段结构对登革热传播的影响机制蚊子种群的阶段结构对登革热的传播有着深远的影响，这种影响主要通过不同阶段蚊子的数量变化、活动能力以及感染病毒的几率等因素来实现。蚊子不同阶段的数量变化直接关系到登革热的传播风险。卵的数量决定了未来蚊子种群的潜在规模，如果环境中存在大量适宜蚊子产卵的场所，如积水的容器、花瓶托盘等，就会导致蚊子卵的大量产生。在适宜的温度和湿度条件下，这些卵迅速孵化，使得幼虫数量急剧增加。幼虫是蚊子生长发育的重要阶段，其数量的多少会影响到后续蛹和成虫的数量。如果幼虫阶段的死亡率较低，大量幼虫成功发育为蛹并最终羽化为成虫，就会导致蚊子成虫种群数量的大幅增长。而成虫是登革热病毒的传播媒介，成虫数量的增加意味着登革热传播的机会增多，病毒更容易在人群中扩散。蚊子不同阶段的活动能力差异也对登革热传播产生影响。幼虫和蛹主要生活在水中，活动范围相对有限，它们与人类的直接接触机会较少，因此在登革热传播中直接作用相对较小。然而，幼虫和蛹的生存环境与成虫的繁殖密切相关，适宜的幼虫和蛹生存环境有助于成虫的大量繁殖。成虫具有飞行能力，能够主动寻找宿主进行吸血。埃及伊蚊和白纹伊蚊等传播登革热的主要蚊种，白天活跃且偏好叮咬人类，它们的活动能力使得登革热病毒能够在不同宿主之间快速传播。成虫的飞行范围和活动频率决定了病毒传播的范围和速度，如果蚊子的活动范围广泛，频繁叮咬不同的人群，就会加速登革热的传播。蚊子不同阶段感染病毒的几率不同，也在很大程度上影响着登革热的传播。虽然幼虫和蛹在自然条件下感染登革病毒的几率相对较低，但如果它们生存的水体受到病毒污染，仍有可能感染病毒。一旦幼虫或蛹感染病毒，可能会影响它们的正常发育，并且在羽化为成虫后，成为具有传播病毒能力的感染源。成虫感染登革病毒的几率与多种因素有关，包括蚊子叮咬感染登革热患者的频率、病毒在蚊子体内的感染和复制效率等。当蚊子叮咬感染登革热病毒的患者后，如果病毒能够在蚊子体内成功复制并扩散到唾液腺，那么这只蚊子在下次叮咬健康人时，就会将病毒传播给新的宿主。感染病毒的成虫数量越多，登革热传播的风险就越高。蚊子种群动态与登革热传播之间存在着紧密的内在联系。蚊子的繁殖、发育和死亡过程受到环境因素、宿主可用性以及自身生物学特性的影响，而这些因素的变化又会反过来影响登革热的传播。例如，气候变化导致的温度和降水模式改变，可能会影响蚊子的繁殖季节和繁殖场所，进而改变蚊子种群的数量和分布。如果温度升高使得蚊子的繁殖周期缩短，或者降水增加提供了更多的积水滋生地，都会导致蚊子种群数量上升，从而增加登革热的传播风险。此外，人类活动如城市化进程、环境卫生状况等也会对蚊子种群动态和登革热传播产生影响。城市化过程中大量的人工建筑和废弃容器为蚊子提供了更多的繁殖场所，而环境卫生不良则有利于蚊子的滋生和生存。了解蚊子种群阶段结构对登革热传播的影响机制，为后续构建登革热传播动力学模型提供了重要的生物学依据，有助于我们更准确地模拟和预测登革热的传播趋势，制定有效的防控策略。三、带有阶段结构蚊子种群的登革热传播动力学模型构建3.1模型假设为了构建更贴合实际的登革热传播动力学模型，我们提出以下一系列假设：种群划分：将蚊子种群依据其生物学特性，细致地划分为卵（E）、幼虫（L）、蛹（P）、易感成虫（S_m）和感染成虫（I_m）五个阶段；同时，把人群划分为易感人群（S_h）、感染人群（I_h）和康复人群（R_h）三类。出生率与死亡率：蚊子种群各阶段均具有特定的出生率和死亡率。其中，成年雌性蚊子的出生率受到环境因素以及自身繁殖能力的影响，我们假设其出生率为常数\beta。各阶段蚊子的死亡率分别为\mu_E（卵）、\mu_L（幼虫）、\mu_P（蛹）、\mu_m（成虫）。对于人群，人的自然出生率和死亡率分别设为\Lambda_h和\mu_h，且假设在研究期间人口总数保持相对稳定，即\Lambda_h=\mu_hN_h，其中N_h=S_h+I_h+R_h为总人口数。发育速率：蚊子从一个阶段发育到下一个阶段的速率是固定的。例如，卵发育为幼虫的速率为k_1，幼虫发育为蛹的速率为k_2，蛹发育为成虫的速率为k_3。这些发育速率受到环境温度、湿度等因素的影响，在本模型中，我们先假设它们为常数，后续可进一步研究环境因素对其的影响。感染与传播：易感蚊子叮咬感染登革热的人后，以概率\beta_{hm}感染病毒，从而转变为感染成虫；感染蚊子叮咬易感人群时，以概率\beta_{mh}将病毒传播给人，使易感人群转变为感染人群。假设病毒传播系数\beta_{hm}和\beta_{mh}为常数，但在实际情况中，它们可能受到蚊子叮咬频率、病毒载量以及人群免疫力等因素的影响，后续可考虑将这些因素纳入模型进行深入分析。人群中的感染人群在经过平均感染期1/\gamma_h后康复，并获得一定的免疫力，康复率为\gamma_h。此外，假设康复人群具有终身免疫，不会再次感染登革热。3.2模型建立基于上述假设，我们利用常微分方程建立带有阶段结构蚊子种群的登革热传播动力学模型，如下所示：\begin{cases}\frac{dE}{dt}=\betaS_m-(\mu_E+k_1)E\\\frac{dL}{dt}=k_1E-(\mu_L+k_2)L\\\frac{dP}{dt}=k_2L-(\mu_P+k_3)P\\\frac{dS_m}{dt}=k_3P-\mu_mS_m-\beta_{hm}\frac{S_mI_h}{N_h}\\\frac{dI_m}{dt}=\beta_{hm}\frac{S_mI_h}{N_h}-\mu_mI_m\\\frac{dS_h}{dt}=\Lambda_h-\mu_hS_h-\beta_{mh}\frac{S_hI_m}{N_h}\\\frac{dI_h}{dt}=\beta_{mh}\frac{S_hI_m}{N_h}-(\mu_h+\gamma_h)I_h\\\frac{dR_h}{dt}=\gamma_hI_h-\mu_hR_h\end{cases}在这个模型中，各变量的含义如下：E(t)：t时刻蚊子卵的数量；L(t)：t时刻蚊子幼虫的数量；P(t)：t时刻蚊子蛹的数量；S_m(t)：t时刻易感蚊子成虫的数量；I_m(t)：t时刻感染蚊子成虫的数量；S_h(t)：t时刻易感人群的数量；I_h(t)：t时刻感染人群的数量；R_h(t)：t时刻康复人群的数量。各参数的含义及取值范围如下：\beta：成年雌性蚊子的出生率，通常根据蚊子的种类和环境条件而定，取值范围一般在[50,300]之间，例如埃及伊蚊在适宜环境下，每只雌性蚊子每次产卵约100-200枚，这里假设\beta=150；\mu_E：蚊子卵的死亡率，受到环境因素如干燥、温度不适宜等影响，取值范围约为[0.1,0.5]，假设\mu_E=0.3；k_1：卵发育为幼虫的速率，与温度密切相关，在适宜温度下，埃及伊蚊卵发育为幼虫约需1-2天，换算为速率约为k_1=0.5（假设时间单位为天）；\mu_L：蚊子幼虫的死亡率，受食物资源、天敌等因素影响，取值范围约为[0.2,0.6]，假设\mu_L=0.4；k_2：幼虫发育为蛹的速率，同样与温度有关，在适宜条件下，埃及伊蚊幼虫发育为蛹约需5-7天，速率约为k_2=0.15；\mu_P：蚊子蛹的死亡率，受环境扰动等因素影响，取值范围约为[0.1,0.3]，假设\mu_P=0.2；k_3：蛹发育为成虫的速率，在适宜条件下，埃及伊蚊蛹发育为成虫约需2-3天，速率约为k_3=0.4；\mu_m：蚊子成虫的死亡率，受环境、天敌和寿命等因素影响，取值范围约为[0.05,0.2]，假设\mu_m=0.1；\beta_{hm}：人传蚊的传播系数，与蚊子叮咬感染人群的频率、病毒载量等有关，取值范围约为[0.1,0.5]，假设\beta_{hm}=0.3；\beta_{mh}：蚊传人传播系数，与蚊子叮咬易感人群的频率、病毒在蚊子体内的复制情况等有关，取值范围约为[0.2,0.6]，假设\beta_{mh}=0.4；\Lambda_h：人的补充率，即自然出生率，根据人口统计数据，假设\Lambda_h=0.01N_h（这里N_h为总人口数）；\mu_h：人的自然死亡率，根据人口统计数据，假设\mu_h=0.01；\gamma_h：人的康复率，根据临床数据，登革热患者的平均康复期约为7-10天，换算为康复率约为\gamma_h=0.1。3.3参数估计参数估计是模型构建中的关键环节，它直接影响模型的准确性和可靠性。为了对上述模型中的参数进行估计，我们广泛收集了相关文献数据以及实际调查数据，并运用统计学方法和最小二乘法等进行分析。首先，对于蚊子种群相关参数，我们参考了大量关于埃及伊蚊和白纹伊蚊的生物学研究文献。例如，通过对不同环境条件下蚊子繁殖和发育的实验数据进行分析，确定了成年雌性蚊子的出生率\beta、各阶段蚊子的死亡率\mu_E、\mu_L、\mu_P、\mu_m以及发育速率k_1、k_2、k_3。同时，我们还结合了当地的实际调查数据，对这些参数进行了进一步的校准和验证。例如，在某地区进行的蚊子种群动态监测中，我们记录了不同时间点蚊子卵、幼虫、蛹和成虫的数量变化，通过拟合这些数据，对参数进行了调整，使其更符合当地的实际情况。对于人群相关参数，我们收集了当地的人口统计数据，以确定人的自然出生率\Lambda_h和自然死亡率\mu_h。同时，通过分析医院的临床病例数据，统计登革热患者的感染期和康复情况，从而估计出人群的感染率\beta_{mh}、康复率\gamma_h以及人传蚊的传播系数\beta_{hm}。例如，在某地区的登革热疫情调查中，我们详细记录了患者的感染时间、症状表现以及康复时间等信息，通过对这些数据的分析，计算出了人群的感染率和康复率。在参数估计过程中，我们运用了最小二乘法来拟合模型与实际数据之间的差异。具体来说，我们定义一个目标函数，该函数表示模型预测值与实际观测值之间的误差平方和。通过调整模型中的参数，使得目标函数达到最小值，从而得到最优的参数估计值。例如，对于蚊子种群数量的预测，我们将模型计算得到的不同阶段蚊子数量与实际监测到的蚊子数量进行比较，通过最小化两者之间的误差平方和，来确定蚊子种群相关参数的估计值。然而，需要注意的是，参数估计存在一定的不确定性。一方面，实际数据的收集可能存在误差，例如蚊子种群数量的监测可能由于采样方法的局限性而不够准确，人群病例数据的统计可能存在漏报或误诊的情况。另一方面，模型本身的假设与实际情况可能存在差异，例如我们假设病毒传播系数为常数，但实际情况中它们可能受到多种因素的动态影响。为了评估参数估计的不确定性，我们采用了蒙特卡洛模拟方法。通过随机生成一系列符合参数取值范围的参数组合，多次运行模型，并分析模型结果的分布情况，我们可以得到参数估计值的置信区间，从而评估参数估计的可靠性。例如，在蒙特卡洛模拟中，我们生成了1000组参数组合，分别运行模型，得到了1000组登革热传播的模拟结果。通过分析这些结果，我们确定了参数估计值的95%置信区间，为模型的应用提供了更可靠的依据。四、模型分析4.1平衡点分析对于我们所构建的带有阶段结构蚊子种群的登革热传播动力学模型，平衡点是指系统处于稳定状态时各变量的值，即当系统达到平衡时，各变量的导数为零。通过求解方程组\frac{dE}{dt}=0，\frac{dL}{dt}=0，\frac{dP}{dt}=0，\frac{dS_m}{dt}=0，\frac{dI_m}{dt}=0，\frac{dS_h}{dt}=0，\frac{dI_h}{dt}=0，\frac{dR_h}{dt}=0，我们可以得到模型的平衡点。4.1.1无病平衡点当I_h=0且I_m=0时，模型存在无病平衡点E_0。此时，方程组可简化为：\begin{cases}\betaS_m-(\mu_E+k_1)E=0\\k_1E-(\mu_L+k_2)L=0\\k_2L-(\mu_P+k_3)P=0\\k_3P-\mu_mS_m=0\\\Lambda_h-\mu_hS_h=0\\\gamma_hI_h-\mu_hR_h=0\end{cases}求解上述方程组，可得无病平衡点E_0的坐标为：E_0=\left(\frac{\beta\Lambda_h}{\mu_h(\mu_E+k_1)},\frac{k_1\beta\Lambda_h}{\mu_h(\mu_E+k_1)(\mu_L+k_2)},\frac{k_1k_2\beta\Lambda_h}{\mu_h(\mu_E+k_1)(\mu_L+k_2)(\mu_P+k_3)},\frac{k_1k_2k_3\beta\Lambda_h}{\mu_h\mu_m(\mu_E+k_1)(\mu_L+k_2)(\mu_P+k_3)},0,\frac{\Lambda_h}{\mu_h},0,0\right)为了分析无病平衡点E_0的稳定性，我们对模型在E_0处进行线性化处理。设x=(E-E_0^1,L-E_0^2,P-E_0^3,S_m-E_0^4,I_m,S_h-E_0^6,I_h,R_h-E_0^8)^T，其中E_0^i表示无病平衡点E_0中第i个变量的值。将模型在E_0处进行泰勒展开，并忽略高阶项，得到线性化后的系统：\frac{dx}{dt}=Ax其中A是线性化后的系数矩阵，其元素可通过对原模型在E_0处求偏导数得到。例如，A_{11}=\frac{\partial(\betaS_m-(\mu_E+k_1)E)}{\partialE}\big|_{E_0}=-(\mu_E+k_1)，A_{14}=\frac{\partial(\betaS_m-(\mu_E+k_1)E)}{\partialS_m}\big|_{E_0}=\beta等。通过计算系数矩阵A的特征值，我们可以判断无病平衡点的稳定性。若所有特征值的实部均小于零，则无病平衡点E_0是局部渐近稳定的，这意味着当系统受到小的扰动时，会逐渐回到无病平衡点，疾病会逐渐消除；若存在特征值的实部大于零，则无病平衡点E_0是不稳定的，疾病有可能爆发并在人群中传播。经计算，系数矩阵A的特征方程为：\begin{vmatrix}-(\mu_E+k_1)-\lambda&0&0&\beta&0&0&0&0\\k_1&-(\mu_L+k_2)-\lambda&0&0&0&0&0&0\\0&k_2&-(\mu_P+k_3)-\lambda&0&0&0&0&0\\0&0&k_3&-\mu_m-\lambda&0&0&0&0\\0&0&0&-\beta_{hm}\frac{I_h}{N_h}\big|_{E_0}&-\mu_m-\lambda&-\beta_{hm}\frac{S_m}{N_h}\big|_{E_0}&0&0\\0&0&0&0&-\beta_{mh}\frac{S_h}{N_h}\big|_{E_0}&-\mu_h-\lambda&-\beta_{mh}\frac{I_m}{N_h}\big|_{E_0}&0\\0&0&0&0&\beta_{mh}\frac{S_h}{N_h}\big|_{E_0}&0&-(\mu_h+\gamma_h)-\lambda&0\\0&0&0&0&0&0&\gamma_h&-\mu_h-\lambda\end{vmatrix}=0由于I_h=0且I_m=0在无病平衡点E_0处，上式可进一步简化。经过一系列计算（具体计算过程见附录），我们得到特征值\lambda_1=-(\mu_E+k_1)，\lambda_2=-(\mu_L+k_2)，\lambda_3=-(\mu_P+k_3)，\lambda_4=-\mu_m，\lambda_5=-\mu_h，以及一个关于\lambda的二次方程：\lambda^2+(\mu_m+\mu_h+\gamma_h)\lambda+\mu_m\mu_h+\mu_m\gamma_h-\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}=0根据二次方程的求根公式\lambda=\frac{-b\pm\sqrt{b^2-4ac}}{2a}，其中a=1，b=\mu_m+\mu_h+\gamma_h，c=\mu_m\mu_h+\mu_m\gamma_h-\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}，我们分析该二次方程的根的情况。当\Delta=b^2-4ac=(\mu_m+\mu_h+\gamma_h)^2-4(\mu_m\mu_h+\mu_m\gamma_h-\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0})\geq0时，两个根为：\lambda_{6,7}=\frac{-(\mu_m+\mu_h+\gamma_h)\pm\sqrt{(\mu_m+\mu_h+\gamma_h)^2-4(\mu_m\mu_h+\mu_m\gamma_h-\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0})}}{2}当\Delta<0时，两个根为共轭复数：\lambda_{6,7}=\frac{-(\mu_m+\mu_h+\gamma_h)\pmi\sqrt{4(\mu_m\mu_h+\mu_m\gamma_h-\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0})-(\mu_m+\mu_h+\gamma_h)^2}}{2}为了使无病平衡点E_0局部渐近稳定，需要\lambda_{6,7}的实部小于零，即：-\frac{\mu_m+\mu_h+\gamma_h}{2}<0这是显然成立的。同时，还需要：\mu_m\mu_h+\mu_m\gamma_h-\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}>0整理可得：\frac{\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}}{\mu_m(\mu_h+\gamma_h)}<1令R_0=\frac{\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}}{\mu_m(\mu_h+\gamma_h)}，R_0即为模型的基本再生数，其生物学意义为在完全易感人群中，一个初始感染者在其感染期内平均能够感染的新个体数。当R_0<1时，无病平衡点E_0是局部渐近稳定的，疾病会逐渐消除；当R_0>1时，无病平衡点E_0是不稳定的，疾病有可能在人群中爆发并持续传播。4.1.2地方病平衡点当I_h\neq0且I_m\neq0时，模型存在地方病平衡点E^*。此时，求解方程组较为复杂，我们可以通过数值方法来近似求解。设地方病平衡点E^*=(E^*,L^*,P^*,S_m^*,I_m^*,S_h^*,I_h^*,R_h^*)，满足：\begin{cases}\betaS_m^*-(\mu_E+k_1)E^*=0\\k_1E^*-(\mu_L+k_2)L^*=0\\k_2L^*-(\mu_P+k_3)P^*=0\\k_3P^*-\mu_mS_m^*-\beta_{hm}\frac{S_m^*I_h^*}{N_h^*}=0\\\beta_{hm}\frac{S_m^*I_h^*}{N_h^*}-\mu_mI_m^*=0\\\Lambda_h-\mu_hS_h^*-\beta_{mh}\frac{S_h^*I_m^*}{N_h^*}=0\\\beta_{mh}\frac{S_h^*I_m^*}{N_h^*}-(\mu_h+\gamma_h)I_h^*=0\\\gamma_hI_h^*-\mu_hR_h^*=0\end{cases}其中N_h^*=S_h^*+I_h^*+R_h^*。对于地方病平衡点E^*的稳定性分析，我们同样对模型在E^*处进行线性化处理，得到线性化后的系数矩阵A^*。通过计算A^*的特征值来判断地方病平衡点的稳定性。若所有特征值的实部均小于零，则地方病平衡点E^*是局部渐近稳定的，这意味着疾病会在一定水平上持续存在；若存在特征值的实部大于零，则地方病平衡点E^*是不稳定的，疾病的传播情况会发生变化。计算A^*的特征值过程较为繁琐，通常需要借助计算机软件进行数值计算。以MATLAB软件为例，我们可以使用eig函数来计算矩阵的特征值。首先，定义线性化后的系数矩阵A^*，然后调用eig函数：A_star=[-(\mu_E+k_1)&0&0&\beta&0&0&0&0;k_1&-(\mu_L+k_2)&0&0&0&0&0&0;0&k_2&-(\mu_P+k_3)&0&0&0&0&0;0&0&k_3&-\mu_m-\beta_{hm}\frac{I_h^*}{N_h^*}&-\beta_{hm}\frac{S_m^*}{N_h^*}&0&0&0;0&0&0&\beta_{hm}\frac{I_h^*}{N_h^*}&-\mu_m&\beta_{hm}\frac{S_m^*}{N_h^*}&0&0;0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);-(\mu_E+k_1)&0&0&\beta&0&0&0&0;k_1&-(\mu_L+k_2)&0&0&0&0&0&0;0&k_2&-(\mu_P+k_3)&0&0&0&0&0;0&0&k_3&-\mu_m-\beta_{hm}\frac{I_h^*}{N_h^*}&-\beta_{hm}\frac{S_m^*}{N_h^*}&0&0&0;0&0&0&\beta_{hm}\frac{I_h^*}{N_h^*}&-\mu_m&\beta_{hm}\frac{S_m^*}{N_h^*}&0&0;0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);k_1&-(\mu_L+k_2)&0&0&0&0&0&0;0&k_2&-(\mu_P+k_3)&0&0&0&0&0;0&0&k_3&-\mu_m-\beta_{hm}\frac{I_h^*}{N_h^*}&-\beta_{hm}\frac{S_m^*}{N_h^*}&0&0&0;0&0&0&\beta_{hm}\frac{I_h^*}{N_h^*}&-\mu_m&\beta_{hm}\frac{S_m^*}{N_h^*}&0&0;0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);0&k_2&-(\mu_P+k_3)&0&0&0&0&0;0&0&k_3&-\mu_m-\beta_{hm}\frac{I_h^*}{N_h^*}&-\beta_{hm}\frac{S_m^*}{N_h^*}&0&0&0;0&0&0&\beta_{hm}\frac{I_h^*}{N_h^*}&-\mu_m&\beta_{hm}\frac{S_m^*}{N_h^*}&0&0;0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);0&0&k_3&-\mu_m-\beta_{hm}\frac{I_h^*}{N_h^*}&-\beta_{hm}\frac{S_m^*}{N_h^*}&0&0&0;0&0&0&\beta_{hm}\frac{I_h^*}{N_h^*}&-\mu_m&\beta_{hm}\frac{S_m^*}{N_h^*}&0&0;0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);0&0&0&\beta_{hm}\frac{I_h^*}{N_h^*}&-\mu_m&\beta_{hm}\frac{S_m^*}{N_h^*}&0&0;0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);0&0&0&0&-\beta_{mh}\frac{S_h^*}{N_h^*}&-\mu_h-\beta_{mh}\frac{I_m^*}{N_h^*}&-\beta_{mh}\frac{S_h^*}{N_h^*}&0;0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);0&0&0&0&\beta_{mh}\frac{S_h^*}{N_h^*}&0&-(\mu_h+\gamma_h)&0;0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);0&0&0&0&0&0&\gamma_h&-\mu_h];eigenvalues=eig(A_star);];eigenvalues=eig(A_star);eigenvalues=eig(A_star);通过上述代码，我们可以得到系数矩阵A^*的特征值eigenvalues。然后，分析这些特征值的实部，判断地方病平衡点E^*的稳定性。若所有特征值的实部均小于零，则地方病平衡点E^*是局部渐近稳定的，疾病会在一定水平上持续存在；若存在特征值的实部大于零，则地方病平衡点E^*是不稳定的，疾病的传播情况会发生变化，可能会出现疫情的爆发或消退等情况。4.2基本再生数计算基本再生数R_0是传染病动力学模型中的一个关键指标，它反映了在完全易感人群中，一个初始感染者在其感染期内平均能够感染的新个体数。对于我们所构建的登革热传播动力学模型，采用下一代矩阵法来计算基本再生数。首先，将模型中的变量分为两类：感染类变量x=(I_m,I_h)^T和非感染类变量y=(E,L,P,S_m,S_h,R_h)^T。根据下一代矩阵法的原理，我们需要先求出模型在无病平衡点E_0处的感染率矩阵F和转移率矩阵V。感染率矩阵F表示新感染的产生速率，其元素F_{ij}表示从第j类感染状态到第i类感染状态的新感染速率。在我们的模型中：F=\begin{pmatrix}0&\beta_{hm}\frac{S_m}{N_h}\big|_{E_0}\\\beta_{mh}\frac{S_h}{N_h}\big|_{E_0}&0\end{pmatrix}其中，F_{12}表示易感蚊子叮咬感染人群后被感染的速率，F_{21}表示感染蚊子叮咬易感人群后导致人群感染的速率。转移率矩阵V表示感染类变量的转移速率，其元素V_{ij}表示从第j类感染状态到第i类感染状态的转移速率。在我们的模型中：V=\begin{pmatrix}\mu_m&0\\0&\mu_h+\gamma_h\end{pmatrix}其中，V_{11}表示感染蚊子的死亡率，V_{22}表示感染人群的死亡率与康复率之和。然后，计算下一代矩阵K=FV^{-1}：V^{-1}=\begin{pmatrix}\frac{1}{\mu_m}&0\\0&\frac{1}{\mu_h+\gamma_h}\end{pmatrix}K=FV^{-1}=\begin{pmatrix}0&\frac{\beta_{hm}\frac{S_m}{N_h}\big|_{E_0}}{\mu_m}\\\frac{\beta_{mh}\frac{S_h}{N_h}\big|_{E_0}}{\mu_h+\gamma_h}&0\end{pmatrix}基本再生数R_0即为下一代矩阵K的谱半径，即R_0=\rho(K)。对于2\times2矩阵K，其谱半径可通过计算矩阵的特征值来得到。矩阵K的特征方程为：\begin{vmatrix}-\lambda&\frac{\beta_{hm}\frac{S_m}{N_h}\big|_{E_0}}{\mu_m}\\\frac{\beta_{mh}\frac{S_h}{N_h}\big|_{E_0}}{\mu_h+\gamma_h}&-\lambda\end{vmatrix}=0即：\lambda^2-\frac{\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}}{\mu_m(\mu_h+\gamma_h)}=0解上述方程可得：\lambda_{1,2}=\pm\sqrt{\frac{\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\big|_{E_0}}{\mu_m(\mu_h+\gamma_h)}}基本再生数R_0为特征值的模的最大值，即：R_0=\sqrt{\frac{\beta_{hm}\beta_{mh}\frac{S_mS_h}{N_h^2}\\##五、数值模拟与结果分析\##\#5.1数值模拟方法为了深入探究所构建的登革热ä¼

播动力学模型的动态行为，本ç

”究选用四阶龙æ

¼-库塔（Runge-Kutta）法进行数值模拟。四阶龙æ

¼-库塔法是一种在工程和科学计算中广泛应用的高精度单步算法，其优势显著。与其他数值方法相比，如简单的欧拉法，四阶龙æ

¼-库塔法在计算精度上有了质的提升。欧拉法仅使用当前点的斜率来预估下一个点的函数值，而四阶龙æ

¼-库塔法通过在区间内多个点上预估斜率值，并采用åŠ

权平均的方式来计算下一个点的函数值，从而有效降低了误差，提高了计算精度。在处理复杂的微分方程系统时，其稳定性和可é

性也得到了广泛验证，能够更准确地模拟系统的动态变化。该方法的具体实施步骤如下：对于一阶常微分方程$\frac{dy}{dt}=f(t,y)$，给定初始条件$y(t_0)=y_0$，在每个时间步长$h$内，通过以下公式计算下一个时间点的函数值$y_{n+1}$：\[\begin{align*}k_1&=hf(t_n,y_n)\\k_2&=hf(t_n+\frac{h}{2},y_n+\frac{k_1}{2})\\k_3&=hf(t_n+\frac{h}{2},y_n+\frac{k_2}{2})\\k_4&=hf(t_n+h,y_n+k_3)\\y_{n+1}&=y_n+\frac{1}{6}(k_1+2k_2+2k_3+k_4)\end{align*}其中，k_1、k_2、k_3和k_4分别表示在不同点上计算得到的斜率值，通过对这四个斜率值进行加权平均，得到更准确的平均斜率，进而计算出下一个时间点的函数值。在对本模型进行数值模拟时，将模型中的八个微分方程分别按照上述四阶龙格-库塔法的步骤进行求解。例如，对于方程\frac{dE}{dt}=\betaS_m-(\mu_E+k_1)E，将\frac{dE}{dt}看作f(t,E)，按照公式计算k_1、k_2、k_3和k_4，进而得到E在下一步的数值。以此类推，对模型中的其他方程进行同样的计算，从而得到整个系统在不同时间点的数值解。数值模拟的初始条件设定如下：假设在初始时刻t=0时，人群中没有感染登革热的个体，即I_h(0)=0，R_h(0)=0，易感人群数量S_h(0)=N_h（假设总人口数N_h=10000）；蚊子种群中，假设初始时没有感染的成虫，即I_m(0)=0，易感成虫数量S_m(0)=100，卵的数量E(0)=50，幼虫数量L(0)=30，蛹的数量P(0)=20。参数取值采用前文在模型构建部分所设定的值，即：\beta=150，表示成年雌性蚊子的出生率；\mu_E=0.3，为蚊子卵的死亡率；k_1=0.5，是卵发育为幼虫的速率；\mu_L=0.4，为蚊子幼虫的死亡率；k_2=0.15，是幼虫发育为蛹的速率；\mu_P=0.2，为蚊子蛹的死亡率；k_3=0.4，是蛹发育为成虫的速率；\mu_m=0.1，为蚊子成虫的死亡率；\beta_{hm}=0.3，是人传蚊的传播系数；\beta_{mh}=0.4，为蚊传人传播系数；\Lambda_h=0.01N_h，为人的补充率；\mu_h=0.01，为人的自然死亡率；\gamma_h=0.1，为人的康复率。5.2模拟结果展示通过数值模拟，我们得到了丰富的结果，这些结果以曲线的形式直观地展示了登革热传播过程中各变量的动态变化。不同参数条件下登革热发病率随时间变化的曲线清晰地呈现了疾病的传播趋势。当保持其他参数不变，仅改变蚊传人传播系数\beta_{mh}时，我们发现随着\beta_{mh}的增大，登革热发病率迅速上升，达到峰值的时间提前，且峰值更高。当\beta_{mh}=0.4时，发病率在第20天左右达到峰值，约为0.05；而当\beta_{mh}增大到0.6时，发病率在第15天左右就达到峰值，峰值约为0.08。这表明蚊传人传播系数对登革热的传播速度和严重程度有着显著影响，传播系数越大，病毒在人群中的传播越迅速，疫情越严重。蚊子各阶段种群数量随时间变化的曲线则揭示了蚊子种群的动态变化过程。蚊子卵的数量在初始阶段迅速增加，这是由于成年雌性蚊子的出生率较高，随着时间的推移，受到卵死亡率和发育为幼虫速率的影响，卵的数量逐渐趋于稳定。幼虫数量在卵孵化后逐渐上升，在第10天左右达到峰值，随后由于发育为蛹以及死亡率的作用，数量开始下降。蛹的数量变化趋势与幼虫相似，在幼虫大量发育为蛹后达到峰值，然后随着蛹羽化为成虫而减少。易感成虫和感染成虫的数量变化与登革热的传播密切相关，感染成虫数量在初期缓慢增加，随着登革热的传播，感染成虫数量迅速上升，在第25天左右达到峰值，随后由于成虫死亡率和病毒传播的动态平衡，数量逐渐稳定。基本再生数随关键参数变化的趋势图也为我们提供了重要信息。当改变人传蚊传播系数\beta_{hm}时，基本再生数R_0呈现出线性增长的趋势。当\beta_{hm}从0.1增加到0.5时，R_0从0.5左右增加到2.5左右，这表明人传蚊传播系数的增大直接导致了疾病在人群中的传播能力增强。而当改变蚊子成虫死亡率\mu_m时，基本再生数R_0呈现出下降的趋势。当\mu_m从0.05增加到0.2时，R_0从2左右下降到1左右，说明提高蚊子成虫死亡率可以有效降低登革热的传播能力。5.3结果分析与讨论通过对数值模拟结果的深入分析，我们可以清晰地探讨蚊子种群阶段结构变化对登革热传播的影响规律。当蚊子幼虫数量增加时，由于更多的幼虫能够发育为成虫，导致蚊子成虫种群数量上升。而成虫是登革热病毒的传播媒介，成虫数量的增加使得登革热传播的机会显著增多。在模拟中，将幼虫的初始数量从30增加到50时，感染成虫数量在后期明显增加，登革热发病率也相应提高，这表明控制幼虫数量对于预防登革热传播至关重要。蚊子成虫寿命延长（即降低成虫死亡率\mu_m）同样会对登革热传播产生重要影响。成虫寿命的延长使得它们有更多的机会叮咬人类，从而增加了病毒传播的几率。在模拟中，当\mu_m从0.1降低到0.05时，感染成虫数量在整个模拟过程中都维持在较高水平，登革热发病率也明显上升，这说明缩短蚊子成虫寿命是控制登革热传播的有效措施之一。针对不同防控措施对登革热传播的控制效果，我们也进行了深入分析。在灭蚊措施方面，通过降低蚊子的出生率（即减小\beta），可以有效减少蚊子种群数量，从而降低登革热的传播风险。在模拟中，将\beta从150降低到100时，蚊子各阶段种群数量都明显减少，感染成虫数量和登革热发病率也大幅下降。在疫苗接种方面，假设疫苗接种能够降低人群的易感性，即减小蚊传人传播系数\beta_{mh}。模拟结果显示，当\beta_{mh}从0.4降低到0.2时，登革热发病率显著降低，疫情得到有效控制。将这些模拟结果与理论分析结果进行对比验证，我们发现两者具有高度的一致性。在理论分析中，我们得出基本再生数R_0是判断登革热传播的关键指标，当R_0<1时，疾病会逐渐消除；当R_0>1时，疾病有可能爆发并持续传播。在数值模拟中，当通过调整参数使得R_0<1时，登革热发病率逐渐降低，最终趋近于零；而当R_0>1时，登革热发病率会上升并维持在一定水平，这充分验证了理论分析的正确性，也为我们基于模型制定科学的登革热防控策略提供了有力支持。六、防控策略研究6.1基于模型的防控策略制定根据前文对登革热传播动力学模型的分析以及数值模拟结果，我们可以制定一系列针对性的防控策略，以有效控制登革热的传播。针对蚊子不同阶段的特点，应采取相应的灭蚊措施。在蚊子的幼虫阶段，清理蚊虫孳生地是关键。由于幼虫主要在水中生长，因此要定期清理各类积水容器，如花盆托盘、花瓶、水桶等，确保无积水残留。对于无法清理的积水，如大型水池、下水道等，可以投放苏云金杆菌以色列亚种（Bti）等生物杀虫剂，这些生物杀虫剂对蚊子幼虫具有特异性毒性，能有效抑制幼虫生长，同时对环境和其他生物的影响较小。据研究表明，在孳生地投放Bti后，蚊子幼虫的密度可降低70%-80%。此外，也可以使用化学杀虫剂，如双硫磷、倍硫磷等，但需注意其对环境的潜在影响，严格按照使用说明进行操作，避免过度使用导致环境污染和蚊子产生抗药性。在蚊子成虫阶段，使用驱蚊剂是减少成虫与人类接触的有效方法之一。驱蚊剂中的有效成分，如避蚊胺（DEET）、驱蚊酯、派卡瑞丁等，能够在皮肤表面形成一层保护膜，干扰蚊子的嗅觉感知，从而避免被蚊子叮咬。研究显示，使用含有10%-30%避蚊胺的驱蚊剂，可提供长达4-6小时的有效防护。同时，安装纱窗、纱门等物理防护设施也能显著减少蚊子进入室内的机会。纱窗、纱门的网眼大小应适中，一般建议网眼直径不超过1毫米，这样既能有效阻挡蚊子进入，又能保证室内通风良好。此外，在室内放置电蚊拍、灭蚊灯等设备，也可以捕杀进入室内的蚊子，降低室内蚊子密度。在疫苗接种策略方面，应合理规划接种方案。首先，要确定最佳接种时间。由于登革热的传播具有季节性特点，在流行季节前1-2个月进行疫苗接种，能够使人群在疫情高发期前获得免疫力。例如，在东南亚地区，登革热流行季节通常为5-10月，那么在3-4月进行疫苗接种较为合适。其次，要明确接种人群。重点接种人群应包括儿童、老年人、医护人员以及在登革热流行地区工作或生活的人群。儿童和老年人由于免疫系统相对较弱，感染登革热后更容易出现重症病例；医护人员在救治患者过程中，感染风险较高；而在流行地区工作或生活的人群，暴露于病毒的机会较多。对于这些人群，应优先安排接种疫苗，以降低感染风险。同时，要根据疫苗的类型和接种程序，确保接种的规范性和有效性。不同类型的登革热疫苗，如减毒活疫苗、灭活疫苗、重组亚单位疫苗等，其接种程序和免疫效果可能存在差异，应严格按照疫苗说明书进行接种。6.2防控策略的成本效益分析对上述提出的防控策略进行成本效益分析，有助于评估不同策略的经济可行性，为实际防控决策提供科学的经济依据。灭蚊措施的实施成本主要包括杀虫剂的采购费用、人工喷洒费用以及清理孳生地所需的人力和物力成本。以化学杀虫剂双硫磷为例，其市场价格约为每升50-100元，对于一个面积为1000平方米的居民区，若每月进行一次喷洒，每次使用量为5升，仅杀虫剂费用每月就需250-500元，再加上人工喷洒费用，每月的灭蚊成本可能在1000-2000元左右。而清理孳生地的成本则主要取决于清理的范围和难度，对于一个中等规模的社区，每月清理孳生地的成本可能在500-1000元左右。然而，这些灭蚊措施带来的效益显著。通过有效控制蚊子种群数量，可大幅减少登革热的传播，从而避免大量的医疗费用支出。根据相关研究，每预防一例登革热病例，可避免医疗费用支出约500-1000元。若一个地区通过灭蚊措施成功预防了100例登革热病例，那么可避免的医疗费用支出就高达5-10万元。疫苗接种策略的成本包括疫苗的研发、生产、运输和接种费用。目前，登革热疫苗的研发成本高昂，一种新型登革热疫苗的研发费用可能高达数亿美元。在生产方面，由于疫苗生产工艺复杂，对质量控制要求严格，导致疫苗的生产成本也较高。例如，某款登革热疫苗的单价可能在50-100元/剂，对于一个10万人的接种人群，仅疫苗采购费用就需500-1000万元，再加上运输、储存和接种费用，总成本可能更高。但是，疫苗接种带来的效益同样可观。通过提高人群的免疫力，减少登革热的发病率，不仅可以避免医疗费用的支出，还能减少因患病导致的劳动力损失。据估算，每接种1000人，可减少登革热发病50-100例，按每人患病导致的劳动力损失为2000元计算，可减少劳动力损失10-20万元。为了更直观地评估不同防控策略的优劣，我们运用成本效益比（C