# Size and structures of disks around very low mass stars in the Taurus star-forming region

N. T. Kurtovic, P. Pinilla, F. Long, M. Benisty, C. F. Manara, A. Natta, I.

Pascucci, L. Ricci, A. Scholz, L. Testi

► To cite this version:

N. T. Kurtovic, P. Pinilla, F. Long, M. Benisty, C. F. Manara, et al. Size and structures of disks around very low mass stars in the Taurus star-forming region. Astronomy and Astrophysics - A&A, 2021, 645, 10.1051/0004-6361/202038983. insu-03636644

HAL Id: insu-03636644

https://insu.hal.science/insu-03636644

Submitted on 12 Apr 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Size and structures of disks around very low mass stars in the Taurus star-forming region*

N. T. Kurtovic<fnref n="1" />, P. Pinilla<fnref n="1" />, F. Long<fnref n="2" />, M. Benisty<fnref n="3" />,<fnref n="4" />, C. F. Manara<fnref n="5" />, A. Natta<fnref n="6" />, I. Pascucci<fnref n="7" />,<fnref n="8" />, L. Ricci<fnref n="9" />, A. Scholz<fnref n="10" />, and L. Testi<fnref n="5" />,<fnref n="11" />

Received 20 July 2020 / Accepted 16 November 2020

## ABSTRACT

**Context.** The discovery of giant planets orbiting very low mass stars (VLMS) and the recent observed substructures in disks around VLMS is challenging planet formation models. Specifically, radial drift of dust particles is a catastrophic barrier in these disks, which prevents the formation of planetesimals and therefore planets.

**Aims.** We aim to estimate if structures, such as cavities, rings, and gaps, are common in disks around VLMS and to test models of structure formation in these disks. We also aim to compare the radial extent of the gas and dust emission in disks around VLMS, which can give us insight about radial drift.

**Methods.** We studied six disks around VLMS in the Taurus star-forming region using ALMA Band 7 (~340 GHz) at a resolution of ~0.1". The targets were selected because of their high disk dust content in their stellar mass regime.

**Results.** Our observations resolve the disk dust continuum in all disks. In addition, we detect the ¹²CO (J = 3–2) emission line in all targets and ¹³CO (J = 3–2) in five of the six sources. The angular resolution allows the detection of dust substructures in three out of the six disks, which we studied by using UV-modeling. Central cavities are observed in the disks around stars MHO 6 (M 5.0) and CIDA 1 (M 4.5), while we have a tentative detection of a multi-ringed disk around J0433. We estimate that a planet mass of ~0.1 M<sub>Jup</sub> or ~0.4 M<sub>Saturn</sub> is required for a single planet to create the first gap in J0433. For the cavities of MHO 6 and CIDA 1, a Saturn-mass planet (~0.3 M<sub>Jup</sub>) is required. The other three disks with no observed structures are the most compact and faintest in our sample, with the radius enclosing 90% of the continuum emission varying between ~13 and 21 au. The emission of ¹²CO and ¹³CO is more extended than the dust continuum emission in all disks of our sample. When using the ¹²CO emission to determine the gas disk extension R<sub>gas</sub>, the ratio of R<sub>gas</sub>/R<sub>dust</sub> in our sample varies from 2.3 to 6.0. One of the disks in our sample, CIDA 7, has the largest R<sub>gas</sub>/R<sub>dust</sub> ratio observed so far, which is consistent with models of radial drift being very efficient around VLMS in the absence of substructures.

**Conclusions.** Given our limited angular resolution, substructures were only directly detected in the most extended disks, which represent 50% of our sample, and there are hints of unresolved structured emission in one of the bright smooth sources. Our observations do not exclude giant planet formation on the substructures observed. A comparison of the size and luminosity of VLMS disks with their counterparts around higher mass stars shows that they follow a similar relation.

**Key words.** accretion, accretion disks – circumstellar matter – stars: pre-main sequence – protoplanetary disks – planets and satellites: formation

## 1. Introduction

For every ten stars that are formed in the Milky Way, around two to five brown dwarfs (BDs) also form (e.g., Scholz et al. 2012; Mužić et al. 2019), and M-dwarfs represent about three-quarters

of all the stars in our galaxy. Exoplanet discoveries show that short-period (<50 days) sub-Neptune planets occur more frequently around M-dwarfs than around FGK stars (Mulders et al. 2015), but also a few giant planets have been discovered around BDs and very low mass stars (VLMS $\le$ 0.1 M<sub>☉</sub>, e.g., Morales et al. 2019). This implies that planets of a large range of masses can form around these objects, although it remains an open question if these massive objects form as binary companions of

<footnote n="1">¹ Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany
e-mail: kurtovic@mpia.de; pinilla@mpia.de</footnote>

<footnote n="2">²Center for Astrophysics | Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA</footnote>

<footnote n="3">³Unidad Mixta Internacional Franco-Chilena de Astronomía, CNRS/INSU UMI 3386, Departamento de Astronomía, Universidad de Chile, Camino El Observatorio 1515, Las Condes, Santiago, Chile</footnote>

<footnote n="4">⁴Univ. Grenoble Alpes, CNRS, IPAG, 38000 Grenoble, France</footnote>

<footnote n="5">⁵European Southern Observatory, Karl-Schwarzschild-Str. 2, 85748 Garching, Germany</footnote>

<footnote n="6">⁶Dublin Institute for Advanced Studies, School of Cosmic Physics, 31 Fitzwilliam Place, Dublin 2, Ireland</footnote>

<footnote n="7">⁷Lunar and Planetary Laboratory, The University of Arizona, Tucson, AZ 85721, USA</footnote>

<footnote n="8">⁸Earths in Other Solar Systems Team, NASA Nexus for Exoplanet System Science, Tucson, Arizona, USA</footnote>

<footnote n="9">⁹Department of Physics and Astronomy, California State University Northridge, 18111 Nordhoff Street, Northridge, CA 91130, USA</footnote>

<footnote n="10">¹⁰SUPA, School of Physics & Astronomy, University of St. Andrews, North Haugh, St. Andrews KY16 9SS, UK</footnote>

<footnote n="11">¹¹INAF - Arcetri, Largo E. Fermi 5, 50125 Firenze, Italy</footnote>

<footnote>* The reduced images and datacubes are only available at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsarc.u-strasbg.fr/viz-bin/cat/J/A+A/645/A139</footnote>
the BDs and very low mass stars (VLMS), or as planets. Current models of planet formation through pebble or planetesimal accretion cannot explain the formation of giant planets around VLMS (Liu et al. 2020).

Observations of BDs and VLMS from the near-infrared to centimeter wavelength reveal the existence of circumstellar disks around these objects (e.g., Luhman 2006; Klein et al. 2003; Scholz et al. 2006, 2007; Ricci et al. 2012, 2017a,b; Daemgen et al. 2016; van der Plas et al. 2016; Sanchis et al. 2020), which are more compact and lower in dust mass when compared to disks around T-Tauri stars (e.g., Pinilla et al. 2017a; Ward-Duong et al. 2018; Hendler et al. 2017, 2020). The typical millimeter fluxes of such disks suggest that they have dust masses of few Earth masses, challenging the formation of giant planets through core or pebble accretion (Liu et al. 2020).

The core accretion scenario for planet formation assumes collisional growth from sub-µm-sized dust particles from the interstellar medium (ISM) to kilometer-sized bodies or planetesimals (e.g. Pollack et al. 1996). The collisions of particles and their dynamics within the disk are regulated by the interaction with the surrounding gas. Different physical processes lead to collisions of particles and their potential growth, such as Brownian motion, turbulence, dust settling, and radial drift (e.g., Brauer et al. 2008). All of these processes have a direct or indirect dependency on the properties of the hosting star, such as the temperature and mass. For instance, from theoretical calculations, settling and radial drift are expected to be more efficient in disks around VLMS and BDs, with BD disks being 15–20% flatter and with radial drift velocities being twice as high or even more in these disks compared to T-Tauri disks (Mulders & Dominik 2012; Pinilla et al. 2013).

With radial drift being a more pronounced problem in disks around BDs and VLMS, it is still unknown how this barrier of planet formation is overcome in these environments where the disks are more compact, colder, and have a lower mass. Millimeter-sized particles have been detected in BD and VLMS disks through measurements of the spectral index (Ricci et al. 2014; Pinilla et al. 2017a), which are only possible to explain when radial drift is significantly reduced by the presence of strong pressure bumps (Pinilla et al. 2013). The presence of pressure bumps produces substructures, such as rings, gaps, spiral arms, and lopsided asymmetries, with a different amplitude, contrasts, and locations depending on the origin of the pressure variations (e.g., Pinilla & Youdin 2017; Andrews 2020). Currently, due to sensitivity limitations, most of our observational knowledge about substructures comes from bright (and probably massive) disks, such as the DSHARP sample (Andrews et al. 2018a). A less biased sample of ALMA observations of disks in the star-formation region of Taurus has demonstrated that at least 33% of disks host substructures at a resolution of 0.1", and the disks that do not have any substructures are compact (dust disk radii lower that ~50 au, Long et al. 2018, 2019). It remains an open question if compact disks are small because they lack pressure bumps or because current observations lack the resolution to detect rings and gaps in these disks (the scale of a radial pressure bump cannot be smaller than one local scale height, e.g., Dullemond et al. 2018).

Pinilla et al. (2018a) have thus far reported the lowest mass star with a resolved large dust cavity (radius ~20 au) in the disk around the M4.5 star CIDA 1. In the context of planets creating such a cavity, it has been shown that a high planet-to-stellar mass ratio is needed to open a gap and trap particles in disks around VLMS, because in these cases, the disk scale height at a given location is higher than in moderate or high mass stars

(Pinilla et al. 2017a; Sinclair et al. 2020). In a typical disk around a VLMS as CIDA 1, at least a Saturn-mass planet is needed to open a gap in the disk. This challenges our current understanding of substructures and the common idea that planets are responsible for their formation, since these disks around VLMS and BDs may not have enough mass to form such massive planets, although they may form from gravitational instability if the disks were much denser in their early stages (e.g., Mercer & Stamatellos 2020).

Based on the previous CIDA 1 observations, we selected a sample of five disks to observe with ALMA at a resolution of 0.1" in the Taurus star-forming region, whose properties are similar to CIDA 1. Specifically, these disks around low mass stars are more massive compared to other disks with hosts in the same stellar regime. These observations included ¹²CO and ¹³CO and aim to estimate how common substructures are around VLMS. In addition, as radial drift is expected to be very efficient in these disks, they are excellent laboratories to search for the difference between the radial extent of the gas and dust in disks, which can be a direct signature of radial drift (e.g., Trapman et al. 2019).

This paper is organized as follows: Sect. 2 summarizes the target selection, the ALMA observations, and the calibration of the data. Section 3 presents the modeling of the data in the visibility plane for the continuum as well as in the image plane for the ¹²CO and ¹³CO emission. In Sect. 4, we discuss our results in the context of different origins for the seen substructures as well as the observed difference between the radial extent of the gas and dust. Section 5 summarizes the conclusions of this paper.

## 2. Target selection and observations

### 2.1. Target selection

For the new ALMA observations, we selected five disks around VLMS in the Taurus star-forming region. Given that this is the first small survey of high angular resolution observations of disks around VLMS, the sample was selected to optimize the chances of finding substructures. The criteria used to select the targets were as follows: (1) the target has been previously observed and detected by either SMA or ALMA in millimeter wavelengths (based on the observations by Andrews et al. 2013 and Ward-Duong et al. 2018, respectively); (2) it has a stellar mass between ~0.1 and 0.2 $M_\odot$; and (3) it has a high disk dust mass compared to the stellar mass. The last condition comes from observations of the $M_{\text{dust}}-M_\star$ relation of transition disks and disks with substructures (Pinilla et al. 2020), which shows that those disks usually have higher dust masses compared to others with stellar hosts of a similar mass. From the list of targets that fulfilled the conditions, we selected the ones with the lowest optical extinction.

### 2.2. Observations

