source: branches/HeuristicLab.TimeSeries/HeuristicLab.Problems.DataAnalysis.Symbolic.TimeSeriesPrognosis/3.4/SingleObjective/SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator.cs @ 7154

#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;
using System.Collections.Generic;
using System.Drawing.Printing;
using System.Linq;
using HeuristicLab.Common;
using HeuristicLab.Core;
using HeuristicLab.Data;
using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;

namespace HeuristicLab.Problems.DataAnalysis.Symbolic.TimeSeriesPrognosis {
  [Item("Mean squared error Evaluator", "Calculates the mean squared error of a symbolic time-series prognosis solution.")]
  [StorableClass]
  public class SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator : SymbolicTimeSeriesPrognosisSingleObjectiveEvaluator {
    [StorableConstructor]
    protected SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator(bool deserializing) : base(deserializing) { }
    protected SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator(SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator original, Cloner cloner)
      : base(original, cloner) {
    }
    public override IDeepCloneable Clone(Cloner cloner) {
      return new SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator(this, cloner);
    }

    public SymbolicTimeSeriesPrognosisSingleObjectiveMeanSquaredErrorEvaluator() : base() { }

    public override bool Maximization { get { return false; } }

    public override IOperation Apply() {
      var solution = SymbolicExpressionTreeParameter.ActualValue;
      IEnumerable<int> rows = GenerateRowsToEvaluate();

      double quality = Calculate(SymbolicDataAnalysisTreeInterpreterParameter.ActualValue,
        solution,
        EstimationLimitsParameter.ActualValue.Lower, EstimationLimitsParameter.ActualValue.Upper,
        ProblemDataParameter.ActualValue,
        rows, HorizonParameter.ActualValue.Value);
      QualityParameter.ActualValue = new DoubleValue(quality);

      return base.Apply();
    }

    public static double Calculate(ISymbolicTimeSeriesPrognosisExpressionTreeInterpreter interpreter, ISymbolicExpressionTree solution, double lowerEstimationLimit, double upperEstimationLimit, ITimeSeriesPrognosisProblemData problemData, IEnumerable<int> rows, int horizon) {
      double[] alpha;
      double[] beta;
      DetermineScalingFactors(solution, problemData, interpreter, rows, out alpha, out beta);
      var scaledSolution = Scale(solution, alpha, beta);
      string[] targetVariables = problemData.TargetVariables.ToArray();
      var meanSquaredErrorCalculators = Enumerable.Range(0, problemData.TargetVariables.Count())
        .Select(i => new OnlineMeanSquaredErrorCalculator()).ToArray();

      var allContinuationsEnumerator = interpreter.GetSymbolicExpressionTreeValues(scaledSolution, problemData.Dataset,
                                                                                  targetVariables,
                                                                                  rows, horizon).GetEnumerator();
      allContinuationsEnumerator.MoveNext();
      // foreach row
      foreach (var row in rows) {
        // foreach horizon
        for (int h = 0; h < horizon; h++) {
          // foreach component
          for (int i = 0; i < meanSquaredErrorCalculators.Length; i++) {
            double e = Math.Min(upperEstimationLimit, Math.Max(lowerEstimationLimit, allContinuationsEnumerator.Current));
            meanSquaredErrorCalculators[i].Add(problemData.Dataset.GetDoubleValue(targetVariables[i], row + h), e);
            if (meanSquaredErrorCalculators[i].ErrorState == OnlineCalculatorError.InvalidValueAdded)
              return double.MaxValue;
            allContinuationsEnumerator.MoveNext();
          }
        }
      }
      var meanCalculator = new OnlineMeanAndVarianceCalculator();
      foreach (var calc in meanSquaredErrorCalculators) {
        if (calc.ErrorState != OnlineCalculatorError.None) return double.MaxValue;
        meanCalculator.Add(calc.MeanSquaredError);
      }
      //int i = 0;
      //foreach (var targetVariable in problemData.TargetVariables) {
      //  var predictedContinuations = allPredictedContinuations.Select(v => v.ElementAt(i));
      //  for (int h = 0; h < horizon; h++) {
      //    OnlineCalculatorError errorState;
      //    meanCalculator.Add(OnlineMeanSquaredErrorCalculator.Calculate(predictedContinuations
      //                                                                    .Select(x => x.ElementAt(h))
      //                                                                    .LimitToRange(lowerEstimationLimit,
      //                                                                                  upperEstimationLimit),
      //                                                                  actualContinuations.Select(x => x.ElementAt(h)),
      //                                                                  out errorState));
      //    if (errorState != OnlineCalculatorError.None) return double.NaN;
      //  }
      //}
      return meanCalculator.MeanErrorState == OnlineCalculatorError.None ? meanCalculator.Mean : double.MaxValue;
    }

