# Expected value and Variance of a stopped random process

Consider the random process where you keep drawing samples from [0,1] uniformly at random as long as the current sample is larger than the previous sample. What are the Expected value and the Variance of the number of samples you draw?

I have tried to think of a discrete case with Wald's Equality with step of 0.1. Summing up p * (1 - p) where p is 0.1, 0.2 ..., equals 0.327, matches my simulation mean=array([0.359285]), variance=array([0.05881588]).

import numpy as np
import matplotlib.pyplot as plt

from scipy import stats

cur_num = None
M = 1000000
li =[]
for i in range(M):

if cur_num is None:
cur_num = np.random.uniform(0,1,1)
tmp = 0

while tmp < cur_num:
tmp = cur_num
cur_num = np.random.uniform(0,1,1)
# print(f" while loop: {cur_num}")

# print(f" for loop: {cur_num}")

li.append(cur_num)

cur_num =None

print(stats.describe(li))

## DescribeResult(nobs=10000000, minmax=(array([1.29312049e-08]), array([0.9996544])),
## mean=array([0.359285]), variance=array([0.05881588]), skewness=array([0.44984084]), kurtosis=array([-0.76854895]))


This gives me the idea of integrate p * (1 - p) dp from 0 to 1. But it gives 1/2 - 1/3, 0.16667, not 0.35

Although Math SE might be a bit more suited for this one, I wanted to give it a try.

The answer relies on the law of total expectation, the law of total variance, and the relationship between Euler's number $$e$$ and series involving the factorial function $$n!$$.

In the first half of the answer, I will present the derivation of $$E(n)$$ and $$Var(n)$$, i.e. mean and variance of the stopping time, in the second half I will show how to derive $$E(x_n)$$ and $$Var(x_n)$$, i.e. the expected value and variance of the stopped process. Finally, in the addendum, I provide the full distribution of the stopped process, $$f(x_n)$$ for $$x_k$$ drawn from a uniform and from a generic distribution.

# Part I: The distribution of the stopping time $$n$$

Given the $$n$$th draw $$x_{n}$$ from the unit interval $$(0,1)$$, we accept the ($$n+1$$)th draw $$x_{n+1}$$ if $$x_n; else we stop at the $$n$$th draw and the process is said to be stopped at $$n$$ with value $$x_n$$. We say that $$n$$ is the stopping time.

Let’s say we are about to draw $$k$$ numbers $$x_1, x_2, \ldots x_k$$ uniformly from $$(0,1)$$ and arrange the results in increasing order, e.g.

\begin{align} x_1&

Then there are $$k!=k\times(k-1)\times(k-2)\ldots\times2\times1$$ possible permutations of the ordering of the draws, each permutation with the same probability $$1/k!$$ - but only one permutation represents the necessary ordering $$x_1 under which the process is not yet stopped at draw $$k$$ (of course, it will be stopped if $$x_k>x_{k+1}$$). Thus, the unconditional probability of $$n\geq k$$, i.e. of drawing at least $$k$$ times, is

\begin{align} P(n\geq k)&=\frac{1}{k!} \\ \Rightarrow P(n=k)&=P(n\geq k)-P(n\geq k+1)\\ &=\frac{1}{k!}-\frac{1}{(k+1)!} \end{align}

### Calculating the expected stopping time

We calculate $$E(n)=\sum_k k\times P(n=k)$$, rearranging the factorials in the meantime:

\begin{align} E(n)&=\sum\limits_{k=1}^\infty \frac{k}{k!}-\frac{k}{(k+1)!}\\ &=\left(\frac{1}{1!}+\frac{2}{2!}+\frac{3}{3!}+\ldots\right)-\left(\frac{1}{2!}+\frac{2}{3!}+\frac{3}{4!}+\ldots\right)\\ &=\frac{1}{1!}+\frac{1}{2!}+\frac{1}{3!}+\ldots\\ &=\sum\limits_{k=1}^{\infty}\frac{1}{k!}\\ &=\sum\limits_{k=0}^{\infty}\frac{1}{k!}-1 \end{align}

In the last row, we recognize the definition of Euler's constant $$e=\sum_k 1/k!=2.718282\ldots$$; hence:

$$E(n)=e-1$$

### Calculating the variance of the stopping time

As $$Var(x)=E(x^2)-E(x)^2$$, we now try to estimate $$E(n^2)$$, again applying some rearranging, and observing relationships with the definition of $$e$$:

