LOF離群因子檢測演算法及python3實現

LOF離群因子檢測演算法及python3實現

來自專欄數據科學與python小記50 人贊了文章

1 相關背景

1.1 異常檢測演算法

隨著數據挖掘技術的快速發展，人們在關注數據整體趨勢的同時，開始越來越關注那些明顯偏離數據整體趨勢的離群數據點，因為這些數據點往往蘊含著更加重要的信息，而處理這些離群數據要依賴於相應的數據挖掘技術。 離群點挖掘的目的是有效的識別出數據集中的異常數據，並且挖掘出數據集中有意義的潛在信息。出現離群點的原因各有不同，其中主要有以下幾種情況：

• 數據來源於異類：如欺詐、入侵、疾病爆發、不同尋常的實驗結果等。這類離群點的產生機制偏離正常軌道，往往是被關注的重點
• 數據變數固有的變化：自然、周期發生等反映數據集的分布特徵，如氣候的突然變化、顧客新型購買模式、基因突變等
• 數據測最或收集誤差：主要是系統誤差、人為的數據讀取偏差、測量設備出現故障或「噪音」干擾
• 隨機或人為誤差：主要原因是實驗平台存在隨機機制，同時人為在數據錄入等過程中可能出現的誤差

離群點檢測具有非常強的實際意義，在相應的應用領域有著廣泛前景。其中工程應用領域主要有以下幾個方面：

• 欺詐檢測：信用卡的不正當行為，如信用卡、社會保障的欺詐行為或者是銀行卡、電話卡的欺詐使用等
• 工業檢測：如計算機網路的非法訪問等
• 活動監控：通過實時檢測手機活躍度或者是股權市場的可疑交易，從而實現檢測移動手機詐騙行為等
• 網路性能：計算機網路性能檢測（穩健性分析），檢測網路堵塞情況等
• 自然生態應用領域：生態系統失調、異常自然氣候的發現等
• 公共服務領域：公共衛生中的異常疾病的爆發、公共安全中的突發事件的發生等

目前，隨著離群點檢測技術的日漸成熟，在未來的發展中將會應用在更多的行業當中，並且能更好地為人類的決策提供指導作用。

離群點檢測的一個目標是從看似雜亂無章的大量數據中挖掘有價值的信息，使這些數據更好地為我們的日常生活所服務。但是現實生活中的數據往往具有成百上千的維度，並且數據量極大，這無疑給目前現有的離群點檢測方法帶來大難題。傳統的離群點檢測方法雖然在各自特定的應用領域裡表現出很好效果，但在高維大數據集中卻不再適用。因此如何把離群點檢測方法有效地應用於大數據、高維度數據，是目前離群點檢測方法的首要目標之一。

1.2 為什麼是異常點檢測

分類學習演算法在訓練時有一個共同的基本假設：不同類別的訓練樣例數目相當。如果不同類別的訓練樣例數目稍有差別，通常影響不大，但若差別很大，則會對學習過程造成困擾。

——周志華《機器學習》

2018 年 Mathorcup 數學建模競賽 B 題——「品牌手機目標用戶的精準營銷」中就出現這樣的問題。原題中，檢測用戶數 1萬+，購買手機用戶 500+。若使用分類演算法，那麼分類器很可能會返回這樣的一個演算法：「所有用戶都不會購買這款手機」，分類的正確率高達96%，但顯然沒有實際意義，因為它並不能預測出任何的正例。

這類問題稱為「類別不平衡」，指不同類別的訓練樣例數目差別很大的情況（上例中，購買的與沒有購買用戶數量差別大）。處理這類問題往往採用「欠採樣」、「過採樣」進行數據處理，但通過這樣的方法，可能會損失原始數據中的信息。因此，從離群點的角度出發，將購買行為視為「異常」，進行離群點挖掘。

1.3 離群點挖掘方法

1.4 LOF 演算法背景

基於密度的離群點檢測方法的關鍵步驟在於給每個數據點都分配一個離散度，其主要思想是：針對給定的數據集，對其中的任意一個數據點，如果在其局部鄰域內的點都很密集，那麼認為此數據點為正常數據點，而離群點則是距離正常數據點最近鄰的點都比較遠的數據點。通常有閾值進行界定距離的遠近。在基於密度的離群點檢測方法中，最具有代表性的方法是局部離群因子檢測方法 (Local Outlier Factor, LOF)。

2 演算法簡介

在眾多的離群點檢測方法中，LOF 方法是一種典型的基於密度的高精度離群點檢測方法。在 LOF 方法中，通過給每個數據點都分配一個依賴於鄰域密度的離群因子 LOF，進而判斷該數據點是否為離群點。若 LOF 1， 則該數據點為離群點；若 LOF 接近於 1，則該數據點為正常數據點。

2.1 距離度量尺度

設對於沒有相同點的樣本集合 ，假設共有 個檢測樣本，數據維數為 ，對於

針對數據集 中的任意兩個數據點 ，定義如下幾種常用距離度量方式

2.1.1 Eucild（歐幾里得）距離：

2.1.2 Hamming（漢明）距離：

注 漢明距離使用在數據傳輸差錯控制編碼裡面，用於度量信息不相同的位數。

取 ，易見 與 中有 2 位數字不相同，因此 與 的漢明距離為 2。

對於數據處理，一種技巧是先對連續數據進行分組，化為分類變數（分組變數），對分類變數可以引入漢明距離進行度量。——沃茲基 · 碩德

2.1.3 Mahalanobis（馬氏）距離：

設樣本集 的協方差矩陣為 ，記其逆矩陣為 若 可逆，對 做 SVD 分解（奇異值分解），得到：

• 若 不可逆，則使用廣義逆矩陣 代替 ，對其求彭羅斯廣義逆，有：

則兩個數據點 的馬氏距離為：

注 馬氏距離表示數據的協方差距離，利用 Cholesky 變換處理不同維度之間的相關性和度量尺度變換的問題，是一種有效計算樣本集之間的相似度的方法。

2.1.4 球面距離

