Stochastic fluctuations of bosonic dark matter

TL;DR

This paper discusses how the stochastic nature of ultralight bosonic dark matter affects experimental constraints, emphasizing the importance of accounting for field fluctuations when interpreting null results in laboratory searches.

Contribution

It highlights the need to incorporate the stochastic properties of bosonic dark matter fields in data analysis, correcting overestimated constraints from previous experiments.

Findings

Constraints were overestimated by factors of 3 to 10 in previous null experiments.

Experiments operating with measurement time shorter than the coherence time do not sample the full field distribution.

Proper interpretation requires modeling the stochastic fluctuations of the ultralight dark matter field.

Abstract

Numerous theories extending beyond the standard model of particle physics predict the existence of bosons that could constitute the dark matter (DM) permeating the universe. In the standard halo model (SHM) of galactic dark matter the velocity distribution of the bosonic DM field defines a characteristic coherence time $\tau_c$. Until recently, laboratory experiments searching for bosonic DM fields have been in the regime where the measurement time $T$ significantly exceeds $\tau_c$, so null results have been interpreted as constraints on the coupling of bosonic DM to standard model particles with a bosonic DM field amplitude $\Phi_0$ fixed by the average local DM density. However, motivated by new theoretical developments, a number of recent searches probe the regime where $T\ll\tau_c$. Here we show that experiments operating in this regime do not sample the full distribution of bosonic DM field amplitudes and therefore it is incorrect to assume a fixed value of $\Phi_0$ when inferring constraints on the coupling strength of bosonic DM to standard model particles. Instead, in order to interpret laboratory measurements (even in the event of a discovery), it is necessary to account for the stochastic nature of such a virialized ultralight field (VULF). The constraints inferred from several previous null experiments searching for ultralight bosonic DM were overestimated by factors ranging from 3 to 10 depending on experimental details, model assumptions, and choice of inference framework.

Figures

Equations

$$s(t)=\gamma \xi \varphi (t)\approx \gamma \xi {\Phi }_0\cos (2\pi f_{\varphi }t+\theta ),$$

$$p({\Phi }_0)=\frac{2{\Phi }_0}{{\Phi }_{\mathrm{DM}}^2}\exp \left(-\frac{{\Phi }_0^2}{{\Phi }_{\mathrm{DM}}^2}\right),$$

$$p(\gamma |D,f_{\varphi },\xi )=\mathcal{C}$$

$${\int }_0^{{\gamma }_{95\%}}p(\gamma |D,f_{\varphi },\xi )d\gamma =0.95.$$

$$p_{\det }(A_s|A)\propto p(A_s)2A\exp \left(-A^2-A_s^2\right)I_0\left(2A{A}_s\right),$$

$$p_{\mathrm{stoch}}(A_s|A)\propto p(A_s)\frac{2A}{(1+A_s^2)}\exp \left(-\frac{A^2}{1+A_s^2}\right).$$

$${\gamma }_{95\%}^{\mathrm{stoch}}\approx 3.0{\gamma }_{95\%}^{\det },$$

Full text

Stochastic fluctuations of bosonic dark matter

Gary P. Centers1,2, John W. Blanchard2, Jan Conrad3, Nataniel L. Figueroa1,2, Antoine Garcon1,2, Alexander V. Gramolin4, Derek F. Jackson Kimball5, Matthew Lawson2,3, Bart Pelssers3, Joseph A. Smiga1,2, Alexander O. Sushkov4, Arne Wickenbrock1,2

  

Dmitry Budker1,2,6,∗

  

Andrei Derevianko7

1Johannes Gutenberg-Universität, Mainz 55128, Germany

2Helmholtz Institute, Mainz 55099, Germany

3Department of Physics, Stockholm University, AlbaNova,10691 Stockholm, Sweden

4Department of Physics, Boston University, Boston, Massachusetts 02215, USA

5Department of Physics, California State University East Bay, Hayward, California 94542-3084,USA

6Department of Physics, University of California, Berkeley, CA 94720-7300,USA

