# \donttest{
#Data Generation
set.seed(123)
m=30
x1=runif(m,0,1)
x2=runif(m,0,1)
x3=runif(m,0,1)
x4=runif(m,0,1)
b0=b1=b2=b3=b4=0.5
u=rnorm(m,0,1)
pi=rgamma(1,1,0.5)
Mu <- exp(b0+b1*x1+b2*x2+b3*x3+b4*x4+u)/(1+exp(b0+b1*x1+b2*x2+b3*x3+b4*x4+u))
A=Mu*pi
B=(1-Mu) * pi
y=rbeta(m,A,B)
y <- ifelse(y<1,y,0.99999999)
y <- ifelse(y>0,y,0.00000001)
MU=A/(A+B)
vardir=A*B/((A+B)^2*(A+B+1))
dataBeta = as.data.frame(cbind(y,x1,x2,x3,x4,vardir))
dataBetaNs=dataBeta
dataBetaNs$y[c(3,14,22,29,30)] <- NA
dataBetaNs$vardir[c(3,14,22,29,30)] <- NA
##Compute Fitted Model
##y ~ x1 +x2
## For data without any nonsampled area
vc = c(1,1,1)
c = c(0,0,0)
formula = y~x1+x2
dat = dataBeta[1:10,]
##Using parameter coef and var.coef
saeHBbeta<-Beta(formula,var.coef=vc,iter.update=10,coef=c,data=dat)
saeHBbeta$Est #Small Area mean Estimates
saeHBbeta$refVar #Random effect variance
saeHBbeta$coefficient #coefficient
#Load Library 'coda' to execute the plot
#autocorr.plot(saeHBbeta$plot[[3]]) # is used to generate ACF Plot
#plot(saeHBbeta$plot[[3]]) # is used to generate Density and trace plot
##Do not use parameter coef and var.coef
saeHBbeta <- Beta(formula,data=dataBeta)
# }
Run the code above in your browser using DataLab