He encontrado un antiguo trabajo en el que necesité crear un gráfico de dispersión entre dos variables numéricas y que en los márgenes de dicho gráfico aparecieran los boxplots de cada variable por separado. Los gráficos están construidos con las funciones gráficas básicas adaptando los márgenes de cada uno con el parámetro fig de la función par() y añadiendo cada gráfico al anterior con el parámetro new = TRUE. La función mtext() sirve para añadir texto en los márgenes de un gráfico y en los boxplots elimino los ejes con axes = TRUE. En el siguiente ejemplo utilizo los datos mtcars:

par(fig = c(0, 0.8, 0, 0.8), new = TRUE)
plot(mtcars$wt, mtcars$mpg, xlab = "Car Weight", ylab = "Miles Per Gallon")
par(fig = c(0, 0.8, 0.55, 1), new = TRUE)
boxplot(mtcars$wt, horizontal = TRUE, axes = FALSE)
par(fig = c(0.65, 1, 0, 0.8), new = TRUE)
boxplot(mtcars$mpg, axes = FALSE)
mtext("Scatterplot with marginal boxplots", side = 3, outer = TRUE, line = -3)

Aunque en este código se podrían cambiar los boxplots marginales por histogramas o densidades, hoy en día existen librerías que simplifican la tarea de tener que programarlos. La librería {ggplot2} tiene una serie de extensiones (https://exts.ggplot2.tidyverse.org/gallery/) que añaden más opciones al ya amplio ecosistema de gráficos creados con dicha librería.

Por ejemplo, la librería {ggExtra} permite añadir gráficos marginales a los scatterplots creados con {ggplot2}:

# install.packages(c('ggplot2', 'ggExtra'))
library(ggplot2)
library(ggExtra)
# Creamos el gráfico de dispersión
p <- ggplot(data = mtcars, aes(x = wt, y = mpg)) + geom_point() + xlab("Car Weight") + ylab("Miles Per Gallon")
# Boxplots marginales
ggMarginal(p, type = 'boxplot')

# Histogramas marginales
ggMarginal(p, type = 'histogram', fill = 'grey', bins = 8)

# Densidades marginales y por categorías
q <- ggplot(data = mtcars, aes(x = wt, y = mpg, colour = factor(cyl))) + geom_point() +
  labs(x = "Car Weight", y = "Miles Per Gallon", colour = 'Cylinders')
ggMarginal(q, type = 'density', groupColour = TRUE, groupFill = TRUE)