A Radial Velocity Survey of Embedded Sources in the Rho Ophiuchi Cluster

Sources were selected based on the spectral index, α, of their SEDs. The majority of our sample have Class I or FS SEDs, with some Class II objects. In addition to this, sources were selected based on their K magnitudes (λ = 2.2 μm) depending on the instrument and telescope used. For NIRSPEC and the Cryogenic Infrared Echelle Spectrograph (CRIRES), this resulted in sources with K ≤ 11; for CSHELL and iSHELL, K ≤ 10. As shown in Figure 1, the sources are concentrated within the L1688 cloud and members of subclusters corresponding to dense cores A (northwest), B (northeast), and E/F (south).

Zoom In Zoom Out Reset image size

Observations were made using CSHELL on NASA's 3 m Infrared Telescope Facility (IRTF) on the summit of Maunakea in Hawaii. Seven first-half nights in 2016 June 23, 29, and 30 and July 10, 11, 13, and 14 were used to observe 10 YSOs and 5 radial velocity standards. Observations from June 29 and 30 were made through some cirrus cloud cover. CSHELL was an infrared echelle spectrograph with a 256 spectral by 160 spatial pixel array (Greene et al. 1993). The 10 wide slit was used providing a spectral resolving power of 21,500. Seven YSOs were observed at a central wavelength of 2.298 μm, three at 2.2935 μm, and one at both. Each setting covered a wavelength range of approximately 50 Å. The 2.298 μm setting was chosen to detect nine prominent absorption lines from the CO v = 2–0 band head. The 2.2935 μm setting was chosen to detect the band head itself in sources where the longer-wavelength absorption lines might be too faint to detect.

In each case, objects were nodded 15'' along the slit in an ABBA pattern. Offset guiding was possible for most, but not all, sources. Telluric standards were taken every time the air mass had changed by approximately 0.2. Data reduction was performed using the Image Reduction and Analysis Facility (IRAF⁴ ). Master flats were made using dark-subtracted frames, and the data were then flat-fielded before being corrected for bad pixels. Due to the nodding, data images were subtracted in their AB pairs, and this subtraction ensured that dark current, as well as the sky background in the data images themselves, was removed. Spectra were then extracted, normalized, and divided by the telluric standards. The final step was to wavelength-calibrate the spectra using 10 telluric methane lines identified from the HITRAN database (Gordon et al. 2017).

Further observations were conducted using the new instrument iSHELL on NASA's 3 m IRTF the following year (Rayner et al. 2016). Three second-half nights in 2017 April 26–28 were used to observe 12 YSOs and 4 radial velocity standards, focusing on a different group of slightly fainter objects than what was possible using CSHELL. iSHELL is a 1.1–5.3 μm cross-dispersed high-resolution echelle spectrograph. For these observations a 075 slit was used, yielding a spectral resolving power of R = 35,000. The K2 setting was selected for each object and covered a wavelength range of 2.09–2.38 μm over 32 orders, ensuring that, once again, the CO v = 2–0 band head, as well as the 3–1 and 4–2 band heads, was observed, plus various atomic absorption lines.

The slit used was only 5'' long, and therefore no nodding was done. Flats, calibration lamps, and telluric standards were observed at regular intervals. Observing efficiency was aided by guiding with an infrared slit viewing camera.

Data reduction was performed using iSHELL Tool, a program developed by Michael Cushing for the specific purpose of reducing iSHELL data based on an earlier program for a cross-dispersed spectrograph (Cushing et al. 2004). The first step in this process was to create master darks and flats. The program would then generate a wavelength calibration solution using Thorium–Argon arc lamp images that were taken as close in time to the data as possible. Following identification of the apertures for spectral extraction, dark subtraction, flat-fielding, background subtraction, spectral extraction, and wavelength calibration would all be performed automatically. The next step was to apply these same steps to the telluric standards, divide out the telluric lines, and clean up any bad pixels. It should be noted that during the telluric division a small wavelength shift was applied to the telluric based on the cross-correlation of the YSO spectrum and the telluric in order to improve telluric division. This implies a small uncertainty in the wavelength calibration solution that varies for each source. The shift was always less than 0.5 pixels and was recorded for the spectrum of each source.

In addition to our own observations, high-resolution infrared spectra were generously provided from the published observations of Doppmann et al. (2005) obtained using NIRSPEC on the 10 m Keck II telescope on Maunakea, Hawaii (McLean et al. 1998). Observations were made on 2000 May 30, 2001 July 7, 8, and 10, and 2003 June 20–21. Fifteen Class I and FS YSOs within the L1688 region were observed, usually having K ≥ 10. The 058 slit was used, providing a spectral resolving power of R = 18,000. Several orders were obtained in the original data, but for this project, only order 33, containing the CO v = 2–0 band head, was analyzed. All data were reduced by Doppmann et al. using IRAF. Images were cleaned of bad pixels and cosmic rays, flat-fielded, and sky-subtracted. Extracted spectra were then wavelength-calibrated and co-added with spectra at the same slit position and similar air masses. Further details can be found in the original paper.

Spectra for seven YSOs were reduced from the European Southern Observatory archive taken with the CRIRES, which is an adaptive-optics-assisted spectrograph mounted on the 8 m VLT UT1 (Antu) located at Paranal Observatory, Chile (Kaufl et al. 2004). Data were from Viana Almeida et al. (2012) (Program ID 081.C-0395(A)), Viana Almeida (2012, unpublished data; Program ID 089.C-0539(A)), and Cottaar (2012, unpublished data; Program ID 089.C-0753(A)). CRIRES is capable of delivering a spectral resolving power of up to 100,000 in the 960–5200 nm wavelength range. Light is projected onto four detectors that are each 1024 × 512 pixels. Viana Almeida et al. chose to cover a wavelength range of 2.2542–2.3047 μm using a 03 slit to achieve a resolution of 60,000. This put the CO v = 2–0 absorption bands on detector 4, but CO bands were placed on detectors 1 and 3 for other sets of observations. Spectra were collected in an ABBA nodding pattern, allowing a very similar reduction process to that which was used for CSHELL as described above, including wavelength calibration using telluric methane lines.

The general procedure was to read in both the data and the synthetic as (x, y) pairs and then bin the synthetic to the same wavelength range and resolution as the data. The likelihood and priors were then constructed. When using this program, however, instead of working with regular probabilities, computations were done with the log of the probabilities in order to ensure that certain values remain positive. The log of the likelihood was then given by

$$\begin{eqnarray}&&\mathrm{ln}[\,p(y| \sigma ,m)]=-\displaystyle \frac{1}{2}\sum \,\left[\displaystyle \frac{{\left({y}_{i}-{m}_{i}\right)}^{2}}{{\sigma }^{2}}+\mathrm{ln}(2\pi {\sigma }^{2})\right]\,,\end{eqnarray} \tag{ 3 }$$

where y_i represents a data point, m_i represents a point on the model, and σ is the value of the data point divided by the approximate S/N. The next step was to construct the priors, which was done by setting up physically reasonable ranges for each parameter informed by physical constraints and by previous estimates in the literature when possible. With this information, the log of the probability was calculated, which was then used by emcee when sampling the distribution. The end result was the relative probability distributions for each of the parameters being fit, as well as their 1σ uncertainties.

The four parameters being fit were the veiling, rotational velocity, radial velocity, and effective temperature. The surface gravity was chosen in advance depending on whether the object in question was a YSO (log g = 3.5) or a radial velocity standard (log g = 4.5). For the program to run correctly, both the synthetic and object spectra needed to be normalized and have the same resolution and number of channels. Dereddening the spectra was not necessary given the narrow wavelength range of each order. For the wavelength range covered by most of the data, the vast majority of lines being fit were CO band heads and absorption lines from the v = 2–0 energy level transition. Physically, the veiling is a measure of the emission from circumstellar gas and dust and is measured through the veiling coefficient, r_k:

$$\begin{eqnarray}&&{F}_{v}=\displaystyle \frac{F+{r}_{k}}{1+{r}_{k}},\end{eqnarray} \tag{ 4 }$$

where F_v is the veiled spectrum and F is the unveiled spectrum. For example, for r_k = 1, the absorption-line depth is half of that compared to r_k = 0. Due to the sensitivity of r_k to the normalization, the veiling coefficient returned by the program was only accurate to first order but important in fitting the other parameters of interest. A simple rotational broadening kernel from the PyAstronomy package called rotBroad() was applied to the synthetic spectra. Unlike r_k, a higher value of v sin i will result in broader and shallower absorption lines but does not change the area of the absorption lines. The key parameter for this study being fit was the radial velocity. This was modeled as a simple fractional pixel shift, assuming that there was a linear variation in the normalized flux between each pixel. Because the shift was measured in pixels, the uncertainties in the radial velocities being measured were directly correlated with the resolution of the instrument. CRIRES and iSHELL data produced the most precise results, with NIRSPEC data being the least precise. An equation describing the modeling of these parameters is given by

