tests/testthat/test-km.R

# Intended to check that km examples results are not changing with version update
# @see http://journal.r-project.org/archive/2011-1/RJournal_2011-1_Wickham.pdf

require("DiceKriging")
require("testthat")

set.seed(1)

test_that.km <- function(model, trend.coef=NULL, covariance.sd2=NULL, covariance.range.val=NULL, covariance.nugget=NULL, covariance.eta=NULL, precision=1e-4) {

    if (!is.null(trend.coef))
        test_that(desc = "Check kriging trend", 
                  expect_true(max(abs((model@trend.coef - trend.coef)/trend.coef)) < precision))

    if (!is.null(covariance.sd2))
        test_that(desc = "Check kriging variance", 
                  expect_true(abs(model@covariance@sd2 - covariance.sd2)/covariance.sd2 < precision))

    if (!is.null(covariance.range.val))
        test_that(desc = "Check kriging range", 
                  expect_true(max(abs((model@covariance@range.val - covariance.range.val)/covariance.range.val)) < precision))

    if (!is.null(covariance.nugget))
        test_that(desc = "Check kriging nugget", 
                  expect_true(max(abs((model@covariance@nugget - covariance.nugget)/covariance.nugget)) < precision))

    if (!is.null(covariance.eta))
        test_that(desc = "Check kriging scaling", 
                  expect_true(max(abs((model@covariance@eta - covariance.eta)/covariance.eta)) < precision))
}

# example(km)

context("Checking km examples: 2D example - Branin-Hoo function")
# ----------------------------------
# A 2D example - Branin-Hoo function
# ----------------------------------

# a 16-points factorial design, and the corresponding response
d <- 2; n <- 16
design.fact <- expand.grid(x1=seq(0,1,length=4), x2=seq(0,1,length=4))
y <- apply(design.fact, 1, branin)

# kriging model 1 : matern5_2 covariance structure, no trend, no nugget effect
m1 <- km(design=design.fact, response=y,
         control=list(trace=FALSE))

test_that.km(m1,trend.coef = 306.5783,covariance.sd2 = 145556.6,covariance.range.val = c(0.8254355,2.0000000))

# kriging model 2 : matern5_2 covariance structure,
#                   linear trend + interactions, no nugget effect
m2 <- km(~.^2, design=design.fact, response=y,
         control=list(trace=FALSE))

test_that.km(m2, trend.coef = c(579.5111, -402.8916, -362.0008, 431.2314), covariance.sd2 = 87350.78, covariance.range.val = c(0.7917705,2.0000000))

m1.loo = leaveOneOut.km(m1,type="UK")

m1.loo.test = list(mean = c(286.993256592263, 61.1933200186092, 13.3103372396603, 6.71932528464657, 165.248905907798, 19.0053295402084,
							27.2225325208522, 7.83789496814171, 59.1350192499509, 36.8594545864432, 90.067237851829, 54.0537533395973,
							27.9241365199806, 97.621906142958, 202.264052547859, 153.096609124748),
				   sd = c(10.772546332554, 5.74056179816984, 5.74056179817999, 10.7725463324928, 4.30907401785216, 2.300089162552,
				   	      2.30008916257098, 4.30907401785217, 4.30907401784204,2.30008916257098, 2.30008916257098, 4.30907401783865,
				   	      10.7725463325132, 5.74056179815971, 5.74056179814956, 10.7725463325064))

test_that(desc="Test leaveOneOut",expect_that(max(abs(m1.loo.test$mean-m1.loo$mean))<1e-6, is_true()))
test_that(desc="Test leaveOneOut",expect_that(max(abs(m1.loo.test$sd-m1.loo$sd))<1e-6, is_true()))

m2.loo = leaveOneOut.km(m2,type="UK")

m2.loo.test = list(mean = c(295.42310365247, 58.7310677087344, 15.7232029614909, 0.759894555079001, 162.564378157371, 19.7760351639429,
						26.4708895401928, 9.63275693530611, 61.6951516883825, 36.1471579561365, 90.7566139900235, 52.5073354395038,
						  20.6033936837633, 99.6133571193964, 200.349297345504, 157.171599505827),
				   sd = c(9.40690817620207, 4.88889794440375, 4.8888979444158, 9.40690817630557, 3.66940746999634, 1.9439593480655,
				   	   1.9439593480619, 3.66940747003794, 3.66940746999794, 1.94395934805049, 1.94395934807647, 3.66940747003195,
				   	   9.40690817621748, 4.88889794440019, 4.88889794440971, 9.40690817632347))

test_that(desc="Test leaveOneOut",expect_that(max(abs(m2.loo.test$mean-m2.loo$mean))<1e-6, is_true()))
test_that(desc="Test leaveOneOut",expect_that(max(abs(m2.loo.test$sd-m2.loo$sd))<1e-6, is_true()))


context("Checking km examples results: 1D example with penalized MLE")
# -------------------------------
# A 1D example with penalized MLE
# -------------------------------

# from Fang K.-T., Li R. and Sudjianto A. (2006), "Design and Modeling for
# Computer Experiments", Chapman & Hall, pages 145-152

n <- 6; d <- 1
x <- seq(from=0, to=10, length=n)
y <- sin(x)
t <- seq(0,10, length=100)

# one should add a small nugget effect, to avoid numerical problems
epsilon <- 1e-3
model <- km(formula<- ~1, design=data.frame(x=x), response=data.frame(y=y),
            covtype="gauss", penalty=list(fun="SCAD", value=3), nugget=epsilon,
            control=list(trace=FALSE))

