Intorno a data.table

Questa nota sarà interessante per coloro che utilizzano la libreria di elaborazione dei dati tabulari per R — data.table, e potrebbe essere contento di vedere la flessibilità della sua applicazione in vari esempi.

Ispirato da un buon esempio colleghi, e sperando che abbiate già letto il suo articolo, propongo di approfondire l'ottimizzazione del codice e delle prestazioni basate su data.table.

Introduzione: da dove proviene data.table?

È meglio iniziare a familiarizzare con la libreria da un po' più lontano, precisamente dalle strutture di dati da cui può essere creato un oggetto data.table (di seguito, DT).

Array

Codice

## arrays ---------

arrmatr <- array(1:20, c(4,5))

class(arrmatr)

typeof(arrmatr)

is.array(arrmatr)

is.matrix(arrmatr)

Una di queste strutture è un array (?base::array). Come in altri linguaggi, gli array qui sono multidimensionali. Tuttavia, è interessante notare che, ad esempio, un array bidimensionale inizia ad ereditare proprietà dalla classe matrice (?base::matrix), mentre un array unidimensionale, il che è altrettanto importante, non eredita da un vettore (?base::vector).

È importante comprendere che il tipo di dati contenuti in un qualche oggetto deve essere controllato con la funzione base::typeof, che restituisce la descrizione interna del tipo secondo R Internals — il protocollo generale del linguaggio, relativo all'originario C.

Un'altra funzione, per determinare la classe di un oggetto, base::class, restituisce per i vettori un tipo vettoriale (che si distingue per il nome da quello interno, ma consente di comprendere anche il tipo di dato).

Elenco

Da un array bidimensionale, ovvero una matrice, è possibile passare a una lista (?base::list).

Codice

## lists ------------------

mylist <- as.list(arrmatr)

is.vector(mylist)

is.list(mylist)

In questo processo accadono diverse cose contemporaneamente:

• Si riduce la seconda dimensione della matrice, quindi otteniamo sia una lista che un vettore.
• La lista, in tal modo, eredita da queste classi. È importante notare che a ogni elemento della lista corrisponderà un (valore scalare) dalla cella della matrice-array.

Grazie al fatto che la lista è anche un vettore, si possono applicare alcune funzioni ai vettori.

DataFrame

Da una lista, matrice o vettore è possibile passare a un DataFrame (?base::data.frame).

Codice

## data.frames ------------

df <- as.data.frame(arrmatr)
df2 <- as.data.frame(mylist)

is.list(df)

df$V6 <- df$V1 + df$V2

Ciò che è interessante è che il DataFrame eredita dalla lista! Le colonne del DataFrame sono celle della lista. Questo sarà importante in seguito, quando utilizzeremo funzioni applicate alle liste.

data.table

Ottenere un DT (?data.table::data.table) può avvenire da un DataFrame, una lista, un vettore o una matrice. Ad esempio, in questo modo (in place).

Codice

## data.tables -----------------------
library(data.table)

data.table::setDT(df)

is.list(df)

is.data.frame(df)

is.data.table(df)

È utile sapere che, proprio come il DataFrame, il DT eredita le proprietà della lista.

DT e memoria

A differenza di tutti gli altri oggetti in R base, i data.table vengono passati per riferimento. Se è necessario fare una copia in una nuova area di memoria, è necessaria una funzione data.table::copy oppure è necessario fare una selezione dal vecchio oggetto.

Codice

df2 <- df

df[V1 == 1, V2 := 999]

data.table::fsetdiff(df, df2)

df2 <- data.table::copy(df)

df[V1 == 2, V2 := 999]

data.table::fsetdiff(df, df2)

Con questa introduzione concludiamo. I data.table rappresentano un'evoluzione delle strutture dati in R, principalmente grazie all'espansione e all'accelerazione delle operazioni effettuate sugli oggetti di classe data.frame. Inoltre, si mantiene l'ereditarietà da altri primitivi.

Alcuni esempi di utilizzo delle proprietà di data.table

Come una lista…