球面距離其實是在歐式距離基礎上進行轉換得到的，並不是一種獨特的距離度量方式，在地理信息轉換中經常使用，本文對此進行詳細介紹。

設 兩點的球面坐標為 。若該球體為地球，則 分別代表緯度和經度。（註：下文的 為 的餘角，便於推導所使用的記號）

如圖所示，連接OA、OB、AB，在 和 中計算：

由於 與 是異面直線， 是它們的公垂線，所成角度經度差為 ，利用異面直線上兩點距離公式：

在 中，由余弦定理：

由於此處的 代表緯度的補角，對其進行轉換：

因此，點 的球面距離為：

此外還有 Chebyshev(切比雪夫)距離、Minkowski(閔科夫斯基)距離、絕對值距離、Lance & Williams 距離，具體問題具體分析，選擇合適的度量方式。

統一使用 表示點 和點 之間的距離。根據定義，易知交換律成立：

2.2 第 距離

定義 為點 的第 距離， ，滿足如下條件

• 在集合中至少存在 個點 ，使得
• 在集合中至多存在 個點 ，使得

簡言之：點 是距離 最近的第 個點

2.3 距離鄰域

定義 設 為點 的第 距離鄰域，滿足：

注 此處的鄰域概念與國內高數教材略有不同（具體的點，而非區間）。該集合中包含所有到點 距離小於點 第 鄰域距離的點。易知有 ，如上圖，點 的第 5 距離鄰域為：

2.4 可達距離

定義 點 到點 的第 可達距離定義為：

注 即點 到點 的第 可達距離至少是點 的第 距離。距離 點最近的 個點，它們到 的可達距離被認為是相當的，且都等於

2.5 局部可達密度

定義 局部可達密度定義為：

注 表示點 的第 鄰域內所有點到 的平均可達距離，位於第 鄰域邊界上的點即使個數大於1，也仍將該範圍內點的個數計為 。如果 和周圍鄰域點是同一簇，那麼可達距離越可能為較小的 ，導致可達距離之和越小，局部可達密度越大。如果 和周圍鄰域點較遠，那麼可達距離可能會取較大值 ，導致可達距離之和越大，局部可達密度越小。

部分資料這裡使用 而不是 。筆者查閱大量資料及數據測試後認為，此處應為 ，否則 會因為過多點在內部同一圓環上（如式13中的 位於同一圓環上）而導致 是一個很小的數，提示此處的密度低，可能為離群值。

此外，本文 2.1 開頭指出「沒有樣本點重合」在這裡也能得到解釋：如果考慮重合樣本點，可能會造成此處的可達密度為 或下文的 為 形式，計算上帶來困擾。

2.6 局部離群因子

注 表示點 的鄰域 內其他點的局部可達密度與點 的局部可達密度之比的平均數。如果這個比值越接近1，說明 的鄰域點密度差不多， 可能和鄰域同屬一簇；如果這個比值小於1，說明 的密度高於其鄰域點密度， 為密集點；如果這個比值大於1，說明 的密度小於其鄰域點密度， 可能是異常點。

3 LOF離群因子檢測演算法python3 實現（第 4 部分進行改進）

此部分先介紹 sklearn 提供的函數，再介紹逐步編程實現方法。

3.1 導入 sklearn 模塊實現

在 python3 中，sklearn 模塊提供了 LOF 離群檢測演算法

3.1.1 導入模塊

import pandas as pdimport matplotlib.pyplot as plt

3.1.2 核心函數

clf = LocalOutlierFactor(n_neighbors=20, algorithm=auto, contamination=0.1, n_jobs=-1)

• n_neighbors = 20：即上文提及的 ，檢測的鄰域點個數超過樣本數則使用所有的樣本進行檢測
• algorithm = auto：使用的求解演算法，使用默認值即可
• contamination = 0.1：範圍為 (0, 0.5)，表示樣本中的異常點比例，默認為 0.1
• n_jobs = -1：並行任務數，設置為-1表示使用所有CPU進行工作
• p = 2：距離度量函數，默認使用歐式距離。（其他距離模型使用較少，這裡不作介紹。具體參考官方文檔）

clf.fit(data)

• 無監督學習，只需要傳入訓練數據data，傳入的數據維度至少是 2 維

clf.kneighbors(data)

• 獲取第 距離鄰域內的每一個點到中心點的距離，並按從小到大排序
• 返回 數組，[距離，樣本索引]

-clf._decision_function(data)

• 獲取每一個樣本點的 LOF 值，該函數範圍 LOF 值的相反數，需要取反號
• clf._decision_function 的輸出方式更為靈活：若使用 clf._predict(data) 函數，則按照原先設置的 contamination 輸出判斷結果（按比例給出判斷結果，異常點返回-1，非異常點返回1）

3.1.3 封裝函數

• localoutlierfactor(data, predict, k)
• 輸入：訓練樣本，測試樣本， 值
• 輸出：每一個測試樣本 LOF 值及對應的第 距離
• plot_lof(result,method)
• 輸入：LOF 值、閾值
• 輸出：以索引為橫坐標的 LOF 值分布圖
• lof(data, predict=None, k=5, method=1, plot = False)
• 輸入：訓練樣本，測試樣本， 值，離群閾值，是否繪製 LOF 圖
• 輸出：離群點、正常點分類情況
• 沒有輸入測試樣本時，默認測試樣本=訓練樣本

