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

Small tweaks to CPM example file #759

Merged
merged 4 commits into from
Jul 16, 2020
Merged

Small tweaks to CPM example file #759

Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
9b6da62
Increase no. points to make closing gap evident
Mv77 Jul 15, 2020
c4d90aa
Drop CGM calibration, use defaults
Mv77 Jul 15, 2020
4e13371
More informative legend
Mv77 Jul 15, 2020
7284cf6
Merge branch 'master' into example_CPM-tweaks
sbenthall Jul 16, 2020
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
81 changes: 17 additions & 64 deletions examples/ConsumptionSaving/example_ConsPortfolioModel.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -22,13 +22,6 @@
MyType.ShareFunc = [MyType.solution[t].ShareFuncAdj for t in range(MyType.T_cycle)]
print('Solving an infinite horizon portfolio choice problem took ' + str(t1-t0) + ' seconds.')

# %%
# Compute the Merton-Samuelson limiting portfolio share when returns are lognormal
MyType.RiskyVar = MyType.RiskyStd**2
MyType.RiskPrem = MyType.RiskyAvg - MyType.Rfree
def RiskyShareMertSamLogNormal(RiskPrem,CRRA,RiskyVar):
return RiskPrem/(CRRA*RiskyVar)

# %%
# Plot the consumption and risky-share functions
print('Consumption function over market resources:')
Expand Down Expand Up @@ -187,66 +180,30 @@ def RiskyShareMertSamLogNormal(RiskPrem,CRRA,RiskyVar):
# %% [markdown]
# The code below tests the mathematical limits of the model.

# %%
import os

# They assume its a polinomial of age. Here are the coefficients
a=-2.170042+2.700381
b1=0.16818
b2=-0.0323371/10
b3=0.0019704/100

time_params = {'Age_born': 0, 'Age_retire': 8, 'Age_death': 9}
t_start = time_params['Age_born']
t_ret = time_params['Age_retire'] # We are currently interpreting this as the last period of work
t_end = time_params['Age_death']

# They assume retirement income is a fraction of labor income in the
# last working period
repl_fac = 0.68212

# Compute average income at each point in (working) life
f = np.arange(t_start, t_ret+1,1)
f = a + b1*f + b2*(f**2) + b3*(f**3)
det_work_inc = np.exp(f)

# Retirement income
det_ret_inc = repl_fac*det_work_inc[-1]*np.ones(t_end - t_ret)

# Get a full vector of the deterministic part of income
det_income = np.concatenate((det_work_inc, det_ret_inc))

# ln Gamma_t+1 = ln f_t+1 - ln f_t
gr_fac = np.exp(np.diff(np.log(det_income)))

# Now we have growth factors for T_end-1 periods.

# Finally define the normalization factor used by CGM, for plots.
# ### IMPORTANT ###
# We adjust this normalization factor for what we believe is a typo in the
# original article. See the REMARK jupyter notebook for details.
norm_factor = det_income * np.exp(0)

# Create a grid of market resources for the plots

# Create a grid of market resources for the plots
mMin = 0 # Minimum ratio of assets to income to plot
mMax = 1e4 # Maximum ratio of assets to income to plot
mMax = 5*1e2 # Maximum ratio of assets to income to plot
mPts = 1000 # Number of points to plot

eevalgrid = np.linspace(0,mMax,mPts) # range of values of assets for the plot

# Number of points that will be used to approximate the risky distribution
risky_count_grid = [5,50]
risky_count_grid = [5,200]
# Plot by ages (time periods) at which to plot. We will use the default
# life-cycle calibration, which has 10 periods.
ages = [2, 4, 6, 8]

# Create a function to compute the Merton-Samuelson limiting portfolio share.
def RiskyShareMertSamLogNormal(RiskPrem,CRRA,RiskyVar):
return RiskPrem/(CRRA*RiskyVar)

# %% Calibration and solution

for rcount in risky_count_grid:

# Create a new dictionary and replace the number of points that
# approximate the risky return distribution

# Create new dictionary
# Create new dictionary copying the default
merton_dict = init_lifecycle.copy()
merton_dict['RiskyCount'] = rcount

Expand All @@ -256,7 +213,7 @@ def RiskyShareMertSamLogNormal(RiskPrem,CRRA,RiskyVar):

# Compute the analytical Merton-Samuelson limiting portfolio share
RiskyVar = agent.RiskyStd**2
RiskPrem = agent.RiskyAvg - agent.Rfree
RiskPrem = agent.RiskyAvg - agent.Rfree
MS_limit = RiskyShareMertSamLogNormal(RiskPrem,
agent.CRRA,
RiskyVar)
Expand All @@ -266,19 +223,15 @@ def RiskyShareMertSamLogNormal(RiskPrem,CRRA,RiskyVar):
agent.updateShareLimit()
NU_limit = agent.ShareLimit

# Plot by ages
ages = [2, 4, 6, 8]
age_born = time_params['Age_born']
plt.figure()
for a in ages:
plt.plot(eevalgrid,
agent.solution[a-age_born]\
.ShareFuncAdj(eevalgrid/
norm_factor[a-age_born]),
label = 'Age = %i' %(a))

plt.axhline(MS_limit, c='k', ls='--', label = 'M&S Limit')
plt.axhline(NU_limit, c='k', ls='-.', label = 'Numer. Limit')
agent.solution[a]\
.ShareFuncAdj(eevalgrid),
label = 't = %i' %(a))

plt.axhline(NU_limit, c='k', ls='-.', label = 'Exact limit as $m\\rightarrow \\infty$.')
plt.axhline(MS_limit, c='k', ls='--', label = 'M&S Limit without returns discretization.')

plt.ylim(0,1.05)
plt.xlim(eevalgrid[0],eevalgrid[-1])
Expand Down