[R-sig-Geo] Trouble with universal kriging

[Roger Bivand]

On Thu, 20 Jun 2013, Edzer Pebesma wrote:

> Dominik, if you want to receive help, please explain why you are
> frustrated, and instead of claims of "won't work", explain in detail
> what you did and what happened, preferably with a well-prepared example
> that people on the list can reproduce without downloading data first.
>
>
> On 06/20/2013 09:57 PM, Dominik Schneider wrote:
>> Hi, I am trying to follow the examples in  Bivand et al 'Applied Spatial
>> Data Analysis' and the examples at
>> http://casoilresource.lawr.ucdavis.edu/drupal/node/442 but I am frustrated.
>> I finally managed an ordinary kriging example to work but with little
>> understanding of why/when it should work.

Well, if you read chapter 8, you should gain sufficient understanding, 
right? The online examples are not intended to be self-explanatory; they 
all work, because they are tested nightly on a fully updated system; the 
online versions have been modified to keep up with successive changes in 
R and the packages used.

>> I posted some data: https://www.dropbox.com/s/rtg7gh67ziai3ec/mar2011.RData
>> Another item I've noticed, if I  define the projection with
>> projection(object)='+proj=lonlat +datum=NAD83'  for any of my spatial
>> objects it won't work.  It only works with an defined projection. any help
>> would be appreciated. Thanks.

Do also pay attention to details:

"+proj=lonlat +datum=NAD83"

is valid PROJ.4, but not recognised in the CRS class, it should be:

"+proj=longlat +datum=NAD83"

All of "latlon, "latlong", "lonlat", and "longlat" are valid PROJ.4, but 
the first two imply northings before eastings, and we decided ten years 
ago to use only "longlat". The documentation is in ?proj4string in sp, but 
the detail is not repeated in ?projection in raster, which does link to 
?"CRS-class" in sp and "spTransform-methods" in rgdal. We will try to make 
sure that use of anything other than "longlat" is an error.

Always include your sessionInfo(), which gives details of your platform 
and the versions of installed packages.

Hope this clarifies,

Roger

>>
>> load('mar2011.RData')
>>
>> vt=variogram(snotel~recon,data=mar2011)
>>
>> vt.fit=fit.variogram(vt,vgm(psill=0.025,nugget=0.015,model='Sph',range=2))
>>
>> #plot(vt,vt.fit)
>>
>> ok<-krige(snotel~1,locations=mar2011,newdata=swe.grid,model=vt.fit)
>>
>> uk<-krige(snotel~mar2011$recon,locations=mar2011,newdata=swe.grid,model=
>> vt.fit)
>>
>> #spplot(ok)
>> Dominik
>>
>> 	[[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-Geo mailing list
>> R-sig-Geo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>>
>
>

-- 
Roger Bivand
Department of Economics, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no