対話型機械学習 (Interactive Machine Learning, IML) について

普段は Computer Graphics や Human-Computer Interaction を勉強しているのですが、最近は機械学習も少し勉強しています。今回は Human-Computer Interaction とも関わりの深いトピックである「対話型機械学習 (Interactive Machine Learning, IML)」について簡単にまとめます。

参考文献

ある機械学習の専門家にこちらの記事をおすすめされました。

Saleema Amershi, Maya Cakmak, W. Bradley Knox, Todd Kulesza:
Power to the People: The Role of Humans in Interactive Machine Learning.
AI Magazine 35(4): 105-120 (2014)
PDF download

この文献は interactive machine learning のサーベイというよりは、この分野を知らない人のための教科書的な導入として書かれているようです。

Interactive Machine Learning とは

そもそも伝統的な機械学習

[...] a powerful tool for transforming data into computational models that can drive user-facing applications.

です。つまり、データがあって、ユーザがいるわけですが、ここで問題となるのはこの "computational models" を構築する作業は肝心のユーザではなく別の "skilled practitioners" (データサイエンティスト?) が行うという点です。したがって、この "practitioners" は実際のユーザと打ち合わせをしながら、適切な特徴量選択やパラメタチューニングを反復的に行っていく必要がありました。


f:id:yuki-koyama:20151010140140p:plain
伝統的な機械学習ではユーザが直接モデル構築を行わないため学習と使用が乖離している。

伝統的な機械学習に対し、モデル構築の過程で直接ユーザと対話することで、より良いユーザエクスペリエンスやより効果的な学習機構を目指すものを interactive machine learning と呼ぶようです。これによって、例えば "practitioners" を挟む間接的なやり方よりも、直接学習と使用のサイクルを回すことができます。


f:id:yuki-koyama:20151010140144p:plain
Interactive Machine Learning ではユーザが直接モデル構築に関与する。

文献の中では、伝統的な機械学習と比較した際の interactive machine learning の特徴として、以下の三点が挙げられています。

  • モデルの更新がより rapid である(ユーザ入力に対し即座にモデルを更新する)
  • モデルの更新がより focused である(モデルの特定の側面だけが更新される)
  • モデルの更新がより incremental である(一回の更新の大きさが小さく、劇的に変わることがない)

これらの条件が満たされることによって、機械学習の専門家でないユーザでも "low-cost trial-and-error or focused experimentation with inputs and outputs" をすることができ、効果的になると述べられています。

事例:画像のセグメンテーションにおける Interactive Machine Learning

Human-computer interaction の文脈において interactive machine learning という言葉が初めて使われたのは以下の論文だそうです:

Jerry Alan Fails and Dan R. Olsen, Jr..
Interactive machine learning.
In Proc. IUI '03. (Best Paper Award)
PDF download

この研究は画像のセグメンテーション(ここでは前景と背景を分ける)において、ユーザが前景(または背景)のピクセル群をブラシストロークによって対話的にマークしていくことで classifier の学習を進めていくという手法が提案されています。

これはちょうど Adobe Photoshop CC 2015 などに搭載されている「クイック選択ツール」の挙動と似ていると思います:

また、この研究で行われたユーザスタディにおいて得られた重要な知見の一つとして、システムの挙動(その時点での classifier によって得られるセグメンテーションの可視化)に応じてユーザが挙動(どこを次にマークするか)を変えた、という点が挙げられます。具体的には、システムがセグメンテーションに失敗した部分をユーザが直す、ということが行われたそうです。この事実は、interactive machine learning では interaction design が学習の効果に大きく寄与するということを示唆していると考えられます。

Interactive Machine Learning に関する学会・会議など

Interactive machine learning の研究は主に CHI, UIST, IUI などの human-computer interaction に関する国際会議で発表されることが多いようです。これは、学習で用いられるモデルそのものよりも、interaction design における知見や human factor の扱いに関する知見が interactive machine learning においては大きな課題であるからだと思われます。

共分散行列を含む多次元のガウス関数の微分

目的:ガウス関数微分

共分散行列 \mathbf{\Sigma} (バンド幅行列) を用いて表された多次元 (多変量) で異方性のガウス関数 (多変量正規分布)微分 (勾配) を考えます。微分したいガウス関数の形は

\displaystyle{
\mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) = \frac{1}{(2\pi)^{\frac{n}{2}}|\mathbf{\Sigma}|^{\frac{1}{2}}}\exp\left( - \frac{1}{2} (\mathbf{x - {\boldsymbol \mu}})^T \mathbf{\Sigma}^{-1} (\mathbf{x - {\boldsymbol \mu}})\right)
}