$$\begin{eqnarray}&&{F}_{\mathrm{model}}[i]=\mathrm{rotBroad}\left(\displaystyle \frac{{F}_{{T}_{\mathrm{avg}},\mathrm{log}g}[i+{dx}]+{r}_{k}}{1+{r}_{k}},v\sin \,i\right),\end{eqnarray} \tag{ 5 }$$

where F_model is the broadened, veiled, and shifted synthetic spectrum; ${F}_{{T}_{\mathrm{avg}},\mathrm{log}g}$ is the unbroadened, unveiled, and unshifted synthetic spectrum corresponding to a certain log g and T_eff; dx is the pixel shift; r_k is the veiling factor, v sin i is the rotational velocity in km s⁻¹; and the rotBroad() function is a convolution determined by v sin i.

As alluded to earlier, the choice was made to vary T_eff continuously within the program while log g was held constant. This was done by averaging together two synthetics from the library in order to estimate a synthetic with an intermediate temperature to the models. The program then used linear interpolation to create weighted averages of these synthetics to fill in the gaps in T_eff and better model the spectra. To estimate the uncertainties associated with assuming a standard log g, separate runs were made for a range of values in T_eff (±500 K) and log g (2.5–4.5 in steps of 0.5) while not letting either be a free parameter in the program. For a given source, the derived v sin i would decrease as log g increased and r_k would decrease as T_eff increased and CO lines weakened. However, these variations had little effect on the pixel shift and derived radial velocity, making the most important parameter from the model fits the most reliable.

In addition to using emcee to model individual sources, it was also used to analyze the entire sample in order to derive values for the mean velocity of the cluster, the intrinsic velocity dispersion, and the binary fraction. This is necessary because single epochs of radial velocity data are likely to have their results affected by the presence of binaries. If a source with a close binary companion happens to be observed in such a way that a component of the binary orbital motion is directed toward or away from Earth, then the observed radial velocity will be offset from the source's true motion through the cluster. The overall effect is to increase the measured radial velocity dispersion of the cluster. For clusters with higher intrinsic velocity dispersions, this effect may be small; however, for clusters with dispersions smaller than a few kilometers per second, the effect can greatly increase the measured dispersion (Kouwenhoven & de Grijs 2008, 2009; Gieles et al. 2010; McConnachie & Côte 2010).

Estimates of the intrinsic velocity dispersion, after corrections for binaries, were made using the program Velbin (Cottaar et al. 2012) in conjunction with emcee. In order to use these two programs together, the following information was required: a set of N radial velocity observations (v_obs,i), each with the associated uncertainties (σ_obs,i) and mass estimates (m_i), and assumptions about the binary period, mass ratio, and eccentricity distributions. Velbin interprets the radial velocity observations as being randomly drawn from a dynamical model describing the intrinsic velocity distribution, which is assumed to be Gaussian. Free parameters for the dynamical model included the mean velocity (v_mean) and the velocity dispersion (v_disp). Then, for a certain subset of stars, a velocity offset v_bin was added to v_obs based on the binary fraction (f_bin). The likelihood function for observing a specific radial velocity for a given star can then be written as

$$\begin{eqnarray}\begin{array}{rcl}{L}_{i}({v}_{\mathrm{obs},i}) & = & (1-{f}_{\mathrm{bin}}){L}_{\mathrm{dyn},i}({v}_{\mathrm{obs},i})\\ & & +\ {f}_{\mathrm{bin}}{\displaystyle \int }_{-\infty }^{+\infty }{L}_{\mathrm{dyn},i}({v}_{\mathrm{obs},i}-{v}_{\mathrm{bin}}){L}_{\mathrm{bin},i}({v}_{\mathrm{bin}}){{dv}}_{\mathrm{bin}}.\end{array}\end{eqnarray} \tag{ 6 }$$

The likelihood function for the entire sample was computed by multiplying the individual likelihood functions. Velbin then determined the log-likelihood function, which was combined with the priors to compute the probability and fed into emcee, which then returned relative probability distributions for how likely the given radial velocity data set is described by v_mean, v_disp, and f_bin.

To measure the radial velocities, all of the spectra obtained were run through the MCMC routine described in Section 3. It was found that to explore adequately the full four-dimensional parameter space (rotational velocity, veiling, radial velocity, and temperature), the program needed 100 walkers to run for 500 steps. The output for each run was then a corner plot of the positions of these walkers along their steps, resulting in relative probability distributions. The corner plots themselves show the probability distribution for each parameter separately along the diagonal, and the off-diagonal plots show correlations between each pair of parameters. Approximately the first 100 steps were excluded from the corner plots, as this was found to be the "burn-in" time for the program to settle on what it thought were the most likely values. This results in the corner plots themselves being approximately Gaussian distributions centered on the most likely values for each parameter. These center values were then determined by finding the value that was higher than 50% of the positions of the walkers and used to create a synthetic to plot over the original data to ensure that the program was returning physically meaningful results. A sample corner plot and overlay of GY 33 taken from an iSHELL spectrum covering the CO v = 2–0 band head are presented in Figures 2 and 3. Parameters derived for YSOs are presented in Table 3, and radial velocities for radial velocity standards are given in Table 2.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As seen in Figure 2, the first two parameters (v sin i and r_k) correspond directly to the physical quantities of the rotational velocity and the veiling. However, dx corresponds to the shift in the spectrum measured in pixels, which is converted to the heliocentric radial velocity given the dispersion of the instrument and corrected for the Earth's motion around the Sun given the date and time of the observation. The last parameter, T_eff, is a measurement of the effective surface temperature of the star, but the actual value of T_eff returned by the program corresponds to the weighted average of two synthetics represented by their array index, as discussed in Section 3. In case of GY 33, the T_eff of 5.8 corresponds to a temperature of about 4260 K. The uncertainties in the fit parameters quoted in Table 3 are from the 16th and 84th percentiles of the marginalized distributions. While the veiling of each source was fit, it was found to be very sensitive to the exact location of the continuum, which is often poorly constrained for these spectra, especially those with low S/Ns or high rotational velocities. For this reason the veiling parameter has been left out of the results in Table 3.

Values for T_eff ranged from 3100 to 5400 K with K or M spectral types, indicating that nearly all of our sources are low-mass YSOs. A large number of sources are fast rotators, with 12 of the 32 objects exhibiting values for v sin i > 30 km s⁻¹. There were two NIRSPEC sources that were outliers in radial velocity: WL 1 and WL 19. We have no explanation for their large negative velocities but note that WL 1 is a known subarcsecond binary and WL 19 is suspected to be a heavily reddened Class III object behind the molecular cloud (see the Appendix). A histogram of these radial velocities with 2 km s⁻¹ bins is shown in Figure 4. A Gaussian fit to the velocities from −12 to +2 km s⁻¹ indicates a velocity dispersion of 2.66 ± 0.70 km s⁻¹. This is significantly higher than what was expected for this cluster based on previous studies. It is also clear from this figure that the distribution itself is marginally Gaussian, making the calculation of a velocity dispersion uncertain. As described in Section 3, the program Velbin was used to analyze these data and attempted to correct for the presence of binaries (see Section 4.3).

Zoom In Zoom Out Reset image size

4.1.1. Error Analysis

4.1.2. Comparisons with Published Data

4.1.3. Radial Velocity Variables

Having established the radial velocity dispersion for 30 sources within the L1688 cloud (Figure 4), the next step was to assess the broadening of the distribution due to the presence of binary companions. The program Velbin was used in order to infer the cluster's intrinsic velocity dispersion and mean velocity assuming a binary fraction (see Section 3.2). In order to use Velbin, the program required a radial velocity, a radial velocity error, and a mass estimate for each source. Masses for YSOs are difficult to estimate, especially for deeply embedded objects with large infrared excesses. Given the log T_eff and distance to L1688, luminosity estimates were made for lightly veiled sources by dereddening the Two Micron All Sky Survey J magnitude (Skrutskie et al. 2006) using the extinction law of Yuan et al. (2013) and intrinsic colors for main-sequence stars and applying a bolometric correction (Pecaut & Mamajek 2013). Mass estimates were made by comparing log T_eff and luminosities to several pre-main-sequence models (D'Antona & Mazzitelli 1997; Palla & Stahler 1999; Siess et al. 1999). Published mass estimates for eight objects from this study ranged from 0.2 to 1.7 M_⊙, and our mass estimates agreed with these to within a factor of two (Ressler & Barsony 2001; Correia et al. 2006; Erickson et al. 2011; Rigliaco et al. 2016; Simon et al. 2017). For more heavily veiled sources, we assumed that the mass scales roughly with log T_eff and assigned masses of 1.5, 1.0, 0.5, or 0.2 M_⊙ accordingly.

