Eigen() warntype

calling LinearAlgebra.eigen() gives warntype. Is it a bug?

import Statistics
import LinearAlgebra
function test(x::Matrix{Float64})
	sigma = Statistics.cov(x)
	eig = LinearAlgebra.eigen(sigma)
	return eig

julia> @code_warntype test(randn(20, 5) )
Body::Union{Eigen{Complex{Float64},Complex{Float64},Array{Complex{Float64},2},Array{Complex{Float64},1}}, Eigen{Float64,Float64,Array{Float64,2},Array{Float64,1}}}
1 ─ %1 = Statistics.cov::Core.Compiler.Const(Statistics.cov, false)
β”‚   %2 = invoke Statistics.:(#cov#37)(1::Int64, true::Bool, %1::Function, _2::Array{Float64,2})::Array{Float64,2}
β”‚   %3 = LinearAlgebra.eigen::Core.Compiler.Const(LinearAlgebra.eigen, false)
β”‚   %4 = invoke LinearAlgebra.:(#eigen#58)(true::Bool, true::Bool, %3::Function, %2::Array{Float64,2})::Union{Eigen{Complex{Float64},Complex{Float64},Array{Complex{Float64},2},Array{Complex{Float64},1}}, Eigen{Float64,Float64,Array{Float64,2},Array{Float64,1}}}
└──      return %4

The eigenvalues and eigenvectors can be real or complex. This is deliberate behaviour (not a bug).

1 Like

okay, even it’s not a bug, it’s a design issue: it gives a warntype, which would slows down runtime.

or, is it a special case that has no runtime impact?

moreover, as far as I know, for positive-semifinite matrix like the covariance, the Eigen values and vectors would always be real. Could we do better?

The amount of time spent on the warning will be tiny compared to the time taken to calculate the eigen-decomposition.

There is no type instability if the matrix is e.g. of type Symmetric.

thanks! so, the function should better be written as:

function test2(x::Matrix{Float64})
	sigma = Statistics.cov(x)
	eig = LinearAlgebra.eigen(LinearAlgebra.Symmetric(sigma) )
	return eig

x = randn(20, 5)

julia> test(x) == test2(x)

@code_warntype test2(x)    # fine
1 Like

Maybe Statistics.cov should automatically return a Symmetric matrix?

1 Like

But probably that would lead to type instability somewhere down the line.

Note that you can test that the result type of the calculation is the same as the type that Julia infers using

using Test
@inferred test2(x)
1 Like

covariance is not only symmetric, it’s also positive-semi-finite (I think). It may be a big change to LinearAlgebra and Statistics to incorporate a special matrix type for it.

Yes, a covariance matrix is at least positive semi-definite.
If the number of linearly independent samples exceeds the dimensionality, it will be positive definite.

If you’re taking the eigendecomposition anyway, you can easily check whether all the eigenvalues are positive.