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

Solve ODE for ETI beliefs #36

Merged
merged 5 commits into from
Feb 25, 2025
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
182 changes: 144 additions & 38 deletions examples/Simulate_all_policies.ipynb
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
},
{
"cell_type": "code",
"execution_count": null,
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
Expand Down Expand Up @@ -617,7 +617,109 @@
},
{
"cell_type": "code",
"execution_count": 16,
"execution_count": 17,
"metadata": {},
"outputs": [],
"source": [
"eti_dict = {\n",
" \"eti_values\": [0.18, 0.106, 0.567, 1.83, 1.9],\n",
" \"knot_points\": [30000, 75000, 250000, 2000000, 10000000]\n",
"}\n",
"# ETI values from Gruber and Saez (2002) (Table 3) and Saez (2004) (Tables 2, 4, 5)\n",
"# Compute MTR schedule under current law\n",
"iot_2023 = iot_user.iot_comparison(\n",
" policies=[{}],\n",
" baseline_policies=[None],\n",
" labels=[\"2023 Law\"],\n",
" years=[2023],\n",
" data=\"CPS\",\n",
" eti=eti_dict\n",
" )\n",
"fig = px.line(\n",
" x=iot_2023.iot[0].df().z,\n",
" y=iot_2023.iot[0].df().mtr\n",
" )\n",
"fig.update_layout(\n",
" template=template,\n",
" xaxis_title=\"Wages and Salaries\",\n",
" yaxis_title=r\"$T'(z)$\",\n",
")\n",
"fig.write_image(\n",
" os.path.join(path, \"MTR_2023.png\"),\n",
" scale=4\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"# Thought experiment: What beliefs about ETI do the candidates need to \n",
"# justify their policies if we assume they are utilitarian? \n",
"# If the elasticities are wildly counterfactual, we can reject this hypothesis\n",
"\n",
"eti_utilitarian = np.zeros((len(iot_all.iot), len(iot_all.iot[0].df().z)))\n",
"for i in range(len(iot_all.iot)):\n",
" eti_utilitarian[i, :] = iot.inverse_optimal_tax.find_eti(\n",
" iot_all.iot[i], g_z = np.ones(len(iot_all.iot[i].df().z))\n",
" )\n",
"\n",
"# Convert the 2D array to a DataFrame for plotting\n",
"eti_df = pd.DataFrame(eti_utilitarian.T, columns=labels)\n",
"\n",
"fig_eti = px.line(\n",
" eti_df,\n",
" x=iot_all.iot[0].df().z,\n",
" y=eti_df.columns,\n",
" labels={\"x\": \"Wages and Salaries\", \"value\": r\"$\\varepsilon$\", \"variable\": \"Candidate\"},\n",
")\n",
"fig_eti.update_layout(\n",
" template=template,\n",
")\n",
"\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"dot\", color=\"blue\"),\n",
" selector=dict(name=\"Obama 2015\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"dot\", color=\"red\"),\n",
" selector=dict(name=\"Romney 2012\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"dash\", color=\"blue\"),\n",
" selector=dict(name=\"Clinton 2016\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"dash\", color=\"red\"),\n",
" selector=dict(name=\"Trump 2016\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"dashdot\", color=\"blue\"),\n",
" selector=dict(name=\"Biden 2020\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"dashdot\", color=\"red\"),\n",
" selector=dict(name=\"Trump 2020\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"solid\", color=\"blue\"),\n",
" selector=dict(name=\"Harris 2024\")\n",
")\n",
"fig_eti.update_traces(\n",
" line=dict(dash=\"solid\", color=\"red\"),\n",
" selector=dict(name=\"Trump 2024\")\n",
")\n",
"fig_eti.write_image(\n",
" os.path.join(path, \"eti_utilitarian.png\"),\n",
" scale=4\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
Expand All @@ -626,6 +728,7 @@
"# one plot with epsilon(z) for each candidate\n",
"# Will pick Trump and Clinton (2016) for example\n",
"\n",
"\n",
"# First, plot just their g(z)\n",
"fig = px.line(\n",
" x=iot_all.iot[2].df().z[10:],\n",
Expand All @@ -649,56 +752,59 @@
" os.path.join(path, \"trump_clinton_g_z_numerical.png\"),\n",
" scale=4\n",
" )\n",
"# Now find the epsilon(z) that would give Trump's policies the same g(z) as Clinton\n",
"eti_beliefs_lw, eti_beliefs_jjz = iot.inverse_optimal_tax.find_eti(iot_all.iot[2], iot_all.iot[3], g_z_type=\"g_z\")\n",
"idx = np.where(np.absolute(eti_beliefs_jjz[1:]) < 10)[0]\n",
"fig2 = px.line(\n",
" x=iot_all.iot[2].df().z[idx],\n",
" y=eti_beliefs_jjz[idx],\n",
" labels={\"x\": \"Wages and Salaries\", \"y\": r\"$\\text{Implied } \\varepsilon$\"},\n",
"\n",
"\n",
"# What ETI does Clinton need to justify Trump's SWW given her policy?\n",
"eti_clinton = iot.inverse_optimal_tax.find_eti(\n",
" iot_all.iot[2], \n",
" g_z=iot_all.iot[3].df().g_z)\n",
"\n",
"fig_eti_clinton = px.line(\n",
" x=iot_all.iot[2].df().z,\n",
" y=eti_clinton,\n",
" labels={\"x\": \"Wages and Salaries\", \"y\": r\"$\\varepsilon$\"},\n",
" )\n",
"fig2.update_layout(\n",
"fig_eti_clinton.update_layout(\n",
" template=template,\n",
")\n",
"fig2.write_image(\n",
" os.path.join(path, \"trump_eti.png\"),\n",
"fig_eti_clinton.update_traces(\n",
" line=dict(dash=\"dash\", color=\"blue\"),\n",
" selector=dict(name=\"Clinton 2016\")\n",
")\n",
"fig_eti_clinton.write_image(\n",
" os.path.join(path, \"eti_clinton.png\"),\n",
" scale=4\n",
" )"
" )\n"
]
},
{
"cell_type": "code",
"execution_count": 17,
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"eti_dict = {\n",
" \"eti_values\": [0.18, 0.106, 0.567, 1.83, 1.9],\n",
" \"knot_points\": [30000, 75000, 250000, 2000000, 10000000]\n",
"}\n",
"# ETI values from Gruber and Saez (2002) (Table 3) and Saez (2004) (Tables 2, 4, 5)\n",
"# Compute MTR schedule under current law\n",
"iot_2023 = iot_user.iot_comparison(\n",
" policies=[{}],\n",
" baseline_policies=[None],\n",
" labels=[\"2023 Law\"],\n",
" years=[2023],\n",
" data=\"CPS\",\n",
" eti=eti_dict\n",
" )\n",
"fig = px.line(\n",
" x=iot_2023.iot[0].df().z,\n",
" y=iot_2023.iot[0].df().mtr\n",
"\n",
"# Converse: What elasticty does Trump need to justify Clinton's weights?\n",
"eti_trump = iot.inverse_optimal_tax.find_eti(\n",
" iot_all.iot[3], \n",
" g_z=iot_all.iot[2].df().g_z)\n",
"\n",
"fig_eti_trump = px.line(\n",
" x=iot_all.iot[2].df().z,\n",
" y=eti_trump,\n",
" labels={\"x\": \"Wages and Salaries\", \"y\": r\"$\\varepsilon$\"},\n",
" )\n",
"fig.update_layout(\n",
"fig_eti_trump.update_layout(\n",
" template=template,\n",
" xaxis_title=\"Wages and Salaries\",\n",
" yaxis_title=r\"$T'(z)$\",\n",
")\n",
"fig.write_image(\n",
" os.path.join(path, \"MTR_2023.png\"),\n",
"fig_eti_trump.update_traces(\n",
" line=dict(dash=\"dash\", color=\"red\"),\n",
" selector=dict(name=\"Trump 2016\")\n",
")\n",
"fig_eti_trump.write_image(\n",
" os.path.join(path, \"eti_trump.png\"),\n",
" scale=4\n",
" )"
" )\n"
]
}
],
Expand All @@ -718,7 +824,7 @@
"name": "python",
"nbconvert_exporter": "python",
"pygments_lexer": "ipython3",
"version": "3.12.3"
"version": "3.11.7"
}
},
"nbformat": 4,
Expand Down
55 changes: 29 additions & 26 deletions iot/inverse_optimal_tax.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,8 @@ class IOT:
parametric, if None, then non-parametric bin weights
mtr_smoother (None or str): method used to smooth our mtr
function, if None, then use bin average mtrs
mtr_smooth_param (scalar): parameter for mtr_smoother
kreg_bw (array_like): bandwidth for kernel regression

