Metadata-Version: 2.4
Name: fasthdfe
Version: 0.1.0
Summary: FastHDFE provides reghdfe and ppmlhdfe to estimate linear and multiplicative models with high-dimensional fixed effects.
Project-URL: Documentation & Evaluation, https://doi.org/10.15495/EPub_UBT_00009520
Author-email: Mario Larch <mario.larch@uni-bayreuth.de>, Mirco Schoenfeld <mirco.schoenfeld@uni-bayreuth.de>
Maintainer-email: Mirco Schoenfeld <mirco.schoenfeld@uni-bayreuth.de>
License-Expression: GPL-3.0-or-later
License-File: LICENSE
Classifier: Development Status :: 4 - Beta
Classifier: Programming Language :: Python
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Programming Language :: Python :: 3.13
Classifier: Programming Language :: Python :: 3.14
Classifier: Programming Language :: Python :: Implementation :: CPython
Classifier: Programming Language :: Python :: Implementation :: PyPy
Requires-Python: >=3.11
Requires-Dist: numpy
Requires-Dist: pandas
Requires-Dist: scipy
Requires-Dist: tabulate
Description-Content-Type: text/markdown

# FastHDFE


This repository provides `FastHDFE`, a Python package that can be used to estimate linear and multiplicative models with high-dimensional fixed effects. 

The core of the package lies in the two functions called `reghdfe` and `ppmlhdfe`. The `reghdfe`-command contains an estimator for linear models, which is optimized for speed and will produce the common output typically expected after a linear regression. Similarly, the `ppmlhdfe`-command contains an estimator for multiplicative models with high-dimensional fixed effects. As we envision the commands to be useful for different fields, the commands allow the choice of the structure of fixed effects as well as the type of standard errors.

This page covers:

* [Requirements](#requirements)
* [Test dataset for code examples](#test-data)
* [The reghdfe command](#reghdfe)
* [The ppmlhdfe command](#ppmlhdfe)
* [Citation details](#citation)

## Requirements

* Numpy
* Scipy
* Pandas
* tabulate

## Test Data

The project uses a small dataset useful for testing. This dataset contains data for 62 exporters, 62 importers and 27
years. Use it like this to execute the example calls:

```
import pandas as pd
df = pd.read_csv('https://www.usitc.gov/data/gravity/example_trade_and_grav_data_small.csv')
df['lntrade_value'] = np.log(cls._df['trade_value'])
```

## reghdfe

`reghdfe` is a command available via our package to estimate linear
models with high-dimensional fixed effects that has the following
structure:

1.  Inputs:

    1.  dependent variable,

    2.  a set of regressors,

    3.  a flexible way to specify the type of fixed effects one wants to
        include,

    4.  options to choose convergence criteria,

    5.  options to choose from various types of standard errors
        (homoskedastic standard errors, heteroskedasticity-robust
        standard errors, clustered standard errors),

    6.  option to enable/disable conjugate gradient acceleration,

    7.  option to save fixed effects.

2.  Outputs:

    1.  coefficient estimates,

    2.  standard errors,

    3.  variance-covariance matrix,

    4.  t-values,

    5.  number of observations,

    6.  model diagnostics (e.g., Aikake Information Criterion, Bayesian
        Information Criterion),

    7.  fitted values, which can be merged with the input data,

    8.  residuals, which can be merged with the input data,

    9.  optional: fixed effects, which can be merged with the input
        data.

3.  Syntax:
    ```
    reghdfe(y,
            x, 
            fixedeffects              = None,
            data                      = None,
            setype                    = "none",
            clustvars                 = None,
            FEgen                     = False,
            verb                      = 0,
            eps_MAP                   = 1e-9,
            minsubiter                = 0,
            maxsubiter                = 1000,
            accel                     = True,
            nt                        = 0,
            print_stata_style_summary = True)
    ```
    

    1.  Parameters:

        -   `y`: single entry list containing the name of the dependent
            variable as a string.

        -   `x`: string list of the names of the explanatory variables.

        -   `fixedeffects`: string list of the names of the fixed
            effects. Variables can be combined by "\#", e.g.
            "state#year" creates state-year fixed effects. Default is
            `None`, which causes estimation with a constant.

        -   `data`: data frame containing all relevant variables.

        -   `setype`: string, type of standard errors. Allowed values
            are "none", "homoskedastic", "heteroskedastic" ("r",
            "robust" and "HC1" are synonyms), "cluster" ("multiway" as
            synonym). Default is `None`.

        -   `clustvars`: string list of cluster variables. Max length 4,
            ignored if setype is not "cluster". Default is `None`.

        -   `FEgen`: bool, if `True`, the observation-wise sum of the
            fixed effects, as well as the different types of fixed
            effects, are returned.

    2.  Additional Parameters:

        -   `accel`: bool, if `True`, MAP employs conjugate gradient
            acceleration technique following
            Hernández-Ramos et al. ( [2011](https://doi.org/10.1080/01630563.2011.591954) ). Default is `True`.

        -   `nt`: integer, sets the number of threads to use. Default is
            `0`, which uses "system settings". Has no effect if the
            backend is Python.

        -   `verb`: int $\in \{0, 1, 2\}$, level of "in-between" output
            displayed (0 is none). Default is `0`.

        -   `eps_MAP`: float $\in (0, \infty)$, convergence criterion of
            the alternating projections algorithm given. Default is
            $\mathtt{10^{-9}}$.

        -   `minsubiter`: minimum number of iterations within the
            demeaning process. Default is 0.

        -   `maxsubiter`: maximum number of iterations within the
            demeaning process. Default is 1000.

        -   `print_stata_style_summary`: bool, regression output table
            is printed. Default is `True`.

    3.  Returns a class object with the following attributes:

        -   `params`: point estimates.

        -   `bse`: standard errors if calculated, "NULL" otherwise.

        -   `covmat`: variance-covariance matrix if calculated, "NULL"
            otherwise.

        -   `covmat_nofix`: variance-covariance matrix prior to the
            "eigenvalue fix". Only present if `eigenvalue_fix` is
            `True`.

        -   `eigenvalue_fix`: `True` only if the covariance matrix was
            not positive semi-definite and the eigenvalue fix had to be
            applied (negative eigenvalues replaced by zero); otherwise
            the attribute is not set.

        -   `cov_type`: string, type of standard errors.

        -   `clustvars`: list, string list with cluster dimensions if
            `cov_type` is "cluster", "NULL" otherwise.

        -   `tvalues`: $t$-values for $H_0 = 0$.

        -   `pvalues`: $p$-values for a two-sided $t$-test with
            $H_0 = 0$.

        -   `nobs`: int, number of observations used for estimation.
            Some observations may have been excluded and do not count in
            here.

        -   `df_model`: number of parameters in $X$ plus the theoretical
            upper bound of fixed effects (collinear fixed effects count
            as estimated parameters here; this may lead to a slight
            over-estimation of small-sample corrected standard errors).

        -   `df_resid`: `nobs` $-$ `df_model`.

        -   `aic`: Aikake Information Criterion.

        -   `bic`: Bayesian Information Criterion.

        -   `llf`: value of the log-likelihood.

        -   `deviance`: sum of the squared deviance residuals, equal to
            the sum of squared residuals.

        -   `pearson_chi2`: Pearson's $\chi^2$-value, equal to the sum
            of squared residuals.

        -   `family_name`: "Gaussian distribution".

        -   `family_link`: "Identity function", ( $g(\cdot)$ ).

        -   `scale`: ML estimate of $\sigma^2$. To obtain the OLS
            estimates, multiply by nobs/(nobs-1).

        -   `method`: "OLS".

        -   `fit_history`: dict with key 'iteration', equal to 1, as
            analytical OLS solution is computed.

        -   `fittedvalues`: array of the fitted values.

        -   `fe_rowsum`: array of the observation-wise sum of the fixed
            effects $\boldsymbol{\eta}-\mathbf{x}'\boldsymbol{\beta}$.

        -   `fixed_effects`: dictionary of all estimates of all fixed
            effects specified.

        -   `resid`: array of residuals.

        -   `msg`: list, information/warning messages during the
            estimation process.

        -   `formula_R`: string, formula as you would type it in R.

        -   `formula_Stata`: string, formula as you would type it in
            Stata.

        -   `model`: GLM class, for compatibility with the gme-package.

    4.  Example:
        ```
        from fasthdfe import reghdfe
        result = reghdfe(y            = ['lntrade_value'],
                         x            = ['LN_DIST','agree_pta','contiguity','common_language'],
                         fixedeffects = ['exp_ind#industry_id#year','imp_ind#industry_id#year','exp_ind#imp_ind#industry_id'],
                         data         = df,
                         setype       = 'cluster',
                         clustvars    = ['exp_ind','imp_ind','year','industry_id'])
        ```

## ppmlhdfe

`ppmlhdfe` is a user-friendly Python command to estimate multiplicative
models with high-dimensional fixed effects with PPML that has the
following structure:

1.  Inputs:

    1.  dependent variable,

    2.  a set of regressors,

    3.  a flexible way to specify the type of fixed effects one wants to
        include,

    4.  options to choose convergence criteria,

    5.  options to choose from various types of standard errors
        (homoskedastic standard errors, heteroskedasticity-robust
        standard errors, clustered standard errors),

    6.  options for initial values of the linear predictor (simple,
        OLS-based, or user-supplied starting values for selected
        coefficients),

    7.  options to detect and drop separated observations
        (fixed-effect-style separation, iterative rectifier, or none,
        mirroring Stata's `separation()` option),

    8.  option to keep or drop singleton observations,

    9.  an offset option for constrained coefficients,

    10. option to enable/disable conjugate gradient acceleration,

    11. option to standardize regressors before estimation for improved
        numerical stability (mirroring Stata's `standardize_data(1)`
        default),

    12. option to save fixed effects.

2.  Outputs:

    1.  coefficient estimates,

    2.  standard errors,

    3.  variance-covariance matrix,

    4.  t-values,

    5.  number of observations,

    6.  diagnostics on the number of observations kept and dropped
        (total, invalid values, singletons/fe-separation, and separation
        detected by the iterative rectifier),

    7.  model diagnostics (e.g., Aikake Information Criterion, Bayesian
        Information Criterion, deviance and pseudo-likelihood),

    8.  fitted values, which can be merged with the input data,

    9.  residuals, which can be merged with the input data,

    10. optional: fixed effects, which can be merged with the input
        data.

3.  Syntax:
    ```
    ppmlhdfe(y,
             x,
             fixedeffects              = None, 
             data                      = None, 
             setype                    = "none", 
             clustvars                 = None,
             off                       = None, 
             eta_0                     = None, 
             FEgen                     = False,
             start                     = "simple", 
             olsguess                  = False,   # olsguess kept for back-compat
             separation                = "ir",
             sep_tol                   = 1e-5, 
             sep_maxiter               = 100, 
             sep_K_penalty             = 1e6,
             keepsingletons            = False,
             nocheck                   = False, 
             verb                      = 0,
             eps_beta                  = 1e-8, 
             miniter                   = 1, 
             maxiter                   = 1000,
             eps_MAP                   = 1e-9, 
             minsubiter                = 0, 
             maxsubiter                = 16000,
             start_inner_tol           = 1e-4,
             standardize_data          = True,
             accel                     = True, 
             nt                        = 0, 
             print_stata_style_summary = True)
    ```

    1.  Parameters:

        -   `y`: single entry list containing the name of the dependent
            variable as string.

        -   `x`: string list of the names of the explanatory variables.

        -   `fixedeffects`: string list of the names of the fixed
            effects. Variables can be combined by "\#", e.g.
            "state#year" creates state-year fixed effects. Default is
            `None`, which causes estimation with a constant.

        -   `data`: data frame containing all relevant variables.

        -   `setype`: string, type of standard errors. Allowed values
            are "none", "homoskedastic", "heteroskedastic" ("r",
            "robust" and "HC1" are synonyms), "cluster" ("multiway" as
            synonym). Default is `None`.

        -   `clustvars`: string list of cluster variables. Max length 4,
            ignored if setype is not "cluster". Default is `None`.

        -   `FEgen`: bool, if `True`, the observation-wise sum of the
            fixed effects (`fe_rowsum`), as well as the different types
            of fixed effects (`fixed_effects`), can be accessed after
            estimation (when the object is saved with "results", then
            with `result.fe_rowsum` and `result.fixed_effects`,
            respectively).

        -   `eta_0`: single entry list containing the string of the
            initial values for the linear predictor $\eta$. If set,
            overrides `start`. Default is `None`.

        -   `start`: strategy for the initial values of $\eta$. Accepts
            one of the following:

            -   `“simple”` (default):
                $\eta = \ln\big((y + \overline{y})/2\big)$, the Stata
                `ppmlhdfe`-style default.

            -   `“ols”`: use OLS on the fixed-effect-demeaned data, then
                take $\eta = \ln(\text{fitted values})$.

            -   `dict`: user-supplied starting values for selected
                regression coefficients,
                e.g. `{“log_distance”:-1.0, “contiguity”: 0.4}`.
                Unspecified variables start at zero.

        -   `olsguess`: bool, deprecated. Kept for backward
            compatibility. If `True`, equivalent to `start=“ols”`.
            Default is `False`.

        -   `off`: single entry list containing the string of the offset
            variable.

        -   `separation`: string, list of strings, bool, or `None`,
            default `“ir”`. Controls detection and removal of separated
            observations, mirroring the `separation()` option of Stata's
            `ppmlhdfe`. Allowed values:

            -   `“fe”`: drop observations whose fixed-effect group has
                $y$ constant across all observations (this is already
                performed unconditionally by the singleton/constant-$y$
                loop; listing it is informational).

            -   `“ir”` / `“relu”`: iterative rectifier; detects
                separation arising from both fixed effects *and*
                regressors ([Correia et al. 2021](https://arxiv.org/pdf/1903.01633.pdf)).

            -   `“none”` / `“off”` / `False` / `None`: disable
                additional separation checks (the implicit "fe"-style
                check still runs).

            -   `“default”` / `“auto”` / `True`: equivalent to `[“ir”]`
                (the default, matching Stata).

            -   list: apply all listed methods, e.g. `[“fe”,“ir”]`.

        -   `sep_tol`: float, threshold above which an observation is
            flagged as a separation certificate by the iterative
            rectifier. Default is $10^{-5}$.

        -   `sep_maxiter`: int, maximum number of outer iterations of
            the iterative rectifier. Default is $100$.

        -   `sep_K_penalty`: float, penalty weight on $y>0$ observations
            in the IR weighted regression. Too large can cause numerical
            issues in the inner weighted MAP; too small allows leakage
            onto $y>0$ observations. Default is $10^{6}$.

        -   `keepsingletons`: bool, default `False` (matching Stata's
            default). If `False`, drop any fixed-effect group with no
            within-group variation in $y$ (singletons, constant-$y$
            groups, and classical fe-style separation). If `True`, only
            drop groups where $y$ is identically zero---the minimum
            required for Poisson identification---and keep singletons
            and constant-positive-$y$ groups. Useful for replicating
            analyses that report the full sample, but singletons carry
            no identifying information and can distort cluster-robust
            standard errors.

        -   `standardize_data`: bool, default `True` (matching Stata's
            `standardize_data(1)`). If `True`, each regressor is divided
            by its sample standard deviation before estimation, and the
            resulting coefficients and variance-covariance matrix are
            rescaled to original units at the end. This substantially
            improves numerical stability when regressors have very
            different magnitudes (a common case in trade gravity, where
            dummies sit next to log distances or large monetary values),
            and prevents the cluster-robust variance computation from
            becoming near-singular when a regressor's identifying
            variation is concentrated in a small subset of observations
            after fixed-effect absorption.

    2.  Additional Parameters:

        -   `accel`: bool, if `True`, MAP employs conjugate gradient
            acceleration technique following
            Hernández-Ramos et al. ( [2011](https://doi.org/10.1080/01630563.2011.591954) ). Default is `True`.

        -   `nt`: integer, sets the number of threads to use. Default is
            `0`, which uses "system settings".

        -   `nocheck`: bool, if `True`, no consistency check is skipped.
            Default is `False`.

        -   `verb`: int $\in \{0, 1, 2\}$, level of output displayed (0
            is none). Default is $0$.

        -   `eps_beta`: float $\in (0, \infty)$, convergence criterion
            of the IRLS loop,
            $||\boldsymbol{\beta}^{(r+1)}-\boldsymbol{\beta}^{(r)}||_2 < \mathtt{eps\_beta}$.
            Convergence is declared only once this criterion is
            satisfied in two consecutive iterations, mirroring Stata's
            convention. Default is $\mathtt{10^{-8}}$, matching Stata
            `ppmlhdfe`'s `tolerance(1e-8)`.

        -   `eps_MAP`: float $\in (0, \infty)$, target convergence
            criterion of the alternating projections algorithm. The
            *effective* inner tolerance applied at convergence is
            $\max(\min(\mathtt{eps\_MAP},\,0.1\cdot\mathtt{eps\_beta}),\,10^{-12})$,
            following Stata `reghdfe`'s rule that the inner solve must
            be at least ten times tighter than the outer IRLS criterion
            but never tighter than machine precision. Default is
            $\mathtt{10^{-9}}$, matching Stata `reghdfe`'s `itol(1e-9)`.

        -   `start_inner_tol`: float $\in (0, \infty)$, starting
            tolerance for the inner partialling-out loop on the first
            IRLS iteration. The inner tolerance is loosened to this
            value early in IRLS (when the outer
            $\boldsymbol{\beta}$-change is large) and tightens
            proportionally to the outer convergence measure across
            iterations, with the floor described under `eps_MAP`. This
            adaptive schedule avoids spending inner iterations on
            precision that the outer loop cannot yet exploit, while
            guaranteeing the required tightness at convergence. Default
            is $\mathtt{10^{-4}}$, matching Stata `reghdfe`'s
            `start_inner_tol(1e-4)`.

        -   `miniter`: minimum number of iterations of the IRLS loop.
            Default is $1$.

        -   `maxiter`: maximum number of iterations of the IRLS loop.
            Default is $1000$, matching Stata `ppmlhdfe`'s
            `maxiter(1000)`.

        -   `minsubiter`: minimum number of iterations within the
            demeaning process. Default is $0$.

        -   `maxsubiter`: maximum number of iterations within the
            demeaning process. Default is $16000$, matching Stata
            `reghdfe`'s `iterations(16000)`.

        -   `print_stata_style_summary`: bool, regression output table
            and Stata-style drop/separation notices are printed. Default
            is `True`.

    3.  Returns a class object with the following attributes:

        -   `params`: point estimates.

        -   `bse`: standard errors if calculated, "NULL" otherwise.

        -   `covmat`: variance-covariance matrix if calculated, "NULL"
            otherwise.

        -   `covmat_nofix`: variance-covariance matrix prior to the
            "eigenvalue fix". Only present if `eigenvalue_fix` is
            `True`.

        -   `eigenvalue_fix`: `True` only if the covariance matrix was
            not positive semi-definite and the eigenvalue fix had to be
            applied (negative eigenvalues replaced by zero); otherwise
            the attribute is not set.

        -   `cov_type`: string, type of standard errors.

        -   `clustvars`: list, string list with cluster dimensions if
            `cov_type` is "cluster", "NULL" otherwise.

        -   `tvalues`: $t$-values for $H_0 = 0$.

        -   `pvalues`: $p$-values for a two-sided $t$-test with
            $H_0 = 0$.

        -   `nobs`: int, number of observations used for estimation.
            Some observations may have been excluded and do not count in
            here.

        -   `df_model`: number of parameters in $X$ plus the theoretical
            upper bound of fixed effects (collinear fixed effects count
            as estimated parameters here; this may lead to a slight
            over-estimation of small-sample corrected standard errors).

        -   `df_resid`: `nobs` $-$ `df_model`.

        -   `aic`: Aikake Information Criterion.

        -   `bic`: Bayesian Information Criterion.

        -   `llf`: value of the log-likelihood.

        -   `deviance`: sum of the squared deviance residuals. For
            $y=0$, the deviance is approximated using L'Hôpital's rule.

        -   `pearson_chi2`: Pearson's $\chi^2$-value.

        -   `family_name`: "Poisson distribution".

        -   `family_link`: "Natural logarithm function", ( $g(\cdot)$ ).

        -   `scale`: 1, per definition of the Poisson distribution.

        -   `method`: "IRLS", Iteratively re-weighted least squares.

        -   `fit_history`: dict with key 'iteration', which is the
            number of iterations until convergence.

        -   `fittedvalues`: array of the fitted values
            $\hat{y}=\mu = g^{-1}(\eta) = \exp(\eta)$.

        -   `fe_rowsum`: array of the observation-wise sum of the fixed
            effects $\boldsymbol{\eta} - \mathbf{x}'\boldsymbol{\beta}$.

        -   `fixed_effects`: dictionary of all estimates of all fixed
            effects specified.

        -   `resid_additive`: array of additive residuals ( $y-\hat{y}$ ).

        -   `resid_multiplicative`: array of multiplicative residuals
            ( $y/\hat{y}$ ).

        -   `msg`: list, information/warning messages during the
            estimation process.

        -   `formula_R`: string, formula as you would type it in R.

        -   `formula_Stata`: string, formula as you would type it in
            Stata.

        -   `model`: GLM class, for compatibility with the gme-package.

        -   `separation_methods`: list of strings, the separation
            methods that were applied (e.g. `[“ir”]`).

        -   `num_separated`: int, number of observations dropped by the
            iterative rectifier. Equal to `N_dropped_ir_sep`.

        -   `keepsingletons`: bool, the value of the `keepsingletons`
            argument passed in.

        -   `N_full_initial`: int, number of observations in the input
            data frame before any drops.

        -   `N_dropped_invalid`: int, number of observations dropped due
            to missing, infinite, or negative values.

        -   `N_dropped_fe_sep`: int, number of observations dropped by
            the singleton / constant-$y$ / fe-style separation loop.

        -   `N_dropped_ir_sep`: int, number of observations dropped by
            the iterative rectifier.

    4.  Example:
        ```
        from fasthdfe import ppmlhdfe
        result = ppmlhdfe(y            = ['trade'],
                          x            = ['LN_DIST','agree_pta','contiguity','common_language'],
                          fixedeffects = ['exp_ind#industry_id#year','imp_ind#industry_id#year','exp_ind#imp_ind#industry_id'],
                          data         = df,
                          eps_beta     = 1e-5,
                          setype       = 'cluster',
                          clustvars    = ['exp_ind', 'imp_ind','year','industry_id'])
        ```


# Citation


```
@misc{ubt_epub9520,
           title = {Fast Estimation of Linear and Poisson Models with High-Dimensional Fixed Effects in Python : The FastHDFE package},
            note = {This is a technical report corresponding to the Python package https://pypi.org/project/fasthdfe},
         address = {Bayreuth, Germany},
           month = {July},
            year = {2026},
        keywords = {high-dimensional fixed effects, Poisson pseudo-maximum-likelihood, gravity models, method of alternating projections, multi-way clustering, Python, econometric software},
          author = {Larch, Mario and Schoenfeld, Mirco and Shikher, Serge},
             url = {https://epub.uni-bayreuth.de/id/eprint/9520/},
             doi = {https://doi.org/10.15495/EPub_UBT_00009520}
}
```
