Modeling and control of Chikungunya with chronic infection

TL;DR

This study models the spread of Chikungunya fever in China, incorporating chronic infection, and identifies the most cost-effective control strategies.

Contribution

The study introduces a novel dynamical model with a chronic infection stage and evaluates cost-effective intervention strategies for Chikungunya.

Findings

The model estimates key transmission parameters and validates the stability of disease-free and endemic equilibria.

Personal protection is identified as the most economically efficient intervention strategy.

Mosquito biting rate and mortality rate are the most influential factors in disease transmission.

Abstract

Recognized globally as a major public health concern in the tropics and subtropics, Chikungunya fever also poses a potential epidemic risk in areas of China such as Guangdong Province, where suitable mosquito vector habitats exist. Based on a Chikungunya fever outbreak in Shunde District, Foshan City, this study develops a dynamical model incorporating a chronic infection stage. We derive R0 and perform a thorough stability analysis of all equilibria. Using daily reported case data from Shunde District, model fitting yields estimates for three key transmission parameters (β, ρ1, ρ2), the total mosquito population (Tv), and the initial number of infected mosquitoes (Iv(0)). Sensitivity analysis identifies that the primary positive and negative parameters on disease transmission are mosquito biting rate (β) and mosquito mortality rate (ϵv), respectively. Accordingly, five types of intervention measures are designed: personal protection, screening and detection, treatment of acute patients, management of chronic cases, and mosquito vector control measures. Based on these findings, we formulate a control framework to optimize intervention strategies. Numerical simulations not only validate the global asymptotic stability of the disease-free equilibrium when R0 < 1 and that of the endemic equilibrium when R0 > 1, but also assess the effectiveness of different control strategies. Strategy A, which emphasizes personal protection, emerges as the most economically efficient option in the cost-effectiveness analysis. It not only effectively interrupts virus transmission but also optimally reduces the burden of chronic cases, thereby offering a scientifically sound and economically feasible approach for public health resource allocation.

Figures

Key components and process of Chikungunya virus transmission.Fig. 1Table 2.1Summary of model parameters.Table 2.1SymbolBiological Meanings*B*_*h*_Susceptible individuals entry rate*B*_*v*_Susceptible mosquitoes entry rate*β*Mosquito biting rate*ρ*_1_Mosquito-to-human infection probability*ρ*_2_Human-to-mosquito infection probability*η*_*h*_Recovery rate for asymptomatic individuals*p*Proportion of cases that become symptomatic*θ*_*h*_Acute infected individuals recovery rate*σ*_*h*_Disease progression rate from acute to chronic phase*δ*_*h*_Chronic infected individuals recovery rate*κ*The ratio

Comparison of data before and after smoothing processing.Fig. 2

Presentation of fitting results.Fig. 3Table 3.5Model parameter estimates with 95 % confidence intervals.Table 3.5ParametersValue95 % CI*β*0.2111[0.0466,0.2411]*ρ*_1_0.1932[0.1715,0.3193]*ρ*_2_0.1118[0.0152,0.2232]*T*_*v*_6071859[1.4562 × 10^6^,7.7513 × 10^6^]*I*_*v*_(0)15536[1.3153 × 10^4^,4.5265 × 10^4^]*R*^2^0.76

Comparison chart between model prediction and actual situation.Fig. 4

Solution Trajectories of system (2.1) when *R*_0_ < 1.Fig. 5

Solution Trajectories of system (2.1) when *R*_0_ > 1.Fig. 6

PRCC of *R*_0_.Fig. 7

Population dynamics of infected subpopulations under constant control.Fig. 8

Strategy AFig. 9

Strategy B.Fig. 10

Strategy C.Fig. 11

Strategy DFig. 12

Full text

Introduction

1

Chikungunya virus (CHIKV), a single-stranded RNA Alphavirus within the Togaviridae family (Tang et al., 2023), is the etiological agent responsible for Chikungunya fever. This arboviral disease is primarily spread by Aedes mosquitoes and has a global distribution. Its cardinal clinical feature is characterized primarily by an acute high fever followed by intense arthralgia (Ortigoza et al., 2020). Chikungunya fever is rarely life-threatening, with a fatality rate of approximately 0.1 % (Pialoux et al., 2007). However, about 30–40 % of patients may develop chronic arthritis that persists for months or even years, resulting in significant loss of workforce capacity (Flandes et al., 2024; Milligan et al., 2019; Puntasecca et al., 2021; Van Aalst et al., 2017). The persisting sequelae and chronic symptoms of the disease place a substantial economic and societal burden on both individuals and communities (Cao et al., 2025). Specifically, in Brazil's private healthcare system, the per capita hospitalization cost for patients with chronic symptoms reaches as high as USD 2,400, accompanied by prolonged work absenteeism (Andrade et al., 2022). In India, Chikungunya outbreaks have resulted in an estimated productivity loss ranging from USD 2.57 million to USD 4.69 million due to absenteeism (Costa et al., 2023). The emergence of Chikungunya poses a serious public health threat, characterized by increased strain on healthcare systems, significant direct economic costs, as well as the deterioration of the socioeconomic order in endemic areas (Simon et al., 2011). Thus, the chronic stage represents a major contributor to the long-term disease burden of Chikungunya, making it imperative to incorporate this phase in relevant models.

Chikungunya fever is primarily distributed in tropical and subtropical regions. However, driven by factors such as climate change, globalization, and urbanization, the risk of importation and re-emergence of the virus into new geographical areas continues to increase (Kang et al., 2024). Since its initial identification in 1952, local transmissions have occurred in over 110 countries across the globe. Southern China, particularly Guangdong Province, is recognized as a high-risk area for both imported and locally transmitted cases of Chikungunya fever, due to the widespread distribution of Aedes mosquito vectors. The high population density and subtropical climate in Guangdong provide conducive conditions for mosquito proliferation, further amplifying transmission risk. In 2010, the first local outbreak in China was reported in Dongguan. Another outbreak occurred in Guangzhou in 2019. Most recently, on July 8, 2025, a locally acquired outbreak triggered by an imported case was detected in Shunde District of Foshan City, Guangdong Province. The epidemic subsequently spread to the entire Foshan municipality, demonstrating the disease's continuing public health threat to southern China.

