The muon $g-2$ and $\alpha(M_Z^2)$: a new data-based analysis

This work presents a complete re-evaluation of the hadronic vacuum polarisation contributions to the anomalous magnetic moment of the muon, $a_{\mu}^{\rm had, \, VP}$ and the hadronic contributions to the effective QED coupling at the mass of the $Z$ boson, $\Delta\alpha_{\rm had}(M_Z^2)$, from the combination of $e^+e^-\rightarrow {\rm hadrons}$ cross section data. Focus has been placed on the development of a new data combination method, which fully incorporates all correlated statistical and systematic uncertainties in a bias free approach. All available $e^+e^-\rightarrow {\rm hadrons}$ cross section data have been analysed and included, where the new data compilation has yielded the full hadronic $R$-ratio and its covariance matrix in the energy range $m_{\pi}\leq\sqrt{s}\leq 11.2$ GeV. Using these combined data and pQCD above that range results in estimates of the hadronic vacuum polarisation contributions to $g-2$ of the muon of $a_{\mu}^{\rm had, \, LO \, VP} = (693.27 \pm 2.46)\times 10^{-10}$ and $a_{\mu}^{\rm had, \, NLO \, VP} = (-9.82 \pm 0.04)\times 10^{-10}$. The new estimate for the Standard Model prediction is found to be $a_{\mu}^{\rm SM} = (11\ 659 \ 182.05 \pm 3.56) \times 10^{-10}$, which is $3.7\sigma$ below the current experimental measurement. The prediction for the five-flavour hadronic contribution to the QED coupling at the $Z$ boson mass is $\Delta\alpha_{\rm had}^{(5)}(M_Z^2)= (276.11 \pm 1.11)\times 10^{-4}$, resulting in $\alpha^{-1}(M_Z^2) = 128.946 \pm 0.015$. Detailed comparisons with results from similar related works are given.

1 Introduction The anomalous magnetic moment of the muon, a µ = (g − 2) µ /2, stands as an enduring test of the Standard Model (SM), where the ∼ 3.5σ (or higher) discrepancy between the experimental measurement a exp µ and the SM prediction a SM µ could be an indication of the existence of new physics beyond the SM. For a exp µ , the value is dominated by the measurements made at the Brookhaven National Laboratory (BNL) [1][2][3], resulting in a world average of [4]  Efforts to improve the experimental estimate at Fermilab (FNAL) [5] and at J-PARC [6] aim to reduce the experimental uncertainty by a factor of four compared to the BNL measurement. It is therefore imperative that the SM prediction is also improved to determine whether the g − 2 discrepancy is well established. The uncertainty of a SM µ is completely dominated by the hadronic contributions, a had µ , attributed to the contributions from the non-perturbative, low energy region of hadronic resonances. The hadronic contributions are divided into the hadronic vacuum polarisation (VP) and hadronic light-by-light (LbL) contributions, which are summed to give a had µ = a had, VP µ + a had, LbL µ . (

1.2)
This analysis, KNT18, is a complete re-evaluation, in line with previous works [7][8][9], of the hadronic vacuum polarisation contributions, a had, VP µ . The hadronic vacuum polarisation contributions can be separated into the leading-order (LO) and higher-order contributions, where the LO and next-to-leading order (NLO) contributions are calculated in this work. 1 These are calculated utilising dispersion integrals and the experimentally measured cross section, where α = α(0), s th = m 2 π , R(s) is the hadronic R-ratio given by R(s) = σ 0 had,γ (s) σ pt (s) ≡ σ 0 had,γ (s) 4πα 2 /(3s) (1.5) and K(s) is a well known kernel function [10,11], given also in equation (45) of [7], but differing by a normalisation factor of m 2 µ /(3s). This kernel function (which behaves as K(s) ∼ m 2 µ /(3s) at low energies), coupled with the factor of 1/s in the integrand of equation (1.4), causes the hadronic vacuum polarisation contributions to be dominated by the low energy domain. At NLO, the data input is identical, with corresponding dispersion integrals and kernel functions [12] (see also [7]). In the previous analysis [9] (denoted as HLMNT11), the LO hadronic vacuum polarisation contributions were found to be This, compared with the experimental measurement, gave a g − 2 discrepancy of 3.3σ. In addition to calculating a had, VP µ , the combination of hadronic cross section data is also used to calculate the running (momentum dependent) QED coupling, α(q 2 ). This, in particular, is then used to determine the effective QED coupling at the Z boson mass, α(M 2 Z ), which is the least precisely known of the three fundamental electro-weak (EW) parameters of the SM (the Fermi constant G F , M Z and α(M 2 Z )) and hinders the accuracy of EW precision fits. Using an identical data input as that used for a had, VP µ , the hadronic contributions to the effective QED coupling are determined from the dispersion relation ∆α (5) had (q 2 ) = − , (1.8) where the superscript (5) indicates the contributions from all quark flavours except the top quark, which is added separately. Together with the leptonic contributions, this is used to determine the running coupling α(q 2 ) = α/ 1 − ∆α had (q 2 ) − ∆α lep (q 2 ) . The structure of this paper is as follows. Section 2 describes all changes to the treatment and combination of hadronic cross section data since [9]. Section 3 gives details concerning the contributions to a had, LO VP µ and ∆α had (M 2 Z ) from the different hadronic final states and energy regions, culminating in updated estimates of both a had, LO VP µ , a had, NLO VP µ and ∆α had (M 2 Z ). These results are then compared with other works and the new prediction of a had, VP µ is combined with all other SM contributions to determine the SM prediction a SM µ , with details given concerning the consequences of this on the existing g − 2 discrepancy and the new result for α(M 2 Z ). Finally, the conclusions of this work, along with discussions of the prospects for g − 2 in the future and potential improvements for a SM µ , are presented in Section 4. . Any new and old data that have not been corrected for VP effects require undressing. However, recent data are more commonly undressed in the experimental analyses already, removing the need to apply a correction to these data sets. This benefits the data combination as new, more precise data undressed of VP effects are dominating the fit for many channels which, in turn, reduces the impact of the extra radiative correction uncertainty which is applied to each channel.
Following the previous analyses, the 'bare' cross section σ 0 had in equation (1.5) is determined to be σ 0 had = C vp σ had , (2.1) where C vp is a multiplicative VP correction. For this purpose, the self-consistent KNT18 vacuum polarisation routine, vp knt v3 0, has been updated via an iteration of the new data input. 2 An in-depth discussion of the determination of the vacuum polarisation and the corresponding routine will be given in [14]. The decisions to undress all data sets previously corrected in [9] are unchanged, with the exception of two measurements of the inclusive hadronic R-ratio by the Crystal Ball collaboration [15,16]. A re-examination of the experimental analyses of these two measurements has shown that not only do they include some treatment of the VP, but that they also both account for this correction with sizeable systematic uncertainties. Applying an extra VP correction here would be an overestimate of the correction and its corresponding uncertainty. Therefore, undressing of these data is not applied here. In addition, it should be noted that in a separate work [17], two measurements by the KLOE collaboration [18,19] in the π + π − channel are now undressed of VP effects using an updated routine [20] compared to the one used previously [21].
The undressing of narrow resonances in the cc and bb regions requires special attention. Importantly, the electronic width of an individual resonance, Γ ee , must be undressed of vac-uum polarisation effects using a parametrisation of the VP where the correction excludes the contribution of that resonance, such that Here, α no res is the effective QED coupling neglecting the contribution of the resonance itself and is given by where ∆α no res (s) is determined from equation (1.8) such that the input R(s) does not include the resonance that is being corrected. To include the resonance would lead to an inconsistent definition of the narrow resonance.

