I have two random variables which are independent and identically distributed, i.e. $\epsilon_{1}, \epsilon_{0} \overset{\text{iid}}{\sim} \text{Gumbel}(\mu,\beta)$:

$$F(\epsilon) = \exp(-\exp(-\frac{\epsilon-\mu}{\beta})),$$

$$f(\epsilon) = \dfrac{1}{\beta}\exp(-\left(\frac{\epsilon-\mu}{\beta}+\exp(-\frac{\epsilon-\mu}{\beta})\right)).$$

I am trying to calculate two quantities:

  1. $$\mathbb{E}_{\epsilon_{1}}\mathbb{E}_{\epsilon_{0}|\epsilon_{1}}\left[c+\epsilon_{1}|c+\epsilon_{1}>\epsilon_{0}\right]$$
  2. $$\mathbb{E}_{\epsilon_{1}}\mathbb{E}_{\epsilon_{0}|\epsilon_{1}}\left[\epsilon_{0}|c+\epsilon_{1}<\epsilon_{0}\right]$$

I get to a point where I need to do integration on something of the form: $e^{e^{x}}$, which seems to not have an integral in closed form. Can anyone help me out with this? Maybe I have done something wrong.

I feel there definitely should be closed form solution. (EDIT: Even if it is not closed form, but there would be software to quickly evaluate the integral [such as Ei(x)], that would be ok I suppose.)


I think with a change of variables, let

$$y =\exp(-\frac{\epsilon_{1}-\mu}{\beta})$$ and

$$\mu-\beta\ln y =\epsilon_{1}$$

This maps to $[0,\;\infty)$ and $\left[0,\;\exp(-\frac{\epsilon_{0}-c-\mu}{\beta})\right] $ respectively.

$|J|=|\dfrac{d\epsilon}{dy}|=\frac{\beta}{y}$. Then under the change of variable, I have boiled (1) down to...

$$\int_{0}^{\infty}\dfrac{1}{1-e^{-x}}\left(\int_{\mu-\beta\ln x-c}^{\infty}\left[c+\mu-\beta\ln y\right]e^{-y}dy\right)e^{-x}dx$$

There might be an algebra mistake but I still cannot solve this integral...

RELATED QUESTION: Expectation of the Maximum of iid Gumbel Variables

  • 1
    $\begingroup$ There's definitely no closed form solution. Why did you feel there must be? $\endgroup$ Commented Feb 20, 2017 at 4:21
  • $\begingroup$ @GordonSmyth How do you know there's no closed form solution? $\endgroup$ Commented May 1, 2019 at 20:27
  • $\begingroup$ How should I interpret the double expectation operators? $\mathbb{E}_{\epsilon_{1}}\mathbb{E}_{\epsilon_{0}|\epsilon_{1}}$ $\endgroup$ Commented Oct 15, 2021 at 8:07

2 Answers 2


Since the parameters $(\mu,\beta)$ of the Gumbel distribution are location and scale, respectively, the problem simplifies into computing $$\mathbb{E}[\epsilon_1|\epsilon_1+c>\epsilon_0]= \frac{\int_{-\infty}^{+\infty} x F(x+c) f(x) \text{d}x}{\int_{-\infty}^{+\infty} F(x+c) f(x) \text{d}x}$$ where $f$ and $F$ are associated with $\mu=0$, $\beta=1$. The denominator is available in closed form \begin{align*} \int_{-\infty}^{+\infty} F(x+c) f(x) \text{d}x &= \int_{-\infty}^{+\infty} \exp\{-\exp[-x-c]\}\exp\{-x\}\exp\{-\exp[-x]\}\text{d}x\\&\stackrel{a=e^{-c}}{=}\int_{-\infty}^{+\infty} \exp\{-(1+a)\exp[-x]\}\exp\{-x\}\text{d}x\\&=\frac{1}{1+a}\left[ \exp\{-(1+a)e^{-x}\}\right]_{-\infty}^{+\infty}\\ &=\frac{1}{1+a} \end{align*} The numerator involves an exponential integral since (according to WolframAlpha integrator) \begin{align*} \int_{-\infty}^{+\infty} x F(x+c) f(x) \text{d}x &= \int_{-\infty}^{+\infty} x \exp\{-(1+a)\exp[-x]\}\exp\{-x\}\text{d}x\\ &\stackrel{z=e^{-x}}{=} \int_{0}^{+\infty} \log(z) \exp\{-(1+a)z\}\text{d}z\\ &= \frac{-1}{1+a}\left[\text{Ei}(-(1+a) z) -\log(z) e^{-(1+a) z}\right]_{0}^{\infty}\\ &= \frac{\gamma+\log(1+a)}{1+a} \end{align*} Hence $$\mathbb{E}[\epsilon_1|\epsilon_1+c>\epsilon_0]=\gamma+\log(1+e^{-c})$$ This result can easily be checked by simulation, since producing a Gumbel variate amounts to transforming a Uniform (0,1) variate, $U$, as $X=-\log\{-\log(U)\}$. Monte Carlo and theoretical means do agree:

This figure demonstrates the adequation of Monte Carlo and theoretical means when $c$ varies from -2 to 2, with logarithmic axes, based on 10⁵ simulations

  • 1
    $\begingroup$ Did you realize epsilon0 is a rv as well? $\endgroup$ Commented Mar 1, 2019 at 2:04
  • 2
    $\begingroup$ @wolfsatthedoor The computation is not $$\mathbb{E}[\epsilon_1]= \frac{\int_{-\infty}^{+\infty} x f(x) \text{d}x}{\int_{-\infty}^{+\infty} f(x) \text{d}x}$$ but instead $$\mathbb{E}[\epsilon_1|\epsilon_1+c>\epsilon_0]= \frac{\int_{-\infty}^{+\infty} x F(x+c) f(x) \text{d}x}{\int_{-\infty}^{+\infty} F(x+c) f(x) \text{d}x}$$ this term $$ \frac{ F(x+c) f(x)}{\int_{-\infty}^{+\infty} F(x+c) f(x) \text{d}x}$$ is the probability density distribution of $\epsilon_1|\epsilon_1+c>\epsilon_0$. $\endgroup$ Commented Oct 15, 2021 at 8:12
  • $\begingroup$ This answer is very helpful, but am I right when I say a=e^c should actually be a=e^-c? Seems like that's what you substituted in the last line $\endgroup$ Commented May 31, 2022 at 17:39

Xi'an computed the answer more directly by evaluating the integrals. We could also get to the answer by arguing that the conditional distribution, when scaled appropriately, is a type 1 Gumbel distribution.

The distribution of $\epsilon_1$ conditional on $\epsilon_1 +c > \epsilon_0$, let's call it $\epsilon_2$, is proportional to the product of a pdf and cdf

$$f_{\epsilon_2}(\epsilon_2) \propto f_{\epsilon_1}(\epsilon_2) \cdot F_{\epsilon_0}(\epsilon_2+c)$$

This will be (with $z =\frac{\epsilon_2-\mu}{\beta}$ and $d = c/\beta$)

$$f_{\epsilon_2}(\epsilon_2) \propto e^{-(z+e^{-z})} \cdot e^{-e^{-z-d}} = e^{-(z+(1+e^{-d})e^{-z})} $$

This is like a type 1 Gumbel distribution or like the more general distribution function

$$f(x) = \frac{b^a}{\Gamma(a)} e^{-(ax+be^{-x})}$$

the mean is for $a=1$ equal to $\mu_{x} = \log (b) +\gamma$ with $\gamma$ the Euler-Mascheroni constant (see below for more details). Then the mean of $z$ is

$$\mu_{z} = \log (1+e^{-c/\beta}) + \gamma$$

and the mean of $\epsilon_2$ is

$$\mu_{\epsilon_2} = \log(1+e^{-c/\beta})\beta + \gamma\beta + \mu$$

Generalized Gumbel distribution and it's mean

This form $e^{-(ax+be^{-x})}$ occurs in Ahuja, J. C., and Stanley W. Nash. "The generalized Gompertz-Verhulst family of distributions." Sankhyā: The Indian Journal of Statistics, Series A (1967): 141-156.

Although the more well known 'Gumbel distribution' is with $a=b=1$, Gumbel studied this type of distribution in 1935 (Les valeurs extrêmes des distributions statistiques) for the distribution of the m-th order observation in a sample, and he used $a=b=m$ with $m$ a positive integer.

The moment generating function is given by Ahuja and Nash:

$$M(t;a,b) = \frac{1}{\Gamma(a)}b^t\Gamma(a-t)$$

and the mean has a closed form expression when $a$ is a positive integer

$$\mu_{X,a,b} = M^\prime(0;a,b) = \log (b) + \left( \gamma - \sum_{k=1}^{a-1} \frac{1}{k} \right)$$

in the case of $a=1$ this is

$$\mu_{X,1,b} = \log (b) + \gamma $$


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.