\begin{align} E(n^2)&=\sum\limits_{k=1}^\infty k^2\times P(n=k)\\ &=\sum\limits_{k=1}^\infty\frac{k^2}{k!}-\frac{k^2}{(k+1)!}\\ &=\left(\frac{1^2}{1!}+\frac{2^2}{2!}+\frac{3^2}{3!}+\ldots\right)-\left(\frac{1^2}{2!}+\frac{2^2}{3!}+\frac{3^2}{4!}+\ldots\right)\\ &=\sum\limits_{k=1}^{\infty}\frac{k^2-(k-1)^2}{k!}\\ &=\sum\limits_{k=1}^{\infty}\frac{2k-1}{k!}\\ &=2\sum\limits_{k=1}^{\infty}\frac{k}{k!}-\left(e-1\right)\\ &=2\sum\limits_{k=0}^{\infty}\frac{1}{k!}-\left(e-1\right)\\ &=e+1 \end{align}

Thus

\begin{align} Var(n)=&E(n^2)-E(n)^2\\ &=e+1-(e-1)^2\\ &=3e-e^2 \end{align}

### A simulation of the stopping time $$n$$

Let's simulate this using R. I choose a cutoff level n=14, the corresponding probability is below 1E-10.

set.seed(42)
nSim <- 1E6
k    <- 14
mu <- exp(1)-1
v  <- 3*exp(1)-exp(2)
u<-sapply(1:nSim,function(i){
x <- rle( diff( runif(k) ) > 0 )
1 + x$$values[1] * x$$lengths[1]
})

print( mean(u)- mu )
print( var(u) - v )


with output

[1] -0.0001078285
[1] -0.0003605153


## Part II: Expectation and variance of $$x_n$$

Here, we will make use of the law of total expectation and the law of total variance.

### Expectation of $$x_n$$

We start with the law of total expectation, applying it to $$x_n$$ as

\begin{align} E(x_n)&=E\left(E\left(x_n|n=k\right)\right)\\ &=\sum\limits_{k=1}^{\infty}E\left(x_n|n=k\right)P(n=k) \end{align}

We already know the distribution of $$n$$, all that remains is calculating $$E(x_n|n=k)$$.

For a (fixed) stopping time $$k$$, it must hold that $$x_1\leq x_2\leq\ldots \leq x_{k-1}\leq x_k > x_{k+1}$$, i.e. the $$k$$th component from $$k+1$$ draws represents the maximum over $$(k+1)$$ draws. The distribution of the maximum is:

\begin{align} F^{max}_{k+1}(x)&\equiv P(\max(x_1,x_2,\ldots,x_k,x_{k+1})\leq x)\\ &=P(x_1\leq x, x_2\leq x, \ldots x_{k-1}\leq x,x_{k+1}\leq x)\\ &=F(x)^{k+1}\\ &=x^{k+1} \end{align}

The corresponding density, $$f^{max}_{k+1}(x)=\partial F^{max}_{k+1}(x)/\partial x$$, is $$f^{max}_{k+1}(x)=(k+1)x^k$$

and the expectation of $$x_n$$ given stopping at $$k$$ is

\begin{align} E(x_n|n=k)&=\int\limits_0^1 xf^{max}_{k+1}(x)\mathrm{d}x\\ &=(k+1)\int\limits_0^1 x^{k+1}\mathrm{d}x\\ &=\frac{k+1}{k+2} \end{align}

We now have both ingredients for the expectation:

\begin{align} E(x_n)&=\sum\limits_{k=1}^{\infty}E\left(x_n|n=k\right)P(n=k)\\ &=\sum\limits_{k=1}^{\infty}\frac{k+1}{k+2}\left(\frac{1}{k!}-\frac{1}{(k+1)!}\right)\\ &=\left(\frac{2/3}{1!}+\frac{3/4}{2!}+\frac{4/5}{3!}+\frac{5/6}{4!}+\ldots\right)-\left(\frac{2/3}{2!}+\frac{3/4}{3!}+\frac{4/5}{4!}+\frac{5/6}{5!}+\ldots\right)\\ &=0.5+\sum\limits_{k=1}^{\infty}\frac{1}{(k+2)!}\\ &=0.5+\sum\limits_{k=1}^{\infty}\frac{1}{k!}-1.5\\ &=\sum\limits_{k=1}^{\infty}\frac{1}{(k)!}-1\\ &=e-2 \end{align} where we have reused the relationships between the inverse factorial sum and Euler's number, above.

