Vektorizácia je spôsob, ako zrýchliť vyčíslenie výrazov, v ktorých sa rovnaká operácia aplikuje veľakrát na rôzne prvky.

Uvažujme napríklad situáciu, kde chceme vypočítať súčet logaritmov 1000 čísel:

x <- abs(rnorm(1000))

Takýto problém vieme riešiť napríklad pomocou for cyklu:

log_sum <- 0
for(i in 1:length(x)) log_sum <- log_sum + log(x[i])

Lepšie je ale aplikovať funkciu sum priamo na vektor x:

log_sum <- sum(log(x))

Úloha: Otestujte rýchlosť oboch implementácií pre \(n=10^7\).

Monte Carlo integrovanie

1. implementácia

monte_carlo <- function(N) {
  hits <- 0
  for (i in seq_len(N)) {
    u1 <- runif(1)
    u2 <- runif(1)
    if (u1 ^ 2 > u2)
      hits <- hits + 1
  }
  return(hits / N)
}

Rýchlosť:

N <- 500000
system.time(monte_carlo(N))
##    user  system elapsed 
##    1.18    0.00    1.19

2. implementácia

monte_carlo_vec <- function(N) sum(runif(N)^2 > runif(N))/N

N <- 500000
system.time(monte_carlo_vec(N))
##    user  system elapsed 
##    0.01    0.00    0.01

Funkcie typu APPLY

Funkcie typu apply sú určené na repetitívnu manipuláciu s podmnožinami dát v maticiach, poliach a dátových rámcoch. Používame ich na to, aby sme sa vyhli, častokrát vnoreným, for cyklom. Funckia, ktorú týmto spôsobom aplikujeme, môže byť napríklad agregačná (priemer, súčet), transformačná alebo iná vektorizovaná funkcia.

K dispozícii máme funkcie apply(), lapply() , sapply(), vapply(), mapply(), rapply() a tapply(). V nasledujúcich príkladoch sa pozrieme na niektoré z nich.

Pri každom príklade si pozrite help danej funkcie a všimnite si, ako sa medzi sebou líšia výstupy jednotlivých funkcií.

1. apply

Ako by sme počítali riadkové alebo stĺpcové súčty matice?

X <- matrix(rnorm(30), nrow=5, ncol=6)
apply(X, 2 ,sum)
## [1] -1.4062491  0.1178578 -1.9623120  0.1005855  1.1808389 -1.7279902

Úloha: Vypočítajte riadkové priemery matice X.

Pre matice máme v R špeciálne funkcie, ktoré nám umožňujú robiť podobné často sa vyskytujúce operácie, napr. colSums alebo rowMeans. Pre polia väčších rozmerov a pre menej bežné funkcie je ale funkcia apply užitočná.

pole <- array(1:(2*3*4), dim = c(2,3,4))
apply(pole, 3, function(x) sum(x^2))
## [1]   91  559 1459 2791
apply(pole, 1:2, var)
##      [,1] [,2] [,3]
## [1,]   60   60   60
## [2,]   60   60   60
apply(pole, c(3,1), max)
##      [,1] [,2]
## [1,]    5    6
## [2,]   11   12
## [3,]   17   18
## [4,]   23   24
apply(pole, c(1,3), max)
##      [,1] [,2] [,3] [,4]
## [1,]    5   11   17   23
## [2,]    6   12   18   24

2.lapply

Funkcia apply aplikovaná na objekty typu list.

l <- list(a = 1:10, b = 11:20)
lapply(l, mean)
## $a
## [1] 5.5
## 
## $b
## [1] 15.5
lapply(l, sum)
## $a
## [1] 55
## 
## $b
## [1] 155

Príklad použitia: ak chceme vykonať rovnakú operáciu na viacerých dátových rámcoch, spojíme ich do listu a použijeme lapply.

data(iris)
data(swiss)
lst <- list(iris=iris[,-5], swiss=swiss[,-1])
lapply(lst, colMeans)
## $iris
## Sepal.Length  Sepal.Width Petal.Length  Petal.Width 
##     5.843333     3.057333     3.758000     1.199333 
## 
## $swiss
##      Agriculture      Examination        Education         Catholic 
##         50.65957         16.48936         10.97872         41.14383 
## Infant.Mortality 
##         19.94255

3.sapply

l <- list(a = 1:10, b = 11:20)
l.mean <- sapply(l, mean)
class(l.mean)
## [1] "numeric"
l.mean[1]
##   a 
## 5.5
l.mean[['a']]
## [1] 5.5

Výstupom sapply je v prípade, že oba listy sú rovnakej veľkosti, matica, nie list, čo môže často byť pohodlnejšie na následnú manipuláciu.

lst <- list(iris=iris[,-5], swiss=swiss[,-(1:2)])
sapply(lst, colMeans)
##                  iris    swiss
## Sepal.Length 5.843333 16.48936
## Sepal.Width  3.057333 10.97872
## Petal.Length 3.758000 41.14383
## Petal.Width  1.199333 19.94255

4.tapply

Funkcia tapply užitočná najmä pre dátové súbory. Vezmime si dáta obsahujúce mená detí v USA v rokoch 1880-2017. Kvôli prehľadnosti budeme uvažovať len roky 2015-2017.

library(babynames)
## Warning: package 'babynames' was built under R version 4.1.3
data(babynames)
bn <- babynames[babynames$year>=2015, ]
tapply(X=bn$n, INDEX=bn$sex, FUN=sum, na.rm=TRUE)
##       F       M 
## 5254610 5633346
tapply(X=bn$n, INDEX=list(bn$sex, bn$year), FUN=sum, na.rm=TRUE)
##      2015    2016    2017
## F 1778883 1763916 1711811
## M 1909804 1889052 1834490

Príklad: dvojitá suma

Uvažujme takúto dvojitú sumu: \[f(n)= \sum_{r=1}^n \sum_{s=1}^r \frac{s^2}{10+4r^3}\]

Ukážeme si, ako túto sumu vypočítať rôznymi (viac alebo menej efektívnymi) spôsobmi.

1.Výpočet pomocou vnorených for cyklov

suma1 <- function(n){
    sum <- 0
    for (r in 1:n){
        for (s in 1:r)
            sum <- sum+s^2/(10+4*r^3)
    }
    sum
}

2.Výpočet pomocou funkcií row a column a súčtu prvkov matice

suma2 <- function(n){
    mat <- matrix(0,ncol=n,nrow=n)
    sum( (col(mat)^2) / (10+4 * row(mat) ^ 3) * (col(mat) <= row(mat)) )
}

3.Výpočet pomocou sapply

suma3 <- function(n){
    f <- function(r){ sum(((1:r)^2 / (10 + 4 * r^3)))}
    sum(sapply(1:n, FUN=f))
}

4. Výpočet pomocou funkcie ‘outer’

suma4 <- function(n)
    sum(outer(1:n, 1:n, FUN = function(r,s){ (s<=r)*(s^2)/(10+4*r^3) }))

Pomocou system.time zistite, ktorá implementácia je najrýchlejšia a ktorá najpomalšia

system.time(suma1(1000))
##    user  system elapsed 
##    0.09    0.00    0.09