STOCHASTIC STELLAR CLUSTER INITIAL MASS FUNCTIONS: MODELS AND IMPACT ON INTEGRATED CLUSTER PARAMETER DETERMINATION

Stellar clusters represent an essential tool to study the stellar populations of galaxies and their evolution. Only stellar clusters can be studied individually out to distances of 100 Mpc or slightly beyond (see, e.g., de Grijs et al. 2003), far beyond distances that allow studies of resolved stellar populations.

Multi-band integrated photometry of stellar clusters is regularly obtained to study large cluster samples. The clusters' spectral energy distributions (SEDs) are compared with model predictions (Leitherer et al. 1999; Anders & Fritze-v. Alvensleben 2003; Bruzual & Charlot 2003; Kotulla et al. 2009; Popescu & Hanson 2009) to provide estimates of the clusters' physical parameters (see, e.g., Anders et al. 2004; de Grijs et al. 2005). The derived parameters are cluster age, mass, metallicity, and foreground extinction. While observed clusters represent populations with finite numbers of stars, the model predictions commonly used are based on fully sampled stellar initial mass functions (IMFs). Stochastic fluctuations around this smooth IMF are consistently taken into account only in the MASSCLEAN models (Popescu & Hanson 2009). Other studies have been performed (Fouesneau & Lançon 2010; Piskunov et al. 2011; Silva-Villa & Larsen 2011), although the stochastic cluster models they obtained are not publicly available.

It is generally believed that a galaxy's field star population originates from stellar clusters, which have either dissolved quickly after their formation or slowly dissolved due to long-term evolution. By combining the stochastic cluster models with a galaxy's star formation history, stochastic effects in field star formation can be taken into account. Non-universal IMFs are found from the analysis of Sloan Digital Sky Survey (SDSS) galaxies, predominantly for faint galaxies (Hoversten & Glazebrook 2008). This effect is related to the Integrated Galactic IMF (IGIMF), which proposes that the IGIMF originates from dissolved stellar clusters (see, e.g., Weidner & Kroupa 2004, 2005).

In this paper we provide a new publicly available model suite of stellar cluster models including the impact of stochastic IMF sampling. This model set is an extension of the GALEV models (Kotulla et al. 2009). The cluster photometry obtained is used to derive physical cluster parameters, i.e., cluster age, mass, metallicity, and foreground extinction. The stochastic models are described in Section 2, where the resulting magnitudes and colors are presented and explained. Section 3 describes the principles of cluster parameter determination based on χ² minimization, and Section 4 shows the results of stochastically sampled clusters analyzed with fully sampled (i.e. non-stochastic) cluster models. The paper presents conclusions in Section 5.

The integrated stellar cluster photometry is derived using the GALEV population synthesis models (Kotulla et al. 2009). Population synthesis models combine stellar evolution data (usually stellar isochrone data), spectral libraries, and an IMF (to weigh the contributions of stars of different masses) to provide an integrated spectrum of the stellar cluster being studied. Convolving the integrated spectrum with filter response curves provides model predictions of integrated magnitudes. Such population synthesis models form the input for analysis and interpretation of observed integrated cluster photometry.

The standard GALEV models use an analytic description of the IMF, similar to most other population synthesis models (i.e., "fully sampled models"). The cluster mass M_cl acts solely as a scaling factor of the integrated magnitude predictions (see Kotulla et al. 2009). The derived cluster colors are independent of cluster mass. The fully sampled approach uses fractional stars to account for the total cluster mass. While this is applicable for the description of idealized clusters, it results in deviations for realistic clusters.

In real stellar clusters, a stellar cluster's IMF is derived by stochastically sampling stars from the analytic IMF description up to the total stellar mass M_cl (i.e., "stochastically sampled models"). In this paper we use stochastically sampled Kroupa (2001) IMFs in the mass range 0.1–100 M_☉ (with dN/dm∝m^−1.3 for stellar masses 0.1 M_☉ ⩽ m ⩽ 0.5 M_☉, dN/dm∝m^−2.3 for stellar masses 0.5 M_☉ < m ⩽ 100.0 M_☉) for a range of total cluster masses (10³–2 × 10⁵ M_☉). Observations have indicated the existence of an upper stellar mass limit (USML; as derived by Pflamm-Altenburg et al. 2007) as a function of cluster mass. For clusters with total masses of 10³ M_☉ we also present models with the appropriate USML of 42.7 M_☉. However, other observations (see, e.g., Calzetti et al. 2010; Andrews et al. 2013) found no strong evidence for a USML as a function of cluster mass exceeding the expectations from stochastic IMF sampling. Calzetti et al. (2010) claim that their observations are only 1σ–2σ away from the USML models, with the observations in between the USML models and models with a cluster-mass-independent stellar mass limit; hence, they cannot unambiguously evaluate the existence of a USML. Therefore, models with and without a USML are provided.

Pflamm-Altenburg et al. (2007) found a maximum stellar mass m_max as a function of stellar cluster mass. Instead, we use their m_max as an upper stellar mass limit, which can be reached, but not exceeded. Therefore, a stochastic distribution of the maximum stellar mass present in a stellar cluster is established, in contrast to Pflamm-Altenburg et al. (2007). The distribution of drawn m_max, based on 2 × 10⁷ test IMFs with 6.96 M_☉ ⩽m_max ⩽ 42.7 M_☉, is shown in Figure 1. In this respect the maximum stellar mass becomes a stochastic quantity, which also applies in general to the model clusters' IMFs.

Zoom In Zoom Out Reset image size

The IMFs of model clusters are populated with individual stars drawn stochastically from the Kroupa (2001) IMF. This process is stopped when the total cluster mass just exceeds the target mass. Due to the piecewise power-law description of the Kroupa (2001) IMF and to simplify the implementation of the IMF sampling for the entire mass range, the stochastic sampling process consists of three steps.

• 1.  
  Determination of the mass fraction (MF) below the Kroupa (2001) IMF break point at 0.5 M_☉.
• 2.  
  Performance of stochastic sampling of stars with M ⩽ 0.5 M_☉, until the total cluster mass reaches MF × M_cl. The last sampled star causes a slight excess over the mass MF × M_cl.
• 3.  
  Stochastic sampling of stars with M > 0.5 M_☉, until the total cluster mass reaches M_cl. The last sampled star causes a slight excess over the total cluster mass M_cl.
• 4.  
  Stars drawn with masses outside the sampled mass range are not considered.

This procedure results in the excess of the total cluster mass by a fraction of a single star (see below), though it resembles the observed IMF accurately. The independent sampling of stars with masses above and below 0.5 M_☉ prevents cluster models with a bias toward an excess of low-mass or high-mass stars. Although the stochastic sampling is performed for the whole stellar mass range, stochastic fluctuations for the most massive stars at a given age dominate the fluctuations in the integrated cluster light. Stochastic fluctuations for low-mass stars have a negligible impact due to the large total number of low-mass stars.

To test the IMF sampling, we designed an alternative sampling method: first a random number⁵ determines whether the stellar mass should be above or below the Kroupa (2001) IMF break point at 0.5 M_☉. A second random number determines the stellar mass using the appropriate relations for the mass range requested by the first random number. Using this sampling technique, clusters with a strong bias toward high-mass or low-mass stars could appear. We created 10⁴ test models with both IMF sampling methods. The resulting magnitudes and colors (i.e., GALEX satellite FUV magnitude to K-band magnitude, using solar metallicity and foreground extinction [E(B − V) = 0.0 mag, 1.0 mag and 2.0 mag] of clusters with 10³ M_☉ and 10⁴ M_☉ without USML) of both IMF sampling methods were evaluated using a K-S test. The K-S test results proved the consistency of the results derived from both IMF sampling methods (with p-values >0.1), except for the H band of 10³ M_☉ clusters. The H band is affected by the impact of stochastic sampling of both red supergiants and asymptotic giant branch (AGB) stars, which are bright but rare stars in such a low-mass cluster. 10⁴ test models might have been a too small sample to evaluate these rare cases in its full extent.

We analyzed the resulting cluster masses for anticipated masses of 10³ M_☉ (without USML), 10⁴ M_☉, and 10⁵ M_☉. The stochastic sampling process leads to a slight excess of the cluster mass, both for low-mass stars (M ⩽ 0.5 M_☉) and for high-mass stars (M > 0.5 M_☉) due to the last sampled star. The size of this excess in cluster mass is dominated by the adopted IMF and is largely independent of the anticipated cluster mass. For high-mass stars a slight tendency for increasing excess with decreasing cluster mass is found: low-mass clusters have a greater chance of being affected by stochastic draws with exceptionally large stellar masses, which cannot be balanced by a large number of further draws. This effect dominates the excess of the total cluster mass. The fraction of low-mass stars is exceeded by at least 0.12 M_☉/0.30 M_☉/0.43 M_☉ for 50%/10%/1% of all clusters. The fraction of high-mass stars is exceeded by at least 1.00 M_☉/14.92 M_☉/58.24 M_☉ for 50%/10%/1% of all clusters. Applying these fractions to the lowest modeled cluster mass of 10³ M_☉ shows that the majority of model clusters experience a mass excess ⩽0.15%, while for 1% of model clusters the mass excess reaches ⩾6% of the intended cluster mass. This fractional mass excess decreases with increasing intended cluster mass. This mass excess is caused by a single star for both mass ranges due to the stochastic sampling process. As the excess is caused by a single star whose mass is determined by the stellar IMF adopted, we do not expect a systematic impact on the modeling results.

