Lowmass galaxies provide a powerful tool with which to investigate departures from the standard cosmological paradigm in models that suppress the abundance of small dark matter structures. One of the simplest metrics that can be used to compare different models is the abundance of satellite galaxies in the Milky Way. Viable dark matter models must produce enough substructure to host the observed number of Galactic satellites. Here we scrutinize the predictions of the neutrino Minimal Standard Model nurm MSM a wellmotivated extension of the Standard Model of particle physics in which the production of sterile neutrino dark matter is resonantly enhanced by a lepton asymmetry in the primordial plasma. This process enables the model to evade current constraints associated with nonresonantly produced dark matter. Independently of assumptions about galaxy formation physics we rule out with at least 95 per cent confidence all parameterizations of the nurm MSM with sterile neutrino rest mass Mrm s leq 1.4 rm keV. Incorporating physically motivated prescriptions of baryonic processes and modelling the effects of reionization strengthen our constraints and we exclude all nurm MSM parameterizations with Mrm s leq 4 rm keV. Unlike other literature our fiducial constraints do not rule out the putative 3.55 keV Xray line if it is indeed produced by the decay of a sterile neutrino however some of the most favoured parameter space is excluded. If the Milky Way satellite count is higher than we assume or if the Milky Way halo is less massive than Mrm MW200  8 times 1011 rm Modot we rule out the nurm MSM as the origin of the 3.55 keV excess. In contrast with other work we find that the constraints from satellite counts are substantially weaker than those reported from Xray nondetections.
['Astrophysics - Cosmology and Nongalactic Astrophysics', 'Astrophysics - Astrophysics of Galaxies']
Constraints on the properties of nuMSM dark matter using the satellite galaxies of the Milky Way
{'Constraints on the properties of 𝜈 MSMdark matter using the satellite galaxies of the Milky Way': 'Oliver Newton , 1 ★ Mark R. Lovell , 2 Carlos S. Frenk , 2 Adrian Jenkins , 2 John C. Helly , 2  \nShaun Cole 2 and Andrew J. Benson 3  \n- 1 Center for Theoretical Physics, Polish Academy of Sciences, al. Lotników 32/46 Warsaw, Poland\n- 2 Institute for Computational Cosmology, Durham University, South Road, Durham, DH1 3LE, UK  \n3 Carnegie Observatories, 813 Santa Barbara Street, Pasadena, CA 91101, U.S.A  \nAccepted XXX. Received YYY; in original form ZZZ', 'ABSTRACT': 'Low-mass galaxies provide a powerful tool with which to investigate departures from the standard cosmological paradigm in models that suppress the abundance of small dark matter structures. One of the simplest metrics that can be used to compare different models is the abundance of satellite galaxies in the Milky Way. Viable dark matter models must produce enough substructure to host the observed number of Galactic satellites. Here, we scrutinize the predictions of the neutrino Minimal Standard Model ( 𝜈 MSM), a well-motivated extension of the Standard Model of particle physics in which the production of sterile neutrino dark matter is resonantly enhanced by a lepton asymmetry in the primordial plasma. This process enables the model to evade current constraints associated with non-resonantly produced dark matter. Independently of assumptions about galaxy formation physics we rule out, with at least 95 per cent confidence, all parameterizations of the 𝜈 MSM with sterile neutrino rest mass, 𝑀 s ≤ 1 . 4 keV. Incorporating physically motivated prescriptions of baryonic processes and modelling the effects of reionization strengthen our constraints, and we exclude all 𝜈 MSM parameterizations with 𝑀 s ≤ 4 keV. Unlike other literature, our fiducial constraints do not rule out the putative 3 . 55 keV X-ray line, if it is indeed produced by the decay of a sterile neutrino; however, some of the most favoured parameter space is excluded. If the Milky Way satellite count is higher than we assume, or if the Milky Way halo is less massive than 𝑀 MW 200 = 8 × 10 11 M ⊙ , we rule out the 𝜈 MSM as the origin of the 3 . 55 keV excess. In contrast with other work, we find that the constraints from satellite counts are substantially weaker than those reported from X-ray non-detections.  \nKey words: dark matter - elementary particles - Galaxy: halo - galaxies: dwarf - galaxies:formation', '1 INTRODUCTION': "Understanding the nature of dark matter is one of the greatest outstanding challenges in contemporary physics. In the prevailing cosmological framework known as Λ + cold dark matter ( Λ CDM)mostof the dark matter is typically composed of weakly interacting, massive fundamental particles (WIMPs) exhibiting negligible thermal velocities at early times (Peebles 1982; Davis et al. 1985; Bardeen et al. 1986). Although Λ CDMsuccessfully predicts and replicates various observable features of the Universe, terrestrial direct detection (e.g. SuperCDMS Collaboration et al. 2015; CRESST Collaboration et al. 2019; DEAP Collaboration et al. 2019; PICO Collaboration et al. 2019; PandaX-4T Collaboration et al. 2021; Aalbers et al. 2023; Abe et al. 2023; Agnes et al. 2023; Aprile et al. 2023), indirect detection (Albert et al. 2017; H.E.S.S. Collaboration et al. 2018; Di Mauro et al. 2019; The Super-Kamiokande Collaboration et al. 2020; Regis et al. 2021; IceCube Collaboration et al. 2022; Di Mauro et al. 2023; McDaniel et al. 2024; Yin et al. 2024, also see Murgia 2020 for a comprehensive review of Fermi Large Area Telescope ( Fermi -LAT) searches for a 𝛾 -ray excess), and particle collider searches (ATLAS Collaboration et al. 2021b,a; CMS Collaboration et al. 2021; Tumasyan et al. 2022b,a; Aad et al. 2023; Tumasyan et al. 2023) have  \nyet to establish that weakly interacting massive dark matter particles exist. This has encouraged closer examination of alternative dark matter models that are well motivated in the particle physics context and which may be compatible with the successes achieved on astrophysical scales by Λ CDM.  \nAstrophysical observations serve as critical benchmarks against which to evaluate various aspects of these alternative models. In particular, low-mass galaxies provide valuable insights into the physics governing the onset of structure formation and its subsequent evolution. It is within this regime that departures from the prevailing paradigm are likely to manifest. Models in which the dark matter is relativistic in the early Universe and undergoes free-streaming on scales that suppress the formation of low-mass haloes may struggle to account for the observed abundance of Milky Way satellite galaxies. This has proven to be an effective means of imposing stringent constraints on several 'warm' dark matter models (WDM, e.g. Kennedy et al. 2014; Lovell et al. 2016; Enzi et al. 2021; Nadler et al. 2021; Newton et al. 2021; Dekker et al. 2022). In this study we adopt this approach, described in detail in Newton et al. (2021), to impose constraints on the parameter space of one well-motivated model: the neutrino Minimal Standard Model ( 𝜈 MSM).  \nThe 𝜈 MSM,aminimalextension of the Standard Model of particle physics, introduces three 'sterile' neutrinos alongside the family of three standard model 'active' neutrinos (Shi & Fuller 1999; Asaka &  \nShaposhnikov 2005; Laine & Shaposhnikov 2008). The sterile neutrinos are produced by the oscillation of the left-chiral active neutrinos into right-chiral states, characterized by the mixing angle, 𝜃 . The mixing of the active and sterile neutrino states produces two distinct sets of mass eigenstates wherein the masses in one set increase as the neutrinos in the other set become less massive (Minkowski 1977; Gell-Mann et al. 1979; Barbieri et al. 1980; Mohapatra & Senjanović 1980; Yanagida 1980). This elegant mechanism offers a plausible explanation for the extremely small masses of the active neutrinos of the Standard Model. The high mass eigenstates are ascribed to the sterile neutrinos: one with a mass O( keV ) that is long-lived and exhibits relativistic behaviour up to nine years after the Big Bang (Lovell 2023), and two others with masses in the GeV scale that are unstable and extremely short-lived (Adhikari et al. 2017; for a detailed review see Abazajian et al. 2012). Their decay into anti-leptons introduces an asymmetry in the lepton abundance (Fukugita & Yanagida 1986) that is partially reprocessed into a baryon asymmetry during baryogenesis. This provides a natural explanation for the observed asymmetry of matter to anti-matter in the Universe (Shaposhnikov 2008; Canetti et al. 2012, 2013). A recent comprehensive description of the production mechanisms underpinning the 𝜈 MSM in an astrophysical context is available in Lovell (2023).  \nAs active neutrinos propagate through a medium their oscillation parameters may diverge from those in vacuo because of interactions with electrons in the medium. This phenomenon is known as the Mikheyev-Smirnov-Wolfenstein (MSW) effect (Wolfenstein 1978; Mikheyev & Smirnov 1985). An analogous mechanism within the primordial plasma affects the mixing of the active and sterile neutrino states in the 𝜈 MSM. This is controlled by the size of the lepton asymmetries in different particle species, which can resonantly enhance the production of sterile neutrinos by increasing the effective mixing angle, 𝜃 e (Laine & Shaposhnikov 2008). Through this process, the 𝜈 MSM may generate sufficient numbers of long-lived keV scale sterile neutrinos to be consistent with the measured abundance of dark matter. Since the sterile neutrinos decouple from the primordial plasma before the onset of Big Bang nucleosynthesis they do not affect the primordial abundance of chemical elements. Consequently, astrophysical measurements of the elemental abundances, which have been used to constrain the properties of thermally produced WDM candidates, cannot be applied to constrain the 𝜈 MSM parameter space (Sabti et al. 2020; An et al. 2022).  \nEfforts to ascertain the nature of the dark matter intensified in response to astrophysical spectroscopic observations revealing a possible spectral line at 3 . 55 keV in systems containing a significant dark matter component (Bulbul et al. 2014; Boyarsky et al. 2014; Cappelluti et al. 2018; Hofmann & Wegg 2019). This has elicited considerable interest, given the absence of a known astrophysical origin for a spectral line at this energy. While charge exchange between sulphur ions and a neutral medium could provide a plausible baryonic explanation (Gu et al. 2015; Shah et al. 2016), it cannot be confirmed with the resolution available in current instruments. 1 Moreover, it may be unable to account for all of the observed flux in any case (Cappelluti et al. 2018). The analysis of the excess is complicated further by the choice of assumptions embedded in the modelling of  \nthe astrophysical background (Dessert et al. 2024b, also see Boyarsky et al. 2019 for a recent review).  \nA non-baryonic explanation for the excess could be provided by the dark sector. Interpreting the line as the decay product of a dark matter particle would place the rest mass of the parent particle within the range expected of a sterile neutrino candidate responsible for the entire dark matter component of the cosmic matter density. Followup observations of the spectral feature in the Galactic centre (Boyarsky et al. 2015; Perez et al. 2017; Hofmann & Wegg 2019), the Perseus cluster (Urban et al. 2015; Franse et al. 2016), and in other galaxy clusters (Bulbul et al. 2016, but see also Hofmann et al. 2016) have been carried out in attempts to improve the significance of the detection and elucidate the origin of the excess, with limited success. These plausible detections are challenging to reconcile with the non-detections reported by other observational studies of stacked galaxy spectra (Anderson et al. 2015), the Galactic centre (RiemerSørensen 2016), the Draco dwarf galaxy (Ruchayskiy et al. 2016; Sonbas et al. 2016), the X-ray background (Sekiya et al. 2016), and Perseus (Aharonian et al. 2017; Tamura et al. 2019). The origin of the 3 . 55 keV excess thus remains unclear, and this emphasizes the need for higher resolution measurements of the soft X-ray spectrum in dark matter-dominated systems. This will be achievable with the JAXA XRISM mission that launched in 2023 September that is now in the final stages of commissioning (Tashiro et al. 2018; Lovell et al. 2019a,b; Terada et al. 2021; Dessert et al. 2024a).  \nIn this work, we impose constraints on the parameter space of resonantly produced 𝜈 MSMmodels by comparing the predicted abundance of low-mass galaxies in Milky Way mass dark matter haloes with the inferred total population of Milky Way satellite galaxies. In Section 2, we summarize the method as applied to the 𝜈 MSM and calibrate it using dark matter-only N -body simulations of 𝜈 MSM models. We present our constraints on the 𝜈 MSM parameter space in Section 3 and compare them with related studies. We discuss our results in Section 4 and present concluding remarks in Section 5.", '2 METHODS': 'Newton et al. (2021) described an approach to constrain the parameter space of thermal relic WDM models by comparing predictions of the abundance of substructure in Milky Way mass dark matter haloes with the total number of satellite galaxies estimated from observations of the Milky Way halo. Here, we apply the same technique to constrain the viable combinations of the sterile neutrino rest mass, 𝑀 s , and the resonantly enhanced effective mixing angle, 𝜃 e , in the 𝜈 MSM. We make predictions of the abundance of substructure in different parameterizations of the 𝜈 MSMusing the Extended Press-Schechter formalism, first computing the 𝜈 MSM momentum distributions and matter power spectra in Section 2.1, then using Extended Press-Schechter theory to predict the subhalo counts; we describe and calibrate the Extended Press-Schechter formalism using gravity only cosmological simulations of Milky Way mass haloes in Section 2.2. In Section 2.3, we summarize the approach we take to evaluate the viability of each parameterization of the model.', '2.1 Computing the 𝜈 MSMsterile neutrino power spectra': 'The primary goal of this study is to predict the abundance of satellites of galaxies like the Milky Way as a function of the dark matter model parameters. The first step is to compute the momentum distribution of the dark matter at an energy of 1 MeV, or approximately one second after the Big Bang, which is then input to a Boltzmann solver code  \nthat calculates the linear matter power spectrum, 𝑃 ( 𝑘 ) . We compute momentumdistributions for an array of 36 values of 𝑀 s and 31 values of the lepton asymmetry parameter, 𝐿 6 , using the code of Laine & Shaposhnikov (2008). The 𝑃 ( 𝑘 ) are then calculated with a modified version of the /c.pc/a.pc/m.pc/b.pc Boltzmann solver (Lewis et al. 2000).  \nOne of the key uncertainties in our results will be the accuracy of the momentum distribution calculation. The Laine & Shaposhnikov (2008) algorithm has been superseded by that of Ghiglieri & Laine (2015), with minor modifications to the free-streaming length. A further study by Venumadhav et al. (2016) produced results that differ significantly from those of Ghiglieri & Laine (2015): the former predicts a cut-off scale some three times larger than the latter at fixed sterile neutrino mass and effective mixing angle (Lovell 2023). We will discuss this discrepancy further in Section 3.3.', '2.2 Calibrating the Extended Press-Schechter formalism with numerical simulations': "In the Extended Press-Schechter formalism the primordial matter density field is filtered using a window function to identify regions that are dense enough to collapse into virialized haloes by the present day (Press & Schechter 1974; Bond et al. 1991; Bower 1991; Lacey & Cole 1993; Parkinson et al. 2008). The abundance of structure at different mass scales can be calculated from this using the conditional halo mass function, d 𝑁 cond / d ln 𝑀 . This describes the fraction of the mass in some halo at time, 𝑡 1 , that was in collapsed progenitors of a given mass, M , at some earlier time, 𝑡 2 . The conditional halo mass function depends on the primordial linear matter power spectrum, and its functional form and the details of its derivation are provided in Schneider (2015).  \nThe subhalo mass function is computed by integrating the conditional halo mass function over the redshift-dependent overdensity collapse threshold of a given progenitor, 𝛿 𝑐 ( 𝑧 ) . It is given by  \nd 𝑁 sub d ln 𝑀 = 1 𝑁 norm ∫ ∞ 𝛿 𝑐 ( 0 ) d 𝑁 cond d ln 𝑀 d 𝛿 𝑐 , (1)  \nwhere 𝑁 sub is the number of subhaloes within a given halo and 𝑁 norm is a normalization constant. The latter term, which is a free parameter that must be calibrated using simulations, corrects the total number count for progenitor subhaloes that exist at multiple redshifts. These structures are counted more than once in this approach (Schneider 2015). We calibrate 𝑁 norm at scales above the power spectrum cut-off to minimise any dependence on the shape of the power spectrum.  \nThe form of the window function can affect the predictions of the mean number of subhaloes in a halo of a given mass. Customarily, using a spherical top-hat function in real space has been favoured because there is a simple relation between the filter volume, the filter scale radius, and the mass enclosed by the virialised object. Somewhat more importantly, this approach has also been shown to produce results that are in good agreement with those of numerical simulations of cold dark matter (Jenkins et al. 2001). However, for models of dark matter with power spectra that are damped at small scales, such as the 𝜈 MSM, this approach predicts an excess of low-mass haloes compared to numerical simulations (Benson et al. 2013; Schneider 2015). This happens because the steep cut-off in the matter power spectrum permits only negligible contributions at high wavenumbers when calculating the variance of the matter density field in the limit of small radii. Meanwhile, in the same regime, the halo mass function is dominated by the derivative of the variance, which depends strongly on the form of the filter function that is used. Consequently, in dark matter models with damped power spectra, the halo mass function calculated when adopting a real space top-hat filter function  \nbecomes independent of the power spectrum, and tends towards a constant value at small radii rather than becoming negligible. Instead adopting a smooth filter in k -space of the form,  \n˜ 𝑊 ( 𝑘 | 𝑅 ) = h 1 + ( 𝑘𝑅 ) ˆ 𝛽 i -1 , (2)  \nwhere k is the wavenumber, R is the filter scale radius, and ˆ 𝛽 > 0 is a free parameter, results in better agreement with simulations (Leo et al. 2018). To achieve this good agreement, the smooth k -space filter introduces a second free parameter, ˆ 𝑐 , that connects the filter mass, M , to R through the mean density, ¯ 𝜌 , as 𝑀 ( 𝑅 ) = 4 𝜋 ¯ 𝜌 ( ˆ 𝑐𝑅 ) 3 / 3. Both ˆ 𝛽 and ˆ 𝑐 , and the normalization constant that we discussed earlier, 𝑁 norm , must be calibrated with simulations. We showed in Newton et al. (2021) that using the Extended Press-Schechter formalism with a smooth k -space filter calibrated using a subset of haloes from the COCOsuite (Bose et al. 2016; Hellwing et al. 2016) produces results in good agreement with the simulations.  \nIn this work we calibrate the Extended Press-Schechter and filter function parameters using gravity only counterparts of the 'V2' and 'V5' sets of Local Group-analogue simulations introduced by Lovell et al. (2017). Each set comprises a cold dark matter volume and two 𝑀 s = 7 keV 𝜈 MSM volumes with sin 2 ( 2 𝜃 e ) ≃ 10 -11 . They were simulated using initial conditions that are identical to the first six simulation volumes of the APOSTLE project (Fattahi et al. 2016; Sawala et al. 2016), except for the differences introduced by changing the properties of the dark matter. All simulations adopt cosmological parameters from the WMAP seventh-year data release (Komatsu et al. 2011): 𝐻 0 = 70 . 4 kms -1 Mpc -1 , Ω M = 0 . 272 , Ω Λ = 0 . 728 , 𝑛 𝑠 = 0 . 967 , 𝜎 8 = 0 . 81.  \nSimulations of dark matter models that impose a cut-off in the primordial power spectrum are susceptible to the formation of 'spurious' haloes. These structures are produced by the artificial fragmentation of filaments caused by resolution-dependent gravitational instabilities that are generated by the discreteness of the simulation particles (Wang & White 2007; Angulo et al. 2013; Lovell et al. 2014). We identify and prune spurious haloes from the merger trees by adopting the procedure introduced by Lovell et al. (2014), which we summarize here. Numerical gravitational instabilities emerge when the gravitational potential is sampled poorly by simulation particles. The onset of these instabilities is associated with a resolution-dependent mass threshold, and haloes that do not surpass this during their evolution are likely to be spurious. In addition, spurious haloes typically form from highly aspherical Lagrangian regions in the initial conditions. These regions are characterized using the particles that compose the halo when it accumulates half of its maximum mass. The shapes of the Lagrangian regions are parametrized by 𝑠 half -max = 𝑐 / 𝑎 , where 𝑐 and 𝑎 are, respectively, the minor and major axes of the diagonalized moment of inertia tensor of the dark matter particles in Lagrangian coordinates. The threshold mass and 𝑠 half -max values for the 𝜈 MSMAPOSTLE are calculated in the appendix of Lovell et al. (2017), who found that almost all spurious haloes can be removed by requiring 𝑀 max > 5 . 35 × 10 7 h -1 M ⊙ and 𝑠 half -max > 0 . 165. We adopt these thresholds to prune the 𝜈 MSM APOSTLE simulations we consider in Fig. 1.  \nWe also correct for other resolution-dependent effects. It has been shown that configuration-space structure finders fail to identify some small haloes that clearly survive to the present day. Subhaloes that are close to the centre of the host halo are particularly likely to be overlooked. Additionally, some small haloes are artificially disrupted by tidal forces due to numerical effects such as small particle number. We recover all prematurely destroyed subhaloes by tracking their most-bound particles through the simulations and determining  \nFigure 1. The total number of dark matter subhaloes within 𝑅 200 as a function of dark matter halo mass, 𝑀 200 . Need to state which nuMSM model parameters these are. The dashed line shows the mean number of subhaloes predicted by the Extended Press-Schechter formalism when adopting the 𝜈 MSM parameters 𝑀 s = 7 . 1 and sin 2 ( 2 𝜃 e ) ≃ 10 -11 . The dark shaded region indicates the associated 68 per cent Poisson scatter, and the light shaded region gives the 68 per cent scatter modelled using a negative binomial distribution (see equation 3). The symbols represent haloes from the 7 keV 𝜈 MSM APOSTLE simulations for sin 2 ( 2 𝜃 e ) ≃ 10 -11 . Unfilled symbols are from subhalo catalogues where the prematurely destroyed subhaloes have not been recovered, and filled symbols show the same haloes after the prematurely destroyed subhaloes have been restored by following the procedure described in Section 2.2. For comparison, we plot as a dotted line the number of subhaloes predicted by the Extended Press-Schechter formalism applied to the 3 . 3 keV thermal relic WDM model from Newton et al. (2021).  \n<!-- image -->  \nwhether they survive to the present day. To do this we use the Simha &Cole(2017) merging scheme implemented in the galaxy formation model /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc (Cole et al. 2000) that recovers the entire population of destroyed substructures, including objects that were disrupted by physical processes. We identify and remove the latter from the recovered population if they satisfy one of the following criteria:  \n- · Thedynamical friction timescale, i.e. the time required for a subhalo to sink to the centre of the host halo, is less than or equal to the time elapsed since the subhalo fell below the resolution limit.\n- · The subhalo passes within the halo tidal disruption radius.  \nIn both cases, the effects of tidal stripping and of interactions between other subhaloes are ignored. This procedure is discussed in more detail in Newton et al. (2021, appendix A).  \nAfter pruning spurious haloes from the merger trees and recovering the prematurely destroyed subhaloes, we obtain excellent agreement between the Extended Press-Schechter predictions and the 𝜈 MSM APOSTLE N -body results using 𝑁 norm = 1 . 4 , ˆ 𝛽 = 4 . 2 , and ˆ 𝑐 = 3 . 9. These values are very similar to the values we found in Newton et al. (2021) when calibrating the Extended Press-Schechter formalism for thermal relic WDM. In Fig. 1 the dashed line shows the mean number of subhaloes predicted by the Extended Press-Schechter formalism within 𝑅 200 of the centre of Milky Way mass dark matter haloes  \nin the 𝑀 s = 7 . 1 keV, sin 2 ( 2 𝜃 e ) ≃ 10 -11 parameterizations of the 𝜈 MSM. Here and throughout the rest of this paper, 𝑀 200 denotes the mass enclosed within 𝑅 200 , the radius within which the mean density equals 200 times the critical density. The shaded regions represent the 68 per cent scatter in 𝑁 sub at fixed halo mass modelled using a negative binomial distribution, which encapsulates the scatter seen in numerical simulations (Boylan-Kolchin et al. 2010). This is given by  \nP ( 𝑁 | 𝑟, 𝑝 ) = Γ ( 𝑁 + 𝑟 ) Γ ( 𝑟 ) Γ ( 𝑁 + 1 ) 𝑝 𝑟 ( 1 -𝑝 ) 𝑁 , (3)  \nwhere Γ ( 𝑥 ) = ( 𝑥 -1 ) ! and N is the number of subhaloes. The parameter, 𝑝 = ⟨ 𝑁 ⟩ / σ 2 , where ⟨ 𝑁 ⟩ and σ 2 = σ 2 Poisson + σ 2 I are, respectively, the mean and the dispersion of the distribution. The parameter, 𝑟 = σ 2 Poisson / σ 2 I , describes the relative contribution of Poisson fluctuations compared to a second distribution that describes the intrinsic variability of the subhalo count within haloes of fixed mass. We model the scatter in the subhalo count by adopting σ I = 0 . 12 ⟨ 𝑁 ⟩ , which is in good agreement with that found in N -body simulations such as COCO (Newton et al. 2021). In Fig. 1, we also overlay the subhalo counts of Milky Way mass haloes in the APOSTLE 𝑀 s = 7 keV, sin 2 ( 2 𝜃 e ) ≃ 10 -11 simulations, both before and after restoring the population of prematurely destroyed subhaloes to the catalogues. The Extended Press-Schechter predictions and model of the scatter in the subhalo counts are in excellent agreement with the latter.  \nAt fixed mass, the Milky Way mass haloes that form in the parameterization of the 𝜈 MSMshowninFig.1(dashedline) are predicted to host more dark matter subhaloes than they would in the 3 . 3 keV thermal relic WDM model (dotted line). Consequently, larger fractions of the 𝜈 MSM systems host at least as many subhaloes as the inferred total Milky Way satellite galaxy population. As Newton et al. (2021) were unable to place constraints on the 3 . 3 keV thermal relic model using this approach, we will be similarly unable to constrain the sin 2 ( 2 𝜃 e ) ≃ 10 -11 parameterizations of the 𝑀 s = 7 . 1 keV 𝜈 MSMmodel. We will discuss our full results in Section 3.1.", '2.3 Calculating model acceptance probability': 'Weobtainourconstraints on the 𝜈 MSMparameterspacebyfollowing the approach described in Newton et al. (2021), which we summarize briefly here. First, we calculate the fraction, 𝑓 v , of Milky Way mass systems that have at least as many satellites 𝑁 𝜈 MSM sat , as the total number of Milky Way satellite galaxies, 𝑁 MW sat , for a given set of 𝜈 MSMparameters. This is given by  \n𝑓 v = ∫ ∞ 0 d 𝑁 MW sat 𝑝 MW GLYPH<16> 𝑁 MW sat GLYPH<17> × ∫ ∞ 𝑁 MW sat d 𝑁 𝜈 MSM sat 𝑝 𝜈 MSM GLYPH<16> 𝑁 𝜈 MSM sat GLYPH<17> , (4)  \nwhere 𝑝 MW and 𝑝 𝜈 MSM are, respectively, the probability distributions of the Milky Way satellite count, and of the number of satellites predicted to be hosted by Milky Way mass haloes at 𝑧 = 0 in a given parameterization of the 𝜈 MSM model. Then, we obtain constraints by computing the intersection of the cumulative distribution of 𝑓 v as a function of Milky Way halo mass at fixed 𝑀 s and 𝜃 e with a five per cent rejection threshold. Parameterizations of the 𝜈 MSM at and below this threshold are ruled out.  \nOur calculation of 𝑓 v in equation (4) accounts for the scatter in the number of dark matter subhaloes at fixed host halo mass and the uncertainty in the size of the total Milky Way satellite galaxy  \npopulation. As both of these quantities depend on the assumed mass of the Galactic halo, the viable fraction is also highly sensitive to it. Estimates of the mass of the Galactic halo are still quite uncertain; indeed, it has been revised downwards in recent years (see e.g. Bird et al. 2022). To account for this uncertainty we calculate 𝑓 v for several Milky Way halo masses spanning the most likely range determined from kinematic analyses of the Galactic potential. Unless stated otherwise the constraints we present in Section 3 are marginalized over the distribution of Milky Way halo masses from Callingham et al. (2019).', '3 CONSTRAINTS ON THE PROPERTIES OF STERILE NEUTRINOS IN THE 𝜈 MSM': 'The constraints on the parameter space of the 𝜈 MSMdepend on the assumed mass of the Milky Way halo, the rest mass of the sterile neutrino dark matter candidate, 𝑀 s , and the resonantly enhanced effective mixing angle of the sterile neutrino, 𝜃 e . We consider Milky Way halo masses in the likely range 𝑀 MW 200 = ( 0 . 5-2 . 0 ) × 10 12 M ⊙ (Callingham et al. 2019; Cautun et al. 2020), and values of sin 2 ( 2 𝜃 e ) between the dark matter under- and over-production bounds if the sterile neutrino comprises all of the dark matter. In this and subsequent sections we assume the Planck Collaboration et al. (2020) cosmology: 𝐻 0 = 67 . 4 ± 0 . 5 kms -1 Mpc -1 , Ω M = 0 . 315 ± 0 . 007 , Ω Λ = 0 . 685 ± 0 . 007 , 𝑛 𝑠 = 0 . 965 ± 0 . 004 , 𝜎 8 = 0 . 811 ± 0 . 006. For each combination of parameters, we compute the fraction of 𝜈 MSM systems that produce at least as many galaxies as the number of Milky Way satellite galaxies inferred from observations using the methodology presented in Newton et al. (2018), as implemented in Newton & Cautun (2018).  \nWe obtain our most robust result, presented in Section 3.1, under the premise that all dark matter subhaloes host a visible galaxy; therefore, we make no assumptions at all about galaxy formation processes. This maximises 𝑓 v and produces extremely robust lower limits on the viable 𝜈 MSM parameter space. In Section 3.2, we use the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc semi-analytic model of galaxy formation to estimate the number of satellite galaxies in each parameterization of the 𝜈 MSM. This allows us to account for astrophysically relevant processes, such as reionization, that can affect the properties of the low-mass galaxy population. In Section 3.3, we consider how our results are affected when we use different estimates of 𝑁 MW sat from that reported by Newton et al. (2018) and compare these results with others in recent literature.', '3.1 Constraints from structure formation considerations': 'In Fig. 2, we present constraints on the parameter space of 𝜈 MSM models obtained by comparing the abundance of dark matter substructure predicted using the Extended Press-Schechter formalism with the Newton et al. (2018) number count of Milky Way satellite galaxies. Parameter combinations of 𝑀 s and 𝜃 e within the shaded regions are ruled out with at least 95 per cent confidence when marginalizing over the uncertainties in the Milky Way halo mass from Callingham et al. (2019). We show recent astrophysical limits from X-ray modelling of M31 (Horiuchi et al. 2014), Galactic bulge non-detections (Roach et al. 2020), and blank sky non-detections (Dessert et al. 2020b; Foster et al. 2021; Roach et al. 2023) using a dark grey curve and shaded region. The viable parameter space decreases significantly at and below 𝑀 s = 7 keV, and we rule out all parameterizations of the 𝜈 MSM with 𝑀 s ≤ 1 . 2 keV. This result is  \nFigure 2. Constraints on the 𝜈 MSM parameter space obtained assuming all subhaloes host a galaxy. Our fiducial structure count result (thick solid curve) is marginalized over the uncertainties in the Callingham et al. (2019) estimate of the Milky Way halo mass. Parameter combinations in the shaded region to the left of this envelope are ruled out with at least 95 per cent confidence. We also show the constraint envelopes we obtain assuming Milky Way halo masses of 𝑀 MW 200 = 0 . 5 × 10 12 M ⊙ (dotted curve) and 𝑀 MW 200 = 2 . 0 × 10 12 M ⊙ (dashed curve). The shaded regions in the upper and lower portions of the parameter space delimit parameter combinations for which the dark matter is over- or under-produced compared with current observational estimates of the composition of the Universe if resonantly produced sterile neutrinos compose all of the dark matter (Asaka et al. 2005; Schneider 2016). The lower limits on the effective mixing angle obtained from X-ray non-detections are shown using the thin solid curve and shading (Neronov et al. 2016; Ng et al. 2019; Roach et al. 2020, 2023). The solid error bar indicates the parameter range favoured by X-ray flux measurements if the dark matter consists of 𝑀 s = 7 . 1 keV decaying sterile neutrinos (Boyarsky et al. 2014; Bulbul et al. 2014). The dash-dot error bar shows the preferred range determined from X-ray line detections in the Galactic bulge (Hofmann &Wegg 2019).  \n<!-- image -->  \nvalid for all Milky Way halo masses we consider up to and including the largest mass of 𝑀 MW 200 = 2 × 10 12 M ⊙ (dashed curve), which produces the least constraining results. When marginalizing over the uncertainty in the Callingham et al. (2019) estimate of the Milky Way halo mass, we rule out all parameterizations of the 𝜈 MSMwith 𝑀 s ≤ 1 . 4 keV. These results are consistent with constraints obtained by theoretical and observational analyses of the phase space density of dark matter within Milky Way dwarf satellite galaxies (Boyarsky et al. 2009).  \nAs we discussed in Section 1, the 3 . 55 keV excess observed emanating from several dark matter-rich astrophysical systems could potentially be evidence of decaying 𝜈 MSM sterile neutrinos with rest masses, 𝑀 s = 7 . 1 keV. If all of the dark matter in the Universe is composed of these 𝜈 MSM neutrinos, X-ray flux observations favour a narrow window of viable parameter space (Boyarsky et al. 2014, 2015; Bulbul et al. 2014; Iakubovskyi et al. 2015; Ruchayskiy et al. 2016). We indicate this in Fig. 2 using a solid error bar. We also show the preferred region of parameter space determined from X-ray line detections in the Galactic bulge (Hofmann & Wegg 2019)  \nusing a dash-dot error bar. Our fiducial results, which we emphasize make no assumptions about galaxy formation physics, constrain the 𝑀 s = 7 . 1 keV 𝜈 MSM parameter space only at small values of 𝜃 e because the amount of substructure produced by other parameterizations is too great. If the mass of the Milky Way is at the lower end of estimates (dotted curve), as suggested by recent studies (e.g. Bird et al. 2022; Necib & Lin 2022; Zhou et al. 2023; Roche et al. 2024), more of the parameter space is excluded, including some parameterizations favoured by Hofmann & Wegg (2019). The limits from X-ray non-detections already exclude the 𝜈 MSMas a viable explanation for the origin of the 3 . 55 keV excess, although they are subject to modelling uncertainty (Dessert et al. 2020b, though see Abazajian 2020; Boyarsky et al. 2020, and Dessert et al. 2020a). Combining this exclusion region with our fiducial satellite count constraints, all parameterizations of the 𝜈 MSM with 𝑀 s = 7 . 1 keV are excluded. This conservative and highly robust result can be extended to more of the parameter space by adopting a physically motivated prescription to populate some of the subhaloes with galaxies. This reduces the subhalo occupation fraction and, consequently, strengthens the constraints on the viable parameter space of favoured models. We explore this in Section 3.2.  \nOur analysis depends on the total number of Milky Way satellite galaxies, which is uncertain. In Appendix A, we present the structure count constraints we obtain when assuming different values of the size of the total Milky Way satellite galaxy population that have been proposed in recent years. We will discuss this within the context of applying a galaxy formation model in Section 3.3.', '3.2 Modelling galaxy formation processes': 'The satellite count constraints presented in Section 3.1 are highly robust and conservative lower limits on the 𝜈 MSMparameter space. However, they ignore the effects of baryon physics on the formation and evolution of the stellar component in low-mass haloes. Although the details of these processes are not fully understood they are expected to have a significant effect on the properties of the satellite galaxy population of the Milky Way. To understand how baryonic physics affects the constraints we explore the parameter space of these processes using the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc semi-analytic model of galaxy formation (Cole et al. 1994, 2000; Lacey et al. 2016). The parameters of this state-of-the-art model are calibrated to reproduce a variety of observational results, including the 𝑧 = 0 luminosity functions of galaxies in the local Universe, the H /i.pc mass function, the distributions of galaxy morphologies and sizes, the 𝑧 = 0 Tully-Fisher relation, and the number counts of sub-millimeter galaxies. The complete set of observations used to calibrate the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc model parameters is described in Lacey et al. (2016, section 4 . 2). Although /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc has been calibrated assuming the Λ CDM model, it is not necessary to recalibrate it for each parameterization of the 𝜈 MSM. This is because /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc is calibrated using massive galaxies far above the mass scales relevant here, and we do not expect that changes in the small-scale dark matter distribution will affect the outcome.  \nAmong the various astrophysical processes that can be explored with /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc, the effect of the reionization of the Universe on lowmass galaxies is of particular interest for this study. During reionization the intergalactic medium is heated by the UV background radiation and is prevented from cooling into low-mass haloes. This inhibits the accretion of fresh supplies of cold gas that are necessary for star formation, thereby curtailing further growth of the stellar component. In /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc this is modelled by assuming that after reionization finishes at redshift, 𝑧 = 𝑧 reion , gas cooling does not take place in haloes with circular velocities, 𝑉 circ < 𝑉 cut . In general, more  \nFigure 3. Constraints on the 𝜈 MSM parameter space obtained assuming our fiducial model of reionization ( 𝑧 reion = 7 and 𝑉 cut = 30 kms -1 ) in the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc galaxy formation model. Our main result marginalizes over the uncertainties in the Callingham et al. (2019) estimate of the Milky Way halo mass (thick solid curve). Parameter combinations in the shaded region to the left of this envelope are ruled out with at least 95 per cent confidence. We also show the constraint envelopes we obtain assuming Milky Way halo masses of 𝑀 MW 200 = 0 . 5 × 10 12 M ⊙ (dotted curve) and 𝑀 MW 200 = 2 . 0 × 10 12 M ⊙ (dashed curve). As in Fig. 2, the parameter space is delimited by upper and lower shaded regions in which the dark matter is over- or under-produced compared with current observational estimates of the composition of the Universe (Asaka et al. 2005; Schneider 2016). Lower limits on the effective mixing angle obtained from X-ray non-detections are shown using the thin solid curve and shading (Neronov et al. 2016; Ng et al. 2019; Roach et al. 2020, 2023). The solid and dot-dashed error bars indicate the parameter ranges favoured by X-ray flux measurements from different systems if all of the dark matter comprises 𝑀 s = 7 . 1 keV sterile neutrinos (Boyarsky et al. 2014; Bulbul et al. 2014; Hofmann & Wegg 2019).  \n<!-- image -->  \nlow-luminosity galaxies are able to form if reionization finishes later, and lower values of the circular velocity cooling threshold, 𝑉 cut , permit those galaxies to become brighter because they are able to continue accreting gas for longer (e.g. Bose et al. 2018). Consequently, we expect to find the most stringent constraints on the 𝜈 MSMin the parameterizations of /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc where reionization finishes at higher redshifts and 𝑉 cut is large. Following the approach of Newton et al. (2021), we explore the parameter space of reionization spanning all combinations of 𝑧 reion = [ 6 , 7 , 8 ] and 𝑉 cut = [ 25 , 30 , 35 ] kms -1 .  \nTo obtain constraints on the 𝜈 MSM, we generate predictions of the number of satellite galaxies brighter than 𝑀 V = 0 in Milky Way mass haloes and compare them with estimates of the total Milky Way satellite complement based on observations (Newton et al. 2018). We adopt the Lacey et al. (2016) version of the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc model, which provides a reasonable description of a large variety of baryonic feedback and evolutionary processes, including those that are particularly relevant for low-mass galaxies such as star formation; gas heating and cooling; supernovae feedback; and the chemical enrichment of stars and gas; and vary the reionization parameters as described above. We generate Monte Carlo merger trees using an improved version of the Parkinson et al. (2008) algorithm, which we calibrate to match as  \nFigure 4. Constraints on the 𝜈 MSM parameter space obtained assuming different parameterizations of reionization in the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc galaxy formation model when marginalizing over the uncertainties in estimates of the Milky Way halo mass from Callingham et al. (2019). Parameter combinations to the left of the envelopes are ruled out with at least 95 per cent confidence. In both panels our fiducial result is indicated by the thick solid curve and shaded region. Left panel : The circular velocity cooling threshold is fixed at 𝑉 cut = 30 kms -1 and we compute constraint envelopes assuming reionization finishes at 𝑧 reion = 6 (dotted curve), 𝑧 reion = 7 (solid curve), and 𝑧 reion = 8 (dashed curve). Parameterizations in which reionization finishes earlier produce more stringent constraints on the 𝜈 MSM parameter space. Right panel : The redshift at which reionization finishes is fixed at 𝑧 reion = 7 and we compute constraint envelopes assuming circular velocity cooling thresholds of 𝑉 cut = 25 kms -1 (dotted curve), 30 km s -1 (solid curve), and 35 km s -1 (dashed curve). Higher cooling thresholds produce more stringent constraints on the 𝜈 MSMparameter space.  \n<!-- image -->  \nclosely as possible the merger trees obtained from recent cold dark matter N -body simulations (A. J. Benson et al., in preparation). In Newton et al. (2021) we noted that the match to the N -body results at the faint end is inexact because the algorithm lacks sufficient free parameters to model both the high- and low-mass ends of the mass function. Here, we try to address this by modifying the Parkinson et al. (2008) algorithm to include an explicit dependence on the linear growth-rate factor and the logarithmic slope of the 𝜎 ( 𝑀 ) relation, which describes the variance of the matter density field as a function of mass. This introduces three new free parameters that we calibrate using N -body simulations. We provide full details of our modifications to the merger tree algorithm in Appendix B. Our alterations significantly improve the agreement between the Monte Carlo satellite galaxy luminosity functions and the N -body results; however, they are still slightly discrepant. Therefore, we adopt the approach described in Newton et al. (2021) and introduce an empirical correction to remap the 𝑀 V values of the Monte Carlo satellite galaxies to new values so that the resulting luminosity functions are consistent with those from the N -body simulations. Although we cannot eliminate the need to apply an empirical correction, the modifications we make to the Parkinson et al. (2008) algorithm reduce the size of the correction we must introduce when compared with that required in Newton et al. (2021). These are typically 10 per cent.  \nTo generate the constraints on the 𝜈 MSM we produce satellite galaxy luminosity functions for 500 realizations of each Milky Way halo mass that we consider and use them to calculate the model acceptance distributions in the manner we described in Section 2.3. Parameterizations of the model for which no more than 5 per cent of the realizations have at least as many satellite galaxies as the obser-  \nions are ruled out. In Fig. 3, the shaded exclusion region shows the constraints we obtain on the 𝜈 MSM parameter space using our fiducial model of reionization with 𝑧 reion = 7 and 𝑉 cut = 30 kms -1 when marginalizing over the uncertainties in the Milky Way halo mass from Callingham et al. (2019). Parameter combinations in this region are ruled out with at least 95 per cent confidence. Once again, we show the constraint envelopes for our fiducial parameterization of reionization assuming Milky Way halo masses at the lower (dotted curve) and upper (dashed curve) ends of the likely range. Independently of Milky Way halo mass, all parameterizations of the 𝜈 MSM with 𝑀 s ≤ 2 . 6 keV are excluded by the Milky Way satellite count constraints with at least 95 per cent confidence (dashed curve). When marginalizing over possible Milky Way halo masses, almost all of the 𝑀 s = 7 . 1 keV 𝜈 MSM parameter space favoured by Hofmann & Wegg (2019) as an explanation for the 3 . 55 keV line, and some of the parameter space favoured by Bulbul et al. (2014) and Boyarsky et al. (2014), is excluded. These constraints complement the reported X-ray non-detection lower limits in closing the 𝑀 s = 7 . 1 keV 𝜈 MSM parameter space. If the Milky Way mass is 𝑀 MW 200 = 0 . 5 × 10 12 M ⊙ , we rule out all parameterizations of the 𝜈 MSM with 𝑀 s ≤ 14 keV. All parameterizations of the 𝜈 MSMwith 𝑀 s = 7 . 1 keV are ruled out when 𝑀 MW 200 ≤ 8 × 10 11 M ⊙ (not shown).  \nWe explore how changing the parameterization of reionization affects the constraints in Fig. 4. The left panel shows the effect of varying the redshift at which reionization finishes while holding the circular velocity cooling threshold fixed at 𝑉 cut = 30 kms -1 . A late-finishing epoch of reionization (dotted curve) enables more galaxies to form in low-mass dark matter haloes and weakens the constraints on the 𝜈 MSMparameter space. This is most pronounced  \nin parameterizations of the 𝜈 MSM with 𝑀 s ≥ 10 keV, where the parameter space is almost completely unconstrained. In the scenario in which reionization finishes earlier (dashed curve) than our fiducial choice the viable parameter space closes slightly, and all parameterizations below 𝑀 s = 5 keV are excluded. The right panel shows the constraints for scenarios in which reionization finishes at 𝑧 reion = 7. Lowering the circular velocity cooling threshold permits a larger fraction of the low-mass galaxy population to accrete fresh supplies of cold gas after the end of reionization. This is then available to fuel further star formation that causes the faintest galaxies to become brighter than 𝑀 V = 0 by 𝑧 = 0. The brightened galaxies populate the satellite galaxy luminosity function of the halo and weaken the constraints on the 𝜈 MSM parameter space (dotted curve). Raising 𝑉 cut has the opposite effect, and a larger fraction of the low-mass galaxies may only form stars from the reservoir of cold gas accreted prior to reionization. This limits how bright they become by 𝑧 = 0, and fewer of them become brighter than 𝑀 V = 0 to populate the luminosity function. Consequently the constraints become stronger in this case (dashed curve). Similarly to the left panel, parameterizations of the 𝜈 MSM with sterile neutrino rest masses larger than 𝑀 s = 10 keV are particularly sensitive to changes to the circular velocity cooling threshold.', '3.3 Dependence on the total number of Milky Way satellite galaxies': 'Recent improvements in the quantity and quality of observations of the Milky Way satellite system have motivated renewed efforts to estimate the size of the total population, 𝑁 MW sat . This number remains highly uncertain because current observations are insufficiently deep to probe as far as the virial radius of the Milky Way halo, and they do not cover the entire sky. As analyses such as ours depend on 𝑁 MW sat it is important to understand how the results are affected when assuming different values. This is especially important when comparing our results with other studies that implement the same technique to impose constraints on the parameter space of the 𝜈 MSM (Nadler et al. 2021). In Fig. 5, we present the constraints derived from our fiducial /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc model ( 𝑧 reion = 7 and 𝑉 cut = 30 kms -1 ) while adopting three different estimates of the size of the Milky Way satellite galaxy population from recent literature (Newton et al. 2018; Nadler et al. 2019, 2020). As in previous sections we have marginalized over the uncertainties in the Callingham et al. (2019) estimate of the Milky Way halo mass. Our constraints become more restrictive as the assumed value of 𝑁 MW sat increases because the 𝜈 MSM models must produce more substructure containing galaxies to remain viable. Adopting 𝑁 MW sat = 134 (dashed curve) as suggested by Nadler et al. (2019), a value that is slightly larger than the Newton et al. (2018) estimate we used for our fiducial results in Sections 3.1 and 3.2, yields little change in the constraints for parameterizations of the 𝜈 MSMwith 𝑀 s = 7 . 1 keV. If the total satellite population is closer to 𝑁 MW sat = 183 (dotted curve), as suggested by Nadler et al. (2020), we find that the viable parameter space of the 𝜈 MSMdecreases significantly, in agreement with Nadler et al. (2021). Moreover, the currently favoured parameter space of plausible 𝜈 MSMmodels corresponding to the putative 3 . 55 keV X-ray excess closes completely.  \nAlongside our constraints in Fig. 5, we plot the constraint envelopes reported by Nadler et al. (2021), who adopt the estimates of the size of the Milky Way satellite galaxy population mentioned above. At fixed 𝑁 MW sat their constraints are more stringent than ours. A comparison between our results is complicated by various factors. To begin with, Nadler et al. (2021) do not calculate constraints on the  \nFigure 5. Constraints on the 𝜈 MSMparameter space from analyses of Milky Way satellite counts. The shaded regions in the upper and lower portions of the parameter space delimit parameter combinations for which the dark matter is over- or under-produced compared with current observational estimates of the composition of the Universe if resonantly produced sterile neutrinos compose all of the dark matter (Asaka et al. 2005; Schneider 2016). Our results adopt our fiducial model of reionization ( 𝑧 reion = 7 and 𝑉 cut = 30 kms -1 ) in the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc galaxy formation model marginalized over the uncertainties in estimates of the Milky Way halo mass from Callingham et al. (2019). Wecompute constraint envelopes assuming a total Milky Way satellite galaxy population with a mean value of 𝑁 MW sat = 124 (thick solid curve), 134 (dashed curve), 183 (dotted curve), and 270 (dot-dashed curve). Parameter combinations within the envelopes are ruled out with at least 95 per cent confidence. For comparison we also show the constraint envelopes from Nadler et al. (2021).  \n<!-- image -->  \n𝜈 MSM directly but rather infer them indirectly from constraints on the mass of a thermal relic neutrino.  \nThere are several other differences between the two analyses. Nadler et al. (2021) adopt the Venumadhav et al. (2016) computation of the momentum distributions, which results in systematically larger cut-offs in the primordial power spectra at fixed sterile neutrino parameters than the Laine & Shaposhnikov (2008) algorithm that we use (Lovell 2023). Nadler et al. make use of subhalo abundance matching and several analytic prescriptions to populate a subset of subhaloes with galaxies in simulations of thermal relic neutrinos. They compare this with their estimate of 𝑁 MW sat , which they determine using a forward modelling framework that accounts for the tidal influence of the LMC on the Milky Way halo, to obtain a constraint on the thermal relic WDM parameter space. They express this result in terms of the half-mode mass, 𝑀 hm , which is the mass corresponding to the half-mode scale, 𝑘 hm , at which the thermal relic power spectrum, 𝑃 WDM , is suppressed by a factor of four relative to the cold dark matter power spectrum, 𝑃 CDM (Viel et al. 2013; Bose et al. 2016), in other words, the scale at which  \n√︄ 𝑃 WDM ( 𝑘 hm ) 𝑃 CDM ( 𝑘 hm ) = 0 . 5 . (5)  \nNadler et al. (2021) remap their thermal relic constraint to the 𝜈 MSM parameter space by comparing 𝑀 hm computed for each parameter-  \nn of the 𝜈 MSM to that of the coldest permitted thermal relic power spectrum, and rule out those parameterizations where 𝑀 hm is greater than that of the thermal relic model.  \nThe indirect approach followed by Nadler et al. (2021) is not directly comparable with our method. However, the stringency of their constraints suggests that their combination of subhalo abundance matching and analytic prescriptions suppresses the formation of lowmass galaxies more aggressively than the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc prescription we adopt. As we discuss in detail in Section 4, this may be connected to the fixed mass threshold Nadler et al. employ to determine whether a galaxy is able to form in a given (sub)halo. Somewhat surprisingly, even when we adopt extremely high values of 𝑁 MW sat we are unable to replicate the stringency of their constraints obtained using smaller values. Some of this discrepancy could be explained by our choice to use the Newton et al. (2018) probability distribution to describe the uncertainty in 𝑁 MW sat , which is broader than that of Nadler et al. (2020) on which the constraints reported in Nadler et al. (2021) are based. Broader uncertainties on 𝑁 MW sat can weaken the constraints by up to 15 per cent (see Newton et al. 2021, section 2 . 4). However, the size of the discrepancy between the results of the two analyses suggests that this is a subdominant effect.  \nWhile a detailed investigation of the sources of disparity between the two studies is beyond the scope of this work, as we will see in Section 4, the methodological disparities in populating subhaloes with galaxies likely play a significant role alongside uncertainties in the particle physics calculation. The results we present here support the view that our approach produces the most conservative constraints on the parameter spaces of alternative dark matter models relative to the standard Λ CDM paradigm. Efforts such as these complement various observation-based analyses such as measurements of the Xray spectra of dark matter-dominated systems, the Lyman 𝛼 forest (Iršič et al. 2024), and gravitational lensing (Enzi et al. 2021; Zelko et al. 2022).', '4 DISCUSSION': "The 𝜈 MSM is an appealing, physically motivated alternative to the Λ CDM paradigm. It provides a plausible explanation for the extremely small masses of the active neutrinos in the Standard Model of particle physics and introduces a mechanism to produce the matteranti-matter asymmetry in the Universe. Through these processes it can produce a dark matter particle candidate in a plausible range of parameter space that appears to be consistent with astrophysical constraints. In certain parameterizations of the 𝜈 MSM the sterile neutrino dark matter behaves relativistically at early times and suppresses the formation of dark matter structure at and below the scale of dwarf galaxies. Such departures from the predictions of Λ CDM could manifest as a suppression of the number of low-mass satellite galaxies around Milky Way mass hosts. This provides a useful metric to evaluate the viability of different parameterizations of the model (this approach is described in detail in Newton et al. 2021). We exploit this to constrain the two-dimensional parameter space of the 𝜈 MSM.  \nThe choices that are adopted to populate dark matter substructure with galaxies play a particularly important role in setting the stringency of the constraints. Scenarios in which the subhalo occupation fraction is high exclude less of the model parameter space because more parameterizations are consistent with the inferred total number of Milky Way satellite galaxies. We rule out a large fraction of the parameter space below 𝑀 s = 7 . 1 keV when assuming that all dark matter substructure hosts galaxies (see Fig. 2); an extreme as-  \nption that is likely to be physically unrealistic. The advantage of this choice is that the results provide robust lower limits that are independent of assumptions about the details of galaxy formation physics. The rest of the parameter space is largely unconstrained; in this regime the most compelling constraints come from observational analyses of X-ray non-detections, although these are subject to modelling uncertainty (Horiuchi et al. 2014; Neronov et al. 2016; Ng et al. 2019; Dessert et al. 2020b; Roach et al. 2020; Foster et al. 2021; Roach et al. 2023).  \nAdopting physically motivated galaxy formation prescriptions, as we do to obtain the results in Fig. 3, lowers the subhalo occupation fraction and rules out large swathes of parameter space below 𝑀 s = 10 keV, and all parameterizations at and below 𝑀 s = 4 keV. In both cases, our fiducial results suggest that Milky Way satellite counts alone are insufficient to exclude all parameterizations of the 𝜈 MSM with 𝑀 s = 7 . 1 keV, which is favoured as an explanation for the 3 . 55 keV excess detected in the X-ray spectra of some galaxies and clusters. In both scenarios, the shapes of the exclusion regions are affected by the non-monotonic evolution of the momentum distribution of the dark matter as a function of 𝜃 e at fixed 𝑀 s . The coldest momentum distributions are typically found at the logarithmic midpoint between the dark matter under- and over-production bounds.  \nThe constraints could be strengthened further if the mass of the Milky Way halo is lower than currently assumed. We determine our fiducial results by marginalizing over the uncertainties in the Callingham et al. (2019) estimate of 𝑀 MW 200 = 1 . 17 + 0 . 21 -0 . 15 × 10 12 M ⊙ . The closest recent comparable analysis to ours assumes 𝑀 MW 200 = 1 . 4 × 10 12 M ⊙ (Nadler et al. 2021). A synthesis of the mass estimates spanning the past two decades favours 𝑀 MW 200 = 8 . 5 × 10 11 M ⊙ (Bird et al. 2022), although this is 50 per cent more massive than the results derived from recent analyses using the extensive Gaia data set of six-dimensional phase space information about the Milky Way stellar halo. Computing our constraints using 𝑀 MW 200 = 0 . 5 × 10 12 M ⊙ improves the stringency of our results, which is especially dramatic when we adopt the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc prescription to populate the subhaloes with galaxies. In this case, we rule out all parameterizations of the 𝜈 MSM at and below 𝑀 s = 14 keV with at least 95 per cent confidence and strongly exclude 𝜈 MSM sterile neutrinos as the sole explanation for the 3 . 55 keV excess. The latter finding holds for nearly all parameterizations of reionization except in the late-finishing scenarios with 𝑧 reion = 6 and 𝑉 cut = 25 kms -1 (not shown). A recent analysis of the Lyman 𝛽 forest suggests that reionization may finish as late as 𝑧 reion = 5 . 5 (Zhu et al. 2022). Our results (see Fig. 4) suggest that the 𝑀 s = 7 . 1 keV parameter space is unlikely to be fully closed by Milky Way satellite counts in this case; however, in combination with the lower limits from X-ray nondetections (also shown in Figs 2 and 3), this parameter space is strongly disfavoured.  \nThe final important consideration folded into satellite count analyses such as ours is the estimated size of the total satellite galaxy population of the Milky Way. As we show in Appendix A and Fig. 5, this can affect the stringency of the constraints significantly. The true size of the Milky Way satellite galaxy population is unknown because most of the virial volume has not been surveyed. In recent years this figure has been updated several times in response to the influx of new data from several wide-area survey programmes and the application of different modelling techniques to account for incompleteness. Now, several estimates exist. When we adopt the value of 𝑁 MW sat = 183 favoured by Nadler et al. (2020), which disagrees with our earlier estimate of 𝑁 MW sat = 124 + 40 -27 (Newton et al. 2018), we rule out all parameterizations at and below 𝑀 s = 7 . 1 keV using our fidu-  \ny Way mass and parameterization of reionization. This is significantly less restrictive than the constraints obtained by Nadler et al. (2021), who rule out all parameterizations of the 𝜈 MSMbelow 𝑀 s = 17 . 5 keV when assuming this number of satellite galaxies. We are unable to account for the discrepancy between our results and theirs by changing the number of satellite galaxies alone.  \nSeveral components of the respective analyses could contribute to a discrepancy of this magnitude: (i) the width of the 𝑁 MW sat uncertainty distribution, (ii) artificially disrupted subhaloes in simulations, (iii) the assumed mass of the Milky Way halo, and (iv) how the subhaloes are populated with galaxies. We now discuss and compare each of these in turn. Newton et al. (2021) demonstrated that the width of the distribution describing the uncertainty in the estimate of 𝑁 MW sat affects the stringency of the constraints, and that broader distributions can weaken them by as much as 15 per cent. We use the Newton et al. (2018) probability distribution to model this uncertainty, which is slightly broader than the one used by Nadler et al. (2021); however, the discrepancy is marginal so this is likely to be a subdominant effect.  \nNewton et al. (2021) also showed that the artificial disruption of low-mass subhaloes caused by limited numerical resolution affects the constraints. Failing to account for this population, some of which could in principle host galaxies, artificially strengthens the constraints in analyses of this type. Nadler et al. (2021) attempt to account for such artificially destroyed subhaloes, finding that they comprise only 10 per cent of the whole subhalo population in gravity only zoom-in simulations of Milky Way mass haloes with a particle mass resolution of 3 × 10 5 h -1 M ⊙ (Mao et al. 2015). This fraction is significantly less than the sevenfold increase in the number of subhaloes suggested by Newton et al. (2021, in particular see fig. 7) for similarly massive host haloes drawn from the COCO simulations, which have a comparable mass resolution. This discrepancy is a potentially significant factor in the disparity between our results and theirs. It also highlights the considerable uncertainty that remains about the proper treatment of prematurely destroyed subhaloes. For example, Benson & Du (2022) found that artificial disruption in the Caterpillar simulations reduces the host's subhalo population by 10 -20 per cent. Most of the destruction is concentrated in the centre of the host, where as much as half of the subhaloes are destroyed by numerical instability. These results are consistent with the estimates used above. A related issue concerns the modelling of the stellar component of prematurely destroyed subhaloes, which is poorly understood at present. It is unclear what fraction of the prematurely destroyed subhaloes would host an observable galaxy, and to what extent their stars would be affected by tidal stripping. It is thus an important systematic in satellite number count analyses.  \nWe mentioned previously that the Milky Way halo mass adopted by Nadler et al. (2021) is 20 per cent larger than our fiducial assumption. As choosing a more massive Milky Way halo mass weakens their constraints slightly, we dismiss this as having only a negligible contribution to the discrepancy between the two sets of results. Finally, as we showed in Sections 3.1 and 3.2, the procedure used to populate dark matter subhaloes with galaxies affects the constraints. Nadler et al. (2021) use subhalo abundance matching to determine what fraction of the subhalo population hosts galaxies, whereas we adopt physically motivated semi-analytic prescriptions implemented in /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc.  \nSubhalo abundance matching in the relevant mass regime has been shown to be inconsistent with the results of high resolution simulations by Sawala et al. (2015). By contrast, /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc provides a more nuanced physical portrayal of the complexities of galaxy formation on small scales, in particular, processes such as star form-  \nation and feedback; reionization; and gas heating and cooling; albeit with some simplifying assumptions. Although /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc is not calibrated to reproduce observational results at the scales of interest here, it nonetheless performs well when compared with observations of the abundance and radial profile of Milky Way satellite galaxies (Bose et al. 2018, 2020) as well as with high-resolution hydrodynamical simulations that self-consistently track the evolution of low-mass dark matter haloes and their baryonic components (e.g. Grand et al. 2021). In comparison, the abundance matching approach suppresses the formation of faint galaxies more aggressively than we find when applying /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc. When combined with the possible under-correction for prematurely destroyed subhaloes, stronger constraints are an inevitable outcome. We believe these two methodological differences are the root causes of the discrepancy between the results, although this can only be confirmed with a more detailed investigation and comparison. An analysis such as this is beyond the scope of this work and we leave it to other studies.  \nOur fiducial results alone do not exclude the 𝜈 MSMas an explanation for the 3 . 55 keV excess; however, they do significantly restrict the parameter space at 𝑀 s = 7 . 1 keV. While X-ray non-detections already suggest that the 𝜈 MSM cannot explain the 3 . 55 keV excess, the lower limits are affected by uncertainties in modelling astrophysical backgrounds that could weaken these conclusions. Further light will be shed on this question by XRISM in the coming months. The Resolve instrument onboard XRISM has sufficient resolution to determine what flux originates from any emission line at 3 . 55 keV compared to any background, and also enough resolution in principle to measure the width of such lines. If these observations rule out all explanations within the realm of known physics, two possibilities remain: the 𝜈 MSMsterile neutrino could account for only some fraction of the cosmic dark matter abundance, or new physics are required to explain its origin.", '5 CONCLUSIONS': 'In this work we impose robust constraints on the parameter space of the 𝜈 MSMby comparing predictions for the number of low-mass galaxies hosted by Milky Way mass 𝜈 MSM haloes against recent estimates of the size of the Milky Way satellite galaxy population. Our results complement the updated constraints determined from X-ray non-detections, and using other astrophysical probes such as gravitational lensing, the Lyman 𝛼 forest, and the modelling of gaps in stellar streams.  \nTo determine our constraints we use the Extended Press-Schechter formalism to produce reasonable estimates of the total abundance of substructure in Milky Way mass haloes (see Fig. 1). We calibrate it to match gravity only simulations corrected to account for subhaloes artificially destroyed by numerical effects. This is important because omitting the prematurely destroyed subhaloes strengthens the constraints artificially (Newton et al. 2021). We adopt two distinct approaches to populate the subhaloes with galaxies. The first approach, in which we make no assumptions about galaxy formation physics and presume that all subhaloes host a galaxy (see Section 3.1), produces the least stringent but most robust constraints. This rules out, with at least 95 per cent confidence, all parameterizations of the 𝜈 MSMwith 𝑀 s ≤ 1 . 2 keV independently of the effective active-sterile neutrino mixing angle and the Milky Way halo mass (see Fig. 2). When adopting the Callingham et al. (2019) estimate of the Milky Way halo mass we exclude all parameterizations with 𝑀 s ≤ 1 . 4 keV. These results are consistent with constraints derived from theoretical and observational analyses of the phase  \nspace density of dark matter within Milky Way dwarf satellite galaxies (Boyarsky et al. 2009). In the latter case, the parameter space with 𝑀 s = 7 . 1 keV is partially constrained at small mixing angles. If, as has been suggested recently, the Milky Way halo mass is at the lower end of estimates, the constraints at this rest mass strengthen somewhat. In combination with the lower limits determined from X-ray non-detections, the parameter space at 𝑀 s = 7 . 1 keV closes completely in the albeit unlikely scenario in which the subhalo occupation fraction is extremely high. Adopting a more realistic estimate of the abundance of satellite galaxies in Milky Way mass haloes strengthens the constraints further.  \nThe second approach, in which we apply the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc semianalytic model of galaxy formation, incorporates physically motivated prescriptions of baryonic processes that suppress the formation of galaxies. In our analysis this renders more parameterizations of the 𝜈 MSM incompatible with estimates of the Milky Way satellite count and tightens the constraint on the viable parameter space significantly. One of the most important physical processes affecting the formation and evolution of low-mass galaxies is the reionization of hydrogen at early times. Using /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc we explore how changing the redshift at which reionization finishes, 𝑧 reion , and the circular velocity threshold, 𝑉 cut , below which galaxy formation in dark matter haloes is suppressed, affects the results (see Section 3.2). Adopting a fiducial parameterization of reionization with 𝑧 reion = 7 and 𝑉 cut = 30 kms -1 , we rule out with at least 95 per cent confidence all parameter combinations of the 𝜈 MSM with 𝑀 s ≤ 4 keV independently of effective mixing angle and Milky Way halo mass (see Fig. 3). The rest of the possible parameter space closes almost completely, and only the parameterizations of the 𝑀 s = 7 . 1 keV 𝜈 MSM with 0 . 4 ≤ sin 2 ( 2 𝜃 e ) / 10 -10 < 10 remain viable. All parameterizations of the 𝜈 MSM with 𝑀 s = 7 . 1 keV are excluded if 𝑀 MW 200 ≤ 8 × 10 11 M ⊙ . Consistent with previous work (Newton et al. 2021), we find that increasing 𝑧 reion and 𝑉 cut strengthens the constraints, and the converse also holds (see Fig. 4).  \nAkey assumption in analyses such as ours is the size of the Milky Way satellite galaxy population, which is unknown because most of the virial volume has not been surveyed yet. Our fiducial analysis adopts the Newton et al. (2018) estimate of the total population; however, several larger estimates have been published since then. Assuming larger values of the Milky Way satellite count rules out larger fractions of the 𝜈 MSM parameter space and strengthens our constraints somewhat (see Fig. 5). Notably, if 𝑁 MW sat = 183 we rule out all parameterizations of the 𝜈 MSM with 𝑀 s ≤ 7 . 1 keV, thus excluding the 𝜈 MSMas the origin of the 3 . 55 keV excess. However, compared to recent literature our constraints are less stringent.  \nWe are unable to reproduce the strongest constraints reported by Nadler et al. (2021) even when we adopt a significantly larger number of Milky Way satellites. We attribute the disparity between our results and theirs to both the different choice of particle physics calculation and to the differing methodologies used to populate the subhaloes with galaxies. The physically motivated semi-analytic prescription in /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc, which we use, contrasts with the subhalo abundance matching adopted in Nadler et al. (2021). The latter approach, which has been shown to be inconsistent with the results of high resolution simulations by Sawala et al. (2015), suppresses the formation of faint galaxies more aggressively than the /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc model. A detailed and systematic comparison of the methodologies used in both studies, while undoubtedly of significant interest, is beyond the scope of this work and we leave this for future studies.  \nThe recent launch of the JAXA XRISM satellite (Terada et al. 2021), an almost like-for-like replacement of the ill-fated Hitomi mission, will considerably improve observational capabilities at soft  \nX-ray energies once commissioning is complete. With higher sensitivity and resolution compared with previous instruments, XRISM is poised to confirm or refute the suggestion that the 3 . 55 keV excess corresponds to a detectable X-ray line. This unconfirmed signal, which currently lacks a clear astrophysical origin, has motivated considerable theoretical and observational efforts to ascertain whether it may originate from the dark sector. A successful confirmation of the signal and subsequent measurements of its width in various astrophysical systems could provide crucial evidence of the nature of the dark matter (Lovell 2023). However, our results imply that it is unlikely to correspond to the sterile neutrino proposed in the 𝜈 MSM framework.', 'ACKNOWLEDGEMENTS': 'ON thanks Mariana Jaber for useful discussions. ON acknowledges support from the Polish National Science Centre under grant UMO-2020/39/B/ST9/03494. ARJ and CSF were supported by STFC grant ST/L00075X/1. CSF acknowledges support from the European Research Council (ERC) through Advanced Investigator grant DMIDAS (GA 786910). MRL, SC, and JCH are supported by STFC grants ST/T000244/1 and ST/X001075/1. This work used the DiRAC@Durham facility managed by the Institute for Computational Cosmology on behalf of the STFC DiRAC HPC Facility ( www.dirac.ac.uk ). The equipment was funded by BEIS capital funding via STFC capital grants ST/K00042X/1, ST/P002293/1, ST/R002371/1 and ST/S002502/1, Durham University, and STFC operations grant ST/R000832/1. DiRAC is part of the National eInfrastructure.  \nSoftware : This research made use of A/s.pc/t.pc/r.pc/o.pc/p.pc/y.pc (The Astropy Collaboration et al. 2013, 2018), /m.pc/a.pc/t.pc/p.pc/l.pc/o.pc/t.pc/l.pc/i.pc/b.pc (Hunter 2007), /n.pc/u.pc/m.pc/p.pc/y.pc (van der Walt et al. 2011; Harris et al. 2020), /p.pc/y.pc/t.pc/h.pc/o.pc/n.pc (Van Rossum & Drake 2009), and /s.pc/c.pc/i.pc/p.pc/y.pc (Jones et al. 2011; Virtanen et al. 2020). This research also made use of the NASA Astrophysics Data System ( http://adsabs.harvard.edu/ ) and the arXiv e-print service ( http://arxiv.org/ ). We thank their developers for maintaining them and making them freely available.', 'DATA AVAILABILITY': 'The data used in this work are available upon reasonable request to the corresponding author. A repository of reduced data and scripts to produce the figures in this manuscript will be made available on GitHub 2 and will be archived in Zenodo.', 'References': '```\nATLAS Collaboration et al., 2021a, Phys. Rev. D, 103, 112006 ATLAS Collaboration et al., 2021b, JHEP, 2021, 226 Aad G., et al., 2023, Phys. Lett. B, 842, 137963 Aalbers J., et al., 2023, PRL, 131, 041002 Abazajian K. N., 2020, arXiv e-prints Abazajian K. N., et al., 2012, preprint ( arXiv:1204.5379 ) Abe K., et al., 2023, Phys. Rev. D, 108, 083022 Adhikari R., et al., 2017, JCAP, 2017, 025 Agnes P., et al., 2023, Phys. Rev. D, 107, 063001 Aharonian F. A., et al., 2017, ApJL, 837, L15 Albert A., et al., 2017, Phys. Lett. B, 769, 249\n```', 'APPENDIX A: DEPENDENCE OF STRUCTURE FORMATION CONSTRAINTS ON THE TOTAL NUMBER OF MILKY WAY SATELLITES': 'The constraints obtained via satellite number counts are sensitive to the size of the Milky Way satellite galaxy population adopted for the analysis. This has not been measured yet because current surveys do not have sufficient depth to survey the entire Milky Way virial volume, nor do they cover the entire sky. In Fig. A1, we present the constraints we obtain by comparing Extended Press-Schechter  \npredictions of the total abundance of substructure in Milky Way mass haloes with four different estimates of the size of the Milky Way satellite galaxy population, 𝑁 MW sat (Newtonetal. 2018; Drlica-Wagner et al. 2020; Nadler et al. 2020). As in Section 3.1, we marginalize over the uncertainties in the Callingham et al. (2019) estimate of the Milky Way halo mass. Our constraints become more restrictive as the assumed size of the Milky Way satellite galaxy population increases because this requires the parameterizations of the 𝜈 MSM model to produce more substructure to remain viable. For completeness we also show our results when adopting 𝑁 MW sat = 270 (dot-dashed line), which was reported by Drlica-Wagner et al. (2020) and was used in the analysis conducted by Dekker et al. (2022). Adopting this value of 𝑁 MW sat produces very stringent constraints on the parameter space; however, Nadler et al. (2020) noted that this estimate is inflated because the satellite distribution was assumed to be isotropic.  \nAlongside our constraints we plot the constraint envelope reported by Dekker et al. (2022), in which they use an Extended PressSchechter framework and adopt the conservative assumption that all subhaloes host galaxies, as we did in Section 3.1. As noted above, their adopted value of 𝑁 MW sat = 270 is now recognised to be too large; however, the overall approach is similar to ours so we would expect our constraints obtained when assuming the same value to be comparable with theirs. As Fig. A1 shows, our constraints are less restrictive, particularly where the effective mixing angle is large. Four aspects of their analysis contribute to their stronger constraints. First, their Extended Press-Schechter framework typically underestimates the number of subhaloes that survive until 𝑧 = 0 by up to 40 per cent compared with the results of N -body simulations from Lovell et al. (2014). This discrepancy rises to 60 per cent in the cold dark matter model (Dekker et al. 2022, appendix C). In part, this may be because they calibrate their Extended Press-Schechter framework by adopting values of the density field filter function parameters from Schneider (2015). The simulations used in that work do not account for subhaloes that were prematurely destroyed by numerical effects. We showed in Newton et al. (2021) that underestimating the abundance of substructure by omitting the population of artificially destroyed subhaloes suppresses 𝑓 v and can artificially strengthen constraints on the parameter space of WDM models.  \nSecondly, their calculation of 𝑓 v does not account fully for sources of uncertainty. The scatter in the subhalo mass function is modelled using a Poisson distribution, which is narrower than the negative binomial distribution favoured by N -body cosmological simulations (Boylan-Kolchin et al. 2010). Additionally, they do not account for the uncertainty in the estimate of the total Milky Way satellite galaxy population. Underestimating or omitting one or both of these sources of uncertainty can produce constraints that are too strict (Newton et al. 2021, section 2 . 4).  \nThird, they assume a Milky Way halo mass, 𝑀 MW 200 = 10 12 M ⊙ , which is 20 per cent less massive than the Callingham et al. (2019) estimate that we use. Adopting smaller values of 𝑀 MW 200 can increase the stringency of the constraints on the parameter space; however, as we showed in Section 3.1, when assuming that all subhaloes are populated by galaxies the variation of the constraints as a function of 𝑀 MW 200 is not significant. Therefore, we think that this is a subdominant contribution to the discrepancy between our results and those of Dekker et al. in Fig. A1.  \nFinally, in common with Nadler et al. (2021), Dekker et al. adopt the 𝜈 MSMmomentumdistributions computed by Venumadhav et al. (2016). They are warmer than the momentum distributions we use, which complicates the comparison of our results with theirs.  \nFigure A1. Structure formation constraints on the 𝜈 MSM parameter space from analyses of Milky Way satellite counts. As for several figures in the main text, the shaded regions in the upper and lower portions of the parameter space delimit parameter combinations for which the dark matter is over- or underproduced if resonantly produced sterile neutrinos compose all of the dark matter (Asaka et al. 2005; Schneider 2016). Our results are marginalized over the uncertainties in estimates of the Milky Way halo mass from Callingham et al. (2019). We compute constraint envelopes assuming a total Milky Way satellite galaxy population with a mean value of 𝑁 MW sat = 124(thick solid line), 134 (dashed line), 183 (dotted line), and 270 (dot-dashed line). Parameter combinations within the envelopes are ruled out with at least 95 per cent confidence. For comparison we show the constraint envelope from Dekker et al. (2022).  \n<!-- image -->', 'APPENDIX B: MODIFICATIONS TO THE MERGER TREE CONSTRUCTION ALGORITHM': "In Newton et al. (2021, appendix B), we noted that the Monte Carlo merger trees generated in /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc using the standard implementation of the Parkinson et al. (2008) algorithm overestimate the abundance of low-mass galaxies by a factor of two. We addressed this by applying an empirical correction to the theoretical galaxy luminosity function in which the satellite galaxy absolute magnitude, 𝑀 V , predicted by /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc applied to the Monte Carlo merger trees is mapped to that predicted by /g.pc/a.pc/l.pc/f.pc/o.pc/r.pc/m.pc applied to the COCO simulation merger trees. This produced results in good agreement with the simulations; however, they are not self-consistent predictions of the Extended Press-Schechter methodology.  \nThe discrepancy between the Monte Carlo and N -body predictions arises partly because the progenitor mass functions generated by the standard Parkinson et al. (2008) algorithm under-predict those producedin cosmological simulations at the low-mass end. This is because the Parkinson et al. (2008) algorithm was calibrated against the Millennium simulation, which cannot probe such low halo masses. The Parkinson et al. (2008) algorithm generates merger trees probabilistically via a sequence of branching events in which the mass, 𝑀 1 , in the halo at time, 𝑡 1 , may be partitioned into two lower-mass progenitors with masses, 𝑀 2 and ( 𝑀 1 -𝑀 2 ) , at an earlier time, 𝑡 2 . A complete merger tree is built up by repeating this process for each progenitor halo at each successive timestep. The size of the timestep, d 𝑡 , is selected to be sufficiently small such that the halo is  \nunlikely to have more than two progenitors above some mass resolution threshold at the earlier time. The rate of progenitor formation is controlled by the binary merger branching rate,  \nd 𝑝 d 𝜔 = ∫ 𝑀 / 2 𝑀 min 𝑀 1 𝑀 ' d 𝑓 ( 𝑀 ' ) d 𝑡 d 𝑆 ( 𝑀 ' ) d 𝑀 ' GLYPH<12> GLYPH<12> GLYPH<12> GLYPH<12> d 𝑡 d 𝜔 GLYPH<12> GLYPH<12> GLYPH<12> GLYPH<12> × 𝐺 GLYPH<2> 𝜔, 𝜎 ( 𝑀 1 ) , 𝜎 GLYPH<0> 𝑀 ' GLYPH<1> GLYPH<3> d 𝑀 ' , (B1)  \nwhere 𝜔 = 𝛿 c , 0 / 𝐷 ( 𝑡 ) is the ratio of the linear theory critical overdensity threshold for gravitational collapse in cold dark matter models, 𝛿 c , 0 , to the linear growth-rate factor, 𝐷 ( 𝑡 ) ; 𝑀 min is the required minimum mass resolution of the merger tree; d 𝑓 / d 𝑡 is the first-crossing rate distribution, the rate at which the excursion set trajectories first cross the overdensity threshold (Bond et al. 1991); 𝜎 ( 𝑀 ) = √︁ 𝑆 ( 𝑀 ) is the square root of the variance, S , of the matter density field; and 𝐺 [ 𝜔, 𝜎 ( 𝑀 ) , 𝜎 ( 𝑀 ' )] is an empirical modification introduced to obtain results consistent with those from cosmological cold dark matter simulations. In the standard implementation of the Parkinson et al. (2008) algorithm G is given by  \n𝐺 [ 𝜔, 𝜎 ( 𝑀 1 ) , 𝜎 ( 𝑀 2 )] = 𝐺 0 GLYPH<20> 𝜎 ( 𝑀 2 ) 𝜎 ( 𝑀 1 ) GLYPH<21> 𝛾 1 GLYPH<20> 𝜔 𝜎 ( 𝑀 1 ) GLYPH<21> 𝛾 2 , (B2)  \nwhich introduces the free parameters, 𝐺 0 , 𝛾 1 , and 𝛾 2 that are calibrated using dark matter-only cosmological simulations (see e.g. Benson et al. 2019, table 1).  \nAs we discussed above, this implementation of the algorithm over-predicts the number of low-mass subhaloes because the empirical modification underestimates the simulated progenitor mass functions. We attempt to address this by introducing an explicit dependence on 𝐷 ( 𝑡 ) and the logarithmic slope of the 𝜎 ( 𝑀 ) relation such that  \n𝐺 [ 𝜔, 𝜎 ( 𝑀 1 ) , 𝜎 ( 𝑀 2 )] = 𝐺 0 GLYPH<20> 𝜎 ( 𝑀 2 ) 𝜎 ( 𝑀 1 ) GLYPH<21> 𝛾 1 GLYPH<20> 𝜔 𝜎 ( 𝑀 1 ) GLYPH<21> 𝛾 2 × ( 1 -GLYPH<20> 𝜎 ( 𝑀 2 ) 𝜎 ( 𝑀 1 ) GLYPH<21> 2 ) 𝛾 3 𝐷 ( 𝑡 2 ) 𝛾 4 GLYPH<12> GLYPH<12> GLYPH<12> GLYPH<12> d ln 𝜎 ( 𝑀 2 ) d ln 𝑀 GLYPH<12> GLYPH<12> GLYPH<12> GLYPH<12> 𝛾 5 . (B3)  \nThis paramterization was chosen to ensure that, if applied to a scalefree cosmology, the merger trees would satisfy the required selfsimilarity constraint (Efstathiou et al. 1988). This introduces a further three new free parameters: 𝛾 3 , 𝛾 4 , and 𝛾 5 that are also determined via calibration to N -body simulations. Throughout this work, we use the latter form of G with 𝐺 0 = 0 . 5024 , 𝛾 1 = 0 . 2493 , 𝛾 2 = -0 . 095846 , 𝛾 3 = 0 . 045166 , 𝛾 4 = 0 . 2342 , and 𝛾 5 = -0 . 1878, based on a preliminary recalibration to recent simulations (A. J. Benson et al., in preparation).  \nThis paper has been typeset from a T E X/L A T E X file prepared by the author."}