GaussianAnomalyDetection.py

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


__author__ = 'herve.schnegg'


class GaussianAnomalyDetection(object):


    m = np.array([])
    sd = np.array([])
    norm_dist = []


    def __init__(self):
        pass


    def fit(self, X):
        self.m = np.mean(X, axis=0)
        self.sd = np.std(X, axis=0)
        if len(X.shape) > 1:
            self.norm_dist = map(scipy.stats.norm, self.m, self.sd)
        else:
            self.norm_dist = scipy.stats.norm(self.m, self.sd)


    def predict(self, X, epsilon=0, bottom_n=0, return_density_estimate=False):
        number_obs = np.array(X).shape[0]
        pred = np.zeros(number_obs)
        
        if len(X.shape) > 1:
            density_estimate = np.apply_along_axis(lambda row: np.prod(map(lambda norm_dist, col: norm_dist.pdf(col), self.norm_dist, row)), 1, X)
        else:
            density_estimate = map(self.norm_dist.pdf, 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