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)
Thanks,
Stephen Collins, MPP | Analyst
Global Strategy | Aon Benfield
[[alternative HTML version deleted]]
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.dalgaard at biostat.ku.dk) FAX: (+45) 35327907