Analytical approaches to the determination of spin-dependent parton distribution functions at NNLO approximation

In this paper, we present {\tt SMKA18} analysis which is a first attempt to extract the set of next-to-next-leading-order (NNLO) spin-dependent parton distribution functions (spin-dependent PDFs) and their uncertainties determined through the Laplace transform technique and Jacobi polynomial approach. Using the Laplace transformations, we present an analytical solution for the spin-dependent Dokshitzer-Gribov-Lipatov-Altarelli-Parisi evolution equations at NNLO approximation. The results are extracted using a wide range of proton $g_1^{p}(x, Q^2)$, neutron $g_1^{n}(x, Q^2)$ and deuteron $g_1^{d}(x, Q^2)$ spin-dependent structure functions dataset including the most recent high-precision measurements from {\tt COMPASS16} experiments at CERN which are playing an increasingly important role in global spin-dependent fits. The careful estimations of uncertainties have been done using the standard 'Hessian error' propagation. We will compare our results with the available spin-dependent inclusive deep inelastic scattering dataset and other results for the spin-dependent PDFs in literature. The results obtained for the spin-dependent PDFs as well as spin-dependent structure functions are clearly explained both in the small and large values of $x$.

From the practical point of view, experiments including nucleon beams at the current and future energy frontiers need most accurate information on the spindependent PDFs to exploit their data. In the absence of improved inclusive deep inelastic scattering (DIS) data, most attention is now turned towards data from the Large Hadron Collider (LHC) experiments. In addition, the nucleon spin structure has been always a fundamental question in high energy physics so it has been extensively studied both in theory and experiment in recent decades. From an experimental point of view, several experiments have been set up to study the longitudinal spin structure of the nucleon, such as COMPASS experiments at CERN [23][24][25], HERMES experiments at DESY [26][27][28], many experiments at Jefferson Laboratory (JLAB) [29][30][31], and PHENIX and STAR experiments at the protonproton (pp) Relativistic Heavy Ion Collider (RHIC) [32][33][34][35]. The majority of this experimental information on * M.Salajegheh@stu.yazd.ac.ir † mmoosavi@yazd.ac.ir ‡ Hamzeh.Khanpour@mail.ipm.ir § Atashbart@gmail.com spin-dependent PDFs come from the neutral-current inclusive and semi-inclusive deep inelastic scattering, DIS and SIDIS, with charged lepton beams ℓ ± and nuclear targets. These inclusive spin-dependent DIS data constrain only the total quark combinations, while the SIDIS data constrain individual quark and anti-quark flavors in the nucleon. In principle, both DIS and SIDIS data are also sensitive to the gluon distribution, however, the constraining power of DIS and SIDIS data on the gluon distribution is rather weak due to the limited kinematical range covered by these data sets.
In our previous study, KTA-I-17 [20], we performed a first analysis of spin-dependent inclusive DIS by taking into account the contributions from target mass corrections (TMCs) and higher twist terms (HT) to the spin-dependent g 1 and g 2 structure functions. In KTA-II-17 [21] we considered the effects of nuclear corrections such as Fermi motions, spin depolarizations, binding and the presence of a non-nucleonic degree of freedom, shadowing and anti-shadowing corrections which are necessary at these kinematics. It is worth mentioning here that, many spin-dependent PDFs analyses impose more stringent cuts on the photon virtuality Q 2 as well as W 2 in order to avoid dealing with the complications associated with the HTs and nuclear corrections. Unfortunately, these restrictions eliminate much of data at the highest x-values. Most of these phenomenological spin-dependent PDFs analyses also utilize standard PDFs fitting technology in which single fits are performed by assuming the basic functional forms for the PDFs. In this approach, the fit parameters are obtained by minimizing the overall χ 2 global . Then, the PDFs errors are typically computed using the Hessian method or Lagrange multiplier or neural network. In the present paper, we construct for the first time a set of spin-dependent PDFs at NNLO approximation using a methodology, the so-called "Laplace transform technique" and "Jacobi polynomials approach", which has been recently used to study the polarized [19] and unpolarized PDFs [41,42].
As in KTA-I-17 [20] and KTA-II-17 [21] analyses, we use the spin-dependent inclusive DIS data whenever available. We include new data sets with high-precision COMPASS [24,25] measurements at CERN on the proton and deuteron. To isolate the impact of recent COMPASS data sets more directly, and to assess the systematics of our new methodology based on Laplace transform technique and Jacobi polynomials approach, we restrict the current analysis to the inclusive DIS data. A full global QCD analysis of all available data including SIDIS, inclusive jet and weak boson production in polarized protonproton collisions, will be presented in our forthcoming studies.
The structure of our paper is as follows: An analytical solution for the NNLO spin-dependent Dokshitzer-Gribov-Lipatov-Altarelli-Parisi(DGLAP)evolution equations is presented in Sec. II. This section also includes the nonsinglet, singlet and gluon solutions in the Laplace space at the NNLO approximation. The theoretical calculations for the leading-twist spin-dependent DIS structure function and Jacobi polynomials approach used in the SMKA18 spin-dependent PDFs fit are summarized in Section III. The dataset used to determine the SMKA18 spin-dependent PDFs are briefly summarized in Sec. IV. Section V, includes the choice of input scale and the SMKA18 parametrization basis. The results of the present spin-dependent PDFs analysis are given in Sec. VI. This section includes a detailed comparison between the present results and available spin-dependent inclusive DIS data. We also present a detailed comparison of our NNLO results with recent results in literature. Finally Sec. VII includes our summary and conclusions. In Appendix A, we present the analytical solutions for the NNLO splitting functions in Laplace space. Appendix B includes the coefficient functions of singlet distributions in the Laplace space.

II. ANALYTICAL SOLUTION FOR THE HELICITY-DEPENDENT DGLAP EVOLUTION EQUATIONS
We perform our analysis in Laplace s-space which has the advantage of significantly shorter fitting time compared to the x-or Mellin-space based analyses. Following this, we describe the novel aspect of our analysis, namely the Laplace technique. It turns out, however, that in this method the computational time can be significantly reduced through the use of Laplace s-space techniques. This is due to that, firstly, the Q 2 evolution equations in Laplace space are the ordinary coupled differential equations which are faster and also simpler to be solved in comparison to the corresponding integro-differential equations in the x-space. Secondly, using this technique it is possible to cast the various multidimensional integrations in terms of precomputed quantities, which can significantly decrease the computational time needed for the observables in the global QCD fits.
Here, we follow the Laplace technique described in our previous unpolarized PDFs analyses [41,42] to present an analytical solution for the coupled NNLO DGLAP evolution equations. This method has been successfully developed and used in variety of QCD analyses, (see Refs.  for clear reviews). We will show that this method can be also used in the spin-dependent case and, hence, one can extract the spin-dependent PDFs inside the nucleon from QCD analysis of spin-dependent inclusive DIS data. Here, we will give the details of our approach and review the method of extracting the spin-dependent PDFs. It must be noted that, first, we will focus on the nonsinglet solution at NNLO approximation in Laplace space and, second, we will present our analytical solution for the singlet and gluon cases.

A.
Nonsinglet solution in Laplace space at the NNLO approximation According to the decoupling of the DGLAP evolution equations in the Laplace method, used in SMKA18 QCD analysis, we convert the solutions into three parts which take the nonsinglet ∆f NS , singlet quark ∆f S and gluon ∆f g distributions. In this section, we present our solution for the nonsinglet sector and left the singlet quark and gluon ones for the next section.
As was mentioned in the I, the DGLAP evolution equations [68][69][70][71] are a set of integro-differential equations which can evolve the PDFs into a desired Q 2 -scale. The nonsinglet sector of the DGLAP evolution equations at NNLO approximation reads: where ∆q NS stands for the nonsinglet spin-dependent PDFs, the symbol ⊗ denotes the convolution and α s (Q 2 ) is the QCD coupling constant. The parameters ] are the nonsinglet spin-dependent Altarelli-Parisi splitting kernels at one, two and three loops corrections, respectively.
Let us now briefly review the method of extracting the spin-dependent PDFs using the analytical solutions of DGLAP evolution equations applying the Laplace transformation technique. Considering the variable definitions as ν ≡ ln(1/x) and w ≡ ln(1/z), one can rewrite the evolution equations (1) in terms of the convolution integrals and also the new variables ν and w. Therefore, one can obtain a simple solution as follows: It is worth noting here that the Q 2 -dependence of the above evolution equations is expressed thorough the variable τ as τ (Q 2 , By defining the Laplace transforms as ∆f NS (s, τ ) ≡ L[∆F NS (ν, τ ); s] and considering the fact that the Laplace transform of convolution factors is simply the ordinary product of the Laplace transform of the factors [66,67], the Laplace transform of Eq. (2) leads to the ordinary first order differential equations in Laplace space s with respect to the τ -variable. Therefore, by working in the Laplace s-space, one can obtain the first order differential equations for the nonsinglet distributions ∆f NS (s, τ ), as A very simplified solution for the above equation is where, up to the NNLO approximation, the Φ NS (s) contains contributions of the splitting functions at the sspace. Up to NNLO approximation, it reads (5) These splitting functions can be calculated from xspace results and presented in Refs. [72,73]. The Q 2dependence variables τ 2 and τ 3 in Eq. (5) are defined as follows and The LO and NLO contributions to the splitting function ∆Φ LO NS (s) and ∆Φ NLO NS (s) are presented in Ref. [43] and the NNLO one, i.e. ∆Φ NNLO NS (s), calculated in this analysis and presented in Appendix A.

B. Singlet and gluon solutions at the NNLO approximation
In this section, we turn to present our solutions for the singlet quark and gluon evolutions in Laplace s-space. The coupled NNLO DGLAP equations, using the convolution symbol ⊗, for the singlet quark ∆F s can be schematically written as The corresponding evolution for the gluon density ∆G at NNLO approximation is given by In the two above equations, ∆P 0 qq , ∆P 0 qg , ∆P 0 gq and ∆P 0 gg are the LO singlet splitting functions. Correspondingly, ∆P 1 qq , ∆P 1 qg , ∆P 1 gq and ∆P 1 gg are the NLO splitting functions, and ∆P 2 qq , ∆P 2 qg , ∆P 2 gq and ∆P 2 gg are the NNLO ones. The α s (Q 2 ) is the NNLO strong coupling constant.
As for the nonsinglet sector presented in Sec. II A, considering new variables τ ≡ 1 and z = e −w , the DGLAP evolutions in Eqs. (8) and (9) can be written in terms of the variables τ , v, and w. Following the Laplace transformation technique, the convolution integrals in DGLAP evolutions can be converted from v space to s space [43,66]. In this respect, one can achieve the following equations for the quark singlet ∆f (s, τ ): and for the gluon ∆g(s, τ ), one has which are expressed in terms of new variable τ .
By introducing the following simplifying notations: one can finally rewrite the singlet distribution (12) as and the gluon density (13) as The solution of Eqs. (17) and (18) can be obtained via an iteration process. We first consider the simple solutions of these equations which are labeled as ∆F 1 (s, U ) and ∆G 1 (s, U ). These are concluded by setting a 1 = 0 and can be given by Their solutions lead to the following equations The denominator D(U, s) is the determinant of the coefficients of ∆F (s, U ) and ∆G(s, U ) in Eqs. (12) and (13), i.e..
We now construct an iterative solution for the Eqs. (17) and (18) to obtain the ∆F and ∆G. To achieve the first iteration solution, we need to change the arguments of ∆F and ∆G in Eq. (20). They can be presented as ∆F (s, U + b 1 ), ∆G(s, U + b 1 ), ∆F (s, U + 2b 1 ) and ∆G(s, U + 2b 1 ), respectively. The numerical values for the unknown parameters in Eq. (16) at this order of iteration are extracted through a fitting procedure so we determined them as: a 1 = 0.0591, b 1 = 7.0038 and a 0 = 0.00498. Now, to complete the calculations and to achieve the second step of the iteration process one can return from U space to τ space using the inverse Laplace transform. This yields the following expressions for the singlet and gluon distributions, The analytical expressions for the k(s, τ ) coefficient functions, up to two step of iteration, are given in Appendix B. To obtain the spin-dependent PDFs as well as the structure functions in x and Q 2 space, we used the numerical Laplace transform algorithms presented in Ref. [45] for the numerical inversion of Laplace transformations and convolutions.

III. LEADING-TWIST SPIN-DEPENDENT DIS STRUCTURE FUNCTION AND JACOBI POLYNOMIALS APPROACH
In the following, in detailed discussions we will describe the basic theoretical issues for a consistent determination of the spin-dependent PDFs from spin-dependent inclusive DIS data. We work in the framework of perturbative QCD at NNLO approximations using the MS scheme for the renormalization and factorization. In the leading twist (twist τ = 2) approximation, the spin-dependent proton structure function g p 1 (s, Q 2 ) in Laplace s-space at NNLO can be expressed as a linear combination of the spin-dependent PDFs and hard-scattering Wilson coeffi-cient functions [18-20, 22, 74], as In the above equation, ∆q v (s, Q 2 ), ∆q s (s, Q 2 ) and ∆g(s, Q 2 ) are the Laplace transforms of spin-dependent valance, sea and gluon densities, respectively. The ∆C (1) q and ∆C (1) g are the NLO spin-dependent quark and gluon hard-scattering Wilson coefficients functions, calculable using the Laplace transform [43]. At NNLO approximation the Wilson coefficients are different for the quarks and antiquarks cases and we use ∆C (2) ns and ∆C (2) s calculated in Ref. [75]. The spin-dependent inclusive DIS is the only polarized process for which the hard-scattering coefficient functions are known up to NNLO.
In Ref. [20], we have shown that the power-suppressed O(1/Q 2 ) corrections to the structure functions as well as the target mass corrections (TMCs) can make important contributions in some kinematic regions. In addition to the TMCs, spin-dependent structure functions in the operator product expansion (OPE) also receive contributions from higher twist terms (HTs) which are associated with matrix elements of multi-quark or quark-gluon operators. In Ref. [20], we have also considered the contributions from the non-perturbative HTs. Since our main aim in this analysis is to study the applicability of Laplace transform to the QCD analysis of spin-dependent structure function, we restrict ourselves to the twist τ = 2 approximation in Eq. (23).
The method employed in the SMKA18 QCD analysis is based on the Jacobi polynomials expansion of the spindependent structure functions. Practical aspects of this method including its major advantages are presented in our previous studies [19,20,41,42] as well as other literature [76][77][78][79][80][81][82][83][84][85][86][87][88]. Since the basic idea of Jacobi polynomials can be found in the mentioned references, we outline a brief review of this method. To illustrate this technique, we consider the case of xg 1 in Eq. (23). The spindependent structure function xg 1 (x, Q 2 ) is given by where, n is the order of expansion terms, N max is the number of polynomials which normally can be set to 7 or 9 and, a n (Q 2 ) are the Jacobi moments. The α and β parameters are fixed as α = 3 and β = 0.5, respectively. In our previous analyses, it is shown that with these selection one can achieve the fastest convergence of the above series [19,20,41,42]. One can conclude form Eq. (24) that the use of Jacobi polynomials for the expansion of the structure functions has an advantage to allow one to factor out the essential part of the x-dependence of the structure functions into a weight function x β (1 − x) α . Also, the Q 2 -dependence is contained in the Jacobi moments.
The Θ α,β n (x) in Eq. (24) are the Jacobi polynomials with the following expansion, where c (n) j (α, β) are the coefficients which are expressed through γ-functions. It is worth mentioning here that the Jacobi polynomials satisfy the following orthogonality relation with the weight function Using the above orthogonality relation, one can relate the spin-dependent proton structure functions with their Laplace s-space moments as follows On the right-hand side of the above equation, the L[xg 1 , s = j + 1] are the Laplace transformation of the structure functions.

IV. OVERVIEW OF SPIN-DEPENDENT INCLUSIVE DIS DATASETS
In this section, we summarize the polarized DIS data sets used in SMKA18 QCD fits. We first review the analyzed data sets used in our work, then we will discuss on the present and future experimental efforts on hardscattering polarized observables which can provide additional constrain on the gluon density and can also be used to separate the polarized quark and antiquark distributions.
The core of all spin-dependent PDFs fits include the spin-dependent DIS data obtained from neutralcurrent inclusive and semi-inclusive deep-inelastic scattering, DIS and SIDIS, with charged lepton beams and nuclear targets at various fixed-target experiments. The set of spin-dependent DIS data are provided in terms of proton, neutron and deuteron spin-dependent structure functions by different experiments at CERN, HERA and JLAB. The polarized fixed-target experiment at HER-MES (at HERA collider) have collected large amounts of data for the proton, neutron and deuteron as well as the heavier targets. The SMKA18 global PDFs analysis uses all available data sets on the inclusive DIS of leptons over proton, neutron and deuteron that pass the required cuts on the invariant final state mass W 2 ≥ 4 GeV 2 and Q 2 ≥ 1 GeV 2 . This includes all of data sets from the EMC, SMC, COMPASS, SLAC, HERMES, and Jefferson Lab Hall A experiments which we also used in our previous global fits [19][20][21]. Our analyzed data sets also include the most recent high-precision data from COMPASS16 experiment at CERN [24,25].  Figure 1: (Color online) Representative kinematic coverage, in the (Q 2 , x) plane, of the most recent neutral current polarized inclusive DIS measurements on proton and deuteron targets reported by the COMPASS16 experiment at CERN [24,25]. These are used as input in SMKA18 global spin-dependent PDFs fit.
A list of published spin-dependent DIS experimental data points used in SMKA18 analysis are described in Table I. For each experiment we provided the x and Q 2 kinematical ranges, the number of data points for each given target as well as the fitted normalization shifts N i obtained from the QCD fit to the data. As is seen from Table I, different combinations of spin-dependent inclusive DIS data sets for the proton, neutron and deuteron obtained by various experiments at CERN, HERA and JLAB are used in SMKA18 analysis. Our analysis does also contain the most recent data by COMPASS16 [24,25] experiment which are not used by other groups. In Refs. [19] and [20], we have shown that the inclusion of the recent COMPASS16 data leads to a reduction in the spin-dependent PDFs uncertainties for the valence and sea quarks as well as in the gluon PDFs uncertainty at the small value of x.
The kinematic coverage of recent COMPASS16 data for the proton and deuteron targets [24,25] is displayed in Fig. 1. These are used in our analysis to set an additional constrain on the spin-dependent PDFs. From this figure, one can conclude that the quantity of high precision data points from COMPASS16 experiments and their kinematic coverage are presently much more limited in the polarized case. Therefore, the polarized PDFs can currently be determined with much less precision than the unpolarized PDFs and only cover the x-range for x ≥ 0.0035. One can expect that the future LHC data will certainly represent important opportunities to improve our knowledge on the spin-dependent PDFs. However, there are many challenges for the spin-dependent PDFs fitters to include such increasingly precise datasets effectively within a spin-dependent PDFs fit. A summary of these data sets can be found in Refs. [2,[19][20][21]105].
In addition to the DIS and SIDIS fixed-target data sets mentioned above, a significant amount of spin-dependent data from longitudinally polarized proton-proton collisions at RHIC has become available recently, which is however, in a limited range of momentum fractions, 0.05 < x < 0.4 [102]. The longitudinal single-spin and double-spin asymmetries for the weak boson productions are sensitive to the flavor decomposition. These can also be used to separate the polarized quark and antiquark distributions [103].
The STAR data sets on the W -boson production at polarized proton-proton collisions at RHIC provide evidence of a positive ∆ū distribution and a negative ∆d distribution [123]. The double-spin asymmetries for the production of jet, di-jet and π 0 are also sensitive to the gluon polarization [104], directly. The kinematic coverage of the spin-dependent data, the quantity of the data points, and the variety of available hard-scattering processes are presently much more limited for the polarized case in comparison with the unpolarized one [105]. Hence, the spin-dependence PDFs can currently be determined with much less precision than the unpolarized PDFs, specially at the small range of x. The kinematic coverage is expected to be significantly improved in the future, with the DIS and SIDIS data from 12 GeV upgrade of Jefferson Lab [106] and future spin-dependent high-energy and high-precision Electron-Ion Collider (EIC) [107]. We should mention here that the QCD analyses of spindependent nuclear-target data requires an accurate account of nuclear corrections. For recent reviews see, Table I: Summary of published spin-dependent DIS experimental data points above Q 2 = 1.0 GeV 2 used in SMKA18 global analysis. For each experiment, it is listed the x and Q 2 ranges, the number of data points for each given target, and the fitted normalization shifts Ni (see the text). for example, Refs. [21,29,[108][109][110][111][112]. Briefly, we have demonstrated that in the spin-dependent case, even with current uncertainties, new and precise data from JLAB, RHIC, EIC and CERN can impose sizable constraints on several important spin-dependent quark combinations. This suggests, the global spin-dependent PDFs analyses should include future data sets in their fits to constrain the gluon density much more. These new data can also constrain some of the less well-known quark combinations, such as the total strangeness.

V. SPIN-DEPENDENT PDFS PARAMETRIZATIONS AND ERRORS
We consider a proton comprised of massless partons with spin-dependent distributions q ± (x, Q 2 ) which carry the momentum fraction of x with a characteristic scale Q. In the present analysis, for the generic parametrization of the spin-dependent PDFs, we take into account the following standard functional form at the initial scale which includes four shape parameters η q , a q , b q , c q , and the normalization coefficients N q . The generic labels q = u v , d v , q, g refer to the partonic flavors of up-valence, down-valence, sea, and gluon, respectively. The normalization constants N q are determined as where the function B(a, b) is the Euler β function.
These are chosen such that η q are the first moments of x∆q(x, Q 2 0 ). Since the present analysis just considers the spin-dependent inclusive DIS data, then we attempt to assume an SU(3) flavor symmetry such that ∆q(x, Therefore, we try to fit only the spin-dependent PDFs of x∆u v , x∆d v , x∆q and x∆g.
One can also consider additional constraints on the moments of the spin-dependent PDFs which could provide by the weak neutron and hyperon decay constants. Hence, the first moments of the spin-dependent valence distribution can be described in terms of the axial charges for octet baryons, F and D, in which measured in hyperon and neutron β decay [113]. These constraints lead to the values of η uv = 0.928 ± 0.014 and η dv = −0.342 ± 0.018 [114]. We fix two valence first moments on their central values. The parameters η q and η g are determined from the QCD fit to inclusive data.
For the spin-dependent quark and gluon distributions in Laplace s-space, we follow our previous PDFs analy-ses [41,42]. As we discussed in Sec. III, in the Jacobi polynomial approach, the DGLAP evolution equations can also be solved in the Laplace s-space. The Laplace transformation of the spin-dependent PDFs ∆q are defined as in Eq. (28), i.e.
The fit parameters are determined by minimizing the χ 2 global ({η k }) function, which are defined as where {η k } is a set of fitted parameters, D i j is the measured value of the observable for the data point j from the experimental data set i. In Eq. (31), T i j ({η k }) is the corresponding theoretical value and σ i j represents the uncorrelated statistical and systematic uncertainties added in quadrature. Another important issue that is addressed in SMKA18 global analysis is the estimation of uncertainties in the extraction of various spin-dependent PDFs, associated with experimental uncertainties. SMKA18 pursues here an approach based on the use of the Hessian method [115]. In the next section, we will show that in the standard single-fit PDF analyses, one often finds that some of the shape parameters in Eq. (28) are not well determined by spin-dependent inclusive data and need to be fixed by hand.

VI. RESULTS AND DISCUSSIONS
In this section, we will present our detailed discussions on the obtained results and compare our extracted results with the others in literature. Having our methodology for fitting the spin-dependent PDFs, which includes the spin-dependent structure functions in Laplace space, the analyzed observables and the SMKA18 spin-dependent input, we are now in a position to present our main results of SMKA18 analysis in this Section.
Before starting our detailed discussions on SMKA18 spin-dependent PDF determinations, we would like to illustrate the present status of spin-dependent PDFs sets at NLO and NNLO QCD approximations and discuss some differences, briefly, which are clearly visible. The available PDF sets at NLO and NNLO approximations are listed in Table. II which include AAC09 [37], BB10 [38], LSS10 [39], NNPDF14 [15], JAM13 [16], JAM14 [17], JAM15 [18] at NLO, and TKAA16 [19], KTA-I-17 [20] and KTA-II-17 [21] at NNLO approximations. In the first column of Table. II, we indicate the name of the group and in the subsequent columns we present the theory accuracy, the data sets and the corresponding references. These spin-dependent PDFs determinations are different in the input parametrizations, the datasets included in the analyses, the details of the QCD analysis such as the theory accuracy and the treatment of HT and TMCs, and finally in the procedure used to determine the spin-dependent PDFs from the data. In the theory setup point of view, the DSSV has developed a method based on "Mellin moments" of the PDFs and JAM collaboration has implemented a new approach called iterative "Monte Carlo" procedure [18,116]. The DSSV14 [40] and NNPDFpol1.1 [15] spin-dependent PDFs update their previous analyses, DSSV08 [36,117] and NNPDFpol1.0 [118] have included almost all available experimental information. The impact of RHIC proton-proton data have been studied by the inclusion of data from double-spin asymmetries for inclusive jet production and π 0 production [119,120] and the data on double-spin asymmetries for high-p T inclusive jet production [119,121,122] and weak boson production [123].
The data from COMPASS16 experiments at LHC [24,25] certainly represent important opportunities to improve the knowledge on spin-dependent PDFs. There are, however, many challenges for the spin-dependent PDFs fitters to include such increasingly precise data effectively within a spin-dependent PDFs fit. These data can be used to further constrain the spin-dependent PDFs and to measure the strong coupling constant α s . These experiments cover a wide range of physical observables. Additional experimental information is expected from ongoing and future experimental efforts. COMPASS16 data sets on polarized proton and deuteron targets have been used in our previous NNLO QCD analyses. As we mentioned, in our present analysis we also use these precise data sets. Hence, our main aim of this paper is to introduce a new method to extract the spin-dependent PDFs; the Laplace transform and Jacobi polynomials approach.
After this brief review on recent efforts on spindependent PDFs analyses, we now present the results obtained for the SMKA18 spin-dependent PDFs at NNLO in the basis of the Laplace transform technique and Jacobi polynomial approach. The SMKA18 fit parameters are presented in Table. III. The parameters c q and c g are set to zero because the current spin-dependent inclusive DIS data can not constrain all the fit parameters. As we have mentioned earlier, despite of the outstanding achievements, at the present the spin-dependent PDFs cannot be well determined in a global QCD analysis with a high accuracy as one has for the unpolarized ones. The experimental data on spin-dependent DIS are relatively span narrow range of x and Q 2 . As a consequence, the quarks, antiquarks, and gluon densities are still affected by large uncertainties. Now, we turn to discuss in more detail , how the results on SMKA18 spin-dependent PDFs depend on the method used for their determination. To illustrate this   dependence, we will compare SMKA18 NNLO set of spindependent PDFs determined by Laplace method with those obtained by KTA-I-17 [20] and TKAA16 [19] using the Mellin space analyses. The differences between the mentioned NNLO spin-dependent PDFs lie in the formalism used in these analyses, Mellin and Laplace space, as well as the treatment of higher twist and TMCs.
The final distributions for the SMKA18 fit are displayed in Fig. 2 as a function of x at fixed Q 2 0 = 1 GeV 2 . For comparison, we have also shown the Mellin space results from recent NNLO analyses of KTA-I-17 [20] (dashed) and TKAA16 [19] (dashed-dotted). In order to see the effect of higher order corrections and its comparison with the NLO analyses, the BB10 [38] and LLS10 [39] results for the spin-dependent PDFs at NLO pQCD accuracy are also shown in Fig. 2. As one can see from this figure, the x∆u v and x∆d v PDFs are the best determined distributions from the spin-dependent inclusive DIS data, with relatively small uncertainty bands. We should stress here that the uncertainties are computed using the "Hessian method". The contributions from the extrapolated regions, x < 0.001 and x > 0.8, the spindependent PDFs are not directly constrained by the inclusive DIS datasets. For the much better determined x∆u v and x∆d v distributions, one can conclude that the shapes and magnitudes from the SMKA18 QCD fit are generally similar to those found in our previous Mellin based analysis [19,20] The strange quark distribution x∆q turns out to be negative. In contrast to the negative sea quark distribution obtained from SMKA18 analysis of inclusive DIS, DSSV09 [36] and LSS10 [39] fits have shown that the inclusion of the semi-inclusive kaon production data in the QCD fit induces a positive value at the Bjorken value of x > 0.05. It is constrained by a combination of Q 2 evolution, the weak baryon decay constants, and the assumption of an SU(3) symmetric sea, ∆q(x, Q 2 ) ≡ ∆u(x, Q 2 ) = ∆d(x, Q 2 ) = ∆s(x, Q 2 ) = ∆s(x, Q 2 ). A much wider error band has been obtained for the sea quark distribution, as is shown in Fig. 2. As a function of x, the shape of the SMKA18 x∆q PDFs is slightly smaller than the one from KTA-I-17 [20] and larger than the one from TKAA16 [19], which are extracted from the Mellin based analysis. Both the BB10 [38] and the LLS10 [39] strange quark densities are larger then all of our NNLO strange quark densities. From Fig. 2, it is also seen that for the gluon distributions there is a much wider error band than all other spin-dependent quark PDFs. This wide uncertainty is expected to be reduced when the pion and jet production data from spin-dependent proton-proton collisions are included in the analysis [18]. The difficulty in constraining the spin-dependent x∆g(x, Q 2 ) distribution is clearly revealed through the spread of gluon distribution from various global PDFs parametrizations. In Fig. 2, the x∆g(x, Q 2 ) PDFs from our previous NNLO global analyses are compared with the SMKA18 results. Note that, the KTA-I-17 and TKAA16 QCD fits have used spindependent inclusive DIS data, similar to the present analysis. As is seen, in our all fits the x∆g PDFs are positive in all range of x. The LLS10 gluon density approaches zero faster than the other curves and the BB10 shows a sign change at the medium value of x. In most of the QCD fits the gluon distributions is positive at large values of x, with a sign change at smaller x values. However, the latest analyses by DSSV14 [40] and JAM15 [18] of the recent high-statistics jet data from RHIC does also give a positive x∆g distribution with no indication of a sign change in the measured x region.
In the rest of this section, we compare the results of our analysis with the analyzed inclusive spin-dependent DIS data. In Fig. 3, we presented a detailed comparison of SMKA18 NNLO pQCD fit results for the spin-dependent structure functions with the analyzed inclusive DIS data for some selected values of Q 2 = 2, 3, 5 and 10 GeV 2 . The results from Mellin space analyses of KTA-I-17 [20] (dashed) and TKAA16 [19] (dashed-dotted) have been also shown for comparison. As one can see, for all these spindependent PDFs sets, the perturbative expansion is stable and the obtained theory computations agree well with the spin-dependent inclusive DIS data for all range of x. The same comparison for the spin-dependent proton and neutron structure functions are displayed in Figs. 4 and 5. Interestingly, the same pattern is also seen for the neutron and deuteron, and the SMKA18 theory predictions are in good agreement with all DIS data both at low and high values of x.
As a short summary, we have performed a new global QCD analysis of spin-dependent PDFs including all available inclusive DIS data from experiments at CERN, SLAC, HERA and JLAB. The analysis is the first QCD fit at NNLO performed using a developed strategy based on Laplace transform technique. One can conclude that our findings for the theory predictions of spin-dependent structure functions are in good agreements with the analyzed inclusive DIS data. These results indicate the validity of the Jacobi polynomials and Laplace technique for the case of spin-dependent QCD analyses.

VII. SUMMARY AND CONCLUSIONS
The detailed study of the spin-dependent PDFs of the nucleon and the spin structure of hadron are active interdisciplinary research field lying at the crossroads of high-energy hadronic and nuclear physics, with important applications in LHC phenomenology. In this paper, our main aim was to maximally utilize the most available spin-dependent inclusive DIS data over the greatest range of kinematics considering the theoretical perturbative QCD in Laplace s-space. The main outcome of this paper is the first quantitative study of the method of Laplace transform technique in the global QCD fits, for the spin-dependent case. While SMKA18 study should and will be improved on a number of aspects, in particular related to the inclusion of new spin-dependent inclusive DIS data or treatment of TMCs and HT effects, we believe that it opens a new method on the determination of spin-dependent structure of the nucleon. Finally, note that the use of our Laplace technique and Jacobi polynomials approach for incorporating NNLO theory accuracy in the analysis has a number of important implications for future practice and can be extended in a number of directions. Overall, very good descriptions of the global spin-dependent inclusive DIS data set has been obtained in SMKA18 QCD fit, over the entire range of x and Q 2 covered by the preferred cuts. In the near future, the data from 12 GeV JLAB experiments will provide stronger constraints on the behavior of spin-dependent PDFs at large x through precise measurements of polarization asymmetries over a greater range of Q 2 and W 2 [16][17][18]. In addition, the inclusion of semi-inclusive DIS asymmetries as well as jet and pion production asymmetries in spin-dependent proton-proton (pp) collisions, will place stronger constraints on the sea quark polarization [18]. Finally, it is pointed out that the very accurate spindependent DIS data in a wide region require a more careful matching of QCD to the data in order to determine the spin-dependent PDFs correctly.