From 33fb4b3cac4f58011b880d22869ffcfe9c9f0b12 Mon Sep 17 00:00:00 2001 From: Phil Elson Date: Mon, 28 Feb 2022 14:04:31 +0100 Subject: [PATCH] Fix the testing to handle the fact that we have monte carlo components, therefore we need to build_model. --- cara/apps/calculator/model_generator.py | 4 +- cara/models.py | 6 +-- cara/monte_carlo/data.py | 41 ++++++++++++------ cara/monte_carlo/models.py | 5 ++- cara/tests/models/test_short_range_model.py | 47 ++++++++++++--------- 5 files changed, 66 insertions(+), 37 deletions(-) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index 19596e0a..f32de6cd 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -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. diff --git a/cara/models.py b/cara/models.py index 7dfae9e7..2032aff2 100644 --- a/cara/models.py +++ b/cara/models.py @@ -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: """ diff --git a/cara/monte_carlo/data.py b/cara/monte_carlo/data.py index 832e1137..a81ea986 100644 --- a/cara/monte_carlo/data.py +++ b/cara/monte_carlo/data.py @@ -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() } diff --git a/cara/monte_carlo/models.py b/cara/monte_carlo/models.py index 6f09023e..f47d3416 100644 --- a/cara/monte_carlo/models.py +++ b/cara/monte_carlo/models.py @@ -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__: diff --git a/cara/tests/models/test_short_range_model.py b/cara/tests/models/test_short_range_model.py index 5198d290..1925177c 100644 --- a/cara/tests/models/test_short_range_model.py +++ b/cara/tests/models/test_short_range_model.py @@ -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 - ) \ No newline at end of file + )