Parameter Estimation for Scarce Stellar Populations

TL;DR

This paper introduces a Bayesian method for estimating stellar population parameters from color magnitude diagrams, effective even for scarce data, and applies it to study star clusters and ultra faint dwarf galaxies.

Contribution

The paper presents a novel Bayesian approach that accurately estimates multiple stellar parameters without data binning, applicable to small stellar populations and complex star formation histories.

Findings

Successfully applied to Magellanic Cloud clusters and ultra faint dwarf galaxies.

Recovered parameters consistent with literature, e.g., age ~13.7 Gyr, low metallicity.

Detected possible multiple populations in Ursa Major I.

Abstract

We present a Bayesian method to determine simultaneously the age, metallicity, distance modulus, and interstellar reddening by dust of any resolved stellar population, by comparing the observed and synthetic color magnitude diagrams on a star by star basis, with no need to bin the data into a carefully selected magnitude grid. We test the method with mock stellar populations, and show that it works correctly even for scarce stellar populations with only one or two hundred stars above the main sequence turn off. If the population is the result of two star formation bursts, we can infer the contribution of each event to the total stellar population. The code works automatically and has already been used to study massive amounts of Magellanic clouds photometric data. In this paper we analyze in detail three Large Magellanic Cloud star clusters and 6 Ultra Faint Dwarf Galaxies. For these galaxies we recover physical parameters in agreement with those quoted in the literature, age $\sim13.7$ Gyr and a very low metallicity $log\,Z\sim-4$. Searching for multiple populations in these galaxies, we find, at a very low significance level, signs of a double stellar population for Ursa Major I: a dominant old population and a younger one which contributes $\sim25$% of the stars, in agreement with independent results from other authors.

Figures

Tables

Table 1. Table 1: Mock galaxy with a single stellar population

Parameter	Input value	Marg Estimate
$t$ [Gyr]	$11.125$	${11.25}_{-2.0}^{+2.75}$
$Log\left(Z\right)$	-4	$-4_{-0.3}^{+0.4}$
$\left(m-M\right)$	21.9	${21.9}_{-0.3}^{+0.4}$
$E\left(V-I\right)$	0.085	$0.075\pm 0.05$

Table 2. Table 2: Mock galaxy with 2 old stellar populations

Parameter	Input value	Marg Estimate
$t_y$ [Gyr]	$8.125$	${8.0}_{-3.2}^{+2}$
$t_o$ [Gyr]	$12.125$	${13.25}_{-2.5}^{+3.8}$
$LogZ$	-4.0	$-4\pm 0.5$
$w_y$	0.5	$0.5\pm 0.35$
$\left(m-M\right)$	21.9	$21.9\pm 0.3$
$E\left(V-I\right)$	0.085	$0.07\pm 0.05$

Table 3. Table 3: Estimates from ℒ N m ​ a ​ r ​ g superscript subscript ℒ 𝑁 𝑚 𝑎 𝑟 𝑔 \mathcal{L}_{N}^{marg} distributions for a mock galaxy with 2 old stellar populations of ages t y = 8.125 subscript 𝑡 𝑦 8.125 t_{y}=8.125 Gyr, t o = 12.125 subscript 𝑡 𝑜 12.125 t_{o}=12.125 Gyr, L ​ o ​ g ​ Z = − 4 𝐿 𝑜 𝑔 𝑍 4 Log\,Z=-4 , ( m − M ) = 21.9 𝑚 𝑀 21.9 (m-M)=21.9 , E ​ ( V − I ) = 0.085 𝐸 𝑉 𝐼 0.085 E(V-I)=0.085 , and varying w y subscript 𝑤 𝑦 w_{y} .

$w_{y,input}$	$w_{y,est}$	$t_{y,est}$ [Gyr]	$t_{o,est}$ [Gyr]	$Log{Z}_{est}$	${\left(m-M\right)}_{est}$	$E{\left(V-I\right)}_{est}$
0.1	${0.25}_{-0.25}^{+0.4}$	${9.25}_{-4.3}^{+3.7}$	${13.75}_{-2.4}^{+4.0}$	$-4_{-0.3}^{+0.5}$	21.9$\pm$0.3	${0.07}_{-0.5}^{+0.7}$
0.3	${0.4}_{-0.3}^{+0.35}$	${8.25}_{-2.5}^{+2.1}$	${13.5}_{-2.0}^{+4.3}$	$-4_{-0.3}^{+0.5}$	21.9$\pm$0.3	$0.07\pm 0.05$
0.5	$0.5\pm 0.35$	${8.0}_{-3.2}^{+2}$	${13.25}_{-2.8}^{+4.5}$	$-4\pm 0.5$	$21.9\pm 0.3$	$0.07\pm 0.05$
0.7	${0.65}_{-0.35}^{+0.3}$	${8.0}_{-2.0}^{+1.5}$	${13.25}_{-2.8}^{+4.5}$	$-4_{-0.3}^{+0.5}$	21.9$\pm$0.25	${0.07}_{-0.6}^{+0.7}$
0.9	${0.75}_{-0.45}^{+0.25}$	${7.75}_{-2.0}^{+1.5}$	${12.0}_{-4.0}^{+5.5}$	$-4\pm 0.5$	21.9$\pm$0.25	0.07$\pm 0.07$

Table 4. Table 4: Estimates from ℒ N m ​ a ​ r ​ g superscript subscript ℒ 𝑁 𝑚 𝑎 𝑟 𝑔 \mathcal{L}_{N}^{marg} distributions for a mock galaxy with 2 young stellar populations of L ​ o ​ g ​ Z = − 2.1 𝐿 𝑜 𝑔 𝑍 2.1 Log\,Z=-2.1 , with t y = 4.05 subscript 𝑡 𝑦 4.05 t_{y}=4.05 Myr, t o = 6.12 subscript 𝑡 𝑜 6.12 t_{o}=6.12 Myr, ( m − M ) = 29.0 𝑚 𝑀 29.0 (m-M)=29.0 , E ​ ( V − I ) = 0.15 𝐸 𝑉 𝐼 0.15 E(V-I)=0.15 , and varying w y subscript 𝑤 𝑦 w_{y} .

$w_{y,input}$	$w_{y,est}$	$t_{y,est}$ [Myr]	$t_{o,est}$ [Myr]	${\left(m-M\right)}_{est}$	$E{\left(V-I\right)}_{est}$
0.1	-	-	${6.0}_{-0.3}^{+1.8}$	29.0$\pm$0.3	0.15$\pm$0.025
0.3	${0.4}_{-0.3}^{+0.4}$	${4.2}_{-0.4}^{+0.9}$	${6.8}_{-1.1}^{+1.0}$	28.9$\pm$0.3	0.145$\pm$0.025
0.5	${0.5}_{-0.35}^{+0.35}$	${4.0}_{-0.2}^{+1.0}$	6.8$\pm$1.1	${29.0}_{-0.4}^{+0.2}$	0.145$\pm$0.025
0.7	${0.6}_{-0.4}^{+0.3}$	${4.2}_{-0.4}^{+0.8}$	${7.0}_{-1.3}^{+0.9}$	${29.1}_{-0.5}^{+0.2}$	0.145$\pm$0.025
0.9	${0.7}_{-0.4}^{+0.3}$	${4.0}_{-0.2}^{+0.9}$	$6.5\pm 1.3$	${29.1}_{-0.4}^{+0.3}$	0.14$\pm$0.025

Table 5. Table 5: Physical parameters of the UFDGs.

