Precision physics with QCD

The four-loop determination of the strong coupling from fully inclusive observables is reviewed. Special attention is given to the low-energy measurement extracted from the hadronic $\tau$ decay width. A recent exhaustive analysis of the ALEPH data, exploring several complementary methodologies with very different sensitivities to inverse power corrections and duality violations, confirms the strong suppression of non-perturbative contributions to $R_\tau$. It gives the value $\alpha_s(m_\tau^2)= 0.328 \pm 0.013$, which implies $\alpha_s(M_Z^2)= 0.1197 \pm 0.0015$. The excellent agreement with the direct measurement at the $Z$ peak, $\alpha_s(M_Z^2)= 0.1196 \pm 0.0030$, provides a beautiful test of asymptotic freedom. Together with the most recent lattice average from FLAG and the NNLO determinations from $e^+e^-$, PDFs and collider data quoted by the PDG, these two inclusive determinations imply a world average value $\alpha_s(M_Z^2)= 0.1180 \pm 0.0010$.


Introduction
All strong interaction phenomena should be described in terms of the strong coupling α s , the single free parameter of Quantum Chromodynamics (QCD). The overwhelming consistency of the many determinations of α s , performed in different processes and at different mass scales provides a beautiful verification of QCD. A good understanding of the uncertainties associated with the different measurements is needed in order to appreciate the significance of this test, which must be then restricted to observables where perturbative techniques are reliable and enough terms in the perturbative expansion are available. The PDG [1] requires a NNLO (or higher) theoretical accuracy. In addition, small non-perturbative corrections are always present, specially at low energies, and one should also worry about the expected asymptotic behaviour of the perturbative series.
The most reliable determinations of α s have been compiled in Refs. [2][3][4][5][6]. I will focus the discussion on the very precise inclusive observables R Z and R τ , which are already known to four loops, i.e., to N 3 LO, and will update the PDG information with the most recent developments, not yet included in the official averages.
(1) The fifth-order coefficient of the β function has been recently computed in Ref. [7] (see also Ref. [8]), which provides a quite precise perturbative control of the scale dependence of α s . In the MS scheme (β 1 and β 2 are scheme independent), the known coefficients are [7,9,10]: The very modest growth of β n with the perturbative order gives rise to a surprisingly smooth power expansion. For n f = 5, for instance, β(α s ) = β 1 a s 1 + 1.26 a s + 1.47 a 2 s + 9.83 a 3 s + 7.88 a 4 s . The scale dependence of α s over a wide range of energies, at different levels of approximation, is shown in figure 1. The 5-loop precision in the β function implies a resummation of N 4 LO logarithmic contributions to the running of α s , i.e., corrections of the form ∆α s (Q 2 ) ∼ α s (µ 2 ) n+5 log n (Q 2 /µ 2 ). Owing to the fast convergence of the β function, the NLO resummation gives already an excellent approximation to the running coupling. The achieved accuracy is quite impressive; the four and five loop corrections are so small that it is difficult to appreciate them in the figure.
The small discontinuities in the plotted curves reflect the crossing of the charm and bottom thresholds where one needs to properly match the different QCD n f effective theories. Since the β n coefficients are functions of n f , the strong coupling depends on the considered number of "active" quark flavours. When a quark is heavy enough to decouple, it is convenient to remove it from the Lagrangian and work with an effective QCD theory which has one quark less and a different value of α s . The matching conditions relating the effective QCD theories with n f and n f − 1 flavours are known to four loops [11,12].

Inclusive observables
Inclusive observables, such as σ(e + e − → hadrons) at high-enough energies, Γ(Z → hadrons) or Γ(W → hadrons), can be accurately predicted with perturbative methods. Since the final hadrons are produced through the vector V µ i j =ψ j γ µ ψ i and axial-vector A µ i j =ψ j γ µ γ 5 ψ i colour-singlet quark currents (i, j = u, d, s . . .), the QCD dynamics is governed by the two-point correlation functions where J = V, A and the superscript L = 0, 1 denotes the angular momentum in the hadronic rest frame. The correlators Π (L) i j,J (q 2 ) are analytic functions of q 2 , in the complex q 2 plane, except along the (physical) positive real axis where their imaginary parts have discontinuities which correspond to the measurable hadronic spectral distributions with the given quantum numbers.
For massless quarks, s Π (0) i j,J (s) = constant (there is a non-perturbative Goldstone-pole contribution to Π (0) i j,A at s = 0, which cancels in Π (0+1) i j,A ). When i j, the two quark currents must necessarily be connected through a quark loop (non-singlet topology), which gives identical contributions to the vector and axial massless correlators: Π(s) ≡ Π (0+1) i j,V (s) = Π (0+1) i j,A (s). They are conveniently parametrized through the Euclidean Adler function (Q 2 = −q 2 and N C = 3 is the number of quark colours) which is known to O(α 4 s ) [13][14][15]: There are additional singlet contributions to the neutral-current correlators (i = j), with each current coupling to a different quark loop. Since gluons have J PC = 1 −− and colour, these topologies start to contribute at O(α 3 s ) and O(α 2 s ), respectively, for the vector and axial-vector currents: The vector-current coefficients are d V 3 = −0.41318 and d V 4 = −5.94225 + 0.191628 n f [16].
The ratio of the electromagnetic e + e − → hadrons and e + e − → µ + µ − cross sections is given by The sum over quark electric charges of different signs strongly suppresses the singlet contribution, which has been included as a small correction to the coefficients F n≥3 . For n f = 5 flavours, one gets F 1 = 1, F 2 = 1.4092, F 3 = −12.805 and F 4 = −80.434 [16]. The perturbative series in Eq. (7) is actually an expansion in powers of α s (µ 2 ) with coefficients containing a polynomial dependence on log (s/µ 2 ). These logarithms are resummed into the running coupling by taking µ 2 = s. Although the physical ratio R e + e − (s) is independent of the renormalization scale µ, the truncated series contains a residual µ-dependence of O(α N+1 s ), where N = 4 is the last included term, which must be taken into account in the theoretical uncertainty. Since non-perturbative corrections are suppressed by Λ 4 /s 2 (the gauge-invariant operators contributing to the current correlators have dimensions D ≥ 4), at high energies one can perform a N 3 LO determination of α s (s). Unfortunately, the experimental uncertainties are large.

Γ(Z → hadrons)
The electroweak neutral current J µ contains vector and axial-vector components, weighted with the corresponding Z couplings. The singlet axial contributions of the two members of a weak isospin doublet cancel each other for equal quark masses because a f = 2I f ; however, the large value of the top mass generates very important singlet axial corrections which start at O(α 2 s ). The ratio of the hadronic and electronic widths of the Z boson involves the QCD series (m b = 0, m t 0) [16]. Taking properly into account the electroweak corrections and QCD contributions suppressed by powers of 17,18], the ratio R Z is included in the global fit to electroweak precision data. This results in a quite accurate value of α s (M 2 Z ) [19]: α This determination assumes the validity of the electroweak Standard Model.

Hadronic decay width of the τ lepton
The hadronic W ± decay width does not provide yet a competitive determination of α s . A much better alternative [20][21][22] is the hadronic τ decay, which proceeds through a virtual W ± boson. The QCD correlation function of two left-handed charged currents receives only non-singlet contributions.
Restricting the analysis to the dominant Cabibbo-allowed decay width, where S EW = 1.0201 ± 0.0003 incorporates the electroweak radiative corrections [23][24][25]. The measured invariant-mass distribution of the final hadrons determines the spectral functions ρ J (s) ≡ 1 π ImΠ (0+1) ud,J (s), shown in figure 2 (the only relevant contribution to the s ImΠ 0 ud,V+A (s) term is the π − final state at s = m 2 π ). Using the analyticity properties of the Π (L) i j,J (s) correlators, the experimental spectral distribution can be related with theoretical QCD predictions through moments of the type [22,27] A where s th is the hadronic mass-squared threshold, ω(s) is any weight function analytic in |s| ≤ s 0 , and the complex integral in the right-hand side (rhs) runs counter-clockwise around the circle |s| = s 0 . For large-enough values of s 0 , the operator product expansion (OPE) can be used to predict the rhs integral as an expansion in inverse powers of s 0 (the D = 0 term contains the perturbative contribution), while the lhs is directly determined by the experimental data. The ratio R τ,V+A in Eq. (10) corresponds to the particular weight ω( Thus, owing to Cauchy's theorem, the contour integral is only sensitive to OPE corrections with D = 6 and 8, which are strongly suppressed by the corresponding powers of the τ mass (there is in addition a further suppression of the D = 6 term because the vector and axial-vector contributions have opposite signs, cancelling to a large extent). Moreover, ω(s) contains a double zero at s = s 0 which heavily suppresses the contribution to the integral from the region near the real axis, where the OPE is not valid. This makes R τ,V+A a very clean observable to measure α s . It is very sensitive to the strong coupling because α s (m 2 τ ) is sizeable, and non-perturbative effects are smaller than the perturbative uncertainties. The availability of good experimental data makes possible to determine the small non-perturbative corrections from the data themselves, using weights with different powers of s which are sensitive to the corresponding power corrections in the OPE [27]. The dominant uncertainty in the α s (m 2 τ ) determination comes from the perturbative error associated with the unknown higher-order corrections to the Adler series in Eq. (4). For a given value of α s , the so-called contour-improved perturbation theory (CIPT) [28,29], which resumms large corrections arising from the long running along the circle s = s 0 , results in a smaller perturbative contribution than the truncated fixed-order perturbation theory (FOPT) approximation [22]. Therefore, CIPT leads to a larger fitted value of α s (m 2 τ ) than FOPT.

Numerical analysis
A detailed reanalysis of the α s (m 2 τ ) determination from τ decay has been recently performed [30], including many consistency checks to assess the potential size of non-perturbative effects. All strategies adopted in previous works have been investigated, studying the stability of the results and trying to uncover any potential hidden weaknesses, and several complementary approaches have been considered. Once their uncertainties are properly estimated, all adopted methodologies result in very consistent values of α s (m 2 τ ). Table 1 summarizes the most reliable determinations. All analyses have been done both in CIPT and FOPT. Within a given approach the perturbative errors have been estimated varying the renormalization scale in the interval µ 2 /s 0 ∈ [0.5 , 2], and taking K 5 = 275 ± 400 as an educated guess of the maximal range of variation of the unknown fifth-order contribution [31]. These two sources of theoretical uncertainty have been combined in quadrature, together with the experimental errors. The different values quoted in the table include, as an additional uncertainty, the variations of the results under various modifications of the fit procedures. The systematic difference between the values obtained with the CIPT and FOPT prescriptions appears clearly manifested in the table. The CIPT and FOPT results have been finally averaged, but adding in quadrature half their difference to the smallest of the CIPT and FOPT errors.
The first determination in table 1 follows the method adopted in the ALEPH analysis of Ref. [26], taking the weights ω kl (x) = (1 − x) 2+k x l (1 + 2x) with (k, l) = {(0, 0), (1, 0), (1, 1), (1, 2), (1, 3)} and s 0 = m 2 τ . With five moments, one can make a global fit of α s (m 2 τ ), the gluon condensate, O 6 and O 8 . To assess possible errors associated with neglected higher-order condensates, a second fit including O 10 has been performed and the variation on the fitted value of the strong coupling has been included as an additional uncertainty. A quite precise value of α s (m 2 τ ) is obtained, in good agreement with Ref. [26]. The extracted condensates have large relative errors exhibiting a very little sensitivity to power corrections. This has been further verified, taking away from the weights the factor (1 + 2x) which eliminates the highest-dimensional condensate contribution to every moment. This gives the  Neglecting all non-perturbative effects, one can determine α s (m 2 τ ) from a single moment. This interesting exercise has been also done in Ref. [30], making 13 separate extractions of the strong coupling with six A (2,m) moments (0 ≤ m ≤ 5), six A (1,m) moments (0 ≤ m ≤ 5) based on the weights ω (1,m) (x) = 1 − x m+1 = (1 − x) m k=0 x k which are only sensitive to O 2(m+2) , and the moment A (0,0) where OPE corrections are absent but it is very exposed to duality-violation effects. In all cases, the resulting determinations of the strong coupling are in agreement with the values in table 1, reflecting the minor numerical role of the neglected non-perturbative corrections.
Non-perturbative contributions should manifest in a distinctive s 0 dependence. Figure 3 shows as function of s 0 the experimental moments A (1,0) (s 0 ) and A (2,0) (s 0 ), in the V, A and 1 2 (V + A) channels, together with their predicted values with α (n f =3) s (m 2 τ ) = 0.329 + 0.020 − 0.018 , neglecting all non-perturbative contributions. A (1,0) (s 0 ), which can only get corrections from O 4 , exhibits a surprisingly good agreement with its pure perturbative prediction. In spite of being only protected by a single pinch factor, the data points above s 0 ∼ 2 GeV 2 closely follow the central values predicted by CIPT. In that energy range non-perturbative contributions appear to be too small to become numerically visible within the much larger perturbative uncertainties covering the shades areas of the figure. The splitting at lower values of s 0 of the V and A moments must be assigned to duality violations, since their D = 4 power corrections are approximately equal. However, these duality-violation effects clearly compensate in V + A, with an impressively flat distribution of the experimental data which does not deviate from the without OPE corrections. A (2,0) (s 0 ) looks slightly more sensitive to non-perturbative contributions and seems to prefer a power correction with different signs for V and A, which cancels to a good extend in V + A. This fits nicely with the expected O 6,V/A contribution, although the merging of the V, A and V + A curves above s 0 ∼ 2.2 GeV 2 suggests a very tiny numerical effect from this source at high invariant masses.
Fitting the s 0 dependence of a single A (2,m) (s 0 ) moment, one can determine the values of α s (m 2 τ ), O 2(m+2) and O 2(m+3) . The sensitivity to power corrections is very bad, as expected, but one finds an amazing stability in the extracted values of α s (m 2 τ ). Including the information from the three lowest moments (m = 0, 1, 2) and the nine energy bins above s 0 = 2.0 GeV 2 , and adding as an additional uncertainty the small fluctuations observed when changing the number of fitted bins, one obtains the values of α s (m 2 τ ) quoted in the fourth line of table 1. Although they are much more sensitive to violations of quark-hadron duality (fitting the s 0 dependence of several consecutive bins, one is using information about the local structure of the spectral function), these results turn out to be in excellent agreement with the more solid determinations in the first three lines of the table. The very flat shape of the V + A hadronic distribution above s 0 = 2.0 GeV 2 implies small duality-violation effects in that region which, moreover, are very efficiently suppressed in the doubly-pinched moments A (2,m) (s 0 ).
The marginal role of power corrections has been also corroborated, making independent α s (m 2 τ ) determinations from seven A (1,m) (s 0 ) (0 ≤ m ≤ 6) and six A (2,m) (s 0 ) (0 ≤ m ≤ 5) moments, as function of s 0 and ignoring all non-perturbative effects. In spite of the fact that these 13 moments get completely different OPE corrections, carrying a broad variety of inverse powers of s 0 , all results exhibit a similar functional dependence on s 0 . The small fluctuations among the different moments stay in all cases well within the much larger perturbative uncertainties shown in figure 3.
Using weights of the type ω (1,m) a (x) = (1− x m+1 ) e −ax , one suppresses potential violations of duality because the exponential factor nullifies the highest invariant-mass region, but paying the price that all condensates contribute to every moment. For a = 0 one recovers the A (1,m) (s 0 ) moments, only affected by O 2(m+2) , while for a 1 the moments become independent of m. Thus, if one neglects all non-perturbative contributions, the OPE corrections should manifest in a larger instability under variations of s 0 than in the a = 0 case. However, with a 0 one gets even more stable results, and the different moments converge very soon when a increases, indicating again that power corrections are not very relevant. From the analysis of seven V + A moments (m = 0, · · · , 6), accepting for each moment all values of α s (m 2 τ ) in the Borel-stable region, and adding as additional theoretical uncertainties the differences among moments and the variations in the region s 0 ∈ [2, 2.8] GeV 2 ), one gets the determination of α s (m 2 τ ) shown in the fifth line of table 1.

Violations of quark-hadron duality
The small differences between the true values of the moments A ω J (s 0 ) and their OPE approximations are known as (global) duality violations. Using analyticity, they can be fomally expressed as [32][33][34][35] ∆A ω,DV with ∆ρ DV V/A (s) the differences between the physical spectral functions and their OPE estimates which, unfortunately, are unknown beyond the experimentally accessed region. Owing to asymptotic freedom, the violations of duality should decrease very fast as s 0 increases. In practice, they are minimized by taking "pinched" weight functions which vanish at s = s 0 and suppress the contributions from the region near the real axis where the OPE is not valid [22,27]. The many tests discussed before clearly indicate that these effects are negligible in the extraction of α s (m 2 τ ) from the V + A distribution.  Instead of using clean moments where duality violations are suppressed, some works focus on observables more sensitive to these uncontrollable effects [36], modelling them with an ansatz for ∆ρ DV J (s) which is fitted to the measured spectral functions. Since the OPE is not valid on the physical cut, one loses theoretical control and gets at best an effective model description with unclear relation with QCD. Let us consider the slightly generalized ansatz (in GeV units) which for λ J = 0 coincides with the model assumed in Ref. [36]. The combination of a dumping exponential with an oscillatory function is expected to describe the fall-off of duality violations at very high energies, but this functional form is completely ad-hoc and difficult to justify at low energies.
Since there are far too many parameters to be fitted to a highly-correlated data set, Ref. [36] concentrates in the moment A (0,0) V (s 0 ) which is very exposed to violations of duality (ω(x) = 1) and does not receive OPE corrections (owing to the tail of the a 1 resonance, the axial channel is not very useful). The model parameters and α s are determined fitting the s 0 dependence for s 0 ≥ŝ 0 = 1.55 GeV 2 . This choice has the largest, but still too small, p-value and gives the smallest α s . However, the p-value falls dramatically when one moves from this point, becoming worse at higherŝ 0 values where the model should work better. The extracted value of α s is very unstable under small modifications of the fit procedure and the fitted ansatz strongly deviates from the data as soon as one moves from the fitted region. This is illustrated in figure 4 which shows the results of this exercise with FOPT, for different values of the power λ V andŝ 0 = 1.55 GeV 2 . The actual uncertainties are much larger than the quoted fit errors; varyingŝ 0 in the range [1.15, 1.75] GeV 2 , with λ V = 0, induces 3σ fluctuations of α s (m 2 τ ). All models reproduce well ρ V (s) in the fitted region (s ≥ 1.55 GeV 2 ), but they fail badly below it. The choice λ V = 0 assumed in Ref. [36] is clearly the worse one. Increasing the power λ V , the ansatz slightly approaches the data below the fitted range, while the exponential parameters δ V and γ V adapt themselves to compensate the growing at high values of s with the net result of a smaller duality-violation correction. The statistical quality of the fit improves also with growing values of λ V , while α s (m 2 τ ) increases approaching the more solid FOPT values in table 1. The strong correlation of the fitted α s (m 2 τ ) with the assumed model should not be a surprise because one is just fitting models to data without any strong theoretical guidance (the OPE is no longer valid), and α s has been converted into one more model parameter. In spite of all caveats, one gets still quite reasonable values of the strong coupling, but they are model dependent and, thus, unreliable.

Updated determination of α s (m τ )
The results shown in table 1 are based on solid theoretical principles (the s 0 -dependence extraction assumes, however, local duality) and exhibit a good stability under small variations of the fit procedures. The overall agreement among determinations extracted under very different assumptions shows their reliability and even indicates that the uncertainties are probably too conservative. Averaging the five determinations, but keeping the smaller uncertainties to account for the large correlations, one finds The same results are obtained irrespective or whether one includes or not in the average the determination from the s 0 dependence of the moments. Averaging the CIPT and FOPT "averages" in table 1, one finally gets α These results nicely agree with the value of the strong coupling extracted from R τ [37]. After evolution up to the scale M Z , the strong coupling decreases to in excellent agreement with the direct measurement at the Z peak in Eq. (9). The comparison of these two determinations, graphically shown in the right panel of figure 1, provides a beautiful test of the predicted QCD running; i.e., a very significant experimental verification of asymptotic freedom: Improvements on the determination of α s (m 2 τ ) from τ decay data would require high-precision measurements of the spectral functions, specially in the higher kinematically-allowed energy bins. Both higher statistics and a good control of experimental systematics are needed, which could be possible at the forthcoming Belle-II experiment. On the theoretical side, one needs an improved understanding of higher-order perturbative corrections. Figure 5 compares the N 3 LO determinations of α s (M 2 Z ) from Z and τ decays with other precise measurements of the strong coupling. Following the PDG criteria [2], only those determinations which are at least of NNLO are taken into account. This includes several event-shape analyses in hadronic final states of e + e − annihilations, and studies of parton distribution functions from deep inelastic scattering and hadron collider data. The numbers quoted in the figure correspond in both cases to the recent PDG compilation [2].

World average value of α s
The PDG includes also in the average the CMS determination from the tt production cross section at √ s = 7 TeV, α s (M 2 Z ) = 0.1151 + 0.0028 − 0.0027 [38], which requires as input a value of the top quark mass (either α s or m t are fitted to the data, but not both). Although there are more recent measurements of this cross section from ATLAS and CMS, at √ s = 7, 8 and 13 TeV, none of them quotes further determinations of α s . Applying the same procedure, these measurements would imply larger values of α s (M 2 Z ) than the one of Ref. [38], which is nevertheless included in the average. The most precise value of α s (M 2 Z ) is obtained from lattice simulations, with a growing number of groups measuring (non-perturbatively) various short-distance quantities, through numerical evaluations of the QCD functional integral, and comparing the results with the corresponding perturbative This number is very close to the 2016 PDG average (0.1181 ± 0.0011), which does not yet include the most recent τ decay and lattice results. The central value has been directly obtained as the weighted average of the six input determinations, while the error has been enlarged applying the PDG prescription, i.e., adjusting all individual uncertainties by a common factor such that χ 2 /dof equals unity. Removing the CMS determination would slightly increase the central value, giving as average α s (M 2 Z ) = 0.1183 ± 0.0011. The overall uncertainty is largely determined by the precise lattice result.