The stochastic IMF's shape follows the underlying analytic Kroupa (2001) IMF. The stochastic realizations of the IMF deviate from the smooth fully sampled IMF, with mass ranges containing more stars and mass ranges lacking stars in comparison with the fully sampled IMF prediction. The deviations are strongest for high-mass stars and increase with decreasing M_cl, both caused by the smaller number of high-mass stars compared to their lower-mass counterparts. The present-day mass function evolves from the IMF due to stellar evolution effects, i.e., the death of stars, and dynamical effects, such as the loss of low-mass stars due to tidal interactions. The impact of dynamical cluster evolution on integrated cluster magnitudes has been studied by Anders et al. (2009). Since effects due to dynamical interactions and stochastic IMF sampling are largely independent of each other, we do not include the former and focus on the latter effect in this paper.

The individual stellar masses from the stochastically sampled IMF are used as input data for the GALEV models. The stellar evolution of these stars is taken from isochrones of Marigo et al. (2008; these isochrones do not include pre-main-sequence models), based solely on the effective temperature, surface gravity, metallicity, and scaled by the bolometric luminosity. Therefore, AGB star magnitude corrections and circumstellar dust are not required. A spectrum from the catalogs of Lejeune et al. (1997, 1998) is assigned to each star. The individual contributions of all stars are summed up to provide the integrated spectrum. Convolving the integrated spectra with a variety of filter functions provides the integrated photometry, for both the fully sampled models and the stochastically sampled models. Models are created for all ages for which isochrones are provided: log(age [yr]) range = 6.6–10.1 (i.e., age range = 4 Myr–12.6 Gyr) with steps of Δlog(age [yr]) = 0.05. In total, predictions for 71 age steps have been calculated. The isochrones for log(age) = 9.20 and 9.25 are strongly dominated by AGB stars (see Marigo et al. 2008) for metallicities −0.125 ⩽ [Fe/H] ⩽ +0.1 dex, causing a significant distortion in the color evolution. These distortions are intrinsic to the isochrones and not related to the chosen time steps (i.e., not related to binning effects of the stellar tracks). Models for the remaining ages are not affected. These ages have been removed to allow a better analysis and presentation of the results.

The time evolution of some regularly used standard colors is presented in Figures 2–5. To indicate the density of model predictions, grayscale colors have been used in logarithmic intervals in Figures 2 and 4. The black solid lines (surrounded by white shading) represent the fully sampled models for comparison. In color evolution plots for 10⁴ M_☉ clusters, four white lines show the time evolution of four individual cluster models. Figures 3 and 5 show the same distributions, based on the quantiles for each input age, representing the median value, 1σ/2σ/3σ equivalent ranges (each represented by the upper and lower quantiles, with the respective ranges in between), and the minimum and maximum value of the distributions. Depending on the cluster mass, even the median color and absolute magnitude change. These changes represent the combined effect of the varying impact of stochasticity as a function of cluster mass and to a lesser degree the finite number of calculated models. The GALEV web page provides the scripts to make similar plots for the colors chosen and downloaded. For some distributions and input ages, the median values are located in sparsely populated regions between two ridges with large numbers of models. However, the median represents the best statistical property of these distributions in combination with the quantiles.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Compared to earlier research, we provide the stochastic cluster models for public access and for a significantly extended parameter range (see Section 2.1). Compared to MASSCLEAN (Popescu & Hanson 2009), the only other available set of stochastic IMF models, we extend the provided metallicity coverage (15 metallicities compared to 2 metallicities provided by MASSCLEAN), magnitudes (184 filters instead of 13), and provide models with foreground extinction E(B − V) for a larger number of models (see Table 2) for nine cluster masses (instead of four). More detailed comparison is presented in Section 2.1.

Figures 2–5 show the main characteristics of the models.

• 1.  
  The distributions of model predictions with stochastic IMF sampling show a wide spread around the predictions with fully sampled IMF description.
• 2.  
  The distributions show substructure based on the stochastic occurrence of a small number of bright stars (e.g., supergiants, AGB stars). Individual bright stars cause certain regions of the model predictions to be far off the main distribution.
• 3.  
  We compared the following quantiles: 1σ for total cluster mass M_cl = 10³ M_☉, 2σ for M_cl = 10⁴ M_☉, and 3σ for M_cl = 10⁵ M_☉. For blue colors (such as U − B and B − V) these quantiles are comparable (with 1σ for 10³ M_☉ being slightly the largest). For red colors (such as V − I and V − K) and ages >100 Myr the 1σ for M_cl = 10³ M_☉ quantiles are more narrow than the other quantiles. The quantiles bluer than the median behave more smoothly than the quantiles redder than the median, as stochastic outliers in the stellar population are more likely to have redder colors (such as red supergiants and AGB stars).
• 4.  
  The individual models, shown for clusters with M_cl = 10⁴ M_☉, indicate strong jumps between successive age steps, tracing the evolution and death of individual stars.
• 5.  
  The spread of model predictions around the median prediction increases with decreasing cluster mass. This effect is also seen in observations (see Figure 5). The cluster masses and ages shown were derived by analyzing the observed SEDs with stochastic MASSCLEAN models (Popescu et al. 2012).
• 6.  
  Observed cluster colors, for derived cluster ages and masses, follow the distributions of our model predictions. Observed U − B colors are slightly bluer than the model predictions. However, the observed U band was not a standard Johnson U-band filter, while the model predictions are for the standard Johnson U band.
• 7.  
  Deviations of stochastic model predictions from fully sampled model predictions in different colors are not directly correlated, but follow widespread distributions as well (see color–color plots in Figure 2, right column).

Each stochastic IMF has been used to create stochastic models for 15 metallicities: Z = 0.0001, 0.0002, 0.0004, 0.0006, 0.0008, 0.001, 0.002, 0.004, 0.006, 0.008, 0.010, 0.015, 0.020, 0.025, 0.030, equivalent to [Fe/H] = −2.3, −2.0, −1.7, −1.5, −1.4, −1.3, −1.0, −0.7, −0.5, −0.4, −0.3, −0.125, 0.0, +0.1, +0.18 dex, respectively. Models are provided in 184 filters (see Table 1) for foreground extinction values E(B − V) = 0.0–2.0 mag (in steps of ΔE(B − V) = 0.2 mag). While magnitudes for the full set of filters given in Table 1 are provided for the "main model set," only magnitudes in the bold-faced filters are provided for additional models. These additional models are only available for zero foreground extinction, i.e., E(B − V) = 0.0 mag. "Main models" and "additional models" are provided for the same isochrones and metallicities, while differences occur in the number of values provided for foreground extinction and magnitudes. The relative number of models varies, depending on cluster mass (see Table 2). All models are available online at http://data.galev.org/models/anders13.

2.1. Comparison with MASSCLEAN Models

The derived integrated photometry of the stochastically sampled models was used to re-determine the physical properties of the stellar cluster (i.e., cluster age, mass, foreground extinction E(B − V), and metallicity [Fe/H]) to study the analysis accuracy of observed cluster photometry. The SED of each stellar cluster was analyzed with fully sampled IMF models, to mimic the approach regularly adopted for observed cluster photometry. For this analysis the program "AnalySED" (Anders et al. 2004) was used.

Briefly, "AnalySED" performs the following steps.

• 1.  
  The cluster SED is compared with all model SEDs.
• 2.  
  Scaling determines the cluster mass.
• 3.  
  χ² determination of the scaled cluster SED with respect to all model SEDs.
• 4.  
  Each model SED is assigned a probability based on the χ² value determined.
• 5.  
  If the total probability <10⁻³⁰ (i.e., all χ² values >70), the fit is considered as "poorly fitting."
• 6.  
  Scaling of the individual probabilities to achieve a total probability of unity.
• 7.  
  The model SED (and its associated cluster parameters) with the highest probability is considered the best fit.
• 8.  
  Models with decreasing probabilities are used to determine the parameter uncertainties, until the summed probability reaches 68.26% (1σ uncertainty range).

The ability to determine a best-fitting model and the cluster parameters of the best fit depend on various factors (see, e.g., Anders et al. 2004; de Meulenaer et al. 2013).

• 1.  
  The combination of filters used.
• 2.  
  The photometric uncertainties.
• 3.  
  The ranges of metallicity and foreground extinction explored in the fitting process.

In this section we explore the accuracy of parameter determination (i.e., cluster age, mass, foreground extinction E(B − V), and metallicity [Fe/H]) of stochastically sampled clusters. In Section 4.1 we will investigate the impact of photometric uncertainties σ and allowed ranges of metallicity and extinction for fitting dust-free input clusters with solar metallicity of 10³ M_☉, 10⁴ M_☉, and 10⁵ M_☉. The results will be used in Section 4.2 to explore the effects for a wide range of cluster parameters and filter combinations.

