pingu_98/cara
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Merge branch 'feature/scientific_model_update_tests' into 'feature/scientific_model_update'

Browse source 
Re-enable pipeline on scientific model updates

See merge request cara/cara!268
This commit is contained in:
Andre Henriques 2021-11-10 17:19:39 +01:00
commit 8ec08c08ca
14 changed files with 253 additions and 155 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

16
.gitlab-ci.yml
View file

@ -11,13 +11,13 @@ variables:
# A full installation of CARA, tested with pytest.
# test_install:
# extends: .acc_py_full_test
test_install:
extends: .acc_py_full_test
# A development installation of CARA tested with pytest.
# test_dev:
# extends: .acc_py_dev_test
test_dev:
extends: .acc_py_dev_test
# A development installation of CARA tested with pytest.
@ -69,10 +69,10 @@ check_openshift_config_prod:
# A development installation of CARA tested with pytest.
# test_dev-39:
# variables:
# PY_VERSION: "3.9"
# extends: .acc_py_dev_test
test_dev-39:
variables:
PY_VERSION: "3.9"
extends: .acc_py_dev_test
.image_builder:

3
cara/apps/calculator/model_generator.py
View file

@ -246,8 +246,9 @@ class FormData:
room=room,
ventilation=self.ventilation(),
infected=self.infected_population(),
evaporation_factor=0.3,
),
exposed=self.exposed_population()
exposed=self.exposed_population(),
)
def build_model(self, sample_size=_DEFAULT_MC_SAMPLE_SIZE) -> models.ExposureModel:

3
cara/apps/expert.py
View file

