# Beam polarization effects on Z-boson pair production at electron-positron
colliders: a full one-loop analysis

Mehmet Demirci<EMAIL_ADDRESS>Department of Physics, University of
Wisconsin–Madison, Madison, Wisconsin 53706, USA  Department of Physics,
Karadeniz Technical University, Trabzon 61080, Turkey  A. Baha Balantekin
<EMAIL_ADDRESS>Department of Physics, University of Wisconsin–Madison,
Madison, Wisconsin 53706, USA

###### Abstract

We present high-precision predictions for Z boson-pair production via
electron-positron collisions by taking into account a full set of one-loop
order scattering amplitudes, i.e., electroweak (EW) corrections together with
soft and hard QED radiation. We provide a detailed numerical evaluation, from
full EW corrections to pure QED corrections, specifically focusing on the
effect of initial beam polarization on production rate. The left-right
asymmetry and angular distributions are also presented. The radiative
corrections are largely affected by initial beam polarizations. We get an
improvement as around three times with completely polarized beams of
$e^{-}_{L}$ and $e^{+}_{R}$. We find that the radiative corrections can
sizably modify the production rate, typically yielding a total relative
correction up to a few tens of a percent. This imply that the full EW
corrections are required for $e^{-}e^{+}\to ZZ$ to match with percent level
accuracy.

Electroweak radiative corrections, Z boson, Beam polarization, Standard Model,
electron-positron collider

##  I Introduction

The Standard Model (SM) [1, 2, 3] has been well established as a self-
consistent gauge theory. Starting with the discovery of the Higgs boson [4, 5]
a decade ago, we have an increasing confidence on the consistency of the
Standard Model (SM) from the LHC experiments [6]. Despite its many successes,
there are remaining questions (such as the hierarchy problem, neutrino masses,
dark matter, the strong CP problem, and the generation of a baryon asymmetry)
which require beyond the SM (BSM) physics. Current and future experiments will
continue not only make sensitivity tests of the SM, but also will search for
BSM physics to meet the challenges listed above.

Proton-proton colliders have significant backgrounds which makes it harder to
explore rare processes. On the other hand, electron-positron colliders have
cleaner backgrounds and offer important opportunities to precisely observe the
interesting phenomena. Various proposals for such colliders exist: The
International Linear Collider (ILC) [7, 8, 9] is one of the proposed linear
colliders, which would be operated in the c.m. energies of
$\sqrt{s}=250-500\mbox{~{}GeV}$ (extendable up to a 1 TeV) with collision
modes of $e^{-}e^{+}$, $e^{-}e^{-}$, $e^{-}\gamma$, and $\gamma\gamma$. At the
ILC, it is also expected to use the electron and positron beams with $80\%$
and $30\%$ polarizations, respectively. Another linear collider is the Compact
Linear Collider (CLIC) [10], which is foreseen to be operated a TeV scale
high-luminosity capacity, in stages at energies of
$\sqrt{s}=380\mbox{~{}GeV}$, $1.5\mbox{~{}TeV}$ and $3\mbox{~{}TeV}$.
Furthermore, there are circular collider projects such as Circular Electron-
Positron Collider (CEPC) [11], and Future Circular Collider (FCC) [12]. The
CEPC is projected with a circumference of 100 km and to have a maximum of
$\sqrt{s}=240\mbox{~{}GeV}$. The FCC project [12], hosted in a 100 km tunnel
at CERN, is considering three collision types: FCC-ee, for electron-positron
collisions, FCC-hh, for hadron-hadron collisions, and FCC-eh, for electron-
hadron collisions. The FCC-ee [13], running at different energies to research
the production of Z, W, H and $t\bar{t}$, is expected to be a Higgs,
electroweak and top factory with high-luminosity. The above lepton colliders
are expected to offer the unique opportunities for precise measurements and
signs of new physics, appearing in the form of small deviations from the SM.

In parallel, high-precision predictions from theory are required to get more
precise measurements for clues of the BSM physics or more sensitivity tests of
the SM. This means that it is necessary to go beyond the lowest-order
calculations for the most production channels. At least, one-loop order
corrections in the production channels must be included to get adequate
accuracy. However, the extent to which higher-order computations beyond one-
loop order will be needed depends largely on the expected experimental
accuracy. In the present study, we focus on the complete one-loop order EW
corrections and beam polarization effects. The beam polarization has the
ability to improve the corresponding signal-to-background ratio along with the
sensitivity of the observables [14, 16, 15]. Furthermore, the CP-violating
couplings can be separated from CP-conserving ones using beam polarization
[16].

The main production reaction of the $Z$ boson pairs at the $e^{-}e^{+}$
colliders are $e^{-}e^{+}\rightarrow ZZ$, which include the Born-level
contributions from electron exchanges via $t$ and $u-$ channels. Since there
are no Yukawa or QCD contributions, the identification of the EW corrections
would be very clean. The process at Born order was studied by [17, 18] a long
ago. Next, one-loop EW radiative corrections were investigated by [19].
However, their results directly depend on the soft cut-off parameter, since
the hard photon correction is not taken into account. In addition to the above
work, in Ref. [20], the high energy behavior of the various helicity
amplitudes was investigated considering all of the MSSM contributions.

Furthermore, neutral gauge boson production via electron-positron collisions
have received a significant experimental [21, 23, 22] and theoretical [24, 25,
26, 27, 28] interest, motivated by the search for anomalous neutral gauge
boson self couplings. In particular, $Z$-pair events provide an opportunity
for the investigation of possible triple neutral gauge boson couplings,
$ZZ\gamma$ and $ZZZ$ [30, 29, 31], not allowed at tree level in the SM.

In the present work, we examine the $Z$ boson pair production via electron-
positron collisions in the SM, including a full set of one-loop order EW
radiative corrections, i.e., the EW corrections together with soft and hard
real photon emission. We provide a detailed numerical discussion with
particular emphasis on the effects of initial beam polarization. Moreover, we
consider the decomposition of the EW corrections into the pure QED photonic
corrections along with the corresponding counterterms and the remaining weak
corrections. The left-right asymmetry and angular distributions are also
presented at the Born and one-loop orders.

The rest of this paper is organized as follows. In Sec. II, we provide some
useful analytical expressions together with the relevant Feynman diagrams and
amplitudes. General shapes of the real photon emission and the virtual
corrections are also discussed. In Sec. III, we give briefly information on
initial beam polarization. In Sec. IV we give a set of input parameters used
in this study. In Sec. V, we present a detailed analysis of numerical results.
A comparison with the results of other approachescis also given. Finally, in
Sec. VI, we present summary and conclusions.

##  II Theoretical Setup for a cross section

We first present our notation regarding the Born amplitudes and one-loop EW
contributions. The relevant production process is expressed as

$e^{+}(p_{1},\sigma_{1})e^{-}(p_{2},\sigma_{2})\rightarrow
Z(k_{1},\lambda_{1})Z(k_{2},\lambda_{2}),$ (1)

where $\sigma_{1}$, $\sigma_{2}$, $\lambda_{1}$ and $\lambda_{2}$ denote the
helicities of initial positron and electron and outgoing $Z$ bosons,
respectively. The helicities take the values $\lambda_{1,2}=0,\pm 1$ and
$\sigma_{1,2}=\pm 1/2$. Polarization vectors of outgoing $Z$-bosons are
denoted by $\varepsilon_{\mu}(k_{1},\lambda_{1})$ and
$\varepsilon_{\nu}(k_{2},\lambda_{2})$. Neglecting the electron mass, the
momenta in the center of mass of the initial state system are given by

$\begin{split}p_{1}&=\frac{\sqrt{s}}{2}(1,0,0,-1),\\\
p_{2}&=\frac{\sqrt{s}}{2}(1,0,0,+1),\\\
k_{1}&=\frac{\sqrt{s}}{2}(1,-\kappa\sin\theta,0,-\kappa\cos\theta),\\\
k_{2}&=\frac{\sqrt{s}}{2}(1,+\kappa\sin\theta,0,+\kappa\cos\theta),\end{split}$
(2)

where

$\kappa=\sqrt{1-\frac{4M_{Z}^{2}}{s}},$ (3)

and $\sqrt{s}$ and $\theta$ denotes the center of mass energy and the
scattering angle, respectively. For further use, we also note the Mandelstam
variables:

$\begin{split}s&=(p_{1}+p_{2})^{2}=(k_{1}+k_{2})^{2},\\\
t&=(k_{1}-p_{1})^{2}=(k_{2}-p_{2})^{2},\\\
u&=(k_{2}-p_{1})^{2}=(k_{1}-p_{2})^{2}.\end{split}$ (4)

The Z-boson polarization vectors are given by

$\begin{split}\varepsilon^{\mu}_{+}(k,\lambda=+1)&=-\frac{1}{\sqrt{2}}(0,1,i,0),\\\
\varepsilon^{\mu}_{0}(k,\lambda=0)&=\frac{1}{M_{Z}}(k_{z},0,0,E),\\\
\varepsilon^{\mu}_{-}(k,\lambda=-1)&=\frac{1}{\sqrt{2}}(0,1,-i,0),\end{split}$
(5)

which obey the condition $\varepsilon^{\mu}k_{\mu}=0$. The sum over physical
polarization states of Z-boson is given by

$\sum_{\lambda}\varepsilon_{\mu}(k,\lambda)\varepsilon^{*}_{\nu}(k,\lambda)=-g_{\mu\nu}+\frac{k_{\mu}k_{\nu}}{M_{Z}^{2}}.$
(6)

The analytical and numerical evaluations are carried out by using the
following tools111 Using the same tools, many calculations have been carried
out with significant results (see, e.g., by one of us [33, 32, 34, 35], and
others [36, 37, 38, 39, 40]).. The Feynman diagrams and amplitudes are created
by using FeynArts [41, 42]. The algebraic evaluation for amplitudes is
provided technically in the same way as defined in Ref. [35]. Then, the
squaring of amplitudes, the simplifying of fermion chains and the numerical
computation are carried out by using FormCalc [43]. The scalar loop integrals
are calculated via LoopTools [43]. The phase-space integrations are calculated
by using the Monte Carlo integration algorithm Vegas, implemented in the CUBA
library [44]. We have checked the cross sections of Born and hard photon
radiation processes against the results obtained by using Whizard [45, 46] and
CalcHEP [47].

###  II.1 Born amplitudes

At lowest order, the process $e^{+}e^{-}\rightarrow ZZ$ includes the leading
contributions from $t$ and $u$-channel electron-exchange diagrams. The
corresponding Feynman diagrams are given in Fig. 1, where a diagram with a
Higgs field coupling to electrons is omitted. It is suppressed by a factor
$m_{e}/\sqrt{s}$ and thus can be neglected.

Figure 1: The lowest order Feynman diagrams contributing to
$e^{-}e^{+}\rightarrow ZZ$.

The explicit expressions are given for the Born matrix elements from these
diagrams as follows:

$\displaystyle\begin{split}{\cal
M}_{1}&=\frac{-ig^{2}}{4\big{[}t-m_{e}^{2}\big{]}c_{W}^{2}}\overline{v}(p_{1},m_{e})\not{\varepsilon}^{*}(k_{1})\big{[}c_{V}^{e}-c_{A}^{e}\gamma_{5}\big{]}\\\
&\times(\not{p}_{2}-\not{k}_{2}+m_{e})\not{\varepsilon}^{*}(k_{2})\big{[}c_{V}^{e}-c_{A}^{e}\gamma_{5}\big{]}u(p_{2},m_{e}),\end{split}$
(7)

$\displaystyle\begin{split}{\cal
M}_{2}&=\frac{-ig^{2}}{4\big{[}u-m_{e}^{2}\big{]}c_{W}^{2}}\overline{v}(p_{1},m_{e})\not{\varepsilon}^{*}(k_{2})\big{[}c_{V}^{e}-c_{A}^{e}\gamma_{5}\big{]}\\\
&\times(\not{p}_{2}-\not{k}_{1}+m_{e})\not{\varepsilon}^{*}(k_{1})\big{[}c_{V}^{e}-c_{A}^{e}\gamma_{5}\big{]}u(p_{2},m_{e}),\end{split}$
(8)

where $g=e/s_{W}$, $s_{W}=\sin\theta_{W}$ and
$c_{W}=\cos\theta_{W}=M_{W}/M_{Z}$. The vector and axial vector couplings are
defined as $c_{V}^{f}=I_{3}^{f}-2Q_{f}s_{W}^{2}$ and $c_{A}^{f}=I_{3}^{f}$ for
the fermion type $f$. For electron, these are $c_{V}^{e}=-1/2+2s_{W}^{2}$ and
$c_{A}^{e}=-1/2$. The Born-level total amplitude is

${\cal M}_{\text{Born}}=\sum_{i=1}^{2}{\cal M}_{i}.$ (9)

