How to generate a smoothed surface for a three dimensional dataset?
On Dec 5, 2013, at 4:02 PM, Jun Shen wrote:
Thanks again, Duncan. Please allow me to ask one more question. Is it possible to generate a contour plot overlaying with the plot3d() plot?
This is not going to give you the kewl rotation capabilities that `rgl` offers but you might be interested in the potential of the 'plot3d' package from Karline Soetaert of NIOZ-Yerseke, The Netherlands. (She also is also part of the team that gave us the very nice ODE package, 'deSolve' and co-wrote the accompanying book. I found it interesting to read that she is a biologist, while co-authors are mathematicians.) http://cran.rstudio.com/web/packages/plot3D/vignettes/plot3D.pdf See Figures 2 and 8 http://cran.rstudio.com/web/packages/plot3D/vignettes/volcano.pdf The second one has the catchy title: "Fifty ways to draw a volcano using package plot3D" You may find that the examples that drew the 4 examples in Figure 7 would give you something useful.
Jun On Thu, Dec 5, 2013 at 11:08 AM, Duncan Murdoch <murdoch.duncan at gmail.com>wrote:
On 05/12/2013 10:33 AM, Jun Shen wrote:
Hi Federico/Duncan/David/Bert, Thank you for your thoughtful comments and it's a great learning experience. I can see the critical point here is to find a right function to make the prediction. So I was thinking to start with "loess". However the predict.loess gave me an error as follows Error in `$<-.data.frame`(`*tmp*`, "z", value = c(0.417071766265867, 0.433916401753023, : replacement has 20 rows, data has 400 Here is the code I tried. Thank you for your help again! Jun ===================================== x<-runif(20) y<-runif(20) z<-runif(20) library(rgl) plot3d(x,y,z) loess(z~x+y,control=loess.control(surface='direct'), span=.5,degree=2)->fit.loess xnew <- seq(min(x), max(x), len=20) ynew <- seq(min(y), max(y), len=20) df <- expand.grid(x = xnew, y = ynew) df$z<-predict(fit.loess,newdata=df)
After the error, use traceback() to find which function called `$<-.data.frame`. It shows that it was your final assignment df$z<-predict(fit.loess,newdata=df) which causes the error, because the predict function returns a matrix. So you can get the plot using surface3d(xnew, ynew, predict(fit.loess,newdata=df), col="gray") You may want aspect3d(1,1,1) afterwards; loess isn't so good at extrapolation. Or you may want to set predictions to NA outside the convex hull of your data. (I'm not sure which function is easiest to calculate that, but here's one way: hullx <- x[chull(x,y)] hully <- y[chull(x,y)] keep <- sp::point.in.polygon(df$x, df$y, hullx, hully) znew <- predict(fit.loess,newdata=df) znew[!keep] <- NA plot3d(x,y,z) surface3d(xnew, ynew, znew, col="gray") aspect3d(1,1,1)
David Winsemius Alameda, CA, USA