Modelling the Periodic Outbreak of Measles in Mainland China

In mainland China, measles infection reached the lowest level in 2012 but resurged again after that with a seasonally fluctuating pattern. To investigate the phenomenon of periodic outbreak and identify the crucial parameters that play in the transmission dynamics of measles, we formulate a mathematical model incorporating periodic transmission rate and asymptomatic infection with waning immunity. We define the basic reproduction number as the threshold value to govern whether measles infection dies out or not. Fitting the reported measles cases from 2013 to 2016 to our proposed model, we estimate the basic reproduction number R0 with immunization to be 1.0077. From numerical simulations, we conclude asymptomatic infection does not cause much new infections and the key parameters affecting the transmission of measles are vaccination rate, transmission rate, and recovery rate, which suggests the public to enhance vaccination and protection measures to reduce effective contacts between susceptible and infective individuals and treat infected individuals timely. To minimize the number of infected individuals at a minimal cost, we formulate an optimal control system to design optimal control strategies. Numerical simulations show the effectiveness of optimal control strategies and recommend us to implement the control strategies as soon as possible. In particular, enhancing vaccination is especially effective in lowering the initial outbreak and making disease recurrence less likely.


Introduction
Measles is an acute, viral infectious disease characterized by high fever, cough, and a maculopapular rash [1]. In 1980, 2.6 million people died from measles. In 2005, the World Health Organization Western Pacific Region set a goal to eliminate measles before 2012. For achieving it, China has taken some immunization measures including continuing a two-dose routine measles vaccine strategy (first and second dose at 8 and 18-23 months) while using supplementary immunization activities (SIAs) [2]. Because of these measures, the number of reported measles cases significantly decreased from approximately 140000 in 2008 to 6678 in 2012. However, in 2013, a resurgence of measles occurred and measles cases obviously increased to 30000 and retained at around 28000 in 2016 with obvious periodic feature. Consequently, measles still remains a serious public health problem in mainland China and how to eliminate it becomes challenging.
Mathematical models play an essential role in investigating the transmission dynamics of the infectious diseases and helping to identifying the key parameters or process that significantly affect disease spread. A number of mathematical models have been formulated to study the transmission and control of measles infection. Hamer in 1906 [3] proposed a model to study the regular occurrence of measles in London. en, a prediction model was developed which successfully predicted the 1997 measles epidemic of New Zealand. Roberts and Tobias [4] extended this model by including age structure to examine optimal timing of vaccine doses and compare the effects of various potential vaccine strategies. In 2005, Pang et al. [5] proposed a SEIR model with vaccination to investigate oscillations of measles and fit the case data of U.S. from 1951 to 1962 to design a control strategy. Considering some developing countries with high measles burden and limited resources, Verguet et al. [6] proposed a model to estimate optimal scheduling of SIAs without second routine vaccination. In contrast, for some countries where measles has been eliminated, Choi et al. [7] estimated susceptibility to measles by age and Mckee et al. [8] investigated the optimal age target for two routine dose without changing vaccination coverage in order to maintain elimination. Furthermore, Bai and Liu [9] fitted the model to the case data in China and Huang et al. [10] investigated the effect of various interventions on measles infection. orrington et al. [11] attempted to quantify the loss of QALY (quality-adjusted life years) due to measles at a population level in England and provided important parameters for future control interventions. In addition, some researchers found critical community size and herd immunity influence the extinction of measles in communities [12][13][14][15][16], and obtained that human mobility has an impact on the periodicity of measles outbreaks [17].
ere are studies showing the evidence of waning immunity in vaccinated individuals [18][19][20], indicating individuals with sufficiently low antibody levels are at risk of mild or subclinical measles infections [21,22]. Although transmission of virus between such individuals with asymptomatic infection has not been demonstrated, serological studies have indicated that this may occur [23,24]. In China, measles cases can still be observed in recent years even if we have initiated enhanced immunization and got high coverage rate of vaccination. Possible reasons may include mature women with low antibody can either be infected by their infected children without showing any symptom and in turn infect others or have babies who have low antibody and can be infected. Hence, the exploration of the influence of asymptomatic infection on the transmission of measles and identification of the key factors or processes that significantly influence the measles outbreak in mainland China remain unclear and fall within the scope of this study.
e main purpose of this paper is to formulate a dynamic model incorporating periodic transmission and asymptomatic infection with waning immunity to identify the key parameters which influence seasonal fluctuation of measles and discuss optimal control measures to minimize the number of infected individuals and costs. In the next section, a periodic model with vaccination and asymptomatic infection due to waning immunity is introduced and then the existence and stability of periodic solution are discussed. By fitting monthly measles data in mainland China from 2013 to 2016, we estimate some unknown parameters and calculate the basic reproduction number. We then numerically evaluate effect of asymptomatic infection and other important parameters on the transmission dynamics of measles. In Section 3, we apply optimal control theory to design optimal control strategies such that the number of infected individuals and costs are minimum. Finally, a discussion is presented in Section 4.

