Hi All,
Sorry to bother you. I'm trying to estimate a set of discrete choice
data in R with mixed logit models where one coefficient is random and
normally distributed. I've searched on the R help archive and don't see
much
information very specific to what I'm doing, so I write the code myself,
which involves simulated maximum likelihood. But it doesn't work, as I
compare the results with Limdep. Could anyone please take a look at it? It's
a bit longer though......Comments are highly appreciated. Thank you very
much!
/*1,000 people, Three alternatives, two variables, their coefficients are
B1 and B2 respectively, B1 is random*/
library(maxLik)
RP<-function(theta,y,X)
{
m1<-theta[1] /* mean
s1<-theta[2] /*standard deviations
b2<-theta[3] /*B2
P<-NULL
b1<-rnorm(500,mean=m1,sd=s1) /*generate 500 random draws for B1
for(m in 0:999)
{
Dm<-X[(1+3*m):(3+3*m),] /*Extract the data for one person
Pn<-NULL
for(n in 1:500)
{
b<-rbind(b1[n],b2)
an<-sum(exp(Dm%*%b))
Pmn<-exp(Dm%*%b)/an /*Under each B1, compute the choice
probabilities
Pn<-cbind(Pn,Pmn)
}
Pm<-rowMeans(Pn) /* The simulated probabilities for one person
P<-rbind(P,Pm) /* Obtain the choice probabilities for all 1,000 people
}
sum(log(P)*(as.numeric(y))) /* Log-likelihood function, where y is the
variable for choices
}
A<-matrix(c(0,1,0),1,3)
B<-0
rp<-maxBFGS(RP,start=c(-0.072,0.8,0.539),y=Y,X=EV,constraints=list(ineqA=A,ineqB=B))
/* Constrained MLE, to ensure that the estimated standard deviation is
positive; the start values are taken from conditional logit estimation
--
Min Chen
Graduate Student
Department of Agricultural, Food, and Resource Economics
208 Cook Hall
Michigan State University
chenmin5@msu.edu / chenmin0905@gmail.com
[[alternative HTML version deleted]]