UFDG	Work	$t$ (Gyr)	$log\left(Z\right)$	$\left(m-M\right)$	$E\left(B-V\right)$
Bootes I	Present	13.25-2.25	$-{3.7}_{-0.5}^{+0.4}$	19.0$\pm 0.35$	0.05$\pm 0.05$
	Brown et al. (2014)	13.3$\pm$0.3$\pm \sim 1$	$-4.1\pm 0.5$	19.11$\pm 0.07\pm 0.01$	0.04$\pm 0.01$
	McConnachie (2012)	-	$-4.32\pm 0.11$	19.11$\pm 0.08$	0.017
Canes Venatici II	Present	13.75-3.25	$-4_{-0.3}^{+0.7}$	20.9$\pm 0.5$	0.05${}_{}{}^{+0.06}{}_{-0.04}{}^{}$
	Brown et al. (2014)	13.6$\pm 0.3\pm \sim 1$	$-4.7\pm 0.6$	21.04$\pm 0.06\pm 0.01$	0.04$\pm 0.01$
	Weisz et al. (2014)	10.0${}_{-1.1}{}^{}{}_{}{}^{+1.5}$	-	-	-
	McConnachie (2012)	-	$-4.0\pm 0.05$	21.02$\pm 0.06$	0.010
Coma Berenices	Present	14.0-3.25	$-{3.7}_{-0.6}^{+0.4}$	18.05$\pm 0.35$	0.04${}_{}{}^{+0.05}{}_{-0.03}{}^{}$
	Brown et al. (2014)	13.9$\pm 0.3\pm \sim 1$	$-4.3\pm 0.5$	18.08$\pm 0.10\pm 0.01$	0.04$\pm 0.01$
	McConnachie (2012)	-	$-4.4\pm 0.05$	18.2$\pm 0.20$	0.017
Hercules	Present	13.75-3.25	$-4_{-0.3}^{+0.7}$	20.5$\pm 0.5$	0.09${}_{}{}^{+0.7}{}_{-0.5}{}^{}$
	Brown et al. (2014)	13.1$\pm 0.3\pm \sim 1$	$-4.4\pm 0.4$	20.92$\pm 0.05\pm 0.01$	0.09$\pm 0.01$
	Weisz et al. (2014)	13.0${}_{-2.8}{}^{}{}_{}{}^{+0.3}$	-	-	-
	McConnachie (2012)	-	$-4.2\pm 0.04$	20.6$\pm 0.2$	0.063
	Musella et al. (2012)	-	$-4.1\pm 0.15$	20.6$\pm 0.2$	0.09$\pm 0.02$
Leo IV	Present	13.75-3.0	$-4_{-0.2}^{+0.4}$	20.9$\pm 0.35$	0.07$\pm 0.02$
	Brown et al. (2014)	13.1$\pm 0.4\pm \sim 1$	$-4.3\pm 0.6$	21.12$\pm 0.07\pm 0.01$	0.08$\pm 0.01$
	Weisz et al. (2014)	11.7${}_{-4.2}{}^{}{}_{}{}^{+1.4}$	-	-	-
	McConnachie (2012)	-	$-4.3\pm 0.07$	20.94$\pm 0.09$	0.026
	Moretti et al. (2009)	-	$-4.1\pm 0.1$	20.94$\pm 0.07$	0.04$\pm 0.01$
Ursa Major I	Present	14.0-2.5	$-{3.7}_{-0.5}^{+0.5}$	19.8$\pm 0.5$	0.05$\pm 0.05$
	Brown et al. (2014)	12.7$\pm 0.3\pm \sim 1$	$-4.4\pm 0.5$	20.10$\pm 0.05\pm 0.01$	0.05$\pm 0.01$
	McConnachie (2012)	-	$-4.0\pm 0.04$	19.93$\pm 0.10$	0.020

Equations

$$ith-\mathrm{isochrone}=h_i[t,Z,m-M,A_{\lambda },E({\lambda }_1-{\lambda }_2)].$$

$$p(h_i|o)=\frac{p(o|h_i)p(h_i)}{p(o)},$$

$$\mathcal{L}(h_i|o)=\frac{p(o|h_i)p(h_i)}{p(o|h_{max})p(h_{max})}.$$

$$\mathcal{L}(h_i|o)=\frac{p(o|h_i)}{p(o|h_{max})}.$$

$$p(o|h_i)=p(s_k|h_i)=\int h_i(x,y)U_k(x-x_k,y-y_k)dxdy,$$

$$U_k(x-x_k,y-y_k)=C_{k,x,y}{e}^{-\left(\frac{(x-x_k{)}^2}{2{\sigma }_{x,k}^2}+\frac{(y-y_k{)}^2}{2{\sigma }_{y,k}^2}\right)},$$

$$C_{k,x,y}=\frac{1}{2\pi {\sigma }_{x,k}{\sigma }_{y,k}}.$$

$$p(s_k|h_i)=\frac{C_{k,x,y}}{n}\sum _{j=1}^n\exp \left(-\frac{(x_j-x_k{)}^2}{2{\sigma }_{x,k}^2}-\frac{(y_j-y_k{)}^2}{2{\sigma }_{y,k}^2}\right),$$

$${\mathcal{L}}_N(h_i|o)=\prod _{k=1}^{N}p(s_k|h_i{)}^{1/N},$$

$$\ln {\mathcal{L}}_N(h_i|o)=\frac{1}{N}\sum _{k=1}^N\ln p(s_k|h_i).$$

$${\mathcal{L}}_N^{marg}({\varphi }_q)=\int \mathcal{L}({\varphi }_q,\nu )d\nu ,$$

$$h_{i,j,w}=w{h}_i(x,y)+(1-w)h_j(x,y),$$

$${\sigma }_x$$

$${\sigma }_y$$

$${\sigma }_V$$

$${\sigma }_I$$

$${\sigma }_{V-I}$$

$$p_{memb}=1-\frac{N_{\ast ,field}}{N_{\ast ,cluster}}\cdot \frac{A_{cluster}}{A_{field}},$$

$$m_{\lambda }=-2.5\times log{f}_{\lambda }-21.1,$$

Full text

††pagerange: Parameter Estimation for Scarce Stellar Populations–Parameter Estimation for Scarce Stellar Populations††pubyear: 2019

1 Introduction

The star formation history (SFH) of resolved stellar populations, like star clusters or dwarf galaxies, is usually derived from their colour-magnitude diagram (CMD). The number of stars in characteristic phases of stellar evolution in the CMD carry information on the age, the metallicity, and the strength of the star formation burst. For example, the position of the Main Sequence Turn Off (MSTO) is sensitive to age, the presence of bright and blue mains sequence (MS) stars is indicative of a young star burst, Horizontal Branch (HB) and RR-Lyrae stars are characteristic of low-metallicity old stellar populations. By fitting synthetic or mock CMDs computed for a wide range of physical parameters to an observed CMD, we can estimate the best age ($t$), metallicity ($Z$), distance modulus $(m-M)$, colour excess $E(\lambda_{1}-\lambda_{2})$, and stellar mass that describe a given stellar population, according to a particular set of stellar evolution models, within the observational errors. Several methods have been developed to find the best match between observed CMDs and theoretical isochrones, ranging from the simple fit by eye, to more refined statistical techniques. Fits by eye are subjective, do not provide confidence intervals for the estimated parameters, and become extremely laborious when the task implies numerous stellar systems.111 For a review on the derivation of SFHs from CMDs see Gallart et al. (2005).

