I often need to do calculations such as the following:
- $A, B, C, D$ are independent random variables distributed uniformly in [0,1].
- $X := A+B$ and $Y:=B+C+D$.
- What is the probability that $X<Y$?
Is there a free mathematical software that can do such calculations easily?
My favorite mathematical software is SageMath, but as far as I know, its abilities in probability theory are not sufficient for such calculations.
$\newcommand{\braces}[1]{\left\lbrace\,{#1}\,\right\rbrace} \newcommand{\bracks}[1]{\left\lbrack\,{#1}\,\right\rbrack} \newcommand{\dd}{\mathrm{d}} \newcommand{\ds}[1]{\displaystyle{#1}} \newcommand{\expo}[1]{\,\mathrm{e}^{#1}\,} \newcommand{\ic}{\mathrm{i}} \newcommand{\mc}[1]{\mathcal{#1}} \newcommand{\mrm}[1]{\mathrm{#1}} \newcommand{\pars}[1]{\left(\,{#1}\,\right)} \newcommand{\partiald}[3][]{\frac{\partial^{#1} #2}{\partial #3^{#1}}} \newcommand{\root}[2][]{\,\sqrt[#1]{\,{#2}\,}\,} \newcommand{\totald}[3][]{\frac{\mathrm{d}^{#1} #2}{\mathrm{d} #3^{#1}}} \newcommand{\verts}[1]{\left\vert\,{#1}\,\right\vert}$
You can do it with a simple $\texttt{javascript}$ code which runs in a terminal with $\texttt{node}$ as
$> \texttt{node random0.js}$
// random0.js Felix Marin
"use strict";
var a = null;
var c = null;
var d = null;
var i = 0;
var ITER = 1000000000; // Total number of iterations.
var total = 0;
while (i < ITER) {
a = Math.random();
c = Math.random();
d = Math.random();
if (a < (c + d)) ++total;
++i;
}
console.log("Result = " + total/ITER);
/*
https://math.stackexchange.com/questions/1959638/mathematical-software-for-calculations-with-random-variables
*/
Result = 0.833326178
The result can changes slightly between 'runs' of the script but they remain close to $\bbox[5px,#eee]{\texttt{0.8333}} \approx 5/6$.
In general, you can follow the $\texttt{LutzL}$ 'graphical answer or perform the following evaluation:
\begin{align} &\int_{0}^{1}\int_{0}^{1}\int_{0}^{1}\int_{0}^{1} \bracks{a + b < b + c + d}\,\dd d\,\dd c\,\dd b\,\dd a = \int_{0}^{1}\int_{0}^{1}\int_{0}^{1}\bracks{d > a - c}\,\dd d\,\dd c\,\dd a \\[5mm] = &\ \int_{0}^{1}\int_{0}^{1}\braces{% \bracks{a - c < 0}\int_{0}^{1}\,\dd d + \bracks{0 < a - c < 1}\int_{a - c}^{1}\,\dd d}\,\dd c\,\dd a \\[5mm] = &\ \int_{0}^{1}\int_{0}^{1}\bracks{c > a}\,\dd c\,\dd a + \int_{0}^{1}\int_{0}^{1}\bracks{a - 1 < c < a}\pars{1 - a + c}\,\dd c\,\dd a \\[5mm] = &\ \int_{0}^{1}\int_{a}^{1}\,\dd c\,\dd a + \int_{0}^{1}\int_{0}^{a}\pars{1 - a + c}\,\dd c\,\dd a = \int_{0}^{1}\pars{1 - a}\,\dd a + \int_{0}^{1}\pars{a - {1 \over 2}\,a^{2}}\,\dd a \\[5mm] = &\ {1 \over 2} + \pars{{1 \over 2} - {1 \over 6}} = \bbox[5px,border:1px groove navy]{5 \over 6} \end{align}
If you would like to do these calculations numerically, I think R is very good (everything free).
Here is a link: https://www.r-project.org/
Here is an editor you may like: https://www.rstudio.com/
A solution in R could look like the following
n <- 1000
A <- runif(n, min = 0, max = 1)
B <- runif(n, min = 0, max = 1)
C <- runif(n, min = 0, max = 1)
D <- runif(n, min = 0, max = 1)
X <- A+B
Y <- B+C+D
T <- (X < Y)
sum(T) / length(T) # I get 0.8295
5/6 # True value is 0.8333...
The true value is given by LutzL in the other post.
Set $A'=1-A$ which is equally a uniform random variable on $[0,1]$, then $X<Y$ is equivalent to $$1<A'+C+D.$$
The probability is the volume of the cube minus one edge pyramid. Which comes up to $$1-\frac16=\frac56.$$
