label shift, exist across domains. Overall, our experimental results highlight
the importance of the label shift correction (IW) step in the presence of
substantial label shift.

| $\rho({\text{Smiling}},{\text{Color\\_Balance}})$ | $P(\text{Smiling}=1)$
---|---|---
Domain | CelebA iv@ | CelebA v@ | CelebA vi@ | CelebA iv@ | CelebA v@ | CelebA vi@
Source Domain 1 | -0.8 | -0.8 | 0.2 | 0.5 | 0.5 | 0.5
Source Domain 2 | -0.6 | -0.6 | 0.4 | 0.5 | 0.5 | 0.5
Source Domain 3 | 0.6 | 0.6 | 0.6 | 0.5 | 0.5 | 0.5
Target Domain | 0.8 | 0.8 | 0.8 | 0.75 | 0.25 | 0.75

Table 5: Settings of CelebA iv@, v@, vi@ under color shift and label shift.
Smiling is the label attribute for prediction. We vary the Color_Balance
feature and the Smiling label distribution across domains. The correlation
between Smiling and the Color_Balance feature remains consistent with Table 3.
The distribution of Smiling is balanced (50% each) in all three source
domains, but becomes unbalanced (either 75% and 25%, or 25% and 75%) in the
target domain.

| CelebA iv@ | CelebA v@ | CelebA vi@
---|---|---|---
DA Algorithm | src_acc | tar_acc | src_acc | tar_acc | src_acc | tar_acc
ERM | 94.4$\pm$0.6 | 94.1$\pm$0.5 | 94.4$\pm$0.6 | 93.9$\pm$0.9 | 94.5$\pm$0.6 | 94.2$\pm$0.5
ERM-Pool | 93.6$\pm$0.2 | 88.7$\pm$0.9 | 93.6$\pm$0.2 | 89.6$\pm$0.9 | 93.8$\pm$0.3 | 93.1$\pm$0.5
DIP | 94.4$\pm$0.4 | 93.5$\pm$1.2 | 94.4$\pm$0.5 | 94.4$\pm$0.9 | 94.4$\pm$0.4 | 93.4$\pm$1.2
DIP-Pool | 93.7$\pm$0.6 | 88.0$\pm$0.6 | 93.6$\pm$0.5 | 89.5$\pm$1.2 | 93.7$\pm$0.5 | 92.9$\pm$0.9
CIP | 93.4$\pm$0.4 | 88.3$\pm$1.0 | 93.4$\pm$0.4 | 89.7$\pm$1.1 | 93.6$\pm$0.2 | 93.3$\pm$0.5
IW-ERM | 92.2$\pm$0.9 | 90.1$\pm$1.0 | 93.2$\pm$0.4 | 90.9$\pm$0.8 | 92.1$\pm$0.7 | 94.4$\pm$0.2
IW-CIP | 91.8$\pm$0.9 | 90.3$\pm$0.5 | 92.7$\pm$1.1 | 90.6$\pm$1.1 | 91.7$\pm$1.3 | 94.3$\pm$0.5
IW-DIP | 93.0$\pm$0.7 | 94.8$\pm$0.3 | 93.5$\pm$0.8 | 95.1$\pm$0.5 | 92.4$\pm$1.4 | 94.7$\pm$0.4
IRM | 93.8$\pm$0.3 | 88.7$\pm$0.9 | 93.5$\pm$0.4 | 89.6$\pm$0.9 | 93.7$\pm$0.2 | 93.1$\pm$0.6
V-REx | 92.5$\pm$0.2 | 89.5$\pm$0.9 | 92.6$\pm$0.3 | 90.2$\pm$1.2 | 93.8$\pm$0.2 | 93.4$\pm$0.6
groupDRO | 93.3$\pm$0.3 | 90.0$\pm$0.8 | 93.3$\pm$0.3 | 91.1$\pm$1.0 | 93.7$\pm$0.3 | 93.0$\pm$0.6

Table 6: Source and target accuracy in CelebA iv@, v@, vi@ under color shift
and label shift. For each setting of the CelebA, DA algorithms are run on 10
different datasets, each generated using 10 different random seeds. The top
four methods are highlighted in bold. The best method is colored in red.

###  5.4 Camelyon17 from WILDS

In this subsection, we conduct numerical experiments on Camelyon17 to
demonstrate the effectiveness of JointDIP on the real data. The Camelyon17
image dataset is one of the WILDS benchmarks (Koh et al., 2021; Sagawa et al.,
2021) for domain adaptation, which comprises whole-slide images (WSIs) of
breast cancer metastases in lymph node sections from five hospitals. The task
involves predicting whether a $96\times 96$ image patch extracted from a WSI
contains any tumor pixels within its central $32\times 32$ area. The labels
are obtained through manual pixel-level annotations of each WSI. Among the
five different hospitals (domains) from which the WSIs are collected, the
first three hospitals are treated as source domains (with 302.4K patches) and
the remaining two serve as validation domain (with 34.9K patches) and target
domain (with 85K patches), respectively. In addition to these labeled
examples, WILDS also provides additional unlabeled WSIs from the same five
hospitals and patches extracted from them, including 600K unlabeled patches
from the target domain. See (Koh et al., 2021; Sagawa et al., 2021) for a
detailed description of the Camelyon17 dataset.

Although the 600K unlabeled target patches are collected from the same
hospital as the 85K labeled target patches, there is no overlap in the WSIs
where these patches are extracted from (Sagawa et al., 2021). Furthermore, the
label distribution of unlabeled target patches is heavily skewed towards
negative, while labeled target patches are sampled in a class-balanced manner
(Sagawa et al., 2021). Consequently, this could create distribution shifts
between the “unlabeled target data” and the “labeled target data” in the
Camelyon17 dataset. In particular, this can lead DA algorithms, which utilize
target covariate information (e.g. DIP or its variants), to identify different
latent features depending on the source of target image patches used by these
algorithms. In the WILDS leaderboard, methods that leverage data from the
target domain utilize these 600K unlabeled patches, rather than treating the
85K labeled target patches as if they were unlabeled. However, if a
distribution shift exists between the “unlabeled target data” and the “labeled
target data”, there is no guarantee that our propose DA algorithms will
improve target performance. Therefore, our experiments examine two approaches
to utilize target data for DA methods: (1) using the 600K unlabeled target
patches, and (2) leveraging the 85K labeled target patches without label
access.

We use the standard DenseNet-121 model architecture (Huang et al., 2017) and
follow the training protocol described in Koh et al. (2021); Sagawa et al.
(2021). We implement CIP, DIP-Pool, and JointDIP-Pool, where CIP only uses
labeled source data, and DIP-Pool and JointDIP-Pool leverage target covariates
(image patches) from the target domain. We find that the pooled versions of
DIP and JointDIP increase reliability of the models by utilizing more samples
from multiple source domains, compared to their single source versions.
Therefore, we choose to implement the pooled version of DIP and JointDIP in
our experiments. Similarly, DANN (Ganin et al., 2016) and CORAL (Sun and
Saenko, 2016), which were initially developed for the single source version,
have been extended in (Sagawa et al., 2021) to take advantage of multiple
source domains. Following the standard submission guidelines in WILDS (Koh et
al., 2021; Sagawa et al., 2021), we do not use any data augmentation in CIP,
but allow DIP-Pool and JointDIP-Pool to utilize the same color augmentation as
described in Koh et al. (2021); Sagawa et al. (2021). Hyperparameters are
chosen via a grid search based on the validation accuracy, where we allow a
larger range of grids compared to previous experiments (see Appendix D.3 for
further details).

DA Algorithm | Target covariates | val accuracy | test accuracy111111The test accuracy is evaluated on 85K labeled target image patches.
---|---|---|---
IRM∗ | - | 86.2$\pm$1.4 | 64.2$\pm$8.1
groupDRO∗ | 85.5$\pm$2.2 | 68.4$\pm$7.3
ERM-Pool∗ | 84.9$\pm$3.1 | 70.3$\pm$6.4
CIP | 87.0$\pm$1.2 | 71.1$\pm$12.5
Fish∗ | 83.9$\pm$1.2 | 74.7$\pm$7.1
DANN∗ | 600K unlabeled target patches | 86.9$\pm$2.2 | 68.4$\pm$9.2
CORAL∗ | 90.4$\pm$0.9 | 77.9$\pm$6.6
DIP-Pool | 91.6$\pm$0.8 | 81.2$\pm$8.2
JointDIP-Pool | 91.5$\pm$0.6 | 82.7$\pm$5.2
DIP-Pool |  85K labeled target patches (target labels are not used) | 91.0$\pm$0.5 | 88.7$\pm$5.7
JointDIP-Pool | 91.7$\pm$0.6 | 91.9$\pm$3.1

Table 7: Accuracy of DA algorithms on Camelyon17 dataset. The first five DA
methods—IRM, GroupDRO, ERM-Pool, CIP, and Fish—only use labeled data from
multiple source domains, whereas DANN, CORAL, DIP-Pool, and JointDIP-Pool
additionally use target image patches from the target domain, but they do not
have access to the target labels. The algorithms are run over 10 trials.
Results for algorithms marked with ∗ are obtained from the WILDS benchmark.

The results are presented in Table 7, where we compare CIP with other DA
algorithms that enforce different invariances across all source domains, such
as IRM (Arjovsky et al., 2019) and Fish (Shi et al., 2021); and we also
compare DIP-Pool and JointDIP-Pool with other domain-invariant algorithms
using different distributional distances, like DANN (Ganin et al., 2016) and
CORAL (Sun and Saenko, 2016). We observe that CIP has a similar performance
compared to ERM-Pool, albeit with a slightly larger variance in test accuracy.
When covariates of 600K unlabeled target patches are utilized in the DIP
matching penalty, DIP-Pool results in an accuracy exceeding $80\%$.
Remarkably, JointDIP-Pool further improves this accuracy by $1.5\%$ while
substantially reducing the variance at the same time. In scenarios where
covariates of 85K labeled target patches are used in the DIP matching penalty,
the accuracy of both methods increases significantly, with JointDIP-Pool
improving DIP-Pool by over $3\%$. This result indicates that jointly matching
the covariates with CICs between source and target can help DIP-based
algorithms to identify the appropriate features to match across domains. In
addition, when comparing the accuracy of two different approaches of using
target covariates in the DIP matching penalty, our results demonstrate that
distribution shifts exist between the unlabeled and labeled target data in
Camelyon17, and leveraging the true target covariates from the 85K labeled
target patches proves more effective than using the target covariates from the
600K unlabeled target patches.

##  6 Related work

As a subfield of transfer learning, domain adaptation (DA), also known as
transductive transfer learning (Redko et al., 2020), aims to develop effective
algorithms when distribution shifts exist across training data (source
domains) and test data (target domain). In DA, it is usually assumed that we
have access to some unlabeled target data, and we aim to build a model
specifically for the target domain. Another related terminology domain
generalization (DG), or out-of-distribution generalization (OOD) (Ben-Tal et
al., 2009), instead assumes that unlabeled target data is unobtainable, and
seeks a model to generalize appropriately for all possible test domains. In
this paper, we do not distinguish these two terms particularly, and simply
consider the DA setting, i.e., we have access to an unlabeled target dataset,
and our goal is to find a classifier that performs well on the target domain.

One line of work in DA focuses on relating source and target domains by
learning a common feature representation across the domains. Pan et al. (2008,
2010) first proposed to find transferable components across domains in a
reproducing kernel Hilbert space, such that data distributions in different
domains are close. Similar designs of creating intermediate representations
across domains were investigated in (Gopalan et al., 2011; Gong et al., 2012),
and later the formal idea of matching probability distributions and extracting
invariant information was introduced in (Baktashmotlagh et al., 2013). They
named their approach domain invariant projection (DIP), which projects data to
a low-dimensional latent space where source and target covariates
distributions are matched. This DIP type of approach was then widely employed
in research on neural networks for DA (Sun and Saenko, 2016; Ganin et al.,
2016; Hoffman et al., 2018b). Specifically, Ganin et al. (2016) first proposed
to use Generative Adversarial Network (GAN) based distributional distances to
learn feature representations. Since then, subsequent research has emerged
that applies different neural-network-based distances (Long et al., 2017;
Courty et al., 2017; Long et al., 2018; Hoffman et al., 2018b; Peng et al.,
2019). While these studies have effectively demonstrated the empirical success
of their methodologies across various image and text datasets, there is still
a lack of understanding about the specific conditions under which these
methods can achieve good target performance.

From a theoretical point of view, the rationale for DIP has been rigorously
established by Ben-David et al. (2006, 2010), where they prove a target risk
bound via Vapnik-Chervonenkis (VC) theory. Further studies have analyzed the
source and target risk difference using different divergence measures (Mansour
et al., 2009; Cortes et al., 2010; Cortes and Mohri, 2011, 2014; Hoffman et
al., 2018a); see the survey by Redko et al. (2020) for a complete review. In
all these works, the target risk bound typically includes three components:
the source risk of the classifier, a divergence term measuring invariance of
the representation, and an optimal joint error term. Then DIP objective can be
viewed as minimizing the sum of the first two terms in order to achieve a low
target risk. However, Johansson et al. (2019); Zhao et al. (2019) argued that
DIP can completely fail in certain scenarios because the joint error term is
not observable and cannot be controlled by DIP. For instance, there can be
label-flipping features that achieve perfect source accuracy and invariance of
representation, but result in poor performance on the target domain.
Addressing how to avoid such cases in DIP has not been adequately explored.

The second class of DA methods arises from exploring invariance solely from
multiple labeled source domains, without using target covariates. One such
category of invariance is known as conditional invariance, which aims to
discover feature representations $\phi(X)$ that are invariant conditional on
the label $Y$ across source distributions. It was first introduced in Gong et
al. (2016), where the authors proposed to find conditional transferable
components after proper location-scale transformations. Later Heinze-Deml and
Meinshausen (2017) applied a related approach by classifying features into
”core” ones, which are conditionally invariant across domains, and ”style”
ones, whose distribution may shift significantly. They sought to construct a
classifier built upon only the core features by imposing the conditional
invariant penalty (CIP). Chen and Bühlmann (2020) further developed a
theoretical framework under structural equation models to analyze the effect
of CIP when the data generation process is anticausal. The emergence of
conditionally invariant features, or CICs, is a natural consequence under
anticausal data generation where the covariates, which are descendants of the
labels, remain unperturbed. Consequently, most tasks tackled by DA methods
utilizing conditional invariance typically involve anticausal problems.

It is also worth noting that under conditional invariance, it is convenient to
study the label shift problem, where the marginal distributions of the label
$Y$ are shifted across domains. This shift in label distribution is common in
many scenarios, for example, the distribution of ArXiv paper categories can be
influenced by changes in research topic trends over time (Wu et al., 2021).
Lipton et al. (2018) first proposed the label shift correction algorithm,
which estimates the amount of label shift by using the conditional invariance
of the covariates and the moment matching equation. To ensure numerically
stability, several variants of the algorithms have been further introduced.
For instance, Azizzadenesheli et al. (2019) formulated an $\ell_{2}$-norm
regularized least squares optimization problem, while Tachet des Combes et al.
(2020) considered a constrained least squares problem and the generalized
label shift assumption. Garg et al. (2020) introduced a maximum likelihood
estimation approach and provided a unified view of previous label shift
correction algorithms. More recently, Chen et al. (2022) has considered a
distributional shift named Sparse Joint Shift (SJS), which allows for both
labels and a few covariates to shift, however, the generalized label shift
considered in Tachet des Combes et al. (2020) and this paper permits
distribution shifts of more general latent feature representations.

Besides conditional invariance, other types of invariance have also attracted
increasing attention. For example, Arjovsky et al. (2019) came up with the
Invariant Risk Minimization (IRM) method to identify causal features via
enforcing the invariance of label $Y$ given the features $\phi(X)$. Following
IRM, analogous invariance has been explored, such as the risk invariance
enforced by variance across domains (Krueger et al., 2021; Xie et al., 2020)
and gradient invariance (Koyama and Yamaguchi, 2020; Shi et al., 2021). While
these approaches were initially motivated to identify invariant features in
causal data models, they have also been applied to anticausal problems.
However, as pointed out by Rosenfeld et al. (2020); Kamath et al. (2021), IRM
may fail to capture the correct invariance, especially when the model is non-
linear. Counterfactual invariance is another type of invariance introduced in
(Veitch et al., 2021; Jiang and Veitch, 2022), where the aim is to seek
representations that are counterfactually invariant to a spurious factor of
variation. (Wang and Veitch, 2022) show that the formulation of learning
counterfactual invariant representations is closely related to different
formulations of invariant representation learning algorithms, depending on the
underlying causal structure of the data. Aside from the robustness achieved by
enforcing specific invariance, a more general framework is distributionally
robust optimization (DRO) (Ben-Tal et al., 2013; Duchi et al., 2016). For
instance, groupDRO (Sagawa et al., 2019) intends to improve worst-group
performance using an online algorithm to update group weights. It should be
noted that related study of invariance and robustness has also appeared in
works from a causal inference point of view (Peters et al., 2016; Meinshausen,
2018; Rothenhäusler et al., 2018; Magliacane et al., 2018).

Other DA algorithms have also been introduced based on different assumptions
relating source and target data, which do not fall into the above classes of
methods relying on invariance. For instance, self training, which originated
from the semi-supervised learning literature (Chapelle et al., 2009),
recursively adapts a classifier to fit pseudolabels predicted by a previous
model using unlabeled data in a new domain (Amini and Gallinari, 2003).
Theoretical properties of self training has been studied recently in (Kumar et
al., 2020; Wei et al., 2020). Data augmentation serves as another beneficial
way of obtaining better performance (Simard et al., 1998; Zhang et al., 2017;
Yao et al., 2022), especially when possible perturbations across domains are
well understood. In addition, a different perspective of DA stems from meta-
learning (Thrun and Pratt, 2012; Finn et al., 2017). This “learning-to-learn”
approach aims to distill knowledge of multiple learning episodes to improve
future learning performance.

##  7 Discussion

In this paper, we highlight three prominent roles of CICs in DA. First, we
propose the IW-CIP algorithm based on CICs, which can solve DA problems beyond
simple covariate shift and label shift with theoretical guarantees.
Specifically, we establish an explicit target risk bound for IW-CIP and
further provide fine-grained analysis for each term in the bound. Second, we
demonstrate the advantage of using CICs to measure the performance of other DA
algorithms. For example, a conditionally invariant classifier built upon CICs
can be applied as an approximate proxy for target labels to refute DA
algorithms with poor target performance. Finally, to solve the label-flipping
issue of DIP, we introduce the JointDIP algorithm which jointly matches the
distribution of DIP features and CICs. Both our theoretical results and
numerical experiments support the benefits of using CICs in DA.

While CICs across multiple source domains can be beneficial in DA, in practice
it may not always be advantageous to use all available source domains for
identifying CICs. Particularly, if not all source domains are related to the
target domain, or if there aren’t any features invariant across all source
domains, it would be more efficient to select a subset of source domains for
learning CICs. However, choosing the source domains that are related to the
target domain or share common features is difficult without access to target
labels. In such settings, one feasible approach is to learn invariant features
across possible combinations of source domains and refer to domain experts to
investigate those features, which would help determine the selection of source
domains. We leave the development of a framework for best source domain
selection for future research.

There are other promising directions to pursue following our work. While we
provide bounds in Section 4.1 to detect failure of DA algorithms, these bounds
may not be tight and potentially can be improved based on other structural
assumptions on the data generation model. Moreover, although in this paper we
mainly focus on combing CICs with other DA algorithms, it is possible that
other forms of invariance can also be applied to improve DA algorithms. For
instance, can we use features obtained by IRM (Arjovsky et al., 2019), V-REx
(Krueger et al., 2021), or other domain generalization algorithms to aid and
enhance existing DA methods such as DIP with target risk guarantees? Finally,
this paper mainly considers the unsupervised domain adaptation scenario where
target labels are unavailable, however, datasets such as Camelyon17 from WILDS
(Koh et al., 2021) may have a small portion of target labels available.
Developing DA algorithms with target risk guarantees in such settings is
therefore of practical interest.

##  Appendix A Illustrative examples of general anticausal model

In this section, we provide several examples for the general anticausal model
and demonstrate how the conditions in Theorem 3(b) and Theorem 3(c) are
satisfied under these models.

###  A.1 Example with CICs and label-flipping features

We present an example of the general anticausal model, as defined in
Definition 4, which incorporates both CICs and label-flipping features.

##### Example 1

