Introduzione ad R (notebook 1)

Statistica II - 750AA

Installare R ed Rstudio

Seguite le istruzioni partendo dal sito di posit: https://posit.co/download/rstudio-desktop/. La pagina web ufficiale del progetto R è invece https://www.r-project.org/ e in particolare CRAN https://cran.r-project.org/.

Al primo avvio, create un progetto in una cartella dedicata ad esempio 750AA (è sempre buona pratica creare progetti per organizzare il proprio flusso di lavoro) e caricate il file .Rmd che potete scaricare da questa pagina cliccando in alto a destra.

Su Rstudio potete eseguire righe di codice cliccando Ctrl-Enter o interi blocchi di codice con Ctrl-Shift-Enter.

Presentiamo ora i comandi e le caratteristiche di base di R – per brevità la presentazione segue quella di https://learnxinyminutes.com/r. Un’ottima referenza caldamente suggerita per usare R per l’analisi dei dati con maggiore consapevolezza è R for Data Science (2e) disponibile alla pagina https://r4ds.hadley.nz/. Una guida per le funzioni più avanzate di R è Advanced R, https://adv-r.hadley.nz/index.html.

R è a tutti gli effetti un linguaggio di programmazione, ma si può usare tranquillamente come una calcolatrice con molte funzioni per la statistica: è facile potenziare R installando pacchetti aggiuntivi con il comando install.packages().

# I commenti si indicano con il simbolo hash #
# installiamo il pacchetto cluster che contiene funzioni per il clustering di dati (lo useremo più avanti nel corso)

install.packages("cluster")

Una volta installato il pacchetto, questo può essere caricato con il comando library().

library("cluster")

Se avete dubbi su un qualsiasi comando, potete cercare nella documentazione ufficiale (tab Help in basso a destra su Rstudio) oppure digitando ?(nome comando) nella console (basso a sinistra).

?hist
# il comando genera istogrammi a partire da un vettore di osservazioni

R contiene già molti dataset (termine tecnico in R è data frame) standard pre-caricati, già a disposizione, che useremo negli esempi (digitare il comando data() per avere un elenco). Anche le informazioni sui data frame principali sono accessibili tramite l’help.

?mtcars

Con il comando head() visualizziamo solo le prime righe di un data frame.

head(mtcars)

Attenzione! per i progetti di esame è richiesto di lavorare su dati nuovi (più avanti vediamo come caricarli). Ma iniziamo dai tipi di dati più semplici.

Con il comando summary() otteniamo informazioni circa la distribuzione delle singole colonne del data frame.

summary(mtcars)

Con il comando plot() otteniamo uno scatterplot tra coppie di colonne (più avanti vediamo meglio).

plot(mtcars$mpg, mtcars$cyl)

Con il comando hist() otteniamo invece un istogramma.

hist(mtcars$mpg)

Un grafico a barre (utile se la variabile osservata è discreta o qualitativa) è invece ottenibile con il comando barplot().

barplot(table(factor(mtcars$cyl)))

Con il comando boxplot() otteniamo un diagramma a scatola e baffi (box and whiskers).

boxplot(mtcars$hp)

Classi di oggetti

In questa sezione presentiamo i tipi di oggetti importanti di R: interi, numeri, caratteri, logici e fattori. Ce ne sono altri, ma questi sono i minimi indispensabili per iniziare. R è abbastanza flessibile sulle classi, ma se avete dubbi potete usare la funzione class() per avere informazioni precise.

?class

Integer

Gli interi (non decimali) sono specificati con L.

5L         
class(5L)

In R, ogni singolo valore, come 5L, è considerato un vettore di lunghezza 1.

length(5L)

Per avere vettori di lunghezza maggiore si può usare la funzione c() (concatena).

c(4L, 5L, 8L, 3L)          
length(c(4L, 5L, 8L, 3L))  
class(c(4L, 5L, 8L, 3L))

Numeric

Un numeric è un numero decimale (a precisione doppia)

5           
class(5)    

Di nuovo, tutto in R è un vettore, quindi possiamo creare un vettore numerico con più di un elemento

c(3, 3, 3, 2, 2, 1)

Possiamo usare anche la notazione scientifica (AeB significa \(A \cdot 10^B\))

5e4         
6.02e23     # numero di Avogadro
1.6e-35     # lunghezza di Planck