Model Formulation.
e underlying structure of the model comprises of eight compartments. We assume that susceptible individuals (S) receive the vaccine and enter into the compartment V 1 at rate of p. en, the vaccinated individuals (V 1 ) progress to the class (V 2 ) with waning of immunity at rate of ω. Susceptible individuals (S) become infected and enter into the latent undetectable class (E) at rate of β(t) (I + ξI v ), where β(t) is the baseline periodic transmission rate and ξ (0 ≤ ξ ≤ 1) quantifies the relative transmissibility of asymptomatic infections. en, the exposed individuals (E) progress to infectious class (I) at rate of σ and recover with rate of c. Note that most individuals with waning immunity may not produce symptoms, once infected again; hence, it is negligible that latent individuals (E) move to asymptomatic infections class (I v ). e individuals in V 1 class are those who are effectively vaccinated with high antibody titer, so we assume they cannot be infected as mentioned in [19,20]. But because the antibody titer is low, individuals located in V 2 may be infected and become exposed individuals E v at rate of ηβ(t) (I + ξI v ) with the modification factor η (0 ≤ η ≤ 1), denoting these individuals having low probability to be infected.
en, these latent individuals may move to symptomatic or asymptomatic infectious class (I) and (I v ) at rate of θσ or (1− θ)σ, where θ (0 ≤ θ ≤ 1) is the proportion of latent individuals becoming symptomatic infections. e recovery rate from I v to R is c v . In general, asymptomatic infections recover faster than symptomatic infection, which means c v > c. A flow diagram of the model is described in Figure 1 and the definitions of variables and parameters are given in Table 1. e specific dynamic model is as follows: where β(t) is nonnegative and periodic function with period T (T > 0) and Λ and μ are recruitment rate and natural death rate per unit time. All parameters are assumed to positive constants.  Proof. Let Hence, lim t⟶ +∞ N(t) � (Λ/μ), which implies for any  It is easy to obtain that system (1) has unique disease-free equilibrium E 0 � (S, V 1 , V 2 , 0, 0, 0, 0, 0), where S � (Λ/p + μ), V 1 � (Λp/(p + μ) (ω + μ)), V 2 � (Λpω/(p + μ)(ω + μ)μ). In the following, the stability of disease-free equilibrium will be discussed. First, we define the reproduction number of system (1) on the basis of [25]. Denote Assume Y(t, s), t ≥ s is the evolution operator of following linear system: at is to say, Let Φ − V (t) be the fundamental solution matrix of linear system (4). Clearly, Φ − V (t) equals to Y(t, 0), t ≥ 0.
Define a linear operator L : where C T is the ordered Banach space of all T-periodic function from R to R 4 with maximum norm ‖·‖ and positive  Figure 1: e flow diagram of the measles transmission.
cone C + T :� ϕ ∈ C T : ϕ(t) ≥ 0, ∀t ∈ R . Suppose that ϕ ∈ C T is the initial distribution of infectious individuals in the periodic environment; then, F(s)ϕ(s) is the distribution of new infections produced by the infected individuals who were introduced at time s. Given t (t ≥ s), then Y(t, s)F(s)ϕ(s) gives the distribution of those infected individuals who were newly infected at time s and remain in the infected compartments. Define the spectral radius of L as the basic reproduction number of system (1), i.e., It can be verified that system (1) satisfies conditions (A1)-(A7) of periodic system [25]; then, the following theorem is established. [25]. en, the following statements are valid: Proof. According to Lemma 1, it is easy to get E 0 is local stable if R 0 < 1. e following will show the global attractive.
By system (1), we obtain dS dt ≤ Λ − pS − μS, By standard comparison theorem [26], we get for any Consider auxiliary system We can write (9) as where X � (Ê, E v ,Î, I v ), and is a solution of (9). Choose t > t 1 and a small value α > 0 such that X(t) ≤ αv(0); then, we can get By standard comparison theorem [26], we get By the last four equations of system (1), considering their limit system, i.e., dS dt And there exists at least one positive T-periodic solution for system (1).
First, it is easy to see that X, X 0 are positively invariant and zX 0 is relatively closed in X. It follows from the Proposition 1 that the solution of system (1) is uniformly and ultimately bounded. us, the semiflow P is point dissipative on R 8 + and P : R 8 + ⟶ R 8 + is compact. By eorem 3.4.8 of [28], P admits a global attractor A. Define It indicates that ( e following shows W s (E 0 ) ∩ X 0 � ∅. First, by the continuity of solutions with respect to initial values, for any In the following, we first claim Hence, from system (1), we obtain Consider auxiliary system For convenience, it can be rewritten as and By Lemma 2.1 of [27], there exists a positive T-periodic function v(t) such that e μ 2 t v(t) is the solution of (18), where By standard comparison theorem [26], we get is is contradiction to Above all, by eorem 3.11 [29], we get P is uniformly persistent with respect to (X 0 , zX 0 ).
Next, prove that for R 0 > 1, there is a positive periodic solution of system (1). It follows from eorem 1.3.6 of [29] that the Poincaré map P has a fixed point X * � (S * (0), Suppose not assuming S * (0) � 0, then the solution S(t) through initial point S * (0) satisfies equation

Mathematical Problems in Engineering
We get (μ+p+β(s)(I(s)+ξI v (s)))ds (S * (0) + t 0 Λe s 0 (μ+p+β(τ)(I(τ)+ξI v (τ)))dτ ds). Since X * is a fixed point of P, we can obtain S(T) � Λe R(t)) be a periodic solution through X * . the uniform persistence of the solution with respect to (X 0 , zX 0 ), i.e., lim m⟶+∞ inf(S(t), In this section, we will numerically analyze model (1) and identify the effects of some parameters we are interested in on the outbreak of measles. According to reference [9], we obtain the recruitment rate Λ � 1589357 per month. Assuming the current average life expectancy of Chinese people is 73 years, then we calculate the natural death rate μ � 0.001142 per month. From [11], we derive that the exposed and infectious periods are about 13.8 and 10 days, which indicates σ � 2.17 month − 1 , c � 3 month − 1 . We let the number of initial infectious individuals be the sum of the numbers of reported measles cases for last 10 days (one infectious period) in December of 2012, that is, I(0) � 300. By fitting model (1) to the monthly reported measles data from 2013 to 2016, we estimate some unknown parameters, which are listed in Table 1. Figure 2 shows a high goodness of fit (R 2 � 0.8378) which demonstrates that the model can capture main trends of measles outbreak events between 2013 and 2016. Note that the recovery rate c v is estimated as 6 month − 1 , that is, 5 days, which means symptomatic patients will take twice the time to recover than asymptomatic infectious individuals. Furthermore, the estimated value of θ is small, demonstrating that most individuals in group E v move into the class I v . In other words, individuals in V 2 class, once infected, mainly become asymptomatic. Using the estimated parameter values, we calculate the basic reproduction number R 0 as 1.0077, which indicates measles will persist.
It is worth noting that two routine immunizations are objected to children under 2 years old in mainland China; hence, the vaccinated individuals in our model actually belong to this age group. Let S 1 be the number of susceptible children younger than 2 years old and p c be the effective coverage rate corresponding these children, while parameter p in the model (1) represents the coverage rate of vaccination for general population; then, we actually have pS � p c S 1 , and in the following simulations, we mainly use p c to investigate the impact of vaccination on disease outbreak and make discussion.
To determine the influence of asymptomatic infection and other important parameters on the measles infection, we investigate how the basic reproduction number and the number of symptomatic infections vary with parameters. It follows from Figures 3(a) and 4(a)-4(c) that reducing baseline transmission rate (a) or increasing coverage rate (p c ) and recovery rate (c) significantly decreases the basic reproduction number and consequently the number of symptomatic infections declines.
is indicates limited contact, extensive vaccination, and timely treatment can effectively mitigate the outbreak and transmission of measles.
en, considering parameters θ, ξ, η relating to asymptomatic infection, we find that decreasing the values of θ, ξ, η has slight impact on reducing the new infections and controlling the transmission of measles except for minutely decreasing the peak size of symptomatic infections, as shown in Figures 3(b)

and 4(d)-4(f ). It indicates that asymptomatic measles infection in population does not cause much new infections.
Further, in order to access the dependence of R 0 on parameter variations, we perform sensitivity analysis using the Latin hypercube sampling (LHS) method and calculate partial rank correlation coefficients (PRCCs) [30,31] for various input parameters against the basic reproduction number. We choose uniform distribution in all parameters as listed in Table 2 since there is limited information on the distribution of each parameter. Figure 5 shows that, in contrast to parameters relative to asymptomatic infection, coverage rate p c , the baseline transmission rate a, and recovery rate c are the first three parameters that mostly influence R 0 .
ese results indicate that reducing transmission rate via limiting contacts (e.g., personal protective and social distancing), treating infected individuals timely, and enhancing vaccination greatly lead to decline in new infections.

Optimal Control Strategies
In this section, we extend the basic model (1) by including some feasible controls and formulate the optimal control system, which is defined by the following model: The changing rate of parameters (k) e control function u 1 (t) represents reducing contact between susceptible and infective individuals on the basis of personal protection (e.g., wearing masks and enhancing awareness about measles) or social distancing (e.g., Fengxiao) measures. e control functions u 2 (t) and u 3 (t) represent the enhanced treatment and vaccination, respectively.
Our goal is to minimize the number of infectious individuals over a time interval [0, T] at a minimal cost. en, the objective function J is defined as where W 1 , W 2 , and W 3 are the weight constants for the control functions. e optimal control strategies can be obtained by finding an optimal pair (U * , Y * ) such that where

Optimality
System. e existence of optimal controls follows from Corollary 4.1 of [32] since the state system satisfies the Lipschitz property with respect to the state variables and is a linear function of U, while the integrand of J is a convex function for u 1 (t), u 2 (t) and u 3 (t).

Theorem 3.
ere exist optimal controls U * (t) and corresponding solutions Y * (t) that minimize J(U(t)) over Ω. In order for the above statement to be true, it is necessary that there exist continuous functions λ i (t) such that with the transversality conditions Furthermore, the optimal control is given as follows:

Mathematical Problems in Engineering
Proof.
e necessary conditions satisfied by optimal controls can be derived from Pontryagin's maximum principles. Let Hamilton function be en, the adjoint equations with transversality satisfy e Hamilton function H is minimized with respect to the controls by differentiating H with respect to u 1 , u 2 , u 3 on the set Ω, respectively, and then gives the following optimal conditions: zH Solving for u * 1 , u * 2 , u * 3 yields By using the bounds 0 ≤ u 1 (t) ≤ b 1 , 0 ≤ u 2 ≤ b 2 , 0 ≤ u 3 ≤ b 3 , we have the properties as eorem 3. □ 3.2. Numerical Simulations. As mentioned in Section 2.3, based on relationship between p c and p, we initially define the control function u 3c (t) as the enhanced vaccination rate for children under 2 years old to replace u 3 (t) for general population and numerically explore the effects of optimal control measures on measles outbreak. e optimality system is solved by using the forward-backward sweep method. With an initial guess for the control variables, the state system and adjoint system are solved forward and backward, separately. e control variables are updated by using optimality condition, and this process is repeated until the system is converged. In the numerical simulations, we use parameter values given in Table 1 To study the effect of optimal control strategies on the transmission of measles infection, we initially plot the time series of the number of symptomatic infected individuals and corresponding optimal functions, as shown in Figure 6. It follows from Figure 6(a) that the optimal control strategies significantly mitigate the outbreak of disease and decrease the number of symptomatic infections. Figures 6(b)-6(d) give the optimal control measures u * 1 (t), u * 2 (t) and u * 3c (t), which monotonically decrease after a period time of implementation with high intensity.
is indicates we should enhance vaccination as soon as possible when measles begin outbreak, while measures such as personal protection and treatment should be implemented as much as possible throughout the whole epidemic.
To further specify the extent of each control measure in reducing the outbreak of measles, we compare the number of symptomatic individuals with only one control measure and without any control measure, as shown in Figure 7. e results illustrate each optimal control measure has great impact on decreasing the measles infections. In particular, the optimal enhanced vaccination for children u * 3c (t), which has limited impact on lowering the first outbreak, results in a low probability of occurrence of repeated outbreaks via decreasing the susceptibility to measles. e timing of implementation of control measures is always a critical issue. erefore, we study the optimal control functions and their influences on the reduction of infected individual under three different starting times: (1) the time when the number of infective individuals is at the lowest level, (2) the time when the number of infections increases but is relatively small, and (3) the time when measles spreads fast with severe level. Here, we can choose the 10th, 13th, and 15th month from January of 2013. Figure 8(a) demonstrates that the relatively late implementation of control measures leads to a great increase in the number of infectious individuals. Comparing optimal control functions u * 1 (t), u * 2 (t), and u * 3c (t) with different starting times of control measure as shown in Figures 8(b)-8(d), we obtained that as the starting time is delayed, the optimal control strategy u * 1 (t) has no obvious changes, while the optimal control function u * 2 (t) must be implemented with maximum level for a long period. Different from u * 1 (t) and u * 2 (t), Figure 8(d) shows that the optimal enhanced vaccination for children u * 3c (t) should be implemented with more maximum level (shown in green curves) and then be implemented with relatively low maximum level (shown in black curves) as the initiation time is delayed.
ese results indicate the earlier the control measures are implemented, the lower the outbreak is and the less (period or intensity) the optimal controls are needed.

Discussion
In mainland China, measles reached its lowest level in 2012 and then resurged again with seasonal transmission since 2013. ese outbreaks have led to an interest in understanding the underlying biological mechanisms and public concern of effective control. To investigate the phenomenon of seasonal outbreak in China and identify the crucial parameters played in the transmission dynamics of measles, we propose a mathematical model with periodic transmission rate. In this system, we extend the existing dynamic model [9,10] by considering asymptomatic infection due to vaccinated populations with low antibodies level being infected. We investigate the threshold dynamics by defining the basic reproduction number. Further, based on the baseline model (1), we propose three types of feasible control measures and formulate a optimal control system to discuss how to design the optimal measures to control epidemic outbreak. We theoretically analyze dynamic behavior of the periodic system by defining the threshold R 0 . We prove the disease-free equilibrium E 0 is globally asymptotically stable for R 0 < 1, meaning the measles can be eliminated. If R 0 > 1, the diseasefree equilibrium E 0 is unstable and system (1) is uniformly persistent, which indicates the disease cannot be effectively eliminated. en, by fitting our proposed model to the reported cases in mainland China from year 2013 to 2016, we estimate some unknown parameters and calculate the basic reproduction number R 0 with immunization as 1.0077. Compared with original basic reproduction number of measles (R 0 � 12.8) [4], this result implies that current vaccination policies have dramatically reduced the spread of disease. Numerical simulations indicate that asymptomatic infection has slight impact on the repeated outbreaks of measles and the crucial factors that significantly affect new infections are vaccination rate, transmission rate, and recovery rate. Further, we propose three types of control measures (limiting contacts between susceptible and infective individuals and enhancing treatment and vaccination) and determine the optimal control strategies to minimize the number of infected individuals at minimum costs. It has been illustrated numerically that the optimal control measures have a very attractive effect on reducing the number of infected individuals and the interventions should be implemented as soon as possible. We further observed that enhancing vaccination is especially effective in lowering the initial outbreak and making disease recurrence less likely.
Note that due to lack of data on measles cases as for different ages, we formulate a measles transmission model without considering age structure. However, the vaccination of measles mostly is implemented among children, and we shall formulate the dynamic model with age structure, conduct cost-effectiveness assessment, and explore the optimal vaccination age in near future.

Data Availability
e data used to support the findings of this study are included within the article.

Conflicts of Interest
e authors declare that they have no conflicts of interest.