def localoutlierfactor(data, predict, k): from sklearn.neighbors import LocalOutlierFactor clf = LocalOutlierFactor(n_neighbors=k + 1, algorithm=auto, contamination=0.1, n_jobs=-1) clf.fit(data) # 記錄 k 鄰域距離 predict[k distances] = clf.kneighbors(predict)[0].max(axis=1) # 記錄 LOF 離群因子，做相反數處理 predict[local outlier factor] = -clf._decision_function(predict.iloc[:, :-1]) return predictdef plot_lof(result, method): import matplotlib.pyplot as plt plt.rcParams[font.sans-serif] = [SimHei] # 用來正常顯示中文標籤 plt.rcParams[axes.unicode_minus] = False # 用來正常顯示負號 plt.figure(figsize=(8, 4)).add_subplot(111) plt.scatter(result[result[local outlier factor] > method].index, result[result[local outlier factor] > method][local outlier factor], c=red, s=50, marker=., alpha=None, label=離群點) plt.scatter(result[result[local outlier factor] <= method].index, result[result[local outlier factor] <= method][local outlier factor], c=black, s=50, marker=., alpha=None, label=正常點) plt.hlines(method, -2, 2 + max(result.index), linestyles=--) plt.xlim(-2, 2 + max(result.index)) plt.title(LOF局部離群點檢測, fontsize=13) plt.ylabel(局部離群因子, fontsize=15) plt.legend() plt.show()def lof(data, predict=None, k=5, method=1, plot=False): import pandas as pd # 判斷是否傳入測試數據，若沒有傳入則測試數據賦值為訓練數據 try: if predict == None: predict = data.copy() except Exception: pass predict = pd.DataFrame(predict) # 計算 LOF 離群因子 predict = localoutlierfactor(data, predict, k) if plot == True: plot_lof(predict, method) # 根據閾值劃分離群點與正常點 outliers = predict[predict[local outlier factor] > method].sort_values(by=local outlier factor) inliers = predict[predict[local outlier factor] <= method].sort_values(by=local outlier factor) return outliers, inliers

3.1.4 實例測試

測試數據：2017年全國大學生數學建模競賽B題數據

測試數據有倆份文件，進行三次測試：沒有輸入測試樣本、輸入測試樣本、測試樣本與訓練樣本互換

測試 1 沒有輸入測試樣本情況：任務密度

數據背景：眾包任務價格制定中，地區任務的密度反映任務的密集程度，從而影響任務的定價，此處不考慮球面距離偏差（即認為是同一個平面上的點），現在需要使用一個合理的指標刻畫任務的密集程度。

import numpy as npimport pandas as pd# 根據文件位置自行修改posi = pd.read_excel(r已結束項目任務數據.xls)lon = np.array(posi["任務gps經度"][:]) # 經度lat = np.array(posi["任務gps 緯度"][:]) # 緯度A = list(zip(lat, lon)) # 按照緯度-經度匹配# 獲取任務密度，取第5鄰域，閾值為2（LOF大於2認為是離群值）outliers1, inliers1 = lof(A, k=5, method = 2)

給定數據中共有835條數據，設置 LOF 閾值為 2，輸出17個離群點信息：

繪製數據點檢測情況分布如下圖所示，其中藍色表示任務分布情況，紅色範圍表示 LOF 值大小：

# 繪圖程序import matplotlib.pyplot as pltfor k in [3,5,10]: plt.figure(k=%d%k) outliers1, inliers1 = lof(A, k=k, method = 2) plt.scatter(np.array(A)[:,0],np.array(A)[:,1],s = 10,c=b,alpha = 0.5) plt.scatter(outliers1[0],outliers1[1],s = 10+outliers1[local outlier factor]*100,c=r,alpha = 0.2) plt.title(k=%d % k)

測試 2 有輸入測試樣本情況：任務對會員的密度

數據背景：眾包任務價格制定中，地區任務的密度反映任務的密集程度、會員密度反映會員的密集程度。而任務對會員的密度則可以用於刻畫任務點周圍會員的密集程度，從而體現任務被完成的相對概率。此時訓練樣本為會員密度，測試樣本為任務密度。

import numpy as npimport pandas as pdposi = pd.read_excel(r已結束項目任務數據.xls)lon = np.array(posi["任務gps經度"][:]) # 經度lat = np.array(posi["任務gps 緯度"][:]) # 緯度A = list(zip(lat, lon)) # 按照緯度-經度匹配posi = pd.read_excel(r會員信息數據.xlsx)lon = np.array(posi["會員位置(GPS)經度"][:]) # 經度lat = np.array(posi["會員位置(GPS)緯度"][:]) # 緯度B = list(zip(lat, lon)) # 按照緯度-經度匹配# 獲取任務對會員密度，取第5鄰域，閾值為2（LOF大於2認為是離群值）outliers2, inliers2 = lof(B, A, k=5, method=2)

給定訓練樣本中共有1877條數據，測試樣本中共有835條數據，設置 LOF 閾值為 2，輸出34個離群點信息：

繪製數據點檢測情況分布如下圖所示，其中藍色表示任務分布情況，綠色表示會員分布情況，紅色範圍表示 LOF 值大小。

# 繪圖程序import matplotlib.pyplot as pltfor k,v in ([1,5],[5,2]): plt.figure(k=%d%k) outliers2, inliers2 = lof(B, A, k=k, method=v) plt.scatter(np.array(A)[:,0],np.array(A)[:,1],s = 10,c=b,alpha = 0.5) plt.scatter(np.array(B)[:,0],np.array(B)[:,1],s = 10,c=green,alpha = 0.3) plt.scatter(outliers2[0],outliers2[1],s = 10+outliers2[local outlier factor]*100,c=r,alpha = 0.2) plt.title(k = %d, method = %g % (k,v))

測試 3 測試樣本與訓練樣本互換：會員對任務的密度

數據背景：眾包任務價格制定中，地區任務的密度反映任務的密集程度、會員密度反映會員的密集程度。而任務對會員的密度則可以用於刻畫會員周圍任務的密集程度，從而體現會員能接到任務的相對概率。此時訓練樣本為任務密度，測試樣本為會員密度。

import numpy as npimport pandas as pdposi = pd.read_excel(r已結束項目任務數據.xls)lon = np.array(posi["任務gps經度"][:]) # 經度lat = np.array(posi["任務gps 緯度"][:]) # 緯度A = list(zip(lat, lon)) # 按照緯度-經度匹配posi = pd.read_excel(r會員信息數據.xlsx)lon = np.array(posi["會員位置(GPS)經度"][:]) # 經度lat = np.array(posi["會員位置(GPS)緯度"][:]) # 緯度B = list(zip(lat, lon)) # 按照緯度-經度匹配# 獲取會員對任務密度，取第5鄰域，閾值為5（LOF大於5認為是離群值）outliers3, inliers3 = lof(A, B, k=5, method=5)

