量子ファイナンスと古典的方法

量子ファイナンス

　量子コンピュータが得意とする分野のうち、ファイナンス分野で活躍が期待されているものがいくつかある。下記(Adam Bouland et al. 2020. "Prospects and challenges of quantum finance"  )では特に３つ挙げられており、それが「モンテカルロ法」「ポートフォリオ最適化」「機械学習」である。本稿ではこのうち前者２つをとりあげる。
　これらの量子Verの詳細までは触れないが、古典的な方法に比べて何が期待できるのか見ていく。

Prospects and challenges of quantum finance

モンテカルロ法

　モンテカルロ法はランダムサンプリングを用いて近似的な解を求める方法のことである。金融の数値解析としては、シミュレーションをたくさんやってみていろんなシナリオの結果を集めよう、という感じである。たとえば６面のサイコロで１の目が出る確率を調べようと思ったら、それを何回も振ってみればいい。我々が良く知るサイコロであればその確率は1/6になるはずで、サイコロを振る回数が多ければ多いほど出目の結果が真の確率に近づく。

　これはオプション価格の期待値計算でも同じように考えることができる。ある単純なヨーロピアンコールオプションの価格を計算したいとき、
現在価格　P
権利行使価格　K
満期の期間　T
リスクフリーレート　r
ボラティリティ　σ
を設定する。これに対し正規乱数を与えて満期時点で権利行使の有無をシミュレートすることで、オプション価格が算出できる。この作業を100回やると100個のペイオフが算出できる。この100個のペイオフの平均（の現在価値）はオプションの真の価格に近いはず、というわけだ。もちろんたくさんやればやるほど精度が上がるのである。
　だが上のような極めて簡単なオプション価格はブラックショールズモデルで計算できる。いわゆる解析解が計算できるのだ（というか上のはまんまBS式の説明である）。モンテカルロ法が力を発揮するのはこのような解析解が求まるものではなく、経路依存やバスケットなどより難しい商品に対してである。ここら辺は専門でないので実務でどうなっているかはわからないが、ポンポン解析解が求まる手法ばっかりではなかろうとは想像に難くない。（さらに詳細は下記など参照）

Monte Carlo methods for option pricing - Wikipedia

　ちなみに手元でBS式とシミュレーションで比較するとこんな感じ(S=100,K=120,σ=20%,r=0.01%,T=3)。
解析解　7.15
10回 2.37
100回 8.19
1000回 6.45
10000回 6.81
100000回 7.16
　チェリーピッキングっぽい感じもするが、試行回数を増やすと確かに解が近づく。

　解析解でなく数値解を求める方法では、例えば勾配法は効率的にパラメータを求めることができるが、それが最適解ではなく局所解である可能性があり、初期値による問題が残る。そういった場合にマルコフ連鎖モンテカルロ法などを用いることで局所解を脱し、最適解に近づきやすくなることがある。
　そもそもモンテカルロ法にはいくつもパターンがあり、サンプリングのアルゴリズムは複数存在する。それっぽいところから多めにサンプリングをとる重点サンプリングや、オプションの文脈で言えばアメリカンオプションの価格を求める最小二乗モンテカルロ法というのがある。

量子的スピードアップ

　さて、モンテカルロ法で誤差を十分小さくするためにはそれなりにサンプル数（計算量）が必要となる。どれくらい必要かというと、具体的にはモンテカルロ積分をチェビシェフの不等式に当てはめるのだが、結果として古典的なモンテカルロ法では$${O(σ^2/ε^2 )}$$の計算量が必要となる。実際に誤差を十分小さくしたい場合、この二乗によるオーダーの増加は致命的となる。例えば誤差0.01に収めたい時、10,000回の施行が必要になったりする。これが量子モンテカルロ法を使うと$${ O(σ/ε)}$$となり、εの二次的なスピードアップが図られる。具体的には量子振幅推定や振幅増幅のアルゴリズムが使われる。

Quantum speedup of Monte Carlo methods

