The first use of the NFXP model was for a bus engine replacement model by Rust. What follows below is a simplified version of Rust's (@Rust1987) model using simulated data.
In the bus engine replacement model, a repairman is trying to minimize their discounted maintenance costs (in an infinite horizon setting). Each period they decide whether to either replace a bus engine or perform maintenance. Maintenance costs that are linear with the mileage acculumated engine. Bus engine mileage is discretized into 10 bins. Each time step the bus engine mileage increases 0-2 bins if they chose to not replace the bus engine. The bus engine replacement cost is unknown to the researcher but can be estimated by observing the repairman.
In [1]:
module Bus
import Distributions
import Optim
import NLsolve
using DynamicDiscreteChoice.General
using DynamicDiscreteChoice.NestedFixedPoint
const D = Distributions
const N = NestedFixedPoint
function call(utility::Utility{Val{:bus}}, state, action)
# maintenance cost normalized to one
mileage = state[1]
replacement_decision = action[1]
replacement_cost = utility.θ[1]
if replacement_decision==:keep
return -mileage
else
return -replacement_cost
end
end
immutable BusTransitionFamily
statespace
actionspace
θ_transition
end
const NullableV = Nullable{Vector{Float64}}
function Base.call(bt::BusTransitionFamily, θ_maybe::NullableV=NullableV())
θ_transition = get(θ_maybe, bt.θ_transition)
push!(θ_transition, 1-sum(θ_transition))
n = length(bt.statespace)
m = length(bt.actionspace)
transit = Array(Vector{Int64}, n, m)
replace = [1, min(2,n), min(3,n)]
for i=1:n
transit[i,1] = [i, min(i+1,n), min(i+2,n)]
transit[i,2] = replace
end
categorical = D.Categorical(θ_transition)
spt = SparseCategoricalTransition(transit, categorical, statespace, actionspace)
return spt
end
function N.stage1_solve(method::NFXPSolver{:bus}, problem)
f(x) = problem(x[1:3]) + x[4]*(sum(x[1:3]) - 1.0)
g = NLsolve.not_in_place(x -> Optim.gradient(f, x))
NLsolve.nlsolve(g, [1/3; 1/3; 1/3; 1.0], xtol=1e-6)
end
N.stage1_n_params(method::NFXPSolver{:bus}) = 2
function N.stage2_solve(method::NFXPSolver{:bus}, problem)
Optim.optimize(x -> problem([x]), 0.0, 10.0)
end
statespace = DiscreteVariableSpace(
:busstate,
[:mileage],
(collect(0.0:9.0),)
)
actionspace = DiscreteVariableSpace(
:busaction,
[:replacement_decision],
([:keep, :replace],)
)
transitionfamily = BusTransitionFamily(statespace, actionspace, [.4, .2])
solver = NFXPSolver{:bus}()
const β = 0.95
const θ_utility = [6.0]
const θ_transition = [0.5, 0.3]
model_family = NFXPFamily(
statespace,
actionspace,
transitionfamily,
UtilityFamily{Val{:bus}}(),
β)
# Create a simulated dataset using the θ_utility and θ_transition parameters
mock_result = mockup(model_family, θ_transition, θ_utility)
srand(100)
data = sample(mock_result, 1, episode_length=25, n_agents=100)
# Fit the model on the simulated data
result = fit(model_family, data, solver)
end
Out[1]:
In [2]:
Bus.result
Out[2]:
In [3]:
vcov = Bus.vcov(Bus.result, Bus.OPG())
stderr = diag(vcov).^0.5
Out[3]:
Since the true value of the replacement cost 6 and the estimated replacement cost was close (5.91) to that the estimator appears to be working. The bus mileage transition probabilities are also close to their true value. The OPG variance covariance matrix shows that the replacement cost parameter is significantly greater than zero.