Fits can also be performed on the number of observed and expected stars inside conveniently chosen colour and magnitude bins in the CMD. The observed and theoretical CMDs, which must include modeling of the photometric errors, are binned in identical fashion and the number of stars in corresponding bins are compared. Minimizing $\chi^{2}$ is not the most convenient parameter estimator for poorly populated CMDs, since to obtain meaningful results, the number of stars in each bin must be statistically significant ($N\geq 10$), a task feasible in richly populated CMDs. Binning methods that overcome this limitation have been proposed, e.g., Mighell (1999), Dolphin (2002), and Aparicio & Hidalgo (2009).

In this paper we develop a Bayesian inference code to estimate simultaneously, from unbinned CMDs, the set of parameters $t,Z,(m-M)$, $E(\lambda_{1}-\lambda_{2})$, and the SFH characterizing single or double stellar populations. A similar technique was used by Tolstoy (1996) and Tolstoy & Saha (1996) to determine the SFH of single stellar populations, and by Hernandez & Valls-Gabaud (2008) and Perren et al. (2015) to determine $t,Z,(m-M)$, and $E(\lambda_{1}-\lambda_{2})$ for the same kind of populations. Our approach differs from theirs in the way that we build the likelihood distribution function for all the parameters under consideration, and in the use of the marginalized probability distribution functions to assess the confidence intervals for the validity of our estimates. More importantly, we show that our method can be used to estimate the parameters of double stellar populations, including the stellar contribution of each burst, and examine the behavior of our results varying the sample size, its limiting magnitude, and the photometric errors in the CMD. Our approach can be easily extended to study multiple ($>2$) populations.

A first aim of this project is to estimate the parameters of a large number of young star clusters in the Magellanic Clouds in an objective and automated manner. In this paper we study three LMC clusters in detail. In Bitsakis et al. (2017, 2018) we used our parameter inference tool to date systematically 4850 clusters in the Largel Magellanic Cloud (LMC) and 1319 clusters in the Small Magellanic Cloud (SMC), respectively. These clusters range in age from 10 to a few hundred Myr and are sparsely populated, containing typically $\sim 100$ stars. Many of them have not been catalogued previously. In this paper we study three LMC clusters in detail, two of them handpicked by inspecting their CMDs and selecting those showing features that denote the presence of double stellar populations. Applying our method, we are able to distinguish the presence of two stellar populations, establish their ages, and the stellar contribution of each burst.

A second goal of this paper is to revisit the derivation of the SFH and other physical parameters of a sample of six Ultra Faint Dwarf Galaxies (UFDGs) observed with the Hubble Space Telescope (HST) by Brown et al. (2012); Brown et al. (2014). UFDGs are interesting because they are thought to be fossils from the first star bursts in the universe. They are indeed old, with ages $\sim 13$ Gyr, show a very low metal content, $Z\sim 0.0001$ (Kirby et al., 2008; Brown et al., 2014), and their SFHs seem to imply that a synchronized global event, such as reionization, deprived them of gas and quenched star formation (Ricotti & Gnedin, 2005). The number of resolved stars in a UFDG is typically $\sim$1000, with only $\sim$ 100 of them above the MSTO. Different methods have been used to estimate their SFHs, from isochrone fitting by eye (Sand et al., 2012) to adequate CMD binning methods (Brown et al., 2014; Weisz et al., 2014). Brown et al. (2014) adopt the maximum likelihoood binning method of Dolphin (2002) (reviewed by Walmswell et al. (2013)), using a binning grid from the MSTO to the top of the Sub Giant Branch (SGB), excluding zones with few stars. As far as we are aware, our analysis of the UFDGs is the first one to use most of the observed stars in the CMD, on a star by star basis without binning the data, exploring the most sensitive areas of the CMD to reduce the confidence intervals of the estimated parameters.

In §2 we describe in detail our implementation of the Bayesian inference approach. In §3 we apply our method to study the CMDs of single and double burst mock populations, recovering their SFHs and physical properties, and measure the capabilities and limitations of the method by performing a series of controlled tests. In §4 we study a sample of three LMC star clusters and six UFDGs and compare our results with previous determinations. Our general conclusions are summarized in §5.

2 BAYESIAN INFERENCE

2.1 Isochrones as PDFs

An isochrone is the theoretical locus in the CMD defined by stars of a given age and metallicity according to a given set of stellar evolutionary tracks. The shape of the isochrone in the CMD depends on the age ($t$) and metallicity ($Z$) of the stellar population, and on the photometric bands in use. The isochrone must be shifted in the vertical (magnitude) axis by the distance modulus $(m-M)$ and the extinction ($A_{\lambda}$) along the line of sight to the stellar system in study. In the horizontal (colour) axis the isochrone is displaced according to the colour excess $E(\lambda_{1}-\lambda_{2})$, which is related to $A_{\lambda}$ by the reddening law. Without lack of generality, we write the $ith$-isochrone selected from a large set as the function

[TABLE]

Stars formed in a single burst are expected to fall in the CMD along the locus corresponding to a specific isochrone. Given our probabilistic approach, we treat $h_{i}$ as a probability distribution function (PDF): the probability of finding a star at a given position in the CMD is then proportional to the number of stars expected at this position according to the assumed stellar initial mass function (IMF).

In this work we use isochrones computed from the Chen et al. (2015) and Marigo et al. (2013) evolutionary tracks by Charlot & Bruzual (in preparation) using the Charlot & Bruzual (1991) isochrone synthesis algorithm. Isochrones are available for 16 values of the metallicity222$Z=0.00005$, 0.0001, 0.0002, 0.0005, 0.001, 0.002, 0.004, 0.006, 0.008, 0.01, 0.014, 0.017, 0.02, 0.03, 0.04, and 0.06 at ages ranging from $10^{4}$ yr to $15$ Gyr in 220 time steps. The BaSeL 3.1 atlas (Westera et al., 2002) is used to derive the photometric properties of the stars. The isochrones are populated stochastically as described by Bruzual (2010) following the Kroupa (2001) universal IMF assuming that stars form initially in the mass range from $m_{L}=0.10$ to $m_{U}=100$ M⊙. Isochrones must be densely populated to resemble as close as possible a continuous PDF in the CMD. By trial and error we found that this is achieved with a population of 1 to 2 million stars, warranting $\sim$ 100,000 stars above the MSTO, where the isochrones tend to be more sparsely populated due to the rapid transit of the stars through the post MS evolutionary stages. This is illustrated in Fig. 1.

2.2 Bayesian procedure

To determine the age and physical properties of the stars in a stellar population, we use a Bayesian approach to infer the stellar population parameters that reproduce the observed CMD. A detailed analytic description on the association of a likelihood value to a given set of parameters is presented in Tolstoy & Saha (1996), and Hernandez et al. (1999), and reviewed in Walmswell et al. (2013).

Let $o$ represent an observed sample of stars in a CMD and $h$ a set of hypotheses (isochrones) which depend on the parameters in Eq. (1) The Bayes rule provides the posterior probability $p(h_{i}|o)$ for each hypothesis given the observation

[TABLE]

where $p(o|h_{i})$ is the probability of obtaining the observation given the $ith$-hypothesis, and $p(h_{i})$ and $p(o)$ are the probability of $h_{i}$ and the observation $o$, respectively. The hypothesis $h_{max}$ that maximizes (2) is known as the estimate and provides the most likely combination of parameters that reproduces the observations. Finding the estimate is equivalent to estimating the free parameters, while the width of the likelihood distribution provides a measure of their confidence intervals. The quantity obtained by dividing Eq. (2) by $p(h_{max}|o)$ is referred to as the likelihood $\mathcal{L}(h_{i}|o)$,

[TABLE]

If, as is our case, we lack prior information on $p(h_{i})$, all $p(h_{i})$ are considered equal, which simplifies Eq. (3) to

[TABLE]

The resulting likelihood distribution, $\mathcal{L}(h_{i}|o)$, is not a PDF properly, but an indicator of the relative probability of occurrence of the hypothesis $h_{i}$.

To compute $p(o|h_{i})$ let us consider a single star $s_{k}$ observed at position $(x_{k},y_{k})$ in the CMD with observational errors $(\sigma_{x,k},\sigma_{y,k})$. The probability that such star arises from a starburst characterized by hypothesis $h_{i}$ (defined by the $ith$ set of parameters), is a function of the position of the star weighted by the errors

[TABLE]

where $U_{k}(x-x_{k},y-y_{k})$ is the error function for the $kth$ star, which we take to be the bivariate Gaussian

[TABLE]

where

[TABLE]

Our isochrones are discrete sets of $n$ stars. If $(x_{j},y_{j})$ is the position in the CMD of the $jth$ star in the isochrone, then Eq. (5) can be expressed as the sum

[TABLE]

where $n$ is the number of stars in the $ith$ isochrone, and $1/n$ corresponds to the normalization constant for this isochrone. From Eq. (8) we see that $p(s_{k}|h_{i})$ is defined if, and only if, $\sigma_{x,k}>0$ and $\sigma_{y,k}>0$, i.e., non zero observational errors must be associated with each star $s_{k}$.

In most cases, observed CMDs contain a large number $N$ of stars. Following Tolstoy & Saha (1996), we write the probability that the $N$ stars in the CMD come from a burst of star formation characterized by the hypothesis $h_{i}$ as

[TABLE]

and in logarithmic form as

[TABLE]

Evaluating the sum in Eq. (10) for the entire set of hypothesis $h_{i}$ we obtain the likelihood distribution function $\mathcal{L}_{N}$.333In practice, we compute $\mathcal{L}_{N}$ from Eq. (10) and divide the resulting distribution by its maximum, resulting in $\mathcal{L}_{N}^{max}=1$. From $\mathcal{L}_{N}$ we obtain its statistical properties, such as its mode, median, maximum, standard deviation, FWHM, skewness, etc. For a single parameter problem, the mode or median and the FWHM of $\mathcal{L}_{N}(\phi)$ provide reasonable estimates of the value of the parameter $\phi$ and its confidence interval $\sigma_{\phi}$. In the case of multi-parameter problems, estimating the value of each parameter and its confidence interval is more cumbersome, and using different schemes may lead to different answers (Basu, 1977; Pawitan, 2013). In this paper we use the marginalized likelihood distribution, $\mathcal{L}_{N}^{marg}(\phi_{q})$, to estimate the value of a parameter $\phi_{q}$. $\mathcal{L}_{N}^{marg}(\phi_{q})$ is commonly used in Bayesian inference studies when there is no information about the nuisance parameters, which are assumed to follow flat probability density functions (Verde, 2010; Pawitan, 2013). $\mathcal{L}_{N}^{marg}$ for $\phi_{q}$ is defined as:

[TABLE]

where $\nu$ refers to the set of nuisance parameters, i.e., all parameters except $\phi_{q}$ (Verde, 2010)444We divide the resulting distribution by its maximum, resulting in $\mathcal{L}^{marg}_{N}(\phi_{q})^{max}=1$.. The mode of the $\mathcal{L}_{N}^{marg}$ distribution leads to the estimate $\phi_{q}^{marg}$, the most likely value of the parameter that integrates all the possibilities of the nuisance parameters. The confidence interval $\sigma^{marg}_{q}$ is obtained from the half maximum of $\mathcal{L}_{N}^{marg}$.

If a problem stellar population arises from a double burst of star formation, we evaluate $\mathcal{L}_{N}$ for all possible linear combinations of pairs $(i,j)$ of isochrones in our hypothesis set, namely

[TABLE]

where the weights $w$, and $1-w$, account for the fraction of stars belonging to the $ith$ and the $jth$ isochrones, respectively. The formalism in Eq. (12) can be extended to multiple stellar populations, i.e, more than two stellar bursts (Dolphin, 1997). Each isochrone in (12) corresponds to a different star formation episode.

2.3 Parameters and priors

From the procedure described above, we estimate four basic parameters characterizing a stellar population $(t,Z,m-M,E(\lambda_{1}-\lambda_{2}))$. For a complete sample of observed stars, we can estimate the total mass of the stellar population. In the case of a double burst population described by Eq. (12), we can recover, in principle, seven parameters ($t_{1},t_{2},Z_{1},Z_{2},m-M,E(\lambda_{1}-\lambda_{2}),w_{1}/w_{2}$).

The number of free parameters to estimate depends on our prior knowledge on them. Using prior information reduces computational time and resources, and may help to break degeneracies in the solution. For instance, prior information about the metallicity of the stellar population, e.g., $Z\pm\sigma_{Z}$, can be adopted to break the age-metallicity degeneracy. The uncertainty $\sigma_{Z}$ is introduced in the analysis by marginalizing over the isochrone sets that cover the range of possible values of $Z$. As a rule, in this work we marginalize over isochrone sets expanding a range of $\pm 5\sigma_{Z}$ (see §4.1).

Prior knowledge on the values of $(m-M)$ and/or $E(\lambda_{1}-\lambda_{2})$, can be adopted to solve the problem only for $t$ and $Z$. In this case, the associated uncertainties, $\sigma_{(m-M)}$ and $\sigma_{E(\lambda_{1}-\lambda_{2})}$, should be added in quadrature to the photometric uncertainties of the stars in the CMD, i.e.,

[TABLE]

where $\sigma_{\lambda_{1}-\lambda_{2}}$ and $\sigma_{\lambda_{2}}$, are the photometric uncertainties on the x and y-axis, respectively (see §4.1). $\sigma_{S_{x}}$ and $\sigma_{S_{y}}$ correspond to the systematic uncertainties.

In absence of previous (or not adequate) measurements of a nuisance parameter, we can estimate its value using our procedure by setting it as a free parameter. Often we have information that can be used to delimit the range of possible values of this parameter, e.g., when analyzing stellar populations in the Milky Way, we constrain the domain of $(m-M)$ to Galactic scales. Delimiting properly the domain of a free parameter is important since the marginalized distributions may depend on the chosen domain, as happens with any prior (Verde, 2010). In our case, when a previous measurement of a nuisance parameter is available (let’s say $\phi=\phi_{0}\pm\sigma_{0}$), and we decide to estimate its value, we choose its domain as: $[\phi_{0}-5\sigma_{0},\,\phi_{0}+5\sigma_{0}]$.

2.4 Advantages of using unbinned CMDs

