RMbubble {RandomFields}

R Documentation

Bubble model for arbitrary areas of scales

Description

A model that allows for arbitray areas of scale applied to an isotropic model, i.e.

C(x,y) = \phi(\|x -y \| / s)

as long as s_x = s_y = s. Here, s_x is the scaling at location x,

The cross-correlations between areas of different scales are given through a modified distance d. Let z_{s} be a finite subset of R^d depending on the scale s. Let w_u be a weight for an auxiliary point u\in z_{s} with \sum_{u \in z_s} w_u = 1. Let \tau_x = s_x^{-2}. Then

d^2(x, y) = \min\{\tau(x), \tau(y)\} \|x - y\|^2 + \sum_{\xi \in_{span(\tau(x), \tau(y))}} \sum_{u \in z_{\xi^{-0.5}}} w_u \|x - u\|^2 \Delta \xi

Here, span(\tau(x), \tau(y)) is the finite set of values s^{-2} that are realized on the locations of interest and \Delta \xi is the difference of two realized and ordered values of the scaling s.

Usage

RMbubble(phi, scaling, z, weight, minscale, barycentre,
         var, scale, Aniso, proj)

Arguments

phi	isotropic submodel
scaling	model that gives the non-stationary scaling s_x
z	matrix of the union of all z_s. The number of rows equals the dimension of the field. If not given, the locations with non-vanishing gradient are taken.
weight	vector of weights w whose length equals the number of columns of z. The points given by z might be weighted.
minscale	vector for partioning z into classes z_s. Its length equals the number of columns of z. The vector values must be descending. See details. If not given then z_s=z for all s. Else see details.
barycentre	logical. If FALSE and z is not given, the reference locations are those with non-vashing gradient. If TRUE then, for each realized value of the scale, the barycentre of the corresponding reference locations is used instead of the reference locations themselves. This leads to higher correlations, but also to highly non-stationary cross-correlation between the areas of different scale. The argument has no effect when z is given. Default: FALSE.
var, scale, Aniso, proj	optional arguments; same meaning for any RMmodel. If not passed, the above covariance function remains unmodified.

Details

minscale gives the minimal scale s value above which the corresponding points z define the set z_s. The validity of the set z_s ends with the next lower value given.

Let minscale = (10, 10, 10, 7, 7, 7, 0.5). Then for some d-dimensional vectors z_1,\ldots, z_7 we have

z_s = \{ z_1, z_2, z_3 \}, s \ge 10

z_s = \{ z_4, z_5, z_5 \}, 7 \ge s < 10

z_s = \{ z_7 \}, s \ge 0.5

Note that, in this case, all realized scaling values must be \ge 0.5. Note further, that the weights for the subset must sum up to one, i.e.

w_1+w_2 +w_3=w_4 + w_5 + w_6 = w_7 = 1.

Value

RMbubble returns an object of class RMmodel.

Note

This model is defined only for grids.

Author(s)

Martin Schlather, schlather@math.uni-mannheim.de, https://www.wim.uni-mannheim.de/schlather/

References

• Bonat, W.H. , Ribeiro, P. Jr. and Schlather, M. (2019) Modelling non-stationarity in scale. In preparation.

See Also

RMSadvanced, RMblend, RMscale

Examples

RFoptions(seed=0) ## *ANY* simulation will have the random seed 0; set
##                   RFoptions(seed=NA) to make them all random again

x <- seq(0,1, if (interactive()) 0.02 else 0.5)
d <- sqrt(rowSums(as.matrix(expand.grid(x-0.5, x-0.5))^2))
d <- matrix(d < 0.25, nc=length(x))
image(d)

scale <- RMcovariate(data=as.double(d) * 2 + 0.5, raw=TRUE)

## two models:
## the frist uses the standard approach for determining the
##           reference point z, which is based on gradients
## the second takes the centre of the ball
model1 <- RMbubble(RMexp(), scaling=scale)
model2 <- RMbubble(RMexp(), scaling=scale, z=c(0.5, 0.5))
model3 <- RMbubble(RMexp(), scaling=scale, barycentre=TRUE) # approx. of model2

## model2 has slightly higher correlations than model1:
C1 <- RFcovmatrix(model1, x, x)
C2 <- RFcovmatrix(model2, x, x)
C3 <- RFcovmatrix(model3, x, x) 
print(range(C2 - C1))
dev.new(); hist(C2 - C1)
print(range(C3 - C2)) # only small differences to C2
print(mean(C3 - C2))
dev.new(); hist(C3 - C2)

plot(z1 <- RFsimulate(model1, x, x))
plot(z2 <- RFsimulate(model2, x, x))
plot(z3 <- RFsimulate(model3, x, x)) # only tiny differences to z2


## in the following we compare the standard bubble model with
## the models RMblend, RMscale and RMS (so, model2 above
## performs even better)
biwm <- RMbiwm(nudiag=c(0.5, 0.5), nured=1, rhored=1, cdiag=c(1, 1), 
                s=c(0.5, 2.5, 0.5))
blend <- RMblend(multi=biwm, blend=RMcovariate(data = as.double(d), raw=TRUE))
plot(zblend <- RFsimulate(blend, x, x)) ## takes a while ...
Cblend <- RFcovmatrix(blend, x, x)

Mscale <- RMscale(RMexp(), scaling = scale, penalty=RMid() / 2)
plot(zscale <- RFsimulate(Mscale, x, x))
Cscale <- RFcovmatrix(Mscale, x, x)

Mscale2 <- RMscale(RMexp(), scaling = scale, penalty=RMid() / 20000)
plot(zscale2 <- RFsimulate(Mscale2, x, x))
Cscale2 <- RFcovmatrix(Mscale2, x, x)

S <- RMexp(scale = scale)
plot(zS <- RFsimulate(S, x, x))
CS <- RFcovmatrix(S, x, x)

print(range(C1 - CS))
print(range(C1 - Cscale))
print(range(C1 - Cscale2))
print(range(C1 - Cblend))
dev.new(); hist(C1-CS)     ## C1 is better
dev.new(); hist(C1-Cscale) ## C1 is better
dev.new(); hist(C1-Cscale2) ## both are equally good. Maybe C1 slightly better
dev.new(); hist(C1-Cblend) ## C1 is better

---