From 26888e294e30bdfde47de8825d7246f1f524b0b0 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Fri, 25 Sep 2026 17:24:10 -0400
Subject: [PATCH 01/15] feat: implement standalone Information-LASSO pathway

---
 causationentropy/core/discovery.py            | 124 ++++++++++++---
 .../tests/test_information_lasso.py           | 149 ++++++++++++++++++
 .../theory/methods/information_lasso.rst      |  31 ++++
 3 files changed, 286 insertions(+), 18 deletions(-)
 create mode 100644 causationentropy/tests/test_information_lasso.py

diff --git a/causationentropy/core/discovery.py b/causationentropy/core/discovery.py
index 3d8b6d7..7418c24 100644
--- a/causationentropy/core/discovery.py
+++ b/causationentropy/core/discovery.py
@@ -10,7 +10,7 @@
 import networkx as nx
 import numpy as np
 import pandas as pd
-from sklearn.linear_model import Lasso, LassoLarsIC
+from sklearn.linear_model import Lasso, LassoCV, LassoLarsIC
 
 from causationentropy.core.information.conditional_mutual_information import (
     conditional_mutual_information,
@@ -249,7 +249,15 @@ def discover_network(
                 bandwidth,
             )
         if method == "information_lasso":
-            S = information_lasso_optimal_causation_entropy(X_lagged, Y, rng)
+            S = information_lasso_optimal_causation_entropy(
+                X_lagged,
+                Y,
+                rng,
+                information=information,
+                metric=metric,
+                k_means=k_means,
+                bandwidth=bandwidth,
+            )
         if method == "lasso":
             S = lasso_optimal_causation_entropy(X_lagged, Y, rng)
         for s in S:
@@ -492,14 +500,40 @@ def alternative_optimal_causation_entropy(
 
 
 def information_lasso_optimal_causation_entropy(
-    X, Y, rng, criterion="bic", max_lambda=100, cross_val=10, information="gaussian"
+    X,
+    Y,
+    rng,
+    criterion="bic",
+    max_lambda=100,
+    cross_val=10,
+    information="gaussian",
+    metric="euclidean",
+    k_means=5,
+    bandwidth="silverman",
 ):
-    """
-    Execute information-theoretic variant of oCSE with LASSO regularization.
+    r"""
+    Select predictors with information-weighted LASSO.
+
+    Each candidate first receives a non-negative information weight from its
+    conditional mutual information with the target (with no additional
+    conditioning set in this standalone pathway):
+
+    .. math::
+
+        w_j = \frac{I(X_j; Y)}{\sum_k I(X_k; Y)}.
+
+    The weighted-LASSO objective is
+
+    .. math::
 
-    This method combines information-theoretic causal discovery with LASSO regularization
-    to handle high-dimensional predictor spaces. The approach balances causal relationship
-    strength with model complexity.
+        \min_{\boldsymbol{\beta}}
+        \frac{1}{2n}\|\mathbf{y}-\mathbf{X}\boldsymbol{\beta}\|_2^2
+        + \lambda \sum_j \frac{1}{w_j}|\beta_j|.
+
+    It is solved with the standard scikit-learn LASSO solvers by scaling
+    feature j by w_j. This is the usual reparameterization of a weighted
+    l1 penalty. Candidates with zero information weight receive an infinite
+    effective penalty and cannot enter the model.
 
     Parameters
     ----------
@@ -508,29 +542,83 @@ def information_lasso_optimal_causation_entropy(
     Y : array-like of shape (T, 1)
         Target variable column.
     rng : numpy.random.Generator
-        Random number generator.
+        Random number generator retained for API consistency. The current
+        weighted-LASSO solvers are deterministic.
     criterion : str, default='bic'
-        Information criterion for model selection ('bic' or 'aic').
+        Information criterion used by LassoLarsIC when there are enough
+        samples relative to predictors.
     max_lambda : int, default=100
-        Maximum number of LASSO iterations.
+        Maximum solver iterations for LassoLarsIC. LassoCV uses at least 1000
+        iterations to avoid premature convergence in high dimensions.
     cross_val : int, default=10
-        Cross-validation folds (currently unused).
+        Number of cross-validation folds used by LassoCV in the
+        high-dimensional regime.
     information : str, default='gaussian'
-        Information measure estimator type.
+        Information estimator passed to conditional_mutual_information.
+    metric : str, default='euclidean'
+        Distance metric used by k-NN information estimators.
+    k_means : int, default=5
+        Number of neighbors used by k-NN information estimators.
+    bandwidth : str or float, default='silverman'
+        Bandwidth used by KDE information estimators.
 
     Returns
     -------
     S : list of int
-        Indices of selected predictor variables.
+        Indices of predictors with non-zero weighted-LASSO coefficients.
 
     Notes
     -----
-    This is a simplified implementation that delegates to LASSO. Future versions
-    will incorporate information-theoretic weighting into the regularization.
+    LassoLarsIC is used when the information criterion is identifiable from
+    the sample size. When predictors are at least as numerous as the available
+    samples, LassoCV uses the existing cross_val parameter to choose the
+    regularization strength instead of falling back to a fixed default penalty.
     """
+    del rng  # Retained in the public signature for consistency with other methods.
+
+    n_features = X.shape[1]
+    information_values = np.empty(n_features, dtype=float)
+
+    for j in range(n_features):
+        value = conditional_mutual_information(
+            X[:, [j]],
+            Y,
+            None,
+            method=information,
+            metric=metric,
+            k=k_means,
+            bandwidth=bandwidth,
+        )
+        if not np.isfinite(value):
+            raise ValueError(
+                "Information-LASSO received a non-finite information weight "
+                f"for predictor {j}."
+            )
+        information_values[j] = max(0.0, value)
+
+    total_information = information_values.sum()
+    if total_information <= 0:
+        return []
+
+    weights = information_values / total_information
+    X_weighted = X * weights.reshape(1, -1)
+
+    if X.shape[0] > n_features + 1:
+        model = LassoLarsIC(criterion=criterion, max_iter=max_lambda).fit(
+            X_weighted, Y.flatten()
+        )
+    else:
+        cv_folds = min(cross_val, X.shape[0])
+        if cv_folds < 2:
+            raise ValueError(
+                "Information-LASSO requires at least two samples for "
+                "high-dimensional cross-validation."
+            )
+        model = LassoCV(cv=cv_folds, max_iter=max(1000, max_lambda)).fit(
+            X_weighted, Y.flatten()
+        )
 
-    # This is a simplified implementation - needs proper information-theoretic weighting
-    return lasso_optimal_causation_entropy(X, Y, rng, criterion, max_lambda, cross_val)
+    return np.flatnonzero(model.coef_ != 0).tolist()
 
 
 def lasso_optimal_causation_entropy(
diff --git a/causationentropy/tests/test_information_lasso.py b/causationentropy/tests/test_information_lasso.py
new file mode 100644
index 0000000..fbd4dbd
--- /dev/null
+++ b/causationentropy/tests/test_information_lasso.py
@@ -0,0 +1,149 @@
+from unittest.mock import MagicMock, patch
+
+import numpy as np
+import pytest
+
+from causationentropy.core.discovery import (
+    discover_network,
+    information_lasso_optimal_causation_entropy,
+)
+
+
+def test_information_lasso_scales_features_by_information_weights():
+    """CMI weights must alter the sparse-selection design matrix."""
+    X = np.arange(10, dtype=float).reshape(5, 2) + 1.0
+    Y = np.arange(5, dtype=float).reshape(-1, 1)
+    rng = np.random.default_rng(0)
+
+    fitted = MagicMock()
+    fitted.fit.return_value = fitted
+    fitted.coef_ = np.array([1.0, 0.0])
+
+    with patch(
+        "causationentropy.core.discovery.conditional_mutual_information",
+        side_effect=[3.0, 1.0],
+    ) as cmi, patch(
+        "causationentropy.core.discovery.LassoLarsIC",
+        return_value=fitted,
+    ):
+        selected = information_lasso_optimal_causation_entropy(
+            X,
+            Y,
+            rng,
+            information="kde",
+            metric="cityblock",
+            k_means=7,
+            bandwidth=0.5,
+        )
+
+    assert selected == [0]
+    np.testing.assert_allclose(
+        fitted.fit.call_args.args[0],
+        X * np.array([0.75, 0.25]),
+    )
+    for call in cmi.call_args_list:
+        assert call.args[2] is None
+        assert call.kwargs == {
+            "method": "kde",
+            "metric": "cityblock",
+            "k": 7,
+            "bandwidth": 0.5,
+        }
+
+
+def test_information_lasso_uses_cross_validation_when_high_dimensional():
+    """p >= n uses the existing cross_val parameter instead of fixed-alpha LASSO."""
+    X = np.arange(20, dtype=float).reshape(4, 5) + 1.0
+    Y = np.arange(4, dtype=float).reshape(-1, 1)
+    rng = np.random.default_rng(0)
+
+    fitted = MagicMock()
+    fitted.fit.return_value = fitted
+    fitted.coef_ = np.array([0.0, 1.0, 0.0, 0.0, -1.0])
+
+    with patch(
+        "causationentropy.core.discovery.conditional_mutual_information",
+        side_effect=[5.0, 4.0, 3.0, 2.0, 1.0],
+    ), patch(
+        "causationentropy.core.discovery.LassoCV",
+        return_value=fitted,
+    ) as lasso_cv:
+        selected = information_lasso_optimal_causation_entropy(
+            X,
+            Y,
+            rng,
+            max_lambda=77,
+            cross_val=10,
+        )
+
+    assert selected == [1, 4]
+    lasso_cv.assert_called_once_with(cv=4, max_iter=1000)
+
+
+def test_information_lasso_returns_empty_when_all_information_is_zero():
+    """No information signal should not fall back to ordinary LASSO."""
+    X = np.ones((6, 3))
+    Y = np.arange(6, dtype=float).reshape(-1, 1)
+
+    with patch(
+        "causationentropy.core.discovery.conditional_mutual_information",
+        return_value=0.0,
+    ), patch("causationentropy.core.discovery.LassoLarsIC") as lars, patch(
+        "causationentropy.core.discovery.LassoCV"
+    ) as lasso_cv:
+        selected = information_lasso_optimal_causation_entropy(
+            X,
+            Y,
+            np.random.default_rng(0),
+        )
+
+    assert selected == []
+    lars.assert_not_called()
+    lasso_cv.assert_not_called()
+
+
+def test_information_lasso_rejects_nonfinite_information_weight():
+    """Estimator failures should surface rather than silently changing selection."""
+    X = np.ones((6, 2))
+    Y = np.arange(6, dtype=float).reshape(-1, 1)
+
+    with patch(
+        "causationentropy.core.discovery.conditional_mutual_information",
+        side_effect=[1.0, np.nan],
+    ):
+        with pytest.raises(ValueError, match="non-finite information weight"):
+            information_lasso_optimal_causation_entropy(
+                X,
+                Y,
+                np.random.default_rng(0),
+            )
+
+
+def test_discover_network_forwards_information_lasso_estimator_options():
+    """The public API estimator settings must reach Information-LASSO."""
+    data = np.arange(36, dtype=float).reshape(12, 3)
+
+    with patch(
+        "causationentropy.core.discovery.information_lasso_optimal_causation_entropy",
+        return_value=[],
+    ) as information_lasso:
+        graph = discover_network(
+            data,
+            method="information_lasso",
+            information="kde",
+            max_lag=1,
+            metric="cityblock",
+            bandwidth=0.4,
+            k_means=9,
+            n_shuffles=2,
+        )
+
+    assert graph.number_of_nodes() == 3
+    assert information_lasso.call_count == 3
+    for call in information_lasso.call_args_list:
+        assert call.kwargs == {
+            "information": "kde",
+            "metric": "cityblock",
+            "k_means": 9,
+            "bandwidth": 0.4,
+        }
diff --git a/docs/source/theory/methods/information_lasso.rst b/docs/source/theory/methods/information_lasso.rst
index d335930..066a1a4 100644
--- a/docs/source/theory/methods/information_lasso.rst
+++ b/docs/source/theory/methods/information_lasso.rst
@@ -175,6 +175,37 @@ Two-Stage Implementation
        selected = np.where(np.abs(beta_original) > 1e-6)[0]
        return selected, beta_original
 
+Current Standalone Implementation
+---------------------------------
+
+The public method="information_lasso" pathway implements the two-stage
+weighted-LASSO approach above.
+
+For each lagged candidate, the standalone pathway computes
+
+.. math::
+
+   w_j = \frac{I(X_j;Y)}{\sum_k I(X_k;Y)}
+
+with no additional conditioning set. The normalized information weights are
+then used as predictor-specific penalty weights. Equivalently, scaling column
+j of the design matrix by w_j lets the existing LASSO solvers optimize
+
+.. math::
+
+   \frac{1}{2n}\|\mathbf{y}-\mathbf{X}\boldsymbol{\beta}\|_2^2
+   + \lambda \sum_j \frac{1}{w_j}|\beta_j|.
+
+Candidates with zero information receive an infinite effective penalty and are
+not selected. When there are enough samples relative to predictors,
+LassoLarsIC selects the regularization strength using AIC or BIC. In the
+high-dimensional regime, LassoCV uses the existing cross_val parameter.
+
+The information estimator selected through discover_network (for example
+Gaussian, KDE, or k-NN) is also used to construct these weights. Iterative
+reweighting and the separate Information-LASSO screening plus exact oCSE
+refinement pathway are not part of this standalone implementation.
+
 Adaptive Implementation
 ----------------------
 

From ac4af8e50216d52713e64f774feb97a2cf946371 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Mon, 28 Sep 2026 20:36:07 -0400
Subject: [PATCH 02/15] feat: add clean information-screened oCSE pathway

---
 causationentropy/core/discovery.py | 77 +++++++++++++++++++++++++++---
 1 file changed, 70 insertions(+), 7 deletions(-)

diff --git a/causationentropy/core/discovery.py b/causationentropy/core/discovery.py
index 7418c24..8c00024 100644
--- a/causationentropy/core/discovery.py
+++ b/causationentropy/core/discovery.py
@@ -5,6 +5,7 @@
 """
 
 import copy
+import math
 from typing import Dict, Tuple, Union
 
 import networkx as nx
@@ -28,6 +29,7 @@ def discover_network(
     bandwidth="silverman",
     k_means: int = 5,
     n_shuffles: int = 200,
+    screen_retention: float = 0.40,
     n_jobs=-1,
     random_state: Union[int, np.random.Generator, None] = 42,
     only_return_significant: bool = True,
@@ -75,6 +77,8 @@ def discover_network(
         - 'standard': Uses initial conditioning set of lagged target variables
         - 'alternative': No initial conditioning set
         - 'information_lasso': Information-theoretic variant with LASSO regularization
+        - 'information_screened': Information-LASSO proposal plus conditional rescue
+          followed by restricted standard oCSE refinement
         - 'lasso': Pure LASSO-based selection
     information : str, default='gaussian'
         Information measure estimator type. Options:
@@ -99,6 +103,10 @@ def discover_network(
     n_shuffles : int, default=200
         Number of permutations for statistical significance testing. Higher values
         provide more accurate p-value estimates but increase computational cost.
+    screen_retention : float, default=0.40
+        Target fraction of external lagged candidates retained by the
+        information_screened method before restricted oCSE refinement.
+        Ignored by other methods.
     n_jobs : int, default=-1
         Number of parallel jobs for computation. -1 uses all available processors.
     random_state : int, numpy.random.Generator, or None, default=42
@@ -174,7 +182,13 @@ def discover_network(
     """
     rng = np.random.default_rng(random_state)
 
-    if method not in ["standard", "alternative", "information_lasso", "lasso"]:
+    if method not in [
+        "standard",
+        "alternative",
+        "information_lasso",
+        "information_screened",
+        "lasso",
+    ]:
         raise NotImplementedError(f"discover_network: method={method} not supported.")
     supported_information_types = ["gaussian", "knn", "kde", "geometric_knn", "poisson"]
     if information not in supported_information_types:
@@ -217,6 +231,8 @@ def discover_network(
         print(f"Estimating edges for node {i} ({var_names[i]})")
 
         Y = Y_all[:, [i]]  # shape: (T - max_lag, 1)
+        output_base = None
+        report_candidate_ids = None
         if method == "standard":
             Z_init = []
             for tau in range(1, max_lag + 1):
@@ -258,6 +274,35 @@ def discover_network(
                 k_means=k_means,
                 bandwidth=bandwidth,
             )
+        if method == "information_screened":
+            Z_init = np.column_stack(
+                [
+                    series[max_lag - tau : T - tau, i]
+                    for tau in range(1, max_lag + 1)
+                ]
+            )
+            candidate_ids = [
+                idx
+                for idx, (source, _lag) in enumerate(feature_names)
+                if source != i
+            ]
+            screened_local = information_screened_optimal_causation_entropy(
+                X_lagged[:, candidate_ids],
+                Y,
+                Z_init,
+                rng,
+                retention=screen_retention,
+                alpha_forward=alpha_forward,
+                alpha_backward=alpha_backward,
+                n_shuffles=n_shuffles,
+                information=information,
+                metric=metric,
+                k_means=k_means,
+                bandwidth=bandwidth,
+            )
+            S = [candidate_ids[int(local_idx)] for local_idx in screened_local]
+            output_base = Z_init
+            report_candidate_ids = candidate_ids
         if method == "lasso":
             S = lasso_optimal_causation_entropy(X_lagged, Y, rng)
         for s in S:
@@ -267,9 +312,17 @@ def discover_network(
             X_predictor = X_lagged[:, [s]]  # predictor at this lag
             Y_target = Y  # target variable
 
-            # Conditioning set: all other selected predictors for this target
+            # Conditioning set: all other selected predictors for this target.
+            # Screened oCSE also keeps the target-history baseline used by the
+            # restricted refinement.
             other_selected = [idx for idx in S if idx != s]
-            Z_cond = X_lagged[:, other_selected] if other_selected else None
+            Z_selected = X_lagged[:, other_selected] if other_selected else None
+            if output_base is None:
+                Z_cond = Z_selected
+            elif Z_selected is None:
+                Z_cond = output_base
+            else:
+                Z_cond = np.hstack((output_base, Z_selected))
 
             # Compute conditional mutual information
             cmi = conditional_mutual_information(
@@ -322,12 +375,23 @@ def discover_network(
             # candidate itself when it is selected).
             selected = set(S)
             n_features = len(feature_names)
-            for cand in range(n_features):
+            candidates_to_report = (
+                report_candidate_ids
+                if report_candidate_ids is not None
+                else range(n_features)
+            )
+            for cand in candidates_to_report:
                 if cand in selected:
                     continue
                 src_var, src_lag = feature_names[cand]
                 X_predictor = X_lagged[:, [cand]]
-                Z_cond = X_lagged[:, S] if S else None
+                Z_selected = X_lagged[:, S] if S else None
+                if output_base is None:
+                    Z_cond = Z_selected
+                elif Z_selected is None:
+                    Z_cond = output_base
+                else:
+                    Z_cond = np.hstack((output_base, Z_selected))
 
                 cmi = conditional_mutual_information(
                     X_predictor,
@@ -620,8 +684,7 @@ def information_lasso_optimal_causation_entropy(
 
     return np.flatnonzero(model.coef_ != 0).tolist()
 
-
-def lasso_optimal_causation_entropy(
+\ndef information_screened_optimal_causation_entropy(\n    X,\n    Y,\n    Z_init,\n    rng,\n    retention=0.40,\n    alpha_forward=0.05,\n    alpha_backward=0.05,\n    n_shuffles=200,\n    information="gaussian",\n    metric="euclidean",\n    k_means=5,\n    bandwidth="silverman",\n):\n    r"""Screen with Information-LASSO plus conditional rescue, then refine.\n\n    X contains only external candidate predictors for one target. Target\n    history is supplied separately through Z_init and remains in the standard\n    oCSE conditioning set during restricted refinement.\n    """\n    if not 0 < retention <= 1:\n        raise ValueError("retention must be in (0, 1].")\n\n    n_features = X.shape[1]\n    if n_features == 0:\n        return []\n\n    endpoint = information_lasso_optimal_causation_entropy(\n        X,\n        Y,\n        rng,\n        information=information,\n        metric=metric,\n        k_means=k_means,\n        bandwidth=bandwidth,\n    )\n    endpoint = [int(j) for j in endpoint]\n    selected = set(endpoint)\n\n    target_size = int(math.ceil(retention * n_features))\n    if len(selected) < n_features:\n        target_size = max(target_size, len(selected) + 1)\n    target_size = min(n_features, target_size)\n\n    if len(selected) < target_size:\n        Z_rescue = X[:, endpoint] if endpoint else None\n        scored = []\n        for j in range(n_features):\n            if j in selected:\n                continue\n            score = conditional_mutual_information(\n                X[:, [j]],\n                Y,\n                Z_rescue,\n                method=information,\n                metric=metric,\n                k=k_means,\n                bandwidth=bandwidth,\n            )\n            if not np.isfinite(score):\n                score = -np.inf\n            scored.append((float(score), int(j)))\n\n        scored.sort(key=lambda item: (-item[0], item[1]))\n        selected.update(\n            j for _score, j in scored[: target_size - len(selected)]\n        )\n\n    screened = sorted(selected)\n    refined_local = standard_optimal_causation_entropy(\n        X[:, screened],\n        Y,\n        Z_init,\n        rng,\n        alpha1=alpha_forward,\n        alpha2=alpha_backward,\n        n_shuffles=n_shuffles,\n        information=information,\n        metric=metric,\n        k_means=k_means,\n        bandwidth=bandwidth,\n    )\n    return [screened[int(local_idx)] for local_idx in refined_local]\n\ndef lasso_optimal_causation_entropy(
     X, Y, rng, criterion="bic", max_lambda=100, cross_val=10
 ):
     r"""

From cb32aa637f171203ac6771a851960a94c6dbd537 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Mon, 28 Sep 2026 20:36:35 -0400
Subject: [PATCH 03/15] test: include information-screened discovery method

---
 causationentropy/tests/test_discovery.py | 8 +++++++-
 1 file changed, 7 insertions(+), 1 deletion(-)

diff --git a/causationentropy/tests/test_discovery.py b/causationentropy/tests/test_discovery.py
index 96c5a65..0168d7d 100644
--- a/causationentropy/tests/test_discovery.py
+++ b/causationentropy/tests/test_discovery.py
@@ -88,7 +88,13 @@ def test_discover_network_valid_methods(self):
         """Test that all valid methods are accepted."""
         data = np.random.normal(0, 1, (30, 2))
 
-        valid_methods = ["standard", "alternative", "information_lasso", "lasso"]
+        valid_methods = [
+            "standard",
+            "alternative",
+            "information_lasso",
+            "information_screened",
+            "lasso",
+        ]
         for method in valid_methods:
             G = discover_network(data, method=method, max_lag=1, n_shuffles=10)
             assert isinstance(G, nx.MultiDiGraph)

From 61d2e67d7e2624a43a128894af371b3d2611de56 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Mon, 28 Sep 2026 20:36:37 -0400
Subject: [PATCH 04/15] test: cover clean information-screened oCSE path

---
 .../tests/test_information_screened_ocse.py   | 104 ++++++++++++++++++
 1 file changed, 104 insertions(+)
 create mode 100644 causationentropy/tests/test_information_screened_ocse.py

diff --git a/causationentropy/tests/test_information_screened_ocse.py b/causationentropy/tests/test_information_screened_ocse.py
new file mode 100644
index 0000000..7b00df9
--- /dev/null
+++ b/causationentropy/tests/test_information_screened_ocse.py
@@ -0,0 +1,104 @@
+import numpy as np
+import pytest
+
+import causationentropy.core.discovery as discovery
+
+
+def test_information_screened_rescue_then_maps_refined_support(monkeypatch):
+    X = np.column_stack(
+        [
+            np.linspace(0.0, 1.0, 20),
+            np.linspace(1.0, 2.0, 20),
+            np.linspace(2.0, 3.0, 20),
+            np.linspace(3.0, 4.0, 20),
+        ]
+    )
+    Y = np.linspace(-1.0, 1.0, 20).reshape(-1, 1)
+    Z_init = np.ones((20, 1))
+
+    monkeypatch.setattr(
+        discovery,
+        "information_lasso_optimal_causation_entropy",
+        lambda *args, **kwargs: [0],
+    )
+
+    def fake_cmi(x, _y, _z, **kwargs):
+        # Candidate 2 is the strongest excluded rescue candidate.
+        first = float(x[0, 0])
+        return {1.0: 1.0, 2.0: 5.0, 3.0: 2.0}.get(first, 0.0)
+
+    monkeypatch.setattr(discovery, "conditional_mutual_information", fake_cmi)
+
+    captured = {}
+
+    def fake_standard(x, _y, z_init, _rng, **kwargs):
+        captured["screened"] = x.copy()
+        captured["z_init"] = z_init.copy()
+        # Return the second column of the screened matrix.
+        return [1]
+
+    monkeypatch.setattr(
+        discovery,
+        "standard_optimal_causation_entropy",
+        fake_standard,
+    )
+
+    selected = discovery.information_screened_optimal_causation_entropy(
+        X,
+        Y,
+        Z_init,
+        np.random.default_rng(0),
+        retention=0.5,
+    )
+
+    assert selected == [2]
+    np.testing.assert_allclose(captured["screened"], X[:, [0, 2]])
+    np.testing.assert_allclose(captured["z_init"], Z_init)
+
+
+def test_information_screened_rejects_invalid_retention():
+    X = np.ones((20, 2))
+    Y = np.ones((20, 1))
+    Z = np.ones((20, 1))
+
+    for value in (0.0, -0.1, 1.1):
+        with pytest.raises(ValueError, match="retention"):
+            discovery.information_screened_optimal_causation_entropy(
+                X,
+                Y,
+                Z,
+                np.random.default_rng(0),
+                retention=value,
+            )
+
+
+def test_discover_information_screened_uses_external_candidates(monkeypatch):
+    captured = []
+
+    def fake_screened(X, Y, Z_init, rng, **kwargs):
+        captured.append((X.shape, Y.shape, Z_init.shape, kwargs["retention"]))
+        return []
+
+    monkeypatch.setattr(
+        discovery,
+        "information_screened_optimal_causation_entropy",
+        fake_screened,
+    )
+
+    data = np.random.default_rng(0).normal(size=(40, 3))
+    graph = discovery.discover_network(
+        data,
+        method="information_screened",
+        max_lag=2,
+        screen_retention=0.35,
+        n_shuffles=3,
+        random_state=0,
+    )
+
+    assert len(graph) == 3
+    assert len(captured) == 3
+    for x_shape, y_shape, z_shape, retention in captured:
+        assert x_shape == (38, 4)  # (n-1) * max_lag external candidates
+        assert y_shape == (38, 1)
+        assert z_shape == (38, 2)
+        assert retention == 0.35

From 195e80445224f5b3f602168a043b4db77f4bced6 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Mon, 28 Sep 2026 20:39:12 -0400
Subject: [PATCH 05/15] fix: restore source newlines in screened-oCSE helper

---
 causationentropy/core/discovery.py | 87 +++++++++++++++++++++++++++++-
 1 file changed, 86 insertions(+), 1 deletion(-)

diff --git a/causationentropy/core/discovery.py b/causationentropy/core/discovery.py
index 8c00024..9433101 100644
--- a/causationentropy/core/discovery.py
+++ b/causationentropy/core/discovery.py
@@ -684,7 +684,92 @@ def information_lasso_optimal_causation_entropy(
 
     return np.flatnonzero(model.coef_ != 0).tolist()
 
-\ndef information_screened_optimal_causation_entropy(\n    X,\n    Y,\n    Z_init,\n    rng,\n    retention=0.40,\n    alpha_forward=0.05,\n    alpha_backward=0.05,\n    n_shuffles=200,\n    information="gaussian",\n    metric="euclidean",\n    k_means=5,\n    bandwidth="silverman",\n):\n    r"""Screen with Information-LASSO plus conditional rescue, then refine.\n\n    X contains only external candidate predictors for one target. Target\n    history is supplied separately through Z_init and remains in the standard\n    oCSE conditioning set during restricted refinement.\n    """\n    if not 0 < retention <= 1:\n        raise ValueError("retention must be in (0, 1].")\n\n    n_features = X.shape[1]\n    if n_features == 0:\n        return []\n\n    endpoint = information_lasso_optimal_causation_entropy(\n        X,\n        Y,\n        rng,\n        information=information,\n        metric=metric,\n        k_means=k_means,\n        bandwidth=bandwidth,\n    )\n    endpoint = [int(j) for j in endpoint]\n    selected = set(endpoint)\n\n    target_size = int(math.ceil(retention * n_features))\n    if len(selected) < n_features:\n        target_size = max(target_size, len(selected) + 1)\n    target_size = min(n_features, target_size)\n\n    if len(selected) < target_size:\n        Z_rescue = X[:, endpoint] if endpoint else None\n        scored = []\n        for j in range(n_features):\n            if j in selected:\n                continue\n            score = conditional_mutual_information(\n                X[:, [j]],\n                Y,\n                Z_rescue,\n                method=information,\n                metric=metric,\n                k=k_means,\n                bandwidth=bandwidth,\n            )\n            if not np.isfinite(score):\n                score = -np.inf\n            scored.append((float(score), int(j)))\n\n        scored.sort(key=lambda item: (-item[0], item[1]))\n        selected.update(\n            j for _score, j in scored[: target_size - len(selected)]\n        )\n\n    screened = sorted(selected)\n    refined_local = standard_optimal_causation_entropy(\n        X[:, screened],\n        Y,\n        Z_init,\n        rng,\n        alpha1=alpha_forward,\n        alpha2=alpha_backward,\n        n_shuffles=n_shuffles,\n        information=information,\n        metric=metric,\n        k_means=k_means,\n        bandwidth=bandwidth,\n    )\n    return [screened[int(local_idx)] for local_idx in refined_local]\n\ndef lasso_optimal_causation_entropy(
+
+def information_screened_optimal_causation_entropy(
+    X,
+    Y,
+    Z_init,
+    rng,
+    retention=0.40,
+    alpha_forward=0.05,
+    alpha_backward=0.05,
+    n_shuffles=200,
+    information="gaussian",
+    metric="euclidean",
+    k_means=5,
+    bandwidth="silverman",
+):
+    r"""Screen with Information-LASSO plus conditional rescue, then refine.
+
+    X contains only external candidate predictors for one target. Target
+    history is supplied separately through Z_init and remains in the standard
+    oCSE conditioning set during restricted refinement.
+    """
+    if not 0 < retention <= 1:
+        raise ValueError("retention must be in (0, 1].")
+
+    n_features = X.shape[1]
+    if n_features == 0:
+        return []
+
+    endpoint = information_lasso_optimal_causation_entropy(
+        X,
+        Y,
+        rng,
+        information=information,
+        metric=metric,
+        k_means=k_means,
+        bandwidth=bandwidth,
+    )
+    endpoint = [int(j) for j in endpoint]
+    selected = set(endpoint)
+
+    target_size = int(math.ceil(retention * n_features))
+    if len(selected) < n_features:
+        target_size = max(target_size, len(selected) + 1)
+    target_size = min(n_features, target_size)
+
+    if len(selected) < target_size:
+        Z_rescue = X[:, endpoint] if endpoint else None
+        scored = []
+        for j in range(n_features):
+            if j in selected:
+                continue
+            score = conditional_mutual_information(
+                X[:, [j]],
+                Y,
+                Z_rescue,
+                method=information,
+                metric=metric,
+                k=k_means,
+                bandwidth=bandwidth,
+            )
+            if not np.isfinite(score):
+                score = -np.inf
+            scored.append((float(score), int(j)))
+
+        scored.sort(key=lambda item: (-item[0], item[1]))
+        selected.update(
+            j for _score, j in scored[: target_size - len(selected)]
+        )
+
+    screened = sorted(selected)
+    refined_local = standard_optimal_causation_entropy(
+        X[:, screened],
+        Y,
+        Z_init,
+        rng,
+        alpha1=alpha_forward,
+        alpha2=alpha_backward,
+        n_shuffles=n_shuffles,
+        information=information,
+        metric=metric,
+        k_means=k_means,
+        bandwidth=bandwidth,
+    )
+    return [screened[int(local_idx)] for local_idx in refined_local]
+
+def lasso_optimal_causation_entropy(
     X, Y, rng, criterion="bic", max_lambda=100, cross_val=10
 ):
     r"""

From 78e06ed15a0b133abf2702ef101235a2887504cc Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Mon, 28 Sep 2026 20:44:09 -0400
Subject: [PATCH 06/15] fix: align screened-oCSE with standard conditioning
 semantics

---
 causationentropy/core/discovery.py | 27 +++++++++++++++++++++++----
 1 file changed, 23 insertions(+), 4 deletions(-)

diff --git a/causationentropy/core/discovery.py b/causationentropy/core/discovery.py
index 9433101..d94b0e0 100644
--- a/causationentropy/core/discovery.py
+++ b/causationentropy/core/discovery.py
@@ -231,6 +231,7 @@ def discover_network(
         print(f"Estimating edges for node {i} ({var_names[i]})")
 
         Y = Y_all[:, [i]]  # shape: (T - max_lag, 1)
+        candidate_ids = list(range(len(feature_names)))
         output_base = None
         report_candidate_ids = None
         if method == "standard":
@@ -238,8 +239,15 @@ def discover_network(
             for tau in range(1, max_lag + 1):
                 Z_init.append(series[max_lag - tau : T - tau, i])  # lagged Y_i
             Z_init = np.column_stack(Z_init)  # shape: (T - max_lag, max_lag)
-            S = standard_optimal_causation_entropy(
-                X_lagged,
+            # Target-history lags are already in Z_init; do not retest them.
+            candidate_ids = [
+                idx
+                for idx, (source, _lag) in enumerate(feature_names)
+                if source != i
+            ]
+            X_candidates = X_lagged[:, candidate_ids]
+            S_local = standard_optimal_causation_entropy(
+                X_candidates,
                 Y,
                 Z_init,
                 rng,
@@ -251,6 +259,9 @@ def discover_network(
                 k_means,
                 bandwidth,
             )
+            S = [candidate_ids[idx] for idx in S_local]
+            output_base = Z_init
+            report_candidate_ids = candidate_ids
         if method == "alternative":
             S = alternative_optimal_causation_entropy(
                 X_lagged,
@@ -496,6 +507,7 @@ def standard_optimal_causation_entropy(
         metric,
         k_means,
         bandwidth,
+        Z_init=Z_init,
     )
 
     return S
@@ -1048,6 +1060,7 @@ def backward(
     metric="euclidean",
     k_means=5,
     bandwidth="silverman",
+    Z_init=None,
 ):
     r"""
     Backward elimination phase of optimal Causation Entropy.
@@ -1101,8 +1114,14 @@ def backward(
     S = copy.deepcopy(S_init)  # working copy
 
     for j in rng.permutation(S_init):
-        # conditioning set Z = S \ {j}
-        Z = X_full[:, [k for k in S if k != j]] if len(S) > 1 else None
+        other_ids = [k for k in S if k != j]
+        Z_selected = X_full[:, other_ids] if other_ids else None
+        if Z_init is None:
+            Z = Z_selected
+        elif Z_selected is None:
+            Z = Z_init
+        else:
+            Z = np.hstack((Z_init, Z_selected))
 
         Xj = X_full[:, [j]]
         cmij = conditional_mutual_information(

From f0cb894b2626d395a7d54852943cb9b7a20e215c Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Mon, 28 Sep 2026 20:44:23 -0400
Subject: [PATCH 07/15] test: cover standard conditioning prerequisite

---
 .../tests/test_information_screened_ocse.py   | 70 +++++++++++++++++++
 1 file changed, 70 insertions(+)

diff --git a/causationentropy/tests/test_information_screened_ocse.py b/causationentropy/tests/test_information_screened_ocse.py
index 7b00df9..67255ef 100644
--- a/causationentropy/tests/test_information_screened_ocse.py
+++ b/causationentropy/tests/test_information_screened_ocse.py
@@ -102,3 +102,73 @@ def fake_screened(X, Y, Z_init, rng, **kwargs):
         assert y_shape == (38, 1)
         assert z_shape == (38, 2)
         assert retention == 0.35
+
+
+def test_backward_preserves_initial_conditioning_for_screened_refinement(monkeypatch):
+    rng = np.random.default_rng(7)
+    X = rng.normal(size=(50, 3))
+    Y = rng.normal(size=(50, 1))
+    Z_init = rng.normal(size=(50, 2))
+    seen = []
+
+    def fake_cmi(_x, _y, z, **kwargs):
+        seen.append(z.copy())
+        return 1.0
+
+    monkeypatch.setattr(discovery, "conditional_mutual_information", fake_cmi)
+    monkeypatch.setattr(
+        discovery,
+        "shuffle_test",
+        lambda *args, **kwargs: {
+            "Threshold": 0.0,
+            "Value": 1.0,
+            "Pass": True,
+            "P_value": 0.0,
+        },
+    )
+
+    selected = discovery.backward(
+        X,
+        Y,
+        [0, 1],
+        np.random.default_rng(11),
+        n_shuffles=3,
+        Z_init=Z_init,
+    )
+
+    assert selected == [0, 1]
+    assert len(seen) == 2
+    for z in seen:
+        np.testing.assert_allclose(z[:, :2], Z_init)
+
+
+def test_standard_discovery_does_not_retest_target_history(monkeypatch):
+    calls = []
+
+    def fake_standard(X, Y, Z_init, rng, *args, **kwargs):
+        calls.append((X.copy(), Z_init.copy()))
+        return []
+
+    monkeypatch.setattr(
+        discovery,
+        "standard_optimal_causation_entropy",
+        fake_standard,
+    )
+
+    data = np.arange(120, dtype=float).reshape(40, 3)
+    discovery.discover_network(
+        data,
+        method="standard",
+        max_lag=2,
+        n_shuffles=3,
+    )
+
+    assert len(calls) == 3
+    for X_candidates, Z_init in calls:
+        assert X_candidates.shape == (38, 4)
+        assert Z_init.shape == (38, 2)
+        for candidate in X_candidates.T:
+            assert not any(
+                np.array_equal(candidate, conditioned)
+                for conditioned in Z_init.T
+            )

From 9d51167eb0448f0bbe7c5713b839f1722d3170b2 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Wed, 30 Sep 2026 22:13:58 -0400
Subject: [PATCH 08/15] fix: condition rescue on target history

---
 causationentropy/core/discovery.py | 10 +++++++---
 1 file changed, 7 insertions(+), 3 deletions(-)

diff --git a/causationentropy/core/discovery.py b/causationentropy/core/discovery.py
index d94b0e0..4062cc0 100644
--- a/causationentropy/core/discovery.py
+++ b/causationentropy/core/discovery.py
@@ -714,8 +714,8 @@ def information_screened_optimal_causation_entropy(
     r"""Screen with Information-LASSO plus conditional rescue, then refine.
 
     X contains only external candidate predictors for one target. Target
-    history is supplied separately through Z_init and remains in the standard
-    oCSE conditioning set during restricted refinement.
+    history is supplied separately through Z_init and is preserved in both the
+    conditional-rescue and restricted-oCSE conditioning sets.
     """
     if not 0 < retention <= 1:
         raise ValueError("retention must be in (0, 1].")
@@ -742,7 +742,11 @@ def information_screened_optimal_causation_entropy(
     target_size = min(n_features, target_size)
 
     if len(selected) < target_size:
-        Z_rescue = X[:, endpoint] if endpoint else None
+        Z_rescue = (
+            np.hstack((Z_init, X[:, endpoint]))
+            if endpoint
+            else Z_init
+        )
         scored = []
         for j in range(n_features):
             if j in selected:

From 3e4e63c6a41a91de20b97e360d9885e601d56eea Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Wed, 30 Sep 2026 22:14:18 -0400
Subject: [PATCH 09/15] test: cover rescue conditioning baseline

---
 .../tests/test_information_screened_ocse.py   | 58 ++++++++++++++++++-
 1 file changed, 57 insertions(+), 1 deletion(-)

diff --git a/causationentropy/tests/test_information_screened_ocse.py b/causationentropy/tests/test_information_screened_ocse.py
index 67255ef..db11fe2 100644
--- a/causationentropy/tests/test_information_screened_ocse.py
+++ b/causationentropy/tests/test_information_screened_ocse.py
@@ -22,7 +22,10 @@ def test_information_screened_rescue_then_maps_refined_support(monkeypatch):
         lambda *args, **kwargs: [0],
     )
 
-    def fake_cmi(x, _y, _z, **kwargs):
+    rescue_conditioning = []
+
+    def fake_cmi(x, _y, z, **kwargs):
+        rescue_conditioning.append(z.copy())
         # Candidate 2 is the strongest excluded rescue candidate.
         first = float(x[0, 0])
         return {1.0: 1.0, 2.0: 5.0, 3.0: 2.0}.get(first, 0.0)
@@ -55,6 +58,59 @@ def fake_standard(x, _y, z_init, _rng, **kwargs):
     np.testing.assert_allclose(captured["screened"], X[:, [0, 2]])
     np.testing.assert_allclose(captured["z_init"], Z_init)
 
+    expected_rescue_conditioning = np.hstack((Z_init, X[:, [0]]))
+    assert len(rescue_conditioning) == 3
+    for z in rescue_conditioning:
+        np.testing.assert_allclose(z, expected_rescue_conditioning)
+
+
+def test_information_screened_rescue_uses_history_when_proposal_is_empty(
+    monkeypatch,
+):
+    X = np.column_stack(
+        [
+            np.linspace(0.0, 1.0, 20),
+            np.linspace(1.0, 2.0, 20),
+        ]
+    )
+    Y = np.linspace(-1.0, 1.0, 20).reshape(-1, 1)
+    Z_init = np.column_stack(
+        [
+            np.linspace(4.0, 5.0, 20),
+            np.linspace(5.0, 6.0, 20),
+        ]
+    )
+    rescue_conditioning = []
+
+    monkeypatch.setattr(
+        discovery,
+        "information_lasso_optimal_causation_entropy",
+        lambda *args, **kwargs: [],
+    )
+
+    def fake_cmi(_x, _y, z, **kwargs):
+        rescue_conditioning.append(z.copy())
+        return 1.0
+
+    monkeypatch.setattr(discovery, "conditional_mutual_information", fake_cmi)
+    monkeypatch.setattr(
+        discovery,
+        "standard_optimal_causation_entropy",
+        lambda *args, **kwargs: [],
+    )
+
+    discovery.information_screened_optimal_causation_entropy(
+        X,
+        Y,
+        Z_init,
+        np.random.default_rng(0),
+        retention=0.5,
+    )
+
+    assert len(rescue_conditioning) == 2
+    for z in rescue_conditioning:
+        np.testing.assert_allclose(z, Z_init)
+
 
 def test_information_screened_rejects_invalid_retention():
     X = np.ones((20, 2))

From ab5c4993d396ec906c5768e6c2a3d2a5b113b89e Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Tue, 13 Oct 2026 19:48:29 -0400
Subject: [PATCH 10/15] docs: define screened oCSE retention contract

---
 .../methods/information_screened_ocse.rst     | 125 ++++++++++++++++++
 1 file changed, 125 insertions(+)
 create mode 100644 docs/source/theory/methods/information_screened_ocse.rst

diff --git a/docs/source/theory/methods/information_screened_ocse.rst b/docs/source/theory/methods/information_screened_ocse.rst
new file mode 100644
index 0000000..bc5cce1
--- /dev/null
+++ b/docs/source/theory/methods/information_screened_ocse.rst
@@ -0,0 +1,125 @@
+==========================
+Information-Screened oCSE
+==========================
+
+``method="information_screened"`` uses Information-LASSO as a proposal stage,
+adds a conditional-information rescue pass, and then runs standard oCSE only
+inside the retained candidate set.
+
+The screen is deliberately conservative. Its purpose is not to estimate the
+final parent set directly; it is to reduce the candidate universe while
+preserving variables needed by the downstream oCSE refinement.
+
+Method
+======
+
+For a target variable :math:`Y`, let:
+
+- :math:`C` be the external lagged-candidate set;
+- :math:`H_Y` be the lagged target-history conditioning set;
+- :math:`S_0` be the Information-LASSO proposal;
+- :math:`\widehat{W}` be the retained working set.
+
+The proposal stage first computes
+
+.. math::
+
+   S_0 = \operatorname{InfoLasso}(C \rightarrow Y).
+
+Candidates outside :math:`S_0` are then ranked by conditional mutual
+information
+
+.. math::
+
+   r_j = \widehat{I}(X_j; Y \mid H_Y \cup S_0).
+
+The target of this rescue pass is a parent-superset working set rather than an
+exact graph estimate. In particular, the screening objective is
+
+.. math::
+
+   \Pr\!\left(\operatorname{Pa}(Y) \subseteq \widehat{W}\right) \rightarrow 1,
+
+while keeping :math:`|\widehat{W}|` substantially smaller than :math:`|C|`.
+The final graph decision remains the responsibility of standard oCSE:
+
+.. math::
+
+   \widehat{\operatorname{Pa}}(Y)
+   = \operatorname{oCSE}(Y; \widehat{W}, H_Y).
+
+Retention Is a Soft Target
+==========================
+
+``screen_retention`` is a target fraction, not a hard cap. Let
+:math:`p = |C|` and
+
+.. math::
+
+   B_0 = \left\lceil \text{screen\_retention} \cdot p \right\rceil.
+
+The implementation never truncates the Information-LASSO proposal. When there
+are candidates outside :math:`S_0`, it also keeps at least one rescue slot so
+that the screened pathway cannot silently collapse back to the proposal-only
+path. The realized target size is therefore
+
+.. math::
+
+   B = \min\!\left(p,\; \max(B_0, |S_0| + 1)\right)
+
+when :math:`|S_0| < p`, and :math:`B=p` otherwise.
+
+This asymmetry is intentional. A false positive at the screening stage only
+increases downstream computation, whereas a false negative removes a candidate
+that restricted oCSE can no longer recover. Consequently, preserving the
+proposal and leaving a rescue channel open takes priority over enforcing an
+exact candidate count.
+
+The resulting working set is
+
+.. math::
+
+   \widehat{W}
+   = S_0 \cup \operatorname{Top}_{B-|S_0|}\{r_j : j \notin S_0\}.
+
+Interpretation
+==============
+
+A useful way to view the method is as a two-level search:
+
+.. code-block:: none
+
+   all external lagged candidates
+              |
+              v
+      Information-LASSO proposal
+              |
+              v
+   conditional rescue given target history + proposal
+              |
+              v
+      conservative working set
+              |
+              v
+       standard oCSE refinement
+              |
+              v
+          causal graph
+
+The screening stage therefore trades a small amount of extra candidate
+retention for protection against proposal-stage misses. The downstream oCSE
+stage still performs the significance-based causal selection.
+
+Current Scope
+=============
+
+The current implementation does not claim a finite-sample sure-screening
+theorem and does not guarantee that every causal parent is retained for every
+finite dataset. The relevant empirical quantity is therefore the
+retention/coverage frontier: how much the candidate universe can be reduced
+before parent or full-oCSE-support coverage begins to degrade.
+
+Likewise, the rescue score is univariate conditional CMI for each excluded
+candidate given :math:`H_Y \cup S_0`. Higher-order synergy searches are outside
+the present scope because they would increase screening cost and change the
+computational objective of this pathway.

From 14beca37ce2fa11518069beb71c25e0de4e0d8a7 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Tue, 13 Oct 2026 19:48:41 -0400
Subject: [PATCH 11/15] test: pin screened retention soft-target behavior

---
 ...information_screened_retention_contract.py | 64 +++++++++++++++++++
 1 file changed, 64 insertions(+)
 create mode 100644 causationentropy/tests/test_information_screened_retention_contract.py

diff --git a/causationentropy/tests/test_information_screened_retention_contract.py b/causationentropy/tests/test_information_screened_retention_contract.py
new file mode 100644
index 0000000..d87bbd8
--- /dev/null
+++ b/causationentropy/tests/test_information_screened_retention_contract.py
@@ -0,0 +1,64 @@
+import numpy as np
+
+import causationentropy.core.discovery as discovery
+
+
+def test_information_screened_retention_is_soft_and_preserves_proposal(monkeypatch):
+    X = np.column_stack(
+        [
+            np.linspace(0.0, 1.0, 24),
+            np.linspace(1.0, 2.0, 24),
+            np.linspace(2.0, 3.0, 24),
+            np.linspace(3.0, 4.0, 24),
+            np.linspace(4.0, 5.0, 24),
+        ]
+    )
+    Y = np.linspace(-1.0, 1.0, 24).reshape(-1, 1)
+    Z_init = np.linspace(10.0, 11.0, 24).reshape(-1, 1)
+
+    # The proposal already exceeds the nominal 40% target (3 > ceil(0.4 * 5)=2).
+    # The screened path must preserve it rather than truncating it.
+    monkeypatch.setattr(
+        discovery,
+        "information_lasso_optimal_causation_entropy",
+        lambda *args, **kwargs: [0, 1, 2],
+    )
+
+    rescue_conditioning = []
+
+    def fake_cmi(x, _y, z, **kwargs):
+        rescue_conditioning.append(z.copy())
+        # Candidate 4 wins the single rescue slot over candidate 3.
+        return 2.0 if np.isclose(x[0, 0], 4.0) else 1.0
+
+    monkeypatch.setattr(discovery, "conditional_mutual_information", fake_cmi)
+
+    captured = {}
+
+    def fake_standard(x, _y, z_init, _rng, **kwargs):
+        captured["screened"] = x.copy()
+        captured["z_init"] = z_init.copy()
+        return []
+
+    monkeypatch.setattr(
+        discovery,
+        "standard_optimal_causation_entropy",
+        fake_standard,
+    )
+
+    selected = discovery.information_screened_optimal_causation_entropy(
+        X,
+        Y,
+        Z_init,
+        np.random.default_rng(0),
+        retention=0.40,
+    )
+
+    assert selected == []
+    np.testing.assert_allclose(captured["screened"], X[:, [0, 1, 2, 4]])
+    np.testing.assert_allclose(captured["z_init"], Z_init)
+
+    expected_conditioning = np.hstack((Z_init, X[:, [0, 1, 2]]))
+    assert len(rescue_conditioning) == 2
+    for z in rescue_conditioning:
+        np.testing.assert_allclose(z, expected_conditioning)

From 32ce85d4b9a03c886aa80b87cadad2742d502f7c Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Tue, 13 Oct 2026 19:48:49 -0400
Subject: [PATCH 12/15] docs: link information-screened oCSE guide

---
 docs/source/theory/index.rst | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/docs/source/theory/index.rst b/docs/source/theory/index.rst
index 65dad02..c243b78 100644
--- a/docs/source/theory/index.rst
+++ b/docs/source/theory/index.rst
@@ -12,6 +12,7 @@ including mathematical foundations, algorithmic details, and methodological comp
    methods/standard_oce
    methods/alternative_oce
    methods/information_lasso
+   methods/information_screened_ocse
    methods/lasso_methods
    information_theory
    statistical_foundations
@@ -65,4 +66,5 @@ Next Steps
 - :doc:`methods/standard_oce`: Detailed explanation of the standard oCSE algorithm  
 - :doc:`methods/alternative_oce`: Alternative oCSE formulation without initial conditioning
 - :doc:`methods/information_lasso`: Information-theoretic variant with LASSO regularization
+- :doc:`methods/information_screened_ocse`: Conservative information screening followed by restricted oCSE
 - :doc:`methods/lasso_methods`: Pure LASSO-based approaches for comparison
\ No newline at end of file

From 2667f17e6ae4d2b8164b227b1ce9b92cd689258e Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Tue, 13 Oct 2026 19:50:29 -0400
Subject: [PATCH 13/15] test: pin screened retention soft-target behavior

---
 .../tests/test_information_screened_ocse.py   | 61 +++++++++++++++++++
 1 file changed, 61 insertions(+)

diff --git a/causationentropy/tests/test_information_screened_ocse.py b/causationentropy/tests/test_information_screened_ocse.py
index db11fe2..414ee7c 100644
--- a/causationentropy/tests/test_information_screened_ocse.py
+++ b/causationentropy/tests/test_information_screened_ocse.py
@@ -112,6 +112,67 @@ def fake_cmi(_x, _y, z, **kwargs):
         np.testing.assert_allclose(z, Z_init)
 
 
+def test_information_screened_retention_is_soft_and_preserves_proposal(monkeypatch):
+    X = np.column_stack(
+        [
+            np.linspace(0.0, 1.0, 24),
+            np.linspace(1.0, 2.0, 24),
+            np.linspace(2.0, 3.0, 24),
+            np.linspace(3.0, 4.0, 24),
+            np.linspace(4.0, 5.0, 24),
+        ]
+    )
+    Y = np.linspace(-1.0, 1.0, 24).reshape(-1, 1)
+    Z_init = np.linspace(10.0, 11.0, 24).reshape(-1, 1)
+
+    # The proposal already exceeds the nominal 40% target (3 > ceil(0.4 * 5)=2).
+    # Preserve the proposal and leave one rescue slot rather than truncating it.
+    monkeypatch.setattr(
+        discovery,
+        "information_lasso_optimal_causation_entropy",
+        lambda *args, **kwargs: [0, 1, 2],
+    )
+
+    rescue_conditioning = []
+
+    def fake_cmi(x, _y, z, **kwargs):
+        rescue_conditioning.append(z.copy())
+        # Candidate 4 wins the single rescue slot over candidate 3.
+        return 2.0 if np.isclose(x[0, 0], 4.0) else 1.0
+
+    monkeypatch.setattr(discovery, "conditional_mutual_information", fake_cmi)
+
+    captured = {}
+
+    def fake_standard(x, _y, z_init, _rng, **kwargs):
+        captured["screened"] = x.copy()
+        captured["z_init"] = z_init.copy()
+        return []
+
+    monkeypatch.setattr(
+        discovery,
+        "standard_optimal_causation_entropy",
+        fake_standard,
+    )
+
+    selected = discovery.information_screened_optimal_causation_entropy(
+        X,
+        Y,
+        Z_init,
+        np.random.default_rng(0),
+        retention=0.40,
+    )
+
+    assert selected == []
+    np.testing.assert_allclose(captured["screened"], X[:, [0, 1, 2, 4]])
+    np.testing.assert_allclose(captured["z_init"], Z_init)
+
+    expected_conditioning = np.hstack((Z_init, X[:, [0, 1, 2]]))
+    assert len(rescue_conditioning) == 2
+    for z in rescue_conditioning:
+        np.testing.assert_allclose(z, expected_conditioning)
+
+
 def test_information_screened_rejects_invalid_retention():
     X = np.ones((20, 2))
     Y = np.ones((20, 1))

From 77773931dd6d825a4b00e1d5004df5d94b3b06ae Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Tue, 13 Oct 2026 19:50:37 -0400
Subject: [PATCH 14/15] test: fold retention contract into screened oCSE
 coverage

---
 ...information_screened_retention_contract.py | 64 -------------------
 1 file changed, 64 deletions(-)
 delete mode 100644 causationentropy/tests/test_information_screened_retention_contract.py

diff --git a/causationentropy/tests/test_information_screened_retention_contract.py b/causationentropy/tests/test_information_screened_retention_contract.py
deleted file mode 100644
index d87bbd8..0000000
--- a/causationentropy/tests/test_information_screened_retention_contract.py
+++ /dev/null
@@ -1,64 +0,0 @@
-import numpy as np
-
-import causationentropy.core.discovery as discovery
-
-
-def test_information_screened_retention_is_soft_and_preserves_proposal(monkeypatch):
-    X = np.column_stack(
-        [
-            np.linspace(0.0, 1.0, 24),
-            np.linspace(1.0, 2.0, 24),
-            np.linspace(2.0, 3.0, 24),
-            np.linspace(3.0, 4.0, 24),
-            np.linspace(4.0, 5.0, 24),
-        ]
-    )
-    Y = np.linspace(-1.0, 1.0, 24).reshape(-1, 1)
-    Z_init = np.linspace(10.0, 11.0, 24).reshape(-1, 1)
-
-    # The proposal already exceeds the nominal 40% target (3 > ceil(0.4 * 5)=2).
-    # The screened path must preserve it rather than truncating it.
-    monkeypatch.setattr(
-        discovery,
-        "information_lasso_optimal_causation_entropy",
-        lambda *args, **kwargs: [0, 1, 2],
-    )
-
-    rescue_conditioning = []
-
-    def fake_cmi(x, _y, z, **kwargs):
-        rescue_conditioning.append(z.copy())
-        # Candidate 4 wins the single rescue slot over candidate 3.
-        return 2.0 if np.isclose(x[0, 0], 4.0) else 1.0
-
-    monkeypatch.setattr(discovery, "conditional_mutual_information", fake_cmi)
-
-    captured = {}
-
-    def fake_standard(x, _y, z_init, _rng, **kwargs):
-        captured["screened"] = x.copy()
-        captured["z_init"] = z_init.copy()
-        return []
-
-    monkeypatch.setattr(
-        discovery,
-        "standard_optimal_causation_entropy",
-        fake_standard,
-    )
-
-    selected = discovery.information_screened_optimal_causation_entropy(
-        X,
-        Y,
-        Z_init,
-        np.random.default_rng(0),
-        retention=0.40,
-    )
-
-    assert selected == []
-    np.testing.assert_allclose(captured["screened"], X[:, [0, 1, 2, 4]])
-    np.testing.assert_allclose(captured["z_init"], Z_init)
-
-    expected_conditioning = np.hstack((Z_init, X[:, [0, 1, 2]]))
-    assert len(rescue_conditioning) == 2
-    for z in rescue_conditioning:
-        np.testing.assert_allclose(z, expected_conditioning)

From 21b97f3041a691e7303c7470ae51594acfb5db42 Mon Sep 17 00:00:00 2001
From: mimosapudical <153621177+mimosapudical@users.noreply.github.com>
Date: Tue, 13 Oct 2026 19:50:50 -0400
Subject: [PATCH 15/15] docs: normalize theory index newline

---
 docs/source/theory/index.rst | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/docs/source/theory/index.rst b/docs/source/theory/index.rst
index c243b78..195a60c 100644
--- a/docs/source/theory/index.rst
+++ b/docs/source/theory/index.rst
@@ -67,4 +67,4 @@ Next Steps
 - :doc:`methods/alternative_oce`: Alternative oCSE formulation without initial conditioning
 - :doc:`methods/information_lasso`: Information-theoretic variant with LASSO regularization
 - :doc:`methods/information_screened_ocse`: Conservative information screening followed by restricted oCSE
-- :doc:`methods/lasso_methods`: Pure LASSO-based approaches for comparison
\ No newline at end of file
+- :doc:`methods/lasso_methods`: Pure LASSO-based approaches for comparison