Merge branch 'feature/reproduction_rate' into 'master'

Implement an ExposureModel.reproduction_number, which is modelled by having precisely one infected person in the concentration model

Closes #98

See merge request cara/cara!88

cara/apps/calculator/report_generator.py

@@ -18,7 +18,7 @@ from .model_generator import FormData
 class RepeatEvents:
     repeats: int
     probability_of_infection: float
-    R0: float
+    expected_new_cases: float
 
 
 def calculate_report_data(model: models.ExposureModel):
@@ -35,7 +35,7 @@ def calculate_report_data(model: models.ExposureModel):
     prob = model.infection_probability()
     er = model.concentration_model.infected.emission_rate_when_present()
     exposed_occupants = model.exposed.number
-    r0 = model.reproduction_rate()
+    expected_new_cases = model.expected_new_cases()
 
     repeated_events = []
     for n in [1, 2, 3, 4, 5, 10, 15, 20]:
@@ -44,7 +44,7 @@ def calculate_report_data(model: models.ExposureModel):
             RepeatEvents(
                 repeats=n,
                 probability_of_infection=repeat_model.infection_probability(),
-                R0=repeat_model.reproduction_rate(),
+                expected_new_cases=repeat_model.expected_new_cases(),
             )
         )
 
@@ -55,7 +55,7 @@ def calculate_report_data(model: models.ExposureModel):
         "prob_inf": prob,
         "emission_rate": er,
         "exposed_occupants": exposed_occupants,
-        "R0": r0,
+        "expected_new_cases": expected_new_cases,
         "scenario_plot_src": embed_figure(plot(times, concentrations)),
         "repeated_events": repeated_events,
     }

cara/apps/calculator/templates/report.html.j2

@@ -129,7 +129,8 @@
 
 	<p class="result_title">Results:</p>
 		<p class="data_text">
-		In this scenario, the estimated probability of one exposed occupant getting infected P(i) is {{ prob_inf | int_format }}% and the expected number of new cases is {{ R0 | float_format }}.
+		In this scenario, the estimated probability of one exposed occupant getting infected P(i) is {{ prob_inf | int_format }}%
+		and the expected number of new cases is {{ expected_new_cases | float_format }}.
 		<p>
 	<p class="data_title">Exposure graph:</p>
 	<img id="scenario_concentration_plot" src="{{ scenario_plot_src }}">
@@ -151,7 +152,7 @@
                 <tr>
                     <td>{{ repeat_event.repeats }}</td>
                     <td>{{ repeat_event.probability_of_infection | int_format }}%</td>
-                    <td style="text-align:right">{{ repeat_event.R0 | float_format }}</td>
+                    <td style="text-align:right">{{ repeat_event.expected_new_cases | float_format }}</td>
                 </tr>
             {% endfor %}
             </tbody>

cara/apps/expert.py

@@ -118,8 +118,12 @@ class WidgetView:
             print(f'Probability of infection: {np.round(P, 0)}%')
 
             print(f'Number of exposed: {model.exposed.number}')
-            R0 = np.round(model.reproduction_rate(), 1)
-            print(f'Number of expected new cases (R0): {R0}')
+
+            new_cases = np.round(model.expected_new_cases(), 1)
+            print(f'Number of expected new cases: {new_cases}')
+
+            R0 = np.round(model.reproduction_number(), 1)
+            print(f'Reproduction number (R0): {R0}')
 
     def _build_widget(self, node):
         self.widget.children += (self._build_room(node.concentration_model.room),)

cara/dataclass_utils.py

@@ -0,0 +1,30 @@
+import dataclasses
+import typing
+
+
+DCInst = typing.TypeVar('T')
+
+
+def nested_replace(obj: DCInst, new_values: typing.Dict[str, typing.Any]) -> DCInst:
+    """
+    Replace an attribute on a dataclass, much like dataclasses.replace, except it
+    supports nested replacement definitions. For example:
+
+    >>> new_obj = nested_replace(obj, {'attr1.sub_attr2.sub_sub_attr3': 4})
+    >>> new_obj.attr1.sub_attr2.sub_sub_attr3
+    4
+
+    """
+    new_inst = obj
+    for name, value in new_values.items():
+        if '.' in name:
+            # Recurse into the desired name and come out with a top-level
+            # dataclass which has been updated appropriately.
+            name, remainder = name.split('.', 1)
+            value = nested_replace(
+                getattr(new_inst, name),
+                {remainder: value}
+            )
+        # We have a plain old name. So set it.
+        new_inst = dataclasses.replace(new_inst, **{name: value})
+    return new_inst

cara/models.py

@@ -1,10 +1,9 @@
+from dataclasses import dataclass
 import functools
 import numpy as np
 import typing
-from abc import abstractmethod
 
-
-from dataclasses import dataclass
+from .dataclass_utils import nested_replace
 
 
 @dataclass(frozen=True)
@@ -138,7 +137,6 @@ class Ventilation:
     def transition_times(self) -> typing.Set[float]:
         return self.active.transition_times()
 
