Quantitative assessment of the risk of airborne transmission of SARS-CoV-2 infection: Prospective and retrospective applications

Graphical abstract


Introduction
The airborne transmission of a virus and the consequent contagion risk assessment is a complex issue that requires multidisciplinary knowledge. It is necessary to understand the characteristics and mechanisms behind the generation of respiratory microdroplets (Ai and Melikov, 2018;Holmgren et al., 2010), the survival of viruses in microdroplets (van Doremalen et al., 2020), the transport of microdroplets and human exposure to them (Ai et al., 2019a), and the airflow patterns that carry microdroplets in buildings (Ai et al., 2019b). Expiratory human activities generate virus-carrying microdroplets that are small enough to remain aloft in air during exhalation, talking, and coughing (Holmgren et al., 2010;Morawska et al., 2009;Morawska and Cao, 2020). Atomization occurs in the respiratory tract, and droplets are expelled at high speed during expiration Morawska, 2006). Toques of liquid originating from different areas of the upper respiratory tract are drawn out from the surface and broken into droplets of different sizes (Hickey and Mansour, 2019). The findings of early investigations (Duguid, 1945;Jennison, 1942;Wells, 1934) served as a foundation for subsequent studies involving temporal and spatial visualization methods using high-speed cameras (Tang et al., 2011), particle image velocimetry ) and, above all, increasingly accurate particle counters , which have facilitated the detailed characterization and quantitation of droplets expelled during various forms of human respiratory exhalation flows.
The issue of the viral load emitted, however, remained difficult to solve. In the past, backward calculation was used to estimate the emission of an infected subject based on retrospective assessments of infectious outbreaks only at the end of an epidemic (Myatt et al., 2008;Rudnick and Milton, 2003;Sze To and Chao, 2010;Wagner et al., 2009). This led to the definition of emission values for each virus regardless of the type of respiratory act and the metabolic activity of the infected subject. Recently, the authors presented an approach to evaluate the viral load emitted by infected individuals with a view to provide new predictive capacities, not currently available . This approach, based on the oral viral load and the infectivity of the virus, takes into account the effect of other parameters such as inhalation rate, type of respiratory activity, and activity level, to estimate the quanta emission rate. This value provides key information for engineers and indoor air quality experts to simulate airborne dispersion of diseases in indoor environments. Indeed, the use of exposure risk models in closed environments (Gammaitoni and Nucci, 1997;Riley et al., 1978) makes it possible to estimate contagion starting from the emission values of a contagious subject.
The overall approach of emission and exposure modelling represents an essential tool to be applied in enclosed spaces, and can support air quality experts and epidemiologists in the management of indoor environments during an epidemic for both prospective and retrospective assessments.
In this paper we apply a novel approach that takes into account the characteristics of the emitting subject, the microenvironment, and the exposed subject to calculate the probability of infection and the individual risk, for both prospective and retrospective assessments of airborne infectious transmission of SARS-CoV-2. In the case of prospective assessment, various exposure scenarios in indoor environments were analyzed in order to assess the influence of risk mitigation parameters. In the case of retrospective assessment, we estimated the probability of infection of two documented outbreaks.

Materials and methods
To evaluate both prospective and retrospective assessments of the airborne transmission of SARS-CoV-2, we used a four-step approach to quantify the probability of infection due to exposure in a microenvironment where a SARS-CoV-2 infected subject is present. The four steps of the proposed approach are: (i) evaluation of the quanta emission rate; (ii) evaluation of the exposure to quanta concentration in the microenvironment; (iii) evaluation of the dose of quanta received by an exposed susceptible subject; and (iv) estimation of the probability of infection on the basis of a dose-response model. The simulations of the probability of airborne transmission of SARS-CoV-2 were performed by applying a Monte Carlo method (Hammersley and Handscomb, 1964) and adopting the infection risk assessment typically implemented to evaluate the transmission dynamics of infectious diseases and to predict the risk of these diseases to the public (Gammaitoni and Nucci, 1997;Riley et al., 1978;Sze To and Chao, 2010).
Once the probability of infection for a given exposure conditions (i. e., as hereinafter explained, for a given quanta emission rate) was obtained, the individual infection risk was calculated, taking into account how likely that condition (i.e. quanta emission rate) can occur. Individual risk can be easily compared to an acceptable risk, i.e. a target reference risk that could be suggested by agencies and regulatory authorities to control the pandemic. In the following sections, the methodologies adopted to evaluate the probability of infection based on the four step approach (Section 2.1) and the individual infection risk (Section 2.2) are described. The application of the proposed approach for prospective and retrospective assessments is described in Sections 2.3 and 2.4.