@ -498,9 +498,10 @@ baseline_model = models.ExposureModel(
presence=models.SpecificInterval(((8., 12.), (13., 17.))),
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Seated'],
expiration=models.Expiration.types['Talking'],
expiration=models.Expiration.types['Speaking'],
host_immunity=0.,
),
evaporation_factor=0.3,
),
exposed=models.Population(
number=10,

39
cara/data/__init__.py
View file

@ -3,7 +3,7 @@ from cara import models
# TODO: The values in this module to be removed and instead use the cara.data.weather functionality.
# average temperature of each month, hour per hour (from midnight to 11 pm)
# Geneva average temperature of each month, hour per hour (from midnight to 11 pm)
Geneva_hourly_temperatures_celsius_per_hour = {
'Jan': [0.2, -0.3, -0.5, -0.9, -1.1, -1.4, -1.5, -1.5, -1.1, 0.1, 1.5,
2.8, 3.8, 4.4, 4.5, 4.4, 4.4, 3.9, 3.1, 2.7, 2.2, 1.7, 1.5, 1.1],
@ -31,6 +31,22 @@ Geneva_hourly_temperatures_celsius_per_hour = {
4.7, 5.2, 5.3, 5.2, 5.2, 4.7, 4.0, 3.7, 3.2, 2.8, 2.6, 2.2]
}
# Toronto average temperature of each month, hour per hour (from midnight to 11 pm)
Toronto_hourly_temperatures_celsius_per_hour = {
"Jan": [ -2.9, -3.0, -3.2, -3.3, -3.3, -3.5, -3.7, -3.8, -3.9, -4.0, -4.1, -4.3, -4.3, -4.3, -4.1, -3.7, -3.2, -2.8, -2.6, -2.3, -2.2, -2.3, -2.6, -2.8],
"Feb": [ -2.4, -2.6, -2.8, -2.9, -2.9, -3.1, -3.3, -3.4, -3.6, -3.8, -3.9, -4.0, -4.3, -4.2, -3.7, -3.2, -2.6, -2.1, -1.7, -1.5, -1.3, -1.4, -1.6, -2.1],
"Mar": [ 1.3, 1.0, 0.7, 0.5, 0.4, 0.1, -0.0, -0.2, -0.4, -0.5, -0.7, -0.8, -0.9, -0.3, 0.4, 1.0, 1.6, 2.0, 2.3, 2.7, 2.7, 2.7, 2.4, 1.9],
"Apr": [ 6.8, 6.5, 6.3, 5.9, 5.7, 5.4, 5.1, 4.9, 4.6, 4.4, 4.2, 4.3, 4.8, 5.5, 6.1, 6.7, 7.1, 7.6, 7.8, 8.1, 8.2, 8.2, 8.0, 7.6 ],
"May": [ 13.0, 12.6, 12.2, 11.8, 11.5, 11.2, 10.8, 10.5, 10.2, 9.9, 9.8, 10.0, 10.9, 11.6, 12.2, 12.7, 13.2, 13.6, 13.9, 14.3, 14.4, 14.3, 14.2, 13.8],
"Jun": [ 18.9, 18.2, 17.8, 17.4, 17.0, 16.6, 16.2, 15.9, 15.6, 15.4, 15.3, 15.8, 16.5, 17.3, 17.9, 18.4, 18.9, 19.4, 19.7, 20.1, 20.3, 20.3, 20.1, 19.7],
"Jul": [ 22.1, 21.4, 20.9, 20.5, 20.0, 19.6, 19.1, 18.9, 18.6, 18.3, 18.1, 18.5, 19.4, 20.3, 21.0, 21.6, 22.2, 22.7, 23.1, 23.4, 23.6, 23.5, 23.3, 22.9],
"Aug": [ 22.0, 21.4, 21.0, 20.7, 20.3, 20.0, 19.6, 19.3, 19.1, 18.8, 18.5, 18.4, 19.4, 20.3, 21.1, 21.7, 22.2, 22.8, 23.1, 23.4, 23.5, 23.3, 23.1, 22.6],
"Sep": [ 18.2, 17.8, 17.4, 17.3, 17.0, 16.6, 16.3, 16.0, 15.8, 15.5, 15.4, 15.0, 15.6, 16.6, 17.5, 18.2, 18.7, 19.2, 19.6, 19.8, 19.7, 19.6, 19.2, 18.6],
"Oct": [ 11.1, 10.9, 10.6, 10.5, 10.2, 10.1, 9.8, 9.7, 9.5, 9.3, 9.2, 9.0, 9.0, 9.7, 10.5, 11.2, 11.7, 12.2, 12.4, 12.6, 12.6, 12.3, 11.8, 11.3],
"Nov": [ 5.3, 5.1, 5.0, 4.7, 4.6, 4.4, 4.3, 4.2, 4.1, 4.0, 3.9, 3.8, 3.7, 4.0, 4.6, 5.2, 5.7, 6.1, 6.2, 6.4, 6.3, 6.0, 5.5, 5.3],
"Dec": [ 0.4, 0.3, 0.2, 0.0, -0.1, -0.2, -0.4, -0.5, -0.6, -0.7, -0.8, -0.8, -0.9, -0.9, -0.6, -0.2, 0.3, 0.7, 0.9, 1.1, 1.1, 0.9, 0.6, 0.5]
}
# Geneva hourly temperatures as piecewise constant function (in Kelvin).
GenevaTemperatures_hourly = {
@ -44,8 +60,27 @@ GenevaTemperatures_hourly = {
}
# Same temperatures on a finer temperature mesh (every 6 minutes).
# Toronto hourly temperatures as piecewise constant function (in Kelvin).
TorontoTemperatures_hourly = {
month: models.PiecewiseConstant(
# NOTE: It is important that the time type is float, not np.float, in
# order to allow hashability (for caching).
tuple(float(time) for time in range(25)),
tuple(273.15 + np.array(temperatures)),
)
for month, temperatures in Toronto_hourly_temperatures_celsius_per_hour.items()
}
# Same Geneva temperatures on a finer temperature mesh (every 6 minutes).
GenevaTemperatures = {
month: GenevaTemperatures_hourly[month].refine(refine_factor=10)
for month, temperatures in Geneva_hourly_temperatures_celsius_per_hour.items()
}
# Same Toronto temperatures on a finer temperature mesh (every 6 minutes).
TorontoTemperatures = {
month: TorontoTemperatures_hourly[month].refine(refine_factor=10)
for month, temperatures in Toronto_hourly_temperatures_celsius_per_hour.items()
}

10
cara/models.py
View file

@ -628,7 +628,7 @@ class MultipleExpiration(_ExpirationBase):
# The correspondence with the BLO coefficients is given.
_ExpirationBase.types = {
'Breathing': Expiration(1.3844), # corresponds to B/L/O coefficients of (1, 0, 0)
'Talking': Expiration(5.8925), # corresponds to B/L/O coefficients of (1, 1, 1)
'Speaking': Expiration(5.8925), # corresponds to B/L/O coefficients of (1, 1, 1)
'Shouting': Expiration(10.0411), # corresponds to B/L/O coefficients of (1, 5, 5)
'Singing': Expiration(10.0411), # corresponds to B/L/O coefficients of (1, 5, 5)
}
@ -788,14 +788,16 @@ class ConcentrationModel:
if (isinstance(self.infected, InfectedPopulation)
and isinstance(self.infected.expiration, Expiration)):
d = self.infected.expiration.diameter * self.evaporation_factor
vg = 1.88e-4 * (d / 2.5)**2 # see CERN-OPEN-2021-04
vg = 1.88e-4 * (d / 2.5)**2
# see https://doi.org/10.1101/2021.10.14.21264988
# (velocity of 1.88e-4 corresponds to diameter of 2.5 microns)
else:
# model is not evaluated for specific values of aerosol
# diameters - we choose a single velocity value
# corresponding to that obtained with a diameter of 2.5 microns
# (geometric average of the breathing expiration distribution,
# taking evaporation into account, see CERN-OPEN-2021-04)
# taking evaporation into account, see
# https://doi.org/10.1101/2021.10.14.21264988)
vg = 1.88e-4
# Height of the emission source to the floor - i.e. mouth/nose (m)
h = 1.5
@ -983,7 +985,7 @@ class ExposureModel:
and isinstance(self.concentration_model.infected.expiration,Expiration)):
# model is not evaluated for specific values of aerosol
# diameters - we choose a single "average" deposition factor,
# as in CERN-OPEN-2021-04.
# as in https://doi.org/10.1101/2021.10.14.21264988.
fdep = 0.6
else:
# deposition factor depends on aerosol particle diameter.

6
cara/monte_carlo/data.py
View file

@ -65,7 +65,7 @@ class BLOmodel:
return result
# From CERN-OPEN-2021-04 and refererences therein
# From https://doi.org/10.1101/2021.10.14.21264988 and refererences therein
activity_distributions = {
'Seated': mc.Activity(LogNormal(-0.6872121723362303, 0.10498338229297108),
LogNormal(-0.6872121723362303, 0.10498338229297108)),
@ -84,7 +84,7 @@ activity_distributions = {
}
# From CERN-OPEN-2021-04 and refererences therein
# From https://doi.org/10.1101/2021.10.14.21264988 and refererences therein
symptomatic_vl_frequencies = LogCustomKernel(
np.array((2.46032, 2.67431, 2.85434, 3.06155, 3.25856, 3.47256, 3.66957, 3.85979, 4.09927, 4.27081,
4.47631, 4.66653, 4.87204, 5.10302, 5.27456, 5.46478, 5.6533, 5.88428, 6.07281, 6.30549,
@ -105,7 +105,7 @@ viable_to_RNA_ratio_distribution = Uniform(0.15, 0.45)
# From discussion with virologists
infectious_dose_distribution = Uniform(10., 100.)
# From CERN-OPEN-2021-04 and refererences therein
# From https://doi.org/10.1101/2021.10.14.21264988 and refererences therein
virus_distributions = {
'SARS_CoV_2': mc.SARSCoV2(
viral_load_in_sputum=symptomatic_vl_frequencies,

4
cara/tests/conftest.py
View file

@ -20,9 +20,11 @@ def baseline_model():
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Light activity'],
known_individual_emission_rate=970 * 50,
host_immunity=0.,
# superspreading event, where ejection factor is fixed based
# on Miller et al. (2020) - 50 represents the infectious dose.
),
evaporation_factor=0.3,
)
return model
@ -36,8 +38,8 @@ def baseline_exposure_model(baseline_model):
presence=baseline_model.infected.presence,
activity=baseline_model.infected.activity,
mask=baseline_model.infected.mask,
host_immunity=0.,
),
fraction_deposited=1.,
)

7
cara/tests/models/test_concentration_model.py
View file

@ -43,10 +43,12 @@ def test_concentration_model_vectorisation(override_params):
viral_load_in_sputum=defaults['viral_load_in_sputum'],
infectious_dose=50.,
viable_to_RNA_ratio = 0.5,
transmissibility_factor=1.0,
),
expiration=models._ExpirationBase.types['Breathing'],
host_immunity=0.,
)
),
evaporation_factor=0.3,
)
concentrations = c_model.concentration(10)
assert isinstance(concentrations, np.ndarray)
@ -67,7 +69,8 @@ def simple_conc_model():
virus=models.Virus.types['SARS_CoV_2'],
expiration=models.Expiration.types['Breathing'],
host_immunity=0.,
)
),
evaporation_factor=0.3,
)

