经典算法笔记：异常检测和推荐系统
本文是吳恩達老師的機器學習課程[1]的筆記和代碼復現部分（聚類、降維）。
作者：黃海廣[2]
備注：筆記和作業（含數據、原始作業文件）、視頻都在github[3]中下載。
我將陸續將課程筆記和課程代碼發布在公眾號“機器學習初學者”，敬請關注。這個是第五部分：異常檢測和推薦系統，是原教程的第9周，包含了筆記和作業代碼（原課程作業是?OCTAVE的，這里是復現的?python?代碼）
已完成部分：
第一部分:回歸
第二部分:邏輯回歸
第三部分:支持向量機
第四部分：無監督學習
本文作業代碼[4]可以下載完整版
筆記的markdown 文件[5]
筆記的pdf 文件[6]
筆記部分-目錄
十五、異常檢測(Anomaly Detection)
15.1 問題的動機
參考文檔: 15 - 1 - Problem Motivation (8 min).mkv
在接下來的一系列視頻中，我將向大家介紹異常檢測(Anomaly detection)問題。這是機器學習算法的一個常見應用。這種算法的一個有趣之處在于：它雖然主要用于非監督學習問題，但從某些角度看，它又類似于一些監督學習問題。
什么是異常檢測呢？為了解釋這個概念，讓我舉一個例子吧：
假想你是一個飛機引擎制造商，當你生產的飛機引擎從生產線上流出時，你需要進行QA(質量控制測試)，而作為這個測試的一部分，你測量了飛機引擎的一些特征變量，比如引擎運轉時產生的熱量，或者引擎的振動等等。
這樣一來，你就有了一個數據集，從到，如果你生產了個引擎的話，你將這些數據繪制成圖表，看起來就是這個樣子：
這里的每個點、每個叉，都是你的無標簽數據。這樣，異常檢測問題可以定義如下：我們假設后來有一天，你有一個新的飛機引擎從生產線上流出，而你的新飛機引擎有特征變量。所謂的異常檢測問題就是：我們希望知道這個新的飛機引擎是否有某種異常，或者說，我們希望判斷這個引擎是否需要進一步測試。因為，如果它看起來像一個正常的引擎，那么我們可以直接將它運送到客戶那里，而不需要進一步的測試。
給定數據集?，我們假使數據集是正常的，我們希望知道新的數據??是不是異常的，即這個測試數據不屬于該組數據的幾率如何。我們所構建的模型應該能根據該測試數據的位置告訴我們其屬于一組數據的可能性?。
上圖中，在藍色圈內的數據屬于該組數據的可能性較高，而越是偏遠的數據，其屬于該組數據的可能性就越低。
這種方法稱為密度估計，表達如下：
欺詐檢測：
用戶的第個活動特征
模型?為我們其屬于一組數據的可能性，通過檢測非正常用戶。
異常檢測主要用來識別欺騙。例如在線采集而來的有關用戶的數據，一個特征向量中可能會包含如：用戶多久登錄一次，訪問過的頁面，在論壇發布的帖子數量，甚至是打字速度等。嘗試根據這些特征構建一個模型，可以用這個模型來識別那些不符合該模式的用戶。
再一個例子是檢測一個數據中心，特征可能包含：內存使用情況，被訪問的磁盤數量，CPU的負載，網絡的通信量等。根據這些特征可以構建一個模型，用來判斷某些計算機是不是有可能出錯了。
15.2 高斯分布
參考視頻: 15 - 2 - Gaussian Distribution (10 min).mkv
在這個視頻中，我將介紹高斯分布，也稱為正態分布。回顧高斯分布的基本知識。
通常如果我們認為變量??符合高斯分布?則其概率密度函數為：??
我們可以利用已有的數據來預測總體中的μ和σ的計算方法如下：?
高斯分布樣例：
注：機器學習中對于方差我們通常只除以而非統計學中的。這里順便提一下，在實際使用中，到底是選擇使用還是其實區別很小，只要你有一個還算大的訓練集，在機器學習領域大部分人更習慣使用這個版本的公式。這兩個版本的公式在理論特性和數學特性上稍有不同，但是在實際使用中，他們的區別甚小，幾乎可以忽略不計。
15.3 算法
參考視頻: 15 - 3 - Algorithm (12 min).mkv
在本節視頻中，我將應用高斯分布開發異常檢測算法。
異常檢測算法：
對于給定的數據集?，我們要針對每一個特征計算??和??的估計值。
一旦我們獲得了平均值和方差的估計值，給定新的一個訓練實例，根據模型計算?：
當時，為異常。
下圖是一個由兩個特征的訓練集，以及特征的分布情況：
下面的三維圖表表示的是密度估計函數，軸為根據兩個特征的值所估計值：
我們選擇一個，將作為我們的判定邊界，當時預測數據為正常數據，否則為異常。
在這段視頻中，我們介紹了如何擬合，也就是?的概率值，以開發出一種異常檢測算法。同時，在這節課中，我們也給出了通過給出的數據集擬合參數，進行參數估計，得到參數??和?，然后檢測新的樣本，確定新樣本是否是異常。
在接下來的課程中，我們將深入研究這一算法，同時更深入地介紹，怎樣讓算法工作地更加有效。
15.4 開發和評價一個異常檢測系統
參考視頻: 15 - 4 - Developing and Evaluating an Anomaly Detection System (13 min). mkv
異常檢測算法是一個非監督學習算法，意味著我們無法根據結果變量??的值來告訴我們數據是否真的是異常的。我們需要另一種方法來幫助檢驗算法是否有效。當我們開發一個異常檢測系統時，我們從帶標記（異常或正常）的數據著手，我們從其中選擇一部分正常數據用于構建訓練集，然后用剩下的正常數據和異常數據混合的數據構成交叉檢驗集和測試集。
例如：我們有 10000 臺正常引擎的數據，有 20 臺異常引擎的數據。 我們這樣分配數據：
6000 臺正常引擎的數據作為訓練集
2000 臺正常引擎和 10 臺異常引擎的數據作為交叉檢驗集
2000 臺正常引擎和 10 臺異常引擎的數據作為測試集
具體的評價方法如下：
根據測試集數據，我們估計特征的平均值和方差并構建函數
對交叉檢驗集，我們嘗試使用不同的值作為閥值，并預測數據是否異常，根據值或者查準率與查全率的比例來選擇?
選出??后，針對測試集進行預測，計算異常檢驗系統的值，或者查準率與查全率之比
15.5 異常檢測與監督學習對比
參考視頻: 15 - 5 - Anomaly Detection vs. Supervised Learning (8 min).mkv
之前我們構建的異常檢測系統也使用了帶標記的數據，與監督學習有些相似，下面的對比有助于選擇采用監督學習還是異常檢測：
兩者比較：
| 非常少量的正向類（異常數據?）, 大量的負向類（） | 同時有大量的正向類和負向類 |
| 許多不同種類的異常，非常難。根據非常 少量的正向類數據來訓練算法。 | 有足夠多的正向類實例，足夠用于訓練 算法，未來遇到的正向類實例可能與訓練集中的非常近似。 |
| 未來遇到的異常可能與已掌握的異常、非常的不同。 | |
| 例如： 欺詐行為檢測 生產（例如飛機引擎）檢測數據中心的計算機運行狀況 | 例如：郵件過濾器 天氣預報 腫瘤分類 |
希望這節課能讓你明白一個學習問題的什么樣的特征，能讓你把這個問題當做是一個異常檢測，或者是一個監督學習的問題。另外，對于很多技術公司可能會遇到的一些問題，通常來說，正樣本的數量很少，甚至有時候是 0，也就是說，出現了太多沒見過的不同的異常類型，那么對于這些問題，通常應該使用的算法就是異常檢測算法。
15.6 選擇特征
參考視頻: 15 - 6 - Choosing What Features to Use (12 min).mkv
對于異常檢測算法，我們使用的特征是至關重要的，下面談談如何選擇特征：
異常檢測假設特征符合高斯分布，如果數據的分布不是高斯分布，異常檢測算法也能夠工作，但是最好還是將數據轉換成高斯分布，例如使用對數函數：，其中??為非負常數； 或者?，為 0-1 之間的一個分數，等方法。(編者注：在python中，通常用np.log1p()函數，就是?，可以避免出現負數結果，反向函數就是np.expm1())
誤差分析：
一個常見的問題是一些異常的數據可能也會有較高的值，因而被算法認為是正常的。這種情況下誤差分析能夠幫助我們，我們可以分析那些被算法錯誤預測為正常的數據，觀察能否找出一些問題。我們可能能從問題中發現我們需要增加一些新的特征，增加這些新特征后獲得的新算法能夠幫助我們更好地進行異常檢測。
異常檢測誤差分析：
我們通常可以通過將一些相關的特征進行組合，來獲得一些新的更好的特征（異常數據的該特征值異常地大或小），例如，在檢測數據中心的計算機狀況的例子中，我們可以用CPU負載與網絡通信量的比例作為一個新的特征，如果該值異常地大，便有可能意味著該服務器是陷入了一些問題中。
在這段視頻中，我們介紹了如何選擇特征，以及對特征進行一些小小的轉換，讓數據更像正態分布，然后再把數據輸入異常檢測算法。同時也介紹了建立特征時，進行的誤差分析方法，來捕捉各種異常的可能。希望你通過這些方法，能夠了解如何選擇好的特征變量，從而幫助你的異常檢測算法，捕捉到各種不同的異常情況。
15.7 多元高斯分布（選修）
參考視頻: 15 - 7 - Multivariate Gaussian Distribution (Optional) (14 min).mkv
假使我們有兩個相關的特征，而且這兩個特征的值域范圍比較寬，這種情況下，一般的高斯分布模型可能不能很好地識別異常數據。其原因在于，一般的高斯分布模型嘗試的是去同時抓住兩個特征的偏差，因此創造出一個比較大的判定邊界。
下圖中是兩個相關特征，洋紅色的線（根據 ε 的不同其范圍可大可小）是一般的高斯分布模型獲得的判定邊界，很明顯綠色的X所代表的數據點很可能是異常值，但是其值卻仍然在正常范圍內。多元高斯分布將創建像圖中藍色曲線所示的判定邊界。
在一般的高斯分布模型中，我們計算??的方法是： 通過分別計算每個特征對應的幾率然后將其累乘起來，在多元高斯分布模型中，我們將構建特征的協方差矩陣，用所有的特征一起來計算?。
我們首先計算所有特征的平均值，然后再計算協方差矩陣：?
注:其中?是一個向量，其每一個單元都是原特征矩陣中一行數據的均值。最后我們計算多元高斯分布的:??其中：
是定矩陣，在 Octave 中用 det(sigma)計算
?是逆矩陣，下面我們來看看協方差矩陣是如何影響模型的：
上圖是 5 個不同的模型，從左往右依次分析：
是一個一般的高斯分布模型
通過協方差矩陣，令特征 1 擁有較小的偏差，同時保持特征 2 的偏差
通過協方差矩陣，令特征 2 擁有較大的偏差，同時保持特征 1 的偏差
通過協方差矩陣，在不改變兩個特征的原有偏差的基礎上，增加兩者之間的正相關性
通過協方差矩陣，在不改變兩個特征的原有偏差的基礎上，增加兩者之間的負相關性
多元高斯分布模型與原高斯分布模型的關系：
可以證明的是，原本的高斯分布模型是多元高斯分布模型的一個子集，即像上圖中的第 1、2、3，3 個例子所示，如果協方差矩陣只在對角線的單位上有非零的值時，即為原本的高斯分布模型了。
原高斯分布模型和多元高斯分布模型的比較：
| 不能捕捉特征之間的相關性 但可以通過將特征進行組合的方法來解決 | 自動捕捉特征之間的相關性 |
| 計算代價低，能適應大規模的特征 | 計算代價較高 訓練集較小時也同樣適用 |
| 必須要有?，不然的話協方差矩陣不可逆的，通常需要??另外特征冗余也會導致協方差矩陣不可逆 |
原高斯分布模型被廣泛使用著，如果特征之間在某種程度上存在相互關聯的情況，我們可以通過構造新新特征的方法來捕捉這些相關性。
如果訓練集不是太大，并且沒有太多的特征，我們可以使用多元高斯分布模型。
15.8 使用多元高斯分布進行異常檢測（可選）
參考視頻: 15 - 8 - Anomaly Detection using the Multivariate Gaussian Distribution (Optional) (14 min).mkv
在我們談到的最后一個視頻，關于多元高斯分布，看到的一些建立的各種分布模型，當你改變參數，?和?。在這段視頻中，讓我們用這些想法，并應用它們制定一個不同的異常檢測算法。
要回顧一下多元高斯分布和多元正態分布：
分布有兩個參數，??和?。其中這一個維向量和??的協方差矩陣，是一種的矩陣。而這里的公式的概率，如按??和參數化?，和你的變量??和?，你可以得到一個范圍的不同分布一樣，你知道的，這些都是三個樣本，那些我們在以前的視頻看過了。
因此，讓我們談談參數擬合或參數估計問題：
我有一組樣本是一個維向量，我想我的樣本來自一個多元高斯分布。我如何嘗試估計我的參數??和??以及標準公式？
估計他們是你設置??是你的訓練樣本的平均值。
?并設置：??這其實只是當我們使用PCA算法時候，有??時寫出來。所以你只需插入上述兩個公式，這會給你你估計的參數??和你估計的參數?。所以，這里給出的數據集是你如何估計??和?。讓我們以這種方法而只需將其插入到異常檢測算法。那么，我們如何把所有這一切共同開發一個異常檢測算法？
首先，我們把我們的訓練集，和我們的擬合模型，我們計算，要知道，設定和描述的一樣。
如圖，該分布在中央最多，越到外面的圈的范圍越小。
并在該點是出路這里的概率非常低。
原始模型與多元高斯模型的關系如圖：
其中：協方差矩陣為：
原始模型和多元高斯分布比較如圖：
十六、推薦系統(Recommender Systems)
16.1 問題形式化
參考視頻: 16 - 1 - Problem Formulation (8 min).mkv
在接下來的視頻中，我想講一下推薦系統。我想講推薦系統有兩個原因：
第一、僅僅因為它是機器學習中的一個重要的應用。在過去幾年，我偶爾訪問硅谷不同的技術公司，我常和工作在這兒致力于機器學習應用的人們聊天，我常問他們，最重要的機器學習的應用是什么，或者，你最想改進的機器學習應用有哪些。我最常聽到的答案是推薦系統。現在，在硅谷有很多團體試圖建立很好的推薦系統。因此，如果你考慮網站像亞馬遜，或網飛公司或易趣，或iTunes Genius，有很多的網站或系統試圖推薦新產品給用戶。如，亞馬遜推薦新書給你，網飛公司試圖推薦新電影給你，等等。這些推薦系統，根據瀏覽你過去買過什么書，或過去評價過什么電影來判斷。這些系統會帶來很大一部分收入，比如為亞馬遜和像網飛這樣的公司。因此，對推薦系統性能的改善，將對這些企業的有實質性和直接的影響。
推薦系統是個有趣的問題，在學術機器學習中因此，我們可以去參加一個學術機器學習會議，推薦系統問題實際上受到很少的關注，或者，至少在學術界它占了很小的份額。但是，如果你看正在發生的事情，許多有能力構建這些系統的科技企業，他們似乎在很多企業中占據很高的優先級。這是我為什么在這節課討論它的原因之一。
我想討論推薦系統地第二個原因是：這個班視頻的最后幾集我想討論機器學習中的一些大思想，并和大家分享。這節課我們也看到了，對機器學習來說，特征是很重要的，你所選擇的特征，將對你學習算法的性能有很大的影響。因此，在機器學習中有一種大思想，它針對一些問題，可能并不是所有的問題，而是一些問題，有算法可以為你自動學習一套好的特征。因此，不要試圖手動設計，而手寫代碼這是目前為止我們常干的。有一些設置，你可以有一個算法，僅僅學習其使用的特征，推薦系統就是類型設置的一個例子。還有很多其它的，但是通過推薦系統，我們將領略一小部分特征學習的思想，至少，你將能夠了解到這方面的一個例子，我認為，機器學習中的大思想也是這樣。因此，讓我們開始討論推薦系統問題形式化。
我們從一個例子開始定義推薦系統的問題。
假使我們是一個電影供應商，我們有 5 部電影和 4 個用戶，我們要求用戶為電影打分。
前三部電影是愛情片，后兩部則是動作片，我們可以看出Alice和Bob似乎更傾向與愛情片， 而 Carol 和 Dave 似乎更傾向與動作片。并且沒有一個用戶給所有的電影都打過分。我們希望構建一個算法來預測他們每個人可能會給他們沒看過的電影打多少分，并以此作為推薦的依據。
下面引入一些標記：
?代表用戶的數量
?代表電影的數量
?如果用戶 j 給電影??評過分則?
?代表用戶??給電影的評分
代表用戶??評過分的電影的總數
16.2 基于內容的推薦系統
參考視頻: 16 - 2 - Content Based Recommendations (15 min).mkv
在一個基于內容的推薦系統算法中，我們假設對于我們希望推薦的東西有一些數據，這些數據是有關這些東西的特征。
在我們的例子中，我們可以假設每部電影都有兩個特征，如代表電影的浪漫程度，?代表電影的動作程度。
則每部電影都有一個特征向量，如是第一部電影的特征向量為[0.9 0]。
下面我們要基于這些特征來構建一個推薦系統算法。 假設我們采用線性回歸模型，我們可以針對每一個用戶都訓練一個線性回歸模型，如是第一個用戶的模型的參數。 于是，我們有：
用戶??的參數向量
電影??的特征向量
對于用戶??和電影?，我們預測評分為：
代價函數
針對用戶?，該線性回歸模型的代價為預測誤差的平方和，加上正則化項：
其中?表示我們只計算那些用戶??評過分的電影。在一般的線性回歸模型中，誤差項和正則項應該都是乘以，在這里我們將去掉。并且我們不對方差項進行正則化處理。
上面的代價函數只是針對一個用戶的，為了學習所有用戶，我們將所有用戶的代價函數求和：
如果我們要用梯度下降法來求解最優解，我們計算代價函數的偏導數后得到梯度下降的更新公式為：
16.3 協同過濾
參考視頻: 16 - 3 - Collaborative Filtering (10 min).mkv
在之前的基于內容的推薦系統中，對于每一部電影，我們都掌握了可用的特征，使用這些特征訓練出了每一個用戶的參數。相反地，如果我們擁有用戶的參數，我們可以學習得出電影的特征。
但是如果我們既沒有用戶的參數，也沒有電影的特征，這兩種方法都不可行了。協同過濾算法可以同時學習這兩者。
我們的優化目標便改為同時針對和進行。
對代價函數求偏導數的結果如下：
注：在協同過濾從算法中，我們通常不使用方差項，如果需要的話，算法會自動學得。 協同過濾算法使用步驟如下：
初始?為一些隨機小值
使用梯度下降算法最小化代價函數
在訓練完算法后，我們預測為用戶??給電影??的評分
通過這個學習過程獲得的特征矩陣包含了有關電影的重要數據，這些數據不總是人能讀懂的，但是我們可以用這些數據作為給用戶推薦電影的依據。
例如，如果一位用戶正在觀看電影?，我們可以尋找另一部電影，依據兩部電影的特征向量之間的距離的大小。
16.4 協同過濾算法
參考視頻: 16 - 4 - Collaborative Filtering Algorithm (9 min).mkv
協同過濾優化目標：
給定，估計：
給定，估計：
同時最小化和：
16.5 向量化：低秩矩陣分解
參考視頻: 16 - 5 - Vectorization_ Low Rank Matrix Factorization (8 min).mkv
在上幾節視頻中，我們談到了協同過濾算法，本節視頻中我將會講到有關該算法的向量化實現，以及說說有關該算法你可以做的其他事情。
舉例子：
當給出一件產品時，你能否找到與之相關的其它產品。
一位用戶最近看上一件產品，有沒有其它相關的產品，你可以推薦給他。
我將要做的是：實現一種選擇的方法，寫出協同過濾算法的預測情況。
我們有關于五部電影的數據集，我將要做的是，將這些用戶的電影評分，進行分組并存到一個矩陣中。
我們有五部電影，以及四位用戶，那么 這個矩陣??就是一個 5 行 4 列的矩陣，它將這些電影的用戶評分數據都存在矩陣里：
| Love at last | 5 | 5 | 0 | 0 |
| Romance forever | 5 | ? | ? | 0 |
| Cute puppies of love | ? | 4 | 0 | ? |
| Nonstop car chases | 0 | 0 | 5 | 4 |
| Swords vs. karate | 0 | 0 | 5 | ? |
推出評分：
找到相關影片：
現在既然你已經對特征參數向量進行了學習，那么我們就會有一個很方便的方法來度量兩部電影之間的相似性。例如說：電影??有一個特征向量，你是否能找到一部不同的電影?，保證兩部電影的特征向量之間的距離和很小，那就能很有力地表明電影和電影??在某種程度上有相似，至少在某種意義上，某些人喜歡電影?，或許更有可能也對電影??感興趣。總結一下，當用戶在看某部電影??的時候，如果你想找 5 部與電影非常相似的電影，為了能給用戶推薦 5 部新電影，你需要做的是找出電影?，在這些不同的電影中與我們要找的電影??的距離最小，這樣你就能給你的用戶推薦幾部不同的電影了。
通過這個方法，希望你能知道，如何進行一個向量化的計算來對所有的用戶和所有的電影進行評分計算。同時希望你也能掌握，通過學習特征參數，來找到相關電影和產品的方法。
16.6 推行工作上的細節：均值歸一化
參考視頻: 16 - 6 - Implementational Detail_ Mean Normalization (9 min).mkv
讓我們來看下面的用戶評分數據：
如果我們新增一個用戶 Eve，并且 Eve 沒有為任何電影評分，那么我們以什么為依據為Eve推薦電影呢？
我們首先需要對結果?矩陣進行均值歸一化處理，將每一個用戶對某一部電影的評分減去所有用戶對該電影評分的平均值：
然后我們利用這個新的??矩陣來訓練算法。 如果我們要用新訓練出的算法來預測評分，則需要將平均值重新加回去，預測，對于Eve，我們的新模型會認為她給每部電影的評分都是該電影的平均分。
代碼部分- 異常檢測和推薦系統
在本練習中，我們將使用高斯模型實現異常檢測算法，并將其應用于檢測網絡上的故障服務器。?我們還將看到如何使用協作過濾構建推薦系統，并將其應用于電影推薦數據集。
Anomaly detection（異常檢測）
我們的第一個任務是使用高斯模型來檢測數據集中未標記的示例是否應被視為異常。?我們有一個簡單的二維數據集開始，以幫助可視化該算法正在做什么。
import numpy as np import pandas as pd import matplotlib.pyplot as plt import seaborn as sb from scipy.io import loadmat data = loadmat('data/ex8data1.mat') X = data['X'] X.shape (307, 2) fig, ax = plt.subplots(figsize=(12,8)) ax.scatter(X[:,0], X[:,1]) plt.show()這是一個非常緊密的聚類，幾個值遠離了聚類。?在這個簡單的例子中，這些可以被認為是異常的。?為了弄清楚，我們正在為數據中的每個特征估計高斯分布。?為此，我們將創建一個返回每個要素的均值和方差的函數。
def estimate_gaussian(X):mu = X.mean(axis=0)sigma = X.var(axis=0)return mu, sigma mu, sigma = estimate_gaussian(X) mu, sigma (array([14.11222578, 14.99771051]), array([1.83263141, 1.70974533]))現在我們有了我們的模型參數，我們需要確定概率閾值，這表明一個樣本應該被認為是一個異常。?為此，我們需要使用一組標記的驗證數據（其中真實異常樣本已被標記），并在給出不同閾值的情況下，對模型的性能進行鑒定。
Xval = data['Xval'] yval = data['yval']Xval.shape, yval.shape ((307, 2), (307, 1))我們還需要一種計算數據點屬于正態分布的概率的方法。?幸運的是 SciPy 有這個內置的方法。
from scipy import stats dist = stats.norm(mu[0], sigma[0]) dist.pdf(15) 0.1935875044615038我們還可以將數組傳遞給概率密度函數，并獲得數據集中每個點的概率密度。
dist.pdf(X[:,0])[0:50] array([0.183842 , 0.20221694, 0.21746136, 0.19778763, 0.20858956,0.21652359, 0.16991291, 0.15123542, 0.1163989 , 0.1594734 ,0.21716057, 0.21760472, 0.20141857, 0.20157497, 0.21711385,0.21758775, 0.21695576, 0.2138258 , 0.21057069, 0.1173018 ,0.20765108, 0.21717452, 0.19510663, 0.21702152, 0.17429399,0.15413455, 0.21000109, 0.20223586, 0.21031898, 0.21313426,0.16158946, 0.2170794 , 0.17825767, 0.17414633, 0.1264951 ,0.19723662, 0.14538809, 0.21766361, 0.21191386, 0.21729442,0.21238912, 0.18799417, 0.21259798, 0.21752767, 0.20616968,0.21520366, 0.1280081 , 0.21768113, 0.21539967, 0.16913173])我們計算并保存給定上述的高斯模型參數的數據集中每個值的概率密度。
p = np.zeros((X.shape[0], X.shape[1])) p[:,0] = stats.norm(mu[0], sigma[0]).pdf(X[:,0]) p[:,1] = stats.norm(mu[1], sigma[1]).pdf(X[:,1])p.shape (307, 2)我們還需要為驗證集（使用相同的模型參數）執行此操作。?我們將使用與真實標簽組合的這些概率來確定將數據點分配為異常的最佳概率閾值。
pval = np.zeros((Xval.shape[0], Xval.shape[1])) pval[:,0] = stats.norm(mu[0], sigma[0]).pdf(Xval[:,0]) pval[:,1] = stats.norm(mu[1], sigma[1]).pdf(Xval[:,1]) pval.shape (307, 2)接下來，我們需要一個函數，找到給定概率密度值和真實標簽的最佳閾值。?為了做到這一點，我們將為不同的 epsilon 值計算 F1 分數。?F1 是真陽性，假陽性和假陰性的數量的函數。?方程式在練習文本中。
def select_threshold(pval, yval):best_epsilon = 0best_f1 = 0f1 = 0step = (pval.max() - pval.min()) / 1000for epsilon in np.arange(pval.min(), pval.max(), step):preds = pval < epsilontp = np.sum(np.logical_and(preds == 1, yval == 1)).astype(float)fp = np.sum(np.logical_and(preds == 1, yval == 0)).astype(float)fn = np.sum(np.logical_and(preds == 0, yval == 1)).astype(float)precision = tp / (tp + fp)recall = tp / (tp + fn)f1 = (2 * precision * recall) / (precision + recall)if f1 > best_f1:best_f1 = f1best_epsilon = epsilonreturn best_epsilon, best_f1 epsilon, f1 = select_threshold(pval, yval) epsilon, f1 (0.009566706005956842, 0.7142857142857143)最后，我們可以將閾值應用于數據集，并可視化結果。
# indexes of the values considered to be outliers outliers = np.where(p < epsilon) outliers (array([300, 301, 301, 303, 303, 304, 306, 306], dtype=int64),array([1, 0, 1, 0, 1, 0, 0, 1], dtype=int64)) fig, ax = plt.subplots(figsize=(12,8)) ax.scatter(X[:,0], X[:,1]) ax.scatter(X[outliers[0],0], X[outliers[0],1], s=50, color='r', marker='o') plt.show()紅點是被標記為異常值的點。?這些看起來很合理。?有一些分離（但沒有被標記）的右上角也可能是一個異常值，但是相當接近。
協同過濾
推薦引擎使用基于項目和用戶的相似性度量來檢查用戶的歷史偏好，以便為用戶可能感興趣的新“事物”提供建議。在本練習中，我們將實現一種稱為協作過濾的特定推薦系統算法，并將其應用于 電影評分的數據集。
我們首先加載并檢查我們將要使用的數據。
data = loadmat('data/ex8_movies.mat')Y 是包含從 1 到 5 的等級的（數量的電影 x 數量的用戶）數組.R 是包含指示用戶是否給電影評分的二進制值的“指示符”數組。?兩者應該具有相同的維度。
Y = data['Y'] R = data['R'] Y.shape, R.shape ((1682, 943), (1682, 943))我們可以通過平均排序 Y 來評估電影的平均評級。
Y[1,np.where(R[1,:]==1)[0]].mean() 3.2061068702290076我們還可以通過將矩陣渲染成圖像來嘗試“可視化”數據。?我們不能從這里收集太多，但它確實給我們了解用戶和電影的相對密度。
fig, ax = plt.subplots(figsize=(12,12)) ax.imshow(Y) ax.set_xlabel('Users') ax.set_ylabel('Movies') fig.tight_layout() plt.show()接下來，我們將實施協同過濾的代價函數。?直覺上，“代價”是指一組電影評級預測偏離真實預測的程度。?代價方程在練習文本中給出。?它基于文本中稱為 X 和 Theta 的兩組參數矩陣。?這些“展開”到“參數”輸入中，以便稍后可以使用 SciPy 的優化包。?請注意，我已經在注釋中包含數組/矩陣形狀（對于我們在本練習中使用的數據），以幫助說明矩陣交互如何工作。
代價函數
def cost(params, Y, R, num_features):Y = np.matrix(Y) # (1682, 943)R = np.matrix(R) # (1682, 943)num_movies = Y.shape[0]num_users = Y.shape[1]# reshape the parameter array into parameter matricesX = np.matrix(np.reshape(params[:num_movies * num_features], (num_movies, num_features))) # (1682, 10)Theta = np.matrix(np.reshape(params[num_movies * num_features:], (num_users, num_features))) # (943, 10)# initializationsJ = 0# compute the costerror = np.multiply((X * Theta.T) - Y, R) # (1682, 943)squared_error = np.power(error, 2) # (1682, 943)J = (1. / 2) * np.sum(squared_error)return J為了測試這一點，我們提供了一組我們可以評估的預訓練參數。?為了保持評估時間的少點，我們將只看一小段數據。
params_data = loadmat('data/ex8_movieParams.mat') X = params_data['X'] Theta = params_data['Theta'] X.shape, Theta.shape ((1682, 10), (943, 10)) users = 4 movies = 5 features = 3X_sub = X[:movies, :features] Theta_sub = Theta[:users, :features] Y_sub = Y[:movies, :users] R_sub = R[:movies, :users]params = np.concatenate((np.ravel(X_sub), np.ravel(Theta_sub)))cost(params, Y_sub, R_sub, features) 22.224603725685675接下來我們需要實現梯度計算。?就像我們在練習 4 中使用神經網絡實現一樣，我們將擴展代價函數來計算梯度。
def cost(params, Y, R, num_features):Y = np.matrix(Y) # (1682, 943)R = np.matrix(R) # (1682, 943)num_movies = Y.shape[0]num_users = Y.shape[1]# reshape the parameter array into parameter matricesX = np.matrix(np.reshape(params[:num_movies * num_features],(num_movies, num_features))) # (1682, 10)Theta = np.matrix(np.reshape(params[num_movies * num_features:],(num_users, num_features))) # (943, 10)# initializationsJ = 0X_grad = np.zeros(X.shape) # (1682, 10)Theta_grad = np.zeros(Theta.shape) # (943, 10)# compute the costerror = np.multiply((X * Theta.T) - Y, R) # (1682, 943)squared_error = np.power(error, 2) # (1682, 943)J = (1. / 2) * np.sum(squared_error)# calculate the gradientsX_grad = error * ThetaTheta_grad = error.T * X# unravel the gradient matrices into a single arraygrad = np.concatenate((np.ravel(X_grad), np.ravel(Theta_grad)))return J, grad J, grad = cost(params, Y_sub, R_sub, features) J, grad (22.224603725685675,array([ -2.52899165, 7.57570308, -1.89979026, -0.56819597,3.35265031, -0.52339845, -0.83240713, 4.91163297,-0.76677878, -0.38358278, 2.26333698, -0.35334048,-0.80378006, 4.74271842, -0.74040871, -10.5680202 ,4.62776019, -7.16004443, -3.05099006, 1.16441367,-3.47410789, 0. , 0. , 0. ,0. , 0. , 0. ]))我們的下一步是在代價和梯度計算中添加正則化。?我們將創建一個最終的正則化版本的功能（請注意，此版本包含一個額外的“學習率”參數，在文本中稱為“lambda”）。
def cost(params, Y, R, num_features, learning_rate):Y = np.matrix(Y) # (1682, 943)R = np.matrix(R) # (1682, 943)num_movies = Y.shape[0]num_users = Y.shape[1]# reshape the parameter array into parameter matricesX = np.matrix(np.reshape(params[:num_movies * num_features], (num_movies, num_features))) # (1682, 10)Theta = np.matrix(np.reshape(params[num_movies * num_features:], (num_users, num_features))) # (943, 10)# initializationsJ = 0X_grad = np.zeros(X.shape) # (1682, 10)Theta_grad = np.zeros(Theta.shape) # (943, 10)# compute the costerror = np.multiply((X * Theta.T) - Y, R) # (1682, 943)squared_error = np.power(error, 2) # (1682, 943)J = (1. / 2) * np.sum(squared_error)# add the cost regularizationJ = J + ((learning_rate / 2) * np.sum(np.power(Theta, 2)))J = J + ((learning_rate / 2) * np.sum(np.power(X, 2)))# calculate the gradients with regularizationX_grad = (error * Theta) + (learning_rate * X)Theta_grad = (error.T * X) + (learning_rate * Theta)# unravel the gradient matrices into a single arraygrad = np.concatenate((np.ravel(X_grad), np.ravel(Theta_grad)))return J, grad J, grad = cost(params, Y_sub, R_sub, features, 1.5) J, grad (31.34405624427422,array([ -0.95596339, 6.97535514, -0.10861109, 0.60308088,2.77421145, 0.25839822, 0.12985616, 4.0898522 ,-0.89247334, 0.29684395, 1.06300933, 0.66738144,0.60252677, 4.90185327, -0.19747928, -10.13985478,2.10136256, -6.76563628, -2.29347024, 0.48244098,-2.99791422, -0.64787484, -0.71820673, 1.27006666,1.09289758, -0.40784086, 0.49026541]))這個結果再次與練習代碼的預期輸出相匹配，所以看起來正則化是正常的。?在我們訓練模型之前，我們有一個最后步驟， 我們的任務是創建自己的電影評分，以便我們可以使用該模型來生成個性化的推薦。?為我們提供一個連接電影索引到其標題的文件。?接著我們將文件加載到字典中。
movie_idx = {} f = open('data/movie_ids.txt',encoding= 'gbk') for line in f:tokens = line.split(' ')tokens[-1] = tokens[-1][:-1]movie_idx[int(tokens[0]) - 1] = ' '.join(tokens[1:]) movie_idx[0] 'Toy Story (1995)'我們將使用練習中提供的評分。
ratings = np.zeros((1682, 1))ratings[0] = 4 ratings[6] = 3 ratings[11] = 5 ratings[53] = 4 ratings[63] = 5 ratings[65] = 3 ratings[68] = 5 ratings[97] = 2 ratings[182] = 4 ratings[225] = 5 ratings[354] = 5print('Rated {0} with {1} stars.'.format(movie_idx[0], str(int(ratings[0])))) print('Rated {0} with {1} stars.'.format(movie_idx[6], str(int(ratings[6])))) print('Rated {0} with {1} stars.'.format(movie_idx[11], str(int(ratings[11])))) print('Rated {0} with {1} stars.'.format(movie_idx[53], str(int(ratings[53])))) print('Rated {0} with {1} stars.'.format(movie_idx[63], str(int(ratings[63])))) print('Rated {0} with {1} stars.'.format(movie_idx[65], str(int(ratings[65])))) print('Rated {0} with {1} stars.'.format(movie_idx[68], str(int(ratings[68])))) print('Rated {0} with {1} stars.'.format(movie_idx[97], str(int(ratings[97])))) print('Rated {0} with {1} stars.'.format(movie_idx[182], str(int(ratings[182])))) print('Rated {0} with {1} stars.'.format(movie_idx[225], str(int(ratings[225])))) print('Rated {0} with {1} stars.'.format(movie_idx[354], str(int(ratings[354])))) Rated Toy Story (1995) with 4 stars. Rated Twelve Monkeys (1995) with 3 stars. Rated Usual Suspects, The (1995) with 5 stars. Rated Outbreak (1995) with 4 stars. Rated Shawshank Redemption, The (1994) with 5 stars. Rated While You Were Sleeping (1995) with 3 stars. Rated Forrest Gump (1994) with 5 stars. Rated Silence of the Lambs, The (1991) with 2 stars. Rated Alien (1979) with 4 stars. Rated Die Hard 2 (1990) with 5 stars. Rated Sphere (1998) with 5 stars.我們可以將自己的評級向量添加到現有數據集中以包含在模型中。
R = data['R'] Y = data['Y'] Y = np.append(Y, ratings, axis=1) R = np.append(R, ratings != 0, axis=1) Y.shape, R.shape, ratings.shape ((1682, 944), (1682, 944), (1682, 1))我們不只是準備訓練協同過濾模型。?我們只需要定義一些變量并對評級進行規一化。
movies = Y.shape[0] # 1682 users = Y.shape[1] # 944 features = 10 learning_rate = 10. X = np.random.random(size=(movies, features)) Theta = np.random.random(size=(users, features)) params = np.concatenate((np.ravel(X), np.ravel(Theta)))X.shape, Theta.shape, params.shape ((1682, 10), (944, 10), (26260,)) Ymean = np.zeros((movies, 1)) Ynorm = np.zeros((movies, users))for i in range(movies):idx = np.where(R[i,:] == 1)[0]Ymean[i] = Y[i,idx].mean()Ynorm[i,idx] = Y[i,idx] - Ymean[i]Ynorm.mean() 5.507036456515984e-19 from scipy.optimize import minimizefmin = minimize(fun=cost, x0=params, args=(Ynorm, R, features, learning_rate),method='CG', jac=True, options={'maxiter': 100}) X = np.matrix(np.reshape(fmin.x[:movies * features], (movies, features))) Theta = np.matrix(np.reshape(fmin.x[movies * features:], (users, features))) X.shape, Theta.shape ((1682, 10), (944, 10))我們訓練好的參數是 X 和 Theta。?我們可以使用這些來為我們添加的用戶創建一些建議。
predictions = X * Theta.T my_preds = predictions[:, -1] + Ymean my_preds.shape (1682, 1) sorted_preds = np.sort(my_preds, axis=0)[::-1] sorted_preds[:10] matrix([[5.00000277],[5.00000236],[5.00000172],[5.00000167],[5.00000018],[4.99999996],[4.99999993],[4.99999963],[4.99999936],[4.99999933]])這給了我們一個排名最高的評級，但我們失去了這些評級的索引。?我們實際上需要使用 argsort 函數來預測評分對應的電影。
idx = np.argsort(my_preds, axis=0)[::-1] idx matrix([[1499],[ 813],[1535],...,[ 313],[ 829],[1431]], dtype=int64) print("Top 10 movie predictions:") for i in range(10):j = int(idx[i])print('Predicted rating of {0} for movie {1}.'.format(str(float(my_preds[j])), movie_idx[j])) Top 10 movie predictions: Predicted rating of 5.0000027697456755 for movie Santa with Muscles (1996). Predicted rating of 5.000002356341384 for movie Great Day in Harlem, A (1994). Predicted rating of 5.000001717626692 for movie Aiqing wansui (1994). Predicted rating of 5.000001669423294 for movie Star Kid (1997). Predicted rating of 5.000000180977937 for movie Prefontaine (1997). Predicted rating of 4.999999960309181 for movie They Made Me a Criminal (1939). Predicted rating of 4.999999928299842 for movie Marlene Dietrich: Shadow and Light (1996) . Predicted rating of 4.999999629264883 for movie Someone Else's America (1995). Predicted rating of 4.9999993570025705 for movie Saint of Fort Washington, The (1993). Predicted rating of 4.999999325074885 for movie Entertaining Angels: The Dorothy Day Story (1996).推薦的電影實際上并不符合練習文本中的內容， 我沒有找到原因在哪里。
參考資料
[1]?機器學習課程:?https://www.coursera.org/course/ml
[2]?黃海廣:?https://github.com/fengdu78
[3]?github:?https://github.com/fengdu78/Coursera-ML-AndrewNg-Notes
[4]?作業代碼:?https://github.com/fengdu78/Coursera-ML-AndrewNg-Notes/blob/master/code/ex8-anomaly%20detection%20and%20recommendation/ML-Exercise8.ipynb
[5]?markdown 文件:?https://github.com/fengdu78/Coursera-ML-AndrewNg-Notes/blob/master/markdown/week9.md
[6]?pdf 文件:?https://github.com/fengdu78/Coursera-ML-AndrewNg-Notes/blob/master/機器學習個人筆記完整版v5.4-A4打印版.pdf
關于本站
“機器學習初學者”公眾號由是黃海廣博士創建，黃博個人知乎粉絲23000+，github排名全球前100名（33000+）。本公眾號致力于人工智能方向的科普性文章，為初學者提供學習路線和基礎資料。原創作品有：吳恩達機器學習個人筆記、吳恩達深度學習筆記等。
往期精彩回顧
那些年做的學術公益-你不是一個人在戰斗
適合初學者入門人工智能的路線及資料下載
吳恩達機器學習課程筆記及資源（github標星12000+，提供百度云鏡像）
吳恩達深度學習筆記及視頻等資源（github標星8500+，提供百度云鏡像）
《統計學習方法》的python代碼實現（github標星7200+）
機器學習的數學精華（在線閱讀版）
備注：加入本站微信群或者qq群，請回復“加群”
加入知識星球（4300+用戶，ID：92416895），請回復“知識星球”
總結
以上是生活随笔為你收集整理的经典算法笔记：异常检测和推荐系统的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： 复现经典：《统计学习方法》第 3 章 k
- 下一篇： 复现经典：《统计学习方法》第1章 统计学