だとします。ただし \mathbf{\Sigma} \in \mathbb{R}^{n \times n} は正値対称な共分散行列 (バンド幅行列) で、\mathbf{x}, {\boldsymbol \mu} \in \mathbb{R}^n です。これはPRMLなどの機械学習の参考書でも見られるガウス分布の表し方です。

共分散行列を \mathbf{\Sigma}、中心を {\boldsymbol \mu} とするような多次元のガウス関数 \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) のイメージ。

ちなみにバンド幅が \sigma \in \mathbb{R} によって表された \mu \in \mathbb{R} を中心とする 1 次元のガウス関数

\displaystyle{
\mathcal{N}(x; \mu, \sigma) = \frac{1}{\sqrt{2\pi}\sigma} \exp\left(- \frac{(x - \mu)^2}{2\sigma^2} \right)
}

微分

などでも紹介されている通り、

\displaystyle{
\frac{d}{dx}\mathcal{N}(x; \mu, \sigma) = - \frac{x - \mu}{\sqrt{2\pi}\sigma^3} \exp\left(- \frac{(x - \mu)^2}{2\sigma^2} \right)
}

となります。今回はこれの n 次元への拡張、さらにバンド幅が等方性ではなく異方性のものを微分します。

準備:指数部にベクトル・行列が含まれるときのベクトルによる微分

\mathbf{x}n \times 1 ベクトル、\mathbf{A}n \times n 対称行列とするとき、

\displaystyle{
\frac{\partial}{\partial \mathbf{x}} \left( \mathbf{x}^T \mathbf{A} \mathbf{x} \right) = 2 \mathbf{A} \mathbf{x}
}

が成り立ちます。この公式自体は様々な記事やウェブサイトで紹介されているので、ウェブ検索をかけるとすぐに見つかります。

この公式と合成関数の微分の公式 (Chain Rule) を用いることで、

\displaystyle{
\frac{\partial}{\partial \mathbf{x}} \left\{ \exp \left( - \frac{1}{2}\mathbf{x}^T \mathbf{A}\mathbf{x} \right) \right\} \\
= \frac{\partial y}{\partial \mathbf{x}} \cdot \frac{\partial}{\partial y} \exp(y) \:\:\: \left( y := - \frac{1}{2}\mathbf{x}^T\mathbf{A}\mathbf{x} \right) \\
= - \frac{1}{2} \cdot \frac{\partial}{\partial \mathbf{x}} \left( \mathbf{x}^T\mathbf{A}\mathbf{x} \right) \cdot \exp(y) \\
= - \frac{1}{2} \cdot 2 \mathbf{A} \mathbf{x} \cdot \exp \left( - \frac{1}{2}\mathbf{x}^T\mathbf{A}\mathbf{x} \right) \\
= - \exp \left( - \frac{1}{2} \mathbf{x}^T \mathbf{A} \mathbf{x} \right) \mathbf{A} \mathbf{x}
}

という関係式を導くことができます。この式は

では公式として紹介されています。ちなみにベクトルによる微分の Chain Rule は

で証明が紹介されています。

導出:式変形と結果

上記で導いた関係式と再び合成関数の微分の公式 (Chain Rule) を用いると

\displaystyle{
\frac{\partial}{\partial \mathbf{x}} \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \\
= \frac{1}{(2\pi)^{\frac{n}{2}}|\mathbf{\Sigma}|^{\frac{1}{2}}} \cdot \frac{\partial}{\partial \mathbf{x}} \exp\left( - \frac{1}{2} (\mathbf{x - {\boldsymbol \mu}})^T \mathbf{\Sigma}^{-1} (\mathbf{x - {\boldsymbol \mu}})\right) \\
= \frac{1}{(2\pi)^{\frac{n}{2}}|\mathbf{\Sigma}|^{\frac{1}{2}}} \cdot \frac{\partial \mathbf{z}}{\partial \mathbf{x}} \cdot \frac{\partial}{\partial \mathbf{z}} \exp\left( - \frac{1}{2} \mathbf{z}^T \mathbf{\Sigma}^{-1} \mathbf{z} \right) \:\:\: \left( \mathbf{z} := \mathbf{x} - {\boldsymbol \mu} \right) \\
= \frac{1}{(2\pi)^{\frac{n}{2}}|\mathbf{\Sigma}|^{\frac{1}{2}}} \cdot \mathbf{I} \cdot \left\{ - \exp\left( - \frac{1}{2} \mathbf{z}^T \mathbf{\Sigma}^{-1} \mathbf{z} \right) \mathbf{\Sigma}^{-1} \mathbf{z} \right\} \\
= - \frac{1}{(2\pi)^{\frac{n}{2}}|\mathbf{\Sigma}|^{\frac{1}{2}}} \exp\left( - \frac{1}{2} (\mathbf{x} - {\boldsymbol \mu})^T \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \right) \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \\
= - \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu})
}

