AILab Creative Researchチームの富樫です。 このブログでは先月末にarxivに投稿された“Item Recommendation from Implicit Feedback”[1]という論文を軸に紹介しつつ、 周辺分野の話題について議論したいと思います。 この論文はitem推薦というタスクにおける手法の各種パラダイムの概観をコンパクトに解説した教科書的内容になっています。 著者はBayesian Personalized Ranking (BPR)[2]を開発したGoogle Research所属のSteffen Rendle氏であり、 長年この分野を開拓してきた権威の一人です。

元論文の内容は元論文を読めばわかることですし、 蛇足かもしれませんが、最近の研究との関連性や議論、個人的な感想などを示すことで、このブログが元論文に対する補足資料のようになることを目指したいと思います。 元論文の参考文献は厳選された機械学習の論文が多いため、情報推薦の既存研究との関連なども示せればと思います。
想定読者は、情報検索・推薦の論文を日頃チェックしている研究者、学生、データサイエンティスト/機械学習エンジニア、の方です。

また、独自の補足を行うため、自分の理解不足により元論文と異なる主張を書く可能性があります。そのため、元論文の内容を正確に知りたい方は先に元論文を通読することをお薦めします。 ブログの内容についてのご指摘、ご質問などあればぜひ遠慮なく富樫までお伝えください。

問題設定と論文の主旨

2.1章冒頭を参考にimplicit feedbackを用いたitem推薦の設定を説明します。
まず使用可能なデータとして、なんらかの固有の嗜好を持ったあるuser \(c \in C\)に対して、 推薦候補となるitemの集合\(I\)の内のどれを過去に好んだか、という情報が得られます。 この過去のuser毎の履歴をimplicit feedbackと呼びます。クリックや視聴履歴のようなものを想定した概念であり、 観測として\(c \in C\)が\(i \in I\)を好んだということを表す組 \((c, i) \in S\)が得られることになります。今、\(c\)それぞれが過去に好んだitem集合を\(I_{c}\)とします。 \(S\)をもとに、\(c\)それぞれに対して、\(I \setminus I_{c}\)を嗜好に沿う順序で並べ替える、というのがitem推薦の目標です。指標としてはcontext毎のRecall@nやnDCG@nの平均を想定することが多いです。
また、Eコマースのサイトで全商品をクリックするユーザが極めて少数であるように、一般のアプリケーションにおいて、大部分のcontext-item対にはimplicit feedbackは観測されません。 つまり、あるcontextに対してimplicit feedbackが観測される回数は一般に限られており(\(|I_c| \ll |I|\))、 \(|S|\ll|C \times I|\)となります。これを、implicit feedbackのsparsityと言います。

ちなみに、\(c\)は協調フィルタリング的にuserと捉えても良いですし、なんらかのクエリなりセッション情報なりを追加で考慮することも可能ですが、 元論文のスコープにおいては、一般化してcontextと呼んでいます。

観測された\(S\)をもとにして、ランキングを推定するために、\(S\)をもとに組\((c, i)\)に対する嗜好スコアのモデル\( \hat{y}(c,i)=\hat{y}(i|c),\,\,\, \hat{y}\colon C \times I \mapsto \mathbb{R} \)を作ることになります。以下、 この\(\hat{y}\)はモデルのパラメータ\(\theta\)によってパラメトライズされていることにします。
元論文の主旨は、モデル\(\hat{y}\)の学習に対して使える主流なパラダイムをいくつか紹介しつつ、現実的なアプリケーションで耐えうる 効率的な学習アルゴリズムの概観を示すことです。

2.4章では、implicit feedbackからの学習という設定の固有の問題について触れています。 Implicit feedbackでは、観測されたある組\((c, i) \in S\)は、\(c\)に対して\(i\)が実際に好まれた(e.g., クリックされた)ということを示すため、 ある\(c\)に対して\(i’ \in I\)が好まれなかった、というデータが得られません。これをone-class problemと呼びます。 また、rerankingのために使われるランキング学習との違いにも触れられています。 特に重要な違いは、rerankingの場合はランキングする対象は全対象のごく一部であり、item推薦のように全対象を並び替える必要はないという点です。

ここでimplicit feedbackを用いたitem推薦の特徴をまとめると、以下のようになります。

• 大量のitem (\(I \setminus I_{c}\))をcontext \(c \in C\)毎に並び替える必要がある
• 得られるimplicit feedbackは好まれたcontext-itemの組しか含んでいない

つまり、item推薦とは呼ばれていますが、これは大規模な候補集合を対象としたランキング学習のタスク一般を扱う話題だと考えることができます。大規模なcontext-item対のそれぞれに対して高価なアノテーションは得られないので、応用のためにはclickなどの自動的に集まるimplicit feedbackを活用する必要があり、それがone-class problemという固有の問題を引き起こしている、という構図です。

内積モデル

