Python如何實(shí)現(xiàn)FM算法

這篇文章主要介紹Python如何實(shí)現(xiàn)FM算法，文中介紹的非常詳細(xì)，具有一定的參考價(jià)值，感興趣的小伙伴們一定要看完！

1. 什么是FM？

FM即Factor Machine，因子分解機(jī)。

2. 為什么需要FM？

1、特征組合是許多機(jī)器學(xué)習(xí)建模過程中遇到的問題，如果對特征直接建模，很有可能會(huì)忽略掉特征與特征之間的關(guān)聯(lián)信息，因此，可以通過構(gòu)建新的交叉特征這一特征組合方式提高模型的效果。

2、高維的稀疏矩陣是實(shí)際工程中常見的問題，并直接會(huì)導(dǎo)致計(jì)算量過大，特征權(quán)值更新緩慢。試想一個(gè)10000*100的表，每一列都有8種元素，經(jīng)過one-hot獨(dú)熱編碼之后，會(huì)產(chǎn)生一個(gè)10000*800的表。因此表中每行元素只有100個(gè)值為1，700個(gè)值為0。

而FM的優(yōu)勢就在于對這兩方面問題的處理。首先是特征組合，通過對兩兩特征組合，引入交叉項(xiàng)特征，提高模型得分；其次是高維災(zāi)難，通過引入隱向量（對參數(shù)矩陣進(jìn)行矩陣分解），完成對特征的參數(shù)估計(jì)。

3. FM用在哪？

我們已經(jīng)知道了FM可以解決特征組合以及高維稀疏矩陣問題，而實(shí)際業(yè)務(wù)場景中，電商、豆瓣等推薦系統(tǒng)的場景是使用最廣的領(lǐng)域，打個(gè)比方，小王只在豆瓣上瀏覽過20部電影，而豆瓣上面有20000部電影，如果構(gòu)建一個(gè)基于小王的電影矩陣，毫無疑問，里面將有199980個(gè)元素全為0。而類似于這樣的問題就可以通過FM來解決。

4. FM長什么樣？

在展示FM算法前，我們先回顧一下最常見的線性表達(dá)式：

其中w0為初始權(quán)值，或者理解為偏置項(xiàng)，wi為每個(gè)特征xi對應(yīng)的權(quán)值?？梢钥吹剑@種線性表達(dá)式只描述了每個(gè)特征與輸出的關(guān)系。

FM的表達(dá)式如下，可觀察到，只是在線性表達(dá)式后面加入了新的交叉項(xiàng)特征及對應(yīng)的權(quán)值。

5. FM交叉項(xiàng)的展開

5.1 尋找交叉項(xiàng)

FM表達(dá)式的求解核心在于對交叉項(xiàng)的求解。下面是很多人用來求解交叉項(xiàng)的展開式，對于第一次接觸FM算法的人來說可能會(huì)有疑惑，不知道公式怎么展開的，接下來筆者會(huì)手動(dòng)推導(dǎo)一遍。

設(shè)有3個(gè)變量（特征）x1 x2 x3，每一個(gè)特征的隱變量分別為v1=(1 2 3)、v2=(4 5 6)、v3=(1 2 1)，即：

設(shè)交叉項(xiàng)所組成的權(quán)矩陣W為對稱矩陣，之所以設(shè)為對稱矩陣是因?yàn)閷ΨQ矩陣有可以用向量乘以向量轉(zhuǎn)置替代的性質(zhì)。
那么W=VVT，即

所以：

實(shí)際上，我們應(yīng)該考慮的交叉項(xiàng)應(yīng)該是排除自身組合的項(xiàng)，即對于x1x1、x2x2、x3x3不認(rèn)為是交叉項(xiàng)，那么真正的交叉項(xiàng)為x1x2、x1x3、x2x1、x2x3、x3x1、x3x2。

去重后，交叉項(xiàng)即x1x2、x1x3、x2x3。這也是公式中1/2出現(xiàn)的原因。

5.2 交叉項(xiàng)權(quán)值轉(zhuǎn)換

對交叉項(xiàng)有了基本了解后，下面將進(jìn)行公式的分解，還是以n=3為例，

所以：

wij可記作或，這取決于vi是1*3 還是3*1 向量。

5.3 交叉項(xiàng)展開式

上面的例子是對3個(gè)特征做的交叉項(xiàng)推導(dǎo)，因此對具有n個(gè)特征，F(xiàn)M的交叉項(xiàng)公式就可推廣為：

我們還可以進(jìn)一步分解：

所以FM算法的交叉項(xiàng)最終可展開為：

5.4隱向量v就是embedding vector?

假設(shè)訓(xùn)練數(shù)據(jù)集dataMatrix的shape為（20000，9)，取其中一行數(shù)據(jù)作為一條樣本i，那么樣本i 的shape為（1，9），同時(shí)假設(shè)隱向量vi的shape為（9，8）（注：8為自定義值，代表embedding vector的長度）

所以5.3小節(jié)中的交叉項(xiàng)可以表示為：

sum((inter_1)^2 - (inter_2)^2)/2

其中：

inter_1 =i*v shape為（1，8）

inter_2 =np.multiply(i)*np.multiply(v) shape為（1，8）

可以看到，樣本i 經(jīng)過交叉項(xiàng)中的計(jì)算后，得到向量shape為（1，8）的inter_1和inter_2。

由于維度變低，所以此計(jì)算過程可以近似認(rèn)為在交叉項(xiàng)中對樣本i 進(jìn)行了embedding vector轉(zhuǎn)換。

故，我們需要對之前的理解進(jìn)行修正：

1. 我們口中的隱向量vi實(shí)際上是一個(gè)向量組，其形狀為（輸入特征One-hot后的長度，自定義長度）；

2. 隱向量vi代表的并不是embedding vector，而是在對輸入進(jìn)行embedding vector的向量組，也可理解為是一個(gè)權(quán)矩陣；

3. 由輸入i*vi得到的向量才是真正的embedding vector。

具體可以結(jié)合第7節(jié)點(diǎn)的代碼實(shí)現(xiàn)進(jìn)行理解。

6. 權(quán)值求解

利用梯度下降法，通過求損失函數(shù)對特征（輸入項(xiàng)）的導(dǎo)數(shù)計(jì)算出梯度，從而更新權(quán)值。設(shè)m為樣本個(gè)數(shù)，θ為權(quán)值。

如果是回歸問題，損失函數(shù)一般是均方誤差（MSE）:

所以回歸問題的損失函數(shù)對權(quán)值的梯度（導(dǎo)數(shù)）為：

如果是二分類問題，損失函數(shù)一般是logit loss：

其中，表示的是階躍函數(shù)Sigmoid。

所以分類問題的損失函數(shù)對權(quán)值的梯度（導(dǎo)數(shù)）為：

相應(yīng)的，對于常數(shù)項(xiàng)、一次項(xiàng)、交叉項(xiàng)的導(dǎo)數(shù)分別為:

7. FM算法的Python實(shí)現(xiàn)

FM算法的Python實(shí)現(xiàn)流程圖如下：

我們需要注意以下四點(diǎn)：

1. 初始化參數(shù)，包括對偏置項(xiàng)權(quán)值w0、一次項(xiàng)權(quán)值w以及交叉項(xiàng)輔助向量的初始化；

2. 定義FM算法；

3. 損失函數(shù)梯度的定義；

4. 利用梯度下降更新參數(shù)。

下面的代碼片段是以上四點(diǎn)的描述，其中的loss并不是二分類的損失loss，而是分類loss的梯度中的一部分：

loss = self.sigmoid(classLabels[x] * p[0, 0]) -1

實(shí)際上，二分類的損失loss的梯度可以表示為：

gradient = (self.sigmoid(classLabels[x] * p[0, 0]) -1)*classLabels[x]*p_derivative

其中 p_derivative 代表常數(shù)項(xiàng)、一次項(xiàng)、交叉項(xiàng)的導(dǎo)數(shù)（詳見本文第6小節(jié)）。

FM算法代碼片段

# 初始化參數(shù)
    w = zeros((n, 1)) # 其中n是特征的個(gè)數(shù)
    w_0 = 0.
    v = normalvariate(0, 0.2) * ones((n, k))
    for it in range(self.iter): # 迭代次數(shù)
      # 對每一個(gè)樣本，優(yōu)化
      for x in range(m):
        # 這邊注意一個(gè)數(shù)學(xué)知識：對應(yīng)點(diǎn)積的地方通常會(huì)有sum，對應(yīng)位置積的地方通常都沒有，詳細(xì)參見矩陣運(yùn)算規(guī)則，本處計(jì)算邏輯在：http://blog.csdn.net/google19890102/article/details/45532745
        # xi·vi,xi與vi的矩陣點(diǎn)積
        inter_1 = dataMatrix[x] * v
        # xi與xi的對應(yīng)位置乘積  與  xi^2與vi^2對應(yīng)位置的乘積  的點(diǎn)積
        inter_2 = multiply(dataMatrix[x], dataMatrix[x]) * multiply(v, v) # multiply對應(yīng)元素相乘
        # 完成交叉項(xiàng),xi*vi*xi*vi - xi^2*vi^2
        interaction = sum(multiply(inter_1, inter_1) - inter_2) / 2.
        # 計(jì)算預(yù)測的輸出
        p = w_0 + dataMatrix[x] * w + interaction
        print('classLabels[x]:',classLabels[x])
        print('預(yù)測的輸出p:', p)
        # 計(jì)算sigmoid(y*pred_y)-1準(zhǔn)確的說不是loss，原作者這邊理解的有問題，只是作為更新w的中間參數(shù)，這邊算出來的是越大越好，而下面卻用了梯度下降而不是梯度上升的算法在
        loss = self.sigmoid(classLabels[x] * p[0, 0]) - 1
        if loss >= -1:
          loss_res = '正方向 '
        else:
          loss_res = '反方向'
        # 更新參數(shù)
        w_0 = w_0 - self.alpha * loss * classLabels[x]
        for i in range(n):
          if dataMatrix[x, i] != 0:
            w[i, 0] = w[i, 0] - self.alpha * loss * classLabels[x] * dataMatrix[x, i]
            for j in range(k):
              v[i, j] = v[i, j] - self.alpha * loss * classLabels[x] * (
                  dataMatrix[x, i] * inter_1[0, j] - v[i, j] * dataMatrix[x, i] * dataMatrix[x, i])

FM算法完整實(shí)現(xiàn)

# -*- coding: utf-8 -*-

from __future__ import division
from math import exp
from numpy import *
from random import normalvariate # 正態(tài)分布
from sklearn import preprocessing
import numpy as np

'''
  data : 數(shù)據(jù)的路徑
  feature_potenital : 潛在分解維度數(shù)
  alpha ： 學(xué)習(xí)速率
  iter ： 迭代次數(shù)
  _w,_w_0,_v ： 拆分子矩陣的weight
  with_col : 是否帶有columns_name
  first_col : 首列有價(jià)值的feature的index
'''


