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

Fix the testing to handle the fact that we have monte carlo components, therefore we need to build_model.

Browse source
This commit is contained in:
Phil Elson 2022-02-28 14:04:31 +01:00
parent f7be248d14
commit 33fb4b3cac
5 changed files with 66 additions and 37 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

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

@ -665,7 +665,7 @@ class FormData:
if self.short_range_interactions else [])
def build_expiration(expiration_definition) -> models._ExpirationBase:
def build_expiration(expiration_definition) -> mc._ExpirationBase:
if isinstance(expiration_definition, str):
return expiration_distributions[expiration_definition]
elif isinstance(expiration_definition, dict):
@ -675,7 +675,7 @@ def build_expiration(expiration_definition) -> models._ExpirationBase:
for exp_type, weight in expiration_definition.items()
], axis=0)
return expiration_distribution(tuple(BLO_factors))
def baseline_raw_form_data():
# Note: This isn't a special "baseline". It can be updated as required.

6
cara/models.py
View file

@ -1076,13 +1076,13 @@ class ConcentrationModel:
@dataclass(frozen=True)
class ShortRangeModel:
#: Short range interactions
presence: typing.List[SpecificInterval]
presence: typing.Tuple[SpecificInterval, ...]
#: Expiration types
expirations: typing.List[Expiration]
expirations: typing.Tuple[_ExpirationBase, ...]
#: The dilution factors for each of the expiratory activity
dilutions: typing.List[_VectorisedFloat]
dilutions: typing.Tuple[_VectorisedFloat, ...]
def _normed_concentration(self, concentration_model: ConcentrationModel, time: float) -> _VectorisedFloat:
"""

41
cara/monte_carlo/data.py
View file

@ -10,6 +10,7 @@ from cara.monte_carlo.sampleable import Normal,LogNormal,LogCustomKernel,CustomK
sqrt2pi = np.sqrt(2.*np.pi)
sqrt2 = np.sqrt(2.)
@dataclass(frozen=True)
class BLOmodel:
"""
@ -65,7 +66,7 @@ class BLOmodel:
return result
# From https://doi.org/10.1101/2021.10.14.21264988 and refererences therein
# From https://doi.org/10.1101/2021.10.14.21264988 and references therein
activity_distributions = {
'Seated': mc.Activity(LogNormal(-0.6872121723362303, 0.10498338229297108),
LogNormal(-0.6872121723362303, 0.10498338229297108)),
@ -84,7 +85,7 @@ activity_distributions = {
}
# From https://doi.org/10.1101/2021.10.14.21264988 and refererences therein
# From https://doi.org/10.1101/2021.10.14.21264988 and references 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,
@ -157,7 +158,10 @@ mask_distributions = {
}
def expiration_distribution(BLO_factors, d_max=30.):
def expiration_distribution(
BLO_factors: typing.Tuple[float, float, float],
d_max=30.,
) -> mc.Expiration:
"""
Returns an Expiration with an aerosol diameter distribution, defined
by the BLO factors (a length-3 tuple).
@ -166,10 +170,15 @@ def expiration_distribution(BLO_factors, d_max=30.):
an historical choice based on previous implementations of the model
(it limits the influence of the O-mode).
"""
dscan = np.linspace(0.1, d_max ,3000)
return mc.Expiration(CustomKernel(dscan,
BLOmodel(BLO_factors).distribution(dscan),kernel_bandwidth=0.1),
cn=BLOmodel(BLO_factors).integrate(0.1, d_max))
dscan = np.linspace(0.1, d_max, 3000)
return mc.Expiration(
CustomKernel(
dscan,
BLOmodel(BLO_factors).distribution(dscan),
kernel_bandwidth=0.1,
),
cn=BLOmodel(BLO_factors).integrate(0.1, d_max),
)
def dilution_factor(activities, distance, D=0.02):
@ -190,8 +199,16 @@ def dilution_factor(activities, distance, D=0.02):
xstar = Cx1*(Q0*u0)**0.25*(tstar + t01)**0.5 - x01
# Dilution factor at the transition point xstar
Sxstar = 2*Cr1*(xstar+x01)/D
factors.append(np.piecewise(distance, [distance < xstar, distance >= xstar],
[lambda distance : 2*Cr1*(distance + x01)/D, lambda distance : Sxstar*(1 + Cr2*(distance - xstar)/Cr1/(xstar + x01))**3]))
factors.append(
np.piecewise(
distance,
[distance < xstar, distance >= xstar],
[
lambda distance: 2*Cr1*(distance + x01)/D,
lambda distance: Sxstar*(1 + Cr2*(distance - xstar)/Cr1/(xstar + x01))**3
]
)
)
return factors
@ -205,11 +222,11 @@ expiration_BLO_factors = {
expiration_distributions = {
exp_type: expiration_distribution(BLO_factors)
for exp_type,BLO_factors in expiration_BLO_factors.items()
for exp_type, BLO_factors in expiration_BLO_factors.items()
}
short_range_expiration_distributions = {
exp_type: expiration_distribution(BLO_factors, d_max=100).build_model(250000)
for exp_type,BLO_factors in expiration_BLO_factors.items()
exp_type: expiration_distribution(BLO_factors, d_max=100)
for exp_type, BLO_factors in expiration_BLO_factors.items()
}

5
cara/monte_carlo/models.py
View file

@ -37,7 +37,7 @@ class MCModelBase(typing.Generic[_ModelType]):
def build_model(self, size: int) -> _ModelType:
"""
Turn this MCModelBase subclass into a cara.models Model instance
Turn this MCModelBase subclass into a cara.model Model instance
from which you can then run the model.
"""
@ -72,6 +72,9 @@ def _build_mc_model(model: _ModelType) -> typing.Type[MCModelBase[_ModelType]]:
elif new_field.type == typing.Tuple[cara.models._ExpirationBase, ...]:
EB = getattr(sys.modules[__name__], "_ExpirationBase")
field_type = typing.Tuple[typing.Union[cara.models._ExpirationBase, EB], ...]
elif new_field.type == typing.Tuple[cara.models.SpecificInterval, ...]:
SI = getattr(sys.modules[__name__], "SpecificInterval")
field_type = typing.Tuple[typing.Union[cara.models.SpecificInterval, SI], ...]
else:
# Check that we don't need to do anything with this type.
for item in new_field.type.__args__:

47
cara/tests/models/test_short_range_model.py
View file

@ -1,30 +1,29 @@
import typing
from unicodedata import decimal
import numpy as np
import pytest
from cara import models
from cara.models import ShortRangeModel
import cara.monte_carlo.models as mc_models
from cara.apps.calculator.model_generator import build_expiration
from cara.monte_carlo.data import dilution_factor, short_range_expiration_distributions
@pytest.fixture
def concentration_model():
return models.ConcentrationModel(
def concentration_model() -> mc_models.ConcentrationModel:
return mc_models.ConcentrationModel(
room=models.Room(volume=75),
ventilation=models.AirChange(
active=models.SpecificInterval(present_times=((8.5, 12.5), (13.5, 17.5))),
air_exch=30.,
),
infected=models.InfectedPopulation(
infected=mc_models.InfectedPopulation(
number=1,
virus=models.Virus.types['SARS_CoV_2'],
presence=models.SpecificInterval(present_times=((8.5, 12.5), (13.5, 17.5))),
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Light activity'],
expiration=build_expiration({'Speaking': 0.33, 'Breathing': 0.67}).build_model(250000),
expiration=build_expiration({'Speaking': 0.33, 'Breathing': 0.67}),
host_immunity=0.,
),
evaporation_factor=0.3,
@ -35,25 +34,32 @@ activities = ['Breathing', 'Speaking', 'Shouting']
@pytest.fixture
def presences():
return [models.SpecificInterval((10.5, 11.0)),
def presences() -> typing.Tuple[models.SpecificInterval, ...]:
return (
models.SpecificInterval((10.5, 11.0)),
models.SpecificInterval((14.5, 15.0)),
models.SpecificInterval((16.5, 17.5)),]
models.SpecificInterval((16.5, 17.5)),
)
@pytest.fixture
def expirations():
return [short_range_expiration_distributions[activity] for activity in activities]
def expirations() -> typing.Tuple[mc_models.Expiration, ...]:
# Monte carlo expirations! So build model.
return tuple(short_range_expiration_distributions[activity] for activity in activities)
@pytest.fixture
def dilutions():
return dilution_factor(activities=activities,
distance=np.random.uniform(0.5, 1.5, 250000))
return dilution_factor(
activities=activities,
distance=np.random.uniform(0.5, 1.5, 250_000),
)
def test_short_range_model_ndarray(concentration_model, presences, expirations, dilutions):
model = ShortRangeModel(presences, expirations, dilutions)
concentration_model = concentration_model.build_model(250_000)
model = mc_models.ShortRangeModel(presences, expirations, dilutions)
model = model.build_model(250_000)
assert isinstance(model._normed_concentration(concentration_model, 10.75), np.ndarray)
assert isinstance(model.short_range_concentration(concentration_model, 14.75), np.ndarray)
assert isinstance(model.normed_exposure_between_bounds(concentration_model, 16.6, 17.7), np.ndarray)
@ -66,13 +72,16 @@ def test_short_range_model_ndarray(concentration_model, presences, expirations,
# [16.75, 1.973757599365769e-05, 19737.57599365769],
]
)
def test_short_range_model(time, expected_sr_normed_concentration, expected_concentration,
concentration_model, presences, expirations, dilutions):
model = ShortRangeModel(presences, expirations, dilutions)
def test_short_range_model(
time, expected_sr_normed_concentration, expected_concentration,
concentration_model, presences, expirations, dilutions,
):
concentration_model = concentration_model.build_model(250_000)
model = mc_models.ShortRangeModel(presences, expirations, dilutions)
model = model.build_model(250_000)
np.testing.assert_almost_equal(
model._normed_concentration(concentration_model, time).mean(), expected_sr_normed_concentration, decimal=0
)
np.testing.assert_almost_equal(
model.short_range_concentration(concentration_model, time).mean(), expected_concentration, decimal=0
)
)