2.2章の話題はRendle氏がいろいろな論文で主張していることで、重要な話題なので紹介したいと思います。 \(\hat{y}(i|c)\)のモデルとして、以下のような内積による表現が可能なモデルを内積モデル(Dot product model)と呼びます。 \[ \hat{y}(i|c) = \langle \phi(c), \psi(i) \rangle \] ここで、\(\phi \colon C \mapsto \mathbb{R}^{d}\)と\(\psi \colon I \mapsto \mathbb{R}^{d}\)はそれぞれ contextとitemの\(d\)次元の特徴へのマッピングです。 協調フィルタリングであれば、Matrix Factorization (MF)やGraph Neural Network (GNN)による実装が主にこの形で表される代表的なモデルです。

モデルがこの形式で表されるということを利用して、色々な高速化のアルゴリズムが考えることが可能になります。
最終的な予測ランキングを作る際には全ての推薦候補をソートする必要があるため、 上記の形で表せないNeural Collaborative Filtering[3]などのmulti-layer perceptron型のモデルは、 test時の全件予測において\(\mathcal{O}(|C||I|)\)の計算コストが生じるため、大規模なアプリケーションでは実行不可能になりえます。 このあたりの議論を性能比較も含めて行った論文がRendle氏からRecSys’20で発表されている[4]ので、興味があればチェックしてみてください。 元論文の6.2章でもこの話題について少し触れられています。

内積モデルを用いることで、item数 \(|I|\)が大きいときのランキングを予測する計算コストを approximate nearest neighbors (ANN)などの確立された技術を使って大幅に削減できます。これは特にオンライン予測が必要なケースにおいてもメリットになります。
単純な協調フィルタリングでは、\(C\), \(I\)はtrain時とtest時に固定なので、事前に予測ランキングを作ってしまえば再計算の必要性があることは少ないですが、context側の特徴マッピング\(\phi(\cdot)\)がリアルタイムに更新される状況や\(c\)が画像やテキストなどの特徴を含んでいる場合などでは、大量のランキングを保持するよりANNなどを使うのが現実的です。リアルタイムに\(\phi(\cdot)\)が変わる状況に関しては6.3章でDynamic user modelという話題で少し触れられています。

主流な損失

Pointwise損失

最初のパラダイムであるpointwise損失について3.1章を参考に説明したいと思います。 pointwise手法では、観測した組に対して\(y=+1\)、観測していない組に対して\(y=0\) or \(y=-1\)という\((c,i) \in C \times I\)それぞれに対するラベルを表す変数を導入します。 これによって単純な二値回帰あるいは分類に帰着することができます。 損失関数は色々なものが使われていますが、特に広く使われているものには二乗誤差(回帰)などがあります。 \[ \begin{align} l(\hat{y}, y) = (\hat{y}- y)^2 \end{align} \] この損失関数を使う場合、未観測itemに対して\(y=0\)として定義します。 損失の和は以下のように書けます \[ L(\theta, S) = \sum_{c \in C} \sum_{i \in I} \alpha (c,i) l(\hat{y}(i|c), y(c,i)) \] ここで、\(\alpha\colon C \times I \mapsto \mathbb{R}\)は組に対する重み関数です。 このアプローチは\(p(y=1|c,i)=\sigma(\hat{y}(i|c))\)なるモデルによってクラス確率を推定したいという雰囲気です。

重み関数は性能と学習効率において、極めて重要な要素です。元論文では 未観測itemに対して\(y=0\)として全対\(C \times I\)上の損失和\(L(\theta, S)\)を最小化するのが単純なpointwise手法ですが、 \(c\)に対して\(I\setminus I_c\)を並べ替えるのが目的なので、未観測itemを全て完全な負例として学習するのは 意味のない予測結果になります。また、one-class problemによって、未観測itemの中には潜在的な正例と負例が含まれています。 そのため、重み関数\(\alpha(c,i)\)は主に未観測itemを重み付けるために使われており、正例\((c,i)\in S\)に対しては\(\alpha(c,i)=1\)とするのが一般的です。 この重み関数の選択によって性能が（というか推定したいものが）大きく変わります。 unbiased recommendationの文脈などでも対応するものが提案されています[5]。 学習効率においては、alternating least squares (ALS)ベースであると一様な重み以外は計算量がcubicなどになる[6]ので、 それを高速化する手法[7]や、SGDを使った場合のEMアルゴリズムの速度の問題を軽減する手法[8]なども提案されています。

それでも全データを使った損失和の要素数は\(\mathcal{O}(|C \times I|)\)となり、大規模な設定では計算コストを削減するための工夫が必要になります。

Pairwise損失

top-\(n\)推薦において長らく主流となっているのは、pairwise手法です。 Recall@\(n\)やnDCG@\(n\)などの指標計算において、\(\hat{y}\)の順序しか最終的に使わないので、順序のみ推定すれば十分です。そこで、ある2つのitemを比べてどちらがより好まれるかを 分類するタスクとして定式化するアプローチです。 損失関数はlogistic回帰型のものが使われることが多いように思います。 \[ \begin{align} l(\hat{y}, y) = \ln(1+\exp(-y\hat{y})) \end{align} \]

