Manuseio de dados com o tidyverse

Aula 04 · Introdução à Linguagem R

Diogo B. Provete

PPG em Ecologia e Conservação · UFMS

A estratégia split-apply-combine

Wickham, H. (2011). The split-apply-combine strategy for data analysis. Journal of Statistical Software 40(1): 1–29. doi:10.18637/jss.v040.i01

Dividir, aplicar, recombinar

sac T anuros 225 linhas G1 B. faber 46 T->G1 split G2 D. minutus 59 T->G2 G3 … mais 3 T->G3 R1 38,44 G1->R1 apply R2 33,25 G2->R2 R3 … G3->R3 F 5 linhas especie · ctmax_med R1->F combine R2->F R3->F

group_by() faz o split, summarise() faz o apply e o combine — e você nunca vê os grupos separados.

O mesmo, em código

anuros %>%
  group_by(especie) %>%                       # split
  summarise(ctmax_med = mean(ctmax, na.rm = TRUE))   # apply + combine

O equivalente no R base seria um for sobre as espécies, um mean() dentro, um objeto crescendo a cada volta e um rbind() no fim. Quatro lugares para errar em vez de um.

É por isso que o artigo de 2011 tem 29 páginas para uma ideia que cabe em duas linhas: a ideia é simples, escrever isso certo à mão é que não era.

O ciclo da ciência de dados

ciclo cluster_ex Explorar IMP Importar ARR Arrumar (tidy) IMP->ARR TRA Transformar ARR->TRA COM Comunicar VIS Visualizar TRA->VIS MOD Modelar VIS->MOD MOD->COM MOD->TRA

Esquema baseado em Wickham, Çetinkaya-Rundel & Grolemund, R for Data Science (2ª ed.). A aula de hoje é a seta Arrumar; a de amanhã, o laço Explorar.

Onde estudar

Instalando e carregando

install.packages("tidyverse")
library(tidyverse)

O tidyverse é um guarda-chuva: carrega dplyr, tidyr, ggplot2, readr, purrr e mais alguns de uma vez. Hoje usamos os dois primeiros.

O dado da semana

De onde ele vem

Bovo, R.P.; Simon, M.N.; Provete, D.B.; Lyra, M.; Navas, C.A.; Andrade, D.V. (2023). Beyond Janzen’s Hypothesis: How Amphibians That Climb Tropical Mountains Respond to Climate Variation. Integrative Organismal Biology 5(1): obad009. doi:10.1093/iob/obad009

225 indivíduos, 5 espécies de anuros, 6 altitudes, 2 serras da Mata Atlântica. Limites térmicos, perda de água, massa e clima de cada sítio.

Dados reais e publicados. É o mesmo arquivo da aula 01 à aula 08.

O que tem dentro

dim(bruto)
[1] 225  26
names(bruto)[1:12]
 [1] "ID"                       "Species"                 
 [3] "Sex"                      "EWL_Ugcm2s1"             
 [5] "WU_Ugcm2s1"               "perc_hidration_after_EWL"
 [7] "Bodymass_g"               "CTmin"                   
 [9] "CTmax"                    "Tbr"                     
[11] "WT"                       "Mountain_Range"

Repare nos nomes: EWL_Ugcm2s1, Bodymass_g, Altitude_m. A unidade está no nome da coluna — ótimo para quem abre o arquivo solto, cansativo para quem vai digitar aquilo cinquenta vezes hoje.

Verbos do dplyr

Os cinco que resolvem quase tudo

  • rename() muda o nome de uma coluna
  • select() escolhe colunas
  • filter() escolhe linhas por um critério — o filtro do Excel
  • mutate() cria colunas novas a partir das que existem
  • summarise() reduz muitas linhas a um valor-resumo
  • arrange() ordena as linhas

select mexe em coluna, filter mexe em linha. Trocar os dois é o erro mais comum da primeira semana.

rename() — tirar a unidade do nome

anuros <- bruto %>%
  rename(especie  = Species,     sexo     = Sex,
         ewl      = EWL_Ugcm2s1, massa    = Bodymass_g,
         ctmin    = CTmin,       ctmax    = CTmax,
         serra    = Mountain_Range, altitude = Altitude_m)
names(anuros)[1:8]
[1] "ID"                       "especie"                 
[3] "sexo"                     "ewl"                     
[5] "WU_Ugcm2s1"               "perc_hidration_after_EWL"
[7] "massa"                    "ctmin"

A unidade não sumiu — mudou de lugar. Ela agora vive no comentário do script e no arquivo de metadados. Encurtar nome sem fazer isso é perder informação.

select() — ficar só com o que interessa

anuros <- anuros %>%
  select(id = ID, especie, sexo, serra, altitude,
         massa, ewl, ctmin, ctmax, bio5 = BIO_5)
head(anuros, 4)
       id           especie sexo        serra altitude  massa    ewl
1  RPB 41 Rhinella icterica Male Serra do Mar      820 218.88 1.9732
2  RPB 42 Rhinella icterica Male Serra do Mar      820 233.12 2.0827
3  RPB 69 Rhinella icterica <NA> Serra do Mar      820  63.65 2.7657
4 RPB 114 Rhinella icterica Male Serra do Mar      820  85.05 1.5170
  ctmin ctmax bio5
1   3.4  39.5 25.3
2   4.2  39.3 25.3
3   3.2  40.6 25.3
4   1.8  39.0 25.3

Repare: select() também renomeia (bio5 = BIO_5). Um verbo a menos.

filter() — escolher linhas

anuros %>% filter(especie == "Boana faber") %>% nrow()
[1] 46
anuros %>% filter(altitude >= 1000, !is.na(ctmax)) %>% nrow()
[1] 106
anuros %>% filter(is.na(sexo)) %>% nrow()
[1] 107

Cento e sete sem sexo registrado, de 225. Não é erro de digitação: nem sempre dá para sexar o animal. Mas é uma variável que vocês não vão conseguir usar como preditora.

Vírgula dentro do filter() significa e. Para “ou”, use |.

mutate() — criar coluna

anuros <- anuros %>%
  mutate(amplitude = ctmax - ctmin,      # amplitude térmica tolerada
         tol_aquec = ctmax - bio5)       # tolerância ao aquecimento
summary(anuros$tol_aquec)
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.     NAs 
   0.60    8.50   10.80   10.69   13.70   15.30      35

O arquivo já trazia essas duas prontas. Recalculamos assim mesmo — não por desconfiança, mas porque agora a definição está escrita no código, e não só num arquivo de metadados que o leitor pode nunca abrir.

O pipe %>%

# sem pipe: lê-se de dentro para fora
head(arrange(filter(anuros, altitude > 1000), ctmax), 3)

# com pipe: lê-se na ordem em que acontece
anuros %>% filter(altitude > 1000) %>% arrange(ctmax) %>% head(3)
       id               especie sexo        serra altitude  massa    ewl
1 RPB 285 Dendropsophus minutus <NA> Serra do Mar     1500 0.5588 1.6490
2 RPB 286 Dendropsophus minutus Male Serra do Mar     1500 0.6915 1.7245
3 RPB 287 Dendropsophus minutus Male Serra do Mar     1500 0.5292 1.7974
  ctmin ctmax bio5 amplitude tol_aquec
1   6.4  28.7 23.8      22.3  4.900001
2   5.1  29.7 23.8      24.6  5.900001
3   3.3  29.7 23.8      26.4  5.900001

Leia o %>% como “e então”. O R moderno também aceita |>, nativo, que faz o mesmo nos casos simples.

group_by() + summarise()

O split-apply-combine do Wickham, em duas linhas:

anuros %>%
  group_by(especie) %>%
  summarise(n = n(),
            ctmax_med = round(mean(ctmax, na.rm = TRUE), 2),
            ctmax_dp  = round(sd(ctmax,   na.rm = TRUE), 2))
# A tibble: 5 × 4
  especie                   n ctmax_med ctmax_dp
  <chr>                 <int>     <dbl>    <dbl>
1 Boana faber              46      38.4     0.78
2 Dendropsophus minutus    59      33.2     1.7 
3 Leptodactylus latrans    31      38.8     0.62
4 Physalaemus cuvieri      39      35.1     1.61
5 Rhinella icterica        50      38.8     0.68

Agora vocês

10 minutos — encadeie os verbos

Escreva um pipeline, do anuros até a resposta, para cada pergunta:

  1. Qual espécie tem o maior CTmax médio?
  2. Quantos indivíduos há por serra, contando só os que têm ctmax medido?
  3. Qual a massa mediana de cada espécie em cada altitude?

Regra: um %>% por linha. Quando der erro, apague a última linha e rode de novo — assim você descobre em qual verbo quebrou, em vez de adivinhar.

A armadilha do .groups

anuros %>%
  group_by(especie, serra) %>%
  summarise(n = n(), .groups = "drop") %>%
  head(4)
# A tibble: 4 × 3
  especie               serra                    n
  <chr>                 <chr>                <int>
1 Boana faber           Serra da Mantiqueira    14
2 Boana faber           Serra do Mar            32
3 Dendropsophus minutus Serra da Mantiqueira    10
4 Dendropsophus minutus Serra do Mar            49

