#First of all we load and attach the data:
data(mldata)
attach(mldata)
#Then we define all the inputs:
Y<-data.frame(measure,age)
X<-data.frame(rep(1,1000),sex)
Z<-data.frame(rep(1,1000))
clus<-data.frame(city)
betap<-matrix(0,2,2)
up<-matrix(0,10,2)
covp<-diag(1,2)
covu<-diag(1,2)
Sp=diag(1,2);
nburn=as.integer(200);
nbetween=as.integer(200);
nimp=as.integer(5);
Sup=diag(1,5);
#And finally we run the imputation function:
imp<-jomo1rancon(Y,X,Z,clus,betap,up,covp, covu,Sp,Sup,nburn,nbetween,nimp)
#Then we run the model on the imputed datsets:
estimates<-rep(0,5)
ses<-rep(0,5)
estimates2<-rep(0,5)
ses2<-rep(0,5)
for (i in 1:5) {
dat<-imp[imp$Imputation==i,]
fit<-lm(measure~age+sex+factor(clus),data=dat)
estimates[i]<-coef(summary(fit))[2,1]
ses[i]<-coef(summary(fit))[2,2]
estimates2[i]<-coef(summary(fit))[3,1]
ses2[i]<-coef(summary(fit))[3,2]
}
#And finally we aggregate results through RUbin's rules with package BaBooN.
#library("BaBooN")
#MI.inference(estimates, ses^2)
#MI.inference(estimates2, ses2^2)
Run the code above in your browser using DataLab