Let $Y\in\mathcal{Y}=\\{1,2\\}$ be a binary random variable and balanced in
both domains, i.e.,
${p^{{(m)}}_{1}}={p^{{(m)}}_{2}}={p^{{({\mathfrak{T}})}}_{1}}={p^{{({\mathfrak{T}})}}_{2}}=1/2$.
Let $X\in{}^{3}$ be a three-dimensional covariate vector with mechanism
function $f^{(1)}$ and $f^{{({\mathfrak{T}})}}$ as given by

$\displaystyle f^{(1)}(1)$ $\displaystyle=\begin{pmatrix}-1\\\ -1\\\
-1\end{pmatrix},\;f^{(1)}(2)=\begin{pmatrix}1\\\ 1\\\ 1\end{pmatrix};$
$\displaystyle f^{{({\mathfrak{T}})}}(1)$ $\displaystyle=\begin{pmatrix}-1\\\
0\\\ 1\end{pmatrix},\;f^{{({\mathfrak{T}})}}(2)=\begin{pmatrix}1\\\ 2\\\
-1\end{pmatrix},$

and where the noise terms $\epsilon^{(1)}$ and $\epsilon^{{({\mathfrak{T}})}}$
are drawn from a Gaussian distribution
$\mathcal{P}_{\epsilon}=\mathcal{N}(0,\mathbb{I}_{3})$. Examining each
coordinate, we can see that the first coordinate ${X}_{[1]}$ maintains its
conditional invariance across source and target, and thus it is a CIC. The
second coordinate ${X}_{[2]}$ is perturbed, but the correlations
$\rho(X^{(1)}_{[2]},Y^{(1)})=\rho(X^{{({\mathfrak{T}})}}_{[2]},Y^{{({\mathfrak{T}})}})=\sqrt{2}/2$
remian unchanged. The third coordinate ${X}_{[3]}$ is also perturbed, and by
Definition 3, it is a label-flipping feature, as
$\rho(X^{(1)}_{[3]},Y^{(1)})=-\rho(X^{{({\mathfrak{T}})}}_{[3]},Y^{{({\mathfrak{T}})}})=\sqrt{2}/2$.

###  A.2 Demonstration of condition in Theorem 3(b)

The condition (35) in Theorem 3(b) requires that for the linear conditionally
invariant feature mapping $\phi_{\textnormal{inv}}$, CICs must have different
conditional means for different labels. In other words, the first-moment of
the CICs should suffice to distinguish between labels. We show that under the
anticausal model given in the example above, this condition is indeed
satisfied.

##### Example 1 Cont.

In Example 1, one possible linear conditionally invariant feature mapping is
$\phi_{\textnormal{inv}}(x)={x}_{[1]}$, i.e., it projects to the first
coordinate of the covariates. Simple algebra shows that

$\displaystyle\mathbb{E}\left[{\phi_{\textnormal{inv}}(X^{(1)})}\ \middle|\
{Y^{(1)}=1}\right]=-1,\text{ and
}\mathbb{E}\left[{\phi_{\textnormal{inv}}(X^{(1)})}\ \middle|\
{Y^{(1)}=2=1}\right],$

which confirms the condition (35).

###  A.3 Demonstration of condition in Theorem 3(c)

While the condition (35) for linear $\phi_{\textnormal{inv}}$ is relatively
straightforward, verifying the full-rank condition of matrix
$C_{\phi_{\textnormal{inv}}}(a)$ for general $\phi_{\textnormal{inv}}$ in
Theorem 3(c) can be more challenging due to the computation of higher-moments.
Here we provide several examples where a nonlinear conditionally invariant
feature mapping $\phi_{\textnormal{inv}}$ meets the full rank condition.

##### Example 1 Cont.

Consider the general anticausal model introduced in Example 1 ($L=2$) and let

$\displaystyle\begin{split}\phi_{\textnormal{inv}}^{k}(x)&={x}_{[1]}^{k},\ \
k=0,1,2,\ldots,\text{ and }\\\
\phi_{\textnormal{inv}}^{e}(x)&=\exp({x}_{[1]}).\end{split}$ (37)

Since ${X}_{[1]}$ is a CIC, we know that
$\phi_{\textnormal{inv}}^{k}({X}_{[1]})$ and
$\phi_{\textnormal{inv}}^{e}({X}_{[1]})$ are also CICs. Take $a=1$ in the
definition of the matrix $C_{\phi_{\textnormal{inv}}}(a)$. From the standard
calculation of moments for Gaussian random variables, we have

$\displaystyle C_{\phi_{\textnormal{inv}}^{k}}(a)$
$\displaystyle=\begin{pmatrix}1&1\\\
\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u-1)^{k}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+1)^{k}}\right]\end{pmatrix},\;\;C_{\phi_{\textnormal{inv}}^{e}}(a)=\begin{pmatrix}1&1\\\
e^{-\frac{1}{2}}&e^{\frac{3}{2}}\end{pmatrix}.$

We conclude that $C_{\phi_{\textnormal{inv}}^{k}}(a)$ is full rank if and only
if $k$ is odd, and $C_{\phi_{\textnormal{inv}}^{e}}(a)$ is a full rank matrix.
Therefore,
$\phi_{\textnormal{inv}}^{2k+1}(x)={x}_{[1]}^{2k+1}(k=0,1,2,\ldots)$ and
$\phi_{\textnormal{inv}}^{e}(x)=\exp({x}_{[1]})$ satisfy the full rank
condition.

##### Example 2

Consider the general anticausal model from Assumption 3 where
$Y\in\mathcal{Y}=\\{1,2,3\\}$ with $L=3$. Let $X\in{}^{2}$ be a two-
dimensional covariate vector with mechanism functions $f^{(1)}$ and
$f^{{({\mathfrak{T}})}}$ defined as follows:

$\displaystyle f^{(1)}(1)$ $\displaystyle=\begin{pmatrix}1\\\
-1\end{pmatrix},\;f^{(1)}(2)=\begin{pmatrix}2\\\ 0\end{pmatrix},\
f^{(1)}(3)=\begin{pmatrix}3\\\ 1\end{pmatrix};$ $\displaystyle
f^{{({\mathfrak{T}})}}(1)$ $\displaystyle=\begin{pmatrix}1\\\
1\end{pmatrix},\;f^{{({\mathfrak{T}})}}(2)=\begin{pmatrix}2\\\
0\end{pmatrix},\ f^{{({\mathfrak{T}})}}(3)=\begin{pmatrix}3\\\
-1\end{pmatrix}.$

That is, the first coordinate ${X}_{[1]}$ is a CIC, while the second
coordinate ${X}_{[2]}$ is perturbed. Assume that the noise distribution
$\mathcal{P}_{\epsilon}=\mathcal{N}(0,\mathbb{I}_{2})$ follows a Gaussian
distribution. Consider the feature mappings $\phi_{\textnormal{inv}}^{2}$ and
$\phi_{\textnormal{inv}}^{e}$ from Eq. (37) as the conditionally invariant
feature mappings, and take $a=1$. Basic algebra shows that

$\displaystyle C_{\phi_{\textnormal{inv}}^{2}}(a)$
$\displaystyle=\begin{pmatrix}1&1&1\\\
\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+1)^{2}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+2)^{2}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+3)^{2}}\right]\\\
\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+1)^{4}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+2)^{4}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{(u+3)^{4}}\right]\end{pmatrix}=\begin{pmatrix}1&1&1\\\
2&5&10\\\ 10&43&138\end{pmatrix},$ $\displaystyle
C_{\phi_{\textnormal{inv}}^{e}}(a)$ $\displaystyle=\begin{pmatrix}1&1&1\\\
\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{e^{u+1}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{e^{u+2}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{e^{u+3}}\right]\\\
\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{e^{2(u+1)}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{e^{2(u+2)}}\right]&\mathbb{E}_{u\sim\mathcal{N}(0,1)}\left[{e^{2(u+3)}}\right]\end{pmatrix}=\begin{pmatrix}1&1&1\\\
e^{\frac{3}{2}}&e^{\frac{5}{2}}&e^{\frac{7}{2}}\\\
e^{4}&e^{6}&e^{8}\end{pmatrix},$

both of which are full-rank, and therefore satisfying the full rank condition
in Theorem 3(c).

##### Example 3

Consider the general anticausal model from Assumption 3 with arbitrary number
of classes $L\geq 2$. Let

$\displaystyle f^{(m)}(y)=\begin{pmatrix}e_{y}\\\
a_{y}^{(m)}\end{pmatrix},f^{{({\mathfrak{T}})}}(y)=\begin{pmatrix}e_{y}\\\
a_{y}^{{({\mathfrak{T}})}}\end{pmatrix},$

where $e_{i}\in{}^{L}$ denotes the zero vector with only the $i$-th element
being $1$, and $a_{y}^{(m)},a_{y}^{{({\mathfrak{T}})}}\in{}^{p-L}$ $(p>L)$ are
arbitrary vectors. That is, the first $L$ covariates are CICs while the last
$p-L$ coordinates can be arbitrarily perturbed. Let

$\displaystyle\phi_{\textnormal{inv}}^{k}(X)=\begin{pmatrix}{X}_{[1]}^{k},{X}_{[2]}^{k},\cdots,{X}_{[L]}^{k}\end{pmatrix}^{\top},\
\ k=1,2,\ldots.$

For a fixed vector $a\in{}^{L}$, verifying
$C_{\phi_{\textnormal{inv}}^{k}}(a)$ is full rank can be technical depending
on the distribution $\mathcal{P}_{\epsilon}$. Instead, we show an illustrative
explanation rather than a rigorous proof. First, let’s consider
$\mathcal{P}_{\epsilon}=\delta_{0}(\cdot)$, the Dirac distribution centered at
0. Then we can compute

$\displaystyle
C_{\phi_{\textnormal{inv}}^{k}}(a)=\begin{pmatrix}1&1&\cdots&1\\\
{a}_{[1]}&{a}_{[2]}&\cdots&{a}_{[L]}\\\ \vdots&\vdots&\ddots&\vdots\\\
{a}_{[1]}^{L-1}&{a}_{[2]}^{L-1}&\cdots&{a}_{[L]}^{L-1}\end{pmatrix},$

which is a Vandermonde matrix and is full rank as long as we take
${a}_{[i]}\neq{a}_{[j]},i\neq j$. Next we consider a case where
$\mathcal{P}_{\epsilon}$ is a distribution concentrated around 0 (e.g. a
Gaussian distribution with very small variance). In this case, the matrix will
be slightly perturbed but is still full rank by the continuity of high-order
moments and matrix determinants.

##  Appendix B Proof of Section 3 and additional results on IW-CIP

In this section we prove our main theorems and propositions given in Section
3, as well as additional results given in Proposition 5 and Lemma 1.
Proposition 5 provides the finite-sample bound for the IW-CIP penalty when it
is defined by the squared mean distance, and Lemma 1 proves the number of
required source domains to validate the conditions in Proposition 3 when the
perturbations are randomly generated.

###  B.1 Proof of Theorem 1.A

Let $w=(w^{(1)},\ldots,w^{(M)})$ where $w^{(m)}$ is the true importance weight
for the $m$-th source domain. For any classifier $h=g\circ\phi$, we have

$\displaystyle\mathcal{R}^{{({\mathfrak{T}})}}({h})-\overline{\mathcal{R}}({h};{w})$
$\displaystyle={\frac{1}{M}}\sum_{m=1}^{M}\left(\mathbb{E}\left[{\mathbf{1}_{{h(X^{{({\mathfrak{T}})}})\neq
Y^{{({\mathfrak{T}})}}}}}\right]-\mathbb{E}\left[{{w}_{[Y^{(m)}]}^{(m)}\cdot\mathbf{1}_{{h(X^{(m)})\neq
Y^{(m)}}}}\right]\right)$
$\displaystyle={\frac{1}{M}}\sum_{m=1}^{M}\sum_{y=1}^{L}\left(\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})\neq
y}\ \middle|\
{Y^{{({\mathfrak{T}})}}=y}\right\\}-\mathbb{P}\left\\{{h(X^{(m)})\neq y}\
\middle|\
{Y^{(m)}=y}\right\\}\right)\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=y}\right\\}$
$\displaystyle={\frac{1}{M}}\sum_{m=1}^{M}\sum_{y=1}^{L}\left(\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=y}\
\middle|\
{Y^{{({\mathfrak{T}})}}=y}\right\\}-\mathbb{P}\left\\{{h(X^{(m)})=y}\
\middle|\
{Y^{(m)}=y}\right\\}\right)\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=y}\right\\},$

where the second step applies
${w}_{[y]}^{(m)}=\frac{\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=y}\right\\}}{\mathbb{P}\left\\{{Y^{(m)}=y}\right\\}}$
and the law of total expectation. By definition of
$D_{\mathcal{G}}\left({\cdot},{\cdot}\right)$, it is then upper bounded by

$\displaystyle\left|\mathcal{R}^{{({\mathfrak{T}})}}({h})-\overline{\mathcal{R}}({h};{w})\right|$
$\displaystyle\leq{\frac{1}{M}}\sum_{m=1}^{M}\sum_{y=1}^{L}D_{\mathcal{G}}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\right)\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=y}\right\\}$
$\displaystyle\leq\max_{\begin{subarray}{c}m=1,\ldots,M,\\\
y=1,\ldots,L\end{subarray}}D_{\mathcal{G}}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\right)\sum_{m=1}^{M}{\frac{1}{M}}\sum_{y=1}^{L}\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=y}\right\\}$
$\displaystyle=\max_{\begin{subarray}{c}m=1,\ldots,M,\\\
y=1,\ldots,L\end{subarray}}D_{\mathcal{G}}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\right)$
$\displaystyle=\Psi_{\mathcal{G},\phi},$ (38)

where we use
$\sum_{y=1}^{L}\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=y}\right\\}=1$ and
$\sum_{m=1}^{M}1/M=1$ in the second step. Furthermore, we can bound

$\displaystyle\left|\mathcal{R}^{(m)}({h};{w^{(m)}})-\mathcal{R}^{(m)}({h};{\widehat{w}^{(m)}})\right|$
$\displaystyle=\left|\mathbb{E}\left[{\left({w}_{[Y^{(m)}]}^{(m)}-{\widehat{w}}_{[Y^{(m)}]}^{(m)}\right)\cdot\mathbf{1}_{{h(X^{(m)})\neq
Y^{(m)}}}}\right]\right|$
$\displaystyle\leq\lVert{w^{(m)}-\widehat{w}^{(m)}}\rVert_{\infty}.$

It follows that

$\displaystyle\left|\overline{\mathcal{R}}({h};{w})-\overline{\mathcal{R}}({h};{\widehat{w}})\right|$
$\displaystyle\leq{\frac{1}{M}}\sum_{m=1}^{M}\left|\mathcal{R}^{(m)}({h};{w^{(m)}})-\mathcal{R}^{(m)}({h};{\widehat{w}^{(m)}})\right|$
$\displaystyle\leq\max_{m=1,\ldots,M}\lVert{w^{(m)}-\widehat{w}^{(m)}}\rVert_{\infty}.$
(39)

The result now follows from combining Eq. (B.1) and Eq. (B.1).

###  B.2 Proof of Theorem 1.B

Denote the empirical average source risk across $M$ environments by

$\displaystyle\widehat{\overline{\mathcal{R}}}({h};{{\rm
w}^{(m)}})=\frac{1}{M}\sum_{m=1}^{M}\widehat{\mathcal{R}}^{(m)}({h};{{\rm
w}^{(m)}}).$

Then we can write

$\displaystyle\quad\mathcal{R}^{{({\mathfrak{T}})}}({\widehat{h}_{\textnormal{IW-
CIP}}})-\mathcal{R}^{{({\mathfrak{T}})}}({h^{\star}})$
$\displaystyle=\underbrace{\mathcal{R}^{{({\mathfrak{T}})}}({\widehat{h}_{\textnormal{IW-
CIP}}})-\overline{\mathcal{R}}({\widehat{h}_{\textnormal{IW-
CIP}}};{w})}_{(A_{1})}+\underbrace{\overline{\mathcal{R}}({\widehat{h}_{\textnormal{IW-
CIP}}};{w})-\widehat{\overline{\mathcal{R}}}({\widehat{h}_{\textnormal{IW-
CIP}}};{w})}_{(B_{1})}+$
$\displaystyle\quad+\underbrace{\widehat{\overline{\mathcal{R}}}({\widehat{h}_{\textnormal{IW-
CIP}}};{w})-\widehat{\overline{\mathcal{R}}}({\widehat{h}_{\textnormal{IW-
CIP}}};{\widehat{w}})}_{(C_{1})}+\underbrace{\widehat{\overline{\mathcal{R}}}({\widehat{h}_{\textnormal{IW-
CIP}}};{\widehat{w}})-\widehat{\overline{\mathcal{R}}}({h^{\star}};{\widehat{w}})}_{(D)}+\underbrace{\widehat{\overline{\mathcal{R}}}({h^{\star}};{\widehat{w}})-\widehat{\overline{\mathcal{R}}}({h^{\star}};{w})}_{(C_{2})}$
$\displaystyle\quad+\underbrace{\widehat{\overline{\mathcal{R}}}({h^{\star}};{w})-\overline{\mathcal{R}}({h^{\star}};{w})}_{(B_{2})}+\underbrace{\overline{\mathcal{R}}({h^{\star}};{w})-\mathcal{R}^{{({\mathfrak{T}})}}({h^{\star}})}_{(A_{2})}.$

##### Term $A_{1}$ and $A_{2}$

From Eq. (B.1) in the proof of Theorem 1.A, we have

$\displaystyle\begin{split}A_{1}&=\mathcal{R}^{{({\mathfrak{T}})}}({\widehat{h}_{\textnormal{IW-
CIP}}})-\overline{\mathcal{R}}({\widehat{h}_{\textnormal{IW-
CIP}}};{w})\leq\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{IW-CIP}}},\\\
A_{2}&=\overline{\mathcal{R}}({h^{\star}};{w})-\mathcal{R}^{{({\mathfrak{T}})}}({h^{\star}})=0.\end{split}$
(40)

The second equation follows from the definition of
$h^{\star}=g^{\star}\circ\phi^{\star}$, i.e., $\phi^{\star}$ is a
conditionally invariant feature mapping and satisfies
$\Psi_{\mathcal{G},\phi^{\star}}=0$.

##### Term $B_{1}$ and $B_{2}$

For any $h=g\circ\phi$, we have
$\left|{w}_{[Y^{(m)}]}^{(m)}\mathbf{1}_{{h(X^{(m)})\neq
Y^{(m)}}}\right|\leq\left\|w^{(m)}\right\|_{\infty}$. Then the generalization
bound based on Rademacher complexity (see e.g. (Wainwright, 2019, Theorem
4.10)) shows that, for the $m$-th source domain, the bound

$\displaystyle\sup_{\begin{subarray}{c}h=g\circ\phi\\\
g\in\mathcal{G},\phi\in\Phi\end{subarray}}\left|\mathcal{R}^{(m)}({h};{w^{(m)}})-\widehat{\mathcal{R}}^{(m)}({h};{w^{(m)}})\right|$
$\displaystyle\leq
2\mathfrak{R}_{n^{(m)},\mathcal{P}^{(m)}}\left({\mathcal{H}(w^{(m)},\mathcal{G},\Phi)}\right)$
$\displaystyle\hskip
50.00008pt+\left\|w^{(m)}\right\|_{\infty}\sqrt{\frac{2\log(M/\delta)}{n^{(m)}}},$

holds with probability at least $1-\delta/M$. Using the union bound over each
$m\in\\{1,\ldots,M\\}$, we obtain that with probability at least $1-\delta$,

$\displaystyle\left|B_{1}\right|+\left|B_{2}\right|$ $\displaystyle\leq
2\sup_{\begin{subarray}{c}h=g\circ\phi,\\\
g\in\mathcal{G},\phi\in\Phi\end{subarray}}{\frac{1}{M}}\sum_{m=1}^{M}\left|\mathcal{R}^{(m)}({h};{w^{(m)}})-\widehat{\mathcal{R}}^{(m)}({h};{w^{(m)}})\right|$
$\displaystyle\leq\max_{m=1,\ldots,M}\left[4\mathfrak{R}_{n^{(m)},\mathcal{P}^{(m)}}\left({\mathcal{H}(w^{(m)},\mathcal{G},\Phi)}\right)+2\left\|w^{(m)}\right\|_{\infty}\sqrt{\frac{2\log(M/\delta)}{n^{(m)}}}\right].$
(41)

##### Term $C_{1}$ and $C_{2}$

For the $m$-th source domain, let us define

$\displaystyle\widehat{\ell}^{(m)}_{j}=\sum_{k=1}^{n}\mathbf{1}_{{Y^{(m)}_{k}=j}}\mathbf{1}_{{h(X^{(m)}_{k})\neq
Y^{(m)}_{k}}}.$