Evaluation of the quanta emission rate: The forward emission approach
Recently,  proposed a forward emission approach to estimate the quanta emission rate of an infectious subject on the basis of the viral load in the sputum and the concentration of droplets expired during different activities. A quantum is defined as the dose of airborne droplet nuclei required to cause infection in 63% of susceptible persons. The quanta emission rate (ER q , quanta h − 1 ) is evaluated as: where c v is the viral load in the sputum (RNA copies mL − 1 ), c i (quanta RNA copies − 1 ) is a conversion factor defined as the ratio between one infectious quantum and the infectious dose expressed in viral RNA copies, IR is the inhalation rate (m 3 h − 1 ), and V d is the droplet volume concentration expelled by the infectious person (mL m − 3 ). The droplet volume concentration V d is a function of the expiratory activities (i.e. breathing, speaking, singing, etc.). Experimental data on the droplet volume emitted are not definitive and the sampling method itself can affect the results due to the rapid dehydration occurring to the large particles emitted (Abbas and Pittet, 2020;Stadnytskyi et al., 2020;Yang and Marr, 2011). In the present paper we derived the droplet volume concentration on the basis of the total volume emitted by a loudspeaking person provided by (Stadnytskyi et al., 2020) as they measured the droplet volume taking into account the droplet dehydration. In particular, they measured an average droplet volume emission rate of 3 × 10 − 2 mL h − 1 (ranging roughly from 1 × 10 − 2 to 5 × 10 − 2 mL h − 1 ). Thus, we have calculated the droplet concentration expelled by the infectious person dividing the volume emission rate by the expiration rate (considered equal to inhalation rate) for standing activity level (Adams, 1993)) obtaining an average droplet concentration of 6 × 10 − 2 mL m − 3 (ranging from 2 × 10 − 2 mL m − 3 to 9 × 10 − 2 mL m − 3 ). Nonetheless, (Stadnytskyi et al., 2020) did not evaluate the emission for other expiratory activities; thus, we evaluated the droplet emission rate for other activity levels by using the inhalation rates provided by (Adams, 1993) (and hereinafter summarized) and for other respiratory activities by adopting the relative volume concentrations of oral breathing (0.03fold) and voiced counting (0.16-fold) with respect to the loud speaking (namely "unmodulated vocalization") provided by .
With reference to the SARS-CoV-2 viral load in the mouth, researchers have recently found c v values in the range 10 3 -10 11 copies mL − 1 , for both symptomatic and asymptomatic persons, which is also variable in the same patient during the course of the disease (Lavezzo et al., 2020;Pan et al., 2020;Rothe et al., 2020;To et al., 2020;Woelfel et al., 2020). In the present paper we have considered an average c v value of 10 7 copies mL − 1 .
The conversion factor, c i (quanta RNA copies − 1 ), i.e. the ratio between one infectious quantum and the infectious dose expressed in viral RNA copies, represents the probability of a pathogen surviving inside the host to initiate the infection; thus c i = 1 implicitly assumes that infection will occur for each pathogen (RNA copy in the case of SARS-CoV-2) received by the exposed people. To convert the infectious dose expressed in viral RNA copies to infectious quanta two parameters should be known: a) the number of infectious virus particles (RNA copies) needed to initiate the infection (c RNA , RNA copies PFU − 1 ), and b) the quanta-to-plaque forming unit (PFU) conversion parameter (c PFU , PFU quanta − 1 ). Then, as expressed in Eq. (1), the conversion factor c i is the product of c RNA − 1 and c PFU − 1 . Actually, there are currently no values available in the scientific literature of c PFU for SARS-CoV-2. (Watanabe et al., 2010) estimated the infectious doses of several coronaviruses on the basis of data sets challenging humans with virus HCoV-229E (known as an agent of human common cold) and animals with SARS-CoV-1. On the basis of those data, we adopted an average c PFU value of 2.1 × 10 2 PFU quanta − 1 . As regards the c RNA parameter, recently an estimate for SARS-CoV-2 virus was provided by (Fears et al., 2020) in an experimental characterization of persistence of SARS-CoV-2 in aerosol suspensions: from their data a c RNA = 1.3 × 10 2 RNA copies PFU − 1 was adopted.
The quanta emission rate calculation was performed for four different emission profiles (which are adopted in the risk evaluations described later) evaluated as a combination of expiratory activities and activity levels: (i) oral breathing during resting; (ii) oral breathing during heavy activity; (ii) speaking during light activity; and (iv) singing (or loudly speaking) during light activity.
Quanta emission rates were calculated using Eq. (1) and applying a Monte Carlo method (Hammersley and Handscomb, 1964) in order to take into account for the possible variation of the input data. To this end, probability density functions characteristics of each parameter were considered. In particular, we considered normal distributions for: (i) the log-transformed c v data (average and standard deviation of log 10 (c v ) equal to 7 and 0.71 log 10 (RNA copies mL − 1 ), respectively); (ii) the droplet volume concentration expelled by the infectious person, V d , (averages and standard deviations of 6 × 10 − 2 mL m − 3 and 2 × 10 − 2 mL m − 3 for unmodulated vocalization, 9 × 10 − 3 mL m − 3 and 2 × 10 − 3 mL m − 3 for voiced counting, and 2 × 10 − 3 mL m − 3 and 5 × 10 − 3 mL m − 4 for oral breathing, respectively), (ii) the c PFU parameter (average and standard deviation equal to 2.1 × 10 2 PFU quanta − 1 and 2.1 × 10 1 PFU quanta − 1 , respectively), and (iv) the c RNA parameter (average and standard deviation equal to 1.3 × 10 2 RNA copies PFU − 1 and 1.3 × 10 1 RNA copies PFU − 1 , respectively). A distribution of quanta emission rates (ER q ), was obtained as a result of application of the Monte Carlo method (Table 2), i.e. the probability density function of ER q (pdf ERq ).

