search

bwlewis / GLM

Notes on generalized linear models
http://bwlewis.github.io/GLM
111 stars 16 forks source link

Cholesky pivoting #3

Open SWotherspoon opened 6 years ago

SWotherspoon commented 6 years ago

In the sparse and incremental variants, the resulting linear system is solved by Cholesky decomposition with pivoting. So $X^T X \beta = X^T y$ is solved essentially with

U <- chol(crossprod(X),pivot=TRUE)
p <- attr(U,"pivot")
backsolve(U,backsolve(U,crossprod(X,y)[p],transpose=TRUE))[p]

But shouldn't the final line be

backsolve(U,backsolve(U,crossprod(X,y)[p],transpose=TRUE))[order(p)]

to invert the original permutation?

bwlewis commented 6 years ago

Yes, indeed, the pivoting has this bug sorry. I am in th emiddle of substantially revising these notes and will include a fix along these in the revised version, which should appear sometime in July 2018. Thanks!

(sorry for the extreme latency)