footnote[page]

# COSTAR: Improved Temporal Counterfactual Estimation with Self-Supervised
Learning

Chuizheng Meng
<EMAIL_ADDRESS>
&Yihe Dong 22footnotemark: 2
<EMAIL_ADDRESS>
&Sercan Ö. Arık 22footnotemark: 2
<EMAIL_ADDRESS>
Yan Liu 11footnotemark: 1 22footnotemark: 2
<EMAIL_ADDRESS>&Tomas Pfister 22footnotemark: 2
<EMAIL_ADDRESS>
University of Southern CaliforniaGoogle Cloud AI Research

###### Abstract

Estimation of temporal counterfactual outcomes from observed history is
crucial for decision-making in many domains such as healthcare and e-commerce,
particularly when randomized controlled trials (RCTs) suffer from high cost or
impracticality. For real-world datasets, modeling time-dependent confounders
is challenging due to complex dynamics, long-range dependencies and both past
treatments and covariates affecting the future outcomes. In this paper, we
introduce COunterfactual Self-supervised TrAnsformeR (COSTAR), a novel
approach that integrates self-supervised learning for improved historical
representations. The proposed framework combines temporal and feature-wise
attention with a component-wise contrastive loss tailored for temporal
treatment outcome observations, yielding superior performance in estimation
accuracy and generalization to out-of-distribution data compared to existing
models, as validated by empirical results on both synthetic and real-world
datasets. ††Code available at https://github.com/google-research/google-
research/tree/master/COSTAR

##  1 Introduction

Accurate estimation of treatment outcomes over time conditioning on the
observed history is a fundamental problem in causal analysis and decision
making in various applications (Mahar et al., 2021; Ye et al., 2023; Wang et
al., 2023). For example, in medical domains, doctors are interested in knowing
how a patient reacts to a treatment or multi-step treatments; in e-commerce,
retailers are concerned about how future sales change if adjusting the price
of an item. While randomized controlled trials (RCTs) are the gold standard
for treatment outcome estimation, most often than not such trials are either
too costly or even impractical to conduct. Therefore, utilizing available
observed data (such as electronic health records (EHRs) and historical sales)
for accurate treatment outcome estimation, has drawn increasing interest in
the community.

Compared to the well-studied i.i.d cases, treatment outcome estimation from
time series observations not only finds more applications in the real world
but also pose significant more challenges, due to the complex dynamics and the
long-range dependencies in time series. Existing works along this endeavors
explore various architectures with improved capacity and training strategies
to alleviate time-dependent confounding. Recurrent marginal structural
networks (RMSNs) (Lim, 2018), counterfactual recurrent networks (CRN) (Bica et
al., 2020), and G-Net (Li et al., 2021) utilize architectures based on
recurrent neural networks. To mitigate time-dependent confounding, they train
proposed models with inverse probability of treatment weighting (IPTW),
treatment invariant representation through gradient reversal, and
G-computation respectively, in addition to the factual estimation loss on
observed data. Causal Transformer (CT) (Melnychuk et al., 2022) further
improves capturing long-range dependencies in the observational data with a
tailored transformer-based architecture and overcomes the temporal confounding
with balanced representations trained through counterfactual domain confusion
loss.

Figure 1: We illustrate the problem of treatment outcome estimation over time
with an example in healthcare. We propose COSTAR as a temporal counterfactual
estimator enhanced with self-supervised learning, inducing transferrability to
both ① cold-start cases from unseen subpopulations and ② counterfactual
outcome estimation.

While existing methods achieve performance gain in empirical evaluation, they
rely on the fully supervised loss of future outcomes to learn representations
of history and thus suffer from its limitations. In many practical
applications, we are confronted with the cold case challenge, where no or
limited observations of testing time series are accessible. For example, after
training the model with historical sales (seen as outcomes) and pricing (seen
as treatments) sequences of products in the food category, we are asked to
estimate the sales of new items in the household category given their pricing,
with no or very limited observations of household items collected beforehand.
Existing methods based on supervised learning have difficulty generalizing to
different domains and handling cold cases in test time.

In this work, we propose a paradigm shift from supervised learning to self-
supervised training for temporal treatment outcome estimation. Our proposed
model, COunterfactual Self-supervised TrAnsformeR (COSTAR), addresses the
aforementioned limitations. To enhance the model capacity, we propose an
encoder architecture composed of alternating temporal and feature-wise
attention, capturing dependencies among both time steps and features. To learn
expressive and transferable representations of the observed history, we refine
the contrastive loss in self-supervised learning to a finer-grained level:
both the entire history and each of the covariate/treatment/outcome components
are contrasted when constructing the loss. Moreover, we view the
counterfactual outcome estimation problem from the unsupervised domain
adaptation (UDA) perspective and provide the theoretical analysis of the error
bound for a counterfactual outcome estimator that gives estimation based on
representations from self-supervised learning.

Our main contributions are summarized as follows:

  1. 1.

We adapt self-supervised learning (SSL) together with component-wise
contrastive losses tailored for temporal observations to learn more expressive
representations of the history in temporal counterfactual outcome estimation.

  2. 2.

We propose a new encoder architecture combining both temporal attention and
feature-wise attention for better modeling the complex temporal dependencies
and feature interactions in observations.

  3. 3.

Empirical results show that our proposed framework outperforms existing
baselines across both synthetic and real-world datasets in both estimation
accuracy and generalization. In addition, we demonstrate that the learned
representations are balanced towards treatments and thus address the temporal
confounding issue.

##  2 Related Work

#### Counterfactual treatment outcome estimation over time.

Early works in counterfactual treatment outcome estimation were first
developed for epidemiology and can be considered under 3 major groups:
G-computation, marginal structural models (MSMs), and structural nested models
(Robins, 1986, 1994; Robins et al., 2000; Robins & Hernan, 2008). One major
shortcoming of these is that they are built on linear models and suffer from
the limited model capacity when facing time series data with complex temporal
dependencies. Follow-up works address the limitation in expressiveness with
Bayesian non-parametric methods (Xu et al., 2016; Soleimani et al., 2017;
Schulam & Saria, 2017) or more expressive deep neural networks (DNNs) such as
recurrent neural networks (RNNs). For example, recurrent marginal structural
networks (RMSNs) (Lim, 2018) replace the linear model in MSM with an RNN-based
architecture for forecasting treatment outcomes. G-Net (Li et al., 2021) also
adopts RNN instead of classical regression models in the g-computation
framework. Inspired by the success of representation learning for domain
adaptation and generalization (Ganin et al., 2016; Tzeng et al., 2015), more
recent works explore learning representations that are both predictive for
outcome estimation and balanced regardless of the treatment bias in training
data. Counterfactual recurrent network (CRN) (Bica et al., 2020) trains an
RNN-based model with both the factual outcome regression loss and the gradient
reversal (Ganin et al., 2016) w.r.t. the treatment prediction loss. The former
loss encourages the learned representations to be predictive of outcomes while
the latter encourages the representations to be homogeneous given different
treatments. The joint training target leads to informative and balanced
representations. With similar motivations, (Melnychuk et al., 2022) replaces
the RNN-based architecture with a Transformer-based (Vaswani et al., 2017) one
along with the domain confusion loss (Tzeng et al., 2015) to learn treatment-
agnostic representations. Given the flexibility of the choice of model
architectures, recent works extend temporal counterfactual outcome estimation
to irregular time series (Seedat et al., 2022; Cao et al., 2023), temporal
point process (Zhang et al., 2022b), and graph-structured spatiotemporal data
(Jiang et al., 2023) with the help of (Kidger et al., 2020) and (Huang et al.,
2020). While existing works claim that both predictive and balanced
representations are critical in accurate counterfactual outcome estimation, we
empirically find that the impact of representation balancing is inconsistent
and marginal. In contrast, improving the expressiveness of representations
brings more robust improvements.

#### Self-supervised learning of time series.

Being widely studied first for computer vision tasks (He et al., 2020; Chen et
al., 2020; Grill et al., 2020; Chen* et al., 2021), self-supervised learning
achieves strong performance with the advantage of not relying on labeled data.
Recent works (Yue et al., 2022; Tonekaboni et al., 2021; Woo et al., 2022;
Zhang et al., 2022a) further generalize and adapt self-supervised learning
methods for time series, including classification, forecasting, and anomaly
detection tasks. However, existing works of counterfactual outcome estimation
over time have neglected self-supervised learning of time series as an
effective way of learning informative representations. Meanwhile, existing
models for self-supervised learning of time series are not tailored for
counterfactual outcome estimation. Hence we propose COSTAR to mitigate the
gap.

##  3 Problem Formulation

Our task is estimating the outcomes of subjects with observed history after
being applied a sequence of treatments from observational data (Lim, 2018;
Bica et al., 2020; Melnychuk et al., 2022). We represent the available
observed dataset as
$\\{\\{{\bm{x}}^{(i)}_{t},{\bm{a}}^{(i)}_{t},{\bm{y}}^{(i)}_{t}\\}_{t=1}^{T^{(i)}},{\bm{v}}^{(i)}\\}_{i=1}^{N}$
of $N$ independently sampled subjects, where $T^{(i)}\in\mathbb{N}^{+}$
denotes the length of the observed history of subject $i$, and
${\bm{x}}^{(i)}_{t}\in\mathbb{R}^{d_{X}}$,
${\bm{a}}^{(i)}_{t}\in\mathbb{R}^{d_{A}}$, and
${\bm{y}}^{(i)}_{t}\in\mathbb{R}^{d_{Y}}$ stand for the observed vector of
covariates, treatments, and outcomes respectively, at time $t$ of subject $i$.
${\bm{v}}^{(i)}\in\mathbb{R}^{d_{V}}$ contain all static features of subject
$i$. We omit the subject index $i$ in the following text for notational
simplicity.

Following the potential outcomes (Splawa-Neyman et al., 1990; Rubin, 1978)
framework extended to time-varying treatments and outcomes (Robins & Hernan,
2008), our target is to estimate
$\mathbb{E}({\bm{y}}_{t+\tau}[\bar{{\bm{a}}}_{t:t+\tau-1}]|\bar{{\bm{H}}}_{t})$
for $\tau\geq 1$, where
$\bar{{\bm{H}}}_{t}=(\bar{{\bm{X}}}_{t},\bar{{\bm{A}}}_{t-1},\bar{{\bm{Y}}}_{t},{\bm{V}})$
is the observed history.
$\bar{{\bm{X}}}_{t}=({\bm{x}}_{1},{\bm{x}}_{2},\dots,{\bm{x}}_{t})$,
$\bar{{\bm{A}}}_{t-1}=({\bm{a}}_{1},{\bm{a}}_{2},\dots,{\bm{a}}_{t-1})$,
$\bar{{\bm{Y}}}_{t}=({\bm{y}}_{1},{\bm{y}}_{2},\dots,{\bm{y}}_{t})$,
${\bm{V}}={\bm{v}}$.
$\bar{{\bm{a}}}_{t:t+\tau-1}=({\bm{a}}_{t},{\bm{a}}_{t+1},\dots,{\bm{a}}_{t+\tau-1})$
is the sequence of the applied treatments in the future $\tau$ discrete time
steps. To ensure the identifiability of treatment effects from observational
data, we take the standard assumptions used in existing works (Bica et al.,
2020; Melnychuk et al., 2022): (1) consistency, (2) positivity and (3)
sequential strong ignorability (See Appendix A).

##  4 Counterfactual Self-Supervised Transformer

We illustrate the detailed design of COSTAR. Our main goal is to learn
representations of observed history sequences that are informative for
counterfactual treatment outcome estimation, which we achieve by tailoring
both the representation encoder architecture and the self-supervised training
loss. On top of the representation learning, we also propose a simple yet
effective decoder for non-autoregressive outcome prediction and demonstrate
improvements in both the accuracy and the speed of multi-step estimation. Fig.
2 overviews the proposed framework.

(a) Encoder architecture.

(b) Self-supervised learning of history.

(c) Non-autoregressive outcome predictor.

Figure 2: Overview of COSTAR. (a) Encoder architecture. The Temporal Attention
Block applies temporal causal attention along the time dimension in parallel
for each feature, while the Feature-wise Attention Block calculates full self-
attention along the feature dimension in all time steps. (b) Self-supervised
learning of the history representations. Positive pairs are generated by
applying random transformations $\mathcal{T}(\cdot)$ on the same sample. We
construct component-wise contrastive losses of historical covariates,
treatments and outcomes in addition to the standard contrastive loss of the
entire sequence. (c) Non-autoregressive outcome predictor architecture.

###  4.1 Encoder architecture

For a given sequence of the observed history
$\bar{{\bm{H}}}_{t}\in\mathbb{R}^{t\times d_{\text{input}}}$ concatenated from
$(\bar{{\bm{X}}}_{t},\bar{{\bm{A}}}_{t-1},\bar{{\bm{Y}}}_{t})$ (we omit static
variables ${\bm{V}}$ here for simplicity and leave the processing details in
Sec. 4.1; $d_{\text{input}}=d_{X}+d_{A}+d_{Y}$), the encoder in Fig. 2(a) maps
the entire history to representations
$\\{{\bm{z}}^{i}_{t}\in\mathbb{R}^{d_{\text{model}}}\\}_{i=1}^{d_{\text{input}}}$
for each feature $f_{i}$. Then, we employ average pooling for feature-wise
representations from corresponding features to get the representations of
covariate, treatment, and outcome components, and all features for the
representation of the entire observed history. Denote the set of covariate,
treatment, and outcome variables as $\mathbf{F}_{X}$, $\mathbf{F}_{A}$, and
$\mathbf{F}_{Y}$ respectively. We have:

$\begin{array}[]{c}{\bm{z}}^{X}_{t}=\text{avg}(\\{{\bm{z}}^{i}_{t}\\}_{f_{i}\in\mathbf{F}_{X}}),~{}{\bm{z}}^{A}_{t}=\text{avg}(\\{{\bm{z}}^{i}_{t}\\}_{f_{i}\in\mathbf{F}_{A}}),~{}{\bm{z}}^{Y}_{t}=\text{avg}(\\{{\bm{z}}^{i}_{t}\\}_{f_{i}\in\mathbf{F}_{Y}}),~{}{\bm{z}}_{t}=\text{avg}(\\{{\bm{z}}^{i}_{t}\\}_{i=1}^{d_{\text{input}}}).\end{array}$
(1)

#### Input feature projection.

Assume the concatenation of time-varying variables
$(\bar{{\bm{X}}}_{t},\bar{{\bm{A}}}_{t-1},\bar{{\bm{Y}}}_{t})$ in history as
$\bar{{\bm{S}}}_{t}\in\mathbb{R}^{t\times d_{S}}$, $d_{S}=d_{X}+d_{A}+d_{Y}$,
and the static variables ${\bm{V}}\in\mathbb{R}^{d_{V}}$. We adopt a linear
transformation
$f_{\text{input}}:\mathbb{R}\rightarrow\mathbb{R}^{d_{\text{model}}}$ to map
$\bar{{\bm{S}}}_{t}$ and ${\bm{V}}$ to the embedding space as
${\bm{E}}^{S}\in\mathbb{R}^{t\times d_{S}\times d_{\text{model}}}$ and
${\bm{E}}^{V}\in\mathbb{R}^{d_{V}\times d_{\text{model}}}$ respectively, where

$\displaystyle{\bm{E}}^{S}[i,j]$
$\displaystyle=f_{\text{input}}(\bar{{\bm{S}}}_{t}[i,j]),~{}1\leq i\leq
t,~{}1\leq j\leq d_{S};$ (2) $\displaystyle{\bm{E}}^{V}[j]$
$\displaystyle=f_{\text{input}}({\bm{V}}[j]),~{}1\leq j\leq d_{V}.$ (3)

