vignettes/statg012-vignette.R

In yijin71/statg012: statg012: Statistical Inference at UCL

## ----set-options, echo=FALSE, cache=FALSE--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
options(width = 999)

## ---- include = FALSE----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
knitr::opts_chunk$set(comment = "#>", collapse = TRUE)

## ----setup, include=FALSE------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
devtools::load_all() # reload all code (after saving them)   or Ctrl-shift-L

## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
library(statg012)
ex <- binombeta(4, 6, 10, 3)

## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
summary(ex)

## ----fig.align='center',fig.height=3, fig.width=3.5----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
plot(ex, lty = 1:2, cex=0.5)

## ----fig.show = "hold", fig.height=3, fig.width=3, results = "hide"------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
uni <- binombeta(1, 1, 10, 3)
plot(uni, lty = 1:2, cex=0.5, main = "Uniform Prior", xlim=c(0,1), ylim=c(0,4))
plot(ex, lty = 1:2, cex=0.5, main = "Beta(4,6) Prior",xlim=c(0,1), ylim=c(0,4))

## ----fig.show = "hold", fig.height=3, fig.width=3, results = "hide"------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
be1 <- binombeta(0.5, 1, 10, 3)
be2 <- binombeta(0.5, 2, 10, 3)
plot(be1, lty = 1:2, cex=0.5, main = "Beta(0.5,1) Prior", xlim=c(0,1), ylim=c(0,4))
plot(be2, lty = 1:2, cex=0.5, main = "Beta(0.5,2) Prior",xlim=c(0,1), ylim=c(0,4))

## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
## generate a sample of 9 from a normal distribution with sd=2
x <- rnorm(9, sd = 2)
## given that the observed sample mean is 20
xx <- x- mean(x) + 20
## find the posterior density
exmp1 <- normnorm(xx, m = 25, s = sqrt(10), sigma = 2)

## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
summary(exmp1)

## ----fig.align='center',fig.height=3, fig.width=3.5----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
plot(exmp1, leg_pos = "right", cex = 0.5)

## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
## generate a sample of 9 from a normal distribution
y <- rnorm(9)
## Given parameters, find the posterior density 
exmp2 <- normnorm(y, m = 1, s = 2, a = 1, b = 1)
summary(exmp2)

## ----fig.show = "hold", fig.height=3, fig.width=3------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
## show the first plot : Prior and Posterior Distribution of mu
 plot(exmp2, which = 1, main = " Distribution of mu", cex = 0.5)

## show the second plot : Prior and Posterior Distribution of tau
 plot(exmp2, which = 2, main = " Distribution of tau", cex = 0.5)

## show the third plot : Prior Contour
plot(exmp2, which = 3, main = "Prior Contour",
       xlim = c(-5,5),ylim = c(0,5), cex = 0.5)

## show the fourth plot : Posterior Contour
plot(exmp2, which = 4, main = "Posterior Contour",
           xlim = c(-5,5), ylim = c(0,5), cex = 0.5)

## show the fifth plot : Prior and Posterior Contour
 plot(exmp2, which = 5, main = "Prior and Posterior Contour",
           xlim = c(-5,5), ylim = c(0,6), cex = 0.5)

## ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
## generate a sample of 16 from a normal distribution with sd=4
y <- rnorm(16, mean = 5.2, sd = 4)
## the observed sample mean is 5.2
yy <- y- mean(y)+5.2
## find the posterior density
exmp1 <- normnorm(yy, m = 4.5, s = sqrt(10), sigma = 4)

## make the hypothesis testing
normaltest(exmp1, 5, 1, 1)