Basic 3D scatter plots

library("scatterplot3d") # load
data(iris)
head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa
scatterplot3d(iris[,1:3])

scatterplot3d(iris[,1:3], angle = 55)

Change the main title and axis labels

scatterplot3d(iris[,1:3],
              main="3D Scatter Plot",
              xlab = "Sepal Length (cm)",
              ylab = "Sepal Width (cm)",
              zlab = "Petal Length (cm)")

# Change the shape and the color of points

scatterplot3d(iris[,1:3], pch = 16, color="steelblue")

Change point shapes by groups

shapes = c(16, 17, 18) 
shapes <- shapes[as.numeric(iris$Species)]
scatterplot3d(iris[,1:3], pch = shapes)

# Change point colors by groups

colors <- c("#999999", "#E69F00", "#56B4E9")
colors <- colors[as.numeric(iris$Species)]
scatterplot3d(iris[,1:3], pch = 16, color=colors)

# Change the global appearance of the graph

Remove the box around the plot

scatterplot3d(iris[,1:3], pch = 16, color = colors,
              grid=TRUE, box=FALSE)

Add grids on the different factes of scatterplot3d graphics:

# 1. Source the function
source('http://www.sthda.com/sthda/RDoc/functions/addgrids3d.r')
# 2. 3D scatter plot
scatterplot3d(iris[, 1:3], pch = 16, grid=FALSE, box=FALSE)
# 3. Add grids
addgrids3d(iris[, 1:3], grid = c("xy", "xz", "yz"))

Add bars

scatterplot3d(iris[,1:3], pch = 16, type="h", 
              color=colors)

# Add legends

Specify the legend position using xyz.convert()

s3d <- scatterplot3d(iris[,1:3], pch = 16, color=colors)
legend(s3d$xyz.convert(7.5, 3, 4.5), legend = levels(iris$Species),
      col =  c("#999999", "#E69F00", "#56B4E9"), pch = 16)

Specify the legend position using keywords

# "right" position
s3d <- scatterplot3d(iris[,1:3], pch = 16, color=colors)
legend("right", legend = levels(iris$Species),
      col =  c("#999999", "#E69F00", "#56B4E9"), pch = 16)

# Use the argument inset
s3d <- scatterplot3d(iris[,1:3], pch = 16, color=colors)
legend("right", legend = levels(iris$Species),
  col = c("#999999", "#E69F00", "#56B4E9"), pch = 16, inset = 0.1)

# "bottom" position
s3d <- scatterplot3d(iris[,1:3], pch = 16, color=colors)
legend("bottom", legend = levels(iris$Species),
      col = c("#999999", "#E69F00", "#56B4E9"), pch = 16)

# Customize the legend position

# Custom point shapes
s3d <- scatterplot3d(iris[,1:3], pch = shapes)
legend("bottom", legend = levels(iris$Species),
       pch = c(16, 17, 18), 
      inset = -0.25, xpd = TRUE, horiz = TRUE)

# Custom colors
s3d <- scatterplot3d(iris[,1:3], pch = 16, color=colors)
legend("bottom", legend = levels(iris$Species),
      col =  c("#999999", "#E69F00", "#56B4E9"), pch = 16, 
      inset = -0.25, xpd = TRUE, horiz = TRUE)

# Custom shapes/colors
s3d <- scatterplot3d(iris[,1:3], pch = shapes, color=colors)
legend("bottom", legend = levels(iris$Species),
      col =  c("#999999", "#E69F00", "#56B4E9"), 
      pch = c(16, 17, 18), 
      inset = -0.25, xpd = TRUE, horiz = TRUE)

# Add point labels

scatterplot3d(iris[,1:3], pch = 16, color=colors)
text(s3d$xyz.convert(iris[, 1:3]), labels = rownames(iris),
     cex= 0.7, col = "steelblue")

# Add regression plane and supplementary points

