Constraints on Lorentz Invariance Violation from Optical Polarimetry of Astrophysical Objects

Theories of quantum gravity suggest that Lorentz invariance, the fundamental symmetry of the Theory of Relativity, may be broken at the Planck energy scale. While any deviation from conventional Physics must be minuscule in particular at attainable energies, this hypothesis motivates ever more sensitive tests of Lorentz symmetry. In the photon sector, astrophysical observations, in particular polarization measurements, are a very powerful tool because tiny deviations from Lorentz invariance will accumulate as photons propagate over cosmological distances. The Standard-Model Extension (SME) provides a theoretical framework in the form of an effective field theory that describes low-energy effects due to a more fundamental quantum gravity theory by adding additional terms to the Standard Model Lagrangian. These terms can be ordered by the mass dimension d of the corresponding operator and lead to a wavelength, polarization, and direction dependent phase velocity of light. Lorentz invariance violation leads to an energy-dependent change of the Stokes vector as photons propagate, which manifests itself as a rotation of the polarization angle in measurements of linear polarization. In this paper, we analyze optical polarization measurements from 63 Active Galactic Nuclei (AGN) and Gamma-ray Bursts (GRBs) to search for Lorentz violating signals. We use both spectropolarimetric measurements, which directly constrain the change of linear polarization angle, as well as broadband spectrally integrated measurements. In the latter, Lorentz invariance violation manifests itself by reducing the observed net polarization fraction. Any observation of non-vanishing linear polarization thus leads to constraints on the magnitude of Lorentz violating effects. We derive the first set limits on each of the 10 individual birefringent coefficients of the minimal SME with d = 4 , with 95% confidence limits on the order of 10−34 on the dimensionless coefficients.


Introduction
Lorentz invariance is the fundamental symmetry of Einstein's theory of Special Relativity. It has been established by many classic experiments, such as Michelson-Morley, Kennedy-Thorndike, and Ives-Stilwell [1][2][3], and tested to great precision by modern experiments [4]. Theories of quantum gravity suggest that there may be minute deviations from Lorentz symmetry, which motivates ever more sensitive tests [5][6][7][8][9][10].
Violations of Lorentz symmetry can lead to an energy-dependent vacuum photon dispersion relation, birefringence as well as anisotropy of the vacuum [11]. All of these effects can be tested with astrophysical observations, which are particularly sensitive because minuscule effects accumulate as photons propagate over very large distances resulting in measurable effects [12]. Vacuum birefringence leads to a wavelength-dependent change of the Stokes parameters, generally resulting in a rotation of the linear polarization angle. Astrophysical tests of Lorentz symmetry include time of flight with the expansion in spin-weighted spherical harmonics s Y jm and mass dimension d: wheren is the direction towards the origin of the photon. For odd d, there are (d − 1) 2 coefficients k (B)jm each. In this paper, we restrict our attention to the birefringent coefficients with d = 4 in the minimal SME, i.e., with a total of 10 coefficients k (4) (E)2m and k (4) (B)2m , which comprise a total of 10 real components since The polarization of an electromagnetic wave is completely described by the four Stokes parameters: intensity I; linear polarization Q and U, where U describes linear polarization at an angle of 45 • relative to Q; and circular polarization V. General elliptical polarization is described by the Stokes vector s = (Q, U, V) T . The polarization of photons with energy E will change as they propagate through a birefringent vacuum: with the birefringence axis ς = (ς 1 , ς 2 , ς 3 ) T . In the CPT-odd case, this axis is aligned with the V axis and as a result, linearly polarized light remains linearly polarized, but the polarization position angle will rotate as light propagates.
In the CPT-even case, the birefringence axis lies in the Q − U plane. Consequently, the Stokes vector will generally rotate out of this plane and linearly polarized light will become elliptically polarized during propagation. However, light with s parallel to ς will remain unaffected. The eigenmode of propagation is described by the polarization angle (following the convention used in Ref. [12]) The observed polarization of light emitted by a source at redshift z can conveniently be calculated in a spin-weighted Stokes basis s = (s (+2) , s (0) , s (−2) ) T with s (0) = V and s (±2) = Q ∓ i U and ς = (ς + , ς 3 , ς − ) T . Then, the observed Stokes vector s is related to the blueshifted Stokes vector s z by [12] s = M z · s z (9) with the Müller matrix M z =    cos 2 Φ z −i sin(2Φ z )e −iξ sin 2 Φ z e −2iξ − i 2 sin(2Φ z )e iξ cos(2Φ z ) i 2 sin(2Φ z )e −iξ sin 2 Φ z e 2iξ i sin(2Φ z )e iξ cos 2 Φ z    (10) and, at d = 4, For convenience, we define the following abbreviations: so that Φ z = Eϑ z (n) (16) and Most astrophysically relevant emission mechanisms are not expected to produce any significant circular polarization. Hence, assuming 100 % linearly polarized light at the source with the polarization angle ψ z , the observer Stokes parameters are Since the data used in this analysis do not contain any information about circular polarization, we do not consider V. The change in Stokes parameters is: The above expressions can be further simplified by realizing that the reference direction for the polarization angle can be chosen freely. Rotating the coordinate system such that ξ = 0, we express all position angles as ψ = ψ − ξ/2 and find All primed quantities are expressed in this rotated frame, while all polarization angles without prime are given in a frame where a polarization angle of 0 corresponds to linear polarization in the north/south direction and 90 • to the east/west direction [25]. Quantities with subscript z refer to the polarization of the source at redshift z, quantities without a subscript to the observer polarization predicted by the SME, and quantities with subscript m, in the following, refer to measured quantities.