4.1. Impact of Fitting Parameters

We first study the special case of dust-free input clusters with solar metallicity ([Fe/H] = 0 and E(B − V) = 0) to illustrate the method and highlight findings applicable to stochastic clusters in general. A more general study of clusters with non-solar metallicity and finite foreground extinction is discussed below in Section 4.2.

First, we investigate which fraction of stochastic models can be fitted at all with the fully-sampled models. A summary is provided in Table 3. We investigated three sets of allowed ranges in metallicity and foreground extinction.

• 1.  
  Range1: Restricting to the input values ([Fe/H] = 0.0 dex and E(B − V) = 0.0 mag).
• 2.  
  Range2: Allowing both extinction and metallicity to vary over a narrow range around the input values (−0.125 ⩽ [Fe/H] ⩽ +0.1 dex and 0.0 ⩽ E(B − V) ⩽ 0.2 mag).
• 3.  
  Range3: Allowing both extinction and metallicity to vary over the full range around the input values (−2.3 ⩽ [Fe/H] ⩽ +0.18 dex and 0.0 ⩽ E(B − V) ⩽ 2.0 mag).

Table 3 indicates three trends.

• 1.  
  The fraction of poorly fitted clusters increases with decreasing cluster mass. The lower the cluster mass is, the stronger the impact of stochastic IMF sampling, relative to the total stellar population in a stellar cluster (with related enhanced changes in the cluster SED). This effect is well known and described, e.g., in Popescu & Hanson (2009) and Silva-Villa & Larsen (2011).
• 2.  
  For a given cluster mass, the allowed fitting ranges in metallicity and extinction strongly influence the number of poorly fitted clusters: the smaller the allowed ranges around the input values are, the higher the number of poorly fitted clusters. The age range with the highest fraction of fitted clusters depends on input cluster mass and fitting constraints. These results suggest large deviations between the original cluster parameters and the fit results, which will be studied below.
• 3.  
  For clusters with a total mass of 10⁴ M_☉ the impact of the assumed photometric uncertainties was studied: The smaller the photometric uncertainty, the more significantly the derived χ² value is affected by SED mismatches caused by the stochastic IMF sampling. This effect leads to an enhanced number of poorly fitted clusters.

Therefore, we investigate the spread of derived cluster parameters around the input values to quantify this effect. The spread is given based on quantiles. For every model the difference between derived parameter and input parameter is determined (for ages and masses the logarithm of the quantity is used). We consider the parameter difference between the 84.135% and the 15.865% quantile as the "quantile 1σ" measure. The median and "quantile 1σ" values for the rederived cluster ages, masses, and foreground extinctions E(B − V) are presented in Table 4 for various cluster masses and fitting range restrictions. In Table 5 the effect of assigned photometric uncertainties is investigated. Median values >0.0 indicate the rederived parameter being larger than the input value.

Table 3. Fraction of Models (in %) Found Poorly Fitted during the AnalySED Fitting, Using UBVRIJHK Magnitudes

Cluster Mass	Fitting Ranges	Age ⩽ 10 Myr	10 Myr < Age ⩽ 100 Myr	100 Myr < Age ⩽ 1 Gyr	Age > 1 Gyr
(M☉)
105	Range1	0.0	0.0	0.0	0.0
105	Range2	0.0	0.0	0.0	0.0
105	Range3	0.0	0.0	0.0	0.0
104	Range1	0.0	0.0	15.1	23.4
104, 0.001 mag	Range2	99.9	99.9	99.9	99.9
104, 0.01 mag	Range2	74.5	62.5	83.3	89.3
104, 0.1 mag	Range2	0.0	0.0	12.0	19.4
104, 0.2 mag	Range2	0.0	0.0	6.0	5.6
104	Range3	0.0	0.0	0.0	0.0
103	Range1	1.6	1.3	20.0	54.8
103	Range2	0.9	1.0	16.2	44.9
103	Range3	0.0	0.0	2.5	0.2

Notes. The stochastic clusters were created with [Fe/H] = 0.0 dex, foreground extinction E(B − V) = 0.0 mag, and an assumed photometric uncertainty of Δmag = 0.1 mag for all magnitudes (unless specified). The percentages in Columns 3–6 represent the fractions of poorly fitted models with input ages in four age ranges. The numbers of models in these age ranges are 12.7%, 28.2%, 28.2%, and 30.9% of the total number of models. Three different fitting ranges were investigated: Range1: [Fe/H] = 0.0 dex and E(B − V) = 0.0 mag (i.e., the true model values); Range2: [Fe/H] = −0.125–+0.1 dex and E(B − V) = 0.0–0.2 mag (i.e., small ranges around the true model values); and Range3: [Fe/H] = −2.3–+0.18 dex and E(B − V) = 0.0–2.0 mag (i.e., wide range in metallicity and extinction).

Download table as:  ASCIITypeset image

Table 4. Cluster Parameter Determination Accuracy Based on AnalySED Fitting, Using UBVRIJHK Magnitudes and Fully Sampled Fitting Models

Cluster Mass	Fitting Ranges	Parameter	Age ⩽ 10 Myr	10 Myr < Age ⩽ 100 Myr	100 Myr < Age ⩽ 1 Gyr	Age > 1 Gyr
(M☉)
105	Range1	Age	0.0 ± 0.10	0.0 ± 0.25	0.0 ± 0.20	0.0 ± 0.20
105	Range2	Age	0.0 ± 0.10	0.0 ± 0.15	0.05 ± 0.15	0.0 ± 0.30
105	Range3	Age	0.0 ± 0.50	−0.05 ± 0.20	−0.05 ± 0.20	0.0 ± 0.45
104	Range1	Age	0.0 ± 0.80	0.0 ± 1.00	−0.1 ± 0.70	−0.15 ± 0.65
104	Range2	Age	0.0 ± 0.30	−0.05 ± 0.80	−0.1 ± 0.65	0.0 ± 0.40
104	Range3	Age	0.15 ± 0.45	−0.05 ± 0.35	−0.1 ± 0.70	0.0 ± 1.15
103	Range1	Age	−0.1 ± 0.35	−0.55 ± 2.20	−0.6 ± 1.15	−0.3 ± 0.95
103	Range2	Age	−0.1 ± 0.35	−0.6 ± 2.15	−0.6 ± 1.20	−0.05 ± 1.05
103	Range3	Age	0.15 ± 0.35	−0.5 ± 0.80	−0.05 ± 2.20	0.15 ± 1.10
105	Range1	Mass	0.0 ± 0.13	0.0 ± 0.21	0.0 ± 0.18	0.0 ± 0.26
105	Range2	Mass	0.03 ± 0.13	0.0 ± 0.16	0.0 ± 0.14	−0.01 ± 0.22
105	Range3	Mass	0.09 ± 0.72	−0.01 ± 0.19	0.0 ± 0.10	−0.01 ± 0.24
104	Range1	Mass	0.04 ± 1.18	−0.02 ± 0.96	−0.08 ± 0.60	−0.15 ± 0.76
104	Range2	Mass	0.05 ± 0.97	−0.03 ± 0.81	−0.07 ± 0.59	−0.05 ± 0.53
104	Range3	Mass	0.27 ± 0.84	0.0 ± 0.50	0.0 ± 0.50	0.04 ± 0.56
103	Range1	Mass	−0.39 ± 1.34	−0.75 ± 2.08	−0.49 ± 1.19	−0.33 ± 1.18
103	Range2	Mass	−0.29 ± 1.24	−0.74 ± 2.09	−0.48 ± 1.22	−0.19 ± 1.06
103	Range3	Mass	0.22 ± 0.57	−0.63 ± 1.13	−0.28 ± 1.38	−0.09 ± 1.69
105	Range2	E(B − V)	0.1 ± 0.20	0.1 ± 0.20	0.1 ± 0.20	0.1 ± 0.20
105	Range3	E(B − V)	0.1 ± 0.25	0.1 ± 0.20	0.15 ± 0.25	0.15 ± 0.30
104	Range2	E(B − V)	0.05 ± 0.20	0.1 ± 0.20	0.2 ± 0.20	0.0 ± 0.20
104	Range3	E(B − V)	0.1 ± 0.40	0.2 ± 0.40	0.3 ± 0.40	0.3 ± 0.70
103	Range2	E(B − V)	0.0 ± 0.15	0.15 ± 0.20	0.2 ± 0.20	0.0 ± 0.15
103	Range3	E(B − V)	0.05 ± 0.10	0.1 ± 0.50	0.35 ± 0.35	0.2 ± 0.80

Notes. The stochastic clusters were created with [Fe/H] = 0.0 dex, foreground extinction E(B − V) = 0.0 mag, and an assumed photometric uncertainty of Δmag = 0.1 mag for all magnitudes. The values in Columns 3–6 represent the median and quantile 1σ of fitted log(ages) with input log(ages), fitted log(mass) with input log(mass), and fitted E(B − V) with input E(B − V) in four input age ranges. Fitting ranges as in Table 3.

Download table as:  ASCIITypeset image

