# Optim interface.jl type instability and efficient exception handling

Hi,

I am currently trying to optimize a simulation that involves repeated parameter estimation, which I do with Optim.jl’s `optimize` with supplied gradients and hessians.

interjace.jl type instability

Looking at the flamegraph of one such simulation run, it seems that there are type instabilities in the internals of optimization. In particular, I get type warnings for:

• interface.jl #optimize#89: line 82
• interface.jl value_gradient!!: line 82
• interface.jl hessian!!: line 95

The flame graph for one such optimization run looks like this, where the intermitted red bars correspond to the type instabilities from the interface.jl calls.

and the code for this is

``````using LinearAlgebra, Optim, ProfileView, Profile, Random, Distributions

n = 200;
p = 11;
beta_0 = randn(11);
beta = randn(p);
X = randn(n,p);
X[:,1] .= 1.0;
y = convert.(Float64, rand.(Binomial.(1, cdf.(Logistic(), X * beta_0))));
beta_init = fill(0.0,11);
buff_X = similar(X);
mus = randn(n);
ws = similar(mus);

function get_logistic_MLE2!(beta::Vector{Float64},X::Matrix{Float64},y::Vector{Float64},beta_init::Vector{Float64},mus::Vector{Float64},ws::Vector{Float64},buff_X::Matrix{Float64}; ϵ::Float64 = 0.0 , kwargs...)
opt = Optim.optimize(β -> logistic_loglikl(β,X,y,mus,ws;ϵ),
(g,β)->logistic_score!(g,β,X,y,mus;ϵ),
(H,β)->logistic_hessian!(H,β,X,mus,ws,buff_X;ϵ),
beta_init;kwargs...)
beta .= opt.minimizer
return Optim.converged(opt) && !any(isnan.(beta)) , Optim.minimum(opt)
end

function logistic_loglikl(beta::Vector{Float64},X::Matrix{Float64},y::Vector{Float64},mus::Vector{Float64},buff_mus::Vector{Float64}; ϵ = 0.0)
mul!(mus,X,beta)
buff_mus .= y .* mus
mus .= 1.0 .+ exp.(mus)
mus .= log.(mus)
if !(ϵ == 0.0)
buff_mus .= (1.0 + ϵ/2) .* buff_mus
mus .= (1.0 + ϵ) .* mus
end
mus .= mus .- buff_mus
return sum(mus)
end

function logistic_score!(g::Vector{Float64},beta::Vector{Float64},X::Matrix{Float64},y::Vector{Float64},mus::Vector{Float64}; ϵ = 0.0)
logistic_mean!(mus,X,beta)
if ϵ == 0.0
mus .= mus .- y
else
mus .= (1.0 + ϵ) .* mus .- (1.0 + ϵ/2) .* y
end
mul!(g, X', mus)
return g
end

function logistic_hessian!(H::Matrix{Float64},beta::Vector{Float64},X::Matrix{Float64},mus::Vector{Float64},ws::Vector{Float64},buff_X::Matrix{Float64}; ϵ = 0.0)
logistic_mean!(mus,X,beta)
ws .= mus .* (1.0 .- mus)
colmult!(buff_X,ws,X)
mul!(H,X',buff_X)
if !(ϵ == 0.0)
H .= (1.0 + ϵ) .* H
end
return H
end

function logistic_score!(g::Vector{Float64},beta::Vector{Float64},X::Matrix{Float64},y::Vector{Float64},mus::Vector{Float64}; ϵ = 0.0)
logistic_mean!(mus,X,beta)
if ϵ == 0.0
mus .= mus .- y
else
mus .= (1.0 + ϵ) .* mus .- (1.0 + ϵ/2) .* y
end
mul!(g, X', mus)
return g
end

function logistic_mean!(mus::Vector{Float64},X::Matrix{Float64},beta::Vector{Float64})
mul!(mus,X,beta,-1.0,false)
mus .= exp.(mus)
mus .= 1.0 ./ (1.0 .+ mus)
end

function colmult!(A::Matrix{Float64},v::Vector{Float64}, B::Matrix{Float64})
@inbounds for j = 1:size(A,2)
@inbounds for i = 1:size(A,1)
A[i,j] = v[i] * B[i,j]
end
end
A
end

Profile.init(10000000, 0.00001)
ProfileView.@profview get_logistic_MLE2!(beta,X,y,beta_init,mus,ws,buff_X)
``````

Is there anything wrong with how I set up my functions that causes this type instability and/or is there anything I can do to get rid of it?

Exception handling

My simulation essentially involves repeatedly drawing some data, estimating some parameters and do some calculations, which involves solving some nonlinear systems of equations whcih I do using NLsolve.jl’s `nlsolve`. With that, there are some data configurations for which estimation will fail, and I will get some error, usually some `LAPACKException, SingularException, NLsolve.IsFiniteException,DomainError`
If estimation fails, all I want to do is set my parameter estimates and other things to NaN, and continue, I do not care about the type of exception.

At the moment I am wrapping all things that can fail in a try-catch block, which works, but appears to be quite costly. Is there a clever/efficient way to simply return an error without throwing it, so that I can get rid of my try-catch statements, e.g. by telling julia to expect either, say a `Vector{Float64}` or some exception as return type.

Thank you in advance for the help.