Convert Optim objective value back to a Float

Specific Domains Optimization (Mathematical)
,
1

I am working with the code in this gist.

In the file test_autochoice_mew2.jl , I have the following pieces of code:

[...]
T = promote_type(promote_type(promote_type(eltype(X),eltype(b)),eltype(Z)),eltype(W))
[...]
ℓ     =  zero(T)
[...]
rs = optimize(td, bstart, LBFGS(; linesearch = LineSearches.BackTracking()), Optim.Options(iterations=100_000,g_tol=1e-8,f_tol=1e-8))
β  = Optim.minimizer(rs)
ℓ  = Optim.minimum(rs)*(-1)

When I try to print like (which is what I’ve labeled ℓ’s output) at the end of the function, I get the following output that means nothing to me:

julia> like
Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(255.55119809155235,-5.033796958153758e-7,-1.0804235078265556e-5,-3.4508112290165815e-6,2.989085825272353e-6,7.117737226280951e-6),Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(-5.033796958153758e-7,66.48776626068965,47.6826072832014,-41.057855033352624,-26.4964747155521,854.6887718388765),Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(-1.0804235078265556e-5,47.6826072832014,47.6826072832014,-26.4964747155521,-26.4964747155521,612.3806230745521),Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(-3.4508112290165815e-6,-41.057855033352624,-26.4964747155521,49.178844548215885,32.485863319462716,-433.9826291198754),Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(2.989085825272353e-6,-26.4964747155521,-26.4964747155521,32.485863319462716,32.485863319462716,-263.0385457700083),Dual{ForwardDiff.Tag{getfield(Main, Symbol("#f#3")){Array{Float64,1},Array{Float64,2},Array{Float64,3},Int64,Array{Float64,1},Int64,Int64},Float64}}(7.117737226280951e-6,854.6887718388763,612.380623074552,-433.9826291198754,-263.0385457700083,12296.248925836317))

I am confident this is coming from the line above where I do a bunch of type promotions. But I don’t know how to get it back to the Float64 that I know like should be.

Does anyone know how I can do that?

2

Are you differentiating through the optimization solver? If not, it sounds strange to have the objective value being a dual number.

2 Likes
3

Thank you for responding. Yes, I’m using ForwardDiff in the optimization:

td = TwiceDifferentiable(f, g!, bstart, autodiff = :forwarddiff)
4

Ah, but you’re using forward diff to differentiate the objective function then, not the optimizer itself. If this is the case, the result should not be a Dual number.

1 Like
5

OK, I’ll have to do some more figuring out what’s going on, then. Thanks for your responsiveness!

6

It’s a scope thing, try this (I only changed the name of the variable which overlap with the variable name of the objective value)


@views @inline function asclogit(bstart::Vector,Y::Array,X::Array,Z::Array,J::Int64,W::Array=ones(length(Y)))

    ## error checking
    @assert ((!isempty(X) || !isempty(Z)) && !isempty(Y))    "You must supply data to the model"
    @assert (ndims(Y)==1 && size(Y,2)==1)                    "Y must be a 1-D Array"
    @assert (minimum(Y)==1 && maximum(Y)==J) "Y should contain integers numbered consecutively from 1 through J"
    if !isempty(X)
        @assert ndims(X)==2          "X must be a 2-dimensional matrix"
        @assert size(X,1)==size(Y,1) "The 1st dimension of X should equal the number of observations in Y"
    end
    if !isempty(Z)
        @assert ndims(Z)==3          "Z must be a 3-dimensional tensor"
        @assert size(Z,1)==size(Y,1) "The 1st dimension of Z should equal the number of observations in Y"
    end

    K1 = size(X,2)
    K2 = size(Z,2)

    function f(b)
        T = promote_type(promote_type(promote_type(eltype(X),eltype(b)),eltype(Z)),eltype(W))
        num   = zeros(T,size(Y))
        dem   = zeros(T,size(Y))
        temp  = zeros(T,size(Y))
        numer =  ones(T,size(Y,1),J)
        P     = zeros(T,size(Y,1),J)
        ℓ     =  zero(T)
        b2 = b[K1*(J-1)+1:K1*(J-1)+K2]
                                             
        for j=1:(J-1)
            temp       .= X*b[(j-1)*K1+1:j*K1] .+ (Z[:,:,j].-Z[:,:,J])*b2
            num        .= (Y.==j).*temp.+num
            dem        .+= exp.(temp)
            numer[:,j] .=  exp.(temp)
        end
        dem.+=1
        P   .=numer./(1 .+ sum(numer;dims=2))

        ℓ = -W'*(num.-log.(dem))
    end

    function g!(G,b)
        T     = promote_type(promote_type(promote_type(eltype(X),eltype(b)),eltype(Z)),eltype(W))
        numer = zeros(T,size(Y,1),J)
        P     = zeros(T,size(Y,1),J)
        numg  = zeros(T,K2)
        demg  = zeros(T,K2)
        b2    = b[K1*(J-1)+1:K1*(J-1)+K2]
                                                                                                                                     
        #G = Array{T}(undef, length(b))
        G .= T(0)
        for j=1:(J-1)
            numer[:,j] .= exp.( X*b[(j-1)*K1+1:j*K1] .+ (Z[:,:,j].-Z[:,:,J])*b2 )
        end
        P   .=numer./(1 .+ sum(numer;dims=2))

        for j=1:(J-1)
            G[(j-1)*K1+1:j*K1] .= -X'*(W.*((Y.==j).-P[:,j]))
        end

        for j=1:(J-1)
            numg .+= -(Z[:,:,j]-Z[:,:,J])'*(W.*(Y.==j))
            demg .+= -(Z[:,:,j]-Z[:,:,J])'*(W.*P[:,j])
        end
        G[K1*(J-1)+1:K1*(J-1)+K2] .= numg.-demg
        G
    end

    td = TwiceDifferentiable(f, g!, bstart, autodiff = :forwarddiff)
    rs = optimize(td, bstart, LBFGS(; linesearch = LineSearches.BackTracking()), Optim.Options(iterations=100_000,g_tol=1e-8,f_tol=1e-8))
    β  = Optim.minimizer(rs)
    ℓℓ  = Optim.minimum(rs)*(-1)
    H  = Optim.hessian!(td, β)
    se = sqrt.(diag(inv(H)))
    q = Array{Float64}(undef, length(β))
    grad = g!(q,β)
    gradzero = g!(q,zeros(size(β)))
    return β,se,ℓℓ,grad,gradzero
end

β,se_β,like,grad,grad0 = asclogit(zeros(size(θ)),Y,X,Zmatlab,3)
7

Thanks! :man_facepalming: This stuff always gets me.