今回は、ファーマフレンチのファクターモデルを使ったパフォーマンス要因分析をPythonを使って行う方法を紹介したいと思います。
ファクターモデルを用いたパフォーマンス要因分析の手法と考え空については以前の記事で扱っているので、まだ、そちらを確認していない場合はあわせて確認してみてください。
データの取得
今回は、FamaFrenchの5ファクターモデルを使って分析を行います。
ポートフォリオはTOPIX100の銘柄の等ウェイト、ベンチマークはTOPIX100銘柄の時価加重ウェイトとします。
from pandas_datareader import famafrench as ff
import yfinance as yf
import pandas as pd
import matplotlib.pyplot as plt
#TOPIX銘柄コードの取得
import requests
# CSVファイルのURL
#https://www.jpx.co.jp/markets/indices/topix/の構成銘柄一覧のcsvを取得。状況に応じてリンクは変更
url = "https://www.jpx.co.jp/automation/markets/indices/topix/files/topixweight_j.csv"
# HTTPリクエストを送信してファイルをダウンロード
response = requests.get(url)
# ステータスコードが200(成功)の場合
if response.status_code == 200:
# ダウンロードしたCSVデータをファイルに保存
with open("topixweight_j.csv", "wb") as f:
f.write(response.content)
# ファイルをPandas DataFrameに読み込む
data_df = pd.read_csv("/content/topixweight_j.csv",encoding='shift-jis')
else:
print(f"Failed to download the file. Status code: {response.status_code}")
#TPX100の銘柄を取り出す
topix_100 = list(data_df[(data_df['ニューインデックス区分']=='TOPIX Core30') | (data_df['ニューインデックス区分']=='TOPIX Large70') ]['コード'].astype(int).astype(str))
ticker_list = [i + '.T' for i in topix_100]
date_s = '2000-01-01'
price = yf.download(ticker_list, start=date_s,interval="1mo")
rtn = price['Adj Close'].pct_change().shift(-1) * 100
rtn.index = rtn.index.to_period()
#FamaFrenchデータ
ff_data = ff.FamaFrenchReader('Japan_5_Factors',start=1960).read()[0]
bm_weight = data_df[(data_df['ニューインデックス区分']=='TOPIX Core30') | (data_df['ニューインデックス区分']=='TOPIX Large70') ][['コード','TOPIXに占める個別銘柄のウエイト']]
bm_weight.columns = ['code','weight']
bm_weight['code'] = bm_weight['code'].astype(int).astype(str) + '.T'
bm_weight = bm_weight.set_index('code')
bm_weight['weight'] = bm_weight['weight'].str.replace('%','').astype(float)/100
bm_weight['weight'] = bm_weight['weight'] / bm_weight['weight'].sum()
port_weight = pd.DataFrame(columns=bm_weight.columns, index=bm_weight.index)
port_weight.loc[:,:] = 1/len(port_weight)
rtn = rtn.clip(upper=100, lower=-100)
rtn = rtn/100
ff_data = ff_data/100
weight_port = pd.DataFrame()
for d in ff_data.index:
port_weight_add = port_weight.copy()
port_weight_add['date'] = d
weight_port = pd.concat([weight_port, port_weight_add])
weight_port = weight_port.reset_index().set_index(['date','code'])
weight_bm = pd.DataFrame()
for d in ff_data.index:
bm_weight_add = bm_weight.copy()
bm_weight_add['date'] = d
weight_bm = pd.concat([weight_bm, bm_weight_add])
weight_bm = weight_bm.reset_index().set_index(['date','code'])
リターンは外れ値処理として、-100%~+100%までに調整しています。(1か月あたり)
また、ウェイトに関しては今回はイメージをつかむためのテストデータであるため、データが存在していない期間についてもその銘柄を保有した前提のウェイトになっています。そのあたりは今回は無視しますが、実際は実運用に使ったウェイトを用いることになるので、今回はテストだと思っていただければと思います。ベンチマークの時価総額加重のウェイトも本来ならその時々で変化しているはずですが、今回は、データ取得の観点から最新のものと過去も同じであったとしています。
今回のポートフォリオウェイトやベンチマークウェイトはあくまで計算過程を理解するためのサンプルと思っていただければと思います。
ファクターエクスポージャーの推定
まずは、各銘柄のファクターエクスポージャーを推定する必要があります。
ファクターエクスポージャーの算出方法は大きく2つあります。一つは、財務データなどから直接求める方法、もう一つはリターンデータの回帰によって求める方法です。
以前の記事で詳細を説明しているのでそちらも併せて確認してみてください。
今回は、必要なデータが少なくて済む後者を選択します。
マルチファクターモデルの考えに基づくと、各銘柄のリターンは、ファクターと残差に分離できます。各ファクターリターンに対する係数が、ファクターエクスポージャーになるので、銘柄リターンをファクターリターンで重回帰分析を行い、ファクターエクスポージャーを計算します。
今回は計算期間を36か月のローリングウィンドウとしました。つまり、直近36ヶ月のデータから係数や切片を推定します。
そのコードが以下になります。
import numpy as np
from tqdm import tqdm
k=0
n=36
df_alpha = pd.DataFrame()
fct_exp = pd.DataFrame()
cols = ['Mkt-RF', 'SMB', 'HML', 'RMW', 'CMA']
for k in tqdm(range(len(rtn)-n-2)):
# データを準備します
y = rtn.iloc[k:k+n, :]
date_s = y.index[0]
date_e = y.index[-1]
X = ff_data.loc[date_s:date_e, cols]
# バイアス項を追加します
X_b = np.c_[np.ones((len(X), 1)), X]
# パラメータを計算します (theta_best = (X^T * X)^-1 * X^T * y)
theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
# 新しいデータに対する予測を行います
i = ff_data.index.get_loc(date_e)
X_pred = ff_data.iloc[[i],:].loc[:,cols]
X_pred_b = np.c_[np.ones((len(X_pred), 1)), X_pred]
y_pred = X_pred.values.dot(theta_best[1:,:])
y_true = rtn.iloc[k+n-1,:].values
y_alpha = y_true - y_pred
y_alpha = pd.DataFrame(y_alpha,columns=rtn.columns, index=[ff_data.index[i]])
fct_exp_d = pd.DataFrame(theta_best,columns=y.columns, index=['alpha'] + cols).T
fct_exp_add = pd.DataFrame(fct_exp_d.unstack()).reset_index()
fct_exp_add.columns = ['factor','code','exp']
fct_exp_add['date'] = date_e
fct_exp = pd.concat([fct_exp, fct_exp_add])
df_alpha = pd.concat([df_alpha, y_alpha])
df_exp = pd.pivot(fct_exp, columns='factor',values='exp', index=['date','code'])
この場合、ファクターリターンや、スペシフィックリターン(固有リターン)は次のようになります。
fct_rtn = ff_data[cols]
srtn = pd.DataFrame(df_alpha.unstack()).reset_index()
srtn.columns = ['code','date','srtn']
srtn = srtn.set_index(['date','code'])
リターンデータの作成
ここまでは共通するデータを作成してきましたが、ここからポートフォリオとベンチマークのそれぞれのリターンを計算していきます。
まずは、ベンチマークとポートフォリオそれぞれの
ファクターリターン、スペシフィックリターンを計算しておきます。
def cal_fctrtn(df_weight, df_exp):
d_exp = df_exp.copy()
d_w = pd.pivot(df_weight.reset_index(),columns='code',index='date')
d_fct_exp = pd.DataFrame()
for factor_name in cols:
d_exp = pd.pivot(df_exp.reset_index(),index='date',columns='code',values=factor_name)
add_df = pd.DataFrame((d_w * d_exp).dropna(how='all').sum(axis=1))
add_df.columns = [factor_name]
d_fct_exp = pd.concat([d_fct_exp, add_df],axis=1)
d_fct_exp = d_fct_exp.shift()
d_fct_rtn = d_fct_exp[cols] * fct_rtn[cols]
return d_fct_rtn, d_fct_exp
frtn_port, fexp_port = cal_fctrtn(weight_port, df_exp)
frtn_bm, fexp_bm = cal_fctrtn(weight_bm, df_exp)
def cal_srtn(df_weight, srtn):
d_w = pd.pivot(df_weight.reset_index(),columns='code',index='date',values='weight').shift()
d_s = pd.pivot(srtn.reset_index(), columns='code',index='date',values='srtn')
return pd.DataFrame((d_w * d_s).dropna(how='all').sum(axis=1),columns=['スペシフィック']).sort_index()
srtn_port = cal_srtn(weight_port, srtn)
srtn_bm = cal_srtn(weight_bm, srtn)
最後に各期間のリターンです。
trtn_port = pd.DataFrame(((pd.pivot(weight_port.reset_index(),columns='code',index='date',values='weight').shift() ) * rtn).dropna(how='all').sum(axis=1),columns=['リターン'])
trtn_bm = pd.DataFrame(((pd.pivot(weight_bm.reset_index(),columns='code',index='date',values='weight').shift() ) * rtn).dropna(how='all').sum(axis=1),columns=['リターン'])
ファクター要因の計算
ここまでのデータを整理します。
def make_rtn_data(trtn_, frtn_, srtn_):
frtn_['ファクター'] = frtn_.sum(axis=1)
rtn_data = pd.merge(trtn_, frtn_, right_index=True, left_index=True)
rtn_data = pd.merge(rtn_data, srtn_, right_index=True, left_index=True)
return rtn_data
rtn_data_port = make_rtn_data(trtn_port, frtn_port, srtn_port)
rtn_data_bm = make_rtn_data(trtn_bm, frtn_bm, srtn_bm)
data = pd.merge(rtn_data_port, rtn_data_bm, suffixes=['_port','_bm'], right_index=True, left_index=True).copy()
data['リターン_act'] = data['リターン_port'] - data['リターン_bm']
data['ファクター_act'] = data['ファクター_port'] - data['ファクター_bm']
data['スペシフィック_act'] = data['スペシフィック_port'] - data['スペシフィック_bm']
for col in cols:
data[f'{col}_act'] = data[f'{col}_port'] - data[f'{col}_bm']
これで各要因のリターンを計算できます。
_portはポートフォリオのリターン、_bmはベンチマークのリターン、_actは超過リターンを意味します。
このように分析を行うことで、自分が作成したポートフォリオのリターン源泉がどのような要因にあるのか理解できるようになります。
例えば、HML_actが大きい場合は、バリューファクターで超過リターンが得られているということで、その期間のバリューファクターのファクターリターンがプラスであれば、ポートフォリオの方がバリュー特性が高く、その期間バリューファクターはプラスだったので、ポートフォリオのパフォーマンスにバリュー特性がプラスに寄与したと考えることができます。
今回のコードはこちらのGoogle Colabで確認できます。