Saltar al contenido

Variación móvil en R

Solución:

Considera el zoo paquete. Por ejemplo filter() da:

> filter(1:100, rep(1/3,3))
Time Series:
Start = 1 
End = 100 
Frequency = 1 
  [1] NA  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
 [26] 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
 [51] 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75
 [76] 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 NA

mientras que rollmean() en zoo da:

> rollmean(1:100, k = 3, na.pad = TRUE)
  [1] NA  2  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
 [26] 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50
 [51] 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75
 [76] 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 NA

que es lo mismo (para una media móvil de 3 puntos en este ejemplo).

Mientras que zoo no tiene un rollsd() o rollvar() tiene rollapply(), que funciona como el apply() funciones para aplicar cualquier función R a la ventana especificada.

> rollapply(1:100, width = 3, FUN = sd, na.pad = TRUE)
  [1] NA  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1
 [26]  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1
 [51]  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1
 [76]  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1  1 NA
Warning message:
In rollapply.zoo(zoo(data), ...) : na.pad argument is deprecated

o en algo más interesante:

> rollapply(vec, width = 3, FUN = sd, na.pad = TRUE)
  [1]        NA 0.3655067 0.8472871 0.5660495 0.3491970 0.4732417 0.9236859
  [8] 0.8075226 1.8725851 1.1930784 0.6329325 1.1412416 0.8430772 0.5808005
 [15] 0.3838545 1.1738170 1.1655400 1.3241700 0.6876834 0.1534157 0.4858477
 [22] 0.9843506 0.6002713 0.6897541 2.0619563 2.5675788 6.3522039 6.0066864
 [29] 6.2618432 5.1704866 2.1360853 2.5602557 1.0408528 1.0316396 4.9441628
 [36] 5.0319314 5.7589716 3.2425000 4.8788158 2.0847286 4.5199291 2.5323486
 [43] 2.1987149 1.8393000 1.2278639 1.5998965 1.5341485 4.4287108 4.4159166
 [50] 4.3224546 3.6959067 4.9826264 5.3134044 8.4084322 9.1249234 7.5506725
 [57] 3.8499136 3.9680487 5.6362296 4.9124095 4.3452706 4.0227141 4.5867559
 [64] 4.7350394 4.3203807 4.4506799 7.2759499 7.6536424 7.8487654 2.0905576
 [71] 4.0056880 5.6209853 1.5551659 1.3615268 2.8469458 2.8323588 1.9848578
 [78] 1.1201124 1.4248380 1.7802571 1.4281773 2.5481935 1.8554451 1.0925410
 [85] 2.1823722 2.2788755 2.4205378 2.0733741 0.7462248 1.3873578 1.4265948
 [92] 0.7212619 0.7425993 1.0696432 2.4520585 3.0555819 3.1000885 1.0945292
 [99] 0.3726928        NA
Warning message:
In rollapply.zoo(zoo(data), ...) : na.pad argument is deprecated

Puede deshacerse de la advertencia utilizando el fill = NA argumento, como en

> rollapply(vec, width = 3, FUN = sd, fill = NA)

El paquete TTR tiene runSD entre otros:

> library(TTR)
> ls("package:TTR", pattern="run*")
 [1] "runCor"    "runCov"    "runMAD"    "runMax"    "runMean"  
 [6] "runMedian" "runMin"    "runSD"     "runSum"    "runVar"

runSD será mucho más rápido que rollapply porque evita hacer muchas llamadas a funciones R. Por ejemplo:

rzoo <- function(x,n) rollapplyr(x, n, sd, fill=NA)
rttr <- function(x,n) runSD(x, n)
library(rbenchmark)
set.seed(21)
x <- rnorm(1e4)
all.equal(rzoo(x,250), rttr(x,250))
# [1] TRUE
benchmark(rzoo(x,250), rttr(x,250))[,1:6]
#           test replications elapsed relative user.self sys.self
# 2 rttr(x, 250)          100    0.58    1.000      0.58     0.00
# 1 rzoo(x, 250)          100   54.53   94.017     53.85     0.06

rollapply en el zoo el paquete toma una función arbitraria. Es diferente a filter en que excluye el NAs por defecto.

Dicho esto, sin embargo, no tiene mucho sentido cargar un paquete para una función que es tan fácil de usar (juego de palabras).

Aquí hay uno que está alineado a la derecha:

my.rollapply <- function(vec, width, FUN) 
    sapply(seq_along(vec), 
           function(i) if (i < width) NA else FUN(vec[i:(i-width+1)]))

set.seed(1)
vec <- sample(1:50, 50)
my.rollapply(vec, 3, sd)
 [1]        NA        NA  7.094599 12.124356 16.522712 18.502252 18.193405  7.234178  8.144528
[10] 14.468356 12.489996  3.055050 20.808652 19.467922 18.009257 18.248288 15.695010  7.505553
[19] 10.066446 11.846237 17.156146  6.557439  5.291503 23.629078 22.590558 21.197484 22.810816
[28] 24.433583 19.502137 16.165808 11.503623 12.288206  9.539392 13.051181 13.527749 19.974984
[37] 19.756855 17.616280 19.347696 18.248288 15.176737  6.082763 10.000000 10.016653  4.509250
[46]  2.645751  1.527525  5.291503 10.598742  6.557439

# rollapply output for comparison
rollapply(vec, width=3, sd, fill=NA, align='right')
 [1]        NA        NA  7.094599 12.124356 16.522712 18.502252 18.193405  7.234178  8.144528
[10] 14.468356 12.489996  3.055050 20.808652 19.467922 18.009257 18.248288 15.695010  7.505553
[19] 10.066446 11.846237 17.156146  6.557439  5.291503 23.629078 22.590558 21.197484 22.810816
[28] 24.433583 19.502137 16.165808 11.503623 12.288206  9.539392 13.051181 13.527749 19.974984
[37] 19.756855 17.616280 19.347696 18.248288 15.176737  6.082763 10.000000 10.016653  4.509250
[46]  2.645751  1.527525  5.291503 10.598742  6.557439
¡Haz clic para puntuar esta entrada!
(Votos: 0 Promedio: 0)


Tags : /

Utiliza Nuestro Buscador

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *