作者celestialgod (天)
看板R_Language
標題Re: [問題] 關於這段程式碼 要如何優化?
時間Sat Mar 12 01:52:25 2016
※ 引述《jackhzt (巴克球)》之銘言：
: [問題類型]:
: 效能諮詢(我想讓R 跑更快)
: [軟體熟悉度]:
: 使用者(已經有用R 做過不少作品)
: [問題敘述]:如何將以下的程式碼跑快一點
: [程式範例]:R中的dist這function 因為想使用不同的計算方式，所以希望可以做到
: 和此function表現差不多的function
: 程式碼可貼於以下網站：
: https://gist.github.com/anonymous/cf844933bb6858936e25
: 希望有提高效率的方法
1. 如果距離函數是常見的，通常建議用dist達成
2. 如果不常見，盡量考慮用矩陣運算求出來
EX: (以歐式距離來說)
library(magrittr)
x <- matrix(rnorm(50), 5, 10)
distMat <- sweep(-x %*% t(x) * 2, 2, rowSums(x^2), '+') %>%
sweep(1, rowSums(x^2), '+')
diag(distMat) <- 0
distMat %<>% sqrt
all.equal(distMat, as.matrix(dist(x)), check.attributes = FALSE) # TRUE
3. 用簡單的平行，我可能會這樣做：
library(magrittr)
library(foreach)
library(doSNOW)
library(plyr)
library(Matrix)
dis <- function(x, y) sum(abs(as.numeric(x)-as.numeric(y)))
x <- matrix(rnorm(50), 5, 10)
allCombinations <- combn(1:nrow(x), 2)
cl <- makeCluster(8, type = "SOCK")
registerDoSNOW(cl)
clusterExport(cl, list = c("dis", "x"))
res <- aaply(allCombinations, 2, function(v){
dis(x[v[1],], x[v[2],])
}, .parallel = TRUE)
stopCluster(cl)
distMat <- sparseMatrix(i = allCombinations[1,],
j = allCombinations[2,], x = res) %>%
rbind(0) %>% as.matrix
distMat[lower.tri(distMat)] <- res
4. 或是乾脆用RcppArmadillo:
(下面是用之前kernel matrix估計的方法，之前有發過更快的方法...)
library(Rcpp)
library(RcppArmadillo)
## For windows user
# library(inline)
# settings <- getPlugin("Rcpp")
# settings$env$PKG_CXXFLAGS <- paste('-fopenmp', settings$env$PKG_CXXFLAGS)
# settings$env$PKG_LIBS <- paste('-fopenmp -lgomp', settings$env$PKG_LIBS)
# do.call(Sys.setenv, settings$env)
sourceCpp(code = '
// [[Rcpp::depends(RcppArmadillo)]]
#include <RcppArmadillo.h>
#include <omp.h>
// [[Rcpp::plugins(openmp)]]
using namespace Rcpp;
using namespace arma;
// [[Rcpp::export]]
NumericMatrix kernelMatrix_cpp(NumericMatrix Xr, NumericMatrix Centerr,
double sigma) {
omp_set_num_threads(omp_get_max_threads());
uword n = Xr.nrow(), b = Centerr.nrow(), row_index, col_index;
mat X(Xr.begin(), n, Xr.ncol(), false);
mat Center(Centerr.begin(), b, Centerr.ncol(), false);
mat KerX(n, b);
#pragma omp parallel private(row_index, col_index)
for (row_index = 0; row_index < n; row_index++)
{
#pragma omp for nowait
for (col_index = 0; col_index < b; col_index++)
{
KerX(row_index, col_index) = exp(sum(square(X.row(row_index)
- Center.row(col_index))) / (-2.0 * sigma * sigma));
}
}
return wrap(KerX);
}')
--
R資料整理套件系列文：
magrittr #1LhSWhpH (R_Language) http://tinyurl.com/1LhSWhpH
data.table #1LhW7Tvj (R_Language) http://tinyurl.com/1LhW7Tvj
dplyr(上) #1LhpJCfB (R_Language) http://tinyurl.com/1LhpJCfB
dplyr(下) #1Lhw8b-s (R_Language)
tidyr #1Liqls1R (R_Language) http://tinyurl.com/1Liqls1R
--
※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 140.109.73.238
※ 文章網址: https://webptt.com/m.aspx?n=bbs/R_Language/M.1457718749.A.B23.html
※ 編輯: celestialgod (140.109.73.238), 03/12/2016 01:55:02
1F：推 jackhzt: 感謝大大 <303/12 02:03
2F：推 jackhzt: 問一下 <anonymous>: ... may be used in an incorrect03/12 03:27
3F：→ jackhzt: context: ‘.fun(piece, ...)出現這是正常的嗎?03/12 03:27
正常
4F：推 jackhzt: 再請教一下 在後面一定要加.parallel = TRUE才有平行運03/12 03:42
對
5F：→ jackhzt: 還有參考以前的文章 使用平行時 如果不使用parApply03/12 03:43
6F：→ jackhzt: 或是其他 parCapply等等的函數 依然會有加速的效果嗎?03/12 03:44
沒有
7F：→ jackhzt: 剛剛看了大大的教學 開始讀了些平行運算的東西03/12 03:45
8F：→ jackhzt: 還是蠻多不懂的地方03/12 03:46
可以回文繼續問
9F：推 jackhzt: 試做3000多條向量，使用snow好像依樣跑不太動...03/12 17:43
po一下你的計算公式，我測試看看
10F：推 jackhzt: 我用的也是測試用的 實際的公式還在推 目前嘗試的是用 03/12 18:24
11F：→ jackhzt: dis = function(x, y)sum(as.numeric(x)!=as.numeric(y)) 03/12 18:24
12F：→ jackhzt: x為 R17的向量 03/12 18:25
13F：→ jackhzt: 都是類別變數 1～4 03/12 18:26
請問什麼是R17？
14F：推 jackhzt: X1=[x11,x12,x13...x17] 03/12 19:34
15F：→ jackhzt: 向量數為 X1,X2,X3...X3000 03/12 19:35
這個可以向量化運算，我寫了一個比較
你可以測試看看雙迴圈跟向量化運算的差異
能向量化的方法絕對不要用for or apply：
http://pastebin.com/y8gC3ap7
16F：推 jackhzt: 謝謝大大 我認真研讀一下 感激涕零~ 03/12 20:03
不客氣，簡單平行還太慢就得想辦法(攤手
※ 編輯: celestialgod (180.218.152.118), 03/12/2016 20:11:27