R中非常缓慢的功能

Question

我试图实现在R. 这里对数似然函数是我用过的功能（我是新来的R）

f <- function(t) 
{ 
s=0 
x=d 
l = dim(x)[1] 
for (i in 1:l) 
    { 
     vector = d[i,] 
     lin_res = t[1] + t[2] * vector[2] + t[3] * vector[3] 
     yi = vector[1] 
     s = s + yi*lin_res - log(1 + exp(lin_res)) 
    } 
return (s[1,1]) 
} 

当d是小矩阵以下数据：

y x1   x2 x3   x4 
1 0 1 0.29944294 5.0 0.71049142 
2 0 2 0.12521669 6.0 0.20554934 
3 1 3 0.97254701 3.0 0.43665094 
4 0 4 0.79952796 1.0 0.64749898 
5 0 5 0.77358425 9.0 0.57564913 
6 0 6 0.09983754 5.0 0.32164782 
7 1 7 0.46133893 10.0 0.86437213 
8 0 8 0.59833493 20.0 0.72545982 
9 0 9 0.80005524 80.0 0.35782812 
10 0 10 0.02979412 115.0 0.76707371 
11 1 11 0.70576655 1.5 0.96908006 
12 0 12 0.67138962 2.0 0.37169164 
13 0 13 0.33446510 8.0 0.23591223 
14 1 14 0.72187427 2.0 0.98578941 
15 0 15 0.28193852 200.0 0.87076869 
16 1 16 0.11258881 3.0 0.05566943 
17 0 17 0.22001868 100.0 0.98197495 
18 1 18 0.54681964 4.0 0.53437931 
19 0 19 0.03336023 5.0 0.26451825 
20 1 20 0.47007378 10.0 0.28463580 

由于某种原因，此功能需要很多时间（运行此功能100次需要约7秒）。

d <- structure(list(y = c(0L, 0L, 1L, 0L, 0L, 0L, 1L, 0L, 0L, 0L, 
1L, 0L, 0L, 1L, 0L, 1L, 0L, 1L, 0L, 1L), x1 = 1:20, x2 = c(0.299442944, 
0.125216695, 0.972547007, 0.799527959, 0.773584254, 0.099837539, 
0.461338927, 0.59833493, 0.800055241, 0.029794123, 0.705766552, 
0.671389622, 0.334465098, 0.721874271, 0.281938515, 0.112588815, 
0.220018683, 0.546819639, 0.033360232, 0.470073781), x3 = c(5, 
6, 3, 1, 9, 5, 10, 20, 80, 115, 1.5, 2, 8, 2, 200, 3, 100, 4, 
5, 10), x4 = c(0.710491422, 0.20554934, 0.436650943, 0.647498983, 
0.575649134, 0.321647815, 0.864372135, 0.725459824, 0.357828117, 
0.767073707, 0.969080057, 0.371691641, 0.23591223, 0.985789413, 
0.870768686, 0.055669431, 0.981974949, 0.534379314, 0.26451825, 
0.284635804)), .Names = c("y", "x1", "x2", "x3", "x4"), class = "data.frame", row.names = c(NA, 
-20L)) 

有人可以帮助我加速这个功能，或明白我做错了什么？

谢谢！

Answer 1

@Andrie：R使用C（和一些地方的Fortran）代码不是C++。

@ user5497：循环速度慢的主要原因是您正在按行访问数据帧。 你的d不是一个矩阵，而是一个数据框，从结构的类参数中可以看出。

看看这个。

F1是你的函数

f2是朱莉娅的功能

f2alt是F2与d由矩阵X替换为F4

F3是F1

F4的字节编译版本与f1相同，但将d转换为变量x中的矩阵并将矢量设置为x [i，]

f5是f4的字节编译版本

f1 <- function(t) 
{ 
s=0 
x=d 
l = dim(x)[1] 
for (i in 1:l) 
    { 
     vector = d[i,] 
     lin_res = t[1] + t[2] * vector[2] + t[3] * vector[3] 
     yi = vector[1] 
     s = s + yi*lin_res - log(1 + exp(lin_res)) 
    } 
return (s[1,1]) 
} 

f2 <- function(t) 
{ 
    lin_res = t[1] + t[2] * d[2] + t[3] * d[3] 
    return (sum(d[1]*lin_res - log(1 + exp(lin_res)))) 
} 

f2alt <- function(t) 
{ 
    x <- as.matrix(d) 
    lin_res = t[1] + t[2] * x[,2] + t[3] * x[,3] 
    return (sum(x[,1]*lin_res - log(1 + exp(lin_res)))) 
} 

library(compiler) 
f3 <- cmpfun(f1) 

f4 <- function(t) 
{ 
s <- 0 
x <- as.matrix(d) 
colnames(x) <- NULL 
l <- dim(x)[1] 
for (i in 1:l) 
    { 
     vector <- x[i,] 
     lin_res <- t[1] + t[2] * vector[2] + t[3] * vector[3] 
     yi <- vector[1] 
     s <- s + yi*lin_res - log(1 + exp(lin_res)) 
    } 
return (s) 
} 

f5 <- cmpfun(f4) 

tstart <- 1:3 

f1(tstart) 
f2(tstart) 
f2alt(tstart) 
f3(tstart) 
f4(tstart) 
f5(tstart) 
all.equal(f1(tstart),f2(tstart)) 
all.equal(f1(tstart),f2alt(tstart)) 
all.equal(f1(tstart),f3(tstart)) 
all.equal(f1(tstart),f4(tstart)) 
all.equal(f1(tstart),f5(tstart)) 

library(rbenchmark) 

benchmark(f1(tstart),f2(tstart),f2alt(tstart),f3(tstart),f4(tstart),f5(tstart),columns=c("test","elapsed","relative")) 

结果是

  test elapsed relative 
1 f1(tstart) 6.912 460.800000 
2 f2(tstart) 0.305 20.333333 
3 f2alt(tstart) 0.015 1.000000 
4 f3(tstart) 6.941 462.733333 
5 f4(tstart) 0.032 2.133333 
6 f5(tstart) 0.024 1.600000 

正如你可以字节编译你的函数几乎有差别看到。 f2很快，但f2alt，f4和f5（字节编译版本的f4）更快，只是因为它们按行访问矩阵而不是数据帧。

f2alt比原来的f2快很多，因为访问矩阵而不是数据框。

警告：我在Mac OS X上使用R-2.15.1补丁，该补丁不接受标准rbenchmark;我使用了一个稍微修改过的版本。

Answer 2

我不确定它是干什么的，但是在使用循环时总是会想如果你可以使用矢量操作。此代码返回相同的值作为功能f：

f2 <- function(t) 
{ 
    lin_res = t[1] + t[2] * d[2] + t[3] * d[3] 
    return (sum(d[1]*lin_res - log(1 + exp(lin_res)))) 
} 

随机数据t：

tt <- cbind(sample(0:100,100,replace=TRUE), sample(0:100,100,replace=TRUE), sample(0:100,100,replace=TRUE)) 

的时间在我的机器：

# original 
ptm <- proc.time() 
for (t in tt) f(t) 
p <- proc.time() - ptm 
print(p) 
# user system elapsed 
# 25.529 0.002 25.533 
# new 
ptm <- proc.time() 
for (t in tt) f2(t) 
p <- proc.time() - ptm 
print(p) 
# user system elapsed 
# 1.612 0.001 1.614