Commit 1a67a243 authored by Matthieu Boileau's avatar Matthieu Boileau
Browse files

Minor updates in srw

parent 375c44d6
%% Cell type:markdown id: tags:
# Course of Marielle Simon: Simple random walk
%% Cell type:markdown id: tags:
First, some python initializations.
%% Cell type:code id: tags:
``` python
%matplotlib inline
from matplotlib import rcParams
from matplotlib import pyplot as plt
from ipywidgets import interact
import ipywidgets as widgets
rcParams['figure.figsize'] = (8., 6.) # Enlarge figure
rcParams['animation.html'] = 'html5' # to render animation in notebook
import srw # Import simple random walk module from srw.py
```
%% Cell type:markdown id: tags:
## Random walk on $\mathbb{Z}^2$
%% Cell type:markdown id: tags:
Create and play a matplotlib animation for a $nstep$-step random walk starting at $(x, y) = (0, 0)$.
%% Cell type:code id: tags:
``` python
walk = srw.Walk2D(nstep=100) # Create a 100-step random walk
anim = walk.generate_animation()
plt.close(anim._fig) # Close the initial figure to display only the animation figure
anim # Now play
```
%% Cell type:markdown id: tags:
Plot entire path for various $nstep$ values.
%% Cell type:code id: tags:
``` python
interact(srw.plot_walk, nstep=widgets.IntSlider(min=50, max=1000, step=50, value=50));
interact(srw.plot_walk, nstep=widgets.IntSlider(min=100, max=1000, step=100, value=50));
```
%% Cell type:markdown id: tags:
### Some quantities as a function of the number of steps
Compute average **final distance** over 1000 random walks.
%% Cell type:code id: tags:
``` python
srw.FinalDistance(nwalk=1000).plot()
```
%% Cell type:markdown id: tags:
Compute average **maximum distance** over 1000 random walks.
%% Cell type:code id: tags:
``` python
srw.MaxDistance(nwalk=1000).plot()
```
%% Cell type:markdown id: tags:
Compute the **number of times** the walk goes back to starting point (average over 1000 random walks).
%% Cell type:code id: tags:
``` python
srw.BackToStart(nwalk=10000, nstepmax=10000).plot()
```
%% Cell type:markdown id: tags:
## Random walk on $\mathbb{Z}$
%% Cell type:markdown id: tags:
Consider the random walk on $\mathbb{Z}$ with $0 < p < 1$, denoted by $(S_n)$. The chain is supposed to start from state 0.
1\. Implement a function `random_walk_z` simulating the behaviour of the random walk for $n_{\max}$ steps, and represent it on a graph. Ensure that the function `random_walk_z`also returns: #
- the number of times the chain is returned to the initial state;
- the largest state reached by the chain.
%% Cell type:code id: tags:
``` python
import numpy as np
from sklearn.utils import check_random_state
from scipy.special import binom
import multiprocessing as mp
mp.set_start_method('spawn', True) # see https://github.com/microsoft/ptvsd/issues/1443
from numba import jit
@jit(nopython=True)
def count_first(item: int, vec: np.ndarray):
"""
Find the index of the first element in the array `vec` equal to the element `item`.
"""
c = 0
for i in range(len(vec)):
if item == vec[i]:
c += 1
return c
def random_walk_z(p, n_max, random_state):
""" Simulate a simple 1D random walk in Z.
:returns:
- Ti (:py:class:`int`) - number of returns to the initial state
- state_max (:py:class:`int`) - farthest state reached by the chain (w.r.t the initial state)
"""
rng = check_random_state(random_state)
Z = 2*rng.binomial(1, p, size=(n_max)) - 1
X = np.empty(shape=(n_max+1), dtype=float)
X[0] = 0
X[1:] = np.cumsum(Z)
Ti = count_first(0, X[1:])
id = np.argmax(np.abs(X))
state_max = X[id]
t = np.arange(0, n_max+1, 1)
plt.plot(t, X)
plt.show()
return Ti, state_max
```
%% Cell type:code id: tags:
``` python
random_walk_z(0.5, 1000, 500)
```
%% Cell type:markdown id: tags:
2\. Assume now that two players $A$ and $B$ play heads or tails, where heads occur with probability $p$. Player $A$ bets $1$ euro on heads at each toss, and $B$ bets $1$ euro on tails. Assume that:
- the initial fortune of $A$ is $a \in \mathbb{N}$;
- the initial fortune of $B$ is $b\in\mathbb{N}$;
- the gain ends when a player is ruined.
Implement a function which returns the empirical frequency of winning for $A$, and compare it with the theoretical probability computed in the lecture.
%% Cell type:code id: tags:
``` python
# Exercise
```
......
......@@ -19,21 +19,19 @@ class Walk2D:
def __init__(self, nstep: int):
self.nstep = nstep
self.x, self.y = self._generate_walk()
self.x, self.y = self._compute_walk(self.nstep)
def _generate_walk(self):
"""
Create a nstep-random walk path on Cartesian grid
"""
x = np.empty(self.nstep)
y = np.empty(self.nstep)
@staticmethod
@jit(nopython=True)
def _compute_walk(nstep: int):
x = np.empty(nstep)
y = np.empty(nstep)
# initial position
x[0] = 0
y[0] = 0
# time loop
for step in range(1, self.nstep):
for step in range(1, nstep):
direction = directions[np.random.randint(4)]
x[step] = x[step - 1] + direction[0]
y[step] = y[step - 1] + direction[1]
......@@ -117,7 +115,7 @@ class NWalk:
nwalk: number of random walks for averaging
nstepmax: maximum final step
num_step: number of nsteps to compute and plot
"""
"""
self.nwalk = nwalk
self.nsteps = np.linspace(1, nstepmax, num=step_num, endpoint=True,
dtype=int)
......@@ -247,8 +245,8 @@ class BackToStart(NWalk):
num_step: number of nsteps to compute and plot
"""
self.nwalk = nwalk
self.nsteps = np.logspace(2, np.log10(nstepmax), num=step_num, endpoint=True,
dtype=int)
self.nsteps = np.logspace(2, np.log10(nstepmax), num=step_num,
endpoint=True, dtype=int)
@staticmethod
@jit(nopython=True)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment