# # This file has functions for computing the UCB of Q(p.0) # using the bootstrap-t method. The function "ucb.boot" # gives the desired UCB, and the other functions are called # by this function. # # See the "example.R" file for an illustration. # # Date: June 24, 2005 ############## BEGIN CODE ############# # # UCB for Q(p.0) using bootstrap-t # ucb.boot <- function(n, p.0, mu.mle, sigma.mle, nboot=1999, alpha=0.05) { est <- logq.estimates(n, mu.mle, sigma.mle, p.0) logq.mle <- est$logq.hat logq.mle.se <- est$logq.hat.se t.resam <- replicate(nboot, t.sim(n, p.0, mu.mle, sigma.mle, logq.mle)) ucb <- exp( logq.mle - sort(t.resam)[alpha*(nboot+1)]*logq.mle.se) return(ucb) } # # simulate a resample of studentized t # t.sim <- function(n, p.0, mu.mle, sigma.mle, logq.mle) { d <- rnorm(n, mean=mu.mle, sd=sigma.mle) mu.hat <- mean(d) sigma.hat <- sqrt(var(d)*(n-1)/n) est <- logq.estimates(n, mu.hat, sigma.hat, p.0) tval <- (est$logq.hat - logq.mle)/est$logq.hat.se return(tval) } # # estimate of logq and its se # logq.estimates <- function(n, mu.hat, sigma.hat, p.0) { q.hat <- sigma.hat*sqrt(qchisq(p.0, df=1, ncp=(mu.hat/sigma.hat)^2)) q.l <- (-q.hat-mu.hat)/sigma.hat q.u <- (q.hat-mu.hat)/sigma.hat c1 <- dnorm(q.u)+dnorm(q.l) c2 <- dnorm(q.u)-dnorm(q.l) c3 <- q.u*dnorm(q.u)-q.l*dnorm(q.l) tau.hat <- sigma.hat*sqrt(c2^2 + 0.5*c3^2)/c1 logq.hat <- log(q.hat) logq.hat.se <- tau.hat/(sqrt(n)*q.hat) return(list(logq.hat=logq.hat, logq.hat.se=logq.hat.se)) } ############# END OF CODE ####################