Evaluation of the exposure to quanta concentration
The second step in evaluating the probability of infection is the calculation of the quanta concentration to which a susceptible subject is exposed. The quanta concentration at time t, n(t), in an indoor environment is based on the quanta mass balance proposed by (Gammaitoni and Nucci, 1997), and can be evaluated as: where IVRR (h − 1 ) represents the infectious virus removal rate in the space investigated, n 0 represents the initial quanta concentration (i.e. at time t = 0), I is the number of infectious subjects, V is the volume of the indoor environment considered, and ER q is the quanta emission rate (quanta h − 1 ) for the specific disease/virus under investigation. The infectious virus removal rate is the sum of three contributions (Yang and Marr, 2011): the air exchange rate (AER) via ventilation, the particle deposition on surfaces (k, e.g. via gravitational settling), and the viral inactivation (λ). The deposition rate was evaluated as the ratio between the settling velocity of super-micrometric particles (roughly 1.0 × 10 − 4 m s − 1 as reported by (Chatoutsidou and Lazaridis, 2019)) and the height of the emission source (1.5 m); thus, k was 0.24 h − 1 . We point out that this is a simplification of the real conditions, in which the droplets emitted span over a wide particle size range, and hence present a wide range of settling velocities. Nonetheless, in basic approaches as the boxmodel here adopted, choosing different particle deposition values for different particle sizes would unnecessarily complicate the calculation; on the contrary, this could be important in more complex approaches, such as CFD ones, allowing the calculation of the path of each droplet. The viral inactivation was evaluated on the basis of the SARS-CoV-2 half-life (1.1 h) detected by (van Doremalen et al., 2020), thus λ was 0.63 h − 1 . The quanta concentration calculation adopted here is based on the following hypotheses: the infectious virus removal rate is constant, the latent period of the disease is longer than the time scale of the model, and the droplets are instantaneously and evenly distributed in the room (Gammaitoni and Nucci, 1997). Thus, since all the other parameters affecting the exposure were considered as constant values, the quanta concentration is a function of t and ER q (n(t, ER q )) and, for each quanta emission rate ER q , a different quanta concentration over time in the environment can be calculated.
In the exposure scenarios tested with the prospective and retrospective approaches, to take the variability of the input parameters into account, the indoor quanta concentration n(t,ER q ) was determined through Eq. (2) by applying a Monte Carlo method that adopted the probability density functions characteristic of quanta emission rates (pdf ERq ). Since the probability density functions of the log-transformed log 10 (ER q ) for the different expiratory activities resulted in a normaldistribution (Shapiro-Wilk test, p < 0.01), the quanta concentration n(t, ER q ) was evaluated by providing a Gaussian distribution of log 10 (ER q ) (average and standard deviation values are summarized in the results section; see Table 2) and then applying a back-transformation from log 10 (ER q ) to ER q . The relative frequency at which a certain quanta concentration occurred for each time step of simulation, i.e. the probability density function of the quanta concentration (pdf n ), was also obtained as result of the Monte Carlo simulations.

Evaluation of the dose of quanta received by an exposed susceptible subject
The dose of quanta received by a susceptible subject exposed to a certain quanta concentration, n(t,ER q ), for a certain time interval, T, can be evaluated by integrating the quanta concentration over time as: It can be concluded from Eq. (3) that the dose of quanta received by a susceptible subject is affected by the inhalation rate (IR) and subsequently by their activity level. As an example, for the same exposure scenario [i.e. identical n(t,ER q ) and T], the dose of quanta received by subjects performing at a light activity level (IR = 1.38 m 3 h − 1 ; e.g. slowly walking) is more than double that received by people just sitting or standing (IR = 0.54 m 3 h − 1 ). For the dose calculation, in the exposure scenarios described in this paper, the IR of the exposed subject was considered as a constant value; thus, once again, the dose of quanta received is a function of the quanta emission rate ER q . In particular, the different dose values for each ER q were obtained applying the Monte Carlo method, i.e. adopting the probability density functions of the dose (pdf D ) resulting from the pdf ERq .

Evaluation of the probability of infection through a dose-response model
The fourth and final step in evaluating the probability of infection is the adoption of a dose-response model. Several dose-response models are available in the scientific literature for assessing the probability of infection of airborne-transmissible pathogens (Rudnick and Milton, 2003;Sze To and Chao, 2010), including deterministic and stochastic models, and threshold and non-threshold models.
The best-suited dose-response models for airborne transmission of pathogens are the stochastic models (Sze To and Chao, 2010). In particular, exponential models have been mostly adopted in previous studies because of their suitability and simplicity (Watanabe et al., 2010). Such models consider the pathogens as discrete bundles (i.e. quanta) distributed in a medium (e.g. saliva/sputum) in a random manner described by the Poisson probability distribution. When the medium is aerosolized, the pathogen distribution in the aerosols, and hence their distribution in the air, also follows the Poisson probability distribution. The complex Poisson summation equations can be simplified in an exponential equation (Haas, 1983;Sze To and Chao, 2010;Watanabe et al., 2010), i.e. the exponential dose-response model, which evaluates the probability of infection, P I (%), of susceptible people as: For a unit dose of quanta (D q = 1), the probability of infection P I is equal to 63%, from which derives the definition of "quantum" as the "amount of infectious material to infect 1-e − 1 (i.e. 63%) of the people in an enclosed space" (Gammaitoni and Nucci, 1997;Wells, 1934).
In the exponential dose-response model, the variation of host sensitivity to the pathogen is not considered. More complex models, such as the Beta-Poisson probability distribution, could take this factor into account (Haas, 1983;Sze To and Chao, 2010;Watanabe et al., 2010); nonetheless, in the present paper the differences in the exposed population in terms of susceptibility to the virus will not be considered.
In our scenarios, since all the other parameters affecting the exposure and dose were considered constant values, the probability of infection is a function of the quanta emission rate, P I (ER q ). In particular, we evaluated the P I (ER q ) of each exposure scenario through Eq. (4), also applying a Monte Carlo method. To this end, the probability density functions of the dose of quanta (pdf D ) obtained as a result of the Monte Carlo simulation on D q (ER q ) were considered; thus, the probability density functions of P I (ER q ) were also obtained (pdf P ).

