Proton scalar dipole polarizabilities from real Compton scattering data, using fixed-t subtracted dispersion relations and the bootstrap method

TL;DR

This paper extracts the proton's scalar dipole polarizabilities from real Compton scattering data using fixed-t dispersion relations and a bootstrap fitting method, providing precise values with systematic error analysis.

Contribution

It introduces a bootstrap-based fitting approach combined with statistical tests to analyze Compton scattering data for proton polarizabilities, ensuring robust error propagation and data consistency.

Findings

Extracted $ ext{α}_{E1}$ and $ ext{β}_{M1}$ with uncertainties

Confirmed data set consistency without excluding points

Provided precise polarizability values with statistical significance

Abstract

We perform a fit of the real Compton scattering (RCS) data below pion-production threshold to extract the electric ($\alpha_{E1}$) and magnetic ($\beta_{M1}$) static scalar dipole polarizabilities of the proton, using fixed-$t$ subtracted dispersion relations and a bootstrap-based fitting technique. The bootstrap method provides a convenient tool to include the effects of the systematic errors on the best values of $\alpha_{E1}$ and $\beta_{M1}$ and to propagate the statistical errors of the model parameters fixed by other measurements. We also implement various statistical tests to investigate the consistency of the available RCS data sets below pion-production threshold and we conclude that there are not strong motivations to exclude any data point from the global set. Our analysis yields $\alpha_{E1} = (12.03^{+0.48}_{-0.54})\times 10^{-4} \text{fm}^3$ and $\beta_{M1} = (1.77^{+0.52}_{-0.54})\times 10^{-4} \text{fm}^3$, with p-value $= 12\%$.

Figures

Tables

Table 1. Table 1 : Results of the fits for the static polarizabilities α E ​ 1 subscript 𝛼 𝐸 1 \alpha_{E1} , β M ​ 1 subscript 𝛽 𝑀 1 \beta_{M1} and γ π subscript 𝛾 𝜋 \gamma_{\pi} using the FULL and TAPS data sets and different fit conditions, together with the corresponding χ ^ 2 / d ​ o ​ f superscript ^ 𝜒 2 𝑑 𝑜 𝑓 \hat{\chi}^{2}/dof and p-values.

FULL data set
fit conditions	${\alpha }_{E1}$	${\beta }_{M1}$	${\gamma }_{\pi }$	${\hat{\chi }}^2/dof$ (p-value)
Fit1	${12.00}_{-0.47}^{+0.41}$	${1.80}_{-0.48}^{+0.46}$	fixed	$1.25$ ($3\%$)
Fit $1^{\prime }$	${12.03}_{-0.54}^{+0.48}$	${1.77}_{-0.54}^{+0.52}$	fixed	$1.25$ ($12\%$)
Fit 2	${11.82}_{-0.91}^{+0.81}$	${1.54}_{-1.00}^{+0.95}$	fixed	$1.26$ ($4\%$)
Fit $2^{\prime }$	${11.86}_{-0.99}^{+0.93}$	${1.54}_{-1.05}^{+0.95}$	fixed	$1.26$ ($13\%$)
Fit 3	${12.08}_{-0.64}^{+0.61}$	${1.71}_{-0.75}^{+0.70}$	${8.52}_{-2.93}^{+2.72}$	$1.26$ ($4\%$)
Fit $3^{\prime }$	${12.12}_{-0.77}^{+0.68}$	${1.68}_{-0.79}^{+0.77}$	${8.59}_{-2.94}^{+2.85}$	$1.26$ ($13\%$)
TAPS data set
fit conditions	${\alpha }_{E1}$	${\beta }_{M1}$	${\gamma }_{\pi }$	${\hat{\chi }}^2/dof$ (p-value)
Fit 1	${11.87}_{-0.54}^{+0.50}$	${1.93}_{-0.56}^{+0.52}$	fixed	$1.32$ ($7\%$)
Fit $1^{\prime }$	${11.82}_{-0.55}^{+0.50}$	${1.98}_{-0.56}^{+0.52}$	fixed	$1.32$ ($7\%$)
Fit 2	${11.62}_{-0.96}^{+0.86}$	${1.57}_{-1.04}^{+1.01}$	fixed	$1.34$ ($8\%$)
Fit $2^{\prime }$	${11.89}_{-1.67}^{+1.50}$	${1.76}_{-2.00}^{+1.83}$	fixed	$1.34$ ($8\%$)
Fit 3	${11.74}_{-0.78}^{+0.68}$	${2.06}_{-0.84}^{+0.77}$	${6.95}_{-3.09}^{+2.81}$	$1.34$ ($9\%$)
Fit $3^{\prime }$	${11.67}_{-0.77}^{+0.68}$	${2.12}_{-0.78}^{+0.76}$	${6.85}_{-3.07}^{+2.83}$	$1.34$ ($9\%$)

Table 2. Table 2 : Results for α E ​ 1 subscript 𝛼 𝐸 1 \alpha_{E1} and β M ​ 1 subscript 𝛽 𝑀 1 \beta_{M1} from the test-fit of the FULL and the SELECTED data set, and taking different values for the leading-order spin polarizabilities: the experimental results from Ref.  Martel  et al. ( 2015 ) and the values predicted in HB χ 𝜒 \chi PT  Griesshammer  et al. ( 2012 ) (results in brackets).

	FULL	SELECTED
${\alpha }_{E1}$	$11.99\pm 0.31\left(11.47\pm 0.30\right)$	$11.02\pm 0.33\left(10.46\pm 0.32\right)$
${\beta }_{M1}$	$1.81\pm 0.31\left(2.33\pm 0.30\right)$	$2.78\pm 0.33\left(3.34\pm 0.32\right)$

Table 3. Table 3 : Angular and energy coverage of the available experimental data on unpolarized cross section for proton RCS.

set label	Ref.	first author	points number	${\theta }_{\text{lab}}$ (∘)	$E_{\gamma }$ (MeV)	symbol
1	Oxley (1958)	Oxley	4	$70-150$	$\simeq 60$
2	Hyman et al. (1959)	Hyman	12	$50,90$	$55-95$	$\bullet$
3	Goldansky et al. (1960)	Goldansky	5	$75-150$	$55-80$
4	Bernardini et al. (1960)	Bernardini	2	$\simeq 135$	$\simeq 140$
5	Pugh et al. (1957)	Pugh	16	$50-135$	$40-120$
6	Baranov et al. (1974, 1975)	Baranov	3	$90,150$	$80-110$
7	Baranov et al. (1974, 1975)	Baranov	4	$90,150$	$80-110$
8	Federspiel et al. (1991)	Federspiel	16	$60,135$	$30-90$
9	Zieger et al. (1992)	Zieger	2	$180$	$100,130$
10	Hallin et al. (1993)	Hallin	13	$45-135$	$130-150$
11	MacGibbon et al. (1995)	MacGibbon	8	$90,135$	$95-145$
12	MacGibbon et al. (1995)	MacGibbon	10	$90,135$	$95-145$
13	Olmos de Leon et al. (2001)	Olmos de Leon	55	$60-155$	$60-150$

Table 4. Table 4 : Correlation coefficients ρ 𝜌 \rho among the fit parameters in the different fitting conditions described in Sect.  IV.2 . The columns 2-4 correspond, from the left to the right, to the correlation coefficients between α E ​ 1 subscript 𝛼 𝐸 1 \alpha_{E1} and β M ​ 1 subscript 𝛽 𝑀 1 \beta_{M1} , α E ​ 1 subscript 𝛼 𝐸 1 \alpha_{E1} and γ π subscript 𝛾 𝜋 \gamma_{\pi} , β M ​ 1 subscript 𝛽 𝑀 1 \beta_{M1} and γ π subscript 𝛾 𝜋 \gamma_{\pi} .

