The stellar masses of galaxies are measured using integrated light via several methods  however few of these methods were designed for lowmass Mstarlesssim108rmModot dwarf galaxies whose properties e.g. stochastic star formation low metallicity pose unique challenges for estimating stellar masses. In this work we quantify the precision and accuracy at which stellar masses of lowmass galaxies can be recovered using UVopticalIR photometry. We use mock observations of 469 lowmass galaxies from a variety of models including both semiempirical models GRUMPY UniverseMachineSAGA and cosmological baryonic zoomin simulations MARVELous Dwarfs and FIRE2 to test literature colorMstarL relations and multiwavelength spectral energy distribution SED mass estimators. We identify a list of best practices for measuring stellar masses of lowmass galaxies from integrated photometry. These include updated prescriptions for stellar mass based on gr color and WISE 3.4 mum luminosity which are less systematically biased than literature calibrations and can recover true stellar masses of lowmass galaxies with sim0.1 dex precision. When using SED fitting to estimate stellar mass we find that the form of the assumed star formation history can induce significant biases parametric SFHs can underestimate stellar mass by as much as sim0.4 dex while nonparametric SFHs recover true stellar masses with insignificant offset 0.03pm0.11 dex. However we also caution that noninformative dust attenuation priors may introduce Mstar uncertainties of up to sim0.6 dex.
['Astrophysics - Astrophysics of Galaxies']
Stellar Mass Calibrations for Local LowMass Galaxies
{'Stellar Mass Calibrations for Local Low-Mass Galaxies': 'Mithi A. C. de los Reyes, 1 Yasmeen Asali, 2 Risa H. Wechsler, 3, 4, 5 Marla Geha, 2 Yao-Yuan Mao, 6 Erin Kado-Fong, 2 Ragadeepika Pucha, 6 William Grant, 1 Pratik J. Gandhi, 2, 7 Viraj Manwadkar, 3, 4, 5 Anna Engelhardt, 8 Ferah Munshi, 8 and Yunchong Wang 3, 4, 5  \n1  \n7 8 Department of Physics, George Mason University, 4400 University Dr, Fairfax, VA 22030, USA  \nDepartment of Physics & Astronomy, Amherst College, 6 East Drive, Amherst, MA 01002, USA 2 Department of Astronomy, Yale University, New Haven, CT 06520, USA 3 Kavli Institute for Particle Astrophysics & Cosmology, Stanford University, 452 Lomita Mall, Stanford, CA 94305, USA 4 Department of Physics, Stanford University, 382 Via Pueblo Mall, Stanford, CA 94305, USA 5 SLAC National Accelerator Laboratory, 2575 Sand Hill Road, Menlo Park, CA 94025, USA 6 Department of Physics and Astronomy, University of Utah, Salt Lake City, UT 84112, USA Department of Physics and Astronomy, University of California - Davis, One Shields Avenue, Davis, CA 95616, USA', 'ABSTRACT': "The stellar masses of galaxies are measured using integrated light via several methods - however, few of these methods were designed for low-mass ( M ⋆ ≲ 10 8 M ⊙ ) 'dwarf' galaxies, whose properties (e.g., stochastic star formation, low metallicity) pose unique challenges for estimating stellar masses. In this work, we quantify the precision and accuracy at which stellar masses of low-mass galaxies can be recovered using UV/optical/IR photometry. We use mock observations of 469 low-mass galaxies from a variety of models, including both semi-empirical models (GRUMPY, UniverseMachine-SAGA) and cosmological baryonic zoom-in simulations (MARVELous Dwarfs and FIRE-2), to test literature colorM ⋆ /L relations and multi-wavelength spectral energy distribution (SED) mass estimators. We identify a list of 'best practices' for measuring stellar masses of low-mass galaxies from integrated photometry. These include updated prescriptions for stellar mass based on g -r color and WISE 3.4 µ m luminosity, which are less systematically biased than literature calibrations and can recover true stellar masses of low-mass galaxies with ∼ 0 . 1 dex precision. When using SED fitting to estimate stellar mass, we find that the form of the assumed star formation history can induce significant biases: parametric SFHs can underestimate stellar mass by as much as ∼ 0 . 4 dex, while non-parametric SFHs recover true stellar masses with insignificant offset ( -0 . 03 ± 0 . 11 dex). However, we also caution that non-informative dust attenuation priors may introduce M ⋆ uncertainties of up to ∼ 0 . 6 dex.", '1. INTRODUCTION': "Much of observational astronomy can be distilled into a single exercise: converting the measured light from astrophysical objects to physical properties. The mass of a galaxy is one of its most fundamental properties, and one that heavily influences not only a galaxy's observable properties, but its evolution over cosmic time. The stellar mass ( M ⋆ ) of a galaxy is a particularly useful property that encodes the galaxy's cumulative star formation, merger, and accretion history. As a result, stellar mass is frequently used as a variable in empirical  \nrelationships, such as the stellar mass-halo mass relation (e.g., Moster et al. 2010; Wechsler & Tinker 2018; Girelli et al. 2020), the star formation main sequence (e.g., Brinchmann et al. 2004; Popesso et al. 2023), and the mass-metallicity relation (e.g., Tremonti et al. 2004; Maiolino & Mannucci 2019), which are in turn used to probe galaxy evolution processes.  \nMeasuring the stellar mass of a galaxy from its integrated light - the only option for most galaxies, which are too distant to resolve individual stars - is a nontrivial exercise. To the first order, a galaxy's stellar mass is correlated with its total luminosity: more stars produce more light. However, determining the exact conversion from light to stellar mass (i.e., the mass-to-light ratio, M ⋆ /L ) is complicated by several factors.  \nFor a single coeval stellar population, the most massive stars contribute the majority of light, while the more numerous lower-mass stars contribute the majority of mass; the ratio of high- to low-mass stars is set at birth by the stellar initial mass function (IMF) and modified by stellar evolution. A galaxy consists of multiple such stellar populations with varying ages and metallicities. The light from these populations is also modified by baryonic components like gas and dust located both within the galaxy and along the observed line of sight. To estimate M ⋆ /L , one can model each of these components - the stellar populations and their evolution, nebular emission, and dust attenuation and emission - to obtain a model for a galaxy's observed spectrum. This is the basis of stellar population synthesis (SPS) modeling (Tinsley 1968), the most common method of estimating a galaxy's stellar mass. SPS models are fit to the galaxy's observed light, which are typically in the form of (1) single fluxes or colors, (2) broadband spectral energy distributions (SEDs), or (3) spectra (for a review, see Conroy 2013). Measurements of M ⋆ therefore depend on all of the assumptions implicit in SPS modeling: the observable properties of a single stellar population with a given age and metallicity (which depends on the form of the IMF and stellar evolution models), the stellar populations that comprise a galaxy (which depends on the galaxy's assumed star formation history or SFH), and the impact of gas and dust (which depends on heavily simplified models).  \nSeveral works have tested the effect of these assumptions on M ⋆ measurements; in particular, quantifying the uncertainties from SED fitting has become especially relevant in the era of JWST, which has begun to obtain SEDs of galaxies in the early universe. For example, Lower et al. (2020) tested the impact of different SFH assumptions when measuring M ⋆ from SED fitting, and Pacifici et al. (2023) tested the effect of using different SED fitting codes. However, nearly all of these tests have focused on relatively massive galaxies: although Lower et al. (2020)'s sample included galaxies down to 10 7 . 6 M ⊙ , the bulk of their sample was in the higher mass range of 10 8 -12 M ⊙ , and Pacifici et al. (2023) did not consider galaxies below 10 8 M ⊙ at all.  \nLow-mass 'dwarf' galaxies ( M ⋆ ≲ 10 8 M ⊙ ) have often been overlooked in these tests, partly because they pose unique problems for SPS modeling. They are thought to have relatively complex, stochastic star formation histories (Emami et al. 2019), making their SEDs and spectra more difficult to fit with SPS models. Furthermore, previous large galaxy surveys have often been incomplete in the low-mass regime, making it natural to prioritize high-mass galaxies.  \nHowever, low-mass galaxies are the most numerous type of galaxy in the universe (Schechter 1976). Due to their small gravitational potentials, they are sensitive to both internal processes like stellar feedback (Collins & Read 2022) and external environmental effects (e.g., Dekel et al. 2003; Farouki & Shapiro 1981; Deason et al. 2014), making them excellent testbeds for cosmological and galaxy evolution models. With the advent of a new generation of extragalactic surveys - e.g., the Dark Energy Spectroscopic Instrument (DESI; DESI Collaboration et al. 2016) and DESI Legacy Surveys (Dey et al. 2019), the Legacy Survey of Space and Time (LSST; Ivezi'c et al. 2019), Euclid (Euclid Collaboration et al. 2024), and Roman (Akeson et al. 2019) - we are beginning to observe unprecedentedly large samples of lowmass galaxies across a range of environments. Accurately characterizing the masses of these galaxies is of critical importance: the low-mass end of the galaxy mass function can constrain cosmological models (Sales et al. 2022), and extending scaling relations to low masses can constrain galaxy physics (e.g., Woo et al. 2008).  \nIn this work, we test a number of stellar mass measurement techniques on simulated low-mass galaxies ( M ⋆ = 10 4 -8 M ⊙ ). In light of several ongoing and upcoming galaxy surveys, we primarily focus on methods that use UV/optical/IR photometry. We identify some best practices for measuring low-mass galaxy masses, and we develop statistical corrections to extend commonly used empirical calibrations to the low-mass regime.  \nThe structure of this paper is as follows. We present the stellar mass calibrations used in this paper in Section 2, and we describe the observed and simulated datasets used to test these stellar mass calibrations in Section 3. We then test the recovery of 'ground-truth' stellar masses from simulated data: using photometric calibrations in Section 4, and using SED fitting in Section 5. We summarize our main results and list 'best practices' for measuring low-mass galaxy masses in Section 6. In Section 7, we discuss potential systematics in our experiment, and we summarize our conclusions in Section 8.", '2. STELLAR MASS CALIBRATIONS': 'Here we describe the stellar mass measurement techniques considered in this work. We emphasize that this is not an exhaustive list, particularly since we only focus on one SED fitting code (Section 2.2). We convert all calibrations to a Chabrier (2003) IMF.', '2.1. Empirical photometric calibrations': "We first consider empirical calibrations based on optical and IR photometry. Such approaches are attractive  \nfrom an observational perspective because measuring a galaxy's flux in one or two photometric bands requires far fewer resources than measuring a multi-band SED over a broad wavelength range, let alone a spectrum. Fortunately, several previous studies have found relatively tight correlations between stellar mass (typically measured from full SED fitting) and combinations of different photometric bands. Historically, these have often been parameterized as relationships between M ⋆ /L and a single color (e.g., Bell & de Jong 2001), but a number of studies have developed stellar mass calibrations based on anywhere from one to three photometric bands. In this work, we use the following calibrations.", '2.1.1. Optical M ⋆ /L calibrations': 'While there are many optical-NIR colorM ⋆ /L relations (e.g., Portinari et al. 2004; Zibetti et al. 2009; Gallazzi & Bell 2009; Taylor et al. 2011; Roediger & Courteau 2015; Schombert et al. 2019; Kouroumpatzakis et al. 2023), we have chosen to focus first on the colorM ⋆ /L relations by Bell et al. (2003) and Into & Portinari (2013) because these have been used in a number of recent low-mass galaxy surveys. For example, the Satellites Around Galactic Analogs survey (SAGA; Geha et al. 2017) and other surveys of low-metallicity/lowmass galaxies (Hsyu et al. 2018) use recalibrations of the Bell et al. (2003) relations to measure stellar masses, while the Exploration of Local VolumE Satellites survey (ELVES; Carlsten et al. 2022; Danieli et al. 2023) use a colorM ⋆ /L relation from Into & Portinari (2013).  \nIn the following equations in this section, all optical magnitudes are AB magnitudes that have been both corrected for Milky Way extinction and k -corrected to z = 0 using g -r color (Chilingarian et al. 2010). The corrected r - and g -band magnitudes are denoted as M r, 0 and M g, 0 . Unless otherwise noted, we assume absolute solar magnitudes from Willmer (2018): 4.65 for r -band and 5.11 for g -band.  \nBell et al. (2003) developed colorM ⋆ /L relations based on optical-NIR SEDs of galaxies with photometry from the Two Micron All Sky Survey (2MASS) and the Sloan Digital Sky Survey (SDSS). We use their g -r relation:  \nlog( M ⋆ / M ⊙ ) = 1 . 461 + 1 . 098(g -r) 0 -0 . 4M r , 0 . (1)  \nThis relation was calibrated using a sample of galaxies with masses down to ∼ 10 8 . 5 M ⊙ . To fit these galaxies, Bell et al. (2003) constructed a grid of SPS models with  \nmetallicities between 10 -2 . 3 < Z/ Z ⊙ < 10 0 . 41 . They assumed exponentially decaying SFHs and did not account for dust attenuation, although they ran some simple tests to quantify the impact of adding star formation bursts and a simple constant dust reddening model.  \nInto & Portinari (2013) also proposed colorM ⋆ /L relations using more detailed asymptotic giant branch models. Here we consider the relation between g -r and M ⋆ /L g used by the ELVES survey (Danieli et al. 2023):  \nlog( M ⋆ / M ⊙ ) = 1 . 332 + 1 . 774(g -r) 0 -0 . 4M g , 0 (2)  \nNote that this equation uses a slightly different absolute solar g -band magnitude of 5.144 (Table 1 of Into & Portinari 2013). Like Equation 1, this relation also assumes an exponentially decaying SFH; however, it was determined from theoretical models, rather than semiempirically from observed photometry.  \nThe majority of colorM ⋆ /L relations in the literature, including Equations 1 and 2 above, were calibrated using samples of relatively massive ( M ⋆ > 10 8 -9 M ⊙ ) galaxies. We therefore also consider empirical relations that are specifically based on low-mass galaxies.  \nHerrmann et al. (2016) measured colorM ⋆ /L relations for 34 nearby dIrr galaxies from Local Irregulars That Trace Luminosity Extremes, The H I Nearby Galaxy Survey (LITTLE THINGS; Hunter et al. 2012). We consider their g -r relation:  \nlog( M ⋆ / M ⊙ ) = 1 . 755 + 0 . 894(g -r) 0 -0 . 4M g , 0 . (3)  \nStellar mass density profiles for these galaxies were measured by modeling galaxy SEDs with a multi-component stellar population model.  \nKlein et al. (2024) proposed a colorM ⋆ /L relation based on low-mass galaxies in the FIRE-2 simulations (see Section 3.1):  \nlog( M ⋆ / M ⊙ ) = 1 . 631 + 1 . 570(g -r) 0 -0 . 4M g , 0 . (4)  \nThis relation is based on 20 simulated galaxies with stellar masses between 10 5 -10 M ⊙ . Klein et al. (2024) used the FIRE studio code (Gurvich 2022; Hopkins et al. 2005) to create post-processed images via ray projection from the star particles, then measured synthetic photometry from these mock images.  \nFinally, Du et al. (2020) did not develop a new photometric calibration, but instead found corrections for literature colorM ⋆ /L relations in order to make them  \ninternally self-consistent when using different photometric bands (McGaugh & Schombert 2014). They used a sample of low surface brightness galaxies to empirically correct several literature relations, including both Bell et al. (2003) and Into & Portinari (2013), listed respectively:  \nlog( M ⋆ / M ⊙ ) = 1 . 461 + 1 . 097( g -r ) 0 -0 . 4 M r, 0 (5)  \nlog( M ⋆ / M ⊙ ) = 1 . 254 + 1 . 530( g -r ) 0 -0 . 4 M r, 0 (6)  \n2.1.2. Infrared (WISE) M ⋆ /L calibrations  \nJarrett et al. (2023) computed updated M ⋆ /L calibrations based on near-infrared photometry from the Widefield Infrared Survey Explorer (WISE). They report two stellar mass calibrations, the first based on W1 (3.4 µ m) flux and the second based on W1 -W2 (3.4 µ m-4.6 µ m) color:  \nlog( M ⋆ / M ⊙ ) = -12 . 62 + 5 . 00 log( L W1 , 0 ) -0 . 44 log( L W1 , 0 ) 2 +0 . 016 log( L W1 , 0 ) 3 (7)  \nlog( M ⋆ / M ⊙ ) = log( L W1 , 0 ) -0 . 376 -1 . 053(W1 -W2) 0 (8)  \nHere, L W1 , 0 = 10 -0 . 4( M W1 , 0 -M W1 , ⊙ ) , where the W1 inband solar value is M W1 , ⊙ = 3 . 24 mag (Jarrett et al. 2013). In Equations 7 and 8, all WISE magnitudes and colors are in Vega magnitudes (to convert from AB to Vega, we follow the prescriptions in Table 1 of Jarrett et al. 2011) and are k -corrected following Equations A1 and A2 in Jarrett et al. (2023). These calibrations were based on galaxies in the Galaxy and Mass Assembly (GAMA) survey with WISE fluxes. This sample contained galaxies with stellar masses down to ∼ 10 6 . 5 M ⊙ , but is likely incomplete below ∼ 10 8 M ⊙ . The GAMA DR4 stellar masses (Driver et al. 2022) were computed by fitting a grid of Bruzual & Charlot (2003) SPS models (for more details, see Taylor et al. 2011) and assumed exponential SFHs, uniform metallicities, and Calzetti et al. (2000) dust attenuation.', '2.2. Prospector': "The empirical relations considered in Section 2.1 use only a few photometric bands to estimate stellar masses. Fitting the multi-wavelength SED of a galaxy takes advantage of photometric information across a much broader wavelength range and is frequently considered the 'gold standard' for obtaining galaxy stellar masses from integrated light. Prospector (Leja et al. 2017; Leja et al. 2019; Johnson et al. 2021) is a Bayesian inference code that estimates galaxy properties by forward modeling galaxy SEDs with the Flexible Stellar Population  \nSynthesis package (FSPS; Conroy et al. 2009; Conroy & Gunn 2010).  \nAs we describe in Section 3, we also use FSPS to produce mock observations of simulated galaxies. Using Prospector to fit these mock observations is therefore an 'apples-to-apples' comparison that directly probes the capabilities of SED fitting rather than the effect of other systematic uncertainties (e.g., stellar model libraries, SED fitting algorithm). We defer tests of these other systematic effects, caused by variations among different SED fitting codes, to future work.", '3. SIMULATED GALAXIES': "We test stellar mass calibrations on mock observations of simulated galaxies, for which the 'true' stellar masses are known. Here we describe the models used and their major assumptions, as well as the method used to produce mock observations from the simulated galaxies.", '3.1. Galaxy models': "3.1.1. GRUMPY  \nThe simplest model we consider is GRUMPY (Galaxy formation with RegUlator Model in PYthon; Kravtsov & Manwadkar 2022), a semi-analytic model of the 'regulator' type. In this model, a low-mass galaxy is simulated using a system of differential equations that track the mass conservation of different baryonic components. The key input to the model is an underlying dark matter halo accretion history, which is taken from the halo mass accretion histories from the ELVIS high-resolution simulation suite (Garrison-Kimmel et al. 2014).  \nThe gas inflow rate is assumed to be proportional to the halo accretion rate with additional factors describing suppression due to UV heating (i.e., from reionization) and the formation of a hot gaseous halo. The gas outflow rate is assumed to be proportional to the star formation rate with the constant of proportionality (the mass-loading factor) parametrized as a function of stellar mass. The SFH is set by the star formation rate, which assumes a constant molecular hydrogen gas depletion time and instantaneous recycling. Only in-situ formed stars are considered in this model; stellar contribution from galaxy mergers is not accounted for in the model. The production and removal of heavy element abundances in the galaxy's ISM are also parameterized, providing a chemical evolution history for the galaxy. Finally, further stochasticity is added to the SFH by introducing a correlated random perturbation in SFR relative to the mean M ⋆ -SFR relation at each timestep (Pan & Kravtsov 2023).  \nIn total, we obtain 69 galaxy SFHs and chemical enrichment histories from GRUMPY (with boosted  \nstochasticity). This sample covers a stellar mass range of 10 4 -11 M ⊙ .", '3.1.2. UniverseMachine-SAGA': "The next model we consider is the UniverseMachine (Behroozi et al. 2019), an empirical galaxy-halo connection model. This model assumes that galaxy SFRs depend on host halo properties (halo mass, assembly history, and redshift). It then fits this relationship by predicting galaxy observables (e.g., stellar mass functions, quenched fractions) from a dark-matter-only simulation, and iteratively comparing the predictions with observations over a wide range of galaxy masses and redshifts simultaneously. Wang et al. (2021) extended the DR1 of the UniverseMachine model 2 to lower-mass galaxies by applying it to a joint set of 45 dark-matter-only zoom-in simulations of isolated Milky Way-mass halos from the Symphony compilation (Nadler et al. 2023).  \nIn this work, we use SFHs predicted by the updated version of UniverseMachine, UM-SAGA (Wang et al. 2024), which is constrained by SAGA satellites (Geha et al. 2017; Mao et al. 2021, 2024; Geha et al. 2024) and SDSS isolated galaxies (Geha et al. 2012) down to M ⋆ ≳ 10 7 M ⊙ . UM-SAGA provides a better match to observed low-mass galaxies properties than UniverseMachine DR1, and is particularly better at predicting the quenched fraction of low-mass galaxies. We apply UMSAGA to one of the 45 Symphony Milky Way-mass halos and obtain 43 satellite galaxy SFHs. Here, we define a 'satellite' galaxy in UM-SAGA as any subhalo within 5 virial radii of the host, with a virial mass of at least 1 . 2 × 10 8 M ⊙ (corresponding to 300 particle masses).  \nOf these 43 simulated galaxies, 35 have stellar masses < 10 7 M ⊙ where the quenched fraction is not well constrained by observations and may lead to additional uncertainties in the SFHs; we include these galaxies in all plots for illustrative purposes but do not include them in our analysis. This model does not track chemical enrichment, so we produce chemical enrichment histories by applying a mass-metallicity relation to the SFHs. For consistency, we obtain this mass-metallicity relation by fitting a simple linear relation 3 to the GRUMPY models (Figure 1):  \nlog ( Z ⋆ Z ⊙ ) = 0 . 43 log ( M ⋆ M ⊙ ) -6 . 24 (9)", '2 https://bitbucket.org/pbehroozi/universemachine/src/main/': '3 This and all other fits in this work are, unless otherwise specified, computed with least-squares minimization using a LevenbergMarquardt algorithm. We bootstrap all fits by performing each fit N = 1000 times, sampling with replacement on each iteration; the reported best-fit parameters are the median estimates.  \nFigure 1. Relationships between stellar mass and stellar metallicity for the simulations that track chemical enrichment: GRUMPY (dark purple points), MARVEL (pink points), and FIRE-2 (yellow points). These models have roughly consistent mass-metallicity relationships, suggesting that any differences between their derived stellar masses are not due to metallicity. The red dashed line indicates the best linear fit to the GRUMPY data; as described in Section 3.1, this fit is used to compute chemical enrichment histories for the UniverseMachine simulated galaxies. The observed stellar mass-metallicity relations from Kirby et al. (2013) (blue solid line indicates the best fit to Local Group dwarf galaxies; shaded regions indicate the rms about the best-fit line) and Gallazzi et al. (2005) (green solid line indicates the running median; shaded regions indicate the 16th and 84th percentiles) are also shown for comparison.  \n<!-- image -->', '3.1.3. MARVELous Dwarfs': "Finally, we use simulated galaxies from two different cosmological baryonic zoom-in simulations. The first is the MARVELous Dwarfs (Bellovary et al. 2021; Munshi et al. 2021; Christensen et al. 2024; Azartash-Namin et al. 2024), a suite of high-resolution zoom-in simulations that use the N-body code ChaNGA (Menon et al. 2015). The MARVELous Dwarfs (hereafter abbreviated as MARVEL) include four simulated volumes of low-mass galaxies in isolated environments using a WMAP3 cosmology (Spergel et al. 2007). These simulations implement gas, initial star, and dark matter particle masses of 1410 M ⊙ , 420 M ⊙ , and 6650 M ⊙ , respectively, and have a force resolution of 60 pc. The high spatial resolution and small particle mass allow MARVEL galaxies to be resolved down to M ∗ ≈ 3 × 10 3 M ⊙ . MARVEL simulates a number of baryonic processes, including: star formation and gas cooling (Christensen et al. 2012), metal line cooling and metal diffusion (Shen et al. 2010), non-equilibrium formation and destruction of molecular hydrogen, photoionization and heating from a cosmological UV field background (Haardt  \n& Madau 2012), and a 'blastwave' model of supernova feedback (Stinson et al. 2006) that temporarily shuts off gas cooling during the momentum-conserving 'snowplow phase' (McKee & Ostriker 1977). The combined processes emulate the energy deposited within the interstellar medium by all processes related to young stars, including UV radiation (see Wise et al. 2012; Agertz et al. 2013). The Amiga Halo Finder (AHF; Knollmann & Knebe 2009) is applied to MARVEL to identify dark matter halos, subhalos, and the baryonic content within. We obtain 30 SFHs and chemical enrichment histories of satellite and field galaxies from MARVEL, with a stellar mass range of 10 5 . 7 -10 9 . 2 M ⊙ .", '3.1.4. FIRE-2': "We also consider the Latte (introduced in Wetzel et al. 2016) and ELVIS on FIRE (introduced in GarrisonKimmel et al. 2019a) suites from the FIRE-2 cosmological baryonic zoom-in simulations, part of the Feedback In Realistic Environments project (public data release introduced in Wetzel et al. 2023). The FIRE-2 simulations are run using Gizmo, a Lagrangian Meshless Finite Mass (MFM) hydrodynamics code (Hopkins 2015), and the FIRE-2 physics model (Hopkins et al. 2018). FIRE-2 models the dense, multiphase interstellar medium (ISM) in galaxies and incorporates physically motivated, metallicity-dependent radiative heating and cooling processes for gas. These include free-free, photoionisation and recombination, Compton, photo-electric and dust collisional, cosmic ray, molecular, metal-line, and fine structure processes. The model tracks 11 element species (H, He, C, N, O, Ne, Mg, Si, S, Ca, Fe) across a temperature range of 10 - 10 10 K. FIRE-2 also includes the following time-resolved stellar feedback processes: core-collapse and white dwarf (Type Ia) supernovae, continuous mass loss, radiation pressure, photoionisation, and photo-electric heating.  \nIn this analysis, we use a sample of 362 low-mass galaxies around MW/M31-mass galaxies from the z = 0 snapshots from the Latte and ELVIS on FIRE suites, and obtain their SFHs and chemical enrichment histories. The simulations have a dark matter mass resolution of 2 -3 . 5 × 10 4 M ⊙ , and gas and star mass resolution of 3500 -7100 M ⊙ , which allows for well-resolved galaxies down to M ∗ ≈ 10 5 M ⊙ . Identifying dark matter halos, and subhalos is done via the ROCKSTAR halo finder (Behroozi et al. 2013). Importantly, this sample of FIRE-2 galaxies has been benchmarked against the stellar mass functions, radial distance distributions, and star-formation histories of low-mass galaxies in the LG (Wetzel et al. 2016; Garrison-Kimmel et al. 2019a,b; Samuel et al. 2020, 2021).  \nIn total, our sample consists of 469 simulated galaxies 4 with a mass range of M ⋆ = 10 4 . 0 -10 . 8 M ⊙ . The majority of our sample (403 galaxies) are low-mass 'dwarf' galaxies, with M ⋆ < 10 8 M ⊙ . The remaining 66 galaxies constitute a relatively small sample that is likely incomplete at the highest masses ( ≳ 10 9 . 5 M ⊙ ); however, this high-mass sample still provides a useful comparison when testing different mass measurement techniques.  \nWhile it is not possible to determine whether our simulated galaxies are truly representative of real low-mass galaxies, we can at least check whether the models produce scaling relations that agree with each other and with observations. As shown in Figure 1, both the FIRE-2 and MARVEL galaxies follow roughly the same mass-metallicity relation as the GRUMPY galaxies (and by extension the UM-SAGA galaxies, which by definition follow the red dotted line in Figure 1, given by Equation 9). Figure 1 also shows the observed stellar mass-metallicity relations measured by Kirby et al. (2013) and Gallazzi et al. (2005). The mass-metallicity relation of our sample is somewhat steeper than that measured for Local Group dwarf galaxies by Kirby et al. (2013), who measured a slope of 0.30 compared to our slope of 0.43, but the two relations agree (within the rms about the Kirby et al. 2013 best-fit line) in the mass range 10 6 -8 M ⊙ , which spans most of our lowmass galaxy sample. Where our sample overlaps with the higher mass range of the Gallazzi et al. (2005) massmetallicity relation, our simulated galaxies and simple linear fit are consistent with the Gallazzi et al. (2005) relation within 1 σ .", '3.2. Mock observations': "To produce mock observations of these simulated galaxies, we use the FSPS package (Conroy et al. 2009; Conroy & Gunn 2010) 5 . FSPS takes galaxy SFHs and chemical enrichment histories as input and generates spectra and photometry. For all mock galaxies, we use MIST isochrones (Choi et al. 2016) and the MILES spectral library (Falc'on-Barroso et al. 2011; Vazdekis et al. 2010, 2015) to model stellar populations and spectra. We assume a Calzetti et al. (2000) attenuation curve, normalized by setting the optical depth at 5500 ˚ A to τ V = 0 . 2. We include dust emission using the default FSPS parameters. Finally, to further mimic realistic conditions, we place all simulated galaxies at a redshift  \nr  \ng  \nFigure 2. The g -r color as a function of stellar mass for the GRUMPY (dark purple circles), UM-SAGA (purple triangles), MARVEL (pink squares), and FIRE-2 (yellow contours) models. UM-SAGA SFHs are not well constrained below log( M ⋆ / M ⊙ ) ∼ 7, so these galaxies are marked with small points; for clarity, FIRE-2 galaxies are represented as density contours rather than points.  \n<!-- image -->  \nTable 1. Photometric filters used to produce mock galaxy SEDs.  \nof z = 0 . 01, consistent with the distances probed by ongoing and upcoming surveys of the nearby universe (e.g., the DESI LOWZ survey aims to target low-mass galaxies out to z < 0 . 03; Darragh-Ford et al. 2023). 6  \nUsing the input SFH and chemical evolution history of a galaxy, FSPS calculates the surviving stellar mass of the galaxy. We use this as 'true' stellar mass of the galaxy, rather than the raw output mass measured from GRUMPY, UM-SAGA, MARVEL, or FIRE. This is because the four models assume a variety of different stel-  \n6 We note that this is not strictly necessary, since many of the tests described in this work assume that the redshift is exactly known. However, we plan to test the effect of unknown redshift on stellar mass measurements in future works.  \nFigure 3. Cumulative SFHs of the simulated galaxies considered in this work. The panels illustrate galaxies of different stellar mass ranges: 10 4 -6 M ⊙ (top), 10 6 -8 M ⊙ (middle), and 10 8 -10 M ⊙ (bottom). UM-SAGA galaxies with M ⋆ < 10 7 M ⊙ are not included in this plot. Solid lines mark the median SFHs of the galaxies from GRUMPY (dark purple), UM-SAGA (purple), MARVEL (pink), and FIRE-2 (yellow). The corresponding filled regions denote the 16thto-84th percentile spread of the SFHs. For comparison, the median SFHs for field dwarf galaxies (blue) and Milky Way satellites (red) measured from resolved color-magnitude diagrams by Weisz et al. (2014) are shown.  \n<!-- image -->  \nlar physics models, including a range of mass loss rates, and these differences can lead to systematic discrepancies in galaxy M ⋆ measurements. We therefore use FSPS to self-consistently convert between the formed SFH and surviving stellar mass, using the same conversion for all four models. With the MIST isochrones (Choi et al. 2016), this corresponds to a typical differ-  \nce of ∼ 40 -50% between the formed and surviving M ⋆ .  \nWe also use FSPS to generate noiseless photometry in bands ranging from the UV through IR, listed in Table 1. Our choice of bands is motivated by the existence of wide-area surveys-GALEX in the UV (Martin et al. 2005; Bianchi et al. 2017), SDSS in the optical (Blanton et al. 2017; Doi et al. 2010), and WISE in the near-IR (Wright et al. 2010)-which, combined, provide multi-wavelength photometry for > 700 , 000 galaxies in the low-redshift universe (Salim et al. 2016). We then add realistic noise to the baseline 'noiseless' photometry in order to compare directly with observations. As a reference, we use the third data release of the Satellites Around Galactic Analogs survey (SAGA DR3; Mao et al. 2024), which identified 378 satellite galaxies around 101 Milky Way-mass analogs between 25-40 Mpc away. SAGA used imaging data from several public datasets-in particular, most of the optical grz photometry came from the DESI Imaging Legacy Surveys (Dey et al. 2019). To estimate observational noise for each galaxy in our sample, we first use a k-d tree to identify the SAGA satellite with the nearest r -band magnitude, g -r color, and r -z color. For each photometric band X , we then assign the photometric uncertainty σ X from the nearest-neighbor SAGA satellite to the simulated galaxy, and we perturb the simulated galaxy's noiseless photometry by adding a noise factor δ X randomly drawn from the Gaussian distribution δ X ∼ N (0 , σ X ).  \nFigure 2 illustrates the results of this FSPS modeling, showing the distribution of the simulated galaxies in ( g -r )-mass space. In general, the galaxies from all four models (particularly when the poorly-constrained < 10 7 M ⊙ UM-SAGA galaxies are removed) occupy a similar region of optical color space. However, each of the models described in Section 3.1 parameterizes galaxy evolution differently, which leads to some minor variations among the mock observations. For example, the GRUMPY galaxies are systematically redder than the UM-SAGA galaxies, particularly in the mass range 10 6 -8 M ⊙ . This is likely related to differences in the SFHs, as illustrated in Figure 3. When comparing galaxies in the mass range 10 6 -8 M ⊙ , the GRUMPY galaxies typically form more stars at earlier times, whereas the UM-SAGA galaxies have more extended SFHs and higher recent star formation rates, leading to bluer colors. We revisit potential implications of these differences on stellar mass calibrations in Section 7.1.", '4. RESULTS: STELLAR MASS RECOVERY USING PHOTOMETRIC CALIBRATIONS': "We now test how well different stellar mass measurement techniques can recover the 'true' stellar masses from mock observations, beginning with empirical photometric calibrations based on optical and IR fluxes.", '4.1. Literature calibrations': "As discussed in Section 2.1, many studies in the literature have developed optical-NIR colorM ⋆ /L relations to estimate galaxy stellar mass, using either samples of observed galaxies (e.g., Bell et al. 2003; Portinari et al. 2004; Schombert 2006; Taylor et al. 2011; Herrmann et al. 2016; Du et al. 2020) or theoretical SEDs (e.g., Zibetti et al. 2009; Gallazzi & Bell 2009; Into & Portinari 2013; Roediger & Courteau 2015; Kouroumpatzakis et al. 2023; Klein et al. 2024). Figure 4 illustrates how well several of these calibrations work on simulated lowmass ( M ⋆, true < 10 8 M ⊙ ) galaxies. Table 2 summarizes their performance, listing the mean and standard deviation of the residuals log( M ⋆, calib /M ⋆, true ) for low-mass galaxies from each model discussed in Section 3, as well as for all models combined.  \nThe optical colorM ⋆ /L relations (top three rows of Figure 4) show two primary features. First, most of the optical calibrations are able to recover M ⋆ of low-mass galaxies ( < 10 8 M ⊙ , left of the vertical red lines in the figure) within ∼ 0 . 2 dex of the true values (dashed horizontal gray lines in the figure). In this low-mass regime, nearly all optical calibrations have average residuals that are consistent with zero within uncertainties, suggesting that they are statistically unbiased. The main exception is the calibration by Klein et al. (2024) (right column, second from top panel in Figure 4), which appears to be overestimate M ⋆ for low-mass galaxies by 0 . 17 ± 0 . 11 dex on average. This is initially somewhat surprising, since Klein et al. (2024) developed this calibration using lowmass galaxies in the FIRE-2 simulation, which also make up the bulk of our sample of mock galaxies. Upon further inspection, this discrepancy is likely due to different assumptions about mass loss prescriptions. Klein et al. (2024) directly measured 'true' stellar masses from counting star particles in simulated FIRE-2 galaxies; in FIRE-2, the typical mass loss (i.e., the percent difference between the formed M ⋆ and the surviving M ⋆ of a galaxy) is ∼ 20 -30%. As described above, we instead obtain M ⋆, true for each simulated galaxy using FSPS with MIST isochrones (Choi et al. 2016), which produce a total mass loss of ∼ 40 -50%. As a result, the surviving M ⋆ of galaxies that Klein et al. (2024) measured directly from the FIRE-2 simulation are higher than the M ⋆, true used in this work, explaining the systematic overestimate. This suggests that uncertain physics in stellar models (i.e., mass loss rates) can lead to uncertainties  \nFigure 4. Comparisons between 'true' stellar masses of model galaxies ( M ⋆, true ) and stellar masses measured from mock observations using different empirical photometric calibrations ( M ⋆, calib ). Different colors and shapes represent mock observations obtained from different galaxy models. UM-SAGA (purple) galaxies with log( M ⋆ / M ⊙ ) < 7 are marked with small unfilled points to denote that their SFHs are not well constrained. The dashed horizontal gray lines mark ± 0 . 2 dex errors, which represent typical systematic uncertainties that are assumed for stellar mass measurements. The horizontal blue lines and shaded bands denote the means and standard deviations of the mass residuals for the low-mass ( M ⋆ < 10 8 M ⊙ ) and high-mass ( M ⋆ > 10 8 M ⊙ samples, separated by the vertical red line.) Note the expanded y -axis limits on the bottom two plots.  \n<!-- image -->  \nof ∼ 0 . 2 dex or more in M ⋆ measured from empirical colorM ⋆ /L relations.  \nAdditionally, all optical relations (top three rows of Figure 4) show clear systematic trends in the residuals log( M ⋆, calib /M ⋆, true ) as a function of stellar mass. As shown by the discontinuous blue horizontal lines  \nin the figure, the residuals for the low-mass galaxies ( < 10 8 M ⊙ , left of the red vertical lines) are consistently ≳ 0 . 15 dex lower than the residuals for high-mass galaxies. The existence of this trend as a function of mass is perhaps unsurprising, since many of these relations have been calibrated based on samples with limited  \nTable 2. Residuals between 'true' and measured masses, log ( M ⋆, calib /M ⋆, true ), for different photometric calibrations. The 'Bands' column lists the various photometric bands used in each calibration, and the bolded band is the reference band (e.g., when calculating M ⋆ /L g , g is the reference band). The values reported are the mean and standard deviation of the residuals for low-mass ( M ⋆, true < 10 8 M ⊙ ) galaxies.  \nmass ranges. It may arise because most of the literature relations treat M ⋆ /L as a linear function of optical color, which inherently assumes that M ⋆ ∝ L ; we will revisit this assumption when we develop new calibrations in the next section.  \nAs shown in the bottom two panels of Figure 4, the literature calibrations based on near-IR WISE bands perform more poorly than the optical calibrations. The Jarrett et al. (2023) W1 calibration (bottom left panel) has a pronounced downward trend in log( M ⋆, calib /M ⋆, true ) as log M ⋆, true decreases below 10 8 M ⊙ , which causes it to underestimate M ⋆ of low-mass galaxies by 0 . 73 dex on average (blue horizontal line on the left side of the plot). Meanwhile, the W1 -W2 calibration (bottom right panel) produces M ⋆ residuals with a large 1 σ dispersion of 1 . 46 dex for low-mass galaxies (blue shaded region on the left side of the plot).  \nThese trends are likely due to a combination of factors. The Jarrett et al. (2023) relation based on W1 is a quadratic function of log( L W1 ), which can produce strong systematic biases when extrapolated. The Jarrett et al. (2023) W1 -W2 calibration, on the other hand, shows dramatic scatter in mass residuals that is almost certainly due to uncertainties in the W2 band. For the low-mass ( < 10 8 M ⊙ ) galaxies in our sample, the assigned uncertainty in W2 magnitude (calculated by finding the nearest-neighbor SAGA galaxy in optical colors; see Section 3.2) is on average 6.2 times larger than the uncertainty in W1 magnitude. For high-mass galaxies, the W2 magnitude is only 3.3 times more uncertain than  \nthe W1 band, which explains why the scatter in the Jarrett et al. (2023) W1 -W2 relation decreases significantly above 10 8 M ⊙ .  \nWe also note that for all galaxies in our sample, the mock photometric data were produced assuming a single dust attenuation model: the Calzetti et al. (2000) attenuation law with τ V = 0 . 2. Changing the dust model may impact the stellar masses recovered by empirical calibrations. We test this by varying the dust attenuation law in two different ways: we either change the form of the dust law by using an SMC-like attenuation law from Gordon et al. (2003), or we change the total amount of dust by increasing the normalization of the Calzetti et al. (2000) law from τ V = 0 . 2 to up to τ V = 1 . 8. In both cases, we find that while altering the dust model may affect the average stellar mass residuals log( M ⋆, calib /M ⋆, true ) (i.e., the blue horizontal lines in Figure 4), these changes are typically smaller than the 1 σ dispersion in the residuals. Furthermore, changing the dust law has the same effect on both low- and high-mass galaxies, so the systematic trends in log( M ⋆, calib /M ⋆, true ) as a function of stellar mass remain unchanged.", '4.2. Updated stellar mass calibrations': 'As discussed in the previous section, many literature optical and NIR calibrations for M ⋆ /L appear to be systematically biased as a function of stellar mass. We can use our mock dataset to define new M ⋆ /L calibrations that avoid these biases.  \n<!-- image -->  \n<!-- image -->  \n<!-- image -->  \nFigure 5. Top panels: linear fits to stellar mass residuals as a function of absolute g -band magnitude (left) and M ⋆ /L W1 as a function of W1 luminosity (right), which are used to produce updated stellar mass calibrations (Equations 11 and 12). Note that the most massive galaxies ( M ⋆ > 10 10 M ⊙ ; open circles) are not included in these fits. Bottom panels: stellar mass residuals as a function of mass for the updated mass calibrations; axes, colors, and symbols are the same as in Figure 4.  \n<!-- image -->  \nFor the optical colors, we begin by attempting to fit M ⋆ /L g as a linear function of g -r color, as has previously been done in the literature; from this best fit, we calculate an initial estimate for stellar mass M ⋆, linear :  \nlog( M ⋆, linear / M ⊙ ) = 1 . 613+1 . 433(g -r) -0 . 4M g , 0 (10)  \nHowever, when we compute the residuals between the true stellar mass and this initial estimate (expressed as log( M ⋆, linear /M ⋆, true )), we find that these residuals are a strong function of g -band luminosity. As shown in the top left panel of Figure 5, we fit the residuals with a quadratic function; we find that this functional form performs better than other forms, including arctangent and piecewise linear functions. Note that we do not include the most massive galaxies ( M ⋆ > 10 10 M ⊙ ) in this fit, since their small sample size might strongly bias the fit. This yields a new calibration for stellar mass:  \nlog( M ⋆ / M ⊙ ) = 1 . 433(g -r) +0 . 00153 M 2 g, 0 -0 . 335 M g, 0 +2 . 072 (11)  \nHere, as before, all magnitudes have been k -corrected using g -r color (Chilingarian et al. 2010), and an abso-  \nlute solar g -band magnitude of 5.11 is assumed (Willmer 2018). The bottom left panel of Figure 5 shows the result of this updated calibration. This stellar mass calibration does not appear to have significant systematic residual trends as a function of mass; the galaxies below and above 10 8 M ⊙ have residuals that are consistent within 1 σ dispersion.  \nBy fitting for M ⋆ as a function of both g -r color and absolute g -band magnitude M g, 0 , we inherently break the assumption that M ⋆ ∝ L g . Future studies examining M ⋆ /L relations may want to consider this multivariable approach, which seems to avoid systematic biases in M ⋆ . However, we caution that given the small size of our simulated sample at low and high masses, this relation is likely to suffer from extrapolation errors outside of the magnitude range -17 < M g, 0 < -5 . 5 (corresponding to an approximate mass range of M ⋆ ∼ 10 5 . 5 -9 M ⊙ ).  \nFor the NIR colors, we focus on the W1 band, since relatively large uncertainties in the W2 band make it less useful as a mass indicator. We plot M ⋆ /L W1 as a function of L W1 and find that it is reasonably well fit by a linear function. This fit, shown in the top right panel  \nof Figure 5, yields an equation for stellar mass:  \nlog( M ⋆ / M ⊙ ) = 0 . 867 log(L W1 ) + 0 . 705 (12)  \nThis equation uses the same equations for L W1 and k -correction as described in the previous section for the Jarrett et al. (2023) calibrations. The result of this updated calibration are shown in the bottom right panel of Figure 5. While this NIR calibration has a slightly larger dispersion in stellar mass residuals ( -0 . 05 ± 0 . 11 dex for low-mass galaxies) than our updated optical calibration, it does not suffer from the same systematic biases as the NIR calibrations considered in the previous section. As before, we note that this calibration might perform poorly outside of the mass range M ⋆ ∼ 10 5 . 5 -9 M ⊙ , but a larger sample of simulated galaxies in this high-mass regime is needed to confirm this.  \nWe briefly comment on the impact of the assumed dust attenuation model on our new photometric calibrations. As with the literature calibrations, we find that if we vary the form or normalization of the dust attenuation law in the simulated galaxies, the average stellar mass residuals vary slightly. However, this effect is again statistically insignificant; regardless of the dust model assumed in the mock photometry, both Equations 11 and 12 produce stellar mass residuals that are approximately constant as a function of stellar mass.', '5. RESULTS: STELLAR MASS RECOVERY USING SED FITTING': "We now consider stellar masses estimated from SED fitting. In a Bayesian framework, the assumptions about prior distributions can have significant effects on the best-fit parameters, and SED fitting is no exception. In this work, we focus on four primary assumptions that may impact the recovery of M ⋆ for low-mass galaxies: (1) the form of the SFH, (2) the dust model, (3) the IMF, and (4) the choice of photometric bands. These choices are particularly impactful in the low-mass regime, since low-mass galaxies have stochastic SFHs that may not be described well by classical SFH models and may lead to stochastic sampling of the IMF. Furthermore, the (sometimes extreme) low-metallicity conditions in lowmass galaxies may affect dust production and evolution in ways that are not captured well by the dust models built for and largely calibrated on more massive galaxies. Many other assumptions are involved in SED modeling (e.g., the stellar model library used, the fitting method), but we defer discussion of these to a later work.  \nThe full priors used in all tests in this section are described in Appendix A. Table 3 summarizes the results of all tests in this section, including the mean and standard deviation of the residuals log( M ⋆, calib /M ⋆, true ) for  \nlow-mass galaxies from each model, as well as for all models combined.  \nAs with any fitting procedure, the best-fit model identified by Prospector may not actually fit the data well. To test the goodness-of-fit, we use the reduced chisquared statistic χ 2 red to compare the true synthetic SED with the SED prediction for the highest probability sample. We identify fits with χ 2 red > 10 as 'poor' fits, marked with x's in all plots in this section. 7 Table 3 also lists the relative fraction of poor fits for each test, since our goal is to identify which prior assumptions most successfully recover M ⋆ ; however, the 'poor' fits are not considered in any further analysis (including the other summary statistics listed in Table 3).", '5.1. Parametric vs. nonparametric SFHs': "We first test the effect of different assumptions about the shape of the SFH on stellar mass measurements. We modify some of the default values in Prospector to accommodate the lower masses and metallicities expected for low-mass galaxies (see Appendix A). For this test, we fix the dust attenuation, dust emission, and IMF to be the same as the input SEDs-we will modify these in the following sections. We also use all 7 photometric bands listed in Table 1, which tests the 'best-case' scenario for lowz galaxies in the SAGA and DESI surveys; we will investigate some of the effects caused by missing photometric coverage later in Section 5.4.  \nThe left column of Figure 6 presents a comparison of Prospector's parametric SFH templates: constant SFH, delayedτ SFH, and delayedτ with an additional instantaneous burst of star formation. We find that the assumption of a constant SFH (top left panel of Figure 6) leads to significant offsets for galaxies with M ⋆ < 10 8 M ⊙ : using a constant SFH model, Prospector can underestimate true stellar masses by up to ∼ 4 dex. This is not unexpected; Lower et al. (2020) also showed that a constant SFH can lead to similarly large offsets in M ⋆ for more massive galaxies. Additionally, we find that a constant SFH produces mostly 'poor' fits: ∼ 97% of the fits produce χ 2 red > 10. This is not surprising, since z ∼ 0 galaxies like the ones in our sample are thought to have quenched over time (e.g., Searle et al. 1973). A constant SFH is more appropriate for high-redshift starforming galaxies that are young enough to have maintained the same SFR over their short lifetimes.  \nFor lowz galaxies whose SFR has likely declined over time, a more realistic parametric SFH model is the  \nFigure 6. Same as Figure 4, but M ⋆, calib represents stellar masses measured using Prospector with various SFH priors. The left column shows the results of different parametric SFH priors: a constant SFH (top), a delayedτ ( ∼ t exp -t/τ ) SFH (middle), and a delayedτ SFH with a δ -burst added (bottom). The right column shows the results of non-parametric SFH priors: a Log M SFH (top), a continuity SFH (middle), and a Dirichlet SFH with a Dirichlet hyperparameter α = 0 . 7 (bottom). Bad fits ( χ 2 red > 10) are marked with x's. Note the expanded y -axis limits on the top plots.  \n<!-- image -->  \ndelayedτ model (left middle panel of Figure 6), in which SFR is described as a delayed declining exponential:  \nSFR( t ) ∝    0 t < T 0 ( t -T 0 ) e -t -T 0 τ t ≥ T 0 (13)  \nwhere T 0 is the delay time of the SFR (or the maximum stellar age of the galaxy) and τ is the e -folding time. For more massive galaxies ( M ⋆ > 10 8 M ⊙ ), the delayedτ model underpredicts the true stellar mass by ∼ 0 . 2 dex, in agreement with Lower et al. (2020). However, for lower-mass galaxies the systematic bias is much less dramatic; on average, after removing poor fits (which constitute only 3% of our sample), stellar masses are only underpredicted by -0 . 05 dex on average, with a dispersion of 0 . 14 dex. As shown in the bottom left panel of Figure 6, increasing the flexibility of the model by adding an instantaneous δ -function burst of star formation (characterized by the time of the burst and the frac-  \ntellar mass produced in the burst) further improves the recovery of M ⋆ for low-mass galaxies (the percentage of poor fits drops even further, to 1.7%) and reduces the dispersion in stellar mass residuals to 0 . 09 dex. Neither model shows a strong trend in residuals as a function of stellar mass below M ⋆ < 10 8 M ⊙ .  \nWe now consider 'nonparametric' SFH models, which are described not with analytic functional forms but with hyperparameters-e.g., some number of time bins (which are each assumed to have constant SFR), the bin widths, and other hyperparameters describing the relationships between the SFRs in the bins. The right column of Figure 6 shows the performance of three of the SFH models implemented in Prospector (Leja et al. 2017). For this test, we fix the number of bins to N bin = 10 for all SFH models (see Section 7.2 for further discussion of this choice). Following the procedure of Lower et al. (2020), the bins are logarithmically spaced  \nexcept for the two youngest bins, which in this work are fixed to span lookback times of 0 -10 Myr and 10 -100 Myr.  \nThe Log M model is the simplest nonparametric model, in which the free parameters are the masses formed in each of the fixed time bins. This SFH model is extremely flexible and can reproduce multiple bursts of star formation and quenching. However, as shown in the top right panel of Figure 6, the Log M model struggles to accurately recover stellar masses. Only 39.5% of the fits are 'good' ( χ 2 red < 10); of these 'good' fits, stellar mass is underestimated by -0 . 24 dex on average, with a dispersion of 0 . 32 dex.  \nNonparametric SFH templates with more constrained functional forms perform better than the Log M model. The Continuity model (right middle panel of Figure 6) fits directly for the difference in SFR between adjacent time bins, but explicitly weights against sharp transitions in SFR. Although low-mass galaxy SFHs are thought to be highly stochastic and bursty, the Continuity model performs remarkably well, with an average stellar mass residual of 0 . 01 dex and a dispersion of 0 . 09 dex. The Dirichlet model (bottom panel), like the Continuity model, also constrains the relationship between SFRs in each of the bins. In this model, the fractional sSFR for each time bin follows a Dirichlet distribution (Betancourt 2012), set by a 'concentration' parameter α : low α values weight toward bursty SFHs, while higher α weights toward smooth SFHs. We assume α = 0 . 7 following Lower et al. (2020), although we discuss the effect of changing α in the following section. The Dirichlet template is comparable to the Continuity model in terms of stellar mass recovery, with an average offset of -0 . 03 dex and a dispersion of 0 . 11 dex. Both the Continuity and the Dirichlet models return the highest number of 'good' fits of all the SFH models (with only 1.28% of all fits being discarded as 'poor' fits), and neither template produces a strong trend in the residuals as a function of stellar mass.  \nThe overall success of the Continuity and Dirichlet SFH models suggests that the additional flexibility provided by non-parametric SFHs improves stellar mass measurements of low-mass galaxies. However, too much flexibility (as in the Log M model) may under-constrain the fit, leading to poor SED fits and unreliable M ⋆ measurements. For the remainder of this work, we use the Dirichlet template as our fiducial SFH model; although the Continuity template performed similarly well, the Dirichlet model provides a more direct comparison to the work by Lower et al. (2020), as discussed in the next section.", '5.1.1. Comparison with literature results': "These results are a direct extension of the work on higher-mass galaxies by Lower et al. (2020), who demonstrated that for galaxies with M ⋆ ≳ 10 8 M ⊙ , assuming simple parametric SFHs can lead to systematic biases in measured stellar mass of up to ∼ 0 . 4 dex. Like this work, Lower et al. (2020) also used FSPS to generate mock photometry of simulated galaxies, then compared the 'true' masses with Prospector's estimates.  \nThere are some differences between the studies-most notably, Lower et al. (2020) used galaxies from a single large-volume hydrodynamic simulation (SIMBA; Dav'e et al. 2019), while we used several different models that all specifically aim to model lower-mass galaxies (including zoom-in simulations that have much finer mass resolution than SIMBA). Rather than using FSPS to apply a simple dust attenuation law while generating mock SEDs, as we have done here, Lower et al. (2020) used the more complex 3D radiative transfer code POWDERDAY to model a dust screen surrounding all stars. Finally, while Lower et al. (2020) considered the effect of applying 3% uncertainties to their mock photometry, we used realistic observational uncertainties that vary as a function of wavelength and for different galaxies.  \nDespite these slight differences in methods, many of our results are qualitatively in agreement. We find that parametric SFHs generally perform worse than nonparametric SFHs, often underestimating the true stellar mass. This is particularly apparent for the Constant SFH template: we find an large stellar mass offset of -1 . 09 ± 0 . 80 dex for low-mass galaxies (when poor fits are included), which is consistent with the offset of -0 . 48 ± 0 . 61 dex reported by Lower et al. (2020). Lower et al. (2020) also find that other parametric SFHs (Delayedτ , Delayedτ +burst) tend to underestimate M ⋆ , which is qualitatively consistent with the high-mass ( M ⋆ > 10 8 M ⊙ ) galaxies in our sample (right of the red vertical lines in the left middle and left bottom panels of Figure 6).  \nHowever, we find that for low-mass galaxies, parametric SFHs can still perform quite well: the Delayedτ and Delayedτ +burst models produce stellar mass residuals that are, on average, consistent with zero within their 1 σ dispersions. Additionally, Lower et al. (2020) find that the non-parametric Continuity SFH template overestimates stellar masses by 0 . 24 ± 0 . 15 dex, while we find that the Continuity template recovers masses that are much closer to the true mass (0 . 01 ± 0 . 09 dex for the low-mass galaxies). This is somewhat surprising: the parametric and Continuity templates prefer smoothly changing SFHs, but low-mass galaxies are expected to have burstier SFHs than massive galaxies.  \nThe exact reason for this discrepancy is unclear. One option is that short-duration bursts of star formation simply do not contribute significantly to the overall stellar mass of a low-mass galaxy. As a result, parametric and smooth SFHs may still recover a reasonably accurate stellar mass even if Prospector's SFH models are unable to capture the full behavior of the SFH. We defer a full test of this hypothesis to future work. In the meantime, we can at least confirm, as Lower et al. (2020) did, that non-parametric SFH models with some constraints (e.g., the Dirichlet and Continuity models) are the best among the Prospector SFH templates tested for stellar mass recovery. Finally, because the Dirichlet and Continuity models performed comparably for our sample, for convenience we choose the Dirichlet model as our fiducial SFH model, since Lower et al. (2020) also adopted the Dirichlet model as their default non-parametric SFH.", '5.1.2. Non-parametric SFH hyperparameters': "As described above, the non-parametric SFH models defined by Prospector are not truly 'non-parametric,' in that they must still be described by hyperparameters such as the number and width of time bins ( N bins ) of constant SFH. Ocvirk et al. (2006) suggest using age bins spaced logarithmically in time, and that no more than eight characteristic star-forming episodes can be recovered with high-quality optical spectra. We therefore test the effect of changing N bins for our fiducial SFH model (Dirichlet SFH with α = 0 . 7). We find that if N bins is too small, many of the fits are poor quality ( χ 2 red > 10); for example, with N bins = 3, 70% of the galaxies have poor fits. However, the 'good' fits are still able to recover stellar masses that are on average consistent with the true stellar masses. The number of poor fits can be reduced by increasing the number of bins, and the results for N bins = 6 , 10 , 12 are essentially the same.  \nWe also test the effect of varying the Dirichlet α parameter, which determines the burstiness of the SFH. Again, we find that if α is too small, the number of poor fits increases - if α = 0 . 2, 41 . 8% of the galaxies have poor fits - but the 'good' fits are still strong predictors of the true stellar mass. Increasing α reduces the number of poor fits, and α = 0 . 7 (the fiducial value) and α = 1 . 0 again return near-identical results.  \nThese tests suggest that there are some minimum values for N bins and α , below which Prospector begins to struggle to fit many low-mass galaxy SFHs. While determining these exact values is beyond the scope of this work (and is likely to depend on the sample of galaxies used - or in this case, the details of the simulations used to model our galaxies), we can at least safely as-  \nume that the fiducial hyperparameters discussed in the previous section N bins = 10 and α = 0 . 7 are reasonable assumptions.", '5.2. Dust prescription': "In the previous sections, we tested SFH assumptions by using idealized models in which the dust attenuation and emission exactly matched the input SEDs. Of course, this is not the case in reality. Having identified a fiducial SFH model, we can now address the impact of dust on our SED fits.  \nDust both attenuates and emits light from a galaxy, generally leading to a reddening of galaxy spectra. While significant work has been done to progress our understanding of both of these effects - see, e.g., Salim & Narayanan (2020) for a recent review of the dust attenuation law - the impact of dust on galaxy observations remains an open question, particularly in low-mass and low-metallicity galaxies. This is both because the physical properties of dust in low-metallicity environments are not well understood (e.g., Nanni et al. 2020; Galliano et al. 2021), and because dust, metallicity, and age have inherently degenerate effects on a galaxy's SED (and this degeneracy is also a complex function of galaxy mass and metallicity; e.g., Nagaraj et al. 2022).  \nIn this work, we do not attempt to choose the most physically correct global dust prescription for low-mass galaxies. Our goal is simply to test the robustness of Prospector's ability to recover M ⋆ of low-mass galaxies, given different assumptions about dust models. To do this, we vary the dust model used to produce the mock observations, then run Prospector with different parameterizations of the dust attenuation model. For these tests, we fix the SFH to our fiducial nonparametric Dirichlet SFH with α = 0 . 7.", '5.2.1. Dust attenuation': "We first test the impact of the assumed dust attenuation model. We consider both (1) the overall normalization of the attenuation law, and (2) the shape of the attenuation law.  \nAs described in Section 3.2, our primary set of mock observations are produced by assuming a Calzetti et al. (2000) attenuation law normalized to τ V = 0 . 2. To test how strongly the normalization of the dust attenuation curve affects stellar mass measurements, we synthesize additional sets of observations using a Calzetti et al. (2000) law normalized to τ V = 1 . 0 and 1 . 8. We then fit these mock SEDs using the Calzetti et al. (2000) law, but where τ V is allowed to vary as a free parameter. This parameterization of the dust attenuation law is a common assumption when fitting the SEDs of low-mass galaxies (e.g., Pandya et al. 2018; Greco et al. 2018),  \nFigure 7. Same as Figure 4, but M ⋆, calib represents stellar masses measured using Prospector with different dust attenuation laws. The rows represent mock observations with varying amounts of dust, produced by using the Calzetti et al. (2000) dust law with different normalization parameters τ V, true : 0.2 (top), 1.0 (middle), and 1.8 (bottom). The left column shows the results when assuming a Calzetti et al. (2000) dust law with the normalization parameter τ V allowed to vary as a free parameter. The right column shows the results when assuming a SMC bar attenuation law (Gordon et al. 2003). Bad fits ( χ 2 red > 10) are marked with x's.  \n<!-- image -->  \nalthough the exact form of the prior on τ V may vary; in this case, we assume a uniform prior for 0 < τ V < 2.  \nThe results of changing the dust normalization are shown in the left column of Figure 7, where the true value of τ V (i.e., the input value used to produce the mock observations) increases from the top row to the bottom row. Two primary features stand out: first, regardless of the true τ V of the mock observations, there is a increase in scatter in the recovered stellar mass residuals below M ⋆, true < 10 8 M ⊙ . This is visibly apparent in the left column of Figure 7; the blue horizontal bands to the left of the red vertical lines are larger than the blue horizontal bands to the right of the vertical lines. This increase in scatter is the direct result of allowing the dust normalization to vary: for low-mass galaxies, the dispersion in stellar mass increases from 0 . 09 dex when the dust prescription exactly matches the input dust law (bottom right panel of Figure 6) to [0 . 19 , 0 . 23 , 0 . 15] dex  \nwhen τ V is allowed to vary (for τ V, true = [0 . 2 , 1 . 0 , 1 . 8], respectively).  \nAdditionally, the average residuals for low-mass galaxies (the blue horizontal bands left of the red vertical lines in Figure 7) appear to be a function of τ V, true . The freeτ V models typically overestimate stellar masses (mean offset of 0 . 21 dex) when τ V, true = 0 . 2 and underestimate stellar masses (mean offset of -0 . 17 dex) when τ V, true = 1 . 8. This makes sense, since we stipulated in our prior assumptions that 0 < τ V < 2. As a result, if τ V, true = 0 . 2 ( τ V, true = 1 . 8), Prospector is more likely to overestimate (underestimate) τ V , and overestimating the total amount of dust will lead to an overestimate (underestimate) of the true stellar mass. To check this hypothesis, Figure 8 plots the residuals in stellar mass as a function of the residuals in τ V (i.e., the difference between the measured and 'true' τ V values). As expected, there is a clear trend between the two residuals:  \nFigure 8. Difference between measured and actual stellar mass log( M ⋆, calib /M ⋆, true ) as a function of the difference between measured dust normalization τ V, calib and actual dust normalization τ V, true . The dashed horizontal gray lines mark ± 0 . 2 dex uncertainties in stellar mass. UM-SAGA galaxies with log( M ⋆ / M ⊙ ) < 7 are marked with small unfilled points to denote that their SFHs are poorly constrained.  \n<!-- image -->  \noverestimating τ V ( τ V, calib -τ V, true > 0) is correlated with overestimating M ⋆ , and underestimating τ V is correlated with underestimating M ⋆ .  \nWe can also test the effect of changing the shape of the dust attenuation law by fitting our mock SEDs with the extinction law measured by Gordon et al. (2003) for the SMC bar, which has a much steeper UV-optical slope than the Calzetti et al. (2000) law. Rather than a single normalization parameter, as was used with the Calzetti et al. (2000) law, we now assume a dust screen geometry in which young stars are more strongly attenuated by their birth clouds, while older stellar populations are less strongly attenuated by diffuse cirrus dust. Prospector uses the Charlot & Fall (2000) parameterization to model how opacity varies as a function of stellar age t . We assume a power-law dust screen for young stars (with ages ≤ 10 7 y), while both young and old stars are attenuated by diffuse dust:  \nτ young = τ 1 ( λ/ 5500 ˚ A) -δ CF (14)  \nτ diff = τ 2 ( A ( λ ) /A ( V )) (15)  \nFollowing Conroy et al. (2010) and Charlot & Fall (2000), we set δ CF = 1 . 0, τ 1 = 1 . 0, and τ 2 = 0 . 3. Here, the diffuse dust attenuation A ( λ ) /A ( V ) is set by the SMC bar extinction law (Table 4 of Gordon et al. 2003).  \nThe results of changing the dust attenuation law are shown in the right column of Figure 7, where the SMC dust law is applied to mock SEDs produced using the Calzetti et al. (2000) law with varying τ V, true . As before, τ V, true traces the overall amount of dust. For galaxies with low amounts of dust ( τ V, true = 0 . 2; top row), the  \nexact form of the dust attenuation law does not significantly change the output; assuming the SMC dust law produces stellar masses that are as accurate as if the 'true' dust law was used (0 . 00 ± 0 . 09 dex, compared to -0 . 03 ± 0 . 09 dex when the dust law is fixed to match the input). However, as the total amount of dust increases, using the incorrect dust law more dramatically impacts the recovered stellar masses. For both τ V, true = 1 . 0 and τ V, true = 1 . 8, the SMC dust law underpredicts the stellar masses of M ⋆ < 10 8 M ⊙ galaxies by ∼ -0 . 3 dex on average. Furthermore, as the true value of τ V increases, the number of 'poor' fits increases, particularly for the most massive galaxies in our sample.", '5.2.2. Dust emission': "We also check the effect of the dust emission model on stellar mass measurements. We again fix the SFH to a Dirichlet prior, this time with the dust attenuation law fixed to be exactly the same as the input SEDs (Calzetti et al. 2000, with τ V = 0 . 2). In this test, we allow the dust emission parameters in the Draine & Li (2007) dust model-the minimum radiation field U min , the fraction of dust heated by starlight γ , and the PAH mass fraction q PAH -to vary. Compared to the same model in which these parameters were fixed (right bottom panel of Figure 6), this additional flexibility improves the fits for some galaxies, particularly those with stellar masses ≳ 10 8 M ⊙ ; the percentage of 'poor' fits decreases to 0.64%. However, for low-mass galaxies with M ⋆, true < 10 8 M ⊙ , varying the dust emission parameters only slightly changes the recovered stellar masses, and the overall mean and standard deviation of the residuals do not change at all.  \nTo summarize, the assumed dust model has a strong impact on both the systematic offset and the scatter in stellar masses recovered from SED fitting, and these effects are particularly apparent for low-mass ( M ⋆ < 10 8 M ⊙ ) galaxies:  \n- · Assumptions about the overall normalization of the attenuation law (i.e., the total amount of dust) affect the recovered stellar mass in two ways.\n- 1. Allowing dust normalization to vary increases the population-level scatter in recovered masses for low-mass galaxies.\n- 2. Any offset in the prior on dust normalization will lead to systematic offsets in the recovered stellar masses, because overestimating (underestimating) the dust normalization will lead to overestimating (underestimating) M ⋆ . Since most low-mass galaxies are likely  \nto have low metallicities and correspondingly low dust masses, a wide prior on the dust normalization is more likely to overestimate the amount of dust and consequently overestimate M ⋆ .  \n- · The shape of the dust attenuation law can also have an effect on the recovered stellar masses, but this only becomes a significant issue for dusty galaxies with τ V > 0 . 2; again, since most low-mass galaxies may have low dust content, this may not be a major issue in this mass regime.\n- · Assumptions about dust emission do not appear to significantly impact the results of SED fitting to UV through near-IR bands, although we note that it may be more relevant when mid- and far-IR photometry is included.", '5.3. Stellar initial mass function': 'We briefly consider the role of the stellar IMF in SED fitting. The IMF directly sets M ⋆ /L , so assuming the wrong IMF should simply add a constant offset to log M ⋆ . Fitting models with a Kroupa (2000) IMF to mock observations produced using a Chabrier (2003) IMF should underestimate stellar masses by ∼ 0 . 03 dex, while fitting a Salpeter (1955) IMF to a Chabrier (2003) IMF should overestimate stellar masses by ∼ 0 . 24 dex. We check this intuition by applying different IMFs in Prospector, assuming the fiducial Dirichlet SFH and dust parameters fixed to match the input SEDs.  \nAs expected, we find that changing the IMF leads to constant shifts in stellar mass residuals (see Table 3): the average offset shifts, but the 1 σ dispersion does not change. The magnitudes of these offsets are also smaller than expected. Compared to the true Chabrier (2003)IMF, assuming a Kroupa (2000) IMF changes the average offset in stellar mass by only 0 . 1 dex, while assuming a Salpeter (1955) IMF increases the stellar mass residuals by 0 . 14 dex on average.', '5.4. Photometric coverage': "In the previous subsections, we considered only the results of fitting to all 7 photometric bands in Table 1. In real lowz galaxy surveys, incomplete photometric coverage is common. We now investigate how well Prospector can recover stellar masses with realistic subsets of the 7 photometric bands. For this test, we run Prospector on our fiducial set of mock observations (produced using a Chabrier 2003 IMF and a Calzetti et al. 2000 attenuation law with τ V = 0 . 2). We aim to mimic a more 'realistic' usage of Prospector, in which the true dust attenuation law is unknown and must be described by  \none or more free parameters. However, some basic assumptions can be made: due to low-mass galaxies' low metallicities, one might expect these galaxies to have relatively little dust - indeed, recent observations suggest that low-mass galaxies in the nearby universe have 0 ≲ τ V ≲ 0 . 4 (e.g., Geha et al. 2024; Greco et al. 2018; Pandya et al. 2018). As discussed in Section 5.2, for galaxies with such low dust content, the exact form of the assumed dust attenuation law does not strongly impact the recovered stellar mass. We therefore assume a Calzetti et al. (2000) attenuation law with a more restricted uniform prior on the normalization parameter τ V : U (0 , 0 . 4).  \nWe then fit Prospector to three different sets of bands: (1) only optical grz bands, (2) optical grz with near-IR coverage (WISE W1 and W2), and (3) optical grz with limited near-UV and near-IR coverage (GALEX NUV and WISE W1). The first case is motivated by the DESI Imaging Legacy Surveys (Dey et al. 2019), an extremely wide-field ( ≈ 14 , 000 deg 2 ) grz survey. The second case is motivated by the fact that while the WISE satellite imaged the entire sky in four near-to-mid-IR bands at 3.4, 4.6, 12, and 22 µ m (W1, W2, W3 and W4; Wright et al. 2010; Cutri & et al. 2012), substantially deeper coverage is available for the two bluest bands, W1 and W2. Finally, the third case is motivated by the SAGA survey (Mao et al. 2024; Geha et al. 2024), which uses GALEX NUV fluxes to measure star formation rates.  \nThe results of this test are shown in Figure 9. As illustrated in the top panel, the optical grz bands alone are able to recover much of the stellar mass information for low-mass galaxies. While the grz bands produce a systematic trend in stellar mass residuals as a function of M ⋆ , similar to the results of the empirical g -r calibrations (Figure 4), on average the low-mass residuals are consistent with zero ( -0 . 13 ± 0 . 13 dex). Adding more bands slightly improves the M ⋆ estimates for low-mass galaxies: the average residuals decrease to -0 . 11 ± 0 . 12 dex if W1 and W2 are added, and to -0 . 05 ± 0 . 09 dex if NUV and W1 are added. The NUV+ grz +W1 test shows no systematic trend in residuals as a function of stellar mass. However, it is worth noting that the inclusion of the W2 band produces 'poor' ( χ 2 red > 10) fits, likely due to the large W2 photometric uncertainties in our sample. We recommend that future studies exercise caution in using uncertain photometric measurements when fitting low-mass galaxy SEDs.  \n6. BEST PRACTICES FOR MEASURING THE STELLAR MASSES OF LOW-MASS GALAXIES  \nTable 3. Residuals between 'true' and measured masses, log ( M ⋆, calib /M ⋆, true ), for different SED fitting tests. The values reported are the mean and standard deviation of the residuals for low-mass ( M ⋆, true < 10 8 M ⊙ ) galaxies. Only 'good' fits ( χ 2 red < 10) are included when computing the mean and 1 σ residuals; the percentage of 'poor' fits ( χ 2 red > 10) out of 645 total fits is also listed.  \nFigure 9. Same as Figure 4, but M ⋆, calib represents stellar masses measured using Prospector with different sets of photometric bands: optical grz only (top), optical grz and near-IR (middle), and optical grz with limited near-UV and near-IR coverage (bottom).  \n<!-- image -->  \nBased on our investigation, we now outline recommendations for measuring the stellar masses of low-mass galaxies (5 < log( M ⋆ / M ⊙ ) < 8).  \n- 1. If multi-wavelength SED fitting is not computationally or observationally feasible, we recommend using Equations 11 or 12 to calculate stellar masses from optical or NIR fluxes, respectively (Figure 5). Based on our simulated sample, these calibrations are less systematically biased as a function of mass than literature relations (Figure 4), and will recover stellar masses with 1 σ errors of ∼ 0 . 1 dex. This is comparable to the best results from multi-wavelength SED fitting, and may significantly outperform SED fitting with certain choices of assumptions (see below); however, we note that these calibrations should only be used in the magnitude range -17 < M g, 0 < -5 . 5 ( M ⋆ ∼ 10 5 . 5 -9 M ⊙ ) to avoid extrapolation errors.\n- 2. If multi-wavelength SED fitting is used, we recommend assuming a non-parametric SFH, rather than a parametric SFH (Figure 6). In a best-case scenario, SED fitting with a non-parametric SFH can produce stellar mass estimates with 1 σ -level dispersions of ∼ 0 . 1 dex.\n- 3. For SED fitting with a non-parametric SFH, we recommend using at least six logarithmically spaced age bins, and applying some restrictions to prevent extremely large variations between SFH bins (i.e., using a Dirichlet SFH template with α ≥ 0 . 7 or a Continuity SFH template).\n- 4. Allowing dust attenuation parameters to vary in an SED fit further increases the dispersion in M ⋆, calib -M ⋆, true (Figure 7) and may potentially introduce uncertainties in individual stellar mass estimates of up to ∼ 0 . 6 dex (Figure 8), depending on how badly the total dust content is overor under-estimated. For low-mass galaxies that are not likely to be very dusty, we therefore recommend using narrow prior distributions (e.g., a Calzetti et al. (2000) dust attenuation law with normalization 0 < τ V < 0 . 4) to avoid significantly overestimating the amount of dust.\n- 5. Some SED fitting assumptions do not appear to strongly impact stellar mass recovery. For example, dust emission parameters do not significantly affect M ⋆ measurements based on optical and near-IR photometry. Changing the assumed stellar IMF may introduce constant offsets to M ⋆ measurements, but these are relatively small ( < 0 . 2 dex) for low-mass galaxies.\n- 6. If attempting to choose a limited subset of photometric bands to observe, the optical ( grz ) bands alone provide the most information about M ⋆ for low-mass galaxies. While adding more bands can further improve stellar mass estimates, any bands with large observational uncertainties (e.g., WISE W2 in this study) may actually increase the uncertainty in M ⋆ .", '7.1. Differences between models': "A fundamental flaw of the tests presented in Sections 4 and 5 is our inherent assumption that simulated galaxies accurately represent real low-mass galaxies. We have attempted to mitigate this by using several different galaxy models that use a variety of approaches and physical assumptions, as described in Section 3. Here  \nwe discuss whether our results depend strongly on these models.  \nFor each of the empirical photometric calibrations discussed in Section 2.1, the mean residuals computed for the low-mass galaxies from each of the four models (GRUMPY, UM-SAGA, MARVEL, FIRE) typically agree within 1 σ (as seen in each row of Table 2). The only exceptions are the UM-SAGA sample, which sometimes has systematically higher stellar mass residuals than the other models (e.g., for the Klein et al. (2024) calibration, the UM-SAGA sample has an average stellar mass residual of 0 . 33 ± 0 . 03 dex, compared to 0 . 18 ± 0 . 08 dex for the GRUMPY sample and 0 . 16 ± 0 . 11 dex for the FIRE sample). This is because only the UM-SAGA galaxies with well-constrained SFHs (and M ⋆ > 10 7 M ⊙ ) are included when calculating average residuals. Not only is this a small sample of galaxies ( N = 7), it also means the UM-SAGA sample is biased to higher stellar masses than the samples from the other models. Because of the upward trend in log( M ⋆, calib /M ⋆, true ) as a function of stellar mass for the optical calibrations (and the Jarrett et al. 2023 W1 relation), the average residuals for the UM-SAGA sample are also higher.  \nA similar effect is apparent when more complex SED fitting is used. As shown in Table 3, the average stellar mass residuals log( M ⋆, calib /M ⋆, true ) are typically highest for the UM-SAGA galaxies and lowest for the GRUMPY galaxies. However, as with the photometric calibrations described above, this difference is not statistically significant (i.e., the average UM-SAGA and GRUMPY residuals are consistent within 1 σ errorbars) for many of the SED fits we performed. For the cases where this difference is > 1 σ , the discrepancy may be due to a number of non-physical factors, including the small sample size and bias to higher M ⋆ of the wellconstrained UM-SAGA sample. Alternatively, as noted in Section 3.2, the UM-SAGA galaxies are systematically bluer than the GRUMPY galaxies at a given mass due to more extended SFHs (Figures 2 and 3). Because more recent bursts of star formation 'outshine' older stellar populations (Narayanan et al. 2024; Haskell et al. 2024), Prospector may overestimate the stellar mass contributed by the older populations in the UMSAGA galaxies with higher recent SFRs. Larger sample sizes are needed to distinguish between these scenarios.  \nHowever, we emphasize that differences between the four models are relatively minor. Qualitatively, the trends in the stellar mass residuals log( M ⋆, calib /M ⋆, true ) for photometric calibrations (Figure 4) and SED fitting (Figures 6-9) are consistent among all four models. This overall agreement implies that our key conclusions - in  \nparticular, the best practices outlined in Section 6 are robust to model-specific assumptions.", '7.2. Additional systematics': 'We now discuss systematic effects that may impact the applicability of our results. Some of these involve physics that have been simplified or ignored in modeling low-mass galaxies, while other systematic effects are due to our experimental setup (i.e., fixing certain variables, choice of SED fitting code).', '7.2.1. Physical factors': "For simplicity, we have ignored a number of physical factors in our modeling of low-mass galaxies.  \nFirst, we have overlooked the potential contribution of active galactic nuclei (AGN) to low-mass galaxies (see, e.g., the recent review of Reines 2022). None of the model galaxies used in this work include AGN, and we fixed the AGN contribution to zero when performing SED fits using Prospector. In real galaxy surveys, the AGN contribution may not be exactly zero. This can lead to systematic offsets in M ⋆ measurements, since ignoring the contribution of an AGN to a galaxy's light can lead to an overestimate of the stellar mass. The magnitude of this effect is unclear and highly dependent on how AGN are modeled: some studies find that AGN in the local universe do not strongly impact M ⋆ measurements from SED fitting (e.g., Leja et al. 2018; Thorne et al. 2022), while others suggest that some AGN models (particularly for high-luminosity AGN) may lead to M ⋆ uncertainties of up to ∼ 0 . 5 dex (e.g., Buchner et al. 2024; Ciesla et al. 2015).  \nMost of this literature has focused on higher-mass ( M ⋆ > 10 8 M ⊙ ) galaxies, so it is even less clear how AGN contributions might affect M ⋆ measurements in the low-mass regime. Fortunately, the fraction of lowmass galaxies in the nearby universe with AGN is expected to be quite low. Recent estimates suggest that the AGN occupation fraction is ≲ 0 . 5% in galaxies with M ⋆ ≲ 10 9 out to redshifts z < 0 . 1 (e.g., Birchall et al. 2020; Mezcua et al. 2018), so the results outlined in this work are likely valid for the vast majority of local lowmass galaxies. However, it is difficult to overstate just how little we know about the effects of AGN on lowmass galaxies. Additional modeling of AGN in low-mass galaxies is therefore needed to better constrain uncertainties on measurements of physical properties.  \nAnother physical phenomenon that may impact lowmass galaxies is quenching from reionization. We briefly consider how different assumptions about reionization might affect the stellar masses of the mock galaxies in our sample. The galaxy models described in Section 3 parameterize reionization in a variety of ways. In  \nthe semi-analytic model GRUMPY, UV heating from reionization is parameterized as a suppression factor on gas inflow, while in the zoom-in simulations MARVEL and FIRE, reionization is directly implemented as a photoionizing UV background (described respectively in Haardt & Madau 2012; Faucher-Gigu'ere et al. 2009). 8 While a full discussion of all the differences among these models is beyond the scope of this paper, we focus on one major difference in the timing of reionization in each model. The redshift of reionization is set to z reioniz = 6 in GRUMPY, z reioniz ∼ 15 in MARVEL, and z reioniz ∼ 10 in FIRE. This may have a direct effect on stellar mass: to first order, reionization occurring later leaves more time to form stars before quenching, leading to higher M ⋆ .  \nTo check this, we would ideally look at the galaxies in our sample with M ⋆ ≲ 10 6 M ⊙ , since these are thought to be the most strongly impacted by reionization (e.g., Wheeler et al. 2019). However, our sample in this mass regime is limited - both GRUMPY and MARVEL suffer from small sample sizes in the 10 4 -6 M ⊙ mass range ( N = 5 and N = 1, respectively, as shown in the top panel of Figure 3), and only a few GRUMPY galaxies have stellar masses < 10 5 M ⊙ (Figure 2) - making it difficult to identify any differences among simulations. A larger sample of model galaxies is needed to fully test whether stellar mass recovery depends on z reioniz . Additionally, we note that all models in this work assume uniform reionization for all galaxies. A potentially more realistic model of 'patchy' reionization (e.g., Pentericci et al. 2014) would imply that low-mass galaxies in denser environments should undergo reionization quenching first; in other words, z reioniz may also be a function of a galaxy's environment. As described in Section 3, nearly all of the model galaxies in our sample are satellite galaxies. Testing the effect of uniform vs. patchy reionization would require not only a larger sample of < 10 6 M ⊙ galaxies, but also a sample that spans a wide range of environments (satellite, field, and isolated).", '7.2.2. Experimental setup': "Other systematic effects may be due not to physical factors, but our choices when setting up our tests of stellar mass measurement techniques.  \nFor example, throughout this work, we assumed the redshifts of our mock galaxies were exactly known. While this is a reasonable assumption for galaxy surveys  \nthat obtain precise spectroscopic redshifts, it may have a larger impact on stellar mass measurements obtained from purely photometric surveys of low-mass galaxies. Over- or underestimating a galaxy's redshift will produce systematic offsets in stellar mass. The magnitude of this effect can be estimated with simple scaling relations: for z ≪ 1, z ∝ d and luminosity L ∝ Fd 2 (where d is luminosity distance and F is observed flux). Assuming M ⋆ ∝ L (which is generally reasonable, as shown by the empirical calibrations in Section 2), this leads to log( M ⋆, calib /M ⋆, true ) ≈ 2 log( z calib /z true ). For a galaxy with z true ≈ 0 . 1, typical photometric redshift errors of ∼ 0 . 03 (Salvato et al. 2019) can produce stellar mass residuals of ∼ 0 . 3 dex. This is an appreciable uncertainty, particularly when compared to the 1 σ dispersions of ∼ 0 . 1 dex that we found for our optical/NIR photometric calibrations and SED fitting (Section 6). 9  \nFurthermore, for most of the SED fitting tests in Section 2.2, redshift was not the only parameter we held fixed. We treated the SFH prior, dust attenuation, dust extinction, and IMF as independent assumptions. By only varying the assumption under investigation and fixing all other parameters, we may have ignored covariant effects. The true uncertainties in stellar mass likely depend on how assumptions about the SFH, dust, and IMF are combined.  \nFinally, as discussed in Section 2.2, we chose Prospector as our SED fitting code because it uses FSPS to create galaxy models. Since our mock observations were produced using FSPS, our experiment is therefore a direct 'apples-to-apples' comparison that tests, to first order, the SED fitting mechanism rather than other systematic effects that may depend on the exact choice of SED fitting code (e.g., stellar model libraries, dust parameterization, SED fitting algorithm). However, when fitting the SEDs of real galaxies, there is no way to select a code that exactly matches the input SED, so it is important to understand the systematic uncertainties introduced by different SED fitting codes. Pacifici et al. (2023) recently investigated the effects of using 14 different SED fitting codes on high-mass ( > 10 8 M ⊙ ) galaxies at z ∼ 1 and z ∼ 3 and found that all codes return stellar mass distributions that are qualitatively similar. We plan to extend this test to the low-mass regime in a future study.", '8. SUMMARY': "As we enter an era of unprecedentedly large extragalactic surveys, we are beginning to observe statistical samples of low-mass galaxies across a range of environments. Accurately estimating the physical properties of these galaxies is critical if we want to understand the physics that drives their evolution. However, many of the methods used to measure galaxy properties from integrated light have primarily been developed for and tested on massive ( > 10 8 M ⊙ ) galaxies. In this study, we have attempted to verify several of these methods for low-mass galaxies.  \nIn particular, we have tested how well different stellar mass measurement techniques can recover the masses of low-mass galaxies from integrated UV/optical/nearIR photometry. Our test sample included 469 simulated galaxies from four different models. Using multiple types of simulations - semi-analytic models, empirical halo-galaxy connection models, and two different hydrodynamic simulations with varying resolutions and prescriptions for feedback - makes our results more robust to variations among different simulations. For each simulated galaxy, we generated mock photometry using the stellar population synthesis code FSPS.  \nWe then tested empirical colorM ⋆ /L calibrations and showed that many literature relations produce systematic trends in stellar mass residuals that are only apparent at the low-mass end. This is particularly true for calibrations based only on near-IR color, which can produce stellar mass errors of > 1 dex for galaxies below 10 8 M ⊙ . We have provided updated prescriptions for stellar mass based on g -r color (Equation 11) and WISE 3.4 µ mluminosity (Equation 12) that reduce these systematic biases. These new calibrations can recover stellar masses with minimal offsets (0 . 04 ± 0 . 08 dex and -0 . 05 ± 0 . 11 dex, respectively).  \nWe also tested different assumptions that go into multi-wavelength SED fitting. In agreement with previous studies of high-mass galaxies (Lower et al. 2020), we found that non-parametric SFH models generally perform better than parametric SFH models when measuring stellar masses of low-mass galaxies. While parametric SFHs can underestimate stellar mass by as much as ∼ 0 . 4 dex, under ideal conditions non-parametric SFH templates can recover M ⋆ of low-mass galaxies with offsets of -0 . 3 ± 0 . 11 dex. Assumptions about dust attenuation parameters introduce larger uncertainties in M ⋆ of low-mass galaxies: over- or under-estimating total dust content can lead to significant ( ∼ 0 . 6 dex) overor under-estimates of stellar mass. The shape of the dust attenuation law may also impact stellar mass estimates, though more work is needed to better understand these effects; in this study we only compared two simple atten-  \ndels with different UV-optical slopes, and we did not consider other features in the attenuation curve (e.g., the UV bump; Salim & Narayanan 2020). Our results are summarized in Section 6, in which we have laid out recommendations for measuring M ⋆ of low-mass galaxies from integrated photometry.  \nThis study represents an initial step towards a full re-evaluation of measurement techniques based on the integrated light of low-mass galaxies. Much remains to be done, including further tests of SED fitting. As discussed in Section 7.2, this work aimed to test the assumptions that go into a single SED fitting code (Prospector, which uses the same SPS models that were used to generate our mock observations). This work also focused on stellar mass as a 'zeroth-order' property that traces a galaxy's integrated history. In the future, we plan to compare other SED and spectral fitting codes, as well as to investigate other galaxy properties-star formation rates and histories, metallicities, AGN properties-that can reveal even more information about the physics underlying low-mass galaxy evolution.  \nMAdlR acknowledges the financial support of the Stanford Science Fellowship while writing this paper. This research has made use of NASA's Astrophysics Data System Bibliographic Services.  \nThere are many communities without whom this work would not have been possible. We acknowledge that this work is rooted in Western scientific practices and is the material product of a long and complex history of settler-colonialism. MAdlR wishes to recognize her status as a settler on the unceded homelands of the Pocumtuc Nation. We hope to work toward a scientific practice guided by pono and a future in which we all honor the land.  \nThis research used data from the SAGA Survey (Satellites Around Galactic Analogs; sagasurvey.org). The SAGA Survey is a galaxy redshift survey with spectroscopic data obtained by the SAGA Survey team with the Anglo-Australian Telescope, MMT Observatory, Palomar Observatory, W. M. Keck Observatory, and the South African Astronomical Observatory (SAAO). The SAGA Survey also made use of many public data sets, including: imaging data from the Sloan Digital Sky Survey (SDSS), the Dark Energy Survey (DES), the GALEX Survey, and the Dark Energy Spectroscopic Instrument (DESI) Legacy Imaging Surveys, which includes the Dark Energy Camera Legacy Survey (DECaLS), the Beijing-Arizona Sky Survey (BASS), and the Mayall z-band Legacy Survey (MzLS); redshift catalogs from SDSS, DESI, the Galaxy And Mass Assem-  \nbly (GAMA) Survey, the Prism Multi-object Survey (PRIMUS), the VIMOS Public Extragalactic Redshift Survey (VIPERS), the WiggleZ Dark Energy Survey (WiggleZ), the 2dF Galaxy Redshift Survey (2dFGRS), the HectoMAP Redshift Survey, the HETDEX Source Catalog, the 6dF Galaxy Survey (6dFGS), the Hectospec Cluster Survey (HeCS), the Australian Dark Energy Survey (OzDES), the 2-degree Field Lensing Survey (2dFLenS), and the Las Campanas Redshift Survey (LCRS); HI data from the Arecibo Legacy Fast ALFA Survey (ALFALFA), the FAST all sky HI Survey (FASHI), and HI Parkes All-Sky Survey (HIPASS); and compiled data from the NASA-Sloan Atlas (NSA),  \nthe Siena Galaxy Atlas (SGA), the HyperLeda database, and the Extragalactic Distance Database (EDD). The SAGA Survey was supported in part by NSF collaborative grants AST-1517148 and AST-1517422 and Heising-Simons Foundation grant 2019-1402. SAGA Survey's full acknowledgments can be found at https:// sagasurvey.org/ack/.  \nSoftware: Prospector (Johnson et al. 2021; Leja et al. 2019), python-FSPS (Johnson et al. 2023), FSPS (Conroy et al. 2009; Conroy & Gunn 2010), Matplotlib (Hunter 2007), Seaborn (Waskom 2021), Astropy (Robitaille et al. 2013), Scipy (Jones et al. 2001)", 'A. PROSPECTOR PRIORS': 'For the sake of reproducibility, the full priors for all Prospector runs used in this study are described in Table 4.', 'REFERENCES': 'Agertz, O., Kravtsov, A. V., Leitner, S. N., & Gnedin, N. Y. 2013, ApJ, 770, 25  \nBruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000  \nAkeson, R., Armus, L., Bachelet, E., et al. 2019, arXiv e-prints, arXiv:1902.05569  \nAzartash-Namin, B., Engelhardt, A., Munshi, F., et al. 2024, ApJ, 970, 40  \nBehroozi, P., Wechsler, R. H., Hearin, A. P., & Conroy, C. 2019, MNRAS, 488, 3143  \n- Behroozi, P. S., Wechsler, R. H., & Wu, H.-Y. 2013, ApJ, 762, 109  \nBell, E. F., & de Jong, R. S. 2001, ApJ, 550, 212  \nBell, E. F., McIntosh, D. H., Katz, N., & Weinberg, M. D.  \n2003, ApJS, 149, 289. https://doi.org/10.1086/378847  \nBellovary, J. M., Hayoune, S., Chafla, K., et al. 2021,  \nMNRAS, 505, 5129  \nBetancourt, M. 2012, in American Institute of Physics Conference Series, Vol. 1443, Bayesian Inference and Maximum Entropy Methods in Science and Engineering: 31st International Workshop on Bayesian Inference and Maximum Entropy Methods in Science and Engineering, ed. P. Goyal, A. Giffin, K. H. Knuth, & E. Vrscay (AIP), 157-164  \nBianchi, L., Shiao, B., & Thilker, D. 2017, ApJS, 230, 24 Birchall, K. L., Watson, M. G., & Aird, J. 2020, MNRAS, 492, 2268  \n- Blanton, M. R., Bershady, M. A., Abolfathi, B., et al. 2017, AJ, 154, 28  \nBrinchmann, J., Charlot, S., White, S. D. M., et al. 2004, MNRAS, 351, 1151  \n- Buchner, J., Starck, H., Salvato, M., et al. 2024, arXiv e-prints, arXiv:2405.19297\n- Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682  \nCarlsten, S. G., Greene, J. E., Beaton, R. L., Danieli, S., & Greco, J. P. 2022, ApJ, 933, 47  \nChabrier, G. 2003, PASP, 115, 763  \nCharlot, S., & Fall, S. M. 2000, ApJ, 539, 718  \n- Chilingarian, I. V., Melchior, A.-L., & Zolotukhin, I. Y. 2010, MNRAS, 405, 1409  \nChoi, J., Dotter, A., Conroy, C., et al. 2016, ApJ, 823, 102 Christensen, C., Quinn, T., Governato, F., et al. 2012, MNRAS, 425, 3058  \nChristensen, C. R., Brooks, A. M., Munshi, F., et al. 2024, ApJ, 961, 236  \n- Ciesla, L., Charmandaris, V., Georgakakis, A., et al. 2015, A&A, 576, A10\n- Collins, M. L. M., & Read, J. I. 2022, Nature Astronomy, 6, 647  \nConroy, C. 2013, ARA&A, 51, 393  \n- Conroy, C., & Gunn, J. E. 2010, FSPS: Flexible Stellar Population Synthesis, Astrophysics Source Code Library, record ascl:1010.043, , , ascl:1010.043  \nConroy, C., Gunn, J. E., & White, M. 2009, ApJ, 699, 486. https://doi.org/10.1088/0004-637X/699/1/486  \nConroy, C., White, M., & Gunn, J. E. 2010, ApJ, 708, 58 Cutri, R. M., & et al. 2012, VizieR Online Data Catalog, II/311.  \nhttps://ui.adsabs.harvard.edu/abs/2012yCat.2311....0C  \nTable 4. Prior distributions used to fit SED models.  \nWetzel, A. R., Hopkins, P. F., Kim, J.-h., et al. 2016, ApJL, 827, L23  \n- Wheeler, C., Hopkins, P. F., Pace, A. B., et al. 2019, MNRAS, 490, 4447\n- Willmer, C. N. A. 2018, ApJS, 236, 47\n- Wise, J. H., Abel, T., Turk, M. J., Norman, M. L., & Smith, B. D. 2012, MNRAS, 427, 311\n- Woo, J., Courteau, S., & Dekel, A. 2008, MNRAS, 390, 1453\n- Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868.\n- https://doi.org/10.1088/0004-6256/140/6/1868\n- Zibetti, S., Charlot, S., & Rix, H.-W. 2009, MNRAS, 400, 1181'}