\[ L(\theta, S) = \sum_{(c,i) \in S} \sum_{j \in I} \alpha (c,i,j) l(\hat{y}(i|c)-\hat{y}(j|c), 1) \] ちなみに、損失\(l(\hat{y}(i|c)-\hat{y}(j|c), 1)\)はスコアの差\(\hat{y}(j|c)-\hat{y}(i|c)\)に対するsoftplus関数になっているので、\(c\)に対する正例\(i\)と負例\(j\)の順序に関する0-1損失の微分できるsurrogateになっています。

\[ \begin{align} \unicode{x1D7D9}\{\hat{y}(j|c) > \hat{y}(i|c)\} &= \unicode{x1D7D9}\{\hat{y}(j|c) – \hat{y}(i|c) > 0 \}\\ &\approx \ln(1+\exp(\hat{y}(i|c)-\hat{y}(j|c))) = l(\hat{y}(i|c)-\hat{y}(j|c), 1) \end{align} \]

0-1損失は順序が間違っている回数を表しています。つまり、この和や期待値はAUCのようなノンパラメトリックな指標をより直接表現しています。 損失和の要素数は\(\mathcal{O}(|S||I|) \supseteq \mathcal{O}(|C||I|)\)となっているのでpointwiseよりも更に計算コストが高い傾向があります。 \(|S|\)はユーザ数x平均click数などになり、ほぼ実行不能な\(|C \times I|\)にcontext平均の観測item数が掛け算されるため、工夫しない限り大規模なデータでの学習は極めて高コストになります。 また、損失のダイナミクスが正例より上位にある負例を学習ステップで得られるかどうかに強く依存するため、特に学習後期において学習が著しく停滞する問題が一般に知られています。

Softmax損失

より素直な定式化として、観測されたimplicit feedbackの密度推定をする、という方法があります。 Softmax損失によるアプローチでは、\(p(i|c)\)を直接モデリングします。 notationがわかりにくいですが、密度推定として定式化した場合にはラベル変数\(y\)を考える必要性がないので、ここにおける\(p(i|c)\)は、 pointwise損失において行った定式化における\(p(i|c,y=+1)\)に相当します。

まずはパラメトリックな密度関数のモデルとして以下のsoftmaxモデルを導入します。 \[ p(i|c) = \frac{\exp(\nu \hat{y}(i|c))}{\sum_{j \in I}\exp(\nu \hat{y}(i|c))} \propto \exp(\nu \hat{y}(i|c)) \] ここで、\(\nu\)はsoftmaxのピーキーさを調整する温度パラメータです。 導入はしましたが、MFなどを使った場合に結局\(\hat{y}\)は埋め込みのノルムによってこれを吸収するので以降無視します。 これを用いて、密度関数の最尤推定を考えると、負の対数尤度として以下の損失和を得ます。 \[ \begin{align} L(\theta, S) &= – \sum_{(c,i)\in S} \ln p(i|c) \\ &= – \sum_{(c,i)\in S} \ln \frac{\exp(\hat{y}(i|c))}{\sum_{j \in I}\exp({\hat{y}(j|c)})} \\ &= – \sum_{(c,i)\in S} \left[ \hat{y}(i|c) – \ln \sum_{j \in I}\exp({\hat{y}(j|c)})\right] \\ &= – \sum_{c \in C}\left[\sum_{i \in I_c} \hat{y}(i|c) – \ln \sum_{j \in I}\exp({\hat{y}(j|c)})\right] \end{align} \] よって、単純に実装した場合(3番目の右辺)の損失の要素数は\(\mathcal{O}(|S||I|)\)でpairwiseと同等ですが、 \(c\)あたりに第二項\(\ln \sum_{j \in I}\exp({\hat{y}(j|c)})\)を一度だけ計算すればいいような場合(4番目の右辺)は\(\mathcal{O}(|C||I|)\)でpointwiseと同等です。 この方法は高速かつ実装も単純、その上でtop-\(n\)性能が良いことが知られており、注目されています。自然言語処理の分野などで広く使われているsoftmaxモデルですが、推薦の分野に限っても割と前から存在はしており、 Collaborative competitive filtering [9]などがあります。 最近はこれを参考にしてVariational autoencoder (VAE)ベースの手法[10]をNetflix researchのチームが提案していました。 Criteo AILabがRecSys’20で発表した論文もこのアプローチをとっています[11]。 しかし、VAEベースの手法は全item数と同数の出力層を持つため、一般にitems数が大きいデータに対してそのままでは(おそらく)適用できません。とはいえ、このような問題はVAEだけの話ではなく、log-partitionと呼ばれる密度関数の正規化項の対数\(\ln \sum_{j \in I}\exp({\hat{y}(j|c)})\)は logの中で和を取っている関係で、intractableであり、これの計算にはitem全件のスコアリングが必要となります。 この点については後でもう少し詳しく説明します。

