Changeset 15322 for branches/MathNetNumerics-Exploration-2789

Timestamp:

08/10/17 17:16:25 (7 years ago)

Author:

gkronber

Message:

#2789 worked on NLR with constraints

Location:

branches/MathNetNumerics-Exploration-2789/HeuristicLab.Algorithms.DataAnalysis.Experimental

Files:

• SymbolicRegressionConstantOptimizationEvaluator.cs (modified) (15 diffs)
• TreeToDiffSharpConverter.cs (modified) (1 diff)

branches/MathNetNumerics-Exploration-2789/HeuristicLab.Algorithms.DataAnalysis.Experimental/SymbolicRegressionConstantOptimizationEvaluator.cs

@@ -33 +33 @@
 using HeuristicLab.Problems.DataAnalysis.Symbolic.Regression;
 using DiffSharp.Interop.Float64;
+using System.Diagnostics;
 
 namespace HeuristicLab.Algorithms.DataAnalysis.Experimental {
@@ -179 +180 @@
 
       // extract inital constants
-      double[] c = new double[initialConstants.Length + 2];
+      double[] c = new double[initialConstants.Length];
       double[] s = new double[c.Length];
       {
-        c[0] = 1.0;
-        c[1] = 0.0;
-        Array.Copy(initialConstants, 0, c, 2, initialConstants.Length);
+        Array.Copy(initialConstants, 0, c, 0, initialConstants.Length);
+
+        // c[0] = 1.0;
+        // c[1] = 0.0;
+        // Array.Copy(initialConstants, 0, c, 2, initialConstants.Length);
 
         // s[0] = 1.0;
@@ -200 +203 @@
 
       IDataset ds = problemData.Dataset;
-      double[,] x = new double[ds.Rows, parameters.Count];
+      double[,] x = new double[rows.Count(), parameters.Count];
       int col = 0;
       int row = 0;
@@ -206 +209 @@
         row = 0;
         // copy training rows
-        foreach (var r in rows) { 
+        foreach (var r in rows) {
           if (ds.VariableHasType<double>(info.variableName)) {
             x[row, col] = ds.GetDoubleValue(info.variableName, r + info.lag);
@@ -222 +225 @@
       // only count constraints
       var constraintRows = Enumerable.Range(0, problemData.Dataset.Rows).Where(rIdx => !string.IsNullOrEmpty(ds.GetStringValue("f(x) constraint-type", rIdx)));
-      double[,] constraints = new double[constraintRows.Count(), parameters.Count + 1]; // +1 for constraint for f(x)
-      string[,] comp = new string[constraintRows.Count(), parameters.Count + 1];
+
+      double[,] constraintX = new double[constraintRows.Count(), parameters.Count]; // inputs for constraint values
+      double[,] constraints = new double[constraintRows.Count(), parameters.Count + 1]; // constraint values +1 column for constraint values for f(x)
+      string[,] comp = new string[constraintRows.Count(), parameters.Count + 1]; // comparison types <= = >=
       int eqConstraints = 0;
       int ieqConstraints = 0;
@@ -233 +238 @@
         var compColName = string.Format("df/d({0}) constraint-type", info.variableName);
 
-        if (ds.VariableNames.Contains(colName)) {
-          foreach (var r in constraintRows) {
+        foreach (var r in constraintRows) {
+          constraintX[row, col] = ds.GetDoubleValue(info.variableName, r);
+          if (ds.VariableNames.Contains(colName)) {
             constraints[row, col] = ds.GetDoubleValue(colName, r);
             comp[row, col] = ds.GetStringValue(compColName, r);
@@ -240 +246 @@
             if (comp[row, col] == "EQ") eqConstraints++;
             else if (comp[row, col] == "LEQ" || comp[row, col] == "GEQ") ieqConstraints++;
-
-            row++;
           }
+
+          row++;
         }
         col++;
@@ -260 +266 @@
       double[] y = ds.GetDoubleValues(problemData.TargetVariable, rows).ToArray();
 
-      alglib.ndimensional_jac jac = CreateJac(x, y, ds, func);
-      double rho = 1000;
+      alglib.ndimensional_jac jac = CreateJac(x, y, constraintX, constraints, comp, func);
+      double rho = 10000;
       int outeriters = 3;
       int updateFreq = 10;
@@ -275 +281 @@
         alglib.minnlcsetscale(state, s);
         alglib.minnlcsetprecexactlowrank(state, updateFreq);
-        // TODO set constraints;
         alglib.minnlcsetnlc(state, eqConstraints, ieqConstraints);
         alglib.minnlcoptimize(state, jac, null, null);
@@ -288 +293 @@
       //  -8  => integrity check failed (e.g. gradient NaN
       if (rep.terminationtype != -7 && rep.terminationtype != -8)
-        UpdateConstants(tree, c.Skip(2).ToArray(), updateVariableWeights);
+        UpdateConstants(tree, c.ToArray(), updateVariableWeights);
+      else {
+        UpdateConstants(tree, Enumerable.Repeat(0.0, c.Length).ToArray(), updateVariableWeights);
+      }
       var quality = SymbolicRegressionSingleObjectivePearsonRSquaredEvaluator.Calculate(interpreter, tree, lowerEstimationLimit, upperEstimationLimit, problemData, rows, applyLinearScaling);
-
-      if (!updateConstantsInTree) UpdateConstants(tree, originalConstants.Skip(2).ToArray(), updateVariableWeights);
+      Console.WriteLine("{0:N4} {1:N4} {2} {3}", originalQuality, quality, state.fi.Skip(1).Where(fii => fii > 0).Count(), rep.terminationtype);
+
+      if (!updateConstantsInTree) UpdateConstants(tree, originalConstants.ToArray(), updateVariableWeights);
       if (
         // originalQuality - quality > 0.001 || 
         double.IsNaN(quality)) {
-        UpdateConstants(tree, originalConstants.Skip(2).ToArray(), updateVariableWeights);
+        UpdateConstants(tree, originalConstants.ToArray(), updateVariableWeights);
         return originalQuality;
       }
@@ -319 +328 @@
 
     private static alglib.ndimensional_jac CreateJac(
-      double[,] x, // x is larger than y 
+      double[,] x, // x is same size as y 
       double[] y, // only targets
-      // double[,] constraints, // df/d(xi), same order as for x, same number of rows as x less columns
-      // string[,] comparison, // {LEQ, GEQ, EQ } same size as constraints
-      IDataset ds,
+      double[,] constraintX, // inputs for constraints 
+      double[,] constraints, // df/d(xi), same size as constraintX, same number of rows as x less columns
+      string[,] comparison, // {LEQ, GEQ, EQ } same size as constraints
       Func<DV, D> func) {
+      Trace.Assert(x.GetLength(0) == y.Length);
+      Trace.Assert(x.GetLength(1) == constraintX.GetLength(1) - 1);
+      Trace.Assert(constraintX.GetLength(0) == constraints.GetLength(0));
+      Trace.Assert(constraintX.GetLength(1) == constraints.GetLength(1));
+      Trace.Assert(constraints.GetLength(0) == comparison.GetLength(0));
+      Trace.Assert(constraints.GetLength(1) == comparison.GetLength(1));
       return (double[] c, double[] fi, double[,] jac, object o) => {
         // objective function is sum of squared errors
@@ -350 +365 @@
 
         int fidx = 1;
-        var constraintRows = Enumerable.Range(0, ds.Rows).Where(rIdx => !string.IsNullOrEmpty(ds.GetStringValue("f(x) constraint-type", rIdx))); 
 
         // eq constraints
-        foreach (var rowIdx in constraintRows) {
-          foreach (var varName in ds.VariableNames.Where(vn => vn.EndsWith("constraint-type"))) {
-            if (ds.GetStringValue(varName, rowIdx) == "EQ") {
+        for (int rowIdx = 0; rowIdx < constraintX.GetLength(0); rowIdx++) {
+          for (var colIdx = 0; colIdx < constraintX.GetLength(1); colIdx++) {
+            if (comparison[rowIdx, colIdx] == "EQ") {
               throw new NotSupportedException();
             }
@@ -361 +375 @@
         }
         // ineq constraints
-        foreach (var rowIdx in constraintRows) {
-XXX 
+        for (int rowIdx = 0; rowIdx < constraintX.GetLength(0); rowIdx++) {
           for (int colIdx = 0; colIdx < constraints.GetLength(1); colIdx++) {
             // there is a constraint value
@@ -369 +382 @@
               : comparison[rowIdx, colIdx] == "GEQ" ? -1.0 : 0.0;
 
+              // copy x_i to the beginning of the function parameters vector
+              for (int cIdx = 0; cIdx < nParams; cIdx++)
+                p[cIdx] = constraintX[rowIdx, cIdx];
+
               // f(x) constraint
               if (colIdx == constraints.GetLength(1) - 1) {
-                // copy x_i to the beginning of the function parameters vector
-                for (int cIdx = 0; cIdx < nParams; cIdx++)
-                  p[cIdx] = x[rowIdx, cIdx];
 
                 fi[fidx] = factor * ((double)(func(p)) - constraints[rowIdx, colIdx]);
@@ -385 +399 @@
                 var g = (double[])AD.Grad(func, p);
                 fi[fidx] = factor * g[colIdx];
+
                 var hess = AD.Hessian(func, p);
                 for (int jacIdx = 0; jacIdx < c.Length; jacIdx++) {
-                  jac[fidx, jacIdx] = factor * (double)hess[nParams + colIdx, nParams + jacIdx]; // skip the elements for the variable values
+                  jac[fidx, jacIdx] = factor * (double)hess[colIdx, nParams + jacIdx]; // skip the elements for the variable values
                 }
                 fidx++;

branches/MathNetNumerics-Exploration-2789/HeuristicLab.Algorithms.DataAnalysis.Experimental/TreeToDiffSharpConverter.cs

@@ -306 +306 @@
       // }
       if (node.Symbol is StartSymbol) {
-        var alpha = Expression.Call(dv, DvIndexer, Expression.Constant(paramIdx++));
-        var beta = Expression.Call(dv, DvIndexer, Expression.Constant(paramIdx++));
-
-        return Expression.Call(d_Add_d, beta,
-          Expression.Call(d_Mul_d, alpha, MakeExpr(node.GetSubtree(0), parameters, dv)));
+        //var alpha = Expression.Call(dv, DvIndexer, Expression.Constant(paramIdx++));
+        //var beta = Expression.Call(dv, DvIndexer, Expression.Constant(paramIdx++));
+
+        // return Expression.Call(d_Add_d, beta,
+        //   Expression.Call(d_Mul_d, alpha, MakeExpr(node.GetSubtree(0), parameters, dv)));
+        return MakeExpr(node.GetSubtree(0), parameters, dv);
       }
       throw new ConversionException();