The individual infection risk and the basic reproduction number
As stated above, for a given exposure scenario (microenvironment, ventilation, inhalation rate of the exposed subject, etc.) the probability of infection can assume different values on the basis of the rate of quanta emitted by the infected subject, P I (ER q ): the lower the quanta emission rates, the lower the probability of infection (since all the other parameters affecting the exposure and the dose were considered to be constant values). Thus, in order to evaluate the individual infection risk (R) of an exposed person for a given exposure scenario, we should know both the probability of infection (P I (ER q )) and the probability of occurrence of each ER q value (P ERq ) which is defined by the abovementioned ER q probability density functions (pdf ERq ). Since the probability of infection (P I (ER q )) and the probability of occurrence P ERq are independent events, the individual infection risk for a given ER q , R(ER q ), can be evaluated as the product of the two terms: where P I (ER q ) is the conditional probability of the infection, given a certain ER q , and P ERq represents the relative frequency of the specific ER q value. Applying the Monte Carlo simulation to the Eq. (5) we obtained the probability density function of R (pdf R ), i.e. the R(ER q ) values. The individual risk (R) of an exposed person, was then evaluated integrating the pdf R for all possible ER q values, i.e. summing up the R(ER q ) values calculated in Eq. (5): The individual risk R also represents the ratio between the number of infection cases (C) and the number of exposed susceptible individuals (S) for a given exposure scenario and taking into account all possible ER q values for the infectious subject under investigation. In retrospective analyses of documented outbreaks, the known C/S ratio is typically defined as the "attack rate".
We point out that R represents the individual risk of being infected for a healthy subject when exposed to a non-zero quanta concentration, i.e. when an infected person was (simultaneously or before) present in the same indoor environment. Indeed, the present model assumes the presence of a contagious individual: thus, to perform a large-scale risk assessment, the probability of running into a contagious person should be considered.
In a conservative application of the proposed approach to estimate and reduce the risk of individuals being together with an infected individual in an indoor environment, the individual infection risk must be less than an acceptable risk.
The US Environmental Protection Agency (EPA) typically uses a target reference risk range of 10 − 4 to 10 − 6 for carcinogens in drinking water (Cotruvo, 1988), which is in line with World Health Organization (WHO) guidelines for drinking water quality, which base guideline values for genotoxic carcinogens on the upper bound estimate of an excess lifetime cancer risk of 10 − 5 (World Health Organization, 2011). If the estimated lifetime cancer risk is lower than 10 − 6 , the risk is considered acceptable, while risks above 10 − 4 are considered unacceptable (Toner, 2008).
The choice of an acceptable contagion risk for SARS-CoV-2 is difficult and certainly questionable. However, considering the mortality rate of SARS-CoV-2, this turns out to be an order of magnitude lower than the corresponding value associated with carcinogenic diseases. For this reason, only for discussion purposes, the value of 10 − 3 is taken as an acceptable risk reference for SARS-CoV-2.
For the purpose of managing an epidemic and keeping the infection under control, it is also important to estimate the basic reproduction number of the infection, R 0 , which is calculated as the ratio between the number of susceptible people infected (C) and the infected subject (I). Thus, R 0 can be easily evaluated by multiplying the individual infection risk, R, by the number of exposed susceptible individuals (S). To control an epidemic, the R 0 value must be <1. Therefore, in addition to estimating an acceptable individual infection risk, it is necessary to specifically verify that, with the crowding expected within the environment, the corresponding value of R 0 is <1. In view of this, adopting an acceptable contagion risk (R max ) of 10 − 3 , would results in a R 0 ≤1 when up to 1000 persons gather in the same closed environment. Thus, the 10 − 3 acceptable risk can be considered a quite effective value in view of containing the pandemic since very few scenarios are characterized by the co-presence of more than 1000 people in a closed environment (some indoor sporting events); moreover, for such large closed environments the hypotheses here adopted (e.g. droplets evenly distributed in the room) are no longer valid and more complex models should be implemented.

