Give each per-variable QRF model its own seed

changelog.d/207.fixed.md

@@ -0,0 +1,2 @@
+Per-variable QRF models now derive distinct seeds, so variables imputed together no longer share one random quantile per row and come out comonotonic. `QRF` also accepts a `seed` argument.
+Derived seeds stay within the supported uint32 range and are used consistently during numeric and classification tuning and target-specific subsampling.

docs/imputation-benchmarking/cross-validation.md

@@ -41,23 +41,23 @@ Returns a dictionary containing separate results for each metric type:
 ```python
 {
     "quantile_loss": {
-        "results": pd.DataFrame,      # rows: ["train", "test"], cols: quantiles (mean across folds)
+        "results": pd.DataFrame,  # rows: ["train", "test"], cols: quantiles (mean across folds)
         "results_std": pd.DataFrame,  # rows: ["train", "test"], cols: quantiles (std across folds)
         "mean_train": float,
         "mean_test": float,
         "std_train": float,
         "std_test": float,
-        "variables": List[str]        # numerical variables evaluated
+        "variables": List[str],  # numerical variables evaluated
     },
     "log_loss": {
-        "results": pd.DataFrame,      # rows: ["train", "test"], cols: quantiles
+        "results": pd.DataFrame,  # rows: ["train", "test"], cols: quantiles
         "results_std": pd.DataFrame,  # rows: ["train", "test"], cols: quantiles (std across folds)
         "mean_train": float,
         "mean_test": float,
         "std_train": float,
         "std_test": float,
-        "variables": List[str]        # categorical variables evaluated
-    }
+        "variables": List[str],  # categorical variables evaluated
+    },
 }
 ```
 
@@ -77,7 +77,7 @@ results = cross_validate_model(
     data=diabetes_df,
     predictors=["age", "sex", "bmi", "bp"],
     imputed_variables=["s1", "s4"],
-    n_splits=5
+    n_splits=5,
 )
 
 # Check performance for numerical variables

docs/imputation-benchmarking/preprocessing.md

@@ -117,10 +117,10 @@ result = autoimpute(
     predictors=["age", "education"],
     imputed_variables=["income", "wealth"],
     preprocessing={
-        "income": "log",       # Log transform (positive values only)
-        "wealth": "asinh",     # Asinh transform (handles zeros/negatives)
-        "age": "normalize"     # Z-score normalization
-    }
+        "income": "log",  # Log transform (positive values only)
+        "wealth": "asinh",  # Asinh transform (handles zeros/negatives)
+        "age": "normalize",  # Z-score normalization
+    },
 )
 ```
 

docs/imputation-benchmarking/visualizations.md

@@ -83,11 +83,7 @@ comparison_viz = method_comparison_results(
 )
 
 # Generate plot
-fig = comparison_viz.plot(
-    title="Method comparison",
-    show_mean=True,
-    plot_type="bar"
-)
+fig = comparison_viz.plot(title="Method comparison", show_mean=True, plot_type="bar")
 fig.show()
 
 # Get summary statistics
@@ -165,7 +161,7 @@ perf_viz = model_performance_results(
     results=cv_results,
     model_name="QRF",
     method_name="Cross-validation",
-    metric="quantile_loss"
+    metric="quantile_loss",
 )
 
 fig = perf_viz.plot(title="QRF performance")

docs/models/imputer/implement-new-model.md

@@ -73,9 +73,7 @@ class NewModelResults(ImputerResults):
 
         except Exception as e:
             self.logger.error(f"Error during Model prediction: {str(e)}")
-            raise RuntimeError(
-                f"Failed to predict with Model: {str(e)}"
-            ) from e
+            raise RuntimeError(f"Failed to predict with Model: {str(e)}") from e
 ```
 
 ## Implementing the main model class

docs/use_cases/index.md

@@ -23,27 +23,27 @@ Before imputation, make sure both datasets have compatible variables. Identify c
 
 ```python
 # Identify common variables
-common_variables = ['age', 'income', 'education', 'marital_status', 'region']
+common_variables = ["age", "income", "education", "marital_status", "region"]
 
 # Ensure variable formats match (example: education coding)
 education_mapping = {
     1: "less_than_hs",
     2: "high_school",
     3: "some_college",
     4: "bachelor",
-    5: "graduate"
+    5: "graduate",
 }
 
 # Apply standardization to both datasets
 for dataset in [scf_data, cps_data]:
-    dataset['education'] = dataset['education'].map(education_mapping)
+    dataset["education"] = dataset["education"].map(education_mapping)
 
     # Convert income to same units (thousands)
-    if 'income' in dataset.columns:
-        dataset['income'] = dataset['income'] / 1000
+    if "income" in dataset.columns:
+        dataset["income"] = dataset["income"] / 1000
 
 # Identify target variable in donor dataset
-target_variable = ['networth']
+target_variable = ["networth"]
 ```
 
 ## Performing imputation

microimpute/models/qrf.py

@@ -10,7 +10,7 @@
 from quantile_forest import RandomForestQuantileRegressor
 from sklearn.ensemble import RandomForestClassifier
 
-from microimpute.config import VALIDATE_CONFIG
+from microimpute.config import RANDOM_STATE, VALIDATE_CONFIG
 from microimpute.models.imputer import Imputer, ImputerResults
 
 try:
@@ -607,6 +607,7 @@ def __init__(
         batch_size: Optional[int] = None,
         cleanup_interval: int = 10,
         max_train_samples: Optional[int] = None,
+        seed: Optional[int] = RANDOM_STATE,
     ) -> None:
         """Initialize the QRF model.
 
@@ -618,8 +619,11 @@ def __init__(
             max_train_samples: If set, subsample X_train to at most this many
                 rows before fitting. Reduces memory and training time while
                 preserving sequential covariance structure.
+            seed: Base random seed. Each imputed variable is given a distinct
+                seed derived from it, so variables imputed together draw
+                independently. Pass None for non-reproducible draws.
         """
-        super().__init__(log_level=log_level)
+        super().__init__(log_level=log_level, seed=seed)
         self.models = {}
         self.log_level = log_level
         self.memory_efficient = memory_efficient