余談ですが、IRGAN [12]もsoftmaxモデルによる密度推定の手法だとみなすことが可能です。 IRGANは、discriminator \(\phi\)とgenerator \(\theta\)に対して、以下のような目的関数を最適化します。 \[ \begin{align} \min_{\theta}\max_{\phi} \sum_{c \in C}\left(\hat{\mathbb{E}}_{p(i|c)}\left[\log D_{\phi}(i)\right] + \mathbb{E}_{q_{\theta}(i|c)}\left[\log (1-D_{\phi}(i))\right]\right) \end{align} \] ここで、\(\hat{\mathbb{E}}_{p(i|c)}\)は\(p(i|c)\)上の経験期待値、サンプル平均です。
今、generatorの分布\(q_{\theta}(c,i)\)は、実はsoftmaxモデルです。 \[ \begin{align} q_{\theta}(i|c) = \frac{\exp(\hat{y}_{\theta}(i|c))}{\sum_{j \in I}\exp_{\theta}(\hat{y}(j|c))} \end{align} \] よって、GANの当初の目的を思い出すと、generatorをrankerにした場合、\(p(i|c)=q_{\theta}(i|c)\)なるsoftmaxモデルによる密度推定だとわかります。 IRGANは\(q_{\theta}(i|c)\)からのサンプリングを必要とするので、単純にはepoch毎などに全件のスコア予測とサンプリングが必要になります。 また、密度関数をdisrcriminatorのために保存するのは\(|C||I|\)の密行列を保持することになるため、工夫が必要になります。
みなせるとは言いましたが、item推薦においてgeneratorをrankerとして使うと性能が悪いので、IRGANのgeneratorはnegative samplerとしてpointwise, pairwiseの手法(discriminator)を補助する目的で主に使われています。

重み関数とサンプリング分布

4.1章冒頭の話題を元に、重み関数とサンプリング分布の関係を説明します。色々な論文を読む際にこのあたりを理解して読むと手法の対応がわかりやすいので、すごく重要な話題だと思います。

今、item数が巨大な場合、一度の学習ステップで全件を使うのは高コストです。 そこで、SGDを使うのが主流化してきています。 ニューラルネットなどのモデルを使う場合にもSGDを使うのが一般的かと思います。 ミニバッチに含まれるデータをどのようなサンプリング分布で選ぶかというのが極めて重要な要素になります。

例えば、以下のようなミニバッチサンプリングのアルゴリズムを考えます。

• 正例\((c,i) \in S\)を一つサンプリング
• 負item \(j \in I \setminus I_c\)を分布\(q(j|c)\)から\(m\)個サンプリング

上のアルゴリズムによって得られたミニバッチに対して、\(\tilde{\alpha}\)なる重み関数を用いてpointwise損失の和を最小化するとき、 SGDを通して最小化される期待損失は、 \[ \begin{align} \hat{\mathbb{E}}_{p(c|y=+1)} \left[\hat{\mathbb{E}}_{p(i|c,y=+1)} \left[\tilde{\alpha} (c,i) l(\hat{y}(i|c), 1)] + m \mathbb{E}_{q(j|c)}\left[\tilde{\alpha} (c,j) l(\hat{y}(j|c), 0)\right] \right]\right] \end{align} \] となります。 よって、期待値の定義から、以下のような重みで損失和を最小化していることになります。 \[ \begin{align} \alpha(c,j) = \begin{cases} \widehat{p(c,i|y=+1)}\tilde{\alpha}(c,j) & \text{if } (c,j) \in S\\ m \widehat{p(c|y=+1)}q(j|c)\tilde{\alpha}(c,j) & \text{if } (c,j) \notin S \end{cases} \end{align} \] ここで\(\widehat{p(c,i|y=+1)}\)と\(\widehat{p(c|y=+1)}\)は経験確率で、 \[ \begin{align} \widehat{p(c,i|y=+1)} = \begin{cases} \frac{1}{|S|} & \text{if } (c,i) \in S \\ 0 & \text{if } (c,i) \notin S \end{cases},\,\,\,\,\,\widehat{p(c|y=+1)} = \frac{|I_c|}{|S|} \end{align} \] になります。 よって、\(1/|S|\)は全サンプルに対して定数なので無視すると、結局以下のような重みを使っていることになります。 \[ \begin{align} \alpha(c,j) = \begin{cases} \tilde{\alpha}(c,j) & \text{if } (c,j) \in S\\ m |I_c|q(j|c)\tilde{\alpha}(c,j) & \text{if } (c,j) \notin S \end{cases} \end{align} \] このように、サンプリング分布の密度関数と重みは、同様の働きをすることがわかります。 このことから、サンプリング手法と重み付け手法は似た手法が多く存在します。 例えばitemの人気度 (\(p(i|y=+1)\)の経験確率)に応じた負例サンプリングやサンプル重みなどがそうです。

Context-independent Sampler/Weights

