From cbb3846eb7df1971669a45532765efd4449c49ba Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Wed, 1 Sep 2021 14:59:13 +0200 Subject: [PATCH 1/4] Added uniform distribution for masks --- cara/apps/calculator/model_generator.py | 7 +++++-- cara/monte_carlo/data.py | 9 ++++++++- cara/monte_carlo/sampleable.py | 12 ++++++++++++ 3 files changed, 25 insertions(+), 3 deletions(-) diff --git a/cara/apps/calculator/model_generator.py b/cara/apps/calculator/model_generator.py index a4b139e7..cafc8077 100644 --- a/cara/apps/calculator/model_generator.py +++ b/cara/apps/calculator/model_generator.py @@ -10,7 +10,7 @@ from cara import models from cara import data import cara.monte_carlo as mc from .. import calculator -from cara.monte_carlo.data import activity_distributions, virus_distributions +from cara.monte_carlo.data import activity_distributions, virus_distributions, mask_distributions LOG = logging.getLogger(__name__) @@ -301,7 +301,10 @@ class FormData: def mask(self) -> models.Mask: # Initializes the mask type if mask wearing is "continuous", otherwise instantiates the mask attribute as # the "No mask"-mask - mask = models.Mask.types[self.mask_type if self.mask_wearing_option == "mask_on" else 'No mask'] + if self.mask_wearing_option == 'mask_on': + mask = mask_distributions[self.mask_type] + else: + mask = models.Mask.types['No mask'] return mask def infected_population(self) -> mc.InfectedPopulation: diff --git a/cara/monte_carlo/data.py b/cara/monte_carlo/data.py index 613b830c..250bd6e6 100644 --- a/cara/monte_carlo/data.py +++ b/cara/monte_carlo/data.py @@ -1,7 +1,7 @@ import numpy as np import cara.monte_carlo as mc -from cara.monte_carlo.sampleable import Normal,LogNormal,LogCustomKernel +from cara.monte_carlo.sampleable import Normal,LogNormal,LogCustomKernel, Uniform # From CERN-OPEN-2021-04 and refererences therein @@ -58,3 +58,10 @@ virus_distributions = { infectious_dose=60/1.6, ), } + + +# From CERN-OPEN-2021-04 and refererences therein +mask_distributions = { + 'Type I': mc.Mask(Uniform(0.25, 0.80)), + 'FFP2': mc.Mask(Uniform(0.83, 0.91)), +} \ No newline at end of file diff --git a/cara/monte_carlo/sampleable.py b/cara/monte_carlo/sampleable.py index 4333c93f..27907e49 100644 --- a/cara/monte_carlo/sampleable.py +++ b/cara/monte_carlo/sampleable.py @@ -28,6 +28,18 @@ class Normal(SampleableDistribution): return np.random.normal(self.mean, self.standard_deviation, size=size) +class Uniform(SampleableDistribution): + """ + Defines a continuous uniform distribution + """ + def __init__(self, low: float, high: float): + self.low = low + self.high = high + + def generate_samples(self, size: int) -> float_array_size_n: + return np.random.uniform(self.low, self.high, size=size) + + class LogNormal(SampleableDistribution): """ Defines a lognormal distribution (i.e. Gaussian distribution vs. the From bdc51d81b8e3315010e7bd9c136ff57957bd2edd Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Thu, 2 Sep 2021 12:11:22 +0200 Subject: [PATCH 2/4] Changed relevant references for mask distributions --- cara/monte_carlo/data.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/cara/monte_carlo/data.py b/cara/monte_carlo/data.py index 250bd6e6..4415efbb 100644 --- a/cara/monte_carlo/data.py +++ b/cara/monte_carlo/data.py @@ -60,7 +60,10 @@ virus_distributions = { } -# From CERN-OPEN-2021-04 and refererences therein +# From: +# https://doi.org/10.1080/02786826.2021.1890687 +# https://doi.org/10.1016/j.jhin.2013.02.007 +# https://doi.org/10.4209/aaqr.2020.08.0531 mask_distributions = { 'Type I': mc.Mask(Uniform(0.25, 0.80)), 'FFP2': mc.Mask(Uniform(0.83, 0.91)), From ab153ca7812fdceda7befb5a5de0bd56663c8f9d Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Thu, 2 Sep 2021 15:23:49 +0200 Subject: [PATCH 3/4] added inward mask efficiency references --- cara/apps/templates/common_text.md.j2 | 3 +++ 1 file changed, 3 insertions(+) diff --git a/cara/apps/templates/common_text.md.j2 b/cara/apps/templates/common_text.md.j2 index 63512dee..ac8829d7 100644 --- a/cara/apps/templates/common_text.md.j2 +++ b/cara/apps/templates/common_text.md.j2 @@ -71,3 +71,6 @@ We wish to thank CERN’s HSE Unit, Beams Department, Experimental Physics Depar [54] Leung, N.H.L et al. Respiratory virus shedding in exhaled breath and efficacy of face masks. Nat Med (2020). 10.1038/s41591-020-0843-2.
[55] Asadi, S., Cappa, C.D., Barreda, S. et al. Efficacy of masks and face coverings in controlling outward aerosol particle emission from expiratory activities. Sci Rep 10, 15665 (2020). https://doi.org/10.1038/s41598-020-72798-7.
[56] Endo A, Abbott S et al. Estimating the overdispersion in COVID-19 transmission using outbreak sizes outside China [version 3; peer review: 2 approved]. Wellcome Open Res 2020, 5:67. doi:10.12688/wellcomeopenres.15842.3.
+[57] Jin Pan, Charbel Harb, Weinan Leng & Linsey C. Marr (2021) Inward and outward effectiveness of cloth masks, a surgical mask, and a face shield, Aerosol Science and Technology, 55:6, 718-733, doi: 10.1080/02786826.2021.1890687.
+[58] C. Makison Booth, M. Clayton, B. Crook, J.M. Gawn, Effectiveness of surgical masks against influenza bioaerosols, Journal of Hospital Infection, Volume 84, Issue 1, 2013, Pages 22-26, https://doi.org/10.1016/j.jhin.2013.02.007.
+[59] Riediker, M., Monn, C. (2021). Simulation of SARS-CoV-2 Aerosol Emissions in the Infected Population and Resulting Airborne Exposures in Different Indoor Scenarios. Aerosol Air Qual. Res. 21, 200531. https://doi.org/10.4209/aaqr.2020.08.0531.
\ No newline at end of file From 250b7c0582ecedf58b3e90eb4d7b3f752dcbdd04 Mon Sep 17 00:00:00 2001 From: Luis Aleixo Date: Tue, 7 Sep 2021 06:43:22 +0200 Subject: [PATCH 4/4] Plot alternative scenarios in report using D3.js, speeding up the response time and improving user interactivity --- cara/apps/calculator/report_generator.py | 77 +----- cara/apps/calculator/static/js/report.js | 226 ++++++++++++++---- .../templates/base/calculator.report.html.j2 | 14 +- 3 files changed, 192 insertions(+), 125 deletions(-) diff --git a/cara/apps/calculator/report_generator.py b/cara/apps/calculator/report_generator.py index d1bff9dc..773481f3 100644 --- a/cara/apps/calculator/report_generator.py +++ b/cara/apps/calculator/report_generator.py @@ -9,11 +9,9 @@ import zlib import loky import jinja2 -import matplotlib -matplotlib.use('agg') -import matplotlib.pyplot as plt import numpy as np import qrcode +import json from cara import models from ... import monte_carlo as mc @@ -162,50 +160,13 @@ def _img2bytes(figure): return img_data -def _figure2bytes(figure): - # Draw the image - img_data = io.BytesIO() - figure.savefig(img_data, format='png', bbox_inches="tight", transparent=True) - return img_data - - def img2base64(img_data) -> str: - plt.close() img_data.seek(0) pic_hash = base64.b64encode(img_data.read()).decode('ascii') # A src suitable for a tag such as f'. return f'data:image/png;base64,{pic_hash}' -def plot(times, concentrations, model: models.ExposureModel): - fig = plt.figure() - ax = fig.add_subplot(1, 1, 1) - datetimes = [datetime(1970, 1, 1) + timedelta(hours=time) for time in times] - ax.plot(datetimes, concentrations, lw=2, color='#1f77b4', label='Mean concentration') - ax.spines['right'].set_visible(False) - ax.spines['top'].set_visible(False) - - ax.set_xlabel('Time of day') - ax.set_ylabel('Mean concentration ($virions/m^{3}$)') - ax.set_title('Mean concentration of virions') - ax.xaxis.set_major_formatter(matplotlib.dates.DateFormatter("%H:%M")) - - # Plot presence of exposed person - for i, (presence_start, presence_finish) in enumerate(model.exposed.presence.boundaries()): - plt.fill_between( - datetimes, concentrations, 0, - where=(np.array(times) > presence_start) & (np.array(times) < presence_finish), - color="#1f77b4", alpha=0.1, - label="Presence of exposed person(s)" if i == 0 else "" - ) - - # Place a legend outside of the axes itself. - ax.legend(bbox_to_anchor=(1.05, 1), loc='upper left') - ax.set_ylim(0) - - return fig - - def minutes_to_time(minutes: int) -> str: minute_string = str(minutes % 60) minute_string = "0" * (2 - len(minute_string)) + minute_string @@ -281,39 +242,7 @@ def manufacture_alternative_scenarios(form: FormData) -> typing.Dict[str, mc.Exp return scenarios -def comparison_plot(scenarios: typing.Dict[str, dict], sample_times: typing.List[float]): - fig = plt.figure() - ax = fig.add_subplot(1, 1, 1) - - dash_styled_scenarios = [ - 'Base scenario with FFP2 masks', - 'Base scenario with HEPA filter', - 'Base scenario with HEPA and FFP2 masks', - ] - - sample_dts = [datetime(1970, 1, 1) + timedelta(hours=time) for time in sample_times] - for name, statistics in scenarios.items(): - concentrations = statistics['concentrations'] - - if name in dash_styled_scenarios: - ax.plot(sample_dts, concentrations, label=name, linestyle='--') - else: - ax.plot(sample_dts, concentrations, label=name, linestyle='-', alpha=0.5) - - # Place a legend outside of the axes itself. - ax.legend(bbox_to_anchor=(1.05, 1), loc='upper left') - ax.spines['right'].set_visible(False) - ax.spines['top'].set_visible(False) - ax.xaxis.set_major_formatter(matplotlib.dates.DateFormatter("%H:%M")) - - ax.set_xlabel('Time of day') - ax.set_ylabel('Mean concentration ($virions/m^{3}$)') - ax.set_title('Mean concentration of virions') - - return fig - - -def scenario_statistics(mc_model: mc.ExposureModel, sample_times: typing.List[float]): +def scenario_statistics(mc_model: mc.ExposureModel, sample_times: np.ndarray): model = mc_model.build_model(size=_DEFAULT_MC_SAMPLE_SIZE) return { 'probability_of_infection': np.mean(model.infection_probability()), @@ -342,7 +271,6 @@ def comparison_report( for (name, model), model_stats in zip(scenarios.items(), results): statistics[name] = model_stats return { - 'plot': img2base64(_figure2bytes(comparison_plot(statistics, sample_times))), 'stats': statistics, } @@ -405,6 +333,7 @@ class ReportGenerator: env.filters['minutes_to_time'] = minutes_to_time env.filters['float_format'] = "{0:.2f}".format env.filters['int_format'] = "{:0.0f}".format + env.filters['JSONify'] = json.dumps return env def render(self, context: dict) -> str: diff --git a/cara/apps/calculator/static/js/report.js b/cara/apps/calculator/static/js/report.js index d1105bec..b815da2f 100644 --- a/cara/apps/calculator/static/js/report.js +++ b/cara/apps/calculator/static/js/report.js @@ -25,56 +25,16 @@ function draw_concentration_plot(svg_id, times, concentrations, exposed_presence yAxis = d3.axisLeft(yRange); // Plot tittle. - vis.append('svg:foreignObject') - .attr("background-color", "transparent") - .attr('width', width) - .attr('height', margins.top) - .style('text-align', 'center') - .html('Mean concentration of virions'); + plot_title(vis, width, margins.top, 'Mean concentration of virions'); // Line representing the mean concentration. - var lineFunc = d3.line() - .defined(d => !isNaN(d.concentration)) - .x(d => xTimeRange(d.time)) - .y(d => yRange(d.concentration)) - .curve(d3.curveBasis); + plot_scenario_data(vis, data, xTimeRange, yRange, '#1f77b4'); - vis.append('svg:path') - .attr('d', lineFunc(data)) - .attr('stroke', '#1f77b4') - .attr('stroke-width', 2) - .attr('fill', 'none'); + // X axis. + plot_x_axis(vis, height, width, margins, xAxis, 'Time of day'); - // X axis declaration. - vis.append('svg:g') - .attr('class', 'x axis') - .attr('transform', 'translate(0,' + (height - margins.bottom) + ')') - .call(xAxis); - - // X axis label. - vis.append('text') - .attr('class', 'x label') - .attr('fill', 'black') - .attr('text-anchor', 'middle') - .attr('x', (width + margins.right) / 2) - .attr('y', height * 0.97) - .text('Time of day') - - // Y axis declaration. - vis.append('svg:g') - .attr('class', 'y axis') - .attr('transform', 'translate(' + margins.left + ',0)') - .call(yAxis); - - // Y axis label. - vis.append('svg:text') - .attr('class', 'y label') - .attr('fill', 'black') - .attr('transform', 'rotate(-90, 0,' + height + ')') - .attr('text-anchor', 'middle') - .attr('x', (height + margins.bottom) / 2) - .attr('y', (height + margins.left) * 0.92) - .text('Mean concentration (virions/m³)'); + // Y axis + plot_y_axis(vis, height, width, margins, yAxis, 'Mean concentration (virions/m³)') // Area representing the presence of exposed person(s). exposed_presence_intervals.forEach(b => { @@ -181,4 +141,178 @@ function draw_concentration_plot(svg_id, times, concentrations, exposed_presence focus.select('#tooltip-time').text('x = ' + time_format(d.hour)); focus.select('#tooltip-concentration').text('y = ' + d.concentration.toFixed(2)); } +} + +// Generate the alternative scenarios plot using d3 library. +// 'alternative_scenarios' is a dictionary with all the alternative scenarios +// 'times' is a list of times for all the scenarios +function draw_alternative_scenarios_plot(svg_id, width, height, alternative_scenarios, times) { + // H:M format + var time_format = d3.timeFormat('%H:%M'); + // D3 array of ten categorical colors represented as RGB hexadecimal strings. + var colors = d3.schemeAccent; + + // Variable for the highest concentration for all the scenarios + var highest_concentration = 0. + + var data_for_scenarios = {} + for (scenario in alternative_scenarios) { + scenario_concentrations = alternative_scenarios[scenario].concentrations + + highest_concentration = Math.max(highest_concentration, Math.max(...scenario_concentrations)) + + var data = [] + times.map((time, index) => data.push({ 'time': time, 'hour': new Date().setHours(Math.trunc(time), (time - Math.trunc(time)) * 60), 'concentration': scenario_concentrations[index] })) + + // Add data into lines dictionary + data_for_scenarios[scenario] = data + } + + // We need one scenario to get the time range + var first_scenario = Object.values(data_for_scenarios)[0] + + var vis = d3.select(svg_id), + width = width, + height = height, + margins = { top: 30, right: 20, bottom: 50, left: 50 }, + + // H:M time format for x axis. + xRange = d3.scaleTime().range([margins.left, width - margins.right]).domain([first_scenario[0].hour, first_scenario[first_scenario.length - 1].hour]), + xTimeRange = d3.scaleLinear().range([margins.left, width - margins.right]).domain([times[0], times[times.length - 1]]), + + yRange = d3.scaleLinear().range([height - margins.bottom, margins.top]).domain([0., highest_concentration]), + + xAxis = d3.axisBottom(xRange).tickFormat(d => time_format(d)), + yAxis = d3.axisLeft(yRange); + + // Plot title. + plot_title(vis, width, margins.top, 'Mean concentration of virions'); + + // Line representing the mean concentration for each scenario. + for (const [scenario_name, data] of Object.entries(data_for_scenarios)) { + var scenario_index = Object.keys(data_for_scenarios).indexOf(scenario_name) + + // Line representing the mean concentration. + plot_scenario_data(vis, data, xTimeRange, yRange, colors[scenario_index]) + + // Legend for the plot elements - lines. + var size = 20 * (scenario_index + 1) + vis.append('rect') + .attr('x', width + 20) + .attr('y', margins.top + size) + .attr('width', 20) + .attr('height', 3) + .style('fill', colors[scenario_index]); + + vis.append('text') + .attr('x', width + 3 * 20) + .attr('y', margins.top + size) + .text(scenario_name) + .style('font-size', '15px') + .attr('alignment-baseline', 'central'); + + } + + // X axis. + plot_x_axis(vis, height, width, margins, xAxis, "Time of day"); + + // Y axis declaration. + vis.append('svg:g') + .attr('class', 'y axis') + .attr('transform', 'translate(' + margins.left + ',0)') + .call(yAxis); + + // Y axis label. + vis.append('svg:text') + .attr('class', 'y label') + .attr('fill', 'black') + .attr('transform', 'rotate(-90, 0,' + height + ')') + .attr('text-anchor', 'middle') + .attr('x', (height + margins.bottom) / 2) + .attr('y', (height + margins.left) * 0.92) + .text('Mean concentration (virions/m³)'); + + // Legend bounding box. + vis.append('rect') + .attr('width', 275) + .attr('height', 25 * (Object.keys(data_for_scenarios).length)) + .attr('x', width * 1.005) + .attr('y', margins.top + 5) + .attr('stroke', 'lightgrey') + .attr('stroke-width', '2') + .attr('rx', '5px') + .attr('ry', '5px') + .attr('stroke-linejoin', 'round') + .attr('fill', 'none'); +} + + +// Functions used to build the plots' components + +function plot_title(vis, width, margin_top, title) { + vis.append('svg:foreignObject') + .attr('width', width) + .attr('height', margin_top) + .attr('fill', 'none') + .append('xhtml:div') + .style('text-align', 'center') + .html(title); + + return vis; +} + +function plot_x_axis(vis, height, width, margins, xAxis, label) { + // X axis declaration + vis.append('svg:g') + .attr('class', 'x axis') + .attr('transform', 'translate(0,' + (height - margins.bottom) + ')') + .call(xAxis); + + // X axis label. + vis.append('text') + .attr('class', 'x label') + .attr('fill', 'black') + .attr('text-anchor', 'middle') + .attr('x', (width + margins.right) / 2) + .attr('y', height * 0.97) + .text(label); + + return vis; +} + +function plot_y_axis(vis, height, width, margins, yAxis, label) { + // Y axis declaration. + vis.append('svg:g') + .attr('class', 'y axis') + .attr('transform', 'translate(' + margins.left + ',0)') + .call(yAxis); + + // Y axis label. + vis.append('svg:text') + .attr('class', 'y label') + .attr('fill', 'black') + .attr('transform', 'rotate(-90, 0,' + height + ')') + .attr('text-anchor', 'middle') + .attr('x', (height + margins.bottom) / 2) + .attr('y', (height + margins.left) * 0.92) + .text(label); + + return vis; + +} + +function plot_scenario_data(vis, data, xTimeRange, yRange, line_color) { + var lineFunc = d3.line() + .defined(d => !isNaN(d.concentration)) + .x(d => xTimeRange(d.time)) + .y(d => yRange(d.concentration)) + .curve(d3.curveBasis); + + vis.append('svg:path') + .attr('d', lineFunc(data)) + .attr("stroke", line_color) + .attr('stroke-width', 2) + .attr('fill', 'none'); + + return vis; } \ No newline at end of file diff --git a/cara/apps/calculator/templates/base/calculator.report.html.j2 b/cara/apps/calculator/templates/base/calculator.report.html.j2 index ad04ffae..146478dc 100644 --- a/cara/apps/calculator/templates/base/calculator.report.html.j2 +++ b/cara/apps/calculator/templates/base/calculator.report.html.j2 @@ -88,9 +88,9 @@

* The results are based on the parameters and assumptions published in the CERN Open Report CERN-OPEN-2021-004.

@@ -108,8 +108,12 @@
- - + + {% block report_scenarios_summary_table %}