Scenario II. Results obtained based by taking twenty times the number of infected and forty times the number of recovered people with respect to the confirmed cases
For our illustrations, we assumed that the number of infected is twenty times the number of the confirmed infected and forty times the number of the confirmed recovered people. Based on this scenario, Fig 9 depicts an estimation of R0 for the period January 16-January 20. Using the first six days from the 11th of January to the 16th of January, R0^ results in 3.2 (90% CI: 2.4-4.0); using the data until January 17, R0^ results in 3.1 (90% CI: 2.5-3.7); using the data until January 18, R0^ results in 3.4 (90% CI: 2.9-3.9); using the data until January 19, R0^ results in 3.9 (90% CI: 3.3-4.5) and using the data until January 20, R0^ results in 4.5 (90% CI: 3.8-5.3).
10.1371/journal.pone.0230405.g009 Fig 9  Scenario II. Estimated values of the basic reproduction number (R0) as computed by least squares using a rolling window with initial date the 11th of January.
The solid line corresponds to the mean value and dashed lines to lower and upper 90% confidence intervals.   It is interesting to note that the above estimation of R0 is close enough to the one reported in other studies (see in the Introduction for a review).
Fig 10 depicts the estimated values of the case fatality (γ^) and case recovery (β^ ratios for the period January 16 to February 10. The confidence intervals are also depicted with dashed lines. Note that the large variation in the estimated values of β^ and γ^ should be attributed to the small size of the data and data uncertainty. This is also reflected in the corresponding confidence intervals. As more data are taken into account, this variation is significantly reduced. Thus,using all the (scaled) data from the 11th of January until the 10th of February, the estimated value of the case fatality ratio γ^ now drops to ∼ 0.147% (90% CI: 0.144%-0.15%) while that of the case recovery ratio is ∼ 0.1 (90% CI: 0.091-0.11). It is interesting also to note that as the available data become more, the estimated case recovery ratio increases slightly (see Fig 10), while the case fatality ratio (in the total population) seems to be stabilized at a rate of ∼ 0.15%.
10.1371/journal.pone.0230405.g010 Fig 10  Scenario II. Estimated values of case fatality (γ^) and case recovery (β^) ratios, as computed by least squares using a rolling window (see in Methodology).
Solid lines correspond to the mean values and dashed lines to lower and upper 90% confidence intervals.   In Figs 11, 12 and 13, we show the coefficients of determination (R2) and the root of mean squared errors (RMSE), for R0^, β^ and γ^, respectively.
10.1371/journal.pone.0230405.g011 Fig 11  Scenario II. Coefficient of determination (R2) and root mean square error (RMSE) resulting from the solution of the linear regression problem with least-squares for the basic reproduction number (R0).
10.1371/journal.pone.0230405.g012 Fig 12  Scenario II. Coefficient of determination (R2) and root mean square error (RMSE) resulting from the solution of the linear regression problem with least-squares for the recovery rate (β^).
10.1371/journal.pone.0230405.g013 Fig 13  Scenario II. Coefficient of determination (R2) and root mean square error (RMSE) resulting from the solution of the linear regression problem with least-squares for the mortality rate (γ^).   The computed values of the “effective” per day mortality and recovery rates of the SIRD model were γ ∼ 0.0005 and β ∼0.16d−1 (corresponding to a recovery period of ∼ 6 d). Note that because of the extremely small number of the data used, the confidence intervals have been disregarded. Instead, for calculating the corresponding lower and upper bounds in our simulations, we have taken intervals of 20% around the expected least squares solutions. Hence, for γ we have taken the interval (0.0004 and 0.0006) and for β, we have taken the interval between (0.13 and 0.19) corresponding to an interval of recovery periods from 5 to 8 days.
Again, we used the SIRD simulator to provide estimation of the infection rate by optimization setting w1 = 1, w2 = 400, w3 = 1 to balance the residuals of deaths with the scaled numbers of the infected and recovered cases. Thus, to find the optimal infection transmission rate, we used the SIRD simulations with β = 0.16d−1, and γ = 0.0005 and as initial conditions one infected, zero recovered, zero deaths on November 16th 2019, and ran until the 10th of February.
The optimal, with respect to the reported confirmed cases from the 11th of January to the 10th of February value of the infected rate (α) was found to be ∼ 0.319(90% CI: 0.318-0.32). This corresponds to a mean value of the basic reproduction number R0^≈2.
Finally, using the derived values of the parameters α, β, γ, we have run the SIRD simulator until the end of February. The simulation results are given in Figs 14, 15 and 16. Solid lines depict the evolution, when using the expected (mean) estimations and dashed lines illustrate the corresponding lower and upper bounds as computed at the limits of the confidence intervals of the estimated parameters.
10.1371/journal.pone.0230405.g014 Fig 14  Scenario II. Simulations until the 29th of February of the cumulative number of infected as obtained using the SIRD model.
Dots correspond to the number of confirmed cases from 16th of Jan to the 10th of February. The initial date of the simulations was the 16th of November with one infected, zero recovered and zero deaths. Solid lines correspond to the dynamics obtained using the estimated expected values of the epidemiological parameters α = 0.319, β = 0.16d−1, γ = 0.0005; dashed lines correspond to the lower and upper bounds derived by performing simulations on the limits of the confidence intervals of the parameters.
10.1371/journal.pone.0230405.g015 Fig 15  Scenario II. Simulations until the 29th of February of the cumulative number of recovered as obtained using the SIRD model.
Dots correspond to the number of confirmed cases from 16th of January to the 10th of February. The initial date of the simulations was the 16th of November, with one infected, zero recovered and zero deaths. Solid lines correspond to the dynamics obtained using the estimated expected values of the epidemiological parameters α = 0.319, β = 0.16d−1, γ = 0.0005; dashed lines correspond to the lower and upper bounds derived by performing simulations on the limits of the confidence intervals of the parameters.
10.1371/journal.pone.0230405.g016 Fig 16  Scenario II. Simulations until the 29th of February of the cumulative number of deaths as obtained using the SIRD model.
Dots correspond to the number of confirmed cases from the 16th of November to the 10th of February. The initial date of the simulations was the 16th of November with zero infected, zero recovered and zero deaths. Solid lines correspond to the dynamics obtained using the estimated expected values of the epidemiological parameters α = 0.319, β = 0.16d−1, γ = 0.0005; dashed lines correspond to the lower and upper bounds derived by performing simulations on the limits of the confidence intervals of the parameters.   Again as Figs 15 and 16 suggest, the forecast of the outbreak at the end of February, through the SIRD model is characterized by high uncertainty. In particular, in Scenario II, by February 29, simulations result in an expected actual number of ∼8m infected cases (corresponding to a ∼13% of the total population) with a lower bound at ∼720,000 and an upper bound at ∼37m cases. Similarly, for the recovered population, simulations result in an expected actual number of ∼4.5m (corresponding to a 8% of the total population), while the lower and upper bounds are at ∼430,000 and ∼23m, respectively. Finally, regarding the deaths, simulations under this scenario result in an average number of ∼14,000, with lower and upper bounds at ∼900 and ∼100,000.
Importantly, under this scenario, the simulations shown in Fig 14 suggest a decline of the outbreak at the end of February. Table 1 summarizes the above results for both scenarios.
10.1371/journal.pone.0230405.t001 Table 1  Model parameters, their computed values and forecasts for the Hubei province under two scenarios: (I) using the exact values of confirmed cases or (II) using estimations for infected and recovered (twenty and forty times the number of confirmed cases, respectively).     We note that the results derived under Scenario II seem to predict a slowdown of the outbreak in Hubei after the end of February.