FULL data set
fit conditions	${\rho }_{{\alpha }_{E1}-{\beta }_{M1}}$	${\rho }_{{\alpha }_{E1}-{\gamma }_{\pi }}$	${\rho }_{{\beta }_{M1}-{\gamma }_{\pi }}$
Fit1	$-0.64$	$--$	$--$
Fit $1^{\prime }$	$-0.72$	$--$	$--$
Fit 2	$0.59$	$--$	$--$
Fit $2^{\prime }$	$0.52$	$--$	$--$
Fit 3	$-0.84$	$0.86$	$-0.88$
Fit $3^{\prime }$	$-0.87$	$0.84$	$-0.86$
TAPS data set
fit conditions	${\rho }_{{\alpha }_{E1}-{\beta }_{M1}}$	${\rho }_{{\alpha }_{E1}-{\gamma }_{\pi }}$	${\rho }_{{\beta }_{M1}-{\gamma }_{\pi }}$
Fit 1	$-0.74$	$--$	$--$
Fit $1^{\prime }$	$-0.74$	$--$	$--$
Fit 2	$0.47$	$--$	$--$
Fit $2^{\prime }$	$0.23$	$--$	$--$
Fit 3	$-0.85$	$0.82$	$-0.84$
Fit $3^{\prime }$	$-0.86$	$0.81$	$-0.83$

Equations

$$\mathrm{Re}[A_i(\nu ,t)]=A_i^B(\nu ,t)+\left[A_i(0,t)-A_i^B(0,t)\right]+\frac{2}{\pi }{\nu }^2\mathcal{P}{\int }_{{\nu }_0}^{+\infty }d{\nu }^{\prime }\frac{{\mathrm{Im}}_s[A_i({\nu }^{\prime },t)]}{{\nu }^{\prime }({{\nu }^{\prime }}^2-{\nu }^2)},$$

$$A_i(0,t)-A_i^B(0,t)$$

$$\begin{array}{rl}{\alpha }_{E1}& =\frac{-a_1-a_3-a_6}{4\pi },\\ {\gamma }_{E1E1}& =\frac{a_2-a_4+2a_5+a_6}{8\pi M_N},\\ {\gamma }_{M1E2}& =\frac{-a_2-a_4-a_6}{8\pi M_N},\end{array}\begin{array}{rl}{\beta }_{M1}& =\frac{a_1-a_3-a_6}{4\pi },\\ {\gamma }_{M1M1}& =\frac{-a_2-a_4-2a_5+a_6}{8\pi M_N},\\ {\gamma }_{E1M2}& =\frac{a_2-a_4-a_6}{8\pi M_N},\end{array}$$

$${\gamma }_0\equiv -{\gamma }_{E1E1}-{\gamma }_{M1M1}-{\gamma }_{E1M2}-{\gamma }_{M1E2},{\gamma }_{\pi }\equiv -{\gamma }_{E1E1}+{\gamma }_{M1M1}-{\gamma }_{E1M2}+{\gamma }_{M1E2}$$

$${\mathcal{S}}_{ij}=(1+{\delta }_{ij})(E_i+{\gamma }_{ij}{\sigma }_i),$$

$${\chi }_{b,j}^2=\sum _{i=1}^{n_{data}}{\left(\frac{{\mathcal{S}}_{ij}-T_i(p)}{{\sigma }_{ij}}\right)}^2,$$

$${\sigma }_{ij}=(1+{\delta }_{ij}){\sigma }_i.$$

$${{\hat{\chi }}_{b,j}}^2=\sum _{i=1}^{n_{data}}{\left(\frac{{\mathcal{S}}_{ij}-T_i({\hat{p}}_j)}{{\sigma }_{ij}}\right)}^2,$$

$$\hat{p}\equiv \frac{1}{N}\sum _{j=1}^N{\hat{p}}_j,{\sigma }_p\equiv {\left(\frac{1}{N-1}\sum _{j=1}^N({\hat{p}}_j-\hat{p}{)}^2\right)}^{1/2}.$$

$${\hat{\chi }}^2=\sum _{i=1}^{n_{data}}{\left(\frac{E_i-T_i(\hat{p})}{{\sigma }_i}\right)}^2.$$

$${\mathcal{M}}_{ij}=(1+{\delta }_{ij})(T_i(\hat{p})+{\gamma }_{ij}{\sigma }_i),$$

$${\chi }_{th,j}^2=\sum _{i=1}^{n_{data}}{\left(\frac{{\mathcal{M}}_{ij}-T_i(p_j)}{{\sigma }_{ij}}\right)}^2.$$

$${{\chi }_{b,j}}^2=\sum _{i=1}^{n_{data}}{\left(\frac{{\mathcal{S}}_{ij}-T_i(p,f)}{{\sigma }_{ij}}\right)}^2.$$

$${\sigma }_{i,TAPS}\to {\left[{\sigma }_{i,TAPS}^2+{\left(\frac{5}{100}{E}_{i,TAPS}\right)}^2\right]}^{1/2}.$$

$$:{\gamma }_{\pi }^{\text{tot}}=-40.9\pm 0.4\pm 2.2,$$

$$:{\gamma }_{\pi }^{\text{tot}}=-35.9\pm 2.3,$$

$${\xi }_i\equiv \frac{E_i-{\hat{T}}_i}{{\sigma }_i},$$

$${\chi }_{set}^2\equiv \frac{1}{n_{set}}\sum _{i=1}^{n_{set}}{\left(\frac{E_i-{\hat{T}}_i}{{\sigma }_i}\right)}^2.$$

$$\text{no rescaling}:$$

$$\sqrt{{\chi }_{set}^2}\text{ rescaling}:$$

$${\alpha }_{E1}=12.03_{-0.53}^{+0.48},{\beta }_{M1}=1.77_{-0.54}^{+0.52},{\hat{\chi }}^2=1.25(\text{p-value}=12\%),{\rho }_{{\alpha }_{E1}-{\beta }_{M1}}=-\mathrm{0.72.}$$

$${\alpha }_{E1}=11.2\pm 0.4,{\beta }_{M1}=2.5\pm \mathrm{0.4.}$$

Full text

I Introduction