Scenarios in the prospective assessment
The proposed four-step approach was applied to different indoor microenvironments by varying the main parameters in order to evaluate the effect of the influencing parameters. In particular, four emission profiles of the infected subject and corresponding profiles of the healthy subjects exposed were chosen. For the sake of simplicity, the simulations were run assuming that the susceptible subjects remained in the microenvironment for the same length of time as the infected subject (i. e. the two subjects enter and leave the environment under test together). Each indoor environment under investigation was tested for three different values of air exchange rate (AER). Table 1 presents a detailed summary of the four different indoor exposure scenarios considered to evaluate the risk of airborne transmission of SARS-CoV-2. Scenario A consists of a hospital room of 100 m 3 where a resting infected patient emits quanta in the room through oral breathing, whereas the exposed susceptible subjects consist of a member of the medical staff (not wearing individual protection devices) in a light exercise activity (scenario A-1) and another patient at rest (scenario A-2). No exceptional events as coughing and sneezing were considered in the evaluation of the quanta emission rate of the infectious person. In scenario B, the infection affects two subjects, both oral breathing during a sports activity in a 300 m 3 gym. Scenario C concerns two subjects (infected and healthy) in light activity while speaking in a generic 300 m 3 office (bank, post office, supermarket, shop, etc.). Finally, scenario D represents an infected subject singing or speaking loudly in an 800 m 3 room with healthy subjects listening at a sedentary activity level. In Table 1 the inhalation rates for both emitting subjects and exposed subjects, as obtained from (Adams, 1993), were provided. Once again, we highlight that the simulations were performed adopting a box-model approach, i. e. considering fully mixing conditions and constant removal rates (e.g. deposition rate not influenced by the droplet size); we are aware that such hypotheses may lead to a simplification of the real conditions for some of the investigated scenarios; nonetheless, they are very helpful in providing a first attempt calculation with respect to more complex approaches (e.g. CFD codes). Indeed, the latter would need exact information in terms of geometries, air ducts positioning, source location within the room, then providing results specific of the simulated environment.

The outbreak in a restaurant in Guangzhou, China
A possible case of airborne transmission was recently documented by (Lu et al., 2020). Here, an index case patient traveled from the Chinese epidemic epicenter, Wuhan, on 23 January 2020 and ate lunch in a restaurant in Guangzhou, China, with his family on 24 January 2020 (family A, 10 people sitting at the same table). Later that day, the index patient experienced onset of fever and cough and SARS-CoV-2 infection was diagnosed. On the following days, nine other people were diagnosed with SARS-CoV-2 infection: four members from family A's table and five other people at two different tables (families B and C). No other customers seated at other tables or waiters were infected.
The restaurant is a 5-floor building without windows; each floor has its own air ventilation system. The third floor dining area, at which the index patient ate lunch, has a floor area of 145 m 2 , with 15 tables arranged with a distance between each table of about 1 m. A total of 91 people (83 customers, 8 staff members) were in the room during that lunch. The exposure time was variable for the customers: those seated at tables close to the index patient had exposure times of 53 min (family B) and 73 min (family C). The ventilation and air conditioning situation is reported (Lu et al., 2020). Five fan coil air-conditioning units are installed in the room and there is no outdoor air supply; thus, the ventilation relies only upon infiltration and natural ventilation. The authors performed computational fluid dynamics analyses and tracer gas decay tests to obtain more information about the possible air-flow pathway in the room, and to determine the air exchange rate expected during that lunch. The analyses performed showed that, due to the particular installation and use of the fan coils, the room can be divided into different air-flow zones, with well-mixed conditions. The air-flow zone involving the table at which the index patient sat also included the two tables at which the other five infected people sat; and covered an area of roughly 45 m 3 . The tracer gas decay tests revealed a low air exchange rate (mostly due to the absence of an outdoor air supply) in the range of 0.56-0.77 h − 1 . Therefore, on the basis of the available information, the retrospective assessment was applied to this outbreak case, through Eqs. (2) and (3), using the following input data: (i) room volume of 45 m 3 ; (ii) documented probability of infection, i.e. attack rate, of 45% (i.e. 5 out of 11 people of families B and C (family A members were excluded as they could easily have been infected through other infection routes); (iii) average exposure time of 1 h; (iv) speaking at a light activity level for all people (both emitting and exposed subjects), and (v) average AER = 0.67 h − 1 .

The outbreak at a choir rehearsal in Skagit (USA)
A further possible case of airborne transmission of SARS-CoV-2 was documented by the USA media (www.latimes.com/world-nation/story/ 2020-03-29/coronavirus-choir-outbreak) (Miller et al., 2020). This case was recorded on 10 March, in Mount Vernon (Skagit County, Washington State, USA). In a 810 m 3 hall, 61 choir members (out of a total of 121 regular members) gathered to rehearse, aware of the practices for the containment of contagion (frequent hand washing and social distancing). None of the members that attended had evident symptoms of SARS-CoV-2 infection. There was hand sanitizer at the front door and members refrained from the usual hugs and handshakes; each person brought their own sheet music. The event lasted from 6:30 pm to 9:00 pm (about 2.5 h). Within few days, 33 of the 61 participants (53%) were diagnosed with SARS-CoV-2 infection, at least three were hospitalized, and two died (Hamner et al., 2020).
As pointed out by (Hamner et al., 2020), the 2.5-hour singing Table 1 Description of the exposure scenarios tested in the prospective assessment.  practice could have provided several opportunities for droplet and fomite transmission (e.g. members sitting close to one another, sharing snacks, and stacking chairs at the end of the practice). Nonetheless, the abovementioned voluntary measures put in place would not support the documented spread of the contagion. On the contrary, the act of singing, itself, might have contributed to transmission through emission of aerosols, which is affected by loudness of vocalization . This is even more relevant considering that attack rate of 53.3% (based on 33 confirmed cases) could represent a conservative estimate, since another 20 probable cases were mentioned by (Hamner et al., 2020). As regard the heating and ventilating system, limited information is available: the Fellowship Hall is heated by a relatively new commercial forced-air furnace with supply and return air grills situated high on a single wall, but it is not known how much fresh air was provided on that evening. During the entire rehearsal no exterior doors were open. We applied a retrospective assessment to the case of the Skagit Valley choir through Eqs. (2) and (3), using the following input data: (i) room volume of 810 m 3 ; (ii) documented probability of infection, i.e. attack rate, equal to 53%; (iii) exposure time of 2.5 h; (iv) singing at a light activity level for all people; and (v) natural ventilation with an AER = 0.5 h − 1 . Table 2 shows the statistics of the quanta emission rates for the four emission profiles considered in Section 2.1. As shown in , there are large differences between the emission profiles. Obviously the lowest values are found under the oral breathing condition during resting (median value of 0.37 quanta h − 1 ), followed by the oral breathing condition during heavy activity as the inhalation rate increases (2.5 quanta h − 1 ), and reaching 5.0 quanta h − 1 for the increase in aerosol emitted during vocalization ) and, finally, peaking during singing/speaking loudly (32 quanta h − 1 ). Indeed, the rate of particle emission during normal human speech is positively correlated with the amplitude of vocalization (Asadi et al., 2019).

