FM（Factorization Machine）因式分解机 与 TensorFlow实现 详解

1，线性回归（Linear Regression）

线性回归，即使用多维空间中的一条直线拟合样本数据，如果样本特征为：

\[x = ({x_1},{x_2},...,{x_n})\]

模型假设函数如下：

\[\hat y = h(w,b) = {w^T}x + b,w = ({w_1},{w_2},...,{w_n})\]

以均方误差为模型损失，模型输入样本为(x⁽¹⁾,y⁽¹⁾),(x⁽²⁾,y⁽²⁾),...,(x^(m),y^(m))，损失函数如下：

\[l(w,b) = \sum\limits_{j = 1}^m {{{({{\hat y}^{(j)}} - {y^{(j)}})}^2}}  = \sum\limits_{j = 1}^m {(\sum\limits_{i = 1}^n {{w_i}x_i^{(j)} + b}  - {y^{(j)}})} \]

2，逻辑回归（Logistic Regression）

线性回归用于预测标记为连续的，尤其是线性或准线性的样本数据，但是有时样本的标记为离散的，最典型的情况莫过于标记为0或者1，此时的模型被称为分类模型。

为了扩展线性回归到分类任务，需要一个函数将(-∞,+∞)映射到(0,1)（或者(-1,1)等任意两个离散值），函数最好连续可导，并且在自变量趋向于-∞时无限趋近于因变量下限，在自变量趋向于+∞时无限趋近于因变量上限。符合条件的函数有不少，比如tanh函数，logistic函数。

如果映射函数采用logistic函数，设模型假设函数为样本预测分类概率，模型假设函数为：

\[\begin{array}{l}
P(y = 1) = h(w,b) = \frac{1}{{1 + {e^{ - {h_{linear}}(w,b)}}}} = \frac{1}{{1 + {e^{ - ({w^T}x + b)}}}}\\
P(y = 0) = 1 - P(y = 1)
\end{array}\]

另外，可以从另一个角度考虑从回归模型扩展为分类模型，令：

\[\ln (\frac{{P(y = 1)}}{{P(y = 0)}}) = \ln (\frac{{h(w,b)}}{{1 - h(w,b)}}) = {h_{linear}}(w,b) = {w^T}x + b\]

同样可以得到映射函数为logistic函数的假设函数：

\[h(w,b) = \frac{1}{{1 + {e^{ - ({w^T}x + b)}}}} = \frac{1}{{1 + {e^{ - (\sum\limits_{i = 1}^n {{w_i}{x_i}}  + b)}}}}\]

使用（负）对数似然函数作为损失函数：

\[\begin{array}{l}
l(w,b) = - \log (\prod\limits_{j = 1}^m {P{{({y^{(j)}} = 1)}^{{y^{(j)}}}}P{{({y^{(j)}} = 0)}^{(1 - {y^{(j)}})}}} )\\
= - \sum\limits_{j = 1}^m {({y^{(j)}}\log P({y^{(j)}} = 1) + (1 - {y^{(j)}})\log P({y^{(j)}} = 0))} \\
= - \sum\limits_{j = 1}^m {({y^{(j)}}\log (g(b + \sum\limits_{i = 1}^n {{w_i}x_i^{(j)}} )) + (1 - {y^{(j)}})\log (1 - g(b + \sum\limits_{i = 1}^n {{w_i}x_i^{(j)}} ))} ,g(z) = \frac{1}{{1 + {e^{ - z}}}}
\end{array}\]

很多场景下，样本包含一些离散的label特征，这些特征很难连续量化，最简单的处理方式就是one-hot，即将离散的每个数映射到多维空间中的一维。离散特征包含多少个可能值，转换后的特征就包含多少维。这样处理后样本特征的维数会变得非常巨大，另外大部分维度的值都是0。

一方面，逻辑回归核心是一个线性模型，因此计算规模随着样本特征数的增长而线性增长，相较其他机器学习模型来说计算量（随特征数增长）的增长率较小；另一方面，逻辑回归假设函数与损失函数中，各个特征对于结果是独立起作用的，因此在样本数足够的前提下，也不会受到大量值为0的特征的干扰。因此特别适合这类场景下的分类问题。

3，因式分解机（Factorization Machine）

逻辑回归最大的优势是简单，最大的劣势也是太简单，没有考虑特征之间的相互关系，需要手动对特征进行预处理。因此就有了包含特征两两交叉项的假设函数：

\[h(w,b) = g(b + \sum\limits_{i = 1}^n {{w_i}{x_i}}  + \sum\limits_{i = 1}^n {\sum\limits_{j = i + 1}^n {{\omega _{ij}}{x_i}} {x_j}} ),g(z) = \frac{1}{{1 + {e^{ - z}}}}\]

此时有一个严重的问题，由于特征维数有可能是非常巨大的，很有可能训练样本中有一些特征组合x_ix_j（x_i、x_j都不为0）是完全不存在的，这样ω_ij就无法得到训练。实际使用模型时，如果出现特征组合x_ix_j，模型就无法正常工作。

为了减小训练计算量，也为了规避上面说的这种情况，我们引入辅助矩阵v，n为特征总维数，k为超参数

 \[v \in {R^{n*k}}\]

然后使用vv^T代替参数矩阵ω，可得

\[{\omega _{ij}} = \sum\limits_{r = 1}^k {{v_{ir}}v_{rj}^T}  = \sum\limits_{r = 1}^k {{v_{ir}}{v_{jr}}} \]

可得假设函数线性部分最后一项可化简为：

\[\begin{array}{l}
\sum\limits_{i = 1}^n {\sum\limits_{j = i + 1}^n {{\omega _{ij}}{x_i}} {x_j}} \\
= \sum\limits_{i = 1}^n {\sum\limits_{j = i + 1}^n {\sum\limits_{r = 1}^k {{v_{ir}}{v_{jr}}{x_i}{x_j}} } } \\
= \sum\limits_{r = 1}^k {(\sum\limits_{i = 1}^n {\sum\limits_{j = i + 1}^n {{v_{ir}}{v_{jr}}{x_i}{x_j})} } } \\
= \frac{1}{2}\sum\limits_{r = 1}^k {(\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{v_{ir}}{v_{jr}}{x_i}{x_j}} } - \sum\limits_{i = 1}^n {{v_{ir}}{v_{ir}}{x_i}{x_i}} )} \\
= \frac{1}{2}\sum\limits_{r = 1}^k {(\sum\limits_{i = 1}^n {{v_{ir}}{x_i}} \sum\limits_{j = 1}^n {{v_{jr}}{x_j}} - \sum\limits_{i = 1}^n {{v_{ir}}{v_{ir}}{x_i}{x_i}} )} \\
= \frac{1}{2}\sum\limits_{r = 1}^k {({{(\sum\limits_{i = 1}^n {{v_{ir}}{x_i}} )}^2} - \sum\limits_{i = 1}^n {{{({v_{ir}}{x_i})}^2}} )}
\end{array}\]

