University of Minnesota, Twin Cities School of Statistics Stat 8701 Rweb

Statistics 8701 (Geyer) Spring 2002 Homework No. 8

The policy is the same as the policy for Homework No. 2.

This assignment is due (initially) Friday, April 19, 2002.

Collateral Reading

Another course I taught fall semester covered the bootstrap. It had web pages with some examples that may (or may not) be helpful (I'll leave that to you do decide). The web pages relevant to the material involved in this assignment are

Problem 1

The ts library in R contains a times series treering (on-line help). It can be view with the following code

library(ts)
data(treering)
plot(treering)

Suppose I want to fit an AR(4) model to this series. For those who don't know, that is a model of the form

X_i = r₁ X_{i - 1} + r₂ X_{i - 2} + r₃ X_{i - 3} + r₄ X_{i - 4} + Z_i

where the Z_i are i. i. d. normal with mean zero (the usual sort of errors encountered in statistics). The following code fits the model

library(ts)
data(treering)
ar(treering, aic = FALSE, order.max = 4)

(if you try without the optional arguments you will see that it uses AIC to choose which AR(k) model it fits and that it choses a different model from the one I did, but that is the subject for another homework).

Assume this time series is stationary and that the estimators of the AR coefficients obey the square root law and determine conficence intervals for the AR coefficients using the subsampling bootstrap with subsample size 100 (if everybody uses 100 our results will all agree).

Problem 2

The file

http://www.stat.umn.edu/geyer/5601/mydata/gamma.txt

contains 200 observations from a gamma distribution. The problem is to estimate the shape parameter using the method of moments estimator

data <- read.table(url("http://www.stat.umn.edu/geyer/5601/mydata/gamma.txt"),
    header = TRUE)
attach(data)
hist(x)
alpha.hat <- mean(x)^2 / var(x)
print(alpha.hat)

The asymptotic variance of this estimator is 2 * alpha * (1 + alpha) / n where alpha is the true unknown parameter value.

The variance stabilizing transformation for alpha.hat and its inverse are given by the following

varstab <- function(alpha) sqrt(2) * asinh(sqrt(alpha))
varstabinv <- function(beta) sinh(beta / sqrt(2))^2

The random variable varstab(alpha.hat) is asymptotically standard normal (mean varstab(alpha), variance one).

Calculate the following confidence intervals for the true unknown shape parameter

The asymptotic confidence interval based on the delta method.
The asymptotic confidence interval based on the variance-stabilizing transformation.
The stupid bootstrap-t confidence interval that uses
```
sdfun <- function(x, nbootsd, theta, ...) return(1)
```
for the standardization.

The (non-stupid) bootstrap-t confidence interval that uses

sdfun <- function(x, nbootsd, theta, ...) {
     alpha.hat <- mean(x)^2 / var(x)
     return(sqrt(2 * alpha.hat * (1 + alpha.hat)))
 }

for the standardization.

The (non-stupid) bootstrap-t confidence interval that does not provide an sdfun and uses the double bootstrap to calculate one.
The (non-stupid) bootstrap-t confidence interval that does uses the double bootstrap to calculate a variance-stabilizing transformation.
The bootstrap percentile interval.
The ABC (or alphabet soup, type I) interval.
The BC_a (or alphabet soup, type II) interval.

In each case use 95% for the confidence level.

The default bootstrap and double bootstrap sample sizes in these functions are too small, in keeping for Efron and Tibshirani's general attitude. As explained on the web page about bootstrap-t intervals for the other class. Use nboott = 1000 and nbootsd = 200 or larger in your calculations. Otherwise the intervals will be so random as to be ungradable.

Also use bootstrap sample sizes of 1000 or more when doing the percentile and BC_a intervals.