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 (base model)",
indivID = "id",
outputDirectory = "output",
panelData = FALSE
)database = read.csv("D:/Projects/20-102/Project Task/Task 4/R code/VO data_base.csv",header=TRUE)
summary(database)## id Choice veh0 veh1
## Min. : 1 Min. :1.000 Min. :0.00000 Min. :0.0000
## 1st Qu.: 2570 1st Qu.:2.000 1st Qu.:0.00000 1st Qu.:0.0000
## Median : 5140 Median :3.000 Median :0.00000 Median :0.0000
## Mean : 5140 Mean :2.942 Mean :0.04145 Mean :0.2525
## 3rd Qu.: 7709 3rd Qu.:4.000 3rd Qu.:0.00000 3rd Qu.:1.0000
## Max. :10278 Max. :4.000 Max. :1.00000 Max. :1.0000
## veh2 veh3 wrktodrv age615todrv
## Min. :0.0000 Min. :0.0000 Min. :0.000 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.500 1st Qu.:0.0000
## Median :0.0000 Median :0.0000 Median :1.000 Median :0.0000
## Mean :0.4284 Mean :0.2777 Mean :0.687 Mean :0.1892
## 3rd Qu.:1.0000 3rd Qu.:1.0000 3rd Qu.:1.000 3rd Qu.:0.0000
## Max. :1.0000 Max. :1.0000 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. :1table(database$Choice)##
## 1 2 3 4
## 426 2595 4403 2854Here, 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_veh1 = 0,
b_wrktodrv_veh0 = 0,
b_615todrv_veh2 = 0,
b_615todrv_veh1 = 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_wrktodrv_veh1 * wrktodrv + b_615todrv_veh1 * age615todrv + 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.00
## Times chosen 426.00 2595.00 4403.00 2854.00
## Percentage chosen overall 4.14 25.25 42.84 27.77
## Percentage chosen when available 4.14 25.25 42.84 27.77
##
## 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_veh1
## 1833.50000 25.50000 -2143.50000 -31.81548
## b_wrktodrv_veh0 b_615todrv_veh2 b_615todrv_veh1 b_1824todrv_veh2
## -1700.64881 569.09464 -124.37202 18.55595
## b_1824todrv_veh1 b_1824todrv_veh0 b_80todrv_veh1 b_inc10_veh2
## -29.96071 -101.21071 67.95417 -45.75000
## b_inc10_veh1 b_inc10_veh0 b_inc1030_veh2 b_inc1030_veh1
## 92.25000 42.25000 11.75000 330.75000
## b_inc1030_veh0 b_inc3060_veh2 b_inc3060_veh1 b_inc3060_veh0
## -150.25000 332.00000 261.00000 -502.00000
## b_inc60100_veh1 b_inc60100_veh0
## -177.75000 -643.75000
## initial value 14248.333444
## iter 2 value 12227.941889
## iter 3 value 12058.141183
## iter 4 value 11973.668443
## iter 5 value 11787.824685
## iter 6 value 11627.504735
## iter 7 value 11466.245417
## iter 8 value 11434.282848
## iter 9 value 11400.840325
## iter 10 value 11375.511280
## iter 11 value 11346.398967
## iter 12 value 11316.210466
## iter 13 value 11300.979497
## iter 14 value 11257.279339
## iter 15 value 11216.255733
## iter 16 value 11183.256269
## iter 17 value 11073.182585
## iter 18 value 10990.647270
## iter 19 value 10968.864338
## iter 20 value 10954.186065
## iter 21 value 10928.740821
## iter 22 value 10918.270067
## iter 23 value 10899.532190
## iter 24 value 10821.465758
## iter 25 value 10819.919001
## iter 26 value 10789.512055
## iter 27 value 10773.261402
## iter 28 value 10765.802585
## iter 29 value 10765.000379
## iter 30 value 10764.949772
## iter 31 value 10764.942332
## iter 32 value 10764.936238
## iter 33 value 10764.934516
## iter 34 value 10764.933895
## iter 35 value 10764.933702
## iter 35 value 10764.933642
## iter 35 value 10764.933640
## final value 10764.933640
## converged
## Additional convergence test using scaled estimation. Parameters will be
## scaled by their current estimates and additional iterations will be
## performed.
## initial value 10764.933640
## iter 1 value 10764.933639
## final value 10764.933639
## converged
## Estimated parameters:
## Estimate
## asc_veh3 0.0000
## asc_veh2 0.3378
## asc_veh1 -1.3170
## asc_veh0 -2.3433
## b_wrktodrv_veh1 0.3571
## b_wrktodrv_veh0 -2.1606
## b_615todrv_veh2 0.3817
## b_615todrv_veh1 -0.1599
## b_1824todrv_veh2 -2.0673
## b_1824todrv_veh1 -2.4203
## b_1824todrv_veh0 -1.7042
## b_80todrv_veh1 0.4839
## b_inc10_veh2 1.2573
## b_inc10_veh1 4.0763
## b_inc10_veh0 5.2541
## b_inc1030_veh2 0.6774
## b_inc1030_veh1 2.8862
## b_inc1030_veh0 3.3790
## b_inc3060_veh2 0.4068
## b_inc3060_veh1 1.7835
## b_inc3060_veh0 1.0693
## b_inc60100_veh1 0.5179
## b_inc60100_veh0 -0.7604
##
## Final LL: -10764.9336
##
## 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: -3.72853
## 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 (base model)
## Model run at : 2023-11-20 12:28:00
## 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) : -12317.2
## LL(final) : -10764.93
## Rho-squared vs equal shares : 0.2445
## Adj.Rho-squared vs equal shares : 0.2429
## Rho-squared vs observed shares : 0.126
## Adj.Rho-squared vs observed shares : 0.1242
## AIC : 21573.87
## BIC : 21733.1
##
## Estimated parameters : 22
## Time taken (hh:mm:ss) : 00:00:5.28
## pre-estimation : 00:00:1.74
## estimation : 00:00:1.15
## post-estimation : 00:00:2.4
## Iterations : 39
## Min abs eigenvalue of Hessian : 3.72853
##
## 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.0000 NA NA NA NA
## asc_veh2 0.3378 0.03128 10.799 0.03061 11.033
## asc_veh1 -1.3170 0.07571 -17.395 0.08059 -16.342
## asc_veh0 -2.3433 0.20022 -11.704 0.23622 -9.920
## b_wrktodrv_veh1 0.3571 0.06427 5.556 0.07028 5.080
## b_wrktodrv_veh0 -2.1606 0.17838 -12.112 0.23294 -9.275
## b_615todrv_veh2 0.3817 0.05882 6.490 0.05686 6.713
## b_615todrv_veh1 -0.1599 0.07630 -2.095 0.08405 -1.902
## b_1824todrv_veh2 -2.0673 0.16023 -12.902 0.17564 -11.770
## b_1824todrv_veh1 -2.4203 0.18646 -12.981 0.23465 -10.314
## b_1824todrv_veh0 -1.7042 0.32380 -5.263 0.41855 -4.072
## b_80todrv_veh1 0.4839 0.08787 5.507 0.09094 5.321
## b_inc10_veh2 1.2573 0.29383 4.279 0.31558 3.984
## b_inc10_veh1 4.0763 0.27892 14.615 0.30222 13.488
## b_inc10_veh0 5.2541 0.33837 15.528 0.36886 14.244
## b_inc1030_veh2 0.6774 0.10383 6.524 0.10831 6.254
## b_inc1030_veh1 2.8862 0.11019 26.193 0.11533 25.026
## b_inc1030_veh0 3.3790 0.22197 15.223 0.24499 13.793
## b_inc3060_veh2 0.4068 0.06558 6.204 0.06637 6.129
## b_inc3060_veh1 1.7835 0.08106 22.003 0.08172 21.823
## b_inc3060_veh0 1.0693 0.24871 4.299 0.26031 4.108
## b_inc60100_veh1 0.5179 0.07229 7.164 0.07193 7.200
## b_inc60100_veh0 -0.7604 0.36574 -2.079 0.36535 -2.081Here, 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_veh1 = Coefficient for ratio of workers to drivers (HH with one 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_615todrv_veh1 = Coefficient for ratio of people age 6-15 to drivers (HH with one vehicle), 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)