MultivariateGaussianAnomalyDetection.py

from __future__ import division
import numpy as np
import scipy.stats
import math


__author__ = 'herve.schnegg'


class MultiGaussianAnomalyDetection(object):


    m = np.array([])
    sigma = np.array([])
    multi_norm_dist = []


    def __init__(self):
        pass


    def fit(self, X):
        self.m = np.mean(X, axis=0)
        self.sigma = np.cov(np.transpose(X))
        self.multi_norm_dist = scipy.stats.multivariate_normal(mean=self.m, cov=self.sigma)


    def predict(self, X, epsilon=0, bottom_n=0, return_density_estimate=False):
        number_obs = np.array(X).shape[0]
        n = np.array(X).shape[1]
        pred = np.zeros(number_obs)

        density_estimate = np.apply_along_axis(lambda row: np.prod(map(lambda col: self.multi_norm_dist.pdf(col), row)), 1, X)

        if epsilon > 0:
            pred[density_estimate < epsilon] = 1 # todo: very slow
        if bottom_n > 0:
            rank = scipy.stats.rankdata(density_estimate, method='ordinal') # todo: very slow
            pred[rank <= bottom_n] = 1

        if return_density_estimate:
            output = np.vstack((density_estimate, pred))
        else:
            output = pred

        output = np.transpose(output)
            
        return output