Problem Set 7

A) Why would we model sires as a random effect?

##This study is looking at average weight gain of calves sired by different bulls.We should model this as a random effect because we are not interested in these bulls in particular, but rather are using these bulls to make an inference about the population of bulls and their calves.

B) Write the factor effects model with a random effect for sire.

#y(ij) = mu + alpha(i) + error(ij)

#in terms of this problem:

#weight gain = intercept(base weight gain) + random effect of sire + random error

C) Fit a mixed model with a random intercept for sire using the lmer function in R

library(lme4)
## Loading required package: Matrix
Weight <- c(1.46, 1.23, 1.12, 1.23, 1.02, 1.15,1.17, 1.08, 1.20, 1.08, 1.01, 0.86, 0.98, 1.06, 1.15, 1.11, 0.83, 0.86, 0.95, 1.10, 1.07, 1.11, 0.89, 1.12)

Sired <- c("A","A","A","A","A","A","B","B","B","B","B","B","C","C","C","C","C","C","D","D","D","D","D","D")



calves <- data.frame(Sired,Weight)


modC <- lmer(Weight~1 +(1|Sired), data = calves)
summary(modC)
## Linear mixed model fit by REML ['lmerMod']
## Formula: Weight ~ 1 + (1 | Sired)
##    Data: calves
## 
## REML criterion at convergence: -23.5
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.6639 -0.5601  0.1081  0.5645  2.3859 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  Sired    (Intercept) 0.005078 0.07126 
##  Residual             0.015945 0.12627 
## Number of obs: 24, groups:  Sired, 4
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)  1.07667    0.04397   24.48

D) Fit a linear model with only an intercept using the lm function.

modD <- lm(Weight~1, data = calves)
summary(modD)
## 
## Call:
## lm(formula = Weight ~ 1, data = calves)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.24667 -0.07417  0.01333  0.07333  0.38333 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.07667    0.02881   37.37   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.1411 on 23 degrees of freedom

E) In order to test the null hypothesis that there is no sire variability in the response, compare the models in part C and part D with the anova function. Report the results. Is it worth including a random effect for sire?

anova(modC,modD)
## refitting model(s) with ML (instead of REML)
## Data: calves
## Models:
## modD: Weight ~ 1
## modC: Weight ~ 1 + (1 | Sired)
##      Df     AIC     BIC logLik deviance  Chisq Chi Df Pr(>Chisq)
## modD  2 -22.898 -20.542 13.449  -26.898                         
## modC  3 -22.095 -18.561 14.047  -28.095 1.1962      1     0.2741

#These results show that there is not a significant difference between the model that includes the random effect of the sire and the model that does not account for the sire. Therefore, it is not worth including the random effect because the added complexity does not significantly improve the prediction of weight gain.

F) What are the variance component estimates for the model with a random effect for sire? (model in part C).

#The variance component for the random effect for sire is 0.005078. This is the amount of variance in weight gain that is explained by the difference in sire. The residual variance component is 0.015945. This is the amount of variance in weight gain that is left unexplained due to random error.

G) Find the 95% confidence intervals for the error standard deviation and the sire to sire standard deviation using the confint function.

confint(modC)
## Computing profile confidence intervals ...
##                  2.5 %    97.5 %
## .sig01      0.00000000 0.1795777
## .sigma      0.09534633 0.1783855
## (Intercept) 0.97994356 1.1733897

#sig01 refers to the confidence interval around the standard deviation of the random effect of sire. The 95% confidence interval is [0.00000000, 0.1795777].

#sigma refers to the confidence interval around the standard deviation of the error. The 95% confidence interval is [0.09534633 0.1783855].