Predicting the evolution and control of the COVID-19 pandemic in Portugal

Coronavirus disease 2019 (COVID-19) is a worldwide pandemic that has been affecting Portugal since 2 March 2020. The Portuguese government has been making efforts to contradict the exponential growth through lockdown, social distancing and the usage of masks. However, these measures have been implemented without controlling the compliance degree and how much is necessary to achieve an effective control. To address this issue, we developed a mathematical model to estimate the strength of Government-Imposed Measures (GIM) and predict the impact of the degree of compliance on the number of infected cases and peak of infection. We estimate the peak to be around 650 thousand infected cases with 53 thousand requiring hospital care by the beginning of May if no measures were taken. The model shows that the population compliance of the GIM was gradual between 30% to 75%, contributing to a significant reduction on the infection peak and mortality. Importantly, our simulations show that the infection burden could have been further reduced if the population followed the GIM immediately after their release on 18 March.


Introduction
Coronavirus disease 2019 (COVID-19) is already considered a world pandemic which is starting to have dramatic effects in Europe, where, as of 27 of March, 265,421 cases have been reported 1,2 . COVID-19 infection in Portugal has been growing exponentially with an average rate of 34±13% new cases per day from 2 March and is far from reaching the peak by the end of March. As of March 27, 4268 infection cases and 76 deaths have been reported 2 . The highest infection burden is found in Porto (317 cases, 7.4%) and in Lisbon (284 cases, 6.7%) but the disease is present throughout the entire country. As in other countries, infection occurs mostly in individuals' with ≥40 years of age (71.9% males; 69.3% females). Death occurs mostly in males (64.5%) all with ≥50 years of age.
Predictive models estimate that the peak of SARS-CoV-2 infection globally would be between mid-April and May, with an estimated total of 48 million people infected 3 . As with most other countries, the Portuguese national health care system cannot deal with the increasing demand of care due to limited ventilators and care units 3 . Therefore, the Portuguese government together with the National Health Directorate (DGS) declared a state of emergency and adopted interventive populational measures through Government-Imposed Measures (GIM) on 18 March 2020 in an attempt to drop the peak of infections even if at the cost of prolonging the infection time. These measures are based on the lockdown of people at home, social distancing and adopting protective antiseptic policies such as the usage of masks. Lockdown was implemented to assure compliance of the population, expect for people that maintain basic services such as medical and food distribution staff.
Most forecasting models are based on the number of cases reported and do not take into account the effects of these government-imposed measures and behavioral change. Thus, accessing the compliance degree and predicting how much is necessary for the control of SARS-CoV2 infection would be a useful tool for fighting COVID-19 pandemic. Recently published mathematical modelling studies of COVID-19 transmission have already provided useful insights that can be used to guide public health measures and resource allocation to better control this pandemic [4][5][6][7] . However, most parameters of statistical models have been estimated with high degree of uncertainty, resulting in predictions with wide intervals of confidence 4,6 . Compartmental models such as susceptible, infected and resistant (SIR) models are deterministic approaches based on solving nonlinear systems of Ordinary Differential Equations (ODE) that have been successful in describing complex dynamics of virus infection in populations, including COVID-19 in several countries [7][8][9][10][11] . Here, we provide a simple SI model that describe the dynamics of transition of COVID-19 in Portugal during the first 21 days, explain the evolution of SARS-CoV-2 infection dynamics up to 19 of August and predicts the degree of compliance of GIM by the Portuguese population.

Methods
Basic transmission dynamics of COVID-19 was modelled using a simple mathematical model based on a system of two ordinary differential equations (ODE) developed specifically for this purpose (Equation 1 and Equation 2). The equations reflect the number of people infected (I) and susceptible (S) to infection per unit of time (dI/ dt and dS/ dt). In this model, we accounted for the reported average time of duration of infection (τ) of 14 days 4,11 . The model was calibrated by adjusting the rate constant (k) to approximate the total infection value reported by the DGS at 17 March. No further fitting was performed in this model. The compliance of GIM by different fractions of the population was modelled through the variation of parameter α in Equation 1 and Equation 2. We considered that these protective measures (GIM) were 97% effective based on recent meta-analysis estimates, accounted through model parameter β 12 . The ODEs were encoded and solved using PLAS software version 1.2.0.120, where a series of simulations were carried scanning various values of the α parameter 13 . Simulations were carried with the initial two cases reported by the DGS and considering only the population of the grand Lisbon and Porto areas (total of 6.5 × 10 6 ) since they represent most of the susceptible population (see Figure 2). For simulations, we used the numerical solver based on the Adams/BDF method, implemented in the LSODA routine of PLAS software. Because a serological screening study made by the Portuguese Nacional Institute of Health (http://www.insa.pt/) found a 6-fold higher infected due to untested asymptomatic exposed to SARS-CoV-2, we have considered this ratio to estimate the reported symptomatic infected by the DGS. Further analysis, computations and plots were conducted using Python 3 in the Jupiter Notebook ipython 7.8.0 programing environment under Anaconda distribution version 4.7.12. Data regarding the daily evolution of number of total infected in Portugal by COVID-19 was collected from the DGS web site (https://covid19.min-saude.pt/ponto-de-situacao-atualem-portugal/) from 2 March to 19 August 2020 (see Source data, Table S1 and Figure S1) 14 . The model is available as Extended data.

Updates from Version 1
As relevant new data on the evolution of COVID-19 pandemic in Portugal give new insights and challenged our model. We have decided to take these data (up to 19 August 2020), improve our model and provide a retrospective insight on the evolution and control of this outbreak in Portugal. We believe that this new version may help to have a better understanding on the evolution of COVID-19 pandemic in Portugal and provide a simple approach to guide on the control and prevention of future outbreaks.

Results and discussion
Simulation of the first 18 days with our model was able to describe the exponential increase of the number of confirmed cases reported by the DGS between 2 and 18 March 2020 ( Figure 1). The predicted peak time for this scenario was 49 days which would be by the 21 of April. This is within the estimated range predicted by statistical modelling of US, Italy and Korea scenarios 3 . Further, the predicted numbers of cases for the end of March if no measures were taken would be around 42,000. This is also in agreement with the number released by the DGS to the social media based on statistical modelling. Thus, the model presented here is consistent with the forecasting made by conven tional models, reinforcing the confidence on our model capacity to generate predictions.
Importantly, our results show that the GIM had an immediate impact on diminishing the exponential increase of the number of infected cases and this depends on the percentage of the population that is in compliance with such measures ( Figure 2). This is evident by the increasing deviation of the reported number of cases relative to the unperturbed simulation (0%) with time. The evolution of the number of cases reported by DGS between 18 and 25 March fit between the simulation curves corresponding to 30% and 40% of model perturbation on parameter α This suggests that the estimated percentage of the Portuguese population that have been start following the GIM was between 30% to 40%. From simulations, we identify other intervals (e.g. 50-60% and 70-75%) that are compatible with the reported data form DGS between April to August 2020, regarding observed peaks of infection and hospitalizations ( Figure 2, Based on the fraction of hospitalized and mortality reported by the DGS on 27 March 2020, together with our model  predictions, we computed several infection indicators for these intervals (Table 1).
Our model analysis indicates that current government-mandated measures together with compliance changes shifted at least two times the expected peak of infections, causing a substantial reduction in the infection numbers ( Figure 2, Table 1). Based on our model, the predicted peak in the number of cases without any interventive measures would be around 650 thousand, whereas current degree of compliance (70-75%) have resulted in a decrease around one half of expected cases, hospitalizations and deaths (  Figure 2). This scenario would result in much less total mortality and hospitalization requirements on peak in comparison to the current trend (Table 1, Figure 2). Meanwhile, percentages >75% comes with the burden of prolonging the time of pandemic control over a year, which can be economically unbearable. Thus, the ideal solutions would be between 70-90% compliance of the GIM. The results obtained during simulations are available as Extended data, Table S2 14 .
Although our model precisely described the exponential curve and explains the shift in the temporal evolution of DGS data, it has limitations that may compromise the exact values of predictions. The fact that we only assume two compartments (susceptible and infected) considering the main populated cities (Lisbon and Porto) as one is huge approximation that neglects regional dynamics. Thus, the model is just an approximation that reflects an average trend and may fail to explain regional observations. In this model we also neglected many important parameters of infection transmission such as age groups, types of social interactions, contact dependent probability, and viral load dependent probability 15 . The inclusion of these parameters would definitely make the model more realistic. However, this data is not available for the Portuguese case and these models require accurate processing of data curation for suitable validation. We have bypassed these limitations by aggregating all of these parameters into one constant, which was fitted to the available data. Overall, the predictions shown here should be taken as semi-quantitative estimates within an upper and lower case-scenario.

Conclusions
In this work we demonstrate the potential of modelling the dynamics of SARS-CoV-2 infection as a useful support tool for predicting the impact of corrective measures as well as estimating the degree of compliance of the GIM by the population. Government-mandated measures on the Portuguese population effectively prevented COVID-19 from reaching dramatic numbers in Portugal but still could be substantially improved to reduce the infection peak. Our estimates and approach may help in guiding additional measures to control the COVID-19 evolution and future epidemies.

Open Peer Review
Notably, the authors now conveyed the idea of estimating the strength of the lockdown in Portugal, during early COVID-19 pandemic. Moreover, the authors improved simulations by taking the estimated number of asymptomatic cases into account. Also, they discussed their predicted number of SARS-CoV-2 infection in the light of newly published data by the Portuguese authorities.
It would be great if the authors could define (biologically or epidemiologically) the parameters A (alpha) and B (beta) (equations 1 and 2). The authors should include, if possible, a citation to the sentence 'This is also in agreement with the number released by the DGS to the social media based on statistical modelling'. The closing punctuation is missing after 'The evolution of the number of cases reported by DGS between 18 and 25 March fit between the simulation curves corresponding to 30% and 40% of model perturbation on parameter α'.
Finally, I advise re-writing the sentence 'For 19 august, the computed total infections and deaths for the 70-75% interval is 30,664 -91,426 and 1,004 -2,995, respectively.' The current version is quite confusing as readers may interpret '70-75% interval' as the confidence interval. It should be clear that these values represent the percentage of lockdown.

Competing Interests:
No competing interests were disclosed.
Reviewer Expertise: Genetics, bioinformatics, epidemiology, medical entomology, Wolbachia, and host-microorganism interactions I confirm that I have read this submission and believe that I have an appropriate level of expertise to confirm that it is of an acceptable scientific standard.

Reviewer Report 14 September 2020
What is the feasible region for the considered model? 2.
Also, simulate the model for a long time that is months, 40 days, etc. 3.
What techniques for numerical simulation have been used I confirm that I have read this submission and believe that I have an appropriate level of expertise to confirm that it is of an acceptable scientific standard, however I have significant reservations, as outlined above.
Author Response 25 Aug 2020

Ricardo Pais, Instituto Universitário Egas Moniz, Caparica, Portugal
We would like to acknowledge the reviewer for the relevant comments and suggestions. We have considered them all and have revised the manuscript accordingly. Comment 2 "the consistently wrote 'transmission dynamics of COVID-19' and not the correct form 'transmission dynamics of SARS-CoV-2" Reply 2 We agree with the suggestion and have changed the manuscript replacing COVID-19 with SARS-CoV-2 (see revised version). Comment 3 "The authors predicted 2 million cases, which I find surprisingly high. Using the optimal isolation percentage (70-75%), the model still predicted over 150,000 cases which are approximately 5-fold higher than the current number of cases in Portugal (~35,000)" Reply 3 Indeed, our predictions are higher than the reported values. This is because only a small fraction of SARS-COV-2-infected individuals are tested for the virus has only a few patients show symptoms. However, our results are in agreement with the recent serologic study conducted by the National Institute of Health (INSA) from Portugal that found that a total of 300.000 people were exposed to the virus, 6-fold higher than the number of reported cases. By the time we conducted the model calibration and analysis, these results were not known, resulting in the deviation between predictions and reported values. Thus, we have corrected the model accounting for the asymptomatic fraction not tested and other recent data. This resulted in novel results which were compared with new data up to August and presented in the revised manuscript. Comment 4 "I would like to see the model's prediction using even higher percentages of isolation. This should be also discussed."

Reply 4
We have simulated the model with higher percentages and discuss the results (see revised manuscript).

Competing Interests:
No competing interests to disclosure.
Author Response 25 Aug 2020

Ricardo Pais, Instituto Universitário Egas Moniz, Caparica, Portugal
We would like to acknowledge the reviewer for the relevant comments and suggestions. We have considered them all and have revised the manuscript accordingly. Comment 1 "Provide the existence of the model" Reply 1 According to F1000 Research rules for supplementary data, we included the model code, data and the analysis code in python as extended on figshare (see revised version references).
Comment 2 "What is the feasible region for the considered model?"

Reply 2
The recommended feasible region for this model is 10 8 > (S+I) > 10 4 where initial S + I = city population.
Comment 3 "Also, simulate the model for a long time that is months, 40 days, etc."

Reply 3
We have performed simulations for 500 days for predictions and 17 days for model calibrations. By the time we submitted the paper no more than 27 days of data was available. Indeed, we only used a 17 days simulation for model calibration to estimate a rate constant in the absence of control. Unfortunately, it is not possible to use more data for model calibration since an uncontrolled lockdown was implemented in Portugal immediately after these 17 days of infection. Using more data would actually result in wrong estimates since many people follow the DGS recommendations but others do not. This is why we have simulated a total of 500 days with multiple % of lockdown scenarios towards estimating how much % of lockdown during the evolution of COVID pandemic. However, we now include relevant reported DGS data (up to August) for contrasting with model simulations and predictions (see revised version). Comment 4 "What techniques for numerical simulation have been used?"

Reply 4
The numerical solver was based on the Adams/BDF method, implemented in the LSODA routine of PLAS software. This is a general-purpose stiff, variable-step and variable-order solver. We add this information in the methods section (see revised version).