6
$\begingroup$

Assume $X$ and $Y$ are iid $N(0,1)$. I am looking for a "neat" expression for $$ P\left(\frac{X+Y}{\sqrt{2}}>c\,\Biggl|\,X<c\right). $$ Related question seem to be discussed here or here, but if they already give the answer to my question, I do not see it.

Simulation suggests it is around 3% for $c$ the 95% normal quantile:

c <- qnorm(0.95)
cprob.num <- rep(NA,50000)

for (i in 1:reps){
  X <- rnorm(1)
  Y <- rnorm(1)
  cprob.num[i] <- (X+Y)/sqrt(2) > c & X<c
}

mean(cprob.num)/0.95 # 0.03117895
$\endgroup$
  • $\begingroup$ With a simple change of variables, you can reduce this to a comparable question for a bivariate Normal distribution $(U,V)$ where $U=X-c$, $V=(X+Y)/\sqrt{2}-c$, $c=0$, and $(U,V)$ are correlated. Then by following the analysis at Using the methods described at stats.stackexchange.com/a/71303/919 or using standard formulas for bivariate normal distributions, the question is reduced to measuring an angle when $c=0$ and otherwise requires numerical integration. That might explain why you have been unable to find a formula. $\endgroup$ – whuber Jun 23 '17 at 13:43
  • $\begingroup$ Thank you. It looks like it will take me a while to digest the linked answer. $\endgroup$ – Christoph Hanck Jun 23 '17 at 14:00
  • $\begingroup$ X and Y have bivariate normal distribution: $X,Y\sim (\begin{bmatrix}0\\0\end{bmatrix},\begin{bmatrix}1&0\\0&1\end{bmatrix})$ However due to constraint $X<c$ we are dealing with truncated multivariate normal distribution. We need to impose linear transformation $A=[\sqrt{2}, \sqrt{2}]$ on it in order to find the result (using R tmvnorm package). Howevere I did not find even characteristic function formula for truncated normal in order to determine linear transformation consequences so it should be not straightforward. However its application could solve the problem. $\endgroup$ – Bogdan Jun 23 '17 at 19:52
  • $\begingroup$ @whuber, I have so far been unsuccessful to turn your help into a solution. In that vein, does Mathematica - see wolfies answer below - overlook some way to produce a closed-form result? $\endgroup$ – Christoph Hanck Jun 24 '17 at 7:57
4
$\begingroup$

Given: $X$ and $Y$ are independent standard Normals with pdf's $\phi(.)$ and cdf's $\Phi(.)$.

Since $X$ and $Y$ are independent, the joint pdf of $\big((X \; \big|\;X<c), \; Y\big)$ is $f(x,y) = {\large\frac{\phi(x)}{\Phi(c)}} \phi(y)$:

enter image description here

where Erf[.] denotes the error function.

Part 1: The pdf of $Z = X+Y \; | \; X<c$

Given $f(x,y)$, consider the transformation $(Z = X+Y, V=Y)$.

If $X <c$ and $Z = X+Y$, then $Z < c + Y$. That is, $Z < c + V$. This dependency is invoked in the following line using the Boole statement. Then the joint pdf of $(Z,V)$, say $g(z,v)$ can be obtained with:

enter image description here

... where I am using the Transform function from mathStatica/Mathematica to automate the nitty-gritties using the Method of Transformations (Jacobian etc).

The pdf of $Z$ that we seek is simply the marginal pdf of $Z$:

enter image description here

... which is our desired closed form solution.

The following diagram plots the pdf of $Z$ (i.e. the sum of 2 independent Normals, conditional on one of them) for six different vales of parameter $c$:

enter image description here


Part 2: Find $P\left(\frac{X+Y}{\sqrt{2}}>c\,\Biggl|\,X<c\right)$

To find $P\left(\frac{X+Y}{\sqrt{2}}>c\,\Biggl|\,X<c\right)$, integrate the above pdfZ over $(\sqrt2 c, \infty)$ wrt $z$.

Alternatively, $P\left(\frac{X+Y}{\sqrt{2}}>c\,\Biggl|\,X<c\right)$ can be obtained directly from the first step by :

enter image description here

... where I am using the Prob function from mathStatica/Mathematica to automate the nitty-gritties. This solution can be written in conventional notation as:

$$\frac{1}{\Phi(c)} \quad \int_{-\infty}^c \phi(x) \; \Phi \left(x-\sqrt{2} c\right) \, dx$$

While the probability does not appear to have a convenient closed-form, it is nevertheless a useful and practical result that is reduced to integrating a single variable. In particular:

a) when $c = 0$, the solution simplifies to $\frac14$

b) for other $c$ values, replace Integrate with NIntegrate for a solution via numerical integration in a single variable, which works very nicely. For instance, here is a plot of the desired probability, as a function of the truncation point $c$:

enter image description here

| cite | improve this answer | |
$\endgroup$
  • 1
    $\begingroup$ Thanks a lot! May I suggest that the notation $(X|X<c, Y)$ is potentially confusing, as it might suggest that $Y$ is part of what is conditioned on? $\endgroup$ – Christoph Hanck Jun 24 '17 at 7:55
  • 1
    $\begingroup$ I don't think the notation "joint pdf of $( X|X<c, Y)$" is ambiguous. (a) it refers to the joint pdf; If the | nested both $X<c$ and $Y$, it would be univariate, not joint. (b) if the alternative was intended, then it would be written: $(X | \{X<c, Y=y\} )$. Either way, if you are suggesting $( (X|X<c), Y)$ is clearer, can certainly do that. $\endgroup$ – wolfies Jun 24 '17 at 14:03
  • 1
    $\begingroup$ You are right that the careful reader has nothing to worry about in your notation. That said, it is maybe just me, but I stumbled over this upon first reading your +1 answer. Your last suggestion avoids this, at of course the cost of more cumbersome notation. Probably it is fine as it is together with this discussion. $\endgroup$ – Christoph Hanck Jun 24 '17 at 14:08
4
$\begingroup$

Sorry for not delivering the details, but $$ \int_{-\infty}^c \phi(x) \; \Phi(x-\sqrt{2} c) \, dx = 2T(c, \sqrt{2}-1) $$ where $T$ is the Owen $T$-function.

This function is available in Mathematica/Wolfram and in the R package OwenQ.

library(OwenQ)
pr <- function(c){
  2*OwenT(c, sqrt(2)-1) / pnorm(c)
}
curve(Vectorize(pr)(x), from=-6, to=6)

enter image description here

Alternatively you can get the Owen $T$-function with the help of the cdf of the noncentral Student distribution:

owenT <- function(h, a) 1/2*(pt(a, 1, h*sqrt(1+a^2)) - pnorm(-h))

But this implementation is not reliable for large values of the noncentrality parameter h*sqrt(1+a^2).

| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ That is very cool. The only other time I have come across the OwenT function was as the CDF of a skewNormal ... which makes me think that the pdf obtained here can likely be characterised as a skew-Normal distribution, with some tweaking. $\endgroup$ – wolfies Jun 25 '17 at 20:10
  • $\begingroup$ Also, the Mma implementation of OwenT appears to have some numerical instability problems - may require sending in a bug report $\endgroup$ – wolfies Jun 25 '17 at 20:15
  • $\begingroup$ @wolfies Where did you observe this numerical instability ? $\endgroup$ – Stéphane Laurent Jun 25 '17 at 20:31
  • $\begingroup$ MyFunc[b_] := With[{c = b}, NIntegrate[(1 + Erf[c - z/2])/ (E^(z^2/4)*(2*Sqrt[Pi]*(1 + Erf[c/Sqrt[2]]))), {z, Sqrt[2]*c, Infinity}]]; Plot[{OwenT[c, 1 - Sqrt[2]], (-2^(-1))* CDF[NormalDistribution[0, 1], c]*MyFunc[c]}, {c, -3.8, -3}, PlotRange->All] $\endgroup$ – wolfies Jun 30 '17 at 18:11
  • $\begingroup$ @wolfies Indeed. The numerical value of OwenT[3.3, Sqrt[2]-1] returned by Mathematica 11.1 is not correct. But the one returned by Wolfram online is the right one (0.0002073998). $\endgroup$ – Stéphane Laurent Jul 8 '17 at 18:47

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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