logo

How to Find Derivatives in Julia 📂Julia

How to Find Derivatives in Julia

Overview1

The package is named Calculus.jl, but it does not support integration.

If automatic differentiation, as discussed in machine learning, is needed, refer to the Zygote.jl package.

Differentiation of Single Variable Function

Derivative function derivative()

It calculates the derivative of $f : \R \to \R$.

  • derivative(f) or derivative(f, :x): Returns the derivative $f^{\prime}$.
  • derivative(f, a): Returns the differential coefficient $f^{\prime}(a)$.
julia> f(x) = 1 + 2x + 3x^2
f (generic function with 1 method)

julia> g(x) = sin(x)
g (generic function with 1 method)

julia> derivative(f)
#1 (generic function with 1 method)

julia> derivative(f, 1)
7.99999999996842

julia> Df = derivative(f)
#1 (generic function with 1 method)

julia> Dg = derivative(g)
#1 (generic function with 1 method)

#f'(x) = 2 + 6x
julia> Df(1)
7.99999999996842

#g'(x) = cos x
julia> Dg(pi)
-0.9999999999441258

Composite functions can also be differentiated.

#f∘g(x) = (2 + 6 sin x)cos x

julia> derivative(f∘g)
#1 (generic function with 1 method)

julia> derivative(f∘g, pi/4)
4.414213562300037

julia> (2+6sin(pi/4))cos(pi/4)
4.414213562373095

Second derivative function second_derivative()

It calculates the second derivative of $f : \R \to \R$.

The function returned by derivative() can take integers as input, but second_derivative() cannot use integer types. Irrational number types are also not allowed.

julia> derivative(f, 1)
7.99999999996842

julia> second_derivative(f, 1)
ERROR: MethodError: no method matching eps(::Type{Int64})
Closest candidates are:
  eps() at float.jl:764
  eps(::AbstractFloat) at float.jl:760
  eps(::Type{Float16}) at float.jl:761

julia> second_derivative(f, 1.)
5.9999956492003905

julia> second_derivative(g, pi)
ERROR: MethodError: no method matching eps(::Type{Irrational{:π}})
Closest candidates are:
  eps() at float.jl:764
  eps(::AbstractFloat) at float.jl:760
  eps(::Type{Float16}) at float.jl:761

julia> second_derivative(g, convert(Float64, pi))
-1.3553766145945872e-7

julia> second_derivative(g, 1pi)
-1.3553766145945872e-7

Differentiation of Multivariable Function

Gradient gradient()

Returns the gradient of $f : \mathbb{R}^{n} \to \mathbb{R}$.

Note that when defining a multivariable function, it should not be defined as a function with multiple real variables. It must be defined as a single-variable function that takes a vector as input. If it’s not a function accepting a vector as input, differentiation can be performed, but the value cannot be calculated. For example, it should not be defined as $f_{1}$, but rather as $f_{2}$.

julia> f₁(x,y,z) = x*y + z^2
f₁ (generic function with 1 method)

julia> Calculus.gradient(f₁)
#2 (generic function with 1 method)

julia> ∇f₁ = Calculus.gradient(f₁)
#2 (generic function with 1 method)

julia> ∇f₁(1,1,1)
ERROR: MethodError: no method matching (::Calculus.var"#2#4"{typeof(f₁), Symbol})(::Int64, ::Int64, ::Int64)

julia> ∇f₁([1,1,1])
ERROR: MethodError: no method matching f₁(::Vector{Float64})

julia> Calculus.gradient(f₁, 1,1,1)
ERROR: MethodError: no method matching gradient(::typeof(f₁), ::Int64, ::Int64, ::Int64)


julia> Calculus.gradient(f₁, [1,1,1])
ERROR: MethodError: no method matching f₁(::Vector{Float64})

julia> f₂(x) = x[1]*x[2] + x[3]^2
f₂ (generic function with 1 method)

julia> Calculus.gradient(f₂, [1,1,1])
3-element Vector{Float64}:
 1.0000000000235538
 1.0000000000235538
 1.9999999999737708

julia> ∇f₂ = Calculus.gradient(f₂)
#2 (generic function with 1 method)

julia> ∇f₂(1,1,1)
ERROR: MethodError: no method matching (::Calculus.var"#2#4"{typeof(f₂), Symbol})(::Int64, ::Int64, ::Int64)

julia> ∇f₂([1,1,1])
3-element Vector{Float64}:
 1.0000000000235538
 1.0000000000235538
 1.9999999999737708

Hessian hessian()

Returns the Hessian of $f : \mathbb{R}^{n} \to \mathbb{R}$.

Similar to second_derivative(), it only accepts Float data types as input. Like gradient(), it can only return values for functions that take vectors as input.

julia> hessian(f₂, [1.,1.,1.])
3×3 Matrix{Float64}:
 3.88008e-7  1.0         0.0
 1.0         3.88008e-7  0.0
 0.0         0.0         2.0

julia> H = hessian(f₂)
#7 (generic function with 1 method)

julia> H([1.,1.,1.])
3×3 Matrix{Float64}:
 3.88008e-7  1.0         0.0
 1.0         3.88008e-7  0.0
 0.0         0.0         2.0

Jacobian jacobian()

Returns the Jacobian of $f : \mathbb{R}^{m} \to \mathbb{R}^{n}$.

Similar to second_derivative(), it only accepts Float data types as input. Like gradient(), it can only return values for functions that take vectors as input.

Unlike other functions, using it like jacobian(f, [x, y, z]) is not possible.