For sample completeness, we included the archival data of CIDA 1 in our list of targets, given that combining these previous archival observations results in a dataset with a similar sensitivity and angular resolution. For consistency with future works, the target properties used in this study (shown in Table 1) were taken from Herczeg & Hillenbrand (2014), with stellar masses derived following the method described in Pascucci et al. (2016), using distances inferred from Gaia DR2 (Gaia Collaboration 2018). The distances in parsecs was calculated as the inverse of the parallax. The spatial distribution of our sample in the Taurus optical extinction map is shown in Fig. A.1.
<table><caption>**Table 1.** Source properties used in this work.</caption><thead><tr><th>2MASS name</th><th>Used name</th><th>Spectral type</th><th>M<sub>⋆</sub> [M<sub>☉</sub>]</th><th>T<sub>eff</sub> [K]</th><th>L<sub>⋆</sub> [L<sub>☉</sub>]</th><th>d [pc]</th></tr></thead><tbody><tr><td>J04141760+2806096</td><td>CIDA 1</td><td>M4.5</td><td>0.19</td><td>3197</td><td>0.20</td><td>135.7</td></tr><tr><td>J04322210+1827426</td><td>MHO 6</td><td>M5.0</td><td>0.17</td><td>3125</td><td>0.06</td><td>141.9</td></tr><tr><td>J04334465+2615005</td><td>J0433</td><td>M5.2</td><td>0.15</td><td>3098</td><td>0.12</td><td>173.3</td></tr><tr><td>J04422101+252034</td><td>CIDA 7</td><td>M5.1</td><td>0.15</td><td>3111</td><td>0.08</td><td>136.2</td></tr><tr><td>J04202555+2700355</td><td>J0420</td><td>M5.25</td><td>0.14</td><td>3091</td><td>0.07</td><td>170.4</td></tr><tr><td>J04155799+2746175</td><td>J0415</td><td>M5.2</td><td>0.15</td><td>3098</td><td>0.05</td><td>135.7</td></tr></tbody></table>

**References.** Spectral type, *T*<sub>eff</sub>, and *L*<sub>☉</sub> comes from Herczeg & Hillenbrand (2014). Stellar masses derived following the Pascucci et al. (2016) method, using distances inferred from Gaia DR2 (Gaia Collaboration 2018).

These targets were observed with ALMA at 0.87 mm (Band 7) as part of the project 2018.1.00310.S (PI: P. Pinilla), during Cycle 6, with the spectral setup configured to observe in four spectral windows centered at 331.3, 333.3, 344.0, and 345.8 GHz, with two of them centered to observe dust continuum emission, and two observing molecular line emission from ¹²CO (J = 3–2) and ¹³CO (J = 3–2). The frequency resolution for ¹²CO was 244.1 kHz per channel; while for ¹³CO and the continuum, it was 976.6 kHz. The most extended antenna configuration used was C43-8, providing an angular resolution of 0.08" at best. Some of our sources had archival data from ALMA project 2012.1.00743.S (PI: G. van der Plas) and 2016.1.01511.S (PI: J. Patience), with observations of ¹²CO and dust continuum, which were also included in the self-calibration and analysis. Archival Band 7 data of CIDA 1 was observed by the project 2015.1.00934.S (PI: L. Ricci) published in Pinilla et al. (2018a), and in this work we combined it with 2016.1.01511.S (PI: J. Patience). A summary of the observation details and the data considered for each target are shown in Table A.1.

The raw datasets were calibrated by applying the ALMA pipeline using the CASA version specified for each project (McMullin et al. 2007). Then, we used CASA v5.6.2 for the subsequent data handling and imaging. We extracted the dust continuum emission of every source by flagging the channels closer than 25 km s⁻<fnref n="1" /> to the targeted molecular lines. To reduce the data volume, we averaged over time (6 s intervals) and channels (with a width of 125 MHz). Before combining all available observations for each source, the centroid spatial position of the emission was determined by fitting a Gaussian using the imfit task, and shifted using fixvis and fixplanets tasks to the centroid of the observations of extended baselines, shown in Table A.2. We also checked for a consistent flux calibration by comparing the amplitude of the emission in different executions. We found a discrepancy of 12% in the fourth compact observation of J0433 (2018-11-24), which we rescaled to match all the others.

In order to boost the signal-to-noise ratio (S/N) of each source, we performed self-calibration of the datasets in two steps: first we combined the compact configuration observations and performed self-calibration. Second, we combined those self-calibrated observations with the extended configuration observations and then self-calibrated again. Phase calibrations were applied until the improvement on the S/N was below ~5%. Only 1 amplitude calibration was applied in each step. The overall S/N improvement was between 1.5 ~ 4.0, depending on the source. The only source where self-calibration was not possible was J0415, because the initial S/N of 9 was too low for improvements. The final continuum images were generated using a Briggs

robust parameter of 0.0 for CIDA 1 and 0.5 for the remaining sources. The image properties are summarized in Table A.2.

All the steps for the dust continuum emission calibrations, including centroid shifting, flux calibration, and self-calibration tables, were then applied to the molecular line emission channels. The continuum emission was subtracted using the uvcontsub task, and the images were generated using a robust parameter of 1. In MHO 6, J0420, and J0415, UV-tapping was applied in order to increase the S/N of the gas images. For MHO 6, we applied a UV-tapping of 0.13" on the ¹²CO; while for CIDA 7 and J0415, we applied a UV-tapping of 0.1" in both molecular line images. All our scripts of self-calibration and imaging are available online<fnref n="1" />.

The final images of the dust continuum, moment 0, and moment 1 of the ¹²CO and ¹³CO (when detected) are shown in Fig. 1, and the velocity channel maps are in Appendix B. The details of the dust continuum and CO images can be found in the appendix as well, summarized in Tables A.2 and A.3, respectively.

# 3. Results

## 3.1. Morphology of the continuum emission

Our sources are spatially compact (radius of 0.1" ~ 0.3"), and both their radial extent and substructures have sizes comparable to the synthesized beam shape. To avoid image reconstruction biases, we fit the deprojected brightness profile of the sources in the UV-plane. The continuum visibilities were extracted from the self-calibrated measurement set, and we used the central frequency of each channel to convert the UV-coordinates to wavelength units. We started by modeling every source with a central Gaussian profile, and then we increased the complexity of the profile if the residuals suggested it. We also guided our parametrization of the profiles based on the best fitting from frank (Jennings et al. 2020), which fits a nonparametric 1D radial brightness profile in the visibilities, using Gaussian processes. For CIDA 7, J0420, and J0415, the function that describes their brightness profile is a centrally peaked Gaussian profile, following an intensity given by:

$$I_g(r) = f_1 \exp \left( -\frac{r^2}{2\sigma_1^2} \right), \quad (1)$$

where $I_g$ is the Gaussian intensity profile of the source as a function of the radius $r$.

<footnote n="1">¹ https://github.com/nicokurtovic/VLMS_ALMA_2018.1.00310.S</footnote>
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** ALMA observations of the disks around VLMS in Taurus. *From left to right:* dust continuum emission, <sup>12</sup>CO moments 0 and 1, and <sup>13</sup>CO moments 0 and 1. All boxes are 3.0" in size. The scale bar represents 20 au, and the white ellipses show the synthesized beam. The centers, beam sizes, measured fluxes, and sensitivity are detailed in Tables A.2 and A.3. Dashed lines show the region used to calculate the radial profiles and to measure R<sub>68</sub> and R<sub>90</sub>. A zoomed-in version of the continuum images can be found in Fig. 2.</figcaption></figure>
For CIDA 1 and MHO 6, we modeled the disk with a radially asymmetric Gaussian ring or a broken Gaussian from hereafter, that is to say the inner and outer width of the ring can differ. This profile is motivated by results of radially asymmetric accumulation of particles in pressure bumps (see e.g., Pinilla et al. 2015, 2017b). Such radially broken Gaussian profiles have been used to describe the morphology of different rings in transition disks and disks with substructures (e.g., Pinilla et al. 2018b; Huang et al. 2020), which is the same model used in Pinilla et al. (2018a) to model CIDA 1. The intensity profile is given by a ring as follows:

$$I_{bg}(r) = \begin{cases} f_1 \exp\left(-\frac{(r-r_1)^2}{2\sigma_1^2}\right) & \text{for } r \le r_1 \\ f_1 \exp\left(-\frac{(r-r_1)^2}{2\sigma_2^2}\right) & \text{for } r > r_1 \end{cases}, \quad (2)$$

where $I_{bg}$ is the broken Gaussian intensity profile as a function of the radius, $r_1$ is the radial location of the ring peak intensity, and $\sigma_{1,2}$ are the Gaussian widths for the inner and outer sides of the ring, respectively.

Finally, for J0433, the profile is the sum of a centrally peaked Gaussian profile and two symmetric Gaussian rings, as suggested by frank. It is also the profile that creates the lowest amount of residuals from our experiments, such as the single Gaussian, Gaussian ring, and broken Gaussian ring. The intensity profile is

$$I_{J0433}(r) = \sum_{i=1}^{3} f_i \exp\left(-\frac{(r-r_i)^2}{2\sigma_i^2}\right), \quad (3)$$

where $r_i=1=0$, so the first Gaussian is peaked at the center.

The visibilities of each profile were computed by combining each model with a spatial offset ($\delta_{RA}$, $\delta_{Dec}$), inclination (inc), and position angle (PA). Therefore, each model has four more free parameters in addition to those that describe the intensity profile. The Fourier transform and the $\chi^2$ calculation were carried out with the galario package (Tazzari et al. 2018). The pixel size used in the models is 1 mas, which is several times smaller than the smallest resolvable scale of the observations. The $\chi^2$ was scaled up by a factor of 2.667 since CASA does not account for the effective channel width, introduced by Hanning smoothing, when it averages the weights during data binning. We adopted a uniform prior probability distribution over a wide parameter range, such that walkers would only be initially restricted by geometric considerations (inc $\in$ [0, 90], PA $\in$ [0, 180], $\sigma \ge 0$).

We used a Monte Carlo Markov chain (MCMC) routine based on the emcee package (Foreman-Mackey et al. 2013) to sample the posterior probability distribution of each parameter space. Furthermore, we ran more than 250 000 steps after converging to find the most likely set of parameters and the error bars, taken from the 16th and 84th percentile.

Our results for each parameter are shown in Table 2, while in Fig. 2 we show the models with and without convolution (right most panel). The visibilities and radial profile were deprojected using the best inclination and position angle. The residual image was generated in CASA using the same parameters and procedure used for the observations, from a measurement set with its visibilities calculated by subtracting the best model from the data.

The radial profile recovered from the UV-modeling was used to measure the dust continuum emission radii ($R_{dust}$) that encloses 68 and 90% of the total flux (the dust $R_{68}$ and $R_{90}$, respectively). We used the different sets of parameters sampled by the walkers to compute their radial profiles, and we found

the 16th and 84th percentile on the $R_{68}$ and $R_{90}$ to calculate the error bars. The results are shown in Table 3 and in Fig. 3. The continuum emission in our sample extends up to 46 au in the biggest disks, MHO 6, and J0433, while we also have the smallest (13 au) and the dimmest (~1 mJy) Taurus disks ever resolved in 0.87 mm emission, CIDA 7, and J0415, respectively.

## 3.2. Rings and cavities in CIDA 1, MHO 6, and J0433

We find evidence of ring structures in three out of six disks in our sample. At the current resolution, CIDA 1 and MHO 6 show a single ring and a cavity, which are well described by a broken Gaussian profile. In CIDA 1, we find the inner side of the Gaussian to be more extended than the outer side, which is similar to the results obtained in Pinilla et al. (2018a). The ratio between the widths of the inner side and the outer side in our analysis is 1.2, and the peak of our ring model is located at 20.8 au.

For MHO 6, we find that a strong asymmetric ring peaked at 10 au, with its inner side having converged to $\sigma=0$ au with a narrow error margin. To ensure that our result was not being affected by numerical biases related to the pixel size, we ran another MCMC with a pixel size of 0.4 mas (about 0.057 au), and we consistently recovered the same result for each parameter. Even though this steep transition could suggest an unresolved inner side of the ring, the best model is also driven by the non-axisymmetric emission of the disk, which is only noticeable when looking at the residuals (see Fig. 2). The contrast of these asymmetries is about 5% of the peak amplitude of the ring, meaning that most of the emission is still well described by a radially axisymmetric ring.

In J0433, our best model finds two rings located at 21 and 32 au, with gaps at 16 and 25 au. The brightness ratio between the first ring and first gap is about 4.4; whereas, the contrast is 1.2 between the second gap and ring. Since the brightness ratio fades at the $1\sigma$ error (as seen in Fig. 3), we cannot exclude that it is an extension of the outer side of the first ring.

## 3.3. Dust disk masses

Assuming the dust continuum emission is optically thin, we estimated the disk dust mass by following Hildebrand (1983):

