source: branches/2886_SymRegGrammarEnumeration/HeuristicLab.Algorithms.DataAnalysis.SymRegGrammarEnumeration/Analysis/BestSolutionAnalyzer.cs @ 16090

#region License Information
/* HeuristicLab
 * Copyright (C) 2002-2018 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.Diagnostics;
using System.Linq;
using HeuristicLab.Analysis;
using HeuristicLab.Common;
using HeuristicLab.Core;
using HeuristicLab.Data;
using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;
using HeuristicLab.Optimization;
using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
using HeuristicLab.Problems.DataAnalysis.Symbolic;
using HeuristicLab.Problems.DataAnalysis.Symbolic.Regression;

namespace HeuristicLab.Algorithms.DataAnalysis.SymRegGrammarEnumeration {
  [Item("Best Solution Analyzer", "Returns the characteristics of the best solution so far.")]
  [StorableClass]
  public class BestSolutionAnalyzer : Item, IGrammarEnumerationAnalyzer {
    public static readonly string BestTrainingQualityResultName = "Best R² (Training)";
    public static readonly string BestTestQualityResultName = "Best R² (Test)";
    public static readonly string BestTrainingModelResultName = "Best model (Training)";
    public static readonly string BestTrainingSolutionResultName = "Best solution (Training)";
    public static readonly string BestComplexityResultName = "Best solution complexity";
    public static readonly string BestSolutions = "Best solutions";
    public static readonly string ParetoFrontResultName = "Pareto Front";
    public static readonly string ParetoFrontAnalysisResultName = "Pareto Front Analysis";
    public static readonly string ParetoFrontSolutionsResultName = "Pareto Front Solutions";

    private static readonly ISymbolicDataAnalysisExpressionTreeInterpreter expressionTreeLinearInterpreter = new SymbolicDataAnalysisExpressionTreeLinearInterpreter();

    public BestSolutionAnalyzer() { }

    [StorableConstructor]
    protected BestSolutionAnalyzer(bool deserializing) : base(deserializing) { }

    protected BestSolutionAnalyzer(BestSolutionAnalyzer original, Cloner cloner) : base(original, cloner) {
    }

    public override IDeepCloneable Clone(Cloner cloner) {
      return new BestSolutionAnalyzer(this, cloner);
    }

    public void Deregister(GrammarEnumerationAlgorithm algorithm) {
      algorithm.DistinctSentenceGenerated -= AlgorithmDistinctSentenceGenerated;
      algorithm.Stopped -= AlgorithmOnStopped;
    }

    public void Register(GrammarEnumerationAlgorithm algorithm) {
      algorithm.DistinctSentenceGenerated += AlgorithmDistinctSentenceGenerated;
      algorithm.Stopped += AlgorithmOnStopped;
    }

    private void AlgorithmOnStopped(object sender, EventArgs eventArgs) {
      var algorithm = (GrammarEnumerationAlgorithm)sender;

      IResult paretoFrontResult;
      if (algorithm.Results.TryGetValue(ParetoFrontAnalysisResultName, out paretoFrontResult)) {
        var plot = (ScatterPlot)paretoFrontResult.Value;

        var solutions = plot.Rows.First().Points.Select(p => (ISymbolicRegressionSolution)p.Tag);

        algorithm.Results.AddOrUpdateResult(ParetoFrontSolutionsResultName, new ItemList<ISymbolicRegressionSolution>(solutions));
      }
    }

    private void AlgorithmDistinctSentenceGenerated(object sender, PhraseAddedEventArgs args) {
      var algorithm = (GrammarEnumerationAlgorithm)sender;
      var sentence = args.NewPhrase;

      var results = algorithm.Results;
      var problemData = algorithm.Problem.ProblemData;

      SymbolicExpressionTree tree = algorithm.Grammar.ParseSymbolicExpressionTree(sentence);
      Debug.Assert(SymbolicRegressionConstantOptimizationEvaluator.CanOptimizeConstants(tree));

      double r2 = algorithm.Evaluator.Evaluate(problemData, tree);
      int rank = GetRank(sentence);

      // Store solution in pareto front
      if (IsParetoOptimal(algorithm, rank, r2)) {
        var model = new SymbolicRegressionModel(problemData.TargetVariable, tree, expressionTreeLinearInterpreter);
        model.Scale(problemData);
        var bestSolution = model.CreateRegressionSolution(problemData);

        AddToParetoFront(algorithm, rank, r2, bestSolution);

        // Store overall best solution
        double bestR2 = results.ContainsKey(BestTrainingQualityResultName)
          ? GetValue<double>(results[BestTrainingQualityResultName].Value)
          : 0.0;
        var bestComplexity = results.ContainsKey(BestComplexityResultName) ? GetValue<int>(results[BestComplexityResultName].Value) : int.MaxValue;
        var complexity = sentence.Complexity;

        if (algorithm.BestTrainingSentence == null || r2 > bestR2 || (r2.IsAlmost(bestR2) && complexity < bestComplexity)) {
          algorithm.BestTrainingSentence = sentence;

          results.AddOrUpdateResult(BestTrainingQualityResultName, new DoubleValue(bestSolution.TrainingRSquared));
          results.AddOrUpdateResult(BestTestQualityResultName, new DoubleValue(bestSolution.TestRSquared));
          results.AddOrUpdateResult(BestTrainingModelResultName, bestSolution.Model);
          results.AddOrUpdateResult(BestTrainingSolutionResultName, bestSolution);
          results.AddOrUpdateResult(BestComplexityResultName, new IntValue(complexity));

          // record best sentence quality & length
          DataTable dt;
          if (!results.ContainsKey(BestSolutions)) {
            var names = new[] { "Quality", "Length", "Complexity", "Timestamp" };
            dt = new DataTable();
            foreach (var name in names) {
              dt.Rows.Add(new DataRow(name) { VisualProperties = { StartIndexZero = true } });
            }
            results.AddOrUpdateResult(BestSolutions, dt);
          }
          dt = (DataTable)results[BestSolutions].Value;
          dt.Rows["Quality"].Values.Add(r2);
          dt.Rows["Length"].Values.Add((double)sentence.Count);
          dt.Rows["Complexity"].Values.Add(complexity);
          dt.Rows["Timestamp"].Values.Add(algorithm.ExecutionTime.TotalMilliseconds / 1000d);
        }
      }

      // stop the algorithm if the best quality was already achieved
      if (r2.IsAlmost(1d)) {
        algorithm.Stop();
      }
    }

    private T GetValue<T>(IItem value) where T : struct {
      var v = value as ValueTypeValue<T>;
      if (v == null)
        throw new ArgumentException(string.Format("Item is not of type {0}", typeof(ValueTypeValue<T>)));
      return v.Value;
    }

    private int GetRank(SymbolList s) {
      return s.Complexity;
    }

    private bool IsParetoOptimal(GrammarEnumerationAlgorithm algorithm, int currRank, double currQuality) {
      if (!algorithm.Results.ContainsKey(ParetoFrontResultName)) return true;

      ItemList<DoubleArray> paretoFront = (ItemList<DoubleArray>)algorithm.Results[ParetoFrontResultName].Value;

      int preceedingRankIndex = -1;
      int lastIndex = paretoFront.Count - 1;
      while (preceedingRankIndex < lastIndex) {
        if (preceedingRankIndex + 1 > currRank)
          break;
        preceedingRankIndex++;
      }

      return preceedingRankIndex < 0 || paretoFront[preceedingRankIndex][1] < currQuality;
    }

    private void AddToParetoFront(GrammarEnumerationAlgorithm algorithm, int currRank, double currQuality, ISymbolicRegressionSolution solution) {
      if (!algorithm.Results.ContainsKey(ParetoFrontResultName)) {
        algorithm.Results.Add(new Result(ParetoFrontResultName, new ItemList<DoubleArray>()));

        var scatterPlot = new ScatterPlot(ParetoFrontAnalysisResultName, ParetoFrontAnalysisResultName);
        algorithm.Results.Add(new Result(ParetoFrontAnalysisResultName, scatterPlot));
        scatterPlot.Rows.Add(new ScatterPlotDataRow());

        scatterPlot.VisualProperties.XAxisTitle = "Complexity";
        scatterPlot.VisualProperties.YAxisTitle = "R²";
        scatterPlot.Rows.First().VisualProperties.PointSize = 10;
      }

      ItemList<DoubleArray> paretoFront = (ItemList<DoubleArray>)algorithm.Results[ParetoFrontResultName].Value;
      ScatterPlotDataRow plot = ((ScatterPlot)algorithm.Results[ParetoFrontAnalysisResultName].Value).Rows.First();

      // Delete solutions with higher rank, which are now dominated.
      int i = 0;
      while (i < paretoFront.Count) {
        if (paretoFront[i][0] >= currRank) { // Go to current rank
          double quality = paretoFront[i][1];
          if (quality <= currQuality) { // If existing solution is worse, delete it
            RemovePoint(plot, currRank);
            paretoFront.RemoveAt(i);
          } else { // Otherwise stop, since following solutions can only be better
            break;
          }
        } else {
          i++;
        }
      }

      paretoFront.Insert(i, new DoubleArray(new double[] { currRank, currQuality }));
      plot.Points.Add(new Point2D<double>(currRank, currQuality, solution));
    }

    private void RemovePoint(ScatterPlotDataRow plot, double rank) {
      plot.Points.RemoveAll(p => p.X.IsAlmost(rank));
    }
  }
}