万字长文,详解推荐系统领域经典模型FM因子分解机

在上一篇文章当中我们剖析了Facebook的著名论文GBDT+LR,虽然这篇paper在业内广受好评,但是毕竟GBDT已经是有些老旧的模型了。今天我们要介绍一个业内使用得更多的模型,它诞生于2010年,原作者是Steffen Rendle。虽然诞生得更早,但是它的活力更强,并且衍生出了多种版本。我们今天剖析的就是这篇2010年最经典的原版论文。

说到推荐、广告的算法模型,几乎很难绕开FM,它是一个非常强的模型。理论简单、推导严谨、实现容易,并且效果不俗。即使是目前仍然在各大厂商当中发挥用场,在一些小厂当中甚至是主力模型。我们初看之前也许还有疑惑,但是相信当大家看完之后,必然会有全新的认识。

同样,这是一篇长文,希望大家耐心读到最后。

简介

FM的全称是Factorization Machines,翻译过来就是因式分解机,这么翻译很拗口,不得其神也不得其形,所以我们一般都不翻译简称为FM模型。

论文开篇就讲了FM模型的种种好处,比如和SVM对比,它能够在稀疏的特征当中仍然拥有很好的表现,并且FM的效率非常高,可以在线性时间内获得结果。并且不像非线性的SVM(核函数),FM并不需要将特征转换成对偶形式,并且模型的参数可以直接训练,也不用借助支持向量或者是其他的方法。

除了SVM之外,FM模型和其他的因式分解模型比如SVD、PITF、FPMC比较起来都有非常明显的优势。这些模型往往只针对一个特定的场景或者是特定的数据集,并且它们的训练以及优化方案都是根据场景定制化的。FM模型不需要定制化就可以实现同样好的效果,可以非常简易地应用在各个预测问题当中。

这段摘要其实没什么太多的内涵,主要就是吹嘘了一通FM。不过作者所言不虚,和他列举的这些模型比较起来,FM的确更加出色。

总结一下,FM模型的优点有如下几点:

  • FM模型的参数支持非常稀疏的特征,而SVM等模型不行
  • FM的时间复杂度为,并且可以直接优化原问题的参数,而不需要依靠支持向量或者是转化成对偶问题解决
  • FM是通用的模型,可以适用于任何实数特征的场景,其他的模型不行

paper在一开始的时候就表明了FM模型的优点,给足了好处,之后才是对FM模型的深入分析,让大家明白这些优点是如何得到的,以及它的工作原理。

稀疏场景

在有监督学习的场景当中,最常见的任务就是给定一个向量x,要求预测一个目标T。如果T是一个实数,那么这就是回归模型,如果T是一个类别的常量,就是一个分类模型。

这些都是机器学习的基础知识,相信大家都了解,然而对于线上排序的功能来说,我们重要的是给商品排序,而不是分类。常规来说我们可以将impression和click看成是两个类别进行分类预测,也可以直接用回归模型预测商品的点击率,根据点击率进行排序。这两种其实本质上是一样的,都是针对商品本身进行学习。还有一种方法是更加侧重商品之间的排序,我们的模型学习的是商品之间的优劣关系。

举个例子,比如我们用向量表示i商品的特征向量,表示j商品的特征向量,那么我们可以将这两个向量合并一起作为输入,进行分类预测。分类的结果表示的是i商品在j商品之前还是反之。这样我们可以通过多次预测,直接得到商品之间的排序关系,而不是根据一个分数进行排序。这样可以在只有正样本的情况下进行训练。这种方法直接训练的模型商品的优劣,业内叫做learning to rank。

然而不论是哪一种做法,都有一个问题绕不开就是特征的稀疏。举个很简单的例子,比如我们对商品的类目进行one-hot处理,在电商场景下商品的类目动辄上万个,那么one-hot之后的结果就是一个长度上万的数组,并且这个数组当中只有一位为1,其他均为0。当然这只是一个例子,除此之外还有许多许多的特征有可能是非常稀疏的。

我们用表示x向量当中非零的元素的个数,用表示所有x样本平均的非零元素的个数,在大多数真实数据场景下,我们都可以得到,这里的n是特征的维数。

真实场景例子

paper当中举了一个真实场景的例子,我们的问题是预测用户对于一部电影的评分。我们先来观察一下下图,下图是样本矩阵。

