From f5ee4c3882c38b8a6acb77e910b4bab5738769cf Mon Sep 17 00:00:00 2001
From: Luis Aleixo <luis.aleixo@cern.ch>
Date: Fri, 6 Aug 2021 11:36:23 +0200
Subject: [PATCH] updated Yan data

---
 .gitignore          |  2 ++
 cara/plot_output.py | 47 ++++++++++++++++++++++++++++++++++++++-------
 2 files changed, 42 insertions(+), 7 deletions(-)

diff --git a/.gitignore b/.gitignore
index 38eee2e6..16187d86 100644
--- a/.gitignore
+++ b/.gitignore
@@ -4,6 +4,8 @@ __pycache__
 *.DS_Store
 *.pyc
 
+data.csv
+
 # Editor stuff
 *.swp
 .idea
diff --git a/cara/plot_output.py b/cara/plot_output.py
index c799f582..7a27f4b0 100644
--- a/cara/plot_output.py
+++ b/cara/plot_output.py
@@ -4,7 +4,7 @@ import matplotlib.pyplot as plt
 import csv
 
 import cara.monte_carlo as mc
-from cara import models,data
+from cara import models, data
 from cara.monte_carlo.data import activity_distributions
 from tqdm import tqdm
 
@@ -32,15 +32,15 @@ for vl in tqdm(viral_loads):
             infected=mc.InfectedPopulation(
                 number=1,
                 virus=models.Virus(
-                    viral_load_in_sputum = 10**vl,
-                    infectious_dose = 50.,
+                    viral_load_in_sputum=10**vl,
+                    infectious_dose=50.,
                 ),
                 presence=mc.SpecificInterval(((0, 2),)),
                 mask=models.Mask.types["No mask"],
                 activity=activity_distributions['Seated'],
                 expiration=models.MultipleExpiration(
                     expirations=(models.Expiration.types['Talking'],
-                                    models.Expiration.types['Breathing']),
+                                 models.Expiration.types['Breathing']),
                     weights=(0.3, 0.7)),
             ),
         ),
@@ -63,14 +63,16 @@ with open('data.csv', 'w', newline='') as csvfile:
     thewriter = csv.DictWriter(csvfile, fieldnames=fieldnames)
     thewriter.writeheader()
     for i, vl in enumerate(viral_loads):
-        thewriter.writerow({'viral load' : 10**vl, 'emission rate' : er_means[i]})
+        thewriter.writerow(
+            {'viral load': 10**vl, 'emission rate': er_means[i]})
 
 ax.plot(viral_loads, er_means)
 ax.fill_between(viral_loads, lower_percentiles,
                 upper_percentiles, alpha=0.2)
 ax.set_yscale('log')
 
-coleman_etal_vl = [8.914378029, 9.140549273, 7.91276252, 8.688206785, 9.366720517,9.172859451, 8.946688207]
+coleman_etal_vl = [8.914378029, 9.140549273, 7.91276252,
+                   8.688206785, 9.366720517, 9.172859451, 8.946688207]
 coleman_etal_er = [127, 455.2, 281.8, 884.2, 448.4, 1100.6, 621]
 plt.scatter(coleman_etal_vl, coleman_etal_er)
 
@@ -82,9 +84,40 @@ yan_vl = [9.347856705]
 yan_er = [45324.55964]
 plt.scatter(yan_vl, yan_er)
 
+# Milton
+boxes = [
+    {
+        'label': "Milton data",
+        'whislo': 0,    # Bottom whisker position
+        'q1': 22,    # First quartile (25th percentile)
+        'med': 220,    # Median         (50th percentile)
+        'q3': 1120,    # Third quartile (75th percentile)
+        'whishi': 260000,    # Top whisker position
+        'fliers': []        # Outliers
+    }
+]
+# `box plot aligned with the viral load value of 5.62325
+ax.bxp(boxes, showfliers=False, positions=[5.62324929])
+
+# Yan
+
+boxes = [
+    {
+        'whislo': 1424.81,    # Bottom whisker position
+        'q1': 8396.78,    # First quartile (25th percentile)
+        'med': 45324.6,    # Median         (50th percentile)
+        'q3': 400054,    # Third quartile (75th percentile)
+        'whishi': 88616200,    # Top whisker position
+        'fliers': []        # Outliers
+    }
+]
+ax.bxp(boxes, showfliers=False, positions=[9.34786])
+
+# box plot aligned with the viral load value of 9.34786
+
 plt.title('Exhaled virions while breathing for 1h', fontsize=14)
 plt.ylabel('RNA copies', fontsize=12)
 plt.xticks(ticks=[i for i in range(2, 13)], labels=[
-            '$10^{' + str(i) + '}$' for i in range(2, 13)])
+    '$10^{' + str(i) + '}$' for i in range(2, 13)])
 plt.xlabel('NP viral load (RNA copies mL$^{-1}$)', fontsize=12)
 plt.show()