爲何與Yes樂隊最類似的歌手是Coda(小田和奏)？

在網易雲音樂中，點入某歌手主頁時會有一欄爲類似歌手：

Yes是1968年成立的著名前衛搖滾樂隊，雲音樂中與其類似的列表中有Rush，然後者則是來自加拿大的一支異常出色的前衛搖滾樂隊，這對於一個喜好前衛搖滾的聽衆來講，固然是一個很不錯的推薦。然而，Yes的推薦列表中的第一位是Coda(小田和奏)，此君是何許人也，此君並不像Rush、Yes同樣是幾十年前便成立並在搖滾史上留下濃墨重彩的前衛搖滾音樂人，那雲音樂又爲什麼會把此君放在Yes類似列表的第一位呢？

稍微看過Yes歌曲的評論便會知道箇中緣由，在當今中國大地上Yes這支樂隊的歌曲獲得曝光，很大程度上是靠《JOJO的奇妙冒險》這部日本動漫中使用了該樂隊的音樂，而云音樂裏的評論亦均是來自此動漫觀光團，而同時這位Coda(小田和奏)的熱門歌曲便是該動漫的片頭片尾曲，如此便爲他們之間的「類似」找到了緣由，該觀光團必然是既聽了Yes的《Roundabout》，又聽了此君的做品，由於該動漫的火爆，使得雲音樂中那些聽過Yes而不看動漫的前衛聽衆相比由於看該動漫而聽Yes的觀光團反而成了少數，這直接越過了音樂內容本質，改變了類似列表的值。

看到了上面的例子，很顯然這種音樂人之間的類似關聯很難僅靠了解音樂的編輯來給出，一方面是音樂人數量巨大，另外一個是可能產生關聯的維度太多，這便須要尋求集體智慧。

這篇文章使用矩陣分解處理Last.fm 360K數據集，利用海量用戶數據來完成音樂人之間類似關係的創建。

加載數據

Last.fm 360K數據集包含了36萬users對30萬artists的收聽狀況，文件大小1.6G，共17559530行記錄，而在現實的工業應用中，世界著名的音樂軟件Spotify在2014年即擁有2400萬users與2000萬songs。

這種數量級，已經不能像以前介紹item-based CF時直接用numpy生成一個這樣shape的用戶物品關係矩陣了，稍做計算17559530/(300000*360000)=0.01%後發現此數據集很是稀疏，僅0.01%的位置上有值，畢竟一個user不可能聽到不少的artists。scipy提供了幾種處理這種稀疏矩陣的方法，可使矩陣中那些有真實值的位置才佔內存，這裏使用csr_matrix。

import time
import pandas as pd
import numpy as np
from scipy.sparse import csr_matrix, diags
from scipy.sparse.linalg import spsolve

df = pd.read_table("~/Desktop/lastfm-dataset-360K/usersha1-artmbid-artname-plays.tsv", 
                         usecols=[0, 2, 3], 
                         names=['user', 'artist', 'plays'],
                         na_filter=False)

# 作user_id與artist_name到user_item_matrix index的映射.
df['user'] = df['user'].astype('category')
df['artist'] = df['artist'].astype('category')

# 根據上一步生成的id到index的映射，及df的id的關係，生成csr_matrix.
plays = csr_matrix((df['plays'].astype(float), 
                   (df['user'].cat.codes, 
                    df['artist'].cat.codes)))

print ('user count ', plays.shape[0])
print ('artist count ', plays.shape[1])
print ('plays matrix memory usage: %d MB.' % (plays.data.nbytes/1024/1024))

def get_row_index_by_user(user):
    for index, i in enumerate(df['user'].cat.categories):
        if i == user:
            return index
    return None

def get_col_index_by_artist(artist):
    for index, i in enumerate(df['artist'].cat.categories):
        if i == artist:
            return index
    return None
    
def get_sparse_matrix_item(i, j):
    return plays.getrow(i).getcol(j).data[0]
複製代碼

user count  358868
artist count  292365
plays matrix memory usage: 133 MB.
複製代碼

加載數據後獲得稀疏矩陣plays，裏面填的值是某user收聽某artist的次數，同時寫了幾個輔助函數，能夠經過user_id，artist_name獲得在矩陣中相應的值。這裏作一個抽樣檢查，head一下文件第5行會看到'00000c289a1829a808ac09c00daf10bc3c4e223b'這位用戶收聽'red hot chili peppers'次數爲691次。

user1_index = get_row_index_by_user('00000c289a1829a808ac09c00daf10bc3c4e223b')
artist1_index = get_col_index_by_artist('red hot chili peppers')
print ('00000c289a1829a808ac09c00daf10bc3c4e223b listened red hot chili peppers count: ', get_sparse_matrix_item(user1_index, artist1_index))
複製代碼

00000c289a1829a808ac09c00daf10bc3c4e223b listened red hot chili peppers count:  691.0
複製代碼

矩陣分解理論基礎

