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]
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@stat.math.ethz.ch
https://stat.ethz.ch/mailman/listinfo/r-sig-geo
