Constraining the stellar populations of ultra-diffuse galaxies in the MATLAS survey using spectral energy distribution fitting

We use spectral energy distribution (SED) fitting to place constraints on the stellar populations of 59 ultra-diffuse galaxies (UDGs) in the low-to-moderate density fields of the MATLAS survey. We use the routine PROSPECTOR , coupled with archival data in the optical from DECaLS, and near-and mid-infrared imaging from WISE , to recover the stellar masses, ages, metallicities and star formation timescales of the UDGs. We find that a subsample of the UDGs lies within the scatter of the mass–metallicity relation (MZR) for local classical dwarfs. However, another subsample is more metal-poor, being consistent with the evolving MZR at high-redshift. We investigate UDG positioning trends in the mass–metallicity plane as a function of surface brightness, effective radius, axis ratio, local volume density, mass-weighted age, star formation timescale, globular cluster (GC) counts and GC specific frequency. We find that our sample of UDGs can be separated into two main classes. Class A: Comprised of UDGs with lower stellar masses, prolonged star formation histories (SFHs), more elongated, inhabiting less dense environments, hosting fewer GCs, younger, consistent with the classical dwarf MZR, and fainter. Class B: UDGs with higher stellar masses, rapid SFHs, rounder, inhabiting the densest of our probed environments, hosting on average the most numerous GC systems, older, consistent with the high-redshift MZR (i.e., consistent with early-quenching), and brighter. The combination of these properties suggests that UDGs of Class A are consistent with a ‘puffed-up dwarf’ formation scenario, while UDGs of Class B seem to be better explained by ‘failed galaxy’ scenarios.

While resembling classical dwarf galaxies in terms of their stellar masses ( ★ < 10 9  ⊙ ) and surface brightnesses ( ,0 > 24 mag arcsec −2 ), they by definition exhibit much larger effective radii ( e > 1.5 kpc).Some of them were found to have populous GC systems and evidence for more massive dark matter halos ( halo > 10 11  ⊙ ) than their classical dwarf counterparts (Beasley et al. 2016;van Dokkum et al. 2019;Forbes et al. 2020;Gannon et al. 2020Gannon et al. , 2022Gannon et al. , 2023;;Toloba et al. 2023;Zaritsky et al. 2023).This mixture of dwarf-and massive-like galaxy properties has led to many different proposed formation scenarios for these galaxies.
A scenario relying on internal processes suggests that UDGs start their lives as classical dwarf galaxies, and through an extended sequence (i.e., repeated episodes over a long period of time) of supernovae feedback (Di Cintio et al. 2017;Chan et al. 2018), they get enlarged to the sizes we observe today.Another scenario suggests that dwarfs with high-spin halos can evolve into UDGs (Amorisco & Loeb 2016;Benavides et al. 2021).We refer to such scenarios as 'puffed-up dwarf' formation scenarios.Other scenarios rely on external processes to form UDGs via e.g., tidal stripping and heating (e.g., Carleton et al. 2019), tidal interaction (e.g., Tremmel et al. 2020), which account for the apparent radial size excess.Some UDGs in the past have been suggested to be consistent with a tidal dwarf galaxy (TDG) origin as well (see e.g., Duc et al. 2014;Buzzo et al. 2023).We refer to this formation scenario throughout the paper as 'tidal UDGs'.
In addition to these, one scenario suggests that UDGs started their lifetimes destined to be large and massive but had their star formation truncated early on.As a consequence, their stellar masses do not increase at the same rate as their halos, resulting in galaxies with dwarf-like stellar masses enclosed by halo masses similar to those of more massive M33-like galaxies (van Dokkum et al. 2015;Peng & Lim 2016;Danieli et al. 2022;Gannon et al. 2023).For the remainder of the manuscript, we refer to this formation scenario as 'failed galaxy' scenario.
One crucial discriminant between these formation scenarios is the stellar populations of UDGs.For example, if they were formed by a puffed-up dwarf scenario, they would overall be expected to have similar stellar populations to those of classical dwarfs.One expectation, for example, is that they may follow the same scaling relations as classical dwarfs, e.g., the mass-metallicity relation (MZR).It is important to mention, however, that different puffed-up dwarf scenarios may lead to very different star formation histories (SFH), i.e., supernovae feedback would lead to a bursty SFH, while high-spin halos or tidal interactions may be more similar to an exponentially decaying SFH, but in both cases the UDGs are expected to preserve the stellar populations of their progenitors.Tidal UDGs are expected to be young, metal-rich, and gas-dominated for their stellar mass (Haslbauer et al. 2019;Duc et al. 2014).On the other hand, if UDGs are failed galaxies, the predictions for their stellar populations are less clear.They are not expected to have stellar populations similar to classical dwarfs (as classical dwarfs were not their progenitors), nor to more massive M33-like galaxies (as they did not evolve to become such).The only stellar population prediction for galaxies formed by this scenario is that they would have old or even ancient stellar populations (≳ 10 Gyr), as they by definition have suffered from early-quenching.For the same reason, they may be expected to have lower metallicities and possibly be alpha enhanced as they have had shorter timescales to form stellar masses comparable to those of classical dwarfs that had prolonged star formation histories (see e.g., Forbes et al. 2020), and also because, as suggested by Danieli et al. (2022), these galaxies could be mainly formed of disrupted GCs, thus being more metal-poor.As a key point in the comparison, the MZR is known to evolve with redshift (Ma et al. 2016;Chartab et al. 2024).Galaxies that follow the MZR at high-redshift have stellar populations consistent with the chemical enrichment up to that point, thus, either galaxies at high-redshift themselves or galaxies that have not had much chemical evolution since then (i.e., early quenched).We thus may expect failed galaxies to follow the MZRs of high-redshift galaxies.
Some of these predictions were recently confirmed by Ferré-Mateu et al. (2023) (hereafter, FM23), who used the largest spectroscopic study of UDGs to date (25) to show that there is a correlation between -enrichment and the star formation history (SFH) of UDGs, further emphasised by the positioning of UDGs in the cluster phase-space diagram.They found that early-infall UDGs are the ones with the fastest star formation episodes and have higher [/Fe], consistent with early quenching scenarios, such as failed galaxy ones.Conversely, UDGs with prolonged star formation histories are the late infallers, having lower [/Fe] ratios, being better explained by puffed-up dwarf scenarios.
Additionally, Buzzo et al. (2022) (hereafter, B22) used SED fitting to show that there is a correlation between the age, metallicity and environment that UDGs reside in.The older ones display the lowest metallicities, are consistent with the MZR at high-redshift ( ∼ 2) and are found in the densest environments.The younger ones (≤ 8 Gyr) follow the classical dwarf MZR, being more metal-rich and are found in less dense environments, such as the field and groups.These findings are in general agreement with the findings of Barbosa et al. (2020), who used SED fitting to study 100 UDGs in the field and to show that they are mostly young (∼ 6 Gyr) and follow the classical dwarf MZR.In their sample, some UDGs were found to be more metal-rich than the classical dwarf MZR ([Fe/H] ≥ −0.5 dex) and extremely young (≤ 1 Gyr), similar to what has been found for some cluster UDGs (Ferré-Mateu et al. 2018;Ruiz-Lara et al. 2018).Their ages and metallicities are consistent with them being tidal UDGs (Collins & Read 2022).
B22 have also shown that the stellar populations of UDGs seem to correlate with their GC richness.Old and metal-poor UDGs have on average the highest number of GCs, while younger and more metalrich ones (consistent with the classical dwarf MZR) have the lowest GC numbers.This trend, however, was not found by FM23 using spectroscopy.Selection effects can be the cause of the differences, since both samples are not fully representative of the population of UDGs.Thus, to understand if this trend holds and whether it can or cannot be extrapolated to other UDGs, it is necessary to test it on a larger sample of UDGs.These would preferably be spread across different environments and span a variety of GC numbers so that selection effects can be diminished and conclusive interpretations made.
In this study, we extend the work of B22 to a sample of 59 UDGs in the MATLAS low-and moderate-density environments (Marleau et al. 2021), to help balancing the sample that was biased towards higher density environments in B22, and to thus start building up a representative sample of UDGs.38 out of these 59 MATLAS UDGs have GC counts from Hubble Space Telescope (HST) imaging.Similar to B22, we employ optical to mid-infrared spectral energy distribution (SED) fitting to explore the stellar populations of these MATLAS UDGs.Buzzo et al. in prep. will combine both samples (i.e.,B22 and this study) with a control sample of classical dwarf galaxies to perform an statistically meaningful and comprehensive study of UDGs spread across environments.
The paper is structured as follows: in Section 2 we present a summary of the UDG sample studied in this work and the data available for each UDG.In Section 3 we describe our methods to obtain the photometry and morphology of the galaxies and our SED fitting methodology.In Section 4 we provide our results.In Section 5 we discuss the implications of our results within the theoretical predictions for UDGs and as compared to the literature.In Section 6 we present the summary and the conclusions of the paper.

DATA SAMPLE
The MATLAS UDGs were identified by the deep optical imaging of the large observing program Mass Assembly of early-Type gaLAxies with their fine Structures (MATLAS) survey (Duc et al. 2014(Duc et al. , 2015;;Duc 2020;Bílek et al. 2020;Poulain et al. 2021;Marleau et al. 2021).The sample of 59 UDGs used in this work was selected by Marleau et al. (2021), using the same UDG criteria as described above (i.e.,  ,0 > 24 mag arcsec −2 ,  e > 1.5 kpc).These UDGs are all in low-to-moderate density environments.There are no UDGs in highdensity environments such as clusters of galaxies in our sample.The environment that the galaxies reside in was calculated using a K-nearest neighbour algorithm.Using the 10 nearest galaxies to the UDGs, the local surface density (Σ 10 ) and local volume density ( 10 ) are calculated as a proxy for the environment.It therefore does not separate central from satellite UDGs.More details about the environment determination can be found in Duc et al. (2014) and Marleau et al. (2021).
For simplicity, throughout this paper we follow the MATLAS convention and refer to the closest massive galaxy to the UDGs as hosts, but that does not necessarily imply that the UDGs are satellites of the massive nearby galaxy.
No selection against star-forming UDGs was made, meaning that our sample has both star forming and quiescent UDGs, differently from the sample of B22 and FM23, who focused only on quiescent UDGs.
As previously mentioned, the identification of the low surface brightness sources was carried out on the imaging fields around massive early-type galaxies (ETG) in the nearby Universe.The UDGs are, thus, classified as being part of the group/in the vicinity of the massive ETG.However, since the UDG definition is based on physical size, it is imperative to understand the real distance to the galaxy to confirm whether or not it is a bonafide UDG.Few galaxies to date have been followed-up spectroscopically to test if they are indeed associated with the group or suffering from projection effects.The largest follow-up was carried out by Heesters et al. (2023) who studied 56 MATLAS dwarfs, out of which 3 were UDGs (MATLAS-585, MATLAS-2019 andMATLAS-2103).They found that 75% of the studied galaxies were at the same distances as their hosts, including the three studied UDGs.In this study, in addition to the SED fitting, we present spectroscopic follow-up of 3 MATLAS UDG candidates with Keck/DEIMOS (MATLAS-342, MATLAS-368 and MATLAS-1059) to recover their redshifts and probe their association with the group/massive neighbour.The spectroscopic observations and results are described in Appendix A. We find that all three UDGs are at the same redshifts as their hosts, in agreement with the findings of Heesters et al. (2023).We therefore consider it reasonable to assume that the rest of the sample is also at the same redshift as their hosts.We further discuss the assumed redshifts for the MATLAS UDGs in Section 4.1.
The description of the imaging data used to perform SED fitting is in Section 2.1.The data used to perform the GC number counts are described in Section 2.2.

Photometric data
In this work, we use data from the optical to mid-infrared to study the stellar population properties of 59 UDGs in the MATLAS survey.
Below we present the data used for each galaxy, along with how the photometry was measured in each band.We note that we tried to recover the photometry in the ultraviolet for the MATLAS UDGs using Galaxy Evolution Explorer (GALEX) imaging.However, all of our measurements turned out to be fainter than the 3 limit of the survey, independently if the imaging was observed with the Deep (DIS, 3 DIS (NUV) = 24.9mag), Medium (MIS, 3 MIS (NUV) = 23.2mag) or All-Sky (AIS, 3 AIS (NUV) = 21.3 mag) survey (Martin et al. 2005;Morrissey et al. 2007).
Optical, near-and mid-IR magnitude measurements are in AB magnitudes and were corrected for Galactic extinction using the two-dimensional dust maps of Schlegel et al. 1998(recalibrated by Schlafly & Finkbeiner 2011) and the extinction law of Calzetti et al. (2000).

DECaLS optical imaging
Although the MATLAS UDGs were identified using Canada-France-Hawaii Telescope (CFHT, surface brightness limit of 28.5 -29 mag arcsec −2 in the -band, Duc et al. 2014;Poulain et al. 2021) data, in this study, we use data from the Dark Energy Camera Legacy Survey (DECaLS, surface brightness limit of 28.5 mag arcsec −2 in the -band, Li et al. 2022) to perform the SED fitting.This is because DECaLS has imaging available in the , ,  and  bands for most galaxies in our sample, while CFHT has only the  and  bands (with a few rare cases where the  band is also available).We caution, nonetheless, that the DECaLS coadded imaging is noisy in the outskirts of galaxies (Li et al. 2022), which may lead to slight differences in the recovered photometry and morphological parameters when compared to CFHT, although consistent within uncertainties.
Differently from B22, in this study we do not use aperture photometry of our galaxies.Interested in obtaining total magnitudes in addition to morphological properties of the galaxies, we use the multi-wavelength galaxy fitting (GALFITM, Häußler et al. 2013;Vika et al. 2013) routine to study the galaxies in the optical.Detailed explanation of the input parameters and configuration of GALFITM are given in Section 3.1.
For all of the galaxies, archival optical coadded data were obtained from the Dark Energy Camera Legacy Survey (DECaLS, Dey et al. 2019).36 out of our 59 UDGs were observed as part of the DECaLS Data Release 10, which includes the , ,  and  bands.The remaining 23 galaxies were observed as part of the DECaLS Data Release 9 and have imaging available in the ,  and  bands.The reduction and calibration of the DECaLS data are described in Dey et al. (2019).
As discussed in B22, the coadded DECaLS data have shallower depths and more uncertain sky subtractions than other optical surveys focused on low surface brightness galaxies.Nonetheless, tests carried out by B22 have shown that the photometry obtained with the DECaLS coadded data are consistent within 1  with the photometry obtained by Lim et al. (2020) and Pandya et al. (2018) using deeper data reduced with an LSB-appropriate pipeline.The photometry recovered from DECaLS is, however, on average 0.1 magnitudes fainter than the ones in the literature.This systematic difference was incorporated into our final magnitudes as described in B22.The final photometry in the optical for all MATLAS UDGs was compared to that obtained by Poulain et al. (2021) using CFHT data and was found to provide similar results.

WISE near-IR and mid-IR imaging
The Wide-field Infrared Survey Explorer (WISE, Wright et al. 2010) is a space telescope that has imaged the entire sky in four filters with effective wavelengths of 3.4, 4.6, 12 and 22 𝜇m (near to midinfrared).For this study, we gathered WISE data for all the galaxies in our sample.These data are a mix of archival ALLWISE data and bespoke data construction and analysis, including custom mosaic construction from WISE single frames.The reduction, calibration and photometric measurement processes are thoroughly described in B22.We use WISE in all of our SED fits as it has been shown by B22 to significantly improve the stellar population results and to help breaking the age-metallicity degeneracy.All photometric measurements can be found in Table B1.

Globular Cluster Numbers
Total GC numbers have been obtained for 38 galaxies in our sample using Hubble Space Telescope (HST)/ACS imaging.The reduction, source detection, GC candidate selection and final GC counts are thoroughly discussed in Marleau et al. subm.Here we use the total (background and completeness corrected) GC numbers of the MAT-LAS UDGs to test the trend between GC-richness and metallicity found by B22 on our more complete sample.In B22 and FM23, UDGs considered to be GC-rich were those with more than 20 GCs, as this value roughly corresponds to a halo mass of 10 11  ⊙ (Burkert & Forbes 2020).UDGs with less than 20 GCs were considered GC-poor.In this study, we do not impose this hard separation at  GC = 20.Instead, we analyse GC numbers in a continuous manner so that evolutionary trends can be identified.

GALFITM
In this study, we use GALFITM (Häußler et al. 2013;Vika et al. 2013) to perform simultaneous galaxy fitting of our optical bands.This technique allows the combination of lower signal-to-noise (S/N) images, such as in the -band, with images with high S/N, and balances them out, resulting in consistent and well-constrained parameters.
We note that GALFITM was only run on the optical images, not on the WISE ones.This is because our methodology (described in Buzzo et al. 2022;Jarrett et al. 2012Jarrett et al. , 2013Jarrett et al. , 2019) ) was already shown to provide reliable total magnitudes for the WISE imaging.
For all galaxies, a single Sérsic model was fitted in the optical bands, using as initial guesses the morphological parameters obtained by Poulain et al. (2021).In our GALFITM setup, the magnitudes were allowed to vary freely amongst bands.Alternatively, the effective radius ( e ), Sérsic index (), axis ratio (/) and position angle (PA) were allowed to vary, with the condition that they are constant in all bands.A discussion on appropriate degrees of freedom to use for each free parameter in GALFITM is provided by Buzzo et al. (2021).
The sky value for each wavelength band was obtained from the image header.To run GALFITM, we created synthetic PSFs of each image using the PSFex routine (Bertin 2002).The parameters used to create such PSFs were obtained using the Source Extractor routine (SExtractor, Bertin 2002, version 2.19.5).We masked all interlopers, including GCs and nuclear star clusters, using the segmentation image output of SExtractor.The routine was run in single-mode in each band.We defined a source detection minimum area of 3 pixels and a threshold of 3 above the background to detect sources.SExtractor was run with an automatic background fitting, with a background cell size of 64 pixels.
For the masking process, we visually inspected every mask to ensure no star forming region was being masked as this could bias the final results to older stellar populations.To do this, we carefully looked into GALEX data on every one of the galaxies to see if any star formation was detected.As previously mentioned, for none of the galaxies UV emission was detected, indicating that the process was not masking any star forming region.
In Fig. 1 we show a GALFITM model for the four optical bands of the galaxy MATLAS-2019 (NGC 5846_UDG1).This example shows the highest mask fraction applied to our sample, as MATLAS-2019 is the galaxy with the most GCs within the sample of MATLAS UDGs.All other galaxies have had none to very few sources masked in front of them, therefore not compromising the final recovered photometry of the galaxies.GALFITM models for the other galaxies are of similar quality to the one shown in Fig. 1.
Results from these fits are given in Table B2.We find that the MATLAS UDGs have a median Sérsic index of  = 0.95 ± 0.14, and a median axis ratio of / = 0.62 ± 0.13.Poulain et al. (2021) also obtained morphological parameters for all of the MATLAS dwarfs.Their results for the MATLAS UDGs show a median Sérsic index of  = 0.79 ± 0.10 and a median axis ratio of / = 0.67 ± 0.14.Thus, their morphological findings for the MATLAS UDGs are consistent with ours.GALFITM, similarly to GALFIT, provides unrealistically small uncertainties.In Appendix C, we discuss some of the implications of using these underestimated uncertainties to the results obtained with SED fitting.
A full and thorough description of the main configuration and models used in Prospector is available in B22.In the latter, Prospector models including dust attenuation provide results closer to those found with spectroscopy for UDGs.Thus, in this study, we only use models that include dust attenuation, assuming the Gordon et al. (2003) attenuation curve.Whether this dust is artificially added to improve the models or the galaxies actually have dust can only be probed using far infrared or ultraviolet deep data, which will be pursued in the future.We note, nonetheless, that tests with Prospector found a dust attenuation of   ∼ 0.1 − 0.2 mag for Milky Way GCs when there should be none, suggesting an artificial addition of dust in Prospector models (Johnson et al. 2021b, fig. 6).
For our fits of all galaxies, we include upper limit fluxes coming from the 12 and 22 m bands from WISE.For cases of extreme low S/N in the -band of DECaLS, this band was also incorporated in the fitting as an upper-limit.For further understanding of how Prospector deals with upper limits, see Appendix A of Sawicki (2012).The inclusion of the upper limits was shown in Buzzo et al. (2022) to be helpful in constraining the amount of dust attenuation found for these galaxies.In the cases where no dust is found, the upper limits do not play an important role in the fit.We assume a delayed exponentially declining star formation history (SFH), as this is consistent with the recovered shape of the SFHs of UDGs by FM23.The delayed model is an extension of the regular exponentially declining SFH.The latter usually assumes that star formation jumps from zero to its maximum value at the time  age and then declines exponentially within a timescale .The delayed model, alternatively, multiplies the exponentially declining SFH by the time since  age .This is capable of removing the discontinuity in the SFH at  age and the condition that star formation can only decline after that point.This results in a more flexible, robust and physical SFH (see Carnall et al. 2018;Leja et al. 2019, for thorough discussions and comparisons between parametric star formation histories).Thus, the form of the chosen SFH is: According to this definition,  age measures the onset of star formation within a galaxy and  is an approximation of the star formation timescale, i.e., how long does a galaxy take to quench after reaching peak star formation.In this study, we are interested in the mass-weighted ages (  ) for the galaxies, which are calculated analytically within Prospector using as input parameters the age since the onset of star formation ( age ) and the star formation timescale ().We note that B22 has incorrectly reported  age as the mass-weighted age, not the true   calculated analytically with Prospector.  is only strongly different from  age in the cases where  is long, which is usually not the case for UDGs.Thus, the bulk of the comparisons between  age from B22 and spectroscopic mass-weighted ages from the literature hold.However, some individual cases may have larger differences and be thus less consistent with spectroscopy than reported.Buzzo et al. (in prep.) is working on combining the dataset in this study with the one from B22, where the photometry for all galaxies will be obtained consistently with GALFITM and the true   will be reported for every galaxy so that we have a homogeneous dataset, where fair comparisons can be made and conclusions can be assessed.
For this study, we use two different configurations in Prospector.
For these two scenarios, we placed linearly uniform priors on our free parameters.The range of the priors in our models are determined by the coverage of the Padova isochrones used in FSPS (Marigo & Girardi 2007;Marigo et al. 2008).These are: log( ★ / ⊙ ) = 6 -10, [M/H] = −2.0 to 0.2 dex,  = 0.1-10 Gyr,  age = 0.1-14 Gyr,   = 0 − 4.344 mag, and redshift  = 0-0.045.This redshift range translates to a luminosity distance range of 0 <   < 200 Mpc.Since all of the galaxies were identified in groups out to 45 Mpc, we assume that this redshift range is representative of the distance of all UDGs.We note that the default stellar masses obtained from Prospector are the total masses formed in the galaxy throughout its life.These were converted to surviving masses (most commonly quoted in the literature) in this work using a builtin function within Prospector.
As mentioned in B22, we emphasize that different prior assumptions do not significantly alter our results or conclusions, with the exception of the prior assumption on the shape of the SFH, as discussed by Webb et al. (2022).