が得られます。また、冒頭で紹介した 1 次元のガウス関数微分の関係式
\displaystyle{
\frac{d}{dx}\mathcal{N}(x; \mu, \sigma) = - \frac{x - \mu}{\sqrt{2\pi}\sigma^3} \exp\left(- \frac{(x - \mu)^2}{2\sigma^2} \right)
}

も、得られた式において \mathbf{\Sigma}\sigma^2 で置き換え n=1 とすることで、確かに成り立っていることが分かります。

おまけ:二階微分(ヘッセ行列)

上記で得られた結果をさらに用いて、ガウス関数ヘッセ行列 (Hessian Matrix) を求めます。

ベクトル \mathbf{x} を引数とするスカラー関数 f(\mathbf{x}) のヘッセ行列 \mathbf{H}

で説明されている通り、

\displaystyle{
\mathbf{H} = \frac{\partial^2}{\partial \mathbf{x} \partial \mathbf{x}^T} f(\mathbf{x}) = \frac{\partial}{\partial \mathbf{x}} \left( \frac{\partial}{\partial \mathbf{x}} f(\mathbf{x}) \right)^T
}

で計算されます。したがって、ガウス関数のヘッセ行列は
\displaystyle{
\mathbf{H} = \frac{\partial^2}{\partial \mathbf{x} \partial \mathbf{x}^T} \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \\
= \frac{\partial}{\partial \mathbf{x}} \left\{ \frac{\partial}{\partial \mathbf{x}} \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \right\}^T \\
= \frac{\partial}{\partial \mathbf{x}} \left\{ - \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \right\}^T \\
= - \left\{ \frac{\partial}{\partial \mathbf{x}} \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \right\} \left( \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \right)^T - \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \left\{ \frac{\partial}{\partial \mathbf{x}} \left( \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \right)^T \right\} \\
= - \left( - \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \right) \left( \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) \right)^T - \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \mathbf{\Sigma}^{-1} \\
= \mathcal{N}(\mathbf{x}; {\boldsymbol \mu}, \mathbf{\Sigma}) \left\{ \mathbf{\Sigma}^{-1} (\mathbf{x} - {\boldsymbol \mu}) (\mathbf{x} - {\boldsymbol \mu})^T \mathbf{\Sigma}^{-1} - \mathbf{\Sigma}^{-1} \right\} 
}

となります。途中、(\mathbf{\Sigma}^{-1})^T = \mathbf{\Sigma}^{-1} が成り立つことと、ベクトルによる微分に関する積の微分の関係式を用いました。一階微分を計算したときと同様に、\mathbf{\Sigma}\sigma^2 で置き換え n=1 とすることで

で示されている結果などと矛盾しないことが分かります。

最後に

この記事には間違いがある可能性があります。ここに書かれている数式を用いる場合には必ず自身で再計算して確認してからにしてください。また、もし間違いを見つけた場合には教えてください。

Multidimensional Scaling (多次元尺度構成法, MDS) の計算方法と実装

Multidimensional Scaling (MDS) について少し勉強したのでメモしておきます。

Multidimensional Scaling とは

はじめに

日本語では「多次元尺度構成法」と呼ばれる統計テクニックの一つです。英語版ウィキペディアの記事が詳細です。

多次元尺度構成法 - Wikipedia
Multidimensional scaling - Wikipedia

MDS には色々な派生がありますが、ここでは Classical MDS (古典的多次元尺度構成法) と呼ばれる最も基礎的で古典的なものについてのみ記述します。

どういうときに使えるか

いくつか要素があり、要素間の距離がなんらかの形で定義できる場合に、その要素を多次元空間に配置する場合に使えます。

具体例

例えば、いま目の前に日本酒が10銘柄(鳳凰美田、花陽浴、...)あるとし、これらを一つずつ味見し「どの銘柄とどの銘柄がどのぐらい似ているか(似ていないか)」という指標によって、以下のように任意の2銘柄間の距離が定義できるとします。


要素間の距離の定義のイメージ*1d(\cdot, \cdot)は要素間の距離を表すとする。
d(農口, 村祐) = 2 d(農口, 鳳凰美田) = 4 d(農口, 花陽浴) = 7 ...
d(十四代, 村祐) = 2 d(十四代, 呼友) = 2 d(十四代, 花陽浴) = 5 ...
... ... ... ...


