source: branches/HeuristicLab.TimeSeries/HeuristicLab.Problems.DataAnalysis/3.4/Implementation/TimeSeriesPrognosis/TimeSeriesPrognosisSolutionBase.cs @ 7160

#region License Information
/* HeuristicLab
 * Copyright (C) 2002-2011 Heuristic and Evolutionary Algorithms Laboratory (HEAL)
 *
 * This file is part of HeuristicLab.
 *
 * HeuristicLab is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * HeuristicLab is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.
 */
#endregion

using System.Collections.Concurrent;
using System.Collections.Generic;
using System.Linq;
using HeuristicLab.Common;
using HeuristicLab.Data;
using HeuristicLab.Optimization;
using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;

namespace HeuristicLab.Problems.DataAnalysis {
  [StorableClass]
  public abstract class TimeSeriesPrognosisSolutionBase : DataAnalysisSolution, ITimeSeriesPrognosisSolution {
    private const string TrainingMeanSquaredErrorResultName = "Mean squared error (training)";
    private const string TestMeanSquaredErrorResultName = "Mean squared error (test)";
    private const string TrainingMeanAbsoluteErrorResultName = "Mean absolute error (training)";
    private const string TestMeanAbsoluteErrorResultName = "Mean absolute error (test)";
    private const string TrainingSquaredCorrelationResultName = "Pearson's R² (training)";
    private const string TestSquaredCorrelationResultName = "Pearson's R² (test)";
    private const string TrainingRelativeErrorResultName = "Average relative error (training)";
    private const string TestRelativeErrorResultName = "Average relative error (test)";
    private const string TrainingNormalizedMeanSquaredErrorResultName = "Normalized mean squared error (training)";
    private const string TestNormalizedMeanSquaredErrorResultName = "Normalized mean squared error (test)";
    private const string TrainingDirectionalSymmetryResultName = "Average directional symmetry (training)";
    private const string TestDirectionalSymmetryResultName = "Average directional symmetry (test)";
    private const string TrainingWeightedDirectionalSymmetryResultName = "Average weighted directional symmetry (training)";
    private const string TestWeightedDirectionalSymmetryResultName = "Average weighted directional symmetry (test)";
    private const string TrainingTheilsUStatisticLastResultName = "Average Theil's U (last) (training)";
    private const string TestTheilsUStatisticLastResultName = "Average Theil's U (last) (test)";
    private const string TrainingTheilsUStatisticMeanResultName = "Average Theil's U (mean) (training)";
    private const string TestTheilsUStatisticMeanResultName = "Average Theil's U (mean) (test)";
    private const string TrainingTheilsUStatisticMAResultName = "Average Theil's U (moving average) (training)";
    private const string TestTheilsUStatisticMAResultName = "Average Theil's U (moving average) (test)";

    public new ITimeSeriesPrognosisModel Model {
      get { return (ITimeSeriesPrognosisModel)base.Model; }
      protected set { base.Model = value; }
    }

    public new ITimeSeriesPrognosisProblemData ProblemData {
      get { return (ITimeSeriesPrognosisProblemData)base.ProblemData; }
      set { base.ProblemData = value; }
    }

    [Storable]
    private int horizon;
    public int Horizon {
      get { return horizon; }
      set {
        if (horizon != value) {
          horizon = value;
          RecalculateResults();
        }
      }
    }

    public abstract IEnumerable<IEnumerable<double>> PrognosedTrainingValues { get; }
    public abstract IEnumerable<IEnumerable<double>> PrognosedTestValues { get; }
    public abstract IEnumerable<IEnumerable<IEnumerable<double>>> GetPrognosedValues(IEnumerable<int> rows, int horizon);

