Systematically Measuring Ultra-diffuse Galaxies (SMUDGes). VI. Nuclear Star Clusters

We present our photometric search for potential nuclear star clusters (NSCs) in ultra-diffuse galaxies (UDGs) as an extension of the SMUDGes catalog. We identify 325 SMUDGes galaxies with NSCs and, from the 144 with existing distance estimates, identify 33 NSC hosts as UDGs (μ 0,g ≥ 24 mag arcsec−2, r e ≥ 1.5 kpc). The SMUDGes with NSCs lie on the galaxy red sequence, satisfy the relationship between NSC and host galaxy stellar masses, have a mean NSC stellar mass fraction of 0.02 but reach as high as 0.1, have NSCs that are displaced from the host center with a standard deviation of 0.10r e , and weakly favor higher-density environments. All of these properties are consistent with previous results from higher surface brightness galaxy samples, allowing for at most a relatively weak dependence of NSC behavior on host galaxy surface brightness.


INTRODUCTION
The origin of massive, compact stellar populations in galaxies, whether those are globular clusters (GCs), massive black holes, or nuclear star clusters (NSCs; Caldwell 1983;Binggeli et al. 1984;Caldwell & Bothun 1987;Bothun & Mould 1988), remains poorly understood.Some models envision their formation through violent, extreme episodes of star formation (e.g., Mihos & Hernquist 1994;Bekki & Couch 2001;Kravtsov & Gnedin 2005;Kruijssen et al. 2012;Renaud 2018), but such episodes may seem somewhat less likely in the low surface brightness galaxies that are the focus here.Other models relate the different population classes to each other, such as those that posit that NSCs form from the infall and merger of GCs (e.g., Tremaine 1976;Gnedin et al. 2014;Sánchez-Salcedo & Lora 2022;Modak et al. 2023) and those where central massive black holes form from the dynamical collapse of NSCs (e.g., Begelman & Rees 1978;Miller & Hamilton 2002;Antonini et al. 2015).Again, processes that may be common and relevant in massive, high surface bright-Corresponding author: Mika Lambert mlambert43@arizona.eduness galaxies, such as dynamical friction, perhaps play a diminished role in the low mass, low surface brightness galaxies.Commonality of features makes it attractive to link these populations into a coherent scenario (e.g., Wehner & Harris 2006;Rossa et al. 2006;Ferrarese et al. 2006;Fahrion et al. 2021).
At least regarding NSCs, basic constraints on any of these scenarios include the rates at which NSCs are found, the relationship between the NSC and host stellar masses and stellar populations, any connection to the host galaxy morphology, and the alignment of the NSC and its host's dynamical center -all specified as a function of the relevant host galaxy properties.One key such property might be the host's central surface brightness, which presumably reflects the degree to which dissipation has concentrated matter toward the galaxy center, where an NSC or a central black hole would reside.
The richness of some of these constraints is already evident, for example, in the NSC occupation fraction (the fraction of galaxies that host an NSC).The occupation fraction varies with galaxy mass, rising and then falling as one proceeds from lower to higher mass galaxies (Neumayer et al. 2020;Hoyer et al. 2021).This behavior is quite distinct from that of the number of GCs in a galaxy, which varies proportionally with the host galaxy mass (Blakeslee et al. 1997;Burkert & Forbes 2020) even to low masses (Forbes et al. 2020;Zaritsky 2022).One might naively have expected a close correspondence between the rate at which GCs and NSCs are found if NSCs are indeed formed from merged GCs, but the observed difference may highlight how certain details of the formation physics, such as the amplitude of dynamical friction, relate to the host galaxy properties.One can test scenarios along these lines and attempt to reproduce the NSC occupation fractions (e.g., Lotz et al. 2001;Capuzzo-Dolcetta & Mastrobuono-Battisti 2009).
The incidence rate of GCs and the NSC occupation fraction are less well determined as a function of host galaxy surface brightness, but there are indications of a dependence of the NSC occupation fraction (Binggeli et al. 2000;Lim et al. 2018).Of course, the incidence rate is only one of various properties to explore.One could imagine that the typical mass of an NSC varies with host mass.In fact, such a trend has been observed (e.g., Balcells et al. 2003;Ferrarese et al. 2006;Turner et al. 2012;Scott & Graham 2013;den Brok et al. 2014), and we are tempted to ask whether there is an analogous relation with the host surface brightness.
Here we begin our exploration of the SMUDGes (Zaritsky et al. 2019(Zaritsky et al. , 2021(Zaritsky et al. , 2022(Zaritsky et al. , 2023) ) set of ultradiffuse galaxy (UDG) candidates to explore the nature of NSCs in low surface brightness galaxies.The value of UDGs to this topic is that they include the most massive, low surface brightness galaxies known (cf.van Dokkum et al. 2015) and thus may help us disentangle the roles of mass vs. surface brightness in shaping NSC properties.The value of the SMUDGes sample is that it is large and spans all environments, enabling us to also explore the possibility of a dependence of NSC properties on environment.This work, in which we focus on the identification of NSCs in SMUDGes galaxies is followed closely by Khim et al. (2023, in prep.) in which we extend the analysis to galaxies of somewhat higher surface brightness, but similar stellar masses, provided by the Sloan Digital Sky Survey (SDSS; Kollmeier et al. 2017) and images from the DESI Legacy Survey (Dey et al. 2019).That work will enable us to place the SMUDGes galaxies in a wider context, without the additional challenge of comparing across studies with disparate image quality and analysis methodology.
Uniform, well-defined criteria are essential for comparisons of NSC properties across host mass, surface brightness, or environment.As we will show, the data characteristics and criteria used to identify NSCs lead to large variations in sample properties.NSCs are broadly defined to be dense and massive stellar agglomerations that reside in or near the centers of galaxies that are brighter than the extrapolated surface brightness pro-file of the inner region of that galaxy (Neumayer et al. 2020).However, consistently-applied, quantitative criteria do not exist for any of these defining characteristics.The purity and completeness of samples therefore varies among studies.Such (potential) differences between studies raise questions about any comparisons one would wish to pursue between, for example, cluster (Lim et al. 2018) and field samples (Carlsten et al. 2022) of galaxies with NSCs.In this particular example, both studies find occupation fractions ∼ 20%, but possible systematic differences leave open the question of whether this agreement is physically meaningful or fortuitous.Other studies have investigated the prevalence of NSCs in both early and late-type galaxies (Côté et al. 2006;Georgiev & Böker 2014) and in dwarf galaxies (Carlsten et al. 2022), but again comparisons among them remain challenging.
As with every study, this one too has its weaknesses and strengths.Among the weaknesses relative to existing NSC studies is that we are working with shallower, lower resolution images than the state-of-the-art (e.g., Hubble Space Telescope images have been used that have a ∼ 3 mag deeper point source magnitude than our images; Lim et al. 2018).Thus, we only probe the bright end of the NSC luminosity function for the majority of our hosts, potentially suffer greater contamination, and obtain more uncertain photometric parameters for the NSCs themselves.On the other hand, the strength of this, and our sister study (Khim et al. 2023, in prep), is that our galaxy sample spans mass, environment in low luminosity, low surface brightness galaxies with consistent classification, thereby simplifying comparisons.In §2, we describe our methodology.In §3, we discuss the results as follows: 1) an NSC classification catalog for the entire SMUDGes sample; 2) our constraints on the relative co-centricity of NSCs and their hosts; 3) NSC properties and their relationship to host galaxy properties and; 4) any relation to the host galaxy environment.We use a standard WMAP9 cosmology (Hinshaw et al. 2013), although the results are insensitive to different choices of cosmological parameters at the level of current uncertainties, and magnitudes are from SDSS/DESI and are thus on the AB system (Oke 1964;Oke & Gunn 1983).

