source: addons/HeuristicLab.FitnessLandscapeAnalysis/HeuristicLab.Analysis.FitnessLandscape/Evolvability/EvolvabilityAnalyzer.cs

using System;
using System.Collections.Generic;
using System.Drawing;
using System.Linq;
using HEAL.Attic;
using HeuristicLab.Common;
using HeuristicLab.Core;
using HeuristicLab.Data;
using HeuristicLab.Operators;
using HeuristicLab.Optimization;
using HeuristicLab.Optimization.Operators;
using HeuristicLab.Parameters;

namespace HeuristicLab.Analysis.FitnessLandscape.Evolvability {

  [Item("Evolvability Analyzer", "Compares two sets of scopes containing the original population and the evolved population generating several evolvability measures.")]
  [StorableType("6B05CC43-161A-45C6-9B08-A4CC2235A96B")]
  public class EvolvabilityAnalyzer : AlgorithmOperator, IAnalyzer {

    public bool EnabledByDefault
    {
      get { return false; }
    }

    #region Parameters
    public LookupParameter<DoubleValue> QualityParameter
    {
      get { return (LookupParameter<DoubleValue>)Parameters["Quality"]; }
    }
    public LookupParameter<IntValue> SampleSizeParameter
    {
      get { return (LookupParameter<IntValue>)Parameters["SampleSize"]; }
    }
    public LookupParameter<HeatMap> EvolvabilityPlotParameter
    {
      get { return (LookupParameter<HeatMap>)Parameters["EvolvabilityPlot"]; }
    }
    public LookupParameter<VariableCollection> ResultsParameter
    {
      get { return (LookupParameter<VariableCollection>)Parameters["Results"]; }
    }
    public LookupParameter<ScatterPlot> FitnessCloudParameter
    {
      get { return (LookupParameter<ScatterPlot>)Parameters["FitnessCloud"]; }
    }
    public LookupParameter<DataTable> EvolvabilityPortraitParameter
    {
      get { return (LookupParameter<DataTable>)Parameters["EvolvabilityPortrait"]; }
    }
    public ValueLookupParameter<IntValue> NrOfBinsParameter
    {
      get { return (ValueLookupParameter<IntValue>)Parameters["NrOfBins"]; }
    }
    public ScopeParameter ScopeParameter
    {
      get { return (ScopeParameter)Parameters["Scope"]; }
    }
    #endregion

    #region Constructors
    public EvolvabilityAnalyzer() {
      Parameters.Add(new LookupParameter<DoubleValue>("Quality", "The qualities of evolved and unevolved solution candidates"));
      Parameters.Add(new LookupParameter<IntValue>("SampleSize", "The number of samples or population size for evolvability analysis"));
      Parameters.Add(new LookupParameter<HeatMap>("EvolvabilityPlot", "Collected evolvability values"));
      Parameters.Add(new LookupParameter<VariableCollection>("Results", "The collection of all results of this algorithm"));
      Parameters.Add(new LookupParameter<ScatterPlot>("FitnessCloud", "A fitness cloud of the original vs. evolved sample."));
      Parameters.Add(new LookupParameter<DataTable>("EvolvabilityPortrait", "Evolvability Portrait as described in Smith et. al (2002). Evol Comput 10:1–34"));
      Parameters.Add(new ValueLookupParameter<IntValue>("NrOfBins", "Nr of bins for evolvability portrait analysis", new IntValue(25)));
      Parameters.Add(new ScopeParameter("Scope", "The current scope"));

      var resultsCollector = new ResultsCollector();
      resultsCollector.CollectedValues.Add(new LookupParameter<HeatMap>(EvolvabilityPlotParameter.Name));
      resultsCollector.CollectedValues.Add(new LookupParameter<ScatterPlot>(FitnessCloudParameter.Name));
      resultsCollector.CollectedValues.Add(new LookupParameter<DataTable>(EvolvabilityPortraitParameter.Name));

      OperatorGraph.InitialOperator = resultsCollector;
      resultsCollector.Successor = null;
    }

    [StorableConstructor]
    protected EvolvabilityAnalyzer(StorableConstructorFlag _) : base(_) { }

    protected EvolvabilityAnalyzer(EvolvabilityAnalyzer original, Cloner cloner)
      : base(original, cloner) {
    }
    #endregion

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

