5 Detailed comparison of WR populations with synthesis models

   
5.1 Procedure

To interpret quantitatively the observational data we use evolutionary synthesis models and proceed essentially as in SGIT00. The following main observational constraints are used:

1.

${\rm H}\beta$ and ${\rm H}\alpha$ equivalent widths. The former is used as a primary age indicator; once $W({\rm H}\beta )$ is reproduced $W({\rm H}\alpha)$ may serve as an independent consistency test for the predicted spectral energy distribution (SED) in the red (cf. SGIT00).

2.

Nebular line intensities. $F({\rm H}\alpha)/F$( ${\rm H}\beta$) determines the extinction of the gas. The use of other line intensities requires detailed photoionization modeling which is beyond the scope of this paper.

3.

Intensities and equivalent widths of the main WR features. The blue bump and C  IV $\lambda$5808 (red bump) serve as main constraints on the WR population. To avoid uncertainties in deblending individual contributions of the blue bump we prefer to use measurements for the entire bump. In contrast to the spectra of metallicity objects our spectra show no evident contamination from nebular lines (e.g. [Fe  III] $\lambda$ 4658, nebular He  II).

To potentially disentangle between various effects (underlying "non-ionizing'' population, loss of photons, differential extinction between gas and stars) it is important to use both equivalent widths and relative $I({\rm WR})/I({\rm H}\beta)$ intensities (cf. Schaerer et al. 1999a).

For the model comparisons we use calculations based on the evolutionary synthesis code of SV98, which in particular includes the most recent calibration of WR line luminosities used to synthesize the WR features, up-to-date stellar tracks, CoStar stellar atmospheres for O stars, pure H-He models for WR stars and Kurucz models for cooler stars (see SV98 for a full description). Except for the improved O star atmospheres used by SV98 the Starburst99 synthesis models (Leitherer et al. 1999) use the same basic input physics. New generation stellar atmosphere models for O and WR stars including a full treatment of non-LTE line blanketing and stellar winds have just now become available for the use in synthesis models (Smith et al. 2002). However, as the quantities of interest here depend only on the total number of Lyman continuum photons which is not altered, the use of these more sophisticated atmosphere models does not affect our results.

It is important to stress that in all cases the high-mass loss stellar tracks of Meynet et al. (1994) are used. It is thought that this adjustment of mass-loss, treated like a free parameter, will become ultimately obsolete when a proper treatment of the various effects of stellar rotation is made in the stellar evolution models. First results tend to indicate that this may indeed be the case (Meynet 1999). The Meynet et al. (1994) tracks are chosen as they reproduce a large number of properties of individual WR stars and WR populations (including especially relative WR/O ratios for a standard Salpeter IMF) in Local Group galaxies (Maeder & Meynet 1994). The use of other tracks (e.g. the "normal'' mass loss tracks) which are known to disagree with these basic constraints on WR and O star populations, would imply a strong inconsistency with the Local Group data.

Figure 8: Observed and predicted equivalent width (left panel) and line intensity with respect to ${\rm H}\beta$ (right panel) as a function of $W({\rm H}\beta )$. Our VLT sample is shown by (black) triangles, the BK02 sample with (red) squares. Typical uncertainties are 5-10% for $W({\rm H}\beta )$, $\le$10% for W(WR bump), and $\sim$0.05 dex in $\log(I($WR)/ $I({\rm H}\beta )$). Model predictions are shown for instantaneous bursts with "standard'' IMFs at Z=0.008 (dashed line), $Z=Z_{\odot }=0.02$ (solid line), and Z=0.04 (long dashed line). Note the overprediction of the WR bump strength in high metallicity models compared with the observations.

The basic model parameters we consider are:

a)

Metallicity. Stellar tracks covering metallicities Z= 0.008, 0.02 (solar), and 0.04.

b)

IMF slope and upper mass cut-off ( $M_{\rm up}$). We adopt a Salpeter IMF (slope $\alpha =2.35$), and $M_{\rm up}=120$ $M_{\odot }$ as our standard model.

c)

Star formation history (SFH). Models for instantaneous bursts (coeval population), extended burst durations (constant SF during period $\Delta t$; in this case age = 0 is defined at the onset of SF, i.e. corresponds to that of the oldest stars present), and constant SF are considered.

d)

Fraction of ionizing Lyman continuum photons ($f_\gamma$). $f_\gamma$ indicates the fraction of ionizing photons absorbed by the gas. Our standard value is $f_\gamma=1$. Values $f_\gamma < 1$ are used to simulate various effects (e.g. dust absorption, photon leakage outside regions, etc.) leading to a reduction of photons available for photoionization.