Within epidemiology, mathematical models have become a core tool for analyzing epidemic patterns and predicting transmission trends, owing to their ability to quantitatively analyze complex transmission processes. A early compartmental model for Chikungunya transmission was constructed by Dumont et al. (2008) based on the outbreak in Réunion Island. Moulay et al. (2011, 2012) later proposed a model incorporating time-varying control. Through theoretical and numerical analysis, they demonstrated the cost-effectiveness advantage of dynamically adjusting control efforts in response to epidemic dynamics. In recent years, the research on Chikungunya fever has attracted extensive efforts from scholars exploring its transmission mechanisms and control strategies from diverse perspectives. Stochastic mathematical modeling of Chikungunya transmission (Cao et al., 2025; Fahimi et al., 2025), has been employed to analyze dynamic behaviors such as disease persistence and extinction through theoretical derivation. Furthermore, several studies have characterized the spread of Chikungunya via fractional-order differential equations (Jain & Chalishajar, 2023; Yangla et al., 2023), systematically investigating the advantages of fractional operators in capturing transmission characteristics and their application in optimizing public health strategies. Some studies have designed optimal intervention strategies by integrating antibody regulation or temperature effects with multiple control measures (El Hajji, 2021; Lusekelo et al., 2023). Another studies have developed models incorporating time-varying parameters, impulse control, nonlocal periodic reaction-diffusion, and periodic models accounting for temperature and rainfall influences, respectively (Liu & Stechlinski, 2015; Liu et al., 2020a, 2020b; Li & Zhao, 2024). These works explore how climatic factors, spatial heterogeneity, and asymptomatic infections influence the transmission mechanisms, epidemic characteristics, and control effectiveness of Chikungunya fever. Feng et al. (2019) and Wang et al. (2022) developed dynamic models from the perspective of CHIKV mutation and analyzed the impact of viral mutations, seasonal temperature variations, and other factors on transmission through real-world case studies and numerical simulations. Additionally, some models have incorporated factors such as chronic infection, asymptomatic infection, or climate influences to study control strategies in specific regions (González-Parra et al., 2019; Helikumi et al., 2022; Trentini et al., 2018; Vázquez-Peña et al., 2023). These studies provide a valuable foundation to comprehend the dynamics of CHIKV transmission and to develop effective, science-based control measures.

However, there is currently limited research focusing on the disease status of chronic infected individuals. In fact, viewing infectious disease dynamics as a multi-state process encompassing asymptomatic infection, a highly contagious acute phase, and a prolonged chronic phase has been established as a key theoretical framework for modeling such complex pathologies. In particular, research (Guo et al., 2012; Guo & Li, 2011; Liu & Li, 2021; Sanz-Lorenzo et al., 2025; Song et al., 2017) demonstrated that accurately characterizing the heterogeneity among individuals across different infection states and their transition dynamics is central to understanding the overall dynamics and long-term impact of an epidemic.

Therefore, although chronic Chikungunya patients are non-infectious, the preceding discussion underscores that their impact on society cannot be overlooked. Neglecting this phase would lead to an underestimation of the pandemic's overall socioeconomic impact and introduce bias into the cost-effectiveness evaluation of public health interventions. Incorporating the chronic stage into model not only expands the analytical scope from mere transmission dynamics to disease burden assessment but also provides critical insights for public health decision-making, thereby guiding rational resource allocation and promoting the design of optimized intervention strategies. Consequently, this paper addresses the Chikungunya fever outbreak in Foshan and aims to develop a model that incorporates chronic infection to simulate the bidirectional human-mosquito transmission process. Designed with a dual purpose, the model aims to replicate CHIKV transmission patterns and quantify the associated disease burden more comprehensively. Following the assessment of system's stability properties, parameter estimation is conducted using daily reported case data from Shunde District of Foshan, identifying key sensitive parameters influencing Chikungunya spread. Based on the findings, relevant control measures are proposed, and numerical simulations along with cost-effectiveness evaluation of intervention strategies are performed.

The paper proceeds as follows. The mathematical model is formulated and its stability is examined in Section 2. The subsequent section addresses the estimation of parameters and presents numerical results of stability properties. Turning to sensitivity analysis, Section 4 identifies the key parameters influencing R0 and disease spread. The optimal control framework is then analyzed in Section 5. Section 6 simulates various control measures and evaluates their economic efficiency. The final section offers a general conclusion and discussion of the study.

Mathematical model and stability analysis

2

Model formulation

2.1

We consider three infection states for human population and categorize them into five compartments: susceptible (Sh), asymptomatic infection (Ah), acute infection (Iac), chronic infection (Ich), and recovered (Rh). The total human population size is denoted by Th. This population experiences a growth rate Bh and a natural mortality rate ϵh. Owing to the low fatality rate of Chikungunya, additional mortality from infection is not considered in this model. The per-bite transmission risk from mosquito to human is ρ1. It is assumed that infected individuals with a proportion p exhibit symptoms and progress to the acute stage, while the remainder become asymptomatic infected individuals. Asymptomatic individuals recover after a certain period, with a recovery rate of ηh. In addition to recovering at rate θh, a fraction of individuals in the acute stage advance into the chronic stage at rate σh, and subsequently recover at rate δh. All recovered individuals acquire permanent immunity and are not susceptible to reinfection.

We classify the mosquito population into two compartments: susceptible (Sv) and infected (Iv), and the total population is denoted by Tv. Bv and ϵv are defined as the birth and mortality rates for the mosquito population, respectively. Since chronic carriers do not carry a sufficient viral load to infect mosquitoes, susceptible mosquitoes typically acquire infection through biting either asymptomatic or acutely infected human hosts. The infection rates from acutely infected and asymptomatic individuals to mosquitoes are given by ρ2 and κρ2, respectively.

