From d83abf9436baed759ac605b0bf5f0f8713590da5 Mon Sep 17 00:00:00 2001 From: Soren Holm Date: Sun, 5 Nov 2017 17:47:30 -0500 Subject: [PATCH 1/5] Solution 1 not working --- solution1.py | 46 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 46 insertions(+) create mode 100644 solution1.py diff --git a/solution1.py b/solution1.py new file mode 100644 index 0000000..52eccd4 --- /dev/null +++ b/solution1.py @@ -0,0 +1,46 @@ +import pandas +import scipy +import scipy.integrate as spint +from plotnine import * + +def ddSim(y,t0,r,K): + N=y[0] + dNdt=r*(1-N/K)*N + return [dNdt] + +N0 = [10] +K = 100 +listr = (-0.1,0.1, 0.4, 0.8, 1.0) +times = range(0,600) +modelOutput=pandas.DataFrame({"time":times,"r1":0, "r2":0, "r3":0, "r4":0, "r5":0}) + +for i in range(len(listr)): + params = (listr[i], K) + modelSim=spint.odeint(func=ddSim,y0=N0,t=times,args=params) + modelOutput.iloc[:,i]=modelSim[:,0] + +ggplot(modelOutput,aes(x="time",y="r1"))+geom_line()+theme_classic() + + +listK = (10,50,100) +r = 0.2 +N0 = [1] +modelOutput=pandas.DataFrame({"time":times,"K1":0, "K2":0, "K3":0}) +for i in range(len(listK)): + params = (r, listK[i]) + modelSim=spint.odeint(func=ddSim,y0=N0,t=times,args=params) + modelOutput.iloc[:,i]=modelSim[:,0] + +ggplot(modelOutput,aes(x="time",y="K1"))+geom_line()+theme_classic() + + +listN0 = (1, 50, 100) +params = (0.1, 50) +modelOutput=pandas.DataFrame({"time":times,"N1":0, "N2":0, "N3":0}) +for i in range(len(listN0)): + modelSim=spint.odeint(func=ddSim,y0=listN0[i],t=times,args=params) + modelOutput.iloc[:,i]=modelSim[:,0] + +ggplot(modelOutput,aes(x="time",y="N1"))+geom_line()+theme_classic() + + From 6628f65b8c6c73e8d4893386b0fc7166b26e05d4 Mon Sep 17 00:00:00 2001 From: Soren Holm Date: Sun, 5 Nov 2017 18:00:22 -0500 Subject: [PATCH 2/5] need help --- solution2.py | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) create mode 100644 solution2.py diff --git a/solution2.py b/solution2.py new file mode 100644 index 0000000..7bff9a4 --- /dev/null +++ b/solution2.py @@ -0,0 +1,26 @@ +import pandas +import scipy +import scipy.integrate as spint +from plotnine import * + +def ddSim(y,t0,b,y): + S=y[0] + I=y[1] + R=y[3] + dSdt= -b * I * S + dIdt = b * I * S - y * I + dRdt = y * I + return [dSdt, dIdt, dRdt] + +listb = (0.0005, 0.005, 0.0001, 0.00005, 0.0001, 0.0002, 0.0001) +listy = (0.05, 0.5, 0.1, 0.1, 0.05, 0.05, 0.06) + +days = range(0,500) +susceptible = 999 +infected = 1 +resistant = 0 +initial = [susceptible, infected, resistant] + +for i in range(len(listb)): + params = (listb[i], listy[i]) + \ No newline at end of file From eba78feef5b9fb24c3a13844dd90800fa4a14613 Mon Sep 17 00:00:00 2001 From: Soren Holm Date: Sun, 5 Nov 2017 18:13:24 -0500 Subject: [PATCH 3/5] n --- solution2.py | 1 + 1 file changed, 1 insertion(+) diff --git a/solution2.py b/solution2.py index 7bff9a4..1bd6049 100644 --- a/solution2.py +++ b/solution2.py @@ -23,4 +23,5 @@ def ddSim(y,t0,b,y): for i in range(len(listb)): params = (listb[i], listy[i]) + modelSim=spint.odeint(func=ddSim,y0=initial,t=times,args=params) \ No newline at end of file From 77da957e32890fa622d203d249a980268aeaf16d Mon Sep 17 00:00:00 2001 From: zoeloh <31715186+zoeloh@users.noreply.github.com> Date: Thu, 9 Nov 2017 14:38:24 -0500 Subject: [PATCH 4/5] Finished problem 2 --- solution2.py | 32 ++++++++++++++++++++++++++------ 1 file changed, 26 insertions(+), 6 deletions(-) diff --git a/solution2.py b/solution2.py index 1bd6049..e612249 100644 --- a/solution2.py +++ b/solution2.py @@ -3,13 +3,13 @@ import scipy.integrate as spint from plotnine import * -def ddSim(y,t0,b,y): +def ddSim(y,t0,b,a): S=y[0] I=y[1] - R=y[3] + R=y[2] dSdt= -b * I * S - dIdt = b * I * S - y * I - dRdt = y * I + dIdt = b * I * S - a * I + dRdt = a * I return [dSdt, dIdt, dRdt] listb = (0.0005, 0.005, 0.0001, 0.00005, 0.0001, 0.0002, 0.0001) @@ -23,5 +23,25 @@ def ddSim(y,t0,b,y): for i in range(len(listb)): params = (listb[i], listy[i]) - modelSim=spint.odeint(func=ddSim,y0=initial,t=times,args=params) - \ No newline at end of file + modelSim=spint.odeint(func=ddSim,y0=initial,t=days,args=params) + modelOutput = pandas.DataFrame({"Days":days,"S":modelSim[:,0], "I":modelSim[:,1], "R":modelSim[:,2]}) + incidence = [] + prevalence = [] + for j in range(0,500): + if j != 0: + incidence.append(modelOutput.loc[j, "I"] - modelOutput.loc[j-1, "I"]) + + prevalence.append(modelOutput.loc[j, "I"] / (modelOutput.loc[j, "I"] + modelOutput.loc[j, "S"] + modelOutput.loc[j, "R"])) + incidence.sort() + prevalence.sort() + + print("B: " + str(listb[i]) + " y: " + str(listy[i])) + print("Max incidence: " + str(incidence[0])) + print("Max prevalence: " + str(prevalence[0])) + pctAffected = (modelOutput.loc[499, "I"] + modelOutput.loc[499, "R"]) / (modelOutput.loc[499, "S"] + modelOutput.loc[499, "I"] + modelOutput.loc[499, "R"]) + print("Percent affected: " + str(pctAffected)) + bscReproduction = listb[i] * (modelOutput.loc[499, "S"] + modelOutput.loc[499, "I"] + modelOutput.loc[499, "R"]) / listy[i] + print("Basic reproduction number: " + str(bscReproduction)) + print("############################################") + + From 761a135b4c1cc8ce81a022d98093c75e8704d0e0 Mon Sep 17 00:00:00 2001 From: zoeloh <31715186+zoeloh@users.noreply.github.com> Date: Fri, 10 Nov 2017 05:38:54 -0500 Subject: [PATCH 5/5] Final --- solution1.py | 12 +++++++----- 1 file changed, 7 insertions(+), 5 deletions(-) diff --git a/solution1.py b/solution1.py index 52eccd4..df6f35f 100644 --- a/solution1.py +++ b/solution1.py @@ -19,7 +19,8 @@ def ddSim(y,t0,r,K): modelSim=spint.odeint(func=ddSim,y0=N0,t=times,args=params) modelOutput.iloc[:,i]=modelSim[:,0] -ggplot(modelOutput,aes(x="time",y="r1"))+geom_line()+theme_classic() +g1=ggplot(modelOutput,aes(x="time",y="r1"))+geom_line()+geom_line(aes(y="R2"))+geom_line(aes(y="R3"))+geom_line(aes(y="R4"))+geom_line(aes(y="R5")) +g1+theme_classic() listK = (10,50,100) @@ -31,8 +32,8 @@ def ddSim(y,t0,r,K): modelSim=spint.odeint(func=ddSim,y0=N0,t=times,args=params) modelOutput.iloc[:,i]=modelSim[:,0] -ggplot(modelOutput,aes(x="time",y="K1"))+geom_line()+theme_classic() - +g2=ggplot(modelOutput,aes(x="time",y="K1"))+geom_line()+geom_line(aes(y="K2"))+geom_line(aes(y="K3")) +g2+theme_classic() listN0 = (1, 50, 100) params = (0.1, 50) @@ -41,6 +42,7 @@ def ddSim(y,t0,r,K): modelSim=spint.odeint(func=ddSim,y0=listN0[i],t=times,args=params) modelOutput.iloc[:,i]=modelSim[:,0] -ggplot(modelOutput,aes(x="time",y="N1"))+geom_line()+theme_classic() - +g3=ggplot(modelOutput,aes(x="time",y="N1"))+geom_line()+geom_line(aes(y="N2"))+geom_line(aes(y="N3")) +g3+theme_classic() +#We dont know how to get all the plots outputted, but the code is there to make them all. We asked stuart but he didn't help us with this specifically