局部异常因子（Local Outlier Factor, LOF）算法详解及实验

#%% 定义函数
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     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()

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