# Optimizing Long-term Value for Auction-Based Recommender Systems via On-
Policy Reinforcement Learning

Ruiyang Xu  Jalaj Bhandari  Dmytro Korenkevych  Fan Liu  Yuchen He  Alex
Nikulkov  Zheqing Zhu

###### Abstract

Auction-based recommender systems are prevalent in online advertising
platforms, but they are typically optimized to allocate recommendation slots
based on immediate expected return metrics, neglecting the downstream effects
of recommendations on user behavior. In this study, we employ reinforcement
learning to optimize for long-term return metrics in an auction-based
recommender system. Utilizing temporal difference learning, a fundamental
reinforcement learning algorithm, we implement an one-step policy improvement
approach that biases the system towards recommendations with higher long-term
user engagement metrics. This optimizes value over long horizons while
maintaining compatibility with the auction framework. Our approach is grounded
in dynamic programming ideas which show that our method provably improves upon
the existing auction-based base policy. Through an online A/B test conducted
on an auction-based recommender system which handles billions of impressions
and users daily, we empirically establish that our proposed method outperforms
the current production system in terms of long-term user engagement metrics.

**footnotetext: These authors contributed equally to this work

##  1 Introduction

With exponential growth of digital information, recommender systems have come
to play a pivotal role in various applications, from personalized movie
recommendations at Netflix to product recommendation over e-commerce
platforms, to help users access relevant/interesting content (Gomez-Uribe and
Hunt, 2015; Smith and Linden, 2017; Lu et al., 2015). Given widespread use,
research and practice of recommender system design has constantly evolved in
the past two decades. These include traditional recommendation strategies,
like content based approaches (Pazzani and Billsus, 2007; Lops et al., 2011)
and collaborative filtering methods (Schafer et al., 2007; Ekstrand et al.,
2011; Shi et al., 2014; Koren et al., 2009) which model user preferences for
different items either using item features or by using past ratings of a user
(or similar users), to modern deep supervised learning based approaches which
personalize user experience by modeling click through rates (Zhang et al.,
2019), with many instances of successful deployment at industry scale. See
(Covington et al., 2016; Cheng et al., 2016; Okura et al., 2017) for example
of prototypes used for video recommendation at YouTube, app recommendations at
Google Play and news recommender systems for Yahoo!.

Auction based systems (Varian, 2007; Evans, 2008) is a critical component for
platforms like Meta, Google, Yahoo! etc. which run large marketplaces such as
online advertising platforms which allocate ad slots based on advertiser bids.
These systems typically use a Vickrey (second price) auction or a Vickrey-
Clarke-Groves (VCG) auction***The VCG auction encompasses the traditional
Vickrey auction as a special case. where the winner pays the next highest
bidder’s bid. An attractive property of the VCG auction is that bidding the
true value is a dominant strategy***Note that for multiple slots, platforms
typically rely on generalized second price (GSP) auctions which truth-telling
is not a dominant strategy (Edelman et al., 2007; Varian, 2007). For our
exposition, we assume a single slot with allocation done using a second price
auction. for all bidders (Edelman et al., 2007; Varian and Harris, 2014),
which simplifies bidder’s decisions. That is, each participant is incentivised
to bid what they value the item to be worth, while the winner gets a discount.
This has led to widespread adoption of the second price auction in large scale
online ad marketplaces (Varian and Harris, 2014).

In addition to bid values, auction based recommender systems also typically
account for “conversion rates” to allocate slots to recommendations. However,
both these metrics only account for the immediate value of showing a
recommendation in a slot. In this way, auction based recommender systems are
designed to optimize myopically for short-term engagement metrics and may not
reflect the impact on long-run user engagement. Thus, optimizing for the long-
term in an auction based system poses a unique challenge due to the auction
mechanism. In this work, we propose a reinforcement learning (RL) approach to
bias an auction based recommender system towards strategies that account for
downstream impact of recommendations on user behavior.

Reinforcement learning provides a mathematical formalism to optimize for long
horizon outcomes and has gained traction through superior performance in
different applications, including arcade games (Mnih et al., 2015), robotics
(Gu et al., 2017; Smith and Linden, 2017) and navigation (Tai et al., 2017).
Many recent works have also shown potential of using RL methods for real-time
recommender systems, including off-policy methods like Q-learning (Zou et al.,
2019; Zheng et al., 2018) and actor-critic algorithms (Chen et al., 2019a,
2022). However, balancing between bidder’s value and long-term user engagement
is one of the key constraints for system designers in applying RL to an
auction based recommender system. Unlike most of the prior work, we cannot
design a system from scratch using RL; instead we must work within a framework
which continues to optimize for bidder’s value as well, at least partially.
Moreover, using off-policy RL algorithms to search for an optimal policy poses
significant challenges due to distribution shift. Offline collected data
rarely satisfies coverage conditions required for off-policy RL algorithms to
find a (near) optimal policy. In this work, we make the following
contributions:

  * •

We use a reinforcement learning approach to bias an auction based recommender
system toward recommendations that improve metrics for long-term user
engagement. Our approach can be understood as one step policy improvement over
the auction based recommender policy.

  * •

Using ideas from classic dynamic programming theory, we make a simple argument
to show that our method provably improves over the base (auction) policy.

  * •

We implement our method in an industrial scale real-time recommender system
serving billions of users daily and empirically show performance improvements
through online A/B testing.

This paper is organized as follows. In section 2, we briefly review prior work
on recommender systems with a focus on reinforcement learning based
approaches. Section 3 instantiates a reinforcement learning setup to optimize
for long-term value in a recommender system problem setting we consider. In
section 4, we describe our online reinforcement learning approach which
provably improves over the baseline auction based recommender systems in terms
of long-run user engagement metrics (we use “conversions” as a proxy metric of
long-run user engagement). We end the paper in section 5 with experimental
results from an online A/B test showing significant improvement in conversion
statistics over an six week period.

######  Remark 1.

Our approach is not specifically tied to the surrogate metrics we use to
quantify long-term user engagement. The surrogate metrics should be looked at
through the lens of reward design.

##  2 Related Work

Many of the existing approaches to personalized recommendations, like content
and collaborative filtering based methods (Lops et al., 2011; Ekstrand et al.,
2011; Shi et al., 2014), sequential recommender systems (Hidasi et al., 2015;
Quadrana et al., 2018; Wang et al., 2019) and deep supervised learning based
methods (Hidasi et al., 2015; Zhang et al., 2019), optimize metrics of
immediate user engagement (e.g. click probabilities) without incorporating the
downstream effect of recommendation policies on metrics of long-run user
engagement. This includes (contextual) multi-arm bandit approaches (Li et al.,
2010), although bandit approaches allow for learning user preferences by
intelligent exploration thereby combating the feedback loops from data
generated by existing recommender systems (Wang et al., 2017; Qin et al.,
2014; Zeng et al., 2016; Zhu and Van Roy, 2023b; Guo et al., 2023).

Given this and the sequential nature of user interactions, reinforcement
learning techniques have emerged as natural alternatives to grapple with the
problem of myopic recommendation strategies more generally (Shani et al.,
2005; Chen et al., 2019a; Ie et al., 2019; Xin et al., 2020; Zhao et al.,
2019). Since RL methods aim to optimize for long-run “value” (estimate of
accrued reward over a long time horizon), the resulting recommendation
strategies account for the efficacy of downstream recommendations in recurring
user interactions. Although many different approaches have been proposed,
including model-based RL methods (Chen et al., 2019b), a lot of recent work
has focused on applying model free RL algorithms to maximize long term
measures of user engagement in different applications; see (Zou et al., 2019;
Zheng et al., 2018; Zhu and Van Roy, 2023a) for examples in e-commerce setting
and a personalizing news article recommendation setting. Similar to these
efforts, we train a deep neural network model to approximate a Q-function of
an auction based recommender policy. However, our focus is on incorporating
long-term metrics in an auction-based recommender system. For this, we
restrict ourselves to biasing the recommendation policy toward actions with
higher Q-values, as opposed to searching for an optimal policy. We use simple
ideas from classical dynamic programming theory to elucidate the theoretical
motivation of our one-step policy improvement approach – modifying the base
policy in this way leads to a policy with improved long-run user engagement
metrics.

##  3 Problem Formulation

In this section, we formalize the long-term value optimization problem of a
recommender system in the framework of reinforcement learning. We later
describe many approximations we make for practical implementation of a large
scale auction based system that is designed around this framework and discuss
the limitations of our choices. A general reinforcement learning framework can
be characterized by an agent and its interactions with the environment Lu et
al. (2021); Hutter (2007)

###  3.1 A model of the environment.

Following Lu et al. (2021), we model the environment by a tuple
$\mathcal{E}=(\mathcal{O},\mathcal{A},\rho)$, where $\mathcal{O}$ and
$\mathcal{A}$ denote the observation set and action set respectively and
$\rho$ prescribes a probability for any observation $o\in\mathcal{O}$. For
simplicity, we assume that time is discrete and take one time period to denote
a 24 hour window for the model we implement. We assume that the set of users
interacting with the recommender system is dynamic and denote $\mathcal{N}(t)$
to be the set of interacting users at time $t$. To avoid overload of notation,
we take the set of all users to be fixed, denoted by $\mathcal{N}$; our
framework as well as implementation can easily account for new users.

  1. 1.

Action space: We take $\mathcal{A}$ to be a fixed***A fixed action space is
assumed to only simplify the exposition in terms of notation. finite set of
all possible items that can be recommended to any user at a given time and let
$a\in\mathcal{A}$ denote one such recommendation. At time step $t$, an action
set $A_{t}\in\mathcal{A}^{|\mathcal{N}(t)|}$ is a collection of recommended
items to different users, $A_{t}=\left(a^{u}_{t}:u\in\mathcal{N}(t)\right)$.

  2. 2.

Observation space: We take $\mathcal{O}$ to be a finite set of observations
about a user, with each observation $o\in\mathcal{O}$ encoding useful
information about user tastes. The observation set,
$O_{t}\in\mathcal{O}^{|\mathcal{N}(t)|}$ is a collection of such observations
about each user interacting with the system,
$O_{t}=(o^{u}_{t}:u\in\mathcal{N}(t))$.

  3. 3.

Outcome space: We take $\mathcal{Y}$ to be a finite set of outcomes, with each
outcome $y\in\mathcal{Y}$ indicating a user’s response to a recommendation.
Similar to the action and observation sets, the outcome set,
$Y_{t}=(y^{u}_{t}:u\in\mathcal{N}(t))\in\mathcal{Y}^{|\mathcal{N}(t)|}$
denotes the collection of outcomes for each user interacting with the system
at time $t$. We later describe how $Y_{t}$ depends on various quantities,
including the action and observation sets, $A_{t}$ and $O_{t}$.

  4. 4.

Observation probability: We let history
$H_{t}=(O_{0},A_{0},Y_{0},\ldots,A_{t-1},Y_{t-1},O_{t})$ denote interactions
with the environment up to time $t$. An agent selects actions $A_{t}$ given
the history $H_{t}$. The probability distribution $\rho(\cdot|H_{t},A_{t})$
determines the next observation conditioned on $H_{t}$ and $A_{t}$. While we
leave the exact form of $\rho$ unspecified, $\rho(\cdot)$ essentially models
the evolution of user tastes as well as their propensity to interact with the
system.

###  3.2 A model of user behavior

Given the environment interface as described above, we model the user response
$y^{u}_{t}$ to a recommended item $a^{u}_{t}$ at time $t$ as a function of its
past interactions with the system***History of a user can be written as
$h^{u}_{t}=\\{(o^{u}_{k},a^{u}_{k}):k\in\mathcal{T}(u)\\}$ where
$\mathcal{T}(u)$ is the set of interaction times of user u with the system up
to time $t$. This can be derived from $H_{t}$ which logs the collection of
observation, action sets up to time $t$. $H_{t}$, the recommended item
$a^{u}_{t}$, user context $x^{u}_{t}$ which encodes demographic information
about the user, as well as some side information $i^{u}_{t}$ which encodes
dynamic user interests summarized by their interaction with the platform,
beyond engagement with the recommender system***For example, this side
information can be derived from organic posts on the platform that users
engage with.. Formally, one can take $g(\cdot)$ to be some unknown fixed
function such that,
$y^{u}_{t}=g(H_{t},a^{u}_{t},x^{u}_{t},i^{u}_{t},\epsilon_{t})$ where
$\epsilon_{t}$ is i.i.d noise capturing idiosyncratic randomness in user
behavior. Throughout we assume that the user contexts as well as side
information $(x_{u},i_{u})$ for each user $u\in\mathcal{N}$ evolves
independently of other users with time. Essentially, this translates to
assuming that there are no network effects.

#### User state representation:

Given our model of the user behavior, we take a user state to include all
information about a user available to the recommender system at the time
period. In particular, we assume vector encoding $z^{u}_{t}=f(H_{t},u)$ to
summarize a user’s interests based on their past interactions with the
recommender system. Taken together with the user context vector $x^{u}_{t}$
and the side information encoding $i^{u}_{t}$, the user state can be
parsimoniously represented as the tuple
$s^{u}_{t}=(z^{u}_{t},x^{u}_{t},i^{u}_{t})$. It is noteworthy that assuming no
network effects, the user states are independent of each other, i.e.
$s^{u_{1}}\perp s^{u_{2}}$ for any pair of users $(u_{1},u_{2})$. Put
differently, a user’s state only depends on its past interactions with the
recommender system - their context and interests evolve exogenously.
Throughout the paper, we denote $s$ to be any generic user state and
$\mathcal{S}$ to denote the state space (collection of all users states). For
simplicity of exposition, we will assume $\mathcal{S}$ to be countable.

###  3.3 Agent design.

Besides the environment, a key component of an RL framework is an agent which
takes actions given observations from the environment. We let
$\pi(\cdot|H_{t})$ denote a policy, a probability distribution over actions
that depends on the history $H_{t}$.

  1. 1.

Agent state: From a computational viewpoint, the dependence of a policy on the
entire history is problematic. We therefore take $S_{t}=\\{s^{u}_{t}:u\in
N(t)\\}$ to represent agent “state”, as a collection of state representations
for all users interacting with the system at time $t$. Note that the agent
state $S_{t}$ encodes all relevant information in history $H_{t}$, including
the observation $O_{t}$. We can write the relationship between observation,
states and history as:

$H_{t}\to S_{t}\to O_{t+1}\to S_{t+1}\,\,\text{or equivalently},\,\,H_{t-1}\to
S_{t-1}\to O_{t}\to S_{t}\to O_{t+1}\to S_{t+1}.$

We use the notation $X\to Y\to Z$, which is standard in information theory, to
denote that random variables $X$ and $Z$ are independent conditioned on Y.
Clearly, the state variable obeys the Markov property,

$\mathbb{P}(S_{t+1}=\cdot\,|\,S_{t},A_{t},\cdots,S_{0},A_{0})=\mathbb{P}(S_{t+1}=\cdot\,|\,S_{t},A_{t}),$

since conditioned on $S_{t}$ and $A_{t}$, $S_{t+1}$ is independent of the
history $H_{t}$ and no network effects are assumed.

  2. 2.

Reward function: We denote $r:\mathcal{A}\times\mathcal{Y}\mapsto\\{0,1\\}$ to
be a reward function which maps an (action, outcome) pair to a scalar,
encoding the agent’s preferences about outcomes. In our setting, rewards are
taken to be binary and $r(a,y)=1$ indicates “conversion”. Essentially, we
associate an “end behavior” for each recommended item $a$. For example, buying
might be the end behavior for some recommendations while getting users to
subscribe might be the end behavior for others. User behavior (as indicated by
the outcome variable) corresponding to the end behavior for an ad results in a
reward of 1; all other user behaviors incur a reward of 0. For action,
observation set $(A_{t},Y_{t})$, we assume rewards of individual users are
additive and take
$R(A_{t},Y_{t})=\sum_{u\in\mathcal{N}(t)}r(a^{u}_{t},y^{u}_{t})$ to be the
scalar reward.

  3. 3.

Constraint function: It is common for recommender systems to have constraints
on which items can be shown to a user at any given time. We denote
$c:\mathcal{A}\times\mathcal{S}\mapsto\bar{\mathcal{A}}$ with
$\bar{\mathcal{A}}\subseteq\mathcal{A}$ to be such a constraint function which
restricts the feasible action set to a subset of all possible actions,
augmented with an option of “no recommendation”.

###  3.4 Reinforcement learning of recommendation policies

While the reward function specifies the agent’s preferences over a single time
step, the goal in reinforcement learning is to optimize for cumulative rewards
over long horizons. As is typical, we take a discounting approach to model
cumulative rewards with the goal to maximize expected cumulative discounted
return,

$\mathbb{E}\Big{[}\sum_{t=0}^{\infty}\gamma^{t}R(A_{t},Y_{t})\Big{]}$ (1)

Thus, a good recommendation policy not only optimizes for immediate rewards
but also for long run conversions by transitioning users to states with a
higher propensity to convert.

As the agent state encodes all relevant information about the history, we
denote the agent policy $\pi_{\text{agent}}(\cdot|S)$ to be a distribution
over actions given the agent state. Furthermore, our assumption on no network
effects (independence of user states across time) and additive reward
functions lets us decompose the policy as

$\pi_{\text{agent}}(\cdot\,|\,S_{t})=\prod_{u\in\mathcal{N}(t)}\pi(\cdot\,|\,s^{u}_{t})$

where $\pi(\cdot|s)$ indicates the action selection policy given any generic
user state $s$. To compare different policies as well as to find the optimal
one, reinforcement learning algorithms use value functions. Formally, the
value function for policy $\pi$ can be defined as

$V_{\pi}(s)=\mathbb{E}_{\pi}\Big{[}\sum_{t=0}^{\infty}\gamma^{t}r(a_{t},y_{t})\,|\,s_{0}=s\Big{]}$
(2)

which measures the cumulative reward under policy $\pi$, from a given user
state $s$. Similarly, we can also define the state-action value function as
the expected return under policy $\pi$ from a user state $s$ and
recommendation $a$,

$Q_{\pi}(s,a)=\mathbb{E}_{\pi}\Big{[}\sum_{t=0}^{\infty}\gamma^{t}r(a_{t},y_{t})\,|\,a_{0}=a,s_{0}=s\Big{]}$
(3)

State-action value functions obey a fixed point equation known as the Bellman
equation, and can equivalently be expressed as,

$Q_{\pi}(s,a)=r(a,y)+\gamma\sum_{s^{\prime}\in\mathcal{S}}P(s^{\prime}|s,a)V_{\pi}(s^{\prime})$
(4)

where $\mathcal{S}$ denotes the set of all states and
$\\{s^{\prime}\\}\in\mathcal{S}$ are the set of successor states from the
state-action pair $(s,a)$. An optimal policy $\pi^{*}$, is defined to be the
one which maximizes the value-to-go from each state,
$\pi^{*}(s)=\operatorname*{arg\,max}_{\pi}V_{\pi}(s)$. The ultimate goal in RL
is to find the optimal policy which maximizes the expected cumulative return
from any user state. However, we take a less ambitious approach of performing
a one-step policy improvement over a given base policy. We outline this below.

##  4 A policy improvement approach

One of our motivation to carefully understand performance of a policy obtained
by a single policy improvement step is tied to business constraints. Auction-
based recommender systems are the main revenue generating work-streams of
advertisement platforms and therefore it is reasonable to not deviate too much
from the given base policy without exhaustive experimentation. Moreover, we
are also constrained by working with real user interaction data - we do not
have access to a simulator of user behavior that we can repeatedly query to
assess the performance of a policy trained with reinforcement learning.
Essentially, from our logged data
$\\{(O_{t},A_{t},Y_{t})\\}_{t\in\mathbb{N}}$, we construct trajectories for
each user,

$\mathcal{D}_{\pi_{\text{base}}}=\\{(o^{u}_{t},s^{u}_{t},a^{u}_{t},r^{u}_{t})\\}_{u\in\mathcal{N},\,t\in\mathcal{T}(u)}$

where $\mathcal{T}(u)$ denotes the set of interaction times of user $u$ with
the system, accounting for the fact that not all users choose to interact with
the system at all times. Here, the subscript $\pi_{\text{base}}$ emphasizes
that this interaction data is collected using the (current) base policy. We
give details about the base policy in the subsection below.

A standard policy iteration step involves approximating the Q-function of the
base policy, $Q_{\pi_{\text{base}}}$ using a RL algorithm and subsequently
solving the single period optimization problem,

$\pi^{+}(s)=\operatorname*{arg\,max}_{a\in\mathcal{A}}\,Q_{\pi_{\text{base}}}(s,a)$
(5)

Notice from (4) that for any policy $\pi$, $Q_{\pi}$ is parameterized by the
immediate action which optimizes for the one step reward $r(a,y)$ and
transition to the next state $s^{\prime}$. Therefore, the optimization problem
in (5) is sometimes referred to as a “single period” optimization problem. A
well known result in dynamic programming shows that a policy $\pi^{+}$ which
simultaneously solves the single step policy iteration problem for all states
is an improved policy (Bertsekas, 2012). That is,

$V_{\pi^{+}}(s)\geq V_{\pi_{\text{base}}}(s)\,\,\,\forall\,s\in\mathcal{S}$

For a problem with finite state and action spaces, a series of policy
improvement steps ultimately lead to an optimal policy geometrically fast
(Bertsekas, 2012). While we can technically obtain an improved policy over the
base recommender policy using (5), we are bound by business constraints of
staying close to the base policy. Nevertheless, we show that it is possible to
find an improved policy by taking a step in the policy improvement direction.

###  4.1 Baseline recommender policy with auction mechanism

The base policy, $\pi_{\text{base}}$ that we consider is a second price
auction based policy. Essentially, for a given user state $s$, bidders provide
a “bid” value $\text{Bid}(s,a)$ for a subset of eligible recommendations
$a\in\mathcal{A}^{b}\subset\bar{\mathcal{A}}$ which represents an estimate of
their utility of recommending item $a$ in user state $s$. As is common in many
practical instantiations of recommender systems, an expected conversion rate
(“eCVR”) factor is multiplied to the bid values to upper bound the expected
revenue of showing an ad in user state $s$. We take
$f(s,a)=\text{Bid}(s,a)\times\text{eCVR}(s,a)$ to denote the “bid-eCVR”
values. Without going too much into details, one can take the deterministic
function $f(\cdot)$ to be a black-box recommendation scoring model. Note
however that this recommendation scoring model is not trained for optimizing
conversions over long time horizons; rather it is trained to be a good
predictor of short run engagement of the users with different items (it takes
immediate conversion rates as input). Another common practice in auction-based
recommender systems is that each bidder typically only bids on a unique subset
of recommendations associated with campaigns. Therefore, it is reasonable to
assume that $\mathcal{A}^{b_{1}}\cap\mathcal{A}^{b_{2}}=\emptyset$ for any two
bidders $b_{1}$ and $b_{2}$.

In this work, we take the agent’s base policy $\pi_{\text{base}}$ to be greedy
with respect to the bid-eCVR values,

$\displaystyle\pi_{\text{base}}(s)=\operatorname*{arg\,max}_{a\in\bar{\mathcal{A}}}\,f(s,a)$
(6)

where $\bar{\mathcal{A}}$ is the set of all eligible recommendations. Note
that our description of the base policy in (6) is not arbitrary; it is a good
approximation to many practical auction based recommender systems. The auction
mechanism is critical to recommender systems and is commonly used across
various monetization products and online marketplaces. Moreover, platforms
invest a lot of engineering resources to build, maintain and fine-tune
conversion rate models. This limits us from re-inventing the wheel from
scratch. Instead, we treat $f(\cdot)$ as a fixed black-box input to the
baseline policy and aim to find an improved policy (in terms of cumulative
“conversions” as defined above) by perturbing $f(\cdot)$ with an estimate of
the expected long run conversions of recommending item $a$ in user state $s$,
$Q_{\pi_{\text{base}}}(s,a)$.

###  4.2 Optimizing for long-run conversions: a one-step improved policy

Our proposed solution is to bias the recommendation scoring model $f$ to
account for long run conversions. Specifically, let $Q_{\pi_{\text{base}}}$
denote the state-action value function of the base policy, with
$Q_{\pi_{\text{base}}}(s,a)$ denoting the expected cumulative conversions
attributable to recommending item $a$ in user state $s$ and following the base
policy thereafter. We take a modified policy to recommend items based on a
weighted score, which is a convex combination of the recommendation score and
estimated long run conversions,

$\pi_{\text{mod}}(s)=\text{SELECT}_{a\in\bar{\mathcal{A}}}\,\left(f(s,a),\hat{Q}_{\pi_{\text{base}}}(s,a)\right)$
(7)

where we define the selection operator as

$\displaystyle\text{SELECT}_{a\in\bar{\mathcal{A}}}\,$
$\displaystyle\left(f(s,a),\hat{Q}_{\pi_{\text{base}}}(s,a)\right):=\operatorname*{arg\,max}_{a\in\bar{\mathcal{A}}}\,\left[(1-\alpha)\cdot
f(s,a)+\alpha\cdot\hat{Q}_{\pi_{\text{base}}}(s,a)\right]$

for a fixed value of $\alpha$. See figure 1 for an illustration of our
approach. Note that when $\alpha=1$, $\bar{\mathcal{A}}=\mathcal{A}$ and exact
Q values are available,

$\pi_{\text{mod}}(s)=\pi^{+}(s)=\operatorname*{arg\,max}_{a\in\mathcal{A}}\,Q_{\pi_{\text{base}}}(s,a)$

and the modified policy $\pi_{\text{mod}}$ mimics the policy iteration step.
When $\alpha=0$, $\pi_{\text{mod}}$ equals the base policy which is greedy
with respect to the recommendation scoring model $f$. Values of
$\alpha\in(0,1)$ interpolate between the two policies, possibly in a non-
linear way and fixing a value of $\alpha$ pins down the trade-off between
optimizing short run engagement versus optimizing for long run conversions. In
practice, we use a value of $\alpha=0.96$ which was chosen empirically after a
few A/B test iterations, taking into account the business constraint of
staying close to the base policy. We discuss more about the choice of $\alpha$
in section 5.

Figure 1: Visualization of an RL approach to improve auction-based recommender
policies. We optimize a weighted combination of recommender score $f(s,a)$,
which is used for item selection in a second price auction, and estimated long
run conversions $\hat{Q}_{\pi_{\text{base}}}(s,a)$.

Next, we show how the modified policy in (7) with exact Q values yields an
improved policy, $V_{\pi_{\text{mod}}}(s)\geq V_{\pi_{\text{base}}}(s)$ for
all user states $s$. That is, the expected cumulative conversions under the
modified policy from each user state improves as compared to the base policy,
with possible strict improvements from a subset of user states. This result is
not unexpected as the base policy is not optimized for long run conversions.
In comparison, the modified policy takes into account expected cumulative
conversions over long horizons, at least partially. Nevertheless, this result
acts as a sanity check and serves as a motivation for our approach. Recall the
definition of value function from (2) and the functional form of the modified
policy given in (7). For our results, we assume that the function $f$ takes
unique values for different actions in each user state, i.e. $f(s,a_{i})\neq
f(s,a_{j})$ for any $(a_{i},a_{j})\in\bar{\mathcal{A}}$. This is reasonable as
we expect bid-eCVR values for different recommendations to be non-unique. In
addition, we also assume that we have access to exact Q values of the base
policy, $Q_{\pi_{\text{base}}}$.

######  Proposition 1 (Policy improvement step).

For any user state $s$, $V_{\pi_{\text{mod}}}(s)\geq
V_{\pi_{\text{base}}}(s)$. Moreover,
$V_{\pi_{\text{mod}}}(s)>V_{\pi_{\text{base}}}(s)$ if
$\pi_{\text{mod}}(s)\neq\pi_{\text{base}}(s)$. That is, the inequality is
strict if recommendations under the modified and base policies differ in a
user state.

The proof follows by a straightforward application of the monotonicity
property of the Bellman operators (Bertsekas, 2012) along with the following
lemma. Note that $Q_{\pi_{\text{base}}}(s,\pi_{\text{mod}}(s))$ denotes the
expected cumulative conversions using the modified policy in user state $s$
and following the base policy thereafter.

######  Lemma 1.

For any user state $s$, $Q_{\pi_{\text{base}}}(s,\pi_{\text{mod}}(s))\geq
V_{\pi_{\text{base}}}(s)$.

###### Proof.

Consider a fixed state $s$ and let $a_{\text{base}}$, $a_{\text{mod}}$ denote
the selected actions under the base and modified policy, with the same
eligible***While comparing action selections under the base and modified
policy, we assume the eligible action set remains the same. This seems
reasonable to compare the outcome of two policies. Note that the constraint
function (which specifies the eligible item set) is part of agent design and
is separate from learning of recommendation policies. item set
$\bar{\mathcal{A}}$

$\displaystyle a_{\text{base}}$
$\displaystyle=\pi_{\text{base}}(s)=\operatorname*{arg\,max}_{a\in\bar{\mathcal{A}}}\,f(s,a)$
(8) $\displaystyle a_{\text{mod}}$
$\displaystyle=\pi_{\text{mod}}(s)=\operatorname*{arg\,max}_{a\in\bar{\mathcal{A}}}\,\left[(1-\alpha)f(s,a)+\alpha
Q_{\pi_{\text{base}}}(s,a)\right]$ (9)

Then, we show that, $Q_{\pi_{\text{base}}}(s,a_{\text{mod}})\geq
Q_{\pi_{\text{base}}}(s,a_{\text{base}})$. Clearly, this holds when
$a_{\text{mod}}=a_{\text{base}}$. Let’s look at the case when the two actions
are different. In this case, first note that
$f(s,a_{\text{base}})>f(s,a_{\text{mod}})$ using the fact that
$a_{\text{base}}=\operatorname*{arg\,max}_{a\in\bar{\mathcal{A}}}\,f(s,a)$ and
that $f(\cdot)$ takes on unique values. However, from (9), it holds that

$\displaystyle(1-\alpha)f(s,a_{\text{mod}})+\alpha Q_{\pi_{\text{base}}}$
$\displaystyle(s,a_{\text{mod}})\geq(1-\alpha)f(s,a_{\text{base}})+\alpha
Q_{\pi_{\text{base}}}(s,a_{\text{base}}).$

Rearranging, we get that

$\displaystyle\alpha(Q_{\pi_{\text{base}}}(s,a_{\text{mod}})-Q_{\pi_{\text{base}}}$
$\displaystyle(s,a_{\text{base}}))\geq(1-\alpha)(f(s,a_{\text{base}})-f(s,a_{\text{mod}}))>0,$

where the last inequality follows as
$f(s,a_{\text{base}})>f(s,a_{\text{mod}})$. Thus, for the case when
$a_{\text{mod}}\neq a_{\text{base}}$, we have

$Q_{\pi_{\text{base}}}(s,a_{\text{mod}})>Q_{\pi_{\text{base}}}(s,a_{\text{base}}).$
(10)

Combining both the cases, we conclude that
$Q_{\pi_{\text{base}}}(s,a_{\text{mod}})\geq
Q_{\pi_{\text{base}}}(s,a_{\text{base}})$. Note that as both
$\pi_{\text{base}}$ and $\pi_{\text{mod}}$ are deterministic policies, we get
our final result,

$\displaystyle
Q_{\pi_{\text{base}}}(s,\pi_{\text{mod}}(s))=Q_{\pi_{\text{base}}}(s,a_{\text{mod}})\geq
Q_{\pi_{\text{base}}}(s,a_{\text{base}})=Q_{\pi_{\text{base}}}(s,\pi_{\text{base}}(s))=V_{\pi_{\text{base}}}(s).$

Here the final equality uses that $Q_{\pi}(s,\pi(s))=V_{\pi}(s)$. This
relationship is well known in RL (see for example (Bertsekas, 2012)) and can
also be derived from equations (2) and (4) above. Note that equation (10)
implies a strict inequality,
$Q_{\pi_{\text{base}}}(s,\pi_{\text{mod}}(s))>V_{\pi_{\text{base}}}(s)$, when
the modified policy makes a different recommendation than the base policy. ∎

#### Proof of Proposition 1:

We use the result in Lemma 1 to show how the modified policy gives an
improvement over the base policy. This essentially follows from using the
standard monotonicity properties of the Bellman operator (Bertsekas, 2012).
Note that as per step rewards are bounded in our setting, $|r(a,y)|\leq 1$,
value functions corresponding to any stationary policy $\pi$ are also bounded.
Let $\mathcal{V}=\\{V_{\pi}:\pi\in\Pi\\}$ be the set of all bounded value
functions corresponding to stationary policies $\pi\in\Pi$. The bellman
operator $T_{\mu}:\mathcal{V}\mapsto\mathcal{V}$ for some policy $\mu\in\Pi$
maps a value functions $V\in\mathcal{V}$ as

$T_{\mu}V(s):=r(\mu(s),y)+\gamma\sum_{s^{\prime}\in\mathcal{S}}P(s^{\prime}|s,\mu(s))V(s^{\prime})\,\,\,\,\forall\,s\in\mathcal{S}.$

Using this and the definition of Q function in Equation (4), it is easy to
check that for any two stationary policies $\pi,\pi^{\prime}$,

$Q_{\pi}(s,\pi^{\prime}(s))=T_{\pi^{\prime}}V_{\pi}(s)$ (11)

With bounded per-step rewards (our setting), Bellman operators are monotone
(Bertsekas, 2012), implying that if $V_{\pi}(s)\geq V_{\pi^{\prime}}(s)$ for
all $s\in\mathcal{S}$, then

$T_{\mu}V_{\pi}(s)\geq
T_{\mu}V_{\pi^{\prime}}(s)\quad\forall\,s\in\mathcal{S}.$ (12)

Another standard result from dynamic programming theory (Bertsekas, 2012)
helps us prove proposition 1.

$V_{\mu}=\lim_{k\to\infty}\,T^{k}_{\mu}V_{\pi}$ (13)

Essentially, this states that repeated application of the Bellman operator
$T_{\mu}(\cdot)$ converges to its corresponding value function $V_{\mu}$.
Using (11), we can express the result of lemma 1 as

$T_{\pi_{\text{mod}}}V_{\pi_{\text{base}}}(s)\geq
V_{\pi_{\text{base}}}(s)\quad\forall\,s\in\mathcal{S}.$ (14)

Applying the monotonicity property to (14) gives,

$T^{2}_{\pi_{\text{mod}}}V_{\pi_{\text{base}}}(s)\geq
T_{\pi_{\text{mod}}}V_{\pi_{\text{base}}}(s)\geq V_{\pi_{\text{base}}}(s)$

Repeatedly applying the monotonicity property shows the desired result,

$\displaystyle
V_{\pi_{\text{mod}}}(s)=\lim_{k\to\infty}T^{k}_{\pi_{\text{mod}}}V_{\pi_{\text{base}}}(s)\geq\ldots\geq
T^{2}_{\pi_{\text{mod}}}V_{\pi_{\text{base}}}(s)\geq
T_{\pi_{\text{mod}}}V_{\pi_{\text{base}}}(s)\geq V_{\pi_{\text{base}}}(s).$

It is noteworthy that the inequality is strict for user states where the
recommendation actions under $\pi_{\text{mod}}$ differ from those under
$\pi_{\text{base}}$. Indeed, proof of lemma 1 shows that the inequality in
(14) is strict for these states.

###  4.3 Implementation using SARSA

To estimate $Q_{\pi_{\text{base}}}$, we train a deep neural network based
model using SARSA, a popular on-policy reinforcement learning algorithm, and
logged data of user interactions
$\\{(s^{u}_{t},a^{u}_{t},r^{u}_{t})\\}_{u\in\mathcal{N},\,t\in\mathcal{T}(u)}$.
We pre-process the trajectory data of each user into tuples,

$d^{u}_{t_{i}}=\left(s^{u}_{t_{i}},a^{u}_{t_{i}},r^{u}_{t_{i}},s^{u}_{t_{i+1}},a^{u}_{t_{i+1}},\tau_{i}\right)\,\,\,\forall\,t_{i}\in\mathcal{T}(u)=\left\\{t_{0},t_{1},\ldots,t_{|\mathcal{T}(u)|}\right\\}$

Recall that $\mathcal{T}(u)$ denote the set of interaction times of a user $u$
with the recommender system. We let $\tau_{i}=t_{i+1}-t_{i}$ denote the
$i^{\text{th}}$ inter-interaction time and append it to each tuple.
Reinforcement learning algorithms like SARSA which do credit assignment by
bootstrapping only require training data in the form of these data tuples,
unlike Monte Carlo methods which require entire user interaction trajectories.
SARSA is a widely used standard reinforcement learning algorithm – however,
for completeness we give some details in Algorithm 1 below.

Let
$\mathcal{D}=\\{d^{u}_{t_{i}}\\}_{u\in\mathcal{N},\,t_{i}\in\mathcal{T}(u)}$
be the collection of all user interaction tuples. In reinforcement learning,
this is usually referred to as the replay buffer. At each training epoch, data
tuples are sampled from $\mathcal{D}$ to train the SARSA algorithm. Note that
the user interaction data is collected under the base policy,
$\pi_{\text{base}}$ and so we can hope to estimate $Q_{\pi_{\text{base}}}$
using an on-policy algorithm like SARSA. We also train our models online,
instead of using a fixed offline data set for training, by updating data set
$\mathcal{D}$ daily with a fresh set of user interaction data. See section 4.4
for more details on how we construct these data tuples from user trajectories
in an online manner. User states which do not have successor interactions are
taken to be terminal and we attribute the immediate reward (conversions) to
the corresponding (state, action) pair. Also, as is typical in the training of
reinforcement learning algorithms with deep neural networks, we use a separate
target network for stability (Mnih et al., 2015). We update the target network
every $k=100$ training steps of the algorithm.

Input: Neural network model $Q_{\theta}$, target network $Q_{\theta^{-}}$,
initial replay buffer $\mathcal{D}$, batch-size $N$, step-size sequence
$\\{\alpha_{k}\\}_{k\in\mathbb{N}}$.;

Initialize: $\theta,\theta^{-}\leftarrow\theta_{0}$.

for  _$t=0,1,\ldots$_ do

Update replay buffer $\mathcal{D}_{t}$ of user transition tuples.

Sample batch
$\mathcal{B}=\\{(s_{j},a_{j},r_{j},s^{\prime}_{j},a^{\prime}_{j},\tau_{j})\\}_{j=1}^{N}$
from $\mathcal{D}_{t}$.;

Compute targets for each tuple:  $y_{j}=\begin{cases}r_{j},&\text{if $s_{j}$
is a terminal user state}\\\
r_{j}+\gamma^{\tau_{j}}Q_{\theta^{-}}(s^{\prime}_{j},a^{\prime}_{j}),&\text{otherwise}\end{cases}$
Compute loss function:
$\ell(\theta_{t})=\sum_{j=1}^{N}(y_{j}-Q_{\theta}(s_{j},a_{j}))^{2}$;

Take a semi-gradient step:
$\theta_{t+1}\leftarrow\theta_{t}-\alpha_{t}\nabla_{\theta}\ell(\theta_{t})$;

if  _$t\,\,\text{mod}\,\,k==0$ :_ then

Update target network: $\theta^{-}\leftarrow\theta_{t}$;

end if

end for

Output: estimate of state action value function $Q_{\theta}$

Algorithm 1 Q-value estimation using SARSA

###  4.4 Data processing

Given that recommender systems serves millions of users every day, storing
interaction data for each user over long time horizons is impossible, even
with distributed storage. Therefore, we take an online training approach –
collecting fresh data every day, integrating it with a fixed length table of
recent interaction data for each user, and using it to generate training
tuples. To describe our data processing pipeline, we introduce some notation
specific to this section. We let $h$ be the effective time horizon, specifying
the number of days we look back to maintain the table of recent user
interactions. Fix a time period $t>h$. Let
$B_{t}=\left(S_{t-h},A_{t-h},R_{t-h},S_{t-h+1},\ldots,S_{t},A_{t},R_{t}\right)$
be a “buffer” table which stores the last $h$ days of interaction data in
terms of the states, recommended items and rewards. Let $\mathcal{N}_{h}(t)$
be the set of users with interaction data in $B_{t}$. We process the agent-
environment interactions at time $t+1$,
$\left(S_{t+1},A_{t+1},R_{t+1}\right)$, into data tuples by dividing user
interactions in three categories:

  1. 1.

New user interactions: Consider a user $u$ which interacts with the system at
time $t+1$ but has not interacted with the system in the last $h$ time
periods, i.e. $u\notin\mathcal{N}_{h}(t)$. Since no predecessor or successor
interactions are associated with this user, we do not use this data to
generate a transition tuple. Information about this interaction remains a part
of the agent-environment interaction, $(S_{t+1},A_{t+1},R_{t+1})$, and is
added to the buffer $B_{t+1}$.

  2. 2.

Active user interactions: Consider a user $u$ which interacts with the system
at time $t+1$, and this user has interacted with the system in the recent
past, i.e. $u\in\mathcal{N}_{h}(t)$. Let $t_{u}\in\\{t-h,t-h+1,\ldots,t\\}$ be
the most recent interaction time of user $u$ before $t+1$. Then, we take
$d^{u}_{t+1}=(s^{u}_{t_{u}},a^{u}_{t_{u}},r^{u}_{t_{u}},s^{u}_{t+1},a^{u}_{t+1})$
to be a transition data tuple and add it to the replay buffer $\mathcal{D}$.

  3. 3.

Inactive user interactions: Consider a user $u$ which interacted with the
system at time $t-h$ but has not interacted with the system since then,
including at time $t+1$. We identify such users to be inactive, and take
$(s^{u}_{t-h},a^{u}_{t-h})$ to be the terminal state action pair such that the
reward $r^{u}_{t-h}$ is attributed to them. That is, the transition tuple
$d^{u}_{t-h}=(s^{u}_{t-h},a^{u}_{t-h},r^{u}_{t-h},s^{*},a^{*},\tau=0)$ is
added to the replay buffer $\mathcal{D}$, where $(s^{*},a^{*})$ are dummy
states and actions such that $Q_{\pi_{\text{base}}}(s^{*},a^{*})=0$.

In some way $h$ can be interpreted as the effective time horizon in our
implementation. As we work in the discounted setting and take a value of
$h=15$ in our experiments, this implicitly implies***This follows as
$(0.8)^{15}\approx 0.03$. We take $\gamma^{15}=0$. a discount factor of
$\gamma\approx 0.8$. Doing data processing this way with a fixed size window
is only an approximation to account for very long term effects of
recommendations, given data size and memory requirements of a real-time
recommender system. We are currently experimenting with a different data
processing pipeline architecture which can allow us to take values of $h$ up
to 60.

#### Updating replay buffer and buffer table:

The new data tuples generated from agent-environment interaction
$(S_{t+1},A_{t+1},R_{t+1})$ are added to the replay buffer $\mathcal{D}$. We
also update the buffer table
$B_{t+1}=(S_{t+1-h},A_{t+1-h},R_{t+1-h},\ldots,S_{t+1},A_{t+1},R_{t+1})$ to
store the agent environment interactions of past $h$ periods at time $t+1$.
This constitutes the “update replay buffer” step in Algorithm 1. See Figure 2
below for an illustration of different types of user interaction tuples
generated by the agent-environment interaction $(S_{t+1},A_{t+1},R_{t+1})$
which are used to update the buffer table $B_{t+1}$.

Figure 2: Data processing pipeline with different user interactions. First two
rows show how new data interaction tuples are created and updated in the
buffer table. The third row illustrates inactive user interaction, where a
dummy state-action is appended to create a data tuple. Fourth row shows that
recent interactions are carried over in the updated buffer table $B_{t+1}$ for
users who did not interact with the system at time $t+1$. The final row
illustrates a new user interaction.

##  5 Experiments

(a)

(b)

Figure 3: Weekly change in average number of conversions and conversion rate
over six weeks of experimentation. Note that change in raw conversion numbers
closely mimics the change in conversion rates as we tried to maintain similar
total impressions over the test and control group. Error bars represent the
standard deviation of change across different bidders.

(c)

Figure 4: Daily change in conversion and impression numbers over the control
group. For most days, our approach (blue curve) results in a significant
increase in conversions over the auction based recommender policy. The red
curve shows the change in number of daily impressions, which approximately
stayed close to 0 for most days.

We conduct an online A/B experiment over a period of six weeks to test our
proposed RL based approach in a web scale recommender system serving billions
of users every day. The base policy in our experiment uses a second price
auction mechanism to recommend items to users. On average, about 2 billion
users interact with our recommender system every day. For our experiments, we
randomly selected 2% (40 million users) to be a part of our test and control
group each. Recall from section 4.4 that our estimated Q values account for
cumulative conversions over a 14 day period.

To stay close to the base policy, we tune the hyperparameter $\alpha$
introduced in section 4.2 to make sure that the contribution of estimated
cumulative conversions, $\alpha Q_{\pi_{\text{base}}}(s,a)$, in the scoring
model $(1-\alpha)f(s,a)+\alpha Q_{\pi_{\text{base}}}(s,a)$ is less than 8%.
That is, we make sure that the fraction $\alpha
Q_{\pi_{\text{base}}}(s,a)/\left((1-\alpha)f(s,a)+\alpha
Q_{\pi_{\text{base}}}(s,a)\right)\leq 0.08$. After hyperparameter tuning or
our experiments, this approximately came to a value of $\alpha=0.96$.

#### Implementation details:

As described in section 3, we model the user state as a tuple of three
features, $s^{u}=(z^{u},x^{u},i^{u})$. The context features $x^{u}$ which
encode user demographics are extracted respecting privacy constraints. A
sequence model is used to generate features $z^{u}$ which summarizes the user
interests by encoding their past interactions with the system while graph
based models pretrained on a large-scale internal dataset (which goes beyond
interactions of users with the recommender system) are used to extract
features $i^{u}$ encoding user preferences. The item embeddings, denoted by
$a$, include textual and image representations along with its content
category. We use a fully connected deep neural network to model the Q-values
$Q_{\pi_{\text{base}}}(s^{u},a)$, by taking state and action representations
as inputs. We retrain our models every 12 hours to adapt to data distribution
shifts and maintain a replay buffer of approximately 300 million data tuples
at anytime to train our models using SARSA.

#### Experimental results:

Figures 4 and 4 summarize the results of our A/B test. In Figure 4, we
aggregate results week-by-week and summarize the overall improvement for both
the number of conversions as well as the conversion rate. On average we see a
lift (change between test and control group) between $4-10\%$ for both metrics
which shows a significant improvement over the base policy. Moreover, Fig. 4
shows that improvements of our approach hold consistently, over most days.
Overall, our results demonstrate the effectiveness of our approach in
optimizing for long run metrics. We remark that change in conversions closely
mimics the change in conversion rates as we tried to maintain a neutral
impression change. This can be seen in Fig. 4 where the change in impressions
for the test group over control is close to zero for all days.

##  6 Conclusion

In this work, we proposed a reinforcement learning based approach to account
for long run user engagement in recommending items to users while still
working within the constraints of an auction based recommender system.
Foundational theory in dynamic programming motivates our method – we show that
our method provably improves over the baseline policy. An online A/B test
conducted with a billion-user scale recommender system demonstrates
significant improvements and illustrates the potential of using reinforcement
learning ideas for recommender systems at scale in production.

## References

  * (1)
  * Bertsekas (2012) Dimitri Bertsekas. 2012\.  _Dynamic programming and optimal control: Volume I_. Vol. 1.  Athena scientific.
  * Chen et al. (2019a) Minmin Chen, Alex Beutel, Paul Covington, Sagar Jain, Francois Belletti, and Ed H Chi. 2019a.  Top-k off-policy correction for a REINFORCE recommender system. In _Proceedings of the Twelfth ACM International Conference on Web Search and Data Mining_. 456–464.
  * Chen et al. (2022) Minmin Chen, Can Xu, Vince Gatto, Devanshu Jain, Aviral Kumar, and Ed Chi. 2022\.  Off-Policy Actor-critic for Recommender Systems. In _Proceedings of the 16th ACM Conference on Recommender Systems_. 338–349.
  * Chen et al. (2019b) Xinshi Chen, Shuang Li, Hui Li, Shaohua Jiang, Yuan Qi, and Le Song. 2019b.  Generative adversarial user model for reinforcement learning based recommendation system. In _International Conference on Machine Learning_. PMLR, 1052–1061.
  * Cheng et al. (2016) Heng-Tze Cheng, Levent Koc, Jeremiah Harmsen, Tal Shaked, Tushar Chandra, Hrishi Aradhye, Glen Anderson, Greg Corrado, Wei Chai, Mustafa Ispir, et al. 2016\.  Wide & deep learning for recommender systems. In _Proceedings of the 1st workshop on deep learning for recommender systems_. 7–10.
  * Covington et al. (2016) Paul Covington, Jay Adams, and Emre Sargin. 2016\.  Deep neural networks for youtube recommendations. In _Proceedings of the 10th ACM conference on recommender systems_. 191–198.
  * Edelman et al. (2007) Benjamin Edelman, Michael Ostrovsky, and Michael Schwarz. 2007\.  Internet advertising and the generalized second-price auction: Selling billions of dollars worth of keywords.  _American economic review_ 97, 1 (2007), 242–259.
  * Ekstrand et al. (2011) Michael D Ekstrand, John T Riedl, Joseph A Konstan, et al. 2011\.  Collaborative filtering recommender systems.  _Foundations and Trends® in Human–Computer Interaction_ 4, 2 (2011), 81–173.
  * Evans (2008) David S Evans. 2008\.  The economics of the online advertising industry.  _Review of network economics_ 7, 3 (2008).
  * Gomez-Uribe and Hunt (2015) Carlos A Gomez-Uribe and Neil Hunt. 2015.  The netflix recommender system: Algorithms, business value, and innovation.  _ACM Transactions on Management Information Systems (TMIS)_ 6, 4 (2015), 1–19.
  * Gu et al. (2017) Shixiang Gu, Ethan Holly, Timothy Lillicrap, and Sergey Levine. 2017\.  Deep reinforcement learning for robotic manipulation with asynchronous off-policy updates. In _2017 IEEE international conference on robotics and automation (ICRA)_. IEEE, 3389–3396.
  * Guo et al. (2023) Hongbo Guo, Ruben Naeff, Alex Nikulkov, and Zheqing Zhu. 2023\.  Evaluating Online Bandit Exploration In Large-Scale Recommender System. In _KDD-23 Workshop on Multi-Armed Bandits and Reinforcement Learning: Advancing Decision Making in E-Commerce and Beyond_.
  * Hidasi et al. (2015) Balázs Hidasi, Alexandros Karatzoglou, Linas Baltrunas, and Domonkos Tikk. 2015.  Session-based recommendations with recurrent neural networks.  _arXiv preprint arXiv:1511.06939_ (2015).
  * Hutter (2007) Marcus Hutter. 2007\.  Universal algorithmic intelligence: A mathematical top→ down approach.  _Artificial general intelligence_ (2007), 227–290.
  * Ie et al. (2019) Eugene Ie, Vihan Jain, Jing Wang, Sanmit Narvekar, Ritesh Agarwal, Rui Wu, Heng-Tze Cheng, Tushar Chandra, and Craig Boutilier. 2019.  SlateQ: A tractable decomposition for reinforcement learning with recommendation sets.  (2019).
  * Koren et al. (2009) Yehuda Koren, Robert Bell, and Chris Volinsky. 2009\.  Matrix factorization techniques for recommender systems.  _Computer_ 42, 8 (2009), 30–37.
  * Li et al. (2010) Lihong Li, Wei Chu, John Langford, and Robert E Schapire. 2010\.  A contextual-bandit approach to personalized news article recommendation. In _Proceedings of the 19th international conference on World wide web_. 661–670.
  * Lops et al. (2011) Pasquale Lops, Marco De Gemmis, and Giovanni Semeraro. 2011\.  Content-based recommender systems: State of the art and trends.  _Recommender systems handbook_ (2011), 73–105.
  * Lu et al. (2015) Jie Lu, Dianshuang Wu, Mingsong Mao, Wei Wang, and Guangquan Zhang. 2015.  Recommender system application developments: a survey.  _Decision support systems_ 74 (2015), 12–32.
  * Lu et al. (2021) Xiuyuan Lu, Benjamin Van Roy, Vikranth Dwaracherla, Morteza Ibrahimi, Ian Osband, and Zheng Wen. 2021\.  Reinforcement learning, bit by bit.  _arXiv preprint arXiv:2103.04047_ (2021).
  * Mnih et al. (2015) Volodymyr Mnih, Koray Kavukcuoglu, David Silver, Andrei A Rusu, Joel Veness, Marc G Bellemare, Alex Graves, Martin Riedmiller, Andreas K Fidjeland, Georg Ostrovski, et al. 2015\.  Human-level control through deep reinforcement learning.  _nature_ 518, 7540 (2015), 529–533.
  * Okura et al. (2017) Shumpei Okura, Yukihiro Tagami, Shingo Ono, and Akira Tajima. 2017\.  Embedding-based news recommendation for millions of users. In _Proceedings of the 23rd ACM SIGKDD international conference on knowledge discovery and data mining_. 1933–1942.
  * Pazzani and Billsus (2007) Michael J Pazzani and Daniel Billsus. 2007.  Content-based recommendation systems.  _The adaptive web: methods and strategies of web personalization_ (2007), 325–341.
  * Qin et al. (2014) Lijing Qin, Shouyuan Chen, and Xiaoyan Zhu. 2014.  Contextual combinatorial bandit and its application on diversified online recommendation. In _Proceedings of the 2014 SIAM International Conference on Data Mining_. SIAM, 461–469.
  * Quadrana et al. (2018) Massimo Quadrana, Paolo Cremonesi, and Dietmar Jannach. 2018\.  Sequence-aware recommender systems.  _ACM Computing Surveys (CSUR)_ 51, 4 (2018), 1–36.
  * Schafer et al. (2007) J Ben Schafer, Dan Frankowski, Jon Herlocker, and Shilad Sen. 2007.  Collaborative filtering recommender systems.  _The adaptive web: methods and strategies of web personalization_ (2007), 291–324.
  * Shani et al. (2005) Guy Shani, David Heckerman, Ronen I Brafman, and Craig Boutilier. 2005.  An MDP-based recommender system.  _Journal of Machine Learning Research_ 6, 9 (2005).
  * Shi et al. (2014) Yue Shi, Martha Larson, and Alan Hanjalic. 2014.  Collaborative filtering beyond the user-item matrix: A survey of the state of the art and future challenges.  _ACM Computing Surveys (CSUR)_ 47, 1 (2014), 1–45.
  * Smith and Linden (2017) Brent Smith and Greg Linden. 2017.  Two decades of recommender systems at Amazon. com.  _Ieee internet computing_ 21, 3 (2017), 12–18.
  * Tai et al. (2017) Lei Tai, Giuseppe Paolo, and Ming Liu. 2017.  Virtual-to-real deep reinforcement learning: Continuous control of mobile robots for mapless navigation. In _2017 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 31–36.
  * Varian (2007) Hal R Varian. 2007\.  Position auctions.  _international Journal of industrial Organization_ 25, 6 (2007), 1163–1178.
  * Varian and Harris (2014) Hal R Varian and Christopher Harris. 2014.  The VCG auction in theory and practice.  _American Economic Review_ 104, 5 (2014), 442–445.
  * Wang et al. (2017) Huazheng Wang, Qingyun Wu, and Hongning Wang. 2017.  Factorization bandits for interactive recommendation. In _Proceedings of the AAAI Conference on Artificial Intelligence_ , Vol. 31.
  * Wang et al. (2019) Shoujin Wang, Liang Hu, Yan Wang, Longbing Cao, Quan Z Sheng, and Mehmet Orgun. 2019\.  Sequential recommender systems: challenges, progress and prospects.  _arXiv preprint arXiv:2001.04830_ (2019).
  * Xin et al. (2020) Xin Xin, Alexandros Karatzoglou, Ioannis Arapakis, and Joemon M Jose. 2020.  Self-supervised reinforcement learning for recommender systems. In _Proceedings of the 43rd International ACM SIGIR conference on research and development in Information Retrieval_. 931–940.
  * Zeng et al. (2016) Chunqiu Zeng, Qing Wang, Shekoofeh Mokhtari, and Tao Li. 2016\.  Online context-aware recommendation with time varying multi-armed bandit. In _Proceedings of the 22nd ACM SIGKDD international conference on Knowledge discovery and data mining_. 2025–2034.
  * Zhang et al. (2019) Shuai Zhang, Lina Yao, Aixin Sun, and Yi Tay. 2019\.  Deep learning based recommender system: A survey and new perspectives.  _ACM computing surveys (CSUR)_ 52, 1 (2019), 1–38.
  * Zhao et al. (2019) Xiangyu Zhao, Long Xia, Jiliang Tang, and Dawei Yin. 2019\.  ” Deep reinforcement learning for search, recommendation, and online advertising: a survey” by Xiangyu Zhao, Long Xia, Jiliang Tang, and Dawei Yin with Martin Vesely as coordinator.  _ACM Sigweb Newsletter_ Spring (2019), 1–15.
  * Zheng et al. (2018) Guanjie Zheng, Fuzheng Zhang, Zihan Zheng, Yang Xiang, Nicholas Jing Yuan, Xing Xie, and Zhenhui Li. 2018.  DRN: A deep reinforcement learning framework for news recommendation. In _Proceedings of the 2018 world wide web conference_. 167–176.
  * Zhu and Van Roy (2023a) Zheqing Zhu and Benjamin Van Roy. 2023a.  Deep Exploration for Recommendation Systems. In _Proceedings of the 17th ACM Conference on Recommender Systems_.
  * Zhu and Van Roy (2023b) Zheqing Zhu and Benjamin Van Roy. 2023b.  Scalable Neural Contextual Bandit for Recommender Systems.  _arXiv preprint arXiv:2306.14834_ (2023).
  * Zou et al. (2019) Lixin Zou, Long Xia, Zhuoye Ding, Jiaxing Song, Weidong Liu, and Dawei Yin. 2019\.  Reinforcement learning to optimize long-term user engagement in recommender systems. In _Proceedings of the 25th ACM SIGKDD International Conference on Knowledge Discovery & Data Mining_. 2810–2818.


[a]D. Bala

# Photon production rate from Transverse-Longitudinal ($T-L$) mesonic
correlator on the lattice

S. Ali  A. Francis  G. Jackson  O. Kaczmarek  T. Ueding

###### Abstract

Thermal photons from the QGP provide important information about the
interaction among plasma constituents. The photon production rate from a
thermally equilibrated system is proportional to the transverse spectral
function $\rho_{T}(\omega=|\vec{k}|,\vec{k})$. One can also calculate the
photon production rate from the difference between $\rho_{T}(\omega,\vec{k})$
(transverse) and $\rho_{L}(\omega,\vec{k})$ (longitudinal) projections, as
$\rho_{L}$ vanishes on the photon point. Because the UV part of
$\rho_{T}-\rho_{L}$ is suppressed, the corresponding Euclidean correlator
receives most of its contribution from the IR part. We calculate the
$T\\!-\\!L$ correlator on $N_{f}=2+1$ flavour HISQ configurations with
$m_{l}=m_{s}/5$ at temperature of about $1.15\,T_{pc}$ (220 MeV). We have used
two ansätze for the spectral function: 1) A polynomial connected to the UV
region consistent with OPE expansion and 2) a hydro-inspired spectral
function. We have also applied the Backus-Gilbert method to estimate the
spectral function. All these different approaches are combined to estimate the
photon production rate.

##  1 Introduction

Lattice QCD determinations of the equation of state (EoS) confirm the
transition from hadronic matter to the quark-gluon plasma (QGP) phase under
increasing temperature [1, 2]. Experimentally this QGP phase has been studied
in high-energy heavy ion collisions at the RHIC and LHC facilities. Photons
and dileptons produced in the plasma are a vital probe of the QGP because the
mean free paths of the photon and the dilepton (i.e. virtual photon) are much
larger than the typical size of the plasma [3]. Therefore photons and
dileptons can transfer information about the plasma from the point at which it
was created and throughout its spacetime evolution.

The photon and dilepton production rate from a thermalized QGP, at temperature
$T$, to leading-order in the electromagnetic coupling $\alpha_{\rm em}$ are
given by [4],

$\frac{d\Gamma_{\gamma}}{d^{3}{\vec{k}}}\;=\;\frac{\alpha_{\rm
em}n_{B}(k)}{\pi^{2}k}\;\Big{(}{\textstyle\sum_{i=1}^{N_{f}}Q_{i}^{2}}\Big{)}\;\rho_{T}(k,\vec{k})\,,$
(1)
$\frac{d\Gamma_{l^{+}l^{-}}}{d{\omega}d^{3}{\vec{k}}}\;\simeq\;\frac{\alpha_{\rm
em}^{2}n_{B}(\omega)}{3\pi^{2}(\omega^{2}-k^{2})}\;\Big{(}{\textstyle\sum_{i=1}^{N_{f}}Q_{i}^{2}}\Big{)}\;\big{(}\,2\rho_{T}(\omega,\vec{k})+\rho_{L}(\omega,\vec{k})\,\big{)}\,,$
(2)

where $n_{B}$ is the Bose distribution and $Q_{i}$ is the charge of a quark of
flavour $i$ in units of the electron charge. All QCD information is encoded in
$\rho_{T}(\omega,\vec{k})$ and $\rho_{L}(\omega,\vec{k})\,$ which denote the
transverse and longitudinal spectral functions at frequencies $\omega$ and
momenta $\vec{k}\,$.

Perturbation theory is an important tool to calculate the spectral function in
the weak coupling limit. One can use naive (NLO) perturbation theory away from
the light cone to calculate the spectral function directly [5]. However, near
the light cone, resummation (LPM) of specific Feynman diagrams is essential to
obtain correct results [6]. In general, one must interpolate between the two
regimes [7]111 These calculations have now been extended to finite baryon
chemical potential [8]..

On the lattice, however, we need to calculate spectral functions from a
lattice correlator by using the following relation,

$G_{E}(\tau,\vec{k})\;=\;\int_{0}^{\infty}\frac{d\omega}{\pi}\;\rho(\omega,\vec{k})\;\frac{\cosh[\omega(1/2T-\tau)]}{\sinh(\omega/2T)}\,.$
(3)

Extraction of the spectral function from noisy lattice data, using the above
equation, is a well-known numerically ill-posed problem. Physically motivated
assumptions are required to extract the spectral function from the lattice
correlator.

Following Ref. [10], we estimate the photon production rate from the
$G_{H}=2(G_{T}-G_{L})$ correlator rather than $G_{T}$. This is because the
spectral function $\rho_{H}=2(\rho_{T}-\rho_{L})$ is identical to $\rho_{T}$
at $\omega=k$ (photon point), and the UV part of $\rho_{H}$ is much suppressed
compared to the IR part of the spectral function. Therefore the resultant
lattice correlator $G_{H}$ is mostly dominated by the important IR part of the
spectral function, in contrast to the $G_{T}$ correlator, which is mostly
dominated by the UV part of $\rho_{T}$. Previously, the photon production rate
has also been estimated from the correlator, which has a large UV contribution
[14, 15]222Recently, the photon production rate has also been estimated from
the $G_{T}$ correlator [9]. In Ref.[11] photon production rate has also been
estimated from imaginary momentum correlators . It is known that the resulting
$T-L$ spectral function also satisfies the following sum rule in the chiral
limit [20],

$\int_{0}^{\infty}d\omega\,\omega\;\rho_{H}(\omega,\vec{k})\;=\;0\,.$ (4)

Our aim is then to compute the $G_{H}$ correlator on the lattice and make
ansatz for the spectral function to fit the lattice data such that it
satisfies the above sum rule. This proceedings contribution is organized as
follows. In the next section, we will present lattice details. In Sec. 3,
reconstruction of the spectral function is performed using two models and the
Backus-Gilbert (BG) method. Sec. 4 contains the final results of the effective
diffusion coefficient (directly related to the photon production rate)
extracted from the spectral functions from Sec. 3.

##  2 Lattice Details

To calculate the correlation functions, we use $N_{f}=2+1$ flavor QCD
gaugefield configurations generated with the HISQ action by the HotQCD
collaboration. These configurations correspond to an unphysical pion mass
$m_{\pi}=315$ MeV with $m_{l}=m_{s}/5$. We used the lattice size $96^{3}\times
32$, corresponding to a temperature of $220$ MeV. On these configurations, we
calculate the correlation function using clover-improved Wilson fermions with
the parameters $\kappa=0.13515$ and tadpole improved tree-level
$c_{sw}=1.34108$. We have used around 1700 configurations for the computation
of these correlators. These parameters for Wilson clover fermion are chosen to
reproduce the above staggered pion mass. The possible spatial momentum is
fixed by the aspect ratio of the lattice, which is, in this case,
$k_{n}/T=2\pi nN_{t}/N_{s}=2\pi n/3$ where $n$ is an integer.

Figure 1: The Euclidean hadronic correlator $G_{H}$ as measured on the
lattice, at the points $\tau
T=\\{\frac{2}{32},\frac{3}{32},\ldots,\frac{16}{32}\\}$, for various momenta
(namely $k_{n}$ for $n=1,2,3$). The continuous lines are the corresponding
estimates from resummed perturbation theory (NLO+LPM) as described in Ref.
[13].

##  3 Lattice correlator and spectral reconstruction

The $G_{H}$ correlator calculated on the lattice is multiplicatively
renormalizable. However, the ratio $G_{H}/2\chi_{q}T$ does not need
renormalization. Here $\chi_{q}$ is the quark number susceptibility. On the
lattice, we can calculate the correlation function $G_{\mu\nu}$ of
electromagnetic current $J_{\mu}$, from which one can easily calculate this
ratio by using the following formula,

$\frac{G_{H}(\tau,\vec{k}=(0,0,k_{z}))}{2\chi_{q}T}=\frac{G_{xx}(\tau,\vec{k})+G_{yy}(\tau,\vec{k})-2(G_{zz}(\tau,\vec{k})-G_{00}(\tau,\vec{k}))}{2\,G_{00}(\tau,\vec{k})}.$
(5)

The renormalized lattice correlator is obtained by multiplying this ratio by
$2\chi_{q}$. We have estimated $\chi_{q}\approx 0.872\,T^{2}$ from the order
$g^{6}\ln(g)$ perturbative calculation [12]. Therefore we are essentially
renormalizing the correlator perturbatively. The resulting perturbatively
renormalized correlator is shown in Fig. 1. In this figure, we have also shown
the NLO+LPM perturbative estimate of the corresponding correlator [13]. The
difference between the perturbative and lattice estimates of the correlator
signifies non-perturbative effects.

As mentioned, spectral reconstruction from lattice using Eq. (3) is an ill-
posed problem. Therefore various possible spectral functions can reproduce the
lattice data and well-motivated assumptions are required in practice for
spectral reconstruction. Here we will present two physically motivated models
of the spectral function that can fit our lattice data.

###  3.1 Model spectral functions

For the first model to be considered, we have used a spectral function which
has the following piecewise polynomial form ($\omega_{0}$ is a parameter of
the model)

$\rho_{H}^{\rm
poly}(\omega)\;=\;\rho_{<}\,\Theta(\omega_{0}-\omega)+\rho_{>}\,\Theta(\omega-\omega_{0})\,,$
(6)

where

$\displaystyle\rho_{<}$
$\displaystyle\equiv\frac{\beta\omega^{3}}{2\omega_{0}^{3}}\left(5-3\frac{\omega^{2}}{\omega_{0}^{2}}\right)-\frac{\gamma\omega^{3}}{2\omega_{0}^{2}}\left(1-\frac{\omega^{2}}{\omega_{0}^{2}}\right)+\frac{\delta\,\omega}{\omega_{0}}\left(1-\frac{\omega^{2}}{\omega_{0}^{2}}\right)^{2}\,,$
(7) $\displaystyle\rho_{>}$
$\displaystyle\equiv-\frac{\beta\,\omega_{0}^{4}}{7\,\omega^{4}}\left(54\frac{\omega_{0}^{4}}{\omega^{4}}\,-\,94\frac{\omega_{0}^{2}}{\omega^{2}}\,+\,33\right)+\frac{\gamma\,\omega_{0}^{5}}{140\,\omega^{4}}\left(-81\frac{\omega_{0}^{4}}{\omega^{4}}\,+\,92\frac{\omega_{0}^{2}}{\omega^{2}}\,-\,11\right)-\frac{16\,\delta\,\omega_{0}^{4}}{35\,\omega^{4}}\left(1-\frac{\omega_{0}^{2}}{\omega^{2}}\right)^{2}\,,$
(8)

and $\Theta$ is the Heaviside step function. The part of the spectral function
for $\omega<\omega_{0}$ is taken from Ref. [15] and has been used based on the
smoothness of the spectral function across the light cone. The second form of
the expansion is motivated by the operator product expansion, according to
which $\rho_{H}(\omega)\sim 1/\omega^{4}$ at large $\omega$ [20].

The various coefficients appearing in Eq. (7) and Eq. (8) have been chosen
such that $\rho_{H}(\omega)$ is smooth at $\omega=\omega_{0}$ and satisfies
the sum rule in Eq. (4). The parameters $\beta$ and $\gamma$ correspond to the
value of the spectral function and its derivative, respectively, at
$\omega_{0}$. The parameter $\delta$ controls the slope of the spectral
function at the origin. We perform the fit of our data with respect to
$\beta$, $\gamma$, and $\delta$, under the constraints $\delta\geq 0$,
$\rho_{H}(k,\vec{k})\geq 0$ and $\frac{\partial G_{H}(\tau)}{\partial\tau}\leq
0$. The first constraint comes from the fact that in hydrodynamics, $\delta$
is related to the diffusion coefficient and, therefore, cannot be negative.
The second constraint amounts to the fact that the photon production rate
cannot be negative.

Figure 2: Mock analysis of perturbative data, showing data points from NLO+LPM
correlators used for the spectral reconstruction, the error on each point is
$0.1$% (left). Red bands are the correlator obtained from the fitted spectral
function, with the width of the band stemming from uncertainty in the fit
parameters. The black curves are the perturbative correlator at all $\tau T$
and thus pass through the mock data points. On the right, we show the
reconstructed spectral function from ‘perturbative data,’ alongside the exact
spectral function as a red dotted curve.

Before fitting the above spectral ansatz to the lattice data, we performed a
mock data analysis wherein we reconstruct the exactly known NLO+LPM spectral
function from the known NLO+LPM correlation function. To perform the fit with
above the spectral function, we have artificially introduced an error into the
perturbative correlator of the following form $\frac{\delta
G_{H}}{G_{H}}=0.001$. For this reconstruction, we used 10 data points between
$\tau T=0.1875$ to $\tau T=0.5$. The resulting spectral function, along with
the fitted correlator, is shown in Fig. 2. In the left panel of Fig. 2, we
depict the perturbative data points used for the spectral reconstruction. The
red band corresponds to the correlator obtained from the fitted spectral
function, and the black curve is the perturbative correlator. We see that we
could predict data points beyond the fitting range with the above-fitted
ansatz. On the right panel of Fig. 2, the reconstructed spectral function for
a few values of the matching point $\omega_{0}$ is shown. The red dotted line
corresponds to the exact spectral function for the correlator. We observe that
the spectral function shows some change with respect to $\omega_{0}$, however
within a systematic uncertainty of $\omega_{0}$ between $\sqrt{k^{2}+(\pi
T)^{2}}$ and $\sqrt{k^{2}+5(\pi T)^{2}}$ we could indeed be able to encompass
the exact spectral function.

Figure 3: (Left) The data points correspond to the lattice correlator
calculated at various momenta. The red band is the correlator from the fitted
spectral function. (Right) The reconstructed spectral function for $k/T=4.19$
along with the perturbative spectral function.

After verifying the above spectral function with the perturbative data, we now
proceed to fit the lattice correlator. The fitting range we have taken is in
$\tau T$ between $0.1875$ and $0.5$, as the short distance part is dominated
by lattice artifacts. The spectral function for various $\omega_{0}$ is shown
in the right panel of Fig. (3). The red shaded band in the right panel
corresponds to the fitted correlator. We see the expected feature that the
spectral function is IR-dominated. For the calculation of photon production
rate, we take the variation with $\omega_{0}$ as a source of systematic error
from this ansatz. We have also shown the corresponding NLO+LPM perturbative
spectral function in the same figure. We see that at $\omega\ll k$, the
perturbative spectral function is suppressed compared to our model spectral
function, whereas near the light cone, the model spectral function is
consistently below the value of the perturbative prediction.

Figure 4: Spectral function for various momenta, as a function of $\omega$. In
the left panel, we show the result from a polynomial ansatz (6), and the right
panel depicts the best fit from a hydro ansatz (9). Both functions have the
correct UV behavior.

We show the spectral function for $\omega_{0}=\sqrt{k^{2}+(\pi T)^{2}}$ for
various momenta in Fig. (4). We observe that the IR part of the spectral
function is suppressed for higher momentum compared to small momentum. This
behavior is expected because the spectral function becomes
$\omega\,\delta(\omega)$ by rotational symmetry at zero momentum. At finite
momentum, this $\delta(\omega)$ function develops a finite width, which
increases with increasing momentum.

The second ansatz is motivated by relaxation-time hydrodynamics [21]. The
spectral function in this model is taken as

$\rho^{\rm
hydro}_{H}(\omega,\vec{k})\;=\;A\frac{\tanh(\omega/T)\,k^{2}\,(1+a^{2}\,k^{2}-2\,a\,b\,\omega^{2}+b^{2}\omega^{2})}{(1+b^{2}\,\omega^{2})(a^{2}\,k^{4}+\omega^{2}-2a\,k^{2}\,b\,\omega^{2}+b^{2}\,\omega^{4})}\,,$
(9)

with parameters $A$, $a$ and $b\,$. At large $\omega$, this spectral function
goes like $1/\omega^{4}$ in accordance with the OPE prediction. The sum rule,
on the other hand, will allow the parameter $a$ to be expressed in terms of
$b$. Then we perform a two-parameter fit of this spectral function with the
lattice data. The resulting spectral function for various momenta has been
shown in the right panel of Fig. 4.

Figure 5: Backus-Gilbert estimated spectral function, (Left) Variation of
spectral function with $\lambda$ for $k/T=4.19$. (Right) Spectral function for
$\lambda=0.01$ for various momenta.

###  3.2 Backus-Gilbert estimate of $\rho_{H}(\omega,\vec{k})$

The third approach that we have used is the Backus-Gilbert (BG) method [16].
BG method has been widely used for the spectral reconstruction from lattice
correlator [17, 18]. In this method, one rewrites the integral transform (3)
in the following form,

$G_{H}(\tau,\vec{k})\;=\;\int_{0}^{\infty}d\omega\,\frac{\rho(\omega,\vec{k})}{f(\omega,\vec{k})}\,K(\omega,\tau)\;;\quad
K(\omega,\tau)\equiv\frac{f(\omega,\vec{k})\cosh[\omega(\beta/2-\tau)]}{\sinh(\omega\beta/2)\pi}\,.$
(10)

The purpose of the function $f(\omega,\vec{k})$ is to put physics information
in the reconstruction process and remove the singularity originating from the
denominator of the kernel.

The BG strategy is to write the spectral function as a linear combination of
the temporal lattice correlator, ($i$ enumerates the measured points)

$\frac{\rho_{H}^{\rm
BG}(\omega,\vec{k})}{f(\omega,\vec{k})}\;=\;\sum_{i}q_{i}(\omega,\vec{k})\,G_{H}(\tau_{i},\vec{k})\\\
\;=\;\int_{0}^{\infty}d\omega^{\prime}\delta(\omega^{\prime},\omega)\frac{\rho_{H}(\omega^{\prime},\vec{k})}{f(\omega^{\prime},\vec{k})}\,,$
(11)

where
$\delta(\omega,\omega^{\prime})=\sum_{i}q_{i}(\omega,\vec{k})K(\omega^{\prime},\tau_{i})$
is the so-called resolution function. The BG estimated spectral function would
be the exact spectral function when
$\delta(\omega,\omega^{\prime})=\delta(\omega-\omega^{\prime})$. Therefore the
goal is to minimize the width of the resolution function as small as possible
by varying $q_{i}$, which can be done by minimizing the following function,

$A(\omega)\;=\;\int
d\bar{\omega}(\omega-\bar{\omega})^{2}\delta(\omega,\bar{\omega})^{2}\;=\;{q^{t}(\omega).W(\omega).q(\omega)}\,,$
(12)

where
$W_{ij}(\omega)\equiv\int_{0}^{\infty}d\bar{\omega}\,(\omega-\bar{\omega})^{2}K(\bar{\omega},\tau_{i})K(\bar{\omega},\tau_{j})$.
Minimization of this function requires the inversion of the matrix $W$.
However, the $W$ matrix is badly conditioned, which causes the minimization to
be unstable. The minimization of the following regularized function can
overcome this,

$F(\omega)=\lambda A(\omega)+(1-\lambda)B(\omega)$ (13)

where $B(\omega)={\rm Var}[\rho^{\rm BG}(\omega)]={q^{t}(\omega).S.q(\omega)}$
with $S$ being the covariance matrix. Various other regularization schemes
have also been used in literature, for example, one such regularization is the
Tikhonov scheme where $S_{ij}=\delta_{ij}$ [19]. The parameter $\lambda$ is
related to the width of the resolution function and as well as the error on
the spectral function. The ‘better’ choices of $f(\omega,\vec{k})$ should also
help reduce the resolution function’s width.

For the calculation of the BG estimated spectral function for our lattice
correlator, we used the following function,

$f(\omega)\;=\;\frac{\tanh(\omega/T)}{1+(\omega/\omega_{0})^{2}+(\omega/\omega_{0})^{4}}$
(14)

where, $\omega_{0}=\sqrt{k^{2}+\pi^{2}T^{2}}$. This functional form is
motivated by the observation that at small $\omega$, the spectral function
should be proportional to $\omega$, whereas, at large $\omega$, the spectral
function should go like $1/\omega^{4}$.

Figure 6: Resolution function obtained from the BG method on the light cone
for various momenta for $\lambda=0.01$.

The BG estimated spectral function is plotted for various values of $\lambda$
in the left panel of Fig. 5. The variation with respect to $\lambda$ is taken
to the systematic error. The spectral function for various momenta for
$\lambda=0.01$ is shown in the right panel of Fig. 5. We observe the same
qualitative trend of the spectral function as from the model fit of the data.
The resolution function at the photon point is plotted in Fig. 6 for all
available momenta. We note that the resolution function is peaked around
$\omega_{0}\simeq|\vec{k}|\,$.

Figure 7: (Left) Effective diffusion coefficient obtained from various
methods. (Right) Effective diffusion coefficients obtained from various models
have been combined to get a systematic error band. In both panels, the
perturbative band is give by varying the scale choice in the running coupling
– see Ref. [13] for details.

##  4 Effective diffusion coefficient

The photon production rate can be rewritten in terms of the so-called
effective diffusion coefficient,

$D_{\rm eff}(k)\;\equiv\;\frac{\rho_{H}(k,\vec{k})}{2\chi_{q}k}\,.$ (15)

The value of $D_{\rm eff}(k)$ estimated from two models, and the BG method is
shown in the left panel of Fig. 7. We see that at small momentum, the model
prediction has considerable uncertainty, whereas, at large momentum, different
methods seem to agree mutually. We also observe that at large momentum,
$D_{\rm eff}$ is smaller than the perturbative estimate. In the right panel,
we combined the model results to get a systematic uncertainty on $D_{\rm
eff}$.

##  5 Conclusion

We estimated the photon production rate from the $T-L$ correlator from
$N_{f}=2+1$ flavor QCD for $m_{\pi}=313$ MeV at $T=220$ MeV. The $T-L$
correlator was calculated using Clover-improved Wilson fermions on HISQ
configurations. Various models have been used for spectral reconstruction from
the lattice correlator. For the first model, we have used a piecewise
polynomial which smoothly connects a quintic function of $\omega$ in the IR
with the asymptotic $1/\omega^{4}$ form in the UV. We apply this model first
to a mock correlator from NLO+LPM resummed perturbation theory, where we know
the exact spectral function. Within systematic uncertainty, we could reproduce
this perturbative spectral function. Then we fit this model to the lattice
correlator, with the result for various momenta shown in the left panel of
Fig. 4. A model based on relaxation time hydrodynamics was similarly
considered, with the resulting spectral function being shown in the right
panel of Fig. 4. We have also used the Backus-Gilbert method to estimate the
spectral function. In all these methods, we found qualitative agreement;
namely, the spectral function is IR-dominated. The effective diffusion
coefficient from all these methods is shown in Fig. 7. We use the spread in
values obtained via the different methods as a proxy for the systematic
uncertainty from reconstructing the spectral function.

##  6 Acknowledgement

We thank Luis Altenkort and Hai-Tao Shu for the production of gauge
configurations and implementation of meson measurements in the QUDA code. This
work is supported by the Deutsche Forschungsgemeinschaft (DFG, German Research
Foundation)-Project number 315477589-TRR 211. For the computations, we used
high-performance facilities from the Bielefeld and Juelich machines. G. J. was
supported by the U.S. Department of Energy (DOE) under grant No. DE-
FG02-00ER41132. A. F. acknowledges support by the Ministry of Science and
Technology Taiwan (MOST) under grant 111-2112-M-A49-018-MY2.

## References

  * [1] A. Bazavov et al, (HotQCD Collaboration), Phys. Rev. D 90, (2014) 094503, [arXiv: 1407.6387[hep-lat]].
  * [2] S. Borsányi, Z. Fodor, C. Hoelbling, S. D.Katz, S. Krieg, Kálmán K. Szabó, Phys. Lett. B 730, (2014) 99 [arXiv: 1309.5258[hep-lat]].
  * [3] A. Adare et al, (PHENIX Collaboration), Phys. Rev. C 91, (2015) 064904 [arXiv: 1405.3940[nucl-ex]].
  * [4] L. D. McLerran and T. Toimela, Phys. Rev. D 31, (1985) 545.
  * [5] M. Laine, JHEP 05, (2013) 083 [arXiv: 1304.0202[hep-ph]].
  * [6] P. B. Arnold, G. D. Moore and L. G. Yaffe, JHEP 12, (2001) 009 [arXiv: hep-ph/0111107[hep-ph]].
  * [7] I. Ghisoiu and M. Laine, JHEP 10, (2014) 083 [arXiv: 1407.7955[hep-ph]].
  * [8] G. Jackson, [arXiv: 2211.09575[hep-ph]].
  * [9] M. Cè, T. Harris, A. Krasniqi, H. B. Meyer, and C. Török, Phys. Rev. D 106, (2022) 054501 [arXiv: 2205.02821[hep-lat]].
  * [10] M. Cè, T. Harris, H. B. Meyer, A. Steinberg, and A. Toniato, Phys. Rev. D 102, (2020) 091501[arXiv: 2001.03368[hep-lat]].
  * [11] C. Török, M. Cè, T. Harris, A. Krasniqi, H. B. Meyer and S. Ruhl, [arXiv: 2212.05622[hep-lat]].
  * [12] A. Vuorinen, Phys. Rev. D 67, (2003) 074032 [arXiv: hep-ph/0212283[hep-ph]].
  * [13] G. Jackson and M. Laine, JHEP 11, (2019) 144 [arXiv: 1910.09567[hep-ph]].
  * [14] H.-T. Ding, A. Francis, O. Kaczmarek, F. Karsch, E. Laermann, and W. Soeldner, Phys. Rev. D 83, (2011) 034504 [arXiv: 1012.4963[hep-lat]].
  * [15] J. Ghiglieri, O. Kaczmarek, M. Laine, and F. Meyer, Phys. Rev. D 94, (2016) 016005 [arXiv: 1604.07544[hep-lat]].
  * [16] G. Backus and F. Gilbert, Geophysical Journal of the Royal Astronomical Society 16, (1968) 169.
  * [17] B. B. Brandt, A. Francis, H. B. Meyer, and D. Robaina, Phys. Rev. D 92, (2015) 094510 [arXiv: 1506.05732[hep-lat]].
  * [18] B. B. Brandt, A. Francis, B. Jäger, and H. B. Meyer, Phys. Rev. D 93, (2016) 054510 [arXiv: 1512.07249[hep-lat]].
  * [19] N. Astrakhantsev, V. Braguta, and A. Kotov, Phys. Rev.D 98, (2018) 054515 [arXiv: 1804.02382[hep-lat]].
  * [20] B. B. Brandt, A. Francis, Tim Harris, Harvey B. Meyer, and A. Steinberg, EPJ Web Conf. 175, (2018) 07044 [arXiv: 1710.07050[hep-lat]].
  * [21] S. Floerchinger, C. Gebhardt and K. Reygers, [arXiv: 2112.12497[nucl-th]].