Stephen Collins wrote:
Say I have a formula Y ~ 1 + X, where X is a categorical variable. A
previous thread showed how to evaluate this model using the mle package
from "stats4" (see below). But, the user had to create the data matrix,
X, including the column of one's for the regression constant. Is there a
way to nest the linear formula in the code below, so the data matrix
doesn't explicitly have to be created by the user?
Y <- c(0,0,1,0,0,1,1,0,0,0,0,1,1,0,1,1,0,1,1,0,1)
X <-
cbind(matrix(1,21,1),matrix(c(-48.5,24.4,82.8,-24.6,-31.6,91.0,52.1,-87.7,-17.0,-51.5,
-90.7,65.5,-44.0,-7.0,51.6,32.4,-61.8,34.0,27.9,-72.9,49.9), 21,1))
log.lo.like <- function(beta,Y,X) {
Fbetax <- 1/(1+exp(-beta%*%t(X)))
loglbeta <- -log(prod(Fbetax^Y*(1-Fbetax)^(1-Y)))
}
#####Using MLE#####
ll <- eval(function(beta0=0,beta1=0)
log.lo.like (c(beta0,beta1),Y,X),
list(X=X,Y=Y))
summary(mle(ll))
####Comparison using glm#####
glm(Y~X-1,family=binomial)
You might get something out of looking at
http://staff.pubhealth.ku.dk/~pd/slides/RSS-sep08.pdf
The code on the last of the three p.12 slides (don't ask) shows how to
go from a model formula to a likelihood function for a GLM in a fairly
automatic way.
--
O__ ---- Peter Dalgaard Ă˜ster Farimagsgade 5, Entr.B
c/ /'_ --- Dept. of Biostatistics PO Box 2099, 1014 Cph. K
(*) \(*) -- University of Copenhagen Denmark Ph: (+45) 35327918
~~~~~~~~~~ - (p.dalga...@biostat.ku.dk) FAX: (+45) 35327907
______________________________________________
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.