給定訓練樣本中共有835條數據，測試樣本中共有1877條數據，設置 LOF 閾值為 5，輸出20個離群點信息：

繪製數據點檢測情況分布如下圖所示，其中藍色表示會員分布情況，綠色表示任務分布情況，紅色範圍表示 LOF 值大小。

# 繪圖程序plt.figure(k=5)outliers3, inliers3 = lof(A, B, k=5, method=5)plt.scatter(np.array(B)[:, 0], np.array(B)[:, 1], s=10, c=b, alpha=0.5)plt.scatter(np.array(A)[:, 0], np.array(A)[:, 1], s=10, c=green, alpha=0.3)plt.scatter(outliers3[0], outliers3[1], s=10 + outliers3[local outlier factor] * 20, c=r, alpha=0.2)plt.title(k = 5, method = 5)

將 method 設置為 0 就能輸出每一個點的 LOF 值，作為密度指標。

3.2 逐步編程實現方法

3.2.1 距離度量尺度

• distances(A, B,model = euclidean)
• 調用 geopy 模塊中的函數進行球面距離計算
• 輸入：A，B兩點坐標，距離模式（euclidean, geo 分別為歐式距離和球面距離）
• 輸出：若A，B為兩個坐標點，則返回坐標點距離；若A，B為矩陣，計算規則如下：

def distances(A, B,model = euclidean): LOF中定義的距離，默認為歐式距離，也提供球面距離 import numpy as np A = np.array(A); B = np.array(B) if model == euclidean: from scipy.spatial.distance import pdist, squareform distance = squareform(pdist(np.vstack([A, B])))[:A.shape[0],A.shape[0]:] if model == geo: from geopy.distance import great_circle distance = np.zeros(A.shape[0]*B.shape[0]).reshape(A.shape[0],B.shape[0]) for i in range(len(A)): for j in range(len(B)): distance[i,j] = great_circle(A[i], B[j]).kilometers if distance.shape == (1,1): return distance[0][0] return distance

3.2.2 鄰域半徑及鄰域點

• k_distance(k, instance_A, instance_B, result, model)
• 調用 3.2.1 的 distances 函數
• 輸入： 值，A，B兩點坐標，pandas.DataFrame 表用於存儲中間信息量，距離模式
• 輸出：pandas.DataFrame 表（距離、鄰域點坐標）

def k_distance(k, instance_A, instance_B, result, model): 計算k距離鄰域半徑及鄰域點 distance_all = distances(instance_B, instance_A, model) # 對 instance_A 中的每一個點進行遍歷 for i,a in enumerate(instance_A): distances = {} distance = distance_all[:,i] # 記錄 instance_B 到 instance_A 每一個點的距離，不重複記錄 for j in range(distance.shape[0]): if distance[j] in distances.keys(): if instance_B[j].tolist() in distances[distance[j]]: pass else: distances[distance[j]].append(instance_B[j].tolist()) else: distances[distance[j]] = [instance_B[j].tolist()] # 距離排序 distances = sorted(distances.items()) if distances[0][0] == 0: distances.remove(distances[0]) neighbours = [];k_sero = 0;k_dist = None # 截取前 k 個點 for dist in distances: k_sero += len(dist[1]) neighbours.extend(dist[1]) k_dist = dist[0] if k_sero >= k: break # 輸出信息 result.loc[str(a.tolist()),k_dist] = k_dist result.loc[str(a.tolist()),neighbours] = str(neighbours) return result

3.2.3 局部可達密度

• local_reachability_density(k,instance_A,instance_B,result, model)
• 調用 3.2.2 的 k_distance 函數
• 輸入： 值，A，B兩點坐標，pandas.DataFrame 表用於存儲中間信息量，距離模式
• 輸出：pandas.DataFrame 表（距離、鄰域點坐標、局部可達密度）

def local_reachability_density(k,instance_A,instance_B,result, model): 局部可達密度 result = k_distance(k, instance_A, instance_B, result, model) # 對 instance_A 中的每一個點進行遍歷 for a in instance_A: # 獲取k_distance中得到的鄰域點坐標，解析為點坐標（字元串 -> 數組 -> 點） try: (k_distance_value, neighbours) = result.loc[str(a.tolist())][k_dist].mean(),eval(result.loc[str(a.tolist())][neighbours]) except Exception: (k_distance_value, neighbours) = result.loc[str(a.tolist())][k_dist].mean(), eval(result.loc[str(a.tolist())][neighbours].values[0]) # 計算局部可達距離 reachability_distances_array = [0]*len(neighbours) for j, neighbour in enumerate(neighbours): reachability_distances_array[j] = max([k_distance_value, distances([a], [neighbour],model)]) sum_reach_dist = sum(reachability_distances_array) # 計算局部可達密度，並儲存結果 result.loc[str(a.tolist()),local_reachability_density] = k / sum_reach_dist return result

3.2.4 局部離群因子

• k_distance(k, instance_A, instance_B, result, model)
• 調用 3.2.3 的 local_reachability_density 函數
• 輸入： 值，A，B兩點坐標，pandas.DataFrame 表用於存儲中間信息量，距離模式
• 輸出：pandas.DataFrame 表（距離、鄰域點坐標、局部可達密度、離群因子）