$$M_{\text{dust}} \approx \frac{d^2 F_{\nu}}{\kappa_{\nu} B_{\nu}(T)}, \quad (4)$$

where d is the distance in parsecs (given in Table 1, taken from Gaia DR2), $\kappa_{\nu}$ is the mass absorption coefficient, and $B_{\nu}(T)$ is the Planck function. For $\kappa_{\nu}$, we assumed the opacity law $\kappa_{\nu} = 2.3 \text{ cm}^2 \text{ g}^{-1} (\nu/230 \text{ GHz})^{0.4}$ (Andrews et al. 2013), while the temperature was assumed constant at $T_{\text{dust}} = 20 \text{ K}$ (e.g., Ansdell et al. 2016; Pinilla et al. 2018a). We estimated the uncertainty of the dust mass by taking the 10% uncertainty in the flux calibration. Our results for each source are compiled in Table A.2 and shown compared with the disk dust mass in the Taurus region in the Fig. 4. The lowest dust mass estimate was obtained for J0415, with only $0.22 \pm 0.02 M_{\oplus}$ (about $2.1 \pm 0.2 M_{\text{Mars}}$), while the most massive dust disk was detected in J0433 with $14.32 \pm 1.43 M_{\oplus}$.

## 3.4. $^{12}$CO and $^{13}$CO emission

### 3.4.1. Lines detection in each source

Both molecular lines were detected in all of our sources, except for J0415, where only $^{12}$CO was observed. As shown in Figs. 1
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Visibility modeling versus observations of our sample. *From left to right:* (1) real part of the visibilities after centering and deprojecting the data versus the best fit model of the continuum data, (2) continuum emission of our sources where the scale bar represents a 10 au distance, (3) model image, (4) residual map (observations minus model), (5) and normalized, azimuthally-averaged radial profile calculated from the beam convolved images in comparison with the model without convolution (purple solid line) and after convolution (red solid line). In the right most plots, the gray scale shows the beam major axis.</figcaption></figure>
<table><caption>**Table 2.** Best parameters from UV-modeling, following Eqs. (1)–(3).</caption><thead><tr><th></th><th>CIDA 1</th><th>MHO 6</th><th>J0433</th><th>CIDA 7</th><th>J0420</th><th>J0415</th><th>Unit</th></tr></thead><tbody><tr><td>δ<sub>RA</sub></td><td>-0.18<sup>+0.09</sup><sub>-0.14</sub></td><td>3.84<sup>+0.15</sup><sub>-0.09</sub></td><td>3.73<sup>+0.07</sup><sub>-0.09</sub></td><td>-1.18<sup>+0.06</sup><sub>-0.06</sub></td><td>2.87<sup>+0.08</sup><sub>-0.09</sub></td><td>-4.19<sup>+4.04</sup><sub>-2.68</sub></td><td>mas</td></tr><tr><td>δ<sub>Dec</sub></td><td>-5.16<sup>+0.11</sup><sub>-0.20</sub></td><td>-2.69<sup>+0.11</sup><sub>-0.07</sub></td><td>-3.91<sup>+0.19</sup><sub>-0.04</sub></td><td>-0.06<sup>+0.06</sup><sub>-0.06</sub></td><td>-4.42<sup>+0.07</sup><sub>-0.10</sub></td><td>0.28<sup>+3.32</sup><sub>-4.24</sub></td><td>mas</td></tr><tr><td>inc</td><td>38.2<sup>+0.15</sup><sub>-0.06</sub></td><td>64.56<sup>+0.06</sup><sub>-0.01</sub></td><td>57.62<sup>+0.01</sup><sub>-0.12</sub></td><td>31.35<sup>+0.30</sup><sub>-0.28</sub></td><td>38.24<sup>+0.26</sup><sub>-0.24</sub></td><td>34.96<sup>+1.21</sup><sub>-28.18</sub></td><td>deg</td></tr><tr><td>PA</td><td>11.2<sup>+0.18</sup><sub>-0.16</sub></td><td>113.55<sup>+0.04</sup><sub>-0.05</sub></td><td>165.26<sup>+0.03</sup><sub>-0.11</sub></td><td>85.95<sup>+0.53</sup><sub>-0.47</sub></td><td>74.83<sup>+0.42</sup><sub>-0.38</sub></td><td>125.55<sup>+22.13</sup><sub>-79.22</sub></td><td>deg</td></tr><tr><td>f<sub>1</sub></td><td>10.25<sup>+0.01</sup><sub>-0.01</sub></td><td>10.45<sup>+0.01</sup><sub>-0.01</sub></td><td>10.76<sup>+0.06</sup><sub>-0.01</sub></td><td>11.00<sup>+0.01</sup><sub>-0.01</sub></td><td>10.63<sup>+0.01</sup><sub>-0.01</sub></td><td>9.45<sup>+0.08</sup><sub>-0.03</sub></td><td>log<sub>10</sub>(Jy, sr<sup>-1</sup>)</td></tr><tr><td>r<sub>1</sub></td><td>153.63<sup>+1.69</sup><sub>-1.88</sub></td><td>68.64<sup>+0.17</sup><sub>-0.03</sub></td><td>—</td><td>—</td><td>—</td><td>—</td><td>mas</td></tr><tr><td>σ<sub>1</sub></td><td>54.35<sup>+1.69</sup><sub>-1.72</sub></td><td>0.001<sup>+0.05</sup><sub>-0.001</sub></td><td>36.46<sup>+0.43</sup><sub>-4.61</sub></td><td>45.05<sup>+0.08</sup><sub>-0.07</sub></td><td>57.69<sup>+0.11</sup><sub>-0.11</sub></td><td>53.09<sup>+0.45</sup><sub>-8.60</sub></td><td>mas</td></tr><tr><td>f<sub>2</sub></td><td>—</td><td>—</td><td>10.39<sup>+0.01</sup><sub>-0.22</sub></td><td>—</td><td>—</td><td>—</td><td>log<sub>10</sub>(Jy sr<sup>-1</sup>)</td></tr><tr><td>r<sub>2</sub></td><td>—</td><td>—</td><td>121.48<sup>+0.58</sup><sub>-5.99</sub></td><td>—</td><td>—</td><td>—</td><td>mas</td></tr><tr><td>σ<sub>2</sub></td><td>41.92<sup>+1.15</sup><sub>-0.97</sub></td><td>129.20<sup>+0.21</sup><sub>-0.08</sub></td><td>8.36<sup>+8.59</sup><sub>-0.12</sub></td><td>—</td><td>—</td><td>—</td><td>mas</td></tr><tr><td>f<sub>3</sub></td><td>—</td><td>—</td><td>10.08<sup>+0.01</sup><sub>-0.01</sub></td><td>—</td><td>—</td><td>—</td><td>log<sub>10</sub>(Jy, sr<sup>-1</sup>)</td></tr><tr><td>r<sub>3</sub></td><td>—</td><td>—</td><td>184.67<sup>+1.77</sup><sub>-3.85</sub></td><td>—</td><td>—</td><td>—</td><td>mas</td></tr><tr><td>σ<sub>3</sub></td><td>—</td><td>—</td><td>60.76<sup>+1.85</sup><sub>-0.81</sub></td><td>—</td><td>—</td><td>—</td><td>mas</td></tr><tr><td>F<sub>0.87mm</sub></td><td>36.08 ± 0.17</td><td>49.05 ± 0.15</td><td>37.52 ± 0.1</td><td>25.52 ± 0.10</td><td>16.50 ± 0.05</td><td>0.95 ± 0.19</td><td>mJy</td></tr></tbody></table>

**Notes.** “mas” stands for milliarcsecond. The resulting $F_{0.87mm}$ of each model is given in the last row (the measured $F_{0.87mm}$ from the data is in Table A.2).

<figure><img src="image_4.png" /></figure>

**Fig. 3.** Normalized azimuthally averaged radial profiles of the dust continuum, $^{12}$CO, and $^{13}$CO of each source. Each curve has a 68% error region. The orange and red solid curves correspond to the average intensity profile for the gas profiles of the deprojected images, while the blue solid line corresponds to the best $\chi^2$ solution from the visibility fit of the dust continuum. The dashed vertical line denotes the position and $1\sigma$ error of the $R_{90}$ radius for the dust (blue) and the gas (red).

and 2, our angular resolution of ~90 mas was just enough to resolve the dust continuum and gas emission of this source, but the limited sensitivity did not allow us to constrain the geometric parameters as in the other systems.

For the two smooth sources, CIDA 7 and J0420, we found strong cloud contamination, so that the western side of the $^{12}$CO emission in J0420 is completely absorbed. The rotational pattern

was, however, recovered in the $^{13}$CO (see Figs. 1 and B.9). On the other hand, in CIDA 7, we observed an extended asymmetric emission in the south region, with a velocity range of at least 0.8 km s$^{-1}$ (it was detected in two velocity channels). The contribution of this emission is not significant to the total flux of the gas emission, and it was not detected in $^{13}$CO. Although the S/N did not allow us to accurately recover geometric parameters
<table><caption>**Table 3.** Continuum and CO radial extension for each source.</caption><thead><tr><th></th><th>R<sub>%</sub></th><th>R<sub>dust</sub> [au]</th><th>R<sub>12CO</sub> [au]<br>(R<sub>gas</sub>)</th><th>R<sub>13CO</sub> [au]</th><th>R<sub>gas</sub> / R<sub>dust</sub></th></tr></thead><tbody><tr><td rowspan="2">CIDA 1</td><td>68</td><td>24.57<sup>+0.03</sup><sub>-0.04</sub></td><td>73.7 ± 15.0</td><td>47.4 ± 6.8</td><td>3.0 ± 0.6</td></tr><tr><td>90</td><td>29.22<sup>+0.07</sup><sub>-0.09</sub></td><td>108.1 ± 15.0</td><td>68.5 ± 6.8</td><td>3.7 ± 0.5</td></tr><tr><td rowspan="2">MHO 6</td><td>68</td><td>34.18<sup>+0.01</sup><sub>-0.05</sub></td><td>130.7 ± 7.0</td><td>78.5 ± 7.1</td><td>3.8 ± 0.2</td></tr><tr><td>90</td><td>46.36<sup>+0.01</sup><sub>-0.07</sub></td><td>191.9 ± 7.0</td><td>113.0 ± 7.1</td><td>4.1 ± 0.2</td></tr><tr><td rowspan="2">J0433</td><td>68</td><td>36.58<sup>+0.11</sup><sub>-0.01</sub></td><td>114.6 ± 11.9</td><td>92.2 ± 11.7</td><td>3.1 ± 0.3</td></tr><tr><td>90</td><td>46.22<sup>+0.08</sup><sub>-0.07</sub></td><td>174.0 ± 11.9</td><td>125.5 ± 11.7</td><td>3.8 ± 0.3</td></tr><tr><td rowspan="2">CIDA 7</td><td>68</td><td>9.26<sup>+0.02</sup><sub>-0.02</sub></td><td>55.1 ± 10.8</td><td>19.6 ± 10.9</td><td>6.0 ± 1.2</td></tr><tr><td>90</td><td>13.16<sup>+0.02</sup><sub>-0.02</sub></td><td>78.9 ± 10.8</td><td></td><td>6.0 ± 0.8</td></tr><tr><td rowspan="2">J0420</td><td>68</td><td>14.84<sup>+0.03</sup><sub>-0.03</sub></td><td>43.1 ± 9.9</td><td>23.3 ± 9.8</td><td>2.9 ± 0.7</td></tr><tr><td>90</td><td>21.09<sup>+0.04</sup><sub>-0.04</sub></td><td>57.7 ± 9.9</td><td>32.2 ± 9.8</td><td>2.7 ± 0.5</td></tr><tr><td rowspan="2">J0415</td><td>68</td><td>10.87<sup>+1.78</sup><sub>-0.15</sub></td><td>26.4 ± 12.1</td><td>—</td><td>2.4 ± 1.2</td></tr><tr><td>90</td><td>15.46<sup>+2.53</sup><sub>-0.22</sub></td><td>34.8 ± 12.1</td><td>—</td><td>2.3 ± 0.9</td></tr></tbody></table>

**Notes.** Uncertainties for continuum comes from the walkers distribution in each MCMC. Gas radii uncertainties were calculated from the synthesized beam radius of each image.

