Careful Attention to Computer Arithmetic
Implement the scheme for Monte Carlo likelihood approximation given on the last slide of the computer arithmetic slides.
For test data use the data
This is a binary Rdata format file. The R commands
load(url("http://www.stat.umn.edu/geyer/8054/hw/hw10/hw10data.rda"))
ls()
should show that three R objects
-
psi, a numeric vector of length p, -
xobs, a numeric vector of length p, and -
x, a numeric matrix having p columns
-
psiis what is denoted ψ on the slides, -
xobsis what is denoted x on the slides, and - the rows of
xare what are denoted xi on the slides.
Task 1
Implement an R function that calculates the Monte Carlo approximation
to the log likelihood, its gradient vector, and its Hessian matrix
given a vector theta, a vector psi, a vector
xobs, and a matrix x.
Your code should make the appropriate input checks (all arguments are
numeric, all have all components is.finite, and the dimensions
are what they are supposed to be).
Your function should return a list with three components
-
value, the log likelihood value (a scalar) -
gradient, the log likelihood gradient (a p vector) -
hessian, the log likelihood Hessian (a p × p matrix)
Name them exactly as specified here. Such a list is what
the trust function wants.
Make theta the first argument to your function.
That is also what the trust function wants.
The Rules
You may implement this in R but
- No loops! Use
%*%(on-line help),apply(on-line help), andsweep(on-line help) instead. - No large objects! No temporary objects larger than
x. Ifxis n × p, this means you can have more n × p matrices and also p × p because we assume n > p. In particular, no n × n.
Task 2
Test the code at several random points not far from psi
by also approximating derivatives by finite differences. This test
code can use loops.
The for small vectors δ you should have
∇2 l(&theta) δ ≈ ∇ l(&theta + &delta) − ∇ l(&theta);
(angle brackets indicate inner product). When δ ranges over the standard basis vectors, this allows approximation of the components of the gradient vector and Hessian matrix.