Our method compares model and observed CMDs on a star by star basis, and does not require a grid of binned data, in contrast to $\chi^{2}$ or Poisson based statistics (Harris & Zaritsky, 2001; Dolphin, 2002), avoiding the dilemma of arbitrarily defining the size, form, and position of the bins in the CMD. Bins should be small enough to resolve sensitive stellar phases on each CMD, but large enough to contain a significant number of stars. Data binning degrades the CMD resolution, results in a lost of valuable information, may introduce unnecessary noise and statistical errors, and compromises the effectiveness of these methods for the study of scarce stellar populations, as pointed out by Walmswell et al. (2013). Badly chosen bins may affect the results of $\chi^{2}$ based fits, and may introduce undesired biases and uncertainties. We avoid rejecting valuable information from stars in poorly filled bins, which allow us to estimate the physical parameters of scarce stellar populations. Using unbinned data is also advantageous in the context of processing automatically and blindly catalogues containing large numbers of CMDs,555 In general, these catalogues are inhomogeneous, and several properties of the CMDs which are critical in designing an appropriate grid of bins, such as the number of stars, their dispersion in position due to the photometric uncertainties, the overall morphology due to different age or reddening, or the possibility that the catalogue mixes data from different observational programs, vary from cluster to cluster. as has been shown by Bitsakis et al. (2017, 2018) using our code on LMC and SMC star cluster catalogues. A brief review on binned and unbinned methods can be found in Walmswell et al. (2013).

2.5 Trimming the sample

Some regions in the CMD are more sensitive to differences in age than others. In general, the brighter the star, the more information it will contribute to unveiling the population parameters under analysis. The differences in the isochrones are more pronounced around the MSTO and the SGB. In contrast, the isochrones overlap towards fainter magnitudes reaching a point where they become indistinguishable (cf. Fig. 1). Therefore, observed stars in the lower part of the MS do not contribute much to our analysis, and since the observational uncertainties increase towards the bottom of the CMD, including lower MS stars only adds noise to our results. In the next section we show, by applying different magnitude cuts to the sample, that including fainter stars does not change the values estimated for the different parameters, only broadens the confidence intervals. Thus, to lower the uncertainties in the parameter estimates, lower MS stars are usually removed from the sample.

3 Mock Galaxies

To test the capabilities and limitations of the parameter estimator algorithm described in §2.2, we generate a series of mock stellar populations and then use the estimator to recover the input values of ($t,Z,m-M,E(\lambda_{1}-\lambda_{2})$) for these synthetic observations. Each mock galaxy is generated as described in §2.1 and its population is left to evolve up to the desired age, that we choose arbitrarily for each experiment. To produce a realistic ($V-I,I$) CMD we add random photometric Gaussian uncertainties to the $V$ and $I$ magnitude of each star on the isochrone. For the mock galaxies discussed in this section we model the photometric errors as a function of apparent magnitude using the relations

[TABLE]

derived from fits to the photometric uncertainties of the UFDGs observed with the HST by Brown et al. (2014).666 In reality, the HST observations were performed in the $F606w$ and $F814w$ bands. For the purpose of the tests discussed in this section, we adopt for the $V$ and $I$ bands the photometric uncertainties characteristic of the $F606w$ and $F814w$ bands, respectively. We allowed for a systematic uncertainty of $\sigma_{S}=0.02$ by adding it in quadrature to the photometric uncertainties before performing the analysis (Eq. 13).

3.1 Mock galaxy with a single stellar population

Fig. 2a shows the ($V-I,I$) CMD of a 11.125 Gyr old, $Log\,Z=-4$ mock simple stellar population at $(m-M)=21.9$ with $E(V-I)=0.085$. To make the test more meaningful, we choose an age that is not in our isochrone grid, adding in equal numbers of stars belonging to the 11 and 11.25 Gyr isochrones777 We use the same combination of ages in all the mock bursts considered in this paper. The resulting mock galaxy contains 1608 stars, only 160 of them above the MSTO, resembling an old dwarf satellite galaxy of the Milky Way. We derive $\mathcal{L}$ using stars brighter than three different magnitude cuts in the CMD, namely $I=25.25$, 26.25, and 27.25. Figs. 2b-e show the resulting marginalized likelihood distributions $\mathcal{L}_{N}^{marg}$. We recover closely the input value of three of the free parameters ($age$, $(m-M)$, and $Z$), while $E(V-I)$ is slightly shifted to lower values. Estimates on $age$, $(m-M)$, and $Z$ do not change by including fainter stars in the analysis, but $\mathcal{L}_{N}^{marg}$ for $E(V-I)$ does change, shifting slightly towards lower reddening values. This effect is mainly due to the reddening-metallicity degeneracy. Fig. 3 shows the 2D $\mathcal{L}_{N}^{marg}$ distributions for all pair of fitted parameters. The partial relative contribution to the $E(V-I)$ $\mathcal{L}_{N}^{marg}$ by stellar populations of different $Z$ can be seen in Fig. 3a. Higher $Z$ values lead to lower reddening values. These figures clearly show that degeneracies may dominate the inference results, producing non-unique solutions.888 A detailed study of the degeneracies in isochrone fitting in the context of a maximum likelihood approach based on Monte Carlo methods can be found in Hernandez & Valls-Gabaud (2008).

Fig. 2 shows that as we include fainter stars, all the $\mathcal{L}_{N}^{marg}$ distributions broaden, indicating that no meaningful information is added by including stars far below the MSTO, as discussed in §2.5. At first sight, this effect may seem puzzling, since one would expect that increasing the sample should result in reduced uncertainties. Here we note that as we include fainter stars, the proportion of stars with large photometric uncertainties increases dramatically. The fainter stars are located in zones of the CMD where the isochrones tend to be identical (cf. Fig. 1); thus, all isochrones are equally able to recreate these faint stars, adding to the same extent to each isochrone’s likelihood while reducing the relative weight from the few bright stars in regions of the CMD where the isochrones clearly separate (around and above the MSTO), leading to wider likelihood distributions. We adopt as our best solution the $\mathcal{L}_{N}^{marg}$ distributions for stars brighter than $I=25.25$, resulting in the parameters listed in Table 1.

3.2 Mock galaxy with two stellar populations

3.2.1 Two old stellar populations

In this section we explore how our parameter estimator handles double stellar populations. We consider a $Log\,Z=-4$ mock galaxy in which $w_{o}=0.5$ of the stars were formed in an initial burst999 As in §3.1, there are no isochrones in our grid for age $t_{o}$ and $t_{y}$. For the older population we add in equal amounts stars from the 12 and 12.25 Gyr isochrones, and for the younger population we add stars from the 8 and 8.25 Gyr isochrones. of age $t_{o}=12.125$ Gyr and $w_{y}=0.5$ of the stars in a younger burst happening at $t_{y}=8.125$ Gyr, seen at a $(m-M)=21.9$ with $E(V-I)=0.085$. Photometric uncertainties were added according to Eq. (14), adding in quadrature a systematic uncertainty of 0.02. The total sample amounts to 1600 stars, the CMD is shown in Fig. 4a. Only the 394 stars above $I=25.5$ were taken in consideration for the analysis. For each burst we estimate its age and weight, along with their common $(m-M)$ and $E(V-I)$, assuming as flat priors $21.4<(m-M)<22.4$, $0.03<E(V-I)<0.13$, and $Log\,Z=\,$-3, -3.3, -3.7, -4, and -4.3. Figs. 4b-d show the resulting $\mathcal{L}_{N}^{marg}$ distributions for $t_{y}$ and $t_{o}$ as 2D and 3D representations. The distribution for $t_{o}$ looks truncated at the high age end in Fig. 4b because the 14.75 Gyr isochrone ($\sim 1$ Gyr older than the universe) is the oldest one in our grid. Figs. 4e-f show the $\mathcal{L}_{N}^{marg}$ distributions for $Log\,Z$, $w_{y}$, $(m-M)$, and $E(V-I)$.