Final state radiative corrections
Final state radiation (FSR) cannot be separated in an unambiguous way in the measured hadronic cross sections. Therefore, while formally of higher order in α, FSR photons have to be taken into account in the definition of the one particle irreducible hadronic blobs and will already be included as part of the leading order hadronic VP contributions. However, depending on the experimental analyses, some amount of real photon FSR may have been removed during the event selection. Adding back these missed contributions is model dependent and not feasible for general hadronic final states. It is therefore necessary to estimate the possible effects and their impact on the accuracy of the data compilations. Similar to the case of VP corrections, an extra radiative correction uncertainty due to FSR effects is then estimated channel by channel. For the important π + π − and K + K − channels, detailed studies have been performed for this analysis. Here, and especially in the limited energy range below 2 GeV, it has been shown that scalar QED provides a good description of photon FSR, see e.g. [22][23][24]. In the past, to estimate possible FSR effects in π + π − and K + K − production, the fully inclusive, order α correction to the cross section, has been used. 3 Here, the function η(s) is given e.g. in [22] and the subscript γ indicates the one photon inclusive cross section. However, experimental cross section measurements, by nature, include all virtual and soft real radiation effects. Therefore, to estimate possibly missing photon FSR, ideally only the effects from (hard) real radiation above/within resolution/cut parameters which are specific for a given experiment or analysis have to be estimated. Whereas in the calculation of the inclusive correction a regularisation of the virtual and real soft contributions is required to obtain the infrared finite result η, the hard real radiation, η hard,real , can be estimated numerically from where m is the mass of the (scalar) particle, Λ is a finite infrared cut-off parameter on the invariant mass of the emitted photon and ρ fin is the radiator function (see appendix B of [23]). σ 0 (e + e -→ K + K -) Figure 1: The effect of final state radiation in the K + K − channel in the φ resonance region. Left panel: the fully inclusive FSR correction η(s). Right panel: hard real radiation η hard,real (s), estimated with accolinearity cuts used in the two analyses [25,26]. The e + e − → K + K − cross section is also plotted for reference.
In the case of the K + K − channel, by far the largest contribution to a µ and ∆α (and their errors) comes from the energy region of the φ peak, where the phase space for real radiation is severely restricted. While the calculation of real radiation accounting for all experimental cuts would be very complicated and beyond the scope of this work, an estimate can easily be made based on equation (2.5), by relating Λ to the cuts in the photon accolinearity given in individual experimental analyses. In Figure 1, the result for the fully inclusive correction η(s) (left panel) is compared to the estimates for the real hard radiation, where η hard,real (s) (right panel) now depends on the accolinearity cuts as given by the two experimental analyses depicted. Clearly, at and around the φ peak, phase space restrictions strongly suppress any hard real radiation and using the inclusive correction would lead to an overestimate of the possible effects. Given the small size of both the possible correction in the φ region and the contribution of the K + K − channel above the φ to both mean value and error of a µ and ∆α, no correction or additional error estimate due to FSR is now applied in the K + K − channel. For the neutral kaon channel, hard photon radiation (which would resolve the quark charges) is similarly suppressed and no FSR correction or additional error are applied in this channel also.
The situation is different in the two pion channel. A study similar to the two kaon channel showed that in principle larger contributions from real radiation of the order of the inclusive correction can arise. However, these contributions are strongly dependent on the cut applied in equation (2.5) and would require a more detailed, measurement-by-measurement analysis, which is beyond the scope of this work. In addition, for many of the data sets the required detailed information is not available. Therefore, in sets which are understood to not include the full FSR corrections, the fully inclusive correction (as shown in Figure 2) is applied to the respective measurements. In the π + π − data combination, recent sets from radiative return, where additional photons are part of the leading order cross section and are studied in detail as part of the analyses, have now become dominant. Therefore, the impact of the fully inclusive FSR correction to older sets is suppressed for both mean value and error in comparison to previous analyses. The procedure to determine the corresponding error estimate is discussed in the next chapter below, whilst resulting numbers are given in Section 3.1, which contains the detailed analysis of the π + π − channel.  Figure 2: The behaviour of the inclusive FSR correction, η(s), for the process e + e − → π + π − . The e + e − → π + π − cross section is also plotted for reference.
For the sub-leading, multi-hadron channels, there are, at present, no equivalent FSR calculations. Depending on the experimental analysis, they are (at least to some extent) estimated by Monte Carlo (MC) simulation and contribute to the systematic uncertainties. 4 However, for many data sets, it is far from clear to which extent FSR effects are included in the systematic errors. Therefore, possible effects are accounted for by applying an additional uncertainty determined as a fraction of the respective contribution.

Estimating extra uncertainty due to radiative corrections
As in [7][8][9], extra uncertainties are estimated whenever additional radiative corrections are applied. This is done, first and foremost, to account for any under-or over-correction that may occur due to a lack of information concerning the treatment of radiative corrections in the experimental analyses. However, these radiative correction uncertainties also account for any possibly incorrect treatment in the analyses, e.g. missed FSR or inconsistent subtraction of VP contributions. This is especially true for older data, where there is very little or even no information at all regarding how the data have been treated.
In each channel, the difference ∆a vp µ between the estimates of a µ with and without additional VP corrections is determined. For the uncertainty due to VP, one third of the shift is taken. This is reduced from one half in [9], safe in the knowledge that the KNT18 VP routine [14], which is determined iteratively in a self-constistent way, is accurate to the level of a permille when correcting the cross section and that newer data sets are commonly undressed of VP effects by the experiment with a modern routine. For the extra uncertainties due to FSR, there are now no contributions from the K + K − and K 0 S K 0 L channels (see the discussion above). For the π + π − channel, the full difference between the estimates of a π + π − µ with and without additional FSR corrections is taken as the FSR uncertainty. 4 In the recent work [27], a version of the MC program Carlomat [28] is discussed, which includes real photon radiation in multi-hadron channels. This could be used to further study the FSR effects.
For all other channels, including the inclusive data combination above 1.937 GeV, a fraction of 1% of the respective cross section is applied as the uncertainty.
The numerical estimates of all additional radiative correction uncertainties are given in the respective sections for the individual channels. The same procedures are applied in the calculation of the contributions to ∆α had .

Data combination
The combination of σ 0 (e + e − → hadrons) data and its corresponding uncertainty has undergone a great deal of scrutiny since [9]. This has been due to the ever increasing amount of data that have become available and due to a better understanding of the treatment of correlated uncertainties that are now a prominent feature of many cross section measurements. In [9], covariance matrices in the π + π − channel from BaBar [29] and KLOE [18,19] for statistical and systematic uncertainties had already been included based on the method of fitted renormalisation factors. In the following, a new method for the combination of σ 0 had data is introduced, which will allow to fully take into account all correlated uncertainties.

Clustering data points
Within each hadronic channel, data points from different experiments are assigned to clusters (see [7][8][9] for details). 5 In this work, the clustering algorithm has been improved. It is universal for all different channels and only differs in the assigned maximum cluster size δ (and δ res , a maximum cluster size applicable at individual narrow resonances). A scan over δ (and δ res if applicable) is performed to determine a suitable clustering configuration which must avoid both over-and under-clustering. Too wide or overpopulated clusters would effectively lead to a rebinning of data points from individual experiments and risk loss of information, while a too narrow clustering would result, in the extreme, in an erratic point-to-point representation of the cross section and no gain in the accuracy. The preferred configuration is then chosen based on the global χ 2 min /d.o.f. and the uncertainty on a had,VP µ , combined with checks by eye of the resulting spectral function.

Fitting data
The previous analyses [7][8][9] employed a non-linear χ 2 -minimisation utilising fitted renormalisation factors as nuisance parameters that represented the energy independent systematic uncertainties. Although a powerful method, recent literature [30] (see also [31]) have highlighted the possibility that an incorrect treatment of multiplicative normalisation uncertainties in a χ 2 minimisation can incur a systematic bias (see chapter 4 of [30]). In addition, although the non-linear χ 2 -minimisation used in [9] was adjusted to include covariance matrices, the method's reliance on fitting energy independent renormalisation factors prevented the use of correlated uncertainties to their full capacity. As recent precise data, specifically radiative return measurements in the π + π − and K + K − channels, have been released with energy dependent uncertainties and non-trivial bin-to-bin correlations for both the statistical and systematic uncertainties, the fit procedure had to be modified to allow the full use of all available correlations in a bias-free 5 The effect that different clusterings have on the χ 2 min /d.o.f. of the fit, the resulting a had,VP µ value and its error has been discussed in detail in [9].
approach. Therefore, instead of fitting renormalisation factors (nuisance parameters), an iterative fit procedure as advocated in [30] has been adopted, which re-initialises the full covariance matrices at each iteration step and, consequently, avoids bias.
Previously, in [7][8][9], a constant cross section had been assumed across the width of each cluster. In this work, the fitted cross section values at the cluster centres are obtained from the iterative χ 2 -minimisation where the cross section is taken to be linear between adjacent cluster centres. This allows for a more stable fit and is consistent with the trapezoidal rule integration utilised for the a had,VP µ and ∆α had integrals. If data points at energies are combined into cluster m, then the weighted average of the cross section value R m and energy E m for the cluster centre are and where R For data points at the borders where no interpolation is possible, R i m is found by linear extrapolation. A covariance matrix is constructed for the combination which contains all the uncertainty and correlation information of all data points. Using the covariance matrix as defined by the data alone could result in bias (see [32,33]). The covariance matrix is therefore redefined at each step of the iteration using the fitted R m values. Convergence of the iteration is observed in this work to occur after only a few steps.
The covariance matrix C i (m) , j (n) is given as the sum of the statistical covariance matrix C stat i (m) , j (n) and the systematic covariance matrix C sys i (m) , j (n) . At each stage of the iteration, it is defined as where the quantities R i,I m and R j,I n are the interpolant cross sections given by either equation (2.10) or (2.11) and I denotes the iteration number of the fit. This is then used as input into the now linear χ 2 -function, where N tot is the total number of contributing data points and C −1 I i (m) , j (n) is simply the inverse of the covariance matrix defined in equation (2.12). Performing the minimisation yields a system of linear equations where, The solution to this yields the cluster centres R m and the covariance matrix V I m, n which describes the correlation between the errors dR m and dR n . As in equations (2.10) and (2.11), the term R j n is to be taken as either where δ denotes the usual Kronecker delta.
In [30], it was stated that for the fit of parton distribution functions, convergence is expected to occur after very few iterations, which is also observed here. The final result includes the total output covariance matrix V (m, n) which, as in [9], is inflated according to the local This is done in order to account for any tensions between the data. The use of the full covariance matrix allows for the inclusion of any-and-all uncertainties and correlations that may exist between the measurements. The flexibility to now make use of fully energy dependent uncertainties ensures that the appropriate influence of the correlations is incorporated into the determination of the cluster centres, R m , with the correct propagation of all experimental errors to the g − 2 uncertainty.
Note that, for all channels, the differences between the old and the new data combination procedures are small and lead to changes of a had,VP µ well within the quoted errors. Importantly, combining data which have only global normalisation uncertainties results in negligible differences between [9] and this work, indicating that previous results were largely unaffected by the potential bias issue.  resonance regions in the π + π − and K + K − channels.

Integration
As in [7][8][9], the data are integrated using a trapezoidal rule integral in order to obtain a had,VP µ and ∆α had . In principle, the use of the trapezoidal rule integral could lead to unreliable results due to the form of the kernel function or at narrow resonances if data are sparse. However, with the current density of cross section measurements, especially in the dominant hadronic channels, the differences between trapezoidal rule integration and any higher order polynomial approximation are consequently small and of no concern. This can be seen in the plots in Figure 3. The calculation of the uncertainty of a had,VP µ and ∆α had has been modified to improve the determination of the error contribution at the integral boundaries. Should the upper or lower integral boundary, E a , with energy E m < E a < E m+1 and cross section value R m < R a < R m+1 , be found by linear interpolation (or extrapolation if it is necessary to extend the integral boundaries), then the output covariance matrix V m, n is interpolated accordingly using the standard error propagation formula where m, n run over all clusters, b is a label for any energy E b and V b, a = V a, b .
3 Determining a had,VP µ and ∆α had (M 2 Z ) The following section summarises the data combination and estimates of a had, LO VP µ and ∆α had (M 2 Z ) from the leading and major sub-leading hadronic final states. All contributions from exclusive hadronic channels are evaluated up to 1.937 GeV, which is the chosen transition point between the sum of exclusive channels and the inclusive R-ratio data in this work. This is discussed in detail in Section 3.6.1. Each contribution to a had, LO VP µ is quoted with its respective statistical (stat) uncertainty, systematic (sys) uncertainty, VP (vp) correction uncertainty and FSR (fsr) correction uncertainty individually. This is followed by the contribution with a total (tot) uncertainty, determined from the individual sources added in quadrature, All results for a µ are given in units 10 −10 . For the contributions to ∆α had (M 2 Z ), only the mean value with the corresponding total uncertainty, is quoted and all results are given in units 10 −4 . In both cases, uncertainties include all available correlations and local χ 2 inflation as discussed in Section 2.2.2. In the following, unless stated explicitly, only those data sets that are new additions since [9], or those that have undergone studies that are mentioned in the text, are discussed and referenced. All other contributing data sets are referenced in [7][8][9].

π + π − channel
The π + π − channel is by far the most important contribution to a had, VP µ , dominating both its mean value and uncertainty. With twenty-six data sets totalling almost one thousand data points, it is also the largest individual data combination. Two new radiative return measurements from the KLOE collaboration [34] and the BESIII collaboration [35] in the important ρ region have greatly improved the estimate of this final state. The new data from the KLOE collaboration (KLOE12) agree well with both other measurements by KLOE (KLOE08 [18] and KLOE10 [19]). The three measurements are, in part, highly correlated, necessitating the construction of full statistical and systematic covariance matrices to be used in any combination of these data. This has been achieved in the separate work [17], where all details regarding the correlations, the resulting combination of the KLOE data and a comparison with other experimental measurements of σ ππ(γ) are presented. These covariance matrices are used here as input in the full π + π − combination in order to fully incorporate the correlation information of the KLOE data as an influence on both the estimate of a π + π − µ and its uncertainty. The BESIII measurement in the ρ resonance region (again with full statistical and systematic covariance matrices) allows for an in-depth comparison of the existing radiative return measurements already contributing to the π + π − channel, namely the three measurements by the KLOE collaboration and the finely binned measurement from the BaBar collaboration [29]. In [9], details were given regarding tension between the KLOE and BaBar measurements, where the BaBar data were considerably higher. As is evident from Figure 4, tension exists between the BaBar data and all other contributing data in the dominant ρ region. When considering this along with the plots of the resulting cross section in Figure 5, it is clear that the new BESIII data agrees well with the KLOE data and the full π + π − combination. Interestingly, however, it is in better agreement with the BaBar data at the peak of the resonance where the cross section is largest. Although BaBar still influences with an increase, the agreement between the other radiative return measurements and the direct scan data largely compensates for this effect. This is demonstrated by Figure 6, with the combination clearly favouring the other measurements. Tension between data sets, however, still exists and is reflected in the local χ 2 error inflation, which results in an ∼ 15% increase in the uncertainty of a π + π − µ . The effect of this energy dependent error inflation is shown in Figure 7, where the difference in using a local scaling of the  error instead of a global one is clearly visible. Tensions arise in particular in the ρ resonance region, where the cross section is large.
The full combination of all π + π − data is found to give 937 GeV] = 502.97 ± 1.14 ± 1.59 ± 0.06 ± 0.14 = 502.97 ± 1.97 (3.3) and BaBar (09) Fit of all π + πdata CMD-2 (03) SND (04) CMD-2 (06) KLOE combination BESIII (15) Global √(χ 2 min /d.o.f.) = 1.30 a µ π + π -(0.6 ≤  √s ≤ 0.9 GeV) = (369.41 ± 1.32) x 10 -10 Figure 6: The relative difference of the radiative return and important direct scan data sets contributing to a π + π − µ and the fit of all data. For comparison, the individual sets have been normalised against the fit and have been plotted in the ρ region. The light green band represents the BaBar data and their errors (statistical and systematic, added in quadrature). The yellow band represents the full data combination which incorporates all correlated statistical and systematic uncertainties. However, the width of the yellow band simply displays the square root of the diagonal elements of the total output covariance matrix of the fit.
Although this value of a π + π − µ stays well within the error estimate of [9], it exhibits a substantial decrease of the mean value. This has been attributed to the new data combination routine which allows for the full use of correlations in the determination of the mean value as well as the uncertainty and the inclusion of the new, precise radiative return data which suppresses the influence of BaBar in the ρ resonance region.
In comparison with equation (3.3), the BaBar data alone in the same energy range give an estimate of a π + π − µ (BaBar data only) = 513.2 ± 3.8. Should all available π + π − data be combined using a simple weighted average as in equation (2.7), which only provides the error weighting to each cluster by its local uncertainty, the estimate for a π + π − µ would be a π + π − µ (Naive weighted average) = 509.1 ± 2.9. In this case, the estimate is strongly pulled up by the fine binning and high statistics of the BaBar data that dominate when no correlations are taken into account for the mean value. This difference of nearly 2σ when comparing to equation (3.3) indicates the importance of fully incorporating all available correlated uncertainties in any combination of the data.
The uncertainty has reduced by approximately one third. Again, this is due to the new, precise radiative return data which further dominate the π + π − fit and the improvement of the overall combination which now fully incorporates the energy dependent correlations. In addition, the radiative corrections uncertainties have reduced since [9], as discussed in detail in Section 2.1. f. in the π + π − channel, which is plotted against the e + e − → π + π − cross section for reference.
3.2 π + π − π 0 channel Since [9], there has only been one new addition to the π + π − π 0 channel [36]. This new data set improves this channel away from resonance, where previously only BaBar data [37] had provided a contribution of notable precision. Compared to [9], an additional change is applied to three separate data scans over the φ resonance in a measurement by CMD-2 [38], where the systematic uncertainties between the three scans are now taken to be fully correlated [39]. These changes, along with the new data combination routine, have resulted in an improved estimate of Figure 8 shows the full integral range of the data for the π + π − π 0 cross section. Figure 9 shows an enlargement of the ω and φ resonance regions in this channel.
this channel is and where there is clear improvement since [9]. The uncertainty contribution from π + π − π 0 π 0 is, however, still relatively large in comparison with its contribution to a had, LO VP µ and ∆α had and requires better new data. The new fit of the data for the bare cross section e + e − → π + π − π 0 π 0 is shown in Figure 11. The fit without the new BaBar data is shown for comparison to highlight the large improvement this new data set provides. Fit of all K + Kdata DM1 (81) DM2 (83) BCF (86) DM2 (87) OLYA (81) CMD (91) CMD-2 (95) SND (00) Scans

KK channels
The K + K − channel now includes a precise and finely binned measurement by the BaBar collaboration, supplemented with full statistical and systematic covariance matrices [45]. This is the first and only example to date of the release of energy dependent, correlated uncertainties outside of the π + π − channel and has an overwhelming influence on the data combination. There is also a new measurement in this channel of the φ resonance by the CMD-3 collaboration [46]. The existing CMD-2 scans in the same region [47] are omitted from this work as they suffer from an overestimation of the trigger efficiency for slow kaons [46] and are awaiting reanalysis [48]. In addition, two new scans by the SND collaboration measured at the tail of the φ and into the continuum are included [49]. The systematic uncertainties of these two scans, along with the  existing two scans by SND [26], are considered to be fully correlated [39]. The combination of the available K + K − data now gives 937 GeV] = 23.03 ± 0.08 ± 0.20 ± 0.03 ± 0.00 = 23.03 ± 0.22 (3.11) and which exhibits an increase of the mean value of more than 1σ from the estimate in [9] attributed to the inclusion of the new BaBar and CMD-3 data. The cross section of the process e + e − → K + K − is displayed in Figure 12. The uncertainty has drastically improved since [9] with much of the change being due to a finer clustering over the φ resonance after the inclusion of the new high statistics BaBar data.
Following from the discussion in Section 2.1.2, that there is now no FSR correction applied to this channel and, therefore, there is no extra radiative correction uncertainty due to FSR. It should also be noted that any FSR correction would result in an increase of a K + K − µ , showing again the strong influence the new data have had in this channel to increase the mean value since [9], where previously an FSR correction was applied.
New data for the K 0 S K 0 L final state is included from the BaBar collaboration above the φ resonance [50] and from the CMD-3 collaboration on the φ [51]. Two existing measurements in this channel [26,52] have multiple data scans of which the systematic uncertainties are now taken to be fully correlated [39]. This combination results in a contribution of 937 GeV] = 13.04 ± 0.05 ± 0.16 ± 0.10 ± 0.00 = 13.04 ± 0.19 (3.13) and ∆α K 0 Again, there is also no additional FSR correction uncertainty applied to this channel, with the reasoning in Section 2.1.2 enforced by the probability of photon emission being highly suppressed for a neutral final state and given the limited phase space. The cross section of the process e + e − → K 0 S K 0 L is displayed in Figure 13.

KKπ and KK2π channels
Since [9], the neutral final state K 0 S K 0 L π 0 has been measured by BaBar [53] and SND [54], completing all modes that contribute to the KKπ final state and removing the reliance on isospin for this channel (other than K 0 S ∼ = K 0 L ). Therefore, the KKπ cross section is now calculated using In [9], the isospin estimate in the same energy range yielded, a KKπ µ (HLMNT11 isospin estimate) = 2.65 ± 0.14. (3.18) This good agreement between the HLMNT11 isospin estimate and the data-based approach in this analysis is also demonstrated in Figure 14.
For KK2π, BaBar have measured the previously missing modes K 0 , such that this contribution is now determined using   Here, again, it is assumed that K 0 S ∼ = K 0 L and, hence,  and examining Figure 15, it is evident that the isospin relations provided a poor estimate of this final state. Using the data, KK2π contributes a much smaller mean value with a greatly reduced uncertainty.

Inclusive R-ratio data
The combination of inclusive hadronic R-ratio data between 1.937 ≤ √ s ≤ 11.2 GeV has three new data additions since [9]. The first of these are the precise BaBar R b data between 10.54 ≤ √ s ≤ 11.20 GeV [56], which must be deconvoluted of the effects from initial state radiation (ISR) and must have the radiative tails of the resonances from the Υ(1S − 4S) states removed (see [57]). These data are then added to the perturbative QCD (pQCD) estimate of R udsc [58] to be included as an accurate data set in the inclusive channel in this region. The inclusion of the BaBar R b data is particularly beneficial as it resolves the resonances of the Υ(5S) and Υ(6S) states, removing the need to estimate these structures using a Breit-Wigner parametrisation as was done in [7][8][9]. These data are shown in Figure 16, together with the previously used incoherent sum of the resonance parametrisations which are clearly very different from the bb cross section as measured by BaBar. Note that apart from the CLEO(98) data point [59] at 10.52 GeV, the CLEO(07) data point [60] at 10.538 GeV and the CUSB data point at 11.09 GeV [61], there are no other data in this bb resonance region.
The other new data are recent precise measurements by the KEDR collaboration: one set between 1.84 ≤ √ s ≤ 3.05 GeV [62] and two scans in the energy range 3.12 ≤ √ s ≤ 3.72 GeV [63]. For the latter, the systematic uncertainties are taken to be fully correlated [39]. The fit of the inclusive data in the range 1.937 ≤ √ s ≤ 3.80 GeV is shown in Figure 17, which demonstrates the good agreement between KEDR and pQCD. In [9], the decision was made to use pQCD in the range 2.6 ≤ √ s ≤ 3.73 GeV, where the quality of inclusive data was poor, with an error inflated according to the percentage errors of the inclusive BES data in this region [64]. With the new KEDR data [62,63], the inclusive data combination is much improved, as shown in Figure 17. In this range, the data combination results in For the larger energy range 1.937 ≤ √ s ≤ 11.2 GeV, the resulting data combination is displayed in Figure 18. As well as the differences observed between the data and pQCD below the charm threshold, the data above it (unchanged since [9]) also show a slight variation from the prediction of pQCD. Considering that with the new, precise KEDR data, the differences between the inclusive data and pQCD are not as large as previously and that this work is aiming at a predominantly data-driven analysis, the contributions in the entire inclusive data region are now estimated using the inclusive data alone. Hence, for this analysis, the contribution from the inclusive data is found to be

Transition region between exclusive and inclusive data
The transition region between the sum of exclusive states and the inclusive R-ratio data is of interest and deserves re-examination. For the sum of exclusive channels, whilst many measurements extend to 2 GeV or beyond, with increasing energy the inclusion of more and more multi-hadronic final states is required to achieve a reliable estimate of the total hadronic cross section. Previously, in [9], the sum of exclusive data was used up to 2 GeV, which defined the transition point between the exclusive sum and the inclusive data combination. In this analysis, the new KEDR data [62] contribute two data points below 2 GeV, extending the lower boundary of the inclusive data down to 1.841 GeV (compared to 2 GeV in [9]) and providing an opportunity to reconsider the previous choice concerning the data input in this region.
From the lower boundary of the KEDR measurement up to 2 GeV, the resulting contributions to a had, LO VP µ from the sum of exclusive states, the inclusive data combination and pQCD are given in Table 1. The integrated values of the inclusive data and pQCD agree within errors. However, the contribution from the sum of exclusive states disagrees with the estimates from both the inclusive data and pQCD, where the sum of exclusive provides a smaller contribution. This is particularly visible in Figure 19, where although the sum of exclusive states agrees with the two inclusive data points below 2 GeV at their respective energies, the combined sum of exclusive states is lower in general. This is largely attributed to the new data for the π + π − π 0 π 0 final state, where Figure 11 shows that these new data result in a clear reduction of the fitted cross section below 2 GeV. 6 Due to this effect, the previous transition point in [9] between the sum of exclusive states and the inclusive data combination at 2 GeV is no longer the preferred choice in this work, where it is clear from Figure 19 that these two different choices for the data input are largely incompatible at this point. A more natural choice for this transition point is now 1.937 GeV, where it can be seen from Figure 19 that all available data choices at this energy are in agreement within errors. This is further substantiated by Table 1, where the value for a had, LO VP µ from the contribution from exclusive states below 1.937 GeV summed 6 Interestingly, as can be seen from Figure 19, the sum of exclusive states is in good agreement with the imprecise and, therefore, unused inclusive hadronic cross section data that exist between 1.43 ≤ √ s ≤ 2.00 GeV. This is in contrast with the findings in the previous analyses [7][8][9], which observed that the inclusive data in this range were lower than the sum of exclusive states.  Figure 19: The energy region between 1.8 ≤ √ s ≤ 2.2 GeV where the fit of inclusive R ratio measurements (light blue band) replaces the sum of exclusive hadronic final states (red band) from 1.937 GeV to 2 GeV. The patterned light blue band and patterned red band show the continuation of the inclusive data combination below 1.937 GeV and the continuation of the exclusive sum above 1.937 GeV, respectively. The recent KEDR data are individually marked and included in the inclusive data fit. The light green band shows the data combination of old inclusive hadronic cross section data that exist between 1.43 ≤ √ s ≤ 2.00 GeV, which were previously discussed in [9] and are not used due to their lack of precision. The estimate from pQCD is included for comparison as a dashed line with an error band which is dominated by the variation of the renormalisation scale µ in the range 1 with the contribution from the inclusive data combination above 1.937 GeV is, within errors, in agreement with the integrated values of all other choices for the data input. Consequently, in this work, this is chosen to be the transition point between the sum of exclusive states and the inclusive R-ratio data.

Total contribution of a had, LO VP
µ and ∆α   (5) had (M 2 Z ) calculated in this analysis. The first column indicates the hadronic final state or individual contribution, the second column gives the respective energy range of the contribution, the third column states the determined value of a had, LO VP µ , the fourth column states the determined value of ∆α (5) had (M 2 Z ) and the last column indicates any new data that have been included since [9]. The last row describes the total contribution obtained from the sum of the individual final states, with the uncertainties added in quadrature.    Full hadronic R ratio π All other states (π + π − π 0 π 0 π 0 ) no η ωηπ 0 ηω π + π − π + π − π + π − (π + π − π 0 π 0 π 0 π 0 ) no η (a) The hadronic R-ratio. Full hadronic R ratio All other states (π + π − π 0 π 0 π 0 ) no η ωηπ 0 ηω π + π − π + π − π + π − (π + π − π 0 π 0 π 0 π 0 ) no η (b) The uncertainty of the hadronic R-ratio. The corresponding result for ∆α (5) where the superscript (5) indicates the contributions from all quark flavours except the top quark, which is added separately. In each case, the errors from the individual channels and sources of uncertainty are added in quadrature to determine the total error. The fractional contributions to the total mean value and uncertainty of both a had, LO VP µ and ∆α (5) had (M 2 Z ) from various energy intervals is shown in Figure 20. Here, the dominance of the energy region below 0.9 GeV for a had, LO VP µ and its uncertainty is clearly evident, with this being predominantly due to the contributions from the π + π − channel. Notably, the pie chart depicting the fractional contributions to the (error) 2 of ∆α (5) had (M 2 Z ) reveals how the uncertainty on this quantity is dominated by the contributions from the radiative correction uncertainties. Mostly, this large error contribution comes from the uncertainty due to possible FSR applied to the combination of inclusive data above 1.937 GeV discussed in Section 2.1.3. This, in particular, highlights the differences in the kernel functions of the respective dispersion integrals for a had, LO VP µ and ∆α (5) had (M 2 Z ), where contributions from higher energies have a larger influence on ∆α (5) had (M 2 Z ) than on a had, LO VP µ . If, instead of a data driven analysis, the region above 1.937 GeV was estimated using pQCD, it would effectively eliminate the impacting radiative correction uncertainties in this region. Figure 21 shows the contributions from all hadronic final states to the hadronic R-ratio and its uncertainty below 1.937 GeV. Here, the individual final states are displayed separately as well as with the resulting total hadronic R-ratio. The full compilation for the hadronic R-ratio is shown in Figure 22. The data vector and corresponding covariance matrix of the hadronic R-ratio in the range m π ≤ √ s ≤ 11.1985 GeV determined in this work are available upon request from the authors.  Figure 23: The normalised difference of the clusters of the π + π − data fit from the KNT18 analysis with respect to those from the HLMNT11 analysis in the range 0.6 ≤ √ s ≤ 0.9 GeV. The width of the yellow band represents the total uncertainties of the clusters of the HLMNT11 analysis. The π + π − cross section is displayed for reference.

Comparison with the HLMNT11 evaluation
To understand further how the changes in the data combination/input have altered the estimate of a had, LO VP µ and its uncertainty, a comparison of the results from this analysis and the previous HLMNT11 evaluation [9] is particularly interesting. Table 3 gives a channel-by-channel comparison of the two works, highlighting the differences in the individual contributions for each channel and the total sum over their respective energy ranges. 7 The largest difference occurs in the π + π − channel, where the mean value in this work is lower by almost 1σ of the HLMNT11 analysis and the uncertainty has reduced by approximately one third. As described in the in-depth discussion of the 2π contribution in Section 3.1, this is largely due to the new, precise and highly correlated radiative return data from KLOE and BESIII and the capability of the new data combination method to utilise the correlations to their full capacity. The global χ 2 min /d.o.f. of the leading and major sub-leading channels in this work are compared to those from the HLMNT11 analysis [9] in Table 4. The reduction of the global χ 2 min /d.o.f. for the π + π − channel further highlights that the data combination for this channel has improved. The energy dependent changes in the resonance region are shown in Figure 23, where it can be seen that, as expected from the comparison of the π + π − results in Table 3, the KNT18 data combination is in good agreement with the HLMNT11 analysis, but sits lower overall. 7 Note that the results for individual contributions to a had, LO VP µ from this work that are listed in Table 3 differ from those given earlier in Section 3 and in Table 2, as for a comparison with HLMNT11 [9], contributions to a had, LO VP µ from exclusive channels are evaluated up to 2 GeV. However, to consistently compare the final results for a had, LO VP µ between the two works, the total KNT18 result given in Table 3 is not determined as the sum of the individual contributions listed above it, but is the final result for a had, LO VP µ calculated in this work using the exclusive channels evaluated up to 1.937 GeV. Summing the KNT18 values listed in Table 3
The K + K − channel shows tension with the HLMNT11 analysis, where the new data in this channel from BaBar [45] and CMD-3 [46] have incurred a large increase in the mean value, whilst also improving the uncertainty despite the small increase in global χ 2 min /d.o.f. This is also the case for the π + π − π + π − channel. Other tensions include the K 0 S K 0 L , ηπ + π − , ηω and ω(→ π 0 γ)π 0 channels, where again, the new, more precise data have resulted in changes outside the quoted HLMNT11 uncertainties. The KK2π channel exhibits a similar change as discussed in Section 3.5. All other channels are in good agreement between the different analyses. It it important to note that this work includes three channels that were not included as part of the HLMNT11 analysis: (ηπ + π − π 0 ) no ω , ηωπ 0 and η → non-purely pionic (npp) KK no φ→KK , where these final states were previously unmeasured by experiment and not estimated through isospin relations. Overall, due to the large reduction in the π + π − channel, it is found that the estimate of a had, LO VP µ has decreased between the HLMNT11 analysis and this work, although this decrease is well within the uncertainty. In total, the uncertainty has been reduced by ∼ 42% with respect to the HLMNT11 analysis.

Comparison with other similar works
The DHMZ group have recently released a new estimate of a had, LO VP µ [78] which, due to a similar data input, is directly comparable with this work and provides insight into how choices with regards to the data combination can affect results. In particular, with the uncertainties of a had, LO VP µ from both the KNT18 and DHMZ17 analyses now less than 0.5%, it is important that these differences are understood in order to quantify the reliability of different approaches and results. In [78], the authors provide a channel-by-channel breakdown of their estimates for the different finals states, which are compared to the respective estimates from this work in Table 5. For the exclusive data channels, the DHMZ group choose to take the contributions from these data up to 1.8 GeV, relying on estimates from pQCD above this (with inflated errors for the pQCD data below the cc threshold). As such, the estimates from this work in Table 5 have been recalculated to mimic the chosen energy regions of the DHMZ analysis and allow for a consistent comparison.
When comparing the total estimate of a had, LO VP µ from the two analyses, the results seem to be in very good agreement. However, as can be seen from Table 5, this masks much larger differences in the estimates from individual channels. The most striking of these is the estimate for the π + π − channel, where there is a tension of slightly more than 1σ between the KNT18 Channel This work (KNT18) DHMZ17 [78] Difference Data based channels ( √ s ≤ 1.8 GeV) π 0 γ (data + ChPT) 4.58 ± 0.10 4.29 ± 0.10 0.29 π + π − (data + ChPT) 503.74 ± 1.96 507.14 ± 2.58 −3.40 π + π − π 0 (data + ChPT) 47 10 . The first column indicates the final state or individual contribution, the second column gives the KNT18 estimate, the third column states the DHMZ17 estimate and the last column gives the difference between the two evaluations. For the final states in this work that have low energy contributions estimated from chiral perturbation theory (see [7]), the contributions from these regions have been added to the contributions from the respective data. and DHMZ17 results. This is unexpected when considering the data input for both analyses are likely to be similar and, therefore, points to marked differences in the way the data are combined. The higher value of the DHMZ17 estimate seems to suggest that their data combination favours the data from the BaBar measurement, with this data set being the only single set that could influence the mean value of the π + π − channel to be as high. This behaviour is similar to the result obtained from combining the π + π − data using only a simple weighted average as discussed in Section 3.1. In turn, this effect is compensated by other major sub-leading final states having larger estimates in this work compared to the DHMZ17 analysis. Specifically, the π + π − π 0 , π + π − π + π − and K 0 S K 0 L estimates are noticeably lower in the DHMZ17 analysis. In addition, there is tension in the region between 1.8 ≤ √ s ≤ 3.7 GeV, where the choice to use data in this region has a higher integrated contribution to a had, LO VP µ than the DHMZ17 estimate from pQCD. This is particularly significant when reconsidering Figure 19, where it was observed that the sum of exclusive states from in the range 1.8 ≤ √ s ≤ 2.0 GeV has a cross section that is lower than the estimate from pQCD. The differences seen in Table 5 above 1.8 GeV are then caused by cross section data below the charm production threshold being higher than pQCD (see Figure 17) and lower than pQCD above it (see Figure 18). It should be noted that the estimate for the ω(→ npp)3π final state from isospin relations, although only a small contribution to a had, LO VP µ , exhibits a significant difference between the two analyses, suggesting a different relation has been used in the DHMZ17 analysis than in this work.
As well as the DHMZ analysis, an updated work by F. Jegerlehner (FJ17) [84] resulted in an estimate of a had, LO VP µ,FJ17 = (688.07±4.14)×10 −10 based on the available e + e − data. Within errors, this result is in agreement with both this work and the DHMZ17 analysis. Interestingly, unlike the comparison with DHMZ17, the two-pion contribution in the energy range 0.316 ≤ √ s ≤ 2.0 GeV is found in the FJ17 analysis to be a π + π − µ,FJ17 = (502.16 ± 2.44) × 10 −10 [85], which is in good agreement with the KNT18 estimate in the same energy range of a π + π − µ,KNT18 = (501.68 ± 1.71) × 10 −10 . However, a more detailed comparison with the estimates of other channels determined in [84,85] is not possible as the FJ17 analysis chooses to estimate certain resonance contributions using available parametrisations [4] instead of using the available data. A comparison of recent and previous evaluations of a had, LO VP µ determined from e + e − → hadrons cross section data is shown in Figure 24, 8 which highlights the agreement between the different works and the improvement in the precision of the respective analyses.

SM prediction of g − 2 of the muon
The Standard Model (SM) prediction of the anomalous magnetic moment of the muon is determined by summing the contributions from all sectors of the SM, such that where a QED µ is the QED contribution, a EW µ is the electro-weak contribution and a had µ are the hadronic contributions due to the hadronic vacuum polarisation and hadronic light-by-light scattering (see equation (1.2)). The QED contributions are known up to and including five-loop accuracy. The five-loop calculation has recently been completed numerically by Kinoshita et al. [86,87]   determined from e + e − → hadrons cross section data. The analyses listed in chronological order are: DEHZ03 [79], HMNT03 [7], DEHZ06 [80], HMNT06 [8], FJ06 [81], DHMZ10 [82], JS11 [83], HLMNT11 [9], FJ17 [84] and DHMZ17 [78]. The prediction from this work is listed as KNT18, which defines the uncertainty band that the other analyses are compared to. where the uncertainties are owing to the uncertainty on the lepton masses, the four-loop contributions, the five-loop contributions and the determination of α using measurements of 87 Rb, respectively. With such a precise determination of a QED µ resulting from a perturbative series that converges extremely well, the QED result seems stable. It should be noted, however, that the four-loop and five-loop contributions rely heavily on numerical integrations and independent checks of these results are crucial. This has been recently accomplished through several different analyses [88][89][90][91][92][93][94], which corroborate the results from Kinoshita and collaborators. Therefore, it is safe to assume that the estimate for the QED contribution is well under control.
The contribution from the EW sector is well known to two-loop accuracy [95][96][97][98][99]. With the mass of the Higgs now known, the updated estimate [100] gives where the knowledge of the Higgs mass has halved the uncertainty of this contribution compared to the estimate used in [9]. Although a relatively small contribution when compared to a QED µ , the uncertainty is not negligible considering the projected experimental accuracy, but is small when compared to the hadronic uncertainties. However, with this contribution known safely to two-loop accuracy, the electroweak estimate is also very well under control.
For the hadronic vacuum polarisation contributions, the leading order and next-to-leading order contributions have been calculated in this work. The LO contribution, from equation (3.28), was found to be a had, LO VP µ = (693.26 ± 2.46) × 10 −10 and the NLO was given in equation (3.29) as a had, NLO VP µ = (−9.82 ± 0.04) × 10 −10 . The calculation of the NNLO hadronic vacuum polarisation contribution has been achieved for the first time in [13] (see also the evaluation in [84]), where contributions from five individual classes of diagrams (each with an independent, respective kernel function) were determined. It was found to be a had, NNLO VP µ = (1.24 ± 0.01) × 10 −10 . This estimate is slightly larger than expected and results in a slight increase to the mean value of a had,VP µ . Summing these, the total contribution to the anomalous magnetic moment from the hadronic vacuum polarisation is a had, VP µ = (684.68 ± 2.42) × 10 −10 . (3.34) It should be noted that the negative NLO contribution results in an anti-correlation between its uncertainty and the uncertainty from the LO contribution, consequently resulting in a slight reduction in the overall uncertainty that has been incorporated into equation (3.34). The hadronic LbL contributions, although small compared to the hadronic vacuum polarisation sector, have, in the past, been determined through model-dependent approaches. These are based on meson exchanges, the large N c limit, ChPT estimates, short distance constraints from the operator product expansion and pQCD. Over time, several different approaches to evaluating a had, LbL µ have been attempted, resulting in good agreement for the leading N c (π 0 exchange) contribution, but differing for sub-leading effects. A commonly quoted determination of the LbL contribution is the 'Glasgow consensus' estimate of a had, LbL µ ('Glasgow consensus') = (10.5 ± 2.6) × 10 −10 [101] (alternatively, see [102][103][104][105]). However, recent works [106][107][108] have re-evaluated the contribution to a had,LbL µ due to axial exchanges, where it has been found that this contribution has, in the past, been overestimated due to an incorrect assumption that the form factors for the axial meson contribution are symmetric under the exchange of two photon momenta [106]. Under this assumption, the determination in [102] previously found the axial vector contribution to be a had, LbL; axial which is the chosen estimate for a had, LbL µ in this work. This result is notably lower than the previously accepted LbL estimates and will incur an overall downward shift on a SM µ . It is, however, still within the original uncertainties when comparing with the original 'Glasgow consensus' estimate. Alternatively, it should be noted that the estimate of a had, LbL µ = (10.2 ± 3.9) × 10 −10 [108,109], which is a result that is independent of the 'Glasgow consensus' estimate, could be employed here. In addition, the recent work [105] has provided an estimate for the next-to-leading order hadronic LbL contribution. It has found a had,NLO−LbL µ = (0.3±0.2)×10 −10 , which does not alter the hadronic LbL contribution significantly, but is taken into account in the full SM prediction given below.
Much work has also been directed at the possibility of a model independent calculation of a had,LbL µ to further consolidate the SM prediction of a µ . One approach involves the measurement of transition form factors by KLOE-2 and BESIII, which can be expected to constrain the   [82], JS11 [83], HLMNT11 [9], FJ17 [84] and DHMZ17 [78]. The prediction from this work is listed as KNT18, which defines the uncertainty band that other analyses are compared to. The current uncertainty on the experimental measurement [1][2][3][4] is given by the light blue band. The light grey band represents the hypothetical situation of the new experimental measurement at Fermilab yielding the same mean value for a exp µ as the BNL measurement, but achieving the projected four-fold improvement in its uncertainty [5].
leading pseudoscalar-pole (π 0 , η, η ) contribution to a precision of approximately 15% [108]. Alternatively, the pion transition form factor (π 0 → γ * γ * ) can be calculated on the lattice for the same purpose [110]. New efforts into the prospects of determining a had,LbL µ using dispersive approaches are also very promising [111][112][113][114][115][116], where the dispersion relations are formulated to calculate either the general hadronic LbL tensor or to calculate a had,LbL µ directly. These approaches will allow for the determination of the hadronic LbL contributions from experimental data and, at the very least, will invoke stringent constraints on future estimates. Lastly, there has been huge progress in developing methods for a direct lattice simulation of a had,LbL µ [110,[117][118][119][120][121][122][123]. With a proof of principle already well established, an estimate of approximately 10% accuracy seems possible in the near future. Considering these developments and the efforts of the Muon g − 2 Theory Initiative [124] to promote the collaborative work of many different groups, the determination of a had,LbL µ on the level of the 'Glasgow consensus' will, at the very least, be consolidated and a reduction of the uncertainty seems highly probable on the timescales of the new g − 2 experiments. where the uncertainty is determined from the uncertainties of the individual SM contributions Analysis ∆α (5) had (M 2 Z ) × 10 4 α −1 (M 2 Z ) DHMZ10 [82] 275.59 ± 1.04 128.952 ± 0.014 HLMNT11 [9] 276.26 ± 1. 38 128.944 ± 0.019 FJ17 [129] 277.38 ± 1. 19 128.919 ± 0.022 DHMZ17 [78] 276.00 ± 0. 94 128.947 ± 0.012 KNT18 [This work] 276.11 ± 1.11 128.946 ± 0.015 Table 6: Comparison of recent and previous evaluations of ∆α (5) had (M 2 Z ) determined from e + e − → hadrons cross section data and the corresponding results for α −1 (M 2 Z ).  Figure 25. In particular, a comparison with the HLMNT11 estimate given in equation (1.7) shows an improvement in the total uncertainty of a SM µ of ∼ 27%. It should be noted that although, as stated above, the DHMZ17 estimate for a had, LO VP µ [78] is lower than the value obtained in this work, the estimate of a SM µ from DHMZ17 is higher than the estimate from this analysis as DHMZ17 choose to use the estimate for the hadronic light-by-light contribution of a had,LbL µ = (10.5 ± 2.6) × 10 −10 [101].

Determination of α(M 2 Z )
From equation (3.30), the five-flavour hadronic contribution to ∆α(M 2 Z ) is found to be ∆α  [127,128], the total value of the QED coupling at the Z boson mass is found in this work to be Here, the top contribution is determined using the top quark mass, m t = 173.1 ± 0.6 GeV and the value for strong coupling constant at the mass of the Z boson, α s (M Z ) = 0.1182(12) [4]. A comparison of these results with other determinations of ∆α (5) had (M 2 Z ) and α −1 (M 2 Z ) is given in Table 6. Note that, as discussed in Section 3.7, the different weighting of the kernel function of the ∆α (5) had (M 2 Z ) dispersion integral compared to the integral for a had, LO VP µ results in contributions from higher energy regions having a larger influence on ∆α (5) had (M 2 Z ) than on a had, LO VP µ . Therefore, the choice to use either the available inclusive data or pQCD above ∼ 2 GeV as discussed in Section 3.9 can have marked differences on the values and errors obtained for ∆α (5) had (M 2 Z ).

Conclusions and future prospects for the muon g − 2
This analysis, KNT18, has completed a full re-evaluation of the hadronic vacuum polarisation contributions to the anomalous magnetic moment of the muon, a had, VP µ , and the hadronic contributions to the effective QED coupling at the Z boson mass, ∆α had (M 2 Z ). These quantities have been determined using the available e + e − → hadrons cross section data as input into corresponding dispersion relations, with an aim to achieve both accurate and reliable results from a predominantly data driven analysis. Since [9], all aspects concerning the radiative corrections of the data and the data combination have been reassessed in this work. Specifically, the data are now combined using an iterative, linear χ 2 -minimisation developed from a method that has been advocated to be free of bias and that has been studied in detail. Importantly, this data combination method allows for the full use of any available correlated statistical and systematic uncertainties, incorporating experimental covariance matrices in the combination in a bias-free approach. These changes, plus the large quantity of new hadronic cross section data, have resulted in improved estimates for nearly all hadronic channels. This is particularly true for the π + π − channel, where the precision of this final state has improved by approximately one third compared to [9], with a π + π − µ from both analyses in very good agreement. Importantly, the inclusion of recently released neutral data in the KKπ and KK2π channels has removed the need to rely on isospin relations to estimate these final states. In addition, new inclusive hadronic R-ratio data from the KEDR collaboration have improved the inclusive data combination. In particular, they have provided the opportunity to reconsider the transition region between the sum of exclusive states and the inclusive data, which has resulted in the transition point being chosen to be at 1.937 GeV in this work, where the different choices for the data input in this point are in agreement within errors. The complete data compilation has yielded the full hadronic R-ratio and its covariance matrix in the energy range m π ≤ √ s ≤ 11.2 GeV. Using these combined data, for the muon g −2, this analysis found a had, LO VP µ = (693.26±2.46)×10 −10 and a had, NLO VP µ = (−9.82 ± 0.04) × 10 −10 . This has resulted in a new estimate for the Standard Model prediction of a SM µ = (11 659 182.04 ± 3.56) × 10 −10 , which deviates from the current experimental measurement by 3.7σ. For the effective QED coupling, the new estimate in this work of ∆α (5) had (M 2 Z ) = (276.11 ± 1.11) × 10 −4 yields an updated value for the QED coupling at the Z boson mass of α −1 (M 2 Z ) = 128.946 ± 0.015. In general, the predictions for a SM µ have been further scrutinised and are well established. In particular, the improvements in the uncertainty on the level discussed in [130] are on track in preparation for the new experimental results from Fermilab and J-PARC. The efforts of the Muon g − 2 Theory Initiative and the groups involved within it show great progress and promise in improving the estimate of a SM µ further. Importantly, it should be noted that there is no indication thus far that the SM prediction does not deviate with the current experimental measurement by more than 3σ. For the hadronic vacuum polarisation contributions, there is scope to further improve the estimates, as new measurements of the π + π − cross section are planned to be released in the near future from BaBar, CMD-3, SND and, possibly, BELLE-2. Also, efforts are currently being made to measure new inclusive R-ratio data by BESIII and the experiments at Novosibirsk (SND, CMD-3, KEDR). This, in particular, will benefit the evaluation of ∆α (5) had (M 2 Z ) which, as discussed in Section 3.11, is more sensitive than a had, LO VP µ to the precision of the data from higher energy regions. In addition, lattice determinations of a had, VP µ are rapidly improving. Recent work that combines data from lattice QCD with those from experimental R-ratio data have already provided extremely accurate results that are in good agreement with the current estimates from the dispersive method [131]. Furthermore, efforts to experimentally measure a had, LO VP µ in the space-like region are progressing [132,133] and would provide an alternative check of the results from the dispersive approach. Should any or all these advances reduce the uncertainty of a had, VP µ even further, the improvement of the hadronic light-by-light estimates will become particularly crucial. This is further driven by the fact that the result for a had, VP µ determined in this analysis is the first estimate of the hadronic vacuum polarisation contribution that is more precise than the currently quoted uncertainties for a had, LbL µ . However, given these developments in improving the Standard Model prediction of a µ and the formidable progress made by the new muon g − 2 experiments at Fermilab and J-PARC, the prospects to either establish the existence of new physics contributing to a µ , or to rule out many new physics scenarios, are highly compelling.