RESULTS
Here we present our results on the stellar population properties of the 59 MATLAS UDGs studied in this work.As explained in Section 3.2, we carried out two model configurations in Prospector, one where the redshift is a free parameter and one with the redshift fixed.In this Section we show our results for both configurations, discuss the differences between them and the implications of using one or the other.
The resulting SEDs from both configurations of every galaxy were inspected to ensure that every photometric point was well modelled and that the shape was as expected.We found a median  2 reduced of 1.41 ± 2.65.The best fits are for MATLAS-2103 and MATLAS-585, both with  2 reduced = 1.01 and the worst is for MATLAS-1589, with  2 reduced = 7.34.Fig. 2 shows a typical fit of an UDG in our sample, with  2 reduced = 1.32.Upper limits were not included in the calculation of the  2 .Figs. 3 and 4 show the distribution of the parameters for all of the galaxies in our sample in both configurations.
It is important to notice that choosing a parametric SFH for the galaxies comes with many limitations.One of such is that the star formation histories of the galaxies may have multiple components, i.e., an old quiescent and another with a recent burst of star formation.In such cases, a model with a parametric SFH that smooths over two distinct star formation episodes will be the most statistically correct, but will not reflect the real stellar populations present in the galaxies.Fitting the galaxies with a non-parametric SFH may be a solution for this problem, but it adds another one: too many free parameters for a limited amount of data, likely resulting in overfitting.To try and avoid this type of problem, we have gathered GALEX data for all of the MATLAS UDGs to be able to correctly fit the star forming regions if the galaxies had any.None of the galaxies in our sample were detected in the FUV or NUV bands of GALEX, and provided only very faint upper limits.We have tested including these UV upper limits on a subsample of our galaxies to understand the effect they would have on the recovered stellar populations.For all galaxies, the results are consistent with the ones obtained without the usage of GALEX, indicating that Prospector is correctly tracing the dominant stellar populations within the galaxies only by using the optical to mid-IR data.

Comparison between models with free and fixed redshifts
In Fig. 2, we show a comparison between the output best-fit SED and covariance matrix for the two configurations of Prospector applied to one of the galaxies in our sample, MATLAS-2019.It can be seen that the model where the redshift is free has a much larger spread in the determination of the stellar mass, while the results with the redshift fixed seem to be better constrained.Apart from the stellar mass, the other parameters do not seem to be strongly affected by the freedom of the redshift.
The same behaviour can be observed if we look into the whole population of galaxies, instead of an individual case.Figs. 3 and 4 provide a comparison of all of the stellar population properties of the MATLAS UDGs obtained with both configurations.As can be seen in these figures, the metallicity, star formation timescale, age and dust attenuation are not strongly affected by the freedom of fitting or lack thereof to fit redshift in the models.Galaxies are found to be on average younger and more metal-rich in the models where the redshift is fixed, but consistent within the uncertainties with the stellar populations found in the models with the redshift free.However, the stellar mass of the galaxies shows a larger difference when comparing models.Overall, galaxies have higher stellar masses in the models where the redshift is a free parameter.That is because in integrated SEDs at non-cosmological distances, the redshift and stellar mass are highly degenerate parameters.
As explained in Section 2, we expect that most UDGs in our sample are at the same redshift as their hosts in the models with the redshift free, if our redshift estimates are reliable.We explore the difference between the redshift of the host and our recovered photometric redshift (  ) in Fig. 5.As can be seen, while there is a show the correlation between parameters.This image can be read and interpreted as a covariance matrix.Columns stand for redshift, stellar mass, metallicity, star formation timescale, mass-weighted age and interstellar diffuse dust extinction.Metallicity, star formation timescale, age and dust attenuation are similar in both configurations.However, the stellar mass is strongly affected by the freedom of the redshift, as these two parameters are highly degenerate.
peak at   = 0, a significant portion (> 70%) of the UDGs have |  |> 0.0067, equivalent to a difference in velocity of |  |> 200 km s −1 ).Thus, they do not have photometric redshifts consistent with those of their hosts.Most of them are in fact found in their background (i.e.,   < 0).Since photometric redshifts at non-cosmological distances are very uncertain, we consider   being different than zero as a confirmation that redshifts with this method cannot be trusted rather than an indication that MATLAS UDGs are not at the redshifts of their hosts.
Going forward, we treat the models with fixed redshift as the best representation of our galaxies.Thus, from now on, we will only discuss the results of this configuration.We reinforce, nonetheless, The bottom-right panel shows the age-metallicity distribution of the 59 MATLAS UDGs in both configurations.Histograms at the top and on the left show the marginal distributions of the data displayed on the age-metallicity plane.We conclude that models with the redshift free (purple) recover on average younger ages and more metal-rich populations for the UDGs.that even if the galaxies are not at this assumed distance, changes in the stellar populations (except stellar mass) between models are small and well within their uncertainties (see Figs. 2, 3, and 4).We carefully checked that the changes in stellar mass do not strongly affect any of our conclusions that follow.Specifically, we have checked that the mass-metallicity bimodality found for the models with fixed redshift and thoroughly discussed in Section 5.1 is also present in the models where the redshift is a free parameter.We have also checked that the change in stellar mass is small and not enough to change the classification of the UDGs according to their positioning on the mass-metallicity plane (i.e., which MZR they follow).