Groverのアルゴリズム

　イメージのために、量子アルゴリズムの代表格で、量子による計算速度向上のお手本でもあるGroverのアルゴリズムを紹介しておく。というか、Groverのアルゴリズムの一般化が量子振幅推定や量子振幅増幅である。
　下の記事でも紹介されているとおり、Groverのアルゴリズムはマーキングと平均値周りの反転を繰り返すことで所望の状態を得る確率を１に近づけていく処理である。

グローバーのアルゴリズムの一般化　Quantum Amplitude Amplification の論文解説 - Qiita

　量子計算においては、測定するまで値がわからない量子ビッドに対して、特定の状態の振幅を増幅させることで所望の状態が観測されやすいようにする。振幅増幅とは、ある空間内で作用する回転演算子を初期ビッドに対して複数回かけ合わせる処理だが、所望の状態がとる確率を最適化する問題は、この回転のなかで最適な角度θを求める問題に置き換えられる。

　角度θとは何かというと、所望の状態$${|Ψ_1⟩}$$とそれ以外の状態$${|Ψ_0⟩}$$で張られる２次元平面において、重ね合わせ状態である$${|Ψ⟩}$$を下記のとおり定義したときに出てくるものである。
$${|Ψ⟩=cos⁡θ |Ψ_0 ⟩+sin⁡θ |Ψ_1 ⟩}$$
　ここで、そもそも$${|Ψ⟩}$$は$${|x⟩}$$の取りうる量子状態の重ね合わせ状態であり、n量子ビッドにアダマール変換を行ったものである。つまり上の式は、もともと下記のように変換されてきたものと言える。

$$
|Ψ⟩=\displaystyle\frac{1}{\sqrt{N}}\sum_{i=0}^{N-1} |x⟩
=\displaystyle\frac{1}{\sqrt{N-1}}\sum_{i=0}^{N-1} |x'⟩+\frac{1}{\sqrt{N}}|x⟩
=cos⁡θ |Ψ_0 ⟩+sin⁡θ |Ψ_1 ⟩
$$

　ここで確率振幅であるsinθが1に近ければ近いほど、量子ビッドを測定した時に所望の状態が得られる確率が高いというわけである。もちろんこのまま測定してしまうと、所望の状態を得る確率は1/Nであり、測定結果は使い物にならない。

　そこで$${|Ψ⟩}$$にいろいろと処理をしていくのであるが、これを細かく説明するとしんどいので、ある演算子Uをm回かけ合わせるということにする。
　すると結果は
$${U^m|Ψ⟩=cos⁡(2m+1)θ |Ψ_0 ⟩+sin⁡(2m+1)θ |Ψ_1 ⟩}$$
と書き表せる。
　ここで$${sin⁡(2m+1)θ =1 ⇔⁡(2m+1)θ=π/2}$$となるのがベスト、というわけである。

　ではこれをmについて解いていく。事前準備として、Nが非常に大きく、そのためsinθがθで近似できるという状態を想定する。この時、
$${sin⁡(2m+1)θ ≈⁡(2m+1)θ≈π/2}$$
であり、もともと$${sinθ=\displaystyle\frac{1}{\sqrt{N}}}$$だったので
$${(2m+1)\displaystyle\frac{1}{\sqrt{N}}≈π/2}$$
と近似できる。これを整理すると
$${m≈\displaystyle\frac{π}{4}({\sqrt{N}-1)}}$$
となる。

　ここでNが非常に大きいとすれば右辺はほぼ$${\sqrt{N}}$$と見なせるので、最適な確率振幅を得るための計算回数はたかだか$${O(\sqrt{N})}$$と言うことができる。

モンテカルロ積分

　話があっちに行ったりこっちに行ったりで申し訳ないが、次にモンテカルロ積分に触れておく。モンテカルロ積分というのは、ある入力xに対して確率p(x)で結果f(x)を得るときに、結果×確率をxで積分して期待値を求めることである。これの試行回数が十分大きければ、離散的なサンプリングの結果を試行回数で割ったものを連続型関数の積分の近似とみなせる、というのがざっくりとしたアイデアである。解析的に期待値を求めるのは難しい一方でサンプリングは可能という条件や、サンプルの生成、試行回数などいろいろと制約はある。
$${ \displaystyle\int_ {}f(x)p(x) dx ≈\frac{1}{N}\sum_{}f(x)}$$

MCMC

　さらに話は外れるが、先ほど話に出たマルコフ連鎖モンテカルロ法にも少し触れておきたい。いわゆるMCMCと略されるやつで（なんとなくかっこいい）、効率的なサンプリングの方法として知られている。マルコフ連鎖とは１個前の状態のみによって次の状態が予測される連鎖のことである。返して言えば、１個前より昔の状態というのは将来予測には役に立たないよ、という前提の世界である。

　もう少し詳しく触れると、ある定常分布を満たすような遷移確率$${Q}$$を設計し、初期状態$${π_0}$$から遷移させていくことで求める定常分布$${π}$$に行きつくようにする。この時初期状態はなんでもよくて、気を使わなくて良いのが利点である。（厳密な定義や条件を省いてますがご勘弁を。。）

　マルコフ連鎖モンテカルロ法の代表的なアルゴリズムであるメトロポリス法では、あるサンプルに対し、提案分布から次の値の候補を生成する。ただしこの候補が採択されるかどうかは目標分布によって変わってくる。ざっくりと言えば目標分布に沿った値なら採択されるし、そうでなければ変な値として捨てられて別の値候補を探す。これの良いところは目標に近いところから多くサンプリングできることと、パラメータによっては局所解を脱することもできることだ。

ポートフォリオ最適化

　本筋に戻って、次にポートフォリオ最適化である。古くはマーコウィッツの平均分散アプローチに始まり、多様なアセット、現実的な制約(投資総額、最小単位)、リスク許容度、時間軸など様々な条件が付帯する問題となっている。投資家が直面する永遠の課題でもあるが、（許容できる範囲で）最もリスクが低くリターンが高い、つまりリスクあたりのリターンが最も効率的で、各種制約条件を満たすポートフォリオを作成しなければならない。ここでの制約条件は完全に満たす必要があるものもあるし、相対的にベターな選択で足りるものもある。しかしこういった多様な条件の最適化は混合整数計画問題（MIP、MILP）、や混合整数二次計画問題（MIQP）と言われ、解くのが非常に難しいのである。

量子アニーリング

　量子計算は重ね合わせ状態や可逆計算に特徴があるが、イジングモデルを利用した量子アニーリングによる組み合わせ最適化が特にこの分野に有望であるらしい（量子回路による方式もあるようですがちょっと置いときます。）。
　イジングモデルというのは、スピン（量子ビッド）が格子状に並んだモデルのことで、磁場と相互作用という特徴を持つ。例えば隣り合うスピンが両方上向きで、相互作用が正であれば、安定した状態と言える。

　量子アニーリングというのはこのイジングモデルの基底状態をより早く求めることができる方法である。量子アニーリングが古典的なシミュレーテッドアニーリングよりも高速である、という主張は下記Denchev et al. (2016)等で代表される。とはいえ限定された問題設定のみで力を発揮することや、マシン性能など制約が様々ある。実際にポートフォリオ最適化問題についても、量子アニーリングマシンへマッピングできることは原理としては証明されている一方、ハードウェアがまだ追いついていないようである。

What is the Computational Value of Finite Range Tunneling?

過学習

　ポートフォリオ最適化で問題になることの一つが、学習データへのオーバーフィッティングであり、これは機械学習にも通じるものである。つまり、今最適なポートフォリオが将来も最適だとは限らなくて、過去の数字を見て最適なポートフォリオを作って寝かしてしまうといつの間にかテキトーなポートフォリオにアンダーパフォームしてしまうのである。つまり常にアップデートが必要だし、不確実性を織り込んだモデルの作成や、バリデーションの工夫が必要となる。

