diff --git a/cara/apps/calculator/__init__.py b/cara/apps/calculator/__init__.py index 8d4e5e45..d7b15698 100644 --- a/cara/apps/calculator/__init__.py +++ b/cara/apps/calculator/__init__.py @@ -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): diff --git a/cara/monte_carlo/data.py b/cara/monte_carlo/data.py index e70a65de..886037cf 100644 --- a/cara/monte_carlo/data.py +++ b/cara/monte_carlo/data.py @@ -6,7 +6,7 @@ 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) @@ -108,8 +108,10 @@ symptomatic_vl_frequencies = LogCustomKernel( # 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 = weibull_min.pdf(viral_load, c=3.47, scale=7.01) -covid_overal_vl_data = Custom(bounds=(2, 10), function=lambda d: np.interp(d, viral_load, frequencies, right=0., left=0.), max_function=0.16) +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, right=0., left=0.), + max_function=0.2) # Derived from data in doi.org/10.1016/j.ijid.2020.09.025 and diff --git a/cara/monte_carlo/sampleable.py b/cara/monte_carlo/sampleable.py index 27907e49..fab478e0 100644 --- a/cara/monte_carlo/sampleable.py +++ b/cara/monte_carlo/sampleable.py @@ -82,6 +82,30 @@ class Custom(SampleableDistribution): return x +class LogCustom(SampleableDistribution): + """ + Defines a distribution which follows a custom curve vs. the the log (in base 10) + of the random variable. Uses a simple algorithm. This is appropriate for a smooth + distribution function (one should know its maximum). + """ + 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