### Setup: Import PyFixest and load data

Source: https://pyfixest.org/tutorials/formula-syntax.html

Imports numpy and pyfixest, then loads the example dataset using pf.get_data().

```python
import numpy as np
import pyfixest as pf
data = pf.get_data()
data.head()
```

--------------------------------

### Load PyFixest and get example data

Source: https://pyfixest.org/tutorials/quantile-regression.html

Sets up the environment by loading the autoreload extension and importing PyFixest, then retrieves the example dataset used in the subsequent examples.

```python
%load_ext autoreload

import pyfixest as pf
data = pf.get_data()
```

--------------------------------

### Quickstart: estimate OLS with fixed effects

Source: https://pyfixest.org/pyfixest.html

Basic PyFixest usage: get sample data, estimate a fixed-effects OLS model, and print the summary. The output shows the estimation details and coefficient table.

```python
import pyfixest as pf

data = pf.get_data()
pf.feols("Y ~ X1 | f1 + f2", data=data).summary()
```

```text
###

Estimation:  OLS
Dep. var.: Y, Fixed effects: f1+f2
Inference:  CRV1
Observations:  997

| Coefficient   |   Estimate |   Std. Error |   t value |   Pr(>|t|) |   2.5% |   97.5% |
|:--------------|-----------:|-------------:|----------:|-----------:|-------:|--------:|
| X1            |     -0.919 |        0.065 |   -14.057 |      0.000 | -1.053 |  -0.786 |
---
RMSE: 1.441   R2: 0.609   R2 Within: 0.2
```

--------------------------------

### Install prek globally

Source: https://pyfixest.org/contributing.html

Install prek globally using pixi, recommended for contributors.

```bash
pixi global install prek
```

--------------------------------

### Set up git hooks with prek

Source: https://pyfixest.org/contributing.html

Run `prek install` once after cloning to install prek as a git hook so checks run automatically on each commit.

```bash
prek install
```

--------------------------------

### Example: create_model_matrix usage

Source: https://pyfixest.org/reference/estimation.formula.model_matrix.create_model_matrix.html

Example usage of create_model_matrix with a parsed formula and data from pyfixest. Imports required modules and creates a model matrix.

```python
import pyfixest as pf
from pyfixest.estimation.formula.model_matrix import create_model_matrix
from pyfixest.estimation.formula.parse import Formula

data = pf.get_data()
formula = Formula.parse("Y ~ X1 + f1 + f2")[0]
model_matrix = create_model_matrix(formula=formula, data=data)
```

--------------------------------

### Install PyFixest from PyPI

Source: https://pyfixest.org/pyfixest.html

Installs the release version of PyFixest from PyPI. Run inside an active virtual environment.

```bash
# inside an active virtual environment
python -m pip install pyfixest
```

--------------------------------

### Install PyFixest development version from GitHub

Source: https://pyfixest.org/pyfixest.html

Installs the development version of PyFixest directly from the GitHub repository.

```bash
python -m pip install git+https://github.com/py-econometrics/pyfixest
```

--------------------------------

### Setup: Import libraries and enable rpy2 conversion

Source: https://pyfixest.org/explanation/compare-fixest-pyfixest.html

Imports pandas, rpy2, and pyfixest, imports R packages (fixest, stats, broom), and sets up a custom converter to keep R list classes intact. Also loads IPython autoreload and fetches example data with pyfixest.

```python
import pandas as pd
import rpy2.robjects as ro
from rpy2.robjects import numpy2ri, pandas2ri
from rpy2.robjects.conversion import Converter
from rpy2.robjects.packages import importr
from rpy2.robjects.vectors import ListVector
from rpy2.rinterface import ListSexpVector

import pyfixest as pf

# Import R packages
fixest = importr("fixest")
stats = importr("stats")
broom = importr("broom")

# Enable automatic pandas <-> R DataFrame conversion (rpy2 >= 3.6)
# Keep R list classes intact (avoid conversion to NamedList).
_list_guard = Converter("pandas2ri-list-guard")

@_list_guard.rpy2py.register(ListSexpVector)
def _keep_r_list_classes(obj):
    return ListVector(obj)

_converter_ctx = (
    ro.default_converter + numpy2ri.converter + pandas2ri.converter + _list_guard
).context()
_converter_ctx.__enter__()

# IPython magic commands for autoreloading
%load_ext autoreload
%autoreload 2

# Get data using pyfixest
data = pf.get_data(model="Feols", N=10_000, seed=99292)
```

--------------------------------

### Install PyFixest pre-release 0.50.0a2

Source: https://pyfixest.org/changelog.html

Commands to install the latest pre-release version of PyFixest to try out new features. The first command installs a specific alpha version, and the second installs the latest pre-release.

```bash
pip install pyfixest==0.50.0a2
pip install --pre pyfixest
```

--------------------------------

### Install PyFixest via pip

Source: https://pyfixest.org/getting-started.html

Install the latest version of PyFixest from PyPI.

```bash
pip install -U pyfixest
```

--------------------------------

### Example usage of wildboottest with different bootstrap types

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.html

Example code showing how to run wildboottest with bootstrap_type '11' and '31', including the output for the '31' case.