很明显上图的左边一大块是特征,右边的Target y表示的预测结果,也就是用户可能对电影做出的评价。这里一共有[1, 2, 3, 4, 5]这5种可能,也就是说这是一个多分类的问题。

接着我们再来看特征,特征一共也有5个部分,其中蓝色的部分表示的用户的one-hot。那么这个数组的长度应该是用户的数量,只有代表当前用户的那一维为1,其他均为0。

红色部分表示电影,同样是电影的one-hot,和用户的one-hot是一样的逻辑。代表当前电影的那一维度为1,其他为0。

之后是黄色的部分,表示的之前用户对于电影的评分行为,维度同样是电影的数量。凡是用户评分过的电影分数大于0,没有评分的等于0。得分的计算逻辑是1除以用户评论过的电影数量。比如第一行当中,第一个用户评价过前三部电影,所以前三部电影每一部分到了的分数。

绿色的部分只有1维,表示的是用户评论电影的时间。计算方法是将记录当中最早的日期作为基数(这里是2009年1月),之后每过一个月,增加1。比如2009年5就可以折算成5。

最后棕色的部分表示的是用户最近一次评论的电影,这同样是一个one-hot的数组,它的维度和电影的数量是一致的。

我们假设用户的数量是U,电影的数量是M,那么最后得到的整个特征的维度数应该是U+3M+1。即使是小众一些的电影评分网站,用户数也至少是以上百万起的,再加上电影的数量,这显然是一个非常庞大的数字。而在这么庞大的维度当中只有少数的一些维度是有值的,其余均为0。

对于这样稀疏的特征矩阵而言,一般的模型是很难保证效果的。为什么FM模型可以置身其外呢?下面,我们就来看看FM模型的原理。

FM模型原理

在我们介绍FM模型的方程之前,我们先来回顾一下线性回归的表达式,其实非常简单,只有一行:

也就是说,W是这样一个n+1维的向量,X是一个n x m的特征矩阵。这里的n是特征的维数,m是样本的数量。所以我们也经常把它写成,不管怎么写,形式都是一样的。

这个式子称为线性回归的原因也很简单,因为它就是一个简单的线性方程,只有一次项。但是一次项有时候效果不好,尤其是在特别稀疏的场景当中,刻画能力不够。我们做特征的时候经常会把两项特征组合起来做成新的组合特征,由于我们这样操作引入了新的特征,找到了新的特征组合,所以能够挖掘出之前无法得到的信息,因此模型也会有更好的效果。

当我们把特征进行二项组合之后,会得到这样的式子:

这里的分别表示两个不同的特征取值,对于n维的特征来说,这样的组合应该一共有种,这个值是,也就意味着我们需要同样数量的权重参数。但是有一个小问题,我们前面已经说过了,由于特征可能非常稀疏,导致n非常大,比如上百万,这里两两特征组合的特征量级大约是n的平方,那么因此带来的参数数量就是一个天文数字。想想看对于一个上百亿甚至更多参数空间的模型来说,我们需要多少训练样本才可以保证完全收敛?这是不可想象的。

既然如此,那么针对性的优化就是必不可少的。FM为人称道的也正是这一点,既引入了特征交叉,又解决了复杂度以及模型参数的问题,真的可以说是非常棒了。

FM解决这个问题的方法非常简单,它不再是简单地为交叉之后的特征对设置参数,而是设置了一种计算特征参数的方法。FM模型引入了新的矩阵V,矩阵V是一个n x k的二维矩阵。这里的k是我们设置的参数,一般不会很大,比如16、32之类。对于特征每一个维度i,我们都可以找到一个,它表示一个长度为k的向量。

于是我们可以用来计算得出上式当中的

也就是说我们用向量的内积来计算得到了就交叉特征的系数,相比于原先量级的参数而言,我们将参数的量级降低到了n x k。由于k是一个常数值,所以可以看成我们的参数数量是。通过这种方式我们把参数的数量降低了一个数量级。

有了V矩阵之后,刚才的公式可以改写成:

我们很容易可以知道,当k这个参数足够大的时候,我们可以保证成立。所以我们引入的参数矩阵V可以看成是对W矩阵做了一个因子分解,这也是FM得名的由来。当然在实际的应用场景当中,我们并不需要设置非常大的K,因为特征矩阵往往非常稀疏,我们可能没有足够多的样本来训练这么大量的参数,并且限制K也可以一定程度上提升FM模型的泛化能力