Returns:
class instance: IOT
Expand Down Expand Up @@ -122,6 +124,8 @@ def compute_mtr_dist(
weight_var (str): name of weight measure from data to use
mtr_smoother (None or str): method used to smooth our mtr
function, if None, then use bin average mtrs
mtr_smooth_param (scalar): parameter for mtr_smoother
kreg_bw (array_like): bandwidth for kernel regression

Returns:
tuple:
Expand Down Expand Up @@ -207,6 +211,7 @@ def compute_income_dist(
weight_var (str): name of weight measure from data to use
dist_type (None or str): type of distribution to use if
parametric, if None, then non-parametric bin weights
kde_bw (array_like): bandwidth for kernel regression

Returns:
tuple:
Expand Down Expand Up @@ -367,7 +372,7 @@ def sw_weights(self):
+ ((self.theta_z * self.eti * self.mtr) / (1 - self.mtr))
+ ((self.eti * self.z * self.mtr_prime) / (1 - self.mtr) ** 2)
)
integral = np.trapz(g_z * self.f, self.z)
integral = np.trapz(g_z * self.f, self.z) # renormalize to integrate to 1
g_z = g_z / integral

# use Lockwood and Weinzierl formula, which should be equivalent but using numerical differentiation
Expand All @@ -386,41 +391,39 @@ def sw_weights(self):
return g_z, g_z_numerical


def find_eti(iot1, iot2, g_z_type="g_z"):
def find_eti(iot, g_z = None, eti_0 = 0.25):
"""
This function solves for the ETI that would result in the
policy represented via MTRs in iot2 be consistent with the
social welfare function inferred from the policies of iot1.
policy represented via MTRs in IOT being consistent with the
social welfare function supplied. It solves a first order
ordinary differential equation.

.. math::
\varepsilon_{z} = \frac{(1-T'(z))}{T'(z)}\frac{(1-F(z))}{zf(z)}\int_{z}^{\infty}\frac{1-g_{\tilde{z}}{1-F(y)}dF(\tilde{z})
\varepsilon'(z)\left[\frac{zT'(z)}{1-T'(z)}\right] + \varepsilon(z)\left[\theta_z \frac{T'(z)}{1-T'(z)} +\frac{zT''(z)}{(1-T'(z))^2}\right]+ (1-g(z))

Args:
iot1 (IOT): IOT class instance representing baseline policy
iot2 (IOT): IOT class instance representing reform policy
g_z_type (str): type of social welfare function to use
Options are:
* 'g_z' for the analytical formula
* 'g_z_numerical' for the numerical approximation
iot (IOT): instance of the I
g_z (None or array_like): vector of social welfare weights
eti_0 (scalar): guess for ETI at z=0

Returns:
eti_beliefs (array-like): vector of ETI beliefs over z
"""
if g_z_type == "g_z":
g_z = iot1.g_z
else:
g_z = iot1.g_z_numerical
# The equation below is a simplication of the above to make the integration easier
eti_beliefs_lw = ((1 - iot2.mtr) / (iot2.z * iot2.f * iot2.mtr)) * (
1 - iot2.F - (g_z.sum() - np.cumsum(g_z))
)
# derivation from JJZ analytical solution that doesn't involved integration
eti_beliefs_jjz = (g_z - 1) / (
(iot2.theta_z * (iot2.mtr / (1 - iot2.mtr)))
+ (iot2.z * (iot2.mtr_prime / (1 - iot2.mtr) ** 2))
)

return eti_beliefs_lw, eti_beliefs_jjz

if g_z is None:
g_z = iot.g_z

# we solve an ODE of the form f'(z) + P(z)f(z) = Q(z)
P_z = 1/iot.z + iot.f_prime/iot.f + iot.mtr_prime/(iot.mtr * (1-iot.mtr))
# integrating factor for ODE: mu(z) * f'(z) + mu(z) * P(z) * f(z) = mu(z) * Q(z)
mu_z = np.exp(np.cumsum(P_z))
Q_z = (g_z - 1) * (1 - iot.mtr) / (iot.mtr * iot.z)
# integrate Q(z) * mu(z), as we integrate both sides of the ODE
int_mu_Q = np.cumsum(mu_z * Q_z)

eti_beliefs = (eti_0 + int_mu_Q) / mu_z

return eti_beliefs


def wm(value, weight):
Expand Down
Loading