Excel con R

1 Riferimento

Questo documento è stato redatto seguendo il libro R for Data Science.

Ci sono diversi pacchetti al di fuori del tidyverse e in questo libro ci si occupa di:

• nycflights13: dati sui voli aerei
• gapminder: dati su sviluppo mondiale
• Lahman: baseball

2 Ricevere aiuto

Quando si cerca aiuto si pensi sempre a fare un reprex con il codice interessanto e soprattutto un dput con le variabili interessate:

dput(mtcars)

## structure(list(mpg = c(21, 21, 22.8, 21.4, 18.7, 18.1, 14.3, 
## 24.4, 22.8, 19.2, 17.8, 16.4, 17.3, 15.2, 10.4, 10.4, 14.7, 32.4, 
## 30.4, 33.9, 21.5, 15.5, 15.2, 13.3, 19.2, 27.3, 26, 30.4, 15.8, 
## 19.7, 15, 21.4), cyl = c(6, 6, 4, 6, 8, 6, 8, 4, 4, 6, 6, 8, 
## 8, 8, 8, 8, 8, 4, 4, 4, 4, 8, 8, 8, 8, 4, 4, 4, 8, 6, 8, 4), 
##     disp = c(160, 160, 108, 258, 360, 225, 360, 146.7, 140.8, 
##     167.6, 167.6, 275.8, 275.8, 275.8, 472, 460, 440, 78.7, 75.7, 
##     71.1, 120.1, 318, 304, 350, 400, 79, 120.3, 95.1, 351, 145, 
##     301, 121), hp = c(110, 110, 93, 110, 175, 105, 245, 62, 95, 
##     123, 123, 180, 180, 180, 205, 215, 230, 66, 52, 65, 97, 150, 
##     150, 245, 175, 66, 91, 113, 264, 175, 335, 109), drat = c(3.9, 
##     3.9, 3.85, 3.08, 3.15, 2.76, 3.21, 3.69, 3.92, 3.92, 3.92, 
##     3.07, 3.07, 3.07, 2.93, 3, 3.23, 4.08, 4.93, 4.22, 3.7, 2.76, 
##     3.15, 3.73, 3.08, 4.08, 4.43, 3.77, 4.22, 3.62, 3.54, 4.11
##     ), wt = c(2.62, 2.875, 2.32, 3.215, 3.44, 3.46, 3.57, 3.19, 
##     3.15, 3.44, 3.44, 4.07, 3.73, 3.78, 5.25, 5.424, 5.345, 2.2, 
##     1.615, 1.835, 2.465, 3.52, 3.435, 3.84, 3.845, 1.935, 2.14, 
##     1.513, 3.17, 2.77, 3.57, 2.78), qsec = c(16.46, 17.02, 18.61, 
##     19.44, 17.02, 20.22, 15.84, 20, 22.9, 18.3, 18.9, 17.4, 17.6, 
##     18, 17.98, 17.82, 17.42, 19.47, 18.52, 19.9, 20.01, 16.87, 
##     17.3, 15.41, 17.05, 18.9, 16.7, 16.9, 14.5, 15.5, 14.6, 18.6
##     ), vs = c(0, 0, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 
##     0, 1, 1, 1, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 1), am = c(1, 
##     1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 
##     0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1), gear = c(4, 4, 4, 3, 
##     3, 3, 3, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 
##     3, 3, 4, 5, 5, 5, 5, 5, 4), carb = c(4, 4, 1, 1, 2, 1, 4, 
##     2, 2, 4, 4, 3, 3, 3, 4, 4, 4, 1, 2, 1, 1, 2, 2, 4, 2, 1, 
##     2, 2, 4, 6, 8, 2)), row.names = c("Mazda RX4", "Mazda RX4 Wag", 
## "Datsun 710", "Hornet 4 Drive", "Hornet Sportabout", "Valiant", 
## "Duster 360", "Merc 240D", "Merc 230", "Merc 280", "Merc 280C", 
## "Merc 450SE", "Merc 450SL", "Merc 450SLC", "Cadillac Fleetwood", 
## "Lincoln Continental", "Chrysler Imperial", "Fiat 128", "Honda Civic", 
## "Toyota Corolla", "Toyota Corona", "Dodge Challenger", "AMC Javelin", 
## "Camaro Z28", "Pontiac Firebird", "Fiat X1-9", "Porsche 914-2", 
## "Lotus Europa", "Ford Pantera L", "Ferrari Dino", "Maserati Bora", 
## "Volvo 142E"), class = "data.frame")

3 Data Visualisation

Attenzione ci sono un sacco di addins per ggplot in questa pagina web

Prima di tutto facciamo qualche grafico con ggplot. Ci chiediamo con il data set di mtcars se le auto con un motore più grande consumano più benzina. Tutta via non vogliamo un si o no bensì che tipo di relazione c’è al crescere del motore sul consumo di benzina.

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy)) +
  labs(x="Car engine in liters", y="Miles x gallon in high")

Con il pacchetto plotly posso fare dei grafici interattivi molto interessanti:

baseDati <- mtcars %>% group_by(cyl) %>% summarise(horsePower = mean(hp))
plotly::ggplotly(
  baseDati %>% ggplot(
  aes(x=cyl,y=horsePower))+
  geom_col(fill="#ba0c2f"))

3.1 Aesthetics

The greatest value of a picture is when it forces us to notice what we never expected to see. — John Tukey

Supponiamo un grafico che metta in relazione il consumo con il motore dell’auto come quello di sopra. Si vedono dei punti nell’angolo a destra in alto che escono fuori dalla linearità. Si può analizzare cosa siano aggiungendo una terza variabile:

ggplot(data = mpg) + 
  geom_point(
    mapping = aes(
      x = displ, 
      y = hwy, 
      #color = class, # Posso usare i colori
      #size = class   # oppure le dimensioni dei punti
      #alpha = class  # oppure la tonalità
      shape = class  # oppure una forma (ce ne sono solo 6!!) 
      )
    )+
  labs(x="Litri motore",y="Consumo autostradale")+
  scale_color_paletteer_d("RColorBrewer::RdGy")

## Warning: The shape palette can deal with a maximum of 6 discrete values because
## more than 6 becomes difficult to discriminate; you have 7. Consider
## specifying shapes manually if you must have them.

## Warning: Removed 62 rows containing missing values (geom_point).

Un modo differente per aggiungere una terza variabile, oltre ai colori, le dimensioni, gli alfa e le forme, potrebbe essere il facet. Sono sottografici che mi servono per mostrare un subset per ogni grafico.

• facet_wrap: utilizza una sola variabile. La variabile deve essere discreta.
• facet_grid: se, invece di una singola variabile, si voglio utilizzare due variabili

Una variabile:

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy)) + 
  facet_wrap(~ class, nrow = 2)

Due variabili:

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy)) + 
  facet_grid( drv ~ cyl)

Si possono fare un facet o per riga o per colonna in questo modo:

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy)) +
  facet_grid(drv ~ .)

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy)) +
  facet_grid(. ~ cyl)

Non tutti i geom accettano gli stessi argomenti. Per esempio gli shape funzionano con i punti ma non con le linee. D’altro canto puoi usare, con le linee, il linetype.

geom_smooth serve per fare una banda di accettazione intorno alla linea principale.

ggplot(data = mpg) + 
  geom_smooth(mapping = aes(x = displ, y = hwy, linetype = drv))

## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

e poi

library(ggpubr)
#> Loading required package: ggplot2
#> Loading required package: magrittr
# Create some data format
# :::::::::::::::::::::::::::::::::::::::::::::::::::
set.seed(1234)
wdata = data.frame(
   sex = factor(rep(c("F", "M"), each=200)),
   weight = c(rnorm(200, 55), rnorm(200, 58)))
head(wdata, 4)

##   sex   weight
## 1   F 53.79293
## 2   F 55.27743
## 3   F 56.08444
## 4   F 52.65430

#>   sex   weight
#> 1   F 53.79293
#> 2   F 55.27743
#> 3   F 56.08444
#> 4   F 52.65430

# Density plot with mean lines and marginal rug
# :::::::::::::::::::::::::::::::::::::::::::::::::::
# Change outline and fill colors by groups ("sex")
# Use custom palette
ggdensity(wdata, x = "weight",
   add = "mean", rug = TRUE,
   color = "sex", fill = "sex",
   palette = c("#00AFBB", "#E7B800"))

data(diamonds)
ggpie3D(data = diamonds, group_key = "cut", count_type = "full", tilt_degrees = -10, label_size=2) + 
  ggtitle("tilt_degrees = -10") + 
  theme(plot.title = element_text(hjust = 0.5))

3.2 Trasformazioni statistiche

Prendiamo in considerazione la base dati diamonds che comprende 53940.

Se noi facciamo un grafico come il seguente:

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut))

vediamo che alcuni grafici come i bar generano una variabile count mentre altri come scatterplot vogliono sia una variabile x sia una y.

I grafici che non richiedono esplicitamente una y fanno implicitamente un conteggio. Per questo motivo posso modificare il grafico di sopra in questo modo:

ggplot(data = diamonds) + 
  stat_count(mapping = aes(x = cut))

Ci sono tre ragioni per le quali posso fare un override della funzione di default:

1. Si vuole fare un override della statistica di default. Nel codice sottostante modifichiamo la statistica di default in identity. In questo modo l’altezza del grafico è esattamente quanto espresso da me.
2. Vuoi fare l’override del mapping di default
3. Vuoi porre attenzione alla trasformazione del codice. Si può usare stat_summary che fa il summarise il valore y per ciascun x.

Caso 1:

demo <- tribble(
  ~cut,         ~freq,
  "Fair",       1610,
  "Good",       4906,
  "Very Good",  12082,
  "Premium",    13791,
  "Ideal",      21551
)

ggplot(data = demo) +
  geom_bar(mapping = aes(x = cut, y = freq), stat = "identity")

Caso 2:

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut, y = stat(prop), group = 1))

Caso 3:

ggplot(data = diamonds) + 
  stat_summary(
    mapping = aes(x = cut, y = depth),
    fun.min = min,
    fun.max = max,
    fun = median
  )

3.3 Aggiustamento della posizione

Ci sono molte altre cose che posso fare con il bar chart. Per esempio puoi colorare o riempire le barre.

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut, colour = cut))

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut, fill = cut))

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut, fill = clarity))

Se non si vuole un chart impilato classico si possono usare una di queste tre opzioni:

• identity

• dodge

• fill

Identity: posiziona ogni oggetto esattamente dove cade

ggplot(data = diamonds, mapping = aes(x = cut, fill = clarity)) + 
  geom_bar(alpha = 1/5, position = "identity")

ggplot(data = diamonds, mapping = aes(x = cut, colour = clarity)) + 
  geom_bar(fill = NA, position = "identity")

fill: funziona come l’impilata ma crea barre tutte della medesima altezza

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut, fill = clarity), position = "fill")

dodge: crea più barre una vicina all’altra

ggplot(data = diamonds) + 
  geom_bar(mapping = aes(x = cut, fill = clarity), position = "dodge")

C’è un altro aggiustamento utile non tanto per i barchart bensì per gli scatter plot che si chiama jitter e serve per ridurre il overplotting. Questo fenomeno viene fuori quando ho tanti punti nella stessa posizione. In questo caso, jitter aggiunge pochissimo rumore sui vari punti in modo da rendere maggiormente visibile la situazione. Si vedano i due grafici sottostanti:

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy))

ggplot(data = mpg) + 
  geom_point(mapping = aes(x = displ, y = hwy), position = "jitter")

Per saperne di più: ?position_dodge, ?position_fill, ?position_identity, ?position_jitter, and ?position_stack

3.4 Sistemi di coordinate

E’ la parte più difficile ma qui introduciamo:

• coord_flip: per invertire x con y

• coord_quickmap: per impostare il ratio alle mappe

• coord_polar: per creare coordinate polari

# Flip
ggplot(data = mpg, mapping = aes(x = class, y = hwy)) + 
  geom_boxplot()

ggplot(data = mpg, mapping = aes(x = class, y = hwy)) + 
  geom_boxplot() +
  coord_flip()

# Quickmap
nz <- map_data("nz")

ggplot(nz, aes(long, lat, group = group)) +
  geom_polygon(fill = "white", colour = "black")

ggplot(nz, aes(long, lat, group = group)) +
  geom_polygon(fill = "white", colour = "black") +
  coord_quickmap()

# Polari
bar <- ggplot(data = diamonds) + 
  geom_bar(
    mapping = aes(x = cut, fill = cut), 
    show.legend = FALSE,
    width = 1
  ) + 
  theme(aspect.ratio = 1) +
  labs(x = NULL, y = NULL)

bar + coord_flip()

bar + coord_polar()

4 Data trasformation

La visualizzazione dei dati è importante ma è raro averli nel formato desiderato. Per modificare i dati, il pacchetto principe è dyplr.

Per far ciò utilizzeremo il dataset denominato flights.

head(flights)

## # A tibble: 6 × 19
##    year month   day dep_time sched_dep…¹ dep_d…² arr_t…³ sched…⁴ arr_d…⁵ carrier
##   <int> <int> <int>    <int>       <int>   <dbl>   <int>   <int>   <dbl> <chr>  
## 1  2013     1     1      517         515       2     830     819      11 UA     
## 2  2013     1     1      533         529       4     850     830      20 UA     
## 3  2013     1     1      542         540       2     923     850      33 AA     
## 4  2013     1     1      544         545      -1    1004    1022     -18 B6     
## 5  2013     1     1      554         600      -6     812     837     -25 DL     
## 6  2013     1     1      554         558      -4     740     728      12 UA     
## # … with 9 more variables: flight <int>, tailnum <chr>, origin <chr>,
## #   dest <chr>, air_time <dbl>, distance <dbl>, hour <dbl>, minute <dbl>,
## #   time_hour <dttm>, and abbreviated variable names ¹​sched_dep_time,
## #   ²​dep_delay, ³​arr_time, ⁴​sched_arr_time, ⁵​arr_delay

flights è un dataframe ma di tipo tibble che è un formato particolarmente utile nel tidyverse. In particolare mi fa vedere solo poche righe e tutte le colonne. Inoltre ogni colonna ha, sotto il nome, il tipo di dato.

In questo capitolo andremo a vedere le 5 più utili funzioni del tidyverse che servono per risolvere la maggior parte dei problemi di manipolazione dati:

