Skip to content

Commit

Permalink
reorganize and rename
Browse files Browse the repository at this point in the history
  • Loading branch information
alanlujan91 committed Feb 22, 2024
1 parent 7eb8e46 commit 018abca
Show file tree
Hide file tree
Showing 9 changed files with 176 additions and 424 deletions.
6 changes: 3 additions & 3 deletions code/do_all.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,7 @@
still run.
"""

from estimark.calibration.options import (
from estimark.options import (
all_replications,
high_resource,
low_resource,
Expand Down Expand Up @@ -137,11 +137,11 @@ def run_replication():
return

if subjective_markets == "2" or subjective_markets == "4":
replication_specs["subjective_stock_market"] = True
replication_specs["subjective_stock"] = True
print("Adding subjective stock market beliefs...")

if subjective_markets == "3" or subjective_markets == "4":
replication_specs["subjective_labor_market"] = True
replication_specs["subjective_labor"] = True
print("Adding subjective labor market beliefs...")

estimate(**replication_specs)
Expand Down
4 changes: 0 additions & 4 deletions code/do_figs.py

This file was deleted.

Original file line number Diff line number Diff line change
Expand Up @@ -3,27 +3,12 @@
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.datasets.life_tables.us_ssa.SSATools import parse_ssa_life_table
from pathlib import Path
from HARK.ConsumptionSaving.ConsIndShockModel import init_lifecycle

# ---------------------------------------------------------------------------------
# Debugging flags
# ---------------------------------------------------------------------------------

show_PermGroFacAgg_error = False
# Error Notes:
# This sets a "quick fix" to the error, AttributeError: 'TempConsumerType' object has no attribute 'PermGroFacAgg'
# If you set this flag to "True" you will see the error. A more thorough fix is to
# fix the place where this error was introduced (Set to "True" and it will appear;
# this was almost certainly introduced when the code was extended to be used in the
# GE setting). An even more thorough solution, which moves beyond the scope of
# fixing this error, is adding unit tests to ID when changes to some code will
# break things elsewhere.
# Note: alternatively, decide that the "init_consumer_objects['PermGroFacAgg'] = 1.0"
# line below fixes it properly ('feature not a bug') and remove all this text.
from HARK.datasets.life_tables.us_ssa.SSATools import parse_ssa_life_table

# ---------------------------------------------------------------------------------
# - Define all of the model parameters for EstimatingMicroDSOPs and ConsumerExamples -
Expand All @@ -32,9 +17,7 @@
exp_nest = 1 # Number of times to "exponentially nest" when constructing a_grid
aXtraMin = 0.001 # Minimum end-of-period "assets above minimum" value
aXtraMax = 100 # Maximum end-of-period "assets above minimum" value
aXtraHuge = None # A very large value of assets to add to the grid, not used
aXtraExtra = None # Some other value of assets to add to the grid, not used
aXtraCount = 200 # Number of points in the grid of "assets above minimum"
aXtraCount = 100 # Number of points in the grid of "assets above minimum"

# Artificial borrowing constraint; imposed minimum level of end-of period assets
BoroCnstArt = 0.0
Expand All @@ -55,65 +38,64 @@
final_age = 90 # Age at which the problem ends (die with certainty)
retirement_age = 65 # Age at which the consumer retires
initial_age = 25 # Age at which the consumer enters the model
TT = final_age - initial_age # Total number of periods in the model
terminal_t = final_age - initial_age # Total number of periods in the model
retirement_t = retirement_age - initial_age - 1

# Initial guess of the coefficient of relative risk aversion during estimation (rho)
CRRA_start = 5.0
init_CRRA = 5.0
# Initial guess of the adjustment to the discount factor during estimation (beth)
DiscFacAdj_start = 0.99
init_DiscFacAdj = 0.99
# Bounds for beth; if violated, objective function returns "penalty value"
DiscFacAdj_bound = [0.0001, 15.0]
bounds_DiscFacAdj = [0.0001, 15.0]
# Bounds for rho; if violated, objective function returns "penalty value"
CRRA_bound = [0.0001, 15.0]
bounds_CRRA = [0.0001, 15.0]

# Income
ss_variances = True
income_spec = CGM_income["HS"]
# 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_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"
DiscFac_timevary = np.genfromtxt(csv_file_path) * 0.0 + 1.0
timevary_DiscFac = np.genfromtxt(csv_file_path) * 0.0 + 1.0 # todo
constant_DiscFac = np.ones_like(timevary_DiscFac)

# Survival probabilities over the lifecycle
liv_prb = parse_ssa_life_table(
female=False, min_age=initial_age, max_age=final_age - 1, cohort=1960
female=False,
min_age=initial_age,
max_age=final_age - 1,
cohort=1960,
)

# 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
empirical_cohort_age_groups = [
[age for age in range(start, start + 5)] for start in range(26, 61, 5)
]


# Three point discrete distribution of initial w
initial_wealth_income_ratio_vals = np.array([0.17, 0.5, 0.83])
init_w_to_y = np.array([0.17, 0.5, 0.83])
# Equiprobable discrete distribution of initial w
initial_wealth_income_ratio_probs = np.array([0.33333, 0.33333, 0.33334])
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 = 31382 # Just an integer to seed the estimation

options = {
"initial_wealth_income_ratio_vals": initial_wealth_income_ratio_vals,
"initial_wealth_income_ratio_probs": initial_wealth_income_ratio_probs,
"init_w_to_y": init_w_to_y,
"prob_w_to_y": prob_w_to_y,
"num_agents": num_agents,
"bootstrap_size": bootstrap_size,
"seed": seed,
"DiscFacAdj_start": DiscFacAdj_start,
"CRRA_start": CRRA_start,
"DiscFacAdj_bound": DiscFacAdj_bound,
"CRRA_bound": CRRA_bound,
"DiscFac_timevary": DiscFac_timevary,
"init_DiscFacAdj": init_DiscFacAdj,
"init_CRRA": init_CRRA,
"bounds_DiscFacAdj": bounds_DiscFacAdj,
"bounds_CRRA": bounds_CRRA,
"timevary_DiscFac": timevary_DiscFac,
}

# -----------------------------------------------------------------------------
Expand All @@ -124,28 +106,28 @@
init_consumer_objects = {
**init_lifecycle,
**{
"CRRA": CRRA_start,
"CRRA": init_CRRA,
"Rfree": Rfree,
"PermGroFac": inc_calib["PermGroFac"],
"PermGroFacAgg": 1.0,
"BoroCnstArt": BoroCnstArt,
"PermShkStd": inc_calib["PermShkStd"],
"PermShkCount": PermShkCount,
"TranShkStd": inc_calib["TranShkStd"],
"TranShkCount": TranShkCount,
"T_cycle": TT,
"T_cycle": terminal_t,
"UnempPrb": UnempPrb,
"UnempPrbRet": UnempPrbRet,
"T_retire": retirement_t,
"T_age": TT,
"T_age": terminal_t,
"IncUnemp": IncUnemp,
"IncUnempRet": IncUnempRet,
"aXtraMin": aXtraMin,
"aXtraMax": aXtraMax,
"aXtraCount": aXtraCount,
"aXtraExtra": [aXtraExtra, aXtraHuge],
"aXtraNestFac": exp_nest,
"LivPrb": liv_prb,
"DiscFac": DiscFac_timevary,
"DiscFac": timevary_DiscFac,
"AgentCount": num_agents,
"seed": seed,
"tax_rate": 0.0,
Expand All @@ -158,24 +140,22 @@
ElnR = 0.020
VlnR = 0.424**2

init_subjective_stock_market = {
TrueElnR = 0.085
TrueVlnR = 0.170**2

init_subjective_stock = {
"Rfree": 1.019, # from Mateo's JMP
"RiskyAvg": np.exp(ElnR + 0.5 * VlnR),
"RiskyStd": np.sqrt(np.exp(2 * ElnR + VlnR) * (np.exp(VlnR) - 1)),
"RiskyAvgTrue": np.exp(TrueElnR + 0.5 * TrueVlnR),
"RiskyStdTrue": np.sqrt(np.exp(2 * TrueElnR + TrueVlnR) * (np.exp(TrueVlnR) - 1)),
}

init_subjective_labor_market = { # from Tao's JMP
init_subjective_labor = { # from Tao's JMP
"TranShkStd": [0.03] * len(inc_calib["TranShkStd"]),
"PermShkStd": [0.03] * len(inc_calib["PermShkStd"]),
}

if show_PermGroFacAgg_error:
pass # do nothing
else:
# print(
# "***NOTE: using a 'quick fix' for an attribute error. See 'Error Notes' in EstimationParameter.py for further discussion.***"
# )
init_consumer_objects["PermGroFacAgg"] = 1.0

if __name__ == "__main__":
print("Sorry, estimation_parameters doesn't actually do anything on its own.")
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -2,54 +2,53 @@
Sets up the SCF data for use in the EstimatingMicroDSOPs estimation.
"""