Sem .groups = "drop", o resultado sai ainda agrupado por especie, e a próxima operação age dentro de cada espécie sem avisar. É a fonte silenciosa de resultado errado mais comum do dplyr.

O desenho amostral, em uma linha

anuros %>% count(serra, altitude)
                 serra altitude  n
1 Serra da Mantiqueira      550 18
2 Serra da Mantiqueira     1600 55
3         Serra do Mar       35 39
4         Serra do Mar      820 36
5         Serra do Mar     1022 54
6         Serra do Mar     1500 23

Leiam com atenção: cada altitude aparece em uma serra só.

Isso não é defeito do arquivo — é como a coleta foi feita. E vai decidir quais modelos são possíveis na aula 06. Guardem.

Verbos do tidyr

O que é um dado tidy

  1. Cada variável é uma coluna
  2. Cada observação é uma linha
  3. Cada unidade observacional é uma tabela

Wickham, H. (2014). Tidy Data. Journal of Statistical Software 59(10): 1–23. doi:10.18637/jss.v059.i10

pivot_longer() — de largo para longo

Nosso arquivo tem ctmin e ctmax em duas colunas. Para o ggplot, muitas vezes queremos uma coluna “limite” e uma coluna “valor”:

longo <- anuros %>%
  select(id, especie, altitude, ctmin, ctmax) %>%
  pivot_longer(cols = c(ctmin, ctmax),
               names_to = "limite", values_to = "temperatura")
head(longo, 4)
# A tibble: 4 × 5
  id     especie           altitude limite temperatura
  <chr>  <chr>                <int> <chr>        <dbl>
1 RPB 41 Rhinella icterica      820 ctmin          3.4
2 RPB 41 Rhinella icterica      820 ctmax         39.5
3 RPB 42 Rhinella icterica      820 ctmin          4.2
4 RPB 42 Rhinella icterica      820 ctmax         39.3

De 225 linhas para 450: cada indivíduo agora ocupa duas.

pivot_wider() — o caminho de volta

longo %>%
  pivot_wider(names_from = limite, values_from = temperatura) %>%
  head(3)
# A tibble: 3 × 5
  id     especie           altitude ctmin     ctmax    
  <chr>  <chr>                <int> <list>    <list>   
1 RPB 41 Rhinella icterica      820 <dbl [1]> <dbl [1]>
2 RPB 42 Rhinella icterica      820 <dbl [1]> <dbl [1]>
3 RPB 69 Rhinella icterica      820 <dbl [1]> <dbl [1]>

“Largo” e “longo” não são certo e errado — são formatos convenientes para propósitos diferentes. Largo para ler na tela; longo para o ggplot e para a maioria dos modelos.

Os outros verbos do tidyr

  • separate() / unite() — dividir e juntar colunas (útil quando alguém gravou "Boana faber" e "macho" na mesma célula)
  • drop_na() — remove as linhas com ausentes
  • replace_na() / fill() — substituem ou propagam

Cuidado com drop_na() sem argumento: ele apaga a linha se qualquer coluna tiver NA. Aqui isso jogaria fora os 107 indivíduos sem sexo, mesmo numa análise que nem usa sexo. Diga em qual coluna: drop_na(ctmax).

Onde isso dá errado de verdade

nrow(anuros)
[1] 225
nrow(tidyr::drop_na(anuros))
[1] 91
nrow(tidyr::drop_na(anuros, ctmax))
[1] 190

O primeiro drop_na() custou mais da metade da amostra para nada. Silêncio total: nenhum aviso, nenhum erro. Só um n menor lá na frente.

Fechando o ciclo

Gravar o resultado

dir.create("dados/processados", showWarnings = FALSE, recursive = TRUE)
write.csv(anuros, "dados/processados/anuros_limpo.csv", row.names = FALSE)
  • O dado bruto nunca é sobrescrito. Fica em dados/, como veio da fonte.
  • O dado processado é descartável. Vive em dados/processados/, que entra no .gitignore — porque o script que o gera está versionado, e isso basta para reconstruí-lo.

Esse arquivo é o que as aulas 05, 06 e 07 vão ler. Rodem o aula-04.R antes delas.

Cheat sheets

Sempre na versão mais recente, direto da fonte:

Posit Software, PBC · CC BY-SA

Exercícios

Com o dado da semana — abram o aula-04.R:

  1. Quantas espécies ocorrem em ambas as serras? Responda com dplyr, não contando na tabela.
  2. Monte um resumo de ewl por serra e sexo. O que fazer com a linha dos NA?
  3. amplitude tem 35 ausentes, os mesmos de ctmax. Mostre com código que são os mesmos indivíduos.

No livro Análises Ecológicas no R, capítulo 5: leiam o tópico 5.6 e façam do 5.7.6 ao 5.8.17.