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.

Sign up for GitHub

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

Jump to bottom

fix default value and change variable name of discrete_approximation #277

Merged
merged 1 commit into from
Jul 23, 2022
Merged

fix default value and change variable name of discrete_approximation #277

Changes from all commits
Commits
Show all changes
1 commit
Select commit
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
38 changes: 19 additions & 19 deletions src/markov/markov_approx.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -687,7 +687,7 @@ Compute a discrete state approximation to a distribution with known moments,
using the maximum entropy procedure proposed in Tanaka and Toda (2013)

```julia
p, lambda_bar, moment_error = discrete_approximation(D, T, TBar, q, lambda0)
p, lambda_bar, moment_error = discrete_approximation(D, T, Tbar, q, lambda0)
```

##### Arguments
Expand All @@ -699,7 +699,7 @@ p, lambda_bar, moment_error = discrete_approximation(D, T, TBar, q, lambda0)
and returns an `L x N` matrix of moments evaluated
at each grid point, where L is the number of moments to be
matched.
- `TBar::AbstractVector` : length `L` vector of moments of the underlying distribution
- `Tbar::AbstractVector` : length `L` vector of moments of the underlying distribution
which should be matched

##### Optional
Expand All @@ -718,8 +718,8 @@ p, lambda_bar, moment_error = discrete_approximation(D, T, TBar, q, lambda0)
discretization minus actual moments) of length `L`

"""
function discrete_approximation(D::AbstractVector, T::Function, TBar::AbstractVector,
q::AbstractVector=ones(N)/N, # Default prior weights
function discrete_approximation(D::AbstractVector, T::Function, Tbar::AbstractVector,
q::AbstractVector=ones(length(D))/length(D), # Default prior weights
lambda0::AbstractVector=zeros(Tbar))

# Input error checking
Expand All @@ -728,15 +728,15 @@ function discrete_approximation(D::AbstractVector, T::Function, TBar::AbstractVe
Tx = T(D)
L, N2 = size(Tx)

if N2 != N || length(TBar) != L || length(lambda0) != L || length(q) != N2
if N2 != N || length(Tbar) != L || length(lambda0) != L || length(q) != N2
error("Dimension mismatch")
end

# Compute maximum entropy discrete distribution
options = Optim.Options(f_tol=1e-16, x_tol=1e-16)
obj(lambda) = entropy_obj(lambda, Tx, TBar, q)
grad!(grad, lambda) = entropy_grad!(grad, lambda, Tx, TBar, q)
hess!(hess, lambda) = entropy_hess!(hess, lambda, Tx, TBar, q)
obj(lambda) = entropy_obj(lambda, Tx, Tbar, q)
grad!(grad, lambda) = entropy_grad!(grad, lambda, Tx, Tbar, q)
hess!(hess, lambda) = entropy_hess!(hess, lambda, Tx, Tbar, q)
res = Optim.optimize(obj, grad!, hess!, lambda0, Optim.Newton(), options)
# Sometimes the algorithm fails to converge if the initial guess is too far
# away from the truth. If this occurs, the program tries an initial guess
Expand All @@ -751,16 +751,16 @@ function discrete_approximation(D::AbstractVector, T::Function, TBar::AbstractVe
lambda_bar = Optim.minimizer(res)
minimum_value = Optim.minimum(res)

Tdiff = Tx .- TBar
Tdiff = Tx .- Tbar
p = (q'.*exp.(lambda_bar'*Tdiff))/minimum_value
grad = similar(lambda0)
grad!(grad, Optim.minimizer(res))
moment_error = grad/minimum_value
return p, lambda_bar, moment_error
end
# discrete_approximation(D::AbstractVector, T::Function, TBar::Real,
# discrete_approximation(D::AbstractVector, T::Function, Tbar::Real,
# q::AbstractVector=ones(N)/N, lambda0::Real=0) =
# discrete_approximation(D, T, [TBar], q, [lambda0])
# discrete_approximation(D, T, [Tbar], q, [lambda0])

"""

Expand Down Expand Up @@ -796,15 +796,15 @@ end
Compute the maximum entropy objective function used in `discrete_approximation`

```julia
obj = entropy_obj(lambda, Tx, TBar, q)
obj = entropy_obj(lambda, Tx, Tbar, q)
```

##### Arguments

- `lambda::AbstractVector` : length `L` vector of values of the dual problem variables
- `Tx::AbstractMatrix` : `L x N` matrix of moments evaluated at the grid points
specified in discrete_approximation
- `TBar::AbstractVector` : length `L` vector of moments of the underlying distribution
- `Tbar::AbstractVector` : length `L` vector of moments of the underlying distribution
which should be matched
- `q::AbstractVector` : length `N` vector of prior weights for each point in the grid.

Expand All @@ -814,9 +814,9 @@ obj = entropy_obj(lambda, Tx, TBar, q)

"""
function entropy_obj(lambda::AbstractVector, Tx::AbstractMatrix,
TBar::AbstractVector, q::AbstractVector)
Tbar::AbstractVector, q::AbstractVector)
# Compute objective function
Tdiff = Tx .- TBar
Tdiff = Tx .- Tbar
temp = q' .* exp.(lambda'*Tdiff)
obj = sum(temp)
return obj
Expand All @@ -831,8 +831,8 @@ Compute gradient of objective function

"""
function entropy_grad!(grad::AbstractVector, lambda::AbstractVector,
Tx::AbstractMatrix, TBar::AbstractVector, q::AbstractVector)
Tdiff = Tx .- TBar
Tx::AbstractMatrix, Tbar::AbstractVector, q::AbstractVector)
Tdiff = Tx .- Tbar
temp = q'.*exp.(lambda'*Tdiff)
temp2 = temp.*Tdiff
grad .= vec(sum(temp2, dims = 2))
Expand All @@ -847,8 +847,8 @@ Compute hessian of objective function

"""
function entropy_hess!(hess::AbstractMatrix, lambda::AbstractVector,
Tx::AbstractMatrix, TBar::AbstractVector, q::AbstractVector)
Tdiff = Tx .- TBar
Tx::AbstractMatrix, Tbar::AbstractVector, q::AbstractVector)
Tdiff = Tx .- Tbar
temp = q'.*exp.(lambda'*Tdiff)
temp2 = temp.*Tdiff
hess .= temp2*Tdiff'
Expand Down