Training and Evaluating Causal Forecasting Models for Time-Series

## Abstract

Abstract Deep learning time-series models are often used to make forecasts that inform downstream decisions.
Since these decisions can differ from those in the training set, there is an implicit requirement that time-series models will generalize outside of their training distribution.
Despite this core requirement, time-series models are typically trained and evaluated on in-distribution predictive tasks.
We extend the orthogonal statistical learning framework to train causal time-series models that generalize better when forecasting the effect of actions outside of their training distribution. To evaluate these models, we leverage Regression Discontinuity Designs popular in economics to construct a test set of causal treatment effects.

## 1 Introduction

Deep learning models have recently seen major improvements, including in time-series modelling Liu et al. (); Nie et al. (); Chen et al. (); Das et al. (); Zhang & Yan ().
Time-series models are often used to make forecasts that inform decisions, such forecasting demand for goods or transportation to optimize pricing, or forecasting health markers to optimize treatment decisions Makridakis et al. (); Nowroozilarki et al. (); Makridakis et al. ().

The decision making procedures that leverage time-series forecast typically optimize for a given outcome, such as revenue or health outcome, and lead to actions that differ from those observed in the training data.
This results in a key implicit requirement on time-series models: that they will generalize outside of the observational distribution.
Despite this requirement, the only complex tasks on which time-series forecasting models are trained and evaluated are in-distribution predictive tasks Liu et al. (); Melnychuk et al. (); Lim (); Nie et al. (); Chen et al. (); Das et al. (); Zhang & Yan (); Nowroozilarki et al. (); Makridakis et al. (; ).

In this work, we develop a training and evaluation procedure for causal forecasting time-series models. Causal models aim to learn the causal relationship between a treatment (e.g., price) and its effect on a specific outcome (e.g., demand), conditioned on other, observational features.
That is, they aim to capture the change in outcome caused by a change in treatment, and not just predict the outcome from an in-distribution observed treatment.
As a result, they generalize better when forecasting the effect of actions outside of their training distribution.

A few causal forecasting models have been proposed,
but they use approaches that do not map to well defined causal effects Bica et al. (; ); Melnychuk et al. (), do not correct for regularization bias Li et al. (); Bica et al. (; ); Melnychuk et al. (); Lim (), or take inspiration from causal frameworks for heuristics Brodersen et al. (); Schultz et al. ().
In this work, we leverage the orthogonal statistical learning framework of Foster & Syrgkanis () to learn a causal forecasting model. We extend this framework to a well-defined time-series causal problem and to high-dimensional treatments, and instantiate it on top of state-of-the-art backbone time-series models to support complex use-cases.

Evaluating causal forecasting models is even more challenging, as outcomes under alternative treatments are by definition unobserved. Existing work thereby rely on simple simulations with known counterfactuals Bica et al. (; ); Li et al. (); Lim (), or in-distribution forecasting performance Lim (); Melnychuk et al. (). The former is often too simplistic to really compare models, while the later can be misleading and inflate the performance of non-causal models. We show details of this effect in §2.1.

To evaluate forecasting models on causal tasks, we take inspiration from Regression Discontinuity Designs (RDDs) popular in econometrics Hahn et al. (); Imbens & Lemieux (), to construct a test set of causal effects estimated under different assumptions than the typical model-fitting assumption. We can then compare the treatment effects estimated by forecasting models using traditional evaluation metrics.
We train and evaluate our causal forecasting models on two relevant tasks: demand forecasting using a proprietary dataset of passenger rail prices and demand over time; and the public health dataset MIMIC-III Johnson et al. ().
We show that, compared to traditional time-series forecasting models and a state-of-the-art causal model, our causal forecasting models predict causal effects that are $36\%$ and $1\%$ closer to those estimated with RDDs, respectively.

## 2 Causal Forecasting with Orthogonal Learning

In the reminder of this paper, we use a demand forecasting task for passenger rails as a running example. For each scheduled train, models aim to predict the daily number of seats sold at a given price using temporal information and other meta-data (see §4.1 for details). The end goal of these forecasts is to set prices to maximize revenue.
In §4, we evaluate our models on this task using a proprietary dataset, as well as on a public health forecasting dataset Johnson et al. ().

We abstract our forecasting task as a dataset of $N$ observed time-series ($N$ trains) indexed by $n$, with time indexed by $t$. Each time-series $n$ consists in stationary context features $S^{n}$, temporal (non-stationary) features measured at each time step $X^{n}_{1:t}$, a treatment at each time step $T^{n}_{1:t}$, and real valued outcomes $Y^{n}_{1:t}$ ($Y^{n}_{t}\in{\mathbb{R}}$). The objective is to predict each outcome $Y^{n}_{t}$ using static, temporal, and treatment features, as well as past outcomes up to a time $\tau<t$ ($Y^{n}_{\tau+1:t-1}$ is unknown).
We thus seek to learn $\hat{Y}^{n}_{t}\triangleq f(S^{n},X^{n}_{1:t},T^{n}_{1:t},Y^{n}_{1:\tau})$, such that $\hat{Y}^{n}_{t}$ is close to $Y^{n}_{t}$. This setting differs from standard time-series forecasting tasks Liu et al. (), but is more relevant to causal tasks (§2).

We next formalize causal effects, show their importance to our task, and introduce orthogonal learning theory that we leverage in our models (§2.1). We then extend orthogonal learning to time-series models (§2.2), and instantiate this theory with deep learning architectures (§2.3).

### 2.1 Background: Causal Treatment Effects and Orthogonal Learning

To formalize the importance of causality in forecasting for optimizing decisions, we first need to introduce treatment effects. In the reminder of this paper, we use the potential outcomes (or Rubin’s) model of causal inference Holland (); Imbens & Rubin (), and flexible data models from the causality Robinson (); Nie & Wager (); Chernozhukov et al. ().

Potential Outcomes represent the values of an outcome under different possible treatments. Consider our demand forecasting setting with two possible prices (treatments) $T_{1}$ and $T_{2}$. For a given time-series $n$ and time step $t$, we define the potential outcomes $Y^{n}_{t}(T_{1})$ and $Y^{n}_{t}(T_{2})$ as the value of the demand $Y^{n}_{t}$ we would observe under each price. Of course, we only ever observe one potential outcome, as $Y^{n}_{t}=\mathds{1}_{T_{1}}Y^{n}_{t}(T_{1})+\mathds{1}_{T_{2}}Y^{n}_{t}(T_{2})$ where $\mathds{1}$ is the indicator function.

Treatment Effects represent the causal effect of a change in treatment (price) on the outcome of interest (demand). Formally, the treatment effect of going from price $T_{1}$ to $T_{2}$, for time-series $n$, at time step $t$, is $Y^{n}_{t}(T_{2})-Y^{n}_{t}(T_{1})$.
The fundamental challenge of causal inference is that we can only ever observe one potential outcome, for the treatment we chose. As a result, we cannot train a supervised treatment effects model. However, we can estimate expected causal effects under unconfoundedness, the main assumption in the literature which we make when learning our models:

###### Assumption 1 (Unconfoundedness) .

Consider a set of possible treatments $\mathcal{T}$. We say that the time-series $(Y_{t}(T_{t}),T_{t},S,X_{t})$ are (conditionally) uncounfounded if, for every time-step $t>\tau$:

$$ $\{Y_{t}(T_{t}),~{}T_{t}\in\mathcal{T}\}\perp\!\!\!\!\perp T_{t}\mid S,~{}X_{1: t},~{}T_{1:\tau},~{}Y_{1:\tau}$ $$

Intuitively, this assumption means that in sub-populations defined by the conditioning variables, the groups receiving each treatment are comparable (have the same potential outcomes distribution).
Unconfoundedness is verified by construction in (conditionally) randomized trials, as treatments are assigned independently of potential outcomes.
It is also verified when context variables $S,~{}X_{1:t},~{}T_{1:\tau},~{}Y_{1:\tau}$ include all pre-treatment variables that are common causes of $T^{n}$ and $Y^{n}$ Pearl (). In our running example, this consists in all variables that influence both demand and pricing decisions, such as the cities of arrival and departure, day of the trip, and more. Under Assumption 1, a well known identifyability result holds (see e.g., Imbens & Rubin (); Pearl ()):

###### Proposition 1 (CATE Identifyablity) .

Under Assumption 1, we can identify the Conditional Average Treatment Effect (CATE) as follows:

$$ CATE $\displaystyle\triangleq\mathbb{E}\big{[}Y_{t}(T_{2})-Y_{t}(T_{1})\mid S,~{}X_{ 1:t},~{}T_{1:\tau},~{}Y_{1:\tau}\big{]}$ $\displaystyle=\mathbb{E}\big{[}Y_{t}\mid T_{t}=T_{2},~{}S,~{}X_{1:t},~{}T_{1: \tau},~{}Y_{1:\tau}\big{]}-\mathbb{E}\big{[}Y_{t}\mid T_{t}=T_{1},~{}S,~{}X_{1 :t},~{}T_{1:\tau},~{}Y_{1:\tau}\big{]}.$ $$

This important result shows that we can estimate the Conditional Average Treatment Effect (CATE) from conditional expectations $\mathbb{E}\big{[}Y_{t}\mid T_{t}=T,~{}S,~{}X_{1:t},~{}T_{1:\tau},~{}Y_{1:\tau}
\big{]}$. This is a familiar quantity to ML practitioners, as a flexible forecasting model $f:T,S,X,T_{1:\tau},Y_{1:\tau}\rightarrow Y$ trained with a mean-squared error will converge to this conditional expectation in the limit of infinite data. Estimating the CATE $\mathbb{E}\big{[}Y_{t}(T_{2})-Y_{t}(T_{1})\mid S,~{}X_{1:t},~{}T_{1:\tau},~{}Y
_{1:\tau}\big{]}$ is then as simple as predicting on the two treatment values $f(T_{t}=T_{2},\ldots)-f(T_{t}=T_{1},\ldots)$.
Unfortunately, practice is more challenging.

The need for causal forecasting arises when fitting complex forecasting models with finite sample sizes. In observational data collected from regular interactions, such as observed demand for passenger rail seats during regular operation, some treatments will be over-represented. For instance, in our data we observe high prices when demand is high, and low prices when demand is low. Of course, we would expect higher prices to cause a decrease in demand. But the correlation in the data is flipped, as the transportation company operators use their experience and domain knowledge about demand to set higher prices when underlying demand is high, thereby increasing revenue.

Figure: (a) TFT
Refer to caption: extracted/5969790/images/background/daily_wrong_gras.png

When trained on finite data with such biased treatment assignments, forecasting models will extrapolate based on feature and treatment values. This introduces bias in treatment effect estimates when models are not flexible enough to represent the true function, or due to regularization (both explicit regularization such as weight decay, or deep learning models’ implicit regularization (Soudry et al., ; Ji et al., )). This effect, called regularization bias, happens even under simple treatment effect functional forms (Chernozhukov et al., ).
In practice, this is exactly what we observe when fitting ML models on our demand forecasting task.
Figure 1 shows predicted treatment effects from a well tuned LightGBM Ke et al. () boosted tree and an improved version of the Temporal Fusion Transformer (TFT) Lim et al. () deep learning model. Both models have comparable in distribution performance, with a small edge for the LGBM model (details in §4). Figures 1(a) and 1(b) show predictions for a specific train: both models predict demand quite accurately using the realized price trajectory (orange line, $T^{n}_{1:t}$ is set to the observed prices). However, when increasing the prices by 30% for all $t$ (green line), both models predict an increase in demand. This is the correlation we observe in the data, of opposite sign to the causal relationship we expect. Figure 1(c) shows that this affects a significant proportion of test set trains: the LGBM and TFT models respectively predict 13% and 28% of positive treatment effects.

Figure: Algorithm 1 Sample Splitting Learning Algorithm from Foster & Syrgkanis ()

Orthogonal Learning (Foster & Syrgkanis, ) is a recent learning theory framework to study the convergence rate of models with nuisance parameters, which can be used to estimate CATE models. The framework is very general: we next describe the core results relevant our work, that we later extend to our setting (§2.2). The results for treatment effect estimation consider a standard binary treatment effect data model (Robinson, ; Nie & Wager, ; Chernozhukov et al., ), with observations $z=(W,Y,T)$ from distribution $\mathcal{D}$, and $T\in\{0,1\}$, such that:

$$ $\displaystyle Y=T\cdot\theta_{0}(W)+f_{0}(W)+\epsilon_{1},$ $\displaystyle\mathbb{E}[\epsilon_{1}|W,T]=0$ (1) $\displaystyle T=e_{0}(W)+\epsilon_{2},$ $\displaystyle\mathbb{E}[\epsilon_{2}|W]=0$ (2) $$

Further define $m_{0}(w)\triangleq\mathbb{E}[Y|W=w]$. We call $g_{0}\triangleq(e_{0},m_{0})$ the nuisance parameters, required to learn the CATE $\theta_{0}$ but not direct objects of interest in this framework. In this model $\theta_{0}(W)$ is the CATE (as it is turned on and off by $T$, see Proposition 1).
This data model is quite flexible, as the CATE depends arbitrarily on the features $W$ through $\theta_{0}$, as do the outcome $Y\in{\mathbb{R}}$ and treatment $T$ through $f_{0},e_{0}$. Note that under Assumption 1, this decomposition is without loss of generality.

Calling $\theta,e,m$ models for $\theta_{0},e_{0},m_{0}$, respectively, we can define the residualized loss (R-loss) $l(\cdot)$, and its population equivalent $L_{\mathcal{D}}(\cdot)$, as:

$$ $l(\theta,m,e;z)=\Big{(}Y-m(W)-\big{(}T-e(W)\big{)}\theta(W)\Big{)}^{2}\ \ ;\ \ L_{\mathcal{D}}(\theta,m,e)=\mathbb{E}_{z\sim\mathcal{D}}[l(\theta,m,e;z)].$ $$

Now consider a model $g$ for $g_{0}$ learned on sample $S$, such that $\|g-g_{0}\|_{L_{2}}\triangleq\sqrt{\mathbb{E}_{z\sim\mathcal{D}}\|g(w)-g_{0}(w
)\|^{2}_{2}}\leq\textrm{Rate}_{\mathcal{D}}(S,\delta)$ w.p. at least $1-\delta$. Since $m_{0}(w)=\mathbb{E}[Y|W=w]$ and $e_{0}(w)=\mathbb{E}[T|W=w]$, these can be models that minimize the mean squared error.

Further consider an estimation algorithm of model $\theta$ for $\theta_{0}$ which, given a model $g$ for nuisance parameters, outputs $\theta\in\Theta$ using sample $S$ (e.g., by minimizing the empirical loss $l(\cdot)$ over $S$) such that $L_{\mathcal{D}}(\theta,g)-L_{\mathcal{D}}(\theta^{*},g)\leq\textrm{Rate}_{
\mathcal{D}}(S,\delta;\theta,g)$, where $\theta^{*}\triangleq\arg\min_{\theta\in\Theta}L_{\mathcal{D}}(\theta,m_{0},e_{
0})$ is the best predictor for $\theta_{0}$ in class $\Theta$.
A low risk $L_{\mathcal{D}}(\theta,g)-L_{\mathcal{D}}(\theta^{*},g)$ does not indicate a good model $\theta$, as the loss $l$ is computed with $g$, and not $g_{0}$ the true nuisance parameter.
However, Foster & Syrgkanis () show (in Theorem 2, and Appendix J.2, with a constant we derive in Appendix B) that:

###### Proposition 2 (Orthogonal Learning for Binary Treatment Effect) .

Consider a class of functions $\Theta$ where $\forall W,\ \forall\theta\in\Theta:\ |\theta(W)|\leq M$.
Fitting $g$ and $\theta$ with the procedures above using Algorithm 1, we have that with probability at least $1-\delta$:

$$ $\displaystyle\mathbb{E}\big{[}\big{(}T-e_{0}(W)\big{)}\big{(}\theta(W)-\theta^ {*}(W)\big{)}\big{]}^{2}$ $\displaystyle\leq L_{\mathcal{D}}(\theta,g_{0})-L_{\mathcal{D}}(\theta^{*},g_{ 0})$ $\displaystyle\leq\textrm{Rate}_{\mathcal{D}}(S_{2},\delta/2;\theta,g)+2(1+M^{2 })\cdot\textrm{Rate}_{\mathcal{D}}(S_{1},\delta/2)^{2}.$ $$

There are three key takeaways. (1) we can get a CATE model $\theta$ close to the best possible predictor in the class (which can be $\theta_{0}$ if $\theta_{0}\in\Theta$) for all values $W$ with treatment variation around the mean. This is similar to the positivity assumption ($\forall T,W:1>{\mathbb{P}}(T|W)>0$) pervasive in causal inference results. (2) the rate at which $\theta$ converges to $\theta^{*}$ under estimated nuisance parameters $g$ is sufficient to bound the error under the true parameters $g_{0}$. (3) the rate of estimation of the nuisance parameters is squared: Algorithm 1 is tolerant of errors in estimating complex nuisance parameters. For instance, a rate of convergence of $1/\sqrt{N}$ for $\theta$, is preserved if $g$ converges as $o(1/N^{\frac{1}{4}})$.

###### Assumption 1 (Unconfoundedness) .

###### Proposition 1 (CATE Identifyablity) .

###### Proposition 2 (Orthogonal Learning for Binary Treatment Effect) .

### 2.2 Orthogonal Learning for Time Series Models

We adapt the orthogonal statistical learning framework described above by making two key extensions: defining daily treatment effects with an observation cut-off; and extending the R-loss to categorical and linear effects with various encodings for $\theta$’s predictions.

Treatment effects in time series. The first step is to formalize the treatment effects to predict, and the context to predict them with. In this paper, we focus on daily treatment effects defined in Proposition 1. While technically past prices are likely to causally impact demand at $t$, domain knowledge tells us that the price at $t$ is the dominant factor. This modelling choice lets us focus on single-day treatments, though supporting sequences of treatments is an interesting avenue for future work. Focusing on a time-step $t$, this implies that $Y=Y_{t}$ in Equation 1 and $T=T_{t}$ in Equation 2.

Defining the relevant feature set ($W$ in Eq. 1, 2) is more tricky. Recall from §2 that our time series consist in static ($S$) and temporal ($X_{1:t}$) features. These features are naturally included ($W\supset(S,X_{1:t})$).
In time-series, results from past time steps, such as past prices and associated demands $(T_{1:t-1},Y_{1:t-1})$, are typically very predictive.
However, both the uncondoundedness Assumption 1 and the orthogonal leaning data model (Eq. 1, 2) implicitly assume that features $W$ follow the observational distribution at training and prediction time.
This assumption is not verified when models are used to influence past decisions, and $T_{<t},Y_{<t}$ cannot be included in $X$.
This is an important point: observational past treatment are very predictive of future demand, but confound the daily treatment effect and are the cause of reversed elasticity predictions we observed (§2.1, Figure 1).
To avoid this confounding, we introduce a cut-off $\tau$ after which treatments can deviate from the distribution. That is, in time-steps $1:\tau$ the model is not used to affect treatments, while for $t>\tau$, the treatment $T_{t}$ can differ from the training distribution. We then set the features as $W=(S,X_{1:t},T_{1:\tau},Y_{1:\tau})\triangleq W_{t}$, which we use as input for our daily treatment effect of $T_{t}$ on $Y_{t}$.

High dimensional effects. The other extension we require is support for higher-dimensional treatments. We consider both categorical and linear treatments and treatment effects. We start with categorical treatments, for which we propose two different encodings. In our demand forecasting dataset, possible prices are discrete, following a large but tractable set of possible values (about thirty possible prices) of size $d\triangleq|\mathcal{T}|$. We extend the approach of Robinson (); Nie & Wager (); Foster & Syrgkanis () (Equations 1, 2) to model each time-step $t$ as follows:

$$ $\displaystyle Y_{t}=T_{t}^{T}\theta_{0}(W_{t})+f_{0}(W_{t})+\epsilon_{1},$ $\displaystyle\mathbb{E}[\epsilon_{1}|W_{t},T_{t}]=0$ (3) $\displaystyle T_{t}=e_{0}(W_{t})+\epsilon_{2},$ $\displaystyle\mathbb{E}[\epsilon_{2}|W_{t}]={\bf 0}_{d}$ (4) $$

where $T_{t}$ is a $d$-dimensional column vector encoding the treatment, $\theta_{0}(W_{t})\in{\mathbb{R}}^{d}$ outputs a column vector treatment effect, and ${\bf 0}_{d}$ is the $d$-dimensional zero vector.
The R-loss becomes:

$$ $l(\theta,{m,e},z_{t})=\Big{(}Y_{t}-m(W_{t})-\big{(}T_{t}-e(W_{t})\big{)}^{T} \theta(W_{t})\Big{)}^{2}.$ (5) $$

We consider two encodings for the categorical treatment $T_{t}$ and treatment effect $\theta(W_{t})$. In the one-hot encoding, $T_{t}$ is a one-hot vector of the treatment. As a result, the $i^{th}$ dimension of the treatment effect model encodes the CATE compared to $f_{0}$. That is, $\theta(W_{t})_{i}=\mathbb{E}\big{[}Y_{t}(T_{i})-f_{0}(W_{t})\ |\ W_{t}\big{]}$.

In the cumulative encoding, the treatment model encodes the treatment effect of incremental price changes, such that the cumulative predictions encode the CATE: $\sum_{j=1}^{i}\theta(W_{t})_{j}=\mathbb{E}\big{[}Y_{t}(T_{i})-f_{0}(W_{t})\ |
\ W_{t}\big{]}$. In this case, $T_{t}$ is a vector of ones in dimensions $1:i$, and zeros after. While we do not do it in this work, the cumulative encoding can be useful to encode constraints, such as fixing the sign of incremental effects of price changes.

Given true categorical treatments and our flexible deep-learning models, we can reasonably assume that $\theta_{0}\in\Theta$, and thus $\theta^{*}=\theta_{0}$. With enough treatment variation, our model will converge to the true CATE.
In practice however, demand forecasting datasets only see local variation around a typical price for any particular value of features $W_{t}$. As a result, a small risk under the R-loss
only ensures an accurate $\theta$ around likely prices. The model is under-constrained further from typical prices, where it will not learn relevant treatment effects.
This is a fundamental limitation of any causal approach, as local variations around a typical price implies a lack of positivity for most dimensions of $\theta(W_{t})$.

In such cases, adding structure to the treatment effect model is practically useful. To this end, we also consider a linear encoding for the treatment effect. Concretely, we use the data model of Equations 3 and 4, where $T$ is a one-dimensional, real treatment value (the price at which we want to predict demand), and the treatment effect $\theta(W_{t})$ is a multiplicative coefficient. This creates a linear effect model, where the linear coefficient depends on $W_{t}$ in a flexible way. As a result, the real treatment effect is probably not in the model class ($\theta_{0}\not\in\Theta$), and we converge to $\theta^{*}$. In exchange, estimations under this model can be much more stable under moderate treatment variability.

All encodings are compatible with the following result:

###### Proposition 3 (Orthogonal Learning for Categorical Treatment Effects) .

Consider a class of functions $\Theta:\mathcal{W}\rightarrow{\mathbb{R}}^{d}$, where $\forall W,\ \forall\theta\in\Theta,\ \forall i\in\{1,d\}:\ |\theta(W)_{i}|\leq
M$.
Fitting $g$ and $\theta$ with the procedures above using Algorithm 1, we have that with probability at least $1-\delta$:

$$ $\displaystyle\mathbb{E}\big{[}\big{(}T-e_{0}(W)\big{)}^{T}\big{(}\theta(W)- \theta^{*}(W)\big{)}\big{]}^{2}$ $\displaystyle\leq L_{\mathcal{D}}(\theta,g_{0})-L_{\mathcal{D}}(\theta^{*},g_{ 0})$ $\displaystyle\leq\textrm{Rate}_{\mathcal{D}}(S_{2},\delta/2;\theta,g)+2(1+dM^{ 2})\cdot\textrm{Rate}_{\mathcal{D}}(S_{1},\delta/2)^{2}.$ $$

###### Proof.

Under the data model of Equations 3, 4, the categorical R-loss from Equation 5 satisfies universal Neyman orthogonality (Assumption 4, proof in Appendix B.2), and has continuous and bounded second directional derivatives (Assumption 5, proof in Appendix B.3) with $\beta=2(1+dM^{2})$. Applying Theorem 2 from Foster & Syrgkanis () concludes the proof.
∎

Since for the linear encoding (and for binary treatments from Proposition 2) $d=1$, the rate of convergence for nuisance parameters is less impactful than for more flexible categorical encodings.

###### Proposition 3 (Orthogonal Learning for Categorical Treatment Effects) .

###### Proof.

### 2.3 Causal Time Series Forecasting Models with Deep Learning

To instantiate the theory described in §2.2, we extend time-series architectures for deep learning, that we use as backbones for models $m,e,\theta$. We then adapt the training procedure to fit Algorithm 1. Finally, we change the prediction procedure to output causal forecasts.

Backbone Time-Series Model Architecture with Orthogonal Learning.
Our time-series backbone takes as input static features $S$, temporal features $X_{1:t}$, and treatments and outcome values before $\tau$, $T_{1:\tau}$ and $Y_{1:\tau}$. It is then used in three models: $e(\cdot)$, which predicts the treatment sequence $T_{>\tau}$; $m(\cdot)$, which predicts expected outcomes (without knowing the treatments); and $\theta(\cdot)$, which predicts a vector of treatment effects (interpreted differently depending on the encoding, see §2.2). Our inputs are non-traditional, and the reason why our main backbone is a modified TFT Lim et al. () (modification details in Appendix C.3). We also experiment with the state-of-the-art iTransformer architecture Liu et al. (), passing static features’ embeddings along temporal ones, and training with randomly sampled $t$ values by truncating the time-series at $t$, as the iTransformer processes entire time-series without a causal structure (i.e., the model “knows the future”). This results in very slow training and prediction, and we can only apply it to one of our datasets.

Fitting a Causal Orthogonal Time-Series Model.
We fit our models following Algorithm 1. We train both models $m(W_{t})$, $e(W_{t})$ as estimators for $\mathbb{E}[Y_{t}|W_{t}]$ and $\mathbb{E}[T_{t}|W_{t}]$ on subset $S_{1}$. $m(\cdot)$ is trained with the mean squared error (MSE) loss. The loss for $e(\cdot)$ depends on the encoding function. We use the cross entropy for the one-hot encoding case, the binary cross-entropy for the cumulative encoding, and the MSE for the linear encoding. In all cases, we train $\theta(\cdot)$ by minimizing the R-loss $l(\cdot,m,e,z_{t})$ from Eq. 5 with data points $z_{t}\in S_{2}$.

Forecasting with a Causal Orthogonal Time-Series Model.
$\theta$ predicts causal changes of outcomes under different treatments.
However, to optimize downstream decisions (e.g. chose prices to maximize revenue) we need a proper forecast for the outcome $Y_{t}$.
Our final estimator combines all models following Eq. 3:

$$ $\hat{Y}_{t}(W_{t})=m(W_{t})+\big{(}T_{t}-e(W_{t})\big{)}^{T}\theta(W_{t}).$ $$

This way under $g_{0}$ we have $\hat{Y}_{t}(W_{t})=f_{0}(W_{t})+e_{0}(W_{t})^{T}\theta(W_{t})+\big{(}T_{t}-e_{
0}(W_{t})\big{)}^{T}\theta(W_{t})=f_{0}(W_{t})+T_{t}^{T}\theta(W_{t})$. When $\theta(\cdot)$ is close to the CATE, our predictions change causally with $T_{t}$.

## 3 Regression Discontinuity Design for Causal Model Evaluation

Once our causal forecasting models are fitted, how can we measure their performance? This is a core challenge in causal models, which fundamentally rely on untestable assumptions Pearl (), such as Assumption 1 or related data models as in Equations 1-4.
In §2.1, we used domain knowledge on the sign of the CATE to argue that non-causal models fail to learn valid causal relationships. However, these are qualitative observations that do not enable comparisons between different models.
To the best of our knowledge, their is currently no reliable method to quantitatively evaluate the quality of causal models, without relying on the same assumptions used in fitting the models.

To address this challenge, we take inspiration from robustness checks common in econometrics, which estimate causal effects of interest using multiple approaches. We scale this approach to evaluate ML models. Specifically, we leverage Regression Discontinuity Designs (RDDs) (§3.1) and changes in treatments over time to estimate a subset of daily CATE values in a time-series. This approach uses a different assumption than Assumption 1 (§3.2). We then compute CATE values for our test set time-series, and use the resulting causal test to evaluate and compare ML models (§3.3).

### 3.1 Background: Regression Discontinuity Designs

Regression Discontinuity Designs (RDD) leverage a continuity assumption to estimate a treatment effect at a cut-off point that triggers a change in treatment Hahn et al. (); Imbens & Lemieux (); Lee & Lemieux ().
We present a slightly non-traditional version that conditions on additional variables, needed for our demand forecast setting.
Formally, RDDs require a variable X with an associated cut-off value $c$ that corresponds to a change in treatment. That is, $T=T_{1}$ when $X<c$, and $T=T_{2}$ when $X\geq c$. RDDs provide identifyability under the following assumption:

###### Assumption 2 (Continuity) .

The potential outcomes’ conditional expectation, $\mathbb{E}[Y(T_{2})|X=x,V]$ and $\mathbb{E}[Y(T_{1})|X=x,V]$ respectively, are both continuous in $x$.

Under assumption 2, we can estimate the CATE at the cutoff $X=c$:

###### Proposition 4 (CATE Identifyablity from RDD) .

Under Assumption 2, we have that:

$$ $\mathbb{E}\big{[}Y(T_{2})-Y(T_{1})\mid X=c,V\big{]}=\lim_{x\to c^{+}}\mathbb{E }\big{[}Y\mid X=x,V\big{]}-\lim_{x\to c^{-}}\mathbb{E}\big{[}Y\mid X=x,V\big{]}$ $$

Practical estimators based on Proposition 4 typically fit linear or polynomial models $g:\mathcal{X},\mathcal{V}\rightarrow\mathcal{Y}$, using an indicator variable for the cut-off to measure the discontinuity at the point of change of treatment $X=c$. This discontinuity captures the CATE at $X=c$. Since we are interested in the limit at $X=c$, one often uses a weight kernel $K(\cdot)$ that decreases the importance of datapoints further from the cut-off $X=c$. Typical weight kernels include the rectangular ”window” kernel of width $h$, or a triangular kernel in which weights decay linearly when moving away from $X=c$, until they reach zero outside of $[c-h,c+h]$.

###### Assumption 2 (Continuity) .

###### Proposition 4 (CATE Identifyablity from RDD) .

### 3.2 Estimating Treatment Effects in Time Series

RDDs are attractive to evaluate our models, as Assumption 2 differs from uncounfoundedness (Assumption 1) which underpins the learning of causal models.
To causally evaluate models, we want to compare a model’s CATE prediction on a given time-series $n$ at time-step $t$, to another estimate of this CATE used as ground truth.
We obtain this ground truth by framing observed changes in treatment on individual time-series as an RDD.
Consider time-series $n$. We call $t^{n}_{i}$ the $i^{th}$ switching time-step at which a price change happens, with $T_{1}\triangleq T_{t^{n}_{i}-1}\neq T_{t^{n}_{i}+1}\triangleq T_{2}$.
Because we use aggregated time-series (e.g., at the daily level), and treatment changes happen in the middle of an observation, $T_{t^{n}_{i}}$ is often ill defined. While this deviates from traditional RDD formulations, in which $\mathbb{E}\big{[}Y\mid X=c,V\big{]}=\mathbb{E}\big{[}Y(T_{2})\mid X=c,V\big{]}$, Prop. 4 and associated estimators still apply.

Another deviation from §3.1 is that our time-steps are discrete, so we technically cannot take the upper- and lower-limits at the $t^{n}_{i}$. We follow Lee & Card ()
and assume that specification errors at $t^{n}_{i}$ (the deviation from $g$’s estimates and the true conditional expectation on each side of the cut-off) are zero in expectation. Formally:

###### Assumption 3 (Unbiased specification errors) .

Consider (continuous at $c$) RDD models $g_{T_{1}}:x\rightarrow g(T_{1},x)$ and $g_{T_{2}}:x\rightarrow g(T_{2},x)$ for time-series $n$.
For $i\in\{1,2\},\ \mathbb{E}[Y_{t^{n}_{i}}(T_{i})-g_{T_{i}}(c)|V]=0$.

In demand forecasting, assumptions 2 and 3 with $V=\emptyset$ are unrealistic, as there are known cyclical patterns of demand based on days of the week. In §3.3, we show how we use $V$ to correct for these effects, making Assumptions 2 and 3 more legitimate.
With these assumptions and estimator, and denoting $w^{n}_{t^{n}_{i}}$ the features we give our orthogonal forecasting model at the cut-off (see §2.2), and $g^{n}$ the associated RDD model, we have the following result:

###### Proposition 5 (RDD for point CATE) .

Consider time-series $n$ and switching-time $t^{n}_{i}$:

$$ $\mathbb{E}\big{[}Y_{t^{n}_{i}}(T_{t^{n}_{i+1}})-Y_{t^{n}_{i}}(T_{t^{n}_{i}})\ |W_{t^{n}_{i}}=w^{n}_{t^{n}_{i}},V\big{]}=\lim_{t\to{t^{n}_{i}}^{+}}\mathbb{E} \big{[}g^{n}(t)\mid~{}V\big{]}-\lim_{t\to{t^{n}_{i}}^{-}}\mathbb{E}\big{[}g^{n }(t)\mid~{}V\big{]}.$ $$

###### Proof.

We start by using the fact that an i.i.d. time-series with a given context $W_{t}$ is an unbiased (one point) estimate of the conditional expectation over the data, that is: $\mathbb{E}\big{[}Y_{t^{n}_{i}}(T_{t^{n}_{i+1}})-Y_{t^{n}_{i}}(T_{t^{n}_{i}})\
|W_{t^{n}_{i}}=w^{n}_{t^{n}_{i}},\mid t=t^{n}_{i}\big{]}=\mathbb{E}\big{[}Y^{n
}_{t^{n}_{i}}(T_{t^{n}_{i+1}})-Y^{n}_{t^{n}_{i}}(T_{t^{n}_{i}})\big{]}=\mathbb
{E}\big{[}Y^{n}_{t}(T_{t^{n}_{i+1}})-Y^{n}_{t}(T_{t^{n}_{i}})\mid t=t^{n}_{i}
\big{]}$.

Invoking Assumption 2 and Proposition 4 yields $\mathbb{E}\big{[}Y^{n}_{t}(T_{t^{n}_{i+1}})-Y^{n}_{t}(T_{t^{n}_{i}})\mid t=t^{
n}_{i}\big{]}=\lim_{t\to{t^{n}_{i}}^{+}}\mathbb{E}\big{[}Y^{n}\mid t,V\big{]}-
\lim_{t\to{t^{n}_{i}}^{-}}\mathbb{E}\big{[}Y^{n}\mid t,V\big{]}$.
Finally, Assumption 3 gives:
$\lim_{t\to{t^{n}_{i}}^{+}}\mathbb{E}\big{[}Y^{n}\mid t,V\big{]}-\lim_{t\to{t^{
n}_{i}}^{-}}\mathbb{E}\big{[}Y^{n}\mid t,V\big{]}=\lim_{t\to{t^{n}_{i}}^{+}}
\mathbb{E}\big{[}g^{n}(t)\mid V\big{]}-\lim_{t\to{t^{n}_{i}}^{-}}\mathbb{E}
\big{[}g^{n}(t)\mid V\big{]}$.
∎

###### Assumption 3 (Unbiased specification errors) .

###### Proposition 5 (RDD for point CATE) .

###### Proof.

### 3.3 Causal Test Sets to Evaluate Causal Forecasting Models

Proposition 5 lets us estimate the CATE on a specific time-series $n$, at a specific time $t^{n}_{i}$, for a specific treatment change $T_{t^{n}_{i}-1}\rightarrow T_{t^{n}_{i}+1}$. We call this quantity $\text{CATE}_{t^{n}_{i}}$. To create our causal test set, we take the in-distribution task test set $S_{\textrm{test}}$, and collect all CATE values that we can estimate from it in $S^{\textrm{CATE}}_{\textrm{test}}$.
We only include switch-times $t^{n}_{i}$ that have at least three time-steps without treatment changes strictly before and after $t^{n}_{i}$, to have enough data for the RDD estimator.

Given enough data, we create a small dataset $D_{i}^{n}$ to fit our RDD model $g$, including all time-steps with constant treatments around $t^{n}_{i}$. Formally, $D_{i}^{n}=\{t:\forall t^{\prime}\in[t,t^{n}_{i}-1],\ T_{t^{\prime}}=T_{t^{n}_{
i}-1}\}\cup\{t:\forall t^{\prime}\in[t^{n}_{i}+1,t],\ T_{t^{\prime}}=T_{t^{n}_
{i}+1}\}$ (note that $t^{n}_{i}\not\in D_{i}^{n}$).
In tasks with known cyclical patterns, such as our demand prediction task with day of the week effects (see §3.2), we create a set of comparable time-series $S^{n}$ (using a rule informed by domain knowledge) on which we fit a linear regression model with time interactions $Y\sim V+V\times t$. For days of the week, this gives us parameters $\alpha_{1,j},\alpha_{2,j}$ for $j\in\{1,\ldots,7\}$. We then compute the residuals $\tilde{Y}_{t}=Y_{t}-\sum_{j=1}^{7}(\alpha_{1,j}+t\alpha_{2,j})$, on which we fit our RDD model.
From early experiments, we noticed that $L_{2}$ regularized linear models perform best for $g$.
We fit such a model with weight kernel $K(.)$ of window size $h$ on dataset $D_{i}^{n}$ using the following specification:

$$ $\tilde{Y}_{t}=\beta_{0}+\beta_{1}(t-t^{n}_{i})+\beta_{2}\mathds{1}_{t>t^{n}_{i }}+\beta_{3}\mathds{1}_{t>t^{n}_{i}}(t-t^{n}_{i})+\epsilon$ $$

Parameter $\beta_{2}$ corresponds to the estimate of $\text{CATE}_{t^{n}_{i}}$ that we put in our causal test set $S^{\textrm{CATE}}_{\textrm{test}}$.
Figure 2 represents the estimation of different $\text{CATE}_{t^{n}_{i}}$ for a given time series, with and without correction.

Figure: (a) RDDs without weekday correction
Refer to caption: extracted/5969790/images/rdd/rdd_fig_linear_small.png

Small sample sizes in RDD datasets make the final causal test set $S^{\textrm{CATE}}_{\textrm{test}}$ noisy. We filter outliers by cutting the 2.5% tails of the distribution on each side. With this test set of causal effect, we can evaluate predictions of causal effects from our models, as described in §2.3. We use traditional metrics such as the root-mean-square error (RMSE) or the mean absolute error (MAE). We can compare such metrics to causal effects estimated with non-causal models, by predicting with those models at different treatments (prices) and subtracting the results to estimate their predicted CATE.

## 4 Evaluation

### 4.1 Demand Forecasting Dataset

Dataset. Each time-series is a sequence of prices and number of sales for one train. Temporal features include the weekday of the timestep and the number of days before departure. Stationary features include the departure and arrival terminals and departure date time. The dataset consists in about 300,000 training time-series, and validation and test sets of around 100,000 each. All the figures and metrics regarding this dataset have been anonymized. Time steps $t$ are normalized to $[-1,0]$ (all time series have the same length).

Experimental setup. We set $\tau=0.33$, and forecast all outcomes $(Y_{t})_{t>\tau}$. Baseline (non-causal) models are trained with the MSE loss, and include a well-tuned LGBM model (used in production) and two deep learning architectures: the state-of-the-art iTransformer Liu et al. (), and a custom version of the Temporal Fusion Transformer (TFT) Lim et al. (). The TFT maps well to our requirements (see §2), and we modify it for performance (details in Appendix C.3).
We tune all models’ hyper-parameters separately and equally, and report metrics by training and testing with 5 random seeds on final hyper-parameters, and reporting the mean and standard-deviation.
The RDD dataset uses a linear weights kernel with $h=14$ time steps.
Keeping all switching times with at least 3 data points strictly before and after the price change retains 70% of observed price changes.

Results. Table 1 shows the performance of several models. We can see that non-causal models perform best in-distribution, and the LGBM is best among them in both RMSE and MAE. This is expected, as the observed prices are very informative in-distribution, since operators use their experience to set prices based on the demand they expect. On causal tasks however, causal models outperform the baselines, with linear models being particularly accurate in terms of RDD RMSE, which is expected in high teatment dimension $d$ (§2.2).
The Causal iTransformer with linear activation has an RDD RMSE 37% better than the baseline TFT, the next best non-causal model, a significant improvement.
Figure 3(a) corresponds to the same time series as figure 1(b)-1(a), and shows that the causal model captures the correct CATE sign (increasing prices decreases the demand forecast).
Figure 3(b) shows that the causal TFT’s CATE distribution is qualitatively much better: model predicts fewer positive CATE effects than baselines on the test set, and those positive CATEs are smaller in magnitude.

**Table 1: Baselines perform better in distribution (RMSE, MAE), while causal models better capture causal effects (RDD metrics).**
| Models | RMSE | MAE | RDD RMSE | RDD MAE |
| --- | --- | --- | --- | --- |
| LGBM | 0.6870 $\pm$ 0.0015 | 0.2877 $\pm$ 0.0005 | 0.2163 $\pm$ 0.0009 | 0.1390 $\pm$ 0.0004 |
| TFT Baseline | 0.7044 $\pm$ 0.0050 | 0.3023 $\pm$ 0.0053 | 0.1910 $\pm$ 0.0014 | 0.1262 $\pm$ 0.0005 |
| Causal TFT Linear | 0.7216 $\pm$ 0.0020 | 0.3104 $\pm$ 0.0047 | 0.1217 $\pm$ 0.0019 | 0.1271 $\pm$ 0.0008 |
| Causal TFT Cumulative | 0.7189 $\pm$ 0.0091 | 0.3067 $\pm$ 0.0019 | 0.1820 $\pm$ 0.0008 | 0.1243 $\pm$ 0.0004 |
| Causal TFT One-hot | 0.7236 $\pm$ 0.0115 | 0.3095 $\pm$ 0.0069 | 0.1874 $\pm$ 0.0007 | 0.1279 $\pm$ 0.0002 |
| Causal iTransformer Linear | 0.7504 $\pm$ 0.0098 | 0.3113 $\pm$ 0.0038 | 0.1197 $\pm$ 0.0014 | 0.1259 $\pm$ 0.0003 |
| Causal iTransformer Cumulative | 0.7377 $\pm$ 0.0070 | 0.3058 $\pm$ 0.0017 | 0.1908 $\pm$ 0.0021 | 0.1282 $\pm$ 0.0009 |
| Causal iTransformer One-hot | 0.7347 $\pm$ 0.0097 | 0.3036 $\pm$ 0.0030 | 0.1930 $\pm$ 0.0023 | 0.1288 $\pm$ 0.0009 |

Figure: (a) Causal TFT
Refer to caption: extracted/5969790/images/eval_WM/Daily_wrong_theta_linear.png

### 4.2 MIMIC Health Prediction Dataset

We also evaluate our approach on a public dataset, and release the code to replicate experiments here: [https://github.com/wiremind/causal-forecasting](https://github.com/wiremind/causal-forecasting).

Dataset. We use MIMIC-extract library Wang et al. () to process MIMIC-III Johnson et al. ().
Every time-series represents data from one patient who stayed between 30 and 60 hours in a critical care unit, a typical setting in causal inference benchmarks Bica et al. (); Lim (); Melnychuk et al. (). Static features include the patient’s gender, age, and ethnicity. Temporal features include 25 vital signals.
Our task is to forecast a patient’s blood pressure, and we aim at estimating the treatment effect of mechanical ventilation and vasopressor drugs.
Since this represents a categorical treatment, our linear encoding does not apply.

Experimental setup. We keep the default variable length time-series, and follow prior art in forecasting blood pressure for the next 5 hours (Lim, ; Melnychuk et al., ). That is, for each $t$ we set $\tau=t-1$ and forecast $[t,t+4]$.
We use the same baseline losses as in §4.1.
We were not able to fit our causal models with the iTransformer backbone, which requires a fixed $\tau$, so we focus on our causal TFT.
We also consider a state-of-the-art causal forecasting baseline, the Causal Transformer Melnychuk et al. (). We slightly modify the code to plug information leakage from the future vitals of hour $t=\tau+1$ (see Appendix C.1 for details).
We evaluate forecasts in distribution as prior work, as well as using our RDD methodology to create a causal test set. Our RDDs use a rectangular weight kernel with $h=5$.
Retaining switching times with at least three points strictly on each side also preserves about 70% of the treatment changes.

Results. We focus on RMSE results. MAE results are in Appendix D.1 and are qualitatively the same.
When forecasting in distribution, Table 2 (up) shows that baseline models outperform causal ones.
The auto-regressive Causal Transformer outperforms our TFT in next-step prediction, but is less accurate when predicting further into the future, which is consistant with prior observations Liu et al. (); Das et al. (); Wu et al. (); Ilbert et al. ().
Once again, baseline models are able to exploit the correlation between treatments and future outcomes, as doctors use their domain knowledge to anticipate ant prevent issues.

On the causal task using our RDD-based CATE estimates though, we can see on Table 2 (down) that the one-hot causal model outperforms the baseline TFT by at least 1.6% in all time-steps, and the Causal Transformer by almost 1% in the first time step.
Contrary to our Causal TFT, the Causal Transformer has access to the treatment sequence between $\tau$ and $t$. Since our RDD methodology only estimates daily treatment effects, the sequence of treatments $\tau:t$ is always in distribution, and the Causal Transformer is able to exploit this signal. However, this means that Causal Transformer CATE predictions for changes in treatment sequences would be confounded, and suffer some the same issues as the TFT baseline.
Finding similar estimators to our RDD approach for treatment sequences to rigorously evaluate these effects is an interesting avenue for future work.

**Table 2: Performance of next five steps forecasts, in distribution and for causal effects (RDD).**
| Models | $t$ = $\tau$+1 | $t$ = $\tau$+2 | $t$ = $\tau$+3 | $t$ = $\tau$+4 | $t$ = $\tau$+5 |
| --- | --- | --- | --- | --- | --- |
| In distribution RMSE |  |  |  |  |  |
| Causal Transformer | 8.823 $\pm$ 0.028 | 9.499 $\pm$ 0.024 | 9.845 $\pm$ 0.035 | 10.098 $\pm$ 0.028 | 10.313 $\pm$ 0.026 |
| TFT Baseline | 8.921 $\pm$ 0.056 | 9.534 $\pm$ 0.057 | 9.834 $\pm$ 0.061 | 10.051 $\pm$ 0.052 | 10.238 $\pm$ 0.043 |
| Causal TFT One-hot | 8.951 $\pm$ 0.023 | 9.596 $\pm$ 0.025 | 9.921 $\pm$ 0.038 | 10.165 $\pm$ 0.028 | 10.371 $\pm$ 0.027 |
| Causal TFT Cumulative | 8.945 $\pm$ 0.026 | 9.608 $\pm$ 0.033 | 9.935 $\pm$ 0.046 | 10.184 $\pm$ 0.035 | 10.390 $\pm$ 0.030 |
| RDD RMSE for CATE prediction |  |  |  |  |  |
| Causal Transformer | 2.887 $\pm$ 0.145 | 2.869 $\pm$ 0.145 | 2.844 $\pm$ 0.141 | 2.833 $\pm$ 0.140 | 2.833 $\pm$ 0.129 |
| TFT Baseline | 2.908 $\pm$ 0.102 | 2.942 $\pm$ 0.103 | 2.969 $\pm$ 0.114 | 3.011 $\pm$ 0.104 | 3.051 $\pm$ 0.096 |
| Causal TFT One-hot | 2.861 $\pm$ 0.094 | 2.861 $\pm$ 0.094 | 2.861 $\pm$ 0.094 | 2.862 $\pm$ 0.094 | 2.862 $\pm$ 0.094 |
| Causal TFT Cumulative | 2.895 $\pm$ 0.050 | 2.895 $\pm$ 0.050 | 2.895 $\pm$ 0.050 | 2.895 $\pm$ 0.050 | 2.895 $\pm$ 0.050 |

## Appendix A Encoding functions

In this section, we expand on the three different encoding functions that we plug into the R-Loss.

### A.1 Linear encoding

A possible encoding function is to pass the treatment as a scalar value in the R-Loss, giving a linear approximation of $\mathbb{E}[Y_{t}(T)|W_{t}]$ with regards to $T$. In this framework, we denote $T_{t}$ as the scalar value of the treatment. Hence the data model reduces to that of Eq. 1, but with a continuous $T$.

We can access the CATE between $T$ and $T^{\prime}$ as:

$$ $\mathbb{E}[Y(T_{t}=T,W_{t})-Y(T_{t}=T^{\prime},W_{t})|W_{t}]=\theta_{0}(W_{t}) (T-T^{\prime})$ $$

This encoding performs best on the the railway dataset, in which treatments are fairly high dimensional and, conditioned on $W_{t}$, we only observe local variation of the treatment around a typical price.
Figure 4 is an example of the estimated values of $\theta_{0}$ for a random time-series. We notice that the effect of prices has little impact at the beginning of the forecast. The effect of price over sales increases as the departure date gets closer. Moreover, we observe that the price has little effect at time step $t=0$ (the day at which the train leaves leaves).

Figure: Figure 4: $\theta_{0}$ values for various time steps for a railway time serie
Refer to caption: extracted/5969790/images/Appendix/encoding_functions/theta_linear.png

### A.2 One-hot encoding

A more intuitive encoding under our categorical treatments is the one-hot encoding. In this case $e_{0}(W_{t})_{i}=P(T_{t}=i|W_{t})$. With $d$ possible treatments, the vector $\theta_{0}(W_{t})$ is:

$$ $\forall i\in[1:k],\theta_{0}(W_{t})_{i}=\mathbb{E}[Y_{t}(T_{t}=i,W_{t})-f_{0}( W_{t})|W_{t}]$ (6) $$

In order to estimate the CATE between $T$ and $T$, we compute the difference between $\theta_{0}(W_{t})_{T}$ and $\theta_{0}(W_{t})_{T^{\prime}}$:

$$ $\mathbb{E}[Y(T_{t}=T,W_{t})-Y(T_{t}=T^{\prime},W_{t})|W_{t}]=\theta_{0}(W_{t}) _{T}-\theta_{0}(W_{t})_{T^{\prime}}$ $$

Thanks to equation 6 and to domain knowledge on the railway dataset, we can make the assumption that the values of $\theta_{0}(W_{t})_{i}$ must decrease when $i$ increases. Figure 5 shows an example of $\theta$’s vector values for a specific time-series and time step. We notice that the values seem to decrease as the treatment increases. Nonetheless, we observe an increase for small treatment values. It is likely caused by the fact that those treatments are extremely rare in our dataset (conditioned on $W_{t}$).

Figure: Figure 5: Example of $\theta$ vector values for the one-hot encoding on the passenger rail dataset.
Refer to caption: extracted/5969790/images/Appendix/encoding_functions/theta_one_hot.png

### A.3 Cumulative encoding

Figure: Figure 6: Example of $\theta_{0}$ vector values for the cumulative encoding
Refer to caption: extracted/5969790/images/Appendix/encoding_functions/theta_cumul.png

The cumulative encoding function $x\mapsto\mathbb{1}_{T_{t}\geq x}$ is another possible encoding function of the treatment in our data model in Eq. 3. The interpretation of the vector $\theta_{0}(W_{t})$ is then:

$$ $\forall i\in[2:k],\theta_{0}(W_{t})_{i}=\mathbb{E}[Y_{t}(T_{t}=i,W_{t})-Y_{t}( T_{t}=i-1,W_{t})|W_{t}]$ $$

Thus the CATE between treatment $T$ and $T^{\prime}$, with $T^{\prime}<T$ equals:

$$ $\mathbb{E}[Y(T_{t}=T,W_{t})-Y(T_{t}=T^{\prime},W_{t})|W_{t}]=\sum_{i=T^{\prime }+1}^{T}\theta_{0}(W_{t})_{i}$ $$

As done in the last subsection, we also can infer some characteristics of the theta vector values thanks to domain knowledge and equation A.3. In that case, the values of $\theta_{0}$ should be negative for all indexes. Figure 6 shows an example of $\theta(W_{t})$. Most of them are indeed negative. Nonetheless, some values associated to small treatment values, and a few very large treatment values, are still positive. Again, this is likely caused by those treatments being rarely seen during training (conditioned on $W_{t}$).

## Appendix B Detailed Theory and Proofs

### B.1 Stating and Verifying Orthogonal Learning Assumptions

Proposition 3 relies on two assumptions that must be verified by the R-loss to ensure convergence through Theorem 2 in Foster & Syrgkanis ().
These two assumptions are Universal Orthogonality and Boundness.
Denote the observations $Z^{n}=(Y_{1:t}^{n},T_{1:t}^{n},X_{1:t}^{n},S^{n})$, sampled from an unknown distribution $\mathcal{D}_{t}$ depending on $t$. We denote $\theta^{*}\in\arg\min_{\theta\in\Theta}L_{\mathcal{D}_{t}}(\theta,m,e)$, and $star(\Theta,\theta)$ the star domain of $\Theta$ at point $\theta$.
Both assumptions require the definition of a Directional Derivative:

###### Definition 1 (Directional Derivative) .

Let $F:\mathcal{F}\rightarrow\mathbb{R}$ be a function from a vector space $\mathcal{F}$, and define the derivative operator at point $f$ with respect to vector $g$ as:

$$ $D_{f}F(f)[g]=\frac{d}{dt}F(f+tg)_{|t=0}$ $$

With this definition, we can state the two assumptions required of the (population level) R-loss (we omit the indexing on $t$ in $\mathcal{D}_{t}$, as the same assumption must old in each $t$ but our definition is identical at each $t$):

###### Assumption 4 (Universal Orthogonality) .

For all $\hat{\theta}\in star(\hat{\Theta},\theta^{*})+star(\hat{\Theta}-\theta^{*})$ :

$$ $D_{g}D_{\theta}L_{D}(\hat{\theta},m,e)[\theta-\theta^{*},g-g_{0}]=0$ $$

###### Assumption 5 (Boundness) .

$D_{g}D_{g}L_{D}$ and $D_{\theta}D_{\theta}D_{g}L_{D}$ are both continuous and there is a constant $\beta$ such as $\forall\theta\in star(\hat{\Theta},\theta^{*})$ and $\forall\tilde{g}\in star(\mathcal{G},g_{0})$:

$$ $|D_{g}D_{g}L_{D}(\theta,\tilde{g})[g-g_{0},g-g_{0}]|\leq\beta||g-g_{0}||_{ \mathcal{G}}^{2}$ $$

###### Definition 1 (Directional Derivative) .

###### Assumption 4 (Universal Orthogonality) .

###### Assumption 5 (Boundness) .

### B.2 Universal Orthogonality

We denote $W_{t}=(S,X_{1:t},Y_{1:\tau})$. We omit the indices $t$ for $T_{t}$, $Y_{t}$ and $W_{t}$. We note $g_{0}=(m_{0},e_{0})$.
Next, we prove the orthogonality of the R-Loss.

$$ $\displaystyle D_{\theta}L_{D}(\bar{\theta},g_{0})[\theta-\theta^{\star}]$ $\displaystyle=-2\mathbb{E}[(T-e_{0}(W)^{T}\big{(}\theta(W)-\theta^{\star}(W) \big{)}\times\big{(}Y-m_{0}(W)-(T-e_{0}(W))^{T}\bar{\theta}(W)\big{)}]$ $$

We first consider the directional derive with respect to $e$:

$$ $\displaystyle D_{e}D_{\theta}L_{D}(\bar{\theta},g_{0})[\theta-\theta^{\star},e -e_{0}]=$ $\displaystyle 2\mathbb{E}[\big{(}Y-m_{0}(W)-(T-e_{0}(W))^{T}\bar{\theta}(W) \big{)}\times\big{(}e(W)-e_{0}(W)\big{)}^{T}\big{(}\theta(W)-\theta^{\star}(W) \big{)}]$ $\displaystyle-2\mathbb{E}[\big{(}e(W)-e_{0}(W)\big{)}^{T}\bar{\theta}(W)\times \big{(}T-e_{0}(W)\big{)}^{T}\big{(}\theta(W)-\theta^{\star}(W)\big{)}]$ $$

We notice that in the first term of this derivative, we have $\forall w\in\mathcal{W}$:

$$ $\displaystyle\mathbb{E}[((Y-m_{0}(W)-(T-e_{0}(W))^{T}\bar{\theta}(W)|W=w]= \mathbb{E}[\epsilon_{1}+\epsilon_{2}^{T}\big{(}\theta_{0}(W)-\bar{\theta}(W) \big{)}|W=w]=0$ $$

Regarding the second term, when developing the dot products we also have that $\forall w\in\mathcal{W}$:

$$ $\displaystyle\mathbb{E}[(e(W)-e_{0}(W))^{T}\bar{\theta}(W)\times(T-e_{0}(W))^{ T}(\theta(W)-\theta^{\star}(W))|W=w]$ $\displaystyle=\sum_{i=1}^{k}\sum_{j=1}^{k}\mathbb{E}[\epsilon_{2j}|W=w]\times \bar{\theta}(w)_{j}\times(e(w)-e_{0}(w))_{i}\times(\theta(w)-\theta_{0}(w))_{j}$ $\displaystyle=0$ $$

We now need to establish the orthogonality for the parameter $m$ of the population risk:

$$ $\displaystyle D_{m}D_{\theta}L_{D}(\bar{\theta},g_{0})[\theta-\theta^{\star},m -m_{0}]$ $\displaystyle=2\mathbb{E}[(T-e_{0}(W))^{T}(\theta(W)-\theta^{\star}(W))\times( m(W)-m_{0}(W))]$ $\displaystyle=2\mathbb{E}[\epsilon_{2}^{T}(\theta(W)-\theta^{\star}(W))\times( m(W)-m_{0}(W))]$ $\displaystyle=0$ $$

with $E[\epsilon_{2}|Z]=0_{d}$.

Hence, the population level version of the R-loss from Eq. 5 verifies assumption 4.

### B.3 Boundness

The second directional derivatives of the population risk are continuous as it is the composition the square function and a linear function, both with continuous gradients.
We now prove the boundedness of the second derivative of the population risk with regard to the nuisance parameters, for the vector case.

Let us denote $l(\gamma,\zeta,z)=[Y-\gamma_{1}-(T-\gamma_{2})^{T}\zeta]^{2}$ with $\gamma_{1}\in\mathbb{R}$ and $(\gamma_{2},\zeta)\in(\mathbb{R^{d}})^{2}$
We have:

i$isequaltothefirstlinemultipliedby$-ζ_i$.Moreovertheuniquenon-zeroeigenvalueofthismatrixis$2*(1+∑_i=1^d ζ_i^2)$.Itseigenvectoris$(1,-ζ_1,…,-ζ_d)$.\par Denote$M$theboundsuchas$∀W, ∀i, —θ(W)_i—¡M$.Thenwhen$ζ_i = θ(W)_i$,wehavethat$2*(1+∑_i=1^d ζ_i^2) ≤2(1+d M^2)$.ThedirectionalderivativeoftheR-lossfunctioninEq.\ref{eq:r-loss-categ}at$g-g_0$isthederivativeofthecompositionofthefunctions\ref{eq:r-loss-categ}and$~g+t(g-g_0)$,withrespectto$t$.Hence,whentakingtwoderivativeswithrespectto$t$,theequationbecome:\begin{aligned} ||D_{g}D_{g}R(\theta,\tilde{g},z)[g-g_{0},g
-g_{0}]||=||(g-g_{0})^{T}\nabla_{\gamma\gamma}l(\gamma,\zeta,z)(g-g_{0})||\end
{aligned}asthehessianof$t ↦t(g-g_0)$equals0.\par Usingthefactthatthelargesteigenvalueof$∇_γγl(γ,ζ,z)$islessthan$2(1+d M^2)$,wehave$——D_gD_g R(θ,~g)[g-g_0,g-g_0] —— ¡ 2(1+d M^2) ——g-g_0——_G^2$.Finally,
wecomputetheexpectationonthisinequalityandweapplythedominatedconvergencetheoremtogetthepopulationrisk
(requiringthat$——g-g_0——_L_2 ¡ ∞$).Toconclude,wegetthefollowinginequality:\begin{gathered}||D_{g}D_{g}L_{D}(
\theta,\tilde{g})[g-g_{0},g-g_{0}]||<2(1+d\times M^{2})||g-g_{0}||_{L_{2}}^{2}
\end{gathered}Thus,assumption\ref{eq:Boundness}isverifiedwith$β= 2(1+d M^2).

## Appendix C Implementation details

### C.1 Causal Transformer

We use the Causal Transformer (Melnychuk et al., ) in order to get a reference model for comparison. We use the public github repository ([https://github.com/Valentyn1997/CausalTransformer/tree/main](https://github.com/Valentyn1997/CausalTransformer/tree/main)). Melnychuk et al. () train an auto regressive model, like us using 25 vital signals, the age, the gender and the ethnicity to estimate the blood pressure of patient.
In that context, the relevant features for predicting blood pressure at time $t$ are all the vitals before time $t$, static features and all treatments up until time $t$. Nonetheless, when evaluating for time $\tau+l$ for $l$ in $[1,5]$, the model had knowledge of vitals at time $\tau+1$, enabling the model to perform very well for $l=1$ and leaking information for all $\tau$ ([https://github.com/Valentyn1997/CausalTransformer/blob/c49a24faa57af966501e241f57a26b528f874a53/src/data/mimic_iii/real_dataset.py#L67](https://github.com/Valentyn1997/CausalTransformer/blob/c49a24faa57af966501e241f57a26b528f874a53/src/data/mimic_iii/real_dataset.py#L67)). This leak manifested as extremely good performance for $l=1$, followed by a large drop of performances at $l\geq 2$. In our evaluation, we fix the leak from the future, and obtain results consistent with other models trained on this dataset during our study.

### C.2 RDD

Algorithm 2 describes how to compute the dataset of RDD values in full details.

Figure: Algorithm 2 RDD Algorithm

### C.3 TFT

We developed a custom version of the Temporal Fusion Transformer Lim et al. (). We focused on this architecture as it supports static features, past and futures features as input, and temporally causal processing of the input, where most architectures focus on predicting future time step from the previous ones Liu et al. (),Zhang & Yan (),Wu et al. (),Nie et al. ().
We apply several changes to the original TFT:

- •
We remove the encoder of the TFT and replace it by a series of encoder blocks. Figure 7(a) presents the architecture of our encoder block. Each encoder block is made of self-attention and residual blocks. Each block divides the length of the input by 2. Hence, at the beginning of the encoder, the model processes high frequency features, and the deeper the encoder is the lower frequency the model is considering.
- •
We use conditional mechanisms that show interesting performances on Image Generation. Static features go through a dense model to create an embedding $s$ gathering all the static information. We then apply successive Cross Attention mechanisms between the token $s$ and the tokens of the temporal features Rombach et al. () to mix temporal and static features. Finally, we apply adaptive group normalization Dhariwal & Nichol () where $s$ is the class embedding.
- •
We replace the masked self-attention decoder by a new decoder block. Figure 7(b) shows the layers used in the decoder block. This block is almost the symmetric of the encoder one. We remove the cross-attention layer in the decoder as it increases the expressivity of the model. We also add skip connections between the encoder the decoder blocks as done in Ronneberger et al. (). Thus, the information of high frequency features are directly transferred to the decoder without passing thought the entire encoder in order to help the expressivity of the model.
- •
We use causal convolution layers to prevent any possible leak of information from future features $X_{t+i}$ to impact $Y_{t}$ prediction. We also apply masking in our self-attention layers as done in Li et al. ().

Figure: (a) Encoder block
Refer to caption: extracted/5969790/images/Appendix/TFT/Encoder.jpg

## Appendix D Additional results

### D.1 Addition MIMIC-II results

Table 3 shows MAE results for in-distribution forecasts. Table 4 shows MAE results for causal effects using our RDD dataset. Results show the same effects as those described in §4.2 for the RMSE and RDD RMSE.

**Table 3: MAE with respect to the ground truth time series, best values are in bold**
| Models | $t$ = $\tau$+1 | $t$ = $\tau$+2 | $t$ = $\tau$+3 | $t$ = $\tau$+4 | $t$ = $\tau$+5 |
| --- | --- | --- | --- | --- | --- |
| Causal Transformer | 6.028 $\pm$ 0.027 | 6.703 $\pm$ 0.027 | 7.048 $\pm$ 0.025 | 7.290 $\pm$ 0.032 | 7.485 $\pm$ 0.040 |
| TFT Baseline | 6.162 $\pm$ 0.044 | 6.745 $\pm$ 0.024 | 7.032 $\pm$ 0.018 | 7.233 $\pm$ 0.022 | 7.402 $\pm$ 0.024 |
| Causal TFT One-hot | 6.170 $\pm$ 0.030 | 6.802 $\pm$ 0.033 | 7.115 $\pm$ 0.033 | 7.337 $\pm$ 0.036 | 7.515 $\pm$ 0.033 |
| Causal TFT Cumulative | 6.147 $\pm$ 0.045 | 6.799 $\pm$ 0.054 | 7.125 $\pm$ 0.044 | 7.359 $\pm$ 0.05 | 7.544 $\pm$ 0.051 |

**Table 4: RDD MAE w.r.t the CATE, best values are in bold**
| Models | $t$ = $\tau$+1 | $t$ = $\tau$+2 | $t$ = $\tau$+3 | $t$ = $\tau$+4 | $t$ = $\tau$+5 |
| --- | --- | --- | --- | --- | --- |
| Causal Transformer | 2.286 $\pm$ 0.096 | 2.271 $\pm$ 0.097 | 2.249 $\pm$ 0.096 | 2.241 $\pm$ 0.094 | 2.242 $\pm$ 0.085 |
| TFT Baseline | 2.301 $\pm$ 0.065 | 2.320 $\pm$ 0.067 | 2.337 $\pm$ 0.071 | 2.368 $\pm$ 0.065 | 2.400 $\pm$ 0.058 |
| Causal TFT One-hot | 2.265 $\pm$ 0.066 | 2.265 $\pm$ 0.066 | 2.265 $\pm$ 0.066 | 2.265 $\pm$ 0.066 | 2.265 $\pm$ 0.066 |
| Causal TFT Cumulative | 2.292 $\pm$ 0.044 | 2.292 $\pm$ 0.044 | 2.292 $\pm$ 0.044 | 2.292 $\pm$ 0.044 | 2.292 $\pm$ 0.044 |

## Appendix E Extended Related Work

Our two background sections §2.1 and §3.1 already discuss the closest related work. In this section, we expand on this discussion and discuss other related papers.

Time-series forecasting has seen a flurry of recent activity. The main task studied in the literature is that of multivariate time series forecasting, with the goal of forecasting $(Y^{n}_{i})_{i\in[t:t+h]}$ from $(Y^{n}_{i})_{i\in[0:t]}$. Another common task is to forecast $(Y^{n}_{i})_{i\in[1:t]}$ from $(X^{n}_{i})_{i\in[1:t]}$. Various approaches have been proposed, though not clear winner has emerged yet. TSMixer Chen et al. () relies only on fully connected linear layers to process entire time-series feature-wise, before mixing features at each time-step. The iTransformer Liu et al. () architecture similarly leverages fully connected layers, but processes intermediate representations with transformer layers channel wise. PatchTST Nie et al. () gathers patches time-steps and features that it treats as tokens in a self-attention mechanism, before reconstructing a time-series shaped output. In all these architectures and others Ilbert et al. (); Wu et al. (), there is no insurance against future information leaking into previous time-steps. While this can be fine in some cases, this is not compatible with our specific task.
Indeed, our causal models predict $Y^{n}_{t+l}$ based on past features $(Y^{n}_{i})_{i\in[0:t]}$, static features $S^{n}$ and temporal features $(X^{n}_{i})_{i\in[0:t+l]}$. A key point in our task is to ensure causality in the forecasting, preventing any leak of information from $(X^{n}_{i})_{i\in[t+1:t+h]}$. Supporting this requirement is not common in recent time series papers, and the reason why we build on the TFT Lim et al. (), which does provide an adequate structure.

Causal time-series models have seen less progress, though some time-series models aim to estimate the CATE.
Bica et al. (); Lim () rely on the same meta-architecture, including a network computing the treatment propensity score and a network computing the outcome (though Bica et al. () tries to learn unobserved confounders, a challenging task). Both models use LSTM or RNN architectures merged with new mechanisms: a propensity network learns to predict propensity scores, which are then used to train the outcome network with the inverse probability of treatment weighting. This loss, mainly used for policy learning Foster & Syrgkanis (), can be high variance and does not full corrects for regularization bias. This is why we decided to train our model with the Residualized Loss (R-Loss) Foster & Syrgkanis (); Nie & Wager (), which also gives us a more flexible approach for treatment encoding.

Brodersen et al. () take a different approach, and use structured time-series models to explicitly model counterfactuals under no-treatment as a synthetic control. This approach is well suited to binary treatments with several observations, to be able to model outcomes without treatment. It is less applicable to our case with high dimensional treatments and complex conditioning variables. However, it could interesting to explore this approach as a multi-time-steps alternative to our RDD causal test set.

With the recent rise of the transformer architecture for sequential models, new techniques have been developed to merge transformers with causal inference Melnychuk et al. (). The authors leverage an adversarial loss Bica et al. () to encourage the model to learn balanced features for the treated and untreated populations and learn the causal link between the treatment and the outcome.
Melnychuk et al. () use a transformer backbone with two small dense layers on top trained with this adversarial loss.
This identification approach differs from classical causal inference frameworks, and it is unclear that it converges to a well defined conditional causal effect.

Regression Discontinuity Designs (RDDs) are studied theoretically, and frequently used empirically, in economics Hahn et al. (); Imbens & Lemieux (); Lee & Lemieux ().
Researchers typically use RDDs to estimate one specific ATE. In contrast, we estimate many ATEs (one for each time series), and combine them with time-series features $W_{t}$ to create a CATE dataset. To the best of our knowledge, this is a novel appraoch to evaluate causal models.
Researchers have used RDDs to estimate ATEs in time-series before though. Hausman & Rapson () surveys typical approaches and pitfalls, and proposes empirical checks to perform on the data to verify if RDDs are applicable. In our setting, we typically do not have enough data around switching times to perform such checks, though it would be interesting to take inspiration from common RDD checks to improve our causal test sets. We leave this direction for future work.

## Peer Reviews

confidence: 4
contribution: 3
initial_rating: 6
presentation: 3
questions: My suggestions would include:
- improving the presentation (especially citations to references and equations)
- improving the way the methodology is presented, so that it does not need to rely on discussing a specific application
- be clear about why we care about improving RDD RMSE - especially here, the TST baseline in table 1 (which does not rely on any causal inference approach) seems to be as competitive as the bespoke causal approaches...
soundness: 3
strengths: It is great to see that the authors focus on both theoretical results and using empirical investigations (for 2 different cases) to explore their approach. It is mostly an extension to existing ideas and methods (e.g., orthogonal learning, for time-series), but still in a way that makes it very valuable since time series learning and forecasting is a very broad field for which any improvement can make a high impact.
summary: The paper concentrates on using orthogonal learning for causal inference and forecasting, for the case of time-series. The paper is well written overall, combining theoretical results (mostly relying on the work of others, e.g., as in the proofs) and empirical investigations. It provides a set of interesting results related to both learning and evaluating forecasting approaches in a causal framework. I find that the paper makes some novel and potentially valuable contribution.
weaknesses: Even though the paper is well written, it is also very compact overall, and difficult to read at stages. Already the motivations could be better developed, as well as the justification for orthogonal learning. Personally, I felt quite frustrated with the readability of the paper. Besides quite a number of points that could be improved in terms of presentation (citation style for references and equations), the flow of the paper is difficult to follow at stages. This is possibly due to page limitations, but still, the authors could have chosen to lighten the introductory parts, in order to leave more room for the novel developments in the paper. Also in terms of the presentation, some of the methodological parts are written as if we only consider a specific application (e.g., mentioning daily prices as treatment at the top of page 5). However, different applications are presented in the empirical investigation part. Maybe there could be some additional effort made with the writing, to make the paper easier to read and to better separate purely methodological considerations and developments from application-specific considerations.

---

confidence: 3
contribution: 3
initial_rating: 5
presentation: 3
questions: 1. How does the causal TFT compare to other models under the original simulation framework used to evaluate Melnychuk 2020, Bica 2020 and Lim 2018?  While more simplistic, it makes it easier to control the degree of bias present in policies, which could widen performance gaps (vs the slim margins seen in the RDD experiments).
2.  Where improvements are seen (i.e. shorter time steps in RDD expts), are these attributable to the use of TFT or new loss function proposed?
soundness: 2
strengths: While many architectures have been proposed for time series forecasting in general, causal temporal models have been relatively understudied by comparison. The paper is also well motivated and clearly presented — with the use of orthogonal learning being an interesting approach, and use of the TFT being sensible to incorporate a wider range of inputs. The outline of the RDD evaluation approach is also useful, providing a way to evaluate causal temporal models on live observational data.
summary: The paper proposes an alternative loss function for training causal deep learning models to estimate conditional average treatment effects (CATE) over time - through the use of orthogonal learning to directly learn the CATE from data under unconfoundedness assumptions. Losses are applied over a customised TFT architecture for time series predictions, and they demonstrate improvements using Regression Discontinuity Design (RDD) experiments.
weaknesses: However, the performance improvements do appear slim, with the causal TFT underperforming the causal transformer on both on all in distribution time-steps and 2/5 of the RDD forecast horizons. Where improvements exist, it is also not immediately clear if they are attributable to the improved loss function or the use of the TFT — with the novelty of the former far out-weighing the latter.

---

confidence: 3
contribution: 3
initial_rating: 6
presentation: 2
questions: Why are more datasets not considered for evaluation?
soundness: 3
strengths: 1. This work has real-world applications and has the potential to improve decision-making processes that depend on forecasts from time series models.
2. The two extensions to adapt the orthogonal statistical learning framework to train causal models—defining daily treatment effects and extending the binary treatment effect data model to categorical and linear treatments—are mathematically justified.
3. Despite the complex approach and new terminologies introduced, the paper is well-structured and relatively easy to follow.
summary: This paper highlights the need for time series models to showcase generalization to out-of-distribution data, where the relationship between treatment and outcome differs from the training data. To this end, the authors propose extending the orthogonal statistical learning framework to train causal forecasting models that can capture changes in outcome based on changes in treatment outside of the training distribution. Due to a lack of out-of-distribution evaluation of existing time series models, the authors propose to evaluate out-of-distribution by creating a test set of causal treatment effects using Regression Discontinuity Designs (RDD).
weaknesses: 1. T_t for treatment at time t, and T_1, T_2 ​ for treatment types can be confusing.
2. The evaluation method using RDD is complex and may generate only a few test set examples.

---

confidence: 2
contribution: 2
initial_rating: 5
presentation: 2
questions: I have outlined some potential confusions and questions in my "Weaknesses" section response.
soundness: 3
strengths: The authors pursue an interesting synthesis of methods from deep learning for time series and statistic tools for treatment effect prediction.  Certainly there are many cases in time series where there might be nuisance variables that lead to unexpected behavior out of distribution, and in cases where there is domain knowledge about these nuisance variables, this method might be helpful.
summary: The authors argue that current deep learning approaches to time series modeling are limited by not explicitly modeling causality. To remedy the situation, the authors propose a framework to convert neural network time series models into causal models by fitting a conditional average treatment effect (CATE) model on a temporally held-out dataset. The authors demonstrate how to use orthogonal learning to construct the CATE model and model high-dimensional treatment using 3 different encodings (one-hot, cumulative, and linear). To evaluate the models, the authors construct ground truth CATEs using regression discontinuity designs (RDDs), which use a temporal cutoff to create examples with and without a specified treatment. The authors then show that while non-causal models perform well in terms of RMSE and MAE in distribution, they perform poorly on recovering the CATEs constructed using RDD and domain knowledge. The authors then conclude that their causal models are superior in settings with clear causal relationships.
weaknesses: Overall I found the paper pretty challenging to parse, for example how exactly the models were constructed and evaluated. I am not a domain expert in casual methods, but I work in machine learning and time series modeling, and I don't think the material was so technical that its complexity was the source of the problem, so sections 2-4 might benefit from some reworking to make the key details more salient.

Beyond style, I take two of the central claims of the paper to be:

(1) Adding causal modeling to neural network time series models helps in practical out-of-distribution scenarios

(2) Estimating CATEs in a time-series setting is novel and necessitates methodological innovation

I don't have enough domain knowledge to fully assess (2), but I am currently unconvinced that the paper demonstrates (1). The evaluation of the model seems a little circular in that the examples are special designed using a similar procedure to the method being evaluated. Perhaps I'm misunderstanding, but I take the evaluation to be whether the causal predictions of the model match a linear RDD analysis. If this is a misunderstanding, a correction might be helpful. If I've understood correctly, when would we be better served by apply the methodology in the paper than just applying the linear RDD analysis directly? This evaluation also seems limited to cases where the causal structure of the problem is known and the changing covariate (e.g. price) is observed. How does this solve prediction in out-of-distribution settings? There is also no explicit evaluation on out-of-distribution data using the method, even though this was one of the stated motivations in the introduction.

For (2), the contribution appears to be mostly the encoding of the treatment and the treatment effect \theta(W_t), because the other tools seem to be borrowed directly from related work. There are some additional proof that those tools can be extended to this setting, but as far as I can tell, the crux of the method is combining a pre-trained time series backbone with existing statistical. Please correct me if this is a misunderstanding. Further from the evaluations that were presented in the results, the method's performance is nearly identical to the baseline for several of the encoding, in terms of the causal evaluation, and its general forecasting performance is worse. What is the practical take-away from this evaluation? What are the settings in which I would want to adopt this method as an alternative to simpler deep learning forecasting methods given that it performs worse in-distribution and there is no explicit out-of-distribution evaluation? In terms of predicting CATE, how can this be used beyond recapitulating the known effects of known treatments?

Meta-review: This paper introduces orthogonal learning for causal time-series forecasting, aiming to assess how interventions impact out-of-distribution outcomes. Reviewers appreciated the integration of theoretical concepts with empirical evaluations on real-world tasks like demand forecasting and healthcare, as well as the innovative use of Regression Discontinuity Designs for evaluating causal effects. However, initial concerns were raised regarding the clarity and presentation, particularly the connection between theory and application, and the limited number of evaluation points due to dependence on specific treatment changes. While the authors’ rebuttal addressed most concerns, the work’s overall significance and novelty remain uncertain for part of the review panel.

Decision: Reject