The median of the absolute value of the difference of two dependent log normal random variables

Question

Assume X and Y have a bivariate lognormal distribution (x,y>0) that is:

$f_{X,Y}(x,y)$=$$\frac{1}{2p\sqrt{1-r^2}xy\sigma_1\sigma_2}exp\{\frac{-1}{2(1-r^2)}[(\frac{ln(x)-\mu_1}{\sigma_1})^2-2r(\frac{ln(x)-\mu_1}{\sigma_1})(\frac{ln(y)-\mu_2}{\sigma_2})+(\frac{ln(y)-\mu_2}{\sigma_2})^2]\}$$
I want to know the median of |x-y|.
first I calculated the density of |x-y| based on changing the variables, as follows:

I took u=|x-y| and v=y then

\begin{cases} x=u+v, & \text{if x>y which result in: u>0 , v>0} \\ x=v−u, & \text{if x<y which result in v>u>0} \end{cases} The Jacobian is one, so the density function of u=|x-y| is:
$$f_{U}(u)=\int_{0}^\infty f_{X,Y}(u+v,v)dv+\int_{u}^\infty f_{X,Y}(v-u,v)dv$$

my question is: how can I calculate the median of the above density using numerical integration(assume for the simplest case where $\mu=\begin{pmatrix} 1 \\ 1 \\ \end{pmatrix} and \sigma= \begin{pmatrix} 1 & 0.95 \\ 0.95 & 1 \\ \end{pmatrix}$)?
I am quite new in R..I really appreciate any help

Answer 1

A simple Monte Carlo approach could work. This method approximates the distribution using a finite set of samples, then computes the median based on the samples. Draw $n$ independent samples of $X$ and $Y$ from the joint distribution: $\{(x_1, y_1), ... (x_n, y_n)\}$. I'll explain how do do this below. This gives samples of $U$: $\{u_1, ..., u_n\}$, where $u_i = |x_i - y_i|$. Take the median of these samples. Accuracy will increase with $n$, which you could easily set in the thousands or millions.

Your problem has special structure that makes it possible to sample $X$ and $Y$ directly from the joint distribution. $\log X$ and $\log Y$ have a bivariate normal distribution from which you can easily sample using a multivariate Gaussian random number generator (using the appropriate mean and covariance matrix). This blog post describes how to do this in R. Once you have samples of $\log X$ and $\log Y$, exponentiate these values to obtain samples of $X$ and $Y$. It's not always possible to sample directly from a given distribution, in which case methods like MCMC algorithms can be used. Luckily, this isn't necessary in your case.