Then $\lVert{\widehat{\ell}^{(m)}}\rVert_{1}\leq n$ and so for any
$h=g\circ\phi$,

$\displaystyle\left|\widehat{\mathcal{R}}^{(m)}({h};{w^{(m)}})-\widehat{\mathcal{R}}^{(m)}({h};{\widehat{w}^{(m)}})\right|$
$\displaystyle=\left|\frac{1}{n}\sum_{k=1}^{n}\left({w}_{[Y^{(m)}_{k}]}^{(m)}-{\widehat{w}}_{[Y^{(m)}_{k}]}^{(m)}\right)\mathbf{1}_{{h(X^{(m)}_{k})\neq
Y^{(m)}_{k}}}\right|$
$\displaystyle=\left|\frac{1}{n}\sum_{j=1}^{L}\left({w}_{[j]}^{(m)}-{\widehat{w}}_{[j]}^{(m)}\right)\widehat{\ell}^{(m)}_{j}\right|$
$\displaystyle\leq\frac{1}{n}\left\|w^{(m)}-\widehat{w}^{(m)}\right\|_{\infty}\left\|\widehat{\ell}^{(m)}\right\|_{1}\leq\left\|w^{(m)}-\widehat{w}^{(m)}\right\|_{\infty}.$

It follows that

$\displaystyle\left|C_{1}\right|+\left|C_{2}\right|$ $\displaystyle\leq
2\sup_{\begin{subarray}{c}h=g\circ\phi,\\\
g\in\mathcal{G},\phi\in\Phi\end{subarray}}{\frac{1}{M}}\sum_{m=1}^{M}\left|\widehat{\mathcal{R}}^{(m)}({h};{w^{(m)}})-\widehat{\mathcal{R}}^{(m)}({h};{\widehat{w}^{(m)}})\right|$
$\displaystyle\leq
2\max_{m=1,\ldots,M}\left\|w^{(m)}-\widehat{w}^{(m)}\right\|_{\infty}.$ (42)

##### Term $D$

By definition of $\widehat{h}_{\textnormal{IW-CIP}}$ from Eq. (3.2), we know
that

$\widehat{\overline{\mathcal{R}}}({\widehat{h}_{\textnormal{IW-
CIP}}};{\widehat{w}})+\widehat{\Lambda}_{\widehat{\phi}_{\textnormal{IW-
CIP}}}\leq\widehat{\overline{\mathcal{R}}}({h^{\star}};{\widehat{w}})+\widehat{\Lambda}_{\phi^{\star}},$

and so rearranging terms,

$\displaystyle
D=\widehat{\overline{\mathcal{R}}}({\widehat{h}_{\textnormal{IW-
CIP}}};{\widehat{w}})-\widehat{\overline{\mathcal{R}}}({h^{\star}};{\widehat{w}})\leq\widehat{\Lambda}_{\phi^{\star}}-\widehat{\Lambda}_{\widehat{\phi}_{\textnormal{IW-
CIP}}}\leq\widehat{\Lambda}_{\phi^{\star}}.$ (43)

The theorem then follows from combining Eq. (40), (B.2), (B.2), and (43).

###  B.3 Proof of Proposition 1

For a family of functions $\mathcal{G}$ mapping $\mathcal{Z}$ to $[0,1]$ and a
distribution $\mathcal{P}$ on $\mathcal{Z}$, let $\widehat{\mathcal{P}}_{n}$
denote the empirical distribution of $n$ $i.i.d.$ samples
$(Z_{1},Z_{2},\ldots,Z_{n})$ from $\mathcal{P}$. The standard Rademacher
complexity theory (see e.g. (Wainwright, 2019, Theorem 4.10)) states that for
any $\delta>0$, with probability at least $1-\delta$, we have

$\displaystyle\sup_{g\in\mathcal{G}}\left|\mathbb{E}_{Z\sim\mathcal{P}}\left[{g(Z)}\right]-\mathbb{E}_{Z\sim\widehat{\mathcal{P}}_{n}}\left[{g(Z)}\right]\right|\leq
2\mathfrak{R}_{n,\mathcal{P}}\left({\mathcal{G}}\right)+\sqrt{\frac{\log(2/\delta)}{2n}}.$

By definition of $\mathcal{G}$-divergence, with probability $1-\delta$ we have

$\displaystyle
D_{\mathcal{G}}\left({\mathcal{P}},{\widehat{\mathcal{P}}_{n}}\right)$
$\displaystyle=\sup_{g\in\mathcal{G}}\max_{y=1,2,\ldots,L}\left|\mathbb{E}_{Z\sim\mathcal{P}}\left[{\mathbf{1}_{{g(Z)=y}}}\right]-\mathbb{E}_{Z\sim\widehat{\mathcal{P}}_{n}}\left[{\mathbf{1}_{{g(Z)=y}}}\right]\right|$
$\displaystyle=\max_{y=1,2,\ldots,L}\sup_{g\in\mathcal{G}}\left|\mathbb{E}_{Z\sim\mathcal{P}}\left[{\mathbf{1}_{{g(Z)=y}}}\right]-\mathbb{E}_{Z\sim\widehat{\mathcal{P}}_{n}}\left[{\mathbf{1}_{{g(Z)=y}}}\right]\right|$
$\displaystyle\leq\max_{y=1,2,\ldots,L}\left(2\mathfrak{R}_{n,\mathcal{P}}\left({\mathcal{G}_{y}}\right)+\sqrt{\frac{\log(2L/\delta)}{2n}}\right),$
(44)

where
$\mathcal{G}_{y}\coloneqq\left\\{\mathbf{1}_{{g(x)=y}},g\in\mathcal{G}\right\\}$.
Then we can bound

$\displaystyle\quad
D_{\mathcal{G}}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}},{\widehat{\mathcal{P}}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}}\right)$
$\displaystyle\overset{(i)}{\leq}D_{\mathcal{G}}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}}\right)+D_{\mathcal{G}}\left({\mathcal{P}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}},{\mathcal{P}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}}\right)$
$\displaystyle\hskip
30.00005pt+D_{\mathcal{G}}\left({\mathcal{P}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}},{\widehat{\mathcal{P}}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}}\right)$
$\displaystyle\overset{(ii)}{=}D_{\mathcal{G}}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}}\right)+D_{\mathcal{G}}\left({\widehat{\mathcal{P}}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}},{\mathcal{P}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}}\right),$
(45)

where in step $(i)$ we apply triangular inequality, and in step $(ii)$ we use
the fact that $\phi^{\star}$ is a conditionally invariant feature mapping
across source domains to cancel the middle term. Applying union bound, we get

$\displaystyle\widehat{\Lambda}_{\phi^{\star}}$
$\displaystyle=\frac{\lambda_{\textnormal{IW-
CIP}}}{LM^{2}}\sum_{y=1}^{L}\sum_{m\neq
m^{\prime}}D_{\mathcal{G}}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}},{\widehat{\mathcal{P}}^{(m^{\prime})}_{{\phi^{\star}(X)}|{Y=y}}}\right)$
$\displaystyle\overset{(i)}{\leq}\frac{2\lambda_{\textnormal{IW-
CIP}}}{LM}\sum_{y=1}^{L}\sum_{m=1}^{M}D_{\mathcal{G}}\left({\mathcal{P}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}},{\widehat{\mathcal{P}}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}}\right)$
$\displaystyle\overset{(ii)}{\leq}\frac{4\lambda_{\textnormal{IW-
CIP}}}{LM}\cdot\sum_{y=1}^{L}\sum_{m=1}^{M}\max_{y^{\prime}=1,2,\ldots,L}\mathfrak{R}_{n^{(m)},\mathcal{P}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}}\left({\mathcal{G}_{y^{\prime}}}\right)+\frac{2\lambda_{\textnormal{IW-
CIP}}}{M}\sum_{m=1}^{M}\sqrt{\frac{\log(2LM/\delta)}{2n^{(m)}}}$
$\displaystyle\leq 2\lambda_{\textnormal{IW-
CIP}}\left(2\max_{\begin{subarray}{c}m=1,2,\ldots,M\\\ y=1,2,\ldots,L\\\
y^{\prime}=1,2,\ldots,L\end{subarray}}\mathfrak{R}_{n^{(m)},\mathcal{P}^{(m)}_{{\phi^{\star}(X)}|{Y=y}}}\left({\mathcal{G}_{y^{\prime}}}\right)+\max_{m=1,2,\ldots,M}\sqrt{\frac{\log(2LM/\delta)}{2n^{(m)}}}\right),$

with probability at least $1-\delta$. Here step $(i)$ uses Eq. (B.3) and step
$(ii)$ uses Eq. (B.3). This completes the proof of the proposition.

###  B.4 Another bound on $\widehat{\Lambda}_{\phi^{\star}}$ under mean
squared distance

In Proposition 1 we established a finite-sample bound for
$\widehat{\Lambda}_{\phi^{\star}}$ when $\mathcal{G}$-divergence is used. If
the IW-CIP penalty is defined by the squared mean distance instead, we provide
another explicit bound for $\widehat{\Lambda}_{\phi^{\star}}$ that converges
to zero as sample size goes to infinity. To show this, we require additional
assumptions on the data generation model and the class of feature mappings
$\Phi$.

###### Proposition 5

Suppose that source and target data are generated under Assumption 1. Let
$\Phi=\\{\phi(x)=Ax+b,A\in{}^{q\times p},b\in{}^{q}\\}$ be the linear class of
feature mappings from p to q $(q\leq p)$, and let the squared mean distance be
the distributional distance in IW-CIP penalty Eq. (19). Suppose $\phi\in\Phi$
is a conditionally invariant feature mapping such that $A\Sigma A^{\top}$ is
non-singular. For any $\delta>0$ satisfying
$\min_{m,y}n^{(m)}{p^{{(m)}}_{y}}\geq 2\log(2ML/\delta)$, we have

$\displaystyle\widehat{\Lambda}_{\phi}$
$\displaystyle\leq\frac{12\lambda_{\textnormal{IW-CIP}}\cdot\lambda_{\rm
max}(A\Sigma{A}^{\top})\cdot\left(p+8\log(2ML/\delta)\right)}{\underset{\begin{subarray}{c}m=1,2,\ldots,M\\\
y=1,2,\ldots,L\end{subarray}}{\min}\
n^{(m)}{p^{{(m)}}_{y}}-2\log(2ML/\delta)},$ (46)

with probability at least $1-\delta$.121212The probability is with respect to
the randomness of
$(X^{(m)}_{k},Y^{(m)}_{k})\stackrel{{\scriptstyle\mathrm{i.i.d.}}}{{\sim}}\mathcal{P}^{(m)}$,
$m=1,\ldots,M,k=1,\ldots,n^{(m)}$; see Eq. (20). Specifically, the statement
holds for the optimal conditionally invariant mapping $\phi^{\star}$ (see Eq.
(5)).

Proof  The exact from of IW-CIP penalty under squared mean distance is

$\displaystyle\widehat{\Lambda}_{\phi}$
$\displaystyle\coloneqq\frac{\lambda_{\textnormal{IW-
CIP}}}{LM^{2}}\sum_{y=1}^{L}\sum_{m\neq
m^{\prime}}\Bigg{\|}\widehat{\mu^{(m)}_{\phi}(y)}-\widehat{\mu^{(m^{\prime})}_{\phi}(y)}\Bigg{\|}_{2}^{2},$

where

$\displaystyle\widehat{\mu^{(m)}_{\phi}(y)}=\frac{1}{n^{(m)}_{y}}\sum_{\begin{subarray}{c}k=1\\\
Y^{(m)}_{k}=y\end{subarray}}^{n^{(m)}}\phi(X^{(m)}_{k}),$

is the empirical mean of $\phi(X)$ on the $m$-th source domain. Since $\phi$
is conditionally invariant, from Eq. (52) we know that for any
$y\in\left\\{1,2,\ldots,L\right\\}$ and for any $m\neq m^{\prime}$,

$\displaystyle A{P_{y}}{v^{{(m)}}_{y}}=A{P_{y}}{v^{{(m^{\prime})}}_{y}}.$

In the finite-sample case, when $Y^{(m)}_{k}=Y^{{({\mathfrak{T}})}}_{k}=y$, we
can write

$\displaystyle\phi(X^{(m)}_{k})=AX^{(m)}_{k}+b$
$\displaystyle=A\left(f^{(1)}(y)+{P_{y}}{v^{{(m)}}_{y}}+\epsilon^{{(m)}}_{k}\right)+b,$
$\displaystyle\phi(X^{{({\mathfrak{T}})}}_{k})=AX^{{({\mathfrak{T}})}}_{k}+b$
$\displaystyle=A\left(f^{(1)}(y)+{P_{y}}{v^{{({\mathfrak{T}})}}_{y}}+\epsilon^{{({\mathfrak{T}})}}_{k}\right)+b,$

and

$\displaystyle\widehat{\mu^{(m)}_{\phi}(y)}-\widehat{\mu^{(m^{\prime})}_{\phi}(y)}$
$\displaystyle=\frac{1}{n^{(m)}_{y}}\sum_{\begin{subarray}{c}k=1\\\
Y^{(m)}_{k}=y\end{subarray}}^{n^{(m)}}\phi(X^{(m)}_{k})-\frac{1}{n^{(m^{\prime})}_{y}}\sum_{\begin{subarray}{c}k=1\\\
Y^{(m^{\prime})}_{k}=y\end{subarray}}^{n^{(m^{\prime})}}\phi(X^{(m^{\prime})}_{k})$
$\displaystyle=A\left(\underbrace{\frac{1}{n^{(m)}_{y}}\sum_{\begin{subarray}{c}k=1\\\
Y^{(m)}_{k}=y\end{subarray}}^{n^{(m)}}\epsilon^{{(m)}}_{k}}_{e^{(m)}_{y}}-\underbrace{\frac{1}{n^{(m^{\prime})}_{y}}\sum_{\begin{subarray}{c}k=1\\\
Y^{(m^{\prime})}_{k}=y\end{subarray}}^{n^{(m^{\prime})}}\epsilon^{(m^{\prime})}_{k}}_{e_{y}^{(m^{\prime})}}\right),$
(47)

where $n^{(m)}_{y}$ is the number of samples with label $y$ in the $m$-th
domain. Since $n_{y}^{(m)}$ follows a multinomial distribution, we conclude
from binomial (multinomial) tail probability bound (e.g. (Chung and Lu, 2006,
Theorem 3.2)) that for all ${t_{y}^{(m)}}>0$,

$\displaystyle\mathbb{P}\left\\{{n^{(m)}_{y}\geq
n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}\right\\}$ $\displaystyle\geq
1-\exp\left(-\frac{n^{(m)}{p^{{(m)}}_{y}}{t_{y}^{(m)}}^{2}}{2}\right).$

When $n^{(m)}_{y}=k>0$ is fixed, we have
$Ae^{(m)}_{y}\sim\mathcal{N}(0,\frac{1}{k}A\Sigma A^{\top})$, and therefore
the standard concentration bound for Gaussian variables (e.g. (Wainwright,
2019, Example 2.11)) implies that for all $t>0$,

$\displaystyle\mathbb{P}\left\\{{\left\|Ae^{(m)}_{y}\right\|_{2}^{2}\leq\frac{1}{k}\lambda_{\rm
max}\left(A\Sigma A^{\top}\right)\cdot p\left(1+t\right)}\ \middle|\
{n^{(m)}_{y}=k}\right\\}\geq 1-e^{-p\min\\{t,t^{2}\\}/8}.$

Next we combine two bounds above into a single bound for
$\left\|Ae^{(m)}_{y}\right\|_{2}^{2}$. Given any $t_{y}^{(m)}>0$ and $t>0$, we
have

$\displaystyle\quad\mathbb{P}\left\\{{\left\|Ae^{(m)}_{y}\right\|_{2}^{2}\leq\frac{\lambda_{\rm
max}\left(A\Sigma A^{\top}\right)\cdot
p\left(1+t\right)}{n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}}\right\\}$
$\displaystyle=\sum_{k=0}^{n^{(m)}}\mathbb{P}\left\\{{\left\|Ae^{(m)}_{y}\right\|_{2}^{2}\leq\frac{\lambda_{\rm
max}\left(A\Sigma A^{\top}\right)\cdot
p\left(1+t\right)}{n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}}\ \middle|\
{n^{(m)}_{y}=k}\right\\}\mathbb{P}\left\\{{n^{(m)}_{y}=k}\right\\}$
$\displaystyle\geq\sum_{k=n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}^{n^{(m)}}\mathbb{P}\left\\{{\left\|Ae^{(m)}_{y}\right\|_{2}^{2}\leq\frac{\lambda_{\rm
max}\left(A\Sigma A^{\top}\right)\cdot
p\left(1+t\right)}{n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}}\ \middle|\
{n^{(m)}_{y}=k}\right\\}\mathbb{P}\left\\{{n^{(m)}_{y}=k}\right\\}$
$\displaystyle\geq\sum_{k=n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}^{n^{(m)}}\mathbb{P}\left\\{{\left\|Ae^{(m)}_{y}\right\|_{2}^{2}\leq\frac{1}{k}\lambda_{\rm
max}\left(A\Sigma A^{\top}\right)\cdot p\left(1+t\right)}\ \middle|\
{n^{(m)}_{y}=k}\right\\}\mathbb{P}\left\\{{n^{(m)}_{y}=k}\right\\}$
$\displaystyle\geq\sum_{k=n^{(m)}{p^{{(m)}}_{y}}(1-{t_{y}^{(m)}})}^{n^{(m)}}\left(1-e^{-p\min\\{t,t^{2}\\}/8}\right)\mathbb{P}\left\\{{n^{(m)}_{y}=k}\right\\}$
$\displaystyle\geq\left(1-e^{-p\min\\{t,t^{2}\\}/8}\right)\left(1-e^{-{n^{(m)}{p^{{(m)}}_{y}}{t_{y}^{(m)}}^{2}}/{2}}\right)$
$\displaystyle\geq
1-e^{-p\min\\{t,t^{2}\\}/8}-e^{-{n^{(m)}{p^{{(m)}}_{y}}{t_{y}^{(m)}}^{2}}/{2}}.$

Take ${t_{y}^{(m)}}=\sqrt{2\log(2ML/\delta)/(n^{(m)}{p^{{(m)}}_{y}})}$ and
$t=8\log(2ML/\delta)/p+\sqrt{8\log(2ML/\delta)/p}$. Then with probability at
least $1-\delta/ML$, we have

$\displaystyle\left\|Ae^{(m)}_{y}\right\|_{2}^{2}$
$\displaystyle\leq\frac{\lambda_{\rm max}\left(A\Sigma
A^{\top}\right)\left(p+8\log(2ML/\delta)+\sqrt{8p\log(2ML/\delta)}\right)}{n^{(m)}{p^{{(m)}}_{y}}-\sqrt{n^{(m)}{p^{{(m)}}_{y}}\cdot
2\log\left(2ML/\delta\right)}}$ $\displaystyle\leq\frac{3\lambda_{\rm
max}\left(A\Sigma
A^{\top}\right)\left(p+8\log(2ML/\delta)\right)}{n^{(m)}{p^{{(m)}}_{y}}-2\log(2ML/\delta)},$

where we apply the inequality $\sqrt{ab}\leq(a+b)/2$ for $a,b>0$ in the last
step. Applying union bound, we get

$\displaystyle\widehat{\Lambda}_{\phi}$
$\displaystyle=\frac{\lambda_{\textnormal{IW-
CIP}}}{LM^{2}}\sum_{y=1}^{L}\sum_{m\neq
m^{\prime}}\Bigg{\|}\widehat{\mu^{(m)}_{\phi}(y)}-\widehat{\mu^{(m^{\prime})}_{\phi}(y)}\Bigg{\|}_{2}^{2}$
$\displaystyle=\frac{\lambda_{\textnormal{IW-
CIP}}}{LM^{2}}\sum_{y=1}^{L}\sum_{m\neq
m^{\prime}}\left\|Ae^{(m)}_{y}-Ae^{(m^{\prime})}_{y}\right\|_{2}^{2}$
$\displaystyle\leq\frac{2\lambda_{\textnormal{IW-
CIP}}}{LM^{2}}\sum_{y=1}^{L}\sum_{m\neq
m^{\prime}}\left(\left\|Ae^{(m)}_{y}\right\|_{2}^{2}+\left\|Ae^{(m^{\prime})}_{y}\right\|_{2}^{2}\right)$
$\displaystyle\leq 4\lambda_{\textnormal{IW-
CIP}}\cdot\underset{\begin{subarray}{c}m=1,2,\ldots,M\\\
y=1,2,\ldots,L\end{subarray}}{\max}\left\|Ae^{(m)}_{y}\right\|_{2}^{2}$
$\displaystyle\leq\frac{12\lambda_{\textnormal{IW-CIP}}\cdot\lambda_{\rm
max}(A\Sigma{A}^{\top})\cdot\left(p+8\log(2ML/\delta)\right)}{\underset{\begin{subarray}{c}m=1,2,\ldots,M\\\
y=1,2,\ldots,L\end{subarray}}{\min}\
n^{(m)}{p^{{(m)}}_{y}}-2\log(2ML/\delta)},$