同逻辑回归一样，令模型假设函数为样本预测分类概率，还是使用（负）对数似然函数为损失函数：

\[\begin{array}{l}
l(b,w,v) = - \log (\prod\limits_{j = 1}^m {P{{({y^{(j)}} = 1)}^{{y^{(j)}}}}P{{({y^{(j)}} = 0)}^{(1 - {y^{(j)}})}}} ) = - \sum\limits_{j = 1}^m {({y^{(j)}}\log (h(w,b){|_{x = {x^{(j)}}}}) + (1 - {y^{(j)}})\log (1 - h(w,b){|_{x = {x^{(j)}}}}))} \\
h(w,b){|_{x = {x^{(j)}}}} = g(b + \sum\limits_{i = 1}^n {{w_i}x_i^{(j)}} + \frac{1}{2}\sum\limits_{r = 1}^k {({{(\sum\limits_{i = 1}^n {{v_{ir}}x_i^{(j)}} )}^2} - \sum\limits_{i = 1}^n {{{({v_{ir}}x_i^{(j)})}^2}} )} )\\
g(z) = \frac{1}{{1 + {e^{ - z}}}}
\end{array}\]

4，关于FM假设函数的讨论

仔细看FM假设函数线性部分最后一项的化简过程，会发现一个问题，这个“化简”过程仅仅是化简了式子的形式，其实并没有减少计算量，反倒是增加了计算量，如果假设函数不是：

\[h(w,b) = g(b + \sum\limits_{i = 1}^n {{w_i}{x_i}}  + \sum\limits_{i = 1}^n {\sum\limits_{j = i + 1}^n {{\omega _{ij}}{x_i}} {x_j}} )\]

而是：

\[h(w,b) = g(b + \sum\limits_{i = 1}^n {{w_i}{x_i}}  + \sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{\omega _{ij}}{x_i}} {x_j}} )\]

表面上看假设函数增加了若干项， 但经过化简后，可得线性部分最后一项为：

\[\sum\limits_{i = 1}^n {\sum\limits_{j = 1}^n {{\omega _{ij}}{x_i}} {x_j}}  = \sum\limits_{i = 1}^n {\sum\limits_{j = i + 1}^n {\sum\limits_{r = 1}^k {{v_{ir}}{v_{jr}}{x_i}{x_j}} } }  = \sum\limits_{r = 1}^k {({{(\sum\limits_{i = 1}^n {{v_{ir}}{x_i}} )}^2})} \]

反倒是更简单了，为什么不用后面式子计算，而要用更复杂的前者呢？

这是由于ω的每一项并不是独立的，我们使用了vv^T代替了ω

\[{\omega _{ij}} = \sum\limits_{r = 1}^k {{v_{ir}}v_{rj}^T}  = \sum\limits_{r = 1}^k {{v_{ir}}{v_{jr}}} \]

即可认为v_i=(v_i1,v_i2,...,v_ik)代表了特征i与其他特征交叉关系的隐藏特征

如果不去掉交叉特征中的x_ix_i项，则等同于v_i不仅隐含了交叉关系，还隐含了x_i的二次项，不符合我们的初衷

5，TensorFlow概述

TensorFlow即tensor+flow，tensor指数据，以任意维的类型化数组为具体形式，flow则指数据的流向

描述flow有三个层次的概念：

• operation：是最基本的单元，每个operation获得0个或多个tensor，产生0个或多个tensor。operation可能为获取一个常量tf.constant()，获取一个变量tf.get_variable()或者进行一个操作，比如加法tf.add()、乘法tf.multiply()。
• graph：定义了数据流，由tf.Graph()初始化（一般情况下省略）。graph由节点和节点间的连线构成，graph中的节点即operation。graph是对数据流向组成的网络结构的定义，本身并不进行任何计算。
• session：定义了graph运行的环境，由tf.Session()初始化。session中可以运行graph中一个或多个operation，可以在运行前通过feed给operation输入数据，也可以通过返回值获取计算结果。

 

概念上讲，图中所有节点间的连线都是tensor，但总需要有一些数据作为数据来源，描述这些数据来源的组件有以下几种类型：

• 常量：用tf.constant()定义，即运行过程中不会改变的数据，一般用来保存一些固有参数。常量可以被认为是一个operation，输入为空，输出为常数（数组）。
• 变量：理论上广义的变量包含所有除常量以外的所有tensor，这里的变量专指tf.Variable()。值得注意的是，在TensorFlow中，小写开头的才是operation，大写开头的都是类，因此tf.Variable()不是一个operation，而是包含了一系列operation，如初始化、赋值的类成员。每个graph运行前必须初始化所有使用到的变量。
• 占位符：用tf.placeholder()定义，是一种特殊类型的operation，不需要初始化，但需要在运行中通过feed接收数据。

6，LR的TensorFlow实现

首先需要定义输入样本的特征维度和接收输入样本特征和样本标签的占位符。

特征维度可以使用普通python变量指定，此时等同于使用tf.constant创建了一个常量

FEATURE_NUM = 8
#shape参数的列表长度代表占位符的维度
#如果shape不为None，则输入数据维度必须等同于shape参数列表长度
#每一维的大小可以不指定（为None），也可以指定，如果指定则输入数据该维度长度必须与shape指定的一致
x = tf.placeholder(tf.float32, shape=[None, FEATURE_NUM])
#等同于x = tf.placeholder(tf.float32, shape=[None, tf.constant(FEATURE_NUM, dtype=tf.int64)])
y = tf.placeholder(tf.int32, shape=[None])

不可以通过使用placeholder指定特征维度

m = tf.placeholder(tf.int32, shape=None)
x = tf.placeholder(tf.float32, shape=[None, m])  #无法执行
y = tf.placeholder(tf.int32, shape=[None])

虽然placeholder初始shape中不能包含其他tensor，但可以根据输入样本的列数动态指定

x = tf.placeholder(tf.float32, shape=[None, None])  #或者可以不指定维数，输入shape=None，此时可以输入任意维数，任意大小的数据
y = tf.placeholder(tf.int32, shape=[None])
m = tf.shape(x)[1]  #此处不能使用x.get_shape()[1].value，get_shape()获取的结果为静态大小，返回结果为[None,None]

接下来要定义放置LR假设函数中w与b的变量，在定义参数w前，需要定义一个w的初始值。对于LR来说，参数初始值为随机数或者全零都没关系，因为LR的损失函数一定是凸函数（求二阶导数可知），但一般情况下还是习惯采用随机数初始化。

#常用的随机函数有random_uniform、random_normal和truncated_normal
#分别是均匀分布、正态分布、被截断的正态分布
weight_init = tf.truncated_normal(shape=[FEATURE_NUM, 1],mean=0.0,stddev=1.0)  #此处shape初始化只能使用常量，不支持使用tensor
weight = tf.Variable(weight_init)  
bais = tf.Variable([0.0])

为了保持习惯用法，将一维向量y扩展为二维列向量，TensorFlow中扩维可以使用tf.expand_dims()或者tf.reshape()：

y_expand = tf.expand_dims(y, axis=1)  #方法1，axis表示在原若干个维度的间隔中第几个位置插入新维度
y_expand = tf.reshape(y, shape=[tf.shape(x)[0],1])  #方法2，shape为输出的维度值
y_expand = tf.reshape(y, shape=[-1,1])  #方法3，shape参数列表中最多只能有一个参数未知，写为-1

可知：

\[x \in {R^{m*n}},y \in {R^{m*1}},w \in {R^{n*1}},b \in {R^0}\]

回头看LR损失函数的定义

