Thanks David for stopping by this thread. I have to admit that the word
"arbitrary" was arbitrarily chosen. As I explained in the second email in
the thread, what I am really interested are some estimates obtained via
other methods that have the same asymptotic distribution as MLE. The
following simple example verifies what Thomas indicated: the vcov matrix
after one iteration is based on the original starting values, rather than
the one-step renewed values:
mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv";)
fit1 <- glm(admit ~ gre + gpa, data = mydata, family = "binomial")
set.seed(669)
beta0 <- jitter(coef(fit1)); beta0
control0 <- glm.control(epsilon=1e20, maxit=1, trace=F)
fit <- glm(admit ~ gre + gpa, data = mydata, family = "binomial",
start=beta0, control=control0)
cbind(beta0, one.step= fit$coef)
V0 <- summary(fit)$"cov.scaled"; V0
X <- as.matrix(cbind(1, mydata$gre, mydata$gpa));
eta <- as.vector(X%*%beta0)
expit <- function(x) (tanh(x/2)+1)/2
p <- expit(eta)
W <- diag(p*(1-p))
V0 <- solve(t(X)%*%W%*%X); V0
As you can see, the two V0s are identical, although sometimes they may be
different sheerly because of floating point rounding errors. I could have
obtained the variance-covariance matrix or Hessian using the known results
for each other type of GLMs. But thought using glm() function could make
the codes more generic. And it is working now.
=============================
Xiaogang Su, Ph.D.
Associate Professor
Department of Mathematical Sciences
University of Texas at El Paso
500 W. University Ave.
El Paso, Texas 79968-0514
x...@utep.edu
xiaogan...@gmail.com
https://sites.google.com/site/xgsu00/
On Wed, Feb 26, 2014 at 3:55 PM, David Winsemius <dwinsem...@comcast.net>wrote:
>
> On Feb 26, 2014, at 12:24 AM, Xiaogang Su wrote:
>
> > Dear All,
> >
> > Does anyone know if there is any way to obtain the variance-covariance
> > matrix for any arbitrarily given estimates with the glm() function?
> >
> > Here is what I really want. Given an arbitrary estimate (i.e., as
> starting
> > points with the start= argument), the glm() function could return the
> > corresponding variance-covariance matrix (or Hessian) and other
> quantities
> > with no Netwon-Raphson iteration?
>
> Not with an arbitrary start argument. You have offered a counter-example.
>
> > This could have been done by setting
> > maxit=0, but very unfortunately that is not an option in glm.control().
> To
> > illustrate the problem,
> >
> > mydata <- read.csv("http://www.ats.ucla.edu/stat/data/binary.csv";)
> > beta0 <- 1:3
> > control0 <- glm.control(epsilon = 1e10, maxit = 0, trace = FALSE)
>
> I admit to curiosity regarding how you interpreted the error message this
> generates:
>
> Error in glm.control(epsilon = 1e+10, maxit = 0, trace = FALSE) :
> maximum number of iterations must be > 0
>
> If you fix that problem you will be faced with another one induced by your
> unrealistic starting values.
>
> --
> David.
>
> > fit <- glm(admit ~ gre + gpa, data = mydata, family = "binomial",
> > start=beta0, control=control0)
> > summary(fit)$"cov.scaled"
> >
> > By the way, I am aware that I could directly compute the covariance
> matrix
> > using the formula. I also know that I could extract the corresponding
> > deviance by using the offset option.
> >
> > Any suggestion is greatly appreicated.
> >
> > Thanks,
> > Xiaogang Su
> >
> > =============================
> >
> >
> > [[alternative HTML version deleted]]
> >
> > ______________________________________________
> > R-help@r-project.org mailing list
> > https://stat.ethz.ch/mailman/listinfo/r-help
> > PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> > and provide commented, minimal, self-contained, reproducible code.
>
> David Winsemius
> Alameda, CA, USA
>
>
[[alternative HTML version deleted]]
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