For arbitrary polarizations of the leptons and bosons, the differential cross
section reads

$\displaystyle\begin{split}\biggl{(}\frac{d\sigma}{d\Omega}\biggr{)}_{\text{Born}}^{e^{+}e^{-}\rightarrow
ZZ}=&\frac{\kappa}{64\pi^{2}s}\sum_{\sigma_{1,2},\lambda_{1,2}}\frac{1}{4}(1+2\sigma_{1}P_{1})(1+2\sigma_{2}P_{2})\\\
&\times|{\cal
M}_{\text{Born}}^{\sigma_{1},\sigma_{2},\lambda_{1},\lambda_{2}}|^{2},\end{split}$
(10)

where $P_{1}$ and $P_{2}$ are the polarization degrees of the incoming
electron and positron. The sum runs over all included boson polarizations.
Following the square of total amplitude and the summation over final
helicities, we have

$\begin{split}\biggl{(}\frac{d\sigma}{dt}\biggr{)}_{\text{Born}}^{e^{-}e^{+}\rightarrow
ZZ}=&\frac{g^{4}}{128\pi^{2}s^{2}}\frac{(c_{V}^{e})^{4}+(c_{A}^{e})^{4}+6(c_{V}^{e})^{2}(c_{A}^{e})^{2}}{c_{W}^{4}}\biggl{[}\frac{u}{t}+\frac{t}{u}+4M_{Z}^{2}\biggl{(}\frac{s}{ut}\biggr{)}-M_{Z}^{4}\biggl{(}\frac{1}{t^{2}}+\frac{1}{u^{2}}\biggr{)}\biggr{]}\end{split}$
(11)

where the electron mass $m_{e}$ is neglected for simplicity. This result is
consistent with the results in the literature [17]. The above coupling term
can be also written in term of chiral coupling constants as
$(c_{V}^{e})^{4}+(c_{A}^{e})^{4}+6(c_{V}^{e})^{2}(c_{A}^{e})^{2}=8(c_{L}^{e})^{4}+8(c_{R}^{e})^{4}$.
Then, the integrated cross section is obtained by

$\sigma_{\text{Born}}^{e^{-}e^{+}\rightarrow
ZZ}=\frac{1}{2!}\int_{t_{\text{min}}}^{t_{\text{max}}}\biggl{(}\frac{d\sigma}{dt}\biggr{)}_{\text{Born}}^{e^{-}e^{+}\rightarrow
ZZ}dt$ (12)

where the lower and upper bounds of the integral are defined as

$t_{\text{max},\text{min}}=M_{Z}^{2}-\frac{1}{2}s\bigl{(}1\pm\kappa\bigr{)}.$
(13)

###  II.2 One-loop EW radiative corrections

The radiative corrections consist of three parts. The first comes from virtual
loop corrections, the second one is soft photon radiation and the last one
corresponds to hard real photon radiation. We now discuss them for
$e^{-}e^{+}\rightarrow ZZ$ below.

####  II.2.1 Virtual corrections

We consider full EW $\mathcal{O}(\alpha)$ contributions at one-loop level for
the process $e^{-}e^{+}\rightarrow ZZ$. The virtual loop corrections for this
process come from three types of diagrams: vertex-type, self-energy, and box-
type. Here, we do not provide their explicit analytical expressions since they
are not particularly illuminating. Instead, we present a complete list of
Feynman diagrams created by the FeynArts. At one-loop order, the process (1)
has a total of 161 one-loop diagrams (112 vertex-type+34 self-energy+ 15 box-
type) as shown in Figs. 3 to 2.

Figure 2: The vertex-correction diagrams contributing to
$e^{-}e^{+}\rightarrow ZZ$. Figure 3: The self-energy diagrams contributing to
$e^{-}e^{+}\rightarrow ZZ$.

We use the following labels on the internal lines: $f$ stands for all fermions
in SM, $G$ is the charged Goldstone boson, and $u_{\pm}$ is the ghosts. The
gauge vector bosons ($\gamma$, $Z$, and $W^{\pm}$) are denoted by wavy lines,
and the Higgs and Goldstone bosons are denoted by the dashed lines. In
diagrams with two arrows on the same lines of the loop, particles move both
clockwise and counterclockwise.

The one-loop diagrams can topologically divided into $s$, $t$, and
$u$-channels with the mediator of gauge bosons ($\gamma$, $Z$, $W^{\pm}$),
Higgs boson ($h^{0}$), and charged Goldstone bosons ($G^{\pm}$). We can also
separate them into QED and Weak corrections. The diagrams of QED virtual
corrections are obtained by the possible virtual photon attachment and fermion
loop insertion to the Born-level diagrams.

First, in Fig. 2, we present the vertex-correction diagrams, which consist of
triangle corrections to $t-$channel electron exchange, bubbles, and triangle
vertices attached to the initial or final state via an mediator of $\gamma$,
$Z$, and $h^{0}$. They can also be divided into three classes. The first comes
from the vertex corrections $ZZA^{*}/Z^{*}$. Here, the off-shell field is
marked by asterisk. These are given with diagrams $(1)$-$(4)$ and $(6)$-$(8)$
in Fig. 2. The second is due to the vertex corrections $eeh^{0*}$ and it is
shown in diagram $(5)$ in Fig. 2. Both types are $s-$channel contributions.
The fermion loop contributions to the vertex of $eeh^{0*}$ are proportional to
the electron mass and hence is suppressed. The third type is the vertex
corrections $Ze\overline{e}^{*}$ in the $t-$ and $u-$ channels. The
$Ze\overline{e}^{*}$ correction is shown in diagrams $(9)-(20)$ of Fig. 2.

Second, in Fig. 3, we show the self-energy diagrams, which consist of all
possible loops of $e^{-}$ or $\nu_{e^{-}}$ with gauge bosons $\gamma$, $Z$, or
$W^{\pm}$ on the electron propagator via $t-$ and $u-$channels. Also, there
are all possible loops of fermions, $W^{\pm}$ and $G^{\pm}$ bosons on the
$\gamma h^{0}$ or $Zh^{0}$ mixing propagator in $s-$channel. Diagrams (4) and
(5) in Fig. 3 are obtained by QED correction, while others are weak
corrections.

Figure 4: The box-type diagrams contributing to $e^{-}e^{+}\rightarrow ZZ$.

Finally, Figure 4 shows the irreducible one-loop diagrams, i.e., the box-type
corrections. They consist of all possible box-loops of $e^{-}$, $\nu_{e^{-}}$,
vector bosons $\gamma$, $Z$, $W^{\pm}$, Higgs and Goldstone bosons. These are
mainly $t$ and $u$-channel contributions. Diagrams (1) and (7) in Fig. 4 are
obtained by QED correction, while others are weak corrections. Diagrams (2),
(4),(8) and (10) are due to neutral current correction and other diagrams are
due to charged current correction.

The total amplitude for the virtual contributions is defined as the summation
of triangle-type, self-energy and box-type corrections:

${\cal\delta M}_{\text{virt}}={\cal M}_{\triangle}+{\cal M}_{\bigcirc}+{\cal
M}_{\Box}.$ (14)

The differential cross section of the virtual corrections is calculated by

$\displaystyle\begin{split}\biggl{(}\frac{d\sigma}{d\Omega}\biggr{)}_{\text{virt}}^{e^{-}e^{+}\rightarrow
ZZ}=&\sum_{\sigma_{1,2},\lambda_{1,2}}\frac{1}{4}(1+2\sigma_{1}P_{1})(1+2\sigma_{2}P_{2})\\\
&\times\frac{\kappa}{64\pi^{2}s}2\text{Re}\bigl{[}{\cal
M}_{\text{Born}}^{*}{\cal\delta M}_{\text{virt}}\bigr{]},\end{split}$ (15)

where $|{\cal\delta M}_{\text{virt}}|^{2}$ is not included because its
contribution is so small that it can be neglected. The one-loop Feynman
diagrams, which form the virtual $\mathcal{O}(\alpha)$ corrections
${\cal\delta M}_{\text{virt}}$, have been calculated in ’t Hooft-Feynman
gauge. We use the on-shell (OS) renormalization scheme (see, Ref. [48]) to fix
all the renormalization constants.

The virtual contributions have ultraviolet (UV) and infrared (IR) divergences.
These divergences can be regularized by extending the dimensions of spinor and
spacetime manifolds to $D=4-2\epsilon$ [49] and adding a photon mass
parameter, respectively.

Figure 5: The counterterm diagrams.

We adopt all Feynman rules of the counterterms (shown in Fig. 5) and of the
renormalization conditions from Ref. [48]. The redefinition of parameters and
fields is carried out in the OS scheme. This turns the Lagrangian into a bare
and counterterms. After applying the renormalization procedure, a UV-finite
virtual contribution is achieved. Although, the soft IR singularity due to
virtual photonic loop corrections still exists. It is regulazed by giving the
photon a fictitious mass, $m_{\gamma}$. The virtual cross section is
independent of the UV regularization parameter
$C_{UV}=1/\epsilon-\gamma_{E}+\log(4\pi)$, but still a function of the IR
regularization parameter $m_{\gamma}$. From the Kinoshita-Lee-Nauenberg
theorem [50, 51]222This was also shown perturbatively in QED by Schwinger
[52]., it is cancelled in the limit $m_{\gamma}\rightarrow 0$ by adding the
real photon corrections. We have checked numerically that our results do not
depend on $m_{\gamma}$. After adding the virtual and real corrections, the
results are still collinear singular. This singularity comes from the initial
state radiation part. We use the phase space slicing method333We have
developed the necessary code for evaluating of collinear contributions and
implemented in FormCalc. [53, 54, 55, 56] to handle the collinear
singularities in the photon radiation off initial state.

####  II.2.2 Real corrections

Real photon emission gives rise to the kinematically different reaction from
$e^{+}e^{-}\rightarrow ZZ$, presented by

$e^{+}(p_{1},\sigma_{1})e^{-}(p_{2},\sigma_{2})\rightarrow
Z(k_{1},\lambda_{1})Z(k_{2},\lambda_{2})\gamma(k_{3},\lambda_{3}),$ (16)

where $k_{3}$ and $\lambda_{3}$ are four-momenta and helicity of the radiated
photon.

Figure 6: The Feynman diagrams for the real photon radiation.

We present the bremsstrahlung Feynman diagrams in Fig. 6, obtained from the
Born diagrams by adding a photon emitted from a lepton line. The lowest-order
of real-photon emission provides an $\mathcal{O}(\alpha)$ correction to
$e^{+}e^{-}\rightarrow ZZ$.

The differential cross section of $e^{+}e^{-}\rightarrow ZZ\gamma$ reads

$\begin{split}\biggl{(}\frac{d\sigma}{d\Omega_{3}}\biggr{)}_{\text{real}}^{e^{+}e^{-}\rightarrow
ZZ\gamma}=&\frac{1}{2s}\sum_{\sigma_{1,2},\lambda_{1,2,3}}\frac{1}{4}(1+2P_{1}\sigma_{1})(1+2P_{2}\sigma_{2})\\\
&\times|{\cal M}_{e^{+}e^{-}\rightarrow
ZZ\gamma}^{\sigma_{1},\sigma_{2},\lambda_{1},\lambda_{2},\lambda_{3}}|^{2}\end{split}$
(17)

with the 3-particle phase-space integral

$\begin{split}\int
d\Omega_{3}=\prod_{i=1}^{3}\int\frac{d^{3}\vec{k}_{i}}{(2\pi)^{3}2k_{i}^{0}}(2\pi)^{4}\delta\biggl{(}p_{1}+p_{2}-\sum_{j=1}^{3}k_{j}\biggr{)}.\end{split}$
(18)

According to the radiated photon energy
$E_{\gamma}=\sqrt{|\overrightarrow{k_{3}}|^{2}+m_{\gamma}^{2}}$, the
bremsstrahlung phase space can be separated into two regions: soft and hard.
Hence, the real photon radiation correction can be organized as

$d\sigma_{\text{real}}^{e^{-}e^{+}\rightarrow
ZZ\gamma}=d\sigma_{\text{soft}}(\Delta_{s})+d\sigma_{\text{hard}}(\Delta_{s})$
(19)

where $\Delta_{s}=\Delta E_{\gamma}/(\sqrt{s}/2)$ is the soft cutoff energy
parameter. The bremsstrahlung photon is named soft when $E_{\gamma}<\Delta
E_{\gamma}=\Delta_{s}\sqrt{s}/2$, whereas it is hard if $E_{\gamma}>\Delta
E_{\gamma}$. The soft photon correction is given by approximation [57, 48]

