Merge remote-tracking branch 'origin/master' into vers

CHANGELOG.md

@@ -5,6 +5,24 @@ All notable changes to this project will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
+## [0.11.13] - 2024-07-28 12:00:00
+
+### Added
+
+- Added three new pension types to the model: (i) defined benefits system, (ii) notional defined contribution system, and (iii) points system.
+
+## [0.11.12] - 2024-07-26 01:00:00
+
+### Bug Fix
+
+- Fixes extrapolation of nested lists of tax function parameters.
+
+## [0.11.11] - 2024-06-24 01:00:00
+
+### Added
+
+- Add new parameters for resource constraint tolerances for steady state and time path solution.
+
 ## [0.11.10] - 2024-06-17 01:00:00
 
 ### Added
@@ -253,6 +271,8 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Any earlier versions of OG-USA can be found in the [`OG-Core`](https://github.com/PSLmodels/OG-Core) repository [release history](https://github.com/PSLmodels/OG-Core/releases) from [v.0.6.4](https://github.com/PSLmodels/OG-Core/releases/tag/v0.6.4) (Jul. 20, 2021) or earlier.
 
 
+[0.11.13]: https://github.com/PSLmodels/OG-Core/compare/v0.11.11...v0.11.13
+[0.11.11]: https://github.com/PSLmodels/OG-Core/compare/v0.11.10...v0.11.11
 [0.11.10]: https://github.com/PSLmodels/OG-Core/compare/v0.11.9...v0.11.10
 [0.11.9]: https://github.com/PSLmodels/OG-Core/compare/v0.11.8...v0.11.9
 [0.11.8]: https://github.com/PSLmodels/OG-Core/compare/v0.11.7...v0.11.8

docs/book/content/api/demographics.rst

@@ -0,0 +1,13 @@
+.. _demographics:
+
+Demographics Functions
+=======================
+
+**demographics.py modules**
+
+ogcore.demographics
+------------------------------------------
+
+.. automodule:: ogcore.demographics
+  :members: get_un_data, get_fert, get_mort, get_pop, pop_rebin, get_imm_rates,
+            immsolve, get_pop_objs, pre_pop_solve

docs/book/content/api/public_api.rst

@@ -12,6 +12,7 @@ There is also a link to the source code for each documented member.
    :maxdepth: 1
 
    aggregates
+   demographics
    elliptical_u_est
    execute
    firm
@@ -22,6 +23,7 @@ There is also a link to the source code for each documented member.
    parameter_plots
    parameter_tables
    parameters
+   pensions
    tax
    txfunc
    utils

docs/book/content/api/tax.rst

@@ -9,6 +9,6 @@ ogcore.tax
 ------------------------------------------
 
 .. automodule:: ogcore.tax
-  :members: replacement_rate_vals, ETR_wealth, MTR_wealth, ETR_income,
-    MTR_income, get_biz_tax, net_taxes, income_tax_liab, pension_amount,
+  :members: ETR_wealth, MTR_wealth, ETR_income,
+    MTR_income, get_biz_tax, net_taxes, income_tax_liab,
     wealth_tax_liab, bequest_tax_liab

docs/book/content/theory/government.md

@@ -317,7 +317,151 @@ Businesses face a linear tax rate $\tau^{b}_{m,t}$, which can vary by industry a
 
 #### Pensions
 
-[TODO: Add description of government pensions and the relevant parameters]
+The `OG-Core` model allows for four different systems for public pensions:
+
+1. U.S.-style social security system
+2. Defined benefit system
+3. Notional defined contribution system
+4. Points system
+
+These can be selected with the `pension_system` parameter.  Accepted values are `US-Style Social Security`, `Defined Benefits`, `Notional Defined Contribution`, `Points System`.  We discuss each of these in turn below.
+
+For all systems, $R$ represents the age at which the individual becomes eligible to receive the government provided retirement benefit.
+
+##### U.S.-style social security system
+
+Under the U.S.-style social security system, households over age $R$ received a pension amount that is a function of their earnings history.  The earings history includes the highest earning `AIME_num_periods` prior to retirement (which OG-Core assumes happens at age $R$).  This history determines the Average Indexed Monthly Earnings (AIME):
+
+  ```{math}
+  :label: eqn:AIME
+  AIME_{j,R,t+R} = {\sum_{u=0}^{AIMEper} w_{t+u}e_{j,u,t+u}n{j,u,t+u}}{ 12 * AIMEper}
+  ```
+
+ The AIME in turn, determines a household's (PIA) based on three rates and brackets:
+
+   ```{math}
+  :label: eqn:PIA
+    PIAbase_{j,R,t+R} =
+      \begin{cases}
+        PIArate_1 \times AIME_{j,R,t+R}, \text{for } AIME_{j,R,t+R} \leq AIMEbkt_1 \\
+        PIArate_2 \times AIME_{j,R,t+R}, \text{for } AIMEbkt_1 < AIME_{j,R,t+R} \leq AIMEbkt_2 \\
+        PIArate_3 \times AIME_{j,R,t+R}, \text{for } AIMEbkt_2 < AIME_{j,R,t+R} \\
+      \end{cases}
+  ```
+
+ The PIA is then capped at a maximum, set by the parameter `PIA_maxpayment`, $PIA_{j,R,t+R} = \max\{PIAbase_{j,R,t+R}, \text{ PIA max payment amount}\}$.
+
+  The replacement rate, $\theta_j$ is then calculated as annual earnings, with an adjustment for the wage rate.:
+
+  ```{math}
+  :label: eqn:theta
+    \theta_{j,R,t+R} = \frac{PIA_{j,R,t+R} \times 12}{factor \times w_{t+R}}
+  ```
+
+  Note that $aIME_{j,R,t+R}$ is a function of each households' earning history, but their choice of earning may depend on their retirement benefit.  Solving this exactly would introduce and additional fixed point problem in both the steady state solution and in the time path solution.  The latter would be extremely computationally taxing.  Therefore, we make the simplification of determining the AIME for each type $j$ using the steady state solution and earnings for type $j$ in the steady state.  This leads to some approximation error, but because $\theta_j$ is adjusted by the current wage rate, this approximation error is minized.
+
+  The pension amount for households under the US-style social security system is then:
+
+  ```{math}
+  :label: eqn:ss_pension
+    pension_{j,s,t} = \theta_j \times w_t \quad \forall s > R
+  ```
+
+
+##### Defined benefit system
+
+The defined benefit system pension amount is given as:
+
+  ```{math}
+  :label: eqn:db_pension
+  pension{j,s,t} = \biggl[\frac{\sum_{s=R-ny}^{R-1}w_{t}e_{j,s,t}n_{j,s,t}}{ny}\biggr]\times Cy \times \alpha_{DB} \quad \forall s > R
+  ```
+
+where:
+  \begin{itemize}
+    \item $ny$ are the number of years over which average earnings are calculated
+    \item $Cy$ are the number of years of contributions.  In our model, there is no exit from the labor force, so workers will contribute for $R$ years, but $Cy$ could be some number less than $R$ if there is a maximum number of years of contributions one can accrue under the DB system.
+    \item $\alpha_{DB}$ is the replacement rate per year of contribution.
+  \end{itemize}
+
+  Given this pension system and the fact that there is only variation in labor supply along the intensive margin (so we don't need to consider changes in $Cy$), the partial derivatives from the household section are given by:
+
+  ```{math}
+  :label: eqn:db_deriv
+    \frac{\partial \theta_{j,u,t+u-s}}{\partial n_{j,s,t}} =
+      \begin{cases}
+        0 , & \text{if}\ s < R - Cy \\
+        w_{t}e_{j,s}\alpha_{DB}\times \frac{Cy}{ny}, & \text{if}\  R - Cy <= s < R  \\
+        0, & \text{if}\ s \geq R \\
+      \end{cases}
+  ```
+
+##### Notional defined contribution system
+
+The pension amount under a notional defined contribution system is given as:
+
+  ```{math}
+  :label: eqn:ndc_pension
+   pension{j,s,t} = \biggl[\sum_{s=E}^{R-1}\tau^{p}_{t}w_{t}e_{j,s,t}n_{j,s,t}(1 + g_{NDC,t})^{R-s-1}\biggr]\delta_{R, t} \quad \forall s > R
+  ```
+
+where:
+
+  * $\tau^p$ is the pension contribution tax rate
+  * $g_{NDC,t}$ the rate of growth applied to contributions.
+    * For example, In the Italian system, $g_{NDC,t}$ is the mean nominal GDP growth rate in the 5 years before seniority
+    * i.e., $g_{NDC,t}=\prod_{j=i}^{R-1}\bar{g}_{j}$
+    * Note, this is not $g_y$. In the SS, it's $(\bar{g}_{y} + \bar{g}_{n})$, and in the transition path equilibrium, it's not a function of exogenous variables since the growth rate of nominal GDP is endogenous.
+  * $\delta_{R, t}$ is the conversion coefficient at time $t$ and its calculation is detailed below.
+
+  ```{math}
+  \delta_{R} = (dir_{R} + ind_{R} - k)^{-1}
+  ```
+
+where $k$ is an adjustment that takes into account the number of payments per year.  In particular, $k=0.5 - (6/13n)$, where $n$ is the number of payments per year.  So if the payments are made monthly, $n=12$ and thus $k=0.4615$.
+
+The $dir_{R, t}$ term is an adjustment to make the payments actuarially fair given mortality risk:
+
+  ```{math}
+  dir_{R, t} = \sum_{u=0}^{E+S-R}\left[\prod_{s=R}^{u}(1-\hat{\rho}_{s, t})\right](1+\hat{g}_{y, t})^{-u}
+  ```
+
+where $\hat{\rho}_{s,t}$ are the mortality tables used in the pension system at time $t$ and $\hat{g}_{y, t}$ is the long run expected nominal GDP growth rate used in the pension system at time $t$.
+
+Finally, $ind_{R}$ is an adjustment for survivor benefits.  Since we model households (and not individuals), we set $ind_{R} = 0$ by default.  This can be changed with the parameter `indR` if one would like to account for the fact that households lose members over time.
+
+Given this pension system, the partial derivatives from the household section are given by:
+
+  ```{math}
+  :label: eqn:ndc_deriv
+  \frac{\partial \theta_{j,u,t+u-s}}{\partial n_{j,s,t}} =
+    \begin{cases}
+      \tau^{p}_{t}w_{t}e_{j,s}(1+g_{NDC,t})^{u - s}\delta_{R,t}, & \text{if}\ s<R-1 \\
+      0, & \text{if}\ s \geq R \\
+    \end{cases}
+  ```
+
+##### Points system
+
+Under a points system, the pension amount is given as:
+
+  ```{math}
+  :label: eqn:ps_pension
+  pension{j,s,t} = \sum_{s=E}^{R-1}w_{t}e_{j,s,t}n_{j,s,t}\times v_{t} \quad \forall s > R
+  ```
+
+where $v_{t}$ is the value of a point at time $t$
+
+Given this pension system, the partial derivatives from the household section are given by:
+
+  ```{math}
+  :label: eqn:ps_deriv
+  \frac{\partial \theta_{j,u,t+u-s}}{\partial n_{j,s,t}} =
+  \begin{cases}
+    0 , & \text{if}\ s < R \\
+    w_{t}e_{j,s}v_{t}, & \text{if}\ s \geq R \\
+  \end{cases}
+ ```
 
 #### Lump sum transfers:
 
@@ -327,6 +471,7 @@ Businesses face a linear tax rate $\tau^{b}_{m,t}$, which can vary by industry a
   :label: EqUnbalGBCtfer
     TR_t = g_{tr,t}\:\alpha_{tr}\: p_t Y_t \quad\forall t
   ```
+
   where total government transfers to households $TR_t$ and GDP ($p_t Y_t$) are in terms of the numeraire good and the term $Y_t$ is in terms of the composite good.
 
   The time dependent multiplier $g_{tr,t}$ in front of the right-hand-side of {eq}`EqUnbalGBCtfer` will equal 1 in most initial periods. It will potentially deviate from 1 in some future periods in order to provide a closure rule that ensures a stable long-run debt-to-GDP ratio. We will discuss the closure rule in Section {ref}`SecUnbalGBCcloseRule`.
@@ -417,8 +562,10 @@ Businesses face a linear tax rate $\tau^{b}_{m,t}$, which can vary by industry a
 
   ```{math}
   :label: EqUnbalGBCgovRev
-    Rev_t &= \underbrace{\sum_{m=1}^M\Bigl[\tau^{corp}_{m,t}\bigl(p_{m,t}Y_{m,t} - w_t L_t\bigr) - \tau^{corp}_{m,t}\delta^\tau_{m,t}K_{m,t} - \tau^{inv}_{m,t}\delta_{M,t}K_{m,t}\Bigr]}_{\text{corporate tax revenue}} \\
-    &\quad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\lambda_j\omega_{s,t}\tau^{etr}_{s,t}\left(x_{j,s,t},y_{j,s,t}\right)\bigl(x_{j,s,t} + y_{j,s,t}\bigr)}_{\text{household tax revenue}} \quad\forall t
+    Rev_t &= \underbrace{\sum_{m=1}^M\Bigl[\tau^{corp}_{m,t}\bigl(p_{m,t}Y_{m,t} - w_t L_t\bigr) - \tau^{corp}_{m,t}\delta^\tau_{m,t}K_{m,t} - \tau^{inv}_{m,t}\delta_{M,t}K_{m,t}\Bigr]}_{\text{corporate income tax revenue}} \\
+    &\quad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\lambda_j\omega_{s,t}\tau^{etr}_{s,t}\left(x_{j,s,t},y_{j,s,t}\right)\bigl(x_{j,s,t} + y_{j,s,t}\bigr)}_{\text{household income tax revenue}} \\
+    &\quad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\sum_{i=1}^I\lambda_j\omega_{s,t}\tau^{c}_{i,t}p{i,t}c_{i,j,s,t}}_{\text{consumption tax revenue}} \\
+    &\quad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\lambda_j\omega_{s,t}\tau^{w}_{t}b_{j,s,t}}_{\text{wealth tax revenue}} \quad\forall t
   ```
 
   where household labor income is defined as $x_{j,s,t}\equiv w_t e_{j,s}n_{j,s,t}$ and capital income $y_{j,s,t}\equiv r_{p,t} b_{j,s,t}$.
@@ -430,7 +577,7 @@ Businesses face a linear tax rate $\tau^{b}_{m,t}$, which can vary by industry a
 
   ```{math}
   :label: EqUnbalGBCbudgConstr
