factrix.multi_factor.bhy_hierarchical ¶

bhy_hierarchical(results: list[EvaluationResult], *, metrics: list[str], group: str, family_size: int | Mapping[tuple[Any, ...], int] | None = None, q: float = 0.05, inactive_policy: InactivePolicy = 'count') -> dict[str, HierarchicalBhyResult]

Two-layer hierarchical BHY screen with selective inference, one per metric.

Each group is represented by the Simes combination of its members' p-values; BHY runs across those representatives (outer layer) and again within each group (inner layer). Flat BHY across the whole input loses group-level interpretability and pays full m-correction even when most groups are dead.

Selection adjustment. The inner layer runs at q · R / G — R groups the outer layer selected out of G — per Benjamini-Bogomolov (2014) (the Yekutieli (2008) hierarchical framework with the selective inner level). A member's adjusted p is max(outer_adj[group], min(1, inner_adj · G / R)), so survivor iff adj_p <= q holds at the supplied q. Because R depends on q, the adjusted p-values are defined at that q: re-screening at another level is a new call, not a threshold change on the returned array.

Why the selective level, not the nominal q at both layers. An earlier version used the nominal q for the inner family and took the max against the outer adjusted p. That bounds the group-level error rate but not member-level FDR — a correctly-passed group can still contribute many false member rejections — and it was measured to break on the library's modal input: small groups of correlated members (a family of four momentum variants is size = 4 with high within-group correlation). At G = 30 groups of 4, five live, one true effect each, realised FDR was 1.09× / 1.18× / 1.14× the nominal 0.10 at rho = 0.5 / 0.7 / 0.9. The selective level restores control there (0.03–0.04) at no power cost in that regime (0.619 vs 0.619).

The cost lands elsewhere: on large groups with dense signal (size = 10, half the members non-null, independent) the nominal-q construction detected more (0.75 vs 0.55) while still controlling FDR. That edge was real, but a procedure whose headline guarantee fails on its typical input cannot keep it. For dense, large families flat :func:bhy is the higher-power alternative that also carries a theorem-backed guarantee.

Dependence assumption — the outer layer is not arbitrary-dependence-valid. BHY itself holds under arbitrary dependence, but each group enters the outer layer through :func:simes_p, and Simes' combination is valid only under independence or positive regression dependence (PRDS). Under negative within-group dependence Simes is anti-conservative, so the composite procedure inherits that limitation despite BHY's own guarantee. Negative cross-regional correlation is not exotic. The simulation evidence below covers positive equicorrelation only (rho in {0.5, 0.7, 0.9}); no negative-dependence cell was measured. For a family with material negative dependence, prefer flat :func:bhy, whose guarantee does not route through Simes.

TestHierarchicalFdrControl pins realised FDR across group counts, sizes, sparsity, selected fraction down to R / G = 0.01, and the small-correlated-group scan above.

Parameters:

Name	Type	Description	Default
`results`	`list[EvaluationResult]`	:class:`EvaluationResult` records. Each is assigned to one group via `result.params[group]` (or via `result.forward_periods` if `group == "forward_periods"`). Within a group, each member is one hypothesis of the inner family; the `(factor, forward_periods, *params)` identifier must be unique across the input.	required
`metrics`	`list[str]`	`list[str]` — one hierarchical screen per metric; return dict keyed by `label`.	required
`group`	`str`	Single key naming the group axis.	required
`family_size`	`int \| Mapping[tuple[Any, ...], int] \| None`	Complete inner-family size per group. A scalar applies to every group; a mapping declares heterogeneous sizes keyed like the returned `family_size`. Declared-but-unsubmitted members are inert non-rejections in both the inner BHY and group Simes value. `None` uses the submitted member count.	`None`
`q`	`float`	Nominal FDR target. The outer layer runs at `q`; the inner layer at the selective `q · R / G`. Must satisfy `0 < q < 1`.	`0.05`
`inactive_policy`	`InactivePolicy`	What `G` and the inner families count. Under `"count"` (default) every declared group is a group: an inactive member stays in its inner family at an inert `p = 1`, and a group made up entirely of inactive members enters the outer layer at a Simes p of 1.0, so `G` is the declared group count. Under `"exclude"` inactive members leave their inner family and a group with neither computed nor declared-unsubmitted members leaves `G`; valid only under an independent or pre-specified activity filter. See :data:`InactivePolicy`.	`'count'`

Returns:

Type	Description
`dict[str, HierarchicalBhyResult]`	`dict[str, HierarchicalBhyResult]` keyed by `label`.

Raises:

Type	Description
`UserInputError`	`group == 'factor'` (it is the identifier); `inactive_policy` not `'count'` / `'exclude'`; only one distinct group value in the input (call `bhy` instead); every result is its own group at `n >= 3`; duplicate `(factor, forward_periods, *params)` identifier; any `_resolve_family` invariant failure.

Warns:

Type	Description
`RuntimeWarning`	More than half of input groups contain a single result — inner BHY on n=1 is a raw cutoff and the outer Simes representative equals that single p-value.

Two-layer false discovery rate (FDR) for factor sets with natural group structure (factor families, regions, sectors). Outer Benjamini-Hochberg-Yekutieli (BHY) on Simes (1986) group representatives, then inner BHY within each selected group at the selective level q · R / G (R groups the outer layer selected of G), per Benjamini-Bogomolov (2014) — the Yekutieli (2008) hierarchical framework with the selection adjustment that makes the member-level guarantee hold.

Measured realised FDR stays under nominal across group counts, sizes, sparsity, selected fractions down to R / G = 0.01, and — the regime that broke the earlier nominal-q variant — small groups of strongly correlated members (0.03–0.04 vs nominal 0.10 at size = 4, rho ≥ 0.5). The selective level costs power only on large groups with dense signal; for those, flat bhy is both higher-power and theorem-backed.

Illustrative

import dataclasses

import factrix as fx
from factrix.metrics import ic

# "Which factor families have signal, and within those, which factors?"
# evaluate() returns dict[str, EvaluationResult]; pull the single result
# and stamp the family label onto it with dataclasses.replace so the
# outer screen can group on params["family"].
def stamp(panel, factor_col, **params):
    res = fx.evaluate(panel, metrics={"ic": ic()}, factor_cols=[factor_col])[factor_col]
    return dataclasses.replace(res, params=params)

results = [
    stamp(panel_mom_1m, "mom_1m", family="momentum"),
    stamp(panel_mom_12m, "mom_12m", family="momentum"),
    stamp(panel_pb, "pb", family="value"),
    stamp(panel_pe, "pe", family="value"),
    # ... + quality, low-vol, etc.
]
survivors = fx.multi_factor.bhy_hierarchical(
    results, metrics=["ic"], group="family", q=0.05
)

Which function fits this question?¶

Same input shape (one result per (factor, condition)), three different claims:

Claim	Survivor unit	Function
"Factor X significant in each condition / universe"	`(factor, condition)` pair	`bhy(expand_over=)`
"Factor X significant in \(\ge k\) of \(m\) conditions"	factor identity	`partial_conjunction`
"Which families have signal, and within those, which factors?"	factor identity (group-then-within)	`bhy_hierarchical`

bhy_hierarchical is the only one of the three that keeps the family-level answer first-class — readers learn both "5 of 8 families showed signal" and "within those, factors A / B / C survived" from a single HierarchicalBhyResult.

How the math works¶

Per group \(g\) with \(m_g\) member p-values:

Compute the group representative

\[ p_{\text{Simes},g} = \min_{k=1,\ldots,m_g} \frac{m_g}{k} \cdot p_{(k),g} \]

where \(p_{(k),g}\) is the \(k\)-th smallest p-value in group \(g\). Simes dominates the Bonferroni representative \(m_g \cdot \min(p)\) and is the Yekutieli 2008 recommended choice.

Outer BHY across the \(G\) group representatives gives \(p_{\text{outer},g}^{\text{adj}}\).
Inner BHY within each group gives \(p_{\text{inner},i}^{\text{adj}}\) for member \(i\) of group \(g(i)\).
Let \(R\) be the number of groups with \(p_{\text{outer},g}^{\text{adj}} \le q\). The cell-level adjusted p folds both layers with the Benjamini-Bogomolov selection scale \(G / R\) on the inner term:

\[ p_i^{\text{adj}} = \max\Bigl(p_{\text{outer},g(i)}^{\text{adj}},\; \min\bigl(1,\; p_{\text{inner},i}^{\text{adj}} \cdot G / R\bigr)\Bigr) \]

so that adj_p[i] <= q is exactly "outer passed at \(q\) and inner passed at \(q \cdot R / G\)". Because \(R\) depends on \(q\), the returned array is defined at the supplied \(q\) — re-screen at another level by calling again, not by re-thresholding.

This preserves the universal survivor[i] iff adj_p[i] <= q duality while encoding the two-layer logic: a cell can fail because its group failed outer, because the cell itself failed inner, or both.

Result output¶

bhy_hierarchical returns a HierarchicalBhyResult per metric — the _FdrResultBase shape (entries / survivors / adj_p / q / family_size) plus a hierarchy-specific field:

Field	Meaning
`survivors`	Surviving `EvaluationResult` records, input order (derived: `adj_p_all <= q`)
`adj_p`	Max-of-layers \(\text{adj}_p\) for the survivors; survivor iff `adj_p <= q`
`q`	The `q` you passed (single target, both layers)
`group`	Context key naming the group axis (a single `str`)
`family_size`	Mapping `(group_value,) -> m_group` for every observed group that entered the outer layer, so `G = len(family_size)`. Each `m_group` includes declared-unsubmitted members. Under the default `inactive_policy="count"`, an all-inactive group enters at a Simes p of 1.0; under `"exclude"`, submitted inactive members leave but declared-unsubmitted members remain. See declared family size.
`family`	Declared / computed / inactive / unsubmitted / adjusted member counts and the `inactive_policy` used

Per-survivor group label: survivor.params[result.group].

to_frame() exports factor | forward_periods | <sorted params> | adj_p | survived for every input hypothesis, not only survivors. The group parameter is therefore retained alongside every other swept knob.

factrix.multi_factor.HierarchicalBhyResult `dataclass` ¶

HierarchicalBhyResult(metric_name: str, entries: list[EvaluationResult], adj_p_all: ndarray, q: float, family_size: Mapping[tuple[Any, ...], int], family: FamilyAccounting, group: str, n_selected_groups: int)

Bases: _FdrResultBase

Two-stage hierarchical BHY survivors for one metric.

Shares metric_name / entries / adj_p_all / q / family_size with :class:_FdrResultBase. adj_p_all is max(outer_adj_p[group], min(1, inner_adj_p[i] · G / R)) aligned with entries so entry[i] survives iff adj_p_all[i] <= q — at the q the screen ran at, since R (groups selected by the outer layer) depends on it; q is shared by both layers; family_size is the per-group inner family size keyed by (group_value,) — covering every group that entered the outer layer, not just the survivors, so G is len(family_size). Under inactive_policy="count" that is every declared group; under "exclude" submitted inactive members leave their inner family and a group with neither a computed nor declared-unsubmitted member is absent from family_size and carries NaN.

Attributes:

Name	Type	Description
`group`	`str`	`params` key naming the group axis.
`n_selected_groups`	`int`	`R` — groups the outer layer passed at `q`. Reported because the inner level is `q · R / G` and `R` is the one term in it not otherwise recoverable from this result (`G` is `len(family_size)`), so without it the adjusted p-values cannot be reproduced by hand. `0` means no group passed and nothing survives.

When not to reach for `bhy_hierarchical`¶

Real intent	Reach for	Why
No natural group structure	`bhy`	The grouping is real or it isn't; faking a group axis trivializes the procedure.
"Factor X passes in every condition"	`partial_conjunction` with `min_pass == m`	Hierarchical is "group-then-within", not "joint across conditions".
Flat BHY split by family for display only	`bhy(expand_over=["family"])`	Independent step-ups per bucket, no group-level inference. Use when you do not need a "this family has signal" answer.
Mixed-sign factors in one bucket	Split the bucket / pre-orthogonalize	Within-group Simes assumes positive regression dependence on a subset (PRDS); structurally opposite factors (e.g. momentum + reversal in one group) can violate it.

Validation summary¶

Trigger	Outcome
`group` shadows an identity field (`factor` / `forward_periods`)	`UserInputError`.
`group` key missing from a result's `params`	`UserInputError`.
Only one distinct group value across input	`UserInputError` — points at `bhy`.
Every result is its own group at \(n \ge 3\) (group axis near-unique)	`UserInputError` — pick a coarser categorical.
Duplicate `(identity, group_value)` partition key	`UserInputError`.
More than half of input groups contain a single result	`RuntimeWarning` — inner BHY on \(n=1\) is a raw cutoff.

Caveats¶

Simes outer representative: not exposed as a kwarg. Dominates Bonferroni-min under PRDS; Edgington-style mean-p has no valid null distribution and is rejected.
PRDS within group: Simes is valid under positive regression dependence — typical for factors within one family that share style exposure. If a group mixes structurally opposite factors (e.g. momentum + reversal in one bucket), the within-group PRDS assumption can fail; split the group or pre-orthogonalize.
Pre-filtered input: bhy_hierarchical assumes the input is the candidate family. If results came from upstream pre-filtering (e.g. top-50 of 500 candidates), the FDR claim does not cover the full screening pipeline — track \(K\) per the experiment-log discipline.
Composed FDR is approximate at exact \(q\): Yekutieli 2008 bounds group-level FDR \(\le q\) and within-group FDR \(\le q\) conditional on group passing; the composed per-hypothesis FDR under PRDS is bounded but not exactly \(q\). Researcher claims should be "FDR-controlled at \(q\) in each layer", not "joint FDR \(= q\)".

References¶

[S1986] Simes, R. J. (1986). An improved Bonferroni procedure for multiple tests of significance. Biometrika, 73(3), 751–754.
[Y2008] Yekutieli, D. (2008). Hierarchical false discovery rate-controlling methodology. JASA, 103(481), 309–316.
[NBER34050] NBER WP 34050 (2025). Hierarchical Multiple Testing in Empirical Asset Pricing.