Possible evidence of QCD axion stars in HSC and OGLE microlensing events

Dark matter in the form of axions is expected to form miniclusters, and their dense regions can harbor compact axion stars. Such axion stars could be discovered by microlensing events. In particular, some candidate events reported by Subaru HSC and OGLE can be explained simultaneously if the axion stars with masses of the order of the Earth mass make up about $\sim27^{+7}_{013}$ percent of dark matter. For QCD axions, this corresponds to the axion mass in the range $10^{-9}-10^{-6}$ eV, which is consistent with the experimental constraints, as well as the cosmological anthropic window of parameters.


Introduction
A natural solution to the strong CP problem is associated with the existence of a light scalar field, the axion [1][2][3]. The axion particles can be produced in the early universe, and they can account for all or part of dark matter [4][5][6]. In addition, string theory predicts a number of axion-like particles, which are not necessarily associated with the strong interactions [7][8][9]. Dark matter in the form of axions is expected to form axion stars [10][11][12][13][14][15][16][17][18][19][20]. Microlensing observations can be used to set upper limits on the abundance of such compact objects [21][22][23].
We study whether microlensing events with short timescale light curves, recently reported in Refs. [24,25], can be explained by axion stars, and we infer the corresponding axion parameters.

Stable solution of axion star
Axions, or more generally axion-like particles, are described in field theory by a real scalar field φ(t, x) with the following potential: where m u 2.2 MeV, m d 4. 7 MeV, and f a are the up and down quark masses and the axion decay constant, respectively. The overall scale of the potential is related to the pion mass and the pion decay constant [26]: Λ 4 = f 2 π m 2 π where f π 92 MeV and m π 135 MeV. In the second line, we approximated the potential by a series expansion around φ = 0 assuming |φ| f a . The first term in the expanded potential acts as an effective mass term m a = m u m d /(m u + m d ) 2 Λ 2 / f a . The second term represents an attractive scalar self-interaction, which, together with self-gravity, can lead to an instability and formation of axion stars [10][11][12][13][14][15][16]. Due to the presence of the factor γ ≡ 1 − 3m u m d /(m u + m d ) 2 in the second term, for which we will use the fiducial value γ ∼ 0.34 in this letter, the field has the effective decay constant as f a = f a γ −1/2 . The fraction of cosmologically produced axions that end up in axion stars is subject of ongoing investigation. The axion stars can be stable, or, if they reach some critical mass, they can collapse emitting a burst of relativistic axions [27][28][29][30]. Among the stable solutions of axion stars, the most massive and compact one has the mass and radius given as [14,16] and Here we adopted, as an example, the axion parameters m a = 10 −8 eV and f a = 10 14 GeV that are consistent with the QCD axion. 1 The best-motivated mass range for the axions extends above and below this mass. If the Peccei-Quinn (PQ) symmetry is broken before inflation, and the axions are produced via the misalignment mechanism, the initial misalignment φ i determines the axion abundance. If φ i / f a ∼ 1, there is a limit m a > 10 −5 eV because the abundance of axion dark matter exceeds the observed values for smaller masses. However, different parts of the universe (or "multiverse"), on the scales much greater than today's horizon, could have different values of φ i [31,32]. Regions with a small misalignment may occupy a small volume of the universe, but they may contain a large fraction of potential observers [31][32][33]. This is particularly intriguing in view of the arguments that stars and planets hosting intelligent life could not have formed in a universe with a dark matter content greater than the observed value [33]. 2 Alternatively, the small axion mass can be reconciled with the dark matter abundance if some additional particle, a diluton decays and produces entropy, diluting the axion density, see, e.g., Refs. [35][36][37].
In the scenario of the PQ symmetry breaking before inflation, there is a potential problem with isocurvature perturbations in excess of the existing observational bounds [38][39][40][41][42][43][44]. This problem can be solved by a relatively low scale of inflation [45][46][47][48]. In addition, if the reheating temperature is sufficiently low and there is significant entropy production, the axion dark matter is consistent with f a 10 13 GeV [45,46,49], and 1 We note that the corresponding axion stars have masses close to the Earth mass ( 3 × 10 −6 M ) and a smaller radius than that of the Earth ( 6400 km). 2 The increased amount of dark matter causes a change in the matter-radiation equality temperature. If the matter comes to dominate the universe too early, the density perturbations grow and become nonlinear before recombination; the baryons and radiation get trapped inside the collapsing halos, and the baryon coupling to photons maintains the Jeans mass at a constant value as the collapse proceeds. As a result, large amounts of coupled baryon-radiation fluid are dragged into the potential wells created by clumps of dark matter, leading to a universe with supermassive black holes, photons, and neutrinos, but without stars and planets [33]. The same arguments apply to some other kinds of dark matter, such as, e.g., moduli [34]. the upper bound can be even higher, up to 10 14 GeV, because of the uncertainties in the modeling of axion strings, their evolution, and the spectrum of emitted axions. Thus we do not impose the cosmological bound of 10 −5 eV on the axion mass.

