# What is the most efficient way of obtaining the orthogonal column space of a matrix?
**URL:** https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212
**Category:** General Usage
**Tags:** linearalgebra
**Created:** [November 28, 2021, 9:44pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212 "2021-11-28T21:44:21Z")
**Posts on this page:** 15
**Page:** 1
### Author:  [@Joris\_Pinkse](https://discourse.julialang.org/u/Joris_Pinkse)
#### Post date: [November 28, 2021, 9:44pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/1 "2021-11-28T21:44:21Z")
What is the most efficient way of obtaining the orthogonal column space of a matrix? (many more rows than columns)
Q, = svd( X; alg = LinearAlgebra.QRIteration() )
gives me what I want, but is there something more efficient?
Relatedly, why isn’t there a command similar to nullspace? (or is there?)
---
### Author:  [@albheim](https://discourse.julialang.org/u/albheim)
#### Post date: [November 28, 2021, 10:24pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/2 "2021-11-28T22:24:33Z")
Maybe this?
```julia
julia> using LinearAlgebra
help?> nullspace
search: nullspace
nullspace(M; atol::Real=0, rtol::Real=atol>0 ? 0 : n*ϵ)
nullspace(M, rtol::Real) = nullspace(M; rtol=rtol) # to be deprecated in Julia 2.0
Computes a basis for the nullspace of M by including the singular vectors
of M whose singular values have magnitudes greater than max(atol,
rtol*σ₁), where σ₁ is M's largest singular value.
By default, the relative tolerance rtol is n*ϵ, where n is the size of
the smallest dimension of M, and ϵ is the eps of the element type of M.
Examples
≡≡≡≡≡≡≡≡≡≡
julia> M = [1 0 0; 0 1 0; 0 0 0]
3×3 Matrix{Int64}:
1 0 0
0 1 0
0 0 0
julia> nullspace(M)
3×1 Matrix{Float64}:
0.0
0.0
1.0
julia> nullspace(M, rtol=3)
3×3 Matrix{Float64}:
0.0 1.0 0.0
1.0 0.0 0.0
0.0 0.0 1.0
julia> nullspace(M, atol=0.95)
3×1 Matrix{Float64}:
0.0
0.0
1.0
```
---
### Author:  [@Joris\_Pinkse](https://discourse.julialang.org/u/Joris_Pinkse)
#### Post date: [November 28, 2021, 10:26pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/3 "2021-11-28T22:26:13Z")
Sorry for being unclear. I know about nullspace; I wanted something similar to that to get the column space.
---
### Author:  [@cvanaret](https://discourse.julialang.org/u/cvanaret)
#### Post date: [November 28, 2021, 10:41pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/4 "2021-11-28T22:41:41Z")
Typically, the QR decomposition uses the Gram-Schmidt process (it computes an orthonormal basis of a set of vectors) to compute Q. So calling Gram-Schmidt would be enough (but computing the R matrix in the QR decomposition really is cheap).
To compute the nullspace, you can apply the QR decomposition to the transposed matrix: [linear algebra - Null-space of a rectangular dense matrix - Computational Science Stack Exchange](https://scicomp.stackexchange.com/a/2511)
---
### Author:  [@albheim](https://discourse.julialang.org/u/albheim)
#### Post date: [November 28, 2021, 11:00pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/5 "2021-11-28T23:00:08Z")
I should also have understood that had I read you post properly 😅
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [November 29, 2021, 1:26am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/6 "2021-11-29T01:26:30Z")
`qr(X).Q`
---
### Author:  [@John\_Gibson](https://discourse.julialang.org/u/John_Gibson)
#### Post date: [November 29, 2021, 1:33am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/7 "2021-11-29T01:33:14Z")
Probably the reduced QR decomposition is what you want: A = QR, with m \times n matrix Q and n \times n matrix R for m \times n matrix A. The columns of Q span the columns of A. The Householder algorithm for computing the QR decomposition is better than Gram-Schmidt, in that it produces a Q that is unitary to machine precision and results in Ax=b solves that are backwards stable. Julia’s `qr` algorithm surely uses Householder.
---
### Author:  [@Joris\_Pinkse](https://discourse.julialang.org/u/Joris_Pinkse)
#### Post date: [November 29, 2021, 1:52am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/8 "2021-11-29T01:52:24Z")
Ah, qr(X).Q [:,1:size(X,2)] is an order of magnitude slower than the svd route, but using @view seems to do the trick, thanks.
---
### Author:  [@Joris\_Pinkse](https://discourse.julialang.org/u/Joris_Pinkse)
#### Post date: [November 29, 2021, 1:55am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/9 "2021-11-29T01:55:51Z")
It does. What I was missing there was that it also appeared to be computing the null space, but the reason qr(X).Q[:,1:size(X,2)] is so much slower is the subsequent indexing.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [November 29, 2021, 2:11am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/10 "2021-11-29T02:11:55Z")
QR factors are stored in an implicit form that makes indexing slow, but they can be quickly multiplied by vectors. To get the first n columns, try multiplying Q by the first n columns of the identity matrix.
---
### Author:  [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)
#### Post date: [November 29, 2021, 2:13am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/11 "2021-11-29T02:13:46Z")
I think I have a PR that fixes this, but it stalled out. I really should revive it though because it’s a dumb performance cliff.
---
### Author:  [@John\_Gibson](https://discourse.julialang.org/u/John_Gibson)
#### Post date: [November 29, 2021, 2:14am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/12 "2021-11-29T02:14:15Z")
Probably another reason for the slowness of `qr(X).Q` is that Householder QR decomp usually does not compute Q explicitly, but rather stores a series of Householder reflector vectors. When you ask for Q, it has to be constructed from the reflector vectors.
EDIT: what @stevengj said!
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [November 29, 2021, 2:21am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/13 "2021-11-29T02:21:13Z")
> [@stevengj](#):
>
> QR factors are stored in an implicit form that makes indexing slow, but they can be quickly multiplied by vectors. To get the first n columns, try multiplying Q by the first n columns of the identity matrix.
In particular, try `qr(X).Q * Matrix(I, size(X)...)`, and you should find that it is much faster than an SVD.
> [@Oscar\_Smith](#):
>
> I think I have a PR that fixes this, but it stalled out. I really should revive it though because it’s a dumb performance cliff.
Yes, it’s pretty common to want to slice a Q matrix and we really should provide a fast algorithm for it.
---
### Author:  [@Joris\_Pinkse](https://discourse.julialang.org/u/Joris_Pinkse)
#### Post date: [November 29, 2021, 2:41am UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/14 "2021-11-29T02:41:23Z")
Thanks, that’s very helpful and I agree.
Having a function like
```julia
colspace(X) = qr(X).Q * Matrix(I, size(X)...)
```
to complement `nullspace` would be welcome, also.
---
### Author:  [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [November 29, 2021, 12:36pm UTC](https://discourse.julialang.org/t/what-is-the-most-efficient-way-of-obtaining-the-orthogonal-column-space-of-a-matrix/72212/15 "2021-11-29T12:36:44Z")
You need a more complicated implementation (and probably the SVD, with tolerance parameters similar to `nullspace`) if you want to handle the case of rank-deficient matrices.