test_that.km(model, trend.coef = -0.5586176, covariance.sd2 = 3.35796, covariance.range.val = 2.417813, covariance.nugget = 0.001)

p <- predict(model, data.frame(x=t), "UK", bias.correct=TRUE)

p.test = list(mean = c(-3.33066907387547e-16, 0.103668583556765,
                    0.204106764415468, 0.301735270357276, 0.395756768311479, 0.485384515282419,
                    0.569851151220041, 0.648417561574452, 0.720381694694425, 0.785087216109229,
                    0.841931880366164, 0.890375501555915, 0.929947405986085, 0.960253254668741,
                    0.980981129345181, 0.991906783620992, 0.992897970333583, 0.983917767397483,
                    0.965026836913487, 0.936384566102604, 0.898249053426139, 0.850975918849646,
                    0.795015933354937, 0.73091147924265, 0.659291869232534, 0.580867568592383,
                    0.496423380244189, 0.406810667751401, 0.312938705039622, 0.215765254418146,
                    0.116286485745949, 0.0155263592432038, -0.0854743976639361, -0.185669583082299,
                    -0.284018942575723, -0.379499577893206, -0.471117200190948, -0.557917087546549,
                    -0.638994610661276, -0.713505196652745, -0.780673608644263, -0.83980242832094,
                    -0.890279639583468, -0.931585223704922, -0.963296689783487, -0.985093478570262,
                    -0.996760192713487, -0.998188621869706, -0.989378546759103, -0.970437321859351,
                    -0.941578251818738, -0.903117791614111, -0.855471614783334, -0.799149607542319,
                    -0.73474985908724, -0.662951729736852, -0.584508088662155, -0.500236821676946,
                    -0.411011716841405, -0.317752841402879, -0.221416527826267, -0.122985089336235,
                    -0.0234563865111159, 0.0761666339406344, 0.174886308980496, 0.271720421685087,
                    0.365712045055562, 0.455938961810148, 0.541522561823235, 0.621636126948619,
                    0.69551242194443, 0.76245051997177, 0.821821801531807, 0.873075076596748,
                    0.915740790936399, 0.9494342891014, 0.97385811805448, 0.988803366903467,
                    0.994150049448424, 0.989866547180333, 0.976008140836365, 0.952714668511005,
                    0.920207357536033, 0.878784885778434, 0.828818735577156, 0.770747910172976,
                    0.705073088118135, 0.632350295735319, 0.553184181193531, 0.468220976160377,
                    0.378141232269578, 0.283652419816319, 0.185481475183269, 0.0843673815421025,
                    -0.0189461355818289, -0.123717720189388, -0.229215232199558,
                    -0.334722260961891, -0.439544244726991, -0.544021110889371),
           sd = c(1.37722860163619e-16, 0.055938047354605, 0.0717239075944272,
                  0.0885333491721034, 0.103479571632164, 0.11558566074085,
                  0.124535231655531, 0.130284178908033, 0.132926953513156,
                  0.132643474051053, 0.129675361680544, 0.124314691041382,
                  0.116899579170262, 0.107815385655198, 0.0975029505547555,
                  0.086477898723425, 0.0753678296196528, 0.0649731906768246,
                  0.0563321548096111, 0.0506662341718909, 0.0489421469571659,
                  0.0511448312539402, 0.0561923009414814, 0.0626840730368677,
                  0.0694846347505854, 0.075822706651813, 0.0812033992665197,
                  0.0853179155081378, 0.0879856320787734, 0.0891211458796282,
                  0.0887162707338793, 0.0868307797971652, 0.0835888470871021,
                  0.0791801846304929, 0.073866202141622, 0.067992148964746,
                  0.0620049698687374, 0.056469588561222, 0.052056108560131,
                  0.0494397154112441, 0.0490738072493008, 0.0509543226071206,
                  0.0546149617347237, 0.0593642281895726, 0.064524637090391,
                  0.0695406284760611, 0.0739914572527713, 0.077571596936056,
                  0.0800691406401005, 0.0813503165577205, 0.0813503165577172,
                  0.0800691406401005, 0.0775715969360595, 0.0739914572527677,
                  0.0695406284760611, 0.064524637090391, 0.059364228189577,
                  0.0546149617347237, 0.050954322607131, 0.0490738072493116,
                  0.0494397154112548, 0.0520561085601361, 0.0564695885612173,
                  0.0620049698687245, 0.067992148964746, 0.0738662021416293,
                  0.0791801846304862, 0.0835888470871022, 0.0868307797971652,
                  0.0887162707338762, 0.0891211458796311, 0.0879856320787612,
                  0.0853179155081378, 0.0812033992665066, 0.07582270665182,
                  0.0694846347505701, 0.0626840730368719, 0.0561923009414766,
                  0.0511448312539402, 0.0489421469571659, 0.0506662341718909,
                  0.0563321548096158, 0.0649731906768328, 0.0753678296196563,
                  0.0864778987234313, 0.097502950554761, 0.107815385655201,
                  0.116899579170262, 0.12431469104138, 0.129675361680544, 0.132643474051059,
                  0.13292695351316, 0.130284178908038, 0.124535231655527, 0.115585660740852,
                  0.103479571632161, 0.0885333491720974, 0.071723907594431,
                  0.055938047354605, 0))
test_that(desc="Test predict",expect_equal(p$mean,p.test$mean))
test_that(desc="Test predict",expect_equal(p$sd,p.test$sd))

context("Checking km examples results: 1D example with known trend and known or unknown covariance parameters")
# ------------------------------------------------------------------------
# A 1D example with known trend and known or unknown covariance parameters
# ------------------------------------------------------------------------