1. filter(): prende le osservazioni in base al loro valore
2. arrange(): riordina i dati
3. select(): prende le variabili in base al nome
4. mutate(): crea nuove variabili
5. summarize(): collassa più valori in uno unico da usare insieme a group_by

4.1 filter()

Se voglio salvare il risultato in una variabile e contemporaneamente visualizzarne l’output posso racchiudere tutta l’istruzione all’iterno di una parentesi.

(feb3 <- flights %>% filter(
  month==1, day==1
))

## # A tibble: 842 × 19
##     year month   day dep_time sched_de…¹ dep_d…² arr_t…³ sched…⁴ arr_d…⁵ carrier
##    <int> <int> <int>    <int>      <int>   <dbl>   <int>   <int>   <dbl> <chr>  
##  1  2013     1     1      517        515       2     830     819      11 UA     
##  2  2013     1     1      533        529       4     850     830      20 UA     
##  3  2013     1     1      542        540       2     923     850      33 AA     
##  4  2013     1     1      544        545      -1    1004    1022     -18 B6     
##  5  2013     1     1      554        600      -6     812     837     -25 DL     
##  6  2013     1     1      554        558      -4     740     728      12 UA     
##  7  2013     1     1      555        600      -5     913     854      19 B6     
##  8  2013     1     1      557        600      -3     709     723     -14 EV     
##  9  2013     1     1      557        600      -3     838     846      -8 B6     
## 10  2013     1     1      558        600      -2     753     745       8 AA     
## # … with 832 more rows, 9 more variables: flight <int>, tailnum <chr>,
## #   origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>, hour <dbl>,
## #   minute <dbl>, time_hour <dttm>, and abbreviated variable names
## #   ¹​sched_dep_time, ²​dep_delay, ³​arr_time, ⁴​sched_arr_time, ⁵​arr_delay

Gli operatori logici sono estremamente utili:

Una cosa particolarmente imporante sono gli NA. Questi vengono esclusi di default e se li voglio, devo includerli esplicitamente come nel seguente caso:

df <- tibble(x = c(1, NA, 3))
filter(df, x > 1)

## # A tibble: 1 × 1
##       x
##   <dbl>
## 1     3

#> # A tibble: 1 x 1
#>       x
#>   <dbl>
#> 1     3
filter(df, is.na(x) | x > 1)

## # A tibble: 2 × 1
##       x
##   <dbl>
## 1    NA
## 2     3

#> # A tibble: 2 x 1
#>       x
#>   <dbl>
#> 1    NA
#> 2     3

4.2 arrange()

Arrange, non seleziona le righe ma, piuttosto le riordina. Per ordinarle dal maggiore al minore si utilizza desc(). I valori NA sono sempre messi al fondo.

arrange(flights,  desc(dep_time), year, day)

## # A tibble: 336,776 × 19
##     year month   day dep_time sched_de…¹ dep_d…² arr_t…³ sched…⁴ arr_d…⁵ carrier
##    <int> <int> <int>    <int>      <int>   <dbl>   <int>   <int>   <dbl> <chr>  
##  1  2013     4     2     2400       2359       1     339     343      -4 B6     
##  2  2013     9     2     2400       2359       1     411     340      31 B6     
##  3  2013     4     4     2400       2355       5     347     345       2 B6     
##  4  2013    12     5     2400       2359       1     427     440     -13 B6     
##  5  2013     2     7     2400       2359       1     432     436      -4 B6     
##  6  2013     2     7     2400       2359       1     443     444      -1 B6     
##  7  2013     7     7     2400       1950     250     107    2130     217 AA     
##  8  2013    12     9     2400       2359       1     432     440      -8 B6     
##  9  2013    12     9     2400       2250      70      59    2356      63 B6     
## 10  2013     8    10     2400       2245      75     110       1      69 B6     
## # … with 336,766 more rows, 9 more variables: flight <int>, tailnum <chr>,
## #   origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>, hour <dbl>,
## #   minute <dbl>, time_hour <dttm>, and abbreviated variable names
## #   ¹​sched_dep_time, ²​dep_delay, ³​arr_time, ⁴​sched_arr_time, ⁵​arr_delay

4.3 select()

Come in SQL serve per selezionare le colonne.

# select(flights, year, month, day) #sel 4 colonne
# select(year:day) #seleziona da year a day (month incluso)
# select(-(year:day)) #seleziona tutto a parte da year a day

# Posso usare anche diverse funzioni quali:
# starts_with("abc"): matches names that begin with “abc”.
# ends_with("xyz"): matches names that end with “xyz”.
# contains("ijk"): matches names that contain “ijk”.
# matches("(.)\\1")

ATTENZIONE Select può anche essere utilizzato per rinominare le variabili. Tuttavia ha l’inconveniente di droppare le variabili non esplicitamente nominate. Al suo posto si può usare rename()

rename(flights, tail_num = tailnum)