data(trees)
head(trees)
##   Girth Height Volume
## 1   8.3     70   10.3
## 2   8.6     65   10.3
## 3   8.8     63   10.2
## 4  10.5     72   16.4
## 5  10.7     81   18.8
## 6  10.8     83   19.7

3D scatter plot with the regression plane:

# 3D scatter plot
s3d <- scatterplot3d(trees, type = "h", color = "blue",
    angle=55, pch = 16)
# Add regression plane
my.lm <- lm(trees$Volume ~ trees$Girth + trees$Height)
s3d$plane3d(my.lm)
# Add supplementary points
s3d$points3d(seq(10, 20, 2), seq(85, 60, -5), seq(60, 10, -10),
    col = "red", type = "h", pch = 8)

# Load plot3D package

library("plot3D")
## Warning: package 'plot3D' was built under R version 4.3.2

Prepare the data

data(iris)
head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa

Scatter plots

3D fancy Scatter plot with small dots on basal plane

# Add small dots on basal plane and on the depth plane
scatter3D_fancy <- function(x, y, z,..., colvar = z)
  {
   panelfirst <- function(pmat) {
      XY <- trans3D(x, y, z = rep(min(z), length(z)), pmat = pmat)
      scatter2D(XY$x, XY$y, colvar = colvar, pch = ".", 
              cex = 2, add = TRUE, colkey = FALSE)
   
      XY <- trans3D(x = rep(min(x), length(x)), y, z, pmat = pmat)
      scatter2D(XY$x, XY$y, colvar = colvar, pch = ".", 
              cex = 2, add = TRUE, colkey = FALSE)
  }
  scatter3D(x, y, z, ..., colvar = colvar, panel.first=panelfirst,
    colkey = list(length = 0.5, width = 0.5, cex.clab = 0.75)) 
}

Regression plane

data(mtcars)
head(mtcars[, 1:6])
##                    mpg cyl disp  hp drat    wt
## Mazda RX4         21.0   6  160 110 3.90 2.620
## Mazda RX4 Wag     21.0   6  160 110 3.90 2.875
## Datsun 710        22.8   4  108  93 3.85 2.320
## Hornet 4 Drive    21.4   6  258 110 3.08 3.215
## Hornet Sportabout 18.7   8  360 175 3.15 3.440
## Valiant           18.1   6  225 105 2.76 3.460
# x, y, z variables
x <- mtcars$wt
y <- mtcars$disp
z <- mtcars$mpg
# Compute the linear regression (z = ax + by + d)
fit <- lm(z ~ x + y)
# predict values on regular xy grid
grid.lines = 26
x.pred <- seq(min(x), max(x), length.out = grid.lines)
y.pred <- seq(min(y), max(y), length.out = grid.lines)
xy <- expand.grid( x = x.pred, y = y.pred)
z.pred <- matrix(predict(fit, newdata = xy), 
                 nrow = grid.lines, ncol = grid.lines)
# fitted points for droplines to surface
fitpoints <- predict(fit)
# scatter plot with regression plane
scatter3D(x, y, z, pch = 18, cex = 2, 
    theta = 20, phi = 20, ticktype = "detailed",
    xlab = "wt", ylab = "disp", zlab = "mpg",  
    surf = list(x = x.pred, y = y.pred, z = z.pred,  
    facets = NA, fit = fitpoints), main = "mtcars")

# text3D: plot 3-dimensionnal texts

data(USArrests)
with(USArrests, text3D(Murder, Assault, Rape, 
  labels = rownames(USArrests), colvar = UrbanPop, 
  col = gg.col(100), theta = 60, phi = 20,
  xlab = "Murder", ylab = "Assault", zlab = "Rape", 
  main = "USA arrests", cex = 0.6, 
  bty = "g", ticktype = "detailed", d = 2,
  clab = c("Urban","Pop"), adj = 0.5, font = 2))

# text3D and scatter3D

# Plot texts
with(USArrests, text3D(Murder, Assault, Rape, 
  labels = rownames(USArrests), colvar = UrbanPop, 
  col = gg.col(100), theta = 60, phi = 20,
  xlab = "Murder", ylab = "Assault", zlab = "Rape", 
  main = "USA arrests", cex = 0.6, 
  bty = "g", ticktype = "detailed", d = 2,
  clab = c("Urban","Pop"), adj = 0.5, font = 2))
# Add points
 with(USArrests, scatter3D(Murder, Assault, Rape - 1, 
    colvar = UrbanPop, col = gg.col(100), 
    type = "h", pch = ".", add = TRUE))

# Zoom near origin: choose suitable ranges
 plotdev(xlim = c(0, 10), ylim = c(40, 150), 
         zlim = c(7, 25))

# display axis ranges
 getplist()[c("xlim","ylim","zlim")] 
## $xlim
## [1]  0.8 17.4
## 
## $ylim
## [1]  45 337
## 
## $zlim
## [1]  7.3 46.0

3D Histogram

data(VADeaths)
#  hist3D and ribbon3D with greyish background, rotated, rescaled,...
hist3D(z = VADeaths, scale = FALSE, expand = 0.01, bty = "g", phi = 20,
        col = "#0072B2", border = "black", shade = 0.2, ltheta = 90,
        space = 0.3, ticktype = "detailed", d = 2)

hist3D (x = 1:5, y = 1:4, z = VADeaths,
        bty = "g", phi = 20,  theta = -60,
        xlab = "", ylab = "", zlab = "", main = "VADeaths",
        col = "#0072B2", border = "black", shade = 0.8,
        ticktype = "detailed", space = 0.15, d = 2, cex.axis = 1e-9)
# Use text3D to label x axis
 text3D(x = 1:5, y = rep(0.5, 5), z = rep(3, 5),
       labels = rownames(VADeaths),
       add = TRUE, adj = 0)
# Use text3D to label y axis
 text3D(x = rep(1, 4),   y = 1:4, z = rep(0, 4),
       labels  = colnames(VADeaths),
       add = TRUE, adj = 1)

scatter2D: 2D scatter plot

# x, y coordinates
set.seed(1234)
x  <- sort(rnorm(10)) 
y  <- runif(10)
# Variable for coloring points
col.v <- sqrt(x^2 + y^2) 

Basic 2D scatter plot:

scatter2D(x, y, colvar = col.v, pch = 16, bty ="n",
          type ="b")

# 2D scatter plot with confidence interval:

# Confidence interval for x variable only
CI <- list()
CI$x <- matrix(nrow = length(x), data = c(rep(0.25, 2*length(x))))
scatter2D(x, y, colvar = col.v, pch = 16, bty ="n", cex = 1.5, 
          CI = CI, type = "b")

# Confidence interval for both x and y variables
CI$y <- matrix (nrow = length(y), data = c(rep(0.05, 2*length(y))))
CI$col <- "black"
scatter2D(x, y, colvar = col.v, pch = 16,  bty ="n", cex = 1.5,
          CI = CI, type ="b")

CI$y[c(2,4,8,10), ] <- NA  # Some points have no CI
CI$x[c(2,4,8,10), ] <- NA  # Some points have no CI
CI$alen <- 0.02            # increase arrow head
scatter2D(x, y, colvar = col.v, pch = 16,  bty ="n", cex = 1.5,
          CI = CI, type ="b")

# text2D

# Only text
with(USArrests, text2D(x = Murder, y = Assault + 5, colvar = Rape, 
     xlab = "Murder", ylab = "Assault", clab = "Rape", 
     main = "USA arrests", labels = rownames(USArrests), cex = 0.6, 
     adj = 0.5, font = 2))

# text with point
 with(USArrests, text2D(x = Murder, y = Assault + 5, colvar = Rape, 
     xlab = "Murder", ylab = "Assault", clab = "Rape", 
     main = "USA arrests", labels = rownames(USArrests), cex = 0.6, 
     adj = 0.5, font = 2))
 with(USArrests, scatter2D(x = Murder, y = Assault, colvar = Rape, 
     pch = 16, add = TRUE, colkey = FALSE))

