On flavourful Easter eggs for New Physics hunger and lepton flavour universality violation.

Within the standard approach of effective field theory of weak interactions for Δ B = 1 transitions, we look for possibly unexpected subtle New Physics effects, here dubbed "flavourful Easter eggs". We perform a Bayesian global fit using the publicly available HEPfit package, taking into account state-of-the-art experimental information concerning these processes, including the suggestive measurements from LHCb of R K and R K ∗ , the latter available only very recently. We parametrise New Physics contributions to b → s transitions in terms of shifts of Wilson coefficients of the electromagnetic dipole and semileptonic operators, assuming CP-conserving effects, but allowing in general for violation of lepton flavour universality. We show how optimistic/conservative hadronic estimates can impact quantitatively the size of New Physics extracted from the fit. With a conservative approach to hadronic uncertainties we find nonzero New Physics contributions to Wilson coefficients at the level of ∼ 3 σ , depending on the model chosen. Furthermore, given the interplay between hadronic contributions and New Physics effects in the leptonic vector current, a scenario with nonstandard leptonic axial currents is comparable to the more widely advocated one with New Physics in the leptonic vector current.

Introduction

Easter eggs nowadays also refer to inside jokes and/or secret messages usually hidden e.g. in computer gaming and hi-tech software. In this work, we take advantage of this terminology to motivate the search for New Physics Beyond the Standard Model in the radiative and in the (semi)leptonic channels of rare B meson decays.

In the decades that have followed the original formulation of flavour mixing [1], the flavour structure of the SM has been experimentally tested and well established. The tremendous progress of the experimental facilities has probed the flavour of the SM to an exquisite level of precision [2], along with the substantial effort on the part of the theoretical community to go well beyond leading order computations [3]. From this perspective of “precision tests”, radiative and (semi)leptonic processes, related at the partonic level to transitions, occupy a special place in probing the SM and its possible extensions in terms of New Physics (NP) models [4, 5].

Firstly, these rare B meson decays belong to the class of flavour-changing neutral current (FCNC) processes, which are well known to be sensitive probes of Physics Beyond the Standard Model (BSM): in fact – within the SM – the flavour structure of the theory allows FCNC to arise only at loop level, as a consequence of the GIM mechanism [6]. This allows for significant room for heavy new degrees of freedom to sizeably contribute to these rare processes.

Secondly, from the experimental side, the study of rare B meson decays offers us some of the most precise measurements amongst the processes. For instance, the measurement of the inclusive branching fraction of is currently performed with a relative uncertainty of a few percent [7-9], while the study of an exclusive mode such as allows for a detailed analysis of the angular distribution of the four final state particles, yielding rich experimental information in terms of angular functions of the dilepton invariant mass, with full kinematic coverage of the latter [10] and – starting from Ref. [11] – also with available experimental correlations among the angular observables.

In B Physics, the recent years have been characterised by the emergence of a striking pattern of anomalies in multiple independent studies of some of these rare transitions [12]. Of particular importance, the measurement of the angular observable [13-16] stands out from all the other ones related to the angular distribution of ; first realised by the LHCb collaboration [17, 18] and later on also by the Belle collaboration [19], the experimental analysis of in the large recoil region of the decay points to a deviation of about with respect to the SM prediction presented in Ref. [20]. The latter, however, suffers from possible hadronic uncertainties which are sometimes even hard to guesstimate [21-24], and this observation has been at the origin of a quite vivid debate in the recent literature about the size of (possibly) known and (yet) unknown QCD power corrections to the amplitude of this process in the infinite mass limit [25-28]. To corroborate even more the cumbersome picture of the “ anomaly”, two new independent measurements of this angular observable (among others) have been recently released by the ATLAS [29] and CMS [30] collaborations, showing, respectively, an appreciable increase and reduction of the tension between data and the SM prediction in Ref. [20], as reported by these experiments.

For the sake of completeness, one should also remark that other smaller tensions have been around, concerning the measurement of differential branching fractions of [31, 32] and [33]. It is worth noting that, while for the latter mode an explanation in terms of hadronic physics may easily be conceivable, the theoretical computation of the former seems to be under control [34].