7Department of Physics, University of Nevada, Reno, Nevada 89557, USA

∗Corresponding Author, Email: budker@uni-mainz.de

**Numerous theories extending beyond the standard model of particle physics predict the existence of bosons Dimopoulos and Giudice (1996); Arkani-Hamed et al. (2000); Taylor and Veneziano (1988); Damour and Polyakov (1994); Peccei and Quinn (1977a, b); Weinberg (1978); Wilczek (1978); Graham et al. (2015); Irastorza and Redondo (2018) that could constitute the dark matter (DM) permeating the universe. In the standard halo model (SHM) of galactic dark matter the velocity distribution of the bosonic DM field defines a characteristic coherence time $\tau_{c}$. Until recently, laboratory experiments searching for bosonic DM fields have been in the regime where the measurement time $T$ significantly exceeds $\tau_{c}$ DePanfilis et al. (1987); Wuensch et al. (1989); Hagmann et al. (1990); Asztalos et al. (2010); Graham and Rajendran (2013); Budker et al. (2014); Brubaker et al. (2017); Caldwell et al. (2017); Miller (2017); Chung (2015); Choi et al. (2017); McAllister et al. (2017); Alesini et al. (2017); Stadnik and Flambaum (2015); Grote and Stadnik (2019), so null results have been interpreted as constraints on the coupling of bosonic DM to standard model particles with a bosonic DM field amplitude $\Phi_{0}$ fixed by the average local DM density. However, motivated by new theoretical developments Marsh (2016); Marsh and Silk (2014); Hu et al. (2000); Hui et al. (2017); Arvanitaki et al. (2010, 2015), a number of recent searches Abel et al. (2017); Garcon et al. (2019); Wu et al. (2019); Terrano et al. (2019); Van Tilburg et al. (2015); Hees et al. (2016); Wcisło et al. (2018) probe the regime where $T\ll\tau_{c}$. Here we show that experiments operating in this regime do not sample the full distribution of bosonic DM field amplitudes and therefore it is incorrect to assume a fixed value of $\Phi_{0}$ when inferring constraints on the coupling strength of bosonic DM to standard model particles. Instead, in order to interpret laboratory measurements (even in the event of a discovery), it is necessary to account for the stochastic nature of such a virialized ultralight field (VULF) Geraci and Derevianko (2016); Derevianko (2018). The constraints inferred from several previous null experiments searching for ultralight bosonic DM were overestimated by factors ranging from 3 to 10 depending on experimental details, model assumptions, and choice of inference framework. **

It has been nearly ninety years since strong evidence of the missing mass we label today as dark matter was revealed Zwicky (1933), and its composition remains one of the most important unanswered questions in physics. There have been many DM candidates proposed and a broad class of them, including scalar (dilatons and moduli Dimopoulos and Giudice (1996); Arkani-Hamed et al. (2000); Taylor and Veneziano (1988); Damour and Polyakov (1994)) and pseudoscalar particles (axions and axion-like particles Peccei and Quinn (1977a, b); Weinberg (1978); Wilczek (1978); Graham et al. (2015); Irastorza and Redondo (2018)), can be treated as an ensemble of identical bosons, with statistical properties of the corresponding fields described by the SHM Kuhlen et al. (2014); Freese et al. (2013). In this work, our model of the resulting bosonic field assumes that the local DM is virialized and neglects non-virialized streams of DM Diemand et al. (2008), Bose-Einstein condensate formation Sikivie and Yang (2009); Davidson (2015); Berges and Jaeckel (2015), and possible small-scale structure such as miniclusters Jackson Kimball et al. (2018); Khlopov et al. (1985). To date it is typical to ignore such DM structure when calculating experimental constraints, and we demonstrate the general weakening of inferred constraints due to the statistical properties of the VULF within this isotropic SHM DM model. We note that including sub-halo structure Knirck et al. (2018); Foster et al. (2018), the formation of which is demonstrated in Refs. Chan et al. (2018); Lin et al. (2018); Veltmaat et al. (2018), can also affect experimental constraints.

During the formation of the Milky Way the DM constituents relax into the gravitational potential and obtain, in the galactic reference frame, a velocity distribution with a characteristic dispersion (virial) velocity $v_{\mathrm{vir}}\approx 10^{-3}c$ and a cut-off determined by the galactic escape velocity. Following Refs. Geraci and Derevianko (2016); Derevianko (2018) we refer to such virialized ultralight fields, $\phi({t,\bm{r}})$, as VULFs, emphasizing their SHM-governed stochastic nature. Neglecting motion of the DM, the field oscillates at the Compton frequency $f_{c}=m_{\phi}c^{2}h^{-1}$. However, there is broadening due to the SHM velocity distribution according to the dispersion relation for massive nonrelativistic bosons: $f_{\phi}=f_{c}+m_{\phi}v^{2}(2h)^{-1}$. The field modes of different frequency and random phase interfere with one another resulting in a net field exhibiting stochastic behavior. The dephasing of the net field can be characterized by the coherence time 111We note that there is some ambiguity in the definition of the coherence time, up to a factor of 2$\pi$, and adopt that which was used in the majority of the literature. See the discussion in the Supplementary Material. $\tau_{c}\equiv\left(f_{c}v_{\mathrm{vir}}^{2}/c^{2}\right)^{-1}$ Schive et al. (2014).

While the stochastic properties of similar fields have been studied before, for example in the contexts of statistical radiophysics, the cosmic microwave background, and stochastic gravitational fields Romano and Cornish (2017), the statistical properties of VULFs have only been explored recently. The 2-point correlation function, $\langle\phi(t,\mathbf{r})\phi(t^{\prime},\mathbf{r}^{\prime})\rangle$, and corresponding frequency-space DM “lineshape” (power spectral density, PSD) were derived in Ref. Derevianko (2018), and rederived in the axion context by the authors of Ref. Foster et al. (2018). While Refs. Derevianko (2018); Foster et al. (2018) explicitly discuss data-analysis implications in the regime of the total observation time $T$ being much larger than the coherence time, $T\gg\tau_{c}$, detailed investigation of the regime $T\ll\tau_{c}$ has been lacking (although we note that Ref. Foster et al. (2018) includes a brief discussion of the change in sensitivity 222 The authors discuss the change of sensitivity due to coherent averaging of the signal in the $T\ll\tau_{c}$ regime, $T^{1/4}\rightarrow T^{1/2}$, in their Appendix E. There is no mention of how the velocity and amplitude distributions would impact the derived limits. in this regime).

Here we focus on this regime, $T\ll\tau_{c}$, characteristic of experiments searching for ultralight (pseudo)scalars with masses $\lesssim 10^{-13}$ eV Abel et al. (2017); Garcon et al. (2019); Wu et al. (2019); Terrano et al. (2019); Van Tilburg et al. (2015); Hees et al. (2016); Wcisło et al. (2018) that have field coherence times $\gtrsim 1~{}\textrm{day}$. This mass range is of significant interest as the lower limit on the mass of an ultralight particle extends to $10^{-22}\,{\rm{eV}}$ and can be further extended if it does not make up all of the DM Marsh (2016). Additionally, there has been recent theoretical motivation for “fuzzy dark matter” in the $10^{-22}-10^{-21}\,{\rm{eV}}$ range Marsh (2016); Marsh and Silk (2014); Hu et al. (2000); Hui et al. (2017) and the so-called string “axiverse” extends to $10^{-33}\,{\rm{eV}}$ Arvanitaki et al. (2010). Similar arguments also apply to dilatons and moduli Arvanitaki et al. (2015).

Figure 1 shows a simulated VULF field, illustrating the amplitude modulation present over several coherence times. At short time scales ($\ll\tau_{c}$) the field coherently oscillates at the Compton frequency, see the inset of Fig. 1, where the amplitude $\Phi_{0}$ is fixed at a single value sampled from its distribution. An unlucky experimentalist could even have near-zero field amplitudes during the course of their measurement.