Fig. 2.1 illustrates the transmission dynamics of Chikungunya fever. All parameters take on non-negative values. Table 2.1 summarizes the complete parameters and their biological interpretations. We can formulate the following ordinary differential equation (ODE) system:

$$[eqn]$$

Fig. 1. Key components and process of Chikungunya virus transmission.Fig. 1. Table 2.1Summary of model parameters.Table 2.1SymbolBiological MeaningsB_hSusceptible individuals entry rateBv_Susceptible mosquitoes entry rateβMosquito biting rateρ_1_Mosquito-to-human infection probabilityρ2_Human-to-mosquito infection probabilityηhRecovery rate for asymptomatic individualspProportion of cases that become symptomaticθhAcute infected individuals recovery rateσhDisease progression rate from acute to chronic phaseδhChronic infected individuals recovery rateκThe ratio of the transmission of asymptomatic infections to that of acute infectionsϵhHuman mortality rateϵv_Mosquito mortality rate

Non-negative initial conditions are assumed for all state variables: Sh(0) ≥ 0, Ah(0) ≥ 0, Iac(0) ≥ 0, Ich(0) ≥ 0, Rh(0) ≥ 0, Sv(0) ≥ 0, Iv(0) ≥ 0. The human and mosquito populations satisfy Th = Sh + Ah + Iac + Ich + Rh, and Tv = Sv + Iv.

Model analysis

2.2

According to the model assumptions, we have: $[eqn]$ , $[eqn]$ . It is straightforward to show that: $[eqn]$ , $[eqn]$ .

Hence, the stability properties of original model is equivalent to that of the limiting system:

$$[eqn]$$

As the variables Ich and Rh are absent from other equations of the system, the subsequent analysis is restricted to the reduced system corresponding to Sh, Ah, Iac, Sv, Iv:

$$[eqn]$$

Let $[eqn]$ , $[eqn]$ , $[eqn]$ , ϵa = ϵh + η, ϵac = ϵh + σh +θh, μ = ϵh, d = ϵv, d1 = ϵac, d2 = ϵa, B1 = Bh, B2 = Bv, q = 1 − p.

And We define the state variables as x = Sh, y = Sv, u = Iac, v = Ah, z = Iv.

Then the reduced system (2.3) takes the form:

$$[eqn]$$

The system's equilibrium E∗ = (x∗, y∗, u∗, v∗, z∗) satisfies

$$[eqn]$$

Solving these equations yields: $[eqn]$ , $[eqn]$ , $[eqn]$ , $[eqn]$ , $[eqn]$ .

Thus, the system has two equilibria:

• •The disease-free equilibrium (z∗ = 0): $[eqn]$ .

Computation of the basic reproduction number, based on the next-generation matrix (Van Den Driessche & Watmough, 2002), leads to the following expression:

$$[eqn]$$

Let $[eqn]$ , where $[eqn]$ , $[eqn]$ .

• •The endemic equilibrium (z∗ ≠ 0): $[eqn]$ , where

$$[eqn]$$

Theorem 2.1When R0 < 1*, the disease-free equilibrium* E0 of system (2.4) is globally asymptotically stable.

Proof. Consider the Jacobian matrix derived from the linearization of system (2.4) at E0:

$$[eqn]$$

Then, $[eqn]$ , where, $[eqn]$ . The condition R0 < 1 ensures that all roots of D(λ) = 0 possess negative real parts, thereby ensuring the local asymptotic stability of E0.

Let P(x) = x − 1 − ln x. Consider the following Lyapunov function:

$$[eqn]$$

where $[eqn]$ , and

$$[eqn]$$ $$[eqn]$$

It follows that $[eqn]$ .

We obtain $[eqn]$ , and $[eqn]$ , if and only if x = x0, y = y0, u = v = 0. The global asymptotic stability of the disease-free equilibrium is established via LaSalle's invariance principle (LaSalle, 1976).Theorem 2.2When R0 > 1*, the endemic equilibrium* E∗ is globally asymptotically stable.

Proof. Consider the Jacobian matrix derived from the linearization of system (2.4) at E∗:

$$[eqn]$$

We denote k = β2u∗ + β3v∗. The characteristic equation is $[eqn]$ , where, $[eqn]$ .

From pβ1x∗z∗ = d1u∗, qβ1x∗z∗ = d2v∗, β2u∗y∗ + β3v∗y∗ = dz∗, $[eqn]$ , we can obtain

$$[eqn]$$

where $[eqn]$ , $[eqn]$ . We have

$$[eqn]$$

Therefore, we obtain:

$$[eqn]$$

where

$$[eqn]$$

An application of the Hurwitz criterion confirms that all roots of $[eqn]$ possess negative real parts, implying that the same holds for all roots of $[eqn]$ . This satisfies the criteria for the local asymptotic stability of E∗.

Next, define a Lyapunov function as follows:

$$[eqn]$$

where $[eqn]$ , $[eqn]$ . Then

$$[eqn]$$ $$[eqn]$$ $$[eqn]$$

Therefore,

$$[eqn]$$

This shows that $[eqn]$ and $[eqn]$ if and only if x = x∗, y = y∗, u = u∗, v = v∗, z = z∗. An application of LaSalle's invariance theorem (LaSalle, 1976) demonstrates the endemic equilibrium is globally asymptotically stable.

Parameter estimation and stability simulation

3

Fixed parameters and data

3.1

The recent Chikungunya fever outbreak in Foshan City remained primarily concentrated in the Shunde District. Therefore, we collect the daily reported number of new Chikungunya cases in Shunde District from July 19 to August 11, 2025 from Foshan Municipal Health Bureau, as shown in Table 3.2. According to information provided by the Foshan Municipal Bureau of Statistics, the total population of Shunde District in 2025 is Th = 3268600, with a mean life expectancy of 83.6 years. From this, the mortality rate can be estimated as $[eqn]$ . Under the assumption of equal birth and death rates, we have Bh = 3268600 × ϵh ≈ 107. Although asymptomatic infected individuals exhibit viremia, their viral load is significantly lower than that of symptomatic infected individuals (Appassakij et al., 2013). Consequently, their efficiency in vector-borne transmission to mosquitoes is reduced and then we estimate that κ = 0.3. The duration of acute infection is typically around 7 days (Rama et al., 2024), hence the removal rate from the acute infection compartment is estimated to be 0.143. Given that approximately 30–40 % of patients may progress to become chronically infected, we set θh = 0.086 and σh = 0.057. Accordingly, we assign the fixed biological parameters based on literature sources, including key epidemiological parameters such as total population size (Th), population growth rate (Bh), and clinical progression rate (σh, θh, ηh, δh). Table 3.3 lists the specific parameter values used in the model.Table 3.2New daily cases in Shunde district.Table 3.2Jul 19Jul 20Jul 21Jul 22Jul 23Jul 24Jul 25Jul 26629368313463383310273310Jul 27Jul 28Jul 29Jul 30Jul 31Aug 1Aug 2Aug 3323357408362312258192158Aug 4Aug 5Aug 6Aug 7Aug 8Aug 9Aug 10Aug 11135122115111102987768Table 3.3Fixed parameters and their values.Table 3.3SymbolValuesReferenceη_h0.167Feng et al., 2019Feng et al. (2019)p0.86Feng et al., 2019Feng et al. (2019)σh0.057See textθh0.086See textδh0.003González-Parra et al., 2019González-Parra et al. (2019)κ0.3See textϵv0.07González-Parra et al., 2019González-Parra et al. (2019)Th3268600See textBh107See textϵh_0.000033See text

Parameter fitting

3.2

This study employs a parameter estimation framework based on the least squares method to fit the daily new case data within the MATLAB environment. The framework is used to estimate three key parameters—the mosquito biting rate (β), the human infection probability (ρ1), and the mosquito infection probability (ρ2), along with the total mosquito population (Tv) and its initial prevalence of infection (Iv(0)). To comprehensively evaluate the model's goodness of fit and predictive capability, the collected daily case data are divided into two segments: the first 17 days for parameter estimation, and the subsequent 7 days for validating the model's predictive performance.

The original epidemic data may be subject to reporting delays and random fluctuations. To improve the quality of model fitting, this study applies the Weighted Moving Average (WMA) method for data smoothing. Specifically, a 5-day sliding window is employed, assigning higher weight to recent data to reduce random noise while preserving essential data characteristics. The weight distribution is generated using an exponential function with decay parameter λ = 1.7. The weight vector ω = [ω1, ω2, ω3, ω4, ω5] is computed as follows:

$$[eqn]$$

For the first five data points (t ≤ 5), a simple arithmetic mean is applied. For t > 5, the smoothed value $[eqn]$ is calculated as:

$$[eqn]$$

where C(t) denotes the observed value at time t.

As shown in Fig. 3.2, the comparison between raw and smoothed data demonstrates that the smoothed curve retains the overall trend of the outbreak development while significantly reducing short-term fluctuations and noise.Fig. 2. Comparison of data before and after smoothing processing.Fig. 2

Based on the characteristics of the early outbreak phase, initial values for the five epidemiological state in humans on July 19th are estimated as follows: susceptible individuals (Sh), asymptomatic infected individuals (Ah), symptomatic infected individuals in the acute stage (Iac), chronically infected individuals (Ich), and recovered individuals (Rh). The assigned values are: Ah = 163, Iac = 1000, Ich = 10, Rh = 151, Sh = Th − Ah − Iac − Ich − Rh. Using the smoothed data, parameter optimization is performed via the nonlinear least squares method starting from day 1 (corresponding to July 1st). The lsqnonlin solver (based on Levenberg–Marquardt algorithm) is employed to minimize the residual sum of squares between simulated and actual values. The model yields a numerical solution for the daily new cases, denoted as $[eqn]$ , computed by:

$$[eqn]$$

resulting in the estimated values $[eqn]$ , $[eqn]$ , …, $[eqn]$ . The objective function is formulated as the sum of squared errors (SSE) across all predicted and observed data points:

$$[eqn]$$

The coefficient of determination R^2^ is used to measure how closely the model's predictions aligned with the observations. It quantifies the proportion of variance explained by the model and is computed as:

$$[eqn]$$

where $[eqn]$ denotes the mean of the observed daily case counts over the fitting period.

To ensure biologically plausible and numerically stable parameter estimation, search bounds are imposed on each parameter during optimization. The specific ranges provided in Table 3.4 are determined through preliminary fitting and practical rationale to enhance optimization efficiency by restricting the search to a realistic parameter space. We employ a bootstrap resampling method to quantify parameter estimation uncertainty by computing 95 % confidence intervals. Based on 1000 resamplings of the model residuals, the confidence intervals are defined by the 2.5th and 97.5th percentiles of the resulting empirical distribution.Table 3.4Parameters for estimation and their range.Table 3.4SymbolValue rangeSymbolValue rangeβ[0.01,0.7]ρ1[0.01,0.7]ρ1[0.01,0.7]Tv[1 × 10^6^,2 × 10^7^]Iv(0)[1,1 × 10^5^]

Fig. 3.3 presents the fitting results, with the estimated parameter values and goodness-of-fit metrics provided in Table 3.5. All parameters have plausible 95 % confidence intervals. The results suggest that the Chikungunya transmission model developed in this study achieves a moderate level of fitting accuracy (R^2^ = 0.76). The maximum residual is approximately 100 cases, suggesting that the model did not fully capture all variations in the epidemic time series. This is common in complex mosquito-borne disease systems, as transmission is influenced by numerous unobserved stochastic factors, such as human mobility and fluctuations in case reporting rates. Nevertheless, the model successfully captures the core dynamic trend of the outbreak and identified intrinsic driving forces of disease spread. The residuals are randomly distributed over time, and smaller residuals in later stages indicate that the model mechanism is relatively reliable during this period. Overall, the model demonstrates satisfactory fitting result, successfully quantifies key transmission parameters, and provides valuable mechanistic insights into the epidemic dynamics. Its limitations also offer clear directions for future model improvements.Fig. 3. Presentation of fitting results.Fig. 3. Table 3.5Model parameter estimates with 95 % confidence intervals.Table 3.5ParametersValue95 % CIβ0.2111[0.0466,0.2411]ρ1_0.1932[0.1715,0.3193]ρ2_0.1118[0.0152,0.2232]Tv6071859[1.4562 × 10^6^,7.7513 × 10^6^]Iv(0)15536[1.3153 × 10^4^,4.5265 × 10^4^]R^2^0.76

As shown in Fig. 3.4, although the actual observed values during the validation period (the last 7 days) mostly fall within the 95 % confidence interval of the model predictions, they consistently cluster near the lower boundary of the prediction interval. This non-random, systematic pattern suggests that while the model captures the overall declining trend of the outbreak, it may have systematically overestimated the intensity of viral transmission during the validation period. We hypothesize that this deviation likely reflects key differences in transmission dynamics between the model calibration period (the first 17 days) and the validation period. A plausible explanation is that around the beginning of the validation period, external mitigating factors begin to take significant effect. These include targeted mosquito control measures implemented after the peak transmission period, which reducing mosquito density, and increased public risk awareness leading to widespread individual protective behaviors (such as repellent use and reducing exposure). Together, these factors effectively suppress the actual transmission level, keeping it consistently lower than the model predictions based on the natural transmission dynamics observed in the earlier phase. Furthermore, according to actual reports, Foshan City activated an emergency response mechanism immediately after the outbreak was detected, focusing on three key interventions: environmental sanitation management, breeding site cleanup, and adult mosquito elimination. A targeted campaign for adult mosquito eradication was launched on July 31, aiming to comprehensively reduce mosquito populations. This real-world context aligns closely with our inferred causes. The systematic deviation of observed values below the predictions provides clear visual and quantitative evidence of the substantial effectiveness of the aforementioned interventions and behavioral changes in the population.Fig. 4. Comparison chart between model prediction and actual situation.Fig. 4

Stability simulation

3.3

This subsection presents numerical simulations that verify the stability of equilibrium points in system (2.1). We employ the ode45 solver in MATLAB to numerically solve this system of ordinary differential equations. The initial values are chosen as follows: Sh(0) = 3267276, Ah(0) = 163, Iac(0) = 1000, Ich(0) = 10, Sv(0) = 6056359, Iv(0) = 15536. For the given initial conditions and parameter set, we calculate R0 = 0.40 < 1. The stability of the equilibrium of system (2.1) can be observed in Fig. 3.5, which depicts the decay of both human and mosquito infections to zero and eventual stabilization. This result validates Theorem 2.1, indicating that disease elimination is the eventual outcome.Fig. 5. Solution Trajectories of system (2.1) when R0 < 1.Fig. 5

Next, we assume climatic conditions become favorable, leading to an increased mosquito reproduction rate and extended survival time. Accordingly, we set ϵv = 0.04, Bv = 2091904, while keeping all other parameters and initial values unchanged. Under this scenario, the basic reproduction number is R0 = 1.01 > 1. We observe from Fig. 3.6 that the equilibrium of system (2.1) is stable, which displays the numbers of asymptomatically infected individuals, acutely infected individuals, and infected mosquitoes stabilize at positive values, indicating disease persistence. This result provides validation for Theorem 2.2.Fig. 6. Solution Trajectories of system (2.1) when R0 > 1.Fig. 6

Compared with relevant Chikungunya studies, our estimate of R0 (0.40) in this study is slightly lower. For instance, González-Parra et al. (2019) calculated a range of 0.74–0.84 for Colombia's 2015 epidemic, while Feng et al. (2019) estimated an R0 of approximately 2.035 for the 2007 Chikungunya outbreak in Italy. In contrast, Li and Zhao (2024) obtained a significantly higher estimate of R0 = 3.6087 in Brazil. We primarily attribute this relatively low estimate to two main factors. First, the model is calibrated using data from the second phase of the outbreak (Zhang et al., 2025), during which intervention measures such as case isolation and public health education have already been implemented in Shunde District, significantly reducing transmission risk. Second, owing to its extensive experience with recurrent Dengue fever epidemics, Guangdong Province has established robust surveillance and control systems, which enabled the timely and effective containment of this Chikungunya outbreak. This finding underscores the critical importance of timely interventions in controlling outbreak dynamics.

Sensitivity analysis

4

Local sensitivity indices

4.1

This section employs local sensitivity methods to investigate the impact of parameter variations on R0 and to determine which parameters with significant influence. Sensitivity indices provide a quantitative measure of local parameter influence and are computed according to Zi (2011):

$$[eqn]$$

The sensitivity indices are defined as follows for each parameter:

$$[eqn]$$

Notably, a number of sensitivity indices exhibit dependence on one or more parameters, while the remaining indices are constants independent of parameter values. Table 4.6 lists their values calculated with the baseline parameters. Among these, the mosquito mortality rate ϵv has the greatest negative effect; conversely, the biting rate β has the strongest positive effect.Table 4.6Sensitivity indicators of R0.Table 4.6ParameterIndex valueParameterIndex valueβ1p0.286Bh−0.5Bv0.5ρ1_0.5ρ2_0.5σh−0.153κ0.016ϵh0.401ϵv−1θh−0.231ηh−0.016

Global sensitivity analysis

4.2

