---
title: "An Aplication to SAE HB ME under Beta Distribution On Sample Data"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{An Aplication to SAE HB ME under Beta Distribution On Sample Data}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
```

## case 1 : Auxiliary variables only contains variable with error in aux variable

## Load Package and Load Data Set
```{r setup}
library(saeHB.ME.beta)
data("dataHBMEbeta")
```
## Fitting HB Model
```{r}
example <- meHBbeta(Y~x1+x2, var.x = c("v.x1","v.x2"),
                   iter.update = 3, iter.mcmc = 10000,
                   thin = 3, burn.in = 1000, data = dataHBMEbeta)
```

## Ekstract Mean Estimation

### Small Area Mean Estimation
```{r}
example$Est
```
### Coefficient Estimation
```{r}
example$coefficient
```
### Random Effect Variance Estimation
```{r}
example$refvar
```
### Extract MSE
```{r}
MSE_HBMEbeta=example$Est$sd^2
```
### Extract RSE
```{r}
RSE_HBMEbeta=sqrt(MSE_HBMEbeta)/example$Est$mean*100
```

## You can compare with Direct Estimation

### Extract Direct  Estimation
```{r}
Y_direct=dataHBMEbeta[,1]
MSE_direct=dataHBMEbeta[,6]
RSE_direct=sqrt(MSE_direct)/Y_direct*100
```
### Comparing Y
```{r}
Y_HBMEbeta=example$Est$mean
Y=as.data.frame(cbind(Y_direct,Y_HBMEbeta))
summary(Y)
```
### Comparing Mean Squared  Error (MSE)
```{r}
MSE=as.data.frame(cbind(MSE_direct,MSE_HBMEbeta))
summary(MSE)
```
### Comparing Relative Standard  Error (RSE) 
```{r}
RSE=as.data.frame(cbind(RSE_direct,RSE_HBMEbeta))
summary(RSE)
```

## case 2: Auxiliary variables contains variable with error and without error

## Fitting HB Model
```{r}
example_mix <- meHBbeta(Y~x1+x2+x3, var.x = c("v.x1","v.x2"),
                   iter.update = 3, iter.mcmc = 10000,
                   thin = 3, burn.in = 1000, data = dataHBMEbeta)
```

## Ekstract Mean Estimation

### Small Area Mean Estimation
```{r}
example_mix$Est
```
### Coefficient Estimation
```{r}
example_mix$coefficient
```
### Random Effect Variance Estimation
```{r}
example_mix$refvar
```
### Extract MSE
```{r}
MSE_HBMEbeta_mix=example_mix$Est$sd^2
```
### Extract RSE
```{r}
RSE_HBMEbeta_mix=sqrt(MSE_HBMEbeta_mix)/example_mix$Est$mean*100
```

## You can compare with Direct Estimation

### Comparing Y
```{r}
Y_HBMEbeta_mix=example_mix$Est$mean
Y_mix=as.data.frame(cbind(Y_direct,Y_HBMEbeta_mix))
summary(Y)
```
### Comparing Mean Squared  Error (MSE)
```{r}
MSE_mix=as.data.frame(cbind(MSE_direct,MSE_HBMEbeta_mix))
summary(MSE_mix)
```
### Comparing Relative Standard  Error (RSE) 
```{r}
RSE_mix=as.data.frame(cbind(RSE_direct,RSE_HBMEbeta_mix))
summary(RSE_mix)
```

###### note : you can use dataHBMEbetaNS for using dataset with non-sampled area