#### Feature positional encoding.

A shared feature projection function among all features is not sufficient to
encode the feature-specific information since the same scalar value represents
different semantics in different features. Meanwhile, feature-specific
information is also critical for modeling the interactions among features.
Therefore we enhance the input embedding with a positional encoding along the
feature dimension. Since the features can be grouped into covariates,
treatments and outcomes and form a hierarchical structure with 2 levels, we
model it with a learnable tree positional encoding (Shiv & Quirk, 2019; Wang
et al., 2021). Denote the lists of covariate, treatment, outcome, and static
features as ${\bm{F}}_{X}$, ${\bm{F}}_{A}$, ${\bm{F}}_{Y}$, and
${\bm{F}}_{V}$. For the $i$-th feature $f_{i}^{\bm{F}}$ in a certain feature
list ${\bm{F}}\in\\{{\bm{F}}_{X},{\bm{F}}_{A},{\bm{F}}_{Y},{\bm{F}}_{V}\\}$,
its positional encoding is:

$\begin{array}[]{l}{\bm{E}}_{\text{fea\\_pos}}(f_{i}^{\bm{F}})={\bm{E}}_{\text{fea}}\cdot\text{Concat}({\bm{e}}_{\bm{F}},{\bm{e}}_{i}),~{}\text{where}~{}{\bm{e}}_{\bm{F}}=\left\\{\begin{array}[]{l}(1,0,0,0)~{}\text{if
${\bm{F}}$ is ${\bm{F}}_{X}$}\\\ (0,1,0,0)~{}\text{if ${\bm{F}}$ is
${\bm{F}}_{A}$}\\\ (0,0,1,0)~{}\text{if ${\bm{F}}$ is ${\bm{F}}_{Y}$}\\\
(0,0,0,1)~{}\text{if ${\bm{F}}$ is
${\bm{F}}_{V}$}\end{array}\right.,\end{array}$ (4)

${\bm{e}}_{i}\in\mathbb{R}^{\max(d_{X},d_{A},d_{Y},d_{V})}$ is a one-hot
vector with only ${\bm{e}}_{i}[i]=1$.
${\bm{E}}_{\text{fea}}\in\mathbb{R}^{d_{\text{model}}\times(4+\max(d_{X},d_{A},d_{Y},d_{V}))}$
are learnable tree embedding weights. After obtaining the stacked feature
positional embeddings
${\bm{E}}_{\text{fea\\_pos}}^{S}\in\mathbb{R}^{d_{S}\times d_{\text{model}}}$,
${\bm{E}}_{\text{fea\\_pos}}^{V}\in\mathbb{R}^{d_{V}\times d_{\text{model}}}$
of time-varying and static features, we broadcast them to the shape of
${\bm{E}}^{S}$ and ${\bm{E}}^{V}$ respectively along the time dimension. The
embedded input is then the sum of input feature projection and feature
positional encoding:

$\displaystyle{\bm{Z}}^{S,(0)}={\bm{E}}^{S}+\text{Broadcast}({\bm{E}}^{S}_{\text{fea\\_pos}}),\quad{\bm{Z}}^{V,(0)}={\bm{E}}^{V}+\text{Broadcast}({\bm{E}}^{V}_{\text{fea\\_pos}}).$
(5)

#### Temporal attention block.

The temporal attention block is designed to capture the temporal dependencies
within each feature. We construct the block based on the self-attention part
in the conventional Transformer decoder (Vaswani et al., 2017). Considering
the importance of relative time interval in modeling treatment outcomes, we
adopt the relative positional encoding (Shaw et al., 2018; Melnychuk et al.,
2022) along the time dimension.

The temporal attention block in the $l$-th layer receives ${\bm{Z}}^{S,(l-1)}$
from the previous layer, reshapes it to $d_{S}$ sequences with lengths $t$,
and passes them through the block in parallel. The outputs
${\bm{Z}}^{S,(l)}_{\text{tmp}}$ have the same shape as ${\bm{Z}}^{S,(l-1)}$.
Since ${\bm{Z}}^{V,(l-1)}$ is static, we only pass it through the point-wise
feed-forward module and get ${\bm{Z}}^{V,(l)}_{\text{tmp}}$

#### Feature-wise attention block.

The feature-wise attention block models interactions among different features.
We reuse the architecture of the conventional Transformer encoder but replace
the positional encoding with the feature positional encoding as described. The
block in the $l$-th layer receives ${\bm{Z}}^{S,(l)}_{\text{tmp}}$ from the
temporal attention block and reshapes it to $t$ sequences, each with a length
$d_{S}$. We broadcast ${\bm{Z}}^{V,(l)}_{\text{tmp}}$ and concatenate it with
each sequence along the feature dimension to enable the attention among both
time-varying and static features. The concatenated $t$ sequences that we apply
attention to are:

${\bm{Z}}^{SV,(l)}_{\text{tmp}}=\text{Concat}({\bm{Z}}^{S,(l)}_{\text{tmp}},\text{Broadcast}({\bm{Z}}^{V,(l)}_{\text{tmp}}))\in\mathbb{R}^{t\times(d_{S}+d_{V})\times
d_{\text{model}}}.$ (6)

We apply full attention across all features and get ${\bm{Z}}^{SV,(l)}$ with
the same shape as ${\bm{Z}}^{SV,(l)}_{\text{tmp}}$. The propagated embeddings
of time-varying features are obtained as:

${\bm{Z}}^{S,(l)}={\bm{Z}}^{SV,(l)}[:,:d_{S},:]\in\mathbb{R}^{t\times
d_{S}\times d_{\text{model}}}.$ (7)

To keep the ${\bm{Z}}^{V,(l)}$ static after feature-wise attention, we only
propagate ${\bm{Z}}^{V,(l)}_{\text{tmp}}$ with full attention among static
features only. The updated embeddings of the static features
${\bm{Z}}^{V,(l)}\in\mathbb{R}^{d_{V}\times d_{\text{model}}}$ have the same
shape as ${\bm{Z}}^{V,(l)}_{\text{tmp}}$.

###  4.2 Self-supervised representation learning of the observed history

We employ pretraining for the encoder in a self-supervised way with the
contrastive learning objectives $\mathcal{L}^{X}_{\text{SSL}}$,
$\mathcal{L}^{A}_{\text{SSL}}$, $\mathcal{L}^{Y}_{\text{SSL}}$,
$\mathcal{L}^{H}_{\text{SSL}}$ for the component representations
${\bm{z}}^{X}_{t},{\bm{z}}^{A}_{t},{\bm{z}}^{Y}_{t}$ and the overall
representation ${\bm{z}}_{t}$ respectively. The overall self-supervised
learning loss is given as:

$\mathcal{L}_{\text{SSL}}=\mathcal{L}^{H}_{\text{SSL}}+(\mathcal{L}^{X}_{\text{SSL}}+\mathcal{L}^{A}_{\text{SSL}}+\mathcal{L}^{Y}_{\text{SSL}})/3.$
(8)

#### Self-supervised training.

We train our proposed encoder to learn representations of the history with a
self-supervised learning framework modified based on MoCo v3 (Chen* et al.,
2021) that achieves the state-of-the-art performance in self-supervised vision
transformer training. Following MoCo v3, we set up our proposed encoder
$f_{\text{enc}}$ as combination of a momentum encoder with the same
architecture and initial weights $f^{\text{mo}}_{\text{enc}}$, and a multi-
layer perceptron (MLP) as the prediction head
$f_{\text{pred\\_head}}:\mathbb{R}^{d_{\text{model}}}\rightarrow\mathbb{R}^{d_{\text{model}}}$.
We first apply the augmentations in (Woo et al., 2022), which includes
scaling, shifting and jittering, on each sample in the input batch
$\\{\bar{{\bm{H}}}^{(i)}_{t}\\}_{i=1}^{B}$ ($B$ is the batch size) and
generates the positive sample pair
$\\{\bar{{\bm{H}}}_{t}^{{}^{\prime}(i)}\\}_{i=1}^{B}$,
$\\{\bar{{\bm{H}}}_{t}^{{}^{\prime\prime}(i)}\\}_{i=1}^{B}$. Their
representations are encoded as follows:

$\begin{array}[]{ll}{\bm{z}}_{t}^{{}^{\prime}X(i)},{\bm{z}}_{t}^{{}^{\prime}A(i)},{\bm{z}}_{t}^{{}^{\prime}Y(i)},{\bm{z}}_{t}^{{}^{\prime}(i)}&=f_{\text{enc}}(\bar{{\bm{H}}}^{{}^{\prime}(i)}_{t}),\\\
{\bm{z}}_{t}^{{}^{\prime\prime}X(i)},{\bm{z}}_{t}^{{}^{\prime\prime}A(i)},{\bm{z}}_{t}^{{}^{\prime\prime}Y(i)},{\bm{z}}_{t}^{{}^{\prime\prime}(i)}&=f_{\text{enc}}(\bar{{\bm{H}}}^{{}^{\prime\prime}(i)}_{t}),\\\
{\bm{z}}_{t}^{{}^{\prime}X,\text{mo}(i)},{\bm{z}}_{t}^{{}^{\prime}A,\text{mo}(i)},{\bm{z}}_{t}^{{}^{\prime}Y,\text{mo}(i)},{\bm{z}}_{t}^{{}^{\prime}\text{mo}(i)}&=f_{\text{enc}}^{\text{mo}}(\bar{{\bm{H}}}^{{}^{\prime}(i)}_{t}),\\\
{\bm{z}}_{t}^{{}^{\prime\prime}X,\text{mo}(i)},{\bm{z}}_{t}^{{}^{\prime\prime}A,\text{mo}(i)},{\bm{z}}_{t}^{{}^{\prime\prime}Y,\text{mo}(i)},{\bm{z}}_{t}^{{}^{\prime\prime}\text{mo}(i)}&=f_{\text{enc}}^{\text{mo}}(\bar{{\bm{H}}}^{{}^{\prime\prime}(i)}_{t}).\end{array}$
(9)

The vanilla MoCo v3 method adopts the InfoNCE contrastive loss (Oord et al.,
2018) as the training objective:

$\displaystyle\mathcal{L}_{\text{SSL}}^{H}=$
$\displaystyle\mathcal{L}_{\text{InfoNCE}}(\\{f_{\text{pred\\_head}}({\bm{z}}_{t}^{{}^{\prime}(i)})\\}_{i=1}^{B},\\{{\bm{z}}_{t}^{{}^{\prime\prime}\text{mo}(i)}\\}_{i=1}^{B})$
(10)
$\displaystyle+\mathcal{L}_{\text{InfoNCE}}(\\{f_{\text{pred\\_head}}({\bm{z}}_{t}^{{}^{\prime\prime}(i)})\\}_{i=1}^{B},\\{{\bm{z}}_{t}^{{}^{\prime}\text{mo}(i)}\\}_{i=1}^{B}),$

where

$\mathcal{L}_{\text{InfoNCE}}(\\{{\bm{q}}^{(i)}\\}_{i=1}^{B},\\{{\bm{k}}^{(i)}\\}_{i=1}^{B})=-\frac{1}{B}\sum_{i=1}^{B}\log\frac{\exp(\cos\left\langle{{\bm{q}}^{(i)}},{{\bm{k}}^{(i)}}\right\rangle)}{\sum_{j=1}^{B}\exp(\cos\left\langle{{\bm{q}}^{(i)}},{{\bm{k}}^{(j)}}\right\rangle)}.$
(11)

#### Component-wise contrastive loss.

In addition to the contrastive loss of the overall representations, we enhance
the training with contrastive losses on each subset of covariates, treatments
and outcomes:

$\displaystyle\mathcal{L}_{\text{SSL}}^{(\cdot)}=$
$\displaystyle\mathcal{L}_{\text{InfoNCE}}(\\{f_{\text{pred\\_head}}({\bm{z}}_{t}^{{}^{\prime}(\cdot)(i)})\\}_{i=1}^{B},\\{{\bm{z}}_{t}^{{}^{\prime\prime}(\cdot),\text{mo}(i)}\\}_{i=1}^{B})$
(12)
$\displaystyle+\mathcal{L}_{\text{InfoNCE}}(\\{f_{\text{pred\\_head}}({\bm{z}}_{t}^{{}^{\prime\prime}(\cdot)(i)})\\}_{i=1}^{B},\\{{\bm{z}}_{t}^{{}^{\prime}(\cdot),\text{mo}(i)}\\}_{i=1}^{B}),$
(13)

where $(\cdot)$ is $X,A,Y$.

###  4.3 Non-autoregressive outcome predictor

The architecture of the proposed predictor model is shown in Fig. 2(c). At the
prediction stage, we first encode the observed history $\bar{{\bm{H}}}_{t}$
with the pretrained encoder, the treatment ${\bm{a}}_{t^{\prime}-1}$ is
modeled right before time $t^{\prime}=t+1,\dots,t+\tau-1$ with a 1x1
convolution layer, and the remaining treatment sequence
$({\bm{a}}_{t},{\bm{a}}_{t+1},\dots,{\bm{a}}_{t^{\prime}-2})$ with a 1D
convolution layer. Then, the concatenated encoding is fed into a multi-layer
perceptron (MLP) to predict future outcomes
$(\hat{{\bm{y}}}_{t+1},\hat{{\bm{y}}}_{t+2},\dots,\hat{{\bm{y}}}_{t+\tau})$.
We jointly train the predictor layers and fine tune the pretrained encoder
with the $L_{2}$ loss of factual outcome estimation weighted for each step:

$\mathcal{L}_{\text{est}}=\sum_{i=1}^{\tau}w_{i}\left\lVert\hat{{\bm{y}}}_{t+i}-{\bm{y}}_{t+i}\right\rVert_{2}^{2},$
(14)

where each $w_{i}$ is a hyperparameter satisfying $\sum_{i=1}^{\tau}w_{i}=1$
and various strategies of setting $w_{i}$s can be selected via validation
errors, which we will discuss in the ablation study (Sec. 5.4).

###  4.4 Unsupervised domain adaptation view of counterfactual outcome
estimation

For the given observed history $\bar{{\bm{H}}}_{t}$, treatment sequence
$\bar{{\bm{a}}}_{t:t+\tau-1}$ to apply, and the outcome
${\bm{y}}_{t+\tau}[{\bm{a}}_{t:t+\tau-1}](\bar{{\bm{H}}}_{t})$ to estimate, we
notice that learning a counterfactual treatment outcome estimator
$f_{{\bm{a}}}(\bar{{\bm{H}}}_{t})=\mathbb{E}({\bm{y}}_{t+\tau}[{\bm{a}}]|\bar{{\bm{H}}}_{t})$
with factual data specifically for a certain treatment sequence ${\bm{a}}$ can
be viewed as an unsupervised domain adaptation (UDA) problem with the
treatment value being discrete – any sample in the factual dataset can be
categorized into one of the two subsets (1)
$\mathcal{S}_{{\bm{a}}}=\\{(\bar{{\bm{H}}}_{t}^{(i)},\bar{{\bm{a}}}_{t:t+\tau-1}^{(i)},{\bm{y}}_{t+\tau}^{(i)})\\}_{\bar{{\bm{a}}}^{(i)}_{t:t+\tau-1}={\bm{a}}}$
and (2)
$\mathcal{T}_{\bar{{\bm{a}}}}=\\{(\bar{{\bm{H}}}_{t}^{(i)},\bar{{\bm{a}}}_{t:t+\tau-1}^{(i)},{\bm{y}}_{t+\tau}^{(i)})\\}_{\bar{{\bm{a}}}^{(i)}_{t:t+\tau-1}\neq{\bm{a}}}$.
With Assumption A.1, we have
${\bm{y}}_{t+\tau}^{(i)}={\bm{y}}_{t+\tau}[{\bm{a}}](\bar{{\bm{H}}}_{t}^{(i)})$
in $\mathcal{S}_{{\bm{a}}}$ and thus ${\bm{y}}_{t+\tau}^{(i)}$ is a label of
$f_{{\bm{a}}}(\bar{{\bm{H}}}_{t})$. This does not hold for
$\mathcal{T}_{{\bm{a}}}$, where
$\bar{{\bm{a}}}^{(i)}_{t:t+\tau-1}\neq{\bm{a}}$. Therefore,
$\mathcal{S}_{{\bm{a}}}$ and $\mathcal{T}_{{\bm{a}}}$ correspond to the
labeled and the unlabeled dataset in UDA. Considering the existence of
treatment bias,
$P_{\mathcal{T}_{{\bm{a}}}}(\bar{{\bm{H}}}_{t})=P(\bar{{\bm{H}}}_{t}|\bar{{\bm{a}}}_{t:t+\tau-1}\neq{\bm{a}})\neq
P(\bar{{\bm{H}}}_{t}|\bar{{\bm{a}}}_{t:t+\tau-1}={\bm{a}})=P_{\mathcal{S}_{{\bm{a}}}}(\bar{{\bm{H}}}_{t})$,
which corresponds to the distribution shift between the labeled source domain
and the unlabeled target domain in UDA. As a natural generalization, we
analyze the upper bound of contrastive learning for counterfactual outcome
estimation based on the transferability analysis of contrastive learning in
UDA (HaoChen et al., 2022):

######  Theorem 4.1 (Upper bound of counterfactual outcome estimator).

Suppose that Assumptions C.2, C.3, and C.4 hold for the set of observed
history $\mathcal{H}$ and its positive-pair graph $G(\mathcal{H},w)$, and the
representation dimension $k\geq 2m$. Let $r$ be a minimizer of the generalized
spectral contrastive loss on factual data and the regression head $f_{\bm{a}}$
be constructed in Alg. 1 with labeled data. We have

$\mathcal{E}_{\mathcal{D}}(f_{\bm{a}})\lesssim
P({\bm{a}})\mathcal{E}_{\mathcal{S}_{\bm{a}}}(f_{\bm{a}})+(1-P({\bm{a}}))\left[\epsilon^{2}+(4B^{2}-\epsilon^{2})\frac{r}{\alpha^{2}\gamma^{4}}\cdot\exp(-\Omega(\frac{\rho\gamma^{2}}{\alpha^{2}}))\right],$
(15)

where $P({\bm{a}})$ is the prior probability of the treatment ${\bm{a}}$ to
apply. $\mathcal{E}_{\mathcal{S}_{\bm{a}}}(f_{\bm{a}})$ is the outcome
estimation error of $f_{\bm{a}}$ in $\mathcal{S}_{{\bm{a}}}$, which can be
optimized with supervised learning. $\alpha$, $r$, $\gamma$, $\rho$ are
parameters in Assumptions C.2, C.3, C.4. $\epsilon$ is the hyperparameter in
Alg. 1. $B$ is the upper bound of outcome and predicted outcome
values111Bounded outcome values can be achieved through normalization.. When
$\gamma\geq\alpha^{1/2}$ and $\rho$ is comparable to $\alpha$,
$\rho\gamma^{2}\gg\alpha^{2}$ and lead to a small upper bound. We provide the
proof in Appendix C.

##  5 Experiments

#### Datasets.

Following the common evaluation setup for counterfactual treatment outcome
estimation overtime in (Lim, 2018; Bica et al., 2020; Melnychuk et al., 2022),
we use datasets from both synthetic simulation and real-world observation in
our experiments. We provide more detailed dataset description in Appendix D.
(1) Tumor growth. Following previous work (Lim, 2018; Bica et al., 2020;
Melnychuk et al., 2022), we run the pharmacokinetic-pharmacodynamic(PK-PD)
tumor growth simulation and generates a fully-synthetic dataset. The PK-PD
simulation (Geng et al., 2017) is a state-of-the-art bio-mathematical model
simulating the combined effect of chemotherapy and radiotherapy on tumor
volumes. (2) Semi-synthetic MIMIC-III. (Melnychuk et al., 2022) constructs a
semi-synthetic dataset by simulating outcomes under endogenous dependencies on
time and exogenous dependencies on observational patient trajectories that are
high dimensional and contain long-range dependencies. We include it in our
evaluation as a more challenging synthetic dataset. (3) M5. The M5 Forecasting
dataset (Makridakis et al., 2022) contains daily sales of Walmart stores
across three US states, along with the metadata of items and stores, as well
as explanatory variables such as price and special events. We transform it to
treatment outcome estimation task with the treatment variable of item price,
and the outcome variable of the sales of items. Covariate variables include
all remaining features. With synthetic data, we report the counterfactual
outcome estimation errors and compare the performance of COSTAR with
baselines. However, with real-world data, the counterfactual outcome cannot be
observed or simulated, thus, we only report the prediction errors of factual
outcome.

#### Feature distribution shifts.

To achieve a comprehensive evaluation of counterfactual outcome estimation
performance, we introduce feature distribution shifts into the datasets. For
each dataset, we split based on static characteristics of subjects into a
subset in the source domain and a subset in the target domain. Each subset is
further divided into train/validation/test sets. We summarize the main
statistics of datasets in Table 5 in the appendix.

#### Baselines.

We select comprehensive methods for estimating counterfactual outcomes over
time as baselines, including MSM (Robins et al., 2000), RMSN (Lim, 2018), CRN
(Bica et al., 2020), G-Net (Li et al., 2021), and Causal Transformer (CT)
(Melnychuk et al., 2022). MSM has difficulty converging when trained with
high-dimensional input in semi-synthetic MIMIC-III and M5 datasets and we thus
only evaluate it for tumor growth. We empirically find that the balanced
representation training losses proposed in CRN and CT do not bring a robust
improvement over their variants trained only with the empirical risk
minimization (ERM) on factual outcomes. Therefore, we also include these
variants, CRN(ERM) and CT(ERM), as baselines.

###  5.1 Zero-shot transfer setup

To showcase cold-start prediction capabilities, in this setup, we focus on the
performance on the target domain after training the model in the source
domain, with distributional difference in features. Results are shown in Table
1. COSTAR demonstrates the state-of-the-art performance in a majority of
horizons across datasets (4/6 in Tumor growth, 6/6 in Semi-synthetic MIMIC-III
and M5). On average, COSTAR decreases the outcome estimation errors by over
6.2%, 22.5% and 26.3% compared to baselines. Results demonstrate the strong
transferability of COSTAR in the zero-shot transfer setting.

Table 1: Results of zero-shot transfer setup for multi-step outcome
estimation. We report the mean $\pm$ standard deviation of Rooted Mean Squared
Errors (RMSEs $\downarrow$) over 5 runs. Bold: the best results. Underline:
the 2nd best results.

Dataset | Method | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---
Tumor growth | MSM | 1.0515$\pm$0.0674 | 0.5048$\pm$0.0591 | 0.7583$\pm$0.0831 | 0.9685$\pm$0.1066 | 1.1561$\pm$0.1243 | 1.3372$\pm$0.1356 | 0.9627$\pm$0.0923 | 6.2%
RMSN | 1.2406$\pm$0.1301 | 1.0914$\pm$0.0346 | 1.1315$\pm$0.0634 | 1.1583$\pm$0.0810 | 1.1674$\pm$0.0913 | 1.1531$\pm$0.0919 | 1.1571$\pm$0.0660 | 22.0%
CRN(ERM) | 1.2924$\pm$0.0772 | 1.1769$\pm$0.1058 | 1.1728$\pm$0.1136 | 1.1906$\pm$0.1106 | 1.1997$\pm$0.1061 | 1.1883$\pm$0.0985 | 1.2035$\pm$0.0901 | 25.0%
CRN | 1.5643$\pm$0.0667 | 10.1978$\pm$14.9899 | 14.5460$\pm$23.3314 | 17.2080$\pm$28.4348 | 18.8712$\pm$31.5756 | 19.9238$\pm$33.5369 | 13.7185$\pm$21.9748 | 93.4%
CT(ERM) | 0.9729$\pm$0.0718 | 1.0217$\pm$0.0292 | 1.1173$\pm$0.0457 | 1.1904$\pm$0.0395 | 1.2359$\pm$0.0618 | 1.2913$\pm$0.0939 | 1.1383$\pm$0.0251 | 20.7%
CT | 1.0272$\pm$0.1077 | 1.1428$\pm$0.2182 | 1.2708$\pm$0.2471 | 1.3608$\pm$0.2681 | 1.4166$\pm$0.2935 | 1.4322$\pm$0.3138 | 1.2751$\pm$0.2326 | 29.2%
G-Net | 1.0492$\pm$0.0529 | 1.0125$\pm$0.0767 | 1.1271$\pm$0.0876 | 1.2153$\pm$0.0777 | 1.2549$\pm$0.0727 | 1.2543$\pm$0.0678 | 1.1522$\pm$0.0537 | 21.7%
COSTAR | 0.8767$\pm$0.0492 | 0.7995$\pm$0.0853 | 0.8282$\pm$0.0801 | 0.9021$\pm$0.1062 | 0.9888$\pm$0.1280 | 1.0210$\pm$0.1168 | 0.9027$\pm$0.0814 | (-)
Semi- synthetic MIMIC-III | RMSN | 0.2551$\pm$0.0303 | 0.6641$\pm$0.1092 | 0.9107$\pm$0.1915 | 1.1217$\pm$0.2916 | 1.2821$\pm$0.3603 | 1.3950$\pm$0.4038 | 0.9381$\pm$0.2210 | 44.7%
CRN(ERM) | 0.2506$\pm$0.0303 | 0.5545$\pm$0.0917 | 0.7581$\pm$0.1112 | 0.9018$\pm$0.1547 | 1.0113$\pm$0.1941 | 1.1068$\pm$0.2324 | 0.7639$\pm$0.1238 | 32.1%
CRN | 0.4041$\pm$0.0537 | 0.8256$\pm$0.1767 | 1.0439$\pm$0.1958 | 1.1807$\pm$0.1725 | 1.3121$\pm$0.2229 | 1.4374$\pm$0.3089 | 1.0340$\pm$0.1606 | 49.8%
CT(ERM) | 0.2762$\pm$0.0804 | 0.5397$\pm$0.1181 | 0.6765$\pm$0.1417 | 0.7728$\pm$0.1636 | 0.8451$\pm$0.1850 | 0.9028$\pm$0.2070 | 0.6688$\pm$0.1472 | 22.5%
CT | 0.3138$\pm$0.0458 | 0.5992$\pm$0.0492 | 0.7576$\pm$0.0694 | 0.8695$\pm$0.0921 | 0.9510$\pm$0.1118 | 1.0128$\pm$0.1274 | 0.7506$\pm$0.0797 | 30.9%
G-Net | 0.5514$\pm$0.1502 | 0.9398$\pm$0.2384 | 1.2461$\pm$0.3321 | 1.4985$\pm$0.4024 | 1.7045$\pm$0.4463 | 1.8731$\pm$0.4660 | 1.3022$\pm$0.3367 | 60.2%
COSTAR | 0.2266$\pm$0.0249 | 0.4501$\pm$0.0893 | 0.5406$\pm$0.0987 | 0.5964$\pm$0.1020 | 0.6344$\pm$0.1040 | 0.6637$\pm$0.1052 | 0.5186$\pm$0.0869 | (-)
M5 | RMSN | 15.1616$\pm$2.0027 | 13.9966$\pm$0.5316 | 13.4899$\pm$1.2632 | 13.5162$\pm$1.7437 | 13.8004$\pm$2.0637 | 14.3366$\pm$2.3891 | 13.8280$\pm$1.5526 | 47.5%
CRN(ERM) | 9.8859$\pm$1.2980 | 20.8199$\pm$3.9049 | 38.2653$\pm$8.9897 | 59.4192$\pm$16.2788 | 82.9515$\pm$26.1928 | 105.8120$\pm$35.5325 | 61.4536$\pm$17.5760 | 88.2%
CRN | 8.1119$\pm$0.3183 | 10.3741$\pm$2.2616 | 12.9356$\pm$3.1588 | 15.4168$\pm$3.8002 | 18.1382$\pm$4.6750 | 21.1337$\pm$5.4694 | 15.5997$\pm$3.7687 | 53.5%
CT(ERM) | 7.1253$\pm$0.5777 | 8.3438$\pm$1.0313 | 9.2014$\pm$1.4146 | 9.9409$\pm$1.7572 | 10.6726$\pm$2.1718 | 11.3597$\pm$2.5966 | 9.9037$\pm$1.7852 | 26.7%
CT | 7.1239$\pm$0.5770 | 8.2939$\pm$0.9702 | 9.1465$\pm$1.3397 | 9.9091$\pm$1.7198 | 10.6311$\pm$2.0328 | 11.3032$\pm$2.4185 | 9.8568$\pm$1.6959 | 26.3%
G-Net | 7.5358$\pm$0.1605 | 8.6077$\pm$0.3166 | 9.7167$\pm$0.4861 | 10.8993$\pm$0.6902 | 12.3477$\pm$0.8940 | 13.8200$\pm$1.1193 | 10.4879$\pm$0.6078 | 30.8%
COSTAR | 6.4054$\pm$0.0547 | 6.9328$\pm$0.0634 | 7.2428$\pm$0.0700 | 7.4585$\pm$0.0580 | 7.7012$\pm$0.0627 | 7.8278$\pm$0.0651 | 7.2614$\pm$0.0609 | (-)

###  5.2 Data-efficient transfer learning setup

Effectively utilizing small amount of target domain data can be important, and
we showcase that it is indeed one of the key strengths of the proposed
approach.

To demonstrate this, for the Tumor Growth Dataset, we fine-tune each method
trained on the source domain with 100 sequences from the target domain. For
the semi-synthetic MIMIC-III and M5 datasets, we set the number of target
domain samples for fine-tuning to be 10% of the number of samples of the
target domain in the original dataset. To achieve a fair comparison, we fine-
tune each method until it reaches the lowest factual outcome estimation error
on a separate validation set in the target domain.

Table 2 compares the performance of all methods in data-efficient transfer
learning setup. For the majority of horizons (4/6 in Tumor growth, 5/6 in
Semi-synthetic MIMIC-III and M5), we observe that COSTAR achieves the state-
of-the-art performance after fine-tuning. Again, COSTAR reduces the outcome
estimation errors by at least 7.8%, 9.9% and 4.11% in the three datasets
respectively.

Table 2: Results of the data-efficient transfer learning setup for multi-step
outcome estimation. We report the mean $\pm$ standard deviation of Rooted Mean
Squared Errors (RMSEs $\downarrow$) over 5 runs. Bold: the best results.
Underline: the 2nd best results.

Dataset | Method | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---
Tumor growth | MSM | 1.0436$\pm$0.0671 | 0.5023$\pm$0.0588 | 0.7475$\pm$0.0829 | 0.9537$\pm$0.1060 | 1.1376$\pm$0.1233 | 1.3146$\pm$0.1338 | 0.9499$\pm$0.0915 | 7.8%
RMSN | 1.1839$\pm$0.0842 | 1.0912$\pm$0.0405 | 1.1215$\pm$0.0593 | 1.1538$\pm$0.0688 | 1.1728$\pm$0.0773 | 1.1740$\pm$0.0830 | 1.1495$\pm$0.0529 | 23.8%
CRN(ERM) | 1.2648$\pm$0.0689 | 1.1740$\pm$0.1015 | 1.1507$\pm$0.1016 | 1.1474$\pm$0.1070 | 1.1414$\pm$0.1041 | 1.1203$\pm$0.0906 | 1.1664$\pm$0.0786 | 24.9%
CRN | 1.4721$\pm$0.0620 | 8.4730$\pm$13.1786 | 13.1474$\pm$21.7125 | 17.2385$\pm$29.3990 | 21.5399$\pm$37.7228 | 25.3391$\pm$45.1447 | 14.5350$\pm$24.5243 | 94.0%
CT(ERM) | 0.8947$\pm$0.0668 | 0.8700$\pm$0.0857 | 0.9507$\pm$0.1309 | 1.0031$\pm$0.1502 | 1.0371$\pm$0.1545 | 1.0668$\pm$0.1565 | 0.9704$\pm$0.1098 | 9.7%
CT | 0.9545$\pm$0.0782 | 0.9494$\pm$0.1597 | 1.0225$\pm$0.1562 | 1.1062$\pm$0.1377 | 1.1455$\pm$0.1192 | 1.1562$\pm$0.0953 | 1.0557$\pm$0.1136 | 17.0%
G-Net | 1.0335$\pm$0.0622 | 1.0154$\pm$0.1100 | 1.1105$\pm$0.1476 | 1.1859$\pm$0.1620 | 1.2257$\pm$0.1693 | 1.2198$\pm$0.1508 | 1.1318$\pm$0.1118 | 22.6%
COSTAR | 0.8654$\pm$0.0328 | 0.7945$\pm$0.0532 | 0.8248$\pm$0.0751 | 0.8754$\pm$0.0987 | 0.9378$\pm$0.1176 | 0.9594$\pm$0.1062 | 0.8762$\pm$0.0720 | (-)
Semi- synthetic MIMIC-III | RMSN | 0.2100$\pm$0.0192 | 0.6084$\pm$0.1114 | 0.7745$\pm$0.1180 | 0.8908$\pm$0.1402 | 0.9776$\pm$0.1505 | 1.0440$\pm$0.1529 | 0.7509$\pm$0.1123 | 31.0%
CRN(ERM) | 0.1946$\pm$0.0158 | 0.4770$\pm$0.0808 | 0.5983$\pm$0.0923 | 0.6786$\pm$0.1004 | 0.7315$\pm$0.1047 | 0.7690$\pm$0.1070 | 0.5748$\pm$0.0823 | 9.9%
CRN | 0.2955$\pm$0.0256 | 0.5051$\pm$0.0748 | 0.6361$\pm$0.0786 | 0.7277$\pm$0.0783 | 0.7919$\pm$0.0764 | 0.8379$\pm$0.0759 | 0.6324$\pm$0.0656 | 18.1%
CT(ERM) | 0.2704$\pm$0.0631 | 0.5347$\pm$0.1061 | 0.6712$\pm$0.1252 | 0.7679$\pm$0.1433 | 0.8402$\pm$0.1607 | 0.8968$\pm$0.1784 | 0.6635$\pm$0.1279 | 21.9%
CT | 0.3105$\pm$0.0459 | 0.5840$\pm$0.0633 | 0.7414$\pm$0.0887 | 0.8530$\pm$0.1157 | 0.9348$\pm$0.1392 | 0.9974$\pm$0.1608 | 0.7368$\pm$0.0971 | 29.7%
G-Net | 0.3814$\pm$0.0556 | 0.6519$\pm$0.0856 | 0.8183$\pm$0.1122 | 0.9413$\pm$0.1365 | 1.0359$\pm$0.1592 | 1.1117$\pm$0.1795 | 0.8234$\pm$0.1191 | 37.1%
COSTAR | 0.2288$\pm$0.0229 | 0.4496$\pm$0.0877 | 0.5393$\pm$0.0962 | 0.5946$\pm$0.0990 | 0.6326$\pm$0.1013 | 0.6626$\pm$0.1026 | 0.5179$\pm$0.0844 | (-)
M5 | RMSN | 13.9705$\pm$0.3867 | 13.6233$\pm$0.8150 | 13.3291$\pm$1.2900 | 13.1984$\pm$1.3892 | 13.0889$\pm$1.2605 | 13.0108$\pm$1.1173 | 13.2501$\pm$1.1696 | 45.92%
CRN(ERM) | 6.3558$\pm$0.0594 | 7.0530$\pm$0.0433 | 7.3452$\pm$0.0447 | 7.5541$\pm$0.0392 | 7.7636$\pm$0.0450 | 7.9247$\pm$0.0561 | 7.5281$\pm$0.0447 | 4.82%
CRN | 6.2868$\pm$0.0471 | 7.0282$\pm$0.0482 | 7.3327$\pm$0.0610 | 7.5378$\pm$0.0521 | 7.7492$\pm$0.0586 | 7.9094$\pm$0.0676 | 7.5115$\pm$0.0572 | 4.60%
CT(ERM) | 6.1720$\pm$0.0354 | 6.9309$\pm$0.0571 | 7.2855$\pm$0.0889 | 7.5418$\pm$0.1191 | 7.7839$\pm$0.1283 | 7.9425$\pm$0.1430 | 7.4969$\pm$0.1058 | 4.42%
CT | 6.2041$\pm$0.0252 | 7.0022$\pm$0.0372 | 7.3675$\pm$0.0513 | 7.6394$\pm$0.0894 | 7.8932$\pm$0.1153 | 8.0701$\pm$0.1456 | 7.5945$\pm$0.0845 | 5.65%
G-Net | 6.7077$\pm$0.1006 | 7.0479$\pm$0.1069 | 7.3872$\pm$0.1349 | 7.6545$\pm$0.1596 | 7.9188$\pm$0.1800 | 8.1186$\pm$0.2058 | 7.4725$\pm$0.1461 | 4.11%
COSTAR | 6.3026$\pm$0.0519 | 6.8364$\pm$0.0560 | 7.1464$\pm$0.0674 | 7.3634$\pm$0.0619 | 7.6058$\pm$0.0640 | 7.7393$\pm$0.0637 | 7.1656$\pm$0.0592 | (-)

###  5.3 Supervised learning setup

Table 3 shows the performance in standard supervised learning setting, with
both train and test data from the source domain. Overall, COSTAR outperforms
other baselines in tumor growth and semi-synthetic MIMIC-III datasets. With
M5, COSTAR also shows comparable performance to the CT(ERM) with a 1.3%
relative difference.

Table 3: Results in standard supervised learning setting, with source and
target datasets coming from the same distribution for multi-step outcome
estimation. We report the mean +- standard deviation of Rooted Mean Squared
Errors (RMSEs) over 5 runs. Bold: the best results. Underline: the 2nd best
results.

Dataset | Method | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---
Tumor growth | MSM | 5.8368$\pm$0.6157 | 2.0400$\pm$0.6719 | 3.0385$\pm$0.9990 | 3.8701$\pm$1.2736 | 4.6173$\pm$1.5246 | 5.3823$\pm$1.7839 | 4.1308$\pm$1.1211 | 12.3%
RMSN | 4.8388$\pm$0.7770 | 5.4447$\pm$1.9202 | 5.9261$\pm$2.1096 | 5.9817$\pm$2.1270 | 5.8705$\pm$2.0544 | 5.5461$\pm$1.8865 | 5.6013$\pm$1.7727 | 35.4%
CRN(ERM) | 5.1601$\pm$0.5222 | 6.0784$\pm$2.3196 | 6.4721$\pm$2.4221 | 6.6142$\pm$2.4206 | 6.5648$\pm$2.3455 | 6.2939$\pm$2.1955 | 6.1972$\pm$2.0226 | 41.6%
CRN | 4.9781$\pm$0.3169 | 11.5230$\pm$9.9577 | 15.2778$\pm$15.9393 | 17.3213$\pm$19.6976 | 18.5698$\pm$22.0477 | 19.2921$\pm$23.5647 | 14.4937$\pm$15.2290 | 75.0%
CT(ERM) | 5.1286$\pm$1.3377 | 5.7262$\pm$2.7601 | 6.5085$\pm$2.9886 | 6.9248$\pm$3.0009 | 7.1971$\pm$2.9346 | 7.2369$\pm$2.7570 | 6.4537$\pm$2.5904 | 43.9%
CT | 6.5485$\pm$1.5221 | 7.5382$\pm$2.8528 | 7.9030$\pm$2.9569 | 7.9828$\pm$2.9332 | 7.8244$\pm$2.8075 | 7.4418$\pm$2.6103 | 7.5398$\pm$2.5976 | 52.0%
G-Net | 3.9371$\pm$0.4023 | 3.7697$\pm$1.1861 | 4.6054$\pm$1.4181 | 4.9730$\pm$1.4773 | 5.0491$\pm$1.4410 | 4.8745$\pm$1.3153 | 4.5348$\pm$1.1778 | 20.1%
COSTAR | 3.7403$\pm$0.3695 | 3.0067$\pm$0.9065 | 3.4619$\pm$1.1557 | 3.8501$\pm$1.3127 | 3.9160$\pm$1.3142 | 3.7525$\pm$1.1493 | 3.6212$\pm$1.0040 | (-)
Semi- synthetic MIMIC-III | RMSN | 0.2107$\pm$0.0261 | 0.5352$\pm$0.0842 | 0.6722$\pm$0.1096 | 0.7669$\pm$0.1203 | 0.8309$\pm$0.1280 | 0.8764$\pm$0.1331 | 0.6487$\pm$0.0976 | 21.7%
CRN(ERM) | 0.1951$\pm$0.0202 | 0.4426$\pm$0.0799 | 0.5530$\pm$0.0859 | 0.6113$\pm$0.0842 | 0.6478$\pm$0.0828 | 0.6708$\pm$0.0819 | 0.5201$\pm$0.0713 | 2.4%
CRN | 0.3276$\pm$0.0301 | 0.5234$\pm$0.0839 | 0.6531$\pm$0.0985 | 0.7234$\pm$0.0985 | 0.7618$\pm$0.0921 | 0.7825$\pm$0.0854 | 0.6286$\pm$0.0801 | 19.2%
CT(ERM) | 0.2130$\pm$0.0164 | 0.4426$\pm$0.0766 | 0.5495$\pm$0.0836 | 0.6191$\pm$0.0851 | 0.6669$\pm$0.0856 | 0.7010$\pm$0.0834 | 0.5320$\pm$0.0695 | 4.6%
CT | 0.2175$\pm$0.0178 | 0.4421$\pm$0.0757 | 0.5458$\pm$0.0854 | 0.6161$\pm$0.0925 | 0.6670$\pm$0.0993 | 0.7047$\pm$0.1040 | 0.5322$\pm$0.0765 | 4.6%
G-Net | 0.3418$\pm$0.0290 | 0.6015$\pm$0.0653 | 0.7542$\pm$0.0758 | 0.8620$\pm$0.0825 | 0.9429$\pm$0.0875 | 1.0035$\pm$0.0915 | 0.7510$\pm$0.0686 | 32.4%
COSTAR | 0.2286$\pm$0.0265 | 0.4417$\pm$0.0876 | 0.5288$\pm$0.0957 | 0.5825$\pm$0.0991 | 0.6190$\pm$0.1018 | 0.6458$\pm$0.1030 | 0.5077$\pm$0.0848 | (-)
M5 | RMSN | 35.7795$\pm$4.3603 | 33.2570$\pm$2.3870 | 33.4138$\pm$4.0678 | 33.4169$\pm$4.4289 | 33.3104$\pm$4.4017 | 33.3819$\pm$4.1602 | 33.7599$\pm$3.9379 | 52.2%
CRN(ERM) | 13.8445$\pm$0.1550 | 15.7926$\pm$0.1278 | 16.7071$\pm$0.2240 | 17.0887$\pm$0.1724 | 17.2709$\pm$0.0923 | 17.9759$\pm$0.0880 | 16.4466$\pm$0.1367 | 1.8%
CRN | 13.5907$\pm$0.0859 | 15.5242$\pm$0.0692 | 16.2694$\pm$0.1157 | 16.7355$\pm$0.0719 | 17.0095$\pm$0.0388 | 17.6874$\pm$0.0558 | 16.1361$\pm$0.0673 | -0.1%
CT(ERM) | 13.4887$\pm$0.1335 | 15.3397$\pm$0.2922 | 16.3415$\pm$0.4096 | 17.0545$\pm$0.5603 | 17.4828$\pm$0.5609 | 18.5832$\pm$0.4930 | 15.9414$\pm$0.3804 | -1.3%
CT | 13.6721$\pm$0.3574 | 15.9384$\pm$1.0910 | 17.2781$\pm$1.6049 | 18.1796$\pm$1.8134 | 18.8805$\pm$2.2159 | 19.2510$\pm$1.9642 | 16.7897$\pm$1.4144 | 3.8%
G-Net | 13.7187$\pm$0.0833 | 14.9851$\pm$0.1205 | 15.9578$\pm$0.1701 | 16.8278$\pm$0.2229 | 17.4833$\pm$0.3111 | 18.1665$\pm$0.3795 | 16.1898$\pm$0.2070 | 0.3%
COSTAR | 14.2556$\pm$0.1792 | 15.6151$\pm$0.2287 | 16.1743$\pm$0.2076 | 16.4791$\pm$0.1276 | 16.9037$\pm$0.1322 | 17.4379$\pm$0.1845 | 16.1443$\pm$0.1742 | (-)

###  5.4 Ablation studies

Table 4: Ablation studies for multi-step outcome estimation. We report the
mean $\pm$ standard deviation of Rooted Mean Squared Errors (RMSEs
$\downarrow$) over 5 runs. Bold: the best results.

Dataset | Component | Choice | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---|---
Semi- synthetic MIMIC-III | COSTAR | 0.2266$\pm$0.0249 | 0.4501$\pm$0.0893 | 0.5406$\pm$0.0987 | 0.5964$\pm$0.1020 | 0.6344$\pm$0.1040 | 0.6637$\pm$0.1052 | 0.5186$\pm$0.0869 | (-)
Encoder |  | w/ VT
---
0.4897$\pm$0.0888 | 0.6161$\pm$0.1139 | 0.6978$\pm$0.1200 | 0.7428$\pm$0.1217 | 0.7705$\pm$0.1196 | 0.7910$\pm$0.1182 | 0.6846$\pm$0.1130 | 24.2%
w/ CT | 0.3519$\pm$0.0584 | 0.4936$\pm$0.0897 | 0.5762$\pm$0.0967 | 0.6279$\pm$0.1017 | 0.6634$\pm$0.1045 | 0.6874$\pm$0.1029 | 0.5667$\pm$0.0919 | 8.5%
| FPE | w/ abs | 0.2981$\pm$0.0444 | 0.4679$\pm$0.0940 | 0.5561$\pm$0.1050 | 0.6091$\pm$0.1079 | 0.6446$\pm$0.1115 | 0.6694$\pm$0.1128 | 0.5409$\pm$0.0951 | 4.1%
| SSL Loss | none | 0.2998$\pm$0.0466 | 0.4718$\pm$0.0905 | 0.5579$\pm$0.1007 | 0.6117$\pm$0.1035 | 0.6460$\pm$0.1060 | 0.6680$\pm$0.1060 | 0.5425$\pm$0.0914 | 4.4%
| w/o comp | 0.2884$\pm$0.0421 | 0.4603$\pm$0.0898 | 0.5475$\pm$0.1032 | 0.6013$\pm$0.1077 | 0.6353$\pm$0.1079 | 0.6610$\pm$0.1075 | 0.5323$\pm$0.0926 | 2.6%
| SupL Loss | w/ uni | 0.2910$\pm$0.0355 | 0.4656$\pm$0.0873 | 0.5547$\pm$0.0998 | 0.6094$\pm$0.1039 | 0.6440$\pm$0.1059 | 0.6681$\pm$0.1070 | 0.5884$\pm$0.1006 | 11.9%
| w/ sq.inv. | 0.1968$\pm$0.0148 | 0.4456$\pm$0.0815 | 0.5415$\pm$0.0906 | 0.6021$\pm$0.0954 | 0.6431$\pm$0.0956 | 0.6761$\pm$0.0948 | 0.5175$\pm$0.0780 | -0.2%
| Decoder | w/ autoreg | 0.2049$\pm$0.0118 | 0.7036$\pm$0.1422 | 1.0234$\pm$0.2214 | 1.2023$\pm$0.2745 | 1.4692$\pm$0.3711 | 1.6577$\pm$0.4959 | 1.0435$\pm$0.2437 | 50.3%
M5 | COSTAR | 6.4054$\pm$0.0547 | 6.9328$\pm$0.0634 | 7.2428$\pm$0.0700 | 7.4585$\pm$0.0580 | 7.7012$\pm$0.0627 | 7.8278$\pm$0.0651 | 7.2614$\pm$0.0609 | (-)
Encoder |  | w/ VT
---
17.8226$\pm$4.7807 | 17.6769$\pm$4.5561 | 17.5937$\pm$4.4219 | 17.5279$\pm$4.2936 | 17.4113$\pm$4.1662 | 17.2706$\pm$4.0453 | 17.5505$\pm$4.3765 | 58.6%
w/ CT | 6.7386$\pm$0.2326 | 7.1911$\pm$0.2474 | 7.4549$\pm$0.2322 | 7.6524$\pm$0.2061 | 7.8488$\pm$0.2021 | 7.9589$\pm$0.2006 | 7.4741$\pm$0.2176 | 2.8%
| FPE | w/ abs | 6.4089$\pm$0.0693 | 6.9648$\pm$0.0617 | 7.2776$\pm$0.0528 | 7.4834$\pm$0.0430 | 7.7214$\pm$0.0421 | 7.8479$\pm$0.0370 | 7.2840$\pm$0.0486 | 0.3%
| SSL Loss | none | 6.4296$\pm$0.1193 | 6.9434$\pm$0.0796 | 7.2548$\pm$0.0699 | 7.4744$\pm$0.0753 | 7.7117$\pm$0.0728 | 7.8429$\pm$0.0817 | 7.2761$\pm$0.0827 | 0.2%
| w/o comp | 6.4637$\pm$0.0926 | 6.9847$\pm$0.0764 | 7.2934$\pm$0.0705 | 7.5093$\pm$0.0661 | 7.7497$\pm$0.0669 | 7.8748$\pm$0.0608 | 7.3126$\pm$0.0715 | 0.7%
| SupL Loss | w/ sq.inv. | 6.3425$\pm$0.0461 | 6.9760$\pm$0.0523 | 7.3170$\pm$0.0538 | 7.5366$\pm$0.0614 | 7.8015$\pm$0.0658 | 7.9439$\pm$0.0736 | 7.3196$\pm$0.0583 | 0.8%
| w/ inv | 6.3575$\pm$0.0473 | 6.9427$\pm$0.0381 | 7.2766$\pm$0.0422 | 7.4932$\pm$0.0447 | 7.7431$\pm$0.0531 | 7.8785$\pm$0.0628 | 7.2819$\pm$0.0464 | 0.3%
| Decoder | w/ autoreg | 6.3572$\pm$0.0621 | $>20$ | $>20$ | $>20$ | $>20$ | $>20$ | $>20$ | $>60\%$

We conduct ablation studies in the zero-shot transfer setup to validate the
design of COSTAR. We choose the feature-rich datasets: semi-synthetic MIMIC-
III and M5 since they contain complex dynamics and thus are more viable for
evaluating components capturing temporal and feature-wise interactions.

Encoder. To validate the impact of our proposed encoder architecture, we
replace the it with two variants: (1) Vanilla Transformer (w/ VT). A vanilla
transformer with temporal causal attention, which takes the history with all
features concatenated as multivariate time series input. (2) CT (w/ CT). The
encoder architecture proposed by Causal Transformer (Melnychuk et al., 2022)
that concatenates features grouped by covariates/treatments/outcomes into 3
subsets first, then applies self-attention/cross-attention among sequences
with each group of features/each pair of feature groups in an alternating way.
Rows “Backbone $|$ w/VT(w/CT)” in Table 4 demonstrate the superior performance
of our proposed encoder architecture. We observe that both methods (CT,
COSTAR) processing features respectively outperform VT that simply
concatenates all features, marking the importance of explicitly modeling
feature interactions. Moreover, the finer-grained modeling of feature
interactions between each pair of features in COSTAR further improves the
estimation performance compared to the coarser modeling of interactions
between feature subsets in CT.

Feature positional encoding (FPE). We replace the tree-based feature
positional encoding with its absolute variant (w/abs): each feature maps to a
separate learnable encoding vector. We observe that the tree-based positional
encoding has gains of 4.1% and 0.3% over the absolute variant in the two
datasets respectively.

Self-supervised loss (SSL). To validate the improvement brought by introducing
self-supervised learning as well as the choice of its training loss, we
compare COSTAR with two variants: (i) none. A model with the same architecture
as COSTAR but trained with factual estimation losses only; and (ii) w/o comp.
with vanilla MoCo v3 training loss in Eq. 10 for self-supervised learning.
Rows “SSL Loss $|$ none(w/o comp)” in Table 4 compare the estimation
performance of the aforementioned choices of self-supervised learning losses
and validate the effectiveness of our component-wise contrastive loss in self-
supervised learning.

Supervised loss (SupL). We consider different choices of the hyperparameter in
the supervised training loss of Eq. 14: (1) w/ uni: a uniform weight with each
$w_{i}=1/\tau$; (2) w/ inv: weights in proportion to the inverse of horizon
$w_{i}=\frac{1/i}{\sum_{j=1}^{\tau}1/j}$; (3) w/ sq.inv.: weights in
proportion to the inverse of the squared horizon $w_{i}\
=\frac{1/i^{2}}{\sum_{j=1}^{\tau}1/j^{2}}$. Both (2) and (3) are designed to
enhance the short-term outcome estimation performance. We select the weights
by validation error for each dataset (w/inv for Semi-synthetic MIMIC-III and
w/uni for M5), and compare it to the other two variants. While the relative
performance order varies across datasets, all variants can outperform the best
baseline results in Table 1.

Decoder. We validate the effectiveness of our non-autoregressive design of the
decoder and compare it with an autoregressive alternative (w/ autoreg) by
including the previous outcome in input features. While results in rows
”Decoder $|$ w/ autoreg” show good performance in very short horizons
($\tau=1$), multistep outcome estimation errors quickly diverges with the
horizon increasing.

###  5.5 Visualization of the Learned Representations

Figure 3: T-SNE visualization of learned representations in Semi-synthetic
MIMIC-III dataset.

Fig. 3 depicts the representations learned for the Semi-synthetic MIMIC-III
dataset after each of the 4 stages: (a) Pretrained: representations after the
self-supervised learning stage of source data. (b) Non-Cold-Start:
representations fine-tuned with factual outcome estimation loss of source
data. (c) Cold-Start: representations of target data when directly applying
the encoder trained in (b). (d) Transfer: representations of target data after
fine-tuned with small amount of target data. We use T-SNE to map each
representation to a 2D space and color each point with values of its upcoming
treatment and outcomes.

As shown in the first row, representations with different types of upcoming
treatments overlap, indicating that the learned representations after each
stage are balanced towards treatments. In the second and the third rows, we
observe clusters of representations corresponding to similar outcome values,
which indicates that the learned representations are informative about the
upcoming outcomes, even including the representations trained only with self-
supervised loss (column “Pretrained”). Such clustered structures also persist
when moving from the source domain data (column “Non-Cold-Start”) to the
target domain data (“Cold-Start”), showing that the learned representations
can generalize to cold-start cases.

##  6 Examples of counterfactual treatment outcome estimation

Figure 4: Examples of counterfactual treatment outcome estimation with semi-
synthetic MIMIC-III data in the zero-shot transfer setting. We plot one of the
two output dimensions for clarity. Each row lists the results of a
counterfactual treatment sequence, while each column shows the estimations of
one method across all treatment sequences tested. In each sub-figure, the
observed historical outcomes are plotted in black solid lines, and the ground
truth counterfactual outcomes in black dash lines. The blue solid lines show
the estimated outcomes.

Fig. 4 qualitatively compare the counterfactual outcome estimation performance
differences between COSTAR and baselines in the zero-shot transfer setting. We
randomly select a sequence from the observed data until time $t=4$ (x-axis),
then apply sequences of treatments sampled uniformly (i.e. no treatment bias)
and simulate the step-wise outcomes for 10 times as the ground truth. We
compare the ground truth of each simulation with all methods tested with semi-
synthetic MIMIC-III dataset. In Fig. 4 we find that the gaps between
estimations and ground truth outcomes are obvious in columns of baseline
results. Instead, they closely match each other in the estimation results (the
rightmost column) given by COSTAR, demonstrating its superior performance.

##  7 Conclusion

In this work, we propose a self-supervised learning framework - Counterfactual
Self-Supervised Transformer - to tackle the challenges associated with
accurately estimating treatment outcomes over time using observed history,
which is a crucial component in areas where randomized controlled trials
(RCTs) are not feasible. By integrating self-supervised learning and the
Transformer-based encoder combining temporal with feature-wise attention,
we’ve achieved notable advances in estimation accuracy and cross-domain
generalization performance.

## References

  * Bica et al. (2020) Ioana Bica, Ahmed M Alaa, James Jordon, and Mihaela van der Schaar.  Estimating counterfactual treatment outcomes over time through adversarially balanced representations.  In _International Conference on Learning Representations_ , 2020.  URL https://openreview.net/forum?id=BJg866NFvB.
  * Cao et al. (2023) Defu Cao, James Enouen, Yujing Wang, Xiangchen Song, Chuizheng Meng, Hao Niu, and Yan Liu.  Estimating treatment effects from irregular time series observations with hidden confounders.  _Proceedings of the AAAI Conference on Artificial Intelligence_ , 37(6):6897–6905, Jun. 2023.  doi: 10.1609/aaai.v37i6.25844.  URL https://ojs.aaai.org/index.php/AAAI/article/view/25844.
  * Chen et al. (2020) Ting Chen, Simon Kornblith, Mohammad Norouzi, and Geoffrey Hinton.  A simple framework for contrastive learning of visual representations.  In _International conference on machine learning_ , pp. 1597–1607. PMLR, 2020.
  * Chen* et al. (2021) Xinlei Chen*, Saining Xie*, and Kaiming He.  An empirical study of training self-supervised vision transformers.  _arXiv preprint arXiv:2104.02057_ , 2021.
  * Farahani et al. (2021) Abolfazl Farahani, Sahar Voghoei, Khaled Rasheed, and Hamid R Arabnia.  A brief review of domain adaptation.  _Advances in data science and information engineering: proceedings from ICDATA 2020 and IKE 2020_ , pp. 877–894, 2021.
  * Ganin et al. (2016) Yaroslav Ganin, Evgeniya Ustinova, Hana Ajakan, Pascal Germain, Hugo Larochelle, François Laviolette, Mario Marchand, and Victor Lempitsky.  Domain-adversarial training of neural networks.  _The journal of machine learning research_ , 17(1):2096–2030, 2016.
  * Geng et al. (2017) Changran Geng, Harald Paganetti, and Clemens Grassberger.  Prediction of treatment response for combined chemo-and radiation therapy for non-small cell lung cancer patients using a bio-mathematical model.  _Scientific reports_ , 7(1):13542, 2017.
  * Grill et al. (2020) Jean-Bastien Grill, Florian Strub, Florent Altché, Corentin Tallec, Pierre Richemond, Elena Buchatskaya, Carl Doersch, Bernardo Avila Pires, Zhaohan Guo, Mohammad Gheshlaghi Azar, et al.  Bootstrap your own latent-a new approach to self-supervised learning.  _Advances in neural information processing systems_ , 33:21271–21284, 2020.
  * HaoChen et al. (2022) Jeff Z. HaoChen, Colin Wei, Ananya Kumar, and Tengyu Ma.  Beyond separability: Analyzing the linear transferability of contrastive representations to related subpopulations.  In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho (eds.), _Advances in Neural Information Processing Systems_ , 2022.  URL https://openreview.net/forum?id=vmjckXzRXmh.
  * He et al. (2020) Kaiming He, Haoqi Fan, Yuxin Wu, Saining Xie, and Ross Girshick.  Momentum contrast for unsupervised visual representation learning.  In _Proceedings of the IEEE/CVF conference on computer vision and pattern recognition_ , pp. 9729–9738, 2020.
  * Huang et al. (2020) Zijie Huang, Yizhou Sun, and Wei Wang.  Learning continuous system dynamics from irregularly-sampled partial observations.  _Advances in Neural Information Processing Systems_ , 33:16177–16187, 2020.
  * Jiang et al. (2023) Song Jiang, Zijie Huang, Xiao Luo, and Yizhou Sun.  Cf-gode: Continuous-time causal inference for multi-agent dynamical systems.  _arXiv preprint arXiv:2306.11216_ , 2023.
  * Kidger et al. (2020) Patrick Kidger, James Morrill, James Foster, and Terry Lyons.  Neural controlled differential equations for irregular time series.  _Advances in Neural Information Processing Systems_ , 33:6696–6707, 2020.
  * Li et al. (2021) Rui Li, Stephanie Hu, Mingyu Lu, Yuria Utsumi, Prithwish Chakraborty, Daby M Sow, Piyush Madan, Jun Li, Mohamed Ghalwash, Zach Shahn, et al.  G-net: a recurrent network approach to g-computation for counterfactual prediction under a dynamic treatment regime.  In _Machine Learning for Health_ , pp. 282–299. PMLR, 2021.
  * Lim (2018) Bryan Lim.  Forecasting treatment responses over time using recurrent marginal structural networks.  _Advances in neural information processing systems_ , 31, 2018.
  * Mahar et al. (2021) Robert K Mahar, Myra B McGuinness, Bibhas Chakraborty, John B Carlin, Maarten J IJzerman, and Julie A Simpson.  A scoping review of studies using observational data to optimise dynamic treatment regimens.  _BMC medical research methodology_ , 21:1–13, 2021.
  * Makridakis et al. (2022) Spyros Makridakis, Evangelos Spiliotis, and Vassilios Assimakopoulos.  M5 accuracy competition: Results, findings, and conclusions.  _International Journal of Forecasting_ , 38(4):1346–1364, 2022.  ISSN 0169-2070.  doi: https://doi.org/10.1016/j.ijforecast.2021.11.013.  URL https://www.sciencedirect.com/science/article/pii/S0169207021001874.  Special Issue: M5 competition.
  * Melnychuk et al. (2022) Valentyn Melnychuk, Dennis Frauen, and Stefan Feuerriegel.  Causal transformer for estimating counterfactual outcomes.  In _International Conference on Machine Learning_ , pp. 15293–15329. PMLR, 2022.
  * Oord et al. (2018) Aaron van den Oord, Yazhe Li, and Oriol Vinyals.  Representation learning with contrastive predictive coding.  _arXiv preprint arXiv:1807.03748_ , 2018.
  * Park et al. (2020) Changhwa Park, Jonghyun Lee, Jaeyoon Yoo, Minhoe Hur, and Sungroh Yoon.  Joint contrastive learning for unsupervised domain adaptation.  _arXiv preprint arXiv:2006.10297_ , 2020.
  * Robins (1986) James Robins.  A new approach to causal inference in mortality studies with a sustained exposure period—application to control of the healthy worker survivor effect.  _Mathematical modelling_ , 7(9-12):1393–1512, 1986.
  * Robins & Hernan (2008) James Robins and Miguel Hernan.  Estimation of the causal effects of time-varying exposures.  _Chapman & Hall/CRC Handbooks of Modern Statistical Methods_, pp. 553–599, 2008.
  * Robins (1994) James M Robins.  Correcting for non-compliance in randomized trials using structural nested mean models.  _Communications in Statistics-Theory and methods_ , 23(8):2379–2412, 1994.
  * Robins et al. (2000) James M Robins, Miguel Angel Hernan, and Babette Brumback.  Marginal structural models and causal inference in epidemiology.  _Epidemiology_ , pp. 550–560, 2000.
  * Rubin (1978) Donald B Rubin.  Bayesian inference for causal effects: The role of randomization.  _The Annals of statistics_ , pp. 34–58, 1978.
  * Sagawa et al. (2021) Shiori Sagawa, Pang Wei Koh, Tony Lee, Irena Gao, Sang Michael Xie, Kendrick Shen, Ananya Kumar, Weihua Hu, Michihiro Yasunaga, Henrik Marklund, et al.  Extending the wilds benchmark for unsupervised adaptation.  _arXiv preprint arXiv:2112.05090_ , 2021.
  * Schulam & Saria (2017) Peter Schulam and Suchi Saria.  Reliable decision support using counterfactual models.  _Advances in neural information processing systems_ , 30, 2017.
  * Seedat et al. (2022) Nabeel Seedat, Fergus Imrie, Alexis Bellot, Zhaozhi Qian, and Mihaela van der Schaar.  Continuous-time modeling of counterfactual outcomes using neural controlled differential equations.  _arXiv preprint arXiv:2206.08311_ , 2022.
  * Shaw et al. (2018) Peter Shaw, Jakob Uszkoreit, and Ashish Vaswani.  Self-attention with relative position representations.  In _Proceedings of the 2018 Conference of the North American Chapter of the Association for Computational Linguistics: Human Language Technologies, Volume 2 (Short Papers)_ , pp. 464–468, 2018.
  * Shiv & Quirk (2019) Vighnesh Shiv and Chris Quirk.  Novel positional encodings to enable tree-based transformers.  _Advances in neural information processing systems_ , 32, 2019.
  * Soleimani et al. (2017) Hossein Soleimani, Adarsh Subbaswamy, and Suchi Saria.  Treatment-response models for counterfactual reasoning with continuous-time, continuous-valued interventions.  _arXiv preprint arXiv:1704.02038_ , 2017.
  * Splawa-Neyman et al. (1990) Jerzy Splawa-Neyman, Dorota M Dabrowska, and Terrence P Speed.  On the application of probability theory to agricultural experiments. essay on principles. section 9.  _Statistical Science_ , pp. 465–472, 1990.
  * Tan et al. (2023) Zhiquan Tan, Yifan Zhang, Jingqin Yang, and Yang Yuan.  Contrastive learning is spectral clustering on similarity graph, 2023.
  * Thota & Leontidis (2021) Mamatha Thota and Georgios Leontidis.  Contrastive domain adaptation.  In _Proceedings of the IEEE/CVF Conference on computer vision and pattern recognition_ , pp. 2209–2218, 2021.
  * Tonekaboni et al. (2021) Sana Tonekaboni, Danny Eytan, and Anna Goldenberg.  Unsupervised representation learning for time series with temporal neighborhood coding.  _arXiv preprint arXiv:2106.00750_ , 2021.
  * Tzeng et al. (2015) Eric Tzeng, Judy Hoffman, Trevor Darrell, and Kate Saenko.  Simultaneous deep transfer across domains and tasks.  In _Proceedings of the IEEE international conference on computer vision_ , pp. 4068–4076, 2015.
  * Vaswani et al. (2017) Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Łukasz Kaiser, and Illia Polosukhin.  Attention is all you need.  _Advances in neural information processing systems_ , 30, 2017.
  * Wang et al. (2023) Chao Wang, Xiaowei Shi, Shuai Xu, Zhe Wang, Zhiqiang Fan, Yan Feng, An You, and Yu Chen.  A multi-stage framework for online bonus allocation based on constrained user intent detection.  In _Proceedings of the 29th ACM SIGKDD Conference on Knowledge Discovery and Data Mining_ , pp. 5028–5038, 2023.
  * Wang et al. (2021) Zhiruo Wang, Haoyu Dong, Ran Jia, Jia Li, Zhiyi Fu, Shi Han, and Dongmei Zhang.  Tuta: Tree-based transformers for generally structured table pre-training.  In _Proceedings of the 27th ACM SIGKDD Conference on Knowledge Discovery & Data Mining_, pp. 1780–1790, 2021.
  * Woo et al. (2022) Gerald Woo, Chenghao Liu, Doyen Sahoo, Akshat Kumar, and Steven Hoi.  CoST: Contrastive learning of disentangled seasonal-trend representations for time series forecasting.  In _International Conference on Learning Representations_ , 2022.  URL https://openreview.net/forum?id=PilZY3omXV2.
  * Xu et al. (2016) Yanbo Xu, Yanxun Xu, and Suchi Saria.  A bayesian nonparametric approach for estimating individualized treatment-response curves.  In _Machine learning for healthcare conference_ , pp. 282–300. PMLR, 2016.
  * Ye et al. (2023) Hangting Ye, Zhining Liu, Wei Cao, Amir M Amiri, Jiang Bian, Yi Chang, Jon D Lurie, Jim Weinstein, and Tie-Yan Liu.  Web-based long-term spine treatment outcome forecasting.  In _Proceedings of the 29th ACM SIGKDD Conference on Knowledge Discovery and Data Mining_ , pp. 3082–3092, 2023.
  * Yue et al. (2022) Zhihan Yue, Yujing Wang, Juanyong Duan, Tianmeng Yang, Congrui Huang, Yunhai Tong, and Bixiong Xu.  Ts2vec: Towards universal representation of time series.  In _Proceedings of the AAAI Conference on Artificial Intelligence_ , volume 36, pp. 8980–8987, 2022.
  * Zhang et al. (2022a) Xiang Zhang, Ziyuan Zhao, Theodoros Tsiligkaridis, and Marinka Zitnik.  Self-supervised contrastive pre-training for time series via time-frequency consistency.  _Advances in Neural Information Processing Systems_ , 35:3988–4003, 2022a.
  * Zhang et al. (2022b) Yizhou Zhang, Defu Cao, and Yan Liu.  Counterfactual neural temporal point process for estimating causal influence of misinformation on social media.  _Advances in Neural Information Processing Systems_ , 35:10643–10655, 2022b.

##  Appendix A Identifiability assumptions

######  Assumption A.1 (Consistency).

The potential outcome of any treatment ${\bm{a}}_{t}$ is always the same as
the factual outcome when a subject is given the treatment ${\bm{a}}_{t}$:
${\bm{y}}_{t+1}[{\bm{a}}_{t}]={\bm{y}}_{t+1}$.

######  Assumption A.2 (Positivity).

If
$P(\bar{{\bm{A}}}_{t-1}=\bar{{\bm{a}}}_{t-1},\bar{{\bm{X}}}_{t}=\bar{{\bm{x}}}_{t})\neq
0$, then
$P({\bm{A}}_{t}={\bm{a}}_{t}|\bar{{\bm{A}}}_{t-1}=\bar{{\bm{a}}}_{t-1},\bar{{\bm{X}}}_{t}=\bar{{\bm{x}}}_{t})>0$
for any $\bar{{\bm{a}}}_{t}$.

######  Assumption A.3 (Sequential strong ignorability).

${\bm{Y}}_{t+1}[{\bm{a}}_{t}]\perp\\!\\!\\!\perp{\bm{A}}_{t}|\bar{{\bm{A}}}_{t-1},\bar{{\bm{X}}}_{t},\forall{\bm{a}}_{t},t$.

##  Appendix B Why use temporally causal attention in COSTAR?

The embeddings of time-varying features
${\bm{Z}}^{S,(L)}\in\mathbb{R}^{t\times(d_{X}+d_{A}+d_{Y})\times
d_{\text{model}}}$ from the final layer is re-organized to stepwise
representations $({\bm{Z}}_{1},{\bm{Z}}_{2},\dots,{\bm{Z}}_{t})$. Each
${\bm{Z}}_{t^{\prime}}=\\{{\bm{z}}_{t^{\prime}}^{i}\in\mathbb{R}^{d_{\text{model}}}\\}_{i=1}^{d_{X}+d_{A}+d_{Y}}$
is further aggregated to
${\bm{z}}^{X}_{t^{\prime}},{\bm{z}}^{A}_{t^{\prime}},{\bm{z}}^{Y}_{t^{\prime}},{\bm{z}}_{t^{\prime}}$
in Equation 1. When the encoder satisfies the temporal causality (i.e.
${\bm{Z}}_{t^{\prime}}$ only depends on $\bar{{\bm{H}}}_{t^{\prime}}$), they
can be seen as a sequence of representations for the observed history
$\bar{{\bm{H}}}_{1},\bar{{\bm{H}}}_{2},\dots,\bar{{\bm{H}}}_{t}$ truncated at
each time step.

When we feed the encoder with factual data in training stages, a major
advantage of encoders satisfying temporal causality is that we can estimate
the outcomes and evaluate the factual estimation losses in every time step of
the input sequence at a single forward pass. Evaluating of counterfactual
data, the encoder only needs to keep ${\bm{Z}}_{t}$ representing the entire
observed history, conditioning on which the predictor rolls out outcome
estimations given counterfactual treatments.

In contrast, feeding the entire history in one pass for training is error-
prone for architectures and can violate the temporal causality (e.g.
transformers with fully temporal attention or frequency-based methods), since
it leaks information of future steps into the representations in previous
steps. Predictors trained with such representations converge quickly to a
trivial model that simply copies future steps as estimations. When we evaluate
counterfactual data where future counterfactual outcomes are no longer
available in input, the performance degenerates. As a result, we have to
explicitly unroll the observed sequence to $t$ truncated sequences and run $t$
forward passes to get the representations and factual errors in every step
when training non-temporally-causal models. This leads to a $\times T$
increase in training time when the batch size remains unchanged due to
hardware restrictions, where $T$ is the maximum length of sequences in
training data. We have $T\geq 50$ in our experiments and find that none of the
architecture violating temporal causality can finish training in a reasonable
time.

##  Appendix C Proof of generalization bound of contrastive learning in
counterfactual outcome estimation

###  C.1 Unsupervised domain adaptation view of counterfactual outcome
estimation

For the given observed history $\bar{{\bm{H}}}_{t}$, treatment sequence
$\bar{{\bm{a}}}_{t:t+\tau-1}$ to apply, and the outcome
${\bm{y}}_{t+\tau}[{\bm{a}}_{t:t+\tau-1}](\bar{{\bm{H}}}_{t})$ to estimate, we
slightly abuse scalar/vector/matrix notations and denote them as $h,a,y[a](h)$
for simplicity. With discrete treatments222For a sequence of discrete
treatments, we can always map it to a single discrete variable with a proper
encoding., we notice that the counterfactual outcome prediction of each type
of treatments from factual data can be viewed as an unsupervised domain
adaptation (UDA) problem:

For a treatment type $a$, we aim at finding a function $f_{a}(h)$ that
specifically estimates $\mathbb{E}(y[a](h))$. Any sample $(h_{i},a_{i},y_{i})$
from the observed dataset $\mathcal{D}_{tr}$ can be categorized into one of
(1) Labeled subset $\mathcal{S}_{a}=\\{(h_{i},a_{i},y_{i})\\}_{a_{i}=a}$ and
(2) Unlabeled subset $\mathcal{T}_{a}=\\{(h_{i},a_{i},y_{i})\\}_{a_{i}\neq
a}$. According to Assumption A.1, for any
$(h_{i},a_{i},y_{i})\in\mathcal{S}_{a}$, we have $y_{i}=y[a](h_{i})$ and thus
$y_{i}$ is a label of $f_{a}(h)$. In contrast, for any
$(h_{i},a_{i},y_{i})\in\mathcal{T}_{a}$, $y_{i}=y[a_{i}](h_{i})$ and
$a_{i}\neq a$, resulting in that $y_{i}$ is no longer a valid label of
$f_{a}(h)$. For simplicity, we omit the treatment symbol $a$ as well as the
$y_{i}$ in $\mathcal{T}_{a}$: $\mathcal{S}_{a}=\\{(h_{i},y_{i})\\}_{a_{i}=a}$,
$\mathcal{T}_{a}=\\{h_{i}\\}_{a_{i}\neq a}$.

Considering the existence of treatment bias, there exists at least a
$a^{\prime}\neq a$ satisfying $P(h|a)\neq P(h|a^{\prime})$, which potentially
leads to $P_{\mathcal{T}_{a}}(h)=P_{\mathcal{D}_{tr}}(h_{i}|a_{i}\neq a)\neq
P_{\mathcal{D}_{tr}}(h_{i}|a_{i}=a)=P_{\mathcal{S}_{a}}(h)$. In counterfactual
estimation, we aim at minimizing the estimation error without treatment bias:

$\mathcal{E}_{\mathcal{D}}(f_{a})=\mathbb{E}_{h\sim P_{\mathcal{D}}(h),y\sim
P_{\mathcal{D}}(y[a]|h)}\ell(f_{a}(h),y).$ (16)

Here, we focus on the case where no covariate and concept shifts happen across
datasets 333The distribution shift between $P_{\mathcal{D}}(h)$ v.s.
$P_{\mathcal{D}_{tr}}(h)$ (covariate shift) and $P_{\mathcal{D}}(y[a]|h)$ v.s.
$P_{\mathcal{D}_{tr}}(y[a]|h)$ (concept shift) can be viewed as the general
covariate shift/concept shift and fit into existing theories of domain
adaptation/generalization (Farahani et al., 2021).:
$P_{\mathcal{D}}(h)=P_{\mathcal{D}_{tr}}(h)=\sum_{a^{\prime}}P(a^{\prime})P_{\mathcal{D}_{tr}}(h|a^{\prime})$,
$P_{\mathcal{D}}(y[a]|h)=P_{\mathcal{D}_{tr}}(y[a]|h)=P_{\mathcal{S}_{a}}(y[a]|h)=P_{\mathcal{T}_{a}}(y[a]|h)$.
Then Eq. 16 becomes:

$\displaystyle\mathcal{E}_{\mathcal{D}}(f_{a})=$
$\displaystyle\sum_{a^{\prime}}P(a^{\prime})\mathbb{E}_{h\sim
P_{\mathcal{D}_{tr}}(h|a^{\prime}),y\sim
P_{\mathcal{D}_{tr}}(y[a]|h)}\ell(f_{a}(h),y)$ (17) $\displaystyle=$
$\displaystyle P(a)\underbrace{\mathbb{E}_{h\sim P_{\mathcal{S}_{a}}(h),y\sim
P_{\mathcal{D}_{tr}}(y[a]|h)}\ell(f_{a}(h),y)}_{\mathcal{E}_{\mathcal{S}_{a}}(f_{a})}+$
$\displaystyle(1-P(a))\underbrace{\mathbb{E}_{h\sim
P_{\mathcal{T}_{a}}(h),y\sim
P_{\mathcal{D}_{tr}}(y[a]|h)}\ell(f_{a}(h),y)}_{\mathcal{E}_{\mathcal{T}_{a}}(f_{a})}.$

So far we can see that for the task of finding an outcome estimator for
treatment type $a$, the counterfactual estimation error is bounded by
estimation errors on both $\mathcal{S}_{a}$ (denoted as
$\mathcal{E}_{\mathcal{S}_{a}}(f_{a})$) and $\mathcal{T}_{a}$ (denoted as
$\mathcal{E}_{\mathcal{T}_{a}}(f_{a})$). Per our analysis above,
$\mathcal{S}_{a}$ and $\mathcal{T}_{a}$ corresponds to the labeled source
domain data and the unlabeled target domain data in UDA problems, where
potential distribution shifts exist between $\mathcal{S}_{a}$ and
$\mathcal{T}_{a}$ due to treatment bias. While
$\mathcal{E}_{\mathcal{S}_{a}}(f_{a})$ can be optimized with supervised
learning using factual data in $\mathcal{S}_{a}$, we cannot directly optimize
$\mathcal{E}_{\mathcal{T}_{a}}(f_{a})$ with labeled data directly.

Recent works (Thota & Leontidis, 2021; Sagawa et al., 2021; Park et al., 2020;
Wang et al., 2021) show that contrastive learning, as an effective self-
supervised representation learning method, demonstrates strong transferability
in UDA and leads to simple state-of-the-art algorithms. Considering the close
connection between counterfactual outcome estimation and UDA, we develop our
model COSTAR based on contrastive learning and our analysis of the
counterfactual estimation error bound from the recent work (HaoChen et al.,
2022), where the authors provide theoretical analysis of the transferability
of contrastive learning in UDA.

###  C.2 Preliminaries

#### Positive pairs.

Pairs of semantically related/similar data samples are positive pairs in
contrastive learning. In contrastive learning, positive pairs are commonly
generated by applying randomized transformation on the same input (He et al.,
2020; Woo et al., 2022).

For exposition simplicity, we assume the set of factual observed history
$\mathcal{H}$ is a finite but large dataset of size $N$. We use $P_{+}$ to
denote the distribution of positive pairs. $P_{+}$ satisfies
$P_{+}(h,h^{\prime})=P_{+}(h^{\prime},h),~{}\forall~{}h,h^{\prime}\in\mathcal{H}$.
$P_{\mathcal{H}}$ denotes the marginal distribution of $P_{+}$:
$P_{\mathcal{H}}(h)=\sum_{h^{\prime}\in\mathcal{H}}P_{+}(h,h^{\prime})$.

#### Positive-pair graph.

Following the definition in (HaoChen et al., 2022), we introduce the postive-
pair graph as a weighted undirected graph $G(\mathcal{H},w)$ with the vertex
set $\mathcal{H}$ and the edge weight $w(h,h^{\prime})=P_{+}(h,h^{\prime})$.
$w(h)=P_{\mathcal{H}}(h)=\sum_{h^{\prime}\in\mathcal{H}}w(h,h^{\prime})$. For
any vertex subset $A$, $w(A)=\sum_{h\in A}w(h)$. For any vertex subsets $A,B$,
$w(A,B)=\sum_{h\in A,h^{\prime}\in B}w(h,h^{\prime})$. For any vertex $h$ and
vetex subset $B$, $w(h,B)=w(\\{h\\},B)$.

#### Generalized spectral contrastive loss.

Let $r:\mathcal{H}\rightarrow\mathbb{R}^{k}$ be a mapping from the input data
to $k$-dimensional features. For the convenience of proof, we consider the
(generalized) spectral contrastive loss proposed in (HaoChen et al., 2022):

$\mathcal{L}_{\sigma}(r)=\mathbb{E}_{(h,h^{+})\sim
P_{+}}\left[\left\lVert{r(h)-r(h^{+})}\right\rVert_{2}^{2}\right]+\sigma\cdot
R(r),$ (18)

where the regularizer is defined as $R(r)=\left\lVert\mathbb{E}_{h\in
P_{\mathcal{H}}}[r(h)r(h)^{T}]-I_{k\times k}\right\rVert_{F}^{2}$ and
$I_{k\times k}$ is the $k$-dimensional identity matrix. Notice that the
InfoNCE loss is more commonly used in empirical study (He et al., 2020; Chen
et al., 2020; Chen* et al., 2021) instead of the spectral contrastive loss,
and their equivalence is still an open problem with some preliminary results
(Tan et al., 2023).

###  C.3 Definitions and assumptions

We reiterate the following definitions and assumptions in (HaoChen et al.,
2022) for self-containment:

######  Definition C.1 (Expansion).

Let $A,B$ be two disjoint subsets of $\mathcal{H}$. We denote the expansion,
max-expansion and min-expansion from $A$ to $B$ as follows:

$\phi(A,B)=\frac{w(A,B)}{w(A)},\quad\bar{\phi}(A,B)=\max_{h\in
A}\frac{w(h,B)}{w(h)},\quad\underline{\phi}(A,B)=\min_{h\in
A}\frac{w(h,B)}{w(h)}.$ (19)

######  Assumption C.2 (Cross-cluster connections).

For some $\alpha\in(0,1)$, we assume that vertices of the positive-pair graph
$G(\mathcal{H},w)$ can be partitioned into $m$ disjoint clusters
$C_{1},\dots,C_{m}$ such that for any $i\in[m]$,

$\bar{\phi}(C_{i},\mathcal{H}\backslash C_{i})\leq\alpha.$ (20)

######  Assumption C.3 (Intra-cluster conductance).

For all $i\in[m]$, assume the conductance of the subgraph restricted to
$C_{i}$ is large, i.e., every subset A of $C_{i}$ with at most half the size
of $C_{i}$ expands to the rest:

$\forall~{}A\subset C_{i}~{}\text{satisfying}~{}w(A)\leq
w(C_{i})/2,~{}\phi(A,C_{i}\backslash A)\geq\gamma.$ (21)

######  Assumption C.4 (Relative expansion).

Let $S$ and $T$ be two disjoint subsets of $\mathcal{H}$, each is formed by
$r$ clusters among $C_{1},C_{2},\dots,C_{m}$ for $r\leq m/2$. Let
$\rho=\min_{i\in[r]}\underline{\phi}(T_{i},S_{i})$ be the minimum min-
expansions from $T_{i}$ to $S_{i}$. For some sufficiently large universal
constant $c$, we assume that $\rho\geq c\cdot\alpha^{2}$ and that

$\rho=\min_{i\in[r]}\underline{\phi}(T_{i},S_{i})\geq c\cdot\max_{i\neq
j}\bar{\phi}(T_{i},S_{j}).$ (22)

###  C.4 Proof of Theorem 4.1

We adapt the preconditioned featurer averaging classifier in (HaoChen et al.,
2022) for regression in our proof:

Algorithm 1 Preconditioned feature averaging (PFA).

Require: Pretrained representation extractor $r$, unlabeled data
$P_{\mathcal{H}}$, source domain labeled data $P_{S}$, target domain test data
$\tilde{h}$, integer $t\in{\mathbb{Z}}^{+}$, outcome discretization
granularity $\epsilon$.

1:Compute the preconditioner matrix $\Sigma=\mathbb{E}_{h\in
P_{\mathcal{H}}}[r(h)r(h)^{T}]$.

2:for every outcome value $y_{i}$ corresponding to the cluster $C_{i},i\in[r]$
do

3: Compute the mean feature of outcome $y_{i}$: $b_{i}=\mathbb{E}_{(h,y)\sim
P_{S}}[\mathbbm{1}[\left\lVert{y-y_{i}}\right\rVert_{2}\leq\epsilon]\cdot
r(h)]$.

4:end for

5:return prediction $y_{i^{*}}$,
$i^{*}=\operatorname*{arg\,max}_{i\in[r]}\langle r(h),\sum^{t-1}b_{i}\rangle$.

For any PFA regressor $f$ constructed with Alg. 1, we can transform it to a
corresponding classifier by defining its 0-1 classification error on the
target domain $T$ as:

$\mathcal{E}_{T}^{01}(f)=\mathbb{E}_{(h,y)\sim
P_{T}}[\mathbbm{1}[\left\lVert{y-f(h)}\right\rVert_{2}>\epsilon]].$ (23)

We can directly apply the main result in (HaoChen et al., 2022) and get an
upper bound of the 0-1 error on the target domain:

######  Theorem C.5 (Upper bound of 0-1 error on the target domain (HaoChen et
al., 2022)).

