Designing multivariate density plot in R

You could start with the persp function to draw the 3 dimensional plot (if you do this from data rather than the formula then you need to use some form of density estimation first, the example plot looks smooth enough that it is probably based on the formula rather than estimated from the data). Then use the return value from persp to project the additional plotting info.

There may also be an option using the rgl package, I seem to remember that it has a way to project a plot onto the axes planes.

Edit

Here is some sample code to get you started. It uses a parametric distribution, but could be adapted to use kde2d from MASS or other ways of estimating the density from data:

x <- seq( -3, 3, length=25 )y <- seq( -3, 3, length=25 )z <- outer( x, y, function(x,y) dnorm(x,0,0.5)*dnorm(y,0,1) )zl <- c(0,4*max(z))## persp plottrmat <- persp(x,y,z, theta=120, zlim=zl, box=FALSE, shade=0.5)## x gridfor( i in seq(-3,3, by=0.5 ) ) {    lines( trans3d( c(i,i), c(-3,-3), zl, trmat ), col='grey' )}for( i in seq(0,zl[2], length=7) ) {    lines( trans3d( c(-3,3), c(-3,-3), c(i,i), trmat ), col='grey' )}## marginal for xlines( trans3d( seq(-3,3,length=100), -3, dnorm(seq(-3,3,length=100),0,.5),     trmat), lwd=2, col='blue' )## y gridfor( i in seq(-3,3, by=0.5 ) ) {    lines( trans3d( c(-3,-3), c(i,i), zl, trmat ), col='grey' )}for( i in seq(0,zl[2], length=7) ) {    lines( trans3d( c(-3,-3), c(-3,3), c(i,i), trmat ), col='grey' )}## marginal for ylines( trans3d( -3, seq(-3,3,length=100), dnorm(seq(-3,3,length=100),0,1),     trmat), lwd=2, col='blue' )