Expectated vs Realized Income Growth in A Standard Life Cycle Model

This notebook uses the income process in Cocco, Gomes & Maenhout (2005) to demonstrate that estimates of a regression of expected income changes on realized income changes are sensitive to the size of transitory shocks.

We first load some tools from the HARK toolkit.

import statsmodels.api as sm
from linearmodels.panel.model import PanelOLS
from HARK.distributions import expected
from HARK.ConsumptionSaving.ConsIndShockModel import (
    IndShockConsumerType,
    init_lifecycle,
)

from HARK.Calibration.Income.IncomeTools import (
    parse_income_spec,
    parse_time_params,
    CGM_income,
)

from HARK.Calibration import parse_ssa_life_table
import pandas as pd
from copy import copy

We now create a population of agents with the income process of Cocco, Gomes & Maenhout (2005), which is implemented as a default calibration in the toolkit.

birth_age = 21
death_age = 66
adjust_infl_to = 1992
income_calib = CGM_income
education = "HS"

# Income specification
income_params = parse_income_spec(
    age_min=birth_age,
    age_max=death_age,
    adjust_infl_to=adjust_infl_to,
    **income_calib[education],
    SabelhausSong=True,
)

# We need survival probabilities only up to death_age-1, because survival
# probability at death_age is 1.
liv_prb = parse_ssa_life_table(
    female=True, cross_sec=True, year=2004, age_min=birth_age, age_max=death_age
)

# Parameters related to the number of periods implied by the calibration
time_params = parse_time_params(age_birth=birth_age, age_death=death_age)

# Update all the new parameters
params = copy(init_lifecycle)
params.update(time_params)
# params.update(dist_params)
params.update(income_params)
params.update(
    {
        "LivPrb": liv_prb,
        "pLogInitStd": 0.0,
        "PermGroFacAgg": 1.0,
        "UnempPrb": 0.0,
        "UnempPrbRet": 0.0,
        "track_vars": ["pLvl", "t_age", "PermShk", "TranShk"],
        "AgentCount": 200,
        "T_sim": 500,
    }
)

Agent = IndShockConsumerType(**params)
Agent.solve()

We simulate a population of agents

# Run the simulations
Agent.initialize_sim()
Agent.simulate();

We assume a standard income process with transitory and permanent shocks: The consumer’s Permanent noncapital income $P$ grows by a predictable factor $Γ$ and is subject to an unpredictable multiplicative shock $𝔼_{t} [ψ_{t + 1}] = 1$ ,

\begin{matrix} P_{t + 1} & = P_{t} Γ_{t + 1} ψ_{t + 1}, \end{matrix}

(1)

and, if the consumer is employed, actual income $Y$ is permanent income multiplied by a transitory shock $𝔼_{t} [θ_{t + 1}] = 1$ ,

\begin{matrix} Y_{t + 1} & = P_{t + 1} θ_{t + 1}, \end{matrix}

(2)

$Γ_{t}$ captures the predictable life cycle profile of income growth (faster when young, slower when old). See our replication of CGM-2005 for a detailed account of how these objects map to CGM’s notation.

Now define $y = \log Y, p = \log P$ and similarly for other variables.

Using this notation, we construct all the necessary inputs to the regressors. The main input is the expected income growth of every agent at every time period, which is given by

\begin{matrix} 𝔼_{t} [Y_{t + 1} / Y_{t}] & = 𝔼_{t} [(\frac{θ_{t + 1} P_{t} Γ_{t + 1} ψ_{t + 1}}{θ_{t} P_{t}})] \\ = (\frac{Γ_{t + 1}}{θ_{t}}) \\ 𝔼_{t} [y_{t + 1} - y_{t}] & = \log Γ_{t + 1} - \log θ_{t} \end{matrix}

(3)

exp = [
    expected(dstn=Agent.IncShkDstn[i], func=lambda x: x["PermShk"] * x["TranShk"])
    for i in range(Agent.T_cycle)
]
exp_df = pd.DataFrame(
    {
        "exp_prod": exp,
        "PermGroFac": Agent.PermGroFac,
        "Age": [x + birth_age for x in range(Agent.T_cycle)],
    }
)

raw_data = {
    "Age": Agent.history["t_age"].T.flatten() + birth_age - 1,
    "pLvl": Agent.history["pLvl"].T.flatten(),
    "PermShk": Agent.history["PermShk"].T.flatten(),
    "TranShk": Agent.history["TranShk"].T.flatten(),
}

Data = pd.DataFrame(raw_data)

# Create an individual id
Data["id"] = (Data["Age"].diff(1) < 0).cumsum()

Data["Y"] = Data.pLvl * Data.TranShk

# Find Et[Yt+1 - Yt]
Data = Data.join(exp_df.set_index("Age"), on="Age", how="left")
Data["ExpIncChange"] = Data["pLvl"] * (
    Data["PermGroFac"] * Data["exp_prod"] - Data["TranShk"]
)

Data["Y_change"] = Data.groupby("id")["Y"].diff(1)

A corresponding version of this relationship can be estimated in simulated data:

𝔼_{t} [Δ y_{i, t + 1}] = γ_{0} + γ_{1} Δ y_{i, t} + f_{i} + ϵ_{i, t}

We now estimate an analogous regression in our simulated population.

Data = Data.set_index(["id", "Age"])

# Create the variables they actually use
Data["ExpBin"] = 0
Data.loc[Data["ExpIncChange"] > 0, "ExpBin"] = 1
Data.loc[Data["ExpIncChange"] < 0, "ExpBin"] = -1

Data["ChangeBin"] = 0
Data.loc[Data["Y_change"] > 0, "ChangeBin"] = 1
Data.loc[Data["Y_change"] < 0, "ChangeBin"] = -1

