Dear all,In continuation of this thread, I've spent some time looking at kriging on the sphere, corrected some bugs, and need further help.
First of all, distances for covariances on the sphere were computed incorrectly as well in gstat, so I hope not too many people have been relying on this--kriging seemed to happen still in some Euclidian mode. The good news is that it seems to work now (gstat 0.9-53, accepted on CRAN). Inverse distance interpolation for spherical data seemed to work already.
Covariances on the sphere now work, but the models present do not include those specially deviced for spherical data. Can anyone provide me with or point me to useful, preferably simple covariance functions that are positive definite on the sphere? The example below seems to work but in certain cases without nugget the interpolation may go crazy. You'll need to update your sp to 0.9-27 to run it (also accepted on CRAN).
Below is an example script. It also needs the new sp and gstat versions.
library(gstat)
library(rgdal)
world = expand.grid(long=seq(-177.5,177.5,5),lat=seq(-87.5,87.5,5))
world.sp = SpatialPixels(SpatialPoints(world,CRS("+proj=longlat")))
plot(world.sp,axes=T)
pts=data.frame(long=runif(100,-180,180),lat=runif(100,-90,90),val=rnorm(100))
coordinates(pts)=~long+lat
proj4string(pts)=CRS("+proj=longlat")
points(pts,col='red')
# inverse distance interpolation on the sphere:
idw.out = idw(val~1,pts,world.sp)
image(idw.out, axes = TRUE, ylim = c(-90,90))
points(pts, pch=3)
idw.spdf = as(idw.out, "SpatialPolygonsDataFrame")
newproj = CRS("+proj=moll")
idw.spdf.moll = spTransform(idw.spdf, newproj)
spplot(idw.spdf.moll, "var1.pred",col.regions=bpy.colors(),col=0,
sp.layout = list(sp.points, spTransform(pts, newproj), col = 'black'))
# kriging on the sphere
kr.out = krige(val~1,pts,world.sp,vgm(1, "Exp", 3000))
idw.spdf.moll$kr = kr.out[[1]]
spplot(idw.spdf.moll, "kr", col.regions=bpy.colors(), col=0,
sp.layout = list(sp.points, spTransform(pts, newproj), col = 'black'))
--
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.springer.com/978-0-387-78170-9 [EMAIL PROTECTED]
<<inline: world.png>>
_______________________________________________ R-sig-Geo mailing list R-sig-Geo@stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-geo
