0

I am running simulation on bivariate longitudinal data. I will describe my code which I hope to reduce variation on csse for group 0. So, there is tow classes Y= (0,1) and for each class I generate linear mixed effects model. There are 4 random parameters (intercept and slopes ) for each group, that is mean the covariance matrix is 4x4 and there are correlated between random effects. Also, the residuals are matrix 2x2 (each variable has error) and there is a correlation . I assume that there are 4 fixed effects parameters for each variable in each group.
Therefore, my code is:

###~~~ 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)

The summaries are:

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    

Just focus on the last summaries, which my question is how to reduce the variation on csse on group 0 Y==0that not includes the range of group 1 . May any one has answer, please? Thanks

R. Saeiti
  • 89
  • 3
  • 11
  • I ma surprised!! that No one has answer!!, – R. Saeiti Oct 08 '15 at 11:23
  • In my humble opinion, it is not very clear what you are trying to do. Ok, you want to reduce the variation on csse for a certain group, but why? and by how much? Anyway, take into account that only 20 people have taken a look to it so far. Maybe for somebody that knows more statistics than me it is very clear what you are asking and you get your answer soon. Good luck. – lrnzcig Oct 10 '15 at 17:34

0 Answers0