Let $X_1$ and $X_2$ be 2 i.i.d. r.v.'s where $\log(X_1),\log(X_2) \sim N(\mu,\sigma)$. I'd like to know the distribution for $X_1 - X_2$.

The best I can do is to take the Taylor series of both and get that the difference is the sum of the difference between two normal r.v's and two chi-squared r.v.'s in addition to the rest of the difference between the rest of the terms. Is there a more straight-forward way to get the distribution of the difference between 2 i.i.d. log-normal r.v.'s?

  • $\begingroup$ Here is a relevant paper. You will find more papers by googling! papers.ssrn.com/sol3/papers.cfm?abstract_id=2064829 $\endgroup$ Commented May 18, 2015 at 15:04
  • 1
    $\begingroup$ I've taken a cursory glance at that paper, and it doesn't seem to answer my question in a satisfying way. They seem to be concerned with numerical approximations to the harder problem of finding the distribution for the sum/difference between correlated lognormal r.v.'s. I was hoping that there would be a simpler answer for the independent case. $\endgroup$
    – frayedchef
    Commented May 18, 2015 at 15:24
  • 2
    $\begingroup$ It might be a simpler answer in the independent case, but not a simple one! The lognormal case is a famously known hard case---the moment-generating function of the lognormal distribution doesnt exist---that is, it doesnt converge on an open interval containing zero. So, you will not find an easy solution. $\endgroup$ Commented May 18, 2015 at 16:42
  • $\begingroup$ I see... So would the approach I outlined above be reasonable? (i.e., if $Y_i = \log(X_i)$, $X_1 - X_2 \approx (Y_1 - Y_2) + (Y_1^2 - Y_2^2)/2 + {} ...$ Do we know anything about the higher order terms, or how to bound them? $\endgroup$
    – frayedchef
    Commented May 18, 2015 at 17:10
  • 1
    $\begingroup$ To illustrate the difficulty---the lognormal mgf is only defined on $(-\infty,0]$. To approximate the difference distribution by saddlepoint methods, we need (K=cumulant gf) $K(s)+K(-s)$, and that sum is only defined in one point, zero. So, doesnt seem to work. Sum or average would be simpler! $\endgroup$ Commented Oct 7, 2015 at 17:20

3 Answers 3


This is a difficult problem. I thought first about using (some approximation of) the moment generating function of the lognormal distribution. That doesn't work, as I will explain. But first some notation:

Let $\phi$ be the standard normal density and $\Phi$ the corresponding cumulative distribution function. We will only analyze the case lognormal distribution $lnN(0,1)$, which has density function $$ f(x)=\frac1{\sqrt{2\pi}x} e^{-\frac12 (\ln x)^2} $$ and cumulative distribution function $$ F(x) =\Phi(\ln x) $$ Suppose $X$ and $Y$ are independent random variables with the above lognormal distribution. We are interested in the distribution of $D=X-Y$, which is a symmetric distribution with mean zero. Let $M(t) = \DeclareMathOperator{\E}{E} \E e^{tX} $ be the moment generating function of $X$. It is defined only for $t\in (-\infty,0]$, so not defined in an open interval containing zero. The moment generating function for $D$ is $M_D(t)=\E e^{t(X-Y)}= \E e^{tX} \E e^{-tY}= M(t)M(-t)$. So, the moment generating function for $D$ is only defined for $t=0$, so not very useful.

That means we will need some more direct approach for finding approximations for the distribution of $D$. Suppose $t\ge 0$, calculate $$ \begin{align} P(D \le t) &= P(X-Y\le t) \\ &= \int_0^\infty P(X-y\le t | Y=y) f(y) \; dy \\ &= \int_0^\infty P(X\le t+y) f(y) \; dy \\ &= \int_0^\infty F(t+y) f(y) \; dy \end{align} $$ (and the case $t<0$ is solved by symmetry, we get $P(D\le t)=1-P(D\le |t|)$).

This expression can be used for numerical integration or as a basis for simulation. First a test:

 integrate(function(y) plnorm(y)*dlnorm(y), lower=0,  upper=+Inf)
  0.5 with absolute error < 2.3e-06

which is clearly correct. Let us wrap this up inside a function:

pDIFF  <-  function(t) {
    d  <-  t
    for (tt in seq(along=t)) {
        if (t[tt] >= 0.0) d[tt] <- integrate(function(y) plnorm(y+t[tt])*dlnorm(y),
                                         lower=0.0,  upper=+Inf)$value else
                          d[tt] <- 1-integrate(function(y) plnorm(y+abs(t[tt]))*dlnorm(y),
                                         lower=0.0, upper=+Inf)$value

> plot(pDIFF,  from=-5,  to=5)

which gives:

cumulative distribution function found by numerical integration

Then we can find the density function by differentiating under the integral sign, obtaining

dDIFF  <-  function(t) {
       d  <- t; t<- abs(t)
       for (tt in seq(along=t)) {
           d[tt]  <-  integrate(function(y) dlnorm(y+t[tt])*dlnorm(y),
                                lower=0.0,  upper=+Inf)$value

which we can test:

> integrate(dDIFF,  lower=-Inf,  upper=+Inf)
0.9999999 with absolute error < 1.3e-05

And plotting the density we get:

plot(dDIFF,  from=-5,  to=5)

density function found by numerical integration

I did also try to get some analytic approximation, but so far didn't succeed, it is not an easy problem. But numerical integration as above, programmed in R is very fast on modern hardware, so is a good alternative which probably should be used much more.

  • 2
    $\begingroup$ Straight numerical convolution may be feasible for some applications $\endgroup$
    – Glen_b
    Commented Sep 18, 2020 at 3:52

This does not strictly answer your question, but wouldn't it be easier to look at the ratio of the $X$ and $Y$? You then simply arrive at

$$ \begin{align} \Pr\left(\frac{X}{Y} \leq t\right) &= \Pr\left(\log\left(\frac{X}{Y}\right) \leq \log(t) \right) \\ &= \Pr(\log(X) - \log(Y) \leq \log(t)) \\ &\sim \mathcal{N}(0, 2 \sigma^2) \end{align}$$

Depending on your application, this may serve your needs.

  • 4
    $\begingroup$ But aren't we looking at X-Y instead of log(X) - log(Y) ? $\endgroup$ Commented Jan 31, 2018 at 12:02
  • 2
    $\begingroup$ Yes, of course. This is just in case somebody would be interested in knowing how two lognormal variables differ from each other without it necessarily needing to be a difference. That's why I also say it doesn't the answer the question. $\endgroup$ Commented Feb 3, 2018 at 21:41

I have a feeling characteristic functions are the best option for handling this question.

For basics, kindly see:https://en.wikipedia.org/wiki/Characteristic_function_(probability_theory)

now, kindly CTRL+F "independen". you'll see the answer

in brief:

  • The characteristic function approach is particularly useful in analysis of linear combinations of independent random variables

  • If X1, ..., Xn are independent random variables, and a1, ..., an are some constants, then the characteristic function of the linear combination of the Xi 's is
      φ a 1 X 1 + ⋯ + a n X n ( t ) = φ X 1 ( a 1 t ) ⋯ φ X n ( a n t ) . {\displaystyle \varphi _{a_{1}X_{1}+\cdots +a_{n}X_{n}}(t)=\varphi _{X_{1}}(a_{1}t)\cdots \varphi _{X_{n}}(a_{n}t).} 
  • One specific case is the sum of two independent random variables X1 and X2 in which case one has

    φ X 1 + X 2 ( t ) = φ X 1 ( t ) ⋅ φ X 2 ( t ) . {\displaystyle \varphi {X{1}+X_{2}}(t)=\varphi {X{1}}(t)\cdot \varphi {X{2}}(t).}

  • Characteristic functions are particularly useful for dealing with linear functions of independent random variables. For example, if X1, X2, ..., Xn is a sequence of independent (and not necessarily identically distributed) random variables, and

    S n = ∑ i = 1 n a i X i , {\displaystyle S_{n}=\sum {i=1}^{n}a{i}X_{i},,!} S_{n}=\sum {i=1}^{n}a{i}X_{i},,!

where the ai are constants, then the characteristic function for Sn is given by

φ S n ( t ) = φ X 1 ( a 1 t ) φ X 2 ( a 2 t ) ⋯ φ X n ( a n t ) {\displaystyle \varphi _{S_{n}}(t)=\varphi _{X_{1}}(a_{1}t)\varphi _{X_{2}}(a_{2}t)\cdots \varphi _{X_{n}}(a_{n}t)\,\!} \varphi _{S_{n}}(t)=\varphi _{X_{1}}(a_{1}t)\varphi _{X_{2}}(a_{2}t)\cdots \varphi _{X_{n}}(a_{n}t)\,\!

In particular, φX+Y(t) = φX(t)φY(t). To see this, write out the definition of characteristic function:

φ X + Y ( t ) = E ⁡ [ e i t ( X + Y ) ] = E ⁡ [ e i t X e i t Y ] = E ⁡ [ e i t X ] E ⁡ [ e i t Y ] = φ X ( t ) φ Y ( t ) {\displaystyle \varphi _{X+Y}(t)=\operatorname {E} \left[e^{it(X+Y)}\right]=\operatorname {E} \left[e^{itX}e^{itY}\right]=\operatorname {E} \left[e^{itX}\right]\operatorname {E} \left[e^{itY}\right]=\varphi _{X}(t)\varphi _{Y}(t)} {\displaystyle \varphi _{X+Y}(t)=\operatorname {E} \left[e^{it(X+Y)}\right]=\operatorname {E} \left[e^{itX}e^{itY}\right]=\operatorname {E} \left[e^{itX}\right]\operatorname {E} \left[e^{itY}\right]=\varphi _{X}(t)\varphi _{Y}(t)}

The independence of X and Y is required to establish the equality of the third and fourth expressions.

  • $\begingroup$ note that lognormal MGF is not always defined;but its CF is. en.wikipedia.org/wiki/Log-normal_distribution $\endgroup$ Commented Sep 18, 2020 at 0:57
  • 1
    $\begingroup$ Have you noticed that the cf for the Lognormal distribution is not analytically tractable? How are you supposing people should work with it -- numerically, perhaps? $\endgroup$
    – whuber
    Commented Sep 18, 2020 at 12:04
  • $\begingroup$ I fail to see the importance of analytical tractability in the context of the OP. $\endgroup$ Commented Sep 18, 2020 at 20:06
  • $\begingroup$ Usually people don't conceive of having a product of unevaluated infinite sums (which is what pertains to the lognormal distribution) as constituting "knowing" a distribution. Suitable mathematical objects are those that, upon analysis, readily yield basic information, of which the most fundamental is the ability to evaluate probabilities (obviously!) and densities (for MLE, for instance). Thus, merely offering a generic method to obtain a characteristic function is rarely seen as a constructive result. $\endgroup$
    – whuber
    Commented Sep 18, 2020 at 20:41
  • $\begingroup$ Indeed, I should have remembered that not everything is for everyone; rarity. $\endgroup$ Commented Sep 18, 2020 at 21:34

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.