Revisão de purrr

Funções map

Vamos começar, como sempre, carregando as bibliotecas que precisaremos para o tutorial. Além de tidyverse, utilizaremos também o pacote purrr, sobre o qual aprenderemos um pouco mais neste tutorial. Instale-o se ainda não tiver feito.

library(tidyverse)
library(purrr)

Geraremos um data frame simples, com três variáveis numéricas, para utilizarmos no início do tutorial. Note no código abaixo que variável “x” terá aproximadamente 5% de missing values:

df <- tibble(x = sample(c(rnorm(95), rep(NA, 5)), 
                            1000, 
                            replace = T),
                 y = rnorm(1000, 42, 3),
                 z = rexp(1000, 1))

Há várias funções para calcularmos a média de todas as colunas de um data frame e em diversos pacotes. Neste tutorial estamos interessados em entender como podemos fazer esse calculo utilizando, em primeiro lugar, for loops, e, a seguir, as funções do pacote purrr que, basicamente, substituem o for loops em R.

Uma possível solução com for loops é:

media_col <- c()

for (i in 1:ncol(df)){
  media_col <- c(media_col, mean(df[[i]], na.rm = T))
}

No exemplo acima, criamos um objeto vazio e a cada iteração do loop “anexamos” a média da i-ésima variável. Uma alternativa ligeiramente mais eficiente é determinar antes do loop o tamanho do objeto que conterá o resultado e a cada i-ésima itearação adicionar a média calculada à i-ésima posição deste objeto:

media_col <- vector("double", ncol(df))

for (i in 1:ncol(df)){
  media_col[[i]] <- mean(df[[i]], na.rm = T)
}

Em ambos os casos temos um bocado de código para uma tarefa bastante simples. A partir de agora, esqueceremos o for loop e aprenderemos um pouco sobre as funções da família “map” que são utlizadas em programação funcional como substitutas de loops.

O uso básico do pacote purrr, para resolver o problema acima é:

map(df, mean)

A função map pode ser aplicada a uma lista, a um data frame (que é uma lista de vetores na vertical e pareados) ou a um vetor e, para cada posição da lista/data frame/vetor aplica uma função. O resultado é uma lista de tais aplicações. Note que o resultado da função acima é uma lista.

As variações da função map permite que o resultado seja outra classe de objeto do que uma lista. map_dbl, por exemplo, retorna um vetor do tipo “double”. Diferentemente de map, as funções da família que especificam a classe de output precisam que a função aplicada a cada elemento também retorne tal classe.

map_dbl(df, mean)

A variável “x” contém, propositalmenete, NAs. Vamos removê-los do cálculo da média. Para adicionar parâmatros à função mean, que aplicamos com map_dbl, basta separá-los com vírgula, em ordem, após a função aplicada. Eles seriam ‘passados’ para a função mean. Para removermos os NAs do cálculo da média fazemos:

map_dbl(df, mean, na.rm = T)

Para calcularmos o desvio padrão:

map_dbl(df, sd, na.rm = T)

Ou o primeiro decil:

map_dbl(df, quantile, probs = 0.1, na.rm = T)

ou o primeiro e último decil:

map_dbl(df, quantile, probs = c(0.1, 0.9), na.rm = T)

Ops! Aqui deu errado. O que aconteceu?

A função map_dbl requer que o output seja um vetor atômico do tipo double. Nesse caso, ao aplicarmos a função quantile e escolhermos 2 decis, temos como resultado um vetor de tamanho 2. O resultado da iteração não pode, por essa razão, ser um vetor. Pode ser, porém, uma lista:

map(df, quantile, probs = c(0.1, 0.9), na.rm = T)

A função map sempre servirá onde alguma outra da família servir. Na dúvida, use-a.

Vamos a outro exemplo. Peguemos o data frame “starwas”, disponível no pacote datasets que venha com R.

data("starwars")

O data frame contém informações sobre cada um dos personagens da série de filmes. Como é de se esperar, algumas variáveis deste dataset não devem ser numéricas. Como investigamos quais variáveis são numéricas?

map_lgl(starwars, is.numeric)

Simples, não? Com map_lgl aplicamos a função is.numeric a todas as variáveis e obtemos como resultado um vetor lógico que contém a informação de quais variávais são numéricas. Poderíamos utilizá-lo, por exemplo, para selecionar apenas as variáveis numéricas do dataset e calcularmos a média:

numerica <- map_lgl(starwars, is.numeric)
map_dbl(starwars[, numerica], mean, na.rm = T)

Algo semelhante pode ser feito com o “typeof”, que retorna o tipo da variável. Este informação, porém, é textual, e por isso utilizamos a função map_chr para aplicá-la a todos os elementos do data frame (ou lista, ou vetor):

map_chr(starwars, typeof)

As funções aplicadas à uma lista (ou data frame ou vetor) não precisam, como vimos no exemplo da função quantiles ter outputs “simples”:

map(df, summary)

Podemos, inclusive, criar funções e aplicá-las com map. Vamos resgatar o data frame “mtcars”:

data(mtcars)

