Halo uncertainties in electron recoil events at direct detection experiments

The dark matter direct detection rates are highly correlated with the phase space distribution of dark matter particles in our galactic neighbourhood. In this paper, we make a systematic study of the impact of astrophysical uncertainties on electron recoil events at the direct detection experiments with Xenon and semiconductor detectors. We find that within the standard halo model there can be up to $ \sim 50\%$ deviation from the fiducial choice in the exclusion bounds from these observational uncertainties. For non-standard halo models, we report a similar deviation from the fiducial standard halo model when fitted with recent cosmological $N$-body simulations while even larger deviations are obtained in case of the observational uncertainties.


Introduction
In the last few decades particulate dark matter (DM) has been probed by its possible scattering with the Standard Model (SM) particles [1][2][3][4][5]. The typical direct detection experiments measure the nuclear recoil of a target material through scattering of ambient DM wind on the surface of the earth [1,2,[6][7][8][9]. While nuclear target experiments are suitable to probe a non-relativistic DM mass at O(100) GeV, however for DM masses in sub-GeV range these looses its sensitivity. This is due to the fact that energy deposited by a sub-GeV non-relativistic DM remains below the threshold of these experiments.
An alternate and novel strategy to search for such light DM is through the DM-electron scattering [10][11][12][13][14]. For an atomic target (e.g Xenon) if DM scatters off the electron on the atomic shell then this may lead to the ionization of electrons. Whereas for a semiconductor target material (e.g. Si, Ge etc.) scattering of DM with electron may transfer an electron from valance band to conduction band. These ionization signals could provide a handle in the search for sub-GeV DM. The boundedness of electron in the target material makes these electron scattering events inelastic in nature. This essentially suggests that the incoming DM particles has to have a sufficient energy to excite these bound electrons. For a semiconductor material the typical energy gap between valance and conduction band is of the order 1 eV whereas for Xe targets the minimum binding energy of a shell is around 10 times larger than that.
Thus for a Xe target materials a relatively light DM can excite an electron if the DM is moving fast. These fast moving DM can only be found near the tail of the galactic DM distribution. Subsequently the DM-electron event rate would be suppressed for such Xe targets as compared to semiconductors targets, implying sensitivity of semiconductor detector below the MeV scale.
The DM-electron scattering rate can be divided into three parts viz. particle physics, atomic physics and astrophysics. The particle physics input depends on the particular model under consideration and determines the hard scattering cross section between the DM and the electron. In this paper we take a model independent approach to estimate the cross sections. The electron ionization form factor constitutes the atomic physics part and depends on the wave function of the scattered electron. For this we have used the result of QEdark [15]. In QEdark, the form factor for Xe targets has been calculated using Hartee-Fock method while for semiconductor materials density functional theory has been utilized.
The local distribution of the ambient DM constitutes the astrophysical part. The Standard Halo Model (SHM) with Maxwell-Boltzmann (MB) distribution truncated at galactic escape velocity is usually assumed for the distribution of DM in the galaxy. The dispersion of the MB distribution is determined by the Sun's circular velocity (v 0 ). The typical choice of these parameters are v 0 = 220 km/s and v esc = 544 km/s, following [16]. However there is still a considerable uncertainty in the measurement of these astrophysical inputs [17]. In this paper we make a systematic study of the impact of these astrophysical uncertainties in DM electron scattering rate.
Cosmological N -body simulations generate a patch of our local universe containing mostly DM particles and in some cases stars and gas to study and compare our local universe including the Milky Way (MW) galaxy and its halo to the present day observations. They usually include the effect of baryons utilizing hydrodynamic simulation. We will refer these collectively as cosmological simulations. These simulations indicate that the SHM may not give an accurate description of a Milky Way-like halo. Modifications of the SHM framework have been introduced to reconcile the astrophysical observations and cosmological simulations [18][19][20][21]. These include the King velocity distribution determines the cut off in the distribution through a self consistent manner [22]. The Double Power Law which can nicely explain the high velocity dependence of double power density profiles like that of NFW [23]. The Tsallis is a theoretical distribution based on Gibbs entropy motivated by Tsallis statistics [24]. The Mao et. al. [25] suggests another distribution showing a strong correlation of particle velocities to their position and characteristic radius of the simulated halo [25]. In this work we have only considered the impact of isotropic velocity distributions on electron scattering events, while leaving the possible impact of anisotropic distributions [26] for a future work.
The paper is organized as follows. In section 2 we briefly review the methodology to calculate rate of the DM-electron scattering for atomic and semiconductor target material. In section 3 we present the variations of DM-electron exclusion limit due to the uncertainty in the SHM parameters. Going beyond, in section 4 we explore the effect of non-standard velocity distribution. In section 5 we present a detailed comparison between the deviations obtained from cosmological simulations and recent astrophysical observations, with reference to the fiducial SHM. Finally, we conclude in section 6.

DM-electron scattering
In this section we will briefly review the scattering of DM particles with electrons that are bound inside the detector material. We will pair it down to the particle physics effects, the atomic physics effects and the astrophysical factors. Let us consider a DM particle (χ) of mass m χ , and initial velocity v scatters of an electron within the target material. Then in the non-relativistic limit the energy conservation of the system implies [10] where q is the momentum transfer by DM and ∆E e is the energy transferred to electron. Note that for a DM mass in sub-GeV scale, the relevant momentum transfer is small compared to the mass of the nucleus, therefore the nuclear recoil piece has been neglected in equation (1) 1 . Following equation (1), the minimum DM speed v min required to transfer an energy ∆E e for a fixed q is given by For an atomic target (e.g. Xenon) the differential DM-electron scattering rate is [11] dR ion where N T denotes the number of atoms in the target. ρ χ is the local DM density. DM-electron reduced mass is denoted by µ χe . σ e stands for DM-electron cross section for a particular momentum transfer q = αm e and F DM is the DM form factor. The ionization form factor and average inverse velocity are represented by f n,l ion and η respectively. Whereas, the differential rate for a semiconductor target (e.g. Ge) can be written as [13] dR crystal where N cell stands for the number of unit cell in a crystal targets [13]. The f crystal denotes the ionization form factor for the crystal. Other parameters have their usual meaning as discussed.
Note that both in equations (3) and (4) F DM (q) takes care of the momentum dependency in DMelectron interaction. This DM form factor and σ e comprise the main particle physics input in equations (3) and (4). Remaining agnostic about any particular model here we have considered three types of interactions between DM and electron. These interactions can be quantified by the DM form factor F DM . The three choices for the F DM are 1, αm e /q, (αm e /q) 2 . Note that F DM = 1 can be induced by an exchange of heavy mediator between DM and electron [27,28], F DM = αm e /q which could arise through electric dipole moment interaction [29], and F DM = (αm e /q) 2 which may be induced by a light mediator [13,28].
The ionization form factor defines the suppression of the event rate to ionize an electron from its bound state to a continuum state of momentum k = √ 2m e E e through q. Hence it depends solely on the target material. Throughout our numerical calculation we have utilized the form factor given in QEdark [15].
Other than these the usual astrophysical inputs are the local DM density ρ χ and average inverse DM speed The f ⊕ (v) in equation (5) is the DM velocity distribution in the detector rest frame. If we assume f (v) as the DM distribution in the galactic frame then the distribution at lab frame can be obtained by where v e = v 0 + v + v ⊕ and v 0 and v are the Sun's circular velocity at local standard rest and Sun's peculiar velocity respectively. The earth velocity in the Solar rest frame is represented by v ⊕ . Note that the variation in v ⊕ with time leads to the familiar annual modulation [12,13] in the DM direct detection rate, has been neglected here. The sun's peculiar motion v = (U , V , W ) = (11.1 ± 1.5, 12.2 ± 2, 7.3 ± 1) is adapted from [30]. We set the earth's rotational velocity following the reference [6].
The mean inverse speed is essentially regulated by the astrophysical parameters discussed above. Therefore any uncertainty in determination of these parameters will have a direct impact in the exclusion limit. In the rest of this paper we systematically study the impact of these astrophysical uncertainties including departure from the MB distribution on the exclusion bounds of direct detection experiments.
The differential rates given in equations (3) and (4) are with respect to the electron recoil energy E e . While the semiconductor detectors are only sensitive to electron and the Xe detector finally detects photo-electron at the PMT's. For the Xe target material this electronic energy is converted into a number of electrons and subsequently to photo-electrons using the prescription of [14]. Note that the Xenon10 experiment [31] which has an exposure of 15 kg-days, sets most stringent limits in most of the region in the parameter space of interest [14], therefore we would only consider the bound from Xenon10 in this paper. The exclusion limit for Xenon10 is obtained by making a conservative assumption that all the observed events arise from DM. Whereas for semiconductor targets the electron-hole yields are obtained following reference [13]. We have considered the one electron threshold to present the bound for semiconductor detectors. For semiconductor detectors, the upper limit on the relevant parameters are presented at the 90% CL (corresponds to 2.3 DM events), assuming an one electron threshold with an exposure of 1 kg-day.

Astrophysical uncertainties within the Standard Halo Model
The most simplified isotropic and isothermal DM distribution in the MW halo is usually described by a MB distribution, with a cut off at the escape velocity of the DM particles [32,33]. The distribution function has the following form where N denotes the normalization constant of the distribution, v 0 is the measure of its velocity dispersion and v esc sets the maximum allowed DM velocity of the distribution. Keeping ourselves within this SHM, in this section we will present the uncertainties in the determination of astrophysical parameters and their implications on the DM-electron scattering events 2 . A discussion about the main observational uncertainties in the SHM parameters are now in order: 1. Local DM density: The typical choice for the local density is 0.3 GeV cm −3 [44]. Recent estimations suggest that it may vary in the range (0.2 − 0.6) GeV cm −3 [17,26,[45][46][47][48][49][50][51][52][53][54]. While others predict slightly different values [50][51][52][53]55]. However note that the differential rate given in equations (3) and (4) scale linearly with the local DM density. Therefore for a change in ρ χ one would expect a proportional vertical shift in the exclusion limits for all the experiments.
Assuming the central value of the local DM density 0.4 GeV cm −3 , we find that there is a maximum 100% relative change due to the aforementioned variation of ρ χ . This change is independent of the DM mass and the target materials used for the detection. Hence we have fixed this to 0.4 GeV cm −3 without considering its variational implications on the exclusion bounds.
2. Circular velocity of the Sun: The local circular velocity of the Sun (v 0 ) with respect to the galactic center is usually assumed to be 220 km/s [9,56,57]. This would be considered as the fiducial choice of the parameter v 0 for rest of the paper. From the orbit of the GD-1 stellar stream, the reference [58] constrained v 0 in the range 221 ± 18 km/s. A similar range of v 0 , namely 225 ± 29 km/s is found to be in consonance with the kinematics of maser [59]. These estimates seem to have around 10% error in the measurement of v 0 . However a more precise assessment of v 0 can be done using the measurement of apparent proper motion of Sgr A * relative to a distant quasar [60,61]. This measurement fixes the total angular velocity of the Sun ((v 0 + V )/R ) in the range 30.24 ± 0.12 km s −1 kpc −1 . On the other hand recently GRAVITY collaboration has estimated the value of R with quite a high accuracy: 8.122±0.031 kpc [62]. Also note that the relevant component of the peculiar velocity V varies in the range 12.24 ± 2.47 [30]. Combining all these observations, the circular velocity of the Sun has been found to be 233 ± 3 km/s [26]. Other recent assessments [63][64][65][66] would also lead to similar result for v 0 . We will explore the impact of the deviations of this from the fiducial value on the exclusion limits from direct detection experiments in the electron scattering events.
As can be seen from the equation (7), v 0 is related to the standard deviation of the distribution. Thus any increment in v 0 would flatten the distribution. Therefore this would make more DM particles available to interact with electrons in the tail region. This will lead to a relatively stronger bound in DM-electron cross section. The effect will be reversed for a decrement in v 0 . Further a change in v 0 would also alter the Galilean boost.
3. Galactic escape velocity: The escape velocity of a massive body in a galaxy is defined by the velocity above which they will no longer remain bound to its gravitational potential. Measurements from the high velocity stars of the RAVE survey determines the v esc in the range 498 − 608 km/s [16] with the median 544 km/s. In the rest of the paper we would consider this as the fiducial choice for v esc . Based on the recent analysis from the velocities of 2850 halo stars from the Gaia velocity survey Data Release-2 [67], the local escape speed has been revised to 580 ± 63 km/s. However it has been argued that this result is sensitive to the prior      (9)) with the fiducial values. The light blue, green and yellow shaded regions express the uncertainties in astrophysical observations for F DM = 1, αm e /q and (αm e /q) 2 respectively. In the lower panel we show the corresponding fractional changes relative to the fiducial choice in the SHM (defined in equation (8)). The uncertainties associated with the recent astrophysical observations are utilized to form light blue shaded bands. Whereas the blue dashed curves denote the relative fractional changes in the Modified Maxwell-Boltzmann distribution for F DM = 1.
chosen for describing high velocity tail of the distribution function. With a prior estimation from simulations, and a more localized sample of 2300 high velocity counter-rotating stars, the escape speed has been obtained to be 528 +24 −25 km/s [68]. This is also in consonance with the previous results. We will be using the central value of the latter as the new central value for the escape velocity of the DM particles within the SHM.
To estimate the relative fractional change in the cross-section we will adopt equation (8) throughout this paper, where σ i e denotes the cross section corresponding to non-fiducial values of the distribution in consideration.
With a decrease in the escape velocity there will be less energetic particles in the halo capable of scattering, hence a larger σ e is required to produce the same number of events at a given experiment. Due to the exponential suppression of the MB distribution near the tail, increasing v esc has a smaller [69][70][71] 646 212.7 APOSTLE DMB [69][70][71] 646 224.1 ARTEMIS DMO [72,73] 521.6 161.4 ARTEMIS DMB [72,73] 521.6 184.3 Table 1: Best fit parameters for the SHM. DMO refers to the DM only simulation. DMB indicates the simulation which takes into account both DM and baryon smooth particle hydrodynamics. effect than decreasing it by the same amount. We find that for similar relative change in v 0 and v esc , the effect due to the change in v 0 is more pronounced. This can be ascribed to the fact that any change in v 0 causes an overall change in the shape of distribution, whereas the v esc determines the distribution's cut-off near the exponentially suppressed tail.
The inverse relation between m χ and v min implies that the required minimum DM velocity is rather close to v esc for light DM. This is because DM particles having lower mass and hence lower kinetic energy can generate the required recoil energy if their minimum velocity becomes closer to the escape velocity of the distribution. And the tail of the distribution is quite sensitive to the choice of astrophysical parameters. Thus the observed fractional change in the exclusion limit is significantly larger for light DM. This is evident from all figures shown in the paper.
In the upper panel of figure 1, the light blue, green and yellow shaded regions represent uncertainties associated with an updated measurements of v esc = 528 +24 −25 km/s and v 0 = 233 ± 6 km/s for the three choices of the DM form factor F DM = 1, αm e /q and (αm e /q) 2 respectively at 95% confidence level. In the lower panel we show the deviation from the fiducial SHM for F DM = 1, the light blue bands arise due the uncertainties related to the measurement of v esc and v 0 . The corresponding deviations for the other choices of F DM lie in a similar range, hence have not been shown in the plots to reduce clutter. Like figure 1, in rest of the paper the variations of a Xenon target are shown in the left panels, whereas similar variations for Germanium and Silicon target materials are depicted in the middle and right panels respectively. In conclusion our study indicates a variation in ∆ between 2% to 50%, due to the astrophysical uncertainties in the exclusion bounds within the SHM. This is over and above the modifications due to the ambiguity in the measurement of the local DM density mentioned earlier.