4.1.1章で触れられているとおり、負例サンプリングの分布\(q(j|c)\)が\(c\)と独立なとき、ミニバッチサンプリングを効率的に行うことができます。 たとえば、人気度に基づくサンプリングは\(c\)と独立な上に、静的な分布なので、\(p(i|y=+1)\)の経験分布からitemを複数個サンプリングするだけで複数のcontextを含むミニバッチ用の負例集合を構築することができます。ここで、ミニバッチ内に含まれるcontextの集合を\(C_B\)、(正例と負例を合わせた)itemの集合を\(I_B\)とすると、context独立なサンプリングにおいては、\(C_B\)の各要素に対して\(I_B\)のほぼ全てを使うことが可能なため、 内積モデルを採用していれば、ミニバッチ内の全てのcontext-item対 (\(C_B \times I_B\))に対するスコアリングを行列積一回で計算することが可能です。 一度の学習ステップで\(|C_B||I_B|\)の対に対する損失を効率的に計算できるため、使用メモリに対してミニバッチ内サンプルを多く確保でき、学習の収束が速くなる傾向があります。 特に、GPUなどを使う場合に、ミニバッチに含まれる全対を省メモリかつ高速に使うことで学習効率の大幅な向上が可能です。

contextに独立なサンプル重みも、ALSなどでよく使われる計算量削減のテクニックが使えます。 この話題は5.1.1章でGramian trickという名前で紹介されています。このブログではグラム行列を使わずに説明してみたいと思います。 さらに、元論文5.1.1章ではcontextとitem両方に独立な定数重みでの例を示していますが、contextかitemの片方に対して独立ならばGramian trickを使うことは可能なので、ここで示したいと思います。 以下のように、二乗誤差を元にしたpointwise損失の和を考えます。 \[ \begin{align} L(\theta, S) &= \sum_{c \in C} \sum_{i \in I} \alpha (c,i) (\hat{y}(i|c) – y(i|c))^2 \\ &= \sum_{c \in C} \sum_{i \in I_c} \alpha (c,i) (1 – \hat{y}(i|c))^2 + \sum_{c \in C} \sum_{i \in I \setminus I_c} \alpha (c,i) (0 – \hat{y}(i|c))^2 \\ &= -2 \sum_{c \in C} \sum_{i \in I_c} \alpha (c,i) \hat{y}(i|c)) + \sum_{c \in C} \sum_{i \in I} \alpha (c,i) \hat{y}(i|c)^2 + const.\\ \end{align} \] このことから、全context-item対を含む第二項が高い計算コストの原因になります。 今、重みがcontextに対して独立なとき (\(\alpha(c,i)=\alpha(i)\))、内積モデルを採用していれば、全対に対する損失は以下のように変形できます。 \[ \begin{align} \sum_{c \in C} \sum_{i \in I} \alpha (c,i) \hat{y}(i|c)^2 &= \sum_{c \in C} \sum_{i \in I} \alpha (i) \langle \phi(c), \psi(i) \rangle^2 \\ &= \sum_{c \in C} \sum_{i \in I} \alpha (i) \left(\sum_{k=1}^{d}\phi(c)_k\psi(i)_k\right)^2 \\ &= \sum_{c \in C} \sum_{i \in I} \alpha (i) \sum_{k=1}^{d}\sum_{l=1}^{d}\phi(c)_k\phi(c)_l\psi(i)_k \psi(i)_l \\ &= \sum_{k=1}^{d}\sum_{l=1}^{d} \left(\sum_{c \in C} \phi(c)_k\phi(c)_l\right) \left(\sum_{i \in I} \alpha (i) \psi(i)_k \psi(i)_l\right) \end{align} \] このように、\(\mathcal{O}(|C||I|d^2)\)の計算量を\(\mathcal{O}((|C|+|I|)d^2)\)に抑えることができます。

さらに、Neural Collaborative Filtering (NCF)の手法で使われている、Generalized Matorix Factorizationに対するGramian trickも存在します。上の説明を前提にするとかなり自明なのですが、元論文には書いていないので示しておこうと思います。 GMFでは、学習可能な\(d\)次元ベクトル\(\bf{h}\)を追加して、以下のようにスコアを作ります。 \[ \begin{align} \hat{y}(i|c) = \langle {\bf h}, \phi(c) \odot \psi(i) \rangle = \sum_{k=1}^d {\bf h}_k \phi(c)_k \psi(i)_k \end{align} \] ここで、\(\odot\)は要素積です。\(\bf{h}\)は次元毎の重み捉えることができます。\[ \begin{align} \sum_{c \in C} \sum_{i \in I} \alpha (c,i) \hat{y}(i|c)^2 &= \sum_{c \in C} \sum_{i \in I} \alpha (i) \langle {\bf h}, \phi(c) \odot \psi(i) \rangle^2 \\ &= \sum_{c \in C} \sum_{i \in I} \alpha (i) \left(\sum_{k=1}^d {\bf h}_k \phi(c)_k \psi(i)_k\right)^2 \\ &= \sum_{c \in C} \sum_{i \in I} \alpha (i) \sum_{k=1}^{d}\sum_{l=1}^{d}{\bf h}_k {\bf h}_l \phi(c)_k\phi(c)_l\psi(i)_k \psi(i)_l \\ &= \sum_{k=1}^{d}\sum_{l=1}^{d} \left( {\bf h}_k {\bf h}_l \right) \left(\sum_{c \in C} \phi(c)_k\phi(c)_l\right) \left(\sum_{i \in I} \alpha (i) \psi(i)_k \psi(i)_l\right) \end{align} \]　このような変形はSGDベースの手法においても巨大なミニバッチを使う場合の学習ステップあたりの計算量を減らすためにも使われています[13]。

Pairwise手法について

4.2章では色々なpairwise手法について分類しています。 今回のブログではここの説明は省略しようと思います。 特に補足することもないのですが、基本は重み関数とサンプリング分布の関係性を理解してから読むと良いかと思います。 元論文に書いていない話題としては、push付きのpairwise損失 (e.g., p-norm push)でしょうか。参考文献として、item推薦に適用した例[14]を示しておきます。

Softmax損失と合わせて使うアルゴリズム

Sampled softmax

個人的にこの話題には興味があるので、 4.3章は少し補足を多めにして説明したいと思います。 先程も触れましたが、対数正規化項 (i.e., log-partition)がsoftmax modelによる密度関数の最尤推定において問題になることを説明します。 密度推定として定式化した場合の、観測\((c,i) \in S\)に対する負の対数尤度は以下でした。 \[ \begin{align} -\ln p(i|c) = – \hat{y}(i|c) + \ln \sum_{j \in I}\exp({\hat{y}(j|c)}) \end{align} \] これの微分を考えると、 \[ \begin{align} \nabla_\theta -\ln p(i|c) &= – \nabla_\theta\hat{y}(i|c) + \frac{1}{Z(\hat{y}(\cdot|c))}\nabla_\theta\sum_{j \in I}\exp({\hat{y}(j|c)})\,\,\,\,(Z(\hat{y}(\cdot|c)) = \sum_{j \in I}\exp({\hat{y}(j|c)})) \\ &= – \nabla_\theta\hat{y}(i|c) + \sum_{j \in I}\frac{\exp({\hat{y}(j|c)})}{Z(\hat{y}(\cdot|c))}\nabla_\theta\hat{y}(j|c) \\ &= – \nabla_\theta\hat{y}(i|c) + \mathbb{E}_{p(j|c)}\left[\nabla_\theta\hat{y}(j|c)\right] \,\,\,\,(p(j|c)=\frac{\exp({\hat{y}(j|c)})}{Z(\hat{y}(\cdot|c))}) \end{align} \] となるので、不偏な勾配を得るための負例のサンプリング分布は\(p(j|c)\)だとわかります。 つまり、全itemに対するスコアリングを避けてサンプリングしようにも、正規化項の対数\(\ln \sum_{j \in I}\exp({\hat{y}(j|c)})\)の勾配を得るために正規化項の計算を含む分布\(p(j|c)\)を構築することになります。 それでは実行不可能なので、より簡単な分布\(q(j|c)\)を用いて重点サンプリングするのがsampled softmaxという手法です。 Sampled softmaxでは、log-partitionの計算時に逆確率で割り引きます。\(m\)個の負例\(j\)を使って、 \[ \begin{align} l(c,i) = – \hat{y}(i|c) + \ln \left(\exp({\hat{y}(i|c)}) + \sum_{l = 1}^{m} \frac{\exp({\hat{y}(j_l|c)})}{m q(j_l|c)}\right) \end{align} \] という損失を使います。\(1/ (m q(j|c))\)の部分が\(\tilde{\alpha}(c,j)\)に相当します。 item推薦では理想的な分布から負例をサンプリングするのは計算コストの都合上難しいことが多く、より単純なサンプリング分布を使った重点サンプリングが基本的な方法になります。 重み関数はその際における重要度重みとして捉えるとわかりやすいのではないかと個人的に思っています。

このときの推定された勾配は、 \[ \begin{align} &- \nabla_\theta\hat{y}(i|c) + \frac{\exp({\hat{y}(i|c)})}{\hat{Z}(\hat{y}(\cdot|c))}\nabla_\theta\hat{y}(i|c) + \sum_{j=1}^{m}\left[\left(\frac{\exp({\hat{y}(j|c)})}{\hat{Z}(\hat{y}(\cdot|c))m q(j|c)}\right)\nabla_\theta\hat{y}(j|c)\right] \end{align} \] ここで、 \[ \hat{Z}(\hat{y}(\cdot|c)) = \exp(\hat{y}(i|c)) + \sum_{l = 1}^{m}\frac{\exp({\hat{y}(j_l|c)})}{m q(j|c)} \] です。このことから、\(m\)個の負例による正規化項の推定値\(\hat{Z}(\hat{y}(\cdot|c))\)が重要度重みの分母に現れてしまうので、一ステップでサンプルを全て使わない限り、学習ステップを増やしても不偏な勾配は得られません。 よって、この定式化である限りバイアスの問題を完全に解決するのは難しいのですが、軽減する方法を次に説明します。
また、元論文でも深入りしていませんし、ブログではあまり説明しませんが、log-partition (or partition)が厄介であるという話題は、energy-based modelや指数分布族の推定の文脈でも色々な方法が最近も提案されています。KL divergenceの双対形式(Donsker-Varadhan variational formula)を使ってlog-partitionを変形できることを利用したりします。

