{'fraction_non_alphanumeric': 0.05235586598792471, 'fraction_numerical': 0.014028649224576773, 'mean_word_length': 4.387276785714286, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Efficiently transporting causal (in)direct effects to new populations under intermediate confounding and with multiple mediators

The same intervention can produce different effects in different sites. Transport mediation estimators can estimate the extent to which such differences can be explained by differences in compositional factors and the mechanisms by which mediating or intermediate variables are produced; however, they are limited to consider a single, binary mediator. We propose novel nonparametric estimators of transported stochastic (in)direct effects that consider multiple, high-dimensional mediators and intermediate variables. They are multiply robust, efficient, asymptotically normal, and can incorporate data-adaptive estimation of nuisance parameters. They can be applied to understand differences in treatment effects across sites and/or to predict treatment effects in a target site based on outcome data in source sites.


Introduction
The same intervention can produce different effects in different populations (e.g., Orr et al., 2003;Miller, 2015;Arnold et al., 2018). Different effects could arise from differences in: i) the distribution of compositional factors that modify aspects of the intervention's effectiveness (e.g., gender, age), ii) probability take-up or degree of adherence to the intervention, iii) the mechanism by which important mediating or intermediate variables are produced, and/or iv) the mechanism by which the outcome is produced in different populations, including differences population-or site-level contextual variables that are predictive of the outcome (Pearl and Bareinboim, 2018).
Transportability has been defined by Pearl and Bareinboim (2018) as the "license to transfer causal information learned in experimental studies to a different environment." Previously, we proposed using the transport graphs of Pearl and Bareinboim (2018) coupled with a transport estimator that predicts effects "transported" to a target population as a tool for quantitatively examining the extent to which differences in effect estimates between sites could be explained by factors i-iii above (Rudolph et al., 2017. In this previous work, we developed an efficient and robust semi-parametric estimator of transported stochastic (also called randomized interventional, see Didelez et al., 2006;VanderWeele et al., 2014) direct and indirect (what we refer to as (in)direct) effects in a target population . Although this previous estimator accounted for the presence of intermediate variables (those affected by treatment/exposure that could affect downstream mediator and outcome variables), it was limited in that it could only consider binary versions of a treatment/exposure variable, intermediate variable, and mediator variable, and assumed that the distribution of the mediator was known . To our knowledge, it is currently the only available estimator for transporting (in)direct effects.
However, many research questions involve continuous and/or multiple mediator variables. Thus, we address this methodologic gap by proposing novel nonparametric estimators of transported stochastic (in)direct effects that allow for multiple, possibly high-dimensional mediators without constraints on their distributions or the intermediate variables.
To motivate this work, we consider a research question from the Moving to Opportunity study (MTO), a multi-site randomized controlled trial conducted by the US Department of Housing and Urban Development, where families living in high-rise public housing were randomized to receive a Section 8 housing voucher that they could use to move to a rental on the private market (Sanbonmatsu et al., 2011). Families were followed up at two subsequent time points with the final time point occurring 10-15 years after randomization. In this study, some unintended harmful effects on children's mental health, substance use, and risk behavior outcomes were documented (Sanbonmatsu et al., 2011), and these overall effects were partially mediated by aspects of the peer and school environments (Rudolph et al., 2018b). However, these unintended harmful effects and their indirect effect components were not universal across sites (Rudolph et al., 2018a(Rudolph et al., , 2020. Our goal is to use the transportability framework and the novel estimators we propose to shed light on possible reasons why the intervention had harmful effects in some sites, particularly in Chicago, but not in others. For example, if we take Chicago as the site we would like to transport to, then we borrow information from the remaining sites to learn the outcome model, we can predict the effect for Chicago, standardizing based on the covariates, intermediate and mediating variables. The utility of borrowing or transporting information across sites applies more broadly than the above MTO example. It applies to questions that seek to: 1) understand differences in treatment, policy, or intervention effects across sites in multi-site trials or cohort studies, or to 2) predict treatment effects in a target site based on outcome data in source sites. This paper is organized as follows. In Section 2, we introduce notation, define the structural causal models we consider, and define and identify the transported stochastic (in)direct effects. In Section 3, we describe the efficient influence function (EIF), including a re-parameterization that allows for estimation with multiple and/or continuously distributed mediators, and derive the robustness properties of the EIF. In Section 4, we describe two efficient estimators for the transported stochastic (in)direct effects, based on the EIF derived in Section 3: an estimator that solves the EIF in one step and a targeted minimum loss-based estimator (TMLE). In Section 5, we present results from a simulation study in which we demonstrate the consistency, efficiency and robustness of the two estimators across various scenarios. In Section 6, we apply the two estimators to estimate the transported indirect effects of housing voucher receipt on subsequent behavioral problems as adolescents among girls in Chicago, operating through aspects of the school environment, borrowing information from the other MTO sites. Section 7 concludes the manuscript.

Notation and definition of (in)direct effects
Let O = (S, W, A, Z, M, SY ) represent the observed data, where S denotes a binary variable indicating membership in the source population (S = 1) or target population (S = 0), W denotes a vector of observed pre-treatment covariates, A denotes a categorical treatment variable, Z denotes an intermediate variable (a mediator-outcome confounder affected by treatment), M denotes a multivariate mediator, and Y denotes a continuous or binary outcome. Let O 1 , . . . , O n denote a sample of n i.i.d. observations of O. Note that the outcome is only observed for the source population/sites, S = 1, but we are interested in estimating effects for the target population/site, S = 0. We formalize the definition of our counterfactual variables using the following nonparametric structural equation model (NPSEM, Pearl, 2009) though equivalent methods may be developed by taking the counterfactual variables as primitives (Rubin, 1974). Assume the datagenerating process satisfies: is a vector of exogenous factors, and the functions f are assumed deterministic but unknown. We use P to denote the distribution of O. We let P be an element of the nonparametric statistical model defined as all continuous densities on O with respect to some dominating measure ν. Let p denote the corresponding probability density function. We denote random variables with capital letters and realizations of those variables with lowercase letters. We define Pf = f (o)dP(o) for a given function f (o). We use the following additional definitions. The function c(a, z, m, w) denotes P(S = 1 | A = a, Z = z, M = m, W = w), g(a | w) denotes P(A = a | W = w, S = 0), e(a | m, w) denotes P(A = a | M = M, W = w, S = 0), q(z | a, w) denotes the density of Z conditional on (A, W, S) = (a, w, 0), r(z | a, m, w) denotes the density of Z conditional on (A, M, W, S) = (a, m, w, 0), b(a, z, m, w) denotes E(Y | A = a, Z = z, M = m, W = w, S = 1), and t denotes P(S = 0). For a random variable X, we let X a denote the counterfactual outcome observed in a hypothetical world in which P(A = a) = 1. For example, we have Z a = f Z (S, W, a, U Z ), M a = f M (S, W, a, Z a , U M ), and Y a = f Y (W, a, Z a , M a , U Y ). Likewise, we let Y a,m = f Y (W, a, Z a , m, U Y ) denote the value of the outcome in a hypothetical world where P(A = a, M = m) = 1.

Transported stochastic (in)direct effects
We define the total effect of A on Y in the target population S = 0 in terms of a contrast between two user-given values a , a ∈ A among those for whom S = 0. The total effect can be decomposed into the natural direct and indirect effect. However, natural direct and indirect effects are not generally identified in the presence of a mediator-outcome confounder affected by treatment (Z, using our notation above) (Avin et al., 2005;Tchetgen and VanderWeele, 2014). Direct and indirect effects may be alternatively defined considering a stochastic intervention on the mediator (Petersen et al., 2006;van der Laan and Petersen, 2008;Zheng and van der Laan, 2012;VanderWeele et al., 2014;Rudolph et al., 2017). Let G a denote a random draw from the conditional distribution of M a conditional on (S, W ). The stochastic indirect effect (also called randomized interventional indirect effect) among those for whom S = 0 can be written: E(Y a ,G a − Y a ,G a | S = 0). This is the effect of A on Y that operates through M . The stochastic direct effect among those for whom S = 0 can be similarly written: E(Y a ,G a − Y a ,G a | S = 0), and is the effect of A on Y that does not operate through M . We focus on identification and estimation of θ = E(Y a ,G a | S = 0). Contrasts of θ under the values of a and a * given in the above definitions correspond to the transported stochastic (in)direct effects. Under the assumptions θ is identified and is equal to Assumption (i) states that, conditional on W , there is no unmeasured confounding of the relation between A and Y ; assumption (ii) states that conditional on W there is no unmeasured confounding of the relation between A and M ; (iii) states that conditional on (A, W, Z) there is no unmeasured confounding of the relation between M and Y ; (iv) states that there is a common outcome model across populations/ sites. It is this last assumption (iv) that allows us to transport or borrow information on the outcome model from other sites. If an alternative data source is available where Y is observed among those for whom S = 0, then the null hypothesis of equivalence between S = 0 and S = 1 can be tested nonparametrically (Luedtke et al., 2019).

Efficient influence function for θ
The efficient influence function (EIF) characterizes the asymptotic behavior of all regular and efficient estimators (Bickel et al., 1997;van der Vaart, 2002). In addition to being locally efficient, estimators constructed using the EIF have advantages of multiple robustness, which means that some components of the data distribution (i.e., nuisance parameters) can be inconsistently estimated while the estimator remains consistent. The multiple robustness property also allows the use data-adaptive machine learning algorithms in estimating nuisance parameters while retaining the ability to compute correct standard errors and confidence intervals. This is due to fact that the asymptotic analysis of the estimators yield second-order bias terms in differences of the nuisance parameters, and therefore allow slow convergence rates (e.g., n −1/4 ) for estimating these nuisance parameters.
Theorem 1 (Efficient influence function). For fixed a , a define (3) The efficient influence function for θ in the nonparametric model M is equal to This theorem makes two important contributions that advance the previous work deriving the EIF for a similar θ, but one that was limited in that it i) assumed that the distribution of M conditional on (A, W, S) was known, and ii) could only consider a single binary M . First, the EIF we derive does not assume that that the distribution of M conditional on (A, W, S) is known, reflected in the D P,M (o) component of the EIF in Equation 4, above. Second, we can overcome the challenge of estimating multivariate or continuous densities on the mediator, M , and intermediate variable, Z, as well as integrals with respect to these densities, if either M or Z is low-dimensional (though it can be multivariate) by using an alternative parameterization of the densities that allows regression methods to be used in estimating the relevant quantities. In the remainder of this work, we assume Z is low-dimensional (e.g., binary, as in our MTO illustrative application), though similar parameterizations may be achieved if M is low-dimensional.
The EIF given in Theorem 1 may be represented in terms of the expressions given in Lemma 1 below, which does not depend on conditional densities or integrals on the mediating variables.
Lemma 1 (Alternative representation of the EIF for univariate Z and multivariate M ). The functions h, u, and v may be parameterized: In the remainder of the paper, we denote η = (c, g, e, q, r, b, u, v) and D P,θ (o) = D η,θ (o). We let η denote an estimator of η, and η 1 denotes the probability limit ofη, which may be different from the true value. We derive the robustness properties of D η,θ (O) in the Supplementary Materials; they are given below in Lemma 2. The behavior of the term PD η 1 ,θ determines the robustness properties of the EIF as an estimating equation. Theorem 1 in the Supplementary Materials, together with the Cauchy-Schwarz inequality shows that PD η 1 ,θ yields a term of the order of: such that consistent estimation of θ is possible under consistent estimation of certain configurations of the parameters in η. The following lemma is a direct consequence.
We note that the cases (b 1 , v 1 , u 1 ) = (b, v, u) and (b 1 , g 1 , u 1 ) = (b, g, u) may be uninteresting if the re-parametrization in Lemma 1 is used to estimate the EIF, because in that case, consistent estimation of u and v will generally require consistent estimation of (b, c, q, r, e) in addition to the outer conditional expectations in Equations (6) and (7).

Estimators
We describe two efficient, robust estimators of θ. In subsection 4.1, we propose an estimator that solves the EIF estimating equation in one step (Pfanzagl and Wefelmeyer, 1985) (which we refer to as a one-step estimator), and in subsection 4.2, we propose a targeted minimum lossbased estimator (TMLE, van der Laan and Rubin, 2006), which is a substitution estimator that also solves the EIF estimating equation, but does it through iterative de-biasing targeted updates to nuisance parameters. We provide the R code to implement the proposed estimators, freely available at https://github.com/kararudolph/transport.
Letθ os andθ tmle denote the estimators defined below in Sections 4.1 and 4.2. Per the theorem below, the two estimators are asymptotically normal and efficient.
The proof of this theorem follows the general proof presented in Appendix 18 of van der Laan and Rose (2011). As a consequence, the variance of the estimators that follow can be estimated as the sample variance of the EIF, withθ and the nuisance parameters estimated as described above. This variance estimate may be used to construct Wald-type confidence intervals.
The Donsker condition of Theorem 2 may be avoided by using cross-fitting (Klaassen, 1987;Zheng and van der Laan, 2011;Chernozhukov et al., 2016) in the estimation procedure. Let V 1 , . . . , V J denote a random partition of the index set {1, . . . , n} into J prediction sets of approximately the same size. That is, In addition, for each j, the associated training sample is given by T j = {1, . . . , n} \ V j . letη j denote the estimator of η, obtained by training the corresponding prediction algorithm using only data in the sample T j . Further, we let j(i) denote the index of the validation set which contains observation i. The one-step and TMLE estimators may be adapted to cross-fitting by substituting all occurrences ofη The third condition of Theorem 2 can be satisfied by many data adaptive algorithms (e.g., lasso (Bickel et al., 2009), regression tress (Wager and Walther, 2015), neural networks (Chen and White, 1999), highly adaptive lasso (HAL)(van der Laan, 2017)); we use HAL in the simulations that follow.

One-step estimator
The one-step estimate of θ is given by the solution to the EIF estimating equation: We first describe how to estimate D η,Y . The regression b(a , z, m, w) can be estimated by fitting a regression of Y on W, A, Z, M among observations with S = 1 and then predicting values of Y setting A = a . The probability t is estimated as the empirical proportion of observations with S = 0 (i.e., in the target population). The regression function c(a , z, m, w) can be estimated by fitting a regression of S on W, A, Z, M and predicting the probability that S = 1 setting A = a . The treatment mechanism g(a | w) for a ∈ {a , a * } can be estimated by fitting a regression of A on (S, W ) and predicting the probability that A = a, setting S = 0. For the motivating example we consider here in which assignment of A is randomized, these can be estimated as the empirical probabilities that A = a and A = a * among those with S = 0. Under the reparameterization in Lemma 1 and in our motivating example, q(z | a, w) can be estimated by fitting a regression of Z on S, A, W and predicting the probability that Z = z setting A = a , S = 0. Likewise, r(z | a, m, w) can be estimated by fitting a regression of Z on S, A, M, W and predicting the probability that Z = z setting A = a , S = 0. The treatment probabilities e(a | m, w) and e(a * | m, w) can be estimated by fitting a regression of A on S, M, W and predicting the probability that A = a and A = a * , respectively, setting S = 0.
We next describe how to estimate D η,Z . For binary Z, the EIF simplifies to be The parameters t, g(a | w) and q(z | a , w) can be estimated as described above. For each z, u(z, a , w) can be estimated by regressing the quantity b(A, Z, M, W ) × h(A, Z, M, W ) on S, A, Z, W and getting predicted values, setting Z = z, A = a , S = 0.
To estimate D η,M , we estimate t, g(a * | w), b(a , z, m, w), q(z | a , w) as described above. The function v(a * , w) can be estimated by marginalizing out Z from b(a , z, M, W ) using q(z | a , W ) as predicted probabilities for each z, and then regressing the resulting quantity on A, W, S, and predicting values setting A = a * , S = 0.

TML estimator
We now describe how to compute a related TML estimator. We assume Y can be bounded in [0, 1], as described previously (Gruber and van der Laan, 2010). Many of the steps are identical to those for the one-step estimator, the differences are in the targeting of b(a , z, m, w), q(z | a , w), and v(a * , w).
Letb(a , z, m, w) be an initial estimate of b(a , z, m, w). We update this initial estimate using covariateĈ Next, letq(z | a , w) be an initial estimate of q(z | a , w). We update this initial estimate using covariateĈ in a logistic regression of Z with logitq(A | A, W )) as an offset, among the subset for which A = a , S = 0. Letˆ q be the MLE fitted coefficient associated with C q (A, W, S). The targeted estimate is given by logitq(z | a , w) = logitq(z | a , w) +ˆ qĈq (a , w).
To potentially improve performance in finite samples, we can move {ĝ(A | W ) ×t} −1 into the weights of a weighted logistic regression model, leavingû(A, 1, W ) −û(A, 0, W ) asĈ q (A, W ). Replacingb andq withb andq, the above steps can be iterated until the score equation is solved up to a factor of ( √ n log(n)) −1 . This iterating process and stopping criterion ensures that the efficient influence function is solved up to n −1/2 and mitigates risk of overfitting.
Next, we marginalize out Z fromb(a , z, M, W ) usingq(z | a , W ) as predicted probabilities for each z, and call the resulting quantity Q. This quantity is then regressed on (A, W ) among units with S = 0 to obtain an estimatorv(A, W ). This estimate is updated using covariatê in a logistic regression of Q with logitv(A, W )) as an offset, among the subset for which A = a * , S = 0. Letˆ v denote the MLE fitted coefficient on C v (A, W, S). The targeted estimate is given by logitṽ(a * , w) = logitv(a * , w) +ˆ vĈv (a * , w).
To potentially improve finite sample performance,Ĉ v (A, W ) may be moved into the weights of a weighted logistic regression model with intercept only. The empirical mean ofṽ(a * , W i ) among those for whom S = 0 is the TMLE estimate. Its variance can be estimated as the sample variance of the estimated EIF, given in Eq 4.

Simulation
We conducted a limited simulation study to examine and compare finite sample performance of these two estimators. We consider the data-generating mechanism (DGM) as follows. All variables are Bernoulli distributed with probabilities given by P (W 1 = 1) = 0.5 P (W 2 = 1 | W 1 ) = 0.4 + 0.2W 1 P (∆ = 1 | W ) = expit(−1 + log(4)W 1 + log (4) This DGM is formulated to align with features of the MTO study we use for the illustrative example. For example, A is randomly assigned and adheres to the exclusion restriction (Angrist et al., 1996), aligned with its role as an instrumental variable. In addition, we consider a modification of the observed data we have considered thus far: ∆ × O = ∆ × (S, W, A, Z, M, SY ), where ∆ is an indicator of selection into the survey sample. We assume the survey sampling weights are known or can be estimated asΓ where Π = P(∆ = 1 | X) and X represents unobserved variables used in the sampling design. Our previous identification result, which can alternatively be written as θ = E[v(a , W ) | S = 0], then becomes where we have added an index ∆ = 1 to emphasize that we are interested in parameters for the population from which the sample was drawn. The EIF is modified to be D P,∆=1 (o) = ΓD P (o), and the estimators of the previous section can be applied by using the weightsΓ i for each subject in the sample.
We consider estimator performance in terms of absolute bias, absolute bias scaled by √ n, influence curve-based standard error relative to the Monte Carlo-based standard error, standard deviation of the estimator relative to the efficiency bound scaled by √ n, mean squared error relative to the efficiency bound scaled by n, and 95% confidence interval (CI) coverage. We run 1,000 simulations for sample sizes N=1,000 and N=10,000. We also consider several model specifications. One in which all nuisance parameters in η are correctly specified, others that misspecify each nuisance parameter one at a time, another in which g(a | w), b(a , z, m, w), q(z | a , w) are correctly specified but the rest are not; and last, correctly specifying b(a , z, m, w), q(z | a , w), v(a * , w) but incorrectly specifying the rest. Under correct specification scenarios, we use HAL van der Laan, 2017) to fit each nuisance parameter. For incorrect specification, we use an intercept-only model. Table 1 shows simulation results for the transported stochastic direct effect, and Table 2 shows simulation results for the transported stochastic indirect effect comparing the one-step and TML estimators under correct specification of all nuisance parameters and various misspecifications. Given the robustness results in Lemma 2, we expect consistent estimates for all specifications in Tables 1 and 2 except when q is misspecified. We see this reflected in the results. We see that when the q model is misspecified, bias is more than an order of magnitude greater than any other specification for the transported stochastic direct effect in Table 1, and also greater, though to a lesser extent for the transported stochastic indirect effect in Table 2. 95% CI coverage using IC-based inference is close to 95% in the correctly specified scenario, but is poor when q is misspecified for the transported stochastic direct effect (Table 1), which is not unexpected given the biased estimates in this scenario. Coverage is less than 95% in other misspecified scenarios for both the transported direct and indirect effects (e.g., 68% when the b model is misspecified for the transported stochastic indirect effect, Table 2). This is not unexpected; the IC may not provide accurate inference when the IC at the estimated distribution using misspecified models does not converge to the IC at the true distribution. For robustness to extend to IC-based inference, further targeting of the nuisance parameters would be necessary that would preserve asymptotic linearity with a known influence curve at the cost of some efficiency. (van der Laan, 2014; Lastly, we note that under the smaller sample size of N=1,000 we see some deterioration in performance, particularly for the indirect effect, which is expected given that the true indirect effect is over five times smaller than the direct effect.

Illustrative Example
We apply the one-step and TML estimators proposed in Section 4 to estimate stochastic indirect effects transported across MTO sites, as described in the Introduction. Specifically, we are interested in the extent to which differences in: a) the distribution of individual-level compositional factors between the sites, b) take-up of the intervention (i.e., using the housing voucher to move), and c) distribution of school environment mediating variables can explain the difference in the indirect effect estimates between MTO sites. For this example, we consider the indirect effect of randomized receipt of a Section 8 housing voucher (A) and subsequent use (Z) on behavioral problems (Y ) (Zill, 1990) through aspects of the school environment (M , i) rank of the schools attended, and ii) whether ever attended a school in the top 50% of rankings, iii) number of schools attended, iv) number of moves since baseline, v) average proportion of students receiving free or reduced lunch, vi) ratio of students to teachers, vii) proportion of schools attended that were Title I, and viii) whether or not the most recent school attended was in the same district as the baseline school) among girls, comparing the Los Angeles (LA) and New York City (NYC) sites (S = 1, N=1,000) to the Chicago site (S = 0, N=600). We do this in order to illustrate our methods: the outcomes in Chicago were actually observed, so we can compared the transported estimate with estimates obtained using Chicago outcome data. Variables W and A were measured at baseline, when the children were 0-10 years old. Mediating variables were measured during the interval between baseline and the final follow-up timepoint 10-15 years later. The outcome was measured at the final follow-up timepoint. We account for a large number of covariates at the child and family levels: child age, race/ethnicity, history of behavioral problems, and gifted/talented status; parental education, marital status, whether or not the parent was under 18 at the birth of the child, employment, receipt of other public benefits, household size, feeling like the neighborhood was unsafe at night, feeling very dissatisfied with the neighborhood, whether or not the family had previously moved more than three times, wanting to move for better schools, whether or not the family had received a Section 8 voucher before, and poverty level of the baseline neighborhood. For this research question, randomization to receive a Section 8 housing voucher is an instrumental variable that affects M and Y through the intermediate variable of using the voucher to move out of public housing and into a rental on the private market (Z). We use the MTO sampling weights as described in Section 5. These weights account for sampling of children within families, changing randomization ratios, and loss to follow-up (Sanbonmatsu et al., 2011). We use data-adaptive methods for fitting the nuisance parameters, using a cross-validated ensemble of machine learning algorithms (Van der Laan et al., 2007), that includes generalized linear models, intercept-only models, and lasso (Tibshirani, 1996) that included all first and second-order predictors. To estimate the observed, non-transported stochastic indirect effects, we use non-transported versions of the one-step and TML estimators developed previously (Díaz et al., 2019). Standard errors are estimated using the sample variance of the influence curve. Figure 1 shows the transported and observed indirect effect estimates and their 95% CIs. Looking at the observed estimates, the indirect pathway from housing voucher receipt and use through the school environment to behavioral problems is protective for girls in LA and NYC, resulting in a reduction in behavioral problems at the final time point. However, the same pathway appears harmful for girls in Chicago, resulting in an increase of behavioral problems. Comparing the transported stochastic indirect effect estimate (one-step estimator: 0.0043, 95% CI: -0.0150, 0.0237, risk difference scale; TMLE: 0.0153, 95% CI: -0.0150, 0.0420) to the observed estimate for girls in Chicago (0.0089, 95% CI: 0.0007, 0.0171), we see that the two are similar even though the outcome data from Chicago was not used in the transported estimates. Thus, by taking the outcome model for LA and NYC and standardizing based on W, A, Z, M in Chicago, the predicted effect for Chicago is close to the observed. In contrast if they were not close to each other, this would suggest that the identification assumptions were not met. In the context of MTO, identification assumption (iv) of a common outcome model is arguably the most tenuous. This assumption would not hold in the presence of any contextual-level effects on the outcome model, such as the local economy, housing market conditions, segregation, etc.

Conclusions
We proposed estimators for transported stochastic direct and indirect effects under intermediate confounding and allowing for multiple, possibly related mediating variables arising from a true, unknown joint distribution. These estimators solve the efficient influence function; one that does so in one step and the other that is a substitution estimator that incorporates a series of targeting steps to optimize the bias-variance trade-off. We derived their multiple robustness properties and examined finite sample performance in a simulation study. Lastly, we applied our proposed estimators to better understand why a particular pathway from a housing intervention through changes in the school environment resulted in an unintended harmful effect on behavioral problems among girls in Chicago, when it led to improvements in behavioral problems among girls in other cities. Figure 1: Stochastic indirect effects estimates of being randomized to the Section 8 voucher group on behavioral problems score in adolescence, 10-15 years later, mediated through features of the school environment, among girls. Estimates and 95% CIs for observed (nontransported) and transported predicted effects. All results were approved for release by the U.S. Census Bureau, authorization numbers CBDRB-FY20-ERD002-023 and CBDRB-FY20-ERD002-024. Proof In this proof we will use Θ(P) to denote a parameter as a functional that maps the distribution P in the model to a real number. We will use P n to denote the empirical distribution of O 1 , . . . , O n . We will assume that the measure ν is discrete so that integrals can be written as sums. It can be checked algebraically that the resulting influence function will also correspond to the influence function of a general measure ν. The true parameter value is thus given by θ = Θ(P) = y,z,m,w y p(y | a , z, m, w, 1)q(z | a , w, 0)p(m | a , w, 0)p(w | 0).
The non-parametric MLE of θ is given by Θ(P n ) = y,z,m,w y P n f y,a ,z,m,w,1 P n f a ,z,m,w,1 P n f z,a ,w,0 P n f a ,w,0 P n f m,a ,w,0 P n f a ,w,0 where we remind the reader of the notation Pf = f dP. Here f y,a,z,m,w,1 = 1(Y = y, A = a, Z = z, M = m, W = w, S = 1), and 1(·) denotes the indicator function. The other functions f are defined analogously. We will use the fact that the efficient influence function in a non-parametric model corresponds with the influence curve of the NPMLE. This is true because the influence curve of any regular estimator is also a gradient, and a non-parametric model has only one gradient. The Delta method shows that ifΘ(P n ) is a substitution estimator such that θ =Θ(P), andΘ(P n ) can be written aŝ Θ (P n f : f ∈ F) for some class of functions F and some mapping Θ , the influence function of Θ(P n ) is equal to Applying this result to (8) with F = {f y,a ,z,m,w,1 , f a ,z,m,w,1 , f z,a ,w,0 , f a ,w,0 , f m,a ,w,0 , f a ,w,0 , f w,0 , f 0 : y, z, m, w} and rearranging terms gives the result of the theorem. The algebraic derivations involved here are lengthy and not particularly illuminating, and are therefore omitted from the proof.

B Lemma 1
Proof This result follows by replacing p(m | a , w) = p(m | a , z, w) g(a | w) g(a | w) q(z | a , w) r(z | a , m, w) e(a | m, w) e(a | m, w) in expression (4) in the main text.

C Additional results
Theorem 3 (Multiple robustness of the EIF). For notational simplicity, in this theorem we omit the dependence of all functions on w. Let P W denote the distribution of W conditional on S = 0. We have Proof For fixed a and a we have + g(a ) g 1 (a ) m 1 (a , z, m)q 1 (z | a )p(m | a )dP W dν(m, z) − v 1 (a ) g(a ) g 1 (a ) − 1 dP W dν(m, z) − v(a )dP W dν(m, z).
Using (9) Putting everything together yields