Cosmological simulations
An alternate way to estimate the value of v 0 and v esc is to fit a model of DM velocity distribution to cosmological simulations. Here we have utilized the results of the APOSTLE [69][70][71] and ARTEMIS [72,73] simulations 3 to find the best-fit values of the velocity distribution parameters. For APOSTLE, two MW-like haloes identified as A1 and A2 in [73] are chosen based on the selection criterion: total halo mass (M 200 ) in the range 5 × 10 11 < M 200 /M < 2 × 10 13 , having rotation curves similar to the observed MW like rotation curves, stellar mass (M ) within the observed stellar mass of the MW i.e in the range 4.5 × 10 10 < M /M < 8.3 × 10 10 and having a substantial stellar disc [71,74]. In this work we have used the A1 halo for all our fits. For the ARTEMIS simulation [72], 42 MW like haloes are  chosen by considering galaxies in the mass range 8 × 10 11 < M 200 /M < 2 × 10 12 . We have utilized the median DM distribution of the 42 MW-like galaxies reported in [73].
Note that these two simulations capture the impact of baryons through hydrodynamic corrections and have relatively high force resolution. This may be contrasted with the older DM only simulations like GHALO [75] or Via Lactea [76] which have large deviations from the SHM. Therefore we do not consider them further in this study. Given these DM distribution of A1 for APOSTLE and the median DM distribution for ARTEMIS, the best fit values of the relevant DM distribution parameters are summarized in table 1. Throughout the paper, we have utilized an inbuilt python function scipy [77] to fit the simulation data with a given distribution. We find that our result matches relatively well with [23,71]. Note that throughout the work DMO implies simulations containing only DM for both APOSTLE and ARTEMIS. Both simulations also have a hydrodynamic version containing both DM and baryons. Simulations that include both DM and baryon effects hereafter will be referred to as DMB. We have employed the values given in table 1, to present the shift in the exclusion bounds as shown in figure 2 for F DM = 1. The light blue curves represent the exclusion bounds for the fiducial choice of v 0 and v esc . The other solid lines represent exclusion limits with the DMB simulations. The dashed lines represent exclusion limits with DMO simulations. The orange and red lines corresponds to APOSTLE and ARTEMIS simulations respectively. For the APOSTLE DMB simulation, we observe a maximum change ∆ ∼ O(30%) in the cross section relative to the fiducial SHM.