## # A tibble: 336,776 × 19
##     year month   day dep_time sched_de…¹ dep_d…² arr_t…³ sched…⁴ arr_d…⁵ carrier
##    <int> <int> <int>    <int>      <int>   <dbl>   <int>   <int>   <dbl> <chr>  
##  1  2013     1     1      517        515       2     830     819      11 UA     
##  2  2013     1     1      533        529       4     850     830      20 UA     
##  3  2013     1     1      542        540       2     923     850      33 AA     
##  4  2013     1     1      544        545      -1    1004    1022     -18 B6     
##  5  2013     1     1      554        600      -6     812     837     -25 DL     
##  6  2013     1     1      554        558      -4     740     728      12 UA     
##  7  2013     1     1      555        600      -5     913     854      19 B6     
##  8  2013     1     1      557        600      -3     709     723     -14 EV     
##  9  2013     1     1      557        600      -3     838     846      -8 B6     
## 10  2013     1     1      558        600      -2     753     745       8 AA     
## # … with 336,766 more rows, 9 more variables: flight <int>, tail_num <chr>,
## #   origin <chr>, dest <chr>, air_time <dbl>, distance <dbl>, hour <dbl>,
## #   minute <dbl>, time_hour <dttm>, and abbreviated variable names
## #   ¹​sched_dep_time, ²​dep_delay, ³​arr_time, ⁴​sched_arr_time, ⁵​arr_delay

Molto utile anche la funzione everything:

select(flights, time_hour, air_time, everything())

## # A tibble: 336,776 × 19
##    time_hour           air_t…¹  year month   day dep_t…² sched…³ dep_d…⁴ arr_t…⁵
##    <dttm>                <dbl> <int> <int> <int>   <int>   <int>   <dbl>   <int>
##  1 2013-01-01 05:00:00     227  2013     1     1     517     515       2     830
##  2 2013-01-01 05:00:00     227  2013     1     1     533     529       4     850
##  3 2013-01-01 05:00:00     160  2013     1     1     542     540       2     923
##  4 2013-01-01 05:00:00     183  2013     1     1     544     545      -1    1004
##  5 2013-01-01 06:00:00     116  2013     1     1     554     600      -6     812
##  6 2013-01-01 05:00:00     150  2013     1     1     554     558      -4     740
##  7 2013-01-01 06:00:00     158  2013     1     1     555     600      -5     913
##  8 2013-01-01 06:00:00      53  2013     1     1     557     600      -3     709
##  9 2013-01-01 06:00:00     140  2013     1     1     557     600      -3     838
## 10 2013-01-01 06:00:00     138  2013     1     1     558     600      -2     753
## # … with 336,766 more rows, 10 more variables: sched_arr_time <int>,
## #   arr_delay <dbl>, carrier <chr>, flight <int>, tailnum <chr>, origin <chr>,
## #   dest <chr>, distance <dbl>, hour <dbl>, minute <dbl>, and abbreviated
## #   variable names ¹​air_time, ²​dep_time, ³​sched_dep_time, ⁴​dep_delay, ⁵​arr_time

4.4 mutate()

Serve per aggiungere nuove variabili alla fine del dataset. Iniziamo con un dataset e l’aggiunta di una nuova colonna:

flights_sml <- select(flights, 
  year:day, 
  ends_with("delay"), 
  distance, 
  air_time,
  hour
) 

#puoi riferirti a var appena create
mutate(flights_sml,
  gain = dep_delay - arr_delay,
  speed = distance / air_time * 60,
  gain_per_hour = gain / hour 
)

## # A tibble: 336,776 × 11
##     year month   day dep_delay arr_d…¹ dista…² air_t…³  hour  gain speed gain_…⁴
##    <int> <int> <int>     <dbl>   <dbl>   <dbl>   <dbl> <dbl> <dbl> <dbl>   <dbl>
##  1  2013     1     1         2      11    1400     227     5    -9  370.  -1.8  
##  2  2013     1     1         4      20    1416     227     5   -16  374.  -3.2  
##  3  2013     1     1         2      33    1089     160     5   -31  408.  -6.2  
##  4  2013     1     1        -1     -18    1576     183     5    17  517.   3.4  
##  5  2013     1     1        -6     -25     762     116     6    19  394.   3.17 
##  6  2013     1     1        -4      12     719     150     5   -16  288.  -3.2  
##  7  2013     1     1        -5      19    1065     158     6   -24  404.  -4    
##  8  2013     1     1        -3     -14     229      53     6    11  259.   1.83 
##  9  2013     1     1        -3      -8     944     140     6     5  405.   0.833
## 10  2013     1     1        -2       8     733     138     6   -10  319.  -1.67 
## # … with 336,766 more rows, and abbreviated variable names ¹​arr_delay,
## #   ²​distance, ³​air_time, ⁴​gain_per_hour

Se al contrario vuoi cestinare tutti gli altri campi e modificarne solo alcuni puoi usare transmutate():

