\begin{Schunk}
\begin{Sinput}
 library(mvtnorm) 
 N <- 100 
 myProbX <- function(x1, x2){
      dmvnorm( cbind(x1,x2), mean = meanX, sigma= varX)
    }
 meanX <- c( 0, 0 )
 varX <- matrix ( c( 1, .3, .3, 1), nrow=2 , byrow = TRUE) 
 x1 <- seq( meanX[1]-2.5*varX[1,1], meanX[1]+2.5*varX[1,1], length.out=N)
 x2 <- seq(  meanX[2]-2.5*varX[2,2], meanX[2]+2.5*varX[2,2] , length.out=N)
 probX <- outer(x1,x2, FUN="myProbX")
 persp(x1, x2, probX, theta = 120, phi = 20, zlim=c(0,.3), zlab="density", xlab="beta_0",ylab=expression(beta_1), col="gray90", ticktype="detailed")
\end{Sinput}
\end{Schunk}