-    D_{t+1} + Rev_t = (1 + r_{gov,t})D_t + G_t + I_{g,t} + TR_t + UBI_t  \quad\forall t
+    D_{t+1} + Rev_t = (1 + r_{gov,t})D_t + G_t + I_{g,t} + Pensions_t + TR_t + UBI_t  \quad\forall t
   ```
 
   where $r_{gov,t}$ is the interest rate paid by the government defined in equation {eq}`EqUnbalGBC_rate_wedge` below, $G_{t}$ is government spending on public goods, $I_{g,t}$ is total government spending on infrastructure investment, $TR_{t}$ are non-pension government transfers, and $UBI_t$ is the total UBI transfer outlays across households in time $t$.

docs/book/content/theory/stationarization.md

@@ -209,14 +209,16 @@ Where $\hat{MTR}^w_{j,s,t} = \left( \frac{h^{w}p_{w}\hat{b}_{j,s,t}}{(\hat{b}_{j
   ```{math}
   :label: EqStnrzGovRev
     \hat{Rev}_t &= \underbrace{\sum_{m=1}^M\Bigl[\tau^{corp}_{m,t}\bigl(p_{m,t}\hat{Y}_{m,t} - \hat{w}_t\hat{L}_t\bigr) - \tau^{corp}_{m,t}\delta^\tau_{m,t}\hat{K}_{m,t} - \tau^{inv}_{m,t}\hat{I}_{m,t}\Bigr]}_{\text{corporate tax revenue}} \\