Kernel based Sampling

\(p(j|c)\)からサンプリングするのは困難ですが、\(q(j|c)\)を\(p(j|c)\)に近づけることでsampled softmaxの損失のバイアスを小さくすることができます。 itemの人気度などを使うことでそれなりに品質の良いサンプリングを行うことはできますが、item数が大きくなると効率が悪いかもしれません。 これを軽減する方法は、密度関数の正規化項\(Z\)を正確かつ低い計算コストで計算できるような分布で\(p(j|c)\)を近似することです。
今、内積モデルにおいて、\(p(i|c) \propto \exp(\langle \phi(c), \psi(i) \rangle)\)でした。 これを\(d\)次元の特徴空間からより大きい\(D\)次元の空間への特徴マッピング\(\pi \colon \mathbb{R}^d \mapsto \mathbb{R}^D\)と線形カーネルだと考えてみます。 \[ \begin{align} p(j|c,\theta) \propto \exp(\langle \phi(c), \psi(j) \rangle) \approx \langle \pi(\phi(c)), \pi(\psi(j)) \rangle \propto q(j|c,\theta) \end{align} \] のようになる都合の良い写像\(\pi\)があるならば、 \[ \begin{align} \sum_{j \in I}\exp(\hat{y}(j|c)) \approx \sum_{j \in I} \langle \pi(\phi(c)), \pi(\psi(j)) \rangle = \left\langle \pi(\phi(c)), \sum_{j \in I}\pi\left(\psi(j) \right) \right\rangle = \left\langle \pi(\phi(c)), {\bf z} \right\rangle \end{align} \] のようになることから、カーネルにpropotionalな分布を考えると、 \[ \begin{align} q(j|c,\theta) =\frac{\langle \phi(c), \psi(j) \rangle}{\sum_{j’ \in I} \langle \phi(c), \psi(j’) \rangle} =\frac{\langle \phi(c), \psi(j) \rangle}{\langle \phi(c), {\bf z} \rangle} \end{align} \] のように密度関数の正規化項をベクトル\({\bf z}=\sum_{j \in I}\pi\left(\psi(j) \right)\)と射影されたcontextベクトルとの内積で計算できることになります。この方法を使うと分割統治法を使って、\(\mathcal{O}(D\log_2 |I|)\)で一つの負例をサンプルできるようです。

\(\pi\)も色々と提案されているので、このブログではQuadratic kernelを用いた例を説明してみたいと思います。Rendle氏らがsampled softmax向けに使うことを提案していました[15]。 Quadratic kernelは \[ \begin{align} K(\phi(c), \psi(i)) &= \alpha \langle \phi(c), \psi(i) \rangle^2 + 1 \\ &= \alpha \left(\sum_{k=1}^d\sum_{l=1}^d \phi(c)_k\phi(c)_l \psi(i)_k\psi(i)_l \right) + 1 \\ &= \alpha \langle \text{vec}(\phi(c) \otimes \phi(c)), \text{vec}(\psi(i) \otimes \psi(i)) \rangle + 1 \end{align} \] と書けます。\(\alpha \in \mathbb{R}^{+}\)はハイパーパラメータです。 \(\otimes\)は外積、\(\text{vec}(\cdot)\)は外積によって生じる\(d \times d\)行列を\(d^2\)次元のベクトルに直す演算です。結果、\(D=d^2\)となります。
展開した形を踏まえて、これに対応する特徴変換として以下があります。 \[ \begin{align} \pi({\bf v}) = \left[ \sqrt{\alpha}\text{vec}({\bf v} \otimes {\bf v}), 1 \right] \end{align} \]
Qudratic kernelを使う気持ちとしては、内積が大きいところではもちろん、大抵のcontext-item対の内積が0まわりに集まっている時に、\(\exp(0)=1=0^2+1\)より、良い近似を与えるからという感じです。
一方、負のスコア\(\langle \phi(c), \psi(i) \rangle < 0\)のケースでQuadratic kernelは近似に失敗します。 そこでabsolute softmaxをモデルに採用することが提案されています。 \[ \begin{align} p(i|c) = \frac{\exp(|\hat{y}(i|c)|)}{\sum_{j \in I} \exp(|\hat{y}(i|c)|)} \end{align} \] このようなモデルで本当に真の分布を十分表現できるのか気になるかと思いますが、 softmaxモデルを含む指数分布族はシフト不変性をもち、ある定数\(\beta \in \mathbb{R}\)でもって、 \(p(i|c) \propto \exp(\hat{y}(i|c))\exp(\beta) = \exp(\hat{y}(i|c)+\beta)\)とも書けます。 つまり、全logitが非負になるよう十分大きい\(\beta\)を考えることで、全てのsoftmaxモデルの解は対応するabsolute softmaxの解を持つことがわかります。