The electric and magnetic static scalar dipole polarizabilities, $\alpha_{E1}$ and $\beta_{M1}$, respectively, are fundamental structure constants of the proton that can be accessed via real Compton scattering (RCS). In the low-energy expansion of the Compton amplitude, they correspond to the leading-order contributions beyond the structure independent terms that describe the scattering process as if the proton were a pointlike particle with anomalous magnetic moment. When approaching the pion-production threshold, also higher-order terms start competing with the scalar dipole polarizabilities. Therefore, one has to resort to reliable theoretical frameworks for extracting the scalar dipole polarizabilities from experimental data. The most accredited theories, which have been used sofar, are fixed-$t$ dispersion relations (DRs), in the unsubtracted L’vov et al. (1997); Babusci et al. (1998); Schumacher (2005) and subtracted Drechsel et al. (1999); Holstein et al. (2000); Pasquini et al. (2007); Drechsel et al. (2003); Pasquini and Vanderhaeghen (2018) formalism, and chiral perturbation theory ($\chi$PT) with explicit nucleons and Delta’s, in the variant of heavy-baryon $\chi$PT (HB$\chi$PT) Bernard et al. (1995); Beane et al. (2003); McGovern et al. (2013) and manifestly covariant Lensky et al. (2015); Lensky and Pascalutsa (2010) $\chi$PT (B$\chi$PT) 111We refer to Kondratyuk and Scholten (2001); Gasparyan et al. (2011) for other theoretical predictions of the scalar and spin polarizabilities not fitted to experimental data.. Based on these theoretical frameworks, extractions of the scalar dipole polarizabilities have been obtained by fitting different data sets for the unpolarized RCS cross section, and adopting a statistical approach based on the conventional $\chi^{2}$-minimization procedure. Recently, a new statistical method has successfully been applied in Ref. Pasquini et al. (2018) to analyze RCS data at low energies and extract values for the energy-dependent scalar dipole dynamical polarizabilities Griesshammer and Hemmert (2002); Hildebrandt et al. (2004). The method is based on the parametric-bootstrap technique, and it is adopted in this work to extract the scalar dipole static polarizabilities, using the updated version of fixed-$t$ subtracted DRs formalism Pasquini and Vanderhaeghen (2018) as theoretical framework. Although the bootstrap method is rarely used in nuclear physics Navarro Pérez and Lei (2018); Pasquini et al. (2018); Navarro Pérez et al. (2014); Nieves and Ruiz Arriola (2000); Bertsch and Bingham (2017), it has high potential and advantages Pastore (2018). In particular, we will show that it allows us to include the systematic errors in the data analysis in a straightforward way and to efficiently reconstruct the probability distributions of the fitted parameters. We will also pay a special attention to discuss the available sets of RCS data below pion-production threshold. Following recent discussions about the possible presence of outliers in the available data sets Krupina et al. (2018); Griesshammer et al. (2012), we perform several tests to judge the data-set consistency.

The manuscript is organized as follows. In Sec. II, we briefly summarize the theoretical framework of fixed-$t$ subtracted DRs. In Sec. III, we describe the main features of the parametric-bootstrap technique, which is applied in Sec. IV to our specific case to fit $\alpha_{E1}$ and $\beta_{M1}$. We perform the fit in different conditions, i.e., switching on/off the effects of the systematic errors, using the constraint of the Baldin’s sum rule for the polarizability sum and including the backward spin polarizability $\gamma_{\pi}$ as additional fit parameter. The consistency of the data set is discussed in Sec. V, where we perform different statistical tests to identify the possible presence of outliers. The results of our analysis are summarized in Sec. VI, in comparison with available extractions of the scalar dipole polarizabilities. Our conclusions are drawn in Sec. VII. In App. A, we give the complete list of the existing data sets of RCS below pion-production threshold, and in App. B we discuss the values of the correlations among the fit parameters in all the different conditions discussed in this work.

II Theoretical framework

We consider RCS off the proton, i.e. $\gamma(q)+P(p)\rightarrow\gamma(q^{\prime})+P(p^{\prime})$, where the variables in brackets denote the four-momenta of the participating particles. The familiar Mandelstam variables are $s=(p+q)^{2}$, $u=(q-p^{\prime})^{2}$ and $t=(q-q^{\prime})^{2}$, and are constrained by $s+u+t=2M^{2}$, with $M$ the proton mass. The RCS amplitude can be described in terms of 6 Lorentz invariant functions $A_{i}(\nu,t)$, which depend on the crossing-symmetric variable $\nu=(s-u)/4M$ and $t$. They are free of kinematical singularities and constraints, and because of the crossing symmetry they obey the relation $A_{i}(\nu,t)=A_{i}(-\nu,t)$. Assuming analyticity, they satisfy the following fixed-$t$ subtracted DRs (with the subtraction point at $\nu=0$) Drechsel et al. (1999); Pasquini et al. (2007)

[TABLE]

where $\nu_{0}$ is the pion-production threshold, and $A_{i}^{B}(\nu,t)$ is the Born term, corresponding to the pole diagrams involving a single nucleon exchanged in $s$- or $u$-channel and $\gamma NN$ vertices taken in the on-shell regime. In Eq. (1), the subtraction functions $\left[A_{i}(0,t)-A_{i}^{B}(0,t)\right]$ can be determined by once-subtracted DRs in the $t$ channel:

[TABLE]

where $A_{i}^{t-pole}(0,t)$ represents the contribution of the poles in the $t$ channel, that amounts to the $\pi^{0}$-pole contribution to the $A_{2}$ amplitude. The subtraction constants $a_{i}\equiv\left[A_{i}(0,0)-A_{i}^{B}(0,0)\right]$ are directly related to the scalar dipole and leading-order spin polarizabilities, i.e.

[TABLE]

with the combination

[TABLE]

defining the forward ($\gamma_{0}$) and backward ($\gamma_{\pi}$) spin polarizabilities. We will consider $\{\gamma_{E1E1},\gamma_{M1M1},\gamma_{0},\gamma_{\pi}\}$ as independent set of spin polarizabilities.

In the actual calculation, the $s$-channel imaginary parts in Eq. (1) are evaluated using the unitarity relation, taking into account the contribution of the $\pi N$ intermediate states from the latest version of the MAID pion-photoproduction amplitudes Drechsel et al. (2007) and approximating the contribution from multipion intermediate channels by the inelastic decay channels of the $\pi N$ resonances, as detailed in Ref. Drechsel et al. (1999). It was found, however, that in the subtracted dispersion relation formalism, the sensitivity to the multipion channels is very small and that subtracted dispersion relations are essentially saturated at $\nu\approx 0.4$ GeV. Furthermore, the $t$-channel imaginary parts in Eq. (2) are calculated using the $\gamma\gamma\rightarrow\pi\pi\rightarrow N\bar{N}$ channel as input for the positive-$t$ cut, while the negative-$t$ cut is strongly suppressed for low values of $t$. The last one can be approximated by the contributions of $\Delta$-resonance and non-resonant $\pi N$ intermediate states in the $s$-channel, which are then extrapolated into the unphysical region at $\nu=0$ by analytical continuation. For more detail in the implementation of the unitarity relations, we refer to the original work Drechsel et al. (1999). Having determined the contributions of the $s$- and $t$-channel integrals, the only remaining unknown are the subtraction constants, i.e. the leading-order static polarizabilities. In principle, all the six leading static polarizabilities can be used as free fit parameters to the Compton observables. However, a simultaneous fit of all them is not feasible at the moment, because of the limited statistics of the available RCS data. In the following, we will limit ourselves to the data sets for unpolarized RCS below pion-production threshold, and consider different variants of fits for two sets of parameters, i.e. $\{\alpha_{E1},\beta_{M1}\}$ or $\{\alpha_{E1},\beta_{M1},\gamma_{\pi}\}$. The remaining constants which do not enter the fit are fixed as described in Sec. IV.

III The fitting method

We consider a generic problem, where we have a model prediction $T(p)$ for an observable, which depends on a set $p$ of parameters, and we want to find the optimal set $\hat{p}$ that better reproduces the available experimental data. We adopt an algorithm based on the parametric bootstrap technique Davidson and Hinkley (1997), i.e., $N$ Monte Carlo replicas of experimental data are produced and a fit of the set $p$ is performed to every bootstrapped data sample. After every cycle $j$, the best values $\hat{p}_{j}$ are stored, to obtain $N$ outcomes of the (unknown) probability distribution of $p$.

In our case we assume that:

every data point is Gaussian distributed with a mean equal to the measured value and a standard deviation given by the experimental (statistical) error; 2. 2.