Suppose that Assumption C.2, Assumption C.3, and Assumption C.4 holds for the
set of observed history $\mathcal{H}$ and its positive-pair graph
$G(\mathcal{H},w)$, and the representation dimension $k\geq 2m$. Let $r$ be a
minimizer of the generalized spectral contrastive loss and the regression head
$f$ be constructed in Alg. 1. We have

$\mathcal{E}^{01}_{T}(f)\lesssim\frac{r}{\alpha^{2}\gamma^{4}}\cdot\exp(-\Omega(\frac{\rho\gamma^{2}}{\alpha^{2}})).$
(24)

######  Lemma C.6 (Relation between the L2 regression error and 0-1
classification error).

Suppose that both $\left\lVert{f(h)}\right\rVert_{2}\leq B$ and
$\left\lVert{y}\right\rVert_{2}\leq B,~{}\epsilon<2B$. The L2 regression error
$\mathcal{E}_{T}(f)$ of the PFA regressor on the target domain $T$ is bounded
by $\mathcal{E}_{T}^{01}(f)$ as:

$\mathcal{E}_{T}(f)\leq\epsilon^{2}+(4B^{2}-\epsilon^{2})\mathcal{E}_{T}^{01}(f).$
(25)

###### Proof.

$\displaystyle\mathcal{E}_{T}(f)=$ $\displaystyle\mathbb{E}_{(h,y)\in
P_{T}}\left\lVert{y-f(h)}\right\rVert_{2}^{2}$ $\displaystyle\leq$
$\displaystyle\sum_{(h,y)\in
T}P(h,y)\left[\mathbbm{1}[\left\lVert{y-f(h)}\right\rVert_{2}>\epsilon]\left\lVert{y-f(h)}\right\rVert_{2}^{2}\right]$
$\displaystyle+\sum_{(h,y)\in
T}P(h,y)(1-\mathbbm{1}[\left\lVert{y-f(h)}\right\rVert_{2}>\epsilon])\epsilon^{2}$
$\displaystyle\leq$ $\displaystyle\sum_{(h,y)\in
T}P(h,y)\left[\mathbbm{1}[\left\lVert{y-f(h)}\right\rVert_{2}>\epsilon]4B^{2}\right]$
$\displaystyle+\sum_{(h,y)\in
T}P(h,y)(1-\mathbbm{1}[\left\lVert{y-f(h)}\right\rVert_{2}>\epsilon])\epsilon^{2}$
$\displaystyle=$
$\displaystyle\epsilon^{2}+(4B^{2}-\epsilon^{2})\mathbb{E}_{(h,y)\in
P_{T}}\mathbbm{1}[\left\lVert{y-f(h)}\right\rVert_{2}>\epsilon]$
$\displaystyle=$
$\displaystyle\epsilon^{2}+(4B^{2}-\epsilon^{2})\mathcal{E}_{T}^{01}(f).$

