# R/data_sim.R In wgeesel: Weighted Generalized Estimating Equations and Model Selection

```data_sim <- function(id, rho,phi, x, beta, x_mis, para, corstr, family,lag_level){

switch(family,
"gaussian"={data <- cont_fun(beta,rho,phi,x,id,corstr)
},
"binary"={data <- binary_fun(beta,rho,x,id,corstr)
},
"poisson"={data <- count_fun(beta,rho,x,id,corstr)
},
stop("Warnings: Invalid type!")
)
complete=F
x_final=list()
prob_miss=0
if(complete==F){
data\$ind <- 1
#data\$prey_mis <- 0
data\$y_first <- NULL;
N=length(unique(id))
for (i in 1:N){
nsubj=max(table(id))
for (k in 2:nsubj){
lagy=NULL
if(lag_level!=0){
for(lag_i in 1:lag_level){
ylagname=paste("ylag",lag_i,sep="")
y_i=data\$response[which(id==i)]
y_new=c(rep(NA,lag_i),y_i[1:(length(y_i)-lag_i)])
lagy=cbind(lagy,y_new)
colnames(lagy)[lag_i]=ylagname
}
}

x_mis_i=cbind(x_mis[which(id==i),],lagy=lagy)
x_final[[i]]=x_mis_i

p <- 1/(1+exp(sum(-x_mis_i[k,]*para,na.rm=T)))
###depends on the first and second observation;
ind_p <- rbinom(1,1,p)
data[data\$id==i,]\$ind[k] <- ind_p
if (data[data\$id==i,]\$ind[k-1]==0){
data[data\$id==i,]\$ind[(k):nsubj] <- rep(0,(nsubj-(k)+1))
}
}
data\$y_first[(i*nsubj-(nsubj-1)):(i*nsubj)]=rep(data[data\$id==i,]\$response[1],nsubj)
}
x_final=do.call("rbind",x_final)
response_mis=ifelse(data\$ind==0,NA,data\$response)
data_final=cbind(data[,c("id","response","ind")],response_mis,x,x_final)
data<-as.data.frame(data_final)
prob_miss<-length(which(data\$ind==0))/length(data\$response)
}
return(list(data=data, prob_miss=prob_miss))
}
```

## Try the wgeesel package in your browser

Any scripts or data that you put into this service are public.

wgeesel documentation built on May 2, 2019, 3:41 a.m.