此外这样做还有一个好处就是有利于模型训练,因为对于有些稀疏的特征组合来说,我们所有的样本当中可能都是空的。比如在刚才的例子当中用户A和电影B的组合,可能用户A在电影B上就没有过任何行为,那么这个数据就是空的,我们也不可能训练出任何参数来。但是引入了V之后,虽然这两项缺失,但是我们仍然可以得到一个参数。因为我们针对用户A和电影B训练出了一个向量参数,我们用这两个向量参数点乘,就得到了这个交叉特征的系数。

这种做法有两种理解方式,一种就是刚才说的,我们对于一些稀疏的组合也可以很好地计算出系数。另外一种理解方式是这其实也是一种embedding的方式,将某一个id映射成向量。所以业内也有使用FM来做embedding的。

复杂度优化

接下来我们来看另外一个很精髓的内容,就是关于FM模型的复杂度优化。我们观察一下刚才上面的式子,不难发现,目前对于预测一条样本的计算复杂度为

n我们前文说过了是一个比较大的数字,那么显然级的计算也是我们不能接受的。所以对它进行优化也是必须的,并且这里的优化非常简单,可以直接通过数学公式的变形推导得到。

这个是它的原式,对于这个式子来说,前面两项的复杂度是,我们可以先忽略,重点来看最后一项。我们要做的就是通过数学公式的变形来对这一项进行化简:

简单来解释一下这个推导过程,第一行我想大家应该都能看懂,第二行也很好理解,其实就是把向量内积展开。第二行到第三行的转化也不难理解,这里有三个,我们提取出的是最里面的,因为是有限项求和,我们调换顺序并不会影响结果。提取出了公因式之后,得到的结果是两个平方项。我们观察一下可以发现,这两个平方项的计算复杂度都是,再加上外面一层的复杂度,整体的复杂度是

这样我们就完成了FM模型预测的优化。

模型训练

我们再来看下模型训练的过程,我们写出变形之后的原式:

针对FM模型我们一样可以使用梯度下降算法来进行优化,既然要使用梯度下降,那么我们就需要写出模型当中所有参数的偏导。

我们可以把参数分成三个部分,第一个部分自然是。第二个部分是,对于其中每一个来说它的系数是。由于样本是固定的,我们要把样本的特征值看成是常数。所以。最后一个部分就是看起来复杂很多,但其实偏导也不难求,我们首先明确这里要优化的参数是,所以我们可以忽略最外层的一层,直接对括号内的进行求导,得出的结果是

我们把这三种综合一下,就得到了:

其中和i是独立的,所以它是可以提前算好的,这样一来对于所有参数项,我们都可以在的时间内计算出它们的梯度。由于我们所有参数的量级是,意味着我们可以在时间内对所有的参数完成梯度下降的更新。

拓展到d维

我们刚才介绍的FM模型专注的是二维特征交叉的情况,如果我们愿意也可以将它拓展到d维特征交叉的情况。这样的话我们的特征可以拟合任意维特征交叉的相互关系。

我们仿照刚才的公式,可以写出FM模型推广到d维的方程:

前面两项都很好理解,我们着重来看第三项。第三项当中包含了从2维到d维交叉特征的情况,我们以d=3为例,那么这一项当中应该包含二维的交叉项以及三维的交叉项,应该是这样的:

这个式子整体上和之前的形式是一样的,我们不难分析出它的复杂度是。当d=2的时候,我们通过一系列变形将它的复杂度优化到了,而当d>2的时候,没有很好的优化方法,而且三重特征的交叉往往没有意义,并且会过于稀疏,所以我们一般情况下只会使用d = 2的情况。

代码实现

上面这段大家了解一下知道有这么回事就可以了,实际上意义不大。最后的时候paper还比较了FM和其他一些经典的模型的效果,比如SVD、SVM等,也没太多价值,因为现在在推荐领域已经几乎没有人直接使用这些模型了。最后我贴一下我用pytorch和TensorFlow两个框架分别实现的FM模型,给大家做一个参考。

Pytorch

import torch
from torch import nn

ndim = len(feature_names)
k = 4

class FM(nn.Module):
    def __init__(self, dim, k):
        super(FM, self).__init__()
        self.dim = dim
        self.k = k
        self.w = nn.Linear(self.dim, 1, bias=True)
        # 初始化V矩阵
        self.v = nn.Parameter(torch.rand(self.dim, self.k) / 100)
        
    def forward(self, x):
        linear = self.w(x)
        # 二次项
        quadradic = 0.5 * torch.sum(torch.pow(torch.mm(x, self.v), 2) - torch.mm(torch.pow(x, 2), torch.pow(self.v, 2)))
        # 套一层sigmoid转成分类模型,也可以不加,就是回归模型
        return torch.sigmoid(linear + quadradic)
    
    