Our parameter estimates closely correspond to the input values listed in Table 2. We remark that the confidence interval for $t_{o}$ is wider than for $t_{y}$, even though in relative terms they are similar ($\sim 23$% for $\mathcal{L}_{N}^{marg}$). Two factors contribute to broaden the confidence interval for $t_{o}$: (1) given the position of the younger starburst in the CMD, its MS stars contaminate the MSTO of the older starburst, while the opposite does not occur (this is the rule for double stellar populations, unless the photometric precision is outstanding), and (2) the MSTO stars of the younger burst are 0.5 mag brighter and have smaller photometric errors than the older burst stars. Both $t_{o}$ and $t_{y}$ have larger confidence intervals as compared to the age of the simple population discussed in §3.1. This is a consequence of the fact that stars from both bursts share common regions in the CMD, contaminating and broadening critical zones in this diagram with respect to the simple population case. Even though the estimated $w_{y}$ is close to the input value, the $\mathcal{L}_{N}^{marg}$ distribution for this parameter in Fig. 4f is wide, flat to a degree, and non-Gaussian, hence the confidence interval for $w_{y}$ is ill defined. The probability of a single stellar population never falls below 60%, hence, these distributions can not be used to rule out this possibility. If the double burst nature of a population is known a priori, our method recovers successfully the age of each burst. The accuracy of the $Log\,Z$, $(m-M)$, and $E(V-I)$ estimates are similar to those in §3.1.

To further explore our parameter estimation power we generated several mock galaxies with identical parameters to the previous one but varying $w_{y}$ as indicated in Table 3, where we also list our results. As expected, the confidence intervals in $t_{o}$ and $t_{y}$ tend to increase as the weight of the concerned burst becomes smaller. This is manly due to the pollution by stars from the more numerous burst. In Table 3 there are two cases in which one of the bursts contributes 90% of the stars. We can consider these two cases as simple populations with a 10% of contamination of stars from the other burst, and compare with the results for the simple population in §3.1. The estimated age for the burst with 90% of the stars has a similar relative uncertainty ($\sim 21$% for $\mathcal{L}_{N}^{marg}$) to that of the simple population. The age of the contaminating bursts (contributing with 10% of the stars) is recovered to some extent in both cases, with a relative uncertainty reaching values as high as $\sim 42\%$. As expected, the accuracy of our age estimates increases with the relative strength of the burst.

3.2.2 Two young stellar populations

Here we repeat the experiment performed in the last part of §3.2 but for a $Log\,Z=-2.1\pm 0.04$ mock galaxy containing two young stellar populations with ages $t_{o}=6.12$ and $t_{y}=4.05$ Myr $(m-M)=29.0$, and $E(V-I)=0.15$. The input weight $w_{y}$ assigned to the younger population is listed in Table 4 together with our results. The galaxy contains 1200 stars but only the 403 stars brighter than $I=25.25$ were used. As priors we use $28.5<(m-M)<29.5$, $0.10<E(V-I)<0.2$. We assume a prior knowledge of $Log\,Z=-2.1\pm{0.05}$, and marginalize over the isochrone sets with $-2.4<Log\,Z<-1.85$, covering the metallicity domain at the $\pm 5\sigma$ level (cf.§2.3). Figs. 5a-e show the CMD, and the ($t_{y},t_{o},m-M,E(V-I),w_{y}$) $\mathcal{L}_{N}^{marg}$ distributions for $w_{y,input}=0.70$. Figs. 6a-b show the CMD, and the $(t_{y},t_{o})$ $\mathcal{L}_{N}^{marg}$ distributions for $w_{y,input}=0.10$. In the latter case only the oldest burst is recovered. The recovery of the parameters for the double young burst galaxy shows similar dependence on $w_{y}$ as the double old burst galaxy discussed in §3.2.1. Again, the $\mathcal{L}_{N}$ distributions for $w_{y}$ are not sharp enough as to discard a single burst nature for this population.

4 Applications

4.1 Star Clusters in the LMC and SMC

One of the advantages of our parameter estimator is that it can be used in a fast and objective manner to characterize resolved stellar groups or clusters in nearby galaxies. Using this procedure Bitsakis et al. (2017, 2018) determined the age of 4850 star clusters in the LMC and 1319 star clusters in the SMC, respectively. The clusters were previously defined by these authors as spatial overdensities in the Magellanic Cloud Photometry Survey,101010 Carried out by Zaritsky et al. (2002, 2004) with the 1m diameter Swope Telescope at Las Campanas Observatory. The data comprise two photometric catalogs of the stars in the LMC and the SMC in the Johnson $BVI$ and Gunn $i$ bands, corrected for internal and Galactic extinction. using a new statistical analysis. The resulting candidate clusters are cleaned from contaminating stars using the following probabilistic algorithm. They draw a box that encloses all the stars defining a candidate cluster, and a similar box in a neighbouring field, and build the CMD of the stars inside each box. Both CMDs are equally binned along their axes, and the number of stars inside each bin is counted. Following Mighell et al. (1996), the membership probability of the stars in a bin of the candidate cluster CMD is given by

[TABLE]

where $N_{*,cluster}$ and $N_{*,field}$ are the number of stars in the cluster and field bins, and ${A_{cluster}}$ and ${A_{field}}$ are the areas of the cluster and field boxes.

To quicken their analysis, Bitsakis et al. (2017, 2018) dated all the star clusters assuming a single stellar population. In this section we study in detail three LMC clusters from Bitsakis et al. (2017): IR1-1959, NUV-1781, and IR1-297. For the first one we use a single stellar population in our analysis, for the other two, a double stellar population. Considering these examples is important since the CMD diagrams of these clusters exhibit more noise and scatter than the HST CMD of UFDGs analyzed in the next section, resulting in a more challenging parameter estimation. Examples of two clusters in the SMC are shown by Bitsakis et al. (2018). For the LMC clusters we have prior knowledge of $Z,(m-M)$, and $E(B-V)$. We adopt $Log\,Z=-2.04\pm 0.03$ and $(m-M)=18.48\pm 0.05$, derived from Cepheids studies in the LMC by Keller & Wood (2006), and (Walker, 2012), respectively. By employing a set of isochrones with $Log\,Z=-2.22,-2.10,-2.00,-2.9,-2.85,-1.85$, we cover the metallicity domain at a $\pm 5\sigma$ level. The uncertainty on $(m-M)$ was added in quadrature to the y-axis uncertainties in the CMD. Since the photometry of the stars has been corrected for extinction, we assume $E(B-V)=0$ and add in quadrature its 0.125 mag uncertainty to the x-axis uncertainties in the CMD (Zaritsky et al., 2004). Finally, an overall systematic uncertainty of 0.02 mag is added in quadrature to both axis. Then, our problem reduces to estimating the age $t$ of the stellar populations (cf.§2.3).

4.1.1 LMC IR1-1959 as a single stellar population

Fig. 7a shows the $(B-V,V)$ CMD for 700 stars in the cluster LMC IR1-1959, $\sim 30$ of them above the MSTO. The $\mathcal{L}^{marg}_{N}$ distributions for $t$ shown in Fig. 7b were obtained using the $V$ magnitude cuts indicated in the Fig. 7a. The three subsamples lead to the same age estimate, but the distributions broaden as we include fainter stars. For the subsample of brighter stars, $t=80^{+18}_{20}$ Myr. The RHS panel in Fig. 7c shows an 80 Myr simple stellar population mock CMD generated as described in §3 next to the observed CMD (left hand side (LHS) panel). The 60, 80, and 100 Myr isochrones are shown in Fig. 7a.

