Proposal for Caching the Within-Preconditioner

Cargo.toml

@@ -15,6 +15,7 @@ rayon   = "1.11.0"
 numpy = "0.26.0"
 thiserror = "2.0.16"
 within = "0.1.0"
+postcard = { version = "1", features = ["use-std"] }
 
 [profile.release]
 opt-level = 3        # Maximize performance

docs/_freeze/how-to/fixed-effects-solvers/execute-results/html.json

@@ -0,0 +1,16 @@
+{
+  "hash": "2e454435a33d9f6487d50ede3e986926",
+  "result": {
+    "engine": "jupyter",
+    "markdown": "---\ntitle: \"Choosing and Reusing Fixed-Effects Solvers\"\ndescription: \"How to configure fixed-effects demeaners, reuse within solvers and preconditioners, and persist reusable preconditioners across sessions.\"\ncategories: [Fixed Effects, Performance]\n---\n\n::: {#661c1809 .cell execution_count=1}\n``` {.python .cell-code}\nimport pickle\nimport tempfile\nfrom pathlib import Path\n\nimport numpy as np\nimport pandas as pd\n\nimport pyfixest as pf\n```\n:::\n\n\n`PyFixest` now accepts object-based demeaner configurations. This makes it easier to choose a backend, tune solver settings, and pass reusable solver objects between regressions.\n\nIf you are deciding whether a fixed-effects problem is likely to be easy or hard for alternating-projections methods, see [When Are Fixed Effects Problems Difficult?](../explanation/difficult-fixed-effects.md).\n\n## Choosing a Demeaner\n\nThe main fixed-effects options are:\n\n| Demeaner | Main idea | Typical use |\n|---|---|---|\n| `pf.NumbaDemeaner(...)` | CPU alternating projections | good default on CPU |\n| `pf.WithinDemeaner(...)` | Krylov solver with Schwarz preconditioning via `within` | hard FE problems and reuse workflows |\n| `pf.ScipyDemeaner(...)` | SciPy sparse iterative solver | CPU sparse linear algebra workflows |\n| `pf.CupyDemeaner(...)` / `pf.JaxDemeaner(...)` | GPU-based demeaning | large GPU workloads |\n\nThe `within` backend is the only one that currently exposes reusable `solver` and `preconditioner` objects.\n\n## Example Data\n\n::: {#08c7114a .cell execution_count=2}\n``` {.python .cell-code}\npanel = pf.get_worker_panel(N_workers=300, N_firms=40, N_years=8, seed=42)\n\nrng = np.random.default_rng(123)\npanel[\"log_bonus\"] = (\n    0.4 * panel[\"worker_fe\"]\n    - 0.2 * panel[\"firm_fe\"]\n    + 0.01 * panel[\"experience\"]\n    + 0.02 * panel[\"tenure\"]\n    + rng.normal(0, 0.2, size=len(panel))\n)\n\nfml_y1 = \"log_wage ~ experience + tenure | worker_id + firm_id\"\nfml_y2 = \"log_bonus ~ experience + tenure | worker_id + firm_id\"\npanel.head()\n```\n\n::: {.cell-output .cell-output-display execution_count=2}\n```{=html}\n<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>worker_id</th>\n      <th>firm_id</th>\n      <th>year</th>\n      <th>female</th>\n      <th>experience</th>\n      <th>tenure</th>\n      <th>log_wage</th>\n      <th>worker_fe</th>\n      <th>firm_fe</th>\n      <th>log_bonus</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>0</th>\n      <td>0</td>\n      <td>16</td>\n      <td>2000</td>\n      <td>0</td>\n      <td>1</td>\n      <td>1</td>\n      <td>0.028921</td>\n      <td>0.152359</td>\n      <td>0.019964</td>\n      <td>-0.110874</td>\n    </tr>\n    <tr>\n      <th>1</th>\n      <td>1</td>\n      <td>34</td>\n      <td>2000</td>\n      <td>1</td>\n      <td>0</td>\n      <td>1</td>\n      <td>-0.526312</td>\n      <td>-0.519992</td>\n      <td>0.009835</td>\n      <td>-0.263521</td>\n    </tr>\n    <tr>\n      <th>2</th>\n      <td>2</td>\n      <td>20</td>\n      <td>2000</td>\n      <td>1</td>\n      <td>3</td>\n      <td>1</td>\n      <td>0.515641</td>\n      <td>0.375226</td>\n      <td>0.234705</td>\n      <td>0.410734</td>\n    </tr>\n    <tr>\n      <th>3</th>\n      <td>3</td>\n      <td>13</td>\n      <td>2000</td>\n      <td>1</td>\n      <td>0</td>\n      <td>1</td>\n      <td>0.528749</td>\n      <td>0.470282</td>\n      <td>0.131698</td>\n      <td>0.220568</td>\n    </tr>\n    <tr>\n      <th>4</th>\n      <td>4</td>\n      <td>31</td>\n      <td>2000</td>\n      <td>0</td>\n      <td>3</td>\n      <td>1</td>\n      <td>-0.337822</td>\n      <td>-0.975518</td>\n      <td>0.191260</td>\n      <td>-0.194413</td>\n    </tr>\n  </tbody>\n</table>\n</div>\n```\n:::\n:::\n\n\n## Pass Demeaners as Dataclasses\n\nInstead of string shorthands, you can pass a configured demeaner directly.\n\n::: {#97c68471 .cell execution_count=3}\n``` {.python .cell-code}\nfit_numba = pf.feols(\n    fml_y1,\n    data=panel,\n    demeaner=pf.NumbaDemeaner(fixef_maxiter=10_000),\n)\n\nfit_within = pf.feols(\n    fml_y1,\n    data=panel,\n    demeaner=pf.WithinDemeaner(\n        krylov_method=\"gmres\",\n        gmres_restart=20,\n        preconditioner_type=\"multiplicative\",\n    ),\n)\n\npd.concat(\n    {\n        \"numba\": fit_numba.coef(),\n        \"within\": fit_within.coef(),\n    },\n    axis=1,\n)\n```\n\n::: {.cell-output .cell-output-display}\n```{=html}\n\n        <div id=\"Yb5WqZ\"></div>\n        <script type=\"text/javascript\" data-lets-plot-script=\"library\">\n            if(!window.letsPlotCallQueue) {\n                window.letsPlotCallQueue = [];\n            };\n            window.letsPlotCall = function(f) {\n                window.letsPlotCallQueue.push(f);\n            };\n            (function() {\n                var script = document.createElement(\"script\");\n                script.type = \"text/javascript\";\n                script.src = \"https://cdn.jsdelivr.net/gh/JetBrains/lets-plot@v4.8.2/js-package/distr/lets-plot.min.js\";\n                script.onload = function() {\n                    window.letsPlotCall = function(f) {f();};\n                    window.letsPlotCallQueue.forEach(function(f) {f();});\n                    window.letsPlotCallQueue = [];\n                    \n                };\n                script.onerror = function(event) {\n                    window.letsPlotCall = function(f) {};    // noop\n                    window.letsPlotCallQueue = [];\n                    var div = document.createElement(\"div\");\n                    div.style.color = 'darkred';\n                    div.textContent = 'Error loading Lets-Plot JS';\n                    document.getElementById(\"Yb5WqZ\").appendChild(div);\n                };\n                var e = document.getElementById(\"Yb5WqZ\");\n                e.appendChild(script);\n            })()\n        </script>\n        \n```\n:::\n\n::: {.cell-output .cell-output-display}\n```{=html}\n\n        <div id=\"BsRUIo\"></div>\n        <script type=\"text/javascript\" data-lets-plot-script=\"library\">\n            if(!window.letsPlotCallQueue) {\n                window.letsPlotCallQueue = [];\n            };\n            window.letsPlotCall = function(f) {\n                window.letsPlotCallQueue.push(f);\n            };\n            (function() {\n                var script = document.createElement(\"script\");\n                script.type = \"text/javascript\";\n                script.src = \"https://cdn.jsdelivr.net/gh/JetBrains/lets-plot@v4.8.2/js-package/distr/lets-plot.min.js\";\n                script.onload = function() {\n                    window.letsPlotCall = function(f) {f();};\n                    window.letsPlotCallQueue.forEach(function(f) {f();});\n                    window.letsPlotCallQueue = [];\n                    \n                };\n                script.onerror = function(event) {\n                    window.letsPlotCall = function(f) {};    // noop\n                    window.letsPlotCallQueue = [];\n                    var div = document.createElement(\"div\");\n                    div.style.color = 'darkred';\n                    div.textContent = 'Error loading Lets-Plot JS';\n                    document.getElementById(\"BsRUIo\").appendChild(div);\n                };\n                var e = document.getElementById(\"BsRUIo\");\n                e.appendChild(script);\n            })()\n        </script>\n        \n```\n:::\n\n::: {.cell-output .cell-output-display execution_count=3}\n```{=html}\n<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>numba</th>\n      <th>within</th>\n    </tr>\n    <tr>\n      <th>Coefficient</th>\n      <th></th>\n      <th></th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>experience</th>\n      <td>0.018400</td>\n      <td>0.018400</td>\n    </tr>\n    <tr>\n      <th>tenure</th>\n      <td>0.016517</td>\n      <td>0.016517</td>\n    </tr>\n  </tbody>\n</table>\n</div>\n```\n:::\n:::\n\n\n## Reuse a Solver Within the Same Session\n\nFor in-process reuse, prefer a `solver`. This is the fast path: the expensive setup is done once, and later solves borrow the stored solver directly.\n\nYou can reuse a solver exposed by a fitted model:\n\n::: {#9d0f1071 .cell execution_count=4}\n``` {.python .cell-code}\nfit1 = pf.feols(\n    fml_y1,\n    data=panel,\n    demeaner=pf.WithinDemeaner(),\n)\n\nfit2 = pf.feols(\n    fml_y2,\n    data=panel,\n    demeaner=pf.WithinDemeaner(solver=fit1.solver_),\n)\n\nfit2.coef()\n```\n\n::: {.cell-output .cell-output-display execution_count=4}\n```\nCoefficient\nexperience    0.011953\ntenure        0.016365\nName: Estimate, dtype: float64\n```\n:::\n:::\n\n\nYou can also build a reusable solver explicitly:\n\n::: {#6e60c279 .cell execution_count=5}\n``` {.python .cell-code}\nsolver = pf.get_solver(\n    fml_y1,\n    data=panel,\n    demeaner=pf.WithinDemeaner(preconditioner_type=\"additive\"),\n)\n\nfit3 = pf.feols(\n    fml_y2,\n    data=panel,\n    demeaner=pf.WithinDemeaner(solver=solver),\n)\n\nfit3.coef()\n```\n\n::: {.cell-output .cell-output-display execution_count=5}\n```\nCoefficient\nexperience    0.011953\ntenure        0.016365\nName: Estimate, dtype: float64\n```\n:::\n:::\n\n\nUse a solver when the sample, weights, and fixed-effects encoding are the same. Solver reuse is intentionally strict.\n\n## Reuse a Preconditioner Across Regressions\n\nIf you want a more portable object, extract a `preconditioner` from a solver and pass that into a later regression:\n\n::: {#8a613a85 .cell execution_count=6}\n``` {.python .cell-code}\npreconditioner_from_fit = fit1.solver_.to_preconditioner()\n\nfit4 = pf.feols(\n    fml_y2,\n    data=panel,\n    demeaner=pf.WithinDemeaner(preconditioner=preconditioner_from_fit),\n)\n\nfit4.coef()\n```\n\n::: {.cell-output .cell-output-display execution_count=6}\n```\nCoefficient\nexperience    0.011953\ntenure        0.016365\nName: Estimate, dtype: float64\n```\n:::\n:::\n\n\nCompared with a solver, a preconditioner is a lower-level object. It only affects convergence speed, not the target problem itself.\n\n## Build a Preconditioner Explicitly\n\nYou can also construct a reusable preconditioner directly:\n\n::: {#39988af0 .cell execution_count=7}\n``` {.python .cell-code}\npreconditioner = pf.get_preconditioner(\n    fml_y1,\n    data=panel,\n    demeaner=pf.WithinDemeaner(preconditioner_type=\"additive\"),\n)\n\nfit5 = pf.feols(\n    fml_y2,\n    data=panel,\n    demeaner=pf.WithinDemeaner(preconditioner=preconditioner),\n)\n\nfit5.coef()\n```\n\n::: {.cell-output .cell-output-display execution_count=7}\n```\nCoefficient\nexperience    0.011953\ntenure        0.016365\nName: Estimate, dtype: float64\n```\n:::\n:::\n\n\nThis is the preferred API when you want a reusable object for persistence or for advanced reuse workflows.\n\n## Persist a Preconditioner Across Sessions\n\nPreconditioners are pickleable, so you can save them and reload them in a later session:\n\n::: {#e1356e5c .cell execution_count=8}\n``` {.python .cell-code}\nwith tempfile.TemporaryDirectory() as tmpdir:\n    path = Path(tmpdir) / \"within_preconditioner.pkl\"\n\n    with path.open(\"wb\") as file:\n        pickle.dump(preconditioner, file)\n\n    with path.open(\"rb\") as file:\n        restored_preconditioner = pickle.load(file)\n\n    print(path)\n\nfit6 = pf.feols(\n    fml_y2,\n    data=panel,\n    demeaner=pf.WithinDemeaner(preconditioner=restored_preconditioner),\n)\n\nfit6.coef()\n```\n\n::: {.cell-output .cell-output-stdout}\n```\n/var/folders/98/c353q4p95v5_fz62gr5c9td80000gn/T/tmpa3kwwj68/within_preconditioner.pkl\n```\n:::\n\n::: {.cell-output .cell-output-display execution_count=8}\n```\nCoefficient\nexperience    0.011953\ntenure        0.016365\nName: Estimate, dtype: float64\n```\n:::\n:::\n\n\nIn practice, replace the temporary file with a durable path.\n\n## Which Reuse Object Should I Use?\n\n- Use `solver` for repeated estimation in the same Python session.\n- Use `preconditioner` when you want a portable, pickleable object.\n- Use `fit.solver_.to_preconditioner()` when you want to turn a fitted model's solver into a persistence artifact.\n- Use `pf.get_preconditioner()` when you want to build the persistence artifact explicitly.\n\nReusable solvers are not supported for `feglm()` and `fepois()`, because those estimators update observation weights across iterations. Reusable preconditioners can still be used there.\n\n",
+    "supporting": [
+      "fixed-effects-solvers_files"
+    ],
+    "filters": [],
+    "includes": {
+      "include-in-header": [
+        "<script src=\"https://cdn.jsdelivr.net/npm/requirejs@2.3.6/require.min.js\" integrity=\"sha384-c9c+LnTbwQ3aujuU7ULEPVvgLs+Fn6fJUvIGTsuu1ZcCf11fiEubah0ttpca4ntM sha384-6V1/AdqZRWk1KAlWbKBlGhN7VG4iE/yAZcO6NZPMF8od0vukrvr0tg4qY6NSrItx\" crossorigin=\"anonymous\"></script>\n<script src=\"https://cdn.jsdelivr.net/npm/jquery@3.5.1/dist/jquery.min.js\" integrity=\"sha384-ZvpUoO/+PpLXR1lu4jmpXWu80pZlYUAfxl5NsBMWOEPSjUn/6Z/hRTt8+pR6L4N2\" crossorigin=\"anonymous\" data-relocate-top=\"true\"></script>\n<script type=\"application/javascript\">define('jquery', [],function() {return window.jQuery;})</script>\n"
+      ]
+    }
+  }
+}

docs/_quarto.yml

@@ -72,6 +72,8 @@ website:
         - text: "All Guides"
           file: how-to/index.qmd
         - text: "---"
+        - text: "Choosing and Reusing Fixed-Effects Solvers"
+          file: how-to/fixed-effects-solvers.qmd
         - text: "Marginal Effects & Hypothesis Testing"
           file: how-to/marginaleffects.qmd
         - text: "Regression Decomposition"