51
cara/tests/models/test_exposure_model.py
View file

@ -17,7 +17,7 @@ class KnownNormedconcentration(models.ConcentrationModel):
which therefore doesn't need other components. Useful for testing.
"""
normed_concentration_function: typing.Callable
normed_concentration_function: typing.Callable = lambda x: 0
def infectious_virus_removal_rate(self, time: float) -> models._VectorisedFloat:
# very large decay constant -> same as constant concentration
@ -41,17 +41,17 @@ populations = [
# A simple scalar population.
models.Population(
10, halftime, models.Mask.types['Type I'],
models.Activity.types['Standing'],
models.Activity.types['Standing'], host_immunity=0.,
),
# A population with some array component for η_inhale.
models.Population(
10, halftime, models.Mask(np.array([0.3, 0.35])),
models.Activity.types['Standing'],
models.Activity.types['Standing'], host_immunity=0.
),
# A population with some array component for inhalation_rate.
models.Population(
10, halftime, models.Mask.types['Type I'],
models.Activity(np.array([0.51, 0.57]), 0.57),
models.Activity(np.array([0.51, 0.57]), 0.57), host_immunity=0.
),
]
@ -64,35 +64,35 @@ def known_concentrations(func):
presence=halftime,
mask=models.Mask.types['Type I'],
activity=models.Activity.types['Standing'],
virus=models.Virus.types['SARS_CoV_2_B117'],
virus=models.Virus.types['SARS_CoV_2_ALPHA'],
expiration=models.Expiration.types['Speaking'],
host_immunity=0.,
)
normed_func = lambda x: func(x) / dummy_infected_population.emission_rate_when_present()
return KnownNormedconcentration(dummy_room, dummy_ventilation,
dummy_infected_population, normed_func)
dummy_infected_population, 0.3, normed_func)
@pytest.mark.parametrize(
"population, cm, f_dep, expected_exposure, expected_probability", [
[populations[1], known_concentrations(lambda t: 36.), 1.,
np.array([432, 432]), np.array([99.6803184113, 99.5181053773])],
"population, cm, expected_exposure, expected_probability", [
[populations[1], known_concentrations(lambda t: 36.),
np.array([432, 432]), np.array([67.9503762594, 65.2366759251])],
[populations[2], known_concentrations(lambda t: 36.), 1.,
np.array([432, 432]), np.array([97.4574432074, 98.3493482895])],
[populations[2], known_concentrations(lambda t: 36.),
np.array([432, 432]), np.array([51.6749232285, 55.6374622042])],
[populations[0], known_concentrations(lambda t: np.array([36., 72.])), 1.,
np.array([432, 864]), np.array([98.3493482895, 99.9727534893])],
[populations[0], known_concentrations(lambda t: np.array([36., 72.])),
np.array([432, 864]), np.array([55.6374622042, 80.3196524031])],
[populations[1], known_concentrations(lambda t: np.array([36., 72.])), 1.,
np.array([432, 864]), np.array([99.6803184113, 99.9976777757])],
[populations[1], known_concentrations(lambda t: np.array([36., 72.])),
np.array([432, 864]), np.array([67.9503762594, 87.9151129926])],
[populations[0], known_concentrations(lambda t: 72.), np.array([0.5, 1.]),
864, np.array([98.3493482895, 99.9727534893])],
[populations[2], known_concentrations(lambda t: np.array([36., 72.])),
np.array([432, 864]), np.array([51.6749232285, 80.3196524031])],
])
def test_exposure_model_ndarray(population, cm, f_dep,
def test_exposure_model_ndarray(population, cm,
expected_exposure, expected_probability):
model = ExposureModel(cm, population, fraction_deposited=f_dep)
model = ExposureModel(cm, population)
np.testing.assert_almost_equal(
model.exposure(), expected_exposure
)
@ -154,7 +154,9 @@ def conc_model():
known_individual_emission_rate=970 * 50,
# superspreading event, where ejection factor is fixed based
# on Miller et al. (2020) - 50 represents the infectious dose.
)
host_immunity=0.,
),
evaporation_factor=0.3,
)
@ -176,9 +178,9 @@ def test_exposure_model_integral_accuracy(exposed_time_interval,
presence_interval = models.SpecificInterval((exposed_time_interval,))
population = models.Population(
10, presence_interval, models.Mask.types['Type I'],
models.Activity.types['Standing'],
models.Activity.types['Standing'], 0.,
)
model = ExposureModel(conc_model, population, fraction_deposited=1.)
model = ExposureModel(conc_model, population)
np.testing.assert_allclose(model.exposure(), expected_exposure)
@ -192,6 +194,7 @@ def test_infectious_dose_vectorisation():
viral_load_in_sputum=1e9,
infectious_dose=np.array([50, 20, 30]),
viable_to_RNA_ratio = 0.5,
transmissibility_factor=1.0,
),
expiration=models.Expiration.types['Speaking'],
host_immunity=0.,
@ -202,9 +205,9 @@ def test_infectious_dose_vectorisation():
presence_interval = models.SpecificInterval(((0., 1.),))
population = models.Population(
10, presence_interval, models.Mask.types['Type I'],
models.Activity.types['Standing'],
models.Activity.types['Standing'], 0.,
)
model = ExposureModel(cm, population, fraction_deposited=1.0)
model = ExposureModel(cm, population) #, fraction_deposited=1.0
inf_probability = model.infection_probability()
assert isinstance(inf_probability, np.ndarray)
assert inf_probability.shape == (3, )

1
cara/tests/test_infected_population.py
View file

@ -33,6 +33,7 @@ def test_infected_population_vectorisation(override_params):
viral_load_in_sputum=defaults['viral_load_in_sputum'],
infectious_dose=50.,
viable_to_RNA_ratio = 0.5,
transmissibility_factor=1.0,
),
expiration=cara.models._ExpirationBase.types['Breathing'],
host_immunity=0.,

12
cara/tests/test_known_quantities.py
View file

@ -73,9 +73,11 @@ def build_model(interval_duration):
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Light activity'],
known_individual_emission_rate=970 * 50,
host_immunity=0.,
# superspreading event, where ejection factor is fixed based
# on Miller et al. (2020) - 50 represents the infectious dose.
),
evaporation_factor=0.3,
)
return model
@ -93,7 +95,7 @@ def test_r0(baseline_exposure_model):
# expected r0 was computed with a trapezoidal integration, using
# a mesh of 100'000 pts per exposed presence interval.
r0 = baseline_exposure_model.reproduction_number()
npt.assert_allclose(r0, 972.880852)
npt.assert_allclose(r0, 776.9419902161412)
def test_periodic_window(baseline_periodic_window, baseline_room):
@ -235,7 +237,9 @@ def build_hourly_dependent_model(
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Light activity'],
known_individual_emission_rate=970 * 50,
host_immunity=0,
),
evaporation_factor=0.3,
)
return model
@ -256,7 +260,9 @@ def build_constant_temp_model(outside_temp, intervals_open=((7.5, 8.5),)):
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Light activity'],
known_individual_emission_rate=970 * 50,
host_immunity=0.,
),
evaporation_factor=0.3,
)
return model
@ -284,7 +290,9 @@ def build_hourly_dependent_model_multipleventilation(month, intervals_open=((7.5
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Light activity'],
known_individual_emission_rate=970 * 50,
host_immunity=0.,
),
evaporation_factor=0.3,
)
return model
@ -368,8 +376,8 @@ def build_exposure_model(concentration_model):
presence=infected.presence,
activity=infected.activity,
mask=infected.mask,
host_immunity=0.,
),
fraction_deposited=1.,
)

2
cara/tests/test_monte_carlo.py
View file

@ -56,6 +56,7 @@ def baseline_mc_model() -> cara.monte_carlo.ConcentrationModel:
expiration=cara.models.Expiration.types['Breathing'],
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc_model
@ -69,6 +70,7 @@ def baseline_mc_exposure_model(baseline_mc_model) -> cara.monte_carlo.ExposureMo
presence=baseline_mc_model.infected.presence,
activity=baseline_mc_model.infected.activity,
mask=baseline_mc_model.infected.mask,
host_immunity=0.,
)
)

248
cara/tests/test_monte_carlo_full_models.py
View file

@ -4,8 +4,7 @@ import pytest
import cara.monte_carlo as mc
from cara import models,data
from cara.monte_carlo.data import activity_distributions, virus_distributions
from cara.monte_carlo.data import expiration_distribution, expiration_distributions
from cara.monte_carlo.data import activity_distributions, virus_distributions, expiration_distributions, infectious_dose_distribution, viable_to_RNA_ratio_distribution
from cara.apps.calculator.model_generator import build_expiration
# TODO: seed better the random number generators
@ -19,85 +18,82 @@ TOLERANCE = 0.05
@pytest.fixture
def shared_office_mc():
"""
Corresponds to the 1st line of Table 5 in CERN-OPEN-2021-04, but
speaking 30% of the time (instead of 1/3)
Corresponds to the 1st line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=50, humidity=0.3),
ventilation=models.MultipleVentilation(
(
room=models.Room(volume=50, humidity=0.5),
ventilation = models.MultipleVentilation(
ventilations=(
models.SlidingWindow(
active=models.PeriodicInterval(period=120, duration=10),
inside_temp=models.PiecewiseConstant((0., 24.), (293,)),
outside_temp=models.PiecewiseConstant((0., 24.), (283,)),
window_height=1.6, opening_length=0.6,
active=models.PeriodicInterval(period=120, duration=120),
inside_temp=models.PiecewiseConstant((0., 24.), (298,)),
outside_temp=data.GenevaTemperatures['Jun'],
window_height=1.6,
opening_length=0.2,
),
models.AirChange(
active=models.SpecificInterval(((0., 24.), )),
air_exch=0.25,
),
),
models.AirChange(active=models.PeriodicInterval(period=120, duration=120), air_exch=0.25),
)
),
infected=mc.InfectedPopulation(
number=1,
virus=virus_distributions['SARS_CoV_2_ALPHA'],
presence=mc.SpecificInterval(((0., 2.), (2.1, 4.), (5., 7.), (7.1, 9.))),
mask=models.Mask(η_inhale=0.3),
presence=mc.SpecificInterval(present_times = ((0, 3.5), (4.5, 9))),
virus=virus_distributions['SARS_CoV_2_DELTA'],
mask=models.Mask.types['No mask'],
activity=activity_distributions['Seated'],
expiration=build_expiration({'Speaking': 0.3, 'Breathing': 0.7}),
expiration=build_expiration({'Speaking': 0.33, 'Breathing': 0.67}),
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
exposed=mc.Population(
number=3,
presence=concentration_mc.infected.presence,
activity=models.Activity.types['Seated'],
mask=concentration_mc.infected.mask,
presence=mc.SpecificInterval(present_times = ((0, 3.5), (4.5, 9))),
activity=activity_distributions['Seated'],
mask=models.Mask.types['No mask'],
host_immunity=0.,
),
)
)
@pytest.fixture
def classroom_mc():
"""
Corresponds to the 2nd line of Table 5 in CERN-OPEN-2021-04
Corresponds to the 2nd line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=160, humidity=0.3),
ventilation=models.MultipleVentilation(
(
ventilation = models.MultipleVentilation(
ventilations=(
models.SlidingWindow(
active=models.PeriodicInterval(period=120, duration=10),
active=models.PeriodicInterval(period=120, duration=120),
inside_temp=models.PiecewiseConstant((0., 24.), (293,)),
outside_temp=models.PiecewiseConstant((0., 24.), (283,)),
window_height=1.6, opening_length=0.6,
outside_temp=data.TorontoTemperatures['Dec'],
window_height=1.6,
opening_length=0.2,
),
models.AirChange(
active=models.SpecificInterval(((0., 24.),)),
air_exch=0.25,
),
),
models.AirChange(active=models.PeriodicInterval(period=120, duration=120), air_exch=0.25),
)
),
infected=mc.InfectedPopulation(
number=1,
presence=models.SpecificInterval(((0, 2), (2.5, 4), (5, 7), (7.5, 9))),
virus=virus_distributions['SARS_CoV_2_ALPHA'],
presence=mc.SpecificInterval(((0., 2.), (2.5, 4.), (5., 7.), (7.5, 9.))),
mask=models.Mask.types['No mask'],
mask=models.Mask.types["No mask"],
activity=activity_distributions['Light activity'],
expiration=expiration_distributions['Speaking'],
expiration=build_expiration('Speaking'),
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
exposed=mc.Population(
number=19,
presence=concentration_mc.infected.presence,
activity=models.Activity.types['Seated'],
mask=concentration_mc.infected.mask,
presence=models.SpecificInterval(((0, 2), (2.5, 4), (5, 7), (7.5, 9))),
activity=activity_distributions['Seated'],
mask=models.Mask.types["No mask"],
host_immunity=0.,
),
)
@ -106,42 +102,118 @@ def classroom_mc():
@pytest.fixture
def ski_cabin_mc():
"""
Corresponds to the 3rd line of Table 5 in CERN-OPEN-2021-04
Corresponds to the 3rd line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=10, humidity=0.5),
ventilation=models.AirChange(
active=models.SpecificInterval(((0., 24.),)),
air_exch=0,
),
room=models.Room(volume=10, humidity=0.3),
ventilation=models.MultipleVentilation(
(models.AirChange(active=models.PeriodicInterval(period=120, duration=120), air_exch=0.0),
models.AirChange(active=models.PeriodicInterval(period=120, duration=120), air_exch=0.25))),
infected=mc.InfectedPopulation(
number=1,
virus=virus_distributions['SARS_CoV_2_ALPHA'],
presence=mc.SpecificInterval(((0., 1/3),)),
mask=models.Mask(η_inhale=0.3),
presence=models.SpecificInterval(((0, 20/60),)),
virus=virus_distributions['SARS_CoV_2_DELTA'],
mask=models.Mask.types['No mask'],
activity=activity_distributions['Moderate activity'],
expiration=expiration_distributions['Speaking'],
expiration=build_expiration('Speaking'),
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
exposed=mc.Population(
number=3,
presence=concentration_mc.infected.presence,
activity=models.Activity.types['Moderate activity'],
mask=concentration_mc.infected.mask,
presence=models.SpecificInterval(((0, 20/60),)),
activity=activity_distributions['Moderate activity'],
mask=models.Mask.types['No mask'],
host_immunity=0.,
),
)
@pytest.fixture
def skagit_chorale_mc():
"""
Corresponds to the 4th line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988,
assuming viral is 10**9 instead of a LogCustomKernel distribution.
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=810, humidity=0.5),
ventilation=models.AirChange(
active=models.PeriodicInterval(period=120, duration=120),
air_exch=0.7),
infected=mc.InfectedPopulation(
number=1,
presence=models.SpecificInterval(((0, 2.5), )),
virus=mc.SARSCoV2(
viral_load_in_sputum=10**9,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=1.,
),
mask=models.Mask.types['No mask'],
activity=activity_distributions['Moderate activity'],
expiration=build_expiration('Shouting'),
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
exposed=mc.Population(
number=60,
presence=models.SpecificInterval(((0, 2.5), )),
activity=activity_distributions['Moderate activity'],
mask=models.Mask.types['No mask'],
host_immunity=0.,
),
)
@pytest.fixture
def bus_ride_mc():
"""
Corresponds to the 5th line of Table 4 in https://doi.org/10.1101/2021.10.14.21264988,
assuming viral is 5*10**8 instead of a LogCustomKernel distribution.
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=45, humidity=0.5),
ventilation=models.AirChange(
active=models.PeriodicInterval(period=120, duration=120),
air_exch=1.25),
infected=mc.InfectedPopulation(
number=1,
presence=models.SpecificInterval(((0, 1.67), )),
virus=mc.SARSCoV2(
viral_load_in_sputum=5*10**8,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=1.,
),
mask=models.Mask.types['No mask'],
activity=activity_distributions['Seated'],
expiration=build_expiration('Speaking'),
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
exposed=mc.Population(
number=67,
presence=models.SpecificInterval(((0, 1.67), )),
activity=activity_distributions['Seated'],
mask=models.Mask.types['No mask'],
host_immunity=0.,
),
)
@pytest.fixture
def gym_mc():
"""
Corresponds to the 4th line of Table 5 in CERN-OPEN-2021-04,
but there the expected number of cases is wrongly reported as 0.56
while it should be 0.63.
Gym model for testing
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=300, humidity=0.5),
@ -158,6 +230,7 @@ def gym_mc():
expiration=expiration_distributions['Breathing'],
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
@ -174,8 +247,7 @@ def gym_mc():
@pytest.fixture
def waiting_room_mc():
"""
Corresponds to the 5th line of Table 5 in CERN-OPEN-2021-04, but
speaking 30% of the time (instead of 20%)
Waiting room model for testing
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=100, humidity=0.5),
@ -192,6 +264,7 @@ def waiting_room_mc():
expiration=build_expiration({'Speaking': 0.3, 'Breathing': 0.7}),
host_immunity=0.,
),
evaporation_factor=0.3,
)
return mc.ExposureModel(
concentration_model=concentration_mc,
@ -205,50 +278,16 @@ def waiting_room_mc():
)
@pytest.fixture
def skagit_chorale_mc():
"""
Corresponds to the 6th line of Table 5 in CERN-OPEN-2021-04, but
infection probability should be 29.8% instead of 32%, and number of
new cases 17.9 instead of 21.
"""
concentration_mc = mc.ConcentrationModel(
room=models.Room(volume=810, humidity=0.5),
ventilation=models.AirChange(
active=models.SpecificInterval(((0,24),)),
air_exch=0.7,
),
infected=mc.InfectedPopulation(
number=1,
virus=virus_distributions['SARS_CoV_2'],
presence=mc.SpecificInterval(((0., 2.5),)),
mask=models.Mask.types["No mask"],
activity=activity_distributions['Light activity'],
expiration=expiration_distribution((5., 5., 5.)),
host_immunity=0.,
),
)
return mc.ExposureModel(
concentration_model=concentration_mc,
exposed=mc.Population(
number=60,
presence=concentration_mc.infected.presence,
activity=models.Activity.types['Moderate activity'],
mask=concentration_mc.infected.mask,
host_immunity=0.,
),
)
@pytest.mark.parametrize(
"mc_model, expected_pi, expected_new_cases, expected_dose, expected_ER",
[
["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, 1968],
["waiting_room_mc", 9.72, 1.36, 34.26, 3534],
["skagit_chorale_mc",29.9, 17.9, 190.0, 141400],
["shared_office_mc", 6.03, 0.18, 24.55, 809],
["classroom_mc", 10.0, 2.0, 79.98, 5624],
["ski_cabin_mc", 17.0, 0.5, 40.25, 7966],
["skagit_chorale_mc",70, 42.5, 241.28, 190422],
["bus_ride_mc", 12.0, 8.0, 63.79, 5419],
["gym_mc", 0.45, 0.13, 0.4852, 1145],
["waiting_room_mc", 1.59, 0.22, 7.23, 737],
]
)
def test_report_models(mc_model, expected_pi, expected_new_cases,
@ -269,10 +308,10 @@ def test_report_models(mc_model, expected_pi, expected_new_cases,
@pytest.mark.parametrize(
"mask_type, month, expected_pi, expected_dose, expected_ER",
[
["No mask", "Jul", 30.0, 405.84, 3894],
["Type I", "Jul", 10.2, 73.38, 702],
["FFP2", "Jul", 4.0, 73.38, 702],
["Type I", "Feb", 4.25, 21.42, 702],
["No mask", "Jul", 10.02, 84.54, 809],
["Type I", "Jul", 1.7, 15.64, 149],
["FFP2", "Jul", 0.51, 15.64, 149],
["Type I", "Feb", 0.57, 4.59, 149],
],
)
def test_small_shared_office_Geneva(mask_type, month, expected_pi,
@ -302,6 +341,7 @@ def test_small_shared_office_Geneva(mask_type, month, expected_pi,
expiration=build_expiration({'Speaking': 0.33, 'Breathing': 0.67}),
host_immunity=0.,
),
evaporation_factor=0.3,
)
exposure_mc = mc.ExposureModel(
concentration_model=concentration_mc,

6
cara/tests/test_predefined_distributions.py
View file

@ -8,7 +8,7 @@ from cara.monte_carlo.data import activity_distributions, virus_distributions
np.random.seed(2000)
# mean & std deviations from CERN-OPEN-2021-04 (Table 4)
# mean & std deviations from https://doi.org/10.1101/2021.10.14.21264988 (Table 3)
# NOTE: a mistake was corrected for the std deviation of the
# "Moderate exercise" case (0.37 in the report, but should be 0.34)
@pytest.mark.parametrize(
@ -30,8 +30,8 @@ def test_activity_distributions(distribution, mean, std):
npt.assert_allclose(activity.inhalation_rate.std(), std, atol=0.01)
# mean & std deviations from CERN-OPEN-2021-04 (Table 4) - with a refined
# precision on the values
# mean & std deviations from https://doi.org/10.1101/2021.10.14.21264988 (Table 3)
# - with a refined precision on the values
@pytest.mark.parametrize(
"distribution, mean, std",[
['SARS_CoV_2', 6.59, 1.74],