def local_outlier_factor(k,instance_A,instance_B,model): 局部離群因子 result = local_reachability_density(k,instance_A,instance_B,pd.DataFrame(index=[str(i.tolist()) for i in instance_A]), model) # 判斷：若測試數據=樣本數據 if np.all(instance_A == instance_B): result_B = result else: result_B = local_reachability_density(k, instance_B, instance_B, k_distance(k, instance_B, instance_B, pd.DataFrame(index=[str(i.tolist()) for i in instance_B]), model), model) for a in instance_A: try: (k_distance_value, neighbours, instance_lrd) = result.loc[str(a.tolist())][k_dist].mean(),np.array(eval(result.loc[str(a.tolist())][neighbours])),result.loc[str(a.tolist())][local_reachability_density].mean() except Exception: (k_distance_value, neighbours, instance_lrd) = result.loc[str(a.tolist())][k_dist].mean(), np.array(eval(result.loc[str(a.tolist())][neighbours].values[0])), result.loc[str(a.tolist())][local_reachability_density].mean() finally: lrd_ratios_array = [0]* len(neighbours) for j,neighbour in enumerate(neighbours): neighbour_lrd = result_B.loc[str(neighbour.tolist())][local_reachability_density].mean() lrd_ratios_array[j] = neighbour_lrd / instance_lrd result.loc[str(a.tolist()), local_outlier_factor] = sum(lrd_ratios_array) / k return result

3.2.5 封裝函數

