Sure, I write the simulation in R since I also need to make calculations to obtain the true parameter. X1 is the dataset I input to SAS. cind = c(0,c2,c2+2,c2+4) ##threshold
varb <- matrix(c(sigma.b0^2,rhob*sigma.b0*sigma.b1,rhob*sigma.b0*sigma.b1,sigma.b1^2),nrow=2)
vare <- matrix(c(sigma.e0^2,rhoe*sigma.e0*sigma.e1,rhoe*sigma.e0*sigma.e1,sigma.e1^2),nrow=2)
k=2 ### 2 timepoints right now
n = 100 ##n samples
set.seed(i)
trt.pure = rbinom(n,1,2/3)
trt = rep(trt.pure,each=k)
time = rep(0:(k-1),times = n)
id = 1:n
X = data.frame(trt = trt, time = time/(k-1), id = rep(id,each=k))
mu1 = model.matrix(~time+trt,X)%*%c(alpha0,alpha1,alpha2)
mu2 = model.matrix(~time+trt,X)%*%c(beta0,beta1,beta2)
set.seed(i)
b <- rmvnorm(n,mean = c(0,0),sigma = varb)
eps <- rmvnorm(n*k,mean = c(0,0),sigma = vare)
y0 = mu1+ rep(b[,1],each=k) + eps[,1]
y1 = mu2+ rep(b[,2],each=k) + eps[,2]
y0 = floor(y0)
y0[y0>upper.limit] = upper.limit
y0[y0<lower.limit] = lower.limit
y1.ord <- dplyr::case_when(y1 < cind[1] ~"0",
y1>=cind[1] & y1< cind[2] ~ "1",
y1>=cind[2] & y1< cind[3] ~ "2",
y1>=cind[3] & y1< cind[4] ~ "3",
y1>=cind[4] ~"4")
X$continuous =y0
X$ordinal = as.numeric(y1.ord)+1
X1 = reshape2::melt(X,id.vars = c("trt","time","id"))
save(X1,file="....")
... View more
