Améliorer la fonction de performance

Question

J'apprends des méthodes pour améliorer la performance d'un code. Pour l'étude que j'ai écrit la fonction qui retourne la statistique descriptive de base comme psych::describe. J'ai essayé différentes versions des boucles et pour le moment c'est tout ce que je pouvais faire.

code:

x <- matrix(rnorm(10*100), nrow=100) # sample data for tests 

descStats <- function(x, na.rm = TRUE, trim = NULL, skew = FALSE, byrow = FALSE, digits =  getOption("digits")) { 
if (!is.matrix(x)) x <- as.matrix(x) 
if(byrow) x <- t(x) 
stats <- c("n", "mean", "se", "sd", "median", "min", "max", "range") # descriptive statistics 
if (skew) { 
    library(moments) 
    stats <- c(stats, "skewness", "kurtosis") 
} 
if (!is.null(trim)) { 
    stats <- append(stats, "trimmed", which(stats == "mean")) 
    trimmed <- function(x) base::mean(x, trim=trim) 
} 
n <- function(x) length(x) 
range <- function(x) max(x) - min(x) 
mean <- function(x) .Internal(mean(x)) # redefined mean function 
sd <- function(x) sqrt(sum((x - mean(x))^2)/(length(x)-1)) # redefined sd function 
se <- function(x) sqrt(sd(x)/length(x)) 
median <- function(x) { # redefined median function 
    n <- length(x) 
    half <- (n + 1L)%/%2L 
    if (n%%2L == 1L) 
    result <- .Internal(sort(x, partial = half))[half] 
    else { 
    result <- mean(.Internal(sort(x, partial = half + 0L:1L))[half + 0L:1L]) 
    } 
} 
describe <- function(x, na.rm=FALSE) { 
    if (na.rm) x <- x[!is.na(x)] 
    result <- vapply(stats, function(fun) eval(call(fun, x)), FUN.VALUE=numeric(1)) 
    return(result) 
} 
out <- t(vapply(seq_len(ncol(x)), function(i) describe(x[,i], na.rm=na.rm), FUN.VALUE=numeric(length(stats)))) 
out <- round(out, digits=digits) 
return(out) 
} 

print(descStats(x)) 
##  n  mean trimmed  se  sd  median  min  max range 
## [1,] 100 0.2524298 0.2763559 0.1024722 1.0500560 0.2842625 -2.905826 3.362598 6.268424 
## [2,] 100 -0.1201740 -0.0627668 0.1027268 1.0552801 -0.0614541 -3.071836 2.247063 5.318899 
## [3,] 100 0.2074781 0.1946393 0.1006384 1.0128089 0.1928790 -2.312749 2.564297 4.877047 
## [4,] 100 0.1088077 0.1127540 0.0935370 0.8749172 0.0864728 -2.757226 2.883687 5.640913 
## [5,] 100 -0.2163515 -0.2147170 0.1064167 1.1324524 -0.2836884 -3.431254 2.950466 6.381720 
## [6,] 100 -0.0324696 -0.0229878 0.0968330 0.9376630 0.0919468 -2.474992 1.860961 4.335953 
## [7,] 100 -0.1497724 -0.1665687 0.1047835 1.0979579 -0.1753578 -2.908781 2.885645 5.794425 
## [8,] 100 -0.0197306 0.0101194 0.1030385 1.0616927 0.0615438 -2.711356 2.506423 5.217779 
## [9,] 100 -0.0346922 -0.0290022 0.1018726 1.0378033 0.0231049 -2.467852 2.528595 4.996447 
## [10,] 100 0.1251403 0.1222156 0.1012441 1.0250359 0.1606492 -2.566209 2.854519 5.420728 

Dans chaque cas, je compare le temps écoulé avec microbenchmaark. Par exemple:

library(microbenchmark) 
bench <- microbenchmark(descStats(x), descStats2(x), times=1000) 
print(bench) 
boxplot(bench, outline=FALSE) 

Quelqu'un peut-il être en mesure d'offrir une version plus efficace ou plus compacte du code?

Mise à jour:

La version finale de cette fonction vous pouvez voir here.

Answer 1

Voici une version qui est ~ 60% plus rapide sur votre objet de test, x et qui calcule toujours l'asymétrie et l'aplatissement (sans le package moments).

descStats2 <- function(x, na.rm = TRUE, trim = 0.1, skew = TRUE, 
    byrow = FALSE, digits = getOption("digits")) { 

    stopifnot(is.matrix(x)) 
    f <- function(x, na.rm) { 
    if(na.rm) x <- x[!is.na(x)] 
    n <- length(x) 
    s <- numeric(11) 
    s[1] <- n 
    s[2] <- sum(x)/n 
    dev <- x-s[2] 
    dev_s2 <- sum(dev^2L) 
    s[3] <- mean.default(x, trim=trim) 
    s[5] <- sqrt(dev_s2/(n-1)) 
    s[4] <- sqrt(s[5]/n) 
    s[6] <- median.default(x) 
    s[7] <- min(x) 
    s[8] <- max(x) 
    s[9] <- s[8]-s[7] 
    s[10] <- (sum(dev^3L)/n)/(dev_s2/n)^1.5 
    s[11] <- n * sum(dev^4L)/(dev_s2^2L) 
    s 
    } 
    if(byrow) x <- t(x) 
    out <- matrix(0.0, ncol(x), 11, dimnames=list(NULL, c("n", "mean", "trimmed", 
    "se", "sd", "median", "min", "max", "range", "skewness", "kurtosis"))) 
    for(i in 1:ncol(x)) 
    out[i,] <- f(x[,i], na.rm=na.rm) 
    round(out, digits=digits) 
} 

La raison pour laquelle ma fonction est plus rapide est parce que j'évite beaucoup d'appels de fonctions inutiles, et je stocke et les valeurs de réutilisation je l'ai déjà calculé (par exemple dev, dev_s2, s[5], etc.).

library(compiler) 
ds2c <- cmpfun(descStats2) 

library(rbenchmark) 
benchmark(descStats(x, skew=TRUE), descStats2(x), ds2c(x), 
    order="relative", replications=1000)[,1:5] 
#      test replications elapsed relative user.self 
# 3     ds2c(x)   1000 3.969 1.000  3.920 
# 2    descStats2(x)   1000 5.115 1.289  4.984 
# 1 descStats(x, skew = TRUE)   1000 8.234 2.075  7.837 
identical(descStats(x, skew=TRUE), descStats2(x)) 
# [1] TRUE