1  #region License Information


2  /* HeuristicLab


3  * Copyright (C) 20022018 Heuristic and Evolutionary Algorithms Laboratory (HEAL)


4  *


5  * This file is part of HeuristicLab.


6  *


7  * HeuristicLab is free software: you can redistribute it and/or modify


8  * it under the terms of the GNU General Public License as published by


9  * the Free Software Foundation, either version 3 of the License, or


10  * (at your option) any later version.


11  *


12  * HeuristicLab is distributed in the hope that it will be useful,


13  * but WITHOUT ANY WARRANTY; without even the implied warranty of


14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the


15  * GNU General Public License for more details.


16  *


17  * You should have received a copy of the GNU General Public License


18  * along with HeuristicLab. If not, see <http://www.gnu.org/licenses/>.


19  */


20  #endregion


21 


22  using System;


23  using System.Collections.Generic;


24  using System.Linq;


25  using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding;


26 


27  namespace HeuristicLab.Problems.DataAnalysis.Symbolic.ConstantsOptimization {


28  public class LMConstantsOptimizer {


29 


30  private LMConstantsOptimizer() { }


31 


32  /// <summary>


33  /// Method to determine whether the numeric constants of the tree can be optimized. This depends primarily on the symbols occuring in the tree.


34  /// </summary>


35  /// <param name="tree">The tree that should be analyzed</param>


36  /// <returns>A flag indicating whether the numeric constants of the tree can be optimized</returns>


37  public static bool CanOptimizeConstants(ISymbolicExpressionTree tree) {


38  return AutoDiffConverter.IsCompatible(tree);


39  }


40 


41  /// <summary>


42  /// Optimizes the numeric constants in a symbolic expression tree in place.


43  /// </summary>


44  /// <param name="tree">The tree for which the constants should be optimized</param>


45  /// <param name="dataset">The dataset containing the data.</param>


46  /// <param name="targetVariable">The target variable name.</param>


47  /// <param name="rows">The rows for which the data should be extracted.</param>


48  /// <param name="applyLinearScaling">A flag to determine whether linear scaling should be applied during the optimization</param>


49  /// <param name="maxIterations">The maximum number of iterations of the LevenbergMarquard algorithm.</param>


50  /// <returns></returns>


51  public static double OptimizeConstants(ISymbolicExpressionTree tree,


52  IDataset dataset, string targetVariable, IEnumerable<int> rows,


53  bool applyLinearScaling, int maxIterations = 10) {


54  if (tree == null) throw new ArgumentNullException("tree");


55  if (dataset == null) throw new ArgumentNullException("dataset");


56  if (!dataset.ContainsVariable(targetVariable)) throw new ArgumentException("The dataset does not contain the provided target variable.");


57 


58  var allVariables = Util.ExtractVariables(tree);


59  var numericNodes = Util.ExtractNumericNodes(tree);


60 


61  AutoDiff.IParametricCompiledTerm term;


62  if (!AutoDiffConverter.TryConvertToAutoDiff(tree, applyLinearScaling, numericNodes, allVariables, out term))


63  throw new NotSupportedException("Could not convert symbolic expression tree to an AutoDiff term due to not supported symbols used in the tree.");


64 


65  //Variables of the symbolic expression tree correspond to parameters in the term


66  //Hence if no parameters are present no variables occur in the tree and the R² = 0


67  if (term.Parameters.Count == 0) return 0.0;


68 


69  var initialConstants = Util.ExtractConstants(numericNodes, applyLinearScaling);


70  double[] constants;


71  double[,] x = Util.ExtractData(dataset, allVariables, rows);


72  double[] y = dataset.GetDoubleValues(targetVariable, rows).ToArray();


73 


74  var result = OptimizeConstants(term, initialConstants, x, y, maxIterations, out constants);


75  if (result > 0.0 && constants.Length != 0)


76  Util.UpdateConstants(numericNodes, constants);


77 


78  return result;


79  }


80 


81  /// <summary>


82  /// Optimizes the numeric coefficents of an AutoDiff Term using the LevenbergMarquard algorithm.


83  /// </summary>


84  /// <param name="term">The AutoDiff term for which the numeric coefficients should be optimized.</param>


85  /// <param name="initialConstants">The starting values for the numeric coefficients.</param>


86  /// <param name="x">The input data for the optimization.</param>


87  /// <param name="y">The target values for the optimization.</param>


88  /// <param name="maxIterations">The maximum number of iterations of the LevenbergMarquard</param>


89  /// <param name="constants">The opitmized constants.</param>


90  /// <param name="LM_IterationCallback">An optional callback for detailed analysis that is called in each algorithm iteration.</param>


91  /// <returns>The R² of the term evaluated on the input data x and the target data y using the optimized constants</returns>


92  public static double OptimizeConstants(AutoDiff.IParametricCompiledTerm term, double[] initialConstants, double[,] x, double[] y,


93  int maxIterations, out double[] constants, Action<double[], double, object> LM_IterationCallback = null) {


94 


95  if (term.Parameters.Count == 0) {


96  constants = new double[0];


97  return 0.0;


98  }


99 


100  var optimizedConstants = (double[])initialConstants.Clone();


101  int numberOfRows = x.GetLength(0);


102  int numberOfColumns = x.GetLength(1);


103  int numberOfConstants = optimizedConstants.Length;


104 


105  alglib.lsfitstate state;


106  alglib.lsfitreport rep;


107  alglib.ndimensional_rep xrep = (p, f, obj) => LM_IterationCallback(p, f, obj);


108  int retVal;


109 


110  try {


111  alglib.lsfitcreatefg(x, y, optimizedConstants, numberOfRows, numberOfColumns, numberOfConstants, cheapfg: false, state: out state);


112  alglib.lsfitsetcond(state, 0.0, 0.0, maxIterations);


113  alglib.lsfitsetxrep(state, LM_IterationCallback != null);


114  alglib.lsfitfit(state, Evaluate, EvaluateGradient, xrep, term);


115  alglib.lsfitresults(state, out retVal, out optimizedConstants, out rep);


116  } catch (ArithmeticException) {


117  constants = new double[0];


118  return double.NaN;


119  } catch (alglib.alglibexception) {


120  constants = new double[0];


121  return double.NaN;


122  }


123 


124  constants = optimizedConstants;


125  return rep.r2;


126  }


127 


128 


129  private static void Evaluate(double[] c, double[] x, ref double fx, object o) {


130  AutoDiff.IParametricCompiledTerm term = (AutoDiff.IParametricCompiledTerm)o;


131  fx = term.Evaluate(c, x);


132  }


133 


134  private static void EvaluateGradient(double[] c, double[] x, ref double fx, double[] grad, object o) {


135  AutoDiff.IParametricCompiledTerm term = (AutoDiff.IParametricCompiledTerm)o;


136  Tuple<double[], double> result = term.Differentiate(c, x);


137  fx = result.Item2;


138  Array.Copy(result.Item1, grad, grad.Length);


139  }


140  }


141  }

