© ESO, 2013

1. Introduction

Circumstellar disks are the natural byproducts of the star formation process that are expected to be the birthplace of planets. While their mass is thought to be dominated by gas with very little amount of dust (1% in mass), it is the dust that dominates the continuum opacity of a circumstellar disk. As the disks reprocess stellar light, they can be characterized best at near- to far-infrared (IR) wavelengths up to millimeter wavelengths. Based on the excess emission in the near-IR, the typical time scale for the disk to dissipate is estimated to be around ~5−6 Myr (e.g., Hernández et al. 2007). Within these few Myr, dust grains are expected to grow from μm-sized pebbles to km-sized planetesimals and eventually to planets before the disk’s material is efficiently dissipated via photo-evaporation, or via grain growth and subsequent radial drift (see Williams & Cieza 2011 for a recent review).

An important question aims at addressing how the planet formation efficiency changes as the function of the environment, especially in the low-mass regime. Numerical simulations of planetary population synthesis, such as the ones led by Alibert et al. (2011) and Mordasini et al. (2012) have shown that the expected planet population should depend on the stellar masses and disk properties. There has been a general effort to better characterize the relationship between the stellar masses and disk properties, which can additionally be related to the formation mechanism of very low-mass stars (Padoan & Nordlund 2002; Reipurth & Clarke 2001; Bayo et al. 2011, 2012). It is therefore of prime importance for better constraining the disk characteristics, as well as dust and gas properties, around stars over a wide range of masses (e.g., Pascucci et al. 2009). Earlier studies based on near-IR Spitzer/IRAC and MIPS 24 μm data (e.g., Szűcs et al. 2010; Mulders & Dominik 2012) investigated the vertical structure of disks as a function of the stellar mass. Szűcs et al. (2010) show that gas disks in hydrostatic equlibrium with perfect mixing of dust and gas are inconsistent with the observations of low-mass (early M-type stars) and very low-mass stars (late M dwarfs and brown dwarfs, hereafter BDs). They find that larger scale height reduction of the dust component is required for BDs than for low-mass stars when compared to the hydrostatic equilibrium solution. Mulders & Dominik (2012) implemented the prescription of dust settling (described in Dullemond & Dominik 2004) to hydrostatic equilibrium models and found similar turbulent mixing strengths (a parameter governing dust settling) and dust disk scale heights over a wide range of stellar masses, at the same emitting regions of the disks.

Supplementing near- and mid-IR observations with far-IR measurements provide additional constraints on the structure of young circumstellar disks, and the Herschel observatory (Pilbratt et al. 2010) has already proven successful at studying these disks (e.g., Cieza et al. 2011; Joergens et al. 2012). Harvey et al. (2012b,a) present the results of their Herschel/PACS Phase 1 Guaranteed Time program with which they observed a significant number of M-type stars and BDs. Thanks to a conservative Bayesian approach to modeling the SEDs of about 50 sources, they find that the median disk mass for their BDs sample may be lower (M_disk ~ 3 × 10^-5 M_⊙) than for disks around more massive T Tauri stars (M_disk ~ 10^-3 − 10^-1 M_⊙). However, they find that disk properties, such as disk flaring and scale height, to be comparable with disks around T Tauri stars. In this paper we seek to complete these efforts by investigating disk properties in the low stellar mass regime, at the midway between T Tauri stars and brown dwarfs. We focus on stars in the M0-M5 spectral type regime in a single star-forming region.

The Chamaeleon-I (Cha-I hereafter) star-forming region is a perfect laboratory to study gas-rich disks around low-mass objects. Given its proximity (160−170 pc, see Luhman et al. 2008 for a detailed review) even M-type stars are bright enough to be observed in the far-IR with Herschel (e.g., Winston et al. 2012). The disk fraction in Cha-I is of about 60%, suggesting a relatively young age (~2 Myr, Hernández et al. 2007 and references therein, ~3 Myr, from the isochronal age distribution, Luhman et al. 2008). In this study we aim for a better characterization of the properties of disks around 17 M-type stars detected with the Herschel/PACS instruments (most likely the brightest M-type stars in Cha-I). Thanks to a careful SED modeling approach, supported with a Bayesian analysis on key parameters, we further discuss the disks properties around low-mass objects. In Sect. 2, we first present the Herschel observations and data processing. The modeling approach is described in Sect. 3 and the results are detailed and discussed in Sects. 4 and 5. We summarize our results and conclude in Sect. 6.

2. Observations and data processing

In the following we first present the Herschel/PACS observations of the Cha-I star-forming region, the data processing, and the sources that were successfully detected in the far-IR.

2.1. Herschel/PACS observations

In this study, we make use of Herschel observations of the Cha-I region from the guaranteed time key programme led by P. André at the three PACS wavelengths (70, 100 and 160 μm, in blue, green and red filters, respectively, Poglitsch et al. 2010). Maps in the green and red filters were produced from OBSIDs 1342224782 and 1342224783 (PACS photometry only), both taken on the 28/07/2011 with durations of 17 461 and 17 721 s, respectively. The map in the blue filter was reduced from OBSIDs 1342213178 and 1342213179, both taken in PACS and SPIRE parallel mode, on the 22/01/2011, with durations of 12 955 and 11 686 s, respectively. Because these observations were performed in parallel mode, their sensitivity is expected to be lower than the green map, and were therefore not used when producing the map at 160 μm. Data processing was performed within the “Herschel Interactive Processing Environment” (HIPE, version 10.0.667) using build-in scripts dedicated to the data reduction of such large maps. This script makes use of highpass filtering to remove the 1/f noise and is best suited for point sources.

The target sample was drawn from previous studies of the Cha-I region by Luhman & Muench (2008) and Luhman et al. (2008). Szűcs et al. (2010) have already compiled a list of sources in Cha-I with near-IR excess indicative of a siginificant amount of circumstellar dust, this list containing 62 M-type stars (the prime focus of this study). Here we only study the sources, that were detected at least at one of the PACS wavelengths (17 sources out of 62), but we will also discuss the detection statistics later on. The map at 100 μm appears to be the most sensitive one compared to other wavelengths. Flux extraction was performed via aperture photometry, with aperture radii of 10′′ for the blue and green maps, and 20′′ for the red map. Fluxes were corrected for the aperture sizes, depending on the radii used (fluxes divided by 0.774, 0.727, and 0.800 for the blue, green, and red filters, respectively). The pixel scales for the blue and green maps are 2′′, and 3′′ for the red map. To estimate the uncertainties, for each sources, we randomly placed 5000 apertures within a circle of a radius of 150 pixels around the source position (radii of 300′′ and 450′′ for both pixel scales). We then build a histogram of the measured fluxes in this area and fitted a Gaussian to the distribution of fluxes. The σ width of the Gaussian was used as the uncertainties on the measured fluxes. The center of the Gaussian fit is a good estimate of the background level, and was consequently subtracted from the source flux. These values were always lower than the derived uncertainties. Table 1 summarizes the Herschel/PACS fluxes and upper limits for M-type stars detected with a signal-to-noise ratio (S/N) higher than 3 at least at one wavelength. Uncertainties on the detections are 1σ, and upper limits are 3σ. A total of 17 sources were successfully detected, with spectral types ranging from M0.5 to M5.75. Young stellar objects with later spectral types, including BDs, were not successfully detected.

2.2. Near- and mid-IR complementary observations

In addition to the 2MASS data (Skrutskie et al. 2006), to construct the SEDs of the sources detected in the far-IR, we gathered photometric observations in the near- and mid-IR from literature: Spitzer/IRAC, MIPS data (Werner et al. 2004; Fazio et al. 2004; Rieke et al. 2004) both from Luhman & Muench (2008) and Luhman et al. (2008) and the WISE catalog (Wright et al. 2010). Table 2 summarizes these fluxes for each of the source. Additionally, several stars in our final sample (15 out of 17) were observed with the IRS spectrograph (Houck et al. 2004) onboard Spitzer. Spectroscopic data were processed using the FEPS pipeline (see, e.g., Bouwman et al. 2008).

The SEDs of the 17 M-type stars detected with PACS are displayed in Figs. 1 and 2. 2MASS and optical data are shown with filled, black circles, WISE measurements are represented with red hexagons and IRAC and MIPS observations are shown with green diamonds. PACS detections are represented by open, black circles, and upper limits by downwards black triangles. Where available, the Spitzer/IRS spectra are represented by a red line. The dashed line represents the stellar photosphere (see Sect. 3).

3. SED modeling: a Bayesian approach

3.1. The RADMC model

We used the radiative transfer code RADMC (Dullemond & Dominik 2004) to model the SEDs of the 17 sources detected at PACS wavelengths. A given model is first described by the stellar photosphere, and we opted for the NextGen atmospheric models (Hauschildt et al. 1999). Table 3 shows the adopted values for the effective temperature T_eff, stellar luminosity (L_bol) and foreground extinction (A_V), taken from the literature (Luhman 2007; Luhman & Muench 2008). We used a surface gravity of log(g) = 3.5, consistent with dwarf stars. For all the stars, we assumed a distance of d_⋆ = 165 pc. The disk itself is then described by the inner and outer radii (r_in and r_out, respectively), the disk mass (M_disk, assuming a gas-to-dust ratio of 100), and a surface density profile (Σ = Σ₀(r/r₀)^α, α < 0). The vertical distribution is assumed to follow a Gaussian profile (∝exp[ − z²/2H(r)²]), where H(r) is defined as H(r) = H₀(r/r₀)^γ where γ is the flaring exponent of the disk. The dust content is fully parametrized by the minimum and maximum grain sizes (s_min and s_max), a grain size distribution (dn(s) ∝ s^pds, p < 0), and the optical constants of astronomical silicates (Draine 2003) combined with the Mie scattering theory to compute the mass absorption coefficients.

Fig. 1 SED of J11071206-7632232 (#5). Grey dashed line represents the photosphere. Filled black circles are 2MASS measurements, green diamonds are Spitzer/IRAC and MIPS observations, red hexagons are the WISE measurements. Spitzer/IRS spectrum is shown in red. PACS observations are in open circles and upper limits are represented as filled black triangles. The most probable model is shown as a solid black line, and the gray lines correspond to models within the 68% confidence interval (see text for details).

Fig. 2 Same as Fig. 1 for the 16 other M-type stars detected at PACS wavelengths.

3.2. The grid of models

Our modeling approach is motivated by the work presented in Harvey et al. (2012a). The goal is not necessarily to achieve the best possible fit to the data, given the degeneracy of certain parameters in SED modeling without spatially resolved observations. Instead, we aim at finding the most probable values for several key parameters, such as M_disk or the scale height of disks around M-type stars. We therefore decided to choose a few relevant parameters, and keep the other constant. Based on the results presented in Harvey et al. (2012a), r_out is chosen to remain constant (100 AU) for all the objects. Determining the value of r_out is a challenging problem, even when millimeter observations are available. Only a handful of studies could constrain the outer radius of disks around low-mass stars and typical values of 15–40 AU were found (e.g., Luhman et al. 2007; Ricci et al. 2013). Several studies investigated whether or not r_out could be inferred from SED modeling of far-IR observations (e.g., Harvey et al. 2012b; Spezzi et al. 2013) and concluded the Herschel observations mostly remained insensitive to this parameter. This could indicate that the far-IR emission arises from regions within a few AUs of the disks and that the parameter r_out mostly has a significant impact on the millimeter emission. In an effort to limit the number of free parameters we opt not to include r_out in the grid of models, an approach in line with the modeling strategies described for instance in Andrews & Williams (2005, 2007); Morrow et al. (2008) and Mathews et al. (2013). On the other end, r_in is a crucial parameter to get the near-IR shape of the SED properly. Even though our study is not focused on the inner disk structure, r_in is a free parameter in the grid. The surface density distribution has an exponent of α =  −1, a value shallower than the Minimum Mass Solar Nebula, but consistent with theoretical estimations (e.g., Bell et al. 1997) and observational results (α ~  −0.9, Andrews et al. 2009). The grain size distribution follows a power-law with a slope in p =  −3.5 (Mathis et al. 1977). The grid of models is therefore computed over the following parameters: r_in, M_dust, H₀, and γ. For H₀, the reference radius is set at r₀ = 100 AU. Each model provides 10 SEDs for different inclinations i between 5° and 86° (0° being a face-on disk).

Even though we do not aim at fitting perfectly the IRS spectra, preliminary tests have shown that the observed emission features (especially at 10 μm) provide valuable information. We computed a first grid of models for the object J11071206-7632232, with s_min = 0.03 μm and s_max = 1 mm for the dust grains. The most probable model (see next paragraph) reproduces the overall shape of the SED but the 10 μm emission feature was severely under-predicted by all the models. We interpreted this as a sign for dust settling, and consequently included a “second” population which contains only small dust grains in the upper layers of the disk, with s_min = 0.01 μm and s_max = 0.1 μm. We parametrized this second population of dust grains via another scale height reference value H_0, small, that can be 0 (no stratification), 1, 2, or 5 AU above the H₀ for the larger grains, therefore mimicking a sedimentation of the larger grains. The small grain component was set to account for 10% of the total dust mass. The values for s_min and s_max for both dust populations were not changed during the fitting process.

Table 4 summarizes the free parameters explored in our grid of models (r_sub is the radius at which the silicate dust grains reach the sublimation temperature of 1500 K). For each source, we therefore computed a total of 25 200 SEDs. To obtain the most probable values for all the parameters aforementioned, we computed probability distributions for all the free parameters. For each model, the goodness of the fit is assessed via a reduced ${\mathit{\chi }}_{\mathrm{r}}^{\mathrm{2}}$. These ${\mathit{\chi }}_{\mathrm{r}}^{\mathrm{2}}$ values are then transformed into probabilities (∝exp[$\mathrm{-}{\mathit{\chi }}_{\mathrm{r}}^{\mathrm{2}}$ /2]), which are then projected onto each dimensions of the parameter space. This approach enables us to get a census on which parameters are well constrained (or not), as well as provide the most probable values. This methodology has already proven to be successfull for such studies (Harvey et al. 2012a; Spezzi et al. 2013).

To be able to compare consistently models for different sources, the goodness of fit has to be estimated in a similar fashion for all 17 sources. Significant photometric variability is observed for several sources at near- and mid-IR wavelengths (e.g., Spitzer/IRAC & MIPS versus WISE fluxes, see J11065906-7718535; Table 2), which is a main issue for SED modeling given that the photometric datapoints are non-simultaneous (e.g., Joergens et al. 2012). Several years have passed between the Spitzer, WISE, and 2MASS observations, and additionally Spitzer/IRAC, MIPS, and IRS data are usually not taken simultaneously. Photometric variability in the Class II phase can be caused for example by variable accretion (e.g., Muzerolle et al. 2009; Joergens et al. 2012) or occultation effects due to the disk (Alencar et al. 2010; Looper et al. 2010) and it affects mostly the near- and mid-IR fluxes and much less the far-IR.

Out of the 17 sources, only two sources have not been observed with the IRS instrument (J10533978-7712338 and J11111083-7641574), while 14 sources have IRS spectra between 5 and 35 μm and one source (J11062554-7633418) was only observed with the Short Low module (5.2−14.5 μm). We therefore opted to use the IRS spectra to construct the ${\mathit{\chi }}_{\mathrm{r}}^{\mathrm{2}}$, to minimize as much as possible the effect of time variability in the near- and mid-IR (i.e., IRAC, MIPS, and WISE measurements are not used when Spitzer/IRS spectra are available). For the 14 sources with a complete IRS spectra, we used five points at 5.5, 7.5, 9.8, 15, and 30 μm (averaged over ±0.1 μm). For J11062554-7633418, we used data from the IRS Short-Low module and used 4 points at 5.5, 7.5, 9.8, and 14 μm, plus the WISE W4 point at 22 μm. For the two sources without IRS data, we used the four WISE points and the IRAC point at 8 μm. For all the sources, the near-IR excess is constrained by the 2MASS K-band photometric point. Finally, the uncertainties were set to 10% of the observed fluxes for all wavelengths but the PACS wavelengths for which we used the uncertainties derived from the observations. The PACS upper limits were not included when computing the ${\mathit{\chi }}_{\mathrm{r}}^{\mathrm{2}}$ values.

4. Results

We detected 17 sources in the far-IR belonging to the young star-forming region Cha-I and successfully model their SEDs. The most probable fits to the SEDs, as defined by models with disk parameters at the peak positions of their respective probability distributions, are displayed as solid black lines in Figs. 1 and 2. The most probable disk parameters are summarized in Table 5. For a given source, the second line of the Table shows the range of validity for the considered parameter, obtained from the probability distribution. Given the distribution P(θ) for the parameter θ between θ_min and θ_max, the range of validity [θ₁,θ₂] is found for ${{}^{\mathrm{\int }}}_{\mathit{\theta }\min }^{\mathit{\theta }\mathrm{1}}\mathit{P}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\hspace{0.17em}\mathrm{d}\mathit{\theta }\mathrm{=}{{}^{\mathrm{\int }}}_{\mathit{\theta }\mathrm{2}}^{\mathit{\theta }\max }\mathit{P}\mathrm{\left(}\mathit{\theta }\mathrm{\right)}\hspace{0.17em}\mathrm{d}\mathit{\theta }\mathrm{=}\mathrm{(}\mathrm{1}\mathrm{-}\mathit{\beta }\mathrm{)}\mathit{/}\mathrm{2}$, where β = 0.68 (see e.g., Pinte et al. 2008). The interval [θ₁,θ₂] corresponds to the 68% confidence interval for each parameters. In Figs. 1 and 2, models that are within this interval are plotted in light gray.

4.1. Edge-on disks and problematic sources

We could hardly find satisfying models for two sources: J10533978-7712338 and J11111083-7641574 (#1 and 16, respectively). Source #16 is most likely an edge-on disk, if we use the luminosity of 0.003 L_⊙, as reported in the literature, the near-IR fluxes can be reproduced regardless of the disk inclination (as long as it is not too high), but the far-IR fluxes would be very much underestimated. The differences between the modeled SED and the far-IR observations could be reduced by setting the inclination to a higher value, and increasing the stellar luminosity to 0.8 L_⊙ to match the near-IR data points. Even though our final model still slightly under-predicts the PACS measurement at 100 μm, our finding compares well with the modeling presented in Robberto et al. (2012), where the authors found i ~ 87.1° for this source.

The main issue in the modeling of source #1 comes from the high fluxes at WISE W4 and MIPS 24 μm wavelengths (59.2 and 55.4 mJy, compared to our modeled fluxes of 26.4 and 28.4 mJy at 22 and 24 μm, respectively). First of all, detailled inspection of the individual WISE W4 images show the source is bright, clearly detected, and isolated from other nearby sources. The overall small reduced χ² in the WISE catalog suggests it is a point source without any significant contamination. This is confirmed by the MIPS 24 point, which is close to the W4 point. However, no disk models could reproduce both the mid-IR datapoint and the very steep decrease toward the PACS measurement (detected at more than 5σ). The object is not located close to any bright region in the Herschel map and does not appear as extended either. Robberto et al. (2012) and Luhman & Muench (2008) suggested the source is seen almost edge-on and is mainly seen in scattered light. In the following statistical analysis both problematic sources will be excluded from the sample.

4.2. Detection statistics and probabilities

The original stellar sample from Szűcs et al. (2010) contained 62 M-type stars, and 17 sources out of these 62 were successfully detected at least at one PACS wavelengths. Table 6 summarizes the number of objects detected and the total number of sources in the stellar sample for different spectral types from M0 to M9. First of all, while it is expected to not detect very low-mass objects, it appears that several M1 stars of the original list are not detected by the PACS instrument (25% detection rate). Close inspection of the 100 μm PACS map reveals that several of these M1 stars are close to bright, extended structures, rendering any detection impossible.

Fig. 3 Detection probabilities at λ = 100 μm, based on the grid of models for the M5 star J11120984-7634366 (#17), as a function of the disk mass for different flaring indices (see text for details).

Based on the numerous SEDs modeled within the grids, we can better investigate the probabilities of detecting a disk for a certain set of parameters. Because the 100 μm map is the most sensitive of all three PACS bands, we considered a 3 σ detection threshold of 60 mJy to be representative of the observations (see Table 1). We then estimated the 100 μm fluxes for all the models of a single source (e.g., source #17 for a M5 star) and computed the fraction of models that would be detected as a function of certain disk parameters. Figure 3 shows the detection probabilities as a function of the disk mass, where different values for the flaring indices γ are color-coded, for the M5 star J11120984-7634366. The fact that for low disk masses the detection probabilities are higher for γ = 1.0 compared to γ = 1.1 is not a numerical issue in the radiative transfer calculations. We tested this against the number of photons in the model, and the resolution of the radial grid. The fact that the detection probabilities become lower for γ = 1.0 at higher masses tends to indicate that different regions in the disks are probed at 100 μm for low disk masses (optical depth effect). Besides this, one can immedialety see that the PACS observations at 100 μm are not sensitive enough to detect disk masses of about 10^-5 M_⊙ around this star. These simple considerations enable us to better understand the detection rates presented in Table 6. All the M0 stars of our stellar sample were detected because no matter how massive, or how flared, the disk is, the chances of detecting it are overall large, while the probabilities of detecting relatively low-mass, reasonably flared disks (e.g., M_disk = 10^-4 M_⊙, γ ~ 1.1) around M5 stars rapidly drop below 50%. According to our modeling results, all the three M5 stars detected in our sample have high flaring indices (γ = 1.3 for J11062554-7633418 and J11120984-7634366, and 1.2 for J11105597-7645325), which should not be interpreted as evidence for higher γ values for disks around M5 stars, in addition to the overall challenge when dealing with low-number statistics.

The stellar sample considered in this study originally aimed at selecting disk-bearing sources (Szűcs et al. 2010), based on their near-IR excesses. However, 15 M5-type stars that are known to harbor a disk were not detected in the Herschel observations. Manoj et al. (2011, and references therein) reported that five M5 stars in Cha I are in multiple systems, some of them having small separations. For instance 0.17′′ separates the primary and secondary of the Cha Ha2 system (Ahmic et al. 2007; Lafrenière et al. 2008) corresponding to about 28 AU at 165 pc. Such small separation can have a strong impact on the disk (see e.g., Bitner et al. 2010), and may explain why the disk is not detected at PACS wavelengths. However, out the of five M5 stars with known multiplicity, two of them are detected with the Herschel data even at small separation: objects #6 (28.9′′ separation, Kraus & Hillenbrand 2007) and #15 (0.13′′ separation Ahmic et al. 2007; Lafrenière et al. 2008). It does not seem (based on five sources) that binarity can solely explain the low detection rate of M5 stars in the Herschel map of Cha-I. Even though the detection probabilities shown in Fig. 3 are derived for a single source, with its own luminosity, this would suggest that the disks around these 15 undetected M5 stars are most likely low-mass (M_disk ≲ 5 × 10^-5 M_⊙ for a flaring index of about γ ~ 1.1), or have low flaring indices (γ ≤ 1.1) for slightly higher disk masses, the two parameters being most likely degenerate.

4.3. Disk masses and disk evolution

As discussed in Sect. 4.2, we most likely detected the most massive, or more flared disks around M3-M5 stars. This detection bias renders the question of disk masses around M-type stars challenging at the sensitivity level of the Herschel maps. Nonetheless, our study complements the current effort towards a better characterization of disks around low-mass stars (e.g., Harvey et al. 2012a; Joergens et al. 2012). Our sample contains slightly more massive objects and fills in the gap for several M0–M3 type stars. We find that disks around such objects have masses in the range 10^-4–10^-3 M_⊙. Based on the detection probabilities as a function of PACS sensitivity at 100 μm, we already concluded that disks around M5 stars (or later spectral types) should either (i) have masses of the order of 10^-5 M_⊙, or below, a result consistent with the study from Harvey et al. (2012a), or (ii) appear flatter than early M-type stars (low γ values). Additionally, we searched for a correlation between the estimated disk masses and the 100 μm fluxes. To estimate the significance (or lack) of the correlation between the two quantities, we used the Kendall τ correlation coefficient (between 1 and − 1, for a perfect correlation and anti-correlation, respectively), and its associated probability P (the lower, the more robust). We obtained a Kendall τ correlation coefficient of 0.23, with a probability P = 0.22. The lack of a correlation implies that the measured far-IR fluxes are only a reliable measure by themselves for disk masses in the regime M_disk ≲ 10^-5 M_⊙ where a larger fraction of the total disk mass lies in regions of smaller optical depths in the far-IR (Harvey et al. 2012a). This result implies that, as opposed to mm observations, fluxes in the far-IR cannot be directly used to measure disk masses around M0-M5 type stars. A careful SED modeling effort is the next best approach to estimate disk masses in the higher mass regime.

Fig. 4 Most probable disk masses Mdisk as a function of near- to far-IR color ([3.6]  −  [100]).

In an attempt to link the most probable disk masses to the evolution of the circumstellar disk, we show in Fig. 4 the disk masses inferred from the modeling as a function of the near- to far-IR color ([3.6]  −  [100]), quantity defined as $\begin{array}{ccc}\mathrm{\left[}\mathit{a}\mathrm{\right]}\mathrm{-}\mathrm{\left[}\mathit{b}\mathrm{\right]}\mathrm{=}\mathrm{-}\frac{\log \mathrm{\left(}{\mathit{\lambda }}_{\mathit{b}}{\mathit{F}}_{\mathit{b}}\mathrm{\right)}\mathrm{-}\log \mathrm{\left(}{\mathit{\lambda }}_{\mathit{a}}{\mathit{F}}_{\mathit{a}}\mathrm{\right)}}{\log \mathrm{\left(}{\mathit{\lambda }}_{\mathit{b}}\mathrm{\right)}\mathrm{-}\log \mathrm{\left(}{\mathit{\lambda }}_{\mathit{a}}\mathrm{\right)}}\mathit{,}& & \end{array}$(1)where F_a,b is the flux density in units of erg s^-1 cm^-2 μm^-1 at wavelength λ_a,b in units of μm. As the near- to far-IR colors become bluer (increasing [3.6] − [100]), the modeled disk masses become lower, a result in line with (sub-)mm surveys of young stellar objects (e.g., Andrews & Williams 2007). Including all our sources, we find a Kendall τ correlation coefficient of − 0.58 with an associated probability P = 1.1 × 10^-3. Consequently, even if the most probable disk masses do not correlate well with the measured Herschel fluxes, they do anti-correlate relatively well with the near- to far-IR shape of the SEDs, a quantity that can be related to the far-IR opacity of the disks (disk evolution or substantial growth of dust particles, Beckwith et al. 1990; Andrews & Williams 2007).

4.4. Disk flaring

Fig. 5 Most probable flaring indices γ as a function of mid- to far-IR color ([11.56]  −  [100]).

Our modeling results point toward a relatively broad distribution of flaring indices (1 ≤ γ ≤ 1.3). We can therefore search for a first order approximation between the SEDs and the flaring indices. Figure 5 shows the mid- to far-IR colors (between the WISE band W3 and the 100 μm flux) as a function of the most probable flaring indices γ. In Fig. 5, sources with γ ≥ 1.2 are shown with open black circles, and γ ≤ 1.15 with red squares. The two sources that are most likely edge-on (#1 and #16) are shown as green diamonds. Excluding the two edge-on sources from the sample (#1 and #16), we find a correlation coefficient τ = 0.73 with the associated probability P = 1.6 × 10^-4, suggesting the two quantities correlate well. We find that a threshold value [11.56] − [100]  ~ 0.5 appears to separate relatively well flat and flared disks.

4.5. A tentative characterization of dust vertical distribution

The Spitzer/IRS observations are probing the optically thin upper layers of the circumstellar disks. Even though we did not aim at fitting the emission features perfectly, our parametrization of the dust stratification (via the parameter H_0, small) has proven to be a good approach to characterize dust settling and reproduce the overall shape of the 10 μm emission feature (especially their peak fluxes). For all the sources but two (sources #6 and 12), a population of small dust grains in the uppermost layers is required to provide a satisfaying match to the IRS data, even though the confidence intervals for H_0, small are fairly large (except for source #9 which has the strongest 10 μm feature in our sample). The IRS spectra of sources #6 and 12 do show a 10 μm emission feature, but the features appear slightly weaker and flatter than for some of the other sources (e.g., sources #9, 11 or 13).

To investigate deeper this possible relation, we computed the ratio between optically thin and thick emission at 10 μm, from the IRS spectra. As a first order approximation, the emission feature at 10 μm indeed arises from optically thin regions of the disk, while the continuum emission originates from optically thick regions. We therefore defined a linear, local continuum (F_ν,cont) underlying the observed 10 μm emission feature (F_ν,obs) between λ₁ = 7.5 and λ₂ = 13.5 μm. The observed ratio between optically thin and thick emission can then be approximated by the quantity ${{}^{\mathrm{\int }}\mathit{\lambda }_{\mathrm{1}}}_{{\mathit{\lambda }}_{\mathrm{2}}}^{}\mathrm{(}{\mathit{F}}_{\mathit{\nu ,}\mathrm{obs}}\mathrm{-}{\mathit{F}}_{\mathit{\nu ,}\mathrm{cont}}\mathrm{)}\hspace{0.17em}\mathrm{d}\mathit{\lambda }\mathit{/}{{}^{\mathrm{\int }}\mathit{\lambda }_{\mathrm{1}}}_{{\mathit{\lambda }}_{\mathrm{2}}}^{}{\mathit{F}}_{\mathit{\nu ,}\mathrm{cont}}\hspace{0.17em}\mathrm{d}\mathit{\lambda }$. We find that this quantity may correlate tentatively with the H_0,small parameter for the most probable fits to the SED (Kendall τ = 0.39, with the associated probability P = 0.05), but there is no sign of correlations with other parameters such as M_disk, H₀ or γ. By adding more small grains in the optically thin regions of the disk, the parametrization of the dust stratification seems to be a good solution to reproduce strong 10 μm emission features as the one of object #9. For the other sources of the sample, the situation is more ambiguous (hence the low correlation Kendall τ coefficient). Computing the value S_peak from the IRS spectra, as in Kessler-Silacci et al. (2006), to quantify the amount of dust processing, we find a median value of 1.7 with a standard deviation of 0.7, indicating “boxy” emission features attributed to μm-sized grains or shallow grain size distributions (Olofsson et al. 2010). These values match the vast majority of previous Spitzer/IRS observations of young circumstellar disks (e.g., Bouwman et al. 2008; Furlan et al. 2009; Olofsson et al. 2009; Manoj et al. 2011; Oliveira et al. 2011). However, while we aimed for reproducing the peak flux of the 10 μm feature, we did not aim for matching their shape (grain size distribution and crystallinity), a question beyond the scope of this paper.

Nonetheless, characterizing dust settling based on SED modeling is a challenging problem. Dullemond & Dominik (2008) have shown that dust sedimentation can have an effect on the 10 μm emission feature when the grain size distribution is bimodal (two different grain sizes), but this effect is severely damped when considering a continuous grain size distribution. Our approach is closer to the bimodal distribution as we do not have a full description of dust settling. Consequently, the question of dust stratification in these disks will have to be further investigated, for instance by observing a flat sub-mm slope confirming the presence of mm-sized dust grains in the disk midplane.

4.6. Sub-mm LABOCA observations

Fig. 6 LABOCA observations of source # 6 (marked with a white “+” symbol), which suggest the measured flux is contaminated by the nearby source. The beam size is displayed in the lower left corner.

The Cha-I star-forming region has been observed in the sub-mm with the bolometer array LABOCA at the APEX telescope (Belloche et al. 2011). Out of the 17 sources detected with Herschel/PACS, three sources were also detected at 870 μm (sources # 6, 9, and 16). The fluxes reported by Belloche et al. (2011) are 94, 132, and 72 mJy for these three sources, respectively. Inspecting the LABOCA map, we find the flux for source # 6 (J11074366-7739411, see Fig. 6) may suffer from contamination by a nearby source. While the close-by source is also detected in the 100 μm PACS map, it is much more compact and does not contaminate the far-IR measurement. Therefore, we estimated the peak flux in the 870 μm map by fitting a Gaussian profile with a full width at half maximum of 21.2′′ (equal to the beam size) on top of the background emission, using the IDL package Starfinder (Diolaiti et al. 2000). We obtain a peak flux of about ~58 mJy (instead of 94 mJy). Our most probable models severely under-predict the observed sub-mm fluxes, with predicted fluxes of 10, 14 and 10 mJy for sources # 6, 9, and 16, respectively.

The first possible explanation for these discrepancies is that we under-estimate the disk masses with our models. Belloche et al. (2011) reported masses of 4.4 × 10^-3, 1 × 10^-2, and 5.5 × 10^-3 M_⊙ for sources # 6, 9 and 16, respectively (assuming F_ν = 58 mJy at λ = 870 μm for source # 6), a factor 5 to 10 above our results. In Sect. 4.3, we already discussed that the measured PACS fluxes do not correlate with the inferred disk masses, indicating that the emission is not fully optically thin at 100 μm. Consequently, the Herschel observations are not direct tracers of the amount of dust in the disks we modeled (hence the effort of modeling the SED from optical to far-IR wavelengths, see also Sect. 5 for further discussion). However, the purpose of this study is to compare relative disk parameters for stars of different spectral types, modeled in a consistent way. Furthermore, the fluxes differences in the sub-mm can also be the consequences of several assumptions that have to be made. First, as demonstrated by Draine (2006), the opacity in the sub-mm wavelength regime highly depends on the maximum grain sizes, the grain size distribution and the dust composition. The dust content is described by the optical constants of astronomical silicates, while carbonaceous grains, icy grains could well be present in the outermost regions of the disk. This would in turn modify the sub-mm slope of our diks models. We choose s_max = 1 mm, which is large enough at PACS wavelengths but may not be the best choice to model observations at 870 μm (following the criterion of s_max ≥ 3λ of Draine 2006). Increasing the maximum grain size to a few mm may slightly increase the sub-mm fluxes of our models. Additionally, far-IR measurements may not be sensitive to cold dust in the outermost regions, suggesting that regions of the disk may remain unseen by Herschel. Finally, optical constants are usually measured at room temperature, while the dust grains in the outermost regions of circumstellar disks can be as cold as a few tens of Kelvin. Coupeaud et al. (2011) have demonstrated that the temperatures of the grains during laboratory experiments have drastic effects on their opacities, especially in the sub-mm. Investigating these effects for three sources in our sample with sub-mm data is out of the scope of this study.

5. Discussion

Fig. 7 Disk mass to stellar mass ratio (left panel) and disk mass (right panel) as a function of stellar mass. On each panel, the shaded area represents log (Mdisk/M⋆) =  −2 ± 1. Vertical gray lines correspond to uncertainties, when available. Data represented with green “×” symbols from Osterloh & Beckwith (1995), open black diamonds from (Natta et al. 2000, and references therein), dark “+” symbols are from Scholz et al. (2006), cyan hexagon from Joergens et al. (2012), open blue stars from Ricci et al. (2012, 2013), orange squares from Harvey et al. (2012a, Cha-I sample only), blue triangles from Spezzi et al. (2013, lower limits), and open red circles from this work. The color-coding in the right panel is the same as in the left panel. Filled symbols correspond to disk masses inferred from modeling of far-IR observations (this work; Joergens et al. 2012; Harvey et al. 2012a; Spezzi et al. 2013).

Several studies have speculated that some disk properties are different for disks around Sun-like stars and very low-mass objects (e.g., Pascucci et al. 2009). The two categories being usually defined for spectral types between G5–K1 to M4.5 for the first group and M5 to M9 for the latter group. For instance, based on Spitzer/IRACs color distributions of a statistically large sample, Szűcs et al. (2010) concluded that disks may appear to be flatter around cooler stars, a consequence of dust settling towards the disk midplane. While the low-number statistics and the biases inherent to the sensitivity of the Herschel/PACS observations prevent us from drawing firm conclusions our findings can still be put into a broader picture by comparing them to other studies.

Mainly, the question of the available mass in disks around low-mass stars has motivated several observational campaigns in the past years, to better assess the efficiency of planet formation in such systems. The underlying motivation being to assess if the disk mass is a function of the central stellar mass. Therefore, to compare the disk masses inferred in our analysis as a function of M_⋆, we gathered results from the literature, for a broader range of stellar masses. First, for our 17 M-type stars, we estimated the masses of the host stars by interpolating isochrones from Baraffe et al. (2003, assuming an age of 2 Myr) over their effective temperatures T_eff (see Table 3). In Fig. 7, the ratio M_disk/M_⋆ as a function of the stellar mass is shown in the left panel and the relation M_disk as a function of M_⋆ is displayed in the right panel. We show results from our analysis and from the literature (Osterloh & Beckwith 1995; Natta et al. 2000; Scholz et al. 2006; Joergens et al. 2012; Ricci et al. 2012, 2013; Harvey et al. 2012a; Spezzi et al. 2013). The stellar masses for the Cha-I sample from Harvey et al. (2012a) were determined from the stellar effective temperatures, in a similar way as for our sample. One should keep in mind that such a representation is affected by detection biases, especially for very low-mass stars. But overall, no clear trends are visible, as there is a significant vertical dispersion at all stellar masses, a conclusion in line with the results of Scholz et al. (2006) and Williams & Cieza (2011). In both panels of Fig. 7 the shaded area correspond to the trend log (M_disk/M_⋆) =  −2 ± 1 observed in (sub-)mm surveys (e.g., Williams & Cieza 2011). It appears that several disk masses inferred from SED modeling of far-IR observations (filled symbols) fall below the trend determined from (sub-)mm observations. This may indicate that optical depth effects are still important in the far-IR and that the modeled disk masses may be lower limits estimates, keeping in mind that far-IR observations may not be sensitive to cold dust in the outermost regions. One can note from Fig. 7 that the dispersion in disk masses is lower for low-mass stars than for sources with solar masses, which may be related to lower far-IR optical depths for disks around low-mass stars. Overall, the absolute values of M_disk seem to be lower for disks around (very) low-mass stars, indicating that massive planets (a few M_Jup) can hardly form in these disks (e.g., Mordasini et al. 2012). However, the relative M_disk/M_⋆ ratio appears to remain constant over a wide range of stellar masses, a valuable input information for planet formation theories. Future ALMA observations will provide more accurate mass determination by probing optically thin emission, with a better sensitivity, and will enable us to investigate deeper the stellar mass dependence.

6. Conclusion

We studied and modeled the SEDs of 17 disk-bearing, low-mass M-type stars (spectral types between M0 and M5.75) detected at Herschel/PACS wavelengths in Cha-I (out of 62 M-type stars). The increased SED coverage at far-IR wavelength, combined with a careful modeling approach, enabled us to discuss further the disk masses, flaring indices and dust settling in regions of the disks responsible for emission up to far-IR wavelengths.

First, we discussed the detection rates of known disk-bearing sources (based on near-IR excess) and concluded the Herschel sensitivity for these observations not to be sufficient to detect disk masses below ~10^-5 M_⊙, or flaring indices below γ ≤ 1.1, for disks around M5 stars (or later spectral types). For sources detected in the PACS map, we find a distribution of disk masses between 10^-4 − 10^-3 M_⊙, which compares well with other studies of similar sources in the far-IR (e.g., Harvey et al. 2012a; Spezzi et al. 2013). We find that the most probable disk masses anti-correlate with the near- to far-IR color of the SED. However, the far-IR fluxes do not correlate well with the derived disk masses, suggesting the contribution of regions with high optical depth, which are often non negligeable even at these wavelengths. Consequently, for several sources in our sample, we find that the most probable masses lie below the trend log (M_disk/M_⋆) =  −2 ± 1 observed in (sub)-mm surveys of young disks. This suggests the disk masses found via modeling of far-IR observations can most likely be considered as lower limits estimates. Future ALMA observations will provide more reliable disk masses estimates (e.g., Ricci et al. 2012, 2013). Finally, we suggest that dust settling may already take place in these young disks (~2 Myr), via the analysis of the 10 μm emission feature, and especially its strength with respect to the underlying continuum. The parametrized stratification enables us to constrain the presence of small dust grains in the uppermost optically thin layers of the disks, while the observations still remain un-sensitive to the grain sizes in the disk midplane. Sub-mm observations are expected to provide further constraints on the dust grain sizes in regions closer to the midplane (e.g., Boehler et al. 2013).

Online material

Acknowledgments

The authors thank the anonymous referee for constructive comments that improved the paper and the discussion about the disk masses. L. Sz. acknowledges support from the Deutsche Forschungsgemeinschaft via SFB project 881 “The Milky Way System” (sub-project B3). The radiative transfer models were calculated on the bwGRiD (http://www.bw-grid.de), member of the German D-Grid initiative, funded by the Ministry for Education and Research (Bundesministerium für Bildung und Forschung) and the Ministry for Science, Research and Arts Baden-Württemberg (Ministerium für Wissenschaft, Forschung und Kunst Baden-Württemberg) and KOLOB cluster at the Heidelberg University. The authors thank Paul M. Harvey and Christophe Pinte for providing the disk masses of their Cha-I sample, and Amelia Bayo for fruitful discussions. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France.

References

All Tables

All Figures

Fig. 1 SED of J11071206-7632232 (#5). Grey dashed line represents the photosphere. Filled black circles are 2MASS measurements, green diamonds are Spitzer/IRAC and MIPS observations, red hexagons are the WISE measurements. Spitzer/IRS spectrum is shown in red. PACS observations are in open circles and upper limits are represented as filled black triangles. The most probable model is shown as a solid black line, and the gray lines correspond to models within the 68% confidence interval (see text for details).

In the text

	Fig. 2 Same as Fig. 1 for the 16 other M-type stars detected at PACS wavelengths.
In the text

	Fig. 3 Detection probabilities at λ = 100 μm, based on the grid of models for the M5 star J11120984-7634366 (#17), as a function of the disk mass for different flaring indices (see text for details).
In the text

	Fig. 4 Most probable disk masses Mdisk as a function of near- to far-IR color ([3.6]  −  [100]).
In the text

	Fig. 5 Most probable flaring indices γ as a function of mid- to far-IR color ([11.56]  −  [100]).
In the text

	Fig. 6 LABOCA observations of source # 6 (marked with a white “+” symbol), which suggest the measured flux is contaminated by the nearby source. The beam size is displayed in the lower left corner.
In the text

Fig. 7 Disk mass to stellar mass ratio (left panel) and disk mass (right panel) as a function of stellar mass. On each panel, the shaded area represents log (Mdisk/M⋆) =  −2 ± 1. Vertical gray lines correspond to uncertainties, when available. Data represented with green “×” symbols from Osterloh & Beckwith (1995), open black diamonds from (Natta et al. 2000, and references therein), dark “+” symbols are from Scholz et al. (2006), cyan hexagon from Joergens et al. (2012), open blue stars from Ricci et al. (2012, 2013), orange squares from Harvey et al. (2012a, Cha-I sample only), blue triangles from Spezzi et al. (2013, lower limits), and open red circles from this work. The color-coding in the right panel is the same as in the left panel. Filled symbols correspond to disk masses inferred from modeling of far-IR observations (this work; Joergens et al. 2012; Harvey et al. 2012a; Spezzi et al. 2013).

In the text