```python
import re
import pyfixest as pf

data = pf.get_data()
fit = pf.feols("Y ~ X1 + X2 | f1", data)

fit.wildboottest(
    param = "X1",
    reps=1000,
    seed = 822
)

fit.wildboottest(
    param = "X1",
    reps=1000,
    seed = 822,
    bootstrap_type = "31"
)
```

```output
param                    X1
t value          -14.843741
Pr(>|t|)                0.0
bootstrap_type           31
inference                HC
impose_null            True
ssc                       1
dtype: object
```

--------------------------------

### Feols.update example output

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.html

Example output from the update method, showing the updated coefficients as a numpy array.

```text
array([ 0.88955689, -0.99519687, -0.17661729])
```

--------------------------------

### Example: Create treatment assignment and outcome plots with panelview()

Source: https://pyfixest.org/reference/did.visualize.panelview.html

This example loads a heterogeneous treatment effects dataset from a URL and creates two panelview plots: one for treatment assignment and one for the outcome variable. It uses subsamp=50 to limit the number of units displayed and sets custom titles.

```python
import pandas as pd
import numpy as np
import pyfixest as pf

url = "https://raw.githubusercontent.com/py-econometrics/pyfixest/master/pyfixest/did/data/df_het.csv"
df_het = pd.read_csv(url)

# Inspect treatment assignment
pf.panelview(
    data = df_het,
    unit = "unit",
    time = "year",
    treat = "treat",
    subsamp = 50,
    title = "Treatment Assignment"
)

# Outcome plot
pf.panelview(
    data = df_het,
    unit = "unit",
    time = "year",
    outcome = "dep_var",
    treat = "treat",
    subsamp = 50,
    title = "Outcome Plot"
)
```

--------------------------------

### get_data

Source: https://pyfixest.org/reference/

Creates a random example data set.

```APIDOC
## get_data

### Description
Create a random example data set.

### Method
Function call

### Endpoint
`pyfixest.utils.utils.get_data`

### Parameters
No parameters documented.

### Request Example
No request example provided.

### Response
#### Success Response
No response schema documented.

#### Response Example
No response example provided.
```

--------------------------------

### Example output of coef()

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.coef.html

Shows the expected output of the coef() method for the example model. The coefficients are displayed as a pandas Series with index 'Coefficient' and values for X1 and X2.

```text
Coefficient
X1   -0.949526
X2   -0.174225
Name: Estimate, dtype: float64
```

--------------------------------

### Import PyFixest and demeaner backends

Source: https://pyfixest.org/how-to/demeaner-backends.html

Imports the main pyfixest package and the LsmrDemeaner and MapDemeaner classes, and loads a sample dataset. This is the setup for using the demeaner backends.

```python
import pyfixest as pf
from pyfixest.demeaners import LsmrDemeaner, MapDemeaner

df = pf.get_data()
```

--------------------------------

### Install PyFixest with plots extra

Source: https://pyfixest.org/pyfixest.html

Installs PyFixest with the 'plots' extra to enable visualization features using the lets-plot backend. matplotlib is included by default, so this extra is optional.

```bash
python -m pip install pyfixest[plots]
```

--------------------------------

### Instantiate MapDemeaner and LsmrDemeaner

Source: https://pyfixest.org/reference/demeaners.BaseDemeaner.html

Creates instances of the default MapDemeaner and LsmrDemeaner backends. This example shows the basic usage for choosing a demeaner backend.

```python
import pyfixest as pf

pf.MapDemeaner(), pf.LsmrDemeaner()
```

--------------------------------

### Install PyFixest with torch extra

Source: https://pyfixest.org/pyfixest.html

Installs PyFixest with the 'torch' extra to enable the LSMR demeaner via PyTorch (CPU and GPU). For CUDA GPU acceleration, a CUDA-enabled torch build is also required.

```bash
python -m pip install pyfixest[torch]
```

--------------------------------

### Basic feols regression examples

Source: https://pyfixest.org/changelog.html

Shows how to load PyFixest data and run three feols regressions with increasing complexity, including fixed effects.

```python
import pyfixest as pf
from pyfixest.report.utils import rename_categoricals

df = pf.get_data()

fit1 = pf.feols("Y ~ X1", data = df)
fit2 = pf.feols("Y ~ X1 + X2", data = df)
fit3 = pf.feols("Y ~ X1 + X2 | f1", data = df)
```

--------------------------------

### Example usage of IV_Diag with feols

Source: https://pyfixest.org/reference/estimation.models.feiv_.Feiv.html

Complete example showing how to simulate data, fit an IV model with feols, run first_stage, and call IV_Diag to compute F statistics. Output shows the unadjusted and effective F statistics.

