mitchell - sampathkumar

EX_9_Script_Q3

@@ -0,0 +1,80 @@
+# Import packages: 
+
+import numpy
+import scipy
+import pandas
+from scipy.optimize import minimize
+from scipy.stats import norm
+from scipy.stats import chi2
+from plotnine import *
+os.listdir('.')
+os.chdir('/Users/sampathkumarbalaji/EX_9/Intro_Biocom_ND_319_Tutorial9')
+
+# Read in data: 
+
+file3=pandas.read_csv("leafDecomp.csv")
+
+# Define custom quadratic function:
+
+def nllike_quadratic(p,obs):
+
+    B0=p[0]
+    B1=p[1]
+    B2=p[2]
+    sigma=p[3]
+
+    expected=B0+B1*obs.Ms+B2*obs.Ms*obs.Ms
+    nll=-1*norm(expected,sigma).logpdf(obs.decomp).sum()
+
+    return nll
+
+# Define custom linear function:
+
+def nllike_linear(p,obs):
+
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0+B1*obs.Ms
+    nll=-1*norm(expected,sigma).logpdf(obs.decomp).sum()
+
+    return nll
+
+# Define custom constant function: 
+
+def nllike_constant(p,obs):
+
+    B0=p[0]
+    sigma=p[1]
+
+    expected=B0
+    nll=-1*norm(expected,sigma).logpdf(obs.decomp).sum()
+
+    return nll
+
+# Set initial guesses for the models:
+
+initialGuess_quadratic=numpy.array([1,1,1,1])
+initialGuess_linear=numpy.array([1,1,1])
+initialGuess_constant=numpy.array([1,1])
+
+# Solve each of the three models:
+
+fitQuad=minimize(nllike_quadratic,initialGuess_quadratic,method="Nelder-Mead",options={'disp': True},args=file3)
+fitLinear=minimize(nllike_linear,initialGuess_linear,method="Nelder-Mead",options={'disp': True},args=file3)
+fitConstant=minimize(nllike_constant,initialGuess_constant,method="Nelder-Mead",options={'disp': True},args=file3)
+
+# Chi squared calculations for the three models:
+
+1-scipy.stats.chi2.cdf(x=-2*(fitQuad.fun-fitConstant.fun),df=2)
+1-scipy.stats.chi2.cdf(x=-2*(fitLinear.fun-fitConstant.fun),df=1)
+1-scipy.stats.chi2.cdf(x=-2*(fitConstant.fun-fitLinear.fun),df=1)
+
+# Print results and interpretation of results:
+
+print("Quadratic Growth Model: " "B0 =",fitQuad.x[0],"B1 =",fitQuad.x[1], "B2 =", fitQuad.x[2], "sigma =", fitQuad.x[3])
+print("Linear Growth Model: " "B0 =",fitLinear.x[0],"B1 =",fitLinear.x[1],"sigma =", fitLinear.x[2])
+print("Constant Growth Model: " "B0 =",fitConstant.x[0],"sigma =",fitConstant.x[1])
+#print("\n")
+#The results suggest that Quadratic Growth Model is best suited. 

Ex9Q1

@@ -0,0 +1,153 @@
+# Import packages, navigate to directory:
+
+import numpy
+import pandas
+from scipy.optimize import minimize
+from scipy.stats import norm
+from scipy.stats import chi2
+from plotnine import *
+os.listdir('.')
+os.chdir('/Users/winghomitchell/Intro_Biocom_ND_319_Tutorial9')
+
+# Load data file:
+
+file=pandas.read_csv("ponzr1.csv")
+file['mutation'] = file['mutation'].map({'WT':0, 'M124K':1, 'V456D':2, 'I213N':3})
+
+# Create custom likelihood function:
+
+def nllike(p,obs):
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0+B1*obs.mutation
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1,1])
+fit=minimize(nllike,initialGuess,method="Nelder-Mead",options={'disp': True},args=file)
+
+print(fit.fun)
+
+# Subset 1: WT & M124K
+
+subset1=file.loc[file.mutation.isin(['0','1']),:]
+
+def nllike(p,obs):
+    B0=p[0]
+    sigma=p[1]
+
+    expected=B0
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1])
+fit1=minimize(nllike,initialGuess,method="Nelder-Mead",options={'disp': True},args=subset1)
+
+def alt(p,obs):
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0+B1*obs.mutation
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1,1])
+fit2=minimize(alt,initialGuess,method="Nelder-Mead",options={'disp': True},args=subset1)
+
+D=2*(fit1.fun-fit2.fun)
+1-chi2.cdf(x=D,df=1)
+sub1ans=1-chi2.cdf(x=D,df=1)
+
+#Subset 2: WT & V456D
+
+subset2=file.loc[file.mutation.isin(['0','2']),:]
+
+def nllike(p,obs):
+    B0=p[0]
+    sigma=p[1]
+
+    expected=B0
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1])
+fit3=minimize(nllike,initialGuess,method="Nelder-Mead",options={'disp': True},args=subset2)
+
+def alt(p,obs):
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0+B1*obs.mutation
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1,1])
+fit4=minimize(alt,initialGuess,method="Nelder-Mead",options={'disp': True},args=subset2)
+
+D=2*(fit3.fun-fit4.fun)
+1-chi2.cdf(x=D,df=1)
+sub2ans=1-chi2.cdf(x=D,df=1)
+
+# Subset 3: WT & I213N
+
+subset3=file.loc[file.mutation.isin(['0','3']),:]
+
+def nllike(p,obs):
+    B0=p[0]
+    sigma=p[1]
+
+    expected=B0
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1])
+fit5=minimize(nllike,initialGuess,method="Nelder-Mead",options={'disp': True},args=subset3)
+
+def alt(p,obs):
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0+B1*obs.mutation
+    nll=-1*norm(expected,sigma).logpdf(obs.ponzr1Counts).sum()
+    return nll
+
+initialGuess=numpy.array([1,1,1])
+fit6=minimize(alt,initialGuess,method="Nelder-Mead",options={'disp': True},args=subset3)
+
+D=2*(fit5.fun-fit6.fun)
+1-chi2.cdf(x=D,df=1)
+sub3ans=1-chi2.cdf(x=D,df=1)
+
+# Print out p-values and interpretation of results for Exercise 9, Question 1. 
+
+print("p = {}" .format(sub1ans))
+
+if sub1ans < 0.05:
+    print ("Reject the null hypothesis. Treatment has an effect.")
+
+else:
+    print ("Fail to reject the null hypothesis. No effect of treatment.")
+
+print("p = {}" .format(sub2ans))
+
+if sub2ans < 0.05:
+    print ("Reject the null hypothesis. Treatment has an effect.")
+
+else:
+    print ("Fail to reject the null hypothesis.")
+
+print("p = {}" .format(sub3ans))
+
+if sub3ans < 0.05:
+    print ("Reject the null hypothesis. Treatment has an effect.")
+
+else:
+    print ("Fail to reject the null hypothesis. No effect of treatment.")
+
+
+

