Aquí tienes un intento de demostración. No tengo una referencia, ni mucha intuición más allá de esta construcción.
Primero considera el caso de una variable $X$ que toma el valor $a<0$ con probabilidad $b/(|a|+b)$ y toma el valor $b>0$ con probabilidad $|a|/(|a|+b).$ Escoge variables aleatorias $Y,Z$ distribuidas uniformemente en $[a,b],$ con $Y=Z+b$ cuando $Z<0$ y $Y=Z+a$ cuando $Z\geq 0.$ Entonces $Y-Z\stackrel{\mathrm d}= X$ como se requiere.
Para el caso general podemos tomar una mezcla, emparejando los eventos $X<0$ y $X>0$ de una manera adecuada.
Observación. Solo necesitamos escribir la ley de $X$ como un promedio ponderado de distribuciones de probabilidad sin media soportadas en dos puntos. Una manera elegante de hacer esto es usando el teorema de Choquet y una caracterización de los puntos extremos del conjunto de variables sin media, por ejemplo, en el artículo de Winkler sobre "Puntos Extremos de Conjuntos de Momentos". El resto de esta respuesta es simplemente una versión más explícita de esto.
Sea $\mu$ la ley de $X.$ Sea $M=\tfrac 1 2 \mathbb E|X|.$ Existe una función boreliana $f^+:[0,M]\to\mathbb R_{>0}$ tal que la imagen inversa de la medida de Lebesgue $f^+_*\lambda$ es la medida $x d\mu(x).$ Por ejemplo, la inversa generalizada de la función $t\mapsto \int_0^t x d\mu(x).$ De manera similar, existe una función boreliana $f^-:[0,M]\to\mathbb R_{<0}$ tal que $f^-_*\lambda$ es la medida $|x| d\mu(x).$ La imagen inversa de $(1/f^+(x))d\lambda(x)$ bajo $f^+$ es $\mu$ restringida a $\mathbb R_{>0}.$ Y de manera similar, la imagen inversa de $(1/|f^-(x)|)d\lambda(x)$ bajo $f^-$ is $\mu$ restringida a $\mathbb R_{<0}.$
Escoge las variables aleatorias independientes $R,S,T$ de la siguiente manera:
- $R$ es una variable aleatoria de Bernoulli con $\mathbb P[R=0]=\mathbb P[X=0],$
- $S$ es una variable aleatoria en $[0,M]$ distribuida según la medida $((1/f^+(x))+(1/|f^-(x)|))d\lambda(x),$
- $T$ está distribuida de forma uniforme en $[0,1].$
Luego define $Y$ y $Z$ de la siguiente manera:
- donde sea que $R=0,$ tomamos $Y=Z=0,$ y omitimos los siguientes pasos.
- $Z$ es $T$ reescalada de $[0,1]$ al intervalo $[f^-(S),f^+(S)].$
- $Y$ es $Z+f^+(S)$ cuando $Z<0$
- $Y$ es $Z+f^-(S)$ cuando $Z\geq 0.$
$Z$ tiene la misma distribución que $Y.$ La variable aleatoria $Y-Z$ tiene la misma probabilidad de ser cero que $X.$ Para cualquier boreliano $U\subseteq\mathbb R_{>0},$ la probabilidad de que $f^+(S)\in U$ y $Y-Z=f^+(S)$ es $\int (1/f^+(x))1_{f^+(x)\in U}d\lambda(x)=\mu(U).$ Por lo tanto, las leyes de $Y-Z$ y $X$ son las mismas en $\mathbb R_{>0}$. Y de manera similar para $\mathbb R_{<0}.$ Así que esta construcción asegura que $Y-Z\stackrel{\mathrm d}= X.$