x0 <- c(0, 0, 0, 0)
y0 <- c(0, 0, 0, 0)
z0 <- c(0, 0, 0, 0)
x1 <- c(0.89, -0.46, 0.99, 0.96)
y1 <- c(0.36,  0.88, 0.02, 0.06)
z1 <- c(-0.28, 0.09, 0.05, 0.24)
cols <- c("#1B9E77", "#D95F02", "#7570B3", "#E7298A")

3D Arrows:

arrows3D(x0, y0, z0, x1, y1, z1, colvar = x1^2, col = cols,
         lwd = 2, d = 3, clab = c("Quality", "score"), 
         main = "Arrows 3D", bty ="g", ticktype = "detailed")
# Add starting point of arrow
points3D(x0, y0, z0, add = TRUE, col="darkred", 
          colkey = FALSE, pch = 19, cex = 1)
# Add labels to the arrows
text3D(x1, y1, z1, c("Sepal.L", "Sepal.W", "Petal.L", "Petal.W"),
       colvar = x1^2, col = cols, add=TRUE, colkey = FALSE)

# 2D arrows:

arrows2D(x0, y0,  x1, y1,  colvar = x1^2, col = cols,
         lwd = 2, clab = c("Quality", "score"), 
          bty ="n", xlim = c(-1, 1), ylim = c(-1, 1))
# Add vertical and horizontal lines at c(0,0)
abline(h =0, v = 0, lty = 2)
# Add starting point of arrow
points2D(x0, y0, add = TRUE, col="darkred", 
          colkey = FALSE, pch = 19, cex = 1)
# Add labels to the arrows
text2D(x1, y1, c("Sepal.L", "Sepal.W", "Petal.L", "Petal.W"),
       colvar = x1^2, col = cols, add=TRUE, colkey = FALSE)

# 3D rectangle

rect3D(x0 = 0, y0 = 0.5, z0 = 0, x1 = 1, z1 = 5, 
       ylim = c(0, 1), bty = "g", facets = TRUE, 
       border = "red", col ="#7570B3", alpha=0.5,
       lwd = 2, phi = 20)

# 2D rectangle:

rect2D(x0 = runif(3), y0 = runif(3), 
       x1 = runif(3), y1 = runif(3), colvar = 1:3, 
       alpha = 0.4, lwd = 2, main = "rect2D")

library("plot3Drgl")
## Warning: package 'plot3Drgl' was built under R version 4.3.2
## Loading required package: rgl
## Warning: package 'rgl' was built under R version 4.3.2
library("car")
## Warning: package 'car' was built under R version 4.3.2
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.3.2

Prepare the data

data(iris)
head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa
rgl.postscript("plot.pdf",fmt="pdf")
## Warning in rgl.postscript("plot.pdf", fmt = "pdf"): Postscript conversion
## failed

Install the RGL package

install.packages("rgl")
## Warning: package 'rgl' is in use and will not be installed

Load the RGL package

library("rgl")

Prepare the data

data(iris)
head(iris)
##   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## 1          5.1         3.5          1.4         0.2  setosa
## 2          4.9         3.0          1.4         0.2  setosa
## 3          4.7         3.2          1.3         0.2  setosa
## 4          4.6         3.1          1.5         0.2  setosa
## 5          5.0         3.6          1.4         0.2  setosa
## 6          5.4         3.9          1.7         0.4  setosa
x <- sep.l <- iris$Sepal.Length
y <- pet.l <- iris$Petal.Length
z <- sep.w <- iris$Sepal.Width

Start and close RGL device

rgl.open(): Opens a new device

rgl.close(): Closes the current device

rgl.clear(): Clears the current device

rgl.cur(): Returns the active device ID

rgl.quit(): Shutdowns the RGL device system

3D scatter plot

Basic graph

rgl.open() # Open a new RGL device
## Warning: 'rgl.open' is deprecated.
## Use 'open3d' instead.
## See help("Deprecated")
rgl.points(x, y, z, color ="lightgray") # Scatter plot
## Warning: 'rgl.points' is deprecated.
## Use 'points3d' instead.
## See help("Deprecated")

