Skip to content
Prev 85675 / 398506 Next

is there a way to visualize 3D normal distributions?

Dear Duncan, Michael, and Ben,

The code to which I referred Michael works by deforming a sphere. Here's the
central function:

ellipsoid <- function(center=c(0, 0, 0), radius=1, shape=diag(3), n=30){
# adapted from the shapes3d demo in the rgl package
  degvec <- seq(0, 2*pi, length=n)
  ecoord2 <- function(p) c(cos(p[1])*sin(p[2]), sin(p[1])*sin(p[2]),
cos(p[2]))
  v <- t(apply(expand.grid(degvec,degvec), 1, ecoord2))
  v <- center + radius * t(v %*% chol(shape))
  v <- rbind(v, rep(1,ncol(v))) 
  e <- expand.grid(1:(n-1), 1:n)
  i1 <- apply(e, 1, function(z) z[1] + n*(z[2] - 1))
  i2 <- i1 + 1
  i3 <- (i1 + n - 1) %% n^2 + 1
  i4 <- (i2 + n - 1) %% n^2 + 1
  i <- rbind(i1, i2, i4, i3)
  qmesh3d(v, i)
  }

Here are the key lines from scatter3d() for the shape and radius (slightly
edited to remove the context):

	dfn <- 3
	dfd <- length(x) - 1
      radius <- sqrt(dfn * qf(level, dfn, dfd))
      ellips <- ellipsoid(center=c(mean(x), mean(y), mean(z)), 
            shape=cov(cbind(x,y,z)), radius=radius)
      shade3d(ellips, col=surface.col[1], alpha=0.1, lit=FALSE)
      wire3d(ellips, col=surface.col[1], lit=FALSE)

Regards,
 John

--------------------------------
John Fox
Department of Sociology
McMaster University
Hamilton, Ontario
Canada L8S 4M4
905-525-9140x23604
http://socserv.mcmaster.ca/jfox 
--------------------------------