Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add inv and derivativeinv to some utility functions #225

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions src/QuantEcon.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
module QuantEcon

import Base: show, isapprox
import Base: show, isapprox, inv
import Statistics: mean, std, var
using Markdown
using FFTW
Expand Down Expand Up @@ -58,7 +58,7 @@ export

# modeltools
AbstractUtility, LogUtility, CRRAUtility, CFEUtility, EllipticalUtility,
derivative,
derivative, derivativeinv, inv,

# gth_solve
gth_solve,
Expand Down
25 changes: 20 additions & 5 deletions src/modeltools/utility.jl
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ abstract type AbstractUtility end
@doc doc"""
Type used to evaluate log utility. Log utility takes the form

u(c) = \log(c)
u(c) = ξ log(c)

Additionally, this code assumes that if c < 1e-10 then

Expand All @@ -26,8 +26,12 @@ LogUtility() = LogUtility(1.0)
c > 1e-10 ? u.ξ*log(c) : u.ξ*(log(1e-10) + 1e10*(c - 1e-10))
derivative(u::LogUtility, c::Float64) =
c > 1e-10 ? u.ξ / c : u.ξ*1e10
inv(u::LogUtility, x::Float64) =
exp.(x / u.ξ)
derivativeinv(u::LogUtility, x::Float64) = u.ξ / x

"""

@doc doc"""
Type used to evaluate CRRA utility. CRRA utility takes the form

u(c) = ξ c^(1 - γ) / (1 - γ)
Expand All @@ -53,17 +57,21 @@ end
c > 1e-10 ?
u.ξ * (c^(1.0 - u.γ) - 1.0) / (1.0 - u.γ) :
u.ξ * ((1e-10^(1.0 - u.γ) - 1.0) / (1.0 - u.γ) + 1e-10^(-u.γ)*(c - 1e-10))
inv(u::CRRAUtility, x::Float64) = (((1 - u.γ)*x) / u.ξ + 1)^(1 / (1 - u.γ))
derivative(u::CRRAUtility, c::Float64) =
c > 1e-10 ? u.ξ * c^(-u.γ) : u.ξ*1e-10^(-u.γ)
derivativeinv(u::CRRAUtility, x::Float64) = (u.ξ / x)^(1.0 / u.γ)


# Labor Utility

"""
@doc doc"""
Type used to evaluate constant Frisch elasticity (CFE) utility. CFE
utility takes the form

v(l) = ξ l^(1 + 1/ϕ) / (1 + 1/ϕ)
v(l) = -ξ l^(1 + 1/ϕ) / (1 + 1/ϕ)

where l is hours worked

Additionally, this code assumes that if l < 1e-10 then

Expand All @@ -88,12 +96,18 @@ end
-u.ξ * (1e-10^(1.0 + 1.0/u.ϕ)/(1.0 + 1.0/u.ϕ) + 1e-10^(1.0/u.ϕ) * (l - 1e-10))
derivative(u::CFEUtility, l::Float64) =
l > 1e-10 ? -u.ξ * l^(1.0/u.ϕ) : -u.ξ * 1e-10^(1.0/u.ϕ)
inv(u::CFEUtility, x::Float64) =
((1 + 1.0/u.ϕ)*x / -u.ξ)^(1.0 / (1.0 + 1.0/u.ϕ))
derivativeinv(u::CFEUtility, x) =
(x / -u.ξ)^u.ϕ


"""
@doc doc"""
Type used to evaluate elliptical utility function. Elliptical utility takes form

v(l) = b (1 - l^μ)^(1 / μ)

where l is hours worked
"""
struct EllipticalUtility <: AbstractUtility
b::Float64
Expand All @@ -107,3 +121,4 @@ EllipticalUtility(;b=0.5223, μ=2.2926) = EllipticalUtility(b, μ)
u.b * (1.0 - l^u.μ)^(1.0 / u.μ)
derivative(u::EllipticalUtility, l::Float64) =
-u.b * (1.0 - l^u.μ)^(1.0/u.μ - 1.0) * l^(u.μ - 1.0)

56 changes: 47 additions & 9 deletions test/test_modeltool.jl
Original file line number Diff line number Diff line change
@@ -1,53 +1,92 @@
@testset "Testing modeltools.jl" begin

@testset "Log Utility" begin
# Test log utility
# Test log utility
u = LogUtility(2.0)

