はじめに
TensorFlowは深層学習の分野で利用が広まっていますが、本当は自動微分機能とGPU資源のシームレスな活用を可能にする数値計算言語です。 事実、機械学習の伝統的な手法やあるいは信号処理、微分方程式のソルバーなどの実装も数多く存在しております。 ただし、活用が大きく広まっているという感じではどうやらありません。
そして下記の記事で取り上げた通り、今年末にTensorFlowは大きく生まれ変わる様子で、注力する分野をある程度絞っていく様子が伺えます。（これまで乱立されていたAPIを統合し、contribは完全に廃止）
TensorFlowが公式にサポートしていくと決めた分野の１つに確率的プログラミングであるTensorFlow Probabilityがあります。 この分野にはPyro、PyMC、PyStanなど、既にPythonで利用できる有用なライブラリが存在しますが、今年、PyMCがTensorFlow Probabilityの肩に乗って開発が進められていく方針となり、今後一大勢力になるのではないかと勝手に思っています。
というわけで、TensorFlow Probabilityの下調べをしておこうという記事になります。 今すぐ活用できる便利な知識が載っているわけではないので、あしからず…。
APIの全体像
APIの設計方針のようなものは明確に打ち出されており、下記の図のようになっています。
一応、低レイヤの方からざっくりと説明します。
Layer 0 : TensorFlow
これは単にTensorFlow ProbabilityがTensorFlowを使って実装されますというだけの話です。TensorFlowは既に大きなコミュニティを持っていますし、クラウドの活用からマイコン・スマホでの動作までサポートしています（あるいはしていく方針です）。 このような強力な土台の上で開発が進められていくわけですから、今後「システムとして実運用していく」というフェイズでも利用していくことができるというわけです。
TensowFlowのtf.linalgの部分に含まれていた場所になります。
Layer 1 : Statistical Building Blocks
ここでは確率分布や確率分布の変数変換など、確率計算の基本的なパーツが取り揃えられています。 このレイヤはあくまで、確率・統計に必要な計算を実装してあるというものですので、直接ココのレイヤを必要とすることは応用上はそう多くはないように思います。
「生TensorFlow Probability」にと言ったところでしょうか。自前のサンプリングアルゴリズムを書いたりする場合には必要になってくると思われます。 混合密度推定をニューラルネットワークで行う場合（混合密度ネットワーク）などのも使えそうです（あるいはその手の手法を開発したい場合）。
tf.ditributionsとtf.contrib.distributionsがtfp.distributionsとしてココに集約されています。
Layer 2 : Model Building
学習が可能な形式の分布を提供してくれるため、確率モデルを組む場合はこのレイヤから利用していくことになるでしょう。 確率・統計に必要な諸々の計算を隠蔽しており、データや潜在変数がどんな確率分布に従っていそうなのかを考えることに集中できます。
tfp.edward2として旧Edwardが入っています。確率モデルを作る時の主役になってくるであろう部分です（ただしtfp.distibutionと大きくは変わらないため、google内部として基本的にはtfp.distributionでいくという記述を見かけた気がします）。
tfp.layersでは確率的なニューラルネットワークの層が提供されています。
tfp.trainable_ditributionsでは確率分布を出力するようなニューラルネットワークの構築を支援する役割を担います。
Layer 3 : Inference techniques
MCMCや変分推論の有名ドコロのアルゴリズムを一通り抑えてあります。このLayerのお陰で、確率モデルを構築したら、あとは推論アルゴリズムを実装済みのものから選べばいいということになります。 仮にモデルに依存する推論アルゴリズムを使いたい場合はLayer1、あるいはLayer0から実装することになるかもしれませんが、多くの場合はその必要はないでしょう。
tfp.mcmcではHMCやランダムウォークMH等の有名なサンプリング手法が予め提供されているほか、自前の遷移核を作るための機能も提供されています（自前の遷移核を作るためにはある程度低Layerの機能も理解する必要があるでしょう）。
tfp.viでは変分推論が、tfp.monte_carloではモンテカルロサンプリングが提供されています。
Layer 4 : Pre-built models + inference
有名な確率モデルが予め準備されています。既存のモデルでデータを学習させるという場合にはこのLayerを利用することになると思われますし、たいていのデータ解析者はこの層で十分な気もしますね。
現状はGLMMがtfp.glmによって提供されている状態ですが、今後は（ベイズ）状態空間モデル等の時系列解析にも活用できるモデルが提供されていく予定のようです。Rでいうところのlme4やBSTSに相当する部分です。
高レベルAPItfp.glmの紹介
提供されている一般化線形モデル（GLM）
GLMを用いれば実際のデータに対してある程度使えるモデルの獲得ができます。（GLMに関して勉強したい方は記事の最後に書籍を載せています）
TensorFlow probabilityではGLMを簡単に利用するためのAPIが提供されているので、これの使い方を見てみましょう。
tfp.glmモジュールの中には下記のクラスが含まれています。
●class Bernoulli
下記の$\sigma(\cdot)$はロジスティックシグモイド関数を表しています。 知らん間によく使っているモデルで、確率変数が2値を取る場合に用いられます（ちなみに、偶然にも（？）リンク関数の逆数を活性化関数としてしまい、損失関数を交差エントロピーに選んだパーセプトロンはこのモデルと等価になります）。
$$ \begin{align} y & \sim \mathcal Bern(\mu) \\ \mu & = \sigma(Wx + b) \end{align} $$
●class BernoulliNormalCDF
ベルヌーイ分布のパラメータを標準正規分布のCDF（累積分布関数）で表現する方法。
$$ \begin{align} y & \sim \mathcal Bern(\mu)\\ \mu & = \int_{-\infty}^{Wx+b} {\mathcal N}(z \mid 0, 1)dz \end{align} $$
●class Poisson
離散の正の確率変数を扱う場合のモデルで、誤差の分布がポアソン分布でリンク関数を対数関数としたモデルです。 個数を数えるような場合（品種ごとのスイカの種の個数）などに使われるモデルです。
$$ \begin{align} y & \sim {\mathcal P} (\mu) \\ \mu & = \exp (Wx + b) \end{align} $$
●class PoissonSoftplus
線形予測子を正の値に補正する関数に$\rm softplus$を利用。
$$ \begin{align} y & \sim {\mathcal P} (\mu) \\ \mu & = {\rm softplus} (Wx + b) \end{align} $$
●class Normal
誤差が正規分布、リンク関数を恒等関数とした一般化線形モデル。いわゆる最小二乗法に相当するモデルです。回帰に使われる単層パーセプトロンに相当します。 何も考えずに回帰を実行していると自動的にコレを用いている場合が多いです。ニューラルネットワークを用いた回帰では、線形予測子$\mu = Wx+b$の部分を非線形で多層のモデルにしています。
$$ \begin{align} y & \sim \mathcal N (\mu, 1) \\ \mu & = Wx + b \end{align} $$
●class NormalReciprocal
線形予測子の逆数を平均としたモデルになります。
$$ \begin{align} y & \sim \mathcal N (\mu, 1) \\ \mu & = \frac{1}{Wx + b} \end{align} $$
●class LogNormal
対数正規分布を用いた一般化線形モデルで、対数正規分布とはこの分布に従う確率変数を対数変換した時に、対応する分布が正規分布になるような分布です。 確率変数が正の連続変数で、正方向に裾歪んでいる場合に用いられます。たとえば機械や部品の故障・破損等のモデルに利用されます。
$$ \begin{align} y & \sim \mathcal \Lambda (\log (\mu)- \frac{\log 2}{2}, \sqrt{\log 2}) \\ \mu & = \exp (Wx + b) \end{align} $$
ちなみに「$\mu$の算出で指数取って、確率分布のパラメータを算出するときに対数取ってるのなんで？」って思うかもしれませんが、正の値を強調したいだけで深い意味は無いと思われます。 実際、次のクラスのように$\mu$の算出を行うときに、正の値にする補正を別の関数で実施する場合もあります。
●class LogNormalSoftplus
線形予測子…略
$$ \begin{align} y & \sim \mathcal \Lambda (\log (\mu)- \frac{\log 2}{2}, \sqrt{\log 2}) \\ \mu & = {\rm softplus} (Wx + b) \end{align} $$
●class GammaExp
ガンマ分布も正の値を取る確率変数を扱うときに利用します。また右方向に歪んだ形であり対数正規分布と同様に信頼性を測る場合に用いられるケースがあります。 ガンマ分布は指数型分布族として理論の中で導かれたものであり、特に共役事前分布として登場したりします。
$$ \begin{align} y & \sim \mathcal \Gamma (1, \frac{1}{\mu}) \\ \mu & = \exp (Wx + b) \end{align} $$
●class GammaSoftplus
線形予測子を正の値に...
$$ \begin{align} y & \sim \mathcal \Gamma (1, \frac{1}{\mu}) \\ \mu & = {\rm softplus}(Wx + b) \end{align} $$
使い方
基本的には上記の既存モデルを下記のように簡単に実行できます。
import tensorflow as tf import tensorflow_probability as tfp # データを準備 train_inputs = ... label_inputs = ... # モデルを準備 model = tfp.glm.Bernoulli() # 訓練 coeffs, linear_response, is_converged, num_iter = tfp.glm.fit( model_matrix=train_inputs, response=train_labels, model=model)
coeffsは線形予測子の係数が格納されたtf.Tensorです。
linear_responseは訓練後の線形予測子による、訓練データに対する予測を格納したtf.Tensorです。
is_convergedはconvergence_criteria_fn引数（デフォルト引数はNone）で指定した基準を満たしているかをtf.Boolで返してきます。（使い方勉強中…。
convergence_criteria_fnはtfp.glm.convergence_criteria_small_relative_norm_weights_changeで生成することが可能な模様）
num_iterは学習ステップ数を返します。
最後に
GLMの場合は勾配法（一次の微分を利用する方法）はもちろんのこと、二次の微分を利用するような最適化手法も利用することができます。この際には、TensorFlow Probabilityユーザーは特に難しいことをする必要はなく、tf.glm.fit_sparse()関数を使うだけです。
学習の関数には正則化項やオフセット項などの設定を簡単に行うための引数も準備されており、すぐに誰でも使える実用性を兼ね備えております。GLMに関して勉強したい人は以下の記事の本がオススメです。
