Title: | Linear Mixed Models - An Introduction with Applications in Veterinary Research |
---|---|
Description: | R Codes and Datasets for Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998). Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute. |
Authors: | Muhammad Yaseen [aut, cre], Luc Duchateau [ctb], Paul Janssen [ctb], John Rowlands [ctb] |
Maintainer: | Muhammad Yaseen <[email protected]> |
License: | GPL-2 |
Version: | 1.0.0 |
Built: | 2024-11-02 02:41:03 UTC |
Source: | https://github.com/myaseen208/vetresearchlmm |
ex121 is.
data(ex121)
data(ex121)
A data.frame
with 40 rows and 4 variables.
herd two treatment 0 and 1
drug unit of observation or observation ID
dose is continuous \& may be assumed Gaussian
PCVDif is the number of "successes"(N and F specify a binomial response)
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex121)
data(ex121)
ex124 is.
data(ex124)
data(ex124)
A data.frame
with 40 rows and 4 variables.
herd two treatment 0 and 1
drug unit of observation or observation ID
dose is continuous \& may be assumed Gaussian
PCVDif is the number of "successes"(N and F specify a binomial response)
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex124)
data(ex124)
ex125 is.
data(ex125)
data(ex125)
A data.frame
with 40 rows and 4 variables.
herd two treatment 0 and 1
drug unit of observation or observation ID
dose is continuous \& may be assumed Gaussian
PCVDif is the number of "successes"(N and F specify a binomial response)
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex125)
data(ex125)
ex127 is.
data(ex127)
data(ex127)
A data.frame
with 40 rows and 4 variables.
herd two treatment 0 and 1
drug unit of observation or observation ID
dose is continuous \& may be assumed Gaussian
PCVDif is the number of "successes"(N and F specify a binomial response)
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex127)
data(ex127)
ex31 is.
data(ex31)
data(ex31)
A data.frame
with 38 rows and 6 variables.
herd Herds of Cattle
animal_id Animal ID
PCV1 Packed Cell Volume (PCV) determined at the time of treatment
PCV2 Packed Cell Volume (PCV) determined at a month later following treatment
dose Dose of Drugs
drug Two drugs against trypanosomosis, Berenil and Samorin, are studied
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex31)
data(ex31)
ex32 is.
data(ex32)
data(ex32)
A data.frame
with 40 rows and 4 variables.
herd two treatment 0 and 1
drug unit of observation or observation ID
dose is continuous \& may be assumed Gaussian
PCVDif is the number of "successes"(N and F specify a binomial response)
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex32)
data(ex32)
ex33 is.
data(ex33)
data(ex33)
A data.frame
with 40 rows and 4 variables.
herd two treatment 0 and 1
drug unit of observation or observation ID
dose is continuous \& may be assumed Gaussian
PCVDif is the number of "successes"(N and F specify a binomial response)
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
data(ex33)
data(ex33)
Examp1.3.2 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 1.3.2 p-16 #------------------------------------------------------------- # PROC GLM DATA=ex124; # CLASS herd dose drug; # MODEL PCVdif=drug herd(drug) dose dose*drug; # RANDOM herd(drug); # RUN; library(lme4) str(ex124) summary(ex124) ex124$herd1 <- factor(ex124$herd) ex124$drug1 <- factor(ex124$drug) ex124$dose1 <- factor(ex124$dose) fm1.1 <- aov( formula = PCVdif ~ drug1 + Error(herd1:drug1) + dose1 + dose1:drug1 , data = ex124 , projections = FALSE , qr = TRUE , contrasts = NULL # , ... ) summary(fm1.1)
#------------------------------------------------------------- ## Example 1.3.2 p-16 #------------------------------------------------------------- # PROC GLM DATA=ex124; # CLASS herd dose drug; # MODEL PCVdif=drug herd(drug) dose dose*drug; # RANDOM herd(drug); # RUN; library(lme4) str(ex124) summary(ex124) ex124$herd1 <- factor(ex124$herd) ex124$drug1 <- factor(ex124$drug) ex124$dose1 <- factor(ex124$dose) fm1.1 <- aov( formula = PCVdif ~ drug1 + Error(herd1:drug1) + dose1 + dose1:drug1 , data = ex124 , projections = FALSE , qr = TRUE , contrasts = NULL # , ... ) summary(fm1.1)
Examp2.4.2.2 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.4.2.2 p-64 #------------------------------------------------------------- # PROC MIXED DATA=ex125 METHOD=ML; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose; # RANDOM region drug*region; # RUN; # # PROC MIXED DATA=ex125 METHOD=REML; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose; # RANDOM region drug*region; # RUN; library(lme4) str(ex125) fm2.4 <- lme4::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = FALSE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.4) anova(fm2.4) fm2.5 <- lme4::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.5) anova(fm2.5) library(lmerTest) fm2.6 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = FALSE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.6) anova(fm2.6) fm2.7 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.7) anova(fm2.7)
#------------------------------------------------------------- ## Example 2.4.2.2 p-64 #------------------------------------------------------------- # PROC MIXED DATA=ex125 METHOD=ML; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose; # RANDOM region drug*region; # RUN; # # PROC MIXED DATA=ex125 METHOD=REML; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose; # RANDOM region drug*region; # RUN; library(lme4) str(ex125) fm2.4 <- lme4::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = FALSE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.4) anova(fm2.4) fm2.5 <- lme4::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.5) anova(fm2.5) library(lmerTest) fm2.6 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = FALSE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.6) anova(fm2.6) fm2.7 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.7) anova(fm2.7)
Examp2.4.3.1 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.4.3.1 p-66 #------------------------------------------------------------- # PROC MIXED DATA=ex127; # CLASS sire; # MODEL ww=; # RANDOM sire/solution; # RUN; library(lme4) str(ex127) fm2.8 <- lme4::lmer( formula = Ww~(1|sire) , data = ex127 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.8) lme4::fixef(fm2.8) lme4::ranef(fm2.8)
#------------------------------------------------------------- ## Example 2.4.3.1 p-66 #------------------------------------------------------------- # PROC MIXED DATA=ex127; # CLASS sire; # MODEL ww=; # RANDOM sire/solution; # RUN; library(lme4) str(ex127) fm2.8 <- lme4::lmer( formula = Ww~(1|sire) , data = ex127 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = NULL , devFunOnly = FALSE # , ... ) summary(fm2.8) lme4::fixef(fm2.8) lme4::ranef(fm2.8)
Examp2.5.1.1 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.5.1.1 p-67 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / solution covb; # RANDOM region drug*region; # RUN; library(lme4) str(ex125) fm2.9 <- lme4::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.9) anova(fm2.9) summary(fm2.9)$vcov
#------------------------------------------------------------- ## Example 2.5.1.1 p-67 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / solution covb; # RANDOM region drug*region; # RUN; library(lme4) str(ex125) fm2.9 <- lme4::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.9) anova(fm2.9) summary(fm2.9)$vcov
Examp2.5.2.1 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.5.2.1 p-68 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / solution covb; # RANDOM region drug*region; # LSMEANS drug*dose; # RUN; library(lmerTest) str(ex125) fm2.10 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.10) anova(fm2.10) summary(fm2.10)$vcov lsmeansLT(model = fm2.10)
#------------------------------------------------------------- ## Example 2.5.2.1 p-68 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / solution covb; # RANDOM region drug*region; # LSMEANS drug*dose; # RUN; library(lmerTest) str(ex125) fm2.10 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.10) anova(fm2.10) summary(fm2.10)$vcov lsmeansLT(model = fm2.10)
Examp2.5.3.1 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.5.3.1 p-70 #------------------------------------------------------------- # PROC GLM DATA=ex125; # CLASS drug dose region; # MODEL pcv=region drug region*drug dose drug*dose; # RANDOM region drug*region; # RUN; # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / ddfm=satterth; # RANDOM region drug*region; # ESTIMATE 'drug dif' drug -1 1 drug*dose -0.5 -0.5 0.5 0.5; # ESTIMATE 'Samorin mean' INTERCEPT 1 drug 0 1 dose 0.5 0.5 # drug*dose 0 0 0.5 0.5; # ESTIMATE 'Samorin HvsL' dose 1 -1 drug*dose 0 0 1 -1; # ESTIMATE 'Samorin high' INTERCEPT 1 drug 0 1 dose 1 0 # drug*dose 0 0 1 0; # RUN; library(lme4) str(ex125) ex125$Region1 <- factor(ex125$Region) fm2.11 <- aov( formula = Pcv ~ Region1 + Drug + Error(Drug:Region1) + dose + dose:Drug , data = ex125 , projections = FALSE , qr = TRUE , contrasts = NULL # , ... ) summary(fm2.11) fm2.12 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.12) anova(object = fm2.12, ddf = "Satterthwaite") library(multcomp) Contrasts1 <- matrix(c( 1, 0.5, 0, 0 , 0, 0, -1, -0.5 , 1, 1, 0, 0 , 0, 1, 0, 0 ) , ncol = 4 , byrow = TRUE , dimnames = list( c("C1", "C2", "C3", "C4") , rownames(summary(fm2.12)$coef) ) ) Contrasts1 summary(glht(fm2.12, linfct=Contrasts1))
#------------------------------------------------------------- ## Example 2.5.3.1 p-70 #------------------------------------------------------------- # PROC GLM DATA=ex125; # CLASS drug dose region; # MODEL pcv=region drug region*drug dose drug*dose; # RANDOM region drug*region; # RUN; # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / ddfm=satterth; # RANDOM region drug*region; # ESTIMATE 'drug dif' drug -1 1 drug*dose -0.5 -0.5 0.5 0.5; # ESTIMATE 'Samorin mean' INTERCEPT 1 drug 0 1 dose 0.5 0.5 # drug*dose 0 0 0.5 0.5; # ESTIMATE 'Samorin HvsL' dose 1 -1 drug*dose 0 0 1 -1; # ESTIMATE 'Samorin high' INTERCEPT 1 drug 0 1 dose 1 0 # drug*dose 0 0 1 0; # RUN; library(lme4) str(ex125) ex125$Region1 <- factor(ex125$Region) fm2.11 <- aov( formula = Pcv ~ Region1 + Drug + Error(Drug:Region1) + dose + dose:Drug , data = ex125 , projections = FALSE , qr = TRUE , contrasts = NULL # , ... ) summary(fm2.11) fm2.12 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.12) anova(object = fm2.12, ddf = "Satterthwaite") library(multcomp) Contrasts1 <- matrix(c( 1, 0.5, 0, 0 , 0, 0, -1, -0.5 , 1, 1, 0, 0 , 0, 1, 0, 0 ) , ncol = 4 , byrow = TRUE , dimnames = list( c("C1", "C2", "C3", "C4") , rownames(summary(fm2.12)$coef) ) ) Contrasts1 summary(glht(fm2.12, linfct=Contrasts1))
Examp2.5.4.1 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.5.4.1 p-74 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / ddfm=satterth; # RANDOM region drug*region; # ESTIMATE 'Samorin mean' INTERCEPT 1 drug 0 1 dose 0.5 0.5 # drug*dose 0 0 0.5 0.5; # RUN; # PROC GLM DATA=ex125; # CLASS drug dose region; # MODEL pcv=region drug region*drug dose drug*dose; # ESTIMATE 'Samorin mean' INTERCEPT 1 drug 0 1 dose 0.5 0.5 # drug*dose 0 0 0.5 0.5; # RUN; library(lme4) str(ex125) ex125$Region1 <- factor(ex125$Region) fm2.13 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.13) library(multcomp) Contrasts2 <- matrix(c( 1, 0.5, 0, 0 ) , ncol = 4 , byrow = TRUE , dimnames = list( c("C5") , rownames(summary(fm2.13)$coef) ) ) Contrasts2 summary(glht(fm2.13, linfct=Contrasts2))
#------------------------------------------------------------- ## Example 2.5.4.1 p-74 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / ddfm=satterth; # RANDOM region drug*region; # ESTIMATE 'Samorin mean' INTERCEPT 1 drug 0 1 dose 0.5 0.5 # drug*dose 0 0 0.5 0.5; # RUN; # PROC GLM DATA=ex125; # CLASS drug dose region; # MODEL pcv=region drug region*drug dose drug*dose; # ESTIMATE 'Samorin mean' INTERCEPT 1 drug 0 1 dose 0.5 0.5 # drug*dose 0 0 0.5 0.5; # RUN; library(lme4) str(ex125) ex125$Region1 <- factor(ex125$Region) fm2.13 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.13) library(multcomp) Contrasts2 <- matrix(c( 1, 0.5, 0, 0 ) , ncol = 4 , byrow = TRUE , dimnames = list( c("C5") , rownames(summary(fm2.13)$coef) ) ) Contrasts2 summary(glht(fm2.13, linfct=Contrasts2))
Examp2.6.1 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 2.6.1 p-76 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / ddfm=satterth; # RANDOM region drug*region; # CONTRAST 'drug dif' drug -1 1 drug*dose -0.5 -0.5 0.5 0.5; # CONTRAST 'all' drug 1 -1 dose 0 0 drug*dose 0.5 0.5 -0.5 -0.5, # drug 0 0 dose 1 -1 drug*dose 0.5 -0.5 0.5 -0.5, # drug 0 0 dose 0 0 drug*dose 0.5 -0.5 -0.5 0.5; # RUN; library(lmerTest) str(ex125) ex125$Region1 <- factor(ex125$Region) fm2.14 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.14) anova(object = fm2.14, ddf = "Satterthwaite") library(multcomp) Contrasts3 <- matrix(c( 0, 0, -1, -0.5 ) , ncol = 4 , byrow = TRUE , dimnames = list( c("C1") , rownames(summary(fm2.14)$coef) ) ) Contrasts3 summary(glht(fm2.14, linfct=Contrasts3)) if(packageVersion("lmerTest") >= "3.0") contest(fm2.14, Contrasts3, joint = FALSE)
#------------------------------------------------------------- ## Example 2.6.1 p-76 #------------------------------------------------------------- # PROC MIXED DATA=ex125; # CLASS drug dose region; # MODEL pcv=drug dose drug*dose / ddfm=satterth; # RANDOM region drug*region; # CONTRAST 'drug dif' drug -1 1 drug*dose -0.5 -0.5 0.5 0.5; # CONTRAST 'all' drug 1 -1 dose 0 0 drug*dose 0.5 0.5 -0.5 -0.5, # drug 0 0 dose 1 -1 drug*dose 0.5 -0.5 0.5 -0.5, # drug 0 0 dose 0 0 drug*dose 0.5 -0.5 -0.5 0.5; # RUN; library(lmerTest) str(ex125) ex125$Region1 <- factor(ex125$Region) fm2.14 <- lmerTest::lmer( formula = Pcv ~ dose*Drug + (1|Region/Drug) , data = ex125 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose = "contr.SAS", Drug = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm2.14) anova(object = fm2.14, ddf = "Satterthwaite") library(multcomp) Contrasts3 <- matrix(c( 0, 0, -1, -0.5 ) , ncol = 4 , byrow = TRUE , dimnames = list( c("C1") , rownames(summary(fm2.14)$coef) ) ) Contrasts3 summary(glht(fm2.14, linfct=Contrasts3)) if(packageVersion("lmerTest") >= "3.0") contest(fm2.14, Contrasts3, joint = FALSE)
Examp3.1 is.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 3.1 Model 1 p-80 #------------------------------------------------------------- # PROC MIXED DATA=ex31; # CLASS drug dose herd; # MODEL PCV2=drug dose(drug)/solution ddfm=satterth; # RANDOM herd(drug); # ESTIMATE 'Mean Samorin' intercept 1 drug 0 1 dose(drug) 0 0 1; # ESTIMATE 'Berenil 2 doses' dose(drug) 1 -1 0; # ESTIMATE 'Ber vs Sam at dose 1' drug 1 -1 dose(drug) 1 0 -1; # CONTRAST 'Mean Samorin' intercept 1 drug 0 1 dose(drug) 0 0 1; # CONTRAST 'Berenil dif 2 doses' dose(drug) 1 -1 0; # CONTRAST 'Ber vs Sam at dose 1' drug 1 -1 dose(drug) 1 0 -l; # CONTRAST 'some difference' drug 1 -1 dose(drug) 0.5 0.5 -1, # drug 0 0 dose(drug) 1 -1 0; # LSMEANS dose(drug); # RUN; library(lmerTest) str(ex31) ex31$drug1 <- factor(ex31$drug) ex31$dose1 <- factor(ex31$dose) ex31$herd1 <- factor(ex31$herd) fm3.1 <- lmerTest::lmer( formula = PCV2 ~ drug1 + dose1:drug1 + (1|herd1:drug1) , data = ex31 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose1 = "contr.SAS", drug1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.1) anova(object = fm3.1, ddf = "Satterthwaite") lsmeansLT(model = fm3.1, test.effs = "dose1:drug1") #------------------------------------------------------------- ## Example 3.1 Model 2 p-84 #------------------------------------------------------------- # PROC MIXED DATA=ex31; # CLASS drug dose herd; # MODEL PCV2=PCV1 drug dose(drug)/solution ddfm=satterth; # RANDOM herd(drug); # RUN; library(lmerTest) str(ex31) ex31$drug1 <- factor(ex31$drug) ex31$dose1 <- factor(ex31$dose) ex31$herd1 <- factor(ex31$herd) fm3.2 <- lmerTest::lmer( formula = PCV2 ~ PCV1 + drug1 + dose1:drug1 + (1|herd1:drug1) , data = ex31 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose1 = "contr.SAS", drug1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.2) anova(object = fm3.2, ddf = "Satterthwaite") lsmeansLT(model = fm3.2, test.effs = "herd1:drug1") #------------------------------------------------------------- ## Example 3.1 Model 3 p-86 #------------------------------------------------------------- # PROC MIXED DATA=ex31; # CLASS drug dose herd; # MODEL PCV2=drug dose(drug) PCV1*dose(drug)/solution ddfm=satterth; # RANDOM herd(drug); # RUN; library(lmerTest) str(ex31) ex31$drug1 <- factor(ex31$drug) ex31$dose1 <- factor(ex31$dose) ex31$herd1 <- factor(ex31$herd) fm3.3 <- lmerTest::lmer( formula = PCV2 ~ drug1 + PCV1*dose1:drug1 + (1|herd1:drug1) , data = ex31 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose1 = "contr.SAS", drug1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.3) anova(object = fm3.3, ddf = "Satterthwaite") lsmeansLT(model = fm3.3, test.effs = "dose1:drug1")
#------------------------------------------------------------- ## Example 3.1 Model 1 p-80 #------------------------------------------------------------- # PROC MIXED DATA=ex31; # CLASS drug dose herd; # MODEL PCV2=drug dose(drug)/solution ddfm=satterth; # RANDOM herd(drug); # ESTIMATE 'Mean Samorin' intercept 1 drug 0 1 dose(drug) 0 0 1; # ESTIMATE 'Berenil 2 doses' dose(drug) 1 -1 0; # ESTIMATE 'Ber vs Sam at dose 1' drug 1 -1 dose(drug) 1 0 -1; # CONTRAST 'Mean Samorin' intercept 1 drug 0 1 dose(drug) 0 0 1; # CONTRAST 'Berenil dif 2 doses' dose(drug) 1 -1 0; # CONTRAST 'Ber vs Sam at dose 1' drug 1 -1 dose(drug) 1 0 -l; # CONTRAST 'some difference' drug 1 -1 dose(drug) 0.5 0.5 -1, # drug 0 0 dose(drug) 1 -1 0; # LSMEANS dose(drug); # RUN; library(lmerTest) str(ex31) ex31$drug1 <- factor(ex31$drug) ex31$dose1 <- factor(ex31$dose) ex31$herd1 <- factor(ex31$herd) fm3.1 <- lmerTest::lmer( formula = PCV2 ~ drug1 + dose1:drug1 + (1|herd1:drug1) , data = ex31 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose1 = "contr.SAS", drug1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.1) anova(object = fm3.1, ddf = "Satterthwaite") lsmeansLT(model = fm3.1, test.effs = "dose1:drug1") #------------------------------------------------------------- ## Example 3.1 Model 2 p-84 #------------------------------------------------------------- # PROC MIXED DATA=ex31; # CLASS drug dose herd; # MODEL PCV2=PCV1 drug dose(drug)/solution ddfm=satterth; # RANDOM herd(drug); # RUN; library(lmerTest) str(ex31) ex31$drug1 <- factor(ex31$drug) ex31$dose1 <- factor(ex31$dose) ex31$herd1 <- factor(ex31$herd) fm3.2 <- lmerTest::lmer( formula = PCV2 ~ PCV1 + drug1 + dose1:drug1 + (1|herd1:drug1) , data = ex31 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose1 = "contr.SAS", drug1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.2) anova(object = fm3.2, ddf = "Satterthwaite") lsmeansLT(model = fm3.2, test.effs = "herd1:drug1") #------------------------------------------------------------- ## Example 3.1 Model 3 p-86 #------------------------------------------------------------- # PROC MIXED DATA=ex31; # CLASS drug dose herd; # MODEL PCV2=drug dose(drug) PCV1*dose(drug)/solution ddfm=satterth; # RANDOM herd(drug); # RUN; library(lmerTest) str(ex31) ex31$drug1 <- factor(ex31$drug) ex31$dose1 <- factor(ex31$dose) ex31$herd1 <- factor(ex31$herd) fm3.3 <- lmerTest::lmer( formula = PCV2 ~ drug1 + PCV1*dose1:drug1 + (1|herd1:drug1) , data = ex31 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(dose1 = "contr.SAS", drug1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.3) anova(object = fm3.3, ddf = "Satterthwaite") lsmeansLT(model = fm3.3, test.effs = "dose1:drug1")
Examp3.2 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 3.3 p-88 #------------------------------------------------------------- # PROC MIXED DATA=ex32; # CLASS sex sire_id breed; # MODEL ww = sex agew breed/SOLUTION DDFM=SATTERTH; # RANDOM sire_id(breed)/SOLUTION; # LSMEANS breed/ADJUST = TUKEY; # RUN; library(lmerTest) str(ex32) ex32$sire_id1 <- factor(ex32$sire_id) ex32$breed1 <- factor(ex32$breed) fm3.4 <- lmerTest::lmer( formula = Ww ~ sex + agew + breed1 + (1|sire_id1:breed1) , data = ex32 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(sex = "contr.SAS", breed1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.4) anova(object = fm3.4, ddf = "Satterthwaite") lsmeansLT(model = fm3.4)
#------------------------------------------------------------- ## Example 3.3 p-88 #------------------------------------------------------------- # PROC MIXED DATA=ex32; # CLASS sex sire_id breed; # MODEL ww = sex agew breed/SOLUTION DDFM=SATTERTH; # RANDOM sire_id(breed)/SOLUTION; # LSMEANS breed/ADJUST = TUKEY; # RUN; library(lmerTest) str(ex32) ex32$sire_id1 <- factor(ex32$sire_id) ex32$breed1 <- factor(ex32$breed) fm3.4 <- lmerTest::lmer( formula = Ww ~ sex + agew + breed1 + (1|sire_id1:breed1) , data = ex32 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(sex = "contr.SAS", breed1 = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.4) anova(object = fm3.4, ddf = "Satterthwaite") lsmeansLT(model = fm3.4)
Examp3.3 is used for inspecting probability distribution and to define a plausible process through linear models and generalized linear models.
Muhammad Yaseen ([email protected])
Duchateau, L. and Janssen, P. and Rowlands, G. J. (1998).Linear Mixed Models. An Introduction with applications in Veterinary Research. International Livestock Research Institute.
#------------------------------------------------------------- ## Example 3.3 Model 1 p-88 #------------------------------------------------------------- # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # RANDOM animal_id(breed)/SOLUTION; # RUN; library(lme4) options(contrasts = c(factor = "contr.SAS", ordered = "contr.poly")) str(ex33) fm3.5 <- lme4::lmer( formula = PCV ~ breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.5) anova(fm3.5) library(lmerTest) fm3.6 <- lmerTest::lmer( formula = PCV ~ breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.6) anova(object = fm3.6, ddf = "Satterthwaite") # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # REPEATED/TYPE=CS SUB = animal_id(breed) R; # RUN; library(nlme) fm3.7 <- nlme::gls( model = PCV ~ breed + breed:time , data = ex33 , correlation = corCompSymm(, form = ~ 1|animal_id/breed) , weights = NULL # , subset = , method = "REML" # c("REML", "ML") , na.action = na.fail , control = list() ) summary(fm3.7) anova(fm3.7) # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = time breed breed*time/SOLUTION; # RANDOM animal_id(breed)/SOLUTION; # RUN; fm3.8 <- lme4::lmer( formula = PCV ~ time + breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.8) anova(fm3.8) fm3.9 <- lmerTest::lmer( formula = PCV ~ time + breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.9) anova(object = fm3.9, ddf = "Satterthwaite", type = 3) # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # REPEATED/TYPE=AR(1) SUBJET = animal_id(breed) R; # RUN; library(nlme) fm3.10 <- nlme::gls( model = PCV ~ breed + breed:time , data = ex33 , correlation = corAR1(, form = ~ 1|animal_id/breed) , weights = NULL # , subset = , method = "REML" # c("REML", "ML") , na.action = na.fail , control = list() ) summary(fm3.10) anova(fm3.10) # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # RANDOM INTERCEPT time/TYPE=UN SUBJET = animal_id(breed) SOLUTION; # RUN; library(nlme) # fm3.11 <- # nlme::gls( # model = PCV ~ breed + breed:time # , data = ex33 # , random = ~1|animal_id/breed # , correlation = corAR1(, form = ~ 1|animal_id/breed) # , weights = NULL # # , subset = # , method = "REML" # c("REML", "ML") # , na.action = na.fail # , control = list() # ) # summary(fm3.11) # anova(fm3.11)
#------------------------------------------------------------- ## Example 3.3 Model 1 p-88 #------------------------------------------------------------- # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # RANDOM animal_id(breed)/SOLUTION; # RUN; library(lme4) options(contrasts = c(factor = "contr.SAS", ordered = "contr.poly")) str(ex33) fm3.5 <- lme4::lmer( formula = PCV ~ breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.5) anova(fm3.5) library(lmerTest) fm3.6 <- lmerTest::lmer( formula = PCV ~ breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.6) anova(object = fm3.6, ddf = "Satterthwaite") # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # REPEATED/TYPE=CS SUB = animal_id(breed) R; # RUN; library(nlme) fm3.7 <- nlme::gls( model = PCV ~ breed + breed:time , data = ex33 , correlation = corCompSymm(, form = ~ 1|animal_id/breed) , weights = NULL # , subset = , method = "REML" # c("REML", "ML") , na.action = na.fail , control = list() ) summary(fm3.7) anova(fm3.7) # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = time breed breed*time/SOLUTION; # RANDOM animal_id(breed)/SOLUTION; # RUN; fm3.8 <- lme4::lmer( formula = PCV ~ time + breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.8) anova(fm3.8) fm3.9 <- lmerTest::lmer( formula = PCV ~ time + breed + breed:time + (1|animal_id:breed) , data = ex33 , REML = TRUE , control = lmerControl() , start = NULL , verbose = 0L # , subset # , weights # , na.action # , offset , contrasts = list(breed = "contr.SAS") , devFunOnly = FALSE # , ... ) summary(fm3.9) anova(object = fm3.9, ddf = "Satterthwaite", type = 3) # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # REPEATED/TYPE=AR(1) SUBJET = animal_id(breed) R; # RUN; library(nlme) fm3.10 <- nlme::gls( model = PCV ~ breed + breed:time , data = ex33 , correlation = corAR1(, form = ~ 1|animal_id/breed) , weights = NULL # , subset = , method = "REML" # c("REML", "ML") , na.action = na.fail , control = list() ) summary(fm3.10) anova(fm3.10) # PROC MIXED DATA=ex33; # CLASS breed animal_id; # MODEL pcv = breed breed*time/SOLUTION; # RANDOM INTERCEPT time/TYPE=UN SUBJET = animal_id(breed) SOLUTION; # RUN; library(nlme) # fm3.11 <- # nlme::gls( # model = PCV ~ breed + breed:time # , data = ex33 # , random = ~1|animal_id/breed # , correlation = corAR1(, form = ~ 1|animal_id/breed) # , weights = NULL # # , subset = # , method = "REML" # c("REML", "ML") # , na.action = na.fail # , control = list() # ) # summary(fm3.11) # anova(fm3.11)