このページは、nAGライブラリのJupyterノートブックExampleの日本語翻訳版です。オリジナルのノートブックはインタラクティブに操作することができます。
nAGライブラリを使用した高速インプライドボラティリティ計算
ヨーロピアン・コール・オプションの価格に対するブラック・ショールズの公式は以下の通りです:
ここで、
金融において重要な問題は、
以下の図に示すように、ボラティリティ・サーフェス(価格と満期までの期間に応じてボラティリティがどのように変化するかを示す3次元プロット)は非常に曲率が高くなる可能性があります。これにより、インプライドボラティリティを正確に計算することは困難な問題となります。
新しいnAGライブラリルーチンを紹介する前に、汎用的な根探索アルゴリズムを使用してインプライドボラティリティを計算する方法を示してみましょう。まず、いくつかのものをインポートする必要があります:
import numpy as np
import time
from naginterfaces.library import rand, roots, specfun
nAGライブラリのランダム数生成器を使用して、いくつかの入力データを生成しましょう:
= 10000 # This is the number of volatilities we will be computing
n = rand.init_nonrepeat(1)
statecomm = rand.dist_uniform(n, 3.9, 5.8, statecomm)
p = rand.dist_uniform(n, 271.5, 272.0, statecomm)
k = rand.dist_uniform(n, 259.0, 271.0, statecomm)
s0 = rand.dist_uniform(n, 0.016, 0.017, statecomm)
t = rand.dist_uniform(n, 0.017, 0.018, statecomm) r
上記の様々な一様分布の限界を選択したのは、入力データが 妥当な値を取るようにするためです。
インプライドボラティリティを計算するために使用できる様々な標準的な根探索技術があります。
一般的な例としてはバイセクション法があります。nAGライブラリのルーチンcontfn_cntin
は、
割線法に基づく汎用的な根探索プログラムです。これはコミュニケーションオブジェクトを介してデータを渡すコールバックを使用します:
def black_scholes(sigma, data):
try:
= specfun.opt_bsm_price('C', [data['k']], data['s0'], [data['t']], sigma, data['r'], 0.0)
price except:
= np.zeros((1,1))
price return price[0, 0] - data['p']
contfn_cntin
operates on scalars, so we need to call the
routine once for every volatility we want to compute. We will time the
computation and count how many errors are caught:
= {}
data = 0
errorcount = time.perf_counter()
tic for i in range(n):
'p'] = p[i]
data['k'] = k[i]
data['s0'] = s0[i]
data['t'] = t[i]
data['r'] = r[i]
data[try:
= roots.contfn_cntin(0.15, black_scholes, 1.0e-14, 0.0, 500, data)
sigma if sigma < 0.0:
+= 1
errorcount except:
+= 1
errorcount = time.perf_counter()
toc print('Using a general purpose root finder:')
print(' Time taken: {0:.3f} seconds'.format(toc-tic))
print(' There were {} failures'.format(errorcount))
Using a general purpose root finder:
Time taken: 51.561 seconds
There were 4 failures
Can a bespoke implied volatility routine do better? Our new routine
at Mark 27.1 is called opt_imp_vol
. We call it as
follows:
sigma, ivalid = specfun.opt_imp_vol('C', p, k, s0, t, r, mode=2)
戻り値の引数 ivalid
は、ボラティリティが計算できなかったデータポイントを記録する配列です。引数
mode
を使用して、使用するアルゴリズムを選択できます -
これについては後ほど詳しく説明しますが、 今は mode=2
を選択します。これにより、Jäckel (2015)
のアルゴリズムが選択されます。これは3次のハウスホルダー反復に基づく非常に精度の高い方法です。
以下は、タイミング計測コードで囲まれた呼び出しです:
= time.perf_counter()
tic = specfun.opt_imp_vol('C', p, k, s0, t, r, mode=2)
sigma, ivalid = time.perf_counter()
toc print('omp_imp_vol with mode = 2 (Jäckel algorithm):')
print(' Time taken: {0:.5f} seconds'.format(toc-tic))
print(' There were {} failures'.format(np.count_nonzero(ivalid)))
omp_imp_vol with mode = 2 (Jäckel algorithm):
Time taken: 0.01415 seconds
There were 0 failures
新しいルーチンは根探索法よりも数桁速く、失敗も報告されていません。nag_roots_contfn_cntinの収束パラメータや反復回数の上限を調整したり、コールバックを通じてデータを渡す方法を改善したりすることはできますが、opt_imp_volの性能に匹敵することは難しいでしょう。
最近、nAGはロンドン大学クイーン・メアリー校の数学者たちとの共同研究を開始しました。彼らはインプライド・ボラティリティを計算するための代替アルゴリズムを開発しています。この新しいアルゴリズム(Glau et al.(2018)に基づく)は、チェビシェフ補間を使用して分岐を除去し、SIMDのパフォーマンスを向上させています。opt_imp_volの呼び出し時にmode=1を設定することで、このアルゴリズムにアクセスできます:
= time.perf_counter()
tic = specfun.opt_imp_vol('C', p, k, s0, t, r, mode=1)
sigma, ivalid = time.perf_counter()
toc print('omp_imp_vol with mode = 1 (Glau algorithm):')
print(' Time taken: {0:.5f} seconds'.format(toc-tic))
print(' There were {} failures'.format(np.count_nonzero(ivalid)))
omp_imp_vol with mode = 1 (Glau algorithm):
Time taken: 0.01770 seconds
There were 0 failures
システムによっては、同程度の精度で、Jäckelアルゴリズムに比べてわずかな速度向上が見られるはずです。数値実験によると、非常に小さな配列(100要素未満)の場合、実際にはJäckelアルゴリズムの方がGlauらのアルゴリズムよりも性能が良いですが、より大きな配列では逆になります。将来的にベクトルユニットが改善されるにつれて、高度にベクトル化可能なGlauらのアプローチの性能も同様に向上すると予想されます。
これまで、倍精度に可能な限り近い相対精度でインプライドボラティリティを計算してきました。しかし、一部のアプリケーションでは、インプライドボラティリティは数桁の精度があれば十分な場合があります。Glauらのアルゴリズムの利点の1つは、より低い精度モードで実行でき、7桁の精度のみを目指すことができることです。これは、opt_imp_volの呼び出し時にmode=0を設定することでアクセスできます。これにより、ルーチンの速度がおおよそ2倍になります:
= time.perf_counter()
tic = specfun.opt_imp_vol('C', p, k, s0, t, r, mode=0)
sigma, ivalid = time.perf_counter()
toc print('omp_imp_vol with mode = 0 (lower accuracy Glau algorithm):')
print(' Time taken: {0:.5f} seconds'.format(toc-tic))
print(' There were {} failures'.format(np.count_nonzero(ivalid)))
omp_imp_vol with mode = 0 (lower accuracy Glau algorithm):
Time taken: 0.01556 seconds
There were 0 failures
以下のグラフは、Intel
Skylakeマシンで収集したタイミングを使用して結果をまとめたものです。Glauらのアルゴリズムは大規模な配列では
Jäckelアルゴリズムを上回りますが、小規模な配列では上回らないことがわかります。なお、汎用の根探索アルゴリズムはopt_imp_vol
よりもはるかに遅いため、ここでは省略しています。
要約すると、nAGの新しい最先端アルゴリズムは、入力配列の長さと必要な精度に応じて3つの異なるモードで実行できます。詳細情報やnAGライブラリへのアクセスについては、以下のURLをご覧ください: https://nag.com/content/nag-library
参考文献
P. Jäckel (2015). Let’s be rational. Wilmott 2015, 40-53.
K. Glau, P. Herold, D. B. Madan, C. Pötz (2019). The Chebyshev method for the implied volatility. Journal of Computational Finance, 23(3).
nAG Library for Pythonのセットアップ
nAG Library for Pythonのインストール手順については、以下のドキュメントをご覧ください: https://www.nag-j.co.jp/naglib/cl/python/install_guide.htm