with probability at least $1-\delta$, proving the result.

###  B.5 Proof of Proposition 2

We have that for any $i\in\left\\{1,2,...,L\right\\}$,

$\displaystyle\mathbb{P}\left\\{{\widehat{h}_{\textnormal{CIP}}(X^{{({\mathfrak{T}})}})=i}\right\\}$
$\displaystyle=\sum_{j=1}^{L}\mathbb{P}\left\\{{\widehat{h}_{\textnormal{CIP}}(X^{{({\mathfrak{T}})}})=i}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j}\right\\}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j}\right\\}$
$\displaystyle=\sum_{j=1}^{L}\left(\mathbb{P}\left\\{{\widehat{h}_{\textnormal{CIP}}(X^{(m)})=i}\
\middle|\
{Y^{(m)}=j}\right\\}+\delta_{i,j}\right)\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j}\right\\},$

where
$\delta_{i,j}=\mathbb{P}\left\\{{\widehat{h}_{\textnormal{CIP}}(X^{{({\mathfrak{T}})}})=i}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j}\right\\}-\mathbb{P}\left\\{{\widehat{h}_{\textnormal{CIP}}(X^{(1)})=i}\
\middle|\ {Y^{(1)}=j}\right\\}$. Writing the above equation in matrix form, we
have

$\displaystyle\mu_{\widehat{h}_{\textnormal{CIP}}}=C_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}w^{(m)}+\Delta\mu^{({\mathfrak{T}})},$
(48)

where $\Delta$ is a matrix with $\delta_{i,j}$ being its $(i,j)$-th element,
and $\mu^{({\mathfrak{T}})}$ denotes the target label distribution. Comparing
with the linear system in Eq. (15), Eq. (48) has an additional term on the
right-hand side because the finite-sample CIP $\widehat{h}_{\textnormal{CIP}}$
does not guarantee exact conditional invariance. Since the importance weight
$\widehat{w}^{(m)}$ is obtained by solving the finite-sample version

$\displaystyle\widehat{\mu}_{\widehat{h}_{\textnormal{CIP}}}=\widehat{C}_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}\widehat{w}^{(m)},$
(49)

combining Eq. (48) and Eq. (49), we have

$\displaystyle\begin{split}&\underbrace{\quad\mu_{\widehat{h}_{\textnormal{CIP}}}-\widehat{\mu}_{\widehat{h}_{\textnormal{CIP}}}}_{(A_{1})}\\\
&=C_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}\left(w^{(m)}-\widehat{w}^{(m)}\right)+\underbrace{\left(C_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}-\widehat{C}_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}\right)}_{(A_{2})}\widehat{w}^{(m)}+\underbrace{\Delta\mu^{\Diamond}}_{(A_{3})}.\end{split}$
(50)

Now we bound $(A_{1}),(A_{2})$, and $(A_{3})$ separately. For bouding
$(A_{1})$ and $(A_{2})$, we closely follow the proof of Theorem 3 in Lipton et
al. (2018).

##### Term ($A_{1}$)

By Hoeffding’s inequality, we have

$\displaystyle\mathbb{P}\left\\{{\left|A_{1,i}\right|>t}\right\\}$
$\displaystyle=\mathbb{P}\left\\{{\left|\frac{1}{n^{({\mathfrak{T}})}}\sum_{k=1}^{n^{({\mathfrak{T}})}}\mathbf{1}_{{\widehat{h}_{\textnormal{CIP}}(\widetilde{X}^{({\mathfrak{T}})}_{k})=i}}-\mathbb{E}\left[{\mathbf{1}_{{\widehat{h}_{\textnormal{CIP}}({X}^{({\mathfrak{T}})})=i}}}\right]\right|>t}\right\\}\leq
2\exp\left\\{-2n^{({\mathfrak{T}})}t^{2}\right\\},$

for all $t>0$. Therefore, for any $t_{1}>0$, we have

$\displaystyle\mathbb{P}\left\\{{\lVert{A_{1}}\rVert_{2}>t_{1}}\right\\}$
$\displaystyle\leq\mathbb{P}\left\\{{\exists
i\in\left\\{1,...,L\right\\},\left|A_{1,i}\right|>\frac{t_{1}}{\sqrt{L}}}\right\\}$
$\displaystyle\leq\sum_{i=1}^{L}\mathbb{P}\left\\{{\left|A_{1,i}\right|>\frac{t_{1}}{\sqrt{L}}}\right\\}$
$\displaystyle\leq
2L\cdot\exp\left\\{-\frac{2n^{({\mathfrak{T}})}t_{1}^{2}}{L}\right\\}.$

##### Term ($A_{2}$)

Note that

$\displaystyle
C_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}-\widehat{C}_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}$
$\displaystyle=\frac{1}{n^{(m)}}\sum_{k=1}^{n^{(m)}}e_{\widehat{h}_{\textnormal{CIP}}(\widetilde{X}^{(m)}_{k})}e_{\widetilde{Y}^{(m)}_{k}}^{\top}-\mathbb{E}\left[{e_{\widehat{h}_{\textnormal{CIP}}(X^{(m)})}{e}_{Y^{(m)}}^{\top}}\right],$

where $e_{i}$ denotes the zero vector with only the $i$-th element being 1.
Let

$\displaystyle
Z_{k}=e_{\widehat{h}_{\textnormal{CIP}}(\widetilde{X}^{(m)}_{k})}e_{\widetilde{Y}^{(m)}_{k}}^{\top}-\mathbb{E}\left[{e_{\widehat{h}_{\textnormal{CIP}}(X^{(m)})}{e}_{Y^{(m)}}^{\top}}\right].$

Then it is straightforward to see that

$\displaystyle\begin{cases}\mathbb{E}\left[{Z_{k}}\right]=0;\\\
\lVert{Z_{k}}\rVert_{2}\leq 2;\\\
\max\left\\{\lVert{\mathbb{E}\left[{Z_{k}Z_{k}^{\top}}\right]}\rVert_{2},\lVert{\mathbb{E}\left[{Z_{k}^{\top}Z_{k}}\right]}\rVert_{2}\right\\}\leq
1.\end{cases}$

By the matrix Bernstein inequality (see e.g. (Tropp et al., 2015, Theorem
6.1.1)), we have

$\displaystyle\mathbb{P}\left\\{{\lVert{A_{2}}\rVert_{2}\geq
t_{2}}\right\\}\leq
2L\exp\left\\{-\frac{n^{(m)}t_{2}^{2}/2}{1+2t_{2}/3}\right\\},$

for all $t_{2}>0$.

##### Term ($A_{3}$)

By definition of $D_{\mathcal{G}}\left({\cdot},{\cdot}\right)$ and
$\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{CIP}}}$, each element of
$A_{3}$ satisfies

$\displaystyle\left|A_{3,i}\right|$
$\displaystyle\leq\sum_{j=1}^{L}\left|\delta_{i,j}\right|\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j}\right\\}$
$\displaystyle\leq\sum_{j=1}^{L}D_{\mathcal{G}}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{\widehat{\phi}_{\textnormal{CIP}}(X)}|{Y=j}}},{\mathcal{P}^{(m)}_{{\widehat{\phi}_{\textnormal{CIP}}(X)}|{Y=j}}}\right)\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j}\right\\}\leq\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{CIP}}}.$

Hence we have

$\displaystyle\lVert{A_{3}}\rVert_{2}=\sqrt{\sum_{i=1}^{L}\left|A_{3,i}\right|^{2}}\leq\sqrt{L}\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{CIP}}}.$

Applying these bounds to Eq. (50) and using the assumption that
$C_{\widehat{h}_{\textnormal{CIP}}}^{{(m)}}$ has condition number
$\kappa^{{(m)}}$ (i.e., the smallest eigenvalue $\geq 1/\kappa^{{(m)}}$
because $1$ is an eigenvalue of the matrix), we get

$\displaystyle\left\|w^{(m)}-\widehat{w}^{(m)}\right\|_{2}$
$\displaystyle\leq\kappa_{m}\left(t_{1}+t_{2}\cdot\left\|\widehat{w}^{(m)}\right\|_{2}+\sqrt{L}\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{CIP}}}\right)$
$\displaystyle\leq\kappa_{m}\left(t_{1}+t_{2}\left(\left\|w^{(m)}-\widehat{w}^{(m)}\right\|_{2}+\left\|w^{(m)}\right\|_{2}\right)+\sqrt{L}\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{CIP}}}\right),$

with probability at least
$1-2L\exp\left(-2n^{({\mathfrak{T}})}t_{1}^{2}/L\right)-2L\exp\left(-3n^{(m)}t_{2}^{2}/(6+4t_{2})\right)$.
Finally, let $t_{1}=\sqrt{L\log(4L/\delta)/2n^{({\mathfrak{T}})}}$ and
$t_{2}=\sqrt{3\log(4L/\delta)/n^{(m)}}$. Using the assumption $n^{(m)}\geq
12\kappa_{m}^{2}\log(4L/\delta)$, we get

$\displaystyle\left\|w^{(m)}-\widehat{w}^{(m)}\right\|_{2}$ $\displaystyle\leq
2\kappa_{m}\left(\sqrt{L}\Psi_{\mathcal{G},\widehat{\phi}_{\textnormal{CIP}}}+\sqrt{\frac{L\log(4L/\delta)}{2n^{({\mathfrak{T}})}}}+\sqrt{\frac{3\log(4L/\delta)}{n^{(m)}}}\left\|w^{(m)}\right\|_{2}\right)$

with probability at least $1-\delta$. Combining the results then proves the
proposition.

###  B.6 Proof of Proposition 3

Using assumption (25), for any positive semidefinite matrix
$S\in{}^{d_{y}\times d_{y}}$, we have

$\displaystyle\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)^{\top}S\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)$
$\displaystyle=\operatorname{trace}\left(\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)^{\top}S^{\frac{1}{2}}S^{\frac{1}{2}}\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)\right)$
$\displaystyle=\operatorname{trace}\left(S^{\frac{1}{2}}\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)^{\top}S^{\frac{1}{2}}\right)$
$\displaystyle\leq\operatorname{trace}\left(S^{\frac{1}{2}}\left(\frac{\zeta_{y}}{M-1}\sum_{m=2}^{M}{v^{{(m)}}_{y}}({v^{{(m)}}_{y}})^{\top}\right)S^{\frac{1}{2}}\right)$
$\displaystyle=\frac{\zeta_{y}}{M-1}\sum_{m=2}^{M}\operatorname{trace}\left(S^{\frac{1}{2}}{v^{{(m)}}_{y}}({v^{{(m)}}_{y}})^{\top}S^{\frac{1}{2}}\right)$
$\displaystyle=\frac{\zeta_{y}}{M-1}\sum_{m=2}^{M}{{v^{{(m)}}_{y}}}^{\top}S{v^{{(m)}}_{y}}.$
(51)

Inequality above uses the fact that
$\operatorname{trace}(A)\leq\operatorname{trace}(B)$ for positive semi-
definite matrix $A,B$ satisfying $A\preceq B$. Given that
$\epsilon^{{(m)}},\epsilon^{{({\mathfrak{T}})}}\stackrel{{\scriptstyle\mathrm{i.i.d.}}}{{\sim}}\mathcal{N}(0,\Sigma)$,
for any $y\in\left\\{1,2,\ldots,L\right\\}$ and
$m\in\left\\{1,2,\ldots,M\right\\}$ we have

$\displaystyle\begin{split}\phi(X^{(m)})\mid
Y^{(m)}=y&\sim\mathcal{N}\left(A\left(f^{(1)}(y)+{P_{y}}{v^{{(m)}}_{y}}\right)+b,A\Sigma
A^{\top}\right),\\\ \phi(X^{{({\mathfrak{T}})}})\mid
Y^{{({\mathfrak{T}})}}=y&\sim\mathcal{N}\left(A\left(f^{(1)}(y)+{P_{y}}{v^{{({\mathfrak{T}})}}_{y}}\right)+b,A\Sigma
A^{\top}\right).\end{split}$ (52)

Note that ${v^{{(1)}}_{y}}=0$ by Assumption 1, so we can further simplify
terms in Eq. (23) as

$\displaystyle{{\Delta}^{(m)}_{\phi}(y)}=AP_{y}\left({v^{{(m)}}_{y}}-{v^{{(1)}}_{y}}\right)=AP_{y}{v^{{(m)}}_{y}}\
\text{ and }\ {\Sigma_{\phi}^{(m)}}(y)=A\Sigma A^{\top}.$

With results above in hand, we can derive

$\displaystyle\quad
D_{\mathcal{G}}^{2}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\right)$
$\displaystyle\leq\sup_{g\in\mathcal{G}}\max_{y^{\prime}=1,2,\ldots,L}\left|\mathbb{E}_{Z\sim\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}}\left[{\mathbf{1}_{{g(Z)=y^{\prime}}}}\right]-\mathbb{E}_{Z\sim\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\left[{\mathbf{1}_{{g(Z)=y^{\prime}}}}\right]\right|$
$\displaystyle\leq\sup_{g\in\mathcal{G}}\sum_{y^{\prime}=1}^{L}\left|\mathbb{E}_{Z\sim\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}}\left[{\mathbf{1}_{{g(Z)=y^{\prime}}}}\right]-\mathbb{E}_{Z\sim\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\left[{\mathbf{1}_{{g(Z)=y^{\prime}}}}\right]\right|$
$\displaystyle=4\sup_{g\in\mathcal{G}}{\rm
TV}^{2}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{g(\phi(X))}|{Y=y}}},{\mathcal{P}^{(m)}_{{g(\phi(X))}|{Y=y}}}\right)$
$\displaystyle\overset{(i)}{\leq}2\sup_{g\in\mathcal{G}}{\rm
KL}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{g(\phi(X))}|{Y=y}}},{\mathcal{P}^{(m)}_{{g(\phi(X))}|{Y=y}}}\right)$
$\displaystyle\overset{(ii)}{\leq}2{\rm
KL}\left({\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi(X)}|{Y=y}}},{\mathcal{P}^{(m)}_{{\phi(X)}|{Y=y}}}\right)$
$\displaystyle\overset{(iii)}{=}\left(A{P_{y}}\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)\right)^{\top}\left(A\Sigma
A^{\top}\right)^{-1}\left(A{P_{y}}\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)\right)$
$\displaystyle\overset{(iv)}{\leq}\frac{\zeta_{y}}{M-1}\sum_{m=2}^{M}\left(A{P_{y}}{v^{{(m)}}_{y}}\right)^{\top}\left(A\Sigma
A^{\top}\right)^{-1}\left(A{P_{y}}{v^{{(m)}}_{y}}\right)$
$\displaystyle\overset{(v)}{=}2\zeta_{y}{\Pi_{\phi}(y)}.$ (53)

Here step $(i)$ applies Pinsker’s inequality and step $(ii)$ applies data
processing inequality. Step $(iii)$ applies the KL divergence formula between
two Gaussians, i.e.,

$\displaystyle{\rm
KL}\left({\mathcal{N}(\mu_{1},\Sigma_{1})},{\mathcal{N}(\mu_{1},\Sigma_{2})}\right)$
$\displaystyle=\frac{1}{2}\left(\log\frac{\left|\Sigma_{2}\right|}{\Sigma_{1}}\right)-p+\text{Tr}\left(\Sigma_{2}^{-1}\Sigma_{1}\right)+(\mu_{1}-\mu_{2})^{\top}\Sigma_{2}^{-1}(\mu_{2}-\mu_{1}).$

Step $(iv)$ uses Eq. (B.6) with $S=P_{y}^{\top}A^{\top}(A\Sigma
A^{\top})^{-1}AP_{y}$ and the last step $(vi)$ follows from definition (24).
Now the proposition follows by applying Eq. (B.6) to the definition of
$\Psi_{\mathcal{G},\phi}$ in Eq. (21).

###  B.7 Number of source domains needed for controlling the deviation from
conditional invariance

To validate condition (25) in Proposition 3, we prove the following lemma
which states that when the perturbations are generated randomly, the required
number of source domains is approximately linear in the dimension of
perturbations.

###### Lemma 1

Under Assumption 1, further assume that perturbations are generated according
to
${v^{{(m)}}_{y}},{v^{{({\mathfrak{T}})}}_{y}}\stackrel{{\scriptstyle\mathrm{i.i.d.}}}{{\sim}}\mathcal{N}(0,\tau_{y}^{2}\mathbb{I}_{d_{y}})$.
For any $\delta\in(0,1)$, if the number of source domains $M$ satisfies

$\displaystyle M\geq
1+r\left(\max_{y}\sqrt{d_{y}}+\sqrt{2\log(2L/\delta)}\right)^{2},$

for some constant $r>1$, then with probability at least $1-\delta$,131313The
probability is with respect to the randomness of perturbations
${v^{{(m)}}_{y}},{v^{{({\mathfrak{T}})}}_{y}}\stackrel{{\scriptstyle\mathrm{i.i.d.}}}{{\sim}}\mathcal{N}(0,\tau_{y}^{2}\mathbb{I}_{d_{y}})$.
condition (25) holds with $\zeta_{y}=c_{r}(M+\log M)$ for some constant
$c_{r}>0$.

Proof  We verify that when
${v^{{(m)}}_{y}},{v^{{({\mathfrak{T}})}}_{y}}\stackrel{{\scriptstyle\mathrm{i.i.d.}}}{{\sim}}\mathcal{N}(0,\tau_{y}^{2}\mathbb{I}_{d_{y}})$,
condition (25) holds with high probability. Sincet
$M\geq\max_{y=1,2,\ldots,L}\\{d_{y}\\}+1$, Theorem 6.1 and Example 2.11 in
Wainwright (2019) show that given $y\in\\{1,2,\ldots,L\\}$, for any
$t_{0}\in(0,1),t_{1},\ldots,t_{M}>0$ we have

$\displaystyle\mathbb{P}\left\\{{\lambda_{\rm
min}\left(\frac{1}{M-1}\sum_{m=2}^{M}{v^{{(m)}}_{y}}{{v^{{(m)}}_{y}}}^{\top}\right)\geq\tau_{y}^{2}\left(1-t_{0}-\sqrt{\frac{d_{y}}{M-1}}\right)}\right\\}$
$\displaystyle\geq 1-e^{-(M-1)t_{0}^{2}/2},$

and

$\displaystyle\mathbb{P}\left\\{{\left\|{v^{{({\mathfrak{T}})}}_{y}}-{v^{{(1)}}_{y}}\right\|_{2}^{2}\leq\tau_{y}^{2}d_{y}\left(1+t_{1}\right)}\right\\}\geq
1-e^{-d_{y}\min\\{t_{1},t_{1}^{2}\\}/8},$
$\displaystyle\mathbb{P}\left\\{{\left\|{v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right\|_{2}^{2}\leq
2\tau_{y}^{2}d_{y}\left(1+t_{m}\right)}\right\\}\geq
1-e^{-d_{y}\min\\{t_{m},t_{m}^{2}\\}/8},\ \ \ \forall
m\in\left\\{2,\cdots,M\right\\}.$

Take $t_{0}=\sqrt{\frac{2}{M-1}\log\frac{2L}{\delta}}$ and
$t_{m}=\frac{8}{d_{y}}\log(\frac{2ML}{\delta})+\sqrt{\frac{8}{d_{y}}\log(\frac{2ML}{\delta})}$
for all $y$ such that $d_{y}>0$. Applying union bound, the following bound
holds for all $y\in\left\\{1,2,\ldots,L\right\\}$ and
$m\in\left\\{1,2,\ldots,M\right\\}$ with probability greater than $1-\delta$:

$\displaystyle\quad\lambda_{\rm
max}\left(\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)^{\top}\right)=\left\|\left({v^{{({\mathfrak{T}})}}_{y}}-{v^{{(m)}}_{y}}\right)\right\|_{2}^{2}$
$\displaystyle\leq
2\tau_{y}^{2}\left(d_{y}+8\log(2ML/\delta)+\sqrt{8d_{y}\log(2ML/\delta)}\right)$
$\displaystyle\leq 3\tau_{y}^{2}\left(d_{y}+8\log(2ML/\delta)\right)$
$\displaystyle\leq\frac{3\left(d_{y}+8\log(2ML/\delta)\right)}{1-\sqrt{{2}\left(\log(2L/{\delta})\right)/(M-1)}-\sqrt{{d_{y}}/(M-1)}}\lambda_{\rm
min}\left(\frac{1}{M-1}\sum_{m=2}^{M}{v^{{(m)}}_{y}}{v^{{(m)}}_{y}}^{\top}\right)$
$\displaystyle\leq\frac{3r}{r-1}\left(d_{y}+8\log\left(\frac{2ML}{\delta}\right)\right)\lambda_{\rm
min}\left(\frac{1}{M-1}\sum_{m=2}^{M}{v^{{(m)}}_{y}}{v^{{(m)}}_{y}}^{\top}\right)$
$\displaystyle\leq c_{r}(M+\log M)\cdot\lambda_{\rm
min}\left(\frac{1}{M-1}\sum_{m=2}^{M}{v^{{(m)}}_{y}}{v^{{(m)}}_{y}}^{\top}\right),$