Quite surprisingly, a possible smoking gun for NP in rare B meson decays already came out in 2014, when the LHCb collaboration presented for the first time the measurement of the ratio of branching fractions [35]:where the subscript refers to the dilepton mass (denoted hereafter ) range going from 1 to 6 GeV. This experimental value shows a deviation of about with respect to the standard theoretical prediction. Indeed, the SM value of in the bin provided by the LHCb collaboration is expected to be equal to unity beyond the percent level of accuracy [36, 37]. In fact, contrary to observables such as , it is important to stress that may be, in general, regarded as insensitive to QCD effects [36]. From the model building point of view, can certainly be considered as quite informative, hinting at a UV completion of the SM where lepton flavour universality violation (LFUV) takes place in the flavour-violating couplings of new heavy degrees of freedom, e.g. leptoquarks and/or gauge bosons [38-67]. Most importantly, the tantalising correlation of this signature of LFUV with the anomaly, suggested by several global analyses [4, 68–72] has triggered different proposals of measurements of such effect in the angular analysis of the channel [73, 74]. Interestingly enough, an analysis from the Belle collaboration aiming at separating the leptonic flavours in [75] shows a consistent deviation from the SM prediction reported in Ref. [20] in the dimuon leptonic final state only. This is compatible with previous experimental findings related only to the mode with muonic final states.

Sitting on similar theoretical grounds to , another intriguing ratio of B decay branching fractions can be measured in the channel: These measurements for the low- bin and the central- one have just been presented by the LHCb collaboration [76], pointing again to a discrepancy of about with respect to the expected SM prediction – again equal to 1 to a very good accuracy for the central- bin and close to 0.9 for the low- one – and yielding more than a deviation when naively combined with the measurement of . Note that with higher degree of braveness (or, depending on the taste of the reader, of unconsciousness), the disagreement of the SM with precision B physics may reach the exciting level of when one naively combines together the single significances coming from ratios, measurements and the minor deviations observed in the other exclusive branching fractions.

Given the excitement of these days for all the above hints of a possible NP discovery in rare B meson decays, in this work we take our first steps towards a positive attitude in the search of a definite BSM pattern aimed at addressing these B anomalies. We perform our study in a model-independent fashion, within the framework of effective field theories for weak interactions [77-79]. In particular, in Sect. 2 we define the setup characterising the whole global analysis, presenting six different benchmark scenarios for NP, together with a discussion as regards two different approaches in the estimate of the hadronic uncertainties that can affect quantitatively our final results. In Sect. 3, we list all the experimental measurements we use to construct the likelihood in our fit, and we discuss in detail our most important findings. The latter are effectively depicted in Figs. 1, 2, 3, 4, 5, 6, and collected in Tables 2, 3, 4, 5 in Appendix A. In Sect. 4 we summarise our conclusions.

Fig. 1

The two NP parameter fit using and . Here and in the following, the left green panel shows the results for the PMD approach and the right red panel shows that for the PDD one. In the 1D distributions we show the 16th, 50th and 84th percentile marked with the dashed lines. In the correlation plots we show the 1, 2 and contours in decreasing degrees of transparency. The blue square and lines identify the values of the NP WCs in the SM limit. The numbers at the bottom left corner of the 2D plots refer to the correlation. We also report IC values for the two approaches (see Eq. (9)). Preferred models are expected to give smaller IC values

Fig. 2

The two NP parameter fit using and . See caption of Fig. 1 for the colour coding and further details

Fig. 3

The three NP parameter fit using , and . See caption of Fig. 1 for the colour coding and further details

Fig. 4

The three NP parameter fit using , and . See caption of Fig. 1 for the colour coding and further details

Fig. 5

The three NP parameter fit using , , and . See caption of Fig. 1 for the colour coding and further details

Fig. 6

The five NP parameter fit using , , , and . See caption of Fig. 1 for the colour coding and further details

Table 2

Results from the fit for WCs and hadronic contributions in the PMD approach. See Sect. 2.1 for details of the six NP scenarios

Par.	(I)	(II)	(III)	(IV)	(V)	(VI)
C7NP	−	−	0.015±0.014	-0.011±0.013	0.003±0.013	0.015±0.014
C9,μNP	-1.58±0.28	-1.53±0.25	-1.66±0.29	−	-0.54±0.17	-1.64±0.29
C9,eNP	-0.10±0.45	−	-0.18±0.46	−	0.09±0.25	-1.6±1.0
C10,μNP	−	0.03±0.16	−	-0.12±0.22	0.54±0.17	0.009±0.200
C10,eNP	−	−	−	-1.22±0.37	-0.09±0.25	-0.91±0.76
|h0(0)|·10^4	2.1±1.2	2.0±1.2	2.2±1.3	1.8±1.2	1.3±1.0	2.0±1.3
|h+(0)|·10^4	0.079±0.067	0.079±0.067	0.076±0.065	0.083±0.069	0.086±0.072	0.076±0.064
|h-(0)|·10^4	0.53±0.19	0.54±0.19	0.52±0.19	0.56±0.20	0.60±0.21	0.52±0.19
|h0(1)|·10^4	0.30±0.23	0.30±0.22	0.30±0.23	0.45±0.26	0.32±0.24	0.28±0.22
|h+(1)|·10^4	0.22±0.20	0.22±0.19	0.22±0.19	0.21±0.19	0.26±0.22	0.22±0.19
|h-(1)|·10^4	0.23±0.19	0.23±0.19	0.23±0.20	0.30±0.21	0.32±0.22	0.23±0.19
|h+(2)|·10^4	0.052±0.045	0.053±0.045	0.052±0.044	0.046±0.042	0.064±0.053	0.050±0.044
|h-(2)|·10^4	0.046±0.038	0.046±0.039	0.046±0.039	0.092±0.050	0.070±0.047	0.045±0.038

Table 3

Results from the fit for WCs and hadronic contributions in the PDD approach. See Sect. 2.1 for details of the six NP scenarios

Par.	(I)	(II)	(III)	(IV)	(V)	(VI)
C7NP	−	−	0.013±0.014	0.008±0.014	0.011±0.014	0.014±0.014
C9,μNP	-1.47±0.63	-1.17±0.46	-1.58±0.64	−	-0.43±0.23	-1.43±0.64
C9,eNP	0.007±0.620	−	-0.08±0.63	−	0.21±0.29	-1.2±1.2
C10,μNP	−	0.26±0.23	−	0.27±0.26	0.43±0.23	0.20±0.25
C10,eNP	−	−	−	-0.86±0.4	-0.21±0.29	-0.60±0.99
|h0(0)|·10^4	2.6±1.6	2.3±1.4	2.6±1.6	1.7±1.3	1.7±1.3	2.6±1.6
|h+(0)|·10^4	0.075±0.066	0.081±0.070	0.077±0.067	0.086±0.075	0.087±0.075	0.077±0.067
|h-(0)|·10^4	0.52±0.21	0.55±0.22	0.52±0.21	0.60±0.23	0.59±0.23	0.53±0.21
|h0(1)|·10^4	0.40±0.32	0.41±0.34	0.39±0.32	0.50±0.36	0.46±0.37	0.40±0.33
|h+(1)|·10^4	0.40±0.29	0.42±0.30	0.40±0.29	0.39±0.29	0.42±0.30	0.41±0.30
|h-(1)|·10^4	0.47±0.35	0.52±0.38	0.48±0.36	0.82±0.46	0.73±0.43	0.50±0.37
|h+(2)|·10^4	0.138±0.087	0.160±0.099	0.131±0.086	0.139±0.094	0.160±0.100	0.145±0.095
|h-(2)|·10^4	0.112±0.085	0.126±0.098	0.111±0.083	0.190±0.100	0.170±0.110	0.124±0.094

Table 4

Experimental results (with symmetrised errors) and results from the fit for key observables in the PMD approach. See Sect. 2.1 for details of the six NP scenarios

Obs.	Exp. value	(I)	(II)	(III)	(IV)	(V)	(VI)
R_K_[1,6]	0.753±0.090	0.722±0.067	0.703±0.047	0.722±0.067	0.781±0.055	0.740±0.061	0.724±0.067
R_K^∗_[0.045,1.1]	0.680±0.093	0.885±0.016	0.881±0.016	0.885±0.016	0.839±0.024	0.858±0.019	0.843±0.030
R_K^∗_[1.1,6]	0.707±0.102	0.803±0.057	0.786±0.049	0.802±0.057	0.713±0.065	0.740±0.060	0.717±0.067
RK^∗L_[1.1,6]	−	0.724±0.067	0.705±0.049	0.723±0.067	0.722±0.065	0.722±0.063	0.723±0.068
RK^∗T_[1.1,6]	−	1.053±0.030	1.046±0.069	1.050±0.030	0.692±0.067	0.794±0.053	0.75±0.25
R_ϕ_[1.1,6]	−	0.782±0.060	0.764±0.048	0.781±0.060	0.720±0.064	0.737±0.061	0.717±0.062
RϕL_[1.1,6]	−	0.723±0.067	0.704±0.049	0.722±0.067	0.728±0.064	0.723±0.063	0.721±0.067
RϕT_[1.1,6]	−	1.044±0.032	1.038±0.066	1.043±0.033	0.691±0.067	0.794±0.053	0.75±0.24
P_5[4,6]LHCb	-0.301±0.160	-0.423±0.060	-0.427±0.061	-0.420±0.059	-0.556±0.063	-0.587±0.053	-0.431±0.061
P_5[6,8]LHCb	-0.505±0.124	-0.602±0.059	-0.607±0.061	-0.598±0.059	-0.677±0.066	-0.704±0.057	-0.61±0.06
P_5[4,6]ATLAS	-0.26±0.39	-0.423±0.060	-0.427±0.061	-0.420±0.059	-0.556±0.063	-0.587±0.053	-0.431±0.061
P_5[4.3,6]CMS	-0.955±0.268	-0.442±0.059	-0.447±0.061	-0.440±0.059	-0.572±0.063	-0.603±0.053	-0.451±0.061
P_5[6,8.68]CMS	-0.660±0.220	-0.623±0.060	-0.626±0.061	-0.618±0.059	-0.688±0.067	-0.711±0.058	-0.63±0.06
P_5[4,8]Belle	-0.025±0.318	-0.520±0.057	-0.524±0.059	-0.516±0.057	-0.620±0.063	-0.649±0.053	-0.528±0.059
P_5,e[4,8]Belle	-0.510±0.272	-0.782±0.054	-0.794±0.039	-0.782±0.054	-0.536±0.063	-0.710±0.044	-0.42±0.23

Table 5

Experimental results (with symmetrised errors) and results from the fit for key observables in the PDD approach. See Sect. 2.1 for details of the six NP scenarios

Obs.	Exp. value	(I)	(II)	(III)	(IV)	(V)	(VI)
R_K_[1,6]	0.753±0.090	0.724±0.068	0.714±0.064	0.723±0.067	0.775±0.057	0.740±0.062	0.719±0.073
R_K^∗_[0.045,1.1]	0.680±0.093	0.883±0.017	0.873±0.016	0.883±0.017	0.842±0.024	0.861±0.019	0.847±0.030
R_K^∗_[1.1,6]	0.707±0.102	0.803±0.059	0.777±0.053	0.801±0.059	0.714±0.066	0.755±0.059	0.722±0.069
RK^∗L_[1.1,6]	−	0.726±0.070	0.719±0.061	0.723±0.070	0.726±0.065	0.729±0.063	0.719±0.073
RK^∗T_[1.1,6]	−	1.050±0.033	0.959±0.083	1.046±0.034	0.688±0.069	0.824±0.049	0.8±0.32
R_ϕ_[1.1,6]	−	0.782±0.062	0.761±0.054	0.779±0.062	0.722±0.065	0.749±0.060	0.718±0.063
RϕL_[1.1,6]	−	0.724±0.070	0.716±0.061	0.721±0.070	0.731±0.064	0.728±0.063	0.716±0.071
RϕT_[1.1,6]	−	1.039±0.036	0.955±0.079	1.036±0.036	0.693±0.068	0.828±0.049	0.79±0.3
P_5[4,6]LHCb	-0.301±0.160	-0.398±0.066	-0.405±0.068	-0.401±0.066	-0.435±0.065	-0.421±0.068	-0.417±0.068
P_5[6,8]LHCb	-0.505±0.124	-0.536±0.080	-0.521±0.083	-0.539±0.080	-0.521±0.080	-0.508±0.082	-0.545±0.083
P_5[4,6]ATLAS	-0.26±0.39	-0.398±0.066	-0.405±0.068	-0.401±0.066	-0.435±0.065	-0.421±0.068	-0.417±0.068
P_5[4.3,6]CMS	-0.955±0.268	-0.414±0.067	-0.420±0.068	-0.417±0.067	-0.447±0.066	-0.434±0.068	-0.433±0.068
P_5[6,8.68]CMS	-0.660±0.220	-0.550±0.083	-0.531±0.086	-0.553±0.083	-0.529±0.084	-0.514±0.085	-0.556±0.086
P_5[4,8]Belle	-0.025±0.318	-0.472±0.070	-0.467±0.072	-0.475±0.070	-0.481±0.070	-0.467±0.072	-0.486±0.072
P_5,e[4,8]Belle	-0.510±0.272	-0.725±0.075	-0.664±0.095	-0.733±0.075	-0.424±0.063	-0.570±0.066	-0.41±0.22

Theoretical framework of the analysis

In this section we present the effective field theory framework at the basis of this work and introduce the benchmark scenarios we focus on for our study of NP effects in rare B decays. We then illustrate the two distinct broad classes of assumptions that characterise our global analysis: the case where we take an optimistic attitude towards the estimate of hadronic uncertainty plaguing the amplitude of both and channels, and a second one where we aim at providing a more conservative approach. All the results in Sect. 3.2 will be classified under these two different setups.

New Physics benchmarks for

Integrating out the heavy degrees of freedom, the resulting effective Hamiltonian of weak interactions for transitions involves the following set of dimension six operators within the SM [80]:where , and we have neglected the chirally suppressed SM dipoles. The effective Hamiltonian can be cast in full generality in the form of a combination of two distinct parts:where, within the SM, the hadronic term involves the first seven operators in Eq. (5):while the second piece includes the electromagnetic dipole and semileptonic operators:with corresponding to the CKM combination for and where are the Wilson coefficients (WCs) encoding the short-distance physics of the theory. All the SM WCs in this work are evolved from the mass scale of the W boson down to  GeV, using state-of-the-art perturbative QCD and QED calculations for the matching conditions [81-83] and the anomalous dimension matrices [83-86] relevant for the processes considered in this analysis.

While a general UV completion of the SM may enter in the effective couplings present in both pieces of Eq. (5), general NP effects in can be phenomenologically parametrised as shifts of the Wilson coefficients of the electromagnetic and semileptonic operators at the typical scale of the processes, . In particular, the most general basis for NP effects in radiative and (semi)leptonic B decays can be enlarged by the presence of scalar, pseudo-scalar and tensorial semileptonic operators, together with right-handed quark currents as the analogue of SM operators [21, 87]. In this work, motivated by previous interesting findings concerning LFUV [69-71] and the measurement of and , we focus on the contributions of NP appearing as shifts of the SM WCs related to the electromagnetic dipole and semileptonic operators with left-handed quark currents only. A comprehensive analysis with different chiral structures as well as a more general effective theory framework will be presented elsewhere [88]. Furthermore, we restrict ourselves to CP-conserving effects, taking NP WCs to be real.

For NP in semileptonic operators we discriminate between couplings to muon and electron fields both in the axial and vector leptonic currents. We characterise our phenomenological analysis for NP through six different benchmark scenarios, studying the impact of combinations of the following NP WCs: We remark that while benchmarks (I) and (II) have been already studied in the literature, none of the other cases has been analysed so far. In particular, NP scenarios (III) and (IV) allow us to study, for the first time, the interesting impact of a NP radiative dipole operator in combination with vector-like and axial-like LFUV effects generated by NP. Most interestingly, scenario (V) allows us to explore the correlation , possibly hinting at a preserving BSM theory. As an additional interesting case to explore, we eventually generalise to simultaneously nonvanishing , , , and in case (VI).

and varied in the range , i.e. adding to the SM two NP parameters;

and varied in the range , adding to the SM again two NP parameters;

and varied in the range , and varied in the range , i.e. a scenario with three NP parameters;

and varied in the range , and varied in the range , i.e. adding again to the SM three NP parameters;

and varied in the range , and varied in the range , i.e. a NP scenario again described by three different parameters.

, , , and varied simultaneously in the respective ranges defined above, i.e. a NP scenario described by five different parameters.

We wish to stress that all of the six benchmarks defined above will be studied for the first time under two different approaches in the estimate of QCD hadronic power corrections, as presented in next section.

Treatment of the hadronic uncertainties

In our previous work [24, 27, 89], we went into considerable detail on the treatment of hadronic contributions in the angular analysis of . Our approach there was to study how large these contributions can be assuming that the LHCb data on branching fractions and angular distributions of these decay modes could be described within the SM. For that purpose we considered four scenarios for the hadronic contributions, with increasing theoretical input from the phenomenological analysis presented in Ref. [90]. The underlying functional form that we used for the hadronic contribution was given bywhere we fitted for the complex, helicity dependent, coefficients , and using the data and the phenomenological model in [90]. Since enters the decay amplitude with an additional factor of with respect to , we drop in our analysis.

In this work we proceed to study the possible existence of NP contributions in semileptonic and radiative decays which requires a re-evaluation of the hadronic uncertainties. For the sake of simplicity, to address hadronic contributions we use the same functional parameterisation as given in Eq. (8). However, we limit ourselves to only two hadronic models. The first, corresponding to the most widely used assumption, relies completely on the phenomenological model in [90] below . The second is a more conservative approach, where we impose the latter only in the large recoil region at GeV while letting the data drive the hadronic contributions in the higher invariant mass region. We will refer to the first approach as phenomenological model driven (PMD) and the second as phenomenologically and data driven (PDD). In our fit we vary the parameters over generous ranges. More detailed discussion of these can be found in [24, 27].

In the present analysis we also need to address modes that were not considered in our previous work, namely , and . The decay has been studied in detail in [34], where the authors show that the hadronic uncertainties are smaller than in . A comparison of the LCSR estimate of the soft gluon contribution and the QCDF estimate of the hard gluon contribution reveals that the soft gluon exchange is subdominant with respect to QCDF hard gluon exchange. Therefore, although in principle the same concerns on the soft gluon contribution we raised for apply in this case, in practice the overall effect of soft gluons can be reasonably neglected. In our computation we therefore only include hard gluon exchange computed using the QCDF formalism in Ref. [91].

The long distance contributions for and follow a similar theoretical derivation to those for and , respectively, barring the fact that the spectator quark in the former is different from that in the latter. No theoretical estimates of power corrections to the infinite mass limit are available for the decays and one has to rely on the ones for the decays to get a handle on the long distance contributions. The spectator quark effects can come through the hard spectator scattering involving matrix elements of , and computable in QCD factorisation [91] which we include in our computation. However, we do not include the sub-leading, and numerically small, QCDF power corrections to spectator scattering involving  [92-94] and contributions to weak spectator scattering involving beyond QCDF computed in LCSR [95-97]. The effect of the difference in all these spectator contributions is expected to be low firstly because they are numerically small and, secondly, because the effect is proportional to the small flavour SU(3) breaking. Different approaches in relating the long distance contributions in the channels to the ones in the channels have been used in the literature [69, 70, 98], which vary in the degree of correlation between the two. While Ref. [70] uses uncorrelated hadronic uncertainties, Refs. [69, 98] have left the two contributions highly correlated noting that the spectator contribution is expected to be numerically small. We take an approach similar to the latter considering the insensitivity of the current data to such effects and use the same value of power corrections in and amplitudes, even though this choice pertains to a quite oversimplifying optimistic attitude. We leave a more detailed analysis of this assumption by relaxing the correlation between the hadronic contributions in the two modes to a future work [88].

Bayesian fit of the dipole and semileptonic operators

Experimental information considered

In this section we discuss the experimental measurements we use in our fit. Note that for the exclusive modes we make use of measurements in the large recoil region only. Our choice harbours on the fact that the QCD long distance effects in the low recoil region are substantially different from the large recoil regime [99-102] and would require a dedicated analysis. For the fit in this study we consider the following experimental information:

For the channel we use the LHCb measurements of CP-averaged angular observables extracted by means of the unbinned maximum likelihood fit, along with the provided correlation matrix [18]. Moreover, we employ the recent results for CP-averaged angular observables from ATLAS [29] and the ones measured by CMS [30, 103].1 Finally, we use the CP-averaged optimised angular observables recently measured by Belle [75]2. Regarding the differential branching fractions, we use the recently updated measurements from LHCb [104] and the ones from CMS [103]. For the channel we consider the LHCb results from [105] and the Belle results from [75]. observable is considered according to the recently presented measurements by LHCb [76] in both the low- and the central- bins; see also Eq. (3). Our theoretical predictions are computed in the helicity basis, whose relevant expressions can be found in [21]; the same framework is employed to study , , and channels. For the latter, we use the full set of form factors extrapolated from the lattice results, along with the provided correlation matrix [106]; for the remaining channels, we use the full set of form factors estimated combining LCSR and lattice results, along with the correlation matrices [107]. For the factorisable and non-factorisable QCD power corrections, we refer to Sect. 2.2.

We include in our analysis the HFAG average for the branching fractions from [2].

We consider the LHCb CP-averaged angular observables and differential branching fractions measurements, along with the provided correlation matrix [33].

We use the LHCb measurement of the branching fraction from [108].

We employ the LHCb measurement of differential branching fraction and from [35].

We use the HFAG average from [2]. We perform our theoretical computation at NNLO in and NLO in , following Ref. [109] and the references therein.

We consider the latest measurement from LHCb [110] and do not consider the measurement from CMS [111], which has the same central value of LHCb, but larger uncertainty. Moreover, we chose not to use results for , since there are only upper bounds for this decay channel so far [110, 111]. Our theoretical predictions include NLO EW corrections, as well as an NNLO QCD correction, following the detailed expressions obtained in Ref. [112].

Results of the global fit

In this section we present the main results of our work. We perform this study using HEPfit  [113] relying on its Markov Chain Monte Carlo-based Bayesian analysis framework implemented with BAT [114]. We fit to the data using 16 real free parameters that characterise the non-factorisable power corrections, as was done in [24], along with the necessary set of NP WCs. We assign to the hadronic parameters and the NP WCs flatly distributed priors in the relevant ranges mentioned in Sect. 2. The remaining parameters used in the fit are listed in Table 1. To better compare different scenarios, we use the Information Criterion [115, 116], defined aswhere is the average of the log-likelihood and is its variance. The second term in Eq. (9) takes into account the effective number of parameters in the model, allowing for a meaningful comparison of models with different number of parameters. Preferred models are expected to give smaller IC values.

Table 1

Parameters used in the analysis. The Gegenbauer parameters and have flat priors with half width reported in the third column. The remaining ones have Gaussian priors. Meson masses, lepton masses, s-quark mass and electroweak couplings are fixed at the PDG value [117]

Parameters	Mean value	Uncertainty	References
α_s(M_Z)	0.1181	0.0009	[117, 118]
μ_W (GeV)	80.385	−
m_t (GeV)	173.34	0.76	[119]
m_c(m_c) (GeV)	1.28	0.02	[120]
m_b(m_b) (GeV)	4.17	0.05	[121]
f_B_s (MeV)	226	5	[122]
f_B_s/f_B_d	1.204	0.016	[122]
ΔΓ_s/Γ_s	0.129	0.009	[2]
λ	0.2250	0.0006	[123, 124]
A	0.829	0.012	[123, 124]
ρ¯	0.132	0.018	[123, 124]
η¯	0.348	0.012	[123, 124]
f_K^∗,|| (MeV)	204	7	[107]
f_K^∗,⊥(1GeV) (MeV)	159	6	[107]
f_ϕ,|| (MeV)	233	4	[107]
f_ϕ,⊥(1GeV) (MeV)	191	4	[107]
λ_B (MeV)	350	150	[125]
a_1(K¯^∗)_⊥,||	0.04	0.03	[126]
a_2(K¯^∗)_⊥,||	0.05	0.1	[127]
a_2(ϕ)_⊥,||	0.23	0.08	[128]
a_1(K)	0.06	0.03	[126]
a_2(K)	0.115	-	[129]

The results for NP WCs for the several cases that we study can be found in Figs. 1, 2, 3, 4, 5 and 6, where the IC value for each model is also reported, and in Tables 2 and 3 in Appendix A. In Tables 4 and 5, we report the results of the fit for observables of interest. We observe that all cases have comparable IC values except cases (IV) and (V), which are disfavoured in the PMD approach while they remain viable in the PDD one. The main difference between the two approaches is that angular observables, in particular , call for NP in in the PMD approach, while they can be accommodated within the SM in the PDD one.

Let us discuss the various cases in more detail. It is important to stress that the evidence of NP in our fit for cases (I)–(V) is always larger than for one of the semileptonic NP WCs used in the analysis, given the need of a source of LFUV primarily from measurements. In particular, we remark that in the PMD scenarios of cases (I) and (II) we get evidence for NP at more than . However, looking at the corresponding PDD scenarios, the NP evidence gets significantly reduced, roughly between and . The reduction in the significance comes from the larger hadronic uncertainties in the PDD approach which weaken the constraining power of the angular observables on the NP WCs.

Concerning case (III), we observe very similar findings to the ones obtained for case (I), since the effective coupling for the radiative dipole operator is well constrained, especially from the inclusive branching fraction.

Regarding case (IV), in which we vary the three NP parameters and , the model comparison between the PDD and PMD realisation of this NP benchmark is quite informative: NP effects in the dipole operator and in the axial semileptonic currents cannot address at the same time ratios and the anomaly in a satisfactory way when we stick to small non-factorisable QCD power corrections; however, this is no longer true when we allow for a more conservative estimate of the hadronic uncertainties. In particular, the tension in the fit coming from the angular analysis of can now be addressed by large QCD effects as those given in Eq. (8), while a at about can successfully describe all the observational hints of LFUV showed by current measurements. This interesting possibility of axial lepton-flavour violating NP is not found in other global analyses [69-72], as it proceeds from the conservative treatment of hadronic uncertainties we proposed in Ref. [24].

Concerning Tables 4 and 5 of Appendix A, we would like to point out the pattern displayed by the transverse ratios and : cases (I)–(III) predict these values to be with a small error, while the remaining cases give different predictions with the central value ranging between and . Therefore, obtaining experimental information on transverse ratios may help in discerning between the different NP scenarios.

We then show results for case (V), in which we vary , , and correlate the semileptonic vector and axial currents according to and . In analogy with case (IV), only within the PDD approach we find for this NP benchmark a fairly good description of data, with compatible with zero at . Again, we are presented with the case where deviations in angular observables are addressed by large QCD power corrections, while LFUV is driven by semielectronic operators. Looking back at Tables 4 and 5, we note that, for this case, as well as for case (IV) and (VI), both transverse and longitudinal muon over electron ratios in the central- bin, namely and , are characterised by similar central values.

We close our presentation with an analysis of case (VI) in which we float simultaneously , , , , and . As can be seen from Fig. 6, current measurements are informative enough to constrain, at the same time, all the NP WCs both in the PMD and PDD approaches. In particular, within the latter case, a nontrivial interplay among NP effects encoded both in and , together with the hadronic contributions reported in Table 3, produces the weakest hint in favour of NP provided by our global analysis – sitting between and level – while allowing for a very good description of the entire data set, similar to the other cases. The power corrections we found are larger than those obtained in Ref. [90], but smaller than those required by the SM fit of  [24]. As discussed in detail in Refs. [27, 89], the size obtained for the power corrections is compatible with the naive power counting relative to the leading amplitude. We stress (once again) that a more optimistic attitude towards the estimate of QCD power corrections (PMD approach) leads to the a much stronger claim in favour of NP, at a statistical significance larger than .

In Tables 2 and 3 we report mean and standard deviation for the NP WCs and absolute values of for all the cases considered in the analysis. It is also relevant to observe that, once we switch on NP effects through in order to attempt at simultaneously explaining observables such as and in the PDD approach we find values for compatible with zero at . Conversely, if we set then a nonvanishing is needed to account for the angular observables, as found in Ref. [24], showing that one cannot disentangle hadronic uncertainties and NP in at present.

Discussion

In this work, we critically examined several BSM scenarios in order to possibly explain the growing pattern of B anomalies, recently enriched by the measurement performed by the LHCb collaboration [76]. We carried out our analysis in an effective field theory framework, describing the non-factorisable power corrections by means of 16 free parameters in our fit along the lines of Ref. [24].

We performed all our fits using two different hadronic models. The first approach, labelled PMD, relies completely on the phenomenological model from Ref. [90] and corresponds to the more widely used choice in the literature. The second one, named PDD, imposes the result of Ref. [90] only at ,3 allowing the data to drive the hadronic contributions in the higher invariant mass region.

Regarding the NP contributions, we analyse six different benchmark scenarios, differentiated by distinct choices of NP WCs employed in the fits. Case (I) allows for and , while case (II) considers the scenario with and ; case (III) studies NP effects coming as , and , and case (IV) is the same as the latter but with instead of ; case (V) studies the possibility described in the third case with and enforced; finally, case (VI) considers the general case with all the five NP WCs being allowed to float independently. Our main results are collected in Figs. 1, 2, 3, 4, 5 and 6 and also reported in Tables 2, 3 and 5.

The comparison of different scenarios using the IC shows that all the considered cases are on the same footing except for cases (IV) and (V). These cases are strongly disfavoured in the PMD approach, as there is no in case (IV) to account for the deviation in , while is constrained by its correlation with and the measured value of BR in case (V).

In fact, from our analysis of radiative and (semi)leptonic B decays we identify two classes of viable NP scenarios:Finally the most general fit we performed, namely case (VI), confirms in the PDD approach that both scenarios above are viable, although a slight preference for is found. More data are needed to assess what kind of NP scenario (if the anomalies persist) is realised in Nature.

The widely studied scenario: from Figs. 1, 2 and 3, we find a remarkable evidence in favour of in the PMD approach. It is indeed nontrivial that a single NP WC can explain all the present anomalies in transitions [4, 68–72]. However, in the more conservative PDD approach, the significance of a nonvanishing drops to about , mainly driven by LFUV.

An alternative scenario with nonvanishing , which emerges in the presence of large hadronic corrections to the infinite mass limit, namely our PDD approach. To our knowledge, a NP electronic axial current has not been studied in the literature, since it does not provide a satisfactory description of the angular observables within the commonly used PMD approach. We think that the present theoretical status of power correction calculations is not robust enough to discard this interesting NP scenario.