利用不同的排队论模型得到limit order book的平均队列长度QL和平均等待时间DLC，与收益率进行相关性分析。
模型假设
使用买一价和卖一价的平均值作为 reference price，用一段时间间隔下reference price 的变化率作为收益率。利用use.total_amount和use.net_amount可以得到买方和卖方的成交单量，分别作为队列的离开与到达。
使用数据
文件/data/intern/3_weeks.pq，每0.25s一个数据点。后面主要用到中频和高频两种取样方法。
中频：每10s的volume求和作为一个点。对于每个点，使用其前30分钟的数据作为队列计算，与后15min的收益率进行相关性检验。
高频：每1s的volume求和作为一个点。对每个点，使用其前60s的数据作为队列计算，与后15s的收益率进行相关性检验。
模型总结
M/M/1模型
假定到达与离开都服从poisson分布，参数记为\lambda, \mu，则有缓冲效用\rho=\frac{\lambda}{\mu}
QL=\frac{\rho}{1-\rho} , DLC=\frac{\rho}{\mu-\lambda}
M/M/1/N
在上模型基础上限定最大队长N
DLC=\rho*\frac{1-(N+1)\rho^N+N\rho^{N+1}}{(1-\rho)(1-\rho^{N+1})}
M/G/1模型
假定离开是独立同分布，其他同M/M/1
QL=\frac{\lambda^2 \overline {X^2}}{2(1-\rho)} , DLC=\frac{\lambda \overline {X^2}}{2(1-\rho)}
G/G/1模型
到达与离开都是独立同分布，没有直接的计算公式，需要先算出队列的概率密度函数再计算。具体到本limit order book的问题，概率分布的公式如下：
\pi_i(n)=\pi_i(0)\prod_{j=1}^n \rho_i(j-1)
\pi_i(0)=(1+\sum_{n=1}^\infty \prod_{j=1}^n\rho_i(j-1))^{-1}
则一般情况下可得：
DLY_i=1-\pi_i(0)
QL_i=\sum_{k=1}^\infty k\pi_i(1+k)
运算代码
M/M/1模型
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import math
from scipy.stats import pearsonr
data = pd.read_parquet('/data/intern/3_weeks.pq')
Pref=(data['ask1']+data['bid1'])/2
data["Pref"]=Pref
use=data.loc[:,['net_amount','total_amount','Pref']]
use['decrease_amount']=(use.total_amount+use.net_amount)/2#买方成交单量 离开
use['increase_amount']=(use.total_amount-use.net_amount)/2#卖方成交单量 到达
use=use.resample('10s').sum()#以10s的数据计算
lambda_=use.increase_amount.rolling('30T').mean()
mu_=use.decrease_amount.rolling('30T').mean()
rate=np.ones(len(use))
for i in range(len(use)):
if i < 60:
rate[i]=(use.Pref[i]-use.Pref[0])/use.Pref[0]
else:
rate[i]=(use.Pref[i-60]+use.Pref[i])/use.Pref[i-60]
use["rate"]=rate
#清除死亡率小于出生率的点
lambda0=[]
mu0=[]
rate0=[]
for i in range(len(mu_)):
if mu_[i]>lambda_[i]:
lambda0.append(lambda_[i])
mu0.append(mu_[i])
rate0.append(rate[i])
rho=np.array(lambda0)/np.array(mu0)
lambda_=np.array(lambda0)
mu_=np.array(mu0)
rate=np.array(rate0)
#平均等待时间W,平均队长N
W=rho/(mu_-lambda_)
N=lambda_*W
N_use=N[180:-60]
rate_use=rate[240:]
pearsonr(N_use,rate_use)
M/G/1模型
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import math
from scipy.stats import pearsonr
data = pd.read_parquet('/data/intern/3_weeks.pq')
Pref=(data['ask1']+data['bid1'])/2
data["Pref"]=Pref
use=data.loc[:,['net_amount','total_amount','Pref']]
use['decrease_amount']=(use.total_amount+use.net_amount)/2#买方成交单量 离开
use['increase_amount']=(use.total_amount-use.net_amount)/2#卖方成交单量 到达
use=use.resample('10s').sum()#以10s的数据计算
lambda_=use.increase_amount.rolling('30T').mean()
mu_=use.decrease_amount.rolling('30T').mean()
rate=np.ones(len(use))
#计算10分钟平均收益率
for i in range(len(use)):
if i < 90:
rate[i]=(use.Pref[i]-use.Pref[0])/use.Pref[0]
else:
rate[i]=(use.Pref[i-90]+use.Pref[i])/use.Pref[i-90]
use["rate"]=rate
key=use.decrease_amount*use.decrease_amount
key=key.rolling('30T').mean()
#清除死亡率小于出生率的点
lambda0=[]
mu0=[]
rate0=[]
for i in range(len(mu_)):
if mu_[i]>lambda_[i]:
lambda0.append(lambda_[i])
mu0.append(mu_[i])
rate0.append(rate[i])
rho=np.array(lambda0)/np.array(mu0)
lambda_=np.array(lambda0)
mu_=np.array(mu0)
rate=np.array(rate0)
#平均等待时间W,平均队长N
W=(lambda_*key)/(2*(1-rho))
N=lambda_*W
N_use=N[180:-60]
rate_use=rate[240:]
pearsonr(N_use,rate_use)
G/G/1模型
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import pearsonr
data = pd.read_parquet('/data/intern/3_weeks_1s.pq')
use=data.loc[:,['net_amount','total_amount']]
use['decrease_amount']=(use.total_amount+use.net_amount)/2#买方成交单量 离开
use['increase_amount']=(use.total_amount-use.net_amount)/2#卖方成交单量 到达
use=use.resample('10s').sum()
use=use[:500]#确定队列长度
Pref=(data['ask1']+data['bid1'])/2
use["Pref"]=Pref
#计算队列
use['total_size']=np.add.accumulate(use.increase_amount)-np.add.accumulate(use.decrease_amount)
#用clip函数处理极值
use.decrease_amount=np.clip(use.decrease_amount.values,0.001,100000)
use.increase_amount=np.clip(use.increase_amount.values,0.001,100000)
use.total_size=use.total_size/4000#size减小倍，这样size的概率密度函数有意义，不然都是0或1
use.total_size=use.total_size.astype(int)#模型需要用整数计算
t=180#30min*6数据点个数
#frho函数，用于确定第i个点的rho（n）的值,
def frho(i,n):
queue=use[i-t:i]
sincrease=0.00001
sdecrease=0.00001
for j in range(len(queue)-1):
if queue.total_size[j] == n:
sincrease+=queue.increase_amount[j]
sdecrease+=queue.decrease_amount[j+1]
return sincrease/sdecrease
use['num']=range(len(use))
times=40#无穷次求和，设置一个大数,观察知rho_times=0，则设为t足够
def fpi(i,n):
a_use1=[]
onlyusethere=[]
for j in range(times):
x=frho(i,j+1)
onlyusethere.append(x)
a_use1=np.add.accumulate(onlyusethere)
t_use=np.multiply.accumulate(a_use1)
pi_0=1/(1+t_use.sum())
if n==0:
return pi_0
else:
return pi_0*t_use[n-1]
QL=[]
for i in range(len(use))[t:len(use)]:
sum=0
for j in range(times):
sum+=j*fpi(i,j)
QL.append(sum)
print(len(QL))
QL
#将QL写入外部文件，便于并行时调用
output = open('1.pk', 'wb')
pickle.dump(QL, output)
output.close()
rate=np.ones(len(use))
#计算10分钟平均收益率
for i in range(len(use)):
if i < 90:
rate[i]=(use.Pref[i]-use.Pref[0])/use.Pref[0]
else:
rate[i]=(use.Pref[i-90]+use.Pref[i])/use.Pref[i-90]
use["rate"]=rate
N_use=QL[:-60]
rate_use=use.rate[240:]
pearsonr(N_use,rate_use)
结果分析
M/M/1
高频：
pearson相关系数：-8.09*10^{-5}
中频：
pearson相关系数：0.0104
M/M/1/N
N 取50时
高频：
pearson相关系数：-0.04734
中频：
pearson相关系数：-0.0378
N 取10时
高频：
pearson相关系数：-0.0059
中频：
pearson相关系数：0.04312
M/G/1
高频：
pearson相关系数：8.66*10^{-4}
中频：
pearson相关系数：0.0015
G/G/1
考虑到中频的相关性明显高于高频情况，只列出中频的结果，并用过去5min收益率、10min收益率、15min收益率作为信号，同未来15min收益率计算相关性进行对比。
计算1.5h的数据（540个数据点）
pearson相关系数：0.44798
前5min收益率与后15min收益率相关性：-0.3894
前10min收益率与后15min收益率相关性：-0.4175
前15min收益率与后15min收益率相关性：-0.5029
计算12小时的数据（4000个数据点）
pearson相关系数：0.08827
前5min收益率与后15min收益率相关性：-0.4396
前10min收益率与后15min收益率相关性：-0.4275
前15min收益率与后15min收益率相关性：-0.4147
计算24小时的数据（8000个数据点）
pearson相关系数：0.05169
前5min收益率与后15min收益率相关性：-0.4582
前10min收益率与后15min收益率相关性：-0.4570
前15min收益率与后15min收益率相关性：-0.4515
数据量变大时，前5min,10min,15min收益率和后15min收益率的相关性系数基本不变，维持较高的值，而GG1模型得出的相关系数持续减小