@@ -695,20 +699,41 @@ def _encode_imputed_variable(
 
         return data
 
+    def _seed_for_variable(self, variable: str) -> Optional[int]:
+        """Derive a distinct seed for one imputed variable.
+
+        Each per-variable model builds its own generator from the seed it is
+        given. Handing every variable the same seed makes them draw the same
+        random quantiles in the same row order, so variables imputed together
+        come out comonotonic regardless of their dependence in the donor. The
+        offset is shared with target-specific subsampling and wraps within
+        sklearn's uint32 seed range.
+        """
+        if self.seed is None:
+            return None
+        if not isinstance(self.seed, (int, np.integer)) or not 0 <= self.seed < 2**32:
+            raise ValueError("seed must be an integer from 0 to 2**32 - 1 or None")
+        try:
+            variable_offset = (self.imputed_variables or []).index(variable)
+        except ValueError:
+            variable_offset = 0
+        return (int(self.seed) + variable_offset) % 2**32
+
     def _create_model_for_variable(self, variable: str, **kwargs) -> Any:
         """Create the appropriate model (classifier or regressor) based on variable type."""
         categorical_targets = getattr(self, "categorical_targets", {})
         boolean_targets = getattr(self, "boolean_targets", {})
+        seed = self._seed_for_variable(variable)
 
         if variable in categorical_targets:
             # Use classifier for categorical targets
-            return _RandomForestClassifierModel(seed=self.seed, logger=self.logger)
+            return _RandomForestClassifierModel(seed=seed, logger=self.logger)
         elif variable in boolean_targets:
             # Use classifier for boolean targets
-            return _RandomForestClassifierModel(seed=self.seed, logger=self.logger)
+            return _RandomForestClassifierModel(seed=seed, logger=self.logger)
         else:
             # Use QRF for numeric targets
-            return _QRFModel(seed=self.seed, logger=self.logger)
+            return _QRFModel(seed=seed, logger=self.logger)
 
     def _fit_model(
         self,
@@ -800,12 +825,7 @@ def _target_fit_data(
             self.max_train_samples is not None
             and len(target_train) > self.max_train_samples
         ):
-            try:
-                variable_offset = (self.imputed_variables or []).index(variable)
-            except ValueError:
-                variable_offset = 0
-            seed = None if self.seed is None else self.seed + variable_offset
-            rng = np.random.default_rng(seed)
+            rng = np.random.default_rng(self._seed_for_variable(variable))
             sel = rng.choice(
                 len(target_train), size=self.max_train_samples, replace=False
             )
@@ -1465,7 +1485,9 @@ def objective(trial: optuna.Trial) -> float:
                         y_val = X_val_fold[var]
 
                         # Create and fit QRF model with trial parameters
-                        model = _QRFModel(seed=self.seed, logger=self.logger)
+                        model = _QRFModel(
+                            seed=self._seed_for_variable(var), logger=self.logger
+                        )
                         model.fit(
                             X_train_augmented[encoded_predictors],
                             X_train_fold[var],
@@ -1623,7 +1645,7 @@ def objective(trial: optuna.Trial) -> float:
 
                         # Create and fit RFC model with trial parameters
                         model = _RandomForestClassifierModel(
-                            seed=self.seed, logger=self.logger
+                            seed=self._seed_for_variable(var), logger=self.logger
                         )
 
                         # Determine variable type and fit appropriately

tests/test_models/test_qrf.py

@@ -1649,3 +1649,159 @@ def test_qrf_fit_predict_with_max_train_samples() -> None:
 
     assert result.shape == (n_test, 1)
     assert not result.isna().any().any()
+
+
+def test_per_variable_models_draw_independently() -> None:
+    """Variables imputed together must not share one random quantile per row.
+
+    Every per-variable model builds its own generator from the seed it is
+    given. When they all receive the same seed they draw the same quantiles in
+    the same row order, so the imputed variables come out comonotonic whatever
+    their dependence in the donor.
+    """
+    rng = np.random.default_rng(7)
+    n = 1500
+    predictors = pd.DataFrame(
+        {"inc": rng.normal(30, 8, n), "age": rng.normal(45, 12, n)}
+    )
+    targets = ["savings", "property_wealth", "corporate_wealth"]
+    train = predictors.copy()
+    for variable in targets:
+        # Shared signal, independent shocks: the conditional dependence is nil.
+        train[variable] = 0.5 * predictors["inc"] + rng.normal(0, 10, n)
+
+    test = pd.DataFrame({"inc": rng.normal(30, 8, 600), "age": rng.normal(45, 12, 600)})
+
+    model = QRF(log_level="WARNING")
+    imputations = model.fit(train, ["inc", "age"], targets).predict(test)
+
+    ranks = imputations[targets].rank().corr()
+    off_diagonal = [
+        abs(ranks.loc[a, b]) for i, a in enumerate(targets) for b in targets[i + 1 :]
+    ]
+    assert max(off_diagonal) < 0.4, (
+        "imputed variables are comonotonic; per-variable models are sharing "
+        f"a seed (rank correlations {off_diagonal})"
+    )
+
+
+def test_seed_is_configurable_and_reproducible() -> None:
+    """QRF should accept a seed, and the same seed should reproduce draws."""
+    rng = np.random.default_rng(3)
+    n = 400
+    train = pd.DataFrame({"x": rng.normal(size=n)})
+    train["y"] = train["x"] + rng.normal(0, 1, n)
+    test = pd.DataFrame({"x": rng.normal(size=120)})
+
+    first = QRF(log_level="WARNING", seed=1234).fit(train, ["x"], ["y"]).predict(test)
+    second = QRF(log_level="WARNING", seed=1234).fit(train, ["x"], ["y"]).predict(test)
+
+    np.testing.assert_allclose(first["y"], second["y"])
+
+
+@pytest.mark.parametrize("target_type", ["numeric", "boolean", "categorical"])
+def test_qrf_max_seed_supports_multiple_targets(target_type: str) -> None:
+    """Every valid sklearn seed must support reproducible multi-target fits."""
+    rng = np.random.default_rng(71)
+    data = pd.DataFrame(
+        {name: rng.normal(size=120) for name in ["x", "first", "second"]}
+    )
+    if target_type == "boolean":
+        data[["first", "second"]] = data[["first", "second"]] > 0
+    elif target_type == "categorical":
+        for name in ["first", "second"]:
+            data[name] = np.where(data[name] > 0, "yes", "no")
+
+    predictions = []
+    for _ in range(2):
+        fitted = QRF(seed=2**32 - 1).fit(
+            data, ["x"], ["first", "second"], n_estimators=12
+        )
+        seeds = [model.seed for model in fitted.models.values()]
+        assert len(set(seeds)) == 2
+        assert all(0 <= seed < 2**32 for seed in seeds)
+        predictions.append(fitted.predict(data[["x"]].iloc[:20]))
+    pd.testing.assert_frame_equal(*predictions)
+    assert not predictions[0].isna().any().any()
+
+
+@pytest.mark.parametrize("seed", [None, 0, 42])
+def test_qrf_child_seeds_preserve_existing_seed_values(seed) -> None:
+    """Ordinary seeds and entropy-based draws retain their established meaning."""
+    model = QRF(seed=seed)
+    model.imputed_variables = ["first", "second"]
+    expected = [None, None] if seed is None else [seed, seed + 1]
+    assert [
+        model._seed_for_variable(name) for name in model.imputed_variables
+    ] == expected
+
+
+@pytest.mark.parametrize("seed", [-1, 2**32, 1.5])
+def test_qrf_invalid_base_seed_is_not_normalized(seed) -> None:
+    """Wrapping child seeds must not silently accept invalid sklearn base seeds."""
+    rng = np.random.default_rng(3)
+    data = pd.DataFrame({"x": rng.normal(size=30), "y": rng.normal(size=30)})
+    with pytest.raises(RuntimeError):
+        QRF(seed=seed).fit(data, ["x"], ["y"], n_estimators=5)
+
+
+@pytest.mark.parametrize("target_type", ["numeric", "boolean"])
+def test_qrf_tuning_uses_distinct_target_seeds(monkeypatch, target_type: str) -> None:
+    """Real Optuna folds must fit the same independent streams as final models."""
+    from microimpute.models.qrf import _RandomForestClassifierModel
+
+    rng = np.random.default_rng(73)
+    data = pd.DataFrame(
+        {name: rng.normal(size=100) for name in ["x", "first", "second"]}
+    )
+    model = QRF(seed=123)
+    model.imputed_variables = ["first", "second"]
+    if target_type == "numeric":
+        internal_model = _QRFModel
+        tune = model._tune_qrf_hyperparameters
+    else:
+        data[["first", "second"]] = data[["first", "second"]] > 0
+        model.boolean_targets = {"first": {}, "second": {}}
+        internal_model = _RandomForestClassifierModel
+        tune = model._tune_rfc_hyperparameters
+
+    fitted_seeds = []
+    original_fit = internal_model.fit
+
+    def record_fit(self, X, y, **kwargs):
+        fitted_seeds.append((y.name, self.seed))
+        return original_fit(self, X, y, **kwargs)
+
+    monkeypatch.setattr(internal_model, "fit", record_fit)
+    tune(data, ["x"], ["first", "second"], n_cv_folds=2, n_trials=1)
+    assert fitted_seeds == [("first", 123), ("second", 124)] * 2
+
+
+def test_qrf_target_subsampling_uses_bounded_child_seed(monkeypatch) -> None:
+    """Filtered training rows and their model use one consistent child seed."""
+    rng = np.random.default_rng(9)
+    data = pd.DataFrame(
+        {name: rng.normal(size=120) for name in ["x", "first", "second"]}
+    )
+    fitted_indices = {}
+    original_fit = _QRFModel.fit
+
+    def record_fit(self, X, y, **kwargs):
+        fitted_indices[y.name] = X.index.to_numpy()
+        return original_fit(self, X, y, **kwargs)
+
+    monkeypatch.setattr(_QRFModel, "fit", record_fit)
+    QRF(seed=2**32 - 1, max_train_samples=50).fit(
+        data,
+        ["x"],
+        ["first", "second"],
+        target_filters={
+            name: np.ones(len(data), dtype=bool) for name in ["first", "second"]
+        },
+        n_estimators=12,
+    )
+    # The second stream wraps to zero at the uint32 boundary.
+    expected_indices = np.random.default_rng(0).choice(
+        len(data), size=50, replace=False
+    )
+    np.testing.assert_array_equal(fitted_indices["second"], expected_indices)