Changeset 15433 for branches/MathNetNumerics-Exploration-2789

Timestamp:

10/27/17 11:10:25 (7 years ago)

Author:

gkronber

Message:

#2789 added the possibility to include interaction terms in GAM

Location:

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

Files:

• GAM.cs (modified) (8 diffs)
• Splines.cs (modified) (1 diff)

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

@@ -47 +47 @@
     private GAM(bool deserializing) : base(deserializing) { }
     [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() {
+    private void AfterDeserialization() {      
     }
 
@@ -62 +62 @@
       Parameters.Add(new ValueParameter<DoubleValue>("Lambda", "Regularization for smoothing splines", new DoubleValue(1.0)));
       Parameters.Add(new ValueParameter<IntValue>("Max iterations", "", new IntValue(100)));
-    }
-
+      Parameters.Add(new ValueParameter<IntValue>("Max interactions", "", new IntValue(1)));
+    }    
 
     protected override void Run(CancellationToken cancellationToken) {
@@ -69 +69 @@
       double lambda = ((IValueParameter<DoubleValue>)Parameters["Lambda"]).Value.Value;
       int maxIters = ((IValueParameter<IntValue>)Parameters["Max iterations"]).Value.Value;
+      int maxInteractions = ((IValueParameter<IntValue>)Parameters["Max interactions"]).Value.Value;
+      if (maxInteractions < 1 || maxInteractions > 5) throw new ArgumentException("Max interactions is outside the valid range [1 .. 5]");
 
       // calculates a GAM model using a linear representation + independent non-linear functions of each variable
@@ -77 +79 @@
       var avgY = y.Average();
       var inputVars = Problem.ProblemData.AllowedInputVariables.ToArray();
-      IRegressionModel[] f = new IRegressionModel[inputVars.Length * 2]; // linear(x) + nonlinear(x)
+      var nTerms = inputVars.Length; // LR
+      for(int i=1;i<=maxInteractions;i++) {
+        nTerms += inputVars.Combinations(i).Count();
+      }
+      IRegressionModel[] f = new IRegressionModel[nTerms]; 
       for(int i=0;i<f.Length;i++) {
         f[i] = new ConstantModel(0.0, problemData.TargetVariable);
@@ -90 +96 @@
       Results.Add(new Result("RMSE", rmseTable));
       rmseRow.Values.Add(CalculateResiduals(problemData, f, -1, avgY, problemData.TrainingIndices).StandardDeviation()); // -1 index to use all predictors
-      rmseRowTest.Values.Add(CalculateResiduals(problemData, f, -1, avgY, problemData.TestIndices).StandardDeviation()); 
+      rmseRowTest.Values.Add(CalculateResiduals(problemData, f, -1, avgY, problemData.TestIndices).StandardDeviation());
+
+      // for analytics
+      double[] rss = new double[f.Length];
+      string[] terms = new string[f.Length];
+      Results.Add(new Result("RSS Values", typeof(DoubleMatrix)));
 
       // until convergence
@@ -99 +110 @@
         foreach (var inputVar in inputVars) {
           var res = CalculateResiduals(problemData, f, j, avgY, problemData.TrainingIndices);
+          rss[j] = res.Variance();
+          terms[j] = inputVar;
           f[j] = RegressLR(problemData, inputVar, res);
           j++;
         }
-        foreach (var inputVar in inputVars) {
-          var res = CalculateResiduals(problemData, f, j, avgY, problemData.TrainingIndices);
-          f[j] = RegressSpline(problemData, inputVar, res, lambda);
-          j++;
+
+        for(int interaction = 1; interaction <= maxInteractions;interaction++) {
+          var selectedVars = HeuristicLab.Common.EnumerableExtensions.Combinations(inputVars, interaction);
+
+          foreach (var element in selectedVars) {
+            var res = CalculateResiduals(problemData, f, j, avgY, problemData.TrainingIndices);
+            rss[j] = res.Variance();
+            terms[j] = string.Format("f({0})", string.Join(",", element));
+            f[j] = RegressSpline(problemData, element.ToArray(), res, lambda);
+            j++;
+          }
         }
 
         rmseRow.Values.Add(CalculateResiduals(problemData, f, -1, avgY, problemData.TrainingIndices).StandardDeviation()); // -1 index to use all predictors
         rmseRowTest.Values.Add(CalculateResiduals(problemData, f, -1, avgY, problemData.TestIndices).StandardDeviation());
+
+        // calculate table with residual contributions of each term
+        var rssTable = new DoubleMatrix(rss.Length, 1, new string[] { "RSS" }, terms);
+        for (int i = 0; i < rss.Length; i++) rssTable[i, 0] = rss[i];
+        Results["RSS Values"].Value = rssTable;
 
         if (cancellationToken.IsCancellationRequested) break;
@@ -117 +142 @@
       model.AverageModelEstimates = false;
       Results.Add(new Result("Ensemble solution", model.CreateRegressionSolution((IRegressionProblemData)problemData.Clone())));
+
     }
 
@@ -158 +184 @@
       } else return new ConstantModel(target.Average(), problemData.TargetVariable);
     }
+    private IRegressionModel RegressSpline(IRegressionProblemData problemData, string[] inputVars, double[] target, double lambda) {
+      if (inputVars.All(problemData.Dataset.VariableHasType<double>)) {
+        var product = problemData.Dataset.GetDoubleValues(inputVars.First(), problemData.TrainingIndices).ToArray();
+        for(int i = 1;i<inputVars.Length;i++) {
+          product = product.Zip(problemData.Dataset.GetDoubleValues(inputVars[i], problemData.TrainingIndices), (pi, vi) => pi * vi).ToArray();
+        }
+        // Umständlich!
+        return Splines.CalculatePenalizedRegressionSpline(
+          product,
+          (double[])target.Clone(), lambda,
+          problemData.TargetVariable, inputVars
+          );
+      } else return new ConstantModel(target.Average(), problemData.TargetVariable);
+    }
+
     private IRegressionModel RegressRF(IRegressionProblemData problemData, string inputVar, double[] target, double lambda) {
       if (problemData.Dataset.VariableHasType<double>(inputVar)) {

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

@@ -763 +763 @@
 
     public IEnumerable<double> GetEstimatedValues(IDataset dataset, IEnumerable<int> rows) {
-      foreach (var x in dataset.GetDoubleValues(VariablesUsedForPrediction.First(), rows)) {
+      var product = dataset.GetDoubleValues(VariablesUsedForPrediction.First(), rows).ToArray();
+      foreach(var factor in VariablesUsedForPrediction.Skip(1)) {
+        product = product.Zip(dataset.GetDoubleValues(factor, rows), (pi, fi) => pi * fi).ToArray();
+      }
+
+      foreach (var x in product) {
         yield return alglib.spline1dcalc(interpolant, x);
       }