```python
import numpy as np
import pandas as pd
from pyfixest.estimation import feols

# Set random seed for reproducibility
np.random.seed(1)

# Number of observations
n = 1000

# Simulate the data
# Instrumental variable
z = np.random.binomial(1, 0.5, size=n)
z2 = np.random.binomial(1, 0.5, size=n)

# Endogenous variable
d = 0.5 * z + 1.5 * z2 + np.random.normal(size=n)

# Control variables
c1 = np.random.normal(size=n)
c2 = np.random.normal(size=n)

# Outcome variable
y = 1.0 + 1.5 * d + 0.8 * c1 + 0.5 * c2 + np.random.normal(size=n)

# Cluster variable
cluster = np.random.randint(1, 50, size=n)
weights = np.random.uniform(1, 3, size=n)

# Create a DataFrame
data = pd.DataFrame({
    'd': d,
    'y': y,
    'z': z,
    'z2': z2,
    'c1': c1,
    'c2': c2,
    'cluster': cluster,
    'weights': weights
})

vcov_detail = "iid"

# Fit OLS model
fit_ols = feols("y ~ 1 + d + c1 + c2", data=data, vcov=vcov_detail)

# Fit IV model
fit_iv = feols("y ~ 1 + c1 + c2 | d ~ z", data=data,
         vcov=vcov_detail,
         weights="weights")
fit_iv.first_stage()
F_stat_pf = fit_iv._f_stat_1st_stage
fit_iv.IV_Diag()
F_stat_eff_pf = fit_iv._eff_F

print("(Unadjusted) F stat :", F_stat_pf)
print("Effective F stat :", F_stat_eff_pf)
```

```text
(Unadjusted) F stat : 52.81535560457512
Effective F stat : 48.61288119858717
```

--------------------------------

### iplot() example with multiple models

Source: https://pyfixest.org/reference/report.coefplot.html

Example usage of iplot() to plot coefficients from three fitted models. It shows how to create models with pyfixest and plot them together.

```python
import pyfixest as pf

df = pf.get_data()
fit1 = pf.feols("Y ~ i(f1)", data = df)
fit2 = pf.feols("Y ~ i(f1) + X2", data = df)
fit3 = pf.feols("Y ~ i(f1) + X2 | f1", data = df)

pf.iplot([fit1, fit2, fit3])
```

--------------------------------

### Setup imports and watermark for PyFixest DID

Source: https://pyfixest.org/tutorials/difference-in-differences.html

Imports required libraries and loads extensions for the notebook. Output shows pandas 2.3.3 and pyfixest 0.40.1.

```python
from importlib import resources

import pandas as pd

import pyfixest as pf
from pyfixest.report.utils import rename_event_study_coefs
from pyfixest.utils.dgps import get_sharkfin

%load_ext watermark
%watermark --iversions
%load_ext autoreload
%autoreload 2
```

```text
pandas  : 2.3.3
pyfixest: 0.40.1
```

--------------------------------

### Import pyfixest and dependencies

Source: https://pyfixest.org/tutorials/poisson-glm.html

Imports numpy, pandas, and pyfixest for the Poisson regression examples.

```python
import numpy as np
import pandas as pd
import pyfixest as pf
```

--------------------------------

### Load Poisson regression dataset

Source: https://pyfixest.org/explanation/compare-fixest-pyfixest.html

Loads the Poisson regression example dataset using pyfixest's get_data function with model='Fepois'.

```python
data = pf.get_data(model="Fepois")
```

--------------------------------

### Frequency weights example with feols and weights_type

Source: https://pyfixest.org/changelog.html

Demonstrates using frequency weights with `pf.feols()` by aggregating data and passing a `weights` column with `weights_type="fweights"`. The example compares an unweighted fit with a weighted fit using `pf.etable()`. The output shows the aggregated data dimensions and the regression table.

```python
import pyfixest as pf

data = pf.get_data(N = 10000, model = "Fepois")
df_weighted = data[["Y", "X1", "f1"]].groupby(["Y", "X1", "f1"]).size().reset_index().rename(columns={0: "count"})
df_weighted["id"] = list(range(df_weighted.shape[0]))

print("Dimension of the aggregated df:", df_weighted.shape)
print(df_weighted.head())

fit = pf.feols(
    "Y ~ X1 | f1",
    data = data
)
fit_weighted = pf.feols(
    "Y ~ X1 | f1",
    data = df_weighted,
    weights = "count",
    weights_type = "fweights"
)
pf.etable([fit, fit_weighted], coef_fmt = "b(se) \n (t) \n (p)")
```

```text
Dimension of the aggregated df: (1278, 5)
     Y   X1   f1  count  id
0  0.0  0.0  0.0     17   0
1  0.0  0.0  1.0     11   1
2  0.0  0.0  2.0     10   2
3  0.0  0.0  3.0     17   3
4  0.0  0.0  4.0     14   4
```

--------------------------------

### Example: Computing SAVI p-values with pyfixest

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.pvalue_savi.html

Example usage of pvalue_savi() on a fitted model. It loads sample data, fits a linear model with HC1 robust standard errors, and computes sequential p-values for each coefficient.

```python
import pyfixest as pf

data = pf.get_data()
fit = pf.feols("Y ~ X1 + X2", data=data, vcov="HC1")
fit.pvalue_savi()
```

```text
Intercept    1.614472e-13
X1           1.879119e-30
X2           4.056323e-13
Name: Pr(>|t|), dtype: float64
```

--------------------------------

### Basic Gelbach decomposition example

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.decompose.html

Example usage of decompose() with a simple decomposition, including creating data, fitting a model, and calling decompose with reps=10 and nthreads=1. Shows how to access results with tidy() and etable().