Iterare sulle righe di un data.frame o di un data.table non è la migliore idea, poiché il codice del ciclo nel linguaggio R è notevolmente più lento C, ma è possibile passare in ciclo sulle colonne, che sono generalmente molto meno numerose. Percorrendo le colonne, ricordiamo che ogni colonna è un elemento di una lista, contenente tipicamente un vettore. E le operazioni sui vettori sono ben vettorizzate nelle funzioni di base del linguaggio. Inoltre, si possono utilizzare operatori di selezione tipici delle liste e dei vettori: `[[`, `$`.

Codice

## operations on data.tables ------------

#using list properties

df$'V1'[1]

df[['V1']]

df[[1]][1]

sapply(df, class)

sapply(df, function(x) sum(is.na(x)))

Vettorizzazione

Se c'è la necessità di scorrere le righe di un grande dataset, la soluzione migliore sarebbe scrivere una funzione con vettorizzazione. Ma se ciò non è possibile, è importante ricordare che il ciclo all'interno è comunque più veloce del ciclo in R, poiché viene eseguito su C.

Proviamo un esempio più grande con 100K righe. Estrarremo la prima lettera dalle parole che compongono la colonna vettore w.

Aggiornato

Codice

library(magrittr)
library(microbenchmark)

## Esempio più grande ----

rown <- 100000

dt %
	.[, d := 1 + b + c + rnorm(nrow(.))]

# vettorizzazione

microbenchmark({
	dt[
		, first_l := unlist(strsplit(w, split = ' ', fixed = TRUE))[1]
		, by = 1:nrow(dt)
	   ]
})

# seconda

first_l_f %
		do.call(rbind, .) %>%
		`[`(,1)
}

dt[, first_l := NULL]

microbenchmark({
	dt[
		, first_l := .(first_l_f(w))
		]
})

# terza

first_l_f2 %
		unlist %>%
		matrix(nrow = 3) %>%
		`[`(1,)
}

dt[, first_l := NULL]

microbenchmark({
	dt[
		, first_l := .(first_l_f2(w))
		]
})

Primo passaggio con iterazione sulle righe:

Unità: millisecondi
espresse minime
{ dt[, `:=`(first_l, unlist(strsplit(w, split = " ", fixed = TRUE))[1]), by = 1:nrow(dt)] } 439.6217
lq media mediana uq max neval
451.9998 460.1593 456.2505 460.9147 621.4042 100

Second run, where vectorization is done by converting the list into a matrix and taking elements with index 1 (this is essentially the vectorization). I stand corrected: vectorization at the function level. strsplit, which can accept a vector as input. It turns out that the procedure for converting a list into a matrix is much heavier than the vectorization itself, but in this case, it is still significantly faster than the non-vectorized version.

Unità: millisecondi
expr min lq mean median uq max neval
{ dt[, `:=`(first_l, .(first_l_f(w)))] } 93.07916 112.1381 161.9267 149.6863 185.9893 442.5199 100

Acceleration by median in 3 times.

Third run, where the conversion scheme to matrix has been changed.

Unità: millisecondi
expr min lq mean median uq max neval
{ dt[, `:=`(first_l, .(first_l_f2(w)))] } 32.60481 34.13679 40.4544 35.57115 42.11975 222.972 100

Acceleration by median in 13 times.

This requires some experimentation; the more you do, the better it will be.

Another example of vectorization, where there is also text, but it is closer to real conditions: different word lengths, different numbers of words. We need to extract the first 3 words. Like this:

Here, the previous function no longer works because the vectors are of different lengths, while we specified the matrix size. Let's redo this by searching through the internet.

Codice

# fourth

rown <- 100000

words <-
	sapply(
		seq_len(rown)
		, function(x){
			nwords <- rbinom(1, 10, 0.5)
			paste(
				sapply(
					seq_len(nwords)
					, function(x){
						paste(sample(letters, rbinom(1, 10, 0.5), replace = T), collapse = '')
					}
				)
				, collapse = ' '
			)
		}
	)

