python实现FM算法

1、通常我们在做逻辑回归或者线性回归的时候一般都是没有考虑特征之间相乘产出的情况（特征交叉）
    假设有3个特征,那么就会有3种特征相乘的组合,如果变量的维度比较少，我们是可以先计算这些特征组合在进行求相关的参数，但是如果遇到维度特别大的时候而且特征又是比较稀疏的时候可以考虑用FM算法

2、FM算法的原理其实很简单，只是利用这个简单的公式  = ，

假设线性回归的模型

, 

那么就可以转化为

，

那么损失函数：

的计算类似

逻辑回归的导数计算类似，只是把损失函数改成交叉墒的公式，可以参考这篇文章https://www.jianshu.com/p/40c7358040c9

3、代码：

import pandas as pd
import numpy as np


# y_hat = w0 + w1*x1 + w2*x2 + w3*x3 + w4*x1*x2 + w5*x1*x3 + w6*x2*x3 =>
# y_hat = w0 + w1*x1 + w2*x2 + w3*x3 + 1/2[(v1*x1 + v2*x2 + v3*x3)^2 -(v1*x1)^2 -(v2*x2)^2 - (v3*x3)^2]
# optimize params:w0,w1,w2,w3,v1,v2,v3

class FM:
    def __init__(self,target,lr=0.1,reg_lambda=0.0001,max_iter=300,bias=True):
        self.__lr = lr
        self.__reg_lambda = reg_lambda
        self.__max_iter = max_iter
        self.__bais = bias
        self.__coef = []
        self.__v_coef = []
        self.__target = target # 线性回归 if reg else 逻辑回归

    @property
    def coef(self):
        return self.__coef

    @property
    def v_coef(self):
        return self.__v_coef

    def args_init(self,columns):
        if self.__bais:
            self._intercept = 0.001
        for _ in columns:
            self.__coef.append(np.random.normal(0,0.1))
            self.__v_coef.append(np.random.normal(0,0.1))

    def sigmoid(self,x):
        ss = np.exp(-np.clip(x,a_min=-1e32,a_max=1e32))
        return 1/(1 + ss)

    @staticmethod
    def calc_ks(y_true, y_pred):
        y_true = y_true.reshape(-1, )
        y_pred = y_pred.reshape(-1, )
        sort_index = np.argsort(y_pred, kind="mergesort")[::-1]
        y_pred = y_pred[sort_index]
        y_true = y_true[sort_index]

        # 获取不同的值
        diff = np.diff(y_pred)
        distinct_value_indices = np.where(np.diff(y_pred))[0]
        threshold_idxs = np.r_[distinct_value_indices, y_pred.size - 1]
        tps = np.cumsum(y_true)[threshold_idxs]
        fps = np.cumsum((1 - y_true))[threshold_idxs]

        threshold = y_pred[threshold_idxs]

        tps = np.r_[0, tps]
        fps = np.r_[0, fps]

        tpr = tps / (tps[-1] + 1e-32)
        fpr = fps / (fps[-1] + 1e-32)
        return max(np.abs(tpr - fpr))

    def predict(self,X):
        bias = np.zeros([X.shape[0],1]) + self._intercept
        W = np.zeros_like(X) + self.__coef
        V = np.zeros_like(X) + self.__v_coef
        W_X = np.sum(W * X,axis=1).reshape(-1,1)
        V_x = (np.square(np.sum(V * X,axis=1)) - np.sum(np.square(V*X),axis=1)).reshape(-1,1) / 2
        if self.__target.startswith('reg'):
            return bias + W_X + V_x
        return self.sigmoid(bias + W_X + V_x)

    def v_opt(self,X,i):
        # i 表示第几列
        V = np.zeros_like(X) + self.__v_coef
        V_opt = np.sum(V * X,axis=1).reshape(-1,1) * X[:,[i]] - np.sum(V[:,[i]] * np.square(X[:,[i]]),axis=1).reshape(-1,1) # 计算关于vi的导数
        return V_opt

    def fit(self,X,y):
        if isinstance(X,np.ndarray):
            cols = ['x%s'%(i+1) for i in range(X.shape[1])]
        elif isinstance(X,pd.core.frame.DataFrame):
            cols = X.columns
        else:
            raise TypeError("input X must be array or dataframe")
        if isinstance(y,pd.core.series.Series):
            y = np.array(y).reshape(-1,1)
        if y.shape == 1:
            y = y.reshape(-1,1)
        self.args_init(cols)
        for _ in range(self.__max_iter):
            r = self.predict(X)
            ks = self.calc_ks(y,r)
            v_opt = {}
            for s in range(X.shape[1]):
                v_opt[s] = self.v_opt(X,s)
            print("ks = ",ks)
            for j in range(len(X)):
                if self.__bais:
                    self._intercept -= float(self.__lr * (r[j] -y[j]))+ self.__reg_lambda * self._intercept
                    for i in range(len(cols)):
                        self.__coef[i] -= float(self.__lr * (r[j] - y[j])*(X[j,[i]])) + \
                                          self.__reg_lambda*self.__coef[i] # 正则
                    for i in range(len(cols)):
                        self.__v_coef[i] -= float(self.__lr * (r[j] - y[j]) * v_opt[i][j]) + \
                                            self.__reg_lambda * self.__v_coef[i]


if __name__ == '__main__':
    from sklearn.datasets import make_moons

    X,y = make_moons()
    fm = FM(target='classify')
    print(type(X))
    fm.fit(X,y)

    from sklearn.linear_model import LogisticRegression

    lr = LogisticRegression()
    xx = X[:,[0]] * X[:,[1]]
    X = np.concatenate([X,xx],axis=1)
    lr.fit(X, y)
    pred = lr.predict_proba(X)[:, [1]]
    print(fm.calc_ks(y, pred))