4.1.2 LMC NUV-1781 as a double stellar population

The CMD of LMC NUV-1781 shown in Fig. 8a resembles that of a cluster with two stellar populations. This cluster has been cleaned of contaminating stars, however, there is still the possibility that remnant field stars may be defining a false secondary population, without being physically associated with the cluster (Cabrera-Ziri et al., 2016). To account for this possibility, we adopt the approach of Bitsakis et al. (2017, 2018). The probability that each star is drawn from an isochrone (Eq. 5) is weighted by the membership probability: $p_{memb}\times p(s_{k}|h_{i})$.

From the $\mathcal{L}_{N}^{marg}$ distribution in Fig. 8b, LMC NUV-1781 contains a young population of age $t_{y}=6_{-5.96}^{+58}$ Myr, and an older population of age $t_{o}=400_{-150}^{+420}$ Myr. The $t_{y}$ confidence interval ranges from 0.04 to 64 Myr. The similarity of the isochrones in this age range, the poorly defined upper MS (due to the scarcity of stars), and the photometric uncertainties, are responsible of the width of this confidence interval. In contrast, the $\mathcal{L}_{N}^{marg}$ distribution for $t_{o}$ in Fig. 8b results in a narrower, in relative terms, confidence interval for $t_{o}$, ranging from 250 to 820 Myr. The presence of RGB stars in the older burst (in contrast with the lack of them in the younger one) provides a second feature that drives our matching isochrone algorithm to obtain a more precise estimate for $t_{o}$ than for $t_{y}$. The $\mathcal{L}_{N}^{marg}$ distribution for $w_{y}$ (Fig. 8c) indicates that the younger burst accounts for $75^{+25}_{-30}\%$ of the stars. This result is consistent with the presence of a single stellar population in this cluster. Higher precision photometry is needed for a more precise age estimate, as well as for supporting the double population hypothesis.. The RHS panel of Fig. 8d shows a mock CMD for this population next to the observed one (LHS).

4.1.3 LMC IR1-297 as a double stellar population

The CMD of LMC IR1-297 is presented in Fig. 8e. The $t_{y}$, $t_{o}$ and $w_{y}$ $\mathcal{L}_{N}^{marg}$ distributions (Fig. 8f, g) were obtained in the same way as for LMC NUV-1781. From these distributions, we estimate $t_{y}=110_{-60}^{+80}$ and $t_{o}=300_{-110}^{+370}$ Myr. The $\mathcal{L}_{N}$ distributions for $w_{y}$ indicate that the younger burst may account for up to $\sim 80^{+20}_{-40}$% of the cluster population. As for NUV-1781, higher precision photometry is required for a more precise age determination and to exclude or validate the double population hypothesis for this cluster. The RHS panel of Fig. 8h shows a mock CMD for this population next to the observed one (LHS).

4.2 UFDGs

In this section we discuss our analysis of six UFDGs, satellites of the Milky Way: Bootes I (Boo I), Canes Venatici II (CVn II), Coma Berenices (Com Ber), Hercules (Her), Leo IV, and Ursa Major I (UMa I). Deep optical images of these galaxies were obtained with the Advanced Camera for Surveys (ACS) on the HST using filters F606w and F814w (GO-12549; PI: T.M. Brown). The photometric data reduction was carried out by Brown et al. (2012); Brown et al. (2014), who also obtained metallicity distributions from spectroscopic observations. The photometry includes stars at least three magnitudes below the MSTO. Contaminating stars and background galaxies were cleaned by rejecting sources with poor photometry, and with no typical stellar profiles on their point spread functions. For a detailed description of the observations and data reduction we refer the reader to Brown et al. (2012); Brown et al. (2014).

The ACS photometry of the identified UFDG stellar sources with photometric uncertainties is provided in the STMAG system as $m_{606}$ and $m_{814}$, where

[TABLE]

with an overall systematic uncertainty of 0.02 mag. This amount is added in quadrature to the photometric uncertainties. We computed a grid of isochrones in the STMAG system spanning from 7 to 14.75 Gyr, in 0.25 Gyr steps, for $log\,Z=-3.3,-3.7,-4$, and $-4.3$, which bracket the value of $log\,Z\sim-4$ for UFDGS from previous works (Brown et al., 2014; Kirby et al., 2008).

For each UFDG we then estimate the parameters $t,Z,(m-M)$, and $E(B-V)$. The extinction correction in STMAG for cool stars was obtained from Bedin et al. (2005).

Fig. 9 shows the CMDs for the six UFDGs.111111 Stars marked in pink in Fig. 9 were removed, in order to clean contaminating outlier sources (such as foreground stars or background galaxies) that are non-modellable points by the isochrone set in consideration. These stars carry a very small $p(s_{k}|h_{i})$, which may be critically risky, in some cases, since a single star with $p(s_{k}|h_{i})=0$ may collapse the whole likelihood during the multiplication process of Eq. (9), as noted by Tolstoy & Saha (1996).

Two different subsampling trimming cuts in $m_{814}$ are used (cf. §2.5). Fig. 10 shows two $\mathcal{L}_{N}^{marg}$ likelihood distributions for each parameter for each galaxy, corresponding to the two cuts in $m_{814}$. The estimated values of the parameters for each UFDG are listed in Table 5, where for comparison we also list previous determinations by other authors. The estimates from the marginalized distributions $\mathcal{L}_{N}^{marg}$ result in an age between 13.25 and 14.25 Gyr, with a relative uncertainty of about 21%. These estimates include only the lower limit for the confidence interval. This is because our isochrones do not extend past 14.75 Gyr, which is not enough to establish the upper limit for the confidence interval. Given the latest estimate of the age of the universe of 13.8 Gyr (Planck Collaboration et al., 2016), this is not a concern for objects as old as this set of UFDGs. Fig. 9 shows the most likely isochrone in the CMD for each UFDG, including the low age confidence boundary isochrone. The high age confidence boundary isochrone is shown only for Boo I. The accuracy of our estimates of $log\,Z,(m-M)$, and $E(m_{606}-m_{814})$ from $\mathcal{L}_{N}^{marg}$ are $\sim 0.5$, $\sim 0.4$ mag, and $\sim 0.04$ mag, respectively. Taking into account the confidence intervals, these results are in line with the values reported by the previous studies listed in Table 5. Our results support a scenario in which the star formation in these UFDGs occurred as early as redshift $z\sim 6$.

Our age estimates are consistent with those of Brown et al. (2014), especially for Boo I, CVn II and Com Ber, for which our estimates are within the reported confidence intervals. Our results are also in good agreement with Weisz et al. (2014) for the case of Her, and Leo IV. For the case of CVn II, there is a discrepancy between the estimated ages by Brown et al. (2014) and Weisz et al. (2014), the latter being significantly younger. Our result supports the age determined by Brown et al. (2014). The younger age for CVn II by Weisz et al. (2014) could be due to their lower quality photometric data, which is remarkably noisier about the MSTO, as already discussed by Brown et al. (2014).

The comparison between the age estimates by Brown et al. (2014) and our work is particularly relevant, considering that both works rely on the same CMDs. In general, the age estimates are consistent, but ours present larger uncertainties. Part of the differences between their results and ours is due to the different methods, with a non-negligible contribution from the different set of isochrones. Whereas our isochrones are computed in steps of 0.25 Gyr, the step in their grid is 0.1 Gyr. Their isochrones follow the proportions of the metallicity distribution from their spectroscopic data, using a grid with a fine step of 0.2 dex in [Fe/H], while our isochrones are limited to discrete published values of $Z$. Moreover, unlike our isochrones, theirs include binary systems, and thus they may be more realistic since binary stars in UFDGS may reach up to 50% (Geha et al., 2013). Even though we do not reach the accuracy of the work by Brown et al. (2014), the goal of our technique is focused on estimating simultaneously a set of multiple parameters using only stellar photometry. Age determinations by procedures like that of Brown et al. (2014) are preceded by several steps, e.g., measuring the metallicity, and estimating $(m-M)$ and $E(m_{606}-m_{814})$. The advantage of our approach is the fast and simultaneous multiple parameter estimation, intended for automated analysis of massive data sets.

4.2.1 Searching for double bursts in UFDGs

We searched for a possible second stellar population in the UFDGs of our sample. From our $\mathcal{L}_{N}$ distributions there is a hint of a possible double population only for UMa I (Fig. 11b,c), but they do not rule out the hypothesis of a single burst. Considering UMa I as a double stellar population, the most likely isochrones describe two possible SGBs in the CMD of this galaxy, as seen in Fig. 11a. The $\mathcal{L}_{N}^{marg}$ distributions lead to ages of $t_{y}=10.0^{+3.5}_{-5.0}$ Gyr for the younger population and $t_{o}=14.75_{-2.0}$ Gyr for the older one, with a relative weight for the younger population of $w_{y}=0.25^{+0.50}_{-0.25}$. Taking into account the uncertainties, there is not conflict between $t_{o}$ and the age of the universe. These estimates are in good agreement with those by Brown et al. (2014), $t_{y}=11.6$ Gyr, $t_{o}=14.1$ Gyr, and $w_{y}=0.55$. Considering the overlap of the $\mathcal{L}$ distributions for $t_{y}$ and $t_{o}$ in Fig. 11a, we cannot rule out a scenario of continuous star formation starting $\sim 13.5$ Gyr ago and lasting for $\sim 1.5$ Gyr until $\sim 12$ Gyr ago. Nevertheless, support for the double population nature of UMa I was found from its metal content by Webster et al. (2015). Webster et al. (2015) compared the observed and the theoretical metallicity distributions of the UFDGs in our sample. For UMa I, they explored two different enrichment scenarios. In the first scenario, an instantaneous burst occurs 14.1 Gyr ago with a extremely poor metallicity. As the stars evolve, they enrich the medium from which a second population forms 2.5 Gyr later (using the parameters determined by Brown et al., 2014). In the second scenario, stars form continuously for $\sim 0.1$ Gyr, enriching the medium as they age. They conclude that the observed distribution of metals in UMa I favours the first scenario with a double burst.

For the other UFDGs in our sample, Boo I, CVn II, Com Ber, Her, and Leo IV we do not find significant evidence of the presence of a second population. Either the estimates of $t_{y}$ and $t_{o}$ are so close within each other inside the confidence interval, or the fraction of stars attributed to one of the bursts is too small ($\leq 15$%). Webster et al. (2015) find that the metallicity distributions of these five galaxies are reproduced better by a single continuous star formation model than by a double population. For Leo IV, their results are not very conclusive due to the scarcity of stars.

5 Conclusions

We present a fast code based on Bayesian inference with the main purpose of estimating the age, metallicity, distance modulus, and color excess of stellar populations from their CMD. We introduce a new method to estimate the ages of double (or multiple) stellar populations, as well as the stellar contribution of each star burst. The code evaluates the likelihood, $\mathcal{L}$, of each combination of parameters, producing a $\mathcal{L}$ distribution. We estimate the parameter values and their confidence intervals from the marginalized likelihood distributions, $\mathcal{L}^{marg}_{N}$. The possibility of evaluating the parameter in this manner represents an alternative to the Monte Carlo approach, which requires a large number of calculations to compare several stochastically populated isochrones with the observed CMD (Dolphin, 2002).

We test our procedure on mock CMDs of simple stellar populations, recovering the full set of input parameters ($age,Z,m-M,E(V-I)$), and examined its performance on the age estimates depending on the number of stars on the CMD. We also explored the effects of the sample selection on the resulting parameter estimates. We trimmed the observed and mock CMD samples by discarding the stars fainter than various cuts in apparent magnitude, and found that keeping only the stars brighter than $\sim 0.5$ mag below the MSTO results in significantly lower uncertainties. Including fainter stars is not inoffensive, enhancing the uncertainty on each of the estimated parameters.

We also study mock CMDs of double stellar populations. We closely recover the correct values of each input parameter for each population, except in the cases of poorly populated bursts and of bursts of very similar age. For the latter we recover the age of the most prominent burst only. The accuracy of the estimated age for each burst depends remarkably on its relative contribution to the number of stars in the double population. The age of the dominating burst is determined with higher accuracy. Although the age estimate of each component of the double populations considered in this work are satisfactory, the resulting $\mathcal{L}$ distributions of the weight of each population tend to be too wide, and non-Gaussian. They entail a low accuracy which do not allow us to discard the possibility of simple stellar populations. To establish the meanigfulness of the pair of estimated ages additional information favouring a double population is required.

We used our tool to analyze in detail three LMC star clusters from the Bitsakis et al. (2017) sample. The CMD of cluster IR1-1959 clearly corresponds visually to a simple stellar population. The CMDs of NUV-1781 and IR1-297 show visual signs of double populations, an hypothesis that was tested in our analysis. Using known priors for $Z$, $(m-M)$, and $E(B-V)$, for the LMC, we estimated the age of IR1-1959, and the age and relative fraction of each stellar population in NUV-1781 and IR1-297. Further work is required in order to establish the reliability of the number of star bursts that have occurred in these clusters.

Finally, we determined $t,Z,(m-M)$, and $E(B-V)$ for 6 UFDGs from their HST ACS $(m_{606}-m_{814},m_{814})$ CMDs. We obtained results consistent with previous works for all the parameters (Table 5). Our results support the notion that UFDGs contain very ancient stars of first or second generation formed in the early universe. The age of these galaxies is older than $\sim 13.5$ $Gyr$ and they are very metal poor, with $log\,Z\sim-4$. We searched for a possible second stellar population in the UFDGs in our sample, and found low signs of a double population for the case of UMa I, in agreement with previous work by other authors, although, our analysis also supports the single stellar burst hypothesis.

6 Acknowledgments

The authors thank Dr. T.M. Brown and his team for sharing the UFDG data used in this work. VRS thanks X. Hernández and M. Cerviño for illuminating discussions during the early phases of this project. The authors acknowledge the thorough reading of the original manuscript by the reviewer, Nicolas Martin, as helpful in reaching a clearer and more complete final version. The research in this paper is part of the PhD thesis of VHS in the Universidad Nacional Autónoma de México graduate program in astrophysics. He thanks the support from the Instituto de Radioastronomía and Astrofísica, its staff, and CONACyT for the scholarship granted. BCS acknowledges financial support through PAPIIT project IA103517 from DGAPA-UNAM. GB acknowledges financial support through PAPIIT projects IG100115 and IG100319 from DGAPA-UNAM. TB acknowledges support from the CONACyT Research Fellowships program. We gratefully acknowledge support from the program for basic research of CONACyT through grant number 252364. The Figures in this article were produced using the TOPCAT graphics viewer and tabular data editor.