Issue arising from QR decomposition changes between 0.6 to 0.7

General Usage
1

I have an implementation of Algorithm 4 from this paper that works on julia 0.6, but is failing on julia 0.7. The crux of the issue can be demonstrated with these two lines of code:

U, T = qr(randn(10, 5))#analogous to line 3 in the paper's algorithm
T \ (transpose(U) * randn(10, 5))#analogous to line 4 in the paper's algorithm

Julia 0.7 complains about a dimension mismatch between T and transpose(U) * randn(10, 5) on the second line. This makes sense since size(T) == (5, 5) and size(transpose(U) * randn(10, 5)) == (10, 5). However, some other operations that ostensibly involve a dimension mismatch do work (e.g., U * T works despite the fact that size(U) == (10, 10) and size(T) == (5, 5)). I have seen this related post, but I haven’t been able to extract a solution from it.

Is there a way to reorganize the operations (or something) so that this will work on julia 0.7? It seems like it should be possible to implement Algorithm 4 from this paper in a few lines of julia code, but I am struggling to do that with julia 0.7 (even though I have such an implementation for julia 0.6). Anyone have ideas?

3

Thanks, but running this on 0.7 doesn’t produce the same result as my original code on 0.6. Is there some reason why we should expect this to work (aside from that the sizes of the objects line up)? Here’s what I came up with for 0.7 that seems to be equivalent to my original code on 0.6:

F = qr(randn(10, 5))
x = F.R \ (Matrix(F.Q)' * randn(10, 5))

This note from the documentation (which unfortunately isn’t in the ?qr documentation for julia 0.7) was the key to figuring this out:

A Q matrix can be converted into a regular matrix with Matrix . This operation returns the “thin” Q factor, i.e., if A is m × n with m>=n , then Matrix(F.Q) yields an m × n matrix with orthonormal columns.

I thought I’d make a note of it here in case someone finds it helpful in the future.

1 Like
4

Same experience here. And the same solution

fact = qr(A); 
p = UpperTriangular(fact.R)\(transpose(Array(fact.Q))*b) # R \ (transpose(Q) * b)
1 Like
5

In Julia 0.6, Q, R = qr(::Matrix) would yield two matrices. In Julia 0.7 and on F = qr(::Matrix) yields a <:Factorization. The Q is kept as QRCompactWYQ. If you want the Matrix version of it, you can do Matrix(F.Q). For additional information do check the documentation.

6

If want to compute R \ (Q'*b) — this is the least-square solution — just do F \ b where F is the factorization object returned by qr.

For example, the following all compute the same thing (up to roundoff errors):

julia> A = randn(10, 5); b = randn(10);

julia> A \ b
5-element Array{Float64,1}:
 -2.1381793415186894
 -1.925894533068401 
 -0.6539752131528996
  0.2931584134939902
 -0.4190898942664864

julia> F = qr(A); F \ b
5-element Array{Float64,1}:
 -2.1381793415186863 
 -1.925894533068399  
 -0.6539752131528995 
  0.29315841349399047
 -0.41908989426648574

julia> F.R \ (Matrix(F.Q)' * b)
5-element Array{Float64,1}:
 -2.1381793415186867 
 -1.9258945330683994 
 -0.6539752131528997 
  0.29315841349399074
 -0.41908989426648574
4 Likes
7

Thanks, Steven. This seems like the ideal solution, and makes me wonder “Why didn’t I think of that?”

8

ERROR: UndefVarError: qr not defined

9

You need julia> using LinearAlgebra. Please do not post on old threads.

2 Likes