　例えば下記Jaydip Sen et al(2021)では業種ごといくつかの企業をピックアップし、効率的フロンティアを作って等分ポートフォリオ、最小分散ポートフォリオ、最適ポートフォリオ(シャープレシオに基づく)に投資しているが、ある時点で最適なポートフォリオがいつまでも勝つわけではないことは明らかである。

Precise Stock Price Prediction for Robust Portfolio Design from Selected Sectors of the Indian Stock Market

階層的リスクパリティ

　ポートフォリオの生成方法としては、先に挙げたものの他にリスクパリティポートフォリオがある。これはすべての資産にリスクを均等に配分したポートフォリオを作る方法である。これをさらに進めた方法が階層的リスクパリティポートフォリオ（HRP）であり、リスクパリティポートフォリオの考え方に、機械学習のクラスタリング手法を適用したものである。一世を風靡したあの「ファイナンス機械学習」でも階層的リスクパリティポートフォリオが優れていることが指摘されている（らしい）。

ファイナンス機械学習―金融市場分析を変える機械学習アルゴリズムの理論と実践

　HRPの経緯としては、そもそもマーコウィッツが開発したCLAアルゴリズムという最適化アルゴリズムがあったのだが、大きな共分散行列を利用することによる不都合があった。具体的には逆行列計算による計算量の増加や、ノイズに弱く結果が不安定化するという点である。HRPは相関行列からユークリッド距離を計算し、距離行列に基づいて、ツリー構造でクラスタリングしていく。これにより逆行列計算の手間がなくなり、ランダムなショックに対しての頑健性も得られるようになる、とされている。

現代ポートフォリオ理論

　量子への繋ぎとして、ポートフォリオ最適化における計算量の問題についてもう少し触れておく。ポートフォリオ最適化とは、要は制約付きのリスク最小化（orリターン最大化）であり、各種制約式と目的関数をラグランジュの未定乗数法で解くことがベーシックな解法といえる。これを解くとなると逆行列計算が必要になる。実際にやればそうなるので簡単な紹介とするが、

目的関数　$${ min \displaystyle\frac{1}{2}w^{T}\sum_{}w}$$（ウェイト×共分散行列の最小化）
制約　$${w^{T}R=R_0}$$  （ウェイト×リターンの合計が要求リターンを満たす）
制約　 $${w^{T}1=1}$$　（ウェイトの合計は1）

　というのを解いてやればよい。今回は簡単のためRは無視する（最小分散ポートフォリオを計算する）が、λについて解いたあとにwについて整理すると
$${w=\displaystyle\frac{∑^{-1}1}{1^{T}∑^{-1}1}}$$
となる。（∑は分散共分散行列）
　分母はスカラーになっており、N次元のベクトルが分散共分散行列によって表現されていることがわかる。これを簡単にして２銘柄だけで考えると、
$${w_1= \displaystyle\frac{σ_2^2-σ_{12}}{σ_1^{2}+σ_2^{2}-2σ_{12}}}$$
$${w_2= \displaystyle\frac{σ_1^2-σ_{12}}{σ_1^{2}+σ_2^{2}-2σ_{12}}}$$
となっていることがわかる。さらに簡単にして相関が無いこととすると、
$${w_1= \displaystyle\frac{σ_2^2}{σ_1^{2}+σ_2^{2}}=1-\frac{σ_1^2}{σ_1^{2}+σ_2^{2}}}$$
$${w_2= \displaystyle\frac{σ_1^2}{σ_1^{2}+σ_2^{2}}=1-\frac{σ_2^2}{σ_1^{2}+σ_2^{2}}}$$
となる。リターンの制約がないもとで分散を最小化するポートフォリオのウェイトは分散の比で決まる、というのは直感的にもそうだろう。が、問題はそこではなく、その算出に逆行列が必要という部分なのである。

逆行列計算

