In the iterative weighted least squares loop of the glm.fit function,
a QR decomposition of the matrix X*sqrt(w) is obtained at each iteration,
using fitted values of the mean response mu from the previous interation.
Here "X" is the design matrix, w = wt*(d(mu)/d(eta))/variance(mu),
"wt" is the vector of weights, and "d(mu)/d(eta)" is the
derivative
of mu wrt the linear predictor eta. Upon convergence, the most recent
QR decomposition of X*sqrt(w) is returned and this is later used by
the function summary.glm to compute standard errors.
The QR decomposition being returned by glm.fit uses the fitted values from
parameter estimates at the (final-1)st iteration. But don't we want
fitted values from parameter estimates at the final iteration?
I tried inserting the following lines upon completion of the iterative
weighted least squares loop:
mu.eta.val <- mu.eta(eta)
w <- sqrt((weights[good] * mu.eta.val[good]^2)/variance(mu)[good])
newQR<-qr(x[good, ] * w)
My standard errors based on newQR$qr agree with the standard
errors I had computed manually from the final fits.
Has anyone else encountered this issue? Thanks in advance
for any clarifications.
Brad
--please do not edit the information below--
Version:
platform = i686-pc-linux-gnu
arch = i686
os = linux-gnu
system = i686, linux-gnu
status =
major = 1
minor = 7.0
year = 2003
month = 04
day = 16
language = R
Search Path:
.GlobalEnv, package:methods, package:ctest, package:mva, package:modreg,
package:nls, package:ts, Autoloads, package:base