Skip to content

estimators

Inputs to construction/optimization: expected returns, covariance, and factor models. Distinct from portpy.metrics.covariance, which decomposes realized risk for a given weight vector against a given matrix - this module produces the matrix itself.

Expected returns

expected_returns 

expected_returns(
    y: DataFrame,
    method: str = "mean_historical",
    periods_per_year: int = TRADING_DAYS_PER_YEAR,
    geometric: bool = False,
    span: int | None = None,
    benchmark: Series | None = None,
    rf: float = 0.0,
    prior: float | Series | None = None,
) -> pd.Series

Estimate an annualized expected-return vector.

Parameters:

Name Type Description Default
y DataFrame

Per-asset periodic returns - the full multi-asset DataFrame (matches Portfolio.asset_returns()), not the portfolio's own weighted series.

 required
method str

One of: - "mean_historical": annualized historical mean return per asset (arithmetic mean * periods_per_year by default; see geometric). - "ewma": exponentially-weighted mean return, annualized. span sets the decay window (default 60 periods) - recent returns matter more than old ones. - "capm_implied": CAPM-implied return, rf + beta_i * (E[market] - rf), using each asset's historical beta against benchmark. - "james_stein": cross-sectional James-Stein shrinkage of the historical mean toward the grand (equal-weighted) mean - a well-documented fix for plain historical means' large estimation error (Jorion, 1986; Efron & Morris, 1975).

 'mean_historical'
periods_per_year int

Periods per year, used to annualize.

 TRADING_DAYS_PER_YEAR
geometric bool

"mean_historical" only - use CAGR-style geometric annualization per asset instead of mean * periods_per_year.

 False
span int | None

EWMA decay in periods, used by "ewma" (default 60 if omitted).

 None
benchmark Series | None

Market/benchmark return series, required by "capm_implied".

 None
rf float

Annual risk-free rate, used by "capm_implied".

 0.0
prior float | Series | None

Optional shrinkage target for "james_stein" - a single annual rate applied to every asset, or a per-asset Series. Defaults to the equal-weighted grand mean of y's own historical means.

 None

Returns:

Type Description
Series

Annualized expected return per asset, indexed by y.columns, named "expected_return".

Raises:

Type Description
ValueError

If method is unrecognized, or benchmark is missing for "capm_implied".
Source code in src/portpy/models/estimators/expected_returns.py 
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
def expected_returns(
    y: pd.DataFrame,
    method: str = "mean_historical",
    periods_per_year: int = TRADING_DAYS_PER_YEAR,
    geometric: bool = False,
    span: int | None = None,
    benchmark: pd.Series | None = None,
    rf: float = 0.0,
    prior: float | pd.Series | None = None,
) -> pd.Series:
    """
    Estimate an annualized expected-return vector.

    Args:
        y: Per-asset periodic returns - the full multi-asset DataFrame (matches
            `Portfolio.asset_returns()`), not the portfolio's own weighted series.
        method: One of:
            - "mean_historical": annualized historical mean return per asset
              (arithmetic mean * periods_per_year by default; see `geometric`).
            - "ewma": exponentially-weighted mean return, annualized. `span`
              sets the decay window (default 60 periods) - recent returns
              matter more than old ones.
            - "capm_implied": CAPM-implied return, `rf + beta_i * (E[market] - rf)`,
              using each asset's historical beta against `benchmark`.
            - "james_stein": cross-sectional James-Stein shrinkage of the
              historical mean toward the grand (equal-weighted) mean - a
              well-documented fix for plain historical means' large estimation
              error (Jorion, 1986; Efron & Morris, 1975).
        periods_per_year: Periods per year, used to annualize.
        geometric: "mean_historical" only - use CAGR-style geometric
            annualization per asset instead of `mean * periods_per_year`.
        span: EWMA decay in periods, used by "ewma" (default 60 if omitted).
        benchmark: Market/benchmark return series, required by "capm_implied".
        rf: Annual risk-free rate, used by "capm_implied".
        prior: Optional shrinkage target for "james_stein" - a single annual
            rate applied to every asset, or a per-asset Series. Defaults to
            the equal-weighted grand mean of `y`'s own historical means.

    Returns:
        Annualized expected return per asset, indexed by `y.columns`, named "expected_return".

    Raises:
        ValueError: If `method` is unrecognized, or `benchmark` is missing for "capm_implied".
    """
    ensure_min_observations(y, 2, "y")

    if method == "mean_historical":
        if geometric:
            values = {col: float(annualized_return(y[col].dropna(), periods_per_year=periods_per_year, geometric=True)) for col in y.columns}
            return pd.Series(values, name="expected_return")
        return (y.mean() * periods_per_year).rename("expected_return")

    if method == "ewma":
        window = 60 if span is None else span
        ewma_per_period = y.ewm(span=window).mean().iloc[-1]
        return (ewma_per_period * periods_per_year).rename("expected_return")

    if method == "capm_implied":
        if benchmark is None:
            raise ValueError("method='capm_implied' requires benchmark=<market return series>.")
        market_ann = float(benchmark.mean() * periods_per_year)
        betas = {}
        for col in y.columns:
            r, b = align_pair(y[col], benchmark, name_a=col, name_b="benchmark")
            cov = float(np.cov(r, b, ddof=1)[0, 1])
            var_b = float(np.var(b, ddof=1))
            betas[col] = cov / var_b if var_b > 1e-30 else np.nan
        beta_s = pd.Series(betas)
        return (rf + beta_s * (market_ann - rf)).rename("expected_return")

    if method == "james_stein":
        mu_periodic = y.mean()
        cov_periodic = y.cov().to_numpy()
        mu = mu_periodic * periods_per_year
        n_assets = len(mu)
        t_obs = len(y)

        if prior is None:
            target_periodic = pd.Series(float(mu_periodic.mean()), index=mu_periodic.index)
        elif isinstance(prior, (int, float)):
            target_periodic = pd.Series(float(prior) / periods_per_year, index=mu_periodic.index)
        else:
            target_periodic = pd.Series(prior).reindex(mu_periodic.index) / periods_per_year

        # Mahalanobis distance (and therefore phi) must be computed at the SAME periodic
        # scale as t_obs, the observation count - Jorion's (1986) shrinkage factor pairs a
        # sample size with a distance measured at that sample's own periodicity. Computing
        # it on annualized mu/cov instead (as an earlier version of this function did) mixes
        # an annualized distance with a periodic T: since mahalanobis computed on annualized
        # inputs scales up by periods_per_year relative to the periodic value, phi collapses
        # toward 0 for any daily-return dataset regardless of the assets' true dispersion -
        # silently defeating the shrinkage entirely at the default periods_per_year=252.
        diff_periodic = (mu_periodic - target_periodic).to_numpy()
        try:
            cov_inv = np.linalg.pinv(cov_periodic)
            mahalanobis = float(diff_periodic @ cov_inv @ diff_periodic)
        except np.linalg.LinAlgError:
            mahalanobis = float(np.sum(diff_periodic**2))

        phi = (n_assets + 2) / ((n_assets + 2) + t_obs * mahalanobis) if mahalanobis > 1e-12 else 1.0
        phi = float(np.clip(phi, 0.0, 1.0))
        target = target_periodic * periods_per_year
        shrunk = phi * target + (1.0 - phi) * mu
        return shrunk.rename("expected_return")

    raise ValueError(f"method must be one of {_METHODS}, got {method!r}.")

Covariance

covariance 

covariance(
    y: DataFrame,
    method: str = "sample",
    span: int | None = None,
    periods_per_year: int = TRADING_DAYS_PER_YEAR,
    annualize: bool = True,
    shrinkage_intensity: float = 0.2,
) -> pd.DataFrame

Estimate an (by default annualized) asset covariance matrix.

Parameters:

Name Type Description Default
y DataFrame

Per-asset periodic returns - the full multi-asset DataFrame (matches Portfolio.asset_returns()).

 required
method str

One of: - "sample": plain sample covariance (y.cov(), ddof=1). Same formula as metrics.covariance.covariance_matrix, just annualized by default here since optimizers expect it paired with an annualized expected_returns. - "ewma": RiskMetrics-style exponentially-weighted covariance - recent co-movements dominate, controlled by span. - "shrinkage": manual-intensity linear shrinkage of the sample matrix toward a diagonal target that keeps each asset's own sample variance and zeroes every cross-correlation, at a user-chosen shrinkage_intensity. This is a different target from the "shrink every variance toward the cross-sectional average" convention some other libraries use for their own "shrinkage" method (e.g. PyPortfolioOpt's shrunk_covariance, which shrinks toward a scaled-identity matrix - constant variance across assets - instead of a diagonal-of-sample matrix) - the two do not, in general, agree at the same intensity, and neither is "more correct" than the other. - "ledoit_wolf": Ledoit-Wolf (2004) automatic-intensity shrinkage toward a scaled-identity target (sklearn.covariance.LedoitWolf) - the shrinkage intensity is estimated from the data, not chosen by hand. - "robust": Minimum Covariance Determinant (sklearn.covariance.MinCovDet) - a high-breakdown-point estimator resistant to outlier periods (e.g. crashes). MinCovDet fits the raw MCD subset covariance and then applies a statistical consistency correction (reweighting using the subset's own outlier threshold) on top of it, so this is not numerically identical to a "raw" MCD covariance computed some other way (e.g. sklearn.covariance.fast_mcd called directly, which some other libraries wrap instead - see the cross-validation notebook for the resulting gap on real data). Needs enough observations relative to the asset count for the MCD subset search to be well-posed; with too few rows relative to columns, MinCovDet raises rather than silently degrading.

 'sample'
span int | None

EWMA decay in periods, used by "ewma" (default 60 if omitted).

 None
periods_per_year int

Periods per year, used when annualize is True.

 TRADING_DAYS_PER_YEAR
annualize bool

If True (default), scale the matrix by periods_per_year (assumes i.i.d., serially uncorrelated returns - same caveat as metrics.covariance.covariance_matrix).

 True
shrinkage_intensity float

Blend weight on the diagonal target for "shrinkage", in [0, 1] (0 = pure sample covariance, 1 = pure diagonal).

 0.2

Returns:

Type Description
DataFrame

The (optionally annualized) covariance matrix, indexed and columned by asset.

Raises:

Type Description
ValueError

If method is unrecognized or shrinkage_intensity is outside [0, 1].
Source code in src/portpy/models/estimators/covariance.py 
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
def covariance(
    y: pd.DataFrame,
    method: str = "sample",
    span: int | None = None,
    periods_per_year: int = TRADING_DAYS_PER_YEAR,
    annualize: bool = True,
    shrinkage_intensity: float = 0.2,
) -> pd.DataFrame:
    """
    Estimate an (by default annualized) asset covariance matrix.

    Args:
        y: Per-asset periodic returns - the full multi-asset DataFrame (matches
            `Portfolio.asset_returns()`).
        method: One of:
            - "sample": plain sample covariance (`y.cov()`, ddof=1). Same
              formula as `metrics.covariance.covariance_matrix`, just annualized
              by default here since optimizers expect it paired with an
              annualized `expected_returns`.
            - "ewma": RiskMetrics-style exponentially-weighted covariance -
              recent co-movements dominate, controlled by `span`.
            - "shrinkage": manual-intensity linear shrinkage of the sample
              matrix toward a diagonal target that keeps each asset's own
              sample variance and zeroes every cross-correlation, at a
              user-chosen `shrinkage_intensity`. This is a different target
              from the "shrink every variance toward the cross-sectional
              average" convention some other libraries use for their own
              "shrinkage" method (e.g. PyPortfolioOpt's `shrunk_covariance`,
              which shrinks toward a *scaled-identity* matrix - constant
              variance across assets - instead of a diagonal-of-sample
              matrix) - the two do not, in general, agree at the same
              intensity, and neither is "more correct" than the other.
            - "ledoit_wolf": Ledoit-Wolf (2004) automatic-intensity shrinkage
              toward a scaled-identity target (`sklearn.covariance.LedoitWolf`) -
              the shrinkage intensity is estimated from the data, not chosen by hand.
            - "robust": Minimum Covariance Determinant (`sklearn.covariance.MinCovDet`) -
              a high-breakdown-point estimator resistant to outlier periods (e.g. crashes).
              `MinCovDet` fits the raw MCD subset covariance and then applies a
              statistical consistency correction (reweighting using the
              subset's own outlier threshold) on top of it, so this is not
              numerically identical to a "raw" MCD covariance computed some
              other way (e.g. `sklearn.covariance.fast_mcd` called directly,
              which some other libraries wrap instead - see the cross-validation
              notebook for the resulting gap on real data). Needs enough
              observations relative to the asset count for the MCD subset
              search to be well-posed; with too few rows relative to columns,
              `MinCovDet` raises rather than silently degrading.
        span: EWMA decay in periods, used by "ewma" (default 60 if omitted).
        periods_per_year: Periods per year, used when `annualize` is True.
        annualize: If True (default), scale the matrix by `periods_per_year`
            (assumes i.i.d., serially uncorrelated returns - same caveat as
            `metrics.covariance.covariance_matrix`).
        shrinkage_intensity: Blend weight on the diagonal target for
            "shrinkage", in [0, 1] (0 = pure sample covariance, 1 = pure diagonal).

    Returns:
        The (optionally annualized) covariance matrix, indexed and columned by asset.

    Raises:
        ValueError: If `method` is unrecognized or `shrinkage_intensity` is
            outside [0, 1].
    """
    ensure_min_observations(y, 2, "y")
    scale = periods_per_year if annualize else 1.0

    if method == "sample":
        return y.cov() * scale

    if method == "ewma":
        window = 60 if span is None else span
        lam = 1.0 - 2.0 / (window + 1.0)
        arr = y.to_numpy(dtype=float)
        arr = arr - np.nanmean(arr, axis=0)
        n, k = arr.shape
        weights = lam ** np.arange(n - 1, -1, -1)
        weights = weights / weights.sum()
        weighted = arr * weights[:, None]
        cov = weighted.T @ arr
        return pd.DataFrame(cov * scale, index=y.columns, columns=y.columns)

    if method == "shrinkage":
        if not (0.0 <= shrinkage_intensity <= 1.0):
            raise ValueError(f"shrinkage_intensity must be in [0, 1], got {shrinkage_intensity}.")
        sample = y.cov()
        target = pd.DataFrame(np.diag(np.diag(sample.to_numpy())), index=sample.index, columns=sample.columns)
        blended = shrinkage_intensity * target + (1.0 - shrinkage_intensity) * sample
        return blended * scale

    if method == "ledoit_wolf":
        from sklearn.covariance import LedoitWolf

        fitted = LedoitWolf().fit(y.dropna().to_numpy())
        return pd.DataFrame(fitted.covariance_ * scale, index=y.columns, columns=y.columns)

    if method == "robust":
        from sklearn.covariance import MinCovDet

        fitted = MinCovDet(random_state=0).fit(y.dropna().to_numpy())
        return pd.DataFrame(fitted.covariance_ * scale, index=y.columns, columns=y.columns)

    raise ValueError(f"method must be one of {_METHODS}, got {method!r}.")

Factor models

factor_models

Factor regressions: CAPM, Fama-French, and general linear/rolling factor models.

Replaces the old experimental metrics/regressions.py module. Everything here is a thin, opinionated wrapper around statsmodels OLS - PortPy adds the excess-return bookkeeping (rf, RF-column handling) and sensible defaults; statsmodels does the actual estimation, so every RegressionResultsWrapper returned here carries the full statsmodels API (.params, .rsquared, .pvalues, .summary(), ...).

linear_regression 

linear_regression(
    y: Series, x: DataFrame
) -> RegressionResultsWrapper

OLS regression of y on x, with an intercept - the general building block every other function here specializes.

Parameters:

Name Type Description Default
y Series

Dependent variable (e.g. asset returns).

 required
x DataFrame

Independent variable(s) (e.g. one or more factor return series).

 required

Returns:

Type Description
RegressionResultsWrapper

A fitted statsmodels RegressionResultsWrapper (.params,
RegressionResultsWrapper

.rsquared, .pvalues, .summary(), ...). The intercept is labeled "const".
Source code in src/portpy/models/estimators/factor_models.py 
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
def linear_regression(y: pd.Series, x: pd.DataFrame) -> RegressionResultsWrapper:
    """
    OLS regression of `y` on `x`, with an intercept - the general building
    block every other function here specializes.

    Args:
        y: Dependent variable (e.g. asset returns).
        x: Independent variable(s) (e.g. one or more factor return series).

    Returns:
        A fitted `statsmodels` `RegressionResultsWrapper` (`.params`,
        `.rsquared`, `.pvalues`, `.summary()`, ...). The intercept is labeled "const".
    """
    if isinstance(x, pd.Series):
        x = x.to_frame(x.name or "x")
    aligned = pd.concat([y.rename("y"), x], axis=1, join="inner").dropna()
    ensure_min_observations(aligned, x.shape[1] + 2, "y/x (aligned)")
    design = sm.add_constant(aligned.drop(columns="y"))
    return sm.OLS(aligned["y"], design).fit()

rolling_regression 

rolling_regression(
    y: Series, x: DataFrame, window: int
) -> pd.DataFrame

Rolling-window OLS of y on x - a rolling alpha and rolling factor loading(s), via statsmodels.regression.rolling.RollingOLS.

Parameters:

Name Type Description Default
y Series

Dependent variable (e.g. asset returns).

 required
x DataFrame

Independent variable(s), same length/index as y.

 required
window int

Rolling window size, in periods.

 required

Returns:

Type Description
DataFrame

A DataFrame indexed like y, with one column "alpha" (the rolling
DataFrame

intercept) plus one column per column of x (its rolling loading).
DataFrame

The first window - 1 rows are NaN.
Source code in src/portpy/models/estimators/factor_models.py 
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
def rolling_regression(y: pd.Series, x: pd.DataFrame, window: int) -> pd.DataFrame:
    """
    Rolling-window OLS of `y` on `x` - a rolling alpha and rolling factor
    loading(s), via `statsmodels.regression.rolling.RollingOLS`.

    Args:
        y: Dependent variable (e.g. asset returns).
        x: Independent variable(s), same length/index as `y`.
        window: Rolling window size, in periods.

    Returns:
        A DataFrame indexed like `y`, with one column "alpha" (the rolling
        intercept) plus one column per column of `x` (its rolling loading).
        The first `window - 1` rows are NaN.
    """
    if isinstance(x, pd.Series):
        x = x.to_frame(x.name or "x")
    aligned = pd.concat([y.rename("y"), x], axis=1, join="inner").dropna()
    ensure_min_observations(aligned, window, "y/x (aligned)")
    design = sm.add_constant(aligned.drop(columns="y"))
    fitted = RollingOLS(aligned["y"], design, window=window, min_nobs=window, missing="drop").fit()
    return fitted.params.rename(columns={"const": "alpha"})

capm 

capm(
    y: Series,
    benchmark: Series,
    rf: float = DEFAULT_RISK_FREE_RATE,
    periods_per_year: int = TRADING_DAYS_PER_YEAR,
) -> RegressionResultsWrapper

Single-factor CAPM regression: (y - rf) = alpha + beta * (benchmark - rf) + e.

Parameters:

Name Type Description Default
y Series

Asset/portfolio return series.

 required
benchmark Series

Market/benchmark return series.

 required
rf float

Annual risk-free rate, de-annualized internally and subtracted from both y and benchmark before regressing.

 DEFAULT_RISK_FREE_RATE
periods_per_year int

Periods per year, used to de-annualize rf.

 TRADING_DAYS_PER_YEAR

Returns:

Type Description
RegressionResultsWrapper

A fitted RegressionResultsWrapper; .params["const"] is the
RegressionResultsWrapper

(periodic-scale) alpha, .params["market"] is beta.
Source code in src/portpy/models/estimators/factor_models.py 
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
def capm(
    y: pd.Series,
    benchmark: pd.Series,
    rf: float = DEFAULT_RISK_FREE_RATE,
    periods_per_year: int = TRADING_DAYS_PER_YEAR,
) -> RegressionResultsWrapper:
    """
    Single-factor CAPM regression: `(y - rf) = alpha + beta * (benchmark - rf) + e`.

    Args:
        y: Asset/portfolio return series.
        benchmark: Market/benchmark return series.
        rf: Annual risk-free rate, de-annualized internally and subtracted
            from both `y` and `benchmark` before regressing.
        periods_per_year: Periods per year, used to de-annualize `rf`.

    Returns:
        A fitted `RegressionResultsWrapper`; `.params["const"]` is the
        (periodic-scale) alpha, `.params["market"]` is beta.
    """
    r, b = align_pair(y, benchmark, name_a="y", name_b="market")
    rf_period = periodic_rate_from_annual(rf, periods_per_year)
    excess_y = r - rf_period
    excess_b = (b - rf_period).rename("market")
    design = sm.add_constant(excess_b)
    return sm.OLS(excess_y, design).fit()

fama_french 

fama_french(
    y: Series, factors: DataFrame, version: int = 3
) -> RegressionResultsWrapper

Multi-factor Fama-French style regression of y on factors.

Parameters:

Name Type Description Default
y Series

Asset/portfolio return series (raw, not pre-adjusted for rf - see the "RF" handling below).

 required
factors DataFrame

Factor return series as columns. If a column named "RF" is present, it's treated as the risk-free rate: y is converted to excess-of-rf returns and "RF" is dropped from the regressors (matching how Ken French's data library distributes factors).

 required
version int

3 or 5 validates that factors has the canonical named Fama-French 3-factor (Mkt-RF, SMB, HML) or 5-factor (+RMW, CMA) columns before regressing. Any other value skips validation and regresses on whatever columns factors provides - useful for extended factor sets (momentum, quality, ...) that don't have one universally agreed column count.

 3

Returns:

Type Description
RegressionResultsWrapper

A fitted RegressionResultsWrapper over factors' columns (minus "RF").

Raises:

Type Description
ValueError

If version is 3 or 5 and factors is missing the corresponding canonical columns.
Source code in src/portpy/models/estimators/factor_models.py 
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
def fama_french(y: pd.Series, factors: pd.DataFrame, version: int = 3) -> RegressionResultsWrapper:
    """
    Multi-factor Fama-French style regression of `y` on `factors`.

    Args:
        y: Asset/portfolio return series (raw, not pre-adjusted for rf - see
            the "RF" handling below).
        factors: Factor return series as columns. If a column named "RF" is
            present, it's treated as the risk-free rate: `y` is converted to
            excess-of-rf returns and "RF" is dropped from the regressors
            (matching how Ken French's data library distributes factors).
        version: 3 or 5 validates that `factors` has the canonical named
            Fama-French 3-factor (`Mkt-RF`, `SMB`, `HML`) or 5-factor
            (+`RMW`, `CMA`) columns before regressing. Any other value skips
            validation and regresses on whatever columns `factors` provides -
            useful for extended factor sets (momentum, quality, ...) that
            don't have one universally agreed column count.

    Returns:
        A fitted `RegressionResultsWrapper` over `factors`' columns (minus "RF").

    Raises:
        ValueError: If `version` is 3 or 5 and `factors` is missing the
            corresponding canonical columns.
    """
    if version in _CANONICAL_FF_FACTORS:
        missing = [c for c in _CANONICAL_FF_FACTORS[version] if c not in factors.columns]
        if missing:
            raise ValueError(
                f"factors is missing canonical Fama-French {version}-factor column(s): {missing}. "
                "Rename your columns to match, or pass a version not in {3, 5} to skip this check "
                "and regress on whatever columns `factors` has."
            )

    factors = factors.copy()
    y_adj = y.rename(y.name or "y")
    if "RF" in factors.columns:
        aligned_rf = pd.concat([y_adj.rename("y"), factors["RF"]], axis=1, join="inner").dropna()
        y_adj = (aligned_rf["y"] - aligned_rf["RF"]).rename(y.name or "y")
        factors = factors.drop(columns="RF")

    aligned = pd.concat([y_adj.rename("y"), factors], axis=1, join="inner").dropna()
    ensure_min_observations(aligned, factors.shape[1] + 2, "y/factors (aligned)")
    design = sm.add_constant(aligned.drop(columns="y"))
    return sm.OLS(aligned["y"], design).fit()

factor_attribution 

factor_attribution(
    y: Series, factor_returns: DataFrame
) -> pd.DataFrame

Decompose y's average per-period return into a contribution from each factor (beta_i * mean(factor_i)) plus an unexplained "alpha" residual.

Parameters:

Name Type Description Default
y Series

Asset/portfolio return series.

 required
factor_returns DataFrame

Factor return series as columns.

 required

Returns:

Type Description
DataFrame

A DataFrame indexed by factor name (plus "alpha"), with columns
DataFrame

"beta" (NaN for the alpha row), "contribution" (return units), and
DataFrame

"pct_of_total" (each contribution's share of y's mean return).
Source code in src/portpy/models/estimators/factor_models.py 
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
def factor_attribution(y: pd.Series, factor_returns: pd.DataFrame) -> pd.DataFrame:
    """
    Decompose `y`'s average per-period return into a contribution from each
    factor (`beta_i * mean(factor_i)`) plus an unexplained "alpha" residual.

    Args:
        y: Asset/portfolio return series.
        factor_returns: Factor return series as columns.

    Returns:
        A DataFrame indexed by factor name (plus "alpha"), with columns
        "beta" (NaN for the alpha row), "contribution" (return units), and
        "pct_of_total" (each contribution's share of `y`'s mean return).
    """
    aligned = pd.concat([y.rename("y"), factor_returns], axis=1, join="inner").dropna()
    ensure_min_observations(aligned, factor_returns.shape[1] + 2, "y/factor_returns (aligned)")
    design = sm.add_constant(aligned.drop(columns="y"))
    model = sm.OLS(aligned["y"], design).fit()

    betas = model.params.drop("const")
    factor_means = aligned[betas.index].mean()
    contributions = betas * factor_means
    alpha = float(model.params["const"])
    total = float(contributions.sum()) + alpha

    out = pd.DataFrame({"beta": betas, "contribution": contributions})
    out.loc["alpha"] = {"beta": np.nan, "contribution": alpha}
    out["pct_of_total"] = out["contribution"] / total if abs(total) > 1e-12 else np.nan
    return out