Changeset 15249 for stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression

Timestamp:

07/15/17 10:29:40 (7 years ago)

Author:

gkronber

Message:

#2699,#2700
merged r14862, r14863, r14911, r14936, r15156, r15157, r15158, r15164, r15169, r15207:15209, r15225, r15227, r15234, r15248 from trunk to stable

Location:

stable

Files:

• . (modified) (1 prop)
• HeuristicLab.Algorithms.DataAnalysis (modified) (1 prop)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression (copied) (copied from trunk/sources/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/CicularKernel.cs (modified) (2 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/GaussianKernel.cs (modified) (2 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/IKernel.cs (modified) (1 diff)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/InverseMultiquadraticKernel.cs (modified) (2 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/KernelBase.cs (modified) (4 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/MultiquadraticKernel.cs (modified) (3 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/PolysplineKernel.cs (modified) (3 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/ThinPlatePolysplineKernel.cs (modified) (3 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelRidgeRegression.cs (modified) (4 diffs)
• HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelRidgeRegressionModel.cs (modified) (11 diffs)

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/CicularKernel.cs

@@ -25 +25 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   [Item("CircularKernel", "A circular kernel function 2*pi*(acos(-d)-d*(1-d²)^(0.5)) where n = ||x-c|| and d = n/beta \n  As described in http://crsouza.com/2010/03/17/kernel-functions-for-machine-learning-applications/")]
   public class CircularKernel : KernelBase {
-
-    #region HLConstructors & Boilerplate
     [StorableConstructor]
     protected CircularKernel(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
+
     protected CircularKernel(CircularKernel original, Cloner cloner) : base(original, cloner) { }
-    public CircularKernel() {
-    }
+
+    public CircularKernel() { }
+
     public override IDeepCloneable Clone(Cloner cloner) {
       return new CircularKernel(this, cloner);
     }
-    #endregion
 
     protected override double Get(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;
@@ -53 +51 @@
     // 4*pi*n^3 / (beta^4 * sqrt(1-n^2/beta^2) 
     protected override double GetGradient(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/GaussianKernel.cs

@@ -27 +27 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   [Item("GaussianKernel", "A kernel function that uses Gaussian function exp(-n²/beta²). As described in http://crsouza.com/2010/03/17/kernel-functions-for-machine-learning-applications/")]
   public class GaussianKernel : KernelBase {
-
-    #region HLConstructors & Boilerplate
     [StorableConstructor]
     protected GaussianKernel(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
+
     protected GaussianKernel(GaussianKernel original, Cloner cloner) : base(original, cloner) { }
+
     public GaussianKernel() {
     }
+
     public override IDeepCloneable Clone(Cloner cloner) {
       return new GaussianKernel(this, cloner);
     }
-    #endregion
 
     protected override double Get(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;
@@ -54 +53 @@
     //2 * n²/b²* 1/b * exp(-n²/b²)
     protected override double GetGradient(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/IKernel.cs

@@ -21 +21 @@
 
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+using System;
+
+namespace HeuristicLab.Algorithms.DataAnalysis {
   public interface IKernel : ICovarianceFunction {
     double? Beta { get; set; } // a kernel parameter
     IDistance Distance { get; set; } // the distance function to use
+
+    event EventHandler BetaChanged;
+    event EventHandler DistanceChanged;
   }
 }

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/InverseMultiquadraticKernel.cs

@@ -25 +25 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   [Item("InverseMultiquadraticKernel", "A kernel function that uses the inverse multi-quadratic function  1 / sqrt(1+||x-c||²/beta²). Similar to http://crsouza.com/2010/03/17/kernel-functions-for-machine-learning-applications/ with beta as a scaling factor.")]
   public class InverseMultiquadraticKernel : KernelBase {
+    private const double C = 1.0;
 
-    private const double C = 1.0;
-    #region HLConstructors & Boilerplate
     [StorableConstructor]
     protected InverseMultiquadraticKernel(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
+
     protected InverseMultiquadraticKernel(InverseMultiquadraticKernel original, Cloner cloner) : base(original, cloner) { }
+
     public InverseMultiquadraticKernel() { }
+
     public override IDeepCloneable Clone(Cloner cloner) {
       return new InverseMultiquadraticKernel(this, cloner);
     }
-    #endregion
 
     protected override double Get(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;
@@ -52 +52 @@
     //n²/(b³(n²/b² + C)^1.5)
     protected override double GetGradient(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/KernelBase.cs

@@ -28 +28 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   public abstract class KernelBase : ParameterizedNamedItem, IKernel {
 
-    #region Parameternames
     private const string DistanceParameterName = "Distance";
-    #endregion
-    #region Parameterproperties
-    public ValueParameter<IDistance> DistanceParameter {
-      get { return Parameters[DistanceParameterName] as ValueParameter<IDistance>; }
+
+    public IValueParameter<IDistance> DistanceParameter {
+      get { return (IValueParameter<IDistance>)Parameters[DistanceParameterName]; }
     }
 
     [Storable]
-    public double? Beta { get; set; }
-    #endregion
-    #region Properties
+    private double? beta;
+    public double? Beta {
+      get { return beta; }
+      set {
+        if (value != beta) {
+          beta = value;
+          RaiseBetaChanged();
+        }
+      }
+    }
+
     public IDistance Distance {
       get { return DistanceParameter.Value; }
-      set { DistanceParameter.Value = value; }
+      set {
+        if (DistanceParameter.Value != value) {
+          DistanceParameter.Value = value;
+        }
+      }
     }
-
-    #endregion
 
     [StorableConstructor]
     protected KernelBase(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
 
     protected KernelBase(KernelBase original, Cloner cloner)
       : base(original, cloner) {
-      Beta = original.Beta;
+      beta = original.beta;
+      RegisterEvents();
     }
 
@@ -64 +71 @@
       Parameters.Add(new ValueParameter<IDistance>(DistanceParameterName, "The distance function used for kernel calculation"));
       DistanceParameter.Value = new EuclideanDistance();
+      RegisterEvents();
+    }
+
+    [StorableHook(HookType.AfterDeserialization)]
+    private void AfterDeserialization() {
+      RegisterEvents();
+    }
+
+    private void RegisterEvents() {
+      DistanceParameter.ValueChanged += (sender, args) => RaiseDistanceChanged();
     }
 
@@ -82 +99 @@
     public ParameterizedCovarianceFunction GetParameterizedCovarianceFunction(double[] p, int[] columnIndices) {
       if (p.Length != GetNumberOfParameters(columnIndices.Length)) throw new ArgumentException("Illegal parametrization");
-      var myClone = (KernelBase)Clone(new Cloner());
+      var myClone = (KernelBase)Clone();
       myClone.SetParameter(p);
       var cov = new ParameterizedCovarianceFunction {
@@ -101 +118 @@
       return dist.Get(r1, r2);
     }
+
+    #region events
+    public event EventHandler BetaChanged;
+    public event EventHandler DistanceChanged;
+
+    protected void RaiseBetaChanged() {
+      var handler = BetaChanged;
+      if (handler != null) handler(this, EventArgs.Empty);
+    }
+
+    protected void RaiseDistanceChanged() {
+      var handler = DistanceChanged;
+      if (handler != null) handler(this, EventArgs.Empty);
+    }
+    #endregion
   }
 }

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/MultiquadraticKernel.cs

@@ -25 +25 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   // conditionally positive definite. (need to add polynomials) see http://num.math.uni-goettingen.de/schaback/teaching/sc.pdf 
@@ -32 +32 @@
 
     private const double C = 1.0;
-    #region HLConstructors & Boilerplate
+
     [StorableConstructor]
     protected MultiquadraticKernel(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
-    protected MultiquadraticKernel(MultiquadraticKernel original, Cloner cloner)
-                : base(original, cloner) { }
 
-    public MultiquadraticKernel() {
-    }
+    protected MultiquadraticKernel(MultiquadraticKernel original, Cloner cloner) : base(original, cloner) { }
+
+    public MultiquadraticKernel() { }
+
     public override IDeepCloneable Clone(Cloner cloner) {
       return new MultiquadraticKernel(this, cloner);
     }
-    #endregion
+
     protected override double Get(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;
@@ -55 +54 @@
     //-n²/(d³*sqrt(C+n²/d²))
     protected override double GetGradient(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/PolysplineKernel.cs

@@ -27 +27 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   // conditionally positive definite. (need to add polynomials) see http://num.math.uni-goettingen.de/schaback/teaching/sc.pdf 
@@ -33 +33 @@
   public class PolysplineKernel : KernelBase {
 
-    #region Parameternames
     private const string DegreeParameterName = "Degree";
-    #endregion
-    #region Parameterproperties
+
     public IFixedValueParameter<DoubleValue> DegreeParameter {
-      get { return Parameters[DegreeParameterName] as IFixedValueParameter<DoubleValue>; }
+      get { return (IFixedValueParameter<DoubleValue>)Parameters[DegreeParameterName]; }
     }
-    #endregion
-    #region Properties
+
     public DoubleValue Degree {
       get { return DegreeParameter.Value; }
     }
-    #endregion
 
-    #region HLConstructors & Boilerplate
     [StorableConstructor]
     protected PolysplineKernel(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
+
     protected PolysplineKernel(PolysplineKernel original, Cloner cloner) : base(original, cloner) { }
+
     public PolysplineKernel() {
       Parameters.Add(new FixedValueParameter<DoubleValue>(DegreeParameterName, "The degree of the kernel. Needs to be greater than zero.", new DoubleValue(1.0)));
     }
+
     public override IDeepCloneable Clone(Cloner cloner) {
       return new PolysplineKernel(this, cloner);
     }
-    #endregion
 
     protected override double Get(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;
@@ -70 +66 @@
     //-degree/beta * (norm/beta)^degree
     protected override double GetGradient(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelFunctions/ThinPlatePolysplineKernel.cs

@@ -27 +27 @@
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   // conditionally positive definite. (need to add polynomials) see http://num.math.uni-goettingen.de/schaback/teaching/sc.pdf 
@@ -33 +33 @@
   public class ThinPlatePolysplineKernel : KernelBase {
 
-    #region Parameternames
     private const string DegreeParameterName = "Degree";
-    #endregion
-    #region Parameterproperties
+
     public IFixedValueParameter<DoubleValue> DegreeParameter {
-      get { return Parameters[DegreeParameterName] as IFixedValueParameter<DoubleValue>; }
+      get { return (IFixedValueParameter<DoubleValue>)Parameters[DegreeParameterName]; }
     }
-    #endregion
-    #region Properties
     public DoubleValue Degree {
       get { return DegreeParameter.Value; }
     }
-    #endregion
 
-    #region HLConstructors & Boilerplate
     [StorableConstructor]
     protected ThinPlatePolysplineKernel(bool deserializing) : base(deserializing) { }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
+
     protected ThinPlatePolysplineKernel(ThinPlatePolysplineKernel original, Cloner cloner) : base(original, cloner) { }
+
     public ThinPlatePolysplineKernel() {
       Parameters.Add(new FixedValueParameter<DoubleValue>(DegreeParameterName, "The degree of the kernel. Needs to be greater than zero.", new DoubleValue(2.0)));
     }
+
     public override IDeepCloneable Clone(Cloner cloner) {
       return new ThinPlatePolysplineKernel(this, cloner);
     }
-    #endregion
 
     protected override double Get(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;
@@ -71 +66 @@
     // (Degree/beta) * (norm/beta)^Degree * log(norm/beta) 
     protected override double GetGradient(double norm) {
+      if (Beta == null) throw new InvalidOperationException("Can not calculate kernel distance gradient while Beta is null");
       var beta = Beta.Value;
       if (Math.Abs(beta) < double.Epsilon) return double.NaN;

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelRidgeRegression.cs

@@ -29 +29 @@
 using HeuristicLab.Parameters;
 using HeuristicLab.Persistence.Default.CompositeSerializers.Storable;
+using HeuristicLab.PluginInfrastructure;
 using HeuristicLab.Problems.DataAnalysis;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [Item("Kernel Ridge Regression", "Kernelized ridge regression e.g. for radial basis function (RBF) regression.")]
   [Creatable(CreatableAttribute.Categories.DataAnalysisRegression, Priority = 100)]
@@ -57 +58 @@
 
     #region parameter properties
-    public ValueParameter<IKernel> KernelParameter {
-      get { return (ValueParameter<IKernel>)Parameters[KernelParameterName]; }
+    public IConstrainedValueParameter<IKernel> KernelParameter {
+      get { return (IConstrainedValueParameter<IKernel>)Parameters[KernelParameterName]; }
     }
 
@@ -102 +103 @@
     public KernelRidgeRegression() {
       Problem = new RegressionProblem();
-      Parameters.Add(new ValueParameter<IKernel>(KernelParameterName, "The kernel", new GaussianKernel()));
+      var values = new ItemSet<IKernel>(ApplicationManager.Manager.GetInstances<IKernel>());
+      Parameters.Add(new ConstrainedValueParameter<IKernel>(KernelParameterName, "The kernel", values, values.OfType<GaussianKernel>().FirstOrDefault()));
       Parameters.Add(new FixedValueParameter<BoolValue>(ScaleInputVariablesParameterName, "Set to true if the input variables should be scaled to the interval [0..1]", new BoolValue(true)));
       Parameters.Add(new FixedValueParameter<DoubleValue>(LambdaParameterName, "The log10-transformed weight for the regularization term lambda [-inf..+inf]. Small values produce more complex models, large values produce models with larger errors. Set to very small value (e.g. -1.0e15) for almost exact approximation", new DoubleValue(-2)));
-      Parameters.Add(new FixedValueParameter<DoubleValue>(BetaParameterName, "The beta parameter for the kernel", new DoubleValue(2)));
+      Parameters.Add(new FixedValueParameter<DoubleValue>(BetaParameterName, "The inverse width of the kernel ]0..+inf]. The distance between points is divided by this value before being plugged into the kernel.", new DoubleValue(2)));
     }
-    [StorableHook(HookType.AfterDeserialization)]
-    private void AfterDeserialization() { }
 
     public override IDeepCloneable Clone(Cloner cloner) {
@@ -125 +125 @@
 
     public static IRegressionSolution CreateRadialBasisRegressionSolution(IRegressionProblemData problemData, ICovarianceFunction kernel, double lambda, bool scaleInputs, out double rmsError, out double looCvRMSE) {
-      var model = new KernelRidgeRegressionModel(problemData.Dataset, problemData.TargetVariable, problemData.AllowedInputVariables, problemData.TrainingIndices, scaleInputs, kernel, lambda);
+      var model = KernelRidgeRegressionModel.Create(problemData.Dataset, problemData.TargetVariable, problemData.AllowedInputVariables, problemData.TrainingIndices, scaleInputs, kernel, lambda);
       rmsError = double.NaN;
       if (problemData.TestIndices.Any()) {

stable/HeuristicLab.Algorithms.DataAnalysis/3.4/KernelRidgeRegression/KernelRidgeRegressionModel.cs

@@ -28 +28 @@
 using HeuristicLab.Problems.DataAnalysis;
 
-namespace HeuristicLab.Algorithms.DataAnalysis.KernelRidgeRegression {
+namespace HeuristicLab.Algorithms.DataAnalysis {
   [StorableClass]
   [Item("KernelRidgeRegressionModel", "A kernel ridge regression model")]
@@ -39 +39 @@
     private readonly string[] allowedInputVariables;
     public string[] AllowedInputVariables {
-      get { return allowedInputVariables; }
+      get { return allowedInputVariables.ToArray(); }
     }
 
@@ -81 +81 @@
       yOffset = original.yOffset;
       yScale = original.yScale;
-      if (original.kernel != null)
-        kernel = cloner.Clone(original.kernel);
+      kernel = original.kernel;
     }
     public override IDeepCloneable Clone(Cloner cloner) {
@@ -88 +87 @@
     }
 
-    public KernelRidgeRegressionModel(IDataset dataset, string targetVariable, IEnumerable<string> allowedInputVariables, IEnumerable<int> rows,
-      bool scaleInputs, ICovarianceFunction kernel, double lambda = 0.1) : base(targetVariable) {
-      if (kernel.GetNumberOfParameters(allowedInputVariables.Count()) > 0) throw new ArgumentException("All parameters in the kernel function must be specified.");
-      name = ItemName;
-      description = ItemDescription;
-      this.allowedInputVariables = allowedInputVariables.ToArray();
+    public static KernelRidgeRegressionModel Create(IDataset dataset, string targetVariable, IEnumerable<string> allowedInputVariables, IEnumerable<int> rows,
+      bool scaleInputs, ICovarianceFunction kernel, double lambda = 0.1) {
       var trainingRows = rows.ToArray();
-      this.kernel = (ICovarianceFunction)kernel.Clone();
-      this.lambda = lambda;
+      var model = new KernelRidgeRegressionModel(dataset, targetVariable, allowedInputVariables, trainingRows, scaleInputs, kernel, lambda);
+
       try {
-        if (scaleInputs)
-          scaling = CreateScaling(dataset, trainingRows);
-        trainX = ExtractData(dataset, trainingRows, scaling);
+        int info;
+        int n = model.trainX.GetLength(0);
+        alglib.densesolverreport denseSolveRep;
+        var gram = BuildGramMatrix(model.trainX, lambda, kernel);
+        var l = new double[n, n];
+        Array.Copy(gram, l, l.Length);
+
+        double[] alpha = new double[n];
+        double[,] invG;
         var y = dataset.GetDoubleValues(targetVariable, trainingRows).ToArray();
-        yOffset = y.Average();
-        yScale = 1.0 / y.StandardDeviation();
         for (int i = 0; i < y.Length; i++) {
-          y[i] -= yOffset;
-          y[i] *= yScale;
-        }
-        int info;
-        int n = trainX.GetLength(0);
-        alglib.densesolverreport denseSolveRep;
-        var gram = BuildGramMatrix(trainX, lambda);
-        var l = new double[n, n]; Array.Copy(gram, l, l.Length);
-
-        double[,] invG;
+          y[i] -= model.yOffset;
+          y[i] *= model.yScale;
+        }
         // cholesky decomposition
         var res = alglib.trfac.spdmatrixcholesky(ref l, n, false);
-        if (res == false) { //throw new ArgumentException("Could not decompose matrix. Is it quadratic symmetric positive definite?");
+        if (res == false) { //try lua decomposition if cholesky faild
           int[] pivots;
           var lua = new double[n, n];
@@ -127 +119 @@
           invG = lua;  // rename
           alglib.rmatrixluinverse(ref invG, pivots, n, out info, out rep);
-          if (info != 1) throw new ArgumentException("Could not invert Gram matrix.");
         } else {
           alglib.spdmatrixcholeskysolve(l, n, false, y, out info, out denseSolveRep, out alpha);
@@ -135 +126 @@
           invG = l;   // rename 
           alglib.spdmatrixcholeskyinverse(ref invG, n, false, out info, out rep);
-          if (info != 1) throw new ArgumentException("Could not invert Gram matrix.");
-        }
+        }
+        if (info != 1) throw new ArgumentException("Could not invert Gram matrix.");
 
         var ssqLooError = 0.0;
@@ -142 +133 @@
           var pred_i = Util.ScalarProd(Util.GetRow(gram, i).ToArray(), alpha);
           var looPred_i = pred_i - alpha[i] / invG[i, i];
-          var error = (y[i] - looPred_i) / yScale;
+          var error = (y[i] - looPred_i) / model.yScale;
           ssqLooError += error * error;
         }
-        LooCvRMSE = Math.Sqrt(ssqLooError / n);
+
+        Array.Copy(alpha, model.alpha, n);
+        model.LooCvRMSE = Math.Sqrt(ssqLooError / n);
       } catch (alglib.alglibexception ae) {
         // wrap exception so that calling code doesn't have to know about alglib implementation
         throw new ArgumentException("There was a problem in the calculation of the kernel ridge regression model", ae);
       }
+      return model;
+    }
+
+    private KernelRidgeRegressionModel(IDataset dataset, string targetVariable, IEnumerable<string> allowedInputVariables, int[] rows,
+      bool scaleInputs, ICovarianceFunction kernel, double lambda = 0.1) : base(targetVariable) {
+      this.allowedInputVariables = allowedInputVariables.ToArray();
+      if (kernel.GetNumberOfParameters(this.allowedInputVariables.Length) > 0) throw new ArgumentException("All parameters in the kernel function must be specified.");
+      name = ItemName;
+      description = ItemDescription;
+
+      this.kernel = (ICovarianceFunction)kernel.Clone();
+      this.lambda = lambda;
+      if (scaleInputs) scaling = CreateScaling(dataset, rows, this.allowedInputVariables);
+      trainX = ExtractData(dataset, rows, this.allowedInputVariables, scaling);
+      var y = dataset.GetDoubleValues(targetVariable, rows).ToArray();
+      yOffset = y.Average();
+      yScale = 1.0 / y.StandardDeviation();
+      alpha = new double[trainX.GetLength(0)];
     }
 
@@ -155 +166 @@
     #region IRegressionModel Members
     public override IEnumerable<double> GetEstimatedValues(IDataset dataset, IEnumerable<int> rows) {
-      var newX = ExtractData(dataset, rows, scaling);
+      var newX = ExtractData(dataset, rows, allowedInputVariables, scaling);
       var dim = newX.GetLength(1);
       var cov = kernel.GetParameterizedCovarianceFunction(new double[0], Enumerable.Range(0, dim).ToArray());
@@ -175 +186 @@
 
     #region helpers
-    private double[,] BuildGramMatrix(double[,] data, double lambda) {
+    private static double[,] BuildGramMatrix(double[,] data, double lambda, ICovarianceFunction kernel) {
       var n = data.GetLength(0);
       var dim = data.GetLength(1);
@@ -190 +201 @@
     }
 
-    private ITransformation<double>[] CreateScaling(IDataset dataset, int[] rows) {
-      var trans = new ITransformation<double>[allowedInputVariables.Length];
+    private static ITransformation<double>[] CreateScaling(IDataset dataset, int[] rows, IReadOnlyCollection<string> allowedInputVariables) {
+      var trans = new ITransformation<double>[allowedInputVariables.Count];
       int i = 0;
       foreach (var variable in allowedInputVariables) {
@@ -205 +216 @@
     }
 
-    private double[,] ExtractData(IDataset dataset, IEnumerable<int> rows, ITransformation<double>[] scaling = null) {
+    private static double[,] ExtractData(IDataset dataset, IEnumerable<int> rows, IReadOnlyCollection<string> allowedInputVariables, ITransformation<double>[] scaling = null) {
       double[][] variables;
       if (scaling != null) {