> ?A> %bjbj 4( .4(\\\\7"Ye,r Qm77mmA\\LAAAmj\\AmAA:,-@)3j
b0R A -A-AxmmmAmmm mmmmmmmmm :
Likelihood 2011, Homework 6
Skink quantitative geneticsAmazing advances in skink husbandry now allow researchers to conduct breeding experiments. One such researcher decided to apply the foursome design to investigate genetic variation in growth rate on different diets. He collected 100 adult males and 100 adult females from a natural population and randomly sorts them into mating groups (2 males and 2 females per group). Each male is mated to each family yielding four clutches per mating group. Six baby skinks were grown from each clutch and measured for growth rate (y) on different diets. The three experimental diets were small insects (S), medium insects (M), and large insects (L). Each clutch had two skinks per diet type. A csv file (which can be read by excel) can be found at:
HYPERLINK "http://phylo.bio.ku.edu/slides/skinks_eat_bugs.csv" http://phylo.bio.ku.edu/slides/skinks_eat_bugs.csv
It contains data from this experiment. A total of 1200 skinks were measured, each of which was adopted into a loving home after completion of the experiment. The model for individual measurements is
Y = (mean given diet) + (mom effect) + (dad effect) + normal error
Here, Y is the dependent variable (growth rate)
Mean given diet = S or M or L depending on diet treatment (fixed effect)
mom effect = value from normal distribution: mean=0 and variance=2g (random effect)
dad effect = value from normal distribution: mean=0 and variance=2g (random effect)
normal error = value from normal distribution: mean=0 and variance=2 (random effect)
Note that because moms and dads were sorted randomly in mating groups, the mom effect is independent of the dad effect for any skink.
Part I. Determine the maximum log-likelihood and MLE for the five parameters given the entire dataset. Using the python template, you will need to populate the mating group specific incidence matrix (X) and V-matrix (V). The program does the rest, including the search.
Part II. Test (A) the null hypothesis that diet treatment has no effect on skink growth and (B) the null hypothesis that there is no genetic variation for growth (2g = 0). For each, calculate the LRT statistic and comparing it to the appropriate critical value from the chi-square distribution.
Part III. Use parametric boostrapping as an alternative to the chi-square distribution for evaluation of LRT in Part II. In each simulation replicate you will need to create a dataset of the same size and character as the real data. For each cyber-skink, this means that you need to create a value for Y given diet treatment and the random effects (see equation above). For this you will need to use the:
normalvariate(mean, sd)
function from pythons random module with an appropriate mean and sd (note that the function does take the standard deviation as the second argument, not the variance)
Part III Addendum (added after class discuss on April 20): Include histograms of the LRT and parameter estimates from your simulations so that you can compare the distribution of LRT values to the chi-square distribution. The template has been modified so that you can it creates two files ParametricBootstrapParam.csv and ParametricBootstrapLRT.csv You will need to modify the loop over simulation replicates to write the appropriate information to these files.
NOTE: repeated running of the script will overwrite the contents of these files. So if you want to preserve the contents of the files, you will need to copy them to a new location before running the script again. The .csv files can be read by excel or R.
A template script can be found at:
HYPERLINK "http://phylo.bio.ku.edu/slides/template_HW6.py.txt" http://phylo.bio.ku.edu/slides/template_HW6.py.txt
Notes about the new template:
1. The script now must be told where to find the data (since we are no longer pasting the data into the script itself). The path to the data file is the first argument.
The easiest thing to do is to:
Download the template script and skinks_eat_bugs.csv to the same directory
Open your terminal (or CMD command prompt on Windows)
Use the cd command to navigate to the directory with the script.
This should show you the help message which lists the order of the parameters in the invocation.
2. In addition to constraining a parameter to a valu89f
()5gijk>@hi.02轹赭hqndh7#h|=5h6W{h6W{H* h6W{H*h6W{h6W{H*h6W{h\hFhn3Ihw}hH0JjhHUhHhHjhHUhHhUtXhkTh~;[h|=hvphf67$ij*
T^w gdH`gdHgd6W{~>?\^_vw 1E]h~*\^_`$%%%%%ļϸ۴ϝϖϋπ|xh2s h|Uhw}hH0JjAhHUhHhHjhHUh7#h7#h7#h7#5h7#hhHOJQJhHhOJQJhHh7#h|=5h|=h6W{hf6H*h6W{hf6H*hqndhf6h7#hf65/_`D$$%%%%%%
&FgdH`gdHgdHe, the template now supports constraining two parameters to have the same value. If you want parameter 0 and parameter 1 to be forced to have the same value, but no other constraints then you would type:
python template_HW6.py.txt skinks_eat_bugs.csv None Equals0 None None None 10
3. This analysis is much slower. It took me 45 minutes to do 500 parametric bootstrap replicates on my laptop.
4. Once again, areas that youll need to modify are near the top and are marked with ALL CAPS.
%%%h|=hlUhlU21h:pqnd/ =!"#$%ADyK3http://phylo.bio.ku.edu/slides/skinks_eat_bugs.csvyKfhttp://phylo.bio.ku.edu/slides/skinks_eat_bugs.csvADyK3http://phylo.bio.ku.edu/slides/template_HW6.py.txtyKfhttp://phylo.bio.ku.edu/slides/template_HW6.py.txt6666666662 0@P`p2( 0@P`p 0@P`p 0@P`p 0@P`p 0@P`p 0@P`p8XV~_HmH nH sH tH @`@NormalCJ_HaJmH sH tH DA DDefault Paragraph FontRiRTable Normal4
l4a(k (No List6U`6H Hyperlink>*B*phFV`FlUFollowedHyperlink>*B*phPK!pO[Content_Types].xmlj0Eжr(]yl#!MB;.n̨̽\A1&ҫ
QWKvUbOX#&1`RT9<l#$>r `С-;c=1g~'}xPiB$IO1Êk9IcLHY<;*v7'aE\h>=^,*8q;^*4?Wq{nԉogAߤ>8f2*<")QHxK
|]Zz)ӁMSm@\&>!7;wP3[EBU`1OC5(F\;ܭqpߡ 69&MDO,ooVM M_ո۹U>7eo >ѨN6}
bvzۜ6?ߜŷiLvm]2SFnHD]rISXO]0 ldC^3شd$s#2.h565!v.chNt9W
dumԙgLStf+]C9P^%AW̯f$Ҽa1Q{B{mqDl
u" f9%k@f?g$p0%ovkrt ֖ ? &6jج="MN=^gUn.SƙjмCR=qb4Y" )yvckcj+#;wb>VD
Xa?p
S4[NS28;Y[,T1|n;+/ʕj\\,E:!
t4.T̡e1
}; [z^pl@ok0e
g@GGHPXNT,مde|*YdT\Y䀰+(T7$ow2缂#G֛ʥ?qNK-/M,WgxFV/FQⷶO&ecx\QLW@H!+{[|{!KAi
`cm2iU|Y+ި [[vxrNE3pmR
=Y04,!&0+WC܃@oOS2'Sٮ05$ɤ]pm3FtGɄ-!y"ӉV
.
`עv,O.%вKasSƭvMz`3{9+e@eՔLy7W_XtlPK!
ѐ'theme/theme/_rels/themeManager.xml.relsM
0wooӺ&݈Э5
6?$Q
,.aic21h:qm@RN;d`o7gK(M&$R(.1r'JЊT8V"AȻHu}|$b{P8g/]QAsم(#L[PK-!pO[Content_Types].xmlPK-!֧6-_rels/.relsPK-!kytheme/theme/themeManager.xmlPK-!!Z!theme/theme/theme1.xmlPK-!
ѐ'(
theme/theme/_rels/themeManager.xml.relsPK]#
(%%%(jGXX8@0( A@
B
S ?
_Hlt164998600@
*-
::::::::.xqC.8^8`o(.
^`hH.
L^ `LhH.
^`hH.
x^x`hH.
HL^H`LhH.
^`hH.
^`hH.
L^`LhH.8^8`o(.
^`hH.
L^ `LhH.
^`hH.
x^x`hH.
HL^H`LhH.
^`hH.
^`hH.
L^`LhH.xqC? ?
$2s kT0eR[f67#f$yr.|=\hFn3IUtX~;[qnd!u6W{clU|HwvpC@%@L**** h@hh @hh(@hH@UnknownGTimes New Roman5Symbol3Arial5MonacoACambria Math"qhp! "! "!243QHP(?|=2!xx,Likelihood 2011, Homework 5biohelpMark Holder
Oh+'0x
4@
LX`hp'Likelihood 2011, Homework 5biohelpNormal.dotmMark Holder6Microsoft Macintosh Word@`4<@@3!
՜.+,D՜.+,Hhp|
'" Likelihood 2011, Homework 5Title 8@_PID_HLINKS'A<8K3http://phylo.bio.ku.edu/slides/template_HW6.py.txtn63http://phylo.bio.ku.edu/slides/skinks_eat_bugs.csv
!"#$%&'()*+,-/012345789:;<=@Root Entry Fw)3BData
1Table
!WordDocument4(SummaryInformation(.DocumentSummaryInformation86CompObj` F Microsoft Word 97-2004 DocumentNB6WWord.Document.8