Astrophysical and Cosmological Probes of Dark Matter

While astrophysical and cosmological probes provide a remarkably precise and consistent picture of the quantity and general properties of dark matter, its fundamental nature remains one of the most significant open questions in physics. Obtaining a more comprehensive understanding of dark matter within the next decade will require overcoming a number of theoretical challenges: the groundwork for these strides is being laid now, yet much remains to be done. Chief among the upcoming challenges is establishing the theoretical foundation needed to harness the full potential of new observables in the astrophysical and cosmological domains, spanning the early Universe to the inner portions of galaxies and the stars therein. Identifying the nature of dark matter will also entail repurposing and implementing a wide range of theoretical techniques from outside the typical toolkit of astrophysics, ranging from effective field theory to the dramatically evolving world of machine learning and artificial-intelligence-based statistical inference. Through this work, the theory frontier will be at the heart of dark matter discoveries in the upcoming decade.


Introduction
Astrophysical and cosmological observations have historically played a critical role in the study of dark matter, underpinning our confidence that there is a missing mass component of the Universe. The evidence that observational measurements provide for dark matter is collected across many length scales. The earliest hints for dark matter arose from its gravitational effects on galaxies, explaining the observed flatness of rotation curves [1][2][3][4]. Gravitational lensing has also detected dark matter surrounding galaxy clusters [5]. On yet larger scales, the cosmic web of large-scale surveys [6], as well as the fluctuations of the cosmic microwave background (CMB) [7], have both been integral in the development of the cold dark matter (CDM) paradigm, where 85% of the Universe's matter budget is dark. A complete theory of particle dark matter 1 will ultimately describe how it interacts with visible matter, as well as whether it interacts with other dark states in its own separate sector. Moreover, any such theory will be successful on both the largest scales of the Universe as well as the smallest (i.e., sub-galactic) scales. Over the next decade, astrophysical and cosmological probes will provide powerful tests of fundamental questions about dark matter, playing a unique and complementary role to the terrestrial dark matter experimental program. This review will focus on five specific theory questions where concrete advancements are anticipated during this time period: • Is the Cold Dark Matter paradigm correct?
In the CDM paradigm, dark matter is collision-less and non-relativistic during structure formation. A natural consequence of this is the prediction of an abundance of low-mass dark matter halos down to ∼ 10 −6 M [8]. Observations that provide information on the matter power spectrum at small scales and various redshifts, therefore, will play a pivotal role in confirming the CDM hypothesis. Evidence of small-scale power Snowmass2021 Theory Frontier: Astrophysical and Cosmological Probes of Dark Matter suppression could, for example, suggest that dark matter is warmer (i.e., not nonrelativistic) during structure formation [e.g. 9], is not collision-less [e.g. 10], is wavelike rather than particle-like [e.g. 11], or underwent non-trivial phase transitions in the early Universe [e.g . 12]. As we will discuss, upcoming astrophysical surveys have the potential to start probing halo masses to much lower values and/or higher redshifts than previously accessible, opening the opportunity of definitively testing the CDM hypothesis.
• Is dark matter production in the early Universe thermal?
The observed relic abundance of dark matter can be explained through a thermal freeze-out mechanism (see [e.g. 13,14] for recent reviews). In this picture, dark matter is kept in thermal equilibrium with the photon bath at high temperatures through weak annihilation processes. Once dark matter becomes non-relativistic, dark matter is still allowed to annihilate, but the reverse process is kinematically forbidden. The continued annihilation of dark matter causes its comoving number density to be Boltzmann suppressed, until it freezes out due to Hubble expansion overcoming the annihilation rate. This process sets the present-day dark matter abundance. Importantly, the predicted abundance is sensitive to the detailed dark matter physics, including its particle mass as well as its specific interactions with the Standard Model. Weakly Interacting Massive Particles (WIMPs) provide a classic example of the freeze-out paradigm. In this case, a O(GeV-TeV) mass particle that is weakly interacting yields the correct relic abundance. As we will demonstrate, upcoming astrophysical surveys will have the opportunity to definitively test key aspects of the WIMP hypothesis by searching for the rare dark matter annihilation and decay products that arise from the same interactions that set its abundance in the early Universe. A combination of improved instruments and the so-far non-observation of WIMPs has also led to the exploration of probing dark matter candidates that are lighter or heavier than the canonical WIMP window, and which often have a non-thermal origin in the early Universe. This broadening of the possible dark matter candidates that one can search for in indirect detection will continue to be driven by the theory community.
• Is dark matter fundamentally wave-like or particle-like?
Model-independent arguments that rely on the phase-space packing of dark matter in galaxies have been used to set generic bounds on its minimum allowed mass. In particular, a fermionic dark matter candidate can have a minimum mass of ∼ keV [15], while a bosonic candidate can have a minimum mass of ∼ 10 −23 eV [11]. Moreover, when the dark matter mass is much less than ∼ eV, its number density in a galaxy is so large that it can effectively be treated as a classical field. Oftentimes referred to as "axions" or "axion-like particles" (ALPs), these ultra-light bosonic states can have distinctive signatures due to their wave-like nature. The QCD axion [16][17][18][19], originally introduced to address the strong CP problem, is a particularly well-motivated dark matter candidate for which there are clear mechanisms for how to generate the correct abundance today [20][21][22][23]. In this framework, the axion mass and coupling are fundamentally related to each other through the symmetry-breaking scale of the theory.

Snowmass2021 Theory Frontier: Astrophysical and Cosmological Probes of Dark Matter
As we show, upcoming searches for astrophysical axions will have the sensitivity reach to probe highly-motivated mass ranges for the QCD axion.
• Is there a dark sector containing other new particles and/or forces?
In a generic and well-motivated theory framework, dark matter can exist in a "dark sector" that communicates with the Standard Model through specific portal interactions. Within the dark sector, there can be multiple new states, as well as new forces that mediate interactions between the dark particles. Recent theory work has demonstrated classes of dark sector models that yield the correct dark matter abundance (see [e.g. 24] for a review), oftentimes for lower dark matter masses than expected for WIMPs. Dark sector models can lead to a rich phenomenology for both astrophysical and terrestrial dark matter searches, as we will discuss. Two properties of the dark sector where upcoming astrophysical surveys will be able to make decisive statements are the presence of self interactions between dark matter particles [25] and new light degrees of freedom.
• How will the development of numerical methods progress dark matter searches?
Given the sheer volume and complexity of data expected from astrophysical surveys in the upcoming decade, the development of effective observational and data analysis strategies is imperative. Novel machine learning and statistical tools will play an important role in maximizing the utility of these datasets. In particular, scalable inference techniques and deep learning methods have the potential to open new dark matter discovery potential across several frontiers. Another critical numerical component to harness the anticipated flood of astrophysical data in the next decade is the further development of cosmological and zoom-in simulations needed to interpret the survey results. We will comment on how such simulations are essential for understanding the implications of particular dark matter models on small-scale structure formation. This list is not intended to be comprehensive, but rather to provide well-motivated examples of areas where fundamental advancements are expected with upcoming astrophysical and cosmological probes. We have divided this white paper into two separate discussions reflecting what we can learn about dark matter from its interactions with visible matter in astrophysical systems (Sec. 2) as well as its early-Universe behavior and its role in the formation of structure (Sec. 3). Each section briefly reviews some of the most promising observational probes for tackling the specific theory questions delineated above. Sec. 4 is dedicated to the exciting advancements expected in applications of statistics and machine learning to astrophysical studies of dark matter. We conclude in Sec. 5.
Complementarity with additional White Papers: We note that there are a number of white papers which contain results complementary to the discussion we provide here. A non-exhaustive list includes Dark Matter Numerical Simulations [26], Data-Driven Cosmology [27], Dark Matter Physics from Halo Measurements [28], Ultra-heavy Particle Dark Matter [29], and Puzzling Excesses and How to Resolve Them [30]. We encourage interested readers to look to these related white papers for further details of how the search for dark matter will proceed in the coming decade.