Velbin was run using a mass ratio and period distribution derived from solar-type field stars (Raghavan et al. 2010; Reggiani & Meyer 2013) and a flat eccentricity distribution. Assuming a binary fraction for Class I and FS YSOs of 0.5 (Connelley et al. 2008) yielded values for v_mean and v_disp of −4.8 ± 0.6 km s⁻¹ and ${2.81}_{-0.49}^{+0.58}$ km s⁻¹, respectively. To account for the possibility that the binary fraction for deeply embedded sources may be close to 1, Velbin was also run setting f_bin = 1, resulting in the largest possible correction for binaries and a v_disp = ${2.39}_{-0.64}^{+0.68}$ km s⁻¹. These results are summarized in Table 4. Because of the large uncertainties in mass estimates, Velbin was run on the same data set, but with the masses set at constant values of 0.5, 1.0, or 1.5 M_⊙. The results of these runs are not shown here, but apart from a marginally smaller dispersion for the highest mass assumed, the dispersions were virtually indistinguishable from those presented in Table 4 for this study.

Table 4.  1D Velocity Dispersions of YSOs in L1688

Type of Study	No. of Sources	SED Class	Velocity Dispersion	Notes	References
			(km s−1)
Relative proper motion (R.A./decl.)	29	Class I/FS	0.84 ± 0.10, 0.71 ± 0.09	Cores A, B, E/F	(1)
Relative proper motion (R.A./decl.)	28	Class II/III	0.88 ± 0.08, 1.28 ± 0.26	Cores A, B, E/F	(1)
Absolute proper motion (R.A.)	68	all Classes	3.30 ± 0.42	Cores E/F	(2)
Absolute proper motion (decl.)	68	all Classes	2.53 ± 0.32	Cores E/F	(2)
Absolute proper motion (R.A.)	100	Class II/III	0.98 ± 0.10	all of L1688	(3)
Absolute proper motion (decl.)	100	Class II/III	1.08 ± 0.14	all of L1688	(3)
Radial velocity	47	Class II/III	1.14 ± 0.35	outskirts of L1688 fbin = 0.56	(4)
Radial velocity	30	Class I/FS	2.66 ± 0.70	Cores A, B, E/F	(5)
Radial velocity	30	Class I/FS	${2.81}_{-0.64}^{+0.68}$	Cores A, B, E/F fbin = 0.5	(5)
Radial velocity	30	Class I/FS	${2.39}_{-0.49}^{+0.58}$	Cores A, B, E/F fbin = 1	(5)

References. (1) Wilking et al. 2015; (2) Ducourant et al. 2017; (3) Brown et al. 2018; (4) Rigliaco et al. 2016; (5) this study.

Download table as:  ASCIITypeset image

Notes regarding individual sources can be found in the Appendix. While the main focus of this project was to derive precise radial velocities from the data, there is a wealth of information that can be obtained from the iSHELL spectra. Brγ emission was detected in the spectra of VSSG 1, SR 24S, GY 235, WLY 2-51, and WLY 2-54. In addition, emission from shocked molecular hydrogen at 2.12 μm was detected in SR 24S, GY 235, and WLY 2-54. Prominent features in the spectrum of WLY 2-54, and present in the spectra of 10 other YSOs, were very narrow interstellar absorption lines from the low R and P branch of the CO v = 2–0 transition that can in principle reveal conditions of the foreground gas. A detailed analysis of these features is beyond the scope of this paper.

