|
| 1 | +import unittest |
| 2 | +import numpy as np |
| 3 | +from sklearn.linear_model import LogisticRegression |
| 4 | +from dte_adj import ( |
| 5 | + SimpleStratifiedDistributionEstimator, |
| 6 | + AdjustedStratifiedDistributionEstimator, |
| 7 | +) |
| 8 | + |
| 9 | + |
| 10 | +def generate_data(n=1000, S=4, d=2, discrete=False): |
| 11 | + d = 20 |
| 12 | + |
| 13 | + Z = np.random.uniform(0, 1, n) |
| 14 | + |
| 15 | + S_i = np.digitize(Z, np.linspace(0, 1, S + 1)[1:-1]) |
| 16 | + |
| 17 | + X = np.random.multivariate_normal(mean=np.zeros(d), cov=np.eye(d), size=n) |
| 18 | + |
| 19 | + W = np.zeros(n, dtype=int) |
| 20 | + unique_strata = np.unique(S_i) |
| 21 | + for s in unique_strata: |
| 22 | + idx = np.where(S_i == s)[0] |
| 23 | + n_s = len(idx) |
| 24 | + W[idx[: n_s // 2]] = 1 |
| 25 | + np.random.shuffle(W[idx]) |
| 26 | + |
| 27 | + b_X = ( |
| 28 | + np.sin(np.pi * X[:, 0] * X[:, 1]) |
| 29 | + + 2 * (X[:, 2] - 0.5) ** 2 |
| 30 | + + X[:, 3] |
| 31 | + + 0.5 * X[:, 4] |
| 32 | + ) |
| 33 | + c_X = 0.1 * (X[:, 0] + np.log(1 + np.exp(X[:, 1]))) |
| 34 | + |
| 35 | + gamma = 0.1 |
| 36 | + u = np.random.normal(0, 1, n) |
| 37 | + |
| 38 | + Y = b_X + c_X * W + gamma * Z + u |
| 39 | + if discrete: |
| 40 | + Y = np.random.poisson(0.2 * np.abs(Y)) |
| 41 | + |
| 42 | + return {"W": W, "X": X, "Z": Z, "Y": Y, "strata": S_i} |
| 43 | + |
| 44 | + |
| 45 | +class TestStratifiedEstimators(unittest.TestCase): |
| 46 | + def setUp(self): |
| 47 | + np.random.seed(42) |
| 48 | + data = generate_data(n=1000, S=4, d=20, discrete=False) |
| 49 | + self.X = data["X"] |
| 50 | + self.W = data["W"] |
| 51 | + self.Y = data["Y"] |
| 52 | + self.strata = data["strata"] |
| 53 | + self.locations = np.linspace(self.Y.min(), self.Y.max(), 20) |
| 54 | + |
| 55 | + def test_simple_stratified_estimator_fit(self): |
| 56 | + estimator = SimpleStratifiedDistributionEstimator() |
| 57 | + result = estimator.fit(self.X, self.W, self.Y, self.strata) |
| 58 | + |
| 59 | + self.assertIsInstance(result, SimpleStratifiedDistributionEstimator) |
| 60 | + self.assertTrue(np.array_equal(estimator.covariates, self.X)) |
| 61 | + self.assertTrue(np.array_equal(estimator.treatment_arms, self.W)) |
| 62 | + self.assertTrue(np.array_equal(estimator.outcomes, self.Y)) |
| 63 | + self.assertTrue(np.array_equal(estimator.strata, self.strata)) |
| 64 | + |
| 65 | + def test_simple_stratified_estimator_predict_dte(self): |
| 66 | + estimator = SimpleStratifiedDistributionEstimator() |
| 67 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 68 | + |
| 69 | + dte, lower_bound, upper_bound = estimator.predict_dte( |
| 70 | + target_treatment_arm=1, |
| 71 | + control_treatment_arm=0, |
| 72 | + locations=self.locations, |
| 73 | + alpha=0.05, |
| 74 | + ) |
| 75 | + |
| 76 | + self.assertEqual(dte.shape, self.locations.shape) |
| 77 | + self.assertEqual(lower_bound.shape, self.locations.shape) |
| 78 | + self.assertEqual(upper_bound.shape, self.locations.shape) |
| 79 | + self.assertTrue(np.all(lower_bound <= dte)) |
| 80 | + self.assertTrue(np.all(dte <= upper_bound)) |
| 81 | + |
| 82 | + def test_simple_stratified_estimator_predict_pte(self): |
| 83 | + estimator = SimpleStratifiedDistributionEstimator() |
| 84 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 85 | + |
| 86 | + pte, lower_bound, upper_bound = estimator.predict_pte( |
| 87 | + target_treatment_arm=1, |
| 88 | + control_treatment_arm=0, |
| 89 | + locations=self.locations, |
| 90 | + alpha=0.05, |
| 91 | + ) |
| 92 | + |
| 93 | + expected_length = len(self.locations) - 1 |
| 94 | + self.assertEqual(pte.shape, (expected_length,)) |
| 95 | + self.assertEqual(lower_bound.shape, (expected_length,)) |
| 96 | + self.assertEqual(upper_bound.shape, (expected_length,)) |
| 97 | + self.assertTrue(np.all(lower_bound <= upper_bound)) |
| 98 | + |
| 99 | + def test_simple_stratified_estimator_predict_qte(self): |
| 100 | + estimator = SimpleStratifiedDistributionEstimator() |
| 101 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 102 | + |
| 103 | + quantiles = np.array([0.25, 0.5, 0.75]) |
| 104 | + qte, lower_bound, upper_bound = estimator.predict_qte( |
| 105 | + target_treatment_arm=1, |
| 106 | + control_treatment_arm=0, |
| 107 | + quantiles=quantiles, |
| 108 | + n_bootstrap=50, |
| 109 | + ) |
| 110 | + |
| 111 | + self.assertEqual(qte.shape, quantiles.shape) |
| 112 | + self.assertEqual(lower_bound.shape, quantiles.shape) |
| 113 | + self.assertEqual(upper_bound.shape, quantiles.shape) |
| 114 | + self.assertTrue(np.all(lower_bound <= upper_bound)) |
| 115 | + |
| 116 | + def test_adjusted_stratified_estimator_fit(self): |
| 117 | + base_model = LogisticRegression(max_iter=1000, random_state=42) |
| 118 | + estimator = AdjustedStratifiedDistributionEstimator(base_model, folds=3) |
| 119 | + result = estimator.fit(self.X, self.W, self.Y, self.strata) |
| 120 | + |
| 121 | + self.assertIsInstance(result, AdjustedStratifiedDistributionEstimator) |
| 122 | + self.assertTrue(np.array_equal(estimator.covariates, self.X)) |
| 123 | + self.assertTrue(np.array_equal(estimator.treatment_arms, self.W)) |
| 124 | + self.assertTrue(np.array_equal(estimator.outcomes, self.Y)) |
| 125 | + self.assertTrue(np.array_equal(estimator.strata, self.strata)) |
| 126 | + self.assertEqual(estimator.folds, 3) |
| 127 | + |
| 128 | + def test_adjusted_stratified_estimator_predict_dte(self): |
| 129 | + base_model = LogisticRegression(max_iter=1000, random_state=42) |
| 130 | + estimator = AdjustedStratifiedDistributionEstimator(base_model, folds=3) |
| 131 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 132 | + |
| 133 | + dte, lower_bound, upper_bound = estimator.predict_dte( |
| 134 | + target_treatment_arm=1, |
| 135 | + control_treatment_arm=0, |
| 136 | + locations=self.locations, |
| 137 | + alpha=0.05, |
| 138 | + variance_type="moment", |
| 139 | + ) |
| 140 | + |
| 141 | + self.assertEqual(dte.shape, self.locations.shape) |
| 142 | + self.assertEqual(lower_bound.shape, self.locations.shape) |
| 143 | + self.assertEqual(upper_bound.shape, self.locations.shape) |
| 144 | + self.assertTrue(np.all(lower_bound <= dte)) |
| 145 | + self.assertTrue(np.all(dte <= upper_bound)) |
| 146 | + |
| 147 | + def test_adjusted_stratified_estimator_predict_pte(self): |
| 148 | + base_model = LogisticRegression(max_iter=1000, random_state=42) |
| 149 | + estimator = AdjustedStratifiedDistributionEstimator(base_model, folds=3) |
| 150 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 151 | + |
| 152 | + pte, lower_bound, upper_bound = estimator.predict_pte( |
| 153 | + target_treatment_arm=1, |
| 154 | + control_treatment_arm=0, |
| 155 | + locations=self.locations, |
| 156 | + alpha=0.05, |
| 157 | + variance_type="moment", |
| 158 | + ) |
| 159 | + |
| 160 | + expected_length = len(self.locations) - 1 |
| 161 | + self.assertEqual(pte.shape, (expected_length,)) |
| 162 | + self.assertEqual(lower_bound.shape, (expected_length,)) |
| 163 | + self.assertEqual(upper_bound.shape, (expected_length,)) |
| 164 | + self.assertTrue(np.all(lower_bound <= upper_bound)) |
| 165 | + |
| 166 | + def test_adjusted_stratified_estimator_predict_qte(self): |
| 167 | + base_model = LogisticRegression(max_iter=1000, random_state=42) |
| 168 | + estimator = AdjustedStratifiedDistributionEstimator(base_model, folds=3) |
| 169 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 170 | + |
| 171 | + quantiles = np.array([0.25, 0.5, 0.75]) |
| 172 | + qte, lower_bound, upper_bound = estimator.predict_qte( |
| 173 | + target_treatment_arm=1, |
| 174 | + control_treatment_arm=0, |
| 175 | + quantiles=quantiles, |
| 176 | + n_bootstrap=50, |
| 177 | + ) |
| 178 | + |
| 179 | + self.assertEqual(qte.shape, quantiles.shape) |
| 180 | + self.assertEqual(lower_bound.shape, quantiles.shape) |
| 181 | + self.assertEqual(upper_bound.shape, quantiles.shape) |
| 182 | + self.assertTrue(np.all(lower_bound <= upper_bound)) |
| 183 | + |
| 184 | + def test_discrete_outcomes(self): |
| 185 | + data = generate_data(n=1000, S=4, d=20, discrete=True) |
| 186 | + |
| 187 | + estimator = SimpleStratifiedDistributionEstimator() |
| 188 | + estimator.fit(data["X"], data["W"], data["Y"], data["strata"]) |
| 189 | + |
| 190 | + locations = np.arange(0, data["Y"].max() + 1) |
| 191 | + dte, lower, upper = estimator.predict_dte(1, 0, locations) |
| 192 | + |
| 193 | + self.assertEqual(dte.shape, locations.shape) |
| 194 | + self.assertTrue(np.all(lower <= upper)) |
| 195 | + |
| 196 | + def test_invalid_input_shapes(self): |
| 197 | + estimator = SimpleStratifiedDistributionEstimator() |
| 198 | + |
| 199 | + X_wrong = self.X[:-10] |
| 200 | + |
| 201 | + with self.assertRaises(ValueError): |
| 202 | + estimator.fit(X_wrong, self.W, self.Y, self.strata) |
| 203 | + |
| 204 | + def test_different_alpha_values(self): |
| 205 | + estimator = SimpleStratifiedDistributionEstimator() |
| 206 | + estimator.fit(self.X, self.W, self.Y, self.strata) |
| 207 | + |
| 208 | + locations = self.locations[:10] |
| 209 | + |
| 210 | + _, lower_005, upper_005 = estimator.predict_dte(1, 0, locations, alpha=0.05) |
| 211 | + _, lower_010, upper_010 = estimator.predict_dte(1, 0, locations, alpha=0.10) |
| 212 | + |
| 213 | + width_005 = upper_005 - lower_005 |
| 214 | + width_010 = upper_010 - lower_010 |
| 215 | + |
| 216 | + self.assertTrue(np.all(width_010 < width_005)) |
0 commit comments