where we apply the lower bound on $M$ in the last two steps.

###  B.8 Proof of Proposition 4

Let $f_{\textnormal{inv}}(1),f_{\textnormal{inv}}(2)\in{}^{p-d}$ denote the
vector of the last $p-d$ coordinates of $f^{(1)}(1),f^{(1)}(2)\in{}^{p}$
respectively, which are conditionally invariant across source and target. We
first prove that the target risks of $h_{\textnormal{oracle}}$ and $h^{\star}$
are given by

$\displaystyle\begin{split}\mathcal{R}^{{({\mathfrak{T}})}}({h_{\textnormal{oracle}}})&=r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}\left(\left\|f^{{({\mathfrak{T}})}}(1)-f^{{({\mathfrak{T}})}}(2)\right\|_{2}\right),\\\
\mathcal{R}^{{({\mathfrak{T}})}}({h^{\star}})&=r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}\Big{(}\left\|f_{\textnormal{inv}}(1)-f_{\textnormal{inv}}(2)\right\|_{2}\Big{)},\end{split}$
(54)

where $r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}(x)$ is a continuous decreasing
function in $x$. The risk of the oracle classifier can be derived as follows.
For simplicity we write $\mu^{({\mathfrak{T}})}_{1}=f^{{({\mathfrak{T}})}}(1)$
and $\mu^{({\mathfrak{T}})}_{2}=f^{{({\mathfrak{T}})}}(2)$. For any
$\phi(x)=\beta^{\top}x+\beta_{0}$ with $\left\|\beta\right\|_{2}=1$, we know
that

$\displaystyle\phi(X^{{({\mathfrak{T}})}})\mid Y^{{({\mathfrak{T}})}}=1$
$\displaystyle\sim\mathcal{N}\left(\beta^{\top}\mu^{({\mathfrak{T}})}_{1}+\beta_{0},\sigma^{2}\right),$
$\displaystyle\phi(X^{{({\mathfrak{T}})}})\mid Y^{{({\mathfrak{T}})}}=2$
$\displaystyle\sim\mathcal{N}\left(\beta^{\top}\mu^{({\mathfrak{T}})}_{2}+\beta_{0},\sigma^{2}\right).$

Then the target risk of $h=g\circ\phi$ can be computed by

$\displaystyle\mathcal{R}^{{({\mathfrak{T}})}}({h})$
$\displaystyle=\mathbb{E}\left[{\mathbf{1}_{{h(X^{{({\mathfrak{T}})}})\neq
Y^{{({\mathfrak{T}})}}}}}\right]$
$\displaystyle=\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=2}\ \middle|\
{Y^{{({\mathfrak{T}})}}=1}\right\\}\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=1}\right\\}$
$\displaystyle\hskip
150.00023pt+\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=1}\ \middle|\
{Y^{{({\mathfrak{T}})}}=2}\right\\}\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=2}\right\\}$
$\displaystyle=\mathbb{P}\left\\{{\phi(X^{{({\mathfrak{T}})}})>0}\ \middle|\
{Y^{{({\mathfrak{T}})}}=1}\right\\}\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=1}\right\\}$
$\displaystyle\hskip
150.00023pt+\mathbb{P}\left\\{{\phi(X^{{({\mathfrak{T}})}})\leq 0}\ \middle|\
{Y^{{({\mathfrak{T}})}}=2}\right\\}\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=2}\right\\}$
$\displaystyle={p^{{({\mathfrak{T}})}}_{1}}\cdot\frac{1}{2}\left(1-{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{1}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)\right)+(1-{p^{{({\mathfrak{T}})}}_{1}})\cdot\frac{1}{2}\left(1+{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{2}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)\right)$
$\displaystyle=\frac{1}{2}\left(1-{p^{{({\mathfrak{T}})}}_{1}}\cdot{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{1}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)+(1-{p^{{({\mathfrak{T}})}}_{1}})\cdot{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{2}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)\right),$
(55)

where ${\rm erf}\left(x\right)=2\int_{0}^{x}e^{-t^{2}}dt/\sqrt{\pi}$ is the
Gauss error function. The optimization problem of $h_{\textnormal{oracle}}$ in
Eq. (3) becomes

$\displaystyle\begin{split}\beta_{\rm orcale},\beta_{0,{\rm
orcale}}=\operatorname*{arg\,min}_{\beta,\beta_{0}}&\hskip
10.00002pt-{p^{{({\mathfrak{T}})}}_{1}}\cdot{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{1}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)+(1-{p^{{({\mathfrak{T}})}}_{1}})\cdot{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{2}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right),\\\
\text{s.t.}&\hskip 10.00002pt\left\|\beta\right\|_{2}=1.\end{split}$ (56)

Setting the gradient of the objective function with respect to $\beta_{0}$ to
$0$, we get

$\displaystyle\beta_{0}=\frac{\sigma^{2}}{\beta^{\top}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2})}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}-\frac{1}{2}\beta^{\top}(\mu^{({\mathfrak{T}})}_{1}+\mu^{({\mathfrak{T}})}_{2}).$

Plugging the expression of $\beta_{0}$ into Eq. (56) and using the method of
Lagrange multipliers, we arrive at

$\displaystyle\frac{\partial}{\partial\beta}\left(-{p^{{({\mathfrak{T}})}}_{1}}\cdot{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{1}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)+(1-{p^{{({\mathfrak{T}})}}_{1}})\cdot{\rm
erf}\left(\frac{\beta^{\top}\mu^{({\mathfrak{T}})}_{2}+\beta_{0}}{\sqrt{2\sigma^{2}}}\right)+\lambda(\left\|\beta\right\|_{2}^{2}-1)\right)=0$
$\displaystyle\Longrightarrow$
$\displaystyle\;\;\;2\lambda\beta={p^{{({\mathfrak{T}})}}_{1}}\frac{2}{\sqrt{\pi}}\exp\left(-\frac{(\beta^{\top}\mu^{({\mathfrak{T}})}_{1}+\beta_{0})^{2}}{2\sigma^{2}}\right)$
$\displaystyle\hskip
100.00015pt\times\left(\frac{1}{2\sqrt{2\sigma^{2}}}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2})-\frac{\sqrt{\sigma^{2}/2}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}}{(\beta^{\top}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}))^{2}}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2})\right)$
$\displaystyle\hskip
20.00003pt-(1-{p^{{({\mathfrak{T}})}}_{1}}){p^{{({\mathfrak{T}})}}_{1}}\frac{2}{\sqrt{\pi}}\exp\left(-\frac{(\beta^{\top}\mu^{({\mathfrak{T}})}_{2}+\beta_{0})^{2}}{2\sigma^{2}}\right)$
$\displaystyle\hskip
100.00015pt\times\left(-\frac{1}{2\sqrt{2\sigma^{2}}}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2})-\frac{\sqrt{\sigma^{2}/2}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}}{(\beta^{\top}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}))^{2}}(\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2})\right).$

That is, $\beta$ and $\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}$
are in the same direction. Using $\left\|\beta\right\|_{2}=1$, we get

$\displaystyle\beta_{\rm oracle}$
$\displaystyle=\frac{\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}}{\left\|\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}\right\|_{2}},$
$\displaystyle\beta_{0,{\rm oracle}}$
$\displaystyle=\frac{\sigma^{2}}{\left\|\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}\right\|_{2}}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}-\frac{\left\|\mu^{({\mathfrak{T}})}_{1}\right\|_{2}^{2}-\left\|\mu^{({\mathfrak{T}})}_{2}\right\|_{2}^{2}}{2\left\|\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}\right\|_{2}}.$

Plugging the above equations into Eq. (55), after some algebra we get the
oracle target risk as

$\displaystyle\mathcal{R}^{{({\mathfrak{T}})}}({h_{\textnormal{oracle}}})$
$\displaystyle=r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}\left(\left\|\mu^{({\mathfrak{T}})}_{1}-\mu^{({\mathfrak{T}})}_{2}\right\|_{2}\right)=r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}\left(\left\|f^{{({\mathfrak{T}})}}(1)-f^{{({\mathfrak{T}})}}(2)\right\|_{2}\right),$

where

$\displaystyle\begin{split}r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}(x)&=\frac{1}{2}\Bigg{(}1-{p^{{({\mathfrak{T}})}}_{1}}{\rm
erf}\left(\frac{x}{2\sqrt{2}\sigma}+\frac{\sigma}{\sqrt{2}x}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)\\\
&\hskip 40.00006pt-(1-{p^{{({\mathfrak{T}})}}_{1}}){\rm
erf}\left(\frac{x}{2\sqrt{2}\sigma}-\frac{\sigma}{\sqrt{2}x}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)\Bigg{)}.\end{split}$
(57)

Next, turning to $h^{\star}$, note that
$span\\{{v^{{(2)}}_{1}},{v^{{(3)}}_{1}},\ldots,{v^{{(M)}}_{1}},{v^{{(2)}}_{2}},{v^{{(3)}}_{2}},\ldots,{v^{{(M)}}_{2}}\\}={}^{d}$,
and so the constraint of conditional invariance in Eq. (5) requires that the
first $d$ coordinates of $\beta^{\star}$ are zero, i.e., $\phi^{\star}$ only
uses the last $p-d$ invariant coordinates. Following exactly the same proof as
above, we have

$\displaystyle\mathcal{R}^{{({\mathfrak{T}})}}({h^{\star}})=r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}\Big{(}\left\|f_{\textnormal{inv}}(1)-f_{\textnormal{inv}}(2)\right\|_{2}\Big{)}.$

Now we verify that $r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}(x)$ is a
decreasing function in $x$, by calculating its gradient.

$\displaystyle\quad\frac{dr_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}(x)}{dx}$
$\displaystyle=-\frac{1}{\sqrt{\pi}}{p^{{({\mathfrak{T}})}}_{1}}\exp\left(-\left(\frac{x}{2\sqrt{2}\sigma}+\frac{\sigma}{\sqrt{2}x}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)^{2}\right)$
$\displaystyle\hskip
200.0003pt\times\left(\frac{1}{2\sqrt{2}\sigma}-\frac{\sigma}{\sqrt{2}x^{2}}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)$
$\displaystyle\hskip
10.00002pt-\frac{1}{\sqrt{\pi}}(1-{p^{{({\mathfrak{T}})}}_{1}})\exp\left(-\left(\frac{x}{2\sqrt{2}\sigma}-\frac{\sigma}{\sqrt{2}x}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)^{2}\right)$
$\displaystyle\hskip
200.0003pt\times\left(\frac{1}{2\sqrt{2}\sigma}+\frac{\sigma}{\sqrt{2}x^{2}}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)$
$\displaystyle=-\sqrt{\frac{{p^{{({\mathfrak{T}})}}_{1}}(1-{p^{{({\mathfrak{T}})}}_{1}})}{2\pi\sigma^{2}}}\exp\left(-\left(\frac{x}{2\sqrt{2}\sigma}\right)^{2}-\left(\frac{\sigma}{\sqrt{2}x}\log\frac{{p^{{({\mathfrak{T}})}}_{1}}}{1-{p^{{({\mathfrak{T}})}}_{1}}}\right)^{2}\right)<0.$

This intermediate result implies that the risk difference between $h^{\star}$
and $h_{\textnormal{oracle}}$ originates from the norm differences between
$f^{{({\mathfrak{T}})}}(1)-f^{{({\mathfrak{T}})}}(2)$ and
$f_{\textnormal{inv}}(1)-f_{\textnormal{inv}}(2)$ only, where the second term
is equivalent to the last $p-d$ dimensions of the first term which corresponds
to CICs in the general anticausal model. Intuitively, when the dimension of
CICs is close to $p$, the discrepancy between the norms of
$f^{{({\mathfrak{T}})}}(1)-f^{{({\mathfrak{T}})}}(2)$ and
$f_{\textnormal{inv}}(1)-f_{\textnormal{inv}}(2)$ becomes insignificant, and
the disparity in the target risks of $h^{\star}$ and $h_{\textnormal{oracle}}$
tends to be negligible. In particular, we establish the explicit bound for
this disparity as follows:

When ${p^{{({\mathfrak{T}})}}_{1}}=1/2$, the explicit formula for
$r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}(x)$ given in Eq.(57) shows that
$r_{{p^{{({\mathfrak{T}})}}_{1}},\sigma}(x)=\left(1-{\rm
erf}\left(x/(2\sqrt{2}\sigma)\right)\right)/2$. Since
${v^{{({\mathfrak{T}})}}_{1}}-{v^{{({\mathfrak{T}})}}_{2}}\sim\mathcal{N}(0,2\tau^{2}\sigma^{2}\mathbb{I}_{d})$,
the standard Gaussian concentration bound shows that

$\displaystyle\mathbb{P}\left\\{{\left\|{v^{{({\mathfrak{T}})}}_{1}}-{v^{{({\mathfrak{T}})}}_{2}}\right\|_{2}\leq\sqrt{2\tau^{2}\sigma^{2}d\left(1+\sqrt{\frac{8}{d}\log\frac{1}{\delta}}+\frac{8}{d}\log\frac{1}{\delta}\right)}}\right\\}\geq
1-\delta.$

By using the target risks given in Eq. (54), we can calculate

$\displaystyle\quad\mathcal{R}^{{({\mathfrak{T}})}}({h^{\star}})-\mathcal{R}^{{({\mathfrak{T}})}}({h_{\textnormal{oracle}}})$
$\displaystyle=\frac{1}{2}\left({\rm
erf}\left(\frac{\left\|f^{{({\mathfrak{T}})}}(1)-f^{{({\mathfrak{T}})}}(2)\right\|_{2}}{2\sqrt{2}\sigma}\right)-{\rm
erf}\left(\frac{\left\|f_{\textnormal{inv}}(1)-f_{\textnormal{inv}}(2)\right\|_{2}}{2\sqrt{2}\sigma}\right)\right)$
$\displaystyle\leq\frac{1}{2}\left({\rm
erf}\left(\frac{\left\|f^{(1)}(1)-f^{(1)}(2)\right\|_{2}+\left\|{v^{{({\mathfrak{T}})}}_{1}}-{v^{{({\mathfrak{T}})}}_{2}}\right\|_{2}}{2\sqrt{2}\sigma}\right)-{\rm
erf}\left(\frac{\left\|f_{\textnormal{inv}}(1)-f_{\textnormal{inv}}(2)\right\|_{2}}{2\sqrt{2}\sigma}\right)\right)$
$\displaystyle\overset{(i)}{=}\frac{1}{2}\left({\rm
erf}\left(\frac{\xi\sigma\sqrt{p}+\left\|{v^{{({\mathfrak{T}})}}_{1}}-{v^{{({\mathfrak{T}})}}_{2}}\right\|_{2}}{2\sqrt{2}\sigma}\right)-{\rm
erf}\left(\frac{\xi\sqrt{p-d}}{2\sqrt{2}}\right)\right)$
$\displaystyle\leq\frac{1}{2}\left({\rm
erf}\left(\frac{\xi\sqrt{p}+\tau\sqrt{3d+24\log(1/\delta)}}{2\sqrt{2}}\right)-{\rm
erf}\left(\frac{\xi\sqrt{p-d}}{2\sqrt{2}}\right)\right)$
$\displaystyle\overset{(ii)}{\leq}\frac{\xi(\sqrt{p}-\sqrt{p-d})+\tau\sqrt{3d+24\log(1/\delta)}}{\sqrt{8\pi}}\exp\left(-\frac{\xi^{2}(p-d)}{8}\right)$
$\displaystyle=\frac{1}{\sqrt{8\pi}}\left(\frac{cd}{\sqrt{p}+\sqrt{p-d}}+\tau\sqrt{3d+24\log\left(\frac{1}{\delta}\right)}\right)\exp\left(-\frac{\xi^{2}(p-d)}{8}\right)$
$\displaystyle\leq\frac{1}{\sqrt{8\pi}}\left({\xi}{\sqrt{d}}+\tau\sqrt{3d+24\log\left(\frac{1}{\delta}\right)}\right)\exp\left(-\frac{\xi^{2}(p-d)}{8}\right)$
$\displaystyle\leq
c_{\xi,\tau}\left(\sqrt{d}+\sqrt{\log(1/\delta)}\right)\exp\left(-\frac{\xi^{2}(p-d)}{8}\right),$

with probability at least $1-\delta$. Here, in step $(i)$, we use the fact
that $f_{\textnormal{inv}}(1)$ and $f_{\textnormal{inv}}(2)$ are the last
$p-d$ invariant coordinates of $f^{(1)}(1)$ and $f^{(1)}(2)$ respectively; and
in step $(ii)$ we use the inequality $\int_{a}^{b}e^{-t^{2}}dt\leq
e^{-a^{2}}(b-a)$ for $0<a<b$. This completes the proof of the proposition.

##  Appendix C Proof of Section 4

In this section we prove our theorems given in Section 4.

###  C.1 Proof of Theorem 2

We prove the following generalized version of the theorem, which does not
require that $\phi_{\textnormal{inv}}$ is exactly conditionally invariant.
That is, for any
$h_{\textnormal{inv}}=g_{\textnormal{inv}}\circ\phi_{\textnormal{inv}}$ where
$g_{\textnormal{inv}}\in\mathcal{G}$ and $\phi_{\textnormal{inv}}\in\Phi$, for
any other classifier $h$ we have

$\left|\mathcal{R}^{{({\mathfrak{T}})}}({h})-\mathcal{R}^{(1)}({h})\right|\leq
2\mathcal{R}^{(1)}({h_{\textnormal{inv}}})+\left|\mathbb{P}\left\\{{h(X^{(1)})\neq
h_{\textnormal{inv}}(X^{(1)})}\right\\}-\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})\neq
h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})}\right\\}\right|\\\
+L\Psi_{\mathcal{G},\phi_{\textnormal{inv}}}.$

If this inequality holds, when $\phi_{\textnormal{inv}}$ is a conditionally
invariant feature mapping across source and target domains, we get
$\Psi_{\mathcal{G},\phi_{\textnormal{inv}}}=0$ and arrive at the original
Theorem 2.

Proof  Define
$\delta_{ij}=\mathbb{P}\left\\{{h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=i}\right\\}-\mathbb{P}\left\\{{h_{\textnormal{inv}}(X^{(1)})=j}\
\middle|\ {Y^{(1)}=i}\right\\}$, then
$\left|\delta_{ij}\right|\leq\Psi_{\mathcal{G},h_{\textnormal{inv}}}$. Because
there is no label shift, for any $i,j\in\left\\{1,2,\ldots,L\right\\}$ we have

$\displaystyle\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}-\mathbb{P}\left\\{{Y^{(1)}=i,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}=\delta_{ij}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i}\right\\}.$

As a result, we get

$\displaystyle\quad\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}-\mathbb{P}\left\\{{h(X^{(1)})=j,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}$
$\displaystyle=\sum_{i=1}^{L}\bigg{[}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}$
$\displaystyle\hskip 40.00006pt-\mathbb{P}\left\\{{h(X^{(1)})=j}\ \middle|\
{Y^{(1)}=i,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}\cdot\mathbb{P}\left\\{{Y^{(1)}=i,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}\bigg{]}$
$\displaystyle=\sum_{i=1}^{L}\Bigg{(}\bigg{[}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}$
$\displaystyle\hskip 50.00008pt-\mathbb{P}\left\\{{h(X^{(1)})=j}\ \middle|\
{Y^{(1)}=i,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}\bigg{]}\cdot\mathbb{P}\left\\{{Y^{(1)}=i,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}$
$\displaystyle\hskip
40.00006pt+\delta_{ij}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i}\right\\}\cdot\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}\Bigg{)}.$
(58)

Similarly we can calculate