    private static ISymbolicExpressionTree Scale(ISymbolicExpressionTree solution, double[] alpha, double[] beta) {
      var clone = (ISymbolicExpressionTree)solution.Clone();
      int n = alpha.Length;
      for (int i = 0; i < n; i++) {
        var parent = clone.Root.GetSubtree(0);
        var rpb = clone.Root.GetSubtree(0).GetSubtree(i);
        var scaledRpb = MakeSum(
          MakeProduct(rpb,
            MakeConstant(beta[i], clone.Root.Grammar), clone.Root.Grammar),
            MakeConstant(alpha[i], clone.Root.Grammar), clone.Root.Grammar);
        parent.RemoveSubtree(i);
        parent.InsertSubtree(i, scaledRpb);
      }
      return clone;
    }

    private static ISymbolicExpressionTreeNode MakeSum(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b, ISymbolicExpressionTreeGrammar grammar) {
      var sum = grammar.Symbols.Where(s => s is Addition).First().CreateTreeNode();
      sum.AddSubtree(a);
      sum.AddSubtree(b);
      return sum;
    }

    private static ISymbolicExpressionTreeNode MakeProduct(ISymbolicExpressionTreeNode a, ISymbolicExpressionTreeNode b, ISymbolicExpressionTreeGrammar grammar) {
      var prod = grammar.Symbols.Where(s => s is Multiplication).First().CreateTreeNode();
      prod.AddSubtree(a);
      prod.AddSubtree(b);
      return prod;
    }

    private static ISymbolicExpressionTreeNode MakeConstant(double c, ISymbolicExpressionTreeGrammar grammar) {
      var node = (ConstantTreeNode)grammar.Symbols.Where(s => s is Constant).First().CreateTreeNode();
      node.Value = c;
      return node;
    }

    private static void DetermineScalingFactors(ISymbolicExpressionTree solution, ITimeSeriesPrognosisProblemData problemData, ISymbolicTimeSeriesPrognosisExpressionTreeInterpreter interpreter, IEnumerable<int> rows, out double[] alpha, out double[] beta) {
      string[] targetVariables = problemData.TargetVariables.ToArray();
      int nComponents = targetVariables.Length;
      alpha = new double[nComponents];
      beta = new double[nComponents];
      var oneStepPredictionsEnumerator = interpreter.GetSymbolicExpressionTreeValues(solution, problemData.Dataset, targetVariables, rows).GetEnumerator();
      var scalingParameterCalculators =
        Enumerable.Repeat(0, nComponents).Select(x => new OnlineLinearScalingParameterCalculator()).ToArray();
      var targetValues = problemData.Dataset.GetVectorEnumerable(targetVariables, rows);
      var targetValueEnumerator = targetValues.GetEnumerator();

      var more = oneStepPredictionsEnumerator.MoveNext() & targetValueEnumerator.MoveNext();
      while (more) {
        for (int i = 0; i < nComponents; i++) {
          scalingParameterCalculators[i].Add(oneStepPredictionsEnumerator.Current, targetValueEnumerator.Current);
          more = oneStepPredictionsEnumerator.MoveNext() & targetValueEnumerator.MoveNext();
        }
      }

      for (int i = 0; i < nComponents; i++) {
        if (scalingParameterCalculators[i].ErrorState == OnlineCalculatorError.None) {
          alpha[i] = scalingParameterCalculators[i].Alpha;
          beta[i] = scalingParameterCalculators[i].Beta;
        } else {
          alpha[i] = 0.0;
          beta[i] = 1.0;
        }
      }
    }

    public override double Evaluate(IExecutionContext context, ISymbolicExpressionTree tree, ITimeSeriesPrognosisProblemData problemData, IEnumerable<int> rows) {
      SymbolicDataAnalysisTreeInterpreterParameter.ExecutionContext = context;
      EstimationLimitsParameter.ExecutionContext = context;
      HorizonParameter.ExecutionContext = context;

      double mse = Calculate(SymbolicDataAnalysisTreeInterpreterParameter.ActualValue, tree, EstimationLimitsParameter.ActualValue.Lower, EstimationLimitsParameter.ActualValue.Upper, problemData, rows, HorizonParameter.ActualValue.Value);

      HorizonParameter.ExecutionContext = null;
      SymbolicDataAnalysisTreeInterpreterParameter.ExecutionContext = null;
      EstimationLimitsParameter.ExecutionContext = null;

      return mse;
    }
  }
}