<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Dust mass compared to stellar mass for the Taurus disks. We note that $M_d$ was computed using distances inferred from Gaia DR2, assuming a 20K disk temperature and using the most updated value of $F_{0.87\text{mm}}$ from Andrews et al. (2013), Tripathi et al. (2017), and Ward-Duong et al. (2018). Sources where $F_{0.87\text{mm}}$ was extrapolated from a $F_{1.3\text{mm}}$ measurement are plotted in faded colors (see Table 2 in Andrews et al. 2013), while direct measurements from SMA or ALMA are plotted in darker colors. Green triangles denote upper limits, and the VLMS studied in this work are plotted as orange stars.</figcaption></figure>

from the gas emission, we can see that for both sources, the PA obtained from the dust continuum emission is in good agreement with the orientation of the major axis in the rotation pattern.

In the sources with detected substructures, both J0433 and CIDA 1 are affected by cloud contamination. For CIDA 1, the side most affected by extinction is the north side, while the same is true for the south side of J0433 (row 1 and 3 of Fig. 1). The least cloud-contaminated source is MHO 6 disk (also seen in Fig. A.1), which is the brightest and more extended disk in CO emission of our sample. Since the cloud contamination is weaker, the gas can be detected continuously in every channel

from both high velocity ends. Given our good spatial detection of both lines, we used the package eddy (Teague 2019) to analyze the Keplerian rotation, which we further describe and discuss in Sects. 3.4.3 and 4.3.

### 3.4.2. Radial gas profiles

The inclinations and position angles obtained from the continuum fitting were used to deproject the distances from the central star in the moment 0 images, and we calculated the azimuthally averaged radial profiles of the $^{12}$CO and $^{13}$CO emission. The vertical structure was neglected in this calculation. In CIDA 1, J0433, and J0420, the profiles were calculated from the side that is less affected by cloud contamination. In CIDA 7, we found that the contribution of the asymmetric emission in the south is negligible, but we nevertheless did not take it into account as we wanted to recover the disk emission profile (see masking in Fig. 1). The gas $R_{68}$ and $R_{90}$ radii that encloses 68 and 90% of the flux are shown in Table 3 and Fig. 3.

As we discuss in Sect. 4.4, the brightness in dust continuum and gas emission of J0415 are much lower than expected from previous SMA observations (Andrews et al. 2013). The low S/N of the detection prevented us from applying self-calibration, and so the sensitivity is two times worse than in the other VLMS disks. Therefore, we calculated the radial profile of the $^{12}$CO emission from an image generated by only considering the channels in the velocity range with line detection, and we did not apply clipping at $3\sigma$.

### 3.4.3. Keplerian rotation of CO emission

Although we recovered the rotational pattern for all the disks in our sample in both CO isotopologues (with the exception of the $^{13}$CO in J0415), the strong cloud contamination and low S/N of our images prevented us from reliably recovering the dynamical mass from our sources without proper modeling. The only system where the cloud contamination does not completely make
the central velocity channels extinct is MHO 6, where we have good spatial detection of both lines. As described in Sect. 2.2, we used CASA to generate the moment 1 image from the velocity maps; in addition, we generated velocity integrated images using the python package bettermoments (Teague & Foreman-Mackey 2018), with a quadratic and Gaussian method, and also by varying the root mean square (RMS) clipping limit. The Keplerian rotation of these images was modeled using the package eddy (Teague 2019), with different models considering the stellar mass ($M_\star$), the central velocity (VLSR), the flaring parameter ($\psi$), the position angle (PA), and the source center ($x_0, y_0$) as fixed or free parameters.

Depending on the free parameters used (e.g., fixing the PA or allowing the fitting of vertical structure) and also on the velocity integrated image used, the stellar mass recovered would vary in the range of 0.16 ~ 0.24 $M_\odot$, which is consistent with a stellar mass derived from evolutionary models. However, it is important to note that all the models and images would leave residuals, which span two times the channel velocity widths and are also strongly structured, as shown in Fig. B.2.

To have a referential value for the stellar host mass, we used CASA5.6 to get the position velocity diagrams (PV diagram from hereafter) along the major axis of each source, based on the position angle obtained from the continuum UV-modeling. The only exception is J0415, whose PV diagram was obtained along the east-west axis. The PV diagrams are shown in Figs. B.3 and B.4 for the sources with continuum substructures and smooth profiles, respectively.

# 4. Discussion

## 4.1. Detection of substructures

We only detected obvious substructures in the brightest disks of our sample (CIDA 1, MHO 6, and J0433), which are also the most radially extended disks in gas and dust emission. The existence of strong dust traps located at larger radial distances from the star is most likely the reason for this observational result, as the dust is allowed to stay for longer timescales in the outer disk, thereby increasing the emitting area. On the other hand, our limited angular resolution only allowed us to detect substructures in the most extended disks.

This direct detection of substructures in 50% of our sample does not represent the occurrence rate of substructures in all disks around VLMS, as our sample is biased towards the brightest disks and our spatial resolution is limited. Therefore, our observations only allowed us to directly detect deviations from a simple Gaussian profile in the disks where the extent of the emission is consistently larger than the synthesized beam size.

CIDA 7 and J0420 are a good example of the limitations that our datasets have when detecting substructures. Even though we are unable to confirm the existence of dust substructures in these systems, the UV-plot in Fig. 2 of CIDA 7 shows some structure which is not described by a Gaussian profile, while the residual image of J0420 suggests that there might be more substructured dust emission that is not described by a single centrally peaked Gaussian in these disks. We lack the resolution and sensitivity to confidently characterize them. Given that these two disks are not totally dust depleted, the expected efficient radial drift must have been counteracted by a dust trapping mechanism. The compactness of these disks suggest that the dust trap is located so close to the star that our resolution did not allow us to detect it. In theory, any pressure bumps cannot be smaller than the local pressure scale height, which implies that if they are located closer to their

star, their radial extent is smaller than our current resolution. CIDA 7 and J0420 are good candidates to be targeted by deep observations with ALMA in the most extended antenna configuration, allowing us to detect substructures of $\le 2$ au in size at the distance of these targets, which is six and ten times smaller than the dust $R_{68}$ of those disks, respectively. Future observations will test if even the very small disks around VLMS are able to generate dust traps, as it is observed in the massive and extended ones.

For J0415, our UV-coverage and sensitivity resolves the emission, and the centrally peaked Gaussian model does not leave any significant residual. Higher sensitivity is needed in order to discern deviations from the Gaussian profile (see Fig. 2).

## 4.2. Origin of dust continuum rings and cavities

All the detected substructures in our sample resemble ring-like emission. MHO 6 is the only disk displaying what could be a hint of non-axisymmetric dust emission in the residuals at our current sensitivity. Rings structures are the most common type of substructure in moderate and massive stars, as shown by surveys such as DSHARP (Andrews et al. 2018a; Huang et al. 2018) and the Taurus survey at 1.3 mm (Long et al. 2018), and now a similar trend is found for the most extended disks around VLMS.

### 4.2.1. Detected structures coincide with possible CO iceline location

The iceline of each volatile marks the location in the disk where that volatile transitions from being mostly gas-phase to being frozen out on dust grains. It is possible that this phenomena could induce ring-like substructures in the dust continuum emission by changing the dust opacity and grain collisional fragmentation and growth rates (e.g., Okuzumi et al. 2016).

To investigate if any of the iceline locations of the major volatiles coincide with the location of our gaps, we calculated the midplane temperature of our disks by following Kenyon & Hartmann (1987), where, for an irradiated flared disk, the temperature was parametrized as:

$$T(r) = T_{\star} \left(\frac{R_{\star}}{r}\right)^{1/2} \phi_{\text{fl}}^{1/4}, \qquad (5)$$

with r is the distance from the star, $T_\star$ and $R_\star$ are the star temperature and radius, respectively, and $\phi_{\text{fl}}$ is the flaring angle, which we assume to be equal to 0.05 (Dullemond & Dominik 2004). The stellar radius $R_\star$ was measured from the stellar luminosity $L_\star$ by assuming black body emission and spherical symmetry ($L_\star = 4\pi R_\star^2 \sigma_{\text{sb}} T_\star^4$, with $\sigma_{\text{sb}}$ the Stefan-Boltzmann constant). As a first approximation, we ignored the contribution of accretion in our calculations of the stellar luminosity (further discussed in Long et al. 2018). If the stellar radius is replaced in Eq. (5), we can obtain the distance $r_i$ at which the temperature $T_i$ is reached, given a star of fixed luminosity $L_\star$, following:

$$r_i = \frac{L_{\star}^{1/2}}{T_i^2} \cdot \left(\frac{\phi_{\text{fl}}}{4\pi\sigma_{\text{sb}}}\right)^{1/2}. \qquad (6)$$

We considered the two coldest icelines presented in Zhang et al. (2015), which are CO with a sublimation temperature of ~20 ~ 28 K (we took the lower limit 20 K from Öberg et al. 2011 for direct comparison with Long et al. 2018), and the $N_2$ iceline which goes from ~12 ~ 15K. The iceline location was calculated for each star using the parameters shown in Table 1, and the results are shown in Fig. 5.
<figure><img src="image_6.png" /><figcaption>Fig. 5. Left panels: temperature profile of the midplane is obtained from Eq. (5) and plotted with a black curve. The blue curve shows the best radial profile obtained from the continuum UV-modeling, which was normalized to the peak emission. Vertical red and light blue shaded regions show the possible locations of the CO (light red) and N₂ (light blue) icelines. Right panel: disk radius is shown versus the square root of stellar luminosity (√L*). The shaded regions show the iceline location for CO (light red) and N₂ (light blue), obtained from Eq. (6). The peak location of the modeled continuum emission rings are shown in squares, while the gap locations of J0433 are shown with dots. The lower line in each target marks half of the theoretical resolution given by the longest baseline for each dataset, as measured from the center, and the upper line marks the continuum R₉₀.</figcaption></figure>

In the right panel of Fig. 5, the positions of the detected substructures are displayed as a function of the root of stellar luminosity. Most of the substructures detected lie within the region where the CO iceline is identified. However, the role that the CO iceline plays in the morphology of our disk is not conclusive from our datasets. Even if we neglect the uncertainties introduced by the stellar parameters and the disk temperature, it is not clear if the structures that were generated by an iceline would result in a ring peaked at the iceline position or in a gap (Pinilla et al. 2017c; van der Marel et al. 2018).

Pinilla et al. (2017c) find that icelines do not strongly change the gas density around the icelines, nor can they carve a dust cavity as the ones we detected in MHO 6 and CIDA 1. Given that the peak of the rings detected in these disks are located very close to the position of the CO iceline, we cannot rule out the possibility that the iceline played a role in triggering the mechanism that is carving the cavity. In J0433, we find both peaks and gaps at the region (or close) where the CO iceline could be located. Given the radial compactness of the disks in our sample, the limited angular resolution, and the wide radial range covered by the CO iceline location, it is also possible that the substructures coincide with the iceline just by chance, in the same way that other surveys of substructures have not found a strong correlation between the position of substructures and iceline location (e.g., Long et al. 2018). Additional deeper observations of the CO isotopologues could allow a better constrained modeling of the temperature of radial and vertical profile of these disks, thus providing better evidence for the role the iceline plays on the substructures detected.

### 4.2.2. CIDA 1

Our modeling of the cavity in CIDA 1 gives similar results to those found in Pinilla et al. (2018a), thus we do not further explore new possibilities for its cavity origin. It is important to note that our finding of the inner side of the ring being wider than the outer side is opposite to what we expected from the physical motivation of using this model, which accounts for the timescale of grain growth from micrometer to millimeter particles in dust traps (Pinilla et al. 2017b, 2018a). One possibility is that there is unresolved emission inside the main ring, and the inner side of the Gaussian is blending with it. Higher angular resolution observations are needed to understand the true nature of the dust distribution of this system.

### 4.2.3. J0433

If the gaps are assumed to be generated by a planet-disk interaction, we can use the width of the gap to estimate the mass of this gap-carving planet, under the assumption that this single planet is located at the position of the gap, and the ring is peaked at the local pressure maxima of the gas. In this scenario, if the physical parameters of the disk are kept constant, a more massive planet would create a wider gap (e.g., Fung et al. 2014; Kanagawa et al. 2015; Rosotti et al. 2016).

For a crude estimate for the mass of the planet in the first J0433 gap, we followed a similar procedure as Long et al. (2018). In this approach, we assumed that the distance between the gap and ring scales with the Hill radius of the planet:

$$R_{\text{Hill}} = r_p (M_p / (3 M_\star))^{1/3}, \quad (7)$$

