diff --git a/EX_9_Script_Q3 b/EX_9_Script_Q3 new file mode 100644 index 0000000..299f994 --- /dev/null +++ b/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. \ No newline at end of file diff --git a/Ex9Q1 b/Ex9Q1 new file mode 100644 index 0000000..601bd97 --- /dev/null +++ b/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.") + + + diff --git a/Ex9Q2 b/Ex9Q2 new file mode 100644 index 0000000..39a1d12 --- /dev/null +++ b/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]) + diff --git a/Ex9Q3 b/Ex9Q3 new file mode 100644 index 0000000..4320865 --- /dev/null +++ b/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...") + \ No newline at end of file