I'm simulating the data underlying a conditional logit model then doing my own MLE estimation using optim
. However, even with an unreasonably large amount of data, the estimator seems to have awful precision. Is this expected, or I am making mistakes in my implementation?
Conditional Logit Model: The conditional logit is a discrete choice model, where a person chooses between multiple alternatives based on the alternatives' covariates. We want to estimate the person's preference for these covariates.
The person chooses the alternative that gives him the highest utility. For example, with 3 choices, the 3 utility functions are:
$$ U_{i1} = \alpha' X_1 + \varepsilon_1 \\ U_{i2} = \alpha' X_2 + \varepsilon_2 \\ U_{i3} = \alpha' X_3 + \varepsilon_3 $$
where $\alpha$ is the preference parameters to be estimated, and $X$ is the alternative's covariates that we observe.
In the conditional logit model, we assume that the random utility components (i.e. the $\varepsilon$) follows a standard Gumbel distribution. With that assumption, the probability of choosing alternative $j$ is:
$$ \Pr(\text{choosing }j) = \Pr(j\text{ offers highest utility}) = \frac{\exp(\alpha'X_j)}{\sum_{\text{all }j} \exp(\alpha'X_j)} $$
The log likelihood is then
$$ LL = \sum_i \left( \alpha'X_{c_i} - \log\left( \sum_j \exp(\alpha' X_j)\right) \right) $$
where $c_i \in \{1, \dots, j\}$ indidates the choice that $i$ makes.
Implementation:
I simulate the data using known preference parameters and estimate them using MLE as below.
num_choices <- 100
xx <- mvtnorm::rmvnorm(num_choices, sigma = diag(2)) # Alternatives' covariates
alpha <- c(0.05, 0.1) # Preference parameters, to be estimated
# Matrix of utilities
mat <- matrix(NA, nrow = 1000, ncol = num_choices)
for (i in 1:num_choices) {
mat[, i] <- sum(alpha * xx[i, ]) + evd::rgumbel(1000)
}
# Choosing the alternative that gives the highest utility
y <- max.col(mat)
# negative log likelihood
cl_nllik <- function(alpha) {
xa <- c(xx %*% alpha)
lse_xa <- log(sum(exp(xa)))
- sum( xa[y] - lse_xa )
}
# MLE estimate -- does NOT produce the true alpha values!
optim(c(0, 0), cl_nllik)