Table 5. Cluster Parameter Determination Accuracy Based on AnalySED Fitting, Using a Variety of Photometric Uncertainties

σ	Parameter	Age ⩽ 10 Myr	10 Myr < Age ⩽ 100 Myr	100 Myr < Age ⩽ 1 Gyr	Age > 1 Gyr
0.001 mag	Age	0.0 ± 0.05	0.0 ± 0.10	0.0 ± 0.00	0.0 ± 0.00
0.01 mag	Age	0.0 ± 0.15	−0.05 ± 0.35	−0.05 ± 0.10	0.05 ± 0.35
0.1 mag	Age	0.0 ± 0.30	−0.05 ± 0.80	−0.1 ± 0.65	0.0 ± 0.40
0.2 mag	Age	0.0 ± 0.30	−0.05 ± 0.85	0.0 ± 0.70	0.0 ± 0.45
0.001 mag	Mass	0.0 ± 0.05	0.01 ± 0.13	−0.01 ± 0.05	0.05 ± 0.06
0.01 mag	Mass	0.0 ± 0.30	−0.02 ± 0.37	−0.01 ± 0.13	0.08 ± 0.31
0.1 mag	Mass	0.05 ± 0.97	−0.03 ± 0.81	−0.07 ± 0.59	−0.05 ± 0.53
0.2 mag	Mass	0.02 ± 0.42	−0.07 ± 0.94	−0.05 ± 0.76	−0.03 ± 1.17
0.001 mag	E(B − V)	0.0 ± 0.05	0.05 ± 0.15	0.05 ± 0.10	0.1 ± 0.15
0.01 mag	E(B − V)	0.1 ± 0.15	0.15 ± 0.20	0.15 ± 0.20	0.15 ± 0.20
0.1 mag	E(B − V)	0.05 ± 0.20	0.1 ± 0.20	0.2 ± 0.20	0.0 ± 0.20
0.2 mag	E(B − V)	0.0 ± 0.20	0.15 ± 0.20	0.2 ± 0.20	0.05 ± 0.20

Notes. The stochastic clusters were created with cluster masses of 10⁴ M_☉, [Fe/H] = 0.0 dex, foreground extinction E(B − V) = 0.0 mag, and a range of photometric uncertainties σ (for all magnitudes). The values in Columns 3–6 represent the median and quantile 1σ of fitted log(ages) with input log(ages), fitted log(mass) with input log(mass), and fitted E(B − V) with input E(B − V) in four age ranges. Fitting ranges are narrow around input values ([Fe/H] = −0.125–+0.1 dex and E(B − V) = 0.0–0.2 mag).

Download table as:  ASCIITypeset image

The following general trends are observed, based on the SEDs of stochastically sampled clusters deviating from the model cluster SEDs using a fully sampled IMF.

• 1.  
  Fitting results are more widely distributed for decreasing cluster mass, and all cluster parameters are affected.
• 2.  
  The deviations do not only correspond to a wide spread of fitting results (i.e., "quantile 1σ" values), but the median of the fitting results also deviates from the input values.
• 3.  
  When fitting parameter ranges are extended, the additional parameter space (i.e., increased ranges for metallicity and foreground extinction) is significantly used for the fitting. However, regarding derived ages and masses, the differences between the three sets of fitting constraints are not very strong and do not appear like a homogeneous trend. Even knowing and using the correct values for metallicity and extinction does not reduce the uncertainty ranges significantly; the stochasticity is the dominating effect.
• 4.  
  Increasing photometric uncertainties allow a larger number of stochastic clusters to be fitted with fully sampled models, although with decreasing fitting accuracy.

We tested the corresponding derived spread when quantified using the Gaussian standard deviation 1σ. In general, the derived Gaussian standard deviations were smaller than the equivalent quantile 1σ. The "quantile 1σ" measure can be applied to much more general cases: Since the derived parameter distributions we studied are significantly non-Gaussian, using a Gaussian standard deviation can lead to misleading results. The "quantile 1σ" can be generally applied to all distributions.

In Figure 6 we show the rederived ages (top row) and masses (bottom row) for a range of cluster masses. All clusters have photometry in the standard Johnson/Cousins/Bessel–Brett system with magnitudes in the UBVRIJHK bands, input values of [Fe/H] = 0.0, and E(B − V) = 0.0, fitted with the narrow parameter range ([Fe/H] = −0.125–+0.1 dex, E(B − V) = 0.0–0.2 mag). From left to right, the clusters have initial masses of 10³ M_☉, 10⁴ M_☉, and 10⁵ M_☉. From the stochastic models and the analytic models two ages have been removed: log(age) = 9.2 and log(age) = 9.25 since for these ages the spectra are dominated by a few very bright AGB stars. For these ages all stochastic models roughly coincide for any color, in contrast to all other ages. In addition, the color evolutions show a strong discontinuity for these ages.

