import pandas as pd
import numpy as np
import calendar
import matplotlib.pyplot as plt
import matplotlib.dates as mdate
from datetime import datetime
from dateutil.parser import parse
from scipy.optimize import minimize
TOLERANCE = 1e-10
# =============================================================================
# part1 :处理数据 ,通过时间区间[time1,time2]中的数据获取日收益率何周收益率
# =============================================================================
start = '20041231'
end = '20181231'
stock = pd.read_excel('C:\\Users\\Administrator\\Desktop\\Risk parity\\test.xlsx',sheet_name='股票')
bond = pd.read_excel('C:\\Users\\Administrator\\Desktop\\Risk parity\\test.xlsx',sheet_name='债券')
#goods = pd.read_excel('C:\\Users\\Administrator\\Desktop\\Risk parity\\stock&bond.xlsx',sheet_name='南华商品')
Name = [stock , bond] # 所要配置的资产
asset = ['stock' , 'bond'] #记录资产组合中资产的名称,方便后期作图
def _deal_data(Name):
## 首先实现表格的连接
if len(Name) == 1 :
data = Name[0]
else:
data = pd.merge(Name[0], Name[1], how='outer',on=['Date'])
# 两个表格直接连接,多表循环
i = 2
while i < len(Name):
data = pd.merge(data, Name[i], how='outer',on=['Date'])
i = i + 1
data = data.sort_values(by = 'Date')
data = data.reset_index(drop = True)
# 返回经时间排序的资产收盘价数据
## 其次删除空值,选取目标日期
data.dropna(axis=0 ,how ='any' ,inplace =True)
data = data.reset_index(drop = True)
data = data[(data['Date'] >= pd.to_datetime(start)) & (data['Date'] <= pd.to_datetime(end))]
data = data.reset_index(drop = True)
dataAll = pd.date_range(start ,end )
dataAll = pd.DataFrame(dataAll)
dataAll.columns = ['Date']
data = pd.merge(dataAll, data, how='outer',on=['Date'])
data = data.fillna(method='pad')
# 返回包含每个日期的收盘价数据,非交易日附之前最近一个交易日数据
data.dropna(axis=0 ,how ='any' ,inplace =True)
data = data.reset_index(drop = True)
data.set_index('Date',inplace=True)
#直接将指标设为日期以方便后续计算
return data
test = _deal_data(Name)
def _day_return(Name):
# 计算收盘价的日收益
data = _deal_data(Name)
name = data.columns
for j in range(data.shape[1]):
data['day_return_'+name[j]] = (data[name[j]]-data[name[j]].shift(1)) / data[name[j]].shift(1)
# 返回dataFrame,前line列为资产收盘价 ,后line列为标的资产日收益率
# 由于是按日历日期计算的,因而由于节假日存在部分零值,为保证波动率的真实性,删除零值
# data =
# data['day_return_中证全债']=4*data['day_return_中证全债']
line = int(0.5*data.shape[1])
return data.iloc[1::,line::]
def _day_return_leverage(Name,leverage):
# 计算收盘价的日收益
# leverage : 表示杠杆的倍数
data = _deal_data(Name)
name = data.columns
for j in range(data.shape[1]):
data['day_return_'+name[j]] = (data[name[j]]-data[name[j]].shift(1)) / data[name[j]].shift(1)
# 返回dataFrame,前line列为资产收盘价 ,后line列为标的资产日收益率
# 由于是按日历日期计算的,因而由于节假日存在部分零值,为保证波动率的真实性,删除零值
data['day_return_中证全债']=leverage*data['day_return_中证全债']-(leverage-1)*0.03/365
line = int(0.5*data.shape[1])
return data.iloc[1::,line::]
def _period_return(data , period):
# 计算收盘价的月收益率
# 输入: period 对应的是 'W':周 ;'M':月 ;'3M':季度 ;'6M':半年
__df= data.resample(period).last()
dataw = pd.concat([data.iloc[0,0::].T,__df.T],axis=1).T
# 返回周末收盘价,同时保留表格第一行数据
name = dataw.columns
for j in range(dataw.shape[1]):
dataw['return_'+name[j]] = (dataw[name[j]]-dataw[name[j]].shift(1)) / dataw[name[j]].shift(1)
# 返回dataFrame,前line列为资产收盘价 ,后line列为标的资产日收益率
line = int(0.5*dataw.shape[1])
return dataw.iloc[1::,line::]
def _calculate_cov(data):
return data.cov().values
return_daily = _day_return(Name) #日收益率
#return_weekly = _weekly_return(Name) # 周收益率
def MaxDrawdown(asset_value):
# '''最大回撤率'''
# np.maximum.accumulate():生成一列当日之前历史最高价值的序列。
# argmax(f(x))是使得 f(x)取得最大值所对应的变量x
i = np.argmax((np.maximum.accumulate(asset_value) - asset_value) / np.maximum.accumulate(asset_value)) # 结束位置
# if i == 0:
# return 0
j = np.argmax(asset_value[:i]) # 开始位置
return (asset_value[j] - asset_value[i]) / (asset_value[j])
def asset_statistics_model(Name,asset0):
#计算资产的各项指标
asset_price_data = _deal_data(Name)
asset_price_data.columns = asset0
A_return = _period_return(asset_price_data , 'M') #'M'表示当前是一月一调仓
A_return_yearly = _period_return(asset_price_data , 'Y')
# A_return_daily = _period_return(asset_price_data , 'D')
total_return = ((A_return + 1).T.values.cumprod(1)).T[-1]-1
total_return = pd.DataFrame(total_return).T
total_return.columns = asset0
annual_return = pow(total_return + 1, 1.0/13) - 1
annual_return_mean = A_return_yearly.mean()
annual_return_mean = pd.DataFrame(annual_return_mean).T
annual_return_mean.columns = asset0
annual_return_std = A_return_yearly.std()
annual_return_std = pd.DataFrame(annual_return_std).T
annual_return_std.columns = asset0
# annual_Volatility = (A_return_daily.std() / A_return_daily.mean())/ np.sqrt(1/365)
yearly_sharpe = annual_return_mean/annual_return_std
# yearly_sharpe = pd.DataFrame(yearly_sharpe).T
yearly_sharpe.columns = asset0
asset_maxDrawdown = []
for i in range(len(Name)):
asset_maxDrawdown.append(MaxDrawdown(asset_price_data[asset0[i]]))
asset_maxDrawdown = pd.DataFrame(asset_maxDrawdown).T
asset_maxDrawdown.columns = asset0
asset_statistics = pd.concat([total_return , annual_return , annual_return_mean , annual_return_std , yearly_sharpe , asset_maxDrawdown])
asset_statistics = asset_statistics.T
asset_statistics.columns = ['total_return' , 'annual_return' , 'annual_return_mean' , 'annual_return_std' , 'yearly_sharpe' , 'asset_maxDrawdown']
return asset_statistics.T
asset_behaviour_statistics = asset_statistics_model(Name , asset)
# =============================================================================
# part2 : Risk parity 二次规划问题求解
# =============================================================================
#########构建求解二次规划问题的方法
#### step1: 构建二次规划问题并优化
def _allocation_risk(weights, covariances):
# 计算给定资产权重下资产组合的风险
# 输入: weights:资产权重,行向量(注意使用np.matrix()将list转化为1*n为行向量) ; covrances:资产的协方差矩阵
portfolio_risk = np.sqrt((weights * covariances * weights.T))[0, 0]
return portfolio_risk
def _assets_risk_contribution_to_allocation_risk(weights, covariances):
# 计算资产组合中单个资产贡献的风险(total risk contribution)
# 输入:与_allocation_risk()函数一致 ; 输出为 N*1的列向量(N表示资产的数量)
portfolio_risk = _allocation_risk(weights, covariances)
assets_risk_contribution = np.multiply(weights.T,covariances*weights.T)/portfolio_risk
return assets_risk_contribution
def _risk_budget_objective_error(weights, args):
# 二次规划问题的目标函数 :weights为规划问题求解的变量 ;args是参数的传递
# 参数1为资产的协方差阵
covariances = args[0]
# 参数2为单个资产对于资产组合的目标贡献率——1/N
assets_risk_budget = args[1]
# convert the weights to a matrix
weights = np.matrix(weights)
# We calculate the risk of the weights distribution
portfolio_risk = _allocation_risk(weights,covariances)
# We calculate the contribution of each asset to the risk of the weights distribution
assets_risk_contribution = _assets_risk_contribution_to_allocation_risk(weights,covariances)
# We calculate the desired contribution of each asset to the risk of the weights distribution
# asset_risk_target 为行向量
assets_risk_target = np.asmatrix(np.multiply(portfolio_risk, assets_risk_budget))
# Error between the desired contribution and the calculated contribution of each asset
error = sum(np.square(assets_risk_contribution - assets_risk_target.T))[0, 0]
# It returns the calculated error ; the aim is minimize error
# error = sum(np.square(assets_risk_contribution[0]-assets_risk_contribution[1]))
return error
def _get_risk_parity_weights(covariances, assets_risk_budget, initial_weights):
# Restrictions to consider in the optimisation: only long positions whose sum equals 100%
constraints = ({'type': 'eq', 'fun': lambda x: np.sum(x) - 1.0},{'type': 'ineq', 'fun': lambda x: x})
# Optimisation process in scipy
optimize_result = minimize(fun=_risk_budget_objective_error,
x0=initial_weights,
args=[covariances, assets_risk_budget],
method='SLSQP',
constraints=constraints,
tol=TOLERANCE,
options={'disp': False})
# Recover the weights from the optimised object
weights = optimize_result.x
# It returns the optimised weights
return weights
#### step2: 确定portfolio中asset_weight
def get_theday_of_last_year(date):
date_parse = parse(date)
year = date_parse.year - 1
month = date_parse.month
day = date_parse.day
if month==2:
#此时涉及润2月的情况(普通闰年:能被4整除但不能被100整除的年份为普通闰年。;世纪闰年:能被400整除的为世纪闰年。)
if (year+1) % 4 == 0:
if (year+1) % 100 == 0:
day = day - 2
else:
day = day - 1
last_date = datetime(year,month,day)
return last_date
def get_theday_of_next_year(date):
date_parse = parse(date)
year = date_parse.year + 1
month = date_parse.month
day = date_parse.day
if month==2:
#此时涉及润2月的情况(普通闰年:能被4整除但不能被100整除的年份为普通闰年。;世纪闰年:能被400整除的为世纪闰年。)
if (year+1) % 4 == 0:
if (year+1) % 100 == 0:
day = day - 2
else:
day = day - 1
next_date = datetime(year,month,day)
return next_date
def get_nmonth_endday(time,month_number):
"""
获取后n个月最后一天的日期
输入 time 以字符串表示的日期
输入 month_number 为后n个月,可正可负
return: 返回 list 日期
需要import datatime,parse,calendar
"""
today = parse(time)
year = today.year
month = today.month
month+=month_number
while month>12:
month -= 12
year+=1
while month <= 0:
month += 12
year -= 1
firstDayWeekDay, monthRange = calendar.monthrange(year, month)
res=datetime(year,month,monthRange).strftime('%Y%m%d')
return res
def _time_interval_build(number):
# number:时间区间所包含的月份 该函数通常用于生成月末,季度末,半年末的日历日期
# 输入:number(int);输出:date_month(DataFrame)
dateAll = pd.date_range(get_theday_of_next_year(start) , end , freq = 'M')
dateAll = pd.DataFrame(dateAll)
date_month = dateAll.iloc[(dateAll.index%number) ==0]
date_month.columns = ['Date']
date_month = date_month.reset_index(drop = True)
return date_month
# =============================================================================
# def _get_nonzero_value(data):
#
# return
# =============================================================================
# =============================================================================
# def get_weights(_return,time):
#
# # 输入:_return: return on aeest ; time: date of rebalance
#
# data = _return.iloc[(_return.index <= time)&(_return.index> get_nmonth_endday(time))]
#
# covariances = data.cov().values*250
# # We calculate the covariance matrixre
# # The desired contribution of each asset to the portfolio risk: we want all asset to contribute equally
# assets_risk_budget = [1/data.shape[1]]*data.shape[1]
# # Initial weights: equally weighted
# init_weights = [1/data.shape[1]]*data.shape[1]
# # Optimisation process of weights
# weights = _get_risk_parity_weights(covariances,assets_risk_budget, init_weights)#两行一列的向量
# # Convert the weights to a pandas Series
# weights = pd.Series(weights, index=data.columns, name='weight')
#
# return weights
# =============================================================================
def get_weights(_return,time):
# 输入:_return: return on aeest ; time: date of rebalance
data = _return.iloc[(_return.index <= time)&(_return.index> get_theday_of_last_year(time))]
#data反映了回看周期 ,即以过去多长时间的数据计算协方差阵
covariances = data.cov().values*250 #*250表示年化
# We calculate the covariance matrixre
# The desired contribution of each asset to the portfolio risk: we want all asset to contribute equally
assets_risk_budget = [1/data.shape[1]]*data.shape[1]
# Initial weights: equally weighted
init_weights = [1/data.shape[1]]*data.shape[1]
# init_weights = [0.2,0.8]
# Optimisation process of weights
weights = _get_risk_parity_weights(covariances,assets_risk_budget, init_weights)#两行一列的向量
# Convert the weights to a pandas Series
weights = pd.Series(weights, index=data.columns, name='weight')
return weights
# =============================================================================
# #def all_portfolio_weight(number):
# # # 返回的是投资组合权重的字典序
# # result={}
# # date_m = _time_interval_build(number)
# #
# #
# # if number == 6:
# # _return=return_weekly
# # else:
# # _return=return_daily
# # # 如果调仓频率为半年,则按照周收益率计算协方差阵,若调仓频率为单月或季度,则按照日收益率计算协方差阵
# #
# # for i in range(date_m.shape[0]-1):
# ## result[date_m.iloc[i+1].values[0]]=pd.DataFrame(get_weights(return_daily , date_m.iloc[i].values[0] , date_m.iloc[i+1].values[0])).iloc[0,0::]
# # result[date_m.iloc[i+1].values[0]]=pd.DataFrame(get_weights(_return , date_m.iloc[i].values[0] , date_m.iloc[i+1].values[0])).T
# #
# # return result
# =============================================================================
def all_portfolio_weight(_return , number):
# 返回的是投资组合权重的字典序
result=[]
result = pd.DataFrame(result)#生成空的dataFrame以方便记录portfolio中asset的权重
date_m = _time_interval_build(number) # date_m表示的是调仓的频率
result = pd.concat([result ,date_m],axis=1)
result.set_index('Date',inplace=True) #result的index指标对应须调仓日
for i in range(len(_return.columns.values)):
result[_return.columns.values[i]]=0.0
# 通过收益率确定对应资产
for i in range(result.shape[0]):
for j in range(result.shape[1]):
result[result.columns.values[j]][date_m.iloc[i].values[0]] = get_weights(_return ,str( date_m.iloc[i].values[0]) )[j]
return result
#weight_weekly_semiannually = all_portfolio_weight(return_weekly , 6)
weight_daily_monthly = all_portfolio_weight(return_daily , 1)
weight_daily_3monthly = all_portfolio_weight(return_daily , 3)
# =============================================================================
# part3 : backtest
# 回测期:200601-201908
# =============================================================================
backtest_start = '20051231'
backtest_end = '20181231'
def plt_weight_change(weight):
## 功能:采用堆叠柱状图反映资产权重变化
#prepare data
weight.columns=asset
weight_value = [] #记录各资产的权重值,以满足堆叠折线图的绘制过程中,stackplot函数对于纵坐标的要求
for i in range(len(asset)):
weight_value.append(weight[asset[i]])
# Draw
# plt.figure(figsize = (12,6) , dpi = 60 )
colors = ['powderblue','steelblue','darkturquoise','coral','darkcyan','gold','teal','pink','y','r','snow','b','k','g','orange','c','bisque','brown','lime','aqua']
plt.stackplot(weight.index, weight_value, labels=asset, colors=colors)
plt.legend()
plt.show()
return
#plt_weight_change(weight_daily_3monthly)
def _rebalance_portfolio_value(Name):
# 输入:资产组合中各资产的权重和日收益率
# 计算资产组合的单位净值以及累计收益
_return = _day_return(Name)
weight = all_portfolio_weight(_return , 1)
asset_value0 = 1
asset_daily_value = []
asset_return = 1
for i in range(weight.shape[0]-1):
asset_value = asset_value0 *weight.iloc[i]
_return_peiod = _return[(_return.index > weight.index[i]) & (_return.index <= weight.index[i+1])]
for j in range(_return_peiod.shape[0]):
asset_value = np.multiply( ( _return_peiod.iloc[j]+1) , asset_value)
asset_daily_value.append(pd.DataFrame(asset_value).T)
asset_return = asset_return * (np.sum(asset_value)/asset_value0)
asset_value0 = np.sum(asset_value)
asset_daily_value = pd.concat(asset_daily_value)
asset_daily_value.columns = asset
return [asset_daily_value,asset_return]
#portfolio_value = _rebalance_portfolio_value (Name)[0]
#test = _rebalance_portfolio_value(Name)
def _rebalance_portfolio_leverage_value(Name,leverage):
# 输入:资产组合中各资产的权重和日收益率
# 计算资产组合的单位净值以及累计收益
_return = _day_return_leverage(Name,leverage)
weight = all_portfolio_weight(_return , 1)
asset_value0 = 1
asset_daily_value = []
asset_return = 1
for i in range(weight.shape[0]-1):
asset_value = asset_value0 *weight.iloc[i]
_return_peiod = _return[(_return.index > weight.index[i]) & (_return.index <= weight.index[i+1])]
for j in range(_return_peiod.shape[0]):
asset_value = np.multiply( ( _return_peiod.iloc[j]+1) , asset_value)
asset_daily_value.append(pd.DataFrame(asset_value).T)
asset_return = asset_return * (np.sum(asset_value)/asset_value0)
asset_value0 = np.sum(asset_value)
asset_daily_value = pd.concat(asset_daily_value)
asset_daily_value.columns = asset
return [asset_daily_value,asset_return]
def portfolio_statistics(Name):
portfolio_value = _rebalance_portfolio_leverage_value(Name,3)[0]
portfolio_value['Asset_sum'] = portfolio_value.apply(lambda x: x.sum(), axis=1)
portfolio = pd.DataFrame(portfolio_value['Asset_sum'])
date = pd.DataFrame( portfolio.index.values)
portfolio = portfolio.reset_index(drop = True)
portfolio = pd.concat([date,portfolio],axis = 1 )
portfolio.columns = ['Date','资产组合']
portfolio_statistics = asset_statistics_model([portfolio],['资产组合'])
return portfolio_statistics
def plt_asset_value(Name):
## 功能:绘制资产以及资产组合的价格变化
# Prepare data
portfolio_value = _rebalance_portfolio_leverage_value (Name,3)[0]
portfolio_value['Asset_sum'] = portfolio_value.apply(lambda x: x.sum(), axis=1) # 确定资产组合的每日价值
data = _deal_data(Name).iloc[(_deal_data(Name).index >= portfolio_value.index[0])&(_deal_data(Name).index <= portfolio_value.index[-1])]
weight = pd.DataFrame(1/data.iloc[0]).T.values
data = np.multiply(weight,data) # 归一化处理
data.columns = asset #data反映的是资产组合中个资产的价值
# Plot
# =============================================================================
# fig,ax1 = plt.subplots(figsize = (10,8))
# ax2 = ax1.twinx() # mirror the ax1
# # 由于股票的变动幅度较大,为更好的反映各个资产价格的变化,将股票价格单独设主坐标轴表示
# ax1.plot(data.index, data[asset[0]] ,'g',label = asset[0])
# ax2.plot(data.index, data[asset[1]],'b')
# ax2.plot(portfolio_value.index, portfolio_value['Asset_sum'],'r')
# plt.show()
# =============================================================================
for i in range(data.shape[1]):
plt.plot(data.index , data[asset[i]] ,label=asset[i])
plt.plot(portfolio_value.index , portfolio_value['Asset_sum'] ,label = 'Portfolio_value')
# Decorations
plt.xlabel('Date')
plt.ylabel('price')
plt.legend()
plt.show()
return
plt_asset_value(Name)
def risk_contribution(Name):
# 根据二次规划的解:计算各项资产对于资产组合的实际贡献程度
# portfolio_value = _rebalance_portfolio_value (Name) # 资产组合中各资产的净值
# _return = _day_return_leverage(Name,3) #资产组合中各资产的日收益率
_return = _day_return_leverage(Name,3)
weight = all_portfolio_weight(_return , 1) #资产组合中各资产的权重
contribution = []
for i in range(weight.shape[0]-1):
_return_peiod = _return[(_return.index > weight.index[i]) & (_return.index <= weight.index[i+1])]
weights = np.matrix(pd.DataFrame(weight.iloc[i]).T)
covariances = np.matrix(_return_peiod.cov())
portfolio_risk = _allocation_risk(weights, covariances)
asset_risk = _assets_risk_contribution_to_allocation_risk(weights, covariances)
asset_contribution = pd.DataFrame(asset_risk/portfolio_risk).T
# asset_contribution= pd.concat( weight.index[i] , asset_contribution)
contribution.append(asset_contribution)
contribution = pd.concat(contribution)
contribution.columns = asset
contribution.index = pd.date_range(weight.index[1] , weight.index[-1] ,freq = 'M')
return contribution
asset_actual_risk_contribution = risk_contribution(Name)