∎

Lemma C.6 connects the L2 error and the 0-1 error. Combining Eq. 17, Eq. 24,
Eq. 25, we immediately get Theorem 4.1.

##  Appendix D Dataset description

Table 5: Statistics of datasets. Dataset | Domain | Property | Seq Length | Train/Validation/Test Seq Num
---|---|---|---|---
Tumor growth | source | $\gamma=10$ | 60 | 10000/1000/1000
target | $\gamma=0$ | 60 | 100/1000/1000
Semi-synthetic MIMIC-III | source | age in [20,45] | 99 | 3704/926/926
target | age$\geq$85 | 99 | 138/347/1737
M5 | source | food items | 50 | 39606/7048/7048
target | household items | 50 | 3623/3512/18005

We summarize the statistics and the way of introducing feature distribution
shifts in Table 5.

#### Tumor growth.

We refer readers to (Bica et al., 2020; Melnychuk et al., 2022) for the
complete descriptions of the pharmacokinetic-pharmacodynamic (PK-PD) model.
Here we focus on how we introduce distribution shifts by adjusting the
treatment bias coefficient $\gamma$.

The volume of tumor after $t$ days of diagnosis is:

$V(t+1)=(1+\rho\log(\frac{K}{V(t)})-\beta_{c}C(t)-(\alpha_{r}d(t)+\beta_{r}d(t)^{2})+e_{t})V(t),$
(26)

