Probing for Lorentz Invariance Violation in Pantheon Plus Dominated Cosmology

The Hubble tension in cosmology is not showing signs of alleviation and thus, it is important to look for alternative approaches to it. One such example would be the eventual detection of a time delay between simultaneously emitted high-energy and low-energy photons in gamma-ray bursts (GRB). This would signal a possible Lorentz Invariance Violation (LIV) and in the case of non-zero quantum gravity time delay, it can be used to study cosmology as well. In this work, we use various astrophysical datasets (BAO, Pantheon Plus and the CMB distance priors), combined with two GRB time delay datasets with their respective models for the {\em {intrinsic time delay}}. Since the intrinsic time delay is considered the largest source of uncertainty in such studies, finding a better model is important. Our results yield as quantum gravity energy bound $E_{QG}\ge 10^{17}$ GeV and $E_{QG}\ge 10^{18}$ GeV respectively. The difference between standard approximation (constant intrinsic lag) and the extended (non-constant) approximations is minimal in most cases we coincide. However, the biggest effect on the results comes from the prior on the parameter $\frac{c}{H_0 r_d}$, emphasizing once again that at current precision, cosmological datasets are the dominant factor in determining the cosmology. We estimate the energies at which cosmology gets significantly affected by the time delay dataset

The search for Lorentz Invariance Violation (LIV) through astrophysical probes has a long history [1,.Some quantum gravity theories predict violations of relativistic symmetries through messenger dispersion (photons, neutrinos, gravitational waves), the detection of which might offer crucial clues for unified theories [35].
There are two possible ways to look for LIV-either locally by dedicated experiments [49], that are so far out of our reach, or alternatively, through cosmological probes.The reason for this is that LIV effects are supposed to be amplified by the distance and also by the energy of the emission.Because of this, astrophysical probes such as gammaray bursts (GRB) are very well suited for such studies.Gamma-ray bursts possess two important qualities for such studies-they can be seen at extreme distances (z max ∼ 9.4 [50]) and at extreme energies (E max ∼ 10 55 erg [51]).Additionally, GRB's emissions have been observed in a very wide energy band, spanning from keV to TeV (for example, GRB 221009A with emission >10 TeV [52,53]).LIV effects are usually measured by the bound of the quantum energy E QG above which they could be observed.A plot of some measurements along with the references can be found in Figure 1.The color legend of the figure refers to the different approximations for the intrinsic time delay used to obtain the points -the blue square uses the standard approximation [22,54], the green color uses the energy fit [55][56][57][58][59], the red one uses the fireball model [60], the brown one uses a variable luminosity [61], and the black circle uses the SME framework [62].The most stringent bounds come from the TeV emissions of GRB 221009A (18 TeV [52]) and GRB 190114C (0.2 TeV [47]).
While the eventual LIV effect would be very small, if it is different from zero, it can also contribute to cosmology studies, providing new datasets independent of the luminosity measurements.Such datasets critically depend on the goodness of the GRB model (affecting the intrinsic time delay) and on the understanding of the propagational systematics of the messenger (affecting the other components of the time delay).Yet, they could help us look at cosmological tensions from another angle.The question of the use of GRBs in cosmology has been studied extensively, for example, using them as standardized candles [65][66][67][68] (and reference therein), through the cosmographic approach [69] and trying to reduce the Hubble tension [70], in combination with BAO, supernovae and quasars and a combined cosmological parameter [71].
In a previous paper [64], we started our investigation on how cosmology affects such constraints.LIV can be constrained either from the different energy bands of a single event or from averaging over multiple events, or both.However, cosmology needs to be taken into account in the estimations.It is particularly important when averaging over multiple GRBs due to the different redshifts employed.In [64], we obtained that the effect of adding cosmology may be significant for certain models and datasets.The biggest unknown in such a study is the intrinsic time delay-i.e., the possibility that the high and low-energy photons were not emitted at the same time.Since we do not have a good enough model of the GRB progenitor to predict the intrinsic time delay, we need to approximate it with a toy model.In [64], we used the standard approximation that assumes a constant (over the energy or the luminosity) intrinsic time delay, common for all the GRBs.
Here, we change that approximation, with two of the most popular other approximations.The first one is the energy-dependent intrinsic time delay introduced in [55,72], which yields as bound for the quantum gravity energy scale E QG,1 ≥ 0.3 × 10 15 GeV from GRB 160625B and 23 more GRBs [63] and E 1 QG ≥ 10 16 GeV [57].The other is the luminosity-dependent approximation for which the previously published result is E 1 QG ≥ 10 15 GeV [61].
In this paper, we use these two extended approximations for the intrinsic time delay and we apply them to two available GRB time delay (TD) datasets, to which we add several robust cosmological datasets.These are the angular baryonic acoustic oscillations (transversal BAO), supernovae type IA (the Pantheon Plus dataset), and the CMB distance prior.To the standard ΛCDM model, we add two dark energy (DE) models and a spatial curvature model (Ω K CDM).We show that at the current level of precision for the LIV parameters, the cosmology is more affected by the priors on the c/H 0 r d parameter than by the intrinsic time delay approximation.We see that indeed the extended approximation has an effect on our measurement of the LIV parameter α, but in most cases it is small.The other LIV parameters remain largely unconstrained.
The paper is organized as follows: in Section 2, we discuss the Theoretical background.In Section 3, we describe our Methods.In Section 4, we elaborate on our Datasets.In Section 5, we discuss the numerical results, and in Section 6, one can find a Discussion of the obtained results.

Theoretical Setup 2.1. Definition of the LIV Time Delay
We consider a modified dispersion relation for photons [19,62] of the type: with E QG -the energy of the QG scale, c is the speed of light, p and E-the momentum and energy of photons and the s ± = ±1 for subluminal or superluminal propagation [62].
Here, we work only with n = 1 and subluminal LIV.Following [23,73,74], the time delay between a low-energy (E low ) and a high-energy (E high ) photon with ∆E = E high − E low in the subluminal case is where H(z) = H 0 E(z) is the Hubble parameter, H 0 -the Hubble constant at z = 0 and E(z)-the equation of state of the Universe.

Different LIV Approximations
The observed time delay is a combination of the quantum-gravity time delay, the intrinsic time delay and different propagational effects [75][76][77].While the propagational effects can be ignored for high-energy photons, to put bounds on the quantum gravity time delay, we need to make an approximation for the intrinsic time delay.Here, we will compare the standard approximation with two of the most popular other approximations.
In this case, for the GRB time delay, one obtains: where a LIV ≡ ∆E/(H 0 E QG ), and If a LIV = 0, there is no LIV, while if a LIV ̸ = 0, there is LIV on energy scales above E QG .The connection between the QG effect and cosmology comes from the term E(z).
In our models, we use the form As in our previous paper [64], we denote a LIV = α and b is replaced with β, to avoid mixing this parameter with b = c/(H 0 r d ).
2. The energy-dependent approximation This approximation has been used in [55,57,58,63] and also in [56] with (α → −α).This popular formula assumes that in the source frame, the intrinsic positive time lag increases with the energy E as a power law: Here, E 0 = 11.34 keV is the median value of the fixed lowest energy band (10-12 keV).The free parameters are τ and α.While this fit is usually associated with studying a single GRB in multiple energy bands, [58] has also applied it to multiple GRBs.Note, from here on, we rename the parameters as follows: α → γ, τ → β.This is in order to differentiate them from the LIV α and the cosmological b.

Intrinsic GRB lag-luminosity relation approximation
This formula has been used in [61] to describe the intrinsic time delay of long GRBs.
Here, L * is arbitrary normalization, and γ and β long are the free parameters.For short GRBs, one takes γ = 0, which recovers the standard approximation.

Cosmology
The equation of state of the Universe for a flat Friedmann-Lemaître-Robertson-Walker metric with the scale parameter a = 1/(1 + z) is: where Ω DE (z) → Ω Λ for ΛCDM.The expansion of the universe is governed by E(z) = H(z)/H 0 , where H(z) := ȧ/a is the Hubble parameter at redshift z.Ω m , Ω K and Ω DE are the fractional matter, spatial and dark energy densities at redshift z = 0. We will consider a few dynamical dark energy models: number of different DE models: Chevallier-Polarski-Linder (CPL, [79][80][81]) and Barboza-Alcaniz (BA [82,83]).The relevant equations for Ω DE (z) for Equation ( 8) can be found in Table 1.For w 0 = −1, w a = 0 one recovers ΛCDM.
Table 1.The DE models we use in this work.The references can be found in the text.
To infer the cosmological parameters, we need to define the angular diameter distance D A : where sinn(x) ≡ sin(x), x, sinh(x) for Ω K < 0, Ω K = 0, Ω K > 0, respectively.It connects to the transversal BAO measurements through the angular scale measurement θ BAO (z): The Pantheon Plus datasets measure the distance modulus µ(z), which is related to the luminosity distance (d where M B is the absolute magnitude.Finally, we include the CMB through the distance priors datapoints provided by [84]: , where l A is the acoustic scale coming from the CMB temperature power spectrum in the transverse direction and R is the "shift parameter" obtained from the CMB temperature spectrum along the line-of-sight direction [85].r s (z * ) is the co-moving sound horizon at redshift photon decoupling z * ( z * ≈ 1089 [86]).

Methods
The idea of the method is to avoid setting priors on H 0 and r d by considering only the quantity b = c/(H 0 r d ).The method has been outlined in [64].For this reason, we just list the different likelihoods here.For BAO we have Here, ⃗ v obs is a vector of the observed points, for BAO corresponding to θ i BAO,obs , ⃗ v model is the theoretical prediction of the model (for BAO -θ BAO,theor (z i )) and σ is the error of each measurement.
For the SN dataset, we marginalize over H 0 and M B .The integrated χ 2 in this case is ( [87][88][89][90]): for where µ i is the observational distance modulus, σ i is its error, and the d L (z) is the luminosity distance, ∆µ = µ i − 5 log 10 [d L (z i )), E is the unit matrix, and C −1 cov is the inverse covariance matrix of the Pantheon Plus dataset as given by [91,92].
Finally, we define the time delay likelihood as Equation ( 12), but here, the quantity we consider is the theoretical time delay ( ⃗ ∆t model ), as defined in Equation ( 5) and its observational value (( ⃗ ∆t obs ), which is provided by the TD dataset.
The final χ 2 is

Datasets
In this work, we use the so-called transversal BAO dataset published by [93], for which the authors claim to be cleaned up from the dependence on the cosmological model and uncorrelated.The CMB distant prior is given by [84].The SN data comes from the Pantheon Plus dataset.It consists of 1701 light curves of 1550 spectroscopically confirmed Type Ia supernovae and their covariances, from which distance modulus measurements have been obtained [91,92,94].
To study the time delays, we use two different time delays (TD) datasets-TD1 provided by [61] and TD2 [59].TD1 uses a combined sample of 49 long and short GRBs observed by Swift, dating between 2005 to 2013.In this dataset [61,95], the time lags have been extracted through a discrete cross-correlation function (CCF) analysis between characteristic rest-frame energy bands of 100-150 keV and 200-250 keV.The redshift for TD1 is z ∈ [0.35, 5.47].TD2 [59] uses 46 short GRBs with measured redshifts at Fixed Energy Bands (15-70 keV and 120-250 keV) gathered between 2004 and 2020 by Swift/BAT or Fermi/GBM.The two datasets have only six common GRBs, which is under 15% of their total number (49 vs. 46 events), which makes them effectively uncorrelated and independent.Because we want to emphasize the effect of cosmology, we prefer to average over multiple GRBs rather than to use measurements from a single GRB [22].
To run the inference, we use a nested sampler, provided by the open-source package Polychord [96] and the package GetDist package [97] for the plots.
We use uniform priors for all quantities:

Results
In [64], we started investigating the effect of cosmology on LIV bounds by studying in the standard approximation (Equation (3)) two databases: the most famous and robust time delay dataset published by Ellis et al. [22] and the one that we currently refer to as TD1.To avoid repetition.In this section, we present only the extended approximations 1 .For completeness, the results for the standard approximation for TD1 and TD2 are presented in Appendix C.
Since our previous work [64] demonstrated a deviation from the expected values for b and Ω m , here, we want to investigate further this observation.To this end, we try two hypotheses for the two datasets-a uniform prior on b and a Gaussian prior focused on the expected value from Planck results, b = 30 ± 1.For the standard approximation, this leads to a very clear higher value of Ω m but a minimal difference for the DE parameters w 0 and w a and Ω K .
The results for the extended models are summarized below. -

Matter density
In Figure 2, one can see the plot of b vs. Ω m for the extended LIV models for the four different cosmological models we consider (and the two priors on b).For both datasets, we see two groups of posteriors-one grouped around the top right corner, which corresponds to the Gaussian priors on b, and one that spans across the whole interval for the uniform prior.In both cases, the mean value and the 95% CL of b are below the expected value of b ∼ 30 for the three models.
For TD1, setting the Gaussian prior for b in this case significantly constrains Ω m and puts it in the limits of [0.29,0.34]at 95% CL.In this case, Ω k CDM gives the highest mean Ω m and CPL gives the lowest one, with ΛCDM in between.For the uniform prior, the posterior for Ω m is much less constrained, with the lowest value coming from Ω K CDM and the posteriors for ΛCDM and CPL largely coinciding.
For the other dataset, TD2, the situation repeats to a great extent, with a little bit higher matter density Ω m ∈ [0.3, 0.34] for the Gaussian prior and b ∈ [27.8, 29, 2] and the lowest value for Ω m is for ΛCDM and the highest for Ω K CDM.For the uniform prior, we have again the lowest value for Ω m and b coming from Ω K CDM and ΛCDM and CPL largely coinciding. -

Spatial curvature
The results for Ω K CDM are presented in Figure 3.We see that our results are not consistent with a flat universe at the 95% CL even for a uniform prior on b.For TD1, we get a bit less constrained contours for Ω K , but more constrained than the flat cases.A rather interesting feature is the huge difference in Ω m between the uniform and the Gaussian cases for both datasets.Our prior on Ω K is relatively small (Ω K ∈ [−0.1, 0.1]), however, due to the large number of parameters of the model. -

DE parameters
The two DE models we consider are shown in Figure 4. Their DE parameters w 0 and w a are largely insensitive to the dataset and their mean values for w 0 are lower than expected for ΛCDM, with w 0 ∼ −0.8.This is a different value from the one obtained [64], which was consistent with ΛCDM at 68% CL.The result in the study [71] without TD on the other hand, shows values about w 0 ∼ −1.1.In the two latter works, we use the Pantheon dataset, while here, we use the Pantheon Plus.To clarify if the difference is due to the new SN dataset or to the new TD datasets, we ran the same experiment with Pantheon instead of Pantheon Plus and it gave w 0 → −1 within 68% CL.As expected, the parameter w 0 is well constrained by the data while w a is not constrained at all.Note, the cosmological results for TD1 and TD2 for the extended models look largely the same.That could be due to similar treatment of the data or other reasons (like adding additional degrees of freedom in the extended models).The standard approximation posteriors presented in Appendix A show a larger dependence on the dataset.Also, the LIV contribution to the fit is quite small due to the smallness of α.As long as it is non-negligible, it still contributes to the cosmology in a minor way, as seen above in b and Ω m .The numerical experiment shows that for the TD1 dataset, the effect of the TD dataset on the cosmology becomes really pronounced at about α ∼ 10 −2 , while for TD2, it becomes pronounced at about α ∼ 10 −3 .This is just an order or two above our current constraints for α.This means that observing an event that would push the bound for E QG lower (or having a GRB intrinsic time delay model that would do so) would also bring the TD measurement on par with the current cosmological probes. -

LIV Parameters
Finally, we are going to discuss the LIV parameters, α, β and γ.They are shown in Figures 5 and 6.In Figure 5, we can see the values for α we obtained (where it should be noted that the lower bound on α is, of course, 0).As a whole, TD1 gives 10 times higher values for α than TD2 (which has been noticed also in our previous paper).We also see that the extended models do not have a conclusive improvement over the standard approximation, even though, for some of the models considered, they lead to significantly smaller α, meaning higher E QG .In Figure 6, we show the parameters β and γ.We note that they are mostly unconstrained in all of the models and while there are differences between the standard and the extended approximations, they are minor.The posteriors for the LIV parameters can be found in Appendix B. −0.5 0.0 0.5 To obtain the values for E QG , one needs to use the formula E QG = ∆E/(H 0 a LIV ).This requires inputting the specific energy band for each dataset (∆E), and choosing a value for H 0 .The error in this case will be Since α = 0 is the low limit of our priors, which gives an upper bound for E QG infinity, we can get only the lower bound of the quantum gravity energy.We take as σ α the 68% CL corresponding to 1σ deviation and we take the minimal and the maximal energies obtained for the different models.The results can be found in Table 2.In it, one needs to remember that the error comes from the incertitude in H 0 and the numerical estimation of α; thus, it is not supposed to be taken as a measurement, but as an estimation of the error we get from this study.

Discussion
We have studied the LIV bounds and the effect of cosmology based on two different datasets TD1 [61,95] and TD2 [59].To allow testing for new approximations for the intrinsic time delay, we have added two extended models-the intrinsic GRB lag-luminosity relation approximation for TD1 and the energy-dependent approximation for TD2 and we have compared them to the standard approximation ("constant" intrinsic time, i.e., not depending on the energy or the luminosity).To study the effect of the cosmological model, we have used additional cosmological datasets including transversal BAO, the Pantheon Plus dataset and the CMB distance priors.These are some of the most robust cosmological datasets; thus, they provide the best opportunity to study the joint effect of cosmology and LIV.We have considered ΛCDN, Ω K CDM, the CPL and the BA dark energy models.
From the results we obtain, we see that the strongest effect is due to the prior on the parameter c H 0 r d and not so much on the approximation for the intrinsic lag.This is because first the LIV effect is expected to be very small, and second, the intrinsic lag parameters are largely unbound from the current data.Instead, we see that the results for α depend a lot on the cosmological model.The lowest bounds of E QG are for Ω K CDM and the highest for ΛCDM G (TD1) and BA (TD2).Surprisingly, despite the new degrees of freedom introduced by the LIV parameters, the Ω K CDM model suggests a closed universe.In terms of LIV energy, we obtain as lowest bound for TD1 E QG > 0.8 × 10 17 GeV and for TD2: E QG > 3.5 × 10 18 GeV.In both cases, the error is significant, regardless of the small error of H 0 .TD1 tends to give 10 times larger values of α than TD2.The effect of the TD datasets on the cosmological parameters becomes noticeable if α is an order higher than the inferred one, meaning lower than the currently estimated E QG .
In conclusion, we see that for the moment, the time delay datasets are not precise enough to constrain the cosmological effects, while the cosmological models have a serious effect on the LIV constraints.For this situation to change, we need to improve the model of the GRB central engine and also to be able to better constrain the propagational effects, which are considered negligible at high energies; however, one needs to remember that not all measurements are made at very high energies, especially for the older GRB collections.Finally, a new and better approximation for the intrinsic time delay could benefit both GRB theoretical models and cosmological studies.

Appendix A. Posteriors of the Models
. Posteriors (95% CL and 68% CL) for the TD1 dataset (top panel) and the TD2 dataset (bottom panel).On the left, we have the SA case with its two parameters α and β, on the right, the EA case with its three parameters α, β and γ.As seen for the plots, only α can be bounded with these datasets.

Appendix C. Tables of the Results
Table A1.Mean values for the TD1 dataset for the standard approximation.

5 ,
−0.5] and w a ∈ [−0.5, 0.5].Since the distance prior is defined at the decoupling epoch (z * ) and the BAO-at drag epoch (z d ), we parametrize the difference between r s (z * ) and r s (z d ) as rat = r * /r d , where the prior for the ratio is rat ∈ [0.9, 1.1] .The LIV priors are α

Figure 2 .
Figure 2. Posteriors (95%CL and 68% CL) for the TD1 dataset (left) and the TD2 dataset (right) for the parameters b and Ω m .

Figure 5 .
Figure 5. Mean values and for the parameter α for the standard approximation (red) and the extended one (black) different models.To the (left) are the results for TD1 and to the (right)from TD2.

Figure 6 .
Figure 6.Mean values and errors (68%CL) for the parameters β and γ for the standard approximation (red) and the extended one (black) different models.To the (left) are the results for TD1 and to the (right)-from TD2.

Figure A2 .
Figure A2.Posteriors for the TD2 dataset for the standard approximation (left panel) and the extended approximation (right panel).On each panel, the uniform prior is on the bottom triangle plot, the Gaussian-on the top one.

Table 2 .
The minimal and maximal energy E QG for the standard approximation (SA) and the extended approximation (EA).

17 Gev E max,SA QG × 10 17 Gev E min,EA QG × 10 17 Gev E max,EA QG × 10 17 Gev
Posteriors for the TD1 dataset for the standard approximation (left panel) and the extended approximation (right panel).On each panel, the uniform prior is on the bottom triangle plot, the Gaussian-on the top one.

Table A2 .
Mean values for the TD1 dataset for the extended approximation.

Table A3 .
Mean values for the TD2 dataset for the standard .

Table A4 .
Mean values for the TD2 dataset for the extended approximation.

Table A5 .
The LIV parameters for TD1: to the left-the standard approximation, to the right-the extended one.

Table A6 .
The LIV parameters for TD2: to the left-the standard approximation, to the right-the extended one.