Possiamo avere anche numeri infinitamente grandi con Inf.

class(Inf)  # "numeric"
class(-Inf) # "numeric"

Un esempio: per ottenere la CDF dalla densità gaussiana (nel punto \(-1\)): \(\Phi(-1) = \int_{-\infty}^{-1}\exp(-x^2/2) dx\)

integrate(dnorm, -Inf, -1)

Confronta comunque con il comando pnorm():

pnorm(-1)

Fare operazioni su un mix di interi e numeric restituisce un altro numeric:

10L + 66L   
53.2 - 4   
2.0 * 2L   
3L / 4     
3 %% 2

Operazioni non valide restituiscono NaN, ossia Not-A-Number.

0 / 0       # NaN
class(NaN)  # "numeric"

Da non confondere con NA, ossia Not-Available per i dati mancanti.

NA
class(NA)

Possiamo fare operazioni su due vettori con lunghezza maggiore di 1, # purché la lunghezza del vettore più grande sia un multiplo intero di quello più piccolo

c(1, 2, 3) + c(1, 2, 3)     # 2 4 6

Poiché un singolo numero è un vettore di lunghezza uno, gli scalari sono applicati elemento per elemento ai vettori.

(4 * c(1, 2, 3) - 2) / 2    # 1 3 5

Ad eccezione degli scalari, attenzione quando eseguite operazioni su vettori con lunghezze diverse. Anche se si può fare, allineare le lunghezze è una pratica migliore e più facile da leggere nella maggior parte dei casi.

c(1, 2, 3, 1, 2, 3) * c(1, 2)               # 1 4 3 2 2 6
c(1, 2, 3, 1, 2, 3) * rep(c(1, 2), 3)   # 1 4 3 2 2 6

Tante funzioni sono già implementate in R di base:

exp(1) # esponenziale a base naturale
log(3) # logaritmo base naturale
log(3, base=10) # logaritmo in base 10
sin(3) # seno (in radianti)
atan(1.3) # arcotangente
sqrt(2) # radice quadrata
2**(1/2) ## esponenziale

Character

Non c’è differenza tra caratteri e stringhe (sequenze di più caratteri) in R:

"Statistica"          
class("Statistica")    # "character"
class("S")          # "character"

Questi sono tutti vettori di caratteri di lunghezza 1. Al solito uno più lungo si può ottenere concatenando.

c("Alice", "Bob", "Carlo", "Davide", "Elisabetta")

length(c("Alice","Bob","Carlo")) # 3

R ha diversi vettori di caratteri incorporati:

letters
LETTERS
month.name

Logic

In R, un oggetto logico è un booleano:

class(TRUE)     # "logical"
class(FALSE)    # "logical"
class(NA)

Espressioni possono essere confrontate con le operazioni booleane == (uguale), !=, (diverso).

TRUE == TRUE    # TRUE
5 == 3
FALSE != FALSE  # FALSE
5!= 3

Si usano |, &, ! per le operazioni logiche di disgiunzione (o) congiunzione (e) e negazione.

TRUE | FALSE    # TRUE
TRUE & FALSE    # FALSE
! (5==3)

Al solito le operazioni si applicano a vettori elemento per elemento:

c(1,2,3,4)==c(1,5,3,2)

c(TRUE, FALSE) & TRUE

Factor

La classe factor (fattore) è per i dati categorici (variabili statistiche discrete). I fattori possono essere ordinati, come ad esempio i giudizi scolastici (ottimo, buono, ecc.) o non ordinati, ad esempio come i colori. La funzione factor trasforma un vettore in un corrispondente vettore di tipo factor. Ad esempio

factor(c("blue", "blue", "green", NA, "blue", "red", "yellow"))

I livelli sono i valori che i dati categoriali possono assumere. Notiamo che i dati mancanti non entrano nei livelli. Per accedere ai livelli usiamo la funzione levels(). La funzione table() permette di ottenere invece una tabella di contingenza contenente le frequenze dei vari livelli.

levels(factor(c("blue", "blue", "green", NA, "blue", "red", "yellow")))
table(factor(c("blue", "blue", "green", NA, "blue", "red", "yellow")))

NULL