Unless stated otherwise our models are calculated assuming an IMF fully sampled over the entire mass range (as in SV98). For some cases we have also done model calculations based on a Monte Carlo sampling of the IMF, in order to quantify the effects of statistical fluctuations due to the finite number of massive stars. We have verified our calculations by comparison with the Monte Carlo models and analytical results of Cerviño et al. (2000, 2002).

Figure 9: Observed and predicted WR bump equivalent width as a function of $W({\rm H}\beta )$. Standard model predictions are shown for instantaneous bursts at $Z=Z_{\odot }=0.02$ (solid line), and Z=0.04 (long dashed line). Left panel: Thick (green) lines with the same styles show models with a standard IMF slope ( $\alpha =2.35$) and upper mass cut-offs of $M_{\rm up}=60$ and 30 $M_{\odot }$ delimiting the singly and doubly shaded regions respectively. Right panel: Thick (green) lines with the same styles show models with a IMF slope of $\alpha =3.3$ and $M_{\rm up}=120$.

   
5.2 Results

A comparison of the observed equivalent widths and relative intensity of the WR bump with standard model predictions at different metallicities is presented in Fig. 8. The following points can be seen from this figure:

• The observed trends of decreasing (increasing) equivalent width (line intensities) of the WR features with decreasing ${\rm H}\beta$ equivalent width agree with the instantaneous burst model predictions over the same $W({\rm H}\beta )$ range. Though shorter, part of the initial phase of increasing (decreasing) $W({\rm WR})$ ( $I({\rm WR})/I({\rm H}\beta)$) corresponding to the onset of the WR rich phase does not seem to be detected.

• Most importantly, essentially all the observed WR features are weaker than the predictions of our standard models for instantaneous bursts at solar or higher metallicity. Very similar observational trends are also found in the metal-rich H  II region samle of BK02 and the metal-rich starbursts of Schaerer et al. (2000). This is in stark contrast with observations at lower metallicity where generally an good agreement is found with short burst models (e.g. Schaerer et al. 1999a; Guseva et al. 2000).

• The above result is found independently from the observed WR equivalent widths and relative line intensities WR bump/ ${\rm H}\beta$. This indicates that the discrepancy between models and observations is not related to the possible presence of an underlying older stellar population, which would dilute (reduce) equivalent widths but not alter the relative line intensities.

The following possibilities (one or a combination thereof) could be invoked to explain the discrepancy between our observations and models:

1.

The metallicities of our HII regions are overestimated. Indeed the observations could be reconciled with burst models with a "standard'' IMF for metallicities $Z \sim (1/2.5$-1) $Z_\odot$, as shown in the left panel of Fig. 8 (short dashed line). However, despite the uncertainties in the O/H determinations (cf. Sect. 3) such low average metallicities seem very implausible.

2.

Extended bursts. Such a scenario has been invoked by SGIT00 for the sample of metal-rich starbursts based on the finding of red supergiant features in their spectra and the fact that these distant objects are mostly nuclear starbursts observed through apertures corresponding to relatively large spatial scales. In this case all observed properties could quite well be fitted with "standard'' solar metallicity models for burst durations $\Delta t \sim 4$-10 Myr. However, in view of the different nature (disk H  II regions) of the present sample, indications of relatively short formation time scales of H  II regions (e.g. Massey et al. 1995), and the lack of direct signatures of older/red supergiant populations (cf. below) it seems quite unjustified to appeal to extended burst to solve the observed discrepancy.

3.

A modified IMF (upper mass cut-off and/or slope). In a plot like Fig. 8, a Salpeter IMF with a lower upper mass cut-off simply implies that the curve plotted here (for $M_{\rm up}=120$ $M_{\odot }$) is joined at lower $W({\rm H}\beta )$ as the WR stars from the most massive stars are absent. This is illustrated for the cases of $M_{\rm up}$ = 30 and 60 $M_{\odot }$ by the shaded domains in Fig. 9. The shape of the predicted WR equivalent width or line intensity remains, however, unchanged. Therefore the observed discrepancy cannot be resolved with an IMF of Salpeter slope and a lower value of $M_{\rm up}$ (see also Sect. 6).

Models with steeper, variable IMF slopes ( $2.35 < \alpha \la 3.3$) and $M_{\rm up}\sim60$-120 $M_{\odot }$ could reproduce most of the objects, with the exception of the lowest $W({\rm H}\beta )$ objects (see Fig. 9). As the least metal-poor objects in our sample are probably of similar nature as young clusters or H  II regions in our Galaxy whose stellar content has been studied in detail, we may presume that their IMF (slope and $M_{\rm up}$) should be similar. Since none of the Galactic regions have shown convincing evidence of a strong deviation of the IMF slope from the Salpeter value (see Massey 1998 and references therein), we think that such a steeper slope is an unlikely explanation.

4.