　さて、逆行列の計算では$${N^3}$$のオーダーが必要であることが知られている。
　線形代数を勉強した方なら、逆行列を求める際の掃き出し法も必ず勉強していると思う。掃き出し法とは行列に対して行基本変形を行い、単位行列まで変形させることで解を確定させる方法である。N×Nの行列を考えると、
まず１行１列成分が1になるように定数倍する。ここでN回の乗算が必要。
次に1行1列以外の１行目成分がゼロになるよう、2列目から1列の定数倍を引く。ここでも定数倍処理がN回と、引き算がN-1回。
それをすべての行に行うとひっくるめてN回。
　合わせると(N+(N-1)×N)×N=$${N^3}$$の計算量が必要となる。より効率的な計算方法（掃き出し法以外の方法）ももちろんあるが、逆行列の計算というのはリソースを食うものなのである。

階層的リスクパリティ ３つの手順

　そこで階層的リスクパリティの方法がある。アルゴリズムの肝は、全体を上から最下層まで２分割していくことである。これは特に量子を使っているわけではなく、古典的な高速化の方法である。
　階層的リスクパリティの手順は主に３段階で、１．階層的クラスタリング、２．準対角化、３．再帰的二分(これは日本語訳が合っているかわからない)である。
　階層的クラスタリングでは、相関行列からユークリッドの距離を算出し、距離が近い順にクラスター化していく。もちろんこれ以外の距離計算やクラスタリング方法でもよろしい。

　まずはある銘柄iとjの距離を下記のように定義する
$${d_{ij}= \displaystyle\sqrt{\frac{1}{2}(1-ρ_{i,j})}}$$
　これで相関が１である自分自身との距離はゼロ、完全に逆相関である-1との距離は1となり、0～1の間をとる距離に変換される。
　ここで銘柄iとjのユークリッド距離は下記のように書ける。
$${\tilde{d_{ij}}= \displaystyle\sqrt{\sum_{n=1}^{N}(d_{ni}-d_{nj})^2}}$$
　例えば３銘柄h,i,jのみの場合を考えると、この意味は直感的にわかりやすい。
$${\tilde{d_{ij}}= \displaystyle\sqrt{(d_{hi}-d_{hj})^2+(d_{ii}-d_{ij})^2+(d_{ji}-d_{jj})^2}}$$
　右辺の第２項、第３項は$${d_{ij}}$$の定義から各々の距離自体であり、同じ値となる。第１項が他の銘柄との相対的な距離と言える。つまり、銘柄iとjの両方とも銘柄hと相関が高い(低い)のであれば小さい値をとるし、片方だけ高くもう片方だけ低いのであればより大きな値となる。極論的には、hとiがほぼ完全相関、hとjがほぼ逆相関であれば最大となる。
　このユークリッドの距離行列を利用して、最近隣法によるクラスタリングを行う。つまり$${\tilde{d_{ij}}}$$の数値が小さいものからクラスターにしていく。

　次に相関行列の準対角化である。と言っても、これはクラスタリングに沿って大きな相関が対角線周りに並ぶようにソートしなおすというだけである。ここで相関行列が対角行列であれば、相関が無いとみなすことができるようになる。つまり前述の
$${w_1= \displaystyle\frac{σ_2^2}{σ_1^{2}+σ_2^{2}}=1-\frac{σ_1^2}{σ_1^{2}+σ_2^{2}}}$$
$${w_2= \displaystyle\frac{σ_1^2}{σ_1^{2}+σ_2^{2}}=1-\frac{σ_2^2}{σ_1^{2}+σ_2^{2}}}$$
と同じ状況が作られるわけである。これからウェイト付けを行っていくのだが、その際に逆行列計算ではなく対角成分の計算だけで済むようになるのが大きな改善である。

