Gradient using ForwardDiff.jl not working

Hi there,

I wrote the following log-likelihood function for which I want to try out ForwardDiff.jl to compute the Jacobian. I have it also by hand but just wanted to try it out and check the differences. However I run into an error that I do not understand how to fix:

using ForwardDiff

struct Observation
    y::Float64
    x1::Float64
    x2::Float64
end

# this function will generate the logit data for Q4 
function simulate_logit_data(; n_obs = 300000)
    # define the DGP 
    X1 = 1.0:0.01:10.0
    X2 = 18:0.01:65.0
    α = -1.0
    θ_1 = 0.5
    θ_2 = 0.02
    # pre-allocate the container for the data 
    data = Array{Observation}(undef, n_obs)
    for i in eachindex(data)
        x1 = rand(X1)
        x2 = rand(X2)
        u = rand()
        y = (α + θ_1*x2 + θ_2*x2 + log(u/(1-u)) > 0) ? 1.0 : 0.0
        data[i] = Observation(y, x1, x2)
    end
    return data
end

function Xθ(θ, obs; f = fieldnames(Observation))
    result::Float64 = θ[1]
    for i in 2:length(f)
        result += θ[i]*getfield(obs, f[i])
    end
    return result
end

# this function will compute the log-likelihood given a parameter space point 
function LL(θ, data)
    ll::Float64 = 0.0
    for i in eachindex(data)
        ll += -log(1+exp(Xθ(θ, data[i]))) + getfield(data[i], :y)*Xθ(θ, data[i])
    end 
    return ll
end

function main()
    data = simulate_logit_data()
    ForwardDiff.gradient(θ -> LL(θ, data), rand(3))
end 

main()

Where the error message reads:

ERROR: LoadError: MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{var"#121#122"{Vector{Observation}}, Float64}, Float64, 3})
The type `Float64` exists, but no method is defined for this combination of argument types when trying to construct it.

Closest candidates are:
  (::Type{T})(::Real, ::RoundingMode) where T<:AbstractFloat
   @ Base rounding.jl:265
  (::Type{T})(::T) where T<:Number
   @ Core boot.jl:900
  Float64(::IrrationalConstants.Logten)
   @ IrrationalConstants ~/.julia/packages/IrrationalConstants/vp5v4/src/macro.jl:112
  ...

Stacktrace:
  [1] convert(::Type{Float64}, x::ForwardDiff.Dual{ForwardDiff.Tag{var"#121#122"{Vector{Observation}}, Float64}, Float64, 3})
    @ Base ./number.jl:7
  [2] Xθ(θ::Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#121#122"{…}, Float64}, Float64, 3}}, obs::Observation; f::Tuple{Symbol, Symbol, Symbol})
    @ Main ~/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:30
  [3] Xθ(θ::Vector{ForwardDiff.Dual{ForwardDiff.Tag{var"#121#122"{Vector{Observation}}, Float64}, Float64, 3}}, obs::Observation)
    @ Main ~/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:29
  [4] LL
    @ ~/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:41 [inlined]
  [5] #121
    @ ~/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:48 [inlined]
  [6] vector_mode_dual_eval!
    @ ~/.julia/packages/ForwardDiff/UBbGT/src/apiutils.jl:24 [inlined]
  [7] vector_mode_gradient
    @ ~/.julia/packages/ForwardDiff/UBbGT/src/gradient.jl:91 [inlined]
  [8] gradient(f::var"#121#122"{…}, x::Vector{…}, cfg::ForwardDiff.GradientConfig{…}, ::Val{…})
    @ ForwardDiff ~/.julia/packages/ForwardDiff/UBbGT/src/gradient.jl:20
  [9] gradient(f::var"#121#122"{Vector{…}}, x::Vector{Float64}, cfg::ForwardDiff.GradientConfig{ForwardDiff.Tag{…}, Float64, 3, Vector{…}})
    @ ForwardDiff ~/.julia/packages/ForwardDiff/UBbGT/src/gradient.jl:17
 [10] gradient(f::var"#121#122"{Vector{Observation}}, x::Vector{Float64})
    @ ForwardDiff ~/.julia/packages/ForwardDiff/UBbGT/src/gradient.jl:17
 [11] main()
    @ Main ~/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:48
 [12] top-level scope
    @ ~/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:51
 [13] include(fname::String)
    @ Main ./sysimg.jl:38
 [14] top-level scope
    @ REPL[73]:1
in expression starting at /Users/miguelborrero/Documents/Psets/StructuralEconometrics/Pset4/Code/test.jl:51
Some type information was truncated. Use `show(err)` to see complete types.

Thanks in advance!

Hi @miguelborrero!
The reason for this error is explained in the limitations of ForwardDiff: for it to work, you must make sure that your functions can accept numbers of generic type, not just Float64. Your code in its current state has several forced conversions to Float64: they do not improve performance and they make autodiff in forward mode impossible.
Here is a version where I made a few changes to:

  • increase generality by adapting to arbitrary types
  • improve clarify and possibly performance by turning Observation into a StaticVector instead of having a loop over its field names (this only helps when Observation has a handful of fields)
using ForwardDiff
using LinearAlgebra
using StaticArrays

# adapt struct to arbitrary field type
struct Observation{T<:Real}
    y::T
    x1::T
    x2::T
end

# efficiently turn Observation into a vector of statically known size
to_vector(o::Observation) = SVector(o.y, o.x1, o.x2)

# this function will generate the logit data for Q4 
function simulate_logit_data(; n_obs=300000)
    # define the DGP 
    X1 = 1.0:0.01:10.0
    X2 = 18:0.01:65.0
    α = -1.0
    θ_1 = 0.5
    θ_2 = 0.02
    # pre-allocate the container for the data, fully specifying its type
    data = Vector{Observation{Float64}}(undef, n_obs)
    for i in eachindex(data)
        x1 = rand(X1)
        x2 = rand(X2)
        u = rand()
        y = (α + θ_1 * x2 + θ_2 * x2 + log(u / (1 - u)) > 0) ? 1.0 : 0.0
        data[i] = Observation(y, x1, x2)
    end
    return data
end

# this function will compute the log-likelihood given a parameter space point 
function LL(θ, data)
    ll = zero(eltype(θ))  # give ll the element type of theta, which may be a Dual number
    for i in eachindex(data)
        obs = data[i]
        Xθ = dot(θ, to_vector(obs))  # precompute this quantity as a dot product
        ll += -log(1 + exp(Xθ)) + obs.y * Xθ
    end
    return ll
end

data = simulate_logit_data()
ForwardDiff.gradient(Base.Fix2(LL, data), @SVector(rand(3)))

Related:

1 Like

Thanks so much! I really appreciate the details in your response. Not only it fixed my issue but I learned the lesson generally.

1 Like