Median stellar populations
The results of the models with the redshift fixed at the redshift of the host galaxy are fully presented in Table B3.The median and median absolute deviations of the stellar populations of the whole sample of MATLAS UDGs are presented below.We find that the MATLAS UDGs have intermediate-to-old ages, with a median mass-weighted In all panels, purple colours stand for the models with the redshift free and green colours are the models with the redshift fixed to the host.We conclude that models with the redshift free (purple) result in slightly more massive UDGs and a much narrower distribution centred at log(  ★ / ⊙ ) ∼ 7.8 than models with the redshift fixed.No significant difference is observed in the star formation timescales and dust attenuation recovered in both configurations.age of   = 7.1±1.8Gyr.The galaxies display an average metal-poor population with [M/H] = −1.2± 0.2 dex.These stellar populations are equivalent to those found for UDGs with spectroscopy by FM23.We find an average  of 1.6 ± 0.7 Gyr.Finally, the median interstellar diffuse dust attenuation coming from the SED fitting of the galaxies in our sample is   = 0.12 ± 0.07 mag, consistent with the expectation that these galaxies should have little-to-no dust component and with the caveat that the non-zero dust content could be a small artificial addition that compensates for inaccuracies in the SED models.It is interesting to notice that although our dust priors extend out to 4.3 mag, the highest   value found is 1.17 mag, highlighting the importance of the inclusion of the WISE upper limits from the 12 and 22m bands to constrain the amount of dust in the galaxies.This power of the near-and mid-infrared bands in constraining the dust amount found had been discussed previously both by Pandya et al. (2018) and B22.Nonetheless, as expected, UDGs with higher amounts of dust show higher degrees of degeneracy between the dust extinction and the metallicity and age.
Trends between the recovered stellar populations and the galaxy morphologies, environments, scaling relations and GC-richnesses are investigated in Section 5.1.

Comparison with the literature
In this Section, we compare our SED fitting results with stellar populations recovered with spectroscopy in the literature for some of the MATLAS UDGs.
MATLAS-2103 was found by Heesters et al. (2023) to have a metallicity of [M/H] = [Fe/H] = −1.69+0.01 −0.22 dex and an age of 11.3 +2.2 −2.5 Gyr.In this study, we find that MATLAS-2103 has a metallicity of [M/H] = −1.40+0.14  −0.13 dex and an age of 9.79 +2.81 −3.84 Gyr.Similar to the case of MATLAS-585, we find ages consistent with those of Heesters et al. (2023), but their stellar populations are more metal-poor than ours.Notably, Heesters et al. (2023) tends to find systematically more metal-poor stellar populations for galaxies if we take this study, FM23 and Müller et al. (2020b) as baselines.
Overall, we show that our results are consistent with those found using spectroscopy.Our SED fitting technique had already been shown to provide consistent results with spectroscopy in the studies of B22 and FM23.In both cases our results were shown to have a metallicity offset of ∼ 0.2 dex, which is well within our uncertainties.We reinforce that while SED fitting results have much larger uncertainties, stellar populations obtained from broad-band photometry are far less restricted than spectroscopic observations in depth and spatial coverage.Additionally, the ability of this technique to gather information from a wide wavelength range makes up for the loss of detailed features provided by spectroscopic data.SED fitting also allows for complete sample studies such as this one, which are prohibitively expensive to do spectroscopically on big telescopes.Apart from the larger uncertainties, another downside of SED fitting is the difficulty in obtaining [/Fe] and dynamical information from the galaxies.Thus, spectral and photometric information are complementary and should be jointly used to better understand UDGs.

DISCUSSION
As mentioned in the Introduction, the stellar populations of UDGs can help understanding their formation histories.This can be made by means of comparing these populations to those of different classes of galaxies to try and identify similarities/differences between them.In Section 5.1, we use scaling relations (i.e., the mass-metallicity relation) to compare UDGs to classical dwarfs (observations) and high-redshift galaxies (simulations) to try and elucidate their origins.In Section 5.2 we discuss if there is evidence for two classes of UDG based on their stellar populations and morphological properties and associate these classes with their most likely formation scenarios.

Mass-metallicity relation: Clues to the origins of UDGs
We explore in this section the positioning of MATLAS UDGs on the stellar mass -metallicity plane, especially in comparison with known mass -metallicity relations (MZR).We investigate trends with different morphological and stellar population parameters, as well as GC-richness.
For reference, in Figs. 6 and 7, we plot two MZRs: Simon (2019) for classical dwarfs, and the simulated MZR from Ma et al. (2016) for evolving sources at a redshift of  = 2.2.To be consistent, we plot total metallicities from Prospector on top of the converted to [M/H] Simon (2019) relation, from its original in [Fe/H] (Vazdekis et al. 2015).A full discussion on how to make this conversion and the caveats that arise in doing so are given in B22.In Fig. 6 we show the distribution of the MATLAS UDGs in the mass-metallicity plane in comparison to three other works in the literature.In blue we show the study of 100 field UDGs using SED fitting from Barbosa et al. (2020).The study of 29 UDGs in a variety of environments also using SED fitting from B22 is shown in magenta.Lastly, the largest spectroscopic study of UDGs to date (25) from FM23 is shown in green.The UDGs in FM23 are in multiple environments, but biased towards clusters.
The MATLAS UDGs, at first glance, seem to form a new MZR, with a metallicity that is anti-correlated with the stellar mass.However, the population shows a bimodality, indicating that instead of an anti-correlation, the MATLAS UDGs split into two samples.One of them is well within the classical dwarf MZR and another has lower metallicities and lies closer to the simulated high- MZR.To make this bimodality clearer, in the bottom panels of Fig, 6, we show the residual MZR (i.e., the difference between the metallicities of the UDG population and the value expected from the MZR).It can be seen that the low-mass metal-rich UDGs are well within the scatter of the classical dwarf MZR, while the high-mass metal-poor UDGs are further down and in fact are more consistent with the MZR at highredshift.This bimodality may be indicative of different formation scenarios for the UDGs populating each mode.
When comparing to the works in the literature, we see that while all have a population of UDGs that follows well the classical dwarf MZR, they all seem to probe different regions of the parameter space.Barbosa et al. (2020), as well as the MATLAS UDGs, are biased towards lower mass UDGs in low density environments, while B22 and FM23 are biased towards more massive cluster UDGs.FM23, in addition, probes a population of more massive, metal-rich, cluster UDGs, although overall their distribution scatters around the MZR relations at log( ★ / ⊙ ) ∼ 8. Barbosa et al. (2020) find many UDGs to be more metal-rich than the classical dwarf MZR.These were found to be the youngest in their sample.FM23 found a similar subpopulation of UDGs more metalrich than the classical dwarf MZR, as well as another subpopulation that is more metal-poor.We recover the population of metal-poor UDGs found by them, but no UDGs in the same mass range that are metal-rich.The explanation for this is that the sample of FM23 is biased towards the cluster environment where tidal effects are more common (Sales et al. 2022) and result in the higher metallicities found.
Finally, the sample covered by B22, although more representative than the current sample in terms of environments, is more constrained in terms of stellar masses, showing only a small range around 8.0 ≤log( ★ / ⊙ ) ≤ 8.5.Similar to the findings of FM23 and the ones of the current study, they found that a subsample of the UDGs follows the classical dwarf MZR, while another population is much more metal-poor and is more consistent with the high- MZR.

