基於RBM的推薦算法

很多人講RBM都要從能量函數講起，由能量最低導出極小化目標函數（你聽說過最常見的建立目標函數的方法可能是最小化平方誤差或者最大化似然函數），然后用梯度下降法求解，得到網絡參數。Introduction to Restricted Boltzmann Machines這篇博客沒有遵循這種套路來講RBM，它直接給RBM網絡權重的訓練方法，講得淺顯易懂，清新脫俗。本文只是對英文版的翻譯。

在基於LFM(Latent Factor Model)的推薦算法一文中我們介紹了用因子分析法來做推薦。比如用戶購買了《推薦系統實踐》、《用Python做數據分析》，背后的隱藏因子是數據挖掘；用戶購買了《一課經濟學》、《郎咸平說》，背后的隱藏因子是經濟學。因子分析法就是找到用戶對各個隱藏因子的喜好程度$U=(uf_1,uf_2,...,uf_n)$，以及商品在各個隱藏因子上的概率分布$I=(if_1,if_2,...,if_n)$，然后兩個向量做乘法即得到用戶對商品的喜好程度。RBM可以理解為一種二值化的因子分析法（當然RBM還有其他的理解方式）。先來一張圖看看RBM的網絡結構。

1.  它只有兩層：可視層和隱藏層，兩層之間是全連接。另有一個偏置單元，跟所有的可視單元和隱藏單元都有連接。可視層之間、隱藏層之間無連接。
2. 所有連接都是雙向的，並且是帶權重的。
3. 所有神經元的狀態只有0和1兩種。
4. RBM是一個隨機網絡，即所有神經元以一個概率值選擇狀態為0還是1。偏置單元是個例外，它總為1。偏置單元用來反應商品固有的受歡迎程度，所謂“固有”就是跟外界無關，反應到RBM網絡里面就是隱藏單元的狀態並不是完全由可視層決定的，也由隱藏層神經元自身固有的一些因素決定，這些固有的因素就由偏置單元來承載。反過來對於可視層也一樣，可視層神經元固有的受歡迎程度由偏置單元來承載。我們在帶偏置的LFM中也說明了偏置的作用。

RBM的運作方式

可視層的神經元用$x_i$表示，隱藏層神經元用$x_j$表示，它們之間的權重用$w_{ij}$表示，可視層神經元個數為$m$，隱藏層神經元個數為$n$。當給定可視層狀態后，用下式更新隱藏層的狀態。

\begin{equation}net_j=\sum_{i=0}^m{x_i{w_{ij}}}\end{equation}$x_0$是偏置單元，總為1

\begin{equation}prob(j)=sigmoid(net_j)=\frac{1}{1+e^{-net_j}}\end{equation}

$x_j$以概率$prob(j)$取1，以概率$1-prob(j)$取0。

$sigmoid$函數關於(0,0.5)這一點中心對稱，$x$為正時$sigmoid(x)>0.5$，$x\to\infty$時$sigmoid(x)\to{1}$。

根據隱藏層求可視層方式雷同，就不寫公式了。

對於推薦系統來說，我們知道用戶購買了哪些商品，將對應的可視層神經元置為1，其他置為0，求出隱藏層狀態，由隱藏層再返回來求可視層狀態，這個時候可視層哪些神經元為1我們就把相應有商品推薦給用戶。

權重學習方法

訓練RBM網絡就是訓練權重$w_{ij}$。首先隨機初始化$w_{ij}$，然后每一次拿一個樣本（即可視層是已知的）經歷下面的步驟。

1. 由可視層的$x_i$算出隱藏層的$x_j$，令$w_{ij}$的正向梯度為\begin{equation}positive(w_{ij})=x_i*{x_j}\end{equation}
2. 由隱藏層$x_j$再來反向計算$x’_i$，注意此時算出的$x’_i$跟原先的$x_i$已經不一樣了，令$w_{ij}$的負向梯度為\begin{equation}negative(w_{ij})=x'_i*{x_j}\end{equation}
3. 更新權重\begin{equation}w_{ij}=w_{ij}+\alpha*(positive(w_{ij})-negative(w_{ij}))\end{equation}

