statspai.structural

structural

Structural estimation methods.

BLPResult

Bases: ResultProtocolMixin

Results from BLP demand estimation.

Attributes:

Name	Type	Description
linear_params	Series	Linear parameter estimates (β, α).
nonlinear_params	Series	Nonlinear parameter estimates (σ, random coefficient std devs).
se_linear	Series	Standard errors for linear parameters.
se_nonlinear	Series	Standard errors for nonlinear parameters.
mean_utility	Series	Estimated mean utility δ for each product-market, indexed like the estimation data.
own_elasticities	Series	Own-price elasticities for each product-market, indexed like the estimation data.
ci_alpha	float	One minus the confidence level of the intervals in summary() (sp.blp(alpha=)).
n_markets	int	Number of markets.
n_products	int	Total number of product-market observations.
gmm_objective	float	Value of the GMM objective at the optimum, (Z'xi)' W (Z'xi) with the step-2 weighting matrix -- N times Hansen's J statistic N * gbar' W gbar (the value pyblp reports).
converged	bool	Whether the outer-loop optimization converged.
model	str	'random_coefficients', 'logit' or 'nested_logit'.
rho	float or None	Nesting parameter of the nested logit (also in nonlinear_params).
instruments	list of str	Excluded instruments the moments were built from.
omitted_instruments	list of dict	Instruments dropped as redundant, {'variable', 'reason'} each.
gmm_steps	int	1 (weighting matrix (Z'Z)^{-1}, the 2SLS weighting) or 2.
hansen_j, hansen_df, hansen_p	(float, int, float)	Hansen's test of the overidentifying restrictions after two-step GMM; None after one step or when the model is just identified.
model_info	dict	How the outer loop went (random coefficients only): optimizer, sigma0, starts (each first-step start with the objective it reached), n_starts_at_best, step1_sigma, step1_objective, n_objective_evaluations, projected_gradient_norm (L-BFGS-B) and sigma_at_bound (standard deviations estimated at zero, whose standard error is NaN). With a random coefficient on price also share_positive_price_coef and n_positive_own_elasticity.

Examples:

>>> import statspai as sp
>>> import numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for m in range(10):  # 10 markets, 5 products each
...     x1 = rng.normal(0, 1, 5)
...     x2 = rng.normal(0, 1, 5)
...     price = rng.uniform(1, 3, 5)
...     delta = x1 + 0.5 * x2 - price + rng.normal(0, 0.2, 5)
...     ex = np.exp(delta)
...     share = ex / (1 + ex.sum())
...     for j in range(5):
...         rows.append({"market_id": m, "product_id": j,
...                      "share": share[j], "price": price[j],
...                      "x1": x1[j], "x2": x2[j]})
>>> df = pd.DataFrame(rows)
>>> res = sp.blp(df, shares="share", prices="price",
...              x_linear=["x1", "x2"], x_random=["x1"],
...              n_draws=50, maxiter=50, seed=0)
>>> type(res).__name__
'BLPResult'
>>> (res.n_markets, res.n_products)
(10, 50)
>>> bool("x1" in res.linear_params.index)
True

summary

summary() -> str

Return a formatted summary of BLP estimation results.

elasticity_matrix

elasticity_matrix(market_id: Optional[Any] = None) -> DataFrame

Return the full own- and cross-price elasticity matrix for a market.

Parameters:

Name	Type	Description	Default
market_id	hashable	Market to return. If None, returns the first market.	None

Returns:

Type	Description
DataFrame	(J_m x J_m) elasticity matrix with product labels.

diversion_ratios

diversion_ratios(market_id: Optional[Any] = None, include_outside: bool = False) -> DataFrame

Diversion ratios for a given market.

D[j, k] = (ds_k/dp_j) / (-ds_j/dp_j): of the buyers product j (row) loses when its price rises, the fraction that switches to product k (column). Each row sums to at most one; the rest goes to the outside good.

Parameters:

Name	Type	Description	Default
market_id	hashable	Market to compute for. If None, uses the first market.	None
include_outside	bool	Append a column "outside" with the fraction that leaves the market, 1 - sum_k D[j, k].	False

Returns:

Type	Description
DataFrame	(J_m, J_m) with a zero diagonal, plus the "outside" column when requested. Rows of products whose own-price derivative is not negative are NaN, with a warning.

Notes

Until October 2026 this method divided cross- by own-price elasticities, which omits the share ratio s_k / s_j; rows could sum to well above one.

markups

markups(market_id: Optional[Any] = None, firm_ids: Any = None) -> Series

Price minus marginal cost under multi-product Bertrand pricing.

Each firm sets the prices of its own products given the others' prices; the first-order conditions give p - mc = -(O * dS')^{-1} s, with dS the share Jacobian and O[j, k] = 1 when j and k belong to the same firm [@berry1995automobile]. Demand is the estimated one; no supply side is estimated.

Parameters:

Name	Type	Description	Default
market_id	hashable	One market (labelled by product) or, if None, every product-market row (indexed like the estimation data).	None
firm_ids	array - like or 'single'	Firm of each row of the estimation data. Defaults to the firm_id column given to sp.blp; 'single' treats every product as its own firm.	None

Returns:

Type	Description
Series	In price units.

marginal_costs

marginal_costs(market_id: Optional[Any] = None, firm_ids: Any = None) -> Series

Marginal costs implied by :meth:markups: mc = p - markup.

to_econometric_results

to_econometric_results() -> EconometricResults

Convert to a standard EconometricResults object.

PathAnalysisResult

Bases: ResultProtocolMixin

Fitted path model.

Attributes:

Name	Type	Description
params	DataFrame	One row per parameter, in lavaan's layout: lhs, op (=~ loading, ~ regression, ~~ (co)variance, ~1 intercept, := defined), rhs, label, est, se, z, pvalue, ci_lower, ci_upper and std_all (the estimate with every variable standardised).
fit	dict	chisq, df, pvalue (the test of the model against the saturated one), baseline_chisq, baseline_df, cfi, tli, rmsea with its 90% interval, srmr, logl, aic, bic and npar; with se='robust' also the Satorra-Bentler chisq_scaled, pvalue_scaled and scaling_factor.
implied_cov, sample_cov	DataFrame	Model-implied and sample covariance matrices of the observed variables (divisor n).
implied_mean	Series or None	Model-implied means, when the model has a mean structure.
latent_cov	DataFrame	Model-implied covariance matrix of the latent variables (empty without any).
factor_scores	DataFrame or None	The expected value of each latent variable given a row's observed values (the regression method), for the rows used in the fit.
vcov	DataFrame	Covariance matrix of the free parameters.
r2	Series	Share of each endogenous variable's variance the model explains (for an indicator, its reliability).
n_obs	int	Rows used. n_dropped counts the rest, and n_patterns the distinct patterns of missing values among the rows used (one with complete data).

Examples:

>>> import numpy as np, pandas as pd
>>> import statspai as sp
>>> rng = np.random.default_rng(0)
>>> x = rng.normal(size=500)
>>> m = 0.5 * x + rng.normal(size=500)
>>> y = 0.4 * m + 0.2 * x + rng.normal(size=500)
>>> df = pd.DataFrame({"x": x, "m": m, "y": y})
>>> fit = sp.path_analysis("m ~ a*x\ny ~ b*m + c*x\nind := a*b", df)
>>> type(fit).__name__, int(fit.fit["df"])
('PathAnalysisResult', 0)
>>> list(fit.params["op"].unique())
['~', '~~', ':=']

effect

effect(name: str) -> Series

The row of a defined (:=) or labelled parameter.

ProductionResult

Bases: EconometricResults

Result object for production function estimation.

Inherits params / std_errors / summary / to_dict from :class:EconometricResults and adds production-function-specific payload:

• coef — input elasticities keyed by input name (e.g. {"l": 0.62, "k": 0.31})
• tfp — firm-time TFP estimates omega_it (in logs); same length as the post-stage-2 working sample
• residuals — i.i.d. shock eta_it from stage 1
• productivity_process — {"rho": float, "sigma": float}: the slope of the linear projection of omega on its lag, and the standard deviation of the innovation. The coefficients of the Markov polynomial itself are in diagnostics["markov_coef"]
• markup — placeholder; populated by :func:statspai.markup

Use .summary() for a Stata-style table or .coef for the raw dict.

Examples:

>>> import numpy as np, pandas as pd
>>> import statspai as sp
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for fid in range(30):
...     k = rng.normal(2.0, 0.5)
...     omega = rng.normal(0.0, 0.3)
...     for yr in range(8):
...         omega = 0.7 * omega + rng.normal(0.0, 0.2)
...         k = 0.9 * k + 0.3 * rng.normal(1.0, 0.3)
...         l = 0.6 * k + omega + rng.normal(1.0, 0.2)
...         m = 0.5 * k + 0.5 * l + omega + rng.normal(0.5, 0.2)
...         y = 0.6 * l + 0.3 * k + omega + rng.normal(0.0, 0.1)
...         rows.append(dict(id=fid, year=yr, y=y, l=l, k=k, m=m))
>>> df = pd.DataFrame(rows)
>>> res = sp.levinsohn_petrin(df, output="y", free="l", state="k", proxy="m")
>>> type(res).__name__
'ProductionResult'
>>> sorted(res.coef.keys())
['k', 'l']
>>> res.method
'lp'

cite

cite(format: str = 'bibtex') -> str

Return the canonical reference string for self.method.

blp_instruments

blp_instruments(data: DataFrame, x: list, market_id: str = 'market_id', firm_id: Optional[str] = None, nest: Optional[str] = None, count: bool = True) -> DataFrame

Characteristics-based instruments for differentiated-products demand.

For every characteristic in x (and the product count), the sums over the other products a product competes with in its market: the firm's own other products and its rivals' products when firm_id is given (Berry, Levinsohn & Pakes 1995, section 5.1), all other products in the market otherwise. With nest the sums over the other products of the same nest are added; they move the within-nest share of a nested logit.

Parameters:

Name	Type	Description	Default
data	DataFrame	One row per product and market.	required
x	list of str	Exogenous product characteristics (not the price).	required
market_id	str	Column identifying markets.	'market_id'
firm_id	str	Column identifying the firm that sells each product. With it the market sum is split into own_<x> (the firm's other products) and rival_<x> (other firms' products).	None
nest	str	Column identifying the nest of each product; adds nest_<x>, the sum over the other products of the same nest and market.	None
count	bool	Also build the sums of a constant, i.e. the number of other products (own_count, rival_count, nest_count).	True

Returns:

Type	Description
DataFrame	One column per instrument, indexed like data. The columns can be passed to sp.blp(instruments=...) or written into an sp.iv formula. A column that does not vary within the data (a rival count when every market has the same number of products, say) is still returned; sp.blp drops redundant instruments itself.

Notes

The instruments are valid when a product's unobserved quality is mean independent of the observed characteristics of all products in its market. They say nothing about whether that holds; an overidentification test after estimation is the usual check.

Examples:

>>> import statspai as sp
>>> import pandas as pd
>>> df = pd.DataFrame({
...     "market_id": [1, 1, 1, 2, 2],
...     "firm": ["a", "a", "b", "a", "b"],
...     "size": [1.0, 2.0, 4.0, 3.0, 5.0],
... })
>>> z = sp.blp_instruments(df, ["size"], firm_id="firm")
>>> z["own_size"].tolist()
[2.0, 1.0, 0.0, 0.0, 0.0]
>>> z["rival_size"].tolist()
[4.0, 4.0, 3.0, 5.0, 3.0]
>>> z["rival_count"].tolist()
[1.0, 1.0, 2.0, 1.0, 1.0]

References

[@berry1995automobile]

ackerberg_caves_frazer

ackerberg_caves_frazer(data: DataFrame, output: str = 'y', free: Sequence[str] | str | None = None, state: Sequence[str] | str | None = None, proxy: str = 'm', panel_id: str = 'id', time: str = 'year', polynomial_degree: int = 3, productivity_degree: int = 3, functional_form: str = 'cobb-douglas', boot_reps: int = 0, seed: Optional[int] = None) -> ProductionResult

Ackerberg-Caves-Frazer (2015) production function estimator.

Corrects the OP / LP "functional dependence" identification problem: when free inputs (labor) are chosen at the same time as the proxy, the labor coefficient is not identified in the stage-1 polynomial. ACF moves all coefficient identification to stage 2, instrumenting free inputs with their lagged values and state inputs at the contemporaneous level (prodest's acf option / prodestACF).

The just-identified moment conditions can have several roots. They are searched from the slopes of the stage-1 polynomial at the sample means and a grid of other starts; the estimate is the root nearest those slopes (so it does not change when a log input is measured in other units, i.e. shifted by a constant), every distinct root is listed in diagnostics["acf_roots"], and a RuntimeWarning names them when there is more than one.

Parameters:

Name	Type	Description	Default
data	DataFrame	Long-form panel.	required
output	str	Log output.	'y'
free	str or list	Free inputs — instrumented with their lag in stage 2.	``["l"]``
state	str or list	State inputs.	``["k"]``
proxy	str	Intermediate input (materials).	``"m"``
panel_id	str	Firm and (numeric) year identifiers.	'id'
time	str	Firm and (numeric) year identifiers.	'id'
polynomial_degree	int	Degree of the stage-1 polynomial in (free, state, proxy).	3
productivity_degree	int	Degree of the productivity polynomial g.	3
functional_form	('cobb-douglas', 'translog')	Translog adds 0.5 * x_j**2 and x_j * x_k terms, instrumented by the same expansion of (lagged free, state).	'cobb-douglas'
boot_reps	int	Firm-cluster bootstrap replications; each replication solves from the full-sample estimate and must reach a root.	0
seed	int	Bootstrap RNG seed.	None

Returns:

Type	Description
ProductionResult

Notes

Requires at least two consecutive time periods per firm so that lagged labor exists.

Examples:

>>> import statspai as sp
>>> import numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for fid in range(30):
...     k = rng.normal(2.0, 0.5)
...     omega = rng.normal(0.0, 0.3)
...     for yr in range(8):
...         omega = 0.7 * omega + rng.normal(0.0, 0.2)
...         k = 0.9 * k + 0.3 * rng.normal(1.0, 0.3)
...         l = 0.6 * k + omega + rng.normal(1.0, 0.2)
...         # The proxy is an exact, invertible function of (k, omega):
...         # noise in it would break the scalar-unobservable condition.
...         m = 0.5 + 0.5 * k + omega
...         y = 0.6 * l + 0.3 * k + omega + rng.normal(0.0, 0.1)
...         rows.append(dict(id=fid, year=yr, y=y, l=l, k=k, m=m))
>>> df = pd.DataFrame(rows)
>>> import warnings
>>> with warnings.catch_warnings():
...     warnings.simplefilter("ignore")
...     res = sp.ackerberg_caves_frazer(
...         df, output="y", free="l", state="k", proxy="m"
...     )
>>> sorted(res.coef.keys())
['k', 'l']
>>> res.method
'acf'

References

Ackerberg, D.A., Caves, K. & Frazer, G. (2015). Identification properties of recent production function estimators. Econometrica, 83(6), 2411-2451.

levinsohn_petrin

levinsohn_petrin(data: DataFrame, output: str = 'y', free: Sequence[str] | str | None = None, state: Sequence[str] | str | None = None, proxy: str = 'm', panel_id: str = 'id', time: str = 'year', polynomial_degree: int = 3, productivity_degree: int = 3, functional_form: str = 'cobb-douglas', boot_reps: int = 0, seed: Optional[int] = None) -> ProductionResult

Levinsohn-Petrin (2003) production function estimator.

Uses an intermediate input (materials / energy) as the productivity proxy, which avoids the OP zero-investment selection problem. Same two steps as :func:olley_pakes (prodest's method(lp) / prodestLP).

Parameters:

Name	Type	Description	Default
data	DataFrame	Long-form panel.	required
output	str	Log output.	'y'
free	str or list	Free inputs.	``["l"]``
state	str or list	State inputs.	``["k"]``
proxy	str	Intermediate input (materials).	``"m"``
panel_id	str	Firm and (numeric) year identifiers.	'id'
time	str	Firm and (numeric) year identifiers.	'id'
polynomial_degree	int	Degree of the stage-1 polynomial in (state, proxy).	3
productivity_degree	int	Degree of the productivity polynomial g.	3
functional_form	'cobb-douglas'	Translog raises NotImplementedError.	'cobb-douglas'
boot_reps	int	Firm-cluster bootstrap replications.	0
seed	int	Bootstrap RNG seed.	None

Returns:

Type	Description
ProductionResult

Notes

The output is value added: the proxy does not enter the production function. The gross-output specification of the original paper, where materials are both the proxy and an input, is refused. In data simulated under the paper's own assumptions its capital and materials coefficients are not pinned down by the moment conditions.

productivity_process["rho"] is the slope of the linear projection of productivity on its lag; the coefficients of the Markov polynomial are in diagnostics["markov_coef"].

Examples:

>>> import statspai as sp
>>> import numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for fid in range(30):
...     k = rng.normal(2.0, 0.5)
...     omega = rng.normal(0.0, 0.3)
...     for yr in range(8):
...         omega = 0.7 * omega + rng.normal(0.0, 0.2)
...         k = 0.9 * k + 0.3 * rng.normal(1.0, 0.3)
...         l = 0.6 * k + omega + rng.normal(1.0, 0.2)
...         # The proxy is an exact, invertible function of (k, omega):
...         # noise in it would break the scalar-unobservable condition.
...         m = 0.5 + 0.5 * k + omega
...         y = 0.6 * l + 0.3 * k + omega + rng.normal(0.0, 0.1)
...         rows.append(dict(id=fid, year=yr, y=y, l=l, k=k, m=m))
>>> df = pd.DataFrame(rows)
>>> res = sp.levinsohn_petrin(df, output="y", free="l", state="k", proxy="m")
>>> sorted(res.coef.keys())
['k', 'l']
>>> res.method
'lp'

References

Levinsohn, J. & Petrin, A. (2003). Estimating production functions using inputs to control for unobservables. Review of Economic Studies, 70(2), 317-341.

olley_pakes

olley_pakes(data: DataFrame, output: str = 'y', free: Sequence[str] | str | None = None, state: Sequence[str] | str | None = None, proxy: str = 'i', panel_id: str = 'id', time: str = 'year', polynomial_degree: int = 3, productivity_degree: int = 3, functional_form: str = 'cobb-douglas', boot_reps: int = 0, seed: Optional[int] = None, drop_zero_proxy: bool = True) -> ProductionResult

Olley-Pakes (1996) production function estimator.

Uses investment as the proxy for unobserved productivity. The free-input coefficients come from the stage-1 regression of output on the free inputs and a polynomial in (state, investment); the state coefficients minimise the sum of squared productivity innovations in stage 2. This is the estimator of Stata's and R's prodest (method(op) / prodestOP).

Parameters:

Name	Type	Description	Default
data	DataFrame	Long-form panel: one row per (firm, year).	required
output	str	Log output column.	``"y"``
free	str or list	Freely chosen inputs (e.g. labor). Multiple are allowed.	``["l"]``
state	str or list	State inputs (capital, predetermined).	``["k"]``
proxy	str	Investment column.	``"i"``
panel_id	str	Firm and year identifiers. time must be numeric; the lag of a year is the firm's value one year earlier (NaN across a gap).	'id'
time	str	Firm and year identifiers. time must be numeric; the lag of a year is the firm's value one year earlier (NaN across a gap).	'id'
polynomial_degree	int	Degree of the stage-1 polynomial in (state, proxy); prodest's Stata default is 3, R prodest uses 2.	3
productivity_degree	int	Degree of the polynomial g in the productivity process (cubic in both prodest implementations).	3
functional_form	'cobb-douglas'	Translog raises NotImplementedError; use :func:ackerberg_caves_frazer.	'cobb-douglas'
boot_reps	int	Firm-cluster bootstrap replications. 0 ⇒ NaN standard errors.	0
seed	int	Bootstrap RNG seed.	None
drop_zero_proxy	bool	Drop rows with non-positive proxy (the OP inversion needs strictly positive investment in levels). Pass False when the proxy is in logs. Dropping period t for a firm also removes period t+1 from stage 2, whose lag is then missing. The number of rows dropped is in diagnostics["n_dropped_nonpositive_proxy"]; negative values among them draw a warning.	True

Returns:

Type	Description
ProductionResult

Examples:

>>> import statspai as sp
>>> import numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for fid in range(30):
...     k = rng.normal(2.0, 0.5)
...     omega = rng.normal(0.0, 0.3)
...     for yr in range(8):
...         omega = 0.7 * omega + rng.normal(0.0, 0.2)
...         k = 0.9 * k + 0.3 * rng.normal(1.0, 0.3)
...         l = 0.6 * k + omega + rng.normal(1.0, 0.2)
...         # The proxy is an exact, invertible function of (k, omega):
...         # noise in it would break the scalar-unobservable condition.
...         inv = 1.0 + 0.4 * k + omega  # investment > 0
...         y = 0.6 * l + 0.3 * k + omega + rng.normal(0.0, 0.1)
...         rows.append(dict(id=fid, year=yr, y=y, l=l, k=k, i=inv))
>>> df = pd.DataFrame(rows)
>>> res = sp.olley_pakes(df, output="y", free="l", state="k", proxy="i")
>>> sorted(res.coef.keys())
['k', 'l']
>>> res.method
'op'

References

Olley, G.S. & Pakes, A. (1996). The dynamics of productivity in the telecommunications equipment industry. Econometrica, 64(6), 1263-1297.

prod_fn

prod_fn(data: DataFrame, output: str = 'y', free: Sequence[str] | str | None = None, state: Sequence[str] | str | None = None, proxy: str | None = None, panel_id: str = 'id', time: str = 'year', method: str = 'acf', polynomial_degree: Optional[int] = None, productivity_degree: Optional[int] = None, functional_form: str = 'cobb-douglas', boot_reps: int = 0, seed: Optional[int] = None, **kwargs: Any) -> ProductionResult

Production function estimation — unified interface.

Parameters:

Name	Type	Description	Default
data	DataFrame	Long panel with one row per (firm, year).	required
output	str	Log output column.	'y'
free	str or list	Free inputs (e.g. labor).	``["l"]``
state	str or list	State / predetermined inputs (capital).	``["k"]``
proxy	str	Productivity proxy. Defaults: "i" for method="op", "m" for all others.	None
panel_id	str	Panel identifier columns.	'id'
time	str	Panel identifier columns.	'id'
method	('op', 'lp', 'acf', 'wrdg')	Estimator. ACF is the modern default (corrects OP/LP identification problem).	'op'
polynomial_degree	int	Degree of the control-function polynomial; None uses the estimator's default (3, as Stata prodest).	None
productivity_degree	int	Degree of the productivity polynomial g; None uses the default (3, cubic, as both prodest implementations for OP / LP / ACF).	None
functional_form	('cobb-douglas', 'translog')	Functional form. Translog adds 0.5 * x_j**2 own-quadratic terms and x_j*x_k cross terms — output elasticities then vary by firm-time and ProductionResult.model_info["elasticities"] carries a per-row DataFrame. Translog is available for method="acf" only (OP, LP and Wooldridge raise NotImplementedError). Its stage-2 instruments are the same polynomial transform of (l_lag, k), so the moment system can be near-singular when state and lagged free inputs are highly correlated and the higher-order coefficients (ll, kk, lk) are much noisier than the linear terms. Use bootstrap SEs (boot_reps>=200) to gauge this.	'cobb-douglas'
boot_reps	int	Firm-cluster bootstrap replications. 0 ⇒ NaN standard errors.	0
seed	int		None

Returns:

Type	Description
ProductionResult	.coef for elasticities, .tfp for log-productivity series, .summary() for a Stata-style table.

Examples:

>>> import statspai as sp
>>> import numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for fid in range(60):
...     omega = rng.normal(0.0, 0.28)
...     k = rng.normal(0.0, 0.5)
...     for t in range(8):
...         omega = 0.7 * omega + rng.normal(0.0, 0.20)
...         l = 0.5 * omega + 0.3 * k + rng.normal(0.0, 0.10)
...         m = 0.8 * omega + 0.5 * k + rng.normal(0.0, 0.05)
...         y = 0.60 * l + 0.35 * k + omega + rng.normal(0.0, 0.10)
...         rows.append({"id": fid, "year": t, "y": y, "l": l, "k": k, "m": m})
...         k = 0.9 * k + 0.1 * rng.normal(0.5, 0.1)
>>> df = pd.DataFrame(rows)
>>> res = sp.prod_fn(df, output="y", free="l", state="k", proxy="m",
...                   panel_id="id", time="year", method="acf", seed=0)
>>> sorted(res.coef)  # output elasticities for free + state inputs
['k', 'l']
>>> # De Loecker-Warzynski markup needs the flexible input among the
>>> # production-function coefficients, so refit with materials as free:
>>> res_m = sp.prod_fn(df, output="y", free=["l", "m"], state="k",
...                    proxy="m", panel_id="id", time="year", method="acf")
>>> mu = sp.markup(
...     res_m, revenue="y", input_cost="m", flexible_input="m"
... )

See Also

olley_pakes, levinsohn_petrin, ackerberg_caves_frazer, wooldridge_prod markup : De Loecker-Warzynski (2012) firm-time markup.

References

Olley & Pakes (1996); Levinsohn & Petrin (2003); Ackerberg, Caves & Frazer (2015); Wooldridge (2009).

wooldridge_prod

wooldridge_prod(data: DataFrame, output: str = 'y', free: Sequence[str] | str | None = None, state: Sequence[str] | str | None = None, proxy: str = 'm', panel_id: str = 'id', time: str = 'year', polynomial_degree: int = 3, productivity_degree: Optional[int] = None, functional_form: str = 'cobb-douglas', boot_reps: int = 0, seed: Optional[int] = None, vce: str = 'cluster', convention: str = 'statspai') -> ProductionResult

Wooldridge (2009) production function estimator.

Estimates the level equation and the productivity-substituted equation jointly; they share the input elasticities and the control-function polynomial h(state, proxy). By default the Markov process omega_t = g(omega_{t-1}) + xi_t has a free polynomial g and the system is estimated by GMM (equation 1 instrumented by current free inputs and h, equation 2 by state, lagged free inputs and lagged h). convention="prodest" reproduces Stata prodest, method(wrdg), which imposes a unit slope (a random walk) and runs one stacked 2SLS. When productivity mean-reverts, that restriction biases the capital elasticity (see the module docstring).

Parameters:

Name	Type	Description	Default
data	DataFrame	Long panel with one row per (firm, year).	required
output	str	Log output column.	``"y"``
free	str or list	Free inputs (labor).	``["l"]``
state	str or list	State inputs (capital).	``["k"]``
proxy	str	Productivity proxy (typically intermediate input).	``"m"``
panel_id	str	Firm and (numeric) year identifiers; lags are calendar lags.	'id'
time	str	Firm and (numeric) year identifiers; lags are calendar lags.	'id'
polynomial_degree	int	Degree of h(state, proxy) (Stata prodest's default).	3
productivity_degree	int	Degree of g (default 3, as OP / LP / ACF). Must be left None with convention="prodest", which fixes g(s) = s.	None
functional_form	'cobb-douglas'	Translog raises NotImplementedError.	'cobb-douglas'
boot_reps	int	If positive, firm-cluster bootstrap standard errors replace the analytic ones.	0
seed	int	Bootstrap RNG seed.	None
vce	('cluster', 'robust', 'unadjusted')	Analytic covariance. "cluster" clusters by firm (G/(G-1)): both equations of a firm-year share eta. "robust" treats firm-years as independent. "unadjusted" (convention="prodest" only) is Stata prodest's, u'u/N over the stacked rows.	"cluster"
convention	('statspai', 'prodest')	Estimator, as described above.	"statspai"

Returns:

Type	Description
ProductionResult

Examples:

>>> import statspai as sp
>>> import numpy as np, pandas as pd
>>> rng = np.random.default_rng(0)
>>> rows = []
>>> for fid in range(60):
...     k = rng.normal(2.0, 0.5)
...     omega = rng.normal(0.0, 0.3)
...     w = rng.normal(0.0, 0.3)
...     for yr in range(8):
...         omega = 0.7 * omega + rng.normal(0.0, 0.2)
...         w = 0.8 * w + rng.normal(0.0, 0.3)
...         l = 0.6 * k + omega - 0.5 * w + rng.normal(1.0, 0.05)
...         m = 0.5 * k + omega + rng.normal(0.5, 0.05)
...         y = 0.6 * l + 0.3 * k + omega + rng.normal(0.0, 0.1)
...         rows.append(dict(id=fid, year=yr, y=y, l=l, k=k, m=m))
...         k = 0.9 * k + 0.3 * rng.normal(1.0, 0.3)
>>> df = pd.DataFrame(rows)
>>> res = sp.wooldridge_prod(df, output="y", free="l", state="k", proxy="m")
>>> sorted(res.coef.keys())
['k', 'l']
>>> res.method
'wrdg'

References

Wooldridge, J.M. (2009). On estimating firm-level production functions using proxy variables to control for unobservables. Economics Letters, 104(3), 112-114.