【笔记】多因子系列报告之一:因子测试框架 (光大)
最近看了一下【光大证券-光大证券多因子系列报告之一:因子测试框架】非常想实现一下,就写了这个笔记记录一下。本篇文章最主要的就是介绍因子有效性的测试方法。本文的框架如下文所示
- 单因子测试具体步骤
1、因子筛选
2、单因子回归模型
3、单因子的有效性检验
4、分层回溯测试
单因子测试具体步骤
样本筛选
- 1、测试样本范围:沪深300指数成分股
- 2、测试样本期:2015-01-01至2017-12-31
- 3、筛选规则:剔除选股日的ST/PT 股票;
- 4、数据清洗:采取去极值法
绝对中位数:
我们将大于的值或小于的值定义为异常值。在对异常值处理时,我们选择将其设为上下限的数值。- 5、因子标准化:Z值标准化
```python获取沪深300指数成分股在测试周期内的因子值
from jqfactor import Factor, calc_factors
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
stock = get_index_stocks('000300.XSHG')
class Hs300Alpha(Factor):
# 设置因子名称
name = 'hs300_alpha'
# 设置获取数据的时间窗口长度
max_window = 10
# 设置依赖的数据
dependencies = ['close']
# 计算因子的函数, 需要返回一个 pandas.Series, index 是股票代码,value 是因子值
def calc(self, data):
# 获取个股的收盘价数据
close = data['close']
# 计算个股近10日收益
stock_return = close.iloc[-1,:]/close.iloc[0,:] -1
# 获取指数(沪深300)的收盘价数据
index_close = self._get_extra_data(securities=['000300.XSHG'], fields=['close'])['close']
# 计算指数的近10日收益
index_return = index_close.iat[-1,0]/index_close.iat[0,0] - 1
# 计算 alpha
alpha = stock_return - index_return
return alpha
factors = calc_factors(stock, [Hs300Alpha()], start_date='2015-01-01', end_date='2017-12-31')
```python
#处理缺失值
data=factors['hs300_alpha']
output=pd.DataFrame()
for i in range(300):
p=sum(data.iloc[:,i].isnull())/len(data.iloc[:,i])
if p<0.2:
data.iloc[:,i].fillna(mean(data.iloc[:,i]))
output[i]=data.iloc[:,i]
output.columns=data.columns[output.columns]
output=output.fillna(mean(data))
#异常值处理
for i in range(len(output.columns)):
MAD=median(abs(output.iloc[:,i]-median(output.iloc[:,i])))
MAX=median(output.iloc[:,i]) 3*1.4826*MAD
MIN=median(output.iloc[:,i])-3*1.4826*MAD
output.iloc[:,i][output.iloc[:,1]>MAX]=MAX
output.iloc[:,i][output.iloc[:,1]<MIN]=MIN
#标准化
for i in range(len(output.columns)):
output.iloc[:,i]=(output.iloc[:,i]-mean(output.iloc[:,i]))/std(output.iloc[:,i])
单因子回归模型
截面回归(Cross-Section Regression)是目前业界较常用于因子测试的方法。我们选择每期针对全体样本做一次回归,回归时因子暴露为已知变量,回归得到每期的一个因子收益值在通过多期回归后我们就可以得到因子值的序列,也就是因子收益率序列,同时可以得到t值序列,也就是因子值与股票收益率相关性的t检验得到的t值。t值序列将有助于我们挑选有效因子,后文中会详细解释t 值的使用方法。
进行截面回归判断每个单因子的收益情况和显著性时,需要特别关注A 股市场中一些显著影响个股收益率的因素,例如行业因素和市值因素。市值因子在过去的很长一段时间内都是A股市场上影响股票收益显著性极高的一个因子,为了能够在单因子测试时得到因子真正收益情况,我们在回归测试时对市值因子也做了剔除。也就是行业市值中性化。
- 获取行业哑变量和市值因子值
- 对测试因子剔除行业和市值因素(用中性化的方法)
- 对剔除后的因子进行RLM回归(Robust Linear Model)
- 得到每一期的因子收益率
#获得行业哑变量矩阵
from jqdata import *
sw=get_industries(name='sw_l1').index
industry=pd.DataFrame(0,columns=output.columns,index=range(0,28))
for i in range(len(sw)):
temp=list(set(output.columns).intersection(set(get_industry_stocks(sw[i]))))
industry.loc[i,temp]=1
#去除市值、行业因素,得到新的因子值
newx=pd.DataFrame()
for i in range(len(output.index)):
m= get_fundamentals(query(valuation.circulating_cap,valuation.code).filter(valuation.code.in_(output.columns)), date=output.index[i])
m.index=np.array(m['code'])
m=m.iloc[:,0]
m=(m-mean(m))/std(m)
x=output.iloc[i,:]
conc=pd.concat([x,m,industry.T],axis=1).fillna(mean(m))
est=sm.OLS(conc.iloc[:,0],conc.iloc[:,1:]).fit()
y_fitted = est.fittedvalues
newx[i]=est.resid
newx=newx.T
newx.index=output.index
newx=newx.iloc[1:,:]
#看图
'''
fig, ax = plt.subplots(figsize=(8,6))
ax.plot(conc.iloc[:,1],conc.iloc[:,0], 'o', label='data')
ax.plot(conc.iloc[:,1], y_fitted, 'r--.',label='OLS')
'''
#将因子值和收益率匹配
output=output.iloc[:-1,:]
df = get_price(list(output.columns), start_date='2017-01-01', end_date='2017-12-31', frequency='daily', fields=['close'])
y=df['close'].diff()/np.array(df['close'])
y=y.iloc[1:,:]
y=y.fillna(mean(y))
y = y.drop((y.index).difference(newx.index))
#做回归 求回归系数也就是因子收益率
f=[0]*len(y.index)
t=[0]*len(y.index)
for i in range(len(y.index)):
rlm_model = sm.RLM(y.iloc[i,:], newx.iloc[i,:], M=sm.robust.norms.HuberT()).fit()
f[i]=float(rlm_model.params)
t[i]=float(rlm_model.tvalues)
'''
#对回归的结果画图
y_fit=rlm_model.fittedvalues
fig, ax = plt.subplots(figsize=(8,6))
ax.plot(newx,y, 'o', label='data')
ax.plot(newx, y_fit, 'r--.',label='OLS')
'''
单因子的有效性检验
因子的有效性及稳定性:
1、因子收益序列的假设检验t值
2、因子收益序列大于0的概率
3、t值绝对值的均值
4、t值绝对值大于等于2的概率因子的有效性及预测能力:
1、IC值的均值
2、IC值的标准差
3、IC大于0的比例
4、IC绝对值大于0.02 的比例
5、IR(IR=IC均值/IC标准差)
#检验结果 #因子收益序列>0的概率 sum(pd.Series(f)>0)/len(f) #t值绝对值的均值---回归假设检验的t值 mean(abs(pd.Series(t))) #t值绝对值大于等于2的概率---回归假设检验的t值 sum(abs(pd.Series(t))>2)/len(t) #计算IC值序列 IC=[0]*len(y.columns) for i in range(len(y.columns)): IC[i]=corrcoef(pd.Series(f).rank(),y.iloc[:,i].rank())[0,1] #计算IC值的均值 mean(IC) #计算IC值的标准差 std(IC) #IC大于0的比例 sum(pd.Series(IC)>0)/len(IC) #IC绝对值大于0.02的比例 sum(pd.Series(IC)>0.02)/len(IC) #IR值 mean(IC)/std(IC)
指标名称 | Hs300Alpha因子 |
---|---|
因子收益序列>0的概率 | 0.45 |
t值绝对值的均值 | 0.99 |
t值绝对值大于等于2的概率 | 0.1 |
IC值的均值 | -0.003 |
计算IC值的标准差 | 0.068 |
IC大于0的比例 | 0.47 |
IC绝对值大于0.02的比例 | 0.38 |
IR值 | -0.05 |
由上表可知Hs300Alpha因子收益序列明没有显著大于0,因子收益状况不佳,同时IC值的均值也只有-0.003,大于0的比例为0.47,可见因子的有效方向,是一个无效因子。(因为我做测试的时候随意挑选了聚宽因子库中的一个因子,所以并没有测试到有效的因子,也懒得再换了)
#图
#因子收益时间序列图
plt.bar(range(len(f)),f)
#因子收益分布直方图
plt.hist(f)
#回归因子收益t值绝对值
plt.bar(range(len(t)),abs(pd.Series(t)))
#因子IC值序列
plt.bar(range(len(IC)),pd.Series(IC))
分层回溯测试
由于单因子回归法所得到的因子收益值序列并不能直观的反应因子在各期的历史收益情况以及单调性,为了同时能够展示所检验因子的单调性,我们将通过分层打分回溯的方法作为补充。
根据因子值的大小,把股票分为5组,查看每组之间的累积收益率,如果区分度大证明该因子效果较为明显。
#做分层回溯法
newy=df['close']/df['close'].iloc[0,:]-1
newy=newy.iloc[1:,:]
newy=newy.fillna(mean(y))
newy = newy.drop((newy.index).difference(newx.index))
fc1=[0]*len(newx.index)
fc2=[0]*len(newx.index)
fc3=[0]*len(newx.index)
fc4=[0]*len(newx.index)
fc5=[0]*len(newx.index)
for i in range(len(newx.index)):
d1=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5)].index
d2=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*2) ].index
d2=(d2).difference(d1)
d3=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*3) ].index
d3=(d3).difference(d2)
d4=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*4) ].index
d4=(d4).difference(d3)
d5=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*5) ].index
d5=(d5).difference(d4)
mean1=mean(newy.iloc[i,:][d1])
mean2=mean(newy.iloc[i,:][d2])
mean3=mean(newy.iloc[i,:][d3])
mean4=mean(newy.iloc[i,:][d4])
mean5=mean(newy.iloc[i,:][d5])
fc1[i]=mean1
fc2[i]=mean2
fc3[i]=mean3
fc4[i]=mean4
fc5[i]=mean5
fcsum=pd.DataFrame(fc1)
fcsum[1]=fc2
fcsum[2]=fc3
fcsum[3]=fc4
fcsum[4]=fc5
#画图
plt.plot(fcsum)
由图中可见,该因子的区分度并不明显,没有明显的区分,所以该因子在测试期内并没有显著的收益。
from jqfactor import Factor, calc_factors
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import statsmodels.api as sm
stock = get_index_stocks('000300.XSHG')
class Hs300Alpha(Factor):
# 设置因子名称
name = 'hs300_alpha'
# 设置获取数据的时间窗口长度
max_window = 10
# 设置依赖的数据
dependencies = ['close']
# 计算因子的函数, 需要返回一个 pandas.Series, index 是股票代码,value 是因子值
def calc(self, data):
# 获取个股的收盘价数据
close = data['close']
# 计算个股近10日收益
stock_return = close.iloc[-1,:]/close.iloc[0,:] -1
# 获取指数(沪深300)的收盘价数据
index_close = self._get_extra_data(securities=['000300.XSHG'], fields=['close'])['close']
# 计算指数的近10日收益
index_return = index_close.iat[-1,0]/index_close.iat[0,0] - 1
# 计算 alpha
alpha = stock_return - index_return
return alpha
factors = calc_factors(stock, [Hs300Alpha()], start_date='2015-01-01', end_date='2017-12-31')
data=factors['hs300_alpha']
#处理缺失值
output=pd.DataFrame()
for i in range(300):
p=sum(data.iloc[:,i].isnull())/len(data.iloc[:,i])
if p<0.2:
data.iloc[:,i].fillna(mean(data.iloc[:,i]))
output[i]=data.iloc[:,i]
output.columns=data.columns[output.columns]
output=output.fillna(mean(data))
#异常值处理
for i in range(len(output.columns)):
MAD=median(abs(output.iloc[:,i]-median(output.iloc[:,i])))
MAX=median(output.iloc[:,i])+3*1.4826*MAD
MIN=median(output.iloc[:,i])-3*1.4826*MAD
output.iloc[:,i][output.iloc[:,1]>MAX]=MAX
output.iloc[:,i][output.iloc[:,1]<MIN]=MIN
#标准化
for i in range(len(output.columns)):
output.iloc[:,i]=(output.iloc[:,i]-mean(output.iloc[:,i]))/std(output.iloc[:,i])
#获得行业哑变量矩阵
from jqdata import *
sw=get_industries(name='sw_l1').index
industry=pd.DataFrame(0,columns=output.columns,index=range(0,28))
for i in range(len(sw)):
temp=list(set(output.columns).intersection(set(get_industry_stocks(sw[i]))))
industry.loc[i,temp]=1
#去除市值、行业因素,得到新的因子值
newx=pd.DataFrame()
for i in range(len(output.index)):
m= get_fundamentals(query(valuation.circulating_cap,valuation.code).filter(valuation.code.in_(output.columns)), date=output.index[i])
m.index=np.array(m['code'])
m=m.iloc[:,0]
m=(m-mean(m))/std(m)
x=output.iloc[i,:]
conc=pd.concat([x,m,industry.T],axis=1).fillna(mean(m))
est=sm.OLS(conc.iloc[:,0],conc.iloc[:,1:]).fit()
y_fitted = est.fittedvalues
newx[i]=est.resid
newx=newx.T
newx.index=output.index
newx=newx.iloc[1:,:]
#看图
'''
fig, ax = plt.subplots(figsize=(8,6))
ax.plot(conc.iloc[:,1],conc.iloc[:,0], 'o', label='data')
ax.plot(conc.iloc[:,1], y_fitted, 'r--.',label='OLS')
'''
#将因子值和y值匹配
output=output.iloc[:-1,:]
df = get_price(list(output.columns), start_date='2017-01-01', end_date='2017-12-31', frequency='daily', fields=['close'])
y=df['close'].diff()/np.array(df['close'])
y=y.iloc[1:,:]
y=y.fillna(mean(y))
y = y.drop((y.index).difference(newx.index))
#做回归 求回归系数
f=[0]*len(y.index)
t=[0]*len(y.index)
for i in range(len(y.index)):
rlm_model = sm.RLM(y.iloc[i,:], newx.iloc[i,:], M=sm.robust.norms.HuberT()).fit()
f[i]=float(rlm_model.params)
t[i]=float(rlm_model.tvalues)
'''
#对回归的结果画图
y_fit=rlm_model.fittedvalues
fig, ax = plt.subplots(figsize=(8,6))
ax.plot(newx,y, 'o', label='data')
ax.plot(newx, y_fit, 'r--.',label='OLS')
'''
#检验结果
#因子收益序列>0的概率
sum(pd.Series(f)>0)/len(f)
#t值绝对值的均值---回归假设检验的t值
mean(abs(pd.Series(t)))
#t值绝对值大于等于2的概率---回归假设检验的t值
sum(abs(pd.Series(t))>2)/len(t)
#计算IC值序列
IC=[0]*len(y.columns)
for i in range(len(y.columns)):
IC[i]=corrcoef(pd.Series(f).rank(),y.iloc[:,i].rank())[0,1]
#计算IC值的均值
mean(IC)
#计算IC值的标准差
std(IC)
#IC大于0的比例
sum(pd.Series(IC)>0)/len(IC)
#IC绝对值大于0.02的比例
sum(pd.Series(IC)>0.02)/len(IC)
#IR值
mean(IC)/std(IC)
#因子收益时间序列图
plt.bar(range(len(f)),f)
#因子收益分布直方图
plt.hist(f)
#回归因子收益t值绝对值
plt.bar(range(len(t)),abs(pd.Series(t)))
#因子IC值序列
plt.bar(range(len(IC)),pd.Series(IC))
#做分层回溯法
newy=df['close']/df['close'].iloc[0,:]-1
newy=newy.iloc[1:,:]
newy=newy.fillna(mean(y))
newy = newy.drop((newy.index).difference(newx.index))
fc1=[0]*len(newx.index)
fc2=[0]*len(newx.index)
fc3=[0]*len(newx.index)
fc4=[0]*len(newx.index)
fc5=[0]*len(newx.index)
for i in range(len(newx.index)):
d1=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5)].index
d2=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*2) ].index
d2=(d2).difference(d1)
d3=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*3) ].index
d3=(d3).difference(d2)
d4=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*4) ].index
d4=(d4).difference(d3)
d5=newx.iloc[i,:][newx.iloc[i,:].rank()<=(len(newx.index)/5*5) ].index
d5=(d5).difference(d4)
mean1=mean(newy.iloc[i,:][d1])
mean2=mean(newy.iloc[i,:][d2])
mean3=mean(newy.iloc[i,:][d3])
mean4=mean(newy.iloc[i,:][d4])
mean5=mean(newy.iloc[i,:][d5])
fc1[i]=mean1
fc2[i]=mean2
fc3[i]=mean3
fc4[i]=mean4
fc5[i]=mean5
fcsum=pd.DataFrame(fc1)
fcsum[1]=fc2
fcsum[2]=fc3
fcsum[3]=fc4
fcsum[4]=fc5
#画图
plt.plot(fcsum)