SOCR ≫ DSPA ≫ Topics ≫

This DSPA section Appendix.3.3 (Brain Surfaces) is part of the DSPA Appendix on visualizaiton of geometric and parametric surfaces. This DSPA Appendix (3) covers the following 3 topics:

1 Brain Surfaces

Often times, we are interested in modeling the boundaries of complex 3D solids using 2D manifolds. One of these interesting examples is reconstructing the mammalian cortical brain surface, e.g., rodents or humans. Below, we demonstrate two examples of loading precomputed, i.e., triangulated, surface models into R/RStudio for rendering and visualization using plotly. You can also see the SOCR Brian Viewer, which provides a more complex WebGL mechanism for visualization and interrogation of complex brain surfaces, models and stereotactic volumes.

1.1 Human Brain

There a many different software tools that can be used to model the human brain cortical surface using 3D sMRI data, e.g., Freesurfer and BrainSuite. Below we demonstrate an example using a Freesurfer surface model using the FS pial surface file format.

rm(list=ls())  # clean the environment
library(plotly)
library(geomorph)
# install.packages("freesurfer")
# https://github.com/joebathelt/ZDHHC9_connectome
# FS Data/Surfaces: http://www.freesurfer.net/pub/data/tutorial_data/long-tutorial/OAS2_0185_MR1.long.OAS2_0185/surf/

library(freesurfer)

fs <- freesurfer_read_surf('https://socr.umich.edu/data/DSPA/surfaces/lh.pial')
## [1] "created by nicks on Tue Aug 14 13:47:54 2007"
## [1] ""
colnames(fs$vertices) <- c('x', 'y', 'z')
# str(fs)
x <- fs$vertices[ , 'x']
y <- fs$vertices[ , 'y']
z <- fs$vertices[ , 'z']
  
mat <- matrix(c(x,y,z), ncol=3, dimnames=list(NULL,c("x","y","z")))

# plot_ly(x = x, y = y, z = z, type = "scatter3d", 
#        mode = "markers", marker = list(size = 3))

# now figure out the colormap for each 2-cell (face) by averaging
# the z-coordinates, mat[ , 3], of all 3 vertices of the face
zmean <- apply(t(fs$faces),MARGIN=2,function(row){mean(mat[row,3])})
# length(zmean); length(mesh$it[1,]) # equal?
# str(zmean); str(mat[,3]); str(mesh) # check structures
 
library(scales)
# color each 2-cell (face), ideally, this would be cortical thickness or another 4D values (e.g., p-value, curvature, etc.)
facecolor = colour_ramp(
  colorRampPalette(c("light blue", "gold"))(5)
)(rescale(x=zmean))

#facecolor = colour_ramp(
#   colorRampPalette(c("navy blue", "gold"))(5)
# )(rescale(x=fs$faces[,4]))
#length(fs$faces[,4])

plot_ly(
  x = x, y = y, z = z,
  i = fs$faces[ , 1], j = fs$faces[ , 2], k = fs$faces[ , 3],
  facecolor = facecolor,
  type = "mesh3d",
  opacity = 0.7,  # change opacity to make see through the surface
  contour=list(show=TRUE, color="#000", width=15)
)

1.2 Rodent Brain

This final example illustrates a simpler model of a rodent brain represented as a PLY surface model.

# library("plotly")
plyFile <- 'https://socr.umich.edu/data/DSPA/surfaces/Mouse_Brain.ply'
dest <- basename(plyFile)
if (!file.exists(dest)) {
  download.file(plyFile, dest)
}

# For ASCII PLY data: 
mesh <- geomorph::read.ply(dest, ShowSpecimen = F)

# For Binary PLY data use: mesh <- vcgPlyRead(dest, updateNormals = TRUE, clean = TRUE)

# see getS3method("shade3d", "mesh3d") for details on how to plot 

# plot point cloud
x <- mesh$vb["xpts",]
y <- mesh$vb["ypts",]
z <- mesh$vb["zpts",]
mat <- matrix(c(x,y,z), ncol=3, dimnames=list(NULL,c("x","y","z")))

# plot_ly(x = x, y = y, z = z, type = "scatter3d", 
#        mode = "markers", marker = list(size = 3))

# now figure out the colormap for each 2-cell (face) by averaging
# the z-coordinates, mat[ , 3], of all 3 vertices of the face
zmean <- apply(t(mesh$it),MARGIN=1,function(row){mean(mat[row,3])})
# length(zmean); length(mesh$it[1,]) # equal?
# str(zmean); str(mat[,3]); str(mesh) # check structures
 
library(scales)
# color each 2-cell (face)
facecolor = colour_ramp(
  colorRampPalette(c("light blue", "gold"))(5)
)(rescale(x=zmean))

plot_ly(
  x = x, y = y, z = -z,  # flip in z-direction to show anatomically
  i = mesh$it[1,]-1, j = mesh$it[2,]-1, k = mesh$it[3,]-1, # JS indices start at zero
  facecolor = facecolor,
  type = "mesh3d",
  opacity = 0.8,
  contour=list(show=TRUE, color="#000", width=15)
)

2 References

SOCR Resource Visitor number Web Analytics SOCR Email