proof

# Heterogeneous Peer Effects in the Linear Threshold Model

Christopher Tran, Elena Zheleva

###### Abstract

The Linear Threshold Model is a widely used model that describes how
information diffuses through a social network. According to this model, an
individual adopts an idea or product after the proportion of their neighbors
who have adopted it reaches a certain threshold. Typical applications of the
Linear Threshold Model assume that thresholds are either the same for all
network nodes or randomly distributed, even though some people may be more
susceptible to peer pressure than others. To address individual-level
differences, we propose causal inference methods for estimating individual
thresholds that can more accurately predict whether and when individuals will
be affected by their peers. We introduce the concept of heterogeneous peer
effects and develop a Structural Causal Model which corresponds to the Linear
Threshold Model and supports heterogeneous peer effect identification and
estimation. We develop two algorithms for individual threshold estimation, one
based on causal trees and one based on causal meta-learners. Our experimental
results on synthetic and real-world datasets show that our proposed models can
better predict individual-level thresholds in the Linear Threshold Model and
thus more precisely predict which nodes will get activated over time.

## Introduction

Social networks play a vital role in spreading information, ideas, and
behaviors. For example, medical and agricultural innovations can spread
through the world (Rogers 2010), and new products can spread via word of mouth
or viral marketing (Kempe, Kleinberg, and Tardos 2003). Emotions such as
happiness (Fowler and Christakis 2008) and hatefulness (Ribeiro et al. 2018)
have also been observed to spread through social networks. These processes,
known as _information diffusion_ , have traditionally been studied in the
social sciences (Granovetter 1978), but more recently have motivated
applications in viral marketing (Kempe, Kleinberg, and Tardos 2003) and
recommender systems (Nikolakopoulos et al. 2019). Many models exist that
capture the diffusion process, including epidemic models (Kermack and
McKendrick 1927), voter models (Clifford and Sudbury 1973), the Independent
Cascade (Kempe, Kleinberg, and Tardos 2003), and the Linear Threshold Model
(LTM) (Granovetter 1978). In this work, we focus on LTM.

According to LTM, an individual is influenced to adopt a product or idea if
the proportion of their friends who have already adopted that product or idea
is above some threshold. LTM is a popular model used in many settings, such as
modeling the spread of ideas (Rogers 2010), predicting diffusion (Soni,
Ramirez, and Eisenstein 2019), and the development of many prominent influence
maximization algorithms (Kempe, Kleinberg, and Tardos 2003; Chen, Yuan, and
Zhang 2010; Goyal, Lu, and Lakshmanan 2011; Li et al. 2018). However, node
thresholds are typically assumed to either be the same for all nodes or
randomly distributed. Moreover, there are no methods for estimating individual
thresholds (Talukder et al. 2019), even though different individuals may have
different susceptibility to social influence. Moreover, LTM has not been
studied through a causal inference lens, even though LTM captures a causal
concept: “how many friends does it take to buy a product before they _cause_
an individual to buy the same product?” In our work, we seek to address these
shortcomings.

We develop two models for estimating individual-level thresholds from data.
Our models reflect real-world scenarios in which some individuals are more
easily influenced than others or differently by specific friends. To address
the variety of characteristics and behaviors of individuals, we propose a new
concept, _heterogeneous peer effects_ (HPEs) in networks, and contrast it with
heterogeneous treatment effects for IID data. While recent work has developed
methods for estimating heterogeneous treatment effects in networks, they only
consider the network structure as a confounder or as a proxy to latent
confounders for individual-level effects and do not estimate heterogeneous
peer effects (Veitch, Wang, and Blei 2019; Guo, Li, and Liu 2020).

To facilitate the estimation of peer effects, we develop a Structural Causal
Model (SCM) that is specific to LTM and encodes _interference_ by _contagion_.
Interference is the influence of treatments or outcomes of peers on an
individual, and contagion is the process of friends’ outcomes influencing an
individual’s outcome (Ogburn and VanderWeele 2014). For example, a user (e.g.,
Angelo in Figure 1) might buy sunglasses (outcome) if their friends already
bought them (contagion). Prior work has studied the role of SCMs in the
presence of interference (Ogburn and VanderWeele 2014; Bhattacharya, Malinsky,
and Shpitser 2020) and the estimation of average peer effects (Arbour, Garant,
and Jensen 2016) but not peer effect heterogeneity. To demonstrate the value
of our node threshold prediction algorithms, we evaluate them on three tasks,
node threshold estimation, activated node prediction, and diffusion size
prediction, using both synthetic and real-world data.

## Related work

Diffusion models have been studied for decades in epidemiology (Kermack and
McKendrick 1927) and opinion dynamics (Holley, Liggett et al. 1975). Two
popular models for social networks are the Linear Threshold Model (LTM)
(Granovetter 1978) and Independent Cascade (IC) Model (Goldenberg, Libai, and
Muller 2001). These models have been used in diffusion prediction (Soni,
Ramirez, and Eisenstein 2019) and influence maximization (Kempe, Kleinberg,
and Tardos 2003). While some work has focused on estimating parameters of IC
models (Saito, Nakano, and Kimura 2008; Bourigault, Lamprier, and Gallinari
2016; Kalimeris et al. 2018), not much work has explored threshold estimation
in the Linear Threshold Model.

A recent survey on influence maximization techniques using LTM highlighted the
lack of threshold estimation models (Talukder et al. 2019). Goyal et al. learn
edge diffusion probabilities but do not learn individual thresholds (Goyal,
Bonchi, and Lakshmanan 2010). Recent work identified a problem in threshold
estimation called the “opacity problem” (Berry et al. 2019). The opacity
problem states that the thresholds estimated from data are upwardly biased and
propose only to use nodes whose thresholds are _precisely measured_. They use
those thresholds in a regression model to estimate thresholds for all nodes.
In contrast, we estimate individual thresholds on “snapshots” of networks.

Our work focuses on defining and identifying heterogeneous effects of peers
through the use of a Structural Causal Model (SCM) (Pearl 2009) to estimate
individual thresholds in the LTM. Shalizi and Thomas explore causal inference
under the presence of homophily and contagion (Shalizi and Thomas 2011). They
address the problem of identification of social effects when influence is
decoupled from group effects (Manski 1993). Ogburn and VanderWeele presented
an extensive discussion of SCMs for interference (Ogburn and VanderWeele
2014). Arbour et al. utilize abstract ground graphs for estimating average
effects in social networks (Arbour, Garant, and Jensen 2016). To the best of
our knowledge, there have been no SCMs developed for the LTM. In addition,
recent work has introduced methods for estimating heterogeneous _treatment_
effects in networks but do not estimate effects of peers (Guo, Li, and Liu
2020).

## Problem Description

To define the problems of individual threshold estimation and conditional
average peer effect estimation, we first present the data, linear threshold
model, and define heterogeneous peer effects.

### Data Model

Figure 1: A toy example illustrating the Linear Threshold Model diffusion
process, where Molly and James are activated. In the first time step, Angelo
will become activated, since his threshold, $\theta_{Angelo}=1/3$.

Let $G=(\mathbf{V},\mathbf{E})$ denote an attributed social network, where
$\mathbf{V}$ is the set of nodes and $\mathbf{E}$ is the set of edges between
nodes. If two nodes $v,u\in\mathbf{V}$ are connected by an edge, then
$(v,u)\in\mathbf{E}$, denotes an edge from $v$ to $u$. Define the set of
neighbors of $v$ to be $N(v)=\\{u:u\in\mathbf{V},(u,v)\in\mathbf{E}\\}$. Each
node $v\in\mathbf{V}$ has an observed $m$-dimensional vector of attributes,
$\mathbf{X}_{v}$, unobserved characteristics, $\mathbf{U}_{v}$, and an outcome
of interest $Y_{v}\in\\{0,1\\}$, which is a binary indicator of whether the
node is activated (e.g., whether an individual has bought new sunglasses). We
define the set of activated nodes at time $t$ to be
$\mathcal{D}_{t}=\\{v:Y_{v}=1\\}$.

### Linear Threshold Model

According to the Linear Threshold Model (LTM), each node $v$ has an activation
threshold $\theta_{v}$. Given an initial set of activated nodes,
$\mathcal{D}_{0}\subseteq\mathbf{V}$, diffusion occurs in discrete steps,
$t=1,2,\dots,T$. In each time step $t$, a node
$v\in\mathbf{V}\setminus\mathcal{D}_{i}$ is activated if the _activation
influence_ , the weighted proportion of its activated neighbors, reaches the
node’s threshold $\theta_{v}$:

$\sum_{u\in N(v)}w_{uv}Y_{u}^{t}\geq\theta_{v},$ (1)

where $w_{uv}$ is the normalized influence weight of neighbor $u$ on $v$.
According to LTM, nodes can only become activated as activated neighbors
increase. In practice, node thresholds are implemented by considering random
or uniform thresholds (Talukder et al. 2019) even though the propensity to be
influenced can vary from individual to individual. Our goal in this work is to
estimate thresholds for all nodes.

######  Problem 1.

(Individual node threshold estimation for LTM) Let $G=(\mathbf{V},\mathbf{E})$
be a graph and $\bm{\theta}=\\{\theta_{v}\mid v\in\bm{V}\\}$ be a set of node
thresholds. The goal is to estimate the thresholds for all nodes,
$\bm{\hat{\theta}}=\\{\hat{\theta}_{v}\mid v\in\mathbf{V}\\}$, such that the
average mean squared error between $\theta$ and $\hat{\theta}$ is minimized:

$\operatorname*{argmin\,\,}_{\hat{\theta}_{v}}\text{MSE}=\frac{1}{|V|}\sum_{v}^{|V|}\Big{(}\theta_{v}-\hat{\theta}_{v}\Big{)}^{2}.$
(2)

Figure 1 shows a toy example of a social network with five individuals. Each
node has their own threshold (e.g., $\theta_{\text{Angelo}}=1/3$ means
Angelo’s threshold is $1/3$). The initial set of activations are the
individuals who have adopted a product (sunglasses), which consists of two
individuals: $\mathcal{D}_{0}$ = {Molly, James}. Assuming equal weights,
Angelo will buy new sunglasses in the first time step since one of his three
friends (James) has already bought them. No one else will buy sunglasses in
subsequent steps since Josh’s threshold is $1$ and Michelle’s threshold is
$2/3$.

### Causal Inference in LTM

Here, we connect causal inference with the Linear Threshold Model under the
potential outcomes framework and define heterogeneous peer effects and their
estimation.

#### Causal inference in networks

A common assumption in estimating causal effects is the stable unit treatment
value assumption (SUTVA) (Rubin 1978), which is the assumption that the
outcome of an individual is independent of the treatment assignment of other
individuals: $Y_{v}(Z_{v})\perp Z_{u},\forall u\neq v\in\bm{V}$. However, this
assumption is violated in network data because of the interaction between
individuals - known as _interference_ (Ogburn and VanderWeele 2014) \- which
can lead to biased causal effect estimations. Interference occurs in the LTM,
where activation depends on the activation of neighbors. Thus, we need to
account for node features _and_ network interference when estimating effects.

Define $\mathbf{A}=\\{Z_{v}\mid v\in\mathbf{V}\\}$ to be the set of treatment
variables for all nodes in the network. Then, the outcome of a node is
dependent on two variables: the individual treatment and the set of treatments
in the network: $Y_{v}(Z_{v},\mathbf{A})$. We can define the average treatment
effect (ATE) and the average peer effect (APE) (Arbour, Garant, and Jensen
2016; Fatemi and Zheleva 2020) as:

$\displaystyle\text{ATE}=E[Y(Z=1)-Y(Z=0)\mid\mathbf{A}=\mathbf{a}],$ (3)
$\displaystyle\text{APE}=E[Y(\mathbf{A}=\mathbf{a})-Y(\mathbf{A}=\mathbf{a}^{\prime})\mid
Z=z].$ (4)

ATE estimates the effect of treatment on a node, keeping other treatment
assignments the same. APE estimates the effect when keeping a node’s treatment
fixed while changing treatments of other nodes.