Statistics of quanta emission rates
The probability density functions of the quanta emission rates (pdf ERq ) were also determined. In particular, the log-transformed ER q values obtained from the Monte Carlo simulations resulted in a normal distribution (Shapiro-Wilk test, p < 0.01). Table 2 shows the average and standard deviation values of the log 10 (ER q ).
We point out that the estimated values present two main uncertainty contributions clearly related to the limited data currently available for the SARS-CoV-2: (i) a still low number of experimental data for the viral load in the mouth, c v , of SARS-CoV-2 infected subjects, (ii) limited information on the infectivity conversion factors for SARS-CoV-2.

Illustrative example of probability of infection and individual risk evaluation
In Fig. 1 an illustrative example of quanta concentration n(t,ER q ), dose of quanta (D q (ER q )), and probability of infection (P I (ER q ) trends as a function of time (here shown for 2 h) and quanta emission rate, ER q , resulting from the Monte Carlo simulation for exposure scenario D (singing exhibition, conference speaker) with an AER = 0.5 h − 1 is shown. In particular, the trends of n(t,ER q ), D q (ER q ), and P I (ER q ) for different percentiles of ER q are reported. The example shows that a person singing/speaking loudly in such a microenvironment and generating a quanta emission rate at 50th percentile (32 quanta h − 1 as reported in Table 2) lead to a n(t,ER q ) value after 2 h equal to 0.019 quanta m − 3 . Such concentrations lead to a corresponding dose of quanta received by the subject exposed for 2 h in a sedentary activity equal to 0.021 quanta, then resulting in a median probability of infection, P I (ER q ), of 2.1%. The figures clearly show that for higher quanta emission rate, the indoor quanta concentrations and the consequent probability of infection can be more than 10-fold than those evaluated for the median emission rate: as an example, for a quanta emission rate at 95th percentile (4.9 × 10 2 quanta h − 1 ), the P I (ER q ) is equal to 27%.
In view of the application of a conservative approach that could be essential to reduce the risk of contagion in indoor environments, using the highest quanta concentration and probability of infection values can be misguiding. Indeed, the probability of occurrence of the ER q leading to such high quanta concentration and probability of infection values is extremely low. Thus, as described in Section 2.2, a proper evaluation of the individual infection risk (R) can be obtained by applying Eq. (5), i.e. multiplying the probability of infection probability of infection (P I (ER q )) by the probability of occurrence of ER q (P ERq ). In Fig. 2 the probability density functions of ER q (pdf ERq ), probability of infection (pdf P ), and individual infection risk (pdf R ) after 2 h of exposure are reported for the illustrative example discussed above (i.e. scenario D, AER = 0.5 h − 1 ). The relative frequency P ERq is here graphed as 99 equally spaced log (ER q ) values (from 1st to 99th percentiles) adopting the log-normal probability density functions of ER q (pdf ERq ) evaluated as described in Section 2.1.1 (and reported in Table 2 as average and standard deviations).
The individual infection risk (R(ER q )) presents a maximum value at the 92nd percentile of ER q of a singer/speaker (i.e. about 330 quanta h − 1 ). The individual risk R for such scenario, i.e. the sum of the R (ER q ) values, resulted equal to 3.8%; whereas, when the AER is increased to 3 h − 1 and 10 h − 1 , the risk decreased to 1.9% and 0.7%, respectively. Due to the similarity of the probability density functions of the four expiration activities resulting from the calculation of the quanta emission rates (log 10 (ER q ) reported in Table 2), the pdf R for all the exposure scenarios tested here were similar to that of the exposure scenario shown in Fig. 2 (i.e. the maximum R(ER q ) values occur in the narrow range of 90th-95th percentile).
Furthermore, as discussed in Section 2.2, the probability density function of the probability of infection (pdf P ) is mostly influenced by the probability density function of the quanta emission rate (pdf ERq ) when moving backwards in the four-step approach; indeed, once the exposure scenario is defined, all the parameters contributing to the calculation of P I (ER q ) (ventilation, room volume, subject activity, etc.) can be considered as constant values. Thus, a simplified estimate of R (instead of the Monte Carlo method) can be carried out. Indeed, using the Eq. (6), the ER q value assumed with a probability of occurrence P ERq = 1 (i.e. considered as a certain emission) which induces a P I (ER q ) equal to the risk R evaluated through the Monte Carlo method can be calculated. Such ER q can be adopted in the four-step calculation using Eqs. (2)-(4) to assess the risk of different scenarios even not running a Monte Carlo simulation. The certain emission value evaluated for the investigated scenarios is the 66th percentile ER q value, which is can be easily calculated on the basis of the averages and standard deviations of the log 10 (ER q ) provided in Table 2. As an example, the 66th percentile ER q values for oral breathing during resting, oral breathing during heavy activity, speaking during light activity, and singing (or loudly speaking) during light activity are equal to 0.72, 4.9, 9.7, and 62 quanta h − 1 , respectively.

Estimate of the individual risk versus exposure time in indoor environments
In Fig. 3 and Table 3 the exposure time-risk relationships are reported for the four exposure scenarios analyzed. Such relationships are essential as they can be used by choosing either the exposure time or the risk as the independent variable. In the first case, knowing the exposure time of the healthy subject in the environment in question, the corresponding individual infection risk can be evaluated and then compared to an acceptable infection risk. In the second case, once an acceptable infection risk has been imposed, the corresponding maximum exposure time value can be easily assessed. The four scenarios are examined assuming an acceptable risk value of 10 − 3 as discussed in Section 2.3.
In scenario A (patient emitting at rest in oral breathing), the maximum exposure time in a hospital room of 100 m 3 for both a medical staff member (scenario A-1) and a patient at rest without infection (scenario A-2) is evaluated. In both cases the maximum exposure time to reach an accepted risk of R max = 10 − 3 increase significantly with the ventilation rate, reaching 30 min and 72 min (scenario A-1), and 54 min and 192 min (scenario A-2) with AER values of 0.5 h − 1 and 10 h − 1 , respectively. If a higher individual risk is accepted, e.g. R max = 10 − 2 , the maximum exposure times increased up to 126 min and 660 min (scenario A-1), and 276 min and 1860 min (scenario A-2) with AER values of 0.5 h − 1 and 10 h − 1 , respectively. Such a higher accepted risk could be adopted, in principle, in view of guaranteeing R 0 ≤ 1 due to the reduced number of patients in a hospital room. Nonetheless, such a longer exposure time could be not adequate as the hospitalization time (and then the exposure time of susceptible individuals) could be likely longer than the maximum exposure times. In exposure scenario B (the gym with infected and healthy subjects during heavy activity with oral breathing), although there is no vocalization in the subject's activity, the high inhalation rate produces considerable ER q values, then increasing the individual risk; thus, in order to guarantee an acceptable infection risk of 10 − 3 the maximum exposure times resulted quite short, i.e. 12 min and 17 min for 0.5 h − 1 and 10 h − 1 , respectively. Once again, the maximum acceptable risk could be increased for specific gyms characterized by a no extremely high number of people, indeed, if a R max = 10 − 2 is accepted, the R 0 is lower than 1 for a maximum of 100 people simultaneously present in the gym for 132 min in the case of mechanical Fig. 2. Probability density functions of ER q (pdf ERq , expressed as P ERq ), probability of infection (pdf P , expressed as P I (ER q )), and individual infection risk (pdf R , expressed as R(ER q )) at t = 120 min for the illustrative example reported in Fig. 1 (exposure scenario D with an air exchange rate of 0.5 h − 1 ). The relative frequency P ERq is here graphed as 99 equally spaced log(ER q ) values (from 1st to 99th percentiles) adopting the log-normal probability density functions of ER q (pdf ERq ).  Table 1. ventilation at 10 h − 1 . On the contrary, for natural ventilation conditions, the maximum exposure time (43 min) appears still too short with respect the typical workout duration in gyms.
Moving to more crowded scenarios, the situation could be even more critical. As an example, for the exposure scenario discussed above (scenario D, AER = 0.5 h − 1 ) the maximum exposure time to reach an accepted risk of R = 10 − 3 is very short (14 min); this is due to the high viral load emitted during singing or speaking loudly leading to high quanta concentrations despite the large volume available. Obviously, the exposure time can increase with higher ventilation rates, e.g. reaching 16 min in the case of mechanical ventilation at 10 h − 1 . Nonetheless, people attending events in conference rooms or auditoriums are expected to remain therein longer than 14 or 21 min. As an example, as shown in the previous section, a 2 h exposure time in the scenario D would result in a risk of 3.8% and 0.7% for AERs of 0.5 h − 1 and 10 h − 1 , respectively. The crowding index of such an indoor environment (800 m 3 ) ranges from 0.75 m 2 (auditorium) to 2 m 2 (conference room) per person (European Committee for Standardisation, 2008); thus, for a room height of 4 m a corresponding floor area of 200 m 2 will be available, resulting in a total number of people simultaneously present in the room (S) ranging from 100 (conference room) to 267 (auditorium). Therefore, after 2 h of exposure in the case of natural ventilation, the value of R 0 will be much larger than 1 for both the auditorium (R 0 = 10.2) and the conference room (R 0 = 3.8). The R 0 values decrease to 2.0 and 0.7 for auditorium and conference room, respectively, in the case of mechanical ventilation at 10 h − 1 . The maximum number of subjects that could attend simultaneously the event in order to guarantee R 0 ≤1 is 26 and 135 for AERs equal to 0.5 h − 1 and 10 h − 1 , respectively. Thus, in the management of the epidemic, reducing the crowding index and exposure time could be essential. Accepting higher R max values would clearly increase the maximum exposure time; indeed, in the case of R max = 10 − 2 , the exposure time values would be 51 min and 165 min, for an AER equal to 0.5 h − 1 and 10 h − 1 , respectively.
In scenario C, the infected subject in light activity speaks in a 300 m 3 environment, along with the healthy subject. The simultaneous reduction of both the quanta emission rate and the volume compared to scenario D makes the maximum exposure times for an acceptable infection risk of 10 − 3 comparable to the previous case (14 min and 20 min for ventilation of 0.5 h − 1 and 10 h − 1 , respectively). Similar conclusions on the crowding index with respect to those resulting from the previous scenarios (D) can be easily derived.
Thus, for all the scenarios investigated, the ventilation conditions strongly influence the risk (or the exposure time) of the exposed subject: this difference increases as the accepted risk increases as shown in the trends presented in Fig. 3. In contrast, if a lower risk was accepted (i.e. 10 − 4 or 10 − 5 ), increasing the air exchange rate is not leading to the significant reduction of the risk, and localized air extraction near the infected person would be more effective.

Retrospective assessment application: The outbreaks in a restaurant in Guangzhou and at a choir rehearsal in Skagit Valley
The retrospective analysis applied to the restaurant in Guangzhou revealed that, under the boundary conditions considered in the simulation (in terms of room volume, ventilation, number of exposed people; see Section 2.4.1), a probability of infection (P I ) after 1 h of exposure equal to the attack rate (45%) can be reached for a quanta emission rate of ER q = 61 quanta h − 1 . Similarly, for the retrospective analysis applied to the Skagit Valley choir, in order to reach an attack rate of 53% after 2.5 h of exposure under the simulation boundary conditions reported in Section 2.4.2, a quanta emission rate of 341 quanta h − 1 is needed. In Fig. 4 the trends of quanta concentration and probability of infection (P I ) evaluated for the two cases adopting the quanta emission rates calculated through the retrospective application of the model (61 and 341 quanta h − 1 for the restaurant in Guangzhou and the choir rehearsal in Skagit Valley, respectively) are shown.
The emission rate obtained for the restaurant in Guangzhou occurs between 93rd and 94th percentile of the probability density function of ER q (pdf ERq ) characteristics of an emitting subject speaking during light exercise; similarly, the emission rate back-calculated for the choir rehearsal in Skagit Valley occurs between 92nd and 93rd percentile of the pdf ERq characteristics of an infected subject while singing. Therefore, such emission rates occur in close proximity of the mode of probability density function of R(ER q ), pdf R , discussed in the previous section. Therefore, for both the analyzed cases in the retrospective analyses, the required ER q values to obtain the documented probability of infection fall perfectly within the possible values of the emission profiles under consideration (i.e. speaking and singing/speaking loudly in light activity reported in Table 2). Moreover, such emission values present the highest probability of occurrence as they are around the 92nd-94th percentile, i. e. at the ER q percentile maximizing the individual infection risk (R (ER q )). In these two cases, an individual risk of < 10 − 3 would be not actually achievable by varying and optimizing the room ventilation (e.g. AER > 1000 h − 1 would be required), and is maybe achievable only by reducing the exposure time of the susceptible subjects and the quanta emission rates, and through advanced ventilation able to remove air exhaled by the infected subject before it is mixed with the room air.
To summarize, the retrospective assessment of the two SARS-CoV-2 outbreaks investigated demonstrate that the documented number of infected people can be explained by means of the airborne transmission route; indeed, the most probable of the expected events (i.e. quanta emission rates) occurred. The approach and consequent calculation reported here clearly highlights that the explanation of such a high number of infected people does not necessarily require the presence of a superspreader in the environment (i.e. an infected person with the highest viral load, c v , and infectious dose, c i ), but rather a co-existence of conditions, including emission and exposure parameters, leading to a highly probable event, which can be defined as a "superspreading event".

Conclusions
The present study proposed an approach aimed at providing a method for quantitative assessment about the airborne transmission risk of SARS-CoV-2. The proposed approach includes the application of the Monte Carlo method both to the estimation of the viral load emission of an infected subject and to the consequent risk assessment.
In the case of prospective assessments great attention must be paid to (i) the situations where specific expiratory activity and/or physical activities are conducted (e.g. subjects singing or speaking aloud, or are engaged in heavy exercising activity) as they can lead to high quanta concentrations and, then, risk, also in large closed environments; (ii) crowded indoor environments and air exchange rates which could lead to basic reproduction number R 0 > 1 also with reduced exposure times. Indeed, the findings of the study revealed that, the exposure times that guarantee an acceptable risk are very limited in typical environments with natural ventilation. In the case of high forced ventilation, the exposure times are longer, but are well below one hour For retrospective assessments, the application of the proposed approach to documented outbreaks highlighted the high probability of occurrence of the events due to the boundary conditions: (i) reduced ventilation and volume for the restaurant in Guangzhou; (ii) high emission of the singing infected subject and high exposure times for the Skagit Valley choir. In both cases, attack rate values can be justified on the basis of boundary conditions rather than the unlikely presence of a superspreader. Airborne transmission represents the main route of contagion in these cases.
The proposed approach is of great relevance as it represents an essential tool to be applied in enclosed space, and it can support public health experts, engineers and epidemiologists in planning exposure times for populations in indoor environments during an epidemic.