source: branches/DataAnalysis/HeuristicLab.Problems.DataAnalysis.MultiVariate/3.3/Evaluators/OnlineMeanMahalanobisDistanceEvaluator.cs @ 4475

using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using HeuristicLab.Core;
using HeuristicLab.Data;
using HeuristicLab.Problems.DataAnalysis.Evaluators;
using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
using HeuristicLab.Parameters;

namespace HeuristicLab.Problems.DataAnalysis.MultiVariate.Evaluators {
  public class OnlineMeanMahalanobisDistanceEvaluator : IMultiVariateOnlineEvaluator {
    private int n;
    private double distance;
    private double[,] covMatrix;
    public double MeanMahalanobisDistance {
      get {
        if (n == 0) throw new InvalidOperationException("no elements");
        else
          return Math.Sqrt(distance) / n;
      }
    }

    public double MeanGeneralizedSquaredInterpointDistance {
      get {
        if (n == 0) throw new InvalidOperationException("no elements");
        else
          return distance / n;
      }
    }

    public OnlineMeanMahalanobisDistanceEvaluator() {
      Reset();
    }

    #region IMultiVariateOnlineEvaluator Members
    public double Value {
      get { return MeanMahalanobisDistance; }
    }

    public void Add(double[] original, double[] estimated) {
      if (original.Length != estimated.Length) throw new ArgumentException("Number of elements of original and estimated doesn't match.");
      if (original.Length != covMatrix.GetLength(0)) throw new ArgumentException("Original and estimated is not compatible with covariance matrix.");

      double[] diff = new double[original.Length];
      for (int i = 0; i < diff.Length; i++) {
        diff[i] = original[i] - estimated[i];
      }
      double[] target = new double[original.Length];

      alglib.ablas.rmatrixmv(covMatrix.GetLength(0), covMatrix.GetLength(1), ref covMatrix, 0, 0, 0, ref diff, 0, ref target, 0);

      double sum = 0.0;
      for (int i = 0; i < original.Length; i++) {
        sum += diff[i] * target[i];
      }
      distance += sum;
      n++;
    }

    public void Reset() {
      n = 0;
      distance = 0.0;
    }

    #endregion

    public void InitializeCovarianceMatrixFromSamples(params IEnumerable<double>[] samples) {
      covMatrix = new double[samples.Length, samples.Length];
      OnlineCovarianceEvaluator covEvaluator = new OnlineCovarianceEvaluator();
      for (int i = 0; i < samples.Length; i++) {
        for (int j = i; j < samples.Length; j++) {
          var xEnumerator = samples[i].GetEnumerator();
          var yEnumerator = samples[j].GetEnumerator();
          covEvaluator.Reset();
          while (xEnumerator.MoveNext() & yEnumerator.MoveNext()) {
            covEvaluator.Add(xEnumerator.Current, yEnumerator.Current);
          }
          if (xEnumerator.MoveNext() | yEnumerator.MoveNext()) {
            throw new ArgumentException("Number of elements must be the same in all enumerations.");
          }
          covMatrix[i, j] = covEvaluator.Covariance;
          covMatrix[j, i] = covEvaluator.Covariance;
        }
      }
      int info = 0;
      alglib.matinv.matinvreport report = new alglib.matinv.matinvreport();
      alglib.matinv.rmatrixinverse(ref covMatrix, covMatrix.GetLength(0), ref info, ref report);
      if (info != 1) throw new InvalidOperationException("Can't invert covariance matrix.");
    }
  }
}