The Data
We begin with the 6,805 visually-confirmed UDG candidates in the SMUDGes catalog (Zaritsky et al. 2023).These candidates were selected to have a low central surface brightness in the g-band, µ 0,g ≥ 24 mag arcsec −2 , and a large effective radius on the sky, r e ≥ 5.3 arcsec.In addition, we now impose a stricter color criterion than in the original work to remove likely background interlopers (∼ 0.2 mag redder than the red sequence) and candidates with unphysical blue colors (0 < g − r < 0.8; which removes 225 galaxies from the sample) and an angular size criterion to remove nearby galaxies that are unlikely to be UDGs and which corresponds to half our extracted image size (r e > 26 arcsec; which removes 38 additional galaxies).We retain 6,542 galaxies to analyze.The angular size cuts help our model fitting (see §2.2) by eliminating images that are either too small relative to our resolution or too large relative to the extracted image.These cuts do not easily translate to criteria on physical size because our candidates span a range of distances.
Importantly, for a study of NSCs in low surface brightness galaxies, the estimate of the central surface brightness of the galaxy used to define the catalog was calculated using Sérsic model fitting where high surface brightness objects, including any potential NSC, were masked (Zaritsky et al. 2022(Zaritsky et al. , 2023)).This aspect of the parent survey, and also whether a limit on the value of the Sérsic n value is imposed (SMUDGes imposes n < 2 in the fitting), can affect the sample selection and care must be taken if comparing results drawn from different catalogs of low surface brightness galaxies.
For our photometric analysis, we extract 200×200 pixel (52.4 × 52.4 arcsec) r-band images of each of the 6,542 candidates from the 9th data release (DR9) of the Legacy Survey (Dey et al. 2019).These cutouts provide a sufficiently large field of view to include adequate background coverage.

Preparing for Model Fitting
To explore the morphology of each UDG candidate further than what was done in the SMUDGes catalog papers, and to determine if there is evidence of an NSC, we use the photometric model fitting software package GALFIT (Peng et al. 2010).We take three preparatory steps before fitting any model.
First, to assess the nature of a possible unresolved source near the center of each UDG candidate, we need the point spread function (PSF) of each image.NSCs in our images are unresolved because their half light radii are typically < 10 pc (Böker et al. 2004;Turner et al. 2012;Georgiev & Böker 2014), which corresponds to an angular size < 1 arcsec for distances > 2 Mpc.For each UDG candidate, we adopt the PSF model provided in the Legacy Survey for the relevant image.
Second, we generate an image of the pixel-by-pixel uncertainties, a σ-image, that is used to assess the likelihoods of the models GALFIT produces.We calculate the σ-image using the inverse-variance image provided by the Legacy Survey that was calculated using the image stack contributing to each pixel.
Finally, we isolate each UDG candidate from any surrounding bright sources that could affect the modeling.This is a complicated, iterative process that we describe in more detail below when discussing our model fitting procedure.