    public override IOperation Apply() {
      var points = GenerateQualityPairs(ScopeParameter.ActualValue);
      GetOrCreateEvolvabilityPlot().Add(points);
      GetOrCreateScatterPlot().Add(points);
      CreateEvolvabilityPortrait(ScopeParameter.ActualValue);
      return base.Apply();
    }

    private IEnumerable<PointF> GenerateQualityPairs(IScope globalScope) {
      foreach (IScope sampleScope in globalScope.SubScopes) {
        double originalQuality = ((DoubleValue)sampleScope.Variables[QualityParameter.ActualName].Value).Value;
        foreach (IScope evolvedScope in sampleScope.SubScopes) {
          double evolvedQuality = ((DoubleValue)evolvedScope.Variables[QualityParameter.ActualName].Value).Value;
          yield return new PointF((float)originalQuality, (float)evolvedQuality);
        }
      }
    }

    private HeatMap GetOrCreateEvolvabilityPlot() {
      HeatMap evolvabilityPlot = EvolvabilityPlotParameter.ActualValue as HeatMap;
      if (evolvabilityPlot == null) {
        evolvabilityPlot = new HeatMap();
        evolvabilityPlot.Name = "Evolvability";
        EvolvabilityPlotParameter.ActualValue = evolvabilityPlot;
      }
      return evolvabilityPlot;
    }

    private ScatterPlot GetOrCreateScatterPlot() {
      ScatterPlot scatterPlot = FitnessCloudParameter.ActualValue as ScatterPlot;
      if (scatterPlot == null) {
        scatterPlot = new ScatterPlot();
        scatterPlot.Name = "Fitness Cloud";
        scatterPlot.XAxisName = "Parent Fitness";
        scatterPlot.YAxisName = "Offspring Fitness";
        FitnessCloudParameter.ActualValue = scatterPlot;
      }
      return scatterPlot;
    }

    private DataTable GetOrCreateEvolvabilityPortrait() {
      DataTable evolvabilityPortrait = EvolvabilityPortraitParameter.ActualValue as DataTable;
      if (evolvabilityPortrait == null) {
        evolvabilityPortrait = new DataTable("Evolvability Portrait");
        evolvabilityPortrait.Rows.Add(new DataRow("E_a", "Probability of non-deleterious mutation"));
        evolvabilityPortrait.Rows.Add(new DataRow("E_b", "Average expected offspring fitness"));
        evolvabilityPortrait.Rows.Add(new DataRow("E_c", "Top offspring fitness"));
        evolvabilityPortrait.Rows.Add(new DataRow("E_d", "Bottom offspring fitness"));
        EvolvabilityPortraitParameter.ActualValue = evolvabilityPortrait;
      }
      return evolvabilityPortrait;
    }

    private double GetQuality(IScope scope) {
      return ((DoubleValue)scope.Variables[QualityParameter.ActualName].Value).Value;
    }

    private void CreateEvolvabilityPortrait(IScope globalScope) {
      DataTable table = GetOrCreateEvolvabilityPortrait();
      var e_a = table.Rows["E_a"].Values;
      var e_b = table.Rows["E_b"].Values;
      var e_c = table.Rows["E_c"].Values;
      var e_d = table.Rows["E_d"].Values;
      var qualities = globalScope.SubScopes.Select(s => GetQuality(s));
      double min = qualities.Min();
      double max = qualities.Max();
      double range = max - min;
      double stepSize = range / NrOfBinsParameter.ActualValue.Value;
      double threshold = min + stepSize;

      int n_not_worse = 0;
      int n = 0;
      double totalOffspringQuality = 0;
      double best = double.MinValue;
      double worst = double.MaxValue;
      foreach (IScope sampleScope in globalScope.SubScopes) {
        double originalQuality = GetQuality(sampleScope);
        foreach (IScope evolvedScope in sampleScope.SubScopes) {
          double evolvedQuality = GetQuality(evolvedScope);
          if (evolvedQuality >= originalQuality)
            n_not_worse++;
          n++;
          totalOffspringQuality += evolvedQuality;
          best = Math.Max(best, evolvedQuality);
          worst = Math.Min(worst, evolvedQuality);
        }
        if (originalQuality > threshold) {
          threshold += stepSize;
          e_a.Add(1.0 * n_not_worse / n);
          e_b.Add(totalOffspringQuality / n);
          e_c.Add(best);
          e_d.Add(worst);
          n_not_worse = 0;
          n = 0;
          totalOffspringQuality = 0;
          best = double.MinValue;
          worst = double.MaxValue;
        }
      }
    }
  }
}