Introduction
Suppose you want to fit a linear regression model
where
is a nx1 vector of response variables,
is an nxp design matrix,
is a px1 coefficient vector and
is a nx1 vector of independent
random variables. Suppose in addition that you cannot observe all of the elements of X directly: instead, you observe V=X+U where the elements of U are themselves independent normal random variables with zero mean. The variance of the (i,j)th element of U is
: note that this does not depend on i. One way to interpret this is: the entries in the jth column of the design matrix are all subject to independent measurement errors, which are normally distributed with variance
.
In this situation, to fit the regression model it is natural to use the usual least-squares procedure with the observed design matrix V in place of the true values X. It can be shown, however, that this produces biased estimates of the coefficient vector β. A corrected estimator can be defined, however, as
where D is a pxp diagonal matrix with jth diagonal element
. The covariance matrix of this estimate is
Standard errors of the estimated parameters can then be calculated as the square roots of the diagonal elements of this covariance matrix. For this question, we will assume that the value of
is known.
Your task
Without using the R command lm() or anything similar, write an R function called lm.adjusted(), to implement the corrected estimator as described above. The arguments to your function should be y, a vector of responses corresponding to Y in the description above; V, a design matrix corresponding to V in the description; sigsq, the known value of
; and tausq, a vector containing the values of the
. The default value fortausq should be a vector of zeroes. You may assume that the matrix V'V-D is invertible.
Your function should return a list containing components betahat (the estimated coefficient vector) and sebeta (a vector of estimated standard errors for the coefficients).