mod = PanelOLS(Data.ExpBin, sm.add_constant(Data.ChangeBin), entity_effects=True)
fe_res = mod.fit()
print(fe_res)

                          PanelOLS Estimation Summary                           
================================================================================
Dep. Variable:                 ExpBin   R-squared:                        0.1958
Estimator:                   PanelOLS   R-squared (Between):             -0.0938
No. Observations:              100000   R-squared (Within):               0.1958
Date:                Thu, Apr 30 2026   R-squared (Overall):              0.1888
Time:                        21:39:45   Log-likelihood                -1.288e+05
Cov. Estimator:            Unadjusted                                           
                                        F-statistic:                   2.376e+04
Entities:                        2435   P-value                           0.0000
Avg Obs:                       41.068   Distribution:                 F(1,97564)
Min Obs:                       2.0000                                           
Max Obs:                       46.000   F-statistic (robust):          2.376e+04
                                        P-value                           0.0000
Time periods:                      45   Distribution:                 F(1,97564)
Avg Obs:                       2222.2                                           
Min Obs:                       1924.0                                           
Max Obs:                       2437.0                                           
                                                                                
                             Parameter Estimates                              
==============================================================================
            Parameter  Std. Err.     T-stat    P-value    Lower CI    Upper CI
------------------------------------------------------------------------------
const          0.1335     0.0028     47.527     0.0000      0.1280      0.1390
ChangeBin     -0.4408     0.0029    -154.14     0.0000     -0.4464     -0.4352
==============================================================================

F-test for Poolability: 1.3343
P-value: 0.0000
Distribution: F(2434,97564)

Included effects: Entity

The estimated ${\hat{γ}}_{1}$ is negative because in usual life-cycle calibrations, transitory shocks are volatile enough that mean reversion of transitory fluctuations is a stronger force than persistent trends in income age-profiles.

However, with less volatile transitory shocks, the regression coefficient would be positive. We demonstrate this by shutting off transitory shocks, simulating another population of agents, and re-running the regression.

params_no_transitory = copy(params)
params_no_transitory.update({"TranShkStd": [0.0] * len(params["TranShkStd"])})

# Create agent
Agent_nt = IndShockConsumerType(**params_no_transitory)
Agent_nt.solve()
# Run the simulations
Agent_nt.initialize_sim()
Agent_nt.simulate();

exp = [
    expected(dstn=Agent_nt.IncShkDstn[i], func=lambda x: x["PermShk"] * x["TranShk"])
    for i in range(Agent_nt.T_cycle)
]
exp_df = pd.DataFrame(
    {
        "exp_prod": exp,
        "PermGroFac": Agent_nt.PermGroFac,
        "Age": [x + birth_age for x in range(Agent.T_cycle)],
    }
)

raw_data = {
    "Age": Agent_nt.history["t_age"].T.flatten() + birth_age - 1,
    "pLvl": Agent_nt.history["pLvl"].T.flatten(),
    "PermShk": Agent_nt.history["PermShk"].T.flatten(),
    "TranShk": Agent_nt.history["TranShk"].T.flatten(),
}

Data = pd.DataFrame(raw_data)

# Create an individual id
Data["id"] = (Data["Age"].diff(1) < 0).cumsum()

Data["Y"] = Data.pLvl * Data.TranShk

# Find Et[Yt+1 - Yt]
Data = Data.join(exp_df.set_index("Age"), on="Age", how="left")
Data["ExpIncChange"] = Data["pLvl"] * (
    Data["PermGroFac"] * Data["exp_prod"] - Data["TranShk"]
)

Data["Y_change"] = Data.groupby("id")["Y"].diff(1)

# Create variables
Data["ExpBin"] = 0
Data.loc[Data["ExpIncChange"] > 0, "ExpBin"] = 1
Data.loc[Data["ExpIncChange"] < 0, "ExpBin"] = -1

Data["ChangeBin"] = 0
Data.loc[Data["Y_change"] > 0, "ChangeBin"] = 1
Data.loc[Data["Y_change"] < 0, "ChangeBin"] = -1

Data = Data.set_index(["id", "Age"])
mod = PanelOLS(Data.ExpBin, sm.add_constant(Data.ChangeBin), entity_effects=True)
fe_res = mod.fit()
print(fe_res)

                          PanelOLS Estimation Summary                           
================================================================================
Dep. Variable:                 ExpBin   R-squared:                        0.0084
Estimator:                   PanelOLS   R-squared (Between):             -0.0257
No. Observations:              100000   R-squared (Within):               0.0084
Date:                Thu, Apr 30 2026   R-squared (Overall):              0.0088
Time:                        21:39:49   Log-likelihood                -1.393e+05
Cov. Estimator:            Unadjusted                                           
                                        F-statistic:                      826.23
Entities:                        2435   P-value                           0.0000
Avg Obs:                       41.068   Distribution:                 F(1,97564)
Min Obs:                       2.0000                                           
Max Obs:                       46.000   F-statistic (robust):             826.23
                                        P-value                           0.0000
Time periods:                      45   Distribution:                 F(1,97564)
Avg Obs:                       2222.2                                           
Min Obs:                       1924.0                                           
Max Obs:                       2437.0                                           
                                                                                
                             Parameter Estimates                              
==============================================================================
            Parameter  Std. Err.     T-stat    P-value    Lower CI    Upper CI
------------------------------------------------------------------------------
const          0.1188     0.0031     38.061     0.0000      0.1127      0.1249
ChangeBin      0.0920     0.0032     28.744     0.0000      0.0857      0.0982
==============================================================================

F-test for Poolability: 1.2509
P-value: 0.0000
Distribution: F(2434,97564)

Included effects: Entity

The estimated ${\hat{γ}}_{1}$ when there are no transitory shocks is positive.