From 38f19cf24386214ed30631e6301b21ee4b691ad2 Mon Sep 17 00:00:00 2001 From: Alan Lujan Date: Thu, 4 Jan 2024 19:36:34 -0500 Subject: [PATCH] move agents to own file --- code/agents.py | 162 +++++++++++++++++++++++++++++++++++++++++++++ code/estimation.py | 153 +++++++----------------------------------- 2 files changed, 186 insertions(+), 129 deletions(-) create mode 100644 code/agents.py diff --git a/code/agents.py b/code/agents.py new file mode 100644 index 0000000..d39e7f9 --- /dev/null +++ b/code/agents.py @@ -0,0 +1,162 @@ +""" +Demonstrates an example estimation of microeconomic dynamic stochastic optimization +problem, as described in Section 9 of Chris Carroll's SolvingMicroDSOPs.pdf notes. +The estimation attempts to match the age-conditional wealth profile of simulated +consumers to the median wealth holdings of seven age groups in the 2004 SCF by +varying only two parameters: the coefficient of relative risk aversion and a scaling +factor for an age-varying sequence of discount factors. The estimation uses a +consumption-saving model with idiosyncratic shocks to permanent and transitory +income as defined in ConsIndShockModel. +""" + +# Parameters for the consumer type and the estimation + +import numpy as np # Numeric Python +from HARK.ConsumptionSaving.ConsBequestModel import ( + BequestWarmGlowConsumerType, + BequestWarmGlowPortfolioType, +) + +# Import modules from core HARK libraries: +# The consumption-saving micro model +from HARK.ConsumptionSaving.ConsIndShockModel import IndShockConsumerType +from HARK.ConsumptionSaving.ConsPortfolioModel import PortfolioConsumerType +from HARK.ConsumptionSaving.ConsWealthPortfolioModel import WealthPortfolioConsumerType +from HARK.core import AgentType + +# Method for sampling from a discrete distribution + +# Estimation methods + + +# Pathnames to the other files: +# Relative directory for primitive parameter files +calibration_dir = "code/calibration/" +# Relative directory for primitive parameter files +tables_dir = "content/tables/" +# Relative directory for primitive parameter files +figures_dir = "content/figures/" +# Relative directory for primitive parameter files +code_dir = "code/" + + +# Set booleans to determine which tasks should be done +# Which agent type to estimate ("IndShock" or "Portfolio") +local_estimation_agent = "IndShock" +local_estimate_model = True # Whether to estimate the model +# Whether to get standard errors via bootstrap +local_compute_standard_errors = False +# Whether to compute a measure of estimates' sensitivity to moments +local_compute_sensitivity = True +# Whether to make a contour map of the objective function +local_make_contour_plot = True + + +# ===================================================== +# Define objects and functions used for the estimation +# ===================================================== + + +class TempConsumerType(AgentType): + def __init__(self, cycles=1, time_flow=True, **kwds): + """ + Make a new consumer type. + + Parameters + ---------- + cycles : int + Number of times the sequence of periods should be solved. + time_flow : boolean + Whether time is currently "flowing" forward for this instance. + + Returns + ------- + None + """ + # Initialize a basic AgentType + super().__init__(cycles=cycles, time_flow=time_flow, **kwds) + # This estimation uses age-varying discount factors as + # estimated by Cagetti (2003), so switch from time_inv to time_vary + self.add_to_time_vary("DiscFac") + self.del_from_time_inv("DiscFac") + + def check_restrictions(self): + return None + + def simBirth(self, which_agents): + """ + Alternate method for simulating initial states for simulated agents, drawing from a finite + distribution. Used to overwrite IndShockConsumerType.simBirth, which uses lognormal distributions. + + Parameters + ---------- + which_agents : np.array(Bool) + Boolean array of size self.AgentCount indicating which agents should be "born". + + Returns + ------- + None + """ + # Get and store states for newly born agents + # Take directly from pre-specified distribution + self.state_now["aNrm"][which_agents] = self.aNrmInit[which_agents] + # No variation in permanent income needed + self.state_now["pLvl"][which_agents] = 1.0 + # How many periods since each agent was born + self.t_age[which_agents] = 0 + # Which period of the cycle each agents is currently in + self.t_cycle[which_agents] = 0 + return None + + +class IndShkLifeCycleConsumerType(TempConsumerType, IndShockConsumerType): + """ + A very lightly edited version of IndShockConsumerType. Uses an alternate method of making new + consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. + """ + + +class PortfolioLifeCycleConsumerType(TempConsumerType, PortfolioConsumerType): + """ + A very lightly edited version of PortfolioConsumerType. Uses an alternate method of making new + consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. + """ + + def post_solve(self): + for solution in self.solution: + solution.cFunc = solution.cFuncAdj + share = solution.ShareFuncAdj + solution.ShareFuncAdj = lambda m: np.clip(share(m), 0.0, 1.0) + + +class BequestWarmGlowLifeCycleConsumerType( + TempConsumerType, BequestWarmGlowConsumerType +): + """ + A very lightly edited version of BequestWarmGlowConsumerType. Uses an alternate method of making new + consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. + """ + + +class BequestWarmGlowLifeCyclePortfolioType( + TempConsumerType, BequestWarmGlowPortfolioType +): + """ + A very lightly edited version of BequestWarmGlowPortfolioType. Uses an alternate method of making new + consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. + """ + + def post_solve(self): + for solution in self.solution: + solution.cFunc = solution.cFuncAdj + share = solution.ShareFuncAdj + solution.ShareFuncAdj = lambda m: np.clip(share(m), 0.0, 1.0) + + +class WealthPortfolioLifeCycleConsumerType( + TempConsumerType, WealthPortfolioConsumerType +): + """ + A very lightly edited version of WealthPortfolioConsumerType. Uses an alternate method of making new + consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. + """ diff --git a/code/estimation.py b/code/estimation.py index e121057..b5ce2e3 100644 --- a/code/estimation.py +++ b/code/estimation.py @@ -13,20 +13,17 @@ import code.calibration.estimation_parameters as parameters import code.calibration.setup_scf_data as data # SCF 2004 data on household wealth import csv +from code.agents import ( + BequestWarmGlowLifeCycleConsumerType, + BequestWarmGlowLifeCyclePortfolioType, + IndShkLifeCycleConsumerType, + PortfolioLifeCycleConsumerType, + WealthPortfolioLifeCycleConsumerType, +) from time import time # Timing utility import matplotlib.pyplot as plt import numpy as np # Numeric Python -from HARK.ConsumptionSaving.ConsBequestModel import ( - BequestWarmGlowConsumerType, - BequestWarmGlowPortfolioType, -) - -# Import modules from core HARK libraries: -# The consumption-saving micro model -from HARK.ConsumptionSaving.ConsIndShockModel import IndShockConsumerType -from HARK.ConsumptionSaving.ConsPortfolioModel import PortfolioConsumerType -from HARK.ConsumptionSaving.ConsWealthPortfolioModel import WealthPortfolioConsumerType # Method for sampling from a discrete distribution from HARK.distribution import DiscreteDistribution @@ -35,6 +32,9 @@ from HARK.estimation import bootstrap_sample_from_data, minimize_nelder_mead from scipy.optimize import approx_fprime +# Import modules from core HARK libraries: +# The consumption-saving micro model + # Pathnames to the other files: # Relative directory for primitive parameter files @@ -64,122 +64,17 @@ # ===================================================== -class TempConsumerType: - def __init__(self, cycles=1, time_flow=True, **kwds): - """ - Make a new consumer type. - - Parameters - ---------- - cycles : int - Number of times the sequence of periods should be solved. - time_flow : boolean - Whether time is currently "flowing" forward for this instance. - - Returns - ------- - None - """ - # Initialize a basic AgentType - super().__init__(cycles=cycles, time_flow=time_flow, **kwds) - # This estimation uses age-varying discount factors as - # estimated by Cagetti (2003), so switch from time_inv to time_vary - self.add_to_time_vary("DiscFac") - self.del_from_time_inv("DiscFac") - - def check_restrictions(self): - return None - - def simBirth(self, which_agents): - """ - Alternate method for simulating initial states for simulated agents, drawing from a finite - distribution. Used to overwrite IndShockConsumerType.simBirth, which uses lognormal distributions. - - Parameters - ---------- - which_agents : np.array(Bool) - Boolean array of size self.AgentCount indicating which agents should be "born". - - Returns - ------- - None - """ - # Get and store states for newly born agents - # Take directly from pre-specified distribution - self.state_now["aNrm"][which_agents] = self.aNrmInit[which_agents] - # No variation in permanent income needed - self.state_now["pLvl"][which_agents] = 1.0 - # How many periods since each agent was born - self.t_age[which_agents] = 0 - # Which period of the cycle each agents is currently in - self.t_cycle[which_agents] = 0 - return None - - -class IndShkLifeCycleConsumerType(TempConsumerType, IndShockConsumerType): - """ - A very lightly edited version of IndShockConsumerType. Uses an alternate method of making new - consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. - """ - - -class PortfolioLifeCycleConsumerType(TempConsumerType, PortfolioConsumerType): - """ - A very lightly edited version of PortfolioConsumerType. Uses an alternate method of making new - consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. - """ - - def post_solve(self): - for solution in self.solution: - solution.cFunc = solution.cFuncAdj - share = solution.ShareFuncAdj - solution.ShareFuncAdj = lambda m: np.clip(share(m), 0.0, 1.0) - - -class BequestWarmGlowLifeCycleConsumerType( - TempConsumerType, BequestWarmGlowConsumerType -): - """ - A very lightly edited version of BequestWarmGlowConsumerType. Uses an alternate method of making new - consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. - """ - - -class BequestWarmGlowLifeCyclePortfolioType( - TempConsumerType, BequestWarmGlowPortfolioType -): - """ - A very lightly edited version of BequestWarmGlowPortfolioType. Uses an alternate method of making new - consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. - """ - - def post_solve(self): - for solution in self.solution: - solution.cFunc = solution.cFuncAdj - share = solution.ShareFuncAdj - solution.ShareFuncAdj = lambda m: np.clip(share(m), 0.0, 1.0) - - -class WealthPortfolioLifeCycleConsumerType( - TempConsumerType, WealthPortfolioConsumerType -): - """ - A very lightly edited version of WealthPortfolioConsumerType. Uses an alternate method of making new - consumers and specifies DiscFac as being age-dependent. Called "temp" because only used here. - """ - - def weighted_median(values, weights): - inds = np.argsort(values) - values = values[inds] - weights = weights[inds] + sorted_indices = np.argsort(values) + sorted_values = values[sorted_indices] + sorted_weights = weights[sorted_indices] - wsum = np.cumsum(inds) - ind = np.where(wsum > wsum[-1] / 2)[0][0] - - median = values[ind] + cumulative_weights = np.cumsum(sorted_weights) + median_index = np.searchsorted( + cumulative_weights, cumulative_weights[-1] / 2, side="right" + ) - return median + return sorted_values[median_index] def get_targeted_moments( @@ -440,7 +335,7 @@ def smmObjectiveFxnReduced(parameters_to_estimate): # Calculate minutes and remaining seconds minutes, seconds = divmod(time_to_estimate, 60) - print(f"Time to estimate: {minutes:.2f} min, {seconds:.2f} sec.") + print(f"Time to estimate: {int(minutes)} min, {int(seconds)} sec.") print(f"Estimated values: DiscFacAdj={model_estimate[0]}, CRRA={model_estimate[1]}") @@ -469,7 +364,7 @@ def compute_standard_errors_bootstrap( t_bootstrap_guess = time_to_estimate * parameters.bootstrap_size minutes, seconds = divmod(t_bootstrap_guess, 60) - print(f"This will take approximately {minutes:.2f} min, {seconds:.2f} sec.") + print(f"This will take approximately {int(minutes)} min, {int(seconds)} sec.") t_start_bootstrap = time() std_errors = calculateStandardErrorsByBootstrap( @@ -484,7 +379,7 @@ def compute_standard_errors_bootstrap( # Calculate minutes and remaining seconds minutes, seconds = divmod(time_to_bootstrap, 60) - print(f"Time to bootstrap: {minutes:.2f} min, {seconds:.2f} sec.") + print(f"Time to bootstrap: {int(minutes)} min, {int(seconds)} sec.") print(f"Standard errors: DiscFacAdj--> {std_errors[0]}, CRRA--> {std_errors[1]}") @@ -602,7 +497,7 @@ def make_contour_plot_obj_func( # Calculate minutes and remaining seconds minutes, seconds = divmod(time_to_contour, 60) - print(f"Time to contour: {minutes:.2f} min, {seconds:.2f} sec.") + print(f"Time to contour: {int(minutes)} min, {int(seconds)} sec.") plt.colorbar(smm_contour) plt.plot(model_estimate[1], model_estimate[0], "*r", ms=15) @@ -773,7 +668,7 @@ def smmObjectiveFxnReduced(parameters_to_estimate): # Calculate minutes and remaining seconds minutes, seconds = divmod(time_to_estimate, 60) - print(f"Time to estimate: {minutes:.2f} min, {seconds:.2f} sec.") + print(f"Time to estimate: {int(minutes)} min, {int(seconds)} sec.") print( f"Estimated values: DiscFacAdj={model_estimate[0]}, CRRA={model_estimate[1]}" @@ -887,7 +782,7 @@ def simulate_moments_reduced(x): # Calculate minutes and remaining seconds minutes, seconds = divmod(time_to_contour, 60) - print(f"Time to contour: {minutes:.2f} min, {seconds:.2f} sec.") + print(f"Time to contour: {int(minutes)} min, {int(seconds)} sec.") fig_sensitivity.savefig(figures_dir + "/AllSensitivity.pdf") fig_sensitivity.savefig(figures_dir + "/AllSensitivity.png")