L`oggetto NULL è un caso strano, si può usare per svuotare un vettore.

class(NULL) # NULL
pappagallo <- c("becco", "piume", "ali", "occhi")
pappagallo

pappagallo <- NULL

pappagallo

Variabili, cicli e funzioni

Una variabile è come una scatola in cui tenere un oggetto (valore) per uso successivo. Questa operazione è detta di assegnazione del valore alla variabile. Una volta che abbiamo variabili, possiamo scrivere cicli (for, while, ecc.) e funzioni.

Variabili

Ci sono molti modi per assegnare un valore a una variabile.

x = 5       # questo è un modo accettato
y <- "1"    # questo è il modo standard suggerito
TRUE -> z   # questo pure funziona

Non ci sono convenzioni ufficiali per i nomi delle variabili. Attenzione! i nomi sono case sensitive (maiuscola/minuscola fa differenza) È vietato fare iniziare il nome con un numero, ed usare nomi di oggetti fondamentali (TRUE, FALSE, NULL, ecc.). Per il resto, è buona pratica usare nomi che siano esplicativi, ma anche non troppo lunghi. Quando introducete una variabile possibilmente aggiungete un commento circa il suo ruolo.

DurataLezione <- 3 #Pascal case
durataLezione <- 3 #lower camelCase
durata_lezione <- 3 #lower_case_with_underscores, aka snake_case

Cicli (loops)

Ci sono i cicli for.

for (i in 1:4) {
    print(i)
}

Il comando A:B definisce il vettore di interi da A a B, estremi inclusi.

(-3):5

10:1

È possibile fare un ciclo su un qualsiasi vettore.

for( animale in c("cane", "gatto", "topo")){ print(animale)}

Ci sono i cicli while, ossia che si ripetono finché la condizione tra parentesi è realizzata (TRUE)

a <- 10
while (a > 4) {
    print(a)
    a <- a - 1
}

In realtà usare i cicli for e while in R è abbastanza lento, e quando le funzioni vengono applicate su interi vettori conviene sfruttare il fatto che R lo fa in automatico.

numeri <- 1:10

# una sola riga per calcolare i quadrati dei numeri da 1 a 10
quadratiVeloci <- numeri**2

# usiamo un ciclo for
quadratiLenti <- NULL
for (i in numeri){
  quadratiLenti <-c( quadratiLenti, i**2)
}

# il risultato è lo stesso
quadratiLenti == quadratiVeloci

If/else

Si può introdurre operazioni condizionali (if, else)

x <- 4
y <- 5

if (x > y) {
    print("x è maggiore di y")
} else if(x==y) {
    print("x è uguale ad y")
} else {
  print("x è minore di y")
}

Funzioni

È possibile definire nuove funzioni di variabili in questo modo:

sommaNuova <- function(x,y) {
    z = x  + y 
    return(z)
}

A questo punto è possibile utilizzarla come qualsiasi altra funzione di R.

sommaNuova(3, 4)

È possibile definire funzioni con alcuni argomenti opzionali (dopo quelli obbligatori).

sommaNuova <- function(x, y=0){
  return(x+y)
}

sommaNuova(3)

Strutture di dati

Abbiamo già visto che R non distingue tra scalari e vettori (purché siano dati dello stesso tipo). Vediamo ora ulteriori strutture come matrici, data frames e liste.

Vettori

Abbiamo già visto i vettori e la funzioni di base c().

vettore <- c(8, 9, 10, 11)
vettore

Ricordate anche il comando : per sequenze di interi. Per generare sequenze più complicate, anche non di interi, si può usare il comando seq().

8:11
seq(8, 11, by=1)
seq(8, 11, by=0.5)

Il comando length() ritorna la dimensione del vettore.

length(vettore)

È spesso utile anche il comando rep per replicare un vettore.

rep(c(1,3), 2)

# crea un vettore di zeri lungo 4

rep(0, 4)

# crea un vettore di TRUE lungo 5

rep(TRUE, 5)

Possiamo chiamare uno o più elementi di un vettore nella posizione \(k\) indicandoli tra parentesi quadrate (usando un altro vettore contenente le posizioni da chiamare). Attenzione! R comincia a contare da \(1\) (altri linguaggi da \(0\)).

vettore[1]
vettore[2]
vettore[c(1,4)]
vettore[1:3]

# ma invece
vettore[5]
vettore[0]

Usando indici negativi invece escludiamo quelle componenti:

vettore[-1]
vettore[-(1:3)]

Inserendo un vettore di valori booleani possiamo filtrare solo le componenti corrispondenti ai valori TRUE.

maschera <- as.logical(vettore %% 2)

maschera

vettore[maschera]

Viceversa, data una condizione che può essere soddisfatta o meno, possiamo estrarre gli indici per cui vale con il comando which.

which(vettore <9 | vettore > 10)

In questo caso mostra gli indici del vettore che corrispondono a valori minori di \(9\) o maggiori di \(10\). Può essere utile il comando which.max() che trova la posizione del valore massimo e simimente which.min() per il minimo.

vettore2 =c(1:10, 15:0)

# il comando max() trova il valore massimo

max(vettore2)

# il comando which.max() trova l'indice corrispondente (argmax, o punto di massimo)
which.max(vettore2)

I comandi `head() e tail() permettono di ottenere le prime o le ultime componenti di un vettore.