where $r_p$ is the location of the planet (the location of the gap coincides with the location of the planet, as the planet is carving the gap), $M_p$ is the mass of the planet, and $M_\star$ is the mass of the host star. If we consider the distance between the gap and ring to be 5 $R_{\text{Hill}}$ (conservative upper limit for the gap width carved by a planet in Dodson-Robinson & Salyk 2011), and considering that our best model gives $r_{\text{ring}} - r_{\text{gap}} = 5$ au for the first ring, then the approximate mass of the planet would be $M_p \sim 0.1 M_{\text{Jup}}$.

This calculation has a large uncertainty depending on the type of simulation. For instance, Dodson-Robinson & Salyk (2011) estimated a maximum of 4 $R_{\text{Hill}}$ between the planet location and the ring position, while Pinilla et al. (2012) estimated 7 $R_{\text{Hill}}$ for planets as massive as Jupiter. Taking these two limits (4 $R_{\text{Hill}}$ and 7 $R_{\text{Hill}}$), the estimated planet mass becomes $M_p = 0.11_{-0.07}^{+0.10} M_{\text{Jup}}$. In this calculation, we did not account for the physical conditions of the disk, such as turbulence or temperature, nor did we consider the minimum mass to open a gap in the disk. Therefore, it just gives us a crude estimate for the mass order of magnitude of a single planet carving the gap. A future analysis should consider higher angular resolution observations to better constrain the gap-ring morphology as well as dedicated hydro-simulations to estimate the planet candidate mass.

### 4.2.4. MHO 6

In MHO 6, the modeling indicates the existence of a central cavity that is slightly asymmetric. Several processes in protoplanetary disks could create this type of an inner cavity, such
as photoevaporation due to stellar irradiation (e.g., Alexander & Armitage 2007; Owen et al. 2012; Owen & Kollmeier 2019), companions (e.g., Price et al. 2018), a planet-disk interaction (e.g., Rice et al. 2006; Zhu et al. 2011), and dead zones (e.g., Flock et al. 2015). In the following, we discuss each one of these scenarios for the MHO 6 cavity.

**Photoevaporation.** The 10 au dust cavity and the low accretion rate ($5 \times 10^{-11} M_{\odot} \text{yr}^{-1}$ estimated from a UV excess in Herczeg & Hillenbrand 2008) are in agreement with the predicted evolution for a very low mass star ($0.1 M_{\odot}$) having its material stripped away by a photoevaporative flow (Owen et al. 2012). However, we do not see clear evidence of a gas depleted cavity in the $^{12}\text{CO}$ emission, although the $^{13}\text{CO}$ emission suggests a decrease in the gas density in the inner region. A follow-up with dedicated VLMS photoevaporation models and higher angular resolution of this target in different molecular lines would be able to characterize the impact of this mechanism in the cavity we detect. From the current observations, this mechanism is still an option for the observed cavity.

**Companion.** Binary companions produce cavities in circumbinary disks, with sizes in the range of 3–5 times the binary semi-major axis, depending on the mass ratio, the eccentricity, and the disk viscosity (Artymowicz & Lubow 1994; Miranda et al. 2017; Ragusa et al. 2017). They also produce quasi-periodic variations in the accretion rate (Muñoz & Lai 2016). If MHO 6 was a circumbinary disk, the location of this companion could not be farther out than ~3 au (~20 mas) from the primary star. However, previous independent observations of MHO 6 have not found any evidence of multiplicity in this system (Briceño et al. 1998; Kraus & Hillenbrand 2007; Herczeg & Hillenbrand 2014, the latest reference identified companions in its survey as late-type as M0.0). A spectroscopic follow-up with radial velocities, given the high inclination of the system, could be useful to constrain the upper mass limit of a companion in the close inner region.

**Embedded planets.** The planet population around very low mass stars spans a wide range of masses, where even giant planets have been detected (e.g., Morales et al. 2019, a $0.46 M_{\text{Jup}}$ around a $0.1 M_{\odot}$ VLMS), meaning that protoplanetary disks around VLMS probably have the potential to create such objects. If we make a simple calculation following the opening gap criterium introduced by Crida et al. (2006), the approximate mass needed for a single planet to open a gap and explain the cavity observed in MHO 6 is 1 or $0.6 M_{\text{Saturn}}$ (0.18 and $0.3 M_{\text{Jup}}$), located at ~7 au for $\alpha = 10^{-3}$ and $\alpha = 10^{-4}$, respectively. If we assume a gas to dust ratio of 100, we estimate a disk mass of ~12.6 $M_{\text{Saturn}}$ or ≈3.8 $M_{\text{Jup}}$ given the dust mass obtained in Sect. 3.3. This means that the cavity opening planet is <10% of the current estimated mass of the disk. This implies that this type of potential single planet could have formed within the disk of MHO 6, although we do not exclude the possibility of multiple planets being responsible of this central cavity.

Another characteristic to take into consideration is the location of the peak brightness in the azimuthally averaged radial profile of $^{13}\text{CO}$ and dust continuum emission. We find hints of the cavity in the $^{13}\text{CO}$ line emission with its peak at 10 au when measured from the image. The same peak appears at 14 au when measured from the dust continuum image. Although these comparisons are strongly biased by the beam shape, it supports the idea that a planet may be responsible for the cavity formation since models predict such segregation between the outer edge of

the gap in gas versus dust (de Juan Ovelar et al. 2013; Facchini et al. 2018). Future higher angular resolution observations could test this hypothesis by better characterizing the CO brightness profile.

**Dead zone.** A dead zone is a low-ionized region at the disk midplane, where the dense environment blocks the high energy radiation, suppressing the magneto-rotation instability and, therefore, the angular momentum transport. It has been shown that the presence of dead zones can open gaps and cavities by forming gas pressure bumps at the outer edge of the dead zones, where dust trapping is efficient (e.g., Regály et al. 2012; Flock et al. 2015). Pinilla et al. (2016) predict that cavities that formed by dead zones alone would have millimeter- and micrometer-sized particles concentrated at the peak of the gas density. As a result, the radial location of the ring in millimeter wavelengths and scatter light would be the same. If we neglect temperature effects, our images suggest that the peak of $^{13}\text{CO}$ is closer to the star than the millimeter peak. Therefore, we expect for the micrometer-sized particles to peak closer to the star as well. In the dead zone scenario, Pinilla et al. (2016) also predict a strong gas depletion in the outer parts of the disk, further out from the dust trap, which is not currently seen in our observations of $^{13}\text{CO}$. If a magnetohydrodynamical wind is included in the models together with a dead zone, the observational diagnostics are very similar to planets. One step forward to try to disentangle between these models is to image this disk in scattered light and search for potential planets in the cavity. However, these disks are too faint to image given the limitations of current telescopes in the optical and near-infrared.

To summarize, the formation of the cavity observed in MHO 6 could be explained by one or a combination of the mechanisms we have discussed above. Several observational efforts can be done to disentangle these possibilities, such as a better estimate of the star accretion rate, a deeper search for planets or companions in its cavity in the optical and infrared wavelengths as well as imaging the disk in this regime, and very deep and higher angular resolution observations of molecular lines.

### 4.3. MHO 6 kinematics

Given that MHO 6 is the only disk with the gas emission detected across the whole velocity range, it was a good candidate for the dynamical mass measurement of its star. After several attempts to model MHO 6 rotation with different images, masks, and a combination of free parameters, we find our range of recovered stellar masses (between ~0.16 and $0.24 M_{\odot}$) to be consistent with the values from pre-main sequence models of MHO 6 (0.09–0.20, $M_{\odot}$ Kraus & Hillenbrand 2009; Herczeg & Hillenbrand 2014; Ward-Duong et al. 2018), but all the kinematic models would leave strong structured residuals with an amplitude spanning two times the channel velocity width.

Although our observations have enough spatial resolution to resolve the vertical structure in some of our disks (as seen in Fig. 1 as a cone-like emission in the $^{12}\text{CO}$ moment 0 of MHO 6 and J0433), the low S/N did not allow us to differentiate between the emission coming from the back side of the disk from that of the front side. The integrated velocity map contains the emission of both sides as if they were the same, and thus parameters such as $M_{\star}$ or $\psi$ from $Z(r, \psi) \propto r^{\psi}$ could not be recovered reliably.

This issue is not an exclusive problem in VLMS disks, but it applies to all the gas measurements with a high angular resolution and poor S/N. Future approaches to accurately recover
$M_*$ should consider more robust methods, such as UV-modeling (similar to the DiskJockey code from Czekala et al. 2015).

MHO 6 is a good candidate to further study the kinematics, substructures, and the physics of planet formation in VLMS. Its brightness allows high S/N observations at high angular resolution with non-prohibitive integration times. Combining the datasets presented in this paper, plus an observation of ≥5 h of time on source with ALMA using long baseline configurations should have enough resolution and sensitivity to precisely characterize the kinematics of the CO isotopologues, as previously done with disks around T-Tauri and Herbig stars (e.g., Pinte et al. 2019; Teague et al. 2019).

## 4.4. Comparison of dust $R_{68}$ and $L_{mm}$ with previous studies

All our sources have observed or estimated 0.87 mm fluxes, and a subset also had their dust continuum $R_{68}$ measured with previous UV-modeling. In the following, we compare and discuss our measurements with previous results, and we include the VLMS in the study of a size-luminosity relation.

**Flux.** J0415 was observed using the Submillimeter Array (SMA) at a wavelength of 1.3 mm, with a reported flux density of 12.6 ± 1.4 mJy in Andrews et al. (2013), where they extrapolated this measurement to 0.87 mm using $F_v \propto v^\alpha$, with $\alpha = 2.4 \pm 0.5$, thus obtaining $F_{0.87\text{mm}} = 32.9 \pm 15.2$ mJy. Although for all our other five sources, the extrapolations from the SMA flux measurements are consistent within the error range, the flux at 0.87 mm received from J0415 is approximately 35 times dimmer than expected, as observed independently by two different ALMA projects. A possible explanation could be a sudden change in millimeter flux due to flares, as it has been observed in low mass stars (e.g., MacGregor et al. 2018). For now, the only ALMA observations of this disk are the ones presented in this work, which are both in Band 7. Future observations at 1.3 mm should solve this discrepancy with the SMA results.

**Dust $R_{68}$.** The radius enclosing 68% of the dust continuum emission ($R_{68}$) was previously estimated for MHO 6 and J0433, using UV-modeling on SMA 340 GHz observations. In Tripathi et al. (2017) and Andrews et al. (2018b), they estimated a value of $R_{68} = 36.9_{-5.7}^{+8.5}$ au for MHO 6 (0.26''$^{+0.06}_{-0.04}$), and $R_{68} = 58.9_{-8.7}^{+6.9}$ au for J0433 (0.34''$^{+0.04}_{-0.05}$), which are consistent with the values measured by this work. Both are slightly overestimated (1.6 times for J0433) compared to the values of this work due to the considerably lower angular resolution of SMA observations compared to ALMA (0.86'' × 0.80'' for MHO 6 and 0.61'' × 0.52'' for J0433).

**$R_{68} - L_{mm}$ relation.** We combined our VLMS measurements with the Taurus sample observed with the SMA at 340GHz from Tripathi et al. (2017) in order to compare this all with the recent analysis by Hendler et al. (2020), where they obtained a relation of $R_{68} \propto L_{mm}^{0.53}$ for the Taurus star-forming region. For comparison purposes with this study, we followed the same approach using the Bayesian linear regression described in Kelly (2007), which was implemented in the python package `linmix` (publicly available in github, see Meyers 2015), to fit a linear relation between $R_{68}$ and $L_{mm}$ following:

$$ \log_{10}(R_{68}) = \alpha + \beta \log_{10}(L_{mm}), \qquad (8) $$

where $\alpha$ and $\beta$ are the regression coefficients. Our best fit was calculated by using the median value of the last 200 000 steps after convergence. We find that including our VLMS sample does not statistically change the previous result. However,

<figure><img src="image_7.png" /><figcaption>**Fig. 6.** Relation between $\log(R_{68})$ and $\log(F_{0.87\text{mm}})$. The points in blue are from Tripathi et al. (2017), with distances scaled to 140 pc using distances inferred from Gaia DR2, similar to Andrews et al. (2018b). The solid line in blue is the linear regression found by Hendler et al. (2020) (abbreviated as H+20), fitting the data from Tripathi et al. (2017). The VLMS from this work are shown in red. The best linear regression fits and their 68% confidence intervals are plotted in red when J0415 is excluded from the fit, and in green when all points are considered.</figcaption></figure>