Dark Matter Interactions with Visible Matter
Historically, a strong motivation for the existence of an interaction between dark and visible matter arises from the simple and compelling cosmologies described by the WIMP miracle or freeze-in production of sterile neutrinos. These scenarios have also long motivated indirect detection searches: the same interactions that generate dark matter could also be occurring today and allow it to decay or annihilate into detectable signatures arising from astrophysical sources. However, in the past decade, a substantial portion-although certainly not all-of the well-motivated parameter space for these models has been excluded (see [e.g. 31,32]). Astrophysical searches for dark matter have broadened in perspective as the theoretical community has realized that the potential mass and interactions dark matter could have are much, much broader. In this section, we will highlight this paradigm shift, demonstrating that while conventional searches continue, ideas to probe significantly heavier and lighter dark matter are appearing and will continue to be developed in the coming years.

X-ray and γ-ray Dark Matter Signatures
If the dark matter of our Universe can decay or annihilate, then one of the most promising channels for determining its particle nature is the detection of high-energy photons. In the past decade, considerable progress has been made in photon-based indirect detection for dark matter with masses O(keV-TeV), which we will predominantly focus on in this section. This improvement has not been solely driven by the experimental observatories: often theoretical insights have led to dramatic leaps forward in our ability to probe dark matter. This bidirectional approach to progress must continue in the coming years.
The fundamental ingredient for indirect searches are observational datasets. A partial summary of the present and future landscape is provided in Fig. 1. At the highest energies, significant progress will be achieved beyond the existing reach of H.E.S.S. [35], HAWC [36], and similar observatories, through a combination of CTA [37], SWGO [38], and LHAASO [39] (the last of which is already in operation, see [e.g. 40]). The combined dark matter discovery potential of these telescopes is significant. CTA has the possibility to discover the higgsino [41], one of the most well-motivated WIMP candidates. However, its ability to do so will depend on whether or not the broad program to understand electroweak effects for TeV scale dark matter can determine its annihilation spectrum and cross section with sufficient accuracy [42][43][44][45][46][47][48][49][50][51]. Observations at these energies further open the path to probing dark matter with masses above the unitarity limit, m χ 100 TeV [52]. The space of models, cosmologies, and production mechanisms for such ultra-heavy dark matter is being actively developed, see [e.g. 29]. Discovering these models requires a detailed understanding of the spectrum of particles that emerge from their annihilation or decay, and how those states propagate to Earth. In recent years, public codes have been developed for the propagation of high-energy states [53][54][55][56][57]. For the spectra, the most widely used approach exploits an analogy with colliders so that Pythia [58][59][60] can be used for the calculation, as in PPPC4DMID [61, 62]. This analogy breaks down at higher energies-indeed, existing LHAASO projections simply end at m χ = 100 TeV due to an absence of theoretical calculations available above that scale [63]. The first steps towards reliable spectra at higher masses has recently been taken in [64], although there remains significant work. The importance of these developments for Snowmass2021 Theory Frontier: Astrophysical and Cosmological Probes of Dark Matter Figure 1: The exposure of existing and upcoming X-ray and γ-ray instruments which can search for the decay or annihilation of dark matter. Significant progress in the coming years is expected at O(TeV-PeV) energies, a full exploitation of which will require theoretical developments in the models and spectra of heavy dark matter. Such work will also complement instruments searching for even higher energy photons, such as PAO [33,34]. At O(keV-GeV) energies, the expected observational progress is far more modest at the level of exposure. Exposure-the product of effective area, E, and observation time, T -partially controls how many photons an instrument will detect on average for a given dark matter flux. We caution that this is just one metric by which instruments can be compared: for certain dark matter searches, the field of view or energy resolution can be critical, and then the improvements made at lower energies will be more substantial. Regardless, improved analysis strategies will be crucial to further enhance the dark matter reach for this lower band. At the top of the figure, we highlight the approximate mass range of several canonical particle dark matter (DM) scenarios; one can roughly associate this mass range with the corresponding energy range of probes, although this connection is only approximate.
In the O(keV-GeV) band, observations must be made from space as the interaction of photons with the atmosphere does not produce a sufficiently detectable signature on the Earth's surface. This sets a fundamental limitation: a 1 m 2 instrument operating for a decade has an exposure of ∼10 12.5 cm 2 s. At keV and GeV energies, this is roughly the sensitivity already achieved by XMM-Newton [ [73], and COMPTEL [74]. In the longer term, instruments such as Athena [75], AMEGO [76], and e-ASTROGAM [77] can improve our sensitivity, but at many energies, the best anticipated datasets are already on disk. Progress will therefore be critically reliant on new insights for how to exploit the data. This is happening on many fronts, including identifying new objects in which to search for dark matter signals, such as newly discovered Milky Way dwarfs [78][79][80][81][82][83], galaxy cata-logs [84,85], or dark substructure [86][87][88][89][90][91][92][93][94][95][96][97][98][99][100][101][102][103][104][105], as well as the development of techniques such as cross correlation with different datasets [106][107][108][109][110][111][112][113], improvements in our modeling of diffuse backgrounds [114][115][116][117], the extension to axion searches as we describe in Sec. 2.4, and the exploitation of the dark matter brightness of the ambient Milky Way [32,[118][119][120]. To expand on a single example, in [120] it was demonstrated that the roughly twenty years of X-ray images collected by XMM-Newton, when combined with the insight that all of these observations occur through a column density of the Milky Way, allowed for a search for dark matter decay that was more than an order of magnitude stronger than previous analyses. The results were strong enough to considerably disfavor the longstanding 3.5 keV line dark matter anomaly [121,122] (although see also [123,124]). More generally, there remains several unexplained dark matter anomalies in this energy window whose resolution will depend on further insights from the theory community; for an extended discussion, we refer to [30].
An additional strategy that has been developed recently considers dark matter that scatters and becomes captured within celestial bodies. The dark matter can then annihilate, generating a signal that depends on the mediator between the dark and visible sectors. If the mediator is short-lived (or insufficiently boosted), the annihilation products will remain within the celestial object, raising its temperature. Instead, a long-lived (or sufficiently boosted) mediator leads to annihilation products outside the body, which can then be searched for by telescopes. The strongest constraints on long-lived or boosted mediator models are due to γ-ray searches, as the γ-ray backgrounds for celestial objects are very low. An excellent candidate is the Sun, and solar γ-ray searches have been performed using both Fermi [125][126][127][128] and HAWC [129], yielding strong constraints on GeV−TeV dark matter [130][131][132]. Optimizing for both proximity and size, the next best celestial body is Jupiter, and Fermi observations have been used to constrain sub-GeV dark matter [133]. More broadly, analogous emission from the full population of brown dwarfs and neutron stars can constrain sub-GeV to TeV dark matter [134] (see also [135]). These searches are inherently multimessenger: solar dark matter searches for neutrinos have been performed with Super-Kamiokande [136], IceCube [137], and ANTARES [138]. Above O(100) GeV, the neutrinos will be attenuated as they exit the Sun, and a long-lived mediator again improves detectability [130,139]. The scenario involving short-lived mediators can be studied with optical and infrared telescopes, including Hubble, JWST, and Roman observations of neutron stars [140][141][142][143][144], white dwarfs [141,145], population III stars [146][147][148][149], and brown dwarfs and exoplanets [150].