dt <- 
	data.table(
		w = words
		, a = sample(letters, rown, replace = T)
		, b = runif(rown, -3, 3)
		, c = runif(rown, -3, 3)
		, e = rnorm(rown)
	) %>%
	.[, d := 1 + b + c + rnorm(nrow(.))]

first_l_f3 <- function(sd, n)
{
	l <- strsplit(sd, split = ' ', fixed = T)
	
	maxl <- max(lengths(l))
	
	sapply(l, "length<-", maxl) %>%
		`[`(n,) %>%
		as.character
}

microbenchmark({
	dt[
		, (paste0('w_', 1:3)) := lapply(1:3, function(x) first_l_f3(w, x))
		]
})

dt[
	, (paste0('w_', 1:3)) := lapply(1:3, function(x) first_l_f3(w, x))
	]

Unità: millisecondi
expr min lq mean median

{ dt[, `:=`((paste0(«w_», 1:3)), strsplit(w, split = " ", fixed = T))] } 851.7623 916.071 1054.5 1035.199
uq max neval
1178.738 1356.816 100

Lo script ha funzionato con una velocità media di 1 secondo. Non male.

Collegati a una catena...

Puoi lavorare con gli oggetti DT utilizzando il chaining. Si presenta come un'applicazione della sintassi delle parentesi a destra, in sostanza, è una dolcezza.

Codice

# chaining

res1 <- dt[a == 'a'][sample(.N, 100)]

res2 <- dt[, .N, a][, N]

res3 <- dt[, coefficients(lm(e ~ d))[1], a][, .(letter = a, coef = V1)]

Scorre attraverso i tubi...

Operazioni simili possono essere eseguite tramite piping, sembra simile, ma funzionalmente è più ricco, poiché puoi utilizzare qualsiasi metodo, non solo DT. Calcoleremo i coefficienti di regressione logistica per i nostri dati sintetici con una serie di filtri su DT.

Codice

# piping

samplpe_b <- dt[a %in% head(letters), sample(b, 1)]

res4 <- 
	dt %>%
	.[a %in% head(letters)] %>%
	.[, 
	  {
	  	dt0 <- .SD[1:100]
	  	
	  	quants <- 
	  		dt0[, c] %>%
	  		quantile(seq(0.1, 1, 0.1), na.rm = T)
	  	
	  	.(q = quants)
	  }
	  , .(cond = b > samplpe_b)
	  ] %>%
	glm(
		cond ~ q -1
		, family = binomial(link = "logit")
		, data = .
	) %>%
	summary %>%
	.[[12]]

Statistica, machine learning e altro all'interno di DT

Puoi utilizzare funzioni lambda, ma a volte è meglio crearle separatamente, descrivere tutto il pipeline di analisi dei dati, e poi avanti: funzionano all'interno di DT. L'esempio è arricchito con tutte le funzionalità sopra menzionate, più alcune cose utili dall'arsenale di DT (come la possibilità di fare riferimento a DT all'interno di DT tramite link, a volte inseriti non in modo sequenziale, ma affinché siano disponibili).

Codice

# function

rm(lm_preds)

lm_preds <- function(
	sd, by, n
)
{
	
	if(
		n < 100 | 
		!by[['a']] %in% head(letters, 4)
	   )
	{
		
		res <-
			list(
				low = NA
				, mean = NA
				, high = NA
				, coefs = NA
			)
		
	} else {

		lmm <- 
			lm(
				d ~ c + b
				, data = sd
			)
		
		preds <- 
			stats::predict.lm(
				lmm
				, sd
				, interval = "prediction"
				)
		
		res <-
			list(
				low = preds[, 2]
				, mean = preds[, 1]
				, high = preds[, 3]
				, coefs = coefficients(lmm)
			)
	}

	res
	
}

res5 <- 
	dt %>%
	.[e < 0] %>%
	.[.[, .I[b > 0]]] %>%
	.[, `:=` (
		low = as.numeric(lm_preds(.SD, .BY, .N)[[1]])
		, mean = as.numeric(lm_preds(.SD, .BY, .N)[[2]])
		, high = as.numeric(lm_preds(.SD, .BY, .N)[[3]])
		, coef_c = as.numeric(lm_preds(.SD, .BY, .N)[[4]][1])
		, coef_b = as.numeric(lm_preds(.SD, .BY, .N)[[4]][2])
		, coef_int = as.numeric(lm_preds(.SD, .BY, .N)[[4]][3])
	)
	, a
	] %>%
	.[!is.na(mean), -'e', with = F]


