※ 引述《ror (回血375)》之铭言：
: [问题叙述]:
: 请简略描述你所要做的事情，或是这个程式的目的
: 各位前辈好
: 想请教一个问题
: 目前我手上有大概几百万(快到千万)组资料 一组资料是两笔内容
: 现在要对每一组资料跑lm()回归后 取出相关系数与P-value 等等资讯
: 但因为要跑的量实在是太大量了 导致速度很慢 已经使用parallel
: 目前写法是 loop-> lm() ->summary-> 取出value
: 想请问是否有方式避免重复宣告lm()物件 ...或是其他可以提升速度的方式
: 感谢 Orz
以下可能跟原PO的问题关系小了一点XD
我提供一个利用RcppArmadillo and RcppParallel做
leave-one-out cross validation (LOOCV) 的例子
十万个样本做LOOCV:
程式：http://pastebin.com/4JK6VSd7
在我的电脑上跑一次花了173.79秒 (平均一次 0.17 milliseconds)
这个应该远比用R的循环快很多...
内存使用的部分也保持在300MB以下，基本上不太会用到太多内存做暂存
给一个简单的比较：
mses = vector('numeric', N)
st = proc.time()
for (i in 1:N){
coef_lm = coef(lm.fit(cbind(1, X[-i,]), y[-i]))
mses[i] = y[i] - c(1, X[i,]) %*% coef_lm
}
proc.time() - st
上面的程式在N = 5000大概花10秒，Rcpp的版本只需要0.3秒
题外话，Rcpp真的值得一学，如果想要直接使用R的BLAS可以考虑RcppArmadillo
如果觉得内建的BLAS太慢可以考用RcppEigen (他本身的BLAS也算是够快的)
不过最简单的方式就是用RRO，或是在linux/mac上也可以考虑openblas
最后再搭配RcppParallel就可以得到很大的加速了XDD
PS: 里面有一段程式是避免RcppParallel跟RRO使用的MKL用到过多的thread
而造成速度的拖累，所以写了一段去减少MKL使用的thread
在有超执行绪(HT)的电脑上，MKL会用2个，其余则是使用1个thread
PS2: RRO是Revolution R Open
最后，原PO的问题也并不复杂，相关系数跟p-value可以分别透过
RcppArmadillo:::cor跟Rcpp:::pt做处理，如果有任何问题，欢迎讨论XD
Rcpp:::pt用法：
http://stackoverflow.com/questions/20144528/how-use-correctly-rcpppt
补充，如果要用K-fold可以用下面的函数去生成cvIndex:
(刚好有现成的，顺便分享XD)
好读版：http://pastebin.com/RhJx3kLs
uvec cvfold_index_f(const uword n, const uword fold)
{
uword fold_n = n / fold, rem = n - fold_n * fold;
uvec tmp_vec = linspace<uvec>(0, fold - 1, fold),
index = vectorise(repmat(tmp_vec, fold_n, 1));
if (rem > 0)
index = join_cols(index, index.head(rem));
// random permutation
int j;
for (int i = 0; i < n; i++)
{
j = as_scalar(randi<uvec>(1, distr_param(i, n - 1)));
index.swap_rows(i, j);
}
return index;
}