TU Wien:Statistik und Wahrscheinlichkeitstheorie UE (Levajkovic)/Übungen 2023W/HW10.5

Simulation of coverage probability

Simulate the coverage probability of the one-sample confidence interval for frequencies - does the confidence interval deliver what it promises? Let ${\textstyle Y_{1},\dots ,Y_{n}}$ be i.i.d. random variables with ${\textstyle Y_{1}\sim bern(p)}$ and ${\textstyle p\in (0,1)}$. Approximate in 10000 simulations the coverage probability of the 95%-confidence interval, i.e., simulate the proportion of coverage events of the parameter p. For that let

(a) ${\textstyle n=45}$ and ${\textstyle p=4/9}$

(b) ${\textstyle n=10}$ and ${\textstyle p=1/10}$

(c) Visualize the simulated relative frequencies from (a) and (b) in two histograms and

comment on the simulated coverage probabilities.

Hint: R-command rbinom and for example ifelse().

Dieses Beispiel ist als solved markiert. Ist dies falsch oder ungenau? Aktualisiere den Lösungsstatus (Details: Vorlage:Beispiel)

Lösungsvorschlag von Friday

--Friday Sa 30 Jan 2021 18:30:42 CET

# Statistics and Probability: HW #12
# Friday 
# Duedate: 18.01.2021

# 6) Simulation of coverage probability
# Simulate the coverage probability of the one-sample confidence interval for 
# frequencies - does the confidence interval deliver what it promises? Let Y1
# ,...,Yn be i.i.d. random variables with Y1 ∼ ber(p) and p ∈ (0,1). 
# Approximate in 10000 simulations the coverage probability of the 
# 95%-confidence interval, i.e., simulate the proportion of coverages of the 
# parameter p. 

plot_sim <- function(n, p) {
  hs<-c()
  inside <- 0
  q <- qnorm(0.025, mean=0, sd=1) *-1
  H <- p
  
  for (i in 1:10000) {
    ys <- rbinom(n, 1, p)
    h <- sum(ys)/n
    hs <- c(hs, h)
    seh <- sqrt((h*(1-h))/n)
    I <- c(h-q*seh,h+q*seh)
    if (H > I[1] && H < I[2]) {
      inside <- inside + 1
    }
  }
  hist(hs, main=sprintf("Histogram n=%d, p=%.2f",n, p), xlab="p-value",col="lightblue")
  inside/10000
}

# 6a) 
# For that let n=45 and p=4/9.
# 6b) 
# For that let n=10 and p=1/10.
# 6c)
# Visualize the simulated relative frequencies from (a) and (b) in a histogram 
# each and comment on the simulated coverage probabilities.
(par(mfrow=c(2,1)))
plot_sim(45, 4/9)
plot_sim(10, 1/10)