Change the background and point colors

rgl.open()# Open a new RGL device
## Warning: 'rgl.open' is deprecated.
## Use 'open3d' instead.
## See help("Deprecated")
rgl.bg(color = "white") # Setup the background color
## Warning: 'rgl.bg' is deprecated.
## Use 'bg3d' instead.
## See help("Deprecated")
rgl.points(x, y, z, color = "blue", size = 5) # Scatter plot
## Warning: 'rgl.points' is deprecated.
## Use 'points3d' instead.
## See help("Deprecated")
rgl.open()# Open a new RGL device
## Warning: 'rgl.open' is deprecated.
## Use 'open3d' instead.
## See help("Deprecated")
rgl.bg(color = "white") # Setup the background color
## Warning: 'rgl.bg' is deprecated.
## Use 'bg3d' instead.
## See help("Deprecated")
rgl.spheres(x, y, z, r = 0.2, color = "grey") 

Scale the data

x1 <- (x - min(x))/(max(x) - min(x))
y1 <- (y - min(y))/(max(y) - min(y))
z1 <- (z - min(z))/(max(z) - min(z))

Use c(-max, max)

lim <- function(x){c(-max(abs(x)), max(abs(x))) * 1.1}

rgl_add_axes(): A custom function to add x, y and z axes

# x, y, z : numeric vectors corresponding to
#  the coordinates of points
# axis.col : axis colors
# xlab, ylab, zlab: axis labels
# show.plane : add axis planes
# show.bbox : add the bounding box decoration
# bbox.col: the bounding box colors. The first color is the
# the background color; the second color is the color of tick marks
rgl_add_axes <- function(x, y, z, axis.col = "grey",
                xlab = "", ylab="", zlab="", show.plane = TRUE, 
                show.bbox = FALSE, bbox.col = c("#333377","black"))
  { 
  
  lim <- function(x){c(-max(abs(x)), max(abs(x))) * 1.1}
  # Add axes
  xlim <- lim(x); ylim <- lim(y); zlim <- lim(z)
  rgl.lines(xlim, c(0, 0), c(0, 0), color = axis.col)
  rgl.lines(c(0, 0), ylim, c(0, 0), color = axis.col)
  rgl.lines(c(0, 0), c(0, 0), zlim, color = axis.col)
  
   # Add a point at the end of each axes to specify the direction
   axes <- rbind(c(xlim[2], 0, 0), c(0, ylim[2], 0), 
                 c(0, 0, zlim[2]))
   rgl.points(axes, color = axis.col, size = 3)
  
  # Add axis labels
  rgl.texts(axes, text = c(xlab, ylab, zlab), color = axis.col,
             adj = c(0.5, -0.8), size = 2)
  
  # Add plane
  if(show.plane) 
    xlim <- xlim/1.1; zlim <- zlim /1.1
    rgl.quads( x = rep(xlim, each = 2), y = c(0, 0, 0, 0),
             z = c(zlim[1], zlim[2], zlim[2], zlim[1]))
  
  # Add bounding box decoration
  if(show.bbox){
    rgl.bbox(color=c(bbox.col[1],bbox.col[2]), alpha = 0.5, 
          emission=bbox.col[1], specular=bbox.col[1], shininess=5, 
          xlen = 3, ylen = 3, zlen = 3) 
  }
}

Change the color of points by groups

#' Get colors for the different levels of 
#' a factor variable
#' 
#' @param groups a factor variable containing the groups
#'  of observations
#' @param colors a vector containing the names of 
#   the default colors to be used
get_colors <- function(groups, group.col = palette()){
  groups <- as.factor(groups)
  ngrps <- length(levels(groups))
  if(ngrps > length(group.col)) 
    group.col <- rep(group.col, ngrps)
  color <- group.col[as.numeric(groups)]
  names(color) <- as.vector(groups)
  return(color)
}