### Variance of $$x_n$$

We make use of the law of total variance

$$Var(x_n)=E(Var(x_n|n=k))+Var(E(x_n|n=k))$$

\begin{align} E(Var(x_n|n=k))&=\sum\limits_{k=1}^{\infty}\left(E(x_n^2|n=k)-E(x_n|n_k)^2\right)P(n=k)\\ &=\sum\limits_{k=1}^{\infty}\left((n+1)\int_{x=0}^1x^2f^{max}_{k+1}(x)\mathrm{d}x-\left(\frac{k+1}{k+2}\right)^2\right)P(n=k)\\ &=\sum\limits_{k=1}^{\infty}\left(\frac{k+1}{k+3}-\left(\frac{k+1}{k+2}\right)^2\right)P(n=k) \end{align}

and stop here. For the second term, we get

\begin{align} Var(E(x_n|n=k))&=E(E(x_n|n=k)^2)-E(E(x_n|n=k))^2\\ &=\sum_{k=1}^{\infty}\left(\frac{k+1}{k+2}\right)^2P(n=k)-\left(\sum_{k=1}^{\infty}\frac{k+1}{k+2}P(n=k)\right)^2\end{align}

Combining both terms yields

\begin{align} Var(x_n)&=E(Var(x_n|n=k))+Var(E(x_n|n=k))\\ &=\sum\limits_{k=1}^{\infty}\frac{k+1}{k+3}P(n=k)-\left(\sum_{k=1}^{\infty}\frac{k+1}{k+2}P(n=k)\right)^2 \end{align}

which, after some more rearranging yields

$$Var(x_n)=2+2e-e^2$$

### A simulation study of $$x_n$$

set.seed(42)
nSim <- 1e6
k    <- 14

u <- sapply(1:nSim,function(i){
x <- runif(k)
z <- rle(diff(x)>0)
x[1 + z$$values[1]*z$$lengths[1]]
})
cat("Variable  : "  ,"\t", "simulation", "\t", "theoretical"    , "\n")
cat("mean(x_n) : "  ,"\t",  mean(u)    , "\t", exp(1)-2         , "\n")
cat("var(x_n)  : "  ,"\t",  var(u)     , "\t", 2+2*exp(1)-exp(2), "\n")


resulting in

Variable  :      simulation      theoretical
mean(x_n) :      0.7182006       0.7182818
var(x_n)  :      0.04752103      0.04750756


We can even produce the distribution of $$x_n$$ in closed form by combining the PMF of $$n$$ and the density $$f^{max}_{k+1}(x)$$:

\begin{align} f(x)&=\sum\limits_{k=1}^{\infty}f(x_n|n=k)P(n=k)\\ &=\sum\limits_{k=1}^{\infty}(k+1)(x)^k\left(\frac{1}{k!}-\frac{1}{(k+1)!}\right)\\ &=\sum\limits_{k=1}^{\infty}(x)^k\frac{1}{(k-1)!}\\ &=x\sum\limits_{k=0}^{\infty}(x)^k\frac{1}{k!}\\ &=xe^{x} \end{align}

By the same line of reasoning we can show that the distribution of $$x_n$$ given draws from a continuous distribution $$\Phi$$ (and corresponding density $$\phi$$) equals

$$f(x)=\phi(x)\Phi(x)e^{\Phi(x)}$$

• Thanks! I only realised now I have to derive E(N), not E(process)! Supposed I have to derive E(process),shouldn't it be (e-1)* E(X_uniform) = (e-1)*(0.5) = 0.859? Why is it not 0.359 as suggested from the simulation? Sep 19, 2021 at 8:42
• I think he also wants the expected value of the most recent U(0.1) given that the process has stopped
– dm63
Sep 19, 2021 at 11:19
• @Kermittfrog Yes! Thank you so much! The proof you provide is very clear. Sep 20, 2021 at 10:49
• The loop is more clear (to me) if written like this: u <- sapply(1:nSim,function(i){ x <- runif(k); x[which(diff(x)<0)[[1]]] }) . It's a bit faster as well. Sep 20, 2021 at 11:28
• Hi @BobJansen thanks! I'll give it a try later; I am always looking for ways to improve speed and readability of my code :) Sep 20, 2021 at 11:35