Expected value of latent utility in logistic regression [duplicate]

Question

I'm looking for an analytical expression for the expected value of the latent utility in a logistic regression.

Setup:

There are two choices indexed by $i \in \{0, 1\}$ with associated utilities $u_i = w_i + \epsilon_i$, with the $\epsilon_i$ i.i.d. Gumbel(0, 1), aka Type I Extreme Value with location 0 and scale 1. The $w_i$ are constants.

(Same setup as http://en.wikipedia.org/wiki/Logistic_regression#As_a_two-way_latent-variable_model)

Let $Z \equiv \max\{u_0, u_1\}$, i.e. $Z$ is the utility you get when choosing the preferred alternative.

Question:

Is there a closed-form expression for the expected value of $Z$? If so, what is it, and how is it derived?

Additional Info:

If you know the answer, skip this section :-)

The $\epsilon_i$ have CDF $F(x) = e^{-e^{-x}}$ and expected value $\gamma \approx 0.5772$, the Euler–Mascheroni constant -- here I'm copying from http://en.wikipedia.org/wiki/Gumbel_distribution, in case this is helpful to anyone reading this question.

I did some googling, and I found something saying $\mathbb{E}[Z] = \ln(\Sigma_i(e^{w_i})) +$ constant. Is this correct?

I wrote some R code to investigate, and it seems like $\mathbb{E}[Z] = \ln(\Sigma_i(e^{w_i})) + \gamma \approx \ln(\Sigma_i(e^{w_i})) + 0.5772$. Is this correct?

Here's the code:

# Binary choice
# w_i + \epsilon_i with \epsilon_i ~ Gumbel
# Is Emax{w_i + \epsilon_i} = log(sum(exp(w_i))) + \gamma?

rgumbel <- function(n) {
  # Slower than Runuran, but avoids installing a new package
  u <- runif(n)
  return(-log(-log(u)))
}
library(Runuran)  # Gumbel
rgumbel <- function(n) {
  distribution <- udgumbel()
  generator <- pinvd.new(distribution)
  return(ur(generator, n))
}

gamma <- 0.57722 # Euler–Mascheroni constant
GetError <- function(w=c(5, 10), n=100000) {
  epsilon <- matrix(nrow=n, ncol=2, data=rgumbel(2*n))
  utility <- pmax(epsilon[, 1] + w[1],
                  epsilon[, 2] + w[2])
  emax <- mean(utility)
  emax.analytical <- log(sum(exp(w))) + gamma  # Correct?
  return(emax - emax.analytical)
}

The GetError function returns the difference between a sample mean and my candidate analytical expression. I played around and tended to see values very close to zero.

Here's my very feeble attempt at getting an analytical expression for $\mathbb{E}[Z]$:

$G(z) \equiv Pr[Z \leq z] = Pr[w_0 + \epsilon_0 \leq z] Pr[w_1 + \epsilon_1 \leq z | w_0 + \epsilon_0 \leq z]$

$G(z) = F(z - w_0) F(z - w_1) = e ^ {-e^{-z + w_0} -e^{-z + w_1}}$

$g(z) \equiv G^\prime(z) = e ^ {-e^{-z + w_0} -e^{-z + w_1}} (e^{-z + w_0} + e^{-z + w_1})$

$\mathbb{E}[Z] = \int z g(z) dz$

So, assuming there are no mistakes up to here, I need to calculate:

$\int z e ^ {-e^{-z + w_0} -e^{-z + w_1}} e^{-z + w_0} dz$ plus a similar second term, but I'm not sure how to deal with this integral.

Many thanks for any help,

Adrian

Answer 1

It's been a while, but since I've run into this recently, I'll post my answer. Letting $i = 1,...,N$, and assuming $\epsilon_i \sim \text{Gumbel}(\mu, 1), \forall i$, the closed form solution for $E[Z]$ is

$$E[Z] = E[\max_i(w_i + \epsilon_i)] = \mu + \gamma + \ln\left(\sum_i^N \exp\left\{w_i \right\} \right), $$

where $\gamma \approx 0.52277$ is the Euler-Mascheroni constant, and $\mu$ is the location parameter of the Gumbel distribution for each $\epsilon_i$. You can find one proof of this here, along with more details.

Answer 2

A partial answer. First, $$ F_{U_i}(u_i) = P(U_i\leq u_i) = P(U_i-w_i\leq u_i-w_i) = P(\epsilon_i\leq u_i-w_i) = e^{-e^{-(u_i - w_i)}} \, . $$ Also, $$ F_Z(z) = P(Z\leq z) = P(\max \{U_0,U_1\} \leq z) = P(U_0\leq z, U_1\leq z) $$ $$ = P(U_0\leq z)P(U_1\leq z) = F_{U_0}(z)F_{U_1}(z) = e^{-\left(e^{-(z - w_0)}+e^{-(z - w_1)}\right)} \, . $$ Therefore, $$ \mathrm{E}[Z] = \int_0^\infty (1 - F_Z(z))dz - \int_{-\infty}^0 F_Z(z)\,dz $$ $$ = \int_0^\infty \left( 1 - e^{-\left(e^{-(z - w_0)}+e^{-(z - w_1)}\right)} \right) dz - \int_{-\infty}^0 e^{-\left(e^{-(z - w_0)}+e^{-(z - w_1)}\right)} dz \, . $$ These integrals can probably be expressed in terms of the Exponential Integral.