Kriging, Setting up custom correlation functions/kernels with more than 1 hyperparameter

I try to use kriging with customized kernels. Unfortunately I run into a problem for kernels with more than 1 hyperparameter. My kernel is stated the following way:

function R = my_first_kernel(x1,x2,theta,~)
exp_sin_sq = @(a,b,c,d) exp(-2*sin(pi/c*abs(a-b)).^2/d^2);
R = exp_sin_sq(x1(:,1),x2(:,1),theta(1),theta(2))
end

If I use theta(1) instead of theta(2) the kernel works. However, in this case I only have one hyperparameter to train, what I don’t want.

I use the kernel by setting the Correlation Family as a handle:

X = 0:0.25:2*pi;
Y = sin(X);
MetaOpts.Type = 'Metamodel';
MetaOpts.MetaType = 'Kriging';
MetaOpts.ExpDesign.Sampling = 'User';
MetaOpts.Corr.Family = @my_first_kernel;
MetaOpts.ExpDesign.X = X;
MetaOpts.ExpDesign.Y = Y;
MetaOpts.Optim.InitialValue = ones(2,1);
test = uq_createModel(MetaOpts);

I saw in the user manual that \theta should be a vector: “The input theta of the function my_eval_R corresponds to the hyperparameters θ and it is expected to be a vector of arbitrary length.” If I state the kernel as in the beginning, I get the error: " Index exceeds the number of array elements (1). " This occurs as the theta handed to the kernel is only a scalar and not a vector. I tried to make theta a vector by setting the optimization initial values to a vector: MetaOpts.Optim.InitialValue = ones(2,1) (also = ones(1,2)) This unfortunately didn’t help. I hope I stated the problem clear enough.

Some system information:

  • UQLab: version 1.4.0
  • Matlab: version R2019b
  • Operating system: windows

Dear Paul

You can make your example work by specifying the bounds of your hyperparameters as a 2\times M matrix. For example:

MetaOpts.Optim.Bounds = [1e-3, 1e-3; 1e4, 1e4];

Best regards
Styfen

Dear Styfen,

thanks for your reply. If I add the bounds to the code above, I get the error Dimension mismatch between Optim.Bounds and Experimental Design! and the model does not initialize. So unfotunately this does not work.

Best regards,
Paul

Dear Paul

When you specify your kernel as a handle at least this part of the code works (I tested it this time :slight_smile: ).

MetaOpts.Optim.Bounds = [1e-3, 1e-3; 1e4, 1e4];
MetaOpts.Corr.Handle = @my_first_kernel;

However, if this example shows the kernel you’re actually going to use, you may run into other problems because your matrix R has an incorrect shape. Its size should be size(x1, 1) by size(x2, 1).

Please let me know in case you encounter other problems.

Best regards
Styfen

1 Like

Dear Styfen,

thank you a lot. This solved the problem.

Best regards,
Paul