fm = FM(ndim, k)
loss_fn = nn.BCELoss()
optimizer = torch.optim.SGD(fm.parameters(), lr=0.005, weight_decay=0.001)
iteration = 0
epochs = 10

for epoch in range(epochs):
    fm.train()
    for X, y in data_iter:
        output = fm(X)
        l = loss_fn(output.squeeze(dim=1), y)
        optimizer.zero_grad()
        l.backward()
        optimizer.step()
        iteration += 1        
        if iteration % 200 == 199:
            with torch.no_grad():
                fm.eval()
                output = fm(X_eva_tensor)
                l = loss_fn(output.squeeze(dim=1), y_eva_tensor)
                acc = ((torch.round(output).long() == y_eva_tensor.view(-11).long()).sum().float().item()) / 1024
                print('Epoch: {}, iteration: {}, loss: {}, acc: {}'.format(epoch, iteration, l.item(), acc))
            fm.train()

TensorFlow

import tensorflow as tf
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import accuracy_score, roc_auc_score
import numpy as np


class FactorizationMachine:

    def __init__(self, n_dim=1, k=4, learning_rate=0.05, epochs=8):
        self._learning_rate = learning_rate
        self._n_dim = n_dim
        self._k = k
        self._epochs = epochs
        self.sess = tf.Session()
        self.x_input = tf.placeholder(shape=[None, self._n_dim], dtype=tf.float32)
        self.y_input = tf.placeholder(shape=[None1], dtype=tf.float32)
        # 初始化W和V
        self.w = tf.Variable(tf.truncated_normal(shape=[self._n_dim, 1], dtype=tf.float32))
        self.V = tf.Variable(tf.truncated_normal(shape=[self._n_dim, self._k], dtype=tf.float32))
        self.b = tf.Variable(tf.truncated_normal(shape=[11]))
        self.linear = tf.add(self.b, tf.matmul(self.x_input, self.w))
        self.quadratic = 1/2 * tf.reduce_sum(tf.square(tf.matmul(self.x_input, self.V)) - tf.matmul(tf.square(self.x_input), tf.square(self.V)), axis=1, keepdims=True)
        self.y_out = self.linear + self.quadratic
        self.y_pred = tf.round(tf.sigmoid(self.y_out))
        self.loss = tf.reduce_mean(tf.nn.sigmoid_cross_entropy_with_logits(logits=self.y_out, labels=self.y_input))
        self.train_op = tf.train.GradientDescentOptimizer(self._learning_rate).minimize(self.loss)
        self.accuracy = tf.reduce_mean(tf.cast(tf.equal(self.y_pred, self.y_input), tf.float32))
        init = tf.global_variables_initializer()
        self.sess.run(init)

    def train(self, X, Y, iterations=1000, batch_size=16, validation_size=0.1):
        x_train, x_test, y_train, y_test = train_test_split(X, Y, test_size=validation_size)
        for epoch in range(self._epochs):
            for i in range(iterations):
                rand_idx = np.random.choice(x_train.shape[0], size=batch_size)
                rand_x = x_train[rand_idx]
                rand_y = y_train[rand_idx]
                self.sess.run(self.train_op, feed_dict={self.x_input: rand_x, self.y_input: rand_y})
                if i % 100 == 99:
                    loss = self.sess.run(self.loss, feed_dict={self.x_input: x_test, self.y_input: y_test})
                    acc = self.sess.run(self.accuracy, feed_dict={self.x_input: x_test, self.y_input: y_test})
                    print('epoch = {}, iteration ={}, loss = {}, accuracy ={}'.format(epoch, i, loss, acc))

    def predict(self, x):
        return self.sess.run(self.y_pred, feed_dict={self.x_input: x})

注:TensorFlow代码使用的1.0版本

FM的paper到这里我们就剖析完了,也给出了代码实现。但是这并没有结束,后续关于FM又有了好几个变种的更新版本。像是AFM、FFM、PFM等等。关于这些paper,将会在后续给大家一一介绍。

今天的文章就到这里,衷心祝愿大家每天都有所收获。如果还喜欢今天的内容的话,请来一个三连支持吧~(点赞、关注、转发

原文链接,求个关注