Indirect Searches with Astrophysical Neutrinos
While searches for dark matter in the electromagnetic spectrum may be more extensively developed, there is no fundamental reason that the first discovery could not happen through a different channel. If that channel is neutrinos, then the coming decade will be particularly exciting.
Already the possibility of dark matter decaying or annihilating to neutrinos can be probed from MeV to PeV masses through a combination of instruments ranging from Borexino to IceCube. For a recent review, see [151]. The detection of O(TeV−PeV) neutrinos at IceCube is particularly tantalizing. While the experimental collaboration has produced limits under the assumption that the observed flux does not originate from dark matter [152,153], a complete understanding of where these neutrinos originate is lacking. (Although in 2017, a ∼ 300 TeV neutrino event was shown to be coincident with a flaring γ-ray blazar [154], which would make this the third extraterrestrial source ever detected in neutrinos after the Sun and Supernova 1987A.) There has been considerable work by theorists to determine whether the IceCube flux could have a component with a dark matter origin, see [e.g. 54, [155][156][157][158][159][160][161][162][163], and the question remains open.
A fundamental aspect of high-energy neutrino searches, and of high-energy indirect detection more generally, is its multimessenger nature: generically, a signal of heavy dark matter annihilation or decay will appear in multiple channels. If PeV scale dark matter decays to neutrinos, there will be a significant probability for the hard neutrinos to emit a W or Z boson and thereby produce additional Standard Model final states, including photons. As a consequence, γ-ray datasets provide important context for dark matter interpretations of the IceCube dataset (see [e.g. 53,118,164]).
This multimessenger strategy will continue to future instruments that will probe neutrinos at ever higher energies. Already, there are tools available that make a partial accounting of the underlying physics of the unbroken Standard Model [64, 165], which build on the observation that electroweak effects are generically relevant for heavy dark matter [61, 62]. As discussed for the corresponding photon signals, work remains to fully understand these processes. The importance of such exploration is emphasized by the wide array of upcoming observatories such as ARIANNA, RNA-G, POEMMA, Grand, IceCube-Gen2, and KM3NET, which have the potential to probe dark matter with masses up to the GUT scale of 10 15 GeV [65, 166-169].

Dark Matter Signals from Charged Cosmic-rays
High-energy cosmic-rays have long been a probe of new physical phenomena in the Milky Way. This is particularly true for antimatter cosmic-rays, a field which PAMELA and AMS-02 have brought into a precision era, providing a challenge to our understanding of the antimatter sources and generating claims of a possible dark matter contribution.
In antiprotons, there are claims of an excess in the AMS-02 data peaking near 10 GeV [170,171]. The excess appears robust to systematic uncertainties on the production cross sections, cosmic-ray injection rates, and the effects of propagation through the interstellar medium and heliosphere [172,173]. The anomaly has a local significance of 3-5 σ and is consistent with the possibility that it is generated by the same dark matter models which could be generating the Galactic center γ-ray excess [170,[172][173][174] (a compelling possibility given the astrophysical uncertainties for the two signals would be uncorrelated, although see [175]). The results at present do not account for the full correlation matrix of the dataset, although attempts to estimate the covariance suggest that it can significantly reduce the significance of the excess [176][177][178]. At present, the AMS-02 collaboration has not released their correlation matrix, which will be critical in establishing or repudiating the antiproton excess. Looking forward, the GAPS experiment will provide an alternative measurement of low-energy antiprotons [179] that will improve the modeling of the propagation of antinuclei in the interstellar medium and heliosphere [180]. More broadly, the combination of high-precision measurements of multiple cosmic-ray species and the observation of cosmic-ray protons and electrons from different time periods will further reduce the astrophysical uncertainties.

Snowmass2021 Theory Frontier: Astrophysical and Cosmological Probes of Dark Matter
For sufficiently massive dark matter, any flux of antiprotons is expected to be accompanied by heavier anti-nuclei cosmic-rays [181][182][183]. While dark matter can only provide at most a small excess over the astrophysical backgrounds in antiprotons, the backgrounds are highly suppressed for more massive anti-nuclei [184][185][186][187][188][189], making them a potential smoking gun signal of new physics. Tantalizingly, the AMS-02 collaboration has presented results of the detection of several antihelium events [190], where essentially no background events were expected. While tentative, if confirmed, this result could revolutionize cosmic-ray and high-energy physics. Heavier anti-nuclei are, however, plagued by a number of uncertainties, particularly as related to their productions (see [e.g. 191, 192]). At present, these uncertainties imply their predicted fluxes can vary by orders of magnitude [181,183]. Future low-energy collider measurements will be instrumental in reducing those uncertainties [191], and the first such measurements have recently been provided [193,194].
Cosmic-ray positrons are another potential dark matter probe. Given that positrons quickly lose their energy as they propagate in the interstellar medium, the sources must be increasingly localized to produce higher energy cosmic-rays. The rising positron fraction measured by PAMELA [195], Fermi [196], and AMS-02 [197] has been widely discussed as a putative signal of dark matter annihilation or, alternatively, of nearby pulsars or supernova remnants (see [e.g. [198][199][200][201][202][203][204][205][206][207][208][209]). Dark matter explanations are particularly challenged by Planck measurements of the CMB temperature and polarization power spectra [210][211][212] and have become increasingly fine-tuned although not entirely ruled out (see [e.g. 213,214]). Regardless, AMS-02 measurements remain a highly sensitive probe of dark matter annihilation [215]; for instance, one can search for spectral features associated with the dark matter mass [216,217], although these must be interpreted carefully as astrophysical sources can also generate such features [218,219]. Existing and upcoming cosmic-ray and electromagnetic observations will be used to develop a deeper understanding of the properties of positron sources [220][221][222][223][224], which will further advance the program of probing dark matter with cosmic-rays.
Radio searches: The mixing of axions and photons is greatly enhanced in the magnetospheres of neutron stars, owing to the large ambient magnetic fields and the resonant  [261] (solid red, labeled "GBT"); shown for comparison is the 10σ discovery limit for SKA in the same frequency band (transparent red). These results are compared with the expected QCD axion parameter space (yellow), two benchmark QCD axion models (brown, dot-dashed), existing constraints from haloscopes (blue) [262][263][264][265][266] and CAST [267] (grey), and indirect searches using magnetic white dwarfs (MWDs) [243,260].
amplification possible due to the ambient plasma [227,261,[268][269][270][271][272][273][274]. The characteristic plasma mass near typical pulsars spans ∼ 0.1 − 100 µeV [275], which roughly corresponds to the range of masses for which axions can simultaneously solve the strong CP problem and constitute dark matter [276], and further to the frequency band of modern radio telescopes. The radio signal is expected to appear as a forest of spectral lines centered about the axion mass, with each line arising from a single neutron star in the Galactic population [261,271]. If dark matter is predominantly in miniclusters rather than smoothly distributed, the events will instead appear as transients spanning hours to weeks [277]. Initial estimates indicate that near-future radio interferometers like the Square Kilometer Array may be capable of discovering the QCD axion [270,277]. Searches using existing infrastructure (including the Effelsberg 100-m telescope, the Green Bank Telescope, and the Very Large Array) are already underway and have leading limits on the axion-photon coupling in the mass range 1 m a 20 GHz [261,278,279] (see Fig. 2 for the most recent analysis).
Recently, there has been significant theoretical progress in our understanding of the radio signal, including a careful treatment of photon refraction, resonant cyclotron absorption, plasma-induced line broadening, anisotropic response of the medium in the photon production process, and general relativistic effects [261,273,274,280]. Yet many open questions remain, including how axions and photons mix in a highly magnetized inhomogeneous plasma, how do charge distributions in active pulsars and magnetars impact the radio flux, what are the properties and distributions of neutron stars in dense dark matter environments, do we expect strong deviations from dipolar magnetic fields (and if so how does this impact the radio signal), how are axions distributed on astrophysical scales (i.e., do they reside in tidally disrupted axion miniclusters, and if so what are the properties of these objects in the Galactic Center), and how can we exploit the spatio-temporal properties of the radio signal to improve analyses. These are questions to be answered in the next decade, and the answers have the potential to establish radio searches as a powerful and robust probe of axion dark matter.
X-ray and γ-ray searches: High-energy photons emitted from astrophysical sources (including galaxies, blazars, supernovae, and quasars) may convert to axions in galacticand cluster-scale magnetic fields (axions could also be produced in stellar cores [281]). The conversion probability depends on the magnetic field strength and configuration along the photon trajectory, as well as the plasma frequency. At sufficiently large photon energies, the photon-to-axion conversion probability becomes O(1) and is thus capable of generating large absorption features in the electromagnetic spectrum. The efficiency of this conversion process decreases at lower energies (at a fixed axion mass), generating small oscillatory features in the observed spectrum. Using this idea, constraints on the axion-photon coupling have been set using X-ray [282][283][284][285][286] and γ-ray [287][288][289][290][291][292] telescopes for masses m 0.1 µeV (excluding couplings g aγγ 10 −12 GeV −1 for m a 10 −11 eV), and they apply regardless of whether or not axions contribute to dark matter. Future progress will be aided by improved high-energy observations and by further understanding of galactic and cluster-scale magnetic fields.

Emission of Dark Sector States from Compact Objects
Standard Model particles in high-density environments, such as stars and supernovae, can emit new weakly-coupled states that might exist beyond the Standard Model. The particles emitted could either be dark matter themselves or, alternatively, part of a broader weaklycoupled "dark sector," which is potentially necessary to endow sub-GeV dark matter with the correct relic abundance. The emission process can result in either observable deviations from the Standard Model predictions on short timescales or long-term global changes to the evolution of the compact object, both of which can be constrained. The theory frontier has long played a crucial role in bridging the gap between astrophysical probes, complex multi-body Standard Model calculations, and inference on new particle properties.
The paradigmatic example of "short-term" constraints comes from the successful explosion of Supernova 1987A and the detection of neutrinos for the predicted ∼ 10-second-long cooling phase [293,294]. The qualitative agreement of this observation with Standard Modelonly numerical simulations [295,296] has been used to constrain the properties of the QCD axion with ever-increasing fidelity and sophistication [253,[297][298][299][300][301][302][303]. The theory frontier is still grappling with these calculations. Upcoming challenges will be centered on application of effective field theory techniques, which promise important changes in expected rates for Standard-Model-only processes (such as nuclear and neutrino matrix elements and scattering rates), as well as for beyond-the-Standard-Model rates.
A successful explosion of Supernova 1987A could also have been inhibited by a large dark sector [303][304][305][306]. Broadly, Supernova 1987A provides the strongest bounds for any number of new particles in the O(1-100 MeV) mass range, being cut off at high masses by Boltzmann suppression from thermal production in the core, which attains temperatures between 30-100 MeV [295,307]. The power of these bounds at larger couplings is generally limited by the existence of a "trapping" limit. At couplings above the trapping limit, the new particles are more tightly coupled than the Standard Model neutrinos and thus are unable to drain the energy the neutrinos were observed to have taken away. Nevertheless, with the increasing fidelity of numerical simulations [308], the region beyond the trapping limit is a clear avenue for future theoretical and numerical investigations.
A second type of stellar constraint pertains most relevantly to particles that are substantially lighter and more weakly coupled than the MeV range in which supernova limits excel. In this mass range, dark matter must generally be non-thermal in order to avoid constraints on the number of new radiation degrees of freedom in the early Universe. Thus, such limits typically focus on bosonic particles such as the QCD axion, axion-like particles, Higgs-portal scalars, and the dark photon. A powerful approach is to require new particle emission to be subdominant to photon emission in stars [306,[309][310][311][312], analogous to the supernova bound wherein the new emission is limited by the neutrino emission. If this rule were to be violated, stellar lifetimes would be unacceptably short or other emission properties would change beyond measured values. A different method of constraining light particle emission from stars has been made possible recently by LIGO-Virgo Collaboration observations of gravitational waves and the corresponding mass census enabled by these observations, which are revealing the characteristics of Pop-III stellar progenitors for the first time [313,314]. One qualitative prediction of Standard Model-only astrophysics is the existence of a "black hole mass gap" formed from these objects at a characteristic mass scale slightly below 50 M [315,316]. New particle emission, gravitational trapping of dark matter, or dark matter coevolution all could change this mass scale [317][318][319][320][321][322]. Theory frontier activities in the stellar domain promise to illuminate new, weakly-coupled particles that are not probed by other mechanisms [323].
A different route by which dark matter could impact a "long-term" observable of compact objects is via the formation of a super-radiant cloud that extracts angular momentum from a central black hole [324][325][326][327][328]. For a sufficiently low-mass dark matter particle, this could lead to detectable changes in the observed black hole spin distribution [239,. Future theoretical explorations will lead to a more comprehensive understanding of the impacts of backreaction of the superradiance on the conditions necessary to support the superradiant instability.

Dark Matter Origins and Structure Formation
Traditional methods of dark matter indirect detection are centered around the idea of detecting the Standard Model byproducts of dark matter interactions within astrophysical systems. New dark matter physics can also be indirectly observed through its impact in the early Universe and on the subsequent formation and evolution of collapsed structures of matter (see Fig. 3). Using such detection methods can cover regions of dark matter parameter space that are complementary to regions probed by traditional methods. Crucially, gravitational indirect detection is necessary to probe certain classes of dark matter theories that require the extreme conditions of the early Universe or dark sector theories that have rich phenomenology but do not couple to known physics except through gravity. In this section, we explore the theoretical progress made in this area, often going hand-in-hand with advancements in numerical simulations, and discuss the importance of continued theory ef- forts to take full advantage of the influx of cosmological and astrophysical data expected over the next decade.

Early-Universe Evolution
The present-day abundance of dark matter may have been established after a period in which dark matter particles were in thermal and chemical equilibrium with some larger thermal reservoir, or it may have been established despite never attaining thermal equilibrium with its environment. The former case, which we refer to as "thermal" dark matter, whether that thermal equilibrium is attained with the Standard Model itself or with a secluded thermal bath, is commonly considered to be UV-insensitive: the initial conditions of the dark matter abundance will be erased by the thermal equilibrium condition. The latter case of "nonthermal" dark matter leads to a large number of observables that can potentially persist from the earliest moments of dark matter genesis. One exception to the rule that thermal dark matter is insensitive to early-Universe microphysics could be if the dark sector underwent a first-order phase transition. This could generate a gravitational wave signal [e.g. [373][374][375][376][377], which would remain thermally decoupled, even if dark matter particles were in thermal equilibrium. In fact, the occurrence of the first-order phase transition requires the dark sector to be at least self-thermalised. The most widely studied example of this kind is a dark sector Higgs mechanism, resulting in a phase transition that can be first order if the sector contains a number of bosonic degrees of freedom [374]. Dark sector confinement may also lead to the generation of a gravitational wave signal, if the global symmetry is large enough at the time of breaking [e.g. [378][379][380][381]. Such scenarios have been studied in the context of axion models [380] and bound state formation such as dark quark nuggets [378]. Theoretical progress is necessary before the gravitational wave phenomenology of first-order phase transitions can be studied reliably: it has been demonstrated that even in perturbative models, two-loop thermal effects must be included to achieve better than O(1) numerical accuracy [382][383][384]. Estimations of the gravitational wave spectrum in non-perturbative hidden sector models have used low-or high-energy effective field theories (see [379] for a comparison), but their validity breaks down in the vicinity of the phase transition. Future theory frontier efforts will be critical in understanding this break down and extending the range of reliable predictions for models of new physics.
Non-thermal dark matter may produce gravitational wave signals, but more generally leads to a multitude of other tests of early-Universe physics. For instance, the evolution of the post-inflationary axion can result in an abundance of axion miniclusters [385][386][387][388][389][390], which can potentially be discovered with dedicated search strategies [391][392][393][394][395][396][397][398]. Accurate predictions for the spectrum and abundance of such miniclusters requires early-Universe simulations of the post-inflationary axion, which are inherently difficult to perform due to a separation of scales and large systematic uncertainties. Therefore, additional large-scale simulations and systematic tests are needed in the future to improve results.
Lastly, we note that the abundances of thermal and non-thermal dark matter are both affected by the evolution of the Universe between the end of inflation and the onset of Big Bang nucleosysthesis (BBN), and we cannot assume that the Universe was radiation dominated throughout this period [399]. If dark matter chemically decouples while the Universe is not radiation dominated, the larger Hubble rate at a given temperature causes an earlier freeze-out and an enhanced relic abundance [400][401][402][403][404]. If the subsequent transition to radiation domination involves the creation of new Standard Matter particles, as occurs after an early matter-dominated era, then the relic abundance of dark matter is diluted [400,[405][406][407][408]]. An early matter-dominated era also enhances dark matter density perturbations on scales that enter the horizon prior to the onset of the final radiation-dominated epoch [409]. Significant progress has been made in understanding the impact this has on the abundance of sub-Earth-mass microhalos [409,410] and how the minimum halo mass depends on the properties of dark matter [411][412][413][414], as well as the properties of the particle responsible for the early matter-dominated era [415][416][417][418]. These microhalos provide a new observational probe of the early Universe; their impact on the dark matter annihilation rate can be constrained using the isotropic gamma-ray background [413,415,419], and they can be detected gravitationally using pulsar timing arrays [420][421][422][423] and observations of stellar microlensing events in galaxy clusters [424,425]. Further study is needed to understand whether the CMB and the 21cm background could further constrain dark matter annihilation in early-forming microhalos.

Early-Universe Structure Formation
The standard cosmological model of ΛCDM describes the large-scale structure of the Universe extremely well. The presence of CDM is crucial through its contribution to the overall energy density and through its density perturbations at early times. Comprising about 84% of the matter content in the Universe [426], the gravitational influence of dark matter is key in the formation of structure. Thus, cosmological observations provide important and unique insight into new dark matter physics that disturb the predictions of CDM.
Light degrees of freedom: A wide variety of dark matter models introduce new light degrees of freedom in the Universe at early times. In particular, dark sectors may contain light or massless force carriers, such as dark photons, that thermalize dark matter at a temperature that generally differs from the Standard Model bath. The introduction of new relativistic species alters the expansion rate of the Universe during radiation domination, in turn affecting CMB anisotropies [e.g. 427] and predictions of the abundances of light elements created during BBN [e.g. 428,429]. Standard contributions to the energy density of radiation during the BBN and CMB eras include photons and neutrinos, assuming their masses are not too large; however, neutrinos decouple from the photon bath at temperatures ∼ 1 MeV, and their energy density contribution is encoded in the parameter N eff . For the three active neutrinos of the Standard Model, recent calculations yield N SM eff = 3.044 [430][431][432][433], and deviations ∆N eff ≡ N eff − N SM eff from this value could imply the existence of non-standard physics. Current CMB and BBN observations constrain N eff to be close to its Standard Model value [e.g. 426]; therefore, in order to incorporate new massless species in a cosmological model, the dark sector temperature has to be lower than that of the photon bath to avoid contributing too much to the relativistic energy budget.
Thermal dark matter itself can contribute to the relativistic energy density during BBN if its mass is 20 MeV. Additionally, if the dark matter relic abundance is set through a standard freeze-out process, dark matter may annihilate into neutrinos or visible particles, which in turn alter weak interaction rates that determine primordial abundances. Dark matter coupled to neutrinos or charged particles generates a positive or negative contribution, respectively, to N eff . Thus, cosmological observations can provide robust bounds on the mass of dark matter, for a given spin and annihilation channel [434,435,[435][436][437][438][439]. In the future, CMB-S4 will obtain a sensitivity to new thermalized, light relics corresponding to ∆N eff < 0.06 at 2σ [440]. Dedicated theory work will be needed to understand the implications of CMB and BBN constraints across a variety of dark sector models.
CMB spectral distortions: Measurements of the CMB energy spectrum provide an opportunity to search for new physics that impacts the thermal history of the Universe. Deviations of the CMB spectrum from a perfect blackbody, referred to as spectral distortions, are sensitive to processes that inject energy into (or extract energy from) the photon-baryon plasma at redshifts z 2 × 10 6 . Current measurements of the CMB spectrum show that it is extremely close to a blackbody with a present-day temperature T 0 = 2.72548 ± 0.00057 K [441,442], with spectral distortions smaller than a few parts in 10 5 [441]. Proposed experimental concepts could probe spectral distortions at least three orders of magnitude smaller [443], thus opening new windows into exotic physics in the early Universe. Several known processes within the standard ΛCDM cosmological model generate spectral distortions [444]. In addition, spectral distortions may be generated through dark matter interactions with Standard Model particles.
Dark matter annihilating or decaying into photons or electrically charged particles would inject energy into the photon-baryon plasma, hence distorting the CMB energy spectrum [445]. CMB anisotropies are significantly more sensitive to s-wave annihilations than spectral distortions [445]; however, spectral distortions could be more sensitive to p-wave annihilations [446], depending on the specifics of the dark matter model. In addition, spectral distortions can constrain decaying particles with lifetimes 10 6 sec τ 10 12 sec, to which CMB anisotropies are insensitive. For such short lifetimes, the decaying particle could only comprise a small fraction of the total dark matter abundance [447].
Alternatively, dark matter particles may extract energy from the photon-baryon plasma if they scatter elastically with photons, electrons, or nuclei [446,448]. Indeed, if dark matter is heavier than ∼ 1 keV, it is non-relativistic by z ∼ 2 × 10 6 and therefore cools down adiabatically faster than the thermalized photon-baryon plasma. Elastic scattering would therefore lead to a systematic transfer of heat from the plasma to the dark matter fluid. This effect is increasingly large for light, thus more abundant, dark matter particles. While current spectral distortion limits only constrain elastically-scattering dark matter particles with masses m χ 100 keV, proposed experiments could extend this sensitivity to ∼ GeV masses. Continued theoretical work is needed to ensure robust predictions [e.g. 449].
CMB anisotropies and dark matter annihilation/decay: Dark matter annihilation or decay via electomagnetic channels injects energy into the photon-baryon plasma, which increases the free-electron fraction x e around and after cosmological recombination at z 1100. The increase in x e delays the last scattering epoch and affects the photon diffusion scale (and hence the damping of CMB anisotropies at small scales). Additionally, an increased x e in the low-redshift tail of recombination suppresses CMB anisotropies on small scales and increases large-scale polarization fluctuations, an effect qualitatively similar to an increase in the reionization optical depth [450].
Through these effects, CMB anisotropies are sensitive to very rare dark matter annihilation or decay processes [212,451], since as little energy as ∼ 1 eV per baryon suffices to significantly alter the ionization history. In particular, even for dark matter produced in a standard freeze-out scenario, residual annihilation at z 1100 may have a significant impact on CMB anisotropies, long after they no longer change the dark matter abundance. Planck data constrain the dark matter s-wave annihilation cross section to σv 3 × 10 −28 (m χ /GeV), ruling out a standard freeze-out production of dark matter for m χ 10 GeV [7]. Planck data also constrain the lifetime of decaying dark matter to τ 10 24 -10 25 sec, depending on the decay channel [452,453], orders of magnitude larger than the age of the Universe. Nearfuture and planned CMB-anisotropy missions could achieve a factor of ∼ 30 improvement in the sensitivity over current Planck limits on dark matter annihilation and decay [454].
More exotic dark matter candidates, such as PHBs, can be probed by CMB anisotropies through similar effects. PBHs can inject energy in the photon-baryon plasma through either Hawking radiation for masses M 10 17 g [453] or accretion-powered radiation for masses M M [455][456][457]. In contrast with the rather clean and well-understood physics involved in dark matter annihilation or decay [212,458], the complex physics of accretion is highly uncertain and much theory work remains to be done to make existing limits more robust.
While the effects described above arise from changes to the average free-electron fraction, inhomogeneous energy injection from dark matter would also lead to spatial fluctuations in x e [459,460]. These fluctuations would induce non-Gaussianities in CMB anisotropies, which could be a complementary avenue to probe energy injection from dark matter.
CMB anisotropies and dark matter scattering: Elastic scattering between dark matter and Standard Model particles in the early Universe can alter the evolution of perturbations, impacting CMB temperature, polarization, and lensing anisotropies. Scattering processes heat the dark matter fluid and induce a drag force from the exchange of momentum. The primary effect of the scattering is inhibiting the clustering capabilities of dark matter, thus washing out structure on a variety of observable scales. As a result, the CMB power spectra experience damping at large multipoles , corresponding to small angular scales on the sky, for models where dark matter decouples from baryons prior to matterradiation equality [461][462][463][464][465][466][467]. For models such as millicharged dark matter, scattering takes place at later times, and CMB anisotropy currently provides some of the best observational bounds on such models [e.g. 464,468,469]. The theoretical developments and numerical implementation into Boltzmann codes, such as CAMB and CLASS, have enabled CMB searches of dark matter scattering with baryons [461][462][463][464][465][466][467], electrons [470,471], photons [472][473][474], and neutrinos [475,476]. The strength of the interaction is a key parameter that controls the amount of power suppression in the CMB primary anisotropy, but including an energy or velocity dependence of the interaction influences the shape of the suppression, potentially allowing a way to distinguish between various scattering models. Even in the case of dark sectors, dark matter scattering with dark radiation [477][478][479] produces features in the CMB power spectra that can be differentiated from other models [480]. Changing N eff can have a similar effect of suppressing the CMB damping tail, but possible degeneracies can be broken using CMB lensing anisotropies [481].
Upcoming ground-based instruments, such as the Simons Observatory [482] and CMB-S4 [440], will measure the CMB with better precision at high and a much higher angular resolution than current experiments, allowing for greater sensitivity to dark matter interactions. In terms of theoretical development, most investigations have focused on thermal relic models, while future investigations of models where dark matter is produced non-thermally will also be of high interest in context of CMB probes. Furthermore, there is a notable synergy between the CMB primary anisotropy and other probes of the structure growth in the Universe, which can be further exploited in self-consistent analyses of CMB with other observables, for specific dark matter models.

21-cm line at high redshifts:
The redshifted 21-cm line of neutral hydrogen presents a unique probe of the post-recombination Universe, prior to the birth of the first stars (cosmic dark ages, z = 30-100) and right after it (cosmic dawn, z = 5-30)-see [e.g. 483,484] for reviews. The corresponding absorption signal imprinted on the CMB backlight during cosmic dawn [485] captures the state of hydrogen gas and various microphysical processes that control it at this early epoch; mapping the absorption signal provides a (3D) view of the Universe at epochs that no other probes can reach. Many dark matter processes can affect the temperature of baryons during cosmic dawn, in turn affecting the 21-cm signal. The experiments targeting this era can use this connection to probe dark matter interactions with the visible sector.

Snowmass2021 Theory Frontier: Astrophysical and Cosmological Probes of Dark Matter
A well-motivated model that can affect the 21-cm signal involves dark matter scattering with baryons through a light mediator [466,468,[486][487][488][489][490][491][492]. For millicharged dark matter, current 21-cm measurements (e.g., from EDGES in the global signal [493] or HERA for fluctuations [494]) exclude millicharges as low as q χ ∼ 10 −6 e [488] (where e is the electron charge), even if less than a percent of dark matter is millicharged. Alternatively, models of dark photons that kinetically mix with Standard Model photons, with mixing parameter ∼ 10 −7 [495], can create a radio background that affects the 21-cm signal [495][496][497]. Additionally, a classical WIMP, for instance, can heat gas through annihilations [498,499], and similar effects occur in case of decaying dark matter and dark-photon dark matter [500,501]; each in turn can affect the thermal properties of baryons during cosmic dawn. Even more exotic models, such as evaporating or accreting PBHs, can produce heating [502,503].
The 21-cm signal during cosmic dawn is sensitive to structure formation because the first galaxies emit photons that produce Wouthuysen-Field coupling [485], critical for the 21-cm absorption feature to arise. Those first galaxies were hosted in halos with masses ∼ 10 6 M [504 -507], which correspond to fluctuations with wavenumbers as large as k ≈ 100-200 Mpc −1 . As a consequence, altering the corresponding scales in the linear matter power spectrum affects the abundance of the first galaxies and the timing of the 21-cm signal [508]. Small scales are often very sensitive to dark matter microphysics and can be used to probe warm dark matter and other beyond-CDM models.
While there are claims that the global (sky-averaged) signal has already been detected [493], the tomographic signal is yet to be measured [494] from interferometers such as HERA [509], LOFAR, MWA, or SKA [510]. However, it is expected that these observatories can yield constraints on the matter power spectrum up to k = 100 Mpc −1 at the ∼ 10% level [508], which can place powerful constraints on specific dark matter candidates such as ETHOS models [511], warm dark matter [512], or fuzzy dark matter [513].
High-redshift galaxy luminosity function: The luminosity function (i.e., counting the number of galaxies versus their luminosity) provides a tracer of the abundance of dark matter halos at the redshifts of measurement. UV luminosity functions are constructed from the HST Ultra Deep Fields (HUDF), where the UV light from the most massive galaxies at z = 4-10 is detected by the visible/IR filters at the HST [514,515]. JWST observations will significantly contribute to constraining the population of early galaxies in the coming decade. There measurements can be used to determine the matter power spectrum during reionization up to k = 10 Mpc −1 [372,516,517] and can provide probes of warm dark matter [518][519][520][521], fuzzy dark matter [522,523], and ETHOS models [524].
Lyman-α forest: The clustering of matter at intermediate redshifts traced by the redshifted forest of Lyman-α absorption lines is a sensitive probe of dark matter physics. Several authors have forward-modeled the Lyman-α forest to place lower limits on the thermal relic warm dark matter mass of O(3-5) keV [525][526][527][528], where the details of the constraints depend on astrophysical assumptions about the temperature, density, and redshift evolution of baryons in the intergalactic medium [529,530]. Other dark matter models such as those with dark matter-baryon interactions, ultra-light axions, and PBHs have also been constrained using similar methods [531][532][533][534], although there are dark matter-related modeling challenges in some cases [e.g. 535]. From a theoretical standpoint, the development of Lyman-α forest emulators [e.g. 536] has contributed to recent advances and will become increasingly impor-tant to enable robust, joint inference of cosmological, astrophysical, and dark matter physics in the coming decade. From an observational standpoint, dark matter analyses have generally been performed using ∼ 10s of high-resolution spectra (e.g., VLT, HIRES/KECK; [525]), ∼ 100s of intermediate-resolution spectra (e.g., XQ-100; [527]), or ∼ 1000s of low-resolution spectra (e.g., SDSS/BOSS; [528]). Ongoing spectroscopic surveys including DESI will significantly enhance the number and redshift coverage of available high-resolution quasar spectra, potentially allowing for percent-level measurements of the Lyman-α flux power spectrum on small scales [537].

Present-Day Structure
Probes of small-scale structure at low redshifts have recently emerged as a key means to test a variety of dark matter properties, including its production mechanism, primordial temperature, self-and Standard Model-interactions, and minimum particle mass. These probes can broadly be categorized according to whether they rely on observations of the baryonic contents of low-mass halos or not. Here, we summarize the current status and future prospects for each of these probes, along with key theoretical considerations. In particular, we emphasize that connecting theoretically motivated predictions for dark matter's gravitational imprints to precise analytic and simulation-based predictions for small-scale structure distributions is a critical area for work over the next decade.
Massive galaxy clusters provide a unique opportunity to stringently test the CDM paradigm of structure formation. Combining strong and weak gravitational lensing detected in high-resolution images of massive clusters has revealed that the dark matter subhalos of cluster galaxies are less massive and less spatially extended compared to those hosting equivalent luminosity field galaxies, indicating that tidal stripping of dark matter is efficient in these dense, violent environments-see review by [538] and, for a critical analysis of the range of lens modeling methodologies, see [539] and [540] for recent developments. Comparison of the derived subhalo mass function from observed cluster lenses with CDM simulations has revealed that while the abundance and mass function of substructures was well reproduced, the radial distribution of subhalos was discrepant [541]. Subhalos are more concentrated in the inner regions of observed clusters than predicted by CDM simulations.
A recent study of Galaxy-Galaxy Strong Lensing (GGSL) in clusters found that observed small-scale cluster substructures (on ∼ 5-10 kpc scales) are more efficient strong lenses than predicted by CDM simulations by more than an order of magnitude [542]. Further theoretical investigation will be needed to evaluate if this large discrepancy arises from hitherto undiagnosed systematic issues within simulations, or if in fact this serves as a hint for deviations from the CDM paradigm. Numerical effects arising from the resolution limits of simulations that lead to artificial subhalo disruption [543] cannot account for the order of magnitude GGSL gap, as they are at most a 20% effect [544]. Importantly, baryonic feedback processes must be carefully investigated as a potential culprit for the discrepancy, as it is well understood that they alter the internal structure of cluster galaxies. This motivates a comparison of simulated galaxy clusters across several independent CDM simulations.
The GGSL discrepancy could potentially be revealing that dark matter might not be collisionless, especially in the extremely dense cluster environments. For elastic self-interactions that are velocity dependent (with large interaction cross sections), halos can undergo gravothermal collapse just as stellar systems. When this occurs, the inner halo develops a negative heat capacity due to outward transfer of energy [545], resulting in a significant enhancement of the concentration of the inner density profile-precisely in the direction needed to address the GGSL discrepancy [546]. The transformation produced by core-collapse motivates the investigation of this particular class of self-interacting models more deeply. Moreover, totally inelastic self-interactions can result in a collapse time-scale up to two orders of magnitude shorter than for the elastic case [547], for smaller interaction cross sections. Further theoretical study is needed to map out the space of self-interacting dark matter models that can potentially account for the GGSL discrepancy.
Dwarf galaxies are the smallest dark matter-dominated baryonic systems in the Universe and form in halos with M halo 10 10 M , down to the galaxy formation threshold of M halo ∼ 10 8 M [548,549]. Thus, the smallest "ultra-faint" dwarf galaxies [550] are a particularly sensitive probe of low-mass halo abundances, which reflect dark matter's smallscale gravitational clustering. To date, ultra-faints have exclusively been detected within the virial radius of the Milky Way as satellite galaxies; recent wide-field photometric surveys, including the Sloan Digital Sky Survey, Pan-STARRS1, and the Dark Energy Survey, have increased the number of known Milky Way satellites to roughly 60 systems (see [551] for a recent census of the Milky Way satellite population). These observations have been used to constrain the warm dark matter particle mass at the level of 6.5 keV, the particle mass of fuzzy dark matter at the level of 2.9 × 10 −21 eV, and the dark matter-proton interaction cross section at the level of 10 −29 cm 2 [552] (also see [548,[553][554][555][556][557]). Over the next decade, observational facilities, including the Vera C. Rubin Observatory and the Nancy Grace Roman Space Telescope, are expected to significantly improve upon current dwarf galaxy discovery power, both within and well beyond the virial radius of the Milky Way [e.g. 558,559]. With these upcoming improvements, continued theoretical developments are key to properly interpreting possible deviations of halo abundances from the CDM expectation.
In addition to their mass abundances, the individual properties and population statistics of dwarf galaxies are also sensitive to dark matter microphysics. For example, dark matter self-interactions can both suppress the inner densities of halos in a mass-dependent fashion [560][561][562][563][564] and eventually drive these systems towards gravothermal core collapse [545,[565][566][567][568][569][570][571][572]. These effects provide a mechanism for explaining the observed diversity of galactic rotation curves [573][574][575][576], although current observations are not yet sufficient to distinguish this scenario from feedback-affected CDM halos [577]. Conservative constraints based on the observed inner-most densities of dwarf galaxies such as Draco [e.g. 578], coupled with constraints from galaxy clusters [e.g. 579], demonstrate that self-interacting dark matter models with velocity-dependent interactions must undergo some degree of gravothermal collapse [580]. This observation may provide a mechanism to explain observations suggesting that the most centrally-dense Milky Way dwarfs also have the smallest orbital pericenters [581]. Improved theoretical modeling of dwarf properties and populations will be needed to harness the full potential of current and upcoming surveys that are amassing information on dwarf galaxies of Milky Way-like systems [582,583], as well as improved stellar kinematic data on the Milky Way's dwarfs from observatories like Gaia [584].
Strong gravitational lensing allows us to detect low-mass halos within lens galaxies and along the line of sight via their effect on the observed multiple images of a background source. This process is purely gravitational and independent of whether these low-mass halos contain any baryons. It thus provides a unique approach to test dark matter models by probing the low-mass end of the halo and subhalo mass functions beyond the local Universe. The detection of low-mass halos with strongly lensed quasars is mainly based on so-called flux ratio anomalies-that is, changes induced to the relative flux of the multiple images [e.g. [585][586][587]. In images of strongly lensed galaxies, a local change of the surface brightness distribution of the data, reflecting a change in the relative position of the images, is the tell-tale signature of the presence of low-mass halos [e.g. [588][589][590][591][592]. This method is often referred to as the gravitational imaging technique.
Flux ratio anomalies and the gravitational imaging technique are complementary approaches that are subject to individual and shared sources of systematic errors [593,594]. The two techniques can also differ in their sensitivity to low-mass haloes. Depending on the size of the background source, flux ratio anomalies typically probe the halo and subhalo mass functions down to masses as low as ∼ 10 7 M and potentially below. The number of available lens systems and the precision of the flux measurement then set the precision with which one can constrain the halo and subhalo mass functions [595]. For example, [596] and [593] constrain the warm dark matter particle mass at 5.2 (5.5) keV at the 95% confidence level using eight (seven) quadruply imaged quasars.
The sensitivity to low-mass halos reached by the gravitational imaging technique is highly dependent on the angular resolution of the observations [597]. At present, only a handful of systems for which Very Long Baseline Interferometry observations are available can probe the halo mass function down to ∼ 10 6 M [598]. From a sample of 20 HST-observed galaxygalaxy lensed systems (sensitivity ∼ 10 10 M ), [594] infer a limit on the warm dark matter particle mass of 1.02 keV using both detections [599] and non-detections [600,601].
At present, the relatively low number of gravitational lens systems with high-enough data quality significantly limits the constraints on dark matter from both flux-ratio anomalies and the gravitational imaging technique. Luckily, ongoing and future surveys with instruments such as Euclid, the Vera Rubin, and SKA will increase the number of known gravitational lens systems by several orders of magnitudes [602][603][604]. This data, coupled with follow-up observations with, e.g., the ELT, TMT, and JWST, is projected to deliver tight constraints on the halo and subhalo mass function [595,605].
Stellar streams are the tidally disrupted remnants of dwarf galaxies and globular clusters. Recently, the combination of astrometric, photometric, and spectroscopic observations has led to the discovery of a plethora of new streams orbiting the Milky Way [606][607][608] and intriguing structure in the density profiles of nearby streams like GD-1 [609] and ATLAS-Aliqa Uma [610]. These substructured streams-and particularly the gaps in stream density profiles and other unexpected, off-stream features-have been modeled to place constraints on the properties of individual perturbers and candidate dark matter subhalos that may have gravitationally perturbed these streams [611], as well as the population statistics of subhalo perturbers in cold and warm dark matter contexts [612,613]. Deeper photometric measurements from upcoming facilities will continue to increase the population of known streams and reveal their fine-grained density structure, with potential sensitivity to subhalos as small as ∼ 10 6 M , regardless of their baryonic content [559,614].

Machine Learning and Statistics
Over the next decade, astrophysical data relevant for the study of the nature of dark matter will increase not only in volume, but also in complexity and detail. Upcoming gammaray and radio observatories like the Cherenkov Telescope Array (CTA) [37] and the Square Kilometer Array (SKA) [615,616] will produce petabytes of data [55]. Fully exploiting this data for scientific purposes will require increasingly detailed and complex physical models, which bring along higher computational costs for simulations, as well as a larger number of uncertain parameters, including those characterizing signal and background systematics. Established algorithms for parameter inference, like Markov Chain Monte Carlo (MCMC), nested sampling [617,618], and Approximate Bayesian Computation (ABC) often require a very large number of simulation runs, which often grow significantly as the number of model parameters increases; see Fig. 4 for an illustration. In many settings, it becomes impractical to compare physical models in their full complexity and detail with the data. As a result, analyses are often limited by the inference tools rather than statistics [619,620]. Modern statistical algorithms based on deep learning and differentiable programming techniques can overcome the limitations of established techniques. Recent reviews can be found in [619,621,622], and we will provide a brief overview of promising techniques and suggest necessary developments here.
Scalable inference techniques: Stochastic variational inference (SVI) [623,624] ap-proaches inference of the posterior as an optimization problem, circumventing the need to sample from high-dimensional parameter spaces. The most commonly employed optimization target is here the evidence lower bound (ELBO) [623]. It can be conveniently optimized using stochastic gradient descent (SGD), provided the physical simulator is fully differentiable with respect to all model parameters. This approach can scale to very high dimensional inference problems. Example applications are deblending starfields [625], disentangling the components of gamma-ray emission [626,627], and strong gravitational lensing [628,629]. One challenge with SVI is that through the mode-seeking nature of the reverse KL divergence [630], it tends to underestimate the posterior variance, potentially leading to over-confident posteriors. An important theoretical development front is the construction of differentiable forward models and simulators [e.g. 631,632], which is easily admitted using modern automatic differentiation tools [633,634], or the use of differentiable surrogate models when this is not feasible [635,636]. The use of differentiable models also allow for inference via gradient-assisted Monte Carlo methods like Hamiltonian Monte Carlo, which have higher sample efficiency and scale better with parameter dimension than traditional Monte Carlo techniques.
Methods based on deep learning: The likelihood function, which is a fundamental input to most established techniques including SVI, can be extremely difficult to compute, due to required marginalization over various unobserved instrumental and physical parameters (see Fig. 4 for an illustration). Simulation-based inference (SBI) methods (see [637,638]) circumvent the evaluation of likelihoods by directly mapping observations and simulations onto summary statistics that are subsequently statistically interpreted. A classical SBI technique is ABC [639]. Various recently developed neural SBI methods use the training of deep neural networks both to generate informative summary statistics as well as performing estimation of posteriors, likelihoods [e.g. 640,641], or likelihood-to-evidence ratios [e.g. 642]. This procedure can require orders of magnitude fewer simulations than established techniques [643,644], see Fig. 4. Neural SBI has been used, for instance, for dark matter substructure inference [628,[645][646][647], dark matter indirect detection with gamma-ray data [648][649][650], and binary microlensing [651]. Fronts where still significant theoretical development is required are neural network architectures tailored to the structure of typical astrophysical data, which can significantly reduce simulation costs for simulation-based algorithms, and simulation-efficient training algorithms. This includes, for example, efforts to develop interpretable and/or explainable architectures for astrophysical data processing and the use of inference algorithms that produce statistically sound results for the purposes of scientific discovery and hypothesis testing [652]. These developments will also help foster community trust in results relying on deep learning methods, which have historically seen reluctance in adoption due to their reputation as black boxes.
In general, deep learning-based techniques enable more information to be extracted from data without requiring the use of simplified data representations and low-dimensional summary observables. Although this has the potential to significantly enhance the sensitivity of astrophysical dark matter searches, it can make typical methods more sensitive to how specific features in the data are modeled. This underscores the need for increased attention to accurately modeling the data.
Infrastructure and workforce: Realizing these goals will require development of both software and human resources. Given the potentially steep theoretical learning curve associated with the aforementioned statistical methods, the development of easy-to-use inference tools [e.g. 644,653] adopting good documentation practices with end practitioners in mind is crucial.
Finally, given the necessity of cross-disciplinary expertise in developing these methods and tools, we recommend, through appropriate hiring practices and promotion options, viable career trajectories at the intersection of statistical methods and astrophysical data analysis for new physics. A concerted effort in this direction has demonstrated success in certain fields of cosmological [654] and collider [655] data analysis. At the training stage, the existence of Ph.D. schools as well as curriculum-based learning of data analysis techniques is encouraged.

Conclusions
We are entering an era which holds the promise of resolving many of the most basic questions we have about dark matter. How was dark matter produced in the early Universe? Is dark matter a single missing piece or part of a broader dark sector? Is dark matter cold, warm, or better thought of as a wave? And above all else: what is dark matter?
A central source for optimism that the answer to these questions may be within our grasp is the upcoming advancements in instruments and observations. Yet, as we have outlined, the role of the theory community is not to simply wait for the experimental program to provide the answers to these questions, but instead to work with, extend, and optimize dark matter searches. Indeed, there is considerable work ahead for theorists to determine the behavior of dark matter and how this would manifest in our observations, and further in the development of techniques such as machine learning that could be required to confidently detect an eventual signal. While there are many challenges to overcome, as we have reviewed, there are clear paths for doing so. Viewed as a whole, there is every reason to be confident that in the coming years we will finally tease apart the mysteries of dark matter and move into a future where rather than wondering what dark matter is, we are instead asking how its particle nature modifies galaxies, cosmology, and possibly even opens a path towards understanding the broader world of physics that exists beyond the Standard Model.