Unlike local approaches, global sensitivity analysis evaluates simultaneous variations across a broader parameter space. This is accomplished by computing the Partial Rank Correlation Coefficient(PRCC) (Blower & Dowlatabadi, 1994; Marino et al., 2008) for each parameter, which are sampled using the Latin Hypercube Sampling (LHS) method.

Fig. 4.7 presents the PRCC analysis result for R0, revealing that parameters β, ρ1, ρ2, p, ϵh, κ and Bv are positively correlated with R0, while Bh, ηh, σh, θh and ϵv are negatively correlated with R0. The analysis identifies the mosquito biting rate (β) and mortality rate (ϵv) as the dominant influences on R0, exhibiting strongly positive and negative effects, respectively. The PRCC values for both parameters close to 1, indicating that they have a highly significant influence on R0.Fig. 7PRCC of R0.Fig. 7

Both local and global sensitivity analyses demonstrate that the model output is most sensitive to the mosquito biting rate (β) and mortality rate (ϵv). This result indicates that reducing the mosquito biting rate and enhancing the mosquito mortality rate are two key intervention strategies for controlling the spread of Chikungunya fever. The findings not only clarify the priority directions for prevention and control measures, but also provide a direct basis for the introduction of time-dependent intervention measures in the subsequent optimal control analysis.

Optimal control

5

In response to the prevalence of Chikungunya fever in Foshan City, Shunde District has established a multi-level prevention and control system centered on “surveillance and early warning—blocking transmission sources—comprehensive vector management”. The system first enhances case screening in medical institutions and dynamic monitoring of mosquito vector density to achieve early warning. It also enforces mosquito-proof isolation of confirmed cases and ensures prompt management of infection hotspots. For vector intervention, a strategy emphasizing “environmental management as the primary approach, supplemented by chemical control” is adopted. This involves mobilizing the public to eliminate breeding sites and integrating the use of biological and chemical larvicides along with adult mosquito control measures. Furthermore, the system promotes public health education, establishes a multi-sectoral joint prevention and control mechanism, and implements risk-based zoning and graded management strategies.

Building upon this framework and considering the results of sensitivity analysis, we identify parameters β and ϵv as the dominant factors influencing the transmission of Chikungunya fever. Accounting for the management needs of chronic cases, this study aims to present a dynamic optimal control framework to optimize the allocation and scheduling of limited public health resources among key intervention measures in Shunde District. Therefore, the following five time-varying control measures are introduced:

• •u1(t): Enhance public awareness campaigns to disseminate knowledge on prevention and control, urging all residents to adopt personal safeguard strategies such as donning long-sleeved garments and employing mosquito netting, and installing window screens to reduce mosquito bites;
• •u2(t): Screen for asymptomatic infections; once identified, isolate such individuals. Since asymptomatic individuals show no apparent adverse symptoms, isolation prevents them from participating in the transmission process;
• •u3(t): For acute infections, adopt interventions including acupoint patches, topical medicinal applications, and traditional Chinese medicine treatments to shorten the duration of the acute phase and alleviate symptoms such as fever, joint pain, and rash;
• •u4(t): Strengthen treatment for chronic infections to mitigate joint symptoms and shorten the chronic phase;
• •u5(t): Implement mosquito eradication campaigns to achieve comprehensive vector elimination.

The following system represents the model incorporating combined intervention strategies:

$$[eqn]$$

The five control measures are denoted by U(t) = (u1(t), u2(t), u3(t), u4(t), u5(t)). The goal of optimal control is to determine a set of optimal controls U∗(t) that minimizes the number of infected individuals (including asymptomatic, acute, and chronic cases), while also minimizing the implementation costs. Accordingly, the corresponding objective function takes the form:

$$[eqn]$$

where the integrand is given by $[eqn]$ . Here K1, K2 and K3 represent the weighting coefficients for the asymptomatic (Ah), acute (Iac), and chronic (Ich) cases, respectively; while Ti(i = 1, …, 5) quantify the implementation costs of the five control measures ui(t).Theorem 5.1For the control model 5.5*, there exists an optimal control* $[eqn]$ such that

$$[eqn]$$

Proof. For an optimal control to exist, the following criteria must be satisfied: (1) According to Lukes’ result (Lukes, 1982), all control and state variables take on non-negative values; (2) The admissible control set U is closed and convex; (3) The system is bounded, ensuring the compactness required for the existence proof. (4) Since the integrand is a quadratic function in the control variables, the objective functional exhibits convexity over the admissible set U; (5) Moreover, owing to the fact that

$$[eqn]$$

We can find constants b > 1, a1 > 0, and a2 > 0 satisfying:

$$[eqn]$$

Theorem 5.2Let $[eqn]$ be the optimal control, $[eqn]$ be the optimal state solution. We formulate the Hamiltonian function:

$$[eqn]$$

There are adjoint variables λ(t) = (λ_1_(t), λ_2_(t), λ_3_(t), λ_4_(t), λ_5_(t), λ_6_(t), λ_7_(t)) satisfying

$$[eqn]$$ $$[eqn]$$

With the transversality conditions given by λ_i_(T) = 0(i = 1, …, 7), the resulting optimality condition for the control is:

$$[eqn]$$

Proof. Following Pontryagin's maximum principle (Pontryagin, 2018) for the defined Hamiltonian, we obtain the adjoint variable vector $[eqn]$ . The adjoint system can be derived through $[eqn]$ . By applying the necessary conditions for optimality, we obtain:

$$[eqn]$$

The optimal control solutions are obtained by solving $[eqn]$ :

$$[eqn]$$

Numerical simulations

6

Control strategy simulation

6.1

To evaluate the static effect of fixed-intensity intervention measures, we assigned three distinct sets of constant values to the control measures U: (1) No control: ui = 0 for i = 1, …, 5; (2)Moderate control: u1 = u2 = u3 = 0.2, u4 = 0.15, u5 = 0.1; (3)Strong control: u1 = u2 = u3 = 0.5, u4 = 0.3, u5 = 0.2. Under varying control intensities, Fig. 6.8 illustrates the comparative trajectories of human infection states and infected mosquito populations.Fig. 8. Population dynamics of infected subpopulations under constant control.Fig. 8

As shown in Fig. 6.8, increasing the control intensity leads to a marked bidirectional suppression effect within the human–mosquito Chikungunya transmission system. In the uncontrolled scenario, asymptomatic and acute human infections attain peaks of 332 and 2252, respectively. Under moderate control intensity, these peaks decline to 167 and 1392, with substantial reductions also observed in the number of chronic infections and infected mosquitoes. With intervention intensity further increased to a strong level, all infected human compartments and infected mosquito continue to decline from the medium control baseline, but the extent of reduction is less pronounced. This indicates that the declining trend does not follow a simple linear relationship with increasing control intensity.

Constant control represents a simplified abstraction of practical disease control measures, with its outcomes serving as a benchmark for evaluating more complex dynamic control strategies. We next consider the optimal control strategy. To explore temporal allocation of dynamically optimized interventions, we employ a forward-backward sweep algorithm that relies on the Runge-Kutta (Lenhart & Workman, 2007) numerical scheme within the MATLAB environment. In the numerical experiments, the control period is set to 50 days, with K1 = 30, K2 = 40, K3 = 20, T1 = 500, T2 = 900, T3 = 1000, T4 = 1500, T5 = 1600. The cost coefficients are assigned based on the complexity and resource intensity of each intervention, with the specific values being hypothetical. This costing structure captures their relative resource use. Public health education, focused on information dissemination, was assigned the lowest value. Asymptomatic screening, requiring testing kits and personnel, received a moderate value. Treatments for acute and chronic patients involve healthcare resources, with chronic treatment incurring greater expense due to its prolonged nature. Comprehensive mosquito spraying is assigned the highest cost, as it demands chemicals, equipment, extensive labor, and accounts for ecological expenses.

First, we consider Strategy A (u1, u4, u5), which emphasizes personal protection to reduce mosquito bites, comprehensive mosquito elimination, and chronic cases management. This strategy primarily focuses on infection prevention measures. As shown in Fig. 6.9, under the implementation of preventive measures such as avoiding mosquito bites and vector elimination, the peaks of asymptomatic and acute infections are significantly reduced, and the decay rate of infected mosquitoes increases. Furthermore, with the addition of chronic case management, the number of chronic infections also decreases substantially. This outcome aligns with the public health principle that prevention is superior to treatment. The final graph in Fig. 6.9 displays the time-varying profiles of the control measures. Both u1 and u5 are deployed at peak intensity for initial periods of 26 and 30 days, respectively, before a gradual reduction. In contrast, the chronic case management measure u4 reaches its peak around day 4, slightly later than the other two controls.Fig. 9. Strategy AFig. 9

We now focus on Strategy B (u1, u2, u3, u4), which emphasizes mosquito-bite prevention, medical testing, and treatment, to evaluate its effectiveness and cost in reducing the spread of Chikungunya (see Fig. 10). Findings confirm that the strategy effectively reduced the peak and overall scale of human infections. However, its suppressive effect on infected mosquitoes is not significant. Additionally, all control measures need to be implemented at maximum intensity for an extended duration, resulting in relatively high costs. This suggests that relying solely on physical protection and medical management, without mosquito elimination measures, still poses certain challenges in controlling the spread of Chikungunya.Fig. 10. Strategy B.Fig. 10

Therefore, we further incorporate comprehensive mosquito elimination into the strategy, denoted as Strategy C (u1, u2, u3, u4, u5). A key benefit of the strategy is the drastic reduction in the high-intensity control period, achieved through its dual impact: a sharp decline in human infections and a markedly accelerated decay of the infected mosquito population, leading to more effective disease containment (see Fig. 11). These findings underscore the importance of mosquito elimination measures in the containment of Chikungunya transmission.Fig. 11. Strategy C.Fig. 11

Similarly, building upon Strategy C, we further adjust the approach by reducing the maximum intensity of comprehensive mosquito elimination by half, resulting in an integrated control Strategy D (u1, u2, u3, u4, u5). This adjustment addresses environmental concerns from overusing control agents, and to evaluate how intervention intensity influences overall efficacy. It can be observed that reducing the mosquito elimination intensity leads to a prolonged duration of other control measures being applied at their maximum strength (see Fig. 12). While this lowers environmental cost, it increases the total control expenditure. Hence, in practical implementation, appropriately raising the intensity of mosquito control could not only reduce overall costs but also rapidly achieve more effective suppression of Chikungunya transmission. Such a strategy will help balance control efficacy, economic cost, and environmental impact.Fig. 12. Strategy DFig. 12

Cost-effectiveness analysis

6.2

When evaluating the practical value of disease prevention and control measures, focusing solely on their efficacy in containing disease spread does not suffice to form a complete basis for decision-making. The implementation of any control measure inevitably requires corresponding resource allocation. Excessively high costs may not only hinder the feasibility of interventions but also divert resources from other public health priorities, thereby compromising overall control effectiveness. Therefore, both the effectiveness and economic efficiency of control measures must be considered. We employ the Incremental Cost-Effectiveness Ratio (ICER) (Okosun et al., 2013) to further quantify investment costs and resulting benefits of different control strategies, which is calculated as follows:

$$[eqn]$$

The strategy cost is computed by the integral: $[eqn]$ , reflecting the cumulative expenditure across all control measures over the entire period. We compute infections averted by the strategy as: $[eqn]$ , where I(t) denotes infections without intervention and $[eqn]$ represents infections under the corresponding control strategy. Table 6.7 summarizes the model simulation results, presenting strategies ranked by ascending infections averted alongside their computed ICER values.Table 6.7Cost-effectiveness comparison of strategies A-DTable 6.7StrategyInfection AvertedCostICERCost of Chonic Case ManagementB43465.1681325.621.8724333.20A51826.8741969.84−4.7122186.76D53700.7672158.4416.1123162.24C56509.8765210.89−2.4722400.98