    #region Results
    public double[] TrainingMeanSquaredError {
      get { return ((DoubleArray)this[TrainingMeanSquaredErrorResultName].Value).ToArray(); }
      private set { this[TrainingMeanSquaredErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TestMeanSquaredError {
      get { return ((DoubleArray)this[TestMeanSquaredErrorResultName].Value).ToArray(); }
      private set { this[TestMeanSquaredErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingMeanAbsoluteError {
      get { return ((DoubleArray)this[TrainingMeanAbsoluteErrorResultName].Value).ToArray(); }
      private set { this[TrainingMeanAbsoluteErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TestMeanAbsoluteError {
      get { return ((DoubleArray)this[TestMeanAbsoluteErrorResultName].Value).ToArray(); }
      private set { this[TestMeanAbsoluteErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingRSquared {
      get { return ((DoubleArray)this[TrainingSquaredCorrelationResultName].Value).ToArray(); }
      private set { this[TrainingSquaredCorrelationResultName].Value = new DoubleArray(value); }
    }
    public double[] TestRSquared {
      get { return ((DoubleArray)this[TestSquaredCorrelationResultName].Value).ToArray(); }
      private set { this[TestSquaredCorrelationResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingRelativeError {
      get { return ((DoubleArray)this[TrainingRelativeErrorResultName].Value).ToArray(); }
      private set { this[TrainingRelativeErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TestRelativeError {
      get { return ((DoubleArray)this[TestRelativeErrorResultName].Value).ToArray(); }
      private set { this[TestRelativeErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingNormalizedMeanSquaredError {
      get { return ((DoubleArray)this[TrainingNormalizedMeanSquaredErrorResultName].Value).ToArray(); }
      private set { this[TrainingNormalizedMeanSquaredErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TestNormalizedMeanSquaredError {
      get { return ((DoubleArray)this[TestNormalizedMeanSquaredErrorResultName].Value).ToArray(); }
      private set { this[TestNormalizedMeanSquaredErrorResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingDirectionalSymmetry {
      get { return ((DoubleArray)this[TrainingDirectionalSymmetryResultName].Value).ToArray(); }
      private set { this[TrainingDirectionalSymmetryResultName].Value = new DoubleArray(value); }
    }
    public double[] TestDirectionalSymmetry {
      get { return ((DoubleArray)this[TestDirectionalSymmetryResultName].Value).ToArray(); }
      private set { this[TestDirectionalSymmetryResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingWeightedDirectionalSymmetry {
      get { return ((DoubleArray)this[TrainingWeightedDirectionalSymmetryResultName].Value).ToArray(); }
      private set { this[TrainingWeightedDirectionalSymmetryResultName].Value = new DoubleArray(value); }
    }
    public double[] TestWeightedDirectionalSymmetry {
      get { return ((DoubleArray)this[TestWeightedDirectionalSymmetryResultName].Value).ToArray(); }
      private set { this[TestWeightedDirectionalSymmetryResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingTheilsUStatisticLast {
      get { return ((DoubleArray)this[TrainingTheilsUStatisticLastResultName].Value).ToArray(); }
      private set { this[TrainingTheilsUStatisticLastResultName].Value = new DoubleArray(value); }
    }
    public double[] TestTheilsUStatisticLast {
      get { return ((DoubleArray)this[TestTheilsUStatisticLastResultName].Value).ToArray(); }
      private set { this[TestTheilsUStatisticLastResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingTheilsUStatisticMean {
      get { return ((DoubleArray)this[TrainingTheilsUStatisticMeanResultName].Value).ToArray(); }
      private set { this[TrainingTheilsUStatisticMeanResultName].Value = new DoubleArray(value); }
    }
    public double[] TestTheilsUStatisticMean {
      get { return ((DoubleArray)this[TestTheilsUStatisticMeanResultName].Value).ToArray(); }
      private set { this[TestTheilsUStatisticMeanResultName].Value = new DoubleArray(value); }
    }
    public double[] TrainingTheilsUStatisticMovingAverage {
      get { return ((DoubleArray)this[TrainingTheilsUStatisticMAResultName].Value).ToArray(); }
      private set { this[TrainingTheilsUStatisticMAResultName].Value = new DoubleArray(value); }
    }
    public double[] TestTheilsUStatisticMovingAverage {
      get { return ((DoubleArray)this[TestTheilsUStatisticMAResultName].Value).ToArray(); }
      private set { this[TestTheilsUStatisticMAResultName].Value = new DoubleArray(value); }
    }
    #endregion

    [StorableConstructor]
    protected TimeSeriesPrognosisSolutionBase(bool deserializing) : base(deserializing) { }
    protected TimeSeriesPrognosisSolutionBase(TimeSeriesPrognosisSolutionBase original, Cloner cloner)
      : base(original, cloner) {
      this.horizon = original.horizon;
    }
    protected TimeSeriesPrognosisSolutionBase(ITimeSeriesPrognosisModel model, ITimeSeriesPrognosisProblemData problemData)
      : base(model, problemData) {
      Add(new Result(TrainingMeanSquaredErrorResultName, "Mean of squared errors of the model on the training partition", new DoubleArray()));
      Add(new Result(TestMeanSquaredErrorResultName, "Mean of squared errors of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingMeanAbsoluteErrorResultName, "Mean of absolute errors of the model on the training partition", new DoubleArray()));
      Add(new Result(TestMeanAbsoluteErrorResultName, "Mean of absolute errors of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingSquaredCorrelationResultName, "Squared Pearson's correlation coefficient of the model output and the actual values on the training partition", new DoubleArray()));
      Add(new Result(TestSquaredCorrelationResultName, "Squared Pearson's correlation coefficient of the model output and the actual values on the test partition", new DoubleArray()));
      Add(new Result(TrainingRelativeErrorResultName, "Average of the relative errors of the model output and the actual values on the training partition", new DoubleArray()));
      Add(new Result(TestRelativeErrorResultName, "Average of the relative errors of the model output and the actual values on the test partition", new DoubleArray()));
      Add(new Result(TrainingNormalizedMeanSquaredErrorResultName, "Normalized mean of squared errors of the model on the training partition", new DoubleArray()));
      Add(new Result(TestNormalizedMeanSquaredErrorResultName, "Normalized mean of squared errors of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingDirectionalSymmetryResultName, "The average directional symmetry of the forecasts of the model on the training partition", new DoubleArray()));
      Add(new Result(TestDirectionalSymmetryResultName, "The average directional symmetry of the forecasts of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingWeightedDirectionalSymmetryResultName, "The average weighted directional symmetry of the forecasts of the model on the training partition", new DoubleArray()));
      Add(new Result(TestWeightedDirectionalSymmetryResultName, "The average weighted directional symmetry of the forecasts of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingTheilsUStatisticLastResultName, "The average Theil's U statistic (reference: previous value) of the forecasts of the model on the training partition", new DoubleArray()));
      Add(new Result(TestTheilsUStatisticLastResultName, "The average Theil's U statistic (reference: previous value) of the forecasts of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingTheilsUStatisticMeanResultName, "The average Theil's U statistic (reference: mean value) of the forecasts of the model on the training partition", new DoubleArray()));
      Add(new Result(TestTheilsUStatisticMeanResultName, "The average Theil's U statistic (reference: mean value) of the forecasts of the model on the test partition", new DoubleArray()));
      Add(new Result(TrainingTheilsUStatisticMAResultName, "The average Theil's U statistic (reference: moving average) of the forecasts of the model on the training partition", new DoubleArray()));
      Add(new Result(TestTheilsUStatisticMAResultName, "The average Theil's U statistic (reference: moving average) of the forecasts of the model on the test partition", new DoubleArray()));
      horizon = 1;
    }

    [StorableHook(HookType.AfterDeserialization)]
    private void AfterDeserialization() {
      if (horizon == 0) horizon = 1;
    }

    protected void CalculateResults() {
      OnlineCalculatorError errorState;
      string[] targetVariables = ProblemData.TargetVariables.ToArray();
      /*
      double[] estimatedTrainingValues = PrognosedTrainingValues.ToArray(); // cache values
      double[] originalTrainingValues = ProblemData.Dataset.GetDoubleValues(ProblemData.TargetVariable, ProblemData.TrainingIndizes).ToArray();
      double[] estimatedTestValues = PrognosedTestValues.ToArray(); // cache values
      double[] originalTestValues = ProblemData.Dataset.GetDoubleValues(ProblemData.TargetVariable, ProblemData.TestIndizes).ToArray();

      double trainingMse = OnlineMeanSquaredErrorCalculator.Calculate(originalTrainingValues, estimatedTrainingValues, out errorState);
      TrainingMeanSquaredError = errorState == OnlineCalculatorError.None ? trainingMse : double.NaN;
      double testMse = OnlineMeanSquaredErrorCalculator.Calculate(originalTestValues, estimatedTestValues, out errorState);
      TestMeanSquaredError = errorState == OnlineCalculatorError.None ? testMse : double.NaN;

      double trainingMae = OnlineMeanAbsoluteErrorCalculator.Calculate(originalTrainingValues, estimatedTrainingValues, out errorState);
      TrainingMeanAbsoluteError = errorState == OnlineCalculatorError.None ? trainingMae : double.NaN;
      double testMae = OnlineMeanAbsoluteErrorCalculator.Calculate(originalTestValues, estimatedTestValues, out errorState);
      TestMeanAbsoluteError = errorState == OnlineCalculatorError.None ? testMae : double.NaN;

      double trainingR2 = OnlinePearsonsRSquaredCalculator.Calculate(originalTrainingValues, estimatedTrainingValues, out errorState);
      TrainingRSquared = errorState == OnlineCalculatorError.None ? trainingR2 : double.NaN;
      double testR2 = OnlinePearsonsRSquaredCalculator.Calculate(originalTestValues, estimatedTestValues, out errorState);
      TestRSquared = errorState == OnlineCalculatorError.None ? testR2 : double.NaN;

      double trainingRelError = OnlineMeanAbsolutePercentageErrorCalculator.Calculate(originalTrainingValues, estimatedTrainingValues, out errorState);
      TrainingRelativeError = errorState == OnlineCalculatorError.None ? trainingRelError : double.NaN;
      double testRelError = OnlineMeanAbsolutePercentageErrorCalculator.Calculate(originalTestValues, estimatedTestValues, out errorState);
      TestRelativeError = errorState == OnlineCalculatorError.None ? testRelError : double.NaN;

      double trainingNmse = OnlineNormalizedMeanSquaredErrorCalculator.Calculate(originalTrainingValues, estimatedTrainingValues, out errorState);
      TrainingNormalizedMeanSquaredError = errorState == OnlineCalculatorError.None ? trainingNmse : double.NaN;
      double testNmse = OnlineNormalizedMeanSquaredErrorCalculator.Calculate(originalTestValues, estimatedTestValues, out errorState);
      TestNormalizedMeanSquaredError = errorState == OnlineCalculatorError.None ? testNmse : double.NaN;
              */

      //double[] trainingDs = new double[targetVariables.Length];
      //double[] testDs = new double[targetVariables.Length];

      //double[] trainingWds = new double[targetVariables.Length];
      //double[] testWds = new double[targetVariables.Length];

      //double[] trainingTheilsU = new double[targetVariables.Length];
      //double[] testTheilsU = new double[targetVariables.Length];

      var trainingDsCalculators = new OnlineDirectionalSymmetryCalculator[targetVariables.Length];
      var testDsCalculators = new OnlineDirectionalSymmetryCalculator[targetVariables.Length];
      var trainingWdsCalculators = new OnlineWeightedDirectionalSymmetryCalculator[targetVariables.Length];
      var testWdsCalculators = new OnlineWeightedDirectionalSymmetryCalculator[targetVariables.Length];
      var trainingTheilsULastCalculators = new OnlineTheilsUStatisticCalculator[targetVariables.Length];
      var testTheilsULastCalculators = new OnlineTheilsUStatisticCalculator[targetVariables.Length];
      var trainingTheilsUMeanCalculators = new OnlineTheilsUStatisticCalculator[targetVariables.Length];
      var testTheilsUMeanCalculators = new OnlineTheilsUStatisticCalculator[targetVariables.Length];
      var trainingTheilsUMovingAverageCalculators = new OnlineTheilsUStatisticCalculator[targetVariables.Length];
      var testTheilsUMovingAverageCalculators = new OnlineTheilsUStatisticCalculator[targetVariables.Length];
      for (int i = 0; i < targetVariables.Length; i++) {
        trainingDsCalculators[i] = new OnlineDirectionalSymmetryCalculator();
        testDsCalculators[i] = new OnlineDirectionalSymmetryCalculator();
        trainingWdsCalculators[i] = new OnlineWeightedDirectionalSymmetryCalculator();
        testWdsCalculators[i] = new OnlineWeightedDirectionalSymmetryCalculator();
        trainingTheilsULastCalculators[i] = new OnlineTheilsUStatisticCalculator();
        testTheilsULastCalculators[i] = new OnlineTheilsUStatisticCalculator();
        trainingTheilsUMeanCalculators[i] = new OnlineTheilsUStatisticCalculator();
        testTheilsUMeanCalculators[i] = new OnlineTheilsUStatisticCalculator();
        trainingTheilsUMovingAverageCalculators[i] = new OnlineTheilsUStatisticCalculator();
        testTheilsUMovingAverageCalculators[i] = new OnlineTheilsUStatisticCalculator();
      }

      var allPrognosedTrainingValues = GetPrognosedValues(ProblemData.TrainingIndizes, horizon).GetEnumerator();
      foreach (var row in ProblemData.TrainingIndizes) {
        if (row + horizon < ProblemData.Dataset.Rows) {
          allPrognosedTrainingValues.MoveNext();
          var prognosedTrainingValues = allPrognosedTrainingValues.Current.SelectMany(x => x).ToArray();
          for (int t = 0; t < targetVariables.Length; t++) {
            var actualContinuation = ProblemData.Dataset.GetDoubleValues(targetVariables[t],
                                                                         Enumerable.Range(row, horizon));
            int maWindow = 10 * horizon;
            var movingAverageContinuation = from h in Enumerable.Range(0, horizon)
                                            select (from r in Enumerable.Range(row + h - maWindow, maWindow - h)
                                                    where r > 0
                                                    select ProblemData.Dataset.GetDoubleValue(targetVariables[t], r)
                                                   ).Average();
            var mean = ProblemData.Dataset.GetDoubleValues(targetVariables[t], ProblemData.TrainingIndizes).Median();
            double startValue = ProblemData.Dataset.GetDoubleValue(targetVariables[t], row - 1);
            var prognosedContinuation = prognosedTrainingValues.Skip(t).TakeEvery(targetVariables.Length);
            trainingDsCalculators[t].Add(startValue, actualContinuation, prognosedContinuation);
            trainingWdsCalculators[t].Add(startValue, actualContinuation, prognosedContinuation);
            trainingTheilsULastCalculators[t].Add(startValue, actualContinuation, prognosedContinuation);
            trainingTheilsUMeanCalculators[t].Add(startValue, actualContinuation.Select(x => mean), actualContinuation, prognosedContinuation);
            trainingTheilsUMovingAverageCalculators[t].Add(startValue, movingAverageContinuation, actualContinuation, prognosedContinuation);
          }
        }
      }
      var allPrognosedTestValues = GetPrognosedValues(ProblemData.TestIndizes, horizon).GetEnumerator();
      foreach (var row in ProblemData.TestIndizes) {
        if (row + horizon < ProblemData.Dataset.Rows) {
          allPrognosedTestValues.MoveNext();
          var prognosedTestValues = allPrognosedTestValues.Current.SelectMany(x => x).ToArray();
          for (int t = 0; t < targetVariables.Length; t++) {
            var actualContinuation = ProblemData.Dataset.GetDoubleValues(targetVariables[t],
                                                                         Enumerable.Range(row, horizon));
            int maWindow = 10 * horizon;
            var movingAverageContinuation = from h in Enumerable.Range(0, horizon)
                                            select (from r in Enumerable.Range(row + h - maWindow, maWindow - h)
                                                    where r > 0
                                                    select ProblemData.Dataset.GetDoubleValue(targetVariables[t], r)
                                                   ).Average();
            // mean must be calculated on the training set!
            var mean = ProblemData.Dataset.GetDoubleValues(targetVariables[t], ProblemData.TrainingIndizes).Median();
            double startValue = ProblemData.Dataset.GetDoubleValue(targetVariables[t], row - 1);
            var prognosedContinuation = prognosedTestValues.Skip(t).TakeEvery(targetVariables.Length);
            testDsCalculators[t].Add(startValue, actualContinuation, prognosedContinuation);
            testWdsCalculators[t].Add(startValue, actualContinuation, prognosedContinuation);
            testTheilsULastCalculators[t].Add(startValue, actualContinuation, prognosedContinuation);
            testTheilsUMeanCalculators[t].Add(startValue, actualContinuation.Select(x => mean), actualContinuation, prognosedContinuation);
            testTheilsUMovingAverageCalculators[t].Add(startValue, movingAverageContinuation, actualContinuation, prognosedContinuation);
          }
        }
      }

      TrainingDirectionalSymmetry = trainingDsCalculators.Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : 0.0)
        .ToArray();
      TestDirectionalSymmetry = testDsCalculators.Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : 0.0)
        .ToArray();
      TrainingWeightedDirectionalSymmetry = trainingWdsCalculators.Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TestWeightedDirectionalSymmetry = testWdsCalculators.Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TrainingTheilsUStatisticLast = trainingTheilsULastCalculators
        .Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TestTheilsUStatisticLast = testTheilsULastCalculators
        .Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TrainingTheilsUStatisticMean = trainingTheilsUMeanCalculators
        .Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TestTheilsUStatisticMean = testTheilsUMeanCalculators
        .Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TrainingTheilsUStatisticMovingAverage = trainingTheilsUMovingAverageCalculators
        .Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
      TestTheilsUStatisticMovingAverage = testTheilsUMovingAverageCalculators
        .Select(c => c.ErrorState == OnlineCalculatorError.None ? c.Value : double.PositiveInfinity)
        .ToArray();
    }
  }
}