-    @abstractmethod
     def air_exchange(self, room: Room, time: float) -> float:
         """
         Returns the rate at which air is being exchanged in the given room
@@ -169,7 +167,6 @@ class MultipleVentilation:
             transitions.update(ventilation.transition_times())
         return transitions
 
-    @abstractmethod
     def air_exchange(self, room: Room, time: float) -> float:
         """
         Returns the rate at which air is being exchanged in the given room
@@ -573,7 +570,24 @@ class ExposureModel:
         # Probability of infection.
         return (1 - np.exp(-inf_aero)) * 100
 
-    def reproduction_rate(self):
+    def expected_new_cases(self):
         prob = self.infection_probability()
         exposed_occupants = self.exposed.number
         return prob * exposed_occupants / 100
+
+    def reproduction_number(self):
+        """
+        The reproduction number can be thought of as the expected number of
+        cases directly generated by one infected case in a population.
+
+        """
+        if self.concentration_model.infected.number == 1:
+            return self.expected_new_cases()
+
+        # Create an equivalent exposure model but with precisely
+        # one infected case.
+        single_exposure_model = nested_replace(
+            self, {'concentration_model.infected.number': 1}
+        )
+
+        return single_exposure_model.expected_new_cases()

cara/tests/conftest.py

@@ -0,0 +1,38 @@
+from cara import models
+
+import pytest
+
+
+@pytest.fixture
+def baseline_model():
+    model = models.ConcentrationModel(
+        room=models.Room(volume=75),
+        ventilation=models.WindowOpening(
+            active=models.PeriodicInterval(period=120, duration=120),
+            inside_temp=models.PiecewiseConstant((0,24),(293,)),
+            outside_temp=models.PiecewiseConstant((0,24),(283,)),
+            cd_b=0.6, window_height=1.6, opening_length=0.6,
+        ),
+        infected=models.InfectedPopulation(
+            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 exercise'],
+            expiration=models.Expiration.types['Unmodulated Vocalization'],
+        ),
+    )
+    return model
+
+
+@pytest.fixture
+def baseline_exposure_model(baseline_model):
+    return models.ExposureModel(
+        baseline_model,
+        exposed=models.Population(
+            number=10,
+            presence=baseline_model.infected.presence,
+            activity=baseline_model.infected.activity,
+            mask=baseline_model.infected.mask,
+        )
+    )

cara/tests/test_dataclass_utils.py

@@ -0,0 +1,27 @@
+import dataclasses
+
+from cara.dataclass_utils import nested_replace
+
+
+@dataclasses.dataclass(frozen=True)
+class Four:
+    four: float
+
+
+@dataclasses.dataclass(frozen=True)
+class Two:
+    three: int
+    four: Four
+
+
+@dataclasses.dataclass(frozen=True)
+class One:
+    one: int
+    two: Two
+
+
+
+def test_nested_replace():
+    inst = One(1, two=Two(3, Four(4)))
+    new_inst = nested_replace(inst, {'two.four': Four(5)})
+    assert new_inst == One(1, two=Two(3, Four(5)))

cara/tests/test_known_quantities.py

@@ -24,41 +24,6 @@ def test_no_mask_emission_rate(baseline_model):
     )
 
 
-@pytest.fixture
-def baseline_model():
-    model = models.ConcentrationModel(
-        room=models.Room(volume=75),
-        ventilation=models.WindowOpening(
-            active=models.PeriodicInterval(period=120, duration=120),
-            inside_temp=models.PiecewiseConstant((0,24),(293,)),
-            outside_temp=models.PiecewiseConstant((0,24),(283,)),
-            cd_b=0.6, window_height=1.6, opening_length=0.6,
-        ),
-        infected=models.InfectedPopulation(
-            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 exercise'],
-            expiration=models.Expiration.types['Unmodulated Vocalization'],
-        ),
-    )
-    return model
-
-
-@pytest.fixture
-def baseline_exposure_model(baseline_model):
-    return models.ExposureModel(
-        baseline_model,
-        exposed=models.Population(
-            number=10,
-            presence=baseline_model.infected.presence,
-            activity=baseline_model.infected.activity,
-            mask=baseline_model.infected.mask,
-        )
-    )
-
-
 @pytest.fixture
 def baseline_periodic_window():
     return models.WindowOpening(

cara/tests/test_model.py

@@ -0,0 +1,12 @@
+import cara.models
+from cara.dataclass_utils import nested_replace
+
+
+def test_exposure_r0(baseline_exposure_model):
+    baseline_n3 = nested_replace(
+        baseline_exposure_model, {'concentration_model.infected.number': 3}
+    )
+    # The number of new cases should be greater if there are more infecteds, but
+    # the reproduction number should be the same (it is a measure of one infected case).
+    assert baseline_n3.expected_new_cases() > baseline_exposure_model.expected_new_cases()
+    assert baseline_n3.reproduction_number() == baseline_exposure_model.reproduction_number()