From 021042ce95226c2df946985cc5cede662a5e6098 Mon Sep 17 00:00:00 2001 From: Soren Holm Date: Sun, 29 Oct 2017 18:31:38 -0400 Subject: [PATCH 1/5] Initial commit --- solution1.py | 16 ++++++++++++++++ solution2.py | 0 solution3.py | 0 3 files changed, 16 insertions(+) create mode 100644 solution1.py create mode 100644 solution2.py create mode 100644 solution3.py diff --git a/solution1.py b/solution1.py new file mode 100644 index 0000000..7b32365 --- /dev/null +++ b/solution1.py @@ -0,0 +1,16 @@ +import numpy +import pandas +from scipy.optimize import minimize +from scipy.stats import norm +from plotnine import * + +def nllike(p, obs): + B0 = p[0] + B1 = p[1] + sigma = p[2] + + expected = B0 + (B1 * obs.x) + nll = -1 * norm(expected, sigma).logpdf(obs.y).sum() + return nll + +data = pandas.read_csv('ponzr1.csv', header = 0, sep = ',') diff --git a/solution2.py b/solution2.py new file mode 100644 index 0000000..e69de29 diff --git a/solution3.py b/solution3.py new file mode 100644 index 0000000..e69de29 From f8bb1304c990427ce440317dec64ff003f5e0a3d Mon Sep 17 00:00:00 2001 From: Zoe Loh Date: Thu, 2 Nov 2017 10:12:13 -0400 Subject: [PATCH 2/5] problem 1 --- solution1.py | 32 +++++++++++++++++++++++++++++--- 1 file changed, 29 insertions(+), 3 deletions(-) diff --git a/solution1.py b/solution1.py index 7b32365..99e0e2a 100644 --- a/solution1.py +++ b/solution1.py @@ -2,15 +2,41 @@ import pandas from scipy.optimize import minimize from scipy.stats import norm +from scipy.stats import chi2 from plotnine import * -def nllike(p, obs): +def linear(p, obs): B0 = p[0] B1 = p[1] sigma = p[2] - expected = B0 + (B1 * obs.x) - nll = -1 * norm(expected, sigma).logpdf(obs.y).sum() + expected = B0 + (B1 * obs.mutation) + nll = -1 * norm(expected, sigma).logpdf(obs.ponzr1Counts).sum() return nll + +def nullH(p, obs): + B0 = p[0] + sigma = p[1] + + + expected = B0 + nll = -1 * norm(expected, sigma).logpdf(obs.ponzr1Counts).sum() + return nll + data = pandas.read_csv('ponzr1.csv', header = 0, sep = ',') +data['mutation'] = data["mutation"].map({'WT' : 0, 'M124K' : 1, 'V456D' : 2, 'I213N' : 3}) + +subset=data.loc[data.mutation.isin(['0','1']),:] + +initialLinGuess = numpy.array([1,1,1]) +linfit = minimize(linear, initialLinGuess, method="Nelder-Mead", options={'disp': True}, args=subset) + +initialNullGuess = numpy.array([1,1]) +nullfit = minimize(nullH, initialNullGuess, method="Nelder-Mead", options={'disp': True}, args=subset) + +linfit.fun +nullfit.fun + +D=2*(nullfit.fun - linfit.fun) +1-chi2.cdf(x=D, df=1) From 90dfc8a8dc968d1b2b363b3d89e417d3445e767c Mon Sep 17 00:00:00 2001 From: Zoe Loh Date: Thu, 2 Nov 2017 10:21:47 -0400 Subject: [PATCH 3/5] Problem 2 --- solution2.py | 22 ++++++++++++++++++++++ 1 file changed, 22 insertions(+) diff --git a/solution2.py b/solution2.py index e69de29..6601d28 100644 --- a/solution2.py +++ b/solution2.py @@ -0,0 +1,22 @@ +import numpy +import pandas +from scipy.optimize import minimize +from scipy.stats import norm +from scipy.stats import chi2 +from plotnine import * + +def monod(p, obs): + UMAX = p[0] + KS = p[1] + sigma = p[2] + + expected = UMAX*(obs.S / (obs.S + KS)) + nll = -1 * norm(expected, sigma).logpdf(obs.u).sum() + return nll + +data = pandas.read_csv('MmarinumGrowth.csv', header = 0, sep = ',') + +initialGuess = numpy.array([1,1,1]) +fit = minimize(monod, initialGuess, method="Nelder-Mead", options={'disp': True}, args=data) + +print(fit.x) \ No newline at end of file From 2ca5b9240fbceb6c58ca83b23ace39a63aef3518 Mon Sep 17 00:00:00 2001 From: Zoe Loh Date: Thu, 2 Nov 2017 10:45:35 -0400 Subject: [PATCH 4/5] problem 3 --- solution3.py | 46 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 46 insertions(+) diff --git a/solution3.py b/solution3.py index e69de29..d1e5b0f 100644 --- a/solution3.py +++ b/solution3.py @@ -0,0 +1,46 @@ +import numpy +import pandas +from scipy.optimize import minimize +from scipy.stats import norm +from scipy.stats import chi2 +from plotnine import * + +def const(p, obs): + a = p[0] + sigma = p[1] + + expected = a + nll = -1 * norm(expected, sigma).logpdf(obs.decomp).sum() + return nll + +def linear(p, obs): + a = p[0] + b = p[1] + sigma = p[2] + + expected = a + b * obs.Ms + nll = -1 * norm(expected, sigma).logpdf(obs.decomp).sum() + return nll + +def hump(p, obs): + a = p[0] + b = p[1] + c = p[2] + sigma = p[3] + + expected = a + b*obs.Ms + c*(obs.Ms) * (obs.Ms) + nll = -1 * norm(expected, sigma).logpdf(obs.decomp).sum() + return nll + +data = pandas.read_csv('leafDecomp.csv', header = 0, sep = ',') + +initialGuess = numpy.array([1,1,1,1]) +constfit = minimize(const, initialGuess, method="Nelder-Mead", options={'disp': True}, args=data) +linfit = minimize(linear, initialGuess, method="Nelder-Mead", options={'disp': True}, args=data) +humpfit = minimize(hump, initialGuess, method="Nelder-Mead", options={'disp': True}, args=data) + + +print(humpfit.x) + +D=2*(const.fun - linfit.fun) +1-chi2.cdf(x=D, df=1) From 0d7033e09b1b32781f77d412c968569023fc9c5c Mon Sep 17 00:00:00 2001 From: Soren Holm Date: Fri, 3 Nov 2017 06:14:21 -0400 Subject: [PATCH 5/5] ready for submission --- solution1.py | 35 ++++++++++++++++++++++++++++++++++- solution3.py | 13 ++++++++++--- 2 files changed, 44 insertions(+), 4 deletions(-) diff --git a/solution1.py b/solution1.py index 99e0e2a..3c1a09f 100644 --- a/solution1.py +++ b/solution1.py @@ -39,4 +39,37 @@ def nullH(p, obs): nullfit.fun D=2*(nullfit.fun - linfit.fun) -1-chi2.cdf(x=D, df=1) +print("M124K P value:") +print(1-chi2.cdf(x=D, df=1)) + +subset=data.loc[data.mutation.isin(['0','2']),:] + +initialLinGuess = numpy.array([1,1,1]) +linfit = minimize(linear, initialLinGuess, method="Nelder-Mead", options={'disp': True}, args=subset) + +initialNullGuess = numpy.array([1,1]) +nullfit = minimize(nullH, initialNullGuess, method="Nelder-Mead", options={'disp': True}, args=subset) + +linfit.fun +nullfit.fun + +D=2*(nullfit.fun - linfit.fun) +print("V456D P value:") +print(1-chi2.cdf(x=D, df=1)) + +subset=data.loc[data.mutation.isin(['0','3']),:] + +initialLinGuess = numpy.array([1,1,1]) +linfit = minimize(linear, initialLinGuess, method="Nelder-Mead", options={'disp': True}, args=subset) + +initialNullGuess = numpy.array([1,1]) +nullfit = minimize(nullH, initialNullGuess, method="Nelder-Mead", options={'disp': True}, args=subset) + +linfit.fun +nullfit.fun + +D=2*(nullfit.fun - linfit.fun) +print("I213N P value:") +print(1-chi2.cdf(x=D, df=1)) + + diff --git a/solution3.py b/solution3.py index d1e5b0f..d6507d5 100644 --- a/solution3.py +++ b/solution3.py @@ -40,7 +40,14 @@ def hump(p, obs): humpfit = minimize(hump, initialGuess, method="Nelder-Mead", options={'disp': True}, args=data) -print(humpfit.x) +D=2*(constfit.fun - linfit.fun) +print("linear over const fit") +print(1-chi2.cdf(x=D, df=1)) -D=2*(const.fun - linfit.fun) -1-chi2.cdf(x=D, df=1) +D=2*(humpfit.fun - linfit.fun) +print("linear over hump fit") +print(1-chi2.cdf(x=D, df=1)) + +D=2*(constfit.fun - humpfit.fun) +print("hump over const fit") +print(1-chi2.cdf(x=D, df=2))