amplitude <- function(x) {
  max(x) - min(x) 
}

map(mtcars, amplitude)

Se quisermos, podemos utilizar funções anônimas, ou seja, que não são geradas como objeto, mas sim criadas no contexto de sua aplicação, o que poupa bastante cógido.

map(mtcars, function(x) max(x)- min(x))

Para tornar o código ainda mais curto e eliminar obviedades, podemos utilizar o símbolo “~” para indicar que haverá na sequência uma função, e “.” em substituição à variável correspondente ao parâmetro da função

map(mtcars, ~ max(.)- min(.))

Para encerrarmos o exemplo, poderíamos utilizar a função map_dbl e obter diretamente um vetor numérico em vez de uma lista:

map_dbl(mtcars, ~ max(.)- min(.))

Finalmente, gostamos de trabalhar com data.frames - eles são bastante intuitivos e podem ser reusados para vários fins subsequentes - tabelas, gráficos, etc. Podemos pedir para o map devolve um data.frame pronto para usar? Sim! Usamos, naturalmente, map_df:

map_df(mtcars, ~ max(.)- min(.))

map_df(mtcars, summary) 

map_df(mtcars, summary) %>% 
  mutate(Statistic=c("Min.", "1st Qu.", "Median", "Mean", "3rd Qu.", "Max.")) %>%
  select(Statistic, everything())

Exercicio

Usando o banco de dados starwars e a função map, idenitificar a pessoa mais alto, mais pesado, e mais jovem. Também identifique a pessoa mais baixo, mais leve, e mais velho. Finalmente, escreva uma função para calcular ambos o máximo e o mínimo ao mesmo tempo, e aplicá-lo com map.

Modelos com purrr

Vamos a um exemplo mais interessante, ainda com “mtcars”. Nosso objetivo será produzir modelos com pequenas variações a partir de um modelo de regressão linear bastante simples (função lm, “linear model”) entre o consumo de combustível de um carro (mpg, “miles per gallon”) e seu peso (wt, “weight”), abaixo apresentado:

reg <- lm(mpg ~ wt , data = mtcars)

O resultado reg da função lm é um objeto bastante complexo da classe lm:

class(reg)
str(reg)

No entanto, trata-se basicamente de uma lista. Podemos extrair seus componentes “navegando” a lista. Por exemplo, para extrairmos os coeficientes gerados no modelo podemos fazer:

reg$coefficients

ou ainda, como nos interessará adiante:

reg[['coefficients']]

ou mesmo utilizando a posição dos coeficientes na estrutura do objeto:

reg[[1]]

Vamos supor que nos interessa agora produzir um modelo seperado para cada categoria de carro em relação ao número de cilindros, que podem ser 4, 6 ou 8. Com a função split, geraremos uma lista de data frames em que cada posição da lista contém um data frame com as observações de uma categoria de número de cilindros.

lista_df <- split(mtcars, mtcars$cyl)

lista_df

Como gerar, sem repetir várias vezes o código, um modelo para cada um dos data frames na lista de data frames (lembrando que as variáveis dos data frames são exatamente as mesmas)?

Com a função map, obviamente:

resultados <- map(lista_df, function(x) lm(mpg ~ wt , data = x))

Se quisermos simplificar o código:

resultados <- map(lista_df, ~ lm(mpg ~ wt , data = .))

Examine o objeto “resultados”

resultados

Note que, para cada categoria de número de cilindros, temos agora um objeto da classe “lm”. Antes, tinhamos uma lista de data frames. Agora, temos uma lista de resultados da aplicação da função lm a data frames.

Mas o que nos interessa não é observar todo o resultado, mas apenas os coeficientes gerados nos modelos. Podemos, assim aplicar novamente a função map ao resultados:

map(resultados, function(x) x[["coefficients"]])

Por termos em cada posição da lista um objeto lm, que é uma lista, podemos simplificar o código acima:

map(resultados, "coefficients")

Vamos guardar os coeficientes em um objeto.

coeficientes <- map(resultados, "coefficients")

Cada modelo tem 2 coeficientes, o intercepto e o coeficiente para a variável “wt”. Procedendo exatamente como acabamos de fazer, podemos extrair apenas os coeficientes de “wt” de nossa lista:

map(coeficientes, "wt")

Por termos apenas vetores atômicos em cada posição da lista, poderíamos ter optado por map_dbl em vez de map:

map_dbl(coeficientes, "wt")

Finalmente, as funções map pode ser utilizadas com “pipe”. O código que produzimos acima pode ser condensado da seguinte maneira:

mtcars %>% 
  split(.$cyl) %>%
  map(~ lm(mpg ~ wt, data = .)) %>%
  map("coefficients") %>%
  map_dbl("wt")

Aproveite agora para examinar a documentação da função map. O exemplo que acabamos de fazer vem de lá:

?map

Vamos criar algo novo agora. Vamos supor que queremos produzir o mesmo modelo não para cada categoria de cilindrada, mas incluindo a cada vez apenas um quintil da variável qsec (tempo que o carro leva para percorrer o primeiro quarto de milha). Ou seja, em primeiro lugar incluiremos apenas 20% carros mais lentos, a seguir os 40% mais lentos, até gerarmos um modelo com todos os carros:

Para gerar os quintis podemos fazer:

quintis <- quantile(mtcars$qsec, probs = c(.2, .4, .6, .8, 1))

Vamos aplicar a função map ao vetor de quintis (e não ao uma lista de data frames, como anteriormente) para gerar o resultado que nos interessa. Exceto pelas duas primeiras linhas de código, as demais são idênticas à anterior:

quintis %>%
  map(~ filter(mtcars, qsec <= .)) %>%
  map(~ lm(mpg ~ wt, data = .)) %>%
  map("coefficients") %>%
  map_dbl("wt")

Note bem que começamos aqui com os quintis e não com o data.frame mtcars como anteriormente. Assim, podemos iterar sobre a variação nos dados que nos interesse.

Mais um exemplo: como aplicamos a mesma regressão para duas variáveis dependentes diferentes? Agora não queremos filtrar o banco de dados mas alterar uma variável no fórmula de regressão cada vez e deixar o banco de dados fixo.

dependentes <- c("mpg","qsec")

dependentes %>%
  map(~lm(as.formula(paste(.x," ~ wt")), data=mtcars)) %>%
  map("coefficients") %>%
  map_dbl("wt")

E se quisemos os “Standard Errors” de cada modelo? Podemos usar a função summary para acessar mais detalhes:

dependentes %>%
  map(~lm(as.formula(paste(.x," ~ wt")), data=mtcars)) %>%
  map(summary) %>%
  set_names(dependentes) %>%
  map("coefficients")

Uma função muito útil aqui para ‘limpar’ os modelos complexos se chama tidy do pacote broom

library(broom)

dependentes %>%
  map(~lm(as.formula(paste(.x," ~ wt")), data=mtcars)) %>%
  map(tidy)

Note que cada um dos nossos modelos venha como um data.frame separado. Mas podemos pedir para um único data.frame com map_df:

models_df <- dependentes %>%
  map(~lm(as.formula(paste(.x," ~ wt")), data=mtcars)) %>%
  set_names(dependentes) %>%
  map_df(tidy, .id="model")

Agora, com um data.frame flexível, temos possibilidades diversas. Por exemplo, podemos reorganizar a tabela para gerar um gráfico das linhas de regressão:

models_df %>%
  select(model, term, estimate) %>%
  spread(key="term", value="estimate") %>%
  rename("Intercept"=`(Intercept)`) %>%
  select(-model) %>%
  ggplot() +
  geom_abline(aes(intercept=Intercept, slope=wt)) +
  xlim(0,10) +
  ylim(0,50) +
  theme_classic()

Exercicio

Usando as funções de map, execute regressões com o banco de dados mtcars para a relação entre ‘hp’ (horsepower) e ‘qsec’ (a variável dependente) separadamente para carros com 3, 4 e 5 ‘gears’. Compare o coeficiente para a variável ‘hp’ para os três tipos de carros.

mtcars %>% 
  split(.$gear) %>%
  map(~lm(qsec~hp, data=.)) %>%
  map("coefficients") %>%
  map_dbl("hp")

Combinando data frames em uma lista

Até agora, vimos basicamente a função map do pacote purrr. Vamos examinar rapidamente outra função bastante útil do pacote: reduce.

Um das chatices de trabalharmos com os dados da MUNIC é o fato das informações estarem espalhadas em diversas planilhas de um mesmo documento.

Para o exemplo, temos 3 planilhas que foram exportadas em formato .txt. Em primeiro lugar, vamos abrí-las:

munic_v_externa <- read.csv("https://raw.githubusercontent.com/leobarone/FLS6397_2018/master/data/munic_v_externa.csv", sep=";")
munic_r_humanos <- read.csv("https://raw.githubusercontent.com/leobarone/FLS6397_2018/master/data/munic_r_humanos.csv", sep=";")
munic_p_urbano <- read.csv2("https://raw.githubusercontent.com/leobarone/FLS6397_2018/master/data/munic_p_urbano.csv")

e juntá-las em uma lista.

lista_munic <- list(munic_v_externa, munic_r_humanos, munic_p_urbano)

Nosso objetivo agora é combiná-las em um único data frame sem precisar fazer múltiplos “joins”. Se estívessemos, por exemplo, trabalhando com dados do TSE, teríamos que combinar 27 data frames em um único e escrever o código das combinações par a par seria pouco inteligente.

Em vez disso, vamos fazer um full_join das 3 planilhas da MUNIC que abrimos pelo id do município (variável A1, comum a todas as planilhas):

munic <- lista_munic %>%
  reduce(full_join, by = "A1")

Bastante simples, não? A função reduce serve para reduzirmos uma lista a um único objeto. Um exemplo mais simples de aplicação da função seria:

lista <- list(c(1, 1, 2, 3), c(5, 8), 13, c(21, 34))

lista %>% 
  map(sum) %>%
  reduce(`+`)

reduce também pode ser utilizado para vetores, e não apenas para listas, como veremos adiante.