Zoom In Zoom Out Reset image size

Earlier research by Silva-Villa & Larsen (2011) found the following impact of stochastic IMF sampling on derived cluster parameters.

• 1.  
  For low-mass clusters with masses ⩽10³ M_☉, they find a wide distribution of wrongly derived ages for stochastic clusters (up to a factor of 100 between input and derived ages).
• 2.  
  For clusters with masses >10³ M_☉, they find a spread of derived ages around the input ages, with the spread decreasing with increasing cluster mass.
• 3.  
  They do not perform stochastic tests to quantify the spread of derived parameters.

These earlier results qualitatively coincide with our results: derived ages and masses show wide distributions around the input values. This spread becomes smaller with increasing cluster mass. However, even for clusters with initial mass = 10⁵ M_☉ we find quantile 1σ uncertainties for log(age) in the range of 0.1–0.3 dex and for log(mass) = 0.13–0.22 dex, depending on input age. Systematic offsets are visible for young ages (e.g., derived ages in the range of 100 Myr–10 Gyr for input ages <30 Myr), i.e., when few bright stars are at the upper main sequence or in rapidly evolving states. For clusters with initial masses of 10³ M_☉, a tendency to underestimate the metallicity is observed, which becomes weaker for clusters with initial masses of 10⁴ M_☉. The impact of input cluster metallicity and foreground extinction on the fitting results will be investigated in the next section. The derived extinction is almost uniformly distributed over the full available range.

4.2. Cluster Parameters and Filter Combinations

Stochastic effects to sample the stellar IMF realistically have not yet been sufficiently taken into account for population synthesis modeling. As we have shown in this paper, the impact on color predictions and parameter determinations using stellar cluster models with a fully sampled IMF is highly significant. The strongest impact on model predictions is seen for parameter ranges and magnitudes/colors, dominated by a small number of bright stars. For young clusters (age ⩽ 15 Myr) such stars are supergiant stars, while for older clusters (100 Myr ⩽ age ⩽ 1 Gyr) the light is dominated by AGB stars. In both cases, the strongest impact is seen in the near-infrared. In addition, the cumulative effects of deviations of the stochastic IMF from the underlying fully sampled IMF cause more general deviations.

The effect on cluster parameter determination is a clear function of cluster mass and parameter fitting techniques: the impact of stochastic IMF sampling increases with a decreasing number of stars dominating the light. In general, this is related to the total cluster mass and the choice of filters used, as filters are affected by varying amounts (e.g., depending on the stellar population and cluster age).

Due to the large number of models provided, classical techniques of fitting observations (such as determining the χ² of all model SEDs compared with the SED of an observed cluster) become unfeasible. For an observed cluster, trying to estimate the cluster parameters that would give rise to such an SED is a challenging problem for a number of reasons. The data do not follow a well-defined distribution, and there may be many modes in the range of the parameters. Our future work will try to tackle this issue by applying the method of Priebe & Marquette (2000). They described a semi-parametric approach to data that builds on the advantages of both mixture model analysis and kernel density estimation. We propose restricting the range of possible parameters through modeling the distributions using this technique, so that an observer could better estimate the range of parameters that would lead to such a cluster.

Our first application to stellar cluster observations are the HST observations of NGC 5253. Analysis of these clusters using model SEDs with a fully sampled IMF has been presented in de Grijs et al. (2013), which indicated cluster masses prone to stochastic IMF sampling. Therefore, we will re-analyze the observational data using the stochastic models presented in this paper (R. de Grijs et al., in preparation).

All model predictions for the range of cluster masses, ages, metallicities [Fe/H], foreground extinction E(B − V), and filters presented are available at http://data.galev.org/models/anders13. The web page also contains the scripts to produce plots of density distributions, as well as quantiles for retrieved magnitudes and colors. Such plots are essential to judge the reliability of cluster analyses for the parameter range and the colors the user chooses.

P.A. and R.d.G. acknowledge funding by the National Natural Science Foundation of China (NSFC, grant No. 11073001). P.A. acknowledges funding by the China Postdoctoral Science Foundation. R.K. gratefully acknowledges financial support from STScI theory grant HST-AR-12840.01-A. Support for Program number HST-AR-12840.01-A was provided by NASA through a grant from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Incorporated, under NASA contract NAS5-26555. P.A. wishes to thank S. Gao and J.-Z. Zhao for technical support and N. Bissantz and E. Silva-Villa for helpful discussions.