# Constraints on the spatially dependent cosmic-ray propagation model from
Bayesian Analysis

Meng-Jie Zhao  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  Key
Laboratory of Particle Astrophysics, Institute of High Energy Physics, Chinese
Academy of Sciences, Beijing 100049, China  Xiao-Jun Bi  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

The energy spectra of primary and secondary cosmic rays (CR) generally harden
at several hundreds of GeV, which can be naturally interpreted by propagation
effects. We adopt a spatially dependent CR propagation model to fit the
spectral hardening, where a slow-diffusion disk (SDD) is assumed near the
Galactic plane. We aim to constrain the propagation parameters with the
Bayesian parameter estimation based on a Markov chain Monte Carlo sampling
algorithm. The latest precise measurements of carbon spectrum and B/C ratio
are adopted in the Bayesian analysis. The $\rm{{}^{10}Be/^{9}Be}$ and Be/B
ratios are also included to break parameter degeneracies. The fitting result
shows that all the parameters are well constrained. Especially, the thickness
of the SDD is limited to 0.4-0.5 kpc above and below the Galactic plane, which
could be the best constraint for the slow-diffusion region among similar
works. The $\bar{p}/p$ ratio and amplitude of CR anisotropy predicted by the
SDD model are consistent with the observations, while the predicted high-
energy electron and positron fluxes are slightly and significantly lower than
the observations, respectively, indicating the necessity of extra sources.

##  I INTRODUCTION

The Galactic cosmic-ray (CR) propagation can be described by the diffusion
process due to the random scattering by MHD waves in the interstellar medium
(ISM). Thus, the properties of the magnetic field turbulence in the ISM
determine the CR diffusion. The turbulence in the Galactic disk is mainly
generated by stellar feedback (such as the supernova explosions), while in the
outer halo, the matter is much rarefied and the turbulence is driven by CRs
themselves [1]. As the turbulence origin and ISM properties are both different
between the Galactic disk and halo, the CR diffusion in the Galaxy is very
likely to be spatially dependent [2, 3].

The spatially dependent diffusion is also supported by observations. The TeV
gamma-ray halos around some middle-aged pulsars indicate that the diffusion
coefficients around these pulsars are more than two orders of magnitude
smaller than the average in the Galaxy [4, 5]. If the slow-diffusion zone is
common in the ISM around Galactic pulsars, the average diffusion coefficient
in the Galactic disk would be significantly suppressed [6]. The spatial
magnetic-energy spectrum in the Galaxy also implies that the magnetic field
turbulence in the Galactic disk is much stronger than in the halo [7], which
means that the diffusion coefficient in the Galactic disk could be
significantly smaller. Besides, the local CR anisotropy predicted by the
standard model [8] is much larger than that observed by multiple experiments
[9], while the assumption of a slower CR diffusion in the Galactic disk could
give it an explanation [10].

The spatially dependent diffusion could account for the well-known spectral
hardening of CRs [11, 12]. Recent experiments such as PAMELA [13], AMS-02 [14,
15, 16], ATIC-2 [17], CALET [18, 19], DAMPE [20, 21], CREAM [22] have all
discovered the spectral hardening at several hundreds of GeV for most primary
and secondary CR nuclei, which cannot be explained by the simplest CR
injection and propagation models. Under the spatially dependent diffusion, the
energy exponent of the diffusion coefficient can also be spatially dependent,
which may explain the spectral hardening. Other possible approaches to
interpreting the spectral hardening include the CR injection reflecting
nonlinear or time-dependent diffusive-shock-acceleration [23], the nonlinear
effects in CR propagation [24], and local anomalies due to nearby sources [25]
or different transport in the Local Bubble [26].

We expect to use a spatially dependent diffusion model to explain anomalies
such as the spectral hardening problem without introducing nearby sources or
spectral breaks in the injection spectra. Our model consists of a slow-
diffusion disk near the Galactic plane and a fast-diffusion halo more extended
in vertical. Variations on diffusion properties can lead to twice spectral
hardening for the secondary nuclei compared with the primaries, which is
consistent with the AMS-02 observations (see Fig. 84 in [27]). We adopt the
Bayesian analysis based on a Markov chain Monte Carlo (MCMC) sampling
algorithm to constrain the model parameters, which is meaningful for depicting
the CR diffusion pattern in the Galaxy. We also notice that the proton and
helium spectra both have a ”knee” around 10TV as recently founded by NUCLEON
[28] and DAMPE [20, 21]. We suppose that other mechanisms give this feature
and only focus on the spectra below this energy.

This paper is organized as follows. In Section II, we introduce our CR
propagation model, the data sets used for analysis, and the method of
parameter inference. In Section III, we present our fitting results in terms
of the parameter posterior probability distributions and the best-fit values.
We discuss the parameter constraints by comparing the best-fit spectra and the
observations. In Section IV, we adopt our spatially dependent propagation
model and the fitting results to predict the $\bar{p}/p$ ratio, the electron
and positron spectra, as well as the anisotropy amplitude of CR nuclei, and
test if they are consistent with these observations. Section V is the
conclusion.

##  II CALCULATIONS

###  II.1 CR propagation model

The propagation equation of Galactic CRs is generally expressed by:

$\displaystyle{\frac{\partial\psi}{\partial t}}=$ $\displaystyle
q(x,p)+\nabla\cdot(D_{xx}\nabla\psi-V_{c}\psi)+{\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
V_{c})\psi]-{\frac{\psi}{\tau_{f}}}-{\frac{\psi}{\tau_{r}}}$

where $\psi$ is the density of CR particles per unit momentum, $q(x,p)$ is the
source term, $D_{xx}$ is the spatial diffusion coefficient, $V_{c}$ is the
convection velocity, $D_{pp}$ is the momentum space diffusion coefficient,
$\dot{p}\equiv dp/dt$ describes ionization and Coulomb losses, $\tau_{f}$ is
the time scales for collisions off gas nuclei, and $\tau_{r}$ is the time
scales for radioactive decay.

Supernova remnants (SNRs) are believed to be the main sources of Galactic CRs,
where charged particles are accelerated by shock waves. As suggested by the
shock acceleration theory, the injection spectrum of primary CRs is assumed to
be a power law as $q\propto R^{-\nu}$, where $R$ is the rigidity of CRs. A
low-energy break $R_{\rm br}$ is needed for all the nuclei to fit the observed
low-energy spectral bumps. The spectral indices below and above the break are
denoted with $\nu_{0}$ and $\nu_{1}$, respectively.

The scattering of CR particles on randomly moving MHD waves leads to
stochastic acceleration, which is described in the transport equation as
diffusion in momentum space $D_{pp}$. Alfvén velocity $V_{a}$ is introduced as
a characteristic velocity of weak propagating in a magnetic field, which is
related to the spatial coefficient $D_{xx}$:

$D_{xx}D_{pp}=\frac{4p^{2}V_{a}^{2}}{3\delta(4-\delta)(4-\delta^{2})\omega}\,.$
(2)

Figure 1: The diffusion coefficient $D$ changes with rigidity $R$ and spacial
vertical position $z$ where we assume ($N$=8, $h$=0.5 kpc).

We introduce a slow-diffusion disk (SDD) model, where the diffusion
coefficient near the Galactic plane is suppressed. As shown in Fig. 1, the SDD
model defines the diffusion coefficient $D_{xx}$ by

$D_{xx}(R,z)=aD_{0}\beta^{\eta}({\frac{R}{R_{0}}})^{b\delta}$ (3a)
$a=1+(\xi-1){\rm exp}[-({\frac{z}{h}})^{N}]$ (3b) $b=1+(\xi_{\delta}-1){\rm
exp}[-({\frac{z}{h}})^{N}]\,,$ (3c)

where $\beta=v/c$ is the particle velocity divided by the speed of light, and
the low-energy random-walk process is shaped by the factor $\beta^{\eta}$.
Here $\eta\neq 1$ is introduced to improve the calculated B/C ratio at low
rigidity to fit the observations. The scale factors $a$ and $b$ define the
spatial variation of the diffusion coefficient. The scale factor $a$ changes
the normalization at the reference rigidity $R_{0}=4$ GV, while $b$ changes
the slope index. The parameter $h$ describes the thickness of this region, and
$N$ describes the smoothness of the variation.

We define the dimensionless rigidity parameter $\rho\equiv R/R_{0}$ so that
the diffusion coefficient in the innermost disk ($z\simeq 0$, also regarded as
local coefficient) can be expressed by $D_{i}=\xi
D_{0}\rho^{\xi_{\delta}\delta}$. The change of slope index from the halo to
the disk can be defined as $\Delta\equiv\delta(1-\xi_{\delta})$. According to
a similar spatially dependent propagation model given by Ref. [29], the
effective height of the slow-diffusion (SD) region can be defined as
$\Lambda(\rho)\equiv h+(L-h)\xi\rho^{-\Delta}$, where $L$ is the overall size
of the Galactic diffusion halo. We can see that the effective height tends to
$h$ at high energies.

We assume $N=\infty$ in the main text, which means that the diffusion has a
two-zone feature. We discuss the effect of $N$ in Appendix B.2. For the two-
zone diffusion scenario, the local CR fluxes and ratios can be approximately
expressed by the following forms as given by Ref. [29]:

$\psi_{\rm{pri}}(0)\propto\rho^{-\nu}\frac{\Lambda(\rho)}{D_{i}}=\frac{\rho^{-\nu-\delta}}{D_{0}}[\frac{h}{\xi}\rho^{\Delta}+(L-h)]\,,$
(4)
$\frac{\psi_{\rm{sec}}(0)}{\psi_{\rm{pri}}(0)}\propto\frac{\Lambda(\rho)}{D_{i}}=\frac{\rho^{-\delta}}{D_{0}}[\frac{h}{\xi}\rho^{\Delta}+(L-h)]\,,$
(5)
$\frac{\psi_{u}(0)}{\psi_{s}(0)}\propto\frac{\sqrt{D_{i}\tau_{r}}}{\Lambda(\rho)}=\frac{\sqrt{\xi
D_{0}\tau_{r}\rho^{\xi_{\delta}\delta}}}{h+(L-h)\xi\rho^{-\Delta}}\,,$ (6)

where $\psi_{\rm{pri}}(0)$ is primary flux,
${\psi_{\rm{sec}}(0)}/{\psi_{\rm{pri}}(0)}$ is secondary to primary flux ratio
and ${\psi_{u}(0)}/{\psi_{s}(0)}$ is unstable to stable flux ratio. It can be
seen from Eq. (4) that the primary CR spectra can be described by the
superposition of a hard component and a soft component. When the particle
rigidity gets larger, the hard component gets dominant, and the spectral index
changes from $(\nu+\delta)$ to $(\nu+\delta\xi_{\delta})$. This feature also
appears in secondary/primary ratios such as B/C. Unstable/stable ratio such as
$\rm{{}^{10}Be/^{9}Be}$ has a different form related with the decay lifetime
$\tau_{r}$ and can help to break the degeneracy between $\Lambda(\rho)$ and
$D_{i}$. The Be/B ratio is quite complex and shows a similar feature to Eq.
(5) at high energies and Eq. (6) at low energies, which is discussed in detail
in Appendix A.

In Appendix B.1, we prove that $\xi_{\delta}$ is required to be very small by
the fitting procedure (See Table 3), which means that the SDD model prefers an
energy-independent $D_{i}$ in the disk. Thus, we fix $\xi_{\delta}=0$ in the
main text to simplify the fitting procedure. In this case, the diffusion
coefficient in the disk at low energies could be larger than that in the outer
halo, which may not be reasonable considering the origin of the ISM
turbulence. Thus, we further add a constraint that the former must be always
smaller than or equal to the latter, thus the scale factors $\xi_{\delta}$ and
$\xi$ are supposed to be equal to 1 below GeV energy.

To solve the propagation equation, we adopt the numerical GALPROP
v56111Current version available at https://galprop.stanford.edu/ [30, 31]. The
information of the interstellar medium (gas, radiation and magnetic fields)
are considered in GALPROP, which makes the calculated results more realistic.
We revise the differencing scheme in the solver by adopting the finite volume
method, which is necessary for the spatially dependent diffusion coefficient
[32].

Solar modulation significantly changes the CR spectra below $\sim$20 GeV. To
account for the solar modulation effect, we adopt a simple force-field
approximation [33], where the strength is described by the solar modulation
potential $\phi$. According to Table 1, all the AMS-02 and ACE-CRIS (except
$\rm{{}^{10}Be/^{9}Be}$ [34]) measurements used in this paper were taken
during the same period (2011/05-2016/05), hence we use a uniform $\phi$ to
modulate carbon, B/C and Be/B. For $\rm{{}^{10}Be/^{9}Be}$ ACE-CRIS data [34]
taken during (1997/08/27-1999/04/09), we use $\phi-(0.1~{}\rm{GV})$ as an
approximation. The $0.1$ GV difference between the two periods is indicated by
the long-term observations of the neutron monitor devices [35].

If we use the default values of cross-section given by GALPROP, there is a
conflict between old statistics of $\rm{{}^{10}Be/^{9}Be}$ [34, 36] and the
newly measured Be/B by AMS-02 [16], as the former predicts a thin diffusion
halo ($\sim 3$ kpc) while the latter predicts a thicker one ($\sim 6$ kpc)
[37, 38, 39]. Since the current uncertainties of cross-section are quite large
($\sim 10-30\%$) [38, 40, 41], a modification of normalization in the entire
energy range $XS$ and the low-energy slope $XS_{\delta}$ of the beryllium
production cross-section [40, 42] can be introduced to reconcile the conflict:

$\displaystyle\sigma=\sigma^{\rm{default}}\cdot
XS\cdot\begin{cases}(\frac{E_{\rm{kin}/n}}{E_{\rm{kin}/n}^{\rm{thresh}}})^{XS_{\delta}}\,,\quad&E_{\rm{kin}/n}<E_{\rm{kin}/n}^{\rm{thresh}}\\\
1\,,\quad&\rm{otherwise}\end{cases}\,.$ (7)

As the cross-section models predict a break around 5 GeV/n energy and a flat
behaviour above it, we choose $E_{\rm{kin}/n}^{\rm{thresh}}=5$ GeV/n.

In summary, the group of free parameters are

$\displaystyle\bm{\theta}=\\{D_{0},\delta,L,V_{a},\eta,\xi,h,A_{c},\nu_{0},\nu_{1},$
$\displaystyle R_{\rm{br}},\phi,XS,XS_{\delta}\\}$

where $D_{0}$, $\delta$, $\eta$, $\xi$, and $h$ are the parameters describing
the diffusion coefficient, $L$ is the half-width of the total diffusive halo,
$V_{a}$ is Alfvén velocity, $A_{c}$ is the abundance of carbon when fixing the
abundance of proton to $1.06*10^{6}$, $\nu_{0}$, $\nu_{1}$, and $R_{\rm{br}}$
are the first and second indices and the break rigidity of overall injection
parameters, respectively, $\phi$ is modulation potential for the AMS-02
measurements, and $XS$ and $XS_{\delta}$ are the modification parameters of
the beryllium production cross-section.

###  II.2 Data sets

According to Ref. [43, 44, 45], the data group of H-He and heavy nucleons
($Z>2$) have different constraints on propagation and injection parameters. We
only use the heavy nuclei data to give a self-consistent constraint on the
propagation model. We assume that all the heavy nuclei share the same
injection parameters ($\nu_{0},\nu_{1},R_{\rm{br}}$) and use the carbon flux,
$\rm{{}^{10}Be/^{9}Be}$ ratio, B/C ratio and Be/B ratio to constrain
parameters. The latter two ratios are mainly decided by the
$\rm(C-N-O)\rightarrow(Be-B)$ series.

Besides the precise measurements of carbon flux, B/C, and Be/B ratios from
AMS-02 [15, 16], other data are also included for better parameter
constraints, which are listed in Table 1. For the carbon flux, we use the
CALET [19], NUCLEON [28], and CREAM-II [46] measurements to cover the multi-
TeV energy region and the ACE-CRIS measurements [47] to cover the MeV energy
region. The low-energy B/C ratio is constrained by the ACE-CRIS data [47]. The
data of $\rm{{}^{10}Be/^{9}Be}$ ratio is taken from ACE-CRIS [34] and ISOMAX
[36].

Table 1: Data Used in This Analysis Experiment | Energy Range | data points | Reference
---|---|---|---
B/C
AMS-02(2011/05-2016/05) | 2-2100 GV | 67 | [16]
ACE-CRIS(2011/05-2016/05) | 0.07-0.17 GeV/n | 6 | [47]
Be/B
AMS-02(2011/05-2016/05) | 2-2100 GV | 67 | [16]
$\mathbf{{}^{10}Be/^{9}Be}$
ISOMAX(1998/08/04-08/05) | 0.5-1.6 GeV/n | 2 | [36]
ACE-CRIS(1997/08/27-1999/04/09) | 0.08-0.14 GeV/n | 3 | [34]
C
NUCLEON(2015/07-2017/06) | 250-17000 GeV/n | 10 | [28]
CREAM-II(2005/12-2006/01) | 85-7500 GeV/n | 9 | [46]
CALET(2015/10-2019/10)*1.27222a multiplication of 1.27 is described in [19] to get aligned with AMS-02 | 10-1700 GeV/n | 22 | [19]
AMS-02(2011/05-2016/05) | 0.4-1200 GeV/n | 68 | [15]
ACE-CRIS(2011/05-2016/05) | 0.06-0.2 GeV/n | 7 | [47]
Voyager1-HET(2012-2015) | 0.02-0.13 GeV/n | 8 | [48]
B
Voyager1-HET(2012-2015) | 0.02-0.11 GeV/n | 8 | [48]
Be
Voyager1-HET(2012-2015) | 0.06-0.1 GeV/n | 2 | [48]

CRs have fully unimpeded access to Voyager 1, free of solar modulation and
local interstellar modulation [49, 50]. Thus, the Voyager 1 data can be
regarded as $\phi=0$ GV. We adopt the carbon, boron, and beryllium fluxes of
Voyager 1 [48] to break the entanglement between $\phi$ and other parameters.
Electron, positron, and antiproton fluxes are modulated differently, which
will be further discussed in Sec. IV.2 and IV.1.

###  II.3 Bayesian inference and MCMC

From the Bayes theorem, the posterior probability distribution of the model
parameters is

$P(\bm{\theta}|{\rm D})=\frac{P({\rm D}|\bm{\theta})P(\bm{\theta})}{P({\rm
D})}\,,$ (8)

where $\rm{D}$ denotes the used data,
$P(\rm{D}|\bm{\theta})=\mathcal{L}(\bm{\theta})$ is the likelihood function,
and $P(\bm{\theta})$ is the prior distribution. The quantity $P(\rm{D})$ in
the denominator of Eq. (8) is the Bayesian evidence, which is a normalizing
constant being independent of the model parameters $\bm{\theta}$ and can be
neglected in parameter inference.

MCMC methods are widely used in Bayesian inference and are powerful to sample
the high-dimensional parameter space for CR propagation models [51, 52, 53,
43]. We use the public code CosmoMC333See https://cosmologist.info/cosmomc/.
as a generic MC sampler to explore parameter space [54, 55], which uses the
Metropolis-Hastings algorithm to generate samples from the posterior
distribution [54]. It also provides tools for analyzing the posterior
distribution and making confidence contour plots [56].

##  III FITTING RESULTS

###  III.1 Posterior distributions of parameters

Figure 2: The 1-D and 2-D distributions of the transport parameters

We first give an expectation about the disentanglement of the important
parameters. As mentioned in Section II.1, the degeneracy between
$\Lambda(\rho)$ and $D_{i}$ can be broken by fitting the data of B/C,
$\rm{{}^{10}Be/^{9}Be}$, and Be/B ratios. Moreover, we have assumed that the
diffusion coefficient in the SD must be smaller or equal to that in the halo,
which means the scale factors $\xi_{\delta}$ and $\xi$ are equal to 1 at low
energies (below $\sim 1$ GeV). According to Eq. (5) and (6), the B/C ratio at
low energies can be approximated by

$\frac{\psi_{\rm{sec}}(0)}{\psi_{\rm{pri}}(0)}\propto\frac{L}{D_{0}\rho^{\delta}}\,,$
(9)

while the $\rm{{}^{10}Be/^{9}Be}$ ratio at low energies can be approximated by

$\frac{\psi_{u}(0)}{\psi_{s}(0)}\propto\frac{\sqrt{D_{0}\rho^{\delta}}}{L}\,.$
(10)

Thus, the degeneracy between $D_{0}$ and $L$ can be broken by the low-energy
data combining Eq. (9) and (10). According to the definitions of $D_{i}$ and
$\Lambda$, the two parameters $\xi$ and $h$ can be successively determined as
long as $D_{0}$ and $L$ are well constrained.

Fig. 2 is the triangle plot of the fitting results, which 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 not shown here for simplicity. As we have expected
above, all the parameters have well-behaved distributions. Weak anti-
correlations in $\delta-V_{a}$ , $\delta-\eta$, and $\delta-D_{0}$ can also be
seen, which are consistent with the standard models (e.g., Fig. 3 in Ref.
[43]).

Table 2: The prior range, best-fit values and posterior 95% range of all parameters in SDD model Parameter | Prior range | Best-fit values | posterior 95% range
---|---|---|---
$D_{0}(10^{28}cm^{2}s^{-1})$ | [0,10.0] | 3.379 | [2.986,4.023]
$\delta$ | [0.2,1.0] | 0.583 | [0.557,0.608]
$L$(kpc) | [1.0,20.0] | 4.743 | [4.323,5.625]
$V_{a}$(km/s) | [0,50] | 19.718 | [17.130,21.706]
$\eta$ | [-3,2] | -1.299 | [-1.518,-1.099]
$\xi$ | [0,4.5] | 1.153 | [0.965,1.277]
$h$(kpc) | [0,2.0] | 0.468 | [0.406,0.515]
$A_{c}(10^{3})$444abundance of proton $A_{p}$ is $1.06*10^{6}$, and the normalization of proton flux at 100 GeV is $4.204*10^{-9}\rm{cm}^{-2}\rm{s}^{-1}\rm{sr}^{-1}\rm{MeV}^{-1}$ | [3.1,3.65] | 3.337 | [3.316,3.377]
$\nu_{0}$ | [0.4,2.0] | 1.266 | [1.076,1.549]
$\nu_{1}$ | [2.2,2.5] | 2.373 | [2.364,2.381]
$R_{\rm{br}}$(GV) | [0,5] | 1.749 | [1.430,2.214]
$\phi$(GV) | [0.5,1.0] | 0.782 | [0.763,0.793]
$XS$ | [0.7,1.1] | 0.973 | [0.968,0.986]
$XS_{\delta}$ | [-0.2,0.2] | 0.0513 | [0.0418,0.0689]
$\chi^{2}_{\rm{min}}/n_{\rm{dof}}$ | - | 167.55/265 | -

The results are also summarized in Table 2, where we list the prior ranges,
best-fit values, and posterior 95% ranges for all the parameters. The total
halo height and SDD thickness are well constrained to
$L=4.743_{-0.420}^{+0.882}$ kpc, $h=0.468_{-0.062}^{+0.047}$ kpc,
respectively. The halo height $L$ is consistent with the earlier findings
using standard models: Ref. [39] got $L\in[3.7,6]$ kpc, Ref. [52] found
$L=4\pm 1$ kpc in a pure diffusion/reacceleration model, and Ref. [57] got
$L=5.4\pm 1.4$ kpc. However, our results are different from another spatially
dependent diffusion mode, Ref. [11]. They obtained a significantly thicker SDD
($h\sim 0.87$ kpc) and larger uncertainty for the halo height, which could be
due to the lack of the current precise measurements555We also notice that they
calculated $\chi^{2}$ by using fewer data on nucleons (above 45 GeV/n) and B/C
ratio (above 2 GeV/n), which could also result in a loose constraint..

The constrained SDD thickness $h$ is one order of magnitude smaller than the
halo height $L$ and a bit thicker than the Galactic disk ($\sim 0.2$ kpc)
where the CR sources concentrates, which may be explained by the convection of
turbulent energy.

The best-fit slope index of diffusion in the halo $\delta=0.583$ is quite
larger than the Kolmogorov type (1/3) and a bit larger than the Iroshnikov-
Kraichnan type (1/2). The normalization scale factor $\xi$ is around 1, which
means that the diffusion coefficient in the disk is close to that in the halo
at the reference rigidity of 4 GV.

The constrained modulation potential $\phi\simeq 0.782$ GV is in agreement
with those found by Ref. [47]. The two parameters modifying the production
cross-section of beryllium, $XS$ and $XS_{\delta}$, are within the experiment
uncertainties of $\sigma_{\rm{norm}}\simeq 0.2$ and $\sigma_{\rm{slope}}\simeq
0.15$ as given in Ref. [39].

We note that the best-fit value $\chi^{2}_{\rm{min}}$ of the SDD model is
167.55 which seems a too good fit compared with the degree of freedom of 265.
The reason may be that we have added the systematic errors of these
measurements in quadrature with the statistical errors to get the total
errors, but have not taken into account the correlation among systematic
uncertainties in the calculation. Covariance matrices may be needed to
properly take into account those data uncertainties [58].

###  III.2 Nucleon fluxes and ratios

The energy spectrum of carbon nuclei is shown in Fig. 3. A clear spectral
hardening can be seen, as predicted by Eq. (4). At lower energies, the
spectrum fits well with the ACE-CRIS measurement, and the LIS also fits well
with the Voyager 1 data, which means that the solar modulation potential is
reasonable. We have also drawn a 95% range band derived from the uncertainties
of parameters to show that the spectrum is strictly constrained.

Figure 3: Carbon spectrum calculated with the best-fit parameters (green solid
line) and the $2\sigma$ confidence interval (green band), compared with the
experimental data, including AMS-02 [15], ACE-CRIS [47], Voyager1 [48], CALET
[19](multiplied by 1.27), CREAM-II[46], and NUCLEON[28]. The LIS of carbon are
drawn with green dotted line. Top: An wide-range comparison between the model
and the experimental data. Bottom: Details around 200 GeV/n where the
hardening appears.

Eq. (5) predicts features of secondary/primary ratios similar to that of the
primary spectrum. We can see from Fig. 4 that our calculation for the B/C
ratio fits well with AMS-02 measurement in the entire energy range. The
calculation also shows a smooth hardening above 100 GeV/n, which could be
confirmed by more precise measurements at higher energies in the future. At
lower energies, the B/C ratio calculated by the SDD model is slightly higher
than the ACE-CRIS measurement. We suppose that a modification on boron cross-
section or/and solar modulation may explain this difference.

Figure 4: B/C ratio calculated with the best-fit parameters (green line) and
the $2\sigma$ confidence interval, compared with the experimental data: AMS-02
[59] and ACE-CRIS [47].

From Fig. 5, the uncertainties reported from ACE-CRIS and ISOMAX are so large
that they cannot give strong constraints on halo height $L$, while more
precise experiments on $\rm{{}^{10}Be/^{9}Be}$ in the future may be required,
as well as additional information from Be/B ratio.

Figure 5: $\rm{{}^{10}Be/^{9}Be}$ ratio calculated with the best-fit
parameters (green line) and the $2\sigma$ confidence interval (green band),
compared with the experimental data: ISOMAX [36] and ACE-CRIS [34].

From Fig. 6, our calculation fits well with the AMS-02 measurement. As
analyzed in Appendix A, the Be/B ratio experiences a transition from the
decay-dominated low-energy region to the diffusion-dominated high-energy
region, which is unique in constraining transport parameters. According to
Ref. [39, 41], the most important region for constraining halo height666To be
exact, effective height $\Lambda$ in SDD model. should be $10\sim 100$ GV
where the Be/B ratio has the strongest dependency on it. At energies below 10
GV, the degeneracy with the modification on cross-section and other factors
would complicate the fitting. We also find that the modification of slope
$XS_{\delta}$ on the cross-section is needed, otherwise, it is hard to
reproduce the feature below 10 GV.

Figure 6: Be/B ratio calculated with the best-fit parameters (green line) and
the $2\sigma$ confidence interval (green band), compared with the experimental
data AMS-02[16].

##  IV PREDICTIONS

AMS-02 has provided observations on antimatter particles like antiprotons and
positrons, which are crucial for predicting dark matter (DM) particles [60,
61]. The antiproton-to-proton ratio [62] and the positron fraction [63]
predicted by the standard models are significantly lower than the AMS-02
measurements at high energies. The excesses could be explained by introducing
DM, while we need to get a proper assessment of antiparticle background
firstly. The spatially dependent propagation model may generate higher
antiparticle fluxes compared with the standard model [11].

Besides, the CR anisotropy predicted by the standard propagation models is
significantly higher than the observations. As the anisotropy is proportional
to the local diffusion coefficient, our SDD model can effectively suppress the
CR anisotropy.

###  IV.1 Antiprotons

To predict the antiproton flux, or $\bar{p}/p$ ratio for the SDD model, we
keep all the parameters in Table 2 fixed to the best-fit values and change the
injection parameters of proton and helium to fit the latest measurements of
AMS-02 [27]. We use the default nuclear scaling routine given by GALPROP to
get the hadronic cross-sections and calculate antiproton flux.

The production mechanism of secondary nucleons like beryllium is different
from that of antiprotons. The nuclear fragmentation keeps the energy per
nucleon of secondary particles the same as that of primary particles, while
the antiproton spectrum is the convolution of the interstellar spectra and the
differential cross-section [64], leading to lower energies compared with the
primary particles ($E_{\bar{p},\rm{max}}\sim\sqrt{E_{p}/2}$). So the hardening
energy of beryllium and antiproton should be around $\sim 200$ GV and $\sim
10$ GV respectively.

The result in Fig. 7 shows that the SDD model can give a good explanation to
the measured $\bar{p}/p$ ratio without introducing an extra source for
antiproton, such as the dark matter annihilation. The predicted hardening
above 10 GV by the propagation effect can explain the antiproton excess.
Furthermore, the injection and cross-section uncertainties may further improve
the fitting result [62]. Here we use a smaller modulation potential
$\phi_{\bar{p}}=0.44$ GV to modulate the low-energy region of antiproton flux,
as the charge of antiproton is opposite from proton and should be effected
differently by solar activities [65].

Figure 7: $\bar{p}/p$ ratio predicted by the best-fit parameters obtained in
Sec. III (green line) and the $2\sigma$ confidence interval (green band),
compared with the experimental data of AMS-02 [27].

###  IV.2 Electrons and Positrons

As AMS-02 collaboration showed [60, 66] that the electron and positron spectra
each could be well described by the sum of two components. Here we use the SDD
model to test if the extra components could be originated from spatially
dependent diffusion.

To calculate electrons and positrons fluxes in the SDD model, we reuse the
fitted proton and helium fluxes from Section IV.1 and choose the Pshirkov-ASS
model [67] in GALPROP to describe the Galactic magnetic field, which plays a
crucial role in lepton energy losses. It has been shown [68, 69] that SNRs
could contribute primary electrons, while a few secondary electrons and
positrons can be produced from the decay of charged pions and kaons created in
collisions of cosmic-ray particles with gas. Besides, electron-positron pairs
can also be produced by PWNe or DM annihilation/decay. To subtract the
possible contributions from PWNe, dark matter and secondary electrons, we
calculate the primary electrons by using a subtracted form
$(\Phi_{e^{-}}-\Phi_{e^{+}})$ between electron and positron fluxes to fit the
injection of electron [70], where we use much recent AMS-02 data [27]. Since
the charge of the electron is opposite from proton and $e^{+}e^{-}$ have
smaller masses than that of the proton, we consider a different $\phi$ to
describe how they are affected by solar modulation.

Figure 8: Positron and primary electron spectra predicted by the best-fit
parameters obtained in Sec. III (green line) and the $2\sigma$ confidence
interval (green band), compared with the experimental data of AMS-02 [27].

In Fig. 8, we can see that the primary electron spectrum fits the data well
below 100 GeV/n. As predicted by Eq. (4), a hard component at high energies
will also give rise to electron spectral hardening, but the energy losses in
the high-energy region are even stronger and steepen the electron spectrum so
rapidly, making the result failed to reproduce the AMS-02 measurement above
100 GeV/n. Extra components above 100 GeV/n may be needed and may be located
within relatively short distances, as electrons have a much shorter lifetime
with strong energy losses. Young and nearby SNRs may produce a much harder
component that gives rise to the excess [71, 72].

We also find that a solar modulation potential larger than 1 GV is required
for positrons to fit the data, which has also been noted by Ref. [73]. The
positron spectrum predicted by the SDD model is harder than that from standard
models since we have introduced a hard component similar to the antiproton
spectrum analyzed in Sec. IV.1, but the overall flux is still significantly
lower than the AMS-02 measurement. The missing flux may come from some extra
contribution of positrons sources, including nearby pulsars (or PWNe) or DM
particles.

###  IV.3 Anisotropy

In the diffusion approximation, the anisotropy is dominated by the radial
streaming of the CR fluxes, and its amplitude $\hat{A}$ is computed as

$\hat{A}=\frac{3D\left|\nabla\psi\right|}{v\psi}\propto\frac{D_{i}}{\Lambda(\rho)}({\rm
local})\,.$ (11)

We note that the anisotropy amplitude has an anti-correlation with the B/C
ratio (Eq. (5)). As the B/C ratio gets harder at higher energies, the
anisotropy amplitude gets softer and its index would change from $\delta$ to
nearly 0, which is different from standard models that assume an unchanged
slope index.

In Fig. 9, we have drawn the anisotropy amplitude together with the
uncertainties given by the fitting procedure in Sec. III, which is consistent
with the current observations (see Ref. [9] and references therein). The grey
dotted line is the anisotropy calculated with the standard propagation model,
which is obviously higher than the observations. Besides, Ref. [29] pointed
out that the anisotropy may be reduced in all energies if one accounts for a
proper radial dependence for the diffusion coefficient.

Figure 9:  CR anisotropy amplitude predicted by best-fit parameters obtained
in Sec. III (black line) and $2\sigma$ confidence interval (green band),
compared with the experimental data given in Ref. [9]. An example of standard
propagation model calculations (grey dotted line) are shown for reference
[47].

The fine structures such as a sudden decrease of the anisotropy amplitude
around $10^{5}$GeV are hard to be explained by the SDD model. There are other
possibilities to further explain these features, for example, the presence of
nearby sources of CRs [74].

##  V SUMMARY

According to the assumption of different diffusion environments in the disk
and the halo, we assume an SDD model for Galactic CR propagation, which has
two different diffusion regions. In the inner region, the diffusion is slow
and the slope index equals zero, while in the outer region, the diffusion
coefficient is similar to the standard models. The SDD model has the advantage
of naturally producing the observed high-energy spectral hardening of both the
primary and secondary CR spectra without assuming a high-energy spectral break
in the injection spectra or introducing nearby sources. By fitting the latest
precise measurement of CR spectra and ratios, the parameters of the SDD model
can be constrained in turn.

We perform a full Bayesian analysis based on an MCMC sampling algorithm to get
the posterior parameter distributions. We use the carbon data as the primary
CR spectrum, B/C ratio as the primary-to-secondary ratio, and
$\rm{{}^{10}Be/^{9}Be}$ ratio as the unstable-to-stable ratio. The Be/B ratio
is also adopted, the low-energy part of which can make up for the
$\rm{{}^{10}Be/^{9}Be}$ data (which is not precise enough) and help to break
parameter degeneracies.

The fitting result shows that all the parameters are well constrained. Our
main finding is the well-constrained thickness of the SD region,
$h=0.468_{-0.062}^{+0.047}$ kpc at 95% confidence, which could be the first
precise estimation on this parameter. The fitted size of this region is a bit
larger than the typical height of the Galactic disk, implying the convection
of turbulent energy in the direction vertical to the Galactic plane. Other
propagation parameters are consistent with those obtained in standard models,
such as the height and diffusion coefficient of the outer diffusive halo.

Based on the fitting result, we also predict the $\bar{p}/p$ ratio,
$e^{-}e^{+}$ fluxes, and the amplitude of CR anisotropy. We find that the
$\bar{p}/p$ ratio fits well with the AMS-02 data, and no extra component is
needed at high energies. The primary $e^{-}$ flux above 100 GeV is lower than
the observation, which is likely due to the spectral fluctuation from nearby
SNRs. The predicted positron flux is significantly lower than the AMS-02 data
in a wide energy range, so primary positron sources are required, such as
pulsars or DM. The anisotropy amplitude predicted by the SDD model fits the
experimental data well in general, which is a remarkable advantage compared
with standard models.

Other indirect observations such as CR-induced diffuse $\gamma$-ray emission
[75] and molecular clouds emission [76] could be used to study the spatial
variation of the diffusion coefficient in the Galaxy. The most prominent
difference between the SDD and standard model is the vertical distribution of
CRs, as the former predicts a more rapid decrease of the CR density within the
SD region. The CR vertical distribution can be traced by the $\gamma$-ray
emission of intermediate-velocity clouds located at various distances away
from the Galactic plane [76]. More precise measurements of the clouds in the
future may give a crucial test to the spatially dependent propagation model.

###### Acknowledgements.

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

## References

  * Erlykin and Wolfendale [2002] A. D. Erlykin and A. W. Wolfendale, J. Phys. G 28, 2329 (2002).
  * Evoli and Yan [2014] C. Evoli and H. Yan, Astrophys. J. 782, 36 (2014), arXiv:1310.5732 [astro-ph.HE] .
  * Evoli  _et al._ [2018] C. Evoli, P. Blasi, G. Morlino, and R. Aloisio, Phys. Rev. Lett. 121, 021102 (2018), arXiv:1806.04153 [astro-ph.HE] .
  * Abeysekara  _et al._ [2017] A. U. Abeysekara _et al._ (HAWC), Science 358, 911 (2017), arXiv:1711.06223 [astro-ph.HE] .
  * Aharonian  _et al._ [2021] F. Aharonian _et al._ (LHAASO), Phys. Rev. Lett. 126, 241103 (2021), arXiv:2106.09396 [astro-ph.HE] .
  * Hooper  _et al._ [2017] D. Hooper, I. Cholis, T. Linden, and K. Fang, Phys. Rev. D 96, 103013 (2017), arXiv:1702.08436 [astro-ph.HE] .
  * Han [2017] J. Han, Annual Review of Astronomy and Astrophysics 55, 111 (2017), https://doi.org/10.1146/annurev-astro-091916-055221 .
  * Blasi and Amato [2012] P. Blasi and E. Amato, JCAP 01, 011, arXiv:1105.4529 [astro-ph.HE] .
  * Ahlers and Mertsch [2017] M. Ahlers and P. Mertsch, Prog. Part. Nucl. Phys. 94, 184 (2017), arXiv:1612.01873 [astro-ph.HE] .
  * Guo  _et al._ [2016] Y.-Q. Guo, Z. Tian, and C. Jin, Astrophys. J. 819, 54 (2016), arXiv:1509.08227 [astro-ph.HE] .
  * Feng  _et al._ [2016] J. Feng, N. Tomassetti, and A. Oliva, Phys. Rev. D 94, 123007 (2016), arXiv:1610.06182 [astro-ph.HE] .
  * Guo and Yuan [2018] Y.-Q. Guo and Q. Yuan, Phys. Rev. D 97, 063008 (2018), arXiv:1801.05904 [astro-ph.HE] .
  * Adriani  _et al._ [2011] O. Adriani _et al._ (PAMELA), Science 332, 69 (2011), arXiv:1103.4055 [astro-ph.HE] .
  * Consolandi [2016] C. Consolandi (AMS), in  _25th European Cosmic Ray Symposium_ (2016) arXiv:1612.08562 [astro-ph.HE] .
  * Aguilar  _et al._ [2017] M. Aguilar _et al._ (AMS), Phys. Rev. Lett. 119, 251101 (2017).
  * Aguilar  _et al._ [2018] M. Aguilar _et al._ (AMS), Phys. Rev. Lett. 120, 021101 (2018).
  * Panov  _et al._ [2009] A. D. Panov _et al._ , Bull. Russ. Acad. Sci. Phys. 73, 564 (2009), arXiv:1101.3246 [astro-ph.HE] .
  * Adriani  _et al._ [2019] O. Adriani _et al._ (CALET), Phys. Rev. Lett. 122, 181102 (2019), arXiv:1905.04229 [astro-ph.HE] .
  * Adriani  _et al._ [2020] O. Adriani _et al._ , Phys. Rev. Lett. 125, 251102 (2020), arXiv:2012.10319 [astro-ph.HE] .
  * An  _et al._ [2019] Q. An _et al._ (DAMPE), Sci. Adv. 5, eaax3793 (2019), arXiv:1909.12860 [astro-ph.HE] .
  * Alemanno  _et al._ [2021] F. Alemanno, Q. An, P. Azzarello, F. Barbato, P. Bernardini, X. Bi, M. Cai, E. Catanzani, J. Chang, D. Chen, and et al., Physical Review Letters 126, 10.1103/physrevlett.126.201102 (2021).
  * Yoon  _et al._ [2017] Y. S. Yoon _et al._ , Astrophys. J. 839, 5 (2017), arXiv:1704.02512 [astro-ph.HE] .
  * Ptuskin  _et al._ [2013] V. Ptuskin, V. Zirakashvili, and E.-S. Seo, Astrophys. J. 763, 47 (2013), arXiv:1212.0381 [astro-ph.HE] .
  * Blasi  _et al._ [2012] P. Blasi, E. Amato, and P. D. Serpico, Phys. Rev. Lett. 109, 061101 (2012), arXiv:1207.3706 [astro-ph.HE] .
  * Thoudam and Horandel [2012] S. Thoudam and J. R. Horandel, Mon. Not. Roy. Astron. Soc. 421, 1209 (2012), arXiv:1112.3020 [astro-ph.HE] .
  * Ohira and Ioka [2010] Y. Ohira and K. Ioka, in  _32nd International Cosmic Ray Conference_ (2010) arXiv:1011.4405 [astro-ph.HE] .
  * Aguilar  _et al._ [2021] M. Aguilar _et al._ (AMS), Phys. Rept. 894, 1 (2021).
  * Gorbunov  _et al._ [2019] N. Gorbunov _et al._ , Adv. Space Res. 64, 2546 (2019), arXiv:1809.05333 [astro-ph.IM] .
  * Tomassetti [2012] N. Tomassetti, Astrophys. J. Lett. 752, L13 (2012), arXiv:1204.4492 [astro-ph.HE] .
  * Strong and Moskalenko [1998] A. W. Strong and I. V. Moskalenko, Astrophys. J. 509, 212 (1998), arXiv:astro-ph/9807150 .
  * Strong  _et al._ [2000] A. W. Strong, I. V. Moskalenko, and O. Reimer, Astrophys. J. 537, 763 (2000), [Erratum: Astrophys.J. 541, 1109 (2000)], arXiv:astro-ph/9811296 .
  * Fang  _et al._ [2018a] K. Fang, X.-J. Bi, P.-F. Yin, and Q. Yuan, Astrophys. J. 863, 30 (2018a), arXiv:1803.02640 [astro-ph.HE] .
  * Gleeson and Axford [1968] L. J. Gleeson and W. I. Axford, Astrophys. J. 154, 1011 (1968).
  * Yanasak  _et al._ [2001] N. E. Yanasak, M. E. Wiedenbeck, R. A. Mewaldt, A. J. Davis, A. C. Cummings, J. S. George, R. A. Leske, E. C. Stone, E. R. Christian, T. T. von Rosenvinge, W. R. Binns, P. L. Hink, and M. H. Israel, Astrophys. J. 563, 768 (2001).
  * Ghelfi  _et al._ [2017] A. Ghelfi, D. Maurin, A. Cheminet, L. Derome, G. Hubert, and F. Melot, Adv. Space Res. 60, 833 (2017), arXiv:1607.01976 [astro-ph.HE] .
  * Hams  _et al._ [2004] T. Hams _et al._ , Astrophys. J. 611, 892 (2004).
  * Evoli  _et al._ [2020] C. Evoli, G. Morlino, P. Blasi, and R. Aloisio, Phys. Rev. D 101, 023013 (2020), arXiv:1910.04113 [astro-ph.HE] .
  * De La Torre Luque  _et al._ [2021] P. De La Torre Luque, M. N. Mazziotta, F. Loparco, F. Gargano, and D. Serini, JCAP 03, 099, arXiv:2101.01547 [astro-ph.HE] .
  * Weinrich  _et al._ [2020a] N. Weinrich, M. Boudaud, L. Derome, Y. Genolini, J. Lavalle, D. Maurin, P. Salati, P. Serpico, and G. Weymann-Despres, Astron. Astrophys. 639, A74 (2020a), arXiv:2004.00441 [astro-ph.HE] .
  * Weinrich  _et al._ [2020b] N. Weinrich, Y. Génolini, M. Boudaud, L. Derome, and D. Maurin, Astron. Astrophys. 639, A131 (2020b), arXiv:2002.11406 [astro-ph.HE] .
  * Tomassetti [2015] N. Tomassetti, Phys. Rev. C 92, 045808 (2015), arXiv:1509.05776 [astro-ph.HE] .
  * Korsmeier and Cuoco [2021] M. Korsmeier and A. Cuoco, Phys. Rev. D 103, 103016 (2021), arXiv:2103.09824 [astro-ph.HE] .
  * Jóhannesson  _et al._ [2016] G. Jóhannesson _et al._ , Astrophys. J. 824, 16 (2016), arXiv:1602.02243 [astro-ph.HE] .
  * Schroer  _et al._ [2021] B. Schroer, C. Evoli, and P. Blasi, Phys. Rev. D 103, 123010 (2021), arXiv:2102.12576 [astro-ph.HE] .
  * Evoli  _et al._ [2019] C. Evoli, R. Aloisio, and P. Blasi, Phys. Rev. D 99, 103023 (2019), arXiv:1904.10220 [astro-ph.HE] .
  * Ahn  _et al._ [2009] H. S. Ahn _et al._ , Astrophys. J. 707, 593 (2009), arXiv:0911.1889 [astro-ph.HE] .
  * Yuan [2019] Q. Yuan, Sci. China Phys. Mech. Astron. 62, 49511 (2019), arXiv:1805.10649 [astro-ph.HE] .
  * 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, Astrophys. J. 831, 18 (2016).
  * Strauss  _et al._ [2013] R. D. Strauss, M. S. Potgieter, S. E. S. Ferreira, H. Fichtner, and K. Scherer, The Astrophysical Journal 765, L18 (2013).
  * Luo  _et al._ [2015] X. Luo, M. Zhang, M. Potgieter, X. Feng, and N. V. Pogorelov, The Astrophysical Journal 808, 82 (2015).
  * Masi [2016] N. Masi, Nuovo Cim. C 39, 282 (2016).
  * Putze  _et al._ [2010] A. Putze, L. Derome, and D. Maurin, Astron. Astrophys. 516, A66 (2010), arXiv:1001.0551 [astro-ph.HE] .
  * 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] .
  * Lewis and Bridle [2002] A. Lewis and S. Bridle, Phys. Rev. D 66, 103511 (2002), arXiv:astro-ph/0205436 .
  * Lewis [2013] A. Lewis, Phys. Rev. D 87, 103529 (2013), arXiv:1304.4473 [astro-ph.CO] .
  * Lewis [2019] A. Lewis, arXiv e-prints , arXiv:1910.13970 (2019), arXiv:1910.13970 [astro-ph.IM] .
  * Trotta  _et al._ [2011] R. Trotta, G. Jóhannesson, I. V. Moskalenko, T. A. Porter, R. R. d. Austri, and A. W. Strong, Astrophys. J. 729, 106 (2011), arXiv:1011.0037 [astro-ph.HE] .
  * Derome  _et al._ [2020] L. Derome, D. Maurin, P. Salati, M. Boudaud, Y. Génolini, and P. Kunzé, PoS ICRC2019, 054 (2020).
  * Aguilar  _et al._ [2016a] M. Aguilar _et al._ (AMS), Phys. Rev. Lett. 117, 231102 (2016a).
  * Aguilar  _et al._ [2019a] M. Aguilar _et al._ (AMS), Phys. Rev. Lett. 122, 041102 (2019a).
  * Aguilar  _et al._ [2016b] M. Aguilar _et al._ (AMS), Phys. Rev. Lett. 117, 091103 (2016b).
  * Giesen  _et al._ [2015] G. Giesen, M. Boudaud, Y. Génolini, V. Poulin, M. Cirelli, P. Salati, and P. D. Serpico, JCAP 09, 023, arXiv:1504.04276 [astro-ph.HE] .
  * Moskalenko and Strong [1998] I. V. Moskalenko and A. W. Strong, Astrophys. J. 493, 694 (1998), arXiv:astro-ph/9710124 .
  * di Mauro  _et al._ [2014] M. di Mauro, F. Donato, A. Goudelis, and P. D. Serpico, Phys. Rev. D 90, 085017 (2014), [Erratum: Phys.Rev.D 98, 049901 (2018)], arXiv:1408.0288 [hep-ph] .
  * Cholis  _et al._ [2016] I. Cholis, D. Hooper, and T. Linden, Phys. Rev. D 93, 043016 (2016), arXiv:1511.01507 [astro-ph.SR] .
  * Aguilar  _et al._ [2019b] M. Aguilar _et al._ (AMS), Phys. Rev. Lett. 122, 101101 (2019b).
  * Pshirkov  _et al._ [2011] M. S. Pshirkov, P. G. Tinyakov, P. P. Kronberg, and K. J. Newton-McGee, Astrophys. J. 738, 192 (2011), arXiv:1103.0814 [astro-ph.GA] .
  * Evoli  _et al._ [2021] C. Evoli, E. Amato, P. Blasi, and R. Aloisio, Phys. Rev. D 103, 083010 (2021), arXiv:2010.11955 [astro-ph.HE] .
  * Di Mauro  _et al._ [2020] M. Di Mauro, F. Donato, and S. Manconi, arXiv e-prints , arXiv:2010.13825 (2020), arXiv:2010.13825 [astro-ph.HE] .
  * Jin  _et al._ [2016] C. Jin, Y. Q. Guo, and H. B. Hu, Chin. Phys. C 40, 015101 (2016), arXiv:1504.06903 [astro-ph.HE] .
  * Fang  _et al._ [2017] K. Fang, B.-B. Wang, X.-J. Bi, S.-J. Lin, and P.-F. Yin, Astrophys. J. 836, 172 (2017), arXiv:1611.10292 [astro-ph.HE] .
  * Fang  _et al._ [2018b] K. Fang, X.-J. Bi, and P.-f. Yin, Mon. Not. Roy. Astron. Soc. 478, 5660 (2018b), arXiv:1706.03745 [astro-ph.HE] .
  * Orlando [2018] E. Orlando, Mon. Not. Roy. Astron. Soc. 475, 2724 (2018), arXiv:1712.07127 [astro-ph.HE] .
  * Yuan  _et al._ [2021] Q. Yuan, B.-Q. Qiao, Y.-Q. Guo, Y.-Z. Fan, and X.-J. Bi, Front. Phys. (Beijing) 16, 24501 (2021), arXiv:2007.01768 [astro-ph.HE] .
  * Ackermann  _et al._ [2012] M. Ackermann _et al._ (Fermi-LAT), Astrophys. J. 750, 3 (2012), arXiv:1202.4039 [astro-ph.HE] .
  * Tibaldo  _et al._ [2015] L. Tibaldo _et al._ , Astrophys. J. 807, 161 (2015), arXiv:1505.04223 [astro-ph.HE] .
  * Maurin  _et al._ [2001] D. Maurin, F. Donato, R. Taillet, and P. Salati, Astrophys. J. 555, 585 (2001), arXiv:astro-ph/0101231 .

##  Appendix A The Be/B grammage

To analyze how the Be/B ratio is related to transport parameters, we firstly
apply resolution on boron and beryllium nucleons respectively as
$\rm{Be=^{10}Be+^{9}Be+^{7}Be}$ and $\rm{B=^{11}B+^{10}B}$. The collision with
gas $(p+^{11}\rm{B)\to(^{10}Be+^{9}Be+^{7}Be})$ and decay process
${}^{10}\rm{Be\to^{10}B}+e^{-}$ also make effects, though not as important as
the main production $\rm{(C-N-O)\to(Be-B)}$.

For a more comprehensive review, one can refer to [77, 37], following which we
can found that stable elements (${}^{10}\rm{B}:I_{a}$) with contribution from
unstable ones (${}^{10}\rm{Be}:I_{b}$) has this form:

$\frac{I_{a}}{X}=\sum_{a^{\prime}>a}\frac{I_{a^{\prime}}\sigma_{a^{\prime}\to
a}}{m}+\frac{I_{b}V_{c}}{\mu v}[\Delta\coth\frac{V_{c}\Delta
L}{2D}-\coth\frac{V_{c}L}{2D}]\,.$ (12)

With diffuse-dominated grammage $X=\mu vL/2D$, decay-dominated grammage
$X_{d}=\mu v\tau/\sqrt{4D\tau}$ and $\Delta=\sqrt{1+4D/V_{c}^{2}\tau}$, when
unstable isotopes decay on a timescale shorter than $4D/V_{c}^{2}$ (usually
below 100GV), $V_{c}\Delta\to\sqrt{4D/\tau}$ ,$V_{c}\to 0$ (since we assumed
few or no galactic wind here), the second term on the right side of Eq. (12)
becomes $I_{b}\sqrt{4D/\tau}/{\mu v}$ and the third term becomes $2I_{b}D/L\mu
v$. Now we can write all isotopes below:

$\begin{cases}{\rm{}^{9}Be,^{7}Be}:\frac{I_{a}}{X}=\sum_{a^{\prime}>a}\frac{I_{a^{\prime}}\sigma_{a^{\prime}\to
a}}{m}\\\
{\rm{}^{10}Be}:\frac{I_{a}}{X_{d}}=\sum_{a^{\prime}>a}\frac{I_{a^{\prime}}\sigma_{a^{\prime}\to
a}}{m}(<100{\rm GeV})\\\
{\rm{}^{11}B}:\frac{I_{a}}{X}+\frac{I_{a}\sigma_{a}}{m}=\sum_{a^{\prime}>a}\frac{I_{a^{\prime}}\sigma_{a^{\prime}\to
a}}{m}\\\
{\rm{}^{10}B}:\frac{I_{a}}{X}=\sum_{a^{\prime}>a}\frac{I_{a^{\prime}}\sigma_{a^{\prime}\to
a}}{m}+\frac{I_{b}}{X_{d}}-\frac{I_{b}}{X}(<100{\rm GeV})\par\end{cases}\\\ $
(13)

By combining them we can further calculate different ratios as

$\begin{cases}{\rm\frac{{}^{9}Be+^{7}Be}{{}^{11}B}}:X\propto\frac{L}{2D}\\\
{\rm\frac{{}^{10}Be}{{}^{11}B}}:X_{d}\propto\frac{\tau}{\sqrt{4D\tau}}(<100{\rm
GeV})\\\
{\rm\frac{{}^{10}Be}{{}^{10}B}}:\cfrac{1}{\frac{X\sum\cfrac{I_{a^{\prime}}\sigma_{a^{\prime}\to
a}}{m}}{X_{d}\sum\cfrac{I_{a^{\prime}}\sigma_{a^{\prime}\to
b}}{m}}+\cfrac{X}{X_{d}}-1}\propto\frac{X_{d}}{X}=\frac{\sqrt{D\tau}}{L}(<100{\rm
GeV})\end{cases}\\\ $ (14)

So the total ratio of Be/B should be a mixture of all ratio forms in Eq. (14),
which shows the B/C-like (L/D) feature at high energy,
$\rm{{}^{10}Be/^{9}Be}$-like ($\sqrt{D}/L$ and $1/\sqrt{D}$) feature at low
energy. The important feature for breaking the degeneracy of $D-L$ is
$\rm{{}^{10}Be/^{9}Be}$-like. As available $\rm{{}^{10}Be/^{9}Be}$ ratio
measurements have large uncertainties, an introduction of precise Be/B ratio
is preferred. It is emphasized in [77] that the
$\rm{{}^{10}Be\to^{10}B}+e^{-}$ channel contributes up to 10% of total Boron
flux and cannot be neglected, but $\rm{{}^{10}Be}$ fluxes make up only 10% of
total beryllium at low energy and there are cross-section uncertainties, which
result in a complicated problem.

##  Appendix B Above 20 GV

###  B.1 specified $h$

Before fitting the free parameters using all measurements according to Section
II.2, we firstly estimate how the SDD model fits the spectral hardening at
high energy, by using AMS-02 carbon and boron fluxes, together with the Be/B
ratio. All experiment points are taken above 20 GV, where the influences of
low energy power-law break, solar modulation, cross-section uncertainties and
non-relativity effect should be lowest.

To make the fitting much simpler we fix $\phi=0.8\rm{GV},\eta=-0.5$, taken
from Yuan’s paper [47] as a reference of standard models. We fix $N=8$ to give
a rapid smoothness from halo to disk, and choose the disk thickness among
specified values $h=\\{0.3,0.5,0.8,1.0\\}$ kpc as we have predicted a strong
degeneracy of $h/\xi$ from Eq. (4). So the free parameters are
$\bm{\theta}=\\{D_{0},\delta,L,V_{a},\xi,\xi_{\delta},A_{c},\nu\\}$, consist
of 6 transport parameters and 2 injection parameters.

Figure 10: The 1-D and 2-D distributions of the transport parameters,
different colors represent different values of specified $h$: 0.3 kpc(green),
0.5 kpc(grey), 0.8 kpc(red), 1.0 kpc(blue)

Figure 11: Nucleon spectra calculated with the best-fit parameters above 20 GV, compared with the experimental data AMS-02 [15, 16], different colors represent different values of specified $h$. Top: carbon. Bottom: boron and beryllium.  Figure 12: CR anisotropy amplitude predicted by best-fit parameters, compared with the experimental data given in [9], different colors represent different values of specified $h$: 0.3 kpc (red), 0.5 kpc (blue), 0.8 kpc (green), 1.0 kpc (black).  Table 3: The best-fit values and posterior 95% range of all parameters in SDD model Parameter | h=0.3 kpc | h=0.5 kpc | h=0.8 kpc | h=1.0 kpc
---|---|---|---|---
$D_{0}(10^{28}\rm{cm^{2}s^{-1}})$ | 4.648 [3.222,5.353] | 4.451 [3.845,4.792] | 4.076 [3.573,5.274] | 3.823 [3.376,5.289]
$\delta$ | 0.683 [0.631,0.707] | 0.700 [0.662,0.729] | 0.708 [0.663,0.733] | 0.730 [0.701,0.784]
$L$(kpc) | 11.414 [4.747,13.703] | 10.574 [7.418,12.144] | 9.446 [6.458,18.076] | 9.240 [6.387,19.176]
$V_{a}$(km/s) | 1.101 [0,13.409] | 2.376 [0,8.862] | 1.251 [0,12.969] | 1.314 [0,9.701]
$\xi$ | 0.2716 [0.1960,0.4251] | 0.6181 [0.4953,0.7230] | 1.3334 [0.9695,1.5242] | 1.760 [0.6905,1.7528]
$\xi_{\delta}$ | 0.002127 [0,0.03082] | 0.02207 [0,0.08207] | 0.02454 [0,0.06832] | 0.03866 [0.01614,0.1612]
$A_{c}(10^{-3})$777abundance of proton $A_{p}$ is $1.06*10^{6}$, and the normalization of proton flux at 100 GeV is $4.204*10^{-9}\rm{cm^{-2}s^{-1}sr^{-1}MeV^{-1}}$ | 3.279 [3.243,3.319] | 3.290 [3.254,3.321] | 3.288 [3.259,3.335] | 3.277 [3.249,3.325]
$\nu$ | 2.352 [2.334,2.368] | 2.359 [2.341,2.376] | 2.363 [2.342,2.379] | 2.360 [2.339,2.380]
$\Delta$888$\Delta=\delta(1-\xi_{\delta})$ | 0.651 | 0.662 | 0.669 | 0.673
$\chi^{2}_{\rm{min}}/n_{\rm{dof}}$ | 47.51/118 | 47.29/118 | 46.78/118 | 47.13/118

Table 3 and Fig. 10 shows the MCMC result of 4 kinds of specified $h$. We
notice that all these $\chi^{2}_{\rm{min}}$s have similar values, thus the
thickness $h$ does not show a strong preference to one of them. As half-
thickness of inner disk gets larger, $\delta,\xi,\xi_{\delta}$ all get
increased and a strong anti-correlation is shown between $\xi$ and
$\xi_{\delta}$. It seems that the Alfvén velocity $V_{a}$ and slope index
scale factor $\xi_{\delta}$ are all converged to zero. The effect of
reacceleration does not significantly affect the hardening and could be
ignored since we are now focused on the energy region above 20 GV. Moreover,
the diffusion coefficient in the disk prefers an energy-independent type as
$\xi_{\delta}\sim 0$. In Fig. 10 we also find that the strong degeneracy of
$D_{0}/L$ worsens the constraining of other parameters.

From Fig. 11 we find that the best-fit values of these specified $h$ have
almost the same results fitting good with AMS-02 measurements, while a few
dispersion appears above $10^{4}$ GV. To estimate the hardening speed of the
slow component $\rho^{\Delta}$, we calculate $\Delta$ in Table 3 and find that
this value keeps nearly unchanged when $h$ become larger, proving that
$\Delta$ should be important in reproducing similar hardening features.

To estimate more differences of these specified $h$, we further predict the
anisotropy amplitude in Fig. 12. As $h$ gets larger, the amplitude in the
entire energy range gets larger (except 1.0 kpc), and the disk thickness
$h\sim 0.5$ kpc fits best with these experiments. Eq. (11) has predicted a
simple relation that $\hat{A}\sim\xi D_{0}/h=\\{4.207,5.502,6.793,6.728\\}$,
which could explain those features.

###  B.2 specified smooth factor $N$

Figure 13: The diffusion coefficient changes with different smooth factor $N$,
where we assume $h=1$ kpc, $\xi=0.05$.

We have introduced the smooth factor $N$ in Eq. (3a), which is used to
describe the change of diffusion coefficient from the innermost disk to the
outer halo. In Fig. 13 we show the spatial dependence of diffusion
coefficient, which changes more rapidly when smooth factor $N$ gets larger.
When $N\sim\infty$, this change becomes a step function around the thickness
$\pm h$ as Ref. [29] used this kind of spatial dependent model.

We choose $h=1.0$ kpc to give a fixed thickness of slow region and choose the
smooth factor among specified values $N=\\{2,4,8\\}$. Other data sets are all
referred to from Appendix. B.1. Fig. 16 shows the fitting results of 3 kinds
of specified $N$. As smooth factor get larger, $\xi,\xi_{\delta}$ all get
increased and the anti-correlation also appears between $\xi$ and
$\xi_{\delta}$. The serious degeneracy of $D_{0}/L$ still worsens the
constraining of other parameters.

From Fig. 14 we find that the best-fit values of these specified $N$ have
almost the same results fitting good with AMS-02 measurements, while large
smooth factor $N$ can give stronger hardening above $10^{3}$ GV. The local CRs
mainly travel from disk to halo and diffuse backward to the solar system. The
real spectra may not be just the superposition of two components like Eq. (4)
and Eq. (5) but include more intermediate states. The presence of intermediate
components can be used to explain how the hardening changes with $N$.

To estimate more differences of these specified $N$, we further predicted the
anisotropy amplitude in Fig. 15. As $N$ gets larger, the amplitude in the
entire energy range gets larger as well, and the smooth factor $N\sim 4$ fits
best with these experiments. Eq. (11) has predicted a simple relation that
$\hat{A}\sim\xi D_{0}/h=\\{3.216,5.737,6.728\\}$, which could explain those
features shown in Fig. 15.

Figure 14: Carbon spectrum calculated with the best-fit parameters, compared
with the experimental data AMS-02 [15] above 20 GV. Different colors represent
different values of specified $N$.  Figure 15: CR anisotropy amplitude
predicted by best-fit parameters, compared with the experimental data given in
[9]. Different colors represent different values of specified $N$: 2(red),
4(blue), 8(green).  Figure 16: The 1-D and 2-D distributions of the transport
parameters. Different colors represent different values of specified $N$:
2(grey), 4(red), 8(blue)