MDS を用いると、これらの距離の情報から、例えば2次元の座標系に銘柄の「マップ」を可視化することができたりします。ここで、各点の座標は MDS によって計算され、その際に各点間の距離は先ほど定義した距離が反映されます。



MDS の使用イメージ*2。MDS によって各点の座標が計算される。

入力と出力

MDS の入力は要素間の距離と、対象とする空間の次元です。一般的には、距離は距離行列 (distance matrix) として表現しておきます。要素数n であった場合、D_{i,j}i 番目と j 番目の要素の距離の二乗*3であるような正方行列 \mathbf{D} \in \mathbb{R}^{n \times n} が入力です。このとき、距離の対称性と非退化性*4から、 \mathbf{D} は対角成分が 0 であるような対称行列となります。

MDS の出力は、各要素の座標 (\mathbf{x}_1, \ldots, \mathbf{x}_n) です。ここで、指定された空間の次元が m だとすると、\mathbf{x}_i \in \mathbb{R}^m です。これらは i 番目の行を \mathbf{x}_i とする行列 \mathbf{X} \in \mathbb{R}^{m \times n} として表現したりします*5

次元圧縮 (Low-Dimensional Embedding)

高次元の点群データから距離行列 \mathbf{D} を計算し、元の次元よりも小さな次元 m を指定し MDS を実行することで、次元圧縮 (Low-Dimensional Embedding) を実現することができたりもします*6

非線形な次元圧縮手法としてよく使われる ISOMAP も MDS を拡張したものです。

MDS の計算方法

行列 \mathbf{X} を求めるためには、まず

\mathbf{K} = \mathbf{X}^T\mathbf{X} \in \mathbb{R}^{n \times n}

を満たす行列 \mathbf{K}*7*8を、入力として与えられた距離行列 \mathbf{D} から計算できるとします。すると、この行列 \mathbf{K}固有値分解 (eigenvalue decomposition) によって
\mathbf{K} = \mathbf{V}\mathbf{\Lambda}\mathbf{V}^T

と分解することができます*9*10。ここで \mathbf{V}固有ベクトルを並べた行列、\mathbf{\Lambda}固有値を対角成分に持つ行列です。これらを用いて、
\mathbf{X} = \mathrm{sqrt}(\mathbf{\Lambda})\mathbf{V}^T

とすることで、所望の \mathbf{X} を計算することができます*11。ここで \mathrm{sqrt}(\cdot) は要素毎に累乗根をとることを意味しています*12。また、このままでは空間の次元は n ですが、これをより小さい次元 m に削減するには、固有値の上位 m 個と、それに対応する固有ベクトルのみを使って前式を計算することで実現できます*13

与えられた距離行列 \mathbf{D} からカーネル行列 \mathbf{K} を計算する方法ですが、これは

\mathbf{K} = - \frac{1}{2} \mathbf{H} \mathbf{D} \mathbf{H}

によって算出するようです。ただし
\mathbf{H} = \mathbf{I} - \frac{1}{n}\mathbf{1}\mathbf{1}^T

は中心化行列*14と呼ばれ、得られる結果の原点を重心に移す効果があるそうです*15*16

MDS の実装

C++ によって MDS を実装した例を以下に示します。行列演算ライブラリとして Eigen を用いています。

#include <Eigen/Core>
#include <vector>
#include <utility>
using namespace std;
using namespace Eigen;

// This function computes low-dimensional embedding by using multi-dimensional scaling (MDS)
// - Input:  A distance (dissimilarity) matrix and a target dimension for embedding
// - Output: A coordinate matrix whose i-th column corresponds to the embedded coordinates of the i-th entry
MatrixXd computeMDS(const MatrixXd &D, unsigned dim)
{
    assert(D.rows() == D.cols());
    assert(D.rows() >= dim);
    const unsigned n = D.rows();
    const MatrixXd H = MatrixXd::Identity(n, n) - (1.0 / static_cast<double>(n)) * VectorXd::Ones(n) * VectorXd::Ones(n).transpose();
    const MatrixXd K = - 0.5 * H * D * H;
    const EigenSolver<MatrixXd> solver(K);
    VectorXd S = solver.eigenvalues().real();
    MatrixXd V = solver.eigenvectors().real();
    extractNLargestEigen(dim, S, V);
    const MatrixXd X = DiagonalMatrix<double, Dynamic>(sqrt(S)) * V.transpose();
    return X;
}

以下にテストコード等を含めた git repository 公開しました:
github.com

*1:これらの銘柄は実在するものですが、その距離として与えられている値には全く根拠はなく、極めて乱数に近いものであると考えてください。

*2:この図はあくまでイメージを掴むためのもので、実際のデータから計算されたものではありません。

