Short range presence and short range activities as part of the backend model

This commit is contained in:
Luis Aleixo 2022-01-21 17:11:05 +01:00
parent cfa2aa858e
commit 1ea01c1274
5 changed files with 128 additions and 4 deletions

View file

@ -3,6 +3,8 @@ import datetime
import html
import logging
import typing
import ast
import json
import numpy as np
@ -11,7 +13,7 @@ from cara import data
import cara.data.weather
import cara.monte_carlo as mc
from .. import calculator
from cara.monte_carlo.data import activity_distributions, virus_distributions, mask_distributions
from cara.monte_carlo.data import activity_distributions, virus_distributions, mask_distributions, initial_concentrations_mouth
from cara.monte_carlo.data import expiration_distribution, expiration_BLO_factors, expiration_distributions
@ -76,6 +78,8 @@ class FormData:
window_width: float
windows_number: int
window_opening_regime: str
short_range_option: str
short_range_interactions: list
#: The default values for undefined fields. Note that the defaults here
#: and the defaults in the html form must not be contradictory.
@ -127,6 +131,8 @@ class FormData:
'windows_frequency': 0.,
'windows_number': 0,
'window_opening_regime': 'windows_open_permanently',
'short_range_option': False,
'short_range_interactions': [],
}
@classmethod
@ -416,6 +422,8 @@ class FormData:
number=infected_occupants,
virus=virus,
presence=self.infected_present_interval(),
short_range_presence=self.short_range_intervals(),
short_range_activities=self.short_range_activities(),
mask=self.mask(),
activity=activity,
expiration=expiration,
@ -448,6 +456,7 @@ class FormData:
exposed = mc.Population(
number=exposed_occupants,
presence=self.exposed_present_interval(),
short_range_presence=self.short_range_interactions,
activity=activity,
mask=self.mask(),
host_immunity=0.,
@ -623,12 +632,30 @@ class FormData:
breaks=self.infected_lunch_break_times() + self.infected_coffee_break_times(),
)
def short_range_intervals(self) -> models.Interval:
if (self.short_range_interactions):
short_range_intervals = []
for interaction in self.short_range_interactions:
start_time = time_string_to_minutes(interaction['start_time'])
duration = float(interaction['duration'])
short_range_intervals.append(models.SpecificInterval((start_time/60, (start_time + duration)/60)))
return short_range_intervals
else:
return []
def exposed_present_interval(self) -> models.Interval:
return self.present_interval(
self.exposed_start, self.exposed_finish,
breaks=self.exposed_lunch_break_times() + self.exposed_coffee_break_times(),
)
def short_range_activities(self):
if (self.short_range_interactions):
return [interaction['activity'] for interaction in self.short_range_interactions]
else:
return []
def build_expiration(expiration_definition) -> models._ExpirationBase:
if isinstance(expiration_definition, str):
@ -642,6 +669,12 @@ def build_expiration(expiration_definition) -> models._ExpirationBase:
return expiration_distribution(tuple(BLO_factors))
def build_concentration_mouth(short_range_def) -> models._ExpirationBase:
expiration_definition = short_range_def[1]['activity']
if isinstance(expiration_definition, str):
return initial_concentrations_mouth[expiration_definition]
def baseline_raw_form_data():
# Note: This isn't a special "baseline". It can be updated as required.
return {
@ -691,7 +724,8 @@ def baseline_raw_form_data():
'window_type': 'window_sliding',
'window_width': '2',
'windows_number': '1',
'window_opening_regime': 'windows_open_permanently'
'window_opening_regime': 'windows_open_permanently',
'short_range_option': 'short_range_no',
}
@ -736,6 +770,15 @@ def time_minutes_to_string(time: int) -> str:
return "{0:0=2d}".format(int(time/60)) + ":" + "{0:0=2d}".format(time%60)
def string_to_list(l: str) -> list:
if (l != []):
return list(ast.literal_eval(l.replace(""", "\"")))
def list_to_string(s: list) -> str:
return json.dumps(s)
def _safe_int_cast(value) -> int:
if isinstance(value, int):
return value
@ -767,3 +810,6 @@ for _field in dataclasses.fields(FormData):
elif _field.type is bool:
_CAST_RULES_FORM_ARG_TO_NATIVE[_field.name] = lambda v: v == '1'
_CAST_RULES_NATIVE_TO_FORM_ARG[_field.name] = int
elif _field.type is list:
_CAST_RULES_FORM_ARG_TO_NATIVE[_field.name] = string_to_list
_CAST_RULES_NATIVE_TO_FORM_ARG[_field.name] = list_to_string

View file

@ -12,8 +12,9 @@ import jinja2
import numpy as np
from cara import models
from cara.monte_carlo.data import initial_concentrations_mouth, dilution_factor
from ... import monte_carlo as mc
from .model_generator import FormData, _DEFAULT_MC_SAMPLE_SIZE
from .model_generator import FormData, _DEFAULT_MC_SAMPLE_SIZE, build_expiration
from ... import dataclass_utils
@ -96,13 +97,45 @@ def interesting_times(model: models.ExposureModel, approx_n_pts=100) -> typing.L
return nice_times
def short_range_interesting_times(model: models.ExposureModel, times: typing.List[float]) -> typing.List[float]:
short_range_times = []
for period in model.concentration_model.infected.short_range_presence:
start, finish = period.boundaries()
short_range_times = short_range_times + [time for time in times if time > start and time < finish]
return short_range_times
def jet_origin_concentrations(model: models.ExposureModel) -> models._ExpirationBase:
return [initial_concentrations_mouth[activity] for activity in model.concentration_model.infected.short_range_activities]
def short_range_initial_concentrations(model: models.ExposureModel, time: float):
dilution = dilution_factor(np.linspace(0.1, 2., 1000))
jet_origin_initial_concentrations = jet_origin_concentrations(model)
for index, interaction in enumerate(model.concentration_model.infected.short_range_presence):
start, finish = interaction.boundaries()
if start <= time <= finish:
expiration = build_expiration(model.concentration_model.infected.short_range_activities[index])
single_exposure_model = dataclass_utils.nested_replace(
model, {'concentration_model.infected.expiration': expiration}
)
concentration = single_exposure_model.concentration_model._normed_concentration(float(time))
jet_origin_concentration = jet_origin_initial_concentrations[index] * model.concentration_model.infected.virus.viral_load_in_sputum
return concentration + ((1/dilution)*(jet_origin_concentration - concentration))
def calculate_report_data(model: models.ExposureModel):
times = interesting_times(model)
times_short_range = short_range_interesting_times(model, times)
concentrations = [
np.array(model.concentration_model.concentration(float(time))).mean()
for time in times
]
short_range_concentrations = [
np.array(short_range_initial_concentrations(model, float(time))).mean()
for time in times_short_range
]
highest_const = max(concentrations)
prob = np.array(model.infection_probability()).mean()
er = np.array(model.concentration_model.infected.emission_rate_when_present()).mean()
@ -118,6 +151,7 @@ def calculate_report_data(model: models.ExposureModel):
"exposed_presence_intervals": [list(interval) for interval in model.exposed.presence.boundaries()],
"cumulative_doses": list(cumulative_doses),
"concentrations": concentrations,
"short_range_concentrations": short_range_concentrations,
"highest_const": highest_const,
"prob_inf": prob,
"emission_rate": er,

View file

@ -496,6 +496,8 @@ baseline_model = models.ExposureModel(
number=1,
virus=models.Virus.types['SARS_CoV_2'],
presence=models.SpecificInterval(((8., 12.), (13., 17.))),
short_range_presence=[],
short_range_activities=[],
mask=models.Mask.types['No mask'],
activity=models.Activity.types['Seated'],
expiration=models.Expiration.types['Speaking'],
@ -506,6 +508,7 @@ baseline_model = models.ExposureModel(
exposed=models.Population(
number=10,
presence=models.SpecificInterval(((8., 12.), (13., 17.))),
short_range_presence=[],
activity=models.Activity.types['Seated'],
mask=models.Mask.types['No mask'],
host_immunity=0.,

View file

@ -678,6 +678,9 @@ class Population:
#: The times in which the people are in the room.
presence: Interval
#: Short range interactions
short_range_presence: list
#: The kind of mask being worn by the people.
mask: Mask
@ -749,6 +752,9 @@ class EmittingPopulation(_PopulationWithVirus):
class InfectedPopulation(_PopulationWithVirus):
#: The type of expiration that is being emitted whilst doing the activity.
expiration: _ExpirationBase
#: The type of expiractory activities in the short range interactions
short_range_activities: list
@method_cache
def fraction_of_infectious_virus(self) -> _VectorisedFloat:

View file

@ -3,6 +3,7 @@ import typing
import numpy as np
from scipy import special as sp
import scipy.integrate
import cara.monte_carlo as mc
from cara.monte_carlo.sampleable import Normal,LogNormal,LogCustomKernel,CustomKernel,Uniform
@ -51,6 +52,9 @@ class BLOmodel:
np.exp(-(np.log(d) - mu) ** 2 / (2 * sigma ** 2))
for A,cn,mu,sigma in zip(self.BLO_factors, self.cn,
self.mu, self.sigma) )
def volume(self, d):
return(np.pi * d**3) / 6
def integrate(self, dmin, dmax):
"""
@ -162,7 +166,7 @@ def expiration_distribution(BLO_factors):
Returns an Expiration with an aerosol diameter distribution, defined
by the BLO factors (a length-3 tuple).
The total concentration of aerosols is computed by integrating
the distribution between 0.1 and 30 microns - these boundaries are
the distribution between 0.1 and D microns - these boundaries are
an historical choice based on previous implementations of the model
(it limits the influence of the O-mode).
"""
@ -172,6 +176,31 @@ def expiration_distribution(BLO_factors):
BLOmodel(BLO_factors).integrate(0.1, 30.))
def dilution_factor(distance, D=0.02):
u0 = 0.6
tstar = 2.0
Cr1 = 0.18
Cr2 = 0.2
Cx1 = 2.4
# The expired flow rate during the expiration period, m^3/s
Q0 = u0 * np.pi/4*D**2
# Parameters in the jet-like stage
x01 = D/2/Cr1
# Time of virtual origin
t01 = (x01/Cx1)**2 * (Q0*u0)**(-0.5)
# The transition point, m
xstar = Cx1*(Q0*u0)**0.25*(tstar + t01)**0.5 - x01
# Dilution factor at the transition point xstar
Sxstar = 2*Cr1*(xstar+x01)/D
return np.mean(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]))
def initial_concentration_mouth(BLO_factors):
value, error = scipy.integrate.quad(lambda d: BLOmodel(BLO_factors).distribution(d) * BLOmodel(BLO_factors).volume(d), 0.1, 1000)
return value * 1e-6 #result in mL/m^3
expiration_BLO_factors = {
'Breathing': (1., 0., 0.),
'Speaking': (1., 1., 1.),
@ -184,3 +213,9 @@ expiration_distributions = {
exp_type: expiration_distribution(BLO_factors)
for exp_type,BLO_factors in expiration_BLO_factors.items()
}
initial_concentrations_mouth = {
exp_type: initial_concentration_mouth(BLO_factors)
for exp_type,BLO_factors in expiration_BLO_factors.items()
}