Modified Maxwell-Boltzmann distribution
In the passing, we briefly discuss the modified version of Maxwell-Boltzmann velocity distribution in the galactic frame, truncated at the galaxy escape velocity v esc . This can be written in the form where the symbols have their usual meaning. In this work β = 1 has been chosen with the desire for an exponential cut-off. The exclusion bounds for this distribution of DM are shown by the dashed lines in figure 1. We use the fiducial values to generate the exclusion curves. We do not observe any significant changes in the exclusion bounds as compared to the standard MB distribution.

Beyond the Standard Halo Model
Several high resolution cosmological simulations suggest that the DM velocity distribution may depart from standard MB distribution, particularly in the high velocity tail [23,78,79]. Due to the sharp cut off at the escape velocity, the SHM over predicts the number of high energetic DM particles that are available for scattering. A possibility in this regard is to look for a non-SHM distribution of DM motivated by cosmological simulations. In what follows, we have considered some of these distributions and their implications on the DM-electron scattering rates in the direct detection experiments.

King Model
For a physical system of finite size, the truncated MB distribution, given in equation (7) may not be a natural solution of collision-less Boltzmann equation. The King distribution [80] is an alternate model which can be formulated self-consistently for a finite size DM halo. In this model instead of the escape velocity, the maximum DM particle velocity v max determines the cut-off criterion (v max < v esc ). It is predicted on the assumption that if a DM particle moves with v max at any position of the halo then it can reach the halo boundary where by construction the density vanishes [81]. This boundary is often called the truncated radius which represents the physical size of a halo. Such a finite size halo provides a more realistic description of galaxies as compared to the isothermal sphere. These so called lowered isothermal model is also preferred by simulations [82]. The distribution function can be written as         In the middle panel we show the exclusion bounds using King distribution with astrophysical parameters obtained from recent observations. The shaded bands represent observational error at the 95% confidence level. The solid light blue and grey lines correspond to the exclusion bounds with the fiducial parameters of the SHM and the King distribution respectively. In the lower panel the fractional changes relative to the fiducial choice (defined in equation (8)) for the APOSTLE and ARTEMIS simulations containing both DM and baryon (DMB) are shown by the solid orange and red lines respectively. Whereas the gray bands represent the deviation due to the recent astrophysical observations. The special feature observed around the DM mass of 50 MeV on the fiducial SHM curve in figures 3a and 3d is originating from the atomic structure of the Xenon atom. This is because around the aforesaid mass the maximum accessible DM energy crosses a threshold to ionize electrons from inner shell. This leads to an increment in the event rate and subsequently the bound becomes tighter. Further, as can be seen from figures 3a and 3d, compared to the fiducial values this special feature   figure 3g is a manifestation of the offset between the two kinks of the SHM and the King distribution.
As indicated in table 2, the best fit values of v 0 for DM only simulations is smaller than the fiducial values. Due to this in the high velocity region, less number of DM particles are available to interact with the electrons. This leads to weaker bounds in DM electron cross section 6 . This has been shown by the dashed lines in figure 3. The best fit values of v 0 with the DMB simulations lie close to fiducial value, implying exclusion limits close to the fiducial exclusion curve. This has been displayed by the orange and red solid lines for APOSTLE and ARTEMIS respectively.
We note that the deviation for King's model lies between 2% to 20% for Xe, Ge and Si target material for the best fit values of the APOSTLE DMB simulation. In contrast to that the best fit values of APOSTLE DMO and ARTEMIS DMO simulations seem to produce larger deviation from the SHM. The corresponding deviation induced by the recent astrophysical observations lie between 5% to 500% for Xe detector and 1% to 70% for both Ge and Si semiconductor detectors.