*3:この部分について明示していない文献も多く私は途中まで勘違いしていたのですが、距離行列は距離の二乗を要素として持つようです。

*4:距離空間 - Wikipedia

*5:文献によっては i 番目の\mathbf{x}_{i}^T とする行列 \mathbf{X} \in \mathbb{R}^{n \times m} として表現したりすることもあります。

*6:この操作は Principal Component Analysis (PCA) と同値?

*7:グラム行列PRML 6.2章)

*8:文献によって \mathbf{B} と表記されることもあります。

*9:\mathbf{K} は対称行列であるため、特異値分解 (SVD) を用いても同様の計算ができます。

*10:\mathbf{D}ユークリッド距離に基づいているとき、まとそのようなときに限り \mathbf{K} は半正定値行列となる?

*11:\mathbf{X}^T\mathbf{X} = (\mathrm{sqrt}(\mathbf{\Lambda})\mathbf{V}^T)^T (\mathrm{sqrt}(\mathbf{\Lambda})\mathbf{V}^T) = \mathbf{V} \mathrm{sqrt}(\mathbf{\Lambda}) \mathrm{sqrt}(\mathbf{\Lambda}) \mathbf{V}^T = \mathbf{V}\mathbf{\Lambda}\mathbf{V}^T

*12:ただし、\mathbf{K} が半正定値行列でない場合には固有値に負の値が現れます。その場合はその固有値を 0 だとみなして以降の計算を行います。これはある種の近似を行っていると考えられるそうです

*13:そのようにして良い根拠はこの資料で説明されています。

*14:文献によって \mathbf{G} と表記されることもあります。

*15:このあたりはよく理解できていませんが、距離を一定に保った点群には剛体移動(平行移動と回転移動)の自由度があるので、そのうち平行移動については重心が原点にくるように制約をかけるということだと思います。

*16:この変換は「ダブル・センタリング変換」と呼ばれる?

日本人研究者によるSIGGRAPH論文

概要

日本人研究者によるSIGGRAPH論文 (2014--2015年分) の情報を以下のページにまとめました。せっかくなので公開します:

SIGGRAPH Papers by Japanese Researchers

YouTube動画の埋め込みや著者の個人ページへのリンクも含んでおり、日本人研究者の活躍ぶりをまとめて把握する上で役に立つかもしれません。

SIGGRAPH論文

当該ページの冒頭にも書きましたが、SIGGRAPHはコンピュータグラフィクスの分野のトップカンファレンスとして知られ、またSIGGRAPHに論文が採択されることは非常に難しいことで有名です。SIGGRAPHに論文が採択されると、分野の研究者だけでなくメディアなどを通じて様々な人の目に触れることになります。

ちなみに、SIGGRAPH AsiaACM Transactions on Graphicsも"uniform high quality"であると明示的に宣言されており、まとめて「SIGGRAPH論文」と表現することも多いです。SIGGRAPHSIGGRAPH Asiaで発表された論文は自動的にACM Transactions on Graphicsに掲載され、またACM Transactions on Graphicsに掲載された論文はSIGGRAPHSIGGRAPH Asiaで発表する権利が与えられます。

また、ここでの論文とは正式にはTechnical Papersと表記されるものであり、Emerging TechnologiesやPostersとは別のものです。

なお、Eigenfactor scoreによると、ACM Transactions on Graphicsはコンピュータグラフィクス分野だけでなくComputer Science, Software Engineering分野においてもトップジャーナルなのだそうです。

日本人研究者

例年、SIGGRAPHに採択される日本人研究者による論文は片手で数えられる程度、0件であることも珍しいことではないようです。例えば、直近だとSIGGRAPH Asia 2014では日本人研究者による論文は採択されなかったようです。ただし今年のSIGGRAPHには非常に多くの論文が日本から採択されたようです。

なお、日本人がSIGGRAPHに論文が採択されると、コンピュータグラフィクスの国内シンポジウムの一つVisual Computing/グラフィクスとCAD合同シンポジウム (VC/GCAD) にて招待講演を依頼されるという制度があります。しかしながら、ACM Transactions on GraphicsやSIGGRAPH Asiaの論文は取りこぼしがあるようで、完全ではありません。

そこで、日本人研究者(海外に在住する人や、外国籍でも日本に深い縁のある人を含む)によるSIGGRAPH論文を、YouTube動画や著者の個人ページへのリンクも含めて集めてみました。収集手法は目grep法を用いたため、精度に問題があるかもしれません。情報に欠けを見つけた場合はご連絡ください。当該ページは、今後順次更新していくかもしれないし、現状のまま放置するかもしれません。