前言
「最佳化」是很酷的觀念，因為現實世界中許多問題，並沒有嚴謹一致的公式解，但可以利用計算機高速運算能力，透過巧妙的演算法，迭代式反覆逼近最佳解，應用領域非常廣。若能多瞭解一點原理，一定可以提昇解決問題的能力。今天從網路上發現一堂手把手的教學課程，就來演練一下整個過程。期望徹底了解之後，後面可以有些延伸應用。
課程連結在此：https://youtu.be/9GA2WlYFeBU?si=0oEjPXGlrk1VrBDm
感謝專業講師 Ryan O'Connell, CFA, FRM 的分享。課程演練
問題描述
對於一籃子的資產種類，如何調整各種類所佔的比例 (投資權重)，使報酬-風險比 (sharpe ratio) 最大。
資料準備
- 從 yahoo finance 取得資料，只取一個欄位「調整後收盤價」，存在一個 dataFrame
import yfinance as yf
import pandas as pd
from datetime import datetime, timedelta
import numpy as np
from scipy.optimize import minimize
import matplotlib.pyplot as plt
import seaborn as sns
# 定義一籃子股票 (台積電，聯發科，聯詠，瑞昱，王品，中美晶，台達電)
tickers = ['2330.TW', '2454.TW', '3034.TW', '2379.TW', '2727.TW', '5483.TWO', '2308.TW']
# 設定七年的期間
end_date = datetime.today()
start_date = end_date - timedelta(days = 7*365)
adj_close_df = pd.DataFrame()
for ticker in tickers:
data = yf.download(ticker, start = start_date, end = end_date)
adj_close_df[ticker] = data['Adj Close']
adj_close_df
初步分析
從收盤價，計算每日報酬率，記得去掉空資料，以免後續計算錯誤。
log_returns = np.log(adj_close_df / adj_close_df.shift(1))
log_returns = log_returns.dropna()
log_returns
有了每日報酬，已經可以得到各資產的歷史表現，以年化指標方式表達。
log_returns.mean()*252
----------------------
2330.TW 0.243656
2454.TW 0.253070
3034.TW 0.291873
2379.TW 0.252808
2727.TW 0.123987
5483.TWO 0.142502
2308.TW 0.182675
dtype: float64
log_returns.std() * np.sqrt(252)
--------------------------------
2330.TW 0.275972
2454.TW 0.378044
3034.TW 0.374934
2379.TW 0.367530
2727.TW 0.354785
5483.TWO 0.420420
2308.TW 0.298477
dtype: float64
接下來計算各標的間的關聯性，以「共變異矩陣」方式表達。
cov_matrix = log_returns.cov()*252
cov_matrix
函式準備
# 計算 portfolio 的年化報酬
def expected_return (weights, log_returns):
return np.sum(log_returns.mean()*weights) * 252
# 計算 portfolio 的標準差，也就是 volatility，或稱波動性
# 重要！！
def standard_deviation (weights, cov_matrix):
variance = weights.T @ cov_matrix @ weights
return np.sqrt(variance)
# 計算 portfolio 的 shapre ratio
def sharpe_ratio (weights, log_returns, cov_matrix, risk_free_rate):
return (expected_return(weights, log_returns) - risk_free_rate) /
standard_deviation(weights, cov_matrix)
# 最佳化的目標函數 (最大化 sharpe，就是最小化 負sharpe)
def neg_sharpe_ratio(weights, log_returns, cov_matrix, risk_free_rate):
return - sharpe_ratio(weights, log_returns, cov_matrix, risk_free_rate)
參數準備
# 執行最佳化時，所需要的重要參數
# 權重加總，必須 = 1
constraints = {'type': 'eq', 'fun': lambda weights: np.sum(weights) - 1}
constraints
-----------
{'type': 'eq', 'fun': <function __main__.<lambda>(weights)>}
# 所有權重，都限制在 0 到 0.5 之間，也就是不准單一標的，佔比超過五成
bounds = [(0, 0.5) for _ in range(len(tickers))]
bounds
------
[(0, 0.5), (0, 0.5), (0, 0.5), (0, 0.5), (0, 0.5), (0, 0.5), (0, 0.5)]
risk_free_rate = 0 # 假設忽略無風險利率
# 權重初始值
initial_weights = np.array([1/len(tickers)] * len(tickers))
initial_weights
---------------
array([0.14285714, 0.14285714, 0.14285714, 0.14285714, 0.14285714,
0.14285714, 0.14285714])
執行最佳化
# 最關鍵的一步，scipy 套件提供 minimize()，未提供 maximize()
optimized_results_sharpe = minimize(neg_sharp_ratio, initial_weights,
args=(log_returns, cov_matrix, risk_free_rate), method='SLSQP',
constraints = constraints, bounds = bounds, options = {'disp': True})
--------------------------------------------------------------------------
Optimization terminated successfully (Exit mode 0)
Current function value: -0.9342487059203365
Iterations: 5
Function evaluations: 40
Gradient evaluations: 5
為了顯示結果，先準備一個函式，讓結果顯示在一張圖表裡：
def show_result(result):
optimal_weights = result.x
optimal_portfolio_return = expected_return(optimal_weights, log_returns)
optimal_portfolio_volatility =
standard_deviation(optimal_weights, cov_matrix)
optimal_sharpe_ratio =
sharpe_ratio(optimal_weights, log_returns, cov_matrix, risk_free_rate)
sns.set_theme(style='whitegrid')
sns.set_color_codes('pastel')
f, ax = plt.subplots(figsize=(7, 5))
title1 = f'Return {optimal_portfolio_return: .4f}'
title2 = f'with Volatility {optimal_portfolio_volatility: .4f}'
title3 = f'and sharpe {optimal_sharpe_ratio: .2f}'
ax.set_title(f'{title1} {title2} {title3}')
sns.barplot(x='weight', y='ticker',
data=pd.DataFrame({'ticker': tickers, 'weight': optimal_weights}),
label='weight', color='b')
以上已經完整掌握所有定義，設定和計算流程。
檢討與變形
最小化波動
以上過程最重要的基本邏輯就是最大化 sharpe，就是報酬 / 波動，但用以過的歷史資料的平均報酬，當作面對未來績效的期望值估計量，明顯存在著大量例外。但本來預測未來就是困難的事，目前也沒有更好的作法，就算專業機構也只能將就使用。但個人研究無需包袱纏身，就大膽來做一個變形，就是「忽略報酬率」！只求「最小化波動」。計算 portfolio 標準差的函式已經存在，就把它直接當成目標函數，執行最佳化作業：
optimized_results_volatility = minimize(standard_deviation, initial_weights,
args=(cov_matrix), method='SLSQP', constraints = constraints,
bounds = bounds, options = {'disp': True})
--------------------------------------------
Optimization terminated successfully (Exit mode 0)
Current function value: 0.2238746037968142
Iterations: 8
Function evaluations: 64
Gradient evaluations: 8
得到結果如下：
雖然 sharpe 降至 0.91，但波動從 24% 降至 22%，不錯，這就是我所要的。歷史波動對未來的預測力，我相信是很有參考價值的，比用單純算術平均的報酬去預測未來有用多了。當波動降到很低，我們就可以審慎的逐步擴大槓桿倍數了。槓桿不是洪水猛獸，要提升絕對績效，這是一定要的！
納入空頭部位
為了更近一步降低波動，「做空」就是必要之惡了，因為再怎樣「長期向上」的市場，也免不了無法預測的「中級回調」！屆時如果部位只有 long-only，免不了大幅縮水，痛苦程度不見得扛得住。簡言之還是降低波動的概念。所以需要改一下 constraint and bounds：
# 允許權重為負值，代表做空，而絕對數字依然一樣，代表整體曝險市值的比例。
bounds = [(-0.6, 0.6) for _ in range(len(tickers))]
constraints = {'type': 'eq', 'fun': lambda weights: np.sum(abs(weights)) - 1}
optimized_results_volatility_long_short =
minimize(standard_deviation, initial_weights, args=(cov_matrix),
method='SLSQP', constraints = constraints, bounds = bounds,
options = {'disp': True})
執行結果與上一次相同，似乎此演算法對 initial weight 的影響很大！怎麼跑都還是維持在正直，看來選股難，選空方標的更難。試著加上「大盤指數」看看，重新抓取資料。此時因為未考慮報酬率，所以除權息影響也就甚微了。同時特別指定此標的的 initial weight 為負值，也就是用做空指數來平衡多頭風險：
# 重新指定 ticker 集合，重新執行抓取資料程序
tickers = ['^TWII', '2330.TW', '2454.TW', '3034.TW', '2379.TW', '2727.TW',
'5483.TWO', '2308.TW']
# 指定初值時，特別把第一個元素變號
initial_weights = np.array([1/len(tickers)] * len(tickers))
initial_weights[0] = - initial_weights[0]
重新執行最佳化，得到很棒棒的「超低風險」(0.03) 的投資組合：
此結果表示，大約配置 50% more 的指數空頭部位，可充分抵銷投資組合的風險。雖然報酬率也被犧牲了，但「歷史績效永遠不等於未來績效」，而我相信「歷史風險大約等於未來風險」。 雖然全球經濟長期向上的事實還是值得相信的，但中期回檔還是很痛苦的，例子不遠，2022 就經歷過了，所以我偏向經常性的 Long-short 配置，在 portfolio 穩定的前提下，加大槓桿與內部波動共舞。
大致完成一回合的探索過程，透過變化問題定義，限制條件，最佳化函數等等，相信會找到新的應用，趣味也將持續。
Newman 2024/11/6
導覽頁：紐曼的技術筆記-索引
