我们如何使用Matlab在以下模型中生成kappa和delta？

Question

我有一个随机模型描述了一个过程(Y)在空间和时间的演变.Ds和Dt是空间中的域(具有x和y轴的2D)和时间(具有t轴的1D).该模型通常称为混合效应模型或变异组件模型

我目前正在开发Y如下:

%# Time parameters
T=1:1:20; % input
nT=numel(T);

%# Grid and model parameters
nRow=100;
nCol=100;


[Grid.Nx,Grid.Ny,Grid.Nt] = meshgrid(1:1:nCol,1:1:nRow,T);

xPower=0.1;
tPower=1;
noisePower=1;
detConstant=1;

deterministic_mu = detConstant.*(((Grid.Nt).^tPower)./((Grid.Nx).^xPower));

beta_s = randn(nRow,nCol); % mean-zero random effect representing location specific variability common to all times

gammaTemp = randn(nT,1);

for t = 1:nT
    gamma_t(:,:,t) = repmat(gammaTemp(t),nRow,nCol); % mean-zero random effect representing time specific variability common to all locations
end

var=0.1;% noise has variance = 0.1
for t=1:nT
    kappa_st(:,:,t) = sqrt(var)*randn(nRow,nCol);
end

for t=1:nT
    Y(:,:,t) = deterministic_mu(:,:,t) + beta_s + gamma_t(:,:,t) + kappa_st(:,:,t);
end

我的问题是:

1. 如何在Y的表达式中产生delta以及kappa和delta的差异？
2. 如果我正确地生成Y,请通过使用Matlab的一些插图帮助解释一下？

如果您需要更多信息/解释,请告诉我.谢谢.

Answer 1

您的实现对 beta、gamma 和 kappa 进行采样，就好像它们是白色的一样（例如，它们在每个 (x,y,t) 处的值是独立的）。术语的描述表明情况并非如此。看起来 delta 应该捕获白噪声，而其他项捕获各自域上的相关性。例如，gamma(t_1) 和gamma(t_1+1) 之间存在非零相关性。

如果您希望将 gamma 建模为具有 gamma(t) 和 gamma(t+1) 之间的方差 var_g 和相关性 cor_g 的平稳高斯马尔可夫过程，您可以使用类似

gamma_t = nan( nT, 1 );
gamma_t(1) = sqrt(var_g)*randn();
K_g = cor_g/var_g;
K_w = sqrt( (1-K_g^2)*var_g );
for t = 2:nT,
   gamma_t(t) = K_g*gamma_t(t-1) + K_w*randn();
end
gamma_t = reshape( gamma_t, [ 1 1 nT ] );

我在上述代码中用于增益 K_g 和 K_w 的公式（以及 gamma_t(1) 的初始化）产生所需的平稳方差 sigma^2_0 和一步协方差 sigma^2_1：

请注意，上面的实现假设稍后您将使用 bsxfun 对项求和来为您执行“repmat”：

Y = bsxfun( @plus, deterministic_mu + kappa_st + delta_st, beta_s );
Y = bsxfun( @plus, Y, gamma_t );

请注意，我尚未测试上述代码，因此您应该通过采样确认它确实会产生相邻样本之间指定方差和协方差的零噪声过程。为了对 beta 进行采样，可以将相同的过程扩展到二维，但原理本质上是相同的。我怀疑 kappa 应该类似地建模为马尔可夫高斯过程，但在所有三个维度中并具有较低的方差来表示 mu、beta 和 gamma 中未捕获的高阶效应。

Delta 应该是零均值平稳白噪声。假设它是具有方差噪声功率的高斯分布，则可以使用以下方式对其进行采样

delta_st = sqrt(noisePower)*randn( [ nRows nCols nT ] );