data points are affected by systematic errors given by different rescaling factors of the data in each subset; 3. 3.

when not explicitly stated otherwise by the experimental groups, every source of systematic error follows an uniform distribution and the published value gives the full estimated interval. If there are more sources, we take the product of such random uniform variables; 4. 4.

the sample in every data subset is independent from the other subsets.

This sampling method can then be written as

[TABLE]

where $S_{ij}$ is a generic bootstrapped point with the index $i$ running over the number of data point ($n_{data}$) and $j$ running over the number of replicas ($N$). $E_{i}$ is the experimental point having an uncertainty $\sigma_{i}$, $\gamma_{ij}$ is the Gaussian normal variable needed for the statistical sampling and $\delta_{ij}$ is a box distributed variable that quantifies the effect of the systematic uncertainties for each data subset independently. Considering a generic subset, labeled with $k$ ($k$ runs from 1 to the number of the different data subsets $n_{set}$) and composed of $n_{k}$ data points, we take $\delta_{ij}\in\mathcal{U}[-\Delta_{k},\Delta_{k}]$ for $i=1,\dots,n_{k}$, where $\pm\Delta_{k}$ is the published systematic error and $\sum_{k=1}^{n_{set}}n_{k}=n_{data}$. If there are $n_{s}$ different and independent sources of systematic uncertainties, $\delta_{ij}$ is the product of all the $n_{s}$ box distributed variables, i.e., $\delta_{ij}=\prod_{f=1}^{n_{s}}\mathcal{U}[-\Delta_{f},\Delta_{f}]$. The systematic sources can be easily excluded from this procedure by just imposing $\delta_{ij}\equiv 0$ in Eq. (5).

The minimization function at the $j^{th}$ iteration is given by

[TABLE]

where

[TABLE]

The minimum in the parameter space can be defined as

[TABLE]

where $\hat{p}_{j}$ are the best values of the fit parameters $p$ at the $j^{th}$ bootstrap cycle.

Repeating this minimization for $N$ cycles, the empirical distribution $\mathcal{P}(\hat{p}_{j})$ of the $\hat{p}_{j}$ random variables gives an estimate of the true probability distribution $\mathcal{P}(p)$ that includes the propagation of both statistical and systematic errors of the experimental data. The best value and the standard deviation of $p$ can be then simply obtained as:

[TABLE]

The goodness of this fit procedure can be estimated in the same way as in the standard case, using the value $\hat{\chi}^{2}$ of the so-called $\chi^{2}$-variable, defined as 222 The link between ${\hat{\chi}_{b,j}}^{2}$ and $\hat{\chi}^{2}$ can be found in Ref. Pedroni et al. .:

[TABLE]

It is worthwhile to notice here that $\hat{\chi}^{2}$ is distributed according to the $\chi^{2}$ distribution only when $\delta_{ij}=0$, i.e. when all the $E_{i}$ are independent random gaussian variables.

Within the bootstrap framework, it is also possible to evaluate the expected theoretical probability distribution associated to ${\hat{\chi}}^{2}$ by replacing $\mathcal{S}_{ij}$ in Eq. (6) with

[TABLE]

and by finding, at each bootstrap cycle, the minimum value $\hat{\chi}_{th,j}^{2}$ of the following function

[TABLE]

After $N$ bootstrap iterations, we are able to empirically reconstruct the probability distribution $\mathcal{P}(\chi_{th}^{2})$ and then to evaluate the final p-value associated to the fit.

It can be easily demonstrated (see Pedroni et al. ) that, when $\delta_{ij}=0$ in Eq. (11), $\mathcal{P}(\chi_{th}^{2})$ coincides with the $\chi^{2}$ distribution, as expected. In any case, we stress that the bootstrap method allows us to obtain a p-value for $\hat{\chi}^{2}$ directly from the evaluated $\mathcal{P}(\chi_{th}^{2})$ distribution, also when systematic errors are taken into account in the fit procedure.

III.1 Uncertainties on additional model parameters

In the most generic situation, the model $T$ may depend on an additional set of parameters $f$ besides the fit parameters $p$, i.e., $T\equiv T(p,f)$. The ${\chi_{b,j}}^{2}$ variable of Eq. (6) is consequently modified as

[TABLE]

Suppose the values of the parameters $f$ are derived from experimental data and are known within an experimental uncertainty $\sigma_{f}$. Within the bootstrap framework, we can easily evaluate how the uncertainties $\sigma_{f}$ affect the values of the fit parameters $\hat{p}$, without using the error-propagation procedure that would require performing numerical derivatives $\partial T/\partial f$. At each $j^{th}$ bootstrap cycle, we can sample the value $f_{j}$ of the model parameters from their known probability distribution, which in the following will be considered to be a Gaussian defined as $\mathcal{G}[f,\sigma^{2}_{f}]$. Then, we can repeat the procedure described above by replacing $T_{i}(\hat{p}_{j})$ with $T_{i}(\hat{p}_{j},f_{j})$ in Eq. (8), and evaluate all the relevant fit parameters.

IV Fit to RCS data

In this section, we apply the fitting method introduced in Sec. III to analyze available RCS data below pion-production threshold. We use fixed-$t$ subtracted DRs for the model predictions, which contain the leading-order static polarizabilities as free parameters, as explained in Sec. II. We discuss two data sets, corresponding to the FULL and TAPS data sets, as described in App. A. Furthermore, we consider different fit conditions, switching on/off the systematic errors and using two sets of free parameters: i) the scalar dipole polarizabilities, with and without the constraint of the Baldin’s sum rule for the polarizability sum $\alpha_{E1}+\beta_{M1}$, and ii) the scalar dipole polarizabilities constrained by the Baldin’s sum rule along with the backward spin polarizability $\gamma_{\pi}$. For the Baldin’s sum rule, we use the weighted average over the available evaluations reported in Ref. Hagelstein et al. (2016), which coincides also with the value used in the fit of Refs. Olmos de Leon et al. (2001); McGovern et al. (2013); Martel et al. (2015), i.e., $\alpha_{E1}+\beta_{M1}=13.8\pm 0.4$. The remaining parameters of fixed-$t$ subtracted DRs are fixed to the experimental values extracted from double polarization RCS Martel et al. (2015), i.e. $\gamma_{E1E1}=-3.5\pm 1.2$ and $\gamma_{M1M1}=3.16\pm 0.85$, and from the GDH experiments Ahrens et al. (2001); Dutz et al. (2003), i.e. $\gamma_{0}=-1.01\pm 0.08\pm 0.10$ 333This value is consistent with the fitting conditions adopted for the extraction of the spin polarizability in Ref. Martel et al. (2015). We note that recent reevaluations Pasquini et al. (2010); Gryniuk et al. (2016) of $\gamma_{0}$ give a slightly smaller central values, with uncertainties consistent with the value used in Ref. Martel et al. (2015).. When the backward spin polarizability is not used as fit parameter, we fixed it to the weighted average of the values extracted at MAMI Schumacher (2005), i.e. $\gamma_{\pi}=-8.0\pm 1.8$. Here and in the following, we used the standard convention to exclude the $t$-channel $\pi^{0}$-pole contribution from the spin polarizabilities. These contributions amount to ${\gamma_{\pi}}^{\pi^{0}-\text{pole}}=-46.7$ Pasquini et al. (2007), $\gamma_{M1M1}^{\pi^{0}-\text{pole}}=-\gamma_{E1E1}^{\pi^{0}-\text{pole}}=\tfrac{1}{4}\gamma_{\pi}^{\pi^{0}-\text{pole}}$, while they vanish in the case of the forward spin polarizability. Finally, for each fitting configuration, we discuss the probability distributions of the fitted parameters and the p-values of the $\hat{\chi}^{2}$ variable. Here and in the following, we use the units of $10^{-4}$ fm3 for the scalar dipole polarizabilities and $10^{-4}$ fm4 for the spin polarizabilities.

