點擊“簡說Python”，選擇“置頂/星標公眾號”

福利干貨，第一時間送達！

本文作者：王圣元

轉載自：王的機器

本文預計閱讀時間30分鐘，實用性高，老表建議先收藏，然后轉發朋友圈，然后吃飯、休閑時慢慢看，反復看，反復記，反復練。

0引言

SciPy 是 Python 里處理科學計算 (scientific computing) 的包，使用它遇到問題可訪問它的官網 (https://www.scipy.org/).?去找答案。 在使用 scipy?之前，需要引進它，語法如下：

import?scipy

這樣你就可以用 scipy?里面所有的內置方法 (build-in methods) 了，比如插值、積分和優化。

numpy.interpolate numpy.integrate numpy.optimize

但是每次寫?scipy?字數有點多，通常我們給 scipy?起個別名 sp，用以下語法，這樣所有出現?scipy?的地方都可以用 sp 替代。

import scipy as sp

?

SciPy 是建立 NumPy 基礎上的，很多關于線性代數的矩陣運算在 NumPy 都能做，因此就不重復在這里講了。此外在〖

本章換一種寫法，我們專門針對科學計算中三個具體問題來介紹 SciPy，它們就是

插值 (interpolation)

積分 (integration)

優化 (optimization)

對于以上每一個知識點我都會介紹一個

• 簡單例子來明晰 SciPy 里各種函數的用法

• 和金融相關的實際例子

  • 計算遠期利率：在零息曲線中插值折現因子

  • 計算期權價格：將期望寫成積分并數值求解

  • 配置資產權重：優化「風險平價」模型權重

1插值

給定一組 (x_i, y_i)，其中 i = 1, 2, ..., n，而且 x_i 是有序的，稱為「標準點」。插值就是對于任何新點?x_new，計算出對應的 y_new。換句話來說，插值就是在標準點之間建立分段函數 (piecewise function)，把它們連起來。這樣給定任意連續 x 值，帶入函數就能計算出任意連續 y 值。

在 SciPy 中有個專門的函數 scipy.interpolate 是用來插值的，首先引進它并記為 spi。?

import scipy.interpolate as spi

簡單例子

用 scipy.interpolate 來插值函數 sin(x) + 0.5x。

基本概念

首先定義 x 和函數 f(x)：

x = np.linspace(-2 * np.pi, 2 * np.pi, 11) f = lambda x: np.sin(x) + 0.5 * x f(x)array([-3.14159265, -1.56221761, -1.29717034, -1.84442231, -1.57937505, 0.? ? ? ? ?,?1.57937505, 1.84442231,1.29717034, 1.56221761, 3.14159265])

接下來介紹 scipy.interpolate 里面兩大殺器：splrep 和 splev。兩個函數名稱都是以 spl 開頭，全稱 spline (樣條)，可以理解這兩個函數都和樣條有關。不同的是，兩個函數名稱以 rep 和 ev 結尾，它們的意思分別是：

• rep：representation 的縮寫，那么?splrep 其實生成的是一個「表示樣條的對象」

• ev：evaluation 的縮寫，那么?splev?其實用于「在樣條上估值」

splrep 和 splev 像是組合拳 (one two punch)

• 前者將 x, y 和插值方式轉換成「樣條對象」tck

• 后者利用它在 x_new 上生成?y_new

一圖勝千言：

接下來仔細分析一下 tck。

tck = spi.splrep( x, f(x), k=1 ) tck

tck 就是樣條對象，以元組形式返回，tck 這名字看起來很奇怪，實際指的是元組 (t, c, k) 里的三元素：

• t - vector of knots (節點)

• c - spline cofficients (系數)

• k - degree of spline (階數)

對照上圖，tck 確實一個元組，包含兩個 array 和一個標量 1，其中

第一個 array 是節點，即標準點，注意到一開始 x 只有 11 個，但現在是 13 個，首尾都往外補了一個和首尾一樣的 x

第二個 array 是系數，注意它就是 y 在尾部補了兩個 0

標量 1 是階數，因為在調用?splrep 時就把 k 設成 1

注：前兩個 array 我只是發現這個規律，但解釋不清楚為什么這樣。它和 matlab 里面的 spline() 的產出不太一樣，希望懂的讀者可以留言區解釋一下。

雖然解釋不清楚前兩個 array，那就把 tck 當成是個黑匣子 (black-box) 直接用了。比如可用 PPoly.from_spline 來查看每個分段函數的系數。

pp = spi.PPoly.from_spline(tck) pp.c.Tarray([[ 1.25682673, -3.14159265],[ 1.25682673, -3.14159265],[ 0.21091791, -1.56221761],[-0.43548928, -1.29717034],[ 0.21091791, -1.84442231],[ 1.25682673, -1.57937505],[ 1.25682673, 0. ],[ 0.21091791, 1.57937505],[-0.43548928, 1.84442231],[ 0.21091791, 1.29717034],[ 1.25682673, 1.56221761],[ 1.25682673, 3.14159265]])

tck 的系數數組里有 13 個元素，而上面數組大小是 (12, 2)，12 表示 12?段，2 表示每段線性函數由 2 個系數確定。

把 x 和 tck 丟進 splev 函數，我們可以插出在 x 點對應的值 iy。

iy = spi.splev( x, tck ) iyarray([-3.14159265,?-1.56221761,?-1.29717034,?-1.84442231,?-1.57937505,?0.? ? ? ? ?,?1.57937505,?1.84442231,1.29717034,?1.56221761,?3.14159265])

用 Matplotlib 來可視化插值的 iy 和原函數 f(x) 發現 iy 都在 f(x) 上。Matplotlib 是之后的課題，現在讀者可忽略其細節。

除了可視化，我們還可以用具體的數值結果來評估插值的效果：

np.allclose(f(x), iy) np.sum((f(x) - iy) ** 2) / len(x)True 0.0

第一行 allclose 的結果都是 True 證明插值和原函數值完全吻合，第二行就是均方誤差 (mean square error, MSE)，0.0 也說明同樣結果。

上面其實做的是在「標準點 x」上插值，那得到的結果當然等于「標準點 y」了。這種插值確實意義不大，但舉這個例子只想讓大家

明晰?splrep 和 splev 是怎么運作的

如何可視化插出來的值和原函數的值

如何用 allclose 來衡量插值和原函數值之間的差異

一旦弄明白了這些基礎，接下來就可以秒懂更實際的例子了。

正規例子

上面在「標準點 x」上插值有點作弊，現在我們在 50 個「非標準點 xd」上線性插值得到 iyd。

xd = np.linspace( 1.0, 3.0, 50 ) iyd = spi.splev( xd, tck ) print( xd, '\n\n', iyd )

密密麻麻的數字啥都看不出來，可視化一下把。

這插得是個什么鬼？紅色插值點在第二段和深青色原函數差別也太遠了吧 (MSE 也顯示有差異)。

np.sum((f(xd) - iyd) ** 2) / len(xd)0.011206417290260647

問題出在哪兒呢？當「標準點 x」不密集時 (只有 11 個)，分段線性函數來擬合 sin(x) + 0.5x 函數當然不會太好啦。那我們試試分段三次樣條函數?(k = 3)。

tck = spi.splrep( x, f(x), k=3 ) iyd = spi.splev( xd, tck )

可視化一下并計算 MSE 看看

np.sum((f(xd) - iyd) ** 2) / len(xd)1.6073247177220542e-05

視覺效果好多了！誤差小多了！

金融例子

用 scipy.interpolate 來插值折現因子來計算遠期利率。

在金融市場中，每個貨幣都有自己相對應的折現曲線，簡單來說，就是在「標準日期」上一組折現因子 (discount factor) 或零息利率 (zero rate)。

那么在「非標準日期」上折現因子或零息利率怎么求呢？插值！

本小節的知識點內容來之〖

知識點

這里面說的插值是分段 (piecewise) 插值！對于線性插值，不是說一條直線擬合上表的 9 個點，這樣也是不可能做到的。但是分段線性插值就可以完美解決這個問題，因為 9 個點，有 8 段，每一段首尾兩個點，可以連一條直線，全部點之間連起來不就是分段線性插值嗎？三種最常見的插值方法

分段常函數

分段線性函數

分段三次樣條函數

首先給出數學符號。給定 N 數據點 (x_i,?f_i), i = 1, 2, …, N，其中?x₁ <?x₂?< ... <?x_N?。我們希望找到一個函數 f(x) 來擬合這 N 個數據點，對于分段函數，因為有 N 個數據點，需要 N -1 段函數。

---

分段常?(piecewise constant) 函數

在這種情況，每一段函數都是一個常數，這種插值方法

• 優點是簡單

• 缺點是在數據點上不連續，更不可導

• 適用于在某些模型的參數 (比如 Heston 模型中的均值回歸率和波動率的波動率) 上插值 (模型參數通常只用常數和分段常函數，但后者比前者能更好的擬合市場數據，因為它有更多自由度)。

• 不適用于曲線和波動率插值

分段常函數不連續，通常稱作?C^-1?函數。

---

分段線性?(piecewise linear) 函數

在這種情況，每一段函數都是一個線性函數，這種插值方法

• 優點是簡單，在數據點上連續，而且形狀保持性很好 (插出的值只和它相鄰兩個數據點有關，別的數據怎么動都不影響它的插值)

• 缺點是在數據點上不可導

• 適用于曲線和波動率插值

• 不適用于在 Hull-White 模型下的曲線插值 (Hull-White 模型需要對曲線求二階導)

分段線性函數連續但是不可導，通常稱作?C^0?函數。

---

分段三次樣條?(piecewise cubic spline) 函數

在這種情況，每一段函數都是一個三次多項式函數，這種插值方法

• 優點是在數據點上可導甚至可導三次 (非常平滑)

• 缺點是有些復雜，而且形狀保持性不好 (插出的值和整個數據點有關，別的數據動以下都會影響它的插值)

• 適用于曲線的插值

分段三次樣條函數連續而且二階可導，通常稱作?C²?函數。

對上面曲線插值有一個概念后，首先用 pandas 讀取數據。Pandas 是下帖內容，你就先把它當成一個可以用字符串來索引或切片的二維數據結構。

import pandas as pd curve = pd.read_excel('CNY zero curve.xlsx') curve

該曲線用于估值日 2019-04-01，上圖第一個點的日期是?2019-04-03，通常稱為即期日，往后的日期分別是從即期日開始往后推 1W, 1M, 2M, 3M, 6M, 9M, 1Y 和 2Y 得到的。

用 Matplotlib 來可視化折現因子和零息利率。

這里用了雙 y 軸來區分折現因子和零息利率，左邊是折現因子，右邊是零息利率，其實通過觀察 y 軸的數值也可以區分出來兩者。

現在實際問題是要計算從起始日 2019-08-05 到終止日 2019-11-05 的 3M 遠期利率，根據其公式 (不推導)：

要計算遠期利率，核心問題就是計算 2019-08-05 和 2019-11-05 兩天的折現因子。為了簡化，我們把這兩天之間的年限差近似定為 0.25?≈ 3個月/12個月。具體步驟如下：

計算曲線上「標準日期」到「估值日」之間的天數差

計算「起始日」和「終止日」到「估值日」之間的天數差

插出「起始日」和「終止日」上的折現因子 (四種方法)

將折現因子帶入公式計算遠期利率

---

第一步：計算曲線上「標準日期」到「估值日」之間的天數差

today = pd.Timestamp('2019-04-01') daydiff = curve['Date'] - today daydiff0 2 days 1 9 days 2 32 days 3 63 days 4 93 days 5 185 days 6 277 days 7 368 days 8 733 days Name: Date, dtype: timedelta64[ns]

上面結果不是數值型變量 (還帶個 days)，用 .dt.days.values 得到相應的 numpy 數組。

d = daydiff.dt.days.values darray([ 2, 9, 32, 63, 93, 185, 277, 368, 733], dtype=int64)

---

第二步：計算「起始日」和「終止日」到「估值日」之間的天數差

import datetime start = datetime.datetime.strptime('2019-08-05','%Y-%m-%d') end = datetime.datetime.strptime('2019-11-05','%Y-%m-%d') d_s = (start - today).days d_e = (end- today).days print( d_s, d_e )126 218

需要引進 datetime 這個庫將字符型日期轉成真正的 date 格式。

---

第三步：插出「起始日」和「終止日」上的折現因子，有多種方法，不同數據商對不同曲線也有不同的設置，常見的四種有：

在折現因子上線性插值

在折現因子上三次樣條插值

在 ln(折現因子) 上線性插值

在零息曲線上線性插值，再計算折現因子

DF 上線性插值tck = spi.splrep( d, curve['Discount Factor'], k=1 ) DF_s = spi.splev( d_s, tck ) DF_e = spi.splev( d_e, tck ) print(?DF_s,?DF_e )0.9909485166188177 0.9828538249018102

splrep 里面 k 設為 1 表示線性插值。

DF 上三次樣條插值tck?=?spi.splrep(?d,?curve['Discount?Factor'],?k=3?) DF_s = spi.splev( d_s, tck ) DF_e = spi.splev( d_e, tck ) print( DF_s, DF_e )0.9909572012597055 0.9827493083500931

splrep 里面 k 設為 3 表示三次樣條插值。

ln(DF) 上線性插值tck = spi.splrep( d, np.log(curve['Discount Factor']), k=1 ) DF_s = np.exp(spi.splev( d_s, tck )) DF_e = np.exp(spi.splev( d_e, tck )) print( DF_s, DF_e )0.9909402218834239 0.9828472942639631

把 ln(DF) 放入 splrep 中，插出來也是 ln 形式，要最終得到折現因子，還需要用 exp 函數還原。

Rate 上線性插值tck = spi.splrep( d, curve['Zero Rate (%)'], k=1 ) r_s = spi.splev( d_s, tck ) r_e = spi.splev( d_e, tck ) DF_s = np.exp(-d_s/365 * r_s/100) DF_e = np.exp(-d_e/365 * r_e/100) print( DF_s, DF_e )0.9921606726777862 0.9843810241053533

插出來的零息利率，需要用以下公式計算出折現因子

????DF = exp( -d/365 × r/100)

d 除以 365?轉換成年限，r 除以 100 是因為 r 單位是 %。

---

第四步：將折現因子帶入公式計算遠期利率

F?=?0.25*(DF_s/DF_e?-?1)?*?100

第三步中四種方法計算出來的遠期利率 (%) 為

DF?上線性插值 - 2.059%

折DF?上三次樣條插值 - 2.088%

ln(DF)?上線性插值 - 2.059%

Rate?上線性插值 - 1.976%

四個遠期利率差別都不大，業界使用較多的是第 3 和 4 種。

2積分

在 SciPy 中有個專門的函數?scipy.integrate?是用來做數值積分的，首先引進它并記為 sci。?

import?scipy.integrate?as?sci

簡單例子

用 scipy.integrate?來對函數 sin(x) + 0.5x 求積分。

首先定義被積函數 f(x)：

f = lambda x: np.sin(x) + 0.5 * x

假設我們想從 x= 0.5 到 9.5 對 f(x) 求積分，可以手推出

在 scipy.integrate?里還有些數值積分的函數：

• fixed_quad：fixed Gaussian quadrature (定點高斯積分)

• quad：adaptive quadrature (自適應積分)

• romberg：Romberg integration (龍貝格積分)

• trapz：用 trapezoidal 法則

• simps：用 Simpson’s 法則

前三個函數 fixed_quad, quad, romberg?的參數是被積函數、下界和上界。代碼如下：

sci.fixed_quad(f, a, b)[0] sci.quad(f, a, b)[0] sci.romberg(f, a, b)24.366995967084602 24.374754718086752 24.374754718086713

對后兩個函數 trapz?和?simps，首先在上下界之間取 n 個點 x_i，再求出對應的函數值 f(x_i)，再把當參數 f(x_i)?和 x_i?傳到函數中。代碼如下：

xi = np.linspace(a, b, 100) sci.trapz( f(xi), xi ) sci.simps( f(xi), xi )24.373463386819818 24.37474011548407

和解析解?24.37475471808675 比較，quad 的結果最精確。一般當被積函數不規則時 (某段函數值激增)，quad (自適應積分) 的結果也是最好。

金融例子

用 scipy.integrate 來以數值積分的形式給歐式期權定價。

注：本節主要將數值積分的用途，因此金融模型上的很多設置我們都用最簡單的，比如常數型的模型參數等等。

股票類的Black-Scholes (BS) 模型下的 SDE 是描述股票價格 (stock price) 的走勢：

其中

????S(t) = 股票價格

????r = 瞬時無風險利率

????σ = S(t)的瞬時波動率

????B(t) = 布朗運動

歐式看漲期權 (call option) 在 BS 模型下的解析解 (closed-form solution) 如下：

編寫一個程序計算 call 的解析解很容易：

這里需要引入 scipy.stats 下的 norm 庫，使用里面 cdf 函數來計算正態分布的累積分布概率。

假設股價 S0 = 100，行權價格 K = 95，利率為 5%，期限為 1 年，波動率為 10%，帶入寫好的 bscall 函數來計算期權的價值。

(S0, K, r, T, sigma) = (100, 95, 0.05, 1, 0.1) bscall( S0, K, r, T, sigma )10.405284289598569

大概記注上面的期權值 10.4053。假設我們推導能力不強或者對于更復雜的期權沒有解析解，只要知道 S_T?的分布，我們可以試著把「期望值」寫成「積分」形式，再用 x = lnS_T?做個轉換，最終可推出下式：

為了求數值積分，我們需要知道 x 是如何分布，也就是推出 x 的密度分布函數 f_X(x)，推導如下 (不是本帖的重點，如無興趣可跳過下框內容)：

給定 S 的隨機微分方程，首先用伊藤公式推出 lnS 的隨機微分方程

在 0 到 T 兩邊求積分，整理得到?S_T?的解。

其中 z 是標準正態分布變量?z ~ N(0, 1)。

用之前的變量轉化 x = lnS_T?得到 x 的解。

顯然 x 是個正態分布，均值為 lnS₀ +(rT - 0.5σ²T)，方差是?σ²T。用 NPDF 代表正態分布 (Normal) 的密度分布函數 (PDF)，可把 call 的價值寫成積分形式，其中

• 被積函數是「支付函數」和「正態分布密度分布函數?」的乘積

• 下界和上界分別是 lnK 和 +∞

最終表達式如下：

跟著「被積函數」的表達式敲代碼

mu = np.log(S0) +(r*T-0.5*sigma**2*T) v = sigma*np.sqrt(T) f = lambda x: np.exp(-r*T) * (np.exp(x)-K) * norm.pdf(x,mu,v)

定義上界和下界

(lb, ub) = (np.log(K), 7)

注意上界不要定義成 +∞。稍微分析下 x = lnS_T，當?S_T=?e^7?≈ 1097 對于 S₀ = 100 已經很大了，因此上界設為 7 比較合理。

看看三個數值積分的結果如何。

sci.quad(f, lb, ub)[0] xi = np.linspace(lb, ub, 1000) sci.trapz( f(xi), xi ) sci.simps( f(xi), xi )

10.405284289598615 10.405170993379011 10.405287100064612

結果和之前的解析解 10.4053 都相當接近。

用數值積分來求解歐式期權的確有點多此一舉 (ovekill)，但很多復雜的產品是沒有解析解的，除了用數值解的「偏微分方法有限差分法」和「蒙特卡洛法」，數值積分也是一種選擇。

3優化

在 SciPy 中有個專門的函數?scipy.optimize?是用來優化的，首先引進它并記為 spo。?

import?scipy.optimize?as?spo

優化問題可分為無約束優化 (unconstrained optimization) 和有約束優化 (constrained optimization)，我們用簡單例子來介紹前者，用金融例子來介紹后者。

簡單例子

用 scipy.optimize 來求出函數

????sin(x) + 0.05x² +?sin(y) + 0.05y²?

的最小值。

首先定義函數

f = lambda x,y: np.sin(x) + 0.05 * x**2 + np.sin(y) + 0.05 * y**2

接著可視化函數

不難看出該函數有多個局部最小值 (local minimum) 和一個全局最小值 (global minimum)。我們目標是求后者，主要步驟如下：

在 (x-y) 定義域上選點，求出函數值 f(x, y)，找出最小值對應的 x* 和 y*

用?x* 和 y* 當初始值，求出函數全局最小值

---

第一步：用蠻力找函數最小值以及對應的參數

之所以用「蠻力」一詞，是因為接下來要用到 brute 函數，而 bruteforce 就是蠻力的意思。首先定義函數 fo (其實就是上面的 f)，只不過 brute 函數要求用一個元組把若干參數集合起來。此外我們添加一個 print() 語句，為了檢查中間產出。

將 x 和 y 在 -10 到 10 以步長為 5 來切片 (回顧切片是前閉后開的，因此切片完得到的是 -10, -5, 0, 5，而不包括 10 這個點)

output = True rranges = ((-10, 10, 5), (-10, 10, 5)) spo.brute(fo, rranges, finish=None)

從上面結果可看出，函數在 (0, 0) 是取最小值 0。真是最小值嗎？我也不知道，但是以 5 為步長是不是太粗糙了些，接下來用 0.1 為步長。這時把 output 設為 False 是因為不想看到打印的內容。

output = False rranges = ((-10, 10, 0.1), (-10, 10, 0.1)) opt1 = spo.brute(fo, rranges, finish=None) opt1array([-1.4, -1.4])fo(opt1)-1.7748994599769203

當步長變小，我們能在更細的網格上計算函數值，這是函數在 (-1.4, -1.4) 取最小值 -1.7749，明顯比函數在 (0, 0) 上的值 0 要小。

---

第二步：把參數當初始值，求函數全局最小值

如果網格足夠密，上面?-1.7749 大概率是全局最小值而 ?(-1.4, -1.4) 是對應的最優解；如果網格不是足夠密，那么以 (-1.4, -1.4) 當初始值，也能很大概率找到全局最小值。

用 fmin 函數，將剛才 opt1 傳進去，并設定 x 和 f 的終止條件 xtol 和 ftol，和最多迭代次數 maxiter 和最多運行函數次數 maxfun。

output = True opt2 = spo.fmin( fo, opt1, xtol=0.001, ftol=0.001, maxiter=15,?maxfun=20?) opt2

此時最優解為 (-1.42702972, -1.42876755)，而對應的函數值為

output = False fo(opt2)-1.7757246992239009

比剛才函數在 (-1.4, -1.4) 取的最小值 -1.7749 又小了一些。

好的初始值對求函數的最優解影響比較大。假設我們無腦的用 (2, 2) 當初始值，看看會發生什么。

output = False opt3 = spo.fmin(fo, (2, 2), maxiter=250) opt3Optimization terminated successfully.Current function value: 0.015826Iterations: 46Function evaluations: 86array([4.2710728 , 4.27106945])

求得函數在 (4.2710728, 4.2710728) 取的最小值 0.015826，是不是錯的太離譜了。

金融例子

用 scipy.optimize 來用「風險平價」模型為資產配置最優權重。

本小節的知識點內容來自〖

投資組合的資產配置是個很重要的課題，投資者為了最大化回報或最小化風險，可以給各種資產配置不同的權重。本節我們看一個很流行的資產配置方法 - 風險平價 (Risk Parity, RP)。但首先我們先來看看它的通用版本，風險預算 (Risk Budgeting, RB)。

知識點風險預算?(RB)?可以基于投資者對資產未來表現?(主要是風險)?的具體看法，或一些通用原則來給資產來分配風險預算，而不是給資產分配權重。下圖畫出兩者的區別。

傳統的?FW?模型把?60%?分給股票而?40%?分給債券，但是這樣的一個投資組合?90%?的風險都來自股票只有?10%?的風險來自債券。那么這個組合更容易出現股票尾部風險?(tail risk)。一個風險更均衡的投資組合應該選擇配置更多債券?(比如?75%)?和更少股票?(比如?25%)，如下圖所示。

RB?模型的思路就是通過分配風險?(上圖的風險比例)?來影響權重?(上圖的資產權重)，通常是給風險低的資產?(如債券)?高風險配額，而風險高的資產?(如股票)?低風險配額。

接下來我們看看?RB?模型的數學公式吧，首先回顧組合風險

對于第?i?個資產，其邊際風險貢獻?(Marginal Risk Contribution, MRC)?是該資產權重?w_i?的微小變化對組合風險?σ_p?所帶來的影響。用數學公式表示就是對?σ_p?求?w_i?的偏導數。

第?i?個資產的總體風險貢獻?(Total Risk Contribution, TRC)?是其?MRC?乘以其權重，顧名思義，這個總體貢獻一方面來自?MRC，一方面來自權重，數學表達式為：

根據?TRC_i?的定義，即第?i?個資產對總體風險的貢獻，可推出它們總和應該等于組合風險?s_p，從數學上也可證實此關系

上式兩邊同時除以?σ_p，并定義風險預算?s_i?為?TRC_i?的占比，可得?s^T1?= 1

由上式看出?s_i?也類似于權重，只不過是風險上的權重，而?w_i?是資產上的權重。下面給出?s_i?和?w_i?之間的關系

在?RB?模型中，股票權重等于風險預算除以貝塔，因此，RB?模型依賴于貝塔的預測質量。歸一化之后的權重等于

事先將一組風險預算分配好，例如?s?= [20%, 30%, 50%]，再數值解下面序列二次規劃?(Sequential Quadratic Programming, SQP)?問題可以得到?RB?模型下的最佳權重

風險平價 (RP) 就是等量的風險預算，即為投資組合中的所有資產分配相等的風險。

知識點

類比?RB，RP?給每個風險配額?s_i?的分配?1/n?的權重，這時組合權重為

同樣可得到?RP?模型下的優化問題?(用?1/n?替代?s_i?)

這是一個有約束 (constrainted) 的優化問題，我們可用?scipy.optimize?里的 minimize 函數來求解 RP 的權重。首先來定義 risk_parity 函數：

該函數的兩個參數 sigma 和 rho 是 n 個資產的波動率向量 (一維數組) 和相關系數矩陣 (二維數組)，其中

• obj 就是用 numpy 把上面目標函數用「匿名函數」的形式表示出來

• 限制條件有兩種形式，等式 (eq) 和不等式 (ineq)，分別用 dict 的形式表達，而限制條件的表達式也用「匿名函數」來表示

最后在 minimize 函數設定好目標函數、初始值、算法、限制條件和終止條件，得到一個 dict 類的結果 w。

兩個資產

先分析簡單的股票和債券兩個資產組合：

• 股票的預期超額回報為?10%，波動率為?20%

• 債券的預期超額回報為?5%，波動率為?10%

• 它們相關系數為 -10%

mu = np.array([0.1, 0.05]) sigma = np.array([0.2, 0.1]) rho = np.array([[1, -0.1], [-0.1, 1]])

運行 risk_partiy 函數

result = risk_parity( sigma, rho ) result? ? fun: 3.26901989274624e-15jac: array([-1.86742928e-07, 1.55459627e-07])message: 'Optimization terminated successfully.'nfev: 22nit: 5njev: 5status: 0success: Truex: array([0.33333332, 0.66666668])

result 是一個字典：

• 'fun' 對應的是目標函數在最優解下的值，非常小接近于零證明找到了最優值。

• 'nfev' 對應的 22 表示目標函數被運行了 22 次

如果只關注最優權重，只需看 ‘x’

result.xarray([0.33333332, 0.66666668])

股票和債券的最優權重為 w* =[33.33%, 66.66%]

三個資產

接著分析股票、債券和信貸三個資產組合：

• 股票的預期超額回報為?10%，波動率為?20%

• 債券的預期超額回報為?5%，波動率為?10%

• 信貸的預期超額回報為?10%，波動率為?15%

• 股票與債券、股票與信貸、債券與信貸的相關系數為?-10%, 30%, -30%

mu = np.array([0.1, 0.05, 0.1]) sigma = np.array([0.2, 0.1, 0.15]) rho = np.array([[1, -0.1, 0.3], [-0.1, 1, -0.3], [0.3, -0.3, 1]])

運行 risk_partiy 函數

result?=?risk_parity(?sigma,?rho?) result.xarray([0.19117648, 0.5147059 , 0.29411762])

股票、債券和信貸的最優權重為 w* = [19.12%, 51.47%, 29.41%]

4總結

本帖只討論用 SciPy 可以實現的三個應用

• 用 scipy.interpolate 來插值?(interpolation)

• 用 scipy.integrate 來積分?(integration)

• 用 scipy.optimize 來優化?(optimization)

學完此貼后，至少你可以解決以下具體金融問題 (training set?)：

在折現曲線上插出折現因子和零息利率

數值積分求解期權價值

優化出風險平價模型的權重

舉一反三一下，你還可以解決新的金融問題 (test set )：

在波動平面上插出波動率

數值積分求解而二維金融衍生品價值

優化出各種資產配置模型的權重 (加各種約束)

我是老表，踏實的人更容易過好生活，本文完。

推薦閱讀：

仔細閱讀下面四篇文章，2小時快速掌握Python基礎知識要點。

完整的Python基礎知識要點

---

/今日留言打卡主題/

仔細想想，你學Python有多久了？你最大的感受是什么？

（參與留言打卡留言字數不少于30字，格式：day xx,balabala,否者視為打卡失敗）

留言打卡獎勵

堅持連續打卡21天，展示進入免費分享資源群堅持連續打卡64天，獲得50以內元型態一贈書本。堅持連續打卡100天，獲得100以內元型態一贈書本。堅持連續打卡101天，免費加入老表的星球知識所有獎品一經兌換，打卡天數也自動清零。

總結

以上是生活随笔為你收集整理的Python小白数据科学教程：SciPy精讲的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。