x <- c(0, 0.4, 0.6, 0.8, 1);
y <- c(-0.3, 0, -0.8, 0.5, 0.9)

theta <- 0.01; sigma <- 3; trend <- c(-1,2)

model <- km(~x, design=data.frame(x=x), response=data.frame(y=y),
            covtype="matern5_2", coef.trend=trend, coef.cov=theta,
            coef.var=sigma^2,
            control=list(trace=FALSE))

# below: if you want to specify trend only, and estimate both theta and sigma:
# model <- km(~x, design=data.frame(x=x), response=data.frame(y=y),
#             covtype="matern5_2", coef.trend=trend, lower=0.2)
# Remark: a lower bound or penalty function is useful here,
#         due to the very small number of design points...

# kriging with gaussian covariance C(x,y)=sigma^2 * exp(-[(x-y)/theta]^2),
#         and linear trend t(x) = -1 + 2x

t <- seq(from=0, to=1, by=0.005)
p <- predict(model, newdata=data.frame(x=t), type="SK")
# beware that type = "SK" for known parameters (default is "UK")

p.test=list(mean = c(-0.3, -0.409945600307312, -0.613204123817726,
                     -0.771785710062141, -0.862937846603047, -0.905542849815739, -0.920593604659762,
                     -0.921829914451626, -0.916656040817311, -0.908661891267655, -0.899474346347788,
                     -0.889796659403621, -0.879922358335064, -0.869970682858201, -0.859989036423894,
                     -0.8499959345099, -0.839998503637175, -0.829999452876737, -0.819999801133533,
                     -0.809999928100776, -0.799999974130126, -0.789999990732495, -0.779999996693293,
                     -0.769999998824458, -0.759999999583497, -0.749999999852888, -0.739999999948188,
                     -0.729999999981801, -0.719999999993623, -0.709999999997771, -0.699999999999222,
                     -0.689999999999729, -0.679999999999906, -0.669999999999967, -0.659999999999989,
                     -0.649999999999996, -0.639999999999999, -0.63, -0.62, -0.61,
                     -0.6, -0.59, -0.58, -0.57, -0.56, -0.549999999999999, -0.539999999999997,
                     -0.529999999999991, -0.519999999999973, -0.509999999999923, -0.499999999999778,
                     -0.489999999999363, -0.479999999998178, -0.4699999999948, -0.459999999985197,
                     -0.449999999957968, -0.439999999880999, -0.429999999664131, -0.419999999055227,
                     -0.409999997352141, -0.399999992608608, -0.389999979457365, -0.37999994318101,
                     -0.369999843679068, -0.359999572467764, -0.3499988384314, -0.339996867549684,
                     -0.329991623673772, -0.319977816667161, -0.309941902686749, -0.299849813242225,
                     -0.28961768321933, -0.27904458309066, -0.267665689843322, -0.254455315617075,
                     -0.237297957090211, -0.212267956172299, -0.17336734573204, -0.115201178233636,
                     -0.044270171516375, 0, -0.024270171516375, -0.0752011782336362,
                     -0.113367345732041, -0.132267956172301, -0.137297957090217, -0.134455315617091,
                     -0.127665689843368, -0.119044583090795, -0.109617683219717, -0.099849813243336,
                     -0.0899419026899333, -0.0799778166762707, -0.0699916236997703,
                     -0.0599968676237009, -0.04999883864156, -0.0399995730627685,
                     -0.0299998453584138, -0.019999947904877, -0.00999999269665746,
                     -2.95655697764231e-08, 0.00999989993468188, 0.0199997168498212,
                     0.0299992187312075, 0.0399978624578226, 0.0499941921990326, 0.0599843377632238,
                     0.069958118374058, 0.0798890833376284, 0.0897095134343814, 0.0992490662113485,
                     0.108088416096727, 0.115222915453328, 0.118328449216618, 0.112276578085377,
                     0.0864897854510566, 0.0213397808614965, -0.113163271339799, -0.34399410883182,
                     -0.638649142418125, -0.8, -0.618649142418125, -0.30399410883182,
                     -0.0531632713397988, 0.101339780861496, 0.186489785451056, 0.232276578085373,
                     0.258328449216605, 0.275222915453288, 0.288088416096611, 0.299249066211015,
                     0.309709513433426, 0.319889083334896, 0.329958118366258, 0.339984337741019,
                     0.349994192135985, 0.359997862279321, 0.369999218227404, 0.379999715432661,
                     0.389999895962894, 0.399999959347342, 0.409999976489389, 0.419999966866637,
                     0.429999920160188, 0.439999785638878, 0.44999941900554, 0.459998433700825,
                     0.469995811810887, 0.479988908324471, 0.48997095134019, 0.499924906620002,
                     0.509808841609278, 0.519522291545196, 0.528832844921614, 0.537227657808521,
                     0.5436489785451, 0.546133978086148, 0.541683672866019, 0.527600589116818,
                     0.507135085758187, 0.5, 0.527135085758188, 0.567600589116818,
                     0.60168367286602, 0.62613397808615, 0.643648978545105, 0.657227657808536,
                     0.668832844921656, 0.679522291545317, 0.689808841609626, 0.699924906621001,
                     0.709970951343056, 0.71998890833267, 0.729995811834286, 0.73999843376744,
                     0.749999419194684, 0.759999786174382, 0.769999921671599, 0.779999971118118,
                     0.789999988404753, 0.799999992608608, 0.809999988404753, 0.819999971118118,
                     0.829999921671599, 0.839999786174382, 0.849999419194684, 0.85999843376744,
                     0.869995811834286, 0.87998890833267, 0.889970951343056, 0.899924906621002,
                     0.909808841609626, 0.919522291545317, 0.928832844921656, 0.937227657808536,
                     0.943648978545105, 0.946133978086149, 0.94168367286602, 0.927600589116818,
                     0.907135085758187, 0.9),
            sd = c(0, 1.6793050315316, 2.55516566296315, 2.87721515634081,
                   2.97102004245279, 2.99394356557226, 2.99884689620796, 2.99979565539401,
                   2.99996576899956, 2.99999451876546, 2.99999915414755, 2.99999987342633,
                   2.99999998154624, 2.99999999736889, 2.99999999963204, 2.9999999999494,
                   2.99999999999315, 2.99999999999908, 2.99999999999988, 2.99999999999998,
                   3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3,
                   3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 2.99999999999998,
                   2.99999999999988, 2.99999999999908, 2.99999999999315, 2.9999999999494,
                   2.99999999963204, 2.99999999736889, 2.99999998154624, 2.99999987342633,
                   2.99999915414755, 2.99999451876546, 2.99996576899956, 2.99979565539401,
                   2.99884689620796, 2.99394356557226, 2.97102004245279, 2.87721515634081,
                   2.55516566296315, 1.6793050315316, 0, 1.6793050315316, 2.55516566296315,
                   2.87721515634081, 2.97102004245279, 2.99394356557226, 2.99884689620796,
                   2.99979565539401, 2.99996576899956, 2.99999451876546, 2.99999915414755,
                   2.99999987342633, 2.99999998154624, 2.99999999736889, 2.99999999963204,
                   2.9999999999494, 2.99999999999315, 2.99999999999908, 2.99999999999988,
                   2.99999999999998, 3, 2.99999999999998, 2.99999999999988, 2.99999999999908,
                   2.99999999999315, 2.9999999999494, 2.99999999963204, 2.99999999736889,
                   2.99999998154624, 2.99999987342633, 2.99999915414755, 2.99999451876546,
                   2.99996576899956, 2.99979565539401, 2.99884689620796, 2.99394356557226,
                   2.97102004245279, 2.87721515634081, 2.55516566296315, 1.6793050315316,
                   0, 1.6793050315316, 2.55516566296315, 2.87721515634081, 2.97102004245279,
                   2.99394356557226, 2.99884689620796, 2.99979565539401, 2.99996576899956,
                   2.99999451876546, 2.99999915414755, 2.99999987342633, 2.99999998154624,
                   2.99999999736889, 2.99999999963204, 2.9999999999494, 2.99999999999315,
                   2.99999999999908, 2.99999999999988, 2.99999999999998, 3, 2.99999999999998,
                   2.99999999999988, 2.99999999999908, 2.99999999999315, 2.9999999999494,
                   2.99999999963204, 2.99999999736889, 2.99999998154624, 2.99999987342633,
                   2.99999915414755, 2.99999451876546, 2.99996576899956, 2.99979565539401,
                   2.99884689620796, 2.99394356557226, 2.97102004245279, 2.87721515634081,
                   2.55516566296315, 1.6793050315316, 0, 1.6793050315316, 2.55516566296315,
                   2.87721515634081, 2.97102004245279, 2.99394356557226, 2.99884689620796,
                   2.99979565539401, 2.99996576899956, 2.99999451876546, 2.99999915414755,
                   2.99999987342633, 2.99999998154624, 2.99999999736889, 2.99999999963204,
                   2.9999999999494, 2.99999999999315, 2.99999999999908, 2.99999999999988,
                   2.99999999999998, 3, 2.99999999999998, 2.99999999999988, 2.99999999999908,
                   2.99999999999315, 2.9999999999494, 2.99999999963204, 2.99999999736889,
                   2.99999998154624, 2.99999987342633, 2.99999915414755, 2.99999451876546,
                   2.99996576899956, 2.99979565539401, 2.99884689620796, 2.99394356557226,
                   2.97102004245279, 2.87721515634081, 2.55516566296315, 1.6793050315316,  0))