$\begin{split}d\sigma_{\text{soft}}=-d\sigma_{\text{Born}}\frac{\alpha
Q^{2}_{\ell}}{2\pi^{2}}\int_{|\overrightarrow{k_{3}}|\leq\Delta
E_{\gamma}}\frac{d^{3}k_{3}}{2E_{\gamma}}\biggl{[}\frac{k_{1}}{k_{1}\cdot
k_{3}}-\frac{k_{2}}{k_{2}\cdot k_{3}}\biggl{]}^{2}\end{split}$ (20)

where $\Delta E_{\gamma}$ satisfies $E_{\gamma}\leq\Delta
E_{\gamma}\ll\sqrt{s}$. Integrating the soft photon phase space in the center-
of-mass system gives

$\begin{split}d\sigma_{\text{soft}}=\delta_{\text{soft}}d\sigma_{\text{Born}}\end{split}$
(21)

with

$\begin{split}\delta_{\text{soft}}=&-\frac{\alpha}{\pi}\biggl{[}2\ln\biggl{(}\frac{2\Delta
E_{\gamma}}{m_{\gamma}}\biggr{)}\biggl{(}1+\ln\biggl{(}\frac{m^{2}_{e}}{s}\biggr{)}\biggr{)}\\\
&+\frac{1}{2}\ln^{2}\biggl{(}\frac{m^{2}_{e}}{s}\biggr{)}+\ln\biggl{(}\frac{m^{2}_{e}}{s}\biggr{)}+\frac{\pi^{2}}{3}\biggr{]}.\end{split}$
(22)

Both the soft and hard photon corrections depend on the soft cutoff parameter
$\Delta_{s}$, but the real correction is independent of this parameter. We
note that the hard photon emission must be taken account into to remove this
dependency. Furthermore, adding the virtual and soft corrections removes the
$m_{\gamma}$ dependency.

In addition to the divergences mentioned above, there also appear mass
singularities as a consequence of the collinear photon emission off massless
particles, so-called collinear divergences. However, the smallness of the
electron mass induces the quasi-collinear IR divergences from the photon
radiation off the electron/positron. To deal with this, we apply the phase
space slicing method. According to this method, the hard bremsstrahlung phase
space can be divided into collinear and finite regions:

$d\sigma_{\text{hard}}^{e^{-}e^{+}\rightarrow
ZZ\gamma}(\Delta_{s})=d\sigma_{\text{coll}}(\Delta_{s},\Delta_{c})+d\sigma_{\text{fin}}(\Delta_{s},\Delta_{c})$
(23)

where $\Delta_{c}$ denotes the angular cut-off parameter. In the hard
collinear region $(E_{\gamma}\geq\Delta
E_{\gamma},\cos\theta_{e\gamma}>1-\Delta_{c})$, the integrand is numerically
unstable, whereas in the hard finite region $(E_{\gamma}\geq\Delta
E_{\gamma},\cos\theta_{e\gamma}\leq 1-\Delta_{c})$, it is finite (numerically
stable). Here, $\cos\theta_{e\gamma}$ is the cosines of angle between the
electron/positron and bremsstrahlung photon. In the stable regions the
integration is carried out numerically, while it is semi-analytically
calculated in the unstable regions [58] by using approximation:

$\begin{split}d\sigma_{\text{coll}}=&\sum_{i=1}^{2}\frac{\alpha}{2\pi}Q^{2}_{i}\int^{1-\Delta_{s}}_{0}dx~{}d\sigma_{\text{Born}}(\sqrt{xs})\\\
&\times\biggl{[}\frac{1+x^{2}}{1-x}\ln\biggl{(}\frac{s\Delta_{c}}{2m_{i}^{2}}\biggr{)}-\frac{2x}{1-x}\biggl{]}\end{split}$
(24)

where the approximation $\Delta_{c}\gg 2m_{e}^{2}/s$ has been taken. To avoid
over-counting in the soft energy region, the integration over all possible
factors $x$ is constrained by the soft cut-off parameter $\Delta_{s}$.

####  II.2.3 Classification of full corrections

As a result, the IR and UV finite EW corrections consist of four parts:

$\begin{split}d\sigma_{\text{NLO}}^{e^{+}e^{-}\rightarrow
ZZ}&=d\sigma_{\text{virt}}(m_{\gamma})+d\sigma_{\text{soft}}(m_{\gamma},\Delta_{s})\\\
&+d\sigma_{\text{coll}}(\Delta_{s},\Delta_{c})+d\sigma_{\text{fin}}(\Delta_{s},\Delta_{c}).\end{split}$
(25)

where they are from virtual corrections, soft photon corrections, collinear
corrections and finite hard photon corrections, respectively.

$\Delta_{s}$ | $\sigma_{\text{1-loop}}/\sigma_{\text{Born}}$
---|---
$10^{-1}$ | $0.9933$
$10^{-2}$ | $0.9904$
$10^{-3}$ | $0.9899$
$10^{-4}$ | $0.9899$
$10^{-5}$ | $0.9899$
$10^{-6}$ | $0.9898$
$10^{-7}$ | $0.9897$
$10^{-8}$ | $0.9897$

$\Delta_{c}$ | $\sigma_{\text{1-loop}}/\sigma_{\text{Born}}$
---|---
$10^{-1}$ | $0.9891$
$10^{-2}$ | $0.9899$
$10^{-3}$ | $0.9899$
$10^{-4}$ | $0.9897$
$10^{-5}$ | $0.9894$
$10^{-6}$ | $0.9892$
$10^{-7}$ | $0.9889$

Figure 7:  Phase space slicing method. The virtual, soft and hard photon
radiation corrections as a function of soft cut-off parameter $\Delta_{s}$
(upper plot) and angular cut-off parameter $\Delta_{c}$ (lower plot) at
$\sqrt{s}=250\mbox{~{}GeV}$.

We have numerically checked that our results are independent of $m_{\gamma}$,
soft cut-off parameter $\Delta_{s}$ and angular cut-off parameter
$\Delta_{c}$. The virtual, soft photon and hard photon emission corrections
are shown as a function of both $\Delta_{s}$ and $\Delta_{c}$ at
$\sqrt{s}=250\mbox{~{}GeV}$ in Fig. 7. It is obvious that the virtual plus
soft and the hard photon emission corrections significantly depend on
$\Delta_{s}$ and $\Delta_{c}$, whereas the total corrections do not depend on
them. Moreover, the relative correction
$\sigma_{\text{1-loop}}/\sigma_{\text{Born}}$ is independent of $\Delta_{s}$
and $\Delta_{c}$ over several orders of magnitude. In order to more clearly
show the cut-off independence, we have also given the values by tables next to
Fig. 7. Finally, our results are also stable over nine digits when varying
$m_{\gamma}$ from $10^{-20}$ GeV to the default value of $1$ GeV. Our
numerical results below have been obtained for $\Delta_{s}=10^{-3}$,
$\Delta_{c}=10^{-2}$ and the range of scattering angles of the final particles
$|\cos\theta|<0.99$.

In order to discuss the origin of the large correction, we also consider the
pure QED corrections $\delta_{\text{QED}}$ as individual in addition to the
full EW corrections $\delta_{\text{Total}}$. The QED corrections consist of
virtual-photon exchange, the corresponding counterterms and real-photon
emission. Hence, we can express the QED relative correction as

$\delta_{\text{QED}}=\delta_{\text{virt, QED}}+\delta_{\text{real}}.$ (26)

The QED-like diagrams consist of only $A$ and $f$ fields as virtual lines. In
this study, the one-loop QED contributions come from the sum of the real-
photonic corrections and contributions of diagrams $(4)$ and $(5)$ of Fig. 3,
$(1)$ and $(7)$ of Fig. 4, and $(1)$ and $(9)-(12)$ of Fig. 2.

The remaining corrections (non-QED) can be called as Weak corrections
$\delta_{\text{Weak}}$, which include the massive gauge bosons $Z^{0}$ and
$W^{\pm}$. The genuine weak relative correction can be obtained from

$\delta_{\text{Weak}}=\delta_{\text{Total}}-\delta_{\text{QED}}.$ (27)

In the above definition of $\delta_{\text{Weak}}$, the Weak corrections are
obtained by subtracting the pure QED corrections from the full EW radiative
corrections. Thus, we can write the full $\mathcal{O}(\alpha)$ EW relative
correction as

$\delta_{\text{Total}}=\delta_{\text{QED}}+\delta_{\text{Weak}}.$ (28)

Overall, we can factorize the full EW corrected cross section into the Born
cross section and the relative corrections. Therefore, the one-loop cross
section $\sigma_{\text{1-loop}}$ becomes

$\begin{split}\sigma_{\text{1-loop}}&=\sigma_{\text{Born}}(1+\delta_{\text{Total}})\\\
&=\sigma_{\text{Born}}(1+\delta_{\text{QED}}+\delta_{\text{Weak}}),\end{split}$
(29)

leading to

$\begin{split}\delta_{\text{X}}&=\frac{\sigma^{\text{X}}_{\text{1-loop}}-\sigma_{\text{Born}}}{\sigma_{\text{Born}}}\end{split}$
(30)

where X can be “QED”, “Weak” and “Total”.

##  III Beam Polarization

Polarization effects are important in $e^{+}e^{-}$ colliders and can be used
to provide significant advantages. This effects could benefit to searches for
new physics with small deviation from the SM predictions in two ways. First,
properly chosen combinations of beam polarization can strengthen the signal
and suppress the background in many cases. Second, it is possible to establish
smart observables that contain beam polarization information. With the above
motives, we analyse the effect of beam polarizations on the production rate of
$e^{+}e^{-}\to ZZ$.

The polarizations with a sign for $e^{-}$ and $e^{+}$ beams are given by [59],
respectively

$\begin{split}P_{e^{-}}&=\frac{(n_{e^{-}_{R}}-n_{e^{-}_{L}})}{(n_{e^{-}_{R}}+n_{e^{-}_{L}})},\,\,\
P_{e^{+}}=\frac{(n_{e^{+}_{R}}-n_{e^{+}_{L}})}{(n_{e^{+}_{R}}+n_{e^{+}_{L}})}\end{split}$
(31)

where $n_{e^{\pm}_{L}}$ and $n_{e^{\pm}_{R}}$ are the number of the left- and
right-handed $e^{\pm}$’s ($e^{\pm}_{L}$ and $e^{\pm}_{R}$) in the beam,
respectively. Here, $P_{e^{-}}$/$P_{e^{+}}$ equal to $+1$ (-1) for the $100\%$
right-(left)-handed polarized $e^{-}/e^{+}$-beams.

If a normalization $n_{e^{\pm}_{L}}+n_{e^{\pm}_{R}}=1$ is applied, the
normalized number of $e^{-}_{R}$’s and $e^{-}_{L}$’s can be obtained as

$\begin{split}n_{e^{-}_{R}}=\frac{1+P_{e^{-}}}{2},\,\,\
n_{e^{-}_{L}}=\frac{1-P_{e^{-}}}{2}.\end{split}$ (32)

Consequently, the cross section for any beam polarizations can be defined by
[60, 59, 16]

$\begin{split}\sigma^{P_{e^{+}}P_{e^{-}}}=&\frac{1}{4}\big{[}(1-P_{e^{+}})(1-P_{e^{-}})\sigma_{LL}\\\
&+(1+P_{e^{+}})(1+P_{e^{-}})\sigma_{RR}\\\
&+(1-P_{e^{+}})(1+P_{e^{-}})\sigma_{LR}\\\
&+(1+P_{e^{+}})(1-P_{e^{-}})\sigma_{RL}\big{]},\end{split}$ (33)

where $\sigma_{LL}$, $\sigma_{RR}$, $\sigma_{LR}$, and $\sigma_{RL}$ indicate
the cross sections with completely polarized beams of the four possible cases.
Namely, $\text{RL},\text{LR},\text{RR}$ and LL stand for
$(P_{e^{+}},P_{e^{-}})=(+1,-1),(-1,+1),(+1,+1),(-1,-1)$, respectively. Figure
8 shows these spin configurations, the corresponding fractions (the fourth
column) and the total spin projections onto the $e^{+}e^{-}$ direction (the
last column). This figure is adapted from Ref. [16].

Figure 8: The longitudinal spin configurations in electron-positron
collisions. The thick arrow denotes the direction of motion of the particle
and the double arrow its spin direction.

Now we give an expression of the left-right asymmetry $A_{LR}$, which has
several advantages such as it is independent of detector efficiency
asymmetries and its measurement has negligible systematic error. It is defined
by

$\displaystyle
A_{LR}=\frac{\sigma(P_{e^{+}}=+1,P_{e^{-}}=-1)-\sigma(P_{e^{+}}=-1,P_{e^{-}}=+1)}{\sigma(P_{e^{+}}=+1,P_{e^{-}}=-1)+\sigma(P_{e^{+}}=-1,P_{e^{-}}=+1)}$

or equivalently

$\begin{split}A_{LR}&=\frac{\sigma_{RL}-\sigma_{LR}}{\sigma_{LR}+\sigma_{RL}}\end{split}$
(34)