Selecting Among Models
To identify NSCs in UDG candidates, we first assess whether there is evidence for a concentration of light beyond what can be described by a Sérsic model with index n ≤ 2. If there is, then we assess the nature of that excess component.Because our targets are extremely diffuse and faint, the GALFIT fitting results are often highly sensitive to the adopted starting parameters for the model fitting, as well as to the presence of a central component, such as a bulge or NSC.To mitigate the impact of these factors on our results, we adopt a twostage approach that we describe below and repeat the fitting multiple times with different adopted starting parameter values in each of the two stages.In the fitting, we use a convolution box set to half the length of the image, a magnitude zero-point of 22.5, originating from the definition of nanomaggies, and a plate scale of 0.262 arcsec pix −1 for the Legacy Survey.We also summarize the procedure in Table 1.
Our goal for the first fitting stage is to obtain a bestfit single Sérsic model for each candidate UDG that is free from the influence of nearby objects or other stellar components (e.g., NSC or bulge).The results from this fit guide the fitting of the more complex models in the second stage.
We fit this single Sérsic model, utilizing distinct image masks for each galaxy.We create these masks (for an example, see Figure 1) using the Source Extractor Python library (SEP, see Barbary 2016, for details) based on (SExtractor; Bertin & Arnouts 1996).We start by subtracting the spatially varying background generated by SEP.We then identify objects defined as groupings of at least 5 adjacent pixels each with a flux that is 1.5σ above the background.We mask these objects except for the galaxy itself.We remove from further consideration the four UDG candidates whose masked regions cover more than 50% of the entire image area (leaving 6538 for study at this point) because these are poorly constrained in the fitting.Furthermore, in many cases, the fitting results would be strongly affected by an existing central NSC or bulge in the host galaxy.Thus, we augment the mask to include any central region of each candidate that contains at least 5 adjacent pixels that are 1.5 times (0.44 mag) brighter than µ = 24 mag arcsec −2 .This step introduces a potential bias against low luminosity NSCs that are not masked because they do not rise to this level and therefore are incorporated into the base Sérsic model.We address sample completeness further below.As we alluded to previously, to mitigate sensitivity to the adopted initial GALFIT fitting parameter values, we repeat the fitting six times using different values.We use combinations of two different effective radii (30 and 50 pix) and three surface brightness values at r e (25, 28, and 31 mag arcsec −2 ).We calculate the reduced chi-squared statistic, χ 2 ν , within a circular region of radius 50 pixels (∼ 13.1 arcsec) centered on the image center and select the model with the smallest χ 2 ν value.The value of 50 pixels ensures that the bulk of the host galaxies falls within the evaluation region.We find that GALFIT occasionally produces a model fit with acceptable χ 2 ν but unusually large final parameter uncertainties, even though other realizations (i.e.different initial parameters) result in fits with typical uncertainties.Because these solutions offer no meaningful constraints on the fit parameters, we only consider models where r e > 2σ re .We settled on this criterion after exploring a range of options, but these odd cases are clearly distinct from the remainder and various criteria could have been adopted without qualitatively affecting the results.In the rare case where all models for a particular galaxy fail this criterion (10 galaxies), we consider the galaxy as having failed our fitting procedure.Additionally, we have 364 galaxies for which all models are rejected with  at least 90% confidence based on their χ 2 ν .These are also categorized as having failed our fitting procedure.
In the second stage, we unmask the center of the host galaxy and all other sources within 0.5r e and compare new fits using five independent model classes for each candidate.These model classes consist of: 1) a single Sérsic profile, to model UDG candidates with no unresolved nuclear source; 2) a pair of Sérsic profiles, to model UDG candidates with an additional extended central source, such as a separate bulge component; 3) a Sérsic profile plus a PSF profile, to model UDG candidates with an unresolved central source; 4) a pair of Sérsic profiles plus a PSF profile, to model UDG candidates with an extended central source, such as a separate bulge component, and an unresolved central source; and 5) a Sérsic profile plus two PSF profiles, to model UDG candidates with two unresolved central sources.In the second and fourth model classes, we refer to the first Sérsic component, intended to model the galaxy as a whole as S 1 and the second, more compact, component, which is intended to model any nuclear excess, as S 2 .For models with an unresolved central source, we refer to that component as N 1,P SF .In the fifth model class, we refer to the PSF component closer to the center of the Sérsic component as N 1,P SF and the further one as N 2,P SF .Among the NSC targets, we have exclude those with a point source magnitude error in either g or r−band larger than 0.2 mags.
As done in the first fitting stage, we fit each model class six times using different initial parameters.The initial parameters we adopt for S 1 are those from the best-fit model in the first stage, except for the central surface brightness.For the initial brightness of S 1 and S 2 we take permutations of the best fit from the first stage and values that are three magnitudes fainter and brighter.For the size of S 2 we adopt a starting size of five pixels.We found no improvement in fitting when varying r e for S 1 , so we simply use the value obtained in the first stage.
We provide GALFIT with the original UDG candidate image, the mask, the PSF, the σ-image, and a constraint file that sets the search range for each of the free parameters.The free parameters are the following: central positions, Sérsic indices (n), effective radii (r e ), magnitudes (m), axis ratios (AR), and position angles (PA) for S 1 and S 2 ; and the position and amplitude for N 1,P SF and N 2,P SF ; and the background level.The constraint file provides initial parameter values and range limits for the chosen model parameters (Peng et al. 2002).We constrain the centers of the various components to lie within a 40 by 40 pixel square centered on the candidate UDG.In some studies of NSCs, the centrality of the source is a critical criterion (e.g., Côté et al. 2006;Georgiev & Böker 2014), but at this stage of our analysis we allow for significantly offset unresolved components.Because our sample consists of diffuse galaxies, we set the upper limit of the Sérsic index to be 2.0 (our UDG candidates have ⟨n⟩ < 1; Zaritsky et al. 2022).On the other hand, we allow the compact Sérsic component (S 2 ) to have n as large as 5.0 because such a component could be morphologically comparable to a bulge.We also set a lower limit on the axis ratio of the Sérsic components of 0.3 to prevent unrealistically elongated models (SMUDGes sample has an axis ratio threshold of 0.34 < b/a; Zaritsky et al. 2023).Finally, we require that the effective radius of each Sérsic component exceed 0.75 times the size of the PSF to ensure differentiation between what we consider to be a resolved component from an unresolved one.
Before comparing the resulting models, we exclude models that can be rejected with at least 90% confidence given their χ 2 ν .Because we are primarily concerned with modeling the center of the UDG, the area within which we evaluate χ 2 ν is now a circular region of radius 0.5r e , centered on the UDG candidate, where r e and the central position come from the single Sérsic model fit from the first fitting stage.Again, we only consider models where r e > 2σ re .If no models survive our χ 2 ν and r e > 2σ re criteria, then that galaxy is considered to have failed our fitting procedure.
Although χ 2 ν values provide a goodness-of-fit measure, they are inappropriate for selecting among models of differing intrinsic complexity.Because a model with greater fitting freedom will naturally fit the data somewhat better, one must account for this additional flexibility when assessing whether there is statistical evidence in favor of the more complex model.The Akaike Information Criterion (AIC; Akaike 1974) is one formulation that incorporates a penalty for models with higher complexity.
We adopt a slight modification of the original AIC formulation that is referred to as the AICc criterion (Sugiura 1978), where p is the number of model parameters and N is the number of data points that are fit, to compare among models.This modification is appropriate when there are a small number of degrees of freedom and AICc will converge to the original AIC criterion as the degrees of freedom increase.Models where S 2 has a large n value, ≥ 4 and so greater than that corresponding to a de Vaucouleurs profile (de Vaucouleurs 1948), often appear visually indistinguishable from models with N 1,P SF .Furthermore, many of these S 2 components have small r e , resulting in their observed morphology being dominated by the convolution with the PSF.For the vast majority of these, the models with N 1,P SF are statistically favored over models with S 2 , but not at the 2σ confidence level.After visual inspection, we decided to reclassify the subset of systems where two-component models are statistically favored over a one-component model at greater than 2σ confidence and where that S 2 component has n ≥ 4, r e < 10 pix (2.62 arcsec), and lies within 5 pixels of the competing N 1,P SF as S 1 + N 1,P SF (see Figure 3).This increases the total number of SMUDGes within which we identify unresolved sources to 842.Nevertheless, the exact division, if one even exists, between a PSF-convolved unresolved source and a steep inner profile is not a simple issue (Côté et al. 2007).
We present a summary of our classifications in Table 2 and an electronic version of a SMUDGes (Zaritsky et al. 2023) line-matched catalog of our classifications as well as an example in Table 3.For each galaxy, we assign the most specific and complex model that is statistically preferred by the AICc criterion at the > 2σ confidence level, except for the S 1 classification which includes all cases where S 1 is favored, even if not at > 2σ confidence, and cases where another more complex model was statistically favored but not at beyond the 2σ confidence level.To provide more detail, in the case where we assign the S 1 +S 2 +N 1,P SF classification it must be statistically preferred over both the S 1 +S 2 model and the S 1 +N 1,P SF model at greater than the 2σ confidence level.In Table 2 we divide the classifications into objects with and without unresolved sources.The latter includes galaxies with indications of possible unresolved sources such as those where models with N 1,P SF are preferred but not at the 2σ level and those for which S 1 +S 2 + N 1,P SF is preferred at the 2σ level.We included the latter in this category because following visual inspection we decided that these are generally highly complex systems that are simply difficult to model and do not necessarily show evidence for NSCs.As a group these are identified as S 1 + N 1,P SF ? in the Table.

