How well can we guess theoretical uncertainties?

Physics Letters B 726 (2013) 266–272

Contents lists available at ScienceDirect

Physics Letters B www.elsevier.com/locate/physletb

How well can we guess theoretical uncertainties? ✩ André David a,∗ , Giampiero Passarino b,c a b c

PH Department, CERN, Switzerland Dipartimento di Fisica Teorica, Università di Torino, Italy INFN, Sezione di Torino, Italy

a r t i c l e

i n f o

Article history: Received 7 July 2013 Received in revised form 1 August 2013 Accepted 9 August 2013 Available online 16 August 2013 Editor: G.F. Giudice

a b s t r a c t The problem of estimating the effect of missing higher orders in perturbation theory is analyzed with emphasis in the application to Higgs production in gluon–gluon fusion. Well-known mathematical methods for an approximated completion of the perturbative series are applied with the goal to not truncate the series, but complete it in a well-defined way, so as to increase the accuracy – if not the precision – of theoretical predictions. The uncertainty arising from the use of the completion procedure is discussed and a recipe for constructing a corresponding probability distribution function is proposed. © 2013 Elsevier B.V. All rights reserved.

1. Introduction In the past 30 years, the commonly accepted way to estimate theoretical uncertainties associated to collider physics observables has been based on the notion of QCD scale variations. We introduce the concept of MHO(U), missing higher order (uncertainty), which is linked to the truncation error in the perturbative expansion. At present, and for some time to come, estimations of observables will be based on a finite number of terms of a series, such that additional information on the behavior of that series should be exploited. Regardless of their precision, truncated calculations are only as accurate as the higher orders that they lack. A more accurate evaluation of the observable may be obtained by estimating the MHO. The issue of precision then becomes more tightly bound to the estimation of the MHOU, taking into account both the uncertainty on the MHO estimation procedure as well as any uncertainties in the terms that have already been calculated. In this Letter, the problem of MHO(U) in Higgs production through gluon–gluon fusion is approached using sequence transformations to improve the rate of convergence of the series and directly estimate the MHO. In Section 2 we discuss foundational issues related to the MHO problem and the applicability of sequence transformations. In Section 3 we summarize existing calculations of Higgs production through gluon–gluon fusion, mapping out the inputs needed to estimate the MHO. Then, in Section 4, we in-

✩ Work supported by MIUR under contract 2001023713_006 and by Compagnia di San Paolo under contract ORTO11TPXK. Corresponding author. E-mail addresses: [email protected] (A. David), [email protected] (G. Passarino).

*

0370-2693/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.physletb.2013.08.025

troduce different types of sequence transformations and discuss in detail their performance in synthetic problems as well as applications to series involving physical observables. The main results are then presented in Section 5 where we apply sequence transformations to the problem of Higgs production through gluon–gluon fusion problem and propose an estimate of MHOU and its probability distribution function (pdf). Finally, in Section 6 we summarize the main arguments and results. 2. MHOU beyond scale uncertainties Consider an observable X ( Q , μ), where Q is the typical scale of the process, and μ ≡ {μ R , μ F } are the renormalization and factorization scales. The traditional procedure to estimate MHOU through scale variations [1] defines

  , X ( Q , ξ μ) , ξ     μ , X ( Q , ξ μ) , X ξ+ ( Q , μ) = max X Q , ξ  

X ξ− ( Q , μ) = min X Q ,

μ

(1)

or variations thereof, see Ref. [2]. Selecting a value for ξ (typically ξ = 2) the prediction is that

X ξ− ( Q , μ) < X ( Q , μ) < X ξ+ ( Q , μ).

(2)

There are several examples in the literature where the ξ = 2 scale uncertainty of the nth-order underestimates the (n + 1)th order calculation. There is also an open and debatable question on how to assign a probability distribution function (pdf) to the MHOU thus obtained [3]. The procedure that is most commonly used is based on a Gaussian (or log-normal) distribution centered at

A. David, G. Passarino / Physics Letters B 726 (2013) 266–272

μ = X c = X ( Q , Q ). This choice of central value is afflicted by the accuracy issues from truncation and there are cases in which the scale has been adapted to match resummation [4,5]. What to use for the standard deviation remains an open problem, though a common ansatz is to use σ = max( X 2+ ( Q , Q ), X 2− ( Q , Q )). Alternatively, it could be assumed that the pdf is a uniform distribution 

P (X) =

1 , X 2+ ( Q , Q )− X 2− ( Q , Q )

X 2− ( Q , Q ) < X < X 2+ ( Q , Q ),

0,

otherwise.

Recently, Cacciari and Houdeau made a proposal to derive the pdf based on a flat (uninformative) Bayesian prior for the MHOU from the scale-variation prescription [2]. More generally, the dependence on scales is only one part of the problem, as the MHO problem is based on how to interpret the relation between an observable O , and a perturbative series

O∼

∞ 

cn g n .

(3)

n =0

The perturbative expansion of Eq. (3) is unlikely to converge [6] (see also Refs. [7–12]) and the asymptotic behavior of the coefficients is expected to be cn ∼ K nα n!/ S n when n → ∞, and where K , α and S are constants [13]. An overview of the mathematical theory of divergent series and interpretation of perturbation series is given in Ref. [14]. The requirement of Eq. (3) (∼) is not a formal one; it has the physical meaning of a smooth transition between the system with interaction and the system without it [15]. Furthermore, Borel and Carleman proved that there are analytic functions corresponding to arbitrary asymptotic power series [16]. Clearly, this is not the proper place to summarize the highorder asymptotics thematics and analyticity conditions, and we refer to the relevant literature, see Ref. [17] and Refs. [15,16,18–21]. We only mention the problem of renormalons: the missing terms in the perturbative expansion are power-suppressed in the hard scale; thus any summation procedure is defined up to power suppressed corrections [22]. Furthermore, it is well known that the leading higher-order behavior in QCD related to (Borel-summable) ultraviolet renormalons can be eliminated by adding irrelevant dimension-six operators to the QCD Lagrangian. The mechanism of these contributions is such that they become important starting at a sufficiently high order in the usual perturbative expansion and are simply not seen to the accuracy of present calculations [23]. We should stress that recoverability of a function by means of its asymptotic series requires “enough” analyticity [15]; for a discussion on the corresponding analyticity conditions and on their failure or fulfillment, in particular in the QCD context, see Ref. [17]. Any work on MHO(U) should face these issues and as we discuss the example of Higgs production via gluon–gluon fusion it also is worth nothing that the authors of Ref. [2] do not make assumptions on the analyticity domain; starting from Eq. (3), they k estimate the remainder R k = O − n=0 cn g n , to be R k ≈ ck+1 g k+1 with ck+1 = max{|c 0 |, . . . , |ck |}. This, in turn, reflects into a width of ck+1 g k+1 for the flat part of the uncertainty pdf. Therefore, the MHO problem and its associated uncertainty can be summarized in one point: how can we make predictions for higher-order perturbative coefficients, whose explicit calculation is cumbersome and time-consuming, while keeping a balance with analyticity? We will not be able to answer general questions (namely to prove uniqueness of our results) and will rather concentrate on predicting higher orders using the well-known concept of “series acceleration” [24–26], i.e., one of a collection of sequence transformations for improving the rate of convergence of a series. If the

267

Table 1 √ Numerical values as derived from Ref. [30] assuming s = 8 TeV. These values are the relevant inputs to an estimation of MHO(U): while traditionally MHOU is estimated from the scale variation of γ3c , the proposed procedure only requires the values in the middle column (μ = M H ). Our notation

Ref. [30]

γ1

1 K gg

γ2

2 K gg

γ3c ± γ3

3 K gg

μ = M H /2

μ = MH

μ = 2M H

11.879 72.254 168.98 ± 30.87

377.20 ± 30.78

681.72 ± 29.93

original series is divergent, the sequence transformation acts as an extrapolation method. In the case of infinite sums that formally diverge, the helpful property of sequence transformations is that they may return a result that can be interpreted as the evaluation of the analytic extension of the series for the sum. The relation between Borel summation (the usual method applied for summing divergent series) and these extrapolation methods was noted for the first time in Refs. [27,28]. Note that the definition of the sum of a factorially divergent series, including those with non-alternating coefficients, is always equivalent to Borel’s definition (see Section 7 of Ref. [14]). 3. Existing calculations of Higgs production via gluon–gluon fusion Let us consider what is presently known of Higgs production via gluon–gluon fusion, i.e., the process gg → H. There have been several attempts to compute approximate N3 LO corrections, see Refs. [29–31]. Here we follow the work of Ref. [30] and define



σgg (τ , M H2 ) = K gg τ , M H2 , αs 0 2 σgg (τ , M H ) =1+

∞ 

n αsn (μ R ) K gg (τ , μ = M H ),

(4)

n =1 2 0 where τ = M H /s, σgg is the LO cross-section, and the K -factor K gg was expanded in powers of αs (μ R ). In Eq. (4) that Nit is understood n when computing the partial sums S N = 1 + n=1 αsn (μ R ) K gg , αs is computed at the highest level, i.e., NLO for S 1 , NNLO for S 2 , etc. n 1 Introducing γn = K gg , the known values are γ1 = K gg (τ , μ = 2 M H ) = 11.879 and γ2 = K gg (τ , μ = M H ) = 72.254. In their recent work, the authors of Ref. [30] computed an approximation √ 3 of αs3 (μ) K gg (μ) at s = 8 TeV for μ = M H /2, M H , and 2M H . Since 3 αs3 (μ) K gg (μ) is only known within a given interval (see Table 1 and discussion after Eq. (4.1) of Ref. [30]) we report in Table 1 the numerical values of γ3 as a central value (γ3c ) and the corresponding uncertainty range (γ3 ). In Table 1 one can immediately see that the approximate cal3 culation of K gg can be varied in two ways: (1) the change in γ3c via scale variation, and (2) the intrinsic uncertainty γ3 due to the approximate nature of the result. While the traditional approach to MHOU estimation considers the effect from scale variation, the procedure that we put forth in later sections combines γ3 with the uncertainty on the estimation of the MHO based on sequence transformations.

4. Sequence transformations The theory of sequence transformations is a well-established branch of numerical mathematics with many applications in science, as described in Refs. [32–34] and Ref. [35]. As an example in connection with the summation of the divergent perturbation expansion of the hydrogen atom in an external magnetic field, the

268

A. David, G. Passarino / Physics Letters B 726 (2013) 266–272

work of Refs. [36,37] introduces a new sequence transformation which uses as input not only the elements of a sequence of partial sums, but also explicit estimates for the truncation errors. Through sequence transformations, slowly convergent and divergent sequences and series can be transformed into sequences and series with hopefully better numerical properties. Thus they are useful for improving convergence. In most situations, a sequence transform does not sum a series exactly; however, in many cases, it correctly predicts some of the unknown terms of the sequence. 4.1. The Levin τ -transform Let us recall the definition of the Levin see nτ -transform, i Refs. [38–40]. Given the partial sums S n = i =0 γi z we define the τ -transform as

k

τkn (β) =

τ

i =i 0 W (n, k, i , β) S n+i k τ i =i 0 W (n, k, i , β)

,

τk ≡ τk0 ≡ τk0 (0)

where i 0 = max{0, n − 1} and W τ (n, k, i , β) = (−1)

(5)



i k (β+n+i )k−1 i

 S n+i −1

,

where ( z)a = ( z + a)/ ( z) is the Pochhammer symbol, and  is the usual forward-difference operator,  S n = S n+1 − S n . The algorithm for estimating the first unknown coefficient is based on the Taylor expansion of τk ; if S 1 , . . . , S k are known, one computes τk − S k = γ k+1 zk+1 + O ( zk+2 ) and γ k+1 is the prediction for γk+1 . Of course, this prediction is not expected to be very reliable for small values of k. Nevertheless, applying τ2 − S 2 = (γ22 /γ1 )z3 + O(z4 ) to the series of Eq. (4), one predicts γ 3 (μ = M H ) = 439.48 which has the correct sign and the right order of magnitude when compared with the results from Ref. [30], 346.42 < γ3 (μ = M H ) < 407.48. 4.1.1. Recursive estimation of unknown coefficients Let us outline our algorithm to improve the convergence of a series. This algorithm can be used with any of the transforms introduced later and to any of the series also discussed in the examples later. We give it below in an explicit form for the Levin τ -transform τk0 (β) and applied to the series of Eq. (4) assuming that the inputs from Table 1 are known: 1.

Use the first 3 terms in Eq. (4), choose γ3 = γ3c ( (5+2β)γ22 −(3+β)γ1 γ3 γ + γ1 γ3 − γ22 ]. derive γ 4 = 3 γ 3γ [2 12+7β+β 2 1 2

μ = M H ), and

2. Construct S 4 assuming γ4 = γ 4 .

3. Derive γ 5 = γ γ 4γ (120 + 72β + 15β 2 + β 3 )−1 , where ϑ = 1 2 3 4γ22 γ3 (6 + 11β + 6β 2 + β 3 ) − 6γ1 γ32 (24 + 26β + 9β 2 + β 3 ) + 4γ1 γ2 γ 4 (60 + 470β + 12β 2 + β 3 ). 4. Construct S 5 assuming γ5 = γ 5 . 5. Repeat the previous steps until τ3 , . . . , τ6 are constructed. 6. Compare the S 3 , . . . , S 6 with the τ3 , . . . , τ6 . 7. Repeat steps 1–6 for γ3 = γ3c + γ3 and γ3 = γ3c − γ3 , always taken at μ = M H . ϑγ

The whole strategy is based on the fact that one can predict the coefficients by constructing an approximant with the known terms of the series (γ0 , . . . , γn ) and expanding the approximant in a Taylor series. The first n terms of this series will exactly agree with those of the original series, while the subsequent terms may be treated as predicted coefficients. Once the series is completed via an algorithm such as the one above, the dependence on μ is removed, and the notion of scale variation with it. This implies that the uncertainty estimation is moved from scale variations to the completion procedure, as discussed in Section 5.3. This procedure represents an extension of the work in Ref. [2].

4.1.2. β -tuning If γ1 , . . . , γ3 are known, the values of γ1 and γ2 can then be γ used to compute β = (1 − 12 γ3 12 )−1 − 2. This value of β is such γ2

that γ 3 = γ3 . With β determined this way, one can then apply the recursive algorithm above. For a discussion on β -tuning of the Levin τ -transform, with applications to predicting new coefficients in the g − 2 of muon and electron, see Ref. [35]. There, given aμ = a(1/2 + 0.7655a + α and tuning β = −0.90935, one 24.05a2 + 125.6a3 ), where a = π derives γ 4 = 513.3 with γ4 expected in the range 433 < γ4 < 713. The β -tuned procedure is used to cross-check results without β -tuning in Section 5.3. 4.2. The Weniger δ -transform A second transform that we have considered in detail in later Sections is the δ -transform introduced by Weniger [39]:

k

δk (β) = i =k 0

W δ (k, i , β) S i

δk ≡ δk (1),

,

W δ (k, i , β) k (β+i ) where W δ (k, i , β) = (−1)i i (β+k)k−1

(6)

i =0

1

γi +1 zi +1

k−1

. Following the same

algorithm described in Section 4.1.1, the predicted γ n values for

δk (1) are γ 4 = 18γ1 γ2 γ 4 ).

γ3

3γ 1 γ 2

(4γ1 γ3 − γ22 ) and γ 5 =

γ4

10γ1 γ2 γ3

(γ22 γ3 − 9γ1 γ32 +

4.3. Other series transformations There are other well-known transforms that we have tested using the algorithm described above:

• Wynn’s ε -algorithm [41], the nonlinear recursive scheme n+1 1 n n n ε− . n 1 = 0, ε0 = S n , εk+1 = εk−1 + n+1 −εk

εk

• Brezinski’s J -algorithm [42], based on the recursive scheme J kn+1 = J kn+1 −

J 0n = S n ,



J kn

=

J kn+1

−



 J kn  J kn+1 2 J kn+1 , n +2 2 n J k  J k −  J kn 2 J kn+1

J kn .

(7) k

• Generalized Levin t-transform [38]: tkn = k (n+i ) where W t (n, k, i ) = (−1)i i  S k−1 . n+i −1

i =i

k0

W t (n,k,i ) S n+i

i =i 0

k

i =i

• Generalized Levin y-transform [43]: ynk = k (n+i ) with W y (n, k, i ) = (−1)i i  S k−2 . n+i −1

k0

W t (n,k,i )

W y (n,k,i ) S n+i

i =i 0

W y (n,k,i )

,

,

4.4. Example applications of sequence transformations To discuss our results, we introduce the following notation: n N S N,n = k=0 γk zk + k=n+1 γ k zk , and τN,n constructed accordingly. For example, τ6,3 = N 6,3 / D 6,3 with

N 6,3 = −

720 S1 − S0

− D 6,3 = −

90 720 S 5,3 − S 4,3 720

S1 − S0

−

S1 +

S 5,3 +

S1 +

90 720 S 5,3 − S 4,3

10 800 S2 − S1

50 400 S3 − S2

30 240 S 6,3 − S 5,3

10 800 S2 − S1

+

S2 −

−

50 400

30 240

.

100 800 S 4,3 − S 3

S 4,3

S 6,3 ,

S3 − S2

S 6,3 − S 5,3

S3 +

S3 +

100 800 S 4,3 − S 3 (8)

A. David, G. Passarino / Physics Letters B 726 (2013) 266–272

Table 2 Actual and predicted coefficients for the series of Eq. (9), which was designed so as to approximately reproduce the values in Table 1 when z = 0.1.

γ1

γ2

γ3

γ4

γ5

γ6

12.620

73.322

259.56

624.24

1076.2

1366.8

γ4

γ5

γ6

624.89

1081.8

1388.9

Table 3 Predictions for the series of Eq. (11). S∞

S6

τ6

S 6,3

τ6,3

1.09773772 –

1.09743700 0.027%

1.09778864 0.005%

1.09705909 0.062%

1.09705234 0.062%

The transformations listed in Section 4.3 have been applied to a suite of test series. Note that for the calculations one could have used readily-available software, e.g. the one described in Ref. [44], Maple [45], or GSL [46]. The suite of test series considered includes:

• We first considered the series S ∞ = (1 + z)ν = 1 +

∞ 

γn zn ,

ν = 12.62,

(9)

n =1

where ν was tuned such that its first 3 coefficients are similar to those of the series in Eq. (4). The sum of the series for z = 0.1 is S ∞ = 3.32947445. Using up to 6 γ -coefficients, shown in Table 2, we derive that the best improvement for the rate of convergence is obtained with the Levin τ -transform of Eq. (5) with β = n = 0:

S ∞ = 3.32947445,

S 6 = 3.32933563,

τ6 = 3.32947445.

(10)

Table 2 also shows the partial results of using our recursive approximations algorithm. Eventually, we obtain τ 6 = 3.32962298, or τ 6 / S ∞ − 1 = 0.0045%. • The goodness of the approximation has also been tested by expanding the hypergeometric function 2 F 1 (n + 1/2, n + 1; n + 3/2; z2 ) for large values of n, with positive results: in all cases convergence is improved.  • Several other examples, (1 + z)1/2 , ln(1 + z), e z , n∞=0 (−1)n × n! zn , Φ(n, z, a), where Φ is the Learch Phi-function, can be found in Ref. [35]. The same work provides examples where higher-order coefficients are estimated, e.g., aμ,e (muon or electron g − 2) and the hadronic ratio R. • Consider now the case of an asymptotic series, e.g.

S∞ =

∞ 

n! zn+1 = e −1/z Ei(1/ z)

(11)

n =0

where the exponential integral is a single-valued function in the plane cut along the negative real axis. However, for z > 0, Ei( z) can be computed to great accuracy using several Chebyshev expansions. Note that the r.h.s. of Eq. (11) is the Borel sum of the series. The approximation returned by the γ n is not of high quality. Nevertheless, the approximation works reasonably well and τ6,3 is not worse that S 6,3 , as shown in Table 3. It has been shown in Ref. [43] that there is a large class of series that have Borel sums that are analytic in the cut-plane and the numerical results of Ref. [47] suggest that Levin–Weniger transforms produce approximations to these Borel sums. Furthermore, in

269

Table 4 Predictions for the series of Eq. (12). nf

S4

τ4,3

4 5

1.145469 1.138618

1.146096 1.140940

Ref. [33], numerical evidence is shown suggesting that the Weniger transform can resum a function with singularities in the Borel plane (but not on the positive axis). • Other relevant examples are: the prediction for the fifth (known) coefficient of the β -function of the Higgs boson coupling, the derivative expansion of QED effective action, and the partition function for zero-dimensional φ 4 theory [33]. • We have also tested the method against some recent calculations like the leptonic contributions to the effective electromagnetic coupling at four-loop order in QED. The coefficients of α /π and αlep are [48]: γ1 = 13.52631(8), γ2 = 14.38553(6), and γ3 = 84.8285(7). The predicted and known results for γ4 are γ 4 = 705.22 and γ4 = 770.76, for a relative difference of γ 4 /γ4 − 1 = −8.5%. • There are cases where the algorithm cannot make a reliable prediction, such as in predicting QCD corrections to the QED β -functions, see Refs. [49,50]. Looking at Eqs. (4.4)–(4.6) of Ref. [49] we see series with sudden jumps of sign in the coefficients; for instance, the series for 5 flavors is



α 2 1.667 + 1.667aS + 2.813a2S − 5.971a3S − 32.336a4S



(12)

with aS = αs /π . Our results, based on τ4,3 are shown in Table 4. Here, neglecting the term O (a4S ) or computing τ4 with an approximated γ 4 gives a difference of the same size. In this case we are considering a series representing a selfenergy that will have a two-particle cut (with the corresponding series of corrections), a three-particle cut (with the corresponding series of corrections), etc. Therefore, at each order in perturbation theory new contributions (i.e. new series) will arise and it is unsafe to make a guess by using only the first 3 orders. However, in this case, using γ 4 as an estimate of the uncertainty in S 3 gives a reasonable result: 0.95 < |τ4 − S 4 |/γ4 a4S < 1.14. Further examples of the performance of different transforms on a number of test sequences can be found in Refs. [51–57]. 5. Application to Higgs production via gluon–gluon fusion For all the examples considered in Section 4.4 we have found that the Levin τ -transform and the Weniger δ -transform provide the fastest convergence. The power of these transformations is due to the fact that the explicit estimates for the truncation error of the series are incorporated into the convergence acceleration. The Levin τ -transform has been shown to work with good accuracy for the prediction of higher-order coefficients of alternating and nonalternating factorially-divergent perturbation series, see Ref. [33]. Arguments supporting the general applicability of Levin transforms to the series of mathematical structures expected from quantum field theory can be found also in Ref. [33]. It should be noted that in the τkn sequence transformation, the superscript n indicates the minimal index occurring in the finite subset of input data, while k, the order of the transformation, is a measure of the complexity for the transformation itself. It is worth noting that τkn requires knowledge of the first n + k partial sums, that is why we limit our considerations to τk ≡ τk0 . The most important question concerns the reliability of the procedure when applied to the series of Eq. (4).

270

A. David, G. Passarino / Physics Letters B 726 (2013) 266–272

Table 5 Effect of QCD scale variation for predicted higher-order terms in the Higgs gluon– −,n +,n ±,n gluon fusion production cross-section. D n = 1 − σ2 /σ2 with the σξ are defined in Eq. (14). In the extrapolation region (n  3) the variation decreases as expected from a reliable estimate of MHO. Dn

γ3c − γ3

γ3c

γ3c + γ3

D2 D3

14.93%

27.02% 16.08%

17.21%

D4 D5 D6

6.59% 2.20% 0.14%

7.99% 3.08% 0.39%

9.68% 4.59% 1.38%

Table 6 Predicted higher-order coefficients in gluon–gluon fusion, computed at

γn γ4 γ5 γ6

Levin-τ

μ = MH .

Weniger-δ

γ3c − γ3

γ3c

γ3c + γ3

γ3c − γ3

γ3c

γ3c + γ3

1437.9 5412.4 18979.0

1806.6 8185.6 35677.0

2214.7 11733.0 61133.0

1512.2 6276.6 25243.0

1860.8 8912.3 41918.0

2244.3 12183.0 65605.0

Table 7 R n = αs γ n+1 /γ n (γn ) for the τ and δ transforms (note that the denominator is not an extrapolation when available). It can be seen that R n is constant to better than 10% in the extrapolation region (n  3) for both transforms. Rn Levin-τ Weniger-δ

R0 1.3280 1.3280

R1 0.6800 0.6800

R2 0.5836 0.5360

R3 0.5354 0.4880

R4 0.5065 0.4640

R5 0.4873 0.4496

5.1. Applicability In motivating the applicability of the procedure, scale variation can be of use. Consider n 0 σgg (τ , μ) = σgg (τ , μ) S n,3 (μ),

S n,3 (μ) = 1 +

3 

αsk (μ)γk (μ) +

k =1 2 where τ = M H /s, Introducing

σξ−,n σξ+,n



√

n 

αsk (μ)γ k (μ)

(13)

k =4

s = 8 TeV, and vary the QCD scales with ξ = 2.

= min σ ( Q , μ/ξ ), σ ( Q , ξ μ) ,  = max σ n ( Q , μ/ξ ), σ n ( Q , ξ μ) n

(15)

and report the result of the calculations in Table 8 where the S,n coefficients needed to construct σgg are based on τ -transform. To understand the comparison one should bear in mind that sequence transforms can also be characterized by the highest coefficient involved: τk requires γ k but δk requires γ k+1 . Therefore, we expect τk and δk−1 to give predictions with comparable quality. The results show that using the Levin τ -transform improves the convergence; indeed n = 3 is already a good approximation with τ ,6 τ ,3 σgg /σgg − 1 = 0.13%. The use of other transforms is compatible with τ6 to within 2%: if we use the Weniger δ -transform of Eq. (6) δ,5 τ ,6 (with β = 1) we obtain σgg /σgg − 1 = 1.38%. Additionally, we have investigated the use of β -tuning, using the Levin τ -transform Eq. (5) with β = 0. To have γ3 = γ3c , we find β = −0.2482, and calculate

τ ,5 σgg (β) = 23.542 pb. This is to

τ ,5

be compared with σgg (β = 0) = 23.105 pb, the difference being within the uncertainty induced by γ3 . Our conclusion is that β -tuning is a procedure to be adopted in those cases where there is a reasonable guess on the value of the next coefficient or on the interval where it is expected. Furthermore, all cases where the β -tuned results are substantially different from β = 0 should be taken with the due caution. Finally, basing the whole procedure on δ -transforms or estimating the coefficients with δk (1) and accelerating the series with τ6 gives consistent results, namely 23.4241 pb (with δ5 ) in the first case and 23.4253 pb in the second. It is worth noting that if any of the transforms predicts at least one extra coefficient of the series, then in principle the whole function is known, which is unlikely to be the case in any physical problem. We can only conclude that a judicious use can make predictions at some relatively good level of accuracy. We also know that all transforms basically differ in the choice of the remainder estimates. A good choice should satisfy the following asymptotic condition [34]: R n = S ∞ω− S n ∼ c, when n → ∞, where ωn is the ren mainder. Levin selects ωn =  S n−1 and, from Table 7, we derive an approximate relation γn+1 αs ≈ K γn , for n > n0 , where K is a constant with K < 1. In this case R n → 1/(1 − K ) for sufficiently large n. 5.3. Discussion of MHOU

n

−,n

0 σggX ,n = σgg (μ = M H ) Xn,3 (μ = M H ), with X ∈ {S, τ , δ}

(14)

+,n

and D n = 1 − σ2 /σ2 we obtain the values reported in Table 5. Comparing the results for D 2 and D 3 it can be seen that the variability due to scale variation is substantially reduced by the inclusion of the N3 LO term. We expect that a reliable estimate of the missing higher orders should follow the trend of further reducing the effect as is the case. The coefficients of the perturbative series, computed with τk at μ = M H , are given in Table 6. The ratio R n = αs γ n+1 /γ n becomes constant to a very good approximation, and is given in Table 7, where δn is defined in Eq. (6). Note that this does not represent a formal proof that there is an upper bound on the remainder but makes plausible the argument in favor of that. 5.2. Numerical results Our strategy for estimating MHO and MHOU can be summarized as follows: we select a scale, μ = M H , for gluon–gluon fusion, and estimate the uncertainty due to higher orders at that scale. This implies that the (scale variation) uncertainty at the chosen scale is part of the uncertainty due to higher orders and should not be counted twice. Therefore, we compare

Given that the sequence transform procedures outlined above provide an estimate for the sum of the full series, when estimating the uncertainty on that quantity we will be deliberately conservative. 5.3.1. Uncertainty due to MHO estimation Given the different nature of the calculations represented by S,3 δ,5 σgg and σgg , it can be expected that, to a very good accuS,3

δ,5

racy, σgg < σgg < σgg . For γ3 = γ3c , this defines the interval [20.13, 23.42] pb that has a relative width of 16.4%. For comparison, the N3 LO calculation for μ = M H and γ3 = γ3c yields σgg = 20.13 pb and traditional QCD scale variations with ξ = 2 leads to the interval [18.90, 21.93] pb that has a relative width of 16.1% (the authors of Ref. [30] quote ±7%). It is worth noting how in our approach the interval is shifted by ≈ +7% with respect to the N3 LO result. This is to be compared to the +17% of N3 LO with respect to NNLO [30]. 5.3.2. Uncertainty due to γ3 We can now discuss how to take into account the uncertainty on γ3 induced by the γ3 (μ = M H ). In line with a simple and conservative approach that can later be refined, we consider all values of γ3 in the interval [γ3c − γ3 , γ3c + γ3 ] as equally likely and take the lowest value of

S,3 δ,5 σgg and the highest value of σgg .

A. David, G. Passarino / Physics Letters B 726 (2013) 266–272

Table 8 Cross-sections obtained using Eq. (15), using μ = M H . For shifted so that rows represent the same order in γ n .

γ3c

3 4 5 6

δ,n σgg , β = 1 is used for the Weniger δ -transform. Note that in the case of the Weniger transform the index is

− γ3

19.89889 21.10181 21.60801 21.80644

γ3c

γ3c

+ γ3

20.12922 21.64063 22.40620 22.77922

20.35954 22.21236 23.30967 23.94883

γ3c

− γ3

21.83017 21.92044 21.91988 21.91988

5.3.3. Result The previous choices lead to an interval with a relative width of 26.01%, shifted by at least +5% with respect to the N3 LO result:



δ,n−1 σgg [pb]

τ ,n σgg [pb]

S,n σgg [pb]

n





S,3 c δ,5 c γ3 − γ3 , σgg γ3 + γ3 σgg ∈ σgg



= [σ− , σ+ ] = [19.89889, 25.07525] pb.

271

(16)

To conclude, our prediction is that the “true” cross-section value is bracketed by the estimations of Eq. (16) as all other transforms fall in that interval. For instance, J 12 from Eq. (7) gives 23.018 pb and τ41 gives 23.244 pb. The advantages of our recipe for estimating MHOU are that the result does not depend on the choice of the parameter expansion (it is based on partial sums) and it takes into account the nature of the coefficients, i.e., that the known terms of the perturbative expansion in gluon–gluon fusion are positive. Starting from the proposal in Eq. (16), the corresponding pdf can be derived following the work of Ref. [2]. 6. Conclusions The flat part of the MHOU pdf has been chosen observing that

σ S,3 is the last known term of the series, that known and predicted coefficients are all positive, and that all transforms “predict” convergence towards a value inside the interval of Eq. (16) and close to σ δ,5 . Therefore, our best guess is the one in Eq. (16) since it would be ambitious to claim that σ τ ,6 or σ δ,5 are the result, with a very small error. One should mention in this regard that there is no proof of the uniqueness of the result reconstructed from its asymptotic series. There is only evidence that all sequence transforms produce a result within a given interval to which we assign an uninformative prior, in the Bayesian sense. It should be mentioned that we have included only the gg-channel. At higher orders we have new channels, new color structures, etc. For instance, the qg-channel contribution is negative; at low orders its contribution is sub-leading but nothing is known at higher orders. This is a general problem that will affect all procedures aimed at estimating MHO(U). Finally, it should be stressed that all re-summation procedures for non-alternating (divergent) series usually fail when the parameter expansion is on a (expected) cut in the complex plane. We support the strategy presented for deriving information on MHO(U) with the following arguments:

• Given the (few) known coefficients in the perturbative expansion, we estimate the next (few) coefficients and the corresponding partial sums by means of sequence transformations. This is the first step towards “reconstructing” the physical observable in Eq. (3). • The use of sequence transformations was tried on a number of test sequences, including several physical observables. • A function can be uniquely determined by its asymptotic expansion if certain conditions are satisfied [19]. • The Borel procedure is a summation method which, under the above conditions, determines uniquely the sum of the series. It

γ3c

γ3c

+ γ3

γ3c − γ3

γ3c

γ3c + γ3

23.07444 23.10458 23.10473 23.10473

24.83209 24.79751 24.79661 24.79658

22.21881 22.21864 22.21864

23.42508 23.42407 23.42405

25.07780 25.07535 25.07525

should be taken into account that there is a large class of series that have Borel sums (analytic in the cut-plane) and there is evidence that Levin–Weniger transforms produce approximations to these Borel sums. This is one of the plausibility arguments supporting our results. • The QCD scale variation uncertainty decreases when we include new (estimated) partial sums. • All known and predicted coefficients are positive and all transforms predict convergence within a narrow interval. • Missing a formal proof of uniqueness, we assume an uninformative prior between the last known partial sum and the (largest) predicted partial sum. The arguments developed in this work support the opinion that perturbation theory up to N3 LO is essential to obtain accurate definition of the theory (MHO) and shed some light on how to formulate consistent procedures for accurate computations (MHOU). We conclude by saying that “new” insights into the properties of perturbative expansions are always important, since computing higher-order corrections is not only cumbersome and costly but also suffers fundamentally from the divergence of the series. The investigation of QCD-scale and renormalization-scheme dependence of a truncated series should not be confused with the attempt to estimate its uncalculated remainder which is the true source of MHOU. Acknowledgements We acknowledge the LHC Higgs Cross Section Working Group where the problem addressed here was posed and as the forum that brought the authors together. A.D. recognizes useful discussions with Kirill Melnikov at the XXV Rencontres de Blois. G.P. acknowledges important discussions with Matteo Cacciari and Stefano Forte. References [1] S. Dittmaier, et al., Handbook of LHC Higgs cross sections: 1. Inclusive observables, arXiv:1101.0593. [2] M. Cacciari, N. Houdeau, Meaningful characterisation of perturbative theoretical uncertainties, J. High Energy Phys. 1109 (2011) 039, http://dx.doi.org/10.1007/ JHEP09(2011)039, arXiv:1105.5152. [3] THUTF meeting: Missing higher orders and PDF uncertainties, https://indico. cern.ch/conferenceDisplay.py?confId=251810, 2013 [online accessed 3 July 2013]. [4] C. Anastasiou, K. Melnikov, Higgs boson production at hadron colliders in NNLO QCD, Nucl. Phys. B 646 (2002) 220–256, http://dx.doi. org/10.1016/S0550-3213(02)00837-4, arXiv:hep-ph/0207004. [5] S. Catani, D. de Florian, M. Grazzini, P. Nason, Soft gluon resummation for Higgs boson production at hadron colliders, J. High Energy Phys. 0307 (2003) 028, arXiv:hep-ph/0306211. [6] B. Simon, Summability methods, the strong asymptotic condition, and unitarity in quantum field theory, Phys. Rev. Lett. 28 (1972) 1145–1146, http://dx.doi.org/10.1103/PhysRevLett.28.1145. [7] C.M. Bender, T.T. Wu, Large order behavior of Perturbation theory, Phys. Rev. Lett. 27 (1971) 461, http://dx.doi.org/10.1103/PhysRevLett.27.461.

272

A. David, G. Passarino / Physics Letters B 726 (2013) 266–272

[8] C.M. Bender, T. Wu, Anharmonic oscillator. 2: A study of perturbation theory in large order, Phys. Rev. D 7 (1973) 1620–1636, http://dx.doi.org/ 10.1103/PhysRevD.7.1620. [9] L. Lipatov, Divergence of perturbation series and pseudoparticles, JETP Lett. 25 (1977) 104–107. [10] A. Bukhvostov, L. Lipatov, High orders of the perturbation theory in scalar electrodynamics, Phys. Lett. B 70 (1977) 48–50, http://dx.doi.org/10.1016/ 0370-2693(77)90341-0. [11] E. Brezin, J.-C. Le Guillou, J. Zinn-Justin, Perturbation theory at large order. 1. The ϕ 2N interaction, Phys. Rev. D 15 (1977) 1544–1557, http://dx.doi.org/ 10.1103/PhysRevD.15.1544. [12] E. Brezin, J.-C. Le Guillou, J. Zinn-Justin, Perturbation theory at large order. 2. Role of the vacuum instability, Phys. Rev. D 15 (1977) 1558–1564, http://dx.doi.org/10.1103/PhysRevD.15.1558. [13] A. Vainshtein, V.I. Zakharov, Ultraviolet renormalon reexamined, Phys. Rev. D 54 (1996) 4039–4048, http://dx.doi.org/10.1103/PhysRevD.54.4039. [14] I. Suslov, Divergent perturbation series, Zh. Eksp. Teor. Fiz. 127 (2005) 1350, http://dx.doi.org/10.1134/1.1995802, arXiv:hep-ph/0510142. [15] J. Fischer, High order behavior and summation methods in perturbative QCD, Acta Phys. Pol. A 27 (1996) 2549–2566, arXiv:hep-ph/9512269. [16] G.H. Hardy, Divergent Series, Clarendon Press, Oxford, 1949, ISBN 978-0-8218-2649-2. [17] G. ’t Hooft, Can we make sense out of quantum chromodynamics? Subnucl. Ser. 15 (1979) 943. [18] V.I. Zakharov, QCD perturbative expansions in large orders, Nucl. Phys. B 385 (1992) 452–480, http://dx.doi.org/10.1016/0550-3213(92)90054-F. [19] A. Sokal, An improvement of Watson’s theorem on Borel summability, J. Math. Phys. 21 (1980) 261–263, http://dx.doi.org/10.1063/1.524408. [20] M. Beneke, V.I. Zakharov, Improving large order perturbative expansions in quantum chromodynamics, Phys. Rev. Lett. 69 (1992) 2472–2474, http://dx.doi. org/10.1103/PhysRevLett.69.2472. [21] J.C. Collins, D.E. Soper, Large order expansion in perturbation theory, Ann. Phys. 112 (1978) 209–234, http://dx.doi.org/10.1016/0003-4916(78)90084-2. [22] M. Beneke, Renormalons, Phys. Rep. 317 (1999) 1–142, http://dx.doi.org/ 10.1016/S0370-1573(98)00130-6, arXiv:hep-ph/9807443. [23] M. Beneke, V.M. Braun, N. Kivel, Large order behavior due to ultraviolet renormalons in QCD, Phys. Lett. B 404 (1997) 315–320, http://dx.doi.org/ 10.1016/S0370-2693(97)00562-5, arXiv:hep-ph/9703389. [24] H. Cohen, F. Villegas, D. Zagier, Convergence acceleration of alternating series, Exp. Math. 9 (1) (2000) 3. [25] C. Brezinski, M. Redivo Zaglia, Extrapolating Methods. Theory and Practice, North-Holland, 1991. [26] C. Brezinski, Acceleration des suites a convergence logarithmique, C. R. Acad. Sci. Paris A 273 (1971) 727–730. [27] A. Sidi, An algorithm for a special case of a generalization of the Richardson extrapolation process, Numer. Math. 38 (1982) 299–307. [28] A. Sidi, Borel Summability and converging factors for some everywhere divergent series, SIAM J. Math. Anal. 17 (5) (1986) 1222–1231. [29] S. Moch, A. Vogt, Higher-order soft corrections to lepton pair and Higgs boson production, Phys. Lett. B 631 (2005) 48–57, http://dx.doi.org/10.1016/ j.physletb.2005.09.061, arXiv:hep-ph/0508265. [30] R.D. Ball, M. Bonvini, S. Forte, S. Marzani, G. Ridolfi, Higgs production in gluon fusion beyond NNLO, arXiv:1303.3590. [31] S. Buehler, A. Lazopoulos, Scale dependence and collinear subtraction terms for Higgs production in gluon fusion at N3 LO, arXiv:1306.2223. [32] E.J. Weniger, Performance of superconvergent perturbation theory, Phys. Rev. A 56 (1997) 5165–5168, http://dx.doi.org/10.1103/PhysRevA.56.5165. [33] U. Jentschura, J. Becher, E. Weniger, G. Soff, Resummation of QED perturbation series by sequence transformations and the prediction of perturbative coefficients, Phys. Rev. Lett. 85 (2000) 2446–2449, http://dx.doi.org/ 10.1103/PhysRevLett.85.2446, arXiv:hep-ph/9911265.

[34] E.J. Weniger, Nonlinear sequence transformations: Computational tools for the acceleration of convergence and the summation of divergent series, arXiv:math/0107080. [35] D. Roy, R. Bhattacharya, Prediction of unknown terms of a sequence and its application to some physical problems, Ann. Phys. 321 (2006) 1483–1523, http://dx.doi.org/10.1016/j.aop.2005.12.010. [36] Skála Cizek Zamastil, New summation technique for rapidly divergent perturbation series. Hydrogen atom in magnetic field, J. Math. Phys. 44 (2003) 962–968. [37] J. Zamastil, Approximate recalculation of the α (Zα )5 contribution to the selfenergy effect on hydrogenic states with a multipole expansion, Ann. Phys. 328 (2013) 139–157, http://dx.doi.org/10.1016/j.aop.2012.09.007. [38] D. Levin, Development of non-linear transformations for improving convergence of sequences, Int. J. Comput. Math. 3 (1973) 371–388. [39] E.J. Weniger, Nonlinear sequence transformations for the acceleration of convergence and the summation of divergent series, Comput. Phys. Rep. 10 (1989) 189–371. [40] E.J. Weniger, Mathematical properties of a new Levin-type sequence transformation introduced by Cizek, Zamastil, and Skala. I. Algebraic theory, J. Math. Phys. 45 (2004) 1209. [41] P. Wynn, On a device for computing the em ( S n ) transformation, Math. Tables Other Aids Comput. 10 (1956) 91–96. [42] C. Brezinski, A general extrapolation algorithm, Numer. Math. 35 (1980) 175–187. [43] A. Sidi, A new method for deriving Pade approximants for some hypergeometric functions, J. Comput. Appl. Math. 7 (1980) 37–40. [44] T. Fessler, W. Ford, D. Smith, HURRY: An acceleration algorithm for scalar sequences and series, ACM Trans. Math. Softw. 9 (1983) 346–354. [45] Maplesoft, Maple Reference: Numerical Summation, http://www.maplesoft. com/support/help/Maple/view.aspx?path=evalf/Sum [online accessed 3 July 2013]. [46] GNU Scientific Library Reference: Series Acceleration, Online, http://www.gnu. org/software/gsl/manual/html_node/Series-Acceleration.html [online accessed 3 July 2013]. [47] W. Ford, D. Smith, Acceleration of linear and logarithmic convergence, SIAM J. Numer. Anal. 16 (1979) 223–240. [48] C. Sturm, Leptonic contributions to the effective electromagnetic coupling at four-loop order in QED, arXiv:1305.0581. [49] P. Baikov, K. Chetyrkin, J. Kuhn, J. Rittinger, Vector correlator in massless QCD at order O (αs4 ) and the QED beta-function at five loop, J. High Energy Phys. 1207 (2012) 017, http://dx.doi.org/10.1007/JHEP07(2012)017, arXiv:1206.1284. [50] A. Kataev, S. Larin, Analytical five-loop expressions for the renormalization group QED β -function in different renormalization schemes, Pisma Zh. Eksp. Teor. Fiz. 96 (2012) 64–67, http://dx.doi.org/10.1134/S0021364012130073, arXiv:1205.2810. [51] D. Roy, R. Bhattacharya, S. Bhowmick, Rational approximants generated by the u-transform, Comput. Phys. Commun. 78 (1993) 29–54. [52] D. Roy, R. Bhattacharya, S. Bhowmick, Rational approximants using Levin– Weniger transforms, Comput. Phys. Commun. 93 (1996) 159–178. [53] D. Roy, R. Bhattacharya, S. Bhowmick, Rational interpolation using Levin– Weniger transforms, Comput. Phys. Commun. 101 (1997) 213–222. [54] D. Roy, R. Bhattacharya, S. Bhowmick, Iterations of convergence accelerating nonlinear transforms, Comput. Phys. Commun. 54 (1989) 31–46. [55] E. Weniger, J. Cizek, Rational approximations for the modified Bessel function of the second kind, Comput. Phys. Commun. 59 (1990) 471–493. [56] E. Weniger, J. Cizek, F. Vinette, Very accurate summation for the infinite coupling limit of the perturbation series expansions of anharmonic oscillators, Phys. Lett. A 156 (1991) 169–174. [57] E. Weniger, Interpolation between sequence transformations, Numer. Algorithms 3 (1992) 447–496.