the inclusion of J0415 changes the steepness of this relation in about $1\sigma$. If we consider J0415 as part of the fitting data, then the relation recovered is $\alpha = 2.09 \pm 0.09$, $\beta = 0.47 \pm 0.08$. If we exclude J0415 from the fitting, the result is $\alpha = 2.19 \pm 0.10$, $\beta = 0.59 \pm 0.10$. Both relations are close to the $1\sigma$ limit of each other, and they also overlap with the previous calculation from Hendler et al. (2020), as shown in Fig. 6. If we check their relations’ intrinsic scatter from the linear regression, we obtain $\sigma_{\text{scatter}} = 0.231 \pm 0.041$ and $\sigma_{\text{scatter}} = 0.217 \pm 0.040$, respectively. Therefore, including J0415 does not significantly increase the scatter.

The difference in the steepness of the recovered relations is not statistically significant from what was previously found. However, it is not completely clear if the same single power law relation between the millimeter luminosity and the size of the disks holds along the whole luminosity range. Given that J0415 is the only source with its size measured in the ~1 mJy brightness range, it is unknown if disks have some mechanisms to remain extended even when they are low in dust content (thus flattening the relation between size and luminosity in the low luminosity regime), or if a J0415-extended dust size is part of the relation scatter that is also observed in bright disks. To understand if J0415 is an outlier and to test if the power law behavior of the $R_{68} - L_{mm}$ relation flattens or holds at the low brightness regime, we need more deep observations at a high angular resolution of disks with $F_{0.87\text{mm}} < 10$ mJy. This could be achieved by observing each source from several tens of minutes to a few hours in ALMA Band 7.

## 4.5. J0415 dust radial extent

Although CIDA 7 has a dust content of at least ~27 times higher than J0415, it remains an open question as to how both can have a similar size (see Table 3). This result could be due to weaker dust traps in J0415 compared to those in CIDA 7, thus CIDA 7 can trap the dust more efficiently.

Alternatively, if the $M_{\text{gas}}$ of the disk is very low, the millimeter grains would have a high Stokes number and radial drift would become negligible. In Pinilla et al. (2017a), they
explore this scenario with a disk of 60 au in radius around a BD of 0.05 $M_\odot$. When $M_{\text{gas}} = 2 \times 10^{-2} M_{\text{Jup}}$, diffusion and drift still depleted the disk from millimeter particles; however, when $M_{\text{gas}} = 2 \times 10^{-3} M_{\text{Jup}}$, the millimeter grains were decoupled from the gas, and they could remain in the disk for longer timescales. For J0415, if we assume a dust to gas ratio of 1/100, we obtain $M_{\text{gas}} \approx 4 \times 10^{-2} M_{\text{Jup}}$, which does not seem to be low enough for a complete dust-gas decouple to occur. Additional observations are needed in order to characterize its gas radial density profile and to test this possibility.

Finally, it has been proposed that dust grains could be grow-
ing in a fractal manner, such that large aggregates would avoid
radial drift by maintaining a low Stokes number (Kataoka et al.
2013). This scenario can be distinguished from the compact mil-
limeter grains by measuring the opacity index $\beta$ between the
1 and 3 mm emission (Kataoka et al. 2014). Those observa-
tions, however, would require a very high sensitivity with enough
angular resolution to spatially resolve the radial profile of J0415.

4.6. $R_{\text{gas}}/R_{\text{dust}}$ ratio

We find that four of our six disks show a ratio between the $R_{90,\text{gas}}$ and $R_{90,\text{dust}}$ that is very close to or above 4.0, which is similar to the ratio measured in CX Tau (Facchini et al. 2019) (ratio of $3.9 \pm 0.5$ at $R_{90}$), as compared in Fig. 7. Our values, however, are conservative measurements of radii ratios for the sources with strong cloud contamination (all the sources in our sample, except for MHO 6), and so they should be considered as a lower limit. As it can be seen in the velocity channel maps in Figs. B.5, B.7, B.8, B.9, and B.10, the channels that are less affected by cloud contamination are the highest velocity channels, which are closer to the star. As soon as the emission gets radially extended in the channels closer to the central velocity, the cloud extinction becomes so high that we cannot see the disk emission anymore. Our integrated flux moment 0 images and the gas emission are then biased towards the compact emission located close to the star, underestimating the $R_{\text{gas}}$ measurement. This effect is particularly strong in J0420 and J0415, as they are the dimmest sources in this work, and J0420 also has the highest extinction (see Fig. A.1).

Our observations show that it is common for bright disks
around VLMS to have a gas radial extension of >3 times
the millimeter dust radial extension, which is consistent with
the efficient radial drift expected for millimeter-sized grains in
VLMS disks (Pinilla et al. 2013; Zhu et al. 2018). In fact, recent
thermochemical modeling including dust evolution by Trapman
et al. (2019, 2020) shows that ratios of $R_{\text{gas}}/R_{\text{dust}} > 4$ cannot
be solely explained by effects of optical depths; additionally, in
those cases, radial drift is required to explain the gas and dust
size difference. Despite the limited sensitivity and cloud con-
tamination, our disks are still very close and even above that
limit. Although similar ratios have been observed in sources with
a moderate stellar mass (~0.5 $M_\odot$, see Ansdell et al. 2018),
as shown in the upper panel of Fig. 7, the direct comparison
between our works is hindered by the data differences, such as in
the sensitivity and angular resolution, and also by the different
approaches used to obtain the gas and dust radii: our continuum
radial extension was calculated from the visibilities rather than
the images, and we did not apply Keplerian masking in the flux
integrated gas images. A more complete sample, which includes
more observations of disks around very low and moderate mass
stars, and also uniform data sensitivity and analysis are required
to confirm if there is a general trend where $R_{\text{gas}}/R_{\text{dust}}$ increases
towards lower stellar masses.

<figure><img src="image_8.png" /></figure>

**Fig. 7.** Upper panel: $R_{\text{gas}}/R_{\text{dust}}$ as a function of the stellar mass of the disks in the Lupus star-forming region (SFR) reported in Ansdell et al. (2018) (abbreviated as A+18), CX Tau (Facchini et al. 2019), and the targets reported in this paper. For CX Tau, we took the $R_{90}$ radii of the dust and gas. The VLMS of this work are shown in red for the sources with a substructure, and they are in yellow for the smooth sources, with the values from Table 3. Ratios from the Lupus SFR by Ansdell et al. (2018) were calculated from $R_{90}$ and a different color was used for disks identified as transition disks (TD). Lower panel: $R_{\text{gas}}/R_{\text{dust}}$ of the same targets, but compared to the dust mass of the disks. The dust mass was calculated from Eq. (4) under the assumption that $T = 20$ K for the whole disk midplane, with distances from Gaia DR2.

According to Trapman et al. (2019, 2020), they expect more
massive disks to have a higher $R_{\text{gas}}/R_{\text{dust}}$ ratio, driven by a larger
observed $R_{\text{gas}}$ due to the greater total CO content, and also
because the higher dust content would produce a more efficient
grain growth and inward radial drift. However, we do not observe
this trend in the lower panel of Fig. 7. Apparently, a decreasing radii ratio is obtained towards higher dust mass disks, which
could be the result of an efficient radial drift in disks around
VLMS, and the linear relation between $\log(M_\star)$–$\log M_{\text{dust}}$. This
supports the idea that smaller disks result from a fast radial drift
(Long et al. 2019) due to their inability to trap dust in the outer
regions, while more massive disks are more capable of creating
dust traps father away from the central star.

CIDA 7 stands out in our sample as having the most extreme $R_{\text{gas}}/R_{\text{dust}}$ ratio observed so far, with the gas being six times more extended than the dust, which is well beyond the ratio of four limit from Trapman et al. (2019). This confirms that radial drift is responsible for the compact size of this source. However, it is not completely dust depleted, so radial drift has been counteracted by another mechanism, or a combination of them. The
southern non-Keplerian emission detected in this system (in the velocity channels 4.4 and 4.8 km s⁻¹ from Figs. B.8 and 1) was masked when measuring the R<sub>90,gas</sub>, so the ratio of six is between the disk rotating gas and the dust size recovered from the model. We were unable to determine the origin of this extended emission in the south, as it is only detected in two different channels and it does not appear to be axisymmetric. A multiwavelength follow-up, with high sensitivity and angular resolution, might be required to understand the nature of this emission, as it could be explained by several different mechanisms, such as winds, outflows, interactions with external companions, an interaction with the surrounding cloud and envelope, among others.

The lowest gas-to-dust size ratio in our sample, measured in J0415, is likely due to the combined effects of lower than expected brightness, a low sensitivity due to our inability to apply self-calibration, and the extinction due to cloud contamination. Although we were unable to confidently recover the gas radius, our observations set a lower limit for its radial extension, and we confirm that in this very low disk mass regime, the gas emission is still more extended than dust emission. A more precise measurement of the gas radius requires a combination of deeper observations of lines less affected by the surrounding cloud and envelope as well as line modeling.

## 5. Conclusions

To understand the process of planet formation in VLMS and compare it with our current knowledge of planet formation in solar type stars, we observed and studied a sample of the brightest disks around VLMS in Taurus, at a 0.1" resolution and at a 0.87 mm wavelength. This sample is composed of CIDA 1, MHO 6, J0433, CIDA 7, J0420, and J0415 (2MASS names in Sect. 2). Here we summarize our main conclusions as follows.

- Detection rate of substructures: millimeter dust substructures were directly detected in only 50% of the targets in our sample. Our results suggest that the detection of substructures in disks around VLMS is limited by angular resolution and sensitivity, since the dust radial extent is very small and these disks are also very faint. Deep, high angular resolution observations over a non-brightness biased sample of VLMS should confirm the ubiquity of substructures in these disks.

- Substructured disks: substructures were detected in CIDA 1, MHO 6, and J0433; with the latest two being new detections. These three disks are the brightest and largest in our sample. They all have axisymmetric ring-like substructures, and only MHO 6 shows a weak asymmetry of amplitude less than 5% of the peak brightness. Both CIDA 1 and MHO 6 show central cavities in their emission. If we assume that a planet-disk interaction is the origin of the MHO 6 cavity, then a Saturn-mass planet (0.3 M<sub>Jup</sub>) is needed (as in the case of CIDA 1 Pinilla et al. 2018a). This planet should be located around 7 au. However, we cannot exclude other mechanisms that can also explain the origin of this cavity, such as multiple planets, a dead zone, a binary companion, or photoevaporation. Our UV-modeling of J0433 suggests that this disk could have two rings located at 21 and 32 au. However, we cannot confirm the separation between both at the 1 σ errorbar. We estimate that a planet of ~0.1 M<sub>Jup</sub> in mass could explain the first gap-ring. The substructures were detected within the region where the CO iceline could be located. Our current datasets lack the necessary S/N and resolution to properly characterize the vertical and radial temperature profile of the CO isotopologues, and so future deeper observations

will be needed to determine if the iceline played any role in triggering or maintaining the substructures observed.

- Smooth disks: the dust disks in CIDA 7, J0420, and J0415 are the less radially extended, less massive disks of our sample. With an angular resolution of 0.1", these disks are well described by a single Gaussian radial profile, which we used to measure their sizes. CIDA 7 is the most compact of them, with an R<sub>90</sub> = 13.16 au, which is similar to the 15.46 au from J0415. Yet, the dust mass estimate suggests that CIDA 7 is ~27 times more massive. In J0420, the residual continuum image shows some structured non-axisymmetric emission with 5σ peaks. However, this emission is very low in contrast to the smooth emission, which is over 300σ at its peak. Higher angular resolution observations are needed to describe the potential substructures in these disks.

- Size-luminosity relation: the disks in our sample follow a similar relation between L<sub>mm</sub> – R<sub>68</sub> as the one found for bright disks in the same star-forming region (see Hendler et al. 2020). However, our single measurement of a disk size in the low luminosity regime (J0415) needs to be complemented with deeper additional observations of other sources with a low stellar mass and low disk brightness. These measurements will help us understand the behavior of the size-luminosity relation across the whole range of disk sizes, enabling us to test if a single power law describes it.