Pairwise損失との関係

Sampled softmax損失とpairwise損失との関係が4.3.3章に短く書かれています。 \(m=1\)としたときのsampled softmax損失は \[ \begin{align} – \hat{y}(i|c) + \ln \left(\exp(\hat{y}(i|c)) + \frac{\exp({\hat{y}(j|c)})}{q(j|c)}\right) &= \ln \frac{\exp({\hat{y}(i|c)}) + \frac{\exp({\hat{y}(j|c)})}{q(j|c)}}{\exp(\hat{y}(i|c))} \\ &= \ln \left(1 + \frac{\exp({\hat{y}(j|c) – \hat{y}(i|c)})}{q(j|c)}\right) \\ &= l(\hat{y}(i|c) – \hat{y}(j|c), 1) \end{align} \] なので、実はpairwise logistic損失と一致します。 ここまでだけだと元論文にも書いてあるので、理想的なsoftmax損失との関係性を見てみます。 \[ \begin{align} &\sum_{(c,i) \in S}\left[- \hat{y}(i|c) + \ln \left(\exp(\hat{y}(i|c)) + \sum_{j \in I, j \neq i}\exp(\hat{y}(j|c))\right)\right] \\ &\propto \hat{\mathbb{E}}_{p(c,i)}\left[ \ln \left(1 + \sum_{j \in I, j \neq i}\exp(\hat{y}(j|c)-\hat{y}(i|c))\right)\right] \\ &= \hat{\mathbb{E}}_{p(c,i)}\left[ \ln \left(1 + (1/p_0) \mathbb{E}_{p_0(j)}\left[\exp({\hat{y}(j|c)-\hat{y}(i|c)})\right]\right)\right] \\ &\leq \hat{\mathbb{E}}_{p(c,i)}\left[ \mathbb{E}_{p_0(j)}\left[\ln \left(1 + \exp({\hat{y}(j|c)-\hat{y}(i|c)} – \ln p_0)\right)\right]\right] \end{align} \] ここで、\(p_0(j)\)は\(I \setminus \{i\}\)上の一様分布で、その密度は\(p_0(j) = p_0 = 1/(|I|-1)\)です。 最後の不等式は\(\ln(1+\cdot)\)の凸性からJensenの不等式を適用しています。
不等式の右辺は、まさにlogistic pairwise損失+一様負例サンプリングをするBPRの期待損失においてlogitを\(- \ln p_0\)だけシフトしたものです。 softmax損失和はマージン付きBPRの期待損失の下界になっていることがわかりました。 ここから思うことは、(sampled) softmaxには元論文の4.2.1章で触れられているBPRと同じ問題があるかもしれないということです。 上で少し触れたように、pairwise logistic損失は順序を間違える回数を意味する0-1損失のsurrogateでした。 \[ \begin{align} \hat{\mathbb{E}}_{p(c,i)}\left[ \mathbb{E}_{p_0(j)}\left[\left\{\hat{y}(j|c) > \hat{y}(i|c)\right\}\right]\right] \approx \hat{\mathbb{E}}_{p(c,i)}\left[ \mathbb{E}_{p_0(j)}\left[\ln \left(1 + \exp({\hat{y}(j|c)-\hat{y}(i|c)})\right]\right)\right] \end{align} \] この左辺は1-AUCという指標と対応していることから、BPRはAUCを最適化していることが知られています。 しかし、AUCは上位に重みを持たない(top-weighted, top-heavyでない)指標であり、top-\(n\)件を評価するitem recommendationにおける典型的な指標であるRecall\(@n\)やnDCG\(@n\)とは興味のあるところが異なります。この性質によって、(特にitem数が多い時)良い性能が出ないことがあります。pairwise損失におけるこの問題に対する解決策は説明を省略した4.2章にまとめられています。 同じ問題がsoftmax損失にあるかはもう少し調べないとわかりませんが、上のような関係性が少なくともあるようです。

少しだけごまかしましたが、普通のBPRの実装だと未観測itemから負例をサンプリングするため、\(p_0(j)\)はあるcontext依存に\(I \setminus I_c\)上で定義します。そこにもsampled softmaxとBPRのズレはあるかもしれません。

まとめ

このブログでは、最近arxivに投稿された素晴らしいレビュー論文“Item Recommendation from Implicit Feedback”[1]の(個人的な興味によって選んだ)要所を、具体例を補足する形で説明しました。 Item recommendation from implicit feedbackという、大規模かつone-classなデータを前提としたランキング学習のタスクにおいて、損失とサンプリング分布/重み関数を主役にして、高いランキング性能をスケールするアルゴリズムで達成するための様々な工夫があることが伝わったなら幸いです。 元論文にはまだまだ色々な話題が書いてありますので、興味がある方はぜひ元論文の方をチェックしてみてください。