Happy soil conservation day!
A paleosol (buried and residual ancient soils) sequence in Cerro Seco, Bogotá.
Some context here · algo de contexto
I happen to be working on a conservation proposal for an area at the southern Bogotá boundary - known as Cerro Seco (google it! 😉) -, which has been afected by sand mining for decades. This is one of the few relicts of subxerophytic ecosystems in the region, and has been little studied -if at all 😤 - due to the socio-environmental conflicts there. We work in close collaboration with the local community through the stages of data collection and management planning. A key aspect here is communication: how to make soil (or paleosol in this case) data more accessible through visualization? Here I am sharing some tricks on how to plot visually appealing profiles using ggplot2. It will require some familiarity with the package. Let’s see the required data structure and packages first!
Trabajo actuamente en una propuesta de conservación para un área en el sur de Bogotá - conocida como Cerro Seco -, la cuál ha sido explotada por décadas por la minería de arena. Es uno de los últimos relictos de ecosistémas subxerofíticos de la región y ha sido muy poco estudiado. Debido a acceso limitado, trabajamos en colaboración con la comunidad local para la recolección de datos y la planeación del uso del terreno. Un aspecto clave aquí es la comunicación: ¿cómo hacer los datos de suelo (o paleosuelos en este caso) mas accesibles a través de la visualización?. Aquí comparto algunos trucos para crear visualizaciones atractivas de perfiles usando ggplot2. Requiere algo del conocimiento de esta librería. Primero, veamos la estructura de los datos y las librerías requeridas!
library(aqp)
library(ggplot2)
library(ggrepel)
library(colorspace)
library(tidyverse)
horizons <- readr::read_csv('https://raw.githubusercontent.com/cmguiob/TCI_CerroSeco_git/main/Datos/Suelos_CS_Horiz.csv')
site <- readr::read_csv('https://raw.githubusercontent.com/cmguiob/TCI_CerroSeco_git/main/Datos/Suelos_CS_Sitio.csv')
#Select four profiles and relevant properties for plot
hz4 <- horizons %>%
dplyr::filter(ID %in% c("CS01", "CS02","CS03","CS04")) %>%
dplyr::select(ID, BASE, TOPE, ESP, HZ, CON_POR, MX_H, MX_V, MX_C, CON_H, CON_V, CON_C )
head(hz4, 10)## # A tibble: 10 x 12
## ID BASE TOPE ESP HZ CON_POR MX_H MX_V MX_C CON_H CON_V CON_C
## <chr> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <dbl> <dbl> <chr> <dbl> <dbl>
## 1 CS01 35 0 35 A 0 10YR 4 2 2.5Y 5 8
## 2 CS01 80 35 45 E 2 10YR 5 1 2.5Y 5 8
## 3 CS01 90 80 10 Bt 80 10YR 4 1 10YR 2 1
## 4 CS01 140 90 50 B/C 10 10YR 5 1 10YR 2 1
## 5 CS01 170 140 30 Bn1 5 10YR 5 1 10YR 2 1
## 6 CS01 177 170 7 Bt1 80 2.5Y 4 1 10YR 2 1
## 7 CS01 250 177 73 C1 5 10YR 7 8 7.5YR 5 8
## 8 CS01 280 250 30 R 0 10YR 8 3 7.5YR 5 8
## 9 CS02 25 0 25 A 0 7.5YR 4 2 10YR 7 8
## 10 CS02 50 25 25 R 0 7.5YR 8 2 7.5YR 5 8For the plot, the relevant variables are profile ID, depth of base and top of horizons (BASE, TOPE), horizon thickness (ESP), horizon name (HZ), percentage of concentrations (CON_POR), and Munsell color coordinates for the matrix and the concentrations (MX_x, CON_x).
Para la gráfica, las variables relevantes son IDdel perfil, profundidad de la base y el tope (BASE, TOPE), espesor de lso horizontes (ESP), nomenclatura del horizonte (HZ), porcentaje de concetraciones (CON_POR), y colores Munsell de la matriz y concentraciones (MX_x, CON_x).
A basic type of graph · un tipo de gráfica básico
Plotting soil profiles is not much different from plotting bar charts, but it takes some tinkering. If you are short of time, functions for plotting soil profiles - and other soil stuff - have been already incorporated in the aqp package by Dylan Beaudette, which has a nice introduction.
So, why bother with ggplot2? Well, it has great compatibility among geom_ family functions, which allow you to create complex plots with aesthetic freedom, it is well structured, extensively documented and transparent - great for debugging! -, and it can be easily incorporated in data workflows with tidyverse functions (spoiler: pipes %>% coming soon).
¿Por qué tomarse la molestia de usar ggplot2 para graficar perfiles, si ya existen funciones para esto? Por ejemplo la librería aqp, que tiene una genial introducción por Dylan Beadette. Pues, ggplot2 permite gran compatibilidad entre las funciones geom_, con lo cual tenemos libertad estética para crear gráficas complejas; es bien estructurado, detalladamente documentado y transparente - genial para depurar errores! - , y puede incorporarse fácilmente en los flujos de trabajo basados en tidyverse (spoiler: pipes %>% a continuación).
# Create color variables | crear variables de color
hz4$RGBmx <- munsell2rgb(hz4$MX_H, hz4$MX_V, hz4$MX_C)
hz4$RGBco <-munsell2rgb(hz4$CON_H,hz4$CON_V , hz4$CON_C)
# Create factor variable with ID and HZ | Crear variable factor con ID y HZ
hz_bdf <- hz4 %>%
dplyr::select(ID, BASE, TOPE, ESP,HZ, CON_POR, RGBmx, RGBco)%>%
dplyr::mutate(ID_HZ = paste(ID, HZ),
ID_HZ2 = factor(ID_HZ, ID_HZ))
head(hz_bdf, 5)## # A tibble: 5 x 10
## ID BASE TOPE ESP HZ CON_POR RGBmx RGBco ID_HZ ID_HZ2
## <chr> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <chr> <chr> <fct>
## 1 CS01 35 0 35 A 0 #6F5F4CFF #9C750FFF CS01 A CS01 A
## 2 CS01 80 35 45 E 2 #83796FFF #9C750FFF CS01 E CS01 E
## 3 CS01 90 80 10 Bt 80 #686057FF #38302AFF CS01 Bt CS01 Bt
## 4 CS01 140 90 50 B/C 10 #83796FFF #38302AFF CS01 B/C CS01 B/C
## 5 CS01 170 140 30 Bn1 5 #83796FFF #38302AFF CS01 Bn1 CS01 Bn1We just transformed the Munsell colors to computer-readable colors (HEX) and stored them in the RGBmx and RGBco variables. The new factor variable ID_HZ2 will be useful to assign fill colors in geom_bar. This is the function that we will tinker all the way down to the final plot.
Acabamos de transformar los colores Munsell a colores legibles por el computador (HEX) y los guardamos en las variables RGBmx y RGBco. La nueva variable tipo factor ID_HZ2 será util para asignar los colores de relleno en geom_bar. Esta es la función que estaremos ajustando hasta el final.
ggplot(hz_bdf, aes(x = ID, y = ESP, fill = ID_HZ2)) +
geom_bar(position="dodge", stat="identity") 
This is the bar chart you are probably most familiar with. Here, each bar represents a horizon and the height its thickness.The bars are grouped by ID, we will keep it that way. You see the default colors are not very soil-like, we will change them next, along with a more useful arrangement of bars: the stack, as position argument.
Este es el diagrama de barras con el que probablemente tengas mas familiaridad. Aquí, cada barra representa un horizote y su altura representa el espesor. Las barras están agrupadas por ID, vamos a mantenerlas así. Como ves los colores por defecto no son de la “paleta de suelos”; vamos a cambiarla a continuación, junto con una configuración más útil de las barras: el stack (apilado), como argumento position.
profiles <- ggplot(hz_bdf, aes(x = ID, y = ESP, fill = ID_HZ2)) +
geom_bar(position="stack", stat="identity") +
# Adds colors from the RGBmx variable
scale_fill_manual(values = hz_bdf$RGBmx,
# Don't plot the fill legend | no grafique leyenda de relleno
guide = "none")
profiles
Now they look pretty much as I remember them 😏. Though you might be wondering: what kind of horizons are those? And, aren’t soil profiles supposed to be dug from the top downwards? We will address these issues next, plus we will modify the bars’ width (aesthetics tip: make some space in your plot!).
Ahora se ven casi como los recuerdo. Sin embargo, quizás te preguntas por el tipo de horizontes y … ¿no se supone que los perfiles de suelo se cavan de la superficie hacia abajo?. Vamos a corregir esto, además de cambiar el ancho de las barras (tip de estética: crea algo de espacio!).
profiles2 <- ggplot(hz_bdf,
aes(x = ID, y = ESP, fill = forcats::fct_rev(ID_HZ2))) +
geom_bar(position="stack", stat="identity", width = 0.4) +
# Adds colors from the RGBmx variable
scale_fill_manual(values = rev(hz_bdf$RGBmx),
# Don't plot the fill legend | no grafique leyenda de relleno
guide = "none") +
# Adds horizon labels | agrega etiquetas a horizontes
ggrepel::geom_text_repel( data = hz_bdf,
aes(y = BASE - (ESP/2), label = HZ),
color = darken(hz_bdf$RGBmx, .2, space = "HCL"),
size = 3,
family = "robotoc",
hjust = 0,
direction = "y",
nudge_x = 0.3,
segment.size = .5,
segment.linetype = "dotted",
segment.square = TRUE,
segment.curvature = 0.1,
segment.angle = 30,
segment.alpha = 0.5,
box.padding = 0.3) +
# Reverse y axis scale | Invierte la escala del eje y
scale_y_reverse(breaks = c(0,100,200,300,400,500),
labels=c("0", "100", "200", "300", "400", "500"))+
scale_x_discrete(position = "top")
profiles2
Adding complexity · agregando complejidad
Still remember the dataset? There is a color variable that we haven´t used yet, the color of concentrations: these are redox or illuvaition coatings. We will create random points that represent these concentrations using the CON_POR: the percentage of concentrations in each horizon, and will encode them using geom_jitter.
Aún recuerdas los datos? Hay una variable de color que aún no hemos usado, el color de las concentraciones: estas son concentraciones tipo redox o recubrimientos de iluviación. Vamos a crear puntos aleatorios que representen estas concentraciones usando CON_POR: el porcentaje de concentraciones por horizonte, y las vamos a graficar usando geom_jitter.
# New data frame with random points | Nuevo data frame con puntos aleatorios
hz_jdf <- hz4 %>%
dplyr::select(ID, BASE, TOPE, ESP,HZ, CON_POR, RGBmx, RGBco)%>%
dplyr::mutate(ID = factor(ID),
ID_HZ = paste(ID, HZ),
ID_HZ2 = factor(ID_HZ, ID_HZ))%>%
dplyr::mutate(CON_POR = ifelse(CON_POR == 0, 1, CON_POR),
n = 5*ESP*CON_POR / 100,
mean = 0.5*(BASE - TOPE),
sd = 0.1*ESP)%>%
dplyr::mutate(samples = purrr::pmap(.[11:13], rnorm))%>%
tidyr::unnest(samples)
head(hz_jdf, 20)## # A tibble: 20 x 14
## ID BASE TOPE ESP HZ CON_POR RGBmx RGBco ID_HZ ID_HZ2 n mean
## <fct> <dbl> <dbl> <dbl> <chr> <dbl> <chr> <chr> <chr> <fct> <dbl> <dbl>
## 1 CS01 35 0 35 A 1 #6F5F4~ #9C75~ CS01~ CS01 A 1.75 17.5
## 2 CS01 80 35 45 E 2 #83796~ #9C75~ CS01~ CS01 E 4.5 22.5
## 3 CS01 80 35 45 E 2 #83796~ #9C75~ CS01~ CS01 E 4.5 22.5
## 4 CS01 80 35 45 E 2 #83796~ #9C75~ CS01~ CS01 E 4.5 22.5
## 5 CS01 80 35 45 E 2 #83796~ #9C75~ CS01~ CS01 E 4.5 22.5
## 6 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 7 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 8 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 9 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 10 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 11 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 12 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 13 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 14 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 15 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 16 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 17 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 18 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 19 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## 20 CS01 90 80 10 Bt 80 #68605~ #3830~ CS01~ CS01 ~ 40 5
## # ... with 2 more variables: sd <dbl>, samples <dbl>
See how the data frame expanded? For each horizon we created some random points with a mean = BASE - TOPE/2, a standard deviation sd= 0.1*ESP and a number of points equal to n = 5*CON_POR*ESP/100 (if you have questions about it 😵, contact me). Ah, the plot! Here it comes.
¿Ves como se expandió el data frame? Para cada horizonte creamos puntos aleatorios con mean = BASE - TOPE/2 , desviación estándar sd= 0.1*ESP y numero de puntos n = 5*CON_POR*ESP/100 (si tienes preguntas de esto, escríbeme). Ah, la gráfica! Ahí viene.
profiles3 <- profiles2 +
#location from where jitter spreads out vertically
geom_jitter(data = hz_jdf, aes(x = ID, y = BASE - (ESP/2)),
width = 0.18,
# how far jitter spreads out to each side
height = hz_jdf$ESP*0.5,
size = 0.3,
col = hz_jdf$RGBco,
shape = 16)+
scale_y_reverse(breaks = c(0,100,200,300,400,500),
labels=c("0", "100", "200", "300", "400", "500\ncm"))+
scale_x_discrete(position = "top")
profiles3
What follows is related to my personal taste. After some adjustments of the theme, here is my final version. You could add some other elements to this simple plot, such as other properties, processes, location, etc. I find the vast possibilities of data-based design pretty exciting!
Lo que sigue tiene que ver con mis gustos personales. Luego de algunos ajustes de theme, esta es mi versión final. A esta visualización simple se pueden agregar elementos que representen otras propiedades, procesos, localización, etc. Las varias posibilidades del diseño basado en datos son emocionantes!