where $K,\rho,\beta_{c},\alpha_{r},\beta_{r}$ are parameters sampled from the
prior distributions defined in (Geng et al., 2017).
$e_{t}\sim\mathcal{N}(0,0.01^{2})$ is the noise term.

PK-PD model constructs time-varying confounding by connecting the probability
of assigning chemotherapy and radiotherapy with the outcome - tumor diameter:

$p_{c}(t)=\sigma(\frac{\gamma_{c}}{D_{\text{max}}}(\bar{D}(t)-\delta_{c})),\quad
p_{r}(t)=\sigma(\frac{\gamma_{r}}{D_{\text{max}}}(\bar{D}(t)-\delta_{r})).$
(27)

$\bar{D}(t)$ is the mean tumor diameter in the past 15 days and
$D_{\text{max}}=13$. $\sigma$ is the sigmoid function.
$\delta_{c}=\delta_{r}=D_{\text{max}}/2$. $\gamma_{c}$ and $\gamma_{r}$
controls the importance of tumor diameter history on treatment assignment,
thus control the strength of time-dependent confounding.

In Tumor growth dataset, we set $\gamma_{c}=\gamma_{r}=\gamma=10$ to generate
data in the source domain, and $\gamma_{c}=\gamma_{r}=\gamma=0$ for the target
domain. As a result, both treatment bias and the data distribution of history
differs between source and target domains.

