Cách nhanh nhất để áp dụng t.test cho mỗi cột của một ma trận lớn là gì?

Question

Giả sử tôi có một ma trận lớn:

M <- matrix(rnorm(1e7),nrow=20) 

Hơn nữa giả sử rằng mỗi cột đại diện cho một mẫu. Giả sử tôi muốn áp dụng t.test() cho mỗi cột, có cách nào để thực hiện việc này nhanh hơn nhiều so với sử dụng apply() không?

apply(M, 2, t.test) 

Phải mất một chút ít hơn 2 phút để chạy các phân tích trên máy tính của tôi:

> system.time(invisible(apply(M, 2, t.test))) 
user system elapsed 
113.513 0.663 113.519 

Answer 1

Nếu bạn có một máy tính đa lõi có một số lợi ích từ việc sử dụng tất cả các lõi, ví dụ sử dụng mclapply.

> library(multicore) 
> M <- matrix(rnorm(40),nrow=20) 
> x1 <- apply(M, 2, t.test) 
> x2 <- mclapply(1:dim(M)[2], function(i) t.test(M[,i])) 
> all.equal(x1, x2) 
[1] "Component 1: Component 9: 1 string mismatch" "Component 2: Component 9: 1 string mismatch" 
# str(x1) and str(x2) show that the difference is immaterial 

Ví dụ nhỏ này cho thấy mọi thứ diễn ra như chúng tôi đã lên kế hoạch. Bây giờ mở rộng quy mô:

> M <- matrix(rnorm(1e7), nrow=20) 
> system.time(invisible(apply(M, 2, t.test))) 
    user system elapsed 
101.346 0.626 101.859 
> system.time(invisible(mclapply(1:dim(M)[2], function(i) t.test(M[,i])))) 
    user system elapsed 
55.049 2.527 43.668 

Điều này đang sử dụng 8 lõi ảo. Số dặm của bạn có thể thay đổi. Không phải là một lợi ích to lớn, nhưng nó đến từ rất ít nỗ lực.

EDIT

Nếu bạn chỉ quan tâm đến những thống kê t chính nó, giải nén các lĩnh vực tương ứng ($statistic) làm cho mọi thứ nhanh hơn một chút, đặc biệt là trong trường hợp đa lõi:

> system.time(invisible(apply(M, 2, function(c) t.test(c)$statistic))) 
    user system elapsed 
80.920 0.437 82.109 
> system.time(invisible(mclapply(1:dim(M)[2], function(i) t.test(M[,i])$statistic))) 
    user system elapsed 
21.246 1.367 24.107 

Hoặc thậm chí nhanh hơn, tính giá trị t trực tiếp

my.t.test <- function(c){ 
    n <- sqrt(length(c)) 
    mean(c)*n/sd(c) 
} 

Sau đó

> system.time(invisible(apply(M, 2, function(c) my.t.test(c)))) 
    user system elapsed 
21.371 0.247 21.532 
> system.time(invisible(mclapply(1:dim(M)[2], function(i) my.t.test(M[,i])))) 
    user system elapsed 
144.161 8.658 6.313 

Answer 2

Bạn có thể làm tốt hơn chức năng colttests từ gói genefilter (trên Bioconductor).

> library(genefilter) 
> M <- matrix(rnorm(40),nrow=20) 
> my.t.test <- function(c){ 
+ n <- sqrt(length(c)) 
+ mean(c)*n/sd(c) 
+ } 
> x1 <- apply(M, 2, function(c) my.t.test(c)) 
> x2 <- colttests(M, gl(1, nrow(M)))[,"statistic"] 
> all.equal(x1, x2) 
[1] TRUE 
> M <- matrix(rnorm(1e7), nrow=20) 
> system.time(invisible(apply(M, 2, function(c) my.t.test(c)))) 
    user system elapsed 
27.386 0.004 27.445 
> system.time(invisible(colttests(M, gl(1, nrow(M)))[,"statistic"])) 
    user system elapsed 
    0.412 0.000 0.414 

Ref: "Đang tính toán hàng ngàn số liệu thống kê kiểm tra đồng thời trong R", SCGN, Vol 18 (1), 2007, http://stat-computing.org/newsletter/issues/scgn-18-1.pdf.