julia> h(x) = [x[1], x[1]*x[2], x[1]*x[3]^2]
h (generic function with 2 methods)

julia> jacobian(h, [1.,1.,1.])
ERROR: MethodError: no method matching jacobian(::typeof(h), ::Vector{Float64})

julia> Jh = jacobian(h)
(::Calculus.var"#g#5"{typeof(h), Symbol}) (generic function with 1 method)

julia> Jh([1.,1.,1.])
3×3 Matrix{Float64}:
 1.0  0.0  0.0
 1.0  1.0  0.0
 1.0  0.0  2.0

Symbolic Differentiation

Symbolic differentiation is also available in the SymEngine.jl package.

differentiate()

Performs symbolic differentiation.

While it returns constant terms and $x$ neatly, for cases like $ax$ or $x^{n}$, it returns in the form of the product rule. For example, differentiating $3x^{2}$ would return it as the product of $3$ and $x^{2}$, seen as $\dfrac{d 3}{dx} x^{2} + 3\dfrac{d x^{2}}{dx}$. Even $x^{2}$ is seen as the product of $1$ and $x^{2}$.

julia> differentiate("1", :x)
0

julia> differentiate("1 + x", :x)
1

julia> differentiate("x^2", :x)
:(2 * 1 * x ^ (2 - 1))

julia> differentiate("x^3", :x)
:(3 * 1 * x ^ (3 - 1))

julia> differentiate("1 + x + x^2", :x)
:(1 + 2 * 1 * x ^ (2 - 1))

julia> differentiate("1 + 2x + 3x^2 + 4x^3", :x)
:((0 * x + 2 * 1) + (0 * x ^ 2 + 3 * (2 * 1 * x ^ (2 - 1))) + (0 * x ^ 3 + 4 * (3 * 1 * x ^ (3 - 1))))

julia> differentiate("x^2 * sin(x) + exp(x) * cos(x)", :x)
:(((2 * 1 * x ^ (2 - 1)) * sin(x) + x ^ 2 * (1 * cos(x))) + ((1 * exp(x)) * cos(x) + exp(x) * (1 * -(sin(x)))))

Characters not specified as symbols are treated as constants, and if more than one symbol is input, it returns the differentiation for each symbol. However, writing "3yx" treats xy as a single variable, so it’s important to explicitly denote multiplication by writing it as 3x*y.

julia> differentiate("1 + x + 3yx + y^2", :x)
:(1 + (0yx + 3 * 0))

julia> differentiate("1 + x + 3y*x + y^2", :x)
:(1 + ((0y + 3 * 0) * x + (3y) * 1))

julia> differentiate("1 + x + 3yx + y^2", [:x, :y])
2-element Vector{Any}:
 :(1 + (0yx + 3 * 0))
 :((0yx + 3 * 0) + 2 * 1 * y ^ (2 - 1))

julia> differentiate("1 + x + 3x*y + y^2", [:x, :y])
2-element Vector{Any}:
 :(1 + ((0 * x + 3 * 1) * y + (3x) * 0))
 :(((0 * x + 3 * 0) * y + (3x) * 1) + 2 * 1 * y ^ (2 - 1))

simplify()

The return of differentiate() lacks readability, but simplify() neatly organizes it. However, it doesn’t perform properly if vectors are used as input.

julia> simplify(differentiate("1 + 2x + 3x^2 + 4x^3", :x))
:(2 + 3 * (2x) + 4 * (3 * x ^ 2))

julia> simplify(differentiate("1 + x + 3x*y + y^2", :x))
:(1 + 3y)

julia> simplify(differentiate("1 + x + 3x*y + y^2", :y))
:(3x + 2y)

julia> simplify(differentiate("1 + x + 3x*y + y^2", [:x, :y]))
2-element Vector{Any}:
 :(1 + ((0 * x + 3 * 1) * y + (3x) * 0))
 :(((0 * x + 3 * 0) * y + (3x) * 1) + 2 * 1 * y ^ (2 - 1))

deparse()

Turns the return of symbolic differentiation into a string.

julia> a = differentiate("1 + x + 3x*y + y^2", :x)
:(1 + ((0 * x + 3 * 1) * y + (3x) * 0))

julia> deparse(a)
"1 + ((0 * x + 3 * 1) * y + (3 * x) * 0)"

julia> deparse(simplify(a))
"1 + 3 * y"

Verification

check_derivative()

One can check how much the derivative obtained by derivative() differs from the actual derivative. It’s implemented for the four types of derivatives excluding jacobian().

  • check_derivative(f, Df, a): Returns the absolute value of derivative(f, a)-Df(a).
julia> f(x) = 1 + x^2
f (generic function with 1 method)

julia> Df(x) = 2x
Df (generic function with 1 method)

julia> Calculus.check_derivative(f, Df, 1)
2.6229241001374248e-11

Application

Polynomials.jl

Although derivative is implemented in Polynomials.jl itself, the derivatives can also be calculated with Calculus.derivative().

julia> p = Polynomial([1,2,4,1])
Polynomial(1 + 2*x + 4*x^2 + x^3)

julia> Polynomials.derivative(p)
Polynomial(2 + 8*x + 3*x^2)

julia> Calculus.derivative(p)
#1 (generic function with 1 method)

julia> Polynomials.derivative(p,2)
Polynomial(8 + 6*x)

julia> Calculus.second_derivative(p)
#6 (generic function with 1 method)

Environment

  • OS: Windows10
  • Version: Julia 1.6.2, Calculus 0.5.1