class fm(object):
  def __init__(self):
    self.data = None
    self.feature_potential = None
    self.alpha = None
    self.iter = None
    self._w = None
    self._w_0 = None
    self.v = None
    self.with_col = None
    self.first_col = None

  def min_max(self, data):
    self.data = data
    min_max_scaler = preprocessing.MinMaxScaler()
    return min_max_scaler.fit_transform(self.data)

  def loadDataSet(self, data, with_col=True, first_col=2):
    # 我就是閑的蛋疼，明明pd.read_table()可以直接度，非要搞這樣的，顯得代碼很長，小數(shù)據(jù)下完全可以直接讀嘛，唉～
    self.first_col = first_col
    dataMat = []
    labelMat = []
    fr = open(data)
    self.with_col = with_col
    if self.with_col:
      N = 0
      for line in fr.readlines():
        # N=1時(shí)干掉列表名
        if N > 0:
          currLine = line.strip().split()
          lineArr = []
          featureNum = len(currLine)
          for i in range(self.first_col, featureNum):
            lineArr.append(float(currLine[i]))
          dataMat.append(lineArr)
          labelMat.append(float(currLine[1]) * 2 - 1)
        N = N + 1
    else:
      for line in fr.readlines():
        currLine = line.strip().split()
        lineArr = []
        featureNum = len(currLine)
        for i in range(2, featureNum):
          lineArr.append(float(currLine[i]))
        dataMat.append(lineArr)
        labelMat.append(float(currLine[1]) * 2 - 1)
    return mat(self.min_max(dataMat)), labelMat

  def sigmoid(self, inx):
    # return 1.0/(1+exp(min(max(-inx,-10),10)))
    return 1.0 / (1 + exp(-inx))

  # 得到對應(yīng)的特征weight的矩陣
  def fit(self, data, feature_potential=8, alpha=0.01, iter=100):
    # alpha是學(xué)習(xí)速率
    self.alpha = alpha
    self.feature_potential = feature_potential
    self.iter = iter
    # dataMatrix用的是mat, classLabels是列表
    dataMatrix, classLabels = self.loadDataSet(data)
    print('dataMatrix:',dataMatrix.shape)
    print('classLabels:',classLabels)
    k = self.feature_potential
    m, n = shape(dataMatrix)
    # 初始化參數(shù)
    w = zeros((n, 1)) # 其中n是特征的個(gè)數(shù)
    w_0 = 0.
    v = normalvariate(0, 0.2) * ones((n, k))
    for it in range(self.iter): # 迭代次數(shù)
      # 對每一個(gè)樣本，優(yōu)化
      for x in range(m):
        # 這邊注意一個(gè)數(shù)學(xué)知識：對應(yīng)點(diǎn)積的地方通常會(huì)有sum，對應(yīng)位置積的地方通常都沒有，詳細(xì)參見矩陣運(yùn)算規(guī)則，本處計(jì)算邏輯在：http://blog.csdn.net/google19890102/article/details/45532745
        # xi·vi,xi與vi的矩陣點(diǎn)積
        inter_1 = dataMatrix[x] * v
        # xi與xi的對應(yīng)位置乘積  與  xi^2與vi^2對應(yīng)位置的乘積  的點(diǎn)積
        inter_2 = multiply(dataMatrix[x], dataMatrix[x]) * multiply(v, v) # multiply對應(yīng)元素相乘
        # 完成交叉項(xiàng),xi*vi*xi*vi - xi^2*vi^2
        interaction = sum(multiply(inter_1, inter_1) - inter_2) / 2.
        # 計(jì)算預(yù)測的輸出
        p = w_0 + dataMatrix[x] * w + interaction
        print('classLabels[x]:',classLabels[x])
        print('預(yù)測的輸出p:', p)
        # 計(jì)算sigmoid(y*pred_y)-1
        loss = self.sigmoid(classLabels[x] * p[0, 0]) - 1
        if loss >= -1:
          loss_res = '正方向 '
        else:
          loss_res = '反方向'
        # 更新參數(shù)
        w_0 = w_0 - self.alpha * loss * classLabels[x]
        for i in range(n):
          if dataMatrix[x, i] != 0:
            w[i, 0] = w[i, 0] - self.alpha * loss * classLabels[x] * dataMatrix[x, i]
            for j in range(k):
              v[i, j] = v[i, j] - self.alpha * loss * classLabels[x] * (
                  dataMatrix[x, i] * inter_1[0, j] - v[i, j] * dataMatrix[x, i] * dataMatrix[x, i])
      print('the no %s times, the loss arrach %s' % (it, loss_res))
    self._w_0, self._w, self._v = w_0, w, v

  def predict(self, X):
    if (self._w_0 == None) or (self._w == None).any() or (self._v == None).any():
      raise NotFittedError("Estimator not fitted, call `fit` first")
    # 類型檢查
    if isinstance(X, np.ndarray):
      pass
    else:
      try:
        X = np.array(X)
      except:
        raise TypeError("numpy.ndarray required for X")
    w_0 = self._w_0
    w = self._w
    v = self._v
    m, n = shape(X)
    result = []
    for x in range(m):
      inter_1 = mat(X[x]) * v
      inter_2 = mat(multiply(X[x], X[x])) * multiply(v, v) # multiply對應(yīng)元素相乘
      # 完成交叉項(xiàng)
      interaction = sum(multiply(inter_1, inter_1) - inter_2) / 2.
      p = w_0 + X[x] * w + interaction # 計(jì)算預(yù)測的輸出
      pre = self.sigmoid(p[0, 0])
      result.append(pre)
    return result

  def getAccuracy(self, data):
    dataMatrix, classLabels = self.loadDataSet(data)
    w_0 = self._w_0
    w = self._w
    v = self._v
    m, n = shape(dataMatrix)
    allItem = 0
    error = 0
    result = []
    for x in range(m):
      allItem += 1
      inter_1 = dataMatrix[x] * v
      inter_2 = multiply(dataMatrix[x], dataMatrix[x]) * multiply(v, v) # multiply對應(yīng)元素相乘
      # 完成交叉項(xiàng)
      interaction = sum(multiply(inter_1, inter_1) - inter_2) / 2.
      p = w_0 + dataMatrix[x] * w + interaction # 計(jì)算預(yù)測的輸出
      pre = self.sigmoid(p[0, 0])
      result.append(pre)
      if pre < 0.5 and classLabels[x] == 1.0:
        error += 1
      elif pre >= 0.5 and classLabels[x] == -1.0:
        error += 1
      else:
        continue
    # print(result)
    value = 1 - float(error) / allItem
    return value


class NotFittedError(Exception):
  """
  Exception class to raise if estimator is used before fitting
  """
  pass


if __name__ == '__main__':
  fm()

以上是“Python如何實(shí)現(xiàn)FM算法”這篇文章的所有內(nèi)容，感謝各位的閱讀！希望分享的內(nèi)容對大家有幫助，更多相關(guān)知識，歡迎關(guān)注億速云行業(yè)資訊頻道！