Dependence with morphology and local environment
The fact that the morphology of galaxies is correlated with the environment they reside in has been known since the 1970's, giving origin to the so-called 'Morphology-Density Relation' (Oemler 1974;Davis & Geller 1976;Dressler 1980;Dressler et al. 1997), first proposed by Edwin Hubble in the 1930's (Hubble & Humason 1931).In a very simplistic explanation, early-type galaxies are found to occupy the densest environments such as the centres of clusters of galaxies, while late-type galaxies are observed in their outskirts and to be prevalent in less dense environments, such as the field.The morphology-density relation was also found to hold for classical dwarf galaxies (Phillipps et al. 1998;Penny et al. 2011) Various works in the literature have looked into the morphological parameters of UDGs (e.g., Rong et al. 2020;Kado-Fong et al. 2020, 2021;Van Nest et al. 2022;Li et al. 2023), such as their ellipticities, which can hold much information about the physical processes involved in forming UDGs.Rong et al. (2020), for example, found that the morphologies of cluster UDGs depend on their luminosities, environments and redshifts.In terms of axis ratios, they showed that there are at least two populations of UDGs, the elongated (/ ∼ 0.4) and the round ones (/ ∼ 0.9).In their study, UDGs at higher redshifts ( ∼ 0.3) are more elongated than their low-redshift counterparts, indicating a possible evolutionary path of UDGs becoming rounder with time.They also show that bright (  > −15.2 mag) UDGs are rounder than the fainter ones (  < −15.2 mag).Finally, they showed that UDGs closer to the centre of clusters are rounder than the ones in the outskirts or outside the cluster virial radius.Their findings apply to both nucleated and non-nucleated UDGs, although it is worth noting that the number of nucleated UDGs increases towards the centre of clusters and these are the roundest sources.A bimodality in the axis ratio, mass, and luminosity of UDGs has also been reported in the Coma and Virgo clusters by Lim et al. (2018) and Lim et al. (2020), respectively.They find that UDGs that are dark matter dominated (M/L > 1000) have relatively rounder shapes (higher /), while UDGs with M/L ∼ 500 are more elongated (low /).They also found that the round UDGs have higher GC specific frequencies than the more elongated ones.
In the top row of Fig. 7, we show the MATLAS UDGs colourcoded by their surface brightnesses, effective radii, axis ratio and local volume density.Similar to the findings of Rong et al. (2020), we also find a correlation between the surface brightness and the mass of the galaxies.Less massive UDGs are, on average, the faintest ones, while more massive UDGs are the brightest.This is expected given the narrow range in surface brightness and size, and the fact that the Sérsic index of the galaxies is close to 1.In addition, a similar correlation exists for the effective radii of the galaxies.Less massive UDGs are smaller, while the more massive ones are bigger.This behaviour is not unique to UDGs, it is in fact expected from the mass-size relation, but it is reassuring to confirm this expectation in our sample of galaxies.
In terms of their axis ratios, we can see a weak correlation with both the mass and metallicity, where more elongated UDGs seem to be positioned in a less massive/ more metal-rich region, while the roundest UDGs occupy a position corresponding to more massive, more metal-poor populations.There is a non-negligible portion of UDGs also with round morphologies located in the upper parts of the plane.One explanation for this is that the ellipticity of galaxies is highly dependent on their inclinations, which could drive oblate/flattened galaxies to look round if they are face-on.It is also worth noticing that more massive edge-on UDGs would likely be pushed over the surface brightness threshold definition for UDGs due to their higher apparent brightness, which can explain why we see such a correlation between the surface brightness, axis ratio and mass of the galaxies.Additionally, less massive galaxies tend to be rounder (e.g., Sánchez-Janssen et al. 2019), however, these faint galaxies are below our mass range, also explaining why we do not have a population of round galaxies at lower masses.
We note the study of Cardona-Barrero et al. (2020), who used the NIHAO simulations to suggest a correlation between the morphology and kinematics of UDGs.They suggested that more elongated UDGs (i.e., smaller axis ratios) are rotationally-supported, while the rounder ones (higher axis ratios) are pressure-dominated.These findings are further supported for classical dwarfs by Pfeffer et al. (subm.).
If we now compare these findings with the panel where the UDGs are coloured according to the density of the environment they reside in (i.e., smaller values of log  10 are for galaxies that are fairly isolated, while higher values signify that galaxies are in more dense environments), we can start to see some forms of the morphologydensity relation appearing.That is, UDGs that show the most elongated morphologies are the ones located in the less dense environments.Alternatively, UDGs that are rounder tend to be in the most dense environments.If we consider that UDGs in denser environments are prone to more mergers and encounters than their isolated counterparts, then having rounder morphologies is a natural consequence of that (Moore et al. 1996(Moore et al. , 1998(Moore et al. , 1999)).We caveat that this interpretation is done for the full sample of galaxies, and individual cases may differ from the ensemble properties.
If we analyse the environment panel on its own, there appears to be a weak trend between the environment and the positioning of the UDGs, where the majority of the UDGs in the most dense of our probed environments are below the classical dwarf MZR.We remind the reader that all of the MATLAS UDGs are, by definition, in low-to-moderate density environments, and thus, the range of densities probed is limited and does not include the most dense environments, such as clusters.Selection effects can be playing a role in our findings and, thus, comparisons with UDGs located in higher density environments are necessary and will be pursued in a follow-up study (Buzzo et al. in prep.).It is interesting, nonetheless, to see that even with such a small range of densities, a correlation between the environment and the stellar populations can still be weakly identified.Although a correlation was identified with the local volume density, we found no correlation between any properties and the projected groupcentric distance of the UDGs.We emphasise that these are the interpretations for the ensemble of MATLAS UDGs, but as can be seen, individual exceptions may exist (e.g., round UDGs at low-density environments or elongated ones at the highest probed densities).We use colour-coding to show the changes in the properties of the MATLAS UDGs across the mass-metallicity plane in Fig. 7, but for an alternative view, we provide in Fig. D1 direct plots of each property against their positioning on the plane to emphasise the discussed trends.

Dependence with age and star formation timescale
We now turn to understanding how the positioning of the UDGs in the stellar mass -metallicity plane depends on their star formation history properties, e.g., their ages and star formation timescales.
Ferré-Mateu et al. ( 2018) have shown that there exists an age dependence with projected clustercentric distance, i.e., UDGs are younger at larger projected clustercentric radii.Alabi et al. (2018) and Kadowaki et al. (2021) also found a colour dependence with the environment, with bluer UDGs residing in lower-density environments than the redder ones.Finally, B22 found, using SED fitting, that UDGs in the field are on average younger and more metal-rich than cluster UDGs.Group UDGs in their work were found to be in a transitional space in the age-metallicity plane.However, these were based on very small number statistics, hence the need for this current study focused mainly on group UDGs.Similar to these recent findings, the ages of the MATLAS UDGs also seem to distill into two main types: younger UDGs following the dwarf MZR and older ones following the high- MZR.We however caution the reader that ages recovered with Prospector are very sensitive to the version of the code used (see Section 3.2), as well as to the uncertainties in the optical photometry (see Appendix C).
One interesting recent finding on the stellar populations of UDGs by FM23 is a correlation between the  enrichment ([/Fe]) and the star formation histories of UDGs.Using a phase-space diagram, which may be particularly related to the environment of the UDGs, they show that the UDGs with the highest [/Fe] are in the region corresponding to early infallers, while UDGs with the lowest [/Fe] tend to be in the outskirts of clusters or in less dense environments.The -enrichment can be associated with the speed at which star formation starts and quenches within galaxies.Thus, this parameter can very well separate early-quenched galaxies (high [/Fe]) from galaxies with prolonged star formation histories (low [/Fe]).
With SED fitting we cannot recover [/Fe].We can, however, look into the star formation timescale of galaxies, i.e., the rate of the exponential decay of star formation, to understand how fast SF was quenched within the galaxy.In the bottom row of Fig. 7, we colour-code our UDGs by their mass-weighted ages (  ) and star formation timescales ().We can see that UDGs that have the fastest star formation episodes (i.e., smallest ) are in the lowest regions of the plane, coinciding with the positioning of the older galaxies in the age panel.They are consistent with the evolving MZR at highredshift, indicative of early-quenching, which is in agreement with the results from FM23, and suggestive of a failed galaxy-like origin.On the other hand, UDGs that show prolonged SFHs (i.e., large ) are the ones that are younger and more metal-rich, aligned with a puffed-up dwarf formation scenario, as also found by FM23.

Dependence with GC richness and specific frequency
One of the properties of UDGs that draws a lot of attention is the large number of GCs that some of them host in spite of their low stellar In parallel, observational work has focused on understanding what are the main differences between GC-poor and GC-rich UDGs, in an attempt to correlate the higher fraction of GCs with other galaxy properties (Forbes et al. 2020;La Marca et al. 2022).
In terms of stellar populations, B22 used SED fitting techniques to find evidence for an inherent difference between GC-rich and GCpoor UDGs.They found that GC-poor UDGs are consistent with the classical dwarf MZR from Simon (2019), having consistently higher metallicities than their GC-rich counterparts.Interestingly, all of the UDGs that are consistent with the evolving MZR at  = 2.2 proposed by Ma et al. (2016) are GC-rich, suggesting a formation scenario relying on early quenching for these galaxies.Following up on this finding, FM23 studied the stellar populations of 25 UDGs, the largest spectroscopic sample to date, using data from Keck/KCWI.They found that GC-poor UDGs consistently follow the classical dwarf MZR.However, they found that the GC-rich UDGs can be found both above and below the Simon (2019) MZR.They suggested that the ones above it are likely tidal UDGs, while the ones below the MZR are consistent with early-quenching (failed galaxy scenarios), similar to what was suggested by B22.
In the bottom row of Fig. 8 we show the MATLAS UDGs colourcoded by their number of GCs.Differently from B22 and FM23, we do not impose the hard separation that GC-rich UDGs are the ones with  GC ≥ 20 and GC-poor UDGs the ones with fewer.Instead, the colour-scheme has been discretized so that evolutionary trends can be identified.In this figure, we see that GC-poor UDGs are overall more metal-rich galaxies, consistent with the classical dwarf MZR, and consistent with the recent findings of Jones et al. ( 2023) that gas-rich, field UDGs host few GCs.Alternatively, we see that all UDGs that follow the evolving MZR at high-redshift are GC-rich.We do not find the population of GC-rich metal-rich cluster UDGs found by FM23, which may be expected as all of our UDGs are in groups/field, and thus less prone to tidal effects (see Fig. 7).
These results are in agreement with the recent findings of Pfeffer et al. (subm.),who used the EAGLE hydrodynamical simulations to find a correlation between the axis ratio and GC-richness of dwarf galaxies.In their study, round (higher /) dwarfs host more GCs than their elongated (low /) counterparts.As a a caveat, the works of Smith et al. (2013) and Smith et al. (2015) have shown that it is not possible to make a galaxy rounder through tidal stripping while also retaining a large number of GCs.If the harassment and stripping are strong enough to significantly alter their shapes, then the bulk of the GC system would be lost as well.In this study, however, we refer to galaxies that are round and GC-rich from their early stages.As a sidenote, the effective radius of galaxies, not only UDGs, has been shown to correlate with the number of GCs by Harris et al. (2013), also in agreement with our findings.
To understand if this trend between the metallicity and GC number is caused solely by the GC-richness or if it is being influenced by the mass of the galaxy, we also calculated the specific frequency per unit stellar mass of the galaxies,   (Zepf & Whitmore 1993).This specific frequency is defined as: (2) The MATLAS UDGs are colour-coded by their specific frequencies in the rightmost panel in the bottom row of Fig. 7.We can see now that the trend between GC-richness and metallicity is not so clear anymore.Although on average metal-poor UDGs have higher specific frequencies, there is a non-negligible population of more metal-rich UDGs that also have high   values.This indicates that the trend found between the metallicity and the GC number may be incorporating correlations with the stellar mass.It is important to note, nonetheless, that the UDGs that follow the high-redshift MZR all show high values of   .This suggests that at least part of the reason that we do not find GC-rich, metal-rich UDGs is our lower stellar mass (on average) than FM23.When controlling for the stellar mass of the galaxies using   , this population becomes apparent.These GC-rich metal-rich galaxies are likely tidal UDGs, as suggested by FM23.This finding is reinforced by the fact that one confirmed tidal UDG in our sample, MATLAS-1824, has a similar metallicity and   as these other GC-rich metal-rich UDGs.
Further supporting these findings, Lim et al. (2018) have studied the specific frequencies (  ) of UDGs in the Coma cluster and also found two different populations of UDGs in terms of their GC content.In their work, UDGs with high   were found to be more dark matter dominated, with mass-to-light (M/L) ratios > 1000, while UDGs with low   have M/L ∼ 500, consistent with the predictions of failed galaxy and puffed-up dwarf scenarios, respectively.Interestingly, a transitional population is also identified that has on average ∼ 10 GCs (10 <   < 200); and these do not follow either the classical dwarf MZR nor the high- MZR (see light blue points in Fig. 8).Seeing a continuous trend in the GC numbers immediately leads to questioning whether GC-poor UDGs could evolve into GC-rich ones (or vice-versa).However, these galaxies, if we put together all of the properties studied so far, seem to be inherently different, not only in terms of their GC systems.GC-poor UDGs are primarily found in the field or low-density groups, they are young, are elongated, and show prolonged star formation histories.GC-rich UDGs, on the other hand, are found in clusters or highdensity groups, are primarily round, have old to ancient ages, and have quenched long ago, having formed the bulk of their stars in early and quick star formation episodes.Given their different properties, an evolutionary path between the two populations seems unlikely and, thus, two (or more) formation scenarios seem to be necessary to Figure 9.The two groups of UDGs identified using the KMeans clustering algorithm.In this figure, the radial axis shows the central value of a given property in each class, while the angular axis shows each property being considered by the clustering algorithm.This polar plot is to be read as the relative difference in the properties of UDGs belonging to the Classes A and B. Class A is shown in green and is comprised of UDGs that have lower stellar masses, prolonged SFHs (i.e., higher star formation timescales), are more elongated (lower /), live in less dense environments (lower log  10 ), host fewer GCs, are younger (smaller   ), have higher  dwarf MZR (i.e., closer to zero and consistent with the classical dwarf MZR), and are fainter (higher  ,0 ).As these properties agree with the predictions of puffed-up dwarf-like formation scenarios, we associate Class A with UDGs having this origin.Class B is shown in pink and includes UDGs with higher stellar masses, rapid SFHs (i.e., smaller ), rounder (high /), living in the densest of our probed environments (higher log  10 ), hosting on average the biggest GC systems, older (higher   ), with smaller  dwarf MZR (i.e., negative values, being thus more metal-poor than the classical dwarf MZR), and brighter (smaller  ,0 ).As these properties are aligned with the predictions of failed galaxy-like formation scenarios, we associate Class B with UDGs having formed through this pathway.

Multiple classes of UDGs?
The results in the previous sections seemed to indicate that there exists at least two populations of UDGs.Driven by the idea that UDGs with different properties were formed by different processes, in what follows we investigate and compare the distribution and possible correlation between their global properties.Our goal is to identify classes of objects with similar properties and possibly associate them with one of the prominent formation scenarios for UDGs.
Assuming the existence of distinct populations of UDGs (Forbes et al. 2020), and in order to better associate each galaxy to one of them, we applied the centroid-based clustering algorithm KMeans (MacQueen et al. 1967) to the data.KMeans has been applied to Astronomy in various contexts and was proven to be a valuable methodology to classify large datasets (see e.g., Sánchez Almeida et al. 2010;Garcia-Dias et al. 2018;Coccato et al. 2022).The KMeans algorithm is an unsupervised learning method and is part of the Scikit-learn li- brary (Pedregosa et al. 2011), responsible for implementing machine learning algorithms in Python.
In our implementation, KMeans groups together galaxies that have similar properties in the multi-parameter space defined in our dataset and flags them accordingly into classes.The initial parameters considered and compared here are stellar mass (log( ★ / ⊙ )), star formation timescale (), mass-weighted age (  ), local volume density (log  10 ), number of GCs ( GC ), GC specific frequency (  ), central surface brightness ( ,0 ), and the residual metallicity between the UDGs and the classical dwarf MZR ( dwarf MZR ).The latter was included instead of metallicity in order to remove the degeneracy between metallicity and stellar mass so both properties could be analysed separately.Positive values of  dwarf MZR mean that UDGs are above the dwarf MZR (i.e. more metal-rich), 0 is when they follow exactly the MZR, and negative values are for the cases where UDGs lie below the MZR (i.e., more metal-poor).The Sérsic index () and effective radius ( e ), were not included in the analysis as they are degenerate with the surface brightness.Amongst our total sample of 59 galaxies, only 38 UDGs have all of these physical quantities available (i.e., 21 UDGs do not have GC counts).Therefore, these 38 are the ones used in our analysis.This is a small sample of UDGs and we caution that results obtained from it need to be further tested in a larger sample of galaxies, as these may not be representative of a diverse populations of UDGs and should not be extrapolated to such.All properties were linearly scaled down to a range between 0 and 1 so that good measures of distance between classes could be assessed for every parameter analysed.The KMeans clustering algorithm is isotropically applied to all parameters in the multiparameter space and, therefore, tends to produce round (rather than elongated) clusters.In this case, if unscaled data is provided (i.e., with unequal variances), then more weight is given to variables with smaller variance.Thus, rescaling the data is a necessary step to allow all properties to be equally examined by the algorithm.
We applied KMeans to the dataset allowing it to freely find the number of classes that best represents our data.We used the silhouette score quantity to measure the goodness of the clustering technique, where a high score means that the classes are well apart from each other and can be clearly distinguished.The score ranges from -1 to 1, where -1 means that sources are being associated with the wrong class, 0 means that classes are indistinguishable, and 1 means that the classes are perfectly separated.Two classes were found as the best number to describe our sample of UDGs, delivering the highest silhouette score (0.4) when compared to other numbers of classes.KMeans allocates every single data point to the closest class centre in the multi-parameter space, which means that all UDGs analysed in this study are allocated to one of the classes.This does not mean that every UDG belonging to a particular class will have all of the properties compatible with that class, but that overall, most of their properties would be consistent.These classes were further assessed using Principal Component Analysis (PCA), which showed that two clusters are able to recover 60% of the variance of the dataset.DBScan, another machine learning clustering algorithm, was also applied to the data and similarly returned that the best number of classes to represent this dataset is two, with a silhouette score of 0.3.
The two classes of UDGs are shown in Fig. 9.This plot shows the relative difference between the properties of the UDGs belonging to each class and provides an intuitive way of understanding which properties provide the biggest separation between classes, and properties that contribute less.From Fig. 9, we can see that the best discriminators of UDGs belonging to different classes are log( ★ / ⊙ ), ,   ,  dwarf MZR , and  ,0 .Properties that provide weaker separations are /, log  10 , and  GC .Several iterations of this plot were analysed to remove properties that did not contribute to the separations, or that were degenerate with others properties already included and were, therefore, not adding anything new to the separation.Properties excluded were  e ,  and PA for showing no significant difference between groups.At this point, it is important to remind the reader that in machine learning algorithms, cleaning the dataset from any properties that are not contributing to the learning process is one of the most important tasks to recover reliable results.One might think that leaving parameters that return a null separation will not make a difference, but these add noise to the clustering process and could pollute and compromise the final recovered results.We emphasise, nonetheless, that the recovered classes are consistent independently of the addition of these parameters, however, the classes are found to have much larger uncertainties as the null parameters contaminate the clear separation provided by the significant ones.On the other hand, it is important to analyse carefully the properties that are providing the biggest separation to ensure these are physically meaningful and not just a result of selection effects or the properties being degenerate/correlated.For example, the surface brightness, effective radius, and number of GCs are known to correlate with the stellar mass.We thus carried out tests excluding every single one of these properties from the clustering algorithm to check whether the classes would change or not.We found that the recovered classes are similar regardless if any of these properties is excluded or included in the process, they however result in larger uncertainties on the class centres.Overall, in all iterations of applying KMeans to our dataset (i.e., with the inclusion/omission of properties) the classes maintained, except for one iteration.The only property capable of changing the identified classes was   , this is because even though all of the GC-rich galaxies have on average high specific frequencies, GC-poor UDGs show a mix of specific frequencies, capable of confusing the clustering algorithm and strongly changing the identified classes.Whether the separation with   is more correct than the one without it needs to be further assessed in a larger sample of galaxies.In this study, we choose to use the results that are stable in the absence of   .As a last comment, we note that if we exclude the GC number from the analysis, we can apply KMeans to the full sample of 59 UDGs (as the GC number is the only property that is not available for all galaxies).Once this is done, the exact same classes are identified, and the variance on the properties is smaller as the classes are comprised of a larger sample of galaxies.This gives us confidence that our clustering methodology is capturing the main properties of our dataset even when constrained to a smaller subsample.
The properties of the two classes of UDGs identified with KMeans as summarised in Table 1 are: • Class A. UDGs in this class have lower stellar masses, prolonged SFHs (i.e., higher star formation timescales), are more elongated (lower /), live in less dense environments (lower log  10 ), host fewer GCs, are younger (smaller   ), have higher  dwarf MZR (i.e., closer to zero and consistent with the classical dwarf MZR), and are fainter (higher  ,0 ).
• Class B. UDGs in this class have higher stellar masses, rapid SFHs (i.e., smaller ), are rounder (high /), live in the densest of our probed environments (higher log  10 ), some host the biggest GC systems, are older (higher   ), have smaller  dwarf MZR (i.e., negative values, being thus more metal-poor than the classical dwarf MZR), and are brighter (smaller  ,0 ).
The centres of the two identified classes and the number of UDGs associated with each class are shown in the upper block of Table 1 for the classification including   .The lower block, on the other hand, shows the centres and number of UDGs following the classification applied to all 59 galaxies (i.e., removing  GC from the analysis).The quoted uncertainties are the standard deviation of each analysed property considering the galaxies that belong to that class, i.e., an approximation of how far apart objects belonging to the same class can be.The cluster centres and uncertainties can also be seen in the histogram distribution of the properties in Fig. 10.Properties that provide the best separation between classes will have more clear bimodal distributions, while the ones that are not strong discriminators will have similar Gaussian distributions.The classes each UDG in our sample are associated in the analyses with and without  GC are given in Table B2.It is reinforcing to see that 37 out of the 38 UDGs included in both analyses (i.e., with and without  GC ) have preserved their classes.Only 1 UDG (MATLAS-984) changed classes when the number of GCs was excluded.We note that changing classes is expected if the classification is being driven by the property that was removed (i.e.,  GC ) and the other properties are in overlapping cluster regions.This is the case for MATLAS-984, which has properties that could fit into any class when  GC is not considered (see Tables B2 and B3).
We emphasise that these classes were obtained using a small sample of UDGs, all in low-density environments.Thus, the centres found may not be representative of a vast population of UDGs with strongly divergent properties.These cluster centres should perhaps be interpreted as a measurement of which properties are important discriminators of UDG formation scenarios and which ones are not.Although there is not a net separation or threshold in the properties of UDGs formed via the two scenarios, and possibly not a one-to-one correspondence between a given class and a formation scenario, we recognise that UDGs in Class A have properties (see Table 1) that are more aligned with being formed through puffed-up dwarf scenarios, i.e., properties similar to those of classical dwarf galaxies.On the other hand, UDGs in Class B share similar properties to those of failed galaxy scenarios (see Table 1), i.e., properties consistent with early-quenching.We note that tidal UDGs, due to low number statistics, was not found as a third class in our dataset, although these types of UDGs have overall different properties than puffed-up dwarfs or failed galaxies.This class is expected to become clearer in larger samples of UDGs, including more massive and cluster ones, where tidal effects are more present.
This clustering algorithm was applied to a small set of galaxies and needs to be further tested, both with a larger number of UDGs and in comparison with classical dwarfs.This way, a more comprehensive comparison between subpopulations can be made.We note, nonetheless, that this separation (in terms of GC number, axis ratio, and mass) has also been reported in the Coma cluster (Lim et al. 2018), so these two classes of UDGs may be a general separation regardless of their environments.A study combining and comparing the MATLAS UDGs, the UDGs studied in B22 and a control sample of classical MATLAS dwarfs is under preparation (Buzzo et al. in prep.) and will focus on deepening these clustering findings.Other properties can also be added to the clustering in the future to enrich the classification, namely [/Fe] (e.g., Ferré-Mateu et al. 2018;Villaume et al. 2022, FM23), dynamical masses and velocity dispersions (e.g., Gannon et al. 2020Gannon et al. , 2021Gannon et al. , 2022Gannon et al. , 2023)), positioning in the phase-space diagram (e.g., Gannon et al. 2022;Forbes et al. 2023, FM23), etc.As a final remark, once the classes are further tested and if they are proven to be well founded, they could be used to predict unknown properties of UDGs.For example, by knowing the stellar populations and morphology of a given UDG, using these clustering results, one would be able to predict if it is GC-poor or GC-rich, in cases where knowing the GC-richness of the galaxy would otherwise be impossible (e.g., at distance higher > 100 Mpc), further showing the importance of such findings.