Astrophysical Polarization Measurements
Essentially, birefringence leads to an energy and ψ z -dependent rotation of the polarization angle and change in linear polarization fraction, as illustrated in Figures 1 and 2. By measuring the linear polarization of photons emitted by distant objects, strong constraints on birefringence can be obtained.
In the analysis presented here, we make use of two kinds of measurements: spectropolarimetric measurements, where polarization fraction and angle are measured as a function of photon energy, and spectrally integrated measurements, where the polarization fraction is measured by integrating over a broad bandwidth determined by a filter in the optical path. Both analyses are based on the results from the previous section, but proceed differently. The goal of this section is to develop statistical measures for each type of observation that allows us to quantify the compatibility of a set of SME coefficients with the observation. The results are then combined in Section 4 into a joint probability function that is used to derive confidence intervals for each individual SME coefficient.  the SME coefficients were set to 0, except for (k (4) (E)21 ) = 2 × 10 −32 . The depolarization is due to averaging over the bandwidth of the filter and the rotation of the polarization angle, as shown for example in Figure 1, as well as the change of linear into circular polarization described by the Müller matrix (Equation (10)). The combination of these two effects leads to the observed "ringing". As described in Section 2, the effect depends on the linear polarization angle ψ z at the source.

Spectropolarimetry
When measuring the polarization angle as a function of energy, ψ(E), we can directly compare the result to the position angle resulting from Equations (19) and (20). Here, we reduce the problem to comparing the change in polarization angle at a given wavelength as predicted by the SME given a set of coefficients to the observed change over an instrument band pass. We start with a linear fit of the measured polarization angle as function of energy, where ψ m (Ē) is the measured polarization angle at the median energyĒ = 2.26 eV of the fit range, and ρ m with the uncertainty σ ρ is the linear rate of change of polarization angle as a function of energy. An example fit is shown in Figure 3.   Table A3.
To compare the measured rate of change ρ m ± σ ρ with the prediction due to a given set of SME parameters, we first find the source polarization angle ψ z from the observed polarization angle ψ m (Ē) at the median energy of the detector band pass. From Equations (23) and (24), we have where Q m (Ē) = cos(2(ψ m (Ē) − ξ/2)) and U m (Ē) = sin(2(ψ m (Ē) − ξ/2)), so that Then, we linearize the predicted change in polarization angle with energy as given by the Stokes parameters in Equations (19) and (20), which is adequate for small changes over the bandwidth: An example of the rotation of the linear polarization angle as a function of one of the SME coefficients is shown in Figure 4. The observer polarization angle in general oscillates around the source polarization angle and the roots ofρ seen in the figure correspond to values of k (4) (E)20 for which this oscillation is extremal atĒ.
This allows us to quantify the compatibility of the measurement with a given set of SME coefficients. The probability to observe a change of polarization angle ρ < |ρ m | assuming a true change |ρ| given by Lorentz violation according to Equation (28) and given the uncertainty of the measurement, σ ρ , is:  (E)20 while keeping all other SME coefficients at 0. The source was assumed to be at redshift z = 1 with a codeclination of θ = 90 • . Results are shown for three different values of the polarization angle at the source, ψ z . For ψ z = 45 • , the Stokes vector rotates out of the Q − U plane, but the linear polarization angle does not change because Q z = 0 (see Equations (26) and (27)).
An example of this probability as a function of one of the SME coefficients is shown in Figure 5. The spikes are due to the roots ofρ, but their width decreases with increasing SME coefficients. A single observation cannot be used to place constraints on the SME coefficients, unless additional assumptions are made. However, combining multiple observations will lead to tight constraints roughly corresponding to the width of the central peak, because the spikes at larger values of the SME coefficients will not line up for different sources. k (E) 20 (4) /10 -32 Figure 5. Probability P(ρ < |ρ m | |ρ|, σ ρ ) according to Equation (29) as a function of the SME coefficient k 20 while keeping all other SME coefficients at 0. A measured change in polarization of ρ m = (0 ± 1) • /eV was assumed. Source distance and direction are the same as in Figure 4, and colors have the same meaning. The roots ofρ in Figure 4 lead to the spikes in the probability function seen here.
The interpretation of the probability P is as follows. If |ρ m | < |ρ|, a strong degree of cancellation of the LIV-induced change of position angle with a source intrinsic change in position angle must have occurred, resulting in a low probability and a strong constraint on the given value ofρ. On the other hand, if |ρ m | > |ρ|, there must be a strong source-intrinsic change of the polarization angle, irrespective ofρ. Hence, P will be large and only a very weak constraint is placed on the given value ofρ. This has deliberately been designed such that no claims of detection of Lorentz invariance violation will be made, because we do not want to make any assumptions about source-intrinsic energy-dependent changes of the polarization angle.
We applied this method to a large sample of publicly available spectropolarimetric measurements of AGN [26] covering observer frame wavelengths between 4000 Å and 7550 Å. From all data published in this archive until 30 March 2016, we selected all sources with redshift z > 0.6, which have at least one observation with a spectrally averaged polarization fraction >10 %. For each source, we chose the measurement that resulted in the largest average polarization fraction. We then fitted each polarization angle measurement as a function (25) of photon energy in the range 1.77 eV to 2.76 eV with a linear function centered at the median energy of this range in order to determine the change of polarization angle ρ m . The resulting dataset was further reduced by removing all sources with ρ m /σ ρ > 3. The final list of measurements is given in Appendix A (Table A3).