我們不去深究為什么正向梯度和負向梯度是這樣一個公式。上述學習方式叫對比散度(contrastive divergence)法。

循環拿樣本去訓練網絡，不停迭代，直到收斂（即$x’_i$和$x_i$很接近）。

實踐中的優化

1. 上面講述中我們拿$x_i$去RBM網絡中返回一次得到$x’_i$，然后就開始計算$negative(w_{ij})$，改進方法是多往返幾次后再計算$negative(w_{ij})$。
2. 計算$positive(w_{ij})$時用$prob(i)*prob(j)$，而非$x_i*{x_j}$。$negative(w_{ij})$同樣。
3. 加正則項，對較大的權重$w_{ij}$進行懲罰。
4. 更新$w_{ij}$時加動量向，即本次前進的方向是本次的梯度與上次迭代中梯度的線性加權。
5. 每次調整權重時使用一批樣本，而非不一個樣本。雖然計算結果是一樣的，但由於numpy對矩陣乘法做了加速優化，比逐個計算向量乘法要快。

rbm.py

# coding=utf-8
__author__ = "orisun"

import numpy as np


class RBM(object):

    def __init__(self, num_visible, num_hidden, learn_rate=0.1, learn_batch=1000):
        self.num_visible = num_visible  # 可視層神經元個數
        self.num_hidden = num_hidden    # 隱藏層神經元個數
        self.learn_rate = learn_rate    # 學習率
        self.learn_batch = learn_batch  # 每次根據多少樣本進行學習

        '''初始化連接權重'''
        self.weights = 0.1 * \
            np.random.randn(self.num_visible,
                            self.num_hidden)  # 依據0.1倍的標准正太分布隨機生成權重
        # 第一行插入全0，即偏置和隱藏層的權重初始化為0
        self.weights = np.insert(self.weights, 0, 0, axis=0)
        # 第一列插入全0，即偏置和可視層的權重初始化為0
        self.weights = np.insert(self.weights, 0, 0, axis=1)

    def _logistic(self, x):
        '''直接使用1.0 / (1.0 + np.exp(-x))容易發警告“RuntimeWarning: overflowencountered in exp”，
           轉換成如下等價形式后算法會更穩定
        '''
        return 0.5 * (1 + np.tanh(0.5 * x))

    def train(self, rating_data, max_steps=1000, eps=1.0e-4):
        '''迭代訓練，得到連接權重
        '''
        for step in xrange(max_steps):  # 迭代訓練多少次
            error = 0.0  # 誤差平方和
            # 每次拿一批樣本還調整權重
            for i in xrange(0, rating_data.shape[0], self.learn_batch):
                num_examples = min(self.learn_batch, rating_data.shape[0] - i)
                data = rating_data[i:i + num_examples, :]
                data = np.insert(data, 0, 1, axis=1)  # 第一列插入全1，即偏置的值初始化為1

                pos_hidden_activations = np.dot(data, self.weights)
                pos_hidden_probs = self._logistic(pos_hidden_activations)
                pos_hidden_states = pos_hidden_probs > np.random.rand(
                    num_examples, self.num_hidden + 1)
                # pos_associations=np.dot(data.T,pos_hidden_states)         #對隱藏層作二值化
                pos_associations = np.dot(
                    data.T, pos_hidden_probs)  # 對隱藏層不作二值化

                neg_visible_activations = np.dot(
                    pos_hidden_states, self.weights.T)
                neg_visible_probs = self._logistic(neg_visible_activations)
                neg_visible_probs[:, 0] = 1  # 強行把偏置的值重置為1
                neg_hidden_activations = np.dot(
                    neg_visible_probs, self.weights)
                neg_hidden_probs = self._logistic(neg_hidden_activations)
                # neg_hidden_states=neg_hidden_probs>np.random.rand(num_examples,self.num_hidden+1)
                # neg_associations=np.dot(neg_visible_probs.T,neg_hidden_states)      #對隱藏層作二值化
                neg_associations = np.dot(
                    neg_visible_probs.T, neg_hidden_probs)  # 對隱藏層不作二值化

                # 更新權重。另外一種嘗試是帶沖量的梯度下降，即本次前進的方向是本次梯度與上一次梯度的線性加權和（這樣的話需要額外保存上一次的梯度）
                self.weights += self.learn_rate * \
                    (pos_associations - neg_associations) / num_examples

                # 計算誤差平方和
                error += np.sum((data - neg_visible_probs)**2)
            if error < eps:  # 所有樣本的誤差平方和低於閾值於終止迭代
                break
            print 'iteration %d, error is %f' % (step, error)

    def getHidden(self, visible_data):
        '''根據輸入層得到隱藏層
           visible_data是一個matrix，每行代表一個樣本
        '''
        num_examples = visible_data.shape[0]
        hidden_states = np.ones((num_examples, self.num_hidden + 1))
        visible_data = np.insert(visible_data, 0, 1, axis=1)  # 第一列插入偏置
        hidden_activations = np.dot(visible_data, self.weights)
        hidden_probs = self._logistic(hidden_activations)
        hidden_states[:, :] = hidden_probs > np.random.rand(
            num_examples, self.num_hidden + 1)
        hidden_states = hidden_states[:, 1:]            # 即首列刪掉，即把偏置去掉
        return hidden_states

    def getVisible(self, hidden_data):
        '''根據隱藏層得到輸入層
           hidden_data是一個matrix，每行代表一個樣本
        '''
        num_examples = hidden_data.shape[0]
        visible_states = np.ones((num_examples, self.num_visible + 1))
        hidden_data = np.insert(hidden_data, 0, 1, axis=1)
        visible_activations = np.dot(hidden_data, self.weights.T)
        visible_probs = self._logistic(visible_activations)
        visible_states[:, :] = visible_probs > np.random.rand(
            num_examples, self.num_visible + 1)
        visible_states = visible_states[:, 1:]
        return visible_states

    def predict(self, visible_data):
        num_examples = visible_data.shape[0]
        hidden_states = np.ones((num_examples, self.num_hidden + 1))
        visible_data = np.insert(visible_data, 0, 1, axis=1)  # 第一列插入偏置
        '''forward'''
        hidden_activations = np.dot(visible_data, self.weights)
        hidden_probs = self._logistic(hidden_activations)
        # hidden_states[:, :] = hidden_probs > np.random.rand(
        #     num_examples, self.num_hidden + 1)
        '''backward'''
        visible_states = np.ones((num_examples, self.num_visible + 1))
        # visible_activations = np.dot(hidden_states, self.weights.T)  #對隱藏層作二值化
        visible_activations = np.dot(hidden_probs, self.weights.T)  # 對隱藏層不作二值化
        visible_probs = self._logistic(visible_activations)  # 直接返回可視層的概率值

        return visible_probs[:, 1:]  # 把第0列(偏置)去掉

if __name__ == '__main__':
    rbm = RBM(num_visible=6, num_hidden=2, learn_rate=0.1, learn_batch=1000)
    rating_data = np.array([[1, 1, 1, 0, 0, 0], [1, 0, 1, 0, 0, 0], [1, 1, 1, 0, 0, 0], [
                           0, 0, 1, 1, 1, 0], [0, 0, 1, 1, 0, 0], [0, 0, 1, 1, 1, 0]])
    rbm.train(rating_data, max_steps=500, eps=1.0e-4)
    print 'weight:\n', rbm.weights
    rating = np.array([[0, 0, 0, 0.9, 0.7, 0]])  # 評分需要做歸一化。該用戶喜歡第四、五項
    hidden_data = rbm.getHidden(rating)
    print 'hidden_data:\n', hidden_data
    visible_data = rbm.getVisible(hidden_data)
    print 'visible_data:\n', visible_data
    predict_data = rbm.predict(rating)
    print '推薦得分:'
    for i, score in enumerate(predict_data[0, :]):
        print i, score  # 第三、四、五項的推薦得分很高，同時用戶已明確表示過喜歡四、五，所以我們把第三項推薦給用戶