# plot

plo <- 
	res5 %>%
	ggplot +
	facet_wrap(~ a) +
	geom_ribbon(
		aes(
			x = c * coef_c + b * coef_b + coef_int
			, ymin = low
			, ymax = high
			, fill = a
		)
		, size = 0.1
		, alpha = 0.1
	) +
	geom_point(
		aes(
			x = c * coef_c + b * coef_b + coef_int
			, y = mean
			, color = a
		)
		, size = 1
	) +
	geom_point(
		aes(
			x = c * coef_c + b * coef_b + coef_int
			, y = d
		)
		, size = 1
		, color = 'black'
	) +
	theme_minimal()

print(plo)

Conclusione

Spero di aver creato un quadro coerente, anche se, ovviamente, non completo, di un oggetto come data.table, partendo dalle sue proprietà legate all'ereditarietà delle classi R e arrivando alle sue caratteristiche uniche e all'ambiente degli elementi tidyverse. Spero che questo possa aiutarvi a studiare e applicare meglio questa libreria per lavoro e divertimento.

Grazie!

Codice completo

Codice

## load libs ----------------

library(data.table)
library(ggplot2)
library(magrittr)
library(microbenchmark)


## arrays ---------

arrmatr <- array(1:20, c(4,5))

class(arrmatr)

typeof(arrmatr)

is.array(arrmatr)

is.matrix(arrmatr)


## lists ------------------

mylist <- as.list(arrmatr)

is.vector(mylist)

is.list(mylist)


## data.frames ------------

df <- as.data.frame(arrmatr)

is.list(df)

df$V6 <- df$V1 + df$V2


## data.tables -----------------------

data.table::setDT(df)

is.list(df)

is.data.frame(df)

is.data.table(df)

df2 <- df

df[V1 == 1, V2 := 999]

data.table::fsetdiff(df, df2)

df2 <- data.table::copy(df)

df[V1 == 2, V2 := 999]

data.table::fsetdiff(df, df2)


## operations on data.tables ------------

#using list properties

df$'V1'[1]

df[['V1']]

df[[1]][1]

sapply(df, class)

sapply(df, function(x) sum(is.na(x)))


## Bigger example ----

rown <- 100000

dt <- 
	data.table(
		w = sapply(seq_len(rown), function(x) paste(sample(letters, 3, replace = T), collapse = ' '))
		, a = sample(letters, rown, replace = T)
		, b = runif(rown, -3, 3)
		, c = runif(rown, -3, 3)
		, e = rnorm(rown)
	) %>%
	.[, d := 1 + b + c + rnorm(nrow(.))]

# vectorization

# zero - for loop

microbenchmark({
	for(i in 1:nrow(dt))
		{
		dt[
			i
			, first_l := unlist(strsplit(w, split = ' ', fixed = T))[1]
		]
	}
})

# first

microbenchmark({
	dt[
		, first_l := unlist(strsplit(w, split = ' ', fixed = T))[1]
		, by = 1:nrow(dt)
	   ]
})

# second

first_l_f <- function(sd)
{
	strsplit(sd, split = ' ', fixed = T) %>%
		do.call(rbind, .) %>%
		`[`(,1)
}

dt[, first_l := NULL]

microbenchmark({
	dt[
		, first_l := .(first_l_f(w))
		]
})

# third

first_l_f2 <- function(sd)
{
	strsplit(sd, split = ' ', fixed = T) %>%
		unlist %>%
		matrix(nrow = 3) %>%
		`[`(1,)
}

dt[, first_l := NULL]

microbenchmark({
	dt[
		, first_l := .(first_l_f2(w))
		]
})

# fourth

rown <- 100000