-    &\qquad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\lambda_j\hat{\omega}_{s,t}\tau^{etr}_{s,t}\left(\hat{x}_{j,s,t},\hat{y}_{j,s,t}\right)\bigl(\hat{x}_{j,s,t} + \hat{y}_{j,s,t}\bigr)}_{\text{household tax revenue}} \quad\forall t
+    &\qquad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\lambda_j\hat{\omega}_{s,t}\tau^{etr}_{s,t}\left(\hat{x}_{j,s,t},\hat{y}_{j,s,t}\right)\bigl(\hat{x}_{j,s,t} + \hat{y}_{j,s,t}\bigr)}_{\text{household tax revenue}} \\
+    &\quad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\sum_{i=1}^I\lambda_j\omega_{s,t}\tau^{c}_{i,t}p{i,t}\hat{c}_{i,j,s,t}}_{\text{consumption tax revenue}} \\
+    &\quad + \underbrace{\sum_{s=E+1}^{E+S}\sum_{j=1}^J\lambda_j\omega_{s,t}\tau^{w}_{t}\hat{b}_{j,s,t}}_{\text{wealth tax revenue}} \quad\forall t
   ```
 
   Every term in the government budget constraint {eq}`EqUnbalGBCbudgConstr` is growing at both the productivity growth rate and the population growth rate, so we stationarize it by dividing both sides by $e^{g_y t}\tilde{N}_t$. We also have to multiply and divide the next period debt term $D_{t+1}$ by $e^{g_y(t+1)}\tilde{N}_{t+1}$, leaving the term $e^{g_y}(1 + \tilde{g}_{n,t+1})$.
 
   ```{math}
   :label: EqStnrzGovBC