CONCLUSIONS
In the current study, we have used the fully Bayesian Monte Carlo Markov Chain inference code Prospector to perform spectral energy distribution fitting on the 59 MATLAS UDGs using data from the optical to the mid-IR.We use Prospector to recover stellar populations using two different configurations, one with the redshift as a free parameter, and one with the redshift fixed to the redshift of the group where the UDG was identified.
We find that the models with the redshift fixed to that of the host group are the most appropriate to study these galaxies.This assumption is based on the results of Heesters et al. (2023), who found that 75% of the MATLAS dwarfs studied by them are at the same distance as their massive hosts.This is supported by our study with Keck/DEIMOS of three MATLAS UDGs, also found to be at the same redshift as their hosts.
Our results indicate that the MATLAS UDGs have intermediateto-old ages, with an average mass-weighted age of 7.1 ± 1.8 Gyr.They are also metal-poor with an average [M/H] of −1.2 ± 0.2 dex.They have an average star formation timescale  of 1.6 ± 0.7 Gyr and are consistent with no dust attenuation, displaying an average of   = 0.12 ± 0.07 mag.
When studying the MATLAS UDGs in the mass-metallicity plane, we find that their distribution is highly bimodal, with the first mode being well explained by the classical dwarf MZR, while the second one is more consistent with the evolving MZR at a high-redshift.
To further investigate the cause of this bimodality, we study the positioning of the UDGs in the mass-metallicity plane according to their surface brightness, effective radii, axis ratios, local volume densities, mass-weighted ages, star formation timescales, number of GCs and GC specific frequencies.We find that UDGs split into two main classes: Class A: They have lower stellar masses, prolonged SFHs, are more elongated, live in less dense environments, host fewer GCs, are younger, are consistent with the classical dwarf MZR, and are fainter.
Class B: UDGs in this class have higher stellar masses, rapid SFHs, are rounder, live in the densest of our probed environments, host on average the biggest GC systems, are older, lie below the classical dwarf MZR, and are brighter.These galaxies are overall better explained by the evolving MZR at high-redshift, i.e., consistent with early-quenching.
This overall picture seems to indicate that UDGs belonging to Class A, combining all of the aforementioned properties, are better explained by puffed-up dwarf formation scenarios.Alternatively, UDGs of Class B, taking into consideration all of the morphological and stellar population properties found, seem to be more aligned with failed galaxy formation scenarios, in agreement with the spectroscopic studies.Further testing of this clustering technique, including a larger sample of UDGs and classical dwarfs for comparison is under preparation (Buzzo et al. in prep.).
This paper provides a continuation of the photometric study of the stellar populations of UDGs across the sky.We demonstrate that SED fitting techniques, coupled with a broad wavelength coverage, are an important approach to statistically understand the origins of UDGs and their impact on the field of galaxy formation and evolution.We fit each spectrum with pPXF (Cappellari 2017), using the MILES stellar library (Vazdekis et al. 2015), assuming a Kroupa initial mass function (Kroupa 2001) and the BaSTI isochrones.In Fig. A1, we show one example of pPXF fit for the galaxy MATLAS-342.Results from the fits of all three galaxies are summarised in Table A1.The results indicate that all three probed UDGs are at the same redshifts as their massive neighbours.
We note that MATLAS-368 was selected to be followed up because it was thought to be nucleated, which would facilitate recovering its velocity because of its brightness.However, the central point source was found to be at  = 0, thus likely being a foreground star.As an alternative, we studied the galaxy light coming from the outer parts of the slit, but this rendered low S/N data.The identification of the H line was possible, resulting in the recessional velocity shown in Table A1.We caution, nonetheless, that the result is based on a single line, which is highly uncertain.
The spectrum of the galaxy that was found to be closer than its host, MATLAS-1059, has very low signal-to-noise (3.8 Å −1 ) due to its low exposure time, resulting in a much more uncertain velocity than the other studied UDGs.Even in this extreme case of low S/N, the recovered velocity is consistent with the one of the host.Additionally, from CFHT and DECaLS imaging, we can see that this UDG seems tidally disrupted by the host, further confirming its proximity to it.These findings show that all three UDGs studied are at the redshifts of their hosts.This is in agreement with the recent study of Heesters et al. (2023) who found, after studying 56 MATLAS dwarfs, that 75% of them were indeed at the same redshifts as their hosts, including the three UDGs present in their sample.Putting together our results and the ones from Heesters et al. (2023), we have a cumulative result that 6 out of 6 MATLAS UDGs followed up with spectroscopy are confirmed as associated with their parent host galaxy.We thus conclude that assuming the redshift of the host for the UDGs that have not been spectroscopically followed-up is a reasonable assumption, bearing in mind that these are based on ensemble properties and may not hold for every individual galaxy.