Ex9Q2

@@ -0,0 +1,33 @@
+# Import packages, navigate to directory:
+
+import numpy
+import pandas
+from scipy.optimize import minimize
+from scipy.stats import norm
+from scipy.stats import chi2
+from plotnine import *
+os.listdir('.')
+os.chdir('/Users/winghomitchell/Intro_Biocom_ND_319_Tutorial9')
+
+# Read in data: 
+
+file2=pandas.read_csv("MmarinumGrowth.csv")
+initialGuess=numpy.array([1,1,1])
+
+# Define custom function:
+
+def nllike_bac(p,obs):
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0*((obs.S)/(obs.S+B1))
+    nll=-1*norm(expected,sigma).logpdf(obs.u).sum()
+
+    return nll
+
+fit7=minimize(nllike_bac,initialGuess,method="Nelder-Mead",options={'disp': True},args=file2)
+
+print("Maximum Growth Rate=",fit7.x[0])
+print("Half Saturation Constant=", fit7.x[1])
+

Ex9Q3

@@ -0,0 +1,80 @@
+# Import packages: 
+
+import numpy
+import pandas
+from scipy.optimize import minimize
+from scipy.stats import norm
+from scipy.stats import chi2
+from plotnine import *
+os.listdir('.')
+os.chdir('/Users/winghomitchell/Intro_Biocom_ND_319_Tutorial9')
+
+# Read in data: 
+
+file3=pandas.read_csv("leafDecomp.csv")
+
+# Define custom quadratic function:
+
+def nllike_quadratic(p,obs):
+
+    B0=p[0]
+    B1=p[1]
+    B2=p[2]
+    sigma=p[3]
+
+    expected=B0+B1*obs.Ms+B2*obs.Ms*obs.Ms
+    nll=-1*norm(expected,sigma).logpdf(obs.decomp).sum()
+
+    return nll
+
+# Define custom linear function:
+
+def nllike_linear(p,obs):
+
+    B0=p[0]
+    B1=p[1]
+    sigma=p[2]
+
+    expected=B0+B1*obs.Ms
+    nll=-1*norm(expected,sigma).logpdf(obs.decomp).sum()
+
+    return nll
+
+# Define custom constant function: 
+
+def nllike_constant(p,obs):
+
+    B0=p[0]
+    sigma=p[1]
+
+    expected=B0
+    nll=-1*norm(expected,sigma).logpdf(obs.decomp).sum()
+
+    return nll
+
+# Set initial guesses for the models:
+
+initialGuess_quadratic=numpy.array([1,1,1,1])
+initialGuess_linear=numpy.array([1,1,1])
+initialGuess_constant=numpy.array([1,1])
+
+# Solve each of the three models:
+
+fitQuad=minimize(nllike_quadratic,initialGuess_quadratic,method="Nelder-Mead",options={'disp': True},args=file3)
+fitLinear=minimize(nllike_linear,initialGuess_linear,method="Nelder-Mead",options={'disp': True},args=file3)
+fitConstant=minimize(nllike_constant,initialGuess_constant,method="Nelder-Mead",options={'disp': True},args=file3)
+
+# Chi squared calculations for the three models:
+
+1-scipy.stats.chi2.cdf(x=-2*(fitQuad.fun-fitConstant.fun),df=2)
+1-scipy.stats.chi2.cdf(x=-2*(fitLinear.fun-fitConstant.fun),df=1)
+1-scipy.stats.chi2.cdf(x=-2*(fitConstant.fun-fitLinear.fun),df=1)
+
+# Print results and interpretation of results:
+
+print("Quadratic Growth Model: " "B0 =",fitQuad.x[0],"B1 =",fitQuad.x[1], "B2 =", fitQuad.x[2], "sigma =", fitQuad.x[3])
+print("Linear Growth Model: " "B0 =",fitLinear.x[0],"B1 =",fitLinear.x[1],"sigma =", fitLinear.x[2])
+print("Constant Growth Model: " "B0 =",fitConstant.x[0],"sigma =",fitConstant.x[1])
+print("\n")
+print("These results suggest...")
+