R is a programming language for statistical computing and graphics. To install R, required installation files can be found in the following link:
Apollo Choice Modeling framework has been employed in running the vehicle ownership model component. We prepared the R notebooks for all models included in the guidebook.
rm(list = ls())
In first step of data processing in R, we need to call required package
library(apollo)
## Apollo 0.2.8
## http://www.ApolloChoiceModelling.com
## See url for a detailed manual, examples and a user forum.
## Sign up to the user forum to receive updates on new releases.
##
## Your version of Apollo is more than six months old.
## Using the latest version will ensure you have all
## current functionality and bug fixes.
## You can update to the latest version by typing:
## install.packages("apollo")
apollo_initialise()
apollo_control = list(
modelName = "MNL_Vehicle Ownership",
modelDescr = "Simple MNL model on vehicle ownership data (random model)",
indivID = "id",
outputDirectory = "output",
panelData = FALSE
)
database = read.csv("D:/Projects/20-102/Project Task/Task 4/R code/VO data_random.csv",header=TRUE)
summary(database)
## id Choice veh0 veh1
## Min. : 1 Min. :1.000 Min. :0.0000 Min. :0.0000
## 1st Qu.: 2570 1st Qu.:2.000 1st Qu.:0.0000 1st Qu.:0.0000
## Median : 5140 Median :3.000 Median :0.0000 Median :0.0000
## Mean : 5140 Mean :2.877 Mean :0.1064 Mean :0.2353
## 3rd Qu.: 7709 3rd Qu.:4.000 3rd Qu.:0.0000 3rd Qu.:0.0000
## Max. :10278 Max. :4.000 Max. :1.0000 Max. :1.0000
## veh2 veh3 wrktodrv age615todrv
## Min. :0.0000 Min. :0.000 Min. :0.000 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:0.000 1st Qu.:0.500 1st Qu.:0.0000
## Median :0.0000 Median :0.000 Median :1.000 Median :0.0000
## Mean :0.3333 Mean :0.325 Mean :0.687 Mean :0.1892
## 3rd Qu.:1.0000 3rd Qu.:1.000 3rd Qu.:1.000 3rd Qu.:0.0000
## Max. :1.0000 Max. :1.000 Max. :4.000 Max. :4.0000
## age1824todrv age80todrv inc_10 inc_1030
## Min. :0.00000 Min. :0.00000 Min. :0.00000 Min. :0.0000
## 1st Qu.:0.00000 1st Qu.:0.00000 1st Qu.:0.00000 1st Qu.:0.0000
## Median :0.00000 Median :0.00000 Median :0.00000 Median :0.0000
## Mean :0.04756 Mean :0.06482 Mean :0.04077 Mean :0.1355
## 3rd Qu.:0.00000 3rd Qu.:0.00000 3rd Qu.:0.00000 3rd Qu.:0.0000
## Max. :2.00000 Max. :4.00000 Max. :1.00000 Max. :1.0000
## inc_3060 inc_60100 uno
## Min. :0.0000 Min. :0.0000 Min. :1
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:1
## Median :0.0000 Median :0.0000 Median :1
## Mean :0.2109 Mean :0.2544 Mean :1
## 3rd Qu.:0.0000 3rd Qu.:1.0000 3rd Qu.:1
## Max. :1.0000 Max. :1.0000 Max. :1
table(database$Choice)
##
## 1 2 3 4
## 1094 2418 3426 3340
Here, Choice = Vehicle ownership level, veh0 = zero vehicle HH, veh1 = HH with one vehicle, veh2= HH with two vehicles, veh3 = HH with three or more vehicles, wrktodrv = Ratio of workers to drivers, age615todrv = Ratio of people age 6-15 to drivers, age1824todrv = Ratio of people age 18-24 to drivers, age80todrv = Ratio of people age 80 or older to drivers, inc_10 = Household income <=$10,000, inc_1030 = $10,000 < Household income <= $30,000, inc_60100 = $60,000 < Household income <= $100,000
Include fixed parameters in estimation process. Vector with names (in quotes) of parameters to be kept fixed at their starting value in apollo_beta, use apollo_beta_fixed = c() if none
apollo_beta=c(asc_veh3 = 0,
asc_veh2 = 0,
asc_veh1 = 0,
asc_veh0 = 0,
b_wrktodrv_veh0 = 0,
b_615todrv_veh2 = 0,
b_1824todrv_veh2 = 0,
b_1824todrv_veh1 = 0,
b_1824todrv_veh0 = 0,
b_80todrv_veh1 = 0,
b_inc10_veh2 = 0,
b_inc10_veh1 = 0,
b_inc10_veh0 = 0,
b_inc1030_veh2 = 0,
b_inc1030_veh1 = 0,
b_inc1030_veh0 = 0,
b_inc3060_veh2 = 0,
b_inc3060_veh1 = 0,
b_inc3060_veh0 = 0,
b_inc60100_veh1 = 0,
b_inc60100_veh0 = 0)
apollo_fixed = c("asc_veh3")
apollo_inputs = apollo_validateInputs()
## All checks on apollo_control completed.
## All checks on database completed.
# Define model and likelihood function
apollo_probabilities=function(apollo_beta, apollo_inputs, functionality="estimate"){
# Attach inputs and detach after function exit
apollo_attach(apollo_beta, apollo_inputs)
on.exit(apollo_detach(apollo_beta, apollo_inputs))
# Create list of probabilities P
P = list()
# List of utilities: these must use the same names as in mnl_settings, order is irrelevant
V = list()
V[["veh3"]] = asc_veh3
V[["veh2"]] = asc_veh2 + b_615todrv_veh2 * age615todrv + b_1824todrv_veh2 * age1824todrv + b_inc10_veh2 * inc_10 + b_inc1030_veh2 * inc_1030 + b_inc3060_veh2 * inc_3060
V[["veh1"]] = asc_veh1 + b_1824todrv_veh1 * age1824todrv + b_80todrv_veh1 * age80todrv + b_inc10_veh1 * inc_10 + b_inc1030_veh1 * inc_1030 + b_inc3060_veh1 * inc_3060 + b_inc60100_veh1 * inc_60100
V[["veh0"]] = asc_veh0 + b_wrktodrv_veh0 * wrktodrv + b_1824todrv_veh0 * age1824todrv + b_inc10_veh0 * inc_10 + b_inc1030_veh0 * inc_1030 + b_inc3060_veh0 * inc_3060 + b_inc60100_veh0 * inc_60100
# Define settings for MNL model component
mnl_settings = list(
alternatives = c(veh0=1, veh1=2, veh2=3, veh3=4),
avail = list(veh0=uno, veh1=uno, veh2=uno, veh3=uno),
choiceVar = Choice,
utilities = V
)
# Compute probabilities using MNL model
P[["model"]] = apollo_mnl(mnl_settings, functionality)
# Prepare and return outputs of function
P = apollo_prepareProb(P, apollo_inputs, functionality)
return(P)
nrow(p)
}
model = apollo_estimate(apollo_beta, apollo_fixed, apollo_probabilities, apollo_inputs)
## Preparing user-defined functions.
##
## Testing likelihood function...
##
## Overview of choices for MNL model component :
## veh0 veh1 veh2 veh3
## Times available 10278.00 10278.00 10278.00 10278.0
## Times chosen 1094.00 2418.00 3426.00 3340.0
## Percentage chosen overall 10.64 23.53 33.33 32.5
## Percentage chosen when available 10.64 23.53 33.33 32.5
##
## Pre-processing likelihood function...
##
## Testing influence of parameters
## Starting main estimation
## Initial function value: -14248.33
## Initial gradient value:
## asc_veh2 asc_veh1 asc_veh0 b_wrktodrv_veh0
## 856.500000 -151.500000 -1475.500000 -1197.832143
## b_615todrv_veh2 b_1824todrv_veh2 b_1824todrv_veh1 b_1824todrv_veh0
## 277.344643 -3.427381 -25.710714 -84.677381
## b_80todrv_veh1 b_inc10_veh2 b_inc10_veh1 b_inc10_veh0
## 37.787500 -24.750000 37.250000 66.250000
## b_inc1030_veh2 b_inc1030_veh1 b_inc1030_veh0 b_inc3060_veh2
## -5.250000 164.750000 -32.250000 143.000000
## b_inc3060_veh1 b_inc3060_veh0 b_inc60100_veh1 b_inc60100_veh0
## 135.000000 -314.000000 -175.750000 -465.750000
## initial value 14248.333444
## iter 2 value 13455.116794
## iter 3 value 13254.844582
## iter 4 value 13211.685064
## iter 5 value 13094.669148
## iter 6 value 12999.749602
## iter 7 value 12938.003327
## iter 8 value 12910.151983
## iter 9 value 12882.476386
## iter 10 value 12857.715221
## iter 11 value 12845.390471
## iter 12 value 12834.974511
## iter 13 value 12825.559600
## iter 14 value 12819.301252
## iter 15 value 12809.730226
## iter 16 value 12793.648454
## iter 17 value 12761.585839
## iter 18 value 12745.817958
## iter 19 value 12734.904398
## iter 20 value 12731.921296
## iter 21 value 12731.770058
## iter 22 value 12731.306044
## iter 23 value 12729.699959
## iter 24 value 12728.184005
## iter 25 value 12717.390090
## iter 26 value 12717.015585
## iter 27 value 12716.976055
## iter 28 value 12716.973412
## iter 28 value 12716.973262
## iter 28 value 12716.973260
## final value 12716.973260
## converged
## Additional convergence test using scaled estimation. Parameters will be
## scaled by their current estimates and additional iterations will be
## performed.
## initial value 12716.973260
## iter 1 value 12716.973259
## final value 12716.973259
## converged
## Estimated parameters:
## Estimate
## asc_veh3 0.00000
## asc_veh2 -0.06776
## asc_veh1 -0.86894
## asc_veh0 -1.68859
## b_wrktodrv_veh0 -0.39984
## b_615todrv_veh2 0.22736
## b_1824todrv_veh2 -1.56179
## b_1824todrv_veh1 -1.57418
## b_1824todrv_veh0 -2.00075
## b_80todrv_veh1 0.32270
## b_inc10_veh2 1.37416
## b_inc10_veh1 2.76353
## b_inc10_veh0 3.91279
## b_inc1030_veh2 0.62826
## b_inc1030_veh1 1.85043
## b_inc1030_veh0 2.39792
## b_inc3060_veh2 0.29759
## b_inc3060_veh1 1.09895
## b_inc3060_veh0 1.13072
## b_inc60100_veh1 0.16799
## b_inc60100_veh0 0.38093
##
## Final LL: -12716.9733
##
## Calculating log-likelihood at equal shares (LL(0)) for applicable
## models...
## Calculating log-likelihood at observed shares from estimation data
## (LL(c)) for applicable models...
## Calculating LL of each model component...
## Calculating other model fit measures
## Computing covariance matrix using analytical gradient.
## 0%....25%....50%....75%..100%
## Negative definite Hessian with maximum eigenvalue: -7.280963
## Computing score matrix...
apollo_modelOutput(model)
## Model run by Surface using Apollo 0.2.8 on R 4.1.3 for Windows.
## www.ApolloChoiceModelling.com
##
## Model name : MNL_Vehicle Ownership
## Model description : Simple MNL model on vehicle ownership data (random model)
## Model run at : 2023-11-20 12:38:25
## Estimation method : bfgs
## Model diagnosis : successful convergence
## Number of individuals : 10278
## Number of rows in database : 10278
## Number of modelled outcomes : 10278
##
## Number of cores used : 1
## Model without mixing
##
## LL(start) : -14248.33
## LL at equal shares, LL(0) : -14248.33
## LL at observed shares, LL(C) : -13467.87
## LL(final) : -12716.97
## Rho-squared vs equal shares : 0.1075
## Adj.Rho-squared vs equal shares : 0.1061
## Rho-squared vs observed shares : 0.0558
## Adj.Rho-squared vs observed shares : 0.0543
## AIC : 25473.95
## BIC : 25618.7
##
## Estimated parameters : 20
## Time taken (hh:mm:ss) : 00:00:5.9
## pre-estimation : 00:00:2.06
## estimation : 00:00:1.21
## post-estimation : 00:00:2.63
## Iterations : 32
## Min abs eigenvalue of Hessian : 7.280963
##
## Unconstrained optimisation.
##
## These outputs have had the scaling used in estimation applied to them.
## Estimates:
## Estimate s.e. t.rat.(0) Rob.s.e. Rob.t.rat.(0)
## asc_veh3 0.00000 NA NA NA NA
## asc_veh2 -0.06776 0.03160 -2.144 0.03128 -2.166
## asc_veh1 -0.86894 0.04790 -18.142 0.04782 -18.172
## asc_veh0 -1.68859 0.09668 -17.466 0.10164 -16.614
## b_wrktodrv_veh0 -0.39984 0.08156 -4.902 0.08861 -4.513
## b_615todrv_veh2 0.22736 0.04932 4.610 0.04853 4.685
## b_1824todrv_veh2 -1.56179 0.16388 -9.530 0.16763 -9.317
## b_1824todrv_veh1 -1.57418 0.17376 -9.060 0.19753 -7.969
## b_1824todrv_veh0 -2.00075 0.24952 -8.018 0.30780 -6.500
## b_80todrv_veh1 0.32270 0.08596 3.754 0.08551 3.774
## b_inc10_veh2 1.37416 0.23205 5.922 0.24028 5.719
## b_inc10_veh1 2.76353 0.22322 12.380 0.23175 11.925
## b_inc10_veh0 3.91279 0.23294 16.797 0.24105 16.233
## b_inc1030_veh2 0.62826 0.09280 6.770 0.09480 6.627
## b_inc1030_veh1 1.85043 0.09505 19.468 0.09692 19.092
## b_inc1030_veh0 2.39792 0.12136 19.759 0.12302 19.493
## b_inc3060_veh2 0.29759 0.06402 4.649 0.06431 4.627
## b_inc3060_veh1 1.09895 0.07395 14.860 0.07399 14.852
## b_inc3060_veh0 1.13072 0.10984 10.294 0.11046 10.236
## b_inc60100_veh1 0.16799 0.06800 2.471 0.06786 2.475
## b_inc60100_veh0 0.38093 0.10698 3.561 0.10706 3.558
Here, asc_veh3 = Alternative specific constant for HH with three or more vehicles, asc_veh2 = Alternative specific constant for HH with two vehicles, asc_veh1 = Alternative specific constant for HH with one vehicle, asc_veh0 = Alternative specific constant for HH with zero vehicle, b_wrktodrv_veh0 = Coefficient for ratio of workers to drivers (zero vehicle HH), b_615todrv_veh2 = Coefficient for ratio of people age 6-15 to drivers (HH with two or more vehicles), b_1824todrv_veh2 = Coefficient for ratio of people age 18-24 to drivers (HH with two or more vehicles), b_1824todrv_veh1 = Coefficient for ratio of people age 18-24 to drivers to drivers (HH with one vehicle), b_1824todrv_veh0 = Coefficient for ratio of people age 18-24 to drivers to drivers (HH with zero vehicle), b_80todrv_veh1 = Coefficient for ratio of people age 80 or older to drivers to drivers (HH with one vehicle), b_inc10_veh2 = Coefficient for household income <=$10,000 (HH with two vehicles), b_inc10_veh1 = Coefficient for household income <=$10,000 (HH with one vehicle), b_inc10_veh0 = Coefficient for household income <=$10,000 (zero vehicle HH), b_inc1030_veh2 = Coefficient for $10,000 < household income <= $30,000 (HH with two vehicles), b_inc1030_veh1 = Coefficient for $10,000 < household income <= $30,000 (HH with one vehicle), b_inc1030_veh0 = Coefficient for $10,000 < household income <= $30,000 (zero vehicle HH), b_inc3060_veh2 = Coefficient for $30,000 < household income <= $60,000 (HH with two vehicles), b_inc3060_veh1 = Coefficient for $30,000 < household income <= $60,000 (HH with one vehicle), b_inc3060_veh0 = Coefficient for $30,000 < household income <= $60,000 (zero vehicle HH), b_nmo_veh2 = Coefficient for NMO availability indicator (urban household indicator) (HH with two vehicles), b_inc60100_veh1 = Coefficient for $60,000 < household income <= $100,000 (HH with one vehicle), b_inc60100_veh0 = Coefficient for $60,000 < household income <= $100,000 (zero vehicle HH)