test_that(desc="Test predict",expect_equal(p$mean,p.test$mean))
test_that(desc="Test predict",expect_equal(p$sd,p.test$sd))


context("Checking km examples results: Kriging with noisy observations (heterogeneous noise variance)")
# --------------------------------------------------------------
# Kriging with noisy observations (heterogeneous noise variance)
# --------------------------------------------------------------

fundet <- function(x){
    return((sin(10*x)/(1+x)+2*cos(5*x)*x^3+0.841)/1.6)
}

level <- 0.5; epsilon <- 0.1
theta <- 1/sqrt(30); p <- 2; n <- 10
x <- seq(0,1, length=n)

# Heteregeneous noise variances: number of Monte Carlo evaluation among
#                                a total budget of 1000 stochastic simulations
MC_numbers <- c(10,50,50,290,25,75,300,10,40,150)
noise.var <- 3/MC_numbers

# Making noisy observations from 'fundet' function (defined above)
y <- fundet(x) + noise.var*rnorm(length(x))

# kriging model definition (no estimation here)
model <- km(y~1, design=data.frame(x=x), response=data.frame(y=y),
            covtype="gauss", coef.trend=0, coef.cov=theta, coef.var=1,
            noise.var=noise.var,
            control=list(trace=FALSE))

# prediction
t <- seq(0,1,by=0.01)
p <- predict.km(model, newdata=data.frame(x=t), type="SK")