words <-
	sapply(
		seq_len(rown)
		, function(x){
			nwords <- rbinom(1, 10, 0.5)
			paste(
				sapply(
					seq_len(nwords)
					, function(x){
						paste(sample(letters, rbinom(1, 10, 0.5), replace = T), collapse = '')
					}
				)
				, collapse = ' '
			)
		}
	)

dt <- 
	data.table(
		w = words
		, a = sample(letters, rown, replace = T)
		, b = runif(rown, -3, 3)
		, c = runif(rown, -3, 3)
		, e = rnorm(rown)
	) %>%
	.[, d := 1 + b + c + rnorm(nrow(.))]

first_l_f3 <- function(sd, n)
{
	l <- strsplit(sd, split = ' ', fixed = T)
	
	maxl <- max(lengths(l))
	
	sapply(l, "length<-", maxl) %>%
		`[`(n,) %>%
		as.character
}

microbenchmark({
	dt[
		, (paste0('w_', 1:3)) := lapply(1:3, function(x) first_l_f3(w, x))
		]
})

dt[
	, (paste0('w_', 1:3)) := lapply(1:3, function(x) first_l_f3(w, x))
	]


# chaining

res1 <- dt[a == 'a'][sample(.N, 100)]

res2 <- dt[, .N, a][, N]

res3 <- dt[, coefficients(lm(e ~ d))[1], a][, .(letter = a, coef = V1)]

# piping

samplpe_b <- dt[a %in% head(letters), sample(b, 1)]

res4 <- 
	dt %>%
	.[a %in% head(letters)] %>%
	.[, 
	  {
	  	dt0 <- .SD[1:100]
	  	
	  	quants <- 
	  		dt0[, c] %>%
	  		quantile(seq(0.1, 1, 0.1), na.rm = T)
	  	
	  	.(q = quants)
	  }
	  , .(cond = b > samplpe_b)
	  ] %>%
	glm(
		cond ~ q -1
		, family = binomial(link = "logit")
		, data = .
	) %>%
	summary %>%
	.[[12]]


# function

rm(lm_preds)

lm_preds <- function(
	sd, by, n
)
{
	
	if(
		n < 100 | 
		!by[['a']] %in% head(letters, 4)
	   )
	{
		
		res <-
			list(
				low = NA
				, mean = NA
				, high = NA
				, coefs = NA
			)
		
	} else {

		lmm <- 
			lm(
				d ~ c + b
				, data = sd
			)
		
		preds <- 
			stats::predict.lm(
				lmm
				, sd
				, interval = "prediction"
				)
		
		res <-
			list(
				low = preds[, 2]
				, mean = preds[, 1]
				, high = preds[, 3]
				, coefs = coefficients(lmm)
			)
	}

	res
	
}

res5 <- 
	dt %>%
	.[e < 0] %>%
	.[.[, .I[b > 0]]] %>%
	.[, `:=` (
		low = as.numeric(lm_preds(.SD, .BY, .N)[[1]])
		, mean = as.numeric(lm_preds(.SD, .BY, .N)[[2]])
		, high = as.numeric(lm_preds(.SD, .BY, .N)[[3]])
		, coef_c = as.numeric(lm_preds(.SD, .BY, .N)[[4]][1])
		, coef_b = as.numeric(lm_preds(.SD, .BY, .N)[[4]][2])
		, coef_int = as.numeric(lm_preds(.SD, .BY, .N)[[4]][3])
	)
	, a
	] %>%
	.[!is.na(mean), -'e', with = F]


# plot

plo <- 
	res5 %>%
	ggplot +
	facet_wrap(~ a) +
	geom_ribbon(
		aes(
			x = c * coef_c + b * coef_b + coef_int
			, ymin = low
			, ymax = high
			, fill = a
		)
		, size = 0.1
		, alpha = 0.1
	) +
	geom_point(
		aes(
			x = c * coef_c + b * coef_b + coef_int
			, y = mean
			, color = a
		)
		, size = 1
	) +
	geom_point(
		aes(
			x = c * coef_c + b * coef_b + coef_int
			, y = d
		)
		, size = 1
		, color = 'black'
	) +
	theme_minimal()

print(plo)

Fonte: habr.com