$\displaystyle\quad\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j,Y^{{({\mathfrak{T}})}}=j}\right\\}-\mathbb{P}\left\\{{h(X^{(1)})=j,Y^{(1)})=j}\right\\}$
$\displaystyle=\sum_{i=1}^{L}\bigg{[}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=i}\right\\}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=i}\right\\}$
$\displaystyle\hskip 40.00006pt-\mathbb{P}\left\\{{h(X^{(1)})=j}\ \middle|\
{Y^{(1)}=j,h_{\textnormal{inv}}(X^{(1)})=i}\right\\}\cdot\mathbb{P}\left\\{{Y^{(1)}=j,h_{\textnormal{inv}}(X^{(1)})=i}\right\\}\bigg{]}$
$\displaystyle=\sum_{i=1}^{L}\Bigg{(}\bigg{[}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=i}\right\\}$
$\displaystyle\hskip 50.00008pt-\mathbb{P}\left\\{{h(X^{(1)})=j}\ \middle|\
{Y^{(1)}=j,h_{\textnormal{inv}}(X^{(1)})=i}\right\\}\bigg{]}\cdot\mathbb{P}\left\\{{Y^{(1)}=j,h_{\textnormal{inv}}(X^{(1)})=i}\right\\}$
$\displaystyle\hskip
40.00006pt+\delta_{ji}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j}\right\\}\cdot\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=i}\right\\}\Bigg{)}.$
(59)

Subtracting Eq. (C.1) from Eq. (C.1), we obtain

$\displaystyle\left[\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}-\mathbb{P}\left\\{{h(X^{(1)})=j,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}\right]$
$\displaystyle\hskip
80.00012pt-\left[\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j,Y^{{({\mathfrak{T}})}}=j}\right\\}-\mathbb{P}\left\\{{h(X^{(1)})=j,Y^{(1)}=j}\right\\}\right]$
$\displaystyle=$
$\displaystyle\bigg{[}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\ {Y^{{({\mathfrak{T}})}}\neq
j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}$
$\displaystyle\hskip 60.00009pt-\mathbb{P}\left\\{{h(X^{(1)})=j}\ \middle|\
{Y^{(1)}\neq j,h_{\textnormal{inv}}(X^{(1)})=j}\right\\}\bigg{]}\cdot B_{j}$
$\displaystyle-\bigg{[}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})\neq
j}\right\\}$ $\displaystyle\hskip 60.00009pt-\mathbb{P}\left\\{{h(X^{(1)})=j}\
\middle|\ {Y^{(1)}=j,h_{\textnormal{inv}}(X^{(1)})\neq
j}\right\\}\bigg{]}\cdot B^{\prime}_{j}$
$\displaystyle\quad+\sum_{i=1}^{L}\bigg{[}\delta_{ij}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i}\right\\}\cdot\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}$
$\displaystyle\hskip
60.00009pt-\delta_{ji}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=j}\right\\}\cdot\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=i}\right\\}\bigg{]}$
(60)

where we write

$\displaystyle B_{j}=\mathbb{P}\left\\{{Y^{(1)}\neq
j,h_{\textnormal{inv}}(X^{(1)})=j}\right\\},\ \
B^{\prime}_{j}=\mathbb{P}\left\\{{Y^{(1)}=j,h_{\textnormal{inv}}(X^{(1)})\neq
j}\right\\}.$

Now it follows that

$\displaystyle\quad\Bigg{|}\left[\mathbb{P}\left\\{{h(X^{(1)})\neq
h_{\textnormal{inv}}(X^{(1)})}\right\\}-\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})\neq
h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})}\right\\}\right]-\left(\mathcal{R}^{(1)}({h})-\mathcal{R}^{{({\mathfrak{T}})}}({h})\right)\Bigg{|}$
$\displaystyle=\Bigg{|}\left[\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})}\right\\}-\mathbb{P}\left\\{{h(X^{(1)})=h_{\textnormal{inv}}(X^{(1)})}\right\\}\right]$
$\displaystyle\hskip
144.54pt-\left[\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=Y^{{({\mathfrak{T}})}}}\right\\}-\mathbb{P}\left\\{{h(X^{(1)})=Y^{(1)}}\right\\}\right]\Bigg{|}$
$\displaystyle\overset{(i)}{\leq}\sum_{j=1}^{L}B_{j}+\sum_{j=1}^{L}B_{j}^{\prime}+\Bigg{|}\sum_{i=1}^{L}\sum_{j=1}^{L}\delta_{ij}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i}\right\\}$
$\displaystyle\hskip
100.00015pt\cdot\bigg{(}\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=j}\
\middle|\
{Y^{{({\mathfrak{T}})}}=i,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=j}\right\\}$
$\displaystyle\hskip
120.00018pt-\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})=i}\ \middle|\
{Y^{{({\mathfrak{T}})}}=j,h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})=i}\right\\}\bigg{)}\Bigg{|}$
$\displaystyle\overset{(ii)}{\leq}2\mathcal{R}^{(1)}({h_{\textnormal{inv}}})+\sum_{i=1}^{L}\sum_{j=1}^{L}\left|\delta_{ij}\right|\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i}\right\\}$
$\displaystyle\overset{(iii)}{\leq}2\mathcal{R}^{(1)}({h_{\textnormal{inv}}})+\sum_{i=1}^{L}\sum_{j=1}^{L}\Psi_{\mathcal{G},\phi_{\textnormal{inv}}}\cdot\mathbb{P}\left\\{{Y^{{({\mathfrak{T}})}}=i}\right\\}$
$\displaystyle=2\mathcal{R}^{(1)}({h_{\textnormal{inv}}})+L\Psi_{\mathcal{G},\phi_{\textnormal{inv}}},$

where in step (i) we sum over $j\in\mathcal{Y}=\\{1,2,\ldots,L\\}$ using Eq.
(C.1), and in step (ii) we use
$\sum_{j=1}^{L}B_{j}=\sum_{j=1}^{L}B^{\prime}_{j}=\mathcal{R}^{(1)}({h_{\textnormal{inv}}})$.
Step (iii) applies
$\left|\delta_{ij}\right|\leq\Psi_{\mathcal{G},\phi_{\textnormal{inv}}}$.
Rearranging terms, we arrive at

$\left|\mathcal{R}^{{({\mathfrak{T}})}}({h})-\mathcal{R}^{(1)}({h})\right|\leq
2\mathcal{R}^{(1)}({h_{\textnormal{inv}}})+\left|\mathbb{P}\left\\{{h(X^{(1)})\neq
h_{\textnormal{inv}}(X^{(1)})}\right\\}-\mathbb{P}\left\\{{h(X^{{({\mathfrak{T}})}})\neq
h_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})}\right\\}\right|\\\
+L\Psi_{\mathcal{G},\phi_{\textnormal{inv}}}.$

###  C.2 Proof of Theorem 3

For notation simplicity, for any $\phi\in\Phi$, we omit the bias term and
write $\phi(x)=Ax$. All proof below holds after adding the bias term. And we
use $u_{1},\ldots,u_{L}$ and $\widetilde{u}_{1},\ldots,\widetilde{u}_{L}$ to
denote $f^{(1)}(1),\ldots,f^{(1)}(L)$ and
$f^{{({\mathfrak{T}})}}(1),\ldots,f^{{({\mathfrak{T}})}}(L)$, respectively.

##### Part (a)

By Assumption 3, for any $A\in\Phi$, the distributions of $AX^{(1)}$ and
$AX^{{({\mathfrak{T}})}}$ are given by the mixture distributions of the form

$\mathcal{P}^{(1)}_{{AX}}=\sum_{j=1}^{L}\frac{1}{L}P_{A\epsilon}(Au_{j}),\ \
\mathcal{P}^{{({\mathfrak{T}})}}_{{AX}}=\sum_{j=1}^{L}\frac{1}{L}P_{A\epsilon}(A\widetilde{u}_{j}),$

where $P_{A\epsilon}(u)$ denotes the distribution of $A\epsilon$ with mean
shifted by a fixed vector $u$. Applying Lemma 2, the constraint in DIP,
$\mathcal{P}^{(1)}_{{AX}}=\mathcal{P}^{{({\mathfrak{T}})}}_{{AX}}$, implies
that we can find a permutation $\pi$ on $\\{1,2,\ldots,L\\}$ such that

$Au_{j}=A\widetilde{u}_{\pi(j)}\text{ for all $j\in\mathcal{Y}$}.$ (61)

Hence any representations learned by DIP should meet the constraint given in
Eq. (61) for some permutation $\pi$. Moreover, letting
$\dim\left(\operatorname{span}\\{u_{1}-\widetilde{u}_{\pi(1)},\ldots,u_{L}-\widetilde{u}_{\pi(L)}\\}\right)=k$,
Lemma 3 shows that any matrix $A$ satisfying Eq. (61) has the row space
orthogonal to
$\operatorname{span}\\{u_{1}-\widetilde{u}_{\pi(1)},\ldots,u_{L}-\widetilde{u}_{\pi(L)}\\}$
and $\operatorname{rank}(A)\leq p-k$.

Now consider all the pairs of $(A_{\pi},\pi)$ that meets (61) and the row
space of $A_{\pi}$ is the orthogonal complement to
$\operatorname{span}\\{u_{1}-\widetilde{u}_{\pi(1)},\ldots,u_{L}-\widetilde{u}_{\pi(L)}\\}$.
By definition of DIP in Eq. (11), we obtain
$\phi_{\textnormal{DIP}}=A_{\pi_{\star}}$, where

$\pi_{\star}=\operatorname*{arg\,min}_{\\{\pi\text{ is
permutation}\\}}\min_{g\in\mathcal{G},A_{\pi}}\mathcal{R}^{(1)}({g\circ
A_{\pi}}).$

We can clearly find examples where $\pi_{\star}\neq\mathbf{I}$, for instance,
as shown in Figure 1. This occurs as long as the label-flipping features offer
lower source risk than the conditionally invariant features while it cannot
generalize to the target domain. In this case, the conditional distributions
of $\phi_{\textnormal{DIP}}(X^{(1)})$ and
$\phi_{\textnormal{DIP}}(X^{{({\mathfrak{T}})}})$, given the labels, are
aligned after the target labels are permuted by $\pi_{\star}$.

##### Part (b)

Next, we consider JointDIP defined in Eq. (32). If $\phi_{\textnormal{inv}}$
is a linear mapping, say $\phi_{\textnormal{inv}}(x)=Bx$ for some matrix
$B\in{}^{r\times p}$, then writing $C=(A^{\top},B^{\top})^{\top}$ to denote
the concatenation of $A$ and $B$, we have

$\displaystyle\mathcal{P}^{(1)}_{{(\phi_{\textnormal{j-DIP}},\phi_{\textnormal{inv}})(X)}}$
$\displaystyle=\mathcal{P}^{(1)}_{{CX}}=\sum_{j=1}^{L}\frac{1}{L}P_{C\epsilon}(Cu_{j}),$
$\displaystyle\mathcal{P}^{{({\mathfrak{T}})}}_{{(\phi_{\textnormal{j-DIP}},\phi_{\textnormal{inv}})(X)}}$
$\displaystyle=\mathcal{P}^{{({\mathfrak{T}})}}_{{CX}}=\sum_{j=1}^{L}\frac{1}{L}P_{C\epsilon}(C\widetilde{u}_{j}).$

Under the JointDIP matching penalty, by Lemma 2 we similarly have

$Cu_{j}=C\widetilde{u}_{\pi(j)}\text{ for all $j$},$ (62)

for some permutation $\pi$. If $\pi\neq\mathbf{I}$, Eq. (62) suggests that
$Bu_{j}=B\widetilde{u}_{\pi(j)}=B\widetilde{u}_{j^{\prime}}$ for some
$j^{\prime}\neq j$, which contradicts the condition (35). Therefore, we must
have $\pi=\mathbf{I}$ (note that it is a valid solution because
$\phi_{\textnormal{inv}}$ is a CIC), which shows that
$\phi_{\textnormal{j-DIP}}$ is always conditionally invariant across source
and target domains.

Moreover, given that there is no label shift and $\phi_{\textnormal{j-DIP}}$
is conditionally invariant across $\mathcal{P}^{(1)}$ and
$\mathcal{P}^{{({\mathfrak{T}})}}$, we have
$\mathcal{R}^{{({\mathfrak{T}})}}({g\circ\phi_{\textnormal{j-DIP}}})=\mathcal{R}^{(1)}({g\circ\phi_{\textnormal{j-DIP}}})$
for all $g\in\mathcal{G}$. Then the optimization objective of population
JointDIP becomes equivalent to minimizing
$\mathcal{R}^{{({\mathfrak{T}})}}({g\circ\phi_{\textnormal{j-DIP}}})$. Since
$\phi_{\textnormal{inv}}^{0}=(\phi_{\textnormal{inv}},\mathbf{0}_{q-r})$ is a
feasible solution to the JointDIP constraint, we conclude that JointDIP is at
least as good as the optimal classifier built upon
$\phi_{\textnormal{inv}}^{0}$, i.e.,

$\displaystyle\mathcal{R}^{{({\mathfrak{T}})}}({h_{\textnormal{j-DIP}}})\leq\min_{g\in\mathcal{G}}\mathcal{R}^{{({\mathfrak{T}})}}({g\circ\phi_{\textnormal{inv}}^{0}}).$

##### Part (c)

In JointDIP, the joint distributions of $AX\in\mathbb{R}^{q}$ and
$\phi_{\textnormal{inv}}(X)\in\mathbb{R}^{r}$ are matched across the source
and target domains. Let

$T(X)=\left(\begin{array}[]{c}AX\\\
\phi_{\textnormal{inv}}(X)\end{array}\right)\in\mathbb{R}^{q+r}.$

The characteristic function of $T(X^{(1)})$ under the source distribution is
then:

$\displaystyle c^{(1)}(t)=c^{(1)}(t_{1},t_{2})$
$\displaystyle=\mathbb{E}\left[{e^{it^{\top}T(X^{(1)})}}\right]$
$\displaystyle=\frac{1}{L}\sum_{j=1}^{L}\mathbb{E}\left[{e^{it^{\top}T(X^{(1)})}}\
\middle|\ {Y^{(1)}=j}\right]$
$\displaystyle=\frac{1}{L}\sum_{j=1}^{L}\mathbb{E}\left[{e^{it_{1}^{\top}AX^{(1)}+it_{2}^{\top}\phi_{\textnormal{inv}}(X^{(1)})}}\
\middle|\ {Y^{(1)}=j}\right]$
$\displaystyle=\frac{1}{L}\sum_{j=1}^{L}e^{it_{1}^{\top}Au_{j}}\cdot\underbrace{\mathbb{E}\left[{e^{it_{1}^{\top}A\epsilon^{(1)}+it_{2}^{\top}\phi_{\textnormal{inv}}(X^{(1)})}}\
\middle|\ {Y^{(1)}=j}\right]}_{:=b_{j}(t_{1},t_{2})},$

where the second step uses the tower property, and the last step follows from
the general anticausal model, where conditional on $Y^{(1)}=j$,
$X^{(1)}=u_{j}+\epsilon^{(m)}$ under the source distribution. Similarly, we
can compute the characteristic function of $T$ under the target distribution:

$\displaystyle c^{{({\mathfrak{T}})}}(t)$
$\displaystyle=c^{{({\mathfrak{T}})}}(t_{1},t_{2})$
$\displaystyle=\frac{1}{L}\sum_{j=1}^{L}e^{it_{1}^{\top}A\widetilde{u}_{j}}\cdot\mathbb{E}\left[{e^{it_{1}^{\top}A\epsilon^{{({\mathfrak{T}})}}+it_{2}^{\top}\phi_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})}}\
\middle|\ {Y^{{({\mathfrak{T}})}}=j}\right]$
$\displaystyle=\frac{1}{L}\sum_{j=1}^{L}e^{it_{1}^{\top}A\widetilde{u}_{j}}\cdot
b_{j}(t_{1},t_{2}).$

Here the last equality holds because $\epsilon^{{({\mathfrak{T}})}}$ and
$\phi_{\textnormal{inv}}(X^{{({\mathfrak{T}})}})$ share the same conditional
distributions under the target distribution as $\epsilon^{(1)}$ and
$\phi_{\textnormal{inv}}(X^{(1)})$ under the source distribution. The JointDIP
matching penalty enforces $c^{(1)}(t)=c^{{({\mathfrak{T}})}}(t)$ for all $t$,
or

$\sum_{j=1}^{L}e^{it_{1}^{\top}Au_{j}}b_{j}(t_{1},t_{2})=\sum_{j=1}^{L}e^{it_{1}^{\top}A\widetilde{u}_{j}}b_{j}(t_{1},t_{2})\text{
for all $t_{1}\in\mathbb{R}^{q},t_{2}\in\mathbb{R}^{r}$}.$

If we take partial derivaties of both sides with respect to $t_{1}$,

$\displaystyle\quad\left(\sum_{j=1}^{L}Au_{j}e^{it_{1}^{\top}Au_{j}}b_{j}(t_{1},t_{2})\right)i+\sum_{j=1}^{L}e^{it_{1}^{\top}Au_{j}}\frac{\partial
b_{j}(t_{1},t_{2})}{\partial t_{1}}$
$\displaystyle=\left(\sum_{j=1}^{L}A\widetilde{u}_{j}e^{it_{1}^{\top}Au_{j}}b_{j}(t_{1},t_{2})\right)i+\sum_{j=1}^{L}e^{it_{1}^{\top}A\widetilde{u}_{j}}\frac{\partial
b_{j}(t_{1},t_{2})}{\partial t_{1}}.$

Setting $t_{1}=0$, it follows that

$\displaystyle\sum_{j=1}^{L}Au_{j}\cdot b_{j}(0,t_{2})$
$\displaystyle=\sum_{j=1}^{L}A\widetilde{u}_{j}\cdot b_{j}(0,t_{2}),$

or equivalently,

$\displaystyle\sum_{j=1}^{L}Au_{j}\cdot\mathbb{E}\left[{e^{it_{2}^{\top}\phi_{\textnormal{inv}}(X^{(1)})}}\
\middle|\ {Y^{(1)}=j}\right]$
$\displaystyle=\sum_{j=1}^{L}A\widetilde{u}_{j}\cdot\mathbb{E}\left[{e^{it_{2}^{\top}\phi_{\textnormal{inv}}(X^{(1)})}}\
\middle|\ {Y^{(1)}=j}\right].$ (63)

For any vector $a\in{}^{r}$, let $t_{2}=ka$ for $k\in\real$. In Eq. (63), take
derivative with respect to $k$ up to $L-1$ times and set $k=0$, then we get

$\displaystyle\begin{cases}\sum_{j=1}^{L}Au_{j}=\sum_{j=1}^{L}A\widetilde{u}_{j},\\\
\sum_{j=1}^{L}Au_{j}\cdot\mathbb{E}\left[{a^{\top}\phi_{\textnormal{inv}}(X^{(1)})}\
\middle|\
{Y^{(1)}=j}\right]=\sum_{j=1}^{L}A\widetilde{u}_{j}\cdot\mathbb{E}\left[{a^{\top}\phi_{\textnormal{inv}}(X^{(1)})}\
\middle|\ {Y^{(1)}=j}\right],\\\
\sum_{j=1}^{L}Au_{j}\cdot\mathbb{E}\left[{\left(a^{\top}\phi_{\textnormal{inv}}(X^{(1)})\right)^{2}}\
\middle|\
{Y^{(1)}=j}\right]=\sum_{j=1}^{L}A\widetilde{u}_{j}\cdot\mathbb{E}\left[{\left(a^{\top}\phi_{\textnormal{inv}}(X^{(1)})\right)^{2}}\
\middle|\ {Y^{(1)}=j}\right],\\\ \cdots\\\
\sum_{j=1}^{L}Au_{j}\cdot\mathbb{E}\left[{\left(a^{\top}\phi_{\textnormal{inv}}(X^{(1)})\right)^{L-1}}\
\middle|\
{Y^{(1)}=j}\right]=\sum_{j=1}^{L}A\widetilde{u}_{j}\cdot\mathbb{E}\left[{\left(a^{\top}\phi_{\textnormal{inv}}(X^{(1)})\right)^{L-1}}\
\middle|\ {Y^{(1)}=j}\right],\end{cases}$

where we use the fact that $\phi_{\textnormal{inv}}$ is a conditionally
invariant feature mapping. Writing the above linear system in a matrix form,
we have

