我在双变量纵向数据上运行模拟 . 我将描述我的代码，我希望减少组0的csse变化 . 因此，有两个类Y =（0,1），并且对于每个类我生成线性混合效应模型 . 每组有4个随机参数（截距和斜率），即协方差矩阵为4x4，随机效应之间存在相关性 . 此外，残差是矩阵2x2（每个变量都有误差）并且存在相关性 . 我假设每组中每个变量有4个固定效果参数 .
因此，我的代码是：
###~~~ simulate two correlated responses variables in longitudinal data ~~~###
## set means of random effects (intercepts and slopes)
m.0 = c(8,8,4,5)
m.1 = c(7,6,3,0)
## correlation matrix of random effects
cor.0 = matrix(c(1.00,0.20,0.55,0.29,
0.20,1.00,0.13,0.18,
0.55,0.13,1.00,0.65,
0.29,0.18,0.65,1.00),nrow=4)
cor.1 = matrix(c(1.00,0.30,0.40,0.56,
0.30,1.00,0.51,0.65,
0.40,0.51,1.00,0.81,
0.56,0.65,0.81,1.00),nrow=4)
### set correlation matrix of residuals
cor.R0 = matrix(c(1.00,0.15,
0.15,1.00),nrow=2)
cor.R1 = matrix(c(1.00,0.42,
0.42,1.00),nrow=2)
## generate covariance matrices
set.seed(1)
sds <- rnorm(4)^2
S0 = cor.0 * sds * rep(sds, each = nrow(cor.0))
sds <- rnorm(4)^2
S1 = cor.1 * sds * rep(sds, each = nrow(cor.1))
sds <- rnorm(2)^2
R0 = cor.R0 * sds * rep(sds, each = nrow(cor.R0))
sds <- rnorm(2)^2
R1 = cor.R1 * sds * rep(sds, each = nrow(cor.R1))
n = 200
## generate intercepte and slpos
library(mixAK)
B0 = rMVN(n,m.0,S0)$x
B1 = rMVN(n,m.1,S1)$x
E0 = rMVN(n,c(0,0),R0)$x
E1 = rMVN(n,c(0,0),R1)$x
Time = rep(c(0,3,6,9),times=n/2)
B0Time = B0[,c(2,4)]*Time[1:n]
B1Time = B1[,c(2,4)]*Time[1:n]
id = rep(1:n/2,each=4)
## set fix intercepts for each variable and groups
B0Fix1 = rep(0.5,n*2) #cs0
B1Fix1 = rep(40,n*2) #cs1
B0Fix2 = rep(40,n*2) #va0
B1Fix2 = rep(20,n*2) #va1
## get the equation
Y1.0 = B0Fix1 + B0[,1] + B0Time[,1] + E0[,1]
Y2.0 = B0Fix2 + B0[,3] + B0Time[,2] + E0[,2]
Y1.1 = B1Fix1 + B1[,1] + B1Time[,1] + E1[,1]
Y2.1 = B1Fix2 + B1[,3] + B1Time[,2] + E1[,2]
csse = c(Y1.0,Y1.1)
vase = c(Y2.0,Y2.1)
Y = as.factor(c(rep(0,400),rep(1,400)))
data.Sim = data.frame(id,Time,csse,vase ,Y)
摘要是：
summary(data.Sim)
summary(data.Sim[data.Sim$Y==1,])
summary(data.Sim[data.Sim$Y==0,])
> summary(data.Sim[data.Sim$Y==1,])
id Time csse vase Y
Min. : 50.50 Min. :0.00 Min. : 15.65 Min. :-47.71 0: 0
1st Qu.: 62.88 1st Qu.:2.25 1st Qu.: 55.83 1st Qu.: 11.38 1:400
Median : 75.25 Median :4.50 Median : 73.44 Median : 23.33
Mean : 75.25 Mean :4.50 Mean : 76.43 Mean : 21.97
3rd Qu.: 87.62 3rd Qu.:6.75 3rd Qu.: 95.21 3rd Qu.: 32.56
Max. :100.00 Max. :9.00 Max. :155.59 Max. : 98.10
> summary(data.Sim[data.Sim$Y==0,])
id Time csse vase Y
Min. : 0.50 Min. :0.00 Min. :-503.75 Min. : 27.17 0:400
1st Qu.:12.88 1st Qu.:2.25 1st Qu.: -18.55 1st Qu.: 51.49 1: 0
Median :25.25 Median :4.50 Median : 10.88 Median : 67.10
Mean :25.25 Mean :4.50 Mean : 49.01 Mean : 67.13
3rd Qu.:37.62 3rd Qu.:6.75 3rd Qu.: 122.14 3rd Qu.: 83.80
Max. :50.00 Max. :9.00 Max. : 564.92 Max. :125.68
只关注最后的摘要，我的问题是如何减少组0 Y==0 上不包括组1范围的csse的变化 . 请问有人有答案吗？谢谢