for a given process. For the process $e^{+}e^{-}\to ZZ$ considered in this
study, it is obtained by

$\begin{split}A_{LR}^{\text{Born}}=\frac{(c_{L}^{e})^{4}-(c_{R}^{e})^{4}}{(c_{L}^{e})^{4}+(c_{R}^{e})^{4}}=0.40979\end{split}$
(35)

at the Born-level.

If $e^{+}e^{-}$ is annihilated into a vector particle, only $\sigma_{LR}$ and
$\sigma_{RL}$ ($J=1$ configurations) have a non-zero contribution. In this
study, the $\sigma_{LL}$ and $\sigma_{RR}$ have the tiny values due to small
electron mass effects, so they can be neglected. In this case, we can also
rewrite the left-right asymmetry as [16]

$\begin{split}A_{LR}=\frac{1}{P_{eff}}A_{LR}^{\text{obs}}=\frac{1}{P_{eff}}\biggl{(}\frac{\sigma_{-+}-\sigma_{+-}}{\sigma_{-+}+\sigma_{+-}}\biggr{)}\end{split}$
(36)

in terms of the effective polarization $P_{eff}$ and the measured left-right
asymmetry $A_{LR}^{\text{obs}}$

$P_{eff}=\frac{(P_{e^{-}}-P_{e^{+}})}{(1-P_{e^{-}}P_{e^{+}})},\,\,A_{LR}^{\text{obs}}=\frac{\sigma_{-+}-\sigma_{+-}}{\sigma_{-+}+\sigma_{+-}},$
(37)

where the corresponding cross-sections are given by

$\begin{split}\sigma_{-+}=&\frac{1}{4}\biggl{[}\bigl{(}1+|P_{e^{+}}||P_{e^{-}}|\bigr{)}\bigl{(}\sigma_{RL}+\sigma_{LR}\bigr{)}\\\
&+\bigl{(}|P_{e^{+}}|+|P_{e^{-}}|\bigr{)}\bigl{(}\sigma_{RL}-\sigma_{LR}\bigr{)}\biggr{]}\\\
\sigma_{+-}=&\frac{1}{4}\biggl{[}\bigl{(}1+|P_{e^{+}}||P_{e^{-}}|\bigr{)}\bigl{(}\sigma_{RL}+\sigma_{LR}\bigr{)}\\\
&-\bigl{(}|P_{e^{+}}|+|P_{e^{-}}|\bigr{)}\bigl{(}\sigma_{RL}-\sigma_{LR}\bigr{)}\biggr{]}.\end{split}$
(38)

Another parameter, the effective luminosity $\mathcal{L}_{eff}$ is given by

$\frac{\mathcal{L}_{eff}}{\mathcal{L}}=\frac{1}{2}(1-P_{e^{-}}P_{e^{+}})$ (39)

which basically reflects the fraction of particles that are interacting.

##  IV Parameter Settings

A set of input parameters must be specified with their corresponding numeric
values, in order to provide consistent higher-order predictions in the SM. We
set the input parameters as follows, see also [61]:

  * •

Mass parameters:

$m_{e}=0.510998928\mbox{~{}MeV}$, | $m_{u}=73.56\mbox{~{}MeV}$
---|---
$m_{\mu}=105.6583715\mbox{~{}MeV}$, |  $m_{d}=73.56\mbox{~{}MeV}$,
$m_{\tau}=1.77682\mbox{~{}GeV}$, |  $m_{s}=95\mbox{~{}MeV}$,
$M_{W}=80.385\mbox{~{}GeV}$, |  $m_{c}=1.275\mbox{~{}GeV}$,
$M_{Z}=91.1876\mbox{~{}GeV}$, |  $m_{b}=4.66\mbox{~{}GeV}$,
$M_{h}=125\mbox{~{}GeV}$, |  $m_{t}=173.21\mbox{~{}GeV}$.

Here, the $u$\- and $d$-quark masses are calculated as a effective parameters
and are specially given for the $M_{Z}$-mass scale via the hadronic
contributions

$\displaystyle\begin{split}\Delta\alpha_{\text{had}}^{(5)}(M_{Z})&=\frac{\alpha}{\pi}\sum_{f=u,c,d,s,b}Q_{f}^{2}\Bigl{(}\ln\frac{M_{Z}^{2}}{m_{f}^{2}}-\frac{5}{3}\Bigr{)}\\\
&\approx 0.027547\,.\end{split}$ (40)

According to [61], $s$-quark mass $m_{s}$ is an estimate of a so-called
”current quark mass” in the $\overline{MS}$ scheme at scale $\mu\approx
2\mbox{~{}GeV}$. $m_{c}\equiv m_{c}(m_{c})$ is “running” mass in the
$\overline{MS}$ scheme and $m_{b}$ is the $\Upsilon(1S)$ bottom quark mass.

  * •

The fine structure constant:

$\alpha(0)=1/137.03599907$ (41)

  * •

The Fermi constant:

$G_{F}=1.1663787(6)\times 10^{-5}\mbox{~{}GeV}^{-2}$ (42)

which is conventionally defined via the muon lifetime.

The renormalization scale $\mu_{R}$ is fixed to center-of-mass energy,
$\sqrt{s}$.

On the other hand, it is important to specify the electromagnetic coupling
$\alpha=e^{2}/(4\pi)$ for the EW $\mathcal{O}(\alpha)$ corrections. For an
obvious choice of $\alpha$, there are two different methods as follows: the
fine-structure constant $\alpha(0)$ in the Thompson limit ($\alpha(0)$ scheme)
and the running electromagnetic coupling $\alpha(Q^{2})$ at any energy scale
$Q$. It is possible to use the value of $\alpha(M^{2}_{Z})\approx 1/129$
($\alpha(M^{2}_{Z})$ scheme), which is calculated by analyzing the
experimental ratio
$R=\sigma(e^{-}e^{+}\rightarrow\text{hadrons})/\sigma(e^{-}e^{+}\rightarrow\mu^{-}\mu^{+})$
[62, 63]. Another choice is a $G_{\mu}$ scheme given by

$\alpha(G_{\mu})=\frac{\sqrt{2}G_{\mu}M^{2}_{W}}{\pi}\biggl{(}1-\frac{M^{2}_{W}}{M^{2}_{Z}}\biggr{)}.$
(43)

An effective value of $\alpha$ in this scheme is obtained as
$\alpha(G_{\mu})\approx 1/132$ from the Fermi constant $G_{\mu}$.

The $G_{\mu}$ scheme provides the possibility to absorb some important
universal corrections associated with the renormalization of the weak mixing
angle into leading order contributions. At next-to-leading order (NLO), the
$\alpha(0)$ and $G_{\mu}$ schemes are related by

$\alpha(G_{\mu})=\alpha(0)\bigl{[}1+\Delta
r^{(1)}\bigr{]}+\mathcal{O}(\alpha^{3})$ (44)

where $\Delta r^{(1)}$ denotes the EW correction to muon decay at NLO [64,
48]. In Sirlin’s relation, the resummation is achieved by the replacement
$(1+\Delta r^{(1)})\rightarrow\frac{1}{(1-\Delta r^{(1)})}$ in Eq. (44). One
then obtains a much closer agreement between the two schemes. This fact also
reveals that the $G_{\mu}$-scheme provides more accurate results above the
$M_{Z}$-mass scale unless one applies the leading log resummation in the
strict NLO $\alpha$-scheme.

Actually, the suitability of the scheme directly relates to the nature of the
considered process. In all cases, a common coupling factor $\alpha^{n}$ should
be used in full gauge-invariant sub-groups, otherwise significant consistency
relations disappear [65].

In this study, we present the results obtained by both $\alpha(0)$ scheme and
$G_{\mu}$ scheme, and discuss the difference due to these choices.

##  V Numerical Results And Discussions

We examine the Z boson pair production in electron-positron collision, by
taking into account a full one-loop EW $\mathcal{O}(\alpha)$ corrections,
including soft and hard QED radiation. We provide the center-of-mass energy
dependence of Born and one-loop cross sections. In order to discuss the origin
of the large correction, we also present the relative QED and Weak corrections
(as defined in Eq. (30)) as a function of the center-of-mass energy. We
include the spin polarization effects of the initial electron and positron
beams in the total cross sections. We discuss how the considered process is
affected by beam polarization. The left-right asymmetry and angular
distributions are also presented.

First of all, for the numerical verification of Born-level and hard photon
bremsstrahlung calculations, we use three different tools
FeynArts$\&$FormCalc, CalcHEP and Whizard. In Table 1, we give the results
obtained by them for Born level process $e^{-}e^{+}\to ZZ$ and hard process
$e^{-}e^{+}\to ZZ\gamma$ at $\sqrt{s}=250$, $500$, and $1000\mbox{~{}GeV}$.

Table 1: The Born and hard photon bremsstrahlung cross sections obtained by FeynArts$\&$FormCalc (FA$\&$FC), Whizard and CalcHEP. $~{}~{}~{}\sqrt{s}$ | $250\mbox{~{}GeV}$ | $500\mbox{~{}GeV}$ | $1000\mbox{~{}GeV}$
---|---|---|---
|  $\sigma^{\text{Born}}(e^{-}e^{+}\to Z^{0}Z^{0})$ [fb]
FA$\&$FC | 991.500(2) | 315.463(9) | 79.592(8)
Whizard | 991.501(5) | 315.464(3) | 79.595(1)
CalcHEP | 991.50 | 315.46 | 79.593
|  $\sigma^{\text{hard}-\gamma}(e^{-}e^{+}\to Z^{0}Z^{0}\gamma)$ [fb]
FA$\&$FC | 144.443(4) | 52.751(9) | 14.852(4)
Whizard | 144.457(1) | 52.739(5) | 14.851(3)
CalcHEP | 144.46 | 52.758 | 14.859

It is found that Born-level results are an excellent agreement. The hard
photon bremsstrahlung results are in good agreement up to five digits.

In Fig. 9, we plot the Born and one-loop cross sections, and the relative
corrections as a function of center-of-mass energy. We also show the proposed
energies of future colliders with the vertical solid lines. The center of mass
energy ranges from 190 GeV to 1.5 TeV in the steps of 10 GeV.

Figure 9: Born and one-loop cross sections as a function of $\sqrt{s}$. The
relative corrections in a percentage is also shown at the bottom panel. The
vertical solid lines show the proposed energies for various future colliders.

Since the colliding energy $\sqrt{s}$ starts near the threshold of $2m_{Z}$,
the Born and one-loop cross sections increase quickly with the opening of the
phase space, reach a maximum value and then fall off rapidly with the
increment of $\sqrt{s}$. This is also the expected behavior. The one-loop
cross section reaches a maximum of 1109.93 fm with
$\delta_{\text{Total}}=-8.45\%$ at $\sqrt{s}=210\mbox{~{}GeV}$ and then
decreases to 40.04 fb with $\delta_{\text{Total}}=+14.17\%$ at $\sqrt{s}=1.5$
TeV. The pure QED corrections make a positive contribution and increase from
$-8.1\%$ to $+31.15\%$ when $\sqrt{s}$ goes from 200 GeV to 1.0 TeV. On the
other hand, the weak corrections make a negative contribution and its relative
correction decreases from $-4.63\%$ to $-17.69\%$. The QED and weak
corrections partially compensate each other, providing relative corrections of
around $-12.73\%$ at the first point and $+13.46\%$ at 1.0 TeV. The EW
radiative corrections significantly increase with the $\sqrt{s}$. This is due
to the presence already at the one-loop order, of large double and single
logarithm terms behaving like $(\alpha/\pi)\ln^{2}s$ and $(\alpha/\pi)\ln s$.

Figure 10: The effect of different individual contributions on cross sections
as a function of $\sqrt{s}$.

These results show that the pure QED and weak corrections are the same order
of magnitude, so that both are equally important. However, the QED correction
make the main contribution to the total EW correction. The relative radiative
correction due to full EW one-loop contributions in the vicinity $\sqrt{s}$
close to the threshold of production ZZ becomes rather large. This effect is
due to the Coulomb singularity in Feynman diagrams, which includes the
instantaneous virtual photon exchange in the loop that has a small spatial
momentum. Overall, the full EW one-loop corrections enhance Born cross
section. While energy-dependent structure of both the QED and weak corrections
is clearly visible, it almost disappears in the total corrections at higher
energies. For proposed colliding energies of the future collider projects:
CEPC (at $\sqrt{s}=240$ GeV), FCCee (at $\sqrt{s}=350$ GeV), CLIC (at
$\sqrt{s}=380$ GeV), and ILC (at $\sqrt{s}=500$ GeV), the unpolarized born-
level cross sections reach 1051.44 fb (with $\delta_{\text{Total}}=-2.20\%$),
588.23 fb (with $\delta_{\text{Total}}=+4.73\%$), 513.82 fb (with
$\delta_{\text{Total}}=+5.70\%$) and 315.46 fb (with
$\delta_{\text{Total}}=+8.66\%$), respectively. The production rate of
$e^{-}e^{+}\to ZZ$ is larger by around one order of magnitude than from the
$\gamma\gamma$-collision mode (see Refs.[66, 67]).

In Fig. 10, we plot the virtual, soft and hard photonic corrections as a
function of $\sqrt{s}$. The virtual corrections decrease from 0.24 to -0.06,
while the soft bremsstrahlung corrections remains nearly constant with a value
of -0.91 with the increment of $\sqrt{s}$. However, hard bremsstrahlung
corrections, i.e., collinear and non-collinear parts, start from their minimum
values and then increase rapidly with the increment of $\sqrt{s}$.
Consequently, the virt+soft and hard bremsstrahlung corrections are partially
canceled, as they are combined into the full EW corrections.

Figure 11: The ratios of cross sections and relative corrections in two
different schemes, $\alpha(0)$ and $G_{\mu}$ schemes, as a function of
$\sqrt{s}$.

Figure 12: Polarized Born and one-loop cross sections as a function of
$\sqrt{s}$.

The soft and hard bremsstrahlung corrections make the negative and positive
dominant contributions to the full EW corrections, respectively.

For $\alpha(0)$ and $G_{\mu}$ schemes, we present the ratios of cross sections
in Fig. 11(a) and the relative corrections in Fig. 11(b) as a function of
$\sqrt{s}$, respectively. In Table 2, we also provide our numerical results
obtained in the $\alpha(0)$ and $G_{\mu}$ schemes for $\sqrt{s}=250$, $500$,
and $1000\mbox{~{}GeV}$. The Born cross section increase by about $7.39\%$ in
$G_{\mu}$ scheme as compared to the $\alpha(0)$ scheme. When $\sqrt{s}$ goes
up from 200 GeV to 1.5 TeV, the one-loop cross section in $G_{\mu}$ scheme is
from $6.58\%$ to $7.69\%$ larger than ones in the $\alpha(0)$ scheme. The QED
relative corrections in $G_{\mu}$ scheme increase by up to about $4\%$ , while
the weak relative corrections decrease by up to about $5\%$, as compared to
the $\alpha(0)$ scheme. However, the difference between total relative
corrections in considered schemes is one percent so small that it can be
considered as a theoretical uncertainty.

Table 2: Born and one-loop cross sections, and relative corrections of $e^{+}e^{-}\to Z^{0}Z^{0}$ in the $\alpha(0)$ scheme and $G_{\mu}$ scheme for various values of $\sqrt{s}$. $\sqrt{s}$ | $250\mbox{~{}GeV}$ | $500\mbox{~{}GeV}$ | $1000\mbox{~{}GeV}$
---|---|---|---
|  $\sigma(e^{+}e^{-}\to Z^{0}Z^{0})$ [fb]
$\sigma_{\text{Born}}^{\alpha(0)}$ | 991.50 | 315.46 | 79.59
$\sigma_{\text{Born}}^{G_{\mu}}$ | 1064.83 | 338.88 | 85.48
$\sigma_{\text{1-loop}}^{\alpha(0)}$ | 981.59 | 342.77 | 90.31
$\sigma_{\text{1-loop}}^{G_{\mu}}$ | 1051.67 | 368.51 | 97.24
|  $\delta(e^{+}e^{-}\to Z^{0}Z^{0})$ [$\%$]
$\delta_{\text{QED}}^{\alpha(0)}$ | +6.35 | +19.50 | +31.15
$\delta_{\text{QED}}^{G_{\mu}}$ | +6.39 | +20.38 | +32.27
$\delta_{\text{Weak}}^{\alpha(0)}$ |  $-$7.35 |  $-$10.85 |  $-$17.69
$\delta_{\text{Weak}}^{G_{\mu}}$ |  $-$7.63 |  $-$11.61 |  $-$18.52
$\delta_{\text{Total}}^{\alpha(0)}$ |  $-$1.00 | +8.65 | +13.46
$\delta_{\text{Total}}^{G_{\mu}}$ |  $-$1.24 | +8.77 | +13.75

Now we examine the initial beam polarization dependence of the Born and the
one-loop cross sections ($\sigma^{P_{e^{+}}P_{e^{-}}}_{\text{Born}}$ and
$\sigma^{P_{e^{+}}P_{e^{-}}}_{\text{1-loop}}$ ) on $\sqrt{s}$. Also, we
present the total relative corrections in order to see effect of polarization
configurations on the EW corrections. Figure 12(a) shows unpolarized and
completely polarized initial beams cases. Here we use the following notation:
$\sigma^{\text{RL}}$ denotes the total cross section with the $100\%$ right-
handed polarized positron ($P_{e^{+}}=+1$), and the $100\%$ left-handed
polarized electron ($P_{e^{-}}=-1$) beams, $({e^{+}_{R}},{e^{-}_{L}})$. Others
can be defined analogously. We note that $\sigma^{\text{LL}}$ and
$\sigma^{\text{RR}}$ are very small ($\mathcal{O}(10^{-10})$ pb) as expected,
and are thus not included here. All curves for polarized and unpolarized cases
increase quickly with the opening of the phase space, reach a maximum value
and then fall off rapidly with the increment of $\sqrt{s}$. Their maximum
values are reached around at $\sqrt{s}\sim 210$ GeV. The one-loop polarized
cross sections $\sigma^{\text{LR}}_{\text{1-loop}}$ and
$\sigma^{\text{RL}}_{\text{1-loop}}$ can be enhanced by about a factor of 1.5
and 2.5, respectively, compared with the unpolarized case. For instance, at
$\sqrt{s}=$ 250 GeV, $\sigma^{\text{RL}}_{\text{Born}}$ reaches a value of
2.80 pb, yielding a total relative correction of about $-10.71\%$. At
$\sqrt{s}=$ 250 GeV, $\sigma^{\text{LR}}_{\text{Born}}$ reaches a value of
1.17 pb, yielding a total relative correction of $+22.20\%$. While
$\sigma^{\text{RL}}_{\text{Born}}$ is around 2 times larger than
$\sigma^{\text{RL}}_{\text{Born}}$, the relative corrections of the former are
smaller than that of the latter. Namely, $({e^{+}_{L}}{e^{-}_{R}})$
polarization case has larger EW radiative corrections than other cases.

Figure 13: (color online). Beam polarization effects on the Born cross section
and on the total relative correction for $\sqrt{s}=250\mbox{~{}GeV}$. The
color heat maps correspond to the values of born cross section (left plot) and
the total EW relative correction (right plot). The asterisk denotes the
unpolarized point $(P_{e^{+}},P_{e^{-}})=(0,0)$.

Figure 14: (color online). Same as in Fig. 13 but for
$\sqrt{s}=500\mbox{~{}GeV}$.

Figure 12(b) shows the results for various polarization degrees of the initial
beams: $(P_{e^{+}},P_{e^{-}})=(0.0,-0.8),(+0.3,-0.8),(+0.6,-0.8)$, proposed by
the future colliders. The cross sections at both Born and one-loop level are
sorted according to various polarization degrees of initial beams as follows:
$\sigma^{\text{RL}}>\sigma^{\text{+0.6,-0.8}}>\sigma^{\text{+0.3,-0.8}}>\sigma^{\text{0.0,-0.8}}>\sigma^{\text{LR}}>\sigma^{\text{UU}}$.
There appears a similar energy dependence behavior in polarized cases as in
the unpolarized case. The one-loop cross section with the polarization degrees
of $(P_{e^{+}},P_{e^{-}})=(+0.6,-0.8)$, i.e.,
$\sigma^{(+0.6,-0.8)}_{\text{1-loop}}$, has a maximum of 2.05 pb, providing a
total relative correction of $-14.59\%$. When $\sqrt{s}$ goes from 240 GeV to
1.0 TeV, the relative corrections vary from about $-11\%$ to $+1.69\%,+0.75\%$
and $+0.14\%$ for $(P_{e^{+}},P_{e^{-}})=(0.0,-0.8),(+0.3,-0.8),(+0.6,-0.8)$,
respectively. We obtain the following results for proposed polarization cases
by the future collider projects: $\sigma^{\text{+0.3,-0.8}}_{\text{1-loop}}=$
1.51 pb with $\delta^{\text{+0.3,-0.8}}_{\text{Total}}=-9.90\%$ at
$\sqrt{s}=250$ GeV (ILC) and $\sigma^{\text{0,-0.8}}_{\text{1-loop}}=$ 0.66 pb
with $\delta^{\text{0,-0.8}}_{\text{Total}}=-3.94\%$ at $\sqrt{s}=380$ GeV
(CLIC), respectively.

In Figs. 13 and 14, the Born cross section and the total relative correction
are also presented in the plane of polarization degrees of the incoming beams
$(P_{e^{+}},P_{e^{-}})$ for $\sqrt{s}=250\mbox{~{}GeV}$ and
$\sqrt{s}=500\mbox{~{}GeV}$, respectively. The $P_{e^{+}}$ and $P_{e^{-}}$
range from $-1$ to $+1$. We also show some values by the contour lines. The
beam polarization dependence of the relative correction is calculated from
$\delta^{P_{e^{+}}P_{e^{-}}}_{\text{Total}}=\bigl{(}\sigma^{P_{e^{+}}P_{e^{-}}}_{\text{1-loop}}/\sigma^{P_{e^{+}}P_{e^{-}}}_{\text{Born}}-1\bigr{)}$
by using Eq. (33). The Born cross section reaches its larger values at the
left top corner ( $0<P_{e^{+}}\leq+1$ and $-1\leq P_{e^{-}}<0$), whereas it
has smaller values at the right top and left bottom corners. As excepted, it
has a maximum value at point $(P_{e^{+}},P_{e^{-}})=(+1,-1)$, namely, $100\%$
left-handed polarized electron and $100\%$ right-handed polarized positron.
The total relative correction reaches positive and larger values at the right
bottom corner ($-1\leq P_{e^{+}}<0$ and $0<P_{e^{-}}\leq+1$). As the values of
polarization degrees approach to $P_{e^{-}}\rightarrow+1$ and
$P_{e^{+}}\rightarrow-1$, the relative correction increases. Particularly,
significant positive corrections are observed in the region below the 0.10
contour line. Together with the larger cross sections, smaller systematic
errors can be expected for the cross-section measurement with the polarized
beams than in the unpolarized case.

Table 3: Born$\&$one-loop cross sections, relative correction and left-right asymmetry for various polarization configurations at $\sqrt{s}=250$. Here, we define the scaling factor as $\text{R}^{P_{e^{+}}P_{e^{-}}}_{X}=\sigma^{P_{e^{+}}P_{e^{-}}}_{\text{X}}/\sigma^{0,0}_{\text{X}}$, where $X$ stands for Born or one-loop level. $(P_{e^{+}},P_{e^{-}})$ | $P_{eff}$ | $L_{eff}/L$ | $A_{LR}^{\text{obs}}(\text{Born})$ | $A_{LR}^{\text{obs}}(\text{1-loop})$ |  $\sigma_{\text{Born}}$ [fb] |  $\sigma_{\text{1-loop}}$ [fb] | $\delta_{\text{Total}}[\%]$ | $\text{R}^{P_{e^{+}}P_{e^{-}}}_{\text{Born}}$ | $\text{R}^{P_{e^{+}}P_{e^{-}}}_{\text{1-loop}}$
---|---|---|---|---|---|---|---|---|---
$(~{}0,~{}0)$ | $0$ | 0.50 | 0 | 0 | 991.50 | 981.59 | $-0.99\%$ | 1.00 | 1.00
$(+1.0,-1.0)$ | $-1.00$ | 1.00 | 0.4098 | 0.2715 | 2795.61 | 2496.35 | $-10.71\%$ | 2.82 | 2.54
$(-1.0,+1.0)$ | $+1.00$ | 1.00 | 0.4098 | 0.2715 | 1170.39 | 1430.25 | $+22.20\%$ | 1.18 | 1.46
$(~{}~{}0.0,-1.0)$ | $-1.00$ | 0.50 | 0.4098 | 0.2715 | 1397.81 | 1248.17 | $-10.70\%$ | 1.41 | 1.27
$(~{}~{}0.0,-0.8)$ | $-0.80$ | 0.50 | 0.3278 | 0.2172 | 1316.54 | 1194.87 | $-9.24\%$ | 1.33 | 1.22
$(~{}~{}0.0,+0.8)$ | $+0.80$ | 0.50 | 0.3278 | 0.2172 | 666.46 | 768.43 | $+15.30\%$ | 0.67 | 0.78
$(+0.3,-0.8)$ | $-0.89$ | 0.62 | 0.3635 | 0.2409 | 1676.40 | 1510.42 | $-9.90\%$ | 1.69 | 1.54
$(-0.3,+0.8)$ | $+0.89$ | 0.62 | 0.3635 | 0.2409 | 782.52 | 924.07 | $+18.09\%$ | 0.79 | 0.94
$(+0.6,-0.8)$ | $-0.95$ | 0.74 | 0.3876 | 0.2568 | 2036.25 | 1825.98 | $-10.33\%$ | 2.05 | 1.86
$(-0.6,+0.8)$ | $+0.95$ | 0.74 | 0.3876 | 0.2568 | 898.59 | 1079.71 | $+20.16\%$ | 0.91 | 1.10

A convenient observable is the left–right asymmetry when applying polarized
beams. This asymmetry is especially important for high precision measurements.
Therefore, we present the left-right asymmetry of the cross sections, which is
defined in Eq. (34), as a function of center-of-mass energy in Fig. 15. For
Born-level, the left-right asymmetry remains a stable of 0.40979 with the
increment of $\sqrt{s}$. This value is exactly the same as the result
calculated theoretically in Eq. (35). On the other hand, at one-loop level,
the left-right asymmetry decrease from 0.2616 to 0.2065, as the center-of-mass
energy goes up from 190 to 1500 GeV. Consequently, the radiative corrections
can conveniently be described by the one-loop induced deviation from its value
at Born-level.

Figure 15: The left-right asymmetry as a function of $\sqrt{s}$.

In Figs. 16(a)-(b), we plot the left-right asymmetry of the cross sections
$A_{LR}^{\text{obs}}$, which is defined in Eq. (36), as a function of
polarization degrees of $(P_{e^{+}},P_{e^{-}})$. Some values are also shown
with contour lines. We note that the results are symmetric with respect to the
exchange of $P_{e^{+}}\leftrightarrow P_{e^{-}}$, since they only appear by
their absolute values in the Eq. (38). It is obvious that the left-right
asymmetry $A_{LR}^{\text{obs}}$ increases as both $|P_{e^{+}}|$ and
$|P_{e^{-}}|$ run from 0 to 1. The maximum values reaches at the point
$(|P_{e^{+}}|,|P_{e^{-}}|)=(1,1)$. These values are 0.40979 and 0.2657 for
Born and one-loop level, respectively. In these figures, we also mark some
special points by asterisk.

Figure 16: (color online). The left-right asymmetry in the plane of
polarization degrees of the incoming beams $(P_{e^{+}},P_{e^{-}})$ for
$\sqrt{s}=500$.

Table 3 presents the numerical values of the polarized cross sections, the
relative corrections and the left-right asymmetry $A_{LR}^{\text{obs}}$ for
various polarization degrees of the initial beams.

Figure 17: Born and one-loop angular distributions for (a)
$\sqrt{s}=250\mbox{~{}GeV}$ and (b) $\sqrt{s}=500\mbox{~{}GeV}$.

The values of effective polarization and effective luminosity, which are
calculated by Eqs. (37) and (39), respectively, are also given. It is clear
that the production rate can be enhanced if $P_{e^{-}}$ and $P_{e^{+}}$ have
negative and positive signs, respectively. In particular, the longitudinally
$(P_{e^{+}},P_{e^{-}})=(+1,-1)$ polarizations of the initial beams
significantly improve the cross-section.

Figures 17(a) and (b) show the differential cross sections at Born and one-
loop level as a function of the angle between the coming electron and the
outing Z boson for $\sqrt{s}=250$ and $500\mbox{~{}GeV}$, respectively. The
angular dependence of the QED, weak and total EW relative corrections is also
shown in the same figures. For Born and one-loop cross sections, the angular
distribution peaks significantly in the (symmetrically) forward and backward
directions. The relative corrections modify somewhat the Born angular
distribution because their influence is larger in the central region. The
total corrections reach their maximum values, when $\cos{\theta}$ goes to the
extreme points -1 or +1 for $\sqrt{s}=250$ and zero point for $\sqrt{s}=500$,
respectively. Namely, the Z bosons are dominantly produced in the forward and
backward directions for $\sqrt{s}=250$, whereas in the central region for
$\sqrt{s}=500$. Therefore, it will be more likely to observe them in these
collision regions. For all values of $\cos{\theta}$, the QED corrections make
a positive contributions, while the weak corrections supply a negative
contributions. The same behavior is observed as in the total cross sections.
The QED and weak contributions are partially offset by each other into EW
contributions. While both corrections show a significant dependence on
$\cos{\theta}$ separately, this decreases for full EW corrections. When
$\cos{\theta}$ runs from 0 to +0.95 or -0.95, at $\sqrt{s}=250\mbox{~{}GeV}$,
the total relative correction $\delta_{\text{Total}}$ increases from $-3.71\%$
to $4.30\%$, while decreases $12.2\%$ to $7.7\%$ at
$\sqrt{s}=500\mbox{~{}GeV}$.

##  VI Summary and Conclusions

For further sensitivity testing of the SM, as well as to look for clues on the
BSM, high precision calculations should be performed. At least, a full set of
one-loop corrections in the production channels must be included to ensure
adequate accuracy. However, the extent to which higher-order computations
beyond one-loop order will be needed depends largely on the expected
experimental accuracy.

In this work, by considering a full set of one-loop EW corrections, we have
investigated the Z-boson pair production at electron-positron collisions. The
UV divergences have been organized by dimensional regularization on the OS
scheme. Furthermore, the IR divergences have been removed by the inclusion of
the soft and hard bremsstrahlung corrections. Also, the collinear divergences
have been corrected by the phase space slicing method. We have verified the
stability of our results against the angular and soft cutoff parameters
$\Delta_{c}$, $\Delta_{s}$, as well as the IR regulator $m_{\gamma}$. We have
also compared the results from Whizard and CalcHEP for the Born and hard
photon bremsstrahlung cross sections with the results of the packages used in
this work, and obtained a very good agreement (up to six digits).

We have carried out the numerical evaluation for the $\alpha(0)$ and $G_{\mu}$
schemes. The differences between these schemes for the QED and weak relative
corrections are about $4\%$ and $5\%$, respectively, whereas this is $1\%$ for
the total relative correction. This can be considered as a theoretical
uncertainty.

Our results show that the Born cross section is commonly increased by the one-
loop EW radiative corrections and the total relative correction is typically
up to about ten percent. The pure QED and weak corrections are the same order
of magnitude, so that both are important for precision calculation. The QED
and weak contributions are partially offset by each other into EW
contributions. However, the QED correction makes the main contribution to the
total EW correction.

Moreover, we have investigated the spin polarization effects of the initial
beams on the total cross sections. We have observed that the radiative
corrections have a large polarization dependence. As a result, an improvement
has been observed by a factor of 2.5 with the $100\%$ right-handed polarized
positron and the $100\%$ left-handed polarized electron beams, compared with
the unpolarized case. The left-right asymmetry and angular distributions have
also presented. The Born and one-loop angular distributions are symmetric and
strongly peaked in the forward and backward directions. For all angles, the
QED corrections make a positive contributions, while the weak corrections
supply a negative contributions. The relative corrections slightly change the
Born angular distribution as their effect is usually larger in the central
region.

In summary, an analysis of the one-loop EW radiative corrections to Z boson-
pair production in electron-positron collisions have been carried out in the
framework of the SM. The initial beam effects on the cross sections have been
discussed in detail. It has been clearly shown that one-loop EW radiative
corrections significantly alter the lowest-order results and should therefore
be fully accounted for a realistic description of experiments at future
collider energies. Our results provide precise predictions for Z-boson pair
production, which can be tested as experiments achieve higher accuracy.

###### Acknowledgements.

The work of M. Demirci was supported by the Scientific and Technological
Research Council of Turkey (TUBITAK) in the framework of 2219-International
Postdoctoral Research Fellowship Program. This work was also supported in part
by the National Science Foundation Grant No. PHY-2108339.

## References

  * [1] S. L. Glashow, Nucl. Phys. 22, 579-588 (1961).
  * [2] S. Weinberg, Phys. Rev. Lett. 19, 1264-1266 (1967).
  * [3] A. Salam, Conf. Proc. C 680519, 367-377 (1968).
  * [4] G. Aad et al. (ATLAS Collaboration), Phys. Lett. B 716, 1 (2012).
  * [5] S. Chatrchyan et al. (CMS Collaboration), Phys. Lett. B 716, 30 (2012).
  * [6] G. Aad et al. (ATLAS and CMS Collaborations), J. High Energ. Phys. 2016, 45 (2016).
  * [7] P. Bambade et al., The International Linear Collider: A Global Project, Report No. DESY 19-037,  arXiv:1903.01629 [hep-ex].
  * [8] H. Baer et al., Physics Chapter of the ILC Detailed Baseline Design Report, Report No. ILC-REPORT-2013-040, arXiv:1306.6352 [hep-ph].
  * [9] H. Aihara et al., The International Linear Collider A Global Project, arXiv:1901.09829 [hep-ex].
  * [10] E. Accomando et al. (CLIC Physics Working Group), Physics at the CLIC Multi-TeV Linear Collider, Report No. CERN-2004-005, arXiv:hep-ph/0412251.
  * [11] The CEPC Study Group, CEPC Conceptual Design Report: Volume 2 - Physics $\&$ Detector, Report No. IHEP-CEPC-DR-2018-02, 2018, arXiv:1811.10545 [hep-ex].
  * [12] M. Benedikt et al., Future Circular Collider - European Strategy Update Documents: The FCC integrated programme (FCC-int), Report No. CERN-ACC-2019-0007, https://cds.cern.ch/record/2653673, Jan, 2019.
  * [13] A. Abada et al. Eur. Phys. J. Spec. Top. 228, 261 (2019).
  * [14] A. A. Pankov, A. V. Tsytrinov, and N. Paver, Phys. Rev. D 73,115005 (2006).
  * [15] P. Osland, A. A. Pankov, and A. V. Tsytrinov, Eur. Phys. J. C 67, 191–204 (2010).
  * [16] G. Moortgat-Pick et al. Phys. Rep. 460, 131-243 (2008).
  * [17] R. W. Brown and K. O. Mikaelian, Phys. Rev. D 19, 922 (1979).
  * [18] K. J. F. Gaemers and G. J. Gounaris, Z. Phys. C - Particles and Fields  1, 259–268 (1979).
  * [19] A. Denner and T. Sack, Nucl. Phys. B 306, 221 (1988).
  * [20] G. J. Gounaris, J. Layssac, and F. M. Renard, Phys. Rev. D 67, 013012 (2003).
  * [21] G. Abbiendi et al. (OPAL Collaboration), Phys. Lett. B 476, 256 (2000).
  * [22] LEP, ALEPH, DELPHI, L3, OPAL, LEP Electroweak Working Group, SLD Electroweak Group, SLD Heavy Flavor Group Collaboration, A Combination of preliminary electroweak measurements and constraints on the standard model, arXiv:hep-ex/0312023 [hep-ex].
  * [23] P. Achard et al. (L3 Collaboration), Phys. Lett. B 572, 133 (2003).
  * [24] S. Jadach, W. Płaczek, B.F.L. Ward, Phys. Rev. D 56, 6939 (1997).
  * [25] G. J. Gounaris, J. Layssac, and F. M. Renard, Phys. Rev. D 61, 073013 (2000).
  * [26] R. Rahaman and R. K. Singh, Eur. Phys. J. C 77, 521 (2017).
  * [27] S. C. İnan and A. V. Kisselev, J. High Energ. Phys. 2021, 121 (2021).
  * [28] S. Spor, E. Gurkanli, M. Köksal, Nucl. Phys. B 979, 115785 (2022).
  * [29] J. Alcaraz, M. A. Falagán, and E. Sánchez, Phys. Rev. D 61, 075006 (2000).
  * [30] G. J. Gounaris, J. Layssac, and F. M. Renard, Phys. Rev. D 62, 073013 (2000).
  * [31] G. Aad et al. (ATLAS Collaboration), J. High Energ. Phys. 2013,128 (2013).
  * [32] M. Demirci and A. I. Ahmadov, Phys. Rev. D 94, 075025 (2016).
  * [33] M. Demirci, Phys. Rev. D 100, 075006 (2019).
  * [34] M. Demirci, Nucl. Phys. B 961, 115235 (2020).
  * [35] M. Demirci and M. F. Mustamin, Phys. Rev. D 103, 113004 (2021).
  * [36] J.J. Su, W.G. Ma, R. Y. Zhang, S.M. Wang, and L. Guo, Phys. Rev. D 78, 016007 (2008).
  * [37] F. Boudjema, Le D. Ninh, S. Hao, and M. M. Weber, Phys. Rev. D 81, 073007 (2010).
  * [38] S. Heinemeyer and C. Schappacher, Eur. Phys. J. C 76, 220 (2016).
  * [39] S. Heinemeyer and C. Schappacher, Eur. Phys. J. C 78, 536 (2018).
  * [40] L. He-Yi, Z. Ren-You , Ma Wen-Gan, J. Yi and L. Xiao-Zhou, Chinese Phys. C 46, 043105 (2022).
  * [41] J. Küblbeck, M. Böhm, and A. Denner, Comput. Phys. Commun. 60, 165 (1990).
  * [42] T. Hahn, Comput. Phys. Commun. 140, 418 (2001).
  * [43] T. Hahn and M. Perez-Victoria, Comput. Phys. Commun. 118, 153 (1999).
  * [44] T. Hahn, Comput. Phys. Commun. 168, 78-95 (2005).
  * [45] W. Kilian, T. Ohl, and J. Reuter, Eur. Phys. J. C 71, 1742 (2011).
  * [46] M. Moretti, T. Ohl, and J. Reuter, O’Mega: An Optimizing matrix element generator, LC-TOOL-2001-040-rev, arXiv:hep-ph/0102195-rev.
  * [47] A. Belyaev, N. D. Christensen, and A. Pukhov, Comput. Phys. Commun. 184, 1729 (2013).
  * [48] A. Denner, Fortschr. Phys. 41, 307 (1993).
  * [49] G. ’t Hooft and M. Veltman, Nucl. Phys. B 44, 189 (1972).
  * [50] T. Kinoshita, J. Math. Phys. 3, 650 (1962).
  * [51] T. D. Lee and M. Nauenberg, Phys. Rev. 133, 1549 (1964).
  * [52] J. Schwinger, Phys. Rev. 76, 790-817 (1949).
  * [53] K. Fabricius, G. Kramer, G. Schierholz and I. Schmitt, Zeit. Phys. C 11, 315 (1982).
  * [54] G. Kramer and B. Lampe, Fortschr. Phys. 37, 161 (1989).
  * [55] H. Baer, J. Ohnemus and J. F. Owens, Phys. Rev. D 40, 2844 (1989).
  * [56] B. W. Harris and J. F. Owens, Phys. Rev. D 65, 094032 (2002).
  * [57] G. ’t Hooft and M. Veltman, Nucl. Phys. B 153, 365 (1979).
  * [58] U. Baur, S. Keller, and D. Wackeroth, Phys. Rev. D 59, 013002 (1998).
  * [59] T. Omori, A Polarized positron beam for linear colliders, 1st ACFA Workshop on Physics Detector at the Linear Collider, KEK-PREPRINT-98-237.
  * [60] K. I. Hikasa, Phys. Rev. D 33, 3203 (1986).
  * [61] P. A. Zyla et al. (Particle Data Group), Prog. Theor. Exp. Phys. 2020, 083C01 (2020).
  * [62] H. Burkhardt and B. Pietrzyk, Phys. Lett. B 356, 398 (1995).
  * [63] S. Eidelman and F. Jegerlehner, Z. Phys. C 67, 585 (1995).
  * [64] A. Sirlin, Phys. Rev. D 22, 971 (1980).
  * [65] A. Denner and S. Dittmaier, Phys. Rep. 864, 1-163 (2020).
  * [66] G. J. Gounaris, J. Layssac, P. I. Porfyriadis, and F. M. Renard, Eur. Phys. J. C 13, 79–97 (2000).
  * [67] D. Bardin, S. Bondarenko, P. Christova, L. Kalinovskaya, W. von Schlippe, and E. Uglov, Phys. Part. Nuclei Lett. 14, 811–816 (2017).


$\lambda_{i}(\boldsymbol{K}_{0})\to\lambda_{i}:=\frac{(1+\theta_{i}^{2})(\gamma_{0}+\theta_{i}^{2})}{\theta_{i}^{2}},\qquad|\boldsymbol{b}_{i}^{\top}\boldsymbol{v}_{i}|^{2}\to
1-\frac{\gamma_{0}(1+\theta_{i}^{2})}{\theta_{i}^{2}(\theta_{i}^{2}+\gamma_{0})}$

where $\boldsymbol{v}_{i}$ is the unit eigenvector of the input Gram matrix
$\boldsymbol{K}_{0}=\boldsymbol{X}^{\top}\boldsymbol{X}$. Thus claims (a) and
(b) of Assumption 2 hold with $r=|\\{i:\theta_{i}>\gamma_{0}^{1/4}\\}|$,
$\mu_{0}=\rho_{\gamma_{0}}^{\textnormal{MP}}$ being the standard Marcenko-
Pastur law, and
$\lambda_{i}=(1+\theta_{i}^{2})(\gamma_{0}+\theta_{i}^{2})/\theta_{i}^{2}$
being the above values.

We note that $\boldsymbol{X}$ is $n^{-1/2+\varepsilon}$-orthonormal for any
$\varepsilon>0$ almost surely for all large $n$, by the given condition
$\max_{i=1}^{r}\|\boldsymbol{b}_{i}\|_{\infty}<n^{-1/2+\varepsilon}$ and the
bounds, for any $\alpha,\beta\in[n]$,

$\displaystyle\|\boldsymbol{x}_{\alpha}\|_{2}$
$\displaystyle=\mathopen{}\mathclose{{}\left\|\boldsymbol{z}_{\alpha}}\right\|_{2}+\sum_{i=1}^{r}O_{\prec}\mathopen{}\mathclose{{}\left(\|\boldsymbol{a}_{i}\|_{2}|\theta_{i}||b_{i,\alpha}|}\right)=\mathopen{}\mathclose{{}\left\|\boldsymbol{z}_{\alpha}}\right\|_{2}+O_{\prec}\mathopen{}\mathclose{{}\left(n^{-1/2+\varepsilon}}\right)=1+O_{\prec}\mathopen{}\mathclose{{}\left(n^{-1/2+\varepsilon}}\right),$
$\displaystyle\boldsymbol{x}_{\alpha}^{\top}\boldsymbol{x}_{\beta}$
$\displaystyle=\boldsymbol{z}_{\alpha}^{\top}\boldsymbol{z}_{\beta}+\sum_{i=1}^{r}O_{\prec}\mathopen{}\mathclose{{}\left(|\theta_{i}|\Big{(}|\boldsymbol{a}_{i}^{\top}\boldsymbol{z}_{\alpha}||b_{i,\alpha}|+|\boldsymbol{a}_{i}^{\top}\boldsymbol{z}_{\beta}||b_{i,\beta}|\Big{)}+\theta_{i}^{2}\|\boldsymbol{a}_{i}\|_{2}^{2}|b_{i,\alpha}b_{i,\beta}|}\right)=O_{\prec}\mathopen{}\mathclose{{}\left(n^{-1/2+\varepsilon}}\right).$

Hence Theorem 3 applies, showing that $\boldsymbol{K}_{\ell}$ has an outlier
eigenvalue corresponding to each input signal $\theta_{i}$ if and only if
$\theta_{i}>\gamma_{0}^{1/4}$ and $i\in\mathcal{I}_{\ell}$. The statement
(2.10) follows from Theorem 3(c) applied with
$\boldsymbol{v}=\boldsymbol{b}_{i}$. ∎

##  Appendix F Proofs for spiked eigenstructure of the trained CK

In this section, we prove Theorem 5. The proof is an application of Theorem 12
as in the one-layer setting of the preceding section, but now reversing the
roles of $\boldsymbol{X}$ and $\boldsymbol{W}$. We abbreviate

$\boldsymbol{W}=\boldsymbol{W}_{\text{trained}},\qquad\boldsymbol{Y}=\frac{1}{\sqrt{N}}\sigma(\tilde{\boldsymbol{X}}\boldsymbol{W}),\qquad\boldsymbol{K}=\boldsymbol{Y}\boldsymbol{Y}^{\top}$

where $\boldsymbol{K}$ is the CK matrix of interest. In contrast to the
preceding section, the theorem requires characterizing the _left_ spike
singular vector of $\boldsymbol{Y}$, and we will do so using Theorem 12(c).

We first recall the following approximation of $\boldsymbol{W}$ from [BES+22].

######  Proposition 35.

Under Assumption 4, set $\eta=\sum_{t=1}^{T}\eta_{t}$, and let
$\theta_{1},\theta_{2}$ be as defined in (2.14). Then

$\displaystyle\|\boldsymbol{W}-\widetilde{}\boldsymbol{W}\|\prec
N^{-1/2}\quad\text{ where
}\quad\widetilde{}\boldsymbol{W}=\boldsymbol{W}_{0}+\frac{b_{\sigma}\eta}{n}\boldsymbol{X}^{\top}\boldsymbol{y}\boldsymbol{a}^{\top}.$
(F.1)

The largest singular value $s_{\max}(\boldsymbol{W})$ falls outside the limit
of its empirical singular value distribution if and only if
$\theta_{1}>\gamma_{0}^{1/4}$, in which case $s_{\max}(\boldsymbol{W})$ and
its unit-norm left singular vector $\boldsymbol{u}(\boldsymbol{W})$ satisfy

$s_{\max}(\boldsymbol{W})\to
s_{1}:=\sqrt{\frac{(1+\theta_{1}^{2})(\gamma_{0}+\theta_{1}^{2})}{\theta_{1}^{2}}},\quad|\boldsymbol{u}(\boldsymbol{W})^{\top}\boldsymbol{\beta}_{*}|^{2}\to\frac{\theta_{2}^{2}}{\theta_{1}^{2}}\mathopen{}\mathclose{{}\left(1-\frac{\gamma_{0}+\theta_{1}^{2}}{\theta_{1}^{2}(\theta_{1}^{2}+1)}}\right)\quad\text{a.s.}$
(F.2)

Proof. Notice that each gradient update matrix $\boldsymbol{G}_{t}$ of (2.12)
takes the form

$\displaystyle\boldsymbol{G}_{t}=\frac{1}{n}\boldsymbol{X}^{\top}\mathopen{}\mathclose{{}\left[\mathopen{}\mathclose{{}\left(\frac{1}{\sqrt{N}}\mathopen{}\mathclose{{}\left(\boldsymbol{y}-\frac{1}{\sqrt{N}}\sigma(\boldsymbol{X}\boldsymbol{W}_{\\!t})\boldsymbol{a}}\right)\boldsymbol{a}^{\top}}\right)\odot\sigma^{\prime}(\boldsymbol{X}\boldsymbol{W}_{\\!t})}\right],$

From the proof of [BES+22, Lemma 16], for each $t=1,\ldots,T$, this matrix
$\boldsymbol{G}_{t}$ satisfies the same rank-one approximation

$\Big{\|}\sqrt{N}\,\boldsymbol{G}_{t}-\frac{b_{\sigma}}{n}\boldsymbol{X}^{\top}\boldsymbol{y}\boldsymbol{a}^{\top}\Big{\|}\prec
N^{-1/2}.$

This implies (F.1) in light of (2.12), and the statements of (F.2) then follow
from [BES+22, Theorem 3]. ∎

We denote the columns of
$\boldsymbol{W}\equiv\boldsymbol{W}_{\text{trained}}\in\mathbb{R}^{d\times N}$
and of the initialization $\boldsymbol{W}_{0}\in\mathbb{R}^{d\times N}$ by

$\boldsymbol{w}_{i}\in\mathbb{R}^{d},\;\boldsymbol{w}_{i,0}\in\mathbb{R}^{d}\text{
for }i\in[N]$

respectively. Fixing a large constant $B>0$ and small constant
$\varepsilon>0$, define the event

$\mathcal{E}(\boldsymbol{W})=\Big{\\{}\|\boldsymbol{W}\|<B,\,|\boldsymbol{w}_{i}^{\top}\boldsymbol{w}_{j}|<n^{-1/2+\varepsilon}\text{
and }|\|\boldsymbol{w}_{i}\|_{2}-1|<n^{-1/2+\varepsilon}\text{ for all }i\neq
j\in[N]\Big{\\}}.$

######  Lemma 36.

Under Assumption 4, for some sufficiently large constant $B>0$ and any fixed
$\varepsilon>0$, $\mathcal{E}(\boldsymbol{W})$ holds almost surely for all
large $n$ and any fixed $T\in\mathbb{N}$.

Proof. By the assumption
$[\boldsymbol{W}_{0}]_{ij}\overset{iid}{\sim}\mathcal{N}(0,1/d)$, it is
immediate to check that $\mathcal{E}(\boldsymbol{W}_{0})$ holds almost surely
for all large $n$. To show that $\mathcal{E}(\boldsymbol{W})$ holds, we apply
the approximation (F.1). Here, under Assumption 4, we have by standard tail
bounds for Gaussian vectors and matrices that

$\mathbf{1}\\{\|\boldsymbol{X}^{\top}\boldsymbol{y}\|>Cn\\}\leq\mathbf{1}\\{\|\boldsymbol{X}\|\cdot(\lambda_{\sigma_{*}}\|\boldsymbol{X}\boldsymbol{\beta}_{*}\|+\|\boldsymbol{\varepsilon}\|_{2})>Cn\\}\prec
0,\qquad\mathbf{1}\\{\|\boldsymbol{a}\|_{2}>C\\}\prec 0$

for a sufficiently large constant $C>0$, and also
$\|\boldsymbol{a}\|_{\infty}\prec N^{-1/2}$. Then this implies

$\displaystyle\max_{i=1}^{N}\|\boldsymbol{w}_{i}-\boldsymbol{w}_{i,0}\|_{2}$
$\displaystyle\leq\max_{i=1}^{N}\|\widetilde{}\boldsymbol{w}_{i}-\boldsymbol{w}_{i,0}\|_{2}+\|\boldsymbol{W}-\widetilde{}\boldsymbol{W}\|\prec
N^{-1/2}$ (F.3)

and $\mathbf{1}\\{\|\boldsymbol{W}-\boldsymbol{W}_{0}\|>C^{\prime}\\}\prec 0$
for a constant $C^{\prime}>0$. Then $\mathcal{E}(\boldsymbol{W})$ also holds
almost surely for all large $n$, as claimed. ∎

Analogous to the argument of Appendix E.1, we may now condition on
$\boldsymbol{W}$, i.e. we assume that $\boldsymbol{W}$ is deterministic and
satisfies $\mathcal{E}(\boldsymbol{W})$ for all large $n$, and we write
$\mathbb{E}$ for the expectation over only the randomness of the new data
$(\tilde{}\boldsymbol{X},\tilde{}\boldsymbol{y})$. Defining

$\boldsymbol{G}=\sqrt{\frac{N}{n}}(\boldsymbol{Y}-\mathbb{E}\boldsymbol{Y})\in\mathbb{R}^{n\times
N},\quad\boldsymbol{u}=\frac{1}{\sqrt{n}}\tilde{}\boldsymbol{y}\in\mathbb{R}^{n}\quad\text{
where
}\quad\boldsymbol{Y}=\frac{1}{\sqrt{N}}\sigma(\tilde{}\boldsymbol{X}\boldsymbol{W}),$
(F.4)

observe that $[\boldsymbol{u},\boldsymbol{G}]\in\mathbb{R}^{n\times(N+1)}$ has
centered i.i.d. rows with respect to the randomness of
$(\tilde{}\boldsymbol{X},\tilde{}\boldsymbol{y})$. We will write
$\mathbb{E}_{\boldsymbol{x}}$ for the expectation with respect to a standard
Gaussian vector
$\boldsymbol{x}\sim\mathcal{N}(\mathbf{0},\boldsymbol{I}_{d})$.

######  Lemma 37.

Suppose $\boldsymbol{W}$ satisfies $\mathcal{E}(\boldsymbol{W})$ for all large
$n$. Then

$\displaystyle\|\mathbb{E}_{\boldsymbol{x}}[\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]\|_{2}\to
0,\qquad\|\mathbb{E}\boldsymbol{Y}\|\to
0,\qquad\mathopen{}\mathclose{{}\left\|{\boldsymbol{\Sigma}-\boldsymbol{\Sigma}_{\mathrm{lin}}}}\right\|\to
0$ (F.5)

where

$\displaystyle\boldsymbol{\Sigma}:=~{}$
$\displaystyle\mathbb{E}_{\boldsymbol{x}}[\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})^{\top}\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]-\mathbb{E}_{\boldsymbol{x}}[\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]^{\top}\mathbb{E}_{\boldsymbol{x}}[\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]$
(F.6) $\displaystyle\boldsymbol{\Sigma}_{\mathrm{lin}}:=~{}$ $\displaystyle
b_{\sigma}^{2}(\boldsymbol{W}^{\top}\boldsymbol{W})+(1-b_{\sigma}^{2})\boldsymbol{I}_{N}.$
(F.7)

Proof. The proof is the same as Lemmas 30 and 31. ∎

Proof of Theorem 5. We condition on $\boldsymbol{W}$ satisfying
$\mathcal{E}(\boldsymbol{W})$ for all large $n$, and we apply Theorem 12(c)
for $[\boldsymbol{u},\boldsymbol{G}]\in\mathbb{R}^{(n+1)\times N}$ (exchanging
$n$ and $N$). It may be checked that Assumption 5 holds for
$[\boldsymbol{u},\boldsymbol{G}]$ by the same argument as in Lemma 34.

By the convergence
$\|\boldsymbol{\Sigma}-\boldsymbol{\Sigma}_{\mathrm{lin}}\|\to 0$ in Lemma 37
and Proposition 35, if $\theta_{1}>\gamma_{0}^{1/4}$, then Assumption 6 holds
for $\boldsymbol{\Sigma}$ with $r=1$ and

$\nu=b_{\sigma}^{2}\otimes\rho^{\textnormal{MP}}_{\gamma_{0}}\oplus(1-b_{\sigma}^{2}),\qquad\lambda_{1}=b_{\sigma}^{2}\frac{(1+\theta_{1}^{2})(\gamma_{0}+\theta_{1}^{2})}{\theta_{1}^{2}}+(1-b_{\sigma}^{2})\notin\mathrm{supp}({\nu}),$

where $\rho_{\gamma_{0}}^{\textnormal{MP}}$ is the standard Marcenko-Pastur
limit for the empirical eigenvalue distribution of
$\boldsymbol{W}^{\top}\boldsymbol{W}$, hence $\nu$ is the limit empirical
eigenvalue distribution of $\boldsymbol{\Sigma}$, and $\lambda_{1}$ is the
limit of $\lambda_{\max}(\boldsymbol{\Sigma})$. If instead
$\theta_{1}\leq\gamma_{0}^{1/4}$, then Assumption 6 holds with $r=0$.

Then Theorem 12(a,c) characterizes the outlier eigenvalue and eigenvector of
$\boldsymbol{G}\boldsymbol{G}^{\top}$, showing:

  * •

$\boldsymbol{G}\boldsymbol{G}^{\top}$ has a spike eigenvalue if and only if
$\theta_{1}>\gamma_{0}^{1/4}$ and $z^{\prime}(-1/\lambda_{1})>0$, where
$z(\cdot)$ is defined by (A.4) with $\gamma=\gamma_{1}$ and the measure $\nu$
given above. In this case,
$\lambda_{\max}(\boldsymbol{G}\boldsymbol{G}^{\top})\to z(-1/\lambda_{1})$
almost surely.

  * •

When $\theta_{1}>\gamma_{0}^{1/4}$ and $z^{\prime}(-1/\lambda_{1})>0$, letting
$\boldsymbol{u}(\boldsymbol{G}),\boldsymbol{v}(\boldsymbol{\Sigma})$ be the
leading unit-norm left singular vector of $\boldsymbol{G}$ and leading unit-
norm eigenvector of $\boldsymbol{\Sigma}$, almost surely

$|\boldsymbol{u}^{\top}\boldsymbol{u}(\boldsymbol{G})|-\frac{\sqrt{z(-1/\lambda_{1})\varphi({-}1/\lambda_{1})}}{\lambda_{1}}\cdot\mathopen{}\mathclose{{}\left|\mathbb{E}_{\boldsymbol{x}}[\sigma_{*}(\boldsymbol{\beta}_{*}^{\top}\boldsymbol{x})\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]\boldsymbol{v}(\boldsymbol{\Sigma})}\right|\to
0.$

where $\varphi(\cdot)$ is defined by (B.6) also with $\gamma=\gamma_{1}$ and
the above measure $\nu$.

By an application of Weyl’s inequality and the Davis-Kahan Theorem as in the
proof of Theorem 3, this implies for
$\boldsymbol{K}=\boldsymbol{Y}\boldsymbol{Y}^{\top}$ that if
$\theta_{1}>\gamma_{0}^{1/4}$ and $z^{\prime}(-1/\lambda_{1})>0$, then its
leading eigenvalue $\lambda_{\max}(\boldsymbol{K})$ and unit eigenvector
$\widehat{}\boldsymbol{u}$ satisfy

$\displaystyle\lambda_{\max}(\boldsymbol{K})\to\gamma_{1}^{-1}z(-1/\lambda_{1}),$
(F.8)
$\displaystyle|\boldsymbol{u}^{\top}\widehat{}\boldsymbol{u}|-\frac{\sqrt{z(-1/\lambda_{1})\varphi({-}1/\lambda_{1})}}{\lambda_{1}}\cdot\mathopen{}\mathclose{{}\left|\mathbb{E}_{\boldsymbol{x}}[\sigma_{*}(\boldsymbol{\beta}_{*}^{\top}\boldsymbol{x})\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]\boldsymbol{v}(\boldsymbol{W})}\right|\to
0,$

where $\boldsymbol{v}(\boldsymbol{W})$ is the leading unit eigenvector of
$\boldsymbol{\Sigma}_{\text{lin}}$, i.e. the leading right singular vector of
$\boldsymbol{W}$. If $\theta_{1}\leq\gamma_{0}^{1/4}$ or
$z^{\prime}(-1/\lambda_{1})\leq 0$, then all eigenvalues of $\boldsymbol{K}$
converge to the support of its limiting empirical eigenvalue law.

Finally, in the case of $\theta_{1}>\gamma_{0}^{1/4}$ and
$z^{\prime}(-1/\lambda_{1})>0$, we may conclude the proof by showing

$\mathopen{}\mathclose{{}\left\|\mathbb{E}_{\boldsymbol{x}}[\sigma_{*}(\boldsymbol{\beta}_{*}^{\top}\boldsymbol{x})\sigma(\boldsymbol{x}^{\top}\boldsymbol{W})]-b_{\sigma}b_{\sigma_{*}}\boldsymbol{\beta}_{*}^{\top}\boldsymbol{W}}\right\|_{2}\to
0\text{ a.s.}$ (F.9)

For each column $i\in[N]$, we have from (F.3) that
$\|\boldsymbol{w}_{i}-\boldsymbol{w}_{i,0}\|_{2}\prec N^{-1/2}$, where
$\boldsymbol{w}_{i,0}\sim\mathcal{N}(0,d^{-1}\boldsymbol{I})$ and
$\boldsymbol{\beta}_{*}$ is deterministic. Hence
$(\boldsymbol{w}_{i},\boldsymbol{\beta}_{*})$ satisfy the approximate
orthonormality conditions $|\|\boldsymbol{w}_{i}\|_{2}-1|\prec N^{-1/2}$,
$\|\boldsymbol{\beta}_{*}\|_{2}-1=0$, and
$|\boldsymbol{w}_{i}^{\top}\boldsymbol{\beta}_{*}|\prec N^{-1/2}$. Then [FW20,
Lemma D.3(a)] implies

$\displaystyle\mathopen{}\mathclose{{}\left|\mathbb{E}_{\boldsymbol{x}}[\sigma_{*}(\boldsymbol{\beta}_{*}^{\top}\boldsymbol{x})\sigma(\boldsymbol{x}^{\top}\boldsymbol{w}_{i})]-b_{\sigma}b_{\sigma_{*}}\boldsymbol{\beta}_{*}^{\top}\boldsymbol{w}_{i}}\right|\prec
N^{-1}.$

(We note that [FW20, Lemma D.3(a)] assumes $\sigma=\sigma_{*}$, but the proof
is identical for $\sigma\neq\sigma_{*}$ both satisfying Assumption 3.)
Applying this to each coordinate $i\in[N]$ yields (F.9). Observe that
$\boldsymbol{\beta}_{*}^{\top}\boldsymbol{W}\boldsymbol{v}(\boldsymbol{W})=s_{\max}(\boldsymbol{W})\cdot\boldsymbol{\beta}_{*}^{\top}\boldsymbol{u}(\boldsymbol{W})$
where $s_{\max}(\boldsymbol{W})$ and $\boldsymbol{u}(\boldsymbol{W})$ are the
leading singular value and left singular vector of $\boldsymbol{W}$, and
recall from the definitions (F.4) that
$\boldsymbol{u}=\frac{1}{\sqrt{n}}\tilde{}\boldsymbol{y}$. Then we can apply
(F.9) and Proposition 35 to (F.8) to conclude that

$\frac{1}{\sqrt{n}}|\tilde{}\boldsymbol{y}^{\top}\widehat{\boldsymbol{u}}|\to
b_{\sigma}b_{\sigma_{*}}\frac{\sqrt{z(-1/\lambda_{1})\varphi(-1/\lambda_{1})}}{\lambda_{1}}\cdot\frac{\theta_{2}\sqrt{(\theta_{1}^{4}-\gamma_{0})(\gamma_{0}+\theta_{1}^{2})}}{\theta_{1}^{3}}>0\text{
a.s.}$

∎