[R-sig-Geo] External Drift Kriging

Paulo Justiniano Ribeiro Jr paulojus at est.ufpr.br
Thu Jan 26 15:00:21 CET 2006


Hi Luke

> krige.conv command. I am trying to krige annual average temperatures
> from 112 met stations located in the Eastern Cascade mountains of
> Oregon. The data are highly non-stationary, due to the presence of a
> linear lapse response between the temperature and the altitude (the
> altitudinal range here is considerable). I have quantified this response
> at a loss of 3.5 degrees centegrade per 1km rise. I have removed this
> trend from the data, and produced a zero mean, gaussian distributed
> stationary process, and from this derived the spatial covarince structure.
>

OK. This is one option indeed, to modelling covariance after getting 
"residuals" from some kind of model to take auto the effect of known 
sources of variation.
And indeed the usual variogram will only make sense for stationary data.

In the max. lik. estimation you can do both together, getting estimates 
for a (linear) model for the mean and for the covariance parameters.

An example call would be something along the lines

> ml <- likfit(foo, trend=~altitude)
provided foo is a geodata object with covariate information on altitude

Later for kriging you would need a call of the type:
> krige.conv(foo, loc=gr, krige=krige.control(obj=ml, trend.d=~altitude, 
trend.l = ~altitude.gr), borders=NULL)

where gr is a matrix with coordinates of prediction points and
altitude.gr is the value of the covariate at those prediction points,

> The problem is that I only want to predict values from a small subset of
> the area which the stations occupy; although some stations are outwith
> the required borders, they are more or less within the range of the
> covariance function. I have an elevation model for the prediction
> locations, but due to computational limitations, cannot extend this to
> cover the full area: thus trend.l is only specified within the borders.

So you can only krige within the borders, unless you extrapolate your 
covariate information somehow.
Kriging with external trend (w/ covariates) requires values of the 
coavriates at all prediction locations.


Hope this helps, and if the problempersists feel free to contact me 
directly

best
P.J.



> All of the stations have elevations attached, and so I am able to
> specify trend.d. To summarise I have the following:
>
> grid                     [a 1km grid which covers the entire domain of
> interest, containing all the met stations]
> lc                        [a two dimensional matrix specifying the
> borders for prediction]
> Normals             [a geodata object with data on the annual station
> temperature average]
> pts                      [a geodata object with data on the station
> elevations]
> cressie_cu           [a variogram model (cressie weights, cubic model),
> derived from the residuals of the linear model average_temp~elevation]
> DEM                  [a geodata object with elevation data at 1km
> spacing for all prediction locations within lc]
>
> I am attempting the following:
>
> pars <- krige.control(obj.m=cressie_cu,
>                               trend.d = trend.spatial(~data, pts),
>                               trend.l = trend.spatial(~data,DEM))
>
> kc <- krige.conv(Normals, loc = grid, borders = lc, krige = pars)
>
> On trying to execute the commands, I get an output which has lots of
> lines and other artifacts running through it, which seems to bear no
> resemblance to the expected output at all. Could you please suggest what
> I am getting wrong?
>
> Many thanks
>
> Luke Spadavecchia
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>

Paulo Justiniano Ribeiro Jr
LEG (Laboratório de Estatística e Geoinformação)
Departamento de Estatística
Universidade Federal do Paraná
Caixa Postal 19.081
CEP 81.531-990
Curitiba, PR  -  Brasil
Tel: (+55) 41 3361 3573
Fax: (+55) 41 3361 3141
e-mail: paulojus at est.ufpr.br
http://www.est.ufpr.br/~paulojus


More information about the R-sig-Geo mailing list