Axa y scalată dublă cu ggplot()

De multe ori m-am trezit nevoit să trasez un singur grafic cu două axe y având scale diferite. De exemplu, acest lucru ar putea fi util pentru reprezentarea datelor despre temperatură și precipitații într-o anumită locație. Din păcate, făcând asta cu ggplot() nu este simplu.

În timp ce căutam soluții posibile, am descoperit că graficele cu axe y duale nu sunt, în general, considerate instrumente bune pentru vizualizarea datelor. Hadley Wickham, autorul ggplot2 și alte câteva pachete R importante, a dat câteva motive întemeiate pentru aceasta (de exemplu, la acest link). Nu intenționez să pun la îndoială aceste argumente generale. Cu toate acestea, ideea cu care nu mă pot descurca ggplot() ceva ce aș putea face cu ușurință cu Excel mă deranjează. La urma urmei, astfel de grafice pot fi destul de utile în anumite circumstanțe specifice, cum ar fi atunci când se afișează date meteo.

Prin urmare, am început să caut o soluție rezonabilă și în cele din urmă mi-am găsit drumul. Aș dori să o împărtășesc în această postare.

Să luăm în considerare setul de date DailyMeteoData.csvcare conține date medii zilnice de temperatură și precipitații timp de doi ani într-o locație din regiunea mea (Umbria, Italia centrală). Să deschidem setul de date și să folosim dplyr pentru a face câteva transformări utile, cum ar fi:

1. conversia șirului de caractere date într-un obiect dată;
2. adăugarea de variabile pentru Anul, Luna și Ziua Anului (DOY).

library(dplyr)
fileName <- "https://www.casaonofri.it/_datasets/DailyMeteoData.csv"
dataMeteo <- read.csv(fileName) |>
  mutate(Date = as.Date(Date, format = "%d/%m/%Y"),
         Year = as.numeric(format(Date, format="%Y")),
         Month = as.numeric(format(Date, format="%m")),
         DOY = as.numeric(format(Date, format="%j")))
head(dataMeteo)

        Date Tavg Rain Year Month DOY
1 2011-01-01  5.9  0.0 2011     1   1
2 2011-01-02  6.2  0.6 2011     1   2
3 2011-01-03  4.5  0.0 2011     1   3
4 2011-01-04 -0.3  0.0 2011     1   4
5 2011-01-05  2.3  0.2 2011     1   5
6 2011-01-06  6.9  0.0 2011     1   6

Temperatura este o variabilă continuă și poate fi reprezentată cu ușurință folosind un grafic cu linii, în timp ce precipitațiile constau din evenimente discrete și, prin urmare, sunt acumulate mai util în perioade precum zece zile sau o lună. Putem folosi dplyr încă o dată pentru a crea un nou set de date privind precipitațiile care să conțină precipitațiile lunare acumulate, care este mai potrivit pentru scopurile noastre. Din motive de simplitate, presupunem că anul este împărțit în 12 luni de lungime egală (aproximativ 30,4 zile) și calculăm DOY-ul corespunzător zilei centrale a fiecărei luni, care va fi folosit ca centru al barei grafice respective.

dataMeteo2 <- dataMeteo |>
  group_by(Year, Month) |>
  summarise(Rain = sum(Rain)) |>
  mutate(DOY = seq(15, 365, by = 365/12))
head(dataMeteo2)

# A tibble: 6 × 4
# Groups:   Year (1)
   Year Month  Rain   DOY
     
1  2011     1  40.4  15  
2  2011     2  32.8  45.4
3  2011     3 113.   75.8
4  2011     4  16.6 106. 
5  2011     5  27.4 137. 
6  2011     6  61.2 167. 

Acum putem produce primul nostru grafic care arată atât precipitațiile, cât și temperatura. În caseta de mai jos, m-am jucat puțin cu bifoanele și etichetele axei x și am lăsat eticheta axa y goală pentru moment.

library(ggplot2)
ggplot() +
  geom_bar(dataMeteo2, mapping = aes(x = DOY, y = Rain), fill = "grey",
           stat = "identity", width = 28) +
  geom_line(dataMeteo, mapping = aes(x = DOY, y = Tavg), col = "blue") +
  scale_x_continuous(breaks = c(365/12, 365/12*4, 365/12*8, 365) - 15, 
                     labels = c("Jan", "Apr", "Aug", "Dec"),
                     name = "") +
  scale_y_continuous(name = "") +
  facet_wrap(~Year) +
  theme_bw()

Graficul anterior nu funcționează foarte bine: linia temperaturii este greu vizibilă deoarece scara sa de măsurare este mult mai mică decât cea a precipitațiilor. Prin urmare, trebuie să „scalăm” variabila de temperatură, astfel încât aceasta să varieze aproximativ de la 150 la 250. Acest lucru va face linia albastră vizibilă clar, fără a interfera cu barele de ploaie. Valorile temperaturii minime și maxime sunt -4,2°C, respectiv 28,9°C; dorim să mapam aceste valori originale la noile valori 150 și, respectiv, 250, așa cum se arată în figura de mai jos.

Figura anterioară ne spune că am putea transforma scala de temperatură inițială utilizând ecuația dreptei care trece prin cele două puncte (-4,2, 150) și (28,9, 250). Datorită a ceea ce ne amintim din cursurile de geometrie, putem calcula panta unei drepte ca:

în timp ce interceptarea este:

Astfel, ecuația de scalare este:

unde este noua scară de temperatură, în timp ce este cel original. Transformarea inversă este:

Acum, suntem gata să trasăm graficul. Transformăm temperatura la noua scară și o diagramăm; în plus, includem a doua axă utilizând sec.axis argument și cel sec_axis() funcția, în care specificăm funcția de transformare înapoi la scara originală.

dataMeteo <- dataMeteo %>% 
  mutate(newTavg = 162.69 + 3.02 * Tavg)

ggplot() +
  geom_bar(dataMeteo2, mapping = aes(x = DOY, y = Rain), fill = "grey",
           stat = "identity", width = 28) +
  geom_line(dataMeteo, mapping = aes(x = DOY, y = newTavg), col = "blue") +
  scale_x_continuous(breaks = c(365/12, 365/12*4, 365/12*8, 365) - 15, 
                     labels = c("Jan", "Apr", "Aug", "Dec"),
                     name = "") +
  scale_y_continuous(name = "Rain (mm)", 
                     sec.axis = sec_axis(~ (. - 162.69)/3.02, 
                                         name = "Daily Temperature (°C)",
                                         breaks = c(-10, 0, 10, 20, 30))) +
  facet_wrap(~Year) +
  theme_bw()

Și am terminat!

library(lubridate)
library(ggh4x)
dataMeteo1 <- dataMeteo |>
  mutate(
    Date = as.Date(Date, format = "%d/%m/%Y"),
    Year = lubridate::year(Date),
    Month = lubridate::month(Date))

dataMeteo2 <- dataMeteo1 |>
  group_by(Year, Month) |>
  summarise(
    Rain = sum(Rain),
    Temp = mean(Tavg, rm.na = TRUE)) |>
  mutate(Mth = factor(month.abb(Month), 
                      levels = c("Jan", "Feb", "Mar", "Apr",
                                 "May", "Jun", "Jul", "Aug",
                                 "Sep", "Oct", "Nov", "Dec")))
dataMeteo3 <- dataMeteo2 |>
  mutate(newTemp = 162.69 + 3.02 * Temp)

ggplot(dataMeteo3, mapping = aes(x = interaction(Mth, Year),
                                 group = 1)) +
  geom_col(aes(y = Rain), fill = "grey") +
  geom_point(aes(y = newTemp), size = 2, col = "blue") +
  geom_smooth(
    method = "gam", formula = y ~ s(x, bs = "cc"),
    aes(x = as.numeric(interaction(Mth, Year)), y = newTemp),
    col = "blue") +
  geom_point(aes(y = newTemp), size = 3, col = "blue", 
             fill = "white", shape = 21, stroke = 1) +
  scale_y_continuous(
    minor_breaks = scales::breaks_width(20),
    name = "Mean Monthly Total Rainfall (mm)",
    sec.axis = sec_axis(~ (. - 162.69) / 3.02,
      name = "Mean Monthly Daily Temperature (°C)",
      breaks = c(-10, 0, 10, 20, 30))) +
  guides(x = "axis_nested",
        y = guide_axis(minor.ticks=TRUE),
        y.sec = guide_axis(minor.ticks=TRUE)) +
  theme_bw(base_size = 18) +
  theme(
     panel.grid.minor = element_blank(),
     axis.title.x = element_blank(),
     axis.text.x = element_text(face = "bold", size = rel(0.5), angle = 90),
     axis.ticks = element_line(colour = "red"),
     ggh4x.axis.nestline.x = element_line(linewidth = 0.6),
     ggh4x.axis.nesttext.x = element_text(colour = "blue", 
                                          face = "bold", 
                                          size = rel(1.0))
  )