Generalized logistic growth modeling of the COVID-19 outbreak in 29 provinces in China and in the rest of the world

The COVID-19 has been successfully contained in China but is spreading all over the world. We calibrate the logistic growth model, the generalized logistic growth model, the generalized growth model and the generalized Richards model to the reported number of infected cases from Jan. 19 to March 10 for the whole of China, 29 provinces in China, four severely affected countries (Japan, South Korea, Iran, Italy) and Europe as a whole. We dissect the development of the epidemics in China and the impact of the drastic control measures, and use this experience from China to analyze the calibration of other countries. The different models provide different future projections including upper and lower bounds of scenarios. We quantitatively document four phases of the outbreak in China with a detailed analysis on the heterogeneous situations across provinces. The extreme containment measures implemented by China were very effective with some instructive variations across provinces. However for other countries, it is almost inevitable to see the continuation of the outbreak in the coming months. Japan and Italy are in serious situations with no short-term end to the outbreak to be expected. There is a significant risk concerning the upcoming July 2020 Summer Olympics in Tokyo. Iran's situation is highly uncertain with unclear and negative future scenarios, while South Korea is approaching the end of the outbreak. Both Europe and the USA are at early stages of the outbreak, posing significant health and economic risks to the world in the absence of serious measures.


Background
Starting from Hubei province in China, the novel coronavirus (COVID -19) has been spreading all over the world, after two months of outbreak in China. Facing uncertainty and irresolution in A lot of efforts have been made in estimating the basic reproduction number R0 and predict the future trajectory of the coronavirus (COVID-2019) outbreak in the first quarter of 2020. In this paper, we focus on using phenomenological models without detailed microscopic foundations, but which have the advantage of allowing simple calibrations to the empirical reported data and providing transparent interpretations. This simple and top-down method can provide straightforward insights regarding the status of the epidemics and future scenarios of the outbreak.
Usually, an epidemic follows an exponential growth at an early stage (following the law of proportional growth), peaks and then the growth rate decays as countermeasures to hinder the transmission of the virus are introduced.
There have been quite an extensive literature reporting statistical analysis and future scenarios based on phenomenological models. Most of previous work use simple exponential growth models and focus on the early growing process [1][2][3][4]. On the other hand, there are also many works arguing that the number of infected people follows a trajectory different from a simple exponential growth [5][6][7][8][9][10][11][12][13][14][15].
In this paper, we employ the logistic growth model, the generalized logistic growth model, the generalized growth model and the generalized Richards model, which have been successfully applied to describe previous epidemics [16][17][18][19][20]. All these models have some limitations and are only applicable in some stages of the outbreak, or when enough data points are available. Thus, we first calibrate different models to the reported number of infected cases in the COVID-19 epidemics from Jan. 19 to March 10 for the whole of China and 29 provinces in mainland China, and then draw some lessons useful to interpret the results of a similar modeling exercise performed on the four countries that are undergoing major outbreaks of this virus: Japan, South Korea, Iran, and Italy. Our analysis dissects the development of the epidemics in China and the impact of the drastic control measures both at the aggregate level and within each province. Borrowing from the experience of China, we made projections on the development of the outbreak in the four key countries and the whole Europe, based on different scenarios provided by the results from different models. Our study employs simple models to quantitatively document the effects of the Chinese containment measures against the COVID-19, and provide informative implications for the coming pandemic.