APPENDIX C: THE EFFECTS OF UNDERESTIMATING UNCERTAINTIES
Galaxy fitting codes are known for severely underestimating the uncertainties of their recovered parameters.In this Appendix, we show the effects on the stellar populations properties obtained with Prospector if the photometric uncertainties from GALFITM were underestimated by 10%, 20%, 50% and 100%.
In Fig. C1, we compare the stellar populations properties of MATLAS-2019 using the nominal uncertainties coming from GALFITM and if they were underestimated.We can see that if the uncertainties were underestimated by 10 or 20%, the recovered stellar populations barely change.However, if they were severely underestimated, such as by 50 or 100%, the changes start to get important.We emphasise, however, that we do not expect the uncertainties to be this highly underestimated.The biggest of these changes is in the mass weighted ages reaching a change of 5 Gyr.The biggest change in metallicity is 0.1 dex, 0.2 dex in mass and 0.02 in   .The minor changes in metallicity and mass give us confidence that even if the uncertainties are underestimated, the major conclusions in this paper (i.e., the bimodality in the mass-metallicity plane) would mantain.

APPENDIX D: ALTERNATIVE VIEW OF FIG. 7
In Fig. 7 we show the MATLAS UDGs overlaid in the massmetallcity plane in an attempt to compare their positioning to known MZRs and to start inferring some of their likely formation histories.We colour-code the MATLAS UDGs by several key parameters and show that these galaxies distill into two main types: some that are consistent with the classical dwarf MZR and some that are consistent with the high-redshift MZR.As the colour gradients may be hard to visualise in some cases, in this appendix we provide an alternative view of Fig. 7, shown in Fig. D1.In this case, the key properties are directly plotted against their distance to both the dwarf MZR and the high-redshift MZR.

Figure 1 .
Figure 1.GALFITM model of MATLAS-2019 (NGC 5846_UDG1).The first column shows the input image.Second column is the mask used in the model.The third column is the model.The fourth column is the residual (input -model).Rows show the modelling in the , ,  and  bands, respectively.