According to the computational results presented in Table 6.7, both ICER(A) and ICER(C) are negative. Comparing Strategy A with Strategy B, and Strategy C with Strategy D, it is evident that Strategy A and Strategy C incur relatively lower costs while averting a larger number of infections, demonstrating superior effectiveness and cost-effectiveness. Therefore, only Strategy A and Strategy C need to be compared in the subsequent analysis.

As clearly shown in Table 6.8, ICER(A) < ICER(C), indicating that Strategy C incurs higher costs while yielding inferior containment effectiveness. It can therefore be concluded that Strategy A achieves the highest cost-effectiveness. In other words, under resource constraints, implementing control measures such as personal protection and mosquito elimination can not only effectively curb the transmission of Chikungunya fever but also reduce the overall cost of control efforts. A comparison of the costs of chronic case management across the various strategies reveals that Strategy A required the lowest expenditure. From the perspective of mitigating long-term chronic disease burden, it also proves to be a highly favorable strategy.Table 6.8Cost-effectiveness comparison of strategies A and C.Table 6.8StrategyInfection AvertedCostICERCost of Chonic Case ManagementA51826.8741969.840.8122186.76C56509.8765210.894.9622400.98

Conclusion and discussions

7

In response to the Chikungunya outbreak in Foshan City, this paper establishes a cross-infection differential equation model between humans and mosquitoes that incorporates a chronic infection compartment in the human population. We analysis the system's dynamics and derive optimal control strategies, yielding the following primary findings:

Firstly, an explicit expression for R0 is obtained. Stability analysis demonstrates that when R0 < 1, the disease-free equilibrium is globally asymptotically stable, leading to disease extinction. Conversely, if R0 > 1, the endemic equilibrium is globally asymptotically stable, indicating that the Chikungunya epidemic may persist. These results are further validated through numerical simulations.

Secondly, the model parameters are fitted using the daily reported Chikungunya case data from Shunde District, yielding a good fit (R^2^ = 0.76). Key transmission parameters are successfully quantified: the mosquito biting rate is estimated as β = 0.2111, the human infection probability per bite as ρ1 = 0.1932, the mosquito infection probability as ρ2 = 0.1118, the total mosquito population as Tv = 6071859, and the initial infected vector count as Iv(0) = 15536. The constructed model accurately captures the actual transmission dynamics of the Chikungunya virus in Shunde District. Although some deviation between the observed data and model predictions is observed during the validation period, this discrepancy is explainable and intuitively demonstrates the success of intervention measures for disease containment. The basic reproduction number, calculated from the fitted parameters, is R0 = 0.40. This indicates a relatively low transmission potential of Chikungunya in Shunde District, which is consistent with the observed downward trend in daily cases. It also reflects the significant impact of timely mosquito control and prevention measures implemented by local authorities since the emergence of the outbreak.

Subsequently, local and global sensitivity analyses reveal that parameters including mosquito biting rate (β) and mosquito mortality rate (ϵv) exert a pronounced influence on R0. This finding underscores that interventions targeting mosquito vectors represent the most efficient approach for controlling Chikungunya outbreaks, thereby providing a theoretical foundation for the formulation of preventive measures.

Finally, incorporating the containment measures implemented in Shunde District, this study introduces five intervention strategies: personal protection, screening and testing, treatment of acute cases, treatment of chronic cases, and comprehensive mosquito vector eradication, to simulate their containment effectiveness against the transmission of Chikungunya. Analysis under both constant control and optimal control frameworks indicates that higher control intensities lead to better suppression of Chikungunya spread, though the relationship is not purely linear. Beyond a certain level of control intensity, the marginal gain in effectiveness diminishes, while the associated costs continue to rise. The optimal control strategy A is dominant in terms of cost-effectiveness, implementing this strategy can effectively reduce the peak number of acute infections by 50.69 % while also mitigating the potential long-term burden of chronic diseases. The results further emphasize that prioritizing preventive measures aimed at curtailing mosquito density and human-vector contact rates can rapidly contain outbreaks, lower the transmission risk of Chikungunya, and optimize short-term cost efficiency.

This study conducts parameter fitting based on the Chikungunya outbreak in Shunde District, Foshan City, thereby enhancing the regional specificity and practical guidance of the conclusions. Through the analysis of the model fitting presented earlier, we have also identified objectives for future work. For instance, the model could be further refined by incorporating seasonal dynamics, integrating local meteorological and vector surveillance data to enable coupled analysis of long-term climate change impacts. Additionally, efforts could focus on parameterizing external suppression factors, such as correlating the mosquito mortality rate (ϵv) with vector control intensity data, to develop a next-generation model that dynamically responds to public health interventions and exhibits enhanced predictive capability. Beyond these, an impulsive control strategy for mosquito elimination could be introduced to quantitatively simulate the effects of Foshan's vector control campaigns.

CRediT authorship contribution statement

Yan Wang: Writing – review & editing, Writing – original draft, Supervision, Software, Methodology, Investigation, Conceptualization. Huan Ma: Software, Data curation, Conceptualization. Qian Yan: Validation, Methodology, Formal analysis. Zhichun Yang: Writing – review & editing, Supervision, Investigation, Conceptualization.

Funding support

This work was supported by the 10.13039/501100001809National Natural Science Foundation of China (Nos. 12101513, 12471153 and 1197108), 10.13039/100010338Chongqing Normal University Foundation Program (No. 24XLB001), Project of Science and Technology Research Program of 10.13039/501100007957Chongqing Municipal Education Commission (No. KJQN202500513), Chongqing Talents Program (No. cstc2024ycjh-bgzxm0046), Higher Education Research Planning Project of China Association of Higher Education (No. 24SX0305), the “Unveiling and Leading” Project of Chongqing Municipal Commission of Economy and Information Technology (No. YJX-2025001001004).

Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.