# Test levels
@test isapprox(u(1.0), 0.0)
@test isapprox(2.0*log(2.5), u(2.5))
@test isapprox(u.(ones(5)), zeros(5))
@test u(2.0) > u(1.0) # Ensure it is increasing
@test isapprox(derivative(u, 1.0), 2.0)
@test isapprox(derivative(u, 2.0), 1.0)

# Test "extrapolation evaluations"
@test isapprox(u(0.5e-10), u(1e-10) + derivative(u, 1e-10)*(0.5e-10 - 1e-10))
@test u(-0.5) < u(-0.1) # Make sure it doesn't fail with negative values

@test isapprox(LogUtility().ξ, 1.0) # Test default constructor
# Test derivatives
@test isapprox(derivative(u, 1.0), 2.0)
@test isapprox(derivative(u, 2.0), 1.0)
@test isapprox(derivative(u, 1.0), (u(1.0 + 1e-4) - u(1.0))/1e-4; atol=1e-3)

# Test extrapolation derivatives
@test isapprox(derivative(u, 0.0), u.ξ*1e10)

# Test inverse of level
@test isapprox(u(inv(u, 1.0)), 1.0)
@test isapprox(u(inv(u, 2.0)), 2.0)

# Test inverse of derivative
@test isapprox(derivative(u, derivativeinv(u, 1.0)), 1.0)
@test isapprox(derivative(u, derivativeinv(u, 2.0)), 2.0)

@test isapprox(LogUtility().ξ, 1.0) # Test default constructor
end

@testset "CRRA Utility" begin
# Test CRRA utility
u = CRRAUtility(2.0)

# Test levels
@test isapprox(u(1.0), 0.0)
@test isapprox((2.5^(-1.0) - 1.0) / (-1.0), u(2.5))
@test isapprox(u.(ones(5)), zeros(5))
@test u(5.0) > u(3.0) # Ensure it is increasing
@test isapprox(derivative(u, 1.0), 1.0)
@test isapprox(derivative(u, 2.0), 0.25)

# Test "extrapolation evaluations"
@test isapprox(u(0.5e-10), u(1e-10) + derivative(u, 1e-10)*(0.5e-10 - 1e-10))
@test u(-0.5) < u(-0.1) # Make sure it doesn't fail with negative values

# Test derivatives
@test isapprox(derivative(u, 1.0), 1.0)
@test isapprox(derivative(u, 2.0), 0.25)

# Test inverse of level
@test isapprox(u(inv(u, 1.0)), 1.0)
@test isapprox(u(inv(u, 0.5)), 0.5)

# Test inverse of derivative
@test isapprox(derivative(u, derivativeinv(u, 1.0)), 1.0)
@test isapprox(derivative(u, derivativeinv(u, 2.0)), 2.0)

@test_throws ErrorException CRRAUtility(1.0) # Test error throwing at γ=1.0
end

@testset "CFE Utility" begin
# Test CFE Utility
v = CFEUtility(2.0)

# Test levels
@test isapprox(v(1.0), -2.0/3.0)
@test isapprox(v(0.5), -v.ξ * 0.5^(1.0 + 1.0/v.ϕ) / (1.0 + 1.0/v.ϕ))
@test v(0.5) > v(0.85) # Ensure it is decreasing
@test isapprox(derivative(v, 0.25), -0.5)
@test isapprox(derivative(v, 1.0), -1.0)

# Test "extrapolation evaluations"
@test isapprox(v(0.5e-10), v(1e-10) + derivative(v, 1e-10)*(0.5e-10 - 1e-10))
@test v(-0.5) > v(-0.1) # Make sure it doesn't fail with negative values

# Test derivatives
@test isapprox(derivative(v, 0.25), -0.5)
@test isapprox(derivative(v, 1.0), -1.0)

# Test inverse of level
@test isapprox(v(inv(v, -0.05)), -0.05)
@test isapprox(v(inv(v, -0.25)), -0.25)

# Test inverse of derivative
@test isapprox(derivative(v, derivativeinv(v, -0.5)), -0.5)
@test isapprox(derivative(v, derivativeinv(v, -0.25)), -0.25)

@test_throws ErrorException CRRAUtility(1.0) # Test error throwing at ϕ=1.0
end

Expand All @@ -63,7 +102,6 @@
# Test default values
@test isapprox(EllipticalUtility().b, 0.5223)
@test isapprox(EllipticalUtility().μ, 2.2926)

end

end
Expand Down