-    e^{g_y}\left(1 + \tilde{g}_{n,t+1}\right)\hat{D}_{t+1} + \hat{Rev}_t = (1 + r_{gov,t})\hat{D}_t + \hat{G}_t + \hat{I}_{g,t} + \hat{TR}_t + \hat{UBI}_t \quad\forall t
+    e^{g_y}\left(1 + \tilde{g}_{n,t+1}\right)\hat{D}_{t+1} + \hat{Rev}_t = (1 + r_{gov,t})\hat{D}_t + \hat{G}_t + \hat{I}_{g,t} + \hat{Pensions}_t + \hat{TR}_t + \hat{UBI}_t \quad\forall t
   ```
 
   The stationarized versions of the rule for total government infrastructure investment spending $I_{g,t}$ in {eq}`EqUnbalGBC_Igt` and the rule for government investment spending in each industry in {eq}`EqUnbalGBC_Igt` are found by dividing both sides of the respective equations by $e^{g_y t}\tilde{N}_t$.
@@ -306,6 +308,105 @@ Where $\hat{MTR}^w_{j,s,t} = \left( \frac{h^{w}p_{w}\hat{b}_{j,s,t}}{(\hat{b}_{j
   \end{split}
   ```
 
+
+  ### Stationarized Pension System Equations
+
+#### Stiationarized Notional Defined Contributions Equations
+
+The stationarized NDC pension amount is given by:
+
+  ```{math}
+  :label: eqn:ndc_amount_stationarized
+  \hat{\theta}_{j,u,t+u-s}=
+    \begin{cases}
+      0, & \text{if}\ u < R \\
+      \biggl[\sum_{s=E}^{R-1}\tau^{p}_{t}\frac{\hat{w}_{t}}{e^{g_{y}(u-s)}}e_{j,s}n_{j,s,t}(1 + g_{NDC,t})^{R-s-1}\biggr]\delta_{R, t}, & \text{otherwise}
+    \end{cases}
+  ```
+
+
+The stationarized derivative of the pension amount it slightly simpler since it involved only current period wages.  We give the derivation first.
+
+The FOC for the choice of labor supply is given by:
+
+  ```{math}
+  \begin{split}
+    \biggl(\frac{1}{1 + \tau^{c}_{j,s,t}}\biggr)&\biggl(w_{t}e_{j,s} - \frac{\partial T_{j,s,t}}{\partial n_{j,s,t}}\biggr)c^{-\sigma}_{j,s,t} + \sum_{u=R}^{E+S}\beta^{u-s}\prod_{v=s}^{u}(1-\rho_{v})\frac{\partial \theta_{j,u,t+u-s}}{\partial n_{j,s,t}}c^{-\sigma}_{j,u,t+u-s}\biggl(\frac{1}{1+\tau^{c}_{j,u,t+u-s}}\biggr) \\
+    & = MDU_l(n_{j,s,t})e^{g_{y}t(1-\sigma)}
+  \end{split}
+  ```
+
+where we now pull the growth factor out of the marginal disutility of labor term to aid in the exposition of the stationarization.  To stationarize this equation, we divide both sides through by $e^{g_{y}t(1-\sigma)}$.
+
+  ```{math}
+  \begin{split}
+    \biggl(\frac{1}{1 + \tau^{c}_{j,s,t}}\biggr)&\biggl(\frac{w_{t}}{e^{g_{y}t}}e_{j,s} - \frac{\partial T_{j,s,t}}{\partial n_{j,s,t}}\biggr)\frac{c^{-\sigma}_{j,s,t}}{e^{g_{y}t(-\sigma)}} + \sum_{u=R}^{E+S}\beta^{u-s}\prod_{v=s}^{u}(1-\rho_{v})\frac{\partial \theta_{j,u,t+u-s}}{\partial n_{j,s,t}e^{g_{y}t}}\frac{c^{-\sigma}_{j,u,t+u-s}}{e^{g_{y}t(-\sigma)}}\biggl(\frac{1}{1+\tau^{c}_{j,u,t+u-s}}\biggr) \\
+    & = MDU_l(n_{j,s,t})\frac{e^{g_{y}t(1-\sigma)}}{e^{g_{y}t(1-\sigma)}}
+  \end{split}
+  ```
+
+Which we can write as:
+
+  ```{math}
+  \begin{split}
+    \biggl(\frac{1}{1 + \tau^{c}_{j,s,t}}\biggr)&\biggl(\hat{w}_{t}e_{j,s} - \frac{\partial T_{j,s,t}}{\partial n_{j,s,t}}\biggr)\hat{c}^{-\sigma}_{j,s,t} + \sum_{u=R}^{E+S}\beta^{u-s}\prod_{v=s}^{u}(1-\rho_{v})\frac{\partial \hat{\theta}_{j,u,t+u-s}}{\partial n_{j,s,t}}\hat{c}^{-\sigma}_{j,u,t+u-s}e^{g_{y}(u-s)(-\sigma)}\biggl(\frac{1}{1+\tau^{c}_{j,u,t+u-s}}\biggr) \\
+    & = MDU_l(n_{j,s,t})
+  \end{split}
+  ```
+
+where $\frac{\partial \hat{\theta}_{j,u,t+u-s}}{\partial n_{j,s,t}}$ is given by:
+
+  ```{math}
+  :label: eqn:ndc_deriv_stationarized
+  \frac{\partial \theta_{j,u,t+u-s}}{\partial n_{j,s,t}} =
+    \begin{cases}
+      \tau^{p}_{t}\hat{w}_{t}e_{j,s}(1+g_{NDC,t})^{u - s}\delta_{R,t}, & \text{if}\ s<R-1 \\
+      0, & \text{if}\ s \geq R \\
+    \end{cases}
+  ```
+
+
+#### Stationarized Defined Benefits Equations
+
+Stationarized pension amount:
+
+  ```{math}
+  :label: eqn:db_amount_staionarized
+  \hat{\theta}_{j,u,t+u-s} = \biggl[\frac{\sum_{s=R-ny}^{R-1}\frac{\hat{w}_{t}}{e^{g_{y}(u-s)}}e_{j,s}n_{j,s,t}}{ny}\biggr]\times Cy \times \alpha_{DB}, \ \ \forall u \geq R
+  ```
+
+Stationarized pension amount derivative:
+
+  ```{math}
+  :label: eqn:db_deriv_stationarized
+    \frac{\partial \hat{\theta}_{j,u,t+u-s}}{\partial n_{j,s,t}} =
+      \begin{cases}
+        0 , & \text{if}\ s < R - Cy \\
+        \hat{w}_{t}e_{j,s}\alpha_{DB}\times \frac{Cy}{ny}, & \text{if}\  R - Cy <= s < R  \\
+        0, & \text{if}\ s \geq R \\
+      \end{cases}
+  ```
+
+#### Stationarized Points System Equations
+
+Stationarized pension amount:
+
+  ```{math}
+  :label: eqn:ps_amount_stationarized
+  \hat{\theta}_{j,u,t+u-s} =\sum_{s=E}^{R-1}\frac{\hat{w}_{t}}{e^{g_{y}(u-s)}}e_{j,s}n_{j,s,t}v_{t}, \ \ \forall u \geq R
+  ```
+
+Stationarized pension amount derivative:
+
+  ```{math}
+  :label: eqn:ps_deriv_stationarized
+  \frac{\partial \hat{\theta}_{j,u,t+u-s}}{\partial n_{j,s,t}} =
+    \begin{cases}
+      \hat{w}_{t}e_{j,s} v_{t}, & \text{if}\ s < R  \\
+      0, & \text{if}\ s \geq R \\
+    \end{cases}
+  ```
+
 (SecStnrzMC)=
 ## Stationarized Market Clearing Equations
 

environment.yml

@@ -8,6 +8,7 @@ dependencies:
 - setuptools
 - scipy>=1.7.1
 - pandas>=1.2.5
+- numba
 - matplotlib
 - dask>=2.30.0
 - dask-core>=2.30.0