# Causal Falsification of Digital Twins

Rob Cornish∗ Department of Statistics, University of Oxford  Muhammad Faaiz
Taufiq∗ Department of Statistics, University of Oxford  Arnaud Doucet
Department of Statistics, University of Oxford  Chris Holmes  Department of
Statistics, University of Oxford  The Alan Turing Institute, UK

###### Abstract

Digital twins hold substantial promise in many applications, but rigorous
procedures for assessing their accuracy are essential for their widespread
deployment in safety-critical settings. By formulating this task within the
framework of causal inference, we show it is not possible to certify that a
twin is “correct” using real-world observational data unless potentially
tenuous assumptions are made about the data-generating process. To avoid these
assumptions, we propose an assessment strategy that instead aims to find cases
where the twin _is not_ correct, and present a general-purpose statistical
procedure for doing so that may be used across a wide variety of applications
and twin models. Our approach yields reliable and actionable information about
the twin under only the assumption of an i.i.d. dataset of real-world
observations, and in particular remains sound even in the presence of
arbitrary unmeasured confounding. We demonstrate the effectiveness of our
methodology via a large-scale case study involving sepsis modelling within the
Pulse Physiology Engine, which we assess using the MIMIC-III dataset of ICU
patients.

**footnotetext: Both authors contributed equally to this work.

##  1 Introduction

Digital twins provide the ability to interact in real-time with a virtual
clone of some real-world process of interest. This modelling approach promises
to radically transform planning and decision-making across a breadth of
fields, since it allows actions to be chosen with a full understanding of the
range of possible outcomes that may result. Early research on digital twins
has considered use-cases including aviation [6], manufacturing [76, 36, 40],
healthcare [70, 12, 10, 11], civil engineering [65], climate modelling [5],
and agriculture [28]. For a comprehensive overview of the surrounding
literature and outlook of this emerging technology, we refer the reader to
[69, 3, 30, 56].

Many applications of digital twins are considered safety-critical, which means
the cost of deploying an inaccurate twin to production is potentially very
high. As such, methodology for assessing the performance of a twin before its
deployment is essential for the safe, widespread adoption of digital twins in
practice [56]. In this work, we consider from first principles the problem of
assessing twin accuracy and propose a concrete, theoretically grounded, and
general-purpose methodology to this end. We focus specifically on the use of
statistical methods that leverage data obtained from the real-world process
the twin is designed to model. Such strategies are increasingly viable for
many applications as datasets grow larger, and offer the promise of lower
overheads compared with alternative strategies that rely for instance on
domain expertise.

A principal design objective for digital twins is to improve decision-making,
with much emphasis placed on their ability to reveal choices of actions that
are expected to produce some desirable outcome in the real world [69, 3, 30,
56]. We therefore adopt the position that, to be considered “accurate”, a twin
must capture the range of outcomes that _would_ occur when certain
interventions are applied to the real-world process of interest, and not
simply the outcomes that _are observed_ to occur as the process evolves on its
own. Accordingly, we formulate twin assessment within the framework of _causal
inference_ [63, 64, 57, 21], which provides a precise mathematical language
for posing interventional or “What if?”-type questions about the system, and
for answering these questions to the extent possible via statistical methods
applied to real-world data.

The primary reason for assessing the accuracy of a twin is to gauge its
reliability and robustness. It is therefore desirable for any assessment
procedure itself to be reliable and robust, and for any conclusions it draws
about the twin to be highly trustworthy. As such, our goal in this paper is to
obtain a methodology that is always _sound_ , even possibly at the expense of
being conservative: we prefer not to draw any conclusion about the accuracy of
the twin at all than to draw some conclusion that is potentially misleading.
To this end, we rely on minimal assumptions about the twin and the real-world
process of interest. In addition to improving robustness in any particular
situation, this also means our resulting methodology is highly generic, and so
may be applied to a wide variety of twins across many application domains.

### Contribution

Our contribution has the following components:

  * •

We show that if interventional information is sought from a twin, such as when
a twin is used for _planning_ , then assessment procedures that attempt to use
observational data to certify that the twin is accurate are unsound in general
unless strong and often tenuous assumptions are made about the data-generating
process.

  * •

We advocate an assessment paradigm based on _falsification_ that aims to find
cases in which the twin is _not_ accurate, for which it is possible to obtain
rigorous methodology under weak and widely applicable assumptions.

  * •

We propose a general-purpose statistical testing procedure for falsifying a
twin that obtains a soundness guarantee under only the assumption of an i.i.d.
dataset of observational trajectories. In particular, the procedure does not
require modelling the dynamics of the real-world process or any internal
implementation details of the twin. Key to our approach is a novel form of a
classical bound from the causal inference literature, which may also be of
interest in other contexts (see Section 5.2).

  * •

We demonstrate the effectiveness our procedure through a large-scale, real-
world case study in which we use the MIMIC-III dataset [29] of ICU patient
trajectories to assess the Pulse Physiology Engine [8], a high-fidelity, open-
source model for human physiology simulation.

### Related work

Various high-level guidelines and workflows have been proposed for the
assessment of digital twins in the literature to-date [18, 33, 55, 12, 34, 56,
14], as well as in the field of computational science more broadly [62, 56].
In some cases, these guidelines have been codified as standards: for example,
the ASME V&V40 Standard [1] provides a risk-based framework for assessing the
credibility of a model from a variety of factors that include source code
quality and the mathematical form of the model [17]. However, a significant
gap still exists between these guidelines and a practical implementation that
could be deployed for real twins, and the need for a rigorous lower-level
framework to enable the systematic assessment of twins has been noted in this
literature [12, 56, 32, 49]. We contribute towards this effort by describing a
precise statistical methodology for twin assessment that can be readily
implemented in practice, and which is accompanied by theoretical guarantees of
robustness that hold under minimal assumptions.

In addition, a variety of concrete assessment procedures have been applied to
individual digital twin models in the literature. For example, the Pulse
Physiology Engine [8], which we consider in our empirical case study, as well
as its predecessor BioGears [51] were both assessed by comparing their outputs
with ad hoc values based either on available medical literature or the
opinions of subject matter experts. Other twins have been assessed by
comparing their outputs with real-world data through a variety of bespoke
numerical schemes [38, 20, 50, 37, 28, 17]. However, to the best of our
knowledge, our paper is the first to analyse twin assessment systematically as
a causal problem, and to propose a general-purpose statistical procedure for
assessing twins that may be applied generically across many applications and
architectures.

### Outline

To better motivate our approach, Section 2 begins by illustrating potential
failure modes that can occur when digital twins are assessed using data
without properly accounting for causal considerations. Our discussion there is
intuitive and high-level. In Sections 3 and 4 we formalise these intuitions
via precise causal models for the twin and the data-generating process, which
we then use to develop our assessment methodology in Sections 5 and 6. Section
7 contains the results of our case study of this methodology.

Figure 1: The discrepancy between observational data and interventional
behaviour: the range of outcomes observed in the data for patients who were
administered the drug (blue) differs from what _would_ be observed if the drug
were administered to the general population (red).

##  2 Motivating example

Suppose that we are interested in using a digital twin to model the effect of
a drug on the weight of patients in a certain population, and we wish to
assess the accuracy of this twin using a dataset. Moreover, assume that this
drug interacts with an enzyme that is only present in part of the population.
Denote by $U\in\\{0,1\\}$ the presence or absence of the enzyme in a patient,
and assume that when $U=1$ the patient’s weight increases after action the
drug is administered, and that when $U=0$ the drug has no effect.
Additionally, suppose that, among the patients whose data we have obtained,
the drug was only prescribed to those for whom $U=1$, perhaps on the basis of
some initial lab reports available to the prescriber. Finally, suppose that
these lab results were not included in the dataset, so that the value of $U$
for each patient is not specified in the data we have available.

In this setup, since the drug was only administered to patients with $U=1$, it
would appear from the data that the drug causes patient weight to increase. A
naive approach to twin assessment would therefore consider a twin to be
accurate if it also replicated this behaviour and simulated an increase in
weight for patients to whom the drug is administered. However, when the drug
is administered to the general population, i.e. regardless of the value of
$U$, we would observe that the drug has no effect on patients for whom $U=0$.
As a result, a twin that always simulates an increase in weight after
administering the drug would not accurately predict the range of possible
outcomes that would occur when the drug is administered to a general member of
the population. Figure 1 illustrates this discrepancy under a toy model for
this scenario.

In the causal inference literature, any unmeasured quantity such as $U$ that
affects both some choice of action and the resulting outcome is referred to as
an _unmeasured confounder_. Beyond the example described here, whenever an
unmeasured confounder exists, a potential discrepancy arises between how the
real-world process was observed to behave in the dataset and how it _would_
behave under certain interventions. In practice, especially for complicated
real-world data, unmeasured confounding can rarely be ruled out as a
possibility altogether [72]. We therefore propose a methodology for assessing
twins using data that is designed to be robust to the possibility of
unmeasured confounders, which we present now.

##  3 Causal formulation

We begin by providing precise causal models for the twin and the real-world
process that the twin is designed to simulate. We do so in the language of
_potential outcomes_ [63, 64], although we note that we could have used the
alternative framework of directed acyclic graphs and structural causal models
[57] (see also [26] for a comparison of these two approaches). To streamline
our presentation, we eschew certain mathematical details in the main text,
especially relating to measure theory; however, the Appendix includes precise
statements of all definitions and results that we mention. Section A of the
Appendix also contains a summary of our notation.

###  3.1 The real-world process

We assume the real-world process operates over a fixed time horizon
$T\in\\{1,2,\ldots\\}$. This simplifies our presentation in what follows, and
it is straightforward to generalise our methodology to variable length time
horizons if needed. For each $t\in\\{0,\ldots,T\\}$, we assume the process
gives rise to an _observation_ at time $t$, which takes values in some set
$\mathcal{X}_{t}$. We also assume that the process can be influenced by some
_action_ taken at each time $t\in\\{1,\ldots,T\\}$. We denote the set of
actions available at time $t$ by $\mathcal{A}_{t}$, which in this work we
assume is always a finite set. For example, in a robotics context, the
outcomes may consist of the readings of all the sensors of the robot, and the
actions may consist of commands that can be input externally. In a medical
context, the outcomes may consist of the vital signs of a patient, and the
actions may consist of possible treatments or interventions. To streamline
notation, we will index these sets using vector notation, so that e.g.
$\mathcal{A}_{1:t}$ denotes the cartesian product
$\mathcal{A}_{1}\times\cdots\times\mathcal{A}_{t}$, and
$a_{1:t}\in\mathcal{A}_{1:t}$ is a choice of
$a_{1}\in\mathcal{A}_{1},\ldots,a_{t}\in\mathcal{A}_{t}$.

We model the dynamics of the real-world process via the longitudinal potential
outcomes framework proposed by Robins [59], which imposes only a weak temporal
structure on the underlying phenomena of interest and so may be applied across
a wide range of applications in practice. In particular, for each
$a_{1:T}\in\mathcal{A}_{1:T}$, we posit the existence of random variables

$X_{0},X_{1}(a_{1}),\ldots,X_{T}(a_{1:T}),$

where $X_{t}(a_{1:t})$ takes values in $\mathcal{X}_{t}$. We will denote this
sequence more concisely as $X_{0:T}(a_{1:T})$. Each $X_{t}(a_{1:t})$ is
referred to as a _potential outcome_. Intuitively, the value of $X_{0}$
represents initial data available before the first action is chosen, while the
remaining values of $X_{1:T}(a_{1:T})$ represent the real-world outcomes that
_would_ occur if actions $a_{1:T}$ were taken successively. The behaviour of
$X_{0:T}(a_{1:T})$ is therefore of fundamental interest to practitioners,
since it determines how actions should be chosen to achieve a sequence of
outcomes with optimal utility.

As random variables, each $X_{t}(a_{1:t})$ may depend on additional randomness
that is not explicitly modelled, and so in particular may be influenced by all
the previous potential outcomes $X_{0:t-1}(a_{1:t-1})$, and possibly other
random quantities. This models a process whose initial state is determined by
external factors, such as when a patient from some population first presents
at a hospital, and where the process then evolves according both to specific
actions chosen from $\mathcal{A}_{1:T}$ as well as additional exogenous
factors. It is clear that this structure applies to a wide range of phenomena
occurring in practice.

###  3.2 The digital twin

We think of the twin as a computational device that, when executed, outputs a
sequence of values intended to simulate a possible future trajectory of the
real-world process when certain actions in $\mathcal{A}_{1:T}$ are chosen,
conditional on some initial data in $\mathcal{X}_{0}$. We allow the twin to
make use of an internal random number generator to produce outputs that vary
stochastically even under fixed inputs. By executing the twin repeatedly, a
user may therefore estimate the range of behaviours that the real-world
process may exhibit under different action sequences, which can in turn be
used to guide planning and decision-making downstream.

To make this precise, we model the output the twin would produce at timestep
$t\in\\{1,\ldots,T\\}$ after receiving initialisation
$x_{0}\in\mathcal{X}_{0}$ and successive inputs $a_{1:t}\in\mathcal{A}_{1:t}$
as the quantity $h_{t}(x_{0},a_{1:t},U_{1:t})$, where $h_{t}$ is a function
taking values in $\mathcal{X}_{t}$, and each $U_{s}$ is some (possibly
multidimensional) random variable. To streamline notation, we will denote

$\widehat{X}_{t}(x_{0},a_{1:t})\coloneqq h_{t}(x_{0},a_{1:t},U_{1:t}),$

which will also refer to as a potential outcome. This means a full trajectory
output by the twin may be written as

$\widehat{X}_{1}(x_{0},a_{1}),\ldots,\widehat{X}_{T}(x_{0},a_{1:T}),$ (1)

or more compactly as $\widehat{X}_{1:T}(x_{0},a_{1:T})$. Conceptually,
$h_{1},\ldots,h_{T}$ constitute the program that executes inside the twin, and
$U_{1:T}$ may be thought of as the collection of all outputs of the internal
random number generator that the twin uses. We will assume

$U_{1:T}\perp\\!\\!\\!\perp(X_{0:T}(a_{1:T}):a_{1:T}\in\mathcal{A}_{1:T}),$

which amounts to saying that this random number generator does not interact
with the real-world process. We will also assume that repeated executions of
the twin give rise to i.i.d. copies of $U_{1:T}$. In turn, this means that,
when given fixed inputs $x_{0}$ and $a_{1:T}$, repeated executions of the twin
produce i.i.d. copies of (1). Otherwise, we make no assumptions about the
precise form or structure of either the $h_{t}$ or the $U_{t}$, which allows
our model to encompass a wide variety of twin implementations that may occur
in practice.

####  3.2.1 Correctness

The purpose of the twin is to predict the future of the real-world process
under different possible sequences of actions. As such, we will assume that
the twin has been designed to be correct in the following precise sense:

######  Definition 1.

The twin is _interventionally correct_ if, for all $x_{0}\in\mathcal{X}_{0}$
and $a_{1:T}\in\mathcal{A}_{1:T}$, the distribution of
$\widehat{X}_{1:T}(x_{0},{a}_{1:T})$ is equal111More generally we could allow
for inequality up to some specified tolerance level. However, our stricter
condition here is simpler, and suffices to convey the main ideas in what
follows. to the conditional distribution of $X_{1:T}(a_{1:T})$ given
$X_{0}=x_{0}$.

Operationally, if this holds, then by repeatedly executing the twin and
applying Monte Carlo techniques, it is possible to approximate arbitrarily
well the conditional distribution of the future of the real-world process
under each possible choice of action sequence. The same can also be shown to
hold when each action at each time $t$ is chosen dynamically on the basis of
previous observations in $\mathcal{X}_{0:t}$.222We refer to Chapter 5 of [72]
for more details about how dynamically chosen interventions can be modelled
within the potential outcomes framework. As a result, an interventionally
correct twin may be used for _planning_ , or in other words may be used to
select a policy for choosing actions that will yield a desirable distribution
over observations at each step.

We emphasise that interventional correctness does not mean the twin will
accurately predict what will happen in an almost sure sense on any specific
trajectory of the real-world process (unless $X_{0:T}(a_{1:T})$ is in fact
deterministic), but only the distribution of outcomes that will be observed
over repeated independent trajectories. However, this is sufficient for many
applications, and appears to be the strongest guarantee possible when dealing
with real-world phenomena whose underlying behaviour is stochastic.

Definition 1 introduces some technical difficulties that arise in the general
case when conditioning on events with probability zero (such as
$\\{X_{0}=x_{0}\\}$ if $X_{0}$ is continuous). In what follows, it is
therefore more convenient to consider an unconditional formulation of
interventional correctness. This is supplied by the following result, which
considers the behaviour of the twin when it is initialised with the (random)
value of $X_{0}$ taken from the real-world process, rather than with a fixed
choice of $x_{0}$. See Section B of the Appendix for a proof.

######  Proposition 1.

The twin is interventionally correct if and only if, for all choices of
$a_{1:T}\in\mathcal{A}_{1:T}$, the distribution of
$(X_{0},\widehat{X}_{1:T}(X_{0},a_{1:T}))$ is equal to the distribution of
$X_{0:T}(a_{1:T})$.

####  3.2.2 Online prediction

Our model here represents a twin at time $t=0$ making predictions about all
future timesteps $t\in\\{1,\ldots,T\\}$. In practice, many twins are designed
to receive new information at each timestep in an online fashion and update
their predictions for subsequent timesteps accordingly [18, 56]. Depending on
how the twin will be deployed, various notions of correctness can be devised
for this online setting. We describe two possibilities in Section C of the
Appendix, and show that these notions of correctness essentially reduce to
correctness in the sense of Definition 1. This motivates our focus on the
setting we have defined above in what follows.

##  4 Data-driven twin assessment

We now consider how to assess the accuracy of a twin. There are many
conceivable methods for doing this, including static analysis of the twin’s
source code and the solicitation of domain expertise, and in practice it seems
most robust to use a combination of different techniques rather than relying
on any single one [1, 56, 17]. However, in this paper, we focus on what we
term _data-driven assessment_ , which we see as an important and frequently-
applicable component of a larger assessment pipeline.333For example, data-
driven assessment falls naturally within the “validation” component of the
Verification, Validation, and Uncertainty Quantification (VVUQ) framework
[62]. That is, we consider the use of statistical methods that rely solely on
data obtained from the real-world process and simulated trajectories obtained
by running the twin. In particular, we do not assume knowledge of the dynamics
of the real-world process, or of the internal implementation details of the
twin.

We would like to obtain a methodology that is general-purpose and can be
applied to many types of twins and datasets with confidence in its validity.
To this end, we consider a deliberately weak model for how our real-world data
are obtained, which we describe below. We show that, without further
assumptions, it is not possible to obtain a data-driven assessment procedure
that can certify that a twin is interventionally correct. We instead propose a
strategy based on _falsifying_ the twin, which we develop into a concrete
statistical testing procedure in the next section.

###  4.1 Data model

We assume access to a dataset of trajectories obtained by observing the
interaction of some behavioural agents with the real-world process. We model
each trajectory as follows. First, we represent the action chosen at time
$t\in\\{1,\ldots,T\\}$ as an $\mathcal{A}_{t}$-valued random variable $A_{t}$.
We then obtain a trajectory in our dataset by recording at each step the
action $A_{t}$ chosen and the observation $X_{t}(A_{1:t})$ corresponding to
this choice of action. As a result, each observed trajectory has the following
form:

$X_{0},A_{1},X_{1}(A_{1}),\ldots,A_{T},X_{T}(A_{1:T}).$ (2)

This corresponds to the usual _stable unit treatment value assumption (SUTVA)_
[4] (also referred to as _consistency_ or _stability_) of the potential
outcomes framework. SUTVA is a standard and mild assumption used to relate the
interventional behaviour of the real-world process under potentially
hypothetical choices of actions $a_{1:t}\in\mathcal{A}_{1:t}$ to the data
actually observed under the actions $A_{1:t}$ actually taken. We then model
our full dataset as a set of i.i.d. copies of (2).

Similarly to the real-world process and the twin, by modelling the actions as
random variables, we implicitly allow their value to depend on additional
random quantities than we have defined here, including quantities that also
influence the values of potential outcomes $X_{t}(a_{1:t})$. This introduces
the possibility that the the data are confounded, which leads to significant
challenges for twin assessment discussed below.

###  4.2 The impossibility of certifying that a twin is correct

A natural high-level strategy for twin assessment, which appears at least in
principle to underlie many approaches used in practice, has the following
structure. First, some hypothesis $\mathcal{H}$ is chosen with the following
property:

If $\mathcal{H}$ is true, then the twin is interventionally correct. (3)

Then, data is used to try to show that $\mathcal{H}$ is true, perhaps up to
some level of confidence to account for the fact that the dataset is finite.
When this is successful, it follows by construction that the twin is
interventionally correct. Assessment procedures with this structure are
appealing because they holds the promise of a strong guarantee about the
accuracy of the twin. Unfortunately, the following fundamental result means it
is not possible to obtain a data-driven procedure in this way that is valid in
the general context we consider:

######  Proposition 2.

If $\mathbb{P}(A_{1:t}\neq a_{1:t})>0$, then the distribution of
$X_{t}(a_{1:t})$ cannot be uniquely determined from the distribution of the
data in (2) without further assumptions.

Proposition 2 is a foundational result in the causal inference literature. The
statement holds due to what has been referred to as the _fundamental problem
of causal inference_ [23], which is namely that the potential outcomes
$X_{t}(a_{1:t})$ are only observed in (2) when $A_{1:t}=a_{1:t}$. As such, the
behaviour of the data in (2) only identifies the behaviour of $X_{t}(a_{1:t})$
on the event $\\{A_{1:t}=a_{1:t}\\}$. On the event $\\{A_{1:t}\neq
a_{1:t}\\}$, $X_{t}(a_{1:t})$ may vary arbitrarily while still giving rise to
the same distribution over datasets. In practice, it is usually very
reasonable to assume that $\\{A_{1:t}\neq a_{1:t}\\}$ occurs with positive
probability, since this is only violated in the extreme case that the agent
always chooses the same action sequence $a_{1:t}$. We include a self-contained
proof of Proposition 2 in Section D of the Appendix for completeness.

Intuitively speaking, the distribution of the data encodes the information
that would be contained in an infinitely large dataset of trajectories.
Proposition 2 therefore imposes a fundamental limit on what can be learned
about the distribution of each $X_{t}(a_{1:t})$ from the data we have assumed.
In particular, it follows that if $\mathcal{H}$ is any hypothesis satisfying
(3), then $\mathcal{H}$ also cannot be determined to be true from the data.
This is because, if we could do so, then we could also determine the
distribution of each $X_{t}(a_{1:t})$, since by Proposition 1 this would be
equal to the distribution of $\widehat{X}_{t}(X_{0},a_{1:t})$. In other words,
we cannot use the data alone to certify that the twin is interventionally
correct.

###  4.3 The assumption of no unmeasured confounding

Proposition 2 is true in the general case, when no additional assumptions
about the data-generating process are made. One way forward is therefore to
introduce assumptions under which the distribution of $X_{0:T}(a_{1:T})$ _can_
be determined from the data. This would mean that it is possible to certify
that the twin is interventionally correct, since, at least in principle, we
could simply compare this distribution with the distribution of
$(X_{0},\widehat{X}_{1:T}(X_{0},a_{1:T}))$ obtained from the twin.

The most common such assumption in the causal inference literature is that the
data are free of _unmeasured confounding_. Informally, this holds when each
action $A_{t}$ is chosen by the behavioural agent solely on the basis of the
information available at time $t$ that is actually recorded in the dataset,
namely $X_{0},A_{1},X_{1}(A_{1}),\ldots,A_{t-1},X_{t-1}(A_{1:t-1})$, as well
as possibly some additional randomness that is independent of the real-world
process, such as the outcome of a coin toss.555This can be made precise via
the _sequential randomisation assumption_ (also known as _sequential
ignorability_) introduced by Robins [59]. Chapter 5 of [72] provides an
overview of this. Unobserved confounding is said to be present whenever this
does not hold, i.e. whenever some unmeasured factor simultaneously influences
both the agent’s choice of action and the observation produced by the real-
world process. For instance, in the motivating example given in Section 2, the
presence or absence of the enzyme encoded by $U$ is an unobserved confounder,
since it affects both the decision to prescribe the drug, as well as the
resulting weight of the patient.

It is reasonable to assume that the data are unconfounded in certain contexts.
For example, in certain situations it may be possible to gather data in a way
that specifically guarantees there is no confounding. Randomised controlled
trials, which ensure that each $A_{t}$ is chosen via a carefully designed
randomisation procedure [39, 53], constitute a widespread example of this
approach. Likewise, it is possible to show666See Section E of the Appendix for
a proof. that the data are unconfounded if each $X_{t}(a_{1:t})$ is a
deterministic function of $X_{0:t-1}(a_{1:t-1})$ and $a_{t}$, which may be
reasonable to assume for example in certain low-level physics or engineering
contexts. However, for stochastic phenomena and for typical datasets, it is
widely acknowledged that the assumption of no unmeasured confounding will
rarely hold, and so assessment procedures based on this assumption may yield
unreliable results in practice [52, 72].

###  4.4 General-purpose assessment via falsification

Our goal in this paper is to obtain an assessment methodology that is general-
purpose, and as such we would like to avoid introducing assumptions such as
unconfoundedness that do not hold in general. To achieve this, borrowing
philosophically from Popper [58], we propose a strategy that replaces the goal
of verifying the interventional correctness of the twin with that of
_falsifying_ it. Specifically, we consider hypotheses $\mathcal{H}$ with the
dual property to (3), namely:

If the twin is interventionally correct, then $\mathcal{H}$ is true. (4)

We will then try to show that each such $\mathcal{H}$ is _false_. Whenever we
are successful, we will thereby have gained some knowledge about a failure
mode of the twin, since by construction the twin can only be correct if
$\mathcal{H}$ is true. In effect, each $\mathcal{H}$ we falsify will
constitute a _reason_ that the twin is not correct, and may suggest concrete
improvements to its design, or may identify cases where its output should not
be trusted.

Importantly, unlike for (3), Proposition 2 does not preclude the possibility
of data-driven assessment procedures based on (4). In particular, as we show
in the next section, there do exist hypotheses $\mathcal{H}$ satisfying (4)
that can in principle be determined to be false from the data alone without
additional assumptions. In this sense, falsification provides a means for
_sound_ data-driven twin assessment, whose results can be relied upon across a
wide range of circumstances. On the other hand, falsification approaches
cannot provide a _complete_ guarantee about the accuracy of a twin: even if we
fail to falsify a large number of $\mathcal{H}$ satisfying (4), we cannot then
infer that the twin is correct. As such, in situations where (for example) it
is reasonable to believe that the data are in fact unconfounded, it is likely
desirable to use this assumption to obtain additional information about the
twin than is possible from falsification alone.

Overall, our approach follows the guidelines for empirical analysis proposed
by econometrician Charles Manski [43], who advocates for researchers always to
establish a “domain of consensus” [46] by analysing their data initially under
minimal assumptions, and only afterwards to consider what follows as
assumptions of various strengths are introduced. To this end, Manski and
others developed a literature on _partial identification_ of probability
distributions [41, 42, 2, 44, 47, 45, 27, 9] in order to provide precise
statistical tools with which this high-level strategy may be implemented.
Similar ideas and motivations have appeared in more recent work: notably, [75]
use partial identification techniques to obtain a more efficient algorithm for
the online learning of optimal dynamic treatment regimes [52], while [25]
propose a meta-algorithm that involves falsifying treatment effect estimators
that can be determined to be invalid using data from a randomised controlled
trial. In the next section, we describe how ideas from this literature give
rise to a family of $\mathcal{H}$ satisfying (4) whose falsification leads to
granular quantitative information about the twin that could be used
downstream. We then provide a statistical testing procedure that can be used
to falsify these hypotheses using a finite dataset. Our resulting methodology
requires no additional assumptions beyond those already described, and so may
constitute a valuable addition to twin assessment pipelines in safety-critical
settings where robustness is essential.

##  5 Hypotheses derived from causal bounds

###  5.1 Hypothesis parameters

We now describe a concrete family of hypotheses that satisfy (4), and which
can in principle be falsified using the data we have assumed. Each hypothesis
will depend on a specific choice of the following parameters:

  * •

A timestep $t\in\\{1,\ldots,T\\}$

  * •

A scalar-valued function $f:\mathcal{X}_{0:t}\to\mathbb{R}$

  * •

A sequence of actions $a_{1:t}\in\mathcal{A}_{1:t}$

  * •

A sequence of subsets $B_{0:t}\subseteq\mathcal{X}_{0:t}$.

To streamline notation, in this section, we will consider these parameters to
be fixed. However, we emphasise that our construction can be instantiated for
many different choices of these parameters, and indeed we will do so below. We
think of $f$ as expressing a specific _outcome_ or eventuality at time $t$ in
terms of the data we have assumed. Accordingly, for each
$a_{1:t}^{\prime}\in\mathcal{A}_{1:t}$, we define new potential outcomes

$Y(a_{1:t}^{\prime})\coloneqq f(X_{0:t}(a_{1:t}^{\prime})).$

For example, if $X_{0:t}(a_{1:t}^{\prime})$ represents a full patient history
at time $t$ after treatments $a_{1:t}^{\prime}$, then $Y(a_{1:t}^{\prime})$
might represent the patient’s heart rate at time $t$ after the same
treatments.

#### Parameter constraints

The hypotheses we consider are only well-defined if $a_{1:t}$ and $B_{0:t}$
satisfy the following positivity conditions:

$\displaystyle\mathbb{P}(X_{0:t}(a_{1:t})\in B_{0:t})$ $\displaystyle>0$ (5)
$\displaystyle\mathbb{P}(\widehat{X}_{1:t}(X_{0},a_{1:t})\in B_{1:t})$
$\displaystyle>0.$ (6)

We will also require that $f$ is such that, for some known values of
$y_{\mathrm{lo}},y_{\mathrm{up}}\in\mathbb{R}$, we have

$\mathbb{P}(y_{\mathrm{lo}}\leq Y(a_{1:t})\leq y_{\mathrm{up}}\mid
X_{0:t}(a_{1:t})\in B_{0:t})=1.$ (7)

In practice, it is easy to ensure that these conditions hold before proceeding
with any further analysis. (6) can be tested directly by running the twin,
while (5) is implied by the stronger condition

$\mathbb{P}(X_{0:t}(A_{1:t})\in B_{0:t},A_{1:t}=a_{1:t})>0,$ (8)

which can also be tested using data since both $X_{0:t}(A_{1:t})$ and
$A_{1:t}$ are observed quantities.777In fact, if
$\mathbb{P}(X_{0:t}(a_{1:t})\in B_{0:t})>0$ but we have
$\mathbb{P}(X_{0:t}(A_{1:t})\in B_{0:t},A_{1:t}=a_{1:t})=0$, then it is easy
to check that the bounds in Theorem 1 are necessarily uninformative. We
describe the practicalities of this in Section 6. Likewise, (7) immediately
follows if $f$ satisfies the stronger condition

$y_{\mathrm{lo}}\leq f(x_{0:t})\leq y_{\mathrm{up}}\qquad\text{for all
$x_{0:t}\in B_{0:t}$.}$ (9)

This holds automatically in certain cases, such as when dealing with a binary
outcome (e.g. to denote survival a patient), or otherwise can be enforced
simply by clipping the value of $f$ to live within
$[y_{\mathrm{lo}},y_{\mathrm{up}}]$. We describe a practical scheme for
choosing $f$ in this way in Section 7.

###  5.2 Causal bounds

Our hypotheses are based on a classical family of bounds from the causal
inference literature whose development began with Manski [42]. We derive a
novel variant of these bounds that is more suited to the longitudinal nature
of our problem. In the structural causal modelling framework [57], a related
bound was recently given by Zhang & Bareinboim as Corollary 1 in [75].
However, this result involves a complicated ratio of unknown quantities that
makes statistical estimation of the bounds difficult, since it is not obvious
how to obtain an unbiased estimator for the ratio term. In contrast, our
proposed causal bounds are considerably simpler, and easily give rise to an
unbiased estimator with attractive theoretical properties. To state this
result, we first define the following random variables:

$\displaystyle N$ $\displaystyle\coloneqq\max\\{0\leq s\leq t\mid
A_{1:s}=a_{1:s}\\}$ (10) $\displaystyle Y_{\mathrm{lo}}$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:t}=a_{1:t})\,Y(A_{1:t})+\mathbbm{1}(A_{1:t}\neq
a_{1:t})\,y_{\mathrm{lo}}$ $\displaystyle Y_{\mathrm{up}}$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:t}=a_{1:t})\,Y(A_{1:t})+\mathbbm{1}(A_{1:t}\neq
a_{1:t})\,y_{\mathrm{up}}.$

We will also denote:

$\displaystyle Q$ $\displaystyle\coloneqq\mathbb{E}[Y(a_{1:t})\mid
X_{0:t}(a_{1:t})\in B_{0:t}]$ (11) $\displaystyle Q_{\mathrm{lo}}$
$\displaystyle\coloneqq\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0:N}(A_{1:N})\in
B_{0:N}]$ $\displaystyle Q_{\mathrm{up}}$
$\displaystyle\coloneqq\mathbb{E}[Y_{\mathrm{up}}\mid X_{0:N}(A_{1:N})\in
B_{0:N}].$

The result is then as follows:

######  Theorem 1.

If (5) and (7) hold, then

$Q_{\mathrm{lo}}\leq Q\leq Q_{\mathrm{up}}.$ (12)

We emphasise that the bounds $Q_{\mathrm{lo}}$ and $Q_{\mathrm{up}}$ here are
expressed simply as (conditional) expectations. This makes their unbiased
estimation straightforward, which we make use of to obtain exact confidence
intervals for both terms in Section 6. A proof of Theorem 1 is given in
Section F of the Appendix.

#### Intuition

The value of $Q$ depends on the distribution of $X_{0:t}(a_{1:t})$, which by
Proposition 2 cannot be determined from the distribution of the data. On the
other hand, both $Q_{\mathrm{lo}}$ and $Q_{\mathrm{up}}$ can be determined,
since the relevant random variables $Y_{\mathrm{lo}}$, $Y_{\mathrm{up}}$, $N$,
and $X_{0:N}(A_{1:N})$ can all be expressed as functions of $X_{0:t}(A_{1:t})$
and $A_{1:t}$. In other words, Theorem 1 bounds the behaviour of an
unobtainable quantity in terms of quantities that we can in fact estimate. At
a high level, this is achieved by replacing the random variable $Y(a_{1:t})$,
whose value is only observed when $A_{1:t}=a_{1:t}$, with $Y_{\mathrm{lo}}$
and $Y_{\mathrm{up}}$, which are equal to $Y(a_{1:t})$ when its value is
observed (i.e. when $A_{1:t}=a_{1:t}$), and which fall back to the worst-case
values of $y_{\mathrm{lo}}$ and $y_{\mathrm{up}}$ otherwise. The full details
of the proof are given in the Appendix. We emphasise that Theorem does not
require adding any new assumptions to our causal model defined earlier, and in
particular remains true under arbitrary degrees of confounding. It therefore
seems reasonable to expect (12) will hold for a very wide variety of datasets
that occur in practice.

The effectiveness of the falsification procedure we describe below depends on
the tightness of the bounds in (12), which in turn is determined by the value
of $\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:N}(A_{1:N})\in B_{0:N})$. This value
is closely related to the classical _propensity score_ in the causal inference
literature [61], and quantifies how likely the action sequence $a_{1:t}$ is to
occur in the data. Intuitively, for very likely action sequences,
$Y(a_{1:t})=Y(A_{1:t})$ holds with very high probability. This means that the
effect of unmeasured confounding on the overall expectation is limited, and
implies that the bounds on $Q$ are tight. The bounds become looser for action
sequences with a low probability of being chosen, since in this case a greater
proportion of interventional outcomes $Y(a_{1:t})$ remain unobserved. In the
extreme case that $A_{1:t}=a_{1:t}$ has zero probability of occurring, (12)
relaxes to the worst-case bounds $[y_{\mathrm{lo}},y_{\mathrm{up}}]$ that are
trivially implied by (7).

#### Optimality of Theorem 1

In Section F of the Appendix we provide several additional results
demonstrating that Theorem 1 cannot be improved upon without introducing
additional assumptions. Specifically, we show that the bounds in Theorem 1 are
sharp when each outcome $Y(a_{1:t}^{\prime})$ is a function of the final
covariates $X_{t}(a_{1:t}^{\prime})$ only (which holds for example throughout
our case study), provided $y_{\mathrm{lo}}$ and $y_{\mathrm{up}}$ are chosen
sensibly. We also consider a more general setting in which each potential
outcome $Y(a_{1:t}^{\prime})$ is not constrained to be a function of
$X_{0:t}(a_{1:t}^{\prime})$ (i.e. it may have an arbitrary joint distribution
with the other random variables in the model), a situation that does not
conform to our model but that may be of interest in other contexts. We show
that Theorem 1 continues to hold in this case, and is in fact always optimal
without further assumptions.

In Section F.4 of the Appendix, we also consider the task of obtaining bounds
on $\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})]$, i.e. where we condition on
the value of $X_{0:t}(a_{1:t})$ rather than on the event
$\\{X_{0:t}(a_{1:t})\in B_{0:t}\\}$. These can be obtained directly from
Theorem 1 if $X_{0:t}(a_{1:t})$ is discrete, or by a simple modification of
the proof of this result if $X_{1:t}(a_{1:t})$ is discrete (but $X_{0}$ is
possibly continuous), which can be used to recover Manski’s original bounds
from [42] as a special case. However, somewhat surprisingly, we show that if
$X_{1:t}(a_{1:t})$ has a continuous component, then it is not possible to
obtain nontrivial bounds on $\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})]$ in
terms of observational quantities without introducing additional assumptions,
and so Theorem 1 is in some sense the best that may be hoped for in this
situation. We show below that this result is nevertheless powerful enough to
obtain useful information about the twin under only the minimal assumptions
that we have made.

###  5.3 Derived hypotheses

The hypotheses we consider are based on the value of the outcome that is
produced by the twin when initialised at $X_{0}$, which we defined for
$a_{1:t}^{\prime}\in\mathcal{A}_{1:t}$ as

$\widehat{Y}(a_{1:t}^{\prime})\coloneqq
f(X_{0},\widehat{X}_{1:t}(X_{0},a_{1:t}^{\prime})).$ (13)

We then define the analogue of $Q$ for the twin as follows:

$\widehat{Q}\coloneqq\mathbb{E}[\widehat{Y}(a_{1:t})\mid X_{0}\in
B_{0},\widehat{X}_{1:t}(X_{0},a_{1:t})\in B_{1:t}].$

It follows from Proposition 1 that if the twin is interventionally correct,
then $\widehat{Q}=Q$. By substituting $\widehat{Q}$ into Theorem 1, we obtain
hypotheses satisfying (4) via the following immediate corollary to Theorem 1:

######  Corollary 1.

Suppose (5), (6), and (7) all hold. If the twin is interventionally correct,
then the following two hypotheses are true:

$\displaystyle\mathcal{H}_{\mathrm{lo}}:Q_{\mathrm{lo}}\leq\widehat{Q}$
$\displaystyle\mathcal{H}_{\mathrm{up}}:\widehat{Q}\leq Q_{\mathrm{up}}.$

We propose a statistical hypothesis testing approach for falsifying
$\mathcal{H}_{\mathrm{lo}}$ and $\mathcal{H}_{\mathrm{up}}$ using a finite
dataset in the next section.

#### Interpreting a falsification

From a practical perspective, falsifying either $\mathcal{H}_{\mathrm{lo}}$ or
$\mathcal{H}_{\mathrm{up}}$ leads to knowledge about a specific failure mode
of the twin with various potential implications downstream. For example, if
$\mathcal{H}_{\mathrm{lo}}$ is false (i.e. if $\widehat{Q}<Q_{\mathrm{lo}}$),
it follows that, among those trajectories for which
$(X_{0},\widehat{X}_{1:t}(X_{0},a_{1:t}))\in B_{0:t}$, the corresponding value
of $\widehat{Y}(a_{1:t})$ will be on average too small. In light of this, a
user might choose not to rely on outputs of the twin produced under these
circumstances, while a developer seeking to improve the twin could focus their
attention on the specific parts of its implementation that give rise to this
behaviour. We illustrate this concretely through our case study in Section 7.

##  6 Statistical methodology

In practice, the values of $Q_{\mathrm{lo}}$, $Q_{\mathrm{up}}$, and
$\widehat{Q}$ will not be known exactly, and must be estimated from data. To
obtain reliable falsifications, it is necessary to control for the finite-
sample error that estimation introduces. The framework of statistical
hypothesis testing, whose explicit aim is to reliably determine when certain
hypotheses are false using data, is particularly well suited for this. In this
section, we provide a hypothesis testing procedure that can be used to falsify
$\mathcal{H}_{\mathrm{lo}}$ and $\mathcal{H}_{\mathrm{up}}$ with precise
control over the probability of drawing false conclusions about the twin.

###  6.1 Testing a single hypothesis

We first describe how to test $\mathcal{H}_{\mathrm{lo}}$ for a fixed choice
of parameters $(t,f,a_{1:t},B_{0:t})$. Our procedure for testing
$\mathcal{H}_{\mathrm{up}}$ is symmetrical, and indeed could be recovered as a
special case of testing $\mathcal{H}_{\mathrm{lo}}$ by replacing $f$ with
$-f$. Further below, we then consider how to test multiple hypotheses across
different choices of these parameters. We will assume throughout that (9)
holds. In practice, as mentioned above, this can easily be achieved by
choosing $f$ to be bounded, and we give a concrete example of such a choice
when we detail our experimental setup.

#### Datasets

As above, we assume access to a finite dataset $\mathcal{D}$ of i.i.d. copies
of (2). We also assume access to a dataset $\widehat{\mathcal{D}}(a_{1:t})$ of
i.i.d. copies of

$X_{0},\widehat{X}_{1}(X_{0},a_{1}),\ldots,\widehat{X}_{t}(X_{0},a_{1:t}).$
(14)

In practice, these copies can be obtained by initialising the twin with some
value $X_{0}$ taken from $\mathcal{D}$ without replacement, supplying inputs
$a_{1:t}$, and recording the output of the twin at each step. If each $X_{0}$
in $\mathcal{D}$ is used to initialise the twin at most once, then the
resulting trajectories in $\widehat{\mathcal{D}}(a_{1:t})$ are guaranteed to
be i.i.d., since we assumed in Section 3 that the potential outcomes
$\widehat{X}_{t}(x_{0},a_{1:t})$ produced by the twin are independent across
runs. We adopt this approach in our case study.

#### Positivity conditions

For $\mathcal{H}_{\mathrm{lo}}$ to be well-defined, both $a_{1:t}$ and
$B_{0:t}$ must satisfy the positivity conditions in (5) and (6). To ensure
this holds, our first step in the testing process is simply to check whether
there at least one trajectory in $\mathcal{D}$ with $A_{1:t}=a_{1:t}$ and
$X_{0:t}(A_{1:t})\in B_{0:t}$, and one trajectory in
$\widehat{\mathcal{D}}(a_{1:t})$ with
$(X_{0},\widehat{X}_{1:t}(X_{0},a_{1:t}))\in B_{0:t}$. If there are not, then
we refuse to reject $\mathcal{H}_{\mathrm{lo}}$ at any significance level
$\alpha$; otherwise, we proceed to test $\mathcal{H}_{\mathrm{lo}}$ as
described in the next subsection. In effect, this means our overall procedure
tests the augmented hypothesis “(5) and (6) implies
$\mathcal{H}_{\mathrm{lo}}$”. It is easily checked that this new hypothesis
has our desired property in (4) for any choice of $a_{1:t}$ and $B_{0:t}$.
Moreover, this sanity check means there is zero probability that we will
proceed to test $\mathcal{H}_{\mathrm{lo}}$ if (5) and (6) do not hold.

#### Testing procedure

We begin by constructing a one-sided lower confidence interval for
$Q_{\mathrm{lo}}$, and a one-sided upper confidence interval for
$\widehat{Q}$. In particular, for each _significance level_ $\alpha\in(0,1)$,
we use $\mathcal{D}$ and $\widehat{\mathcal{D}}(a_{1:t})$ to obtain
$R_{\mathrm{lo}}^{\alpha}$ and $\widehat{R}^{\alpha}$ such that

$\displaystyle\mathbb{P}(Q_{\mathrm{lo}}\geq R_{\mathrm{lo}}^{\alpha})\geq
1-\frac{\alpha}{2}$
$\displaystyle\mathbb{P}(\widehat{Q}\leq\widehat{R}^{\alpha})\geq
1-\frac{\alpha}{2}.$ (15)

We will also ensure that these are nested in the sense that
$R_{\mathrm{lo}}^{\alpha}\leq R_{\mathrm{lo}}^{\alpha^{\prime}}$ and
$\widehat{R}^{\alpha^{\prime}}\leq\widehat{R}^{\alpha}$ if
$\alpha\leq\alpha^{\prime}$. We describe two methods for obtaining
$R_{\mathrm{lo}}^{\alpha}$ and $\widehat{R}^{\alpha}$ satisfying these
conditions below.

From these confidence intervals, we obtain a test for the hypothesis
$\mathcal{H}_{\mathrm{lo}}$ that rejects when
$\widehat{R}^{\alpha}<R_{\mathrm{lo}}^{\alpha}$. By a union bound argument
given in Section G of the Appendix, it is straightforward to show that this
test has the correct size $\alpha$. In other words, there is probability at
most $\alpha$ of inferring something incorrect about the twin via this
procedure. Our nestedness assumption also gives rise to a $p$-value obtained
as the smallest value of $\alpha$ for which this test rejects, i.e.

$p_{\textup{lo}}\coloneqq\inf\\{\alpha\in(0,1)\mid\text{$\widehat{R}^{\alpha}<R_{\mathrm{lo}}^{\alpha}$}\\},$
(16)

or $1$ if this set is empty, which occurs if the test does not reject for any
value of $\alpha$. The $p$-value provides a measure of the _degree_ to which
the data contradict the hypothesis, with lower values indicating greater
evidence in favour of rejection.

#### Confidence intervals

We now consider how to obtain confidence intervals for $Q_{\mathrm{lo}}$ and
$\widehat{Q}$ satisfying the desired conditions above. To this end, observe
that both quantities are (conditional) expectations involving random variables
that can be computed from $\mathcal{D}$ or $\widehat{\mathcal{D}}(a_{1:t})$.
This allows both to be estimated unbiasedly, which in turn can be used to
derive confidence intervals via standard techniques. For example, consider the
subset of trajectories in $\widehat{\mathcal{D}}(a_{1:t})$ with
$(X_{0},\widehat{X}_{t}(X_{0},a_{1:t}))\in B_{0:t}$. For each such trajectory,
we obtain a corresponding value of $\widehat{Y}(a_{1:t})$ from (13). It is
then easily seen that each $\widehat{Y}(a_{1:t})$ obtained in this way is
i.i.d. and has expectation $\widehat{Q}$. Similarly, for $Q_{\mathrm{lo}}$, we
extract the subset of trajectories in $\mathcal{D}$ for which
$X_{0:N}(A_{1:N})\in B_{0:N}$ holds, where $N$ is computed for each trajectory
via (1). The values of $Y_{\mathrm{lo}}$ obtained from each such trajectory
are then i.i.d. and have expectation $Q_{\mathrm{lo}}$.

At this point, our problem now reduces to that of constructing a confidence
interval for the expectation of a random variable using i.i.d. copies of it.
Various techniques exist for this, and we consider two possibilities in our
case study. The first leverages the fact that $\widehat{Y}(a_{1:t})$ and
$Y_{\mathrm{lo}}$ are bounded in $[y_{\mathrm{lo}},y_{\mathrm{up}}]$, which
gives rise to $R_{\mathrm{lo}}^{\alpha}$ and $\widehat{R}^{\alpha}$ via an
application of Hoeffding’s inequality. This approach has the appealing
property that (15) holds exactly, although may be overly conservative as we
find in our case study. To address this, we also consider confidence intervals
obtained using the bootstrap [16, 71, 15, 22]. Unlike Hoeffding’s inequality,
bootstrapping is only exact in the limit as the size of
$\widehat{\mathcal{D}}(a_{1:t})$ and $\mathcal{D}(a_{1:t})$ approaches
infinity, although several variants of the procedure exist that are first or
second-order accurate [19], which constitutes acceptable coverage in many
applications. Moreover, as we show in our case study, bootstrapping produces
more informative confidence intervals than our exact approach, and so in
practice is more effective when used to falsify the twin. Section G of the
Appendix describes both methods in greater detail.

###  6.2 Testing multiple hypotheses

The tests just described control the probability of falsely rejecting any
individual hypothesis. In practice, when considering multiple hypotheses, it
is often desirable to control the rate of false rejections across all tests
performed. For this we can apply standard techniques for multiple testing. In
our case study, we use the Holm-Bonferroni [24] method, which controls the
family-wise error rate across all tests. If it is instead preferred to control
the false discovery rate, then this could easily be replaced by (for example)
the Benjamini–Yekutieli [7] procedure. Both are easily carried out using the
collection of $p$-values obtained from all hypotheses tested across all
choices of parameters $(t,f,a_{1:t},B_{0:t})$.888This includes hypotheses that
were not rejected because (5) and (6) could not be verified, which
automatically give rise to a $p$-value of $1$. Moreover, both procedures are
valid under arbitrary dependencies among the $p$-values, which means the same
data in $\mathcal{D}$ and $\widehat{\mathcal{D}}(a_{1:t})$ may be reused
across any number of individual tests performed.

###  6.3 Choosing hypotheses from data

In practice, rather than specifying the parameters of each hypothesis
manually, it may be convenient to do so automatically in terms of some data.
However, unless care is taken, this can introduce spurious dependencies
between the tests chosen and the subsequent outcomes of those tests,
invalidating the statistical guarantees of our procedure above. A
straightforward and standard way to avoid this possibility is by _sample
splitting_ [13], which selects hypotheses to test using an independent, held-
out dataset that is subsequently discarded. In other words, data used to
define hypotheses is not reused as part of $\mathcal{D}$ or
$\widehat{\mathcal{D}}(a_{1:t})$ when actually testing any these hypotheses.
Subject to this restriction, our methodology above remains statistically valid
when parameters are chosen depending on data in this way. We selected each
$f$, $a_{1:t}$, and $B_{0:t}$ using held-out data in our case study, and
describe our procedure for doing so below.

##  7 Case Study: Pulse Physiology Engine

We used our proposed methodology to assess sepsis modelling within the Pulse
Physiology Engine [8], a large-scale, open-source model for human physiology
simulation. Pulse generates synthetic trajectories of physiological quantities
including vital signs like heart rate and respiratory rate, as well as lab
values like blood concentration of potassium and sodium, for patients with
medical conditions such as sepsis, COPD, ARDS, etc. To initialise a patient
trajectory, users provide a set of parameters including sex, initial weight,
heart rate, systolic and diastolic blood pressure, and respiratory rate. At
runtime, users may then perform actions such as administering drugs, fluids,
and nutrients, and having the patient exercise. We assessed Pulse using our
methodology and identified a large number of specific instances in which its
output it is inaccurate. We describe our experimental procedure and results
now. Additional details are given in Section H of the Appendix, and code to
reproduce our experiments can be found at
github.com/faaizT/CausalTwinAssessment.

###  7.1 Observational data

We used the MIMIC-III dataset [29] as the observational dataset for this case
study. MIMIC-III is a database of longitudinal health data taken from over
40,000 patients who stayed in the critical care units of the Beth Israel
Deaconess Medical Center between 2001 and 2012. Since our goal was to assess
the accuracy of Pulse for sepsis modelling, we first extracted the subset of
patients conforming to the sepsis-3 criteria [67]. We did so using the same
procedure as [35], who also consider the MIMIC-III dataset. We detail this
procedure in Section H.1 of the Appendix. In total, we obtained 11,677 sepsis
patient trajectories after preprocessing. Before proceeding further, we
randomly partitioned the full dataset into a held-out subset $\mathcal{D}_{0}$
of 583 trajectories (5% of the total), which we used to choose our hypotheses
of interest, and a subset $\mathcal{D}$ of the remaining 11094 trajectories,
which we used then to test these hypotheses.

#### Observation and action spaces

We considered hourly observations of each patient over the the first four
hours of their ICU stay, i.e. $T=4$. We defined the observation spaces
$\mathcal{X}_{0:T}$ using a total of 17 features included in our extracted
MIMIC trajectories. These features including static demographic quantities
(e.g. patient age), as well as patient vitals and laboratory data taken hourly
over the first 4 hours of their ICU stay (i.e. $T=4$). See Section H.2 of the
Appendix for a complete list. Following [35], we considered actions based on
the administration of intravenous fluids and vasopressors, which both play a
primary role in the treatment of sepsis in clinical practice. Since these are
recorded as continuous quantities in MIMIC, in order to obtain a finite action
space suitable for our framework, we discretized their values via the same
procedure as [35], which we detail in H.3 of the Appendix. This gave us action
spaces $\mathcal{A}_{1}=\cdots=\mathcal{A}_{4}$, each with $25$ distinct
actions.

#### Hypothesis parameters

We used $\mathcal{D}_{0}$ to choose a collection of hypothesis parameters
$(t,f,a_{1:t},B_{0:t})$, each of which we then used to define an
$\mathcal{H}_{\mathrm{lo}}$ and $\mathcal{H}_{\mathrm{up}}$ to test.
Specifically, for each of the 14 different physiological quantities listed in
Table 1, we did the following. For each $t\in\\{0,\ldots,4\\}$, we obtained 16
choices of $B_{t}$ by discretizing $\mathcal{X}_{t}$ into 16 bins based on
sex, age, and the physiological quantity of interest at time $t$. We detail
this discretization procedure in Section H.4 of the Appendix. We then
collected all triples $(t,a_{1:t},B_{0:t})$ with $t\in\\{1,\ldots,4\\}$ and
$a_{1:t}\in\mathcal{A}_{1:t}$ that had at least one corresponding trajectory
$(X_{0:T}(A_{1:T}),A_{1:T})\in\mathcal{D}_{0}$, i.e. with $A_{1:t}=a_{1:t}$
and $X_{0:t}(A_{1:T})\in B_{0:t}$. For each triple $(t,a_{1:t},B_{0:t})$, we
set $y_{\mathrm{lo}}$ and $y_{\mathrm{up}}$ to be the .2 and .8 quantiles of
the physiological quantity of interest at time $t$ across all trajectories in
$\mathcal{D}_{0}$ corresponding to that triple. We then defined $f$ for that
triple as the function that extracts the physiological quantity of interest
from $\mathcal{X}_{t}$ and clips its value to between $y_{\mathrm{lo}}$ and
$y_{\mathrm{up}}$, which ensured that (7) held by construction. Overall, we
obtained 721 distinct choices of $(t,f,a_{1:t},B_{0:t})$ in this way. Each of
these produced two hypotheses $\mathcal{H}_{\mathrm{lo}}$ and
$\mathcal{H}_{\mathrm{up}}$, which gave a total of 1,442 hypotheses overall.

###  7.2 Digital twin trajectories

We now describe how we generated twin data to test the chosen hypotheses. For
each action sequence $a_{1:t}$ selected according to the procedure above, we
obtained the dataset $\widehat{\mathcal{D}}(a_{1:t})$ defined in Section 6.
Specifically, to obtain i.i.d. copies of (14), we sampled $X_{0}$ without
replacement from $\mathcal{D}$, and used this to initialise a twin trajectory.
At each hour $t^{\prime}\in\\{1,\ldots,t\\}$ in simulation time, we
administered a dose of intravenous fluids and vasopressors equal to the median
dose value observed in the bin corresponding to $a_{t^{\prime}}$ in the held-
out dataset $\mathcal{D}_{0}$. Additionally, after each hour of simulation
time, we recorded observations corresponding to the same 19 features as we
extracted from MIMIC. We generated 26,115 simulated trajectories in total,
which took approximately 48 hours running on 100 parallel CPU threads.

###  7.3 Results

We tested the chosen hypotheses using our methodology from Section 6. Here we
report the results when using Hoeffding’s inequality to obtain confidence
intervals for $Q_{\mathrm{lo}}$, $Q_{\mathrm{up}}$, and $\widehat{Q}$. We also
tried confidence intervals obtained via bootstrapping, and obtained similar
(if less conservative) results that we report in Section H of the Appendix. We
used the Holm-Bonferroni method to adjust for multiple testing, with our
overall family-wise error rate taken to be 0.05.

Physiological quantity | # Rejections | # Hypotheses
---|---|---
Chloride Blood Concentration (Chloride) | 24 | 94
Sodium Blood Concentration (Sodium) | 21 | 94
Potassium Blood Concentration (Potassium) | 13 | 94
Skin Temperature (Temp) | 10 | 86
Calcium Blood Concentration (Calcium) | 5 | 88
Glucose Blood Concentration (Glucose) | 5 | 96
Arterial CO2 Pressure (paCO2) | 3 | 70
Bicarbonate Blood Concentration (HCO3) | 2 | 90
Systolic Arterial Pressure (SysBP) | 2 | 154
Arterial O2 Pressure (paO2) | 0 | 78
Arterial pH (Arterial_pH) | 0 | 80
Diastolic Arterial Pressure (DiaBP) | 0 | 72
Mean Arterial Pressure (MeanBP) | 0 | 92
Respiration Rate (RR) | 0 | 172
Heart Rate (HR) | 0 | 162
Table 1: Total hypotheses and rejections per physiological quantity

#### Hypothesis rejections

Table 1 shows the number of hypotheses that were rejected for each
physiological quantity of interest. We obtained rejections for 10 different
quantities, the majority of which were for sodium, chloride, and potassium
levels, and skin temperature. We may therefore infer that, at a high level,
Pulse does not simulate these quantities accurately for the population of
sepsis patients we consider. This appears of interest in a variety of
downstream settings: for example, a developer could use this information as a
natural starting point from which to improve the accuracy of Pulse, while a
practitioner using Pulse to inform their decision-making may wish to rely less
on these outputs as a result.

####  $p$-value plots

To obtain more granular information about the failure modes of the twin just
identified, we examined the $p$-values obtained for each hypothesis
$\mathcal{H}_{\mathrm{lo}}$ and $\mathcal{H}_{\mathrm{up}}$ tested, which we
denote here by $p_{\textup{lo}}$ and $p_{\textup{up}}$. Figure 2 shows the
distributions of $-\log_{10}{p_{\textup{lo}}}$ and
$-\log_{10}{p_{\textup{up}}}$ that we obtained for all physiological
quantities for which some hypothesis was rejected. (The remaining $p$-values
are shown in Figure 7 in the Appendix.) Notably, in each row, one distribution
is always tightly concentrated at $-\log_{10}p=0$ (i.e. $p=1$). This means
that, for all physiological outcomes of interest, there was either very little
evidence in favour of rejecting _any_ $\mathcal{H}_{\mathrm{lo}}$, or very
little in favour of rejecting any $\mathcal{H}_{\mathrm{up}}$. In other words,
across configurations of $(t,f,a_{1:t},B_{0:t})$ that were rejected, the twin
consistently either underestimated or overestimated each quantity on average.
For example, Pulse consistently underestimated chloride blood concentration
and skin temperature, while it consistently overestimated sodium and glucose
blood concentration levels. Like Table 1, this information appears of interest
and actionable in a variety of downstream tasks.

Figure 2: Distributions of $-\log_{10}{p_{\textup{lo}}}$ and
$-\log_{10}{p_{\textup{up}}}$ across hypotheses, grouped by physiological
quantity. Higher values indicate greater evidence in favour of rejection.

###  7.4 Pitfalls of naive assessment

A naive approach to twin assessment involves simply comparing the output of
the twin with observational data directly. We now show that, unlike our
methodology, the results produced in this way can be misleading. In Figure 3,
for two different choices of $(a_{1:4},B_{1:4})$, we plot estimates of
$\widehat{Q}_{t}$ and $Q^{\textup{obs}}_{t}$ for $t\in\\{1,\ldots,4\\}$, where

$\displaystyle\widehat{Q}_{t}$
$\displaystyle\coloneqq\mathbb{E}[\widehat{Y}(a_{1:t})\mid\widehat{X}_{0:t}(a_{1:t})\in
B_{0:t}],$ $\displaystyle Q^{\textup{obs}}_{t}$
$\displaystyle\coloneqq\mathbb{E}[Y(A_{1:t})\mid X_{0:t}(A_{1:t})\in
B_{0:t},A_{1:t}=a_{1:t}].$

Each plot also shows one-sided 95% confidence intervals on $Q_{\mathrm{lo}}$
and $Q_{\mathrm{up}}$ at each $t\in\\{1,\ldots,4\\}$ obtained from Hoeffding’s
inequality. Directly comparing the estimates of $\widehat{Q}_{t}$ and
$Q^{\textup{obs}}_{t}$ would suggest that the twin is comparatively more
accurate for the right-hand plot, as these estimates are closer to one another
in that case. However, our methodology falsifies the output of the twin for
the right-hand plot at $t=1$, as can be seen from the fact that confidence
interval for $\widehat{Q}_{1}$ lies entirely above the one-sided confidence
interval for $Q_{\mathrm{up}}$ at that timestep. On the other hand, the output
of the twin is not falsified at any of the timesteps shown for the left-hand
plot, so that the twin may in fact be accurate for these $(a_{1:4},B_{1:4})$,
contrary to what a naive assessment strategy would suggest. Our methodology
provides a principled means for twin assessment that avoids drawing
potentially misleading inferences of this kind.

A similar phenomenon appears in Figure 4, which, for two choices of $B_{0:t}$
and $a_{1:t}$, shows histograms of raw glucose values obtained from the
observational data conditional on $A_{1:t}=a_{1:t}$ and $X_{0:t}(A_{1:t})\in
B_{0:t}$, and from the twin conditional on $\widehat{X}_{0:t}(a_{1:t})\in
B_{0:t}$. (Note that these raw values differ from $Y(A_{1:t})$ and
$\widehat{Y}(a_{1:t})$ since they are not clipped to lie between
$y_{\mathrm{lo}}$ and $y_{\mathrm{up}}$.) Below each histogram we also show
95% confidence intervals for $Q_{\mathrm{lo}}$ and $\widehat{Q}$ obtained from
Hoeffding’s inequality. While Figures 4(a) and 4(b) appear visually very
similar, the inferences produced by our testing procedure are different: the
hypothesis corresponding to the right-hand plot is rejected, since there is no
overlap between the confidence intervals underneath, while the hypothesis
corresponding to the left-hand plot is not. This was not an isolated case and
several other examples of this phenomenon are shown in Figure 9 in the
Appendix. This demonstrates that our methodology does not simply rely on
direct comparison of observed and simulated trajectories, but also accounts
for the possibility of arbitrary confounding in the data.

Figure 3: Estimates and 95% confidence intervals for $\widehat{Q}_{t}$ and
$Q^{\textup{obs}}_{t}$ at each $t\in\\{1,\ldots,4\\}$ for two choices of
$(B_{0:4},a_{1:4})$, where $\widehat{Y}(a_{1:t})$ and $Y(a_{1:t})$ correspond
to bicarbonate blood concentration. The dashed lines indicate one-sided lower
and upper 95% confidence intervals for $Q_{\mathrm{lo}}$ and $Q_{\mathrm{up}}$
respectively.

(a) Not rejected

(b) Rejected

Figure 4: Raw glucose values from the observational data and twin for two
choices of $(B_{0:t},a_{1:t})$, with confidence intervals for $\widehat{Q}$
and $Q_{\mathrm{lo}}$ shown below. The horizontal axes are truncated to the
.025 and .975 quantiles of the observational data to make visually clearer
whether or not the confidence intervals overlap.

##  8 Discussion

We have advocated for a causal approach to digital twin assessment, and have
presented a statistical procedure for doing so that obtains rigorous
theoretical guarantees under minimal assumptions. We now highlight the key
limitations of our approach. Importantly, our methodology implicitly assumes
that the interventional distribution of the real-world process does not differ
between when our dataset was obtained and when the twin is deployed to
production. If the conditional distribution of $X_{1:T}(a_{1:T})$ given
$X_{0}$ changes at deployment time, then so too does the set of twins that are
interventionally correct, and if this change is significant enough, our
assessment procedure may yield misleading results. This may occur in cases
where $X_{0}$ encodes previous actions that were taken before $t=0$, which is
useful in certain decision-support cases (see Section C of the Appendix). In
such cases, distribution shift can arise if the policy used to determine these
historical actions changes at deployment time. Distribution shift in this
sense is a separate issue to unobserved confounding, and arises in a wide
variety of statistical problems beyond ours.

Additionally, the procedure we used in our case study to choose the hypothesis
parameters $B_{0:t}$ was ad hoc and somewhat crude. For scalability, it would
likely be necessary to obtain $B_{0:t}$ via a more automated procedure. It may
also be desirable to choose $B_{0:t}$ dynamically in light of previous
hypotheses tested, zooming in to regions containing possible failure modes to
obtain increasingly granular information about the twin. We see opportunities
here for using machine learning techniques, but leave developing this further
to future work.

Various other extensions and improvements appear possible. For example, it is
possible to replace our one-sided confidence intervals for $Q_{\mathrm{lo}}$,
$Q_{\mathrm{up}}$, and $\widehat{Q}$ with two-sided ones, and thereby to
obtain a procedure that may yield more precise information about the twin than
we obtain by rejecting one of $\mathcal{H}_{\mathrm{lo}}$ or
$\mathcal{H}_{\mathrm{up}}$. We outline this at a high level in Section G.3 of
the Appendix. It also seems possible to leverage ideas from the literature on
partial identification [45] to obtain greater statistical efficiency, for
example by building on the line of work initiated by [27] for obtaining more
informative confidence intervals. Beyond this, it may be useful in some
contexts to consider additional assumptions that lead to less conservative
assessment results. For example, various methods for _sensitivity analysis_
have been proposed that model the _degree_ to which the actions of the
behavioural agent are confounded [60, 68, 74, 54]. This can yield tighter
bounds on $Q$ than are implied by Theorem 1, albeit at the expense of less
robustness if these assumptions are violated.

## Acknowledgements

RC, AD, and CC were supported by the Engineering and Physical Sciences
Research Council (EPSRC) through the Bayes4Health programme [Grant number
EP/R018561/1]. MFT acknowledges PhD funding from Google DeepMind. The authors
are also highly appreciative of troubleshooting and development assistance
provided by the Pulse team.

## References

  * [1] AMSE.  Assessing Credibility of Computational Modeling through Verification and Validation: Application to Medical Devices.  AMSE, 2018.
  * [2] A. Balke and J. Pearl.  Bounds on treatment effects from studies with imperfect compliance.  Journal of the American Statistical Association, 92(439):1171–1176, 1997.
  * [3] B. R. Barricelli, E. Casiraghi, and D. Fogli.  A survey on digital twin: Definitions, characteristics, applications, and design implications.  IEEE access, 7:167653–167671, 2019.
  * [4] D. Basu.  Randomization analysis of experimental data: The fisher randomization test.  Journal of the American Statistical Association, 75(371):575–582, 1980.
  * [5] P. Bauer, B. Stevens, and W. Hazeleger.  A digital twin of earth for the green transition.  Nature Climate Change, 11(2):80–83, 2021.
  * [6] N. Bellinger, E. J. Tuegel, A. R. Ingraffea, T. G. Eason, and S. M. Spottswood.  Reengineering aircraft structural life prediction using a digital twin.  International Journal of Aerospace Engineering, 2011:154798, 2011\.
  * [7] Y. Benjamini and D. Yekutieli.  The control of the false discovery rate in multiple testing under dependency.  The Annals of Statistics, 29(4):1165 – 1188, 2001.
  * [8] A. Bray, J. B. Webb, A. Enquobahrie, J. Vicory, J. Heneghan, R. Hubal, S. TerMaath, P. Asare, and R. B. Clipp.  Pulse physiology engine: an open-source software platform for computational modeling of human medical simulation.  SN Comprehensive Clinical Medicine, 1(5):362–377, 2019.
  * [9] V. Chernozhukov, H. Hong, and E. Tamer.  Estimation and confidence regions for parameter sets in econometric models.  Econometrica, 75(5):1243–1284, 2007.
  * [10] G. Coorey, G. A. Figtree, D. F. Fletcher, and J. Redfern.  The health digital twin: advancing precision cardiovascular medicine.  Nature Reviews Cardiology, 18(12):803–804, 2021.
  * [11] G. Coorey, G. A. Figtree, D. F. Fletcher, V. J. Snelson, S. T. Vernon, D. Winlaw, S. M. Grieve, A. McEwan, J. Y. H. Yang, P. Qian, et al.  The health digital twin to tackle cardiovascular disease—a review of an emerging interdisciplinary field.  NPJ Digital Medicine, 5(1):1–12, 2022.
  * [12] J. Corral-Acero, F. Margara, M. Marciniak, C. Rodero, F. Loncaric, Y. Feng, A. Gilbert, J. F. Fernandes, H. A. Bukhari, A. Wajdan, et al.  The ’digital twin’ to enable the vision of precision cardiology.  European Heart Journal, 41(48):4556–4564, 2020.
  * [13] D. R. Cox.  A note on data-splitting for the evaluation of significance levels.  Biometrika, 62(2):441–444, 1975.
  * [14] U. R. Dahmen, T. Osterloh, and H.-J. Roßmann.  Verification and validation of digital twins and virtual testbeds.  International journal of advances in engineering sciences and applied mathematics, 11(1):47–64, 2022.
  * [15] A. C. Davison and D. V. Hinkley.  Bootstrap Methods and their Application.  Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 1997.
  * [16] B. Efron.  Bootstrap methods: Another look at the jackknife.  The Annals of Statistics, 7(1):1 – 26, 1979.
  * [17] S. Galappaththige, R. A. Gray, C. M. Costa, S. Niederer, and P. Pathmanathan.  Credibility assessment of patient-specific computational modeling using patient-specific cardiac modeling as an exemplar.  PLoS computational biology, 18(10):e1010541, 2022.
  * [18] M. Grieves and J. Vickers.  Digital twin: Mitigating unpredictable, undesirable emergent behavior in complex systems.  In Transdisciplinary Perspectives on Complex Systems, pages 85–113. Springer, 2017.
  * [19] P. Hall.  Theoretical comparison of bootstrap confidence intervals.  The Annals of Statistics, 16(3):927 – 953, 1988.
  * [20] A. Hemmler, B. Lutz, G. Kalender, C. Reeps, and M. W. Gee.  Patient-specific in silico endovascular repair of abdominal aortic aneurysms: application and validation.  Biomechanics and Modeling in Mechanobiology, 18(4):983–1004, 2019\.
  * [21] M. A. Hernán and J. M. Robins.  Causal Inference: What If.  Chapman and Hall/CRC, Boca Raton, 2020.
  * [22] T. C. Hesterberg.  What teachers should know about the bootstrap: Resampling in the undergraduate statistics curriculum.  The American Statistician, 69(4):371–386, 2015.  PMID: 27019512.
  * [23] P. W. Holland.  Statistics and causal inference.  Journal of the American Statistical Association, 81(396):945–960, 1986.
  * [24] S. Holm.  A simple sequentially rejective multiple test procedure.  Scandinavian Journal of Statistics, 6(2):65–70, 1979.
  * [25] Z. Hussain, M. Oberst, M.-C. Shih, and D. Sontag.  Falsification before extrapolation in causal effect estimation.  In Advances in Neural Information Processing Systems, September 2022\.
  * [26] G. W. Imbens.  Potential outcome and directed acyclic graph approaches to causality: Relevance for empirical practice in economics.  Journal of Economic Literature, 58(4):1129–79, 2020.
  * [27] G. W. Imbens and C. F. Manski.  Confidence intervals for partially identified parameters.  Econometrica, 72(6):1845–1857, 2004.
  * [28] M. Jans-Singh, K. Leeming, R. Choudhary, and M. Girolami.  Digital twin of an urban-integrated hydroponic farm.  Data-Centric Engineering, 1:e20, 2020.
  * [29] A. E. W. Johnson, T. J. Pollard, L. Shen, L.-w. H. Lehman, M. Feng, M. Ghassemi, B. Moody, P. Szolovits, L. Anthony Celi, and R. G. Mark.  Mimic-iii, a freely accessible critical care database.  Scientific Data, 3(1):160035, 2016.
  * [30] D. Jones, C. Snider, A. Nassehi, J. Yon, and B. Hicks.  Characterising the digital twin: A systematic literature review.  CIRP Journal of Manufacturing Science and Technology, 29:36–52, 2020.
  * [31] O. Kallenberg and O. Kallenberg.  Foundations of Modern Probability, volume 2.  Springer, 1997.
  * [32] M. G. Kapteyn, J. V. Pretorius, and K. E. Willcox.  A probabilistic graphical model foundation for enabling predictive digital twins at scale.  Nature Computational Science, 1(5):337–347, 2021.
  * [33] A. Khan, M. Dahl, P. Falkman, and M. Fabian.  Digital twin for legacy systems: Simulation model testing and validation.  In 2018 IEEE 14th International Conference on Automation Science and Engineering (CASE), pages 421–426. IEEE, 2018.
  * [34] B. Kochunas and X. Huan.  Digital twin concepts with uncertainty for nuclear power applications.  Energies, 14(14):4235, 2021.
  * [35] M. Komorowski, L. A. Celi, O. Badawi, A. C. Gordon, and A. A. Faisal.  The artificial intelligence clinician learns optimal treatment strategies for sepsis in intensive care.  Nature Medicine, 24(11):1716–1720, 2018.
  * [36] W. Kritzinger, M. Karner, G. Traar, J. Henjes, and W. Sihn.  Digital twin in manufacturing: A categorical literature review and classification.  IFAC-PapersOnLine, 51(11):1016–1022, 2018.
  * [37] A. Lal, G. Li, E. Cubro, S. Chalmers, H. Li, V. Herasevich, Y. Dong, B. Pickering, K. Oguz, and O. Gajic.  Development and verification of a digital twin patient model to predict treatment response in sepsis.  Critical Care Medicine, 49:611–611, 01 2021.
  * [38] I. Larrabide, M. Kim, L. Augsburger, M. C. Villa-Uriol, D. Rüfenacht, and A. F. Frangi.  Fast virtual deployment of self-expandable stents: method and in vitro evaluation for intracranial aneurysmal stenting.  Medical Image Analysis, 16(3):721—730, April 2012.
  * [39] P. W. Lavori and R. Dawson.  Dynamic treatment regimes: practical design considerations.  Clinical trials, 1(1):9–20, 2004.
  * [40] Y. Lu, C. Liu, I. Kevin, K. Wang, H. Huang, and X. Xu.  Digital twin-driven smart manufacturing: Connotation, reference model, applications and research issues.  Robotics and Computer-Integrated Manufacturing, 61:101837, 2020\.
  * [41] C. F. Manski.  Anatomy of the selection problem.  The Journal of Human Resources, 24(3):343–360, 1989.
  * [42] C. F. Manski.  Nonparametric bounds on treatment effects.  The American Economic Review, 80(2):319–323, 1990.
  * [43] C. F. Manski.  Identification Problems in the Social Sciences.  Harvard University Press, 1995.
  * [44] C. F. Manski.  Monotone treatment response.  Econometrica, 65(6):1311–1334, 1997.
  * [45] C. F. Manski.  Partial Identification of Probability Distributions, volume 5.  Springer, 2003.
  * [46] C. F. Manski.  Identification for Prediction and Decision.  Harvard University Press, 2009.
  * [47] C. F. Manski and J. V. Pepper.  Monotone instrumental variables: With an application to the returns to schooling.  Econometrica, 68(4):997–1010, 2000.
  * [48] R. Marcus, E. Peritz, and K. R. Gabriel.  On closed testing procedures with special reference to ordered analysis of variance.  Biometrika, 63(3):655–660, 1976.
  * [49] J. Masison, J. Beezley, Y. Mei, H. A. L. Ribeiro, A. C. Knapp, L. Sordo Vieira, B. Adhikari, Y. Scindia, M. Grauer, B. Helba, et al.  A modular computational framework for medical digital twins.  Proceedings of the National Academy of Sciences, 118(20):e2024287118, 2021.
  * [50] O. Mazumder, D. Roy, S. Bhattacharya, A. Sinha, and A. Pal.  Synthetic ppg generation from haemodynamic model with baroreflex autoregulation: a digital twin of cardiovascular system.  In 2019 41st Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), pages 5024–5029, 2019.
  * [51] McDaniel, Matthew, and A. Baird.  A full-body model of burn pathophysiology and treatment using the biogears engine.  41st Annual International Conference of the IEEE Engineering in Medicine and Biology Society (EMBC), 2019.
  * [52] S. A. Murphy.  Optimal dynamic treatment regimes.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 65(2):331–355, 2003.
  * [53] S. A. Murphy.  An experimental design for the development of adaptive treatment strategies.  Statistics in Medicine, 24(10):1455–1481, 2005.
  * [54] H. Namkoong, R. Keramati, S. Yadlowsky, and E. Brunskill.  Off-policy policy evaluation for sequential decisions under unobserved confounding.  In H. Larochelle, M. Ranzato, R. Hadsell, M. Balcan, and H. Lin, editors, Advances in Neural Information Processing Systems, volume 33, pages 18819–18831. Curran Associates, Inc., 2020.
  * [55] S. A. Niederer, Y. Aboelkassem, C. D. Cantwell, C. Corrado, S. Coveney, E. M. Cherry, T. Delhaas, F. H. Fenton, A. V. Panfilov, P. Pathmanathan, G. Plank, M. Riabiz, C. H. Roney, R. W. dos Santos, and L. Wang.  Creation and application of virtual patient cohorts of heart models.  Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences, 378(2173):20190558, 2020.
  * [56] S. A. Niederer, M. S. Sacks, M. Girolami, and K. Willcox.  Scaling digital twins from the artisanal to the industrial.  Nature Computational Science, 1(5):313–320, 2021.
  * [57] J. Pearl.  Causality.  Cambridge University Press, 2 edition, 2009.
  * [58] K. Popper.  The Logic of Scientific Discovery.  Routledge, 2005.
  * [59] J. Robins.  A new approach to causal inference in mortality studies with a sustained exposure period—application to control of the healthy worker survivor effect.  Mathematical Modelling, 7(9-12):1393–1512, 1986.
  * [60] P. R. Rosenbaum.  Observational Studies.  Springer, New York, NY, 2002.
  * [61] P. R. Rosenbaum and D. B. Rubin.  The central role of the propensity score in observational studies for causal effects.  Biometrika, 70(1):41–55, 1983.
  * [62] C. J. Roy and W. L. Oberkampf.  A comprehensive framework for verification, validation, and uncertainty quantification in scientific computing.  Computer Methods in Applied Mechanics and Engineering, 200(25):2131–2144, 2011.
  * [63] D. B. Rubin.  Estimating causal effects of treatments in randomized and nonrandomized studies.  Journal of Educational Psychology, 66:688–701, 1974.
  * [64] D. B. Rubin.  Causal inference using potential outcomes.  Journal of the American Statistical Association, 100(469):322–331, 2005.
  * [65] R. Sacks, I. Brilakis, E. Pikas, H. S. Xie, and M. Girolami.  Construction with digital twin information systems.  Data-Centric Engineering, 1, 2020.
  * [66] C. W. Seymour, V. X. Liu, T. J. Iwashyna, F. M. Brunkhorst, T. D. Rea, A. Scherag, G. Rubenfeld, J. M. Kahn, M. Shankar-Hari, M. Singer, C. S. Deutschman, G. J. Escobar, and D. C. Angus.  Assessment of Clinical Criteria for Sepsis: For the Third International Consensus Definitions for Sepsis and Septic Shock (Sepsis-3).  JAMA, 315(8):762–774, 02 2016.
  * [67] M. Singer, C. S. Deutschman, C. W. Seymour, M. Shankar-Hari, D. Annane, M. Bauer, R. Bellomo, G. R. Bernard, J.-D. Chiche, C. M. Coopersmith, R. S. Hotchkiss, M. M. Levy, J. C. Marshall, G. S. Martin, S. M. Opal, G. D. Rubenfeld, T. van der Poll, J.-L. Vincent, and D. C. Angus.  The Third International Consensus Definitions for Sepsis and Septic Shock (Sepsis-3).  JAMA, 315(8):801–810, 02 2016.
  * [68] Z. Tan.  A distributional approach for causal inference using propensity scores.  Journal of the American Statistical Association, 101(476):1619–1637, 2006.
  * [69] F. Tao, H. Zhang, A. Liu, and A. Y. Nee.  Digital twin in industry: State-of-the-art.  IEEE Transactions on Industrial Informatics, 15(4):2405–2415, 2018\.
  * [70] C. Thuemmler and C. Bai.  Health 4.0: Application of Industry 4.0 Design Principles in Future Asthma Management, pages 23–37.  Springer International Publishing, Cham, 2017.
  * [71] R. J. Tibshirani and B. Efron.  An introduction to the bootstrap.  Monographs on statistics and applied probability, 57:1–436, 1993\.
  * [72] A. A. Tsiatis, M. Davidian, S. T. Holloway, and E. B. Laber.  Dynamic treatment regimes: Statistical methods for precision medicine.  Chapman and Hall/CRC, 2019.
  * [73] D. Williams.  Probability with Martingales.  Cambridge University Press, 1991.
  * [74] S. Yadlowsky, H. Namkoong, S. Basu, J. Duchi, and L. Tian.  Bounds on the conditional and average treatment effect with unobserved confounding factors.  The Annals of Statistics, 50(5):2587 – 2615, 2022.
  * [75] J. Zhang and E. Bareinboim.  Near-optimal reinforcement learning in dynamic treatment regimes.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
  * [76] R. Y. Zhong, X. Xu, E. Klotz, and S. T. Newman.  Intelligent manufacturing in the context of industry 4.0: A review.  Engineering, 3(5):616–630, 2017.

##  Appendix A Notation

$z_{t:t^{\prime}}$ |  The sequence of elements $(z_{t},\ldots,z_{t^{\prime}})$ (or the empty sequence when $t>t^{\prime}$)
---|---
$\mathcal{Z}_{t:t^{\prime}}$ (where each $\mathcal{Z}_{i}$ is a set) |  The cartesian product $\mathcal{Z}_{t}\times\cdots\times\mathcal{Z}_{t^{\prime}}$ (or the empty set when $t>t^{\prime}$)
$Z_{t:t^{\prime}}(a_{1:t^{\prime}})$ |  The sequence of potential outcomes $Z_{t}(a_{1:t}),\ldots,Z_{t^{\prime}}(a_{1:t^{\prime}})$ (or the empty sequence when $t>t^{\prime}$)
$Z\overset{\mathrm{d}}{=}Z^{\prime}$ |  The distributions of random variables $Z$ and $Z^{\prime}$ are identical
$Z\overset{\mathrm{a.s.}}{=}Z^{\prime}$ |  The random variables $Z$ and $Z^{\prime}$ are almost surely equal, i.e. $\mathbb{P}(Z=Z^{\prime})=1$
$Z\perp\\!\\!\\!\perp Z^{\prime}$ |  The random variables $Z$ and $Z^{\prime}$ are independent
$Z\perp\\!\\!\\!\perp Z^{\prime}\mid Z^{\prime\prime}$ |  The random variables $Z$ and $Z^{\prime}$ are conditionally independent given the random variable $Z^{\prime\prime}$
$\mathbbm{1}(E)$ |  Indicator function of some event $E$

###  A.1 Conditioning

In most places throughout this paper, we only require the elementary
definition of conditioning. That is, for an event $E$ with $\mathbb{P}(E)>0$,
the conditional probability of any other event $E^{\prime}$ is defined to be

$\mathbb{P}(E^{\prime}\mid E)\coloneqq\frac{\mathbb{P}(E^{\prime}\cap
E)}{\mathbb{P}(E)},$

the conditional distribution of any $\mathcal{Z}$-valued random variable $Z$
given $E$ is defined to be the distribution $\mathrm{Law}[Z\mid E]$ such that,
for measurable $A\subseteq\mathcal{Z}$, we have

$\mathrm{Law}[Z\mid E](A)\coloneqq\mathbb{P}(Z\in A\mid E),$

and the conditional expectation of any real-valued random variable
$Z^{\prime}$ is defined to be

$\mathbb{E}[Z^{\prime}\mid
E]\coloneqq\frac{\mathbb{E}[Z^{\prime}\,\mathbbm{1}(E)]}{\mathbb{P}(E)}.$

However, in some places in this Supplement (in particular in Section B and
Part F.4 of Section F), we make use of the general definition that allows
conditioning on the value of a random variable. This reduces to the elementary
definition when conditioning on a discrete random variable, but encompasses
other cases such as conditioning on continuous quantities as well. See e.g.
Chapter 9 of [73] or Chapter 6 of [31] for a review of this topic.

##  Appendix B Proof of Proposition 1 (unconditional form of interventional
correctness)

In this section we prove Proposition 1 from the main text. To account for
measure-theoretic technicalities, we first clarify slightly our earlier
definition of interventional correctness. We will assume here that
$\mathcal{X}_{0:T}$ is Borel, which means that $X_{0:T}(a_{1:T})$ admits a
regular conditional distribution given $X_{0}$ [31, Theorem 6.3]. By
interventional correctness, we then mean that, for all
$a_{1:T}\in\mathcal{A}_{1:T}$, the map
$(x_{0},B_{1:T})\mapsto\mathrm{Law}[\widehat{X}_{1:T}(x_{0},a_{1:T})](B_{1:T})$
is a version of this conditional distribution, i.e. it is a Markov kernel such
that

$\mathrm{Law}[\widehat{X}_{1:T}(x_{0},a_{1:T})]=\mathrm{Law}[X_{1:T}(a_{1:T})\mid
X_{0}=x_{0}]\qquad\text{for $\mathrm{Law}[X_{0}]$-almost all
$x_{0}\in\mathcal{X}_{0}$}.$ (17)

The result is then as follows.

######  Proposition 1.

The twin is interventionally correct if and only if, for all
$a_{1:T}\in\mathcal{A}_{1:T}$, it holds that

$\mathrm{Law}[X_{0},\widehat{X}_{1:T}(X_{0},a_{1:T})]=\mathrm{Law}[X_{0:T}(a_{1:T})]$
(18)

###### Proof.

Fix any choice of $a_{1:T}\in\mathcal{A}_{1:T}$. Since both
$\mathrm{Law}[X_{0},\widehat{X}_{1:T}(X_{0},a_{1:T})]$ and
$\mathrm{Law}[X_{0:T}(a_{1:T})]$ have the same $\mathcal{X}_{0}$-marginal,
namely $\mathrm{Law}[X_{0}]$, (18) holds if and only if

$\mathrm{Law}[\widehat{X}_{1:T}(X_{0},a_{1:T})\mid
X_{0}=x_{0}]=\mathrm{Law}[X_{1:T}(a_{1:T})\mid X_{0}=x_{0}]\qquad\text{for
$\mathrm{Law}[X_{0}]$-almost all $x_{0}\in\mathcal{X}_{0}$}.$ (19)

But now, our definition of $\widehat{X}_{1:T}(x_{0},a_{1:T})$ in terms of
$h_{t}$ and $U_{1:t}$ means we can write

$\widehat{X}_{1:T}(X_{0},a_{1:T})=\bar{h}(X_{0},a_{1:T},U_{1:T}),$

where
$\bar{h}(x_{0},a_{1:T},u_{1:T})\coloneqq(h_{1}(x_{0},a_{1},u_{1}),\ldots,h_{T}(x_{0},a_{1:T},u_{1:T}))$.
For all $x_{0}\in\mathcal{X}_{0}$ and measurable
$B_{1:T}\subseteq\mathcal{X}_{1:T}$, we then have

$\displaystyle\mathrm{Law}[\widehat{X}_{1:T}(x_{0},a_{1:T})](B_{1:T})$
$\displaystyle=\mathbb{E}[\mathbbm{1}(\bar{h}(x_{0},a_{1:T},U_{1:T})\in
B_{1:T})]$ $\displaystyle=\int\mathbbm{1}(\bar{h}(x_{0},a_{1:T},u_{1:T})\in
B_{1:T})\,\mathrm{Law}[U_{1:T}](\mathrm{d}u_{1:T}).$

It is standard to show that the right-hand side is a Markov kernel in $x_{0}$
and $B_{1:T}$. Moreover, for any measurable $B_{0}\subseteq\mathcal{X}_{0}$,
we have

$\displaystyle\int_{B_{0}}\mathrm{Law}[\widehat{X}_{1:T}(x_{0},a_{1:T})](B_{1:T})\,\mathrm{Law}[X_{0}](\mathrm{d}x_{0})$
$\displaystyle\qquad=\int_{B_{0}}\left[\int\mathbbm{1}(\bar{h}(x_{0},a_{1:T},u_{1:T})\in
B_{1:T})\,\mathrm{Law}[U_{1:T}](\mathrm{d}u_{1:T})\right]\,\mathrm{Law}[X_{0}](\mathrm{d}x_{0})$
$\displaystyle\qquad=\int\mathbbm{1}(x_{0}\in
B_{0},\bar{h}(x_{0},a_{1:T},u_{1:T})\in
B_{1:T})\,\mathrm{Law}[X_{0},U_{1:T}](\mathrm{d}x_{0},\mathrm{d}u_{1:T})$
$\displaystyle\qquad=\mathrm{Law}[X_{0},\widehat{X}_{1:T}(X_{0},a_{1:T})](B_{0:T}),$

where the second step follows because $X_{0}\perp\\!\\!\\!\perp U_{1:T}$. It
therefore follows that
$(x_{0},B_{1:T})\mapsto\mathrm{Law}[\widehat{X}_{1:T}(x_{0},a_{1:T})](B_{1:T})$
is a regular conditional distribution of $\widehat{X}_{1:T}(X_{0},a_{1:T})$
given $X_{0}$, i.e.

$\mathrm{Law}[\widehat{X}_{1:T}(x_{0},a_{1:T})]=\mathrm{Law}[\widehat{X}_{1:T}(X_{0},a_{1:T})\mid
X_{0}=x_{0}]\qquad\text{for $\mathrm{Law}[X_{0}]$-almost all
$x_{0}\in\mathcal{X}_{0}$.}$

Substituting this into (19), we see that (18) holds if and only if (17) does.
The result now follows since $a_{1:T}$ was arbitrary. ∎

##  Appendix C Online prediction

A distinguishing feature of many digital twins is their ability to integrate
real-time information obtained from sensors in their environment [3]. It is
therefore relevant to consider a setting in which a twin is used repeatedly to
make a sequence of predictions over time, each time taking all previous
information into account. One way to formalise this is to instantiate our
model for the twin at each timestep. For example, we could represent the
predictions made by the twin at $t=0$ after observing initial covariates
$x_{0}$ as potential outcomes $\widehat{X}^{1}_{1:T}(x_{0},a_{1:T})$,
similarly as we did in the main text. We could then represent the predictions
made by the twin after some action $a_{1}$ is taken and an additional
observation $x_{1}$ is made via potential outcomes
$\widehat{X}^{2}_{2:T}((x_{0:1},a_{1}),a_{2:T})$. More generally, for
$t\in\\{1,\ldots,T\\}$, we could introduce potential outcomes
$\widehat{X}^{t}_{t:T}((x_{0:t-1},a_{1:t-1}),a_{t:T})$ to represent the
predictions that the twin would make at time $t$ after the observations
$x_{0:t-1}$ are made and the actions $a_{1:t-1}$ are taken.

###  C.1 Correctness in the online setting

This extended model requires a new definition of correctness than our
Definition 1 from the main text. A natural approach is to say that the twin is
correct in this new setting if

$\mathrm{Law}[\widehat{X}_{t:T}^{t}((x_{0:t-1},a_{1:t-1}),a_{t:T})]=\mathrm{Law}[X_{t:T}(a_{1:T})\mid
X_{0:t-1}(a_{1:t-1})=x_{0:t-1}]\\\ \text{for all $t\in\\{1,\ldots,T\\}$,
$a_{1:T}\in\mathcal{A}_{1:T}$, and almost all
$x_{0:t-1}\in\mathcal{X}_{0:t-1}$}.$ (20)

A twin with this property would at each step be able to accurately simulate
the future in light of previous information, use this to choose a next action
to take, observe the result of doing so, and then repeat.

It is possible to show that (20) holds if and only if:

  1. 1.

$\mathrm{Law}[\widehat{X}^{1}_{1:T}(x_{0},a_{1:T})]=\mathrm{Law}[X_{1:T}(a_{1:T})\mid
X_{0}=x_{0}]$ for all $t\in\\{1,\ldots,T\\}$, $a_{1:t}\in\mathcal{A}_{1:t}$
and almost all $x_{0}\in\mathcal{X}_{0}$; and

  2. 2.

$\mathrm{Law}[\widehat{X}^{t}_{t:T}((x_{0:t-1},a_{1:t-1}),a_{t:T})]=\mathrm{Law}[\widehat{X}^{s}_{t:T}((x_{0:s-1},a_{1:s-1}),a_{s:T})\mid\widehat{X}^{s}_{s:t-1}((x_{0:s-1},a_{1:s-1}),a_{s:t-1})=x_{s:t-1}]$
for all $t\in\\{1,\ldots,T\\}$, $s\in\\{1,\ldots,t-1\\}$,
$a_{1:T}\in\mathcal{A}_{1:T}$, and almost all
$x_{0:t-1}\in\mathcal{X}_{0:t-1}$.

The first condition here says that $\widehat{X}^{1}_{1:T}(x_{0},a_{1:T})$ must
be interventionally correct in the sense of Definition 1 from the main text.
The second condition says that the predictions made by the twin across
different timesteps must be internally consistent with each other insofar as
their conditional distributions must align. This holds automatically in many
circumstances, such as if the predictions of the twin are obtained from a
Bayesian model (for example), and otherwise could be checked numerically given
the ability to run simulations from the twin, without the need to obtain data
or refer to the real-world process in any way. As such, the problem of
assessing the correctness of the twin in this new sense primarily reduces to
the problem of assessing the correctness of
$\widehat{X}^{1}_{1:T}(x_{0},a_{1:T})$ in the sense of Definition 1 in the
main text, which motivates our focus on that condition.

###  C.2 Alternative notions of online correctness

An important and interesting subtlety arises in this context that is worth
noting. In general it does not follow that a twin correct in the sense of (20)
satisfies

$\mathrm{Law}[\widehat{X}_{t:T}^{t}((x_{0:t-1},a_{1:t-1}),a_{t:T})]\
=\mathrm{Law}[X_{t:T}(A_{1:t-1},a_{t:T})\mid
X_{0:t-1}(A_{1:t-1})=x_{0:t-1},A_{1:t-1}=a_{1:t-1}]\\\ \text{for all
$a_{1:T}\in\mathcal{A}_{1:T}$, and almost all
$x_{0:t-1}\in\mathcal{X}_{0:t-1}$},$ (21)

since in general it does not hold that

$\mathrm{Law}[X_{t:T}(a_{1:T})\mid
X_{0:t-1}(a_{1:t-1})=x_{0:t-1}]=\mathrm{Law}[X_{t:T}(A_{1:t-1},a_{t:T})\mid
X_{0:t-1}(A_{1:t-1})=x_{0:t-1},A_{1:t-1}=a_{1:t-1}]\\\ \text{for all
$a_{1:T}\in\mathcal{A}_{1:T}$, and almost all
$x_{0:t-1}\in\mathcal{X}_{0:t-1}$}$

unless the actions $A_{1:t-1}$ are unconfounded. (Here $A_{1:T}$ denotes the
actions of a behavioural agent; see Section 4 of the main text.) In other
words, a twin that is correct in the sense of (20) will make accurate
predictions at time $t$ when every action taken before time $t$ was
unconfounded (as occurs for example when the twin is directly in control of
the decision-making process), but in general not when certain taken actions
before time $t$ were chosen by a behavioural agent with access to more context
than is available to the twin (as may occur for example when the twin is used
as a decision-support tool). However, should it be desirable, our framework
could be extended to encompass the alternative condition in (21) by
relabelling the observed history $(X_{0:t-1}(A_{1:t-1}),A_{1:t-1})$ as
$X_{0}$, and then assessing the correctness of
$\widehat{X}_{t:T}^{t}((x_{0:t-1},a_{1:t-1}),a_{t:T})$ in the sense of
Definition 1 from the main text.

Overall, the “right” notion of correctness in this online setting is to some
extent a design choice. We believe our causal approach to twin assessment
provides a useful framework for formulating and reasoning about these
possibilities, and consider the investigation of assessment strategies for
additional usage regimes to be an interesting direction for future work.

##  Appendix D Proof of Proposition 2 (interventional distributions are not
identifiable)

It is well-known in the causal inference literature that the interventional
behaviour of the real-world process cannot be uniquely determined from
observational data. For completeness, we now provide a self-contained proof of
this result. Our statement here is lengthier than Proposition 2 in the main
text in order to clarify what we mean by “uniquely determined”: intuitively,
the idea is that there always exist distinct families of potential outcomes
whose interventional behaviours differ and yet give rise to the same
observational data. We also account for the corner case that
$\mathcal{X}_{0:T}$ is trivial (e.g. if each $\mathcal{X}_{0:T}$ is a
singleton set), which is of no consequence in practice.

######  Proposition 2.

Suppose we have $t\in\\{1,\ldots,T\\}$ such that the underlying
$\sigma$-algebra of $\mathcal{X}_{t}$ is nontrivial, and that for some
$a_{1:t}\in\mathcal{A}_{1:t}$ it holds that $\mathbb{P}(A_{1:t}\neq
a_{1:t})>0$. Then there exist $\mathcal{X}_{0:T}$-valued potential outcomes
$(\tilde{X}_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$
defined on the same probability space as
$(X_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$ such that

$\tilde{X}_{t}(a_{1:t})\overset{\mathrm{d}}{\neq}X_{t}(a_{1:t})$

but for which

$(\tilde{X}_{0:T}(A_{1:T}),A_{1:T})\overset{\mathrm{a.s.}}{=}(X_{0:T}(A_{1:T}),A_{1:T}).$
(22)

###### Proof.

The underlying idea is to define, for each $s\in\\{0,\ldots,T\\}$ and
$a_{1:s}^{\prime}\in\mathcal{A}_{1:s}$,

$\tilde{X}_{s}(a_{1:s}^{\prime})\coloneqq\begin{cases}\mathbbm{1}(A_{1:t}=a_{1:t})\,X_{t}(a_{1:t})+\mathbbm{1}(A_{1:t}\neq
a_{1:t})\,Z&\text{if $s=t$ and $a_{1:s}^{\prime}=a_{1:t}$}\\\
X_{s}(a_{1:s}^{\prime})&\text{otherwise},\end{cases}$

where $Z$ is some $\mathcal{X}_{t}$-valued random variable with
$Z\perp\\!\\!\\!\perp A_{1:t}$ and with
$\mathrm{Law}[Z]\neq\mathrm{Law}[X_{t}(a_{1:t})\mid A_{1:t}\neq a_{1:t}]$.
(Note that our assumption that $\mathbb{P}(A_{1:t}\neq a_{1:t})>0$ means this
conditional distribution is well-defined.) It is then easily checked that (22)
holds, but

$\displaystyle\mathrm{Law}[\tilde{X}_{t}(a_{1:t})]$
$\displaystyle=\mathrm{Law}[\tilde{X}_{t}(a_{1:t})\mid
A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t})+\mathrm{Law}[\tilde{X}_{t}(a_{1:t})\mid
A_{1:t}\neq a_{1:t}]\,\mathbb{P}(A_{1:t}\neq a_{1:t})$
$\displaystyle=\mathrm{Law}[X_{t}(a_{1:t})\mid
A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t})+\mathrm{Law}[Z]\,\mathbb{P}(A_{1:t}\neq
a_{1:t})$ $\displaystyle\neq\mathrm{Law}[X_{t}(a_{1:t})\mid
A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t})+\mathrm{Law}[X_{t}(a_{1:t})\mid
A_{1:t}\neq a_{1:t}]\,\mathbb{P}(A_{1:t}\neq a_{1:t})$
$\displaystyle=\mathrm{Law}[X_{t}(a_{1:t})],$

from which the result follows.

To obtain $Z$, we use our assumption that the underlying $\sigma$-algebra of
$\mathcal{X}_{t}$ is nontrivial. This means there exists a measurable subset
$B\subseteq\mathcal{X}_{t}$ with $B\neq\emptyset$ and $B\neq\mathcal{X}_{t}$.
We may therefore choose $x\in B$ and $x^{\prime}\not\in B$, and define

$Z\coloneqq\begin{cases}x&\text{if $\mathbb{P}(X_{t}(a_{1:t})\in B\mid
A_{1:t}\neq a_{1:t})=0$}\\\ x^{\prime}&\text{otherwise}.\end{cases}$

In either case, certainly $Z\perp\\!\\!\\!\perp A_{1:t}$ since $Z$ is
constant, and likewise

$\mathbb{P}(Z\in B)\neq\mathbb{P}(X_{t}(a_{1:t})\in B\mid A_{1:t}\neq
a_{1:t}),$

which means $\mathrm{Law}[Z]\neq\mathrm{Law}[X_{t}(a_{1:t})\mid A_{1:t}\neq
a_{1:t}]$. ∎

##  Appendix E Deterministic potential outcomes are unconfounded

In this section we expand on our earlier claim that, if the real-world process
is deterministic, then the observational data is unconfounded. We first make
this claim precise. By “deterministic”, we mean that there exist measurable
functions $g_{t}$ for $t\in\\{1,\ldots,T\\}$ such that

$X_{t}(a_{1:t})\overset{\mathrm{a.s.}}{=}g_{t}(X_{0:t-1}(a_{1:t-1}),a_{1:t})\qquad\text{for
all $t\in\\{1,\ldots,T\\}$ and $a_{1:t}\in\mathcal{A}_{1:t}$.}$ (23)

By “unconfounded”, we mean that the _sequential randomisation assumption
(SRA)_ introduced by Robins [59] holds, i.e.

$(X_{s}(a_{1:s}):s\in\\{1,\ldots,T\\},a_{1:s}\in\mathcal{A}_{1:s})\perp\\!\\!\\!\perp
A_{t}\mid X_{0:t-1}(A_{1:t-1}),A_{1:t-1}\qquad\text{for all
$t\in\\{1,\ldots,T\\}$},$ (24)

where $\perp\\!\\!\\!\perp$ denotes conditional independence. Intuitively,
this says that, apart from the historical observations
$(X_{0:t-1}(A_{1:t-1}),A_{1:t-1})$, any additional factors that influence the
agent’s choice of action $A_{t}$ are independent of the behaviour of the real-
world process. The SRA provides a standard formulation of the notion of
unconfoundedness in longitudinal settings such as ours (see [72, Chapter 5]
for a review).

It is now a standard exercise to show that (23) implies (24). We include a
proof below for completeness. Key to this is the following straightforward
Lemma.

######  Lemma A.

Suppose $U$ and $V$ are random variables such that, for some measurable
function $g$, it holds that $U\overset{\mathrm{a.s.}}{=}g(V)$. Then, for any
other random variable $W$, we have

$U\perp\\!\\!\\!\perp W\mid V.$

###### Proof.

By standard properties of conditional expectations, for any measurable sets
$S_{1}$ and $S_{2}$, we have almost surely

$\displaystyle\mathbb{P}(U\in S_{1},W\in S_{2}\mid V)$
$\displaystyle=\mathbb{E}[\mathbbm{1}(g(V)\in S_{1})\,\mathbbm{1}(W\in
S_{2})\mid V]$ $\displaystyle=\mathbbm{1}(g(V)\in
S_{1})\,\mathbb{E}[\mathbbm{1}(W\in S_{2})\mid V]$
$\displaystyle=\mathbb{E}[\mathbbm{1}(U\in S_{1})\mid V]\,\mathbb{P}(W\in
S_{2}\mid V)$ $\displaystyle=\mathbb{P}(U\in S_{1}\mid V)\,\mathbb{P}(W\in
S_{2}\mid V),$

which gives the result. ∎

It is now easy to see that (23) implies (24). Indeed, by recursive
substitution, it is straightforward to show that there exist measurable
functions $h_{t}$ for $t\in\\{1,\ldots,T\\}$ such that

$X_{t}(a_{1:t})\overset{\mathrm{a.s.}}{=}h_{t}(X_{0},a_{1:t})\qquad\text{for
all $t\in\\{1,\ldots,T\\}$ and $a_{1:t}\in\mathcal{A}_{1:t}$},$

and so

$(X_{s}(a_{1:s}):s\in\\{1,\ldots,T\\},a_{1:s}\in\mathcal{A}_{1:s})=(h_{t}(X_{0},a_{1:s}):s\in\\{1,\ldots,T\\},a_{1:s}\in\mathcal{A}_{1:s}).$

The right-hand side is now seen to be a measurable function of $X_{0}$ and
hence certainly of $(X_{0:t-1}(A_{1:t-1}),A_{1:t-1})$, so that the result
follows by Lemma A.

##  Appendix F Causal bounds

###  F.1 Preliminaries

We now provide a self-contained proof of the causal bounds provided in Theorem
1 in the main text, as well as additional results demonstrating the optimality
of these bounds. Throughout this section, we slightly generalise our model by
relaxing our assumption from the main text that

$Y(a_{1:t}^{\prime})\overset{\mathrm{a.s.}}{=}f(X_{0:t}(a_{1:t}^{\prime})).$
(25)

Instead, we allow $(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$
to be arbitrary real-valued random variables defined jointly on the same
probability space as
$(X_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$ and
$A_{1:T}$, a setup that encompasses (25) as a special case. At various points,
we also discuss considerations that arise when (25) is known to hold.

###  F.2 Proof of Theorem 1

We now give a proof of Theorem 1 from the main text.

######  Theorem 1.

Given some choice of $t\in\\{1,\ldots,T\\}$ and $a_{1:t}\in\mathcal{A}_{1:t}$,
suppose we have measurable $B_{0:t}\subseteq\mathcal{X}_{0:t}$ such that

$\mathbb{P}(X_{0:t}(a_{1:t})\in B_{0:t})>0,$ (26)

and $y_{\mathrm{lo}},y_{\mathrm{up}}\in\mathbb{R}$ such that

$\mathbb{P}(y_{\mathrm{lo}}\leq Y(a_{1:t})\leq y_{\mathrm{up}}\mid
X_{0:t}(a_{1:t})\in B_{0:t})=1.$ (27)

Then, defining

$\displaystyle N$ $\displaystyle\coloneqq\max\\{0\leq s\leq t\mid
A_{1:s}=a_{1:s}\\}$ $\displaystyle Y_{\mathrm{lo}}$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:t}=a_{1:t})\,Y(A_{1:t})+\mathbbm{1}(A_{1:t}\neq
a_{1:t})\,y_{\mathrm{lo}}$ $\displaystyle Y_{\mathrm{up}}$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:t}=a_{1:t})\,Y(A_{1:t})+\mathbbm{1}(A_{1:t}\neq
a_{1:t})\,y_{\mathrm{up}},$

it follows that

$\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0:N}(A_{1:N})\in
B_{0:N}]\leq\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in
B_{0:t}]\leq\mathbb{E}[Y_{\mathrm{up}}\mid X_{0:N}(A_{1:N})\in B_{0:N}].$

###### Proof.

It is easily checked that

$\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in
B_{0:t}]=\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in
B_{0:t},A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:t}(a_{1:t})\in
B_{0:t})\\\ +\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in B_{0:t},A_{1:t}\neq
a_{1:t}]\,\mathbb{P}(A_{1:t}\neq a_{1:t}\mid X_{0:t}(a_{1:t})\in B_{0:t}).$
(28)

If $\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:t}(a_{1:t})\in B_{0:t})>0$, then

$\displaystyle\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in
B_{0:t},A_{1:t}=a_{1:t}]$ $\displaystyle=\mathbb{E}[Y(A_{1:t})\mid
X_{0:t}(A_{1:t})\in B_{0:t},A_{1:t}=a_{1:t}]$
$\displaystyle=\mathbb{E}[Y(A_{1:N})\mid X_{0:N}(A_{1:N})\in
B_{0:N},A_{1:t}=a_{1:t}],$

where the second step follows because $\mathbb{P}(N=t\mid A_{1:t}=a_{1:t})=1$.
Similarly, if $\mathbb{P}(A_{1:t}\neq a_{1:t}\mid X_{0:t}(a_{1:t})\in
B_{0:t})>0$, then (27) implies

$\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in B_{0:t},A_{1:t}\neq
a_{1:t}]\geq y_{\mathrm{lo}}.$

Substituting these results into (28), we obtain

$\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in
B_{0:t}]\geq\mathbb{E}[Y(A_{1:t})\mid X_{0:N}(A_{1:N})\in
B_{0:N},A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:t}(a_{1:t})\in
B_{0:t})\\\ +y_{\mathrm{lo}}\,\mathbb{P}(A_{1:t}\neq a_{1:t}\mid
X_{0:t}(a_{1:t})\in B_{0:t}).$ (29)

Now observe that the right-hand side of (29) is a convex combination with
mixture weights $\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:t}(a_{1:t})\in B_{0:t})$
and $\mathbb{P}(A_{1:t}\neq a_{1:t}\mid X_{0:t}(a_{1:t})\in B_{0:t})$. We can
bound

$\displaystyle\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:t}(a_{1:t})\in B_{0:t})$
$\displaystyle=\frac{\mathbb{P}(X_{0:t}(a_{1:t})\in
B_{0:t},A_{1:t}=a_{1:t})}{\mathbb{P}(X_{0:t}(a_{1:t})\in B_{0:t})}$
$\displaystyle\geq\frac{\mathbb{P}(X_{0:t}(a_{1:t})\in
B_{0:t},A_{1:t}=a_{1:t})}{\mathbb{P}(X_{0:N}(a_{1:N})\in B_{0:N})}$
$\displaystyle=\frac{\mathbb{P}(X_{0:N}(A_{1:N})\in
B_{0:N},A_{1:t}=a_{1:t})}{\mathbb{P}(X_{0:N}(A_{1:N})\in B_{0:N})}$
$\displaystyle=\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:N}(A_{1:N})\in B_{0:N}),$
(30)

where the inequality holds because $t\geq N$ almost surely, and the second
equality holds because the definition of $N$ means

$X_{0:N}(a_{1:N})\overset{\mathrm{a.s.}}{=}X_{0:N}(A_{1:N}).$

As such, we can bound the convex combination in (29) from below by replacing
its mixture weights with $\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:N}(A_{1:N})\in
B_{0:N})$ and $\mathbb{P}(A_{1:t}\neq a_{1:t}\mid X_{0:N}(a_{1:N})\in
B_{0:N})$, which shifts weight from the $\mathbb{E}[Y(A_{1:t})\mid
A_{1:t}=a_{1:t},X_{0:N}(A_{1:N})\in B_{0:N}]$ term onto the $y_{\mathrm{lo}}$
term. This yields

$\displaystyle\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in B_{0:t}]$
$\displaystyle\qquad\qquad\geq\mathbb{E}[Y(A_{1:t})\mid X_{0:N}(A_{1:N})\in
B_{0:N},A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t}\mid X_{0:N}(a_{1:N})\in
B_{0:N})$
$\displaystyle\qquad\qquad\qquad+y_{\mathrm{lo}}\,\mathbb{P}(A_{1:t}\neq
a_{1:t}\mid X_{0:N}(a_{1:N})\in B_{0:N})$
$\displaystyle\qquad\qquad=\mathbb{E}[Y(A_{1:t})\,\mathbbm{1}(A_{1:t}=a_{1:t})+y_{\mathrm{lo}}\,\mathbbm{1}(A_{1:t}\neq
a_{1:t})\mid X_{0:N}(A_{1:N})\in B_{0:N}]$
$\displaystyle\qquad\qquad=\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0:N}(A_{1:N})\in
B_{0:N}].$

∎

###  F.3 Sharpness of bounds

The following result shows that Theorem 1 is sharp in the sense that it is not
possible to obtain tighter bounds on $\mathbb{E}[Y(a_{1:t})\mid
X_{0:t}(a_{1:t})\in B_{0:t}]$ in terms of observational quantities. Roughly
speaking, for any model consisting of
$(X_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$,
$(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$, and $A_{1:T}$,
it is always possible to define alternative models whose observational
behaviour is equivalent but for which the bounds in Theorem 1 are achieved.

######  Proposition C.

Assume the same setup as in Theorem 1. There always exists potential outcomes
$(\tilde{X}_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$ and
$(\tilde{Y}(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$ defined
on the same probability space as
$(X_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$,
$(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$, and $A_{1:T}$
that satisfy

$(\tilde{X}_{0:T}(A_{1:T}),\tilde{Y}(A_{1:t}),A_{1:T})\overset{\mathrm{a.s.}}{=}(\tilde{X}_{0:T}(A_{1:T}),Y(A_{1:t})A_{1:T}),$

but for which

$\mathbb{E}[\tilde{Y}(a_{1:t})\mid\tilde{X}_{0:t}(a_{1:t})\in
B_{0:t}]=\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0:N}(A_{1:N})\in B_{0:N}].$

A similar statement is true for the upper bound
$\mathbb{E}[Y_{\mathrm{up}}\mid X_{0:N}(A_{1:N})\in B_{0:N}]$.

###### Proof.

We consider the case of the lower bound; the case of the upper bound is
analogous. Choose $x_{1:T}\in B_{1:T}$ arbitrarily. (Certainly some choice is
always possible, since each $B_{s}$ has positive measure and is therefore
nonempty.) Define

$\displaystyle\tilde{X}_{0}$ $\displaystyle\coloneqq X_{0}$
$\displaystyle\tilde{X}_{s}(a_{1:s}^{\prime})$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:s}=a_{1:s}^{\prime})\,X_{s}(a_{1:s}^{\prime})+\mathbbm{1}(A_{1:s}\neq
a_{1:s}^{\prime})\,x_{s}\qquad\text{for each $s\in\\{0,\ldots,T\\}$ and
$a_{1:s}^{\prime}\in\mathcal{A}_{1:s}$},$

and similarly let

$\tilde{Y}(a_{1:t}^{\prime})=\mathbbm{1}(A_{1:t}=a_{1:t}^{\prime})\,Y(a_{1:t}^{\prime})+\mathbbm{1}(A_{1:t}\neq
a_{1:t}^{\prime})\,y_{\mathrm{lo}}\qquad\text{for all
$a_{1:t}^{\prime}\in\mathcal{A}_{1:t}$}.$

It is easy to check that

$(\tilde{X}_{0:T}(A_{1:T}),\tilde{Y}(A_{1:t}),A_{1:T})\overset{\mathrm{a.s.}}{=}(X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}).$

But now, we have directly

$\tilde{Y}(a_{1:t})=Y_{\mathrm{lo}}.$

Moreover, it is easily checked from the definition of $N$ and
$\tilde{X}_{0:t}(a_{1:t})$ that

$\tilde{X}_{0:t}(a_{1:t})\overset{\mathrm{a.s.}}{=}(X_{0:N}(A_{1:N}),x_{N+1:t}),$

so that

$\displaystyle\mathbbm{1}(\tilde{X}_{0:t}(a_{1:t})\in B_{0:t})$
$\displaystyle\overset{\mathrm{a.s.}}{=}\mathbbm{1}(\tilde{X}_{0:N}(a_{1:N})\in
B_{0:N},x_{N+1:t}\in B_{N+1:t})$
$\displaystyle\overset{\mathrm{a.s.}}{=}\mathbbm{1}(X_{0:N}(A_{1:N})\in
B_{0:N})$

since each $x_{s}\in B_{s}$. Consequently,

$\displaystyle\mathbb{E}[\tilde{Y}(a_{1:t})\mid\tilde{X}_{0:t}(a_{1:t})\in
B_{0:t}]$
$\displaystyle=\mathbb{E}[Y_{\mathrm{lo}}\mid\tilde{X}_{0:t}(a_{1:t})\in
B_{0:t}]$ $\displaystyle=\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0:N}(A_{1:N})\in
B_{0:N}],$

which gives the result. ∎

####  F.3.1 Our case

The previous discussion considered
$(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$ to be arbitrary
potential outcomes defined jointly on the same probability space as
$(X_{0:t}(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$. In
contrast, our methodology assumes the particular form

$Y(a_{1:t}^{\prime})\overset{\mathrm{a.s.}}{=}f(X_{0:t}(a_{1:t}^{\prime})),$

for some known measurable function $f:\mathcal{X}_{0:t}\to\mathbb{R}$.
Proposition C does not directly imply that the bounds in Theorem 1 are sharp
under this additional assumption; however, this remains true for many cases of
interest in practice. In particular, if $f$ depends only on the final
observation space $\mathcal{X}_{t}$ (which is true for example throughout our
case study), i.e. we have

$Y(a_{1:t}^{\prime})\overset{\mathrm{a.s.}}{=}f(X_{t}(a_{1:t}^{\prime})),$

then it still holds that the bounds can be achieved, provided the worst-case
values are chosen sensibly as

$y_{\mathrm{lo}}=\min_{x_{t}\in B_{t}}f(x_{t})\qquad
y_{\mathrm{up}}=\max_{x_{t}\in B_{t}}f(x_{t}).$

This follows straightforwardly by modifying the proof of Proposition C to
define $x_{t}$ as either the minimiser or maximiser of $f$ on $B_{t}$.

###  F.4 Bounds on the conditional expectation given specific covariate values

Theorem 1 provides a bound on $\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})\in
B_{0:t}]$, i.e. the conditional expectation given the event
$\\{X_{0:t}(a_{1:t})\in B_{0:t}\\}$, which is assumed to have positive
probability. We consider here the prospect of obtaining bounds on
$\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})]$, i.e. the conditional
expectation given the value of $X_{0:t}(a_{1:t})$. For falsification purposes,
this would provide a means for using a sufficiently large observational
dataset to determine that the behaviour of the twin is incorrect when it
outputs specific values of $\widehat{X}_{0:t}(a_{1:t})$, rather than just that
it is incorrect on average across all runs that output values
$\widehat{X}_{0:t}(a_{1:t})\in B_{0:t}$.

#### Discrete covariates

When $X_{0:t}(a_{1:t})$ is discrete, Theorem 1 directly yields a bound on
$\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})]$, i.e. a conditional expectation
given the value of $X_{0:t}(a_{1:t})$, which can be seen by taking
$B_{0:t}=\\{x_{0:t}\\}$ to be any singleton set with
$\mathbb{P}(X_{0:t}(a_{1:t})=x_{0:t})>0$. Since the constants
$y_{\mathrm{lo}},y_{\mathrm{up}}\in\mathbb{R}$ in Theorem 1 were allowed to
depend on $B_{0:t}$, and hence here on each choice of $x_{0:t}$, it is natural
to think of these now as measurable functions
$y_{\mathrm{lo}},y_{\mathrm{up}}:\mathcal{X}_{0:t}\to\mathbb{R}$ that are
required to satisfy

$y_{\mathrm{lo}}(X_{0:t}(a_{1:t}))\leq Y(a_{1:t})\leq
y_{\mathrm{up}}(X_{0:t}(a_{1:t}))\qquad\text{almost surely}.$ (31)

In other words, when $X_{0:t}(a_{1:t})$ is discrete, Theorem 1 may be thought
of as providing functions $\Psi_{\mathrm{lo}}$ and $\Psi_{\mathrm{up}}$ that
map observational distributions to measurable functions
$\mathcal{X}_{0:t}\to\mathbb{R}$ such that, for all models consisting of
$(X_{0:t}(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$,
$(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$, and $A_{1:t}$
satisfying (31), we have almost surely

$\Psi_{\mathrm{lo}}(\mathrm{Law}[X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}])(X_{0:t}(a_{1:t}))\leq\mathbb{E}[Y(a_{1:t})\mid
X_{0:t}(a_{1:t})]\\\
\leq\Psi_{\mathrm{up}}(\mathrm{Law}[X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}])(X_{0:t}(a_{1:t})).$
(32)

####  F.4.1 Continuous initial covariates

When $\mathbb{P}(X_{1:t}(a_{1:t})\in B_{1:t})>0$, a fairly straightforward
modification of the proof of Theorem 1 yields bounds of the following form:

$\displaystyle\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0},X_{1:N}(A_{1:N})\in
B_{1:N}]$ $\displaystyle\leq\mathbb{E}[Y(a_{1:t})\mid
X_{0},X_{1:t}(a_{1:t})\in B_{1:t}]$
$\displaystyle\qquad\qquad\leq\mathbb{E}[Y_{\mathrm{up}}\mid
X_{0},X_{1:N}(A_{1:N})\in B_{1:N}]\qquad\text{almost surely}.$ (33)

In particular, this holds regardless of whether or not $X_{0}$ is discrete. In
turn, if $X_{1:t}(a_{1:t})$ is discrete, then by a similar argument as was
given in the previous subsection, this yields almost sure bounds on
$\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})]$ of the form in (32), provided
(31) holds. Alternatively, by taking $B_{1:t}\coloneqq\mathcal{X}_{1:t}$, (33)
yields bounds of the form

$\mathbb{E}[Y_{\mathrm{lo}}\mid X_{0}]\leq\mathbb{E}[Y(a_{1:t})\mid
X_{0}]\leq\mathbb{E}[Y_{\mathrm{up}}\mid X_{0}].$

If the action sequence $a_{1:t}$ is thought of as a single choice of an action
from the extended action space $\mathcal{A}_{1:t}$, then this recovers the
bounds originally proposed by Manski [42], which allowed conditioning on
potentially continuous pre-treatment covariates corresponding to our $X_{0}$.

####  F.4.2 Continuous subsequent covariates

When the distribution of $X_{1:t}(a_{1:t})$ is continuous, Theorem 1 does not
directly give rise to bounds of the form in (32) via the argument made for
discrete $X_{1:t}(a_{1:t})$, since then
$\mathbb{P}(X_{1:t}(a_{1:t})=x_{1:t})=0$ for all
$x_{1:t}\in\mathcal{X}_{1:t}$. As such, it is natural to ask whether such
bounds can be obtained through a more general strategy. Unfortunately, and
somewhat surprisingly, it turns out that no nontrivial bounds of the form in
(32) exist apart from certain corner cases, such as if
$\mathbb{P}(A_{1:t}=a_{1:t})=1$ (in which case $\mathbb{E}[Y(a_{1:t})\mid
X_{0:t}(a_{1:t})]\overset{\mathrm{a.s.}}{=}\mathbb{E}[Y(A_{1:t})\mid
X_{0:t}(A_{1:t})]$ is always identifiable). Precisely, we have the following
result:

######  Theorem B.

Given $t\in\\{1,\ldots,T\\}$, $a_{1:t}\in\mathcal{A}_{1:t}$, and measurable
functions $y_{\mathrm{lo}},y_{\mathrm{up}}:\mathcal{X}_{0:t}\to\mathbb{R}$,
suppose we have $\Psi_{\mathrm{lo}}$ and $\Psi_{\mathrm{up}}$ that satisfy
(32) for all models consisting of
$(X_{0:T}(a_{1:T}^{\prime}):a_{1:T}^{\prime}\in\mathcal{A}_{1:T})$,
$(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$, and $A_{1:T}$
that satisfy (31). Given a model satisfying (31), if for some
$s\in\\{1,\ldots,t\\}$ we have

$\displaystyle\mathbb{P}(X_{s}(A_{1:s})=x_{s})$
$\displaystyle=0\qquad\text{for all $x_{s}\in\mathcal{X}_{s}$},$ (34)
$\displaystyle\mathbb{P}(A_{1:s}\neq a_{1:s}\mid X_{0:s-1}(A_{1:s-1}))$
$\displaystyle>0\qquad\text{almost surely},$ (35)

then the following must also hold:

$\displaystyle\Psi_{\mathrm{lo}}(\mathrm{Law}[X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}])(X_{0:t}(a_{1:t}))$
$\displaystyle\leq y_{\mathrm{lo}}(X_{0:t}(a_{1:t}))\qquad\text{almost
surely}$ (36)
$\displaystyle\Psi_{\mathrm{up}}(\mathrm{Law}[X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}])(X_{0:t}(a_{1:t}))$
$\displaystyle\geq y_{\mathrm{up}}(X_{0:t}(a_{1:t}))\qquad\text{almost
surely}.$ (37)

A proof is given further below. The condition in (34) holds if
$X_{s}(A_{1:s})$ is continuous, or more generally if it has some continuous
component. The condition in (35) constitutes a mild positivity assumption that
in practice can be expected always to hold. In particular,
$\mathbb{P}(A_{1:t}\neq a_{1:t}\mid X_{0:t-1}(A_{1:t-1}))$ may become
arbitrarily small with arbitrarily high probability, provided that it remains
nonzero. At a high level, this result says that, under these conditions, the
observational distribution cannot be used to obtain any new information about
the behaviour of $\mathbb{E}[Y(a_{1:t})\mid X_{0:t}(a_{1:t})]$ than is already
contained in our assumption in (31), which always trivially implies

$y_{\mathrm{lo}}(X_{0:t}(a_{1:t}))\leq\mathbb{E}[Y(a_{1:t})\mid
X_{0:t}(a_{1:t})]\leq y_{\mathrm{up}}(X_{0:t}(a_{1:t}))\qquad\text{almost
surely}$

regardless of the observational distribution.

With minor modification to the proof, (35) could be replaced by the (in
general disjoint) condition that

$\mathbb{P}(A_{1:s}\neq a_{1:s}\mid X_{0:s-1}(a_{1:s-1}))>0\qquad\text{almost
surely}.$ (38)

However, this alternative condition is somewhat unsatisfying since it is in
general unidentifiable, and so there exist models for which it holds but for
which this could not be determined solely from the observational distribution.
More generally, and possibly at the expense of a more complicated argument, we
conjecture that Theorem B remains true if

$\mathbb{P}(A_{1:s}\neq a_{1:s}\mid X_{0})>0\qquad\text{almost surely},$

which is weaker than both (35) and (38).

We emphasise that Theorem B does not hold in general if (34) is only true for
$s=0$, which is to be expected in light of (33). Intuitively, the phenomenon
underlying Theorem B arises when conditioning on (in general) unobserved
quantities that are not discrete. As such, roughly speaking, the phenomenon
does not arise when only $X_{0}$ is continuous, since this quantity is always
observed. In the next section, we illustrate this point via a toy example.

##### Illustrative example

To gain intuition for the phenomenon underlying Theorem B, consider a
simplified model consisting of $\mathcal{X}$-valued potential outcomes
$(X(a^{\prime}):a\in\mathcal{A})$, $\mathbb{R}$-valued potential outcomes
$(Y(a^{\prime}):a\in\mathcal{A})$, and an $\mathcal{A}$-valued random variable
$A$ representing the choice of action. (This constitutes a special case of our
setup with $T=1$ and $\mathcal{X}_{0}$ taken to be a singleton set.) Suppose
moreover that the following conditions hold:

$\displaystyle\mathbb{P}(X(A)=x)$ $\displaystyle=0\qquad\text{for all
$x\in\mathcal{X}$}$ $\displaystyle\mathbb{P}(A=a)$ $\displaystyle<1.$

We then have

$\mathbb{E}[Y(a)\mid X(a)]\overset{\mathrm{a.s.}}{=}\mathbb{E}[Y(A)\mid
X(A),A=a]\,\mathbb{P}(A=a\mid X(a))+\mathbb{E}[Y(a)\mid X(a),A\neq
a]\,\mathbb{P}(A\neq a\mid X(a)).$ (39)

But now, since the behaviour of $X(a)$ is only observed on $\\{A=a\\}$, for
any given value of $x\in\mathcal{X}$, we cannot rule out the possibility that

$X(a)=\mathbbm{1}(A=a)\,X(A)+\mathbbm{1}(A\neq a)\,x\qquad\text{almost
surely}.$

In turn, since $\mathbb{P}(A=a)>0$, this would imply $\mathbb{P}(X(a)=x)>0$,
and, since $\mathbb{P}(X(A)=x)=0$, that $\mathbb{P}(A=a\mid X(a)=x)=0$. From
(39), this would yield

$\mathbb{E}[Y(a)\mid X(a)=x]=\mathbb{E}[Y(a)\mid X(a)=x,A\neq a].$

But now, since the behaviour of $Y(a)$ is unobserved on $\\{A\neq a\\}$,
intuitively speaking, the observational distribution does not provide any
information about the value of the right-hand side, and therefore about the
behaviour of $\mathbb{E}[Y(a)\mid X(a)]$ more generally since
$x\in\mathcal{X}$ was arbitrary.

##### Our case

The previous discussion considered
$(Y(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$ to be arbitrary
potential outcomes defined jointly on the same probability space as
$(X_{0:t}(a_{1:t}^{\prime}):a_{1:t}^{\prime}\in\mathcal{A}_{1:t})$. In
contrast, our methodology assumes the particular form

$Y(a_{1:t})\overset{\mathrm{a.s.}}{=}f(X_{0:t}(a_{1:t})),$

which means $\mathbb{E}[Y(a_{1:t})\mid
X_{0:t}(a_{1:t})]\overset{\mathrm{a.s.}}{=}f(X_{0:t}(a_{1:t}))$ is known
trivially. In this context, an alternative quantity to consider is
$\mathbb{E}[Y(a_{1:t})\mid X_{0:r}(a_{1:r})]$ with $r\in\\{0,\ldots,t-1\\}$,
which in general will be unknown and therefore still interesting to bound. In
the discrete case, Theorem 1 yields a bound on this quantity obtained by
taking $B_{0:r}=\\{x_{0:r}\\}$ with $\mathbb{P}(X_{0:r}(a_{1:r})=x_{0:r})>0$
and then $B_{r+1:t}\coloneqq\mathcal{X}_{r+1:t}$, and with $y_{\mathrm{lo}}$
and $y_{\mathrm{up}}$ in (31) now replaced with

$\min_{x_{t}\in\mathcal{X}_{t}}f(x_{t})\qquad\text{and}\qquad\max_{x_{t}\in\mathcal{X}_{t}}f(x_{t})$

respectively. However, in the continuous case, the same issues described above
continue to apply in many cases of interest. For example, when $f$ is a
function of $\mathcal{X}_{t}$ only (as is the case for example throughout our
case study), then if (34) and (35) now hold for some $s\in\\{1,\ldots,r\\}$,
without further assumptions, the most informative almost sure lower bound on
$\mathbb{E}[Y(a_{1:t})\mid X_{0:s}(a_{1:s})]$ is

$\min_{x_{t}\in\mathcal{X}_{t}}f(x_{t}),$

which is already known trivially. Roughly, this follows by modifying the proof
of Theorem B so that $\tilde{X}_{t}(a_{1:t}^{\prime})$ becomes

$\tilde{X}_{t}(a_{1:t}^{\prime})\coloneqq\mathbbm{1}(A_{1:r}=a_{1:r}^{\prime})\,X_{r}(a_{1:r}^{\prime})+\mathbbm{1}(A_{1:r}\neq
a_{1:r}^{\prime})\,\operatorname*{arg\,min}_{x_{t}\in\mathcal{X}_{t}}f(x_{t}),$

and $\tilde{Y}(a_{1:t}^{\prime})$ becomes

$\tilde{Y}(a_{1:t}^{\prime})\coloneqq f(\tilde{X}_{t}(a_{1:t})).$

The remainder of the argument is then unchanged. An analogous result holds for
the upper bound also.

##### Proof of Theorem B

###### Proof.

We will show (36); an analogous argument can be used to show (37). To this
end, choose $x_{1:T}\in\mathcal{X}_{1:T}$ arbitrarily, and define new
potential outcomes

$\displaystyle\tilde{X}_{0}$ $\displaystyle\coloneqq X_{0}$
$\displaystyle\tilde{X}_{r}(a_{1:r}^{\prime})$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:r}=a_{1:r}^{\prime})\,X_{r}(a_{1:r}^{\prime})+\mathbbm{1}(A_{1:r}\neq
a_{1:r}^{\prime})\,x_{r}\qquad\text{for $r\in\\{1,\ldots,T\\}$ and
$a_{1:r}^{\prime}\in\mathcal{A}_{1:r}$}.$

Similarly, define

$\displaystyle\tilde{Y}(a_{1:t}^{\prime})$
$\displaystyle\coloneqq\mathbbm{1}(A_{1:t}=a_{1:t}^{\prime})\,Y(a_{1:t}^{\prime})+\mathbbm{1}(A_{1:t}\neq
a_{1:t}^{\prime})\,y_{\mathrm{lo}}(\tilde{X}_{0:t}(a_{1:t}^{\prime}))\qquad\text{for
all $a_{1:t}^{\prime}\in\mathcal{A}_{1:t}$}.$

It follows that

$(\tilde{X}_{0:T}(A_{1:T}),\tilde{Y}(A_{1:t}),A_{1:T})\overset{\mathrm{a.s.}}{=}(X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}),$

and as such

$g_{\mathrm{lo}}\coloneqq\Psi_{\mathrm{lo}}(\mathrm{Law}[X_{0:T}(A_{1:T}),Y(A_{1:t}),A_{1:T}])=\Psi_{\mathrm{lo}}(\mathrm{Law}[\tilde{X}_{0:T}(A_{1:T}),\tilde{Y}(A_{1:t}),A_{1:T}]).$

Moreover, it is easily checked that

$y_{\mathrm{lo}}(\tilde{X}_{0:t}(a_{1:t}))\leq\tilde{Y}(a_{1:t})\leq
y_{\mathrm{up}}(\tilde{X}_{0:t}(a_{1:t}))\qquad\text{almost surely},$

and as such, from our assumed property of $\Psi_{\mathrm{lo}}$, we have almost
surely

$\displaystyle g_{\mathrm{lo}}(\tilde{X}_{0:t}(a_{1:t}))$
$\displaystyle\leq\mathbb{E}[\tilde{Y}(a_{1:t})\mid\tilde{X}_{0:t}(a_{1:t})]$
$\displaystyle=\begin{multlined}\mathbb{E}[\tilde{Y}(A_{1:t})\mid\tilde{X}_{0:t}(A_{1:t}),A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t}\mid\tilde{X}_{0:t}(a_{1:t}))\\\
+\underbrace{\mathbb{E}[\tilde{Y}(a_{1:t})\mid\tilde{X}_{0:t}(a_{1:t}),A_{1:t}\neq
a_{1:t}]}_{=y_{\mathrm{lo}}(\tilde{X}_{0:t}(a_{1:t}))}\,\mathbb{P}(A_{1:t}\neq
a_{1:t}\mid\tilde{X}_{0:t}(a_{1:t})).\end{multlined}\mathbb{E}[\tilde{Y}(A_{1:t})\mid\tilde{X}_{0:t}(A_{1:t}),A_{1:t}=a_{1:t}]\,\mathbb{P}(A_{1:t}=a_{1:t}\mid\tilde{X}_{0:t}(a_{1:t}))\\\
+\underbrace{\mathbb{E}[\tilde{Y}(a_{1:t})\mid\tilde{X}_{0:t}(a_{1:t}),A_{1:t}\neq
a_{1:t}]}_{=y_{\mathrm{lo}}(\tilde{X}_{0:t}(a_{1:t}))}\,\mathbb{P}(A_{1:t}\neq
a_{1:t}\mid\tilde{X}_{0:t}(a_{1:t})).$ (42)

Now, almost surely,

$\displaystyle\mathbbm{1}(A_{1:s}\neq
a_{1:s})\,\mathbb{P}(A_{1:t}=a_{1:t}\mid\tilde{X}_{0:t}(a_{1:t}))$
$\displaystyle=\mathbbm{1}(A_{1:s}\neq
a_{1:s},\tilde{X}_{s}(a_{1:s})=x_{s})\,\mathbb{P}(A_{1:t}=a_{1:t}\mid\tilde{X}_{0:t}(a_{1:t}))$
$\displaystyle=\mathbbm{1}(A_{1:s}\neq
a_{1:s})\,\mathbb{E}[\mathbbm{1}(A_{1:t}=a_{1:t},\tilde{X}_{s}(a_{1:s})=x_{s})\mid\tilde{X}_{0:t}(a_{1:t})]$
$\displaystyle=\mathbbm{1}(A_{1:s}\neq
a_{1:s})\,\mathbb{E}[\mathbbm{1}(A_{1:t}=a_{1:t},X_{s}(A_{1:s})=x_{s})\mid\tilde{X}_{0:t}(a_{1:t})]$
$\displaystyle=0,$

where the last step follows by our assumption in (34). Combining this with
(42), we get

$\mathbbm{1}(A_{1:s}\neq
a_{1:s})\,g_{\mathrm{lo}}(\tilde{X}_{0:s-1}(a_{1:s-1}),x_{s:t})\,\leq\mathbbm{1}(A_{1:s}\neq
a_{1:s})\,y_{\mathrm{lo}}(\tilde{X}_{0:s-1}(a_{1:s-1}),x_{s:t}).$ (43)

But now

$\mathbbm{1}(A_{1:s}\neq
a_{1:s})=\sum_{r=1}^{s}\mathbbm{1}(A_{1:r-1}=a_{1:r-1},A_{r}\neq a_{r}),$

and moreover

$\displaystyle\mathbbm{1}(A_{1:r-1}=a_{1:r-1},A_{r}\neq a_{r})$
$\displaystyle\qquad\overset{\mathrm{a.s.}}{=}\mathbbm{1}(A_{1:r-1}=a_{1:r-1},A_{r}\neq
a_{r},\tilde{X}_{0:r-1}(a_{1:r-1})=X_{0:r-1}(a_{1:r-1}),\tilde{X}_{r:t}(a_{1:t})=x_{r:t}).$

As such, (43) becomes

$\displaystyle\sum_{r=1}^{s}\mathbbm{1}(A_{1:r-1}=a_{1:r-1},A_{r}\neq
a_{r})\,g_{\mathrm{lo}}(X_{0:r-1}(a_{1:r-1}),x_{r:t})$
$\displaystyle\qquad\leq\sum_{r=1}^{s}\mathbbm{1}(A_{1:r-1}=a_{1:r-1},A_{r}\neq
a_{r})\,y_{\mathrm{lo}}(X_{0:r-1}(a_{1:r-1}),x_{r:t})$ (44)

Since $x_{1:s-1}$ was arbitrary, for $r\in\\{1,\ldots,s-1\\}$, we can