IV.1 Handling the experimental and model errors

We apply the method described above using $N=10000$ bootstrap replicas. Within this framework and following the method outlined in Sec. III.1, we take into account the uncertainties of the model parameters on the values of the polarizabilities not treated as free parameters in the fit procedure. In particular, we take $\gamma_{0}\in\mathcal{G}[-1.01,0.13^{2}]$444The uncertainty value $0.13^{2}$ is the sum of the squares of the statistical and systematic errors., $\gamma_{E1E1}\in\mathcal{G}[-3.5,1.2^{2}]$ and $\gamma_{M1M1}\in\mathcal{G}[3.16,0.85^{2}]$. When keeping fixed the backward spin polarizability, we propagate the error of $\gamma_{\pi}$ using $\gamma_{\pi}\in\mathcal{G}[8.0,1.8^{2}]$. Furthermore, the Baldin’s sum rule constraint is implemented using $\alpha_{E1}+\beta_{M1}\in\mathcal{G}[13.8,0.4^{2}]$. The uncertainties on the fitted $\alpha_{E1}$ and $\beta_{M1}$ thus automatically include the propagation of the errors of the spin polarizabilities and the Baldin’s sum rule. The statistical and systematic uncertainties of the experimental data are taken into account as described in Sec. III, except for the TAPS data points Olmos de Leon et al. (2001). As discussed in Ref. Griesshammer et al. (2012), they are affected by a $5\%$ point-to-point systematic error, and, accordingly, the statistical error of each point is modified as follows

[TABLE]

IV.2 Results

We discuss in this section the results of the fit, performed under several configurations:

• •

Fit 1: with Baldin’s sum rule, and systematic errors excluded: $\alpha_{E1}-\beta_{M1}$ as free parameter;

• •

Fit $1^{\prime}$: with Baldin’s sum rule, and systematic errors included: $\alpha_{E1}-\beta_{M1}$ as free parameter;

• •

Fit 2: without Baldin’s sum rule, and systematic errors excluded: $\alpha_{E1}$ and $\beta_{M1}$ as free parameters;

• •

Fit $2^{\prime}$: without Baldin’s sum rule, and systematic errors included: $\alpha_{E1}$ and $\beta_{M1}$ as free parameters;

• •

Fit $3$: with Baldin’s sum rule, and systematic errors excluded: $\alpha_{E1}-\beta_{M1}$ and $\gamma_{\pi}$ as free parameters;

• •

Fit $3^{\prime}$: with Baldin’s sum rule, and systematic errors included: $\alpha_{E1}-\beta_{M1}$ and $\gamma_{\pi}$ as free parameters.

All these different fits are performed using both the FULL and TAPS data sets. The corresponding results are summarized in Table 1 and shown in Figs. 1-3. In all the cases, the probability distributions of the fit parameters are very similar to Gaussian functions.

A few comments are in order:

• •

the values of the fitted $\alpha_{E1}$ and $\beta_{M1}$ depend on the choice of the data set, but are all consistent within the uncertainties;

• •

the sum of the values of $\alpha_{E1}$ and $\beta_{M1}$ from the Baldin-unconstrained fit is well compatible, within the fit errors, with the Baldin’s sum rule value

• •

the inclusion of systematic errors does not change the central values of the fitted parameters, but increases their uncertainties. This effect is mostly visible for the TAPS data set fitted in the Fit 2 and Fit $2^{\prime}$ conditions, while it is reduced for the FULL data set, where the effects of the systematic errors in the different subsets are, at least partially, compensated (see Figs. 1-3).

• •

when the systematic errors are taken into account, the central values of the $\hat{\chi}^{2}/dof$ do not change. However, the corresponding p-values significantly change for the FULL data set, since higher values of $\hat{\chi}^{2}/dof$ are more likely to occur. This effect is clearly visible from the cumulative distribution functions (CDFs) of $\hat{\chi}^{2}$ shown in Figs. 4 and 5. When we fit a single data set, as in the case of the TAPS data set, the systematic error becomes a common scale factor for all the data points and it does not change the p-value. Therefore, the main effect of the systematic-error propagation is to increase the statistical errors on the fitted parameters (see Fig. 3);

• •

the fitted values of $\gamma_{\pi}$ in the Fit $3^{\prime}$ conditions and with the additional contribution from the $\pi^{0}$-pole, i.e. $\gamma_{\pi}^{\text{tot}}=\gamma_{\pi}+\gamma_{\pi}^{\pi^{0}-\text{pole}},$ are in very good agreement with the values extracted within the fixed-$t$ unsubtracted DR analysis Schumacher (2005); Wolf et al. (2001); Camen et al. (2002); Olmos de Leon et al. (2001):

[TABLE]

The results of Refs. Wolf et al. (2001); Camen et al. (2002) are obtained using data above the pion-production threshold, while the result of Ref. Olmos de Leon et al. (2001) is extracted from the complete TAPS data set, ranging up to photon energies of $165$ MeV.

From all the results we conclude that the inclusion of the systematic errors in the fitting procedure is very important since, in this case, the p-value associated to the $\hat{\chi}^{2}$ changes significantly (see Fig. 4), while the uncertainty on the fitted parameters changes in a much less pronounced way. This behavior can be observed only thanks to the bootstrap method, since it is not possible to compute the correct p-values without resorting to the Monte Carlo replicas, once included the systematic errors.

As mentioned before, the $\hat{\chi}^{2}$ parameter obtained in the different fitting conditions is never distributed like a $\chi^{2}$ probability function. This effect is due to the correlation present among all the points in each subset, and is clearly visible from the CDFs shown in Fig. 4.

V Consistency checks on the available data set

The scientific community has not reach so far a common agreement on the definition of the data set of proton RCS below pion-production threshold Pasquini et al. (2018); Griesshammer et al. (2012); Krupina et al. (2018). As pointed out in Ref. Krupina et al. (2018) and in Sec. IV.2 of this work, the values obtained from a fit of $\alpha_{E1}$ and $\beta_{M1}$ strongly depend on the choice of the data set. In this section, we apply a few basic statistical tests to investigate the consistency of the data set and the possible occurrence of outliers. We will discuss the FULL data set, the TAPS data set and the SELECTED data set, as an example of selection from the complete data set, all below the pion-production threshold as detailed in App. A.

In this section, we will use the standard minimization technique and we will not take systematic errors into account, in order to work in a well-established fitting condition and investigate the pure statistical features of the experimental data. We set the Fit 1 condition (i.e., we use $\alpha_{E1}-\beta_{M1}$ as free parameter and we neglect the systematic errors) and we use the conventional $\chi^{2}$ of Eq. (10) as minimization function, without implementing the bootstrap procedure. Therefore, the errors on the fixed spin polarizabilities are not included in the results of the tests, while the uncertainties of $\alpha_{E1}+\beta_{M1}$ affect the electric and magnetic polarizabilities errors as $\epsilon_{\alpha_{E1},\beta_{M1}}=\sqrt{(\epsilon_{\alpha_{E1}+\beta_{M1}}^{2}+\epsilon_{\alpha_{E1}-\beta_{M1}}^{2})}/2$. We will refer to this fitting configuration as test-fit. The result of the test-fit applied to the FULL data set leads to the best values $\alpha_{E1}=11.99\pm 0.31$ and $\beta_{M1}=1.81\pm 0.31$, which are almost identical to the values given in Table 1, with the Fit 1 condition applied to the FULL data set. The tiny difference in the central values and the different statistical errors are due to the propagation of the uncertainties of those polarizabilities that are not treated as free parameters in the bootstrap fit.

V.1 The Jackknife resampling

A possible strategy to discuss the consistency of the data set is the Jackknife, a resampling technique that can be considered as a particular case of the non-parametric bootstrap technique. Given a data set $D=\{d_{i}\},i=1,\dots,n$, composed by $n$ points, we can define $n$ data subsets by removing one datum at a time, i.e., $D_{k}=D\setminus\{d_{k}\}$, where $k=1,\dots,n$. We then fit the model $T(p)$ to every $D_{k}$ data set, obtaining a best value of the parameters $\hat{p}_{k}$ for each set. From the $n$-tuple of $\hat{p}_{k}$, we can compute the average $\hat{p}_{Jack}$ and its sample standard deviation $\sigma_{Jack}$. An outlier $k$ is expected to give a result far from the average value, i.e., $|\frac{\hat{p}_{k}-\hat{p}_{Jack}}{\sigma_{Jack}}|\gg 1$. Instead, if there are no evident outliers, we expect that all the variables $\hat{p}_{k}$ follow, at least approximately, Gaussian confidence levels Hinkley (1977). In this way, we can identify possible deviations of a data subset from the other ones.

We apply the Jackknife to the FULL, TAPS and SELECTED data sets: the best values of $\alpha_{E1}$ and $\beta_{M1}$ versus the index $k$ of the excluded point in each subset are plotted in Fig. 6. In the case of the FULL data set, we note that the statistical fluctuations are well in agreement with the expected Gaussian confidence levels ($\sim 95\%$ of the occurrences within the $2\sigma$ range). We can then conclude that there is no clear evidence of outliers.

In the case of the SELECTED data set, we obtain very similar results, with less pronounced fluctuations ($\sim 98\%$ of the occurrences within the $2\sigma$ range). This does not necessarily implies that there is an improvement in the data set. Instead, this behavior may simply reflect the fact that the data points excluded from the set are not ”close enough to” the model predictions.

The same test applied to the TAPS data set shows a clear dependence of the values of $\alpha_{E1}$ and $\beta_{M1}$ on the scattering angle. This feature is due to the fact that the data are ordered by increasing scattering angles and that the sensitivity of the unpolarized RCS cross section to $\alpha_{E1}-\beta_{M1}$ is higher in the backward scattering region: when a single datum is removed in the backward region, the value of $\beta_{M1}$ decreases and $\alpha_{E1}$ increases.

V.2 Residual analysis

In order to cross-check the stability of the FULL data set, we performed the analysis of the residuals, defined as

[TABLE]

where $E_{i}$ is the $i^{th}$ experimental datum with the uncertainty $\sigma_{i}$, and $\hat{T}_{i}$ is the model prediction obtained with the best values of the fitted parameters. If the model is able to correctly describe the experimental datum, the value $E_{i}$ can be considered a possible outcome of the probability distribution of $T_{i}$. In this case, the variable $\xi_{i}$ of Eq. (16) is Gaussian distributed as $\mathcal{G}[0,1]$.

The residual analyses for the FULL and SELECTED data sets are shown in Fig. 7, together with the q-q plots, representing the CDF($\xi_{i}$) vs CDF($z$), with $z$ a Gaussian distributed variable according to $\mathcal{N}[0,1]$. The variable $\xi_{i}$ has mean value and standard deviation in fairly good agreement with the expectations. In the case of the SELECTED data set, we observe again less pronounced statistical fluctuations mostly due to the exclusions of the subsets 1 Oxley (1958) and 7 Baranov et al. (1975). This is also shown by the fact that the CDF($\xi_{i}$) for the SELECTED data set approaches the maximum value of unity faster than in the case of the FULL data set.

More precisely, the FULL data set shows three points lying outside the $3\sigma$ range: this configuration can happen with a $1\%$ probability, assuming only Gaussian random choice. Such a value, even if not extremely low, points to possible outliers inside the data set. On the other hand, the SELECTED data set has only 2 points outside the $2\sigma$ range; in this case the associated probability is $\simeq 3\%$. Such a low value seems to indicate that at least some of the discarded points were not outliers.

The FULL and SELECTED data sets have a very similar statistical significance and this ambiguity can be only resolved with new sets of precise and accurate data, especially in the in the backward angular region where the sensitivity to $\alpha_{E1}-\beta_{M1}$ is higher.

All these conclusions are in agreement with the results obtained with the bootstrap method and shown in Table 1. Under the Fit 1 condition, taking only into account the statistical (Gaussian) errors, the p-value of the fit is about $3\%$, which is very close to the $1\%$ occurrence probability. However, when the systematic errors are included in the fitting procedure (Fit 1′), the statistical significance strongly increases ($12\%$), thus indicating that the occurrence probability of the FULL data set is higher, when taking properly into account all the data error sources. Moreover, since there is not a clearly identified source of possible experimental problems that could affect the data discarded in the SELECTED data base, we prefer not to exclude any point to keep the highest sensitivity to $\alpha_{E1}-\beta_{M1}$. Instead, we propose to treat the suspicious points with the approach outlined below.

V.3 The $\chi^{2}$ per set

Given a data set composed by subsets with $n_{set}$ points, we can define for each subset the following variable

[TABLE]

If the model $\hat{T}_{i}$ is able to well describe the data, all the $\chi^{2}_{set}$ values should be fairly close to one. Viceversa, if $\chi^{2}_{set}\gg 1$, we cannot automatically deduce that a data subset should be excluded. This parameter is evaluated using the particular model used in the fit procedure and a bias may be introduced by using large values of $\chi^{2}_{set}$ as criterion to exclude data sets.

We applied this kind of analysis to the FULL data set and the results are shown in Fig. 8. We can notice that most of the subsets have $\chi^{2}_{set}\approx 1$, while the subsets 1 Oxley (1958) and 7 Baranov et al. (1975) give higher $\chi^{2}_{set}$ values. As mentioned before, these subsets are indeed excluded in the definition of the SELECTED data set. However, both data sets have only 4 points each and with such a small number of points we can not exclude the occurrence of pure statistical fluctuations, as mentioned in the previous section.

An alternative method, first suggested in Birge (1932), is to rescale the statistical errors of the points of each data subset by a factor $\sqrt{\chi^{2}_{set}}$ and to repeat again the fit procedure (see also Behnke et al. (2013)). This relies on the assumption that a large $\chi^{2}_{set}$ value indicates underestimated measurement uncertainties that should be equally attributed to all the points of a given subset. We then obtain new values for the fitted parameters ${\hat{p}}^{\prime}$ with the minimum of the $\chi^{2}$ function equal to 1, by construction. This strategy is again model dependent, but it can be used as an indication for the identification of outliers. If there are no data subsets that behave as outliers and then could determine very different values for the fitted parameters, we would expect that ${\hat{p}}^{\prime}\simeq\hat{p}$.

In our case, the values of the fitted parameters obtained from the FULL data set with and without rescaling of the statistical errors are consistent within the (large) fit errors, i.e.

[TABLE]

If we exclude from the fit the subsets 1 and 7, without rescaling the errors, we obtain the value $\alpha_{E1}-\beta_{M1}=9.01\pm 0.50$, which is very similar to the result in Eq. (19) obtained with the rescaling method. As a matter of fact, the rescaling method is equivalent to reduce the impact of the data points in the subsets 1 and 7, since these points are weighted by their relatively high $\sqrt{\chi^{2}_{set}}$ factor. However, the rescaling does not lead to a significant improvement in the accuracy of the extraction, since the error bars in Eqs. (18) and (19) are very similar. Furthermore, the difference between the central values in Eqs. (18) and (19) can be related to the angular distribution of the experimental data of sets 1 and 7, which is mainly in the backward scattering region, where the sensitivity to $\alpha_{E1}-\beta_{M1}$ is higher.

Given all these findings, we once more conclude that there is no clear evidence that these sets are outliers that should be excluded from the fit.

V.4 Behavior of the minimization function

In order to investigate the effect on the fit results of the exclusion of some data points, we examined the behavior of the minimization function versus the values of the fit parameters $\alpha_{E1}$ and $\beta_{M1}$.

The results for the FULL data set (with and without rescaling the statistical error by a factor $\sqrt{\chi^{2}_{set}}$), for the TAPS data set and for the SELECTED data set Griesshammer et al. (2012) are shown in Fig. 9.

When outliers are discarded from the fit, we would expect a significant reduction of the minimum of the $\chi^{2}$ function as well as a more pronounced convexity, corresponding to smaller errors for the fitted parameters. In the case of the SELECTED data set, we indeed observe that the minimum value of the reduced $\chi^{2}$ function is closer to 1, but the shape of the minimization function is the same as in the case of the FULL data set, i.e., the errors on the fitted parameters remain ultimately the same.

This simple analysis gives another additional hint that there is no clear evidence of the presence of outliers and no strong enough motivation for the exclusion of some data points from the FULL data set.

V.5 Summary of the tests

All the previous consistency tests led us to the conclusion that there are no strong motivations for the exclusion of any data point from the global RCS data set below pion-production threshold, even though we observed significant deviations for a few data points at the backward scattering angles. The residual analysis lets us to conclude that the FULL and the SELECTED data sets have almost the same statistical significance and also the fit errors are basically the same in both cases.

An alternative approach to handle the suspicious points is to rescale their statistical errors by a factor $\sqrt{\chi^{2}_{set}}$ rather than exclude them from the data set. Also the bootstrap fitting technique is useful in these cases to check if the inclusion of systematic errors in the fitting procedure improves the significance of the obtained results.

Since most of these points are located in the backward scattering region, where the RCS unpolarized cross section has the larger sensitivity to $\alpha_{E1}-\beta_{M1}$, their exclusion may lead to biased results.

We conclude that the main reason of the sizeable uncertainties that are present at the moment in the extraction of the scalar polarizabilities and especially of $\beta_{M1}$ are mainly due to the intrinsic limitations (poor accuracy and scarcity) of the data set at our disposal.

VI Available extractions of RCS scalar dipole static polarizabilities

In Fig. 10, we collect the available results for the extraction of the scalar dipole static polarizabilities from RCS at low energies. The red solid curve show the results from this work, obtained from the bootstrap-based fit using the FULL data set with the constraint of the Baldin’s sum rule and taking into account the effects of the systematic errors of the experimental data and the propagation of the statistical errors of the fixed polarizabilities $\alpha_{E1}+\beta_{M1}$, $\gamma_{0}$, $\gamma_{\pi}$, $\gamma_{E1E1}$ and $\gamma_{M1M1}$ (Fit $1^{\prime}$ conditions). Within our fitting technique, we are able to evaluate the correlation coefficient $\rho_{\alpha_{E1}-\beta_{M1}}$ among $\alpha_{E1}$ and $\beta_{M1}$: this determines the ellipse-shape in Fig. 10. All the other correlation terms are given in App. B. Numerically, we obtain the following best values with the 68% confidence-level error bar

[TABLE]

These results are in very good agreement with the ones obtained using a traditional $\chi^{2}$ fitting procedure in a fixed-$t$ subtracted DRs framework Drechsel et al. (2003). The experimental fits shown by black curves have been obtained within unsubtracted DRs MacGibbon et al. (1995); Federspiel et al. (1991); Olmos de Leon et al. (2001). The light-green band shows the experimental constraint on the difference $\alpha_{E1}-\beta_{M1}$ from Zieger et al. Zieger et al. (1992). The green solid curve shows the B$\chi$PT predictions of Ref. Lensky et al. (2015). The blue solid curve corresponds to the $68\%$ ellipse of the Baldin constrained fit of Ref. McGovern et al. (2013); Griesshammer et al. (2016), using the SELECTED data set and the HB$\chi$PT framework. These results are in excellent agreement also with the fit within B$\chi$PT of Ref. Lensky and McGovern (2014). We also show the latest value from PDG Patrignani et al. (2016) (solid black disk):

[TABLE]

They differ from the 2012 and earlier editions by inclusion of the data fit analysis within HB$\chi$PT McGovern et al. (2013).

We note that there is a discrepancy between the values obtained in the framework of effective field theories McGovern et al. (2013); Griesshammer et al. (2012); Lensky et al. (2015) and the results obtained using DRs, even if they are compatible within the $2\sigma$-range.

In order to shed some light on the origin of the difference between the results from the extraction within HB$\chi$PT and fixed-$t$ subtracted DRs, we performed some test-fits, in the condition described in Sec. V, using fixed-$t$ subtracted DRs with input from the central values of HB$\chi$PT predictions for the spin polarizabilities. The results for the leading-order spin polarizabilities in HB$\chi$PT read McGovern et al. (2013); Griesshammer et al. (2016) $\gamma_{E1E1}=-1.1\pm 1.9$, $\gamma_{M1M1}=2.2\pm 0.5\text{(stat)}\pm 0.6$, $\gamma_{0}=-2.6\pm 0.5\text{(stat)}\pm 1.8$, and $\gamma_{\pi}=5.6\pm 0.5\text{(stat) }\pm 1.8$, and are quite different from the experimental values used in our DR analysis. On top of that, we noticed a different evaluation for the $\pi^{0}$-pole contribution calculated in Ref. McGovern et al. (2013), which is $-45.9$ for $\gamma_{\pi}$. In Table 2, we compare the test-fit values for $\alpha_{E1}$ and $\beta_{M1}$ in the case we use the results of the spin polarizabilities and the $\pi^{0}$pole from the experimental extraction Martel et al. (2015) or the corresponding values from HB$\chi$PT McGovern et al. (2013); Griesshammer et al. (2016), with the $\pi^{0}$-pole contribution reported in McGovern et al. (2013) (results in brackets). This analysis has been performed for both the FULL and SELECTED data sets, in order to investigate the dependence of the results not only on the values of the spin polarizabilities, but also on the choice of the data set (see Ref. Krupina et al. (2018) for a more comprehensive discussion). If we focus on the central values of $\beta_{M1}$, we notice that the different input for the spin polarizabilities affects the results by 20-30$\%$, while the choice of the data set leads to a 40-50$\%$ increase. It is certainly too simplistic to estimate the model dependence of the two extractions with the different values of the spin polarizabilities. However, in the energy range below pion production threshold, this gives a rather good indication of the main effects due to the model dependence.

The results for the RCS differential cross section obtained with the values of Eq. (20) for the scalar dipole polarizabilities and the experimental values of Ref. Martel et al. (2015) for the leading-order spin polarizabilities are shown in Fig. 11 as a function of the lab photon energy $E_{\gamma}$ and the lab scattering angle $\theta_{\text{lab}}$, in comparison with the experimental data of the FULL data set. The grey bands correspond to the 1-$\sigma$ error range, computed in the bootstrap framework. For each values of $E_{\gamma}$ and $\theta_{\text{lab}}$, we calculate the differential cross section $d\sigma/d\Omega$ as function of the best values of $\alpha_{E1}$ and $\beta_{M1}$ obtained at every bootstrap cycle. We then have $N=10000$ values for $d\sigma/d\Omega$, from which we can reconstruct its probability distribution and the $68\%$ confidence level range.

In Fig. 12, we also show the differential cross section at forward angle $\theta_{\mathrm{lab}}=0^{o}$ from our analysis (red band) in comparison with the results obtained with the empirical forward RCS amplitudes of Refs. Gryniuk et al. (2016, 2015) (blue band). The last ones are evaluated from dispersive sum rules, using as input the total photoabsorption cross sections fitted to the available experimental data. In particular, we used the empirical amplitudes from the fit I of Refs. Gryniuk et al. (2016, 2015), that are tabulated from $E_{\gamma}=50$ MeV and correspond to $\alpha_{E1}+\beta_{M1}=14.29\pm 0.27$ and $\gamma_{0}=-0.929\pm 0.105$. We observe a remarkable agreement between the two analysis.

VII Conclusions

We performed a fit of the electric $\alpha_{E1}$ and magnetic $\beta_{M1}$ polarizabilities to the proton RCS unpolarized cross section data below pion-production threshold, using subtracted fixed-$t$ DRs and a bootstrap-based statistical analysis. Within the subtracted DR formalism, all the leading-order static polarizabilities enter as subtraction constants to be fitted to the data. However, due to the limited statistic of the RCS data, a simultaneous fit of all of them is not achievable at the moment. We then have restricted ourselves to fit the sets $\{\alpha_{E1},\beta_{M1}\}$ or $\{\alpha_{E1},\beta_{M1},\gamma_{\pi}\}$, which mainly affect the unpolarized RCS cross section below pion-production threshold. The remaining spin polarizabilities have been fixed to the available experimental information Martel et al. (2015); Ahrens et al. (2001); Dutz et al. (2003); Schumacher (2005). Furthermore, we consider different fit conditions, switching on/off the systematic errors and with/without the constraint of the Baldin’s sum rule for the polarizability sum $\alpha_{E1}+\beta_{M1}$.

We summarized the main features of the parametric-bootstrap method, in particular the advantages of taking into account both the effect of the systematic errors of the experimental data and the propagation of the statistical errors of the polarizability values not treated as free parameters in the fit procedure.

We showed that the inclusion of the sources of systematic errors in the data analysis changes significantly the expected theoretical probability distribution of the final $\hat{\chi}^{2}/d.o.f.$ variable and we were able to give realistic p-values for every fitting condition. We also presented a critical discussion of the data set consistency. We showed some simple but meaningful tests, which led us to conclude that there are no strong motivations for the exclusion of any data point from the global RCS data set below pion-production threshold. Even if we observed sizeable deviations between our fit model and two data subsets, there is not a clearly identified source of possible experimental problems for these data. Therefore, instead of excluding them from the fit, we discussed the possibility to handle them with a suitable rescaling factor of the statistical error bar. Also the bootstrap fitting technique showed to be useful in these cases to check if the inclusion of systematic errors in the fitting procedure improves the significance of the fit results.

The bootstrap fit using fixed-$t$ subtracted DRs and the global RCS data set below pion-production threshold yields $\alpha_{E1}=(12.03^{+0.48}_{-0.54})\times 10^{-4}\text{fm}^{3}$ and $\beta_{M1}=(1.77^{+0.52}_{-0.54})\times 10^{-4}\text{fm}^{3}$, with p-value $=12\%$. The results are in agreement with previous analysis obtained with different variants of DRs and the traditional $\chi^{2}$ fitting procedure. They differ from the extractions using the $\chi$PT frameworks, even if they are compatible within the $2\sigma$ range. This discrepancy can be traced back to the different data sets used in the analyses and, partially, also to the different theoretical estimates of the higher-order contributions beyond the scalar dipole polarizabilities to the RCS cross section.

Future measurements planned by the A2 collaboration at MAMI below pion-production threshold Sokhoyan et al. (2017); Downie and et al. (2016) hold the promise to improve the accuracy and the statistic of the available data set and will help to extract with better precision the values of the proton scalar dipole polarizabilities.

VIII Acknowledgments

We are grateful to V. Bertone and A. Rotondi for a careful reading of the manuscript and useful comments. We thank D. Phillips for stimulating discussions and useful suggestions on the fitting procedure, and H. Griesshammer, J. Mc Govern and V. Lensky for the help on the correct representation of the results of $\chi$PT.

Appendix A Data sets

In Table 3, we list all the available data sets for RCS in the energy range below pion production threshold ($\sim 150$ MeV in lab frame). For the sets Oxley (1958); Hyman et al. (1959); Goldansky et al. (1960) and Pugh et al. (1957), we use the Baranov data-selection Baranov et al. (2001). Furthermore, as done also in Ref. McGovern et al. (2013); Griesshammer et al. (2012), we discard the data from Table I in the Hallin paper Hallin et al. (1993), because it is not clear if they are really independent from the data given in Table II of the same work. The data sets used in our analysis are:

• •

FULL, which includes all the available data sets below pion-production threshold listed in Table 3, for a total of 150 data points.

• •

SELECTED, which is based on the data selection proposed in Ref. McGovern et al. (2013); Griesshammer et al. (2012), corresponding to the FULL data set except for the data from Ref. Oxley (1958); Bernardini et al. (1960); Baranov et al. (1975), a single point ($\theta_{\text{lab}}=133^{\circ}$, $E_{\gamma}=108$ MeV) from Ref. Olmos de Leon et al. (2001) and a single point ($\theta_{\text{lab}}=135^{\circ}$, $E_{\gamma}=44$ MeV) from Ref. Federspiel et al. (1991), for a total of 137 data points.

• •

TAPS, which is the most comprehensive available subset with 55 data points below pion-production threshold Olmos de Leon et al. (2001).

The sets 6 and 7 from Ref. Baranov et al. (1974, 1975) are from the same experimental measurements, but they differ for the values of the systematic errors. The same for the sets 11 and 12 from Ref. MacGibbon et al. (1995).

Appendix B Correlation coefficients among fit parameters

In the bootstrap framework, the correlation coefficients $\rho$ among the fit parameters are obtained from the reconstructed probability distribution in the parameters space. In Table 4, we list these coefficients for all the different fitting conditions used in this work.

In the Baldin-constrained fits, we do not obtain $\rho_{\alpha_{E1}-\beta_{M1}}=-1$, due to the fact that $\alpha_{E1}+\beta_{M1}$ is not fixed to its central value, but is sampled within its uncertainty with a Gaussian distribution, as explained in Sec. IV.1. This behavior was already observed in the extraction of the scalar dipole dynamical polarizabilities in Ref. Pasquini et al. (2018). We also note a large and negative (positive) correlation between $\gamma_{\pi}$ and $\beta_{M1}$ ($\alpha_{E1}$). This behavior is mainly a consequence of low sensitivity of the existing data to the $\gamma_{\pi}$ polarizability.