2017-04-13 115 views
0

我正尝试使用nlme软件包在R中使用重复测量(MMRM)模型拟合混合模型。nlme:使用CSH协方差模型拟合混合模型

数据的结构如下: 每个患者属于三个组(grp)之一,并被分配到一个治疗组(trt)。 患者结果(y)在6次访问(访问)期间测量。

我想在不同访问中使用具有异构差异的复合对称模型(如SAS的PROC MIXED的CSH类型,https://support.sas.com/documentation/cdl/en/statug/63347/HTML/default/viewer.htm#statug_mixed_sect020.htm)。

为此,我使用lme中的相关性参数将相关结构设置为CS(corCompSymm)和权重参数,因此方差是访问的函数。

我也尝试添加访问corCompSymm本身的窗体参数。

我有这个问题:看起来我得到了相同的结果,不管我是否在lme调用中设置了权重参数(换句话说,看起来我得到CS模型而不是CSH模型)。

执行下面的代码,您会注意到模型参数估计的协方差矩阵的对角线无论使用什么模型都是相同的,这表明权重参数被忽略。

remove(list = objects()) 
library(nlme) 

set.seed(55) 

npatients  = 200; 
nvisits  = 6; 

#--- 
# Generate some data: 
subject_table = data.frame(subject = sprintf("S%03d", 1:npatients), 
          trt  = sample(x = c("P", "D"),  replace = T, size = npatients), 
          grp  = sample(x = c("A", "B", "C"), replace = T, size = npatients)) 
subject_table = merge(subject_table, 
         data.frame(visit.number = 1:6)) 
subject_table = transform(subject_table, 
          visit = sprintf("V%02d", visit.number), 
          y  = rnorm(nrow(subject_table), mean = 0, sd = visit.number^2)) 
subject_table = transform(subject_table, 
          visit = factor(visit), 
          subject = factor(subject, ordered = T, levels =  sort(unique(as.character(subject)))), 
          grp  = factor(grp), 
          trt  = factor(trt)) 
#--- 
# Fit MMRM model to data using nlme 
cs_model  = lme(y ~ trt*visit*grp,        # fixed  effects 
        random  = ~1|subject,      # random effects 
        data  = subject_table,     # data 
        correlation = corCompSymm(form=~1|subject))  # CS correlation matrix within patient 

csh_model_v1 = lme(y ~ trt*visit*grp,        # fixed effects 
        random  = ~1|subject,      # random effects 
        data  = subject_table,     # data 
        weights  = varIdent(~1|visit),    # different "weight" within each visit (I think) 
        correlation = corCompSymm(form=~1|subject))  # CS correlation matrix within patient 

csh_model_v2 = lme(y ~ trt*visit*grp,        # fixed effects 
        random  = ~1|subject,      # random effects 
        data  = subject_table,     # data 
        weights  = varIdent(~visit|subject),   # different "weight" within each visit (I think) 
        correlation = corCompSymm(form=~1|subject))  # CS correlation matrix within patient 

csh_model_v3 = lme(y ~ trt*visit*grp,        # fixed effects 
        random  = ~1|subject,      # random effects 
        data  = subject_table,     # data 
        correlation = corCompSymm(form=~visit|subject)) # CS correlation matrix within patient 

diag(vcov(cs_model)) 
diag(vcov(csh_model_v1)) 
diag(vcov(csh_model_v2)) 
diag(vcov(csh_model_v3)) 

问题: 如何获取NLME以适应不同的访问不同的变化的参数?

回答

0

经过几次死路,看来问题在于确保在调用varIdent时设置了正确的参数。

做正确的做法似乎是:

csh_model_right = lme(y ~ trt*visit*grp,       # fixed effects 
        random  = ~1|subject,     # random effects 
        data  = subject_table,    # data 
        weights  = varIdent(form=~1|visit),  # different "weight" within each visit (I know) 
        correlation = corCompSymm(),    # CS correlation matrix within subject per random statement above 
        control  = lme.control) 

它看起来是一样的,但是请注意,传递给varIdent参数被明确认定为“形式”。如果这种解释有任何其他的方式,我预计会发生崩溃,但我错了。