head(1:1000)
tail(1:1000, n=10)

Indicatori statistici

Tante funzioni di statistica, in particolare descrittiva, sono già implementate di base.

mean(vettore2)
var(vettore2)
sd(vettore2)
median(vettore2)
summary(vettore2)
quantile(vettore2, 1/4)
quantile(vettore2, 3/4)
quantile(vettore2, .95)

Per trovare la moda ad esempio di un vettore di fattori (variabile statistica discreta) basta usare il comando which.max() in combinazione con table().

colori <- factor(c("black", "red", "blue", "blue", "green", NA, "blue", "red", "yellow"))
table(colori)

moda <- which.max(table(colori))
moda

Notate che il \(2\) non è la frequenza (che sarebbe \(3\)) bensì la posizione della moda nel vettore dei livelli.

Matrici

Possiamo creare una matrice con entrate tutte dello stesso tipo (spesso numeric) cambiando la forma di un vettore, con il comando matrix().

mat <- matrix(c(1, 2, 3, 4, 5, 6), nrow = 3, ncol = 2,)
mat

Diversamente da un vettore, la classe di una matrice è sempre matrix (anche se non contiene numeri).

class(mat)

Per accedere alla componente di riga \(i\) e colonna \(j\) (ricordare che si conta da \(1\)) si scrive mat[i,j]. Per chiedere l’intero vettore della riga \(i\) invece mat[i,], mentre per la colonna \(j\) si scrive mat[,j].

# prima riga
mat[1, ]       
# tutte le colonne tranne la prima
 mat[, -1]  
#  prima e terza riga, seconda colonna
mat[c(1,3), 2]      

Operazioni di base su matrici sono già disponibili senza pacchetti aggiuntivi.

 
# matrice trasposta
t(mat)

# prodotto di matrici (righe per colonne)

mat %*% t(mat)

# attenzione a non confondere usando solo * (prodotto componente per componente)

mat * mat

I comandi cbind() e rbind() uniscono vettori (o matrici) tra loro, per colonne o per righe rispettivamente. Attenzione! le matrici devono comunque avere entrate delle stesse classi (se non lo sono vengono convertite).

matEstesa <- cbind(mat, c(1,3,10))

matEstesa

matEstesaZeri <- rbind(c(0,0,0), matEstesa)

matEstesaZeri

Data frame

Quando una tabella ha colonne di classi diverse, si usa un data frame. Questa struttura è estremamente utile per rappresentare dati (una osservazione per ciascuna riga) di cui si osservano caratteristiche multiple (una caratteristica per ciascuna colonna). La funzione per definire un data frame partendo da vettori (tutti con la stessa lunghezza) è data.frame.

altezza <- c(148, 170, NA, 179, 190, 168, 181, 158, 166)

# creiamo un data frame con i dati di altezza e colore 

osservazioni <- data.frame(altezza, colori)

class(osservazioni)

Con il comando head() possiamo visualizzare le prime righe. Altri comandi permettono di ottenerne le dimensioni

head(osservazioni)
nrow(osservazioni)
ncol(osservazioni)  
dim(osservazioni)   

I nomi delle colonne possono essere visualizzati e pure cambiati con la funzione colnames()

colnames(osservazioni)
colnames(osservazioni) <- c("height", "colors")
head(osservazioni)

Si può accedere alle singole colonne con la sintassi nomeDataFrame$nomeColonna. Ad esempio:

osservazioni$height

# la colonna selezionata è ora un vettore

class(osservazioni$height)

Se si vuole selezionare alcune colonne mantenendo la struttura di data frame, si può indicarne ad esempio il numero.

head(osservazioni[1])

class(osservazioni[1])

Ci sono diverse estensioni della struttura dei data frames. Una sono le data tables, per gestire grandi quantità di dati (installare pacchetto data.table). Un’altra sono le tibble, per funzioni più intuitive e codice più pulito e interpretabile (installare la suite tidyverse https://www.tidyverse.org/).

install.packages("tidyverse")

Liste

Infine, R ha liste di vettori (anche di lunghezze e classi diverse).

lista <- list(tempo = 1:40)
lista$prezzo = rnorm(10, 1,3)

class(lista)

lista$tempo[1:4]
lista$prezzo[1:4]

Le liste non sono efficienti per lavorare con grandi quantità di osservazioni (usare data frames o data tables), ma è bene sapere che alcune funzioni restituiscono una lista.

Caricare e salvare dati (comandi di base)

Le funzioni di input/output possono essere complicate dal fatto che ci sono molti formati per salvare i dati. Se si usano formati standard come .csv (comma separated values) o .tsv (tab separated values) si può usare la funzione di base read.csv() con la sintassi nomeDataFrame <- read.csv("nomeFile.csv"). L’importante è che il file si trovi nella cartella di lavoro corrente (working directory).

# per sapere la cartella di lavoro (di default quella del progetto R su cui state lavorando)

getwd()

# per modificare la cartella di lavoro usare setwd()

?setwd

# se il file non è separato da virgola ma da altri simboli, consultare l'help

?read.csv

Per salvare un data frame come file .csv basta usare il comando write.csv(nomeDataFrame, "nomeFile.csv"). Vedremo esempi più avanti (anche su come caricare altri formati ad esempio Excel .xlsx)

Plot

Il comando di base per le funzioni grafiche è plot(). Questo produce un diagramma a nuvola di punti (scatterplot).

head(iris)

plot(iris$Sepal.Length, iris$Sepal.Width)

Il comando hist() produce un istogramma.

hist(iris$Sepal.Length)

Il comando boxplot() produce un diagramma a scatola.

boxplot(iris$Sepal.Length)

Questi plot tuttavia sono piuttosto di base, il pacchetto ggplot2 (suite tidyverse) permette funzionalità più avanzate e comandi facilmente intepretabili.

Carichiamo il pacchetto (se non è installato usare install.packages())

library("ggplot2")

Visualizziamo i plot di prima aggiungendo colori in base alle specie (ultima colonna). Scatterplot:

ggplot(data=iris, aes(x=Sepal.Length, y=Sepal.Width, colour=Species)) +
  geom_point() +
  xlab("Lunghezza sepali") +
  ylab("Lunghezza petali") +
  ggtitle("Dataset Iris") 

Istogramma:

ggplot(data=iris, aes(x=Sepal.Length, fill=Species)) +
  geom_histogram(bins=10)+
  xlab("Lunghezza sepali") +
  ylab("frequenza assoluta") +
  ggtitle("Dataset Iris") 

Possiamo rappresentare anche con una densità continua:

ggplot(data=iris, aes(x=Sepal.Length, fill=Species)) +
  geom_density(position="stack")+
  xlab("Lunghezza sepali") +
  ylab("frequenza assoluta") +
  ggtitle("Dataset Iris") 

Possiamo anche rappresentare tre istogrammi (uno per specie) uno accanto all’altro.

ggplot(data=iris, aes(x=Sepal.Length, fill=Species)) +
  geom_histogram(bins=10, position=position_dodge())+
  xlab("Lunghezza sepali") +
  ylab("frequenza assoluta") +
  ggtitle("Dataset Iris") 

E pure le tre densità:

ggplot(data=iris, aes(x=Sepal.Length, fill=Species)) +
  geom_density(alpha=0.4)+
  xlab("Lunghezza sepali") +
  ylab("frequenza assoluta") +
  ggtitle("Dataset Iris") 

Boxplot:

ggplot(iris, aes(x=Sepal.Length, fill=Species))+
  geom_boxplot() +
  xlab("Lunghezza sepali") +
  ylab("") +
  ggtitle("Dataset Iris")