　というわけで最後に再帰的二分である。すでにクラスタリングされている各集合について、それぞれ重みをつけていく。順番としては上から全体を二分割していき(これは資産数で２分割したり、デンドログラムに基づく方法もある)、最終的にクラスターがなくなったら終了である。ここで、二分された部分のウェイトは、各部分の分散から
$${w_i=\displaystyle\frac{diag[V_i]^{-1}}{tr[diag[V_i]^{-1}]}}$$
と計算される。さらに二分された部分の分散は
$${σ_i^2=w_i^{T}V_iw_i}$$
と表すことができる。この分散を使って、以下のようにαを定義する。
$${α_1=1-\displaystyle\frac{V_1}{V_1+V_2}    ;    α_2=1-α_1}$$
　このαを使ってウェイトを再評価する。
　初期状態ではすべてのアセットのウェイトは$${w_i=1 }$$であり、ウェイトの合計は銘柄数Nに等しい。αによる再評価は各銘柄のウェイトを１からだんだん小さくしていき、最終的にウェイトの合計を１にする作業である。この作業をクラスターの分割ができなくなるまで続けると、実際にウェイトの和は1となる。この割り当てのオーダーについては、２分割してαをかける操作をn回するだけという意味で $${O(log_2n)}$$となる。

　申し訳程度に最後に量子の話を付け加えるが、下記では階層的リスクパリティについて、量子アニーラーによって計算可能な形で定式化しており、通常の階層的リスクパリティ方式よりも優れた結果を得られたとしている。将来的に量子アニーラーの性能が向上することで、より早い最適化計算が可能となることに期待している。

Quantum-Inspired Hierarchical Risk Parity | 1QBit

まとめ

　機械学習についてはまとめる前に力尽きたのでここまで。そもそも量子計算周りはほぼ初心者だったので、初等的な部分に時間がかかって深堀するには至らなかった。
　古典的な手法と新しい手法をどちらも理解するというのは大変なことである。

参考文献

Adam Bouland, Wim van Dam, Hamed Joorati, Iordanis Kerenidis, Anupam Prakash "Prospects and challenges of quantum finance", arXiv:2011.06492.
Elham Alipour, Clemens Adolphs, Arman Zaribafiyan, and Maxwell Rounds "Quantum-inspired hierarchical risk parity "  1QB Information Technologies, 2016
Marcos López de Prado "Building Diversified Portfolios that Outperform Out-of-Sample" Journal of Portfolio Management, 2016
Ashley Montanaro "Quantum speedup of Monte Carlo methods" 
arXiv:1504.06987
Lov K. Grover "A fast quantum mechanical algorithm for database search" arXiv:quant-ph/9605043
Gilles Brassard , Peter Hoyer , Michele Mosca , Alain Tapp  "Quantum Amplitude Amplification and Estimation" arXiv:quant-ph/0005055
Ashwin Kumar R S, Geetha Joseph, Kaushik Muthukrishnan, Koushik Tulasi , Praveen Varukolu "Precise Stock Price Prediction for Robust Portfolio Design from Selected Sectors of the Indian Stock Market" arXiv:2201.05570.
Vasil S. Denchev, Sergio Boixo, Sergei V. Isakov, Nan Ding, Ryan Babbush, Vadim Smelyanskiy, John Martinis, Hartmut Neven "What is the Computational Value of Finite Range Tunneling?" arXiv:1512.02206

宇野 隼平  量子コンピュータを用いた高速数値積分 みずほ情報総研技報 Vol.10 2019
森谷 博之　階層的リスクパリティ―理想的なポートフォリオ構築への道― 中央大学 企業研究 第34号 2018
杉田 洋　モンテカルロ積分のためのサンプリング法　生産技術振興協会 生産と技術 第56巻 第1号 2004
https://qiita.com/koke-saka/items/16b20edfeca49e1049cd
https://hudsonthames.org/an-introduction-to-the-hierarchical-risk-parity-algorithm/
https://whyitsso.net/physics/quantum_mechanics/AA_AE_Grover.html
https://en.wikipedia.org/wiki/Monte_Carlo_methods_for_option_pricing