Estimation and inference for causal spillover effects in egocentric-network randomized trials in the presence of network membership misclassification

Summary

1. INTRODUCTION

Causal inference is often conducted under the potential outcomes framework, which usually rules out interference, or spillover, among individuals, that is, it assumes that an individual’s potential outcomes depend only on their own exposure and not that of others (Rubin 1980; Hernán and Robins 2020; Imbens and Rubin 2015). However, in many settings, interference may be present. For example, in the HIV Prevention Trials Network (HPTN) 037 study (Latkin et al. 2009), the intervention is assigned to networks of people who inject drugs, and a participant’s engagement in HIV prevention behaviors can be affected by their peers, which can then affect their transmission risk. Understanding social interactions that drive disease transmission is key for strengthening interventions. In particular, peer education interventions, for which individuals are trained on prevention of risky behaviors and encouraged to disseminate knowledge and behavioral change in their peer networks, have been shown to effectively increase HIV knowledge and reduce risk behaviors among both those who receive training and network members (Cai et al. 2008; Aroke et al. 2023).

A common study design for assessing peer-based interventions, and frequently used in HIV risk behavior research, is the egocentric-network randomized trial (ENRT), in which index participants are randomly assigned to the intervention, and data on behavioral change is collected on both the index and members of their egocentric network (e.g. sexual or drug–injection partners) (Friedman et al. 2008; Latkin et al. 2009; Davey-Rothwell et al. 2011). ENRTs are able to assess both the direct effect of receiving the intervention on indexes, as well as the spillover effect on non-index members. We refer to the latter as the Average Spillover Effect (ASpE). In the presence of interference, most causal inference methods for the estimation of direct and spillover effects rely on the specification of an interference set, among whom interference is possible, for each individual. Assumptions on the extent of interference are required to define these interference sets. One common assumption is partial interference, typically assumed for cluster randomized trials (CRTs) either using the potential outcomes framework (Sobel 2006; Hudgens and Halloran 2008) or non-parametric structural equations with directed acyclic graphs (Balzer et al. 2019; Benitez et al. 2021), where interference is possible within randomization clusters and not between them. On the other hand, under network interference, interference sets are defined by network connections and typically restricted to the first degree, known as the neighborhood interference (Aronow et al. 2017; Forastiere et al. 2021; Lee et al. 2023). Additionally, it is typically assumed that an individual’s potential outcomes depend on not only their own exposure, but also a function of the exposures of those in their interference set, commonly referred to as an “exposure mapping function” (Aronow et al. 2017; Forastiere et al. 2021). Networks are usually assumed to be correctly measured; however, investigators may not have accurate information for identifying true network ties. In fact, it has been shown that improved methods are needed for ascertaining social connections in network-based studies (Rudolph et al. 2017; Young et al. 2018). When interference sets are misspecified, we will show that the ASpE estimated by usual approaches is biased. Bias correction methods are crucial to produce valid impact assessment.

Mismeasured exposures have been extensively studied in the statistical literature. For instance, the matrix and inverse matrix methods have long been available to non-parametrically correct for misclassified proportions and functions of proportions, when misclassification parameters can be estimated (Barron 1977; Marshall 1990). In causal inference, misclassified spillover exposures can arise when the exposure mapping or interference set is misspecified. Under misspecified exposure mapping, Aronow et al. (2017) showed that the Horvitz–Thompson estimator, a conventional inverse probability weighting estimator, is unbiased for the overall exposure effect, and Sävje (2024) showed that this estimator is consistent as long as the misspecification process is not systematic. Interference sets are misspecified when the extent of interference is wrongly assumed. Leung (2022) addresses wrongly assumed interference based on K-degree neighbors and proposes methods for “approximate neighborhood interference,” under which exposures of network members further from the index have smaller effects on the index.

On the other hand, we address the misspecification of interference sets due to network mismeasurement, which can result in mismeasured spillover exposures even if the extent of interference was correctly specified. A few papers have studied this issue. Bhattacharya et al. (2020) explored learning network structures using a score-based model selection algorithm. Egami (2021) showed that the average network-specific spillover effect is biased when interference may occur through multiple partially observed networks, and proposed sensitivity analysis methods to address this concern. Hardy et al. (2019) used mixture models to model the distribution of the latent true exposure conditions, and estimated causal effects using the Expectation–Maximization algorithm. Zhang (2020) used 2 network proxies observed for each individual to identify true networks, one of which is assumed to be an instrumental variable, and the other containing only one type of measurement error. While these publications aim to first estimate the true networks in order the estimate causal effects, which can be computationally intensive, we pursue a different approach where we augment the study design to include a validation study in which true networks are ascertained alongside the observed, error-prone ones, allowing us to estimate the network misclassification process and bias-correct the ASpE (Spiegelman et al. 2000; Carroll et al. 2006). To our knowledge, our novel approach has not been applied previously.