p.test = list(mean=c(0.539334747148277, 0.586647044674522, 0.633654722392378, 0.679929639445603,
                     0.725032146520035, 0.76851701239259, 0.809939706096152, 0.848862926253087,
                     0.884863259581319, 0.917537843577757, 0.946510904204198, 0.971440038235026,
                     0.992022111884504, 1.00799865243506, 1.01916061776519, 1.02535243976294,
                     1.02647525135751, 1.02248922297386, 1.01341495221093, 0.999333869999283,
                     0.980387646904218, 0.956776604070301, 0.9287571540107, 0.896638316492772,
                     0.860777373649944, 0.821574745689538, 0.779468183753871, 0.734926389283163,
                     0.688442179356868, 0.640525324773073, 0.59169519196919, 0.542473321284192,
                     0.493376072589807, 0.444907465129216, 0.397552331718294, 0.351769898564123,
                     0.307987891159104, 0.266597254366536, 0.227947561289564, 0.19234317117567,
                     0.16004018180712, 0.131244206895013, 0.106108994229466, 0.0847358860015263,
                     0.067174109021187, 0.0534218696840629, 0.0434282166177865, 0.0370956230593633,
                     0.0342832312338992, 0.0348106923531448, 0.0384625283391518, 0.04499293500192,
                     0.0541309411535644, 0.0655858340206116, 0.0790527583232632, 0.0942183945394055,
                     0.110766621187819, 0.128384066487297, 0.146765456523336, 0.1656186701334,
                     0.184669415155653, 0.203665446515514, 0.222380253873431, 0.24061615522411,
                     0.258206742887227, 0.275018639686927, 0.290952535661456, 0.305943489205345,
                     0.319960490905191, 0.333005303219061, 0.345110604257797, 0.356337478906972,
                     0.366772315006854, 0.376523175895532, 0.385715732926623, 0.394488852219889,
                     0.402989938540434, 0.411370145522012, 0.419779565200525, 0.428362510821178,
                     0.437253005021073, 0.446570580747562, 0.456416494719598, 0.46687044303276,
                     0.477987855893409, 0.489797833769091, 0.502301770857379, 0.51547269415915,
                     0.52925532809377, 0.543566876039111, 0.558298491954252, 0.573317397873672,
                     0.588469587048248, 0.603583038305923, 0.618471355209542, 0.632937734124466,
                     0.64677915861476, 0.659790713814218, 0.671769913623502, 0.682520935728924,
                     0.691858664392208),
              sd=c(0.377029458451379, 0.353477635741769, 0.331048780750455, 0.309901789954737,
                   0.290194209894752, 0.272078140711739, 0.255693335726301, 0.241156971388702,
                   0.228550243358871, 0.217903121026768, 0.209180023038072, 0.202270201518728,
                   0.196986396115609, 0.193073449055664, 0.190225693933524, 0.188109475678877,
                   0.186386324483874, 0.184733202947482, 0.182858030705805, 0.180510369041955,
                   0.17748818331955, 0.173641976171823, 0.168877543352864, 0.163158408232265,
                   0.156508809270597, 0.149018011275447, 0.140846648943365, 0.132235611420845,
                   0.123517178517844, 0.11512572777864, 0.107599737842161, 0.101557254845469,
                   0.0976194088289012, 0.0962714353566698, 0.0977066069260894, 0.101752687919136,
                   0.107939608904669, 0.115653145670172, 0.124274628914087, 0.133258903590641,
                   0.142159563724345, 0.150626545790194, 0.158394033615692, 0.165267084194709,
                   0.171109726773692, 0.175834882207011, 0.179395709724203, 0.181777922609756,
                   0.18299276787915, 0.183070560424202, 0.18205483806674, 0.179997341432579,
                   0.176954116156102, 0.172983084635437, 0.168143445517905, 0.162497243795756,
                   0.156113433919499, 0.149074759545596, 0.141487817200238, 0.133496739715259,
                   0.125300892938463, 0.117176389739067, 0.109499017316081, 0.102760375180372,
                   0.0975580180680758, 0.0945294254198149, 0.094211568198508, 0.0968685832496973,
                   0.102397003846006, 0.110380425208506, 0.12023962124095, 0.131371888401085,
                   0.143228876118868, 0.155342329216061, 0.167323289645481, 0.178852522758449,
                   0.189670189557165, 0.199567202683516, 0.20837840113671, 0.215977000068187,
                   0.22226973231672, 0.22719226116034, 0.230704645851273, 0.232786831193859,
                   0.233434297578214, 0.232654148383716, 0.230462031132744, 0.226880394359816,
                   0.221938688196329, 0.21567625058695, 0.208148830747982, 0.199440062098813,
                   0.189679808332715, 0.179072247643501, 0.167937686154268, 0.156772333879503,
                   0.146325713173994, 0.137675771860785, 0.132232680285778, 0.13154331895153,
                   0.136834814016462))
test_that(desc="Test predict",expect_equal(p$mean,p.test$mean))
test_that(desc="Test predict",expect_equal(p$sd,p.test$sd))

context("Checking km examples results: Kriging with non-linear scaling on Xiong et al.'s function")
# ----------------------------------------------------------
# Kriging with non-linear scaling on Xiong et al.'s function
# ----------------------------------------------------------

f11_xiong <- function(x){
    return( sin(30*(x - 0.9)^4)*cos(2*(x - 0.9)) + (x - 0.9)/2)
}

t <- seq(0,1,,300)
f <- f11_xiong(t)

doe <- data.frame(x=seq(0,1,,20))
resp <- f11_xiong(doe)