On these short time scales the DM signal $s(t)$ exhibits a harmonic signature,

[TABLE]

where $\gamma$ is the coupling strength to a standard-model field and $\theta$ is an unknown phase. Details of the particular experiment are accounted for by the factor $\xi$. In this regime the amplitude $\Phi_{0}$ is unknown and yields a time-averaged energy density $\langle\phi(t)^{2}\rangle_{T\ll\tau_{c}}=\Phi_{0}^{2}/2$. However, for times much longer than $\tau_{c}$ the energy density approaches the ensemble average determined by $\langle\Phi_{0}^{2}\rangle=\Phi_{\mathrm{DM}}^{2}$. This field oscillation amplitude is estimated by assuming that the average energy density in the bosonic field is equal to the local DM energy density $\rho_{\mathrm{DM}}\approx 0.4\,{\rm{GeV/cm^{3}}}$, and thus $\Phi_{\mathrm{DM}}=\hbar(m_{\phi}c)^{-1}\sqrt{2\rho_{\mathrm{DM}}}$.

The oscillation amplitude sampled at a particular time for a duration $\ll\tau_{c}$ is not simply $\Phi_{\mathrm{DM}}$, but rather a random variable whose sampling probability is described by a distribution characterizing the stochastic nature of the VULF. Until recently, most experimental searches have been in the $m_{\phi}\gg 10^{-13}\,{\rm{eV}}$ regime with short coherence times $\tau_{c}\ll 1\,{\rm{day}}$. However, for smaller boson masses it becomes impractical to sample over many coherence times: for example, $\tau_{c}\gtrsim 1\,{\rm{year}}$ for $m_{\phi}\lesssim 10^{-16}\,{\rm{eV}}$. Assuming the value $\Phi_{0}=\Phi_{{\rm{DM}}}$ neglects the stochastic nature of the bosonic dark matter field Abel et al. (2017); Garcon et al. (2019); Wu et al. (2019); Terrano et al. (2019); Van Tilburg et al. (2015); Hees et al. (2016); Wcisło et al. (2018).

The net field $\phi(t)$ is a sum of different field modes with random phases. The oscillation amplitude, $\Phi_{0}$, results from the interference of these randomly phased oscillating fields. This can be visualized as arising from a random walk in the complex plane, described by a Rayleigh distribution Foster et al. (2018)

[TABLE]

analogous to that of chaotic (thermal) light Loudon (1983). This distribution implies that $\approx$ 63% of all amplitude realizations will be below the r.m.s. value $\Phi_{\mathrm{DM}}$. Equation (2) Foster et al. (2018) is typically represented in its exponential form Knirck et al. (2018) (see Supplementary Material), and is well sampled in the $T\gg\tau_{c}$ regime. However, this stochastic behavior should not be ignored in the opposite limit.

We refer to the conventional approach assuming $\Phi_{0}=\Phi_{{\rm{DM}}}$ as deterministic and approaches that account for the VULF amplitude fluctuations as stochastic. To compare these two approaches we choose a Bayesian framework and calculate the numerical factor affecting coupling constraints, allowing us to illustrate the effect on exclusion plots of previous deterministic constraints Abel et al. (2017); Garcon et al. (2019); Wu et al. (2019); Terrano et al. (2019); Van Tilburg et al. (2015); Hees et al. (2016); Wcisło et al. (2018). It is important to emphasize that different frameworks to interpret experimental data than presented here can change the magnitude of this numerical factor Protassov et al. (2002); Cowan et al. (2011); Conrad (2015); Tanabashi et al. (2018), see Supplementary Material for a detailed discussion. In any case, accounting for this stochastic nature will generically relax existing constraints as we show below.

*Establishing constraints on coupling strength — * We follow the Bayesian framework Gregory (2010) (see application to VULFs in Ref. Derevianko (2018)) to determine constraints on the coupling-strength parameter $\gamma$. Bayesian inference requires prior information on the parameter of interest to derive its respective posterior probability distribution, in contrast to purely likelihood-based inference methods. The central quantity of interest in our case is the posterior distribution for possible values of $\gamma$, derived from Bayes theorem,

