# Efficient IRLS

**URL:** <https://discourse.julialang.org/t/efficient-irls/9936>\
**Category:** Statistics\
**Tags:** regression\
**Created:** [March 24, 2018, 5:01am UTC](https://discourse.julialang.org/t/efficient-irls/9936 "2018-03-24T05:01:05Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![Nosferican](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nosferican/32/9275_2.png) [@Nosferican](https://discourse.julialang.org/u/Nosferican)\
**Post date:** [March 24, 2018, 5:01am UTC](https://discourse.julialang.org/t/efficient-irls/9936/1 "2018-03-24T05:01:05Z")

</div>

I was inspecting [GLM.jl](https://github.com/JuliaStats/GLM.jl) IRLS implementation, but haven’t figured out from the code how many matrix decompositions it performs.

For any given, `A::AbstractMatrix{<:Real}`, `b::AbstractVector{<:Real}`, and `wts::AbstractWeights `,  
for solving the linear system: Ax = b with weights w, x = (A^{\top}WA)^{-1}A^{\top}Wb.  
Assuming I compute the `F::QRCompactWY = qrfact(sqrt.(w) .* A)` and solve `x = F \ (sqrt.(w) .* y)`, for new weights w\_{2}, how can I solve the new system relying on my already computed `F::QRCompactWY`? Do I need to compute a new factorization every iteration? This [study](http://repositori.uji.es/xmlui/bitstream/handle/10234/164827/Belloch_Jos_preprint.pdf?sequence=1) suggests the best would be a Q-less decomposition.

---

<div class="post-metadata">

**Author:** ![mohamed82008](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mohamed82008/32/18171_2.png) [@mohamed82008](https://discourse.julialang.org/u/mohamed82008)\
**Post date:** [March 24, 2018, 9:36am UTC](https://discourse.julialang.org/t/efficient-irls/9936/2 "2018-03-24T09:36:13Z")

</div>

I am not sure if there is a way to re-use the QR factorization since stretching the coordinates of orthogonal vectors by different amounts distorts the orthogonality. But how big is A^TWA? If it is small, you can just factorize/invert it, possibly using [StaticArrays.jl](https://github.com/JuliaArrays/StaticArrays.jl) if it is of fixed size. If it is large and the weights are not aggressively changing then consider [IterativeSolvers.jl](https://github.com/JuliaMath/IterativeSolvers.jl) from the second iteration onwards since the previous solution would be a good starting point for the next one. One advantage of iterative solvers is that you don’t have to multiply out the matrices, you can just define a linear operator that applies the matrices one by one to a vector at the cost of some loss of accuracy.

---

<div class="post-metadata">

**Author:** ![Nosferican](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nosferican/32/9275_2.png) [@Nosferican](https://discourse.julialang.org/u/Nosferican)\
**Post date:** [April 4, 2018, 10:58pm UTC](https://discourse.julialang.org/t/efficient-irls/9936/3 "2018-04-04T22:58:27Z")

</div>

I found this cool [post](https://bwlewis.github.io/GLM/) with a few variants that achieves efficient IRLS.