- Evidence of efficient radial drift: when considering the dust and gas radii as the location where 90% of the emission is enclosed, four out of six disks in our sample show a ratio between R<sub>gas</sub>/R<sub>dust</sub> above 3.5. This is expected for disks where radial drift is depleting the dust, suggesting that radial drift is more efficient in VLMS than in moderate or high mass stars. However, the analysis of the sizes of disks in Lupus by Ansdell et al. (2018) also allowed them to find moderate stellar mass sources with similar disk size ratios as those observed in VLMS, which suggest that we need more observations to confirm this trend. Our comparison with Ansdell et al. (2018) does not account for the differences in data acquisition and R<sub>90</sub> calculation, which should be considered in future works, aiming for a uniform analysis of the extension of disks between different star-forming regions. The most extreme case of high R<sub>gas</sub>/R<sub>dust</sub> in our VLMS is observed in CIDA 7, with a value of six. This very high R<sub>gas</sub>/R<sub>dust</sub> ratio suggests that strong radial drift is at play (Trapman et al. 2019), raising the question about how this disk remains massive in dust.

Our observations do not exclude giant planet formation as an explanation for the substructures detected. Disks around VLMS follow similar trends as those that have been observed in disks around higher mass stars, based on our sample of bright disks. A future confirmation of a deviation from current correlations of physical parameters will require the recalculation of fluxes and sizes of the Taurus disks by using more sensitive and higher angular resolution observations from ALMA in a larger and more complete sample.

**Acknowledgements.** We thank Timmy Delage and Reema Joshi for their very helpful comments. We also thank Sean M. Andrews for the insightful discussions, which significantly improved this manuscript. N.K. and P.P. acknowledge support provided by the Alexander von Humboldt Foundation in the framework of the Sofja Kovalevskaja Award endowed by the Federal Ministry of Education and Research. This Project has received funding from the European Union's Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie grant agreement No 823823 (DUSTBUSTERS). This work was partly supported by the Deutsche Forschungs-Gemeinschaft (DFG, German Research Foundation) - Ref no. FOR 2634/1 TE 1024/1-1. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2012.1.00743.S,
ADS/JAO.ALMA#2015.1.00934.S, ADS/JAO.ALMA#2016.1.01511.S, and
ADS/JAO.ALMA#2018.1.00310.S. ALMA is a partnership of ESO (representing
its member states), NSF (USA) and NINS (Japan), together with NRC (Canada),
MOST and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation
with the Republic of Chile. The Joint ALMA Observatory is operated by ESO,
AUI/NRAO and NAOJ.

References

Alexander, R. D., & Armitage, P. J. 2007, MNRAS, 375, 500

Andrews, S. M. 2020, ARA&A, 58, 483

Andrews, S. M., Rosenfeld, K. A., Kraus, A. L., & Wilner, D. J. 2013, ApJ, 771, 129

Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018a, ApJ, 869, L41

Andrews, S. M., Terrell, M., Tripathi, A., et al. 2018b, ApJ, 865, 157

Ansdell, M., Williams, J. P., van der Marel, N., et al. 2016, ApJ, 828, 46

Ansdell, M., Williams, J. P., Trapman, L., et al. 2018, ApJ, 859, 21

Artymowicz, P., & Lubow, S. H. 1994, ApJ, 421, 651

Brauer, F., Dullemond, C. P., & Henning, T. 2008, A&A, 480, 859

Briceño, C., Hartmann, L., Stauffer, J., & Martín, E. 1998, AJ, 115, 2074

Crida, A., Morbidelli, A., & Masset, F. 2006, Icarus, 181, 587

Czekala, I., Andrews, S. M., Jensen, E. L. N., et al. 2015, ApJ, 806, 154

Daemgen, S., Natta, A., Scholz, A., et al. 2016, A&A, 594, A83

de Juan Ovelar, M., Min, M., Dominik, C., et al. 2013, A&A, 560, A111

Dodson-Robinson, S. E., & Salyk, C. 2011, ApJ, 738, 131

Dullemond, C. P., & Dominik, C. 2004, A&A, 417, 159

Dullemond, C. P., Birnstiel, T., Huang, J., et al. 2018, ApJ, 869, L46

Facchini, S., Pinilla, P., van Dishoeck, E. F., & de Juan Ovelar, M. 2018, A&A, 612, A104

Facchini, S., van Dishoeck, E. F., Manara, C. F., et al. 2019, A&A, 626, L2

Flock, M., Ruge, J. P., Dzyurkevich, N., et al. 2015, A&A, 574, A68

Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306

Fung, J., Shi, J.-M., & Chiang, E. 2014, ApJ, 782, 88

Gaia Collaboration (Brown, A. G. A., et al.) 2018, A&A, 616, A1

Hendler, N. P., Mulders, G. D., Pascucci, I., et al. 2017, ApJ, 841, 116

Hendler, N., Pascucci, I., Pinilla, P., et al. 2020, ApJ, 895, 126

Herczeg, G. J., & Hillenbrand, L. A. 2008, ApJ, 681, 594

Herczeg, G. J., & Hillenbrand, L. A. 2014, ApJ, 786, 97

Hildebrand, R. H. 1983, QJRAS, 24, 267

Huang, J., Andrews, S. M., Dullemond, C. P., et al. 2018, ApJ, 869, L42

Huang, J., Andrews, S. M., Dullemond, C. P., et al. 2020, ApJ, 891, 48

Jennings, J., Booth, R. A., Tazzari, M., Rosotti, G. P., & Clarke, C. J. 2020, MNRAS, 495, 3209

Kanagawa, K. D., Muto, T., Tanaka, H., et al. 2015, ApJ, 806, L15

Kataoka, A., Tanaka, H., Okuzumi, S., & Wada, K. 2013, A&A, 557, L4

Kataoka, A., Okuzumi, S., Tanaka, H., & Nomura, H. 2014, A&A, 568, A42

Kelly, B. C. 2007, ApJ, 665, 1489

Kenyon, S. J., & Hartmann, L. 1987, ApJ, 323, 714

Klein, R., Apai, D., Pascucci, I., Henning, T., & Waters, L. B. F. M. 2003, ApJ, 593, L57

Kraus, A. L., & Hillenbrand, L. A. 2007, ApJ, 662, 413

Kraus, A. L., & Hillenbrand, L. A. 2009, ApJ, 704, 531

Liu, B., Lambrechts, M., Johansen, A., Pascucci, I., & Henning, T. 2020, A&A, 638, A88

Long, F., Pinilla, P., Herczeg, G. J., et al. 2018, ApJ, 869, 17

Long, F., Herczeg, G. J., Harsono, D., et al. 2019, ApJ, 882, 49

Luhman, K. L. 2006, ApJ, 645, 676

MacGregor, M. A., Weinberger, A. J., Wilner, D. J., Kowalski, A. F., & Cranmer, S. R. 2018, ApJ, 855, L2

McMullin, J. P., Waters, B., Schiebel, D., Young, W., & Golap, K. 2007, ASP Conf. Ser., 376, 127

Mercer, A., & Stamatellos, D. 2020, A&A, 633, A116

Meyers, J. 2015, linmix, https://github.com/jmeyers314/linmix

Miranda, R., Muñoz, D. J., & Lai, D. 2017, MNRAS, 466, 1170

Morales, J. C., Mustill, A. J., Ribas, I., et al. 2019, Science, 365, 1441

Muñoz, D. J., & Lai, D. 2016, ApJ, 827, 43

Mulders, G. D., & Dominik, C. 2012, A&A, 539, A9

Mulders, G. D., Pascucci, I., & Apai, D. 2015, ApJ, 814, 130

Mužić, K., Scholz, A., Peña Ramírez, K., et al. 2019, ApJ, 881, 79

Öberg, K. I., Murray-Clay, R., & Bergin, E. A. 2011, ApJ, 743, L16

Okuzumi, S., Momose, M., Sirono, S.-i., Kobayashi, H., & Tanaka, H. 2016, ApJ, 821, 82

Owen, J. E., & Kollmeier, J. A. 2019, MNRAS, 487, 3702

Owen, J. E., Clarke, C. J., & Ercolano, B. 2012, MNRAS, 422, 1880

Pascucci, I., Testi, L., Herczeg, G. J., et al. 2016, ApJ, 831, 125

Pinilla, P., & Youdin, A. 2017, Astrophys. Space Sci. Lib., 445, 91

Pinilla, P., Benisty, M., & Birnstiel, T. 2012, A&A, 545, A81

Pinilla, P., Birnstiel, T., Benisty, M., et al. 2013, A&A, 554, A95

Pinilla, P., van der Marel, N., Pérez, L. M., et al. 2015, A&A, 584, A16

Pinilla, P., Flock, M., Ovelar, M. d. J., & Birnstiel, T. 2016, A&A, 596, A81

Pinilla, P., Quiroga-Nuñez, L. H., Benisty, M., et al. 2017a, ApJ, 846, 70

Pinilla, P., Pérez, L. M., Andrews, S., et al. 2017b, ApJ, 839, 99

Pinilla, P., Pohl, A., Stammler, S. M., & Birnstiel, T. 2017c, ApJ, 845, 68

Pinilla, P., Natta, A., Manara, C. F., et al. 2018a, A&A, 615, A95

Pinilla, P., Tazzari, M., Pascucci, I., et al. 2018b, ApJ, 859, 32

Pinilla, P., Pascucci, I., & Marino, S. 2020, A&A, 635, A105

Pinte, C., van der Plas, G., Ménard, F., et al. 2019, Nat. Astron., 41, 1109

Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62

Price, D. J., Cuello, N., Pinte, C., et al. 2018, MNRAS, 477, 1270

Ragusa, E., Dipierro, G., Lodato, G., Laibe, G., & Price, D. J. 2017, MNRAS, 464, 1449

Regály, Z., Juhász, A., Sándor, Z., & Dullemond, C. P. 2012, MNRAS, 419, 1701

Ricci, L., Testi, L., Natta, A., Scholz, A., & de Gregorio-Monsalvo, I. 2012, ApJ, 761, L20

Ricci, L., Testi, L., Natta, A., et al. 2014, ApJ, 791, 20

Ricci, L., Cazzoletti, P., Czekala, I., et al. 2017a, AJ, 154, 24

Ricci, L., Rome, H., Pinilla, P., et al. 2017b, ApJ, 846, 19

Rice, W. K. M., Armitage, P. J., Wood, K., & Lodato, G. 2006, MNRAS, 373, 1619

Rosotti, G. P., Juhasz, A., Booth, R. A., & Clarke, C. J. 2016, MNRAS, 459, 2790

Sanchis, E., Testi, L., Natta, A., et al. 2020, A&A, 633, A114

Schlafly, E. F., Green, G., Finkbeiner, D. P., et al. 2014, ApJ, 789, 15

Scholz, A., Jayawardhana, R., & Wood, K. 2006, ApJ, 645, 1498

Scholz, A., Jayawardhana, R., Wood, K., et al. 2007, ApJ, 660, 1517

Scholz, A., Muzic, K., Geers, V., et al. 2012, ApJ, 744, 6

Sinclair, C. A., Rosotti, G. P., Juhasz, A., & Clarke, C. J. 2020, MNRAS, 493, 3535

Tazzari, M., Beaujean, F., & Testi, L. 2018, MNRAS, 476, 4527

Teague, R. 2019, J. Open Source Softw., 4, 1220

Teague, R., & Foreman-Mackey, D. 2018, Res. Notes Am. Astron. Soc., 2, 173

Teague, R., Bae, J., & Bergin, E. A. 2019, Nature, 574, 378

Trapman, L., Facchini, S., Hogerheijde, M. R., van Dishoeck, E. F., & Bruderer, S. 2019, A&A, 629, A79

Trapman, L., Ansdell, M., Hogerheijde, M. R., et al. 2020, A&A, 638, A38

Tripathi, A., Andrews, S. M., Birnstiel, T., & Wilner, D. J. 2017, ApJ, 845, 44

van der Marel, N., Williams, J. P., & Bruderer, S. 2018, ApJ, 867, L14

van der Plas, G., Ménard, F., Ward-Duong, K., et al. 2016, ApJ, 819, 102

Ward-Duong, K., Patience, J., Bulger, J., et al. 2018, AJ, 155, 54

Zhang, K., Blake, G. A., & Bergin, E. A. 2015, ApJ, 806, L7

Zhu, Z., Nelson, R. P., Hartmann, L., Espaillat, C., & Calvet, N. 2011, ApJ, 729, 47

Zhu, Z., Andrews, S. M., & Isella, A. 2018, MNRAS, 479, 1850
## Appendix A: Summary of ALMA observations