knots <- list(  c(0,0.5,1) )
m <- km(design=doe, response=resp, scaling=TRUE, gr=TRUE,
        knots=knots, covtype="matern5_2",  coef.var=1, coef.trend=0,
        control=list(trace=FALSE))

test_that.km(m, covariance.eta = c(17.6113829, 2.4169448,0.8873958), precision=1e-4)

p <- predict(m, data.frame(x=t), "UK", bias.correct=TRUE)

p.test = list(mean=c(-0.6181866493757, -0.619384327994964, -0.617544523595288, -0.612923071362958,
                     -0.605845139548477, -0.596671063973597, -0.5857719188513, -0.573512506951456,
                     -0.560239960336162, -0.546276545346517, -0.531915583696315, -0.517419652595045,
                     -0.503020424530901, -0.48891966269035, -0.475291009762516, -0.462282303092248,
                     -0.450018223374948, -0.438603069295669, -0.42812274959666, -0.418645515567079,
                     -0.410222490091608, -0.402888224309673, -0.396661205967063, -0.391544260046136,
                     -0.387524798637093, -0.384574889326519, -0.382651120532267, -0.381694248922229,
                     -0.381628618885048, -0.382361347422901, -0.383781270163437, -0.385757645719527,
                     -0.388138619983493, -0.390750112138151, -0.39339879496282, -0.395878591204417,
                     -0.397976974292908, -0.399480772708841, -0.400181725829701, -0.399882006665258,
                     -0.39839989913914, -0.395575795258101, -0.39127865957273, -0.3854130938907,
                     -0.377927123515424, -0.368820816728299, -0.358155841284965, -0.346066054910513,
                     -0.33276910000895, -0.318572838493337, -0.303859639167353, -0.289060874908515,
                     -0.274635057702686, -0.261049657449351, -0.24876576473141, -0.238224905538517,
                     -0.229837438656032, -0.223972068333405, -0.220946089591973, -0.221016054175952,
                     -0.224368604223206, -0.231111270339583, -0.24126307264834, -0.254744799000975,
                     -0.271369340122068, -0.290841799346952, -0.312780978672547, -0.336744665077429,
                     -0.362250734912443, -0.388794612006908, -0.415863768733658, -0.44294984162595,
                     -0.46955883415162, -0.495219795913324, -0.519492297338891, -0.541972959729411,
                     -0.562301250572315, -0.58016471181062, -0.595303753045532, -0.607516110964843,
                     -0.616660528928325, -0.622654732051651, -0.625466489130159, -0.62510481364984,
                     -0.621612795849773, -0.615061819013151, -0.605546931866299, -0.593183182715417,
                     -0.578102749616398, -0.560452725215194, -0.54039343557363, -0.518097189866798,
                     -0.493747372773389, -0.467537804093884, -0.439672300960193, -0.410364383850582,
                     -0.379836606390905, -0.348317278525024, -0.316035481422044, -0.283216733342154,
                     -0.250079473593566, -0.216832238706492, -0.18367142428487, -0.150779543024492,
                     -0.118323903822031, -0.0864556491104871, -0.0553090978625722,
                     -0.025001350370222, 0.00436788182654946, 0.0327162514325693,
                     0.0599785818530151, 0.0861065136495609, 0.111067500215288, 0.134842073141954,
                     0.157420631276143, 0.17880074278872, 0.198984948559774, 0.217978999387435,
                     0.235790468125003, 0.252427685264059, 0.267898952858084, 0.282211997182547,
                     0.295373625258048, 0.3073895544446, 0.318264387829763, 0.328001711173836,
                     0.336604289799667, 0.34407438300404, 0.350414777597394, 0.355630963113317,
                     0.359733138904088, 0.362737619456597, 0.364667717790052, 0.365554190300026,
                     0.365435316966647, 0.364356682516446, 0.36237071675713, 0.35953604575857,
                     0.355916699693483, 0.351581217897373, 0.346601686947497, 0.341052743241555,
                     0.335010567551571, 0.328551867513987, 0.321752661075919, 0.314686950737844,
                     0.307425655285841, 0.300035827577346, 0.292580112854382, 0.285116407162208,
                     0.27766204250694, 0.270041862952395, 0.262242107025742, 0.254283613101952,
                     0.246187224420236, 0.237973694135407, 0.22966359627162, 0.22127724230822,
                     0.212834603129049, 0.204355235204216, 0.195858207012237, 0.187362026833072,
                     0.178884575211622, 0.17044304226518, 0.162053869619148, 0.153732696771108,
                     0.145494311686322, 0.137352605442916, 0.129320530751226, 0.121410064183481,
                     0.113632171952961, 0.105996779097326, 0.0985127419209676, 0.0911878235616712,
                     0.084028673156181, 0.0770408295665373, 0.0702287931929301, 0.0635961260874189,
                     0.0571455481314747, 0.0508790277534277, 0.0447978674345578, 0.0389027842294168,
                     0.0331939855252174, 0.0276712402481287, 0.0223339457174205, 0.0171811903426216,
                     0.012211812341782, 0.00742445466112723, 0.00281761625927527,
                     -0.00161030008117329, -0.00586094412986754, -0.00993599177390887,
                     -0.0138371466425015, -0.0175661490456186, -0.0211247843624052,
                     -0.0245148907297168, -0.0277383660603413, -0.0307971744181088,
                     -0.0336933517779194, -0.0364290111936353, -0.0390063474041107,
                     -0.0414276408938766, -0.043695261437668, -0.0458116711464816,
                     -0.0477794270378376, -0.0496011831495988, -0.0512796917672645,
                     -0.0528178004990926, -0.0542184449995801, -0.0554846412842739,
                     -0.0566194783879771, -0.0576261113543447, -0.0585077545441052,
                     -0.0592676752434355, -0.0599091875622562, -0.0604356466079271,
                     -0.0608504429197739, -0.0611569971574735, -0.0613587550243052,
                     -0.0614591824230813, -0.0614617608275543, -0.0613699828459608,
                     -0.0611873466596634, -0.0609173445435782, -0.0605634478750595,
                     -0.0601290929642404, -0.0596176680799324, -0.0590325016143571,
                     -0.0583768513313831, -0.0576538946484541, -0.0568667199002712,
                     -0.056018318538126, -0.0551115782169636, -0.0541492767282164,
                     -0.053134076733538, -0.0520685212602773, -0.0509550299185329,
                     -0.0497958959771331, -0.0485932893203502, -0.0473492739893414,
                     -0.0460658305532157, -0.0447448768083503, -0.0433882866015425,
                     -0.0419979068691594, -0.0405755729807539, -0.0391231224764197,
                     -0.037642407279874, -0.0361353044671106, -0.0346037256676119,
                     -0.0330496251735891, -0.0314750068279022, -0.029881929758972,
                     -0.0282725130298932, -0.0266489379356994, -0.0250134308116316,
                     -0.0233682159255298, -0.0217154651520595, -0.0200572529111883,
                     -0.0183955162365097, -0.0167320196916362, -0.0150683248582279,
                     -0.0134057641334375, -0.0117454185846687, -0.0100880996189366,
                     -0.00843433422987598, -0.00678435360090385, -0.00513808484810206,
                     -0.00349514569129573, -0.00185484187708767, -0.00021617168983519,
                     0.00142213417059717, 0.00306154084187579, 0.00470363940656769,
                     0.00635008412999384, 0.00800253644346121, 0.0096626152789829,
                     0.0113318533696893, 0.0130116591535249, 0.0147032839285056, 0.0164077939210321,
                     0.01812604694847, 0.0198586733609143, 0.0216060609695321, 0.0233683436753124,
                     0.0251453936053517, 0.026936821616966, 0.0287420016740626, 0.0305601077090647,
                     0.0323901483514226, 0.0342309984794159, 0.0360814277847514, 0.0379401265308031,
                     0.039805728677887, 0.0416768325421688, 0.0435520191471147, 0.0454298684185662,
                     0.0473089733708269, 0.0491879524196559, 0.0510654599590658, 0.052940195323104
                     ),
              sd=c(5.40521927259871e-09, 0.0278955268448099, 0.0527585896172703,
                   0.0739907998166573, 0.0911593329653399, 0.103979113019168, 0.112298629160489,
                   0.116088760303684, 0.115434249987665, 0.110527531524034, 0.101664573762998,
                   0.0892423550839467, 0.0737575014939018, 0.0558055602221966, 0.0360803469228653,
                   0.0153728564251457, 0.00543056383264634, 0.0253804824448798,
                   0.0436243232213521, 0.0594847396849461, 0.0724484332714464, 0.0821529857215723,
                   0.0883749441803885, 0.0910191504001311, 0.090109354682045, 0.0857800864787522,
                   0.0782696560953018, 0.0679140679135962, 0.0551415471531945, 0.0404673284002107,
                   0.0244883430773698, 0.00787751923542965, 0.00862301047745113,
                   0.0242523811374775, 0.0383578795645493, 0.0504259263091857, 0.060073235236072,
                   0.0670380051926555, 0.0711714614319289, 0.0724298026008374, 0.0708666151632932,
                   0.0666257608974258, 0.0599346743115707, 0.0510979442963852, 0.0404910063063968,
                   0.0285537455708464, 0.0157838237623598, 0.00272962748301085,
                   0.0100194860268291, 0.0218917845732854, 0.0324154436456769, 0.041226072210972,
                   0.0480604358327453, 0.0527500887865998, 0.0552149089868492, 0.0554566080881815,
                   0.0535522815305141, 0.0496480244197844, 0.0439525957426021, 0.036731076282206,
                   0.0282984402728045, 0.0190129534423633, 0.00926933292299627,
                   0.000508314633666543, 0.00987887094446615, 0.0184394549981615,
                   0.025868805193392, 0.0319263863521959, 0.0364476944644532, 0.0393392568502533,
                   0.0405733622311219, 0.0401825961301444, 0.0382542477369581, 0.0349246297860518,
                   0.0303733273987832, 0.0248173713629215, 0.0185053195971253, 0.0117112320949676,
                   0.00472854851496484, 0.00213608466720373, 0.00857703536284724,
                   0.0143320779778923, 0.01920010213466, 0.0230382212933763, 0.0257579258195494,
                   0.0273209044938779, 0.027734593873248, 0.0270475309664019, 0.0253445760290491,
                   0.0227420566204376, 0.0193828688312309, 0.0154315601428472, 0.0110694129370792,
                   0.00648955128236172, 0.00189211244663229, 0.00252062580429939,
                   0.00655821986761304, 0.0100680349620087, 0.0129397490195468,
                   0.0151023575824548, 0.0165208562061587, 0.0171926877066967, 0.0171440201894679,
                   0.0164259252930565, 0.0151105180215727, 0.0132871086353029, 0.0110584070205879,
                   0.00853681248148238, 0.00584081822587245, 0.0030915617872107,
                   0.000409564033737249, 0.00208876593896904, 0.00430216449017573,
                   0.00615645492601835, 0.00760391452504312, 0.00862060533720077,
                   0.00920357645564576, 0.00936798525194681, 0.00914419573867131,
                   0.00857490906518254, 0.00771237319010031, 0.00661571005253904,
                   0.00534839078003959, 0.00397588338780098, 0.00256349360766343,
                   0.00117441883154567, 0.000131960472590607, 0.00130221158917693,
                   0.00229518025109844, 0.00308525654371616, 0.00366061098970859,
                   0.00402110530459799, 0.00417623390728292, 0.00414313350870577,
                   0.00394469727146971, 0.00360782388655352, 0.00316182448568462,
                   0.00263700322359383, 0.00206342123657317, 0.00146984871518563,
                   0.000882906106068735, 0.000326393111302364, 0.00017919945124684,
                   0.000617502861537928, 0.000978176647015097, 0.00125681052512315,
                   0.00145375833407683, 0.00157304170827566, 0.00162136080526537,
                   0.00160722109959803, 0.00153974554393112, 0.00142436837949148,
                   0.00126633761399362, 0.00107215305790654, 0.000849243596628756,
                   0.000605830445054434, 0.000350781121539972, 9.34557230022058e-05,
                   0.000156461894606401, 0.000389530359206717, 0.000597578015332369,
                   0.000773965170886398, 0.000913601489404074, 0.00101294403454977,
                   0.00106997595903501, 0.00108416608705489, 0.00105641026773893,
                   0.000988955608436533, 0.000885308725116187, 0.00075012910692986,
                   0.000589108661666367, 0.000408838589325305, 0.000216665045473663,
                   2.05358786258385e-05, 0.000171195663440373, 0.00035049470042588,
                   0.000510462170714594, 0.00064547147590016, 0.000751196250130842,
                   0.00082461791854646, 0.000864013024818227, 0.000868921373466751,
                   0.000840096400322934, 0.000779439245878901, 0.000689917932191898,
                   0.000575472935835753, 0.000440910388231675, 0.000291784190593147,
                   0.00013426865201414, 2.49759352367094e-05, 0.000179029580407467,
                   0.000321463534618553, 0.000446854479379327, 0.000550845208077309,
                   0.000630162898882065, 0.000682619492371508, 0.000707094717834278,
                   0.00070350298582261, 0.000672745558617544, 0.000616649392677178,
                   0.000537893952064388, 0.000439927185715771, 0.000326871812311767,
                   0.000203423136934822, 7.47399457217233e-05, 5.36703955013989e-05,
                   0.000176223765368053, 0.000287924529942511, 0.000384636393809221,
                   0.000463108413248862, 0.000520981923362827, 0.00055678272003298,
                   0.000569899274228118, 0.000560548166050523, 0.0005297280036614,
                   0.000479163051962901, 0.00041123769638931, 0.000328922784868723,
                   0.000235694860014172, 0.000135449388533949, 3.24094344400267e-05,
                   6.89750596878059e-05, 0.000164340771316577, 0.000249916380224616,
                   0.000322641253834587, 0.000380175530727422, 0.000420898010528653,
                   0.000443892364955128, 0.000448922521742606, 0.00043639830101167,
                   0.000407332408300389, 0.000363289826996956, 0.000306330563905879,
                   0.000238946643657319, 0.000163994227643533, 8.46218434857493e-05,
                   4.19605289148908e-06, 7.37936649095409e-05, 0.000146060736710977,
                   0.000209846922736215, 0.000262963820935719, 0.000303793912147134,
                   0.000331282275677673, 0.000344919445923973, 0.000344716229990787,
                   0.000331171444458688, 0.000305233497624261, 0.000268256691983734,
                   0.000221953047406635, 0.000168340390248118, 0.000109687473029649,
                   4.84569757950173e-05, 1.27524710410294e-05, 7.13034155784107e-05,
                   0.000124811915473892, 0.000171320676715324, 0.00020930911251436,
                   0.000237689136521643, 0.000255794183156826, 0.000263361980089882,
                   0.00026051181346646, 0.000247717087779459, 0.000225773936214625,
                   0.000195766598363813, 0.000159030220130916, 0.000117111705794583,
                   7.17292720886261e-05, 2.47314354002609e-05, 2.19442339162618e-05,
                   6.63678680666483e-05, 0.000106831375282275, 0.000141922537176468,
                   0.000170527257494718, 0.000191826601721701, 0.000205289795982779,
                   0.000210663613661126, 0.000207958686742575, 0.000197433290365047,
                   0.00017957513853237, 0.000155081669579709, 0.000124839293810563,
                   8.99020365281582e-05, 5.14700095128775e-05, 1.08682044601612e-05,
                   3.04760387232468e-05, 7.11069314639068e-05, 0.000109651244740311,
                   0.000144843205375932, 0.000175529086439696, 0.000200670149077007,
                   0.000219344018510228, 0.000230744672357331, 0.000234181207482184,
                   0.000229075585637918, 0.000214959535503335, 0.000191470767953811,
                   0.000158348663688514, 0.000115429547853075, 6.26416716073123e-05,
                   5.4052192726006e-09))
test_that(desc="Test predict",expect_that(max(abs(p$mean-p.test$mean))<1e-6, is_true()))
test_that(desc="Test predict",expect_that(max(abs(p$sd-p.test$sd))<1e-6, is_true()))

Try the DiceKriging package in your browser

Any scripts or data that you put into this service are public.

DiceKriging documentation built on Feb. 24, 2021, 1:07 a.m.