三角図(三角グラフ)は、
$ x + y + z = 1 , x ≧ 0, y ≧ 0, z ≧ 0 $
を満たすような $ x, y , z $ を2次元の三角形で表す図である。
三角図の例や見方については一旦おいておき、python-ternary を使って、関数を用いた三角図のヒートマップを実装してみようと思う。
!pip install python-ternary
ここで、今回取り上げるディリクレ分布についての紹介を。
ディリクレ分布は、各次元が正であり総和が1となるような $ K $ 次元ベクトルが従う分布である。
数式は以下。
$$
Dir(\textbf{π}|\textbf{α}) = C_D(\textbf{α}) \prod_{k=1}^{K} {\pi_k}^{\alpha_k - 1}
$$
ここで、$ \bf{α} $ は $ K $ 次元ベクトルのパラメータであり、各要素 $ \alpha_k $ は正の実数値である。
定数項は、
$$
C_D(\textbf{α}) = \frac{Γ(\sum_{k=1}^{K} \alpha_k)}{\prod_{k=1}^{K} Γ(\alpha_k)}
$$
$ \pi_k $ ($ \textbf{π} $ のk番目の要素)の期待値は、
$$
\frac{\alpha_k}{\sum_{i=1}^{K} \alpha_i}
$$
である。
$ K = 3 $ のとき、ディリクレ分布は3つの次元それぞれが正であり総和が1となる3次元ベクトルを生成してくれる。
python-ternaryを使う上ではディリクレ分布についてはあまり深く知る必要はない。
とりあえずディリクレ分布を実装してみよう。
先に言っておくと、下記の関数は引数を二つとることに注意。
import numpy as np
from scipy.special import gamma
def Dirichlet(pi_vec, alpha_vec=[1, 2, 3]):
# πベクトル、αベクトルをnumpyのndarrayとして扱う
pi_vec = np.array(pi_vec)
alpha_vec = np.array(alpha_vec)
# Πは、各インデックスに対する値を格納したnp.ndarrayにnp.prod(・)を適用することで実装
main = np.prod(pi_vec ** (alpha_vec - 1))
C_D = gamma(np.sum(alpha_vec)) / np.prod(gamma(alpha_vec))
return main * C_D
上記のコードでは、$ K = 3 $ を想定して引数のデフォルト値を用意した。
早速やってみよう。
全体のコードから断片を拾い上げ、最後にまとめたコードを紹介する。
まずは必要なライブラリのインポート。python-ternaryのインポートは import ternary でいい。
$ \textbf{α} $ は固定で決める。
from functools import partial
import matplotlib.pyplot as plt
import matplotlib
# 三角グラフ描画用のライブラリをインポート
import ternary
alpha = (2, 4, 6)
ternary.figureで土台をつくる。これはなんだかMatplotlibと似たような感じ。
ternary.figureの引数にscaleがあるが、これは分割数を意味している。
戻り値でfigureとtaxを取得しているが、このtaxは、Matplotlibのaxと似たような感じで使える。
# 分割数。
scale = 100
figure, tax = ternary.figure(scale=scale)
# 表示させる画像のサイズを指定
figure.set_size_inches(10, 8)
今回はtaxのheatmapfというメソッドを利用する。関数を使ってヒートマップを描画しよう、という感じ。
このheatmapfのfunc引数に渡された関数には、順次三つ組みのタプル ―例えば、(0.15, 0.28, 0.57)ー が引き渡され、各座標における値が計算される。
今回は自作のDirichlet関数を使用しているが、Dirichlet関数には引数が二つある。引数alpha_vecは固定したいので、Python標準ライブラリfunctoolsのpartialを用いて、部分適用された関数を作成する。
heatmapfのboundary引数は、境界での値を計算するかどうか。今回取り上げたディリクレ分布が返すベクトル $ \textbf{π} $ は各次元共に正である条件があるので、境界の座標 ―例えば、(0.2, 0.8, 0.0)― では計算しないでほしいのでFalseにした。
heatmapfのstyle引数は、"triangular" "dual-triangular" "hexagonal" のいずれかである。試してみてほしい。
gridlinesは三角図の内部に描かれる直線のことで、各辺に平行なこれが何本も描かれることで、特定変数についての等高線を得ることもできる。
# funcには、総和1のベクトルがscaleの分割数に合わせて渡される。
# funcには関数を渡す必要がある。今回は引数pi_vec, alpha_vecのうち、alpha_vecのみ先に指定したうえで関数として渡したいので、Python標準ライブラリfunctoolsのpartialを用いて、部分適用された関数を作り、func引数
tax.heatmapf(func=partial(Dirichlet, alpha_vec=alpha), boundary=False, style="triangular", cmap='rainbow')
tax.boundary(linewidth=2.0)
tax.gridlines(color="black", multiple=10, linewidth=1)
ここは、ディリクレ分布の期待値を計算しているところ。
本編とは少しずれるかも。
alpha = np.array(alpha)
alpha_sum = alpha.sum()
# tupleに戻すのは、このまま表示させるとそれっぽく見えるから。
Exp = tuple([np.around((alpha_i / alpha_sum), decimals=2) for alpha_i in alpha])
ここは割と見たまんま。
下、右、左に軸のラベルを表示させる。
offsetは、三角図の辺からの距離のことをいっているみたい。
# Set Axis labels and Title
fontsize = 25
offset = 0.14
tax.set_title(title=f'α = {tuple(alpha)}, Exp = {Exp}', loc='center', fontsize=fontsize)
tax.bottom_axis_label("$\\alpha_1$", fontsize=fontsize, offset=offset)
tax.right_axis_label("$\\alpha_2$", fontsize=fontsize, offset=offset)
tax.left_axis_label("$ \\alpha_3 $", fontsize=fontsize, offset=offset)
tax.ticksで設定しているのは、数字の表示。
今回は、実装の初めにscale=100とした。
ここで初めて言う事実がひとつある。
heatmapfで三角図を表示させた場合、各辺ごとに1～scaleまでの整数が表示される。今回であれば、1～100の整数だ。
だが、あくまでもこれはscale(100)分割したうちの何番目かを示しているだけで、あくまでも座標値の総和は1である。
辺ごとに表示される整数を1/scaleすれば座標値が得られるという認識でいてくれればいい。
tax.ticksのmultiple引数を設定することで、値をとびとびにできる。
1～100までの100個の数字を表示させるのはだれしも嫌だろう。
multipleに設定した数だけ値がとびとびになる。
tax.clear_matplotlib_ticks() は公式 (https://github.com/marcharper/python-ternary) にて、MatplotlibのデフォルトのAxexを取り除く、との目的で使われている。今回は外しても何も変わらなかった。
tax.get_axes().axis('off') をしないと、三角図を囲う四角い枠が表示されることになり、ダサい。
# multiple、デフォルトでは1みたい。
tax.ticks(axis='lbr', linewidth=1, multiple=10)
tax.clear_matplotlib_ticks()
tax.get_axes().axis('off')
tax.show()
というわけで、上記のコードをまとめたのが以下である。
from functools import partial
import matplotlib.pyplot as plt
import matplotlib
# 三角グラフ描画用のライブラリをインポート
import ternary
alpha = (2, 4, 6)
# 分割数。
scale = 100
# taxはmatplotlibのaxのようなものだと思ってもらえればたぶん大丈夫
figure, tax = ternary.figure(scale=scale)
# 表示させる画像のサイズを指定
figure.set_size_inches(10, 8)
# funcには、総和1のベクトルがscaleの分割数に合わせて渡される。
# funcには関数を渡す必要がある。今回は引数pi_vec, alpha_vecのうち、alpha_vecのみ先に指定したうえで関数として渡したいので、Python標準ライブラリfunctoolsのpartialを用いて、部分適用された関数を作り、func引数
tax.heatmapf(func=partial(Dirichlet, alpha_vec=alpha), boundary=False, style="triangular", cmap='rainbow')
tax.boundary(linewidth=2.0)
tax.gridlines(color="black", multiple=10, linewidth=1)
#tax.gridlines(color="black", multiple=100, linewidth=1)
# (r,p)の組と、期待値('Exp')を表示
alpha = np.array(alpha)
# np.array()に対しても、内法表記使えるみたいだね。便利。
alpha_sum = alpha.sum()
Exp = tuple([np.around((alpha_i / alpha_sum), decimals=2) for alpha_i in alpha])
# 結局αはtupleに戻すんだけどね。表示用に。
# Set Axis labels and Title
fontsize = 25
offset = 0.14
tax.set_title(title=f'α = {tuple(alpha)}, Exp = {Exp}', loc='center', fontsize=fontsize)
tax.bottom_axis_label("$\\alpha_1$", fontsize=fontsize, offset=offset)
tax.right_axis_label("$\\alpha_2$", fontsize=fontsize, offset=offset)
tax.left_axis_label("$ \\alpha_3 $", fontsize=fontsize, offset=offset)
# multiple、デフォルトでは1みたい。
tax.ticks(axis='lbr', linewidth=1, multiple=10)
#tax.clear_matplotlib_ticks()
tax.get_axes().axis('off')
tax.show()
python-ternaryでは、3次元の座標データから三角図の散布図を作成することもできるし、色々なメソッドが用意されているらしい。
興味のある方はぜひ使ってみてね。