Double Power Law
The Double Power Law (DPL) distribution is an isotropic distribution of DM which has been obtained empirically [23]. The DPL very well describes the empirical matter distributions such as NFW, Hernquist, etc [23]. The matter distribution of a generic empirical double power law has the form where r s is the characteristic radius, ρ s is the characteristic density, α and γ determines the slope of the density profile at small and large radii respectively. For (α, γ) = (1, 3) this reduces to the NFW profile [84] while for (α, γ) = (1, 4) this reproduces Hernquist profile [85]. The DPL velocity where the symbols have their usual meaning. Cosmological simulation which take into account the quasi-static equilibrium nature of the virialised objects [32] and it's formation history attributed from hierarchical merging, smooth accretion and violent relaxations [86] favour such distribution [23]. Unlike the SHM, the velocity distribution in DPL smoothly goes to zero at the escape velocity. Thus it differs in the high velocity tail region from the SHM and predicts lower number of DM particles near the tail of the velocity distribution. For k → 0 it reduces to the SHM and for k = 1 it tends to the King distribution. The best fit parameters from cosmological simulations, that have been used in this work are given in table 2. We note that our best fit values is matched with [23].
Note that the DPL distribution can be viewed as a generalization of the King model, described in section 4.1. Therefore the impact on electron recoil event rates here would be similar to what has been discussed in section 4.1. The difference here is in the power index k. Any change in the numerical value of k would proportionally change the event rate and subsequently the direct detection limits. This can also be understood by comparing the relative fractional change in cross section depicted in lower panel of the figures 3 and 4. Like figure 3, in the upper panel of figure 4 the orange and red coloured lines correspond to the APOSTLE and ARTEMIS simulations respectively. Whereas in the middle panel, the solid grey lines correspond to the exclusion bounds for a representative k = 1.5 with the fiducial parameters. The shaded bands represent the uncertainties associated with the recent astrophysical observations of v 0 = 233 ± 6 km/s and v esc = 528 +24 −25 km/s for k = 1.5. We have also shown the fiducial SHM exclusion bounds by the light blue solid lines for reference.
We find that for APOSTLE DMB simulation with DPL distribution, the deviation obtained ranges between 5% to 20% for Xe and 4% to 14% for both Ge and Si semiconductor detectors. This deviation is fairly flat in the DM mass range of interest. The corresponding deviation induced by the recent astrophysical observation ranges between 2% to 100% for Xe detector, between 1% to 21% for Ge and 1% to 26% for Si semiconductor detectors.

Tsallis
The Tsallis distribution is explicitly derived through a factorization approximation of the Tsallis statistics [24] which is a generalisation of Boltzmann-Gibbs entropy. The distribution is widely used in high energy collisions [87], Bose-Einstein condensation [88], black-body radiation, neutron star [89], early universe cosmology [90] and superconductivity [91]. The velocity distribution function goes by the form where the symbols have their usual meaning. For this distribution with q < 1 the escape velocity is determined by the relation v 2 esc = v 2 0 /(1 − q). This inherent cut off criterion makes this distribution appealing as compared to the SHM. While for q > 1 escape velocity still remains a somewhat arbitrary parameter. In q → 1 limit the Tsallis distribution reduces to the Gaussian form of the SHM. Further from equation (13) it is evident that this distribution predicts a continuous and smooth fall near the tail, favoured by cosmological simulations [78,79,82]. In particular, it has been argued in reference [82] that the Tsallis distribution seems to fit better with Milky Way like simulations including Baryonic physics.
With the best fit values provided in table 2, the bound on σ e considering Tsallis as the distribution for DM is presented in the upper panel of figure 5. In the middle panel, the solid grey curve corresponds to the Tsallis distribution when fitted with the fiducial value. The grey bands correspond to the Tsallis distribution when fitted to the recent astrophysical observations of v 0 = 233 ± 6 km/s and v esc = 528 +24 −25 km/s. The deviation obtained for APOSTLE DMB simulation ranges in between 12% to 33% for the three set of target materials considered here. The corresponding deviations induced by the recent astrophysical measurement are between 3% to 290% for Xe targets, 1% to 90% for Ge and 1% to 100% for Si semiconductor detectors.