_Contagion_ is a case of _interference_ , where peer outcomes in the previous
time step affect individual outcomes. In the LTM, the outcome of peers in the
previous time step are captured by $\mathbf{A}=\\{Y_{v}^{t-1}\mid
v\in\mathbf{V}\\}$, which is the treatment variable of interest. Neighborhood-
level variables (e.g. $\mathbf{A}$) are _relational_ and individual-level
variables (e.g. $\mathbf{X},Y,Z$) are _propositional_ variables.

#### Heterogeneous peer effects

Heterogeneous treatment effects (HTEs) occur when subpopulations react
differently to treatment. HTE is expressed through the conditional average
treatment effect (CATE). CATE is defined as the expected difference in
potential outcomes, conditional on an individual’s features:
$\tau(\mathbf{x})\equiv E[Y(z)-Y(z^{\prime})\mid\mathbf{X}=\mathbf{x}]$, which
results in personalized effect estimations. In LTM, this is the expected
change in activation if the node was activated vs. not activated in the
previous time point with features $\mathbf{X}_{v}=\mathbf{x}$. So far, we have
ignored the influence of neighbors. Next, we describe how to incorporate it
into the HTE estimation. First, we define CATE under interference as:

$\displaystyle\tau(\mathbf{x})=E[Y(Z=1)-Y(Z=0)|\mathbf{X}=\mathbf{x},\mathbf{A}=\mathbf{a}].$
(5)

Different individuals may have different peer effects based on personal
characteristics, such as personality and susceptibility to influence. To
capture this idea, we define _conditional average peer effect (CAPE)_ :

$\rho(\mathbf{x})=E[Y(\mathbf{A}=\mathbf{a})-Y(\mathbf{A}=\mathbf{a}^{\prime})\mid\mathbf{X}=\mathbf{x},Z=z]$
(6)

_Heterogeneous peer effects (HPEs)_ occur when there are at least two possible
assignments of $\mathbf{X}$, $\mathbf{x}_{i}$ and $\mathbf{x}_{j}$, for which
$\rho(\mathbf{x}_{i})\neq\rho(\mathbf{x}_{j})$. Since not all nodes in the
network affect every other node, we use neighbors of nodes for contagion. We
can define $\mathbf{A}_{v}=\\{Y_{u}^{t-1}\mid u\in N(v)\\}$, the set of
outcomes of neighbors of $v$ in the previous time step, which affects the
outcome and APE in the current time step.

######  Problem 2.

(Conditional average peer effect estimation) Let $G=(\mathbf{V},\mathbf{E})$
be a graph. Let $\hat{\rho}(\mathbf{X})$ be an estimate for the CAPE. The goal
is to minimize the expected error between the estimate and its true value:

$\min_{\hat{\rho}}E[(\rho(\mathbf{X})-\hat{\rho}(\mathbf{X}))^{2}].$ (7)

In general, it is not trivial to estimate CAPE since relational variables are
sets, and $\mathbf{A}_{v}$ for $v\in\mathbf{V}$ are different sizes. To
circumvent this problem, we model relational variables through aggregate
functions, such as the mean or variance (Arbour, Garant, and Jensen 2016).
Next, we discuss the link between CAPE estimation and threshold estimation.

## HPE and Threshold Estimation

We begin by linking CAPE estimation to threshold estimation. Second, we
introduce a causal model to identify the CAPE. Then, we map the CAPE
estimation to a problem of estimating triggers for heterogeneous effects,
which allows the estimation of node thresholds for LTM.

### Causal model for LTM

In LTM, we define the relational variable at time $t$ as the _activation
influence_ on the left side of Equation 1:

$I_{v}^{t}=\sum_{u\in N(v)}w_{uv}Y_{u}^{t-1}.$ (8)

We define the treatment variable, $Z_{v}$, as the previous outcome:
$Z_{v}=Y_{v}^{t-1}$ and the CAPE is defined with the aggregation function in
(8):

$\displaystyle\rho(\mathbf{x})=$ $\displaystyle
E[Y(I^{t}=i^{t})-Y(I^{t}=i^{t^{\prime}})\mid\mathbf{X}=\mathbf{x},Z=z],$ (9)

where $i^{t}$ and $i^{t^{\prime}}$ correspond to two different values of
$\mathbf{A}$. In order to identify this effect, we need to account for correct
_adjustment variables_. Hence, we now introduce a causal model for the LTM
process to identify CAPE for LTM.

$\mathbf{U}_{v}$$\mathbf{X}_{v}$$Y_{v}^{0}$$Y_{v}^{1}$$Y_{v}^{2}$$Y_{v}^{T}$$\theta_{v}$$\mathbf{X}_{N}$$I_{v}^{0}$$I_{v}^{1}$$I_{v}^{2}$$I_{v}^{T}$

Figure 2: A causal model of peer effects for a node $v$ that reflects the LTM
process.

To estimate CAPE in LTM, we develop a Structural Causal Model (SCM) of
diffusion that corresponds to the LTM process. SCMs are graphical models that
encode cause-effect relationships between variables and allow for the
identification of effects (Pearl 2009).

Figure 2 shows the causal model of diffusion we develop to capture the LTM
process. Here, $\mathbf{X}_{N}$ represents a _set_ of neighbor features, and
$I_{v}^{t}$ is the random variable of _activation influence_ on the outcomes
of $v$’s neighbors. An arrow from one variable to another denotes a cause-
effect relationship. Contagion from node $v$’s neighbors are captured through
the arrows from $I_{v}^{t}$ to $Y_{v}^{t+1}$.

In order to _identify_ the effect of activation influence on a node’s own
activation (e.g., the effect of $I_{v}^{1}$ on $Y_{v}^{2}$, marked in blue in
Figure 2), we need to find the correct _adjustment set_ that blocks all back-
door paths, based on the back-door criterion (Pearl 2009; Arbour, Garant, and
Jensen 2016). In our example, three arrows are going into the treatment
$I_{v}^{1}$: $Y_{v}^{0}$, $I_{v}^{0}$ and $\mathbf{X}_{N}$, opening potential
backdoor paths. Based on this, a sufficient set to block all backdoor paths is
the set of $\mathbf{X}_{v}$, the previous outcome $Y_{v}^{1}$, and the
threshold $\theta_{v}$. In the CAPE formula defined in Equation (6), this is
the heterogeneous subgroup $\mathbf{X}=\mathbf{x}$ that we need to consider.
The special case is the initial peer effect of $I_{v}^{0}$ on $Y_{v}^{1}$,
where there are no backdoor paths from $I_{v}^{0}$ to $Y_{v}^{1}$. In this
case, we can estimate the effect of $I_{v}^{0}$ without adjustment.

To match the characteristics of LTM, we define a functional form of the
outcome of node $v$. The outcome $Y_{v}^{t+1}$ can be represented as a
function of $\big{(}\mathbf{X}_{v},Y_{v}^{t},I_{v}^{t},\theta_{v}\big{)}$. We
can define the functional form of $Y_{v}^{t+1}$ for $t\geq 1$ as:

$Y_{v}^{t+1}(\mathbf{X}_{v},Y_{v}^{t},I_{v}^{t},\theta_{v})=\begin{cases}Y_{v}^{t}&\text{if
$Y_{v}^{t}=1$,}\\\ \mathbbm{1}\big{[}I_{v}^{t}\geq\theta_{v}\big{]}&\text{if
$Y_{v}^{t}=0$.}\end{cases}$

This form correctly captures the LTM process: a node stays activated if
already activated. Otherwise, it activates based on its threshold and the
activations of neighbors. theore

### Individual threshold estimation for LTM

To estimate CAPE in LTM, we estimate the threshold $\theta$ assuming that
$I^{t}$ is known. We first map the problem of estimating thresholds to the
problem of estimating triggers (Tran and Zheleva 2019) in the context of
heterogeneous peer effects. This mapping allows us to adapt causal trees to
the problem of estimating thresholds and opens avenues for further algorithm
development.

Since the peer influence, $I_{v}^{t}$, is a continuous value, it cannot be
modeled using traditional causal inference methods which consider binary
treatment values. We reformulate the peer effect as an effect of a binary
treatment: the expected difference when a node’s activation influence,
$I_{v}^{t}$, is above and below its threshold, $\theta_{v}$. Then our goal is
transformed into estimating the correct node threshold that correctly
identifies when a node is above and below its actual threshold. To do this, we
map the threshold estimation problem to the problem of estimating triggers for
heterogeneous effects (Tran and Zheleva 2019).

A _trigger_ is defined as the minimum amount of treatment necessary to change
an outcome (Tran and Zheleva 2019). Some examples of triggers are the minimum
number of days a patient needs to take a medicine to be cured or a minimum
discount needed for a customer to buy a product. For the problem of threshold
estimation, our causal question is: “what is the minimum number of activated
neighbors that can cause a node with certain attributes to become activated”?
Define a trigger, $\hat{\theta}$, with sets of potential outcomes above and
below the trigger:

$\displaystyle\mathcal{Y}^{t}(I^{t}\geq\hat{\theta})=\\{Y_{v}^{t}\mid
I_{v}^{t}\geq\hat{\theta}\\}$ (10)
$\displaystyle\mathcal{Y}^{t}(I^{t}<\hat{\theta})=\\{Y_{v}^{t}\mid
I_{v}^{t}<\hat{\theta}\\}$ (11)

Then, we can define the average peer effect with a global trigger,
$\hat{\theta}$ as:
$E[\mathcal{Y}^{t}(I^{t}\geq\hat{\theta})-\mathcal{Y}^{t}(I^{t}<\hat{\theta})]$,
which is the case when a global threshold is defined for LTM.

Since we are interested in individual node thresholds, we estimate
_heterogeneous triggers_ rather than global triggers. We can now define CAPE
with a trigger:

$\mathcal{P}(\mathbf{x})=E[\mathcal{Y}(I^{t}\geq\hat{\theta})-\mathcal{Y}(I^{t}<\hat{\theta})\mid\mathbf{X}=\mathbf{x},Z=z].$
(12)

Estimating CAPE with a trigger translates to estimating the effect of
influence crossing the trigger. The implication of mapping the node threshold
estimation to a trigger-based heterogeneous treatment effect estimation
problem is that an accurate trigger is the best estimate of the node
threshold.

######  Theorem 1.

For any node $v\in\bm{V}$, let $\theta_{v}$ be $v$’s true threshold. Then
$\hat{\theta}_{v}$ that maximizes the CAPE with a trigger

$\operatorname*{argmax\,\,}_{\hat{\theta}_{v}}E[\mathcal{Y}(I^{t}_{v}\geq\hat{\theta}_{v})-\mathcal{Y}(I^{t}_{v}<\hat{\theta}_{v})\mid\mathbf{X}_{v},Z_{v}],$
(13)

provides the best estimate of the node threshold, $\theta_{v}$.

The proof is in the Appendix111https://github.com/edgeslab/hpe-ltm.

#### Trigger-based Causal Trees

Adapting trigger-based causal tree methods to find individual thresholds is
straightforward. Causal trees work similarly to decision trees in that they
greedily split using a partition function. The main difference between causal
trees and decision trees is that the goal of a causal tree is to estimate CATE
for different populations of individuals rather than to predict a label.
Causal trees work by maximizing a partition function based on the difference
of mean in outcomes while keeping low variance outcomes (Athey and Imbens
2016). The intuition is that finding splits that maximize the difference in
means finds the heterogeneity in effect.

Tran and Zheleva developed a trigger-based causal tree method called CT-HV
(Tran and Zheleva 2019). CT-HV works by introducing a validation set for
generalizing CATE estimations and reducing variance in outcomes. In order to
learn triggers, an additional search is done at each tree split to find the
trigger that maximizes the effect estimation for that split. We can utilize
CT-HV for threshold estimation through CAPE with a trigger. Instead of an
individual-level treatment, we use the activation influence, $I^{t}$, to
estimate triggers. Details of their algorithm are in the Appendix.

#### ST-Learner

Now we present a novel causal _meta-learner_ for estimating triggers for
heterogeneous effects. A causal meta-learner uses the outputs of base learners
for estimating effects. A base learner can be any regression or classification
method, such as Linear Regression or Decision Trees. Künzel et al. name two
commonly used meta-learners in causal inference literature, S- and T-learners,
and develop a meta-learner called X-Learner (Künzel et al. 2019). However,
none of them consider the problem of trigger-based HTE and are not directly
suitable for solving the problem of threshold estimation. Following the idea
of causal meta-learners, we propose a novel algorithm, _ST-Learner_ 222The “S”
refers to a _single_ learner and “T” is for _triggers_ , following the naming
scheme..

The goal of our ST-Learner is to learn the trigger that maximizes the effect.
A base learner is trained to predict the outcome $(Y_{v})$ given all node
features $(\mathbf{X}_{v})$ and the activation influence $(I_{v})$:
$E[Y_{v}|\mathbf{X}_{v},I_{v}]$. Once a learner has been built on the training
data (i.e., known activations), the next step is to estimate the trigger for
each node. The estimation considers different possible values of $I_{v}$ and
consists of two steps: outcome prediction and trigger estimation.

Let $\bm{\beta}=\\{\beta_{1},\beta_{2},\dots,\beta_{n}\\}$ be all the possible
treatment values (i.e., activation influence levels $I_{v})$ in the data.
Define $\Theta=\\{r_{1},r_{2},\dots,r_{m}\\}$, $m\leq n$, to be a set of
triggers. We use $r_{i}$ to refer to any potential trigger, while $\theta_{v}$
is a _node’s_ individual trigger. In practice, we can consider all potential
values of treatment in the training data to be possible triggers, or we can
consider a discretization of these values ($m\leq n$). For any node, we can
compute $n$ predictions for all possible treatment values
$\beta_{k}\in\bm{\beta}$: ${E}[Y_{v}|\mathbf{X}_{v},I=\beta_{k}]$.

With the predictions based on different activation influence levels, we can
estimate CAPE with a trigger. Let $\Theta_{i}^{1}=\\{\beta_{k}:\beta_{k}\geq
r_{i}\\}$ and $\Theta_{i}^{0}=\\{\beta_{j}:\beta_{j}<r_{i}\\}$ be the set of
treatment values above and below some trigger $r_{i}$. For an arbitrary
trigger $r_{i}$, we get the trigger estimation for ST-Learner:

$\displaystyle\operatorname*{argmax\,\,}_{r_{i}}\hat{\mathcal{P}}_{s}(\mathbf{x})$
$\displaystyle=\frac{1}{|\Theta_{i}^{1}|}\sum_{\beta_{k}\geq
r_{i}}{E}[Y_{v}|\mathbf{X}_{v}=\mathbf{x},I=\beta_{k}]$
$\displaystyle-\frac{1}{|\Theta_{i}^{0}|}\sum_{\beta_{j}<r_{i}}{E}[Y_{v}|\mathbf{X}_{v}=\mathbf{x},I=\beta_{j}].$

We use the base learner to estimate the outcomes above and below the trigger
by taking average outcomes. This results in a time complexity of $O(|B|\cdot
L)$ where $O(L)$ is the complexity of prediction for a base learner.
Pseudocode for ST-Learner is available in the Appendix.

## Experiments

We study the effectiveness of our proposed methods using four synthetic
network generation models and three real-world datasets. We consider three
tasks: threshold prediction, activated node prediction, and diffusion size
prediction.

### Experimental setup

Our main goal in this work is to estimate node thresholds (task 1) and
understand their role in downstream tasks, such as activated node prediction
(task 2) and diffusion size prediction (task 3). Comparing LTM to other
diffusion models (e.g., cascade models) is outside the scope of this work. To
compare the effectiveness of threshold prediction on these three tasks, we
estimate the node thresholds given the data at a snapshot (e.g., end of Jan
2016). A snapshot is the current structure and activations at network time
$t$. Each model predicts thresholds based on individual snapshots of data. For
example, models learn thresholds using the network snapshot at the end of Jan
2016 ($\mathcal{D}_{0}$) and estimate node thresholds for all inactivated
nodes in the test data. In this way, we also assess how having snapshots with
more activations, and thus more training data, impacts predictions.

(a) Erdos-Renyi

(b) Pref. attachment

(c) Forest fire

(d) Watts-Strogatz

Figure 3: MSE of threshold prediction for the _Linear_ setups over different
parameters.

We use four baseline threshold estimation methods, two of which are based on a
recent LTM survey (Talukder et al. 2019). _Heuristic Expected_ computes the
same threshold for all nodes in the network. _Heuristic Individual_ estimates
a range of values and samples individual node thresholds randomly from that
range. We also employ a baseline called _Random_ , that assign node thresholds
uniformly random from 0 to 1. In addition, we use a more practical baseline
for threshold estimation using _Linear Regression_. To get labels for Linear
Regression, we take all activated nodes in the network and estimate their
threshold by computing the proportion of activated neighbors. All baseline
methods focus on node threshold prediction for LTM.

We compare _Causal Tree_ and _ST-Learner_ to the baseline methods. For _ST-
Learner_ , we use two base learners, Linear Regression (ST-LR) and Decision
Trees (ST-DT). For neighbor influences, we use degree centrality:
$w_{uv}=1/|N(v)|$, where $|N(v)|$ is the degree of node $v$.

#### Task 1: node threshold prediction

For the first task, the goal is to accurately predict individual node
thresholds. Given a snapshot of the data at time $t$, we estimate a threshold
for each node. To evaluate the effectiveness of each model, we use the average
mean squared error of predicted thresholds across all snapshots.

$\text{MSE}=\frac{1}{T\cdot|V|}\sum_{t=1}^{T}\sum_{v}^{|V|}\Big{(}\theta_{v}-\hat{\theta}_{v}^{(t)}\Big{)}^{2}.$
(14)

This task is evaluated based on synthetic data because the true thresholds are
unknown in real-world data.

#### Task 2: activated node prediction

In the second task, we are interested in how well each model can predict which
specific nodes activate at each time step. To do this, we use the predicted
node thresholds and simulate a diffusion according to LTM. Let
$\mathcal{D}_{t}$ be the set of activated nodes at time $t$ and
$\hat{\mathcal{D}}_{t}$ be the set of predicted activated nodes at time $t$.
We compute the average Jaccard index on all sets as:

$J=\frac{1}{T}\sum_{t=0}^{T}\frac{|\mathcal{D}_{t}\cap\hat{\mathcal{D}}_{t}|}{|\mathcal{D}_{t}\cup\hat{\mathcal{D}}_{t}|}.$
(15)

The Jaccard index shows how accurate a model’s estimated thresholds are for
predicting activated nodes.

#### Task 3: diffusion size prediction

The final task is diffusion size prediction. Using the simulated diffusion, we
plot the number of predicted activations and compare them to the true number
of activations at each time step. This task shows how node threshold
prediction impacts diffusion size prediction. We refer to the diffusion size
as the _reach_.

### Datasets

Table 1: Our models predict the specific nodes that will activate most
accurately based on highest average Jaccard index.

| Linear setup | Quadrant setup | Real-world data
---|---|---|---
Threshold prediction method |  Erdos-Renyi |  Pref. Attach. |  Forest Fire |  Watts-Strogatz |  Erdos-Renyi |  Pref. Attach. |  Forest Fire |  Watts-Strogatz |  Hateful 2016 |  Hateful 2017 |  Cannabis |  Higgs
Random |  0.6981 |  0.7986 |  0.9094 |  0.8082 |  0.9830 |  0.8319 |  0.9187 |  0.7877 |  0.2450 |  0.2262 |  0.2033 |  0.6287
Heuristic Expected |  0.6268 |  0.7415 |  0.9649 |  0.7557 |  0.9783 |  0.7736 |  0.9630 |  0.7082 |  0.2608 |  0.2300 |  0.2433 |  0.6429
Heuristic Individual |  0.9127 |  0.8728 |  0.9103 |  0.8398 |  0.9785 |  0.7801 |  0.9692 |  0.7196 |  0.2617 |  0.2301 |  0.2691 |  0.6436
Linear Regression |  0.9082 |  0.8589 |  0.9454 |  0.8306 |  0.9891 |  0.8532 |  0.9494 |  0.7940 |  0.1522 |  0.2297 |  0.3674 |  0.6449
Causal Tree |  0.9599 |  0.9357 |  0.9509 |  0.8383 |  0.9978 |  0.9559 |  0.9638 |  0.8542 |  0.7207 |  0.5793 |  0.7830 |  0.9427
ST-Learner |  0.9658 |  0.9471 |  0.9622 |  0.8618 |  0.9972 |  0.9508 |  0.9526 |  0.8059 |  0.6877 |  0.5785 |  0.7564 |  0.9457
ST-Learner (Tree) |  0.9813 |  0.9777 |  0.9636 |  0.8669 |  0.9983 |  0.9801 |  0.9663 |  0.8672 |  0.7128 |  0.5988 |  0.7811 |  0.9734

We evaluate task 1 under multiple synthetic network scenarios and tasks 2 and
3 with three real-world datasets.

Synthetic datasets. We study four graph generation models: Erdos-Renyi (Erdős
and Rényi 1960), preferential attachment (Barabási and Albert 1999), forest
fire (Leskovec, Kleinberg, and Faloutsos 2005), and Watts-Strogatz (Watts and
Strogatz 1998) and two threshold generation models. For each set of {network
generation model, network parameter, threshold generation model}, we run 10
simulations and report average results. We set the number of nodes to 1000,
and for each node, we randomly generate 100 node attributes from a Gaussian,
$N(0,1)$.

In Erdos-Renyi, we vary the probability of edge creation $p$ from $0.05$ to
$0.5$. In preferential attachment, the number of new attachments $k$ from $1$
to $50$. For forest fire networks, we fix the backward probability of an edge
to $0.1$ and vary the forward probability of an edge $f$ from $0.05$ to $0.5$.
For the Watts-Strogatz networks, we fix the probability of rewiring an edge to
$0.1$, and vary nearest neighbors $k$ from $1$ to $50$.

In the first threshold generation, we use a random linear regression model
with $10$ of $100$ attributes and normalize the output to be the user
threshold, called the _Linear_ setup. In the second threshold generation, we
use $2$ features to separate thresholds into four quadrants, where each
quadrant receives a single threshold uniformly from $U(0,1)$, called the
_Quadrant_ setup. 50 nodes are randomly activated, and diffusion events are
generated based on LTM for $8$ time steps.

Hateful Users. The Hateful Users dataset is a retweet network from Twitter,
with 200 most recent tweets for each user (Ribeiro et al. 2018). Each user is
represented by the average _Empath_ category based on their tweets (Fast,
Chen, and Bernstein 2016). Empath captures a wide variety of topics such as
violence, fear, and warmth. A sample of users were selected to be annotated as
_hateful_ or _not hateful_ and the rest were predicted based on the history of
tweets using the methodology in (Ribeiro et al. 2018).

We estimate node thresholds for how “hatefulness” spreads through the network,
where being activated means you change from _not hateful_ to _hateful_. We
consider a month to month diffusion: how “hatefulness” diffuses on a month to
month basis. We look at two time periods: Jan 2016 to Dec 2016 and Jan 2017 to
Oct 2017.

Cannabis. The Cannabis dataset is a follower network. The dataset covers all
users who tweet about both cannabis and the e-cigarette Juul. From the users
who tweet about Juul, we identify those users who also tweet about cannabis or
marijuana. Empath categories are used as attributes. We estimate thresholds
for how cannabis tweets spread. Activation means an individual tweets a
cannabis related tweet. We consider the period between Jan 2017 to Dec 2017.

Higgs Boson. This dataset is based on the announcement of the Higgs-boson like
particle at CERN on July 4, 2012. The dataset was collected between July 1 and
July 8 of 2012 and is a follower network on Twitter (De Domenico et al. 2013).
We estimate node thresholds for the mention of the Higgs-boson discovery.
Activation means a user tweets about the Higgs-boson discovery. There are no
node attributes, so we construct features based on the graph. We use degree
centrality, both in-degree and out-degree, and counts of user and neighborhood
tweets. We consider hourly diffusions from July 4, 12:00am until July 8,
12:00am.

(a) Hateful Users

(b) Cannabis

(c) Higgs

Figure 4:  Comparison of diffusion size prediction on three real world
datasets. Our models have the closest estimation of reach over longer time
periods whereas the baselines incorrectly predict diffusion saturation in the
early stages.

### Task 1 results: node threshold prediction

Since real-world datasets do not have true node thresholds, we present results
on the synthetic datasets for this task. Figure 3 shows the MSE for each
threshold estimation method varied across network generation models for the
_Linear_ threshold generation method. _Quadrant_ dataset results are available
in the Appendix and show ST-DT performing the best overall. We see that the
node threshold estimators perform better across setups compared to baseline
methods that do not consider node threshold prediction. Our methods also
perform noticeably better than the Linear Regression baseline in forest fire,
and Watts-Strogatz. All individualized threshold prediction models perform
noticeably better than models that do not consider node threshold prediction.

From the results, ST-DT performs the best overall on the _Linear_ setup. Just
estimating thresholds using a regression method, like Linear Regression, does
not always perform better than other baselines. For example, Linear Regression
has very high error in Watts-Strogatz networks for smaller number of
neighbors. ST-DT also performs the best overall on the nonlinear threshold in
the _Quadrant_ setup, followed by Causal Tree.

### Task 2 results: activated node prediction

For task 2, we show results for both synthetic and real-world datasets. Table
1 shows the average Jaccard index for all datasets. Overall, our models are
able to achieve the best average Jaccard index across all datasets which is
especially pronounced in the real-world datasets. ST-DT achieves the highest
average Jaccard index in all synthetic datasets. In the real-world datasets,
we see that our models are able to achieve significantly higher average
Jaccard index compared to the baseline methods. One potential reason for the
discrepancy between the synthetic and real-world datasets is that synthetic
datasets have generated features and a specific linear or quadrant threshold
function. On the other hand, our formulation and estimation using individual
level features can better capture correct node thresholds for LTM for more
accurate prediction in real-world scenarios.

### Task 3 results: diffusion size prediction

We show results on diffusion size prediction on one network snapshot for real-
world datasets in Figure 4. A good model should have a curve that is close to
the true diffusion curve (in blue). Other snapshots can be found in the
Appendix and are consistent with the results here.

Hateful Users dataset. Figure 4(a) shows the diffusion size predictions for
the Hateful Users dataset from Jan 2016 to Dec 2016. Our models initially
overestimate reach, but are able to predict close to the final reach. With
more training data they estimate reach more accurately. In this dataset,
Linear Regression significantly overestimates diffusion prediction in all
snapshots. In the 2017 dataset we notice the same trends as in the 2016
dataset, so we omit the plots.

Cannabis dataset. Figure 4(b) shows diffusion size predictions on the Cannabis
dataset from Jan 2017 to Dec 2017. In contrast to the Hateful Users dataset,
all models predict reach that saturates after a number of time steps. One
reason for this could be the sparsity of edges. Additionally, different parts
of the network may not necessarily interact with each other, which results in
disjoint subgraphs.

Higgs dataset. Figure 4(c) shows the diffusion size predictions for the Higgs
dataset. Specifically, we start at noon, July 4th and increase the starting
snapshots by 12 hours each time. We see that the ST-DT performs the best
overall among our methods, and our methods perform significantly better than
baselines. Additionally, the baselines significantly overestimate the reach
predictions.

## Conclusion

In this work, we proposed a causal inference approach for estimating node
thresholds in the Linear Threshold Model (LTM). We defined a new concept of
heterogeneous peer effect estimation, and developed a structural causal model
for LTM to identify and estimate peer effects. We developed a new meta-
learner, the ST-Learner, and adapted trigger-based causal trees to solve the
threshold estimation problem through heterogeneous peer effects. Our results
on synthetic and real-world datasets showed our models are able to estimate
individualized thresholds from data better than baseline methods, and can
produce more accurate sets of activated nodes and diffusion size predictions
in the context of LTM, especially for real-world data. A fruitful avenue of
research would be to develop models that estimate edge influence weights, or
to jointly learn edge influence weights and thresholds through a causal
inference lens.

## Acknowledgments

This material is based on research sponsored in part by the Defense Advanced
Research Projects Agency (DAPRA) under contract numbers HR00111990114 and
HR001121C0168 and the National Science Foundation under grant No. 2047899. The
views and conclusions contained herein are those of the authors and should not
be interpreted as necessarily representing the official policies, either
expressed or implied, of DARPA or the U.S. Government. The U.S. Government is
authorized to reproduce and distribute reprints for governmental purposes
notwithstanding any copyright annotation therein.

## References

  * Arbour, Garant, and Jensen (2016) Arbour, D.; Garant, D.; and Jensen, D. 2016.  Inferring network effects from observational data.  In _Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining_ , 715–724.
  * Athey and Imbens (2016) Athey, S.; and Imbens, G. 2016.  Recursive partitioning for heterogeneous causal effects.  _Proceedings of the National Academy of Sciences_ , 113(27).
  * Barabási and Albert (1999) Barabási, A.-L.; and Albert, R. 1999.  Emergence of scaling in random networks.  _science_ , 286(5439): 509–512.
  * Berry et al. (2019) Berry, G.; Cameron, C. J.; Park, P.; and Macy, M. 2019.  The opacity problem in social contagion.  _Social Networks_ , 56: 93–101.
  * Bhattacharya, Malinsky, and Shpitser (2020) Bhattacharya, R.; Malinsky, D.; and Shpitser, I. 2020.  Causal inference under interference and network uncertainty.  In _Uncertainty in Artificial Intelligence_ , 1028–1038. PMLR.
  * Bourigault, Lamprier, and Gallinari (2016) Bourigault, S.; Lamprier, S.; and Gallinari, P. 2016.  Representation Learning for Information Diffusion Through Social Networks: An Embedded Cascade Model.  In _Proceedings of the Ninth ACM International Conference on Web Search and Data Mining_ , WSDM ’16, 573–582. New York, NY, USA: ACM.  ISBN 978-1-4503-3716-8.
  * Chen, Yuan, and Zhang (2010) Chen, W.; Yuan, Y.; and Zhang, L. 2010.  Scalable influence maximization in social networks under the linear threshold model.  In _2010 IEEE international conference on data mining_ , 88–97. IEEE, IEEE.
  * Clifford and Sudbury (1973) Clifford, P.; and Sudbury, A. 1973.  A model for spatial conflict.  Oxford University Press.
  * De Domenico et al. (2013) De Domenico, M.; Lima, A.; Mougel, P.; and Musolesi, M. 2013.  The anatomy of a scientific rumor.  _Scientific reports_ , 3: 2980.
  * Erdős and Rényi (1960) Erdős, P.; and Rényi, A. 1960.  On the evolution of random graphs.  _Publ. Math. Inst. Hung. Acad. Sci_ , 5(1): 17–60.
  * Fast, Chen, and Bernstein (2016) Fast, E.; Chen, B.; and Bernstein, M. S. 2016.  Empath: Understanding topic signals in large-scale text.  In _Proceedings of the 2016 CHI Conference on Human Factors in Computing Systems_ , 4647–4657. ACM.
  * Fatemi and Zheleva (2020) Fatemi, Z.; and Zheleva, E. 2020.  Minimizing Interference and Selection Bias in Network Experiment Design.  In _Proceedings of the International AAAI Conference on Web and Social Media_ , volume 14, 176–186.
  * Fowler and Christakis (2008) Fowler, J. H.; and Christakis, N. A. 2008.  Dynamic spread of happiness in a large social network: longitudinal analysis over 20 years in the Framingham Heart Study.  _Bmj_ , 337.
  * Goldenberg, Libai, and Muller (2001) Goldenberg, J.; Libai, B.; and Muller, E. 2001.  Talk of the network: A complex systems look at the underlying process of word-of-mouth.  _Marketing letters_ , 12(3): 211–223.
  * Goyal, Bonchi, and Lakshmanan (2010) Goyal, A.; Bonchi, F.; and Lakshmanan, L. V. 2010.  Learning influence probabilities in social networks.  In _Proceedings of the third ACM international conference on Web search and data mining_ , 241–250. ACM.
  * Goyal, Lu, and Lakshmanan (2011) Goyal, A.; Lu, W.; and Lakshmanan, L. V. 2011.  Simpath: An efficient algorithm for influence maximization under the linear threshold model.  In _2011 IEEE 11th international conference on data mining_ , 211–220. IEEE.
  * Granovetter (1978) Granovetter, M. 1978.  Threshold models of collective behavior.  _American journal of sociology_ , 83(6): 1420–1443.
  * Guo, Li, and Liu (2020) Guo, R.; Li, J.; and Liu, H. 2020.  Learning individual causal effects from networked observational data.  In _Proceedings of the 13th International Conference on Web Search and Data Mining_ , 232–240.
  * Holley, Liggett et al. (1975) Holley, R. A.; Liggett, T. M.; et al. 1975.  Ergodic theorems for weakly interacting infinite systems and the voter model.  _The annals of probability_ , 3(4): 643–663.
  * Kalimeris et al. (2018) Kalimeris, D.; Singer, Y.; Subbian, K.; and Weinsberg, U. 2018.  Learning diffusion using hyperparameters.  In _International Conference on Machine Learning_ , 2425–2433.
  * Kempe, Kleinberg, and Tardos (2003) Kempe, D.; Kleinberg, J.; and Tardos, É. 2003.  Maximizing the spread of influence through a social network.  In _Proceedings of the ninth ACM SIGKDD international conference on Knowledge discovery and data mining_ , 137–146. ACM.
  * Kermack and McKendrick (1927) Kermack, W. O.; and McKendrick, A. G. 1927.  A contribution to the mathematical theory of epidemics.  _Proceedings of the royal society of london. Series A, Containing papers of a mathematical and physical character_.
  * Künzel et al. (2019) Künzel, S. R.; Sekhon, J. S.; Bickel, P. J.; and Yu, B. 2019.  Metalearners for estimating heterogeneous treatment effects using machine learning.  _Proceedings of the national academy of sciences_ , 116(10): 4156–4165.
  * Leskovec, Kleinberg, and Faloutsos (2005) Leskovec, J.; Kleinberg, J.; and Faloutsos, C. 2005.  Graphs over time: densification laws, shrinking diameters and possible explanations.  In _Proceedings of the eleventh ACM SIGKDD international conference on Knowledge discovery in data mining_ , 177–187. ACM.
  * Li et al. (2018) Li, Y.; Fan, J.; Wang, Y.; and Tan, K.-L. 2018.  Influence maximization on social graphs: A survey.  _IEEE Transactions on Knowledge and Data Engineering_ , 30(10): 1852–1872.
  * Manski (1993) Manski, C. F. 1993.  Identification of endogenous social effects: The reflection problem.  _The review of economic studies_ , 60(3): 531–542.
  * Nikolakopoulos et al. (2019) Nikolakopoulos, A. N.; Berberidis, D.; Karypis, G.; and Giannakis, G. B. 2019.  Personalized Diffusions for Top-n Recommendation.  In _Proceedings of the 13th ACM Conference on Recommender Systems_. ACM.
  * Ogburn and VanderWeele (2014) Ogburn, E. L.; and VanderWeele, T. J. 2014.  Causal diagrams for interference.  _Statistical science_ , 29(4): 559–578.
  * Pearl (2009) Pearl, J. 2009.  _Causality_.  Cambridge university press.
  * Ribeiro et al. (2018) Ribeiro, M. H.; Calais, P. H.; Santos, Y. A.; Almeida, V. A.; and Meira Jr, W. 2018\.  Characterizing and detecting hateful users on twitter.  In _Twelfth International AAAI Conference on Web and Social Media_.
  * Rogers (2010) Rogers, E. M. 2010.  _Diffusion of innovations_.  Simon and Schuster.
  * Rubin (1978) Rubin, D. B. 1978.  Bayesian inference for causal effects: The role of randomization.  _The Annals of statistics_ , 34–58.
  * Saito, Nakano, and Kimura (2008) Saito, K.; Nakano, R.; and Kimura, M. 2008.  Prediction of information diffusion probabilities for independent cascade model.  In _International conference on knowledge-based and intelligent information and engineering systems_ , 67–75. Springer.
  * Shalizi and Thomas (2011) Shalizi, C. R.; and Thomas, A. C. 2011.  Homophily and contagion are generically confounded in observational social network studies.  _Sociological methods & research_, 40(2): 211–239.
  * Soni, Ramirez, and Eisenstein (2019) Soni, S.; Ramirez, S. L.; and Eisenstein, J. J. 2019.  Detecting Social Influence in Event Cascades by Comparing Discriminative Rankers.  In _The 2019 ACM SIGKDD Workshop on Causal Discovery_ , 78–99.
  * Talukder et al. (2019) Talukder, A.; Alam, M. G. R.; Tran, N. H.; Niyato, D.; Park, G. H.; and Hong, C. S. 2019.  Threshold Estimation Models for Linear Threshold-Based Influential User Mining in Social Networks.  _IEEE Access_ , 7: 105441–105461.
  * Tran and Zheleva (2019) Tran, C.; and Zheleva, E. 2019.  Learning Triggers for Heterogeneous Treatment Effects.  In _Proceedings of the AAAI Conference on Artificial Intelligence_.
  * Veitch, Wang, and Blei (2019) Veitch, V.; Wang, Y.; and Blei, D. 2019.  Using Embeddings to Correct for Unobserved Confounding in Networks.  In _Advances in Neural Information Processing Systems 32_ , 13792–13802.
  * Watts and Strogatz (1998) Watts, D. J.; and Strogatz, S. H. 1998.  Collective dynamics of ‘small-world’networks.  _nature_ , 393(6684): 440.

## Appendix

## Individual threshold estimation for LTM

The implication of mapping the node threshold estimation to a trigger-based
heterogeneous treatment effect estimation problem is that an accurate trigger
correctly estimates the node threshold.

######  Theorem 1.

For any node $v\in\bm{V}$, let $\theta_{v}$ be $v$’s true threshold. Then
$\hat{\theta}_{v}$ that maximizes the CAPE with a trigger in eq (12):

$\operatorname*{argmax\,\,}_{\hat{\theta}_{v}}E[\mathcal{Y}(I^{t}_{v}\geq\hat{\theta}_{v})-\mathcal{Y}(I^{t}_{v}<\hat{\theta}_{v})\mid\mathbf{X}_{v},Z_{v}],$
(16)

provides the best estimate of the node threshold, $\theta_{v}$.

###### Proof.

Suppose $|N(v)|=n$ and $w_{uv}$ be an arbitrary weight such that
$\sum_{u}w_{uv}=1$. Let the true threshold of node $v$ be $\theta_{v}$. Define
$\mathbf{A_{v}}=\mathbf{a_{i}}$ to be an assignment of activations of
neighbors of $v$ and $W_{v}(\alpha_{i})$ to be the outcome for the activated
set $\mathbf{a_{i}}$:

$\displaystyle W_{v}(\alpha)$ $\displaystyle=\begin{cases}0&\text{if
$I_{v}(\mathbf{a_{i}})<\theta_{v}$},\\\ 1&\text{if
$I_{v}(\mathbf{a_{i}})\geq\theta_{v}$}.\end{cases}$ (17)

Where $I_{v}$ is the activation influence. Let the set of potential outcomes,
$\mathcal{Y}_{v}$, below and above the estimated trigger $\hat{\theta_{v}}$
be:

$\displaystyle\mathcal{Y}_{v}(I_{v}<\hat{\theta}_{v})$
$\displaystyle=\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})<\hat{\theta_{v}}\\}$ (18)
$\displaystyle\mathcal{Y}_{v}(I_{v}\geq\hat{\theta}_{v})$
$\displaystyle=\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})\geq\hat{\theta}_{v}\\}.$ (19)

Let $N_{0}$ and $N_{1}$ be the size of the sets
$\mathcal{Y}_{v}(I_{v}<\hat{\theta}_{v})$ and
$\mathcal{Y}_{v}(I_{v}\geq\hat{\theta}_{v})$, respectively. Define the set of
outcomes when node $v$ is not activated and activated with the true threshold
as $\mathcal{Z}_{v}(0)$ and $\mathcal{Z}_{v}(1)$:

$\displaystyle\mathcal{Z}_{v}(0)$ $\displaystyle=\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})<\theta_{v}\\},$ (20) $\displaystyle\mathcal{Z}_{v}(1)$
$\displaystyle=\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})\geq\theta_{v}\\}.$ (21)

Let $N_{\mathcal{Z}_{0}}$ and $N_{\mathcal{Z}_{1}}$ be the cardinalities of
set $\mathcal{Z}_{v}(0)$ and $\mathcal{Z}_{v}(1)$, respectively. The
difference between $\mathcal{Y}_{v}$ and $\mathcal{Z}_{v}$ is that
$\mathcal{Y}_{v}$ contains the true potential activations based on the
_estimated_ threshold and $\mathcal{Z}_{v}$ contains the true potential
activations based on the _true_ threshold. We define expectation over the sets
$\mathcal{Y}_{v}$ and $\mathcal{Z}_{v}$ as the mean over the set. There are
three cases for the estimated threshold, $\hat{\theta}_{v}$:

######  Case 1 ($\hat{\theta}_{v}<\theta_{v}$):

We compute the expected outcomes below and above the estimated trigger
$\hat{\theta}_{v}$:

$\displaystyle E[\mathcal{Y}_{v}(I_{v}<\hat{\theta}_{v})]$
$\displaystyle=E[\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})<\hat{\theta}_{v}\\}]$
$\displaystyle=\frac{1}{N_{0}}\sum_{i}W_{v}(\mathbf{a_{i}})=0,$ (22)

since there are no times when $W_{v}$ will be $1$ (activated) since all
outcomes below the estimated threshold $\hat{\theta}_{v}$ are also below the
true threshold $\theta_{v}$. The expected value above the estimated threshold
is:

$\displaystyle
E[\mathcal{Y}_{v}(I_{v}\geq\hat{\theta}_{v})]=E[\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})\geq\hat{\theta}_{v}\\}]$
$\displaystyle=\frac{1}{N_{1}}\sum_{i}W_{v}(\mathbf{a_{i}})=\frac{1}{N_{1}}\sum_{i}W_{v}(\mathbf{a_{i}})+\frac{1}{N_{1}}\sum_{j}W_{v}(\mathbf{a_{j}})$
$\displaystyle=\frac{1}{N_{1}}\cdot
N_{\mathcal{Z}_{1}}=\frac{N_{\mathcal{Z}_{1}}}{N_{1}}.$ (23)

The effect is the difference above and below the estimated trigger:

$\displaystyle
E[\mathcal{Y}_{v}(I_{v}\geq\hat{\theta_{v}})-\mathcal{Y}_{v}(I_{v}<\hat{\theta_{v}})]$
$\displaystyle=\frac{N_{\mathcal{Z}_{1}}}{N_{1}}-0=\frac{N_{\mathcal{Z}_{1}}}{N_{1}}<1.$
(24)

Note that $N_{\mathcal{Z}_{1}}<N_{1}$ since $N_{1}$ has cases above the
estimated threshold, which is smaller than the true threshold.

######  Case 2 ($\hat{\theta}_{v}>\theta_{v}$):

$\displaystyle E[\mathcal{Y}_{v}(I_{v}<\hat{\theta}_{v})]$
$\displaystyle=E[\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})<\hat{\theta}_{v}\\}]$
$\displaystyle=\frac{1}{N_{0}}\sum_{i}W_{v}(\mathbf{a_{i}})+\frac{1}{N_{0}}\sum_{j}W_{v}(\mathbf{a_{i}})$
$\displaystyle=0+\frac{1}{N_{0}}\cdot(N_{0}-N_{\mathcal{Z}_{0}})$
$\displaystyle=\frac{N_{0}-N_{\mathcal{Z}_{0}}}{N_{0}}.$ (25) $\displaystyle
E[\mathcal{Y}_{v}(I_{v}\geq\hat{\theta}_{v})]$
$\displaystyle=E[\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})\geq\hat{\theta}_{v}\\}]$
$\displaystyle=\frac{1}{N_{1}}\cdot N_{1}=1.$ (26) $\displaystyle
E[\mathcal{Y}_{v}(I_{v}\geq\hat{\theta_{v}})$
$\displaystyle-\mathcal{Y}_{v}(I_{v}<\hat{\theta_{v}})]=1-\frac{N_{0}-N_{\mathcal{Z}_{0}}}{N_{0}}<1.$
(27)

######  Case 3 ($\hat{\theta}_{v}=\theta_{v}$):

$\displaystyle E[\mathcal{Y}_{v}(I_{v}<\hat{\theta}_{v})]$
$\displaystyle=E[\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})<\hat{\theta}_{v}\\}]$
$\displaystyle=\frac{1}{N_{0}}\sum_{i}W_{v}(\mathbf{a_{i}})=0.$ (28)
$\displaystyle E[\mathcal{Y}_{v}(I_{v}\geq\hat{\theta}_{v})]$
$\displaystyle=E[\\{W_{v}(\mathbf{a_{i}})\mid
I_{v}(\mathbf{a_{i}})\geq\hat{\theta}_{v}\\}]$
$\displaystyle=\frac{1}{N_{1}}\sum_{i}W_{v}(\mathbf{a_{i}})=1.$ (29)
$\displaystyle E[\mathcal{Y}_{v}(I_{v}\geq\hat{\theta_{v}})$
$\displaystyle-\mathcal{Y}_{v}(I_{v}<\hat{\theta_{v}})]=1-0=1.$ (30)

The maximum causal effect only occurs when $\hat{\theta}_{v}=\theta_{v}$.
Therefore if we can estimate the trigger that maximizes CAPE we find the true
node threshold. As a result, we can use any trigger-based HTE estimation
method to estimate node thresholds. ∎

### Threshold estimation algorithms

#### Trigger-based causal trees

Here we describe the Causal Tree (CT-HV) algorithm proposed in (Tran and
Zheleva 2019). Let $\mathbf{X}^{\ell}$ be the features in partition $\ell$,
which represents a node in the tree which contains $N_{\ell}$ samples, and let
$\hat{\mu}_{1}(\ell)$ and $\hat{\mu}_{0}(\ell)$ be the mean of outcomes when
treated and non-treated. The estimate of CATE of any partition is
$\hat{\tau_{c}}(\mathbf{X}^{\ell})=\hat{\mu}_{1}(\ell)-\hat{\mu}_{0}(\ell)$.
Given a partition $\ell$ that needs to be partitioned further into two
children $\ell_{1},\ell_{2}$, the causal tree algorithm finds the split on
features that maximizes the weighted CATE in each child:

$\operatorname*{max\,\,}_{\ell_{1},\ell_{2}}N_{\ell_{1}}\cdot\hat{\tau_{c}}(\mathbf{X}^{\ell_{1}})+N_{\ell_{2}}\cdot\hat{\tau_{c}}(\mathbf{X}^{\ell_{2}}).$
(31)

We refer to the two quantities
$N_{\ell_{1}}\cdot\hat{\tau_{c}}(\mathbf{X}^{\ell_{1}})$ and
$N_{\ell_{2}}\cdot\hat{\tau_{c}}(\mathbf{X}^{\ell_{2}})$ as partition measures
for $\ell_{1}$ and $\ell_{2}$. In our work, we adapt the CT-HV algorithm (Tran
and Zheleva 2019) for the problem of node threshold estimation.

In order to learn triggers, an additional search is done at each split to find
the trigger that maximizes the effect estimation in each split. Let $F(\ell)$
represent the partition measure for CT-HV and let $\hat{\theta_{v}}$ be some
estimated trigger for partition $\ell$. Define $M_{1}(\ell,\hat{\theta_{v}})$
and $M_{0}(\ell,\hat{\theta_{v}})$ to be the expected mean outcomes when above
and below the trigger $\hat{\theta_{v}}$, and $F(\ell,\hat{\theta_{v}})$ be
the partition measure with trigger $\hat{\theta_{v}}$:
$F(\ell,\hat{\theta_{v}})=M_{1}(\ell,\hat{\theta_{v}})-M_{0}(\ell,\hat{\theta_{v}})=E[\mathcal{Y}_{\ell}(I\geq\hat{\theta_{v}})-\mathcal{Y}_{\ell}(I<\hat{\theta_{v}})].$
Then to find a trigger in partition $\ell$, we find the trigger that maximizes
the partition measure: $F(\ell,\hat{\theta_{v}})$. This results in a unique
trigger in each partition.

### ST-Learner

Pseudocode is provided for ST-Learner in Algorithm 1. Both the training and
prediction subroutine are included. In the train subroutine, we train a base
learner, $f$ to predict the outcome $Y$ using node features $\mathbf{X}$ and
the activation influence $I$. To predict a trigger for a new test example, we
compute and score the predicted outcome for all treatment values found in the
dataset. Then, for each potential trigger to consider, we compute an average
above and below that trigger, and find the trigger that maximizes the
difference in outcomes. The additional complexity added by the search for
triggers is on the order of potential treatment values $|B|$. For example,
suppose a base learner of Linear Regression requires $O(d)$ time for
prediction, then ST-Learner with Linear Regression will take $O(d\cdot|B|)$
time.

Algorithm 1 ST-Learner

1:Training set $(\mathbf{X},Y,I)$, test example $\mathbf{x}^{te}$, base
learner $f$, possible treatments $\bm{\beta}$, possible triggers $\Theta$

2:The causal effect estimate for $\mathbf{x}^{te}$

3:function Train($\mathbf{X},Y,I$)

4: Fit $f$ to predict $Y$, $f$: $(\mathbf{X},I)\rightarrow Y$

5: $\implies f(\mathbf{X},I)=E[Y\mid\mathbf{X},I]$

6:function Predict($\mathbf{x}^{te}$)

7: for each $\beta_{i}$ in $\bm{\beta}$ do

8: Compute and store $f(\mathbf{x}^{te},\beta_{i})$

9: max_trigger $=0$, max_effect $=0$

10: for each $r_{i}$ in $\Theta$ do

11: Compute $t_{1}=\frac{1}{|\Theta_{i}^{1}|}\sum_{\beta_{k}\geq
r_{i}}f(\mathbf{x}^{te},\beta_{k})$

12: Compute
$t_{0}=\frac{1}{|\Theta_{i}^{0}|}\sum_{\beta_{j}<r_{i}}f(\mathbf{x}^{te},\beta_{j})$

13: $e=t_{1}-t_{0}$

14: if $e>$ max_effect

15: max_trigger $=r_{i}$, max_effect $=e$

16: return max_trigger, max_effect

### Additional results

#### Task 1: node threshold prediction

(a) Erdos-Renyi

(b) Pref. attachment

(c) Forest fire

(d) Watts-Strogatz

(e) Erdos-Renyi

(f) Pref. attachment

(g) Forest fire

(h) Watts-Strogatz

Figure 5: MSE of threshold prediction for the _Linear_ (top row) and
_Quadrant_ (bottom row) setups over different parameters.

Figure 5 shows the change in MSE in the _Linear_ and _Quadrant_ setup. From
the results, we see that Causal Tree in general performs better on the
nonlinear threshold in the _Quadrant_ dataset, compared to the linearly
generated threshold.

### Task 3: diffusion size prediction

(a)

(b)

(c)

(d)

(e)

(f)

(g)

(h)

(i)

(j)

(k)

(l)

Figure 6:  Comparison of diffusion size prediction on three real world
datasets. Our models have the closest estimation of reach over longer time
periods whereas the baselines incorrectly predict diffusion saturation in the
early stages.

Figure 6 shows additional figures for the real-world datasets. Each figure
starts at a different network snapshot, where a snapshot is the current
structure and activations at time $t$. For example, Figure 6(b) learns on data
starting at the snapshot at the end of March 2016. A later snapshot means
there are more activations in the dataset, and thus more training data. As the
training data size increases, our methods can achieve better diffusion size
prediction.

Hateful Users dataset. Figures 6(a)-6(d) show the diffusion predictions for
the Hateful Users dataset from Jan 2016 to Dec 2016, with varying starting
points (Jan, Mar, May, Jul). We show diffusion simulations at four snapshots
before Dec 2016. Our models initially overestimate the diffusion size
prediction The first thing we notice is that our models initially
overestimates the reach the diffusion prediction, but is able to predict close
to the final diffusion amount. Additionally, we see that as we get more time
steps, our models obtains more accurate reach estimates. For example, in
Figures 6(a) and 6(b), where we start with information in January and March,
our models overestimates the prediction in the beginning. With more
information, the reach estimates are better, as in Figures 6(c) (starting in
May), and 6(d) (starting in July). In the 2017 dataset, we notice the same
trends as in the 2016 dataset: our models overestimates the reach prediction
in the beginning, but predicts close to the final true amount, so we omit the
plots.

Cannabis dataset. Figures 6(e)-6(h) show results for threshold estimation on
the Cannabis dataset from Jan 2017 to Dec 2017. Contrast to the Hateful Users
dataset, all models predict diffusion that saturate after a number of time
steps. Our models predict saturation closer to the final diffusion amount with
more training data.

Higgs dataset. Figure 6(i) shows the diffusion size predictions for the Higgs
dataset. Specifically, we start at noon, July 4th and increase the starting
snapshots by 12 hours each time. We see that the ST-Learner performs slightly
better than Causal Tree and both outperform the baselines by a significant
margin. Additionally, the baselines significantly overestimate the reach
predictions.


* [23] M. Jung, W. Saad, M. Debbah, and C. S. Hong, “On the optimality of reconfigurable intelligent surfaces (RISs): Passive beamforming, modulation, and resource allocation,” _IEEE Trans. Wireless Commun._ , vol. 20, no. 7, pp. 4347–4363, 2021.
  * [24] W. Mei and R. Zhang, “Performance analysis and user association optimization for wireless network aided by multiple intelligent reflecting surfaces,” _IEEE Trans. Wireless Commun._ , vol. 69, no. 9, pp. 6296–6312, 2021.
  * [25] H. Sun, Z. Zhao, H. Cheng, J. Lyu, X. Wang, Y. Zhang, and T. Q. Quek, “IRS-assisted RF-powered IoT networks: System modeling and performance analysis,” _IEEE Trans. Wireless Commun._ , 2023.
  * [26] Y. Fang, S. Atapattu, H. Inaltekin, and J. Evans, “Optimum reconfigurable intelligent surface selection for wireless networks,” _IEEE Trans. Commun._ , 2022.
  * [27] W. Khawaja, O. Ozdemir, Y. Yapici, F. Erden, and I. Guvenc, “Coverage enhancement for NLOS mmWave links using passive reflectors,” _IEEE Open J. Commun. Soc._ , vol. 1, pp. 263–281, 2020.
  * [28] T. Noh and J. Choi, “Cell-free MIMO systems powered by intelligent reflecting surfaces,” _IEEE Commun. Lett._ , 2022.
  * [29] E. Bjornson, O. Ozdogan, and E. G. Larsson, “Intelligent reflecting surface versus decode-and-forward: How large surfaces are needed to beat relaying?” _IEEE Wireless Commun. Lett._ , vol. 9, no. 2, pp. 244–248, Feb. 2020.
  * [30] Y. Liu, X. Liu, X. Mu, T. Hou, J. Xu, M. Di Renzo, and N. Al-Dhahir, “Reconfigurable intelligent surfaces: Principles and opportunities,” _IEEE Commun. Surv. Tutorials_ , vol. 23, no. 3, pp. 1546–1577, 2021.
  * [31] C. Pan, H. Ren, K. Wang, J. F. Kolb, M. Elkashlan, M. Chen, M. Di Renzo, Y. Hao, J. Wang, A. L. Swindlehurst _et al._ , “Reconfigurable intelligent surfaces for 6G systems: Principles, applications, and research directions,” _IEEE Commun. Mag._ , vol. 59, no. 6, pp. 14–20, 2021.
  * [32] S. W. H. Shah, A. N. Mian, S. Mumtaz, A. Al-Dulaimi, I. Chih-Lin, and J. Crowcroft, “Statistical QoS analysis of reconfigurable intelligent surface-assisted D2D communication,” _IEEE Trans. Veh. Technol._ , vol. 71, no. 7, pp. 7343–7358, 2022.
  * [33] W. Tang, X. Chen, M. Z. Chen, J. Y. Dai, Y. Han, S. Jin, Q. Cheng, G. Y. Li, and T. J. Cui, “On channel reciprocity in reconfigurable intelligent surface assisted wireless networks,” _IEEE Wireless Commun._ , vol. 28, no. 6, pp. 94–101, 2021.
  * [34] M. Cui, Z. Wu, Y. Lu, X. Wei, and L. Dai, “Near-field MIMO communications for 6g: Fundamentals, challenges, potentials, and future directions,” _IEEE Communications Magazine_ , vol. 61, no. 1, pp. 40–46, Jan. 2023. [Online]. Available: https://doi.org/10.1109/mcom.004.2200136
  * [35] M. Najafi, V. Jamali, R. Schober, and H. V. Poor, “Physics-based modeling and scalable optimization of large intelligent reflecting surfaces,” _IEEE Trans. Commun._ , vol. 69, no. 4, pp. 2673–2691, 2020.
  * [36] K. Dovelos, S. D. Assimonis, H. Quoc Ngo, B. Bellalta, and M. Matthaiou, “Intelligent reflecting surfaces at terahertz bands: Channel modeling and analysis,” in _ICC Workshops_. IEEE, Jun. 2021, pp. 1–6.
  * [37] Y. U. Ozcan, O. Ozdemir, and G. K. Kurt, “Reconfigurable intelligent surfaces for the connectivity of autonomous vehicles,” _IEEE Trans. Veh. Technol._ , vol. 70, no. 3, pp. 2508–2513, 2021.
  * [38] Ö. Özdogan, E. Björnson, and E. G. Larsson, “Intelligent reflecting surfaces: Physics, propagation, and pathloss modeling,” _IEEE Wireless Commun. Lett._ , vol. 9, no. 5, pp. 581–585, 2019.
  * [39] X. Hu, C. Zhong, and Z. Zhang, “Angle-domain intelligent reflecting surface systems: Design and analysis,” _IEEE Trans. Commun._ , vol. 69, no. 6, pp. 4202–4215, 2021.
  * [40] T. L. Nguyen, T. N. Do, G. Kaddoum, D. B. da Costa, and Z. J. Haas, “Channel characterization for RIS-aided terahertz communications: A stochastic approach,” _IEEE Wireless Commun. Lett._ , vol. 11, no. 9, pp. 1890–1894, 2022.
  * [41] S. W. Ellingson, “Path loss in reconfigurable intelligent surface-enabled channels,” in _2021 IEEE PIMRC_. IEEE, 2021, pp. 829–835.
  * [42] W. Ni, X. Liu, Y. Liu, H. Tian, and Y. Chen, “Resource allocation for multi-cell IRS-aided NOMA networks,” _IEEE Trans. Wireless Commun._ , vol. 20, no. 7, pp. 4253–4268, 2021.
  * [43] C. A. Balanis, _Advanced engineering electromagnetics_. John Wiley & Sons, 2012.
  * [44] J. Cui, Y. Liu, Z. Ding, P. Fan, and A. Nallanathan, “Optimal user scheduling and power allocation for millimeter wave NOMA systems,” _IEEE Trans. Wireless Commun._ , vol. 17, no. 3, pp. 1502–1517, 2017.
  * [45] A. Alizadeh and M. Vu, “Load balancing user association in millimeter wave MIMO networks,” _IEEE Trans. Wireless Commun._ , vol. 18, no. 6, pp. 2932–2945, 2019.
  * [46] S. Sun, F. Yang, and J. Song, “Sum rate maximization for intelligent reflecting surface-aided visible light communications,” _IEEE Communications Letters_ , vol. 25, no. 11, pp. 3619–3623, 2021.
  * [47] M. N. Jung, C. Kirches, and S. Sager, “On perspective functions and vanishing constraints in mixed-integer nonlinear optimal control,” in _Facets of combinatorial optimization_. Springer, 2013, pp. 387–417.
  * [48] S. Burer and A. N. Letchford, “Non-convex mixed-integer nonlinear programming: A survey,” _Surveys in Operations Research and Management Science_ , vol. 17, no. 2, pp. 97–106, 2012.
  * [49] F. Fang, H. Zhang, J. Cheng, and V. C. Leung, “Energy-efficient resource allocation for downlink non-orthogonal multiple access network,” _IEEE Trans. Commun._ , vol. 64, no. 9, pp. 3722–3732, 2016.
  * [50] Z. Wang, L. Liu, and S. Cui, “Channel estimation for intelligent reflecting surface assisted multiuser communications: Framework, algorithms, and analysis,” _IEEE Trans. Wireless Commun._ , vol. 19, no. 10, pp. 6607–6620, 2020.
  * [51] D. Gale and L. S. Shapley, “College admissions and the stability of marriage,” _The American Mathematical Monthly_ , vol. 69, no. 1, pp. 9–15, 1962.
  * [52] S. Shamaei, S. Bayat, and A. M. A. Hemmatyar, “Interference management in d2d-enabled heterogeneous cellular networks using matching theory,” _IEEE Trans. Mob. Comput._ , vol. 18, no. 9, pp. 2091–2102, 2018.
  * [53] H. Zhang, Y. Duan, K. Long, and V. C. Leung, “Energy efficient resource allocation in terahertz downlink NOMA systems,” _IEEE Trans. Commun._ , vol. 69, no. 2, pp. 1375–1384, 2020.
  * [54] Y. Xu and S. Mao, “User association in massive MIMO hetnets,” _IEEE Syst. J._ , vol. 11, no. 1, pp. 7–19, 2015.
  * [55] M. Rahim, G. Kaddoum, and T. N. Do, “Joint devices and IRSs association for terahertz communications in industrial IoT networks,” _IEEE Trans. Green Commun. Networking_ , pp. 1–1, 2023.
  * [56] C. A. Balanis, _Antenna theory: analysis and design_. Wiley-Interscience, 2005.
  * [57] S. Tarboush, H. Sarieddeen, H. Chen, M. H. Loukil, H. Jemaa, M.-S. Alouini, and T. Y. Al-Naffouri, “TeraMIMO: A channel simulator for wideband ultra-massive MIMO Terahertz communications,” _IEEE Trans. Veh. Technol._ , vol. 70, no. 12, pp. 12 325–12 341, Dec. 2021.


[a]Xuan Li

# Silicon detector R$\&$D for the future Electron-Ion Collider

###### Abstract

The high-luminosity high-energy Electron-Ion Collider (EIC) to be built at
Brookhaven National Laboratory (BNL) will provide a clean environment to study
several fundamental questions in the high energy and nuclear physics fields. A
high granularity and low material budget silicon vertex and tracking detector
is required to provide precise measurements of primary and displaced vertex
and track reconstruction with good momentum and spatial resolutions. The
reference design of the EIC silicon vertex and tracking detector includes the
Monolithic Active Pixel Sensor (MAPS) based vertex and tracking subsystem and
the AC-Coupled Low Gain Avalanche Diode (AC-LGAD) based outer tracker, and it
has the track reconstruction capability in the pseudorapidity region of -3.5
to 3.5 with full azimuthal coverage. Further detector geometry optimization
with a new magnet based on the EIC project detector reference design are being
performed by the newly formed ePIC collaboration. The latest ePIC vertex and
tracking detector geometry and its performance evaluated in simulation will be
presented. Details of the EIC silicon vertex and tracking detector R$\&$D,
which include the proposed detector technologies, prototype sensor
characterization and detector mechanical design will be discussed as well.

##  1 Introduction to the Electron-Ion Collider (EIC)

The future Electron-Ion Collider (EIC), which will be built at Brookhaven
National Laboratory, has received the CD1 approval to define its funding scope
in 2021 by the US Department of Energy (DOE). The EIC will operate
electron+proton ($e+p$) and electron+nucleus ($e+A$) collisions with the
electron beam energy at 2.5-18 GeV and the proton/nucleus beam energy at 41,
100-275 GeV. The Instantaneous luminosities of the EIC are around
$10^{33-34}~{}\rm{cm}^{-2}\rm{s}^{-1}$ and the bunch crossing rate is around
10 ns [1]. The EIC project detector, which locates at the 6 o’clock location
of the accelerator ring, is scheduled to start construction in 2025. The
reference design of the EIC project detector [2], which was selected in early
2022, consists of high granularity vertex, tracking, particle identification,
electromagnetic calorimeter and hadronic calorimeter subsystems. Figure 1
shows the bird-view of the EIC accelerator apparatus (left) and the EIC
project detector reference design (right). Following the EIC project detector
reference design decision, a new collaboration: ePIC has been formed to lead
the EIC project detector design optimization towards the CD3A approval, which
is scheduled in 2024.

Figure 1: The EIC accelerator apparatus (left) has been designed to host the
EIC project detector at the 6 o’clock position and potentially a second
detector at the 8 o’clock location. The detector layout of the EIC project
detector reference design is shown in the right.

##  2 EIC Project Detector Silicon Tracking Subsystem Design and Performance

The vertex and tracking subsystem of the EIC project detector reference design
utilizes the next generation 65 nm Monolithic Active Pixel Sensor (MAPS) [3],
Micro-Pattern Gaseous Detector (MPGD) [4] and AC coupled Low Gain Avalanche
Diode (AC-LGAD) [5] technologies. The EIC project detector design optimization
is underway by the ePIC collaboration considering the detector integration
needs and performance improvements. Two iterations of geometry optimizations
have been performed for the ePIC MAPS vertex and tracking detector design. The
updated ePIC MAPS vertex and tracking detector covers a relatively larger
active detector area with a larger number of channels than the EIC project
detector reference design. The geometry parameters of the current ePIC MAPS
vertex and tracking detector are shown in Table 1, Table 2 and Table 3.

The ePIC barrel MAPS vertex and tracking subsystem consists of 3 vertex layers
based on bent sensors and 2 sagitta layers of stitched sensor staves. The
radii of the MAPS vertex layers vary from 3.6 cm to 12.0 cm and the radii of
the MAPS sagitta layers are at 27 cm and 42 cm. The lengths of the barrel MAPS
vertex and tracking layers are from 27 cm to 84 cm. This design aims to obtain
low material budgets and good tracking momentum and spatial resolutions listed
in the EIC yellow report [1]. The ePIC MAPS hadron endcap detector (geometry
parameters listed in Table 2) is composed of 5 disks with expanded z coverage
from 25 cm to 135 cm. Due to the constraints of the EIC beam pipe dimensions
and the allocated detector space limitation, the inner radii of the ePIC MAPS
hadron endcap disks vary from 3.67 cm to 7.01 cm and the outer radii vary from
23 cm to 43 cm. Compared to the EIC project detector reference design, the
ePIC MAPS electron endcap detector (geometry parameters listed in Table 3)
includes one additional disk and the longitudinal and radial parameters of all
five disks have been updated to make full usage of the allocated space. The
pixel size of the current ePIC MAPS vertex and tracking detector design is
assumed at 10 $\mu\rm{m}$, and is open fors further tuning depends on the
progress of the sensor design. The material budgets of the ePIC MAPS vertex
detector design can reach as low as 0.05$\%~{}X/X_{0}$ per layer and the
material budgets are expected to not exceed 0.55$\%~{}X/X_{0}$ per layer for
the outer most MAPS batter layer.

Table 1: ePIC MAPS barrel vertex and tracking detector geometry Index |  $r$ (cm) |  $z_{\rm min}$ (cm) |  $z_{\rm max}$ (cm) | Pixel Pitch ($\mu\rm{m}$) | Material Budget ($X/X_{0}$)
---|---|---|---|---|---
1 | 3.6 | -13.5 | 13.5 | 10 | 0.05$\%$
2 | 4.8 | -13.5 | 13.5 | 10 | 0.05$\%$
3 | 12 | -13.5 | 13.5 | 10 | 0.05$\%$
4 | 27 | -27 | 27 | 10 | 0.25$\%$
5 | 42 | -42 | 42 | 10 | 0.55$\%$
Table 2: ePIC MAPS hadron endcap tracking detector geometry Index |  $z$ (cm) |  $r_{\rm in}$ (cm) |  $r_{\rm out}$ (cm) | Pixel Pitch ($\mu\rm{m}$) | Material Budget ($X/X_{0}$)
---|---|---|---|---|---
1 | 25 | 3.676 | 23 | 10 | 0.24$\%$
2 | 45 | 3.676 | 43 | 10 | 0.24$\%$
3 | 70 | 3.842 | 43 | 10 | 0.24$\%$
4 | 100 | 5.443 | 43 | 10 | 0.24$\%$
5 | 135 | 7.014 | 43 | 10 | 0.24$\%$
Table 3: ePIC MAPS electron endcap tracking detector geometry Index |  $z$ (cm) |  $r_{\rm in}$ (cm) |  $r_{\rm out}$ (cm) | Pixel Pitch ($\mu\rm{m}$) | Material Budget ($X/X_{0}$)
---|---|---|---|---|---
1 | -25 | 3.676 | 23 | 10 | 0.24$\%$
2 | -45 | 3.676 | 43 | 10 | 0.24$\%$
3 | -65 | 3.676 | 43 | 10 | 0.24$\%$
4 | -90 | 4.006 | 43 | 10 | 0.24$\%$
5 | -115 | 4.635 | 43 | 10 | 0.24$\%$
Table 4: ePIC AC-LGAD outer tracker detector geometry Region |  $z$ locations (cm) |  $r$ locations (cm) | pixel size (mm) | Material Budget ($X/X_{0}$)
---|---|---|---|---
Barrel | -120 to 120 | 64.6 | $0.5\times 1.0$ |  $\sim$ 1$\%$
h-endcap | 192 | 8.5 to 67 | $0.5\times 0.5$ |  $\sim$ 7$\%$

The current ePIC tracking detector includes the AC-LGAD outer tracker
subsystem in the barrel and hadron endcap region as well. The geometry
parameters of the AC-LGAD outer tracker are listed in Table 4. The pixel size
of the barrel AC-LGAD tracker is 0.5 mm by 1.0 mm and the hadron endcap AC-
LGAD tracker pixel size is 0.5 mm by 0.5 mm. The barrel AC-LGAD tracker
consists of around 2.4 million pixels to cover around 10.9 $\rm{m}^{2}$ active
area. The hadron endcap AC-LGAD tracker is composed of around 8.8 million
channels with around 2.22 $\rm{m}^{2}$ active area. The barrel AC-LGAD tracker
is expected to achieve around 30 $\mu\rm{m}$ spatial resolution in the
$r\varphi$ plane and around 30 ps timing resolution. The hadron endcap AC-LGAD
tracker is designed to obtain around 30 $\mu\rm{m}$ spatial resolution in the
$xy$ plane and around 25 ps timing resolution.

Figure 2: Tracking performance of the current ePIC detector design in the
pseudorapidity regions of $-3.5<\eta<-2.5$, $-1.0<\eta<1.0$ and
$2.5<\eta<3.5$. The track momentum dependent momentum resolution (top row) and
transverse momentum dependent transverse Distance of Closest Approach
($\rm{DCA}_{\rm{2D}}$) resolution (bottom row) are evaluated with the new 1.7
T ePIC magnet. The EIC yellow report tracking requirements in the respective
pseudorapidity regions are highlighted in blue dashed lines.

To evaluate the ePIC detector tracking performance, the active detector
volume, the related mechanical structure and the estimated service parts of
the current ePIC detector design have been implemented in GEANT4 [6]
simulation. Figure 2 shows the tracking momentum and spatial resolutions of
the current ePIC detector design in GEANT4 simulation. With the optimized ePIC
detector design, the track momentum dependent momentum resolutions in the
pseudorapidity regions of $-1<\eta<1$ and $2.5<\eta<3.5$ meet the EIC yellow
report detector requirements as shown in the top row of Fig. 2. Although the
tracking momentum resolution of the current ePIC detector design in the
pseudorapidity region of $-3.5<\eta<-2.5$ is worse than the EIC yellow report
detector requirements, further studies to optimize the integrated detector
subsystems and reduce the detector material budgets in the electron endcap
region are underway. The transverse momentum dependent transverse Distance of
Closest Approach ($\rm{DCA}_{\rm{2D}}$) resolution in most pseudorapidity
regions meets the EIC yellow report requirements. As the EIC backgrounds such
as the beam gas and synchrotron radiation backgrounds, are still under study
in simulation, additional tuning of the ePIC tracking detector geometry is
expected to provide sufficient number of hits for tracking pattern recognition
studies. The ePIC detector tracking performance will be updated accordingly
with new simulation samples which will include the EIC backgrounds.

##  3 EIC Silicon Detector Technology and R$\&$D Status

The EIC silicon R$\&$D supported by the EIC project detector program primarily
focuses on the EIC MAPS and AC-LGAD sensor developments, vertex and tracking
detector mechanical design and readout options. The 65 nm MAPS sensor R$\&$D
for the EIC silicon vertex detector aligns well with the ALICE ITS3
developments [3]. Better than 4 $\mu\rm{m}$ spatial resolution has been
achieved by the ALICE ITS3 prototype sensor beam tests with 15 $\mu\rm{m}$
pixel size and $10^{15}~{}\rm{n}_{eq}\rm{cm}^{-2}$ irradiation dose [7, 8].
Additional developments towards 0.05$\%~{}X/X_{0}$ material budgets per layer,
on the order of 100 ns timing resolution and low fake-hit rate are ongoing.
The EIC flat MAPS prototype sensors for the sagitta layers and endap disks are
still under design to obtain similar features of the ALICE ITS3 technology.
Figure 3 shows the mechanical design of the ePIC MAPS vertex and tracking
detector current design (see geometry parameters in Table 1, Table 2 and Table
3) with the evaluated support structure, power and data cables and cooling
system. Other components such as the readout modules are under developments.

Figure 3: The mechanical design of the ePIC MAPS vertex and tracking detector
according to the current optimized ePIC detector geometry (details in Table 1,
Table 2 and Table 3). As highlighted inside the red box, the current MAPS disk
layout is based on stitched flat MAPS sensors.

The ePIC AC-LGAD outer tracker will use two types of sensors: 0.5 mm by 0.5 mm
pixel sensor and 0.5 mm by 1.0 mm strip sensor. As shown in the left panel of
Figure 4, new AC-LGAD pixel and strip prototype sensors have been produced at
BNL and Hamamatsu Photonics K.K. (HPK) company. The AC-LGAD strip sensors with
different strip lengths have been characterized with the 120 GeV proton beam
at FNAL (see the setup configuration in the right panel of Figure 4). Around
30 $\mu\rm{m}$ spatial resolution and better than 30 ps timing resolution per
hit have been achieved from the FNAL AC-LGAD strip sensor beam tests, which
meets the EIC AC-LGAD detector design requirement.

Figure 4: EIC AC-LGAD prototype sensor production at BNL and beam test setup
at FNAL. The left panel shows the EIC AC-LGAD pixel and strip sensors have
been diced on a silicon wafer at BNL. The configuration of AC-LGAD prototype
sensor beam tests at the FNAL test beam facility is shown in the right panel.
Figure 5: The carrier board of the AC-LGAD prototype pixel sensor and the
3-layer AC-LGAD sensor telescope ${}^{90}\rm{Sr}$ source test configuration is
shown in the left panel. The MALTA prototype sensor carrier board and the
bench test configuration of a single MALTA sensor is shown in the right panel.

Other advanced silicon prototype sensors have been studied in parallel for the
EIC detectors. Through supports by the Los Alamos National Laboratory (LANL)
Laboratory Directed Research $\&$ Development (LDRD), the LANL EIC team has
setup the test bench and performed characterizations of LGAD, AC-LGAD and
Depleted MAPS (DMAPS) prototype sensors in collaboration with other
institutions. Figure 5 shows the bench setup of the 3-layer AC-LGAD pixel
sensor telescope ${}^{90}\rm{Sr}$ source tests and the DMAPS MALTA [9] sensor
bench tests. The pixel size of the AC-LGAD prototype pixel sensor is 500
$\mu\rm{m}$ and the MALTA prototype sensor pixel size is 36.4 $\mu\rm{m}$. The
AC-LGAD prototype sensor contains a 4 by 4 pixel array and the MALTA prototype
sensor is composed of a 512 by 512 pixel matrix. The AC-LGAD prototype sensor
bench test system consists of a scintillator trigger unit to trigger on
${}^{90}\rm{Sr}$ beta decay electron events, low and high power supplies to
operate the sensors and their carrier boards, a 2.5 GHz oscilloscope to
evaluate the analog signals from the sensor, daisy-chained 16-channel CAEN
1730s digitizers to perform fast digitization, and a Data Acquisition (DAQ)
computer to collect the digitized data. For the AC-LGAD sensor data
digitization, one 16-channel CAEN 1730s module is connected to one AC-LGAD
sensor carrier board, which is wire-bounded to 16 AC-LGAD pixels. The MALTA
prototype sensor test system is composed of a XILINX KC705 module for data
processing, a scintillator trigger unit to provide the ${}^{90}\rm{Sr}$ beta
decay electron event trigger, a XILINX KC705 module for parallel trigger
synchronization and transmission, and a DAQ computer to provide the MALTA
sensor operation, online monitoring and data collection.

Figure 6: The AC-LGAD prototype pixel sensor ${}^{90}\rm{Sr}$ source test
results. The accumulated digitized signals of a AC-LGAD pixel (left) and the
average pulse shape of this AC-LGAD pixel (middle) in the ${}^{90}\rm{Sr}$
source tests. The event display of the reconstructed electron track from
clustered hits collected from the three-layer AC-LGAD telescope
${}^{90}\rm{Sr}$ tests.

In the AC-LGAD sensor ${}^{90}\rm{Sr}$ source tests, the amplitude of the AC-
LGAD pixel analog signal is around -500 mV and the pulse Full Width Half
Maximum (FWHM) is around 1 ns. The left panel of Figure 6 shows the
accumulated digitized signals of a single AC-LGAD pixel collected from the
${}^{90}\rm{Sr}$ source tests and the average pulse shape is shown in the
middle panel of Figure 6. The gain uniformity of all 16 pixels in the AC-LGAD
sensor has been verified. A three-layer AC-LGAD prototype sensor telescope has
been assembled and tested with the ${}^{90}\rm{Sr}$ beta decay electrons. For
the AC-LGAD telescope tests, synchronized triggers are distributed to the
connected CAEN 1730s digitizers and a offline analysis framework has been
developed to reconstruct electron tracks from clustered hits located in
different AC-LGAD sensors. The right panel of Figure 6 illustrates one event
display of a reconstructed electron track (highlighted in red) in the AC-LGAD
telescope ${}^{90}\rm{Sr}$ source tests. The radiation tolerance of the AC-
LGAD prototype sensors has been studied in 2022 with the 500 MeV proton beam
at the LANL LANSCE facility and the accumulated radiation doses are from
$10^{13}-10^{16}~{}\rm{n}_{eq}\rm{cm}^{-2}$. Ongoing data analysis will
characterize and compare the AC-LGAD sensor performance at different radiation
doses. The MALTA sensor bench tests include threshold and noise scan and hit
occupancy studies with and without a ${}^{90}\rm{Sr}$ source [10]. A new
version of the MALTA sensor: MALTA2 [11] has been recently produced and
characterized with beam tests at CERN. Better than 2 ns timing resolution has
been obtained by the MALTA2 sensor. The threshold over noise ratio of the
MALTA2 sensor is better than 10 even with $3\times
10^{15}~{}\rm{n}_{eq}\rm{cm}^{-2}$ radiation dose.

##  4 Summary and Outlook

The EIC project detector design is evolving rapidly led by the newly formed
ePIC collaboration. The MAPS and AC-LGAD based silicon detector subsystems
play an essential role in providing precise vertex determination, track
reconstruction and particle identification in a broad kinematic coverage for
the proposed EIC physics measurements. Good progresses have been achieved for
the ePIC MAPS and AC-LGAD detector design optimization and related R$\&$D,
which are in good synergy with other major detector development projects,
which have similar timeline and schedules, such as the High Luminosity LHC
detector upgrade. The current ePIC detector design will deliver an integrated
silicon detector that can obtain better than 10 $\mu\rm{m}$ tracking spatial
resolution, sub-percent tracking momentum resolution and better than 30 ps
timing resolution. Ongoing EIC silicon detector R$\&$D will provide new sensor
concepts, updated mechanical design and options for the readout architecture
for the EIC detector construction, which is scheduled to start in 2025.

## References

  * [1] R. Abdul Khalek et al., _Science Requirements and Detector Concepts for the Electron-Ion Collider: EIC Yellow Report_ , arXiv: 2103.05419 (2021).
  * [2] J. K. Adkins et al., _Design of the ECCE Detector for the Electron Ion Collider_ , arXiv: 2209.02590 (2022).
  * [3] ALICE collaboration, _Letter of intent for an ALICE ITS upgrade in LS3_ , CERN-LHCC-2019-018, LHCC-I-034 (2019).
  * [4] F. Sauli and A. Sharma, _Micropattern Gaseous Detectors_ , Annu. Rev. Nucl. 49 (1999) 341-388.
  * [5] G. Giacomini et al., _Fabrication and performance of AC-coupled LGADs_ , INST 14 (2019) P09004.
  * [6] S. Agostinelli, et al., _Geant4 - a simulation toolkit_ , Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 506 (3) (2003) 250-303.
  * [7] M. Buckland, ALICE collaboration, _Development of the ITS3: A bent-silicon vertex detector for ALICE in the LHC Run 4_ , Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 1039 (2022) 166875.
  * [8] A. Kluge, ALICE collaboration, _ALICE - ITS3 - A bent, wafer-scale CMOS detector_ , Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment 1041 (2022) 167315.
  * [9] H. Pernegger et al., _First tests of a novel radiation hard CMOS sensor process for Depleted Monolithic Active Pixel Sensors_ , JINST 12 (2020) P06008.
  * [10] X. Li et al., _Forward silicon vertex/tracking detector design and R $\&$D for the future Electron-Ion Collider_, PoS PANIC2021 (2022) 084.
  * [11] M. van Rijnbach et al., _Radiation hardness and timing performance in MALTA monolithic pixel sensors in TowerJazz 180 nm_ , JINST 17 (2022) C04034.
  * [12] G. Gustavino et al., _Timing performance of radiation hard MALTA monolithic Pixel sensors_ , arXiv: 2209.14676 (2022).


problem settings, while also providing data scientists with a clear
understanding of the practical significance and open challenges in various
application domains. Moreover, the paper highlights the current status of
research in this area, shedding light on future trends and potential avenues
for interdisciplinary collaboration. As the field of LLM domain specialization
continues to evolve, this survey serves as a valuable resource for researchers
and practitioners, fostering further advancements and innovations in the
application of artificial intelligence across diverse domains.