I've been experimenting with the rgl library to draw 3d plots of some linear models. The package works really nicely, but I have a question on the math behind how one draws a plane (using planes3d) that shows the 'fit' (predicted values ) of the linear model, for a model with 2 parameters.
The R package documentation says that the arguments -- a, b, c, and d -- of a call to planes3d must be chosen so they are the coordinates of the normal to the plane (given by the model). And there is an example on page 33 of this document: http://cran.r-project.org/web/packages/rgl/rgl.pdf, which I have slightly modified to use the data from the mtcars package (see R code below).
In the listing I create the model from the mtcars data at line 3, then I plot the data points at line 6. Next I plot the plane in lines 9 - 15. The a,b,c,d coordinates are chosen per instructions in the documentation. In particular,'c' is set to be -1, which is currently just a 'magic' value to me. I don't understand why -1 is chosen for 'c'.
But it works out. In line 18 - 26 I define vectors for the normal to the plane, and for a point lying within the plane. I take the dot product of these two and add the intercept (offset) and I get zero. It all works out as expected. But I am hoping someone can explain the theory behind >why< it works out.
1 library(rgl) 2 attach (mtcars) 3 model <- lm(mpg ~ cyl + disp) 4 5 # Plot data points 6 plot3d(cyl, disp, mpg, type="s", col="red", size=1) 7 8 9 # Plot the plane of points that make up the predicted values of the model 10 coefs <- coef(model) 11 a <- coefs["cyl"] 12 b <- coefs["disp"] 13 c <- -1 14 d <- coefs["(Intercept)"] 15 planes3d(a, b, c, d, alpha=0.5) 16 17 18 # define a normal vector to the plane 19 # 20 normal = c(a,b,c) 21 22 # define a vector within the plane 23 # 24 newdata <- data.frame(cyl = 5, disp = 200) 25 predicted = predict(model, newdata) 26 pointOnPlane = c(5, 200, predicted) 27 28 29 30 # Verify that the dot products of the normal and the vector that defines pointOnPlane is 0 31 # (after we subtract the offset 'd', which came from the intercept of the model). 32 33 # This will print as 'TRUE' 34 sum(pointOnPlane * normal) + d == 0
epilogue: Got a great answer from Andre Silva, below. Thanks, Andre ! makes sense now ;^)