A summary of the observation log is provided in Table A.1. The difference for configurations “Compact” and “Extended” denotes the stage at which we started its self-calibration based on the spatial extension of the antenna array, see Sect. 2.2. The properties of the dust continuum emission images and gas emission datacube are shown in Tables A.2 and A.3, respectively. In Fig. A.1 we show the reddening map from Schlafly et al. (2014) of the Taurus star-forming region, with the location of our sources. The reddening values for the sources are 0.61 for CIDA 1, 0.48 for MHO 6, 0.51 for J0433, 0.50 for CIDA 7, 1.38 for J0420, and 0.44 for J0415.

<figure><img src="image_9.png" /><figcaption>Fig. A.1. Spatial distribution of the Taurus sample studied in this work. The background is an extinction map compiled by Schlafly et al. (2014).</figcaption></figure>

<table><caption>**Table A.1.** Summary of ALMA observations.</caption><thead><tr><th>Source</th><th>Program ID</th><th>Obs. date</th><th>Exp. time (min)</th><th>N° Antennas</th><th>Baselines (m)</th><th>Configuration</th></tr></thead><tbody><tr><td rowspan="2">CIDA 1</td><td>2015.1.00934.S</td><td>2016-08-12</td><td>47.68</td><td>38</td><td>15-1462</td><td>Compact</td></tr><tr><td>2016.1.01511.S</td><td>2017-07-06</td><td>4.23</td><td>42</td><td>17-2647</td><td>Compact</td></tr><tr><td rowspan="3">MHO 6</td><td>2012.1.00743.S</td><td>2013-11-19</td><td>4.66</td><td>28</td><td>17-1284</td><td>Compact</td></tr><tr><td></td><td>2014-07-27</td><td>2.71</td><td>33</td><td>24-820</td><td>Compact</td></tr><tr><td>2018.1.00310.S</td><td>2019-08-21</td><td>47.15</td><td>45</td><td>41-3189</td><td>Extended</td></tr><tr><td rowspan="7">J0433</td><td>2012.1.00743.S</td><td>2013-11-17</td><td>2.80</td><td>29</td><td>17-1284</td><td>Compact</td></tr><tr><td></td><td>2014-07-27</td><td>2.27</td><td>33</td><td>24-820</td><td>Compact</td></tr><tr><td></td><td>2018.1.00310.S</td><td>2018-11-20</td><td>12.98</td><td>47</td><td>15-1398</td><td>Compact</td></tr><tr><td></td><td>2018-11-24</td><td>12.98</td><td>46</td><td>15-1261</td><td>Compact</td></tr><tr><td></td><td>2019-08-13</td><td>26.61</td><td>43</td><td>41-3144</td><td>Extended</td></tr><tr><td></td><td>2019-08-13</td><td>10.45</td><td>43</td><td>41-3144</td><td>Extended</td></tr><tr><td></td><td>2018.1.00310.S</td><td>2018-11-13</td><td>12.48</td><td>47</td><td>15-1398</td><td>Compact</td></tr><tr><td rowspan="2">CIDA 7</td><td></td><td>2019-08-24</td><td>26.04</td><td>47</td><td>41-3396</td><td>Compact</td></tr><tr><td></td><td>2019-08-24</td><td>26.04</td><td>47</td><td>41-3638</td><td>Extended</td></tr><tr><td rowspan="8">J0420</td><td>2012.1.00743.S</td><td>2013-11-19</td><td>4.66</td><td>26</td><td>17-1284</td><td>Compact</td></tr><tr><td></td><td>2014-07-27</td><td>2.27</td><td>33</td><td>24-820</td><td>Compact</td></tr><tr><td></td><td>2018-10-30</td><td>12.98</td><td>48</td><td>15-1398</td><td>Compact</td></tr><tr><td></td><td>2018-11-13</td><td>12.98</td><td>47</td><td>15-1398</td><td>Compact</td></tr><tr><td></td><td>2019-08-12</td><td>27.22</td><td>47</td><td>41-3638</td><td>Extended</td></tr><tr><td></td><td>2019-08-13</td><td>27.22</td><td>43</td><td>41-3144</td><td>Extended</td></tr><tr><td></td><td>2019-09-18</td><td>27.25</td><td>41</td><td>15-3638</td><td>Extended</td></tr><tr><td></td><td>2016.1.01511.S</td><td>2015-09-20</td><td>4.8</td><td>42</td><td>15-3189</td><td>Extended</td></tr><tr><td></td><td>2018.1.00310.S</td><td>2019-09-20</td><td>27.88</td><td>45</td><td>15-3189</td><td>Extended</td></tr></tbody></table>

<table><caption>**Table A.2.** Continuum imaging summary.</caption><thead><tr><th>Source</th><th>RA (ICRS)</th><th>Dec (ICRS)</th><th>Beam (mas × mas, deg)</th><th>Peak I<sub>v</sub> (mJy beam<sup>-1</sup>)</th><th>RMS noise (mJy beam<sup>-1</sup>)</th><th>F<sub>0.87mm</sub> (mJy)</th><th>M<sub>dust</sub> (M<sub>☉</sub>)</th></tr></thead><tbody><tr><td>CIDA 1</td><td>04:14:17.620</td><td>+28:06:09.289</td><td>150 × 111, -24</td><td>5.48</td><td>0.064</td><td>36.3 ± 0.1</td><td>8.40 ± 0.84</td></tr><tr><td>MHO 6</td><td>04:32:22.128</td><td>+18:27:42.286</td><td>104 × 72, 35</td><td>4.62</td><td>0.025</td><td>48.40 ± 0.15</td><td>12.25 ± 1.23</td></tr><tr><td>J0433</td><td>04:33:44.670</td><td>+26:15:00.080</td><td>128 × 86, -8</td><td>6.21</td><td>0.024</td><td>37.93 ± 0.12</td><td>14.32 ± 1.43</td></tr><tr><td>CIDA 7</td><td>04:42:21.022</td><td>+25:20:33.996</td><td>101 × 80, 10</td><td>11.53</td><td>0.029</td><td>25.96 ± 0.10</td><td>6.05 ± 0.61</td></tr><tr><td>J0420</td><td>04:20:25.581</td><td>+27:00:35.242</td><td>113 × 82, 22</td><td>6.45</td><td>0.021</td><td>17.35 ± 0.05</td><td>6.33 ± 0.63</td></tr><tr><td>J0415</td><td>04:15:58.016</td><td>+27:46:16.811</td><td>161 × 78, -36</td><td>0.44</td><td>0.052</td><td>0.96 ± 0.29</td><td>0.22 ± 0.02</td></tr></tbody></table>
<table><caption>**Table A.3.** $^{12}$CO and $^{13}$CO imaging summary.</caption><thead><tr><th>Source</th><th>Ch. width<br>(km s<sup>-1</sup>)</th><th>Detection<br>LSRK<br>(km s<sup>-1</sup>)</th><th>Beam<br>(mas × mas,<br>deg)</th><th>Peak I<sub>v</sub><br>(mJy beam<sup>-1</sup><br>km s<sup>-1</sup>)</th><th>RMS noise<br>(mJy beam<sup>-1</sup>)</th><th>Flux<br>(Jy)</th></tr></thead><tbody><tr><td rowspan="2">CIDA 1 <sup>12</sup>CO<br><sup>13</sup>CO</td><td>0.20</td><td>1.7–11.1</td><td>268 × 173, 1</td><td>91.2</td><td>4.0</td><td>0.725</td></tr><tr><td>1.00</td><td>2.0–9.1</td><td>118 × 83, 36</td><td>21.6</td><td>1.3</td><td>0.215</td></tr><tr><td rowspan="2">MHO 6 <sup>12</sup>CO<br><sup>13</sup>CO</td><td>0.25</td><td>0.65–10.1</td><td>115 × 81, 35</td><td>33.8</td><td>2.0</td><td>3.179</td></tr><tr><td>0.90</td><td>2.0–9.1</td><td>118 × 83, 36</td><td>21.6</td><td>1.3</td><td>1.061</td></tr><tr><td rowspan="2">J0433 <sup>12</sup>CO<br><sup>13</sup>CO</td><td>0.35</td><td>0.7–11.2</td><td>165 × 110, −7</td><td>72.7</td><td>2.2</td><td>1.241</td></tr><tr><td>1.0</td><td>1.0–11.0</td><td>159 × 110, −8</td><td>22.7</td><td>1.3</td><td>0.618</td></tr><tr><td rowspan="2">CIDA 7 <sup>12</sup>CO<br><sup>13</sup>CO</td><td>0.4</td><td>0.8–8.4</td><td>168 × 148, −1</td><td>61.9</td><td>2.6</td><td>1.065</td></tr><tr><td>1.0</td><td>3.0–8.0</td><td>170 × 149, 0</td><td>22.8</td><td>1.6</td><td>0.160</td></tr><tr><td rowspan="2">J0420 <sup>12</sup>CO<br><sup>13</sup>CO</td><td>0.45</td><td>7.75–10.9</td><td>135 × 98, 23</td><td>35.1</td><td>1.5</td><td>0.224</td></tr><tr><td>0.90</td><td>5.5–10.9</td><td>133 × 96, 22</td><td>8.2</td><td>1.0</td><td>0.038</td></tr><tr><td rowspan="2">J0415 <sup>12</sup>CO<br><sup>13</sup>CO</td><td>0.8</td><td>4.4–8.4</td><td>216 × 142, −38</td><td>32.1</td><td>3.8</td><td>0.020</td></tr><tr><td>1.0</td><td>—</td><td>222 × 141, −40</td><td>—</td><td>2.6</td><td>&lt;0.004</td></tr></tbody></table>

**Appendix B: Spectral profiles, eddy fitting, and channel maps**

<figure><img src="image_10.png" /><figcaption>**Fig. B.1.** Spectral profile of the $^{12}$CO and $^{13}$CO in our sample. The gray line denotes the 0 Jy level.</figcaption></figure>

The spectral profiles of $^{12}$CO and $^{13}$CO for each source, shown in Fig. B.1, were obtained with the CASA 5.6.2 software by placing an elliptical mask in the region where emission was detected. The errorbars represent the standard deviation of the noise, measured from the channels without emission.

The moment images shown in Fig. B.2 were also computed with CASA 5.6.2 by clipping the emission in each channel at the 3 σ level. For the Keplerian fitting with eddy (Teague 2019), we fixed the inclination to be one recovered from the continuum UV-modeling. The fitting was performed inside an elliptical masked region, with the same inclination of the source, extending up to 0.9" in the $^{12}$CO emission, and up to 0.6" in the $^{13}$CO. Each image was modeled separately.

The models and residuals shown in Fig. B.2 allowed the center of the source ($x_0, y_0$), the position angle (PA), the mass of the star ($M_*$), the central velocity (VLSR), and the vertical structure ($Z(r) = z_0 \cdot r^{\phi}$) to vary. The central velocity was consistent between the two lines, and PA was also consistent with continuum. For the stellar mass, $^{12}$CO would give a value of $M_* = 0.195 \pm 0.02 M_\odot$, and $M_* = 0.19 \pm 0.1 M_\odot$ from $^{13}$CO. Fixing some of these parameters to simplify the model, or using images obtained with bettermoments, would change $M_*$ between 0.16~0.24 $M_\odot$, but the structure of residuals would remain the same.
<figure><img src="image_11.png" /><figcaption>Fig. B.2. Keplerian fitting of $^{12}$CO and $^{13}$CO. The grid shows the best surface recovered, extending up to the distance of the mask for the fitting.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. B.3. Position velocity diagrams for the sources where substruc-
tures are detected in the continuum emission. Lines follow the Keplerian
velocities for central mass objects of 0.1, 0.2, and 0.3 M⊙ (solid, dashed,
and dotted, respectively).</figcaption></figure>

<caption>Fig. B.4. Same as Fig. B.3, but for the sources without a substructure.</caption>
<figure><img src="image_13.png" /><figcaption>**Fig. B.5.** Images of CIDA 1 ¹²CO channel maps are shown in the *top*, and ¹³CO channel maps are shown in the *bottom*. Each square is 3.0" in width and height, centered at the samed ¹³CO channel maps are shown in the bottom position as the dust continuum image. The contour level traces the 5σ emission in the continuum image. The scale bar in the *lower left panel* is 20 au in size, and the beam size is found in the lower left corner of the same panel. Central velocities in km s⁻¹ of each channel are given in the upper left corner of each panel.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>Fig. B.6. Same as Fig. B.5, but for MHO 6.</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>Fig. B.7. Same as Fig. B.5, but for J0433.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>Fig. B.8. Same as Fig. B.5, but for CIDA 7.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>Fig. B.9. Same as Fig. B.5, but for J0420.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Fig. B.10. Same as Fig. B.5, but for J0415.</figcaption></figure>