From 9ed65a2b044448483dfc6b3de34047705cb0e3b2 Mon Sep 17 00:00:00 2001 From: Nicolas Mounet Date: Mon, 13 Sep 2021 12:00:06 +0200 Subject: [PATCH 1/5] Model optimisation: avoid computing before first presence time --- cara/models.py | 14 +++++++++++++- cara/tests/models/test_concentration_model.py | 4 ++++ 2 files changed, 17 insertions(+), 1 deletion(-) diff --git a/cara/models.py b/cara/models.py index 53cf52bc..86d2acb9 100644 --- a/cara/models.py +++ b/cara/models.py @@ -776,6 +776,14 @@ class ConcentrationModel: state_change_times.update(self.ventilation.transition_times()) return sorted(state_change_times) + @method_cache + def _first_presence_time(self) -> float: + """ + First presence time. Before that, the concentration is zero. + + """ + return self.infected.presence.boundaries()[0][0] + def last_state_change(self, time: float) -> float: """ Find the most recent/previous state change. @@ -823,7 +831,9 @@ class ConcentrationModel: Note that time is not vectorised. You can only pass a single float to this method. """ - if time == 0: + # The model always starts at t=0, but we avoid running concentration calculations + # before the first presence as an optimisation. + if time <= self._first_presence_time(): return 0.0 next_state_change_time = self._next_state_change(time) IVRR = self.infectious_virus_removal_rate(next_state_change_time) @@ -841,6 +851,8 @@ class ConcentrationModel: """ Get the integrated concentration dose between the times start and stop. """ + if stop <= self._first_presence_time(): + return 0.0 state_change_times = self.state_change_times() req_start, req_stop = start, stop total_concentration = 0. diff --git a/cara/tests/models/test_concentration_model.py b/cara/tests/models/test_concentration_model.py index 46e7e88e..9373d4a9 100644 --- a/cara/tests/models/test_concentration_model.py +++ b/cara/tests/models/test_concentration_model.py @@ -121,6 +121,10 @@ def test_next_state_change_time_out_of_range(simple_conc_model: models.Concentra simple_conc_model._next_state_change(3.1) +def test_first_presence_time(simple_conc_model): + assert simple_conc_model._first_presence_time() == 0.5 + + def test_integrated_concentration(simple_conc_model): c1 = simple_conc_model.integrated_concentration(0, 2) c2 = simple_conc_model.integrated_concentration(0, 1) From fdca7da53a0b8cd5866d688307681e8a11007356 Mon Sep 17 00:00:00 2001 From: jdevine Date: Mon, 13 Sep 2021 13:50:12 +0200 Subject: [PATCH 2/5] Changed default strain in calculator to Delta. Also applies retroactively to QR codes missing delta. The rational for the retrospective change is that old simulations (wild type) re-run need to be re-run with the current dominant strain. --- cara/apps/calculator/model_generator.py | 2 +- cara/apps/calculator/templates/calculator.form.html.j2 | 6 +++--- cara/apps/calculator/templates/userguide.html.j2 | 4 ++-- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index 06c29229..61f4e9fd 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -117,7 +117,7 @@ class FormData: 'simulation_name': _NO_DEFAULT, 'total_people': _NO_DEFAULT, 'ventilation_type': 'no_ventilation', - 'virus_type': 'SARS_CoV_2', + 'virus_type': 'SARS_CoV_2_B16172', 'volume_type': _NO_DEFAULT, 'window_type': 'window_sliding', 'window_height': 0., diff --git a/cara/apps/calculator/templates/calculator.form.html.j2 b/cara/apps/calculator/templates/calculator.form.html.j2 index c119d4e5..97209c98 100644 --- a/cara/apps/calculator/templates/calculator.form.html.j2 +++ b/cara/apps/calculator/templates/calculator.form.html.j2 @@ -62,10 +62,10 @@ v{{ calculator_version }} Please sen
@@ -344,13 +344,13 @@ v{{ calculator_version }} Please sen Quick Guide:
This tool simulates the long range airborne spread SARS-CoV-2 virus in a finite volume and estimates the risk of COVID-19 infection. It is based on current scientific data and can be used to compare the effectiveness of different mitigation measures.
Virus data:
- SARS-CoV-2 covers typical strains of the virus and three variants of concern (VOC):
+ SARS-CoV-2 covers the original wild strain of the virus and three variants of concern (VOC):
  • Alpha (also known as B.1.1.7, first identified in UK, Dec 2020),
  • Gamma (also known as P.1, first identified in Brazil/Japan, Jan 2021).
  • Delta (also known as B.1.617.2, first identified in India, Oct 2020).
