Imke Mayer
email: imke.mayer@ehess.fr

Julie Josse
email: julie.josse@inria.fr

Traumabase Group

Generalizing treatment effects with incomplete covariates

Keywords: Causal effect transportability, missing values, external validity, data integration, multiple imputation, evidence-based medicine

We focus on the problem of generalizing a causal effect estimated on a randomized controlled trial (RCT) to a target population described by a set of covariates from observational data. Available methods such as inverse propensity weighting are not designed to handle missing values, which are however common in both data sources. In addition to coupling the assumptions for causal effect identifiability and for the missing values mechanism and to defining appropriate estimation strategies, one difficulty is to consider the specific structure of the data with two sources and treatment and outcome only available in the RCT. We propose and compare three multiple imputation strategies (separate imputation, joint imputation with fixed effect, joint imputation without source information), as well as a technique that uses estimators that can handle missing values directly without imputing them. These methods are assessed in an extensive simulation study, showing the empirical superiority of fixed effect multiple imputation followed with any complete data generalizing estimators. This work is motivated by the analysis of a large registry of over 20,000 major trauma patients and a RCT studying the effect of tranexamic acid administration on mortality. The analysis illustrates how the missing values handling can impact the conclusion about the effect generalized from the RCT to the target population.

Introduction

Observational and clinical trial data can provide different perspectives when evaluating an intervention or a medical treatment. Combining the information gathered from experimental and observational data is a promising avenue for medical research, because the knowledge that can be acquired from integrative analyses would not be possible from any single-source analysis alone. Such integrative analyses can be used for example to predict a treatment effect on a specific target population using the one estimated from the RCT, to validate methods, especially applied to observational data by emulating a trial, to better estimate heterogeneous effects (which generally cannot be estimated from experimental data due to underpowered studies). Here, we are interested in the former case, where the experimental data or randomized controlled trial (RCT) is considered as a biased sample of a target population and we would like to estimate the treatment effect on the target population represented by an observational study. More precisely, the effect is estimated on a RCT composed of covariates, treatment and outcome while the observational study is composed only of covariates. There exists a multitude of methods to generalize a treatment effect, for detailed reviews we refer to [START_REF] Colnet | Causal inference methods for combining randomized trials and observational studies: a review[END_REF] and [START_REF] Degtiar | A review of generalizability and transportability[END_REF]. But all these methods do not consider the problem of missing data which is ubiquitous in data analysis practice [START_REF] Josse | Introduction to the special section on missing data[END_REF]. We emphasize here that we focus on incomplete covariates in both studies, while the outcome and treatment (available in the RCT) are assumed to be fully observed. In certain cases, naive approaches such as complete-case analysis can yield unbiased treatment effect estimates [START_REF] Bartlett | Asymptotically unbiased estimation of exposure odds ratios in complete records logistic regression[END_REF]; however, in many settings, especially for observational data, the estimations are known to be biased since the complete-case observations are generally not a representative subsample of the population of interest [START_REF] Little | Statistical Analysis with Missing Data[END_REF]. Furthermore, this approach is sometimes not even possible since in high-dimensional settings the probability of having complete observations decreases rapidly [START_REF] Zhu | High-dimensional principal component analysis with heterogeneous missingness[END_REF]. There exists a multitude of methods to handle missing values [START_REF] Little | Statistical Analysis with Missing Data[END_REF][START_REF] Van Buuren | Flexible Imputation of Missing Data[END_REF][START_REF] Mayer | R-miss-tastic: a unified platform for missing values methods and workflows[END_REF], such as maximum likelihood estimation or multiple imputations with the aim of estimating as well as possible a parameter and its variance. These methods make assumptions on the mechanism that generated the missing values. More recent works also consider the question of supervised learning with missing values which is a different problem from statistical inference of model parameters [START_REF] Josse | On the consistency of supervised learning with missing values[END_REF][START_REF] Le Morvan | Neumiss networks: differential programming for supervised learning with missing values[END_REF]. In the context of causal inference there exist several recent works that address missing data [START_REF] Mattei | Estimating and using propensity score in presence of missing background data: an application to assess the impact of childbearing on wellbeing[END_REF][START_REF] Seaman | Inverse probability weighting with missing predictors of treatment assignment or missingness[END_REF][START_REF] Yang | Causal inference with confounders missing not at random[END_REF][START_REF] Kallus | Causal inference with noisy and missing covariates via matrix factorization[END_REF][START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF]. One difficulty in the context of causal inference, is that one has to couple identifiability assumptions for the causal parameter with assumptions on the missing values mechanisms to define an appropriate strategy. According to [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF], identifiability of the causal effect with missing values is ensured by adapting the causal inference assumptions to the missing values setting with an unconfoundedness despite missing values (UDM) assumption [START_REF] Rosenbaum | Reducing bias in observational studies using subclassification on the propensity score[END_REF]. However, these works only consider the case of a single dataset-or potentially multiple datasets with the same data distribution, i.e., sampled from the same population of interest-and do not treat the case of transporting or generalizing a treatment effect from an RCT to a target distribution defined through an observational dataset; the RCT representing a "distorted" population due to sampling or selection bias (generally defined via eligibility criteria). In practice, the observed distributions between the observational data and the RCT do not only differ due to the selection bias but may also differ in terms of missing values patterns. Another field that explicitly studies data integration and missing values therein is meta-analysis. [START_REF] Burgess | Combining multiple imputation and meta-analysis with individual participant data[END_REF] study a multiple imputation approach in the context of meta-analysis where missing values can occur in different data sources.

In this paper we propose to revisit the standard identifiability assumptions and estimators for generalizing a treatment effect from one to the other under the perspective of identifying the impact of missing values and suggesting appropriate strategies to handle missing values according to the different assumptions, such as the missing values mechanisms. Indeed, to the best of our knowledge, there exists no directly applicable method for this setting nor guidelines for a correct handling of missing values in either of or both data sources. Our main contributions are:

• we define several multiple imputation strategies adapted for the aforementioned data integration problem in Section 3;

• we propose alternative identifiability assumptions which can be seen as an an extension of the unconfoundedness despite missingness (UDM) assumption from the observational data case [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF] and suggest adapted estimators in Section 4;

• we assess the performance of the proposed estimators and naive complete case estimators in an extensive simulation study in Section 5;

• we present the results of the real-world data analysis which has motivated this work in Section 6.

For simplicity, we will assume that the relevant covariates that allow for adjustment of the sampling bias are observed both in the RCT and the observational dataset, and we focus on the case of different missing values mechanisms for these covariates. Before diving into the statistical aspects of the problem of treatment effect generalization, we briefly introduce a medical question about major trauma patients that has motivated this work.

Major trauma denotes injuries that endanger the life or the functional integrity of a person. The World Health Organization (WHO) has recently shown that major trauma including roadtraffic accidents, interpersonal violence, falls, etc. remains a world-wide public health challenge and major source of mortality and handicap [START_REF] Leigh | Global,regional,and national disability-adjusted life years(dalys) for 359 diseases and injuries and healthy life expectancy(hale) for 195 countries and territories,1990-2017: a systematic analysis for the global burden of disease study 2017[END_REF]. We focus on trauma patients suffering from a traumatic brain injury (TBI). TBI is a sudden damage to the brain caused by a blow or jolt to the head and can lead to intracranial bleeding that can be observed on a computed tomography (CT) scan. Ongoing intracranial bleeding can lead to raised intracranial pressure, brain herniation, and death. Over 10 million people pass away or are hospitalized worldwide because of TBI each year [START_REF] Dewan | CRASH-3: Tranexamic acid for the treatment of significant traumatic brain injury: study protocol for an international randomized, double-blind, placebo-controlled trial[END_REF]. Tranexamic acid (TXA) is an antifibrinolytic agent that limits excessive bleeding, commonly given to surgical patients. Previous clinical trials showed that TXA decreases mortality in patients with traumatic extracranial bleeding [START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF]. Such a result raises the possibility that it might also be effective in TBI, because intracranial hemorrhage is common in TBI patients. Therefore the question here is to assess the potential decrease of mortality in patients with intracranial bleeding when using TXA.

To answer this question, we have at disposal two data sources:"CRASH-2", a multi-center international RCT, and "Traumabase", an observational national registry. The details about these data are provided in Section 6, as well as the analysis results for generalizing the treatment effect from the CRASH-2 study to the Traumabase registry.

Background and notations

We first introduce the notations, standard assumptions and estimators in the complete case, following the line of [START_REF] Colnet | Causal inference methods for combining randomized trials and observational studies: a review[END_REF].

Notations

The general case consists in assuming that each patient from the two populations is fully characterized by a random tuple (X, Y (0), Y (1), A, S), where X ∈ X is a p-dimensional vector of covariates, A denotes the binary treatment assignment, Y (a) is the potential outcome for treatment level a, a ∈ {0, 1} and S indicates RCT participation1 . The data considered is composed of n + m independent random tuples: (X i , Y i (0), Y i (1), A i , S i ) i=1,...,n+m . For simplicity of exposition, we introduce an additional indicator variable Q that allows to distinguish between observations from the RCT (indexed by the set Set R {1, . . . , n}) and from the observational cohort (indexed by

Set O {n + 1, . . . , n + m}). Thus Q i R if i ∈ Set R O if i ∈ Set O . Note that (i) P (Q = R|S = 1) = 1
, and we assume that conditionally on S = 0, Q is independent of all other variables considered in this setting. In this paper, we consider the case where for each RCT sample i such that

Q i = R, we observe (X i , A i , Y i , S i = 1
), while for observational data i such that Q i = O, we consider that we only observe the covariates X i , more precisely:

• the RCT samples i = 1, . . . , n are identically distributed according to P(X, Y (0), Y (1), A, S | S = 1),

• and the observational data samples i = n+1, . . . , n+m are identically distributed following P(X, S).

An illustration of the data we consider in this work is provided in Figure 1. We define the conditional average treatment effect (CATE):

Q Covariates Treatment Outcome under A=0 Outcome under A=1 X1 X2 X3 A Y (0) Y (1) 1 R 1.1 20 5.4 1 23.4 24.1 . . . R . . . . . . . . . . . . n -1 R -6 45 8.3 0 26.3 27.6 n R 0 15 6.2 1 28.1 23.5 n + 1 O -2 52 7.1 NA NA NA n + 2 O -1 35 2.4 NA NA NA . . . O . . . NA NA NA n + m O -2 22 3.4 NA NA NA Q Covariates Treatment Outcome under A X1 X2 X3 A Y 1 R 1.1 20 5.4 1 24.1 . . . R . . . . . . . . . n -1 R -6 45 8.3 0 26.3 n R 0 15 6.2 1 23.5 n + 1 O -2 52 7.1 NA NA n + 2 O -1 35 2.4 NA NA . . . O . . . NA NA n + m O -2 22 3.4 NA NA
∀x ∈ X , τ (x) = E[Y (1) -Y (0) | X = x] , (1)
and the population average treatment effect (ATE):

τ = E[Y (1) -Y (0)] = E[τ (X)] ,
while we define the RCT (or sample) average treatment effect as

τ 1 = E[Y (1) -Y (0) | S = 1].
We denote by µ a (x) and µ a,1 (x) the conditional response surfaces under treatment a ∈ {0, 1} in the general and in the RCT population, respectively:

µ a (x) = E [Y (a) | X = x] , µ a,1 (x) = E [Y (a) | X = x, S = 1] ,
and by π S (x) the selection score:

π S (x) = P (S = 1 | X = x) .
Note that π S (x) is the probability of being eligible for selection in the RCT and of being willing to participate given covariate values x. It is different from the probability that an individual with covariates x, known to be in the study (RCT or observational population), is selected in the RCT:

π S (x) = π R (x), π R (x) = P (∃i, Q i = R, X i = x | ∃i ∈ Set R ∪ Set O , X i = x) .
We similarly note

π O (x) = P (∃i, Q i = O, X i = x | ∃i ∈ Set R ∪ Set O , X i = x) = 1 -π R (x) .
Finally, we denote by α(x) the conditional odds that an individual with covariates x is in the RCT or in the observational cohort:

2 α(x) = P (i ∈ R | ∃i ∈ Set R ∪ Set O , X i = x) P (i ∈ O | ∃i ∈ Set R ∪ Set O , X i = x) = π R (x) π O (x) = π R (x) 1 -π R (x)
.

This quantity arises in several approaches that have been proposed to generalize a treatment effect from an RCT to a target population, see for example [START_REF] Westreich | Transportability of trial results using inverse odds of sampling weights[END_REF]. As we will see in the following, this conditional odds is identifiable under certain assumptions and can be used instead of the selection score π S to generalize a treatment effect. Indeed the latter is only identifiable in the case of a nested trial design [START_REF] Dahabreh | Study designs for extending causal inferences from a randomized trial to a target population[END_REF] which we do not cover in this work.

Assumptions for identifiability of the ATE on the target population in the full data case

The main identifiability assumptions that allow for generalizing an ATE from the RCT onto a target population are as follows:

• Internal validity of the RCT:

-Consistency of potential outcomes Y = A Y (1) + (1 -A) Y (0)
-Treatment randomization Y (a) ⊥ ⊥ (A, X) | S = 1 for all X and a = 0, 1.

• Generalizability of the RCT to the target population:

-Ignorability on trial participation

{Y (0), Y (1)} ⊥ ⊥ S | X. (2)
-Positivity of trial participation ∃ c such that for all x, P (π S (x) ≥ c > 0) = 1, and 0 < P (A = a | X = x, S = 1) < 1 for all a and for all x such that P (S = 1 | X = x) > 0.

(3) Assumption (2) implies that we have measured all variables related to the trial eligibility indicator S that are treatment effect modifiers. In other words, participation in the RCT is randomized within levels of X. This is analogous to the ignorability assumption on treatment assignment in causal inference with observational data. With Assumption (3) we require adequate overlap of the covariate distribution between the trial sample and the target population, as well as between the treatment groups in the trial.

Estimators in the full data case

The covariate distribution of the RCT sample is generally different from that of the target population; therefore, τ 1 is different from τ , and an estimator based solely on the RCT is biased for the ATE of interest τ . Under the previous identifiability assumptions, different estimators are available to estimate the ATE τ : one proposes to reweight the RCT sample so that it "resembles" the target population with respect to the observed shifted covariates and treatment effect modifiers.

Another proposes to model the conditional outcomes with and without treatment in the RCT, and then to apply the model to the target population of interest. Doubly robust approaches are combinations of the former two, improving the robustness and efficacy. For simplicity, we assume a standard RCT with a constant propensity score of 0.5 for all individuals.

Inverse probability of sampling weighting (IPSW). This estimator is defined as the weighted difference of average outcomes between the treated and control groups in the trial. The observations are weighted by the inverse odds 1/α(x) = π O (x)/π R (x) to account for the shift of the covariate distribution from the RCT sample to the target population. The IPSW estimator can be written as follows:

τ IP SW = 2 m n i=1 Y i (2A i -1) α(X i ) , (4)
where α is an estimate of α. The IPSW estimator is consistent as soon as α is consistently estimated by α.

Conditional outcome-based estimation. It fits models of the conditional response surfaces among trial participants. Applying these models to the covariates of the observational data, gives the corresponding expected outcome [START_REF] Robins | A new approach to causal inference in mortality studies with a sustained exposure period-application to control of the healthy worker survivor effect[END_REF]. This outcome-model-based estimator, also called g-formula estimator, is then defined as:

τ CO = 1 m n+m i=n+1 ( µ 1,1 (X i ) -µ 0,1 (X i )) , (5)
where µ a,1 (X i ) is an estimator of µ a,1 (X i ). If the model is correctly specified, then the estimator is consistent.

Doubly robust estimators

The selection score and outcome models used in the first two estimators can be combined to form an augmented IPSW estimator (AIPSW):

τ AIPSW = 2 m n i=1 1 α(X i ) [A i {Y i -µ 1,1 (X i )} -(1 -A i ) {Y i -µ 0,1 (X i )}] + 1 m m+n i=n+1 { µ 1,1 (X i ) -µ 0,1 (X i )} .
It is doubly robust, i.e., consistent and asymptotically normal when either one of the two models for π R and µ a,1 (X) (a = 0, 1) is consistent.

Calibration weighting. It is well known that IPSW is likely to be unstable as soon as some of the estimated odds are very small. To resolve the instability of IPSW calibration weighting [START_REF] Dong | Integrative analysis of randomized clinical trials with real world evidence studies[END_REF] have been proposed. They calibrate, i.e. weight subjects in the RCT sample in such a way that afterwards, the covariates are balanced between the RCT sample and the target population: usually the balance is enforced on the first and second moments of the covariates X 1 , . . . , X p such that the weighting mean and variance for each variable in the RCT match the ones in the observational data. More precisely, in order to calibrate, they assign an entropybalancing weight ω i to each subject i in the RCT sample obtained by solving an optimization problem:

min ω1,...,ωn n i=1 ω i log ω i , (6)
subject to ω i ≥ 0, for all i,

n i=1 ω i = 1, n i=1 ω i g(X i ) = g, (the balancing constraint)
where )] from the observational sample. The balancing constraint calibrates the covariate distribution of the RCT sample to the target population in terms of g(X). The objective function in (6) is the negative entropy of the calibration weights; thus, minimizing this criterion ensures that the empirical distribution of calibration weights are not too far away from the uniform, such that it minimizes the variability due to heterogeneous weights. Based on the calibration weights, the CW estimator is then defined as

g = m -1 m+n i=n+1 g(X i ) is a consistent estimator of E[g(X
τ CW = 2 n i=1 ω i Y i (2A i -1). ( 7
)
This estimator is doubly robust in that it is a consistent estimator for τ if either the selection score follows a log-linear model, or if the CATE (1) is linear in the calibration constraint.

Missing values mechanisms

In the taxonomy proposed by [START_REF] Rubin | Inference and missing data[END_REF], an ignorable missingness mechanism is either missing completely at random (MCAR) or missing at random (MAR). The former means that the missingness mechanism is independent of the data, whereas the latter states that the missingness only depends on the observed values. In a nutshell, ignorable missingness means that the missing data mechanism can be "ignored" when doing inference for the parameter in the data likelihood function since the two are separable in the full data likelihood function, whereas non-ignorable missingness complicates analyses more significantly.

More formally, we denote the response pattern of the i-th sample as

R i ∈ {0, 1} p such that R ij = 1 if X ij is observed and R ij = 0 otherwise.
We model 1 -R i as a random vector and its (conditional) distribution is known as the missing values mechanism. Additionally, for our problem, for all response patterns r and X = (X obs(r) , X mis(r) ) the partition of the data in realized observed and missing values given a specific realization of the pattern,

(M CAR) ∀ r, P (R = r|X, S, Q, A, Y ) = P (R = r) (8) (M AR) ∀ r, P (R = r|X, S, Q, A, Y ) = P (R = r|X obs(r) , S, Q, A, Y ) (9)
If the missingness mechanism is non-ignorable, it is qualified as missing not at random (MNAR) and it formally states that the mechanism does not satisfy (8) or ( 9), in other words the missingness is allowed to depend on the missing values themselves.

In Figure 2 we provide an example of observed incomplete data and recall the underlying data for our problem. 3 Multiple imputation

Q Covariates Treatment Outcome under A=0 Outcome under A=1 X1 X2 X3 A Y (0) Y (1) 1 R 1.1 20 5.4 1 23.4 24.1 . . . R . . . . . . . . . . . . n -1 R -6 45 8.3 0 26.3 27.6 n R 0 15 6.2 1 28.1 23.5 n + 1 O -2 52 7.1 NA NA NA n + 2 O -1 35 2.4 NA NA NA . . . O . . . NA NA NA n + m O -2 22 3.4 NA NA NA Q Covariates Treatment Outcome under A X * 1 X * 2 X * 3 A Y 1 R 1.1 20 NA 1 24.1 . . . R . . . . . . . . . n -1 R -6 NA 8.3 0 26.3 n R 0 15 6.2 1 23.5 n + 1 O -2 52 NA NA NA n + 2 O -1 NA 2.4 NA NA . . . O . . . NA NA n + m O NA NA 3.4 NA NA

General concept

Multiple imputation (MI) is one of the most powerful approaches to estimate parameters and their variance from incomplete data [START_REF] Little | Statistical Analysis with Missing Data[END_REF][START_REF] Kim | Statistical methods for handling incomplete Data[END_REF][START_REF] Schafer | Analysis of incomplete multivariate data[END_REF]. In a nutshell, for a single dataset, it consists in generating M plausible values for each missing entry, which leads to M completed datasets. Then, an analysis is performed on each imputed data set m = 1, . . . , M , to get an estimate for the parameter of interest, say θ as θ m and an estimate of its variance V m ( θ m ) and the results are combined using [START_REF] Rubin | Multiple Imputation for Nonresponse in Surveys[END_REF]'s rules to get correct inference with missing values, namely confidence intervals with the appropriate coverage.

Adapted multiple imputation for multiple data sources with different data design

For our problem, there are multiple possibilities to derive a multiple imputation strategy to generalize a treatment effect. This is due to the multi-source structure of the data and the fact that there is an additional complication due to the number of variables are not the same in the RCT and in the observational study. Indeed, we assume that the observational study does not include treatment and outcome but only covariates. We stress again that we assume missing values only occur in the covariates of both data. We suggest and describe three strategies to tackle this problem:

1. Within-study multiple imputation, see also Figure 3a: (a) Multiple imputation of the RCT: Impute M times the covariates of the RCT using

(X i , A i , Y i ) i: Qi=R .
(b) Multiple imputation of the observational data: Impute M times the covariates of the observational data using only the covariates (X i ) i: Qi=O .

(c) Create M × M complete tables by concatenating all possible combinations of imputed RCT and observational data. Estimate the treatment effect on every combination using any complete case estimator as in Section 2 and aggregate these estimations using Rubin's rules.

2. Ad-hoc joint covariates multiple imputation, ignoring the group variable, see also Figure 3b:

(a) Impute M times the joint datasets (the concatenation of the covariates from the RCT and the ones from the observational study) with covariates X and ignore the "source" indicator variable Q during the imputation.

(b) Concatenate the outcome Y and treatment A for each imputed RCT.

(c) Compute the M treatment effect estimators using any complete case estimator as in Section 2 and aggregate them using Rubin's rules.

3. Joint covariates multiple imputation, modeling the group variable as a fixed effect, i.e. keep the indicator variable Q indicating the corresponding "group"/"source" during the imputation, see also Figure 3c:

(a) Impute M times the joint datasets with covariates X and model the source indicator variable Q as a fixed effect during the imputation.

(b) Concatenate the outcome Y and treatment A for each imputed RCT.

(c) Compute the M treatment effect estimators using any complete case estimator as in Section 2 and aggregate them using Rubin's rules.

The first strategy has the advantage that it takes into account the outcome Y and treatment A which are dependent variables of the covariates X, when imputing the covariates of the RCT as suggested by [START_REF] Leyrat | Propensity score analysis with partially observed covariates: How should multiple imputation be used? Statistical methods in medical research[END_REF]; [START_REF] Seaman | Inverse probability weighting with missing predictors of treatment assignment or missingness[END_REF]; [START_REF] Mattei | Estimating and using propensity score in presence of missing background data: an application to assess the impact of childbearing on wellbeing[END_REF] as it models the entire joint distribution. The other strategies only consider the covariates X. The second strategy solely relies on the relationships between the covariates X and the assumption that these are stable across the data sources, i.e., Cov(X) = Cov(X|S = 1). The third strategy allows to additionally take into differences between both sources when imputing by modeling the source as a fixed effect. Strategy 3 can thus be seen as a fixed effect method, where the variable Q is included as a variable in the imputation model which allows, e.g., in case of multiple imputation with conditional regression models, to impute according to an analysis of covariance model. A drawback of this approach is that it generally inflates the between-group variability [START_REF] Andridge | Quantifying the impact of fixed effects modeling of clusters in multiple imputation for cluster randomized trials[END_REF].

In the case of several observational data sources, a fourth strategy could be considered, namely a multilevel multiple imputation approach, adding random effects into the model for each (observational) data source. Indeed, multilevel multiple imputation is specifically designed for cases of hierarchical or clustered observations, allows for a random intercept between groups (or, in our case, data sources) (van Buuren, 2018; [START_REF] Burgess | Combining multiple imputation and meta-analysis with individual participant data[END_REF][START_REF] Audigier | Multiple imputation for multilevel data with continuous and binary variables[END_REF]. It could be also appropriate when the observational data consists of patient records coming from different hospitals. It is indeed useful to encode the clustered structure of these records to account for between-hospital variability in terms of patient population and treatment practices. For a broader overview of multilevel imputation, we refer to [START_REF] Audigier | Multiple imputation for multilevel data with continuous and binary variables[END_REF]; van Buuren (2018).

All three strategies can be implemented easily using the R package mice (van Buuren, 2018) which uses conditional models such as linear and logistic regressions to perform multiple imputation.

The sketched fourth strategy could be implemented using the micemd R package [START_REF] Audigier | Multiple imputation for multilevel data with continuous and binary variables[END_REF]. Multiple imputation is suited if the missing values are ignorable as described in Section 2.4 and if the identifiability assumptions of the ATE τ in the full data case are met, namely Assumptions (2) and (3). Note that for the considered imputation strategies, increasing the number of imputations does not improve (significantly) the final result. There exists no clear rule about the number of multiple imputations to achieve good performances, however an accepted rule of thumb is to choose the number of imputations to be similar to the percentage of incomplete cases [START_REF] Hippel | How to impute interactions, squares, and other transformed variables[END_REF] or to the average percentage of missing data (van Buuren, 2018).

! ! ! ! ! ! ! ! ! ! " ! ! ! ! ! ! ! ! ! ! m n ! " !" !! !" !" ! !" !# !" "! !" "" ! !" "# ! ! ! ! !" #! !" #" ! !" ## Original observational data Original RCT Imputation 1 Imputation 2 Imputation M Imputation 1 Imputation 2 Imputation M ! ! ! " ! " !" # $ % & ! ' ' ! "# ( $ #%& $ "%& (
m n " 1 " 2 … " M Original observational data Original RCT Imputation 1 Imputation 2 … " = 1 % & "#$ % ' " " Q Imputation M Q Q Q (c)

Missing incorporated in attributes under adapted ignorability assumption

Multiple imputation makes classical identifiability assumptions of the causal effect and ignorability assumptions of the missing values mechanism [START_REF] Seaman | Inverse probability weighting with missing predictors of treatment assignment or missingness[END_REF][START_REF] Leyrat | Propensity score analysis with partially observed covariates: How should multiple imputation be used? Statistical methods in medical research[END_REF]. An alternative to handle missing covariates values consists in modifying the identifiability assumptions so that they directly handle missing values but do not necessarily require assumptions on the missing values mechanism. This can be seen as an advantage as it possibly allows for MNAR data, but the new identifiability assumptions may be more difficult to satisfy than in the full data case. More precisely, we extend the work of [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF] who adapt the unconfoundedness assumption to the incomplete covariates in order to identify the (average) treatment effect in the observational data case. In the following we lay out the modified assumptions for generalizability of the RCT to the target population we make to achieve a similar identifiability in the case of missing values in the RCT and the observational data. 3 We then explain how to generalize treatment effects under these assumptions. First, we introduce an additional notation, required for the following approach: The matrix of observed covariates can be written with

X * i X i R i + NA (1 -R i )
, with the element-wise multiplication and 1 the matrix filled with 1, so that X * i takes its value in the half discrete space

X * × 1≤j≤|X | {X j ∪ {NA}}.
And, similar to D'Agostino and [START_REF] D'agostino | Estimating and using propensity scores with partially missing data[END_REF]; [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF], we define the generalized conditional response surfaces µ * a and µ * a,1 as follows:

µ * a (x * ) = E(Y (a) | X * = x * ), µ * a,1 (x * ) = E(Y (a) | X * = x * , S = 1), (10)
The possibility to infer causal effects and to generalize the effect(s) from the RCT to another (broader) target population in the presence of missing data, i.e., using observations (X * i , A i , Y i , S i = 1) i=1, ..., n and (X * i ) i=n+1, ..., n+m , depends on the following additional assumptions on the joint law of (

X i , A i , Y i (0), Y i (1), S i = 1, R i ) i=1, ..., n and (X i , R i ) i=n+1, ..., n+m .
• Ignorability on trial participation, conditionally independent selection (CIS)

Assumption (2) and {Y (0), Y (1)} ⊥ ⊥ S | X * . ( 11
)
• Positivity of trial participation Assumption (3) and ∃ c * such that for all x * , P (π S (x * ) ≥ c > 0) = 1, and 0 < P (A = a | X * = x * , S = 1) < 1 for all a and for all x * such that

P (S = 1 | X * = x * ) > 0. ( 12
)
Assumption ( 11) means that being eligible to the RCT does not affect the potential outcomes conditionally on the covariates X * . The intuition behind these additions is to assume that instead of requiring conditional independence of potential outcomes and trial eligibility conditionally on all covariates, we only require conditional independence conditionally on the observed information, meaning the observed values and the pattern of missing values. Taking the example given in Figure 2, for observation 1, only X 1 and X 2 and the fact that X 3 is unobserved are decisive for trial eligibility, while for observation 2, only X 1 and X 3 and the fact that X 2 is missing decide upon trial eligibility, etc.

A possible scenario could be the existence of a list of sufficient but not necessary trial eligibility criteria. In other words, one could imagine a "check list" of L conditions and it is necessary to fulfill at least l < L of these to be eligible. Similar to the UDM assumption in [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF], Assumption (11) can be replaced by two sufficient assumptions: Assumption (2) and S ⊥ ⊥ X | X * , thus the term conditionally independent selection.

Generalized parameters and estimators

Since the conditional response surfaces µ * a and µ * a,1 now depend on X * rather than X and thus depends on the pattern of missing values, the estimators that involve an estimation of these quantities require an adaptation. Analogously to the generalized response surfaces, on can also define the generalized conditional odds α * , following the same logic.

µ * a (x * ) = E [Y (a) | X * = x * ] , µ * a,1 (x * ) = E [Y (a) | X * = x * , S = 1] . (13)
The resulting estimators are then formed analogously to the estimators in the full data case, by substituting the corresponding nuisance or intermediary parameters with their generalized counterparts. More explicitly, the outcome-model-based estimator defined by (5) becomes

τ CO, * = 1 m n+m i=n+1 µ * 1,1 (X * i ) -µ * 0,1 (X * i ) . ( 14
)
Fitting the new nuisance or intermediary parameters (13) is not straightforward, since these require to fit a separate regression model for each possible pattern r of missing values. For example, if we have three incomplete covariates X 1 , X 2 , X 3 , this means we need to fit a separate regression on {X 1 , X 2 , X 3 }, on {X 1 , X 2 }, on {X 2 , X 3 }, on {X 1 }, etc. We can see from this example that this is not possible in moderate and high dimensions with classical regression methods. This is why we propose to use random forests with a splitting criterion adapted to missing data. Indeed, as noted already by [START_REF] Athey | Generalized random forests[END_REF]; [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF], many modern machine learning methods, including tree ensembles and neural networks, can be adapted to this context and thus readily handle missing data and enable direct fitting of the generalized models above [START_REF] Josse | On the consistency of supervised learning with missing values[END_REF].

Nonparametric estimation. As an example of such a modern nonparametric estimation approach, we propose to estimate the generalized parameters via random forests [START_REF] Breiman | Random forests[END_REF][START_REF] Athey | Generalized random forests[END_REF], with missing data handled using the missing incorporated in attributes (MIA) method of [START_REF] Twala | Good methods for coping with missing data in decision trees[END_REF]. The resulting IPSW, CO and AIPSW estimators will be denoted by τ IP SW M IA , τ CO M IA , and τ AIP SW M IA respectively. In random trees, the MIA approach extends the classical splitting rules such that missing values are incorporated in the splitting criterion. More specifically, consider splitting on the j-th at-tribute and assume that for some individuals, the value of X j is missing, MIA treats the missing values as a separate category or code and considers the following splits:

• {i : X ij ≤ t or X ij is missing} vs. {i : X ij > t} • {i : X ij ≤ t} vs. {i : X ij > t or X ij is missing} • {X ij is missing} vs. {X ij is observed},
for some threshold t. The MIA approach does not seek to model why some features are unobserved; instead, it simply tries to use information about missingness to make the best possible splits for modeling the desired outcome. Thus the MIA strategy works with arbitrary missingness mechanisms and does not require the missing data to follow a specific mechanism. This MIA approach for (generalized) random forests is implemented in the R package grf [START_REF] Tibshirani | grf: Generalized Random Forests[END_REF] which is also used in the simulation part of this work presented in Section 5.

Parametric alternative. Parametric estimation is however possible in the case of logistic and linear regression models. This is based on work by [START_REF] Jiang | Logistic regression with missing covariates-parameter estimation, model selection and prediction within a joint-modeling framework[END_REF] and [START_REF] Schafer | Analysis of incomplete multivariate data[END_REF] for logistic and linear regressions with missing covariates. The functions µ * and α * that take in incomplete covariates x * are estimated via EM [START_REF] Dempster | Maximum likelihood from incomplete data via the em algorithm[END_REF]. The resulting IPSW, CO and AIPSW estimators will be denoted by τ IP SW EM , τ CO EM , and τ AIP SW EM respectively. The details of this approach are given in the Appendix A. However, a major limitation of this approach is that, in addition to the modified identifiability assumptions ( 11) and ( 12), in order to justify the use of the EM algorithm, one typically needs to make further assumptions on the missing value mechanism; in particular, this approach assumes the MAR mechanism (9). In other words, although we did not require the missing at random assumption to identify τ , this assumption is used for consistent parametric estimation of the generalized conditional models α * and µ * a,1 .

Simulations

We conduct a detailed simulation study to assess the performance of the previously introduced estimators to handle missing values. This controlled study allows to quantify the impact of different missing values mechanisms and identifiability assumptions on the final estimate for the ATE τ .

Data generation

Classical assumptions for identifiability in treatment effect generalization and on missing values mechanisms

We consider the selection model generated as a logistic model as follows

logit {π S (X)} = -2.5 -0.5X 1 -0.3X 2 -0.5X 3 -0.4X 4 , (15)
where every X is drawn from a multivariate normal distribution with mean 1 and covariance

matrix Σ such that Σ ij = 1 if i = j 0.6 if i = j
to have correlated covariates. The outcome is generated according to the linear model below such that X 1 is a treatment effect modifier and the true ATE τ is set to 27.4.

Y (a) = -100 + 27.4aX 1 + 13.7X 2 + 13.7X 3 + 13.7X 4 + with ∼ N (0, 1). ( 16)

We do not modify the treatment assignment mechanism since by assumption it is independent of the rest and in standard RCT design it is constant for all individuals. Missing values in the covariates are generated as follows, defining different models for the response indicator R (see Section 2.4):

• Missing values can occur in all four covariates.

• Proportion of missing values in each incomplete covariate: 20%.

• The missing values mechanism can be one out of the following and is implemented in the produce NA function proposed by [START_REF] Mayer | R-miss-tastic: a unified platform for missing values methods and workflows[END_REF]:

-MCAR where the probability to have missing values does not depend on any variable:

P (R i• = r|X i , S i , Q i , Y i , A i ) = P (R i• = r) (17)
We choose P (R ij = 0) = 0.2 for all i ∈ {1, . . . , n + m}, and j ∈ {1, 2, 3, 4}. -MAR:

P (R i• = r|X i , S i , Q i , Y i , A i ) = f (X obs(r) ), (18)
for some function f : X → [0, 1], for all i ∈ {1, . . . , n + m}. For example, missing values in X 1 are introduced for observation i using a logistic model on X 2 , X 3 , X 4 , assuming these three variables are observed for observation i. -MNAR:

P (R ij = 0|X i , S i , Q i , Y i , A i ) = g(X ij ), (19)
for some function g : X j → [0, 1], for all j and i ∈ {1, . . . , n + m}. In this study, we use a self-masking MNAR mechanism, i.e., the missingness of a variable depends on its value alone. More precisely, we use an upper quantile censorship approach. The quantile level q is chosen such that when missing values are generated on the qquantile at random, the requested proportion of missing values is achieved. For more details about this chosen approach, we refer to the documentation of the produce NA function4 .

We assume that the trial selection, randomization and potential outcomes are completely independent from the missing values. The simulation design is summarized by Algorithm 1.

Algorithm 1: Steps for simulation design under the standard assumption Result: Joint data table X * of RCT and observational data and additional variables A, Y for the RCT. 1 Sample N >> n observations X 1 , . . . , X N from the target population P(X); 2 Sample S according to the logistic model (15) on X; 3 Keep the {S = 1} indexed observations X R ← X {i: Si=1} as the RCT; 4 Sample A according to a Bernoulli distribution B(0.5) (coin flip); 5 Sample Y according to the linear model ( 16) on X R ; 6 Sample m observations X O from the target population P(X) as the observational data; 7 Concatenate the datasets: X ← [X T R , X T O ] T and append the indicator Q to the data (X ← [X , Q]); 8 Sample missing values for the n + m observations according to either (17), ( 18) or (19); For ease of comparison with the alternative simulation scenario, we provide the expression of the joint distribution over (X, R, A, Y, S, Q) factorized according to the underlying generative model:

p(X, R, A, Y, S, Q) ∝ p(R|S, Q, X, A, Y )p(S, Q, Y |X, A)p(A|X)p(X) ∝ p(R|S, Q, X, A, Y )p(Q|S)p(S|X, A)p(Y |X, A)p(A|X)p(X) (20)
In the MCAR case (8), this factorization simplifies as follows:

p(X, R, A, Y, S, Q) ∝ p(R)p(Q|S)p(S, Y |X, A)p(A|X)p(X).

Modified assumptions on treatment effect generalization

We also generate data according where the CIS assumption ( 11) is met. In such scenarios, we expect that the methods described in section 4 will work best.

The main difference with the previous setting of simulation, lies in the definition of the selection model, the outcome model remaining unchanged. To relate this setting to the previous one, we begin by giving the expression of the factorization of the joint distribution under these assumptions on the data generating process.

p(X, R, A, Y, S, Q) ∝ p(S, Q|X, R, A, Y )p(R|X, A, Y )p(Y |X, A)p(A|X)p(X) ∝ p(Q|S)p(S|X, R, A, Y )p(R|X)p(Y |X, A)p(A|X)p(X) (21)
Note that we can see from this factorization, the difference w.r.t. the previous case induced by the modified transportability assumptions, and in particular the modified ignorability (11) which induces a (conditional) dependence between R and {S, Q} by assumption.

In order to simulate data under the CIS assumption (11), we need to modify the definition of π S such that it becomes pattern-dependent.

logit {π S (X)} = -2.5 -0.5X 1 R 1 -0.3X 2 R 2 -0.5X 3 R 3 -0.4X 4 R 4 , (22)
The simulation design and the adapted CIS ignorability assumption is summarized in Algorithm 2 and is different from the one described in Algorithm 1.

Algorithm 2: Steps for simulation design under the CIS assumption.

Result: Joint data table X of RCT and observational data and additional variables A, Y for the RCT. 1 Sample N >> n observations X 1 , . . . , X N from the target population P(X); 2 Sample missing values for the N observations according to either (17), ( 18) or (19); 3 Sample S according to pattern-dependent logistic model on X; 4 Keep the {S = 1} indexed observations X R ← X {i: Si=1} as the RCT; 5 Sample A according to a Bernoulli distribution B(0.5) (coin flip); 6 Sample Y according to the linear model ( 16) on X R ; 7 Sample m observations X O from the target population as the observational data; 8 Sample missing values for the m observations X O using the same mechanism as before but possibly with different proportions;

9 Concatenate X R and X O (X ← [X T R , X T O ] T )
, and append the indicator Q to the data (X ← [X , Q]);

Estimation methods

We consider different scenarios of data generating processes by varying the type of missing values (MCAR, MAR, MNAR), ignorability assumption (standard or CIS), as well as the number of observations. We compare the following methods to handle missing values (the following acronyms are identical to the method labels used in Figures 45• Full data: we apply the standard full data estimators from Section 2 on the full data before introducing missing values (this would serve as a reference).

• Complete cases (CC): we apply the standard full data estimators from Section 2 on the complete observations extracted from the incomplete data (by deleting observations with missing values).5

• EM (see Section 4.1): we use EM to fit logistic and linear regression models of α * and µ * a,1

on the incomplete data using the R package misaem [START_REF] Jiang | Logistic regression with missing covariates-parameter estimation, model selection and prediction within a joint-modeling framework[END_REF].

• MIA (see Section 4.1): we use generalized random forests with MIA splitting criterion to estimate the generalized models α * and µ * a,1 on the incomplete data, using the R package grf [START_REF] Athey | Generalized random forests[END_REF].

• Multiple imputation (MI, Section 3): we apply the standard full data estimators from Section 2 on the imputed data (5-10 imputations obtained using the R package mice) where we use either within-study multiple imputation (WI-MI), or ad-hoc multiple imputation (AH-MI), or fixed effect multiple imputation (FE-MI).

We do not assess the random effect multiple imputation (i.e., the joint covariates multiple imputation with a 2-level model accounting for the multiple data sources, sketched Strategy 4 from Section 3) in this simulation study because this does not correspond to our motivating data example from the introduction. Note that for the EM and MIA approach, we only compute the IPSW, CO and AIPSW estimators (see Section 2.3) since the calibration weighting estimator in its current form is not applicable on incomplete data and future work is required to adapt this estimator to incomplete data.

Results

Due to the large number of different scenarios we consider in this simulation study, we first provide a summary of the expected behaviors of the different estimators in various cases before we report results of our experiments.

Summary of expected and empirical results

In Table 1 we summarize the expected results in terms of consistency of the different approaches, depending on the mechanisms that generate the data. For example, the MIA approach being suited for the modified ignorability assumption CIS (11), we expect it to perform well under CIS, independently of the missingness mechanism. In contrast, we do not expect it to be consistent under the standard ignorability assumption (2), whereas multiple imputation should be consistent under this standard assumption, provided the missingness is either MCAR or MAR. Some simulations results are in agreement with what is expected but there are also some gaps between the expected and empirical behavior in the chosen simulation settings. Indeed, in Figures 4 and5 we report the empirical bias with at 95% Monte Carlo confidence interval (based on a Monte Carlo standard error) of the estimated generalized ATE τ relative to the true value τ = 27.4, using n sim = 100 repetitions of each scenario. We define the empirical bias and its Monte Carlo standard error respectively by

B τ = Bias(τ ) = 1 n sim nsim i=1 τi -τ SE( B τ ) = 1 n sim (n sim -1) nsim i=1 (τ i -τ ) 2 ,
where τ = nsim i=1 τi , following [START_REF] Morris | Using simulation studies to evaluate statistical methods[END_REF].

• As expected, in Figure 4 we note that the full-data estimations are unbiased in all scenarios under the standard ignorability assumption (2). Under the CIS assumption, only the fulldata estimators that (partly) rely on the outcome model, namely CO, AIPSW, and CW are unbiased, whereas the parametric IPSW estimator fails under CIS. Surprisingly, the nonparametric full-data IPSW estimator recovers the true value in the MAR and MNAR case (see Figure 5).

• The complete case estimations are, unsurprisingly, unbiased only in the MCAR case.

• The behavior of the EM estimations is as expected: all estimators are unbiased under CIS under MCAR and MAR. However the AIPSW estimator performs better than the IPSW and the CO. In the MNAR case, the algorithm fails to converge.

• The MIA estimations overall have either small or no bias under CIS, especially the AIPSW estimator. Furthermore, under the CIS assumption, the AIPSW estimator always performs at least as well as the simple estimators (IPSW and CO). Under the standard ignorability assumption, the behavior is heterogeneous and tends toward biased results for all missingness mechanisms.

• Surprisingly, the within-study MI IPSW estimator is biased in all cases except the CIS+MCAR, and CIS+MNAR cases. This behavior is not expected and it remains to be investigated whether this is due to the simulation design.

• The joint fixed effect MI estimator (FE-MI) comes closer to the expected behavior of the multiple imputation approach than the ad-hoc MI (AH-MI) and the within-study MI (WI-MI) estimator as it has small or no bias under the standard ignorability and ignorable missingness (I+MCAR and I+MAR), but all three fail in the MNAR case (as expected).

Note that for the full data case, the choice of the estimator, namely parametric (in our case, generalized linear models) or non-parametric (here generalized random forests), has an impact on the bias, especially for the single-model estimators IPSW and CO. This is not surprising given the linear specification of the conditional odds and outcome models from ( 15) and ( 16) and the rather slow convergence of the chosen non-parametric method, random forest (grf), for linear models. 6In view of the results, if one has to recommend a method, it is preferable, from the present empirical evidence, to choose the MIA-AIPSW estimator or the joint multiple imputation coupled with the CW estimator. The MIA approach is simple to use in R thanks to the grf package which directly handles incomplete variables with the MIA splitting criterion for random forests. However, in terms of computational costs, this approach can be more expensive due to (automatic) parameter tuning. The joint multiple imputation approaches are easy to implement (with the mice package) but the relative computational running time is of similar magnitude as the MIA approach.

Impact of different proportions of missing values in the RCT and observational data

It is common that the RCT presents significantly less missing values than the observational study due to a more systematic monitoring of the data collection process. This invokes the question of how the above studied methods behave in the case of unbalanced proportions of missing values or different missing values mechanisms in the different data sources.

Extending the previous simulation study by such a case, we summarize in Figure 6 the performance of the different estimators under different scenarios of varying proportions of missing values in the RCT and the observational data when the data is MAR given S (or equivalently, we say it is MCAR in each data set). Note that this implies a slightly different factorization of the joint distribution over all variables than the one in the standard MCAR case (8):

p(X, R, A, Y, S) ∝ p(R|S, X, A, Y )p(S, Y |X, A)p(A|X)p(X) ∝ p(R|S, X, A, Y )p(S|X, A)p(Y |X, A)p(A|X)p(X) ∝ p(R|S)p(S|X, A)p(Y |X, A)p(A|X)p(X)
As expected, the complete case estimators are unbiased in this special case since conditionally on S, the data is MCAR. The doubly robust multiple imputation estimators can cope with very different proportions of missing values, either 10%, 50% in RCT and observational data respectively, or 5% and 22% respectively. These results are supporting our claim that the previous results and methods apply as well to the likely case of different proportions of missing values in the two studies. Indeed the following data analysis in Section 6 is an example of this case.

Application on critical care data

In this part, we come back to the medical question introduced in the beginning of this work about the potential effect of tranexamic acid (TXA) on mortality in patients with intracranial bleeding. We recall that, in order to answer this question, we have at disposal two data sources:

(1) CRASH-2, a multi-center international RCT, (2) Traumabase, an observational national registry.

A detailed data analysis of the observational registry to address the above medical question has been conducted by [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF]. We thus refer to this previous analysis for a detailed description of the observational registry as well as their findings. In a nutshell, leveraging only the observational registry does not provide evidence towards a beneficial (or detrimental) effect of TXA on trauma patients with TBI in terms of head-injury-related mortality.

We will first recall a summary of the findings of the original CRASH-2 study [START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF] before turning to focus on how the handling of missing values in the RCT and the observational registry impacts the final estimations of the population average treatment effect.

Findings of the CRASH-2 RCT

The CRASH-2 (Clinical Randomisation of Antifibrinolytic in Significant Haemorrhage) trial enrolled 20,211 patients in 274 hospitals in 40 countries between May 2005and 2009[START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF]. The aim of this trial was to study the effect of tranexamic acid in adult trauma patients with ongoing significant hemorrhage or at risk of significant hemorrhage, within 8 hours of injury (inclusion criteria), except those for whom antifibrinolytic agents were thought to be clearly indicated or clearly counter-indicated (exclusion criteria).7 . More precisely, eligible patients were defined as trauma patients within 8 hours of the injury, of age at least 16 years:

• with ongoing significant hemorrhage (systolic blood pressure less than 90 mmHg and/or heart rate more than 110 beats per minute)

• or who are considered to be at risk of significant hemorrhage.

The inclusion criteria and other baseline regressors are summarized in the graph of Figure 7. The results of the CRASH-2 study are reported in [START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF] and show a beneficial effect of TXA on the trial population for the primary outcome of interest (all-cause 28 day death).

Integration of the CRASH-2 trial and the Traumabase registry

In the following, we discuss common variables definition, outcome, treatment, and designs in order to leverage both sources of information. We recall the causal question of interest: "What is the effect of the TXA on brain-injury death on patients suffering from TBI?" This part is important for the harmonization of the study protocol.

Treatment exposure. The treatment protocol of CRASH-2 frames the timing and mean of administration precisely (a first dose given by intravenous injection shortly after randomization, i.e., within 8 hours of the accident, and a maintenance dose given afterwards [START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF]). The Traumabase study being a retrospective analysis, this level of granularity concerning TXA is unfortunately not available. Neither the exact timing, nor the type of administration are specified for patients who received the drug. However, the expert committee agreed that the assumption of treatment within 3 hours of the accident is very likely since this drug is administered in pre-hospital phase or within the first 30 minutes at the hospital [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF].

Outcome of interest. The CRASH-2 trial defined primary outcome as any-cause death in hospital within 28 days of injury. This outcome is also available in the Traumabase.

Covariates accounting for trial eligibility. For the CRASH-2 trial, four criteria determined inclusion: age (patients of at least 16 years old were eligible), ongoing or risk of significant hemorrhage (defined as systolic blood pressure below 90 mmHg or heart rate above 110 beats per minute, or clinicians evaluation of a risk), within 8 hours of injury and absence of a clear indication or counter-indication of antifibrinolytic agents. The necessary variables are also available in the Traumabase, either exactly or in form of proxies, which allows the estimation of the trial inclusion model on the combined data.

Additional covariates. Note that other covariates are (partially) available in both data sets, while not responsible of trial inclusion according to CRASH-2 investigators. But as this could still be covariates moderating the outcome and treatment effect, we include them in the outcome models used in the CO and AIPSW estimators ( 5) and ( 6). According to the two data sets, we could add three of them: sex (binary), type of injury (categorical, 1 =blunt, 2 = penetrating, 3 = blunt and penetrating), and initial Glasgow coma scale (numeric, integers from 3 to 15). Note that this three covariates are all mentioned in the baseline of CRASH-2 results [START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF], arguing that they should impact the outcome. The variables central capillary refill time and respiratory rate are also mentioned but are not available in the Traumabase, we thus omit them from this joint study.

Missing values. First, note that the RCT contain almost no missing values, whereas the variables for determining eligibility in the observational data contains important fractions of missing values, as shown in Table 8, while the sample sizes of the data sets are similar, see Table 3. While the MCAR assumption is admitted for the RCT [START_REF] Shakur-Still | Effects of tranexamic acid on death, vascular occlusive events, and blood transfusion in trauma patients with significant haemorrhage (CRASH-2): A randomised, placebo-controlled trial[END_REF], the missing values in the observational Traumabase are more complex and, according to the medical experts monitoring the collection process, partly non-ignorable. For example, the pre-hospital systolic blood pressure (SBP) is likely to be missing for patients with severe ongoing bleeding. Since the latter is informed in the hemorrhage risk variable, we could admit the missing values in the SBP variable as being MAR. A similar reasoning can be applied for the delay between the accident and treatment administration. However, there remains uncertainty as to whether the observed variables allow to fully explain the missingness in this variable.

Distribution shift. There are different ways of assessing the shift between the distributions of the two studies, e.g., by univariate comparisons. We provide a simplified comparison of the means of the covariates between the treatment groups of the two studies in Figure 9. This graph illustrates again the fundamental difference between the two studies, namely the treatment bias in the observational study and the balanced treatment groups in the RCT, but also a covariate shift between the two studies. To provide an example, the average patient age in the RCT is 7-9 years below the average age in the observational study; and there are only 16% of female patients in the RCT while there are over 20% in the observational study. Figure 9: Distributional shift and difference in terms of univariate means of the trial inclusion criteria (red: group mean greater than overall mean, blue: group mean less than overall mean, white: no significant difference with overall mean, numeric values: group mean (resp. proportion for binary variables). Graph obtained with the catdes function of the FactoMineR package [START_REF] Lê | FactoMineR: A package for multivariate analysis[END_REF].

Ignorability. Due to the design of the CRASH-2 study, namely the eligibility criteria which all need to be informed to decide upon trial eligibility, the modified ignorability assumption CIS (11) is less plausible to hold in this case and we rather consider the standard assumptions (2) and (3) to be satisfied by the CRASH-2 and Traumabase studies. 8 The additional assumptions concerning the missing values mechanism(s) have been outlined above and we consider that the assumptions for the multiple imputation strategy to be applicable are sufficiently plausible in this real-world example. The distribution of the estimated selection scores are given in Figures 10a (logistic regression via EM), 10b (generalized random forest with MIA), and 10c (logistic regression on joint fixed effect multiple imputations, MI). We notice that the scores obtained using EM and MI are similar and suggest that the positivity assumption is satisfied since we observe a good degree of overlap between the distributions of the scores for the two data sets. The scores estimated via MIA however concentrate around 0 and 1 for the observational and RCT observations respectively, suggesting poor overlap under this model. In Appendix B, we provide further comparisons of the estimated selection scores, pointing towards the multiple imputation strategy as the suited approach in this case. These results provide an additional argument in favor of the multiple imputation strategy which appears to be more adapted to the handling of missing values in this analysis; in particular we will apply a joint fixed effect multiple imputation strategy since it outperforms the other multiple imputation strategies in the simulation study of Section 5.

Final results when generalizing the ATE from the CRASH-2 trial to the observational data population

We now apply the estimators presented in this work and implemented first for the simulation study of Section 5. The confidence intervals for the corresponding point estimators are computed via non-parametric stratified bootstrap [START_REF] Efron | An introduction to the bootstrap[END_REF] using 100 bootstrap samples (using stratified sampling to preserve the study-specific sample sizes). We additionally report two consistent ATE estimators from the solely CRASH-2 data:

• Difference in mean: the difference in mean estimator (classical estimator for RCT);

• Difference in condmean ols the difference in conditional means where the we assume linear-logistic outcome models for Y (1) and Y (0) (estimator for RCT with smaller variance).

The former only involves treatment assignment A and outcome Y and thus requires no additional handling of the incomplete covariates; the latter is obtained using an EM algorithm for logistic regression with ignorable missing values in the covariates [START_REF] Jiang | Logistic regression with missing covariates-parameter estimation, model selection and prediction within a joint-modeling framework[END_REF]. And we present the AIPW estimators [START_REF] Robins | Estimation of regression coefficients when some regressors are not always observed[END_REF] for the observational study applied solely on the Traumabase data. For details about the derivation and properties of these estimators applied on incomplete observations we refer to [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF]:

• Generalized random forest after multiple imputation (MI AIPW),

• Generalized random forest using MIA splitting criterion (MIA AIPW).

Since AIPW combined with either missing incorporated in attributes (MIA) or multiple imputation (MI) is recommended by [START_REF] Mayer | Doubly robust treatment effect estimation with missing attributes[END_REF] when analyzing observational data, these are the estimators kept in this analysis.

When summarizing the results from the separate analyses on the RCT and the observational data respectively and the results from the joint analysis of both studies, we observe on Figure 11 a discrepancy between the different results. The only approach with consistent estimations throughout all estimators is the EM approach that points towards a beneficial effect of TXA on all-cause mortality. The joint fixed effect multiple imputation IPSW and AIPSW estimators (MI IPSW and MI AIPSW) also conclude on a beneficial effect of the treatment. However, the calibration weighting estimator (MI CW) does not find a significant effect. The large confidence intervals could be partly explained by the measurement noise in the administration delay variable in the Traumabase9 : contrary to the RCT, the Traumabase does not encode the exact delay of treatment administration, but is defined by a noisy proxy (delay between accident and admission to the resuscitation bay). However there exists evidence that administration delay is a treatment modifier for TXA and that only early administration has a beneficial effect [START_REF] Hijazi | Endogenous plasminogen activators mediate progressive intracerebral hemorrhage after traumatic brain injury in mice[END_REF]CRASH-2 Collaborators et al., 2011). This remark and the discrepancies between the findings, especially the different conclusions of the joint fixed effect multiple imputation estimators call for additional attempts to further refine the administration delay proxy variable in the Traumabase and potentially for additional analyses with supplementary data such as the CRASH-3 study [START_REF] Cap | Crash-3: a win for patients with traumatic brain injury[END_REF] which describes another slightly different TBI patient population.

Conclusion

In this work we have shown that missing values in multiple data sources require additional assumptions either on the mechanism that generated the missing values or on the data generating processes of both data sources to ensure generalizability (or transportability) of the treatment effect. We have proposed several estimators that are suited for generalizing a treatment effect from an RCT to a different target population described by observational covariate data. Which of the proposed methods is preferable depends on the underlying identifiability assumptions for generalizing the treatment effect. If the identifiability assumptions on the full data case are kept and only assumptions on the missingness mechanism are added, we recommend a joint multiple imputation that models the data source as fixed effect as long as the missing values are ignorable.

If the identifiability assumptions are altered to account for informative missing values implicated in the selection process, then estimators involving generalized conditional models are suited. These estimators rely on strong assumptions about the form of ignorability that is imposed, and we note that in many common examples in medicine or epidemiology it does not appear to be the most plausible one. The approach(es) that only imply additional assumptions on the missing values mechanism seem both closer to real application contexts and methodologically more feasible. Indeed, the presented data analysis on a treatment administered in critical care falls into this latter case and the results obtained with the proposed estimators from the two different presented approaches do not come to the same conclusion. This illustrates again the importance of choosing adequate assumptions for identifiability with missing covariate values and corresponding estimation strategies. On the methodological side, for all considered approaches and simulation settings, the question of varying missing values mechanisms across data sources remains to be addressed in more detail. Note as well that we have focused on missing values in both data sets, but the recommendations extend to the case where there are only missing values in the observational data and not the RCT due to different levels of systematic data collection.

Finally, the problem of incomplete observations addressed here is different from the problem of inconsistent variable sets between an RCT and an observational dataset, e.g., one variable is completely missing in one set, which is a challenging problem of a different kind that is left for future work. We note that this latter problem may cause issues of identifiability and is thus more related to the problem of unobserved confounding and a recent work proposes sensitivity analysis to address this issue [START_REF] Colnet | Generalizing a causal effect: sensitivity analysis and missing covariates[END_REF].

A Details on the estimation methods with missing values A.1 Prediction on new incomplete observations with parametric model

As mentioned in Section 4, it is possible to predict the outcome y for new incomplete observations, using the regression model estimated via EM, by marginalizing over the distribution of missing data given the observed. More formally, in the logistic regression case, using a Monte Carlo approach and maximum a posteriori estimator, it is possible to predict the response y for a new observation x i as follows:

1. Sample x For the linear case, the prediction proceeds in two steps:

1. Imputation of the new observation using the estimated variance-covariance matrix of the covariates Σ. xnew mis = -Σ -1 mis,mis

B Details on the critical care management application

We have seen in Section 6 that the scores obtained using EM and MI suggest that the positivity assumption is satisfied, while the scores estimated via MIA however concentrate around 0 and 1 for the observational and RCT observations respectively, suggesting poor overlap under this model. This can be also observed in the scatter plots in Figure 12, where we color the observations according to the values of the systolic blood pressure (SBP, a variable with an important fraction of missing values in the observational data, see Figure 8). We notice that MIA attributes a very low selection score to all observations with missing SBP value and thus selects the response indicator of the SBP variable to partly predict trial eligibility. This method has been studied more extensively in a regression framework and not in a classification framework and here we find that it predicts a class according to the presence of missing values and this is not necessarily what we intend when applying this method.

Figure 1 :

 1 Figure 1: Example of data structure in the full data problem setting. Left: complete underlying data with potential outcomes in the RCT. Right: observed data with factual outcomes.

Figure 2 :

 2 Figure 2: Example of data structure in the incomplete data problem setting. Left: complete underlying data with potential outcomes in the RCT. Right: observed incomplete data with factual outcomes.

  Joint covariates multiple imputation, modeling the source variable as fixed effect.

Figure 3 :

 3 Figure 3: Schematic illustrations of different multiple imputation strategies.

Figure 4 :Figure 5 :

 45 Figure 4: Empirical bias of generalizing ATE estimators under the standard ignorability assumption, 95% Monte Carlo confidence intervals, n = 1000.

Figure 6 :

 6 Figure 6: Bias estimates of generalized ATE under standard ignorability where missing values are "study-wise MCAR" (≡ MAR given S). For n ∈ {1000, 5000}, 20 repetitions. Case A={m=10 × n, RCT=10% NA, Obs=50% NA}; case B={m = 10 × n, RCT=5% NA, Obs=22% NA}.

Figure 7 :

 7 Figure 7: Causal graph of CRASH-2 trial representing treatment, outcome, inclusion criteria with S and other predictors of outcome (Figure generated using the Causal Fusion software by Bareinboim and Pearl (2016)).

Figure 8 :

 8 Figure 8: Percentages of missing values in each covariate for the Traumabase and CRASH-2 RCT.

  EM. Obtained with the misaem R package. MIA. Obtained with the grf R package. Joint fixed effect MI. Obtained with the mice R package.

Figure 10 :

 10 Figure 10: Estimated densities of the fitted selection scores.

Figure 11 :

 11 Figure 11: Separate and joint ATE estimators and 95% confidence intervals computed on the Traumabase (observational data set; purple), on the CRASH-2 trial (RCT; cyan), and generalized from CRASH-2 to the Traumabase target population (red). Number of variables used for adjustment in each context is given in the legend. The confidence intervals obtained on the observational data set and on the joint data sets are obtained via nonparametric bootstrap.

  ≤ s ≤ S ∼ p (x mis | x obs ) 2. Predict the response y by maximum a posteriori ŷ = arg max y p (y | x obs ) = arg max y p(y | x)p (x mis | x obs ) dx mis = arg max y E p z mis |x obs p(y | x)

  of response y using the imputed observation [x new obs , xnew mis ].

  x-axis: EM; y-axis: MIA. x-axis: EM; y-axis: joint fixed effect MI.

Figure 12 :

 12 Figure 12: Scatter plots of different estimated selection scores. The point color is set according to the systolic blood pressure (SBP) covariate values (missing SBP values are indicated by gray points).

Table 1 :

 1 Expected behavior under different assumptions about the data generating process and used estimation approach. Color code: blue=no bias; red=bias. This Table1is based on the Table2that summarizes the required assumptions for each method.

		Full data	Complete cases	EM	MIA	Multiple imputation
		CIS
	MCAR	(Assmpt. 11)
		I
		(Assmpt. 2)
		CIS
	MAR	(Assmpt. 11))
		I
		(Assmpt. 2)
		CIS
	MNAR	(Assmpt. 11)
		I
		(Assmpt. 2)

Table 2 :

 2 Methods for handling incomplete observations in treatment effect generalization and their assumptions on the underlying data generating process. (indicates cases that can be handled by a method, whereas marks cases where a method is not applicable in theory; () indicates cases without theoretical guarantees but with good empirical performance.)

	Covariates	Missingness	Identifiability of generalized τ	Models for (S, Y )
	multiva-riate normal	other MCAR MAR MNAR	I (≡ (2) & (3))	CIS (≡ (11) & (12))	Generalized linear models	Non-parametric models
	CC
	EM
	MIA
	MI						()

Table 3 :

 3 Sample sizes for the two studies.

		m	n	#treated #all-cause 28d death
	Traumabase 8248	-	686	1648
	CRASH-2	-	3727	1866	1176

This indicator S comprises several components: eligibility, selection into the trial and willingness to participate.

In the statistical and econometric literature, this term is sometimes also referred to as conditional odds ratio even though, by definition, it is an odds and not a ratio of odds.

Note that the assumptions for internal validity of the RCT are the same as in the full data case.

https://rmisstastic.netlify.app/how-to/generate/missSimul.pdf

Note that this approach is the most common default option in many implementations.

The random forest approach would require a lot of data to estimate linear regression functions; random forests are however known for their good performance in the presence of non-linearities and high order interaction terms[START_REF] Breiman | Random forests[END_REF].

Extract from the study protocol available at https://www.thelancet.com/protocol-reviews/05PRT-1.

The CIS assumption may be plausible in a context where trial eligibility (and participation) are defined via a set of sufficient conditions that are not all necessary conditions. More specifically, if we consider five covariates X 1 , . . . , X 5 and assume that there exist different possibilities to be included in the trial, e.g., a condition on X 1 , X 2 and X 3 regardless of X 4 and X 5 and another alternative condition on X 2 and X 5 . In the case of such a design, the CIS assumption is potentially a suited assumption to generalize a treatment effect onto another population.

Another direction to explore could be to trim the scores such as to avoid extreme weights. This could potentially reduce the size of the confidence intervals but it is known to generally provide a biased result since trimming of the weights induces an implicit change in the definition of the target population[START_REF] Li | Balancing covariates via propensity score weighting[END_REF]

Acknowledgements We would like to thank Shu Yang for fruitful discussions and her valuable feedback on our work. We thank Tobias Gauss, Jean-Denis Moyer and François-Xavier Ageron for their medical insights and interpretation of our data analysis; finally we thank the CRASH-2 trial investigators for sharing the trial data with us.

Funding IM was supported by a EHESS PhD fellowship and a Google PhD fellowship. The funding institutions had no role in the design and conduct of the study; collection, management, analysis, and interpretation of the data; preparation, review, or approval of the manuscript; and decision to submit the manuscript for publication.