False Positives
Among the 842 candidate UDGs for which the S 1 + N 1,P SF or S 1 + N 1,P SF + N 2,P SF models are preferred with ≥ 2σ confidence many show large offsets between S 1 and N 1,P SF (for example, 44 have offsets larger than 20 pixels or 5.2 arcsec).This result raises the question of what constitutes a nuclear star cluster.
Our candidates span a range of projected offsets from the center of the host galaxy (Figure 4).The offset distribution, plotted in normalized, distance independent terms of r p /r e , where r p is the projected separation between S 1 and N 1,P SF and r e is the effective radius of the host galaxy (as measured in our initial fitting pass), appears to have two components: a concentrated central one and a far more extended one.The distribution is somewhat difficult to interpret because the positional offsets we allowed for GALFIT translate to different cutoffs in r p /r e for each galaxy.Nevertheless, it suggests that our sample consists of populations we might call 'true NSCs' and 'contamination'.The contaminating population may be a combination of sources physically associated with the galaxy, such as non-nuclear star clusters and star-forming regions, and physically unassoci- Note-Classifications for the SMUDGes candidates.The second column describes the classification outcomes for the SMUDGes that satisfy 0 < g − r < 0.8 and re < 26 arcsec based on the original catalog (Zaritsky et al. 2023) and are not heavily masked.The third column presents the numbers of those with identified unresolved sources that lie within a normalized projected separation from the S1 component that is < 0.10, where the normalization is done using re as measured from the initial fitting of the single Sérsic model.The fourth column refers back to the second column and identifies the number of galaxies with estimated distances in the original catalog.The fifth column refers back to the objects in the fourth column and identifies the number of those that satisfy the projected separation criterion.
ated sources, such as foreground stars and background galaxies.Critically, however, the details of how these populations are differentiated will impact certain questions related to the possibility of off-center nuclear star clusters, whether that be because they have stalled in their inward migration or been jostled off-center by a dynamical event.
To better understand the NSC and background populations, so that we can optimize our selection of NSCs, we explore a set of models for the radial distribution of all candidate NSCs in Figure 4.In those models, we describe the NSC population alternatively as an empirically-motivated projected 1-D Gaussian distribu-tion in normalized projected radial offsets, r p /r e , a 2-D Gaussian distribution in r p /r e on the sky, or a 2-D exponential in r p /r e .The contaminating population we describe alternatively as a Sérsic distribution with the parameters of the host galaxy, assuming the contaminants come primarily from the galaxy itself, or a uniform distribution on the sky, assuming the contaminants are principally either foreground or background sources.Finally, we also account for our radial completeness by evaluating the fraction of the sample at each radius for which our selection criteria would have allowed us to find an NSC.The radial completeness is affected both by our criteria that we only search for unresolved sources within 0.5r e and allow GALFIT to explore positional offsets within a box of 40 by 40 pixels.
We evaluate the parameters for each of these model combinations using a Bayesian approach and the Markov Chain Monte Carlo (MCMC) ensemble sampler called EMCEE (Foreman-Mackey et al. 2013a).We model the distribution using different combinations of either a 1-D Gaussian or exponential primary distribution and a Sérsic profile or uniform background secondary distribution.We adopt uniform priors for the amplitude and standard deviation of the central population and for the amplitude of the second component.We find a slight preference for the models that describe the distribution as a 1-D Gaussian central component plus Sérsic-distributed contamination, but the exponentialdistributed central component plus uniform background is nearly as good a fit (Figure 4).As such, we conclude that we cannot distinguish whether the contaminating population is primarily within the host galaxy or unassociated on the basis of this data and fitting.This question will be reexamined in future work where we explore the extended population in greater detail.The contaminating population is potentially physically interesting because it might include clusters that are otherwise similar to NSCs but found at large radii.
To distinguish true NSCs from contamination, the choice of model becomes irrelevant because using either model we conclude that the contamination in our recovered NSC sample is 15% (the percentage that we have set as our target) when we reject candidates with r p /r e > 0.10.Setting 0.10 as an upper limit on the normalized radial offset, we retain 325 UDG candidates with NSCs, and reject 517.Examples of our final NSC sample are presented in Figure 5 and this is the sample we present as NSC hosting SMUDGes.This severe selection demonstrates directly the impact of any imposed radial selection on the overall population.Consider that a slightly more permissive criterion of r p /r e < 0.2 results in ∼50% more candidates (463), although a larger .The distribution of normalized radial offsets between S1 and N1,P SF , rp/re.Solid blue lines represent the data, the blue dotted lines show the models for the central Gaussian and background that together, modulated by the radial completeness, combine to produce the red dashed line.In the left panel, the background is assumed to follow the S1 profile, while in the right panel, it is assumed to be randomly distributed on the sky.The models fit the data indistinguishably well, demonstrating that we cannot differentiate between the two background scenarios at the current time and that our conclusions are insensitive to this choice.
fraction of these will come from the contaminating population.
To highlight again how comparing among studies is fraught, we note that Poulain et al. (2021) accept the brightest unresolved source out to a radial offset of 0.5r e as an NSC.After an initial reading of our results one would conclude that their sample is dominated by contaminants.However, their use of superior imaging data (both deeper and of higher angular resolution) might allow them to reject contamination far better than we can in our data.Without carefully examining both data sets and redoing parts of the analysis similarly, it is not possible to reach definitive conclusions regarding a comparison between our two studies.
Our strict radial selection cut ensures relatively high purity (defined to be 85%) but excludes true NSCs that lie beyond this radial cut.We use our best-fit model for this component (a 1-D Gaussian with σ = 0.10) to calculate that 35% of this component lies outside of r p /r e = 0.10.From the 325 candidate UDGs with NSCs at projected offsets less than r p /r e = 0.10, we calculate that 276 are true NSCs after correcting for the 15% contamination and that this corresponds to a total population of UDG candidates with NSCs of 340 after correcting for systems that lie at r p /r e > 0.10.