Our methods are applicable to ENRTs, where network misclassification is present. In Section 2, we present methods to bias-correct the point and interval estimates of the ASpE developed under a main study/validation study design, where the main study contains data on the outcome and observed networks of the study population, and the validation study contains measurements of the true and observed networks. Section 3 presents a simulation study used to investigate finite sample properties of our methods. Section 4 presents an application of our methods to the HPTN 037 study, where a contamination study allowed us to determine the participants’ true exposures.

2. MATERIALS AND METHODS

2.1. Notation, assumptions, and causal estimand

In an ENRT, a sample of |$ K $|  index participants are recruited and randomly assigned to the intervention, each nominating a set of network members. Throughout this article, we refer to those who are non-index participants as network members, and we define a (egocentric) network as the index and their set of network members. Let |$ ik $| be the |$ i $|-th participant in the |$ k $|-th egocentric network, with |$ i=1,\dots, n_{k} $| and |$ k=1,\dots, K $|⁠, for |$ N=\sum_{k=1}^{K}n_{k} $| study participants. For ease of notation, let |$ i=1 $| correspond to the index and |$ i\gt1 $| correspond to network members for all |$ k $|⁠. In addition, let |$ \mathcal{N}_{ik} $| be participant |$ ik $|’s network neighborhood, comprising participants who share a network link with |$ ik $|⁠. Note that in an ENRT, only indexes are asked to nominate their network members, and network members are typically not asked to delineate any of their network members (e.g. Latkin et al. 2009; Davey-Rothwell et al. 2011). Therefore, we only observe |$ \mathcal{N}_{1k} $| for all |$ k $|⁠, while |$ \mathcal{N}_{ik} $| with |$ i\gt1 $| is only partially observed, in that only the link with the index is observed, even though |$ \mathcal{N}_{ik} $| may include network members in the same or a different egocentric network, as well as out-of-study individuals (Fig. 1a).

Fig. 1.

Subfigure 1(a) depicts 2 egocentric networks, |$ k $| and |$ k^{\prime} $|⁠, with index participants |$ 1k $| and |$ 1k^{\prime} $|⁠, and their respective network members |$ ik $| and |$ ik^{\prime} $|⁠, |$ i\gt1 $|⁠. |$ 1k $| is randomized to the intervention, represented by a gray solid circle, while all participants not receiving the intervention are represented by white solid circles. Dashed circles represent out-of-sample individuals. Solid lines are the observed network connections between the index and their network members, while dashed lines are network connections that may exist but are unobserved. Subfigure 1(b) depicts the same egocentric networks in the presence of network misclassification, where solid lines are correctly measured network links, long dashed lines represent network links that were observed but not true, while dash-dot lines represent network links that were not observed but are true.

Open in new tabDownload slide

Let |$ \mathbf{R}_{K\times 1} $| be the intervention allocation vector, where |$ R_{k}=1 $| if network |$ k $| is randomized to intervention with probability |$ P_{R} $|⁠. Let |$ A_{ik} $| be the individual exposure status, where |$ A_{ik}=1 $| if |$ ik $| receives the intervention and |$ 0 $| otherwise. In an ENRT, |$ A_{1k}=1 $| if index |$ 1k $| receives the intervention, and |$ A_{ik}=0 $| for other study participants. Throughout this article, the exposure vector |$ \mathbf{A} $| is assumed to be perfectly measured under perfect compliance.

Let |$ Y_{ik} $| be |$ ik $|’s observed outcome, assumed to be binary throughout this paper. To define potential outcomes, we first make the neighborhood interference assumption (Aronow et al. 2017; Forastiere et al. 2021; Lee et al. 2023). Let |$ Y_{ik}(\mathbf{A}) $| be |$ ik $|’s potential outcome under |$ \mathbf{A} $|⁠. Let |$ \mathbf{A}_{\mathcal{N}_{ik}} $| be the exposure vector in |$ ik $|’s network neighborhood, |$ \mathcal{N}_{ik} $|⁠, and |$ g(\cdot):\{0,1\}^{|\mathcal{N}_{ik}|}\rightarrow\mathcal{G} $| the exposure mapping function applied to |$ \mathbf{A}_{\mathcal{N}_{ik}} $| with codomain |$ \mathcal{G}\subseteq\mathbb{R} $|⁠. We formally define the neighborhood interference assumption as follows.

Under Assumption 1, |$ ik $|’s potential outcomes depend on their own exposure and a function of the exposures of those in |$ \mathcal{N}_{ik} $|⁠, but not those outside of |$ \mathcal{N}_{ik} $|⁠. We refer to |$ G_{ik}=g(\mathbf{A}_{\mathcal{N}_{ik}}) $| as |$ ik $|’s spillover exposure, and denote |$ ik $|’s potential outcome by |$ Y_{ik}(a, g) $|⁠, which is the outcome that would be observed under |$ A_{ik}=a $| and |$ G_{ik}=g $|⁠. To define |$ G_{ik} $|⁠, we additionally assume the following:

Under Assumption 2, we assume that in addition to the individual exposure, |$ ik $|’s potential outcome depends on the overall number of neighbors who received the intervention, regardless of who among them is actually treated (Hudgens and Halloran 2008; Forastiere et al. 2021).

We further assume non-overlapping networks, that is, indexes are not connected among themselves, and network members can only be connected to one index. Let |$ \mathbf{M}_{N\times K} $| be the membership matrix, where |$ M_{ik, k}=1 $| if |$ ik $| is a network member of index |$ 1k $|⁠, and |$ M_{ik, k^{\prime}}=0 $| for all |$ k\neq k^{\prime} $|⁠.

Here, we allow for a network member to not belong to any network, that is, |$ \sum_{k^{\prime}}M_{ik, k^{\prime}}=0 $| for |$ i\gt1 $|⁠, which we will further discuss in later sections. Under Assumption 3, indexes are not connected among themselves, and since |$ A_{ik}=0 $| for all |$ i\gt1 $|⁠, |$ G_{1k}=0 $| for all |$ k $|⁠. For network members, although |$ \mathcal{N}_{ik} $| is only partially observed, since |$ A_{ik}=0 $| for all |$ i\gt1 $|⁠, by design no out-of-study individual can receive the intervention, and network members are only connected to one index under Assumption 3, |$ G_{ik} $| for |$ i\gt1 $| depends solely on the intervention received by their index and can be written as |$ G_{ik}=A_{1k} $|⁠. As such, |$ \mathbf{G}=\mathbf{M}\mathbf{R} $|⁠, and the set of possible potential outcomes for |$ ik $| are |$ Y_{ik}(a, g)\in\{Y_{ik}(1,0),Y_{ik}(0,1),Y_{ik}(0,0)\} $|⁠, as illustrated in Fig. 1a. Although Assumptions 1 and 2 are not empirically verifiable, in many ENRTs, Assumption 3 is ensured by the protocol and/or can be empirically verified for the observed networks. For example, as part of the inclusion-exclusion criteria in HPTN 037, a participant could not have been enrolled in more than one network. However, if the non-overlapping assumption is not guaranteed by the protocol and network members nominated by indexes may not be uniquely identified for anonymity purposes, a network member can, in theory, also be present in the data as an index or as a network member of a different index. In this case, the non-overlapping assumption must be made such that |$ G_{ik} $| may be correctly specified to be a binary variable, which is needed in the methods we propose in Section 2.5.

Our causal estimand, the ASpE, is the comparison of potential outcomes of network members if, possibly contrary to fact, all were exposed to the intervention through an index, versus not, while not receiving the intervention themselves. We consider the ASpE (⁠|$ \delta $|⁠) as a risk difference (RD) and a risk ratio (RR), defined as the following, where the expectation is taken over the sampling distribution of potential outcomes of network members (Imbens and Rubin 2015):

We make identification assumptions necessary for valid causal inference of the ASpE:

where Assumption 4 states that the spillover exposure for network members is probabilistic, and Assumption 5 states that the intervention assignment does not depend on potential outcomes, both of which hold by the ENRT design given that the intervention is randomized to index participants and therefore |$ 0 \lt Pr(G_{ik}=1) \lt 1 $| for |$ i\gt1 $| under the exposure mapping defined in Assumption 2. Given both assumptions, the ASpE can be identified as follows:

where |$ \mathbb{E}[Y_{ik}|A_{ik}=0, G_{ik}=1] $| does not condition on |$ i\gt1 $| as only those |$ i\gt1 $| can have the exposure status |$ A_{ik}=0 $| and |$ G_{ik}=1 $|⁠. Then, an unbiased estimator of the right side quantities would be unbiased for |$ \delta_{RD} $| and |$ \delta_{RR} $|⁠.

Unlike |$ A_{ik} $|⁠, |$ G_{ik} $| can be misclassified due to network misclassification. Let |$ \mathbf{M^{*}}_{N\times K} $| be the observed membership matrix, where |$ M^{*}_{1k, k^{\prime}}=0 $| for all |$ k, k^{\prime} $|⁠; |$ M^{*}_{ik, k^{\prime}}=0 $| for all |$ k\neq k^{\prime} $| and |$ i\gt1 $|⁠; and |$ M^{*}_{ik, k^{\prime}}=1 $| for |$ k=k^{\prime} $| and |$ i\gt1 $|⁠. Then, we define |$ P_{M}=Pr(M_{ik, k^{\prime}}^{*}=1|M_{ik, k^{\prime}}=1) $|⁠, |$ k=k^{\prime} $|⁠, as the probability of being classified into the correct network. Throughout this paper, |$ P_{M} $| is assumed to be independent of covariates; however, this assumption may be relaxed, and |$ P_{M} $| may be represented by a model that includes covariates. Let |$ \mathbf{G^{*}}=\mathbf{M^{*}}\mathbf{R} $| denote the observed spillover exposure vector. An illustration of possible misclassification scenarios is presented in Fig. 1b. In the absence of |$ G_{ik} $| where only |$ G^{*}_{ik} $| is observed, we would instead estimate the following quantities:

The relationship between |$ \delta $| and |$ \delta^{*} $| is derived in Section 2.3 under the assumption of non-differential misclassification, which assumes that the measurement error process is independent of potential outcomes conditional on the true spillover exposure, that is, formally,

For a network member, if the gain or loss of connections to an index is independent of their potential outcome and therefore their spillover effect, then Assumption 6 holds. For example, the misclassification process is non-differential if the risk behaviors that a network member would have if their index received the training did not depend on the likelihood of falling out touch with the index with whom they enrolled. On the other hand, if an index has the tendency to not correctly recall their network members and this is associated with ineffective conversations with them, or when one’s tendency to be wrongly considered as a friend by others is associated with their susceptibility to their friends’ influence, then the spillover effect would depend on the measurement error in |$ G_{ik} $| and Assumption 6 would be violated.

2.2. Main study: estimation of the ASpE under misclassified spillover exposure

Under Assumption 5 (unconfounded intervention), the ASpE can be estimated using sample average estimators. Using Table 1, if |$ G_{ik} $| were observed, the ASpE would be estimated as follows:

where |$ a=\mathbb{I}(\sum Y_{ik}=1, G_{ik}=1) $|⁠, |$ b=\mathbb{I}(\sum Y_{ik}=1, G_{ik}=0) $|⁠, |$ n_{1}=\mathbb{I}(G_{ik}=1) $|⁠, and |$ n_{0}=\mathbb{I} (G_{ik}=0) $|⁠. These estimators are unbiased for the quantities on the right-hand side (RHS) of the identifying formulas (2.3) and (2.4), and therefore, under Assumption 5, are unbiased for |$ \delta_{RD} $| and |$ \delta_{RR} $|⁠. When instead |$ G_{ik} $| is not observed, estimates must rely on statistics given |$ \mathbf{M}^{*} $|⁠, such as the upper script letters in Table 1. In this case, the ASpE estimators would be as follows:

where |$ \mathcal{A}=\mathbb{I}(\sum Y_{ik}=1, G^{*}_{ik}=1) $|⁠, |$ \mathcal{B}=\mathbb{I}(\sum Y_{ik}=1, G^{*}_{ik}=0) $|⁠, |$ N_{1}=\mathbb{I}(G^{*}_{ik}=1) $|⁠, and |$ N_{0}=\mathbb{I}(G^{*}_{ik}=0) $|⁠. These estimators are unbiased for the quantities on the RHS of the identifying formulas in (2.5) and (2.6), and therefore, under Assumption 5, are unbiased for |$ \delta_{RD}^{*} $| and |$ \delta_{RR}^{*} $|⁠.

2.3. Bias analysis of |$ \hat{\delta}_{RD}^{*} $| and |$ \hat{\delta}_{RR}^{*} $|

Let the outcome rate for |$ i\gt1 $| be |$ P_{Y0}\delta_{RR}^{G_{ik}} $|⁠, where |$ P_{Y0} $| is the baseline outcome rate. Under |$ G^{*}_{ik} $|⁠, |$ \hat{\delta}^{*}_{RD} $| and |$ \hat{\delta}^{*}_{RR} $| are biased. For simplicity, here we report formulas assuming that each network member is connected to one index in the true network under Assumption 3, that is, |$ \sum_{k^{\prime}}M_{ik, k^{\prime}}=1 $| for |$ i\gt1 $|⁠, and we relax this Assumption in Section S2 to allow for network members to not be connected to any index in the true network. Under Assumption 6 (non-differential misclassification), the bias with respect to the true effect, are:

where the bias is zero when |$ P_{M}=1 $|⁠, and towards the null when |$ P_{M} \lt 1 $|⁠, as shown in Section S2, along with formulas for the relative bias.

2.4. Validation study: estimation of network misclassification parameters

When |$ G^{*}_{ik}\neq G_{ik} $| as a result of network misclassification, the estimated ASpE is biased. However, when data on true networks or spillover exposures are available in a subpopulation, we can extend methods developed under a main study/validation study approach to bias-correct the ASpE (Barron 1977; Marshall 1990; Spiegelman et al. 2000). In the validation study, the relationship between |$ G_{ik} $| and |$ G^{*}_{ik} $| is modeled. Here, we estimate the sensitivity, |$ \theta=Pr(G_{ik}^{*}=1|G_{ik}=1) $|⁠, and specificity, |$ \phi=Pr(G_{ik}^{*}=0|G_{ik}=0) $|⁠, of spillover exposure classification among network members, to be used to bias-correct the ASpE estimated in the main study. We can define |$ \theta $| and |$ \phi $| as functions of |$ P_{M} $| and |$ P_{R} $|⁠. Intuitively, if |$ ik $| were exposed to a treated index, they could either truly belong to their assigned intervention network |$ k $|⁠, or a different network |$ k^{\prime} $| that is also an intervention network. An analogous argument can be made if |$ ik $| were unexposed to a treated index. Therefore, for |$ i\gt1 $|⁠, |$ \theta $|⁠, and |$ \phi $| are defined as follows:

where it is easily seen that |$ \theta=\phi $| under |$ P_{R}=0.5 $| as is often the case. Derivations are provided in Section S1, with extensions to allow network members to not be connected to any index in the true network. In this case, the derivation of |$ \theta $| is not affected; however, the derivation of |$ \phi $| is extended to additionally included the probability that |$ G_{ik}=0 $| because |$ ik $| does not belong to any network, and as a result, |$ \phi \gt\theta $| even under |$ P_{R}=0.5 $|⁠.

Additionally, if data are only available for true spillover exposures but not true networks, and |$ P_{M} $| cannot be estimated, assuming that network misclassification does not vary by covariates, |$ \theta $| and |$ \phi $| among |$ i\gt1 $| can also be estimated as follows:

where |$ \mathbb{I}(\cdot) $| is an indicator function. Besides |$ \theta $| and |$ \phi $|⁠, other misclassification parameters can also be estimated in the validation study, such as the positive and negative predictive values, |$ PPV=Pr(G_{ik}=1|G^{*}_{ik}=1) $| and |$ NPV=Pr(G_{ik}=0|G^{*}_{ik}=0) $|⁠, respectively.

In order to bias-correct the ASpE estimated in the main study using parameters estimated from the validation study, we assume that the misclassification parameters estimated in the validation study are generalizable to the main study had they been observed, that is:

2.5. Methods for bias correction

2.5.1. Matrix method estimator

Because the ASpE under |$ G^{*}_{ik} $| is estimated with misclassified counts, one approach to bias correction is to bias-correct these counts by estimating their true counterparts. Given the binary nature of |$ G_{ik} $|⁠, we extend the matrix method (Barron 1977), developed under Assumption 6 (non-differential misclassification), and correctly measured outcomes, such that |$ \mathcal{A}+\mathcal{B}=a\,+\,b $|⁠. Given |$ \hat{\theta} $| and |$ \hat{\phi} $| from the validation study, the bias-corrected estimators are given by:

Derivations are provided in Section S3.

Sampling variance estimation is provided by the multivariate delta method previously given by Greenland (1988); however, this would not account for the clustering of outcomes within networks, as may be expected in ENRTs. For example, in HPTN 037, the estimated intracluster correlation coefficient (ICC) was 0.15. Therefore, to take into account this correlation, we propose 2 approaches for variance adjustment. First, the variance can be inflated by the design effect, |$ 1+(\bar{m}-1)\hat{ICC} $|⁠, where |$ \bar{m} $| is the average network size (Donner and Klar 2000). Second, the variance can be estimated by network bootstrapping, where networks are resampled as a whole to maintain the within-network correlation structure (Davison and Hinkley 1997).

Lastly, for good finite-sample performance of these estimators, a few constraints on the values of |$ \theta $| and |$ \phi $| need to be considered, which are described in Section S4.

2.5.2. Inverse matrix method estimator

If the requirements needed for the matrix method estimator are not met, the ASpE estimate will be inaccurate, and the inverse matrix method (Greenland 1988; Marshall 1990), which does not impose these restrictions, may be considered. Under Assumption 6 (non-differential misclassification), neither method is uniformly more efficient than the other; however, the inverse matrix method is more efficient under differential misclassification, while the matrix method estimator approaches the efficiency of the inverse matrix method estimator when the validation study is small (Morrissey and Spiegelman 1999), as is the case in our illustrative example, HPTN 037. The inverse matrix method corrects for bias using |$ PPV $| and |$ NPV $| estimated separately for |$ Y_{ik}=1 $| and |$ Y_{ik}=0 $|⁠, with the following estimators:

where the subscript for |$ PPV $| and |$ NPV $| indicates the estimation among |$ Y_{ik}=1 $| or 0. Derivations are provided in Section S5. Variance estimation of the inverse matrix method estimators can also be achieved as described in Section 2.5.1.

2.5.3. Maximum likelihood estimator

While the matrix and inverse matrix estimators are easily implementable, there is no clear way to directly incorporate the effect of clustering. As such, we also propose a maximum likelihood method for clustered data by constructing an observed data likelihood using a mixed effects model, where the variance of the ASpE can be estimated by the inverse of the observed information matrix. Let |$ Pr(Y_{ik}=1)=P_{Y0}+\delta_{RD}G_{ik}+b_{k} $| be a model that estimates |$ \delta_{RD} $|⁠, where |$ b_{k}\sim N(0,\sigma^{2}_{B}) $| is the network random effect, where the limits of integration over |$ b_{k} $| shown below are imposed to ensure that probabilities are |$ \in(0,1) $|⁠. Note that the random effects are with respect to the sampling mechanism given by sampling the index participants together with their reported network members. For the estimation of |$ \delta_{RR} $|⁠, we replace the linear model with a log-binomial model. Let |$ V_{ik} $| be an indicator for being in the internal validation study, and |$ \boldsymbol{\omega}=(P_{Y0},\delta_{RD},\sigma^{2}_{B},PPV, NPV)^{T} $|⁠. Following Spiegelman et al. (2000) and Zhou et al. (2020), assuming no clustering of the misclassification process, the likelihood function is given by:

3. SIMULATION STUDY

3.1. Design

Based on a real-world ENRT, HPTN 037 (Latkin et al. 2009), we conducted a simulation study under scenarios of varying values of |$ P_{Y0} $|⁠, |$ \delta_{RR} $|⁠, |$ P_{M} $|⁠, and |$ P_{R} $| to assess trends in bias and confirm that empirical results align with those derived analytically. The simulation sample consisted of |$ K=1000 $| with |$ n_{k}=3 $| for all |$ k $|⁠, for a total of |$ 1000 $| index participants and |$ 2000 $| network members, such that each network consisted of 2 network members, which was the average number of network members per network in HPTN 037. A validation study was randomly sampled from the main study with |$ n_{v}=1000 $|⁠, among whom |$ \hat{\theta} $| and |$ \hat{\phi} $| were estimated. Under |$ P_{M}\in $| (0.1, 0.25, 0.5, 0.75, 0.9), |$ P_{Y0}\in $| (0.1, 0.25, 0.5, 0.75, 0.9), |$ \delta_{RR}\in $| (0.25, 0.75, 1.25, 3, 5), and |$ P_{R}\in $| (0.2, 0.5, 0.8), a total of 255 scenarios were studied, excluding 120 scenarios with combinations of |$ P_{Y0} $| and |$ \delta_{RR} $| that gave |$ P_{Y0}\delta_{RR} \gt 1 $|⁠. Outcomes were generated non-parametrically, first with probability |$ P_{Y0}\delta_{RR}^{G_{ik}} $|⁠, then to be correlated within networks to assess the robustness of the estimators to such correlations, with |$ Pr(Y_{ik}=1)\sim P_{Y0}\delta_{RR}^{G_{ik}}\times\exp(b_{k}) $| where |$ b_{k}\sim N(0,\sigma_{b}^{2}) $| were network random effects, |$ \sigma^{2}_{b}=\frac{\sigma^{2}_{e}ICC}{(1-ICC)} $| for |$ ICC\in $| (0.10, 0.25), and |$ \sigma_{e}^{2}= $| is the residual variance set to 1. Each network member was correctly classified with probability |$ P_{M} $|⁠, and if misclassified, each network except the true network had equal probability of being the misclassified network. Each scenario was simulated 2000 times, where the bias-corrected ASpE estimates were taken as empirical averages across all replications and compared to the true values. All simulations were repeated for |$ K=300 $| and |$ n_{k}=10 $| for all |$ k $| to assess any differences in results for larger network sizes.

3.2. Results

The bias and relative bias of |$ \hat{\delta}^{*}_{RD} $| and |$ \hat{\delta}^{*}_{RR} $| by |$ 1-P_{M} $|⁠, |$ P_{Y0} $|⁠, |$ \delta_{RR} $|⁠, and |$ P_{R} $|⁠, are presented in Tables S1–S4, respectively, where bias trends may be observed. The empirical bias and relative bias of the misclassified ASpE were virtually identical to the analytical ones, and obtained 95% coverage on average. As expected from the bias formulas given in Section 2.3, the ASpE estimated under misclassified networks was biased towards the null: when the true RR, |$ \delta_{RR} $|⁠, was <1, or equivalently, when |$ \delta_{RD} $| was negative, positive bias was observed as the misclassified risk was closer to the null; conversely, when |$ \delta_{RR} $| was >1, negative bias was observed. The bias of |$ \hat{\delta}^{*}_{RD} $| increased with |$ 1-P_{M} $|⁠, |$ P_{Y0} $|⁠, and |$ \delta_{RR} $| away from the null and did not depend on |$ P_{R} $|⁠. Similarly, the bias of |$ \hat{\delta}^{*}_{RR} $| increased with |$ 1-P_{M} $| and |$ \delta_{RR} $| away from the null, but decreased slightly with increased |$ P_{R} $|⁠, while did not depend on |$ P_{Y0} $|⁠.

We next compared the bias, standard error (SE), and coverage of the bias-corrected ASpE estimated using the matrix and inverse matrix estimators for select scenarios for uncorrelated outcomes, as presented in Table 2. The bias-corrected estimates were on average unbiased using either method, and the empirical and analytical SEs were virtually identical. Neither method performed uniformly better than the other across all scenarios, as discussed in Section 2.5.2; however, in certain scenarios where the matrix method failed due to unmet constraints mentioned in Section 2.5, the inverse matrix method performed much better.

When outcomes were correlated, the empirical SEs of the matrix method estimates were slightly larger than the analytical ones in some scenarios; however, these differences were not remarkable for the majority, suggesting robustness of our variance estimators to such correlations in ENRTs. For |$ \hat{\delta}^{*}_{RD} $|⁠, the empirical SEs were greater than their analytical counterparts in 3.92% and 8.24% of the scenarios when ICC = 0.10 and 0.25, respectively, compared to 1.17% when ICC = 0, while for |$ \hat{\delta}^{*}_{RR} $|⁠, this was true for 9.41% and 20.8% of the scenarios, compared to 4.71% when ICC = 0. Examples of such scenarios are presented in Table 3, with comparisons of empirical SEs to analytical ones that were unadjusted, inflated by the design effect, obtained from network bootstrapping, and using the likelihood-based method. The SEs obtained under all adjustment methods achieved good comparability to the empirical estimates.

Lastly, we re-ran our simulations with |$ K=300 $| and |$ n_{k}=10 $| (presented in Tables S5 and S6). There were no notable differences in bias trends, as expected since the bias formulas did not depend on sample size, whereas the inclusion of more network members in the study slightly decreased the SEs of the bias-corrected estimates. When |$ n_{k}=10 $| and outcomes were correlated, the clustering effect was slightly larger, also expected given the larger network size. Compared to results for |$ n_{k}=3 $|⁠, the empirical SEs for |$ \hat{\delta}^{*}_{RD} $| were greater than the analytical ones in 34.1% and 54.9% of the scenarios when ICC = 0.10 and 0.25, respectively. Network bootstrapping SEs were most comparable to the empirical values, while the design effect SEs tended to be over-inflated, as the differences between the empirical and unadjusted analytical SEs were not large enough compared to the design effect under |$ n_{k}=10 $|⁠.

4. ILLUSTRATIVE EXAMPLE: HPTN 037

HPTN 037 was a phase III ENRT conducted in the Philadelphia, PA, and Chiang Mai, Thailand, that assessed the impact of a peer education program on HIV risk reduction behaviors among people who inject drugs (Latkin et al. 2009). The study recruited 414 index participants who were asked to recruit at least one drug or sex network member to form their network neighborhoods. Both indexes and network members provided informed consent as part of the study eligibility criteria, and indexes randomized to the intervention received peer-educator training to encourage injection and sexual risk reduction with their network members. The primary study found that a decrease in reported injection risk behaviors between baseline and follow-up in both arms, with the intervention arm having a larger decrease. Here, we included only participants at the Philadelphia site and focused on a composite outcome defined as the occurrence of any reported HIV risk behavior, such as sharing injection equipment, at the 12-month visit, as did Buchanan et al. (2018) and Buchanan et al. (2024).

Here, |$ G^{*}_{ik} $| was given by the randomization networks, while |$ G_{ik} $| by an internal exposure contamination survey (Simmons et al. 2015; Aroke et al. 2023). At the 6-month follow-up, each participant was asked to recall terminology specific to the training program unlikely to be known outside of it. If a network member from recalled any of the 5 intervention-associated terms, we assumed that they were exposed to the intervention through a treated index such that |$ G_{ik}=1 $|⁠; likewise, if none of the 5 terms were recalled, we assume that there was no exposure from an index and |$ G_{ik}=0 $|⁠. Participants were also asked to recall a positive control term, to which everyone was exposed, as well as 3 negative control terms unrelated to the intervention that should not have been recalled. To mitigate the possibility that a network member was exposed to the intervention but could not recall the terms specifically, or as well as one who may not have been exposed but said they remembered the terms due to social desirability, only network members who recalled the positive control term and none of the negative control terms were included in the validation study for greater accuracy in the modeling of the misclassification process.

At the 12-month follow-up, the main study included 269 network members from 184 networks followed between December 2002 to July 2006, with 38 network members from 35 networks also included in the validation study. |$ \hat{\theta} $| and |$ \hat{\phi} $| were 0.60 and 0.79, respectively, such that 60% of network members were observed to be in an intervention network given that they were, while 79% of network members were observed to be from a control network given that they were either really connected to a control index, did not belong to any network, or were friends with a treated index but did not converse with them about the intervention and were therefore excluded from the set of network members of that index. Some degree of network misclassification was expected, as not all network members recruited by a treated index may have talked to them about the intervention, and some network members recruited by a control index may have become friends with a treated index and learned about the intervention.

In this study, the observed networks were non-overlapping by design, and we assume that the non-overlapping networks assumption was also satisfied for the true networks. There was no evidence of differential misclassification, as |$ \hat{\theta} $| and |$ \hat{\phi} $| estimated among |$ Y_{ik}=1 $| and 0 were not significantly different (P-value = 0.92). The conditions described in Section S4 for using the matrix method were also satisfied. The inverse matrix method may not be appropriate here because of the small validation size, such that |$ \hat{PPV}_{y} $| and |$ \hat{NPV}_{y} $| for |$ y\in\{1,0\} $| may be inaccurate and adversely affect the accuracy of the bias-corrected estimates as described in Section 2.5.2. Lastly, we assume that given the stringent criteria we set on the validation sample, the validation data provided us with an accurate measure of the spillover exposure, such that we may bias-corrected the spillover effects under the assumptions we invoked in this article.

Table 4 presents |$ \hat{\delta}_{RD} $| and |$ \hat{\delta}_{RR} $| estimated under the randomization networks, compared to those bias-corrected using the matrix method estimator and the maximum likelihood estimator. The ASpE was considerably biased towards the null due to substantial network misclassification: there was |$ 66\% $| bias in the RD, and |$ 124\% $| bias in the RR. We reported several confidence intervals. First, without considering clustering, we estimated the variance using the multivariate delta method, which when accounting for the estimation of |$ \theta $| and |$ \phi $| using a small validation study sample size, was substantially inflated. To take clustering into account, we inflated the variance by the design effect with an estimated ICC of 0.15 under a mixed effects model, performed network bootstrapping with 1000 bootstrap samples, and obtained the likelihood-based estimates. Given the small validation study sample size, confidence intervals obtained from bootstrap and likelihood-based variances were likely more justifiable alternatives in this study, and results showed a significant spillover effect of the intervention on the reduction of risk behaviors.

5. DISCUSSION

The impact of public health interventions, such as their diffusion and societal benefit, can be increased by leveraging behavioral spillover; however, the assessment of the spillover effect can be challenging because there is often substantial uncertainty about true social connections. We analyzed the bias of the ASpE when the recorded networks in an ENRT are misclassified, and developed bias-corrected estimators. Both our simulation study and illustrative example showed that the ASpE estimate is biased towards the null under non-differential misclassification of the networks, confirming our bias formulas. This finding is in line with past literature that showed the estimated effect is always biased towards the null when the outcome is perfectly classified and non-differential misclassification occurs on a binary exposure (Wacholder et al. 1995).

In this article, we assumed that networks were egocentric and non-overlapping, such that (i) network members can only belong to one network; and (ii) index participants cannot interfere with one another. We may relax these assumptions for the observed and true networks if information on the possible overlap between networks were available in the main and validation studies. To relax (i), network members may be influenced by more than one index. By re-defining the spillover exposure to be the connection to at least one treated index, |$ G_{ik}=\mathbb{I}(\sum_{j\in\mathcal{N}_{ik}}{A_{jk}} \gt 0) $|⁠, our methods still apply. On the other hand, if we keep the exposure mapping function as is, |$ G_{ik} $| is now a discrete variable, and our methods would need to be extended. To relax (ii), indexes can be exposed to the intervention from other indexes, and the potential outcome |$ Y_{ik}(1,1) $|⁠, which was previously undefined, is now identifiable given we can now have |$ G_{1k}=1 $|⁠. Other causal effects, such as the total effect of having both the individual and spillover exposures versus neither, |$ \mathbb{E}[Y_{ik}(1,1)]-\mathbb{E}[Y_{ik}(0,0)] $|⁠, can now be estimated.

Our methods have a few limitations. Our estimators rely on the availability of validation data such that the extent of the misclassification can be determined. In the absence of data on true networks, sensitivity analysis may be conducted to assess the robustness of the bias-corrected results (Egami 2021; Weinstein and Nevo 2023). Here, we focused on methods for which both the exposure and outcome are binary, and no covariate adjustment was considered. While it is reasonable to assume that both observed and unobserved covariates are balanced on average under randomization, covariate adjustment may still be needed. A likelihood-based approach to estimate the ASpE and misclassification parameters, as detailed in Section 2.5, can include covariates, although this may be more computationally intensive. It is also possible that the misclassification parameters may depend on covariates, and future research will develop estimators for this situation. Furthermore, our methods rely on the neighborhood interference assumption; however, it is also possible that network members may influence peers in other networks or even reinforce the knowledge of their index (Latkin et al. 2013). Our methods are currently being extended to the setting of CRTs where interference sets are typically defined by randomization clusters, but where interactions can exist across these clusters. We will propose estimators to bias-correct causal effects using approaches that will allow for non-binary spillover exposures and outcomes, as well as covariate adjustment to account for intervention non-compliance.

SUPPLEMENTARY MATERIAL

Supplementary material is available at Biostatistics Journal online.

FUNDING

This study was funded by the National Institutes of Health [1R01MH134715-01].

CONFLICT OF INTEREST

None declared.

DATA AVAILABILITY

Analysis code applied to a simulated dataset similar to the design of HPTN 037 can be found at https://github.com/arielchao/ENRT. The HPTN 037 study datasets are publicly available and can be requested from the Statistical Center for HIV/AIDS Research and Prevention through https://atlas.scharp.org/cpas.

REFERENCES