flights_trmt <- transmute(flights,
  gain = dep_delay - arr_delay,
  hours = air_time / 60,
  gain_per_hour = gain / hours
)

Ci sono molteplici funzioni che possono essere utili in questa funzione e ne vediamo alcune.

• Operatori aritmetici + - * / ^: sono vettorizzati con la regola che se finisce un vettore, ricomincia il calcolo
• Aritmetica modulare %/% (la parte intera di una divisione) %% (resto di una divisione – non ho capito bene!!)
• Logaritmici: log() log2() log10()
• Offset: lead() lag() permettono di fare riferimento alla parte iniziale o terminale. Utile per trovare le differenze incrementali (x-lag(x)) oppure quando cambia il valore (x!=lag(x))
• Cumulativi: cumsum cumprod cummax cummean
• Ranking: min_rank() row_number(), dense_rank(), percent_rank(), cume_dist(), ntile()

Esempio offset:

vettore <- c(1,12,33,44,11,67)
lag(vettore)

## [1] NA  1 12 33 44 11

lead(vettore)

## [1] 12 33 44 11 67 NA

4.5 summarise()

Serve per collassare un data.frame in un’unica riga:

summarise(flights, media = mean(dep_delay, na.rm=TRUE))

## # A tibble: 1 × 1
##   media
##   <dbl>
## 1  12.6

flights %>% 
  group_by(origin) %>% 
  summarise( media = mean(dep_delay, na.rm=TRUE))

## # A tibble: 3 × 2
##   origin media
##   <chr>  <dbl>
## 1 EWR     15.1
## 2 JFK     12.1
## 3 LGA     10.3

Tutti questi comandi possono essere messi in pipe.

flights %>% 
  group_by(dest) %>% 
  summarise(
    count = n(),
    dist = mean(distance, na.rm = TRUE),
    delay = mean(arr_delay, na.rm = TRUE)
  ) %>% 
  filter(count > 20, dest != "HNL") %>% 
  ggplot(mapping = aes(x = dist, y = delay)) +
  geom_point(aes(size = count), alpha = 1/3) +
  geom_smooth(se = FALSE)

## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

oppure così

# Molto interessante mettere in ascisse il count 
# ed in ordinate il valore tipo il numero di rate
flights %>% 
  filter(!is.na(dep_delay), !is.na(arr_delay)) %>%  
  group_by(tailnum) %>% 
  summarise(n=n(), delay = mean(arr_delay)) %>% 
  ggplot(mapping = aes(x = n, y = delay)) + 
  geom_point(alpha = 1/10)

Vediamo un esempio con il dataset Lahman:

batting <- tibble(Lahman::Batting)
batters <- batting %>% 
  group_by(playerID) %>% 
  summarise(
    ba = sum(H, na.rm = TRUE) / sum(AB, na.rm = TRUE),
    ab = sum(AB, na.rm = TRUE)
  ) %>% filter(ab > 100) 
   
batters %>% 
  ggplot(mapping = aes(x = ab, y = ba)) +
    geom_point() + 
    geom_smooth(se = FALSE)

## `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'

#> `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'

Ci sono diverse istruzioni molto interessanti che viaggiano insieme al summarise:

• Misure di posizione: mean(), median(), sum(), first(x), nth(x, 2), last(x), range() (range ti dà il minimo e massimo valore del vettore)
• Misure di spread: sd(), IQR(), mad() (IQR è interquartile)
• Misure di rank: min(), quantile(x,0.25), max()
• Misure di conteggio n(), n_distinct()
• Misure di conteggio logico: interessante! Quando utilizzo delle condizioni delle misure, i TRUE vengono convertiti in 1 mentre i FALSE in 0. Quindi posso fare il conteggio di sum(x<10) oppure mean(y==0).

flights %>% 
  group_by(year, month) %>% 
  summarise(numero =n(), n_early = sum(dep_time < 1000))

## `summarise()` has grouped output by 'year'. You can override using the
## `.groups` argument.

## # A tibble: 12 × 4
## # Groups:   year [1]
##     year month numero n_early
##    <int> <int>  <int>   <int>
##  1  2013     1  27004      NA
##  2  2013     2  24951      NA
##  3  2013     3  28834      NA
##  4  2013     4  28330      NA
##  5  2013     5  28796      NA
##  6  2013     6  28243      NA
##  7  2013     7  29425      NA
##  8  2013     8  29327      NA
##  9  2013     9  27574      NA
## 10  2013    10  28889      NA
## 11  2013    11  27268      NA
## 12  2013    12  28135      NA

4.6 ungroup()

Utile se necessiti di togliere il raggruppamento e ritornare alla situazione precedente.

5 Esplorare la Analisi dei dati

In questo capitolo vediamo come utilizzare la Data Visualisation e la Data transformation per esplorare i dati in maniera sistematica, un’attività che gli statistici chiamano Exploring Data Analysis o EDA.

