R多層面分析案例

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

We first import the tab-delimited file, rat_pup.dat, assumed to be located in the C:\temp directory, into a data frame object in R named ratpup. The h = T argument in the read.table() function indicates that the first row (the header) in the rat_pup.dat file contains variable names. After reading the data, we “attach” the ratpup data frame to R’s working memory so that the columns (i.e., variables) in the data

frame can be easily accessed as separate objects. Note that we show the “>” prompt for each command as it would appear in R, but this prompt is not typed as part of the commands.
> ratpup <- read.table("c:\\temp\\rat_pup.dat", h = T)
> attach(ratpup)
To facilitate comparisons with the analyses performed using the other software procedures, we recode the variable SEX into SEX1, which is an indicator variable for females, so that males will be the reference group:
> ratpup$sex1[sex == "Female"] <- 1
> ratpup$sex1[sex == "Male"] <- 0Step 1: Fit a model with a “loaded” mean structure (Model 3.1). We first load the nlme package,* so that the lme() function will be available for model fitting:
> library(nlme)
We fit the initial LMM, Model 3.1, to the Rat Pup data using the lme() function:
> # Model 3.1.
> model3.1.fit <- lme(weight ~ treatment + sex1 + litsize +treatment*sex1, random = ~1 | litter, ratpup, method = "REML")

We next explain each part of the syntax used for the lme() function:

• model3.1.fit is the name of the object that will contain the results of the fitted model.
• The first argument of the function, weight ~ treatment + sex1 + litsize+ treatment*sex1, defines a model formula. The response variable, WEIGHT, and the terms that have associated fixed effects in the model (TREATMENT,SEX1, LITSIZE, and the TREATMENT×SEX1 interaction), are listed.
• The factor() function is not necessary for the categorical variable TREATMENT, because the original treatment variable has string values High, Low, and Control, and will therefore be considered as a factor automatically. We also do not need to declare SEX1 as a factor, because it is an indicator variable having only values of 0 and 1.
• The second argument, random = ~1 | litter, includes a random effect for each level of LITTER in the model. These random effects will be associated with the intercept, as indicated by ~1.
• The third argument of the function, ratpup, indicates the name of the data frame object to be used in the analysis.
• The final argument, method = "REML", specifies that the default REML estimation method is to be used.
  

By default, the lme() function treats the lowest level (alphabetically or numerically) of a categorical fixed factor as the reference category. This means that “Control” will be the reference category of TREATMENT because “Control” is the lowest level of treatment alphabetically.
We obtain estimates from the model fit by using the summary() function:
> summary(model3.1.fit)
Additional results for the fit of Model 3.1 can be obtained by using other functions in conjunction with the model3.1.fit object. For example, we can obtain F-tests for the fixed effects in the model by using the anova() function:
> anova(model3.1.fit)
The anova() function performs a series of Type I (or sequential) F-tests for the fixed effects in the model, each of which are conditional on the preceding terms in the model specification. For example, the F-test for SEX1 is conditional on the TREATMENT effects, but the F-test for TREATMENT is not conditional on the SEX1 effect. The random.effects() function can be used to display the EBLUPs for the random litter effects:
> # Display the random effects (EBLUPs) from the model.
> random.effects(model3.1.fit)
R01.pdf (1.88 MB)

2014-6-3 20:31:00 上传

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享1 收藏1 回帖

关键词：分析案例 fixed effect conditional comparisons Additional function separate reading objects memory

挺详细，还是不容易做到

原来在这儿