\[l(w,b) =  - \sum\limits_{j = 1}^m {({y^{(j)}}\log (g(b + \sum\limits_{i = 1}^n {{w_i}x_i^{(j)}} )) + (1 - {y^{(j)}})\log (1 - g(b + \sum\limits_{i = 1}^n {{w_i}x_i^{(j)}} ))} ,g(z) = \frac{1}{{1 + {e^{ - z}}}}\]

在TensorFlow中，“+”等同于tf.add()，“-”等同于tf.subtract()，“*”等同于tf.multiply()，“/”等同于tf.div()，所有这些操作是“逐项操作”：

• 如果两个操作数都为标量，则结果为标量
• 如果一个操作数为标量，另一个操作数为向量，则标量会分别与向量每一个元素逐个操作，得到结果
• 如果两个操作数为相同维度，每个维度大小相同的向量，则结果为向量每两个对应位置上的元素的操作得到的结果
• 如果两个操作数维度不相同，如果向量维度数为a、b，并且a>b，则b的各维度大小需要与a的高维度相对应，低维向量在高维向量的低维展开，操作得到结果
• 除上面几种情况，调用不合法

因此，在逻辑回归的损失函数中，我们无需关心“b+”、“1-”等操作，这些标量自然会展开到正确的维度。同样，tf.log()、tf.sigmoid()函数也会对输入向量每一个元素逐个操作，不会改变向量的维度和对应值的位置。

对于式子中的求和Σ，在TensorFlow中有两种处理方式：

• 一种是矩阵乘法tf.matmul()，如矩阵a、b、c关系为c=tf.matmul(a,b)，则可得：

\[{c_{ij}} = \sum\limits_k {{a_{ik}}{b_{kj}}} \]

• 另一种是矩阵压缩函数tf.reduce_sum()（求和）、tf.reduce_mean()（求平均）等函数。这些函数除了第一个参数为要操作的矩阵外，还有几个很有用的参数：axis是要压缩的维度。对于一个二维矩阵来说，压缩维度为0意味着压缩行，生成列向量，压缩维度为1意味着压缩列，生成行向量，axis=None意味着压缩成标量；keepdim指压缩后是否保持原有维度，如果keepdims=True，操作矩阵被压缩的维度长度会保持，但长度一定是1

在损失函数中，先看最外层求和项（以j求和），可以看做y与log函数的结果的逐项相乘求和，因此log函数的输出向量形状必须与y或者y的转置一致。而log函数的输出向量形状取决于内部求和项（以i求和）。可以通过求和式得到合适的运算形式为x*w，也可以通过形状来推理合适的运算形式，即：

\[tf.matmul(x,w) \in {R^{m*1}}\]

然后经过log函数得到的结果与y逐项相乘，再求压缩平均（此处之所以不是求和而是求平均，是为了避免样本数量对损失函数结果产生影响），可以得到假设函数、（负）似然函数、损失函数的运算式：

y_float = tf.to_float(y_expand)
hypothesis = tf.sigmoid(tf.matmul(x, weight) + bais)
likelyhood = -(y_float*tf.log(hypothesis) + (1.0-y_float)*(tf.log(1.0-hypothesis)))
loss = tf.reduce_mean(likelyhood, axis=0)

还可以使用TensorFlow的损失函数计算loss

TensorFlow的损失函数主要有以下几个（下述所有公式假设样本数只有一个，多个的话需要取平均）：

• log_loss：即交叉熵计算，设损失函数l，假设函数h，样本标签y，二分类损失公式如下：

\[l =  - y*\log (h) - (1 - y)\log (1 - h)\]

        多分类损失公式如下，注意可以有多个分类（分类数为c）为正确分类：

\[l =  - \sum\limits_i^c {({y_i}\log ({h_i}) - (1 - {y_i})\log (1 - {h_i}))} \]

        注意不是下式，这里跟一些其他库有区别：

\[l =  - \sum\limits_i^c {{y_i}\log ({h_i})} \]

• sigmoid_cross_entropy：sigmoid+交叉熵，与log_loss的差别在于先对输入预测值求sigmoid函数，然后再计算交叉熵。sigmoid函数如下：

\[g(z) = \frac{1}{{1 + {e^{ - z}}}}\]

• softmax_cross_entropy：softmax+交叉熵，设分类总数为c，softmax函数如下：

\[g({z_i}) = \frac{{{e^{{z_i}}}}}{{\sum\limits_j^c {{e^{{z_j}}}} }}\]

• sparse_softmax_cross_entropy：计算方式与softmax_cross_entropy，但输入样本标签为整数标签，而不是one-hot形式，这是唯一一个输入两个参数各维度大小不同的损失函数

 

• mean_squared_error：均方误差损失函数，损失公式如下：

\[l = {(y - h)^2}\]

• huber_loss：由于均方误差对误差比较大的点相当敏感，为了控制异常点对模型训练的影响，huber loss对误差比较小的样本点使用二次项计算损失，对误差比较大的样本点使用线性方程计算损失，

\[l = \left\{ \begin{array}{l}
\frac{1}{2}{(y - h)^2},|y - h| < = d\\
\frac{1}{2}{d^2} + d(|y - h| - d),|y - h| > d
\end{array} \right.\]

• hinge_loss：hinge loss不太在意预测的有多准确，而主要比较预测正确类别与错误类别的间隔，设损失函数l，假设函数h，样本标签y，二分类损失公式如下：

\[l = \left\{ \begin{array}{l}
\max (0,1 - h),y = 1\\
\max (0,1 + h),y = 0
\end{array} \right.\]

        设正确分类为i，分类总数为c，多分类损失公式如下：

\[l = \sum\limits_{j \ne i}^c {\max (0,1 + {h_j} - {h_i})} \]

在此处，使用tf.losses.log_loss()或者tf.losses.sigmoid_cross_entropy()可以得到与上面运算式完全相同的结果，即：

hypothesis = tf.sigmoid(tf.matmul(x, weight) + bais)
loss = tf.losses.log_loss(y_expand, hypothesis)

或者

hypothesis = tf.matmul(x, weight) + bais
loss = tf.losses.sigmoid_cross_entropy(y_expand, hypothesis)

然后利用TensorFlow的自动微分功能，即优化类之一来定义模型的优化过程：

LEARNING_RATE = 0.02
optimizer = tf.train.GradientDescentOptimizer(LEARNING_RATE)
training_op = optimizer.minimize(loss)

TensorFlow的优化类主要有以下几个：

• GradientDescentOptimizer：最普通的批量梯度下降，令学习速率为η，t代表本次迭代，t+1代表下次迭代，则梯度迭代公式如下：

\[{\theta _{t + 1}} = {\theta _t} - \eta \frac{{\partial l(\theta )}}{{\partial \theta }}\]

• AdagradOptimizer：进行参数迭代的同时记录了每个参数每次迭代的梯度的平方和，下次迭代时梯度与累积平方和的平方根成反比。这样会对低频的参数做较大的更新，对高频的参数做较小的更新，对于稀疏数据表现的更好；但是由于学习速率越来越小，有可能没有到达最低点学习速率就变得很慢了，难以收敛。令s为梯度累积平方和，ε为极小量，t代表本次迭代，t-1代表上次迭代，t+1代表下次迭代，梯度迭代公式如下：

\[{s_t} = {s_{t - 1}} + {(\frac{{\partial l(\theta )}}{{\partial \theta }})^2},{\theta _{t + 1}} = {\theta _t} - \frac{\eta }{{\sqrt {{s_t} + \varepsilon } }}\frac{{\partial l(\theta )}}{{\partial \theta }}\]

• RMSPropOptimizer：为解决AdagradOptimizer后期更新速率过慢的问题，RMSprop使用加权累积平方和替换累积平方和。令m代表梯度加权累积平方和，ε为极小量，β为权重，t代表本次迭代，t-1代表上次迭代，t+1代表下次迭代，梯度迭代公式如下：

\[{m_t} = \beta {m_{t - 1}} + (1 - \beta ){(\frac{{\partial l(\theta )}}{{\partial \theta }})^2},{\theta _{t + 1}} = {\theta _t} - \frac{\eta }{{\sqrt {{m_t} + \varepsilon } }}\frac{{\partial l(\theta )}}{{\partial \theta }}\]

• MomentumOptimizer：多一个必须参数“动量速率”，每次迭代时参考前一次迭代的“动量”，在迭代中方向不变的维度做较大的更新，迭代中方向反复改变的维度做较小的更新。适用于在不同维度梯度差距很大的情况，更新不会在小梯度方向反复震荡。令γ代表动量速率，t代表本次迭代，t-1代表上次迭代，t+1代表下次迭代，梯度迭代公式如下：

\[{v_t} = \gamma {v_{t - 1}} + \eta \frac{{\partial l(\theta )}}{{\partial \theta }},{\theta _{t + 1}} = {\theta _t} - {v_t}\]

•  AdamOptimizer：综合了MomentumOptimizer和RMSPropOptimizer，既包含动量（一次项）部分也包含衰减（两次项）部分。令f代表动量，g代表衰减，β₁、β₂为权重，t代表本次迭代，t-1代表上次迭代，t+1代表下次迭代：

\[\begin{array}{l}
{f_t} = {\beta _2}{f_{t - 1}} + (1 - {\beta _2})\frac{{\partial l(\theta )}}{{\partial \theta }}\\
{g_t} = {\beta _1}{g_{t - 1}} + (1 - {\beta _1}){(\frac{{\partial l(\theta )}}{{\partial \theta }})^2}
\end{array}\]

        注意此处第一次迭代的梯度值和速率值偏小，为了校正偏差，计算校正后的动量f′和衰减g′，再计算迭代公式：

\[\begin{array}{l}
{{f'}_t} = \frac{{{f_t}}}{{1 - {\beta _1}}}\\
{{g'}_t} = \frac{{{g_t}}}{{1 - {\beta _2}}}\\
{\theta _{t + 1}} = {\theta _t} - \frac{\eta }{{\sqrt {{{g'}_t} + \varepsilon } }}{{f'}_t}
\end{array}\]

最后还需要一个预测操作和准确率判断操作：

THRESHOLD = 0.5  
predictions = tf.sign(hypothesis-THRESHOLD)  #符号函数，判断向量对应位置的符号，输出对应位置为-1、0或1
labels = tf.sign(y_float-THRESHOLD)
corrections = tf.equal(predictions, labels)  #比较向量对应位置的两个值，相等则输出对应位置为True
accuracy = tf.reduce_mean(tf.cast(corrections, tf.float32))

也可以使用：

THRESHOLD = 0.5
predictions = tf.to_int32(hypothesis-THRESHOLD)
corrections = tf.equal(predictions, y_expand)
accuracy = tf.reduce_mean(tf.cast(corrections, tf.float32))

完整代码（包括测试运行代码）如下：

tf.reset_default_graph()  #清空Graph

FEATURE_NUM = 8  #特征数量
with tf.name_scope("input"):
    x = tf.placeholder(tf.float32, shape=[None, FEATURE_NUM])
    y = tf.placeholder(tf.int32, shape=[None])

with tf.name_scope("lr"):
    weight_init = tf.truncated_normal(shape=[FEATURE_NUM, 1],mean=0.0,stddev=1.0)
    weight = tf.Variable(weight_init)  
    bais = tf.Variable([0.0])
    
    y_expand = tf.expand_dims(y, axis=1)
    
    hypothesis = tf.sigmoid(tf.matmul(x, weight) + bais)

with tf.name_scope("loss"):
    y_float = tf.to_float(y_expand)
    likelyhood = -(y_float*tf.log(hypothesis) + (1.0-y_float)*(tf.log(1.0-hypothesis)))
    loss = tf.reduce_mean(likelyhood, axis=0)

LEARNING_RATE = 0.02  #学习速率
with tf.name_scope("train"):
    optimizer = tf.train.GradientDescentOptimizer(LEARNING_RATE)
    training_op = optimizer.minimize(loss)
    
THRESHOLD = 0.5  #判断门限
with tf.name_scope("eval"):
    predictions = tf.sign(hypothesis-THRESHOLD)
    labels = tf.sign(y_float-THRESHOLD)
    corrections = tf.equal(predictions, labels)
    accuracy = tf.reduce_mean(tf.cast(corrections, tf.float32))

init = tf.global_variables_initializer()  #初始化所有变量

EPOCH = 10  #迭代次数
with tf.Session() as sess:
    sess.run(init)
    for i in range(EPOCH):
        _training_op, _loss = sess.run([training_op, loss], feed_dict={x: np.random.rand(10,8), y: np.random.randint(2,size=10)})
        _accuracy = sess.run([accuracy], feed_dict={x: np.random.rand(5,8), y: np.random.randint(2,size=5)})
        print "epoch:", i, _loss, _accuracy

7，FM的TensorFlow实现

同样，先定义变量占位符：

FEATURE_NUM = 8  #特征数量
x = tf.placeholder(tf.float32, shape=[None, FEATURE_NUM])
y = tf.placeholder(tf.int32, shape=[None])
y_exband = tf.expand_dims(y, axis=1)

然后实现假设函数：

\[\begin{array}{l}
l(b,w,v) = - \sum\limits_{j = 1}^m {({y^{(j)}}\log (h(w,b){|_{x = {x^{(j)}}}}) + (1 - {y^{(j)}})\log (1 - h(w,b){|_{x = {x^{(j)}}}}))} \\
h(w,b){|_{x = {x^{(j)}}}} = g(b + \sum\limits_{i = 1}^n {{w_i}x_i^{(j)}} + \frac{1}{2}\sum\limits_{r = 1}^k {({{(\sum\limits_{i = 1}^n {{v_{ir}}x_i^{(j)}} )}^2} - \sum\limits_{i = 1}^n {{{({v_{ir}}x_i^{(j)})}^2}} )} )
\end{array}\]

先定义模型权重：

HIDDEN_N = 5
bais = tf.Variable([0.0])
weight = tf.Variable(tf.random_normal([FEATURE_N, 1], 0.0, 1.0))
weight_mix = tf.Variable(tf.random_normal([FEATURE_N, HIDDEN_N], 0.0, 1.0))

一次项和零次项比较好计算，跟LR一样，但二次项就比较麻烦了。这里的计算过程，一些广为流传的博客上代码都是错的。

首先可知：

\[x \in {R^{m*n}},v \in {R^{n*k}},y \in {R^{m*1}}\]

先尝试矩阵乘法降维，可得：

\[xv \in {R^{m*k}},{(xv)_{ab}} = \sum\limits_{i = 1}^n {x_i^{(a)}{v_{ib}}} \]

可见正好与二次项的第一项形式一致，可以通过以下方式计算第一项，一个m*k的矩阵，得到结果后可以通过先压缩第1维再压缩第0维，得到损失函数结果（0维标量）

x_weight_mix = tf.matmul(x, weight_mix)
x_weight_mix_square = tf.square(x_weight_mix)

但观察二次项的第二项，就没这么简单了，去除维度无关的加法、减法、log等操作的影响，第二项等同于计算如下式子：

\[\sum\limits_{j = 1}^m {\sum\limits_{r = 1}^k {\sum\limits_{i = 1}^n {{{(x_i^{(j)}{v_{ir}})}^2}} } } \]

即需要得到每个j、r、i的组合，再平方，再依次求和，可知这必须是一个三维矩阵才能完成运算。因此，如果能设计一个三维矩阵d满足如下就条件，就好计算了：

\[{d_{jri}} = x_i^{(j)}{v_{ir}},d \in {R^{m*k*n}}\]

由于三维矩阵d的每项不涉及求和，只有两个数相乘，因此最直接的方式是扩展x、v的维度与d相同，然后直接逐项相乘，为方便计算，调整三维矩阵d的定义：

\[{d_{jir}} = x_i^{(j)}{v_{ir}},d \in {R^{m*n*k}}\]

这样，x在第2维扩展，v在第0为扩展，逐项相乘结果即为d：

\[\begin{array}{l}
{d_{jir}} = {x_{(new),}}_{ir}^{(j)}{v_{(new)}}_{,jir},d \in {R^{m*n*k}}\\
{x_{(new),}}_{ir}^{(j)} = x_i^{(j)},{x_{(new)}} \in {R^{m*n*k}}\\
{v_{(new)}}_{,jir} = {v_{ir}},{v_{(new)}} \in {R^{m*n*k}}
\end{array}\]

这里需要用到tf.tile()和tf.reshape()，tf.tile()的功能是将矩阵的某一维（或者n维）的数据重复数次，tf.reshape()可以将两维向量变为三维矩阵。这里需要注意，如果先tile再reshape，其实是很容易出错的，例如将x的第2维扩展k次，可得：

\[{x_{(tile)}} \in {R^{m*(n*k)}}\]

但是，到底x_(tile)经过“合适的”reshape升维后是m*n*k维矩阵还是m*k*n维矩阵，是一个难以直接判断的问题，因此，最好先reshape（或者expand_dims），再tile。二次项的第二项计算方式如下：

SAMPLE_N = tf.shape(x)[0]
x_ = tf.reshape(x, [SAMPLE_N,FEATURE_N,1]) #SAMPLE_N*FEATURE_N*1
x__ = tf.tile(x_, [1,1,HIDDEN_N]) #SAMPLE_N*FEATURE_N*HIDDEN_N
w_ = tf.reshape(weight_mix, [1,FEATURE_N,HIDDEN_N]) #1*FEATURE_N*HIDDEN_N
w__ = tf.tile(x_, [SAMPLE_N,1,1]) #SAMPLE_N*FEATURE_N*HIDDEN_N
embeddings = tf.multiply(x__, w__) #SAMPLE_N*FEATURE_N*HIDDEN_N
embeddings_square = tf.square(embeddings) #SAMPLE_N*FEATURE_N*HIDDEN_N
embeddings_square_sum = tf.reduce_sum(embeddings_square, 1) #SAMPLE_N*HIDDEN_N

此时，可以发现，构造出embeddings矩阵后，二次项的第一项也可以通过embeddings来计算：

embeddings_sum = tf.reduce_sum(embeddings, 1) #SAMPLE_N*HIDDENT_N
embeddings_sum_square = tf.square(embeddings_sum) #SAMPLE_N*HIDDENT_N

然后将损失函数的零次项、一次项、二次项组合起来，可得损失函数：

z = bais + tf.matmul(X, weight) + 1.0/2.0*tf.reduce_sum(tf.subtract(embeddings_sum_square,embeddings_square_sum), 1, keepdims=True)
hypothesis  = tf.sigmoid(z)

这里可以发现，z的计算公式中包含大量二次项的求和，绝对值很可能会比较大，此时再经过sigmoid函数，由于浮点数表达精度有限，结果可能会近似等于1或等于0，造成模型无法更新，因此此处需要使用tf.clip_by_value()对z进行截断，防止超出精度：

z = bais + tf.matmul(x, weight) + 1.0/2.0*tf.reduce_sum(tf.subtract(embeddings_sum_square,embeddings_square_sum), 1, keepdims=True)
z_ = tf.clip_by_value(z,-4.0,4.0)
hypothesis  = tf.sigmoid(z_)

接下来优化过程与判断过程与LR基本一样，可得完整代码（包括测试运行代码）如下：

tf.reset_default_graph()  #清空Graph

FEATURE_NUM = 8  #特征数量
with tf.name_scope("input"):
    x = tf.placeholder(tf.float32, shape=[None, FEATURE_NUM])
    y = tf.placeholder(tf.int32, shape=[None])
    y_expand = tf.expand_dims(y, axis=1)

HIDDEN_NUM = 5  #隐藏特征维度
with tf.name_scope("fm"):
    bais = tf.Variable([0.0])
    weight = tf.Variable(tf.random_normal([FEATURE_NUM, 1], 0.0, 1.0))
    weight_mix = tf.Variable(tf.random_normal([FEATURE_NUM, HIDDEN_NUM], 0.0, 1.0))
    
    SAMPLE_NUM = tf.shape(x)[0]  #获取样本数
    
    x_ = tf.reshape(x, [SAMPLE_NUM,FEATURE_NUM,1]) #SAMPLE_NUM*FEATURE_NUM*1
    x__ = tf.tile(x_, [1,1,HIDDEN_NUM]) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
    w_ = tf.reshape(weight_mix, [1,FEATURE_NUM,HIDDEN_NUM]) #1*FEATURE_NUM*HIDDEN_NUM
    w__ = tf.tile(w_, [SAMPLE_NUM,1,1]) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
    
    embeddings = tf.multiply(x__, w__) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
    embeddings_sum = tf.reduce_sum(embeddings, 1) #SAMPLE_NUM*HIDDEN_NUM
    embeddings_sum_square = tf.square(embeddings_sum) #SAMPLE_NUM*HIDDEN_NUM
    embeddings_square = tf.square(embeddings) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
    embeddings_square_sum = tf.reduce_sum(embeddings_square, 1) #SAMPLE_NUM*HIDDEN_NUM
    
    z = bais + tf.matmul(x, weight) + 1.0/2.0*tf.reduce_sum(tf.subtract(embeddings_sum_square,embeddings_square_sum), 1, keepdims=True)
    z_ = tf.clip_by_value(z,-4.0,4.0)
    hypothesis  = tf.sigmoid(z_)

with tf.name_scope("loss"):
    loss = tf.losses.log_loss(y_expand, hypothesis)

LEARNING_RATE = 0.02  #学习速率
with tf.name_scope("train"):
    optimizer = tf.train.GradientDescentOptimizer(LEARNING_RATE)
    training_op = optimizer.minimize(loss)
    
THRESHOLD = 0.5  #判断门限
with tf.name_scope("eval"):
    predictions = tf.to_int32(hypothesis-THRESHOLD)
    corrections = tf.equal(predictions, y_expand)
    accuracy = tf.reduce_mean(tf.cast(corrections, tf.float32))

init = tf.global_variables_initializer()  #初始化所有变量

EPOCH = 10  #迭代次数
with tf.Session() as sess:
    sess.run(init)
    for i in range(EPOCH):
        _training_op, _loss = sess.run([training_op, loss], feed_dict={x: np.random.rand(10,8), y: np.random.randint(2,size=10)})
        _accuracy = sess.run([accuracy], feed_dict={x: np.random.rand(5,8), y: np.random.randint(2,size=5)})
        print "epoch:", i, _loss, _accuracy

但是这时可以发现，模型本身已经比较复杂了，其中包含了若干变量，在嵌入工程，或者作为模型中的一个模块嵌入网络中时，很有可能会造成变量名称的混乱，并且也不利于生成明晰的网络结构，因此，这里考虑将FM构建网络部分包裹成函数或者类：

#FM模型
class FmModel(object):
    def __init__(self, x, y, feature_num, hidden_num):
        self.x = x
        self.y = y
        self.feature_num = feature_num  #获取特征数，这个值要建Variable，所以不能动态获取
        self.sample_num = tf.shape(x)[0]  #获取样本数
        self.hidden_num = hidden_num  #获取隐藏特征维度
        
        self.bais = tf.Variable([0.0])
        self.weight = tf.Variable(tf.random_normal([self.feature_num, 1], 0.0, 1.0))
        self.weight_mix = tf.Variable(tf.random_normal([self.feature_num, self.hidden_num], 0.0, 1.0))
        
        x_ = tf.reshape(self.x, [self.sample_num,self.feature_num,1]) #SAMPLE_NUM*FEATURE_NUM*1
        x__ = tf.tile(x_, [1,1,self.hidden_num]) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
        w_ = tf.reshape(self.weight_mix, [1,self.feature_num,self.hidden_num]) #1*FEATURE_NUM*HIDDEN_NUM
        w__ = tf.tile(w_, [self.sample_num,1,1]) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
        
        embeddings = tf.multiply(x__, w__) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
        embeddings_sum = tf.reduce_sum(embeddings, 1) #SAMPLE_NUM*HIDDEN_NUM
        embeddings_sum_square = tf.square(embeddings_sum) #SAMPLE_NUM*HIDDEN_NUM
        embeddings_square = tf.square(embeddings) #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
        embeddings_square_sum = tf.reduce_sum(embeddings_square, 1) #SAMPLE_NUM*HIDDEN_NUM
        
        z = self.bais + tf.matmul(self.x, self.weight) + 1.0/2.0*tf.reduce_sum(tf.subtract(embeddings_sum_square,embeddings_square_sum), 1, keepdims=True)
        z_ = tf.clip_by_value(z,-4.0,4.0)
        
        self.hypothesis  = tf.sigmoid(z_)
        y_expand = tf.expand_dims(self.y, axis=1)
        self.loss = tf.losses.log_loss(y_expand, self.hypothesis)

这里需要注意，无论使用类也好，函数也好，其中很多变量从全局变量变为了局部变量，函数执行完变量就被销毁了。但只要执行完一次，网络就已经建立起来了，变量的生存期并不影响网络节点的生存期。这时又涉及到另一个问题：如何在Session中想要得到这些类或函数中定义的网络节点的值。

TensorFlow中会话运行网络节点的方式有两种：

• 使用网络建立过程中，仍处于生存期内的python变量。这样做有两个局限性，一种情况是python变量生存期外就无法使用，另一情况是如果模型加载自已经训练好的网络，没有对应节点的python变量。
• 使用网络中的节点名，赋给python变量。TensorFlow中的Tensor都可以在定义时添加一个附加参数name，可以通过tf.get_default_graph().get_tensor_by_name()来获取网络中任何的节点。

比如上例中如果要获取embeddings变量的值，可以如此操作：

#FM模型
class FmModel(object):
    def __init__(self, x, y, feature_num, hidden_num):
        ...
        embeddings = tf.multiply(x__, w__, name="embeddings") #SAMPLE_NUM*FEATURE_NUM*HIDDEN_NUM
        ...
        
HIDDEN_NUM = 5  #隐藏特征维度
with tf.name_scope("fm"):
    fm = FmModel(x, y, FEATURE_NUM, HIDDEN_NUM)
    
embeddings = tf.get_default_graph().get_tensor_by_name("fm/embeddings:0")  #注意name_scope与后面的冒号和序号

使用FmModel类的完整代码（包括测试运行代码）如下：

tf.reset_default_graph()  #清空Graph

FEATURE_NUM = 8  #特征数量
with tf.name_scope("input"):
    x = tf.placeholder(tf.float32, shape=[None, FEATURE_NUM])
    y = tf.placeholder(tf.int32, shape=[None])

HIDDEN_NUM = 5  #隐藏特征维度
with tf.name_scope("fm"):
    fm = FmModel(x, y, FEATURE_NUM, HIDDEN_NUM)

LEARNING_RATE = 0.02  #学习速率
with tf.name_scope("train"):
    optimizer = tf.train.GradientDescentOptimizer(LEARNING_RATE)
    training_op = optimizer.minimize(fm.loss)
    
THRESHOLD = 0.5  #判断门限
with tf.name_scope("eval"):
    y_expand = tf.expand_dims(y, axis=1)
    predictions = tf.to_int32(fm.hypothesis-THRESHOLD)
    corrections = tf.equal(predictions, fm.y_expand)
    accuracy = tf.reduce_mean(tf.cast(corrections, tf.float32))

init = tf.global_variables_initializer()  #初始化所有变量

EPOCH = 10  #迭代次数
with tf.Session() as sess:
    sess.run(init)
    for i in range(EPOCH):
        _training_op, _loss = sess.run([training_op, fm.loss], feed_dict={x: np.random.rand(10,8), y: np.random.randint(2,size=10)})
        _accuracy = sess.run([accuracy], feed_dict={x: np.random.rand(5,8), y: np.random.randint(2,size=5)})
        print "epoch:", i, _loss, _accuracy

7，稀疏FM的TensorFlow实现

上述TensorFlow实现中使用了不少二维甚至三维的矩阵运算，在生产环境中，使用FM模型的样本数据往往是维数巨大并且稀疏的，这时直接使用矩阵乘法会占用大量的系统资源、浪费大量的运算时间。

针对稀疏矩阵乘法，TensorFlow有一种专门的方式简化运算过程。此时要用到tf.nn.embedding_lookup()函数，这个函数的输入一个参数矩阵w和一个稀疏索引矩阵i，作用是在在参数矩阵w的第一维中查找索引矩阵i中各值对应的“行”，抽取出索引值对应的“一行”参数，放到索引矩阵中对应位置，组合为有效特征对应的参数矩阵。可知索引矩阵的值必须小于w第一维的大小，结果维数等于w的维数与i的维数之和减一。

由于有效特征矩阵会参与矩阵运算，矩阵运算中不能出现每行长度不同的情况，所以这里有一点局限性，就是每个样本的有效特征数必须是一样的，如果多了需要截取，如果少了需要补零。

针对稀疏特征的FM完整代码（包括测试运行代码）如下：

#FM模型
class FmModel(object):
    def __init__(self, i, x, y, feature_num, valid_num, hidden_num):
        self.i = i
        self.x = x
        self.y = y
        self.feature_num = feature_num  #获取特征数，这个值要建Variable，所以不能动态获取
        self.valid_num = valid_num  #获取有效特征数，这个值要建Variable，所以不能动态获取
        self.sample_num = tf.shape(x)[0]  #获取样本数
        self.hidden_num = hidden_num  #获取隐藏特征维度
        
        self.bais = tf.Variable([0.0])
        self.weight = tf.Variable(tf.random_normal([self.feature_num, 1], 0.0, 1.0))
        self.weight_mix = tf.Variable(tf.random_normal([self.feature_num, self.hidden_num], 0.0, 1.0))
        
        x_ = tf.reshape(self.x, [self.sample_num,self.valid_num,1]) #SAMPLE_NUM*VALID_NUM*1
        w_ = tf.nn.embedding_lookup(self.weight, self.i) #SAMPLE_NUM*VALID_NUM*1
        
        expressings = tf.multiply(x_, w_) #SAMPLE_NUM*VALID_NUM*1
        expressings_reduce = tf.reshape(self.x, [self.sample_num,self.valid_num]) #SAMPLE_NUM*VALID_NUM
        
        x__ = tf.tile(x_, [1,1,self.hidden_num]) #SAMPLE_NUM*VALID_NUM*HIDDEN_NUM
        w__ = tf.nn.embedding_lookup(self.weight_mix, self.i) #SAMPLE_NUM*VALID_NUM*HIDDEN_NUM
        
        embeddings = tf.multiply(x__, w__) #SAMPLE_NUM*VALID_NUM*HIDDEN_NUM
        embeddings_sum = tf.reduce_sum(embeddings, 1) #SAMPLE_NUM*HIDDEN_NUM
        embeddings_sum_square = tf.square(embeddings_sum) #SAMPLE_NUM*HIDDEN_NUM
        embeddings_square = tf.square(embeddings) #SAMPLE_NUM*VALID_NUM*HIDDEN_NUM
        embeddings_square_sum = tf.reduce_sum(embeddings_square, 1) #SAMPLE_NUM*HIDDEN_NUM
        
        z = self.bais + \
            tf.reduce_sum(expressings_reduce, 1, keepdims=True) + \
            1.0/2.0*tf.reduce_sum(tf.subtract(embeddings_sum_square,embeddings_square_sum), 1, keepdims=True)
        z_ = tf.clip_by_value(z,-4.0,4.0)
        
        self.hypothesis  = tf.sigmoid(z_)
        
        self.y_expand = tf.expand_dims(self.y, axis=1)
       
        self.loss = tf.losses.log_loss(self.y_expand, self.hypothesis)

tf.reset_default_graph()  #清空Graph

VALID_NUM = 8  #有效特征数量
with tf.name_scope("input"):
    i = tf.placeholder(tf.int32, shape=[None, VALID_NUM])
    x = tf.placeholder(tf.float32, shape=[None, VALID_NUM])
    y = tf.placeholder(tf.int32, shape=[None])
    y_expand = tf.expand_dims(y, axis=1)

FEATURE_NUM = 20  #特征数量
HIDDEN_NUM = 5  #隐藏特征维度
with tf.name_scope("fm"):
    fm = FmModel(i, x, y, FEATURE_NUM, VALID_NUM, HIDDEN_NUM)

LEARNING_RATE = 0.02  #学习速率
with tf.name_scope("train"):
    optimizer = tf.train.GradientDescentOptimizer(LEARNING_RATE)
    training_op = optimizer.minimize(fm.loss)
    
THRESHOLD = 0.5  #判断门限
with tf.name_scope("eval"):
    predictions = tf.to_int32(fm.hypothesis-THRESHOLD)
    corrections = tf.equal(predictions, fm.y_expand)
    accuracy = tf.reduce_mean(tf.cast(corrections, tf.float32))

init = tf.global_variables_initializer()  #初始化所有变量

EPOCH = 10  #迭代次数
with tf.Session() as sess:
    sess.run(init)
    for epoch in range(EPOCH):
        _training_op, _loss = sess.run([training_op, fm.loss],
            feed_dict={i: np.array([np.random.choice(20, 8) for cnt in range(10)]), x: np.random.rand(10,8), y: np.random.randint(2,size=10)})
        _accuracy = sess.run([accuracy],
            feed_dict={i: np.array([np.random.choice(20, 8) for cnt in range(5)]), x: np.random.rand(5,8), y: np.random.randint(2,size=5)})
        print "epoch:", epoch, _loss, _accuracy

 

 

 

 

参考文献：

https://www.cnblogs.com/pinard/p/6370127.html