1. Genera domande sui tuoi dati,
2. Cerca risposte con trasformazioni, visualizzazioni, modelli
3. Usa ciò che ottieni per raffinare le domande.

Non si tratta di un processo formale ma casuale e sempre più raffinato. Quando ti fai una domanda, vai ad esplorare un determinato ambito del dataset.

Le domande che devono sempre guidare sono:

1. che tipo di variazione è occorsa alle mie variabili?
2. che tipo di correlazione è occorsa alle mie variabili?

E’ ineteressante fare un taglio della base dati per cluster tipo. Altra cosa geom_freqpoly è come geom_histogram ma al posto degli istogrammi fa apparire delle rette:

diamonds  %>% 
  group_by(cut_width(carat, 0.5)) %>% 
  summarise(numero = n())

## # A tibble: 11 × 2
##    `cut_width(carat, 0.5)` numero
##    <fct>                    <int>
##  1 [-0.25,0.25]               785
##  2 (0.25,0.75]              29498
##  3 (0.75,1.25]              15977
##  4 (1.25,1.75]               5313
##  5 (1.75,2.25]               2002
##  6 (2.25,2.75]                322
##  7 (2.75,3.25]                 32
##  8 (3.25,3.75]                  5
##  9 (3.75,4.25]                  4
## 10 (4.25,4.75]                  1
## 11 (4.75,5.25]                  1

smaller <- diamonds %>% 
  filter(carat < 3)
  
ggplot(data = smaller, mapping = aes(x = carat)) +
  geom_histogram(binwidth = 0.1)

ggplot(data = smaller, mapping = aes(x = carat, colour = cut)) +
  geom_freqpoly(binwidth = 0.1)

5.1 Missing values

Nel caso si ottengano risultati inaspettati è possibile sostituire le righe con valori strani con degli NA e un ifelse.

diamonds2 <- diamonds %>% 
  mutate(y = ifelse(y<3 | y>20, NA, y))

Esiste anche il case_when() particolarmente utile all’interno di mutate().

Se si creano dei grafici con degli NA, ggplot non li visualizza nel grafico ma fornisce un warning che può essere eliminato con la clausola na.rm = TRUE.

5.2 Covarianza

Se la varianza indica la variazione all’interno di una variabile, la covarianza indica le variazioni tra variabili ovvero la tendenza di 2 o più variabili di cambiare in una direzione.

5.2.1 Una variabile categorica ed una continua

E’ importante capire che è molto differente una variabile categorica e una continua. E’ molto comune voler analizzare una variabile continua differenziata da una categorica: in tal caso geom_freqpoly non è così utile perchè sulle ascisse c’è il conteggio. Per questo motivo se un gruppo è molto più piccolo di un altro è difficile vederne la sagoma. Per esempio:

ggplot(data = diamonds, mapping = aes(x = price)) + 
  geom_freqpoly(mapping = aes(colour = cut), binwidth = 500)

è difficile vedere le differenze perchè Fair ha un numero piccolo di osservazioni:

ggplot(diamonds) + 
  geom_bar(mapping = aes(x = cut))

In questi casi possiamo cambiare l’asse y ed invece del count, possiamo mettere in questo caso, la densità che è il conteggio a somma 100. IMPORTANTE:

ggplot(
  data = diamonds, 
  mapping = aes(x = price, y = ..density..)) + 
  geom_freqpoly(mapping = aes(colour = cut), binwidth = 500)

Una visualizzazione alternativa sono i boxplot:

• Ogni box parte dal 25esimo al 75esimo percentile ovvero la distanza conosciuta come intervallo interquartile o IQR e nel mezzo c’è una linea che determina la mediana dando la sensazione della distribuzione osservano la mediana se è in mezzo o più spostata da un lato.
• Punti visuali che dicono se ci sono distribuzioni ad oltre 1.5 volte dal interquartile
• Una linea che parte da ogni lato del IQR

Spiegazione grafica dei tre punti sopra

Co questo grafico si veda la situazione:

ggplot(data = diamonds, mapping = aes(x = cut, y = price)) +
  geom_boxplot()

Questo grafico dà meno informazioni sulla distribuzione ma è compatto e più facile da comparare avvallando la testi che il prezzo medio dei diamanti più pregiati è mediamente più basso di quelli meno pregiati.

