text:
# 1 Introduction

## 1 Introduction

One of the most attractive features of the Minimal Supersymmetric extension of the Standard Model (MSSM) , is the fact that it provides a mechanism for breaking radiatively the electroweak gauge $`SU(2)_L\times U(1)_Y`$ symmetry down to $`U(1)_{\mathrm{EM}}`$. It was first shown  that a supersymmetry (SUSY) breaking term for the gluino can induce an effective potential which spontaneously breaks the electroweak symmetry. At the same time, a mechanism relying on the renormalization group evolution from a grand unification (GUT) scale $`M_{GUT}`$ down to the weak scale was proposed . In this framework, at the scale $`M_{GUT}`$, the parameters entering the scalar potential of the MSSM obey simple boundary conditions dictated by the underlying theory of SUSY breaking, and the electroweak symmetry is unbroken. When the parameters are evolved down to the weak scale by means of the MSSM renormalization group equations (RGE), which amounts to resumming the leading logarithmic corrections to all orders, the soft SUSY–breaking mass $`m_{H_2}^2`$ is driven towards negative values, due to corrections controlled by the top Yukawa coupling $`h_t`$. This helps to destabilize the origin in field space, so that the Higgs fields acquire non–vanishing vacuum expectation values (VEVs) and the electroweak symmetry is spontaneously broken. Although the studies in Refs.  where differing on the initial boundary conditions, the result was one: the radiative electroweak symmetry breaking (REWSB) takes place if the top Yukawa coupling is large, such that $`60\mathrm{GeV}\stackrel{<}{}m_t\stackrel{<}{}200\mathrm{GeV}`$, with the upper bound coming from the requirement that $`h_t`$ remains in the perturbative range up to the GUT scale. It could be a coincidence that the top quark is found at the Tevatron to have mass around 175 GeV, but certainly this is consistent with the REWSB mechanism in the MSSM.

In the RGE–improved potential of the MSSM employed at tree level, the VEVs of the Higgs fields, and the occurrence of spontaneous symmetry breaking itself, depend critically on the renormalization scale at which the parameters entering the potential are computed; an inappropriate choice of that scale can lead to results that are even qualitatively wrong. In fact, the electroweak symmetry is either broken or unbroken, independently of the renormalization scale choice, and the critical behavior described above is just an artifact of the tree–level approximation. The correct way of determining the ground state of the theory is to minimize the Coleman–Weinberg effective potential , i.e. the tree–level potential plus a correction coming from the sum of all the one–loop diagrams with zero–momentum external lines. Since Refs.  this procedure has become standard in the renormalization group analyses of the MSSM (for early examples see Refs. ).

The effective potential is also a useful tool for computing the leading corrections to the MSSM Higgs masses, both at the one loop  and the two loop  level111Other two–loop computations of the MSSM Higgs masses have been performed in the renormalization group  and diagrammatic  approaches., in the approximation of zero external momentum. The leading one–loop corrections are $`𝒪(\alpha _t)`$, i.e. they are controlled by the top Yukawa coupling $`h_t\sqrt{4\pi \alpha _t}`$. For stop masses of $`𝒪(1\mathrm{TeV})`$, such corrections increase by 40–60 GeV the mass $`m_h`$ of the lightest Higgs boson (which at tree level must be lighter than $`m_Z`$), allowing it to escape the direct searches at LEP. Also, the leading two–loop corrections have sizeable effects: the $`𝒪(\alpha _t\alpha _s)`$ corrections, controlled by the strong gauge coupling $`g_s\sqrt{4\pi \alpha _s}`$, typically reduce $`m_h`$ by 15–20 GeV, whereas the $`𝒪(\alpha _t^2)`$ ones may increase it by up to 7–8 GeV.

Motivated by the relevance of the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ two–loop corrections in the case of the Higgs masses, we study in this paper the effect of the same corrections on the electroweak symmetry breaking conditions. The contributions to the two–loop MSSM effective potential that are relevant to the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections have been discussed in Refs. , and a complete computation of the two–loop effective potential has been presented in Ref. . However, practical studies of REWSB usually require explicit formulae for the Higgs tadpole diagrams, i.e. the first derivatives of the effective potential with respect to the Higgs fields. Such formulae are presently available at the one–loop order , but they have not been presented so far at the two–loop order. Using the techniques developed in Ref. , we compute in this paper explicit and compact analytical expressions for the two–loop $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ part of the tadpoles. As a byproduct from our $`𝒪(\alpha _t\alpha _s)`$ corrections we obtain also the $`𝒪(\alpha _b\alpha _s)`$ corrections, that are relevant for large values of $`\mathrm{tan}\beta `$. Once we assume that the electroweak symmetry is indeed broken, giving rise to the observed value of the $`Z`$ boson mass, the corrections to the tadpoles translate into $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections to the values of $`\mu `$, the Higgs mass term in the superpotential, and $`m_A`$, the running mass of the $`A`$ boson. We discuss the effect of our two–loop corrections in the framework of gravity mediated SUSY breaking , also denoted as minimal supergravity (mSUGRA), referring in particular to various “benchmark” scenarios suggested at Snowmass . We find that the inclusion of the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections significantly improves the renormalization scale dependence of the results, and that partial cancellations occur between the $`𝒪(\alpha _t\alpha _s)`$ corrections and the $`𝒪(\alpha _t^2)`$ ones. Our corrections are also required for consistency in the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ two–loop computation of the MSSM Higgs masses, if the input parameters are computed via renormalization group evolution from a set of high energy boundary conditions.

The paper is organized as follows: in section 2 we recall the basic concepts of radiative electroweak symmetry breaking, and introduce some notation which will be used in the rest of the paper; in section 3 we describe the main features of our $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ computation of the two–loop tadpoles; in section 4 we discuss the numerical effect of our corrections, and we show how they improve the dependence of $`\mu `$ and $`m_A`$ on the renormalization scale at which the effective potential is minimized; section 5 contains our conclusions. In addition, we present in the appendix A some useful formulae for the integrals entering the two–loop effective potential, and in the appendix B the explicit analytical formulae for the $`𝒪(\alpha _t\alpha _s)`$ part of the corrections. The formulae for the $`𝒪(\alpha _t^2)`$ part are indeed rather long, thus we make them available, upon request, in the form of a computer code 222E–mail: slavich@mppmu.mpg.de.

## 2 Radiative electroweak symmetry breaking

We start our discussion from the tree–level scalar potential of the MSSM, that reads, keeping only the dependence on the neutral Higgs fields $`H_1^0`$ and $`H_2^0`$:

$$V_0=\mathrm{\Lambda }+m_1^2|H_1^0|^2+m_2^2|H_2^0|^2+m_3^2(H_1^0H_2^0+\mathrm{h}.\mathrm{c}.)+\frac{g^2+g^{\mathrm{\hspace{0.17em}2}}}{8}(|H_1^0|^2|H_2^0|^2)^2,$$
(1)
where: $`\mathrm{\Lambda }`$ is a field–independent vacuum energy; $`m_1^2=m_{H_1}^2+\mu ^2`$, $`m_2^2=m_{H_2}^2+\mu ^2`$ (we assume $`\mu `$ to be real, neglecting all possible CP–violating phases); $`m_{H_1}^2`$, $`m_{H_2}^2`$ and $`m_3^2`$ are soft SUSY–breaking masses; $`g`$ and $`g^{}`$ are the $`SU(2)_L`$ and $`U(1)_Y`$ gauge couplings, respectively. At the classical level, the mass parameters entering $`V_0`$ must satisfy the following conditions:

$$m_1^2+m_2^2\mathrm{\hspace{0.33em}2}|m_3^2|,m_1^2m_2^2m_3^4.$$
(2)
The first condition guarantees that the potential is bounded from below; the second condition destabilizes the origin in field space, making sure that the neutral components of the Higgs fields acquire non–vanishing VEVs $`H_1^0v_1/\sqrt{2}`$ and $`H_2^0v_2/\sqrt{2}`$. It is not restrictive to choose $`m_3^2`$ real and negative333Our conventions differ by a sign in the parameters $`\mu `$ and $`m_3^2`$ with respect to those used in the second paper of Ref.., so that $`v_1`$ and $`v_2`$ are real and positive, and the neutral Higgs fields can be decomposed into their VEVs plus their CP–even and CP–odd fluctuations as $`H_i^0=(v_i+S_i+iP_i)/\sqrt{2}`$.

Since the parameters entering $`V_0`$ are taken as ”running” ones (i.e., they vary with the renormalization scale), also the validity of the conditions in Eq. (2) depends on the scale, as well as the numerical values of $`v_1`$ and $`v_2`$. As discussed in Ref. , the minimization of the tree–level potential may lead to grossly inaccurate results, unless the renormalization scale is chosen in such a way that the radiative corrections to the scalar potential are small. To obtain the correct results, one should rather minimize the effective potential $`V_{\mathrm{eff}}`$, defined as:

$$V_{\mathrm{eff}}=V_0+\mathrm{\Delta }V,$$
(3)
where $`\mathrm{\Delta }V`$ contains the radiative corrections to the scalar potential $`V_0`$. The minimization conditions for $`V_{\mathrm{eff}}`$ can be written as:

$`{\displaystyle \frac{1}{v_1}}{\displaystyle \frac{V_{\mathrm{eff}}}{S_1}}|_{\mathrm{min}}`$ $`=`$ $`m_{H_1}^2+\mu ^2+{\displaystyle \frac{g^2+g^{\mathrm{\hspace{0.17em}2}}}{4}}(v_1^2v_2^2)+m_3^2{\displaystyle \frac{v_2}{v_1}}+\mathrm{\Sigma }_1=\mathrm{\hspace{0.33em}0},`$ (4)
$`{\displaystyle \frac{1}{v_2}}{\displaystyle \frac{V_{\mathrm{eff}}}{S_2}}|_{\mathrm{min}}`$ $`=`$ $`m_{H_2}^2+\mu ^2+{\displaystyle \frac{g^2+g^{\mathrm{\hspace{0.17em}2}}}{4}}(v_2^2v_1^2)+m_3^2{\displaystyle \frac{v_1}{v_2}}+\mathrm{\Sigma }_2=\mathrm{\hspace{0.33em}0},`$ (5)
where the “tadpoles” $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$ are defined as:

$$\mathrm{\Sigma }_i\frac{1}{v_i}\frac{\mathrm{\Delta }V}{S_i}|_{\mathrm{min}}.$$
(6)
In principle, a renormalization group study of the MSSM should start from some large scale $`M_{GUT}`$, where the input parameters have a simple structure dictated by the underlying theory of SUSY breaking, and the electroweak symmetry is unbroken. The parameters are then evolved, by means of appropriate renormalization group equations, down to some lower scale, where the electroweak symmetry breaking occurs and the VEVs $`v_1`$ and $`v_2`$ can be obtained by solving Eqs. (4)–(5). A set of high–energy input parameters is then acceptable if it leads to the correct value of the squared running mass for the $`Z`$ boson, $`m_Z^2=(g^2+g^{\mathrm{\hspace{0.17em}2}})(v_1^2+v_2^2)/4`$. However, in most practical applications of the renormalization group procedure, it is more convenient to assume that there is successful electroweak symmetry breaking, and trade two of the high–energy input parameters for $`v_1`$ and $`v_2`$ (or, equivalently, for $`v^2v_1^2+v_2^2`$ and $`\mathrm{tan}\beta v_2/v_1`$). Eqs. (4)–(5) can thus be rephrased into the following conditions among the parameters at the weak scale:

$`\mu ^2`$ $`=`$ $`{\displaystyle \frac{m_Z^2}{2}}+{\displaystyle \frac{m_{H_1}^2+\mathrm{\Sigma }_1(m_{H_2}^2+\mathrm{\Sigma }_2)\mathrm{tan}^2\beta }{\mathrm{tan}^2\beta 1}},`$ (7)
$`m_3^2`$ $`=`$ $`{\displaystyle \frac{1}{2}}\mathrm{sin}2\beta \left(m_{H_1}^2+m_{H_2}^2+2\mu ^2+\mathrm{\Sigma }_1+\mathrm{\Sigma }_2\right),`$ (8)
i.e., the terms proportional to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$ in the above equations can be viewed as the radiative corrections to the values of $`\mu ^2`$ and $`m_3^2`$ obtained from the requirement of successful electroweak symmetry breaking. We recall that $`m_3^2`$ is related to the squared running mass for the $`A`$ boson through $`m_A^2=2m_3^2/\mathrm{sin}2\beta `$. Notice also that the sign of $`\mu `$ is not fixed by Eq. (7), and it must be supplemented as an additional input quantity. If the right side of Eq. (7) is such that $`\mu ^2`$ is negative, then our choice of input parameters is inconsistent, and the electroweak symmetry fails to be broken. We remark in passing that the choice of the input parameters is constrained by further requirements: it must lead to a spectrum of physical masses for the MSSM superpartners and Higgs bosons compatible with the present experimental lower bounds, and such that the lightest supersymmetric particle (LSP) is electrically neutral; it must satisfy phenomenological constraints coming from radiative B-meson decays, muon anomalous magnetic moment and cosmological relic density; finally, it must guarantee that the MSSM scalar potential is bounded from below and does not lead to charge and color breaking minima. However, a detailed study of the (theoretically and experimentally) allowed regions in the MSSM parameter space goes beyond the scope of this paper, and will not be pursued in the following.

The full one–loop corrections to the REWSB conditions have been extensively discussed in the literature  in the framework of the mSUGRA scenario. The dominant one–loop contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$ come to a large extent from the top/stop (and, for large $`\mathrm{tan}\beta `$, bottom/sbottom and tau/stau) diagrams. The contributions of the diagrams involving charginos and neutralinos can also be sizeable and comparable to the top/stop ones in some regions of the parameter space, while the Higgs and gauge bosons and the first two generations of (s)quarks and (s)leptons give only subdominant corrections.

In the following sections we provide explicit analytical formulae for the two–loop top/stop contributions to the tadpoles $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$, resulting into $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections to $`\mu ^2`$ and to the running mass $`m_A^2`$. In analogy with the case of the Higgs masses, we expect such corrections to be the leading ones at the two–loop level, giving rise to sizeable effects at least in some regions of the MSSM parameter space. In addition, most public codes that compute the MSSM mass spectrum from a set of unified parameters at the scale $`M_{GUT}`$, such as SuSpect , SoftSusy , SPheno  and FeynSSG , include a two–loop $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ computation of the Higgs masses,444Another widely used public code, Isajet 7.58 , relies on a one–loop effective potential computation of the Higgs masses. but employ one–loop results for the tadpoles (see also Ref.  for a recent discussion). Since $`m_A^2`$ enters the tree–level mass matrix of the CP–even Higgs bosons, the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections to $`m_A^2`$ should be included in those codes for consistency.

## 3 Computation of the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections

We shall now describe our two–loop, $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ computation of the tadpoles $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$, involving the first derivatives of $`\mathrm{\Delta }V`$ with respect to the CP–even parts of the neutral Higgs fields \[see Eq. (6)\]. The computation is consistently performed by setting to zero all the gauge couplings but $`g_s`$ and by keeping $`h_t`$ as the only non–vanishing Yukawa coupling (with a slight abuse of language, in the following we will refer to this approximation as to the gaugeless limit). In this limit, the tree–level (field–dependent) spectrum of the MSSM simplifies considerably: gauginos and Higgsinos do not mix; the charged and neutral Higgsinos combine into Dirac spinors $`\stackrel{~}{h}^0`$ and $`\stackrel{~}{h}^\pm `$ with degenerate mass $`\mu `$; the gaugino masses coincide with the soft SUSY–breaking parameters $`M_A(A=1,2,3)`$ (among them, only the gluino mass $`m_{\stackrel{~}{g}}=M_3`$ is relevant to our calculation); the only massive Standard Model (SM) fermion is the top quark; all other fermions and gauge bosons have vanishing masses; besides the top squarks, the only sfermion with non–vanishing couplings is the bottom squark $`\stackrel{~}{b}_L`$; the lighter CP–even Higgs boson, $`h`$, is massless, and the same is true for the Goldstone bosons $`G`$ and $`G^\pm `$; all the remaining Higgs states, $`(H,A,H^\pm )`$, have degenerate mass eigenvalues $`m_A^2`$. The tree–level mixing angle in the CP–even sector is just $`\alpha =\beta \pi /2`$.

To begin with, we address the renormalization of the effective potential. In the loop expansion, the correction to the effective potential can be decomposed as $`\mathrm{\Delta }V=V_1\mathrm{}+V_2\mathrm{}+\mathrm{}`$, where the ellipsis stand for higher loops. Using the Landau gauge and dimensional reduction  in $`d=42ϵ`$ dimensions, and including for later convenience also terms that vanish when $`ϵ0`$, the unrenormalized one–loop effective potential reads:

$$V_1\mathrm{}=\frac{1}{64\pi ^2}\mathrm{Str}^4\left[\frac{1}{ϵ}+\frac{3}{2}\mathrm{ln}\frac{^2}{Q^2}+ϵ\left(\frac{7}{4}\frac{3}{2}\mathrm{ln}\frac{^2}{Q^2}+\frac{1}{2}\mathrm{ln}^2\frac{^2}{Q^2}+\frac{\pi ^2}{12}\right)\right],$$
(9)
where $`^2`$ is the matrix of the field–dependent squared masses, and the supertrace of a generic function $`f(^2)`$ is defined as a sum over the eigenvalues $`m_i^2`$:

$$\mathrm{Str}f(^2)=\underset{i}{}(1)^{2s_i}(2s_i+1)f(m_i^2),$$
(10)
where $`s_i`$ is the spin of the corresponding particles. In Eq. (9), $`Q^2=4\pi \mu ^2e^{\gamma _E}`$, i.e. the finite terms that are removed together with $`1/ϵ`$ in the modified subtraction schemes have been reabsorbed in the renormalization scale (the same convention will be adopted in the following). In the gaugeless limit described above, only the top and stop contributions to $`V_1\mathrm{}`$ are relevant, giving rise to $`𝒪(\alpha _t)`$ contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$.

The Feynman diagrams that contribute to the two–loop effective potential $`V_2\mathrm{}`$ and give rise to $`𝒪(\alpha _t\alpha _s)`$ contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$ are shown in Fig. 1, while the diagrams relevant to the $`𝒪(\alpha _t^2)`$ contributions are shown in Fig. 2. The corresponding analytical formulae for $`V_2\mathrm{}`$ can be found e.g. in the last paper of Ref.  555The formulae of Ref.  are obtained for vanishing CP–odd fields. While requiring some modifications  for the computation of the CP–even Higgs boson masses in terms of the physical $`m_A`$, those formulae can be used as they stand for the purposes of the present analysis.. These formulae involve two basic integrals, $`I(m_1^2,m_2^2,m_3^2)`$ and $`J(m_1^2,m_2^2)`$, that have been evaluated with different methods in Refs. . Explicit expressions for $`I`$ and $`J`$ in the formalism of Ref.  are presented in the appendix A.

To carry out the renormalization of $`V_{\mathrm{eff}}`$, at the two–loop order and in the $`\overline{\mathrm{DR}}`$ scheme , we start from the unrenormalized effective potential, written in terms of generic bare parameters $`x_i`$. Then, we expand the parameters as $`x_i=x_i^{\overline{\mathrm{DR}}}+\delta x_i`$, where $`\delta x_i`$ are purely divergent quantities, so that all the poles in $`1/ϵ`$ and $`1/ϵ^2`$ are cancelled. After taking the limit $`ϵ0`$, the non vanishing part of the renormalized effective potential is:

$$V_{\mathrm{eff}}=V_0(x_i^{\overline{\mathrm{DR}}})+V_1\mathrm{}^{(1)}(x_i^{\overline{\mathrm{DR}}})+V_2\mathrm{}^{(1)}+\frac{V_1\mathrm{}^{(ϵ)}}{x_i}\delta ^1\mathrm{}x_i,$$
(11)
where $`V_1\mathrm{}^{(1)}`$ and $`V_2\mathrm{}^{(1)}`$ denote the finite parts of the one–loop and two–loop effective potential, respectively, $`V_1\mathrm{}^{(ϵ)}`$ denotes the terms proportional to $`ϵ`$ in Eq. (9), and $`\delta ^1\mathrm{}x_i`$ is the coefficient of $`1/ϵ`$ in the one–loop part of the generic counterterm (notice that we need to compute explicitly only the one–loop counterterms for the top and stop masses). In Eq. (11), $`V_0`$ and $`V_1\mathrm{}^{(1)}`$ are expressed in terms of $`\overline{\mathrm{DR}}`$–renormalized parameters, while the renormalization of the parameters entering the two–loop part is irrelevant, amounting to a higher–order (i.e., three–loop) effect. In summary, it is possible to define the renormalized two–loop effective potential as:

$$\widehat{V}_2\mathrm{}=V_2\mathrm{}^{(1)}+\frac{V_1\mathrm{}^{(ϵ)}}{x_i}\delta ^1\mathrm{}x_i.$$
(12)
We have checked that $`\widehat{V}_2\mathrm{}`$ corresponds to the finite part of the potential obtained by replacing the integrals $`I`$ and $`J`$ in $`V_2\mathrm{}`$ with the “subtracted” integrals $`\widehat{I}`$ and $`\widehat{J}`$, first introduced in Ref. . More precisely, this is true only up to terms that give a null contribution to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$, unless we include in $`V_{\mathrm{eff}}`$ also diagrams that do not depend on the Higgs fields (such as, e.g., the one–loop diagram involving gluinos and the two–loop diagram involving gluinos and gluons).

Compact analytical formulae for the derivatives of the renormalized effective potential in the gaugeless limit can be obtained with a procedure similar to that of Ref. . The relevant field–dependent quantities are the top mass $`m_t`$, the stop masses $`m_{\stackrel{~}{t}_1}^2`$ and $`m_{\stackrel{~}{t}_2}^2`$, and the stop mixing angle $`\theta _{\stackrel{~}{t}}`$ (the top and stop phases $`\phi `$ and $`\stackrel{~}{\phi }`$, introduced in  to take into account the dependence on the CP–odd part of the Higgs fields, do not enter the computation of $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$). At the minimum of the effective potential, the parameters in the stop sector are related by:

$$\mathrm{sin}2\theta _{\stackrel{~}{t}}=\frac{2m_t(A_t+\mu \mathrm{cot}\beta )}{m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2},$$
(13)
where $`A_t`$ is the soft SUSY–breaking trilinear coupling of the stops \[notice that Eq. (13) defines our convention for the sign of $`\mu `$\]. We will use in the following the shortcuts $`c_{2\theta }\mathrm{cos}2\theta _{\stackrel{~}{t}}`$ and $`s_{2\theta }\mathrm{sin}2\theta _{\stackrel{~}{t}}`$. After a straightforward application of the chain rule for the derivatives of the effective potential, we get:

$`v_1^2\mathrm{\Sigma }_1`$ $`=`$ $`m_t\mu \mathrm{cot}\beta s_{2\theta }F,`$ (14)
$`v_2^2\mathrm{\Sigma }_2`$ $`=`$ $`m_tA_ts_{2\theta }F+2m_t^2G.`$ (15)
The functions $`F`$ and $`G`$ are combinations of the derivatives of $`\mathrm{\Delta }V`$ with respect to the field–dependent parameters, computed at the minimum of the effective potential:

$`F`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Delta }V}{m_{\stackrel{~}{t}_1}^2}}{\displaystyle \frac{\mathrm{\Delta }V}{m_{\stackrel{~}{t}_2}^2}}{\displaystyle \frac{2c_{2\theta }}{m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2}}{\displaystyle \frac{\mathrm{\Delta }V}{c_{2\theta }}},`$ (16)
$`G`$ $`=`$ $`{\displaystyle \frac{\mathrm{\Delta }V}{m_t^2}}+{\displaystyle \frac{\mathrm{\Delta }V}{m_{\stackrel{~}{t}_1}^2}}+{\displaystyle \frac{\mathrm{\Delta }V}{m_{\stackrel{~}{t}_2}^2}}.`$ (17)
The one–loop parts of $`F`$ and $`G`$, giving rise to $`𝒪(\alpha _t)`$ contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$, are easily computed from the derivatives of $`V_1\mathrm{}^{(1)}`$. In units of $`N_c/(16\pi ^2)`$, where $`N_c=3`$ is a color factor, they read:

$`F^1\mathrm{}`$ $`=`$ $`m_{\stackrel{~}{t}_1}^2\left(\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}1\right)m_{\stackrel{~}{t}_2}^2\left(\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_2}^2}{Q^2}}1\right),`$ (18)
$`G^1\mathrm{}`$ $`=`$ $`m_{\stackrel{~}{t}_1}^2\left(\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}1\right)+m_{\stackrel{~}{t}_2}^2\left(\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_2}^2}{Q^2}}1\right)2m_t^2\left(\mathrm{ln}{\displaystyle \frac{m_t^2}{Q^2}}1\right).`$ (19)
Although a naive, brute–force computation of the derivatives of the renormalized two–loop potential $`\widehat{V}_2\mathrm{}`$ presents no major conceptual difficulties, the number of terms involved blows up quickly, giving rise to very long and complicated analytical expressions. However, in the spirit of Ref. , it is possible to obtain recursive relations for the derivatives of the integral $`I(m_1^2,m_2^2,m_3^2)`$ with respect to the internal masses (see the appendix A), that simplify considerably the results. In this way, we obtained compact analytical formulae for the two–loop parts of $`F`$ and $`G`$, giving rise to $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$. As a non trivial check of the correctness of our computation, we have verified that the quantities $`\mu ^2`$ and $`m_3^2`$ defined in Eqs. (7)–(8) obey the appropriate two–loop RGE , specialized to the gaugeless limit:

$`{\displaystyle \frac{\mu ^2}{\mathrm{ln}Q^2}}`$ $`=`$ $`{\displaystyle \frac{3\alpha _t}{4\pi }}\mu ^2+{\displaystyle \frac{\alpha _t\alpha _s}{\pi ^2}}\mu ^2{\displaystyle \frac{9\alpha _t^2}{16\pi ^2}}\mu ^2,`$ (20)
$`{\displaystyle \frac{m_3^2}{\mathrm{ln}Q^2}}`$ $`=`$ $`{\displaystyle \frac{3\alpha _t}{4\pi }}\left({\displaystyle \frac{m_3^2}{2}}+A_t\mu \right)+{\displaystyle \frac{\alpha _t\alpha _s}{\pi ^2}}\left({\displaystyle \frac{m_3^2}{2}}+A_t\mu m_{\stackrel{~}{g}}\mu \right){\displaystyle \frac{9\alpha _t^2}{16\pi ^2}}\left({\displaystyle \frac{m_3^2}{2}}+2A_t\mu \right).`$
We have also checked explicitly that, in the special supersymmetric limit in which all the soft–SUSY breaking parameters as well as $`\mu `$ are set to zero, such that $`m_{\stackrel{~}{t}_1}=m_{\stackrel{~}{t}_2}=m_t`$, the two–loop parts of $`F`$ and $`G`$ are indeed vanishing.

For illustrative purposes, we present in the appendix B the explicit formulae for the $`𝒪(\alpha _t\alpha _s)`$ corrections, valid for arbitrary values of the relevant MSSM parameters. The corresponding formulae for the $`𝒪(\alpha _t^2)`$ corrections are indeed rather long, thus we make them available, upon request, in the form of a Fortran code.

The computation described above allows us to obtain also the two–loop $`𝒪(\alpha _b\alpha _s)`$ corrections 666In contrast, the $`𝒪(\alpha _b^2+\alpha _t\alpha _b)`$ corrections would require a dedicated computation. induced by the bottom/sbottom sector ($`h_b\sqrt{4\pi \alpha _b}`$ being the bottom Yukawa coupling), that can be relevant for large values of $`\mathrm{tan}\beta `$. To this purpose, the substitutions $`tb,v_1v_2`$ (i.e., $`\mathrm{tan}\beta \mathrm{cot}\beta `$) and $`\mathrm{\Sigma }_1\mathrm{\Sigma }_2`$ must be performed in the corresponding formulae for the $`𝒪(\alpha _t\alpha _s)`$ part of the corrections. The complications relative to the on–shell definition of the sbottom parameters, discussed in Ref. , do not arise in this case since we are working in the $`\overline{\mathrm{DR}}`$ renormalization scheme. However, the $`\mathrm{tan}\beta `$–enhanced threshold corrections  to the relation between $`h_b`$ and the bottom mass $`m_b`$ must be resummed to all orders  in a redefinition of the bottom Yukawa coupling (see e.g. Ref.  for the details).

## 4 Numerical results

In this section we discuss the numerical effect of our $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ two–loop corrections on the minimization conditions of the MSSM effective potential.

For definiteness, we work in the mSUGRA scenario, in which the MSSM Lagrangian at the large scale $`M_{GUT}`$ contains only five independent mass parameters: a common soft SUSY-breaking scalar mass $`m_0`$, a common soft gaugino mass $`m_{1/2}`$, a common soft trilinear term $`A_0`$, the superpotential Higgs mixing parameter $`\mu _0`$ and its soft SUSY-breaking counterpart $`B_0`$ (the subscript “0” denotes the fact that the parameters are computed at the boundary scale). The soft Higgs mixing parameter $`B`$ has the dimensions of a squared mass, and is defined in such a way that in the low–energy Higgs potential of Eq. (1) it coincides with $`m_3^2`$ (to avoid confusion, we will refer to $`B`$ as to $`m_3^2`$ from now on). As anticipated in section 2, rather than providing input values for all the five mass parameters at the GUT scale, we assume that the electroweak symmetry is successfully broken at the weak scale, and we trade $`\mu _0`$ and $`m_3^2(M_{GUT})`$ for the weak scale input parameters $`v`$ and $`\mathrm{tan}\beta `$.

Before discussing our results, it is useful to describe in some detail the numerical procedure for the renormalization group evolution of the MSSM parameters. We start by defining the $`SU(3)\times SU(2)_L\times U(1)_Y`$ gauge couplings at the weak scale (which we identify with the pole $`Z`$ boson mass, $`M_Z=91.187`$ GeV), from the running weak mixing angle $`\widehat{s}^{\mathrm{\hspace{0.17em}2}}=0.2315`$, the electromagnetic coupling $`\widehat{\alpha }_{EM}=1/127.9`$ and the strong coupling $`\alpha _s=0.119`$. The electroweak symmetry breaking parameter, $`v^2=v_1^2+v_2^2`$, is defined in terms the muon decay constant according to the relation $`v=(\sqrt{2}G_\mu )^{1/2}=246.218`$ GeV, and then translated to the $`\overline{\mathrm{DR}}`$ scheme by means of the formulae of Ref. . In addition, $`\mathrm{tan}\beta =v_2/v_1`$ is taken as an input parameter at the weak scale, allowing us to determine $`v_1`$ and $`v_2`$. The Yukawa couplings of the light SM fermions are obtained from the corresponding masses at the scale $`Q=2`$ GeV, and then evolved up to $`M_Z`$ by means of the two–loop SM RGEs . In the case of the bottom coupling, the $`\mathrm{tan}\beta `$–enhanced threshold corrections  to the relation between $`h_b`$ and $`m_b`$ are included at the scale $`M_Z`$ according to the formulae of Ref. . Finally, the top Yukawa coupling is also defined at the scale $`M_Z`$ through the relation $`h_t=\sqrt{2}m_t/v_2`$, where $`m_t`$ is the $`\overline{\mathrm{DR}}`$ mass for the top quark, obtained from the pole mass $`M_t=174.3`$ GeV by means of the formulae of Ref. .

The evolution of the MSSM parameters from the weak scale to the GUT scale, and back, is performed by means of the full one–loop MSSM RGEs. However, for consistency with our two–loop $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ analysis of the REWSB conditions, we supplement the RGEs with the two–loop strong and top–Yukawa contributions as given in Ref. . We also make use of the two–loop RGEs in the Landau gauge for the VEVs, following Refs. . As a first step, we evolve the gauge and Yukawa couplings from the scale $`M_Z`$ up to the scale where the $`SU(2)_L\times U(1)_Y`$ gauge couplings $`g_1`$ and $`g_2`$ meet, that we identify with $`M_{GUT}`$ (we do not force the strong gauge coupling $`g_3`$ to meet $`g_1`$ and $`g_2`$). At this scale, which turns out to be of the order of $`10^{16}`$ GeV, we set the input boundary conditions for the soft SUSY–breaking masses $`m_0,m_{1/2}`$ and $`A_0`$. The values of $`\mu _0`$ and $`m_3^2(M_{GUT})`$, to be later determined from the REWSB conditions, are provisionally set to zero. At this point, we start an iterative procedure: first we run the MSSM parameters from $`M_{GUT}`$ down to some scale $`Q_{\mathrm{min}}`$, of the order of the weak scale, where the values of $`\mu ^2(Q_{\mathrm{min}})`$ and $`m_3^2(Q_{\mathrm{min}})`$ are computed through Eqs. (7) and (8), with the sign of $`\mu `$ supplied as an extra input parameter. Then we run all the parameters, including $`\mu `$ and $`m_3^2`$, down to the scale $`M_Z`$, which we regard as the end point of the RGE evolution. At this scale, we compute the threshold corrections to the top and bottom Yukawa couplings, using the newly obtained values of the relevant MSSM parameters. Finally, we run all the parameters back to the scale $`M_{GUT}`$, where the resulting values for $`\mu `$ and $`m_3^2`$ are taken as new guesses for the corresponding boundary conditions. We iterate the procedure until convergence is reached, i.e. the values of $`\mu ^2(Q_{\mathrm{min}})`$ and $`m_3^2(Q_{\mathrm{min}})`$ obtained from the RGE evolution of $`\mu _0`$ and $`m_3^2(M_{GUT})`$ coincide with those obtained from the minimization conditions (7) and (8). If however $`\mu ^2`$ turns out to be negative, then our choice of input parameters is inconsistent (i.e., it does not lead to successful REWSB) and must be discarded.

In order to discuss the effect of our $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections to the REWSB conditions, we show in Figs. 37 the values of $`|\mu (M_Z)|`$ and $`m_A(M_Z)`$, the latter obtained through $`m_A^2(M_Z)=2m_3^2(M_Z)/\mathrm{sin}2\beta `$, as functions of the minimization scale $`Q_{\mathrm{min}}`$. Stability of the results with respect to moderate changes in $`Q_{\mathrm{min}}`$ (which should anyway lie in the weak range, i.e. between $`M_Z`$ and a few TeV) indicates that the higher–order corrections not included in the computation of $`\mu `$ and $`m_A`$ are small and can be safely neglected. We will see that in general the inclusion of our $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections significantly improves the scale–dependence of the results.

In the choice of the input parameters, we refer to the so–called Snowmass Points , which represent typical “benchmark” scenarios that are commonly investigated in the phenomenological analyses of the mSUGRA parameter space. In particular, Figs. 36 correspond to:

$$\begin{array}{cccccc}\mathrm{SPS}\mathrm{\hspace{0.33em}1}\mathrm{a}:\hfill & m_0=100\mathrm{GeV},& m_{1/2}=250\mathrm{GeV},& A_0=100\mathrm{GeV},& \mathrm{tan}\beta =10,& \mu <0,\\ \mathrm{SPS}\mathrm{\hspace{0.33em}3}:\hfill & m_0=90\mathrm{GeV},& m_{1/2}=400\mathrm{GeV},& A_0=0,& \mathrm{tan}\beta =10,& \mu <0,\\ \mathrm{SPS}\mathrm{\hspace{0.33em}4}:\hfill & m_0=400\mathrm{GeV},& m_{1/2}=300\mathrm{GeV},& A_0=0,& \mathrm{tan}\beta =50,& \mu <0,\\ \mathrm{SPS}\mathrm{\hspace{0.33em}5}:\hfill & m_0=150\mathrm{GeV},& m_{1/2}=300\mathrm{GeV},& A_0=1\mathrm{TeV},& \mathrm{tan}\beta =5,& \mu <0,\end{array}$$

respectively. Notice that our convention for the sign of $`\mu `$ differs from the one in Ref. , where a discussion on the characteristics of the various scenarios can be found. A further scenario, denoted as “Focus point” (SPS2) and characterized by a common scalar mass, $`m_0=1450`$ GeV, much larger than the common fermion mass, $`m_{1/2}=300`$ GeV, has also been suggested in Ref. . However, we found that in this scenario the results for $`\mu `$ and $`m_A`$, including the qualitative effect of the various corrections and the occurrence of REWSB itself, depend dramatically on very small adjustments of the input value for the top pole mass (e.g. $`|\mu (M_Z)|`$ varies roughly between 400 and 100 GeV if $`M_t`$ is varied between 174.3 and 175 GeV). The extreme sensitivity of the SPS2 scenario on the input top mass has already been discussed in Ref. . Since this scenario appears to lead to unstable results, we will not consider it further in this work. However, in order to investigate the situation in which the common scalar mass is considerably larger than the common fermion mass, we show in Fig. 7 a scenario with $`m_0=1`$ TeV and $`m_{1/2}=300`$ GeV (the other parameters being chosen as $`A_0=0,\mathrm{tan}\beta =10`$ and $`\mu <0`$). We have checked that this “Large $`m_0`$” scenario is not unreasonably sensitive to small variations in the input top mass.

In all the plots of Figs. 37 the minimization scale varies in the range $`M_Z<Q_{\mathrm{min}}<2Q_{}`$, where $`Q_{}(m_0^2+4m_{1/2}^2)^{1/2}`$ is a scale roughly comparable with the squark masses. The dotted curves in the upper and lower panels of each figure represent $`|\mu (M_Z)|`$ and $`m_A(M_Z)`$, respectively, as obtained by including in the minimization conditions only the one–loop $`𝒪(\alpha _t)`$ top/stop contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$; the dashed lines include instead the full one–loop computation of $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$; the dot–dashed lines include in addition the two–loop $`𝒪(\alpha _t\alpha _s)`$ contributions; finally, the solid lines include our full two–loop result, i.e. the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ contributions to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$. In Fig. 5, corresponding to the “large $`\mathrm{tan}\beta `$” (SPS4) scenario with $`\mathrm{tan}\beta =50`$, we show also the effect of the $`𝒪(\alpha _b\alpha _s)`$ corrections, obtained from the $`𝒪(\alpha _t\alpha _s)`$ ones as described at the end of the previous section. The dot–dot–dashed lines in Fig. 5 include the $`𝒪(\alpha _t\alpha _s+\alpha _b\alpha _s)`$ contributions, while the solid lines represent the full $`𝒪(\alpha _t\alpha _s+\alpha _b\alpha _s+\alpha _t^2)`$ result. The effect of the $`𝒪(\alpha _b\alpha _s)`$ corrections is indeed negligible in the other SPS scenarios, where $`\mathrm{tan}\beta `$ takes on more moderate values. In any case we include the $`𝒪(\alpha _b\alpha _s)`$ corrections in all the scenarios we investigate here.

We see from Figs. 37 that the inclusion of the two–loop $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections improves the dependence of $`|\mu (M_Z)|`$ and $`m_A(M_Z)`$ on the minimization scale with respect to the full one–loop result. The effect is particularly manifest in the case of $`\mu `$, where, in all the scenarios, the two–loop corrected result appears to depend only very weakly (within 3 GeV at most) on $`Q_{\mathrm{min}}`$, while the one–loop result tends to decrease for both small and large $`Q_{\mathrm{min}}`$. In the case of $`m_A`$, the improvement is less striking: although the two–loop corrected result has in general a better scale dependence than its one–loop counterpart, especially for increasing values of $`Q_{\mathrm{min}}`$, a small residual scale dependence is visible in most plots when $`Q_{\mathrm{min}}`$ gets close to $`M_Z`$. In any scenario the residual uncertainty on $`m_A`$ is never larger than 10 GeV (the latter case occurring in the “large $`\mathrm{tan}\beta `$” scenario, where $`m_A(M_Z)`$ is around 1200 GeV) and might be due to the corrections that we neglect in our two–loop computation of the tadpoles, i.e. those controlled by the electroweak gauge couplings and, for large $`\mathrm{tan}\beta `$, those of $`𝒪(\alpha _b^2+\alpha _t\alpha _b)`$. From the small residual scale dependence visible in Figs. 37, we estimate that the effect of the neglected two–loop and higher–order corrections on the parameters $`\mu `$ and $`m_A`$ should be at most of 1%.

Other interesting observations can be drawn from Figs. 37: first of all, at the one–loop level, the inclusion of the top/stop contributions only is in general not a good approximation of the full result (this has been already observed e.g. in Ref. ). Moreover, at the two–loop level, a significant compensation occurs between the $`𝒪(\alpha _t\alpha _s)`$ and the $`𝒪(\alpha _t^2)`$ contributions to the tadpoles, thus including only the former may lead to rather inaccurate predictions (as shown by the dot–dashed curves). This partial cancellation between the $`𝒪(\alpha _t\alpha _s)`$ and $`𝒪(\alpha _t^2)`$ corrections is similar to the one occurring in the case of the Higgs masses . Finally, it is worth noticing that in the “large $`\mathrm{tan}\beta `$” scenario of Fig. 5 the inclusion of the $`𝒪(\alpha _b\alpha _s)`$ corrections to the tadpoles has a sizeable effect on the minimization scale dependence of $`m_A(M_Z)`$ (see the difference between the dot–dashed and dot–dot–dashed curves), while it does not affect significantly that of $`|\mu (M_Z)|`$.

It is clear from the above discussion that the numerical effect on $`|\mu (M_Z)|`$ and $`m_A(M_Z)`$ of the two–loop corrections to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$ depends critically on the choice of the minimization scale. In all the plots we find a range of values of $`Q_{\mathrm{min}}`$, usually in the vicinity of $`Q_{}`$, for which the one–loop and two–loop curves are close to each other, implying that the effect of the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections is small. On the other hand, for values of $`Q_{\mathrm{min}}`$ far from this optimal choice, the omission of the two–loop corrections can lead to an error of several (possibly, tenths of) GeV, especially in the case of $`\mu `$. Thus, the proper inclusion the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ two–loop tadpole corrections on the top of the full one–loop ones allows us to obtain more precise and reliable results for $`\mu `$ and $`m_A`$, i.e. results that do not depend on a preconceived choice of the minimization scale for the MSSM effective potential.

## 5 Conclusions

In this paper we presented explicit and general results for the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ two–loop corrections to the minimization conditions of the MSSM effective potential, which translate into corrections to the $`\overline{\mathrm{DR}}`$–renormalized parameters $`\mu `$ and $`m_A`$. We discussed the numerical impact of our corrections in some representative scenarios of gravity–mediated SUSY breaking, and we found that the inclusion of the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ corrections significantly improves the renormalization scale dependence of the results. Due to partial cancellations between the $`𝒪(\alpha _t\alpha _s)`$ and $`𝒪(\alpha _t^2)`$ corrections, including only the former may lead to inaccurate results. Our corrections are also required for consistency in the $`𝒪(\alpha _t\alpha _s+\alpha _t^2)`$ two–loop computations of the MSSM Higgs masses, if the parameters entering the tree–level Higgs mass matrix are computed via renormalization group evolution from a set of high energy boundary conditions.

A complete study of the electroweak symmetry breaking at the two–loop level would require also the knowledge of the corrections that are neglected in our gaugeless limit, among which the most relevant are controlled by the bottom Yukawa coupling and the electroweak gauge couplings. Concerning the corrections controlled by $`\alpha _b`$, they can be numerically non–negligible only for large values of $`\mathrm{tan}\beta `$. As discussed at the end of section 3, the formulae for the $`𝒪(\alpha _b\alpha _s)`$ corrections can be obtained by performing simple substitutions in their $`𝒪(\alpha _t\alpha _s)`$ counterparts. On the other hand, the $`𝒪(\alpha _b^2+\alpha _t\alpha _b)`$ corrections, which in some cases might be as relevant as the $`𝒪(\alpha _b\alpha _s)`$ ones, cannot be obtained in a straightforward way from the presently computed corrections and would require further work.

In Ref.  a complete two–loop computation of the MSSM effective potential is presented, including also the terms controlled by the electroweak gauge couplings that are neglected in our analysis. The inclusion of such terms improves further the scale dependence of the parameters $`\mu `$ and $`m_A`$, that in Ref.  are determined through a numerical minimization of the effective potential. However, the explicit analytical formulae for the two–loop tadpoles $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$, which are usually needed for practical applications, would be quite involved in the general case, and have not been presented so far.

In conclusion, our work should lead to a more precise and reliable determination of the MSSM parameters at the weak scale, once the boundary conditions are provided at some larger scale according to the underlying theory of SUSY breaking. The fact that, among its many attractive features, the MSSM provides a natural mechanism for breaking radiatively the electroweak symmetry, with the heavy top quark mass nicely falling in the required range777As an example, we find that in the “typical” mSUGRA (SPS1a) scenario the acceptable range for the top quark mass is $`80\mathrm{GeV}<M_t<215\mathrm{GeV}`$., seems to indicate the MSSM as the most viable theory for physics at the weak scale. However, only the forthcoming experimental results from the Tevatron and the LHC will tell us if this is indeed the case.

## Acknowledgments

We thank G. Degrassi, S. P. Martin, W. Porod, K. Tamvakis, F. Zwirner and especially M. Drees for useful comments and discussions. A. D. also thanks J. R. Espinosa and R. J. Zhang for discussions in the early stages of the project. This work was partially carried out at the Physics Department of the University of Bonn, and it was partially supported by the European Programmes HPRN-CT-2000-00148 (Across the Energy Frontier) and HPRN-CT-2000-00149 (Collider Physics).

## Appendix A: Two–loop integrals

We give here explicit expressions for the momentum integrals that appear in the two-loop part of the effective potential. The basic integrals in $`d=42ϵ`$ dimensions are:

$`{\displaystyle \frac{1}{(16\pi ^2)^2}}J(x,y)`$ $``$ $`{\displaystyle \frac{\mu ^{2(4d)}}{(2\pi )^{2d}}}{\displaystyle \frac{\mathrm{d}^dp\mathrm{d}^dq}{[p^2x][q^2y]}},`$ (A1)
$`{\displaystyle \frac{1}{(16\pi ^2)^2}}I(x,y,z)`$ $``$ $`{\displaystyle \frac{\mu ^{2(4d)}}{(2\pi )^{2d}}}{\displaystyle \frac{\mathrm{d}^dp\mathrm{d}^dq}{[p^2x][q^2y][(pq)^2z]}}.`$ (A2)
Following Ref. , the functions $`J(x,y)`$ and $`I(x,y,z)`$ defined in Eqs. (A1)–(A2) are:

$`J(x,y)`$ $`=`$ $`{\displaystyle \frac{xy}{ϵ^2}}{\displaystyle \frac{xy}{ϵ}}\left(\overline{\mathrm{ln}}x+\overline{\mathrm{ln}}y2\right)xy\left[2\overline{\mathrm{ln}}x+2\overline{\mathrm{ln}}y{\displaystyle \frac{1}{2}}\overline{\mathrm{ln}}^{\mathrm{\hspace{0.17em}2}}xy\left(3+{\displaystyle \frac{\pi ^2}{6}}\right)\right],`$ (A3)
$`I(x,y,z)`$ $`=`$ $`{\displaystyle \frac{x+y+z}{2ϵ^2}}+{\displaystyle \frac{1}{ϵ}}\left[x\overline{\mathrm{ln}}x+y\overline{\mathrm{ln}}y+z\overline{\mathrm{ln}}z{\displaystyle \frac{3}{2}}(x+y+z)\right]`$ (A4)
$`+{\displaystyle \frac{1}{2}}(x\overline{\mathrm{ln}}y\overline{\mathrm{ln}}z+y\overline{\mathrm{ln}}x\overline{\mathrm{ln}}z+z\overline{\mathrm{ln}}x\overline{\mathrm{ln}}y){\displaystyle \frac{x+y+z}{2}}(7+\pi ^2/6)`$
$`{\displaystyle \frac{1}{2}}(x\overline{\mathrm{ln}}x+y\overline{\mathrm{ln}}y+z\overline{\mathrm{ln}}z)(\overline{\mathrm{ln}}x+\overline{\mathrm{ln}}y+\overline{\mathrm{ln}}z6){\displaystyle \frac{\mathrm{\Delta }(x,y,z)}{2z}}\mathrm{\Phi }(x,y,z).`$
In the above formulae, $`\overline{\mathrm{ln}}x`$ stands for $`\mathrm{ln}(x/Q^2)`$, where $`Q^2=4\pi \mu ^2e^{\gamma _E}`$ is the renormalization scale ($`\gamma _E`$ is the Euler constant). The functions $`\mathrm{\Delta }`$ and $`\mathrm{\Phi }`$ are respectively:

$`\mathrm{\Delta }(x,y,z)`$ $`=`$ $`x^2+y^2+z^22(xy+xz+yz),`$ (A5)
$`\mathrm{\Phi }(x,y,z)`$ $`=`$ $`{\displaystyle \frac{1}{\lambda }}[2\mathrm{ln}x_+\mathrm{ln}x_{}\mathrm{ln}u\mathrm{ln}v2(\mathrm{Li}_2(x_+)+\mathrm{Li}_2(x_{}))+{\displaystyle \frac{\pi ^2}{3}}],`$ (A6)
where $`\mathrm{Li}_2(z)=_0^zdt\left[\mathrm{ln}(1t)/t\right]`$ is the dilogarithm function and the auxiliary (complex) variables are:

$$u=\frac{x}{z},v=\frac{y}{z},\lambda =\sqrt{(1uv)^24uv},x_\pm =\frac{1}{2}\left[1\pm (uv)\lambda \right].$$
(A7)
The definition (A6) is valid for the case $`x/z<1`$ and $`y/z<1`$. The other branches of $`\mathrm{\Phi }`$ can be obtained using the symmetry properties:

$$\mathrm{\Phi }(x,y,z)=\mathrm{\Phi }(y,x,z),x\mathrm{\Phi }(x,y,z)=z\mathrm{\Phi }(z,y,x).$$
(A8)
Finally, the following recursive relation for the derivatives of $`\mathrm{\Phi }`$ proved very useful 888We thank G. Degrassi for explanations on how to derive Eq. (A9). for obtaining compact analytical results:

$$\mathrm{\Delta }(x,y,z)\frac{\mathrm{\Phi }(x,y,z)}{x}=(y+zx)\mathrm{\Phi }(x,y,z)+\frac{z}{x}\left[(yz)\mathrm{ln}\frac{z}{y}+x\left(\mathrm{ln}\frac{x}{y}+\mathrm{ln}\frac{x}{z}\right)\right].$$
(A9)
The derivatives of $`\mathrm{\Phi }`$ with respect to $`y`$ and $`z`$ can be obtained from the above equation with the help of the symmetry properties of Eq. (A8).

## Appendix B: Explicit formulae for the $`𝒪(\alpha _t\alpha _s)`$ corrections

We present here explicit expressions for the two–loop part of the functions $`F`$ and $`G`$, giving rise to $`𝒪(\alpha _t\alpha _s)`$ corrections to $`\mathrm{\Sigma }_1`$ and $`\mathrm{\Sigma }_2`$. These formulae are valid when the parameters entering the one–loop parts of $`F`$ and $`G`$ are expressed in the $`\overline{\mathrm{DR}}`$ scheme and computed at the renormalization scale $`Q`$. In units of $`g_s^2C_FN_c/(16\pi ^2)^2`$, where $`C_F=4/3`$ and $`N_c=3`$ are color factors, $`F^2\mathrm{}`$ and $`G^2\mathrm{}`$ read:

$`F^2\mathrm{}`$ $`=`$ $`{\displaystyle \frac{4m_{\stackrel{~}{g}}m_t}{s_{2\theta }}}(1+4c_{2\theta }^2)\left[2(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)+{\displaystyle \frac{4m_{\stackrel{~}{g}}m_t}{s_{2\theta }}}\right]\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{g}}^2}{Q^2}}\mathrm{ln}{\displaystyle \frac{m_t^2}{Q^2}}`$ (B1)
$`2(4s_{2\theta }^2)(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)+{\displaystyle \frac{4m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2s_{2\theta }^2(m_{\stackrel{~}{t}_1}^2+m_{\stackrel{~}{t}_2}^2)^2}{m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_2}^2}{Q^2}}`$
$`+`$ $`\{[4(m_{\stackrel{~}{g}}^2+m_t^2+2m_{\stackrel{~}{t}_1}^2)s_{2\theta }^2(3m_{\stackrel{~}{t}_1}^2+m_{\stackrel{~}{t}_2}^2){\displaystyle \frac{16c_{2\theta }^2m_{\stackrel{~}{g}}m_tm_{\stackrel{~}{t}_1}^2}{s_{2\theta }(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)}}4s_{2\theta }m_{\stackrel{~}{g}}m_t]\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}`$
$`+{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2}}\left[s_{2\theta }^2(m_{\stackrel{~}{t}_1}^2+m_{\stackrel{~}{t}_2}^2)2(2m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)\right]\mathrm{ln}^2{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}`$
$`+2\left[m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{g}}^2m_t^2+m_{\stackrel{~}{g}}m_ts_{2\theta }+{\displaystyle \frac{2c_{2\theta }^2m_{\stackrel{~}{g}}m_tm_{\stackrel{~}{t}_1}^2}{s_{2\theta }(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)}}\right]\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{g}}^2m_t^2}{Q^4}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}`$
$`+{\displaystyle \frac{4m_{\stackrel{~}{g}}m_tc_{2\theta }^2(m_t^2m_{\stackrel{~}{g}}^2)}{s_{2\theta }(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)}}\mathrm{ln}{\displaystyle \frac{m_t^2}{m_{\stackrel{~}{g}}^2}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}+[{\displaystyle \frac{4m_{\stackrel{~}{g}}^2m_t^2+2\mathrm{\Delta }(m_{\stackrel{~}{g}}^2,m_t^2,m_{\stackrel{~}{t}_1}^2)}{m_{\stackrel{~}{t}_1}^2}}`$
$`{\displaystyle \frac{2m_{\stackrel{~}{g}}m_ts_{2\theta }}{m_{\stackrel{~}{t}_1}^2}}(m_{\stackrel{~}{g}}^2+m_t^2m_{\stackrel{~}{t}_1}^2)+{\displaystyle \frac{4c_{2\theta }^2m_{\stackrel{~}{g}}m_t\mathrm{\Delta }(m_{\stackrel{~}{g}}^2,m_t^2,m_{\stackrel{~}{t}_1}^2)}{s_{2\theta }m_{\stackrel{~}{t}_1}^2(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)}}]\mathrm{\Phi }(m_{\stackrel{~}{g}}^2,m_t^2,m_{\stackrel{~}{t}_1}^2)`$
$`(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2,s_{2\theta }s_{2\theta })\},`$
$`G^2\mathrm{}`$ $`=`$ $`{\displaystyle \frac{5m_{\stackrel{~}{g}}s_{2\theta }}{m_t}}(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)10(m_{\stackrel{~}{t}_1}^2+m_{\stackrel{~}{t}_2}^22m_t^2)4m_{\stackrel{~}{g}}^2+12m_t^2\left(\mathrm{ln}^2{\displaystyle \frac{m_t^2}{Q^2}}2\mathrm{ln}{\displaystyle \frac{m_t^2}{Q^2}}\right)`$ (B2)
$`+\left[4m_{\stackrel{~}{g}}^2{\displaystyle \frac{m_{\stackrel{~}{g}}s_{2\theta }}{m_t}}(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2)\right]\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{g}}^2}{Q^2}}\mathrm{ln}{\displaystyle \frac{m_t^2}{Q^2}}+s_{2\theta }^2(m_{\stackrel{~}{t}_1}^2+m_{\stackrel{~}{t}_2}^2)\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_2}^2}{Q^2}}`$
$`+`$ $`\{[4(m_{\stackrel{~}{g}}^2+m_t^2+2m_{\stackrel{~}{t}_1}^2)+s_{2\theta }^2(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2){\displaystyle \frac{4m_{\stackrel{~}{g}}s_{2\theta }}{m_t}}(m_t^2+m_{\stackrel{~}{t}_1}^2)]\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}`$
$`+\left[{\displaystyle \frac{m_{\stackrel{~}{g}}s_{2\theta }}{m_t}}(5m_t^2m_{\stackrel{~}{g}}^2+m_{\stackrel{~}{t}_1}^2)2(m_{\stackrel{~}{g}}^2+2m_t^2)\right]\mathrm{ln}{\displaystyle \frac{m_t^2}{Q^2}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}`$
$`+\left[{\displaystyle \frac{m_{\stackrel{~}{g}}s_{2\theta }}{m_t}}(m_{\stackrel{~}{g}}^2m_t^2+m_{\stackrel{~}{t}_1}^2)2m_{\stackrel{~}{g}}^2\right]\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{g}}^2}{Q^2}}\mathrm{ln}{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}(2+s_{2\theta }^2)m_{\stackrel{~}{t}_1}^2\mathrm{ln}^2{\displaystyle \frac{m_{\stackrel{~}{t}_1}^2}{Q^2}}`$
$`+\left[2{\displaystyle \frac{m_{\stackrel{~}{g}}^2}{m_{\stackrel{~}{t}_1}^2}}(m_{\stackrel{~}{g}}^2+m_t^2m_{\stackrel{~}{t}_1}^22m_{\stackrel{~}{g}}m_ts_{2\theta })+{\displaystyle \frac{m_{\stackrel{~}{g}}s_{2\theta }}{m_tm_{\stackrel{~}{t}_1}^2}}\mathrm{\Delta }(m_{\stackrel{~}{g}}^2,m_t^2,m_{\stackrel{~}{t}_1}^2)\right]\mathrm{\Phi }(m_{\stackrel{~}{g}}^2,m_t^2,m_{\stackrel{~}{t}_1}^2)`$
$`+(m_{\stackrel{~}{t}_1}^2m_{\stackrel{~}{t}_2}^2,s_{2\theta }s_{2\theta })\}.`$
The functions $`\mathrm{\Delta }(x,y,z)`$ and $`\mathrm{\Phi }(x,y,z)`$ appearing in $`F^2\mathrm{}`$ and $`G^2\mathrm{}`$ are defined in Eqs. (A5)–(A6). The parameter $`s_{2\theta }`$ is defined in Eq. (13).


format: text

text:
# On-Line Learning Theory of Soft Committee Machines with Correlated Hidden Units – Steepest Gradient Descent and Natural Gradient Descent –

## I Introduction

One of the biggest problems of neural network learning is the plateau of the learning curve. Considering the gradient learning method and its generalization error, this plateau is mainly caused by the saddle structure of the error function. The permutation symmetry prevents the identification of the hidden units in multilayer perceptrons if they have the same weight vectors, and produces this saddle structure rf:1 ; rf:2 . In the learning scenario of a teacher and a student network, the saddle is thought to be affected by the strength of the correlation of the hidden units in the teacher network, which may be closely related to the length of the plateau. More specifically, in the conventional gradient descent (GD), the weight vectors in the student network are known to approach the saddle before reaching their final states rf:2 . Since the saddle is located between the weight vectors of the teacher hidden units, their stronger correlation is supposed to force the student weight vectors closer to the saddle, resulting in a longer plateau.

Natural gradient descent (NGD), however, may be able to avoid the saddle because it can update the network parameters to the optimal direction in the Riemannian space rf:3 . NGD is a fairly general method for effectively adjusting the parameters of stochastic models, but its validity in multilayer perceptrons is uncertain because of three intrinsic problems: 1) NGD needs prior knowledge of the input distribution to calculate the Fisher information matrix, 2) NGD is unstable around the singular points of the Fisher information matrix, 3) matrix inversion is time consuming, which might be critical especially in real-time learning. The method proposed by Yang and Amarirf:4  can be used to calculate NGD efficiently in the case of a large input dimension in multilayer perceptrons. Also, the adaptive method can be used to approximate the inverse of the Fisher information matrix asymptotically without prior knowledge or matrix inversion rf:5 . In this paper, we discuss the problem of singularity; since the saddle is one of the singular points, how NGD works around there is one of our main topics.

On-line learning is one of the most popular forms of training. Analysis of the network dynamics in on-line learning is much easier than for batch learning because the state of the network and the learning samples are independent of each other. In this framework, the statistical mechanics method proposed by Saad and Solla can be used to analyze the GD dynamics exactly at the large limit of the input dimension rf:2 . Rattray and Saad extended this technique to NGD and reported that it works efficiently in multilayer perceptrons rf:6 . In this paper, we also use this method and contrast the dynamics for GD and NGD, focusing on the corrupted saddle structure under a strong correlation of the hidden units in the teacher network.

## II Model

Soft committee machines (Fig. 1) are considered where the teacher network has $`M`$ hidden units while the student has $`K`$ units. To apply NGD, Gaussian noise $`n𝒩(0,\sigma ^2)`$ is added to the output of the student;

$`\zeta `$ $``$ $`f{}_{𝑩}{}^{}(\xi ),f{}_{𝑩}{}^{}(𝝃){\displaystyle \underset{k=1}{\overset{M}{}}}g(𝑩_k^T𝝃),`$ (1)
$`\zeta ^{}`$ $``$ $`f{}_{𝑱}{}^{}(\xi )+n,f{}_{𝑱}{}^{}(𝝃){\displaystyle \underset{k=1}{\overset{K}{}}}g(𝑱_k^T𝝃),`$ (2)
where $`𝝃^N`$ denotes the input vector while $`𝑩_i^N`$ and $`𝑱_i^N`$ are the $`i`$th weight vectors of the teacher and the student networks, respectively. Here, $`T`$ means the transposition while $`g`$ is an activation function.

The joint probability distribution of the input $`𝝃`$ and the output $`\zeta ^{}`$ of the student network is given by

$`p{}_{𝑱}{}^{}(𝝃,\zeta ^{})p(𝝃)p{}_{𝑱}{}^{}(\zeta ^{}|𝝃),`$ (3)
$`p{}_{𝑱}{}^{}(\zeta ^{}|𝝃){\displaystyle \frac{1}{\sqrt{2\pi \sigma ^2}}}\mathrm{exp}({\displaystyle \frac{\{\zeta ^{}f{}_{𝑱}{}^{}(𝝃)\}^2}{2\sigma ^2}}).`$ (4)
The parameter vector of (3), $`𝑱[𝑱_1^T,𝑱_2^T,\mathrm{},𝑱_K^T]^T^{KN}`$ , is updated iteratively to approximate the joint probability distribution of the input $`𝝃`$ and the output $`\zeta `$ of the teacher network,

$`p(𝝃,\zeta )p(𝝃)\delta (\zeta f{}_{𝑩}{}^{}(𝝃)),`$ (5)
where $`\delta `$ is the delta function. The loss function for a given set of a learning sample $`\{𝝃,\zeta \}`$, defined using the logarithmic loss of the conditional probability distribution of (4), is

$`ϵ{}_{𝑱}{}^{}(𝝃,\zeta )\mathrm{ln}p{}_{𝑱}{}^{}(\zeta |𝝃)+c_0={\displaystyle \frac{1}{2\sigma ^2}}\{\zeta f{}_{𝑱}{}^{}(\xi )\}^2,`$ (6)
where $`c_0\mathrm{ln}\sqrt{2\pi \sigma ^2}`$ is constant. The generalization error is then defined as the expected loss:

$`ϵ_g(𝑱)ϵ{}_{𝑱}{}^{}(𝝃,\zeta )_\{{}_{𝝃}{}^{}{}_{,\zeta \}}{}^{}.`$ (7)
The definitions of (6) can be written, by applying (1) and (2), as

$`ϵ{}_{𝑱}{}^{}(𝝃,\zeta )=ϵ{}_{𝑱}{}^{}(𝝃){\displaystyle \frac{1}{2\sigma ^2}}\{f{}_{𝑩}{}^{}(𝝃)f{}_{𝑱}{}^{}(𝝃)\}^2.`$ (8)
We consider on-line learning in this paper, where the parameter vector $`𝑱`$ is updated for each set of an independently given sample $`\{𝝃,\zeta \}`$. The updating rule, the differential of $`𝑱`$, for GD is defined with a learning rate $`\eta `$ as

$`\mathrm{\Delta }𝑱={\displaystyle \frac{\eta }{N}}{}_{𝑱}{}^{}ϵ{}_{𝑱}{}^{}(𝝃,\zeta ),`$ (9)
where

$`{}_{𝑱}{}^{}{}_{i}{}^{}ϵ{}_{𝑱}{}^{}(𝝃,\zeta )={\displaystyle \frac{1}{\sigma ^2}}g^{}(𝑱_i^T𝝃)\{f{}_{𝑩}{}^{}(𝝃)f{}_{𝑱}{}^{}(𝝃)\}𝝃,`$ (10)
where $`g^{}`$ denotes the derivative of $`g`$. One for NGD is also defined as

$`\mathrm{\Delta }𝑱={\displaystyle \frac{\eta }{N}}𝑮^1{}_{𝑱}{}^{}ϵ{}_{𝑱}{}^{}(𝝃,\zeta ),`$ (11)
where $`𝑮`$ denotes the Fisher information matrix of the parameter vector $`𝑱`$:

$`𝑮[{}_{𝑱}{}^{}\mathrm{ln}p{}_{𝑱}{}^{}(𝝃,\zeta ^{})][{}_{𝑱}{}^{}\mathrm{ln}p{}_{𝑱}{}^{}(𝝃,\zeta ^{})]^T_\{{}_{𝝃}{}^{}{}_{,\zeta ^{}\}}{}^{}.`$ (12)
The $`𝑮`$ can be written, in block form, as

$`𝑮=[\begin{array}{ccc}𝑮_{1,1}& \mathrm{}& 𝑮_{1,K}\\ \mathrm{}& \mathrm{}& \mathrm{}\\ 𝑮_{K,1}& \mathrm{}& 𝑮_{K,K}\end{array}],`$ (16)
$`𝑮_{ij}={\displaystyle \frac{1}{\sigma ^2}}g^{}(𝑱_i^T𝝃)g^{}(𝑱_j^T𝝃)𝝃𝝃^T_\{{}_{𝝃}{}^{}{}_{\}}{}^{}.`$ (17)
In the case of the standard multivariate normal distribution input, $`𝝃𝒩(\mathrm{𝟎},𝑰)`$, the inverse of the Fisher information matrix is also given by

$`𝑮^1=[\begin{array}{ccc}𝑮_{1,1}^1& \mathrm{}& 𝑮_{1,K}^1\\ \mathrm{}& \mathrm{}& \mathrm{}\\ 𝑮_{K,1}^1& \mathrm{}& 𝑮_{K,K}^1\end{array}],`$ (21)
$`𝑮_{ij}^1=\sigma ^2\{\theta _{ij}𝑰+𝑱^{}𝚯_{ij}𝑱^{}{}_{}{}^{T}\},`$ (22)
where $`𝑱^{}[𝑱_1,\mathrm{},𝑱_K]`$ is a $`N`$ by $`K`$ matrix, while $`\theta _{ij}`$ is a scalar and $`𝚯_{ij}`$ is a $`K`$ by $`K`$ matrix rf:4 .

## III Theory

### III.1 Order parameters and generalization error

At the thermodynamics limit, the limit of $`N\mathrm{}`$, the dynamics of the network can be analyzed using statistical mechanics. Here, the order parameters that represent the correlations of the weight vectors are used instead of the $`N`$-dimensional vectors $`𝝃`$, $`𝑩_i`$, and $`𝑱_i`$. To make the present paper self-contained, we briefly summarize the derivation of the order parameter equations of the soft committee machine rf:2 ; rf:6 .

From here on, the input vector is assumed to obey a $`N`$-dimensional multivariate Gaussian noise with zero mean and a unit covariance matrix: $`𝝃𝒩(\mathrm{𝟎},𝑰)`$. The correlation between the input and each weight vector, denoted by $`x_i𝑱_i^T𝝃`$ and $`y_i𝑩_i^T𝝃`$, is then distributed as a normal distribution; $`x_i𝒩(0,𝑱_i^T𝑱_i)`$ and $`y_i𝒩(0,𝑩_i^T𝑩_i)`$, while each covariance of them is given by $`x_ix_j_\{{}_{𝝃}{}^{}{}_{\}}{}^{}=𝑱_i^T𝑱_j`$, $`x_iy_j_\{{}_{𝝃}{}^{}{}_{\}}{}^{}=𝑱_i^T𝑩_j`$, and $`y_iy_j_\{{}_{𝝃}{}^{}{}_{\}}{}^{}=𝑩_i^T𝑩_j`$. Therefore, a new vector, defined as

$`𝒛[x_1,\mathrm{},x_K,y_1,\mathrm{},y_M]^T^{K+M},`$ (23)
is distributed as a multivariate normal distribution $`𝒩(\mathrm{𝟎},𝑪)`$:

$`p(𝒛)={\displaystyle \frac{1}{\sqrt{\{2\pi \}^{K+M}|𝑪|}}}\mathrm{exp}\left({\displaystyle \frac{1}{2}}𝒛^T𝑪^1𝒛\right),`$ (24)
where $`𝑪`$ is the variance-covariance matrix:

$`𝑪[\begin{array}{cc}𝑸& 𝑹\\ 𝑹^T& 𝑻\end{array}]`$ (27)
with

$`𝑸𝑱^{}{}_{}{}^{T}𝑱_{}^{}`$ $`=`$ $`[\begin{array}{ccc}Q_{1,1}& \mathrm{}& Q_{1,K}\\ \mathrm{}& \mathrm{}& \mathrm{}\\ Q_{K,1}& \mathrm{}& Q_{K,K}\end{array}],`$ (31)
$`𝑹𝑱^{}{}_{}{}^{T}𝑩_{}^{}`$ $`=`$ $`[\begin{array}{ccc}R_{1,1}& \mathrm{}& R_{1,M}\\ \mathrm{}& \mathrm{}& \mathrm{}\\ R_{K,1}& \mathrm{}& R_{K,M}\end{array}],`$ (35)
$`𝑻𝑩^{}{}_{}{}^{T}𝑩_{}^{}`$ $`=`$ $`[\begin{array}{ccc}T_{1,1}& \mathrm{}& T_{1,M}\\ \mathrm{}& \mathrm{}& \mathrm{}\\ T_{M,1}& \mathrm{}& T_{M,M}\end{array}],`$ (39)
and $`𝑱^{}[𝑱_1\mathrm{}𝑱_K]`$, $`𝑩^{}[𝑩_1\mathrm{}𝑩_M]`$. Here, $`𝑸`$ and $`𝑹`$ are the order parameters of this system.

Using these order parameters, the generalization error in (7), $`ϵ_g(𝑱)ϵ{}_{𝑱}{}^{}(𝝃,\zeta )_\{_{,\zeta \}}^{}{}_{𝝃}{}^{}`$, can be calculated by

$`ϵ_g(𝑱)`$ $`=`$ $`{\displaystyle 𝑑𝒛p(𝒛)\frac{1}{2\sigma ^2}\left\{\underset{k=1}{\overset{M}{}}g(y_k)\underset{k=1}{\overset{K}{}}g(x_k)\right\}^2}.`$ (40)
If we define the activation function $`g`$ as $`g(x)=\mathrm{erf}(x/\sqrt{2})`$ from here on, the generalization error is given by

$`ϵ_g(𝑱)`$ $`=`$ $`{\displaystyle \frac{1}{\pi \sigma ^2}}\{2{\displaystyle \underset{i=1}{\overset{K}{}}}{\displaystyle \underset{j=1}{\overset{M}{}}}\mathrm{arcsin}{\displaystyle \frac{R_{ij}}{\sqrt{Q_{ii}+1}\sqrt{T_{jj}+1}}}`$ (41)
$`+`$ $`{\displaystyle \underset{i,j=1}{\overset{K}{}}}\mathrm{arcsin}{\displaystyle \frac{Q_{ij}}{\sqrt{Q_{ii}+1}\sqrt{Q_{jj}+1}}}`$
$`+`$ $`{\displaystyle \underset{i,j=1}{\overset{M}{}}}\mathrm{arcsin}{\displaystyle \frac{T_{ij}}{\sqrt{T_{ii}+1}\sqrt{T_{jj}+1}}}\},`$
which depends on only the order parameters.

### III.2 Dynamics of the order parameters

Here we substitute the dynamics of the order parameters for those of the system. First, we can replace the updating rule (9) with

$`\mathrm{\Delta }𝑱_i={\displaystyle \frac{\eta }{\sigma ^2N}}\delta _i𝝃,`$ (42)
where

$`\delta _ig^{}(x_i)\left\{{\displaystyle \underset{k=1}{\overset{M}{}}}g(y_k){\displaystyle \underset{k=1}{\overset{K}{}}}g(x_k)\right\}.`$ (43)
Thus, the updating rule of the order parameters is given by

$`\mathrm{\Delta }R_{ij}`$ $`=`$ $`[𝑱_i+\mathrm{\Delta }𝑱_i]^T𝑩_j𝑱_i^T𝑩_j`$ (44)
$`=`$ $`{\displaystyle \frac{\eta }{\sigma ^2N}}\delta _iy_j,`$
and

$`\mathrm{\Delta }Q_{ij}`$ $`=`$ $`[𝑱_i+\mathrm{\Delta }𝑱_i]^T[𝑱_j+\mathrm{\Delta }𝑱_j]𝑱_i^T𝑱_j`$ (45)
$`=`$ $`{\displaystyle \frac{\eta }{\sigma ^2N}}\{\delta _ix_j+\delta _jx_i\}+{\displaystyle \frac{\eta ^2}{\sigma ^4N^2}}\delta _i\delta _j𝝃^T𝝃.`$
Here we introduce the time $`\alpha `$; a short period, $`\mathrm{\Delta }\alpha =1/N`$, is defined to be consumed for each learning iteration. At the large limit of $`N`$, the differential dynamics of $`R_{ij}`$ and $`Q_{ij}`$ are calculated as

$`{\displaystyle \frac{R_{ij}}{\alpha }}`$ $`=`$ $`\underset{\mathrm{\Delta }\alpha 0}{lim}{\displaystyle \frac{\mathrm{\Delta }R_{ij}}{\mathrm{\Delta }\alpha }}=\underset{N\mathrm{}}{lim}N\mathrm{\Delta }R_{ij}={\displaystyle \frac{\eta }{\sigma ^2}}\delta _iy_j_\{_{\}}^{}{}_{𝒛}{}^{}`$ (46)
$`=`$ $`{\displaystyle \frac{\eta }{\sigma ^2}}\psi _{ij}`$
and

$`{\displaystyle \frac{Q_{ij}}{\alpha }}={\displaystyle \frac{\eta }{\sigma ^2}}\{\varphi _{ij}+\varphi _{ji}\}+{\displaystyle \frac{\eta ^2}{\sigma ^4}}\upsilon _{ij},`$ (47)
where $`𝝃^T𝝃N`$ is applied, while the new variables are defined as

$`\psi _{ij}\delta _iy_j_\{{}_{𝒛}{}^{}{}_{\}}{}^{},\varphi _{ij}\delta _ix_j_\{{}_{𝒛}{}^{}{}_{\}}{}^{},\upsilon _{ij}\delta _i\delta _j_\{{}_{𝒛}{}^{}{}_{\}}{}^{},`$ (48)
The dynamics for NGD can be provided in the same way:

$`\mathrm{\Delta }𝑱_i={\displaystyle \frac{\eta }{\sigma ^2N}}{\displaystyle \underset{k=1}{\overset{K}{}}}\delta _k𝑮_{ik}^1𝝃.`$ (49)
Thus, the dynamics of the order parameters are

$`{\displaystyle \frac{R_{ij}}{\alpha }}=\eta {\displaystyle \underset{k=1}{\overset{K}{}}}\{\theta _{ik}\psi _{kj}+\varphi _k𝚯_{ki}R_j\},`$ (50)
$`{\displaystyle \frac{Q_{ij}}{\alpha }}`$ $`=`$ $`\eta {\displaystyle \underset{k=1}{\overset{K}{}}}\{\theta _{ik}\varphi _{kj}+\theta _{jk}\varphi _{ki}+\varphi _k𝚯_{ik}^TQ_j`$ (51)
$`+\varphi _k𝚯_{jk}^TQ_i\}+\eta ^2{\displaystyle }_{k,l=1}^K\theta _{ik}\theta _{jl}\upsilon _{kl},`$
where $`\varphi _k`$ denotes the $`k`$th row of the matrix $`\{\varphi _{ij}\}_{i,j=1,\mathrm{},K}`$, while $`R_j`$ denotes the $`j`$th column of the matrix $`𝑹`$, and so on rf:6 .

## IV Numerical Results

In this section, we discuss how the learning dynamics depend on the correlation of teacher weight vectors $`T_{1,2}`$. The results are also contrasted between GD and NGD.

We set the number of the hidden units and the lengths of the teacher weight vectors as follows:

$`K=M=2,T_{1,1}=T_{2,2}.`$ (52)
We also restrict the initial conditions to

$`Q_{1,1}=Q_{2,2},R_{1,1}=R_{2,2},R_{1,2}=R_{2,1}.`$ (53)
Because of the symmetry of the system, these restrictions are preserved throughout the learning. Specifically, we use

$`𝑸=[\begin{array}{cc}1& 0\\ 0& 1\end{array}],𝑹=[\begin{array}{cc}10^2& 0\\ 0& 10^2\end{array}].`$ (58)
Therefore, we have four free parameters $`Q_{1,1}`$, $`Q_{1,2}`$, $`R_{1,1}`$, and $`R_{1,2}`$ in this system. Note that $`𝑸`$ and $`𝑻`$ are always symmetric matrices from the definitions of (31) and (39). Other parameters are set as $`T_{1,2}=1`$, $`\eta =10^2`$, and $`\sigma ^2=5\times 10^2`$. Various values for $`T_{1,2}`$ are employed to examine the influence of the correlation of the teacher hidden units. We sometimes use $`\kappa \mathrm{arccos}\frac{T_{1,2}}{T_{1,1}}`$, the angle of the teacher weight vectors, instead of $`T_{1,2}`$.

In this case, $`\theta _{ij}`$ and $`𝚯_{ij}`$ in the inverse of the Fisher information matrix (21) can be simplified as

$`\theta _{1,1}`$ $`=`$ $`\theta _{2,2}=c\sqrt{a},\theta _{1,2}=\theta _{2,1}=c\sqrt{b},`$
$`𝚯_{1,1}`$ $`=`$ $`d\sqrt{a}[\begin{array}{cc}2a\{ab\}bQ_{1,1}& bQ_{1,2}\\ bQ_{1,2}& bQ_{1,1}\end{array}],`$ (61)
$`𝚯_{2,2}`$ $`=`$ $`d\sqrt{a}[\begin{array}{cc}bQ_{1,1}& bQ_{1,2}\\ bQ_{1,2}& 2a\{ab\}bQ_{1,1}\end{array}],`$ (64)
$`𝚯_{1,2}`$ $`=`$ $`𝚯_{2,1}`$ (67)
$`=`$ $`d\sqrt{b}[\begin{array}{cc}a\{Q_{1,1}+1\}b^2& aQ_{1,2}\\ aQ_{1,2}& a\{Q_{1,1}+1\}b^2\end{array}],`$
where

$`a\{Q_{1,1}+1\}^2Q_{1,2}^2,b2Q_{1,1}+1,`$
$`c{\displaystyle \frac{\pi }{2}}{\displaystyle \frac{\sqrt{ab}}{ab}},d{\displaystyle \frac{c}{a^2b^2}}.`$ (68)
Here we summarize the order of each variable to $`N`$. Since the length of the input vector $`𝝃`$ is $`O(\sqrt{N})`$, $`x_i`$ and $`y_i`$ are $`O(1)`$. This guarantees that the arguments of the activation function $`g`$ are $`O(1)`$. Therefore, the lengths of the weight vectors, $`\sqrt{Q_{ii}}`$ and $`\sqrt{T_{ii}}`$, are $`O(1)`$. If the direction of the initial $`𝑱_i`$ is chosen randomly, the size of $`R_{ii}`$, the correlation between $`𝑱_i`$ and $`𝑩_i`$, is $`O(1/\sqrt{N})`$. The initial numerical values in (58) are defined according to these sizes.

Figure 3 shows the time evolution of the generalization error. In the GD (Fig. 3a), the plateau was greatly prolonged as the correlation of the teacher weight vectors rose. In NGD (Fig. 3b), almost no plateau occurred at any $`T_{1,2}`$ if $`\eta `$ was set small enough relative to the initial $`R_{1,1}`$, and the generalization error was exponentially decreased. The plateau periods of Fig. 1a were measured and are shown in Fig. 3c, where we defined a plateau as occurring if $`\frac{\mathrm{ln}ϵ_g}{\alpha }>0.0005`$. The order of the plateau lengths was about $`O(\kappa ^3)`$ in GD.

Figure 3 shows the trajectories of the order parameters $`(R_{1,1},R_{1,2})`$ and $`(R_{2,1},R_{2,2})`$. Because of the symmetry, the latter plots are mirror images of the former. As $`R_{1,1}`$ is the correlation between the first student and the corresponding teacher, the initial value is almost $`0`$ and the goal is $`1`$; $`R_{1,2}`$ is the correlation between the first student and the not corresponding teacher, and the initial value is almost $`0`$ and the goal is $`T_{1,2}`$. Therefore, the target location of the plots are $`(1,T_{1,2})`$ and $`(T_{1,2},1)`$, respectively (shown as $`\mathrm{}`$). The other order parameters $`Q_{1,1}`$ and $`Q_{1,2}`$ are not shown. In the case of GD (Fig. 3a), the plots start at $`\mathrm{}`$, turn back at $`\mathrm{}`$, then approach $``$ (the saddle, as explained in the next section), and finally reach $`\mathrm{}`$. Actually, the parameters never pass through the same place again because $`Q_{1,1}`$ and $`Q_{1,2}`$ are updated. In the case of NGD (Fig. 3b), the plots start at $`\mathrm{}`$ and reach $`\mathrm{}`$ while avoiding $``$.

We performed a numerical simulation to confirm the dynamics at the above thermodynamics limit. The input dimension was $`N=1000`$, the teacher weight vectors were set as

$`𝑩_1=[\begin{array}{c}1\\ 0\\ 0\\ \mathrm{}\\ 0\end{array}],𝑩_2=[\begin{array}{c}\mathrm{cos}\kappa \\ \mathrm{sin}\kappa \\ 0\\ \mathrm{}\\ 0\end{array}],`$ (79)
and every initial $`𝑱_i`$ was randomly and independently chosen from $`𝒩(\mathrm{𝟎},𝑰/N)`$ for each try. Thus, the order parameters $`𝑸`$ and $`𝑹`$ were no longer limited by the restriction of (53). The learning was performed using these real weight vectors and the original equations: (9) for GD and (11) for NGD. Figures 5 and 5 show the time evolution of the generalization error and the trajectories of the order parameters in the same manner as Figs. 3 and 3, respectively. Both figures support the statistical dynamics well, which suggests the constraint of (53) is a rather minor problem and the system retains most of its generality even with that restriction.

## V Saddle

Here, we discuss why NGD is so effective even with a strong correlation between teacher hidden units. We consider the dynamics around the saddle of the generalization error under the conditions of (52) and (53). This point, where all the differentials of the order parameters are zero and the Hessian matrix is not positive definite nor negative definite, is shown as $``$ in Figs. 3 and 5:

$`Q_{1,1}=Q_{1,2}`$ $`=`$ $`{\displaystyle \frac{T_{1,1}+T_{1,2}}{T_{1,1}T_{1,2}+2}},`$
$`R_{1,1}=R_{1,2}`$ $`=`$ $`{\displaystyle \frac{T_{1,1}+T_{1,2}}{\sqrt{2\{T_{1,1}T_{1,2}+2\}}}}.`$ (80)
This saddle is a special point because 1) it corresponds to the goal both in the case of $`T_{1,1}=T_{1,2}`$ (the teacher is a smaller network: $`f{}_{𝑩}{}^{}(𝝃)=2g(𝑩_1^T𝝃)`$) and in the case that the student is a smaller network: $`f{}_{𝑱}{}^{}(𝝃)=2g(𝑱_1^T𝝃)`$, 2) in GD, the plateau occurs around it, and in NGD the student vectors avoid it, 3) it coincides with one of the singular points of the Fisher information matrix since $`Q_{1,1}=Q_{1,2}`$.

We simplify the situation as shown in Fig. 6; the two student weight vectors belong to the plane made by the two teacher weight vectors. This simplification is useful because we are now interested in how fast the student vectors leave this point for the goals. The correlations are re-parameterized by $`\kappa `$ and $`\lambda `$ as

$`T_{1,2}=T_{1,1}\mathrm{cos}\kappa ,Q_{1,2}=Q_{1,1}\mathrm{cos}\lambda ,`$
$`R_{1,1}=\sqrt{Q_{1,1}T_{1,1}}\mathrm{cos}{\displaystyle \frac{\kappa \lambda }{2}},`$
$`R_{1,2}=\sqrt{Q_{1,1}T_{1,1}}\mathrm{cos}{\displaystyle \frac{\kappa +\lambda }{2}}.`$ (81)
Now, we have only two free parameters $`Q_{1,1}`$ and $`\lambda `$. Since the first derivative of $`\lambda `$ can be written with $`Q_{1,1}`$ and $`Q_{1,2}`$ as

$`{\displaystyle \frac{\lambda }{\alpha }}={\displaystyle \frac{}{\alpha }}\mathrm{arccos}{\displaystyle \frac{Q_{1,2}}{Q_{1,1}}}={\displaystyle \frac{Q_{1,2}\frac{Q_{1,1}}{\alpha }Q_{1,1}\frac{Q_{1,2}}{\alpha }}{Q_{1,1}\sqrt{Q_{1,1}^2Q_{1,2}^2}}},`$ (82)
we can formulate the angular velocity of $`\lambda `$ at $`0<\lambda 1`$. The term $`\eta ^2`$ included in $`\frac{Q_{ij}}{\alpha }`$ can be ignored if the learning rate $`\eta `$ is set small enough.

The angular velocity for GD is

$`{\displaystyle \frac{\lambda }{\alpha }}=c_1\lambda \mathrm{sin}^2\kappa ,`$ (83)
where $`c_1\frac{4\eta }{\pi \sigma ^2}T_{1,1}\{T_{1,1}\{1\mathrm{cos}\kappa \}+2\}^{\frac{1}{2}}\{T_{1,1}\{3+\mathrm{cos}\kappa \}+2\}^{\frac{3}{2}}`$. We notice that the order of $`c_1`$ is not greatly changed by $`\kappa `$. The velocity converges to zero in the first order of $`\lambda `$. Moreover, it decreases as $`\kappa `$ decreases. Therefore, this equation supports the simulation results showing that the plateau is prolonged as the teacher correlation rises. The angular velocity for NGD is

$`{\displaystyle \frac{\lambda }{\alpha }}=c_2{\displaystyle \frac{1}{\lambda }}\mathrm{tan}^2{\displaystyle \frac{\kappa }{2}},`$ (84)
where $`c_22\eta `$. This velocity diverges to infinity as $`\lambda `$ goes to zero. Although it decreases as $`\kappa `$ decreases, this effect would be canceled by $`\lambda ^1`$ near the saddle. Therefore, this equation means that the student weight vectors are repelled by the saddle. In addition, this also supports the simulation results showing that the student weight vectors avoid the saddle and that the plateau does not occur even in the case of strongly correlated teacher hidden units.

## VI Conclusion

We have studied the on-line learning of soft committee machines under correlated teacher hidden units. The plateau in GD is largely prolonged at about $`O(\kappa ^3)`$ as the correlation of the teacher weight vectors rises, but almost no plateau occurs in NGD with a low learning rate $`\eta `$ and this does not depend on the correlation. Our analytical results for around the saddle reveal that the NGD avoided the saddle, even though the strong correlation of the teacher weight vectors forced the student weight vectors close to the saddle where the Fisher information matrix is singular.