![Expected value $\mathbb{E}[X-Y \mid X>Y, Y< c]$ for independent variables $X,Y$?](https://stats.stackexchange.com/Content/Sites/stats/Img/apple-touch-icon@2.png?v=344f57aa10cc)
There are several questions here. Before addressing them, let's look at the situation.
The gray and cyan dots in the figure display five (independently generated) datasets corresponding to five sets of parameters $[x_1,x_2]$ and $[y_1,y_2=c].$
The region $\mathcal E$ in question, which statisticians know as an event, is defined by the conditions $X\gt Y$ and $Y\lt c$ (along with the givens, $x_1\le X\le x_2$ and $y_1\le Y \le y_2.$ A simpler description is obtained upon replacing $y_2$ and $c$ by the lesser of those two values, so for now on I will use $y_2=c$ interchangeably.
The expectation of $X-Y,$ conditional on $\mathcal E,$ is readily found from the conditional expectation of the vector $(X,Y)$--and that is what the physicists know as the center of mass of a region. Everywhere you read "mass" you may substitute "total probability" or, in the case of these uniform 2D regions, "area." The concepts are mathematically identical and so, therefore, are the rules for working with them.
1. What is a correct formula?
Provided $\mathcal E$ has finite probability--which means finite area when the probability is uniform, as it is here--the center of mass can be expressed (abstractly) as a ratio of integrals,
$$\operatorname{CM}(\mathcal E) = E[(X,Y)\mid \mathcal E] = \frac{\int_{\mathcal E} (x,y)\,\mathrm dx \mathrm dy}{\int_{\mathcal E}\, \mathrm dx \mathrm dy}.$$
Because $X-Y$ is a linear function of $(X,Y);$ namely, $(X,Y)\cdot (1,-1)^\prime,$ the conditional expectation of $X-Y$ is obtained by applying the same linear function to this formula:
$$E[X-Y\mid \mathcal E] = \frac{\int_{\mathcal E} (x,y)\,\mathrm dx \mathrm dy}{\int_{\mathcal E}\, \mathrm dx \mathrm dy} \cdot (1, -1)^\prime =\frac{\int_{\mathcal E} x - y\,\mathrm dx \mathrm dy}{\int_{\mathcal E}\, \mathrm dx \mathrm dy}.$$
2. Is the formula in the question correct?
As we can see, the areas in question (the denominator of this formula) are given by $(x_2-x_1)(c-y_1)$ only in example (C), where the event is a rectangle. Comparing this to the expression given in the question (strongly) suggests that formula is incorrect.
3. What is a practical, insightful formula or algorithm?
I will limit this discussion to mathematics that is over two thousand years old.
We can exploit basic properties of centers of mass (CM) for regions of uniform (constant) density. These are elementary and well-known:
The CM of a right triangle is located at one-third of its height. More generally, the CM of any triangle determined by points $A,B,C$ is $(A+B+C)/3.$ This can be demonstrated in many ways and often is proven in a first course in Euclidean geometry.
The CM of a rectangle is its center. This follows from the observation that the center of mass must be located along any line of (reflection) symmetry. When opposite corners of the rectangle are points $A$ and $B,$ the CM is $(A+B)/2.$
The CM of two points $A_1$ and $A_2$ with corresponding masses $m_1$ and $m_2$ is their mass-weighted average (provided the total mass $m_1+m_2$ is nonzero), $$\frac{m_1A_1 + m_2A_2}{m_1+m_2}$$
The CM of two regions can be found by moving all the mass of each region to its center of mass and applying (3).
The foregoing properties imply (through simple algebra or by considering signed measures) that removing a region is equivalent to including a region of negative mass. (This is a nice application of a concept of "negative probability.") Specifically, when we remove $\mathcal F\subset \mathcal E$ of mass $m_F$ from $\mathcal E$ of mass $m_E,$ we may treat the remaining region $\mathcal E\setminus \mathcal F$ as consisting of the region $\mathcal E$ with mass $m_E$ together with the region $\mathcal F$ of mass $-m_F.$
Because we need to know masses, let's add these facts about areas in the plane:
The area of a rectangle determined by vertices $(x_1,y_1)$ and $(x_2,v)$ is $(x_2-x_1)(v-y_1).$
The area of a right triangle whose legs have lengths $a$ is $a^2/2.$
As you can see from the simulations in the figures, the problem amounts to finding the CM of a rectangle (shown with gray points in examples A, B, and C) and a triangle (cyan points) from which, possibly, a triangle has been removed (as in examples B and E). The foregoing rules permit us to express these regions as a combination of triangles (possibly with negative masses) and, optionally, a rectangle. As the illustration suggests, there are never more than three terms.
Software implementation
I found it convenient to use coordinates in which
The vertex $(x_2, y_2)$ (shown as yellow circles in the figure) is the origin.
The original $x$ coordinates are $0$ and $1.$ (Divisions by $x_2-x_1=1$ disappear, simplifying the calculations.)
The $y$ coordinates are all positive.
To illustrate how simple this solution is, the following R function implements the basic rules for computing CMs of triangles and rectangles and combining them, along with the area formulas, for the (limited) triangles and rectangles that can appear in this problem. It changes coordinates, applies these rules, and changes the coordinates back again, returning the CM and area that it has calculated.
The CMs as computed by this function are plotted as red triangles in the figure. They look correct--and longer simulations (involving $10^7$ sample points) bear out their accuracy, serving to verify the correctness of this solution.
CM <- function(x1, x2, y1, y2) {
cm.Triangle <- function(y) list(CM = y * c(1/3, 2/3), Area = y^2 / 2)
cm.Rectangle <- function(y) list(CM = c(1, y + 1)/2, Area = pmax(y - 1, 0))
cm <- function(X, Y, w = c(1, 1)) {
a = X$A * w[1] + Y$A * w[2]
if (abs(a) <= 2^(-52))
list(CM = c(0, 0), Area = 0)
else
list(CM = (X$CM * X$A * w[1] + Y$CM * Y$A * w[2]) / a, Area = a)
}
f <- function(y) cm(cm.Triangle(pmin(y, 1)), cm.Rectangle(y))
#
# Standardize coordinates.
#
s <- x2 - x1
eta <- c(pmax(0, x2 - y2), x2 - y1) / s
#
# Compute.
#
Z <- cm(f(eta[1]), f(eta[2]), w = c(-1, 1))
#
# Convert to original coordinates.
#
list(CM = (c(1, 1) - Z$CM) * s + x1, Area = Z$Area * s^2)
}
This solution has one limitation: it will not correctly compute the CM when the area of $\mathcal E$ is zero. That either can be treated separately or you can thicken the region by adding a positive infinitesimal to $y_2.$
Appendix
Because somebody always asks, here is the R code to generate the figure. It shows how the simulations are performed and how to use the CM function.
Examples <- list(
A = c(0, 1, -1/2, 3/2) * 2,
B = c(0, 1, -1/2, 3/5) * 2,
C = c(0, 1, -2/3, -1/3) * 3,
D = c(0, 1, 1/5, 1) * 5 - 1,
E = c(0, 1, 1/4, 3/4) * 4 - 4
)
n <- 1e4
d <- 1 + ceiling(log10(n) / 2) # Digits for printing estimates
set.seed(17)
args <- par(mfrow = c(1, length(Examples)), mar = c(2, 3, 3, 2))
for (name.p in names(Examples)) {
p <- Examples[[name.p]]
#
# Compute.
#
x1 <- p[1]; x2 <- p[2]; y1 <- p[3]; y2 <- p[4]
with(CM(x1, x2, y1, y2), {r <<- CM; a <<- Area})
#
# Simulate.
#
xlim <- range(c(x1, x2))
ylim <- range(c(y1, y2, x2))
N <- ceiling(n * diff(range(xlim)) * diff(range(ylim)) / a)
x <- runif(N, x1, x2)
y <- runif(N, y1, y2)
i <- which(x > y) # The condition
#
# Compare.
#
print(round(cbind(x = c(Calculation = r[1], Simulation = mean(x[i])),
y = c(Calculation = r[2], Simulation = mean(y[i]))), d))
#
# Plot.
#
j <- head(i, 1500 * n / N) # Don't plot too many points
plot(x[j], y[j], xlim = xlim, ylim = ylim, asp = 1,
main = name.p, bty = "n", xaxt = "n", yaxt = "n",
col = ifelse(y[j] - x1 < 0, hsv(0, 0, 0, 0.25), hsv(0.45, 0.9, 0.5, 0.25)),
cex = 0.75, xlab = "", ylab = "")
#
# Apply decorations and labels.
#
abline(h = c(y1, y2), v = c(x1, x2), lty = 3)
abline(0:1)
points(r[1], r[2], pch = 24, cex = 1.5, bg = "Red")
points(x2, x2, pch = 21, cex = 1.25, bg = "Yellow")
mtext(c(bquote(.(signif(x1, 2))), bquote(.(signif(x2, 2)))),
side = 1, line = 0.5, at = c(x1, x2))
mtext(c(bquote(.(signif(y1, 2))), bquote(.(signif(y2, 2)))), las = 1,
side = 2, line = 0.5, at = c(y1, y2))
}
par(args)