+
+
+
+
+
+
+
+
+
+
+The first column consists of the weights of chicks and the second column consists of the feed type (either horsebean, linseed, soybean, sunflower, meatmeal, or casein). The two columns are separated by a comma, so the file is a "csv".¶
+
+
+
+
+
+
+
+
+Import all the necessary packages
+ +
+
+
+
+
+In [36]:
+
+
+
+
+
+
+#Load packages
+import numpy
+import pandas
+from plotnine import *
+
+
+
+
+
+
+
+
+Load the csv file by pandas.read, separator = ,
+ +
+
+
+
+
+In [37]:
+
+
+
+
+
+
+#Load file
+chicks=pandas.read_csv("chickwts.txt", header=0, sep=",")
+
+
+
+
+
+
+
+
+Generating a plot that summarizes the mean weights of chicks when fed different feed types¶
+
+
+
+
+
+In [38]:
+
+
+
+
+
+
+#Make a bar graph that summarizes all data
+d=ggplot(chicks)+theme_classic()+xlab("feed")+ylab("weight")
+d+geom_bar(aes(x="factor(feed)",y="weight"),stat="summary",fun_y=numpy.mean)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+Out[38]:
+
+
+
+
+
+
+
+
+
+
+
+
+
+Hnull : chick weight when fed soybean = chick weight when fed sunflower
+Halternative: chick weight when fed soybean ≠chick weight when fed sunflower
+ +
+
+
+
+
+
+
+
+Testing the null hypthothesis using likelihood ratio test¶
+
+
+
+
+
+In [39]:
+
+
+
+
+
+
+#Import packages
+import numpy
+import pandas
+from scipy.optimize import minimize
+from scipy.stats import norm
+from plotnine import *
+
+
+
+
+
+In [40]:
+
+
+
+
+
+
+#Subsetting data so it only includes sunflower & soybean data points
+SFvsSB=chicks.loc[chicks.feed.isin(['sunflower', 'soybean']),:]
+
+
+
+
+
+In [41]:
+
+
+
+
+
+
+#Creating a new dataframe to change x column to 0's and 1's
+#For example, replaces "sunflower" with "0" and "soybean" with "1"
+SFvsSBdf=pandas.DataFrame({'y':SFvsSB.weight, 'x':0})
+SFvsSBdf.loc[SFvsSB.feed=='sunflower', 'x']=1
+
+
+
+
+
+In [42]:
+
+
+
+
+
+
+#Null hypothesis likelihood ratio equation
+def nllike(p,obs):
+ B0=p[0]
+ sigma=p[1]
+ expected=B0
+ nll=-1*norm(expected,sigma).logpdf(obs.y).sum()
+ return nll
+
+
+
+
+
+In [43]:
+
+
+
+
+
+
+#Alternative hypothesis likelihood ratio equation
+def nllike2(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
+
+
+
+
+
+In [44]:
+
+
+
+
+
+
+#Estimating parameters by minimizing the nll
+initialVals1=numpy.array([1,1,1])
+
+
+
+
+
+In [45]:
+
+
+
+
+
+
+#Fit null and alternative hypotheseses to Nelder-Mead model
+fitNull=minimize(nllike,initialVals1, method="Nelder-Mead",options={'disp': True}, args=SFvsSBdf)
+fitAlt=minimize(nllike2,initialVals1, method="Nelder-Mead",options={'disp': True}, args=SFvsSBdf)
+print(fitNull.x)
+print(fitAlt.x)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+In [46]:
+
+
+
+
+
+
+#Calculating p-value
+from scipy.stats import chi2
+D=(2*(fitNull.fun-fitAlt.fun))
+feedanswer=(1-chi2.cdf(x=D,df=1))
+print('weight p value')
+print(feedanswer)
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+Because p=0.00023, which is less than 0.01, the alternative hypothesis is true¶
+
+
+
+
+
+
+
+
+