• lof(k,instance_A,instance_B,k_means=False,n_clusters=False,k_means_pass=3,method=1,model = euclidean
• 調用 3.2.4 的 k_distance 函數
• 輸入： 值，A，B兩點坐標，是否使用聚類演算法，聚類簇數，跳過聚類樣本數低於3的簇，閾值，距離模式
• 輸出：離群點信息、正常點信息、LOF 分布圖

# 函數中的 k_means將在第4部分介紹def lof(k, instance_A, instance_B, k_means=False, n_clusters=False, k_means_pass=3, method=1, model=euclidean): A作為定點，B作為動點 import numpy as np instance_A = np.array(instance_A); instance_B = np.array(instance_B) if np.all(instance_A == instance_B): if k_means == True: if n_clusters == True: n_clusters = elbow_method(instance_A, maxtest=10) instance_A = kmeans(instance_A, n_clusters, k_means_pass) instance_B = instance_A.copy() result = local_outlier_factor(k, instance_A, instance_B, model) outliers = result[result[local_outlier_factor] > method].sort_values(by=local_outlier_factor, ascending=False) inliers = result[result[local_outlier_factor] <= method].sort_values(by=local_outlier_factor, ascending=True) plot_lof(result, method) # 該函數見3.1.3 return outliers, inliers

4 LOF 演算法相關思考及改進

4.1 一維數據可以使用 LOF 嗎？

本文 3.1 中 sklearn 模塊提供的 LOF 方法進行訓練時會進行數據類型判斷，若數據類型為list、tuple、numpy.array 則要求傳入數據的維度至少是 2 維。實際上要篩選 1 維數據中的離群點，直接在坐標系中繪製出圖像進行閾值選取判斷也很方便。但此情形下若要使用 LOF 演算法，可以為數據添加虛擬維度，並賦相同的值：

# data 是原先的1維數據，通過下面的方法轉換為2維數據data = list(zip(data, np.zeros_like(data)))

此外，也可以通過將數據轉化為 pandas.DataFrame 形式避免上述問題：

data = pd.DataFrame(data)

4.2 多大的 LOF 才是離群值？

LOF 計算結果對於多大的值定義為離群值沒有明確的規定。在一個平穩數據集中，可能 1.1 已經是一個異常值，而在另一個具有強烈數據波動的數據集中，即使 LOF 值為 2 可能仍是一個正常值。由於方法的局限性，數據集中的異常值界定可能存在差異，筆者認為可以使用統計分布方法作為參考，再結合數據情況最終確定閾值。

基於統計分布的閾值劃分

將 LOF 異常值分數歸一化到 [0, 1] 區間，運用統計方法進行劃分下面提供使用箱型圖進行界定的方法，根據異常輸出情況參考選取。

• box(data, legend=True)
• 輸入：LOF 異常分數值，在箱型圖中繪製異常數據（設置為False不繪製）
• 輸出：異常識別情況

def box(data, legend=True): import matplotlib.pyplot as plt import pandas as pd plt.rcParams[axes.unicode_minus] = False plt.rcParams[font.sans-serif] = [SimHei] plt.style.use("ggplot") plt.figure() # 如果不是DataFrame格式，先進行轉化 if type(data) != pd.core.frame.DataFrame: data = pd.DataFrame(data) p = data.boxplot(return_type=dict) warming = pd.DataFrame() y = p[fliers][0].get_ydata() y.sort() for i in range(len(y)): if legend == True: plt.text(1, y[i] - 1, y[i], fontsize=10, color=black, ha=right) if y[i] < data.mean()[0]: form = 低 else: form = 高 warming = warming.append(pd.Series([y[i], 偏 + form]).T, ignore_index=True) print(warming) plt.show()

box 函數可以插入封裝函數 lof 中，傳入 data = predict[local outlier factor] 實現；也可以先隨機指定一個初始閾值，（輸出的離群點、正常點分別命名為outliers, inliers）再輸入：

box(outliers[local outlier factor].tolist()+inliers[local outlier factor].tolist(), legend=True)

此時，交互控制台中輸出情況如下左圖所示，箱型圖如下右圖所示。輸出情況提示我們從數據分布的角度上，可以將 1.4 作為離群識別閾值，但實際上取 7 更為合適（從 2 到 7 間有明顯的斷層，而上文中設定為 5 是經過多次試驗後選取的數值）。

4.3 數據維度過大，還能使用 LOF 演算法嗎？

數據維度過大一方面會增大量綱的影響，另一方面增大計算難度。此時直接使用距離度量的表達形式不合理，並有人為放大較為分散數據影響的風險。一種處理方式是採用馬氏距離作為距離度量的方式（去量綱化）。另一種處理方式，參考隨機森林的決策思想，可以考慮在多個維度投影上使用 LOF 並將結果結合起來，以提高高維數據的檢測質量。

集成學習

集成學習：通過構建並結合多個學習器來完成學習任務，通常可以獲得比單一學習器更顯著優越的泛化性能。

——周志華《機器學習》

數據檢測中進行使用的數據應該是有意義的數據，這就需要進行簡單的特徵篩選，否則無論多麼「離群」的樣本，可能也沒有多大的實際意義。根據集成學習的思想，需要將數據按維度拆分，對於同類型的數據，這裡假設你已經做好了規約處理（如位置坐標可以放在一起作為一個特徵「距離」進行考慮），並且數據的維度大於 1，否則使用 4.1 中的數據變換及一般形式 LOF 即可處理。

4.3.1 投票表決模式

投票表決模式認為每一個維度的數據都是同等重要，單獨為每個維度數據設置 LOF 閾值並進行比對，樣本的 LOF 值超過閾值則異常票數積 1 分，最終超過票數閾值的樣本認為是離群樣本。

• localoutlierfactor(data, predict, k)
• 輸入：訓練樣本，測試樣本， 值
• 輸出：每一個測試樣本 LOF 值及對應的第 距離
• plot_lof(result,method)
• 輸入：LOF 值、閾值
• 輸出：以索引為橫坐標的 LOF 值分布圖
• ensemble_lof(data, predict=None, k=5, groups=[], method=1, vote_method = auto)
• 輸入：訓練樣本，測試樣本， 值，組合特徵索引，離群閾值，票數閾值
• 組合特徵索引（列位置 - 1）：如第一列數據與第二列數據作為同類型數據，則傳入groups = [[0, 1]]
• 離群閾值：每一個特徵的離群閾值，缺少位數則以 1 替代
• 票數閾值：正常點離群因子得票數上限
• 輸出：離群點、正常點分類情況
• 沒有輸入測試樣本時，默認測試樣本=訓練樣本

def localoutlierfactor(data, predict, k, group_str): from sklearn.neighbors import LocalOutlierFactor clf = LocalOutlierFactor(n_neighbors=k + 1, algorithm=auto, contamination=0.1, n_jobs=-1) clf.fit(data) # 記錄 LOF 離群因子，做相反數處理 predict[local outlier factor %s % group_str] = -clf._decision_function(predict.iloc[:, eval(group_str)]) return predictdef plot_lof(result, method, group_str): import matplotlib.pyplot as plt plt.rcParams[font.sans-serif] = [SimHei] # 用來正常顯示中文標籤 plt.rcParams[axes.unicode_minus] = False # 用來正常顯示負號 plt.figure(local outlier factor %s % group_str) try: plt.scatter(result[result > method].index, result[result > method], c=red, s=50, marker=., alpha=None, label=離群點) except Exception: pass try: plt.scatter(result[result <= method].index, result[result <= method], c=black, s=50, marker=., alpha=None, label=正常點) except Exception: pass plt.hlines(method, -2, 2 + max(result.index), linestyles=--) plt.xlim(-2, 2 + max(result.index)) plt.title(LOF局部離群點檢測, fontsize=13) plt.ylabel(局部離群因子, fontsize=15) plt.legend() plt.show()def ensemble_lof(data, predict=None, k=5, groups=[], method=1, vote_method = auto): import pandas as pd import numpy as np # 判斷是否傳入測試數據，若沒有傳入則測試數據賦值為訓練數據 try: if predict == None: predict = data.copy() except Exception: pass data = pd.DataFrame(data); predict = pd.DataFrame(predict) # 數據標籤分組，默認獨立自成一組 for i in range(data.shape[1]): if i not in pd.DataFrame(groups).values: groups += [[i]] # 擴充閾值列表 if type(method) != list: method = [method] method += [1] * (len(groups) - 1) else: method += [1] * (len(groups) - len(method)) vote = np.zeros(len(predict)) # 計算LOF離群因子並根據閾值進行票數統計 for i in range(len(groups)): predict = localoutlierfactor(pd.DataFrame(data).iloc[:, groups[i]], predict, k, str(groups[i])) plot_lof(predict.iloc[:, -1], method[i], str(groups[i])) vote += predict.iloc[:, -1] > method[i] # 根據票數閾值劃分離群點與正常點 predict[vote] = vote if vote_method == auto: vote_method = len(groups)/2 outliers = predict[vote > vote_method].sort_values(by=vote) inliers = predict[vote <= vote_method].sort_values(by=vote) return outliers, inliers

測試 4 仍然使用測試3的情況進行分析，此時將經度、緯度設置為獨立的特徵，分別對兩個維度數據進行識別（儘管單獨的緯度、經度數據沒有太大的實際意義）

import numpy as npimport pandas as pdposi = pd.read_excel(r已結束項目任務數據.xls)lon = np.array(posi["任務gps經度"][:]) # 經度lat = np.array(posi["任務gps 緯度"][:]) # 緯度A = list(zip(lat, lon)) # 按照緯度-經度匹配posi = pd.read_excel(r會員信息數據.xlsx)lon = np.array(posi["會員位置(GPS)經度"][:]) # 經度lat = np.array(posi["會員位置(GPS)緯度"][:]) # 緯度B = list(zip(lat, lon)) # 按照緯度-經度匹配# 獲取會員對任務密度，取第5鄰域，閾值分別為 1.5，2，得票數超過 1 的認為是異常點 outliers4, inliers4 = ensemble_lof(A, B, k=5, method=[1.5,2], vote_method = 1)# 繪圖程序plt.figure(投票集成 LOF 模式)plt.scatter(np.array(B)[:, 0], np.array(B)[:, 1], s=10, c=b, alpha=0.5)plt.scatter(np.array(A)[:, 0], np.array(A)[:, 1], s=10, c=green, alpha=0.3)plt.scatter(outliers4[0], outliers4[1], s=10 + 1000, c=r, alpha=0.2)plt.title(k = 5, method = [1.5, 2])

4.3.2 LOF 異常分數加權模式

異常分數加權模式則是對各維度數據的 LOF 值進行加權，獲取最終的 LOF 得分作為整體數據的 LOF 得分。權重可以認為是特徵的重要程度，也可以認為是數據分布的相對離散程度，若視為後面一種情形，可以根據熵權法進行設定，關於熵權法的介紹詳見筆者另一篇博文。

式中 表示第 個數據的第 維度 LOF 異常分數值。

• localoutlierfactor(data, predict, k)
• 輸入：訓練樣本，測試樣本， 值
• 輸出：每一個測試樣本 LOF 值及對應的第 距離
• plot_lof(result,method)
• 輸入：LOF 值、閾值
• 輸出：以索引為橫坐標的 LOF 值分布圖
• ensemble_lof(data, predict=None, k=5, groups=[], method=2, weight=1)
• 輸入：訓練樣本，測試樣本， 值，組合特徵索引，離群閾值，特徵權重
• 組合特徵索引（列位置 - 1）：如第一列數據與第二列數據作為同類型數據，則傳入groups = [[0, 1]]
• 離群閾值：加權 LOF 離群閾值，默認為 2
• 特徵權重：為每個維度的 LOF 設置權重，缺少位數則以 1 代替
• 輸出：離群點、正常點分類情況
• 沒有輸入測試樣本時，默認測試樣本=訓練樣本

def localoutlierfactor(data, predict, k, group_str): from sklearn.neighbors import LocalOutlierFactor clf = LocalOutlierFactor(n_neighbors=k + 1, algorithm=auto, contamination=0.1, n_jobs=-1) clf.fit(data) # 記錄 LOF 離群因子，做相反數處理 predict[local outlier factor %s % group_str] = -clf._decision_function(predict.iloc[:, eval(group_str)]) return predictdef plot_lof(result, method): import matplotlib.pyplot as plt plt.rcParams[font.sans-serif] = [SimHei] # 用來正常顯示中文標籤 plt.rcParams[axes.unicode_minus] = False # 用來正常顯示負號 plt.scatter(result[result > method].index, result[result > method], c=red, s=50, marker=., alpha=None, label=離群點) plt.scatter(result[result <= method].index, result[result <= method], c=black, s=50, marker=., alpha=None, label=正常點) plt.hlines(method, -2, 2 + max(result.index), linestyles=--) plt.xlim(-2, 2 + max(result.index)) plt.title(LOF局部離群點檢測, fontsize=13) plt.ylabel(局部離群因子, fontsize=15) plt.legend() plt.show()def ensemble_lof(data, predict=None, k=5, groups=[], method=auto, weight=1): import pandas as pd # 判斷是否傳入測試數據，若沒有傳入則測試數據賦值為訓練數據 try: if predict == None: predict = data except Exception: pass data = pd.DataFrame(data); predict = pd.DataFrame(predict) # 數據標籤分組，默認獨立自成一組 for i in range(data.shape[1]): if i not in pd.DataFrame(groups).values: groups += [[i]] # 擴充權值列表 if type(weight) != list: weight = [weight] weight += [1] * (len(groups) - 1) else: weight += [1] * (len(groups) - len(weight)) predict[local outlier factor] = 0 # 計算LOF離群因子並根據特徵權重計算加權LOF得分 for i in range(len(groups)): predict = localoutlierfactor(pd.DataFrame(data).iloc[:, groups[i]], predict, k, str(groups[i])) predict[local outlier factor] += predict.iloc[:, -1] * weight[i] if method == auto: method = sum(weight) plot_lof(predict[local outlier factor], method) # 根據離群閾值劃分離群點與正常點 outliers = predict[predict[local outlier factor] > method].sort_values(by=local outlier factor) inliers = predict[predict[local outlier factor] <= method].sort_values(by=local outlier factor) return outliers, inliers

測試 5 仍然使用測試3的情況進行分析，此時將經度、緯度設置為獨立的特徵，分別對兩個維度數據進行識別（儘管單獨的緯度、經度數據似乎沒有太大的實際意義）

import numpy as npimport pandas as pdposi = pd.read_excel(r已結束項目任務數據.xls)lon = np.array(posi["任務gps經度"][:]) # 經度lat = np.array(posi["任務gps 緯度"][:]) # 緯度A = list(zip(lat, lon)) # 按照緯度-經度匹配posi = pd.read_excel(r會員信息數據.xlsx)lon = np.array(posi["會員位置(GPS)經度"][:]) # 經度lat = np.array(posi["會員位置(GPS)緯度"][:]) # 緯度B = list(zip(lat, lon)) # 按照緯度-經度匹配# 獲取會員對任務密度，取第5鄰域，閾值為 100，權重分別為5，1 outliers5, inliers5 = ensemble_lof(A, B, k=5, method=100,weight = [5,1])# 繪圖程序plt.figure(LOF 異常分數加權模式)plt.scatter(np.array(B)[:, 0], np.array(B)[:, 1], s=10, c=b, alpha=0.5)plt.scatter(np.array(A)[:, 0], np.array(A)[:, 1], s=10, c=green, alpha=0.3)plt.scatter(outliers5[0], outliers5[1], s=10 + outliers5[local outlier factor], c=r, alpha=0.2)plt.title(k = 5, method = 100)

4.3.3 混合模式

混合模式適用於數據中有些特徵同等重要，有些特徵有重要性區別的情況，即對 4.3.1、4.3.2 情形綜合進行考慮。同等重要的數據將使用投票表決模式，重要程度不同的數據使用加權模式並根據閾值轉換為投票表決模式。程序上只需將兩部分混合使用即可，本文在此不做展示。

4.4 數據量太大，演算法執行效率過低，有什麼改進方法嗎？

LOF 演算法在檢測離群點的過程中，遍歷整個數據集以計算每個點的 LOF 值，使得演算法運算速度慢。同時，由於數據正常點的數量一般遠遠多於離群點的數量，而 LOF 方法通過比較所有數據點的 LOF 值判斷離群程度，產生了大量沒必要的計算。因此，通過對原始數據進行修剪可以有效提高 LOF 方法的計算效率。此外，實踐過程中也發現通過數據集修剪後，可以大幅度減少數據誤判為離群點的幾率。這種基於聚類修剪得離群點檢測方法稱為 CLOF (Cluster-Based Local Outlier Factor) 演算法。

基於 K-Means 的 CLOF 演算法

在應用 LOF 演算法前，先用 K-Means 聚類演算法，將原始數據聚成 n_clusters 簇。對其中的每一簇，計算簇的中心 ，求出該簇中所有點到該中心的平均距離並記為該簇的半徑 。對該類中所有點，若該點到簇中心的距離大於等於 則將其放入「離群點候選集」 ，最後對 中的數據使用 LOF 演算法計算離群因子。

設第 個簇中的點的個數為 ，點集為

中心和半徑的計算公式如下：

如何確定最佳的 n_clusters —— 肘部法則

K-Means演算法通過指定聚類簇數n_clusters及隨機生成聚類中心，對最靠近他們的對象進行迭代歸類，逐次更新各聚類中心的值，直到最好的聚類效果（代價函數值最小）。

對於n_clusters的選取將直接影響演算法的聚類效果。肘部法則將不同的n_clusters值的成本函數值刻畫出來，隨著n_clusters增大，每個簇包含的樣本數會減少，樣本離其中心更接近，代價函數會減小。但隨著n_clusters繼續增大，代價函數的改善程度不斷下降（在圖像中，代價函數曲線趨於平穩）。n_clusters值增大過程中，代價函數改善程度最大的位置對應的n_clusters就是肘部，使用此n_clusters一般可以取得不錯的效果。但肘部法則的使用僅僅是從代價水平進行考慮，有時候還需結合實際考慮。

由於離群值樣本數量一般較少，如果聚類出來的簇中樣本量太少（如 1-4 個，但其他簇有成百上千個樣本），則這種聚類簇不應進行修剪。

定義代價函數：

• elbow_method(data,maxtest = 11)
• 輸入：需要修剪的數據集，最大測試範圍（默認為11）
• 輸出：代價曲線圖

def elbow_method(data,maxtest = 11): 使用肘部法則確定n_clusters值 from sklearn.cluster import KMeans from scipy.spatial.distance import cdist import numpy as np import matplotlib.pyplot as plt plt.rcParams[font.sans-serif] = [SimHei] plt.rcParams[axes.unicode_minus] = False ax = plt.figure(figsize=(8,4)).add_subplot(111) N_test = range(1, maxtest) # 代價函數值列表 meandistortions = [] for n_clusters in N_test: model = KMeans(n_clusters=n_clusters).fit(data) # 計算代價函數值 meandistortions.append(sum(np.min(cdist(data, model.cluster_centers_, euclidean), axis=1)) / len(data)) plt.plot(N_test, meandistortions, bx-,alpha = 0.4) plt.xlabel(k) plt.ylabel(代價函數,fontsize= 12) plt.title(用肘部法則來確定最佳的n_clusters值,fontsize= 12) ax.spines[top].set_visible(False) ax.spines[right].set_visible(False) ax.set_xticks(np.arange(0, maxtest, 1)) plt.show()

• kmeans(data, n_clusters, m):
• 輸入：需要修剪的數據集，聚類簇數，最小修剪簇應包含的樣本點數
• 輸出：修剪後的數據集

def kmeans(data, n_clusters, m): 使用K-Means演算法修剪數據集 from sklearn.cluster import KMeans import numpy as np data_select = [] model = KMeans(n_clusters=n_clusters).fit(data) centeroids = model.cluster_centers_ label_pred = model.labels_ import collections for k, v in collections.Counter(label_pred).items(): if v < m: data_select += np.array(data)[label_pred == k].tolist() else: distance = np.sqrt(((np.array(data)[label_pred == k] - centeroids[k]) ** 2).sum(axis=1)) R = distance.mean() data_select += np.array(data)[label_pred == k][distance >= R].reshape(-1, np.array(data).shape[-1]).tolist() return np.array(data_select)

測試 6 對B數據集進行修剪分析，B數據集共有 1877 條數據

# 肘部法則確定最佳修剪集elbow_method(B,maxtest = 11)

# 根據上面設定的 n_clusters 為3，最小樣本量設置為3B_cut = kmeans(B, n_clusters = 3, m = 3)

執行上述程序，原先包含 1877 條數據的 B 數據集修剪為含有 719 條數據的較小的數據集。使用 LOF 演算法進行離群檢測，檢測結果如下：

# 獲取會員分布密度，取第10鄰域，閾值為3（LOF大於3認為是離群值）outliers6, inliers6 = lof(B_cut, k=10, method=3)# 繪圖程序plt.figure(CLOF 離群因子檢測)plt.scatter(np.array(B)[:, 0], np.array(B)[:, 1], s=10, c=b, alpha=0.5)plt.scatter(outliers6[0], outliers6[1], s=10 + outliers6[local outlier factor]*10, c=r, alpha=0.2)plt.title(k = 10, method = 3)

修剪後的數據集 LOF 意義不再那麼明顯，但離群點的 LOF 仍然會是較大的值，並且 選取越大的值，判別效果越明顯。

4.5 如何提高識別精度

合理增大 值能顯著提高識別精度。但 值的增加會帶來不必要的計算、影響演算法執行效率，也正因此本文所用的 值都取較小。合理選取 與閾值將是 LOF 成功與否的關鍵。

5 後記

本文內容主要參考演算法原文及筆者學習經驗進行總結。在異常識別領域，LOF 演算法和 Isolation Forest 演算法已經被指出是性能最優、識別效果最好的演算法。對於常用的人群密度（或其他）的刻畫，LOF異常分數值不失為一種「高端」方法（參考文獻[3]），相比傳統方法，更具有成熟的理論支撐。

6 參考文獻

[1] 維基百科. Moore–Penrose inverse[EB/OL]. [2018-6-2]. https://en.wikipedia.org/wiki/Moore%E2%80%93Penrose_inverse

[2] Breunig M M, Kriegel H P, Ng R T, et al. LOF: identifying density-based local outliers[C]// ACM SIGMOD International Conference on Management of Data. ACM, 2000:93-104.

[3] 董天文，潘偉堤，戚銘珈，張曉敏. 「拍照賺錢」的任務定價. 教育部中國大學生在線網站[J/OL]. [2017-11-13]. http://upload.univs.cn/2017/1113/1510570109509.pdf

---

作者：張柳彬 仇禮鴻

如有疑問，請聯繫QQ：965579168

轉載請聲明出處

推薦閱讀：