Skip to content

Covariate Adaptive Randomization Estimators

This page documents estimators that work with stratified experimental designs, particularly for covariate-adaptive randomization (CAR) within strata.

These estimators are designed to handle stratified block randomization where participants are grouped into strata based on baseline covariates before treatment assignment. The key methodological contribution is leveraging additional covariates beyond strata indicators using machine learning methods to enhance the precision of distributional treatment effect estimates.

Byambadalai et al. (2025)1 propose a flexible distribution regression framework that achieves the semiparametric efficiency bound for distributional treatment effects under CAR, demonstrating that regression-adjusted estimators can optimally utilize covariate information in stratified designs.

SimpleStratifiedDistributionEstimator

Bases: DistributionEstimatorBase

A class is for estimating the empirical distribution function and computing the Distributional parameters for CAR.

Initializes the DistributionFunctionMixin.

Returns: DistributionFunctionMixin: An instance of the estimator.

Source code in dte_adj/base.py
def __init__(self):
    """
    Initializes the DistributionFunctionMixin.

    Returns:
        DistributionFunctionMixin: An instance of the estimator.
    """
    self.covariates = None
    self.outcomes = None
    self.treatment_arms = None
    self.last_locations = None

predict_dte

predict_dte(
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]

Compute Distribution Treatment Effects (DTE) based on the estimator for the distribution function.

The DTE measures the difference in cumulative distribution functions between treatment groups at specified locations. It quantifies how treatment affects the probability of observing outcomes below each threshold.

Args: target_treatment_arm (int): The index of the treatment arm of the treatment group. control_treatment_arm (int): The index of the treatment arm of the control group. locations (np.ndarray, optional): Scalar values to be used for computing the cumulative distribution. If None, evenly-spaced locations spanning the observed outcome range are generated automatically. The number of points is determined from data size and distribution via np.histogram_bin_edges(outcomes, bins='auto'). The actual array used is stored on self.last_locations. alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05. variance_type (str, optional): Variance type to be used to compute confidence intervals. Available values are "moment", "simple", and "uniform". Defaults to "moment". n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing: - Expected DTEs (np.ndarray): Treatment effect estimates at each location - Lower bounds (np.ndarray): Lower confidence interval bounds - Upper bounds (np.ndarray): Upper confidence interval bounds

Example:

import numpy as np
from dte_adj import SimpleDistributionEstimator

# Generate sample data
X = np.random.randn(1000, 5)
D = np.random.binomial(1, 0.5, 1000)
Y = X[:, 0] + 2 * D + np.random.randn(1000)

# Fit estimator
estimator = SimpleDistributionEstimator()
estimator.fit(X, D, Y)

# Compute DTE
locations = np.linspace(Y.min(), Y.max(), 20)
dte, lower, upper = estimator.predict_dte(
    target_treatment_arm=1,
    control_treatment_arm=0,
    locations=locations,
    variance_type="moment"
)

print(f"DTE shape: {dte.shape}")  # Should match locations.shape
print(f"Average DTE: {dte.mean():.3f}")

