[R] Cross-variograms

[Edzer J. Pebesma]

Jacques, provided that X and Y are colocated (i.e., have
exactly the same observation locations), you get the
cross variogram right; the definition of this cross
variogram is however:
gamma(h)= E[(X(s)-X(s+h))*(Y(s)-Y(s+h))]
also, where you select:
cv <- v$gamma[1:14]
you may be better off using the more general
v$gamma[v$id == "X.Y"]
Best regards,
--
Edzer
__
[EMAIL PROTECTED] mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

[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

__
[EMAIL PROTECTED] mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html