Условные симуляции кригинга/гауссова процесса в Matlab

Question

Я хотел бы выполнить условное моделирование для моделей Gaussian process (GP) в Matlab. Я нашел учебник Мартина Коларжа (http://mrmartin.net/?p=223).

sigma_f = 1.1251; %parameter of the squared exponential kernel 
l = 0.90441; %parameter of the squared exponential kernel 
kernel_function = @(x,x2) sigma_f^2*exp((x-x2)^2/(-2*l^2)); 

%This is one of many popular kernel functions, the squared exponential 
%kernel. It favors smooth functions. (Here, it is defined here as an anonymous 
%function handle) 

% we can also define an error function, which models the observation noise 
sigma_n = 0.1; %known noise on observed data 
error_function = @(x,x2) sigma_n^2*(x==x2); 
%this is just iid gaussian noise with mean 0 and variance sigma_n^2s 

%kernel functions can be added together. Here, we add the error kernel to 
%the squared exponential kernel) 
k = @(x,x2) kernel_function(x,x2)+error_function(x,x2); 

X_o = [-1.5 -1 -0.75 -0.4 -0.3 0]'; 
Y_o = [-1.6 -1.3 -0.5 0 0.3 0.6]'; 

prediction_x=-2:0.01:1; 

K = zeros(length(X_o)); 
for i=1:length(X_o) 
    for j=1:length(X_o) 
     K(i,j)=k(X_o(i),X_o(j)); 
    end 
end 

%% Demo #5.2 Sample from the Gaussian Process posterior 
clearvars -except k prediction_x K X_o Y_o 

%We can also sample from this posterior, the same way as we sampled before: 
K_ss=zeros(length(prediction_x),length(prediction_x)); 
for i=1:length(prediction_x) 
    for j=i:length(prediction_x)%We only calculate the top half of the matrix. This an unnecessary speedup trick 
     K_ss(i,j)=k(prediction_x(i),prediction_x(j)); 
    end 
end 
K_ss=K_ss+triu(K_ss,1)'; % We can use the upper half of the matrix and copy it to the 

K_s=zeros(length(prediction_x),length(X_o)); 
for i=1:length(prediction_x) 
    for j=1:length(X_o) 
     K_s(i,j)=k(prediction_x(i),X_o(j)); 
    end 
end 

[V,D]=eig(K_ss-K_s/K*K_s'); 
A=real(V*(D.^(1/2))); 

for i=1:7 
    standard_random_vector = randn(length(A),1); 
    gaussian_process_sample(:,i) = A * standard_random_vector+K_s/K*Y_o; 
end 
hold on 
plot(prediction_x,real(gaussian_process_sample)) 

set(plot(X_o,Y_o,'r.'),'MarkerSize',20) 

Учебное пособие генерирует условное моделирование с использованием метода прямого моделирования, основанного на разложении матрицы ковариации. Я понимаю, что существует несколько способов генерации условных симуляций, которые могут быть лучше, когда число симуляционных точек велико, например, обусловливание Кригинга с использованием локальной окрестности. Я нашел информацию о нескольких методах в J.-P. Chilès and P. Delfiner, «Глава 7 - Условные симуляции», в «Геостатистике: моделирование пространственной неопределенности», второе издание, John Wiley & Sons, Inc., 2012, pp. 478-628.

Есть ли существующий набор инструментов Matlab, который можно использовать для условного моделирования? Я знаю DACE, GPML и mGstat (http://mgstat.sourceforge.net/). Я считаю, что только mGstat предлагает возможность выполнять условное моделирование. Однако mGstat также, по-видимому, ограничивается только трехмерными моделями, и меня интересуют модели с более высоким размерностью.

Может ли кто-нибудь предложить какие-либо советы по началу выполнения условных симуляций с помощью существующего инструментария, такого как GPML?

================================================================================================================= ====================== EDIT

Я нашел несколько Matlab инструментарии: СТК, ScalaGauss, ooDACE

оказывается СТК способный к условному моделированию с использованием разложения ковариационной матрицы. Однако ограничивается умеренным числом (возможно, нескольких тысяч?) Точек моделирования из-за факторизации Cholesky.

Answer 1

Я использовал набор инструментов STK и я рекомендую его для других:

http://kriging.sourceforge.net/htmldoc/

Я обнаружил, что, если вам нужно условное моделирование при большом количестве точек, то вы могли бы рассмотреть генерируя условное моделирование в точках в большой конструкции эксперимента (DoE), а затем просто полагаясь на среднее предсказание, обусловленное этим DoE.