pyfector

pyfector: counterfactual estimators for panel data in Python.

A Python implementation of counterfactual estimators for panel data, following Liu, Wang and Xu (2024). See the R fect package for the R implementation.

Usage
    import pyfector

    result = pyfector.fect(
        data=df,
        Y="outcome", D="treat",
        index=("unit", "year"),
        method="ife",
        r=(0, 5),
        se=True,
        seed=42,
    )
    result.summary()
    result.plot()
API modules

pyfector.diagnostics provides pre-trend tests, holdout/refit placebo tests, and post-treatment sensitivity diagnostics. pyfector.inference provides bootstrap, jackknife, and permutation inference. pyfector.panel, pyfector.estimators, and pyfector.cv describe panel construction, estimation, and cross-validation.

 1"""
 2pyfector: counterfactual estimators for panel data in Python.
 3
 4A Python implementation of counterfactual estimators for panel data,
 5following Liu, Wang and Xu (2024). See the
 6[R fect package](https://github.com/xuyiqing/fect) for the R implementation.
 7
 8Usage
 9-----
10
11```python
12
13    import pyfector
14
15    result = pyfector.fect(
16        data=df,
17        Y="outcome", D="treat",
18        index=("unit", "year"),
19        method="ife",
20        r=(0, 5),
21        se=True,
22        seed=42,
23    )
24    result.summary()
25    result.plot()
26```
27
28API modules
29-----------
30`pyfector.diagnostics` provides pre-trend tests, holdout/refit placebo tests,
31and post-treatment sensitivity diagnostics. `pyfector.inference` provides
32bootstrap, jackknife, and permutation inference. `pyfector.panel`,
33`pyfector.estimators`, and `pyfector.cv` describe panel construction,
34estimation, and cross-validation.
35
36"""
37
38__version__ = "0.2.3"
39
40from .fect import fect, FectResult
41from .backend import set_device, get_device
42from .diagnostics import (
43    run_diagnostics,
44    Diagnostics,
45    DiagnosticResult,
46    TostResult,
47    PretrendFResult,
48    EquivFResult,
49    PlaceboResult,
50    CarryoverResult,
51    LooResult,
52)
53from .plotting import plot
54
55__all__ = [
56    "fect",
57    "FectResult",
58    "set_device",
59    "get_device",
60    "run_diagnostics",
61    "Diagnostics",
62    "DiagnosticResult",
63    "TostResult",
64    "PretrendFResult",
65    "EquivFResult",
66    "PlaceboResult",
67    "CarryoverResult",
68    "LooResult",
69    "plot",
70]
def fect( data, Y: str, D: str, index: tuple[str, str], X: list[str] | None = None, W: str | None = None, group: str | None = None, method: Literal['fe', 'ife', 'mc', 'cfe', 'both'] = 'ife', force: Literal['none', 'unit', 'time', 'two-way'] = 'two-way', r: int | tuple[int, int] = 0, lam: float | None = None, nlambda: int = 10, lambda_candidates: list[float] | numpy.ndarray | None = None, CV: bool = True, k: int = 10, cv_prop: float = 0.1, cv_nobs: int = 3, cv_treat: bool = True, cv_donut: int = 0, criterion: str = 'mspe', cv_rule: Literal['min', 'onepct'] = 'min', se: bool = False, vartype: Literal['bootstrap', 'jackknife'] = 'bootstrap', nboots: int = 200, alpha: float = 0.05, tol: float = 1e-07, max_iter: int = 5000, min_T0: int = 1, min_T0_strict: bool = False, max_missing: float = 1.0, normalize: bool = False, Z: list[str] | None = None, Q: list[str] | None = None, device: Literal['cpu', 'gpu'] = 'cpu', n_jobs: int | None = -1, seed: int | None = None, diagnostics: Union[Literal['none', 'full'], list[str]] = 'none', diagnostics_options: dict | None = None, *, cluster: str | None = None) -> FectResult: 
233def fect(
234    data,
235    Y: str,
236    D: str,
237    index: tuple[str, str],
238    X: list[str] | None = None,
239    W: str | None = None,
240    group: str | None = None,
241    method: Literal["fe", "ife", "mc", "cfe", "both"] = "ife",
242    force: Literal["none", "unit", "time", "two-way"] = "two-way",
243    r: int | tuple[int, int] = 0,
244    lam: float | None = None,
245    nlambda: int = 10,
246    lambda_candidates: list[float] | np.ndarray | None = None,
247    CV: bool = True,
248    k: int = 10,
249    cv_prop: float = 0.1,
250    cv_nobs: int = 3,
251    cv_treat: bool = True,
252    cv_donut: int = 0,
253    criterion: str = "mspe",
254    cv_rule: Literal["min", "onepct"] = "min",
255    se: bool = False,
256    vartype: Literal["bootstrap", "jackknife"] = "bootstrap",
257    nboots: int = 200,
258    alpha: float = 0.05,
259    tol: float = 1e-7,
260    max_iter: int = 5000,
261    min_T0: int = 1,
262    min_T0_strict: bool = False,
263    max_missing: float = 1.0,
264    normalize: bool = False,
265    # CFE-specific
266    Z: list[str] | None = None,
267    Q: list[str] | None = None,
268    # Performance
269    device: Literal["cpu", "gpu"] = "cpu",
270    n_jobs: int | None = -1,
271    seed: int | None = None,
272    # Diagnostics (optional; run at fit time and attached to result)
273    diagnostics: Literal["none", "full"] | list[str] = "none",
274    diagnostics_options: dict | None = None,
275    *,
276    cluster: str | None = None,
277) -> FectResult:
278    """Estimate counterfactual treatment effects for panel data.
279
280    This is the main Python entry point for the counterfactual estimator
281    workflow.  Where the paper and the historical R package differ,
282    pyfector defaults to the paper's statistical definition and exposes
283    R-package-style behavior through explicit options.
284
285    Missing outcome policy
286    ----------------------
287    pyfector distinguishes raw missing outcomes from counterfactual
288    missingness caused by treatment.  Observed untreated cells
289    (``D == 0`` and non-missing ``Y``) fit the response surface.  Observed
290    treated cells (``D == 1`` and non-missing ``Y``) contribute to ATT as
291    ``Y - Y_ct``.  If a treated outcome is missing in the input data, the
292    model can still produce a counterfactual ``Y_ct`` for that cell, but
293    the cell is not counted in ``att_avg`` or ``att_on`` because the
294    treated potential outcome was not observed.
295
296    By default, ``min_T0`` is enforced only for treated and reversal
297    units.  Sparse controls are retained if they have at least one
298    observed outcome, because they may still inform the low-rank response
299    surface.  Set ``min_T0_strict=True`` to require controls to satisfy
300    ``min_T0`` too, matching the more conservative R fect sparse-panel
301    behavior.
302
303    Parameters
304    ----------
305    data : polars.DataFrame, pandas.DataFrame
306        Long-format panel data.
307    Y, D : str
308        Column names for outcome and binary treatment indicator.
309    index : (str, str)
310        Column names for (unit_id, time_period).
311    X : list of str, optional
312        Time-varying covariates.
313    W : str, optional
314        Observation weight column.
315    group : str, optional
316        Reserved for grouped estimation. Currently raises
317        ``NotImplementedError`` when supplied.
318    method : {"fe", "ife", "mc", "cfe", "both"}
319        Estimation method.
320    force : {"none", "unit", "time", "two-way"}
321        Fixed effects specification.
322    r : int or (int, int)
323        Number of factors.  If tuple, CV selects from range.
324    lam : float, optional
325        Nuclear norm penalty for MC.  If None with CV=True, auto-selected.
326    nlambda : int
327        Number of automatically generated lambda candidates for MC CV.
328    lambda_candidates : array-like, optional
329        Explicit non-negative lambda candidates for MC CV. When supplied,
330        ``nlambda`` is ignored.
331    CV : bool
332        If True, cross-validate over ``r`` for IFE when ``r`` is a tuple,
333        or over ``lam`` for MC when ``lam`` is None.
334    k : int
335        Number of CV folds.
336    cv_prop : float
337        Fraction of eligible observed control cells masked per CV fold.
338    cv_nobs : int
339        Number of consecutive within-unit observations to mask as a block.
340    cv_treat : bool
341        If True, restrict CV masks to pre-treatment cells of ever-treated
342        units. If False, use all observed control cells.
343    cv_donut : int
344        Exclude this many periods around treatment onset from CV evaluation.
345    criterion : {"mspe", "gmspe", "mad"}
346        Cross-validation loss.
347    cv_rule : {"min", "onepct"}
348        CV selection rule. ``"min"`` chooses the strict minimum-score
349        candidate and is the paper-faithful default. ``"onepct"`` chooses
350        the simplest candidate within 1% of the best score (lower ``r`` for
351        IFE, higher ``lam`` for MC).
352    se : bool
353        Compute standard errors via bootstrap/jackknife.
354    vartype : {"bootstrap", "jackknife"}
355        Inference method when ``se=True``.
356    cluster : str, optional
357        Column of cluster labels, constant within each unit. Bootstrap
358        inference resamples whole clusters in both the main fit and the
359        placebo refit. Requires at least two retained clusters.
360    nboots : int
361        Number of bootstrap replications, at least 2. Ignored for jackknife.
362    alpha : float
363        Significance level for confidence intervals and tests.
364    tol : float
365        EM convergence tolerance for final point estimation.
366    max_iter : int
367        Maximum EM iterations.
368    min_T0 : int
369        Minimum untreated/pre-treatment observed periods. By default this is
370        enforced only for treated and treatment-reversal units.
371    min_T0_strict : bool
372        If True, enforce ``min_T0`` on all units, including controls. This
373        matches R fect's conservative handling of sparse control rows.
374    max_missing : float
375        Maximum missing-outcome fraction per unit, in ``[0, 1]``. Units with
376        no observed outcomes are always dropped, regardless of this threshold,
377        because they provide neither fitting information nor observed treated
378        effects.
379    normalize : bool
380        If True, estimate on an outcome standardized by its observed standard
381        deviation, then transform effects back to the original scale.
382    Z, Q : list of str, optional
383        Reserved CFE interaction arguments. Currently raise
384        ``NotImplementedError`` when supplied.
385    device : {"cpu", "gpu"}
386        Compute device.
387    n_jobs : int, optional
388        Parallel workers for CV and bootstrap. ``-1`` or ``None`` uses
389        all available CPUs.
390    seed : int, optional
391        Random seed for full reproducibility.
392
393    Returns
394    -------
395    FectResult
396        Treatment effects, fitted counterfactuals, and requested inference
397        and diagnostics.
398    """
399    # Set device
400    set_device(device)
401    xp = get_backend()
402    n_jobs = _resolve_n_jobs(n_jobs)
403    if device == "gpu":
404        n_jobs = 1
405
406    if group is not None:
407        raise NotImplementedError("The `group` argument is not implemented yet.")
408    if Z is not None or Q is not None:
409        raise NotImplementedError("The `Z` and `Q` CFE interaction arguments are not implemented yet.")
410    if criterion not in {"mspe", "gmspe", "mad"}:
411        raise ValueError("criterion must be 'mspe', 'gmspe', or 'mad'")
412    if cv_rule not in {"min", "onepct"}:
413        raise ValueError("cv_rule must be 'min' or 'onepct'")
414    if min_T0 < 0:
415        raise ValueError("min_T0 must be non-negative")
416    if cluster is not None and vartype != "bootstrap":
417        raise ValueError("cluster is supported only with vartype='bootstrap'")
418    if se and vartype == "bootstrap" and (
419        isinstance(nboots, (bool, np.bool_))
420        or not isinstance(nboots, (int, np.integer)) or nboots < 2
421    ):
422        raise ValueError("nboots must be an integer of at least 2")
423    if not 0.0 <= max_missing <= 1.0:
424        raise ValueError("max_missing must be between 0 and 1")
425
426    # Validate diagnostics request before doing any expensive estimation
427    # so users don't wait for a 30-min MC fit only to find their config
428    # is wrong.
429    requested_diag = validate_diagnostics_request(
430        diagnostics, diagnostics_options, se,
431    )
432
433    # Map force string to int
434    force_map = {"none": 0, "unit": 1, "time": 2, "two-way": 3}
435    force_int = force_map[force]
436
437    # Prepare panel data
438    panel = prepare_panel(
439        data, Y=Y, D=D, index=index, X=X, W=W,
440        group=group, min_T0=min_T0, min_T0_strict=min_T0_strict,
441        max_missing=max_missing,
442        cluster=cluster,
443    )
444
445    if not np.any((panel.D > 0) & (panel.I > 0)):
446        raise ValueError(
447            "No observed treated units remain after filtering. "
448            "Check treatment, min_T0, and missing-outcome filters."
449        )
450    if se and panel.cluster is not None and len(np.unique(panel.cluster)) < 2:
451        raise ValueError("cluster bootstrap requires at least two clusters after filtering")
452
453    # Move to device
454    Y_mat = to_device(panel.Y)
455    D_mat = to_device(panel.D)
456    I_mat = to_device(panel.I)
457    II_mat = to_device(panel.II)
458    X_mat = to_device(panel.X) if panel.X is not None else None
459    W_mat = to_device(panel.W) if panel.W is not None else None
460
461    # Normalize
462    norm_factor = 1.0
463    if normalize:
464        sd_y = float(xp.std(Y_mat[I_mat > 0]))
465        if sd_y > 0:
466            Y_mat = Y_mat / sd_y
467            norm_factor = sd_y
468
469    # Initial fit
470    Y0, beta0 = initial_fit(Y_mat, X_mat, II_mat, force_int)
471
472    # Determine r and lambda
473    r_cv = None
474    lambda_cv = None
475    cv_result = None
476
477    if method == "ife":
478        if isinstance(r, tuple) and CV:
479            cv_result = cv_ife(
480                Y_mat, Y0, X_mat, I_mat, II_mat, D_mat, W_mat, beta0,
481                force=force_int, r_range=r, k=k, cv_prop=cv_prop,
482                cv_nobs=cv_nobs, cv_treat=cv_treat, cv_donut=cv_donut,
483                criterion=criterion, cv_rule=cv_rule,
484                tol=tol, max_iter=max_iter,
485                n_jobs=n_jobs, seed=seed,
486            )
487            r_cv = cv_result.best_r
488        else:
489            r_cv = r if isinstance(r, int) else r[0]
490
491    elif method == "mc":
492        if lam is None and CV:
493            cv_result = cv_mc(
494                Y_mat, Y0, X_mat, I_mat, II_mat, D_mat, W_mat, beta0,
495                force=force_int, lambda_candidates=lambda_candidates,
496                nlambda=nlambda, k=k, cv_prop=cv_prop,
497                cv_nobs=cv_nobs, cv_treat=cv_treat, cv_donut=cv_donut,
498                criterion=criterion, cv_rule=cv_rule,
499                tol=tol, max_iter=max_iter,
500                n_jobs=n_jobs, seed=seed,
501            )
502            lambda_cv = cv_result.best_lambda
503        else:
504            lambda_cv = lam if lam is not None else 0.0
505
506    elif method == "fe":
507        r_cv = 0
508
509    elif method == "cfe":
510        r_cv = r if isinstance(r, int) else r[0]
511
512    # Point estimation
513    if method in ("fe", "ife"):
514        est = estimate_ife(
515            Y_mat, Y0, X_mat, II_mat, W_mat, beta0,
516            r=r_cv, force=force_int, tol=tol, max_iter=max_iter,
517        )
518    elif method == "mc":
519        est = estimate_mc(
520            Y_mat, Y0, X_mat, II_mat, W_mat, beta0,
521            lam=lambda_cv, force=force_int, tol=tol, max_iter=max_iter,
522        )
523    elif method == "cfe":
524        est = estimate_cfe(
525            Y_mat, Y0, X_mat, II_mat, W_mat, beta0,
526            r=r_cv, force=force_int, tol=tol, max_iter=max_iter,
527        )
528    elif method == "both":
529        # Run both IFE and MC, return IFE results with MC comparison
530        if isinstance(r, tuple) and CV:
531            cv_result = cv_ife(
532                Y_mat, Y0, X_mat, I_mat, II_mat, D_mat, W_mat, beta0,
533                force=force_int, r_range=r, k=k, cv_prop=cv_prop,
534                cv_nobs=cv_nobs, cv_treat=cv_treat, cv_donut=cv_donut,
535                criterion=criterion, cv_rule=cv_rule,
536                tol=tol, max_iter=max_iter,
537                n_jobs=n_jobs, seed=seed,
538            )
539            r_cv = cv_result.best_r
540        else:
541            r_cv = r if isinstance(r, int) else r[0]
542        est = estimate_ife(
543            Y_mat, Y0, X_mat, II_mat, W_mat, beta0,
544            r=r_cv, force=force_int, tol=tol, max_iter=max_iter,
545        )
546    else:
547        raise ValueError(f"Unknown method: {method}")
548
549    # Compute effects
550    eff = Y_mat - est.fit
551    Y_ct = est.fit
552
553    # Additive-FE baseline residual variance (Liu et al. 2024 sigma2.fect).
554    # For method="fe" the main estimator IS the additive-FE pass, so reuse
555    # est.sigma2. For ife/mc/cfe/both, run an extra r=0 IFE pass on the
556    # same panel with the user's requested FE structure.
557    if method == "fe":
558        sigma2_fect_value = float(est.sigma2)
559    else:
560        est_fect = estimate_ife(
561            Y_mat, Y0, X_mat, II_mat, W_mat, beta0,
562            r=0, force=force_int, tol=tol, max_iter=max_iter,
563        )
564        sigma2_fect_value = float(est_fect.sigma2)
565
566    # Denormalize
567    if normalize and norm_factor != 1.0:
568        eff = eff * norm_factor
569        Y_ct = Y_ct * norm_factor
570        Y_mat = Y_mat * norm_factor
571        if est.beta is not None:
572            est = est._replace(beta=est.beta * norm_factor)
573        sigma2_fect_value *= norm_factor ** 2
574
575    # ATT computation
576    T_on = to_device(panel.T_on)
577    att_avg, att_on, time_on, count_on, att_avg_unit = _compute_effects(
578        to_numpy(eff), to_numpy(D_mat), to_numpy(panel.T_on), to_numpy(I_mat),
579    )
580
581    # Build result
582    result = FectResult(
583        method=method,
584        r_cv=r_cv,
585        lambda_cv=lambda_cv,
586        att_avg=att_avg,
587        att_avg_unit=att_avg_unit,
588        att_on=att_on,
589        time_on=time_on,
590        count_on=count_on,
591        beta=to_numpy(est.beta) if est.beta is not None else None,
592        covariate_names=panel.covariate_names,
593        mu=est.mu,
594        alpha=to_numpy(est.alpha) if est.alpha is not None else None,
595        xi=to_numpy(est.xi) if est.xi is not None else None,
596        factors=to_numpy(est.factors) if est.factors is not None else None,
597        loadings=to_numpy(est.loadings) if est.loadings is not None else None,
598        Y_ct=to_numpy(Y_ct),
599        eff=to_numpy(eff),
600        residuals=to_numpy(est.residuals),
601        sigma2=est.sigma2,
602        sigma2_fect=sigma2_fect_value,
603        IC=est.IC,
604        PC=est.PC,
605        niter=est.niter,
606        converged=est.converged,
607        cv_result=cv_result,
608        panel=panel,
609        fit_options={
610            "force": force,
611            "force_int": force_int,
612            "tol": tol,
613            "max_iter": max_iter,
614            "normalize": normalize,
615            "norm_factor": norm_factor,
616            "vartype": vartype,
617            "nboots": nboots,
618            "n_jobs": n_jobs,
619            "min_T0": min_T0,
620            "min_T0_strict": min_T0_strict,
621            "cluster": cluster,
622        },
623        seed=seed,
624    )
625
626    # Inference
627    if se:
628        result.inference = _run_inference(
629            result, panel, Y_mat, X_mat, W_mat, beta0, Y0,
630            method=method, r_cv=r_cv, lambda_cv=lambda_cv,
631            force_int=force_int, tol=tol, max_iter=max_iter,
632            vartype=vartype, nboots=nboots, alpha=alpha,
633            n_jobs=n_jobs, seed=seed, normalize=normalize,
634            norm_factor=norm_factor,
635        )
636
637    # Run requested diagnostics at fit time. requested_diag is None when
638    # diagnostics="none". Validation already enforced se=True and
639    # required-config presence.
640    if requested_diag is not None:
641        opts = dict(diagnostics_options or {})
642        if "loo" in requested_diag:
643            opts["loo"] = True
644        else:
645            opts.setdefault("loo", False)
646        result.diagnostics = _run_diagnostics(
647            result, _requested=requested_diag,
648            _strict=isinstance(diagnostics, list), **opts,
649        )
650
651    return result

Estimate counterfactual treatment effects for panel data.

This is the main Python entry point for the counterfactual estimator workflow. Where the paper and the historical R package differ, pyfector defaults to the paper's statistical definition and exposes R-package-style behavior through explicit options.

Missing outcome policy

pyfector distinguishes raw missing outcomes from counterfactual missingness caused by treatment. Observed untreated cells (D == 0 and non-missing Y) fit the response surface. Observed treated cells (D == 1 and non-missing Y) contribute to ATT as Y - Y_ct. If a treated outcome is missing in the input data, the model can still produce a counterfactual Y_ct for that cell, but the cell is not counted in att_avg or att_on because the treated potential outcome was not observed.

By default, min_T0 is enforced only for treated and reversal units. Sparse controls are retained if they have at least one observed outcome, because they may still inform the low-rank response surface. Set min_T0_strict=True to require controls to satisfy min_T0 too, matching the more conservative R fect sparse-panel behavior.

Parameters
  • data (polars.DataFrame, pandas.DataFrame): Long-format panel data.
  • Y, D (str): Column names for outcome and binary treatment indicator.
  • index ((str, str)): Column names for (unit_id, time_period).
  • X (list of str, optional): Time-varying covariates.
  • W (str, optional): Observation weight column.
  • group (str, optional): Reserved for grouped estimation. Currently raises NotImplementedError when supplied.
  • method ({"fe", "ife", "mc", "cfe", "both"}): Estimation method.
  • force ({"none", "unit", "time", "two-way"}): Fixed effects specification.
  • r (int or (int, int)): Number of factors. If tuple, CV selects from range.
  • lam (float, optional): Nuclear norm penalty for MC. If None with CV=True, auto-selected.
  • nlambda (int): Number of automatically generated lambda candidates for MC CV.
  • lambda_candidates (array-like, optional): Explicit non-negative lambda candidates for MC CV. When supplied, nlambda is ignored.
  • CV (bool): If True, cross-validate over r for IFE when r is a tuple, or over lam for MC when lam is None.
  • k (int): Number of CV folds.
  • cv_prop (float): Fraction of eligible observed control cells masked per CV fold.
  • cv_nobs (int): Number of consecutive within-unit observations to mask as a block.
  • cv_treat (bool): If True, restrict CV masks to pre-treatment cells of ever-treated units. If False, use all observed control cells.
  • cv_donut (int): Exclude this many periods around treatment onset from CV evaluation.
  • criterion ({"mspe", "gmspe", "mad"}): Cross-validation loss.
  • cv_rule ({"min", "onepct"}): CV selection rule. "min" chooses the strict minimum-score candidate and is the paper-faithful default. "onepct" chooses the simplest candidate within 1% of the best score (lower r for IFE, higher lam for MC).
  • se (bool): Compute standard errors via bootstrap/jackknife.
  • vartype ({"bootstrap", "jackknife"}): Inference method when se=True.
  • cluster (str, optional): Column of cluster labels, constant within each unit. Bootstrap inference resamples whole clusters in both the main fit and the placebo refit. Requires at least two retained clusters.
  • nboots (int): Number of bootstrap replications, at least 2. Ignored for jackknife.
  • alpha (float): Significance level for confidence intervals and tests.
  • tol (float): EM convergence tolerance for final point estimation.
  • max_iter (int): Maximum EM iterations.
  • min_T0 (int): Minimum untreated/pre-treatment observed periods. By default this is enforced only for treated and treatment-reversal units.
  • min_T0_strict (bool): If True, enforce min_T0 on all units, including controls. This matches R fect's conservative handling of sparse control rows.
  • max_missing (float): Maximum missing-outcome fraction per unit, in [0, 1]. Units with no observed outcomes are always dropped, regardless of this threshold, because they provide neither fitting information nor observed treated effects.
  • normalize (bool): If True, estimate on an outcome standardized by its observed standard deviation, then transform effects back to the original scale.
  • Z, Q (list of str, optional): Reserved CFE interaction arguments. Currently raise NotImplementedError when supplied.
  • device ({"cpu", "gpu"}): Compute device.
  • n_jobs (int, optional): Parallel workers for CV and bootstrap. -1 or None uses all available CPUs.
  • seed (int, optional): Random seed for full reproducibility.
Returns
  • FectResult: Treatment effects, fitted counterfactuals, and requested inference and diagnostics.
@dataclass
class FectResult: 
 74@dataclass
 75class FectResult:
 76    """Container for all fect estimation results."""
 77    # Method info
 78    method: str
 79    """Estimator used for the point fit."""
 80    r_cv: int | None = None
 81    """Selected IFE rank, or the fixed rank supplied by the caller."""
 82    lambda_cv: float | None = None
 83    """Selected or supplied MC penalty."""
 84
 85    # Point estimates
 86    att_avg: float = 0.0
 87    """ATT averaged over observed treated outcome cells."""
 88    att_avg_unit: float = 0.0
 89    """Mean of unit-specific ATTs over units with observed treated outcomes."""
 90
 91    # Dynamic effects
 92    att_on: np.ndarray | None = None
 93    """Mean counterfactual gap at each relative event time, including pre-treatment periods."""
 94    time_on: np.ndarray | None = None
 95    """Relative event-time grid; treatment onset is zero."""
 96    count_on: np.ndarray | None = None
 97    """Number of observed ever-treated units contributing at each event time."""
 98
 99    # Exit effects (treatment reversal)
100    att_off: np.ndarray | None = None
101    """Post-reversal event-study effects when populated."""
102    time_off: np.ndarray | None = None
103    """Post-reversal event-time grid when populated."""
104
105    # Coefficients
106    beta: np.ndarray | None = None
107    """Estimated covariate coefficients."""
108    covariate_names: list[str] = field(default_factory=list)
109    """Column names corresponding to beta."""
110
111    # Fixed effects
112    mu: float = 0.0
113    """Estimated intercept."""
114    alpha: np.ndarray | None = None   # unit FE
115    """Estimated unit fixed effects."""
116    xi: np.ndarray | None = None      # time FE
117    """Estimated time fixed effects."""
118    factors: np.ndarray | None = None
119    """Estimated latent time factors."""
120    loadings: np.ndarray | None = None
121    """Estimated unit loadings on the latent factors."""
122
123    # Counterfactual and effects matrices
124    Y_ct: np.ndarray | None = None    # T×N counterfactual
125    """Fitted untreated counterfactual outcomes, shaped (T, N)."""
126    eff: np.ndarray | None = None     # T×N treatment effects
127    """Outcome minus fitted counterfactual, shaped (T, N)."""
128    residuals: np.ndarray | None = None
129    """Masked residuals on the untreated fitting sample."""
130
131    # Model fit
132    sigma2: float = 0.0
133    """Chosen estimator residual variance."""
134    sigma2_fect: float = 0.0   # additive-FE baseline residual variance
135    """Additive fixed-effect baseline residual variance used for the TOST bound."""
136    IC: float = 0.0
137    """Estimator information criterion."""
138    PC: float = 0.0
139    """Estimator panel criterion."""
140    rmse: float = 0.0
141    """Reserved field; defaults to 0.0."""
142    niter: int = 0
143    """Number of EM iterations."""
144    converged: bool = False
145    """Whether the estimator reached its convergence tolerance."""
146
147    # Inference
148    inference: InferenceResult | None = None
149    """Bootstrap or jackknife inference results; see pyfector.inference.InferenceResult."""
150
151    # Diagnostics (populated when fect(..., diagnostics="full" | list))
152    diagnostics: Diagnostics | None = None
153    """Requested diagnostic tests and their resolved options."""
154
155    # CV
156    cv_result: CVResult | None = None
157    """Cross-validation candidates and scores; see pyfector.cv.CVResult."""
158
159    # Panel metadata
160    panel: PanelData | None = None
161    """Filtered panel matrices and unit metadata; see pyfector.panel.PanelData."""
162    fit_options: dict[str, Any] = field(default_factory=dict)
163    """Fit settings used by diagnostic refits."""
164
165    # Reproducibility
166    seed: int | None = None
167    """Seed used for cross-validation and bootstrap draws."""
168
169    def summary(self) -> str:
170        """Print summary table of results."""
171        lines = []
172        lines.append(f"pyfector estimation results")
173        lines.append(f"{'='*60}")
174        lines.append(f"Method: {self.method}")
175        if self.r_cv is not None:
176            lines.append(f"Number of factors (CV): {self.r_cv}")
177        if self.lambda_cv is not None:
178            lines.append(f"Lambda (CV): {self.lambda_cv:.6f}")
179        lines.append(f"Converged: {self.converged} (iter={self.niter})")
180        lines.append(f"Sigma^2: {self.sigma2:.6f}")
181        lines.append(f"Sigma^2_fect (FE baseline): {self.sigma2_fect:.6f}")
182        lines.append(f"")
183        lines.append(f"ATT (average): {self.att_avg:.6f}")
184        if self.inference is not None:
185            inf = self.inference
186            lines.append(f"  SE:     {inf.att_avg_se:.6f}")
187            lines.append(f"  CI:     [{inf.att_avg_ci[0]:.6f}, {inf.att_avg_ci[1]:.6f}]")
188            lines.append(f"  p-val:  {inf.att_avg_pval:.4f}")
189
190        if self.beta is not None and len(self.beta) > 0:
191            lines.append(f"")
192            lines.append(f"Coefficients:")
193            for i, name in enumerate(self.covariate_names):
194                lines.append(f"  {name}: {self.beta[i]:.6f}")
195
196        if self.att_on is not None and self.time_on is not None:
197            lines.append(f"")
198            lines.append(f"Dynamic effects (ATT by relative time):")
199            lines.append(f"  {'Time':>6s}  {'ATT':>10s}  {'Count':>6s}", )
200            for i, t in enumerate(self.time_on):
201                count = self.count_on[i] if self.count_on is not None else ""
202                att = self.att_on[i]
203                if self.inference is not None:
204                    se = self.inference.att_on_se[i]
205                    lines.append(f"  {t:>6.0f}  {att:>10.4f}  ({se:.4f})  {count}")
206                else:
207                    lines.append(f"  {t:>6.0f}  {att:>10.4f}  {count}")
208
209        lines.append(f"{'='*60}")
210        if self.panel is not None:
211            lines.append(f"N={self.panel.N}, T={self.panel.T}")
212        if self.seed is not None:
213            lines.append(f"Seed: {self.seed}")
214        if self.diagnostics is not None:
215            lines.append("")
216            lines.append(self.diagnostics.summary())
217        return "\n".join(lines)
218
219    def __repr__(self):
220        return self.summary()
221
222    def plot(self, kind="gap", **kwargs):
223        """Plot results. Shortcut for ``pyfector.plot(self, kind, ...)``."""
224        from .plotting import plot as _plot
225        return _plot(self, kind=kind, **kwargs)
226
227    def diagnose(self, **kwargs):
228        """Run diagnostic tests. Shortcut for ``pyfector.run_diagnostics(self, ...)``."""
229        from .diagnostics import run_diagnostics
230        return run_diagnostics(self, **kwargs)

Container for all fect estimation results.

FectResult( method: str, r_cv: int | None = None, lambda_cv: float | None = None, att_avg: float = 0.0, att_avg_unit: float = 0.0, att_on: numpy.ndarray | None = None, time_on: numpy.ndarray | None = None, count_on: numpy.ndarray | None = None, att_off: numpy.ndarray | None = None, time_off: numpy.ndarray | None = None, beta: numpy.ndarray | None = None, covariate_names: list[str] = <factory>, mu: float = 0.0, alpha: numpy.ndarray | None = None, xi: numpy.ndarray | None = None, factors: numpy.ndarray | None = None, loadings: numpy.ndarray | None = None, Y_ct: numpy.ndarray | None = None, eff: numpy.ndarray | None = None, residuals: numpy.ndarray | None = None, sigma2: float = 0.0, sigma2_fect: float = 0.0, IC: float = 0.0, PC: float = 0.0, rmse: float = 0.0, niter: int = 0, converged: bool = False, inference: pyfector.inference.InferenceResult | None = None, diagnostics: Diagnostics | None = None, cv_result: pyfector.cv.CVResult | None = None, panel: pyfector.panel.PanelData | None = None, fit_options: dict[str, typing.Any] = <factory>, seed: int | None = None)
method: str

Estimator used for the point fit.

r_cv: int | None = None

Selected IFE rank, or the fixed rank supplied by the caller.

lambda_cv: float | None = None

Selected or supplied MC penalty.

att_avg: float = 0.0

ATT averaged over observed treated outcome cells.

att_avg_unit: float = 0.0

Mean of unit-specific ATTs over units with observed treated outcomes.

att_on: numpy.ndarray | None = None

Mean counterfactual gap at each relative event time, including pre-treatment periods.

time_on: numpy.ndarray | None = None

Relative event-time grid; treatment onset is zero.

count_on: numpy.ndarray | None = None

Number of observed ever-treated units contributing at each event time.

att_off: numpy.ndarray | None = None

Post-reversal event-study effects when populated.

time_off: numpy.ndarray | None = None

Post-reversal event-time grid when populated.

beta: numpy.ndarray | None = None

Estimated covariate coefficients.

covariate_names: list[str]

Column names corresponding to beta.

mu: float = 0.0

Estimated intercept.

alpha: numpy.ndarray | None = None

Estimated unit fixed effects.

xi: numpy.ndarray | None = None

Estimated time fixed effects.

factors: numpy.ndarray | None = None

Estimated latent time factors.

loadings: numpy.ndarray | None = None

Estimated unit loadings on the latent factors.

Y_ct: numpy.ndarray | None = None

Fitted untreated counterfactual outcomes, shaped (T, N).

eff: numpy.ndarray | None = None

Outcome minus fitted counterfactual, shaped (T, N).

residuals: numpy.ndarray | None = None

Masked residuals on the untreated fitting sample.

sigma2: float = 0.0

Chosen estimator residual variance.

sigma2_fect: float = 0.0

Additive fixed-effect baseline residual variance used for the TOST bound.

IC: float = 0.0

Estimator information criterion.

PC: float = 0.0

Estimator panel criterion.

rmse: float = 0.0

Reserved field; defaults to 0.0.

niter: int = 0

Number of EM iterations.

converged: bool = False

Whether the estimator reached its convergence tolerance.

inference: pyfector.inference.InferenceResult | None = None

Bootstrap or jackknife inference results; see pyfector.inference.InferenceResult.

diagnostics: Diagnostics | None = None

Requested diagnostic tests and their resolved options.

cv_result: pyfector.cv.CVResult | None = None

Cross-validation candidates and scores; see pyfector.cv.CVResult.

panel: pyfector.panel.PanelData | None = None

Filtered panel matrices and unit metadata; see pyfector.panel.PanelData.

fit_options: dict[str, typing.Any]

Fit settings used by diagnostic refits.

seed: int | None = None

Seed used for cross-validation and bootstrap draws.

def summary(self) -> str: 
169    def summary(self) -> str:
170        """Print summary table of results."""
171        lines = []
172        lines.append(f"pyfector estimation results")
173        lines.append(f"{'='*60}")
174        lines.append(f"Method: {self.method}")
175        if self.r_cv is not None:
176            lines.append(f"Number of factors (CV): {self.r_cv}")
177        if self.lambda_cv is not None:
178            lines.append(f"Lambda (CV): {self.lambda_cv:.6f}")
179        lines.append(f"Converged: {self.converged} (iter={self.niter})")
180        lines.append(f"Sigma^2: {self.sigma2:.6f}")
181        lines.append(f"Sigma^2_fect (FE baseline): {self.sigma2_fect:.6f}")
182        lines.append(f"")
183        lines.append(f"ATT (average): {self.att_avg:.6f}")
184        if self.inference is not None:
185            inf = self.inference
186            lines.append(f"  SE:     {inf.att_avg_se:.6f}")
187            lines.append(f"  CI:     [{inf.att_avg_ci[0]:.6f}, {inf.att_avg_ci[1]:.6f}]")
188            lines.append(f"  p-val:  {inf.att_avg_pval:.4f}")
189
190        if self.beta is not None and len(self.beta) > 0:
191            lines.append(f"")
192            lines.append(f"Coefficients:")
193            for i, name in enumerate(self.covariate_names):
194                lines.append(f"  {name}: {self.beta[i]:.6f}")
195
196        if self.att_on is not None and self.time_on is not None:
197            lines.append(f"")
198            lines.append(f"Dynamic effects (ATT by relative time):")
199            lines.append(f"  {'Time':>6s}  {'ATT':>10s}  {'Count':>6s}", )
200            for i, t in enumerate(self.time_on):
201                count = self.count_on[i] if self.count_on is not None else ""
202                att = self.att_on[i]
203                if self.inference is not None:
204                    se = self.inference.att_on_se[i]
205                    lines.append(f"  {t:>6.0f}  {att:>10.4f}  ({se:.4f})  {count}")
206                else:
207                    lines.append(f"  {t:>6.0f}  {att:>10.4f}  {count}")
208
209        lines.append(f"{'='*60}")
210        if self.panel is not None:
211            lines.append(f"N={self.panel.N}, T={self.panel.T}")
212        if self.seed is not None:
213            lines.append(f"Seed: {self.seed}")
214        if self.diagnostics is not None:
215            lines.append("")
216            lines.append(self.diagnostics.summary())
217        return "\n".join(lines)

Print summary table of results.

def plot(self, kind='gap', **kwargs): 
222    def plot(self, kind="gap", **kwargs):
223        """Plot results. Shortcut for ``pyfector.plot(self, kind, ...)``."""
224        from .plotting import plot as _plot
225        return _plot(self, kind=kind, **kwargs)

Plot results. Shortcut for pyfector.plot(self, kind, ...).

def diagnose(self, **kwargs): 
227    def diagnose(self, **kwargs):
228        """Run diagnostic tests. Shortcut for ``pyfector.run_diagnostics(self, ...)``."""
229        from .diagnostics import run_diagnostics
230        return run_diagnostics(self, **kwargs)

Run diagnostic tests. Shortcut for pyfector.run_diagnostics(self, ...).

def set_device(device: Literal['cpu', 'gpu']) -> None: 
35def set_device(device: Literal["cpu", "gpu"]) -> None:
36    """Set the compute device globally."""
37    global _DEVICE
38    if device == "gpu" and not _check_cupy():
39        raise ImportError(
40            "CuPy is required for GPU support. Install it with: "
41            "pip install cupy-cuda12x  (adjust for your CUDA version)"
42        )
43    _DEVICE = device

Set the compute device globally.

def get_device() -> Literal['cpu', 'gpu']: 
46def get_device() -> Literal["cpu", "gpu"]:
47    """Return the current device."""
48    return _DEVICE

Return the current device.

def run_diagnostics( result, f_threshold: float = 0.5, tost_threshold: float = 0.36, placebo_period: tuple[int, int] | None = None, carryover_period: tuple[int, int] | None = None, loo: bool = False, alpha: float = 0.05, *, _requested: list[str] | None = None, _strict: bool | None = None) -> Diagnostics: 
369def run_diagnostics(
370    result,
371    f_threshold: float = 0.5,
372    tost_threshold: float = 0.36,
373    placebo_period: tuple[int, int] | None = None,
374    carryover_period: tuple[int, int] | None = None,
375    loo: bool = False,
376    alpha: float = 0.05,
377    *,
378    _requested: list[str] | None = None,
379    _strict: bool | None = None,
380) -> Diagnostics:
381    """Run diagnostic tests on a FectResult.
382
383    Parameters
384    ----------
385    result : FectResult
386        Must have inference results (``se=True``).
387    f_threshold : float
388        Noncentrality per treated unit for the equivalence F-test:
389        ``ncp = N * f_threshold``. The default coefficient is 0.5; in the
390        paper's notation this is ``kappa**2``. Small p-values support
391        equivalence. ``N`` is the largest ``count_on`` in the tested
392        pre-periods and also determines both F-tests' scale and degrees
393        of freedom.
394    tost_threshold : float
395        Equivalence bound for TOST. The literal ``0.36`` (default or
396        explicit) triggers Liu et al. (2024)'s scale-aware bound,
397        ``0.36 * sqrt(result.sigma2_fect)``. Any other positive float
398        is taken as an absolute outcome-scale bound. ``None`` is
399        invalid.
400    placebo_period : (start, end), optional
401        Relative pre-treatment window for the holdout/refit placebo test.
402        Observed cells with ``start <= time_on <= end`` and ``time_on < 0``
403        are removed from the fitting mask, the model is refit using the
404        selected rank/lambda configuration, and the placebo ATT is
405        computed from those withheld cells. For example, ``(-3, 0)``
406        withholds relative periods ``-3, -2, -1``.
407        Treated and reversal units must retain at least ``max(min_T0, 1)``
408        untreated observations after holdout. Bootstrap draws use the
409        fit's cluster column when supplied.
410    carryover_period : (start, end), optional
411        Relative time window for carryover test.
412    loo : bool
413        If True, run leave-one-out post-period sensitivity.
414    alpha : float
415        Significance cutoff used for ``TostResult.all_pass``.
416
417    Returns
418    -------
419    Diagnostics
420        Available tests and their resolved options. Unavailable configured
421        tests raise for an explicit list and warn in best-effort mode.
422    """
423    from scipy import stats
424
425    diag = Diagnostics()
426
427    requested = set(_requested) if _requested is not None else _VALID_DIAG_NAMES.copy()
428    strict = _requested is not None if _strict is None else _strict
429    expected = requested.copy()
430    unavailable_reasons = {}
431    if _requested is None:
432        if placebo_period is None:
433            expected.discard("placebo")
434        if carryover_period is None:
435            expected.discard("carryover")
436        if not loo:
437            expected.discard("loo")
438
439    def report_unavailable():
440        missing = sorted((expected - set(diag.tests)) | (expected & set(unavailable_reasons)))
441        if missing:
442            descriptions = [
443                f"{name} ({unavailable_reasons[name]})" if name in unavailable_reasons else name
444                for name in missing
445            ]
446            message = f"Requested diagnostics unavailable: {', '.join(descriptions)}"
447            if strict:
448                raise ValueError(message)
449            warnings.warn(message, RuntimeWarning, stacklevel=3)
450
451    if result.inference is None:
452        report_unavailable()
453        return diag
454
455    inf = result.inference
456    time_on = result.time_on
457    att_on = result.att_on
458
459    if time_on is None or att_on is None:
460        report_unavailable()
461        return diag
462
463    pre_mask = time_on < 0
464
465    pre_idx = np.where(pre_mask)[0]
466    k = len(pre_idx)
467    att_pre = att_on[pre_idx]
468
469    if placebo_period is not None:
470        placebo_period = _validate_placebo_period(placebo_period)
471
472    needs_tost_threshold = bool({"tost", "placebo"} & requested)
473    threshold_abs = float("nan")
474    threshold_source: str | None = None
475    sigma2_fect_used: float | None = None
476    if needs_tost_threshold:
477        _validate_tost_threshold_value(tost_threshold)
478        try:
479            threshold_abs, threshold_source, sigma2_fect_used = _resolve_tost_threshold(
480                result, tost_threshold
481            )
482        except ValueError as exc:
483            if strict or _strict is None:
484                raise
485            for name in requested & {"tost", "placebo"}:
486                unavailable_reasons[name] = str(exc)
487
488    # Pre-trend / equivalence F-tests share the bootstrap covariance pass.
489    if (
490        k > 0 and ("pretrend_f" in requested or "equiv_f" in requested)
491        and inf.att_on_boot is not None
492    ):
493        boot_pre_all = inf.att_on_boot[pre_idx, :]
494        valid_boot = np.all(np.isfinite(boot_pre_all), axis=0)
495        boot_pre = boot_pre_all[:, valid_boot]
496        n_boot = boot_pre.shape[1]
497        if n_boot <= k:
498            boot_pre = None
499
500    else:
501        boot_pre = None
502
503    if boot_pre is not None:
504        counts = getattr(result, "count_on", None)
505        if counts is None:
506            reason = "observed treated-unit counts (count_on) are required"
507            n_treated = 0
508        else:
509            pre_counts = np.asarray(counts)[pre_idx]
510            valid_counts = np.all(np.isfinite(pre_counts)) and np.all(pre_counts >= 0)
511            n_treated = int(np.max(pre_counts)) if valid_counts else 0
512            reason = "insufficient observed treated units for the number of tested pre-periods"
513        if n_treated <= k:
514            for name in requested & {"pretrend_f", "equiv_f"}:
515                unavailable_reasons[name] = reason
516            boot_pre = None
517
518    if boot_pre is not None:
519        S = np.cov(boot_pre)
520        if k == 1:
521            S = np.asarray(S).reshape(1, 1)
522
523        try:
524            cond = np.linalg.cond(S)
525            S_inv = np.linalg.pinv(S) if cond > 1e12 else np.linalg.inv(S)
526            F_raw = float(att_pre @ S_inv @ att_pre)
527            scale = (n_treated - k) / ((n_treated - 1) * k)
528            F_stat = F_raw * scale
529
530            if np.isfinite(F_stat) and F_stat >= 0:
531                if "pretrend_f" in requested:
532                    diag.pretrend_f = PretrendFResult(
533                        f_stat=F_stat,
534                        p_value=float(stats.f.sf(F_stat, k, n_treated - k)),
535                        df1=k,
536                        df2=n_treated - k,
537                    )
538                if "equiv_f" in requested:
539                    ncp = n_treated * f_threshold
540                    diag.equiv_f = EquivFResult(
541                        p_value=float(stats.ncf.cdf(F_stat, k, n_treated - k, ncp)),
542                        f_threshold=f_threshold,
543                        n_treated=n_treated,
544                        ncp=ncp,
545                    )
546        except np.linalg.LinAlgError:
547            pass
548
549    # Per-period TOST.
550    if k > 0 and "tost" in requested and "tost" not in unavailable_reasons and inf.att_on_se is not None:
551        se_pre = inf.att_on_se[pre_idx]
552
553        tost_pvals = np.full(k, np.nan)
554        for i in range(k):
555            if se_pre[i] > 0:
556                if inf.att_on_boot is not None:
557                    n_boot_i = int(np.isfinite(inf.att_on_boot[pre_idx[i], :]).sum())
558                else:
559                    n_boot_i = 200
560                df = max(n_boot_i - 1, 1)
561                t_upper = (att_pre[i] - threshold_abs) / se_pre[i]
562                t_lower = (att_pre[i] + threshold_abs) / se_pre[i]
563                p_upper = float(stats.t.cdf(t_upper, df))
564                p_lower = float(1 - stats.t.cdf(t_lower, df))
565                tost_pvals[i] = max(p_upper, p_lower)
566
567        finite = tost_pvals[np.isfinite(tost_pvals)]
568        max_p = float(finite.max()) if finite.size else float("nan")
569        all_pass = bool(finite.size == k and (finite < alpha).all())
570        if finite.size < k:
571            unavailable_reasons["tost"] = "missing usable standard errors for pre-treatment periods"
572        if finite.size:
573            diag.tost = TostResult(
574                pvals=tost_pvals,
575                periods=time_on[pre_idx].copy(),
576                threshold=threshold_abs,
577                threshold_source=threshold_source,
578                sigma2_fect=sigma2_fect_used,
579                max_pval=max_p,
580                all_pass=all_pass,
581            )
582
583    # Holdout/refit placebo test. This follows R fect's placeboTest path:
584    # selected pre-treatment cells are removed from II before fitting.
585    if "placebo" in requested and "placebo" not in unavailable_reasons and placebo_period is not None:
586        try:
587            diag.placebo = _run_placebo_test(
588                result,
589                placebo_period=placebo_period,
590                threshold_abs=threshold_abs,
591            )
592        except (ValueError, RuntimeError) as exc:
593            if strict or _strict is None:
594                raise
595            unavailable_reasons["placebo"] = str(exc)
596
597    # Carryover (estimate only; bootstrap p deferred).
598    if (
599        "carryover" in requested
600        and carryover_period is not None
601        and hasattr(result, "att_off")
602        and result.att_off is not None
603    ):
604        c_start, c_end = carryover_period
605        time_off = getattr(result, "time_off", None)
606        if time_off is not None:
607            off_mask = (time_off >= c_start) & (time_off <= c_end)
608            if np.any(off_mask):
609                diag.carryover = CarryoverResult(
610                    estimate=float(np.mean(result.att_off[off_mask])),
611                    period=(int(c_start), int(c_end)),
612                )
613
614    # Leave-one-out post-period.
615    if "loo" in requested and loo:
616        post_mask = time_on >= 0
617        post_idx = np.where(post_mask)[0]
618        if len(post_idx) > 1:
619            full_att = result.att_avg
620            loo_atts: list[float] = []
621            loo_periods: list[float] = []
622            for drop_i in post_idx:
623                remaining = np.delete(post_idx, np.where(post_idx == drop_i))
624                if len(remaining) > 0:
625                    loo_atts.append(float(np.mean(att_on[remaining])))
626                    loo_periods.append(float(time_on[drop_i]))
627            atts_arr = np.array(loo_atts)
628            diag.loo = LooResult(
629                atts=atts_arr,
630                periods=np.array(loo_periods),
631                max_change=float(np.max(np.abs(atts_arr - full_att))),
632            )
633
634    diag.options = {
635        "tost_threshold": threshold_abs if needs_tost_threshold else None,
636        "tost_threshold_source": threshold_source,
637        "f_threshold": f_threshold,
638        "alpha": alpha,
639        "placebo_period": placebo_period,
640        "carryover_period": carryover_period,
641        "loo": loo,
642        "sigma2_fect": sigma2_fect_used if sigma2_fect_used is not None
643                       else getattr(result, "sigma2_fect", None),
644    }
645    try:
646        from . import __version__ as _ver
647        diag.options["pyfector_version"] = _ver
648    except Exception:
649        pass
650
651    report_unavailable()
652    return diag

Run diagnostic tests on a FectResult.

Parameters
  • result (FectResult): Must have inference results (se=True).
  • f_threshold (float): Noncentrality per treated unit for the equivalence F-test: ncp = N * f_threshold. The default coefficient is 0.5; in the paper's notation this is kappa**2. Small p-values support equivalence. N is the largest count_on in the tested pre-periods and also determines both F-tests' scale and degrees of freedom.
  • tost_threshold (float): Equivalence bound for TOST. The literal 0.36 (default or explicit) triggers Liu et al. (2024)'s scale-aware bound, 0.36 * sqrt(result.sigma2_fect). Any other positive float is taken as an absolute outcome-scale bound. None is invalid.
  • placebo_period ((start, end), optional): Relative pre-treatment window for the holdout/refit placebo test. Observed cells with start <= time_on <= end and time_on < 0 are removed from the fitting mask, the model is refit using the selected rank/lambda configuration, and the placebo ATT is computed from those withheld cells. For example, (-3, 0) withholds relative periods -3, -2, -1. Treated and reversal units must retain at least max(min_T0, 1) untreated observations after holdout. Bootstrap draws use the fit's cluster column when supplied.
  • carryover_period ((start, end), optional): Relative time window for carryover test.
  • loo (bool): If True, run leave-one-out post-period sensitivity.
  • alpha (float): Significance cutoff used for TostResult.all_pass.
Returns
  • Diagnostics: Available tests and their resolved options. Unavailable configured tests raise for an explicit list and warn in best-effort mode.
@dataclass
class Diagnostics: 
172@dataclass
173class Diagnostics:
174    """Slim-safe registry of diagnostic test results.
175
176    ``tests`` is the future-proof extension point. Built-in tests expose
177    convenience properties for stable user code, but the container itself
178    does not need a new field every time pyfector adds a diagnostic.
179    """
180    options: dict = field(default_factory=dict)
181    tests: dict[str, Any] = field(default_factory=dict)
182
183    @property
184    def available(self) -> tuple[str, ...]:
185        """Names of populated diagnostics, with built-ins first."""
186        known = [name for name in _DIAGNOSTIC_SUMMARY_ORDER if name in self.tests]
187        extra = sorted(name for name in self.tests if name not in set(known))
188        return tuple(known + extra)
189
190    def get(self, name: str, default: Any = None) -> Any:
191        """Return a diagnostic by name, or ``default`` if absent."""
192        return self.tests.get(name, default)
193
194    def set_test(self, name: str, value: Any | None) -> None:
195        """Set or remove a diagnostic result by name."""
196        if not name or not isinstance(name, str):
197            raise ValueError("diagnostic name must be a non-empty string")
198        if value is None:
199            self.tests.pop(name, None)
200        else:
201            self.tests[name] = value
202
203    def set(self, name: str, value: Any | None) -> None:
204        """Compatibility shortcut for :meth:`set_test`."""
205        self.set_test(name, value)
206
207    def __contains__(self, name: str) -> bool:
208        return name in self.tests
209
210    def __getitem__(self, name: str) -> Any:
211        return self.tests[name]
212
213    def _typed(self, name: str, typ: type) -> Any | None:
214        value = self.tests.get(name)
215        return value if isinstance(value, typ) else None
216
217    @property
218    def tost(self) -> TostResult | None:
219        return self._typed("tost", TostResult)
220
221    @tost.setter
222    def tost(self, value: TostResult | None) -> None:
223        self.set_test("tost", value)
224
225    @property
226    def pretrend_f(self) -> PretrendFResult | None:
227        return self._typed("pretrend_f", PretrendFResult)
228
229    @pretrend_f.setter
230    def pretrend_f(self, value: PretrendFResult | None) -> None:
231        self.set_test("pretrend_f", value)
232
233    @property
234    def equiv_f(self) -> EquivFResult | None:
235        return self._typed("equiv_f", EquivFResult)
236
237    @equiv_f.setter
238    def equiv_f(self, value: EquivFResult | None) -> None:
239        self.set_test("equiv_f", value)
240
241    @property
242    def placebo(self) -> PlaceboResult | None:
243        return self._typed("placebo", PlaceboResult)
244
245    @placebo.setter
246    def placebo(self, value: PlaceboResult | None) -> None:
247        self.set_test("placebo", value)
248
249    @property
250    def carryover(self) -> CarryoverResult | None:
251        return self._typed("carryover", CarryoverResult)
252
253    @carryover.setter
254    def carryover(self, value: CarryoverResult | None) -> None:
255        self.set_test("carryover", value)
256
257    @property
258    def loo(self) -> LooResult | None:
259        return self._typed("loo", LooResult)
260
261    @loo.setter
262    def loo(self, value: LooResult | None) -> None:
263        self.set_test("loo", value)
264
265    def summary(self) -> str:
266        lines = ["Diagnostic Tests", "=" * 50]
267        for name in self.available:
268            sub = self.tests[name]
269            if sub is not None:
270                if hasattr(sub, "summary"):
271                    lines.append(sub.summary())
272                else:
273                    lines.append(f"{name}: {sub!r}")
274        return "\n".join(lines)

Slim-safe registry of diagnostic test results.

tests is the future-proof extension point. Built-in tests expose convenience properties for stable user code, but the container itself does not need a new field every time pyfector adds a diagnostic.

Diagnostics(options: dict = <factory>, tests: dict[str, typing.Any] = <factory>)
options: dict
tests: dict[str, typing.Any]
available: tuple[str, ...] 
183    @property
184    def available(self) -> tuple[str, ...]:
185        """Names of populated diagnostics, with built-ins first."""
186        known = [name for name in _DIAGNOSTIC_SUMMARY_ORDER if name in self.tests]
187        extra = sorted(name for name in self.tests if name not in set(known))
188        return tuple(known + extra)

Names of populated diagnostics, with built-ins first.

def get(self, name: str, default: Any = None) -> Any: 
190    def get(self, name: str, default: Any = None) -> Any:
191        """Return a diagnostic by name, or ``default`` if absent."""
192        return self.tests.get(name, default)

Return a diagnostic by name, or default if absent.

def set_test(self, name: str, value: typing.Any | None) -> None: 
194    def set_test(self, name: str, value: Any | None) -> None:
195        """Set or remove a diagnostic result by name."""
196        if not name or not isinstance(name, str):
197            raise ValueError("diagnostic name must be a non-empty string")
198        if value is None:
199            self.tests.pop(name, None)
200        else:
201            self.tests[name] = value

Set or remove a diagnostic result by name.

def set(self, name: str, value: typing.Any | None) -> None: 
203    def set(self, name: str, value: Any | None) -> None:
204        """Compatibility shortcut for :meth:`set_test`."""
205        self.set_test(name, value)

Compatibility shortcut for set_test().

tost: TostResult | None 
217    @property
218    def tost(self) -> TostResult | None:
219        return self._typed("tost", TostResult)
pretrend_f: PretrendFResult | None 
225    @property
226    def pretrend_f(self) -> PretrendFResult | None:
227        return self._typed("pretrend_f", PretrendFResult)
equiv_f: EquivFResult | None 
233    @property
234    def equiv_f(self) -> EquivFResult | None:
235        return self._typed("equiv_f", EquivFResult)
placebo: PlaceboResult | None 
241    @property
242    def placebo(self) -> PlaceboResult | None:
243        return self._typed("placebo", PlaceboResult)
carryover: CarryoverResult | None 
249    @property
250    def carryover(self) -> CarryoverResult | None:
251        return self._typed("carryover", CarryoverResult)
loo: LooResult | None 
257    @property
258    def loo(self) -> LooResult | None:
259        return self._typed("loo", LooResult)
def summary(self) -> str: 
265    def summary(self) -> str:
266        lines = ["Diagnostic Tests", "=" * 50]
267        for name in self.available:
268            sub = self.tests[name]
269            if sub is not None:
270                if hasattr(sub, "summary"):
271                    lines.append(sub.summary())
272                else:
273                    lines.append(f"{name}: {sub!r}")
274        return "\n".join(lines)
DiagnosticResult = <class 'Diagnostics'>
@dataclass(frozen=True)
class TostResult: 
54@dataclass(frozen=True)
55class TostResult:
56    """Per-period TOST equivalence test on pre-treatment coefficients."""
57    pvals: np.ndarray
58    periods: np.ndarray
59    threshold: float
60    threshold_source: str
61    sigma2_fect: float | None
62    max_pval: float
63    all_pass: bool
64
65    def summary(self) -> str:
66        lines = [
67            f"TOST per pre-period (threshold={self.threshold:.4f}, "
68            f"source={self.threshold_source}):"
69        ]
70        for i, t in enumerate(self.periods):
71            p = float(self.pvals[i])
72            status = "PASS" if p < 0.05 else "fail"
73            lines.append(f"  t={t:+.0f}: p={p:.4f} [{status}]")
74        lines.append(
75            f"  max p={self.max_pval:.4f}  all_pass={self.all_pass}"
76        )
77        return "\n".join(lines)

Per-period TOST equivalence test on pre-treatment coefficients.

TostResult( pvals: numpy.ndarray, periods: numpy.ndarray, threshold: float, threshold_source: str, sigma2_fect: float | None, max_pval: float, all_pass: bool)
pvals: numpy.ndarray
periods: numpy.ndarray
threshold: float
threshold_source: str
sigma2_fect: float | None
max_pval: float
all_pass: bool
def summary(self) -> str: 
65    def summary(self) -> str:
66        lines = [
67            f"TOST per pre-period (threshold={self.threshold:.4f}, "
68            f"source={self.threshold_source}):"
69        ]
70        for i, t in enumerate(self.periods):
71            p = float(self.pvals[i])
72            status = "PASS" if p < 0.05 else "fail"
73            lines.append(f"  t={t:+.0f}: p={p:.4f} [{status}]")
74        lines.append(
75            f"  max p={self.max_pval:.4f}  all_pass={self.all_pass}"
76        )
77        return "\n".join(lines)
@dataclass(frozen=True)
class PretrendFResult: 
80@dataclass(frozen=True)
81class PretrendFResult:
82    """Joint F-test for all pre-treatment ATTs equal to zero."""
83    f_stat: float
84    p_value: float
85    df1: int
86    df2: int
87
88    def summary(self) -> str:
89        return (
90            f"Pre-trend F-test:\n"
91            f"  F({self.df1},{self.df2}) = {self.f_stat:.4f}, "
92            f"p = {self.p_value:.4f}"
93        )

Joint F-test for all pre-treatment ATTs equal to zero.

PretrendFResult(f_stat: float, p_value: float, df1: int, df2: int)
f_stat: float
p_value: float
df1: int
df2: int
def summary(self) -> str: 
88    def summary(self) -> str:
89        return (
90            f"Pre-trend F-test:\n"
91            f"  F({self.df1},{self.df2}) = {self.f_stat:.4f}, "
92            f"p = {self.p_value:.4f}"
93        )
@dataclass(frozen=True)
class EquivFResult: 
 96@dataclass(frozen=True)
 97class EquivFResult:
 98    """Equivalence F-test (non-central F) for pre-trend bound."""
 99    p_value: float
100    f_threshold: float
101    n_treated: int | None = None
102    ncp: float | None = None
103
104    def summary(self) -> str:
105        return (
106            f"Equivalence F-test:\n"
107            f"  p = {self.p_value:.4f} (f_threshold={self.f_threshold})"
108        )

Equivalence F-test (non-central F) for pre-trend bound.

EquivFResult( p_value: float, f_threshold: float, n_treated: int | None = None, ncp: float | None = None)
p_value: float
f_threshold: float
n_treated: int | None = None
ncp: float | None = None
def summary(self) -> str: 
104    def summary(self) -> str:
105        return (
106            f"Equivalence F-test:\n"
107            f"  p = {self.p_value:.4f} (f_threshold={self.f_threshold})"
108        )
@dataclass(frozen=True)
class PlaceboResult: 
111@dataclass(frozen=True)
112class PlaceboResult:
113    """Holdout/refit placebo ATT with bootstrap and TOST p-values.
114
115    ``n_obs`` counts withheld cells from eligible units. ``n_dropped_units``
116    counts units excluded after holdout; ``n_boot`` and ``n_boot_failed``
117    count valid and failed placebo bootstrap replications.
118    """
119    estimate: float
120    se: float
121    p_value: float
122    equiv_p_value: float
123    period: tuple[int, int]
124    n_obs: int
125    n_boot: int
126    n_dropped_units: int = 0
127    n_boot_failed: int = 0
128
129    def summary(self) -> str:
130        return (
131            f"Placebo test (holdout/refit, window {self.period}):\n"
132            f"  ATT={self.estimate:.4f} (SE={self.se:.4f}) "
133            f"p={self.p_value:.4f} equiv_p={self.equiv_p_value:.4f}"
134        )

Holdout/refit placebo ATT with bootstrap and TOST p-values.

n_obs counts withheld cells from eligible units. n_dropped_units counts units excluded after holdout; n_boot and n_boot_failed count valid and failed placebo bootstrap replications.

PlaceboResult( estimate: float, se: float, p_value: float, equiv_p_value: float, period: tuple[int, int], n_obs: int, n_boot: int, n_dropped_units: int = 0, n_boot_failed: int = 0)
estimate: float
se: float
p_value: float
equiv_p_value: float
period: tuple[int, int]
n_obs: int
n_boot: int
n_dropped_units: int = 0
n_boot_failed: int = 0
def summary(self) -> str: 
129    def summary(self) -> str:
130        return (
131            f"Placebo test (holdout/refit, window {self.period}):\n"
132            f"  ATT={self.estimate:.4f} (SE={self.se:.4f}) "
133            f"p={self.p_value:.4f} equiv_p={self.equiv_p_value:.4f}"
134        )
@dataclass(frozen=True)
class CarryoverResult: 
137@dataclass(frozen=True)
138class CarryoverResult:
139    """Carryover-window mean ATT (post-reversal). SE/p deferred."""
140    estimate: float
141    period: tuple[int, int]
142
143    def summary(self) -> str:
144        return (
145            f"Carryover test (window {self.period}):\n"
146            f"  ATT={self.estimate:.4f}"
147        )

Carryover-window mean ATT (post-reversal). SE/p deferred.

CarryoverResult(estimate: float, period: tuple[int, int])
estimate: float
period: tuple[int, int]
def summary(self) -> str: 
143    def summary(self) -> str:
144        return (
145            f"Carryover test (window {self.period}):\n"
146            f"  ATT={self.estimate:.4f}"
147        )
@dataclass(frozen=True)
class LooResult: 
150@dataclass(frozen=True)
151class LooResult:
152    """Post-treatment sensitivity from dropping each post-treatment period.
153
154    This is not a leave-one-pre-period-out no-pretrend test.
155    """
156    atts: np.ndarray
157    periods: np.ndarray
158    max_change: float
159
160    def summary(self) -> str:
161        return (
162            f"Leave-one-out:\n"
163            f"  Max ATT change = {self.max_change:.6f}"
164        )

Post-treatment sensitivity from dropping each post-treatment period.

This is not a leave-one-pre-period-out no-pretrend test.

LooResult(atts: numpy.ndarray, periods: numpy.ndarray, max_change: float)
atts: numpy.ndarray
periods: numpy.ndarray
max_change: float
def summary(self) -> str: 
160    def summary(self) -> str:
161        return (
162            f"Leave-one-out:\n"
163            f"  Max ATT change = {self.max_change:.6f}"
164        )
def plot( result, kind: Literal['gap', 'status', 'factors', 'counterfactual', 'equiv', 'calendar'] = 'gap', units: list | None = None, show_ci: bool = True, title: str | None = None, figsize: tuple[float, float] = (10, 6), ax=None, **kwargs): 
20def plot(
21    result,
22    kind: Literal["gap", "status", "factors", "counterfactual", "equiv", "calendar"] = "gap",
23    units: list | None = None,
24    show_ci: bool = True,
25    title: str | None = None,
26    figsize: tuple[float, float] = (10, 6),
27    ax=None,
28    **kwargs,
29):
30    """Plot fect results.
31
32    Parameters
33    ----------
34    result : FectResult
35        Output from ``pyfector.fect()``.
36    kind : str
37        Plot type.
38    units : list, optional
39        Unit IDs for counterfactual plot.
40    show_ci : bool
41        Show confidence intervals (requires ``se=True`` in estimation).
42    """
43    import matplotlib.pyplot as plt
44
45    if ax is None:
46        fig, ax = plt.subplots(figsize=figsize)
47    else:
48        fig = ax.get_figure()
49
50    if kind == "gap":
51        _plot_gap(result, ax, show_ci, **kwargs)
52    elif kind == "status":
53        _plot_status(result, ax, **kwargs)
54    elif kind == "factors":
55        _plot_factors(result, ax, **kwargs)
56    elif kind == "counterfactual":
57        _plot_counterfactual(result, ax, units, **kwargs)
58    elif kind == "equiv":
59        _plot_equiv(result, ax, **kwargs)
60    elif kind == "calendar":
61        _plot_calendar(result, ax, show_ci, **kwargs)
62    else:
63        raise ValueError(f"Unknown plot kind: {kind}")
64
65    if title:
66        ax.set_title(title)
67
68    fig.tight_layout()
69    return fig

Plot fect results.

Parameters
  • result (FectResult): Output from pyfector.fect.
  • kind (str): Plot type.
  • units (list, optional): Unit IDs for counterfactual plot.
  • show_ci (bool): Show confidence intervals (requires se=True in estimation).