False Negatives
We classify as single component 2979 candidate UDGs.Among this set, there may be some for which the NSCs fall below our detection limit and these are, as such, false negatives.To estimate our limiting mag-nitude, we randomly select 100 galaxies that are preferred by a single Sérsic model and insert artificial point sources of varying brightness ranging from 18 to 30 magnitudes, in intervals of 0.04 magnitudes, at the centers of the images.We model these point sources using a Gaussian of 1.5 pixel width that extends out to 20 pixels and run our pipeline.When we do not recover the point source, there are two failure modes: (1) the bestfit model is not the Sérsic+PSF model; (2) the best-fit model is the Sérsic+PSF model but has a confidence level of less than 2σ (∆AICc <11.83).For each of the 100 galaxies, we set the magnitude of the brightest inserted point source that we fail to recover as the detection limit for that galaxy.In Figure 6 we present the set of these detection limits as a function of the host galaxy central surface brightness.The detection limits lie mostly between magnitudes of 23 to 25.We find that the limits correlate weakly with central surface brightness (confidence level 99.2% and Spearman Rank correlation coefficient of 0.26).This result matches our intuition that it should be more difficult to detect an NSC in a galaxy that is itself intrinsically brighter in its center.However, given the large scatter in detection limits about this mean trend (Figure 6), there must be other factors at play as well and we neglect the subtle, but real, surface brightness dependence in our subsequent qualitative discussion of completeness.
We now examine the consequences of using shallower, lower resolution data by comparing our classification results to those of Lim et al. (2018) for Coma galaxies using Hubble Space Telescope images.Lim et al. (2018) presented 44 UDGs in their Table 1, 26 of which match SMUDGes.The majority of the unmatched galaxies are either of smaller angular extent than the SMUDGes limit or at a surface brightness where SMUDGes is incomplete.From the matched 26, Lim et al. (2018) concluded from visual inspection that five contain NSCs.Among these, we find an NSC in only one of those, but the normalized separation (r p /r e = 0.226) is greater than our separation criterion (r p /r e < 0.10).We attribute our failure to identify NSCs in the other four as the result of the difference between our magnitude limit of ∼24 mag and theirs of ∼ 27.4 mag (Lim et al. 2018).Lim et al. (2018) do not present photometry for their NSCs, so we cannot confirm that these are indeed fainter than our detection limit, but as we show in Figure 7, our detection limit is likely to exclude nearly all NSCs by the time we are considering galaxies at the distance of the Coma cluster.This incompleteness also helps to explain why our occupation fraction (the fraction of candidates that we identify to host NSCs), which is globally ∼ 0.05 (340/6538) for our sample, is so much lower than the  ≳ 0.20 that is commonly found (e.g., den Brok et al. 2014;Lim et al. 2018;Eigenthaler et al. 2018;Carlsten et al. 2022), although it may also reflect that our sample is not as constrained to high density environments and the occupation fraction is measured to be higher in high density environments (Lim et al. 2018;Sánchez-Janssen et al. 2019;Poulain et al. 2021;Carlsten et al. 2022).

RESULTS
A principal product of this study is an NSC-related classification for each galaxy in the SMUDGes catalog and a measurement of the NSC properties when there is one.We present in Table 3 the first five lines of the full catalog, available electronically, where the objects are row matched with the Zaritsky et al. (2023) catalog.We provide a summary of the classifications in Table 2 Note-We present here the first five lines of the full catalog, which is available in electronic form.Magnitudes refer to the S1 and N1,P SF components of the best fit model, while colors and effective radii refer to those derived from the initial single Sérsic fit, which are found to be more stable.Entries of −99.00 signal invalid values corresponding to SMUDGes galaxies for which an unresolved source is not identified.We retain these SMUDGes to maintain line-by-line matching with the Zaritsky et al. ( 2023) catalog.The numerical values in the Class column correspond, in order, to the classification categories in Table 2 ( S1 ≡ 1, S1+S2 ≡ 2, S1+N1,P SF ?≡ 3, S1+N1,P SF ≡ 4, S1+S2+N1,P SF ≡ 5, failed fitting≡ 6) and galaxies that were excluded prior to fitting (≡ 7).Note that some S1+N1,P SF and S1+S2+N1,P SF targets have photometric errors larger than 0.2 mags.These were classified as S1 ≡ 1 in Table 2 and subsequently excluded from our analysis.
by listing the number of objects in each of our classification classes, the number of which have an inferred unresolved source that meets our r p /r e < 0.10 criterion, the number in each class with distance estimates from the original SMUDGes catalog (Zaritsky et al. 2023), and the number of those with distance estimates that have an inferred unresolved source component that meets our r p /r e criterion.

NSC Positional Offsets
We have defined NSCs as the centrally located subpopulation of unresolved sources coincident with our UDG candidates.The degree to which NSCs are truly found at the dynamical center of their host galaxy is somewhat difficult to address because NSCs are often required, by definition, to be at the galaxy's center (e.g., Böker et al. 2002;Côté et al. 2006;Neumayer et al. 2011).For specific examples, we cite Côté et al. (2007) which sets an offset upper limit of 0.02r e (about 20 pc for the typical distance in their sample) and Neumayer et al. (2020) who propose an offset limit of 50 pc.In contrast, Poulain et al. (2021), who allow for larger offsets, find NSCs out to 0.58r e .Binggeli et al. (2000) identify some with even larger offsets.Offsets are potentially interesting to measure because they may be caused by dynamical interactions (Bellovary et al. 2018) and could help us measure the shape of the gravitational potential (Miller & Smith 1992;Taga & Iye 1998).
Although we too have imposed such a requirement, based on the distribution of r p /r e , we do measure differences in the degree of alignment between S 1 and N 1,P SF .Careful examination of Figure 4 shows that the observed r p /r e distribution does not peak at zero separation and we have measured a dispersion of 0.10.To assess whether these findings reflect physical scatter in the co-centricity of the NSC and host or are simply the results of measurement errors always leading to positive offsets, we compare the measured distribution to the expected scatter arising simply from our observational uncertainties.When we adopt the GALFIT positional uncertainties in S 1 and N 1,P SF for the individual observed systems and use those to randomly generate S 1 and N 1,P SF pairs, we find a distribution of r p /r e that also does not peak at zero but which has a dispersion of only 0.02.From this result, we conclude that the observed scatter does constitute evidence of physical offsets.Our measured dispersion is in excellent agreement with the median offset measured in the MATLAS dwarf galaxy sample of 0.10r e (Poulain et al. 2021).Because the Poulain et al. (2021) sample has galaxies with a somewhat brighter central surface brightness than SMUDGes1 , this agreement may be in conflict with findings of increasing offsets in lower surface brightness hosts (Binggeli et al. 2000;Barazza et al. 2003), although, as noted by Poulain et al. (2021) this measurement is complicated by the likely larger uncertainties in the measurements of the centers of lower surface brightness galaxies.
A different complication is that it is not clear which component is the better tracer of the dynamical center.It is possible that the photometric center of S 1 is precisely measured but that it does not reflect the dynamical center.All we can conclude is that the photometric center of S 1 and the position of N 1,P SF show greater scatter than accounted for by the observational centroiding uncertainties.In fact, the entire dark matter halo may be offset leading to strong lopsidedness in the center of the galaxy (Prasad & Jog 2017).

Nuclear and Non-nuclear Stellar Clusters: Host Galaxy Colors
A related question is whether the unresolved sources associated with the more extended population include a substantial number of stellar clusters that are otherwise indistinguishable from those in the NSC population.Such a population would have important repercussions on models of NSC formation.
To work with a set of unresolved sources that are mostly independent of what we have identified as the NSC population, we define a non-nuclear class in an analogous way as we did for the nuclear class, again focusing on purity.We use the model of the r p /r e distribution but this time search for a lower limit on r p /r e that ensures ≤ 15% contamination of the non-nuclear population by the nuclear population.By setting that lower limit to be 0.38, we find 134 UDG candidates with unresolved sources that we define to be a 'clean' sample of non-nuclear unresolved sources.We find that the distributions in color differ markedly (Figure 8).NSC hosts have g − r ∼ 0.6, with a modest dispersion ∼ 0.05, indicating that these hosts are predominantly on the red sequence.The number of NSC hosts that are much bluer (g − r < 0.5) is only 17 (out of the total of 325) and is smaller than the anticipated level of contamination (i.e.15% or 30 galaxies), suggesting that NSCs, as we have defined them, might only be found in red sequence hosts galaxies in our sample.In contrast, the hosts of the non-nuclear unresolved-source population are distributed broadly in color.This flatter color distribution, and its similarity to the color distribution of the entire SMUDGes sample, further suggests that a significant fraction of this population is indeed contamination and not physically associated with the hosts.
The relative deficit of NSCs in our blue hosts raises the question of whether NSCs are exclusively a red galaxy phenomenon or reflect a selection bias that is working against us.Previous studies of late-type galaxies have found large occupation fractions (Georgiev & Böker 2014), albeit not in hosts of this low surface brightness.We identify various challenges in identifying NSCs in blue SMUDGes.The hosts are more likely to be irregular and difficult to model, leading both to failures in the fitting and to greater statistical noise in the centroiding.The latter can both complicate our procedure for assessing the significance of an additional unresolved source and corrupt our measurement of the offset.Nevertheless, SMUDGes with color 0.4 < g − r < 0.5 are not highly irregular and yet these already host relatively few NSCs, as can be seen by comparing the various distributions presented in Figure 8.

NSC properties
We now focus on the internal properties of NSCs in SMUDGes.We have two distinct populations to discuss.First, we draw inferences from our sample of 325 NSCs for which we do not have distance estimates, to discuss distance-independent aspects of the population.Their hosts are all galaxies of low central surface brightness, although it is likely that many will not satisfy the UDG physical size criterion.Second, we select only those hosts for which we do have distance estimates (see Zaritsky et al. 2022Zaritsky et al. , 2023, for a discussion of distance estimation), and then either focus on the physical properties of the 144 such systems or select only the 33 of those that meet the UDG size criterion (r e > 1.5 kpc).
The luminosity (or corresponding stellar mass for similar mass-to-light ratios) of an NSC scales with that of the host galaxy (Neumayer et al. 2020).We find this broadly holds for our sample (Figure 9).We adopt a stellar M/L of 1.8 M ⊙ /L ⊙ for our transformation to stellar masses.The masses we derive, typically between 10 5 and 10 7 M ⊙ are within the previously measured range for NSCs (Neumayer et al. 2020), although on the lower end as expected given that the stellar masses of our host sample is also lower than average.The data are mostly consistent with the previously published relationship (Neumayer et al. 2020), although a steeper proportionality relation appears to fit the data better.However, recall that we are incomplete at lower NSC masses for more distant, hence typically larger and more massive hosts, which could help fill in the distribution in the lower right portion of the diagram.Furthermore, distance errors, and we expect that a significant fraction of the estimated redshifts (∼ 30 %) are incorrect (Zaritsky et al. 2023), will tend to have a preference for scattering objects to larger distances and proportionally higher masses along both axes.We conclude that our findings are consistent with the previous relation obtained with higher surface brightness galaxies.As such, we conclude that we find no systemic difference in the relationship between NSC and host galaxy stellar masses for low surface brightness galaxies relative to what was previously found for the more general NSC host population.
We find that the stellar mass fraction in NSCs can reach close to 0.1 in the most extreme objects, and is typically about 1/50th (Figure 9).This measurement is distance independent.At the upper end of the mass fractions, our findings are consistent with those of Binggeli et al. (2000) for Virgo dwarfs and mostly consistent with those of Poulain et al. (2021) for the MATLAS dwarfs.The latter do find objects with mass fractions >0.1, with one reaching 0.45, although they find only 6 in their sample of 508 with mass fraction > 0.2.For the lower mass fractions our sample is incomplete, so we do not discuss that side of the distribution.Our typical value of 0.02 is in excellent agreement with their median value of 0.017, again suggesting no strong dependence on host surface brightness.
The objects with the largest mass ratios appear to pose interesting constraints on formation models as it is difficult to envision how nearly 10% of a galaxy's stellar mass could end up as an NSC either in the globular cluster infall model or in the central star formation model.In Figure 10 we present the five galaxies and the five UDGs with the most extreme NSC mass fractions in our sample.In our entire sample, where we can use galaxies without estimated distances because we are considering only the ratio of stellar masses, we do find that the most extreme galaxies, such as those shown in the upper row of the Figure, have slightly more than 10% of their stel-lar mass in their NSCs.For the most extreme UDGs that percentage drops to less than 3%.Given the small number of confirmed UDGs in the sample, we do not yet know if this represents a real difference between UDGs and other galaxies or if we have simply not yet sampled the high fraction tail sufficiently well.In closing, we note that the NSC in one of the UDGs (SMDG0430339-052909) has a measured color that is anomalously red and so we suspect that this is an unresolved background galaxy rather than an NSC.Having one contaminant among the ten galaxies shown in the Figure is consistent with our anticipated 15% contamination rate.

Environment
The hosts of NSCs in the SMUDGes sample are almost exclusively on the red sequence (Figure 8).Given the connection between UDG color and environment (Prole et al. 2018;Kadowaki et al. 2021), this finding may point to an environmental dependence of the occupation fraction, in the sense as others have found (Lim et al. 2018;Sánchez-Janssen et al. 2019;Poulain et al. 2021;Carlsten et al. 2022).In Figure 11 we compare the distribution of UDGs with NSCs to those without NSCs in the environmental parameter space defined by the estimated velocity dispersion, σ EN V , and the richness, N EN V , of the hosting overdensity that are provided by Zaritsky et al. (2023).Those parameters are known to be highly uncertain measures of environment, so Zaritsky et al. (2023) suggest using them in combination to define poor and rich environments (the lower left and upper right quadrants in Figure 11, respectively).We are limited in this comparison by the small number of UDGs with NSCs (33), but the visual impression from the Figure is that UDGs with NSCs tend to higher values of σ EN V than the overall sample.Although the sense of this behavior is as previously observed (e.g., Lim et al. 2018;Sánchez-Janssen et al. 2019;Poulain et al. 2021), statistical tests (comparison of means, medians, KS test) do not yield a statistically significant detection of a difference in the two populations presented here.
Confirming this preliminary result is important because it would extend the relation between NSC occupation fraction and environment to that of groups and the field and to UDGs.This will require spectroscopic redshifts of many more SMUDGes galaxies to provide distances and enable us to convert angular measurements to physical ones.

SUMMARY
We present the results of our photometric search for potential nuclear star clusters (NSCs) hosted by the ultra-diffuse galaxy candidates in the SMUDGes catalog (Zaritsky et al. 2023).Using r-band images from DR9 of the Legacy Survey (Dey et al. 2019), we develop an algorithm to statistically determine if additional photometric components beyond a single Sérsic model within 0.5r e are needed, and then if among those components there is one that is best modeled as an unresolved object.We find that slightly over half of the SMUDGes sample does show evidence for additional components.Among those, we identify 1059 for which we find with 90% confidence the need for an unresolved source.
The distribution in projected radius of these unresolved sources shows a peak near zero separation and a second more extended component.We explore models and quantify the nature of the two components, attributing only the central component to NSCs.We use our models of the radial distribution to define a maximum projected separation for our defined NSC sample that ensures an 85% pure sample of candidate NSCs (0.10r e ).
We explore our incompleteness using simulations and establish that we are significantly incomplete due to our relatively shallow imaging (NSC magnitude limits ranging from ∼ 23 to 25 mag depending on the specifics of the host galaxy and surroundings).Nevertheless, we are able to confirm with confidence 325 SMUDGes galaxies with NSCs, 144 of which also have estimated distances provided by Zaritsky et al. (2023).Among those 144, we confirm 33 as UDGs with NSCs.
Despite our identification of NSCs as a population that is closely aligned with the center of their host galaxy, we find that the observed scatter of positional offsets between NSCs and their hosts is greater than expected from measurement errors alone.Our estimate of the dispersion in offsets (0.10r e ) is in good agreement with the median offset measured from MATLAS dwarf galaxies (Poulain et al. 2021).Such offsets could be used to constrain formation scenarios and models of the host's gravitational potential (Miller & Smith 1992;Taga & Iye 1998;Bellovary et al. 2018).
We find that our sample of NSCs is hosted almost exclusively by galaxies on the red sequence.The number of NSCs found in bluer hosts is consistent with the expected level of contamination.This result may reflect the color-environment trend identified for such galaxies (Prole et al. 2018;Kadowaki et al. 2021) and the greater NSC occupation fraction in denser environments (Lim et al. 2018;Sánchez-Janssen et al. 2019;Poulain et al. 2021;Carlsten et al. 2022).Unresolved sources away from the nucleus are found in hosts whose color distribution matches that of the SMUDGes sources overall, suggesting no physical connection.We discuss some possible selection biases against our finding NSCs in blue hosts.
Despite our focus on low surface brightness galaxies, and UDGs in particular, we find that the NSCs in our sample fall on the NSC-host galaxy stellar mass relationship found previously (Neumayer et al. 2020) from higher surface brightness objects.There is potentially a deviation for the largest objects in our sample (the UDGs) but our strong selection effects and the currently necessity of comparing across studies rather than within a single study limits our conclusions in this regard.The typical NSC in our sample contains about 0.02 of the total stellar mass of the host galaxy, although the most extreme objects reach a fraction of nearly 0.1.These results are in agreement with previous studies of NSCs in somewhat different galaxy samples (Binggeli et al. 2000;Poulain et al. 2021).
Finally, we search for possible environmental effects in the NSC population.Despite this being principally a field sample, we find a suggestion that NSCs are more likely in UDGs in higher density environments, in agreement with previous results (Lim et al. 2018;Sánchez-Janssen et al. 2019;Poulain et al. 2021;Carlsten et al. 2022).However, quantitative analysis of this trend in our data does not yet yield statistically significant results.Increasing the sample of UDGs with NSCs, by having spectroscopic redshifts of a larger fraction of the SMUDGes NSC sample, will allow us to further assess this possibility in the future.
The SMUDGes catalog and the NSC extension catalog provided here enable ongoing work on the nature of NSCs and their hosts, extending the latter to low surface brightness, physically large galaxies.This sample will benefit greatly from complementary high angular resolution imaging to come from surveys carried out with the Euclid and Nancy Grace Roman telescopes, spectroscopy from highly multiplexed surveys such as DESI, and from continued dedicated follow-up observations.

Figure 1 .
Figure 1.An example UDG candidate with a plausible NSC (SMDG0004118+163159, the second galaxy with what we will classify as an NSC in our right ascension ordered catalog) with several nearby contaminating sources is presented in the left panel.The mask in the right panel highlights most of the visible nearby sources with the central object shaded in gray indicating that it is unmasked for the final model fitting, but would masked in the initial model fitting (see text for details).
Note-Columns indicate morphological componentsused in the six different model combinations.The initial parameters used for of the second stage are those of the best-fit model from the first stage.We adopt two-or three-component models only when they are significantly favored over the one-or twocomponent models, respectively.Similarly, we only adopt models with unresolved sources when they are significantly favored over those with resolved components.A more comprehensive description can be found in Section 2.3.

Figure 2 .
Figure 2. Model and residual images for SMDG0004118+163159.Left to right: the original image, single Sérsic model, residual image using the single Sérsic model, Sérsic + PSF model, and residual image using the Sérsic + PSF model.The need for a central unresolved source is evident when comparing the residual images.

Figure 3 .
Figure 3. Distribution of n and re for S2 components when two components are statistically favored over a single component, but no statistically significant difference exists between competing multi-component models.The box highlights the parameter region where, after visual examination, we reclassify systems as having unresolved sources.

Figure 4
Figure4.The distribution of normalized radial offsets between S1 and N1,P SF , rp/re.Solid blue lines represent the data, the blue dotted lines show the models for the central Gaussian and background that together, modulated by the radial completeness, combine to produce the red dashed line.In the left panel, the background is assumed to follow the S1 profile, while in the right panel, it is assumed to be randomly distributed on the sky.The models fit the data indistinguishably well, demonstrating that we cannot differentiate between the two background scenarios at the current time and that our conclusions are insensitive to this choice.

Figure 5 .
Figure 5. Mosaic of the first 21 SMUDGes sources in Right Ascension that have an identified NSC.Images span slightly over an arcmin on a side with North at the top and East to the left.These are drawn from the Legacy Surveys on-line viewer (https://www.legacysurvey.org/viewer).Object labels are included at the top of each panel.

Figure 6 .
Figure 6.NSC detection limits in the r band determined for 100 Single Sérsic SMUDGes images as a function of r−band central surface brightness, µ0,r.Points represent the brightest simulated point source that was not recovered as an NSC by our procedure.Horizontal lines mark where we reach the corresponding incompleteness percentage.The three labeled values correspond to limiting magnitudes of rP SF of 23.2, 24.1, and 24.9 mag.

Figure 7 .
Figure 7. Absolute magnitude of NSCs, for SMUDGes with estimated distances, across redshift.Systems that are UDGs (re ≥ 1.5 kpc; filled blue circles) and non-UDGs (re < 1.5 kpc; open red circles) are plotted.The curves designate completeness limits corresponding to the three incompleteness fractions shown in Figure 6.We also include a vertical line that indicates the redshift beyond which any object in the SMUDGes catalog would satisfy the UDGs size criterion, and the redshifts corresponding to the Virgo and Coma clusters.

Figure 8 .
Figure 8.The distribution of host galaxy colors, (g −r)UDG from Table 3, and apparent magnitudes for systems with unresolved sources with rp/re < 0.10 (NSCs; filled circles) and rp/re > 0.38 (non-NSCs; open circles).The right panel shows the normalized distributions in color of these two populations (filled histogram represents NSCs; unfilled blue, solid line represents non-NSCs), as well as that for the full SMUDGes sample (represented by the dashed red line).

Figure 9 .
Figure 9. NSC stellar mass vs. host galaxy stellar mass.Stellar masses calculated assuming a fixed stellar mass-tolight ratio of 1.8 M⊙/L⊙.Solid line and shaded region represent the relationship and its uncertainty presented by Neumayer et al. (2020).The heavy dotted line represents a normalized 1:1 relation where MNSC = 0.019MHOST .The thinner dotted line corresponds to MNSC = 0.1MHOST and represents the upper limit on the NSC mass fraction.Red stars represent UDGs.Blue dots represent the galaxies that do not meet the UDG physical size criterion.

Figure 10 .
Figure 10.Mosaic of galaxies with extreme NSC mass fractions.The top row presents the five galaxies in our sample with the largest NSC mass fraction.The bottom row presents the five UDGs (re ≥ 1.5 kpc) with the largest NSC mass fraction.Images span slightly over an arcmin on a side with North at the top and East to the left.These are drawn from the Legacy Surveys on-line viewer.Object labels are included at the top of each panel.

Figure 11 .
Figure 11.Environmental properties of UDGs with and without NSCs.We compare the local measures of environment, σENV and NENV , measures of the local velocity dispersion and richness provided by Zaritsky et al. (2023), for UDGs with NSCs (red stars) and UDGs without NSCs (smoothed blue distribution).The dotted lines mark the regions set by Zaritsky et al. (2023) to differentiate low (NENV < 15 and σENV < 500) and high (NENV > 15 and σENV > 500) density environments.

Table 1 .
Model Summary to calculate the confidence level corresponding to any specific value of ∆AICc.For our situation, a 2σ confidence level corresponds to ∆AICc = 11.83.We adopt this threshold to assess if the best-fitting twocomponent models are significantly preferred over the best-fitting one-component model.If it is, then we use this same threshold to assess if the model with S 1 + N 1,P SF or S 1 + N 1,P SF + N 2,P SF is significantly preferred over that with S 1 + S 2 .Figure2shows an example of a model and residual images created by GAL-FIT for the single Sérsic model and the S 1 + N 1,P SF model for SMDG0004118+163159.At this stage, we identify 757 SMUDGes for which the S 1 + N 1,P SF or S 1 + N 1,P SF + N 2,P SF model is preferred over other competing models with greater than 2σ confidence.
The model with the smaller AICc value is the statistically preferred model although the greater the difference (∆AICc) the greater the confidence with which one can discriminate among them.Because AICc values describe likelihoods and are distributed like χ 2 , we are able

Table 2 .
Classification Summary