[R-sig-Geo] sim kriging

[Edzer Pebesma]

On 11/11/2011 02:25 AM, Javier Leon Patino wrote:

> 1. I would like to remove a quadratic trend from my elevation dataset before ordinary kriging. Maybe using surf.ls?
> 
> Would the following work (Z is elevation in file ahd)?:
> 
> depth.ok <- krige(ahd~Z, shp, mask, depth.vgm)

this would do universal kriging (external kriging, regression kriging)
with a first order linear trend in Z. It means that the trend is not
removed as such, but is modelled as a first order linear function in Z,
and that depth.vgm should be the variogram of residuals wrt this function.

depth.ok2 = krige(ahd~Z+I(Z^2), shp, mask, depth.vgm)

would do this wrt to a second order (quadratic) trend, but depth.vgm
should again be the variogram of residuals from such a model.

specifying such models in calls to variogram() indeed gives you the
variogram of residuals (computed, by default, using OLS).

> 
> 
> 
> 2. How do I write simulation kriging to ascii? I am using:
>  
> depth.ok <- krige(ahd~1, shp, mask, depth.vgm, nsim=2, nmax=10)
> for (i in 1:length(depth.ok)) {
> 	write.asciigrid(depth.ok[i], c(paste("DEMair",as.character(i),".asc",sep="")))
>  }
> 

length(depth.ok) gives you the number of grid cells depth.ok has, use
ncol(depth.ok) to obtain the number of "layers", or attribute variables,
it has.

> Ascii grids seem fine but I get the following error:
> 
>> Error in gridded(x) : 
>>
>> error in evaluating the argument 'obj' in selecting a method for 
>> function 'gridded': Error in `[.data.frame`(x at data, , i, drop = FALSE) 
>> :
>>
>> undefined columns selected
> 
> 
> Thank you very much for your help. These two issues have stopped my work for more than a week now!
> 
> Cheers,
> 
> Javier
> 
> 
> Dr Javier Leon
> Post-doctoral Research Fellow
> 
> Global Change Institute/
> Centre for Spatial Environmental Research (CSER)
>  
> The University of Queensland
> Gehrmann Laboratories (60) 
> The University of Queensland 
> St Lucia QLD 4072 Australia
> Ph: (+61 7) 3346 7074   
> Fax: (+61 7) 3346 3299
> email: j.leonpatino at uq.edu.au
> http://www.gpem.uq.edu.au/javier-leon-patino
> 
> 
> 


-- 
Edzer Pebesma
Institute for Geoinformatics (ifgi), University of Münster
Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251
8333081, Fax: +49 251 8339763  http://ifgi.uni-muenster.de
http://www.52north.org/geostatistics      e.pebesma at wwu.de