levyrandomstep_2D

import random
import numpy as np
import matplotlib.pyplot as plt
from datetime import datetime
from scipy.stats import levy

random.seed(datetime.now())

"""pick an initial random direction with uniform distribution"""
directions = ['+1', '-1']
direction_idx = random.randint(0, 1)
direction = directions[direction_idx]
print('initial random direction: '+direction+'\n')

# custom parameters
loc = 0
scale = 1
switch_threshold = 4
levy_max_extension = 1000
#

plot_resolution = 50000
x = np.linspace(0, levy_max_extension, plot_resolution)
y = levy.pdf(x, loc, scale)
plt.axvspan(switch_threshold, levy_max_extension, facecolor='b', alpha=0.05)

# plot normal distribution (loc, scale)
plt.plot(x, 1/(scale * np.sqrt(2 * np.pi)) *
    np.exp( - (x - switch_threshold)**2 / (2 * scale**2) ),
    linewidth=2, color='g', alpha=0.4)

# plot levy distribution
plt.plot(x, y, 'b-', lw=2, alpha=0.6)

# plot threshold
plt.axvline(x=switch_threshold, color='red')

# rescale axis
plt.axis([0,10,0,0.5])

for i in range(50):
    """direction inversion using a random pick from levy distribution"""
    r = levy.rvs(loc, scale)
    n = np.random.normal(switch_threshold, scale)

    plt.plot(r, 0.5, 'b+')
    plt.plot(n, 0.0, 'g+')

    # use > to invert direction on tail
    # use < to invert direction on head
    if r > switch_threshold:
        direction_idx = 1 if direction_idx==0 else 0
        print(str(r)+'\tinverted direction => '+directions[direction_idx])
    else:
        print(str(r)+'\t'+directions[direction_idx])

    print('\tnormal:', n, '\n')

plt.show()