On March 1 on ABC's This Week, US Secretary of Health and Human Services Alex Azar told host George Stephanopoulos that asymptomatic spread of COVID-19 is "not the major driver" of the spread of the new coronavirus, a sentiment corroborated on the CDC's website. Meanwhile other researchers disagreed, for example Dr. Osterholm, director of the Center for Infectious Disease Research and Policy at the University of Minnesota and Dr. Schaffner of the Vanderbilt University Medical School and adviser to the CDC (see the CNN article, Cohen (2020) ). At the time, there was clearly no unimodal opinion. What does the data say? There is a clear need for generating actionable intelligence on the nature of a pandemic from early and minimal data.
One approach is scenario analysis. Recently, Chinazzi et al. (2020) used the Global Epidemic and Mobility Model (GLEAM) to perform infection scenario analyses on China COVID-19 data, using a networked meta-population model based on transportation hubs. A similar model for the US was reported in Wilson (2020) where their web-app predicts from 150,000 to 1.4 million infected cases by April 30, depending on the intervention level. Such complex scenario analysis requires user input such as infection sites and properties of the contagion. However, a majority of infection sites may be hidden, especially if asymptomatic transmission is significant, and the contagion parameters are unknown and must be deduced, perhaps using domain expertise.
Data driven approaches are powerful. A long range regression analysis of COVID-19 out to 2025 on US data using immunological, epidemiological and seasonal effects is given in Kissler et al. (2020) , which predicts recurrent outbreaks. We also follow a pure data-driven machine learning approach to understand early dynamics of COVID-19 on the first 54 days of US confirmed infection reports (downloadable from the European Center For Disease Control). We address the challenge of real-time data-intelligence from early data. Our data-approach is simple, requires minimal data or human input and generates fast actionable insights. Is asymptomatic spread significant? Our data-driven analysis says yes, emphatically. What are the properties of the contagion, for example is it unusually virulent? Our analysis says no.
Early data is both a curse and a blessing. The curse is that "early" implies not much information, so quantitative models must be simple and robust to be identifiable from the data. The blessing is that early data is a peek at the pandemic as it roams free, unchecked by public health protocols, for to learn the true intentions of the lion, you must follow the beast on the savanna, not in the zoo. As we will see, much can be learned from early data and these insights early in the game, that are not based on opinion but data, are crucial to public health governance.
We analyze daily confirmed COVID-19 from January 21 to March 14. Our model is time stamped, hence we can provably test it moving forward. Our model predictions up to March 15 are in Figure 1 . Red circles are cumulative confirmed infections, the blue curve is our learned model and the blue circle is our prediction for March 15, 665 new confirmed infections. Qualitatively the model captures the data. The blue model curve is, in some sense, an optimal fit to the data, with an asymptomatic infectious force if left unchecked of 30% new cases per day (doubling in about 2.6 days) and a virulence of γ = 0.14% (1 or 2 in a thousand conversions from mild infection to serious). Not all serious cases are fatal. Such early intelligence has public health significance: stay calm, the virulence is not astonishingly high, but the spread rate is high, hence one must curtail that spread, or else a health-system can get overloaded.
Our goal is a simple robust model. The majority of disease models are based on aggregation of individuals according to disease status, such as SI, SIR, SIS Kermack and McKendrick (1927) ; Bailey (1957) ; Anderson and May (1992) . We use a similar model by considering a mild infected state which can transition to a serious state. Early data allows us to make simplifying assumptions. In the early phase, when public health protocols have not kicked in, a confirmed infection is selfreported. That is, you decide to get tested. Why? Because the condition became serious. This is important. A confirmed case is a transition from mild infection to serious. This is not true later when, for example, public health protocols may institute randomized testing. At time t let there be C(t) confirmed cases and correspondingly M (t) mild unreported asymptomatic infections. The new confirmed cases at time t correspond to mild infections at some earlier time t − k which have now transitioned to serious and hence got self-reported. Let us suppose that a fraction γ of those mild cases transitioned to serious,
Another advantage of early dynamics is that we may approximate the growth from each infection as independent and exponential, according to the infectious force of the disease. So,
Here, the second term is the loss of mild cases that transitioned to serious, the third term is the remaining cases that don't transition to serious recovering at some later time r and the fourth term accounts for new infections from confirmed cases. We will further simplify and assume that confirmed cases are fully quarantined and so α = 0 and recovery to a non-infectious state occurs late enough to not factor into early dynamics. Our simplified model is:
We set t = 0 at the first confirmed infection (Jan 21 in USA). There are four parameters.
β, the asymptomatic infectious force parameter governing the exponential spreading γ, the virulence, the fraction of mild cases that become serious some time later k, the lag time for mild asymptomatic infection to become serious (incubation time) M 0 , The number of unconfirmed mild infections at time 0
We assess the quality of the parameters β, γ, k, M 0 by how well they reproduce the observed dynamics in Figure 1 , using an error-measure, see for example Abu-Mostafa et al. (2012) . We used a combination of root-mean-squared-error and root-mean-squared-percentage-error between observed dynamics and model predictions. The model parameters which optimally fit the data are β * = 1.30 γ * = 0.0014 k * = 10 days
The asymptomatic infectious force, β, is very high, and corresponds to a doubling time of 2.6 days. On the other hand, the virulence at 0.14% is small, comparable to a standard flu. That is promising. The incubation period of 10 days seems in line with physician observations. The data analysis predicts that when the first confirmed case appeared, there were 4 other infections in the USA. The parameters β * , γ * and M * 0 are new knowledge, gained with relative ease by calibrating a simple robust model to the early dynamics. But, these optimal parameters are not the whole story, when it comes to prediction.
Let us walk through the process of obtaining the optimal parameters (β * , γ * , k * , M * 0 ). First we determined the incubation period k and the number of initial unseen mild infections M 0 . We exhaustively searched over (k, M 0 ) and for each pair of values, we used gradient descent with 4 million iterations to identify the optimal β and γ for that specific k, M 0 . We show the quality of fit for various (k, M 0 ) in Figure 2 (a). The deep-blue region contains essentially equivalent models within 0.5% of the optimal fit, our (user defined) error tolerance. The deep-blue region shows the extent to which the model is ill-identified by the data. Indeed, all these deep-blue models equally fit the data which results in a range of predictions. For robustness, we pick the white square in the middle of the deep-blue region, but note here that it is only one of the models which are consistent period and initial infections. Blue is better fit, red is worse fit. The deep-blue region corresponds to comparable models, within 0.5% of the optimal fit. The white square is the model chosen, (k * , M * 0 ) = (10, 4) which has optimal fit and also is "robust" being in the middle of the deep-blue regions. (b) Model fit for the chosen k * and M * 0 . Again, blue is a better fit and the deep-blue region is an equivalence class of optimal models. We choose the "robust" model corresponding to the white square in the middle of the deep-blue region, (β * , γ * ) = (1.30, 0.0014) The deep-blue regions represent uncertainty.
with the data. In making predictions, we should consider all equivalent models to get a range for the predictions that are all equally consistent with the data. Similarly, in Figure 2 (b), we fix k * and M * 0 to their optimal robust values and show the uncertainty with respect to β and γ (the deep-blue region). Again, we pick the white square in the middle of the deep-blue region of equivalent models with respect to the data. Hence we arrive at our optimal parameters in Equationeq:param-pot.
In the supplementary material we give a cross-sectional study across countries. The different countries have different cultures, social networks, mobility, demographics, as well as different times at which the first infection was reported (the "delay"). The infectious force β has significant variability, and we study how β statistically depends on a number of factors. We find:
• The larger the delay, the larger is β. The more that has been witnessed, the faster the spread.
That seems unusual, but is strongly indicated by the data. • Population density at the infection site has a strong positive effect but the country's population density does not. • There is faster spread in countries with more people under the poverty level defined as the percentage of people living on less than $5.50 per day. • Median age has a strong effect. Spread is faster in younger populations. The youth are more mobile and perhaps also more carefree. • Wealth and per-capita income have a slight negative effect. Spread is slower in richer countries, perhaps due to some combination of more risk-aversion, higher levels of education and less Figure 3 .
reliance on public transportation. Whatever the cause, it does have an impact, but relatively smaller than the other effects.
As mentioned in the previous section, to get honest estimates, we must consider all models which are equally consistent with the data, indicated by the deep-blue regions in Figure 2 . We arrive at the following ranges for our estimates, based on a 0.5% error tolerance. Asymptomatic infectionis strong, around 30%, however it is not that virulent, converting at most 1.2% of mild infections to serious. There is significant uncertainty in the lag, from 1 up to 13 days, and there could have been up to 12 unknown initial infections when the first case was confirmed.
As of March 14, the model suggests on the order of a million unknown infections. We now come to predictions beyond March 14, see Table 1 and Figure 3 . Data beyond March 14 was not used to calibrate the model. You may complain about the large uncertainty-range for each prediction. Unfortunately, that is what the data tells us. There is uncertainty, it's a fact of life. There are many models equally consistent with the data, resulting in a range of predictions. The model says that if the pandemic is left to its business unchecked, we should expect a total of about 29 thousand new confirmed infections by March 24, and the range is from 28 to 125 thousand.
We can learn a lot from the early dynamics of COVID-19. It's infectious force, virulence, incubation period, unseen infections and predictions of new confirmed cases. All this, albeit within error tolerances, from a simple model and little data. And, the information is useful for planning and health-system preparedness. A side benefit of the predictions in Table 1 , if the model is to be trusted, is as a benchmark against which to evaluate public health interventions. If moving forward, observed new infections are low compared to the data in Table 1 , it means the interventions are working by most likely reducing β. The US instuted broad and aggressive social distancing protocols starting on or before March 12. If our "lag" k = 10 is accurate, and the social distancing protocols have been effective at reducing β, the observed new infections should drop below the predictions starting on or before March 22. Without such a target to compare with, it would be hard to evaluate intervention protocols. Our approach is simple, general and works with coarse, aggregated data. We end with some limitations. First, the independent evolution of infection sites only applies to early dynamics. Hence, when the model infections increase beyond some point, and the pandemic begins to saturate the population, a more sophisticated network model that captures the dependence between infection sites would be needed Balcan et al. . While we did present an optimal model, it should be taken with a grain of salt because many models are nearly equivalent, resulting in prediction uncertainty. The model and the interpretation of its parameters will change once public health protocols kick in. The model may have to be re-calibrated (for example if β decreases) and the parameters may have to be reinterpreted (for example γ is a virulence only in the self-reporting setting, not for random testing). It is also possible to build a more general model with an early phase β 1 and a latter phase β 2 (after social distancing). Lastly, the model was learned on USA data. The learned model parameters may not be appropriate for another society. The virulence could be globally invariant, but it could also depend on genetic and demographic factors like age, as well as what "serious" means for the culture -that is when do you get yourself checked. In a high-strung society, you expect high virulence since the threshold for serious is lower. One certainly expects the infectious force to depend on the underlying social network and interaction patterns between people, which can vary drastically from one society to another and depending on interventions. Hence, one should calibrate the model to country specific data to gain more useful insights.
Discussion. Early dynamics allows us to learn useful properties of the pandemic. Later dynamics may be contaminated by human intervention, which renders the data not as useful without a more general model. Are our parameters correct? We are in a unique position to test our predictions because our model was time-stamped as version 1 of the preprint Magdon-Ismail (2020). The model has provably not changed, and we just added in new test data as it arrives (Figure 3) . The model predictions in Figure 3 are in no way forward looking, data snooped or overfitted. Training used only data before 03/14. The model is holding up.
The Lag, k, and Public Policy. The lag k is important for public policy. Our lag, k = 10, says the confirmed infections should fall away from the predictions on about 03/22, because agressive social distancing began on around 03/12. The lag is important due to how public policy can be driven by human psychology as opposed to quantitative models. March 20 saw a 40% increase in new infections compared to the previous day. This led to a significant heightening of vigilance and calls for "shelter in place" protocols. But, by March 20, we still do not know the effectiveness of social distancing protocols from March 12th have been. Hence, such calls to vigilance can be premature. More importantly, suppose strict shelter in place protocols are enacted on March 21 and infections fall the next day. The Human tendency is to associate the improvement to the March 21 strict protocols. But, if the model is correct, the March 21 protocols have nothing to do with the improvement on March 22. Such lags are present in traditional machine learning, for example the delayed reward in reinforcement learning settings. Credit assignment to prior actions in the face of delayed reward is a notoriously difficult task, and this remains so with humans in the loop. Knowledge of the lag helps to assign credit appropriately to prior actions, and the public health setting is no exception.
We perform our analysis on early dynamics data available from the ECDC giving infection numbers starting from December 31, 2019 ECDC (2020). We use data for 20 countries selected qualitatively because they appear to have reasonably efficient testing procedures for self-reported cases. We iniclude China data for completeness, even though China dynamics since December 31 are not early dynamics. We show these countries below, together with some demographic data which might determine spread dynamics: Start City; Start and End dates; Delay to first infection in days; Population Density; Median Age ; Wealth as defined by adjusted net national income per capita (World-Bank, 2017); Average Income (World-Data-Info, 2015); Poverty Level (Wikipedia, Source: World Bank, 2020) defined by living on less than $5.50 per day; Population Density around the first infection site for the country.
Recall the model,
For fixed k, M 0 , we must perform a gradient descent to optimize β, γ. Unfortunately, the dependence on β is exponential and hence very sensitive. So if the starting point is not chosen carefully, the optimization gets stuck in a very flat region, and many millions of iterations are needed to converge. Hence it is prudent to choose the starting conditions carefully. To do so, we need to analyze the recursion. First, we observe that the recursion for M (t) is a standalone k-th order recurrence. Hence, we can guess a solution M (t) = M 0 φ t−1 , for t ≥ 1, which requires
We do a perturbation analysis in γ → 0. At γ = 0, φ = β, so we set φ = β + , to get (β + ) k − β(β + ) k−1 + γ = 0, which to first order in is solved by ≈ −γ/β k−1 and so
where φ = β(1 − γ/β k ). Given this approximation, we can solve for S(t). For t > k,
We can independently control two parameters φ and γ. We use this to match the observed S(t) at two time points. Since the growth is exponential, we match the end time, S(T ) and some time τ in the middle, for example τ = 3T /4 . Let ∆ T = (S(T ) − S(1))/M 0 and ∆ τ = (S(τ ) − S(1))/M 0 . Then,
Dividing gives ∆ T /∆ τ = (φ T −k − 1)/(φ τ −k − 1) ≈ φ T −τ , because φ > 1. Let us consider the equation κ = (φ r − 1)/(φ s − 1), which gives φ r − κφ s + κ − 1 = 0, or more generally φ r − κφ s + ρ = 0, where r > s > 1 and κ > ρ 1. This means φ > 1. When ρ = 0, we have φ r−s = κ, so we do a perturbation analysis with φ r−s = κ + , and our perturbation parameter is . Then, φ r = (κ + )φ s and plugging into the equation gives
Solving for gives ≈ −ρ(r − s)/((r − s)κ s/(r−s) + s), which gives
For our setting, r = T − k, s = τ − k, κ = ∆ T /∆ τ and ρ = κ − 1. Finally, since φ is approximate, we may not be able to satisfy both equations ineq:app:1, hence we can instead minimize the mean squared error, which gives
We now need to get β which satisfies φ = β(1 − γ/β k ). Again, we do a perturbation analysis, omitting the details, to obtain
If one wishes, a fixed point iteration starting at the above will quickly approach a solution to φ = β(1 − γ/β k ). We show the approximate fit on the US data ( Figure 4) . We show the optimal fit, the initial fit using the parameters constructed fromeq:app:bet andeq:app:gam. The parameters and fit error are β γ fit error optimal 1.306 0.0013282 1.6521 Equationseq:app:bet andeq:app:gam 1.3055 0.0013423 1.6619 The approximate fit works pretty well, and is certainly good enough to initialize an optimization.
From the public health perspective, perhaps the most important parameter is β, since actions can be taken to mitigate the spread by reducing β, whereas γ, k and M 0 are somewhat givens for the country. We show the fits in Table 3 . As you can see, there is much variability in β.
We perform a simple statistical analysis to test if β can be explained by any of the country parameters in Table 2 . We include the delay as a global explanatory variable, which would account for a global increase in vigilence as time passes and awareness of the pandemic increases. One expects β to decrease with the delay. A table of correlations of β with the various parameters is shown below. For our analysis we use the best case β, although similar results follow from the optimal β. As expected, there is a very significant correlation of β with delay, but in the opposite direction.
• The larger the delay, the larger is β. The more a country has observed, the faster the spread in that country. That seems unusual but seems strongly indicated by the data. • Population density at the infection site has a strong positive effect but the country's population density does not.
• There is faster spread in poorer countries.
• Median age has a strong effect. Spread is faster in younger countries. The youth are more mobile and perhaps also more carefree.
• There is a slight negative effect from wealth and per-capita income. Spread is slower in richer countries. Perhaps this is due to more risk-aversion, perhaps higher levels of education, perhaps less use of public transportation. Whatever the cause, it does have an impact, but relatively smaller than the other effects.
We now use regularized regression to perform a linear model fit to explain β. To make the weight magnitudes meaningful, we normalize the data. We use a leave-one-out cross validation to select the optimal regularization parameter (which happens to be 20). The optimal regularized fit with this regularization parameter gives a new feature X = w 1 ·(Delay)+w 2 ·(Pop-Den)+w 3 ·(Age)+w 4 ·(Wealth)+w 5 ·(Income)+w 6 ·(Poverty)+w 7 ·(Pop-Den-Init)
The learned weights and the their ranges which yield a cross-validation error within 10% of optimal are shown in the table. Figure 5 : Optimally constructed feature to predict β, within a cross-valudation setting to select the regularization parameter. The gray region is the range of the predicted value for each country. The R 2 = 0.57, so the cross-validation based optimal linear feature captures 57% of the residual variance. The predictions of β using this feature are also shown in Figure 5 . As we observed from the correlations, Delay, Poverty and Population Density at the initial infection site have strong positive weights. Age has a strong negative weight. Wealth and income have weak negative effects, but non-zero. The population density of the country as a whole seems to have no effect, with a weight range that includes 0.