# Reanalysis of the Systematic Uncertainties in Cosmic-Ray Antiproton Flux

Xing-Jian Lv<EMAIL_ADDRESS>Key Laboratory of Particle Astrophysics,
Institute of High Energy Physics, Chinese Academy of Sciences, Beijing 100049,
China  University of Chinese Academy of Sciences, Beijing 100049, China  Xiao-
Jun Bi<EMAIL_ADDRESS>Key Laboratory of Particle Astrophysics, Institute of
High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China
University of Chinese Academy of Sciences, Beijing 100049, China  Kun Fang
<EMAIL_ADDRESS>Key Laboratory of Particle Astrophysics, Institute of High
Energy Physics, Chinese Academy of Sciences, Beijing 100049, China  Peng-Fei
Yin<EMAIL_ADDRESS>Key Laboratory of Particle Astrophysics, Institute of
High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China  Meng-
Jie Zhao<EMAIL_ADDRESS>Key Laboratory of Particle Astrophysics, Institute
of High Energy Physics, Chinese Academy of Sciences, Beijing 100049, China
University of Chinese Academy of Sciences, Beijing 100049, China

###### Abstract

Recent studies on cosmic rays (CRs) have reported the possibility of an excess
in the antiproton flux around $10-20$ GeV. However, the associated systematic
uncertainties have impeded the interpretation of these findings. In this
study, we conduct a global Bayesian analysis to constrain the propagation
parameters and evaluate the CR antiproton spectrum, while comprehensively
accounting for uncertainties associated with interstellar CR propagation,
production cross sections for antiprotons and other secondaries, and the
charge and energy dependent effects of solar modulation. We establish that the
most recent AMS-02 $\bar{p}$ spectrum is in agreement with a pure secondary
origin. Based on this, we establish upper limits on dark matter (DM)
annihilation. We also determine that the AMS-02 data favors the empirical
hadronic interaction models over phenomenological ones. Finally, we find that
the latest AMS-02 antiproton data from 2011 to 2018 disfavors the antiproton
excess at $\mathcal{O}$(10) GeV and the corresponding DM interpretation that
can simultaneously account for the Galactic Center excess in the gamma-ray
observation.

##  I INTRODUCTION

Antimatter investigations in cosmic rays (CRs) have critical implications for
the fields of astrophysics and dark matter (DM) indirect detection. Over the
past few years, significant advancements have been made in the measurement of
antimatter particles in CRs. Of particular note is the Alpha Magnetic
Spectrometer (AMS-02) [1], which was launched in 2011 and has yielded
extremely precise CR measurements. The characteristics of CR antimatter
particles have been extensively analyzed quantitatively in the literature [2],
drawing on the findings of previous experiments such as PAMELA [3] and Fermi-
LAT [4], in conjunction with those from AMS-02.

Several recent studies have reported an overabundance of antiprotons in the
flux at energies between 10 and 20 GeV, which cannot be explained by their
production as secondary particles alone [5, 6, 7, 8, 9, 10]. These studies
proposed that the excess might be attributed to the annihilation or decay of
some weakly interacting massive particles (WIMPs), which are candidates of DM.
The proposed WIMPs have a mass in the range of $50\sim 100$ GeV and a
thermally-averaged annihilation cross section $\langle\sigma v\rangle$ near
the so-called nature value of $\sim 3\times 10^{-26}$ cm3s-1, and
predominantly annihilate to some hadronic final states. Additionally, several
authors suggested that these WIMPs can explain the Galactic Center excess
(GCE) observed in gamma-ray emissions by Fermi-LAT [11, 12, 13], and are in
line with other direct or indirect DM search results [14, 15]. On the other
hand, many recent studies also argued that the antiproton flux is consistent
with a pure secondary origin [16, 17, 18, 19, 20, 21, 22].

At the crux of the present controversy lies the intricacy of accurately
predicting the secondary antiproton flux. These antiprotons are generated via
the interaction between high energy CRs and the interstellar gas, with their
flux being determined by the process of CR propagation and its hadronic
interaction with the interstellar gas [23]. However, both of these factors
have considerable uncertainties. In this study, we present a systematic
analysis of the secondary antiproton flux, wherein we deliberately consider
the aforementioned sources of uncertainty. Then we can derive constraints on
the cross section of DM annihilation based on the predicted $\bar{p}$ flux and
the AMS-02 data.

The propagation of CRs within the Galaxy is a complex process that entails
diverse phenomena including diffusion, energy loss, convection, and
reacceleration. Stable primary-to-secondary ratios, like B/C and B/O, are the
main probes for studying CR propagation. Recent analysis utilizing high-
precision data provided by the AMS-02 collaboration has indicated that a
modified diffusion-reacceleration model (hereafter referred to as DR2 model)
is preferred [24, 25]. However, the present measurements of CR data are
insufficient to pinpoint the various propagation parameters as the
uncertainties involved in secondary production cross sections can bias these
parameters [26, 27, 28, 29]. The solar modulation effect on CRs induces
further uncertainties in propagation as it significantly affects low-energy CR
spectra, while its precise quantification remains challenging [26, 30].

We employ the numerical tool GALPROP v56111Current version available at
https://galprop.stanford.edu/ [31, 32] to solve the CR propagation equation.
The propagation parameters are constrained by the Bayesian analysis, where the
latest data on B/C and B/O are used. A Markov Chain Monte Carlo (MCMC)
algorithm [33] is applied to explore multidimensional parameter space with
high efficiency. The uncertainties in cross sections for Boron production are
fully taken into account to eliminate potential bias in the resulting
propagation parameters. A charge-sign dependent solar modulation model [34] is
adopted to describe the different modulation processes for $p$ and $\bar{p}$.

Another important source of uncertainties in predicting secondary antiproton
flux is the $\bar{p}$ production cross section. Fortunately, its uncertainty
has been reduced owing to the advancements in collider experiments [35, 36,
37]. We adopt the cross sections derived in Ref. [38] that incorporate the
most recent data. The uncertainties in $\bar{p}$ production are fully
propagated into the final flux. Additionally, we compare the performance of
the latest model with some older models that do not incorporate the latest
data, aiming to investigate the potential differences between them.

The paper is organized as follows. In Sec. II we detail the propagation model
and methodology employed to derive the propagation and source parameters which
are the basis of the calculation of the background antiproton flux. Our setup
of DM annihilation is also introduced in this section. In Sec. III we show the
resulting secondary antiproton spectrum under various settings and its
implications on DM annihilation. Finally, we summarize our findings and
provide insightful discussions in Sec. IV.

##  II METHOD

###  II.1 CR propagation model

Following acceleration within sources, galactic CRs undergo injection and
diffusion in the interstellar medium, subsequently experiencing various
propagation effects before reaching the Earth. It is commonly assumed that the
propagation of CRs takes place within a cylindrical halo with a vertical
extent of half height $L$, beyond which CRs are able to escape freely. This
process can be described mathematically with a propagation equation [39]

$\displaystyle\frac{\partial\psi}{\partial t}=$ $\displaystyle
Q(\mathbf{x},p)+\nabla\cdot\left(D_{xx}\nabla\psi-\mathbf{V}_{c}\psi\right)+\frac{\partial}{\partial
p}[p^{2}D_{pp}\frac{\partial}{\partial p}(\frac{\psi}{p^{2}})]$ (1)
$\displaystyle-\frac{\partial}{\partial
p}[\dot{p}\psi-\frac{p}{3}(\nabla\cdot\mathbf{V}_{c})\psi]-\frac{\psi}{\tau_{f}}-\frac{\psi}{\tau_{r}}\;,$

where $Q(x,p)$ denotes the CR source term, $\psi=\psi(x,p,t)$ is the CR
density per momentum interval, $\dot{p}\equiv dp/dt$ represents the momentum
loss rate, and the time scales $\tau_{f}$ and $\tau_{r}$ characterize
fragmentation processes and radioactive decays, respectively.

Supernova remnants are widely recognized as the main sources of Galactic CRs,
wherein charged particles experience acceleration through shock waves.
Consequently, the spatial distribution of the CR sources follows that of
supernova remnants [40]. According to the shock acceleration theory [41], the
injection spectrum of primary CRs follows a power-law distribution that is
dependent on the rigidity of CRs, $q\propto R^{-\nu}$. To account for the
observed structures in CR spectra, a three-piece broken power law is adopted
to describe the rigidity dependence of the injection spectrum [42].

The spatial diffusion coefficient $D_{xx}$ can be parametrized as [43]

$D_{xx}=D_{0}\beta^{\eta}\left(R/R_{0}\right)^{\delta}$ (2)

where $R\equiv pc/Ze$ is the rigidity of the CR particle, $\beta$ is the CR
particle velocity in units of the speed of light $c$, and $R_{0}$ and $D_{0}$
represent the reference rigidity and normalization parameter, respectively.
Although the theoretical prediction for the slope of the diffusion coefficient
$\delta$ is expected to be 1/3 for a Kolmogorov spectrum of interstellar
turbulence or 1/2 for a Kraichnan cascade, it is generally treated as a free
parameter to fit the data. Additionally, the slope of the velocity $\eta$
reflects the possibility that the turbulence dissipation could alter the
diffusion coefficient at low velocities [24].

The process of reacceleration of CRs caused by collisions with interstellar
weak hydrodynamic waves can be characterized by the diffusion in momentum
space, which is determined by the momentum diffusion coefficient $D_{pp}$.
This coefficient is interrelated with the spatial diffusion coefficient
$D_{xx}$ through [44]

$D_{pp}D_{xx}=\frac{4p^{2}v_{A}^{2}}{3\delta\left(4-\delta^{2}\right)(4-\delta)\omega}\;,$
(3)

where $v_{A}$ and $\omega$ denote the Alfvén velocity and the ratio of
magnetohydrodynamic wave energy density to magnetic field energy density,
respectively. Since $D_{pp}\propto v^{2}_{A}/\omega$, we can set $\omega=1$
without loss of generality. Furthermore, previous research has indicated that
convection is not favored by the observed CR spectrum [25]. Thus, for the
purposes of this study, $\mathbf{V}_{c}$ is set to $0$ throughout the
remainder of the article.

The CR spectra are greatly influenced by solar modulation below $\sim 20$ GeV.
The conventional approach to address this phenomenon involves utilizing the
force-field approximation, with the solar modulation potential $\phi$ employed
to quantify the associated strength [45]. However, the use of this approach in
situations involving particles of opposite charges is deemed inadequate due to
the presence of drift effects [30]. In Ref. [34], a modified technique was
proposed, incorporating a simple amendment to the solar modulation potential’s
structure to account for this challenge

$\phi^{\pm}(t,R)=\phi_{0}(t)+\phi_{1}^{\pm}(t)\mathcal{F}\left(R/R_{0}\right)\;,$
(4)

where
$\mathcal{F}=(1+\left(R/R_{0}\right)^{2})/(\beta\left(R/R_{0}\right)^{3})$ and
$R_{0}=0.5$ GV is the reference rigidity. The second term on the right-hand
side of the equation accounts for the augmented energy dissipation encountered
by particles whose charge polarity is mismatched with the current sheet. For
the positive (negative) polarity phase, the values of
$\phi^{+}_{1}$($\phi^{-}_{1}$) are assigned as 0.

In line with previous studies [18, 22], we fix the value of $\phi^{+}_{1}$ to
0, while treating the parameters $\phi_{0}$ and $\phi^{-}_{1}$ as time-
averaged values during the experimental process, which are obtained by fitting
the CR data. We further allow for the possibility that $\phi^{+}_{0}$ and
$\phi^{-}_{0}$ may have distinct values, as particles with opposite charges
traverse different regions of the solar system [46]. The results show that
they only differ by $\mathcal{O}(0.1\mathrm{GV})$.

###  II.2 Analysis setup

We employ a two-step approach, as previously proposed in Ref. [5, 7], to
compute the antiproton spectrum. The first step entails utilizing an MCMC
technique based on Bayesian inference to derive the propagation parameters
$D_{0}$, $\delta$, $v_{A}$, and $\eta$, and their corresponding covariance
matrix from the secondary to primary ratios B/C and B/O. In consideration of
the strong dependence of the DM induced signature flux on the propagation halo
height $L$ and the conflict between the Be/B and
$\ce{{}^{10}Be}/\ce{{}^{9}Be}$ data on $L$ [47, 48, 49], $L$ is not fitted in
this study. Instead, two benchmark values $L=3.44$ and 7.17 Kpc are selected,
which are inferred from the combined B/C+$\ce{{}^{10}Be}/\ce{{}^{9}Be}$ and
combined B/C+$\ce{{}^{10}Be}/\ce{{}^{9}Be}$+Be/B fits, respectively [50]. The
B/C and B/O ratios are rather insensitive to the source spectra of their
progenitors and the choice of the Fisk potential $\phi$, thus enabling us to
maintain a fixed value for the primary injection spectrum and $\phi$ for the
sake of simplicity. A preliminary fit yields the values of these parameters,
which are tabulated in Table 1.

To incorporate the uncertainties in the Boron production cross section, a data
driven method from earlier research is used [51, 26, 50], whereby
$\sigma_{\mathrm{prod.}}$ of 7 different channels are renormalized to match
the available data points through overall rescaling. These central values and
dispersions are included as priors during the fit and can be found in table
III of Ref. [50]. Consequently, the uncertainties in the production cross
section are naturally propagated into the transportation parameters.

In the second step, we determine the injection spectrum of protons and Helium
by performing a combined MCMC fit of the proton, Helium, and antiproton
spectra. Additionally, we conduct an ancillary fitting procedure that excluded
antiprotons, to ascertain whether the results obtained are consistent.
Furthermore, we allow the propagation parameters to vary to account for the
possibility of non-universality of CR transportation [52, 50]. Multivariate
Gaussian priors on propagation parameters are utilized to incorporate the
central values and covariance matrix of the propagation parameters from the
first step.

The production cross section of the antiproton is an essential factor in
calculating their spectrum, and we utilize the parametrization proposed in
Ref. [53] and improved in Ref. [38] (hereafter referred to as Winkler +
Korsmeier), which account for the effects of isospin violation, hyperon
contribution, and scaling violation. To incorporate the uncertainties in the
antiproton production cross section, we introduce an overall rescaling factor
$XS_{\bar{p}}$ with a Gaussian prior of $\mathcal{N}(1.0,0.08)$ to reproduce
the uncertainty band shown in Fig. 5 of Ref. [38].

Table 1: Injection spectrum of primary particles (except protons and Helium) and Fisk potential of heavy elements adopted in this study $\nu_{1}$ | $\nu_{2}$ | $\nu_{3}$ |  $R_{\mathrm{br1}}$(GV) |  $R_{\mathrm{br2}}$(GV) |  $\phi$(GV)
---|---|---|---|---|---
1.23 | 2.38 | 2.24 | 2.1 | 222.8 | 0.788

Regarding data selection, we employ the most up-to-date CR data available to
us. Specifically, the B/C and B/O data are sourced from AMS-02 observations
spanning the 2011-2016 period [54]. We further utilize proton flux data
collected by Voyager (2012-2015) [55], AMS-02 (2011-2018) [54], and DAMPE
(2016-2018) [56] to derive the proton injection spectral parameters. For
Helium, the AMS-02 (2011-2018) [54] and Voyager (2012-2015) [55] data are
used. Notably, AMS-02 proton and Helium data falling below 6 GV are excluded
from the analysis since they are below the antiproton production threshold
[23]. Finally, the $\bar{p}$ data are taken from AMS-02 observations during
the 2011-2018 period [54].

In this study, the MCMC procedure is implemented using the publicly available
codebase, cobaya222https://cobaya.readthedocs.io/en/latest/ [33, 57, 58]. The
Metropolis-Hastings algorithm is employed to generate samples from the
posterior distribution, which is a key step in Bayesian statistical inference.
In addition to the MCMC sampler, we utilize the Python package
GetDist333https://getdist.readthedocs.io/en/latest/ [59] which offers a suite
of analytical tools that facilitate the interpretation of the posterior
distribution, including methods for estimating credible intervals and
producing confidence contour plots.

In addition to the default setup, we perform several alternative analyses to
ensure the robustness of our results. One such approach involved examining the
use of different models for antiproton production cross sections. These models
can be categorized as either phenomenological or empirical. The former are
based on Monte Carlo generators and are often used to model hard-scattering
processes in colliders, which are not necessarily well-suited to CR antiproton
production dominated by soft production. Conversely, the empirical models are
derived from the fits to measurement results with certain analytical
parametrization. The model utilized in our default analysis belongs to the
second category and was fitted to the latest NA61 experiment [36] in the $pp$
channel and the first ever determination of $p$He data by the LHCb experiment
[37]. We also examine the use of the default hadronic model embedded in
GALPROP, which uses the parametrization proposed by Tan & Ng [60]. In
addition, we explore two phenomenological models that are best suited to CR
studies according to Ref. [20], QSJET II-04m [61] and LHC 1.99 [62], both of
which are tuned to low-energy data.

Moreover, we have conducted a quantitative assessment of the influence that
the latest antiproton data from AMS-02 2011-2018 [54] exerts when compared to
the previously employed data of AMS-02 2011-2015 [63], which remains the most
widely utilized data, even in contemporary publications. Our analysis reveals
that the utilization of the new data manifests a nuanced yet significant
impact on the ultimate outcomes, as demonstrated in Sec. III.4.

(a)  (b)

Figure 1: The 1-D and 2-D posterior distributions of the transport parameters
for different elements. The left and right panels show the results for the
halo height $L=7.17$ and 3.44 kpc, respectively. The fitting outcomes of step
one employing B/C and B/O are represented by the green line. The blue and
orange lines represent the results of step two utilizing the proton, Helium,
and antiproton data, and the outcomes obtained from the proton and Helium data
alone, respectively.

###  II.3 Computation of a WIMP contribution

Finally, we conduct an MCMC analysis to investigate the potential annihilation
of DM particles, which are WIMPs in this case, into hadronic states that
generate $\bar{p}$. The purpose of this analysis is to estimate the mass and
$\langle\sigma v\rangle$ of DM, that can account for the possible observed
inconsistencies when assuming a purely secondary origin of antiprotons.
Specifically, we focus on the $\chi\chi\rightarrow b\bar{b}$ annihilation
channel, as other quark final state channels yield similar contributions to
the production of antiprotons [12, 19]. The corresponding source term is given
by

$q_{\bar{p}}^{\mathrm{DM}}=\frac{\langle\sigma
v\rangle}{2m_{\chi}^{2}}\frac{\mathrm{d}N}{\mathrm{~{}d}E}\rho(\bm{x})\;,$ (5)

where the factor 1/2 corresponds to the DM particle being scalar or Majorana
fermion, $m_{\chi}$ is the mass of the DM particle, $\langle\sigma v\rangle$
is the thermally averaged DM annihilation cross section, and
$\mathrm{d}N/\mathrm{d}E$ is the antiproton production spectrum per
annihilation, which are given by PPPC 4 DM ID [64]. The DM density profile
$\rho(\bm{x})$ is taken to be the Navarro-Frenk-White distribution [65]

$\rho(r)=\frac{\rho_{s}}{\left(r/r_{s}\right)\left(1+r/r_{s}\right)^{2}}\;,$
(6)

with $r_{s}=20$ kpc and $\rho_{s}=0.35$ GeV. This choice of parameters
corresponds to a local DM density of 0.4 GeV $\text{cm}^{-3}$, which is
consistent with the recent constraints from the Galactic rotation curve [66,
67]. Alternative density profiles, such as the Einasto [68, 69] or Burkert
[70, 71] profiles, are not taken into account in our study, due to their
negligible impact on the final antiproton flux, as evidenced in Ref. [6, 72].

##  III RESULTS & DISCUSSION

###  III.1 Results of fiducial analysis

Table 2 displays the outcome of the fitting process in step one. Our
calculation for the B/C and B/O ratios fit well with the AMS-02 measurement,
as evident by the small $\chi^{2}$ statistic. It is noteworthy that the
optimal scaling factors for the various cross sections associated with
secondary production are confined within the $2\sigma$ threshold established
by the cross section data, with all factors exhibiting a variation of less
than 5%. As a result, the impact of secondary production cross sections on the
propagation parameters is negligible and does not significantly alter the
central values of propagation parameters. Their primary effect on the
propagation parameters is the slight broadening of their error bars.

Table 2:  Posterior mean and 95% credible uncertainties of the model parameters and $\chi^{2}$ value Parameter | L=7.17 kpc | L=3.44kpc
---|---|---
$D_{xx}(10^{28}\text{cm}^{2}\text{s}^{-1})$ | $6.03^{+0.50}_{-0.45}$ | $3.30^{+0.27}_{-0.27}$
$\delta$ | $0.440^{+0.021}_{-0.022}$ | $0.436^{+0.022}_{-0.022}$
$v_{A}$(km/s) | $21^{+3}_{-3}$ | $23^{+4}_{-3}$
$\eta$ | $-0.48^{+0.25}_{-0.25}$ | $-0.53^{+0.20}_{-0.20}$
$XS_{\ce{{}^{12}C}\to\ce{{}^{11}B}}$ | $0.961^{+0.033}_{-0.034}$ | $0.959^{+0.037}_{-0.036}$
$XS_{\ce{{}^{12}C}\to\ce{{}^{11}C}\to\ce{{}^{11}B}}$ | $1.013^{+0.015}_{-0.018}$ | $1.013^{+0.015}_{-0.015}$
$XS_{\ce{{}^{12}C}\to\ce{{}^{10}Be}}$ | $1.023^{+0.063}_{-0.065}$ | $1.019^{+0.061}_{-0.058}$
$XS_{\ce{{}^{16}O}\to\ce{{}^{10}Be}(\ce{{}^{10}C})}$ | $1.041^{+0.063}_{-0.063}$ | $1.044^{+0.059}_{-0.058}$
$XS_{\ce{{}^{16}O}\to\ce{{}^{11}Be}(\ce{{}^{11}C})}$ | $0.973^{+0.053}_{-0.055}$ | $0.979^{+0.050}_{-0.050}$
$XS_{\ce{{}^{14}N}\to\ce{{}^{11}Be}(\ce{{}^{11}C})}$ | $1.026^{+0.057}_{-0.053}$ | $1.028^{+0.056}_{-0.054}$
$XS_{\ce{{}^{15}N}\to\ce{{}^{11}Be}}$ | $1.009^{+0.090}_{-0.086}$ | $1.008^{+0.088}_{-0.093}$
$\chi^{2}/\text{d.o.f.}$ | $96.56/123$ | $98.3/123$

Fig. 1 is the triangle plot of the fitting results for step two and shows the
1D marginalized posterior probability density functions of the parameters and
2D contour plots of 68% and 95% credible regions for all the combinations. The
injection parameters are omitted for simplicity. It is evident that most of
the propagation parameters are consistent with those obtained during step one,
except for $\eta$. Notably, the outcomes attained in the presence or absence
of $\bar{p}$ are almost indistinguishable, differing only in a negligible
discrepancy in the width of their respective posterior distributions. This
similarity suggests that the impact of $\bar{p}$ on the overall results is
minimal, indicating the internal consistency of the study.

The first two columns of Table 3 present a summary of the results, listing the
mean values and posterior 95% ranges. The propagation parameters with and
without antiproton exhibit full consistency within uncertainties for both $L$
values. The halo diffusion slope indexes $\delta$ fall within the range of the
Kolmogorov type with $\delta=1/3$ and the Kraichnan type with $\delta=1/2$.
The fitted scale factors $c_{\bar{p}}$ are within the experiment’s uncertainty
of $\sigma_{\mathrm{norm.}}\simeq 0.08$ as given in Ref. [38]. Regarding solar
modulation, $\phi^{-}_{0}$ is smaller than $\phi^{-}_{0}$ by
$\mathcal{O}(0.1)$ GV. Confirmation of this result must wait for the monthly
antiproton flux’s release by the AMS-02 collaboration.

It is worth mentioning that the value of
$\chi^{2}_{\mathrm{min}}/\mathrm{DoF}$ is significantly lower than 1. This may
be due to that we employ a quadratic combination of the AMS-02 statistical and
systematic errors, thereby neglecting the correlations present within the
AMS-02 systematic errors. Nevertheless, it has been previously suggested that
the inclusion of these correlations, as proposed in Refs. [16, 17], would
result in even flatter residues, which does not alter the conclusions drawn in
this study. However, future investigations that account for error correlations
are contingent upon the AMS-02 collaboration publishing their covariance
matrix [73].

The $\bar{p}$ spectrum, calculated using the propagation parameters determined
in step two, is depicted in Fig. 2. The spectrum is in excellent agreement
with the corresponding AMS-02 data in most energy ranges. However, a slight
discrepancy is observed at the low energy end of the spectrum, which may be
attributed to the over-simplified solar modulation model. Furthermore, a few
data points at the high energy end exceed our calculation. Nevertheless, the
large error bars associated with these data points prevent drawing any
meaningful conclusions.

The results obtained above lead to two significant conclusions. Firstly, the
exceptional precision of the proton and Helium measurements has the potential
to alleviate the degeneracy between the source and propagation parameters to
some extent, which is a crucial finding as it suggests the possibility of
determining propagation parameters solely based on primary particles.
Therefore, this indicates that the resulting propagation parameters remain
unaffected by the uncertainties arising from hadronic interactions. Secondly,
the antiproton spectrum observed from AMS-02 aligns entirely with a secondary
astrophysical origin, suggesting that no further primary antiproton sources
are required to explain the AMS-02 data.

Figure 2: Best-fit model predictions of the antiproton fluxes, compared with
the AMS-02 data [54]. The top and bottom panels show the results for the thick
halo with $L=7.17$ kpc and thin halo with $L=3.44$ kpc, respectively.

###  III.2 Impact of hadronic interaction models

In this subsection, various hadronic models are tested to determine their
efficacy in antiproton production cross section. Three other hadronic models,
namely LHC 1.99 [62], QSJET II-04 m [61], and Tan & Ng [60], are considered
for $L=7.17$ kpc here. We firstly plot the predictions of different hadronic
interaction models in Fig. 3, utilizing the parameters established from the
fit to the protons and Helium data. It is clear from the figure that the two
phenomenological models, LHC 1.99 and QSJET II-04m, are unable to accurately
replicate the $\bar{p}$ data. In contrast, the two empirical models, Winkler +
Korsmeier and Tan & Ng, produce comparable results that align with the data
within experimental error in most energy ranges. It should be noted that the
GALPROP default Tan & Ng model appears to underestimate the $\bar{p}$ flux at
the high energy end, with notable discrepancies compared to the last several
data points.

Subsequently, various MCMC fits are performed utilizing these models. The
results are presented in Figure 4. Our examination of Figure 3 is reinforced
by the analysis, which demonstrates that the phenomenological models are less
effective in comparison to the empirical ones. Notably, the propagation
parameters with and without antiproton are no longer entirely consistent, with
the greatest discrepancy occurring when the QSJET II-04m model is employed.
The results of the fits are presented in the final column of Table 3. It is
evident that the $\chi^{2}$ statistic significantly worsens when compared to
that of the empirical models. Furthermore, the two phenomenological models
display a scale factor of $c_{\bar{p}}$ approximately 20%, which is
considerably higher than the 5% scale factor observed in the empirical models.
Moreover, it is worth noting that the QSJET II-04m model requires a solar
modulation potential greater than 1GV, which could be unphysical.

Concerning the comparison between the two empirical models, the latest Winkler
+ Korsmeier model that incorporates the collider data is slightly favored over
the Tan & Ng model. This can be seen from the fact that the propagation
parameters derived with the Tan & Ng model demonstrate a slight discrepancy
from the fit involving only the proton and Helium data, necessitating a larger
$c_{\bar{p}}$ of 5% in contrast to 2%. Nevertheless, these divergences are all
contained within the 2$\sigma$ uncertainty limits, thus we cannot assert
definitively that Winkler + Korsmeier is the preferred model.

Table 3: The mean values and posterior 95% ranges of key parameters in different scenarios. For the results of the alternative hadronic models listed in the last three columns, $L$ is taken to be 7.1.7 kpc. | Standard Analysis |  | w/o $\bar{p}$ |  | Hadronic Models
---|---|---|---|---|---
|  $L=7.17$ kpc |  $L=3.44$ kpc |  |  $L=7.17$ kpc |  $L=3.44$ kpc |  | Tan & Ng | EPOS 1.99 | QGSJETII-04m
$c_{\bar{p}}$ | $1.019^{+0.036}_{-0.034}$ | $1.002^{+0.033}_{-0.032}$ |  | - | - |  | $1.050^{+0.034}_{-0.033}$ | $0.824^{+0.030}_{-0.028}$ | $1.209^{+0.040}_{-0.039}$
$\phi^{+}_{0}$(GV) | $0.672^{+0.057}_{-0.053}$ | $0.735^{+0.054}_{-0.058}$ |  | $0.667^{+0.061}_{-0.065}$ | $0.736^{+0.058}_{-0.053}$ |  | $0.677^{+0.047}_{-0.047}$ | $0.652^{+0.053}_{-0.053}$ | $0.708^{+0.064}_{-0.064}$
$\phi^{-}_{0}$(GV) | $0.54^{+0.10}_{-0.11}$ | $0.49^{+0.10}_{-0.11}$ |  | - | - |  | $0.627^{+0.097}_{-0.10}$ | $0.654^{+0.094}_{-0.10}$ | $1.20^{+0.10}_{-0.10}$
$\phi^{-}_{1}$(GV) | $0.27^{+0.25}_{-0.27}$ | $<0.348$ |  | - | - |  | $0.43^{+0.29}_{-0.27}$ | $<0.203$ | $0.98^{+0.32}_{-0.30}$
$\chi^{2}$/DoF | 73.6/183 | 65.1/183 |  | 49.2/128 | 43.1/128 |  | 70.0/183 | 93.4/183 | 99.3/183
$D_{xx}(10^{28}\text{cm}^{2}\text{s}^{-1})$ | $5.96^{+0.42}_{-0.38}$ | $3.35^{+0.21}_{-0.21}$ |  | $5.94^{+0.47}_{-0.44}$ | $3.37^{+0.22}_{-0.21}$ |  | $6.02^{+0.36}_{-0.36}$ | $5.91^{+0.38}_{-0.36}$ | $6.03^{+0.43}_{-0.48}$
$\delta$ | $0.446^{+0.016}_{-0.017}$ | $0.434^{+0.016}_{-0.016}$ |  | $0.447^{+0.018}_{-0.019}$ | $0.432^{+0.017}_{-0.017}$ |  | $0.443^{+0.015}_{-0.014}$ | $0.451^{+0.015}_{-0.016}$ | $0.438^{+0.019}_{-0.019}$
$v_{A}$(km/s) | $23^{+3}_{-2}$ | $25.2^{+2.5}_{-2.6}$ |  | $22^{+3}_{-3}$ | $25.3^{+2.7}_{-2.5}$ |  | $23.0^{+2.2}_{-2.3}$ | $22^{+3}_{-3}$ | $24^{+3}_{-3}$
$\eta$ | $-0.82^{+0.12}_{-0.12}$ | $-0.75^{+0.11}_{-0.11}$ |  | $-0.83^{+0.16}_{-0.17}$ | $-0.74^{+0.13}_{-0.12}$ |  | $-0.81^{+0.12}_{-0.12}$ | $-0.85^{+0.11}_{-0.13}$ | $-0.78^{+0.13}_{-0.14}$

Figure 3: The $\bar{p}$ flux expected from different hadronic models for
$L=7.17$ kpc, in comparison with the AMS-02 measurement [54]. Four lines
represent the results for the models of Winkler + Korsmeier (orange), LHC 1.99
(purple), QSJET II-04 m (red), and Tan & Ng (green), respectively. The
parameters used here give the best fit to the proton and Helium data. The
solar modulation potentials are set to be $\phi^{-}_{0}=0.56$GV,
$\phi^{-}_{1}=0.40$GV.

In conclusion, our findings show that the two empirical models are strongly
preferred over the two phenomenological models, which are already superior to
other phenomenological models in reproducing collider data in the low energy
region. This finding also highlights the potential for CR data to be utilized
as a means for selecting the optimal hadronic models. Furthermore, with a
deeper understanding of the propagation mechanism and access to more accurate
CR data, the possibility of determining nuclear production cross section via
CRs is becoming increasingly plausible.

Figure 4: Posterior distribution of propagation parameters using different
hadronic models for halo height $L=7.17$ kpc. The solid line represents the
fitting outcomes obtained exclusively from the proton and Helium data, while
the dashed lines depict the fitting results with antiproton data, which are
calculated using the models proposed by LHC 1.99 (blue), QSJET II-04 m
(yellow), and Tan & Ng (gray).

###  III.3 Testing possible antiproton production from DM

The AMS-02 antiproton data exhibits a high degree of agreement with our
calculations, thereby resulting in a robust constraint on DM annihilation. We
utilize the standard Bayes analysis protocol to ascertain the corresponding
95% upper limit for a specific value of $m_{\chi}$ through the following
equation [74]

$\frac{\int_{0}^{\langle\sigma
v\rangle_{95}}\mathcal{P}(x)\mathrm{d}x}{\int_{0}^{\infty}\mathcal{P}(x)\mathrm{d}x}=0.95\;.$
(7)

The results are shown in Fig. 5. The limits derived from the investigation are
demonstrated to be typically stronger in comparison to the outcomes of gamma-
ray observations of dwarf galaxies [75]. However, it should be noted that the
effect of varying halo sizes is not a mere shift of the exclusion limit as
naively expected. Instead, the line shape of the limit undergoes moderate
changes. Nonetheless, the qualitative conclusion that the DM annihilation
cross section increases as the halo thickness decreases remains valid.
Additionally, it is important to mention that the constraints established in
this study may be subject to a constant factor reduction (increase) if the
local density of DM is proven to be higher (lower) [76].

Table 4: DM parameters, $\chi^{2}$ values, number of free fit parameters, and significance for the best fits. The total $\chi^{2}$ values refer to the fits with (without) DM. Parameter | L=7.17 kpc | L=3.44kpc
---|---|---
$m_{\chi}$(GeV) | $132$ | $166$
$\langle\sigma v\rangle(10^{-26}\text{cm}^{3}\text{s}^{-1}$) | $0.98$ | $1.71$
$\chi^{2}_{\text{tot}}$(200 points) | $65.1(73.6)$ | $59.9(65.1)$
No. of fit param. | $19(17)$ | $19(17)$
$\Delta\chi^{2}$ | $8.5$ | $5.2$
local sig. | $2.2\sigma$ | $1.4\sigma$

Table 4 presents the properties of DM as inferred from fitting, with the local
significance of the signal determined using the likelihood ratio test. It
should be noted that the global significance of the DM signal is lower than
the provided local significance due to the impact of the look-elsewhere effect
(LEE) [77]. Because the local evidence for a DM signal is statistically
insignificant, we refrain from conducting a computationally demanding
evaluation of its global significance. Nonetheless, according to Ref. [17],
the global significance of the signal was estimated to be approximately
$0.8\sim 1\sigma$ lower than its local significance.

Figure 5: The 95% upper limits on the DM annihilation rate $\langle\sigma
v\rangle$ derived from the AMS-02 $\bar{p}$ data. The top panel is for
$L=7.17$ kpc, and the bottom panel is for $L=3.44$ kpc. For comparison, we
also show the upper limits from the Fermi-LAT observation of dwarf galaxies
[75] as the red lines. The dashed gray curve corresponds to the thermal relic
cross section from Ref. [78]. The green cross is the best-fit DM signal.

###  III.4 Impact of new AMS-02 data

The derived DM mass accounting for the potential antiproton excess, which
falls within the range of 130-160 GeV, is significantly greater than that
observed in prior investigations [7, 8, 79] and is incompatible with the
gamma-ray Galactic Center excess [11]. The study attributes this result to the
use of the most recent AMS-02 antiproton data acquired between 2011 and 2018
[54]. Fig. 6 illustrates the comparison of the old and new antiproton data
sets, demonstrating that while both data sets exhibit consistency within their
respective margins of error, there exist significant discrepancies in the
central values of data points between 10 and 20 GV. Specifically, the bump
structure evident in the old data, which has previously been a subject of
controversy in the context of the antiproton excess, is absent in the new
data.

Figure 6: The subtle difference between the AMS-02 data taken during 2011-2015
[63] and 2011-2018 [54] in the energy range of 10-20 GV. The solid green line
is the best-fit spectrum to the data taken during 2011-2018.

We further perform an MCMC fit using the old data. As shown in Fig. 7, the
results of the MCMC analysis align with the observations described above.
Specifically, the estimated mass of the hypothetical DM signal is 76.7 GeV,
which is in agreement with the findings of previous studies conducted using
the old data.

The results of the current investigation disfavor the feasibility of utilizing
a signal DM model to explain both the GCE and the antiproton excess. The
antiproton excess, whose existence is contested based on the study findings,
would necessitate a DM mass that exceeds the range compatible with the GCE.
This same limitation also pertains to efforts to address the W-boson mass
anomaly while simultaneously accounting for the antiproton excess, as
previously discussed in Ref [10]. This discovery is in qualitative accordance
with the findings presented in Ref. [22].

Figure 7: Same as Fig. 5 for $L=7.17$ kpc, but obtained with the old AMS-02
$\bar{p}$ data [63].

##  IV SUMMARY

The current research endeavor involves the evaluation of the CR antiproton
spectrum, which is then compared with the most recent AMS-02 results. Given
the absence of any discernible excess in the $\bar{p}$ spectrum, the study
subsequently establishes upper limits on the annihilation cross section of DM.

In this study, the DR2 propagation model is selected due to its preference for
the latest secondary-to-primary ratios from the AMS-02 experiment. Using the
latest AMS-02 B/C and B/O data and incorporating uncertainties in the Boron
production cross section, the constraints on the propagation parameters are
updated. The injection spectrum for proton and Helium is then determined via a
combined fit to the AMS-02 proton, Helium, and antiproton data. The fitting
process incorporates a charge-sign dependent solar modulation model and the
most current $\bar{p}$ production cross section model. Alternative $\bar{p}$
production models that do not incorporate the most recent collider data are
also considered, and the effect of using the latest AMS-02 $\bar{p}$ data is
evaluated.

The most important findings of this work are summarized as follows:

  * •

The $\bar{p}$ spectrum is consistent with a secondary origin, and therefore no
contribution from DM is necessary.

  * •

The AMS-02 proton and Helium data are of sufficient precision to break the
degeneracy between injection and propagation.

  * •

Accounting for the various systematic uncertainties, the derived upper bound
on the DM annihilation cross section, calculated using the $\bar{p}$ data,
could be more stringent than the corresponding constraint obtained from the
Fermi observation of the dwarf galaxies.

  * •

The utilization of the most recent AMS-02 $\bar{p}$ data has a notable impact.
Specifically, the plausibility of a self-consistent explanation for the GCE
and the antiproton excess is no longer tenable upon the utilization of the
recent $\bar{p}$ data.

  * •

The precision of the AMS-02 CR data is sufficient to enable the identification
of the most appropriate cross section models: the empirical hadronic
interaction models are superior to the phenomenological ones.

In the future, the modeling of the transport of CRs in the Galaxy can be
improved through the incorporation of more precise models, such as the
spatial-dependent diffusion [80]. Utilizing various numerical tools to model
solar modulation [81, 82] can also improve our understanding of its impact on
the $\bar{p}$ flux. Ultimately, these models will facilitate a more
comprehensive understanding of the mechanisms driving the $\bar{p}$ flux,
which will enable a more accurate understanding of the nature and properties
of DM.

###### Acknowledgements.

This work is supported by the National Natural Science Foundation of China
under the grants No. U1738209, No. U2031110, and No. 12105292.

## References

  * [1] AMS-02 — The Alpha Magnetic Spectrometer Experiment, https://ams02.space/.
  * Blum  _et al._ [2017] K. Blum, R. Sato, and E. Waxman, arXiv e-prints 10.48550/arXiv.1709.06507 (2017), arxiv:1709.06507 [astro-ph.HE] .
  * [3] PAMELA experiment web page, https://pamela-web.web.roma2.infn.it/.
  * [4] The Fermi Gamma-ray Space Telescope, https://fermi.gsfc.nasa.gov/.
  * Cui  _et al._ [2017] M.-Y. Cui, Q. Yuan, Y.-L. S. Tsai, and Y.-Z. Fan, Physical Review Letters 118, 191101 (2017), arxiv:1610.03840 [astro-ph.HE] .
  * Cuoco  _et al._ [2017a] A. Cuoco, M. Krämer, and M. Korsmeier, Physical Review Letters 118, 191102 (2017a), arxiv:1610.03071 [astro-ph.HE] .
  * Cui  _et al._ [2018] M.-Y. Cui, X. Pan, Q. Yuan, Y.-Z. Fan, and H.-S. Zong, Journal of Cosmology and Astroparticle Physics 06 (06), 024, arxiv:1803.02163 [astro-ph.HE] .
  * Cuoco  _et al._ [2019] A. Cuoco, J. Heisig, L. Klamt, M. Korsmeier, and M. Krämer, Physical Review D 99, 103014 (2019), arxiv:1903.01472 [astro-ph.HE] .
  * Cholis  _et al._ [2019] I. Cholis, T. Linden, and D. Hooper, Physical Review D 99, 103026 (2019), arxiv:1903.02549 [astro-ph.HE] .
  * Zhu  _et al._ [2022] C.-R. Zhu, M.-Y. Cui, Z.-Q. Xia, Z.-H. Yu, X. Huang, Q. Yuan, and Y.-Z. Fan, Phys. Rev. Lett. 129, 231101 (2022), arxiv:2204.03767 [astro-ph.HE] .
  * Di Mauro and Winkler [2021] M. Di Mauro and M. W. Winkler, Physical Review D 103, 123005 (2021), arxiv:2101.11027 [astro-ph.HE] .
  * Cuoco  _et al._ [2017b] A. Cuoco, J. Heisig, M. Korsmeier, and M. Krämer, Journal of Cosmology and Astroparticle Physics 10 (10), 053, arxiv:1704.08258 [astro-ph.HE] .
  * Di Mauro [2021] M. Di Mauro, Physical Review D 103, 063029 (2021), arxiv:2101.04694 [astro-ph.HE] .
  * Leane [2020] R. K. Leane, arXiv e-prints (2020), arxiv:2006.00513 [hep-ph] .
  * Hooper  _et al._ [2020] D. Hooper, R. K. Leane, Y.-D. Tsai, S. Wegsman, and S. J. Witte, Journal of High Energy Physics 07, 163 (2020), arxiv:1912.08821 [hep-ph] .
  * Boudaud  _et al._ [2020] M. Boudaud, Y. Génolini, L. Derome, J. Lavalle, D. Maurin, P. Salati, and P. D. Serpico, Physical Review Research 2, 023022 (2020), arxiv:1906.07119 [astro-ph.HE] .
  * Heisig  _et al._ [2020] J. Heisig, M. Korsmeier, and M. W. Winkler, Physical Review Research 2, 043017 (2020), arxiv:2005.04237 [astro-ph.HE] .
  * Reinert and Winkler [2018] A. Reinert and M. W. Winkler, JCAP 2018 (01), 055, arxiv:1712.00002 [astro-ph.HE] .
  * Calore  _et al._ [2022] F. Calore, M. Cirelli, L. Derome, Y. Genolini, D. Maurin, P. Salati, and P. D. Serpico, SciPost Phys. 12, 163 (2022).
  * Lin  _et al._ [2017] S.-J. Lin, X.-J. Bi, J. Feng, P.-F. Yin, and Z.-H. Yu, Physical Review D 96, 123010 (2017), arxiv:1612.04001 [astro-ph.HE] .
  * Lin  _et al._ [2019] S.-J. Lin, X.-J. Bi, and P.-F. Yin, Physical Review D 100, 103014 (2019), arxiv:1903.09545 [astro-ph.HE] .
  * Luque [2021] P. D. L. T. Luque, JCAP 2021 (11), 018, arxiv:2107.06863 [astro-ph.HE] .
  * Donato  _et al._ [2001] F. Donato, D. Maurin, P. Salati, A. Barrau, G. Boudoul, and R. Taillet, Astrophys. J. 563, 172 (2001), arxiv:astro-ph/0103150 .
  * Yuan  _et al._ [2017] Q. Yuan, S.-J. Lin, K. Fang, and X.-J. Bi, Phys. Rev. D 95, 083007 (2017), arXiv:1701.06149 [astro-ph.HE] .
  * Yuan [2018] Q. Yuan, Science China Physics, Mechanics & Astronomy 62, 49511 (2018), arxiv:1805.10649 [astro-ph.HE] .
  * Tomassetti [2017] N. Tomassetti, Physical Review D 96, 103005 (2017), arxiv:1707.06917 [astro-ph, physics:hep-ph] .
  * Génolini  _et al._ [2018] Y. Génolini, D. Maurin, I. V. Moskalenko, and M. Unger, Physical Review C 98, 034611 (2018), arxiv:1803.04686 [astro-ph.HE] .
  * Weinrich  _et al._ [2020a] N. Weinrich, Y. Génolini, M. Boudaud, L. Derome, and D. Maurin, Astronomy & Astrophysics 639, A131 (2020a), arxiv:2002.11406 [astro-ph.HE] .
  * Korsmeier and Cuoco [2021] M. Korsmeier and A. Cuoco, Physical Review D 103, 103016 (2021), arxiv:2103.09824 [astro-ph.HE] .
  * Potgieter [2013] M. S. Potgieter, Living Reviews in Solar Physics 10, 3 (2013), arxiv:1306.4421 [physics.space-ph] .
  * Strong and Moskalenko [1998] A. W. Strong and I. V. Moskalenko, The Astrophysical Journal 509, 212 (1998), arxiv:astro-ph/9807150 .
  * Strong  _et al._ [2000] A. W. Strong, I. V. Moskalenko, and O. Reimer, The Astrophysical Journal 537, 763 (2000), arxiv:astro-ph/9811296 .
  * Lewis and Bridle [2002] A. Lewis and S. Bridle, Phys. Rev. D66, 103511 (2002), arXiv:astro-ph/0205436 [astro-ph] .
  * Cholis  _et al._ [2016] I. Cholis, D. Hooper, and T. Linden, Physical Review D 93, 043016 (2016), arxiv:1511.01507 [astro-ph.SR] .
  * Baatar  _et al._ [2013] B. Baatar, G. Barr, J. Bartke, _et al._ (NA49), Eur. Phys. J. C 73, 2364 (2013), arxiv:1207.6520 [hep-ex] .
  * Collaboration [2017] N. Collaboration (NA61/SHINE), The European Physical Journal C 77, 671 (2017), arxiv:1705.02467 [nucl-ex] .
  * Aaij  _et al._ [2018] R. Aaij _et al._ (LHCb), Phys. Rev. Lett. 121, 222001 (2018), arxiv:1808.06127 [hep-ex] .
  * Korsmeier  _et al._ [2018] M. Korsmeier, F. Donato, and M. Di Mauro, Physical Review D 97, 103019 (2018), arxiv:1802.03030 [astro-ph.HE] .
  * Strong  _et al._ [2007] A. W. Strong, I. V. Moskalenko, and V. S. Ptuskin, Annual Review of Nuclear and Particle Science 57, 285 (2007), arxiv:astro-ph/0701517 .
  * Trotta  _et al._ [2011] R. Trotta, G. Johannesson, I. V. Moskalenko, T. A. Porter, R. R. de Austri, and A. W. Strong, The Astrophysical Journal 729, 106 (2011), arxiv:1011.0037 [astro-ph.HE] .
  * Achterberg  _et al._ [2001] A. Achterberg, Y. A. Gallant, J. G. Kirk, and A. W. Guthmann, Monthly Notices of the Royal Astronomical Society 328, 393 (2001), arxiv:astro-ph/0107530 .
  * Ahn  _et al._ [2010] H. S. Ahn, P. Allison, M. G. Bagliesi, J. J. Beatty, G. Bigongiari, J. T. Childers, N. B. Conklin, S. Coutu, M. A. DuVernois, O. Ganel, J. H. Han, J. A. Jeon, K. C. Kim, M. H. Lee, L. Lutz, P. Maestro, A. Malinin, P. S. Marrocchesi, S. Minnick, S. I. Mognet, J. Nam, S. Nam, S. L. Nutter, I. H. Park, N. H. Park, E. S. Seo, R. Sina, J. Wu, J. Yang, Y. S. Yoon, R. Zei, and S. Y. Zinn, The Astrophysical Journal 714, L89 (2010), arxiv:1004.1123 [astro-ph.HE] .
  * Maurin  _et al._ [2010] D. Maurin, A. Putze, and L. Derome, Astronomy and Astrophysics 516, A67 (2010), arxiv:1001.0553 [astro-ph.HE] .
  * Drury and Strong [2017] L. O. Drury and A. W. Strong, Astronomy & Astrophysics 597, A117 (2017), arxiv:1608.04227 [astro-ph.HE] .
  * Gleeson and Axford [1968] L. J. Gleeson and W. I. Axford, The Astrophysical Journal 154, 1011 (1968).
  * Strauss  _et al._ [2012] R. D. Strauss, M. S. Potgieter, I. Büsching, and A. Kopp, Astrophysics and Space Science 339, 223 (2012).
  * Weinrich  _et al._ [2020b] N. Weinrich, M. Boudaud, L. Derome, Y. Génolini, J. Lavalle, D. Maurin, P. Salati, P. Serpico, and G. Weymann-Despres, Astronomy & Astrophysics 639, A74 (2020b), arxiv:2004.00441 [astro-ph.HE] .
  * Luque  _et al._ [2021] P. d. l. T. Luque, M. N. Mazziotta, F. Loparco, F. Gargano, and D. Serini, Journal of Cosmology and Astroparticle Physics 03 (03), 099, arxiv:2101.01547 [astro-ph.HE] .
  * Evoli  _et al._ [2020] C. Evoli, G. Morlino, P. Blasi, and R. Aloisio, Physical Review D 101, 023013 (2020), arxiv:1910.04113 [astro-ph.HE] .
  * Zhao  _et al._ [2023] M.-J. Zhao, X.-J. Bi, and K. Fang, Phys. Rev. D 107, 063020 (2023), arxiv:2209.03799 [astro-ph.HE] .
  * Tomassetti [2015] N. Tomassetti, Physical Review C 92, 045808 (2015), arxiv:1509.05776 [astro-ph.HE] .
  * Jóhannesson  _et al._ [2016] G. Jóhannesson, R. R. de Austri, A. C. Vincent, I. V. Moskalenko, E. Orlando, T. A. Porter, A. W. Strong, R. Trotta, F. Feroz, P. Graff, and M. P. Hobson, The Astrophysical Journal 824, 16 (2016), arxiv:1602.02243 [astro-ph.HE] .
  * Winkler [2017] M. W. Winkler, Journal of Cosmology and Astroparticle Physics 02 (02), 048, arxiv:1701.04866 [hep-ph] .
  * Aguilar  _et al._ [2021] M. Aguilar _et al._ (AMS), Physics Reports The Alpha Magnetic Spectrometer (AMS) on the International Space Station: Part II - Results from the First Seven Years, 894, 1 (2021).
  * Cummings  _et al._ [2016] A. C. Cummings, E. C. Stone, B. C. Heikkila, N. Lal, W. R. Webber, G. Jóhannesson, I. V. Moskalenko, E. Orlando, and T. A. Porter, The Astrophysical Journal 831, 18 (2016).
  * An  _et al._ [2019] Q. An _et al._ (DAMPE), Science Advances 5, eaax3793 (2019), arxiv:1909.12860 [astro-ph.HE] .
  * Lewis [2013] A. Lewis, Phys. Rev. D87, 103529 (2013), arXiv:1304.4473 [astro-ph.CO] .
  * Torrado and Lewis [2021] J. Torrado and A. Lewis, Journal of Cosmology and Astroparticle Physics 2021 (05), 057, arxiv:2005.05290 [astro-ph.IM] .
  * Lewis [2019] A. Lewis, arXiv e-prints , arXiv:1910.13970 (2019), arXiv:1910.13970 [astro-ph.IM] .
  * Tan and Ng [1983] L. C. Tan and L. K. Ng, Journal of Physics G: Nuclear Physics 9, 1289 (1983).
  * Kachelriess  _et al._ [2015] M. Kachelriess, I. V. Moskalenko, and S. S. Ostapchenko, The Astrophysical Journal 803, 54 (2015), arxiv:1502.04158 [astro-ph.HE] .
  * Pierog and Werner [2009] T. Pierog and K. Werner, Nuclear Physics B - Proceedings Supplements 196, 102 (2009), arxiv:0905.1198 [hep-ph] .
  * Aguilar  _et al._ [2016] M. Aguilar _et al._ (AMS), Physical Review Letters 117, 091103 (2016).
  * Cirelli  _et al._ [2011] M. Cirelli, G. Corcella, A. Hektor, G. Hütsi, M. Kadastik, P. Panci, M. Raidal, F. Sala, and A. Strumia, Journal of Cosmology and Astroparticle Physics 03 (03), 051, arxiv:1012.4515 [hep-ph] .
  * Navarro  _et al._ [1997] J. F. Navarro, C. S. Frenk, and S. D. M. White, The Astrophysical Journal 490, 493 (1997), arxiv:astro-ph/9611107 .
  * Karukes  _et al._ [2019] E. V. Karukes, M. Benito, F. Iocco, R. Trotta, and A. Geringer-Sameth, Journal of Cosmology and Astroparticle Physics 09 (09), 046, arxiv:1901.02463 [astro-ph.GA] .
  * Benito  _et al._ [2019] M. Benito, A. Cuoco, and F. Iocco, Journal of Cosmology and Astroparticle Physics 03 (03), 033, arxiv:1901.02460 [astro-ph.GA] .
  * Merritt  _et al._ [2006] D. Merritt, A. W. Graham, B. Moore, J. Diemand, and B. Terzic, The Astronomical Journal 132, 2685 (2006), arxiv:astro-ph/0509417 .
  * Navarro  _et al._ [2010] J. F. Navarro, A. Ludlow, V. Springel, J. Wang, M. Vogelsberger, S. D. M. White, A. Jenkins, C. S. Frenk, and A. Helmi, Monthly Notices of the Royal Astronomical Society 402, 21 (2010), arxiv:0810.1522 [astro-ph] .
  * Burkert [1996] A. Burkert, The Astrophysical Journal 171, 175 (1996), arxiv:astro-ph/9504041 .
  * Salucci and Burkert [2000] P. Salucci and A. Burkert, The Astrophysical Journal 537, L9 (2000), arxiv:astro-ph/0004397 .
  * Cuoco  _et al._ [2018] A. Cuoco, J. Heisig, M. Korsmeier, and M. Krämer, JCAP 2018 (04), 004, arxiv:1711.05274 [hep-ph] .
  * Derome  _et al._ [2019] L. Derome, D. Maurin, P. Salati, M. Boudaud, Y. Génolini, and P. Kunzé, Astronomy & Astrophysics 627, A158 (2019), arxiv:1904.08210 [astro-ph.HE] .
  * Brooks  _et al._ [2011] S. Brooks, A. Gelman, G. Jones, and X.-L. Meng, eds.,  _Handbook of Markov Chain Monte Carlo_ (Chapman and Hall/CRC, New York, 2011).
  * Collaborations [2017] F.-L. Collaborations (Fermi-LAT, DES), The Astrophysical Journal 834, 110 (2017), arxiv:1611.03184 [astro-ph.HE] .
  * de Salas and Widmark [2021] P. F. de Salas and A. Widmark, Rept. Prog. Phys. 84, 104901 (2021), arxiv:2012.11477 [astro-ph.GA] .
  * Gross and Vitells [2010] E. Gross and O. Vitells, The European Physical Journal C 70, 525 (2010), arxiv:1005.1891 [physics.data-an] .
  * Steigman  _et al._ [2012] G. Steigman, B. Dasgupta, and J. F. Beacom, Physical Review D 86, 023506 (2012), arxiv:1204.3622 [hep-ph] .
  * Heisig [2020] J. Heisig, Mod. Phys. Lett. A 36, 2130003 (2020), arxiv:2012.03956 [astro-ph.HE] .
  * Zhao  _et al._ [2021] M.-J. Zhao, K. Fang, and X.-J. Bi, Physical Review D 104, 123001 (2021), arxiv:2109.04112 [astro-ph.HE] .
  * Boschini  _et al._ [2018] M. Boschini, S. Della Torre, M. Gervasi, G. La Vacca, and P. Rancoita, Advances in Space Research 62, 2859 (2018), arxiv:1704.03733 [astro-ph.SR] .
  * Kappl [2016] R. Kappl, Computer Physics Communications 207, 386 (2016), arxiv:1511.07875 [astro-ph.SR] .