ATTENZIONE: il cut è un fattore ma in altri casi potrei voler ordinare la mia variabile con `reorder``.

Per esempio questo grafico:

ggplot(data = mpg, mapping = aes(x = class, y = hwy)) +
  geom_boxplot()

posso voler riordinare la classe per la mediana del hwy che sarebbero le ordinate:

ggplot(data = mpg) +
  geom_boxplot(mapping = aes(x = reorder(class, hwy, FUN = median), y = hwy))

Nel caso di nomi di variabili lunghi posso invertire gli assi con coord_flip.

5.2.2 Due variabili categoriche

Per visualizzare la correlazione tra due categoriche puoi contare il numero di osservazioni per ogni combinazione. Per esempio puoi usare geom_count.

ggplot(data = diamonds) +
  geom_count(mapping = aes(x = cut, y = color))

Un altro approccio può essere una semplice tabella:

diamonds %>% 
  count(color, cut)

## # A tibble: 35 × 3
##    color cut           n
##    <ord> <ord>     <int>
##  1 D     Fair        163
##  2 D     Good        662
##  3 D     Very Good  1513
##  4 D     Premium    1603
##  5 D     Ideal      2834
##  6 E     Fair        224
##  7 E     Good        933
##  8 E     Very Good  2400
##  9 E     Premium    2337
## 10 E     Ideal      3903
## # … with 25 more rows

Oppure un geom_tile

diamonds %>% 
  count(color, cut) %>%  
  ggplot(mapping = aes(x = color, y = cut)) +
    geom_tile(mapping = aes(fill = n))

5.2.3 Due variabili continue

Uno dei migliori grafici è lo scatterplot. Il problema è che quando il numero di punti aumenta tanto, inizia a diventare confusionario: in tal caso posso usare la trasparenza (con alpha=1/100) oppure bin.

geom_bin2() e geom_hex dividono il piano in due contenitori e poi li usano un colore per visualizzare quanti punti appaiono in ciascun contenitore: il primo con i rettangoli, il secondo con gli esagoni.

ggplot(data = smaller) +
  geom_bin2d(mapping = aes(x = carat, y = price))

library("hexbin")

## Warning: il pacchetto 'hexbin' è stato creato con R versione 4.2.1

ggplot(data = smaller) +
  geom_hex(mapping = aes(x = carat, y = price))

Un’altra tecnica è far diventare una variabile continua come se fosse una variabile categorica e poi sfruttare una delle visualizzazioni viste sopra:

ggplot(data = smaller, mapping = aes(x = carat, y = price)) + 
  geom_boxplot(mapping = aes(group = cut_width(carat, 0.1)))

IMPORTANTE: cut_width(x, width) divide x in contenitori di larghezza width.

Un altro approccio è visualizzare lo stesso numero di osservazioni per contenitore. Anzichè usare cut_width si userà cut_number:

ggplot(data = smaller, mapping = aes(x = carat, y = price)) + 
  geom_boxplot(mapping = aes(group = cut_number(carat, 20)))

5.3 Patterns e Modelli

Se esiste un pattern ovvero una relazione fra variabili bisogna chiedersi:

• Si tratta di una coincidenza?
• Come posso descrivere la relazione del pattern?
• Quanto è forte la relazione del pattern?
• La relazione cambia se guardo un sottogruppo del pattern?

Se si pensa la variazione di una variabile come fonte di incertezza, la covarianza è un fenomeno che riduce tale incertezza. Posso fare delle previsioni sulla seconda variabile al cambiare della prima. I modelli sono strumenti per estrarre pattern dai dati. Per esempio: nel dataset diamonds è difficile trovare una relazione tra cut, price perchè cut, carat e price, carat sono molto forti. Si può usare un modello per rimuovere la forte relazione price, carat e quindi esplorare ciò che rimane.

Posso quindi impostare un modello che preveda il prezzo dai carati e calcoli i residui (differenza tra il valore predetto ed attuale). I residui mi danno una vista del prezzo del diamante una volta che l’effetto dei carati è stato rimosso:

library(modelr)

## Warning: il pacchetto 'modelr' è stato creato con R versione 4.2.1

mod <- lm(log(price) ~ log(carat), data = diamonds)

diamonds2 <- diamonds %>% 
  add_residuals(mod) %>% 
  mutate(resid = exp(resid))

ggplot(data = diamonds2) + 
  geom_point(mapping = aes(x = carat, y = resid))

Una volta rimossa la forte relazione tra prezzo e carati, si può vedere la relazione tra taglio e prezzo: in relazione alla dimensione, migliore qualità è più caro.

ggplot(data = diamonds2) + 
  geom_boxplot(mapping = aes(x = cut, y = resid))

6 Wrangle

Data wrangling sarebbe l’arte di rendere i dati facilmente utilizzabili per la visualizzazione e la modellazione. E’ una parte estremamente importante senza la quale non si riuscirebbe a lavorare con i propri dati e si divide in tre parti:

• tibbles In italiano significa bocconcino. Vediamo in cosa sono diversi dai dataframe e come costruirli In autonomia.
• data import Vedremo come importare i dati dal disco e ci focalizziamo sul formato plaintext Ciò non di meno andremo a vedere anche altri tipi di formati
• tidy data Vedremo un modo consistente per immagazzinare i dati. Questo ci permetterà in futuro di visualizzarli, manipolarli, modellarli.

Ci orienteremo in particolare in alcuni tipologie di dati: Dati relazionali, stringhe, date, fattori.