Gravitational microlensing
When an axion star and a background star are aligned along the line-of-sight direction of an observer, the star is multiply imaged by strong lensing [50]. The multiple images are usually not resolved by a telescope, and an observer can identify this lensing event from a time-varying brightness of the same star, which forms a characteristic light curve, because the source star, lens (here an axion star), and an observer have the relative motion. This is the so-called microlensing. The length scale characterizing the cross section of microlensing is the Einstein radius, given by where D ≡ d l d ls /d s (d l , d s and d ls are distances to lens, to source and between lens and source from an observer, respectively). If separation between a background star and a lens is smaller than R E on the sky, the microlensing occurs. Comparing Eqs. (3) and (4) manifests R E R AS for an axion star with M AS = 1.2 × 10 −6 M corresponding to the axion parameters taken in the equations, meaning that such an axion star is sufficiently compact to cause a microlensing event.
In this paper we assume that axion stars are formed in the early universe and constitute some mass fraction of dark matter (DM) that exists in the Milky Way (MW) (and also the Andromeda galaxy). A shape of the mass function of axion stars, even if formed, depends on details of physics inherent in the formation process of axion stars in the early universe. Generally speaking, axion stars with masses smaller than the mass scale (Eq. 2) can be formed, and the mass function can be extended to smaller mass scales. In this paper, for simplicity, we assume a monotonic (delta-function) mass function of a single mass scale M AS for axion stars, and that the mass fraction of axion stars to total DM in the MW is parameterized by f AS . The result with a monochromatic mass function in this paper can be always re-interpreted into general mass functions by convolving our result with a mass function of interest [51].

Microlensing constraints on axion parameters
To constrain parameters of axion stars, we use the published results of microlensing observations in the Subaru HSC data [24] and in the Optical Gravitational Lensing Experiment (OGLE) [52].
Niikura et al. [24] used the dense-cadence (2 min sampling rate), 7-hour long observation data of the Andromeda galaxy (M31) with Subaru HSC, and reported one possible microlensing event with short timescale of ∼ 1 hour, implying a lens with moon mass scale (∼ 10 −8 M ) [also see 53]. The direction of M31 is through the halo region of the MW, and we assume that the HSC microlensing event is due to an axion star that has the same spatial and velocity distributions in the MW and M31 halo regions as those implied from the DM halo model [see Ref. 24, for the details].
Ref. [52] used the publicly-available 2622 microlensing events, obtained from the 5-year OGLE observation of ∼ 5×10 7 stars in the Galactic bulge region, to discuss the populations of lensing objects. The OGLE events revealed a distinct population of 6 ultra-short timescale microlensing events ([0.1, 0.3] day scales) compared to the other majority of events that have a smooth distribution of microlensing timescales and are explained by standard populations of lensing objects such as brown dwarfs, main-sequence stars and white dwarfs in the standard model of the MW [also see 52, 54-57]. Niikura et al. [25] discussed that these ultra-short timescale events can be explained by compact objects such as primordial black holes that have the same spatial and velocity distributions as predicted in the MW halo model, if the compact objects are in the range of Earth mass scale (∼ 10 −6 M ). Since the OGLE observation is towards the Galactic bulge, the compact objects might be due to free-floating planets that are in the Galactic disk and bulge regions. Even in this case, theory needs to explain why free-floating planets have a preferred mass scale, i.e. Earth mass scale, rather than a smooth mass distribution. In this paper we use 1, 2, 2 and 1 microlensing events in each of 4 bins of the timescales t E = [0.1, 0.3] days as given in Fig. 5 of Niikura et al. [25], assuming that the events are due to axion stars. Here axion stars are distributed within the entire MW halo region extending to ∼ 200 kpc in radius from the Galactic center, and the microlensing occurs if an axion star comes across the Galactic disk and bulge regions and passes through in front of a background star in the Galactic bulge region along the line-of-sight direction of an observer.
Assuming that axions stars follow the spatial and velocity distributions predicted by the standard DM halo where t obs is a duration of the microlensing observation, and N s is the number of monitored stars used in the microlensing observation. dΓ/dt E is the differential event rate giving the expected number of microlensing events of a given timescale t E per unit observation time and per a single source star; the dimension is [events sec −2 ].
(t E ) is the detection efficiency of microlensing that quantifies the probability that a microlensing event of timescale t E is successfully recovered (detected) by the observation, which was estimated in Refs. [24] and [52] for the Subaru HSC and OGLE data, respectively. For the halo models of the MW and M31, we assume a Navarro-Frenk-White model [58] given in Ref. [59] that reproduces the observations such as the rotation curve. For axion stars we employ the two parameters, M AS and where L is the likelihood and Π is the prior. For the likelihood, we can safely assume that the observed number of microlensing events in a given timescale bin follows the Poisson distribution, because different microlensing events are safely considered independent due to the smallness of the microlensing optical depth for each of source stars: where N i is the number of observed microlensing events in the i-th timescale bin (t E,i ), and N exp,i is the expectation value predicted from the model assuming the abundance of axion stars in the MW and M31 halo regions (Eq. 5). As described above, we consider one bin and four bins of the microlensing light curve timescales for the HSC and OGLE events, respectively. We employ flat priors in the logarithmic scale for f AS and M AS . The contours in Fig. 1 shows the 68% and 95% credible regions on M AS and f AS , obtained from either of the HSC or OGLE short timescale microlensing events or the joint constraints, respectively. The marginalized constraints on f AS and M AS are summarized in Table 1 for each posterior. To be more precise, log f AS = −0.57 +0. 16 −0.18 and log[M AS /M ] = −5.52 +0.15 −0.10 for the HSC and OGLE events, respectively. It is intriguing that the two independent microlensing events of HSC and OGLE are consistent with each other, if the axion stars have an Earth mass scale (10 −6 M ) and make up about a few tens percent for the mass fraction to the total DM. (A similar conclusion can be applied to primordial black holes [25]. ) We also note that a slight difference between the HSC and OGLE parameter ranges in Fig. 1 does not necessarily imply a "tension" of these observations. The two sets of observations are sensitive to different mass ranges due to the difference in the observational time scales. The single event from HSC favors the mass scale M AS ∼ [10 −9 , 10 −6 ]M , but this detection has no implications for the OGLE sensitivity range, M AS ∼ [10 −6 , 10 −3 ]M where the other six events are detected. If one assumes a monochromatic mass function, as we have done for simplicity, then the overlap region is favored by the combination of HSC and OGLE data. However, many models predict an extended mass function, as e.g. Ref. [66]. A broader mass function can simultaneously be consistent with the HSC and OGLE events as the two observations probe different parts of the mass function.
We then discuss implications of the results in Fig. 1 on axion parameters. By marginalizing the posterior distribution P(M AS , f AS ) over f AS , we can obtain the projected posterior distribution of M AS . We then convert the credible interval of M AS to the posterior region of axion parameters using Eq. (2). Here we consider, for the axion parameters, m a and g aγγ , where g aγγ is the axion and photon coupling constant. We used the relation g aγγ = α EM /(2π f a )(E/N − 1.92), where α EM = 1/137 is the fine-structure constant of electromagnetic interaction, and E and N are anomaly coefficients depending on the QCD model. We use KSVZ model, where E/N = 0, as a benchmark model. We note that an uncertainty in E/N alters g aγγ by up to a factor of 10 [60]. Fig. 2 shows that the Subaru and OGLE ultra-shorttimescale microlensing events are consistent with an interesting range of axion parameters, including the QCD axions in the mass range m a ∼ 10 −8 eV. It is also interesting to note that the degeneracy direction of the microlensing constraint is quite complementary to existing constraints from terrestrial experiments. The region, denoted by "7R AS > R E (M AS )", indicates that the axion star is not compact enough to be regarded as a point lens object [23]. We also note that the axion mass range obtained in this paper is complementary to the direct detection of axion-like particle stars probed by optical magnetometer, e.g. the Global Network of Optical Magnetometers to search for Exotic physics (GNOME) col-   Fig. 1 using Eq. (2). Color scheme is similar to that of Fig. 1. Blue color bands are the allowed regions for the joint constraints of the HSC and OGLE results. Green line and yellow band denote a range of the QCD axion parameters [e.g. 60], displaying an overlapping region with the microlensing allowed regions around m a ∼ 10 −8 eV. Note that axion stars with masses lighter than the scale indicated by the microlensing constraints can be stable and exist (see text for details). The shaded region, labeled as "7R AS > R E (M AS )", denote the range of axion parameters where the axion star is not compact enough to be regarded as a point lens object [23]. Some bounds may be altered if a large fraction of the axion dark matter resides in axion stars. For comparison the existing upper limits on the coupling constant by the terrestrial experiments and previous works are shown in gray: SN1987A [61], HESS [62], Haloscopes [63], and CAST [64,65]. We assume all the galactic dark matter is in the form of diffuse axion in the Haloscope limits rather than shifting them by our estimation 1 − f AS . laboration, which has sensitivity around m a 10 −4 eV [67,68]. In this letter, we focus on the most massive and compact axion star solution in Eqs. (2) and (3) among possible stable solutions. Depending on the detail of the axion star formation scenario, lighter and more dispersed axion stars can be formed. In this case, our constraints in (m a , g aγγ ) are shifted toward the lighter axion mass.

Conclusion
Microlensing events in the HSC and OGLE data are consistent with ∼ 27 +7 −13 percent of dark matter in the form of axion stars. If the QCD axions make up all or most of dark matter, it is likely that ∼ 27 +7 −13 percent of the axions have condensed into axion stars [10][11][12][13][14][15]. For the QCD axion, the detected microlensing events correspond to the masses in the range 10 −9 − 10 −6 eV, which is consistent with the experimental bounds [69]. Further microlensing observations such as those using the Sub-aru HSC and the upcoming Vera C. Rubin Observatory telescope will enable us to test the possible contribution of axion stars to dark matter.

Acknowledgments
We thank Mark Hertzberg and Kaz Kohri for useful discussions. The work of A.