$\displaystyle\underbrace{\begin{pmatrix}1&1&\cdots&1\\\
m_{1}^{1}(a)&m_{2}^{1}(a)&\cdots&m_{L}^{1}(a)\\\
\vdots&\vdots&\ddots&\vdots\\\
m_{1}^{L-1}(a)&m_{2}^{L-1}(a)&\cdots&m_{L}^{L-1}(a)\end{pmatrix}}_{=C_{\phi_{\textnormal{inv}}}(a)}\begin{pmatrix}({Au_{1}}-{A\widetilde{u}_{1}})^{\top}\\\
({Au_{2}}-{A\widetilde{u}_{2}})^{\top}\\\ \vdots\\\
({Au_{L}}-{A\widetilde{u}_{L}})^{\top}\end{pmatrix}=0,\ \
\forall\ell\in\\{1,2,\ldots,L\\},$

where $m_{j}^{l}(a)$ =
$\mathbb{E}\left[{\left(a^{\top}\phi_{\textnormal{inv}}(X^{(1)})\right)^{l}}\
\middle|\ {Y^{(1)}=j}\right]$. By assumption that
$C_{\phi_{\textnormal{inv}}}(a)$ is full rank, we get
$Au_{j}=A\widetilde{u}_{j}$ for all $j\in\\{1,\ldots,L\\}$, i.e.,

$\displaystyle\mathcal{P}^{(1)}_{{\phi_{\textnormal{j-DIP}}(X)}|{Y=y}}=\mathcal{P}^{{({\mathfrak{T}})}}_{{\phi_{\textnormal{j-DIP}}(X)}|{Y=y}},\
\ \forall y\in\\{1,2,\ldots,L\\},$

proving the result.

###### Lemma 2

Consider the mixture distributions $P_{1},P_{2}$ given by

$P_{1}=\sum_{j=1}^{L}w_{j}P(u_{j})\text{ and
}P_{2}=\sum_{j=1}^{L}\widetilde{w}_{j}P(\widetilde{u}_{j}),$

where $w_{j},\widetilde{w}_{j}$ are the mixing weights such that
$\sum_{j}w_{j}=\sum_{j}\widetilde{w}_{j}=1$ and $w_{j},\widetilde{w}_{j}\geq
0$, and $u_{j},\widetilde{u}_{j}\in\mathbb{R}^{p}$ are the fixed vectors for
the $j$-th components, and $P(u_{j}),P(\widetilde{u}_{j})$ denote the
distribution $P$ location shifted by $u_{j},\widetilde{u}_{j}$ respectively.
If $P_{1}=P_{2}$, then there exists a permutation $\pi$ such that
$u_{j}=\widetilde{u}_{\pi j}$ for all $j$.

Proof Let $\nu=\sum_{j=1}^{L}w_{j}\delta_{u_{j}}$ and
$\widetilde{\nu}=\sum_{j=1}^{L}\widetilde{w}_{j}\delta_{\widetilde{u}_{j}}$ be
the mixing distributions of $P_{1}$ and $P_{2}$ respectively. Let $X\sim
P_{1}$ and $\widetilde{X}\sim P_{2}$. Then we can write $X=U+Z$ and
$\widetilde{X}=\widetilde{U}+\widetilde{Z}$ where

$\begin{cases}U\sim\nu,Z\sim P\textnormal{ independent of $U$};\\\
\widetilde{U}\sim\widetilde{\nu},\widetilde{Z}\sim P\textnormal{ independent
of $\widetilde{U}$}.\end{cases}$

The characteristic functions of $X$ and $\widetilde{X}$ can be expressed as

$\varphi_{X}(t)=\varphi_{U}(t)\varphi_{Z}(t),\ \
\varphi_{\widetilde{X}}(t)=\varphi_{\widetilde{U}}(t)\varphi_{\widetilde{Z}}(t).$

Since $\varphi_{X}(t)=\varphi_{\widetilde{X}}(t)$ and
$\varphi_{Z}(t)=\varphi_{\widetilde{Z}}(t)$, it follows that
$\varphi_{U}(t)=\varphi_{\widetilde{U}}(t)$, or equivalently
$\nu=\widetilde{\nu}$.

Next let $S=\\{u_{1},\ldots,u_{L}\\}$ and
$\widetilde{S}=\\{\widetilde{u}_{1},\ldots,\widetilde{u}_{L}\\}$ to denote the
support sets of $\nu$ and $\widetilde{\nu}$. We can easily extend Lemma 2 in
Wu and Yang (2020) to the high-dimensional setting to show that

$d_{H}(S,\widetilde{S})\leq L\cdot W_{1}(\nu,\widetilde{\nu})=0,$

where $d_{H}(S,\widetilde{S})$ is the Hausdorff distance between $S$ and
$\widetilde{S}$, defined as

$d_{H}(S,\widetilde{S})=\max\left\\{\sup_{x\in
S}\inf_{\widetilde{x}\in\widetilde{S}}\lVert{x-\widetilde{x}}\rVert_{1},\sup_{\widetilde{x}\in\widetilde{S}}\inf_{x\in
S}\lVert{x-\widetilde{x}}\rVert_{1}\right\\},$

and where $W_{1}$ denotes the $1$-Wasserstein distance. By definition of
$d_{H}(S,\widetilde{S})$, this proves that we can find a permutation $\pi$
such that $u_{j}=\widetilde{u}_{\pi j}$ for all $j$, thus proving the lemma.

###### Lemma 3

Suppose that $\\{u_{1},\ldots,u_{L}\\}$ and
$\\{\widetilde{u}_{1},\ldots,\widetilde{u}_{L}\\}\subseteq\mathbb{R}^{p}$ are
collections of vectors such that

$\dim\left(\operatorname{span}\\{u_{1}-\widetilde{u}_{1},\ldots,u_{L}-\widetilde{u}_{L}\\}\right)=k<p.$

Then for any $A\in\mathbb{R}^{q\times p}$ satisfying
$Au_{j}=A\widetilde{u}_{j}$ for all $j=1,\ldots,L$, we have
$\operatorname{rank}(A)\leq p-k$. The equality is achieved if and only if the
rows of $A$ span
$\operatorname{span}\\{u_{1}-\widetilde{u}_{1},\ldots,u_{L}-\widetilde{u}_{L}\\}^{\perp}$.

Proof This lemma is a simple consequence of linear algebra. By the condition
of the lemma, $A(u_{j}-\widetilde{u}_{j})=0$ for all $j=1,\ldots,L$. This
implies that the row space of $A$ is orthogonal to the subspace spanned by
$\\{u_{1}-\widetilde{u}_{1},\ldots,u_{L}-\widetilde{u}_{L}\\}$, i.e.,

$\textnormal{row}({A})\perp\operatorname{span}\\{u_{1}-\widetilde{u}_{1},\ldots,u_{L}-\widetilde{u}_{L}\\}.$

In particular, $\dim(\textnormal{row}({A}))\leq p-k$, with equality achieved
if and only if $\textnormal{row}({A})$ is the orthogonal complement of
$\operatorname{span}\\{u_{1}-\widetilde{u}_{1},\ldots,u_{L}-\widetilde{u}_{L}\\}$.
By the rank theorem, $\operatorname{rank}(A)=\dim(\textnormal{row}({A}))\leq
p-k$, proving the lemma.

##  Appendix D Experiment details

This section provides details on various DA algorithms, training procedures,
network architectures, and hyperparameter selection strategies that we
employed in our numerical experiments (Section 5).

###  D.1 DA algorithms in numerical experiments

Table 8 provides a summary of the algorithms that we used in our numerical
experiments. Both mean distance and maximum mean discrepancy (MMD) were used
in DIP and CIP, but for JointDIP, we only considered MMD distance, because
joint matching with mean distance is essentially the same as the original DIP.
In addition to the DA algorithms introduced in this work, we also include IRM
(Arjovsky et al., 2019), V-REx (Krueger et al., 2021), and groupDRO (Sagawa et
al., 2019) in our experiments.

DA Algorithm |  $\widehat{w}^{(m)}$ Estimation | $\phi_{\textnormal{inv}}$ | Minimization Objective | $\mathfrak{D}\left({\cdot},{\cdot}\right)$
---|---|---|---|---
Tar | - | - | $\widehat{\mathcal{R}}^{{({\mathfrak{T}})}}({g\circ\phi})$ | -
ERM | - | - | $\widehat{\mathcal{R}}^{(1)}({g\circ\phi})$ | -
ERM-Pool | - | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\widehat{\mathcal{R}}^{(m)}({g\circ\phi})$ | -
DIP | - | - | $\widehat{\mathcal{R}}^{{({\mathfrak{T}})}}({g\circ\phi})+\lambda_{\textnormal{DIP}}\cdot\mathfrak{D}\left({\widehat{\mathcal{P}}^{(1)}_{{\phi(X)}}},{\widehat{\mathcal{P}}^{{({\mathfrak{T}})}}_{{\phi(X)}}}\right)$ | mean, MMD
DIP-Pool | - | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\left\\{\widehat{\mathcal{R}}^{(m)}({g\circ\phi})+\lambda_{\textnormal{DIP}}\cdot\mathfrak{D}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi(X)}}},{\widehat{\mathcal{P}}^{{({\mathfrak{T}})}}_{{\phi(X)}}}\right)\right\\}$ | mean, MMD
CIP | - | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\widehat{\mathcal{R}}^{(m)}({g\circ\phi})+\frac{\lambda_{\textnormal{CIP}}}{LM^{2}}\cdot\sum\limits_{y=1}^{L}\sum\limits_{m\neq m^{\prime}}\mathfrak{D}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi(X)}|{Y=y}}},{\widehat{\mathcal{P}}^{(m^{\prime})}_{{\phi(X)}|{Y=y}}}\right)$ | mean, MMD
IW-ERM | ERM-Pool | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\widehat{\mathcal{R}}^{(m)}({g\circ\phi};{\widehat{w}^{(m)}})$ | -
IW-CIP | CIP | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\widehat{\mathcal{R}}^{(m)}({g\circ\phi};{\widehat{w}^{(m)}})+\frac{\lambda_{\textnormal{CIP}}}{LM^{2}}\cdot\sum\limits_{y=1}^{L}\sum\limits_{m\neq m^{\prime}}\mathfrak{D}\left({\widehat{\mathcal{P}}^{(m)}_{{\phi(X)}|{Y=y}}},{\widehat{\mathcal{P}}^{(m^{\prime})}_{{\phi(X)}|{Y=y}}}\right)$ | mean, MMD
IW-DIP | CIP | - | $\widehat{\mathcal{R}}^{(1)}({g\circ\phi};{\widehat{w}^{(1)}})+\lambda_{\textnormal{DIP}}\cdot\mathfrak{D}\left({\widehat{\mathcal{P}}^{(1)}_{{\widehat{w}^{(1)}(Y)\phi(X)}}},{\widehat{\mathcal{P}}^{{({\mathfrak{T}})}}_{{\phi(X)}}}\right)$ | mean, MMD
JointDIP | - | $\widehat{\phi}_{\textnormal{CIP}}$ | $\widehat{\mathcal{R}}^{(1)}({g\circ\phi})+\lambda_{\textnormal{j-DIP}}\cdot\mathfrak{D}\left({\widehat{\mathcal{P}}^{(1)}_{{(\phi,\phi_{\textnormal{inv}})(X)}}},{\widehat{\mathcal{P}}^{{({\mathfrak{T}})}}_{{(\phi,\phi_{\textnormal{inv}})(X)}}}\right)$ | MMD
JointDIP-Pool | - | $\widehat{\phi}_{\textnormal{CIP}}$ | $\frac{1}{M}\sum\limits_{m=1}^{M}\left\\{\widehat{\mathcal{R}}^{(m)}({g\circ\phi})+\lambda_{\textnormal{j-DIP}}\cdot\mathfrak{D}\left({\widehat{\mathcal{P}}^{(m)}_{{(\phi,\phi_{\textnormal{inv}})(X)}}},{\widehat{\mathcal{P}}^{{({\mathfrak{T}})}}_{{(\phi,\phi_{\textnormal{inv}})(X)}}}\right)\right\\}$ | MMD
IW-JointDIP | CIP | $\widehat{\phi}_{\textnormal{CIP}}$ | $\widehat{\mathcal{R}}^{(1)}({g\circ\phi};{\widehat{w}^{(1)}})+\lambda_{\textnormal{j-DIP}}\cdot\mathfrak{D}\left({\widehat{\mathcal{P}}^{(1)}_{{\widehat{w}^{(1)}(Y)\cdot(\phi,\phi_{\textnormal{inv}})(X)}}},{\widehat{\mathcal{P}}^{{({\mathfrak{T}})}}_{{(\phi,\phi_{\textnormal{inv}})(X)}}}\right)$ | MMD
IRM | - | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\left\\{\widehat{\mathcal{R}}^{(m)}({g\circ\phi})+\lambda_{{\rm IRM}}\cdot\left\|\nabla_{w\mid w=1.0}\widehat{\mathcal{R}}^{(m)}({w\cdot g\circ\phi})\right\|_{2}^{2}\right\\}$ | -
V-REx | - | - | $\frac{1}{M}\sum\limits_{m=1}^{M}\widehat{\mathcal{R}}^{(m)}({g\circ\phi})+\lambda_{{\rm VREx}}\cdot\text{Var}\left(\widehat{\mathcal{R}}^{(m)}({g\circ\phi})\right)$ | -
groupDRO | - | - | $\sup_{\sum_{m=1}^{M}q_{m}=1,q_{m}\geq 0}\sum_{m=1}^{M}q_{m}\widehat{\mathcal{R}}^{(m)}({g\circ\phi})$ | -

Table 8: DA algorithms used in numerical experiments.

###  D.2 Training details and network architectures

Except for Camelyon17, the datasets we used in Section 5 are either generated
or perturbed synthetically. To avoid any potential bias that can result from
repeatedly running DA algorithms on a single, fixed synthetic dataset, for
each classification task, we use 10 different seeds to construct 10 different
datasets, and each DA algorithm is run once per seed.

Table 9 shows details on the training parameters and model architectures used
in all the experiments in Section 5. For multi-stage DA algorithms such as
JointDIP (where CIP is run before solving JointDIP optimization problem), the
number of epochs is the same in each training stage. As for the feature layer
used for CIP and DIP matching penalty, we utilize the last layer of the neural
network in experiments on SCM, MNIST, and CelebA; for experiments on
Camelyon17, we use the flattened output from CNN in DenseNet-121 following the
practice in Koh et al. (2021).

Classification Task | Epochs | Optimizer | Batch Size | NN Architecture | Feature Layer |  $\mathfrak{D}\left({\cdot},{\cdot}\right)$141414This column only applies to CIP and DIP. We always use MMD distance when implementing JointDIP.
---|---|---|---|---|---|---
SCM | 50 | Adam(lr=1e-2) | 100 | Linear Model | Last layer | mean
SCM_binary151515Ths setting is used for the additional experiments in Appendix E. | 50 | Adam(lr=1e-2) | 100 |  Linear$(10,10)$ or Linear$(10,10)$-ReLU Linear$(10,2)$ | Last layer | mean
MNIST | 20 | Adam(lr=1e-3) | 256 |  Conv2d$(1,20,5,1)-$ReLU MaxPool2d$(2)$ Conv2d$(20,50,5,1)-$ReLU MaxPool2d$(2)$ Linear$(800,500)-$ReLU Linear$(500,2)$ | Last layer | MMD
CelebA | 10 | Adam(lr=1e-3) | 64 |  Conv2d$(3,16,5,1)-$ReLU MaxPool2d$(2)$ Conv2d$(16,32,5,1,1)-$ReLU MaxPool2d$(2)$ Conv2d$(32,64,5,1)-$ReLU MaxPool2d$(2)$ Linear$(960,256)-$ReLU Linear$(256,2)$ | Last layer | MMD
Camelyon17 | 5 |  SGD(lr=1e-3, weight_decay=0.01, momentum=0.9) | 32 |  DenseNet-121 $($pretrained$=$False$)$ |  Flattened CNN output | MMD

Table 9: Training details of each classification task.

###  D.3 Hyperparameter selection

In our experiments on SCM, MNIST, and CelebA, instead of searching for
different hyperparameters for each seed, which would be time-consuming, we
choose the same hyperparameters across all seeds. In Appendix D.4, we compare
between two hyperparameter selection strategies, i.e., using different
hyperparameters versus using the same hyperparameters across seeds on SCMs,
and find that these two strategies yield similar results.

In order to select hyperparameters which allow for the highest average
accuracy, we use a simple grid search. Three different ways of hyperparameter
selection are presented in Gulrajani and Lopez-Paz (2020): training-domain
validation set, leave-one-out cross validation, and test-domain validation set
(oracle). The training-domain validation set approach is not suitable for us
as it assumes the target domain is very similar to the training domains, which
does not apply to our case. The leave-one-out cross validation approach is
cumbersome and would also suffer if the validation domain is significantly
different from the target domain. The test-domain validation set approach has
often been adopted in several domain adaptation papers (Ganin et al., 2016;
Arjovsky et al., 2019; Krueger et al., 2021) when an appropriate validation
set is unavailable. However, this approach contradicts the domain adaptation
principle where most target labels should be unavailable. Therefore, we
restrict access to just 10% of labeled samples from the target domain in each
seed, and check the accuracy of these labeled samples only at the final
checkpoint during hyperparameter selection.

Note that Camelyon17 provides a validation set. Hence, we adhere to the
standard submission guidelines from WILDS Koh et al. (2021); Sagawa et al.
(2021) for this dataset, selecting hyperparameters based on accuracy scores on
the validation set.

Algorithm | Hyperparameter | Grid Search Space
---|---|---
DIP (mean / MMD, Pool) | $\lambda_{\textnormal{DIP}}^{*}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
CIP (mean / MMD) | $\lambda_{{\rm CIP}}^{*}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
IW-DIP (mean / MMD, Pool) | $\lambda_{{\rm CIP}}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
$\lambda_{\textnormal{IW-DIP}}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
IW-CIP (mean / MMD) | $\lambda_{{\rm CIP}}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
$\lambda_{\textnormal{IW-CIP}}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
JointDIP (Pool) | $\lambda_{{\rm CIP}}^{*}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
$\lambda_{{\rm JointDIP}}^{*}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
IW-JointDIP | $\lambda_{{\rm CIP}}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
$\lambda_{{\rm JointDIP}}$ | $10^{\left\\{-2,-1,0,1,2\right\\}}$
IRM | $\lambda_{{\rm IRM}}$ | $10^{\left\\{-1,0,1,2,3,4\right\\}}$
iterations annealing | $0,10,100,1000,3000$
V-REx | $\lambda_{\rm{VREx}}$ | $10^{\left\\{-1,0,1,2,3,4\right\\}}$
iterations annealing | $0,10,100,1000,3000$
groupDRO | eta | $10^{\left\\{-2,-1,0,1\right\\}}$

Table 10: Hyperparameter search space. In experiments on camelyon17, the
search space of hyperparameters marked with ∗ is defined by
$10^{\left\\{-3,-2,-1,0,1,2\right\\}}$.

###  D.4 Different hyperparameter vs. same hyperparameter across all seeds

We carry out a comparison of two strategies for selecting hyperparameters on
SCMs: (1) using the same hyperparameter across all seeds, and (2) using
different hyperparameters for different seeds. We observe that in most cases,
the second approach tends to give a better performance. However, this
difference typically doesn’t affect the order of the top-performing
algorithms; the best performing algorithm tends to remain the same regardless
of the hyperparameter selection strategies. For the sake of simplicity and
time efficiency, we therefore adopt the first approach in our final results.