固然使用基於領域的協同過濾也能夠作這個事情，但已經再也不合適，首先向量維度越高相互之間計算類似度越慢，更重要的是，決定用戶是否喜歡聽某歌手的緣由多是多方面的，但毫不是要到30萬這個數量級這麼多的方面來決定的，這裏面並非一個線性的關係。

矩陣分解是對上述用戶物品關係矩陣的降維打擊，爲user和item各生成一個低維的隱因子，來表達各自的特徵：

將原矩陣分解成兩個矩陣，兩矩陣的乘積儘量與原矩陣相等，越接近越好。對於電影評分這樣的顯示反饋來說，損失函數爲：

而對於聽音樂這樣的隱示反饋，損失函數須要添加一個置信度：

其中即爲置信度，與收聽次數有關：

爲收聽次數，爲超參數，經驗默認值40，爲正則化係數，用來防止過擬合。

對上述損失函數的優化方式爲ALS的變體Weighted-ALS(加權最小交替二乘法)，將item固定，對損失函數求導可得用戶向量爲：

上式中爲這次計算中被固定的矩陣，爲以用戶u的置信度向量生成的對角矩陣，是用戶u的收聽向量，在用戶聽過的artist的位置填入1，其餘位置填0。

矩陣分解代碼

有了上面的公式，很容易使用numpy表示這一計算過程。

def weighted_alternating_least_squares(plays, factors, alpha=40, regularization=0.1, iterations=20):
    Cui = (plays * alpha).astype('double') # 這裏爲了保持矩陣的稀疏，暫不加1，在後面的過程當中補上；
    users, items = Cui.shape

    X = np.random.rand(users, factors) * 0.01
    Y = np.random.rand(items, factors) * 0.01

    Ciu = Cui.T.tocsr()
    for iteration in range(iterations):
        least_squares(Cui, X, Y, regularization)
        least_squares(Ciu, Y, X, regularization)

    return X, Y

def least_squares(Cui, X, Y, regularization):
    users, factors = X.shape
    for u in range(users):
        conf = Cui[u,:].toarray() ；
        pref = conf.copy()
        conf = conf + 1 # 此時僅取出一行，補上1
        pref[pref != 0] = 1
        Cu = diags(conf, [0])
        A = Y.T.dot(Cu).dot(Y) + regularization * np.eye(factors)
        b = Y.T.dot(Cu).dot(pref.T)
        X[u] = spsolve(A, b) 
複製代碼

依然是數據量的問題，上述代碼基本是無法用的。大矩陣各類點乘，空間和時間消費都難以承受。優化方法是利用數學，顯然等於，因爲當用戶u沒有收聽artist i時，爲1，使得很是稀疏，同時p(u)亦很是稀疏，計算時將矩陣點乘拆開，僅取出有值的部分循環相加便可，如此可極大程度提高計算速度並對內存的需求。

implicit的做者在blog中實現了這一過程，裏面全是線性代數：

def nonzeros(m, row):
    for index in range(m.indptr[row], m.indptr[row+1]):
        yield m.indices[index], m.data[index]

def least_squares(Cui, X, Y, regularization):
    users, factors = X.shape
    YtY = Y.T.dot(Y)

    for u in range(users):
        if u % 10000 == 0 and u > 0:
            print (u)
        # accumulate YtCuY + regularization * I in A
        A = YtY + regularization * np.eye(factors)

        # accumulate YtCuPu in b
        b = np.zeros(factors)

        for i, confidence in nonzeros(Cui, u):
            factor = Y[i]
            A += (confidence - 1) * np.outer(factor, factor)
            b += confidence * factor

        # Xu = (YtCuY + regularization * I)^-1 (YtCuPu)
        X[u] = np.linalg.solve(A, b)
複製代碼

user_factors, artist_factors = weighted_alternating_least_squares(plays, 50)
複製代碼

儘管如此，這個分解依然在mac上跑了三個多小時，固然還有更快的方法，implicit庫中，做者加了C++代碼，能夠多線程，而且可使用GPU加速，10幾分鐘就能夠跑下來360K數據集。

對於置信度的計算方法，論文裏還給了另外一個效果更好的公式：，本文提到的第一個置信度計算方式無疑會致使the beatles問題，the beatles太火了你們都聽過，致使會出現一些莫名其秒的類似結果，implicit做者推薦使用bm25算法來計算此weight，能夠消除這種影響，使結果趨於更加合理的方向。

from implicit.nearest_neighbours import bm25_weight
from implicit.als import AlternatingLeastSquares

model = AlternatingLeastSquares(factors=50, regularization=0.01, iterations = 50)
model.fit(bm25_weight(plays.T.tocsr()))

user_factors = model.user_factors
artist_factors = model.item_factors
複製代碼

100%|██████████| 50.0/50 [20:15<00:00, 23.36s/it]
複製代碼

獲取類似artists

生成兩矩陣後，表示用戶u和artist i的低維稠密向量分別爲user_factors[u]與artist_factors[i]，它們維數相同。在顯示反饋中可用來作用戶u對物品i的評分預測，兩向量求點積便可；同時，對於兩個artists的隱因子向量artist_factors[i1]與artist_factors[i2]，依然可使用餘弦類似度公式來計算二者之間的類似度。

