double spline working ok

iot/example.py

@@ -1,5 +1,5 @@
 # %%
-from iot_user import iot_comparison
+from iot.iot_user import iot_comparison
 
 # %%
 iot1 = iot_comparison(
@@ -10,17 +10,20 @@
     baseline_policies=[None, None],
     labels=["2017 Law", "Biden 2020"],
     years=[2017, 2020],
+    mtr_smoother="spline",
 )
 
-# iot2 = iot_comparison(
-#     policies=[
-#         "https://raw.githubusercontent.com/PSLmodels/examples/main/psl_examples/taxcalc/2017_law.json",
-#         "https://raw.githubusercontent.com/PSLmodels/examples/main/psl_examples/taxcalc/Biden2020.json",
-#     ],
-#     labels=["2017 Law", "Biden 2020"],
-#     years=[2017, 2020],
-#     inc_elast=2,
-# )
+iot2 = iot_comparison(
+    policies=[
+        "https://raw.githubusercontent.com/PSLmodels/examples/main/psl_examples/taxcalc/2017_law.json",
+        "https://raw.githubusercontent.com/PSLmodels/examples/main/psl_examples/taxcalc/Biden2020.json",
+    ],
+    baseline_policies=[None, None],
+    labels=["2017 Law", "Biden 2020"],
+    years=[2017, 2020],
+    inc_elast=2,
+    mtr_smoother="poly",
+)
 
 # %%
 iot1.iot[-1].df().head()
@@ -38,3 +41,5 @@
 thetaplot1.show()
 gzplot1.show()
 gzplot2.show()
+
+# %%

iot/inverse_optimal_tax.py

@@ -4,6 +4,9 @@
 from scipy.interpolate import UnivariateSpline
 from statsmodels.nonparametric.kernel_regression import KernelReg
 from sklearn.kernel_ridge import KernelRidge
+from sklearn.linear_model import LinearRegression
+from sklearn.preprocessing import PolynomialFeatures
+import plotly.express as px
 
 
 class IOT:
@@ -46,7 +49,7 @@ def __init__(
         dist_type="kde_full",
         kde_bw=None,
         mtr_smoother="spline",
-        mtr_smooth_param=4,
+        mtr_smooth_param=3,
     ):
         # keep the original data intact
         self.data_original = data.copy()
@@ -116,26 +119,59 @@ def compute_mtr_dist(
                 * mtr_prime (array_like): rate of change in marginal tax rates
                     for each income bin
         """
-        # get rid of duplicate values of income
-        data = data.groupby(income_measure).mean().reset_index()
+
         if mtr_smoother == "spline":
+            # get rid of duplicate values of income
+            data = data[[income_measure, weight_var, "mtr"]].groupby(income_measure).mean().reset_index()
             spl = UnivariateSpline(
                 data[income_measure],
                 data["mtr"],
                 w=data[weight_var],
                 k=mtr_smooth_param,
                 s=len(data[weight_var]) * 6300
             )
+            mtr = spl(self.z)
+            # try running through another spline to see what happens
+            spl = UnivariateSpline(
+                self.z,
+                mtr,
+                w=None,
+                k=mtr_smooth_param,
+                s=len(mtr)
+            )
+            mtr = spl(self.z)
         elif mtr_smoother == "kr":
+            # get rid of duplicate values of income
+            data = data[[income_measure, weight_var, "mtr"]].groupby(income_measure).mean().reset_index()
             krr = KernelRidge(alpha=1.0)
             krr.fit(
                 data[income_measure].values.reshape(-1, 1),
                 data["mtr"].values,
                 sample_weight=data[weight_var].values,
             )
             mtr = krr.predict(self.z)
+        # elif mtr_smoother == "kreg":
+        #     mtr_function = KernelReg(
+        #         data_group.values,
+        #         self.z,
+        #         var_type="c",
+        #         reg_type="ll",
+        #         bw=[mtr_smooth_param],
+        #         ckertype="gaussian",
+        #         defaults=None,
+        #     )
+        #     mtr, _ = mtr_function.fit(self.z)
+        elif mtr_smoother == "poly":
+            poly = PolynomialFeatures(degree=20, include_bias=False)
+            poly_features = poly.fit_transform(data[income_measure].values.reshape(-1,1))
+            poly_reg_model = LinearRegression()
+            poly_reg_model.fit(poly_features, data["mtr"], data[weight_var])
+            z_poly = poly.fit_transform(self.z.reshape(-1, 1))
+            mtr = poly_reg_model.predict(z_poly)
+            px.line(x=self.z, y=mtr).show()
         else:
-            pass
+            print('Please enter a value mtr_smoother method')
+            assert False
         mtr_prime = np.gradient(mtr, edge_order=2)
 
         return mtr, mtr_prime
@@ -181,8 +217,8 @@ def compute_income_dist(
                 ).values
                 / data[weight_var].sum()
             ).sum()
-            print("mu = ", mu)
-            print("sigma = ", np.sqrt(sigmasq))
+            # print("mu = ", mu)
+            # print("sigma = ", np.sqrt(sigmasq))
             f = st.lognorm.pdf(z_line, s=(sigmasq) ** 0.5, scale=np.exp(mu))
         elif dist_type == "kde_full":
             # uses the original full data for kde estimation
@@ -201,6 +237,7 @@ def compute_income_dist(
             )
             f = f_function.pdf(z_line)
         else:
+            print('Please enter a valid value for dist_type')
             assert False
         # normalize f
         f = f / np.sum(f)

iot/iot_user.py

@@ -53,7 +53,7 @@ def __init__(
         data="CPS",
         compare_default=True,
         mtr_wrt="e00200p",
-        income_measure="expanded_income",
+        income_measure="e00200",
         weight_var="s006",
         inc_elast=0.25,
         bandwidth=1000,