Figure 2 .
Figure 2. Prospector best-fit SED and covariance matrix for MATLAS-2019.In all panels, purple colours stand for the models with the redshift free ( ≠  host ) and green colours are the models with the redshift fixed to the host ( =  host ).Top-right Best-fit spectra of MATLAS-2019 in both configurations.Red points are the observed magnitudes of the galaxy in the g, r, i, z, W1, W2, W3 and W4 bands, respectively.Yellow curves show the DECaLS transmission curves.Blue curves are the WISE transmission curves.The remaining panels show the Monte Carlo Markov Chain (MCMC) corner plot comparing the posterior distribution of the two configurations.The first panel in each column shows the 1D posterior distribution of the fitted parameter, while the remaining panelsshow the correlation between parameters.This image can be read and interpreted as a covariance matrix.Columns stand for redshift, stellar mass, metallicity, star formation timescale, mass-weighted age and interstellar diffuse dust extinction.Metallicity, star formation timescale, age and dust attenuation are similar in both configurations.However, the stellar mass is strongly affected by the freedom of the redshift, as these two parameters are highly degenerate.

Figure 3 .
Figure3.Age and metallicity of the 59 MATLAS UDGs recovered with Prospector using the two different configurations discussed in Section 3.2.In all panels, purple colours stand for the models with the redshift free ( ≠  host ) and green colours are the models with the redshift fixed to the host ( =  host ).The bottom-right panel shows the age-metallicity distribution of the 59 MATLAS UDGs in both configurations.Histograms at the top and on the left show the marginal distributions of the data displayed on the age-metallicity plane.We conclude that models with the redshift free (purple) recover on average younger ages and more metal-rich populations for the UDGs.

Figure 4 .
Figure 4. Stellar mass, star formation timescale and dust attenuation of the 59 MATLAS UDGs recovered with Prospector using the two different configurations discussed in Section 3.2.In all panels, purple colours stand for the models with the redshift free and green colours are the models with the redshift fixed to the host.We conclude that models with the redshift free (purple) result in slightly more massive UDGs and a much narrower distribution centred at log(  ★ / ⊙ ) ∼ 7.8 than models with the redshift fixed.No significant difference is observed in the star formation timescales and dust attenuation recovered in both configurations.

Figure 5 .
Figure 5. Accuracy in recovering photometric redshifts as a function of redshift.  measures the difference between the spectroscopic redshifts of the candidate hosts and the recovered photometric redshifts of the UDGs.Marginal histogram shows the distribution of   .Positive values mean that UDGs were found to be in the foreground of their hosts.Negative values mean that UDGs were found in the background of their hosts.From previous findings (Heesters et al. 2023), the expectation is that at least 75% of the dwarf galaxies in the MATLAS survey are at the same redshifts as their hosts.Our results show a much larger portion of the UDGs at different redshifts (> 70% of UDGs with |   |> 0.0067, equivalent to a difference in velocity of |   |> 200 km s −1 ), indicating that photometric redshifts are not reliable at non-cosmological distances.
2019 has been spectroscopically studied by Müller et al. (2020b), Heesters et al. (2023) and FM23.Müller et al. (2020b) found, using VLT/MUSE data, that MATLAS-2019 has a metallicity of [M/H] = [Fe/H] = −1.33 +0.19 −0.01 , showing that they behave in a similar way to their massive counterparts in terms of the environments they inhabit.Noteworthy, Poulain et al. (2021) also found a morphology-density relation for the dwarfs on the MATLAS fields comparable to that seen in dense clusters, and Li et al. (2023) found a similar behaviour for UDGs surrounding Milky Way analogs.

Figure 6 .
Figure 6.Stellar mass-metallicity distribution of UDGs.Top panel: MAT-LAS UDGs are shown in red.Field UDGs from Barbosa et al. (2020) are shown in blue.UDGs in different environments from the spectroscopic sample of FM23 are shown in green.UDGs in a variety of environments from B22 are shown in magenta.The Simon (2019) MZR for classical Local Group dwarf galaxies is shown with the black dashed line.The dash-dotted line is the evolving MZR at redshift  = 2.2 from Ma et al. (2016).Middle panel: Difference in metallicity between the MATLAS UDGs and the Simon (2019) classical dwarf MZR.UDGs at  MZR = 0 follow perfectly the MZR, if  MZR > 0, then the UDG lies above the MZR, if  MZR < 0, then the UDG lies below the MZR.Bottom panel: Difference in metallicity between the MATLAS UDGs and the Ma et al. (2016) MZR at  = 2.2.MATLAS UDGs show a bimodality in the mass-metallicity distribution, with the first mode being consistent with the classical dwarf MZR and the second being more consistent with the high-redshift MZR.

Figure 7 .
Figure 7. Distribution of the physical properties of the MATLAS UDGs in the stellar mass-metallicity plane.Plus signs show the results obtained from SED fitting with PROSPECTOR for the 59 UDGs in our sample, with the 2 confirmed tidal UDGs in our sample being highlighted with star symbols (MATLAS-1824 and MATLAS-478).Top-row: UDGs are colour-coded by surface brightness, effective radii, axis ratio and local environment (Marleau et al. 2021), respectively.Bottom-row: UDGs are colour-coded by their mass-weighted ages, star formation timescales, number of GCs (Marleau et al. subm.) and GC specific frequencies (per unit stellar mass), respectively.The Simon (2019) MZR for classical Local Group dwarf galaxies is shown with the dashed black line.The dash-dotted black line is the evolving MZR at redshift  = 2.2 from Ma et al. (2016).An alternative view of this figure is shown in Appendix D1.These results indicate a possible bimodality in the MATLAS UDGs population: Population i) UDGs are fainter, smaller, elongated, live in less dense environments, are younger, have prolonged SFHs, have a smaller number of GCs and on average smaller specific frequencies.Population ii) UDGs in this population are brighter, bigger, rounder, live in denser environments, are older, have shorter SFHs, host on average more GCs and have on average larger specific frequencies.

Figure 8 .
Figure 8. Difference in metallicity (  MZR ) between known MZRs and the MATLAS UDGs colour-coded by their number of GCs.Description is the same as Fig. 6.GC-poor MATLAS UDGs follow the classical dwarf MZR from Simon (2019).On another hand, all UDGs that follow the Ma et al. (2016) MZR at high-redshift are GC-rich.

Figure 10 .
Figure 10.The two groups of UDGs identified using the KMeans clustering algorithm.Distribution of UDGs in the parameter space defined by stellar mass (log(  ★ / ⊙ )), star formation timescale (), axis ratio (/), local volume density (log  10 ), number of GCs ( GC ), mass-weighted age (  ),  dwarf MZR (i.e., level of consistency between the UDGs and the classical dwarf MZR), and central surface brightness ( ,0 ).Green points indicate galaxies that belong to class A (associated to the "puffed-up dwarf" scenario), pink points indicate galaxies that belong to class B (associated to the "failed galaxy" scenario).The histograms are smoothed using a Kernel Density Estimation (KDE) technique and indicate the distribution of a given parameter; green and pink histograms indicate those of class A and B, respectively.

Figure A1 .
Figure A1.pPXF fitting of MATLAS-342, one MATLAS UDG with available spectra from Keck/DEIMOS.Black curve shows the spectrum of each source, while red curve shows the best-fit.Green dots are the residuals from the fit.Grey regions show the masked areas in the fit.Legend box within the plot shows the recovered radial velocity of the UDG ( = 2142 ± 31 km/s) and the radial velocity of its host galaxy ( = 2153 km/s).Results indicate that this UDG is at the same redshift as its host.
GALFITM results and physical properties of MATLAS UDGs B3 Stellar population properties from Prospector

Figure C1 .
Figure C1.Differences in the recovered stellar populations of MATLAS-2019 if the photometric uncertainty in the optical obtained from GALFITM was underestimated by 10% (red), 20% (blue), 50% (green) and 100% (orange).Left: Differences in metallicity and mass-weighted ages obtained.Right: Differences in stellar mass and dust attenuation.We conclude that the photometry being underestimated by 50% or 100% can severely change the ages of the galaxies, reaching a change of 5 Gyr.The maximum change in metallicity is 0.11 dex, in stellar mass is 0.2 dex and in   is 0.02.Due to the small changes in mass and metallicity, the main conclusions of the paper regarding the bimodality in the mass-metallicity plane remain even if the uncertainties are found to be extremely underestimated.

Figure D1 .
Figure D1.Alternative view of Fig. 7. Distribution of the physical properties of the MATLAS UDGs in the stellar mass-metallicity plane.Plus signs show the results obtained from SED fitting with PROSPECTOR for the 59 UDGs in our sample.The Simon (2019) MZR for classical Local Group dwarf galaxies is shown with the dashed black line.The dash-dotted line is the evolving MZR at  = 2.2 from Ma et al. (2016).In all panels, the properties are plotted against distance to the dwarf (left-hand side y-axis) and high-redshift MZRs (right-hand side y-axis).Top-row: Surface brightness (left) and effective radii (right).Second row: Axis ratio and local environment(Marleau et al. 2021), respectively.Third row: Mass-weighted ages and star formation timescales.Fourth row: Number of GCs(Marleau et al. subm.)  and GC specific frequencies (per unit stellar mass), respectively.Results are as described in the caption of Fig.7 and Heesters et al. (2023)(S/N = 12.4 pix −1 ).Gyr and a stellar mass of log( ★ / ⊙ ) = 8.01 +0.03 −0.05 for MATLAS-2019, for a distance of 25.2 Mpc (i.e., the distance to the host NGC 5846).Our study is consistent within the uncertainties with both Müller et al. (2020b) and FM23.The results, nonetheless, in terms of age and metallicity, are closer to the ones found by

Table 1 .
Mean values of UDG classes obtained with the KMeans clustering results.The first block is the classification done using  GC (only 38 galaxies).The second block is the classification without  GC (all 59 MATLAS UDGs).

Table B2 .
GALFITM morphological parameters and physical properties of the MATLAS UDGs.