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.

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

Transition Matrix Methods #1155

Merged
merged 2 commits into from
Oct 4, 2022
Merged

Conversation

wdu9
Copy link
Collaborator

@wdu9 wdu9 commented Jul 21, 2022

Please ensure your pull request adheres to the following guidelines:

  • [x ] Tests for new functionality/models or Tests to reproduce the bug-fix in code.
  • [x ] Updated documentation of features that add new functionality.
  • [x ] Update CHANGELOG.md with major/minor changes.

@codecov-commenter
Copy link

codecov-commenter commented Jul 21, 2022

Codecov Report

Base: 73.62% // Head: 72.83% // Decreases project coverage by -0.78% ⚠️

Coverage data is based on head (59c16d3) compared to base (e6c936b).
Patch coverage: 32.74% of modified lines in pull request are covered.

Additional details and impacted files
@@            Coverage Diff             @@
##           master    #1155      +/-   ##
==========================================
- Coverage   73.62%   72.83%   -0.79%     
==========================================
  Files          72       72              
  Lines       11509    11735     +226     
==========================================
+ Hits         8473     8547      +74     
- Misses       3036     3188     +152     
Impacted Files Coverage Δ
HARK/utilities.py 29.63% <10.71%> (-5.23%) ⬇️
HARK/ConsumptionSaving/ConsIndShockModel.py 79.27% <36.36%> (-5.99%) ⬇️
...nsumptionSaving/tests/test_IndShockConsumerType.py 76.73% <100.00%> (+1.57%) ⬆️

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

☔ View full report at Codecov.
📢 Do you have feedback about the report comment? Let us know in this issue.

@llorracc
Copy link
Collaborator

llorracc commented Jul 21, 2022 via email

@wdu9
Copy link
Collaborator Author

wdu9 commented Jul 22, 2022

Will, If you think this is ready to merge, could you add a tag to that effect to it (we have a predefined tag for that).

This is ready to be merged.

I think we discussed this before but I don't think I am able to add tags.

@alanlujan91 alanlujan91 added Ready-To-Merge Has been reviewed and when branch is updated and checks pass it should be merged automerge Status: Review Needed TSC: Review Needed Item escalated to the TSC. and removed Ready-To-Merge Has been reviewed and when branch is updated and checks pass it should be merged labels Jul 22, 2022
@AMonninger
Copy link
Collaborator

Happy to review this

@alanlujan91 alanlujan91 self-requested a review July 22, 2022 18:52
@alanlujan91 alanlujan91 removed the Ready-To-Merge Has been reviewed and when branch is updated and checks pass it should be merged label Jul 22, 2022
#Dist_pGrid is taken to cover most of the ergodic distribution
p_variance = self.PermShkStd[0]**2 #set variance of permanent income shocks
max_p = max_p_fac*(p_variance/(1-self.LivPrb[0]))**0.5 #Maximum Permanent income value
one_sided_grid = make_grid_exp_mult(1.05+1e-3, np.exp(max_p), num_points, 2)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

min_p should maybe be a parameter to this model instead of hardcoded

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Agreed, good catch.

@AMonninger
Copy link
Collaborator

For calc_ergodic_dist: (based on numerical summer school lecture)
As we are only interested in one eigenvalue, we can transform the problem:

  1. b = vector(Nx1) with element (1,1) = 1, rest zeros
  2. transition matrix (T): replace first row with zeros and element (1,1)=1
  3. solve: s = T^(-1)*b
  4. Normalize: s = s/sum(s)

Should be faster then searching for the eigenvalues with sp.linalg.eigs

Comment on lines 59 to 60
jump_to_grid,
jump_to_grid_fast,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

More descriptive names would be useful for extensions/alternative uses. In the past, we've used _fast to refer to numba-fied methods that do the same as their reference methods, but here, both are using numba.jit.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok got it.

Comment on lines 2402 to 2410
for i in range(m_density):
axtra_shifted = np.delete(aXtra_Grid,-1)
axtra_shifted = np.insert(axtra_shifted, 0,1.00000000e-04)
dist_betw_pts = aXtra_Grid - axtra_shifted
dist_betw_pts_half = dist_betw_pts/2
new_A_grid = axtra_shifted + dist_betw_pts_half
aXtra_Grid = np.concatenate((aXtra_Grid,new_A_grid))
aXtra_Grid = np.sort(aXtra_Grid)
self.dist_mGrid = aXtra_Grid
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is this "thickening" the grid of the distribution of cash-on-hand? what is particular about this construction that is different from make_grid_exp_mult? If this is doing something that can't be accomplished by make_grid_exp_mult maybe it should be its own utility method.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It does thicken if that means to increase the density. Chris had asked me to do this a while back but I'm not sure whethe rmake_grid_exp_mult can thicken, I feel as though it just concentrates points near the bottom end.

elif self.cycles > 1:
raise Exception('define_distribution_grid requires cycles = 0 or cycles = 1')

elif self.T_cycle != 0:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So this condition is (cycles == 1 and T_cycle !=0)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is T_cycle < 1 ever?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

yes this condition is meant as (cycles == 1 and T_cycle !=0).

T_cycle is never less than 1 I believe

else:
m_points = num_pointsM

if self.cycles == 0:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This could be
if (self.cycles == 0) or (self.cycles == 1 and self.T_cycle !=0):
to avoid code repetition I think

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The T_cycle !=0 problem is exclusive of the cycles=0 problem so it is not unfortunately.

LivPrb = self.LivPrb[0] # Update probability of staying alive

#New borns have this distribution (assumes start with no assets and permanent income=1)
NewBornDist = jump_to_grid(tran_shks,np.ones_like(tran_shks),shk_prbs,dist_mGrid,dist_pGrid)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is this a probability distribution? Would it be useful to use our DiscreteDistribution constructor here?

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is a function that takes in a grid, a set of values that may lie between the gridpoints, and finally the associated probabilities of each the values and returns the probability mass of each gridpoint. I don't think itd be useful to use DiscreteDistribution here

elif self.cycles > 1:
raise Exception('calc_transition_matrix requires cycles = 0 or cycles = 1')

elif self.T_cycle!= 0 :
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there a model that has cycles==1 and T_cycle==0? I'm thinking this should never happen, and if it does, it's probably a bug.

Copy link
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I'm not sure, that is why I was so specific on the conditional cases.

Copy link
Member

@alanlujan91 alanlujan91 left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Left a few in-line comments. Overall looks great!

Comments about new tools or using DiscreteDistribution can probably be addressed later.

@wdu9 wdu9 closed this Sep 9, 2022
@wdu9 wdu9 reopened this Sep 9, 2022
@llorracc
Copy link
Collaborator

@wdu9 ,

Glad it is now passing tests. Will ask @alanlujan91 to make a final review before merging. (Thanks to @AMonninger for earlier comments; maybe he could take a final look with an eye to how it might work for his durables problem; but @alanlujan91 should be the decider about merging).

@wdu, one thing I'd ask you to add to the chain is some insight about why it was failing before and is now working. What did you change? If it's another instance of our tests just defaulting to a degree of precision that doesn't make sense, that would be a useful insight. If it's something else, that would also be useful to know.

@wdu9
Copy link
Collaborator Author

wdu9 commented Sep 19, 2022

For calc_ergodic_dist: (based on numerical summer school lecture) As we are only interested in one eigenvalue, we can transform the problem:

1. b = vector(Nx1) with element (1,1) = 1, rest zeros

2. transition matrix (T): replace first row with zeros and element (1,1)=1

3. solve: s = T^(-1)*b

4. Normalize: s = s/sum(s)

Should be faster then searching for the eigenvalues with sp.linalg.eigs

I had this in the code however sometimes it would fail for unknown reasons. It was indeed twice as fast at very large grid dimensions.

@wdu9
Copy link
Collaborator Author

wdu9 commented Sep 19, 2022

I reduced the tolerance of the self.assertalmostequal function to resolve the original problem where the self assertion was off by a difference of 10e-6.
This minuscule difference was due to the randomness in the sp.linalg.eigs function under the calc_ergodic_dist method.
Although, I eliminated the randomness in my own code and the test was passing in my directory, it was failing in the tests here.

@AMonninger
Copy link
Collaborator

Hi Will, some additional comments:

  1. Can you make the linear algebra described above as an option?
  2. For transition matrices in a lifecycle model, it saves a lot of RAM when you make it sparse as well. Eg.
  • TranMatrix_sparse = sp.csr_matrix(TranMatrix)
  • self.tran_matrix.append(TranMatrix_sparse)
  • Helps for many grid points or large number of periods. If you use the transitionmatrices afterwards, you can use: MODEL.tran_matrix[i].A

@wdu9 wdu9 closed this Sep 30, 2022
@wdu9 wdu9 reopened this Sep 30, 2022
@alanlujan91 alanlujan91 merged commit dc6cbc6 into econ-ark:master Oct 4, 2022
@sbenthall sbenthall added this to the 0.13.0 milestone Jan 4, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

6 participants