Changeset 17200 for branches/2994-AutoDiffForIntervals/HeuristicLab.Problems.DataAnalysis.Regression.Symbolic.Extensions

Timestamp:

08/12/19 10:50:02 (5 years ago)

Author:

gkronber

Message:

#2994 implemented checking of autodiff with numeric differentials and fixed a bug in IntervalEvaluator.

File:

• branches/2994-AutoDiffForIntervals/HeuristicLab.Problems.DataAnalysis.Regression.Symbolic.Extensions/ConstrainedNLSInternal.cs (modified) (5 diffs)

branches/2994-AutoDiffForIntervals/HeuristicLab.Problems.DataAnalysis.Regression.Symbolic.Extensions/ConstrainedNLSInternal.cs

@@ -1 +1 @@
 using System;
 using System.Collections.Generic;
+using System.Diagnostics;
 using System.Linq;
 using System.Runtime.InteropServices;
@@ -193 +194 @@
       bestSolution = x;
       bestError = minf;
-      
+
       if (res < 0) {
         throw new InvalidOperationException($"NLOpt failed {res} {NLOpt.nlopt_get_errmsg(nlopt)}");
@@ -232 +233 @@
         // average
         for (int j = 0; j < grad.Length; j++) { grad[j] /= target.Length; }
+
+        #region check gradient
+        var eval = new VectorEvaluator();
+        for (int i = 0; i < dim; i++) {
+          // make two additional evaluations
+          var xForNumDiff = (double[])curX.Clone();
+          double delta = Math.Abs(xForNumDiff[i] * 1e-5);
+          xForNumDiff[i] += delta;
+          UpdateThetaValues(xForNumDiff);
+          var evalHigh = eval.Evaluate(preparedTree, problemData.Dataset, trainingRows);
+          var mseHigh = MSE(target, evalHigh);
+          xForNumDiff[i] = curX[i] - delta;
+          UpdateThetaValues(xForNumDiff);
+          var evalLow = eval.Evaluate(preparedTree, problemData.Dataset, trainingRows);
+          var mseLow = MSE(target, evalLow);
+
+          var numericDiff = (mseHigh - mseLow) / (2 * delta);
+          var autoDiff = grad[i];
+          if ((Math.Abs(autoDiff) < 1e-10 && Math.Abs(numericDiff) > 1e-2)
+            || (Math.Abs(autoDiff) >= 1e-10 && Math.Abs((numericDiff - autoDiff) / numericDiff) > 1e-2))
+            throw new InvalidProgramException();
+        }
+        #endregion
       } else {
         var eval = new VectorEvaluator();
         var prediction = eval.Evaluate(preparedTree, problemData.Dataset, trainingRows);
 
-        // calc sum of squared errors and gradient
+        // calc sum of squared errors
         sse = 0.0;
         for (int i = 0; i < target.Length; i++) {
@@ -249 +273 @@
       if (double.IsNaN(sse)) return double.MaxValue;
       return sse / target.Length;
+    }
+
+    private double MSE(double[] a, double[] b) {
+      Trace.Assert(a.Length == b.Length);
+      var sse = 0.0;
+      for (int i = 0; i < a.Length; i++) sse += (a[i] - b[i]) * (a[i] - b[i]);
+      return sse / a.Length;
     }
 
@@ -280 +311 @@
       // we transformed this to a constraint c(x) <= 0, so only the upper bound is relevant for us
       if (grad != null) for (int j = 0; j < grad.Length; j++) { grad[j] = upperGradient[j]; }
+
+      #region check gradient
+      for (int i = 0; i < dim; i++) {
+        // make two additional evaluations
+        var xForNumDiff = (double[])curX.Clone();
+        double delta = Math.Abs(xForNumDiff[i] * 1e-5);
+        xForNumDiff[i] += delta;
+        UpdateThetaValues(xForNumDiff);
+        var evalHigh = intervalEvaluator.Evaluate(constraintData.Tree, dataIntervals, constraintData.ParameterNodes,
+          out double[] unusedLowerGradientHigh, out double[] unusedUpperGradientHigh);
+
+        xForNumDiff[i] = curX[i] - delta;
+        UpdateThetaValues(xForNumDiff);
+        var evalLow = intervalEvaluator.Evaluate(constraintData.Tree, dataIntervals, constraintData.ParameterNodes,
+          out double[] unusedLowerGradientLow, out double[] unusedUpperGradientLow);
+        
+        var numericDiff = (evalHigh.UpperBound - evalLow.UpperBound) / (2 * delta);
+        var autoDiff = grad[i];
+        // XXX GRADIENT FOR UPPER BOUND FOR FUNCTION IS ZERO FOR THE FIRST SET OF VARIABLES? GRADIENT FOR ADDITION INCORRECT?!
+        if ((Math.Abs(autoDiff) < 1e-10 && Math.Abs(numericDiff) > 1e-2)
+          || (Math.Abs(autoDiff) >= 1e-10 && Math.Abs((numericDiff - autoDiff) / numericDiff) > 1e-2))
+          throw new InvalidProgramException();
+      }
+      #endregion
+
+
       UpdateConstraintViolations(constraintData.Idx, interval.UpperBound);
       if (double.IsNaN(interval.UpperBound)) return double.MaxValue;