moderndid.npiv#

moderndid.npiv(data=None, yname=None, xname=None, wname=None, y=None, x=None, w=None, x_eval=None, x_grid=None, alpha=0.05, basis='tensor', biters=99, j_x_degree=3, j_x_segments=None, k_w_degree=4, k_w_segments=None, k_w_smooth=2, knots='uniform', ucb_h=True, ucb_deriv=True, deriv_index=1, deriv_order=1, check_is_fullrank=False, w_min=None, w_max=None, x_min=None, x_max=None, seed=None)[source]#

Estimate nonparametric instrumental variables model with uniform confidence bands.

Estimates the structural function \(h_0\) and its derivatives in the nonparametric IV model

\[\mathbb{E}[Y - h_0(X) \mid W] = 0 \quad \text{(a.s.)}\]

where \(Y\) is a scalar outcome, \(X\) is a (possibly endogenous) regressor vector, and \(W\) is a vector of instrumental variables. The function is approximated by a B-spline sieve \(h_0(x) \approx (\psi^J(x))' c_J\) and coefficients are estimated by two-stage least squares using \(K\) B-spline basis functions of \(W\) as instruments

\[\hat{c}_J = (\boldsymbol{\Psi}_J' \mathbf{P}_K \boldsymbol{\Psi}_J)^{-} \boldsymbol{\Psi}_J' \mathbf{P}_K \mathbf{Y},\]

where \(\mathbf{P}_K = \mathbf{B}_K (\mathbf{B}_K' \mathbf{B}_K)^{-} \mathbf{B}_K'\) projects onto the instrument space. Function and derivative estimates are then given by

\[\hat{h}_J(x) = (\psi^J(x))' \hat{c}_J, \quad \partial^a \hat{h}_J(x) = (\partial^a \psi^J(x))' \hat{c}_J.\]

When j_x_segments is None, a bootstrap implementation of Lepski’s method selects the sieve dimension \(\tilde{J}\) that adapts to the unknown smoothness of \(h_0\) and instrument strength, achieving the minimax sup-norm convergence rate for both \(h_0\) and its derivatives.

The adaptive CCK procedure then constructs honest uniform confidence bands that guarantee coverage uniformly over a class of data-generating processes. When a fixed j_x_segments is supplied, the standard undersmoothing approach of [1] is used instead.

Parameters:
dataDataFrame, optional

Input data. Accepts any object implementing the Arrow PyCapsule Interface (__arrow_c_stream__), including polars, pandas, pyarrow Table, and cudf DataFrames. When provided, yname, xname, and wname are required.

ynamestr, optional

Name of the outcome column in data.

xnamestr or list of str, optional

Name(s) of the endogenous regressor column(s) in data.

wnamestr or list of str, optional

Name(s) of the instrumental variable column(s) in data.

ynumpy.ndarray of shape (n,), optional

Outcome variable. Required when data is not provided.

xnumpy.ndarray of shape (n,) or (n, p_x), optional

Endogenous regressors. Automatically promoted to 2-d if needed. Required when data is not provided.

wnumpy.ndarray of shape (n,) or (n, p_w), optional

Instrumental variables. Requires \(K \geq J\). Required when data is not provided.

x_evalnumpy.ndarray of shape (m, p_x), optional

Points at which to evaluate \(\hat{h}\) and its derivatives. If None, evaluates at the sample points x.

x_gridnumpy.ndarray, optional

Alias for x_eval. Ignored when x_eval is provided.

alphafloat, default=0.05

Significance level for \(100(1-\alpha)\%\) confidence bands.

basis{“tensor”, “additive”, “glp”}, default=”tensor”

Multivariate basis construction for \(X\):

  • "tensor": Full tensor product of univariate B-splines.

  • "additive": Sum of univariate B-splines (additive model).

  • "glp": Generalized linear product (hierarchical interactions).

bitersint, default=99

Number of multiplier bootstrap draws for critical value computation. Each draw generates i.i.d. \(N(0,1)\) weights \((\varpi_i)_{i=1}^n\) to form bootstrap sup-\(t\) statistics.

j_x_degreeint, default=3

Degree of B-spline basis for \(X\) (order \(r = \text{degree} + 1\)). For UCBs of first derivatives, degree \(\geq 2\) is required; for second derivatives, \(\geq 3\).

j_x_segmentsint, optional

Number of segments for the \(X\) basis, determining sieve dimension \(J\). When None, the data-driven Lepski procedure selects \(\tilde{J}\) adaptively. Supplying a fixed value triggers the undersmoothing UCB approach.

k_w_degreeint, default=4

Degree of B-spline basis for \(W\). Defaults to j_x_degree + 1 because the reduced form \(\mathbb{E}[h_0(X) \mid W]\) is smoother than \(h_0\).

k_w_segmentsint, optional

Number of segments for the instrument basis. When None, chosen proportionally to j_x_segments via the resolution-level mapping \(l_w = \lceil (l + q) \, d / d_w \rceil\), where \(q\) is controlled by k_w_smooth.

k_w_smoothint, default=2

Controls the resolution gap \(q\) between the \(X\) and \(W\) bases in the data-driven procedure. Larger values yield more instrument basis functions relative to the \(X\) basis.

knots{“uniform”, “quantiles”}, default=”uniform”

Knot placement strategy:

  • "uniform": Equally spaced knots on the support.

  • "quantiles": Knots at empirical quantiles of the data.

ucb_hbool, default=True

Compute uniform confidence bands for \(\hat{h}\).

ucb_derivbool, default=True

Compute uniform confidence bands for \(\partial^a \hat{h}\).

deriv_indexint, default=1

Which component of \(X\) to differentiate with respect to (1-based indexing).

deriv_orderint, default=1

Order \(|a|\) of the derivative (1 = first, 2 = second, etc.).

check_is_fullrankbool, default=False

Verify that the basis matrices \(\boldsymbol{\Psi}_J\) and \(\mathbf{B}_K\) have full column rank before estimation.

w_min, w_maxfloat, optional

Override the support bounds for \(W\). Defaults to data range.

x_min, x_maxfloat, optional

Override the support bounds for \(X\). Defaults to data range.

seedint, optional

Random seed for bootstrap reproducibility.

Returns:
NPIVResult

Named tuple with the following fields:

  • h – Estimated \(\hat{h}_J(x)\) at evaluation points.

  • deriv – Estimated \(\partial^a \hat{h}_J(x)\).

  • h_lower, h_upper – Lower/upper UCB for \(h_0\).

  • h_lower_deriv, h_upper_deriv – Lower/upper UCB for \(\partial^a h_0\).

  • beta – Sieve coefficient vector \(\hat{c}_J\).

  • asy_se – Pointwise asymptotic standard errors \(\hat{\sigma}_J(x)\).

  • deriv_asy_se – Pointwise asymptotic standard errors \(\hat{\sigma}_J^a(x)\) for derivatives.

  • cv, cv_deriv – Bootstrap critical values \(z_{1-\alpha}^*\) used for band construction.

  • residuals – TSLS residuals \(\hat{u}_{i,J} = Y_i - \hat{h}_J(X_i)\).

  • j_x_degree, j_x_segments – Basis parameters for \(X\) (segments may differ from input when data-driven).

  • k_w_degree, k_w_segments – Basis parameters for \(W\).

  • args – Diagnostic dictionary. When data-driven selection is used, includes j_x_seg, k_w_seg, j_hat_max, theta_star, and other selection diagnostics.

See also

npiv_est

Core sieve TSLS estimation (no confidence bands).

compute_ucb

Multiplier bootstrap confidence band construction.

npiv_choose_j

Data-driven sieve dimension selection.

References

[1]

Chen, X., & Christensen, T. M. (2018). Optimal sup-norm rates and uniform inference on nonlinear functionals of nonparametric IV regression. Quantitative Economics, 9(1), 39-84.

[2]

Chen, X., Christensen, T. M., & Kankanala, S. (2024). Adaptive estimation and uniform confidence bands for nonparametric structural functions and elasticities. Review of Economic Studies. https://arxiv.org/abs/2107.11869.

[3]

Newey, W. K., & Powell, J. L. (2003). Instrumental variable estimation of nonparametric models. Econometrica, 71(5), 1565-1578.

Examples

The Engel dataset contains household expenditure shares and income measures for 1655 households. We estimate a nonparametric Engel curve relating food share (food) to log-expenditure (logexp), using log-wages (logwages) as an instrument for potentially endogenous expenditure:

In [1]: import numpy as np
   ...: from moderndid import npiv, load_engel
   ...: 
   ...: df = load_engel()
   ...: df.head()
   ...: 
Out[1]: 
shape: (5, 10)
┌──────────┬──────────┬──────────┬──────────┬───┬──────────┬──────────┬──────────┬───────┐
│ food     ┆ catering ┆ alcohol  ┆ fuel     ┆ … ┆ leisure  ┆ logexp   ┆ logwages ┆ nkids │
│ ---      ┆ ---      ┆ ---      ┆ ---      ┆   ┆ ---      ┆ ---      ┆ ---      ┆ ---   │
│ f64      ┆ f64      ┆ f64      ┆ f64      ┆   ┆ f64      ┆ f64      ┆ f64      ┆ i64   │
╞══════════╪══════════╪══════════╪══════════╪═══╪══════════╪══════════╪══════════╪═══════╡
│ 0.28026  ┆ 0.0      ┆ 0.0      ┆ 0.0      ┆ … ┆ 0.202545 ┆ 3.609024 ┆ 5.013565 ┆ 0     │
│ 0.379358 ┆ 0.092832 ┆ 0.0      ┆ 0.190952 ┆ … ┆ 0.033882 ┆ 3.933002 ┆ 2.71866  ┆ 0     │
│ 0.226277 ┆ 0.01666  ┆ 0.133963 ┆ 0.143409 ┆ … ┆ 0.119021 ┆ 4.064315 ┆ 3.881564 ┆ 0     │
│ 0.167698 ┆ 0.074122 ┆ 0.0      ┆ 0.245357 ┆ … ┆ 0.107083 ┆ 4.130275 ┆ 4.900374 ┆ 0     │
│ 0.343115 ┆ 0.042245 ┆ 0.0      ┆ 0.166295 ┆ … ┆ 0.026279 ┆ 4.259548 ┆ 5.564099 ┆ 0     │
└──────────┴──────────┴──────────┴──────────┴───┴──────────┴──────────┴──────────┴───────┘

Estimate the structural function with 5 B-spline segments and 95% uniform confidence bands. The output is an NPIVResult containing the estimated function, confidence bands, derivatives, and diagnostics:

In [2]: result = npiv(
   ...:     data=df,
   ...:     yname="food",
   ...:     xname="logexp",
   ...:     wname="logwages",
   ...:     j_x_segments=5,
   ...:     biters=500,
   ...:     seed=42,
   ...: )
   ...: print(f"Estimates at {len(result.h)} points")
   ...: print(f"h[:5] = {result.h[:5]}")
   ...: print(f"95% UCB critical value: {result.cv:.3f}")
   ...: print(f"Basis: degree={result.j_x_degree}, segments={result.j_x_segments}")
   ...: 
Estimates at 1655 points
h[:5] = [ 9.2505101   0.85939871 -0.26872492 -0.51409879 -0.58251832]
95% UCB critical value: 2.893
Basis: degree=3, segments=5

The derivative of the Engel curve (the marginal propensity to spend on food) is estimated simultaneously:

In [3]: print(f"deriv[:5] = {result.deriv[:5]}")
   ...: print(f"Derivative UCB critical value: {result.cv_deriv:.3f}")
   ...: 
deriv[:5] = [-41.98831557 -12.49571971  -5.12728049  -2.42403742   0.93826829]
Derivative UCB critical value: 2.858

You can also pass numpy arrays directly instead of a DataFrame:

In [4]: rng = np.random.default_rng(0)
   ...: n = 200
   ...: w = rng.uniform(0, 1, (n, 1))
   ...: x = w + 0.2 * rng.normal(0, 1, (n, 1))
   ...: y = np.sin(2 * np.pi * x).ravel() + 0.1 * rng.normal(0, 1, n)
   ...: result = npiv(y=y, x=x, w=w, j_x_segments=4, biters=500, seed=0)
   ...: result.h.shape
   ...: 
Out[4]: (200,)