Polarimetry Integrated over a Broad Bandwidth
When integrating over the bandpass of a broadband polarimeter, vacuum birefringence will lead to a reduction of the observed polarization compared to the polarization at the source due to the rotation of the polarization angle. We derive the largest possibly observable linear polarization fraction for an instrument with energy-dependent detection efficiency T(E), assuming that the emitted light is 100 % linearly polarized with an energy-independent polarization angle ψ z . We then quantify the compatibility of this maximum possible polarization with measured polarization fractions and angles. Any energy dependence of ψ z will lead to an additional reduction of measured polarization, making this a conservative approach. While it is in principle possible that birefringence and source-intrinsic effects cancel, this is very unlikely, in particular when observing multiple astrophysical sources.
We start by computing the effective Stokes parameters of a measurement for a given set of SME coefficients and a 100 % polarized source. Additivity of Stokes parameters allows us to integrate them over the detector bandpass, and where we use the definitions in Equations (23) and (24) and introduce the instrument-dependent normalization constant and the instrument-dependent function These integrals must be computed numerically, since T(E) is typically measured for each individual instrument. The advantage of formulating the problem in this way, however, is that N solely depends on the instrument being used, and F (ϑ) can be tabulated for efficient evaluation. The integrals in Equations (32) and (33) were calculated in the range 1.2 eV to 2.8 eV. All source properties (distance and direction) and SME coefficients are combined into the single instrument-independent parameter ϑ.
In this analysis, we used data from various optical telescopes employing a variety of filters. The filter transmission curves used in this analysis are shown in Figure 6. Table 1 lists the resulting normalization constants N and Figure 7 shows the tabulated functions F (ϑ). With those definitions, the maximum observable polarization for a 100 % linearly polarized source is where we used Q 2 z + U 2 z = 1 in the last step. The corresponding observed polarization angle is where the sign in the denominator is chosen to match the sign of Q z .
HOWPol R-band Figure 6. Filter transmission as a function of photon energy for the instruments and filters used in this analysis [27][28][29][30][31]. The function F (ϑ) determines the reduction of the observable linear polarization. It shows only a minor dependence on the exact shape of the transmission curve (compare for example the FORS1 R-band and the FORS2 R Special filters), but a clear dependence on the energy range covered by the filter can be seen: for small values of ϑ, the V-band filters lead to larger values of F (ϑ) than the R-band filters, and hence stronger sensitivity to Lorentz invariance violation. It is easy to show that in the limit of large SME coefficients lim ϑ→∞ F (ϑ) = 1 (36) and it is obvious from Figure 7 that F (ϑ) oscillates around this value as ϑ increases. In these cases, i.e., for F (ϑ) = 1, one finds This result implies that only certain effective polarization angles Ψ can be observed in case ϑ is large. Hence, the observed polarization angle itself already places constraints on the SME coefficients.
Given a set of SME coefficients, source distance L z and directionn, we use Equations (34) and (35) to find the largest possible polarization fraction Π max given a measured polarization angle Ψ m . It is important to note that this does not imply the assumption of 100 % polarization at the source, but simply reflects the fact that we do not want to make any assumptions about the astrophysics of the source. In realistic models, the polarization of optical emission from blazars is expected to be at most on the order of 20-30 % [32,33]. Using such models as input, significantly tighter constraints would be possible. However, optical polarization of blazars is highly variable (see, e.g., [34]), and in case of GRB afterglows significantly higher degrees of polarization are possible [35]. To be conservative and to avoid systematic uncertainties, we decided against using detailed source models, and follow the approach used in previous work [19,20].
We numerically solve Equation (35) for U z by requiring Ψ = Ψ m , where Ψ m = Ψ m − ξ/2 is the measured polarization angle in the rotated frame. Figure 8 illustrates this for different values of F (ϑ) as function of the observed polarization angle Ψ m . The result is then used to calculate Π max from Equation (34), shown as a function of Ψ m and F (ϑ) in Figure 9. These results were also tabulated in the range 0.01 ≤ F (ϑ) ≤ 1 in steps of δF (ϑ) = 0.005 and δΨ m = 1 • for fast lookup at later stages of the analysis. Values of Π max for arbitrary F (ϑ) and Ψ m can be found from this table using bilinear interpolation. Figure 10 shows an example of the maximum theoretically possible net polarization for GRB 091208B as a function of one of the SME coefficients.   The probability to observe a polarization Π given a true polarizationΠ, can be written as [36,37]: where I 0 is the modified Bessel function of order zero, i 0 (x) = exp(−|x|)I 0 (x), and N is related to the statistical quality, e.g., the number of photons detected in a photon counting experiment. Use of the scaled modified Bessel function i 0 (x) is advantageous for numerical implementation. Expectation value and standard deviation of Π are then where I 1 is the modified Bessel function of order 1. For each polarization measurement Π m ± σ m , we determine N by numerically solvingσ Π = σ m for N assumingΠ = Π m . This allows us to calculate the cumulative probability by numeric integration of Equation (38): This result quantifies the probability that a given set of SME coefficients is compatible with the measurement. Typical examples of the integral as a function of the upper limit Π max are shown in Figure 11. In practice, we replace P(Π|Π, N) from Equation (38) with a simple Gaussian distribution with mean Π m and standard deviation σ m if Π m /σ m > 10 due to numerical issues when evaluating the Bessel functions for large values of N. The error in this case is P Gauss (x < 0|µ > 10σ) < 7.7 × 10 −24 . An example of this cumulative distribution as a function of one of the SME coefficients is shown in Figure 12.  We applied this method to a large sample of polarization measurements of AGN [38,39]. From this catalog, we selected 36 sources for which polarization was measured with at least 5σ significance. Furthermore, we included optical polarization measurements of eight GRBs [40][41][42][43][44][45][46][47][48][49][50][51]. This selection of measurements is identical to the data used in Ref. [20] with the exception that we were unable to use data from GRB 090102 because no polarization angle was published [52]. Details about all astrophysical sources and measurements used here are shown in Tables A1 and A2 in Appendix A.

Results
The goal of this work is to obtain constraints on the individual coefficients k where the P i (k (4) (E)2m , k(B)2m (4) ) are the probabilities that each measurement is compatible with the given set of SME coefficients calculated as described in the previous section (Equations (29) and (41)). The product goes over all measurements of both types, listed in Tables A1-A3 in Appendix A.
To calculate constraints on each SME coefficient, we evaluate the probability in Equation (42) and marginalize over all other coefficients. Using the Metropolis-Hastings Markov Chain Monte Carlo method [53] with a normal distribution with σ = 6 × 10 −35 for each coefficient as the proposal distribution g(δk (4) (E)2m , δk (4) (B)2m ), we sample the combined distribution. At each step t, we randomly choose a new set of parameters k (E,B)2m according to g and calculate the acceptance ratio α = P(k (B)2m,t ) according to Equation (42). The step is accepted if u ≤ α for a uniform random number u ∈ [0, 1], otherwise we set k In this way, we iteratively construct a set of coefficients k (4) (E,B)2m that can be shown to be distributed according to the probability density in Equation (42).
The choice of proposal distribution g leads to a step acceptance rate of about 13 %. Starting from all coefficients set to 0, we discarded the first 10,000 steps to reduce the dependence of the result on the initial set of coefficients, and then took 10 7 steps recording each set of selected SME coefficients.
The distribution of these sets of coefficients is proportional to P(k (4) (E)2m , k (4) (B)2m ). From the distribution of values of each individual coefficient, we then find the 5th and 95th percentile as lower and upper limits. The resulting distributions are shown in Figure A1 in Appendix B and the corresponding upper and lower limits on the SME coefficients are listed in Table 2. The constraints turned out to be symmetrical around 0 within at least two-digit precision.
We also derived constraints from each of the two methods individually as a cross check. When using only the spectrally integrated measurements described in Section 3.2, the upper limits on the coefficients are about a factor 100 larger than those in Table 2. Consistently, we found that the results obtained when only using the spectropolarimetric results from Section 3.1 differ from the final ones on the order of 1 %.
By combining multiple spectropolarimetric measurements, the probability spikes in Figure 5 are eliminated as their location depends on source distance, direction, and measured polarization angle. The "ringing" observed in the probability distribution derived from each spectrally integrated measurement seen in Figure 12 cancels to some degree when combining multiple spectrally integrated measurements. Furthermore, the spectropolarimetric measurements are significantly more constraining than the spectrally integrated results. Combined, these effects result in the relatively smooth probability distribution of each SME coefficient shown in Figure A1. Figure A2 in Appendix B shows the correlation between the different coefficients. There is significant correlation among most pairs of coefficients most likely due to the non-uniform sky coverage. Therefore, the best way to improve on these results is to improve the sky coverage particularly with spectropolarimetric measurements, which at this point were not available for most of the Southern sky. Table 2. Limits at the 95 % confidence level on all independent SME parameters k

Discussion
Using optical polarization measurements of 63 AGN and GRBs, we searched for signals of Lorentz invariance violation. We derived 95 % confidence level limits on each of the 10 coefficients with mass dimension d = 4 of the minimal SME on the order of 10 −34 , as listed in Table 2. The results summarized in Table 2 are the first constraints on the individual SME parameters in this sector. Furthermore, the results are based on highly significant polarization measurements with well quantified uncertainties. In comparison, the limits published in [19] only constrain combinations of SME parameters, and are based on measurements of the polarization of the X-ray emission from GRBs with large systematic uncertainties.
While there may not be a theory of quantum gravity currently under active consideration that predicts a form of Lorentz invariance violation described by these coefficients, there is no generally accepted quantum gravity theory, either. A broad search for possible effects is, thus, warranted as neither the form of potential Lorentz invariance violation, nor its degree is known a priori. Given the tight constraints on the SME coefficients of d = 4 presented here and in prior work [4], any significant violation of Lorentz symmetry at the Planck scale must be strongly suppressed at lower energies. If a hierarchy of scales exists, suppression with (M low /M Planck ) n is possible, where M low is the lower mass scale and n is some power. SUSY-breaking is one scenario that could provide such a hierarchy [54]. Given extensive direct searches for deviations from the Standard Model of particle physics at energies up to about 1 TeV [55], M low 1 TeV must be assumed. A model with M low = 10 TeV and n = 2 would lead to a suppression of the order of 10 −30 . In this scenario, Lorentz invariance violating effects of order unity at the Planck scale due to operators of mass dimension d = 4 can be ruled out by the constraints presented here. For a similar model to be viable and to result in significant Lorentz invariance violating effects at the Planck scale, one must assume that either n > 2 or M low 10 TeV. While optical polarization measurements are a powerful tool to constrain coefficients of d = 4 and d = 5 [20], at higher mass dimension polarization measurements at higher energy will be significantly stronger due to the E d−3 dependence of the polarization signature. The Imaging X-ray Polarimetry Explorer (IXPE, [56]) to be launched in early 2021 will provide highly significant X-ray polarization measurements of a large number of astrophysical objects and one can expect that a similar analysis with those data will significantly improve on the results presented here. Future Compton gamma-ray telescopes, such as the proposed All-sky Medium-Energy Gamma-ray Observatory (AMEGO [57]) or the Compton Spectrometer and Imager (COSI [58]) will be sensitive to gamma-ray polarization up to MeV energies. These future instruments will significantly enhance our capability to search for Lorentz invariance violation in the photon sector.

Conflicts of Interest:
The author declares no conflict of interest.

Appendix A. Astrophysical Sources
Tables A1-A3 list all astrophysical sources used in this analysis, where the spectropolarimetric results are given in Table A3.   [51] † † No uncertainty of the position angle was specified in the paper. We are able to use this measurement because we always find the position angle the source to exactly match the observed position angle, irrespective of its uncertainty. Table A3. Sources selected from the Steward Observatory spectropolarimetric monitoring project [26]. The second column lists the highest observed polarization fraction during cycles 1-7 of the monitoring program. Coordinates have been obtained from the SIMBAD database [59]. Individual references are given for the red shifts. The last two columns give the position angle at the median energy of the linear fit and the differential change in position angle.   Figure A1 shows the probability distributions for all 10 SME coefficients with d = 4 derived in this analysis. Figure A2 shows the correlations between coefficients. Re(k Figure A2. Correlations between SME coefficients found this analysis.