```python
import re
import pyfixest as pf
from pyfixest.utils.dgps import gelbach_data

data = gelbach_data(nobs = 1000)
fit = pf.feols("y ~ x1 + x21 + x22 + x23", data=data)

# simple decomposition
gb = fit.decompose(decomp_var = "x1", reps = 10, nthreads = 1)
type(gb)

gb.tidy()
gb = fit.decompose(decomp_var = "x1", reps = 10, nthreads = 1, x1_vars = ["x21"])
# combine covariates
gb = fit.decompose(decomp_var = "x1", reps = 10, nthreads = 1, combine_covariates = {"g1": ["x21", "x22"], "g2": ["x23"]})
# supress inference
gb = fit.decompose(decomp_var = "x1", reps = 10, nthreads = 1, combine_covariates = {"g1": ["x21", "x22"], "g2": ["x23"]}, only_coef = True)
# print results
gb.etable()

# group covariates via regex
res = fit.decompose(decomp_var="x1", combine_covariates={"g1": re.compile("x2[1-2]"), "g2": re.compile("x23")})
```

--------------------------------

### Minimal IV example with fixed effects and summary output

Source: https://pyfixest.org/tutorials/formula-syntax.html

Estimates an IV model with fixed effects and prints the summary. The output shows the estimation type, coefficients, and inference details.

```python
fit_iv = pf.feols("Y ~ X2 | f1 + f2 | X1 ~ Z1", data=data)
fit_iv.summary()
```

```text
###

Estimation:  IV
Dep. var.: Y, Fixed effects: f1+f2
sample: None = all
Inference:  iid
Observations:  997

| Coefficient   |   Estimate |   Std. Error |   t value |   Pr(>|t|) |   2.5% |   97.5% |
|:--------------|-----------:|-------------:|----------:|-----------:|-------:|--------:|
| X2            |     -0.174 |        0.015 |   -11.701 |      0.000 | -0.204 |  -0.145 |
| X1            |     -1.050 |        0.089 |   -11.793 |      0.000 | -1.225 |  -0.875 |
---
```

--------------------------------

### iplot() with rename_models using formula strings

Source: https://pyfixest.org/reference/report.coefplot.html

Example of iplot() with rename_models using the formula strings as keys to rename models.

```python
pf.iplot(
    models = [fit1, fit2, fit3],
    rename_models = {
        "Y~i(f1)": "Model 1",
        "Y~i(f1)+X2": "Model 2",
        "Y~i(f1)+X2|f1": "Model 3"
    },
)
```

--------------------------------

### Predict Values with predict()

Source: https://pyfixest.org/reference/estimation.api.feols.feols.html

Use predict() to get predicted values for the fitted model. The example shows the first five predictions.

```python
fit = pf.feols("Y ~ X1 + X2 | f1 + f2", data)
fit.predict()[0:5]
```

```text
array([ 3.06336634, -0.69574135, -0.91240433, -0.46370256, -1.67331152])
```

--------------------------------

### Build the documentation

Source: https://pyfixest.org/contributing.html

Commands to build, render, and preview the documentation.

```bash
pixi run docs-build
pixi run docs-render
pixi run docs-preview
```

--------------------------------

### Example: Using pvalue() on a fitted feols model

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.pvalue.html

Fits a linear regression with fixed effects using pyfixest's feols and then calls pvalue() to get the coefficient p-values. The output is a pandas Series with the p-values for each coefficient.

```python
import pyfixest as pf

fit = pf.feols("Y ~ X1 + X2 | f1", pf.get_data())
fit.pvalue()
```

```text
Coefficient
X1    0.0
X2    0.0
Name: Pr(>|t|), dtype: float64
```

--------------------------------

### Fit feols model and call tidy()

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.tidy.html

Fits a fixed effects OLS model using pyfixest's feols() and calls tidy() to get a tidy DataFrame of regression results. The example uses the built-in dataset from get_data().

```python
import pyfixest as pf

fit = pf.feols("Y ~ X1 + X2 | f1", pf.get_data())
fit.tidy()
```

--------------------------------

### Import pyfixest and generate data with get_data()

Source: https://pyfixest.org/reference/utils.utils.get_data.html

Basic usage example: import pyfixest, call get_data() to create a DataFrame, and display the first few rows with head(). The resulting data includes columns Y, Y2, X1, X2, f1, f2, f3, group_id, Z1, Z2, and weights, with some missing values.

```python
import pyfixest as pf

data = pf.get_data()
data.head()
```

--------------------------------

### Example: Using se() to get standard errors

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.se.html

Fits a fixed effects model using pyfixest's feols() and then calls se() to retrieve the coefficient standard errors. The output shows the resulting pandas Series with standard errors for each coefficient.

```python
import pyfixest as pf

fit = pf.feols("Y ~ X1 + X2 | f1", pf.get_data())
fit.se()
```

```output
Coefficient
X1    0.066373
X2    0.017596
Name: Std. Error, dtype: float64
```

--------------------------------

### Create PDF from LaTeX table using pylatex

Source: https://pyfixest.org/tutorials/regression-tables.html

Defines make_pdf() to compile a LaTeX table into a PDF using pylatex. It adds required packages (booktabs, threeparttable, makecell, tabularx) and generates the PDF. The example shows a conditional run block and a FileLink to the generated PDF; requires a TeX installation.

```python
## Use pylatex to create a tex file with the table


def make_pdf(tab, file):
    "Create a PDF document with tex table."
    doc = pl.Document()
    doc.packages.append(pl.Package("booktabs"))
    doc.packages.append(pl.Package("threeparttable"))
    doc.packages.append(pl.Package("makecell"))
    doc.packages.append(pl.Package("tabularx"))

    with (
        doc.create(pl.Section("A PyFixest LateX Table")),
        doc.create(pl.Table(position="htbp")) as table,
    ):
        table.append(pl.NoEscape(tab))

    doc.generate_pdf(file, clean_tex=False)


# Compile latex to pdf & display a button with the hyperlink to the pdf
# requires tex installation
run = False
if run:
    make_pdf(tab, "latexdocs/SampleTableDoc")
display(FileLink("latexdocs/SampleTableDoc.pdf"))
```

--------------------------------

### Simulate backend migration certification with SAVI confidence sequences

Source: https://pyfixest.org/how-to/anytime-valid-inference.html

Runs a simulation for the Kotlin backend migration scenario. It defines parameters, simulates experiment paths, finds the first day the confidence sequence or confidence interval is within ±1 pp, and computes median certification days and certification rates. Requires prior definitions of ALPHA, N_REPRESENTATIVE_PATHS, pf, np, and helper functions _experiment_path and _representative_seed.

```python
DELTA = 0.01
BACKEND_USERS = 2_000
BACKEND_MAX_DAYS = 30
BACKEND_FIRST_SEED = 50_000
BACKEND_MIXTURE_PRECISION = pf.optimal_mixture_precision(
    nobs=2 * BACKEND_USERS * BACKEND_MAX_DAYS,
    number_of_coefficients=2,
    alpha=ALPHA,
)

backend_cert_days = np.zeros(N_REPRESENTATIVE_PATHS, dtype=int)
backend_ci_cert_days = np.zeros(N_REPRESENTATIVE_PATHS, dtype=int)
for simulation in range(N_REPRESENTATIVE_PATHS):
    candidate = _experiment_path(
        p_control=0.10,
        p_treatment=0.10,
        seed=BACKEND_FIRST_SEED + simulation,
        with_intervals=True,
        n_days=BACKEND_MAX_DAYS,
        users_per_day=BACKEND_USERS,
        mixture_precision=BACKEND_MIXTURE_PRECISION,
    )
    candidate_certified = (
        (candidate["cs_lower"] > -DELTA)
        & (candidate["cs_upper"] < DELTA)
    )
    if candidate_certified.any():
        backend_cert_days[simulation] = int(
            candidate.loc[candidate_certified, "day"].iloc[0]
        )
    candidate_ci_certified = (
        (candidate["ci_lower"] > -DELTA)
        & (candidate["ci_upper"] < DELTA)
    )
    if candidate_ci_certified.any():
        backend_ci_cert_days[simulation] = int(
            candidate.loc[candidate_ci_certified, "day"].iloc[0]
        )

median_cs_cert = int(np.median(backend_cert_days[backend_cert_days > 0]))
median_ci_cert = int(np.median(backend_ci_cert_days[backend_ci_cert_days > 0]))
cs_cert_rate = (backend_cert_days > 0).mean()
ci_cert_rate = (backend_ci_cert_days > 0).mean()

backend_seed = _representative_seed(
    backend_cert_days,
    BACKEND_FIRST_SEED,
)

backend_full = _experiment_path(
    p_control=0.10,
    p_treatment=0.10,
    seed=backend_seed,
    with_intervals=True,
    n_days=BACKEND_MAX_DAYS,
    users_per_day=BACKEND_USERS,
    mixture_precision=BACKEND_MIXTURE_PRECISION,
).set_index("day", drop=False)

certified = (backend_full["cs_lower"] > -DELTA) & (backend_full["cs_upper"] < DELTA)
cert_day = int(backend_full.loc[certified, "day"].iloc[0])
cert_cs = (
    f"[{backend_full.loc[cert_day, 'cs_lower']:.2%}, "
    f"{backend_full.loc[cert_day, 'cs_upper']:.2%}]"
)
ci_certified = (backend_full["ci_lower"] > -DELTA) & (backend_full["ci_upper"] < DELTA)
ci_cert_day = int(backend_full.loc[ci_certified, "day"].iloc[0])

backend_path = backend_full.loc[backend_full["day"] <= cert_day]
b_range = [
    float(backend_path["cs_lower"].min()) - 0.005,
    float(backend_path["cs_upper"].max()) + 0.005,
]
```

--------------------------------

### Quantile Regression with multi_method = 'cfm1' and 'cfm2'

Source: https://pyfixest.org/tutorials/quantile-regression.html

This example shows how to run a quantile regression process using two different algorithms (CFM1 and CFM2) via the `multi_method` argument. The first method (`cfm1`) is exactly identical to running separate regressions per quantile, while the second method (`cfm2`) is only asymptotically identical. The code sets up the formula, method, seed, and quantiles, then fits two models and displays the results.

```R
fml = "Y~X1"
method = "pfn"
seed = 929
quantiles = [0.1, 0.5, 0.9]

fit_multi1 = pf.quantreg(
  fml = fml,
  data = data,
  method = method,
  multi_method = "cfm1",  # this is algorithm 2 in CFM, the 1rst algorithm for the full qr process
  seed = seed,
  quantile = quantiles,
)

fit_multi2 = pf.quantreg(
  fml = fml,
  data = data,
  method = method,
  multi_method = "cfm2",  # this is algorithm 3 in CFM, the 2nd algorithm for the full qr process
  seed = seed,
  quantile = quantiles
)

pf.etable(fit_multi1.to_list() +  fit_multi2.to_list())
```

--------------------------------

### Example: computing confidence intervals with feols

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.confint.html

Example usage of confint() on a fitted feols model. It demonstrates regular, simultaneous, and SAVI confidence intervals. The output table shows the resulting confidence intervals for the SAVI example.

```python
from pyfixest.utils import get_data
from pyfixest.estimation import feols

data = get_data()
fit = feols("Y ~ C(f1)", data=data)
fit.confint(alpha=0.10).head()
fit.confint(alpha=0.10, inference_type="simult", reps=9999).head()

savi_fit = feols("Y ~ X1 + X2", data=data, vcov="hetero")
savi_fit.confint(alpha=0.10, inference_type="savi").head()
```

```text
| | 5.0% | 95.0% |
| --- | --- | --- |
| Intercept | 0.522146 | 1.255412 |
| X1 | -1.264565 | -0.721307 |
| X2 | -0.250161 | -0.102524 |
```

--------------------------------

### Generate test dataset with pyfixest

Source: https://pyfixest.org/how-to/stata-2-pyfixest.html

Use pyfixest's get_data() function to create a sample dataset for testing. This is the first step before running regressions.

```python
import pyfixest as pf
df = pf.get_data()
```

--------------------------------

### Compare 'fn' and 'pfn' solver methods

Source: https://pyfixest.org/tutorials/quantile-regression.html

Runs quantile regressions with the default 'fn' solver and the 'pfn' solver with pre-processing, setting a seed for reproducibility, then displays both results side by side with etable().

```python
%%capture
fit_fn = pf.quantreg(
  fml = "Y ~ X1",
  method = "fn",     # standard frisch newton interior point solver
  data = data,
)
fit_pfn = pf.quantreg(
  fml = "Y ~ X1",
  method = "pfn",   # standard frisch newton interior point solver with pre-processing
  seed = 92,         # set a seed for reproducibility
  data = data,
)

pf.etable([fit_fn, fit_pfn])
```

--------------------------------

### Generate LaTeX output with pf.etable(type='tex')

Source: https://pyfixest.org/tutorials/regression-tables.html

Use type='tex' to return a LaTeX table. The example sets digits=2. The output requires the LaTeX packages booktabs, threeparttable, makecell, and tabularx in your document. Use file_name= to save the table to a .tex file.

```python
# LaTex output (include latex packages booktabs, threeparttable, makecell, and tabularx in your document):
tab = pf.etable(
    [fit1, fit2, fit3, fit4, fit5, fit6],
    signif_code=[0.01, 0.05, 0.1],
    digits=2,
    type="tex",
)
```

--------------------------------

### List available pixi tasks

Source: https://pyfixest.org/contributing.html

Run `pixi task list` to see all available tasks and their descriptions.

```bash
pixi task list
```

--------------------------------

### Combining Different Estimation Methods with multi_method

Source: https://pyfixest.org/tutorials/quantile-regression.html

This example demonstrates how to combine different estimation methods (`fn` and `pfn`) with the same `multi_method` (`cfm1`). It fits two models with different seeds and quantiles, then displays the results.

```R
fit_multi2a = pf.quantreg(
  fml = "Y~X1",
  data = data,
  method = "fn",
  multi_method = "cfm1",
  seed = 233,
  quantile = [0.25, 0.75]
)

fit_multi2b = pf.quantreg(
  fml = "Y~X1",
  data = data,
  method = "pfn",
  multi_method = "cfm1",
  seed = 233,
  quantile = [0.25, 0.75]
)

pf.etable(fit_multi2a.to_list() +  fit_multi2b.to_list())
```

--------------------------------

### Example: i() with continuous interaction in feols formula

Source: https://pyfixest.org/reference/estimation.formula.transforms.factor_interaction.factor_interaction.html

Example of using the i() operator to interact a categorical variable with a continuous variable, giving group-specific slopes.

```python
pf.feols("Y ~ i(f1, X2, ref=0)", data).tidy().head()
```

--------------------------------

### Load sample heterogeneous DID dataset

Source: https://pyfixest.org/reference/did.estimation.did2s.html

Imports pandas, numpy, and pyfixest, then reads the heterogeneous treatment effects dataset from the pyfixest GitHub repository and displays the first few rows.

```python
import pandas as pd
import numpy as np
import pyfixest as pf

url = "https://raw.githubusercontent.com/py-econometrics/pyfixest/master/pyfixest/did/data/df_het.csv"
df_het = pd.read_csv(url)
df_het.head()
```

--------------------------------

### did2s() and feols() example for difference-in-differences

Source: https://pyfixest.org/changelog.html

Example from the PyFixest 0.10.8 changelog showing how to use the new did2s() function for two-stage difference-in-differences, and comparing it with a TWFE model using feols(). The example loads a dataset from a URL, fits both models, and plots the results with iplot(). Note the i_ref1 argument is used to set reference levels for the i() syntax.

```python
from pyfixest.did.did import did2s
from pyfixest.estimation import feols
from pyfixest.visualize import iplot
import pandas as pd
import numpy as np

df_het = pd.read_csv("https://raw.githubusercontent.com/py-econometrics/pyfixest/master/pyfixest/did/data/df_het.csv")

fit = did2s(
    df_het,
    yname = "dep_var",
    first_stage = "~ 0 | state + year",
    second_stage = "~i(rel_year)",
    treatment = "treat",
    cluster = "state",
    i_ref1 = [-1.0, np.inf],
)

fit_twfe = feols(
    "dep_var ~ i(rel_year) | state + year",
    df_het,
    i_ref1 = [-1.0, np.inf]
)

iplot([fit, fit_twfe], coord_flip=False, figsize = (900, 400), title = "TWFE vs DID2S")
```

--------------------------------

### LsmrDemeaner constructor with default parameters

Source: https://pyfixest.org/reference/demeaners.LsmrDemeaner.html

Shows the constructor signature with all parameters and their default values. The 'within' backend collapses fixef_atol and fixef_btol to max(fixef_atol, fixef_btol); the 'torch' backend uses both independently. The 'local_size' field only applies to backend='within'. The 'precision', 'device', and 'warn_on_cpu_fallback' fields are only relevant for the 'torch' backend.

```python
LsmrDemeaner(
    fixef_maxiter=1000,
    backend='within',
    precision='float64',
    device='auto',
    fixef_atol=1e-08,
    fixef_btol=1e-08,
    warn_on_cpu_fallback=True,
    preconditioner='auto',
    local_size=None,
)
```

--------------------------------

### Wald test example with pyfixest

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.wald_test.html

Example of using wald_test with a custom R matrix and q vector, specifying a chi-squared distribution. The code also accesses the _f_statistic and _p_value attributes and prints them.

```python
import numpy as np
import pandas as pd
import pyfixest as pf

data = pf.get_data()
fit = pf.feols("Y ~ X1 + X2| f1", data, vcov={"CRV1": "f1"}, ssc=pf.ssc(k_adj=False))

R = np.array([[1,-1]] )
q = np.array([0.0])

# Wald test
fit.wald_test(R=R, q=q, distribution = "chi2")
f_stat = fit._f_statistic
p_stat = fit._p_value

print(f"Python f_stat: {f_stat}")
print(f"Python p_stat: {p_stat}")
```

```output
Python f_stat: 126.40650474043491
Python p_stat: 2.505309282814043e-29
```

--------------------------------

### Example: i() with reference level in feols formula

Source: https://pyfixest.org/reference/estimation.formula.transforms.factor_interaction.factor_interaction.html

Example of using the i() operator in a feols formula to expand a categorical variable into indicators, dropping a reference level. This is commonly used for event study specifications.

```python
import pyfixest as pf

data = pf.get_data()

fit = pf.feols("Y ~ i(f1, ref=0)", data)
fit.tidy().head()
```

--------------------------------

### Pickle and reuse a within preconditioner across sessions

Source: https://pyfixest.org/how-to/demeaner-backends.html

Builds a preconditioner, pickles it to a temp file, loads it in a new session, and reuses it in a fit. Includes a warning about dimension and structure validation.

```python
import pickle
import tempfile, os

# session 1: build and save
fit = pf.feols("Y ~ X1 | f1 + f2", data=df, demeaner=LsmrDemeaner())
tmp = os.path.join(tempfile.gettempdir(), "precond.pkl")
with open(tmp, "wb") as f:
    pickle.dump(fit.preconditioner, f)

# session 2: load and reuse
with open(tmp, "rb") as f:
    pre = pickle.load(f)

fit2 = pf.feols(
    "Y ~ X1 | f1 + f2",
    data=df,
    demeaner=LsmrDemeaner(preconditioner=pre),
)

pf.etable([fit, fit2], digits = 6)
```

--------------------------------

### Feols.update usage example

Source: https://pyfixest.org/reference/estimation.models.feols_.Feols.html

Example of updating a fitted model with a new observation. Fits on all but the last observation, then adds the last observation using the update method. X_new must include the intercept column.

```python
import numpy as np
import pyfixest as pf

data = pf.get_data().dropna()
fit = pf.feols("Y ~ X1 + X2", data.iloc[:-1])

last = data.iloc[[-1]]
X_new = np.column_stack(
    [np.ones(1), last["X1"].to_numpy(), last["X2"].to_numpy()]
)
y_new = last["Y"].to_numpy()

fit.update(X_new, y_new)
```

--------------------------------

### NYC-Taxi performance regression example with pyfixest

Source: https://pyfixest.org/changelog.html

Shows a full workflow for loading NYC taxi data with DuckDB, converting columns to categorical/int types, and running a fixed-effects OLS regression with `pf.feols()`. Includes a mock regression for JIT compilation and uses `copy_data=False` and `store_data=False` to save time. The output indicates the regression takes about 22 seconds.

```python
%load_ext autoreload
%autoreload 2

import duckdb
import time
import numpy as np
import pyfixest as pf

# %%
nyc = duckdb.sql(
    '''
    FROM 'C:/Users/alexa/Documents/nyc-taxi/**/*.parquet'
    SELECT
        tip_amount, trip_distance, passenger_count,
        vendor_id, payment_type, dropoff_at,
        dayofweek(dropoff_at) AS dofw
    WHERE year = 2012 AND month <= 3
    '''
    ).df()

# convert dowf, vendor_id, payment_type to categorical
tic = time.time()
nyc["dofw"] = nyc["dofw"].astype(int)
nyc["vendor_id"] = nyc["vendor_id"].astype("category")
nyc["payment_type"] = nyc["payment_type"].astype("category")
print(f"""
    I am convering columns of type 'objects' to 'categories' and 'int'data types outside
    of the regression, hence I am cheating a bit. This saves {np.round(time.time() - tic)} seconds.
    """
)
#    I am convering columns of type 'objects' to 'categories' and 'int'data types outside
#    of the regression, hence I am cheating a bit. This saves 7.0 seconds.

run = True
if run:

    # mock regression for JIT compilation
    fit = pf.feols(
        fml = "tip_amount ~ trip_distance + passenger_count | vendor_id + payment_type + dofw",
        data = nyc.iloc[1:10_000],
        copy_data = False,
        store_data = False
        )

    import time
    tic = time.time()
    fit = pf.feols(
        fml = "tip_amount ~ trip_distance + passenger_count | vendor_id + payment_type + dofw",
        data = nyc,
        copy_data = False, # saves a few seconds
        store_data = False # saves a few second
        )
    passed = time.time() - tic
    print(f"Passed time is {np.round(passed)}.")
    # Passed time is 22.
```

--------------------------------

### Import PyFixest

Source: https://pyfixest.org/tutorials/instrumental-variables.html

Import the PyFixest library to begin using its estimation functions.

```python
import pyfixest as pf
```

--------------------------------

### Feols.get_performance

Source: https://pyfixest.org/reference/

Get Goodness-of-Fit measures.

```APIDOC
## Feols.get_performance

### Description
Get Goodness-of-Fit measures.

### Method
Method call on a fitted model object.

### Endpoint
Feols.get_performance()

### Parameters
None documented.

### Response
Returns goodness-of-fit measures.
```

--------------------------------

### etable() for sw() models

Source: https://pyfixest.org/reference/estimation.api.feols.feols.html

Displays the comparison table for the two models estimated with sw() syntax. The table shows coefficients, fixed effects, and fit statistics for each model.

```python
pf.etable(fit)
```

--------------------------------

### Feglm usage example with logit family

Source: https://pyfixest.org/reference/estimation.models.feglm_.Feglm.html

Fits a fixed-effects logit model using pf.feglm and prints the tidy coefficient table. The example uses synthetic data from pf.get_data() and creates a binary outcome Y_bin.

```python
import numpy as np
import pyfixest as pf

data = pf.get_data()
data["Y_bin"] = np.where(data["Y"] > 0, 1, 0)

fit = pf.feglm("Y_bin ~ X1 + X2 | f1", data, family="logit")
fit.tidy()
```

--------------------------------

### Example: demeaning data with pyfixest

Source: https://pyfixest.org/reference/estimation.demean.html

Example usage of demean with pyfixest and formulaic. It creates a model matrix, converts fixed effects to integers, and calls demean with weights. The output shows the coefficient from a least squares regression and the pyfixest feols result.

```python
import numpy as np
import pyfixest as pf
from pyfixest.utils.dgps import get_blw
from pyfixest.estimation import demean
from formulaic import model_matrix

fml = "y ~ treat | state + year"

data = get_blw()
data.head()

Y, rhs = model_matrix(fml, data)
X = rhs[0].drop(columns="Intercept")
fe = rhs[1].drop(columns="Intercept")
YX = np.concatenate([Y, X], axis=1)

# to numpy
Y = Y.to_numpy()
X = X.to_numpy()
YX = np.concatenate([Y, X], axis=1)
fe = fe.to_numpy().astype(int)  # demean requires fixed effects as ints!

YX_demeaned, success = demean(YX, fe, weights = np.ones(YX.shape[0]))
Y_demeaned = YX_demeaned[:, 0]
X_demeaned = YX_demeaned[:, 1:]

print(np.linalg.lstsq(X_demeaned, Y_demeaned, rcond=None)[0])
print(pf.feols(fml, data).coef())
```

```text
[-6.73780555]
Coefficient
treat   -6.737806
Name: Estimate, dtype: float64
```

--------------------------------

### Run all linting hooks manually

Source: https://pyfixest.org/contributing.html

Run `pixi run lint` to execute all linting and formatting hooks manually.

```bash
pixi run lint
```