Changeset 16329 for branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs

Timestamp:

11/27/18 08:21:13 (6 years ago)

Author:

gkronber

Message:

#2925: made several extensions in relation to blood glucose prediction

• added sqr function,
• added outputs for latent variables (linechart and model),
• added optimization of initial values for latent variables (for each episode separately)
• TODO: test with CVODES (so far only our own integration scheme has been tested)

File:

• branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs (modified) (14 diffs)

branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs

@@ -259 +259 @@
       }
       // init params randomly from Unif(-1e-5, 1e-5)
-      var theta = paramNodes.Select(_ => random.NextDouble() * 2.0e-2 - 1.0e-2).ToArray();
+      // theta contains parameter values for trees and then the initial values for latent variables (a separate vector for each episode)
+      // inital values for latent variables are also optimized
+      var theta = new double[paramNodes.Count + latentVariables.Length * episodes.Count()];
+      for (int i = 0; i < theta.Length; i++)
+        theta[i] = random.NextDouble() * 2.0e-2 - 1.0e-2;
 
       optTheta = new double[0];
@@ -355 +359 @@
       int r = 0;
       foreach (var y_pred in predicted) {
-        for (int c = 0; c < y_pred.Length; c++) {
+        // y_pred contains the predicted values for target variables first and then predicted values for latent variables
+        for (int c = 0; c < targetVariables.Length; c++) {
 
           var y_pred_f = y_pred[c].Item1;
@@ -413 +418 @@
         }
 
-        // only for actual target values
+        // only for target values
         var trainingRows = TrainingEpisodes.SelectMany(e => Enumerable.Range(e.Start, e.End - e.Start));
         for (int colIdx = 0; colIdx < targetVars.Length; colIdx++) {
@@ -450 +455 @@
                                    OdeSolver,
                                    NumericIntegrationSteps).ToArray();
-        // only for actual target values
+        // for target values and latent variables
         var trainingRows = TrainingEpisodes.SelectMany(e => Enumerable.Range(e.Start, e.End - e.Start));
-        for (int colIdx = 0; colIdx < targetVars.Length; colIdx++) {
-          var targetVar = targetVars[colIdx];
-          var trainingDataTable = new DataTable(targetVar + " prediction (training)");
-          var actualValuesRow = new DataRow(targetVar, "The values of " + targetVar, problemData.Dataset.GetDoubleValues(targetVar, trainingRows));
-          var predictedValuesRow = new DataRow(targetVar + " pred.", "Predicted values for " + targetVar, trainingPrediction.Select(arr => arr[colIdx].Item1).ToArray());
-          trainingDataTable.Rows.Add(actualValuesRow);
-          trainingDataTable.Rows.Add(predictedValuesRow);
-          trainingList.Add(trainingDataTable);
+        for (int colIdx = 0; colIdx < trees.Length; colIdx++) {
+          // is target variable
+          if (colIdx < targetVars.Length) {
+            var targetVar = targetVars[colIdx];
+            var trainingDataTable = new DataTable(targetVar + " prediction (training)");
+            var actualValuesRow = new DataRow(targetVar, "The values of " + targetVar, problemData.Dataset.GetDoubleValues(targetVar, trainingRows));
+            var predictedValuesRow = new DataRow(targetVar + " pred.", "Predicted values for " + targetVar, trainingPrediction.Select(arr => arr[colIdx].Item1).ToArray());
+            trainingDataTable.Rows.Add(actualValuesRow);
+            trainingDataTable.Rows.Add(predictedValuesRow);
+            trainingList.Add(trainingDataTable);
+          } else {
+            var latentVar = latentVariables[colIdx - targetVars.Length];
+            var trainingDataTable = new DataTable(latentVar + " prediction (training)");
+            var predictedValuesRow = new DataRow(latentVar + " pred.", "Predicted values for " + latentVar, trainingPrediction.Select(arr => arr[colIdx].Item1).ToArray());
+            var emptyRow = new DataRow(latentVar);
+            trainingDataTable.Rows.Add(emptyRow);
+            trainingDataTable.Rows.Add(predictedValuesRow);
+            trainingList.Add(trainingDataTable);
+          }
         }
         // TODO: DRY for training and test
@@ -475 +491 @@
          NumericIntegrationSteps).ToArray();
 
-        for (int colIdx = 0; colIdx < targetVars.Length; colIdx++) {
-          var targetVar = targetVars[colIdx];
-          var testDataTable = new DataTable(targetVar + " prediction (test)");
-          var actualValuesRow = new DataRow(targetVar, "The values of " + targetVar, problemData.Dataset.GetDoubleValues(targetVar, testRows));
-          var predictedValuesRow = new DataRow(targetVar + " pred.", "Predicted values for " + targetVar, testPrediction.Select(arr => arr[colIdx].Item1).ToArray());
-          testDataTable.Rows.Add(actualValuesRow);
-          testDataTable.Rows.Add(predictedValuesRow);
-          testList.Add(testDataTable);
+        for (int colIdx = 0; colIdx < trees.Length; colIdx++) {
+          // is target variable
+          if (colIdx < targetVars.Length) {
+            var targetVar = targetVars[colIdx];
+            var testDataTable = new DataTable(targetVar + " prediction (test)");
+            var actualValuesRow = new DataRow(targetVar, "The values of " + targetVar, problemData.Dataset.GetDoubleValues(targetVar, testRows));
+            var predictedValuesRow = new DataRow(targetVar + " pred.", "Predicted values for " + targetVar, testPrediction.Select(arr => arr[colIdx].Item1).ToArray());
+            testDataTable.Rows.Add(actualValuesRow);
+            testDataTable.Rows.Add(predictedValuesRow);
+            testList.Add(testDataTable);
+
+          } else {
+            var latentVar = latentVariables[colIdx - targetVars.Length];
+            var testDataTable = new DataTable(latentVar + " prediction (test)");
+            var predictedValuesRow = new DataRow(latentVar + " pred.", "Predicted values for " + latentVar, testPrediction.Select(arr => arr[colIdx].Item1).ToArray());
+            var emptyRow = new DataRow(latentVar);
+            testDataTable.Rows.Add(emptyRow);
+            testDataTable.Rows.Add(predictedValuesRow);
+            testList.Add(testDataTable);
+          }
         }
 
@@ -492 +520 @@
 
         int nextParIdx = 0;
-        foreach (var tup in targetVars.Zip(trees, Tuple.Create)) {
-          var targetVarName = tup.Item1;
-          var tree = tup.Item2;
+        for (int idx = 0; idx < trees.Length; idx++) {
+          var varName = string.Empty;
+          if (idx < targetVars.Length) {
+            varName = targetVars[idx];
+          } else {
+            varName = latentVariables[idx - targetVars.Length];
+          }
+          var tree = trees[idx];
 
           // when we reference HeuristicLab.Problems.DataAnalysis.Symbolic we can translate symbols
-          var shownTree = new SymbolicExpressionTree(TranslateTreeNode(tree.Root, optTheta, ref nextParIdx));
+          var shownTree = new SymbolicExpressionTree(TranslateTreeNode(tree.Root, optTheta.ToArray(),
+            ref nextParIdx));
 
           // var shownTree = (SymbolicExpressionTree)tree.Clone();
@@ -514 +548 @@
           // }
 
-          var origTreeVar = new HeuristicLab.Core.Variable(targetVarName + "(original)");
+          var origTreeVar = new HeuristicLab.Core.Variable(varName + "(original)");
           origTreeVar.Value = (ISymbolicExpressionTree)tree.Clone();
           models.Add(origTreeVar);
-          var simplifiedTreeVar = new HeuristicLab.Core.Variable(targetVarName + "(simplified)");
+          var simplifiedTreeVar = new HeuristicLab.Core.Variable(varName + "(simplified)");
           simplifiedTreeVar.Value = TreeSimplifier.Simplify(shownTree);
           models.Add(simplifiedTreeVar);
@@ -541 +575 @@
 
     private static IEnumerable<Tuple<double, Vector>[]> Integrate(
-      ISymbolicExpressionTree[] trees, IDataset dataset, string[] inputVariables, string[] targetVariables, string[] latentVariables, IEnumerable<IntRange> episodes,
+      ISymbolicExpressionTree[] trees, IDataset dataset,
+      string[] inputVariables, string[] targetVariables, string[] latentVariables, IEnumerable<IntRange> episodes,
       double[] parameterValues,
       string odeSolver, int numericIntegrationSteps = 100) {
 
       // TODO: numericIntegrationSteps is only relevant for the HeuristicLab solver
-
+      var episodeIdx = 0;
       foreach (var episode in episodes) {
         var rows = Enumerable.Range(episode.Start, episode.End - episode.Start);
-        // return first value as stored in the dataset
-        yield return targetVariables
-          .Select(targetVar => Tuple.Create(dataset.GetDoubleValue(targetVar, rows.First()), Vector.Zero))
-          .ToArray();
 
         // integrate forward starting with known values for the target in t0
@@ -565 +596 @@
         }
         // add value entries for latent variables which are also integrated
+        // initial values are at the end of the parameter vector
+        // separete initial values for each episode
+        var initialValueIdx = parameterValues.Length - episodes.Count() * latentVariables.Length + episodeIdx * latentVariables.Length;
         foreach (var latentVar in latentVariables) {
-          variableValues.Add(latentVar, Tuple.Create(0.0, Vector.Zero)); // we don't have observations for latent variables -> assume zero as starting value TODO
-        }
+          var arr = new double[parameterValues.Length]; // backing array
+          arr[initialValueIdx] = 1.0;
+          var g = new Vector(arr);
+          variableValues.Add(latentVar,
+            Tuple.Create(parameterValues[initialValueIdx], g)); // we don't have observations for latent variables therefore we optimize the initial value for each episode 
+          initialValueIdx++;
+        }
+
         var calculatedVariables = targetVariables.Concat(latentVariables).ToArray(); // TODO: must conincide with the order of trees in the encoding
+
+        // return first value as stored in the dataset
+        yield return calculatedVariables
+          .Select(calcVarName => variableValues[calcVarName])
+          .ToArray();
 
         var prevT = rows.First(); // TODO: here we should use a variable for t if it is available. Right now we assume equidistant measurements.
@@ -582 +627 @@
           //if (variableValues.Count == targetVariables.Length) {
           // only return the target variables for calculation of errors
-          var res = targetVariables
+          var res = calculatedVariables
             .Select(targetVar => variableValues[targetVar])
             .ToArray();
@@ -597 +642 @@
           }
         }
+        episodeIdx++;
       }
     }
@@ -1012 +1058 @@
             );
           }
+        case "sqr": {
+            var x = InterpretRec(node.GetSubtree(0), nodeValues);
+            return Tuple.Create(
+              x.Item1 * x.Item1,
+              2.0 * x.Item1 * x.Item2
+            );
+          }
         default: {
             return nodeValues[node];  // value and gradient for constants and variables must be set by the caller
@@ -1106 +1159 @@
       l.Add(new StringValue("sin").AsReadOnly());
       l.Add(new StringValue("cos").AsReadOnly());
+      l.Add(new StringValue("sqr").AsReadOnly());
       return l.AsReadOnly();
     }
@@ -1195 +1249 @@
       } else if (n.Symbol.Name == "cos") {
         translatedNode = new Cosine().CreateTreeNode();
+      } else if (n.Symbol.Name == "sqr") {
+        translatedNode = new Square().CreateTreeNode();
       } else if (IsConstantNode(n)) {
         var constNode = (ConstantTreeNode)new Constant().CreateTreeNode();