$\rm{[Python]}$$\rm{[Tableau]}$ 台股 $\rm{Fama\;French\; 5 \;Factor}$¶
$\rm{周育呈}$ e11559560532@gmail.com¶
前言¶
根據過去使用 R 語言的經驗，我還是覺得從做中學是熟悉一個工具最快的方式，所以想找個計算過程不過於簡單且未來用的到的東西來練習 Python 跟 Tableau，決定來計算台股的 Fama Factor，之前在寫碩論時有用 R 計算過 Fama Factor 的經驗，順便可以利用過去經驗熟悉與 dplyr 對應的 pandas 方法。
所使用的公司資料皆來自於 FinMind，包括調整股價、公司市值、資產負債表、綜合損益表，原本調整股價打算用 yfinance 抓來的資料，但在算報酬時發現有些股票的交易日沒資料，加上其他資料都從 FinMind 上抓，所以最後決定所有資料都使用 FinMind 資料庫的。
FinMind 是真的方便，但資料品質沒有到很高，不論是資產負債表、綜合損益表和調整股價，財務報表內提供的資料沒有很多，尤其是調整股價只到小數點第2位，在計算報酬率的時候跟真實資料會有些差距，但臨時需要股市資料的時候 FinMind 真的迅速方便，感謝 FinMind 提供開源資料庫。
因子計算方式¶
計算方式皆採 Kenneth R. French 網站的方法，即 Fama, E. F., & French, K. R. (2015). A five-factor asset pricing model. Journal of Financial Economics, 116(1), 1–22. 中的內容，基於資料的侷限性計算方式稍微不同，以下分別說明：
- MKTPR (Market Premium)
加權指數報酬率減去無風險利率，無風險利率使用臺灣十年期國債債券收益率，無風險利率資料來源為 Investing.com。 - HML (High Minus Low)
依照樣本公司的市值中位數分為2組 (s1,s2)，並依照帳面市值比 (Book to Market ratio) 的30與70百分位數分為3組 (bm1,bm2,bm3) ，共2x3=6組投資組合，以大公司價值股減去大公司成長股的報酬率計算溢酬 (premuim)，小公司同上，兩溢酬平均數即為HML，即： $$\tt{HML=\frac{(s2bm3 - s2bm1)+(s1bm3 - s1bm1)}{2}}$$ - RMW (Robust Minus Weak)
依照樣本公司的市值中位數分為2組 (s1,s2)，並依照營業獲利 (Operating Profitability) 的30與70百分位數分為3組 (op1,op2,op3)，共2x3=6組投資組合，以大公司高獲利減去大公司低獲利的報酬率計算溢酬 (premuim)，小公司同上，兩溢酬平均數即為RMW，即： $$\tt{RMW=\frac{(s2op3 - s2op1)+(s1op3 - s1op1)}{2}}$$ - CMA (Conservative Minus Aggressive)
依照樣本公司的市值中位數分為2組 (s1,s2)，並依照總資產成長率 (Investment) 的30與70百分位數分為3組 (inv1,inv2,inv3)，共2x3=6組投資組合，以大公司低投資減去大公司高投資的報酬率計算溢酬 (premuim)，小公司同上，兩溢酬平均數即為CMA，即： $$\tt{CMA=\frac{(s2inv1 - s2inv3)+(s1inv1 - s1inv3)}{2}}$$ - SMB (Small Minus Big)
SMB分為3因子模型 $\tt SMB_3$ 與5因子模型 $\tt SMB_5$，以各因子的小公司減去大公司的報酬率計算溢酬(premuim)，所有溢酬平均數即為SMB，即： $$\tt{SMB_{(BM)}=\frac{(s1bm1-s2bm1)+(s1bm2-s2bm2)+(s1bm3-s2bm3)}{3}} $$ $$\tt{SMB_{(OP)}=\frac{(s1op1-s2op1)+(s1op2-s2op2)+(s1op3-s2op3)}{3}} $$ $$\tt{SMB_{(INV)}=\frac{(s1inv1-s2inv1)+(s1inv2-s2inv2)+(s1inv3-s2inv3)}{3}} $$ 其中 $\tt SMB_3$ 與 $\tt SMB_5$ 分別為： $$\tt{SMB_{3}=SMB_{(BM)}} $$ $$\tt{SMB_{5}=\frac{SMB_{(BM)}+SMB_{(OP)}+SMB_{(INV)}}{3}} $$
補充
其中營業獲利 (Operating Profitability) 依照文獻定義為營收 (revenues) 減去銷貨成本 (cost of goods sold) 減去營業費用 (selling, general, and administrative expenses) 再減去利息支出 (interest expense) 最後除以股東權益 (book equity)，而 FinMind 的綜合損益表沒有包含完整的利息支出且每期損益表資料有的項目不完全相同(詳見附錄)，故在計算營業獲利時忽略利息支出。
投資組合分組方式與文獻相同，以每年7月重新分組，分組依據以去年資料為基準，公司市值採去年年底市值，其餘公司資料採去年年報，以 2024 年 7 月分組為例，依照 2023 年年底公司市值與 2023 年年財報作為分組依據，該分組組別持續到 2025 年 6 月，2025 年 7 月再以 2024 年資料重新分組，即分組頻率為每年 1 次。
所有因子皆以獨立分組，並採市值加權 (Value-Weighted) 計算投資組合報酬率，樣本為台股上市上櫃公司並剔除金融業，樣本期間為 2014 年 7 月至 2025 年 3 月 (因資產負債表資料始於 2013 年)。
除了上述的 Fama 五因子，其他因子中最知名的就是動能因子，忘記在哪看到有人說除非有特殊需求，否則業界檢驗投組是否有 alpha 基本上是用 Carhart 四因子模型，不知道真的假的總之順便計算動能因子。
- UMD (Up Minus Down) 依照樣本公司的市值中位數分為 2 組 (s1,s2)，並依照報酬動能 (Momentum) 的 30 與 70 百分位數分為 3 組 (mom1,mom2,mom3)，共 2x3=6 組投資組合，以大公司高動能減去大公司低動能的報酬率計算溢酬 (premuim)，小公司同上，兩溢酬平均數即為 UMD，即： $$\tt{UMD=\frac{(s2mom3 - s2mom1)+(s1mom3 - s1mom1)}{2}}$$
補充
報酬動能 (Momentum) 的定義為：過去12個月到過去2個月的累積報酬，不同於上述因子，動能因子採每月重新分組，且動能因子的市值分組 (s1,s2) 也是每月分組，採上個月月底市值。
詳細方法可參考 Kenneth R. French 動能因子，或文獻 Carhart, M. M. (1997). On Persistence in Mutual Fund Performance. The Journal of Finance, 52(1), 57–82.，這裡以 Fama 他們的方式計算，即考慮規模效應。
$\rm{Fama}$ $\rm{factor}$¶
上述 FinMind 的資料我下載完後，都先轉成 csv 並存在 MySQL 裡，故先把資料從資料庫抓出來再整理。
In [4]:
import pandas as pd
import numpy as np
import mysql.connector
from sqlalchemy import create_engine
engine = create_engine("mysql+mysqlconnector://root:password@localhost/finmind_tw")
FF_MC 為公司市值的 DataFrame(df)，FinMind 原始資料為日資料共3個 columns。
先依照日期排序並把 $\tt{date}$ 字串的後 5 個字元統一成年 '12-01'，最後透過 agg 抓出每間公司每年的最後一筆市值資料作為當年市值，順便同時改個 colname。
In [5]:
FF_MC = pd.read_sql_query(
sql = 'SELECT * FROM marketcap',
con = engine).assign(date = lambda x: x.date.astype('str'))
FF_MC.head(5)
Out[5]:
| date | stock_id | market_value | |
|---|---|---|---|
| 0 | 2004-02-12 | 1101 | 47617612848 |
| 1 | 2004-02-13 | 1101 | 50927928180 |
| 2 | 2004-02-16 | 1101 | 53474324589 |
| 3 | 2004-02-17 | 1101 | 53474324589 |
| 4 | 2004-02-18 | 1101 | 51691847103 |
In [6]:
FF_MC = FF_MC.sort_values(['stock_id','date']).assign(
year = lambda x: x.date.str.slice_replace(start=-5,repl ='12-01')).\
groupby(['stock_id','year'],as_index = False).market_value.\
agg([('MC', lambda x: x.tail(n = 1) )]).loc[lambda x: x.MC != 0]
FF_AstEqt 為資產負債表的 df，FinMind 原始資料為長表共有 5 個 columns，只從資料庫取出資產負債表中的總資產和股東權益。
先把 df 轉成寬表，再把 $\tt{date}$ 的日期統一改成 '01' 並更改 colnames 為 $\tt{year}$，最後只取出 $\tt{year}$ 含'12-01'的 row，即只取出公司的年報。
In [7]:
FF_AstEqt = pd.read_sql_query(
sql = "SELECT * FROM balancesheet WHERE type IN ('Equity','TotalAssets')",
con = engine).assign(date = lambda x: x.date.astype('str'),
stock_id = lambda x: x.stock_id.astype('str'))
FF_AstEqt.head(5)
Out[7]:
| date | stock_id | type | value | origin_name | |
|---|---|---|---|---|---|
| 0 | 2012-03-31 | 1101 | TotalAssets | 2.710318e+11 | 資產總額\r |
| 1 | 2012-03-31 | 1101 | Equity | 1.414830e+11 | 權益總額\r |
| 2 | 2012-06-30 | 1101 | TotalAssets | 2.734726e+11 | 資產總額\r |
| 3 | 2012-06-30 | 1101 | Equity | 1.334564e+11 | 權益總額\r |
| 4 | 2012-09-30 | 1101 | TotalAssets | 2.687747e+11 | 資產總額\r |
In [8]:
FF_AstEqt = FF_AstEqt.pivot(index = ['date','stock_id'],
columns = 'type',
values = 'value').reset_index()#.assign(date = lambda x: x.date.astype('str'))
FF_AstEqt = FF_AstEqt.assign(
date = lambda x: x.date.str.slice_replace(start = -2,repl = '01')).rename(
columns = {'date':'year'}).loc[lambda x: x.Equity > 0]
FF_AstEqt = FF_AstEqt.loc[lambda x: x.year.str.contains('12-01')]
FF_OP 為綜合損益表的 df，FinMind 原始資料為長表共有 5 個 columns，只從資料庫取出綜合損益表中的營收、銷貨成本與營業費用。
先把 df 轉成寬表後計算營業獲利，再新增1個年 column 用於最後加總同間公司同年的年營業獲利。
In [9]:
FF_OP = pd.read_sql_query(
sql = "SELECT * FROM incomestatement WHERE type IN ('Revenue','CostOfGoodsSold','OperatingExpenses')",
con = engine).assign(date = lambda x: x.date.astype('str'))
FF_OP.head(5)
Out[9]:
| date | stock_id | type | value | origin_name | |
|---|---|---|---|---|---|
| 0 | 1990-03-31 | 1101 | CostOfGoodsSold | 2.503496e+09 | 營業成本\r |
| 1 | 1990-03-31 | 1101 | OperatingExpenses | 1.127620e+08 | 營業費用\r |
| 2 | 1990-06-30 | 1101 | CostOfGoodsSold | 2.559507e+09 | 營業成本\r |
| 3 | 1990-06-30 | 1101 | OperatingExpenses | 1.260250e+08 | 營業費用\r |
| 4 | 1990-09-30 | 1101 | CostOfGoodsSold | 2.668818e+09 | 營業成本\r |
In [10]:
FF_OP = FF_OP.pivot(index = ['date','stock_id'],
columns = 'type',
values = 'value').reset_index().dropna().assign(
OP = lambda x: x.Revenue - x.CostOfGoodsSold - x.OperatingExpenses)[['date','stock_id','OP']]
#FF_OP = FF_OP.assign(date = lambda x: x.date.astype('str'))
FF_OP = FF_OP.assign(
year = FF_OP.date.str.slice_replace(start = -5,repl = '12-01') )
FF_OP = FF_OP.groupby(['stock_id','year'],as_index = False).OP.agg([('yOP', lambda x: x.sum())])
把公司財報資料清理完後合併起來並去除金融業，公司所屬產業資料來源為公開資訊觀測站，
剔除完金融業公司後再計算對應因子的公司資料，包含：帳市值比 (BMratio)、營業獲利除以股東權益 (OPtoE)、總資產成長率 (Inv)。
In [11]:
FF_material = pd.merge(FF_AstEqt,FF_MC, on=['stock_id','year'], how='left').\
merge(FF_OP, on=['stock_id','year'], how='left')
In [12]:
TSE_company = pd.read_csv('C:/Users/e1155_l2c4ye3/Desktop/famafactor_py/TSE_company.csv')
OTC_company = pd.read_csv('C:/Users/e1155_l2c4ye3/Desktop/famafactor_py/OTC_company.csv')
bank_code = pd.concat([TSE_company,OTC_company])
bank_code = bank_code[['公司代號','產業類別']].rename(columns = {'公司代號':'stock_id','產業類別' : 'Industry'})
bank_code = bank_code[bank_code.Industry.str.contains('金融')].drop(columns = ['Industry']).stock_id.astype('str').to_list()
FF_material = FF_material[~FF_material.stock_id.isin(bank_code)]
FF_material = FF_material.sort_values(['stock_id','year']).dropna()
In [13]:
FF_material = FF_material.assign(
BMratio = FF_material.Equity / FF_material.MC,
OPtoE = FF_material.yOP/FF_material.Equity,
Inv = FF_material.groupby('stock_id').TotalAssets.transform( lambda x: np.log(x) - np.log(x.shift(1)) ) ).dropna()
最後依照因子對應的資料做分組，並保留需用到的 columns，其中 $\tt{MC}$ 用於後續投資組合的加權報酬。
FF_material 為每間公司每年在不同因子中對應的組別。
In [14]:
FF_material = FF_material.assign(
g_size = FF_material.groupby('year').\
MC.transform(
lambda x: pd.qcut(x,q = [0,0.5,1],
labels=['s1','s2']) ),
g_bm = FF_material.groupby('year').\
BMratio.transform(
lambda x: pd.qcut(x, q = [0,0.3,0.7,1],
labels=['bm1','bm2','bm3']) ),
g_op = FF_material.groupby('year').\
OPtoE.transform(
lambda x: pd.qcut(x, q = [0,0.3,0.7,1],
labels=['op1','op2','op3']) ),
g_inv = FF_material.groupby('year').\
Inv.transform(lambda x: pd.qcut(x, q = [0,0.3,0.7,1],
labels=['inv1','inv2','inv3']) ))
FF_material = FF_material[['year','stock_id','MC','g_size','g_bm','g_op','g_inv']]
FF_material.head(5)
Out[14]:
| year | stock_id | MC | g_size | g_bm | g_op | g_inv | |
|---|---|---|---|---|---|---|---|
| 1502 | 2013-12-01 | 1101 | 1.707631e+11 | s2 | bm2 | op3 | inv2 |
| 3070 | 2014-12-01 | 1101 | 1.602404e+11 | s2 | bm3 | op3 | inv2 |
| 4723 | 2015-12-01 | 1101 | 1.007964e+11 | s2 | bm3 | op2 | inv2 |
| 6406 | 2016-12-01 | 1101 | 1.297800e+11 | s2 | bm3 | op2 | inv1 |
| 8125 | 2017-12-01 | 1101 | 1.547853e+11 | s2 | bm2 | op2 | inv2 |
接著準備把 FF_material 與公司月報酬資料合併，由於報酬是月資料而 FF_material 是年資料，故先製作一個用以對應日期的 df$\tt -Corresponding$ 用來連結 FF_material 和月報酬資料，
Corresponding 的內容為報酬日期對應分組依據的財報日期，$\tt{ret\_date}$ 為月報酬日期、$\tt{year}$ 為財報日期。
有了連結用的 $\tt Corresponding$ 後就把 FF_material 與 月報酬資料 FF_ret 合併在一起，資料前處理到這完成，後續開始計算各因子的加權報酬投資組合。
In [15]:
Corresponding = pd.DataFrame({
'ret_date':pd.Series(pd.date_range(
start = '2001-07-01',end='2026-06-01', freq = 'MS')).\
dt.date.astype('str'),
# 因為 year 有重複值, 他的serires.index 也是重複值:0,0,0,1,1,1....,所以要reset_index才可以組合成df
'year':pd.Series(pd.date_range(
start = '2000-12-01',end='2024-12-01', freq = '12MS')).\
repeat(12).dt.date.astype('str').reset_index(drop = True)
})
Corresponding.head(13)
# 以下面資料為例，月報酬日期2001-07~2002-06對應的財報日期(也就是FF_material的year)為2000-12，然後在7月重新分組。
Out[15]:
| ret_date | year | |
|---|---|---|
| 0 | 2001-07-01 | 2000-12-01 |
| 1 | 2001-08-01 | 2000-12-01 |
| 2 | 2001-09-01 | 2000-12-01 |
| 3 | 2001-10-01 | 2000-12-01 |
| 4 | 2001-11-01 | 2000-12-01 |
| 5 | 2001-12-01 | 2000-12-01 |
| 6 | 2002-01-01 | 2000-12-01 |
| 7 | 2002-02-01 | 2000-12-01 |
| 8 | 2002-03-01 | 2000-12-01 |
| 9 | 2002-04-01 | 2000-12-01 |
| 10 | 2002-05-01 | 2000-12-01 |
| 11 | 2002-06-01 | 2000-12-01 |
| 12 | 2002-07-01 | 2001-12-01 |
In [16]:
FF_ret = pd.read_sql_query(
sql = 'SELECT stock_id,`month` as ret_date,ret FROM monthret_to202503',
con = engine).assign(ret_date = lambda x: x.ret_date.astype('str'))
FF_ret = pd.merge(FF_ret,Corresponding,on='ret_date',how='left',validate='m:1').dropna()
In [17]:
factor_cal = pd.merge(FF_ret,FF_material,on=['stock_id','year'],how='left').dropna()
factor_cal.head(5)
Out[17]:
| stock_id | ret_date | ret | year | MC | g_size | g_bm | g_op | g_inv | |
|---|---|---|---|---|---|---|---|---|---|
| 156 | 1101 | 2014-07-01 | 0.044147 | 2013-12-01 | 1.707631e+11 | s2 | bm2 | op3 | inv2 |
| 157 | 1101 | 2014-08-01 | 0.059382 | 2013-12-01 | 1.707631e+11 | s2 | bm2 | op3 | inv2 |
| 158 | 1101 | 2014-09-01 | -0.046188 | 2013-12-01 | 1.707631e+11 | s2 | bm2 | op3 | inv2 |
| 159 | 1101 | 2014-10-01 | 0.025388 | 2013-12-01 | 1.707631e+11 | s2 | bm2 | op3 | inv2 |
| 160 | 1101 | 2014-11-01 | -0.021550 | 2013-12-01 | 1.707631e+11 | s2 | bm2 | op3 | inv2 |
投資組合依照不同因子分開計算，計算步驟如下：
- 先計算出各組別每個月份的投組市值 (GMktCap)
- 以個股市值除以投組市值得到該月個股在投組中的市值占比 (weight)
- 以個股報酬率乘上個股在投組中的市值占比得到個股在投組中的加權報酬率 (weight_ret)
- 把同個投組中的所有個股加權報酬率加總即得到每月該投組的加權報酬率 (Portret)
- 把構成投組的兩組別名稱字串貼在一起 (Port)，並把多餘的投組 column 移除即完成該因子所有投組的計算
In [18]:
HML_data = factor_cal.assign(
GMktCap = factor_cal.groupby(['ret_date','g_size','g_bm'],as_index = False,observed=True).MC.transform(lambda x: x.sum() ),
weight = lambda x: x.MC / x.GMktCap,
weight_ret = lambda x: x.weight * x.ret).\
groupby(['ret_date','g_size','g_bm'],as_index = False,observed=True).weight_ret.agg([('Portret', lambda x: x.sum() )]).\
assign(Port = lambda x: x.g_size.astype('str') +''+ x.g_bm.astype('str')).drop(columns = ['g_size','g_bm'])
RMW_data = factor_cal.assign(
GMktCap = factor_cal.groupby(['ret_date','g_size','g_op'],as_index = False,observed=True).MC.transform(lambda x: x.sum() ),
weight = lambda x: x.MC / x.GMktCap,
weight_ret = lambda x: x.weight * x.ret).\
groupby(['ret_date','g_size','g_op'],as_index = False,observed=True).weight_ret.agg([('Portret', lambda x: x.sum() )]).\
assign(Port = lambda x: x.g_size.astype('str') +''+ x.g_op.astype('str')).drop(columns = ['g_size','g_op'])
CMA_data = factor_cal.assign(
GMktCap = factor_cal.groupby(['ret_date','g_size','g_inv'],as_index = False,observed=True).MC.transform(lambda x: x.sum() ),
weight = lambda x: x.MC / x.GMktCap,
weight_ret = lambda x: x.weight * x.ret).\
groupby(['ret_date','g_size','g_inv'],as_index = False,observed=True).weight_ret.agg([('Portret', lambda x: x.sum() )]).\
assign(Port = lambda x: x.g_size.astype('str') +''+ x.g_inv.astype('str')).drop(columns = ['g_size','g_inv'])
計算完的各因子投資組合的 df 是長表，為了方便計算先把所有 df 合併後轉成寬表，再透過 column 直接計算因子數值。
In [19]:
FF5_long = pd.concat([HML_data, RMW_data,CMA_data])
FF5_wide = FF5_long.pivot(index='ret_date',columns='Port',values='Portret').reset_index().\
assign(HML = lambda x:(x.s1bm3+x.s2bm3-x.s1bm1-x.s2bm1)/2,
RMW = lambda x:(x.s1op3+x.s2op3-x.s1op1-x.s2op1)/2,
CMA = lambda x:(x.s1inv1+x.s2inv1-x.s1inv3-x.s2inv3)/2,
SMB3 = lambda x:((x.s1bm1+x.s1bm2+x.s1bm3)-
(x.s2bm1+x.s2bm2+x.s2bm3))/3,
SMB5 = lambda x:((x.s1bm1+x.s1bm2+x.s1bm3+
x.s1op1+x.s1op2+x.s1op3+
x.s1inv1+x.s1inv2+x.s1inv3)-
(x.s2bm1+x.s2bm2+x.s2bm3+
x.s2op1+x.s2op2+x.s2op3+
x.s2inv1+x.s2inv2+x.s2inv3))/9)
最後計算市場因子 MKTPR，資料包括加權指數與無風險利率資料。
加權指數與臺灣十年期國債債券報酬率皆為日資料，加權指數取每月最後一筆數值計算月報酬，債券報酬率一樣取每月最後一筆最為當月無風險利率。
In [20]:
TAIEX_index = pd.read_sql_query(
sql = 'SELECT * FROM taiex',
con = engine).assign(date = lambda x: x.date.astype('str'))
print(TAIEX_index.head(5))
#
TW10yearbond = pd.read_csv("C:/Users/e1155_l2c4ye3/Desktop/Py_practice/臺灣十年期國債債券報酬率歷史數據.csv")
print(TW10yearbond.head(5))
date stock_id price
0 2003-01-02 TAIEX 4524.92
1 2003-01-03 TAIEX 4626.36
2 2003-01-06 TAIEX 4690.26
3 2003-01-07 TAIEX 4701.48
4 2003-01-08 TAIEX 4837.35
日期 收市 開市 高 低 升跌（%）
0 2025/4/8 1.525 1.524 1.525 1.524 -1.55%
1 2025/4/7 1.549 1.610 1.610 1.549 -3.79%
2 2025/4/2 1.610 1.610 1.610 1.610 0.00%
3 2025/4/1 1.610 1.610 1.610 1.610 0.00%
4 2025/3/31 1.610 1.630 1.630 1.610 -1.23%
In [21]:
TAIEX_mret = TAIEX_index.sort_values(['date']).assign(
ret_date = lambda x: x.date.str.slice_replace(start = -2, repl = '01')).groupby('ret_date').tail(1).assign(
ret = lambda x:(x.price-x.price.shift(1))/x.price)[['ret_date','ret']]
#
RFrate_month = TW10yearbond[['日期','收市']].rename(columns = {'日期':'date','收市':'RFrate'}).\
assign(
date = lambda x: pd.to_datetime(x.date,format='%Y/%m/%d').dt.strftime('%Y-%m-%d'),
ret_date = lambda x: x.date.str.slice_replace(start = -2,repl = '01')).\
sort_values('date').groupby('ret_date',as_index = False).tail(1).\
assign(RFrate = lambda x: x.RFrate/(12*100))[['ret_date','RFrate']]#年利率轉月利率,去除百分比
MKTPR_wide = pd.merge(TAIEX_mret,RFrate_month,on='ret_date',how='left').dropna()
MKTPR_wide = MKTPR_wide.assign(
MKTPR = lambda x: x.ret - x.RFrate)[['ret_date','MKTPR']]
#
FF5_wide = pd.merge(FF5_wide, MKTPR_wide,on='ret_date',how='left')
FF5_wide.head(5)
Out[21]:
| ret_date | s1bm1 | s1bm2 | s1bm3 | s1inv1 | s1inv2 | s1inv3 | s1op1 | s1op2 | s1op3 | ... | s2inv3 | s2op1 | s2op2 | s2op3 | HML | RMW | CMA | SMB3 | SMB5 | MKTPR | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 2014-07-01 | -0.026396 | -0.013806 | 0.011818 | -0.016122 | -0.003636 | -0.005209 | -0.007376 | -0.001083 | -0.025156 | ... | -0.006655 | -0.003731 | 0.011595 | -0.003502 | 0.035731 | -0.008775 | 0.001537 | -0.019694 | -0.014931 | 0.007662 |
| 1 | 2014-08-01 | -0.036342 | -0.015890 | 0.001389 | -0.020354 | -0.009311 | -0.017363 | -0.030037 | -0.002780 | -0.017835 | ... | 0.022984 | -0.025665 | 0.006190 | 0.027434 | 0.022899 | 0.032650 | -0.007127 | -0.035823 | -0.028553 | 0.020245 |
| 2 | 2014-09-01 | -0.029662 | -0.019099 | -0.032586 | -0.031595 | -0.026148 | -0.018325 | -0.027614 | -0.019153 | -0.038611 | ... | -0.047579 | -0.026829 | -0.047688 | -0.048998 | -0.009584 | -0.016583 | -0.014112 | 0.023747 | 0.020627 | -0.052961 |
| 3 | 2014-10-01 | -0.070788 | -0.067014 | -0.056872 | -0.072808 | -0.064529 | -0.052790 | -0.067993 | -0.064383 | -0.059546 | ... | 0.008625 | -0.042162 | -0.021088 | 0.011748 | 0.001793 | 0.031178 | -0.022846 | -0.058993 | -0.053728 | -0.000398 |
| 4 | 2014-11-01 | 0.013760 | 0.004767 | -0.004190 | 0.002906 | 0.000281 | 0.011716 | -0.004152 | 0.006470 | 0.011183 | ... | 0.037532 | -0.000578 | 0.008069 | 0.034810 | -0.021055 | 0.025361 | -0.015345 | -0.012022 | -0.012091 | 0.021794 |
5 rows × 25 columns
到這就計算完 Fama-French 五因子了。
下面繼續計算動能因子，同上述因子計算方式提到，動能因子是每月重新分組，故跟五因子分開來獨立計算。
$\rm{Momentum}$ $\rm{factor}$¶
In [22]:
FF_mom = pd.read_sql_query(
sql = 'SELECT stock_id,month AS ret_date, ret FROM monthret_to202503',
con = engine).assign(ret_date = lambda x: x.ret_date.astype('str')).sort_values(['stock_id','ret_date']).dropna()
計算動能因子只需要個股月報酬資料，載入完資料後依照每個公司計算每個月的報酬動能 $\tt mom$，先把前 12 個月到前 2 個月的報酬率加 1，再把這 11 筆數值相乘即為當月的 $\tt mom$ (累積報酬理論上最後要再減 1，但後續 $\sf mom$ 只要用來分組所以沒差)。
$\tt mom$ 計算完後即可分組，同樣以 30 與 70 百分位數分成 3 組(mom1,mom2,mom3)。
In [23]:
FF_mom = FF_mom.assign(mom = FF_mom.groupby('stock_id',as_index = False).ret.transform(
lambda x: pd.Series([1+x.shift(i) for i in range(2,13)]).prod(skipna=False)) ).dropna()
FF_mom = FF_mom.assign(
g_mom = FF_mom.groupby('ret_date').mom.transform(
lambda x: pd.qcut(x, q = [0,0.3,0.7,1],labels=['mom1','mom2','mom3']) ))
接著處理動能因子的市值分組，同樣以公司每月月底市值為當月市值，再依照每月樣本市值中位數分為2組(s1,s2)。
In [24]:
MOM_mc = pd.read_sql_query(
sql = 'SELECT date, stock_id, market_value as MC FROM marketcap',
con = engine).assign(date = lambda x: x.date.astype('str'))
MOM_mc = MOM_mc.sort_values(['stock_id','date']).assign(
mc_month = lambda x: x.date.str.slice_replace(start = -2,repl = '01')).\
groupby(['stock_id','mc_month']).tail(1)[['stock_id','mc_month','MC']]
MOM_mc = MOM_mc.assign(g_size = MOM_mc.groupby('mc_month').MC.transform(
lambda x: pd.qcut(x, q = [0,0.5,1],labels=['s1','s2']) ))
動能因子同上述因子需要一個 $\tt Corresponding$ 把分組資料與報酬資料合併，
由於 t 月的 $\tt mom$ 數值在 t-1 月就可以獲得，所以 $\tt mom$ 與報酬率使用同個日期 $\tt{ret\_date}$，
而市值是以上個月資料為準，故 $\tt mc\_month$ 的日期早 $\tt{ret\_date}$ 1 個月。
In [25]:
Corresponding = pd.DataFrame({
'ret_date':pd.Series(pd.date_range(start = '1980-02-01',end='2026-02-01', freq = 'MS')).dt.date.astype('str'),
'mc_month':pd.Series(pd.date_range(start = '1980-01-01',end='2026-01-01', freq = 'MS')).dt.date.astype('str')
})
UMD_material = Corresponding.\
merge(FF_mom,on = ['ret_date'],how = 'left').\
merge(MOM_mc,on = ['stock_id','mc_month'],how = 'left').dropna()
資料合併後同上述因子，以市值加權計算各個投組的報酬率，計算完後一樣為長表，轉成寬表後計算 $\tt UMD$ 即完動能因子所有計算。
In [26]:
UMD_long = UMD_material.assign(
GMktCap = UMD_material.groupby(['ret_date','g_size','g_mom'],observed=True).MC.transform(lambda x: x.sum()),
weight = lambda x: x.MC / x.GMktCap,
weight_ret = lambda x: x.ret * x.weight).\
groupby(['ret_date','g_size','g_mom'],as_index = False,observed=True).weight_ret.agg([('Portret',lambda x: x.sum())]).\
assign(Port = lambda x: x.g_size.astype('str') +'_'+ x.g_mom.astype('str')).drop(columns = ['g_size','g_mom'])
UMD_wide = UMD_long.pivot(index = ['ret_date'],columns = 'Port',values = 'Portret').\
assign(UMD = lambda x:(x.s1_mom3+x.s2_mom3-x.s1_mom1-x.s2_mom1)/2)
把 $\tt{FF5\_wide}$ 和 $\tt{UMD\_wide}$ 合併順便統一日期範圍，再把計算完因子的寬表轉回長表就完成計算了，最後清除上述計算過程產生的變數，留下最後計算完的的因子資料。
In [27]:
FF6_wide = FF5_wide.merge(UMD_wide,on = 'ret_date')
FF6_long = FF6_wide.melt(id_vars = 'ret_date',var_name = 'Port',value_name = 'Portret')
In [28]:
ls = %who_ls
vtr = [x for x in ls if x not in ['FF6_wide','FF6_long','pd','engine'] ]
for x in vtr:
del globals()[x]
del ls,x,vtr
$\rm{seaborn}$ 視覺化¶
透過 $\tt{seaborn}$ 把六因子的報酬趨勢圖畫出來，看看因子之間有沒有較明顯的關係。
先把六因子從所有投組中取出，順便把數值放大一百倍跟把日期 column 轉成時間格式。
In [29]:
plotdf = FF6_long.loc[lambda x: x.Port.isin(['MKTPR','SMB5','HML','RMW','CMA','UMD'])]
plotdf = plotdf.assign(ret_date = lambda x: pd.to_datetime(x.ret_date),
Portret = lambda x: x.Portret*100)
In [30]:
import matplotlib.pyplot as plt
import seaborn as sns
import matplotlib.dates as mdates
plt.rcParams['font.sans-serif'] = ['Microsoft JhengHei'] # 修改中文字體
plt.rcParams['axes.unicode_minus'] = False # 顯示負號
In [31]:
fig, ax = plt.subplots(figsize = (15,5))
# theme
plt.style.use('fivethirtyeight')
sns.lineplot(data=plotdf,
x='ret_date',y='Portret',
hue='Port',
hue_order= ['MKTPR','SMB5','HML','RMW','CMA'],
alpha=0.5,
palette = ['#404040','#2eb82e','#ff0000','#3377ff','#cc6600'],
linewidth = 1.5,
ax= ax)
# grid
plt.grid(axis= 'both',linestyle = '--')
# x axis label format
plt.xticks(ticks=plotdf.ret_date,rotation = 25)
plt.gca().xaxis.set_major_formatter(mdates.DateFormatter("%Y-%m"))
plt.gca().xaxis.set_major_locator(mdates.MonthLocator(interval= 6))
# x axis range
plt.xlim(pd.to_datetime(plotdf.ret_date.min())-pd.DateOffset(months = 2),
pd.to_datetime(plotdf.ret_date.max())+pd.DateOffset(months = 2))
# legend format
plt.legend(title="",
ncol = 6,
loc='center',
bbox_to_anchor=(0.45, 1.05),
fontsize = 14,
frameon = False,
alignment = 'left',
columnspacing = 1.2,
handletextpad = 0)
# title format
plt.title(label = 'Taiwan Stock Market Fama-French 5 Factors Monthly Return',
fontdict= {'fontsize': 20},
loc= 'center',
y = 1.1)
# second title (subtitle) format
plt.title(label = 'Period: '+
plotdf.ret_date.min().strftime(format = "%Y-%m-%d") + ' ~ ' +
plotdf.ret_date.max().strftime(format = "%Y-%m-%d"),
fontdict= {'fontsize': 14},
loc= 'right',
color = '#4d4d4d')
# axis title
plt.xlabel(xlabel = "")
plt.ylabel(ylabel = "Return (%)")
plt.show()
5條線還是太雜亂了，但 HML 與 RMW 看起來似乎有反向關係，準備用 $\tt{Tableau}$ 做個儀錶板，但需再用長表整理一些數據給 $\tt{Tableau}$ 使用。
Tableau 資料前處理¶
老實說我覺得 $\tt{matplotlib}$ 和 $\tt{seaborn}$ 比起 $\tt{ggplot2}$ 來的難用，比起用很多不同段 code 完成一張圖，我還是喜歡一整串 code 直接完成一張圖，對我來說可讀性較高使用起來也較上手，且總覺得用 Google 找 $\tt{ggplot2}$ 的資訊比 $\tt{matplotlib}$ 來的多，$\tt{matplotlib}$ 用起來更像是 R 語言原生的繪圖函數，我是沒有花時間學 R 語言原生的繪圖函數而是直接進到 $\tt{ggplot2}$，但這都不重要了~
因為有 $\tt{Tableau}$ 可用，打算用 $\tt{Tableau}$ 做一個有四張圖表的儀錶板，左上為報酬樹狀圖-TreeMap、右上為報酬分布圖-Histogram、左下為報酬趨勢圖-Lineplot、右下為累積報酬圖-Lineplot，示意圖如下。
最後先把上述提到的圖表用到的資料在 python 中整理完再匯出成 excel 檔準備匯入 $\tt{Tableau}$ 使用。
In [32]:
fig, ax = plt.subplots()
plt.hlines(y = 0.5, xmin = 0, xmax = 1, linestyles = 'solid')
plt.vlines(x = 0.5, ymin = 0, ymax = 1, linestyles = 'solid')
plt.text(x = 0.10, y = 0.70, s = ' Return\nTreemap',fontsize=15)
plt.text(x = 0.65, y = 0.70, s = ' Return\nDistribution',fontsize=15)
plt.text(x = 0.15, y = 0.15, s = 'Month\nReturn',fontsize=15)
plt.text(x = 0.65, y = 0.15, s = 'Cumulative\n Return',fontsize=15)
plt.title(label = 'Fama-French Factors Dashboard for the Taiwan Stock Market')
plt.axis('off')
plt.show()
首先計算各投組的累積報酬並以 100 為起始數值，順便把報酬率放大 100 倍把單位變為 %。
以投資組合的最後一個月的累積報酬率做排名，並把投組排名合併到每個日期。
在累積報酬率的圖表中想做出「不論選哪些投資組合，市場因子(MKTPR) 都會出現在圖表上當作基準可用於比較」的效果，因此做出 MKT_cumret 每個日期都對應到市場因子的累積報酬數值。
In [33]:
import pandas as pd
#Tableau_df = FF6_long.copy()
Tableau_df = FF6_long.sort_values(['Port','ret_date']).assign(
cumret = lambda y: y.groupby(['Port'],as_index = False).Portret.transform(lambda x: ((x+1).cumprod())*100 ),
Portret = lambda x: x.Portret * 100)
rank_dict = Tableau_df.sort_values(['ret_date']).groupby('Port',as_index = False).tail(1).assign(
Portrank = lambda x: x.cumret.rank(ascending = False))[['Port','Portrank']]
rank_dict= dict(zip(rank_dict.Port, rank_dict.Portrank))
Tableau_df = Tableau_df.assign(Portrank = lambda x: x.Port.map(rank_dict))
#
MKT_dict = Tableau_df.loc[lambda x: x.Port == 'MKTPR'][['ret_date','cumret']]
MKT_dict = dict(zip(MKT_dict.ret_date,MKT_dict.cumret))
Tableau_df = Tableau_df.assign(MKT_cumret = lambda x: x.ret_date.map(MKT_dict))
#Tableau_df.to_excel('Tableau_df.xlsx')
最後做出來的儀錶板放在我的 $\rm{Tableau\;Public}$，
要查看不同投資組合的圖表可透過點左上 Treemap Chart 做選擇，
如果要多選右下累積報酬圖表的投資組合，可透過右下圖表的選擇清單或者用 Ctrl 多選左上的 Treemap Chart。
如果透過 Ctrl 多選左上的 Treemap Chart，左下的月報酬趨勢圖會變成所選之投資組合的平均報酬趨勢圖，右上也會變成所選投資組合的報酬分布。
同上面提到的，如果同時選 HML 與 RMW，可看到報酬趨勢圖的波動相比於單選兩個別因子來的小，報酬會在更接近 0 的附近震盪，即市場的價值溢酬與獲利溢酬普遍有反向關係，當然這是建立在每年重新分組的前提下。
Appendix¶
其實計算個股在所屬投資組合裡的市值占比 (即上述各因子 df 中的 $\tt{weight}$) 可在 FF_material 就計算完，不用在計算各因子時每年重新分組後都重複計算 12 遍，但我把 code 寫完才發現，我是先在 Spyder 寫完 code 再來寫這篇 Jupyterlab 轉 html，要改的話需要改很多地方所以也懶得改了，反正資料量不大也沒啥必要，故寫在這紀錄一下。
FinMind 的損益表資料為長表，可看到 營業收入 以下的項目數量急遽減少，基本上以下的項目沒辦法完整使用，且連 營業收入 都比 營業成本 少了 3 萬筆了。
In [34]:
incomestatement_items = pd.read_sql_query(
sql = 'SELECT count(origin_name) AS counts,origin_name from incomestatement GROUP BY origin_name ORDER BY counts DESC',
con = engine)
print(incomestatement_items.to_string())
counts origin_name 0 153153 本期淨利(淨損)\r 1 146710 營業費用\r 2 144440 營業成本\r 3 117895 稅前淨利(淨損)\r 4 111480 營業收入\r 5 79789 基本每股盈餘（元）\r 6 78869 本期綜合損益總額\r 7 78584 繼續營業單位本期淨利（淨損）\r 8 77767 營業利益（損失）\r 9 77535 營業外收入及支出\r 10 77487 所得稅費用（利益）\r 11 72020 其他綜合損益（淨額）\r 12 67370 非常損益\r 13 62959 營業毛利（毛損）淨額\r 14 58308 淨利（淨損）歸屬於母公司業主\r 15 48122 營業毛利(毛損)\r 16 39419 淨利（淨損）歸屬於非控制權益\r 17 38796 停業單位損益\r 18 35382 所得稅(費用)利益\r 19 34776 繼續營業單位稅前淨利(淨損)\r 20 34043 會計原則變動累積影響數\r 21 34028 營業外收入及利益\r 22 34028 營業外費用及損失\r 23 33877 合併前非屬共同控制股權損益\r 24 33804 營業淨利(淨損)\r 25 33804 基本每股盈餘\r 26 33804 繼續營業單位淨利(淨損)\r 27 33804 聯屬公司間已實現利益\r 28 33804 聯屬公司間未實現利益\r 29 33334 所得稅(利益)\r 30 33298 營業利益\r 31 33205 調整項目\r 32 33163 停業部門損益\r 33 33095 營業外收入\r 34 33095 營業毛利\r 35 33095 營業外支出\r 36 33095 稅前純益\r 37 33095 累積影響數\r 38 33095 稅後純益\r 39 33089 每股稅後盈餘(元)\r 40 13792 綜合損益總額歸屬於母公司業主\r 41 9458 綜合損益總額歸屬於非控制權益\r 42 7551 其他收益及費損淨額\r 43 4756 未實現銷貨（損）益\r 44 4280 已實現銷貨（損）益\r 45 1812 本期其他綜合損益(稅後淨額)\r 46 1805 利息淨收益\r 47 1397 綜合損益總額歸屬於共同控制下前手權益\r 48 1344 淨利（淨損）歸屬於共同控制下前手權益\r 49 1299 利息以外淨收益\r 50 1169 本期稅後淨利（淨損）\r 51 1032 淨收益\r 52 902 繼續營業單位稅前損益\r 53 819 本期綜合損益總額(稅後)\r 54 686 呆帳費用、承諾及保證責任準備提存\r 55 636 非常損益(稅後)\r 56 636 會計原則變動之累積影響數(稅後)\r 57 636 呆帳費用\r 58 636 停業單位損益(稅後)\r 59 506 利息以外淨損益\r 60 506 繼續營業單位本期稅後淨利（淨損）\r 61 506 其他綜合損益（稅後）\r 62 500 營業收入淨額\r 63 466 保險負債準備淨變動\r 64 464 營業外損益\r 65 464 支出及費用\r 66 464 收益\r 67 453 其他綜合損益(稅後淨額)\r 68 441 繼續營業單位稅前純益(純損)\r 69 406 呆帳費用及保證責任準備提存(各項提存)\r 70 374 所得稅利益(費用)\r 71 369 合併總損益\r 72 369 繼續營業單位稅後合併淨利(淨損)\r 73 369 繼續營業單位稅前合併淨利(淨損)\r 74 369 收回(提存)各項保險責任準備淨額\r 75 369 合併總損益歸屬予母公司股東\r 76 369 合併總損益歸屬予少數股權\r 77 343 每股盈餘\r 78 285 繼續營業單位本期純益（純損）\r 79 280 呆帳費用及保證責任準備提存\r 80 267 手續費淨收益\r 81 267 繼續營業單位稅後淨利(淨損)\r 82 239 本期損益\r 83 232 會計原則變動之累積影響數\r 84 231 生物資產當期公允價值減出售成本之變動利益（損失）\r 85 207 收入\r 86 199 其他綜合損益\r 87 199 支出\r 88 181 原始認列生物資產及農產品之利益（損失）\r 89 156 繼續營業單位稅後純益(純損)\r 90 152 少數股權淨利(損失)\r 91 143 淨利（損）歸屬於母公司業主\r 92 136 繼續營業部門稅前淨利(淨損)\r 93 136 繼續營業部門稅後淨利(淨損)\r 94 117 合併前非屬共同控制股權綜合損益淨額\r 95 101 淨利（損）歸屬於非控制權益\r 96 12 淨利（損）歸屬於共同控制下前手權益\r 97 8 費損\r