### Multinomial Logit model estimation
### Create the data array from the existing dataset
Data <- read.csv("C:/Users/1/Desktop/R sample code/R sample code/MNL/Data_Clean_English.csv",header=TRUE)
## Count the number of data rows
hh <- nrow(Data)
## Set the initial values of the parameters (the number in parenthesis corresponds to the number of parameters to be estimated)
b0 <- numeric(5)
## Recall that we estimate parameters by maximizing the log-likelihood function, hence we first need to define that function
##### Define the log-likelihood function of the logit model#####
fr <- function(x) {
### declare the parameters###
## Alternative specific constants
b1 <- x[1]
b2 <- x[2]
b3 <- x[3]
b4 <- x[4]
## Travel time to destination
d1 <- x[5]
## declare the log-likelihood variable, set value to 0
LL = 0
### For this choice problem, we consider the following 5 modes???
## calculate the utility function: :introduce the desired explanatory variables in the function
train <- Data$ModeAvailableTrain*exp(d1*Data$TotalTimeTrain/100 +b1*matrix(1,nrow =hh,ncol=1))
bus <- Data$ModeAvailableBus *exp(d1*Data$TotalTimeBus/100 +b2*matrix(1,nrow =hh,ncol=1))
car <- Data$ModeAvailableCar *exp(d1*Data$TimeCar/100 +b3*matrix(1,nrow =hh,ncol=1))
bike <- Data$ModeAvailableBike *exp(d1*Data$TimeBike/100 +b4*matrix(1,nrow =hh,ncol=1))
walk <- Data$ModeAvailableWalk *exp(d1*Data$TimeWalk/100 )
### Calculate the choice probabilities
## calculate the Inclusive Value (the denominator of the choice probabilities equation)
deno <- (car + train + bus + bike + walk)
## Calculate indiviudal choice probabilities
Ptrain <- Data$ModeAvailableTrain*(train / deno)
Pbus <- Data$ModeAvailableTrain *(bus / deno)
Pcar <- Data$ModeAvailableCar *(car / deno)
Pbike <- Data$ModeAvailableBike *(bike / deno)
Pwalk <- Data$ModeAvailableWalk *(walk / deno)
## Avoid problems stemming from choice probabilities becoming zero.
Ptrain <- (Ptrain!=0)*Ptrain + (Ptrain==0)
Pbus <- (Pbus!=0)*Pbus + (Pbus==0)
Pcar <- (Pcar !=0)*Pcar + (Pcar ==0)
Pbike <- (Pbike !=0)*Pbike + (Pbike ==0)
Pwalk <- (Pwalk!=0)*Pwalk + (Pwalk ==0)
## Choice results
Ctrain <- Data$MainModeENG =="Rail"
Cbus <- Data$MainModeENG =="Bus"
Ccar <- Data$MainModeENG =="Car"
Cbike <- Data$MainModeENG =="Bicycle"
Cwalk <- Data$MainModeENG =="Walk"
## Calculate the Log-likelihood function
LL <- colSums(Ctrain*log(Ptrain) + Cbus*log(Pbus) +
Ccar *log(Pcar) + Cbike *log(Pbike) +Cwalk *log(Pwalk))
}
##### Maximize the Log-likelihood function#####
##Parameter optimization
res<-optim(b0,fr,method="Nelder-Mead",hessian=TRUE,control=list(fnscale=-1))
## Parmeter estimation銆丠essian matrix calculation
b <- res$par
hhh <- res$hessian
运行Parameter optimization那里报错
res<-optim(b0,fr,method="Nelder-Mead",hessian=TRUE,control=list(fnscale=-1))
Error in optim(b0, fr, method = "Nelder-Mead", hessian = TRUE, control = list(fnscale = -1)) :
function cannot be evaluated at initial parameters
