diff --git a/code/estimark/agents.py b/code/estimark/agents.py index 1693663..0864c35 100644 --- a/code/estimark/agents.py +++ b/code/estimark/agents.py @@ -23,7 +23,6 @@ class TempConsumerType(AgentType): - def check_restrictions(self): return None diff --git a/code/estimark/estimation.py b/code/estimark/estimation.py index fffda22..27bd8c7 100644 --- a/code/estimark/estimation.py +++ b/code/estimark/estimation.py @@ -17,9 +17,6 @@ import numpy as np import pandas as pd -# Method for sampling from a discrete distribution -from HARK.distribution import DiscreteDistribution - # Estimation methods # TODO: use estimagic from HARK.estimation import bootstrap_sample_from_data, minimize_nelder_mead @@ -38,10 +35,12 @@ # Parameters for the consumer type and the estimation from estimark.parameters import ( age_mapping, + bootstrap_options, init_consumer_objects, + init_params_options, init_subjective_labor, init_subjective_stock, - options, + minimize_options, sim_mapping, ) @@ -57,54 +56,55 @@ figures_dir = "content/figures/" Path(figures_dir).mkdir(parents=True, exist_ok=True) - # ===================================================== # Define objects and functions used for the estimation # ===================================================== -def make_estimation_agent(agent_name, subjective_stock=False, subjective_labor=False): - if agent_name == "IndShock": +def make_agent( + init_agent_name, + subjective_stock=False, + subjective_labor=False, +): + if init_agent_name == "IndShock": agent_type = IndShkLifeCycleConsumerType - elif agent_name == "Portfolio": + elif init_agent_name == "Portfolio": agent_type = PortfolioLifeCycleConsumerType - elif agent_name == "WarmGlow": + elif init_agent_name == "WarmGlow": agent_type = BequestWarmGlowLifeCycleConsumerType - elif agent_name == "WarmGlowPortfolio": + elif init_agent_name == "WarmGlowPortfolio": agent_type = BequestWarmGlowLifeCyclePortfolioType - elif agent_name == "WealthPortfolio": + elif init_agent_name == "WealthPortfolio": agent_type = WealthPortfolioLifeCycleConsumerType + local_consumer_objects = init_consumer_objects.copy() + agent_name = init_agent_name + if subjective_stock or subjective_labor: agent_name += "Sub" if subjective_stock: agent_name += "(Stock)" - init_consumer_objects.update(init_subjective_stock) + local_consumer_objects.update(init_subjective_stock) if subjective_labor: agent_name += "(Labor)" - init_consumer_objects.update(init_subjective_labor) + local_consumer_objects.update(init_subjective_labor) agent_name += "Market" # Make a lifecycle consumer to be used for estimation, including simulated # shocks (plus an initial distribution of wealth) # Make a TempConsumerType for estimation - estimation_agent = agent_type(**init_consumer_objects) + agent = agent_type(**local_consumer_objects) # Set the number of periods to simulate - estimation_agent.T_sim = estimation_agent.T_cycle + 1 + agent.T_sim = agent.T_cycle + 1 # Choose to track bank balances as wealth track_vars = ["bNrm"] if "Portfolio" in agent_name: track_vars += ["Share"] - estimation_agent.track_vars = track_vars - # Draw initial assets for each consumer - estimation_agent.aNrmInit = DiscreteDistribution( - options["prob_w_to_y"], - options["init_w_to_y"], - seed=options["seed"], - ).draw(N=options["num_agents"]) - estimation_agent.make_shock_history() + agent.track_vars = track_vars - return estimation_agent, agent_name + agent.name = agent_name + + return agent def weighted_median(values, weights): @@ -120,7 +120,13 @@ def weighted_median(values, weights): return median -def get_targeted_moments(data, variable, weights=None, groups=None, mapping=None): +def get_targeted_moments( + data, + variable, + weights=None, + groups=None, + mapping=None, +): # Common variables that don't depend on whether weights are None or not data_variable = data[variable] data_groups = data[groups] @@ -142,13 +148,11 @@ def get_targeted_moments(data, variable, weights=None, groups=None, mapping=None ) else: print(f"Warning: Group {key} does not have any data.") - # Uncomment the following line if you want to assign a default value - # target_moments[key] = 0 return target_moments -def get_initial_guess(agent_name): +def get_initial_guess(agent_name, init_DiscFac=0.99, init_CRRA=2.0): # start from previous estimation results if available csv_file_path = f"{tables_dir}{agent_name}_estimate_results.csv" @@ -157,13 +161,13 @@ def get_initial_guess(agent_name): res = pd.read_csv(csv_file_path, header=None) initial_guess = res.iloc[:2, 1].astype(float).tolist() except (FileNotFoundError, IndexError): - initial_guess = [options["init_DiscFac"], options["init_CRRA"]] + initial_guess = [init_DiscFac, init_CRRA] return initial_guess # Define the objective function for the simulated method of moments estimation -def simulate_moments(params, agent, agent_name): +def simulate_moments(params, agent): """A quick check to make sure that the parameter values are within bounds. Far flung falues of DiscFac or CRRA might cause an error during solution or simulation, so the objective function doesn't even bother with them. @@ -174,38 +178,29 @@ def simulate_moments(params, agent, agent_name): agent.DiscFac = DiscFac agent.CRRA = CRRA - if "(Stock)" in agent_name and "Portfolio" in agent_name: + # Homothetic + # if hasattr(agent, "BeqCRRA"): + # agent.BeqCRRA = CRRA + + # ensure subjective beliefs are used for solution + if "(Stock)" in agent.name and "Portfolio" in agent.name: agent.RiskyAvg = init_subjective_stock["RiskyAvg"] agent.RiskyStd = init_subjective_stock["RiskyStd"] agent.update_RiskyDstn() - if "(Labor)" in agent_name: - agent.TranShkStd = init_subjective_labor["TranShkStd"] - agent.PermShkStd = init_subjective_labor["PermShkStd"] - agent.update_income_process() - - # if hasattr(agent, "BeqCRRA"): - # agent.BeqCRRA = [CRRA] * len(options["timevary_DiscFac"]) # Solve the model for these parameters, then simulate wealth data agent.solve() # Solve the microeconomic model - # "Unpack" the consumption function for convenient access - # agent.unpack("cFunc") # simulate with true parameters (override subjective beliefs) - if "(Stock)" in agent_name and "Portfolio" in agent_name: + if "(Stock)" in agent.name and "Portfolio" in agent.name: agent.RiskyAvg = agent.RiskyAvgTrue agent.RiskyStd = agent.RiskyStdTrue agent.update_RiskyDstn() - # TODO: use perceived for simulation too - if "(Labor)" in agent_name: - agent.TranShkStd = init_consumer_objects["TranShkStd"] - agent.PermShkStd = init_consumer_objects["PermShkStd"] - agent.update_income_process() - max_sim_age = agent.T_cycle + 1 # Initialize the simulation by clearing histories, resetting initial values agent.initialize_sim() + # agent.make_shock_history() agent.simulate(max_sim_age) # Simulate histories of consumption and wealth # Take "wealth" to mean bank balances before receiving labor income sim_w_history = agent.history["bNrm"] @@ -216,7 +211,7 @@ def simulate_moments(params, agent, agent_name): for key, cohort_idx in sim_mapping.items(): sim_moments[key] = np.median(sim_w_history[cohort_idx]) - if "Portfolio" in agent_name: + if "Portfolio" in agent.name: sim_share_history = agent.history["Share"] suffix = "_port" for key, cohort_idx in sim_mapping.items(): @@ -225,7 +220,7 @@ def simulate_moments(params, agent, agent_name): return sim_moments -def smm_obj_func(params, agent, moments, agent_name): +def smm_obj_func(params, agent, moments): """The objective function for the SMM estimation. Given values of discount factor adjuster DiscFac, coeffecient of relative risk aversion CRRA, a base consumer agent type, empirical data, and calibrated parameters, this function calculates @@ -272,8 +267,10 @@ def smm_obj_func(params, agent, moments, agent_name): median wealth-to-permanent-income ratio in the simulation. """ - sim_moments = simulate_moments(params, agent, agent_name) + sim_moments = simulate_moments(params, agent) + # TODO: make sure all keys in moments have a corresponding + # key in sim_moments, raise an error if not common_moments = list(set(sim_moments) & set(moments)) errors = np.array([sim_moments[key] - moments[key] for key in common_moments]) @@ -284,10 +281,9 @@ def smm_obj_func(params, agent, moments, agent_name): # Define the bootstrap procedure def calculate_se_bootstrap( - initial_estimate, - N, agent, - agent_name, + initial_estimate, + n_draws=50, seed=0, verbose=False, ): @@ -316,11 +312,11 @@ def calculate_se_bootstrap( # Generate a list of seeds for generating bootstrap samples RNG = np.random.default_rng(seed) - seed_list = RNG.integers(2**31 - 1, size=N) + seed_list = RNG.integers(2**31 - 1, size=n_draws) # Estimate the model N times, recording each set of estimated parameters estimate_list = [] - for n in range(N): + for n in range(n_draws): t_start = time() # Bootstrap a new dataset by resampling from the original data @@ -342,7 +338,6 @@ def smm_obj_func_bootstrap(params): params, agent=agent, moments=bootstrap_moments, - agent_name=agent_name, ) # Estimate the model with the bootstrap data and add to list of estimates @@ -353,7 +348,7 @@ def smm_obj_func_bootstrap(params): # Report progress of the bootstrap if verbose: print( - f"Finished bootstrap estimation #{n + 1} of {N} in {t_now - t_start} seconds ({t_now - t_0} cumulative)", + f"Finished bootstrap estimation #{n + 1} of {n_draws} in {t_now - t_start} seconds ({t_now - t_0} cumulative)", ) # Calculate the standard errors for each parameter @@ -369,7 +364,12 @@ def smm_obj_func_bootstrap(params): # ================================================================= -def do_estimate_model(agent_name, estimation_agent, target_moments, initial_guess): +def do_estimate_model( + agent, + target_moments, + initial_guess, + minimize_options=None, +): print("----------------------------------------------------------------------") print( f"Now estimating the model using Nelder-Mead from an initial guess of {initial_guess}...", @@ -384,25 +384,12 @@ def smm_obj_func_redux(params): """ return smm_obj_func( params, - agent=estimation_agent, + agent=agent, moments=target_moments, - agent_name=agent_name, ) t_start_estimate = time() - res = em.minimize( - smm_obj_func_redux, - initial_guess, - algorithm="scipy_neldermead", - upper_bounds=np.array( - [options["bounds_DiscFac"][1], options["bounds_CRRA"][1]], - ), - lower_bounds=np.array( - [options["bounds_DiscFac"][0], options["bounds_CRRA"][0]], - ), - # multistart=True, - error_handling="continue", - ) + res = em.minimize(smm_obj_func_redux, initial_guess, **minimize_options) t_end_estimate = time() time_to_estimate = t_end_estimate - t_start_estimate @@ -411,11 +398,11 @@ def smm_obj_func_redux(params): # Calculate minutes and remaining seconds minutes, seconds = divmod(time_to_estimate, 60) print(f"Time to estimate: {int(minutes)} min, {int(seconds)} sec.") - print(f"Estimated model: {agent_name}") + print(f"Estimated model: {agent.name}") print(f"Estimated values: DiscFac={model_estimate[0]}, CRRA={model_estimate[1]}") # Create the simple estimate table - estimate_results_file = tables_dir + agent_name + "_estimate_results.csv" + estimate_results_file = tables_dir + agent.name + "_estimate_results.csv" with open(estimate_results_file, "w") as f: writer = csv.writer(f) @@ -437,29 +424,29 @@ def smm_obj_func_redux(params): def do_compute_se_boostrap( - agent_name, - estimation_agent, + agent, model_estimate, time_to_estimate, + bootstrap_size=50, + seed=0, ): # Estimate the model: print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~") print( - f"Computing standard errors using {options['bootstrap_size']} bootstrap replications.", + f"Computing standard errors using {bootstrap_size} bootstrap replications.", ) print("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~") - t_bootstrap_guess = time_to_estimate * options["bootstrap_size"] + t_bootstrap_guess = time_to_estimate * bootstrap_size minutes, seconds = divmod(t_bootstrap_guess, 60) print(f"This will take approximately {int(minutes)} min, {int(seconds)} sec.") t_start_bootstrap = time() std_errors = calculate_se_bootstrap( model_estimate, - N=options["bootstrap_size"], - agent=estimation_agent, - agent_name=agent_name, - seed=options["seed"], + n_draws=bootstrap_size, + agent=agent, + seed=seed, verbose=True, ) t_end_bootstrap = time() @@ -472,7 +459,7 @@ def do_compute_se_boostrap( print(f"Standard errors: DiscFac--> {std_errors[0]}, CRRA--> {std_errors[1]}") # Create the simple bootstrap table - bootstrap_results_file = tables_dir + agent_name + "_bootstrap_results.csv" + bootstrap_results_file = tables_dir + agent.name + "_bootstrap_results.csv" with open(bootstrap_results_file, "w") as f: writer = csv.writer(f) @@ -489,14 +476,14 @@ def do_compute_se_boostrap( ) -def do_compute_sensitivity(agent_name, estimation_agent, model_estimate, initial_guess): +def do_compute_sensitivity(agent, model_estimate, initial_guess): print("``````````````````````````````````````````````````````````````````````") print("Computing sensitivity measure.") print("``````````````````````````````````````````````````````````````````````") # Find the Jacobian of the function that simulates moments def simulate_moments_redux(params): - return simulate_moments(params, agent=estimation_agent, agent_name=agent_name) + return simulate_moments(params, agent=agent) n_moments = len(scf_mapping) jac = np.array( @@ -530,14 +517,14 @@ def simulate_moments_redux(params): axs[1].set_ylabel("Sensitivity") axs[1].set_xlabel("Median W/Y Ratio") - plt.savefig(figures_dir + agent_name + "Sensitivity.pdf") - plt.savefig(figures_dir + agent_name + "Sensitivity.png") - plt.savefig(figures_dir + agent_name + "Sensitivity.svg") + plt.savefig(figures_dir + agent.name + "Sensitivity.pdf") + plt.savefig(figures_dir + agent.name + "Sensitivity.png") + plt.savefig(figures_dir + agent.name + "Sensitivity.svg") plt.show() -def do_make_contour_plot(agent_name, estimation_agent, model_estimate, target_moments): +def do_make_contour_plot(agent, model_estimate, target_moments): print("``````````````````````````````````````````````````````````````````````") print("Creating the contour plot.") print("``````````````````````````````````````````````````````````````````````") @@ -559,9 +546,8 @@ def do_make_contour_plot(agent_name, estimation_agent, model_estimate, target_mo CRRA = CRRA_list[k] smm_obj_levels[j, k] = smm_obj_func( np.array([DiscFac, CRRA]), - agent=estimation_agent, + agent=agent, moments=target_moments, - agent_name=agent_name, ) smm_contour = plt.contourf(CRRA_mesh, DiscFac_mesh, smm_obj_levels, level_count) t_end_contour = time() @@ -575,14 +561,14 @@ def do_make_contour_plot(agent_name, estimation_agent, model_estimate, target_mo plt.plot(model_estimate[1], model_estimate[0], "*r", ms=15) plt.xlabel(r"coefficient of relative risk aversion $\rho$", fontsize=14) plt.ylabel(r"discount factor adjustment $\beth$", fontsize=14) - plt.savefig(figures_dir + agent_name + "SMMcontour.pdf") - plt.savefig(figures_dir + agent_name + "SMMcontour.png") - plt.savefig(figures_dir + agent_name + "SMMcontour.svg") + plt.savefig(figures_dir + agent.name + "SMMcontour.pdf") + plt.savefig(figures_dir + agent.name + "SMMcontour.png") + plt.savefig(figures_dir + agent.name + "SMMcontour.svg") plt.show() def estimate( - agent_name, + init_agent_name, estimate_model=True, compute_se_bootstrap=False, compute_sensitivity=False, @@ -608,8 +594,8 @@ def estimate( None """ - estimation_agent, agent_name = make_estimation_agent( - agent_name=agent_name, + agent = make_agent( + init_agent_name=init_agent_name, subjective_stock=subjective_stock, subjective_labor=subjective_labor, ) @@ -622,7 +608,7 @@ def estimate( mapping=age_mapping, ) - if "Portfolio" in agent_name: + if "Portfolio" in agent.name: share_moments = get_targeted_moments( data=snp_data, variable="share", @@ -634,31 +620,30 @@ def estimate( for key, value in share_moments.items(): target_moments[key + suffix] = value - initial_guess = get_initial_guess(agent_name) + initial_guess = get_initial_guess(agent.name, **init_params_options) # Estimate the model using Nelder-Mead if estimate_model: model_estimate, time_to_estimate = do_estimate_model( - agent_name, - estimation_agent, + agent, target_moments, initial_guess, + minimize_options=minimize_options, ) # Compute standard errors by bootstrap if compute_se_bootstrap: do_compute_se_boostrap( - agent_name, - estimation_agent, + agent, model_estimate, time_to_estimate, + **bootstrap_options, ) # Compute sensitivity measure if compute_sensitivity: do_compute_sensitivity( - agent_name, - estimation_agent, + agent, model_estimate, initial_guess, ) @@ -666,8 +651,7 @@ def estimate( # Make a contour plot of the objective function if make_contour_plot: do_make_contour_plot( - agent_name, - estimation_agent, + agent, model_estimate, target_moments, ) @@ -689,7 +673,7 @@ def estimate( local_subjective_labor = False estimate( - agent_name=local_agent_name, + init_agent_name=local_agent_name, estimate_model=local_estimate_model, compute_se_bootstrap=local_compute_se_bootstrap, compute_sensitivity=local_compute_sensitivity, diff --git a/code/estimark/parameters.py b/code/estimark/parameters.py index c9fa756..c4138b5 100644 --- a/code/estimark/parameters.py +++ b/code/estimark/parameters.py @@ -2,13 +2,14 @@ model. The empirical data is stored in a separate csv file and is loaded in setup_scf_data. """ -from pathlib import Path - import numpy as np from HARK.Calibration.Income.IncomeTools import CGM_income, parse_income_spec from HARK.ConsumptionSaving.ConsIndShockModel import init_lifecycle from HARK.datasets.life_tables.us_ssa.SSATools import parse_ssa_life_table +# Method for sampling from a discrete distribution +from HARK.distribution import DiscreteDistribution + # --------------------------------------------------------------------------------- # - Define all of the model parameters for EstimatingMicroDSOPs and ConsumerExamples - # --------------------------------------------------------------------------------- @@ -57,37 +58,12 @@ income_spec = CGM_income["HS"] # College? # Replace retirement age income_spec["age_ret"] = retirement_age -inc_calib = parse_income_spec( - age_min=initial_age, - age_max=final_age, - **income_spec, - SabelhausSong=ss_variances, -) - -# # Age-varying discount factors over the lifecycle, lifted from Cagetti (2003) -# # Get the directory containing the current file and construct the full path to the CSV file -# csv_file_path = Path(__file__).resolve().parent / ".." / "data" / "Cagetti2003.csv" -# # timevary_DiscFac = np.genfromtxt(csv_file_path) * 0.0 + 1.0 # TODO -# # constant_DiscFac = np.ones_like(timevary_DiscFac) -# timevary_DiscFac = np.ones_like(inc_calib["PermShkStd"]) - -# Survival probabilities over the lifecycle -liv_prb = parse_ssa_life_table( - female=False, - min_age=initial_age, - max_age=final_age - 1, - cohort=1960, -) - # Three point discrete distribution of initial w init_w_to_y = np.array([0.17, 0.5, 0.83]) # Equiprobable discrete distribution of initial w prob_w_to_y = np.array([0.33333, 0.33333, 0.33334]) num_agents = 10000 # Number of agents to simulate -bootstrap_size = 50 # Number of re-estimations to do during bootstrap -seed = 1132023 # Just an integer to seed the estimation - # Age groups for the estimation: calculate average wealth-to-permanent income ratio # for consumers within each of these age groups, compare actual to simulated data @@ -109,17 +85,51 @@ remove_ages_from_scf = np.arange(61, 71) # remove retirement ages 61-70 remove_ages_from_snp = np.arange(71) # only match ages 71 and older +# Bootstrap options +bootstrap_size = 50 # Number of re-estimations to do during bootstrap +seed = 1132023 # Just an integer to seed the estimation + +params_to_estimate = ["CRRA", "DiscFac"] +init_params = [init_CRRA, init_DiscFac] +init_bounds = [bounds_CRRA, bounds_DiscFac] -options = { - "init_w_to_y": init_w_to_y, - "prob_w_to_y": prob_w_to_y, - "num_agents": num_agents, +inc_calib = parse_income_spec( + age_min=initial_age, + age_max=final_age, + **income_spec, + SabelhausSong=ss_variances, +) + +# Survival probabilities over the lifecycle +liv_prb = parse_ssa_life_table( + female=False, + min_age=initial_age, + max_age=final_age - 1, + cohort=1960, +) + +aNrmInit = DiscreteDistribution( + prob_w_to_y, + init_w_to_y, + seed=seed, +).draw(N=num_agents) + +init_params_options = { + "init_CRRA": init_CRRA, + "init_DiscFac": init_DiscFac, +} + +bootstrap_options = { "bootstrap_size": bootstrap_size, "seed": seed, - "init_DiscFac": init_DiscFac, - "init_CRRA": init_CRRA, - "bounds_DiscFac": bounds_DiscFac, - "bounds_CRRA": bounds_CRRA, +} + +minimize_options = { + "algorithm": "scipy_neldermead", + "upper_bounds": np.array([bounds_DiscFac[1], bounds_CRRA[1]]), + "lower_bounds": np.array([bounds_DiscFac[0], bounds_CRRA[0]]), + "multistart": True, + "error_handling": "continue", } # ----------------------------------------------------------------------------- @@ -156,6 +166,7 @@ "tax_rate": 0.0, "vFuncBool": vFuncBool, "CubicBool": CubicBool, + "aNrmInit": aNrmInit, } # from Mateo's JMP for College Educated @@ -178,7 +189,6 @@ "PermShkStd": [0.03] * len(inc_calib["PermShkStd"]), } - if __name__ == "__main__": print("Sorry, estimation_parameters doesn't actually do anything on its own.") print("This module is imported by estimation, providing calibrated ") diff --git a/code/estimark/scf.py b/code/estimark/scf.py index eca1d37..4945cf2 100644 --- a/code/estimark/scf.py +++ b/code/estimark/scf.py @@ -1,5 +1,4 @@ -"""Sets up the SCF data for use in the EstimatingMicroDSOPs estimation. -""" +"""Sets up the SCF data for use in the EstimatingMicroDSOPs estimation.""" from pathlib import Path diff --git a/code/estimark/snp.py b/code/estimark/snp.py index 175803b..1e487ec 100644 --- a/code/estimark/snp.py +++ b/code/estimark/snp.py @@ -1,5 +1,4 @@ -"""Sets up the S&P data for use in the EstimatingMicroDSOPs estimation. -""" +"""Sets up the S&P data for use in the EstimatingMicroDSOPs estimation.""" from pathlib import Path diff --git a/code/notebooks/testing_notebook.ipynb b/code/notebooks/testing_notebook.ipynb new file mode 100644 index 0000000..e0a67ee --- /dev/null +++ b/code/notebooks/testing_notebook.ipynb @@ -0,0 +1,104 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [], + "source": [ + "from estimark.agents import IndShkLifeCycleConsumerType\n", + "from estimark.parameters import (\n", + " init_consumer_objects,\n", + " init_subjective_labor,\n", + " init_subjective_stock,\n", + ")" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [], + "source": [ + "agent = IndShkLifeCycleConsumerType(**init_consumer_objects)" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [], + "source": [ + "pre = agent.PermShkDstn[0].atoms" + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": {}, + "outputs": [], + "source": [ + "agent.PermShkStd = init_subjective_labor[\"PermShkStd\"]\n", + "agent.update_income_process()" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": {}, + "outputs": [], + "source": [ + "post = agent.PermShkDstn[0].atoms" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "array([[-0.29788637, -0.18268043, -0.10486472, -0.0302585 , 0.05262022,\n", + " 0.16142052, 0.40164928]])" + ] + }, + "execution_count": 16, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "pre - post" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "estimatingmicrodsops", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.10.13" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +} diff --git a/code/run_all.py b/code/run_all.py index 12d93e4..b0f9565 100644 --- a/code/run_all.py +++ b/code/run_all.py @@ -24,20 +24,20 @@ def run_replication(): replication_specs = {} if which_model == "1" or which_model == "": - replication_specs["agent_name"] = "IndShock" + replication_specs["init_agent_name"] = "IndShock" elif which_model == "2": - replication_specs["agent_name"] = "Portfolio" + replication_specs["init_agent_name"] = "Portfolio" elif which_model == "3": - replication_specs["agent_name"] = "WarmGlow" + replication_specs["init_agent_name"] = "WarmGlow" elif which_model == "4": - replication_specs["agent_name"] = "WarmGlowPortfolio" + replication_specs["init_agent_name"] = "WarmGlowPortfolio" elif which_model == "5": - replication_specs["agent_name"] = "WealthPortfolio" + replication_specs["init_agent_name"] = "WealthPortfolio" else: print("Invalid model choice.") return - print("Model: ", replication_specs["agent_name"]) + print("Model: ", replication_specs["init_agent_name"]) if which_replication == "q": return diff --git a/code/tests.py b/code/tests.py index ae6eed1..c09474b 100644 --- a/code/tests.py +++ b/code/tests.py @@ -14,19 +14,19 @@ def test_medium_resource(): def test_portfolio_low_resource(): print("Running medium-resource replication...") - estimate(**low_resource, agent_name="Portfolio") + estimate(**low_resource, init_agent_name="Portfolio") def test_warmglow_low_resource(): print("Running medium-resource replication...") - estimate(**low_resource, agent_name="WarmGlow") + estimate(**low_resource, init_agent_name="WarmGlow") def test_warmglowportfolio_low_resource(): print("Running medium-resource replication...") - estimate(**low_resource, agent_name="WarmGlowPortfolio") + estimate(**low_resource, init_agent_name="WarmGlowPortfolio") def test_wealthportfolio_low_resource(): print("Running medium-resource replication...") - estimate(**low_resource, agent_name="WealthPortfolio") + estimate(**low_resource, init_agent_name="WealthPortfolio") diff --git a/content/tables/IndShockSub(Labor)Market_estimate_results.csv b/content/tables/IndShockSub(Labor)Market_estimate_results.csv index 7d82142..eebce10 100644 --- a/content/tables/IndShockSub(Labor)Market_estimate_results.csv +++ b/content/tables/IndShockSub(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.9610311710007147 +DiscFac,0.9610311704660466 CRRA,20.0 -time_to_estimate,174.0795018672943 -params,"[0.9610311710007147, 20.0]" -criterion,114.86969241191717 -start_criterion,114.86969241191717 -start_params,"[0.9610311710007147, 20.0]" -algorithm,scipy_neldermead +time_to_estimate,609.466682434082 +params,"[0.9610311704660466, 20.0]" +criterion,114.86969240385115 +start_criterion,115.32980579195647 +start_params,"[0.977286257563632, 17.73372158221288]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,59 +n_criterion_evaluations,190 n_derivative_evaluations, -n_iterations,26 +n_iterations,97 diff --git a/content/tables/IndShockSub(Stock)(Labor)Market_estimate_results.csv b/content/tables/IndShockSub(Stock)(Labor)Market_estimate_results.csv index ed4bd91..001ce13 100644 --- a/content/tables/IndShockSub(Stock)(Labor)Market_estimate_results.csv +++ b/content/tables/IndShockSub(Stock)(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.9672837526561269 +DiscFac,0.967283752656127 CRRA,20.0 -time_to_estimate,180.29548740386963 -params,"[0.9672837526561269, 20.0]" -criterion,115.43440076933648 +time_to_estimate,510.584103345871 +params,"[0.967283752656127, 20.0]" +criterion,115.43440076933652 start_criterion,115.43440076933648 start_params,"[0.9672837526561269, 20.0]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,61 +n_criterion_evaluations,155 n_derivative_evaluations, -n_iterations,28 +n_iterations,79 diff --git a/content/tables/IndShockSub(Stock)Market_estimate_results.csv b/content/tables/IndShockSub(Stock)Market_estimate_results.csv index 836728e..3cfd855 100644 --- a/content/tables/IndShockSub(Stock)Market_estimate_results.csv +++ b/content/tables/IndShockSub(Stock)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.946635883013375 -CRRA,3.4163818390123857 -time_to_estimate,274.5633170604706 -params,"[0.946635883013375, 3.4163818390123857]" -criterion,263.5745293428449 -start_criterion,263.57452934480165 -start_params,"[0.9466358833023842, 3.4163818274425566]" -algorithm,scipy_neldermead +DiscFac,0.9475779695403024 +CRRA,3.3786216820593937 +time_to_estimate,908.9846773147583 +params,"[0.9475779695403024, 3.3786216820593937]" +criterion,263.5489922316287 +start_criterion,263.57452934284487 +start_params,"[0.946635883013375, 3.416381839012386]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,141 +n_criterion_evaluations,312 n_derivative_evaluations, -n_iterations,74 +n_iterations,162 diff --git a/content/tables/IndShock_estimate_results.csv b/content/tables/IndShock_estimate_results.csv index c002fbf..fd94e2f 100644 --- a/content/tables/IndShock_estimate_results.csv +++ b/content/tables/IndShock_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.9354653321855162 +DiscFac,0.9354653321855163 CRRA,3.6415024998929013 -time_to_estimate,263.70834851264954 -params,"[0.9354653321855162, 3.6415024998929013]" -criterion,249.01488060992747 -start_criterion,249.01488060992747 -start_params,"[0.9354653321855162, 3.6415024998929013]" -algorithm,scipy_neldermead +time_to_estimate,912.5639476776123 +params,"[0.9354653321855163, 3.6415024998929013]" +criterion,249.01488060992742 +start_criterion,249.01488060992742 +start_params,"[0.9354653321855163, 3.6415024998929013]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,127 +n_criterion_evaluations,304 n_derivative_evaluations, -n_iterations,64 +n_iterations,154 diff --git a/content/tables/PortfolioSub(Labor)Market_estimate_results.csv b/content/tables/PortfolioSub(Labor)Market_estimate_results.csv index be0754e..f6fa8d1 100644 --- a/content/tables/PortfolioSub(Labor)Market_estimate_results.csv +++ b/content/tables/PortfolioSub(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,1.0154424685349777 -CRRA,9.806412073325795 -time_to_estimate,355.2289125919342 -params,"[1.0154424685349777, 9.806412073325795]" -criterion,65.97027813743722 -start_criterion,65.97027813743722 +DiscFac,0.97142958573396 +CRRA,13.81956746082277 +time_to_estimate,2032.6699523925781 +params,"[0.97142958573396, 13.81956746082277]" +criterion,126.8840660215122 +start_criterion,66.43894799279691 start_params,"[1.0154424685349777, 9.806412073325795]" -algorithm,scipy_neldermead +algorithm,multistart_pounders direction,minimize n_free,2 -message,Optimization terminated successfully. +message,Identical model used in successive iterations. success,True -n_criterion_evaluations,133 +n_criterion_evaluations, n_derivative_evaluations, -n_iterations,69 +n_iterations,864 diff --git a/content/tables/PortfolioSub(Stock)(Labor)Market_estimate_results.csv b/content/tables/PortfolioSub(Stock)(Labor)Market_estimate_results.csv index ba595d1..14a8c23 100644 --- a/content/tables/PortfolioSub(Stock)(Labor)Market_estimate_results.csv +++ b/content/tables/PortfolioSub(Stock)(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,1.0215310271927214 -CRRA,6.065962292352818 -time_to_estimate,402.04997158050537 -params,"[1.0215310271927214, 6.065962292352818]" -criterion,92.7176131910088 -start_criterion,92.71761323086302 -start_params,"[1.021531027302729, 6.06596228373667]" -algorithm,scipy_neldermead +DiscFac,0.9967824545259933 +CRRA,9.902526713453287 +time_to_estimate,899.9833602905273 +params,"[0.9967824545259933, 9.902526713453287]" +criterion,97.23941740958084 +start_criterion,97.26028276149013 +start_params,"[0.9973873863798245, 9.817651707796575]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,160 +n_criterion_evaluations,303 n_derivative_evaluations, -n_iterations,70 +n_iterations,151 diff --git a/content/tables/PortfolioSub(Stock)Market_estimate_results.csv b/content/tables/PortfolioSub(Stock)Market_estimate_results.csv index 1ef8cfc..122dfd9 100644 --- a/content/tables/PortfolioSub(Stock)Market_estimate_results.csv +++ b/content/tables/PortfolioSub(Stock)Market_estimate_results.csv @@ -1,15 +1,15 @@ DiscFac,0.5 -CRRA,16.02486136184392 -time_to_estimate,168.6807029247284 -params,"[0.5, 16.02486136184392]" -criterion,276.8831690170802 -start_criterion,276.8831690170802 +CRRA,15.774479359921266 +time_to_estimate,606.3084945678711 +params,"[0.5, 15.774479359921266]" +criterion,277.68974198948183 +start_criterion,277.847051758136 start_params,"[0.5, 16.02486136184392]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,82 +n_criterion_evaluations,182 n_derivative_evaluations, -n_iterations,42 +n_iterations,96 diff --git a/content/tables/WarmGlowPortfolioSub(Stock)(Labor)Market_estimate_results.csv b/content/tables/WarmGlowPortfolioSub(Stock)(Labor)Market_estimate_results.csv index 4dba77f..fd26ec5 100644 --- a/content/tables/WarmGlowPortfolioSub(Stock)(Labor)Market_estimate_results.csv +++ b/content/tables/WarmGlowPortfolioSub(Stock)(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,1.015454607302058 -CRRA,11.349728320202992 -time_to_estimate,324.64174580574036 -params,"[1.015454607302058, 11.349728320202992]" -criterion,72.61134217197446 -start_criterion,72.61134217197446 -start_params,"[1.015454607302058, 11.349728320202992]" -algorithm,scipy_neldermead +DiscFac,0.9970513630476063 +CRRA,10.081230274824245 +time_to_estimate,921.4927318096161 +params,"[0.9970513630476063, 10.081230274824245]" +criterion,97.26569375298486 +start_criterion,97.34158517443015 +start_params,"[0.9954560516566981, 10.380978428407523]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,117 +n_criterion_evaluations,315 n_derivative_evaluations, -n_iterations,59 +n_iterations,157 diff --git a/content/tables/WarmGlowPortfolioSub(Stock)Market_estimate_results.csv b/content/tables/WarmGlowPortfolioSub(Stock)Market_estimate_results.csv index b022e97..b85a51b 100644 --- a/content/tables/WarmGlowPortfolioSub(Stock)Market_estimate_results.csv +++ b/content/tables/WarmGlowPortfolioSub(Stock)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.8401378419885737 -CRRA,6.9449399419378 -time_to_estimate,246.84836077690125 -params,"[0.8401378419885737, 6.9449399419378]" -criterion,262.6929560482342 -start_criterion,262.6929560482342 +DiscFac,0.9398281943218568 +CRRA,2.739125419778687 +time_to_estimate,1047.3171365261078 +params,"[0.9398281943218568, 2.739125419778687]" +criterion,242.88587325369724 +start_criterion,nan start_params,"[0.8401378419885737, 6.9449399419378]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,123 +n_criterion_evaluations,385 n_derivative_evaluations, -n_iterations,56 +n_iterations,199 diff --git a/content/tables/WarmGlowSub(Labor)Market_estimate_results.csv b/content/tables/WarmGlowSub(Labor)Market_estimate_results.csv index b6be7d5..868334d 100644 --- a/content/tables/WarmGlowSub(Labor)Market_estimate_results.csv +++ b/content/tables/WarmGlowSub(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.997603902620452 -CRRA,4.095522481064101 -time_to_estimate,401.39938735961914 -params,"[0.997603902620452, 4.095522481064101]" -criterion,132.55119356283873 -start_criterion,132.55119356289453 -start_params,"[0.9976039026204292, 4.095522481063449]" -algorithm,scipy_neldermead +DiscFac,1.0036720534088468 +CRRA,4.010479090085397 +time_to_estimate,1034.5716614723206 +params,"[1.0036720534088468, 4.010479090085397]" +criterion,127.04771996041589 +start_criterion,132.55119356283873 +start_params,"[0.997603902620452, 4.095522481064101]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,167 +n_criterion_evaluations,406 n_derivative_evaluations, -n_iterations,42 +n_iterations,157 diff --git a/content/tables/WarmGlowSub(Stock)(Labor)Market_estimate_results.csv b/content/tables/WarmGlowSub(Stock)(Labor)Market_estimate_results.csv index f63e9ac..2982208 100644 --- a/content/tables/WarmGlowSub(Stock)(Labor)Market_estimate_results.csv +++ b/content/tables/WarmGlowSub(Stock)(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,1.0079244332693773 -CRRA,4.734637245990841 -time_to_estimate,398.02410554885864 -params,"[1.0079244332693773, 4.734637245990841]" -criterion,135.75860439946902 -start_criterion,135.7586043995409 -start_params,"[1.0079244332693316, 4.734637245990195]" -algorithm,scipy_neldermead +DiscFac,1.0154158636549848 +CRRA,4.625581972362632 +time_to_estimate,1081.8968892097473 +params,"[1.0154158636549848, 4.625581972362632]" +criterion,130.421394767837 +start_criterion,135.75860439946902 +start_params,"[1.0079244332693773, 4.734637245990841]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,163 +n_criterion_evaluations,436 n_derivative_evaluations, -n_iterations,41 +n_iterations,182 diff --git a/content/tables/WarmGlowSub(Stock)Market_estimate_results.csv b/content/tables/WarmGlowSub(Stock)Market_estimate_results.csv index 52b40d7..d862d2c 100644 --- a/content/tables/WarmGlowSub(Stock)Market_estimate_results.csv +++ b/content/tables/WarmGlowSub(Stock)Market_estimate_results.csv @@ -1,15 +1,15 @@ DiscFac,0.9505752102756507 CRRA,3.1219574140853 -time_to_estimate,262.18266248703003 +time_to_estimate,887.5714428424835 params,"[0.9505752102756507, 3.1219574140853]" criterion,242.64184346164961 start_criterion,242.64184346164961 start_params,"[0.9505752102756507, 3.1219574140853]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,129 +n_criterion_evaluations,294 n_derivative_evaluations, -n_iterations,67 +n_iterations,149 diff --git a/content/tables/WarmGlow_estimate_results.csv b/content/tables/WarmGlow_estimate_results.csv index 5e1182f..27f324b 100644 --- a/content/tables/WarmGlow_estimate_results.csv +++ b/content/tables/WarmGlow_estimate_results.csv @@ -1,11 +1,11 @@ -DiscFac,0.9387761860256366 +DiscFac,0.9387761860256367 CRRA,3.3382252591815313 -time_to_estimate,264.959814786911 -params,"[0.9387761860256366, 3.3382252591815313]" -criterion,227.81638720311375 +time_to_estimate,473.51883602142334 +params,"[0.9387761860256367, 3.3382252591815313]" +criterion,227.81638720311378 start_criterion,227.81638720311375 start_params,"[0.9387761860256366, 3.3382252591815313]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. diff --git a/content/tables/WealthPortfolioSub(Labor)Market_estimate_results.csv b/content/tables/WealthPortfolioSub(Labor)Market_estimate_results.csv index 40189af..17e8e85 100644 --- a/content/tables/WealthPortfolioSub(Labor)Market_estimate_results.csv +++ b/content/tables/WealthPortfolioSub(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ DiscFac,0.914060173596014 CRRA,5.7238865307868485 -time_to_estimate,441.63514614105225 +time_to_estimate,1074.5020213127136 params,"[0.914060173596014, 5.7238865307868485]" criterion,110.24776730302995 start_criterion,110.24776730302995 start_params,"[0.914060173596014, 5.7238865307868485]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,137 +n_criterion_evaluations,287 n_derivative_evaluations, -n_iterations,76 +n_iterations,150 diff --git a/content/tables/WealthPortfolioSub(Stock)(Labor)Market_estimate_results.csv b/content/tables/WealthPortfolioSub(Stock)(Labor)Market_estimate_results.csv index 2b0365b..ba24084 100644 --- a/content/tables/WealthPortfolioSub(Stock)(Labor)Market_estimate_results.csv +++ b/content/tables/WealthPortfolioSub(Stock)(Labor)Market_estimate_results.csv @@ -1,15 +1,15 @@ DiscFac,0.8650319312758133 CRRA,16.869456069412628 -time_to_estimate,429.3529198169708 +time_to_estimate,1116.6552760601044 params,"[0.8650319312758133, 16.869456069412628]" criterion,110.45040833862312 -start_criterion,110.4504083396946 -start_params,"[0.8650319315162985, 16.86945605399145]" -algorithm,scipy_neldermead +start_criterion,110.45040833862312 +start_params,"[0.8650319312758133, 16.869456069412628]" +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,129 +n_criterion_evaluations,306 n_derivative_evaluations, -n_iterations,69 +n_iterations,159 diff --git a/content/tables/WealthPortfolioSub(Stock)Market_estimate_results.csv b/content/tables/WealthPortfolioSub(Stock)Market_estimate_results.csv index e7702b5..721c0ae 100644 --- a/content/tables/WealthPortfolioSub(Stock)Market_estimate_results.csv +++ b/content/tables/WealthPortfolioSub(Stock)Market_estimate_results.csv @@ -1,15 +1,15 @@ -DiscFac,0.8701946659837727 -CRRA,4.247835330643327 -time_to_estimate,393.2902042865753 -params,"[0.8701946659837727, 4.247835330643327]" -criterion,163.29576612377568 +DiscFac,0.8706337679058813 +CRRA,4.185375372200809 +time_to_estimate,1142.2072291374207 +params,"[0.8706337679058813, 4.185375372200809]" +criterion,163.29053671548135 start_criterion,163.29576612377568 start_params,"[0.8701946659837727, 4.247835330643327]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,158 +n_criterion_evaluations,321 n_derivative_evaluations, -n_iterations,78 +n_iterations,160 diff --git a/content/tables/WealthPortfolio_estimate_results.csv b/content/tables/WealthPortfolio_estimate_results.csv index 7c05b40..0b4b2a6 100644 --- a/content/tables/WealthPortfolio_estimate_results.csv +++ b/content/tables/WealthPortfolio_estimate_results.csv @@ -1,15 +1,15 @@ DiscFac,0.8519952862716866 CRRA,4.404726074797133 -time_to_estimate,359.2156708240509 +time_to_estimate,1066.6434442996979 params,"[0.8519952862716866, 4.404726074797133]" criterion,156.78926007862455 start_criterion,156.78926007862455 start_params,"[0.8519952862716866, 4.404726074797133]" -algorithm,scipy_neldermead +algorithm,multistart_scipy_neldermead direction,minimize n_free,2 message,Optimization terminated successfully. success,True -n_criterion_evaluations,138 +n_criterion_evaluations,287 n_derivative_evaluations, -n_iterations,75 +n_iterations,152