corr_matrix.RdComputes the power exponential or Matern correlation matrix for a set of n design points in d-dimensional input region and a vector of d correlation hyper-parameters (beta).
corr_matrix(X, beta, corr = list(type = "exponential", power = 1.95))the (n x d) design matrix
a (d x 1) vector of correlation hyper-parameters in
\((-\infty, \infty)\)
a list that specifies the type of correlation function
along with the smoothness parameter. The default corresponds to power
exponential correlation with smoothness parameter "power=1.95". One
can specify a different power (between 1.0 and 2.0) for the power
exponential, or use the Matern correlation function, specified as
corr=list(type = "matern", nu=(2*k+1)/2), where \(k \in
\{0,1,2,...\}\)
The (n x n) correlation matrix, R, for the design matrix
(X) and the hyper-parameters (beta).
The power exponential correlation function is given by
\(R_{ij} = \prod_{k=1}^{d} \exp({-10^{\beta_k}|x_{ik}-x_{jk}|^{power}})\).
The Matern correlation function given by Santner, Williams, and Notz (2003) is
\(R_{ij} = \prod_{k=1}^{d} \frac{1}{\Gamma(\nu)2^{\nu-1}}(2\sqrt{\nu}|x_{ik} - \)\(x_{jk}|10^{\beta_k})^\nu \kappa_{\nu}(2\sqrt{\nu}|x_{ik} - \)\(x_{jk}|10^{\beta_k})\),
where \(\kappa_{\nu}\) is the modified Bessel function of
order \(\nu\).
Both Matern and power exponential correlation functions use the new
\(\beta\) parametrization of hyper-parameters given by \(\theta_k =
10^{\beta_k}\) for easier likelihood optimization. That is, beta is a
log scale parameter (see MacDonald et al. (2015)).
MacDonald, K.B., Ranjan, P. and Chipman, H. (2015).
GPfit: An R Package for Fitting a Gaussian Process Model to
Deterministic Simulator Outputs.
Journal of Statistical Software, 64(12), 1-23.
http://www.jstatsoft.org/v64/i12/
Ranjan, P., Haynes, R., and Karsten, R. (2011).
A Computationally Stable
Approach to Gaussian Process Interpolation of
Deterministic Computer Simulation Data,
Technometrics, 53(4), 366 - 378.
Santner, T.J., Williams, B., and Notz, W. (2003), The design and analysis of computer experiments, Springer Verlag, New York.
## 1D Example - 1
n = 5
d = 1
set.seed(3)
library(lhs)
x = maximinLHS(n,d)
beta = rnorm(1)
corr_matrix(x,beta)
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 1.0000000 0.9896462 0.9971117 0.9725082 0.9987826
#> [2,] 0.9896462 1.0000000 0.9975013 0.9954172 0.9819401
#> [3,] 0.9971117 0.9975013 1.0000000 0.9866792 0.9924232
#> [4,] 0.9725082 0.9954172 0.9866792 1.0000000 0.9609520
#> [5,] 0.9987826 0.9819401 0.9924232 0.9609520 1.0000000
## 1D Example - 2
beta = rnorm(1)
corr_matrix(x,beta,corr = list(type = "matern"))
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 1.000000e+00 1.030361e-08 2.858752e-04 4.463807e-15 1.096842e-02
#> [2,] 1.030361e-08 1.000000e+00 6.061287e-04 1.620249e-05 6.549658e-12
#> [3,] 2.858752e-04 6.061287e-04 1.000000e+00 4.870739e-10 2.600229e-07
#> [4,] 4.463807e-15 1.620249e-05 4.870739e-10 1.000000e+00 2.344908e-18
#> [5,] 1.096842e-02 6.549658e-12 2.600229e-07 2.344908e-18 1.000000e+00
## 2D example - 1
n = 10
d = 2
set.seed(2)
library(lhs)
x = maximinLHS(n,d)
beta = rnorm(2)
corr_matrix(x, beta,
corr = list(type = "exponential", power = 2))
#> [,1] [,2] [,3] [,4] [,5] [,6]
#> [1,] 1.000000000 0.014022024 0.34714595 0.89060933 0.3276537 0.0076817079
#> [2,] 0.014022024 1.000000000 0.33303062 0.03388843 0.2581943 0.9112261104
#> [3,] 0.347145948 0.333030622 1.00000000 0.53402810 0.9135157 0.2208700162
#> [4,] 0.890609331 0.033888427 0.53402810 1.00000000 0.5679417 0.0172374680
#> [5,] 0.327653733 0.258194311 0.91351568 0.56794168 1.0000000 0.1434039387
#> [6,] 0.007681708 0.911226110 0.22087002 0.01723747 0.1434039 1.0000000000
#> [7,] 0.737847362 0.061653019 0.66793462 0.95424505 0.7415377 0.0310328673
#> [8,] 0.225684376 0.452321621 0.87260109 0.32410337 0.6468200 0.3753906863
#> [9,] 0.702845281 0.001125376 0.08851206 0.55462602 0.1035367 0.0004641268
#> [10,] 0.042484708 0.896650481 0.57042147 0.09564464 0.4969425 0.6991987724
#> [,7] [,8] [,9] [,10]
#> [1,] 0.73784736 0.22568438 0.7028452806 0.042484708
#> [2,] 0.06165302 0.45232162 0.0011253760 0.896650481
#> [3,] 0.66793462 0.87260109 0.0885120621 0.570421469
#> [4,] 0.95424505 0.32410337 0.5546260208 0.095644643
#> [5,] 0.74153771 0.64682002 0.1035367021 0.496942549
#> [6,] 0.03103287 0.37539069 0.0004641268 0.699198772
#> [7,] 1.00000000 0.40594252 0.3982684429 0.160378309
#> [8,] 0.40594252 1.00000000 0.0398888440 0.633357372
#> [9,] 0.39826844 0.03988884 1.0000000000 0.005047204
#> [10,] 0.16037831 0.63335737 0.0050472037 1.000000000
## 2D example - 2
beta = rnorm(2)
R = corr_matrix(x,beta,corr = list(type = "matern", nu = 5/2))
print(R)
#> [,1] [,2] [,3] [,4] [,5]
#> [1,] 1.000000e+00 1.001810e-51 3.574760e-26 6.447775e-09 6.380465e-31
#> [2,] 1.001810e-51 1.000000e+00 1.142756e-22 1.106351e-40 1.962757e-25
#> [3,] 3.574760e-26 1.142756e-22 1.000000e+00 7.952601e-16 1.567619e-07
#> [4,] 6.447775e-09 1.106351e-40 7.952601e-16 1.000000e+00 9.375523e-20
#> [5,] 6.380465e-31 1.962757e-25 1.567619e-07 9.375523e-20 1.000000e+00
#> [6,] 1.396997e-49 7.213687e-08 3.766691e-26 4.364914e-42 2.977643e-35
#> [7,] 1.369537e-15 2.738439e-34 3.051641e-13 3.782901e-05 6.137748e-13
#> [8,] 4.240675e-25 4.111048e-24 8.674054e-10 9.628599e-25 4.657027e-19
#> [9,] 2.322358e-13 2.593713e-60 9.913060e-35 7.847741e-17 3.304532e-39
#> [10,] 7.600152e-47 9.558113e-07 1.004926e-17 1.146816e-35 4.815743e-17
#> [,6] [,7] [,8] [,9] [,10]
#> [1,] 1.396997e-49 1.369537e-15 4.240675e-25 2.322358e-13 7.600152e-47
#> [2,] 7.213687e-08 2.738439e-34 4.111048e-24 2.593713e-60 9.558113e-07
#> [3,] 3.766691e-26 3.051641e-13 8.674054e-10 9.913060e-35 1.004926e-17
#> [4,] 4.364914e-42 3.782901e-05 9.628599e-25 7.847741e-17 1.146816e-35
#> [5,] 2.977643e-35 6.137748e-13 4.657027e-19 3.304532e-39 4.815743e-17
#> [6,] 1.000000e+00 1.074790e-41 5.263661e-22 1.168921e-58 4.870246e-16
#> [7,] 1.074790e-41 1.000000e+00 7.147910e-25 1.293539e-23 7.052656e-29
#> [8,] 5.263661e-22 7.147910e-25 1.000000e+00 9.063823e-38 1.702791e-19
#> [9,] 1.168921e-58 1.293539e-23 9.063823e-38 1.000000e+00 2.393528e-55
#> [10,] 4.870246e-16 7.052656e-29 1.702791e-19 2.393528e-55 1.000000e+00