Mao et. al.
Mao et. al. [92] postulates an empirical model for the velocity distribution of DM having a wider peak and a steeper tail in comparison to the MB distribution. The Mao distribution function is given by where the symbols have their usual meaning. It is favoured by the simulations that have taken into account the sequence of mergers, violent-relaxation and accretion in the simulated halos [93]. Unlike other variants of the SHM, this empirical model is not based on a Gaussian but rather on an exponential distribution function having a power law cut-off at the binding energy or the equivalent escape velocity. We fit the distribution given in (14) with two different simulation namely APOSTLE and ARTEMIS with the corresponding best fit values given in table 2.        Assuming Mao as the velocity distribution for DM, in the upper panel of figure 6 we have presented the shifts in the exclusion limits for the best fit value given in table 2. In the middle panel the grey curves correspond to the exclusion bounds for a representative value of p = 2.7, while keeping v 0 and v esc at the fiducial choice. In the lower panel, with p = 2.7, the grey shaded band emerge due to the uncertainties associated with the recent astrophysical observations of v 0 = 233 ± 6 km/s and v esc = 528 +24 −25 km/s.
From the upper panel of figure 6 it is clear that compared to the fiducial SHM, the red curves set weaker bounds on DM-electron scattering cross section. This is mainly due to the large difference between the best fit values of v 0 and the fiducial values of the SHM. For the best fit DMO simulations, we observe a significant deviation from the SHM, compared to DMB simulations. Further the best fit APOSTLE DMB simulation which takes into account the baryonic contribution, shows a variation between 10% to 20% for all the three considered targets. This signifies the correlation between the

