
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Example 1 - Exploratory Beta Model
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

source("Simulation Study 1/bsub_Sim1Aux.R")

out <- NULL
for(p0 in c(5,10,20)){
  for(n0 in c(200,500,1000,5000)){
    out <- rbind(out,outCppE1(n = n0, p = p0, plot = T))
  }
}

colnames(out) <- c("p","n","AvMSEm0","AvMSEm1","AvMSEs0","AvMSEs1",
                   "AvABm0","AvABm1","AvABs0","AvABs1",
                   "AvCOm0","AvCOm1","AvCOs0","AvCOs1","AvCT (m)")

# Table 5 (Section 4.1)
# ~~~~~~~~~~~~~~~~~~~~~
# xtable::xtable(out[,-c(1:2)], digits = 4)
round(out,4)

# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Example 2 - Confirmatory SN + Bernoulli model
# ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

source("Simulation Study 2/bsub_Sim2Aux.R")

out <- NULL
for(n0 in c(500,1000,3000)){
  out <- rbind(out,outCppE2(n = n0, plot = T))
}

colnames(out) <- c("n","AvMSEpi","AvMSEmu","AvMSEsg","AvMSEnu","AvMSErh",
                   "AvABpi","AvABmu","AvABsg","AvABnu","AvABrh",
                   "AvCOpi","AvCOmu","AvCOsg","AvCOnu","AvCOrh","AvCT (m)")
# round(out,4)

forpi = grep("^[Av].*pi$",colnames(out),ignore.case = T)
formu = grep("^[Av].*mu$",colnames(out),ignore.case = T)
forsg = grep("^[Av].*sg$",colnames(out),ignore.case = T)
fornu = grep("^[Av].*nu$",colnames(out),ignore.case = T)
forrh = grep("^[Av].*rh$",colnames(out),ignore.case = T)

# Table 6 (Section 4.2)
# ~~~~~~~~~~~~~~~~~~~~~
xtable::xtable(
  rbind(
    out[3,forpi],
    out[3,formu],
    out[3,forsg],
    out[3,fornu],
    out[3,forrh])
  , digits = 4)
