stochastic process homework

homework1

import numpy as np
import matplotlib.pyplot as plt
def random_walk_steps(n_steps):
    """
    生成n_steps步随机行走的变量
    """
    X =np.random.rand(n_steps)
    X[X>=0.5] = 1.0
    X[X<0.5] = -1.0
    return X
def generate_random_walk(steps,n):
    """
    steps 存放了每一步随机行走的情况
    生成 Y,注意到这里steps 和 n其实是两个独立的变量
    """
    Y0 = np.array([0])
    Y = np.concatenate((Y0,np.cumsum(steps)/np.sqrt(n)))
    return Y
def brownian_motion_simulate(T,dt):
    """
    函数把[0,T]这个时间段分成[0,dt],[dt,2dt]... 然后用正态分布对每个时间段进行模拟再求和
    """
    nsteps = int(np.ceil(T/dt))+1 # ceil是向上取整的意思
    steps = np.random.randn(nsteps)*np.sqrt(dt) #每一个时间步都是 N(0,dt)的随机变量
    steps[0] = 0.
    bmotion = np.cumsum(steps)
    return bmotion
T = 80.
dt = 0.001
nsteps = int(np.ceil(T/dt))+1
NPATH = 100
SAMPLES=np.zeros((NPATH,nsteps))
for ii in range(NPATH):
    SAMPLES[ii,:]=brownian_motion_simulate(T,dt)
## 作图
plt.figure(figsize=(7,7))
TSPAN = np.arange(0,nsteps*dt,dt)
plt.plot(np.tile(TSPAN,(NPATH,1)).T,SAMPLES.T,alpha=0.4)
plt.xlabel('Time',fontsize=14)
plt.ylabel(r'$W_t$',fontsize=14)
plt.show()
FirstOrderDifference=np.zeros((NPATH,nsteps))
for ii in range(NPATH):
    for i in range(1,nsteps):
        FirstOrderDifference[ii,i]= FirstOrderDifference[ii,i-1]+abs(SAMPLES[ii,i]-SAMPLES[ii,i-1])
plt.figure(figsize=(7,7))
TSPAN = np.arange(0,nsteps*dt,dt)
plt.plot(np.tile(TSPAN,(NPATH,1)).T,FirstOrderDifference.T,alpha=0.4)
plt.xlabel('Time',fontsize=14)
plt.ylabel(r'$S^{A,1} _t$',fontsize=14)
plt.show()
SecondOrderDifference=np.zeros((NPATH,nsteps))
for ii in range(NPATH):
    for i in range(1,nsteps):
        SecondOrderDifference[ii,i]=SecondOrderDifference[ii,i-1]+ abs(SAMPLES[ii,i]-SAMPLES[ii,i-1])**2
plt.figure(figsize=(7,7))
TSPAN = np.arange(0,nsteps*dt,dt)
plt.plot(np.tile(TSPAN,(NPATH,1)).T,SecondOrderDifference.T,alpha=0.4)
plt.xlabel('Time',fontsize=14)
plt.ylabel(r'$S^{W,W} _t$',fontsize=14)
plt.show()

homework2

import numpy as np
import matplotlib.pyplot as plt
# 参数
a, v, sigma,x0= 0.02, 0.3, 1.0,0


def Euler(T, m, n):
    # n: sample path的数目
    # m: 时间离散的步数
    # T: 仿真的终止时间

    h = T / m  # 时间步长

    Z = np.random.normal(0.0, 1.0, (n, m))  # 标准正态分布变量

    X = np.zeros((n, m + 1))  # 存储n条路径的矩阵

    X[:, 0] = x0  # 初始化

    for i in range(m):
        # 迭代：矩阵化实现
        X[:, i + 1] = X[:, i] + a*(v-X[:, i ])*h+sigma*np.sqrt(h)*Z[:,i]
    return X
T, m, n = 1.0, 1000, 10
X_Euler = Euler(T, m, n)
t_grid = np.linspace(0,T,m+1)

# 展示Euler Scheme的结果
for i in range(n):
    plt.plot(t_grid, X_Euler[i,:])

plt.xlabel('$t$')
plt.ylabel('$X_t$')
plt.title('Euler Scheme for BS')
plt.savefig('Euler_Scheme_BS.jpg')
plt.show()

homework3

import numpy as np
import matplotlib.pyplot as plt
# 参数
a, v, sigma,x0= 0.02, 0.3, 1.0,0
s0,sigma1,sigma2,thteta=0,0.01,0.01,0
def Euler(T, m, n):
    # n: sample path的数目
    # m: 时间离散的步数
    # T: 仿真的终止时间

    h = T / m  # 时间步长

    Z1 = np.random.normal(0.0, 1.0, (n, m))  # 标准正态分布变量
    Z2 = np.random.normal(0.0, 1.0, (n, m))  # 标准正态分布变量
    S=np.zeros((n, m + 1))
    X = np.zeros((n, m + 1))  # 存储n条路径的矩阵
    S[:,0]=s0
    X[:, 0] = x0  # 初始化

    for i in range(m):
        # 迭代：矩阵化实现
        X[:, i + 1] = X[:, i] + a*(v-X[:, i ])*h+sigma*np.sqrt(h)*Z1[:,i]
        S[:,i+1]=sigma1*np.sqrt(h)*Z1[:,i]+sigma2*np.sqrt(h)*Z2[:,i]+thteta*(X[:,i]-S[:,i])
    return X,S
T, m, n = 1.0, 1000, 10
X_Euler,S_Euler = Euler(T, m, n)
t_grid = np.linspace(0,T,m+1)

# 展示Euler Scheme的结果
for i in range(n):
    plt.plot(t_grid, X_Euler[i,:])

plt.xlabel('$t$')
plt.ylabel('$X_t$')
plt.title('Euler Scheme for BS')
plt.savefig('Euler_Scheme_BS.jpg')
plt.show()

for i in range(n):
    plt.plot(t_grid, S_Euler[i,:])

plt.xlabel('$t$')
plt.ylabel('$S_t$')
plt.title('Euler Scheme for BS')
plt.savefig('Euler_Scheme_BS.jpg')
plt.show()