Comparison of deviation: observational and cosmological simulation
A discussion about the effects of best fit parameters for the various cosmological simulations on the derived exclusion bounds and its comparison with observational estimation of astrophysical parameters is now in order.
Comparing the bounds from various simulations, it can be inferred that DM only simulations show large deviation with respect to SHM fiducial as compared to simulation which takes into account the baryonic effects. Further among the two simulations, considered in the paper the APOSTLE DMB closely resemble the fiducial SHM scenario. For instance, the King, DPL, Tsallis, and Mao distributions show a maximum deviation of O(10%) relative to the fiducial SHM. We typically observe an order of magnitude deviation from the fiducial SHM for DM only counterpart of the considered simulations. This indicates that these deviations are correlated with the underlying dynamics of the simulations. Interestingly, for these non-SHM models, the recent astrophysical measurements indicate a maximum O(100%) and O(10%) deviations from fiducial SHM for Xe and semiconductor detectors respectively.
In table 3, we have provided a comparison of the percentage deviations from the fiducial SHM for sophisticated cosmological simulations and for the recent observational values of v esc and v 0 . As evident from the table, we conclude that for most of the halo distributions it is the updated astrophysical measurements that gives the most significant deviations in the exclusion bounds which are traditionally represented using the fiducial choice. This can imply significant reinterpretation of the conclusion drawn from direct detection experiments.

Conclusions
The non-observation of DM in the typical nuclear recoil direct detection experiments and the inability of GeV scale cold DM to address certain small scale structure formation issues have increased the interest in sub-GeV scale DM. An elegant avenue to probe such light DM is to consider the scattering of DM with electron in the direct detection experiments. Bounds on DM electron cross section are typically presented assuming the SHM for the DM distribution in our galaxy, with a fiducial choice for v 0 = 220 km/s and v esc = 544 km/s. However, recent progress in the measurement of these parameters shows a deviation from these values. In this paper we have systematically investigated the effects of uncertainties associated with the determination of these astrophysical quantities on the exclusion limits of DM electron cross section. We consider the uncertainties within the SHM and empirical models of DM distribution beyond the SHM that have been motivated by recent high resolution cosmological simulations.
We find that the exclusion bounds are expectedly sensitive to the population of DM particles in the high velocity tails of the distributions. Within the SHM the velocity distribution is assumed to be MB like. The tail shape is controlled by the Sun's circular velocity (v 0 ) and to a lesser extent by the escape velocity (v esc ). We find that within SHM for contact interaction between DM and electron, these uncertainties imply a 2% to 50% change in the event rates in the three target materials that have been considered. Further, inclusion of uncertainties in the local DM density leads to additional change in the exclusion bounds.
Going beyond the SHM, we have considered the simulation motivated King's model, Double Power Law, Mao, and Tsallis distributions. Relative to SHM all these non-SHM models fall smoothly near the high velocity tail, predicting less number particles in the region. This causes a reduction in the event rate. Therefore for the same set of astrophysical parameters, the non-SHM models seem to provide weaker bounds as compared to their SHM counterparts. Further for these models, depending on the fitted parameters we find that for most of the region in the parameter space the fractional changes in the cross section could vary substantially. Interestingly the amount of deviations observed from cosmological simulations fit can be traced to their treatment of baryonic content. In these models the deviation from the fiducial choice ranges between 1% − 30% for the APOSTLE hydrodynamic simulation which includes baryons. Whereas for these non-SHM models the uncertainties associated with the recent astrophysical observations lead to a maximum O(100%) change in the exclusion limit.
Note added: While this paper was under preparation a related work [94] appeared in the arXiv. Their treatment of halo uncertainties for semiconductor targets with Tsallis and Mao distribution is complementary with our discussion in section 4.

A A brief detail of the cosmological simulations
In this appendix we briefly summarize the cosmological simulations used in this paper. We give a brief and generic outline of the APOSTLE and ARTEMIS simulations along with their relevant resolutions.

A.1 APOSTLE
The APOSTLE (A Project Of Simulating The Local Environment) project comprises of a suite of highresolution cosmological simulation of 12 Milky way like environments, selected from a volumes of the ΛCDM universe. Each contain a pair of halo with a virial mass in the range 5×10 11 M to 2.5×10 12 M , with median values of 1.4 × 10 12 M for the more massive and 0.9 × 10 12 M for the less massive haloes. The project uses a combination of a hydrodynamical SPH implementation named Anarchy superimposed on the Tree-PM SPH Gadget code. The numerical resolution of the hydrodynamical simulation reaches up to a mass of 10 4 M per gas particle, using the sub-grid developed for the EAGLE project [95]. Each of the main galaxy haloes contain more than 2 × 10 7 particles. The high resolution initial conditions were generated using second-order Lagrangian perturbation theory. The Apostle simulation suite has three different resolution levels for the primordial gas (DM) particle masses of approximately 1.0(5.0) × 10 4 M , 1.2(5.9) × 10 5 M and 1.5(7.5) × 10 6 M respectively. Maximum gravitational softening lengths of 134 pc, 307 pc and 711 pc. In this work we use the DM distribution data of the APOSTLE IR halo marked A1, which fits the selection criteria as discussed in section 3.1 having a stellar mass of 4.88 × 10 10 M and a virial mass of 1.64 × 10 12 M .

A.2 ARTEMIS
The ARTEMIS (Assembly of high-ResoluTion Eagle-simulations of MIlky Way-type galaxieS) is a suite of high-resolution cosmological hydrodynamical simulations of Milky Way like haloes [72]. It uses the 'zoom-in' technique to simulate the halos at a high resolution with DM and by incorporating baryons with DM using SPH. The initial conditions were generated using the MUSIC code [96]. Haloes were selected from a base periodic box is 25 Mpc h −1 on one side with 256 3 particles. The initial conditions were generated at a redshift of 127 using a transfer function computed using the CAMB Boltzmann code cite for a flat ΛCDM WMAP cosmology was run to z = 0 using the Tree-PM SPH Gadget-3 code. Haloes were selected with masses in the range 8 ×