#### Semi-synthetic MIMIC-III.

We split the semi-synthetic MIMIC-III dataset introduced in (Melnychuk et al.,
2022) by ages of patients to the source/target domain. More specifically, we
generate simulation data from patients with ages falling in $[20,45]$ as the
source domain data and simulation based on patients with ages over 85 as the
target domain data.

#### M5.

We adapt the M5 forecasting dataset
(https://www.kaggle.com/competitions/m5-forecasting-accuracy) for treatment
effect estimation over time. In M5, we select the item pricing as treatment,
its sales as outcome and all other features as covariates. We aggregate the
item sales by week to reduce the sequence length to the same level as the
other two datasets for the convenience of evaluation. We also discretize the
continuous pricing by mapping $(p_{t,i}-p_{0,i})/p_{0,i}$ to buckets divided
by its 20-quantiles, where $p_{t,i},p_{0,i}$ are the prices of item $i$ at
time $t$ and at its initial sale.

To introduce the feature distribution shift, we select 5000 items in the food
category as the source domain data and another 5000 items in the household
category as the target domain data.

##  Appendix E Baselines

#### Baseline implementation.

We reuse the implementation in (Melnychuk et al., 2022) for evaluating all the
baselines, including: MSM (Robins et al., 2000), RMSN (Lim, 2018), CRN (Bica
et al., 2020), G-Net (Li et al., 2021), and Causal Transformer (CT) (Melnychuk
et al., 2022).

#### Hyperparameter tuning.

For all baselines, we follow the ranges of hyperparameter tuning in (Melnychuk
et al., 2022) and select the hyperparameters with the lowest factual outcome
estimation error on the validation set from the source domain. For each method
and each dataset, the same set of hyperparameters are used in the zero-shot
transfer/data-efficient transfer/standard supervised learning settings. The
detailed hyperparameters used for baselines and COSTAR are listed in the
configuration files in our code repository. Here we list the main
hyperparameters for reference.

MSM. There is no tuneable hyperparameter in MSM.

RMSN. We list the hyperparameters of RMSN in Table 6.

Table 6: RMSN hyperparameters. |  | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---|---
Propensity Treatment | RNN Hidden Units | 8 | 6 | 44
| Dropout | 0.5 | 0.1 | 0.4
| Layer Num | 1 | 2 | 1
| Max Gradient Norm | 1.0 | 0.5 | 2.0
| Batch Size | 128 | 256 | 128
| Learning Rate | 0.01 | 0.01 | 0.001
Propensity History | RNN Hidden Units | 24 | 74 | 92
| Dropout | 0.1 | 0.5 | 0.5
| Layer Num | 1 | 2 | 2
| Max Gradient Norm | 2.0 | 1.0 | 0.5
| Batch Size | 128 | 64 | 128
| Learning Rate | 0.01 | 0.001 | 0.01
Encoder | RNN Hidden Units | 24 | 74 | 46
| Dropout | 0.1 | 0.1 | 0.1
| Layer Num | 1 | 1 | 2
| Max Gradient Norm | 0.5 | 0.5 | 0.5
| Batch Size | 64 | 1024 | 128
| Learning Rate | 0.01 | 0.001 | 0.0001
Decoder | RNN Hidden Units | 48 | 196 | 45
| Dropout | 0.1 | 0.1 | 0.1
| Layer Num | 1 | 1 | 1
| Max Gradient Norm | 0.5 | 0.5 | 4.0
| Batch Size | 256 | 1024 | 256
| Learning Rate | 0.0001 | 0.0001 | 0.0001

CRN(ERM). See Table 7.

Table 7: CRN(ERM) hyperparameters. |  | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---|---
Encoder | RNN Hidden Units | 24 | 74 | 46
| Balancing Representation Size | 18 | 74 | 46
| FC Hidden Units | 18 | 37 | 46
| Layer Num | 1 | 1 | 2
| Dropout | 0.1 | 0.1 | 0.1
| Batch Size | 256 | 64 | 128
| Learning Rate | 0.01 | 0.001 | 0.001
Decoder | RNN Hidden Units | 18 | 74 | 46
| Balancing Representation Size | 6 | 98 | 90
| FC Hidden Units | 6 | 98 | 22
| Layer Num | 1 | 2 | 1
| Dropout | 0.1 | 0.1 | 0.1
| Batch Size | 256 | 256 | 256
| Learning Rate | 0.001 | 0.0001 | 0.0001

CRN. See Table 8.

Table 8: CRN hyperparameters. |  | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---|---
Encoder | RNN Hidden Units | 18 | 74 | 46
| Balancing Representation Size | 3 | 74 | 46
| FC Hidden Units | 12 | 37 | 46
| Layer Num | 1 | 1 | 2
| Dropout | 0.2 | 0.1 | 0.1
| Batch Size | 256 | 64 | 128
| Learning Rate | 0.001 | 0.001 | 0.001
Decoder | RNN Hidden Units | 3 | 74 | 46
| Balancing Representation Size | 3 | 98 | 90
| FC Hidden Units | 3 | 98 | 22
| Layer Num | 1 | 2 | 1
| Dropout | 0.2 | 0.1 | 0.1
| Batch Size | 256 | 256 | 256
| Learning Rate | 0.001 | 0.0001 | 0.0001

CT(ERM). See Table 9.

Table 9: CT(ERM) hyperparameters. | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---
Transformer Hidden Units | 12 | 24 | 24
Balancing Representation Size | 2 | 88 | 94
FC Hidden Units | 12 | 44 | 47
Layer Num | 1 | 1 | 2
Head Num | 2 | 3 | 2
Max Relative Position | 15 | 20 | 30
Dropout | 0.1 | 0.1 | 0.1
Batch Size | 64 | 64 | 64
Learning Rate | 0.001 | 0.01 | 0.001

CT. See Table 10.

Table 10: CT hyperparameters. | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---
Transformer Hidden Units | 16 | 24 | 24
Balancing Representation Size | 16 | 88 | 94
FC Hidden Units | 16 | 44 | 47
Layer Num | 1 | 1 | 2
Head Num | 2 | 3 | 2
Max Relative Position | 15 | 20 | 30
Dropout | 0.2 | 0.1 | 0.1
Batch Size | 64 | 64 | 64
Learning Rate | 0.001 | 0.01 | 0.001

G-Net. See Table 11.

Table 11: G-Net hyperparameters. | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---
RNN Hidden Units | 24 | 148 | 144
FC Hidden Units | 48 | 74 | 72
Dropout | 0.1 | 0.1 | 0.1
Layer Num | 1 | 1 | 2
Batch Size | 128 | 256 | 256
Learning Rate | 0.001 | 0.01 | 0.001

COSTAR. See Table 12.

Table 12: COSTAR hyperparameters. |  | Tumor growth | Semi-synthetic MIMIC-III | M5
---|---|---|---|---
Encoder | Transformer Hidden Units | 24 | 36 | 36
| Encoder Momentum | 0.99 | 0.99 | 0.99
| Temperature | 1.0 | 1.0 | 1.0
| Layer Num | 1 | 1 | 2
| Head Num | 2 | 3 | 2
| Dropout | 0.1 | 0.1 | 0.1
| Batch Size | 64 | 64 | 64
| Learning Rate | 0.001 | 0.001 | 0.001
Decoder | Hidden Units | 128 | 128 | 128
| Batch Size | 32 | 32 | 32
| Learning Rate | 0.001 | 0.001 | 0.001

#### Comparison of numbers of model parameters.

In Table 13 we list the number of trainable parameters in each baseline as
well as COSTAR in the experiments of each dataset.

Table 13: Number of trainable parameters. #trainable params | Tumor growth | semi-synthetic MIMIC-III | M5
---|---|---|---
MSM |  $<$1K | (-) | (-)
RMSN | 18.8K | 387K | 213K
CRN(ERM) | 6.5K | 164K | 78K
CRN | 2.3K | 164K | 78K
CT(ERM) | 5.2K | 45K | 80.3K
CT | 9.4K | 45K | 80.3K
G-Net | 3.4K | 151K | 323K
COSTAR | 20.7K | 43.6K | 77.5K

##  Appendix F A detailed look into CRN results with large variances.

In Table 1, 2, 3 we notice that the results of CRN models have large variances
on the Tumor growth dataset. Here we take a closer look into the results of
each repeated run and find that a part of the runs suffer from training
instability with the selected seeds. In Table 14, 15, 16, we use the 3 runs
with the lowest average estimation errors as the updated performance of the
CRN baseline (denoted as CRN∗), whose variances now are close to those of
other baselines. Still, we observe that a stronger CRN baseline does not
change our conclusions.

Table 14: Results of zero-shot transfer setup for multi-step outcome
estimation. We report the mean $\pm$ standard deviation of Rooted Mean Squared
Errors (RMSEs $\downarrow$) over 5 runs. Bold: the best results. Underline:
the 2nd best results.

Dataset | Method | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---
Tumor growth | MSM | 1.0515$\pm$0.0674 | 0.5048$\pm$0.0591 | 0.7583$\pm$0.0831 | 0.9685$\pm$0.1066 | 1.1561$\pm$0.1243 | 1.3372$\pm$0.1356 | 0.9627$\pm$0.0923 | 6.2%
RMSN | 1.2406$\pm$0.1301 | 1.0914$\pm$0.0346 | 1.1315$\pm$0.0634 | 1.1583$\pm$0.0810 | 1.1674$\pm$0.0913 | 1.1531$\pm$0.0919 | 1.1571$\pm$0.0660 | 22.0%
CRN(ERM) | 1.2924$\pm$0.0772 | 1.1769$\pm$0.1058 | 1.1728$\pm$0.1136 | 1.1906$\pm$0.1106 | 1.1997$\pm$0.1061 | 1.1883$\pm$0.0985 | 1.2035$\pm$0.0901 | 25.0%
CRN∗ | 1.6047$\pm$0.0487 | 2.0846$\pm$0.1665 | 2.0963$\pm$0.1274 | 2.1574$\pm$0.1255 | 2.2609$\pm$0.0569 | 2.3704$\pm$0.1368 | 2.0957$\pm$0.0514 | 68.9%
CT(ERM) | 0.9729$\pm$0.0718 | 1.0217$\pm$0.0292 | 1.1173$\pm$0.0457 | 1.1904$\pm$0.0395 | 1.2359$\pm$0.0618 | 1.2913$\pm$0.0939 | 1.1383$\pm$0.0251 | 20.7%
CT | 1.0272$\pm$0.1077 | 1.1428$\pm$0.2182 | 1.2708$\pm$0.2471 | 1.3608$\pm$0.2681 | 1.4166$\pm$0.2935 | 1.4322$\pm$0.3138 | 1.2751$\pm$0.2326 | 29.2%
G-Net | 1.0492$\pm$0.0529 | 1.0125$\pm$0.0767 | 1.1271$\pm$0.0876 | 1.2153$\pm$0.0777 | 1.2549$\pm$0.0727 | 1.2543$\pm$0.0678 | 1.1522$\pm$0.0537 | 21.7%
COSTAR | 0.8767$\pm$0.0492 | 0.7995$\pm$0.0853 | 0.8282$\pm$0.0801 | 0.9021$\pm$0.1062 | 0.9888$\pm$0.1280 | 1.0210$\pm$0.1168 | 0.9027$\pm$0.0814 | (-)

Table 15: Results of the data-efficient transfer learning setup for multi-step
outcome estimation. We report the mean $\pm$ standard deviation of Rooted Mean
Squared Errors (RMSEs $\downarrow$) over 5 runs. Bold: the best results.
Underline: the 2nd best results.

Dataset | Method | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---
Tumor growth | MSM | 1.0436$\pm$0.0671 | 0.5023$\pm$0.0588 | 0.7475$\pm$0.0829 | 0.9537$\pm$0.1060 | 1.1376$\pm$0.1233 | 1.3146$\pm$0.1338 | 0.9499$\pm$0.0915 | 7.8%
RMSN | 1.1839$\pm$0.0842 | 1.0912$\pm$0.0405 | 1.1215$\pm$0.0593 | 1.1538$\pm$0.0688 | 1.1728$\pm$0.0773 | 1.1740$\pm$0.0830 | 1.1495$\pm$0.0529 | 23.8%
CRN(ERM) | 1.2648$\pm$0.0689 | 1.1740$\pm$0.1015 | 1.1507$\pm$0.1016 | 1.1474$\pm$0.1070 | 1.1414$\pm$0.1041 | 1.1203$\pm$0.0906 | 1.1664$\pm$0.0786 | 24.9%
CRN∗ | 1.5019$\pm$0.0587 | 1.5362$\pm$0.0248 | 1.7824$\pm$0.1060 | 1.9842$\pm$0.1707 | 2.1317$\pm$0.2431 | 2.2546$\pm$0.3414 | 1.8651$\pm$0.1318 | 53.0%
CT(ERM) | 0.8947$\pm$0.0668 | 0.8700$\pm$0.0857 | 0.9507$\pm$0.1309 | 1.0031$\pm$0.1502 | 1.0371$\pm$0.1545 | 1.0668$\pm$0.1565 | 0.9704$\pm$0.1098 | 9.7%
CT | 0.9545$\pm$0.0782 | 0.9494$\pm$0.1597 | 1.0225$\pm$0.1562 | 1.1062$\pm$0.1377 | 1.1455$\pm$0.1192 | 1.1562$\pm$0.0953 | 1.0557$\pm$0.1136 | 17.0%
G-Net | 1.0335$\pm$0.0622 | 1.0154$\pm$0.1100 | 1.1105$\pm$0.1476 | 1.1859$\pm$0.1620 | 1.2257$\pm$0.1693 | 1.2198$\pm$0.1508 | 1.1318$\pm$0.1118 | 22.6%
COSTAR | 0.8654$\pm$0.0328 | 0.7945$\pm$0.0532 | 0.8248$\pm$0.0751 | 0.8754$\pm$0.0987 | 0.9378$\pm$0.1176 | 0.9594$\pm$0.1062 | 0.8762$\pm$0.0720 | (-)

Table 16: Results in standard supervised learning setting, with source and
target datasets coming from the same distribution for multi-step outcome
estimation. We report the mean +- standard deviation of Rooted Mean Squared
Errors (RMSEs) over 5 runs. Bold: the best results. Underline: the 2nd best
results.

Dataset | Method | $\tau=1$ | $\tau=2$ | $\tau=3$ | $\tau=4$ | $\tau=5$ | $\tau=6$ | Avg | Gain(%)
---|---|---|---|---|---|---|---|---|---
Tumor growth | MSM | 5.8368$\pm$0.6157 | 2.0400$\pm$0.6719 | 3.0385$\pm$0.9990 | 3.8701$\pm$1.2736 | 4.6173$\pm$1.5246 | 5.3823$\pm$1.7839 | 4.1308$\pm$1.1211 | 12.3%
RMSN | 4.8388$\pm$0.7770 | 5.4447$\pm$1.9202 | 5.9261$\pm$2.1096 | 5.9817$\pm$2.1270 | 5.8705$\pm$2.0544 | 5.5461$\pm$1.8865 | 5.6013$\pm$1.7727 | 35.4%
CRN(ERM) | 5.1601$\pm$0.5222 | 6.0784$\pm$2.3196 | 6.4721$\pm$2.4221 | 6.6142$\pm$2.4206 | 6.5648$\pm$2.3455 | 6.2939$\pm$2.1955 | 6.1972$\pm$2.0226 | 41.6%
CRN∗ | 4.8130$\pm$0.2296 | 6.3126$\pm$2.9523 | 6.6993$\pm$3.8805 | 6.7520$\pm$3.8551 | 6.8386$\pm$3.5630 | 6.8852$\pm$3.1150 | 6.3834$\pm$2.8863 | 43.3%
CT(ERM) | 5.1286$\pm$1.3377 | 5.7262$\pm$2.7601 | 6.5085$\pm$2.9886 | 6.9248$\pm$3.0009 | 7.1971$\pm$2.9346 | 7.2369$\pm$2.7570 | 6.4537$\pm$2.5904 | 43.9%
CT | 6.5485$\pm$1.5221 | 7.5382$\pm$2.8528 | 7.9030$\pm$2.9569 | 7.9828$\pm$2.9332 | 7.8244$\pm$2.8075 | 7.4418$\pm$2.6103 | 7.5398$\pm$2.5976 | 52.0%
G-Net | 3.9371$\pm$0.4023 | 3.7697$\pm$1.1861 | 4.6054$\pm$1.4181 | 4.9730$\pm$1.4773 | 5.0491$\pm$1.4410 | 4.8745$\pm$1.3153 | 4.5348$\pm$1.1778 | 20.1%
COSTAR | 3.7403$\pm$0.3695 | 3.0067$\pm$0.9065 | 3.4619$\pm$1.1557 | 3.8501$\pm$1.3127 | 3.9160$\pm$1.3142 | 3.7525$\pm$1.1493 | 3.6212$\pm$1.0040 | (-)