Brownian motion–derived process: Rayleigh Distribution

Proceso derivado del movimiento Browniano

Definición:

Sean \(\{B_1(t), t\geq 0 \}\) y \(\{B_2(t), t\geq 0 \}\) dos movimientos Brownianos estándar independientes. Y sea \[ \begin{eqnarray} R(t) = \sqrt{ B_1(t)^2 + B_2(t)^2}, \ \ \ t \geq 0 \end{eqnarray} \] donde \(R(t)\) representa la distancia al origen de un movimiento Browniano bidimensional. \(R(t)\) tiene distribución Rayleigh \[ \begin{eqnarray} f(x) & = &\frac{x}{t} e^{-\frac{x^2}{2t}}, \ \ \ x>0 \\ & = & 0, \ \ \ d.o.m. \end{eqnarray} \] con valor esperado

\[\mathbb{E}[X(t)] = \sqrt{\frac{\pi t}{2}}\]

y varianza

\[\mathbb{V}ar(X(t)) = \frac{4-\pi}{2} t\]

Movimiento Browniano con deriva y volatilidad en una dimensión

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import matplotlib.cm as cm
import matplotlib.patches as patches
sns.set_theme()

# Function para simular una Movimiento Aleatorio 
class StochasticSimulation():
  
    def __init__(self, T=None, steps=100, n_times=1):
        self.steps = steps
        self.n_times = n_times 
        # Si T no se proporciona usar T = steps
        if T is None:
            self.T = steps
        else:
            self.T = T

        self.dt = self.T/ steps
        self.t = np.linspace(0, self.T, steps+1)
             
    # Movimiento Browniano
    def BrownianMotion(self):
        # creamos arreglo
        bm = np.zeros( (self.n_times, self.steps+1))
        # inicializa simulación
        for j in range(bm.shape[1]-1):
            for i in range(bm.shape[0]):
                bm[i][j+1] = bm[i][j] + np.sqrt(self.dt) * np.random.normal(size=1)
    
        df = pd.DataFrame(bm)  
        return df

En este ejemplo, primero se definen los parámetros de la simulación: se tomarán \(𝑁=10,000\) pasos, lo que permitirá generar una trayectoria muy detallada del proceso, y se simulará \(𝑀=1\) única trayectoria.

A continuación, se crea una instancia de la clase \(StochasticSimulation()\), pasando como argumentos el número de pasos y el número de trayectorias a simular. Esta clase se encarga de manejar los cálculos necesarios para generar procesos estocásticos.

Finalmente, se llama al método \(Brownian()\) com parámetros \(\mu\) y \(\sigma\) de la instancia, el cual genera la trayectoria simulada del movimiento Browniano. El resultado se almacena en df, un DataFrame que contiene los valores de la trayectoria a lo largo de los pasos de tiempo definidos.

Solución: Mil trayectorias del movimiento Browniano

# Parámetros
N = 10000        # número de pasos
M = 1000            # número de trayectorias simuladas

# Instanciamos las clases para dos procesos independientes
sim1 = StochasticSimulation(steps=N, n_times=M)
sim2 = StochasticSimulation(steps=N, n_times=M)
# generamos trajectorias de B1 y B2
B1 = sim1.BrownianMotion()
B2 = sim2.BrownianMotion()
# generamos trajectorias de Rt
Rt = np.sqrt(B1**2 + B2**2)
# Copia de Rt
df = Rt.copy()

# Colormap 
cmap = cm.get_cmap("Greys", M*2)
# lista de t
t = sim1.t
std_teorica = np.sqrt(((4 -np.pi)/2)*t)
# Media teórica 
mean_B = np.sqrt((t*np.pi)/2)
# Intervalos de confianza 95% 
conf_upper =  mean_B + 1.96 * std_teorica
#conf_lower =  mean_B - 1.96 * std_teorica
# Gráfico
plt.figure(figsize=(14,6))
for i in range(M):
    plt.plot(t, df.iloc[i,:].to_list(), color=cmap(i), linewidth=0.5)
plt.plot(t, conf_upper, 'blue', linewidth=0.3)
#plt.plot(t, conf_lower, 'blue', linewidth=0.3)
#plt.fill_between(t, conf_lower, conf_upper, color='lightblue', alpha=0.2, label='Intervalo de confianza 95% (teórico)')
plt.title(f"Simulación de R(t) [{M} trayectorias]")
plt.xlabel("Tiempo t")
plt.ylabel("R(t)")
plt.grid(True)
plt.legend()
plt.show()

Distribución al tiempo \(t\)

Vamos a hacer cortes verticales a distintos tiempos \(t\) y hechar un vistazo a las distribuciones. Para cada \(t\) se puede verificar que tanto la media como la varianza se hacercan a los correspondientes valores teóricos definidos al inicio.

# histogramas para distintos valores de t
valores_t = [10, 100, 1000, 5000]

# Crear figura con subplots 4x4
fig, axes = plt.subplots(2, 2, figsize=(14, 10))
axes = axes.flatten()

for i, step in enumerate(valores_t):
    data = df[step]
    mean = data.mean()
    var = data.var()

    sns.histplot(data, bins=30, kde=True, ax=axes[i])
    axes[i].set_title(f"Distribución al tiempo {step}")
    axes[i].set_xlabel('Posición')
    axes[i].set_ylabel('Frecuencia')
    axes[i].legend([f'Media = {mean:.2f}, Varianza = {var:.2f}'])

# Eliminar los subplots vacíos 
for j in range(len(valores_t), len(axes)):
    fig.delaxes(axes[j])

plt.tight_layout()
plt.show()

Ejemplo 1: para t=100

El valor teórico es:

\[\mathbb{E}[R(100)] = \sqrt{\frac{\pi t}{2}} = \sqrt{\frac{\pi *100}{2}} = 12.53\] \[\mathbb{V}ar(R(100)) = \sqrt{\frac{4- \pi}{2}}*t = \sqrt{\frac{4- \pi}{2}}*100 = 42.92\]

Simulación de un Movimiento Browniano en dos dimensiones

# Colormap (gradiente de gris)
cmap = cm.get_cmap("Greys", M*2)
centro = (0, 0)
circulo = patches.Circle(centro, 400, edgecolor='grey', facecolor='black', alpha=0.8, lw=1, ls='--')

# Crear la figura y los ejes
fig, ax = plt.subplots(figsize=(14, 10))
ax.add_patch(circulo)
# Brownian Motion
for i in range(df.shape[0] - 1):
    if i < 800:
        ax.plot(df.iloc[i][:], df.iloc[i + 1][:], linewidth=0.2, color=cmap(i))  # Normaliza el valor de i
    else:
        continue
plt.title('Simulación de R(t): {} steps and {} paths'.format(N, M))
plt.xlabel('R(t)')
plt.ylabel('R(t)')
plt.gca().set_xticklabels([])  # quitar etiquetas eje X
plt.gca().set_yticklabels([])  # quitar etiquetas eje Y
plt.grid(True)
plt.show()

Referencias

Ubbo F Wiersema (2008). Brownian Motion Calculus. John Wiley & Sons Ltd.
Sheldon M. Ross(1995). Stochastic Processes. Wiley Series in Probability and Statistics: Probability and Statistics. John Wiley & Sons.
LefebvreMario(2007). Applied Stochastics Processes. Springer, New York.