{'fraction_non_alphanumeric': 0.03629876448559549, 'fraction_numerical': 0.023029516704508597, 'mean_word_length': 4.5945232596502805, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 12, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Identifying epidemic spreading dynamics of COVID-19 by pseudocoevolutionary simulated annealing optimizers

At the end of 2019, a new coronavirus (COVID-19) epidemic has triggered global public health concern. Here, a model integrating the daily intercity migration network, which constructed from real-world migration records and the Susceptible–Exposed–Infected–Removed model, is utilized to predict the epidemic spreading of the COVID-19 in more than 300 cities in China. However, the model has more than 1800 unknown parameters, which is a challenging task to estimate all unknown parameters from historical data within a reasonable computation time. In this article, we proposed a pseudocoevolutionary simulated annealing (SA) algorithm for identifying these unknown parameters. The large volume of unknown parameters of this model is optimized through three procedures co-adapted SA-based optimization processes, respectively. Our results confirm that the proposed method is both efficient and robust. Then, we use the identified model to predict the trends of the epidemic spreading of the COVID-19 in these cities. We find that the number of infections in most cities in China has reached their peak from February 29, 2020, to March 15, 2020. For most cities outside Hubei province, the total number of infected individuals would be less than 100, while for most cities in Hubei province (exclude Wuhan), the total number of infected individuals would be less than 3000.


Introduction
Infectious diseases have been raging in the world many times in history. For example, the Black Death (also known as the Pestilence) in the fourteenth century lasted for 30 years in Europe, with more than 25 million deaths, accounting for about 1/3 of the European population at that time [2]. In 1918, Spanish influenza, which initially outbreaks within the U.S. military, eventually swept the world, infected nearly 600 million people, and caused about 40-50 million death [3]. In 2003, SARS, with a fatality rate of 11%, spread from Guangdong Province to the whole country, bringing huge losses to China's national economy [4]. In 2009, the H1N1 flu outbreak spread to 214 countries and regions, causing 1, 220 deaths in a few months [5]. In 2014, the outbreak of the Ebola epidemic resulted in 28,637 infections and 11,315 deaths. At the end of 2019, a highly contagious disease, which is caused by infection of the SARS-CoV-2 virus and named the 2019 Coronavirus Disease (COVID- 19), broke out and caused millions of infections [6,7]. The spreading trends of the COVID-19, including when the peaks would occur, how many people would eventually be infected, the final infection rate of the population of each city, and which cities would run the risk of being out of control, are the core questions that need to be answered.
Scholars from various disciplines have participated in the research on epidemic transmission and control. The epidemic spreading model can be traced back to the This paper is an updated version of a preprint uploaded on Feb 20, 2020, to medRxiv.org (doi: 2020.02. 18 [8], while the most classic epidemic spreading model is the Susceptible-Infected-Removed (SIR) model proposed by Kermack and McKendrich in 1927 [9]. Based on the SIR model, scholars found that there exists an epidemic threshold depended on the infection and recovery rate. If the infection rate is greater than a threshold, epidemics will spread on a large scale in the population. In 1932, Kermack and McKendrich established the Susceptible-Infected-Susceptible (SIS) model [10], which is similar to the SIR model, except that infected individuals would return to a susceptible state instead of an immune state after recovery. In 1992, for infectious diseases with a limited immune period, J. Mena-Lorca and et al. improved and proposed a more complex Susceptible-Infected-Removed-Susceptible (SIRS) model [11]. Inspired by these pioneering studies, a lot of efforts have been devoted to investigating epidemic spreading under various circumstances [12,13]. One of the basic problems in theoretical epidemiology is to study the epidemic threshold of classic epidemic spreading models, e.g., SIR, SIS, and SEIR, in various networks, such as scale-free network and its extension [12,13], complex heterogeneous networks [14], the real-world Oregon graph [15,16], adaptive networks [17], and the complete worldwide air travel network [18]. Previous literature has unveiled that the epidemic threshold is highly related to the spectral radius of the adjacent matrix of the network [19,20]. Other researches focused on developing epidemic spreading models in different scenarios and the temporal evolution [21,22]. However, to our best knowledge, few of them focus on developing models to describe and predict the dynamic of real epidemic spreading cases. The migration of individuals, especially intercity migration, plays a core role in the spread of SARS-CoV-2 [7,18]. Wuhan is a metropolis with a population of more than 11 million and is also one of the transportation hubs in China. During the Spring Festival in 2019, 5 million people set off from Wuhan to other cities in China. Large-scale migration greatly enhanced the spread of COVID-19 from outbreak areas to other cities in China. Therefore, the intercity migration data is also an important indicator for describing and predicting the spread of the virus. Here, daily intercity migration data for 367 cities in China are collected and utilized to construct intercity migration networks. Further, a model established by combining complex network theory and the classic SEIR model can be used to describe how COVID-19 spreads from Wuhan to other cities in China. This dynamic model has more than 1800 unknown parameters to be determined from the historical data. The inference of model parameter values from rare historical time-course data can be reformulated as an optimization problem and is still one of the most challenging tasks [23,24]. Evolutionary algorithms have outstanding performance in solving nonlinear optimization problem [25,26]. Hence, it is worth developing evolutionary algorithms, which should be robust against noise, efficient in computation, and flexible enough to meet different constraints for estimating these 1800 unknown parameters. In this article, a novel pseudocoevolutionary simulated annealing algorithm is proposed to solve this problem [27]. Results show that the proposed algorithm successfully identified optimal parameter sets of this epidemic spreading model. Also, the model can fit the number of infected, recovered individuals, and the death toll of each city with a minor error.
Based on the model, we find that migration control was extremely effective in controlling the spread of the epidemic. If the government continues strict migration control, the infections numbers of most cities in China would peak between mid-February to early March 2020. The peak number of infections in most cities is smaller than 100, while the proportion of infected individuals in each cities population is smaller than 0.01%. However, if the epidemic spreading is out of control, it would infect about 1% of the population in Hubei province, while infecting about 0.3% population outside Hubei province. The peak number of infections in most cities would come at the end of April 2020. Evidence shows that China has controlled the spreading of COVID-19.
The main contributions of this study are as follows: • First, we integrate daily intercity migration data and traditional SEIR model to develop an extended SEIR model; • A novel pseudocoevolutionary simulated annealing algorithm is proposed. Additionally, we compared the estimation result with simulated annealing, particle swarm optimization, and pattern search algorithms. Results show that the proposed algorithm provides the best results; • The pandemic situation of China has been investigated.
Results show that this technique can accurately reflect the spread of COVID-19. Study shows that migration control is extremely effective in controlling the spread of the epidemic.
The rest of this paper is organized as follows. Section 2 reviews related works. The description of data is introduced in Sect. refsec:dataspsdescription, including official released confirmed cases, recovered cases, death toll, and intercity migration data. In Sect. 4, the SEIR-migration model and pseudocoevolutionary simulated annealing algorithm are introduced. Then, the experimental design and results are shown in Sect. 5. Finally, conclusion and future work are presented in Sect. 6.

Literature review
Complex network theory is a powerful tool for researchers to study epidemic spreading. With the development of complex network theory, a large amount of work has investigated the effect of the structure of complex networks (such as degree of relevance, clustering coefficients, community structure, hierarchical structure, and edge weights) on the propagation properties of infectious diseases (such as spreading rate, scale of propagation, and epidemic threshold) [13,28]. In 2001, Pastor-Satorras et al. studied the propagation model on scale-free networks and proved that scale-free networks are weak against infectious diseases and can maintain spreading at any small infection rate [12,13]. In the same year, May and Lloyd investigate the effect of network scale on spreading behavior of scalefree networks and pointed out that there are positive propagation thresholds for limited scale-free networks [29].
In 2002, Newman investigated the SIR propagation model in a scale-free network and proved that there is an epidemic threshold when the cutoff value of the degree of nodes in the network is relatively small [30].  [33]. Epidemic threshold b c plays an important role in epidemic spreading. For a large-scale system, when the infection rate b [ b c , the proportion of infected individuals will reach a limited proportion. Otherwise, if b\b c , the proportion of infected people will reduce to almost zero. Therefore, to control the outbreak of epidemics, reducing the infection rate is one of the effective ways. Research shows that frequent absorption, wearing a mask, and disconnecting from the infectious individuals will reduce the probability of infection and then effectively control or slow the outbreak of an epidemic [34]. Additionally, studies have shown that quarantining, closing schools, and restricting individuals from attending public events can make people's contact networks sparse and reduce infection rate [35]. In 2011, Jin et al. developed an epidemiological model of influenza A, demonstrating that an immunization strategy targeting specific populations with given connectivity can greatly reduce epidemic spreading [36]. Two years later, Guo et al. introduced a continuoustime adaptive susceptible-infectious-susceptible (ASIS) model, proving that the adaption of the topology can inhibit infection [37]. In the same year, Peng et al. investigate several epidemiological models, including susceptibility, infection, and incomplete vaccination segment models, on the Watts-Strogatz small world, Barabasi-Albert scalefree, and random scale-free networks, for analyzing the epidemic threshold and infection rate [38]. More information can be found in the review papers [39]. However, to our understanding, most of the literature about epidemiological models analyzing the propagation in a network is based on analytical methods and large-scale simulations, without the support of real-world data.
A common approach for explaining and analyzing realworld phenomena is to establish epidemiological models based on real-world observation data. These epidemiological models are always nonlinear dynamical models with high parameter dimension, which is often presented as a set of ordinary differential equations (ODEs) or discrete-time equations containing a large volume of unknown parameters [40]. The identification of a large volume of unknown parameters from historical observation data is critical for judging the performance of an epidemiological model. The method for identifying unknown parameters can be classified as ''reverse engineering techniques,'' which usually formulate the problem of parameter identification into a nonlinear optimization problem that minimizes an objective function representing the fitness of the model with respect to the observation data [23,24]. The identification results are highly dependent on the optimization algorithm. Due to the simplicity and ease of use, evolutionary algorithms are widely utilized to identify unknown parameters of nonlinear dynamic models [41][42][43][44]. However, for nonlinear dynamical systems with high parameter dimension, the objective function is always complicated and has tremendous local minima. Parameter estimation algorithms face a high possibility of converging at local optima but not global minima [45]. Additionally, one parameter identification trial always required tens of hours computation time [46,47]. Therefore, the parameter estimation problem is still a challenging task and even a bottleneck for nonlinear dynamical models with high parameter dimension. Evolutionary algorithms have been extensively used in nonlinear optimization problems and shown that can provide satisfying results [25][26][27]48]. In this article, a new pseudocoevolutionary algorithm is proposed to solve this hard engineering problem, and the more detail information about hard engineering can be found in [49][50][51][52].

Official data of COVID-19 cases
Testing is the only way to know whether a susceptible individual is infected or not. At present, there exist two kinds of tests (techniques) for testing COVID-19. One kind of tests (techniques) checks the presence of the COVID-19 virus, aiming to establish whether an individual is currently infected. The other kind of tests examines the presence of antibodies, which can figure out whether an individual has been infected in the past, even this individual has recovered and not carried COVID-19 virus now. A summary of the current state of testing technologies associated with their implementations can be found in [53]. Now, the most common way to perform a COVID-19 test is adopted by detecting the viral RNA through polymerase chain reaction (PCR). In this work, the official number of infected cases only contains individual who has a positive COVID-19 testing result. In China, Wuhan is received the first confirmed case of COVID-19 infection on December 8, 2019 [6]. Most other cities in China released data of COVID-19 infections around January 20, 2020. The data of COVID-19 infections, recovery, and death toll used in the study were derived from official data released by the National Health Commission of China. Hubei Province was the epicenter of the epidemic in China. Most of the infections occurred in Hubei Province (as shown in Fig.1), while the number of other provinces is relatively small. In this study, one of our aims is to develop an epidemiological model that accurately describes how the number of infections, recovery, and death toll change over time in various cities.

Intercity travel data
COVID-19 mainly spread through human-to-human transmission. In this case, the intercity migration of infected and exposed individuals has become the main driving force for COVID-19 to spread from one city to another. Chinese New Year (mid-January to early February 2020) is the most important holiday for Chinese people. In 2020, the period around during the Spring Festival holiday in 2019 is approximately from mid-January to early February 2020. Wuhan, as one of the most important transportation hubs in China and the world, is one of the cities with the largest flow of entry and exit around the Spring Festival. China's Ministry of Transport estimates that Wuhan has about 5 million trips, while China as a whole has about 3 billion trips during the Spring Festival holiday. We have collected daily intercity travel data for 08Dec19 21Dec19 03Jan20 16Jan20 29Jan20 11Feb20 24Feb20 367 cities in China. The data provide the intensity of population migration and also indicate the strength of the population in and out of various cities. Based on these data, we can develop the migration networks (shown in Fig. 2). After the outbreak of COVID-19, the Chinese government has rapidly restricted the intercity migration since January 23, so the strength of intercity migration has dramatically reduced since January 23. Figure 3 shows the total inflow/ outflow of travelers of 6 metropolis in Chinese. Note that after migration control, the strength of the intercity migration of Wuhan has almost reduced to zero. The control measure effectively reduced the speed of virus transmission and ultimately successfully controlled the further spread of the virus. Base on this data, we can construct the migration matrix, which is given as where K ¼ 367 is the number of the cities, and m ij ðtÞ is the migrant volume from city i to city j at time t. Migration matrix M thus effectively describes the network of cities with human movement constituting the links of the network. Figure 3 plots the daily total inflow and outflow migration strengths of Wuhan, showing the abrupt decrease in migrant strength after the city shut down all inflow and outflow traffic from February 01, 2020.
4 The SEIR-migration model and pseudocoevolutionary simulated annealing algorithm First, we will give a brief description of human contact networks with multiple sub-networks representing a city or administrative regions for epidemic spreading propagation. A real human contact network consists of multiple subnetworks, just as a country consisting of many cities, towns, and villages. Here, we consider a human contact network G ¼ ðV; EÞ contains K sub-networks fG 1 ; G 2 ; . . .; G M g. V stands for the set of nodes, and E is the set of edges. Here, each node represents an individual. If two individuals/nodes v i and v j have contacts, there will be a link e i;j between them, otherwise, no connection (shown in Fig. 4a). Note that nodes in the same sub-network have plenty and strong connections with network neighbors, which results in a highly clustered sub-network. However, nodes belonging to different sub-networks have less and weak connections. In this work, a sub-network can be treated as a city, while nodes in a sub-network stand for citizens. The K sub-network (city) forms a huge contact network of a country (shown in Fig. 4b).
In the classic SIR model, each individual can be in three different states: infected (I), susceptible (S), and recovered (R). In an epidemic spreading case, infected individuals (I) can infect susceptible individuals (S) through human contact network, while an infected individual can be cured and turn into recovery state (R). Once an infectious disease starts to spread in a certain sub-network, due to the denseness of the sub-network nodes and the short distance between neighboring nodes, the epidemic will outbreak within the sub-network in a short time (G 1 in Fig. 4b).
With the epidemic spread, plenty of susceptible nodes in the network transit into infected nodes. Some infected nodes have connections with nodes in other sub-network. Then, the virus can spread from one sub-network to another sub-network through the nodes connecting two sub-networks and eventually spread to the entire network (G 2 in Fig. 4b).

SEIR-Migration model
Studies reveal that the median incubation period of COVID-19 to be 5.6 days (95% CI 4.8-6.3) [54]. Exposed individuals can also infect other individuals during the incubation period. Each node in human contact networks may assume one of four possible states in the epidemic spreading process, namely susceptible (S), exposed (E), infected (I), and recovered/removed (R). For sub-network (city) j, the number of nodes in the four states S j ðtÞ, E j ðtÞ, I j ðtÞ, and R j ðtÞ, at time t. 01Jan20 10Jan20 19Jan20 28Jan20 06Feb20 15Feb20 25Feb20  Here, m ij ðtÞ represents the volume of individuals moving from city i to city j at time t. Here, we assume the population of city j is P j . Then, the number of infected individuals moving from city i to city j is DI in ij ðtÞ ¼ Also, the volume of individual migrating from city j is P N i¼1 m ji . Then, the number of infected individuals moving out of city j is DI out j ðtÞ ¼ Thus, the dynamic change of the number of infected cases in city j at time t is given by where DI j ðtÞ ¼ I j ðt þ 1Þ À I j ðtÞ. Moreover, if city i has a population of P i and the eventual percentage of infection is where N s i is the eventual number of infections. Similarly, the dynamic changes of infected, susceptible, recovered individuals, and the population of a city can be obtained. Then, we have the modified SEIR model with consideration of human migration dynamics as follows: where DE j ðtÞ ¼ E j ðt þ 1Þ À E j ðtÞ, DS j ðtÞ ¼ S j ðt þ 1Þ ÀS j ðtÞ, DR j ðtÞ ¼ R j ðt þ 1Þ À R j ðtÞ, DN s j ðtÞ ¼ N s j ðt þ 1Þ ÀN s j ðtÞ, and DP j ðtÞ ¼ P j ðt þ 1Þ À P j ðtÞ. The physical meaning of each parameter of model (6) is presented in Table 1, while a detailed description of the model is given in [1]. Note that we assume the recovered individuals are assumed to stay in the city j.

Parameter identification problem
Model (6) has a large volume of unknown parameter. The parameter estimation problem can be transformed into a nonlinear optimization problem (NLP). The purpose of optimization is to find a set of suitable parameters to make the estimated growth trajectory that matches historical data. Here, we define I j;0 ¼ I j ðt 0 Þ and E j;0 ¼ E j ðt 0 Þ, which are the initial number of infected and exposed individuals in city j, respectively. Note hat Wuhan is the epicenter of the COVID-19 pandemic in China, with Hubei province being the region immediately surrounding it. Therefore, it is reasonable to assume that Wuhan has initially infected individuals.
Then, I WH ðt 0 Þ; E WH ðt 0 Þ 6 ¼ 0, and I j ðt 0 Þ; E j ðt 0 Þ ¼ 0 for all other cities, where I WH ðt 0 Þ and E WH ðt 0 Þ represent the initial number of infected and exposed individuals in Wuhan, respectively. For city j, there exist a set of unknown parameters, i.e., Then, the unknown parameter set is H ¼ fI WH;0 ; E WH;0 ; k I ; h 1 ; h 2 ; . . .; h K g. Totally, there exist 5K þ 3 unknown parameters, where K is the number of cities. Thus, an enormous effort of computation is required to estimate the suitable parameters.
Let X j ðtÞ be the extended state vector, i.e., X j ðtÞ ¼ ½S T j ðtÞ E T j ðtÞ I T j ðtÞ R T j ðtÞ P T j ðtÞ ðN s j ðtÞÞ T T , then, we define: Model (6) can be reformulated as where f(x) is the right side of (6), and H is the set of unknown parameters. Note that DX j ðtÞ ¼ Xðt i þ 1Þ À Xðt i Þ, then, Eq. (9) can be reformulated as: Finally, the parameter estimation problem can be formulated as the following constrained nonlinear optimization problem: whereÎðt i H j Þ represents the estimated number of infected individuals at time t i with parameter set H and initial condition fI WH;0 ; E WH;0 g. w ij stands for the weighted coefficient. The unknown parameter set is bounded between H L and H U . In this work, an inverse approach is taken to find the unknown parameters and states by solving (11). The rate at which the infected individuals infect the susceptible individuals in city j a j :

Proposed pseudocoevolutionary simulated annealing algorithm
The rate at which the exposed individuals infect the susceptible individuals in city j j j : The rate at which exposed individuals become infected in city j c j : The recovery rate in city j k I : The possibility of an infected individual moving from one city to another d j : The eventual percentage of infections in city j I j;0 : The initial number of infected individuals in city j E j;0 : The initial number of individuals in city j Algorithm 1 Sub-algorithm for calculating time series I j (t i ) Input: Parameter set and initial number of infected and exposed individuals of each city Θ = {I W H,0 , E W H,0 , k I , θ 1 , θ 2 , · · · , θ K }, where θ j = {α j , β j , γ j , κ j , δ j } and j = 1, 2, · · · K; Output: A time series {I(t 0 ), I(t 1 ), · · · , I(t N )}; Initialization  Calculate X(t i + 1) from model (6), namely, Since I j (t i ) ∈ X j (t i ) and X j (t i ) ∈ X(t i ), we can derive I j (t i ) from X(t i ). 4: end for 5: return I(t 0 ), I(t 1 ), , I(t N ) Algorithm 2 Sub-algorithm for estimating the optimal parameter setΘ * Input: The initial parameters set of all the K cities Θ 0 and the set of index (Φ = {a 1 , a 2 , · · · , a n }, 1 ≤ a j ≤ K) of the city, of which the parameters need to be optimized. Output: Optimal parameter setΘ * ; Initialisation : 1: Initialize temperature T and set Θ current = Θ 0 Compute the value of objective function from Algorithm-1 and objective function (11) with parameter set Θ current Initialize temperature T and random starting point where k rand is a randomly generated real number between 0 and 1. LOOP Process 2: for i iter = 0 to i max do 3: i iter = i iter + 1, temp iter = 0, Θ previous = Θ current cost previous = cost current 4: while temp iter ≤ n rep do 5: temp iter = temp iter + 1. Keep θ j = θ j,0 .

end if
Compute the value of objective function (11) and derive δ δ = cost current − cost previous . 11: if δ < 0 then 12: Accept new parameter set. Accept new parameter set with probability exp(−δ/T ). 15: end if 16: end while 17: T = α * T, (0 < α < 1). 18: end for 19: returnΘ * Note that the cumulative number of infections varies widely by cities. For example, Wuhan, the city with the largest number of infected people in China, is infected by more than 50,000 people. However, in some small cities, only a dozen people are infected, or even no one has been infected yet. Therefore, the weighted coefficient of the objective function 11 should be carefully adopted. In this work, w ij is defined as follows: If the city is Wuhan, we have a ¼ 5000; otherwise, a ¼ 1.
In our model, each city has a set of unknown variables h j ¼ fa j ; b j ; c j ; j j ; d j g, which controls the size of infected population, spreading rate, and the death rate.
Note that this model has 5K þ 3 parameters, namely the proposed model has a high-dimension unknown parameter set, which should be optimized. The search space for these optimization problems may be highly nonlinear and contain many local minima. Evolutionary algorithms have been extensively used in nonlinear optimization [25,26,48]. In this article, we proposed a new pseudocoevolutionary algorithm to solve this inverse engineering problem. This main procedure tunes all the parameters, while the other two processes tune part of the parameters. The parameter estimation problem is separated into three co-adapted SA-based optimization processes. This main procedure tunes all the parameters, while the other two processes tune part of the parameters. Algorithm 3 Pesudocoevolutionary algorithm for estimating the optimal parameter set Θ * Input: The set of unknown parameters and initial number of infected and exposed individuals of each city Θ = {I W H,0 , E W H,0 , k I , θ 1 , θ 2 , · · · , θ K }, where θ j = {α j , β j , γ j , κ j , δ j } and j = 1, 2, · · · K; Output: Optimal parameter set Θ * ; Initialisation : Initialize temperature T , and random starting point The index of adopted cities is Φ = 1, 2, · · · , K. 1: Apply Algorithm-2 to optimize parameter set to achieveΘ * . Θ 0 ←Θ * LOOP Process 2: for i = 0 to M do 3: for j = 1 to K do 4: Set the model parameter as Θ 0 and apply Algorithm-1 to derive

5:
Using evaluation criteria (13) to derive RMSPE j for each city.

6:
Find the index of the M 1 largest RMSPE j and set Φ = {a 1 , a 2 , · · · , a M 1 }, where a i represents the index of the i-th largest RMSPE j .

9:
Apply Algorithm-2 to optimize parameter set to achieveΘ * . 1. In the main procedure, we tune the all the 5K þ 3 parameters. Then, we adopt root mean square percentage error (RMSPE) as follows: Here, RMSPE j is utilized as the criterion to measure the difference between the real daily infection data and the estimated infected individuals generated by this extended SEIR-migration model with an optimal parameter set H; 2. In this process, we find the index of the M 1 largest RMSPE j and only tune the parameter sets of the corresponding cities. 3. In order to avoid that the parameters of some cities have not been adopted and adjusted individually during 01Jan20 29Jan20 27Feb20 27Mar20 25Apr20 24May20 22Jun20 the whole identification process, we randomly select M 2 cities and adjust their parameters.
The whole optimization procedure is summarized in Algorithm-3, while Algorithm-2 is utilized for searching optimal parameters in subspace in each step.

Experimental results
The National Health Committee of China has published data on the spread of the COVID-19 epidemic from January 20, 2020. We use these historical data for parameter estimation of the SEIR-migration model (6). The pseudocoevolutionary simulated annealing algorithm, as described in Sect. 4.3, is adopted to find the optimal parameter set of this model. Since the parameters of the model are all estimated from historical data through reverse engineering techniques, the accuracy and integrity of the data are essential. In the early stage of COVID-19, people do not know much about this virus, and the diagnostic techniques are limited. Therefore, during the early stage of the outbreak of COVID-19, the historical data of Wuhan City are possible to deviate from the true value in a wide range. Hence, we reduce the weighting coefficients corresponding to Wuhan data in the objective function. In addition, after the outbreak of the epidemic, Chinese government has adopted effective quarantine measures and promoted epidemic prevention knowledge. These measures can effectively reduce the infection rate in each city. Therefore, the parameters in the model should be time-varying. Nevertheless, to simplify the calculation, we assume that these parameters are constant over the whole process. We also applied traditional simulated annealing, particle swarm optimization, genetic algorithm, and pattern search method to estimate the parameters of the model. However, these methods cannot provide a satisfied result or cannot even converge in an acceptable computation time (such as one day), while the proposed method can converge to the global optima in two hours.
This model can estimate the daily number of infected, exposed, and recovered individuals in all 367 cities. Due to space limitations, we only show the results of 17 cities in Fig. 5. Assuming that the migration control measures, infection rate, and recovery rate will remain unchanged for a period of time in the future, this model can provide a prediction of the amount of actively infected individuals in each city, as shown in 5. Results clearly show that the number of actively infected individuals has reached or peaked in most cities in China from February 20 to March 15, 2020. The prediction results show that there will be few new confirmed cases from early March, and the number of actively infected individuals will gradually decrease. The pseudocoevolutionary simulated annealing algorithm successfully finds the optimal parameter set: a ¼ 0:5631 AE 0:0161, b ¼ 0:2742 AE 0:0178, c ¼ 0:0509AE 0:0007, j ¼ 0:1310 AE 0:0035. Figure 6a shows the estimated peak number of infected individual in each province, while Fig. 6b reveals the estimated total number of infected individuals in each province. Results show that the number of infections in most provinces is smaller than 1000. It is reasonable to assume that our knowledge of the COVID-19 will gradually increase. Then, the medical treatment will improve and the recovery rate will increase each day. We assume the recovery rate will increase 0.0005 each day, namely every 20 days, the number of daily recovered individuals increases by 1% of the total number of infected individuals. Then, most cities in China will have almost zero infections before July 2020. Therefore, we can claim that China has already controlled the spreading of COVID-19 in China.

Conclusion
The novel Coronavirus Disease 2019 (COVID-19) epidemic has caused 75,204 confirmed cases and 2,006 death utile February 19, 2020, which triggers global public health concern. Peak prediction informing social and non- Esitimated total number of infections (b) Fig. 6 a Peak number of infections in each province; b estimated total number of infected individuals eventually infected in a province pharmaceutical prevention interventions is illuminating but remains difficult to achieve accuracy. Intercity migration is one of the essential factors in the spread of the disease. We construct migration networks of 367 cities in China. An SEIR-migration model with more than 1800 unknown parameters is utilized to model the spreading of COVID-19 in the 367 cities in China. We proposed a pseudocoevolutionary simulated annealing algorithm to identify these unknown parameters from historical data of the number of infected, recovered, and death toll. From this model, we can achieve all the essential information about the epidemic spreading, including infection rates, recovery rates, and eventual percentage of the infected population for 367 cities in China. The main conclusion of our study is that the COVID-19 epidemic spreading would peak between mid-February to early March 2020, with about 0:8%, less than 0:1%, and less than 0:01% of the population eventually infected in Wuhan, Hubei Province, and the rest of China, respectively. Results indicate that the COVID-19 epidemic has been controlled. This work provides a method for estimating the proportion of infected people. However, only seroprevalence studies may actually estimate the proportion of infected individuals.