Changeset 8484 for trunk/sources/HeuristicLab.Algorithms.DataAnalysis/3.4/GaussianProcess/CovariancePeriodic.cs

Timestamp:

08/14/12 13:25:17 (12 years ago)

Author:

gkronber

Message:

#1902 changed interface for covariance functions to improve readability, fixed several bugs in the covariance functions and in the line chart for Gaussian process models.

File:

• trunk/sources/HeuristicLab.Algorithms.DataAnalysis/3.4/GaussianProcess/CovariancePeriodic.cs (modified) (6 diffs)

trunk/sources/HeuristicLab.Algorithms.DataAnalysis/3.4/GaussianProcess/CovariancePeriodic.cs

@@ -21 +21 @@
 
 using System;
+using System.Collections.Generic;
 using HeuristicLab.Common;
 using HeuristicLab.Core;
@@ -29 +30 @@
   [Item(Name = "CovariancePeriodic", Description = "Periodic covariance function for Gaussian processes.")]
   public class CovariancePeriodic : Item, ICovarianceFunction {
-    [Storable]
-    private double[,] x;
-    [Storable]
-    private double[,] xt;
     [Storable]
     private double sf2;
@@ -43 +40 @@
     public double Period { get { return p; } }
 
-    private bool symmetric;
-
-    private double[,] sd;
     public int GetNumberOfParameters(int numberOfVariables) {
       return 3;
@@ -53 +47 @@
     protected CovariancePeriodic(CovariancePeriodic original, Cloner cloner)
       : base(original, cloner) {
-      if (original.x != null) {
-        x = new double[original.x.GetLength(0), original.x.GetLength(1)];
-        Array.Copy(original.x, x, x.Length);
-        xt = new double[original.xt.GetLength(0), original.xt.GetLength(1)];
-        Array.Copy(original.xt, xt, xt.Length);
-      }
       sf2 = original.sf2;
       l = original.l;
       p = original.p;
-      symmetric = original.symmetric;
     }
     public CovariancePeriodic()
@@ -77 +64 @@
       this.p = Math.Exp(hyp[1]);
       this.sf2 = Math.Exp(2 * hyp[2]);
-      // sf2 = Math.Min(10E6, sf2); // upper limit for the scale
-
-      sd = null;
-    }
-    public void SetData(double[,] x) {
-      SetData(x, x);
-      this.symmetric = true;
     }
 
-    public void SetData(double[,] x, double[,] xt) {
-      this.x = x;
-      this.xt = xt;
-      this.symmetric = false;
-
-      sd = null;
-    }
-
-    public double GetCovariance(int i, int j) {
-      if (sd == null) CalculateSquaredDistances();
-      double k = sd[i, j];
+    public double GetCovariance(double[,] x, int i, int j) {
+      double k = i == j ? 0.0 : GetDistance(x, x, i, j);
       k = Math.PI * k / p;
       k = Math.Sin(k) / l;
@@ -104 +75 @@
     }
 
-    public double GetGradient(int i, int j, int k) {
-      double v = Math.PI * sd[i, j] / p;
-      switch (k) {
-        case 0: {
-            double newK = Math.Sin(v) / l;
-            newK = newK * newK;
-            return 4 * sf2 * Math.Exp(-2 * newK) * newK;
-          }
-        case 1: {
-            double r = Math.Sin(v) / l;
-            return 4 * sf2 / l * Math.Exp(-2 * r * r) * r * Math.Cos(v) * v;
-          }
-        case 2: {
-            double newK = Math.Sin(v) / l;
-            newK = newK * newK;
-            return 2 * sf2 * Math.Exp(-2 * newK);
-
-          }
-        default: {
-            throw new ArgumentException("CovariancePeriodic only has three hyperparameters.", "k");
-          }
-      }
+    public IEnumerable<double> GetGradient(double[,] x, int i, int j) {
+      double v = i == j ? 0.0 : Math.PI * GetDistance(x, x, i, j) / p;
+      double gradient = Math.Sin(v) / l;
+      gradient *= gradient;
+      yield return 4.0 * sf2 * Math.Exp(-2.0 * gradient) * gradient;
+      double r = Math.Sin(v) / l;
+      yield return 4.0 * sf2 / l * Math.Exp(-2 * r * r) * r * Math.Cos(v) * v;
+      yield return 2.0 * sf2 * Math.Exp(-2 * gradient);
     }
 
-    private void CalculateSquaredDistances() {
-      if (x.GetLength(1) != xt.GetLength(1)) throw new InvalidOperationException();
-      int rows = x.GetLength(0);
-      int cols = xt.GetLength(0);
-      sd = new double[rows, cols];
+    public double GetCrossCovariance(double[,] x, double[,] xt, int i, int j) {
+      double k = GetDistance(x, xt, i, j);
+      k = Math.PI * k / p;
+      k = Math.Sin(k) / l;
+      k = k * k;
 
-      if (symmetric) {
-        for (int i = 0; i < rows; i++) {
-          for (int j = i; j < cols; j++) {
-            sd[i, j] = Math.Sqrt(Util.SqrDist(Util.GetRow(x, i), Util.GetRow(x, j)));
-            sd[j, i] = sd[i, j];
-          }
-        }
-      } else {
-        for (int i = 0; i < rows; i++) {
-          for (int j = 0; j < cols; j++) {
-            sd[i, j] = Math.Sqrt(Util.SqrDist(Util.GetRow(x, i), Util.GetRow(xt, j)));
-          }
-        }
-      }
+      return sf2 * Math.Exp(-2.0 * k);
+    }
+
+    private double GetDistance(double[,] x, double[,] xt, int i, int j) {
+      return Math.Sqrt(Util.SqrDist(Util.GetRow(x, i), Util.GetRow(xt, j)));
     }
   }