Data
Confirmed cases: In this report, we focus on the daily data of confirmed cases in provinces in mainland China. We exclude the epicenter province, Hubei, which had a significant issue of underreporting at the early stage and also data inconsistency during mid-Feb due to a change of classification guidelines. For the provinces other than Hubei, the data is consistent except for one special event on Feb 20 concerning the data coming from several prisons. The data source is the national and provincial heath commission. For international data, the source is WHO. Note that the cases of the Diamond Princess cruise are excluded from Japan, following the WHO standard.
Data adjustment: On Feb 20, for the first time, infected cases in the Chinese prison system were reported, including 271 cases from Hubei, 207 cases from Shandong, 34 cases from Zhejiang.
These cases were concealed before because the prison system was not within the coverage of each provincial health commission system. Given that the prison system is relatively independent and the cases are limited, We remove these cases in our data for the modelling analysis to ensure consistency.
Migration data: the population travels from Hubei and Wuhan to other provinces from Jan 1 st to Jan 23 rd are retrieved from the Baidu Migration Map (http://qianxi.baidu.com).

Method
We use the generalized Richards model (GRM), which is an extension of the original Richards growth model [21]. With three free parameters, the Richards growth model has been fitted to a range of logistic-type epidemic curves [16]. The generalized Richards model is defined by the differential equation: where C(t) represents the cumulative number of cases at time t, r is the growth rate at the early stage, and K is the final epidemic size. ∈ [0,1] is a parameter that allows the model to capture different growth profiles including the constant incidence ( = 0), sub-exponential growth (0 < < 1 ) and exponential growth ( = 1 ). The exponent α measures the deviation from the symmetric s-shaped dynamics of the simple logistic curve. The model recovers the original Richards model for = 1, and reduces to the generalized logistic model [17] for α = 1 and =

1.
For the calibration, we use the standard Levenberg-Marquardt algorithm to solve the non-linear least square optimization. For the GRM, the initial point C(0) ≔ C 0 is fixed at the empirical value, and there are 4 remaining parameters ( , , , ) to be determined. For the fitting of the standard logistic growth function ( = 1 and α = 1), we free the initial point C 0 and allow it being one of the 3 parameters to be optimized, as the early stage growth does not follow a logistic growth.
To estimate the uncertainty of our model estimates, we use a bootstrap approach with a negative binomial error structure NB( , 2 ), where and 2 are the mean and variance of the distribution, estimated from the empirical data.

Analysis at the aggregate level of mainland China (excluding Hubei)
As of March 10 th , 2020, there are in total 13172 infected cases reported in the 30 provinces in mainland China outside Hubei. The initially impressive rising statistics has given place to a tapering associated with the limited capacity for transmission, exogenous control measure, and so on. In Figure 1, the trajectory of the total confirmed cases, the daily increase of confirmed cases, and the daily growth rate of confirmed cases in whole China excluding Hubei province are   presented. The fits with the generalized Richards model and with the standard logistic growth   model are shown in red and blue lines respectively in the upper, middle and lower left panel, with the data up to March 1 st , 2020. In the lower left panel of Figure 1, the daily empirical growth rate of the confirmed cases is plotted in log scale against time. We can observe two exponential decay regimes of the growth rate with two different decay parameters before and after This gives where 0 is a constant of integration determined from matching this asymptotic solution with the non-asymptotic dynamics far from the asymptote. Thus, the leading behavior of the growth rate at long times is  30, 25, 20, 15, 10 and 5 days before Feb 22 were presented. For each of the six data sets, we generated 500 simulations of ( ) based on the best fit parameters using parametric bootstrap with a negative binomial error structure, as in prior studies [17]. Each of these 500 simulations constitutes a plausible scenario for the daily number of new confirmed cases, which is compatible with the data and GRM. The dispersion among these 500 scenarios provides a measure of stability of the fits and their range of values gives an estimation of the confidence intervals. The first conclusion is the non-surprising large range of scenarios obtained when using data before the maximum, which however encompass the realized data. We observe a tendency for early scenarios to predict a much faster and larger number of new cases than observed, which could be expected in the absence of strong containment control. With more data, the scenarios become more accurate, especially when using realized data after the peak, and probably account now well for the impact of the containment measures that modified the dynamics of the epidemic spreading.

Analysis at the provincial level (29 provinces) of mainland China (excluding Hubei)
As of March 1, 2020, the daily increase of the number of confirmed cases in China excluding Hubei province has decreased to less than 10 cases per day. The preceding one-month extreme quarantine measures thus seems to have been very effective from an aggregate perspective. At this time, it is worthwhile to take a closer look at the provincial level to study the effectiveness of  Phase I (Jan 19 -Jan 24, 6 days): early stage outbreak. The data mainly reflects the situation before Jan 20, when no measures were implemented, or they were of limited scope. On Jan 19, Guangdong became the first province to declare a confirmed case outside Hubei in mainland China [22]. On Jan 20, with the speech of President Xi, all provinces started to react. As of Jan 24, 28 provinces reported confirmed cases with daily growth rates of confirmed cases ranging from 50% to more than 100%. In this period, all provinces continued to implement their strict measures, striving to bring the epidemics to an end. The growth rate of the number of confirmed cases declined exponentially with similar rates as in Phase II, pushing down the growth rate from 10% to 1%. In phase III, all provinces have passed the peak of the incidence curve, which allows us to obtain precise scenarios for the dynamics of the end of the outbreak from the model fits (

Quantification of the initial reactions and ramping up of control measures
On Jan 19, Guangdong was the first province to report a confirmed infected patient outside Hubei.
On Jan 20, 14 provinces reported their own first case. During Jan 21-23, another 14 provinces reported their first cases. If we determine the peak of the outbreak from the 5 days moving average of the incidence curve, then there are 15 provinces taking 7-11 days from their first case to their peak, 9 provinces taking 12-15 days, and 6 provinces taking more than 15 days. If we define the end of the outbreak as the day when the 5 days moving average of the growth rate becomes smaller than 1%, then 7 provinces spent 8-12 days from the peak to the end, 7 provinces spent 13-16 days, 13 provinces spent 17-20 days, and 2 provinces spent 21-22 days. For the six provinces that have the longest duration from the start of their outbreak to the peak (more than 15 days), it took 8-13 days for them to see the end of the outbreak (Figure 3). This means that these 6 provinces were able to control the local transmissions of the imported cases quite well, so that the secondary transmissions were limited. In contrast, 20 provinces took 28-31 days from the start to the end of the outbreak. Thus, those provinces that seem to have responded sluggishly during the early phase of the epidemics seem to have ramped up aggressively their countermeasures to achieve good results.

Quality of fits with the logistic equation as proxy to early under-reporting
Both the generalized Richards model and the logistic growth model can capture most of the dynamics at the cumulative level, daily increase level (1 st derivative) and the daily growth rate level (2 nd derivative).
However, likely due to the potential underreporting at the beginning, the logistic growth curve fails to capture the growth dynamics at the early stage in most provinces. The logistic growth provides a better quality of fit after Feb 1 st in general, as we can see from the lower right panel in the figure for each province (see supplementary material).

Diagnostic of the efficiency of control measures from the exponential decay of the growth rate of infected cases
The 10 most infected provinces (Guangdong, Henan, Zhejiang, Hunan, Anhui, Jiangxi, Shandong, Jiangsu, Chongqing, Sichuan) have done quite well in controlling the transmission, as indicated by the fact that their daily growth rates follow well-defined exponential decays, with a 2 larger than 90%. This exponential decay continued for all ten cities until the situation was completely under control during Feb 15-18, when the daily incidence was at near zero or a single-digit number.
Eight out of these ten provinces have a decreasing parameter of the exponential decay of the growth rate ranging from 0.142 to 0.173, similar to what is observed at the national average level (0.157).
Note that this exponential decay can be inferred from the generalized Richards model, as we noted in section 4.1.

Zhejiang and Henan exemplary developments
Zhejiang and Henan are the 2 nd and 3 rd most infected provinces but have the fastest decaying speed of the incidence growth rate (decay parameter for Zhejiang: 0.223, Hunan: 0.186) among the most infected provinces. This is consistent with the fast and strong control measures enforced by both provincial governments, which have been praised a lot on Chinese social network [23,24]. implemented among all provinces.

Heterogeneity of the development of the epidemic and responses across various provinces
Less infected provinces exhibit a larger variance in the decaying process of the growth rate.
However, we also see good examples like Shanghai, Fujian and Shanxi, which were able to reduce the growth rate consistently with a low variance. These provinces benefited from experience obtained in the fight against the 2003 SARS outbreak or enjoy richer local medical resources [7].
This enabled the government to identify as many infected/suspected cases as possible in order to contain continuously the local transmissions. Bad examples include Heilongjiang, Jilin, Tianjin, Gansu, which is consistent with the analysis of [7].
Most provinces have a small parameter of of the GRM (see equation (1)) and an exponent α Heilongjiang has been criticized a lot for its high infected cases and death rate (2.7%), given that it is far from Hubei and does not have a large number of migrating people from Hubei.

Initial and total confirmed numbers of infected cases correlated with travel index
The initial value C 0 of the logistic equation could be used as an indicator of the early number of cases, reflecting the level of early contamination from Hubei province as the epicenter of the outbreak. To support this proposition, the upper panel of Figure 4 plots the estimated C 0 as a function of the migration index from Hubei & Wuhan to each province. The migration index is calculated as equal to 25% of the population migrating from Hubei (excluding Wuhan) plus 75% of the population migrating from Wuhan, given that Wuhan was the epicenter and the risks from the Hubei region excluding Wuhan are lower. One can observe a clear positive correlation between the estimated C 0 and the migration index. The lower panel of figure 4 shows an even stronger correlation between the total number of cases recorded on March 6 st and the travel index, expressing that a strong start of the epidemics predicts a larger number of cases, which is augmented by infections resulting from migrations out of the epidemic epicenter.

Analysis of the recent development in Japan, South Korea, Iran and Italy and the world
Although it is clear that the virus in China has been successfully controlled up to the time of finalizing this paper, it has spread to 117 territories in the world as of March 10, with Italy, South Korea, Japan, and Iran as the four major hotspots. The highly connected European continent is at the beginning phase of the outbreak. In this section, we will analyze the situation in the four countries that have recently reported increasing numbers, and provide forecasts based on several models.

Models
Most of the countries analyzed here are still at an early stage of the outbreak, i.e., the cumulative number of cases is not yet reaching the inflection point. In this regime, the generalized Richards model, which has an additional parameter α describing the asymmetry between growth and decay of the incidence curve, is not sufficiently constrained by the limited available data. Therefore, we consider the following three simpler models with fewer parameters: -Logistic growing model: -Generalized Logistic model (GLM): -Generalized growth model (GGM): The generalized growth model (6) allows for a sub-exponential growth of the outbreak in the early stage (for p < 1), but cannot describe the decay of the incidence rate. It thus serves as a rough upper limit, obtained by assuming that the outbreak continues to grow following the same process as in the past. The generalized Logistic model (5) and Logistic growth model (4) both assume a logistic decay of the growth rate as the total number of confirmed cases increases. Note that they are nested in the generalized Richards model (1), by fixing respectively = 1 and = 1; p = 1 .
Both Logistic type models tend to under-estimate the total size of the infected population at the early stage, so they should provide lower bounds. The generalized Logistic model allows for an early sub-exponential growth and can better describe the possible asymmetry between growth and decay dynamics. We now calibrate these three models to the data of each of the four countries: Italy, South Korea, Japan, and Iran. The fitting procedure and simulation of confidence intervals are the same as for the generalized Richards model presented for the Chinese data sets in previous sections. However, it took China two weeks to reduce the growth rate from 10% to 1% under extreme containment measures. Given that we do not expect the same level of measures can be implemented in Japan sooner or later, the future scenario for Japan is highly uncertain, and will depend on whether the government decides to increase the containment measures. Without stronger and fast measures introduced in Japan, we see a significant risk concerning the upcoming July 2020 Summer Olympics in Tokyo. The fitted data is plotted every three days, and the error bar is the 95% confidence interval extracted from the 500 simulations assuming a negative binomial error structure. The fitting is based on data since Feb 12, 2020. Figure 6. Same as figure 5 for South Korea. The fitting is based on data since Feb 18, 2020.

South Korea
In contrast to Japan, the analysis of South Korea's data indicates that the transmission is under control, although the total number of confirmed cases is much larger than Japan. This is likely due to the much larger efforts to test the population, with a ratio of close to 4 tests per 1000 inhabitants, compared with 0.066 tests per 1000 inhabitants for Japan. Figure 6 shows that the dynamics of the infected population is clearly not an exponential growth, as the growth rate has decreased significant from 40% at the end of Feb to 1.  could be expected to be reached within 20 days. However, whether the daily incidence can decay symmetrically to the increase phase will depend largely on the strategy of the government. In China, it was mainly the extreme measures of quarantines that helped most provinces' incidences to decay fast. The South Korean government has been actively tracing all possible cases and performing more than 10,000 tests per day as of May 8 [26], trying to cut off all the potential transmission sources.

Italy
The situation in Italy is quite severe at the moment of writing. As of March 10, Italy has the largest number of confirmed cases (9172) and deaths (463) among all the countries except China. The outbreak in Italy was mainly due to a few outbreak clusters, such as Lombardy and Veneto. Similar to Japan, the growth rate of the confirmed cases at Italy has been fluctuating between 20% and 30% for a week, which is still far from the inflection point.
The generalized logistic model provides the best scenario, i.e. predicting that the inflection point will be reached in 10 days. This scenario will already lead to a total cumulative number of infected until the end of the outbreak.
In contrast, a pessimistic scenario provided by the Logistic model predicts that there could be up to 2 million infected cases adding up until the end of May.
In any case, the situation in Italy is predicted to become worse than at the epicenter province Hubei   Figure 9 shows the number of confirmed cases per million people in the four hotspot countries, in comparison with mainland China excluding Hubei, the epicenter Hubei province, Europe, Singapore, Hong Kong, and the USA. Although Hubei has the highest infection ratio above 0.1%, South Korea and Italy have become the most infected country in the world. The high statistics in South Korea is likely as mentioned above due to the largest number of performed tests. We should indeed stress that the actual cumulative number of infected cases is probably misleading as it is likely that the incidence is an increasing function of testing.

Europe and the world
Since it is delicate to interpret the absolute number of infected cases, we propose to rather focus on the dynamics, which is more informative as we have presented above. For South Korea, we have shown that it is approaching the ceiling in the total number of infected cases and we do not expect this number to grow significantly in the coming future. If we assume a generalized logistic model for the behavior of Italy, the estimated final fraction of population infected will be 0.15% (95% CI: [0.03%, 0.30%]), which will be much higher than Hubei. Nevertheless, the estimates from the generalized growth model and the generalized logistic model can provide a negative and positive scenario respectively. The generalized growth model assumes that the future growth will continue following the same path as the past weeks, which could be expected as the outbreak is just starting. This scenario will bring the total number of confirmed

Conclusion
In this paper, we have calibrated the logistic growth model, the generalized logistic growth model, We quantified the initial reactions and ramping up of control measures on the dynamics of the epidemics and unearth an inverse relationship between the number of days from peak to the quasiend and the duration from start to the peak of the epidemic among the 29 analyzed Chinese provinces. We identified the dynamic signatures of the exemplary developments in Zhejiang and Henan provinces and the heterogeneity of the development of the epidemic and responses across various other provinces. We found a strong correlation between the initial and total confirmed numbers of infected cases and travel index quantifying the mobility between provinces.
For countries that are in the middle of the outbreak, we used some useful experience from China to interpret the calibration results from different models, and made future scenario projections. For Japan, the situation is found similar to phase II-III in mainland China and the estimated total number of cases as of March  which will be higher than Hubei, the epicenter in China of the epidemic. Our statistical analysis showed that the epidemic in Europe is in the early exponential regime and this is likely to continue for a while. This negative but probable scenario will provide 114867 people infected in Europe in 10 days, corresponding to 0.015% European population. The generalized logistic model provides a lower (optimistic) bound, leading to a prediction of a total cumulative number of 164219 people infected, corresponding to 0.02% (95% CI: [0.005%, 0.044%]) of the European population.
However, for this forecast to come true, all European countries need to coordinate and fight together to avoid another major outbreak like in Italy. The situation in the United States is also at the very beginning stage, and may pose significant health and economic risks to the world in absence of serious measures.