- Choose variant according to local area prevalence, e.g. for Geneva + Modify the default (Delta) as necessary, according to local area prevalence e.g. for Geneva or Ain (France).
Ventilation data:
    diff --git a/cara/apps/calculator/templates/userguide.html.j2 b/cara/apps/calculator/templates/userguide.html.j2 index ede8c8e3..ac0e69ca 100644 --- a/cara/apps/calculator/templates/userguide.html.j2 +++ b/cara/apps/calculator/templates/userguide.html.j2 @@ -56,12 +56,12 @@ Changing this setting alters the properties of the virus which are used for the This has a significant effect on the probability of infection. The choices are:

      -
    • SARS-CoV-2 (nominal strain), covering typical strains and varaints which are not of concern from an epidemiologic point of view of the virus;
      • +
    • SARS-CoV-2 (nominal strain), covering typical strains and variants which are not of concern from an epidemiologic point of view of the virus;
    • SARS-CoV-2 (Alpha VOC), first identified in the UK at the end of 2020 which is found to be approximately 1.5x more transmissible compared to the non-VOCs;
    • SARS-CoV-2 (Gamma VOC), first identified in Brazil in January 2021 which is found to be approximately 2.2x more transmissible compared to the non-VOCs.
    • SARS-CoV-2 (Delta VOC), first identified in India towards the end of 2020 which is found to be approximately 60% more transmissible compared to the ALPHA VOC.
    -

    The user can base their choice according to the prevalence of the different variants in the local area. Access to this information can be found here:

    +

    The user can modify the selected variant from the default (Delta VOC), according to the prevalence of the different variants in the local area. Access to this information can be found here:

    • Geneva: https://www.covid19.admin.ch/fr/epidemiologic/virus-variants?detGeo=GE
    • Ain (France): https://www.santepubliquefrance.fr/dossiers/coronavirus-covid-19/covid-19-cartographie-des-variants-en-france-donnees-par-region-et-par-departement
      • From 8fbe7856b98cd908934ee9798b1ac06039fb40ae Mon Sep 17 00:00:00 2001 From: jdevine Date: Tue, 14 Sep 2021 14:57:48 +0200 Subject: [PATCH 3/5] Addressed comments from Andre. --- cara/apps/calculator/__init__.py | 2 +- cara/apps/calculator/templates/calculator.form.html.j2 | 6 +++--- cara/apps/calculator/templates/userguide.html.j2 | 2 +- 3 files changed, 5 insertions(+), 5 deletions(-) diff --git a/cara/apps/calculator/__init__.py b/cara/apps/calculator/__init__.py index 1688d494..a60d5f74 100644 --- a/cara/apps/calculator/__init__.py +++ b/cara/apps/calculator/__init__.py @@ -33,7 +33,7 @@ from .user import AuthenticatedUser, AnonymousUser # calculator version. If the calculator needs to make breaking changes (e.g. change # form attributes) then it can also increase its MAJOR version without needing to # increase the overall CARA version (found at ``cara.__version__``). -__version__ = "3.0.0" +__version__ = "3.0.1" class BaseRequestHandler(RequestHandler): diff --git a/cara/apps/calculator/templates/calculator.form.html.j2 b/cara/apps/calculator/templates/calculator.form.html.j2 index 97209c98..65391454 100644 --- a/cara/apps/calculator/templates/calculator.form.html.j2 +++ b/cara/apps/calculator/templates/calculator.form.html.j2 @@ -62,10 +62,10 @@ v{{ calculator_version }} Please sen
      @@ -344,13 +344,13 @@ v{{ calculator_version }} Please sen Quick Guide:
      This tool simulates the long range airborne spread SARS-CoV-2 virus in a finite volume and estimates the risk of COVID-19 infection. It is based on current scientific data and can be used to compare the effectiveness of different mitigation measures.
      Virus data:
      - SARS-CoV-2 covers the original wild strain of the virus and three variants of concern (VOC):
      + SARS-CoV-2 covers the original "wild type" strain of the virus and three variants of concern (VOC):
      • Alpha (also known as B.1.1.7, first identified in UK, Dec 2020),
      • Gamma (also known as P.1, first identified in Brazil/Japan, Jan 2021).
      • Delta (also known as B.1.617.2, first identified in India, Oct 2020).
      - Modify the default (Delta) as necessary, according to local area prevalence e.g. for Geneva + Modify the default as necessary, according to local area prevalence e.g. for Geneva or Ain (France).
      Ventilation data:
        diff --git a/cara/apps/calculator/templates/userguide.html.j2 b/cara/apps/calculator/templates/userguide.html.j2 index ac0e69ca..f51cba61 100644 --- a/cara/apps/calculator/templates/userguide.html.j2 +++ b/cara/apps/calculator/templates/userguide.html.j2 @@ -61,7 +61,7 @@ The choices are:

      • SARS-CoV-2 (Gamma VOC), first identified in Brazil in January 2021 which is found to be approximately 2.2x more transmissible compared to the non-VOCs.
      • SARS-CoV-2 (Delta VOC), first identified in India towards the end of 2020 which is found to be approximately 60% more transmissible compared to the ALPHA VOC.
      -

      The user can modify the selected variant from the default (Delta VOC), according to the prevalence of the different variants in the local area. Access to this information can be found here:

      +

      The user can modify the selected variant from the default, according to the prevalence of the different variants in the local area. Access to this information can be found here:

      • Geneva: https://www.covid19.admin.ch/fr/epidemiologic/virus-variants?detGeo=GE
      • Ain (France): https://www.santepubliquefrance.fr/dossiers/coronavirus-covid-19/covid-19-cartographie-des-variants-en-france-donnees-par-region-et-par-departement
        • From 92f9bce7e5c2de2de65ad13ea3d96139bb32c2b7 Mon Sep 17 00:00:00 2001 From: Nicolas Mounet Date: Tue, 14 Sep 2021 16:11:44 +0200 Subject: [PATCH 4/5] Normalization of the computation of the concentration --- cara/models.py | 94 ++++++++++++++-------- cara/tests/models/test_exposure_model.py | 16 ++-- cara/tests/test_monte_carlo_full_models.py | 2 +- 3 files changed, 70 insertions(+), 42 deletions(-) diff --git a/cara/models.py b/cara/models.py index 86d2acb9..98245997 100644 --- a/cara/models.py +++ b/cara/models.py @@ -701,11 +701,11 @@ class InfectedPopulation(Population): elif np.isinf(aerosols): ER = 970 * self.virus.infectious_dose - return ER + return ER * self.number - def individual_emission_rate(self, time) -> _VectorisedFloat: + def emission_rate(self, time) -> _VectorisedFloat: """ - The emission rate of a single individual in the population. + The emission rate of the population. """ # Note: The original model avoids time dependence on the emission rate @@ -721,13 +721,6 @@ class InfectedPopulation(Population): return self.emission_rate_when_present() - def emission_rate(self, time) -> _VectorisedFloat: - """ - The emission rate of the entire population. - - """ - return self.individual_emission_rate(time) * self.number - @dataclass(frozen=True) class ConcentrationModel: @@ -753,16 +746,22 @@ class ConcentrationModel: ) @method_cache - def _concentration_limit(self, time: float) -> _VectorisedFloat: + def _normed_concentration_limit(self, time: float) -> _VectorisedFloat: """ Provides a constant that represents the theoretical asymptotic value reached by the concentration when time goes to infinity, if all parameters were to stay time-independent. + This is normalized by the emission rate, the latter acting as a + multiplicative constant factor for the concentration model that + can be put back in front of the concentration after the time + dependence has been solved for. """ + if not self.infected.person_present(time): + return 0. V = self.room.volume IVRR = self.infectious_virus_removal_rate(time) - return (self.infected.emission_rate(time)) / (IVRR * V) + return 1. / (IVRR * V) @method_cache def state_change_times(self) -> typing.List[float]: @@ -815,15 +814,16 @@ class ConcentrationModel: ) @method_cache - def _concentration_cached(self, time: float) -> _VectorisedFloat: - # A cached version of the concentration method. Use this method if you - # expect that there may be multiple concentration calculations for the - # same time (e.g. at state change times). - return self.concentration(time) + def _normed_concentration_cached(self, time: float) -> _VectorisedFloat: + # A cached version of the _normed_concentration method. Use this + # method if you expect that there may be multiple concentration + # calculations for the same time (e.g. at state change times). + return self._normed_concentration(time) - def concentration(self, time: float) -> _VectorisedFloat: + def _normed_concentration(self, time: float) -> _VectorisedFloat: """ - Virus exposure concentration, as a function of time. + Virus exposure concentration, as a function of time, and + normalized by the emission rate. The formulas used here assume that all parameters (ventilation, emission rate) are constant between two state changes - only the value of these parameters at the next state change, are used. @@ -837,25 +837,36 @@ class ConcentrationModel: return 0.0 next_state_change_time = self._next_state_change(time) IVRR = self.infectious_virus_removal_rate(next_state_change_time) - concentration_limit = self._concentration_limit(next_state_change_time) + conc_limit = self._normed_concentration_limit(next_state_change_time) t_last_state_change = self.last_state_change(time) - concentration_at_last_state_change = self._concentration_cached(t_last_state_change) + conc_at_last_state_change = self._normed_concentration_cached(t_last_state_change) delta_time = time - t_last_state_change fac = np.exp(-IVRR * delta_time) - return concentration_limit * (1 - fac) + concentration_at_last_state_change * fac + return conc_limit * (1 - fac) + conc_at_last_state_change * fac + + def concentration(self, time: float) -> _VectorisedFloat: + """ + Virus exposure concentration, as a function of time. + + Note that time is not vectorised. You can only pass a single float + to this method. + """ + return (self._normed_concentration(time) * + self.infected.emission_rate_when_present()) @method_cache - def integrated_concentration(self, start: float, stop: float) -> _VectorisedFloat: + def normed_integrated_concentration(self, start: float, stop: float) -> _VectorisedFloat: """ - Get the integrated concentration dose between the times start and stop. + Get the integrated concentration dose between the times start and stop, + normalized by the emission rate. """ if stop <= self._first_presence_time(): return 0.0 state_change_times = self.state_change_times() req_start, req_stop = start, stop - total_concentration = 0. + total_normed_concentration = 0. for interval_start, interval_stop in zip(state_change_times[:-1], state_change_times[1:]): if req_start > interval_stop or req_stop < interval_start: continue @@ -863,17 +874,24 @@ class ConcentrationModel: start = max([interval_start, req_start]) stop = min([interval_stop, req_stop]) - conc_start = self._concentration_cached(start) + conc_start = self._normed_concentration_cached(start) next_conc_state = self._next_state_change(stop) - conc_limit = self._concentration_limit(next_conc_state) + conc_limit = self._normed_concentration_limit(next_conc_state) IVRR = self.infectious_virus_removal_rate(next_conc_state) delta_time = stop - start - total_concentration += ( + total_normed_concentration += ( conc_limit * delta_time + (conc_limit - conc_start) * (np.exp(-IVRR*delta_time)-1) / IVRR ) - return total_concentration + return total_normed_concentration + + def integrated_concentration(self, start: float, stop: float) -> _VectorisedFloat: + """ + Get the integrated concentration dose between the times start and stop. + """ + return (self.normed_integrated_concentration(start, stop) * + self.infected.emission_rate_when_present()) @dataclass(frozen=True) @@ -890,14 +908,22 @@ class ExposureModel: #: The fraction of viruses actually deposited in the respiratory tract fraction_deposited: _VectorisedFloat = 0.6 - def exposure(self) -> _VectorisedFloat: - """The number of virus per meter^3.""" - exposure = 0.0 + def _normed_exposure(self) -> _VectorisedFloat: + """ + The number of virus per meter^3, normalized by the emission rate + of the infected population. + """ + normed_exposure = 0.0 for start, stop in self.exposed.presence.boundaries(): - exposure += self.concentration_model.integrated_concentration(start, stop) + normed_exposure += self.concentration_model.normed_integrated_concentration(start, stop) - return exposure * self.repeats + return normed_exposure * self.repeats + + def exposure(self) -> _VectorisedFloat: + """The number of virus per meter^3.""" + return (self._normed_exposure() * + self.concentration_model.infected.emission_rate_when_present()) def infection_probability(self) -> _VectorisedFloat: exposure = self.exposure() diff --git a/cara/tests/models/test_exposure_model.py b/cara/tests/models/test_exposure_model.py index cf8cfb6e..11a135a1 100644 --- a/cara/tests/models/test_exposure_model.py +++ b/cara/tests/models/test_exposure_model.py @@ -11,20 +11,20 @@ from cara.dataclass_utils import replace @dataclass(frozen=True) -class KnownConcentrations(models.ConcentrationModel): +class KnownNormedconcentration(models.ConcentrationModel): """ A ConcentrationModel which is based on pre-known exposure concentrations and which therefore doesn't need other components. Useful for testing. """ - concentration_function: typing.Callable + normed_concentration_function: typing.Callable def infectious_virus_removal_rate(self, time: float) -> models._VectorisedFloat: # very large decay constant -> same as constant concentration return 1.e50 - def _concentration_limit(self, time: float) -> models._VectorisedFloat: - return self.concentration_function(time) + def _normed_concentration_limit(self, time: float) -> models._VectorisedFloat: + return self.normed_concentration_function(time) def state_change_times(self): return [0., 24.] @@ -32,8 +32,8 @@ class KnownConcentrations(models.ConcentrationModel): def _next_state_change(self, time: float): return 24. - def concentration(self, time: float) -> models._VectorisedFloat: # noqa - return self.concentration_function(time) + def _normed_concentration(self, time: float) -> models._VectorisedFloat: # noqa + return self.normed_concentration_function(time) halftime = models.PeriodicInterval(120, 60) @@ -67,7 +67,9 @@ def known_concentrations(func): virus=models.Virus.types['SARS_CoV_2_B117'], expiration=models.Expiration.types['Talking'] ) - return KnownConcentrations(dummy_room, dummy_ventilation, dummy_infected_population, func) + normed_func = lambda x: func(x) / dummy_infected_population.emission_rate_when_present() + return KnownNormedconcentration(dummy_room, dummy_ventilation, + dummy_infected_population, normed_func) @pytest.mark.parametrize( diff --git a/cara/tests/test_monte_carlo_full_models.py b/cara/tests/test_monte_carlo_full_models.py index 5cf88198..08bd2dad 100644 --- a/cara/tests/test_monte_carlo_full_models.py +++ b/cara/tests/test_monte_carlo_full_models.py @@ -238,7 +238,7 @@ def skagit_chorale_mc(): ["shared_office_mc", 10.7, 0.32, 57.24, 654], ["classroom_mc", 36.1, 6.85, 780.0, 28464], ["ski_cabin_mc", 16.3, 0.49, 35.94, 7404], - ["gym_mc", 2.25, 0.63, 0.7842, 984], + ["gym_mc", 2.25, 0.63, 0.7842, 1968], ["waiting_room_mc", 9.72, 1.36, 34.26, 3534], ["skagit_chorale_mc",29.9, 17.9, 190.0, 141400], ] From 5d4dc4b91d7b3c22bc30ded340e7a65f3cd2f3aa Mon Sep 17 00:00:00 2001 From: jdevine Date: Tue, 14 Sep 2021 16:22:49 +0200 Subject: [PATCH 5/5] reverted model_generator default virus type --- cara/apps/calculator/model_generator.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index 61f4e9fd..06c29229 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -117,7 +117,7 @@ class FormData: 'simulation_name': _NO_DEFAULT, 'total_people': _NO_DEFAULT, 'ventilation_type': 'no_ventilation', - 'virus_type': 'SARS_CoV_2_B16172', + 'virus_type': 'SARS_CoV_2', 'volume_type': _NO_DEFAULT, 'window_type': 'window_sliding', 'window_height': 0.,