2017-05-26 14 views
3

を使用して、より速いGLM推定を達成しようとしていますが、それはなぜより遅いのですか?なぜ `speedglm`が` glm`より遅いのですか?

set.seed(0) 
n=1e3 
p=1e3 
x=matrix(runif(n*p),nrow=n) 
y=sample(0:1,n,replace = T) 

ptm <- proc.time() 
fit=glm(y~x,family=binomial()) 
print(proc.time() - ptm) 
# user system elapsed 
# 10.71 0.07 10.78 

library(speedglm) 
ptm <- proc.time() 
fit=speedglm(y~x,family=binomial()) 
print(proc.time() - ptm) 
# user system elapsed 
# 15.11 0.12 15.25 
+0

@李哲源ZheyuanLiコメントのおかげで。どうして?私は$ O(p^2)$ – hxd1011

+1

のような文書を読んで、あなたが働いている文脈を説明すべきです。できるだけ速く(妥当な可能性のある)小さなGLM適合をしたいなら、 'glm.fit()'を直接使うことを考えてみてください。 –

+0

@BenBolker私は3百万行のロジスティック回帰を試しています。 〜1000列、それが異なるパッケージでどのくらい速く実行されるのかを見たい。 – hxd1011

答えて

7

speedglmglm上の効率は、それがp * pマトリックスにn * pモデル行列を減少させる方法です。ただし、n = pがある場合、効果的な削減はありません。本当に確認したいのはn >> pケースです。


フィッシャースコアリングの各繰り返しで、より多くの洞察力の計算の複雑さ。 n * p行列に対してQR分解を用い

glmspeedglmp * p行列のQR分解に続いn * p行列の行列クロス積を形成しながら、np^2 + (4/3)p^3 FLOPを含む、2np^2 - (2/3)p^3 FLOPをとります。したがって、n >> pのように、の計算量は、glmの半分に過ぎません。さらに、speedglmによって使用されるブロッキング、キャッシング戦略は、コンピュータハードウェアをより良く使用し、高性能を提供します。

あなたはn = pを持っている場合は、すぐにglm(4/3)p^3 FLOPかかりますが、p^3 + (4/3)p^3 FLOPをspeedglmとり、より高価ことを参照してください!実際、この場合、行列の積は剪断オーバーヘッドになります!

+1

詳しいFLOPのための偉大な答え! – hxd1011

+0

もう1つ質問できますか?私が1e8行と1e4列を正則化ロジスティック回帰と呼んでいたら、どのパッケージを使うべきですか? – hxd1011

関連する問題