3

I have 2 variables x and y which are Cartesian coordinates at [0,1], and z is the value of a (continuous) variable at these coordinates. The z vector has some important outliers

x<-sample(seq(0,1,0.001), replace=F)
y<-sample(seq(0,1,0.001), replace=F)
z<-runif(1001,min=0,max=1)
z[100]<-8;z[400]<-16;z[800]<-4

These outliers I would like to emphasize when presenting these data in a filled.contour

I have used until now

library(akima)
a<-interp(x,y,z)
filled.contour(a$x,a$y,a$z)

But I am not happy with this linear interpolation. For example (the outliers do not show up correctly).

I am thinking what I need is a some kind nearest neighbor "spatial" smoothing of z (based on the x,y location). Can anyone help or pinpoint to data/examples/packages/code that could help me? I would prefer a base R solution but if ggplot2 or lattice can do my job it would be fine. Any other idea/proposal of better visualization would be also welcomed.

ECII
  • 10,297
  • 18
  • 80
  • 121
  • 1
    You could try spline interpolation, e.g., `a<-interp(x,y,z,linear=FALSE,xo=seq(min(x), max(x), length = 200),yo=seq(min(y), max(y), length = 200)); filled.contour(a$x,a$y,a$z,color.palette=terrain.colors)`. – Roland Dec 02 '12 at 10:53
  • I dont know how to understand the key values (z range) when I use spline. In this example linear=T give a range of [0,14] as in the data, but linear=F gives a range of (-100,100) – ECII Dec 02 '12 at 11:05

2 Answers2

7

Check out @Paul Heimstra's answer here. They suggest thin plate spline interpolation from the fields package. If kriging is right for you, the package automap may be useful.

If you're concerned with the interpolation of geographic values, it's worth mentioning DG Rossiter's online/free course on Geostatistics and Open-source statistical computing.

Community
  • 1
  • 1
metasequoia
  • 7,014
  • 5
  • 41
  • 54
4

Thin plate splines, or the Tps() function, from the fields package is a your ticket. The model results must then be placed into the predictSurface() function from the same package to generate smoothed spatial data. Just a few lines of code is all you need:

test.spline <- Tps(data.frame(x,y), z)
new.grid <- predictSurface(test.spline, nx = 200, ny = 200)
image(new.grid)

The nature of the data will define how fine a grid is appropriate. You control this with nx and ny arguments in predictSurface(). In this example, increased grid resolution has no impact. With real spatial data, smoothing with a high-res surface can be dramatic. Check function help for more details.

Brad Horn
  • 649
  • 6
  • 12