Source code in dte_adj/base.py
def predict_dte(
    self,
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[np.ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
    """
    Compute Distribution Treatment Effects (DTE) based on the estimator for the distribution function.

    The DTE measures the difference in cumulative distribution functions between treatment groups
    at specified locations. It quantifies how treatment affects the probability of observing
    outcomes below each threshold.

    Args:
        target_treatment_arm (int): The index of the treatment arm of the treatment group.
        control_treatment_arm (int): The index of the treatment arm of the control group.
        locations (np.ndarray, optional): Scalar values to be used for computing the cumulative
            distribution. If None, evenly-spaced locations spanning the observed outcome range
            are generated automatically. The number of points is determined from data size and
            distribution via ``np.histogram_bin_edges(outcomes, bins='auto')``. The actual array
            used is stored on ``self.last_locations``.
        alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05.
        variance_type (str, optional): Variance type to be used to compute confidence intervals.
            Available values are "moment", "simple", and "uniform". Defaults to "moment".
        n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing:
            - Expected DTEs (np.ndarray): Treatment effect estimates at each location
            - Lower bounds (np.ndarray): Lower confidence interval bounds
            - Upper bounds (np.ndarray): Upper confidence interval bounds

    Example:
        ```python
        import numpy as np
        from dte_adj import SimpleDistributionEstimator

        # Generate sample data
        X = np.random.randn(1000, 5)
        D = np.random.binomial(1, 0.5, 1000)
        Y = X[:, 0] + 2 * D + np.random.randn(1000)

        # Fit estimator
        estimator = SimpleDistributionEstimator()
        estimator.fit(X, D, Y)

        # Compute DTE
        locations = np.linspace(Y.min(), Y.max(), 20)
        dte, lower, upper = estimator.predict_dte(
            target_treatment_arm=1,
            control_treatment_arm=0,
            locations=locations,
            variance_type="moment"
        )

        print(f"DTE shape: {dte.shape}")  # Should match locations.shape
        print(f"Average DTE: {dte.mean():.3f}")
        ```
    """
    if locations is None:
        locations = _infer_default_locations(self.outcomes, for_intervals=False)
    self.last_locations = locations
    return self._compute_dtes(
        target_treatment_arm,
        control_treatment_arm,
        locations,
        alpha,
        variance_type,
        n_bootstrap,
        display_progress,
    )

predict_pte

predict_pte(
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]

Compute Probability Treatment Effects (PTE) based on the estimator for the distribution function.

The PTE measures the difference in probability mass between treatment groups for intervals defined by consecutive location pairs. It quantifies how treatment affects the probability of observing outcomes within specific ranges.

Args: target_treatment_arm (int): The index of the treatment arm of the treatment group. control_treatment_arm (int): The index of the treatment arm of the control group. locations (np.ndarray, optional): Scalar values defining interval boundaries for probability computation. For each interval (locations[i], locations[i+1]], the PTE is computed. If None, boundaries spanning the observed outcome range are generated automatically with the left endpoint placed just below outcomes.min() so that minimum-valued samples fall inside the first interval. The number of boundaries is determined from data size and distribution via np.histogram_bin_edges(outcomes, bins='auto'). The actual array used is stored on self.last_locations. alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05. variance_type (str, optional): Variance type to be used to compute confidence intervals. Available values are "moment", "simple", and "uniform". Defaults to "moment". n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing: - Expected PTEs (np.ndarray): Treatment effect estimates for each interval, shape (len(locations)-1,) - Lower bounds (np.ndarray): Lower confidence interval bounds - Upper bounds (np.ndarray): Upper confidence interval bounds

Example:

import numpy as np
from dte_adj import SimpleDistributionEstimator

# Generate sample data
X = np.random.randn(1000, 5)
D = np.random.binomial(1, 0.5, 1000)
Y = X[:, 0] + 2 * D + np.random.randn(1000)

# Fit estimator
estimator = SimpleDistributionEstimator()
estimator.fit(X, D, Y)

# Define interval boundaries
locations = np.array([-2, -1, 0, 1, 2])  # Creates intervals: (-2,-1], (-1,0], (0,1], (1,2]

# Compute PTE
pte, lower, upper = estimator.predict_pte(
    target_treatment_arm=1,
    control_treatment_arm=0,
    locations=locations,
    variance_type="moment"
)

print(f"PTE shape: {pte.shape}")  # Should be (4,) for 4 intervals
print(f"Interval effects: {pte}")

Source code in dte_adj/base.py
def predict_pte(
    self,
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[np.ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
    """
    Compute Probability Treatment Effects (PTE) based on the estimator for the distribution function.

    The PTE measures the difference in probability mass between treatment groups for intervals
    defined by consecutive location pairs. It quantifies how treatment affects the probability
    of observing outcomes within specific ranges.

    Args:
        target_treatment_arm (int): The index of the treatment arm of the treatment group.
        control_treatment_arm (int): The index of the treatment arm of the control group.
        locations (np.ndarray, optional): Scalar values defining interval boundaries for
            probability computation. For each interval (locations[i], locations[i+1]], the PTE
            is computed. If None, boundaries spanning the observed outcome range are generated
            automatically with the left endpoint placed just below ``outcomes.min()`` so that
            minimum-valued samples fall inside the first interval. The number of boundaries is
            determined from data size and distribution via
            ``np.histogram_bin_edges(outcomes, bins='auto')``. The actual array used is stored
            on ``self.last_locations``.
        alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05.
        variance_type (str, optional): Variance type to be used to compute confidence intervals.
            Available values are "moment", "simple", and "uniform". Defaults to "moment".
        n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing:
            - Expected PTEs (np.ndarray): Treatment effect estimates for each interval,
              shape (len(locations)-1,)
            - Lower bounds (np.ndarray): Lower confidence interval bounds
            - Upper bounds (np.ndarray): Upper confidence interval bounds

    Example:
        ```python
        import numpy as np
        from dte_adj import SimpleDistributionEstimator

        # Generate sample data
        X = np.random.randn(1000, 5)
        D = np.random.binomial(1, 0.5, 1000)
        Y = X[:, 0] + 2 * D + np.random.randn(1000)

        # Fit estimator
        estimator = SimpleDistributionEstimator()
        estimator.fit(X, D, Y)

        # Define interval boundaries
        locations = np.array([-2, -1, 0, 1, 2])  # Creates intervals: (-2,-1], (-1,0], (0,1], (1,2]

        # Compute PTE
        pte, lower, upper = estimator.predict_pte(
            target_treatment_arm=1,
            control_treatment_arm=0,
            locations=locations,
            variance_type="moment"
        )

        print(f"PTE shape: {pte.shape}")  # Should be (4,) for 4 intervals
        print(f"Interval effects: {pte}")
        ```
    """
    if locations is None:
        locations = _infer_default_locations(self.outcomes, for_intervals=True)
    self.last_locations = locations
    return self._compute_ptes(
        target_treatment_arm,
        control_treatment_arm,
        locations,
        alpha,
        variance_type,
        n_bootstrap,
        display_progress,
    )

predict_qte

predict_qte(
    target_treatment_arm: int,
    control_treatment_arm: int,
    quantiles: Optional[ndarray] = None,
    alpha: float = 0.05,
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]

Compute Quantile Treatment Effects (QTE) based on the estimator for the distribution function.

The QTE measures the difference in quantiles between treatment groups, providing insights into how treatment affects different parts of the outcome distribution. For stratified estimators, the computation properly accounts for strata.

Variance is estimated by stratified bootstrap: indices are resampled with replacement within each stratum independently, which preserves per-stratum sample sizes and reflects the covariate-adaptive randomization (CAR) design. For estimators without strata (single stratum), this degenerates to a plain bootstrap.

Args: target_treatment_arm (int): The index of the treatment arm of the treatment group. control_treatment_arm (int): The index of the treatment arm of the control group. quantiles (np.ndarray, optional): Quantiles used for QTE. Defaults to [0.1, 0.2, ..., 0.9]. alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05. n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing: - Expected QTEs (np.ndarray): Treatment effect estimates at each quantile - Lower bounds (np.ndarray): Lower confidence interval bounds - Upper bounds (np.ndarray): Upper confidence interval bounds

Example:

import numpy as np
from dte_adj import SimpleStratifiedDistributionEstimator

# Generate stratified sample data
X = np.random.randn(1000, 5)
strata = np.random.choice([0, 1, 2], size=1000)
D = np.random.binomial(1, 0.5, 1000)
Y = X[:, 0] + 2 * D + 0.5 * strata + np.random.randn(1000)

# Fit stratified estimator
estimator = SimpleStratifiedDistributionEstimator()
estimator.fit(X, D, Y, strata)

# Compute QTE at specific quantiles
quantiles = np.array([0.25, 0.5, 0.75])  # 25th, 50th, 75th percentiles
qte, lower, upper = estimator.predict_qte(
    target_treatment_arm=1,
    control_treatment_arm=0,
    quantiles=quantiles,
    n_bootstrap=100
)

print(f"QTE at quantiles {quantiles}: {qte}")
print(f"Median effect (50th percentile): {qte[1]:.3f}")

Source code in dte_adj/base.py
def predict_qte(
    self,
    target_treatment_arm: int,
    control_treatment_arm: int,
    quantiles: Optional[np.ndarray] = None,
    alpha: float = 0.05,
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
    """
    Compute Quantile Treatment Effects (QTE) based on the estimator for the distribution function.

    The QTE measures the difference in quantiles between treatment groups, providing insights
    into how treatment affects different parts of the outcome distribution. For stratified
    estimators, the computation properly accounts for strata.

    Variance is estimated by stratified bootstrap: indices are resampled with replacement
    within each stratum independently, which preserves per-stratum sample sizes and reflects
    the covariate-adaptive randomization (CAR) design. For estimators without strata
    (single stratum), this degenerates to a plain bootstrap.

    Args:
        target_treatment_arm (int): The index of the treatment arm of the treatment group.
        control_treatment_arm (int): The index of the treatment arm of the control group.
        quantiles (np.ndarray, optional): Quantiles used for QTE. Defaults to [0.1, 0.2, ..., 0.9].
        alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05.
        n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing:
            - Expected QTEs (np.ndarray): Treatment effect estimates at each quantile
            - Lower bounds (np.ndarray): Lower confidence interval bounds
            - Upper bounds (np.ndarray): Upper confidence interval bounds

    Example:
        ```python
        import numpy as np
        from dte_adj import SimpleStratifiedDistributionEstimator

        # Generate stratified sample data
        X = np.random.randn(1000, 5)
        strata = np.random.choice([0, 1, 2], size=1000)
        D = np.random.binomial(1, 0.5, 1000)
        Y = X[:, 0] + 2 * D + 0.5 * strata + np.random.randn(1000)

        # Fit stratified estimator
        estimator = SimpleStratifiedDistributionEstimator()
        estimator.fit(X, D, Y, strata)

        # Compute QTE at specific quantiles
        quantiles = np.array([0.25, 0.5, 0.75])  # 25th, 50th, 75th percentiles
        qte, lower, upper = estimator.predict_qte(
            target_treatment_arm=1,
            control_treatment_arm=0,
            quantiles=quantiles,
            n_bootstrap=100
        )

        print(f"QTE at quantiles {quantiles}: {qte}")
        print(f"Median effect (50th percentile): {qte[1]:.3f}")
        ```
    """
    if quantiles is None:
        quantiles = np.arange(1, 10) / 10
    if np.any((quantiles <= 0) | (quantiles >= 1)):
        raise ValueError("quantiles must be in the open interval (0, 1)")

    qte = self._compute_qtes(
        target_treatment_arm,
        control_treatment_arm,
        quantiles,
        self.covariates,
        self.treatment_arms,
        self.outcomes,
        self.strata,
    )

    # Precompute stratum indices for stratified bootstrap.
    # When there is a single stratum this is equivalent to plain bootstrap.
    unique_strata = np.unique(self.strata)
    strata_indices = [np.where(self.strata == s)[0] for s in unique_strata]

    qtes = np.zeros((n_bootstrap, qte.shape[0]))
    bootstrap_iter = range(n_bootstrap)
    if display_progress:
        bootstrap_iter = tqdm(bootstrap_iter, desc="Bootstrap QTE")
    for b in bootstrap_iter:
        bootstrap_indexes = np.concatenate(
            [
                np.random.choice(idx, size=len(idx), replace=True)
                for idx in strata_indices
            ]
        )

        qtes[b] = self._compute_qtes(
            target_treatment_arm,
            control_treatment_arm,
            quantiles,
            self.covariates[bootstrap_indexes],
            self.treatment_arms[bootstrap_indexes],
            self.outcomes[bootstrap_indexes],
            self.strata[bootstrap_indexes],
        )

    qte_var = qtes.var(axis=0)

    qte_lower = qte + norm.ppf(alpha / 2) * np.sqrt(qte_var)
    qte_upper = qte + norm.ppf(1 - alpha / 2) * np.sqrt(qte_var)

    return qte, qte_lower, qte_upper

predict

predict(
    treatment_arm: int,
    locations: ndarray,
    display_progress: bool = True,
) -> np.ndarray

Compute cumulative distribution values.

Args: treatment_arm (int): The index of the treatment arm. locations (np.ndarray): Scalar values to be used for computing the cumulative distribution. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: np.ndarray: Estimated cumulative distribution values for the input.

Source code in dte_adj/base.py
def predict(
    self, treatment_arm: int, locations: np.ndarray, display_progress: bool = True
) -> np.ndarray:
    """
    Compute cumulative distribution values.

    Args:
        treatment_arm (int): The index of the treatment arm.
        locations (np.ndarray): Scalar values to be used for computing the cumulative distribution.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        np.ndarray: Estimated cumulative distribution values for the input.
    """
    if self.outcomes is None:
        raise ValueError(
            "This estimator has not been trained yet. Please call fit first"
        )

    if treatment_arm not in self.treatment_arms:
        raise ValueError(
            f"This target treatment arm was not included in the training data: {treatment_arm}"
        )

    return self._compute_cumulative_distribution(
        treatment_arm,
        locations,
        self.covariates,
        self.treatment_arms,
        self.outcomes,
        display_progress=display_progress,
    )[0]

fit

fit(
    covariates: ArrayLike,
    treatment_arms: ArrayLike,
    outcomes: ArrayLike,
    strata: ArrayLike,
) -> DistributionEstimatorBase

Train the DistributionEstimatorBase.

Args: covariates: Pre-treatment covariates. treatment_arms: The index of the treatment arm. outcomes: Scalar-valued observed outcome. strata: Stratum indicators.

Returns: DistributionEstimatorBase: The fitted estimator.

Source code in dte_adj/stratified.py
def fit(
    self,
    covariates: ArrayLike,
    treatment_arms: ArrayLike,
    outcomes: ArrayLike,
    strata: ArrayLike,
) -> DistributionEstimatorBase:
    """
    Train the DistributionEstimatorBase.

    Args:
        covariates: Pre-treatment covariates.
        treatment_arms: The index of the treatment arm.
        outcomes: Scalar-valued observed outcome.
        strata: Stratum indicators.

    Returns:
        DistributionEstimatorBase: The fitted estimator.
    """
    covariates = _convert_to_ndarray(covariates)
    treatment_arms = _convert_to_ndarray(treatment_arms)
    outcomes = _convert_to_ndarray(outcomes)
    strata = _convert_to_ndarray(strata)

    if covariates.shape[0] != treatment_arms.shape[0]:
        raise ValueError("The shape of covariates and treatment_arm should be same")

    if covariates.shape[0] != outcomes.shape[0]:
        raise ValueError("The shape of covariates and outcome should be same")

    self.covariates = covariates
    self.treatment_arms = treatment_arms
    self.outcomes = outcomes
    self.strata = strata

    return self

AdjustedStratifiedDistributionEstimator

Bases: DistributionEstimatorBase

A class is for estimating the adjusted distribution function and computing the Distributional parameters for CAR.

Initializes the AdjustedDistributionEstimator.

Args: base_model (scikit-learn estimator): The base model implementing used for conditional distribution function estimators. The model should implement fit(data, targets) and predict_proba(data). folds (int): The number of folds for cross-fitting. is_multi_task(bool): Whether to use multi-task learning. If True, your base model needs to support multi-task prediction (n_samples, n_features) -> (n_samples, n_targets).

Returns: AdjustedDistributionEstimator: An instance of the estimator.

Source code in dte_adj/stratified.py
def __init__(self, base_model: Any, folds=3, is_multi_task=False):
    """
    Initializes the AdjustedDistributionEstimator.

    Args:
        base_model (scikit-learn estimator): The base model implementing used for conditional distribution function estimators. The model should implement fit(data, targets) and predict_proba(data).
        folds (int): The number of folds for cross-fitting.
        is_multi_task(bool): Whether to use multi-task learning. If True, your base model needs to support multi-task prediction (n_samples, n_features) -> (n_samples, n_targets).

    Returns:
        AdjustedDistributionEstimator: An instance of the estimator.
    """
    if (not hasattr(base_model, "predict")) and (
        not hasattr(base_model, "predict_proba")
    ):
        raise ValueError(
            "Base model should implement either predict_proba or predict"
        )
    self.base_model = base_model
    self.folds = folds
    self.is_multi_task = is_multi_task
    super().__init__()

predict_dte

predict_dte(
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]

Compute Distribution Treatment Effects (DTE) based on the estimator for the distribution function.

The DTE measures the difference in cumulative distribution functions between treatment groups at specified locations. It quantifies how treatment affects the probability of observing outcomes below each threshold.

Args: target_treatment_arm (int): The index of the treatment arm of the treatment group. control_treatment_arm (int): The index of the treatment arm of the control group. locations (np.ndarray, optional): Scalar values to be used for computing the cumulative distribution. If None, evenly-spaced locations spanning the observed outcome range are generated automatically. The number of points is determined from data size and distribution via np.histogram_bin_edges(outcomes, bins='auto'). The actual array used is stored on self.last_locations. alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05. variance_type (str, optional): Variance type to be used to compute confidence intervals. Available values are "moment", "simple", and "uniform". Defaults to "moment". n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing: - Expected DTEs (np.ndarray): Treatment effect estimates at each location - Lower bounds (np.ndarray): Lower confidence interval bounds - Upper bounds (np.ndarray): Upper confidence interval bounds

Example:

import numpy as np
from dte_adj import SimpleDistributionEstimator

# Generate sample data
X = np.random.randn(1000, 5)
D = np.random.binomial(1, 0.5, 1000)
Y = X[:, 0] + 2 * D + np.random.randn(1000)

# Fit estimator
estimator = SimpleDistributionEstimator()
estimator.fit(X, D, Y)

# Compute DTE
locations = np.linspace(Y.min(), Y.max(), 20)
dte, lower, upper = estimator.predict_dte(
    target_treatment_arm=1,
    control_treatment_arm=0,
    locations=locations,
    variance_type="moment"
)

print(f"DTE shape: {dte.shape}")  # Should match locations.shape
print(f"Average DTE: {dte.mean():.3f}")

Source code in dte_adj/base.py
def predict_dte(
    self,
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[np.ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
    """
    Compute Distribution Treatment Effects (DTE) based on the estimator for the distribution function.

    The DTE measures the difference in cumulative distribution functions between treatment groups
    at specified locations. It quantifies how treatment affects the probability of observing
    outcomes below each threshold.

    Args:
        target_treatment_arm (int): The index of the treatment arm of the treatment group.
        control_treatment_arm (int): The index of the treatment arm of the control group.
        locations (np.ndarray, optional): Scalar values to be used for computing the cumulative
            distribution. If None, evenly-spaced locations spanning the observed outcome range
            are generated automatically. The number of points is determined from data size and
            distribution via ``np.histogram_bin_edges(outcomes, bins='auto')``. The actual array
            used is stored on ``self.last_locations``.
        alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05.
        variance_type (str, optional): Variance type to be used to compute confidence intervals.
            Available values are "moment", "simple", and "uniform". Defaults to "moment".
        n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing:
            - Expected DTEs (np.ndarray): Treatment effect estimates at each location
            - Lower bounds (np.ndarray): Lower confidence interval bounds
            - Upper bounds (np.ndarray): Upper confidence interval bounds

    Example:
        ```python
        import numpy as np
        from dte_adj import SimpleDistributionEstimator

        # Generate sample data
        X = np.random.randn(1000, 5)
        D = np.random.binomial(1, 0.5, 1000)
        Y = X[:, 0] + 2 * D + np.random.randn(1000)

        # Fit estimator
        estimator = SimpleDistributionEstimator()
        estimator.fit(X, D, Y)

        # Compute DTE
        locations = np.linspace(Y.min(), Y.max(), 20)
        dte, lower, upper = estimator.predict_dte(
            target_treatment_arm=1,
            control_treatment_arm=0,
            locations=locations,
            variance_type="moment"
        )

        print(f"DTE shape: {dte.shape}")  # Should match locations.shape
        print(f"Average DTE: {dte.mean():.3f}")
        ```
    """
    if locations is None:
        locations = _infer_default_locations(self.outcomes, for_intervals=False)
    self.last_locations = locations
    return self._compute_dtes(
        target_treatment_arm,
        control_treatment_arm,
        locations,
        alpha,
        variance_type,
        n_bootstrap,
        display_progress,
    )

predict_pte

predict_pte(
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]

Compute Probability Treatment Effects (PTE) based on the estimator for the distribution function.

The PTE measures the difference in probability mass between treatment groups for intervals defined by consecutive location pairs. It quantifies how treatment affects the probability of observing outcomes within specific ranges.

Args: target_treatment_arm (int): The index of the treatment arm of the treatment group. control_treatment_arm (int): The index of the treatment arm of the control group. locations (np.ndarray, optional): Scalar values defining interval boundaries for probability computation. For each interval (locations[i], locations[i+1]], the PTE is computed. If None, boundaries spanning the observed outcome range are generated automatically with the left endpoint placed just below outcomes.min() so that minimum-valued samples fall inside the first interval. The number of boundaries is determined from data size and distribution via np.histogram_bin_edges(outcomes, bins='auto'). The actual array used is stored on self.last_locations. alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05. variance_type (str, optional): Variance type to be used to compute confidence intervals. Available values are "moment", "simple", and "uniform". Defaults to "moment". n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing: - Expected PTEs (np.ndarray): Treatment effect estimates for each interval, shape (len(locations)-1,) - Lower bounds (np.ndarray): Lower confidence interval bounds - Upper bounds (np.ndarray): Upper confidence interval bounds

Example:

import numpy as np
from dte_adj import SimpleDistributionEstimator

# Generate sample data
X = np.random.randn(1000, 5)
D = np.random.binomial(1, 0.5, 1000)
Y = X[:, 0] + 2 * D + np.random.randn(1000)

# Fit estimator
estimator = SimpleDistributionEstimator()
estimator.fit(X, D, Y)

# Define interval boundaries
locations = np.array([-2, -1, 0, 1, 2])  # Creates intervals: (-2,-1], (-1,0], (0,1], (1,2]

# Compute PTE
pte, lower, upper = estimator.predict_pte(
    target_treatment_arm=1,
    control_treatment_arm=0,
    locations=locations,
    variance_type="moment"
)

print(f"PTE shape: {pte.shape}")  # Should be (4,) for 4 intervals
print(f"Interval effects: {pte}")

Source code in dte_adj/base.py
def predict_pte(
    self,
    target_treatment_arm: int,
    control_treatment_arm: int,
    locations: Optional[np.ndarray] = None,
    alpha: float = 0.05,
    variance_type="moment",
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
    """
    Compute Probability Treatment Effects (PTE) based on the estimator for the distribution function.

    The PTE measures the difference in probability mass between treatment groups for intervals
    defined by consecutive location pairs. It quantifies how treatment affects the probability
    of observing outcomes within specific ranges.

    Args:
        target_treatment_arm (int): The index of the treatment arm of the treatment group.
        control_treatment_arm (int): The index of the treatment arm of the control group.
        locations (np.ndarray, optional): Scalar values defining interval boundaries for
            probability computation. For each interval (locations[i], locations[i+1]], the PTE
            is computed. If None, boundaries spanning the observed outcome range are generated
            automatically with the left endpoint placed just below ``outcomes.min()`` so that
            minimum-valued samples fall inside the first interval. The number of boundaries is
            determined from data size and distribution via
            ``np.histogram_bin_edges(outcomes, bins='auto')``. The actual array used is stored
            on ``self.last_locations``.
        alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05.
        variance_type (str, optional): Variance type to be used to compute confidence intervals.
            Available values are "moment", "simple", and "uniform". Defaults to "moment".
        n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing:
            - Expected PTEs (np.ndarray): Treatment effect estimates for each interval,
              shape (len(locations)-1,)
            - Lower bounds (np.ndarray): Lower confidence interval bounds
            - Upper bounds (np.ndarray): Upper confidence interval bounds

    Example:
        ```python
        import numpy as np
        from dte_adj import SimpleDistributionEstimator

        # Generate sample data
        X = np.random.randn(1000, 5)
        D = np.random.binomial(1, 0.5, 1000)
        Y = X[:, 0] + 2 * D + np.random.randn(1000)

        # Fit estimator
        estimator = SimpleDistributionEstimator()
        estimator.fit(X, D, Y)

        # Define interval boundaries
        locations = np.array([-2, -1, 0, 1, 2])  # Creates intervals: (-2,-1], (-1,0], (0,1], (1,2]

        # Compute PTE
        pte, lower, upper = estimator.predict_pte(
            target_treatment_arm=1,
            control_treatment_arm=0,
            locations=locations,
            variance_type="moment"
        )

        print(f"PTE shape: {pte.shape}")  # Should be (4,) for 4 intervals
        print(f"Interval effects: {pte}")
        ```
    """
    if locations is None:
        locations = _infer_default_locations(self.outcomes, for_intervals=True)
    self.last_locations = locations
    return self._compute_ptes(
        target_treatment_arm,
        control_treatment_arm,
        locations,
        alpha,
        variance_type,
        n_bootstrap,
        display_progress,
    )

predict_qte

predict_qte(
    target_treatment_arm: int,
    control_treatment_arm: int,
    quantiles: Optional[ndarray] = None,
    alpha: float = 0.05,
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]

Compute Quantile Treatment Effects (QTE) based on the estimator for the distribution function.

The QTE measures the difference in quantiles between treatment groups, providing insights into how treatment affects different parts of the outcome distribution. For stratified estimators, the computation properly accounts for strata.

Variance is estimated by stratified bootstrap: indices are resampled with replacement within each stratum independently, which preserves per-stratum sample sizes and reflects the covariate-adaptive randomization (CAR) design. For estimators without strata (single stratum), this degenerates to a plain bootstrap.

Args: target_treatment_arm (int): The index of the treatment arm of the treatment group. control_treatment_arm (int): The index of the treatment arm of the control group. quantiles (np.ndarray, optional): Quantiles used for QTE. Defaults to [0.1, 0.2, ..., 0.9]. alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05. n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing: - Expected QTEs (np.ndarray): Treatment effect estimates at each quantile - Lower bounds (np.ndarray): Lower confidence interval bounds - Upper bounds (np.ndarray): Upper confidence interval bounds

Example:

import numpy as np
from dte_adj import SimpleStratifiedDistributionEstimator

# Generate stratified sample data
X = np.random.randn(1000, 5)
strata = np.random.choice([0, 1, 2], size=1000)
D = np.random.binomial(1, 0.5, 1000)
Y = X[:, 0] + 2 * D + 0.5 * strata + np.random.randn(1000)

# Fit stratified estimator
estimator = SimpleStratifiedDistributionEstimator()
estimator.fit(X, D, Y, strata)

# Compute QTE at specific quantiles
quantiles = np.array([0.25, 0.5, 0.75])  # 25th, 50th, 75th percentiles
qte, lower, upper = estimator.predict_qte(
    target_treatment_arm=1,
    control_treatment_arm=0,
    quantiles=quantiles,
    n_bootstrap=100
)

print(f"QTE at quantiles {quantiles}: {qte}")
print(f"Median effect (50th percentile): {qte[1]:.3f}")

Source code in dte_adj/base.py
def predict_qte(
    self,
    target_treatment_arm: int,
    control_treatment_arm: int,
    quantiles: Optional[np.ndarray] = None,
    alpha: float = 0.05,
    n_bootstrap=500,
    display_progress: bool = True,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
    """
    Compute Quantile Treatment Effects (QTE) based on the estimator for the distribution function.

    The QTE measures the difference in quantiles between treatment groups, providing insights
    into how treatment affects different parts of the outcome distribution. For stratified
    estimators, the computation properly accounts for strata.

    Variance is estimated by stratified bootstrap: indices are resampled with replacement
    within each stratum independently, which preserves per-stratum sample sizes and reflects
    the covariate-adaptive randomization (CAR) design. For estimators without strata
    (single stratum), this degenerates to a plain bootstrap.

    Args:
        target_treatment_arm (int): The index of the treatment arm of the treatment group.
        control_treatment_arm (int): The index of the treatment arm of the control group.
        quantiles (np.ndarray, optional): Quantiles used for QTE. Defaults to [0.1, 0.2, ..., 0.9].
        alpha (float, optional): Significance level of the confidence bound. Defaults to 0.05.
        n_bootstrap (int, optional): Number of bootstrap samples. Defaults to 500.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        Tuple[np.ndarray, np.ndarray, np.ndarray]: A tuple containing:
            - Expected QTEs (np.ndarray): Treatment effect estimates at each quantile
            - Lower bounds (np.ndarray): Lower confidence interval bounds
            - Upper bounds (np.ndarray): Upper confidence interval bounds

    Example:
        ```python
        import numpy as np
        from dte_adj import SimpleStratifiedDistributionEstimator

        # Generate stratified sample data
        X = np.random.randn(1000, 5)
        strata = np.random.choice([0, 1, 2], size=1000)
        D = np.random.binomial(1, 0.5, 1000)
        Y = X[:, 0] + 2 * D + 0.5 * strata + np.random.randn(1000)

        # Fit stratified estimator
        estimator = SimpleStratifiedDistributionEstimator()
        estimator.fit(X, D, Y, strata)

        # Compute QTE at specific quantiles
        quantiles = np.array([0.25, 0.5, 0.75])  # 25th, 50th, 75th percentiles
        qte, lower, upper = estimator.predict_qte(
            target_treatment_arm=1,
            control_treatment_arm=0,
            quantiles=quantiles,
            n_bootstrap=100
        )

        print(f"QTE at quantiles {quantiles}: {qte}")
        print(f"Median effect (50th percentile): {qte[1]:.3f}")
        ```
    """
    if quantiles is None:
        quantiles = np.arange(1, 10) / 10
    if np.any((quantiles <= 0) | (quantiles >= 1)):
        raise ValueError("quantiles must be in the open interval (0, 1)")

    qte = self._compute_qtes(
        target_treatment_arm,
        control_treatment_arm,
        quantiles,
        self.covariates,
        self.treatment_arms,
        self.outcomes,
        self.strata,
    )

    # Precompute stratum indices for stratified bootstrap.
    # When there is a single stratum this is equivalent to plain bootstrap.
    unique_strata = np.unique(self.strata)
    strata_indices = [np.where(self.strata == s)[0] for s in unique_strata]

    qtes = np.zeros((n_bootstrap, qte.shape[0]))
    bootstrap_iter = range(n_bootstrap)
    if display_progress:
        bootstrap_iter = tqdm(bootstrap_iter, desc="Bootstrap QTE")
    for b in bootstrap_iter:
        bootstrap_indexes = np.concatenate(
            [
                np.random.choice(idx, size=len(idx), replace=True)
                for idx in strata_indices
            ]
        )

        qtes[b] = self._compute_qtes(
            target_treatment_arm,
            control_treatment_arm,
            quantiles,
            self.covariates[bootstrap_indexes],
            self.treatment_arms[bootstrap_indexes],
            self.outcomes[bootstrap_indexes],
            self.strata[bootstrap_indexes],
        )

    qte_var = qtes.var(axis=0)

    qte_lower = qte + norm.ppf(alpha / 2) * np.sqrt(qte_var)
    qte_upper = qte + norm.ppf(1 - alpha / 2) * np.sqrt(qte_var)

    return qte, qte_lower, qte_upper

predict

predict(
    treatment_arm: int,
    locations: ndarray,
    display_progress: bool = True,
) -> np.ndarray

Compute cumulative distribution values.

Args: treatment_arm (int): The index of the treatment arm. locations (np.ndarray): Scalar values to be used for computing the cumulative distribution. display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

Returns: np.ndarray: Estimated cumulative distribution values for the input.

Source code in dte_adj/base.py
def predict(
    self, treatment_arm: int, locations: np.ndarray, display_progress: bool = True
) -> np.ndarray:
    """
    Compute cumulative distribution values.

    Args:
        treatment_arm (int): The index of the treatment arm.
        locations (np.ndarray): Scalar values to be used for computing the cumulative distribution.
        display_progress (bool, optional): Whether to display a progress bar. Defaults to True.

    Returns:
        np.ndarray: Estimated cumulative distribution values for the input.
    """
    if self.outcomes is None:
        raise ValueError(
            "This estimator has not been trained yet. Please call fit first"
        )

    if treatment_arm not in self.treatment_arms:
        raise ValueError(
            f"This target treatment arm was not included in the training data: {treatment_arm}"
        )

    return self._compute_cumulative_distribution(
        treatment_arm,
        locations,
        self.covariates,
        self.treatment_arms,
        self.outcomes,
        display_progress=display_progress,
    )[0]

fit

fit(
    covariates: ArrayLike,
    treatment_arms: ArrayLike,
    outcomes: ArrayLike,
    strata: ArrayLike,
) -> DistributionEstimatorBase

Train the DistributionEstimatorBase.

Args: covariates: Pre-treatment covariates. treatment_arms: The index of the treatment arm. outcomes: Scalar-valued observed outcome. strata: Stratum indicators.

Returns: DistributionEstimatorBase: The fitted estimator.

Source code in dte_adj/stratified.py
def fit(
    self,
    covariates: ArrayLike,
    treatment_arms: ArrayLike,
    outcomes: ArrayLike,
    strata: ArrayLike,
) -> DistributionEstimatorBase:
    """
    Train the DistributionEstimatorBase.

    Args:
        covariates: Pre-treatment covariates.
        treatment_arms: The index of the treatment arm.
        outcomes: Scalar-valued observed outcome.
        strata: Stratum indicators.

    Returns:
        DistributionEstimatorBase: The fitted estimator.
    """
    covariates = _convert_to_ndarray(covariates)
    treatment_arms = _convert_to_ndarray(treatment_arms)
    outcomes = _convert_to_ndarray(outcomes)
    strata = _convert_to_ndarray(strata)

    if covariates.shape[0] != treatment_arms.shape[0]:
        raise ValueError("The shape of covariates and treatment_arm should be same")

    if covariates.shape[0] != outcomes.shape[0]:
        raise ValueError("The shape of covariates and outcome should be same")

    self.covariates = covariates
    self.treatment_arms = treatment_arms
    self.outcomes = outcomes
    self.strata = strata

    return self

  1. Byambadalai, U., Hirata, T., Oka, T., & Yasui, S. (2025). On Efficient Estimation of Distributional Treatment Effects under Covariate-Adaptive Randomization. arXiv preprint arXiv:2506.05945