Skip to content

Supporting higher derivatives of the log of the gamma function #196

Closed
@fredRos

Description

@fredRos

The highest derivative of lgamma is currently limited to the first derivative of trigamma(x) = polygamma(1,x), probably because it has a special name. But of course the derivative of polygamma(n, x) is polygamma(n+1,x) and they are all available. I ran into this problem when computing the Hessian of a function involving trigamma. Here is a minimal example

julia> harr(x::Vector) = 1/2 * log(x[1]*trigamma(x[1])-1) - log(x[2])
harr (generic function with 1 method)

julia> ForwardDiff.hessian(harr, [1.5, 60.])
ERROR: MethodError: Cannot `convert` an object of type ForwardDiff.Dual{2,Float64} to an object of type Float64
This may have arisen from a call to the constructor Float64(...),
since type constructors fall back to convert methods.
 in harr(::Array{ForwardDiff.Dual{2,ForwardDiff.Dual{2,Float64}},1}) at ./REPL[15]:1
 in vector_mode_gradient at /home/beaujean/.julia/v0.5/ForwardDiff/src/gradient.jl:60 [inlined]
 in gradient(::#harr, ::Array{ForwardDiff.Dual{2,Float64},1}, ::ForwardDiff.GradientConfig{2,ForwardDiff.Dual{2,Float64},Array{ForwardDiff.Dual{2,ForwardDiff.Dual{2,Float64}},1}}) at /home/beaujean/.julia/v0.5/ForwardDiff/src/gradient.jl:7
 in vector_mode_jacobian(::ForwardDiff.##33#34{#harr,ForwardDiff.HessianConfig{2,Float64,Array{ForwardDiff.Dual{2,Float64},1},ForwardDiff.Dual{2,Float64},Array{ForwardDiff.Dual{2,ForwardDiff.Dual{2,Float64}},1}}}, ::Array{Float64,1}, ::ForwardDiff.JacobianConfig{2,Float64,Array{ForwardDiff.Dual{2,Float64},1}}) at /home/beaujean/.julia/v0.5/ForwardDiff/src/jacobian.jl:75
 in jacobian(::ForwardDiff.##33#34{#harr,ForwardDiff.HessianConfig{2,Float64,Array{ForwardDiff.Dual{2,Float64},1},ForwardDiff.Dual{2,Float64},Array{ForwardDiff.Dual{2,ForwardDiff.Dual{2,Float64}},1}}}, ::Array{Float64,1}, ::ForwardDiff.JacobianConfig{2,Float64,Array{ForwardDiff.Dual{2,Float64},1}}) at /home/beaujean/.julia/v0.5/ForwardDiff/src/jacobian.jl:7
 in hessian(::#harr, ::Array{Float64,1}) at /home/beaujean/.julia/v0.5/ForwardDiff/src/hessian.jl:6

I checked the source code of Calculus and will post an issue there to add this functionality. I'm wondering why the generic support for polygamma hasn't been added. Perhaps it can't pick up the first argument to polygamma?

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions