局部异常因子(Local Outlier Factor, LOF)算法详解及实验
局部异常因子(Local Outlier Factor, LOF)通过计算样本点的局部相对密度来衡量这个样本点的异常情况,可以算是一类无监督学习算法。下面首先对算法的进行介绍,然后进行实验。
LOF算法
下面介绍LOF算法的每个概念,以样本点集合中的样本点$P$为例。下面的概念名称中都加了一个k-,实际上部分名称原文没有加,但是感觉这样更严谨一些。
k-邻近距离(k-distance):样本点$P$与其最近的第$k$个样本点之间的距离,表示为$d_k(P)$。其中距离可以用各种方式度量,通常使用欧氏距离。
k-距离邻域:以$P$为圆心,$d_k(P)$为半径的邻域,表示为$N_k(P)$。
k-可达距离:$P$到某个样本点$O$的k-可达距离,取$d_k(O)$或$P$与$O$之间距离的较大值,表示为
$reach\_dist_k(P,O)=\max\{d_k(O),d(P,O)\}$
也就是说,如果$P$在$N_k(O)$内部,$reach\_dist_k(P,O)$取$O$的k-邻近距离$d_k(O)$,在$N_k(O)$外部则取$P$与$O$之间距离$d(P,O)$。需要注意k-可达距离不是对称的。
k-局部可达密度(local reachability density, lrd):$P$的k-局部可达密度表示为
$\displaystyle lrd_k(P)=\left(\frac{\sum\limits_{O\in N_k(P)}reach\_dist_k(P,O)}{|N_k(P)|}\right)^{-1}$
括号内,分子计算了$P$到其k-距离邻域内所有样本点$O$的k-可达距离之和,然后除以$P$的k-距离邻域内部的样本点个数进行平均。再加一个倒数,表示为密度,即$P$到每个点的平均距离越小,密度越大。可以推理出,如果$P$在所有$O$的k-邻域内部,其局部可达密度即为
$\displaystyle lrd_k(P)=\left(\frac{\sum\limits_{O\in N_k(P)}d_k(O)}{|N_k(P)|}\right)^{-1}$
可以看出,如果$P$是一个离群点,那么它不太可能存在于$N_k(P)$中各点的k-距离邻域内,从而导致其局部可达密度偏小;如果$P$不是离群点,其局部可达密度最大取为上式。
实际上我有点奇怪为什么要用一个最大值来将距离作一个限制,也就是使用k-可达距离,而不是直接使用距离,即定义局部密度为下式
$\displaystyle ld_k(P)=\left(\frac{\sum\limits_{O\in N_k(P)}d(O,P)}{|N_k(P)|}\right)^{-1}$
k-局部异常因子(Local Outlier Factor, LOF):$P$的k-局部异常因子表示为
$\displaystyle LOF_k(P)=\frac{\frac{1}{|N_k(P)|}\sum\limits_{O\in N_k(P)}lrd_k(O)}{lrd_k(P)}$
从直觉上理解:当$LOF_k(P)\le1$时,表明$P$处密度比其周围点的平均密度大或相当,则$P$是内点;当$LOF_k(P)>1$时,表明$P$处密度比其周围点小,可以判别为离群(异常)点。
实验
LOF算法实现
实验设置样本维度为2以便可视化。由于样本点只包含连续值,实验默认设置$|N_k(P)|=k$。设置$k=5$,并将阈值设为2,即将LOF大于2的样本点视作异常。函数定义、抽样、计算以及可视化代码如下。
#%% 定义函数 import torch import matplotlib.pyplot as plt #计算所有样本点[N, M]之间的距离,得到[N, N] def get_dists(points:torch.Tensor): x = torch.sum(points**2, 1).reshape(-1, 1) y = torch.sum(points**2, 1).reshape(1, -1) dists = x + y - 2 * torch.mm(points, points.permute(1,0)) #数值计算问题,防止对角线小于0 dists = dists - torch.diag_embed(torch.diag(dists)) return torch.sqrt(dists) #计算所有样本的LOF def get_LOFs(dists:torch.Tensor, k): #距离排序,获取所有样本点的k-临近距离 sorted_dists, sorted_inds = torch.sort(dists, 1, descending=False) k_dists = sorted_dists[:, k] neighbor_inds = sorted_inds[:, 1:k+1].reshape(-1) neighbor_k_dists = k_dists[neighbor_inds].reshape(-1, k) neighbor_k_reach_dists = torch.max(neighbor_k_dists, sorted_dists[:, 1:k+1]) lrds = (neighbor_k_reach_dists.sum(1)/k)**-1 neighbor_lrds = lrds[neighbor_inds].reshape(-1, k) LOFs = neighbor_lrds.sum(1)/k/lrds return LOFs #%% 随机生成聚集点和异常点 from torch.distributions import MultivariateNormal torch.manual_seed(0) crowd_mu_covs = [ [[0.0, 0.0], [[1.0, 0.0], [0.0, 2.0]], 10], [[-10.0, -1.0], [[2.0, 0.8], [0.8, 2.0]], 10], [[-10.0, -20.0], [[5.0, -2], [-2, 3.0]], 50], [[5.0, -10.0], [[5.0, -2], [-2, 3.0]], 50], [[-3.0, -10.0], [[5.0, -2], [-2, 3.0]], 50], [[-13.0, -10.0], [[0.3, -0.1], [-0.1, 0.5]], 10], [[-4.0, -10.0], [[0.3, -0.1], [-0.1, 0.1]], 100], ]#正态分布点的均值和方差 outliers = [[5, 5.], [3, 4], [5, -3], [4, -30], [-2, -35], [-10, -35]] #异常点 points = [] for i in crowd_mu_covs: mu = torch.tensor(i[0]) cov = torch.tensor(i[1]) ps = MultivariateNormal(mu, cov).sample([i[2]]).to('cuda') points.append(ps) for o in outliers: points.append(torch.tensor([o]).to('cuda')) points = torch.cat(points) #%% 等高线图 k, threshold = 5, 2 x = torch.arange(-16, 11, 0.5) y = torch.arange(-36, 6, 0.5) X, Y = torch.meshgrid(x, y) Z = torch.zeros_like(X).to('cuda') M = torch.stack([X,Y]).permute(1,2,0).to('cuda') for i in range(len(x)): for j in range(len(y)): ps = torch.cat([points, M[i,j:j+1]]) dists = get_dists(ps) LOFs = get_LOFs(dists, k) Z[i,j] = LOFs[-1] plt.contourf(X.cpu().numpy(),Y.cpu().numpy(),Z.cpu().numpy()) # 点 dists = get_dists(points) LOFs = get_LOFs(dists, k) for i, p in enumerate(points.cpu().numpy()): shape, color = '.', 'black' if len(points) - i <= len(outliers): shape = '^' plt.annotate("%.2f"%LOFs[i].cpu().numpy(), (p[0], p[1])) if LOFs[i] > threshold: color = 'red' plt.annotate("%.2f"%LOFs[i].cpu().numpy(), (p[0], p[1]), color='blue') plt.plot(p[0], p[1], shape, color=color) plt.show()
实验可视化结果如下图所示,其中红色点表示被标为异常的点,三角点表示实验设置的真实异常点。
距离代替局部可达距离
根据前面的疑问,用距离代替局部可达距离进行相应实验。仅在get_LOFs函数处做了相关改动,并将阈值threshold改为2.5。get_LOFs函数修改如下。
def get_LOFs(dists:torch.Tensor, k): #距离排序,获取所有样本点的k-临近距离 sorted_dists, sorted_inds = torch.sort(dists, 1, descending=False) densities = (sorted_dists[:, 1:k+1].sum(1)/k)**-1 neighbor_inds = sorted_inds[:, 1:k+1].reshape(-1) neighbor_densities = densities[neighbor_inds].reshape(-1, k) LOFs = neighbor_densities.sum(1)/k/densities return LOFs
可视化结果如下图所示
效果看起来和原始LOF差不多。至于用局部可达距离的理论解释,本文不再作深究,欢迎前来讨论。