3

I've created a nice plot using scatter3d() and Rcmdr. That plot contains two nice surface smooths. Now I'd like to add to this plot one more surface, the truth (i.e. the surface defined by the function generating my observations minus the noise component).

Here is my code so far:

library(car) 
set.seed(1)

n <- 200 # number of observations (x,y,z) to be generated
sd <- 0.3 # standard deviation for error term
x <- runif(n) # generate x component
y <- runif(n) # generate y component

r <- sqrt(x^2+y^2) # used to compute z values

z_t <- sin(x^2+3*y^2)/(0.1+r^2) + (x^2+5*y^2)*exp(1-r^2)/2 # calculate values of true regression function

z <-  z_t + rnorm(n, sd = sd) # overlay normally distrbuted 'noise' 
dm <- data.frame(x=x, y=y, z=z) # data frame containing (x,y,z) observations
dm_t <- data.frame(x=x,y=y, z=z_t) # data frame containing (x,y) observations and the corresponding value of the *true* regression function

# Create 3D scatterplot of:
# - Observations (this includes 'noise')
# - Surface given by Additive Model fit 
# - Surface given by bivariate smoother fit

scatter3d(dm$x, dm$y, dm$z, fit=c("smooth","additive"), bg="white", 
          axis.scales=TRUE, grid=TRUE, ellipsoid=FALSE, xlab="x", ylab="z", zlab="y")

The solution given in another thread is to then define a function:

my_surface <- function(f, n=10, ...) { 
  ranges <- rgl:::.getRanges()
  x <- seq(ranges$xlim[1], ranges$xlim[2], length=n)
  y <- seq(ranges$ylim[1], ranges$ylim[2], length=n)
  z <- outer(x,y,f)
  surface3d(x, y, z, ...)
}

f <- function(x, y)
  sin(x^2+3*y^2)/(0.1+r^2) + (x^2+5*y^2)*exp(1-r^2)/2

my_surface(f, alpha=0.2)

This however yields an error, saying (translated from German since this is my system language, I apologize):

Error in outer(x, y, f) : 
  Dimension [Product 100] does not match the length of the object [200]

I then tried an alternative approach:

x <- seq(0,1,length=20)
y <- x
z <- outer(x,y,f)
surface3d(x,y,z)

This does add a surface to my plot but it doesn't look right at all (i.e. the observations are not even close to it). Here's what the supposed true surface looks like (this is obviously wrong): Thanks!

The wrong "true" regression surface

I think the problem may in fact be scaling. Here I created a couple of points that sit on the plane z = x+y. Then I proceeded to try to plot that plane using my method above:

library(car)

n <- 50
x <- runif(n)
y <- runif(n)
z <- x+y

scatter3d(x,y,z, surface = FALSE)

f <- function(x,y) 
  x + y 

x_grid <- seq(0,1, length=20)
y_grid <- x_grid 
z_grid <- outer(x_grid, y_grid, f)
surface3d(x_grid, y_grid, z_grid)

This gives me the following plot:

The scaling issue

Maybe one of you can help me out with this?

JasonMArcher
  • 14,195
  • 22
  • 56
  • 52
user2249626
  • 453
  • 2
  • 6
  • 15
  • 3
    You asked a simillar question, for which you received a good answer, but you have *never* voted on, or accepted an answer. This will make people unwilling to help you because you give nothing back to the community. Read the [**Help**](http://stackoverflow.com/help) and the [**About**](http://stackoverflow.com/about) page to see how SO is *supposed* to work. – Simon O'Hanlon Aug 14 '13 at 13:20
  • @HaskellElephant I never said anything about *my* choices. I simply pointed out that persons may be unwilling to help when they see a user has never pointed. Where do you draw that conclusion from? That aside, at the point at which I left the comment the question included no code and no reproducible example. Now I see the OP has left a great reproducible example, and has a better sense of the community because they have accepted an answer (which by their own admission solved their problem). So my comment enriched the community and the OPs experience of it. More involvement = better. – Simon O'Hanlon Aug 14 '13 at 14:16
  • 1
    +1 for taking the time to edit this question. This is now a great reproducible example./ – Simon O'Hanlon Aug 14 '13 at 14:17

1 Answers1

0

The scatter3d function in car rescales data before plotting it, which makes it incompatible with essentially all rgl plotting functions, including surface3d.

You can get a plot something like what you want by using all rgl functions, e.g. plot3d(x, y, z) in place of scatter3d, but of course it will have rgl-style axes rather than car-style axes.

user2554330
  • 37,248
  • 4
  • 43
  • 90