I have a 1-D random walker. It starts at a value of x at time t=0. It ends with a value of y at time t=1. Along the way, it hits a low of w, and a high of z. I know nothing about its path except the endpoints and the min/max values achieved.

Suppose I want to simulate the behavior of this random walker. So at any given time t between 0 and 1, I need to be able to sample from a realistic distribution with the described endpoint and extremes. What would be an appropriate distribution for this task?

Properties of distribution function f(t)

  • f(0)=x
  • f(1)=y
  • The central tendency should be weighted towards x at values close to t=0, and should lean more and more towards y as t approaches 1.
  • The variance should scale according to w and z, should be greatest at t=0.5, and should go to zero at t=0 and t=1.
  • The the tails of the probability distributions should go to zero at w and z, without abrupt discontinuity (i.e. it shouldn't be truncated).

Possible options?

  • A beta distribution visually looks like it has the right properties, but I'm not entirely sure how to interpret its parameters as a random walk scenario.
  • I understand that Rayleigh distributions describe random walks, but I'm not sure how one might incorporate the concept of knowing the bounds and final value of the walk.

After a couple of days hammering away at this, I think I have a workable solution, reparameterizing a beta distribution. The key insight (from this answer) was to solve for the beta distribution parameters $\alpha$ and $\beta$ as functions of the mean $\mu$ and variance $\sigma^2$. Then, I just had to decide how $\mu$ and $\sigma^2$ should vary with $t$ on the interval $[0,1]$.

As a reminder:

$w$ = min boundary, $x$ = starting value, $y$ = final value, $z$ = max boundary

So then from $t=0$ to $t=1$, $\mu$ can travel linearly from $x$ to $y$. The simplest form of this would be:

$$\mu = x (1 - t) + yt$$

Next is the variance and standard deviation. At $t=0$ and $t=1$, the function value is exactly known (equivalent to the starting or end point mean), so $\sigma=0$. The simplest version of this would be:


In the middle ($t=0.5$), where the time is furthest from the two endpoints, we have the least information about the location of the random walker, so then the distribution should be approximately normal, spread across the entire interval. Thus, the min and max should each be about 3 standard deviations from the mean at that time: $$z-w=6*\sigma(t=0.5)$$ rearranges to: $$\sigma(t=0.5)=(\frac{z-w}6)$$ Solve the previous equation for $C$ in this scenario: $$\sigma(t=0.5)=C*0.5(1-0.5)$$ $$C=4*\sigma(t=0.5)$$ $$C=4(\frac{z-w}6)$$ Finally: $$\sigma=4(\frac{z-w}6)*t(1-t)$$ But the beta distribution is only defined on the interval $[0,1]$, so we can just set $w$ and $z$ to $0$ and $1$ respectively, and rescale everything later: $$\sigma=4(\frac{1}6)*t(1-t)$$ $$\sigma=\frac{2}3t(1-t)$$

Hence, keeping in mind that $x$ and $y$ are predefined (and should be between $w$ and $z$ on the $[0,1]$ interval, to be rescaled later), these 4 equations fully determine the $\alpha$ and $\beta$ parameters of the beta distribution based on $t$, with all of the desired behaviors and constraints of a 1-D random walk from $x$ to $y$: $$\mu = x (1 - t) + yt$$ $$\sigma=\frac{2}3t(1-t)$$ $$\alpha=\left(\frac{1-\mu}{\sigma^2}-\frac{1}{\mu}\right)\mu^2$$ $$\beta=\alpha\left(\frac{1}{\mu}-1\right)$$ Plugging in some values of $x$, $y$, and $t$, the plots of these beta distributions give a very convincing simulation!

  • $\begingroup$ Actually, some edge cases (e.g. when x=y=0.2) don't look great, because when the mean is very close to one of the boundaries, too much probability density is squeezed into that small space. The random walker shouldn't actually be spending equal time on either side of x and y if x and y are skewed away from the center of the distribution. I think I can fix this by reformulating the equations in terms of the mode and variance, instead of the mean and variance. $\endgroup$ Feb 16 at 19:21

A friend helpfully suggested that this is similar to a Brownian Bridge, and I was able to use that to derive an analytical solution.

Conditions: $w$ = min boundary, $x$ = starting value, $y$ = final value, $z$ = max boundary

1. Incorporating the x and y conditions using a Brownian Bridge

The generalized Brownian Bridge is a specific Wiener Process, defined as a time-dependent normal distribution: $$B(t)=N(\mu_t,\sigma^2_t)$$ that is tethered at either end by the conditions: $$B(t_1)=x,\; B(t_2)=y$$ The mean and variance for this distribution are: $$\mu_t=x+\frac{t-t_1}{t_2-t_1}(y-x)$$ $$\sigma^2_t=\frac{(t_2-t)(t-t_1)}{t_2-t_1}$$ When our time interval is $[0,1]$, these reduce to: $$\mu_t=x(1-t) + yt$$ $$\sigma^2_t=t(1-t)$$ But then we have to condition our distribution on the observed min and max values, $w$ and $z$.

2. Min/Max Conditional

The total distribution of points observed across the entire $t\in[0,1]$, $B_T$ will also be normal (Brownian Bridges are proven to be normally distributed at all scales). We want to find the normal distribution that would be expected to have $w$ and $z$ as its min and max observed values.

It can be proven using the joint cumulative distribution functions (CDFs) that $B_T$ has a mean $\mu_T=\frac{1}2(w+z)$ halfway between the endpoints, and a standard deviation $\sigma_T=k(n)*(z-w)$ that depends on the sample size $n$.

For example, a well-known result derived from the normal CDFs is that 99.73% of the normal probability distribution lies within $3\sigma$ on either side of the mean. If we assume that $w$ and $z$ are the only two values observed outside this range, then $0.0027 = \frac{2}n$, giving $n=741$. Since this is the case where $2*3\sigma_T=z-w$, we can use this in our previous standard deviation equation to find $k(741)=\frac{1}6$.

Generalizing, where $p$ is the tail of the normal distribution containing a single extremum, and $z_p$ is the quantile: $$2(1-p) = \frac{2}n$$ $$p=\frac{n-1}n$$ $$2*z_p\sigma_T=z-w$$ $$\sigma_T=\frac{z-w}{2z_p}$$ We can use the quantile function of the normal distribution to find $\sigma_T$ as a function of $n$: $$z_p=\sqrt2\text{erf}^{-1}(2p-1)$$ $$z_p=\sqrt2\text{erf}^{-1}\left(\frac{n-2}n\right)$$ $$\sigma_T=\frac{z-w}{2\sqrt2\text{erf}^{-1}\left(\frac{n-2}n\right)}$$ Or, we can just use quantile tables that give us $z_p$ for values of $p\approx p(n)$.

3. Bringing it all together

Our solution will be a particular Wiener Process distribution that behaves like the Brownian Bridge, but also adds up over time to the total distribution described in Part 2:

$$W(t)\sim B(t)$$ $$B_T=\int_{0}^{1}W(t)dt$$

A helpful property is that the product of two independent normal distributions is another normal distribution. So, we can decompose our desired model into the joint probability of $B(t)$ and an independent, time-invariant normal distribution $N_C$: $$W(t)\ = N_CB(t)$$ $$B_T=N_C\int_{0}^{1}B(t)dt$$ The time-integral of a Brownian or Wiener process is another specific normal distribution: $$B_T=N_CN(0,\tfrac13)$$ In general, the product of two normal distributions is the normal distribution with mean and variance as: $$\mu=\frac{\mu_1\sigma^2_2+\mu_2\sigma^2_1}{\sigma^2_1+\sigma^2_2}$$ $$\sigma^2=\frac{\sigma^2_1\sigma^2_2}{\sigma^2_1+\sigma^2_2}$$ So then $N_C$ is fully specified by: $$\sigma^2_T=\frac{\sigma^2_C\tfrac13}{\sigma^2_C+\tfrac13}$$ $$\sigma^2_C=\frac{\sigma^2_T}{1-3\sigma^2_T}$$ and $$\mu_T=\frac{\tfrac13\mu_C+0}{\sigma^2_C+\tfrac13}$$ $$\mu_C=\mu_T(3\sigma^2_C+1)$$ $$\mu_C=\frac{\mu_T}{1-3\sigma^2_T}$$

Finally, we can use the same property to fully specify our model $W(t)=N_CB(t)=N(\mu_W,\sigma_W^2)$ in terms of the values solved for in Parts 1 and 2: $$\mu_W=\frac{\mu_C\sigma^2_t+\mu_t\sigma^2_C}{\sigma^2_C+\sigma^2_t}$$ $$\sigma_W^2=\frac{\sigma^2_C\sigma^2_t}{\sigma^2_C+\sigma^2_t}$$


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.