| SCM i@ | SCM ii@ | SCM iii@ | SCM iv@
---|---|---|---|---
Mean shift | Y | Y | Y | Y
CICs | N | Y | Y | Y
Label shift | N | Y | N | Y
Label-flipping features | N | N | Y | Y
Hyperparameter | Different | Same | Different | Same | Different | Same | Different | Same
DIP-mean | 88.2$\pm$1.9 | 87.6$\pm$1.5 | 73.3$\pm$10.7 | 62.0$\pm$2.9 | 36.3$\pm$11.8 | 34.5$\pm$14.9 | 35.6$\pm$14.0 | 35.3$\pm$14.6
DIP-MMD | 86.2$\pm$2.2 | 86.6$\pm$2.2 | 72.0$\pm$16.2 | 59.9$\pm$3.3 | 77.5$\pm$7.9 | 63.3$\pm$31.7 | 74.0$\pm$4.2 | 60.2$\pm$30.1
DIP-Pool-mean | 86.2$\pm$2.8 | 86.4$\pm$2.2 | 70.8$\pm$10.1 | 60.1$\pm$3.1 | 81.9$\pm$1.8 | 82.0$\pm$1.1 | 81.9$\pm$3.8 | 82.3$\pm$3.8
DIP-Pool-MMD | 86.2$\pm$1.9 | 85.5$\pm$3.1 | 74.2$\pm$13.5 | 61.3$\pm$9.3 | 82.5$\pm$1.4 | 83.3$\pm$0.8 | 82.2$\pm$4.1 | 82.5$\pm$3.8
CIP-mean | 56.6$\pm$11.8 | 55.9$\pm$12.0 | 79.3$\pm$6.4 | 75.7$\pm$6.5 | 82.2$\pm$1.5 | 81.8$\pm$1.3 | 83.2$\pm$1.8 | 82.1$\pm$1.2
CIP-MMD | 63.2$\pm$14.1 | 60.4$\pm$14.8 | 84.1$\pm$9.6 | 73.5$\pm$7.0 | 81.7$\pm$2.7 | 81.7$\pm$2.6 | 83.6$\pm$2.8 | 82.7$\pm$1.6
IW-CIP-mean | 54.1$\pm$11.8 | 54.0$\pm$11.8 | 90.8$\pm$1.4 | 90.4$\pm$0.8 | 80.3$\pm$4.2 | 81.2$\pm$4.1 | 83.1$\pm$2.3 | 83.8$\pm$2.2
IW-CIP-MMD | 64.4$\pm$12.7 | 56.8$\pm$13.6 | 89.5$\pm$4.3 | 90.4$\pm$0.8 | 80.9$\pm$3.5 | 80.5$\pm$6.3 | 84.1$\pm$1.9 | 83.5$\pm$2.8
IW-DIP-mean | 54.3$\pm$11.7 | 54.3$\pm$11.7 | 92.7$\pm$2.1 | 92.1$\pm$2.7 | 37.6$\pm$12.9 | 37.2$\pm$14.9 | 66.6$\pm$5.2 | 64.2$\pm$7.3
IW-DIP-MMD | 75.2$\pm$16.2 | 68.0$\pm$19.4 | 91.5$\pm$2.9 | 90.1$\pm$1.1 | 70.9$\pm$13.8 | 65.1$\pm$23.6 | 81.5$\pm$2.1 | 80.1$\pm$2.6
JointDIP | 85.8$\pm$2.2 | 86.8$\pm$1.9 | 77.7$\pm$12.5 | 70.6$\pm$6.2 | 85.3$\pm$2.8 | 85.4$\pm$2.1 | 81.6$\pm$3.9 | 82.8$\pm$1.9
IW-JointDIP | 72.7$\pm$17.3 | 68.0$\pm$19.4 | 91.3$\pm$2.8 | 90.0$\pm$1.3 | 82.9$\pm$3.3 | 82.9$\pm$8.1 | 85.4$\pm$1.9 | 85.1$\pm$3.7
IRM | 66.1$\pm$12.1 | 56.7$\pm$10.4 | 86.7$\pm$5.4 | 71.9$\pm$17.9 | 82.0$\pm$1.9 | 80.2$\pm$1.8 | 84.0$\pm$2.2 | 83.7$\pm$3.3
V-REx | 63.0$\pm$13.1 | 55.6$\pm$11.5 | 72.7$\pm$26.3 | 62.9$\pm$25.7 | 80.8$\pm$7.6 | 80.4$\pm$7.6 | 84.0$\pm$4.0 | 83.8$\pm$3.7
groupDRO | 54.4$\pm$10.0 | 54.4$\pm$10.0 | 65.2$\pm$26.6 | 64.2$\pm$25.7 | 82.5$\pm$6.7 | 81.3$\pm$7.4 | 84.7$\pm$3.8 | 84.3$\pm$3.2

Table 11: Target accuracy in linear SCMs under different hyperparameter tuning
strategies.

##  Appendix E Additional Experiments

In this section, we provide additional experiments on SCMs and MNIST to
investigate the effect of the number of source domains on the generalization
performance of CIP and its variants.

###  E.1 Effect of number of source domains in linear SCMs

To show that a sufficient number of source domains is necessary for finding
CICs, we conduct experiments on linear SCMs to compare the results of
utilizing fewer versus more source domains. Table 12 reports the risk
difference
$\widehat{\mathcal{R}}^{{({\mathfrak{T}})}}({h})-\widehat{\overline{\mathcal{R}}}({h};{\widehat{w}})$
between the source and target domains for four different CIP-based methods. We
observe that the risk difference tends to be closer to 0 for larger values of
$M$ in SCM ii@, iii@, and iv@, whereas smaller values of $M$ lead to a larger
risk difference. This finding suggests that CIP and IW-CIP may underperform
when the number of source domains $M$ is small, because they can only find
feature representations that are conditionally invariant across the source
domains, but with limited domains, it becomes challenging to obtain feature
representations that are conditionally invariant across both source and target
domains. Note that in SCM i@, the number of domains does not make a
significant difference, because CICs do not exist in SCM<EMAIL_ADDRESS>
| SCM i@ | SCM ii@ ($\times 10^{-2}$) | SCM iii@ ($\times 10^{-2}$) | SCM iv@ ($\times 10^{-2}$)
---|---|---|---|---
DA Algorithm | $M=2$ | $M=3$ | $M=2$ | $M=11$ | $M=2$ | $M=11$ | $M=2$ | $M=11$
CIP-mean | $1.7\pm 2.0$ | $2.1\pm 1.4$ | $11.3\pm 57.9$ | $-0.1\pm 0.9$ | $22.2\pm 26.3$ | $1.3\pm 3.6$ | $14.2\pm 25.3$ | $0.1\pm 0.9$
CIP-MMD | $1.7\pm 1.9$ | $1.3\pm 1.1$ | $-5.9\pm 7.6$ | $3.6\pm 9.6$ | $6.2\pm 7.1$ | $2.1\pm 5.2$ | $2.0\pm 4.7$ | $0.3\pm 2.1$
IW-CIP-mean | $1.9\pm 2.2$ | $2.0\pm 0.8$ | $89.0\pm 120.0$ | $-0.6\pm 10.9$ | $41.1\pm 50.1$ | $2.6\pm 6.9$ | $25.6\pm 29.9$ | $5.1\pm 6.7$
IW-CIP-MMD | $2.0\pm 2.3$ | $1.7\pm 1.0$ | $41.7\pm 60.0$ | $-4.7\pm 9.0$ | $30.2\pm 31.6$ | $6.6\pm 16.8$ | $11.4\pm 20.8$ | $5.3\pm 9.6$

Table 12: Risk difference
$\widehat{\mathcal{R}}^{{({\mathfrak{T}})}}({h})-\widehat{\overline{\mathcal{R}}}({h};{\widehat{w}})$
between source and target in linear SCMs. These differences are evaluated
under varying number of source domains $M$.

The minimum number of source domains necessary for finding CICs depends on the
types of interventions and the hypothesis class in use. Proposition 3 and
Lemma 1 show that the number of source domains needs to be greater than the
dimension of the perturbations, which is the case for linear SCMs above. In
this context, we also consider another SCM under binary interventions, as
described below. The purpose is to demonstrate that the required number of
source domains that enables successful finding of CICs via CIP can vary when
the hypothesis class changes.

  * •

SCM binary: The data generation model is

$\displaystyle Y^{(m)}$ $\displaystyle=\text{Bernoulli}(0.5),$ $\displaystyle
X^{(m)}_{[1:5]}$
$\displaystyle=0.2(Y^{(m)}-0.5)\cdot\mathds{1}_{5}+0.4\cdot\mathcal{N}(0,\mathbb{I}_{5}),$
$\displaystyle X^{(m)}_{[6:10]}$
$\displaystyle=0.2(Y^{(m)}-0.5)\cdot\mathds{1}_{5}+0.4\cdot\mathcal{N}(0,\mathbb{I}_{5})+A^{(m)},$

where $A^{(m)}$ is randomly sampled from $\\{a\in{}^{5}\mid
a_{i}\in\\{-1,1\\}\\}$ ($1\leq m\leq M$) without replacement in source
domains, and the target domain suffers larger intervention with
$A^{({\mathfrak{T}})}=(2,2,2,2,2)^{\top}$.

We use two different models in this experiment: a linear model and a two-layer
neural network. For the linear model, as shown in Proposition 3 and Lemma 1,
the number of source domains should be at least on the order of the dimension
of intervention plus one, that is, $M\geq 6$. We vary the number of source
domains $M$ from $2$ to $7$ and evaluate the performance of CIP-mean and CIP-
MMD in Table 13 (there is no label shift between source and target domains).
We find that when $M\geq 6$ the target accuracy and risk are close to those in
the source domains. We also observe that using more source domains than the
theoretically required minimum number of source domains can further improve
performance in practice.

For the two-layer neural network, we vary $M$ from $2^{1}$ to $2^{5}$. In
Table 14, we observe that when $M=32$, the target accuracy and target risk are
close to the source accuracy and source risk. This suggests that if the
hypothesis class is nonlinear and complex, CIP may need more domains to
accurately identify CICs. In this case, $M=2^{5}=32$ source domains covering
all possible binary interventions are required.

| CIP-mean | CIP-MMD
---|---|---
#source domains | src_acc | tar_acc | risk_diff | src_acc | tar_acc | risk_diff
$M=2$ | $85.1\pm 1.5$ | $50.6\pm 1.6$ | $1.8\pm 1.0$ | $84.5\pm 1.0$ | $51.7\pm 3.4$ | $2.4\pm 1.4$
$M=3$ | $84.7\pm 2.6$ | $54.3\pm 7.9$ | $9.1\pm 6.3\times 10^{-1}$ | $83.1\pm 1.7$ | $58.4\pm 9.9$ | $1.3\pm 1.2$
$M=4$ | $82.4\pm 2.3$ | $62.3\pm 11.2$ | $4.1\pm 5.1\times 10^{-1}$ | $81.9\pm 1.5$ | $67.9\pm 11.5$ | $5.2\pm 7.0\times 10^{-1}$
$M=5$ | $81.9\pm 2.0$ | $63.4\pm 12.2$ | $2.9\pm 4.3\times 10^{-1}$ | $81.7\pm 1.5$ | $69.8\pm 11.4$ | $3.8\pm 5.6\times 10^{-1}$
$M=6$ | $80.6\pm 1.5$ | $72.2\pm 10.8$ | $7.7\pm 13.2\times 10^{-2}$ | $81.3\pm 0.9$ | $73.0\pm 8.7$ | $2.5\pm 4.7\times 10^{-1}$
$M=7$ | $80.0\pm 0.5$ | $78.6\pm 1.6$ | $6.7\pm 9.9\times 10^{-3}$ | $81.1\pm 0.7$ | $73.9\pm 5.0$ | $1.1\pm 1.1\times 10^{-1}$

Table 13: Performance of CIP under a linear model on SCM-binary. The number of
source domains varies from $M=2$ to $M=7$.

| CIP-mean | CIP-MMD
---|---|---
#source domains | src_acc | tar_acc | risk_diff | src_acc | tar_acc | risk_diff
$M=2$ | $85.3\pm 1.7$ | $50.8\pm 1.7$ | $2.8\pm 1.8$ | $86.5\pm 1.9$ | $50.6\pm 1.6$ | $2.7\pm 1.1$
$M=4$ | $84.1\pm 1.3$ | $55.5\pm 8.7$ | $9.3\pm 7.0\times 10^{-1}$ | $85.2\pm 1.3$ | $53.5\pm 8.4$ | $2.1\pm 1.0$
$M=8$ | $82.5\pm 1.4$ | $63.8\pm 12.1$ | $3.1\pm 4.2\times 10^{-1}$ | $83.3\pm 1.4$ | $59.7\pm 11.3$ | $9.5\pm 7.2\times 10^{-1}$
$M=16$ | $80.8\pm 0.3$ | $76.4\pm 4.3$ | $2.7\pm 2.3\times 10^{-2}$ | $81.1\pm 0.5$ | $72.1\pm 7.4$ | $1.5\pm 1.8\times 10^{-1}$
$M=32$ | $80.5\pm 0.4$ | $78.4\pm 1.2$ | $1.6\pm 1.1\times 10^{-2}$ | $80.6\pm 0.3$ | $78.7\pm 1.2$ | $2.6\pm 1.6\times 10^{-2}$

Table 14: Performance of CIP under a two-layer neural network on SCM-binary.
The number of source domains varies from $M=2^{1}$ to $M=2^{5}$.

###  E.2 Effect of number of source domains in MNIST

Next we compare the performance of CIP-based methods on MNIST with varying
number of source domains ($M=3$ versus $M=5$). For $M=3$, we use the middle
three domains with rotation angles of $-15^{\circ},0^{\circ}$, and
$15^{\circ}$ from Section 5.2; for $M=5$, we use all source domains mentioned
in our experiments from Section 5.2. Table 15 shows the risk difference
between source and target domains for four different CIP-based methods.
Similar to our findings in the SCM experiments, the results demonstrate that
using a higher number of source domains helps in identifying CICs and reducing
the risk difference between source and target.

| MNIST i@ ($\times 10^{-2}$) | MNIST ii@ ($\times 10^{-4}$) | MNIST iii@ ($\times 10^{-4}$) | MNIST iv@ ($\times 10^{-4}$)
---|---|---|---|---
DA Algorithm | $M=3$ | $M=5$ | $M=3$ | $M=5$ | $M=3$ | $M=5$ | $M=3$ | $M=5$
CIP-mean | $43.4\pm 8.9$ | $16.0\pm 1.9$ | $31.5\pm 3.4$ | $17.1\pm 5.1$ | $45.0\pm 9.6$ | $26.8\pm 3.8$ | $44.4\pm 8.7$ | $33.7\pm 5.7$
CIP-MMD | $60.1\pm 4.1$ | $17.9\pm 4.9$ | $65.3\pm 13.6$ | $15.6\pm 3.2$ | $47.5\pm 6.8$ | $35.2\pm 4.1$ | $50.1\pm 10.6$ | $36.0\pm 6.3$
IW-CIP-mean | $58.2\pm 7.5$ | $15.7\pm 1.4$ | $33.3\pm 4.9$ | $11.5\pm 2.0$ | $40.3\pm 5.5$ | $27.5\pm 4.4$ | $43.5\pm 6.5$ | $27.6\pm 2.9$
IW-CIP-MMD | $62.9\pm 5.6$ | $20.3\pm 6.9$ | $62.5\pm 13.3$ | $14.1\pm 2.2$ | $42.7\pm 5.9$ | $36.3\pm 4.3$ | $46.5\pm 5.6$ | $35.3\pm 5.9$

Table 15: Risk difference
$\widehat{\mathcal{R}}^{{({\mathfrak{T}})}}({h})-\widehat{\overline{\mathcal{R}}}({h};{\widehat{w}})$
between source and target in MNIST experiments. These differences are
evaluated under varying number of source domains $M$.

## References

  * Amini and Gallinari (2003) M. Amini and P. Gallinari.  Semi-supervised learning with an explicit label-error model for misclassified data.  In _Proceedings of the 18th IJCAI_ , pages 555–560. Citeseer, 2003\.
  * Arjovsky et al. (2019) M. Arjovsky, L. Bottou, I. Gulrajani, and D. Lopez-Paz.  Invariant risk minimization.  _arXiv preprint arXiv:1907.02893_ , 2019.
  * Azizzadenesheli et al. (2019) K. Azizzadenesheli, A. Liu, F. Yang, and A. Anandkumar.  Regularized learning for domain adaptation under label shifts.  _arXiv preprint arXiv:1903.09734_ , 2019.
  * Baktashmotlagh et al. (2013) M. Baktashmotlagh, M. T. Harandi, B. C. Lovell, and M. Salzmann.  Unsupervised domain adaptation by domain invariant projection.  In _Proceedings of the IEEE International Conference on Computer Vision_ , pages 769–776, 2013.
  * Ben-David et al. (2006) S. Ben-David, J. Blitzer, K. Crammer, and F. Pereira.  Analysis of representations for domain adaptation.  _Advances in neural information processing systems_ , 19, 2006.
  * Ben-David et al. (2010) S. Ben-David, J. Blitzer, K. Crammer, A. Kulesza, F. Pereira, and J. W. Vaughan.  A theory of learning from different domains.  _Machine learning_ , 79(1-2):151–175, 2010.
  * Ben-Tal et al. (2009) A. Ben-Tal, L. El Ghaoui, and A. Nemirovski.  _Robust optimization_ , volume 28.  Princeton university press, 2009.
  * Ben-Tal et al. (2013) A. Ben-Tal, D. Den Hertog, A. De Waegenaere, B. Melenberg, and G. Rennen.  Robust solutions of optimization problems affected by uncertain probabilities.  _Management Science_ , 59(2):341–357, 2013.
  * Bühlmann (2020) P. Bühlmann.  Invariance, causality and robustness.  2020\.
  * Cadrin-Chênevert (2022) A. Cadrin-Chênevert.  Moving from imagenet to radimagenet for improved transfer learning and generalizability.  _Radiology: Artificial Intelligence_ , 4(5):e220126, 2022.
  * Chapelle et al. (2009) O. Chapelle, B. Scholkopf, and A. Zien.  Semi-supervised learning (chapelle, o. et al., eds.; 2006)[book reviews].  _IEEE Transactions on Neural Networks_ , 20(3):542–542, 2009.
  * Chen et al. (2022) L. Chen, M. Zaharia, and J. Zou.  Estimating and explaining model performance when both covariates and labels shift.  _arXiv preprint arXiv:2209.08436_ , 2022.
  * Chen and Bühlmann (2020) Y. Chen and P. Bühlmann.  Domain adaptation under structural causal models.  _arXiv preprint arXiv:2010.15764_ , 2020.
  * Chung and Lu (2006) F. Chung and L. Lu.  Concentration inequalities and martingale inequalities: a survey.  _Internet mathematics_ , 3(1):79–127, 2006.
  * Cortes and Mohri (2011) C. Cortes and M. Mohri.  Domain adaptation in regression.  In _Algorithmic Learning Theory: 22nd International Conference, ALT 2011, Espoo, Finland, October 5-7, 2011. Proceedings 22_ , pages 308–323. Springer, 2011.
  * Cortes and Mohri (2014) C. Cortes and M. Mohri.  Domain adaptation and sample bias correction theory and algorithm for regression.  _Theoretical Computer Science_ , 519:103–126, 2014.
  * Cortes et al. (2010) C. Cortes, Y. Mansour, and M. Mohri.  Learning bounds for importance weighting.  _Advances in neural information processing systems_ , 23, 2010.
  * Courty et al. (2017) N. Courty, R. Flamary, A. Habrard, and A. Rakotomamonjy.  Joint distribution optimal transportation for domain adaptation.  _Advances in Neural Information Processing Systems_ , 30, 2017.
  * Deng et al. (2009) J. Deng, W. Dong, R. Socher, L.-J. Li, K. Li, and L. Fei-Fei.  Imagenet: A large-scale hierarchical image database.  In _2009 IEEE conference on computer vision and pattern recognition_ , pages 248–255. Ieee, 2009.
  * Duchi et al. (2016) J. Duchi, P. Glynn, and H. Namkoong.  Statistics of robust optimization: A generalized empirical likelihood approach.  _arXiv preprint arXiv:1610.03425_ , 2016.
  * Finn et al. (2017) C. Finn, P. Abbeel, and S. Levine.  Model-agnostic meta-learning for fast adaptation of deep networks.  In _International conference on machine learning_ , pages 1126–1135. PMLR, 2017.
  * Ganin et al. (2016) Y. Ganin, E. Ustinova, H. Ajakan, P. Germain, H. Larochelle, F. Laviolette, M. Marchand, and V. Lempitsky.  Domain-adversarial training of neural networks.  _The Journal of Machine Learning Research_ , 17(1):2096–2030, 2016.
  * Garg et al. (2020) S. Garg, Y. Wu, S. Balakrishnan, and Z. C. Lipton.  A unified view of label shift estimation.  _arXiv preprint arXiv:2003.07554_ , 2020.
  * Gong et al. (2012) B. Gong, Y. Shi, F. Sha, and K. Grauman.  Geodesic flow kernel for unsupervised domain adaptation.  In _2012 IEEE conference on computer vision and pattern recognition_ , pages 2066–2073. IEEE, 2012.
  * Gong et al. (2016) M. Gong, K. Zhang, T. Liu, D. Tao, C. Glymour, and B. Schölkopf.  Domain adaptation with conditional transferable components.  In _International conference on machine learning_ , pages 2839–2848. PMLR, 2016.
  * Gopalan et al. (2011) R. Gopalan, R. Li, and R. Chellappa.  Domain adaptation for object recognition: An unsupervised approach.  In _2011 international conference on computer vision_ , pages 999–1006. IEEE, 2011.
  * Gretton et al. (2012) A. Gretton, K. M. Borgwardt, M. J. Rasch, B. Schölkopf, and A. Smola.  A kernel two-sample test.  _The Journal of Machine Learning Research_ , 13(1):723–773, 2012.