[TABLE]

The left-hand side of the equation is the posterior distribution for $\gamma$, where $D$ represents the data, and the Compton frequency $f_{\phi}$ is a model parameter. $\mathcal{C}$ is the normalization constant, and the likelihood $\mathcal{L}(\cdots)$ is the probability of obtaining the data $D$ given that the model and prior information, such as those provided by the SHM, are true. The integral on the right-hand side accounts for (marginalizes over) the unknown VULF amplitude $\Phi_{0}$, which we assume follows the Rayleigh distribution described by Eq. (2). For the choice of prior $p(\gamma,\Phi_{0})$ we use what is known as an objective prior Kass and Wasserman (1996): the Berger-Bernardo reference prior 333This approach is equivalent to starting with the marginal likelihood $\int d\Phi_{0}p(\Phi_{0})\mathcal{L(\cdots)}$ and using Jefferey’s prior to calculate the posterior Bernardo (1979). See details in the Supplementary Material. Berger and Bernardo (1992). Results from Bayesian inference are sensitive to the choice of prior Berger and Bernardo (1992), and we find better agreement with frequentist based approaches when using an objective prior rather than a uniform prior $p(\gamma)=1$ (see Supplementary Material).

It is important to note that experiments searching for couplings of VULFs to fermion spins (axion “wind” searches) are sensitive not only to the amplitude of the bosonic filed but also to the relative velocity between the laboratory and the VULF, which stochastically varies on a time scale $\approx\tau_{c}$ Graham and Rajendran (2013). The signal due to the axion wind is proportional to the projection of this stochastically varying velocity onto the sensitive axis of the experiment. Accounting for the stochastic nature of the relative velocity increases the uncertainty of the derived coupling strength for a given measurement. Axion-wind experiments can also utilize the daily modulation of this projection, due to rotation of the Earth, to search for signals with an oscillation period much longer than the measurement time $T\ll 1/f_{\phi}$. The unknown initial phase $\theta$ of the VULF sets the amplitude of this daily oscillation and also needs to be marginalized over. We describe how we account for stochastic variations of velocity and daily modulation in the relevant experiments Abel et al. (2017); Garcon et al. (2019); Wu et al. (2019); Terrano et al. (2019) in the Supplementary Material and focus solely on the stochastic variations of the amplitude, $\Phi_{0}$, here.

Using the posterior distribution, $p(\gamma|D,f_{\phi},\xi)$, one can set constraints on the coupling strength $\gamma$. Such a constraint at the commonly employed 95% confidence level (CL), $\gamma_{95\%}$, is given by

[TABLE]

The posteriors in both the deterministic and stochastic treatments are derived in the Supplementary Material. In short, the two posteriors differ due to the marginalization over $\Phi_{0}$ for the stochastic case, see the integral of Eq. (3). Assuming white noise of variance $\sigma^{2}$ and that the data are in terms of excess amplitude $A$ (observed Fourier amplitude divided by expected noise, an analog to the excess power statistic) we can derive the posterior for excess signal amplitude $A_{s}$. The posteriors are

[TABLE]

Here $A_{s}\equiv\gamma\times\xi\Phi_{\mathrm{DM}}\sqrt{N}/(2\sigma)$, $I_{0}(x)$ is the modified Bessel function of the first kind, and $p(A_{s})$ is effectively the prior on $\gamma$. In Fig. 2 we plot the normalized posteriors assuming $A$ at the 95% detection threshold $A^{dt}=\sqrt{-\ln{(1-0.95)}}$ and using Berger-Bernardo reference priors for $p(A_{s})$; we compare other choices of prior in the Supplementary Material. The derivation relies on the discrete Fourier transform for a uniform sampling grid of $N$ points and the assumptions of the uniform grid and white noise can be relaxed Derevianko (2018).

