我需要减少 R 中的 for 循环吗?如何减少?
Do I need to reduce for-loop in R and how?
我是 R 编码的新手,但我听说 R 中的循环比 Python 或 C 等其他语言慢得多。那么在 R 中编码时是否需要减少循环?
具体来说,在这个模拟代码中,我该如何提高自己的编码水平?
library(moments)
n <- c(5:20)
m <- c(1:10000)
skew <- c()
kurt <- c()
for(num in n){
beta1 <- c()
beta2 <- c()
for(i in m){
set.seed(num * 10000 + i)
x <- rnorm(num, mean = 0, sd = 1)
beta1 <- c(beta1, skewness(x))
beta2 <- c(beta2, kurtosis(x) - 3)
}
skew <- c(skew, quantile(beta1, probs = c(0, 0.01, 0.1, 0.2, 0.5, 0.8, 0.9, 0.99, 1)))
kurt <- c(kurt, quantile(beta2, probs = c(0, 0.01, 0.1, 0.2, 0.5, 0.8, 0.9, 0.99, 1)))
}
在 R 中不使用 for
循环的一个主要优点是利用其向量化。因此,在 Python 或 C 等语言中,您可以为向量的每个元素编写向量计算代码,而在 R 中,您可以方便地一次为整个向量编写计算代码(请参阅下面的编辑),还可以通过实际使用来减少计算时间快速底层 C、Fortran 等函数。
我会把你想为单个样本量进行的所有计算放入函数 statFUN
并将其放入 lapply
以遍历样本量向量 n
.
对于分位数,我们可以使用 apply
或我推荐的 matrixStats::rowQuantiles
,因为它更快。
set.seed()
在 运行 和 lapply
之前只需要一次,所有 res
结果都可以用那个种子重现。
n <- 5:20 ## different sample sizes
m <- 1e4 ## number of replications in each iteration
probs <- c(0, 0.01, 0.1, 0.2, 0.5, 0.8, 0.9, 0.99, 1)
library(moments)
library(matrixStats)
statFUN <- function(i, num) {
r <- replicate(i, {
x <- rnorm(num, mean=0, sd=1)
c(kurt=kurtosis(x) - 3, skew=skewness(x))
})
# t(apply(r, 1, quantile, probs=probs)) ## using base R
rowQuantiles(r, probs=probs) ## using matrixStats
}
set.seed(42)
res <- lapply(n, statFUN, m)
结果
res
ult 是每个样本大小的峰度分位数和偏度分位数的列表。
res
# [[1]]
# 0% 1% 10% 20%
# kurt -0.04710729 -0.04658709 -0.04190536 -0.03670343
# skew -0.03045563 -0.02969417 -0.02284104 -0.01522645
# 50% 80% 90% 99%
# kurt -0.03388803 -0.006250622 1.068998e-03 0.007656657
# skew -0.01028591 -0.006132523 -5.883157e-05 0.005407491
# 100%
# kurt 0.008388619
# skew 0.006014860
#
# [[2]]
# 0% 1% 10% 20%
# kurt -0.09089922 -0.08859363 -0.06784329 -0.04478737
# skew -0.03252828 -0.03165837 -0.02382918 -0.01513009
# 50% 80% 90% 99%
# kurt -0.023634727 -0.005277533 0.01038904 0.02448896
# skew 0.003433589 0.017711708 0.01947178 0.02105585
# 100%
# kurt 0.02605562
# skew 0.02123186
#
# [...]
哪里
length(res)
# [1] 16
编辑
这里有一个小例子可以更好地说明 R 中向量化的实际含义。虽然在大多数编程语言中,两个向量的加法是按元素编码的,但在 R 中,向量的加法可以直接编码(即在矢量化的方式)。
a <- 1:9
b <- rev(a)
## element wise addition of vectors a and b
s1 <- c()
for (i in seq(a)) {
s1[i] <- a[i] + b[i]
}
s1
# [1] 10 10 10 10 10 10 10 10 10
## direct addition of vectors a and b (i.e. vectorized)
s2 <- a + b
s2
# [1] 10 10 10 10 10 10 10 10 10
我们可以研究 *apply
系列而不是 for
循环。但是,大多数情况下仍然隐藏着 for 循环。 (要查看功能代码类型,例如 lapply
不带括号或任何内容。)
您可能想阅读例如那些很棒的问答:
- Grouping functions (tapply, by, aggregate) and the *apply family
- Why are loops slow in R?
- How to see the source code of R .Internal or .Primitive function?
注:向量化其实只是R的语言特性,所谓的“向量化函数”内部往往使用C、Fortran等代码,在其中你还是在最后找到 for 循环,但是使用一种更快的语言。例如,当我们使用 sum()
.
时调用 source code of summary.c
我是 R 编码的新手,但我听说 R 中的循环比 Python 或 C 等其他语言慢得多。那么在 R 中编码时是否需要减少循环?
具体来说,在这个模拟代码中,我该如何提高自己的编码水平?
library(moments)
n <- c(5:20)
m <- c(1:10000)
skew <- c()
kurt <- c()
for(num in n){
beta1 <- c()
beta2 <- c()
for(i in m){
set.seed(num * 10000 + i)
x <- rnorm(num, mean = 0, sd = 1)
beta1 <- c(beta1, skewness(x))
beta2 <- c(beta2, kurtosis(x) - 3)
}
skew <- c(skew, quantile(beta1, probs = c(0, 0.01, 0.1, 0.2, 0.5, 0.8, 0.9, 0.99, 1)))
kurt <- c(kurt, quantile(beta2, probs = c(0, 0.01, 0.1, 0.2, 0.5, 0.8, 0.9, 0.99, 1)))
}
在 R 中不使用 for
循环的一个主要优点是利用其向量化。因此,在 Python 或 C 等语言中,您可以为向量的每个元素编写向量计算代码,而在 R 中,您可以方便地一次为整个向量编写计算代码(请参阅下面的编辑),还可以通过实际使用来减少计算时间快速底层 C、Fortran 等函数。
我会把你想为单个样本量进行的所有计算放入函数 statFUN
并将其放入 lapply
以遍历样本量向量 n
.
对于分位数,我们可以使用 apply
或我推荐的 matrixStats::rowQuantiles
,因为它更快。
set.seed()
在 运行 和 lapply
之前只需要一次,所有 res
结果都可以用那个种子重现。
n <- 5:20 ## different sample sizes
m <- 1e4 ## number of replications in each iteration
probs <- c(0, 0.01, 0.1, 0.2, 0.5, 0.8, 0.9, 0.99, 1)
library(moments)
library(matrixStats)
statFUN <- function(i, num) {
r <- replicate(i, {
x <- rnorm(num, mean=0, sd=1)
c(kurt=kurtosis(x) - 3, skew=skewness(x))
})
# t(apply(r, 1, quantile, probs=probs)) ## using base R
rowQuantiles(r, probs=probs) ## using matrixStats
}
set.seed(42)
res <- lapply(n, statFUN, m)
结果
res
ult 是每个样本大小的峰度分位数和偏度分位数的列表。
res
# [[1]]
# 0% 1% 10% 20%
# kurt -0.04710729 -0.04658709 -0.04190536 -0.03670343
# skew -0.03045563 -0.02969417 -0.02284104 -0.01522645
# 50% 80% 90% 99%
# kurt -0.03388803 -0.006250622 1.068998e-03 0.007656657
# skew -0.01028591 -0.006132523 -5.883157e-05 0.005407491
# 100%
# kurt 0.008388619
# skew 0.006014860
#
# [[2]]
# 0% 1% 10% 20%
# kurt -0.09089922 -0.08859363 -0.06784329 -0.04478737
# skew -0.03252828 -0.03165837 -0.02382918 -0.01513009
# 50% 80% 90% 99%
# kurt -0.023634727 -0.005277533 0.01038904 0.02448896
# skew 0.003433589 0.017711708 0.01947178 0.02105585
# 100%
# kurt 0.02605562
# skew 0.02123186
#
# [...]
哪里
length(res)
# [1] 16
编辑
这里有一个小例子可以更好地说明 R 中向量化的实际含义。虽然在大多数编程语言中,两个向量的加法是按元素编码的,但在 R 中,向量的加法可以直接编码(即在矢量化的方式)。
a <- 1:9
b <- rev(a)
## element wise addition of vectors a and b
s1 <- c()
for (i in seq(a)) {
s1[i] <- a[i] + b[i]
}
s1
# [1] 10 10 10 10 10 10 10 10 10
## direct addition of vectors a and b (i.e. vectorized)
s2 <- a + b
s2
# [1] 10 10 10 10 10 10 10 10 10
我们可以研究 *apply
系列而不是 for
循环。但是,大多数情况下仍然隐藏着 for 循环。 (要查看功能代码类型,例如 lapply
不带括号或任何内容。)
您可能想阅读例如那些很棒的问答:
- Grouping functions (tapply, by, aggregate) and the *apply family
- Why are loops slow in R?
- How to see the source code of R .Internal or .Primitive function?
注:向量化其实只是R的语言特性,所谓的“向量化函数”内部往往使用C、Fortran等代码,在其中你还是在最后找到 for 循环,但是使用一种更快的语言。例如,当我们使用 sum()
.
summary.c