:
LEAVE UNSET 2023 LEAVE UNSET LEAVE UNSET Machine Learning for Health (ML4H)
2023

# Adaptive Interventions with User-Defined Goals for Health Behavior Change

Aishwarya Mandyam<EMAIL_ADDRESS>
Stanford University These authors contributed equally Matthew
Jörke11footnotemark: 1<EMAIL_ADDRESS>
Stanford University  Barbara E. Engelhardt<EMAIL_ADDRESS>
Gladstone Institutes  Stanford University  Emma Brunskill
<EMAIL_ADDRESS>
Stanford University

###### Abstract

Physical inactivity remains a major public health concern, having associations
with adverse health outcomes such as cardiovascular disease and type-2
diabetes. Mobile health applications present a promising avenue for low-cost,
scalable physical activity promotion, yet often suffer from small effect sizes
and low adherence rates, particularly in comparison to human coaching. Goal-
setting is a critical component of health coaching that has been underutilized
in adaptive algorithms for mobile health interventions. This paper introduces
a modification to the Thompson sampling algorithm that places emphasis on
individualized goal-setting by optimizing personalized reward functions. As a
step towards supporting goal-setting, this paper offers a balanced approach
that can leverage shared structure while optimizing individual preferences and
goals. We prove that our modification incurs only a constant penalty on the
cumulative regret while preserving the sample complexity benefits of data
sharing. In a physical activity simulator, we demonstrate that our algorithm
achieves substantial improvements in cumulative regret compared to baselines
that do not share data or do not optimize for individualized rewards.

###### keywords:

Mobile health, contextual bandits, health behavior change

##  1 Introduction

Physical inactivity is a leading risk factor for global mortality, associated
with adverse health outcomes such as cardiovascular disease, cancer, and
type-2 diabetes. Over a quarter of the world’s population is physically
inactive, failing to meet the WHO’s global recommendations for weekly physical
activity (World Health Organization, 2022).

Researchers are increasingly exploring mobile health (mHealth) applications as
a low-cost, scalable, and accessible approach to motivate health behavior
change (Domin et al., 2021; Hicks et al., 2023). Within machine learning and
statistics, there is a burgeoning interest in personalizing health behavior
change interventions by applying adaptive experimentation or reinforcement
learning algorithms to automatically discover which interventions work best
for different individuals across diverse contexts (Klasnja et al., 2019; Baek
et al., 2023; Mintz et al., 2019; Ruggeri et al., 2023). Despite the potential
for low-cost, personalized, contextually tailored interventions to promote
positive health outcomes, mHealth interventions are known to suffer from small
effect sizes (Yang and Van Stee, 2019) and low adherence (Yang et al., 2020),
particularly when compared to human health coaches (McEwan et al., 2016).

A key component of effective health coaching is goal-setting (Olsen and
Nesbitt, 2010; Wolever et al., 2013; Epton et al., 2017), which is one of the
most common strategies for encouraging physical activity behavior change
(Howlett et al., 2019). Goal-setting theory highlights that effective goals
can focus attention towards goal-related activities and lead to greater effort
and persistence (Locke and Latham, 2002). While goal-setting strategies have
been explored in mHealth and personal informatics tools (Ekhtiar et al.,
2023), it is largely absent from algorithms for adaptive, personalized
interventions. Instead, prior work typically optimizes for some shared,
measureable outcome (e.g., step count), implicitly assuming that each user
wants to maximize this quantity as much as possible. Not only does this
approach neglect the positive psychological effects of goal setting on long-
term motivation, it also ignores the diversity of people’s physical activity
goals (Epton et al., 2017). Moreover, when individuals receive feedback on
misaligned or overly ambitious goals, it can lead to abandonment or
habituation (Locke and Latham, 2002; Peng et al., 2021).

One approach to optimizing over individualized goals is to train an
independent policy for each person in a cohort, using their own goal as a
reward function. However, this approach does not share data across individuals
and thus requires more samples to learn an optimal intervention assignment
rule. Prior work on data sharing across users has been shown to be useful, but
existing work does not allow users to have distinct preferences (Zhou and
Brunskill, 2016; Tomkins et al., 2020). Meanwhile, existing work on optimizing
individual preferences focuses on efficiently learning complex, unknown
preferences for a single user (Roijers et al., 2021) or adapting to an
individual’s preference changes over time (Hariri et al., 2015). In contrast,
we assume that preferences are known for each user and the objective is to
quickly maximize individualized rewards across all users.

In summary, the contributions of this work to the existing mHealth literature
are as follows:

  1. 1.

We propose a modification to the Thompson sampling algorithm for linear
contextual bandits (Agrawal and Goyal, 2013) that can optimize for
individualized reward functions (i.e., goals) while maintaining the ability to
share information across individuals. We formalize user goals as a Lipschitz
continuous function of a shared outcome variable (e.g., step count) that has
common structure across individuals.

  2. 2.

We prove that our modification incurs only a Lipschtiz constant penalty on the
cumulative regret while preserving the sample complexity benefits of data
sharing.

  3. 3.

We apply this approach in a modified version of a previously proposed
simulator for physical activity (Yao et al., 2021) and demonstrate that it
outperforms policies trained separately for each user and policies that
optimize for a misspecified reward that ignores personal goals.

##  2 Methods

Multi-Objective Multi-User Thompson Sampling Input: number of participants
$P$, number of timesteps $T$, prior
$p(\boldsymbol{\theta})\sim\mathcal{N}(0,\lambda^{-1}I)$, user-specific
utilities $U_{i,t,k}$ and weights $w_{i,t,k}$, $\forall i\in[P],k\in[K]$

For all $m\in[M]$, let $p_{m}=p(\boldsymbol{\theta})$, $\mu_{m}=\mathbf{0}$,
$b_{m}=\mathbf{0}$, and $\Sigma=\lambda I$ for  _$t=1,\dots,T$_ do

for  _$i=1,\dots,P$_ do

Sample parameters from posterior $\tilde{\boldsymbol{\theta}}_{m}\sim p_{m}$
and compute outcomes:
$[\tilde{\mathbf{y}}(s_{i,t},a)]_{m}=\boldsymbol{\phi}(s_{i,t},a)^{\top}\tilde{\boldsymbol{\theta}}_{m}$
Choose action to maximize weighted utility using known $U_{i,t,k}$ and
$w_{i,t,k}$:
$a_{i,t}=\mathop{\mathrm{arg\,max}}_{a}\sum_{k=1}^{K}w_{i,t,k}U_{i,t,k}\left([\tilde{\mathbf{y}}(s_{i,t},a)]_{m_{k}},\boldsymbol{\phi}(s_{i,t},a)\right)$
Observe the true outcome $\mathbf{y}_{i,t}$ and update the posterior:
$\Sigma\leftarrow\Sigma+\boldsymbol{\phi}_{i,t}\boldsymbol{\phi}_{i,t}^{\top},\quad\forall
m\in[M]:b_{m}\leftarrow
b_{m}+\boldsymbol{\phi}_{i,t}[\mathbf{y}_{i,t}]_{m},\quad\mu_{m}\leftarrow\Sigma^{-1}b_{m},\quad
p_{m}\leftarrow\mathcal{N}(\mu_{m},\Sigma^{-1})$

end for

end for

We present an augmentation to the popular Thompson sampling algorithm for
linear contextual bandits (Agrawal and Goyal, 2013) that enables optimization
of user-specific goals while sharing data across users. At each timestep, the
algorithm samples potential outcomes for each action from the posterior.
Rather than directly optimizing for a shared outcome, our augmentation
propagates outcomes through a known, user-specific reward function, modeled as
a weighted sum of utilities. We formalize our setting as a Bayesian contextual
bandit environment in which a policy selects actions for $i\in[P]$
participants over $t\in[T]$ rounds. At each time $t$, each individual $i$’s
context is sampled i.i.d. from a known distribution $s_{i,t}\sim\rho$. A
policy is tasked with learning to select an action $a_{i,t}\in[A]$ as a
function of $s_{i,t}$ to maximize reward $r_{i,t}(s_{i,t},a_{i,t})$, described
below (Eq. 2). Let
$\boldsymbol{\phi}_{i,t}\triangleq\boldsymbol{\phi}(s_{i,t},a_{i,t})\in\mathbb{R}^{d}$
be a known state-action featurization. For example, $\boldsymbol{\phi}_{i,t}$
may contain demographic information about user $i$ (e.g., age, gender, pre-
existing medical conditions), action-specific information (e.g., indicator
variables for a given treatment or treatment interaction terms), and other
contextual information (e.g., time, day of the week, location, number of
previous interventions). Let
$\mathbf{y}_{i,t}\triangleq\mathbf{y}(s_{i,t},a_{i,t})\in\mathbb{R}^{M}$ be a
vector containing $M$ potentially unknown, but measurable quantities about a
user. For example, $\mathbf{y}_{i,t}$ can represent step count, sleep
duration, heart rate, or blood pressure. We assume each entry $m$ in
$\mathbf{y}_{i,t}$ follows a linear model,111While we assume linearity for
Theorem 2.1, Algorithm 2 is amenable to non-linear extensions so long as
$\mathbf{y}$ admits an efficiently computable posterior given the history.
$[\mathbf{y}_{i,t}]_{m}=\boldsymbol{\phi}_{i,t}^{\top}\boldsymbol{\theta}_{m}+\varepsilon$
(1) where $\varepsilon$ is independent Gaussian noise. Note that
$\boldsymbol{\theta}_{m}$ is specific to the outcome and is shared across
users. In our Bayesian setting, we assume that each $\boldsymbol{\theta}_{m}$
is drawn from a known prior $p(\boldsymbol{\theta})$. Instead of optimizing
directly for the outcome $\mathbf{y}_{i,t}$ as is standard in prior work, our
algorithm selects actions to maximize a weighted sum of $K$ user-specific
utilities. Each utility $U_{i,t,k}$ can encode a user’s goal (e.g., walking
10,000 steps per day) or preferences over interventions (e.g., penalizing
excessive notifications), which may be context-specific (e.g., ignoring
interventions on weekends). Each utility is associated with a known weight
$w_{i,t,k}$ satisfying $\sum_{k=1}^{K}w_{i,t,k}=1$. These weights are either
fixed by the algorithm designer or can be learned in advance, using preference
elicitation. Specifically, upon taking action $a_{i,t}\in[A]$ in state
$s_{i,t}$, the environment reveals reward  $\displaystyle
r_{i,t}(s_{i,t},a_{i,t})$ $\displaystyle=\sum_{k=1}^{K}w_{i,t,k}\cdot
U_{i,t,k}([\mathbf{y}_{i,t}]_{m_{k}},\boldsymbol{\phi}_{i,t}).$ (2) We assume
that each utility function $U_{i,t,k}$ depends only on the $m_{k}$-th entry of
$\mathbf{y}_{i,t}$ and is $L$-Lipschitz in $[\mathbf{y}_{i,t}]_{m_{k}}$, i.e.,
$\forall y,|U(y,\phi)-U(y^{\prime},\phi)|\leq L|y-y^{\prime}|$. While this
reward function is user-specific and can vary across time, it is a known and
deterministic function of $\mathbf{y}_{i,t}$, which is modeled as a stationary
function of $\boldsymbol{\phi}_{i,t}$ that is shared across users. In
practice, one might offer a range of potential utilities to a user and they
could select the weights over the utilities they care most about. Users could
also update utilities or weights over time if their priorities change. Since
the parameters corresponding to all utility functions and the weights for each
user are fully known, the objective is to efficiently learn
$\boldsymbol{\theta}_{m}$, which is unknown. Our procedure is described in
Algorithm 2. We analyze performance with respect to Bayesian cumulative regret
(BR) as a function of $N=P\cdot T$, the total number of samples observed by
the algorithm.

$\displaystyle\text{BR}(N)=\mathbb{E}\left[\sum_{t=1}^{T}\sum_{i=1}^{P}\max_{a\in[A]}r_{i,t}(s_{i,t},a)-r_{i,t}(s_{i,t},a_{i,t})\right]$

(3) The expectation is taken both over the prior $p(\boldsymbol{\theta})$ and
the policy trajectories. In Theorem 2.1 below, we demonstrate that by
propagating uncertainty through these known utilities, regret is only scaled
by a Lipschitz constant with a logarithmic penalty on the number of outcomes
$M$.

######  Theorem 2.1.

After running Algorithm 2 for $N=P\cdot T$ samples, where $P$ is the number of
participants in the cohort and $T$ is the time horizon, we achieve Bayesian
cumulative regret on the order of

$\text{BR}(N)\leq O(Ld\sqrt{N\log(NM)\log(N/d)}).$ (4)

The proof is provided in Appendix A and follows from a modification of
Lattimore and Szepesvári (2020) Theorem 36.4, which establishes an
$O(d\sqrt{N\log(N)\log(N/d)}$ regret bound for standard Thompson sampling with
linear contextual bandits.

##  3 Experiments

Figure 1: We compare step count and reward for two policies: one policy (blue,
dashed) always sends notifications, and another policy (red, solid) sends
notifications only when a user’s step count drops below their desired goal
(green dashed line). A policy that always sends a notification achieves higher
step count, but lower reward due to notification burden. Shaded area is $\pm
1\sigma$ across $P=20$ participants.

Figure 2: We compare our approach to traditional Thompson sampling optimizing
for step count, a random policy, and separate policies for each user. Our
approach generates less cumulative regret because it is selective about
sending notifications and is optimizing for the correct element, the reward.
Shaded area is $\pm 1\sigma$ across 100 trials.

We now demonstrate the results of our algorithm using a simulator inspired by
mobile health intervention studies such as HeartSteps (Klasnja et al., 2018).
The HeartSteps study aims to improve users’ physical activity outcomes by
adaptively sending notifications encouraging them to walk. We draw inspiration
from prior simulators designed to model the HeartSteps study (Yao et al.,
2021; Liao et al., 2015), with some modifications specific to our setting.

First, we model treatment effect heterogeneity and non-stationarity by
positing two cluster of users: active and inactive users. Empirically, studies
have found that some users consistently respond positively to notifications,
whereas others’ response declines over time (Liao et al., 2015). The MyHeart
Counts study (Shcherbina et al., 2019), which contains a similar setup to
HeartSteps without adaptive interventions, also describes clusters of users
with varying baseline activity levels and treatment responses. Next, we adopt
an autoregressive structure to model these individual variations. This allows
for an intervention at timestep $t$ to impact future timesteps $t^{\prime}>t$
while preserving a bandit model, allowing us model habit formation (Hagger,
2019; Peng et al., 2021). Lastly, we keep track of the number of notifications
that a user has received over the course of the study in order to model
notification burden. Wang et al. (2021) notes that too many notifications
negatively impact on a user’s engagement with the study, and we want to
explicitly model the tradeoff between sending notifications and increasing
step count.

To construct our simulator (described in full detail in Appendix B), we model
one outcome $y_{i,t}$ (step count) as a function of $\boldsymbol{\phi}_{i,t}$,
which contains information such as the previous day’s step count $y_{i,t-1}$
and group specific treatment effects. We define two utility functions
$U_{i,t,0}$, a piecewise linear function that models step count goals, and
$U_{i,t,1}$, which penalizes higher numbers of notifications. Inspired by
habituation and recovery dynamics evidenced in prior work, we choose a
quadratic function to represent notification burden (Mintz et al., 2019;
Bertsimas et al., 2022). We model group differences by assuming that inactive
users have a lower baseline activity level and treatment effect than active
users, with individual variation resulting from the auto-regressive structure
in each population. We assume that active users place a higher weight on their
step count goal and a lower weight on the notification goal, and vice versa
for inactive users.

Now, we evaluate our proposed algorithm in the step count simulator. We first
note that optimizing for the outcome, step count, is not equivalent to
optimizing for reward, which is a combination of utility functions. In Figure
2, we see that it is possible to increase step count by always sending
notifications, but that this eventually reduces reward due to our notification
penalty, meaning that a user would not achieve their personal goal of
achieving their step count with minimal notification burden.

Next, we present results comparing four choices for algorithms: a random
policy, Thompson sampling optimizing for $y_{i,t}$, Thompson sampling
optimizing for $r_{i,t}$ (Algorithm 2), and one which learns indepdent
policies for each user without sharing data. We find that Thompson sampling
optimzing for $r_{i,t}$ generates the lowest cumulative regret across our
horizon ($T=100$, $P=20$) when cumulative regret is measured with respect to
reward (Eq. 3), because it can correctly balance the number of notifications
sent while optimizing for a user’s goals. Additionally, we find that we get a
speedup in cumulative regret when we share data, indicating that sharing data
allows us to learn the optimal policy for each user more quickly.

##  4 Discussion

In this work, we introduce a modification to Thompson sampling that allows us
to effectively optimize for user-specific goals that only introduces a
constant Lipschitz penalty on cumulative regret. In simulation, we demonstrate
that our approach achieves a cumulative regret speedup in comparison to
policies that are trained separately for each user and policies that do not
account for user-specific utilities.

There are a number of exciting directions for future work. We are actively
developing semi-synthetic simulators based on empirical study data to evaluate
our algorithm under more realistic conditions. This may require us to relax
our linearity assumptions, e.g., using Gaussian processes (Chowdhury and
Gopalan, 2017) or semi-parametric models (Greenewald et al., 2017;
Krishnamurthy et al., 2018), and explicitly account for treatment effect
heterogeneity, e.g., using hierarchical (Hong et al., 2022) or mixed-effects
models (Tomkins et al., 2020). Moreover, our bandit model accounts for non-
stationarity and habituation by including time-dependent treatment effect
terms in the context and penalizing excessive notifications in the reward.
While this bandit model has sample complexity benefits, it may be more
appropriate to model long-term treatment effect dynamics using an MDP model,
which could also be learned using Thompson sampling (Osband et al., 2013).

Towards expanding our theoretical results, our algorithm and regret bounds
assume that each action is assigned sequentially (i.e., the model is updated
between assigning treatment to user $i$ and $i+1$), though in practice
treatment would likely be assigned in parallel. Further, it may be possible to
achieve similar regret guarantees for non-Lipschitz utility functions, e.g.,
indicator functions with a fininte number of discontinuities. It may also
prove useful to derive simple regret bounds to investigate whether the
presence of user-specific goals can change exploration strategies in a pure
exploration setting.

We thank Jonathan Lee, Jiayu Yao, and Scott Fleming for feedback. AM was
supported in part by a Stanford Engineering Fellowship. This work was
supported in part by NSF grant #2112926.

## References

  * Agrawal and Goyal (2013) Shipra Agrawal and Navin Goyal.  Thompson sampling for contextual bandits with linear payoffs.  In _International conference on machine learning_ , pages 127–135. PMLR, 2013.
  * Baek et al. (2023) Jackie Baek, Justin J. Boutilier, Vivek F. Farias, Jonas Oddur Jonasson, and Erez Yoeli.  Policy optimization for personalized interventions in behavioral health, 2023.
  * Bertsimas et al. (2022) Dimitris Bertsimas, Predrag Klasnja, Susan Murphy, and Liangyuan Na.  Data-driven interpretable policy construction for personalized mobile health.  In _2022 IEEE International Conference on Digital Health (ICDH)_ , pages 13–22, 2022.  10.1109/ICDH55609.2022.00010.
  * Chowdhury and Gopalan (2017) Sayak Ray Chowdhury and Aditya Gopalan.  On kernelized multi-armed bandits.  In _International Conference on Machine Learning_ , pages 844–853. PMLR, 2017.
  * Domin et al. (2021) Alex Domin, Donna Spruijt-Metz, Daniel Theisen, Yacine Ouzzahra, Claus Vögele, et al.  Smartphone-based interventions for physical activity promotion: scoping review of the evidence over the last 10 years.  _JMIR mHealth and uHealth_ , 9(7):e24308, 2021\.
  * Ekhtiar et al. (2023) Tina Ekhtiar, Armağan Karahanoğlu, Rúben Gouveia, and Geke Ludden.  Goals for goal setting: A scoping review on personal informatics.  2023\.
  * Epton et al. (2017) Tracy Epton, Sinead Currie, and Christopher J Armitage.  Unique effects of setting goals on behavior change: Systematic review and meta-analysis.  _Journal of consulting and clinical psychology_ , 85(12):1182, 2017.
  * Greenewald et al. (2017) Kristjan Greenewald, Ambuj Tewari, Susan Murphy, and Predag Klasnja.  Action centered contextual bandits.  _Advances in neural information processing systems_ , 30, 2017.
  * Hagger (2019) Martin S Hagger.  Habit and physical activity: Theoretical advances, practical implications, and agenda for future research.  _Psychology of Sport and Exercise_ , 42:118–129, 2019.
  * Hariri et al. (2015) Negar Hariri, Bamshad Mobasher, and R. Burke.  Adapting to user preference changes in interactive recommendation.  In _International Joint Conference on Artificial Intelligence_ , 2015\.
  * Hicks et al. (2023) Jennifer L Hicks, Melissa A Boswell, Tim Althoff, Alia J Crum, Joy P Ku, James A Landay, Paula ML Moya, Elizabeth L Murnane, Michael P Snyder, Abby C King, et al.  Leveraging mobile technology for public health promotion: A multidisciplinary perspective.  _Annual Review of Public Health_ , 44:131–150, 2023.
  * Hong et al. (2022) Joey Hong, Branislav Kveton, Manzil Zaheer, and Mohammad Ghavamzadeh.  Hierarchical bayesian bandits, 2022.
  * Howlett et al. (2019) Neil Howlett, Daksha Trivedi, Nicholas A Troop, and Angel Marie Chater.  Are physical activity interventions for healthy inactive adults effective in promoting behavior change and maintenance, and which behavior change techniques are effective? a systematic review and meta-analysis.  _Translational behavioral medicine_ , 9(1):147–157, 2019.
  * Klasnja et al. (2018) Predrag Klasnja, Shawna Smith, Nicholas Seewald, Andy Jinseok Lee, Kelly Hall, Brook Luers, Eric Hekler, and Susan Murphy.  Efficacy of contextually tailored suggestions for physical activity: A micro-randomized optimization trial of heartsteps.  _Annals of Behavioral Medicine_ , 53, 09 2018.  10.1093/abm/kay067.
  * Klasnja et al. (2019) Predrag Klasnja, Shawna Smith, Nicholas J Seewald, Andy Lee, Kelly Hall, Brook Luers, Eric B Hekler, and Susan A Murphy.  Efficacy of contextually tailored suggestions for physical activity: a micro-randomized optimization trial of heartsteps.  _Annals of Behavioral Medicine_ , 53(6):573–582, 2019.
  * Krishnamurthy et al. (2018) Akshay Krishnamurthy, Zhiwei Steven Wu, and Vasilis Syrgkanis.  Semiparametric contextual bandits.  In _International Conference on Machine Learning_ , pages 2776–2785. PMLR, 2018.
  * Lattimore and Szepesvári (2020) Tor Lattimore and Csaba Szepesvári.  _Bandit algorithms_.  Cambridge University Press, 2020.
  * Liao et al. (2015) Peng Liao, Predrag Klasnja, Ambuj Tewari, and Susan Murphy.  Sample size calculations for micro-randomized trials in mhealth.  _Statistics in medicine_ , 35, 12 2015.  10.1002/sim.6847.
  * Locke and Latham (2002) Edwin A Locke and Gary P Latham.  Building a practically useful theory of goal setting and task motivation: A 35-year odyssey.  _American psychologist_ , 57(9):705, 2002.
  * McEwan et al. (2016) Desmond McEwan, Samantha M Harden, Bruno D Zumbo, Benjamin D Sylvester, Megan Kaulius, Geralyn R Ruissen, A Justine Dowd, and Mark R Beauchamp.  The effectiveness of multi-component goal setting interventions for changing physical activity behaviour: a systematic review and meta-analysis.  _Health psychology review_ , 10(1):67–88, 2016\.
  * Mintz et al. (2019) Yonatan Mintz, Anil Aswani, Philip Kaminsky, Elena Flowers, and Yoshimi Fukuoka.  Non-stationary bandits with habituation and recovery dynamics, 2019.
  * Olsen and Nesbitt (2010) Jeanette M Olsen and Bonnie J Nesbitt.  Health coaching to improve healthy lifestyle behaviors: an integrative review.  _American journal of health promotion_ , 25(1):e1–e12, 2010.
  * Osband et al. (2013) Ian Osband, Daniel Russo, and Benjamin Van Roy.  (more) efficient reinforcement learning via posterior sampling.  _Advances in Neural Information Processing Systems_ , 26, 2013.
  * Peng et al. (2021) Wei Peng, Lin Li, Anastasia Kononova, Shelia Cotten, Kendra Kamp, Marie Bowen, et al.  Habit formation in wearable activity tracker use among older adults: qualitative study.  _JMIR mHealth and uHealth_ , 9(1):e22488, 2021\.
  * Roijers et al. (2021) Diederik Roijers, Luisa Zintgraf, Pieter Libin, Mathieu Reymond, Eugenio Bargiacchi, and Ann Nowe.  _Interactive Multi-objective Reinforcement Learning in Multi-armed Bandits with Gaussian Process Utility Models_ , pages 463–478.  02 2021.  ISBN 978-3-030-67663-6.  10.1007/978-3-030-67664-3_28.
  * Ruggeri et al. (2023) Kai Ruggeri, Amel Benzerga, Sanne Verra, and Tomas Folke.  A behavioral approach to personalizing public health.  _Behavioural Public Policy_ , 7(2):457–469, 2023\.  10.1017/bpp.2020.31.
  * Shcherbina et al. (2019) Anna Shcherbina, Steven G Hershman, Laura Lazzeroni, Abby C King, Jack W O’Sullivan, Eric Hekler, Yasbanoo Moayedi, Aleksandra Pavlovic, Daryl Waggott, Abhinav Sharma, Alan Yeung, Jeffrey W Christle, Matthew T Wheeler, Michael V McConnell, Robert A Harrington, and Euan A Ashley.  The effect of digital physical activity interventions on daily step count: a randomised controlled crossover substudy of the myheart counts cardiovascular health study.  _The Lancet Digital Health_ , 1(7):e344–e352, 2019.  ISSN 2589-7500.  https://doi.org/10.1016/S2589-7500(19)30129-3.
  * Tomkins et al. (2020) Sabina Tomkins, Peng Liao, Predrag Klasnja, and Susan Murphy.  Intelligentpooling: Practical thompson sampling for mhealth, 2020.
  * Wang et al. (2021) Shihan Wang, Chao Zhang, Ben Kröse, and Herke Hoof.  Optimizing adaptive notifications in mobile health interventions systems: Reinforcement learning from a data-driven behavioral simulator.  _Journal of Medical Systems_ , 45, 12 2021.  10.1007/s10916-021-01773-0.
  * Wolever et al. (2013) Ruth Q Wolever, Leigh Ann Simmons, Gary A Sforzo, Diana Dill, Miranda Kaye, Elizabeth M Bechard, Mary Elaine Southard, Mary Kennedy, Justine Vosloo, and Nancy Yang.  A systematic review of the literature on health and wellness coaching: defining a key behavioral intervention in healthcare.  _Global advances in health and medicine_ , 2(4):38–57, 2013.
  * World Health Organization (2022) World Health Organization.  Physical activity, 2022.  URL https://www.who.int/en/news-room/fact-sheets/detail/physical-activity.
  * Yang and Van Stee (2019) Qinghua Yang and Stephanie K Van Stee.  The comparative effectiveness of mobile phone interventions in improving health outcomes: meta-analytic review.  _JMIR mHealth and uHealth_ , 7(4):e11244, 2019\.
  * Yang et al. (2020) Xiaotian Yang, Lin Ma, Xi Zhao, and Atreyi Kankanhalli.  Factors influencing user’s adherence to physical activity applications: A scoping literature review and future directions.  _International Journal of Medical Informatics_ , 134:104039, 2020.
  * Yao et al. (2021) Jiayu Yao, Emma Brunskill, Weiwei Pan, Susan Murphy, and Finale Doshi-Velez.  Power constrained bandits.  In _Proceedings of the 6th Machine Learning for Healthcare Conference_ , pages 209–259, 2021.
  * Zhou and Brunskill (2016) Li Zhou and Emma Brunskill.  Latent contextual bandits and their application to personalized recommendations for new users.  In _Proceedings of the Twenty-Fifth International Joint Conference on Artificial Intelligence_ , pages 3646–3653, 2016.

##  Appendix A Bayesian Regret Proof

The proof structure largely follows Lattimore and Szepesvári (2020) Theorem
36.4, with modifications to account for our Lipschitz utility transformations.

###### Theorem 1

After running Algorithm 2 for $N=P\cdot T$ samples, where $P$ is the number of
participants in the cohort and $T$ is the time horizon, we achieves Bayesian
cumulative regret on the order of

$\text{BR}(N)\leq O(Ld\sqrt{N\log(NM)\log(N/d)})$

Proof: For simplicity, we define $n=Pt+i$ to convert the nested index $(i,t)$
to a single index. Using this notation, recall that

$\displaystyle r_{n}(s,a)$ $\displaystyle=\sum_{k=1}^{K}w_{n,k}\cdot
U_{n,k}\Big{(}[\mathbf{y}(s,a)]_{m_{k}},\boldsymbol{\phi}(s,a)\Big{)},\quad\text{where
}[\mathbf{y}(s,a)]_{m}=\boldsymbol{\theta}_{m}^{\top}\boldsymbol{\phi}(s,a)+\varepsilon,\quad\sum_{k=1}^{K}w_{n,k}=1$
(5)

Writing $y_{n,m}\triangleq[\mathbf{y}(s_{n},a_{n})]_{m}$ and
$\boldsymbol{\phi}_{n}\triangleq\boldsymbol{\phi}(s_{n},a_{n})$, define

$\displaystyle\hat{\boldsymbol{\theta}}_{n,m}$
$\displaystyle=\Sigma_{n}^{-1}\Phi_{n}^{\top}Y_{n,m},\quad\Phi_{n}=\begin{bmatrix}\rule[2.15277pt]{10.76385pt}{0.5pt}&\boldsymbol{\phi}_{1}^{\top}&\rule[2.15277pt]{10.76385pt}{0.5pt}\\\
&\vdots&\\\
\rule[2.15277pt]{10.76385pt}{0.5pt}&\boldsymbol{\phi}_{n}^{\top}&\rule[2.15277pt]{10.76385pt}{0.5pt}\\\
\end{bmatrix},\quad Y_{n,m}=\begin{bmatrix}y_{1,m}\\\ \vdots\\\
y_{n,m}\end{bmatrix},\quad\Sigma_{n}=\Phi_{n}^{\top}\Phi_{n}+\lambda I_{d}$
(6)

Assume that $\forall m,\|\boldsymbol{\theta}_{m}\|_{2}\leq S$,
$\sup_{s,a}\|\boldsymbol{\phi}(s,a)\|_{2}\leq Q$, $\forall
m,\sup_{s,a}|\boldsymbol{\theta}_{m}^{\top}\boldsymbol{\phi}(s,a)|\leq 1$. Let
each $U_{n,k}$ be $L$-Lipschitz in $y$:
$|U(y,\boldsymbol{\phi})-U(y^{\prime},\boldsymbol{\phi})|\leq
L|y-y^{\prime}|$. Define the upper confidence bound

$\displaystyle\text{UCB}_{n}(s,a)$ $\displaystyle=\sum_{k=1}^{K}w_{n,k}\cdot
U_{n,k}\Big{(}\hat{y}_{n,m_{k}}(s,a),\boldsymbol{\phi}(s,a)\Big{)}+\beta
L\|\boldsymbol{\phi}(s,a)\|_{\Sigma_{n-1}^{-1}}$ (7)
$\displaystyle\hat{y}_{n,m}(s,a)$
$\displaystyle=\boldsymbol{\phi}(s,a)^{\top}\hat{\boldsymbol{\theta}}_{n-1,m}$
(8) $\displaystyle\beta$
$\displaystyle=\sqrt{\lambda}S+\sqrt{2\log\left(NM\right)+d\log\left(1+\frac{NQ^{2}}{\lambda
d}\right)}$ (9)

By Lattimore and Szepesvári (2020) Theorem 20.5 and a union bound over
$m\in[M]$, it holds that

$\displaystyle P\left(\exists m,\exists n\leq
N:\|\hat{\boldsymbol{\theta}}_{n-1,m}-\boldsymbol{\theta}_{m}\|_{\Sigma_{n-1}}\geq\beta\right)\leq\frac{1}{N}$
(10)

Let $E_{n}$ denote the event that $\forall
m,\|\hat{\boldsymbol{\theta}}_{n-1,m}-\boldsymbol{\theta}_{m}\|_{\Sigma_{n-1}}\leq\beta$
and let $E=\bigcap_{n=1}^{N}E_{n}$. Define
$H_{n}=\\{(\boldsymbol{\phi}_{n^{\prime}},\mathbf{y}_{n^{\prime}})\\}_{n^{\prime}=1}^{n}$
to be the history up until time $n$. We have that

$\displaystyle\text{BR}(N)$
$\displaystyle=\mathbb{E}\left[\sum_{n=1}^{N}r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\right]$
(11)
$\displaystyle=\mathbb{E}\left[\mathbbm{1}\\{E^{c}\\}\sum_{n=1}^{N}r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\right]+\mathbb{E}\left[\mathbbm{1}\\{E\\}\sum_{n=1}^{N}r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\right]$
(12)
$\displaystyle=2LN\cdot\mathbb{E}\left[\mathbbm{1}\\{E^{c}\\}\right]+\mathbb{E}\left[\mathbbm{1}\\{E\\}\sum_{n=1}^{N}r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\right]$
(13) $\displaystyle\leq
2L+\mathbb{E}\left[\sum_{n=1}^{N}\mathbb{E}\Big{[}\mathbbm{1}\\{E_{n}\\}\left(r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\right)\
|\ H_{n-1}\Big{]}\right]$ (14)

Note that the first term in Eq. 12 is bounded by $2LN$ since
$|r_{n}(s,a)-r_{n}(s^{\prime},a^{\prime})|\leq
L\sum_{k}w_{n,k}|y_{m_{k}}(s,a)-y_{m_{k}}(s^{\prime},a^{\prime})|\leq 2L$.

Define $a_{s_{n}}^{\star}=\mathop{\mathrm{arg\,max}}_{a}r_{n}(s_{n},a)$ to be
the optimal action in state $s_{n}$. Noting that
$P(a_{s_{n}}^{\star}=\cdot|H_{n})=P(a_{n}=\cdot|H_{n})$, we have that

$\displaystyle\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\Big{(}r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\Big{)}\
|\ H_{n-1}\right]$ (15)
$\displaystyle=\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\Big{(}r_{n}(s_{n},a_{s_{n}}^{*})-\text{UCB}_{n}(s_{n},a_{n})+\text{UCB}_{n}(s_{n},a_{n})-r_{n}(s_{n},a_{n})\Big{)}\
|\ H_{n-1}\right]$ (16)
$\displaystyle=\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\Big{(}r_{n}(s_{n},a_{s_{n}}^{*})-\text{UCB}_{n}(s_{n},a_{s_{n}}^{*})+\text{UCB}_{n}(s_{n},a_{n})-r_{n}(s_{n},a_{n})\Big{)}\
|\ H_{n-1}\right]$ (17)
$\displaystyle\leq\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\Big{(}\text{UCB}_{n}(s_{n},a_{n})-r_{n}(s_{n},a_{n})\Big{)}\
|\ H_{n-1}\right]$ (18)

The inequality in Eq. 18 follows from the fact that $\text{UCB}_{n}(s,a)\geq
r_{n}(s,a)$ given $E_{n}$. Continuing, we find that

$\displaystyle\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\Big{(}\text{UCB}_{n}(s_{n},a_{n})-r_{n}(s_{n},a_{n})\Big{)}\
|\ H_{n-1}\right]$ (19)
$\displaystyle\leq\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\left(\sum_{k=1}^{K}w_{n,k}\left(U_{n,k}(\hat{y}_{n,m_{k}},\boldsymbol{\phi}_{n})-U_{n,k}(y_{n,m_{k}},\boldsymbol{\phi}_{n})\right)+\beta
L\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\right)\ |\ H_{n-1}\right]$ (20)
$\displaystyle\leq\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\left(\sum_{k=1}^{K}w_{n,k}L|\hat{y}_{n,m_{k}}-y_{n,m_{k}}|+\beta
L\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\right)\ |\ H_{n-1}\right]$ (21)
$\displaystyle=\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\left(\sum_{k=1}^{K}w_{n,k}L|\langle\boldsymbol{\phi}_{n},\hat{\boldsymbol{\theta}}_{n-1,m_{k}}-\boldsymbol{\theta}_{m_{k}}\rangle|+\beta
L\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\right)\ |\ H_{n-1}\right]$ (22)
$\displaystyle\leq\mathbb{E}\left[\mathbbm{1}\\{E_{n}\\}\left(\sum_{k=1}^{K}w_{n,k}L\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\|\hat{\boldsymbol{\theta}}_{n-1,m_{k}}-\boldsymbol{\theta}_{m_{k}}\|_{\Sigma_{n-1}}+\beta
L\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\right)\ |\ H_{n-1}\right]$ (23)
$\displaystyle\leq\mathbb{E}\left[2L\beta\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\
|\ H_{n-1}\right]$ (24)

Plugging the above into Eq. 14 combined with
$\mathbb{E}\left[r_{n}(s_{n},a_{s_{n}}^{*})-r_{n}(s_{n},a_{n})\right]\leq 2L$,
we conclude that

$\displaystyle\text{BR}(N)$ $\displaystyle\leq
2L+2L\beta\cdot\mathbb{E}\left[\sum_{n=1}^{N}\min\left(1,\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}\right)\right]$
(25) $\displaystyle\leq
2L+2L\beta\cdot\sqrt{N\cdot\mathbb{E}\left[\sum_{n=1}^{N}\min\left(1,\|\boldsymbol{\phi}_{n}\|_{\Sigma_{n-1}^{-1}}^{2}\right)\right]}$
(26) $\displaystyle\leq
2L+2L\beta\cdot\sqrt{2Nd\log\left(1+\frac{NQ^{2}}{\lambda d}\right)}$ (27)
$\displaystyle=2L\left(1+\left(\sqrt{\lambda}S+\sqrt{2\log\left(NM\right)+d\log\left(1+\frac{NQ^{2}}{\lambda
d}\right)}\right)\sqrt{2Nd\log\left(1+\frac{NQ^{2}}{\lambda d}\right)}\right)$
(28)

Eq. 26 follows from Cauchy-Schwarz and Eq. 27 follows from the elliptical
potential lemma (Lattimore and Szepesvári (2020) Lemma 19.4).

In conclusion, under the assumption that $Q,S,\lambda$ are $O(1)$, we have
that the Bayesian regret is of order $O(Ld\sqrt{N\log(NM)\log(N/d)})$

##  Appendix B Step Count Simulator

Here we discuss the details of the simulator. We model one outcome $(M=1)$,
step count, and assume the step count outcome follows a linear model,

$y_{i,t}=\boldsymbol{\phi}_{i,t}^{\top}\boldsymbol{\theta}^{\star}+\varepsilon$
(29)

where $\varepsilon\sim\mathcal{N}(0,0.15^{2})$. The featurization
$\boldsymbol{\phi}_{i,t}=[1,y_{i,t-1},t_{a},t_{i},nn]$ contains two components
that represent the treatment effect for each of the two groups ($t_{a}$ for
active users and $t_{i}$ for inactive users) and one counter that accumulates
the number of notifications ($nn$). The treatment effect is represented as a
decaying sigmoid of the number of notifications. For the active group it is

$t_{a}(nn)=1-\sigma\left(\left(nn-\frac{T}{1.9}\right)\cdot\frac{10}{T}\right)$
(30)

and for the inactive group it is

$t_{i}(nn)=1-\sigma\left(\left(nn-\frac{T}{1.3}\right)\cdot\frac{4}{T}\right)$
(31)

We set $\boldsymbol{\theta}^{\star}=[-0.02,0.999,0.3,0.21,0]$, which ignores
the number of notifications that have been sent to a user.

Figure 3: Comparing treatment effect between active and inactive groups.

The first utility function upweights step count and is given by

$U_{i,t,0}(y_{i,t},\boldsymbol{\phi}_{i,t})=\begin{cases}\alpha_{1}\cdot
y_{i,t}&y_{i,t}\leq\text{step goal}\\\ \alpha_{2}\cdot
y_{i,t}+\gamma&y_{i,t}>\text{step goal}\end{cases}$ (32)

where $\alpha_{1}=0.005,\alpha_{2}=0.001$, $\gamma=0.42$ if the user is an
active user, and $\gamma=0.3$ if the user is an inactive user. The step count
goal for active users is $10,000$ steps, and for inactive users it is $5,600$
steps. This utility is a piecewise linear function with a high slope before a
user achieves their step count goal, and a much smaller slope afterwards. The
second utility function is given by

$U_{i,t,1}(y_{i,t},\boldsymbol{\phi}_{i,t})=-\beta((nn)^{2})$ (33)

where $\beta=0.00003$. This utility decreases quadratically as the number of
notifications increases.


#  Long Range Azimuthal Correlation, Entanglement and Bell Inequality
Violation by Spinning Gluons at the LHC

Yuxun Guo<EMAIL_ADDRESS>Nuclear Science Division, Lawrence Berkeley
National Laboratory, Berkeley, CA 94720, USA  Xiaohui Liu<EMAIL_ADDRESS>Center of Advanced Quantum Studies, School of Physics and Astronomy, Beijing
Normal University, Beijing, 100875, China  Key Laboratory of Multi-scale Spin
Physics, Ministry of Education, Beijing Normal University, Beijing 100875,
China  Feng Yuan<EMAIL_ADDRESS>Nuclear Science Division, Lawrence Berkeley
National Laboratory, Berkeley, CA 94720, USA  Institute for Theoretical
Physics, Universität Tübingen, Auf der Morgenstelle 14, D-72076 Tübingen,
Germany  Hua Xing Zhu<EMAIL_ADDRESS>School of Physics, Peking University,
Beijing 100871, China  Center for High Energy Physics, Peking University,
Beijing 100871, China

###### Abstract

We apply the recently developed concept of the nucleon energy-energy
correlator (NEEC) for the gluon sector to investigate the long-range azimuthal
angular correlations in proton-proton collisions at the LHC. The spinning
gluon in these collisions will introduce a significant nonzero $\cos(2\phi)$
asymmetries in both Higgs Boson and top quark pair productions. The genesis of
the $\cos(2\phi)$ correlation lies in the intricate quantum entanglement.
Owing to the substantial $\cos(2\phi)$ effect, the NEEC observable in Higgs
Boson and $t{\bar{t}}$ production emerges as a pivotal avenue for delving into
quantum entanglement and scrutinizing the Bell inequality at high-energy
colliders.

Introduction. Long range correlation in particle productions in proton-proton
($pp$) collisions at the LHC has attracted great attention in the last decade
with tremendous efforts from both experiment and theory sides [1, 2, 3, 4]. In
this paper, we investigate this physics from a different perspective, applying
the nucleon energy-energy correlator (NEEC) [5, 6, 7] at the LHC. We will show
that the spinning gluon distribution in this framework [7] leads to sizable
$\cos(2\phi)$ azimuthal asymmetries in forward-backward energy correlators in
$pp$ collisions. These long range $\cos(2\phi)$ asymmetries are signatures of
the quantum entanglement, thereby providing the first test of the Bell
Inequality [8, 9] within the entangled gluon system. Pursuing such a test in
the Standard Model (SM) of particle physics at high energy colliders has been
very active in recent years [10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21,
22, 23].

The NEEC was introduced in [5] as a new method to explore the nucleon
structures. It employs an asymptotic energy flow operator $\hat{{\cal
E}}(\theta_{a})$ which measures energy deposits in the detector at a fixed
angle $\theta_{a}$ relative to the nucleon incoming beam direction in collider
experiments. Previous studies mainly focused on the deep inelastic scattering
(DIS) [5, 6, 24, 7] which will be explored at the future electron-ion collider
(EIC) [25, 26, 27]. In the following, we will study the NEEC observables in
$pp$ collisions. The comparison between these two collision systems will
provide an opportunity to test the universality of the NEECs. Meanwhile, the
novel phenomena unveiled below will stimulate further experiment
investigations and help decipher the origin of nearside ridge in $pp$
collisions.

Figure 1: Nucleon energy-energy correlator measurements in proton-proton
collisions at the LHC. Energy deposits in the forward directions of both
incoming hadron beams with polar angles $\theta_{a,b}$ and azimuthal angles
$\phi_{a,b}$ represented by $\vec{n}_{a,b}$, respectively.

To investigate the NEEC at the LHC, we propose to measure the energy deposits
along the beam directions of incoming hadrons with polar angles $\theta_{a,b}$
and azimuthal angles $\phi_{a,b}$, respectively, see the illustration in Fig.
1. The hard partonic scattering produces, e.g., the Higgs Boson, top quark
pair. The experiment can be carried out by a coincidence measurement between
the forward/backward energy flows and the hard interactions in the central.
Because $\theta_{a}$ and $\theta_{b}$ are both small and in opposite
directions, their rapidity difference will be large, for which we refer to as
a long range correlation. Meanwhile, we will show that different processes
lead to different $\cos(2\phi)$ asymmetries. Especially, we find that the
asymmetries in Higgs Boson and top quark pair productions are quite sizable
but with opposite signs. Therefore, a detailed study of these correlations
will open a new avenue for precision SM physics.

In the following, we focus on the gluon NEEC [7],

$\displaystyle f^{\alpha\beta}_{g,{\rm
EEC}}(x,\vec{n}_{a})=\int\frac{dy^{-}}{2\pi
xP^{+}}e^{-ixP^{+}\frac{y^{-}}{2}}$ $\displaystyle\hskip
17.22217pt\times\langle P|{\cal F}^{+\alpha}\left(y^{-}\right){\cal
L}^{\dagger}[\bm{\infty},y^{-}]\hat{{\cal E}}({\vec{n}}_{a}){\cal
L}[\bm{\infty},0]{\cal F}^{+\beta}(0)|P\rangle\,$
$\displaystyle=\left(-{g_{T}^{\alpha\beta}}/{2}\right)f_{g,{\rm
EEC}}+h_{T}^{\alpha\beta}d_{g,{\rm EEC}}\,,$ (1)

for the proton moving in $+\hat{z}$ direction with momentum $P$, where ${\cal
F}$ is the gauge field strength tensor, and ${\cal L}$ is the gauge link. We
have kept the azimuthal dependence of the energy flow direction
$n_{a}^{\alpha}=(1,\sin\theta_{a}\cos\phi_{a},\sin\theta_{a}\sin\phi_{a},\cos\theta_{a})$.
To parameterize the spinning gluon distribution, we introduced two projection
tensors:
$g_{T}^{\alpha\beta}=g^{\alpha\beta}-(P^{\alpha}{\bar{n}}^{\beta}+{\bar{n}}^{\alpha}P^{\beta})/{{\bar{n}}\cdot
P}$ and
$h_{T}^{\alpha\beta}=n_{a,T}^{\alpha}n_{a,T}^{\beta}/|n_{a,T}^{2}|+g_{T}^{\alpha\beta}/2$,
with ${\bar{n}}\cdot P=P^{0}+P^{z}\equiv P^{+}$ and
$n_{a,T}^{\alpha}=(0,\vec{n}_{a},0)$ is the transverse component of
$n_{a}^{\alpha}$. These two tensors help to define the normal gluon NEEC
$f_{g,\rm EEC}$ and the spinning gluon NEEC $d_{g,\rm EEC}$, respectively.
Similarly, we can define the gluon NEECs for the proton moving in $-\hat{z}$
direction with momentum ${P}^{\prime}$ and energy flow direction
$n_{b}^{\alpha}=(1,\sin\theta_{b}\cos\phi_{b},\sin\theta_{b}\sin\phi_{b},\cos\theta_{b})$.
The spinning gluon NEEC $d_{g,{\rm EEC}}(\theta_{a})$ originates from the
interference between different helicity states. To generate a long range
correlation between $\vec{n}_{a}$ and $\vec{n}_{b}$, we need to couple two
$d_{g,{\rm EEC}}(\theta)$ from both incoming protons, resulting into a
$\cos(2\phi)$ asymmetry, where $\phi=\phi_{a}-\phi_{b}$.

The spinning gluon distributions of the nucleon has also been studied in the
literature under different context. In the generalized parton distribution
(GPD) framework [28, 29, 30, 31], the spinning gluon GPD, also called
helicity-flip gluon GPD, predicts a $\cos(2\phi)$ asymmetry in the exclusive
processes [32, 33, 34]. Meanwhile, in the transverse momentum dependent (TMD)
formalism, the spinning gluon distribution, referred as the linearly polarized
gluon distribution, leads to a $\cos(2\phi)$ asymmetry in the associated TMD
processes [32, 33, 34, 35, 36, 37, 38, 39, 40, 41]. More recently, the
$\cos(2\phi)$ asymmetry has also been discussed in the context of jet
substructure [42, 43, 44, 45]. The comparison of these measurements will help
us understand the QCD dynamics associated with the spinning gluon.

NEEC for Higgs Boson and top quark pair processes at the LHC. The
factorization for NEEC in $pp$ collisions is similar to that for the DIS
processes [7]. As shown in Fig. 1, we measure the energy flows in $2$
arbitrary pixels on the calorimeter located at
${\vec{n}}_{a}=\sin\theta_{a}(\cos\phi_{a},\sin\phi_{a})$ and
${\vec{n}}_{b}=\sin\theta_{b}(\cos\phi_{b},\sin\phi_{b})$. The polar angles
are measured with respect to the $z$-axis, i.e., the particular rapidities,
and the azimuthal angles are measured from the transverse plane perpendicular
to the beam direction. We require each of the two pixels much closer to one of
the hadron beams. Therefore, these two particles are in opposite directions,
forward/backward in the Lab frame, e.g., $\theta_{a}\to 0$ and
$\theta_{b}\to\pi$. The generic cross section measurement takes the following
form

$\displaystyle\Sigma(Q^{2};\theta_{a,b},\phi)$ $\displaystyle=$
$\displaystyle\sum_{ij}\int
d\sigma(Q^{2})\frac{E_{i}}{E_{P}}\frac{E_{j}}{E_{P}}{\cal
F}(\phi;\vec{n}_{a,b})$ (2)
$\displaystyle~{}~{}~{}~{}~{}~{}\times\delta(\vec{n}_{a}-\vec{n}_{i})\delta(\vec{n}_{b}-\vec{n}_{j})\,,$

where ${\cal F}(\phi;\vec{n}_{a,b})$ imposes the phase space measurement to
construct $\phi$. In the above equation, $d\sigma(Q)$ represents partonic
scattering cross section. Following previous examples, the factorization
formula can be written as,

$\displaystyle\Sigma(Q^{2};\theta_{a,b},\phi)$ $\displaystyle=\int
d\Omega\left\\{x_{a}f_{g,{\rm
EEC}}\left(x_{a},\theta_{a}^{2}\right)x_{b}f_{g,{\rm
EEC}}\left(x_{b},\theta_{b}^{2}\right)\hat{\sigma}_{0}\right.$ (3)
$\displaystyle\left.+x_{a}d_{g,{\rm
EEC}}\left(x_{a},\theta_{a}^{2}\right)x_{b}d_{g,{\rm
EEC}}\left(x_{b},\theta_{b}^{2}\right)\hat{\sigma}_{2}(Q^{2})\cos(2\phi)\right\\}\
,$

where $Q^{2}=x_{a}x_{b}S_{pp}$ with $S_{pp}$ the center of mass energy
squared, $d\Omega$ represents additional phase space integral.
$\hat{\sigma}_{0,2}$ are partonic cross sections calculable perturbatively.
The $\cos(2\phi)$ term $\hat{\sigma}_{2}$ comes from the interference between
double helicity-flip amplitudes where both incoming gluons have the same
helicity as illustrated in Fig. 2.

Figure 2: Long range $\cos(2\phi)$ asymmetry comes from the interference
between double helicity-flip amplitudes in the partonic scattering processes.

Of course, the final results of $\cos(2\phi)$ asymmetries also depend on the
NEEC gluon distributions. When $P^{+}\theta_{a}\gg\Lambda_{\rm QCD}$, they can
be computed from perturbative QCD with collinear splitting contributions,

$\displaystyle f_{g,{\rm
EEC}}(x,\theta_{a}^{2})=\frac{\alpha_{s}}{2\pi}\frac{1}{\theta_{a}^{2}}\int_{x}^{1}\frac{dz}{z}\frac{x(1-z)}{z}$
$\displaystyle\hskip 43.05542pt\times\left[{\cal
P}_{g/q}(z)f_{q}\left(\frac{x}{z}\right)+{\cal
P}_{g/g}(z)f_{g}\left(\frac{x}{z}\right)\right]\,,$ (4) $\displaystyle
d_{g,{\rm
EEC}}(x,\theta_{a}^{2})=\frac{\alpha_{s}}{2\pi}\frac{1}{\theta_{a}^{2}}\int_{x}^{1}\frac{dz}{z}\frac{x(1-z)}{z}$
$\displaystyle\hskip
43.05542pt\times\frac{2(1-z)}{z}\left[C_{F}f_{q}\left(\frac{x}{z}\right)+C_{A}f_{g}\left(\frac{x}{z}\right)\right]\,,$
(5)

where ${\cal P}_{g/q}$ and ${\cal P}_{g/g}$ are usual collinear splitting
kernels. Additional DGLAP resummation will modify the power behavior, for
which we expect a similar effect for both $f_{g,\rm EEC}$ and $d_{g,\rm EEC}$
[42]. In the following, we will apply the above leading order results to
demonstrate the azimuthal asymmetries for the Higgs Boson production and top
quark pair production in $pp$ collisions at the LHC.

To study the spinning gluon effect at the LHC, the simplest process is the
Higgs Boson production. Similar to the TMD case calculated before [46, 47],
the Higgs Boson can couple to the spinning gluons directly, and at the leading
order

$\hat{\sigma}_{2}=\hat{\sigma}_{0}=\pi g_{\phi}^{2}/64\ ,$ (6)

where $g_{\phi}$ represents the coupling between the Higgs Boson and the gluon
fields in the effective theory ${\cal L}_{eff}=-\frac{1}{4}g_{\phi}\Phi
F^{a}_{\mu\nu}F^{a\mu\nu}$ [48]. The above shows that the $\cos(2\phi)$
asymmetry for Higgs production is positive and can reach a sizable value
depending on the ratio between $d_{g,\rm EEC}$ and $f_{g\rm EEC}$. A similar
$\cos(2\phi)$ asymmetry has also been found for Higgs plus two jets
production, where $\phi$ is the azimuthal angle between the two jets [49]. In
the common kinematics, the physics behind these two $\cos(2\phi)$ is the same,
originating from the spinning gluon.

On the other hand, for the top quark pair production, $\hat{\sigma}_{2}$ is
different from $\hat{\sigma}_{0}$,

$\displaystyle\hat{\sigma}_{0}$ $\displaystyle=$
$\displaystyle\frac{\alpha_{s}^{2}\pi}{\hat{s}^{2}}\left[\frac{1}{6}\frac{1}{\hat{t}_{1}\hat{u}_{1}}-\frac{3}{8}\frac{1}{\hat{s}^{2}}\right]\left[\hat{t}_{1}^{2}+\hat{u}_{1}^{2}+4m_{t}^{2}\hat{s}-\frac{4m_{t}^{4}\hat{s}^{2}}{\hat{t}_{1}\hat{u}_{1}}\right]$
$\displaystyle\hat{\sigma}_{2}$ $\displaystyle=$
$\displaystyle\frac{\alpha_{s}^{2}\pi}{\hat{s}^{2}}\left[\frac{3}{8}\frac{1}{\hat{s}^{2}}-\frac{1}{6}\frac{1}{\hat{t}_{1}\hat{u}_{1}}\right]\frac{2m_{t}^{4}\hat{s}^{2}}{\hat{t}_{1}\hat{u}_{1}}\
,$ (7)

for the dominant $gg\to t\bar{t}$ channel, where
$\hat{t}_{1}=\hat{t}-m_{t}^{2}$ and $\hat{u}_{1}=\hat{u}-m_{t}^{2}$,
$\hat{s}$, $\hat{t}$ and $\hat{u}$ are usual Mandelstam variables. Contrary to
the Higgs case, the $\cos(2\phi)$ asymmetry for top quark pair production is
negative. Interestingly, the asymmetry will reach the maximum value when the
pair are close to the threshold where $\hat{s}=4m_{t}^{2}$. The opposite sign
of $\cos(2\phi)$ asymmetries between these two processes is due to their
difference in parity: Higgs Boson is parity even while the threshold top quark
pair is parity odd. As a result, they couple to the incoming gluon helicity
states with a different sign between $|++\rangle$ and $|--\rangle$, which
leads to the opposite $\cos(2\phi)$ asymmetries.

Figure 3: Long range $\cos(2\phi)$ azimuthal asymmetries associated with Higgs
boson production and top quark pair threshold production as functions of their
rapidity $y$.

In Fig. 3, we show the $\cos(2\phi)$ asymmetries as functions of rapidity in
Higgs boson production and threshold top quark pair production. From this
plot, we find that both asymmetries are quite sizable at mid-rapidity. They
decrease with rapidity, which reflects $x$-dependence of the spinning gluon
and the normal gluon distributions as described in Eqs. (4,5). The experiment
measurements of these asymmetries will provide important constraints on the
gluon spinning effects.

Figure 4: $\cos(2\phi)$ azimuthal asymmetries in the NEEC observable
associated with top quark pair production as functions of the rapidity
difference between the pair $\Delta y$ at different $p_{T}$.

As shown in Fig. 4, the $\cos(2\phi)$ asymmetry also depends on the transverse
momentum and the rapidity difference between the pair $\Delta
y=y_{t}-y_{\bar{t}}$ for the top quark pair. We have also computed two-photon
production through the gluon-gluon fusion process by applying the amplitudes
derived in the literature [50, 51, 52] and the $\cos(2\phi)$ asymmetry is
smaller as compared to that of Higgs Boson production with an opposite sign.

These results demonstrate that the $\cos(2\phi)$ asymmetries can provide a
strong case to study the spinning gluon physics at the LHC. More importantly,
this shall open a new avenue to study precision physics in the SM. It may also
lead to a unique probe of new physics beyond the SM. Especially, the
asymmetries crucially depend on the couplings between the gluon fields with
different helicities and the Higgs Boson, which have been argued to be
sensitive to the new physics beyond the SM and similar studies on the TMD
related observable have been investigated in Refs. [47, 53].

Although the above results are based on the leading order calculations, we
expect that higher order corrections will not modify the large $\cos(2\phi)$
asymmetries for the above processes. For example, in Higgs production at one-
loop order, the soft gluon radiation leads to the same double logarithmic
contributions for both $\sigma_{2}$ and $\sigma_{0}$. As a result, the
$\cos(2\phi)$ will remain the same after resummation. Studies on azimuthal
asymmetry between the jets in Higgs plus two jets production found mild
dependence on both higher $\alpha_{s}$ order corrections [54] and parton
showers [55]. Therefore, we anticipate this attribute to persist for NEEC.

Quantum entanglement and test of Bell inequality. The $\cos(2\phi)$
correlation can be interpreted as a signature of entanglement. In an
experiment, what is being measured are the real particles that hit the forward
detectors. Although these forward-moving particles never come into contact,
they remain entangled in their helicities. The physics picture is as follows.
Two pairs of entangled real particles and virtual gluons are created though
the splitting of the incoming partons. At this point, the two real particles
are separable entities. The virtual gluons will participate the partonic hard
process, while the real particles will travel towards the forward detectors at
opposite ends of the beam with large momentum $E\sim P_{z}\gg P_{t}\sim
E\theta$. Once the hard process entangles the virtual gluons, it can be
demonstrated that the two real particles become entangled instantaneously.
Notably, the two real particles are always separated by space-like intervals,
since $\Delta s^{2}=\Delta t^{2}-\Delta r_{t}^{2}-\Delta
z^{2}\sim\frac{1}{E^{2}}-\frac{1}{E^{2}\theta^{2}}-\frac{1}{P_{z}^{2}}<0$, due
to $\theta\ll 1$ in the forward limit. This implies the causal disconnection
between the detected particles, and any non-trivial correlation between them
should signify their entanglement. This remarkable property also allows us to
conduct localized measurements when testing the Bell inequality.

This observation provides a basis for testing Bell’s theorem [8] through the
$\cos(2\phi)$ correlation. Leveraging the NEEC in Eq. (Long Range Azimuthal
Correlation, Entanglement and Bell Inequality Violation by Spinning Gluons at
the LHC), one can formulate the Bell observable

$\displaystyle
S(\phi_{a},\phi_{b})\equiv\frac{\Sigma(\phi_{a},\phi_{b})+\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})-\Sigma(\phi_{a}^{\prime},\phi_{b})-\Sigma(\phi_{a},\phi_{b}^{\prime})}{\Sigma(\phi_{a},\phi_{b})+\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})+\Sigma(\phi_{a}^{\prime},\phi_{b})+\Sigma(\phi_{a},\phi_{b}^{\prime})}$
(8)

where $\phi_{a}$ and $\phi_{b}$ are azimuthal angles of the energy flow
directed towards the detector, measured with respect to arbitrary reference
vectors $r_{a,b}$. $\phi^{\prime}=\phi+\frac{\pi}{2}$ and can be regarded as
one measures the azimuthal angles with the reference vectors perpendicular to
$r_{a,b}$. For appropriate choices of $\phi_{a,b}$, $\tilde{\phi}_{a,b}$, the
Clauser-Horne-Shimony-Holt (CHSH) inequality [9], an equivalent version of the
Bell’s original inequality,

$\displaystyle
B\equiv|S(\phi_{a},\phi_{b})-S(\phi_{a},\tilde{\phi}_{b})+S(\tilde{\phi}_{a},\phi_{b})+S(\tilde{\phi}_{a},\tilde{\phi}_{b})|\leq
2\,.$

can potentially be violated. The maximum violation of the CHSH inequality for
any quantum state is given by the Tsirelson’s bound, $B_{\rm
max}=2\sqrt{2}\approx 2.828$ [56]. A proof of Eq. (Long Range Azimuthal
Correlation, Entanglement and Bell Inequality Violation by Spinning Gluons at
the LHC) can be found in the supplementary material.

Fig. 5 demonstrates the concept by measuring the CHSH inequality in Eq. (Long
Range Azimuthal Correlation, Entanglement and Bell Inequality Violation by
Spinning Gluons at the LHC) using NEEC factorization in Eq. (Long Range
Azimuthal Correlation, Entanglement and Bell Inequality Violation by Spinning
Gluons at the LHC). We choose $\phi_{a}=0$, $\phi_{b}=\frac{\pi}{8}$,
$\tilde{\phi}_{a}=\frac{\pi}{4}$ and $\tilde{\phi}_{b}=\frac{3\pi}{8}$.
Violation of the CHSH inequality is observed for the Higgs rapidity $y_{\rm
Higgs}<0.5$. We note that the significance can be dramatically improved by
quark jet tagging as manifest from Fig. 5 where the CHSH inequality violation
is observed for both Higgs and $t{\bar{t}}$ threshold production. We also
check that increasing the machine energy leads to a more significant
violation, reaching $B\approx 2.36$ for $y_{\rm Higgs}=0$ at
$\sqrt{S_{pp}}=33\>{\rm TeV}$ without jet tagging, as the entanglement between
the detected forward-moving particles intensifies near small $x$ values.

Figure 5: Violation of the CHSH inequality in Higgs (red) and top pair (blue)
production at the LHC. Quark jet tagging (dashed lines) substantially enhances
the significance.

Conclusion. In summary, we studied the long-range azimuthal angular
correlations in NEEC measurements in $pp$ collisions at the LHC. For a number
of processes, we found significant large $\cos(2\phi)$ asymmetries. The
comparison between these and future studies at the EIC will provide an
important test of the universality of the NEEC distribution functions. Because
of large asymmetries in these processes at the LHC, we emphasize that this
will also open a new avenue to study precision physics for the Standard Model
(SM), in particular, through comparison between Higgs Boson production, top
quark pair production, and two-photon production.

The connection between the $\cos(2\phi)$ correlation and the entanglement
makes the long-range correlation in NEEC a promising alternative for testing
the Bell Inequality, one of the most fundamental properties of quantum
mechanics. We demonstrate the feasibility of this approach using Higgs and
threshold $t{\bar{t}}$ production at the LHC in which violation of the Bell
inequality is significant when we perform quark jet tagging. Compared to the
other collider-based tests discussed in the literature [10, 11, 12, 13, 14,
15, 16, 17, 18, 19, 20, 21, 22, 23], the long-range correlation in NEEC
enables, for the first time, a test of this fundamental quantum property in
confined quantities like gluons. Our method benefits from the NEEC
factorization theorem, ensuring the test remains local, thus closing the major
potential loophole [22] present in LHC-based tests. Moreover, unlike previous
proposals that often require reconstructing the full kinematics which is
usually challenging at the LHC, the NEEC measurement only requires determining
the azimuthal angles of the energy flow deposit at the forward detectors,
making it more practical for experimental implementation.

Looking ahead, extending this research to other QCD processes, including
multi-jet production, and heavy quarkonium production, will be interesting to
follow. Additionally, recent investigations [57, 58, 59, 60, 61, 62, 63, 64,
65, 66, 67] have indicated that the quantum entanglement may bring novel
perspectives into nuclear and particle physics. We thus anticipate our work
may spark similar endeavors in unraveling the nucleon structures using the
quantum information properties. These studies will promise to yield deeper
insights into the effects of spinning gluons, complement our current
understanding, and potentially reveal new physics beyond the SM.

###### Acknowledgements.

Acknowledgement. We thank Meng Xiao for discussions on jet tagging. This work
is supported by the Natural Science Foundation of China under contract No.
12175016 (X. L.), the Office of Science of the U.S. Department of Energy under
Contract No. DE-AC02-05CH11231 and under the umbrella of the Quark-Gluon
Tomography (QGT) Topical Collaboration with Award DE-SC0023646 (Y.G and F.Y.),
the Startup Grant of Peking University, and the Asian Young Scientist
Fellowship (H. X. Z.).

## References

  * [1] K. Dusling, W. Li, and B. Schenke,  Int. J. Mod. Phys. E 25, 1630002 (2016), 1509.07939.
  * [2] C. Loizides,  Nucl. Phys. A 956, 200 (2016), 1602.09138.
  * [3] M. Strickland,  Nucl. Phys. A 982, 92 (2019), 1807.07191.
  * [4] J. L. Nagle and W. A. Zajc,  Ann. Rev. Nucl. Part. Sci. 68, 211 (2018), 1801.03477.
  * [5] X. Liu and H. X. Zhu,  Phys. Rev. Lett. 130, 091901 (2023), 2209.02080.
  * [6] H. Cao, X. Liu, and H. X. Zhu,  Phys. Rev. D 107, 114008 (2023), 2303.01530.
  * [7] X. L. Li, X. Liu, F. Yuan, and H. X. Zhu,  Phys. Rev. D 108, L091502 (2023), 2308.10942.
  * [8] J. S. Bell,  Physics Physique Fizika 1, 195 (1964).
  * [9] J. F. Clauser, M. A. Horne, A. Shimony, and R. A. Holt,  Phys. Rev. Lett. 23, 880 (1969).
  * [10] M. Fabbrichesi, R. Floreanini, and G. Panizzo,  Phys. Rev. Lett. 127, 161801 (2021), 2102.11883.
  * [11] C. Severi, C. D. E. Boschi, F. Maltoni, and M. Sioli,  Eur. Phys. J. C 82, 285 (2022), 2110.10112.
  * [12] A. J. Barr,  Phys. Lett. B 825, 136866 (2022), 2106.01377.
  * [13] J. A. Aguilar-Saavedra and J. A. Casas,  Eur. Phys. J. C 82, 666 (2022), 2205.00542.
  * [14] J. A. Aguilar-Saavedra,  Phys. Rev. D 107, 076016 (2023), 2209.14033.
  * [15] R. Ashby-Pickering, A. J. Barr, and A. Wierzchucka,  JHEP 05, 020 (2023), 2209.13990.
  * [16] M. Fabbrichesi, R. Floreanini, E. Gabrielli, and L. Marzola,  Eur. Phys. J. C 83, 823 (2023), 2302.00683.
  * [17] J. A. Aguilar-Saavedra,  Phys. Rev. D 108, 076025 (2023), 2307.06991.
  * [18] Q. Bi, Q.-H. Cao, K. Cheng, and H. Zhang,  Phys. Rev. D 109, 036022 (2024), 2307.14895.
  * [19] T. Han, M. Low, and T. A. Wu,  (2023), 2310.17696.
  * [20] K. Ma and T. Li,  (2023), 2309.08103.
  * [21] A. Bernal, P. Caban, and J. Rembieliński,  (2024), 2405.16525.
  * [22] A. J. Barr, M. Fabbrichesi, R. Floreanini, E. Gabrielli, and L. Marzola,  (2024), 2402.07972.
  * [23] F. Maltoni, C. Severi, S. Tentori, and E. Vryonidou,  (2024), 2404.08049.
  * [24] H.-Y. Liu, X. Liu, J.-C. Pan, F. Yuan, and H. X. Zhu,  Phys. Rev. Lett. 130, 181901 (2023), 2301.01788.
  * [25] A. Accardi et al.,  Eur. Phys. J. A 52, 268 (2016), 1212.1701.
  * [26] R. Abdul Khalek et al.,  Nucl. Phys. A 1026, 122447 (2022), 2103.05419.
  * [27] Proceedings, Probing Nucleons and Nuclei in High Energy Collisions: Dedicated to the Physics of the Electron Ion Collider: Seattle (WA), United States, October 1 - November 16, 2018, WSP, 2020, 2002.12333.
  * [28] X.-D. Ji,  Phys. Rev. Lett. 78, 610 (1997), hep-ph/9603249.
  * [29] D. Müller, D. Robaschik, B. Geyer, F. M. Dittes, and J. Hořejši,  Fortsch. Phys. 42, 101 (1994), hep-ph/9812448.
  * [30] X.-D. Ji,  Phys. Rev. D 55, 7114 (1997), hep-ph/9609381.
  * [31] A. V. Radyushkin,  Phys. Rev. D 56, 5524 (1997), hep-ph/9704207.
  * [32] M. Diehl, T. Gousset, B. Pire, and J. P. Ralston,  Phys. Lett. B 411, 193 (1997), hep-ph/9706344.
  * [33] P. Hoodbhoy and X.-D. Ji,  Phys. Rev. D 58, 054006 (1998), hep-ph/9801369.
  * [34] A. V. Belitsky and D. Mueller,  Phys. Lett. B 486, 369 (2000), hep-ph/0005028.
  * [35] D. Boer, S. J. Brodsky, P. J. Mulders, and C. Pisano,  Phys. Rev. Lett. 106, 132001 (2011), 1011.4225.
  * [36] A. Metz and J. Zhou,  Phys. Rev. D 84, 051503 (2011), 1105.1991.
  * [37] C. Pisano, D. Boer, S. J. Brodsky, M. G. A. Buffing, and P. J. Mulders,  JHEP 10, 024 (2013), 1307.3417.
  * [38] Y. Hatta, B.-W. Xiao, F. Yuan, and J. Zhou,  Phys. Rev. Lett. 126, 142001 (2021), 2010.10774.
  * [39] Y. Hatta, B.-W. Xiao, F. Yuan, and J. Zhou,  Phys. Rev. D 104, 054037 (2021), 2106.05307.
  * [40] R. Esha, Z.-B. Kang, K. Lee, D. Y. Shao, and F. Zhao,  (2022).
  * [41] P. Caucal, F. Salazar, B. Schenke, T. Stebel, and R. Venugopalan,  (2023), 2308.00022.
  * [42] H. Chen, I. Moult, and H. X. Zhu,  Phys. Rev. Lett. 126, 112003 (2021), 2011.02492.
  * [43] H. Chen, I. Moult, and H. X. Zhu,  JHEP 08, 233 (2022), 2104.00009.
  * [44] A. Karlberg, G. P. Salam, L. Scyboz, and R. Verheyen,  Eur. Phys. J. C 81, 681 (2021), 2103.16526.
  * [45] Z. Yu and C. P. Yuan,  Phys. Rev. Lett. 129, 112001 (2022), 2110.11539.
  * [46] P. Sun, B.-W. Xiao, and F. Yuan,  Phys. Rev. D 84, 094005 (2011), 1109.1354.
  * [47] D. Boer, W. J. den Dunnen, C. Pisano, M. Schlegel, and W. Vogelsang,  Phys. Rev. Lett. 108, 032002 (2012), 1109.1444.
  * [48] S. Dawson,  Nucl. Phys. B 359, 283 (1991).
  * [49] T. Plehn, D. L. Rainwater, and D. Zeppenfeld,  Phys. Rev. Lett. 88, 051801 (2002), hep-ph/0105325.
  * [50] D. A. Dicus and S. S. D. Willenbrock,  Phys. Rev. D 37, 1801 (1988).
  * [51] Z. Bern, A. De Freitas, and L. J. Dixon,  JHEP 09, 037 (2001), hep-ph/0109078.
  * [52] J.-W. Qiu, M. Schlegel, and W. Vogelsang,  Phys. Rev. Lett. 107, 062001 (2011), 1103.3861.
  * [53] D. Boer, W. J. den Dunnen, C. Pisano, and M. Schlegel,  Phys. Rev. Lett. 111, 032002 (2013), 1304.2654.
  * [54] J. M. Campbell, R. K. Ellis, and G. Zanderighi,  JHEP 10, 028 (2006), hep-ph/0608194.
  * [55] V. Del Duca et al.,  JHEP 10, 016 (2006), hep-ph/0608158.
  * [56] B. S. Cirelson,  Lett. Math. Phys. 4, 93 (1980).
  * [57] D. E. Kharzeev and E. M. Levin,  Phys. Rev. D 95, 114008 (2017), 1702.03489.
  * [58] S. R. Beane, D. B. Kaplan, N. Klco, and M. J. Savage,  Phys. Rev. Lett. 122, 102001 (2019), 1812.03138.
  * [59] Z. Tu, D. E. Kharzeev, and T. Ullrich,  Phys. Rev. Lett. 124, 062001 (2020), 1904.11974.
  * [60] N. Mueller, A. Tarasov, and R. Venugopalan,  Phys. Rev. D 102, 016007 (2020), 1908.07051.
  * [61] W. Gong, G. Parida, Z. Tu, and R. Venugopalan,  Phys. Rev. D 106, L031501 (2022), 2107.13007.
  * [62] N. Armesto, F. Dominguez, A. Kovner, M. Lublinsky, and V. Skokov,  JHEP 05, 025 (2019), 1901.08080.
  * [63] D. E. Kharzeev,  Phil. Trans. A. Math. Phys. Eng. Sci. 380, 20210063 (2021), 2108.08792.
  * [64] I. Low and T. Mehen,  Phys. Rev. D 104, 074014 (2021), 2104.10835.
  * [65] M. Carena, I. Low, C. E. M. Wagner, and M.-L. Xiao,  Phys. Rev. D 109, L051901 (2024), 2307.08112.
  * [66] K. Sakurai and M. Spannowsky,  Phys. Rev. Lett. 132, 151602 (2024), 2310.01477.
  * [67] M. Hentschinski, D. E. Kharzeev, K. Kutak, and Z. Tu,  Phys. Rev. Lett. 131, 241901 (2023), 2305.03069.

## Supplemental Material: Proof of the CHSH inequality in Eq. (9)

For the convenience of the readers, we provide a self-contained proof of the
CHSH inequality in Eq. (Long Range Azimuthal Correlation, Entanglement and
Bell Inequality Violation by Spinning Gluons at the LHC). First, we note that
the two-particle NEEC defined in Eq. (8) satisfies

$\displaystyle|S(\phi_{a},\phi_{b})|=$ $\displaystyle\
\frac{|\Sigma(\phi_{a},\phi_{b})+\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})-\Sigma(\phi_{a}^{\prime},\phi_{b})-\Sigma(\phi_{a},\phi_{b}^{\prime})|}{|\Sigma(\phi_{a},\phi_{b})+\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})+\Sigma(\phi_{a}^{\prime},\phi_{b})+\Sigma(\phi_{a},\phi_{b}^{\prime})|}$
$\displaystyle=$ $\displaystyle\
\frac{|\Sigma(\phi_{a},\phi_{b})+\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})-\Sigma(\phi_{a}^{\prime},\phi_{b})-\Sigma(\phi_{a},\phi_{b}^{\prime})|}{|\Sigma(\phi_{a},\phi_{b})|+|\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})|+|\Sigma(\phi_{a}^{\prime},\phi_{b})|+|\Sigma(\phi_{a},\phi_{b}^{\prime})|}$
$\displaystyle\leq$ $\displaystyle\
\frac{|\Sigma(\phi_{a},\phi_{b})|+|\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})|+|-\Sigma(\phi_{a}^{\prime},\phi_{b})|+|-\Sigma(\phi_{a},\phi_{b}^{\prime})|}{|\Sigma(\phi_{a},\phi_{b})|+|\Sigma(\phi_{a}^{\prime},\phi_{b}^{\prime})|+|\Sigma(\phi_{a}^{\prime},\phi_{b})|+|\Sigma(\phi_{a},\phi_{b}^{\prime})|}$
$\displaystyle=$ $\displaystyle\ 1\,.$ (10)

Therefore, $-1\leq S(\phi_{a},\phi_{b})\leq 1$. For the particular two-
particle NEEC generated by the spinning gluon effects, where
$\Sigma(\phi_{a},\phi_{b})=A_{0}+A_{2}\cos(2(\phi_{a}-\phi_{b}))$ and
$S(\phi_{a},\phi_{b})=\frac{A_{2}}{A_{0}}\cos(2(\phi_{a}-\phi_{b}))$ for
$\phi_{a,b}^{\prime}=\phi_{a,b}+\frac{\pi}{2}$, the inequality can be
saturated when $A_{2}=\pm A_{0}$. The key assumption in the proof of Eq. (Long
Range Azimuthal Correlation, Entanglement and Bell Inequality Violation by
Spinning Gluons at the LHC) is that the two-particle correlation is induced
through a set of hidden variables, collectively denoted as $\lambda$,

$S(\phi_{a},\phi_{b})=\int
d\lambda\,\rho(\lambda)S_{\lambda}(\phi_{a},\phi_{b})=\int
d\lambda\,\rho(\lambda)R_{a}(\phi_{a},\lambda)R_{b}(\phi_{b},\lambda)\,,$ (11)

where the last equation holds for local measurement, i.e., the measurement
events are separated by a space-like interval and the result of a measurement
on particle $a$ be unaffected by operations on the distant particle $b$, and
vice versa. Here,

$\displaystyle R_{a}(\phi_{a},\lambda)=\int
d\phi_{b}S_{\lambda}(\phi_{a},\phi_{b})\,,\qquad R_{b}(\phi_{a},\lambda)=\int
d\phi_{a}S_{\lambda}(\phi_{a},\phi_{b})\,,\quad$ (12)

where we note that $R_{b}(\phi_{b},\lambda)\neq R_{a}(\phi_{a},\lambda)$ in
general, and

$\int\rho(\lambda)d\lambda=1\,,\qquad\rho(\lambda)\geq
0\,,\qquad|R_{a,b}(\phi,\lambda)|\leq 1\,.$ (13)

We proceed by first prove that

$|S(\phi_{a},\phi_{b})-S(\phi_{a},\tilde{\phi}_{b})|+|S(\tilde{\phi}_{a},\phi_{b})+S(\tilde{\phi}_{a},\tilde{\phi}_{b})|\leq
2\,.$ (14)

Using (11), we have

$\displaystyle\
|S(\phi_{a},\phi_{b})-S(\phi_{a},\tilde{\phi}_{b})|+|S(\tilde{\phi}_{a},\phi_{b})+S(\tilde{\phi}_{a},\tilde{\phi}_{b})|$
$\displaystyle=$ $\displaystyle\ \int
d\lambda\,\rho(\lambda)\Big{(}|R_{a}(\phi_{a},\lambda)||R_{b}(\phi_{b},\lambda)-R_{b}(\tilde{\phi}_{b},\lambda)|+|R_{a}(\tilde{\phi}_{a},\lambda)||R_{b}(\phi_{b},\lambda)+R_{b}(\tilde{\phi}_{b},\lambda)|\Big{)}$
$\displaystyle\leq$ $\displaystyle\ \int
d\lambda\,\rho(\lambda)\Big{(}|R_{b}(\phi_{b},\lambda)-R_{b}(\tilde{\phi}_{b},\lambda)|+|R_{b}(\phi_{b},\lambda)+R_{b}(\tilde{\phi}_{b},\lambda)|\Big{)}\,.$
(15)

Without loss of generality, suppose that
$R_{b}(\phi_{b},\lambda)\geq|R_{b}(\tilde{\phi}_{b},\lambda)|\geq 0$. Then

$|R_{b}(\phi_{b},\lambda)-R_{b}(\tilde{\phi}_{b},\lambda)|+|R_{b}(\phi_{b},\lambda)+R_{b}(\tilde{\phi}_{b},\lambda)|=R_{b}(\phi_{b},\lambda)-R_{b}(\tilde{\phi}_{b},\lambda)+R_{b}(\phi_{b},\lambda)+R_{b}(\tilde{\phi}_{b},\lambda)=2R_{b}(\phi_{b},\lambda)\leq
2\,.$ (16)

Then we have

$|S(\phi_{a},\phi_{b})-S(\phi_{a},\tilde{\phi}_{b})|+|S(\tilde{\phi}_{a},\phi_{b})+S(\tilde{\phi}_{a},\tilde{\phi}_{b})|\leq
2\int\rho(\lambda)d\lambda=2\,,$ (17)

and (14) is proved. From (11), we can apply the triangular inequality to
obtain

$|S(\phi_{a},\phi_{b})-S(\phi_{a},\tilde{\phi}_{b})+S(\tilde{\phi}_{a},\phi_{b})+S(\tilde{\phi}_{a},\tilde{\phi}_{b})|\leq|S(\phi_{a},\phi_{b})-S(\phi_{a},\tilde{\phi}_{b})|+|S(\tilde{\phi}_{a},\phi_{b})+S(\tilde{\phi}_{a},\tilde{\phi}_{b})|\leq
2\,.$ (18)

This completes the proof of (Long Range Azimuthal Correlation, Entanglement
and Bell Inequality Violation by Spinning Gluons at the LHC).