diff --git a/cara/models.py b/cara/models.py index 98245997..296c260c 100644 --- a/cara/models.py +++ b/cara/models.py @@ -620,7 +620,6 @@ _ExpirationBase.types = { 'Talking': Expiration((1., 1., 1.)), 'Shouting': Expiration((1., 5., 5.)), 'Singing': Expiration((1., 5., 5.)), - 'Superspreading event': Expiration((np.inf, 0., 0.)), } @@ -669,44 +668,21 @@ class Population: @dataclass(frozen=True) -class InfectedPopulation(Population): +class _PopulationWithVirus(Population): #: The virus with which the population is infected. virus: Virus - #: The type of expiration that is being emitted whilst doing the activity. - expiration: _ExpirationBase - @method_cache def emission_rate_when_present(self) -> _VectorisedFloat: """ - The emission rate if the infected population is present. - - Note that the rate is not currently time-dependent. - + The emission rate if the infected population is present + (in virions / h). It should not be a function of time. """ - # Emission Rate (virions / h) - # Note on units: exhalation rate is in m^3/h, aerosols in mL/cm^3 - # and viral load in virus/mL -> 1e6 conversion factor - aerosols = self.expiration.aerosols(self.mask) - - ER = (self.virus.viral_load_in_sputum * - self.activity.exhalation_rate * - 10 ** 6 * - aerosols) - - # For superspreading event, where ejection_factor is infinite we fix the ER - # based on Miller et al. (2020). - if isinstance(aerosols, np.ndarray): - ER[np.isinf(aerosols)] = 970 * self.virus.infectious_dose - elif np.isinf(aerosols): - ER = 970 * self.virus.infectious_dose - - return ER * self.number + raise NotImplementedError("Subclass must implement") def emission_rate(self, time) -> _VectorisedFloat: """ - The emission rate of the population. - + The emission rate of the population vs time. """ # Note: The original model avoids time dependence on the emission rate # at the cost of implementing a piecewise (on time) concentration function. @@ -722,11 +698,49 @@ class InfectedPopulation(Population): return self.emission_rate_when_present() +@dataclass(frozen=True) +class EmittingPopulation(_PopulationWithVirus): + #: The emission rate of a single individual, in virions / h. + known_individual_emission_rate: float + + @method_cache + def emission_rate_when_present(self) -> _VectorisedFloat: + """ + The emission rate if the infected population is present. + """ + return self.known_individual_emission_rate * self.number + + +@dataclass(frozen=True) +class InfectedPopulation(_PopulationWithVirus): + #: The type of expiration that is being emitted whilst doing the activity. + expiration: _ExpirationBase + + @method_cache + def emission_rate_when_present(self) -> _VectorisedFloat: + """ + The emission rate if the infected population is present. + Note that the rate is not currently time-dependent. + """ + # Emission Rate (virions / h) + # Note on units: exhalation rate is in m^3/h, aerosols in mL/cm^3 + # and viral load in virus/mL -> 1e6 conversion factor + + aerosols = self.expiration.aerosols(self.mask) + + ER = (self.virus.viral_load_in_sputum * + self.activity.exhalation_rate * + 10 ** 6 * + aerosols) + + return ER * self.number + + @dataclass(frozen=True) class ConcentrationModel: room: Room ventilation: _VentilationBase - infected: InfectedPopulation + infected: _PopulationWithVirus @property def virus(self): diff --git a/cara/tests/conftest.py b/cara/tests/conftest.py index 9499dfcf..e828e0cd 100644 --- a/cara/tests/conftest.py +++ b/cara/tests/conftest.py @@ -13,13 +13,15 @@ def baseline_model(): active=models.SpecificInterval(((0., 24.), )), air_exch=30., ), - infected=models.InfectedPopulation( + infected=models.EmittingPopulation( number=1, virus=models.Virus.types['SARS_CoV_2'], presence=models.SpecificInterval(((0., 4.), (5., 8.))), mask=models.Mask.types['No mask'], activity=models.Activity.types['Light activity'], - expiration=models.Expiration.types['Superspreading event'], + known_individual_emission_rate=970 * 50, + # superspreading event, where ejection factor is fixed based + # on Miller et al. (2020) - 50 represents the infectious dose. ), ) return model diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index 11a135a1..ac82dae3 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -144,13 +144,15 @@ def conc_model(): return models.ConcentrationModel( models.Room(25), models.AirChange(always, 5), - models.InfectedPopulation( + models.EmittingPopulation( number=1, presence=interesting_times, mask=models.Mask.types['No mask'], activity=models.Activity.types['Seated'], virus=models.Virus.types['SARS_CoV_2'], - expiration=models.Expiration.types['Superspreading event'], + known_individual_emission_rate=970 * 50, + # superspreading event, where ejection factor is fixed based + # on Miller et al. (2020) - 50 represents the infectious dose. ) ) diff --git a/cara/tests/test_known_quantities.py b/cara/tests/test_known_quantities.py index 89c0268d..d1db42c2 100644 --- a/cara/tests/test_known_quantities.py +++ b/cara/tests/test_known_quantities.py @@ -66,13 +66,15 @@ def build_model(interval_duration): active=models.PeriodicInterval(period=120, duration=interval_duration), q_air_mech=500., ), - infected=models.InfectedPopulation( + infected=models.EmittingPopulation( number=1, virus=models.Virus.types['SARS_CoV_2'], presence=models.SpecificInterval(((0., 4.), (5., 8.))), mask=models.Mask.types['No mask'], activity=models.Activity.types['Light activity'], - expiration=models.Expiration.types['Superspreading event'], + known_individual_emission_rate=970 * 50, + # superspreading event, where ejection factor is fixed based + # on Miller et al. (2020) - 50 represents the infectious dose. ), ) return model @@ -226,13 +228,13 @@ def build_hourly_dependent_model( outside_temp=outside_temp, window_height=1.6, opening_length=0.6, ), - infected=models.InfectedPopulation( + infected=models.EmittingPopulation( number=1, virus=models.Virus.types['SARS_CoV_2'], presence=models.SpecificInterval(intervals_presence_infected), mask=models.Mask.types['No mask'], activity=models.Activity.types['Light activity'], - expiration=models.Expiration.types['Superspreading event'], + known_individual_emission_rate=970 * 50, ), ) return model @@ -247,13 +249,13 @@ def build_constant_temp_model(outside_temp, intervals_open=((7.5, 8.5),)): outside_temp=models.PiecewiseConstant((0., 24.), (outside_temp,)), window_height=1.6, opening_length=0.6, ), - infected=models.InfectedPopulation( + infected=models.EmittingPopulation( number=1, virus=models.Virus.types['SARS_CoV_2'], presence=models.SpecificInterval(((0., 4.), (5., 7.5))), mask=models.Mask.types['No mask'], activity=models.Activity.types['Light activity'], - expiration=models.Expiration.types['Superspreading event'], + known_individual_emission_rate=970 * 50, ), ) return model @@ -275,13 +277,13 @@ def build_hourly_dependent_model_multipleventilation(month, intervals_open=((7.5 model = models.ConcentrationModel( room=models.Room(volume=75), ventilation=vent, - infected=models.InfectedPopulation( + infected=models.EmittingPopulation( number=1, virus=models.Virus.types['SARS_CoV_2'], presence=models.SpecificInterval(((0., 4.), (5., 7.5))), mask=models.Mask.types['No mask'], activity=models.Activity.types['Light activity'], - expiration=models.Expiration.types['Superspreading event'], + known_individual_emission_rate=970 * 50, ), ) return model