The inability of Velbin to derive an intrinsic radial velocity dispersion for our sample that was significantly smaller than that observed could arise for several reasons. Velbin makes its calculations by assuming that the intrinsic dispersion is Gaussian in nature, and as can be seen from the histogram of velocities (Figure 4), the radial velocity distribution is only vaguely Gaussian in shape. To assess the Gaussian nature of the radial velocity distribution, we applied the Shapiro–Wilk test of normality to this sample with the two outliers between −25 and −30 km s⁻¹ removed. The test yielded a p-value = 0.054, which for a significance level of α = 0.05 means that we cannot reject the null hypothesis that the radial velocities were drawn from a normal distribution, albeit only by a small margin. It is also possible that the assumed mass and period distributions of binaries provided by Velbin for main-sequence stars do not match up with less evolved YSOs. Perhaps a greater limitation of Velbin is that it does not consider any triple- or higher-order systems, which could be common among pre-main-sequence objects. Indeed, there are a few known triple systems in this sample (e.g., WL 4, WL 20, WLY 2-44) and a possible quadruple system (SR 24). For example, Stassun et al. (2014) found that 6 of 13 pre-main-sequence eclipsing binaries they studied were, in fact, triple systems, implying that an overall fraction of triple systems for young stars could be on the order of 0.25.

Velocity dispersions from recent proper-motion surveys in L1688, as well as the Gaia DR2 data release, can be compared with those derived by this study. One advantage of proper-motion studies is that the effect of binaries on the dispersion is minimized since, over decade-long baselines, any effects due to binaries would either be averaged out over several periods or be too small to have a significant impact on the measured motion. Two proper-motion surveys were conducted in the infrared and hence not biased against less evolved YSOs. Ducourant et al. (2017) observed dense cores E and F in L1688 with an astrometric K_s filter over a 02 × 02 area over a 9+ yr time baseline. These data were complemented with infrared images from publicly available databases that encompassed a larger area. Absolute proper motions were obtained for 68 high-probability members with an accuracy of ∼1.2 km s⁻¹; about 40% of the sample exhibited Class I or FS SEDs. Gaussian fits to these data yielded 1D velocity dispersions for these objects in R.A. and decl. of 3.3 ± 0.4 km s⁻¹ and 2.9 ± 0.3 km s⁻¹, respectively. The authors note that the lack of a near-infrared reference frame and the inhomogeneous distribution of stars may introduce systematic errors that are unaccounted for. However, on face value, the 1D velocity dispersions are consistent with those from our radial velocity study.

Wilking et al. (2015) monitored four fields in three subclusters (dense cores A, B, and E/F) in the λ = 2.1 μm filter over a 12 yr period. In the absence of an absolute reference frame of background stars, relative proper motions were derived for 65 YSOs with median uncertainties in R.A. and decl. of 0.39 km s⁻¹ and 0.46 km s⁻¹, respectively. Assuming that each field had the same mean proper motion, 1D velocity dispersions in R.A. and decl. were found to be 0.88 ± 0.10 km s⁻¹ and 0.96 ± 0.09 km s⁻¹, respectively. Furthermore, the sample could be split up nearly evenly into two groups by SED class: 29 Class I and FS objects and 28 Class II/III objects. As shown in Table 4, these groups each had relative 1D velocity dispersions derived from their R.A. and decl. motions, which were also consistent with ∼1 km s⁻¹. To reconcile these results with the 1D radial velocity dispersion of this study and that of Ducourant et al. (2017), the assumption that the average proper motion in each subcluster is the same would not be correct.

One hundred sources were extracted from the Gaia DR2 data release covering a square degree centered on L1688 with distances in the range of 137 ± 5 pc. The median uncertainties in the absolute proper motions in R.A. and decl. were 0.15 km s⁻¹ and 0.10 km s⁻¹, respectively. The areal extent of the sample overlaps well with the Rigliaco et al. radial velocity survey, and as it is an optical survey, it is also heavily biased toward low-extinction Class III YSOs. Not surprisingly, the 1D velocity dispersions of ∼1 km s⁻¹ are consistent with that derived by the Rigliaco et al. survey.

As summarized in Table 4, the results of radial velocity and absolute proper-motion surveys appear to show differences in the 1D velocity dispersions between samples of Class I/FS YSOs compared to Class II/III dominated samples at the 2σ level. These samples also differ depending on whether they are concentrated toward the dense molecular cores or encompass a much larger areal extent that includes the low-extinction outskirts or surface layers of L1688. Future infrared proper-motion surveys that can establish a reliable absolute frame of reference will be able to determine whether these differences are of statistical significance. Below we explore several possibilities for the apparent higher velocity dispersions of less evolved YSOs.

5.3.1. Underestimating the Binary Fraction

5.3.2. Supervirial State