Incorrect stellar evolution models and/or "calibration data'' Although the adopted tracks (Meynet et al. 1994) compare fairly well with various observations, several failures of the non-rotating stellar models are also known (see e.g. Maeder 1999). However, the used tracks have essentially been calibrated/adjusted to fit the observed WR/O ratio in various regions of our Galaxy and Local Group objects which are though to be at equilibrium, i.e. showing relative populations corresponding to constant star formation (see compilation in Maeder & Meynet 1994). The relative WR/O star ratio is the one most directly related to our (time resolved) observables. As this calibration yields a fairly good agreement over a large metallicity range ( $1/10 \la Z/Z_{\odot}\ \la 2$) there seems little room for changes in the tracks which could reduce the predicted WR bump by the required factor of $\sim$2 without violating the WR/O constraints in the Local Group.

One could argue that the calibration data, the observed WR/O number ratio at solar metallicity and above could be incorrect due to possible incompleteness or biases in the stellar counts (see e.g. related discussions in Massey & Johnson 1998). However, to reconcile our WR observations in H  II regions with the corresponding counts for our Galaxy and M31 would require a downward revision of the relative WR/O ratio by up to a factor of 2, which seems highly unlikely.

5.

Uncertainties in synthesis of the WR bump. Presently the calculation of observables related to WR stars is simply done in the following way in evolutionary synthesis models. The different emission line strengths are computed by multiplying the predicted number of WR stars (grouped in different types and/or subtypes) with their average line luminosity as derived from observations of a sample of WR stars (see SV98). Interestingly the intrinsic line luminosity of the strongest line of the WR bump, He  II $\lambda$4686, shows a rather large scatter, namely $L_{\rm 4686}=(1.6 \pm 1.5) \times 10^{36}$ erg s^-1 in the Galactic and LMC WNL calibration sample of SV98 with a possible increase of $L_{\rm 4686}$ with the stellar bolometric luminosity L (see Fig. 1 of SV98). Such a luminosity dependence of $L_{\rm 4686}$ with L could in fact (partly or fully) explain the observed discrepancy as we will now show.

Splitting the WNL calibration in two domains with luminosities above/below $\log L/L_{\odot}= 6$, SV98 found average line luminosities $L_{\rm 4686} = 5.6 \times 10^{35}$ ( $\log L < 6$) and $L_{\rm 4686} = 3.1 \times 10^{36}$ ( $\log L > 6$). Replacing in the synthesis models the overall average for WNL stars by these quantities leads to an important reduction of $W({\rm WR})$ in solar metallicity bursts with $W({\rm H}\beta )$ $\la 60$-70 Å, as shown in Fig. 10. At larger ${\rm H}\beta$ equivalent widths (corresponding to ages $\la$4-5 Myr for Z=0.02) the WR bump predictions are less modified, since a) WC stars contribute more importantly to the bump and b) only the youngest bursts with very high $W({\rm H}\beta )$ are dominated by very luminous WNL stars showing thus larger $W({\rm WR})$.

Figure 10: Observed and predicted WR bump equivalent width as a function of $W({\rm H}\beta )$. Standard model predictions are shown for instantaneous bursts at $Z=Z_{\odot }=0.02$ (solid line), and Z=0.04 (long dashed line). Thick (green) lines with the same styles show the predictions with the modified $L_{\rm 4686}$ calibration for WNL stars leading to an important reduction of the WR bump for $\log W({\rm H}\beta)\protect\la 1.8$, due to the lower average luminosity of WNL stars in bursts with ages $\protect\ga$4-5 Myr (for Z= 0.02).

The last option (5) seems the most likely explanation to explain the surprisingly low WR equivalent widths and intensities in our sample of metal-rich H  II regions. Implications on earlier studies of WR galaxies are briefly discussed in Sect. 5.3.

In contrast, the following hypothesis or effects altering observed equivalent widths and/or relative line intensities cannot be the cause of the discrepancy:

• Underlying (old) populations diluting the $W({\rm WR})$ measurements. The line intensities are unaffected by underlying populations and show the same discrepancy as $W({\rm WR})$ (cf. above). Furthermore, inspection of our spectra show little or no indication of signatures of older stellar populations.
• Differential extinction between gas and stars as frequently observed in H  II galaxies and starbursts (cf. Fanelli et al. 1988; Calzetti 1997; Mas-Hesse & Kunth 1999, SGIT00). If present such an effect alters $W({\rm H}\beta )$ and $I({\rm WR})/I({\rm H}\beta)$ but not $W({\rm WR})$. To bring the observations to agreement with standard models would require that the gas suffered a lower extinction than the stars (implying thus lower corrected $W({\rm H}\beta )$ and larger $I({\rm WR})/I({\rm H}\beta)$), opposite to what is found empirically.
• Leakage of ionising photons outside the observed regions, dust absorption of ionising photons, or other mechanisms reducing the fraction $f_\gamma$ of Lyman continuum photons used for photoionisation. Correcting the observations for such an effect would increase the observed $W({\rm H}\beta )$ and decrease $I({\rm WR})/I({\rm H}\beta)$, which does not reconcile observations and theory (see Fig. 8).
• Effects due to varying seeing conditions and limited slit widths could lead to a loss of nebular emission in the aperture or a fraction of the stellar light. The former was just discussed ("leakage''). If the WR distribution were systematically more extended than the other stars responsible of the continuum, this effect could lead to reduced WR bump equivalent widths. No such trend between $W({\rm WR})$ and the seeing conditions is found.

• Stochastical fluctuations of the IMF due to small number statistics of the massive stars, as invoked by Bresolin & Kennicutt (2002). Although indeed expected to introduce some scatter (see Cerviño et al. 2000, 2002) this effect cannot explain the discrepancy for several reasons. First, the sample discussed here is sufficiently large and shows clear observational trends with relatively small scatter. In addition, no systematic trend of $W({\rm WR})$ (and $I({\rm WR})/I({\rm H}\beta)$) with absolute scale, such as measured by the continuum luminosity or ${\rm H}\alpha$ luminosity, is observed in our sample as shown in Fig. 11. Second, Monte Carlo simulations we have undertaken for cluster sizes corresponding roughly to the observed average continuum luminosity of $<\log L_{\rm cont}> \sim \!36.8 \pm 0.4$ predict a typical relative scatter of $\sigma(W({\rm WR bump}))/W({\rm WR})\sim 25$% (cf. Cerviño et al. 2002) - too small to account for the discrepancy - and no significant bias towards lower values as illustrated in Fig. 12.
• The use of different stellar tracks, such as e.g. the Geneva tracks with standard mass loss tracks adopted by Bresolin & Kennicutt (2002) for the Monte Carlo models, which do not reproduce basic constraints of massive stars and stellar populations in Local Group objects cannot be a solution as they would lead to important inconsistencies (cf. discussion in Sect. 5.1).

Figure 11: Equivalent widths of the WR bump as a function of the monochromatic continuum luminosity at ${\rm H}\beta$ in erg s-1 Å-1 (top panel), and as a function of the ${\rm H}\alpha$ luminosity for the WR regions of our sample. The mean and dispersion (1 $\sigma$) of $\log L_{\rm cont}({\rm H}\beta)$ is plotted in the top panel.

Figure 12: Same as Figs. 8 (left) and 10 showing the comparison between model predictions using a fully sampled, analytical IMF (dotted line, black) and the predicted mean $W({\rm WR})$ and $1 \sigma$ scatter (solid line, blue) for a burst of scale/mass corresponding to the average observed continuum luminosity. Observations are shown using the same symbols and in the previous figures. Note, the deviation of the mean values for $\log W({\rm H}\beta)\sim 1.9$-2.2 is due to a numerical artifact. The comparison shows that no significant bias is expected and that the scatter is too small to resolve the discrepancy with observations.

   
5.3 Discussion

In Sect. 5.2 we have argued that, compared to the normal prescription used in our SV98 synthesis models, a different prescription should preferrably be adopted to predict more accurately the He  II $\lambda$4686 emission from WN stars. As several earlier studies including ours (e.g. Schaerer 1996, 1999; Schaerer et al. 1999a; Guseva et al. 2000, SGIT00) are based on the use of the simple average He  II $\lambda$4686 line luminosity of SV98 for WNL stars, it is important to assess if or to what extent the use of a luminosity dependent prescription would affect the results from previous studies.

To verify this we have recomputed several sets of models for sub-solar metallicities. The maxima of the WR bump intensity and $W({\rm WR})$ (cf. Fig. 6) are only slightly modified (increased at $12 + \log({\rm O/H})\la 8.5$, and decreased above) and lead to a somewhat smaller increase with O/H, improving the agreement with the observations. For metallicities $Z \la$ 1/2 $Z_\odot$ the predicted WR bump is found to be larger at all ages (as the bulk of WN stars are of high luminosity), whereas for higher metallicities both larger/smaller WR bump strengths are predicted depending on the burst age ( $W({\rm H}\beta )$), as for the cases shown in Fig. 10. These changes improve the comparison with observations at low Z (see e.g. Fig. 7 of Guseva et al. 2000). No clear statement can be made for intermediate metallicities. A better understanding of the dependence of the WR emission lines on the stellar parameters appears necessary to improve the accuracy of the predictions of WR features in evolutionary synthesis models. The impact of newly available stellar evolution models including the effects of rotation on interior mixing and mass loss on massive star populations remains also to be explored.