Examination of Eqs. (5), (6) and Fig. 2 reveals that the fat-tailed stochastic posterior is much broader than the Gaussian-like deterministic posterior. It is clear that for the stochastic posterior, the integration must extend considerably further into the tail, leading to larger values of $\gamma_{95\%}$ and thereby to weaker constraints, $\gamma_{95\%}^{\mathrm{stoch}}>\gamma_{95\%}^{\mathrm{det}}$. Explicit evaluation of Eq. (4) with the derived posteriors results in a relation between the constraints

[TABLE]

where the numerical value of the correction factor depends on CL and assumed value of $A$ (the factor increases for higher CL and decreases for smaller $A$).

This correction factor becomes $\approx 10$ when derived using a uniform prior, as discussed in the Supplementary Material. However, the result obtained with the uniform prior is not invariant under a change of variables (e.g. from excess amplitude to power). Additionally, using the objective prior yields better agreement with frequentist-based results of a factor $\approx 2.7$. For the pseudoscalar coupling, the additional stochastic parameters (field velocity and phase) increase this factor up to $\approx 8.4$ as shown in the Supplementary Material.

Ultralight DM candidates are theoretically well motivated and an increasing number of experiments are searching for them. Most of the experiments with published constraints thus far are haloscopes, sensitive to the local galactic DM and affected by Eq. (7). However, experiments that measure axions generated from a source, helioscopes or new-force searches, for example, do not fall under the assumptions made here. We illustrate how the existing constraints have been affected in Fig. 3 and provide more detailed exclusion plots for the axion-nucleon coupling $g_{aNN}$ Abel et al. (2017); Wu et al. (2019); Garcon et al. (2019); Vasilakis et al. (2009) and for dilaton couplings Hees et al. (2016); Van Tilburg et al. (2015); Wcisło et al. (2018) in the Supplementary Material.

Figure 3 shows that a few previously published constraints for the axion-nucleon coupling Abel et al. (2017); Wu et al. (2019); Garcon et al. (2019) no longer constrain new parameter space with respect to the new force constraint at $\approx 10^{-4}$ GeV*-1* Vasilakis et al. (2009).

*Conclusion – * To interpret the results of an experiment searching for bosonic DM in the regime of measurement times smaller than the coherence time, stochastic properties of the net field must be taken into account. An accurate description accounts for the Rayleigh-distributed amplitude $\Phi_{0}$, where the variation is induced by the random phases of individual virialized fields. Accounting for this stochastic nature yields a correction factor of $\approx 2.7-10$, relaxing existing experimental bosonic DM constraints in this regime. In the event of a bosonic DM discovery, the stochastic properties of the field would result in increased uncertainty in the determination of coupling strength or local average energy density in this regime.

It is important to note that observational knowledge of the local distribution of DM can constrain stochastic behavior of the amplitude (energy density). The smallest features observed so far are on the order of $\approx 10$ pc Iocco et al. (2015) (corresponding to a $m_{\phi}\approx 10^{-21}$ eV coherence length), however the analysis in Ref. Iocco et al. (2015) performs radial averages which would remove the stochastic variation discussed in this paper.

*Data Availability – * All conclusions made in this paper can be reproduced using the information presented in the manuscript and/or Supplementary Material. Additional information is available upon reasonable request to the corresponding author. For access to the experimental data presented here please contact the corresponding authors of the respective papers.

*Acknowledgments – * We thank Eric Adelberger and William A. Terrano for pointing out the need to account for the unknown phase in the CASPEr-ZULF Comagnetometer analysis. We thank Kent Irwin, Marina Gil Sendra, and Martin Engler for helpful discussions and suggestions. We thank M. Zawada, N. A. Leefer, and A. Hees for providing raw data for the published deterministic constraints. We also thank Jelle Aalbers for helpful discussions and expert advice on the blueice inference framework. Jan Conrad appreciates the support by the Knut and Alice Wallenberg Foundation. This project has received funding from the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation programme (grant agreement No 695405). We acknowledge the partial support of the U.S. National Science Foundation, the Simons and Heising-Simons Foundations, and the DFG Reinhart Koselleck project.