On interpolating anomalous dimension of twist-two operators with general spins

We study non-perturbative interpolating functions to probe the physics of anomalous dimensions associated with twist-two operators in ${\cal N}=4$ SYM of finite and infinite spin. Compared to previous studies, the novel result of this paper is to introduce single multivariate functions of both coupling $g$ and spin $j$ to approximate such anomalous dimensions. We provide an unified framework to study such operators in interim ranges of the parameters which so far has eluded previous results. Explicitly, we consider twist-two anomalous dimensions in two distinct scenarios using interpolating functions. For the large $N$ case, we stick to simple Pad\'{e} approximants and its generalizations . For the finite $N$ case, ${\cal N}=4$ SYM is expected to be S-dual invariant, hence the observables are expected be modular invariant. To probe the finite $N$ physics, we take into account the non-planar and instanton contributions by constructing modular invariant interpolating functions to approximate the cusp and twist-two anomalous dimensions. We also consider interpolating functions for the twist-four operators and study level crossing phenomenon between the twist-two and twist-four operators.


Introduction
As a close cousin of quantum chromodynamics (QCD), N = 4 super Yang-Mills (SYM) theory has drawn significant attentions in the past twenty years. Unlike QCD which is analytically under control only in the perturbative ultraviolet (UV) regime, in N = 4 SYM, both weak and strong coupling behavior can be understood quantitively, thanks to the AdS/CFT correspondence [1,2,3]. In the 't Hooft large N limit [4], the integrability makes it even possible to study various physical observables non-perturbatively, see e.g. [5] for a review. This makes N = 4 SYM an ideal model to study the non-perturbative dynamics of a four dimensional interacting gauge theory, which hopefully may lead to better understanding in generic gauge theories, including realistic QCD.
A particularly interesting class of observables are the anomalous dimensions of twisttwo Wilson operators, which will be the central objects of our study. In N = 4 SYM, the simplest twist-two operator consists of two scalars written as tr(φD j φ), which has classical dimension ∆ 0 = j + 2 and spin j, and the twist is defined by the difference ∆ 0 − j = 2. The classical dimension receives quantum corrections which is the anomalous dimension γ(g, j), where g is the gauge coupling. In QCD, twist-two operators appear in the operator product expansion description of deep inelastic scattering processes, and they play an important role for determining the parton distribution functions, see e.g. [6].
At least in this UV regime where the twist-two operators contribute to high energy QCD processes, one expects that the two theories, QCD and N = 4 SYM, have many features in common. Let us point out two interesting connections between the two theories. Firstly, in the large spin limit, the twist-two anomalous dimensions in both theories develops a logarithmic behaviour with respect to the large spin dependence, and the coefficient is given by the cusp anomalous dimension (CAD) γ cusp [7,8], i.e. γ(g, j) ∼ Γ cusp (g) log(j) , for j 1 and j g . (1.1) Another intriguing connection is the so-called "maximal transcendentality principle" which was first observed in [9,10]. It says that the anomalous dimensions of twist-two operators in N = 4 SYM can be obtained from the maximally transcendental part of the QCD results [11]. (Here, "transcendentality" refers to transcendentality degree which characterizes the "complexity" of mathematical numbers or functions, for example, the Riemann zeta value ζ n and the polylogarithm Li n have degree n.) This is a conjecture expected to be true perturbatively to all orders. Further evidences of these correspondences were also found for other observables such as amplitudes and Wilson lines [12,13,14,15,16]. These little understood mystical correspondence between N = 4 SYM and QCD indicates there is a deep connection between the two theories.
While the non-perturbative QCD is an extremely hard question to address, the goal of the present paper is to study non-perturbative interpolating functions to approximate the twist-two anomalous dimensions in N = 4 SYM. We hope such a study might help us to probe the real physics in QCD, as also indicated by the aforementioned connections. Similar studies of using interpolating functions have been considered in [17,18,19,20,21,22]. Compared to previous studies, the main new aspect of this work is that under a single multivariate function we can analytically consider the anomalous dimension as a function of both the coupling g and spin j. This provides a unified framework for many previous results. For example, the dependence on spin allows us to incorporate the information of the cusp anomalous dimension as the large spin limit of the constructed functions. Explicitly, we will consider twist-two anomalous dimensions in two distinct scenarios, using two distinct classes of interpolating functions.
The first scenario is the planar large N limit, where N = 4 SYM is integrable and one can employ the tools of AdS/CFT duality to explore non-trivial gauge invariant objects. This case is relatively well understood due the significant progress made in past years, and in principle, non-perturbative results can be obtained via integrability techniques [5]. In N = 4 SYM, one could think of twist-two operators tr(φD j φ) as adding large number of derivative impurities to the protected half-BPS vacuum tr(φ 2 ), which forms the backbone of various integrability related investigations for such operators and their dual string states. While such a closed form for twist-two anomalous dimensions for arbitrary spin is not yet known, perturbative expansion in various limits are known to higher orders. These data allows us to construct relatively reliable interpolating functions, from which one can study the non-perturbative properties in regimes unreachable herein-before. One interesting feature is that the anomalous dimensions in the small spin (j g) and large spin (j g) limits have very different behaviours. In this work we will use the framework of interpolating functions to explicitly capture such disparate physical phenomena and discuss the consequences.
The second scenario is the case of finite N , which physically is more closely related to realistic QCD (where N = 3). This case is much harder to study, because the theory is no longer integrable and much less data are available. On the other hand, an important new ingredient in this case is the S-duality property of N = 4 SYM [23]. (Note that we do not expect S-duality for twist-two anomalous dimensions in the large N limit.) This suggests use of interpolating functions that are invariant understand SL(2, Z) modular transformations providing strong constraints on the result. Furthermore, since we will be using a basis of modular invariant Eisenstein series (instead of polynomials) for our interpolating functions, both the instanton corrections as well as the non-planar data can be incorporated in our unified framework. Compared to the large N case, a subtle new physical feature expected to appear at finite N is the level crossing phenomenon, see e.g. [24]. We will try to address this issue based on the interpolating function. Since the available finite N data is limited, the main goal here is to have a qualitative or even speculative understanding of the physical picture. Hopefully this can provide new insight for (or be tested by) further studies in this regime. Besides, we will also encounter several technical challenges, such as getting the correct coupling powers in the strong coupling expansion and encoding the non-planar and instanton contributions in the modular invariant functions, which we will explicitly address in our construction.
The structure of the paper is planned as follows. We first review the construction of various classes of interpolating functions in section 2. Then we explicitly compile the data available in the literature for cusp and finite spin twist-two anomalous dimensions in section 3, which will serve as input constraints for the interpolating functions in the subsequent sections. In section 4, we construct the interpolating function for the twist-two anomalous dimensions in the planar limit. We elaborate on both large and small spin physics here based on construction of two-point approximants. In section 5, we turn to the theory with finite value of N . We first discuss the interpolating function of cusp anomalous dimension with S-duality, which is related to the large spin physics. Then we discuss the interpolating function for the anomalous dimensions with finite spin. We discuss the results and present an outlook for our constructions in section 6. Finally, we give several appendices covering the technical details the construction of interpolating functions.

Interpolating functions: Construction
Before going into the details of our interpolating functions for the anomalous dimensions of twist-two operators, let us first briefly discuss some well known techniques employed to study non-perturbative interpolations via approximants. In literature there exists a variety of sophisticated techniques to resum perturbative expansions (both strong coupling and weak coupling) to get non-perturbative answers in quantum mechanics and quantum field theory. If the perturbative series of the theory under consideration is Borel summable [25], powerful techniques can be employed to resum the perturbative series. So it should be possible to obtain a great deal of insight into the behaviour of the resummed function associated to an observable at any finite coupling with some additional information about the strong coupling expansion along with the perturbative series. For us, the physical observable would be the anomalous dimensions of twist-two operators.
Although it is always theoretically possible to write a non-perturbative function which encodes the perturbative expansions of such anomalous dimensions, in general it becomes quite hard to encode extra symmetries such as those implied by duality in N = 4 SYM on the functions themselves. It would be interesting to find appropriate resummation methods which are not only compatible with the asymptotic behaviour of perturbation theory, but also clearly encodes the features of duality symmetry. It should be noted that such dualities are extremely valuable to constrain the structure of the anomalous dimensions. As a starting point to discuss such non-perturbative duality invariant results, one could explicitly construct interpolating functions which are inherently invariant under such dualities.
Motivated by this, in this work we construct interpolating functions to study the nonperturbative properties of twist-two anomalous dimension in N = 4 SYM both at finite N and large N limits. In the former case where it is expected that anomalous dimensions will be invariant under the action of the full SL(2, Z) modular group, we construct interpolating functions that are invariant under modular transformations.
In the later case, i.e., in the large N limit where the observables are not invariant under the symmetries imposed by S-duality, 1 we would construct simple interpolating functions consistent with the weak coupling and the strong coupling expansions.
To this effect, we will make use of distinct classes of interpolating functions constructed in the literature [17,18,19,20,22].
Broadly speaking, we will be employing two different approaches to approximate the anomalous dimensions for leading twist two operators both at finite and large spin limits.
• Interpolating functions without S-duality

• S-duality invariant interpolating functions
In what follows we will briefly describe various classes of such interpolating functions which will be important to us and build up to the machinery we extensively use in the later sections.

Interpolating functions without S-duality
As mentioned before, in the large N limit it is expected that the observables i.e, the anomalous dimensions of twist-two operators in N = 4 SYM are not invariant under the full modular group. We provide the methodology of constructing such interpolating functions which are consistent with both the weak coupling and the strong coupling expansions without invoking any modular invariance.
The idea of an interpolating function is to broadly know the perturbative results at strong and weak coupling and match the results from both ends. Suppose we want to interpolate a function F (g) which has the weak coupling expansion F (Nw) w (g) up to (a + N w )th order around g = 0 and strong coupling expansion F Then for a consistent interpolation we expect that the Taylor expansion of the interpolating function will match these two expansions around the weak and strong coupling, In terms of these expansions, we would like to construct smooth interpolating function which coincides with the small-g and large-g expansions up to the given orders.

Padé approximant:
A simple possibility for an interpolating function is the two-point Padé approximant. Let us construct the Padé approximant P m,n (g) for the function F (g), with m ≤ N w and n ≤ N s . The formal definition of the Padé approximant interpolating function for b − a ∈ Z is given by The coefficients c k and d k in (2.4) can be fixed such that power series expansions of the Padé approximation, P m,n (g) around g = 0 and g = ∞ agrees with the weak coupling expansion and strong coupling expansions up to the given order in the perturbation up to O(g a+m+1 ) and O(g b−n−1 ) respectively. Thus by construction the Padé approximant would satisfy F (g) − P m,n (g) = O(g a+m+1 , g b−n−1 ).
For this construction we need Padé approximants are widely used to construct non-perturbative answers for perturbatively known functions in diverse areas of physics. However in some situations the denominator in (2.4) runs into zeroes in physically interesting regions and poles show up in the total function. This situation signals limitation of approximation by the Padé, and except for cases where F (g) itself has poles, it would become necessary to investigate the radius of convergence of the Padé approximant.

Fractional Power of Polynomial (FPP):
In [17] the author constructed a new type of interpolating function, which we refer to as the fractional power of polynomial method (FPP), having the following skeleton structure, (2.7) We can determine the coefficients c k and d k in a similar fashion as explained in the case of Padé approximant. By construction, the FPP would satisfy The FPP does not have constraint on the parameters such as (2.6) in the Padé approximant. But as in the previous case, these functions are also not free from running into non-analytic regions. One can encounter cases where the polynomial itself can become negative in physically important regions. Consequently, when the power (b − a)/(m + n + 1) is not integer, the FPP takes complex value and signals a breakdown of the approximation.
Although we will use this FPP in its current form only to approximate functions which are not expected to be invariant under duality symmetries, it is interesting to note that this method has been applied to capture S-duality improvements of string perturbation theory [17,26]. 2

Fractional Power of Rational function (FPR):
Based on the earlier two constructions, one can construct a more general class of interpolating function [20] with the basic structure as follows where the parameters are We can easily deduce that the Padé and FPP both are special cases of the above interpolating function, by taking the following special limits of (2.8): • If 2 + 1 = a − b for a − b ∈ Z and m + n is even, (2.8) reduces to Padé approximant.
The coefficients c k and d k could be determined in a similar way as explained in the Padé and FPP. We refer to this interpolating function as "fractional power of rational function method" (FPR). To construct this interpolating function we require that for m + n : odd , with ∈ Z. (2.11) In analogy to the case of FPP, when the rational function in the parenthesis has poles or takes negative values for non-integer α, we cannot trust approximation of the function F (g) by the FPR.

Interpolating functions with S-duality
In the previous section we have introduced a class of interpolating functions which are not inherently invariant under any duality transformations. As mentioned earlier, one can use, for example FPP, to approximate functions invariant under S-duality [17,26]. The strategy in these cases is to demand S-duality invariance for the whole function, which in turn gives strong constraints on the coefficients of the polynomials. However, this is not enough to handle functions with all possible non-perturbative effects (like instantons for example) since these contributions cannot be taken care of using polynomials. In our approach, we will be using a better guiding principle to ascertain S-duality in the function. Now to study the non-perturbative properties of the anomalous dimensions at finite N case, it is expected that S-duality would play an important role. This would imply that observables should transform appropriately under the full modular transformations of coupling parameters which is a combination of S-and T -transformations: where τ is the complex gauge coupling . (2.14) In the case of finite N , we will mostly focus on the anomalous dimensions of leadingtwist operators and construct interpolating functions which satisfy full S-duality invariance. To construct such a function, the basic philosophy is to choose inherently modular invariant building blocks, instead of polynomials in g, as expansion basis. One such natural choice is the real or non-holomorphic Eisenstein series 3 Since the Eisenstein series are invariant under the duality transformation (2.12), by construction interpolating functions constructed out of Eisenstein series as the basic building blocks are invariant under the full S-duality. The non-holomorphic Eisenstein series E s (τ ) has the following expansion for small argument (see e.g. section 5.3 of [27]) where σ s (k) is the divisor function defined by 4 σ s (k) = d|k d s . In terms of (g, θ), it is written as Here the Γ is the usual gamma function and K is the modified bessel function. Below we present a brief account of the interpolating function methods developed using such Eisenstein series in [19] and [22]. We note an important feature of the Eisenstein's E s (τ ), it contains both a perturbative and a non-perturbative part. The non-perturbative part of the Eisenstein series contains power of q = e 2πiτ , which would play a crucial role in reproducing the correct instanton contributions, as we will use later for the cusp anomalous dimension.

FPP-like interpolating function involving Eisenstein series:
In [19] the following type of interpolating function has been constructed where the coefficient c k 's are determined such that expansion ofF (s) m around g = 0 agrees with the weak coupling expansion of the anomalous dimension, F (τ ) up to O(g m+1 ). One could see that the above function has structural similarities with FPP like interpolating functions (2.7), where instead of a polynomial in g, the Eisenstein series has been used. Thus an appropriate choice of c k correctly gives the weak coupling expansion of F (τ ). However, since this function is not actually constrained by strong coupling expansion, it is only natural to consider generalisations of it where strong coupling data has significant role to play.

FPR-like duality invariant interpolating function:
In section 5.1 we will construct FPR-like duality invariant interpolating functions to study the cusp anomalous dimensions and the anomalous dimensions for finite spin operators at finite N . Here, we give a brief methodology of building such generalised function and will provide a more detailed construction procedures in the later section. The FPR-like duality invariant interpolating functions constructed in [22] has the following structure: where we can determine the coefficients c k and d k (except d 1 ) such that expansion ofF (s,α) m around g = 0 agrees with F (τ ) up to O(g m+1 ). 5 Matching at O(g) leads to Now, the interpolating function is invariant under the scaling c k , d k → λc k , λd k , so without any loss of generality we can take Matching at other orders leads to the constraint and hence we find the relation between the parameters (2.23) 5 Note that m should be m ≥ 2 since we need two coefficients at least for this interpolating function.
We also require p, q ∈ Z ≥1 , which implies α = 1 2 for m : odd 1 2 +1 for m : even , with ∈ Z. (2.24) The above interpolating function can be further constrained if we impose the matching of coefficients from either the strong coupling or the finite N results. Let us discuss the constraints imposed by strong coupling results with O( 1 N 2 ) corrections on the interpolating functions. We start with 't Hooft expansion of the interpolating function: 6 Then we determine the yet unknown coefficient d 1 to satisfy is the result in the supergravity limit given by (3.23). Imposing matching of other orders leads us to p + q − 2 = m, (2.27) and therefore we get We also require p, q ∈ Z ≥1 , which constrains α as α = 1 2 +1 for m : odd 1 2 for m : even , with ∈ Z. (2.29) The constraints coming from the supergravity, where the data is of the same form as above, will be important to construct the interpolating function for the finite spin twisttwo operators. Note that there are clearly three parameters (m, s, α) driving the interpolating function. There could be infinite choices for this set of parameters, leading to infinite number of possible interpolating functions. This "landscape problem of interpolating functions" was studied earlier in [20]. It is a priori unclear which set of values of (m, s, α) would give us the best approximation. We will briefly discuss the procedure to choose optimal values of (m, s, α) in our construction, see also [20,22].
• Choice of m.
By construction the interpolating function should reproduce the correct weak coupling expansion up to the given m-th loop order. The best choice for m depends on the details of the weak coupling expansion. One of the important criteria would be the convergence property of the weak coupling expansion. Suppose the weak coupling expansion is convergent series then we can take m as large as possible otherwise we have to judiciously select m.
• Choice of s. Since the weak coupling expansion of the function of interest only contains positive integer powers of g, however in principle the interpolating function (2.19) can contain fractional power of g for any arbitrary value of s. In order to guarantee absence of such fractional powers, we should take 7 2s ∈ Z.
As discussed in [22], most of the constructions has little dependence on s and thus we get an infinite class of interpolating functions with extremely close numerical values but different structures.
The parameter α determines the type of branch cuts of the interpolating functions. Hence a correct choice of α would be related to the analytic properties of the interpolating functions.
We will provide more details of the choice of such parameters in section 5 where we construct such functions explicitly. Finally, let us mention an issue of matching strong coupling expansion which would require a further generalization of the above function. An important feature of the above function is that in the limit λ 1 (planar limit), since the non-perturbative part of the Eisenstein series is suppressed (2.17), the coefficients c k and d k are determined only in terms of the perturbative part i.e. O(g −s ) and O(g s−1 ) terms. Furthermore, in the regime discussed in [22], at large s only the O(g −s ) part is relevant and in the planar limit the total function have a schematic form 8 where we havec (2.31) In the large λ limit,c k andd k are essentially O(1). Furthermore, the function in-principle can't generate any fractional powers in λ in the strong coupling limit. In order to take into account any fractional powers of λ in the strong coupling limit we have to further generalize this function, which we defer for a details discussion in section 5.1.

Data on anomalous dimensions
There has been many well known investigations on twist-two operators in the relevant literature. In this section, we will very briefly review some aspects of these investigations 7 Note that s can be non-integer and Es(τ ) has a pole at s = 1 hence we take s > 1. 8 We would expect it to be O(N −(s+k) ) but since we normalize ds+q = 1, we multiply each coefficient by and mention the main results of such works, which will serve as effective input data in our work. We summarize our convention for the coupling constants here: where the g Y M is the Yang-Mills coupling constant and λ is the 't Hooft coupling. The 't Hooft large N limit is taken by N → ∞, keeping λ org constant. For the study of finite N case, we also introduce

Results on cusp anomalous dimension
Cusp anomalous dimension is an important observable that governs the universal scaling behavior of various gauge invariant quantities. As we already mentioned, it governs the large spin scaling behavior of twist-two anomalous dimensions. By definition it is also the anomalous dimension of Wilson loop with a light-like cusp singularity [7,8]. Furthermore, it provides the leading infrared divergences of on-shell amplitudes and is an essential ingredient in constructing amplitudes, such as in [30]. In AdS/CFT correspondence, it is related to dual description in terms of spinning strings [31] or cusped minimal surface [32] in the AdS background.
In the planar limit, in principle we can find the weak and strong coupling expansions up to any loop order for cusp anomalous dimension using the BES formula [33]. This is an integral equation derived from all-loop Bethe Ansatz equation with a mathematically complicated kernel structure. At weak coupling, the planar expansion for cusp anomalous dimension has been obtained up to four loops from rigorous perturbative analysis [34,35,36]. From semiclassical computations in string theory the strong coupling expansion has been explicitly computed up to two loops from the analysis of quantum string sigma model in AdS [37,38].
On the other hand, non-planar corrections to quantities like cusp anomalous dimension is so far hard to compute within a framework like AdS/CF T , where these corrections correspond to string loop corrections. Also, the power of integrability fails here since it's role beyond planar limit is yet to be fully uncovered. Recently, progress in computing the non-planar corrections to cusp anomalous dimension has been made via a numerical calculation of Sudakov form factor [39,40], where the non-planar part enters into the result first in the fourth loop order in weak coupling.
Let us summarize the result of cusp anomalous dimension up to four loops: where the the non-planar four-loop cusp anomalous dimension, Γ np , is given by 9 Γ np ∼ −2400.
(3.4) 9 We use the central value of the non-planar result [39]. Note that the definition of cusp anomalous dimension is different from that in [39] by an overall factor 2.
At strong coupling, we quote the result [37] Γ cusp,s = 2g − 3 log 2 2π As mentioned before, higher order data can be extracted from the BES equation recursively. For example, in the next orders of weak and strong coupling expansions, one can read off, Γ (5) cusp,w = 16 where K is the Catalan constant. This allows us in principle to predict the planar cusp anomalous dimension perturbatively to arbitrary loop order [41].
One important motivation of this paper is to consider truly non-perturbative corrections to the anomalous dimensions. In a recent work [42], leading instanton contribution to the light-like cusp anomalous dimension has been computed. We briefly review this below. In general quantum corrections to four point functions of half-BPS operators have the following form in the weak coupling limit, where the first term is the perturbative part and the other term is the non-perturbative correction and u, v are the two cross ratios. In [42] it is mentioned that in large spin limit the leading instanton contribution scales as This result is valid for the SU(2) gauge group. Remember that in the light-like limit the cross ratios u, v → 0 and this reproduces the contribution for the cusp anomalous dimension, as have been quoted above. For the one instanton correction, the contribution has the following form Now we may generalize this instanton corrections to the SU (N ) group as follows. The non-perturbative correction appears in the non-planar sector so there is a explicit factor of 1/2 2 for the SU (2) case. We can rewrite (3.9) 10 as follows (3.10) In the large N limit, the result for the leading instanton correction gets multiplied by a N dependent factor [42] due to contribution from all bosonic and fermionic modes arising 10 Here we have assumed that the leading instanton corrections starts at 1/N 2 order [43].
from embedding the SU (2) instanton into SU (N ). The appearance of such a factor has been worked out in details in [44] and in our case it has a form For large N , one can see that the total factor before the instanton contribution then has a factor O(N −3/2 ).

Data for finite spin anomalous dimensions
Next, we consider the anomalous dimension of twist-two operator with generic spins, O j = Tr (φD j φ). The scaling dimension of these operators can be written schematically as, where γ(g) is the anomalous dimension. Using conventional methods, one can study these functions only up to first few loops in the strong/weak coupling regimes. With integrability techniques, employing the TBA or Y-system [45], one gets analytically very complicated integral equations. With the advent of techniques associated to Quantum Spectral Curve (QSC) [46,47,48,49], in principle one can extract such data at any value of the spin and coupling constant. However, such an explicit computation in particular for generic spin dependence still seems to be absent in the present literature. As mentioned earlier, the results for anomalous dimensions in N = 4 SYM can be extracted from QCD results by isolating the maximally transcendental part perturbatively at different orders in the coupling constant [9,50,51]. The internal symmetries of N = 4 SYM makes it evident that the basic building block of such anomalous dimension are sums of the form i 1 j i , where i is the level of transcendentality, indicating that these anomalous dimensions are polynomials in Riemann Zeta value or it multi-index generalisations. The basis for the results is formed from these harmonic sums defined as follows, Due to this remarkable structure, we can define the anomalous dimension for finite spin twist-two operators using a basis of harmonic sums. We can now write the anomalous dimension up to three loops as whereg has been defined as before. This γ is the anomalous dimension defined at finite values of j.
Using the above assumptions, the leading order (LO) and the next-to-leading order (NLO) anomalous dimensions for twist-two operators were found in [52,53]. The three loop expressions were obtained in [10] by extracting the most complicated contributions from the three loop non-singlet anomalous dimensions in QCD [11]. To keep our considerations simpler, we will use the expression for first few loop data as follows, 11 The j → ∞ results are important here, especially for matching to the cusp anomalous dimension. From each of these terms we will get a log(j) contribution in the j → ∞ limit, the coefficients of which can be exactly matched to the ones written in equation (3.3). So, in the weak coupling case, the finite spin data smoothly connects to the infinite spin case and gives rise to the log behaviour.
On the other hand, we are also interested in the coefficients of the strong coupling expansion of the anomalous dimension at finite spin values. In planar limit, this can be (partially) computed from AdS/CFT prescription and in principle via integrability. Let us call the strong coupling anomalous dimension G(g), For the twist-operators in the SL(2) sector, there are analytical predictions for the first four coefficients of (3.18). The coefficients for twist-two operators can be written as functions of j and is known to take the following form via Quantum Spectral Curve calculation [47]: Our conventions are the same as of [10]. Higher loop data for twist two operators have been explored in [54,55,56,57]. We note that the non-planar corrections were computed for spin 2, 4, 6 and 8 [58,59,60], but for generic spin the results are still not available.
The first two coefficients in (3.19) can be determined either from Basso's slope function [61] or from semi-classical computations in string theory [62,63,64]. The next two coefficients were determined by matching the O(j 2 ) term of the small spin expansion with classical and semi-classical results [47].
In the small spin regime the data can be compared with the dispersion relation of a small-spin circular string moving near the centre of AdS. The other regime of the strong coupling data, i.e. at j → ∞ limit is given by the folded string dispersion relation. It is important to note here that unlike weak coupling, it is not possible to reproduce the log(j) behaviour from the strong coupling data (3.19)-(3.21) at large spin limit. This would indicate a non-trivial "phase-transition" -like physics in large spin regime when we go to the dual AdS string picture. We will come back to this picture in section 4.3.
At finite N , another important ingredient that will enter our calculations at strong coupling is the double-trace twist-four operator from operator mixing effect. At planar large 't Hooft coupling, the anomalous dimension of the the twist-two operator O M (3.18) grows without any bound [2], and mixing effect with double-trace operators kicks in. We will discuss more in section 5 when constructing interpolating functions at finite N . Here, let us introduce the double-trace operators which has a schematic form where tr(φ (i φ j) ) is the symmetric traceless part of tr(φ i φ j ) and a chiral primary operator belonging to 20 representation of SU (4) R . There are known results for the anomalous dimension of such operators from supergravity [65,66,67] computation in AdS 5 × S 5 . Below we note down the the anomalous dimension for spin j double trace operator to O(1/N 4 ) order as considered in [68] and [69] (3.23) Here we have defined anomalous dimension of the double trace operator as γ(j) = ∆ s − (j + 2), where ∆ s is the scaling dimension of the double trace operator.

The planar case: Finite spin twist-two operators
In this section, we would try to understand the physics of anomalous dimension of "short" (small spin) and "long" (large spin) twist-two operators from the interpolating function point of view in the large N limit. In what follows, we would try to address this issue using simple Padé type approximants as toy models. Later we will propose betterment over such simple approximants and try to see whether this improves the interpolation. As we will see, the underlying proposed physical picture is fairly independent of the construction itself. The 't Hooft coupling is denoted byg or equivalently λ, and we again recall our convention: Figure 1: A visualisation of the data in theg-j plane as presented above.

Results from Padé approximant
In the last subsection, we discussed various regimes of available data for anomalous dimensions of twist-two operators. Let us list them in one place in the following simplistic way, A. Smallg and small j (can also connect to large j): Data given in terms of Harmonic Sums as function of j as in (3.15) D. Largeg and large j: Data given by cusp anomalous dimension (coefficient of log(j)) in the strong coupling expansion (3.5), also from AdS/CFT via folded strings.
E. Largeg and small j (which can also connect to large j): Data from results of string computations which take care of small spins, i.e. "short strings".
In short, we now know input data for four corners in ag-j parameter space but we don't know the expressions for anomalous dimensions in the interim regime of parameters. In principle, one might be able to construct a master interpolating function along both g and j direction to extract data at any point in this space, but that turns out to be a difficult job. In what follows, we would try to take a different route. By constructing two interpolating functions in the coupling constant, one for the small j region and one for the large j region, we would try to piecewise fit the total parameter space. It is quite evident that the intersection line for such two functions would signify a change of physical regimes.
We start with constructing an interpolating function in the large j region. In this case, we can approximate the anomalous dimension by Γ cusp (g) log j, and it suffices to construct an approximant for only Γ cusp (g). Using the data given in (3.1), we can construct the following two-point minimal Padé approximant in the tune of eq. (2.4), 12 In figure 2(a) , we show the smooth interpolation resulting from this approximant. We next move on to construct another Padé approximant for the small j region using the anomalous data up to three loops (weak coupling as in (3.15)) and QSC data up to four loops (strong coupling (3.18)). Let us call this function G 8/7 for simplicity whose explicit is given in the Appendix A. In figure 2(b), we compare the interpolating function with the perturbative data for j = 5 to showcase the efficiency of interpolation.
As we have seen earlier, in general the small j, smallg data (Section A in figure 1) can also be interpolated to the large j region via the property of the Harmonic Sums. But our function G 8/7 also includes strong coupling data from semi-classical and QSC computations, so we can't a priori assume that the whole function will also be valid into the large j region. We can however expect that up to finite values of j, the approximant works quite well. However, for larger values of j andg there is a competition between the two regimes of data as a dependence in order of limits kicks in. Physically, it seems plausible that in j g region the log j dominates and ing j region the finite spin behaviour (starting with √ j) dominates. We will discuss this more in section 4.3. Now, we have two different models to fit two different j regimes in theg-j parameter space, i.e. Γ (G 6/5 ) cusp log(j) and G 8/7 . We can plot them together and discuss the implications. 12 Here, only the first four loops in the weak coupling and first two loops in the strong coupling side has been taken into account. For related construction, see for example [70].
It is indeed notable that there is a well defined 'sharp' transition region between the two functions, given by the equation In figure 3(a), we show the real solution of j = j(g) from this equation (apart fromg = 0). A comparison with figure 1 can assure the reader that this line roughly separates regions (A+B) and (C+D), while maintaining a sharp transition along a j ∝g line for higher values ofg. One can see from the figure that this curve splits the parameter space into two, and it can be checked that in large j, in the left region value of G 8/7 dominates while in right region Γ The line that appears in the small j region is, however, not trustworthy.
We also plot the whole parameter space in figure 3(b) to show clearly the features of two dimensional surfaces corresponding to the functions and their intersection region. In the figure, it is evident that the two regions are separated by the almost linear transition region i.e. the transition occurs along a j ∝g curve, although in the plot we can mostly see this line in large enough j andg region.
We must note here that by our construction G 8/7 should be well suited for small j region for allg and Γ (G 6/5 ) cusp log(j) should be valid for predominantly large j region. But it seems that due to the logarithmic term, the latter takes over in the parameter space sooner than expected and hence the intersection becomes important. In figure 3(b), we can clearly see that around the intersection, our expectation of log j dominating in j g region and √ j dominating in theg j region is fulfilled by the subdominant branches of the two-dimensional surfaces, i.e. the ones drawn in solid colour. Now due to the extra line in the small j region, it seems Γ (G 6/5 ) cusp (g) log(j) is the right choice in the region with low j but highg (analogous to section B in figure 1), but this appears to be counter-intuitive. We conclude that in small j region, our approximation may not be trustable. Before we go into the physical implications of the above constructions, let us first see how we can in principle improve this interpolating function and show that j ∝g transition behaviour is universal even with improvements.

Improving the interpolating functions
A notable issue with construction of interpolating functions for a physical object is that there could be various ways one could improve the behaviour of such a function. 13 For example, one could take data up to more loops or even construct a different approximant that gives correct expansion at two ends. We have already encountered a number of different interpolating functions that could do the trick. However, our skeleton construction can 13 In general, the structure of large spin expansion of the cusp anomalous dimension follows the expression be dubbed sufficient if, with various improvements, the physical properties captured by the functions do not change much. In what follows we will briefly talk about possible improvements to the construction we presented in the last section.

Improvement via FPP construction
For the sake of completeness we should try to construct better interpolating functions for our small spin and large spin data. One example could be the Fractional Power Polynomial (FPP) type of interpolating function, which has been briefly introduced in section 2.1. We remind the reader again the structure of such function. For our case with small spin data, we can see that the parameters of FPP can be fixed as following, b − a = −3, m + n + 1 = 15. (4.4) Then the total interpolating function for the small spin case will have the form where we write an FPP for cusp anomalous dimension as The problem with such a function is that the denominator has a fractional power, so depending on the coefficients it can run into complex values, which will signal the breakdown of our approximation. One can in principle solve the denominator for negative values to find the breakdown parameters. Careful investigation shows that such breakdown of the approximations in F small occur here not before j ∼ 300 which is a improvement over our model with Padé approximants. A plot of the two functions covering theg-j plane is shown in figure 4. One may note the qualitative similarities between this picture and figure 3(b).

Improvement via nonlinear variable transformations
Another problem with both Padé and FPP type approximants lie in the appearance of anomalous powers in the expansion along weak or strong coupling points. This becomes particularly clear as we expand our Padé approximants or FPP's beyond the orders up to which we have known data. For example, if one fixes the expansion of such a function uptõ g 8 , and tries to expand upto higher orders, anomalous powers likeg 9 can also creep in.
These unwanted powers can be in principle managed by adding "offset" interpolating functions to cancel them up to arbitrary orders, however it does not actually provide a permanent solution. A tentatively better way is to tweak the variables in the interpolating functions so that it only spews out the right powers in both weak and strong coupling expansion. As an example let us consider the variable change 14 Then, we can consider an interpolating function of the form the expansion of which only gives terms of order {g 2 ,g 4 ,g 6 ...} in the weak coupling and {g 1/2 ,g −1/2 ,g −3/2 ...} in the strong coupling expansion, i.e exactly reproduces our small spin data. Similarly, the following function can also be considered, This gives {g 2 ,g 4 ,g 6 ...} in the weak coupling and {g,g 0 ,g −1 ...} in the strong coupling expansion, i.e exactly reproduces our cusp data. One should note that, similar problem with generation of right powers will reappear in the case of finite N (section 5) too, but will be dealt with in a different way.

Towards an interpretation of the "transition region"
As we mentioned in the last section, in the weak coupling, the twist-two anomalous dimension has a distinct log(j) scaling behaviour in the large spin limit to all perturbative orders. On the other hand, the story at strong coupling has two different cases, making it more subtle. When the spin is very large the log(j) scaling is still true. However, when we have small spins, the scaling becomes γ ∼ √ j, see (3.19). Intriguingly, the dependence on the 't Hooft coupling is also different in these two cases. The former case has √ λ dependence at the leading order, while the latter is λ 1/4 dependent. To summarize, we have γ(j, λ) ∼ √ λ log(j), Large spins, (4.11) γ(j, λ) ∼ λ 1/4 √ j, Small spins. (4.12) Obviously, there is an order of limits issue when taking large λ and large j together. The region where such a transition from one behaviour to another happens can be characterized by the intersection line from eq. (4.2) obtained in the previous section. Below we will try to give a physical account for the dependence of this two limits and comment on the physics of the transition. The best way to describe this physics seems to be in the realm of classical strings in AdS/CF T , where large λ finds a natural existence. The general expansion of string energy calculated from such prescription then appears as a dispersion relation between conserved charges from the sigma model.
On the classical string side, the appearance of λ 1/4 scaling is not that surprising, even though one would expect that leading order energy contribution from string sigma model would start at √ λ. Since spin of string roughly counters the contracting effect due to the string tension (T ∼ √ λ), strings with smaller spins actually correspond to smaller lengths for strings with a centre of mass near AdS centre (ρ = 0). In [72,73], it has been argued that one could then define an effective parameter for the case of the spinning folded string as J = j/ √ λ, which actually measures this interplay of spin and string tension and remains fixed even when j → ∞ and λ → ∞.
So, in general there is a case of order of limits involved in this energy expansion. Formally, a string with J 1 (j √ λ) can still be described to be in the 'short string' phase. Following the computation of [72,73] for energies of such string states, we can write the desired dispersion relation E = aJ 1/2 + bJ 3/2 + cJ 5/2 + ... where evidently A(j) goes as √ j. This expression can easily be checked on physical grounds to be low lying excited string states and always behave as m 2 ∼ 1 α , i.e., they have the flat space Regge trajectory behaviour and hence the dimension of dual operators go as λ 1/4 in the leading order. 15 One has to include subleading corrections starting from O(λ −1/4 ) due to curvature of the target space.
On the other hand, it is well known that to capture the physics of cusp anomalous dimensions we need to take the large spin limit, where j √ λ 1 [31], and the dispersion relation clearly gives rise to the √ λ log j term in the leading order. So, it can be clearly seen that the parameter space is divided into two regions along the line j ∼ √ λ. But we must remind ourselves that this works only for large enough value of the coupling, where the string description is still valid. It seems that at least in the large enough coupling region, the transition region between our two interpolating functions are in tune with this description. The situation has been summarised in figure 5. Comparing with figure 3(b) and figure 4, we could see that the solid colour surfaces, i.e. the subdominant ones are the physically trustable regions.
Physically, these two regions on the parameter space signify two ends of the AdS spinning string spectrum, the 'long' string that almost touches the boundary and the 'short' string which is is not stretched much compared to the radius of curvature of AdS and stays near the centre. Both of these pictures are defined in the strong coupling regime, but as we mentioned earlier, the physics depends on the order of limits between the coupling 15 We may also recall that λ = R 4 α 2 ∼ R 4 l 4 s , where R is the AdS radius, and ls is the string length. For large λ and for small spin states, it is equivalent to think that the AdS radius (which we set to be one) is much large comparable to the string length. Thus we can approximate the close string as moving in the flat background which gives the well-known Regge trajectory behavior between mass and spin, m ∼ √ j/α . and the spin. The two regimes of classical strings are discriminated by the appearance of "cusps" or "folds" on the string profile. As we increase the spin and go towards j ∼ √ λ, the string starts to get more and more stretched along the AdS radius, and proper length of this highly excited string starts becoming comparable to AdS radius.
A nice way to see this would be to consider the profile and conserved charges associated to spinning strings in the AdS 3 ∈ AdS 5 spaces. For example the profile of a string spinning with angular velocity ω is given by 16 where ρ is the radial direction in the AdS 3 and σ is the spacelike worldsheet coordinate. We parameterize the AdS coordinates in the form t ∼ τ, φ ∼ ωτ, ρ = ρ(σ) and consider a spinning folded closed string whose center lies at rest at ρ = 0. The extent of the string ρ 0 along AdS radial direction is roughly given by coth ρ 0 = ω. It is evident that for a closed string with ω 2 > 1 and sinh 2 ρ < 1 ω 2 −1 , we can have two cases: 1. ω → ∞: No derivative discontinuities appear on ρ(σ), i.e. no cusps appear on the string, which corresponds to the "short" string phase.
2. ω → 1: ρ(σ) develops a derivative discontinuity at particular point(s), i.e cusps begin to appear on the string, which corresponds to the "long" string phase.
For actual ω = 1, the cusps on the string touches the AdS boundary. So, in general the information of both the short and long strings can be explored from the spinning string 16 See [74] for a recent take on this issue. setup at large coupling. The Noether charges associated to the string can be calculated using (4.15) and is given by [ where E and K are the usual complete Elliptic integrals. One can see here that the charges implicitly depend on each other, and hence the dispersion relation is given by E = E(j). A way to look at the anomalous dimension for any j is to consider the quantity γ = π 2 √ λ (E−j). Although it is quite complicated to evaluate this quantity for different j and λ, one can try to do so numerically. In figure 5 we plot this quantity from the equations above. We see that no such sharp transition occurs in this case as in the last, and the small j and large j regimes interpolate smoothly into each other. From the interpolating function point of view, if one could properly construct a master interpolating function in both j and g then it is expected that the physics of this transition would become more clear and one could predict data at any non-trivial point on the surface which could later be checked against analytics.

Probing the finite N case: Modular invariant interpolating functions
So far in this work, we have focused mainly on the planar (large N ) limit where N = 4 SYM is expected to be integrable and there exists a well defined AdS/CFT dictionary to study various gauge-invariant observables. However, moving away from the large N , N = 4 SYM is a much harder model to study yet interesting to explore due to close relation to realistic QCD (where N =3). For finite N , it is expected that N = 4 SYM possesses S-duality [23] which connects weak and strong coupling regimes of the theory. It is a well known fact for gauge invariant observables in N = 4 SYM that one has to include all non-planar corrections to restore full S-duality to the result. In particular, the instanton contributions are very important in the S-duality context, but they are supposed to be exponentially small at large N , see e.g. [42]. The main challenge in studying S-duality invariant object is that one has to look at truly non-planar and non-perturbative data, which is hard to find in many cases.
With this new symmetry at hand, we would attempt to take a step towards constructing interpolating functions on the same line of modular invariant interpolating forms presented in [22] to approximate the anomalous dimensions of twist-two operators. As before, we pay special attention to the dependence on the spins. Below we will construct the modular invariant interpolating function first for the cusp anomalous dimension and then for the twist-two anomalous dimension with general spin dependence.
We recall our convention

Cusp anomalous dimension
As we mentioned before, the cusp anomalous dimension is related to the UV singularity of a cusped Wilson loop [8], and since a Wilson loop is expected to have a S-duality completion, it is natural to expect that there is a non-perturbative definition for the cusp anomalous dimension that also satisfies the modular invariance. We will thus construct the interpolating function by using Eisenstein series as the building block as introduced in section 2. Concretely, we will use the FPR-type S-duality invariant interpolating functions as eq. (2.19) for our construction. Let us explain the physical data we will use for the construction.
• In the weak coupling expansion, we will use the perturbative data up to four loops given in (3.3)-(3.4), including the non-planar result that appears first at the fourth loop.
• At strong coupling, under the S-duality it is expected that a Wilson loop is related to a 't Hooft loop, and in this regard one should use corresponding quantity associated to the cusped 't Hooft loop, a construction which seems to be still missing in the literature. Instead, here we assume that the cusp anomalous dimension has 'weak' mixing effect, and we will (naively) use the holographic data of cusp anomalous dimension at strong coupling in our construction.
It would be very interesting to study this quantity, for example, by considering the D1 string coupled to cusped 't Hooft loop at AdS boundary [75,76].
• Finally, another important contribution that we will take into account is the instanton contribution [42].

Resolving power issue at strong coupling
For the FPR type interpolating function (2.19), there is a problem of matching the correct power at strong coupling as explained around (2.30)-(2.31); the constraints on the param-eters of the function does not allow for fractional powers of λ. 17 Note that the change of variable methods used in section 4.2 will not work here, since the new non-linear variable as in eq. (4.8) won't meet the requirement of modular invariance. Instead, we would have to take an alternative approach.
To take into account the failure of (2.19) to encode the fractional power in the strong coupling expansion, let us consider the following generalization of the usual construction of F (s,α) m by introducing where c 1 ∼ O(N q−t−1 ) and the only difference from (2.19) is the presence of a new integer parameter t in the first term of the numerator. The planar limit of this interpolating function for large s is given by For a negative t and large λ, the leading order for the above expression is O(λ −αt ). We can now match any fractional power in the strong coupling expansion of the form O(λ c ), where c is any fractional power by solving −α t = c . and so on. So essentially, our problem with fractional powers is solved by tweaking the previous construction and yet it is consistent with the weak coupling data (3.3) and the strong coupling data (3.5). As is evident from the above discussion, we now need to fix αt = −1/2 to get the right strong coupling expansion. But, there is a caveat as the subsequent order in the expansion of the function is O(λ −1/2 ), rather than O(1). To reproduce all the correct powers of λ in (3.5) we would need to construct a interpolating function which is a linear combination of the basic interpolating functions. The total interpolating function that has the correct properties to generate strong coupling powers like {λ 1/2 , λ 0 , λ −1/2 ....} is where all the coefficients including w 1 , w 2 are fixed by respectively matching to the O(λ 1 2 ) and O(1) coefficients in the strong coupling data. Essentially, while the first function generates terms of the order {λ 1/2 , λ −1/2 ....}, the second function complements with the 17 We will shift our convention of coupling constant to the t' Hooft coupling λ, to study the finite N scenario. This will be essential in defining the holographic limit in a systematic way. orders {λ 0 , λ −1 ....} and together they can explicitly match the strong coupling power series. For a detailed account on how to fix the coefficients of the function, the reader is directed to Appendix C.
Note that, in general the linear combination of functions could involve different number of basic functions depending on which powers we want to have in the strong coupling expansion, i.e. the general schematic ansatz would be

Non-planar corrections
To take into account the non-planar piece of data (3.4), let us consider weak coupling data up to four loops and strong coupling up to two loops. As explained before, our minimal ansatz then should have a form of (5.5), i.e. we will need a linear combination of two interpolating functions. We should note that in the large N limit the coefficients for a fixed order of g or λ in the weak or strong coupling expansions has O 1 N 2 , O 1 N 4 · · · nonplanar corrections which induces a similar non-planar contribution to the solved coefficients of the interpolating function. Hence, the planar parts of the interpolating function receives no corrections from the non-planar data and the construction of the previous section goes through with systematic non-planar corrections.
However, the effect of adding such a term shows up starting from O(λ 4 ) and subsequent orders in weak coupling expansion receive non-planar corrections, so does the strong coupling expansion. 19 From the weak coupling expansion (3.3), we can see that the numerical values of the planar and non-planar terms are of the same order. As a check of the construction, we can predict the data for the 5th-loop O The γ (5) cusp,w predicted from our construction can be checked against the result from the BES equation(3.6) and our result is within 35% error bar. We can also predict the next order in the strong coupling limit and we found that our result is only 0.6% off from the result obtained from BES equation (3.6) 21 18 For example, in [22] we have studied the anomalous dimension of the Konishi constrained on the strong coupling by supergravity. However if we assume that Konishi operator doesn't level cross on the strong coupling to double trace operator and used the naive gauge theory result at the strong coupling which starts at O(λ 1/4 ), one would have to start with n = 3. 19 Here, we would like to stress the fact that strong coupling limit of the interpolating function is first the N → ∞ limit and then the λ → ∞ limit. It allows us to drop the exponentially suppressed O e − 8nπ 2 N λ terms arising from the non-perturbative parts of the Eisenstein series and focus only on the the perturbative expansion. 20 We mention that the 1/N 4 correction appears only from six loops, where higher order group theory invariants enter, see e.g. Appendix B of [77]. 21 For completeness, let us give an example of parameter counting here. Lets start with combination of two I5's as in 5.5, each of which contains five unknown coefficients. After matching the weak coupling data up to four loops, three unknown variables remain, of which two can be solved using the two loop strong  Figure 7: In the above figures the blue and green lines respectively denote the modular invariant interpolating functions for the cusp anomalous dimension (Γ cusp (λ)) with zero and non-zero value of the non-planar term Γ np . The non-planar correction decreases with increasing N and hence the maximum difference is visible at N = 2 and g = 1. The same is plotted for N = 20 zooming around g = 1, which shows a decrease in the difference with increasing N .

Adding Instanton corrections
The next hurdle in the process is to include instanton corrections at the weak coupling expansion of the function. Such a correction would first occur at O(g 4 e − 2π g ) (3.8) . Since the weak coupling expansion does not have any fractional powers of λ, we could start with a single interpolating function which gives the right powers up to four loop expansion and solve the coefficients accordingly. For example, taking the FPR-like function (see (2.19)) To incorporate all N powers and g powers at the right places it is better to follow with a construction discussed previously in this section and use a linear combination of interpolating functions (I . (5.8) The choice of t = 0 simply points to the fact that we do not constrain our function from the strong coupling data and as we discussed such functions for this choice are nothing but the FPR-like duality invariant functions as in (2.19). Here the first function simply generates the weak coupling perturbative terms, the second one generates the O(g 4 e − 2π g ) with right power of N and the third one compensates the extra terms generated in the process, for details and exact forms we refer the readers to Appendix D. It is the combined effect of coupling data. The remaining one variable should be thought of as an adjustable free parameter and we can fix it by demanding that the absolute relative error in predicting the weak coupling planar O(λ 5 ) and strong coupling O λ −1/2 terms is minimised. the last two functions that enables us to get the right instanton contribution at the weak coupling. 22

Finite spin twist-two anomalous dimensions
In this subsection we consider the twist-two anomalous dimension with finite spin dependence. There is an essential difference between large N and finite N cases. In the planar limit the single-trace twist-two operators have well-defined anomalous dimension to all range of 't Hooft coupling. On the other hand, at finite N , we have to be cautious in our program at strong coupling. For Konishi, it is clearly known that while the weak coupling anomalous dimension starts at 2 + O(λ) and grows to O(λ 1/4 ) in the strong coupling, the leading twist operator in strong coupling are the double-trace twist-four operator as (3.22), which is protected from corrections in the planar limit. This indicates that their anomalous dimensions will cross each other at certain finite coupling. The operator mixing effect would indicate that Konishi should not be self S-dual since the strong coupling description changes the operator altogether. The behaviour for higher spin operators could be more complicated. Thus for modular invariant interpolating function at finite N , one expects that there is an operator mixing behaviour such that there is change of operators from weak to strong coupling, such as from single-trace twist-two operators to double-trace twist-four operators.
The above consideration leads us to the following choice of physical data for constructing the modular invariant interpolating function.
• At weak coupling, the operator corresponds to (single-trace) twist-two operators, and the anomalous dimensions are given by the gauge theory computation as (3.14). We use the perturbative twist-two result up to three loops from (3.15). Note that the four-loop non-planar corrections have been computed for twist-two operators up to spin 8.
• At strong coupling, we will assume the dominant contribution comes from the doubletrace twist-four operators of the schematic form [tr(φ 2 )D j tr(φ 2 )](x), and we use the results from supergravity dual picture as quoted in eq. (3.23). It is indeed true that one could have other possible tower of multi trace operator eigenstates in the strong coupling side, but their contributions are expected to be suppressed in O( 1 N ). We explicitly focus on this set of double trace operator as we hope this approximation may provide some qualitative picture for the physics. 23 • The non-zero instanton correction is known for Konishi operator [43]. For operators with higher spin (j > 2) the leading instanton correction vanishes [78]. In our study of 22 Although the numerical value of the instanton correction are exponentially small in the large N limit, nevertheless it plays a crucial role in restoring the S-duality and studying AdS/CFT beyond the planar limit [42]. In principle one could also constrain the interpolating function with the non-planar and strong coupling data over this basic construction including the instanton corrections, but we find that the leading instanton sector gives very small correction to the interpolation (see Appendix D for details). 23 In general for a multi trace operator with twist 2t and spin j, one could have (t − 1) degenerate double trace operators of the form Ot n D j Ot with n = t − 2, t − 3...0. For our case, we restrict to twist-four double trace operators. Figure 8: g − j plot for the interpolating function (5.9) for N = 5 for θ = 0 and τ = i/g. One can notice here, since the supergravity data falls like O( 1 j 2 ) for a fixed N , the shape of the curve almost does not change after some initial values of j. The interpolating function has peak value at g = 1 which corresponds to one of the duality invariant point τ = τ S .
interpolating functions for operators with finite spin we will not consider the instanton corrections.
Note that in this picture, we have neglected the mixing effect between different spins, by assuming the mixing between different spin is weak. Thus we study the modular invariant function as a function of both spin and gauge coupling and expect this provides an qualitative approximation to the true physics.
From three loops in weak coupling and the supergravity result, we can see that while the weak coupling starts at O(g), there is no g dependence in the strong coupling. So following the construction we presented in section 5.1, we can again restrict ourselves to interpolating functions with t = 0, i.e. which reduces to FPR-like functions of (2.19). Note here, that in contrast to cases in the previous subsection, we have O( 1 N 2 ) corrections in the strong coupling side. Here, our discussion presented in section 2.2 on O( 1 N 2 ) correction comes into play. For this case the interpolating function does not need to match non-trivial powers on both sides and is simply given by where the coefficients c k and d k are functions of spin j and the details can be found in appendix E. Features of this function, namely the presence of extremas in τ plane, has been plotted in figures 9 and 10 for simplicity. One could actually see here that the peak values of the function on the whole τ plane appears respectively at two points τ = i and τ = e iπ/3 , which turn out to be two special values of τ in (2.13) invariant under S-transformation and (T · S)-transformation.

On level crossing
We can consider the phenomena of level crossing between the single-trace and the double- Figure 9: g−j plot for the interpolating function (5.9) for N = 5 for θ = π and τ = 1/2+i/g. The interpolating function has peak value at g = 2 √ 3 which corresponds to one of the duality invariant point τ = τ T S . trace operators of the same spin using our formalism. This in particular means that as we increase the coupling g the dimension of the leading twist operator increase while the dimension of the subleading twist operator decreases. Therefore, it is possible that for some finite value of the coupling the dimension of both the leading and subleading operators becomes equal and they cross over.
As an attempt to probe this, we consider two interpolating functions, for purely singletrace operators and purely double-trace operators respectively, and we study their crossing.
At first, consider the function which gives out weak and strong coupling data for twist-four double trace operators in the finite N case. The strong coupling data for such an operator with leading order corrections in O( 1 N ) can be simply written from the supergravity approximation, ∆ DT,sugra (j, N ) = 2 − 96 (j + 1)(j + 6)N 2 . (5.10) In the weak coupling, since not much is known for leading order correction to anomalous dimensions of double trace operators of general spin, 24 we will naively just consider the bare dimension of such operators and construct a simple interpolating function which reproduces the correct leading order results The other case we need is a function defined in finite N that mimics the behaviour of finite spin twist-two operators at both ends of the spectrum. We had already constructed such a function in section 4 using Padé approximants for large N case. As further naive approximation, one could demand that such a function could be valid also for smaller N as in large N one would not expect any crossing. We can compare the anomalous dimension of the single-trace and double-trace operators, as shown in figure 10. We observe that the anomalous dimensions of the twist-two and twist-four operators do cross-over at some finite value of the coupling g (red and the blue dashed curve in figure 10) which marks the onset of level crossing in this case. Hence, in such a physical crossing over between the dimensions, the interpolating functions are reliable only up to the crossing region as the dimensions themselves change their behaviour in the vicinity of such point.
Our modular invariant interpolating function (5.9), on the other hand, explicitly takes single-trace twist-two anomalous dimension in weak coupling to double-trace twist-four anomalous dimension in strong coupling, thus avoiding the crossing of anomalous dimensions. This non-crossing of the dimensions is consistent with the Wigner-von Neumann no-crossing rule, where the dimensions of the new eigenstate (due to operator mixing) would repel each other. The mechanism for such non-crossing is precisely due to the non-planar corrections as studied in [24].
Let us consider the spin dependence of the crossing points. At large spin the anomalous dimension the crossing is expected to happen at lower values of the coupling. In large spin the anomalous dimensions of leading twist operators grows logarithmically as γ(j) ∼ Γ cusp,w log(j) . the dimensions of the leading and the subleading twist operators will cross over. As we increase the spin of the operator the cross-over will happen at lower values of the coupling [19]. In figure 11 we have plotted the crossing point for large spin limit and fitted with approximate fitting function We observe that as we increase the spin of the operator the crossing happens at lower values of coupling and the coefficient of the approximate fitting function is close to the 2π N . For small spin (j ∼ coupling) the cross over between the anomalous dimension of the single trace operator (3.18) and the double trace operator in the strong coupling(supergravity limit) occurs at the strong coupling. As we increase the spin of the operator the cross over shifts towards the weak coupling end.
In figure 12 we have presented plot of the crossing point for the anomalous dimension of the leading and the subleading twist operators as a function of j for finite but smaller spin. From figure 12 we observe that as we increase the spin of the operators the crossing happens for lower value of the coupling. Hence, for very large spin of the operator the crossing occurs in the weak coupling and there is maximal mixing of the operator in the large spin limit.

Upper bounds on anomalous dimensions and conformal bootstrap
One immediate application of the finite N interpolating function could in principle be in constraining the upper bound on the anomalous dimensions of leading twist operators for any arbitrary spin j. The N = 4 superconformal bootstrap approach [81] has obtained some upper bounds on the dimensions of the unprotected leading twist operators by studying the four-point function (see also [82,83,84]) where O I 20 is a superconformal primary scalar operator of dimension two in energymomentum tensor multiplets transforming as 20 representation in SU (4) R . The N = 4 superconformal symmetry allows us to describe the four-point function in terms of the N = 4 superconformal block [65,85,86,87]. It is conjectured in [81] that the bound on the anomalous dimension γ(j) has a global maximum at one (or both) of the duality invariant points τ = τ S = i or τ = τ T S = e iπ/3 .
To compare our interpolating functions with results from N = 4 superconformal bootstrap, we have to investigate where the interpolating function takes its maximal value as a function of τ . We can actually expect that the extremal value of the interpolating function is given at either of the duality invariant points in the τ direction since the building blocks of the interpolating functions have extremas at precisely these duality invariant fixed points. 25 From figure 8 we observe that the interpolating function has a peak at g = 1, which is one of the duality invariant points τ s .
We also show in figure 9 that for θ = π the interpolating function has maximum at g = 2 √ 3 , which corresponds to the other duality invariant point τ = τ T S . To predict an upper bound on the anomalous dimension it is important to find the global maxima of the interpolating function at either of this two duality invariant points. In figure 13 we have plotted the interpolating function at these two duality invariant fixed points for some arbitrary values of j at a fixed N = 2 and observe that the value of the interpolating function at τ = τ T S always stays larger than the value at τ = τ S . It seems to be a generic feature for other values of N .
With the observations that our interpolating function for the anomalous dimensions of leading-twist operators takes a maximum value at one of the duality invariant point τ = τ T S , we could give a conservative prediction for the maximal value saturated by the anomalous dimension arising from the conformal bootstrap. In figure 13 we present this maximal value of the anomalous dimension at τ = τ T S (orange) and τ = τ S (blue) 25 Global maxima of Es(τ ) is given by τ = τT S . for arbitrary j (up to j=100). It should be noted that we just present a very crude approximation here based only on the maximal value of our interpolating function, which can be compared with rigorous results from superconformal bootstrap data when they become available.

Discussion
In this paper, we focus on one of the simplest classes of observables in N = 4 SYM, namely, the anomalous dimension of twist-two operators. We study their non-perturbative completion via interpolating function method, paying special attention to the dependence on the general spin parameter. These observables in large N theory have been extensively studied in various aspects and is relatively well understood thanks to the AdS/CFT correspondence and integrability. Based on the results both in the weak coupling and strong coupling regime, we construct simple interpolating functions with generic spin and coupling dependence. The interpolating function as a function of spin allows us to encode the cusp anomalous dimension as the large spin limit. Here a particular interesting aspect is the intriguing 'transition' between small and large spin at strong coupling which is discussed in some details in section 4.3. When spin is small, the transition can be described within supergravity where massive string excitations effectively decouple. In the large spin limit, the stringy effects play an important role. This shows that there is non-trivial rich physics even in the planar limit.
In the case of finite N , which is more close to the realistic QCD, the physical picture becomes significantly more complicated, mainly due to the operator mixing effect. The study in this case is more on the qualitative side. We apply the constraint from modular invariance, for which we use Eisenstein series as building blocks in our interpolating Weak coupling Strong coupling Figure 14: From the flux tube picture, at strong coupling, the energy is large enough to generate pairs of charged particles, so that a single trace operator at weak coupling splits into multiple smaller single traces.
function. We first consider the cusp anomalous dimension and construct the interpolating function by taking into account the four-loop non-planar result and the instanton contribution for the first time. Here we also solve a few technical challenges of the construction, such as correctly reproducing the strong coupling expansion and encoding the instanton contribution. We also provide a prediction to five-loop non-planar result in (5.6) based on the interpolating function.
We then focus on twist-two anomalous dimension with a finite spin parameter. Since both the data at finite N is very limited and the operator mixing may in principle be rather complicated, we have to make some assumptions to simplify the picture. Concretely, we consider a spin dependent modular invariant function such that at weak coupling it is given by the single-trace twist-two operators tr(φD J φ), while at strong coupling, it is dominated by the double trace twist-four operators tr(φ 2 )D J tr(φ 2 ). We make this approximation by assuming that the mixing with other operators such as higher twists (traces) are sub dominant.
Let us comment on the relation between cusp anomalous dimension and finite spin operators. In planar limit, the large spin scaling behaviour as Γ cusp log j is true in both weak and strong coupling regime. However, in the finite N theory at strong coupling, one may expect such log j scaling is broken. In particular, in our modular invariant function the strong coupling expansion corresponds to double trace operators, which indeed have no log(j) scaling in the large spin limit. Physically, according to the flux tube picture discussed in [88], the energy at strong coupling is large enough to generate pairs of colorcharged particles, so that a single trace operator (at weak coupling) splits into multi-traces (at weak coupling), as illustrated in figure 14.
As we have a connection to QCD in mind, we might hope the discussion in N = 4 SYM provide certain qualitative picture of the physics of generic gauge theories. Indeed, one may make an analogy with QCD: in the perturbative UV regime, fundamental degrees of freedom are partons (gluons or quarks) carrying color charges and a color singlet is given by a single trace operator; while in the IR non-perturbative regime, the fundamental degrees of freedom are themselves color singlets (i.e. single traces), such as pions or baryons in the effective Chiral Lagrangian theory for QCD. This is consistent with a single trace to multi-trace transition from weak to strong coupling that we discuss in N = 4 SYM.
We have also analysed the phenomena of level-crossing between the leading twisttwo and the sub-leading twist-four operators. A rigorous one loop computation of the anomalous dimension of finite spin subleading twist operators might give us more insight into the phenomena of level-crossing. Finally, we mention that similar studies as of this paper may be applied to more general observables, such as operators with θ-angle dependence or the OPE coefficients, as well as in other theories such as ABJM. We hope to address these questions in the future.

Acknowledgments
The work of AB, ST and GY is supported by the Chinese Academy of Sciences (CAS) Hundred- Talent  A. Construction of G 8/7 In this appendix, we give the explicit construction of the finite spin interpolating Padé approximant as discussed in section (4.1). The function in question has a form The coefficients of the approximant are then given by,

B. Construction of FPP's for finite spin
We constructed FPP's to improve upon our discussion of Padé approximants in large N case in section 4.2. Here we will give the explicit FPP functions for finite spin case, which read

D. Adding instanton corrections to cusp anomalous dimension
Here we discuss systematic inclusion of the instanton corrections to cusp anomalous dimension. Note that the Eisenstein (2.17) series used in the construction of the interpolating functions has in-built within itself an infinite series of non-perturbative terms taking the form e − 2πn g (where n ∈ Z) from the weak coupling expansion of the modified Bessel (K) terms. From section 3.1, eq. (3.9) we see that we have to include such a correction at O g 4 e − 2π g .
To begin with, from the construction of interpolating function mentioned in section (2.2) and the weak coupling expansion (3.3), we notice that if an interpolating function with parameters (α, s) generates first weak coupling term at g β , the we get the first instanton term having a form O g min(p,q)+β+s e − 2π g . As min(p, q) is at least 1, the choice for s is very limited and thereby we lose the possible infinite class of modular invariant interpolating functions parametrized by s. Even if we construct a valid interpolating function say with α = 1 3 , m = 4 and s = 2 where p = 1 and q = 4, the structure of the interpolating function is as (c 1 X 1 g + · · · ) + c 2 X 1 g 4 + · · · e − 2π g + (· · · ) e − 4π g + · · · , (D.1) where c 1 & c 2 are some numerical constants and X 1 is an unsolved coefficient always shared by O(g) and O g 4 e − 2π g terms. Note that from weak coupling data X 1 ∼ N but from the instanton data X 1 ∼ N − 3 2 and hence it is impossible to consolidate both trends in a single interpolating function. As an alternative, we could have two interpolating functions I 1 & I 2 with I 1 carrying a large s such that the power of g multiplying O(e − 2π g ) is large and I 2 a minimal interpolating function used just to reproduce the O g 4 e − 2π g term.
To construct I 2 we only use the O g 4 e − 2π g data with α = 1, m = 2 & s = 2. The structure of the interpolating function turns out as follows, where c s are some constants. Firstly, note that we have generated a set of extra weak coupling terms, but they are largely suppressed at reasonably large N compared to the weak coupling data (3.3). Secondly, at O e . Though we have no data to match at this order but if we had we could easily incorporate it either by creating another interpolating function with s = 3 starting at O g 5 N 1 2 e − 2π g or by allowing more unsolved coefficients in I 2 . Now, to remove the extra weak coupling terms from eq. (D.2) up to some order in g we could consider it as input weak coupling data for yet another interpolating function (I 3 ) with a large value of s and finally subtract it from I 2 .
As an concrete example of construction and comparison with previous constructions, lets take the extra weak coupling terms from I 2 up to O g 28 . As instanton is non-perturbative, we could expect that including instanton data would be insignificant and this is clearly reflected in Table 1 where the differences between the interpolating function with instanton correction (I 1 +(I 2 −I 3 )) and without the correction (I 1 ) are analysed for the critical points of the difference function at N = 2, N = 20 and N = 200. It seems instanton correction is far less significant than the strong coupling corrections or the non-planar corrections to the modular invariant interpolating functions.
The planar interpolating functions (for Γ np = 0) and only the weak coupling data have been used in the construction of the instanton corrected modular invariant interpolating  Table 1: The difference between including the instanton correction and not including it at some critical values of g. We would like to highlight that even at the worst case of N = 2 & g = 1 the difference is O 10 −6 .
functions. The explicit forms of the functions are as follows The Eisenstein series (E s+k ) has an infinite sum over modified Bessel (K) which has been truncated to include the first twenty sums.