Loading the libraries and the dataset

# initialisation
graphics.off()
rm(list = ls())

# global parameters
flag.verbose = FALSE
opers = EStatOption_fromKeys(c("NUM", "MINI", "MAXI", "MEAN", "STDV"))
set.seed(25854)

Initial data set

# grid
nx  = c(100, 100)
dx  = 1/nx
x0  = c(0., 0.)
grd = DbGrid_create(nx = nx, dx = dx, x0 = x0)

# neighborhood
neigh = NeighMoving_create(nmaxi = 50, nmini = 10, radius = 1.0)

# model
rho = 0.7
c0  = matrix(c(1.0, rho, rho, 1), nrow = 2, ncol = 2)
mod = Model_createFromParam(type = ECov_EXPONENTIAL(), range = 0.1, sills = c0)
err = mod$setMeans(c(100, 0))
# Initial simulation
err = simtub(dbin = NULL, dbout = grd, model = mod, nbsimu = 1, 
             namconv = NamingConvention("Y"))

err = grd$setName("Y.1", "Z1")
err = grd$setName("Y.2", "Z2")

knitr::kable(dbStatisticsMono(grd, names = "Z*", opers = opers, flagIso = FALSE)$toTL(),
             digits = 3, 
             caption = "Statistics on the initial simulations")
Statistics on the initial simulations
Number Minimum Maximum Mean St. Dev.
Z1 10000 96.239 103.799 100.056 1.034
Z2 10000 -3.874 3.877 -0.001 1.047 
print(paste0("Correlation = ", round(cor(grd["Z1"], grd["Z2"]), 3)))
## [1] "Correlation = 0.727"

Deriving a heterotopic data base

np = 200

dat = Db_createSamplingDb(db = grd, number = np, names = c("x1", "x2", "Z1", "Z2"), 
                          seed = 123)

# discarding 1/4 of the initial data on Z1
sel = sort(sample(1:np, np/4, replace = FALSE))
Z1 <- dat["Z1"]
Z1[sel] <- NaN
dat["Z1.heterotopic"] <- Z1

# discarding 1/5 of the initial data on Z2
sel = sort(sample(1:np, np/5, replace = FALSE))
Z2 <- dat["Z2"]
Z2[sel] <- NaN
dat["Z2.heterotopic"] <- Z2

knitr::kable(dbStatisticsMono(dat, names = "Z*.heterotopic", opers = opers, 
                              flagIso = FALSE)$toTL(), digits = 3, 
             caption = "Statistics on the heterotopic data set")
Statistics on the heterotopic data set
Number Minimum Maximum Mean St. Dev.
Z1.heterotopic 150 97.695 102.961 100.087 1.02
Z2.heterotopic 160 -2.507 2.602 0.073 1.08 
knitr::kable(dbStatisticsMono(dat, names = "Z*.heterotopic", opers = opers, 
                              flagIso = TRUE)$toTL(), digits = 3, 
             caption = "Statistics on the homotopic sub set")
Statistics on the homotopic sub set
Number Minimum Maximum Mean St. Dev.
Z1.heterotopic 125 97.695 102.961 100.128 1.028
Z2.heterotopic 125 -2.124 2.602 0.127 1.057

Conditional simulations

- with homotopic data set

err = dat$setLocators(names = paste(c("Z1", "Z2"), sep = "."), 
                      locatorType = ELoc_Z(), cleanSameLocator = TRUE)
err = simtub(dbin = dat, dbout = grd, model = mod, neigh = neigh, nbsimu = 1, 
             seed = 2597,
             namconv = NamingConvention("homotopic"))

# statistics
knitr::kable(
  rbind(
    dbStatisticsMono(grd, names = paste0("homotopic*"), opers = opers, 
                     flagIso = FALSE)$toTL()
  ), digits = 4, caption = "Statistics on simulated data"
)
Statistics on simulated data
Number Minimum Maximum Mean St. Dev.
homotopic.Z1 10000 95.8759 103.4401 100.0228 0.9996
homotopic.Z2 10000 -3.2586 3.7649 0.0872 1.0332 
# histogram and base map
p1 = plot.init() + plot.hist(grd, name = "homotopic.Z1")
p2 = plot.init() + plot.hist(grd, name = "homotopic.Z2")
p3 = plot.init(asp=1) + plot.raster(grd, name = "homotopic.Z1", palette = "Spectral", 
                             flagLegend = TRUE, legendName = "Z1")
p4 = plot.init(asp=1) + plot.raster(grd, name = "homotopic.Z2", palette = "Spectral", 
                             flagLegend = TRUE, legendName = "Z2")
ggarrange(p1, p2, p3, p4, nrow = 2, ncol = 2)
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.

- with heterotopic data set

err = dat$setLocators(names = paste(c("Z1", "Z2"), "heterotopic", sep = "."), 
                      locatorType = ELoc_Z(), cleanSameLocator = TRUE)
err = simtub(dbin = dat, dbout = grd, model = mod, neigh = neigh, nbsimu = 1, 
             seed = 2597,
             namconv = NamingConvention("heterotopic"))

# statistics
knitr::kable(
  rbind(
    dbStatisticsMono(grd, names = paste0("heterotopic*"), opers = opers, 
                     flagIso = FALSE)$toTL()
  ), digits = 4, caption = "Statistics on simulated data"
)
Statistics on simulated data
Number Minimum Maximum Mean St. Dev.
heterotopic.Z1.heterotopic 10000 96.0065 103.6438 100.0153 1.0215
heterotopic.Z2.heterotopic 10000 -3.3246 3.5293 0.0841 1.0305 
# histogram and base map
p1 = plot.init() + plot.hist(grd, name = "heterotopic.Z1.heterotopic")
p2 = plot.init() + plot.hist(grd, name = "heterotopic.Z2.heterotopic")
p3 = plot.init(asp=1) + plot.raster(grd, name = "heterotopic.Z1.heterotopic", 
                             palette = "Spectral", 
                             flagLegend = TRUE, legendName = "Z1")
p4 = plot.init(asp=1) + plot.raster(grd, name = "heterotopic.Z2.heterotopic", 
                             palette = "Spectral", 
                             flagLegend = TRUE, legendName = "Z2")
ggarrange(p1, p2, p3, p4, nrow = 2, ncol = 2)
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.
## `stat_bin()` using `bins = 30`. Pick better value `binwidth`.