using System; using HeuristicLab.Optimization; using HEAL.Attic; using HeuristicLab.Common; using System.Threading; using HeuristicLab.Core; using HeuristicLab.Data; using HeuristicLab.Parameters; using System.Linq; using HeuristicLab.Encodings.SymbolicExpressionTreeEncoding; using HeuristicLab.Analysis; namespace HeuristicLab.Problems.DataAnalysis.Symbolic.Regression { [StorableType("676B237C-DD9C-4F24-B64F-D44B0FA1F6A6")] [Creatable(CreatableAttribute.Categories.DataAnalysisRegression, Priority = 120)] [Item(Name = "ConstrainedNLS", Description = "Non-linear Regression with non-linear constraints")] public class ConstrainedNLS : BasicAlgorithm { public static readonly string IterationsParameterName = "Iterations"; public static readonly string SolverParameterName = "Solver"; public static readonly string ModelStructureParameterName = "Model structure"; public IFixedValueParameter IterationsParameter { get { return (IFixedValueParameter)Parameters[IterationsParameterName]; } } public IFixedValueParameter ModelStructureParameter { get { return (IFixedValueParameter)Parameters[ModelStructureParameterName]; } } public IConstrainedValueParameter SolverParameter { get { return (IConstrainedValueParameter)Parameters[SolverParameterName]; } } public IFixedValueParameter FuncToleranceRelParameter { get { return (IFixedValueParameter)Parameters["FuncToleranceRel"]; } } public IFixedValueParameter FuncToleranceAbsParameter { get { return (IFixedValueParameter)Parameters["FuncToleranceAbs"]; } } public IFixedValueParameter MaxTimeParameter { get { return (IFixedValueParameter)Parameters["MaxTime"]; } } public IFixedValueParameter CheckGradientParameter { get { return (IFixedValueParameter)Parameters["CheckGradient"]; } } public int Iterations { get { return IterationsParameter.Value.Value; } set { IterationsParameter.Value.Value = value; } } public StringValue Solver { get { return SolverParameter.Value; } set { throw new NotImplementedException(); } } public string ModelStructure { get { return ModelStructureParameter.Value.Value; } set { ModelStructureParameter.Value.Value = value; } } public bool CheckGradient { get { return CheckGradientParameter.Value.Value; } set { CheckGradientParameter.Value.Value = value; } } public double FuncToleranceRel { get { return FuncToleranceRelParameter.Value.Value; } set { FuncToleranceRelParameter.Value.Value = value; } } public double FuncToleranceAbs { get { return FuncToleranceAbsParameter.Value.Value; } set { FuncToleranceAbsParameter.Value.Value = value; } } public double MaxTime { get { return MaxTimeParameter.Value.Value; } set { MaxTimeParameter.Value.Value = value; } } public ConstrainedNLS() : base() { Problem = new RegressionProblem(); Parameters.Add(new FixedValueParameter(ModelStructureParameterName, "The function for which the parameters must be fit (only numeric constants are tuned).", new StringValue("1.0 * x*x + 0.0"))); Parameters.Add(new FixedValueParameter(IterationsParameterName, "Determines how many iterations should be calculated while optimizing the constant of a symbolic expression tree (0 indicates other or default stopping criterion).", new IntValue(10))); var validSolvers = new ItemSet(new[] { "MMA", "COBYLA", "CCSAQ", "ISRES", "DIRECT_G", "NLOPT_GN_DIRECT_L", "NLOPT_GN_DIRECT_L_RAND", "NLOPT_GN_ORIG_DIRECT", "NLOPT_GN_ORIG_DIRECT_L", "NLOPT_GD_STOGO", "NLOPT_GD_STOGO_RAND", "NLOPT_LD_LBFGS_NOCEDAL", "NLOPT_LD_LBFGS", "NLOPT_LN_PRAXIS", "NLOPT_LD_VAR1", "NLOPT_LD_VAR2", "NLOPT_LD_TNEWTON", "NLOPT_LD_TNEWTON_RESTART", "NLOPT_LD_TNEWTON_PRECOND", "NLOPT_LD_TNEWTON_PRECOND_RESTART", "NLOPT_GN_CRS2_LM", "NLOPT_GN_MLSL", "NLOPT_GD_MLSL", "NLOPT_GN_MLSL_LDS", "NLOPT_GD_MLSL_LDS", "NLOPT_LN_NEWUOA", "NLOPT_LN_NEWUOA_BOUND", "NLOPT_LN_NELDERMEAD", "NLOPT_LN_SBPLX", "NLOPT_LN_AUGLAG", "NLOPT_LD_AUGLAG", "NLOPT_LN_BOBYQA", "NLOPT_AUGLAG", "NLOPT_LD_SLSQP", "NLOPT_LD_CCSAQ", "NLOPT_GN_ESCH", "NLOPT_GN_AGS", }.Select(s => new StringValue(s).AsReadOnly())); Parameters.Add(new ConstrainedValueParameter(SolverParameterName, "The solver algorithm", validSolvers, validSolvers.First())); Parameters.Add(new FixedValueParameter("FuncToleranceRel", new DoubleValue(0))); Parameters.Add(new FixedValueParameter("FuncToleranceAbs", new DoubleValue(0))); Parameters.Add(new FixedValueParameter("MaxTime", new DoubleValue(10))); Parameters.Add(new FixedValueParameter("CheckGradient", "Flag to indicate whether the gradient should be checked using numeric approximation", new BoolValue(false))); CheckGradientParameter.Hidden = true; } public ConstrainedNLS(ConstrainedNLS original, Cloner cloner) : base(original, cloner) { } [StorableHook(HookType.AfterDeserialization)] public void AfterDeserializationHook() { if (!Parameters.ContainsKey("CheckGradient")) { Parameters.Add(new FixedValueParameter("CheckGradient", "Flag to indicate whether the gradient should be checked using numeric approximation", new BoolValue(false))); CheckGradientParameter.Hidden = true; } } [StorableConstructor] protected ConstrainedNLS(StorableConstructorFlag _) : base(_) { } public override bool SupportsPause => false; public override IDeepCloneable Clone(Cloner cloner) { return new ConstrainedNLS(this, cloner); } protected override void Run(CancellationToken cancellationToken) { var parser = new InfixExpressionParser(); var tree = parser.Parse(ModelStructure); var problem = (IRegressionProblem)Problem; #region prepare results var functionEvaluations = new IntValue(0); Results.AddOrUpdateResult("Evaluations", functionEvaluations); var bestError = new DoubleValue(double.MaxValue); var curError = new DoubleValue(double.MaxValue); Results.AddOrUpdateResult("Best error", bestError); Results.AddOrUpdateResult("Current error", curError); Results.AddOrUpdateResult("Tree", tree); var qualitiesTable = new DataTable("Qualities"); var curQualityRow = new DataRow("Current Quality"); var bestQualityRow = new DataRow("Best Quality"); qualitiesTable.Rows.Add(bestQualityRow); qualitiesTable.Rows.Add(curQualityRow); Results.AddOrUpdateResult("Qualities", qualitiesTable); var curConstraintValue = new DoubleValue(0); Results.AddOrUpdateResult("Current Constraint Value", curConstraintValue); var curConstraintIdx = new IntValue(0); Results.AddOrUpdateResult("Current Constraint Index", curConstraintIdx); var curConstraintRow = new DataRow("Constraint Value"); var constraintsTable = new DataTable("Constraints"); constraintsTable.Rows.Add(curConstraintRow); Results.AddOrUpdateResult("Constraints", constraintsTable); #endregion var state = new ConstrainedNLSInternal(Solver.Value, tree, Iterations, problem.ProblemData, FuncToleranceRel, FuncToleranceAbs, MaxTime); if (CheckGradient) state.CheckGradient = true; int idx = 0; var formatter = new InfixExpressionFormatter(); var constraintDescriptions = state.ConstraintDescriptions.ToArray(); foreach(var constraintTree in state.constraintTrees) { // HACK to remove parameter nodes which occurr multiple times var reparsedTree = parser.Parse(formatter.Format(constraintTree)); Results.AddOrUpdateResult($"{constraintDescriptions[idx++]}", reparsedTree); } // we use a listener model here to get state from the solver state.FunctionEvaluated += State_FunctionEvaluated; state.ConstraintEvaluated += State_ConstraintEvaluated; state.Optimize(); bestError.Value = state.BestError; curQualityRow.Values.Add(state.CurError); bestQualityRow.Values.Add(bestError.Value); Results.AddOrUpdateResult("Best solution", CreateSolution((ISymbolicExpressionTree)state.BestTree.Clone(), problem.ProblemData)); var bestConstraintValues = new DoubleArray(state.BestConstraintValues); bestConstraintValues.ElementNames = constraintDescriptions; Results.AddOrUpdateResult("Best solution constraint values", bestConstraintValues); // local function void State_FunctionEvaluated() { if (cancellationToken.IsCancellationRequested) state.RequestStop(); functionEvaluations.Value++; bestError.Value = state.BestError; curError.Value = state.CurError; curQualityRow.Values.Add(state.CurError); bestQualityRow.Values.Add(bestError.Value); } // local function void State_ConstraintEvaluated(int constraintIdx, double value) { curConstraintIdx.Value = constraintIdx; curConstraintValue.Value = value; curConstraintRow.Values.Add(value); } } private static ISymbolicRegressionSolution CreateSolution(ISymbolicExpressionTree tree, IRegressionProblemData problemData) { var model = new SymbolicRegressionModel(problemData.TargetVariable, tree, new SymbolicDataAnalysisExpressionTreeLinearInterpreter()); // model.Scale(problemData); // var sol = model.CreateRegressionSolution((IRegressionProblemData)problemData.Clone()); // model.CreateRegressionSolution produces a new ProblemData and recalculates intervals var sol = new SymbolicRegressionSolution(model, (IRegressionProblemData)problemData.Clone()); // ==> set variable ranges to same range as in original problemData // foreach(var interval in problemData.VariableRanges.GetIntervals()) { // sol.ProblemData.VariableRanges.SetInterval // } return sol; } } }