Python Homework 4: Arcsine Laws

The purpose of this python homework is to explore the so-called Arcsine laws numerically. The arcsine laws are a number of fascinating results for random walks. They relate path properties of the simple symmetric random walk to the arcsine distribution.

A random variable XX on [0,1][0,1] is arcsine-distributed if the cumulative distribution function is given by

[Xx]=2arcsin(x)P[Xx]=2arcsin(x) for all 0x10x1 and the probability density function is given by fX(x)=1x(1x)fX(x)=1x(1x) on (0,1)(0,1).

Given a simple symmetric random walk (Sn)n0(Sn)n0 with S0=0S0=0, we define the following random variables:

The total number of periods from 00 to 2N2N the random walk spends above zero:

C2N:=|{n{1,,2N}:Sn>0}|.C2N:=|{n{1,,2N}:Sn>0}|.

The time of the last visit to 00 before time 2N2N:

L2N:=max{0n2N:Sn=0}.L2N:=max{0n2N:Sn=0}.

The time when the random walk reaches its unique maximum value between time 00 and 2N2N:

M2N:=argmax{Sn:0n2N}(this notation means that SM2N=max{Sn:0n2N}.M2N:=argmax{Sn:0n2N}(this notation means that SM2N=max{Sn:0n2N}.

As usual, we start with loading some packages:

In [1]: import numpy as np import math import matplotlib.pyplot as plt from scipy.stats import arcsine # Allows to render plots directly within the notebook %matplotlib inline

Problem 1 (1 points)

Write a Python code to simulate one path of a simple symmetric random walk with 2N2N time steps (i.e. lengthOfRandomWalk = 2N) starting at S0=0S0=0.

The function randomWalk should return an array of length 2N with the entire path of the random walk.

In [6]: import random def randomwalk(n): random.randint(0,1) possiblemoves = [-1,1] #set from -1 to 1 because it will extend all possiblites on the integer line simplepath = [] simplepath.append(0) for i in range(2*int(n)): simplepath.append(int(simplepath[int(i)])+int(possiblemoves[random.randint(0,1)])) #computing the number of steps in the random path return simplepath n = input("Enter number of steps (n) : ") print(randomwalk(n)) Enter number of steps (n) : 5 [0, 1, 0, 1, 0, 1, 2, 1, 2, 3, 2]

Problem 2 (6 points)

Write three Python functions which take a path of a random walk as an input and then compute the values of the random variables C2N,L2N,M2NC2N,L2N,M2N, respectively, as defined above:

In [2]: ### Insert your code here for total time spend above zero def timeAboveZero(path): ## Write your own code here return result In [3]: ### Insert your code here for last time the random walk visited zero def TimeOfLastVisitOfZero(path): ## Write your own code here return result In [ ]: ### Insert your code here for the time when the random walk hits its maximum def timeHitMax(randomWalk): ## Write your own code here return result

Problem 3 (3 points)

Random Walk theory provides us with the following nice (anf maybe counterintuitive) asymptotic statements:

limN[C2N2Nx]=2arcsin(x)limNP[C2N2Nx]=2arcsin(x)limN[L2N2Nx]=2arcsin(x)limNP[L2N2Nx]=2arcsin(x)limN[M2N2Nx]=2arcsin(x)limNP[M2N2Nx]=2arcsin(x)

We say that the random variables C2N/2N,L2N/2N,M2N/2NC2N/2N,L2N/2N,M2N/2N converge in distribution to the Arcsine Distribution.

The interesting property about the Arcsine distribution is that its density (see its formula above) is U-shaped on (0,1)(0,1). In other words, if XX is arcsine-distributed on (0,1)(0,1), the probabilty that XX takes very small values near 0 or very large values near 1 is rather high, but the probability for taking values around, say, 0.5, is low.

In [11]: x = np.linspace(arcsine.ppf(0.05), arcsine.ppf(0.95), 100) plt.title("Density of the arcsine distribution") plt.plot(x, arcsine.pdf(x), linewidth=2, color='b') Out[11]: []

For 2N=10002N=1000 sample 10,000 realisations of each of the random variables C2N/2N,L2N/2N,C2N/2N,L2N/2N, and M2N/2NM2N/2N, respectively. Display a normalized histogram for all three simulations, along with the probability density function of the arcsine distribution, to check the above facts numerically!

In [ ]: ### Implement your Simulations here In [15]: ### Complete the plot commands accordingly c = arcsine.rvs(size=10000) # This has to be replaced by the simulated values for C_2N plt.figure(figsize=(10,5)) plt.title("Normalized histogram for 10000 realisations of $C_{1000}$") plt.hist(c, bins='auto', normed='True') plt.plot(x, arcsine.pdf(x), linewidth=2, color='r', label="true arcsine density") plt.legend() l = arcsine.rvs(size=10000) # This has to be replaced by the corresponding simulated values for L_2N plt.figure(figsize=(10,5)) plt.title("Normalized histogram for 10000 realisations of $L_{1000}$") plt.hist(l, bins='auto', normed='True') plt.plot(x, arcsine.pdf(x), linewidth=2, color='r', label="true arcsine density") plt.legend() m = arcsine.rvs(size=10000) # This has to be replaced by the corresponding simulated values for M_2N plt.figure(figsize=(10,5)) plt.title("Normalized histogram for 10000 realisations of $M_{1000}$") plt.hist(m, bins='auto', normed='True') plt.plot(x, arcsine.pdf(x), linewidth=2, color='r', label="true arcsine density") plt.legend() Out[15]: