# What does cosmology tell us about the mass of thermal–relic dark matter?

Rui An  Vera Gluscevic  Department of Physics and Astronomy, University of
Southern California, Los Angeles, CA, 90089, USA  Erminia Calabrese  School of
Physics and Astronomy, Cardiff University, The Parade, Cardiff, Wales CF24
3AA, UK  J. Colin Hill  Department of Physics, Columbia University, New York,
NY 10027, USA  Center for Computational Astrophysics, Flatiron Institute, New
York, NY 10010, USA

###### Abstract

The presence of light thermally coupled dark matter affects early expansion
history and production of light elements during the Big Bang Nucleosynthesis.
Specifically, dark matter that annihilates into Standard Model particles can
modify the effective number of light species in the universe
$N_{\mathrm{eff}}$, as well as the abundance of light elements created buring
BBN. These quantities in turn affect the cosmic microwave background (CMB)
anisotropy. We present the first joint analysis of small-scale temperature and
polarization CMB anisotropy from Atacama Cosmology Telescope (ACT) and South
Pole Telescope (SPT), together with Planck data and the recent primordial
abundance measurements of helium and deuterium to place comprehensive bounds
on the mass of light thermal–relic dark matter. We consider a range of models,
including dark matter that couples to photons and Standard-Model neutrinos. We
discuss the sensitivity of the inferred mass bounds on measurements of
$N_{\mathrm{eff}}$, primordial element abundances and the baryon density, and
quantify the sensitivity of our results to a possible existence of additional
relativistic species. We find that the combination of ACT, SPT, and Planck
generally leads to the most stringent mass constraint for dark matter that
couples to neutrinos, improving the lower limit by 40%–80%, with respect to
previous Planck analyses. On the other hand, the addition of ACT and SPT leads
to a slightly weaker bound on electromagnetically coupled particles, due to a
shift in the preferred values of $Y_{\mathrm{p}}$ and $N_{\mathrm{eff}}$
driven by the ground based experiments. In most scenarios, the combination of
CMB data has a higher constraining power than the primordial abundance
measurements alone, with the best results achieved when all data are combined.
Combining all CMB measurements with primordial abundance measurements, we rule
out masses below $\sim$4 MeV at 95% confidence, for all models. We show that
allowing for new relativistic species can weaken the mass bounds for dark
matter that couples to photons by up to an order of magnitude or more.
Finally, we discuss the reach of the next generation of the CMB experiments in
terms of probing the mass of the thermal relic dark matter.

##  I Introduction

A variety of astrophysical and cosmological observations indicate that a
significant fraction of the matter in our Universe is composed of dark matter
(DM) [1, 2, 3, 4, 5]. Despite many decades of dedicated searching, the nature
of DM remains a mystery and exploring its essence is one of the most
challenging tasks for fundamental physics today. Some of the most compelling
candidate models invoke DM that is in thermal equilibrium with the
Standard–Model (SM) particles in the early universe [6, 7, 8, 9]. Among the
thermal–relic models are the Weakly Interacting Massive Particles (WIMPs), the
main focus of many current direct detection experiments [10, 11]. While WIMPs
were originally considered to have masses in the GeV–TeV range, current null
results from direct detection searches have inspired compelling WIMP-like and
other models for light thermal–relic DM, with sub-GeV masses [12, 13, 14]. In
this work, we consider the most general observable consequences of light
thermal–relic DM, including the effects on the expansion history and the Big
Bang Nucleosynthesis (BBN), and infer comprehensive bounds on its particle
mass $m_{\chi}$, using all available data.

BBN took place in the very early universe, beginning a fraction of a second
after the Big Bang, and ending tens of minutes later. During this time, masses
and interactions of particles present in the primordial plasma have defined
the expansion history and thermal history, defining the abundances of light
elements created during BBN. In particular, if DM particles are in thermal
contact with the rest of the plasma during BBN (as is the case for majority of
the popular WIMP models in current literature [15, 16, 17, 18, 19, 20, 21, 22,
23, 24, 25]), and if their masses are around $\sim$1–20 MeV, they become non-
relativistic right around the time of BBN. The resulting DM annihilation into
SM particles can affect the early expansion history and the production of
light elements [26, 27, 28], including the abundance of helium–4,
$Y_{\mathrm{p}}$, and deuterium, $Y_{\mathrm{D}}\equiv(D/H)\times 10^{5}$.
Annihilation products can additionally alter the radiation content in the
universe, changing the effective number of light species, $N_{\mathrm{eff}}$.
These effects are captured by the cosmic microwave background (CMB) anisotropy
and can be directly traced by measurements of primordial abundances of light
elements in Lyman-$\alpha$ forest systems.

We focus on inferring the lower bound on $m_{\chi}$ from both these data sets,
for a minimal set of modeling assumptions, requiring only that DM is in
thermal equilibrium with the rest of the universe prior to BBN. In this
context, we consider two general scenarios: one in which DM couples
electromagnetically to the SM, and another in which it only couples to the SM
neutrinos. Within each scenario, we consider four types of DM spin statistic:
real scalar, complex scalar, Majorana Fermion, and Dirac Fermion. Notably, we
do not attempt to relate the cosmological effects of $m_{\chi}$ to the
strength of the coupling to the SM; provided that the coupling is sufficient
to ensure thermal equilibrium prior to BBN, the effects of $m_{\chi}$ are
independent on details of the interaction model. This choice means that we do
not provide a connection between the relic abundance of DM and its mass, nor
do we consider late–time DM annihilation constraints—these are additionally
dependent on the strength and type of the interaction at hand, and we leave
such considerations for future work.

We note that the CMB constraints on the mass of light DM have been previously
derived using Planck data and primordial element abundance measurements (e. g.
see [18, 29, 20, 21, 23, 30]). In this work, in addition to the latest Planck
measurements [5], we include the most recent public releases of the CMB small-
scale measurements from the Atacama Cosmology Telescope (ACT) [31, 32] and
South Pole Telescope (SPT) [33], as well as the primordial element abundance
measurements from [34]. ACT and SPT in particular provide higher resolution
measurements of the polarization anisotropy, increasing the sensitivity to
small angular scales, and thus enabling complementary constraints on
cosmological parameters to Planck data. Since a primary effect of $m_{\chi}$
on the CMB is through $N_{\mathrm{eff}}$, we further extend the analyses by
allowing for additional light degrees of freedom to exist in the universe. We
quantify the impact of this additional freedom on the inferred mass bounds.

We find that the DM mass bounds can have a notable dependence on the choice of
the data set and assumptions about other light degrees of freedom. In
particular, the inclusion of ACT and SPT data leads to slightly different
preferred values of $Y_{\mathrm{p}}$ and $N_{\mathrm{eff}}$, as compared to
Planck alone, resulting in an improved mass limit for neutrino coupled DM and
a slightly less stringent bound on electromagnetically coupled DM mass.
Combining all CMB measurements with the measurements of helium and deuterium
abundance, we infer the lower mass limit of $\sim$4 MeV at $95\%$ confidence
level (CL), regardless of the model. We further show that allowing for the
existence of new relativistic species, parameterized by the effictive number
of additional (neutrino–like) degrees of freedom $\Delta N_{\nu}$, the mass
bounds are significantly less stringent for electromagnetically coupled
models. We discuss the implications of possible systematic effects on the
reported mass inference. Finally, we show that the future–generation
ground–based CMB measurements from the Simons Observatory [35] and CMB-S4 [36,
37] will saturate the precision of CMB mass bounds.

The organization of the paper is as follows. In Section II, we briefly review
the most relevant physics behind the standard BBN (SBBN) model, and discuss
the alterations of SBBN that include light thermal–relic DM. In Section III,
we quantify the impact of light DM on the CMB power spectrum. In Section IV,
we present the current observational data sets, our analysis method, and the
resulting constraints on DM mass. In Section V, we show the projected
sensitivity of the upcoming CMB experiments. Finally, we summarize and discuss
our findings in Sec. VI.

##  II Effects on primordial abundances

Within the standard models of particle physics and cosmology, the universe was
radiation-dominated during BBN. At the beginning of this process, it contained
electrons and positrons (denoted here as $e^{\pm}$), photons ($\gamma$), three
neutrino species ($\nu$), a small number of protons ($p$) and neutrons ($n$),
and DM particles—all in thermal equilibrium with each other. The SM particles
were initially tightly coupled via electromagnetic and weak interactions, but
as the universe expanded and cooled, the rate of weak interactions dropped,
resulting in neutrino decoupling at the temperature of $\sim$few MeV [38, 39,
40]. Later on, at $T_{\gamma}$$\sim$$m_{e}$, the rest of the plasma received
an energy injection from electron-positron annihilation; within the standard
cosmological model, this process heats photons with respect to neutrinos,
leading to the present-day neutrino–to–photon temperature ratio of
$(T_{\nu}/T_{\gamma})_{0}\simeq(4/11)^{1/3}\simeq 1.4$ 111This result is under
the assumption of instantaneous neutrino decoupling.. Soon after neutrino
decoupling, at $\sim$0.7 MeV, the weak interactions could no longer sustain
neutron–proton chemical equilibrium; the neutron number density began to
plummet due to neutron decay, up until the onset of BBN. Once the universe
cooled down enough for the formation of helium–4, deuterium, helium–3, and
lithium to take place, the leftover free neutrons were captured (and
preserved) within light nuclei. The time elapsed until the universe expanded
and cooled enough for nuclei to form is thus the main parameter determining
$Y_{\mathrm{p}}$ and $Y_{\mathrm{D}}$.

The rate of expansion is in turn controlled by the content of the universe,
including the amount of DM. If DM is heavy enough, it is non-relativistic
during BBN and thus contributes negligibly to the overall energy density at
the time of BBN. However, if DM is as light as $\sim$20 MeV or less, it may
transition to being non-relativistic during BBN, affecting this process. In
this case, DM can no longer be efficiently produced within the thermal bath
during BBN, and its annihilation into other species takes place, transferring
entropy into the rest of the plasma. Depending on the value of $m_{\chi}$,
spin statistic of the specific DM species, and whether it couples to photons,
electron–positrons, or neutrinos (the only other radiation species around at
the time), the entropy transfer can alter the rate of expansion (and cooling)
by different amount, compared to the standard cosmological scenario.
Regardless of the modeling specifics, one of the main effects of the presence
of light thermal–relic particles is to modify the time at which
proton–to–neutron conversion and various other nuclear processes freeze
out—and therefore change the values of $Y_{\mathrm{p}}$ and $Y_{\mathrm{D}}$
as compared to the standard BBN (SBBN); for a detailed review of the BBN
calculation, see Ref. [41].

Figure 1: Primordial abundance of helium-4, as a function of DM mass
$m_{\chi}$, obtained using AlterBBN code. We assume only the standard neutrino
species ($\Delta N_{\nu}=0$) and use $\Omega_{b}h^{2}=0.0224$ [5]. The
abundances are shown for two sets of models: with DM coupling
electromagnetically (top panel) and to the SM neutrinos (bottom panel). The
dashed black line is for the standard BBN scenario, where DM has a large mass
$m_{\chi}\gtrsim 20$ MeV. The grey bands show the current measurement
uncertainty on $Y_{\mathrm{p}}$ from Ref. [34].

Figure 2: Same as Fig. 1, except for primordial abundance of deuterium
$Y_{\mathrm{D}}=(D/H)\times 10^{5}$.

In the presence of light thermal–relic DM, the primordial abundance of light
elements, including $Y_{\mathrm{p}}$ and $Y_{\mathrm{D}}$, depends on three
cosmological parameters: baryon–to–photon number density ratio $\eta$ (or,
equivalently, the present-day baryon energy density $\Omega_{b}h^{2}$), the
effective number of new relativistic species beyond the SM neutrinos, $\Delta
N_{\nu}$ (where $\Delta N_{\nu}=0$ in SBBN), and DM mass $m_{\chi}$. To
quantify the impact of DM mass on the synthesis of the light elements, we use
the publicly-available AlterBBN code [42, 43]. AlterBBN enables high-accuracy
predictions for the abundance of light elements generated during BBN in
different cosmological scenarios, including the presence of light DM
particles222Apart from AlterBBN, there are other public BBN codes, such as the
Fortran program PArthENoPE [44, 45] and the Mathematica code PRIMAT [46]. We
have verified that the specific choice of BBN code does not significantly
affect our reported results.. In Figs. 1 and 2, we show the predicted
primordial abundances in cosmologies with light thermally-coupled DM, as a
function of $m_{\chi}$, for the two DM scenarios considered in this work,
obtained using AlterBBN; these Figures also show the current best bounds on
these quantities, obtained from Lyman-$\alpha$ forest measurements and
presented in Ref. [34].

From Fig. 1, we see that DM with $m_{\chi}\gtrsim 20$ MeV has no significant
effect on $Y_{\mathrm{p}}$, because it annihilates before the weak–interaction
freezeout. In contrast, lower–mass DM annihilates after neutrino decoupling
and can significantly alter the BBN yields. Beyond this global trend, there
are additional subtleties captured in the BBN calculation shown here. For
example, in the electromagnetically coupled DM scenario, the expansion rate
can slow down for a narrow range of masses, increasing the time available for
neutron decay, but also leading to a slower conversion rate between neutrons
and protons, due to the relative decrease in neutrino energy density, as
compared to photons [26, 20]. These two effects nearly cancel out, leading to
a very slightly reduced relic abundance in the range $2\text{ MeV}\lesssim
m_{\chi}\lesssim 20$ MeV, as compared to SBBN (as shown in the top panel of
Fig. 1). Lighter DM particles simply contribute to the expansion rate and
result in a higher value of $Y_{\mathrm{p}}$. In the neutrino coupled DM
scenario, neutrino–to–photon temperature ratio is increased compared to SBBN,
speeding up the rate of weak interactions, and driving up the value of
$Y_{\mathrm{p}}$ at all masses below $\sim$20 MeV [21].

From Fig. 2, we see that $Y_{\mathrm{D}}$ is altered relative to SBBN only for
light DM with $m_{\chi}\lesssim 20$ MeV. The electromagnetically coupled DM
particles reduce the neutrino–to–photon temperature ratio, resulting in slower
expansion rate at a fixed temperature $T_{\gamma}$, leaving more time for
deuterium destruction, and reducing its relic abundance. Conversely, the
increased neutrino–to–photon temperature ratio in the neutrino coupled DM
scenario leads to a higher value of $Y_{\mathrm{D}}$.

##  III Effects on CMB

Light thermal–relic particles can leave imprints on the CMB through two
distinct effects. First, as discussed in Sec. II, they can alter the process
of BBN and change $Y_{\mathrm{p}}$, in turn affecting the temperature and
polarization anisitropy. In particular, since helium–4 recombines at a lower
temperature than hydrogen, increasing $Y_{\mathrm{p}}$ (and fixing all other
cosmological parameters) reduces the number of free electrons present in the
universe at a particular fixed temperature. Thus, the mean free path of
photons increases, resulting in increased diffusion length and damping of the
small–scale anisotropy in the CMB [47]. As a result, enhancing
$Y_{\mathrm{p}}$ due to the presence of light thermal–relic DM during BBN
generically leads to a suppression of the CMB power spectra. We modify the
standard Boltzmann code CLASS [48, 49] in order to include the effects of
$m_{\chi}$, as described in more detail below; using the output of our
modified version of CLASS, we illustrate the effect on the temperature
anisotropy in Fig. 3, for a fixed total $N_{\mathrm{eff}}=3.046$.

Figure 3: The ratio of the CMB temperature power spectrum for a cosmology
featuring light electromagnetically coupled Dirac Fermion DM, with a given
mass $m_{\chi}$, as compared to the $\Lambda$CDM scenario (with all other
parameters kept at their best-fit Planck values). Suppression of power at
small angular scales occurs due to an increased primordial yield of helium in
light–DM cosmologies. Smaller particle masses lead to a more prominent
suppression, as they alter the expansion rate during BBN more significantly.
For these purposes, masses above $\gtrsim 20$ MeV behave as standard (heavy)
CDM.

In addition to the effect on $Y_{\mathrm{p}}$, light thermal relics can affect
the overall budget of radiation, quantified by the $N_{\mathrm{eff}}$
parameter, which controls the expansion rate at later epochs and affects the
CMB anisotropy. For example, DM annihilating into photons after neutrino
decoupling can heat photons relative to the decoupled neutrinos, reducing the
present–day ratio of neutrino–to–photon temperature (assuming that the
present–day temperature of the CMB is kept fixed), corresponding to a reduced
value of $N_{\mathrm{eff}}$ [19]. Conversely, DM annihilating to neutrinos
heats up neutrinos relative to photons, increasing $N_{\mathrm{eff}}$.

Figure 4: The number of relativistic degrees of freedom $N_{\mathrm{eff,w}}$
in a cosmology featuring light thermal–relic DM, as a function of $m_{\chi}$
in the electromagnetically coupled (dashed lines) and neutrino coupled (solid
lines) scenarios. Since neutrino coupled DM particles increase the temperature
of neutrinos relative to photons, they also increase the effective number of
relativistic species, while the opposite is true for DM that annihilates into
photons. The curves are obtained using AlterBBN code, assuming only the
standard three neutrino species are present.

In standard cosmology, with only the thee SM neutrino species and SBBN,
$N_{\mathrm{eff}}=3.046$. On the other hand, in light–DM cosmologies, we can
define (following Ref. [19])

$N_{\mathrm{eff}}(m_{\chi},\Delta N_{\nu})\equiv
N_{\mathrm{eff,w}}(m_{\chi})(1+\Delta N_{\nu}/3)$ (1)

to be the total number of relativistic degrees of freedom, where

$N_{\mathrm{eff,w}}\equiv
3\left[\frac{11}{4}\left(\frac{T_{\nu}}{T_{\gamma}}\right)_{0}^{3}\right]^{\frac{4}{3}}$
(2)

includes the contribution from the SM neutrinos, as well as contribution from
light DM, and the subscript zero denotes the present time. $\Delta N_{\nu}$ is
the effective number of any additional new (neutrino–like) relativistic
species that may be present in the universe. We note that $N_{\mathrm{eff,w}}$
depends on the DM mass $m_{\chi}$, but differs for different particle models,
as illustrated in Fig. 4, using AlterBBN code (for $\Delta N_{\nu}=0$). We see
that particles with $m_{\chi}\gtrsim 20$ MeV have no significant effect on
$N_{\mathrm{eff,w}}$, because their annihilation is complete prior to the
decoupling of the SM neutrinos, and therefore does not affect the
neutrino–to–photon temperature, giving $N_{\mathrm{eff,w}}=3.046$.
Interestingly, presence of lighter DM particles coupled to photons can conceal
a possible existence of additional relativistic species, by offsetting their
effect on $\Delta N_{\nu}$; we consider this possibility in Sec. IV.

In order to account for the effects of $m_{\chi}$ in cosmologies with
thermal–relic DM, both through its effects on the BBN yields and the effects
on $N_{\mathrm{eff}}$, we modify CLASS code [48, 49]. In the standard version
of this code, the value of $Y_{\mathrm{p}}$ assumes a standard BBN scenario,
where $Y_{\mathrm{p}}$ depends only on $\Omega_{b}h^{2}$ and $\Delta N_{\nu}$
(for a fixed present-day CMB temperature). In the modified version of CLASS,
we use an updated prediction for $Y_{\mathrm{p}}$, as a function of
$\Omega_{b}h^{2}$, $\Delta N_{\nu}$, and $m_{\chi}$, obtained from AlterBBN.
Fig. 3 illustrates the output of our modified CLASS code, including both
effects discussed in this Section (note that we fix the effective number of
neutrino at $N_{\mathrm{eff}}=3.046$ in this Figure).

##  IV Current constraints

We consider the most recent CMB measurements from ACT, SPT, and Planck, as
well measurements of primordial element abundances from spectroscopic
measurements of high–redshift absorption systems, in order to infer
$m_{\chi}$. We now describe the data sets, our analysis method, and the
constraints we obtain.

###  IV.1 Data

ACT DR4: We use the multifrequency TT, EE, and TE power spectra from ACT Data
Release 4 (DR4) [31, 32], implemented within the actpollite_dr4 likelihood
333https://github.com/ACTCollaboration/pyactlike. These data products are
derived from a four–year survey, with the power spectrum measurements
reconstructed from the deepest 5400 $\text{deg}^{2}$ of the sky. This data
provide high resolution measurements of the polarization anisotropy,
complementing the data from Planck. In the ACT likelihood, the covariance of
the foreground–marginalized CMB power spectra already includes the effects of
noise, foreground uncertainty, beam and calibration uncertainties, with one
nuisance parameter $y_{p}$ included to marginalize over an overall
polarization efficiency; we allow this variable to vary in a range centred at
1. ACT alone cannot constrain the optical depth to reionization $\tau$, as it
is mainly determined by low–$\ell$ polarization power spectra; thus, when
analyzing ACT data alone, we assume a Gaussian prior on $\tau=0.065\pm 0.015$,
following Ref. [31].

SPT–3G: We use the most recent publicly–available measurements of EE and TE
power spectra from the SPT–3G survey [33]. These power spectra are obtained
from the observations of a $1500\text{deg}^{2}$ region of the sky, taken over
four months (in 2018) at three frequency bands centered on 95, 150, and 220
GHz. The SPT likelihood is publicly
available444https://pole.uchicago.edu/public/data/dutcher21 and includes the
effects of the aberration due to relative motion with respect to the CMB rest
frame [50], super-s-ample lensing [51], polarized foregrounds, uncertainty in
the calibration of the bandpowers, and uncertainty in the beam measurements.
The priors on many of these terms are listed in the Table III of Ref. [33].
Similar to ACT, SPT has the advantage of high resolution and high sensitivity
specifically in polarization measurements on small scales. Since the
low–$\ell$ polarization are not probed by this experiment, we adopt a
Planck–based Gaussian prior of $\tau=0.0543\pm 0.007$555Widening the prior to
$\tau=0.065\pm 0.015$ has no significant effect on cosmological parameter
constraints [5] when analyzing SPT data alone, following Ref. [33].

Planck 2018: We use both low–$\ell$ and high–$\ell$ multifrequency power
spectra TT, TE, and EE from Planck PR3 (2018) [5], available through the
Legacy Archive for Microwave Background Data Analysis
(LAMBDA666https://lambda.gsfc.nasa.gov/). We rely only on the lite (pre-
marginalized) likelihood for high-$\ell$ TTTEEE, for computational efficiency;
using the full likelihood does not lead to appreciable changes in the inferred
parameter values for models we consider here [5]. Planck data provides lower
noise level at large scales, as compared to the ground-based observations,
enabling complementary constraints on cosmological parameters.

Planck+ACT+SPT: When combining ACT DR4 and SPT–3G with Planck 2018 data, the
covariance between SPT and Planck is expected to be negligible since the SPT
observation region is only a small fraction of the Planck field. Following
Ref. [31], we minimize covariance between ACT and Planck by discarding large
scales ($\ell<1800$) in TT (and with no cuts in TE and EE). The correlations
between SPT and ACT can be ignored because they observe different parts of the
sky. Note that we rely on Planck low-$\ell$ EE to constrain $\tau$ in this
case.

Primordial element abundances: In addition to the CMB anisotropy, we also
consider the measurements of the primordial helium–4 and deuterium abundances,
derived from Lyman–$\alpha$ forest observations and reported in Ref. [34],
where $Y_{\mathrm{p}}=0.245\pm{0.003}$ at 68% CL, and the bound on deuterium
is $Y_{\mathrm{D}}=2.547\pm{0.025}$. Besides the observational uncertainties,
we also account for theoretical uncertainty related to uncertainties on
various nuclear reaction rates,
$\sigma_{Y_{\mathrm{p}}^{\mathrm{th}}}=1.6\times 10^{-4}$ and
$\sigma_{Y_{\mathrm{D}}^{\mathrm{th}}}=0.03$ [41, 43], as discussed in the
following.

###  IV.2 Method

The first step in our analysis is to determine the values of
$N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ that are consistent with CMB data,
without reference to $m_{\chi}$. To do this, we carry out a series of Markov
Chain Monte Carlo (MCMC) runs within the Cobaya sampling framework [52, 53],
individually and jointly applied to ACT, SPT, and Planck data. We utilize mcmc
sampler, and employ the convergence criterion $R-1=0.01$, where R is the
Gelman–Rubin threshold [54]. In each MCMC run, we sample the posterior
distributions of the six standard cosmological parameters (baryon density
$\Omega_{b}h^{2}$, DM density $\Omega_{c}h^{2}$, acoustic scale $100\theta$,
reionization optical depth $\tau$, scalar spectral index $n_{s}$, and
amplitude of the scalar perturbations $A_{s}$), with the addition of
$Y_{\mathrm{p}}$ and $N_{\mathrm{eff}}$. For each parameter, we employ broad
priors listed in Table 1. Importantly, we allow $Y_{\mathrm{p}}$ to vary in an
unconstrained way, in order to determine values that are consistent with CMB
observations, regardless of BBN predictions. Results of the MCMC runs are
shown as 68% and 95% CL contours in Fig. 5, and are consistent with previous
analyses [5, 31, 55].

Table 1: Priors on cosmological parameters used in the likelihood analysis. Parameter |  Prior
---|---
$\Omega_{b}h^{2}$ |  [0.005, 0.1]
$\Omega_{c}h^{2}$ |  [0.001, 0.99]
$100\theta$ |  [0.5, 10]
$\tau$ 777For ACT DR4 data alone, we use a Gaussian prior of $\tau=0.065\pm 0.015$. For SPT–3G data alone, we use a Gaussian prior of $\tau=0.0543\pm 0.007$. |  [0.01, 0.8]
$n_{s}$ |  [0.8, 1.2]
$\text{log}(10^{10}A_{s})$ |  [1.61, 3.91]
$Y_{\mathrm{p}}$ |  [0.1, 0.5]
$N_{\mathrm{eff}}$ |  [1, 5]

Once we obtain the full posterior probability distributions using CMB data,
the second step we take is to relate the preferred values of
$\Omega_{b}h^{2}$, $N_{\mathrm{eff}}$, and $Y_{\mathrm{p}}$ to the
corresponding allowed values of $m_{\chi}$ and $\Delta N_{\nu}$ and obtain
their bounds under each DM model. To do this, we approximate our posterior
distributions as Gaussian and employ the following chi–squared statistic to
determine the relevant CL intervals

$\chi_{\rm{CMB}}^{2}=\left(\textbf{X}-\textbf{X}_{\rm{obs}}\right)\textbf{Cov}^{-1}\left(\textbf{X}-\textbf{X}_{\rm{obs}}\right)^{T},$
(3)

where $\textbf{X}=\\{\Omega_{b}h^{2},N_{\mathrm{eff}},Y_{\mathrm{p}}\\}$ is a
vector of three relevant cosmological parameters; $\textbf{X}_{\rm{obs}}$ is
the inferred maximum–likelihood value of the three parameters; and Cov is the
covariance matrix in the 3d parameter space, derived from the sample chains.
Using AlterBBN, we recover the mapping $N_{\mathrm{eff}}(m_{\chi},\Delta
N_{\nu})$ and $Y_{\mathrm{p}}(m_{\chi},\Delta N_{\nu},\Omega_{b}h^{2})$ for
each DM model of interest. Finally, using this mapping, we identify the
appropriate values of $\chi_{\rm{CMB}}^{2}$ with 95% CL intervals for
$m_{\chi}$ and $\Delta N_{\nu}$.

To add the primordial abundance measurements from Lyman–$\alpha$ forest
spectra requires additionally imposing the constraints on $Y_{\mathrm{p}}$ and
$Y_{\mathrm{D}}$. To do this, we compute the primordial–abundance chi–squared
statistic $\chi_{\rm{PA}}^{2}$, as

$\chi_{\rm{PA}}^{2}=\frac{\left[Y_{\mathrm{p}}-Y_{\mathrm{p}}^{\rm
obs}\right]^{2}}{{\sigma_{Y_{\mathrm{p}}^{\rm
th}}}^{2}+{\sigma_{Y_{\mathrm{p}}^{\rm
obs}}}^{2}}+\frac{\left[Y_{\mathrm{D}}-Y_{\mathrm{D}}^{\rm
obs}\right]^{2}}{{\sigma_{Y_{\mathrm{D}}^{\rm
th}}}^{2}+{\sigma_{Y_{\mathrm{D}}^{\rm obs}}}^{2}},$ (4)

which includes both observational and theoretical uncertainties. To infer
combined bounds from CMB and primordial–abundance measurements, we evaluate
the sum of the relevant statistics,
$\chi_{\rm{tot}}^{2}\equiv\chi_{\rm{CMB}}^{2}+\chi_{\rm{PA}}^{2}$.

###  IV.3 CMB–only bounds

Figure 5: 95 % CL bounds are shown as lighter colored areas in the
$N_{\mathrm{eff}}$–$Y_{\mathrm{p}}$ plane, as derived from 2d marginalized
posterior distribution, obtained from different CMB data sets, denoted in the
legend. The horizontal band represents the helium–4 abundance measurement from
Ref. [34]. Theoretical predictions are shown as solid and dashed lines, for
the electromagnetically and neutrino coupled DM, respectively, both for
$\Delta N_{\nu}=0$ (no other relativistic species are allowed); DM spin
statistic is denoted in the legend. All theoretical curves converge to a
single point that corresponds to standard BBN with $m_{\chi}\gtrsim 20$ MeV.
Along each theoretical curve, $m_{\chi}$ varies from the this point to
$\lesssim$0.01 MeV at the upper left (solid) and upper right (dashed) end of
the lines.

Figure 6: The constraints on $N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ derived
from CMB data are the same as in Fig. 5. The theoretical curves are shown for
the Dirac Fermion DM that couples electromagnetically (top panel) or to
neutrinos (bottom panel). Theoretical predictions are shown as black lines for
different values of $\Delta N_{\nu}$. Along each line, $m_{\chi}$ varies from
$\gtrsim$20 MeV to $\lesssim$0.01 MeV, as labeled in the plots.

Figs. 5 and 6 show CMB–only bounds in the $N_{\mathrm{eff}}$–$Y_{\mathrm{p}}$
plane, with darker and lighter shaded areas representing 68% and 95% CL
regions, respectively. The bounds are obtained by varying the full set of
eight cosmological parameters shown in Table 1. The $N_{\mathrm{eff}}$
best–fit value is higher for SPT data than for Planck, which in turn is higher
than that for ACT. However, at the 2$\sigma$ level, all three data sets are
consistent with each other. We note that the bounds in these plots are not
dependent on DM model, and no DM parameters are fit to the data at this stage;
rather, the bounds on $N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ are used to
obtain posterior distribution for $m_{\chi}$ in the second stage of the
analysis. For illustration, the Figures also show theoretical relationship
between $m_{\chi}$, $\Delta N_{\nu}$, and $Y_{\mathrm{p}}$; this relationship
is illustrated in Figs. 1 and 4, and Eq. (1) (for $\Omega_{b}h^{2}=0.02237$).

In Fig. 5, we assume no other relativistic degrees of freedom and fix $\Delta
N_{\nu}=0$. We then illustrate the theoretical relationship between
$N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ for DM coupling either
electromagnetically (solid lines) or to SM neutrinos (dashed lines), while the
spin statistic of DM particles is denoted in the legend. Mass $m_{\chi}$
varies along each curve from $\lesssim$0.01 MeV (at the far left and far right
ends of the theoretical lines) to $\gtrsim$20 MeV (at the center of the panel,
where the dashed and the solid lines converge to a point). In the former, the
effects on $N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ saturate at $\sim$0.01 MeV
since lighter particles annihilate after the end of BBN, while the latter case
represents the SBBN scenario where DM is too massive to affect BBN through
either effect.

In Fig. 6, we focus on one example of DM spin statistic (for a Dirac Fermion)
and explore the effect of $\Delta N_{\nu}$. The two panels correspond to
electromagnetically coupled (top panel) and neutrino coupled (bottom panel)
particles. The black curves show the theoretical relationship between
$N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ for several values of $\Delta
N_{\nu}$, with $m_{\chi}$ varying along each curve from $\lesssim$0.01 MeV (at
the red dashed line) to $\gtrsim$20 MeV (at the green dashed line). We note
that any new relativistic degrees of freedom can only contribute positive
values of $\Delta N_{\nu}$, while negative values are included here for
illustration purposes only, and do not correspond to known physical models.

In Figs. 5 and 6, the allowed values of $m_{\chi}$ correspond to the parameter
space where the theoretical curve (for a given model and a given choice of
$\Delta N_{\nu}$) is consistent with the shaded region of the parameter space,
preferred by the CMB data; the values of $m_{\chi}$ outside these regions are
excluded by our CMB analyses. Precise bounds on $m_{\chi}$ obtained using our
$\chi^{2}$ procedure are listed in Table 2, both for the case where $\Delta
N_{\nu}=0$ is fixed, and for the case where $\Delta N_{\nu}$ is varied as a
free parameter (and marginalized over).

Table 2: 95% CL lower bounds on DM mass $m_{\chi}$ (in units of MeV), inferred from CMB data only, both in the absence of other relativistic species ($\Delta N_{\nu}=0$) and allowing for their presence ($\Delta N_{\nu}\geq 0$, where $\Delta N_{\nu}$ is marginalized over). ‘$-$’ indicates that there are no lower bounds.  | $\Delta N_{\nu}=0$ | $\Delta N_{\nu}\geq 0$
---|---|---
| ACT | SPT | Planck |  Planck+ACT+SPT | ACT | SPT | Planck |  Planck+ACT+SPT
Electromagnetically coupled DM
Real scalar | - | 1.94 | - | - | - | - | - | -
Complex scalar | - | 4.99 | 4.96 | 3.61 | - | - | 0.06 | 0.21
Majorana Fermion | - | 4.91 | 4.85 | 3.35 | - | - | - | 0.18
Dirac Fermion | - | 7.72 | 7.65 | 6.45 | - | 0.21 | 0.30 | 0.45
Neutrino coupled DM
Real scalar | 6.69 | - | 6.57 | 11.26 | 6.92 | - | 6.72 | 11.38
Complex scalar | 9.24 | - | 8.81 | 13.36 | 9.49 | - | 9.14 | 13.54
Majorana Fermion | 9.21 | - | 8.78 | 13.29 | 9.42 | - | 9.06 | 13.43
Dirac Fermion | 11.53 | 1.17 | 11.02 | 15.65 | 11.72 | 1.36 | 11.39 | 15.94
Table 3: Same as Table 2, except for either primordial abundance analysis only (denoted as $Y_{\mathrm{p}}+Y_{\mathrm{D}}$), or a joint analysis of all CMB data (Planck+ACT+SPT) and the primordial abundance measurements.  | $\Delta N_{\nu}=0$ | $\Delta N_{\nu}\geq 0$
---|---|---
| $Y_{\mathrm{p}}+Y_{\mathrm{D}}$ |  $Y_{\mathrm{p}}$+CMB |  $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+CMB | $Y_{\mathrm{p}}+Y_{\mathrm{D}}$ |  $Y_{\mathrm{p}}$+CMB |  $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+CMB
Electromagnetically coupled DM
Real scalar | 0.34 | 0.46 | 4.16 | 0.34 | 0.43 | 0.62
Complex scalar | 0.46 | 3.51 | 7.94 | 0.48 | 1.78 | 4.25
Majorana Fermion | 0.46 | 3.19 | 7.76 | 0.47 | 1.52 | 4.18
Dirac Fermion | 0.72 | 6.42 | 10.99 | 0.77 | 4.17 | 8.96
Neutrino coupled DM
Real scalar | 1.12 | 11.76 | 8.06 | 1.25 | 13.91 | 9.24
Complex scalar | 2.85 | 13.64 | 11.03 | 4.56 | 15.81 | 12.92
Majorana Fermion | 2.53 | 13.50 | 10.79 | 3.92 | 15.31 | 12.34
Dirac Fermion | 4.81 | 16.14 | 13.63 | 5.29 | 18.34 | 15.09
Table 4: The projected lower bounds on DM mass $m_{\chi}$ at 95% CL, in units of MeV, for Simons Observatory (SO) and CMB-S4 (S4) (both including data from Planck), and their combination with the current primordial abundance measurements. Results are shown for the absence of other relativistic species ($\Delta N_{\nu}=0$) and allowing for their presence ($\Delta N_{\nu}\geq 0$, where $\Delta N_{\nu}$ is marginalized over). We set SBBN values of $N_{\mathrm{eff}}$ and $Y_{p}$ as the fiducial models in this analysis.  | $\Delta N_{\nu}=0$ | $\Delta N_{\nu}\geq 0$
---|---|---
| SO | S4 |  $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+SO |  $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+S4 | SO | S4 |  $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+SO |  $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+S4
Electromagnetically coupled DM
Real scalar | 9.28 | 10.02 | 9.75 | 10.32 | 0.76 | 0.97 | 2.89 | 3.79
Complex scalar | 11.72 | 12.33 | 12.40 | 12.94 | 1.13 | 2.10 | 7.61 | 7.94
Majorana Fermion | 11.63 | 12.29 | 12.34 | 12.84 | 1.07 | 1.97 | 7.44 | 7.69
Dirac Fermion | 13.70 | 14.42 | 14.59 | 15.01 | 2.39 | 5.29 | 10.73 | 10.98
Neutrino coupled DM
Real scalar | 9.74 | 10.55 | 8.05 | 9.86 | 10.39 | 11.23 | 10.62 | 11.81
Complex scalar | 12.06 | 12.83 | 11.01 | 12.14 | 12.71 | 13.44 | 13.12 | 13.95
Majorana Fermion | 12.02 | 12.74 | 10.76 | 12.08 | 12.65 | 13.34 | 12.95 | 13.81
Dirac Fermion | 14.27 | 15.03 | 13.06 | 14.01 | 14.91 | 15.80 | 15.77 | 16.62

We note that the CMB bound on $m_{\chi}$ is primarily driven by the
measurement of $N_{\mathrm{eff}}$. In turn, the $N_{\mathrm{eff}}$ is
primarily constrained by Planck, with a small contribution from the current
ACT and SPT data. However, the small shifts in the values of
$N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$, preferred by different CMB
experiments, can notably affect the resulting CMB bounds on $m_{\chi}$.

For example, for neutrino coupled DM, ACT has a similar constraining power as
Planck, and the combination of ACT, SPT, and Planck strengthens the Planck
bounds by a factor of $1.4$–$1.8$, depending on the DM spin statistic. This
leads to the lower limits between 11 MeV and 16 MeV, at 95% CL, reported in
Table 2. Conversely, SPT has a similar constraining power as Planck for
electromagnetically coupled DM, but the addition of ACT and SPT weakens the
bound by up to 30% for complex scalars and Fermions. For the case of a real
scalar, neither Planck nor the combination of all CMB data can rule out any DM
masses, and SPT only rules out masses below $\sim$2 MeV.

Assuming the presence of additional relativistic degrees of freedom allows for
$\Delta N_{\nu}>0$; this additional freedom can offset the reduction in the
value of $N_{\mathrm{eff}}$ caused by electromagnetically coupled DM, putting
small $m_{\chi}$ in agreement with CMB measurements (as discussed in Sec. II).
When we fit for $\Delta N_{\nu}>0$ as a free parameter in our CMB likelihood
analysis, Planck, ACT, SPT and their combination are consistent with a far
broader range of DM masses, down to hundreds of keV. Finally, letting $\Delta
N_{\nu}>0$ vary as a free parameter generally has comparably smaller effect on
the inference of neutrino–coupled–DM mass. Details are listed in Table 2.

###  IV.4 Bounds from primordial abundances only

We next consider only the measurements of primordial helium–4 and deuterium
abundances, inferred from the spectra of absorption systems from Ref. [34],
and shown as grey bands in Figs. 1 and 2. Employing the chi–squared method
described earlier in this Section, we obtain results presented in Table 3
(columns denoted as “$Y_{\mathrm{p}}+Y_{\mathrm{D}}$”). We find that the
primordial abundance measurements alone typically place a significantly weaker
lower bound on $m_{\chi}$ than the CMB data, giving $m_{\chi}\gtrsim 0.3$ MeV,
for electromagnetically coupled DM, and $m_{\chi}\gtrsim 1.1$ MeV, for
neutrino coupled DM. Letting $\Delta N_{\nu}>0$ vary as a free parameter
generally has a negligible effect on the primordial–abundance bounds on
electromagnetically coupled DM mass, while the bounds on neutrino coupled DM
mass are improved.

We further find a small difference in the preferred value of the baryon
density $\Omega_{b}h^{2}$ between the primordial abundance measurements and
CMB results, shown in Fig. 7, consistent with the results of Ref. [25, 56].
Specifically, measurements of primordial element abundances prefer a lower
value of $\Omega_{b}h^{2}$, driven by the measurement of $Y_{\mathrm{D}}$
(since $Y_{\mathrm{p}}$ has very little sensitivity to $\Omega_{b}h^{2}$).
$\Omega_{b}h^{2}$ is degenerate with low values of $m_{\chi}$, which will
affect the mass bounds reported here. However, since the values of
$\Omega_{b}h^{2}$ measured from the CMB and the primordial abundances are
consistent at the $3\sigma$ level, we proceed to combine all the available
data in the following.

Figure 7: The $68\%$ and $95\%$ CL contours on the DM mass and baryon density,
for electromagnetically coupled (upper panel) and neutrino coupled (lower
panel) Dirac Fermion. The contours are derived using primordial abundance
measurements of helium–4 and deuterium only (green), CMB data only (gray;
including the measurements from Planck, ACT, and SPT), and from a joint
analysis of the CMB and primordial abundance data (brown).

###  IV.5 Joint bounds

Figure 8: The $68\%$ and $95\%$ CL contours on the DM mass and the effective
number of additional relativistic species, for electromagnetically coupled
(upper panel) and neutrino coupled (lower panel) Dirac Fermion DM. These
results are derived using primordial abundance measurements of helium–4 and
deuterium only (green), CMB data only (grey; including the measurements from
Planck, ACT, and SPT), and using joint analyses of CMB data and the primordial
abundance measurements (with helium–4 only: blue, with both helium–4 and
deuterium: brown). Note that $\Omega_{b}h^{2}$ is marginalized over.

Figure 9: The current and future $68\%$ and $95\%$ CL contours on the DM mass
and the effective number of additional relativistic species, for
electromagnetically coupled (upper panel) and neutrino coupled (lower panel)
Dirac Fermion DM. The contours are derived assuming joint analyses of current
primordial abundance measurements with CMB data from Planck (blue), Simons
Observatory (grey), and CMB-S4 (red), allowing for the presence of equivalent
neutrinos $\Delta N_{\nu}\geq 0$. Note that $\Omega_{b}h^{2}$ is marginalized
over.

When no additional relativistic species are allowed, the addition of
$Y_{\mathrm{p}}$ measurement to the combination of the CMB data excludes
masses below $\sim$0.5 MeV for electromagnetically coupled real scalar, which
is unconstrained by CMB alone. For all other models where the mass bounds from
CMB are present, the addition of $Y_{\mathrm{p}}$ leads to only a small
improvement. The joint CMB + $Y_{\mathrm{p}}$ limits on $m_{\chi}$ for all DM
models are presented in Table 3.

We next combine the measurement of $Y_{\mathrm{D}}$ with $Y_{\mathrm{p}}$ \+
CMB data and show the resulting 95% lower limits on $m_{\chi}$ in Table 3. If
$\Delta N_{\nu}=0$ is kept fixed, we see that the mass bound improves for
electromagnetically coupled DM, excluding masses below $\sim$4 MeV at 95% CL,
for all DM spin statistics; conversely, the bound relaxes for neutrino coupled
DM, excluding masses below $\sim$8 MeV at 95% CL, for all DM spin statistics.

Allowing $\Delta N_{\nu}$ to vary as a free parameter leads to the bounds
illustrated in Fig. 8 and in the rest of Table 3. From the Figure, we see that
there is a degeneracy between $m_{\chi}$ and $\Delta N_{\nu}$ at
$m_{\chi}\lesssim 20$ MeV, for electromagnetically coupled particles,
regardless of the spin statistic. This means that smaller masses are
consistent with the data, if new relativistic species are allowed in the
cosmological model. In contrast, for neutrino coupled DM, the mass bounds are
not very sensitive to the choice of $\Delta N_{\nu}$.

##  V Future constraints

We now quantify the degree to which measurements of small–scale CMB
polarization and temperature anisotropy can improve current bounds on
$m_{\chi}$. For this purpose, we consider two future ground–based CMB
experiments: Simons Observatory (SO) [35] and CMB-S4 [36, 37]. We next use the
forecasted covariance matrices for the three key parameters
$\\{\Omega_{b}h^{2},N_{\mathrm{eff}},Y_{\mathrm{p}}\\}$ reported in Ref. [25]
for the two CMB experiments888This reference has only focused on a subset of
models discussed in this work., computed assuming that both SO and CMB-S4 are
combined with the current Planck data. The forecasts assume a Gaussian
likelihood, such that

$\chi_{\rm{CMB}}^{2}=\left(\textbf{X}-\textbf{X}_{\mathrm{Fiducial}}\right)\textbf{Cov}^{-1}\left(\textbf{X}-\textbf{X}_{\mathrm{Fiducial}}\right)^{T},$
(5)

where we assume that the future best-fit cosmology corresponds to standard
BBN,

$\textbf{X}_{\mathrm{Fiducial}}=\\{\Omega_{b}h^{2},N_{\mathrm{eff}},Y_{\mathrm{p}}\\}=\\{0.02236,\,3.046,\,0.2472\\}$
(6)

and

$\displaystyle\textbf{Cov}=\,\left(\begin{array}[]{ccc}\sigma_{1}^{2}&\sigma_{1}\sigma_{2}\rho_{12}&\sigma_{1}\sigma_{3}\rho_{13}\\\
\sigma_{1}\sigma_{2}\rho_{12}&\sigma_{2}^{2}&\sigma_{2}\sigma_{3}\rho_{23}\\\
\sigma_{1}\sigma_{3}\rho_{13}&\sigma_{2}\sigma_{3}\rho_{23}&\sigma_{3}^{2}\\\
\end{array}\right)\,,$ (10)

where $(\sigma_{1},\sigma_{2},\sigma_{3})=(0.000073,0.11,0.0066)$ and
$(\rho_{12},\rho_{13},\rho_{23})=(0.072,0.33,-0.86)$ for SO,
$(\sigma_{1},\sigma_{2},\sigma_{3})=(0.000047,0.081,0.0043)$ and
$(\rho_{12},\rho_{13},\rho_{23})=(0.25,0.22,-0.84)$ for CMB-S4. We then
convert these covariances into 95% CL intervals for DM mass $m_{\chi}$, and
list the results in Table 4. Fig. 9 shows the projected constraints on
$m_{\chi}$ and $N_{\mathrm{eff}}$ for Dirac Fermion, when future CMB
measurements are combined with the current measurements of primordial element
abundances. We find that the addition of primordial element abundance
measurements has minimal effect on the projected bounds from future CMB data,
even in the case where $\Delta N_{\nu}$ is a free parameter.

However, if future CMB experiments return the best–fit values of
$\\{\Omega_{b}h^{2},N_{\mathrm{eff}},Y_{\mathrm{p}}\\}$ that are inconsistent
with the standard BBN, they might have sufficient accuracy to detect
thermal–relic DM mass. For example, if the best–fit parameters for SO are the
current best–fit values obtained from a combined analysis of all CMB
experiments (Planck+ACT+SPT), such measurements will have constraining power
to exclude neutrino coupled DM at $>2\sigma$ significance. Equivalently, if
thermal–relic DM has mass just below the detection threshold from the current
CMB$+Y_{\mathrm{p}}+Y_{\mathrm{D}}$ analyses (around 10 MeV for a Dirac
Fermion), future CMB experiments will be able to favor light thermal–relic DM
at 2–3$\sigma$ confidence and return measurement of its mass with a few MeV
accuracy.

##  VI Summary and discussion

Figure 10: Summary of current cosmological bounds on thermal–relic DM mass,
for the example of a Dirac Fermion, where the values for the mass consistent
with a given data set (at 95% confidence) correspond to the colored regions.
The purple shaded regions are constraints derived in Sec. IV and the forecasts
for joint analyses of $Y_{\mathrm{p}}+Y_{\mathrm{D}}$+Planck with Simons
Observatory (SO) and CMB-S4 (S4) are shown as shades of green. The light
shaded regions result from analyses that allow additional relativistic species
in the universe, and dark shaded regions are from analyses that assume no
additional light species. The top and bottom panels represent the case of
electromagnetically coupled and neutrino coupled Dirac Fermion DM,
respectively.

We have combined Planck 2018 data with the most recent public releases of
ground–based measurements of the small–scale CMB temperature and polarization
anisotropy from ACT DR4 and SPT–3G, and with the measurement of the primordial
abundances of helium–4, $Y_{\mathrm{p}}$, and deuterium, $Y_{\mathrm{D}}$. We
derived the most general bounds on the mass of a light thermally–coupled DM
particle that couples to either SM neutrinos or photons, affecting the
expansion history and the process of BBN. To quantify the effects of mass
$m_{\chi}$ on CMB anisotropy, our likelihood analyses have included
$N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ as free parameters, in addition to the
standard cosmological parameters. We converted the resulting posterior
probability distribution on the relevant cosmological parameters to the
posterior probability distribution on $m_{\chi}$, inferring the lower mass
limits reported in Tables 2 and 3. We also explored DM mass measurements that
might be possible to achieve in the future, with SO and CMB-S4, with key
results reported in Table 4. Importantly, we quantified how the presence of
additional relativistic species may affect the inference of DM mass, and found
a moderate impact on electromagnetically coupled DM models, and a negligible
impact on neutrino coupled models. A summary of the bounds obtained from
various combinations of data is illustrated in Fig. 10, for the representative
case of a Dirac Fermion.

###  VI.1 Key results

For DM coupled to neutrinos, combining ACT and SPT with Planck improves the
lower limit on $m_{\chi}$ by 40%–80%, depending on DM spin statistic. The
strongest lower limit on the mass is around 15 MeV, at 95%CL, regardless of
whether $\Delta N_{\nu}$ is fixed to 0 or allowed to vary. The CMB is
significantly more constraining than the $Y_{\mathrm{p}}$ and $Y_{\mathrm{D}}$
measurements alone.

For electromagnetically coupled DM, the addition of ACT and SPT to Planck
weakens previous bounds by $\lesssim 30$% for $\Delta N_{\nu}=0$. This occurs
mainly because the joint CMB data prefer a lower helium abundance
$Y_{\mathrm{p}}$ and a slightly higher value of $N_{\mathrm{eff}}$ than Planck
alone. We emphasize that these shifts have small statistical significance, and
the standard BBN is still within 2$\sigma$ of the best–fit value for the
combined CMB measurements. Allowing additional light species also
significantly weakens CMB bounds on electromagnetically coupled models,
resulting in the change of the lower limit on the mass from a few MeV to
hundreds of keV. In this case, the mass of a Majorana Fermion is only
constrained by the CMB when all three CMB data sets are combined. Finally,
when $\Delta N_{\nu}>0$, the primordial abundance measurements and the CMB are
highly complementary in their constraining power, and the mass bounds are
improved by a factor of a few when the two sets of data are combined; the
lower limit on the mass is between $\sim 600$ keV and $\sim 9$ MeV, depending
on DM spin statistic.

###  VI.2 Robustness of the current measurements

We highlight the small differences in the preferred values of cosmological
parameters found between different data sets considered in this study, and
discuss their effects on the current bounds on thermal–relic DM. First, as
discussed in Sec. IV, the addition of $Y_{\mathrm{D}}$ measurement to the CMB
worsens most of the CMB–only bounds on neutrino coupled DM because of the mild
tension between the $\Omega_{b}h^{2}$ values inferred from the two data sets
[25, 56]. This tension is currently at $\sim 1.8\sigma$ level, but future
high–precision spectral characterization of damped Lyman–$\alpha$ systems at
high redshift can lead to further improvement in the measurement of
$Y_{\mathrm{D}}$ [57], while the next-generation CMB measurements from SO and
CMB–S4 will significantly decrease the uncertainty on the measurement of the
baryon fraction as well. If the tension persists and increases in significance
with future data, neutrino coupled DM model will be insufficient to describe
all data, and new physics may be required to explain this discrepancy and
infer bounds on DM mass accurately.

Secondly, as discussed in Sec. IV, the mass bounds for all DM models shift
when ground–based CMB measurements are combined with Planck, because different
CMB experiments prefer slightly different values of $N_{\mathrm{eff}}$.
Notably, ACT is only consistent with the standard BBN within 2$\sigma$,
preferring a slightly lower $N_{\mathrm{eff}}$, as discussed in Ref. [31]; at
the same time, SPT prefers a larger value than the Planck best–fit. This mild
discrepancy can be a statistical fluctuation, or a consequence of a systematic
effect. Preliminary checks we performed in our analysis do not point directly
to either, and a more detailed investigation is needed to understand this
difference. Notably, the preference towards lower value seen in ACT is not
obviously driven by low multipoles in EE and TE power spectra that drive
preference for other beyond–$\Lambda$CDM models [58]. Future data sets from
ACT and those from SO and CMB–S4 will be able shed more light on this issue
and recover a more robust bound on DM mass.

###  VI.3 The potential of future CMB measurements

Key science goals of Simons Observatory and CMB-S4 include high–precision
measurements of $N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$. As discussed in Sec.
V, if these results confirm the standard BBN scenario, the bounds on
thermal–relic DM will improve by several MeV, as compared to current analyses
for most models999Due to the shift in the preferred values of
$N_{\mathrm{eff}}$ and $Y_{p}$, the combination of all CMB data can provide
tighter constraints than the forecast bounds for neutrino coupled DM in some
scenarios., excluding most of the mass range that could affect the process of
BBN.

Alternatively, in the scenario where future CMB measurements converge around
the current best–fit values of $N_{\mathrm{eff}}$ and $Y_{\mathrm{p}}$ from
Planck+ACT+SPT, both SO and CMB–S4 will have sufficient sensitivity to exclude
standard BBN at several $\sigma$; if the DM mass is just above the current
upper limit, these experiments will be able to measure it with a few–MeV
accuracy.

In the end, it is important to note that the effects of neutrinophilic and
electrophilic DM on cosmological observables can mutually cancel, and may be
obscured by the presence of other new physics. If reality is not well
described by standard cosmology and standard BBN, many of the mass bounds
derived in this and similar studies will require revision.

## Acknowledgements

The authors acknowledge helpful conversations with Jo Dunkley, related to
testing the results of ACT analysis. RA and VG acknowledge the support from
NASA through the Astrophysics Theory Program, Award Number 21-ATP21-0135. EC
acknowledges support from the STFC Ernest Rutherford Fellowship ST/M004856/2,
STFC Consolidated Grant ST/S00033X/1and from the European Research Council
(ERC) under the European Union’s Horizon 2020 research and innovation
programme (Grant agreement No. 849169). JCH acknowledges support from NSF
grant AST-2108536. The Flatiron Institute is supported by the Simons
Foundation.

## References

  * [1] F. Zwicky, “Republication of: The redshift of extragalactic nebulae,” _General Relativity and Gravitation_ , vol. 41, no. 1, pp. 207–224, Jan. 2009\.
  * [2] V. C. Rubin and J. Ford, W. Kent, “Rotation of the Andromeda Nebula from a Spectroscopic Survey of Emission Regions,” _ApJ_ , vol. 159, p. 379, Feb. 1970.
  * [3] A. Refregier, “Weak Gravitational Lensing by Large-Scale Structure,” _ARA &A_, vol. 41, pp. 645–668, Jan. 2003.
  * [4] D. Clowe, M. Bradač, A. H. Gonzalez, M. Markevitch, S. W. Randall, C. Jones, and D. Zaritsky, “A Direct Empirical Proof of the Existence of Dark Matter,” _ApJL_ , vol. 648, no. 2, pp. L109–L113, Sep. 2006.
  * [5] Planck Collaboration, N. Aghanim, Y. Akrami, M. Ashdown _et al._ , “Planck 2018 results. VI. Cosmological parameters,” _A &A_, vol. 641, p. A6, Sep. 2020.
  * [6] Y. B. Zel’dovich, “Magnetic Model of the Universe,” _Soviet Journal of Experimental and Theoretical Physics_ , vol. 21, p. 656, Sep. 1965.
  * [7] Y. B. Zeldovic, L. B. Okun, and S. B. Pikelner, “Quarks, astrophysical and physico-chemical aspects,” _Physics Letters_ , vol. 17, no. 2, pp. 164–166, Jul. 1965.
  * [8] H.-Y. Chiu, “Symmetry Between Particle and Antiparticle Populations in the Universe,” _PhRvL_ , vol. 17, no. 13, pp. 712–714, Sep. 1966.
  * [9] G. Bertone, D. Hooper, and J. Silk, “Particle dark matter: evidence, candidates and constraints,” _PhR_ , vol. 405, no. 5-6, pp. 279–390, Jan. 2005.
  * [10] G. Arcadi, M. Dutra, P. Ghosh, M. Lindner, Y. Mambrini, M. Pierre, S. Profumo, and F. S. Queiroz, “The waning of the WIMP? A review of models, searches, and constraints,” _European Physical Journal C_ , vol. 78, no. 3, p. 203, Mar. 2018.
  * [11] L. Roszkowski, E. M. Sessolo, and S. Trojanowski, “WIMP dark matter candidates and searches—current status and future prospects,” _Reports on Progress in Physics_ , vol. 81, no. 6, p. 066201, Jun. 2018.
  * [12] D. Hooper and K. M. Zurek, “Natural supersymmetric model with MeV dark matter,” _PhRvD_ , vol. 77, no. 8, p. 087302, Apr. 2008.
  * [13] J. L. Feng and J. Kumar, “Dark-Matter Particles without Weak-Scale Masses or Weak Interactions,” _PhRvL_ , vol. 101, no. 23, p. 231301, Dec. 2008.
  * [14] D. Hooper, M. Kaplinghat, L. E. Strigari, and K. M. Zurek, “MeV dark matter and small scale structure,” _PhRvD_ , vol. 76, no. 10, p. 103515, Nov. 2007.
  * [15] G. Jungman, M. Kamionkowski, and K. Griest, “Supersymmetric dark matter,” _PhR_ , vol. 267, pp. 195–373, Mar. 1996.
  * [16] C. B. hm, M. J. Dolan, and C. McCabe, “Increasing Neff with particles in thermal equilibrium with neutrinos,” _JCAP_ , vol. 2012, no. 12, p. 027, Dec. 2012.
  * [17] L. Bergström, “Dark matter evidence, particle physics candidates and detection methods,” _Annalen der Physik_ , vol. 524, no. 9-10, pp. 479–496, Oct. 2012.
  * [18] C. Bœhm, M. J. Dolan, and C. McCabe, “A lower bound on the mass of cold thermal dark matter from Planck,” _JCAP_ , vol. 2013, no. 8, p. 041, Aug. 2013.
  * [19] G. Steigman, “Equivalent neutrinos, light WIMPs, and the chimera of dark radiation,” _PhRvD_ , vol. 87, no. 10, p. 103517, May 2013.
  * [20] K. M. Nollett and G. Steigman, “BBN and the CMB constrain light, electromagnetically coupled WIMPs,” _PhRvD_ , vol. 89, no. 8, p. 083508, Apr. 2014.
  * [21] ——, “BBN and the CMB constrain neutrino coupled light WIMPs,” _PhRvD_ , vol. 91, no. 8, p. 083505, Apr. 2015.
  * [22] R. J. Wilkinson, A. C. Vincent, C. BÅ‘hm, and C. McCabe, “Ruling out the light weakly interacting massive particle explanation of the Galactic 511 keV line,” _PhRvD_ , vol. 94, no. 10, p. 103525, Nov. 2016\.
  * [23] M. Escudero, “Neutrino decoupling beyond the Standard Model: CMB constraints on the Dark Matter mass with a fast and precise Neff evaluation,” _JCAP_ , vol. 2019, no. 2, p. 007, Feb. 2019.
  * [24] P. F. Depta, M. Hufnagel, K. Schmidt-Hoberg, and S. Wild, “BBN constraints on the annihilation of MeV-scale dark matter,” _JCAP_ , vol. 2019, no. 4, p. 029, Apr. 2019.
  * [25] N. Sabti, J. Alvey, M. Escudero, M. Fairbairn, and D. Blas, “Refined bounds on MeV-scale thermal dark sectors from BBN and the CMB,” _JCAP_ , vol. 2020, no. 1, p. 004, Jan. 2020.
  * [26] E. W. Kolb, M. S. Turner, and T. P. Walker, “Effect of interacting particles on primordial nucleosynthesis,” _PhRvD_ , vol. 34, no. 8, pp. 2197–2205, Oct. 1986.
  * [27] P. D. Serpico and G. G. Raffelt, “MeV-mass dark matter and primordial nucleosynthesis,” _PhRvD_ , vol. 70, no. 4, p. 043526, Aug. 2004.
  * [28] C. Boehm, T. A. Enßlin, and J. Silk, “Can annihilating dark matter be lighter than a few GeVs?” _Journal of Physics G Nuclear Physics_ , vol. 30, no. 3, pp. 279–285, Mar. 2004.
  * [29] G. Steigman and K. M. Nollett, “Light WIMPs, equivalent neutrinos, BBN, and the CMB,” _MmSAI_ , vol. 85, p. 175, Jan. 2014.
  * [30] C. Giovanetti, M. Lisanti, H. Liu, and J. T. Ruderman, “Joint CMB and BBN Constraints on Light Dark Sectors with Dark Radiation,” _arXiv e-prints_ , p. arXiv:2109.03246, Sep. 2021.
  * [31] S. Aiola, E. Calabrese, L. Maurin _et al._ , “The Atacama Cosmology Telescope: DR4 maps and cosmological parameters,” _JCAP_ , vol. 2020, no. 12, p. 047, Dec. 2020.
  * [32] S. K. Choi, M. Hasselfield, S.-P. P. Ho _et al._ , “The Atacama Cosmology Telescope: a measurement of the Cosmic Microwave Background power spectra at 98 and 150 GHz,” _JCAP_ , vol. 2020, no. 12, p. 045, Dec. 2020\.
  * [33] D. Dutcher, L. Balkenhol, P. A. R. Ade _et al._ , “Measurements of the E -mode polarization and temperature-E -mode correlation of the CMB from SPT-3G 2018 data,” _PhRvD_ , vol. 104, no. 2, p. 022003, Jul. 2021.
  * [34] Particle Data Group, P. A. Zyla, , and others., “Review of Particle Physics,” _Progress of Theoretical and Experimental Physics_ , vol. 2020, no. 8, p. 083C01, Aug. 2020.
  * [35] P. Ade, J. Aguirre, Z. Ahmed _et al._ , “The Simons Observatory: science goals and forecasts,” _JCAP_ , vol. 2019, no. 2, p. 056, Feb. 2019\.
  * [36] K. N. Abazajian, P. Adshead, Z. Ahmed _et al._ , “CMB-S4 Science Book, First Edition,” _arXiv e-prints_ , p. arXiv:1610.02743, Oct. 2016\.
  * [37] K. Abazajian, G. Addison, P. Adshead _et al._ , “CMB-S4 Science Case, Reference Design, and Project Plan,” _arXiv e-prints_ , p. arXiv:1907.04473, Jul. 2019.
  * [38] K. Enqvist, K. Kainulainen, and V. Semikoz, “Neutrino annihilation in hot plasma,” _Nuclear Physics B_ , vol. 374, no. 2, pp. 392–404, Apr. 1992\.
  * [39] A. D. Dolgov, “Neutrinos in cosmology,” _PhR_ , vol. 370, no. 4-5, pp. 333–535, Nov. 2002.
  * [40] S. Hannestad, “Oscillation effects on neutrino decoupling in the early universe,” _PhRvD_ , vol. 65, no. 8, p. 083006, Apr. 2002.
  * [41] E. S. Jenssen, “New alterbbn: A code for big bang nucleosynthesis with light dark matter,” 2016.
  * [42] A. Arbey, “AlterBBN: A program for calculating the BBN abundances of the elements in alternative cosmologies,” _Computer Physics Communications_ , vol. 183, no. 8, pp. 1822–1831, Aug. 2012.
  * [43] A. Arbey, J. Auffinger, K. P. Hickerson, and E. S. Jenssen, “AlterBBN v2: A public code for calculating Big-Bang nucleosynthesis constraints in alternative cosmologies,” _arXiv e-prints_ , p. arXiv:1806.11095, Jun. 2018\.
  * [44] O. Pisanti, A. Cirillo, S. Esposito, F. Iocco, G. Mangano, G. Miele, and P. D. Serpico, “PArthENoPE: Public algorithm evaluating the nucleosynthesis of primordial elements,” _Computer Physics Communications_ , vol. 178, no. 12, pp. 956–971, Jun. 2008.
  * [45] R. Consiglio, P. F. de Salas, G. Mangano, G. Miele, S. Pastor, and O. Pisanti, “PArthENoPE reloaded,” _Computer Physics Communications_ , vol. 233, pp. 237–242, Dec. 2018.
  * [46] C. Pitrou, A. Coc, J.-P. Uzan, and E. Vangioni, “Precision big bang nucleosynthesis with improved Helium-4 predictions,” _PhR_ , vol. 754, pp. 1–66, Sep. 2018.
  * [47] R. Trotta and S. H. Hansen, “Constraining the helium abundance with CMB data,” _PhRvD_ , vol. 69, no. 2, p. 023509, Jan. 2004.
  * [48] J. Lesgourgues, “The Cosmic Linear Anisotropy Solving System (CLASS) I: Overview,” _arXiv e-prints_ , p. arXiv:1104.2932, Apr. 2011.
  * [49] D. Blas, J. Lesgourgues, and T. Tram, “The Cosmic Linear Anisotropy Solving System (CLASS). Part II: Approximation schemes,” _JCAP_ , vol. 2011, no. 7, p. 034, Jul. 2011.
  * [50] D. Jeong, J. Chluba, L. Dai, M. Kamionkowski, and X. Wang, “Effect of aberration on partial-sky measurements of the cosmic microwave background temperature power spectrum,” _PhRvD_ , vol. 89, no. 2, p. 023003, Jan. 2014\.
  * [51] A. Manzotti, W. Hu, and A. Benoit-Lévy, “Super-sample CMB lensing,” _PhRvD_ , vol. 90, no. 2, p. 023003, Jul. 2014.
  * [52] J. Torrado and A. Lewis, “Cobaya: code for Bayesian analysis of hierarchical physical models,” _JCAP_ , vol. 2021, no. 5, p. 057, May 2021\.
  * [53] ——, “Cobaya: Bayesian analysis in cosmology,” p. ascl:1910.019, Oct. 2019\.
  * [54] A. Gelman and D. B. Rubin, “Inference from Iterative Simulation Using Multiple Sequences,” _Statistical Science_ , vol. 7, pp. 457–472, Jan. 1992\.
  * [55] L. Balkenhol, D. Dutcher, P. A. R. Ade _et al._ , “Constraints on $\Lambda$CDM extensions from the SPT-3G 2018 E E and T E power spectra,” _PhRvD_ , vol. 104, no. 8, p. 083509, Oct. 2021.
  * [56] C. Pitrou, A. Coc, J.-P. Uzan, and E. Vangioni, “A new tension in the cosmological model from primordial deuterium?” _MNRAS_ , vol. 502, no. 2, pp. 2474–2481, Apr. 2021.
  * [57] R. Maiolino, M. Haehnelt, M. T. Murphy _et al._ , “A Community Science Case for E-ELT HIRES,” _arXiv e-prints_ , p. arXiv:1310.3163, Oct. 2013.
  * [58] J. C. Hill, E. Calabrese, S. Aiola, and others., “The Atacama Cosmology Telescope: Constraints on Pre-Recombination Early Dark Energy,” _arXiv e-prints_ , p. arXiv:2109.04451, Sep. 2021.


# Leveraging Multi-Method Evaluation for
Multi-Stakeholder Settings

Christine Bauer  Johannes Kepler University Linz, Austria
<EMAIL_ADDRESS>and  Eva Zangerle  University of Innsbruck, Austria
<EMAIL_ADDRESS>
###### Abstract.

In this paper, we focus on recommendation settings with multiple stakeholders
with possibly varying goals and interests, and argue that a single evaluation
method or measure is not able to evaluate all relevant aspects in such a
complex setting. We reason that employing a multi-method evaluation, where
multiple evaluation methods or measures are combined and integrated, allows
for getting a richer picture and prevents blind spots in the evaluation
outcome.

recommender systems, evaluation, multi-methods, multi-stakeholder, user
experience

††conference: ImpactRS ’19: 1st Workshop on the Impact of Recommender Systems;
September 19, 2019; Copenhagen, Denmark
Copyright ©2019 for this paper by its authors. Use permitted under Creative
Commons License Attribution 4.0 International (CC BY 4.0).
††price: ††copyright: none

##  1\. Introduction

In recommender systems (RS) research, we observe a strong focus on advancing
systems such that they accurately predict items that an individual user may be
interested in. The approach of evaluating an RS is thereby largely focused on
system-centric methods and metrics (e.g., recall and precision in
leave-$n$-out analyses (Knijnenburg and Willemsen, 2015)). By employing such
an evaluation approach and aiming at optimizing these metrics, the following
crucial components in the ecosystem are neglected (Burke, 2017; Ekstrand and
Willemsen, 2016): (i) multiple stakeholders are embedded in the ecosystem, but
current research largely considers merely the role of the end consumer; (ii)
the stakeholders typically have diverging interests and objectives for an RS;
however, accurately predicting a user’s interests is the predominant focus in
current RS research; and (iii) with taking a mainly accuracy-driven, system-
centric approach to evaluation, many aspects that determine a user’s
experience with an RS are not considered (Knijnenburg and Willemsen, 2015).
This results in an incomplete picture of user experience, leaving “blind
spots” that are not captured in the quality evaluation of an RS. Although
studies (Azaria et al., 2013) could show that a lower accuracy rate may
increase the business utility (e.g., revenue) without any significant drop in
user satisfaction, the objectives and interests of stakeholders other than the
user are typically not the focus of research in academic settings in the RS
community.

In this paper, we call for considering the multiple stakeholders in RS
evaluation and postulate that only taking a multi-method evaluation approach
allows for capturing and assessing the various interests, objectives, and
experiences of these very stakeholders; thus, contributing to eliminating the
blind spots in RS evaluation.

To illustrate the opportunities of multi-method evaluations in evaluating an
RS from multiple stakeholders’ perspectives, we exemplarily focus on the music
domain. In Section 2, we outline related work. Section 3 points out the
multitude and diversity of stakeholders to be considered in the evaluation of
music RS. In Section 4, we exemplify for selected stakeholders how the
integration of multiple evaluation methods contributes to eliminating blind
spots in RS evaluation while balancing the stakeholders’ interests.

##  2\. Related Work

The idea to combine different research methods is not a new one. The concept
of mixed methods research (Creswell, 2003), for instance, combines
quantitative and qualitative research approaches. It has been termed the third
methodological paradigm, with quantitative and qualitative methods
representing the first and second paradigm respectively (Teddlie and
Tashakkori, 2009). Yet, it seems that mixed methods research appears to
attract considerable interest but is rarely brought into practice (Ågerfalk,
2013). From a practical point of view, the reasons for the low adoption of
evaluations leveraging multiple methods are manifold, including higher costs,
higher complexity, wider skill requirements compared to adopting one method
only (Celik et al., 2018).

For RS research, Gunawardana and Shani (2015) point out that there is an
extensive number of aspects that may be considered when assessing the
performance of a recommendation algorithm. Indeed, already early research on
RS pointed towards the wide variety of metrics available for system-centric RS
evaluation, including classification metrics, predictive metrics, coverage
metrics, confidence metrics, and learning rate metrics (Said et al., 2012). As
accuracy-driven evaluation has been shown not to be able to capture all the
aspects that are relevant for user satisfaction (Konstan and Riedl, 2012),
more user-relevant metrics and measures have been introduced and considered
over time (Kaminskas and Bridge, 2016) (so-called “quality factors beyond
accuracy” (Jannach et al., 2016)). This wider range of objectives includes
qualities such as novelty, serendipity (Herlocker et al., 2004), or diversity
(Kaminskas and Bridge, 2016).

Kohavi et al. (2013) stress the importance of applying multiple metrics also
in the field of A/B testing and online experiments, pointing out that
different metrics reflect different concerns. For A/B testing in RS research,
Ekstrand and Willemsen (2016) emphasize the need to include methods and
metrics that go beyond the typical A/B behavior metrics. They argue that the
currently dominating RS evaluation based on implicit feedback and A/B testing
(they refer to this combination as “behaviorism”) is often very limited in its
ability to explain why users acted in a particular way. They emphasize that
experiments need to be thoroughly grounded in theory and point to the
advantages of collecting subjective responses from users which may help to
explain their behavior.

Jannach and Adomavicius (2017) point out that academic research in the field
of RS tends to focus on the consumer’s perspective with the goal to maximize
the consumer’s utility (measured in terms of the most accurate items for a
user), while maximizing the provider’s utility (e.g., in terms of profit)
appears to be neglected. While industry research on RS will naturally build
around the provider’s perspective, publications in this area are scarce
(Zanker et al., 2019).

##  3\. Digital Music Stakeholders

Various stakeholders are involved in the digital music value chain
(Abdollahpouri and Essinger, 2017). From songwriters who create songs; to
performers (e.g., (solo) artists, bands); to music producers who take a broad
role in the production of a track; to record companies, including the three
major labels; to music platform providers with huge repositories of music
tracks, acting at the interface to music consumers; and hundreds of millions
of music consumers with different music preferences and various objectives for
using RS (e.g., discovering previously unknown items, rediscovering items not
having listened to in a while ); to society at large with its social,
economic, and political objectives and needs.

Some stakeholders focus on user experience, where the goal is to propose “the
right music, to the right user, at the right moment” (Laplante, 2014). Other
stakeholders have business-oriented utility functions (Abdollahpouri and
Essinger, 2017). For instance, artists will most likely want to have their own
songs recommended to consumers. While some artists may be fine with any of
their songs being recommended, others may prefer to increase the playcount of
a particular song (e.g., to reach the top charts, which would open an
opportunity to draw an even broader audience; or some song may generate higher
revenues than others due to contract rules). Achieving additional 1,000
playcounts will not get apparent for highly popular artists with yearly
playcounts of several billions, but could be an important milestone for a
comparatively less popular (local) artist.

##  4\. Balancing Stakeholder Interests in Evaluation

In the following, we aim to make the case for multi-method evaluations that
contribute to identifying the strong and weak spots of a music RS for the
stakeholders involved, where we focus on the users’ and artists’ perspectives
in this section.

From a user perspective, recommendations that are adequate in terms of system-
centric measures—e.g., the predictive accuracy of recommendation algorithms—do
not necessarily meet a user’s expectations (Konstan and Riedl, 2012). User-
centric evaluation methods, in contrast, involve users who interact with an RS
(Knijnenburg and Willemsen, 2015) to gather user feedback (Beel et al., 2013)
either implicitly or explicitly (depending on the concrete evaluation design).
Such methods measure a user’s perceived quality of the RS at the time of
recommendation, e.g., by established questionnaires (Pu et al., 2011). Still,
relying solely on user-centric methods does not reveal the accuracy of the
recommendations, because, given the vast amount of items, users are not able
to judge whether a given recommendation was indeed the most relevant one (Beel
et al., 2013).

Measuring accuracy does not capture the recommendations’ usefulness for users,
because higher accuracy scores do not necessarily imply higher user
satisfaction (McNee et al., 2006). For instance, a user’s most favorite song
is an accurate prediction; still, repeating the same song five times is,
though accurate, likely not a satisfying experience. Hence, we argue that for
evaluating the user’s perspective of a RS—the user being only one of the many
stakeholders involved—multiple evaluation methods and measures are required.
This may include combining a set of different measures (ranging from recall
and precision to serendipity, list diversity or novelty) or integrating
different evaluation methods (ranging from leave-n-out offline experiments to
user studies and A/B testing). Furthermore, although A/B-testing using user’s
implicit feedback is effective for testing the impact of different algorithms
or designs on user behavior—and is, thus, frequently considered the “golden
standard” for recommender evaluation—, it has limited ability in explaining
why users acted in a particular way (Ekstrand and Willemsen, 2016). Additional
information (e.g., users’ subjective responses) is necessary to allow for
explaining behavior.

In short, sticking to a single evaluation method narrows our view on the RS,
literally having blinders on, while devising and evaluating RS. We can borrow
from social and behavioral sciences, where, e.g., mixed-methods research
combines quantitative and qualitative evaluations using different designs
(Creswell, 2003). Creswell’s proposed designs include—among others—the
convergent parallel design and the sequential design. In the convergent
parallel design, two evaluation methods are first applied in parallel, and
finally integrated into a single interpretation. The sequential design uses
sequential timing, employing the methods in distinct phases. The second phase
of the study, using the second method, is designed such that it follows from
the results of the first phase. Depending on the research goal and the
concrete choice of methods, researchers may either interpret how the second
phase’s results help to explain the initial results (explanatory design) or
they build on the exploratory results of the first phase to subsequently
employ a different method (in the second phase) to test or generalize the
initial findings (exploratory design). For instance, Kamehkhosh and Jannach
(2017) showed that in the field of music RS, the results of a conducted
offline evaluation could be reproduced with online studies assessing the
users’ perceived quality of recommendations. Similarly, for the Recommender
Systems Challenge 2017, participants firstly evaluated their prototypes in
offline evaluations, before actually deploying them and evaluating them in the
live system utilizing A/B tests (Abel et al., 2017), showing that many of RS
who performed well in offline evaluations were able to repeat this in online
experiments. However, some of the devised RS also performed substantially
worse in online experiments—highlighting a shortcoming that was not revealed
by evaluating from solely an offline perspective. Along the same lines,
Ekstrand and Willemsen (Ekstrand and Willemsen, 2016) state that utilizing
behaviorism for evaluation purposes (e.g., through A/B tests) is not
sufficient to understand _why_ users act in a particular way and, for
instance, like a particular recommendation.

While academic research in the field of RS tends to focus on maximizing the
users’ utility, some authors (e.g., (Jannach and Adomavicius, 2017)) emphasize
the importance of profit (or value) maximization. Profit maximization may not
only be a goal for platform providers, but also for artists who are the
content providers for music platforms. From an artist’s perspective, a good RS
recommends the respective artist’s songs sufficiently frequently, which may
ultimately lead to playcounts, likes, purchases, profit maximization, etc.
Evaluating for profit may, though, leave blind spots. For example, depending
on the chosen strategy, an artist may want to emphasize other values such as
expanding the audience (thus, reaching new listeners) or increasing the
listening or purchase volume within the current fan base. Hence, metrics such
as number of unique listeners per artist, the sum of playcounts over all songs
of an artist, and metrics such as profit-per-audience type may be valuable for
RS optimization and need to be considered in the RS evaluation strategy.
Accordingly, evaluation efforts need to elicit and integrate the artists’
goals and preferences need to be elicited and integrated into the evaluation
efforts. While evaluation on a per-artist-basis might be interesting for the
individual artists (e.g., for a comparison between platforms and their
integrated RS), it may not be adequate for an overall RS evaluation. Still, an
RS needs to be evaluated for its ability to serve the various strategies and
for revealing potential tendencies towards the one or other strategy. As the
targeted strategy might correlate with artist characteristics (e.g., top-of-
the-top artists vs. “long tail” artist; early career vs. come-back phase vs.
long-term career; mainstream artists vs. niche genres), it might be in the
society’s interest to evaluate for and ensure a balance.

Having given these examples, we emphasize that, due to interdependencies
between the RS and the various stakeholders’ actions, the entire RS ecosystem
has to be taken into account in the evaluation. For instance, low accuracy of
recommendations and low user experience are not likely to continuously
increase profits for the platform provider and all kinds of artists; high
accuracy does not automatically imply high user experience and may not
contribute to profit maximization.

##  5\. Conclusions

In this position paper, we exemplarily focused on the digital music ecosystem
to illustrate that multiple stakeholders are impacted by music RS, and
discussed the opportunities of multi-method evaluations to consider the
multiple stakeholders’ perspectives. We emphasize that—irrespective of the
application domain—there are always multiple stakeholders involved in
recommendation settings. Hence, there are always multiple—and possibly
diverging—perspectives and goals of these very stakeholders which need to be
considered in evaluating an RS. Consequently, multiple evaluation methods and
criteria have to be combined and possibly also weighted.

_Multi-method evaluations_ allow for gathering a richer and more integrated
picture of the quality of a RS and contributes to understanding the various
phenomena involved in a multi-stakeholder setting, for which one method in
isolation would be insufficient (Venkatesh et al., 2013).

###### Acknowledgements.

This research is supported by the Austrian Science Fund (FWF): V579.

## References

  * (1)
  * Abdollahpouri and Essinger (2017) Himan Abdollahpouri and Steve Essinger. 2017.  Multiple stakeholders in music recommender systems. In _1st International Workshop on Value-Aware and Multistakeholder Recommendation at RecSys 2017_ _(VAMS ’17)_.  arXiv:1708.00120
  * Abel et al. (2017) Fabian Abel, Yashar Deldjoo, Mehdi Elahi, and Daniel Kohlsdorf. 2017.  Recsys challenge 2017: Offline and online evaluation. In _Proceedings of the Eleventh ACM Conference on Recommender Systems_. ACM, New York, NY, USA, 372–373.
  * Ågerfalk (2013) Pär J Ågerfalk. 2013\.  Embracing diversity through mixed methods research.  _European Journal of Information Systems_ 22, 3 (2013), 251–256.  https://doi.org/10.1057/ejis.2013.6
  * Azaria et al. (2013) Amos Azaria, Avinatan Hassidim, Sarit Kraus, Adi Eshkol, Ofer Weintraub, and Irit Netanely. 2013\.  Movie Recommender System for Profit Maximization. In _Proceedings of the 7th ACM Conference on Recommender Systems_ _(RecSys ’13)_. ACM, New York, NY, USA, 121–128.  https://doi.org/10.1145/2507157.2507162
  * Beel et al. (2013) Joeran Beel, Marcel Genzmehr, Stefan Langer, Andreas Nürnberger, and Bela Gipp. 2013\.  A Comparative Analysis of Offline and Online Evaluations and Discussion of Research Paper Recommender System Evaluation. In _Proceedings of the International Workshop on Reproducibility and Replication in Recommender Systems Evaluation_ _(RepSys ’13)_. ACM, New York, NY, USA, 7–14.  https://doi.org/10.1145/2532508.2532511
  * Burke (2017) Robin Burke. 2017\.  Multisided fairness for recommendation. In _4th Workshop on Fairness, Accountability, and Transparency in Machine Learning_ _(FAT/ML ’17)_.  arXiv:1707.00093
  * Celik et al. (2018) Ilknur Celik, Ilaria Torre, Frosina Koceva, Christine Bauer, Eva Zangerle, and Bart Knijnenburg. 2018\.  UMAP 2018 Intelligent User-Adapted Interfaces: Design and Multi-Modal Evaluation (IUadaptMe). In _Adjunct Publication of the 26th Conference on User Modeling, Adaptation and Personalization_ _(UMAP ’18)_. ACM, New York, NY, USA, 137–139.  https://doi.org/10.1145/3213586.3226202
  * Creswell (2003) John W. Creswell. 2003\.  _Research design: qualitative, quantitative, and mixed methods approaches_ (2nd ed.).  Sage Publications, Thousand Oaks, CA, USA.
  * Ekstrand and Willemsen (2016) Michael D. Ekstrand and Martijn C. Willemsen. 2016.  Behaviorism is not enough: better recommendations through listening to users. In _Proceedings of the 10th ACM Conference on Recommender Systems_ _(RecSys ’16)_. ACM, New York, NY, USA, 221–224.  https://doi.org/10.1145/2959100.2959179
  * Gunawardana and Shani (2015) Asela Gunawardana and Guy Shani. 2015.  Evaluating Recommender Systems.  In _Recommender Systems Handbook_ (2nd ed.), Francesco Ricci, Lior Rokach, and Bracha Shapira (Eds.). Springer, Boston, MA, USA, 265–308.  https://doi.org/10.1007/978-1-4899-7637-6_8
  * Herlocker et al. (2004) Jonathan L. Herlocker, Joseph A. Konstan, Loren G. Terveen, and John T. Riedl. 2004.  Evaluating Collaborative Filtering Recommender Systems.  _ACM Transaction on Information Systems_ 22, 1 (Jan. 2004), 5–53.  https://doi.org/10.1145/963770.963772
  * Jannach and Adomavicius (2017) Dietmar Jannach and Gediminas Adomavicius. 2017.  Price and profit awareness in recommender systems. In _1st International Workshop on Value-Aware and Multistakeholder Recommendation at RecSys 2017_ _(VAMS ’17)_.  arXiv:1707.08029
  * Jannach et al. (2016) Dietmar Jannach, Paul Resnick, Alexander Tuzhilin, and Markus Zanker. 2016.  Recommender Systems — Beyond Matrix Completion.  _Commun. ACM_ 59, 11 (2016), 94–102.  https://doi.org/10.1145/2891406
  * Kamehkhosh and Jannach (2017) Iman Kamehkhosh and Dietmar Jannach. 2017.  User Perception of Next-Track Music Recommendations. In _Proceedings of the 25th Conference on User Modeling, Adaptation and Personalization_ _(UMAP ’17)_. ACM, New York, NY, USA, 113–121.  https://doi.org/10.1145/3079628.3079668
  * Kaminskas and Bridge (2016) Marius Kaminskas and Derek Bridge. 2016.  Diversity, Serendipity, Novelty, and Coverage: A Survey and Empirical Analysis of Beyond-Accuracy Objectives in Recommender Systems.  _ACM Transactions on Interactive Intelligent Systems_ 7, 1, Article 2 (2016), 42 pages.  https://doi.org/10.1145/2926720
  * Knijnenburg and Willemsen (2015) Bart P. Knijnenburg and Martijn C. Willemsen. 2015.  Evaluating Recommender Systems with User Experiments.  In _Recommender Systems Handbook_ (2nd ed.), Francesco Ricci, Lior Rokach, and Bracha Shapira (Eds.). Springer, Boston, MA, USA, 309–352.  https://doi.org/10.1007/978-1-4899-7637-6_9
  * Kohavi et al. (2013) Ron Kohavi, Alex Deng, Brian Frasca, Toby Walker, Ya Xu, and Nils Pohlmann. 2013\.  Online Controlled Experiments at Large Scale. In _Proceedings of the 19th ACM SIGKDD International Conference on Knowledge Discovery and Data Mining_ _(KDD ’13)_. ACM, New York, NY, USA, 1168–1176.  https://doi.org/10.1145/2487575.2488217
  * Konstan and Riedl (2012) Joseph A. Konstan and John Riedl. 2012.  Recommender systems: from algorithms to user experience.  _User Modeling and User-Adapted Interaction_ 22, 1 (2012), 101–123.  https://doi.org/10.1007/s11257-011-9112-x
  * Laplante (2014) Audrey Laplante. 2014\.  Improving music recommender systems: what can we learn from research on music tags?. In _15th International Society for Music Information Retrieval Conference_ _(ISMIR ’14)_. International Society for Music Information Retrieval, 451–456.
  * McNee et al. (2006) Sean M. McNee, John Riedl, and Joseph A. Konstan. 2006\.  Being Accurate is Not Enough: How Accuracy Metrics Have Hurt Recommender Systems. In _CHI ’06 Extended Abstracts on Human Factors in Computing Systems_ _(CHI EA ’06)_. ACM, New York, NY, USA, 1097–1101.  https://doi.org/10.1145/1125451.1125659
  * Pu et al. (2011) Pearl Pu, Li Chen, and Rong Hu. 2011.  A User-centric Evaluation Framework for Recommender Systems. In _Proceedings of the 5th ACM Conference on Recommender Systems_ _(RecSys ’11)_. ACM, New York, NY, USA, 157–164.  https://doi.org/10.1145/2043932.2043962
  * Said et al. (2012) Alan Said, Domonkos Tikk, Klara Stumpf, Yue Shi, Martha Larson, and Paolo Cremonesi. 2012\.  Recommender Systems Evaluation: A 3D Benchmark. In _Proceedings of the Workshop on Recommendation Utility Evaluation: Beyond RMSE_ _(RUE ’12)_ , Vol. 910. CEUR Workshop Proceedings, 21–23.  http://ceur-ws.org/Vol-910/
  * Teddlie and Tashakkori (2009) Charles Teddlie and Abbas Tashakkori. 2009.  _Foundations of mixed methods research: Integrating quantitative and qualitative approaches in the social and behavioral sciences_.  Sage Publications, Thousand Oaks, CA, USA.
  * Venkatesh et al. (2013) Viswanath Venkatesh, Susan A. Brown, and Hillol Bala. 2013\.  Bridging the qualitative-quantitative divide: Guidelines for conducting mixed methods research in information systems.  _MIS Quarterly_ 37, 1 (2013), 21–54.
  * Zanker et al. (2019) Markus Zanker, Laurens Rook, and Dietmar Jannach. 2019\.  Measuring the impact of online personalisation: Past, present and future.  _International Journal of Human-Computer Studies_ (2019).  https://doi.org/10.1016/j.ijhcs.2019.06.006