Statistics 8701 (Geyer, Spring 2003) Homework No. 2

Policy

This assignment will be handed in as a regular paper assignment, although you can hand in code by giving a URL where it can be found on the web or a location in your home directory on the School of Statistics linux boxes.

For the programming problems, you should hand in enough so that it is clear exactly what you did.

The rules for grading this stuff is that only perfect is passing. If I find anything wrong, you have to debug the problem and hand in a fix (for that problem or part of a problem). And the process will be repeated until everything is correct. Thus everyone will eventually do every assignment perfectly (if they keep working at it).

You are not in competition with each other for grades. There's no reason why everyone shouldn't get an A. If you don't get A, that's because you either failed to complete some assignments or took some huge number of iterations to complete some.

This assignment is due (initially) Friday, February 7, 2002.

Problem 1

Write a C function myrunif and a corresponding R function does the same as the R function runif. That is, your function should have the prototype

void myrunif(double *x, int *n)

and each call should set x to be a vector of n[0] uniform random numbers, so that, after you have dynamically loaded myrunif.so, the R commands

n <- 7
invisible(runif(1))
save.seed <- .Random.seed
runif(n)
print(.Random.seed)
.Random.seed <- save.seed
.C("myrunif", x = double(n), n = as.integer(n))$x
print(.Random.seed)

should produce exactly the same results.

Section 5.3 in the Writing R Extensions book tells how to do this.

Problem 2

Write a C function rar that simulates an AR(k) time series. This reimpliments part of what the function arima.sim function in the ts package of R does.

The result should be a time series having the form

Xi = ρ1 Xi - 1 + . . . + ρk Xi - k + ei

where the ei are i. i. d. Normal(0, &sigma2) random variables.

Your R function should have the prototype

rar(rho, sigma, n, nburn, x)

where

Note that nburn is unnnecessary if x is supplied and vice versa, but we allow both to be supplied and do what the description describes when both are supplied.

For a first pass you can write this entirely in R, though we will eventually write the inner loop in C, so you can just go straight to C if you like.

Problem 3

Write your own R package with help page for your rar function.

Some helpful R commands are the following. Suppose you have a directory rar that is a subdirectory of the current working directory and that is your package (it has a file DESCRIPTION and subdirectories R, man, and src). Then

cd rar # that is, go into the package directory
R CMD Rdindex man > INDEX
cd .. # and go back up

creates the INDEX file so you don't have to do it yourself.

Then

R CMD check rar

checks your package, installs it in the directory rar.Rcheck so you can just use the package right there, say with

R
library(rar, lib.loc = "rar.Rcheck")

and so forth.

You can also do a regular install elsewhere. For example,

R CMD INSTALL -l ~/library rar

installs the package in ~/library (which is the library subdirectory of your home directory). Then use it with

Problem 4

Show that with the q functions (qbinom, qpois) defined the way they are in R code like

qbinom(runif(nsim), n, p)

simulates random variates having the desired distribution.

Generally, for any distribution there exists a definition of the quantile function Q (what the R q functions calculate) such that Q(U) is a random variate having this distribution where U is Uniform(0, 1).