from estimark.calibration.estimation_parameters import (
empirical_cohort_age_groups,
initial_age,
)
from pathlib import Path

import numpy as np # Numerical Python
import pandas as pd


from pathlib import Path
from estimark.calibration.parameters import initial_age

# Get the directory containing the current file and construct the full path to the CSV file
csv_file_path = Path(__file__).resolve().parent / ".." / "data" / "SCFdata.csv"

# Read the CSV file
scf_data = pd.read_csv(csv_file_path)
scf_full_data = pd.read_csv(csv_file_path)


# Keep only observations with normal incomes > 0,
# otherwise wealth/income is not well defined
scf_data = scf_data[scf_data.norminc > 0.0]
scf_full_data = scf_full_data[scf_full_data.norminc > 0.0]

# 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
age_groups = [list(range(start, start + 5)) for start in range(26, 61, 5)]

# Initialize empty lists for the data
w_to_y_data = [] # Ratio of wealth to permanent income
empirical_weights = [] # Weighting for this observation
empirical_groups = [] # Which age group this observation belongs to (1-7)
scf_data = [] # Ratio of wealth to permanent income
scf_weights = [] # Weighting for this observation
scf_groups = [] # Which age group this observation belongs to (1-7)

# Only extract the data required from the SCF dataset
search_ages = [ages[-1] for ages in empirical_cohort_age_groups]
search_ages = [ages[-1] for ages in age_groups]
for idx, age in enumerate(search_ages):
age_data = scf_data[scf_data.age_group.str.contains(f"{age}]")]
w_to_y_data.append(age_data.wealth_income_ratio.values)
empirical_weights.append(age_data.weight.values)
age_data = scf_full_data[scf_full_data.age_group.str.contains(f"{age}]")]
scf_data.append(age_data.wealth_income_ratio.values)
scf_weights.append(age_data.weight.values)
# create a group id 1-7 for every observation
empirical_groups.append([idx + 1] * len(age_data))
scf_groups.append([idx + 1] * len(age_data))

# Convert SCF data to numpy's array format for easier math
w_to_y_data = np.concatenate(w_to_y_data)
empirical_weights = np.concatenate(empirical_weights)
empirical_groups = np.concatenate(empirical_groups)
scf_data = np.concatenate(scf_data)
scf_weights = np.concatenate(scf_weights)
scf_groups = np.concatenate(scf_groups)

# Generate a single array of SCF data, useful for resampling for bootstrap
scf_data_array = np.array([w_to_y_data, empirical_groups, empirical_weights]).T
scf_array = np.array([scf_data, scf_groups, scf_weights]).T

# Generate a mapping between the real ages in the groups and the indices of simulated data
simulation_map_cohorts_to_age_indices = []
for ages in empirical_cohort_age_groups:
simulation_map_cohorts_to_age_indices.append(np.array(ages) - initial_age)
scf_mapping = []
for ages in age_groups:
scf_mapping.append(np.array(ages) - initial_age)

if __name__ == "__main__":
print("Sorry, setup_scf_data doesn't actually do anything on its own.")
Expand Down
Loading

0 comments on commit 018abca

Please sign in to comment.