問題依然存在，要找到與artist i最類似的幾個artists，須要遍歷30萬隱因子向量計算並排序，這個代價依然是巨大的。既然這裏如今有這個問題，那麼Spotify確定早就遇到了，它曾經的推薦組技術帶頭人Erik Bernhardsson使用C++與Python寫了一個用於解決這種近臨搜索問題的庫annoy，使用起來就是創建一個索引把向量加入進去，搜索時拿着索引去搜很快能獲得結果，Approximate Nearest Neighbors，這個Approximate是指在時間與類似程度準確度上進行了取捨。另外類似的工具還有nmslib與facebook開源的faiss。

from annoy import AnnoyIndex
import random

artist_nn_index = AnnoyIndex(50)
for i in range(artist_factors.shape[0]):
    artist_nn_index.add_item(i, artist_factors[i])

artist_nn_index.build(25)
複製代碼

def get_similar_artists(artist, n = 20):
    similar_artist_list = list()
    for i in artist_nn_index.get_nns_by_item(artist, n):
        similar_artist_list.append(df['artist'].cat.categories[i])
    return similar_artist_list

yes = get_col_index_by_artist('yes')
the_clash = get_col_index_by_artist('the clash')
the_smiths = get_col_index_by_artist('the smiths')
pink_floyd = get_col_index_by_artist('pink floyd')
blur = get_col_index_by_artist('blur')

print ('yes similar artists:\n', get_similar_artists(yes))
print ('----------')
print ('the_clash similar artists:\n', get_similar_artists(the_clash))
print ('----------')
print ('the_smiths similar artists:\n', get_similar_artists(the_smiths))
print ('----------')
print ('pink_floyd similar artists:\n', get_similar_artists(pink_floyd))
print ('----------')
print ('blur similar artists:\n', get_similar_artists(blur))
複製代碼

yes similar artists:
 ['yes', 'emerson, lake & palmer', 'genesis', 'rush', 'jethro tull', 'king crimson', 'gentle giant', 'the moody blues', 'the alan parsons project', 'camel', 'kansas', 'david gilmour', 'focus', 'supertramp', 'jeff beck', 'roger waters', 'peter gabriel', 'steely dan', 'marillion', 'van der graaf generator']
----------
the_clash similar artists:
 ['the clash', 'ramones', 'pixies', 'iggy pop', 'david bowie', 'the pogues', 'the specials', 'the smiths', 'the rolling stones', 'the cure', 'the white stripes', 'lou reed', 'violent femmes', 'the velvet underground', 'johnny cash', 'joy division', 'beastie boys', 'the kinks', 'nirvana', 'misfits']
----------
the_smiths similar artists:
 ['the smiths', 'morrissey', 'the cure', 'joy division', 'david bowie', 'pixies', 'echo & the bunnymen', 'the clash', 'the jesus and mary chain', 'pulp', 'interpol', 'the velvet underground', 'sonic youth', 'arcade fire', 'elliott smith', 'radiohead', 'my bloody valentine', 'blur', 'lou reed', 'nick drake']
----------
pink_floyd similar artists:
 ['pink floyd', 'led zeppelin', 'the doors', 'jimi hendrix', 'queen', 'jethro tull', 'the police', 'jefferson airplane', 'the jimi hendrix experience', 'the who', 'creedence clearwater revival', 'the rolling stones', 'nirvana', 'dire straits', 'deep purple', 'genesis', 'santana', 'david gilmour', 'john lennon', 'pearl jam']
----------
blur similar artists:
 ['blur', 'franz ferdinand', 'supergrass', 'pulp', 'the dandy warhols', 'the verve', 'manic street preachers', 'the white stripes', 'oasis', 'the stone roses', 'beck', 'primal scream', 'arctic monkeys', 'the coral', 'kasabian', 'david bowie', 'eels', 'the smiths', 'kaiser chiefs', 'new order']
複製代碼

獲得的結果使人欣喜。

yes的類似列表裏有king crimson，rush，emerson, lake & palmer，genesis，都是前衛得不行的樂隊，在前衛搖滾全家福上很容易找到他們的身影；
the clash的列表中ramones，pixies，iggy pop都是朋克樂隊，joy division與the clash都是1976年成立，前者是後朋的先驅，看到the velvet underground真的笑出聲，畢竟「每一位朋克、後朋克和先鋒流行藝術家在過去的30年中都欠下了‘地下絲絨’樂隊一筆靈感的債務，哪怕只是受到了間接的影響。」；
the smiths是上個世紀80年代英國獨立搖滾的表明，列表中第一位是主唱莫老師，而後有同時期的joy division，echo & the bunnymen，都是獨立搖滾的表明；
pink floyd的列表也沒必要多說，關鍵字70年代，同一時期歷史評價極高的幾支樂隊都在推薦之列；
最後看到blur的列表裏面有pulp和oasis，就放心了。

參考

《Collaborative Filtering for Implicit Feedback Datasets》
Finding Similar Music using Matrix Factorization
A Gentle Introduction to Recommender Systems with Implicit Feedback