Eu fiz duas novas versões. Na primeira eu simplifiquei o código como eu faria em um programa meu e na segunda eu executei o loop em paralelo com o pacote furrr
. Não encontrei um jeito de otimizar mais que isso… O programa parecer ser O(bn), não tem muito jeito.
# Argumentos utilizados para testar os algoritmos
data <- cbind(
sample(cars$speed, 10000, replace = TRUE),
sample(cars$dist, 10000, replace = TRUE)
)
bw <- seq(10)
# ORIGINAL -----------------------------------------------------------------
### cross-validation function
cv <- function(data, used.function, bandwidth.grid)
{
n <- nrow(data)
mse <- matrix(, nrow=length(bandwidth.grid), ncol=2)
for (b in 1:length(bandwidth.grid)){
cv.value <- rep(0, n-1)
for (i in 1:(n-1)){ # *
new.data <- data[-i,]
funcargs <- list(reg.x=new.data[,1],reg.y=new.data[,2],x=data[i,1], h = bandwidth.grid[b])
cv.value[i] <- do.call(used.function, funcargs)
}
mse[b,] <- c(bandwidth.grid[b], 1/n*sum((new.data[,2]-cv.value)^2)) # *
}
## MSE
colnames(mse) <- c("bandwidth", "mse")
return(mse)
}
### kernel estimator usind nadaraya-watson:
fcn1 <- function(reg.x, reg.y, x, h){
return(ksmooth(reg.x, reg.y, x.point = x, kernel = "normal", bandwidth = h)$y)
}
# MAP ----------------------------------------------------------------------
fcn2 <- function(reg_x, reg_y, x, h) {
ksmooth(reg_x, reg_y, "normal", h, x.points = x)$y
}
make_cv_value <- function(i, fun, data, b) {
fun(data[-i, 1], data[-i, 2], data[i, 1], b)
}
cv2 <- function(data, fun, bandwidth_grid) {
n <- nrow(data)
mse <- c()
for (b in bandwidth_grid) {
cv_value <- purrr::map_dbl(seq_len(n - 1), make_cv_value, fun, data, b)
mse <- append(mse, 1 / n * sum((data[-n + 1, 2] - cv_value)^2))
}
mse <- cbind(seq_along(mse), mse)
colnames(mse) <- c("bandwidth", "mse")
return(mse)
}
# PARALELO -----------------------------------------------------------------
library(furrr)
plan("multicore", workers = 4)
cv3 <- function(data, fun, bandwidth_grid) {
n <- nrow(data)
mse <- c()
for (b in bandwidth_grid) {
cv_value <- future_map_dbl(seq_len(n - 1), make_cv_value, fun, data, b)
mse <- append(mse, 1 / n * sum((data[-n + 1, 2] - cv_value)^2))
}
mse <- cbind(seq_along(mse), mse)
colnames(mse) <- c("bandwidth", "mse")
return(mse)
}
# BENCHMARKS ---------------------------------------------------------------
microbenchmark::microbenchmark(
`for` = cv(data, fcn1, bw),
`map` = cv2(data, fcn2, bw),
`par` = cv3(data, fcn2, bw),
times = 10
)
#> Unit: seconds
#> expr min lq mean median uq max neval
#> for 26.466603 26.738418 27.013609 26.996009 27.334605 27.679403 10
#> map 23.987672 24.193279 24.411272 24.301028 24.686168 25.075567 10
#> par 7.743104 7.812796 7.929976 7.928873 8.021716 8.136829 10
Created on 2021-10-27 by the reprex package (v2.0.1)
Só mais duas coisas:
-
Eu marquei duas linhas da sua versão com *
. Eu não sei nada sobre esse algoritmo, mas me parece que na primeira delas você quis dizer 1:n
ao invés de 1:(n-1)
e, na segunda, data
ao invés de new.data
. Eu fiz as minhas versões baseadas na sua, mas estranhei você usar data[-n + 1, 2]
na hora de calcular mse
; por que usar a base sem a penúltima linha?
-
Eu sempre gosto de dar dicas de como as pessoas podem melhorar seus códigos. Tenho 4 sugestões para você baseadas nesse programa:
-
Não use attach()
. Ele é uma armadilha que polui o seu ambiente com objetos invisíveis! Prefira sempre usar o $
.
-
Evite colocar .
no nome das variáveis. Essa é uma prática antiga e desencorajada hoje em dia por causa das classes S3.
-
Não precisa usar do.call()
quando você passa uma função como argumento. Ela continua funcionando do mesmo jeito (veja como eu fiz no meu código).
-
Prefira nomes mais esclarecedores para as suas funções. É muito fácil se confundir quando as funções não indicam nada do que elas fazem.