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

Merge branch 'feature/new_viral_data' into 'master'

Browse source 
Added weibull_min dist to data.py on virus viral_load_in_sputum

See merge request cara/cara!339
This commit is contained in:
Andre Henriques 2022-06-01 08:58:57 +02:00
commit 9ebb521769
5 changed files with 66 additions and 23 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

2
cara/apps/calculator/__init__.py
View file

@ -33,7 +33,7 @@ from .user import AuthenticatedUser, AnonymousUser
# calculator version. If the calculator needs to make breaking changes (e.g. change
# form attributes) then it can also increase its MAJOR version without needing to
# increase the overall CARA version (found at ``cara.__version__``).
__version__ = "4.1.2"
__version__ = "4.2"
class BaseRequestHandler(RequestHandler):

29
cara/monte_carlo/data.py
View file

@ -3,9 +3,10 @@ import typing
import numpy as np
from scipy import special as sp
from scipy.stats import weibull_min
import cara.monte_carlo as mc
from cara.monte_carlo.sampleable import LogNormal,LogCustomKernel,CustomKernel,Uniform, Custom
from cara.monte_carlo.sampleable import LogCustom, LogNormal,LogCustomKernel,CustomKernel,Uniform, Custom
sqrt2pi = np.sqrt(2.*np.pi)
@ -101,47 +102,61 @@ symptomatic_vl_frequencies = LogCustomKernel(
kernel_bandwidth=0.1
)
# Weibull distribution with a shape factor of 3.47 and a scale factor of 7.01.
# From https://elifesciences.org/articles/65774 and first line of the figure in
# https://iiif.elifesciences.org/lax:65774%2Felife-65774-fig4-figsupp3-v2.tif/full/1500,/0/default.jpg
viral_load = np.linspace(weibull_min.ppf(0.01, c=3.47, scale=7.01),
weibull_min.ppf(0.99, c=3.47, scale=7.01), 30)
frequencies_pdf = weibull_min.pdf(viral_load, c=3.47, scale=7.01)
covid_overal_vl_data = LogCustom(bounds=(2, 10),
function=lambda d: np.interp(d, viral_load, frequencies_pdf, left=0., right=0.),
max_function=0.2)
# Derived from data in doi.org/10.1016/j.ijid.2020.09.025 and
# https://iosh.com/media/8432/aerosol-infection-risk-hospital-patient-care-full-report.pdf (page 60)
viable_to_RNA_ratio_distribution = Uniform(0.01, 0.6)
# From discussion with virologists
infectious_dose_distribution = Uniform(10., 100.)
# 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,
viral_load_in_sputum=covid_overal_vl_data,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=1.,
),
'SARS_CoV_2_ALPHA': mc.SARSCoV2(
viral_load_in_sputum=symptomatic_vl_frequencies,
viral_load_in_sputum=covid_overal_vl_data,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=0.78,
),
'SARS_CoV_2_BETA': mc.SARSCoV2(
viral_load_in_sputum=symptomatic_vl_frequencies,
viral_load_in_sputum=covid_overal_vl_data,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=0.8,
),
'SARS_CoV_2_GAMMA': mc.SARSCoV2(
viral_load_in_sputum=symptomatic_vl_frequencies,
viral_load_in_sputum=covid_overal_vl_data,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=0.72,
),
'SARS_CoV_2_DELTA': mc.SARSCoV2(
viral_load_in_sputum=symptomatic_vl_frequencies,
viral_load_in_sputum=covid_overal_vl_data,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=0.51,
),
'SARS_CoV_2_OMICRON': mc.SARSCoV2(
viral_load_in_sputum=symptomatic_vl_frequencies,
viral_load_in_sputum=covid_overal_vl_data,
infectious_dose=infectious_dose_distribution,
viable_to_RNA_ratio=viable_to_RNA_ratio_distribution,
transmissibility_factor=0.2,

30
cara/monte_carlo/sampleable.py
View file

@ -62,7 +62,9 @@ class Custom(SampleableDistribution):
"""
Defines a distribution which follows a custom curve vs. the random
variable. Uses a simple algorithm. This is appropriate for a smooth
distribution function (one should know its maximum).
distribution function.
Note: in max_function, a value slightly above the maximum of the distribution
function should be provided.
"""
def __init__(self, bounds: typing.Tuple[float, float],
function: typing.Callable, max_function: float):
@ -82,6 +84,32 @@ class Custom(SampleableDistribution):
return x
class LogCustom(SampleableDistribution):
"""
Defines a distribution which follows a custom curve vs. the log (in base 10)
of the random variable. Uses a simple algorithm. This is appropriate for a smooth
distribution function.
Note: in max_function, a value slightly above the maximum of the distribution
function should be provided.
"""
def __init__(self, bounds: typing.Tuple[float, float],
function: typing.Callable, max_function: float):
self.bounds = bounds
self.function = function
self.max_function = max_function
def generate_samples(self, size: int) -> float_array_size_n:
fvalue = np.random.uniform(0,self.max_function,size)
x = np.random.uniform(*self.bounds,size)
invalid = np.where(fvalue>self.function(x))[0]
while len(invalid)>0:
fvalue[invalid] = np.random.uniform(0,self.max_function,len(invalid))
x[invalid] = np.random.uniform(*self.bounds,len(invalid))
invalid = np.where(fvalue>self.function(x))[0]
return 10 ** x
class CustomKernel(SampleableDistribution):
"""
Defines a distribution which follows a custom curve vs. the

22
cara/tests/test_monte_carlo_full_models.py
View file

@ -308,17 +308,17 @@ def waiting_room_mc():
)
@retry()
@retry(tries=10)
@pytest.mark.parametrize(
"mc_model, expected_pi, expected_new_cases, expected_dose, expected_ER",
[
["shared_office_mc", 5.55, 0.17, 2.699, 809],
["classroom_mc", 9.58, 1.82, 9.034, 5624],
["ski_cabin_mc", 16.0, 0.47, 17.315, 7966],
["shared_office_mc", 5.38, 0.16, 3.350, 1056],
["classroom_mc", 8.21, 1.56, 11.356, 7416],
["ski_cabin_mc", 12.92, 0.39, 21.796, 10231],
["skagit_chorale_mc",61.01, 36.53, 84.730, 190422],
["bus_ride_mc", 10.59, 7.06, 6.65, 5419],
["gym_mc", 0.43, 0.12, 0.197, 1145],
["waiting_room_mc", 1.34, 0.18, 0.670, 737],
["bus_ride_mc", 10.59, 7.06, 6.650, 5419],
["gym_mc", 0.52, 0.14, 0.249, 1450],
["waiting_room_mc", 1.53, 0.21, 0.844, 929],
]
)
def test_report_models(mc_model, expected_pi, expected_new_cases,
@ -339,10 +339,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", 8.46, 8.113, 809],
["Type I", "Jul", 1.44, 0.727, 149],
["FFP2", "Jul", 0.43, 0.197, 149],
["Type I", "Feb", 0.54, 0.253, 149],
["No mask", "Jul", 7.689, 10.050, 1034.435],
["Type I", "Jul", 1.663, 0.938, 193.52],
["FFP2", "Jul", 0.523, 0.253, 193.52],
["Type I", "Feb", 0.659, 0.325, 193.52],
],
)
def test_small_shared_office_Geneva(mask_type, month, expected_pi,

6
cara/tests/test_predefined_distributions.py
View file

@ -34,11 +34,11 @@ def test_activity_distributions(distribution, mean, std):
# - with a refined precision on the values
@pytest.mark.parametrize(
"distribution, mean, std",[
['SARS_CoV_2', 6.59, 1.74],
['SARS_CoV_2', 6.22, 1.80],
['SARS_CoV_2_ALPHA', 6.59, 1.74],
['SARS_CoV_2_ALPHA', 6.22, 1.80],
['SARS_CoV_2_GAMMA', 6.59, 1.74],
['SARS_CoV_2_GAMMA', 6.22, 1.80],
]
)
def test_viral_load_logdistribution(distribution, mean, std):