[R] Cross-variograms

[Jacques VESLOT]

Dear All,

In order to appreciate the spatial link between to variables: X (numeric)
and Y (factor), measured on a grid at the same locations (i.e. a
geodata-like object with a two-column dataframe as $data), I want to use a
so-called cross-variogram based on such a formula:

h(r)= E[(X(0)-X(h))*(Y(0)-Y(h))] = D(0) - D(h) - D(-h)  where D(h) is the
covariance between X(0) and Y(h),

and then draw an 95-envelop estimated by sampling one dataset (Y) by
translation.

GeoR and geoRglm doesn't include cross-variograms, thus I tried with
variogram() in gstat. But, despite help pages and gstat website, I
encountered difficulties in understanding how to deal with this function,
and notably with the formula argument.

I tried to create the gstat object I need and get the related
cross-variogram (cv), just as follows:

g <- gstat(id="X", formula=X~1, locations=~x+y, data = mat)		## mat is a
4-column df with variables X and Y
g <- gstat(g, id="Y", formula=Y~1, locations=~x+y, data = mat)	## and
(x,y)-coords
v <- variogram(g, cutoff=15)
cv <- v$gamma[1:14]

But I am definitely not sure of the way to deal with formula argument, since
I didn't understand its role here, and how to create a gstat object with two
variables.

Could you please let me know if this is the right way to operate (and what
formula argument means here) ?


Best regards,

Jacques VESLOT
CIRAD