# Small and Large Dust Cavities in Disks around mid-M Stars in Taurus

Yangfan Shi (施杨帆) Kavli Institute for Astronomy and Astrophysics, Peking
University, Beijing 100871, China  Department of Astronomy, Peking University,
Beijing 100871, China  European Southern Observatory, Karl-Schwarzschild-Str.
2, D-85748 Garching bei München, Germany  Feng Long (龙凤) NASA Hubble
Fellowship Program Sagan Fellow  Lunar and Planetary Laboratory, University of
Arizona, Tucson, AZ 85721, USA  Gregory J. Herczeg (沈雷歌) Kavli Institute for
Astronomy and Astrophysics, Peking University, Beijing 100871, China
Department of Astronomy, Peking University, Beijing 100871, China  Daniel
Harsono Institute of Astronomy, Department of Physics, National Tsing Hua
University, Hsinchu, Taiwan  Yao Liu Purple Mountain Observatory & Key
Laboratory for Radio Astronomy, Chinese Academy of Sciences, Nanjing 210023,
China  Paola Pinilla Mullard Space Science Laboratory, University College
London, Holmbury St Mary, Dorking, Surrey RH5 6NT, UK  Enrico Ragusa
Dipartimento di Matematica, Università degli Studi di Milano, Via Saldini 50,
20133, Milano, Italy  Univ Lyon, Univ Lyon1, Ens de Lyon, CNRS, Centre de
Recherche Astrophysique de Lyon UMR5574, F-69230 Saint-Genis-Laval, France
Doug Johnstone NRC Herzberg Astronomy and Astrophysics, 5071 West Saanich Rd,
Victoria, BC, V9E 2E7, Canada  Department of Physics and Astronomy, University
of Victoria, Victoria, BC, V8P 5C2, Canada  Xue-Ning Bai Institute for
Advanced Study and Department of Astronomy, Tsinghua University, 100084,
Beijing, China  Ilaria Pascucci Lunar and Planetary Laboratory, the University
of Arizona, Tucson, AZ 85721, USA  Carlo F. Manara European Southern
Observatory, Karl-Schwarzschild-Str. 2, D-85748 Garching bei München, Germany
Gijs D. Mulders Facultad de Ingeniería y Ciencias, Universidad Adolfo Ibáñez,
Av. Diagonal las Torres 2640, Peñalolén, Santiago, Chile  Millennium Institute
for Astrophysics, Chile  Lucas A. Cieza Instituto de Estudios Astrofísicos,
Facultad de Ingeniería y Ciencias, Universidad Diego Portales, Av. Ejército
441, Santiago, Chile

###### Abstract

High-angular resolution imaging by ALMA has revealed the near-universality and
diversity of substructures in protoplanetary disks. However, disks around
M-type pre-main-sequence stars are still poorly sampled, despite the
prevalence of M-dwarfs in the galaxy. Here we present high-resolution
($\sim$50 mas, 8 au) ALMA Band 6 observations of six disks around mid-M stars
in Taurus. We detect dust continuum emission in all six disks, 12CO in five
disks, and 13CO line in two disks. The size ratios between gas and dust disks
range from 1.6 to 5.1. The ratio of about 5 for 2M0436 and 2M0450 indicates
efficient dust radial drift. Four disks show rings and cavities and two disks
are smooth. The cavity sizes occupy a wide range: 60 au for 2M0412, and
$\sim$10 au for 2M0434, 2M0436 and 2M0508. Detailed visibility modeling
indicates that small cavities of 1.7 and 5.7 au may hide in the two smooth
disks 2M0450 and CIDA 12. We perform radiative transfer fitting of the
infrared SEDs to constrain the cavity sizes, finding that micron-sized dust
grains may have smaller cavities than millimeter grains. Planet-disk
interactions are the preferred explanation to produce the large 60 au cavity,
while other physics could be responsible for the three $\sim$10 au cavities
under current observations and theories. Currently, disks around mid-to-late M
stars in Taurus show a higher detection frequency of cavities than earlier
type stars, although a more complete sample is needed to evaluate any
dependence of substructure on stellar mass.

Protoplanetary disks (1300); Planetary-disk interactions (2204); Planetary
system formation (1257)

††facilities: ALMA††software: CASA (CASA Team et al., 2022), GALARIO (Tazzari
et al., 2018), emcee (Foreman-Mackey et al., 2013), RADMC-3D (Dullemond et
al., 2012)

##  1 Introduction

The discovery to date of over 5000 exoplanets reveals that planetary systems
occupy a wide parameter space in architecture (Zhu & Dong 2021). In planetary
systems around low-mass stars (e.g., lower than 0.4 M⊙), small planets occur
more frequently than those around solar-mass stars (e.g., Mulders et al. 2015;
Hardegree-Ullman et al. 2019), while giant planets around those low-mass stars
are rare (but have a non-zero occurrence rate, e.g., Bryant et al., 2023).

The existence of a few giant planets around mid-to-late (later than M3) M
stars has challenged the core accretion planet formation theory (e.g., Morales
et al. 2019; Stefánsson et al. 2023). Planet population synthesis simulations
often fail to form any giant planet around these very low-mass stars, in both
the planetesimal accretion case (e.g., Miguel et al. 2020) and the pebble
accretion case (e.g., Liu et al. 2019; Mulders et al. 2021; Burn et al. 2021;
Schlecker et al. 2022). Special circumstances, such as a reduction in the
type-I migration velocities by a factor of 10, are needed to form planetary
cores more massive than 10 $M_{\oplus}$ (Burn et al. 2021; Schlecker et al.
2022), which then leads to the runaway gas accretion that enables giant planet
formation (Pollack et al. 1996). Besides core accretion, direct collapse in
gravitationally unstable disks may also form giant planets at early stages
(e.g., Boss 1997; Mercer & Stamatellos 2020; Longarini et al. 2023; Boss &
Kanodia 2023), when those class 0/I prostellar disks own higher disk mass than
their later stage counterparts (e.g., Tychoniec et al. 2020; Tobin et al.
2020). Possible evidence of gravitational instability has been found in some
recent ALMA observations, even though they targeted higher mass stars than
mid-to-late M stars focused in this work (e.g., Pérez et al., 2016; Paneque-
Carreño et al., 2021; Weber et al., 2023).

In the past decade, high angular-resolution observations with the Atacama
Large Millimeter/submillimeter Array (ALMA) have revealed that substructures
are common in protoplanetary disks. These substructures are mostly seen as
rings and gaps (e.g., ALMA Partnership et al. 2015; Andrews et al. 2018b; Long
et al. 2018; Cieza et al. 2021), asymmetric substructures like arcs and spiral
arms are also found but relatively rare (e.g., van der Marel et al. 2013;
Huang et al. 2018a; Dong et al. 2018). These substructures are expected to be
induced by various disk physical processes and planet-disk interactions (see
reviews by Andrews 2020; Bae et al. 2023). The prevalence of rings suggests
the presence of pressure bumps that halt the inward radial drift of the dust
(Pinilla et al. 2012).

Thus far these conclusions are obtained mostly from observations of solar-type
stars and Herbig stars. For disks around M stars, the dust drift should be
faster (Pinilla et al., 2013), so their disks may need stronger pressure bumps
to trap dust particles in order to explain their millimeter emission at ages
of a few Myr (e.g., van der Plas et al. 2016; Sanchis et al. 2020). However,
only little is known about substructures in disks around mid-to-late M stars.
Most structured disks around M-dwarf stars show rings and gaps (González-
Ruilova et al. 2020; Kurtovic et al. 2021; Pinilla et al. 2021; Cieza et al.
2021; van der Marel et al. 2022), where the gaps are mostly central cavities
or large gaps surrounding a compact inner disk (Pinilla 2022). Only one disk
around a mid M-dwarf shows a clear asymmetric ring, with properties similar to
asymmetries found around T Tauri and Herbig AeBe stars (Hashimoto et al.
2021). Detailed hydrodynamical simulations toward the disk around the mid-M-
dwarf CIDA 1 suggest that a planet with minimum mass of $\sim 1.4M_{\rm Jup}$
is needed to carve out the observed cavities present in 0.9 mm and 2.1 mm
continuum images (Curone et al. 2022), challenging the core accretion planet
formation theory around very low mass stars.

Following the pilot studies above, additional deep high-resolution
observations on disks around M stars (or very low mass stars) are needed to
evaluate the dust morphology and test disk physical processes in extreme
situations, following similar strategies applied to solar-mass stars (e.g.,
Andrews et al., 2018b; Cieza et al., 2019; Long et al., 2019). This paper
presents observations from a program designed to obtain high resolution
observations ($\sim 0.05^{{}^{\prime\prime}}$) of 16 mid-M star disks in the
Taurus star-forming region. Of the proposed targets, we have obtained
observations of six disks, including the double-ringed disk of 2MASS
J04124068+2438157, which was presented in Long et al. (2023). Here we present
all six observed disks (including 2M0412 from Long et al. 2023) to study
properties of substructures in M star disks and their possible origins.

The paper is organized as follows: in Section 2, we describe our target
selection, ALMA observations, data reduction and calibration. In Section 3, we
present our modeling method in the visibility plane and the modeling results,
as well as CO gas measurement. Section 4 discusses the global properties of
the six disks and possible origins for their detected substructures. In
Section 5 we summarize our main findings.

##  2 ALMA Program and Observations

###  2.1 Sample Selection and Properties

We analyze disks observed in the ALMA project 2019.1.00566.S (PI: G. Herczeg),
which aimed at studying the dependence of dust substructure properties on
stellar mass by targeting 16 disks around M3-M5 stars selected from Taurus.
The sample selection started from Luhman et al. (2017) with disk
identification based on excess emission in WISE W2 and W3 bands. We excluded
known binaries (e.g., Kraus et al. 2011; Daemgen et al. 2015) and sources with
$A_{V}>3$ mag or $J-$band brightness 1 mag fainter than the median to avoid
edge-on disks and embedded objects. Furthermore, targets with declination
$>$26 deg are excluded for visibility purposes. Like Long et al. (2019), the
selection ignored the millimeter brightness of disks. Six out of the 16
proposed targets were observed with the C43-9/10 configuration and Band 6
receivers between 2021-08-27 and 2021-09-27.

For these six sources, their stellar properties (listed in Table 2.1) were re-
evaluated based on optical spectra when available. For 2M0450 and 2M0508, the
spectral type and extinction are obtained from Luhman (2018), with uncertain
accretion properties. We obtained flux-calibrated low-resolution optical
spectra for 2M0412, 2M0434, and 2M0436 using UH88/SNIFS (3200–10000 Å, Lantz
et al. 2004) and for 2M0507 (CIDA12) using Palomar-Hale 5m/DBSP (3200–8700 Å
with a gap from 5600–6250 Å, Oke & Gunn 1982). The SNIFS data reduction is
described by Guo et al. (2018), while the DBSP data reduction is described by
Herczeg & Hillenbrand (2014). For these four targets, the spectral type,
extinction, and accretion rate are calculated from a simultaneous fit,
following approaches described in Herczeg & Hillenbrand (2008) and Herczeg &
Hillenbrand (2014) and using a plane-parallel slab model for accretion
developed by Valenti et al. (1993). The properties of 2M0412 were presented in
Long et al. (2023).

The spectral type for each star is converted to temperature using the
relationship of Herczeg & Hillenbrand (2014). The luminosity is then
calculated from the 2MASS $J$-band magnitude (Cutri et al., 2003) and assuming
zero $J$-band veiling, $J$-band bolometric corrections from Pecaut & Mamajek
(2013), and $J$-band extinction from the $A_{J}/A_{V}$ ratio developed by Wang
& Chen (2019). All properties are calculated from distances obtained from
inverting the parallax in Gaia Collaboration et al. (2023). The masses are
then estimated using the Somers et al. (2020) evolutionary tracks for pre-main
sequence stars with 51% spots coverage. These tracks lead to higher masses
than standard evolutionary tracks (see discussion for 2M0412 in Long et al.,
2023).

Table 1: Properties of Host Stars

2MASS $|$ Short Name | D | SpT | $A_{V}$ | $J$ | Lphot | $M_{*}$ | $R_{*}$ | Lacc | Macc | Refs
---|---|---|---|---|---|---|---|---|---|---
| [pc] |  | [mag] |  | [L⊙] | [M⊙] | [R⊙] | [L⊙] | [M⊙ yr-1] |
J04124068+2438157 $|$ 2M0412 | 148.7 | M4.3 | 0.84 | 11.151 | 0.126 | 0.30 | 1.21 | 0.00204 | 3.2e-10 | L23
J04343128+1722201 $|$ 2M0434 | 145.7 | M4.3 | 0.30 | 11.205 | 0.115 | 0.30 | 1.16 | 0.00168 | 2.6e-10 | –
J04360131+1726120 $|$ 2M0436 | 144.5 | M2.7 | 1.08 | 11.105 | 0.146 | 0.58 | 1.07 | 0.0062 | 4.5e-10 | –
J04504003+1619460 $|$ 2M0450 | 144.4 | M4.75 | 0.0 | 11.725 | 0.054 | 0.24 | 0.84 | – | – | L18
J05075496+2500156 $|$ CIDA 12 | 170.0 | M3.7 | 0.5 | 11.415 | 0.13 | 0.39 | 1.13 | 0.0028 | 3.2e-10 | HH14
J05080709+2427123 $|$ 2M0508 | 170.7 | M3.5 | 0.0 | 11.396 | 0.11 | 0.44 | 1.03 | – | – | L18

References. — L23: Long et al. (2023); L18: Luhman (2018); HH14: Herczeg &
Hillenbrand (2014).

###  2.2 ALMA Observations

For our ALMA observations, the receivers were configured into four spectral
windows, with two centered at 217.0 and 234.4 GHz for dust continuum with
bandwidths of 1.875 GHz, and the other two centered at 220.0 and 230.5 GHz
targeting 12CO, 13CO, and C18O $J=2-1$ lines. Line channel intervals are
spaced at 0.244 MHz ($\sim$ 0.3 km/s). Table 2 shows the detailed ALMA
observation log.

The raw visibilities were pipeline calibrated using the specified CASA
versions (6.1.1 for 2M0412, CIDA 12 & 2M0508; 6.2.1 for the rest) for each
object that can be found in the QA2 report (CASA Team et al., 2022). We
identified residual features of atmospheric absorption correction around
channel 500 in the 234.4 GHz spectral window for all six targets and flagged
corresponding channels 400-600. For dust continuum imaging, we flagged the
channels within 25 $\rm km\ s^{-1}$ of the rest frame of CO lines and then
averaged the dataset with a channel width of 125 MHz, which is the recommended
maximum bandwidth to avoid bandwidth smearing for ALMA Band 6.

These observations were intended to provide snapshots of these disks111Young
M-stars can happen to be active in millimeter wavelength in time scale of
minutes (e.g., MacGregor et al. 2018; Mairs et al. 2019), which would be an
issue for disk analysis especially for compact sources like 2M0450 and CIDA
12. After performing time series measurement of the fluxes for 2M0450 and CIDA
12, no significant flux variations are found during our observations.. Due to
the short on-source time ($\sim$15 min), the peak signal-to-noise (S/N) for
each source ranges from 5–20. We attempted one round of phase-only self-
calibration on sources 2M0412, 2M0434, and 2M0450. For 2M0436 and 2M0508, the
emission morphology was largely altered after self-calibration, so we did not
apply the solutions. CIDA 12 has too low S/N to improve from self-calibration.
The self-calibration solutions are not applied to CO line emission channels,
since the improvement is negligible.

The continuum images were generated using the ‘tclean’ task with multiscale
imaging, with the Briggs ‘robust’ weighting parameters being 0.5 for 2M0434,
2M0436, 2M0450 and 2M0508, 1 for 2M0412 and 2 for CIDA 12. The adopted
weighting parameters are compromises between resolution and sensitivity. The
uv coverage of observations produced elongated beams (aspect ratio $\sim$2).
For better visualization, we performed uv-tapering to get more circular beams.
Images with original beams and uv-tapered beams are shown in Figure 1. The
beam sizes are $0.03^{\prime\prime}\sim 0.07^{\prime\prime}$ and the RMS noise
is 30-40 mJy/beam, with detailed values summarized in Table 2.

For line images, after subtracting the continuum emission using ‘uvcontsub’
task we applied uv-tapering to enhance the SNR of the line images. The beams
are tapered to $0.2^{\prime\prime}$ for 2M0412 and 2M0434 and
$0.1^{\prime\prime}$ for the other four targets. We detected (or marginally
detected) 12CO for all six targets and 13CO for 2M0412 and 2M0434. The channel
maps in velocity and moment maps of detected 12CO and 13CO are shown in
Appendix A.

Table 2: Summary of ALMA observations

Source | Obs. Date | $N_{\rm ant}$ | Baselines | On-source Time | Mean PWV | Peak $I_{\nu}$ (taper) | RMS Noise (taper) | Beam (taper)
---|---|---|---|---|---|---|---|---
|  |  | [m] | [min] | [mm] | [mJy/beam] | [mJy/beam] | [mas$\times$mas, deg]
2M0412 | 2021-08-27 | 38 | 92.1-10803.3 | 16.80 | 0.3 | 0.19 (0.26) | 0.032 (0.039) | 74$\times$35, 37 (76$\times$73, 32)
2M0434 | 2021-09-27 | 45 | 70.1-14361.8 | 15.29 | 1.0 | 0.47 (0.69) | 0.039 (0.041) | 45$\times$24, -47 (43$\times$41, -41)
2M0436 | 2021-09-27 | 45 | 70.1-14361.8 | 15.39 | 1.0 | 0.31 (0.42) | 0.037 (0.039) | 44$\times$24, -46 (43$\times$42, 38)
2M0450 | 2021-09-28 | 47 | 70.1-14361.8 | 15.05 | 0.7 | 0.23 (0.27) | 0.029 (0.029) | 36$\times$26, -40 (36$\times$36, -35)
CIDA 12 | 2021-09-13 | 37 | 178.3-12594.5 | 16.80 | 0.6 | 0.19 (0.30) | 0.031 (0.036) | 72$\times$33, 4 (76$\times$73, -9)
2M0508 | 2021-09-13 | 37 | 178.3-12594.5 | 16.87 | 0.6 | 0.38 (0.63) | 0.035 (0.039) | 65$\times$25, 4 (63$\times$60, 36)

Note. — The numbers in parenthese in last three columns are corresponding
values for uv-tapered images.

Figure 1: The top panels show ALMA continuum images of disks around mid-M
stars in Taurus at 1.3mm. The box size for 2M0412 is $2.4^{\prime\prime}$ and
the rest are $1.2^{\prime\prime}$. The white scale bars represent 20 au. The
bottom panels are for the same targets with beams tapered to a more circular
shape for better visualization. The beam sizes and rms noise of each image are
detailed in Table 2.

##  3 Results

From our ALMA observations, four disks in the sample show rings and cavities
in their dust emission, while the other two show compact smooth emission. In
this section, we characterize the dust emission through fitting the directly
observable visibilities and present the measurement of dust and gas disk
properties.

###  3.1 Visibility Fitting of Dust Morphology

The peak fluxes of our disks range from 0.19 mJy/beam to 0.47 mJy/beam, which
corresponds to peak SNRs about 6 to 12 from the non-uv-tapered images in
Figure 1. Dust emission from these disks is faint, as expected for low-mass M
stars. To quantify the dust emission morphology without biases introduced by
the image reconstruction, we fit the brightness profiles in the visibility
domain. The model profiles are selected based on eye identification of disk
substructures. For 2M0450 and CIDA 12, the function to describe their
deprojected brightness is a centrally peaked axis-symmetric Gaussian profile
(GP):

$I_{\rm GP}(r)=10^{f}\exp{\left(-\frac{r^{2}}{2\sigma^{2}}\right)}.$ (1)

where $I_{\rm GP}(r)$ is the Gaussian brightness profile as a function of disk
radius $r$, with $10^{f}$ and $\sigma$ being the amplitude and Gaussian width.
For the other four disks showing cavities and rings, we modeled their
morphologies with a radially asymmetric Gaussian ring whose inner and outer
width can differ:

$I_{\rm AGR}=\begin{cases}10^{f}\exp{\left(-\frac{(r-r_{\rm
peak})^{2}}{2\sigma^{2}_{1}}\right)}\ (r\leq r_{\rm peak}).\\\
10^{f}\exp{\left(-\frac{(r-r_{\rm peak})^{2}}{2\sigma^{2}_{2}}\right)}\
(r>r_{\rm peak}).\end{cases}$ (2)

where $10^{f}$, $r_{\rm peak}$, and $\sigma_{i}$ are the amplitude, peak
location and ring width for the inner and outer side. A radially asymmetric
Gaussian ring has been previously used to describe the rings in transition
disks (e.g., Pinilla et al. 2017b, 2018; Huang et al. 2020; Kurtovic et al.
2021). This model fitting is motivated by the accumulation of dust particles
trapped in radial pressure bumps (Pinilla et al. 2017b). Dust particles in the
outer disk are expected to take a longer time to grow and drift radially
toward the pressure bump. Hence, the external width of the ring is expected to
be larger than the internal width (e.g., Figure 4 in Pinilla et al. 2015).

Using the adopted radial profiles, we first generate corresponding face-on
images. Each model image is projected with an inclination ($i$) and a position
angle (PA) and combined with a spatial offset ($\delta_{\rm RA},\delta_{\rm
Dec}$) which adds another four parameters. We generate each model image with a
pixel size of 3 mas, which is far smaller than our synthesized beam ($\sim$ 40
mas). We then use galario (Tazzari et al., 2018) to transform each model image
into complex visibilities sampled at the same $(u,v)$ points as those in the
observations and calculate the $\chi^{2}=\sum_{k}w_{k}|V_{\rm
obs}(u_{k},v_{k})-V_{\rm mod}(u_{k},v_{k})|^{2}$, where $w_{k}$ is the
observed visibility weight, the likelihood function is then calculated as
$L\propto\exp(-\chi^{2}/2)$. We assumed a uniform prior probability
distribution for the parameters above (for inclination the probability
$\propto\sin i$ to get a uniform disk orientation).

We sample the posterior probability distribution using a Markov Chain Monte
Carlo (MCMC) method (emcee, Foreman-Mackey et al. 2013) with 50 walkers and
10000 steps. After running MCMC, we get posterior samples after 10 times the
integrated autocorrelation time to ensure stationary posterior samples.
Results for each model parameter are shown in Table 3, with the reference
values adopted as the median values of posterior samples and their
uncertainties estimated from the 16th and 84th percentiles. The comparisons
between observations and models in the visibility plane, the image plane, and
deprojected azimuthally averaged radial profiles are shown in Figure 2. The
brightness profiles of each model are shown in Figure 3. All best-fit model
images show reasonable matches with observations, with no residuals above
$5\sigma$.

Figure 2: Best-fit visibility fitting results versus observations. Columns
from left to right: (1) the real part of deprojected and binned visibilities
from observations and best-fit models, (2) observational continuum images, (3)
model images convolved with the same beam as observations, (4) residual images
with white contours at $-3\sigma$ and green contours at $3\sigma$, overlaid
orange contours are from continuum images at $3\sigma$ level (5) azimuthally
averaged radial profiles for observations and convolved model images, light
gray shaded regions represent the standard deviation at each radial bin
divided by the square root of the number of beams in the radial bin.

###  3.2 Rings and cavities

####  3.2.1 Disks with substructures in images

As mentioned in Section 3.1, rings are modeled as radially asymmetric Gaussian
rings. Our results show that the rings generally have wider outer width than
the inner, while 2M0436 shows the opposite. Here we summarize the peak
locations and Gaussian $\sigma$ widths from the best-fit models.

2M0412 has the largest rings and cavity in the sample, with the inner ring
peak at 61.7 au and the outer ring peak at 113.9 au. The inner ring has an
inner width of 0.9 au and outer width of 17.5 au, while the outer ring has 5.5
au inner side and 9.5 au outer side. For the ratios between the ring outer
width and inner width, best-fit model gives 19.8 for the inner ring and 1.7
for the outer ring. The outer ring is brighter than the inner ring with the
outer peak intensity 2.1 times as large as the inner peak. The brightness
ratio between the first ring and the exterior gap is around 6.8.

For 2M0434, we find a ring peak at 9.8 au, followed by radially extended
emission beyond the ring. The ring has 1.5 au inner width and 8.2 au outer
width. A Gaussian profile in Section 3.1 fits the outer broad emission well.
For 2M0508, the best model finds a ring peak at 6.8 au, where its inner width
is 0.8 au and outer width is 7.1 au. In 2M0436, the best-fit ring peaks at
12.2 au, with an outer width being 1.0 au and a larger inner width being 4.7
au showing different ring width behavior.

Some potential asymmetries appear in the images. For example, the northern
part of 2M0434 seems to be brighter than its southern part as shown by the
residual map and a bright spot appears at the west-north of 2M0436, which are
both at $3\sigma$ levels. The west part of 2M0508 is brighter though the
significance is less than $3\sigma$ level. Deeper observations are necessary
to check whether above potential asymmetries exist or not. In this work, we
focus on only axis-symmetric substructures.

####  3.2.2 Potential cavities around 2M0450 and CIDA 12

Gaussian profiles (Equation 1) can only describe continuous disks with a
monotonically decreasing radial profile. The Nuker profile, which was first
introduced to describe the brightness profiles of elliptical galaxies (Lauer
et al. 1995), can also reproduce ring-like emission with a central depression.
The Nuker profile is formulated as:

$I_{\nu}(r)=10^{f}\left(\frac{r}{r_{t}}\right)^{-\gamma}\left[1+\left(\frac{r}{r_{t}}\right)^{\alpha}\right]^{(\gamma-\beta)/\alpha}.$
(3)

where $r_{t}$ is the transition radius, $\alpha$ is the transition index,
$\beta$ is the outer disk index, and $\gamma$ is the inner disk index. The
Nuker profile behaves as $r^{-\gamma}$ at $r\ll r_{t}$ and $r^{-\beta}$ at
$r\gg r_{t}$. The transition index $\alpha$ determines how smooth/sharp the
transition is between these two asymptotic behaviors (see Figure 2 in Tripathi
et al. 2017). With $\beta>0,\ \gamma<0$ the Nuker profile can describe a ring-
like morphology.

Appendix C presents the detailed modeling with Nuker profiles on 2M0450 and
CIDA 12. Interestingly, the best-fit Nuker profiles for the two smooth-
appearing disks show depleted inner emission with emission peak at $\sim$ 1.7
au for 2M0450 and 5.7 au for CIDA 12 (see Figure 16).

However, the Nuker profile has more parameters than the Gaussian profile,
which results in a larger Bayesian Information Criterion222The Bayesian
Information Criterion (BIC) is a criterion for the preference of models. It is
defined as ${\rm BIC}=k\ln(n)-2\ln(\hat{L})$, where $k$ is the number of model
parameters, $n$ is the number of data points, and $\hat{L}$ is the maximized
value of the likelihood function. value. Hence, we still treat 2M0450 and CIDA
12 as non-structured disks in this work. Future higher angular resolution
observations are needed to check their potential small cavities.

###  3.3 Potential unresolved central emission in 2M0436

2M0436 shows a ring with a wider inner width, which contradicts with models of
dust trapping in pressure bump. One possibility is that an external object has
truncated the outer disk, which we expect to be substellar since no stellar
companions have been detected (Kraus et al. 2011; Kraus & Hillenbrand 2012).
Constraining the mass of the substellar object given current knowledge would
be highly uncertain and beyond the scope of this work. Other possibilities are
the presence of central compact emission blending with the ring emission and
the blending of multiple unresolved rings. Compact emission at the center is
also hinted at because the model intensity of 2M0436 slightly underpredicts
the emission of the very inner disk (right column in Figure 2), though this
difference is not statistically significant. An example is CIDA 1, where a
model considering only one ring results in a ring with wider inner width
(Kurtovic et al. 2021), while higher-angular resolution observations reveal a
compact inner disk and the inclusion of an inner disk recovers a ring with
wider outer width (Pinilla et al. 2021).

In Appendix B, we show the results of modeling 2M0436 with a ring and a
central point. The flux of the added central point emission is
$0.04^{+0.03}_{-0.03}$ mJy (corner plot statistics in Figure 15). The new
radial profile increases the inner emission and matches the observations
better than a ring-only model, while the ring still displays a wider inner
width after including a point emission. However, the result of modeling with
central point emission could be limited by current observation, deeper
sensitivity and higher angular resolution are necessary for future
observations.

###  3.4 Millimeter Flux, Dust Disk Mass and Size

The millimeter fluxes, masses and sizes of the dust disk are not explicit
parameters in the models and are inferred from fitting results. Since our
models are axis-symmetric, the total fluxes are obtained through
$F_{\nu}=f_{\nu}(\infty)$, where $f_{\nu}(r)$ is the flux within radius $r$:

$f_{\nu}(r)=2\pi\int_{0}^{r}I_{\nu}(r^{\prime})r^{\prime}{\rm d}r^{\prime}.$
(4)

The dust disk sizes are defined to be a radius that encloses 68% and 90% of
total fluxes, consistent with choices from previous studies (see, for example,
the review by Miotello et al., 2023). The continuum fluxes, sizes, and their
uncertainties are then computed from the last 5000 steps of the chains, as
shown in Table 3. Most adopted models reproduce lower fluxes than those seen
at short baselines (left column in Figure 2), likely due to the fact that
either large scale emission is resolved out by high angular resolution beams
or faint emission buried by the noise, so that we lack extra components in our
modeling which could account for the flux difference between models and
observations. Future observations with deeper sensitivity and shorter
baselines are needed to fix this issue.

For CIDA 12, Nuker profiles reproduce millimeter flux better than that from
Gaussian profiles. The visibility near zero-spacing baselines indicates its
flux above 1 mJy and it is consistent with $1.16^{+0.09}_{-0.09}$ mJy reported
in Akeson et al. (2019), specifically the Gaussian profile gives
$0.60^{+0.06}_{-0.06}$ mJy, while the Nuker profile recovers a flux of
$1.20^{+0.98}_{-0.51}$ mJy.

Under the assumption that the dust disk is optically thin at millimeter
wavelengths, the dust disk mass $M_{\rm dust}$ is estimated following
Hildebrand (1983), as

$M_{\rm dust}\approx\frac{D^{2}F_{\nu}}{\kappa_{\nu}B_{\nu}(T_{\rm dust})}.$
(5)

where $D$ is distance from Gaia DR3 (Gaia Collaboration et al., 2023),
$\kappa_{\nu}$ is dust opacity for which we adopt a power law of
$\kappa_{\nu}=2.3\ {\rm cm}^{2}{\rm g}^{-1}(\nu/230\ {\rm GHz})^{0.4}$
(Andrews et al., 2013), and $B_{\nu}$ is the Planck function where the $T_{\rm
dust}$ is assumed at $20\ {\rm K}$ (e.g., Ansdell et al. 2016).

We assume 10% uncertainty in the source flux calibration (Diaz Trigo et al.
2019, see also Francis et al. 2020) to calculate the dust mass uncertainty.
Table 3 summarizes the estimated dust masses. 2M0434 has the most massive disk
in our sample, with a dust mass $19.31^{+1.93}_{-1.93}M_{\oplus}$, while CIDA
12 has the lowest dust mass with only $0.51^{+0.05}_{-0.05}M_{\oplus}$. We
note that these dust masses serve as lower limits, since part of the disk
could be optically thick in reality (see Section 4.1 for a brief calculation).
The assumed dust temperature would also affect the estimated dust mass: for
example, the dust mass would be higher/lower by $\sim 40\%/50\%$ if the actual
averaged dust temperature is 15 K/30 K.

Table 3: Dust Disk Model Parameters from uv Modeling

Source Name | 2M0412 | 2M0434 | 2M0436 | 2M0450 | CIDA 12 | 2M0508 | Unit
---|---|---|---|---|---|---|---
Model | 2AGR | GP&AGR | AGR | GP | GP | AGR |
$f_{1}$ |  $8.63^{+0.03}_{-0.03}$ | $9.83^{+0.01}_{-0.01}$ |  $9.98^{+0.03}_{-0.03}$ | $10.45^{+0.05}_{-0.05}$ |  $9.54^{+0.09}_{-0.08}$ | $9.93^{+0.03}_{-0.03}$ |  $\log_{10}$ Jy/sr  |  |  |
$r_{\rm peak,1}$ |  $415.10^{+8.71}_{-6.41}$ | -.- |  $84.68^{+6.05}_{-7.35}$ | -.- | -.- | $39.80^{+6.38}_{-4.54}$ | mas  |  |  |
$\sigma_{11}$ |  $5.93^{+6.51}_{-4.21}$ | $268.38^{+3.04}_{-3.02}$ |  $32.90^{+5.84}_{-6.23}$ | $16.89^{+1.44}_{-1.28}$ |  $54.06^{+11.11}_{-11.15}$ | $4.69^{+5.94}_{-3.47}$ | mas  |  |  |
$\sigma_{12}$ |  $117.71^{+13.71}_{-11.45}$ | -.- |  $6.88^{+5.34}_{-4.63}$ | -.- | -.- | $41.91^{+3.61}_{-4.01}$ | mas  |  |  |
$f_{2}$ |  $8.96^{+0.02}_{-0.02}$ | $9.95^{+0.03}_{-0.03}$ | -.- | -.- | -.- | -.- |  $\log_{10}$ Jy/sr  |  |  |
$r_{\rm peak,2}$ |  $766.27^{+7.44}_{-14.94}$ | $67.41^{+8.91}_{-7.60}$ | -.- | -.- | -.- | -.- | mas  |  |  |
$\sigma_{21}$ |  $37.05^{+7.31}_{-12.02}$ | $10.15^{+7.94}_{-6.71}$ | -.- | -.- | -.- | -.- | mas  |  |  |
$\sigma_{22}$ |  $63.74^{+8.73}_{-5.68}$ | $56.22^{+4.81}_{-5.52}$ | -.- | -.- | -.- | -.- | mas  |  |  |
Inc |  $15.90^{+0.70}_{-0.82}$ | $68.54^{+0.23}_{-0.23}$ |  $53.46^{+1.41}_{-1.49}$ | $43.26^{+7.84}_{-10.37}$ |  $65.53^{+7.43}_{-12.65}$ | $54.34^{+1.67}_{-1.78}$ | deg  |  |  |
PA |  $124.14^{+2.81}_{-2.06}$ | $96.67^{+0.24}_{-0.24}$ |  $46.96^{+2.03}_{-2.05}$ | $63.01^{+12.43}_{-12.87}$ |  $169.94^{+5.22}_{-6.59}$ | $98.90^{+2.10}_{-2.10}$ | deg  |  |  |
$\Delta_{\rm RA}$ |  $2.77^{+1.86}_{-2.06}$ | $10.89^{+1.03}_{-1.05}$ |  $7.86^{+1.12}_{-1.15}$ | $-3.13^{+0.96}_{-0.97}$ |  $9.05^{+3.21}_{-3.30}$ | $12.36^{+1.38}_{-1.40}$ | mas  |  |  |
$\Delta_{\rm Dec}$ |  $-0.69^{+2.11}_{-2.12}$ | $10.37^{+0.50}_{-0.49}$ |  $5.76^{+1.11}_{-1.11}$ | $-8.98^{+0.91}_{-0.92}$ |  $-2.15^{+6.02}_{-6.18}$ | $5.77^{+1.10}_{-1.12}$ | mas  |  |  |
$F_{\rm 1.3mm}$ |  $17.58^{+0.33}_{-0.32}$ | $30.73^{+0.26}_{-0.24}$ |  $2.71^{+0.05}_{-0.05}$ | $0.85^{+0.03}_{-0.03}$ |  $0.60^{+0.06}_{-0.06}$ | $3.02^{+0.09}_{-0.09}$ | mJy  |  |  |
$M_{\rm dust}$ |  $11.50^{+1.15}_{-1.15}$ | $19.31^{+1.93}_{-1.93}$ |  $1.68^{+0.17}_{-0.17}$ | $0.53^{+0.05}_{-0.05}$ |  $0.51^{+0.05}_{-0.05}$ | $2.61^{+0.26}_{-0.26}$ |  $M_{\oplus}$ |  |  |
$R_{68\%}$ |  $118.0^{+0.5}_{-0.3}$ | $55.1^{+0.5}_{-0.5}$ |  $11.8^{+0.2}_{-0.2}$ | $3.6^{+0.3}_{-0.3}$ |  $13.8^{+2.8}_{-2.9}$ | $15.6^{+0.5}_{-0.5}$ | au  |  |  |
$R_{90\%}$ |  $126.0^{+0.9}_{-0.6}$ | $81.2^{+0.8}_{-0.8}$ |  $13.2^{+0.4}_{-0.3}$ | $5.1^{+0.5}_{-0.4}$ |  $19.6^{+4.0}_{-4.1}$ | $20.4^{+0.8}_{-0.8}$ | au  |  |  |

Note. — Model row: AGR for asymmetric Gaussian ring, GP for Gaussian profile.

###  3.5 CO emission

Gas disk observations and comparisons to dust disks are helpful to understand
disk evolution as a whole. To compute the integrated fluxes of CO lines,
extraction regions need to be determined. We first generate the cumulative
flux distributions and determine the radii where CO emission no longer
increases. Then we draw circular regions as large as these radii, and project
the regions with the same inclination and position angle as those of dust
disks. CO fluxes are then retrieved from these elliptical regions. The
uncertainties of the integrated flux are measured in non-emission channels
with the same velocity range as those used to generate moment-0 images, as the
standard deviation of integrated fluxes estimated from 1000 randomly
distributed elliptical regions with the same shape as the above extraction
regions. The fluxes and uncertainties of 12CO, 13CO and C18O $J=2-1$ lines are
listed in Table 4, as well as the $3\sigma$ upper limits for the non-
detections.

We detect (SNRs of integrated fluxes above $3\sigma$) 12CO in five disks and
13CO in 2M0412 & 2M0434. 2M0450 has a marginal detection of 12CO with SNR of
$2.7\sigma$. The gas disk sizes are measured the same way as dust disk sizes,
adopting the same inclination and position angles from the fits to the dust
continuum emission. 2M0450 has the smallest gas disk about 26 au while 2M0412
has the largest gas disk about 199 au. When compared to dust disk sizes, the
ratio ranges from 1.6 (2M0412) to 5.1 (2M0450). Since the largest angular
scale recovered by our observations is only $\sim 0.6^{\prime\prime}$, the gas
emission at large radii is likely resolved out; cloud absorption could also
result in underestimation of disk sizes (e.g., Long et al. 2022); by 10%-30%
with a typical cloud emission line width of 1-2 $\rm km\ s^{-1}$, we treat the
gas radii $R_{\rm gas,90\%}$ as lower limits. The measured gas disk sizes and
size ratios compared to dust disks are summarized in Table 4.

Figure 3 shows the azimuthally averaged radial profiles of gas emission and
the best-fit continuum model intensity profiles. CO emission around 2M0412
suffers from severe cloud absorption, so we show the profile extracted from
the moment-8 map333Moment-8 map represents the maximum value along the
spectrum. instead. Among the six disks, 2M0412 shows a drop in emission inside
the dust cavity, with peak emission located at $\sim$ 36 au for 12CO and 54 au
for 13CO. The cloud absorption may affect these values. A tentative drop of
12CO is seen around 2M0508, with the peak location overlapping with the dust
cavity size. No other clear evidence of gas emission depletion is found in the
other four disks, under current $\sim 0.1^{\prime\prime}$
($0.2^{\prime\prime}$ for 2M0434) resolution for gas images.

Figure 3: Normalized radial profiles of dust continuum, 12CO and 13CO (if detected). The thick orange curves are the best-fit model profiles of dust continuum (i.e. unconvolved); the thin light orange lines correspond to 1000 randomly selected model intensities from posterior sample (Section 3.1). The gray lines are radial profiles of dust continuum from observations. The blue and green curves correspond to 12CO and 13CO azimuthally averaged radial profiles retrieved from moment-0 images, except for 2M0412 the moment-8 images are used due to severe cloud contamination. The color-filled regions represent $1\sigma$ uncertainty for each radial profile. The sticks at the top of each panel represent the radius enclosing $90\%$ emission of each profile, orange for dust continuum, blue for 12CO and green for 13CO. Table 4: Measurement of gas emission from CO isotopes Molecule | Integrated Flux | Peak Intensity | Velocity Range | Chan. Width | Chan. RMS | Beam | $R_{90\%}$ | $R_{\rm gas}/R_{\rm dust}$
---|---|---|---|---|---|---|---|---
| [Jy $\rm km\ s^{-1}$] | [mJy beam$\rm{}^{-1}km\ s^{-1}$] | [$\rm km\ s^{-1}$] | [$\rm km\ s^{-1}$] | [$\rm mJy\ beam^{-1}$] | [mas$\times$mas, deg] | [au] |
2M0412
12CO | $1.08\pm 0.06$ | 71.4 | 4.0-8.0 | 0.4 | 5.5 | 215$\times$209, -3 |  $\geq$198.9 |  $\geq$1.6
13CO | $0.54\pm 0.06$ | 39.2 | 4.0-8.0 | 0.4 | 5.8 | 220$\times$212, -2 |  $\geq$187.2 |  $\geq$1.5
C18O | $<0.15$ | - | 4.0-8.0 | 0.4 | 4.2 | 222$\times$211, 5 | - | -
2M0434
12CO | $1.21\pm 0.11$ | 130.5 | 2.0-10.0 | 0.4 | 6.9 | 209$\times$175, -67 |  $\geq$153.9 |  $\geq$1.9
13CO | $0.61\pm 0.09$ | 74.2 | 2.0-10.0 | 1.0 | 4.4 | 170$\times$159, -57 |  $\geq$155.9 |  $\geq$1.9
C18O | $<0.21$ | - | 2.0-10.0 | 1.0 | 3.4 | 169$\times$159, -55 | - | -
2M0436
12CO | $0.38\pm 0.07$ | 56.1 | 1.4-9.4 | 0.8 | 3.5 | 112$\times$99, -51 |  $\geq$64.4 |  $\geq$4.9
13CO | $<0.24$ | - | 1.4-9.4 | 0.8 | 3.8 | 113$\times$100, -55 | - | -
C18O | $<0.17$ | - | 1.4-9.4 | 0.8 | 2.9 | 113$\times$101, -54 | - | -
2M0450
12CO | $0.10\pm 0.04$ | 30.3 | 2.0-12.0 | 1.0 | 2.6 | 110$\times$97, -77 |  $\geq$25.8 |  $\geq$5.1
13CO | $<0.11$ | - | 2.0-12.0 | 1.0 | 2.9 | 112$\times$98, -82 | - | -
C18O | $<0.10$ | - | 2.0-12.0 | 1.0 | 2.2 | 111$\times$99, -85 | - | -
CIDA 12
12CO | $0.24\pm 0.05$ | 57.6 | 2.0-12.0 | 1.0 | 3.0 | 101$\times$94, -37 |  $\geq$53.7 |  $\geq$2.7
13CO | $<0.15$ | - | 2.0-12.0 | 1.0 | 3.3 | 102$\times$96, -40 | - | -
C18O | $<0.12$ | - | 2.0-12.0 | 1.0 | 2.5 | 101$\times$97, -35 | - | -
2M0508
12CO | $0.19\pm 0.06$ | 39.6 | 2.0-12.0 | 1.0 | 2.9 | 100$\times$93, -39 |  $\geq$56.3 |  $\geq$2.8
13CO | $<0.20$ | - | 2.0-12.0 | 1.0 | 3.3 | 101$\times$96, -43 | - | -
C18O | $<0.19$ | - | 2.0-12.0 | 1.0 | 2.5 | 101$\times$97, -38 | - | -

##  4 Discussion

Among the six disks, the four brighter disks show cavity and ring
substructures and the two faintest disks appear as smooth emission. The
cavities can be large to 60 au while the most common sizes are about 10 au.
For the two smooth disks, visibility fitting suggests the possibility of them
harboring cavities down to 2 au. Comparing the gas disks with dust disks, two
disks show large size ratios of about 5. With these results, we organise our
discussions as follows: In Section 4.1, we discuss the disk global properties
under disk evolution context. Section 4.2 focuses on the detected
substructures at millimeter wavelength and SED hints for cavities for micron-
sized dust disk. At last, we discuss the mechanisms that are possibly
responsible for the observed rings and cavities in Section 4.3.

Figure 4: Scaling relations for Taurus disks. Left: Millimeter luminosity at
1.3 mm versus stellar mass. Open orange stars are structured disks and solid
orange stars are smooth disks in this work. Blue circles are for disks in Long
et al. (2019), with open ones for structured and filled ones for smooth disks.
Open green circles represent other structured disks in literature (see
references in text). The gray dots (detection) and triangles (upper limits)
are other Taurus disks (Manara et al. 2023). The light pink line is the
scaling relation fitted to disks with large cavities ($>$20 au) from Pinilla
et al. (2020). Central: Thumbnails of each target’s continuum image, the box
sizes are varied so that the targets are shown clearly. Right: Dust disk size
$R_{\rm dust,68\%}$ versus millimeter luminosity. The orange stars are the
same as left panel, blue circles are from Long et al. (2019) with open and
filled for structured and smooth disks respectively, and gray dots are from
Tripathi et al. (2017); Kurtovic et al. (2021). The scaling relation is from
Hendler et al. (2020). The gray dots and the scaling relation are extrapolated
from observations at ALMA Band 7 ($\sim$ 0.9mm).

###  4.1 Context with disk evolution

Putting disk properties together with the disk population can help constrain
global disk evolution (Manara et al., 2023). Previous studies have shown that
the millimeter luminosity of disks scales with stellar mass (e.g., Andrews et
al., 2013; Ansdell et al., 2016, 2017; Barenfeld et al., 2016; Pascucci et
al., 2016). The left panel of Figure 4 compares the six disks in this work
with other Taurus disks with detected substructures in the $L_{\rm mm}-M_{*}$
plane (Long et al. 2018; Hashimoto et al. 2021; Kurtovic et al. 2021;
Yamaguchi et al. 2021; Jennings et al. 2022; Zhang et al. 2023; review by Bae
et al. 2023), also with other Taurus disks from Manara et al. (2023). The
properties of the six disks follow the correlation for the whole Taurus
population (Andrews et al., 2013). The four structured disks are brighter than
the two smooth disks. It could be that the two smooth disks simply formed with
less dust, while it could also be the lack of substructures to trap the dust
and prevent dust from drifting toward the star. However, if the potential
small cavities shown in Section 3.2.2 are indeed present in 2M0450 and CIDA
12, either their trapping efficiency is lower or the formation time of dust
traps is later than the brighter disks.

If disks are optically thick, the low $L_{\rm mm}$ might not correspond to a
low dust mass. To have a simple estimate of the optical depth for the two
smooth disks, we estimate the disk midplane temperature following Huang et al.
(2018b):

$T_{\rm mid}(r)=\left(\frac{\phi_{\rm fl}L_{*}}{8\pi r^{2}\sigma_{\rm
SB}}\right)^{1/4}.$ (6)

where $\sigma_{\rm SB}$ is the Stefan-Boltzmann constant, $\phi_{\rm fl}$ is
the disk flaring angle. A conservative $\phi_{\rm fl}=0.02$ is used (Huang et
al. 2018b). Then the optical depth $\tau_{\nu}$ is calculated using:

$I_{\nu}(r)=B_{\nu}(T_{\rm mid}(r))(1-\exp(-\tau_{\nu}(r))).$ (7)

where $B_{\nu}$ is the Planck’s law for black-body radiation. The
conservatively estimated peak $\tau$ is $\sim 0.5$ for 2M0450 and $0.15$ for
CIDA 12 so that 2M0450 is partially optically thick and CIDA 12 is nearly
optically thin. We expect the above arguments regarding dust trapping to hold
for CIDA 12 at least.

The $L_{\rm mm}-M_{*}$ relation has been suggested to be flatter for disks
with inner cavities (Pinilla et al., 2018, 2020, pink line in left panel of
Figure 4). While those studies lack data for the very low mass stars and disks
with small cavities ($<$20 au). Furthermore, the flatter relation is not
observed in any individual star-forming region. The two fainter disks with
small inner cavities (2M0436 and 2M0508) deviate significantly from this
flatter relation. Conservative estimates of their optical depth yield maximum
values of 0.6 for 2M0436 and 0.5 for 2M0508, suggesting the two disks are
partially optically thick. Boulder formation could also lower the observable
millimeter flux (Pinilla et al. 2020). A combination of optical thickness and
boulder formation in these disks might explain their deviation.

In addition to stellar mass, the millimeter luminosity also scales with disk
size (Tripathi et al., 2017; Andrews et al., 2018a). The right panel of Figure
4 compares the six disks in this work with other Taurus disks (Tripathi et
al., 2017; Long et al., 2019; Kurtovic et al., 2021) in the $R_{\rm
eff}-L_{\rm mm}$ plane with the best-fit scaling relation $R_{\rm eff}\propto
L_{\rm mm}^{0.53\pm 0.12}$ from Hendler et al. (2020).444Fluxes of sources in
Kurtovic et al. (2021) and the scaling relation in Hendler et al. (2020) are
extrapolated from 0.9 mm to 1.3 mm using a spectral index of 2.2 (Andrews,
2020), the slopes of the scaling relation at these two wavelength are found to
be identical (Tazzari et al., 2021). Numerical simulations have shown that
disks fall along the observed relation $L_{\rm mm}\propto R^{2}_{\rm eff}$ if
disks are in the radial drift dominated regime, while for disks with strong
substructures they follow a relation of $L_{\rm mm}\propto R^{5/4}_{\rm eff}$
(Rosotti et al., 2019; Zormpas et al., 2022). The two relations intersect at
the top end of the size-luminosity relation (i.e. for bright and large disks,
Zormpas et al. 2022).

Most of our six disks follow the observed relation $L_{\rm mm}\propto
R^{2}_{\rm eff}$ except 2M0412, whose large size stands out at its luminosity
regime, suggesting that a pressure bump was built early at large disk radii of
2M0412 (Long et al., 2023). 2M0436, 2M0450, CIDA 12 and 2M0508 fall at the
faint end of $L_{\rm mm}\propto R^{2}_{\rm eff}$ relation hence are likely
drift dominated. For the structured 2M0436 and 2M0508, either their
substructures are too weak to retain the dust in their rings or their rings
formed late when the disks have already evolved to the faint end of the
$L_{\rm mm}\propto R^{2}_{\rm eff}$ relation.

The gas disk is universally found to be more extended than the dust disk
(e.g., Barenfeld et al. 2017; Ansdell et al. 2018; Long et al. 2022). In our
six disks, we find that the two largest dust disks (2M0412 and 2M0434) have a
$R_{\rm 90,gas}$/$R_{\rm 90,dust}$ of $\sim$ 1.5 and 1.9 respectively. CIDA 12
and 2M0508 have a ratio of $\sim$ 2.7. The two smallest dust disks (2M0436 and
2M0450) have the most extreme ratios around 5, which are also among the
largest values for disk population (Long et al. 2022). Early formation of
pressure bumps in the two largest dust disks 2M0412 and 2M0434 may explain
their lower size ratios. 2M0412 show a large ring at large radii of 114 au,
while 2M0434’s ring peaks at about 10 au which is far smaller than its size of
81 au, indicating its extended emission outside the ring holds unresolved
substructures. Trapman et al. (2019) suggests that disks with $R_{90,\rm
gas}/R_{90,\rm dust}>4$ can only be explained with dust evolution and radial
drift. This is likely what 2M0436 and 2M0450 have undergone due to their such
extreme size ratios, which is also consistent with more efficient dust radial
drift expected in lower mass stars (Pinilla et al., 2013). These size ratios
should be considered as lower limits since they are affected by cloud
absorption and the observations lack short-spacing data (see Figures in
Appendix A), hence the gas emission is probably not recovered as well as dust
emission.

###  4.2 Substructures in disks around mid-to-late M stars

Current studies of substructures mostly come from disks around solar-like
stars (e.g., Andrews et al. 2018b; Long et al. 2018; Cieza et al. 2021). van
der Marel & Mulders (2021) found that substructures are less common around M
stars, which was mostly based on studies with intermediate spatial resolution
($\sim$ 25 au in radius) and biased towards larger disks, while detecting
substructures around M stars needs higher resolution since they are smaller.
Currently only a few disks around M stars (especially mid-to-late type) have
been imaged at high spatial resolution. In Taurus, Kurtovic et al. (2021)
surveyed at a resolution of $\sim 0.1^{\prime\prime}$ (14 au) for a sample of
six disks around M4-M5 stars and found two disks with a cavity and one ringed
disk. The sample in Kurtovic et al. (2021) is biased to the brightest disks in
the corresponding stellar regime hence favorable for substructure detection.
Hashimoto et al. (2021) reported an asymmetric dust ring around a cavity
around the M4.5 star ZZ Tau IRS with resolution $\sim 0.2^{\prime\prime}$.
Other studies toward Lupus and Ophiuchus regions have also detected a few
structured disks around mid-to-late M stars with resolutions down to $\sim 4$
au (González-Ruilova et al., 2020; Cieza et al., 2021; van der Marel et al.,
2022).

Cavity+ring seem to be the most common substructure around low mass stars. In
our sample, substructures are detected in four (2M0412, 2M0434, 2M0436, and
2M0508) of six disks. The substructures we identify are all cavities
surrounded by rings, with varying numbers in individual disks. Combined with
structured disks around very low-mass stars collected in Pinilla (2022), the
cavity frequency seems higher than that around solar-like stars. Taking disks
in Taurus as an example (Hashimoto et al. 2021; Kurtovic et al. 2021; this
work), high spatial resolution imaging revealed 7 disks with cavities out of
13 disks around M3-M6 stars. In comparison to stars earlier than M3 in Taurus,
only 4 disks out of 32 disks analyzed by Long et al. (2019) show cavities
(5/32 considering the reanalysis in Zhang et al. 2023). However, the sample
around mid-to-late M stars used here might be more biased than that in Long et
al. (2019). Based on SED modeling (see below and Figure 5), two disks out of
the six disks in this work have cavities $>1$ au, which is higher than the
transition disk fraction $\sim 8\%$ within the full Class II disks (van der
Marel et al. 2016, similar fraction is obtained when limiting to M3-M6 stars).
Furthermore, disks in this work have been imaged with $\sim
0.05^{\prime\prime}$ compared to $0.1^{\prime\prime}$ in Long et al. (2019). A
complete sample around mid-to-late M stars and higher resolution imaging of
the inner regions of disks around solar-like stars are needed to evaluate any
difference in their cavity occurrence rate in millimeter wavelength.

Gap+ring pairs are only detected in the large double ring disk 2M0412. 2M0434
has a ring located at about 10 au while its dust disk size is about 80 au, the
very extended emission in between may hold shallow rings and gaps inside.
Rings around 2M0436 and 2M0508 have radial widths comparable to the beam size
hence are not resolved – at higher resolution they might be resolved into
multiple narrower rings (e.g., Facchini et al. 2020; Pérez et al. 2020). It is
also possible that fainter rings at large radii are not recovered with current
sensitivity. Observations with a spatial resolution that is smaller than the
pressure scale height and with deeper sensitivity are needed to have complete
characterization of substructure type and occurrence rate in disks around mid-
to-late M star.

Smooth disks 2M0450 and CIDA 12 are candidates to host small cavities (Section
3.2.2). Though their visibility profiles do not show clear null points (Figure
2), small cavities below $5$ au could remain hidden at our current angular
resolution, especially for CIDA 12, which has an inclination of $\sim
65^{\circ}$ (see analysis of high inclination disks by van der Marel et al.,
2022).

ALMA revealed a central depression of large grains in four structured disks in
our sample. The SED complements ALMA imaging by providing information on small
grains, since deficits of near infrared emission indicate an inner disk
clearing of small dust grains. Figure 5 shows the collected photometry and the
modelled SEDs with modeling processes given in Appendix D. In short, the
modeling uses a continuous disk component with a cavity sized $R_{\rm\mu
m,in}$ and an outer radius adopting $R_{\rm 90\%,dust}$ retrieved from
millimeter images (Table 3). The goal is not to have a stringent constraint on
disk parameters since the available infrared photometry are limited, but to
have a rough sense of the cavity sizes of small dust grains and how they
compare to the millimeter-sized grains observed by ALMA. Since we lack
knowledge of other parameters describing disk density profiles (e.g., power-
law index of surface density and pressure scale height), the cavity size is
degenerated with other parameters, including the power index for dust mass
density and the pressure scale height. Through experiments, the constrained
cavity sizes of small dust grains can differ by a few au under different
combinations of other parameters. However, the uncertainty on cavity sizes
does not affect the discussions below.

For the four disks with cavities in millimeter images, three disks show
evidence of inner dust clearing from their SEDs. 2M0412’s SED is consistent
with $R_{\rm\mu m,in}\sim$ 0.9 au for small dust grains, which is far smaller
than its $\sim 60$ au millimeter cavity. 2M0434 and 2M0508 are consistent with
$R_{\rm\mu m,in}\sim$ 4 au cavities. This test of simple SED fitting suggests
that micron-sized dust grains are not as depleted as mm-sized grains in the
cavity, consistent with a morphology that is also seen around larger and more
massive disks with both scattered light and ALMA images (Villenave et al.,
2019). $R_{\rm\mu m,in}$ for 2M0436 is around 0.03 au, consistent with dust
sublimation radius, so there is no evidence for clearing of small grains.
There is evidence of an inner disk around 2M0436 both from its brightness
profile (Figure 4) and visibility fitting result (Section 3.1). Models with
the same inclination as that from 2M0436’s ALMA image have a bit of difficulty
in reproducing the high flux at near infrared. These suggest that 2M0436 may
host a puffed up or misaligned inner disk (e.g., Dullemond et al. 2001; van
der Marel et al. 2018, 2022). For the two smooth disks, 2M0450 show no clear
small grain cavity with $R_{\rm\mu m,in}$ consistent with dust sublimation
radius, while CIDA 12 show evidence of a 0.3 au cavity. Further evaluation of
these potential cavities will require follow-up at much higher spatial
resolution and sensitivity. Scattered light imaging would likely require the
next-generation of ground-based telescopes, since these disks are too faint to
observe with current ground-based telescopes, and the cavities are too small
for JWST imaging.

Figure 5: Spectral energy distributions from observations (black dots) and
radiative transfer modeling (blue lines). The red lines are the stellar
photospheric emission applied with extinction. The orange dots are measurement
at 1.3mm in this work, for CIDA 12 the open black circle represents the 1.3mm
flux from Akeson et al. (2019). $R_{\rm in}$ at upper right of each panel are
cavity sizes from radiative transfer modeling for micron-sized dust grains.

###  4.3 Origins of small cavities

Rings are the most common substructures in our sample, similar to surveys on
disks around early-type stars as DSHARP (Andrews et al. 2018b; Huang et al.
2018b) and the Taurus survey (Long et al. 2018). With the current resolution
and sensitivity of our observations, the detected rings all appear axis-
symmetric, without significant azimuthally asymmetric features. Most rings
encircle cavities, except for the outer ring of 2M0412 (see Long et al. 2023
for a detailed analysis of the origins of substructures of the 2M0412 disk).
In short, planet-disk interaction is preferred over other mechanisms in
2M0412, with Saturn-mass planets capable of shaping the observed disk
morphology. However, a combination of dead zones and photoevaporation cannot
be ruled out. Below we focus on discussion of possible origins555We exclude
the dead zone outer edge as a possible origin of radial dust traps. In the up-
to-date picture where all the three non-ideal MHD effects are considered, the
level of MRI turbulence is damped smoothly in the outer disk (no abrupt
change). As a result, a dead zone outer edge will not develop (see, e.g.,
Figure 2 in Bai 2016). for substructures of 2M0434, 2M0436, and 2M0508, which
show small cavities around 10 au.

##### Condensation fronts.

Major volatiles in protoplanetary disks freeze onto dust grains from the gas
phase as the disk temperature decreases outward. Across these regions, which
are referred as condensation fronts or icelines, the dust opacity and critical
fragmentation velocity are expected to change and further ring/gap
substructures can be produced (e.g., Zhang et al. 2015; Okuzumi et al. 2016).

We consider an irradiated flared disk with a disk midplane temperature
formulated as Equation 6, where the flaring angle $\phi_{\rm fl}$ is assumed
to be 0.02 (e.g., Huang et al. 2018b; Dullemond et al. 2018). The iceline
location of a volatile can be calculated from a stellar luminosity and the
freezing point of the molecule. Since icelines of abundant species like $\rm
H_{2}O$ and $\rm NH_{3}$ are close to the host star ($<$ 2 au) and are not
resolvable by our observations, we focus on species with lower condensation
temperatures in Zhang et al. (2015), which are clathrate-hydrated CO and $\rm
N_{2}$ (41-46 K), CO (23-28 K) and $\rm N_{2}$ (12-15 K).

The comparison between iceline locations and brightness profiles from best-fit
models is shown in Figure 6. Interestingly, the peak emission radii are all
close to the CO icelines. However, inferring the iceline radius is a
challenging task. Besides the uncertainties in the estimation of disk midplane
temperature (e.g., Liu 2021), icelines can be thermally unstable and
dynamically evolve on time scales from 1000 to 10000 years (Owen, 2020). In
simulations, icelines have not been found to carve out cavities (see Pinilla
et al. 2017a for a simulation around a Herbig AeBe star). Deeper observations
of CO isotopologues are needed to better constrain disk mid-plane temperature,
and future simulations tailored to M stars could together help determine the
contribution of icelines to the cavities/rings detected.

Figure 6: Left panels (1)-(3): disk midplane temperature profiles calculated
using Eq. 6 are shown as black curves. Orange curves show the normalized
radial profile of best-fit models from visibility fitting. Shaded regions are
the iceline locations for $\rm N_{2}$ (light green), CO (light blue) and
clathrate-hydrated CO and $\rm N_{2}$ (light purple). Right panel: Disk radius
versus the square root of stellar luminosity. The peak model intensity
locations are shown for 2M0434 (cross), 2M0436 (plus) and 2M0508 (circle). The
shaded regions for icelines are the same as the left panels.

##### Photoevaporation and MHD wind.

When the mass loss rate by photoevaporation surpasses the accretion rate
through the disk, a gas and dust gap will open in the disk. Recent new
photoevaporation models by Picogna et al. (2019) predict transition disks with
accretion rates $\leq 10^{-9}M_{\odot}\ {\rm yr}^{-1}$ and maximum cavity
sizes about 30 au when surrounding solar-type stars. The lower mass loss rate
at large radii around very low-mass stars (Picogna et al., 2021) will likely
lead to a smaller maximu cavity size, which we expect to be above 10 au. This
is because the old photoevaporation models around $0.1M_{\odot}$ stars by Owen
et al. (2012) have cavities $\leq 10$ au and new mass loss profiles by Picogna
et al. (2019) are more efficient at removing material at larger disc radii,
which will lead to larger cavities. We therefore expect the maximum cavity
size around very low mass stars will fall between 10 and 30 au. 2M0434’s 10 au
cavity, 2M0436’s 12 au cavity and their accretion rates about a few $\times
10^{-10}M_{\odot}\ {\rm yr}^{-1}$ fall in the predicted parameter space.

Further evaluating the role of photoevaporation in the creation of cavities
will require additional observations. Possible discriminators include high
angular resolution imaging of the gas – in this paper, we have limited
sensitivity that allows us to only tentatively identify a cavity around one
source, 2M0508. In addition, empirical measurements of photoevaporation rates
through lines of O I and Ne II (Pascucci et al., 2011) but challenging to
interpret (e.g., Banzatti et al., 2019; Fang et al., 2023b; Rab et al., 2023),
given the possibility that these lines may form in a magnetothermal wind
(e.g., Wang et al., 2019).

Furthermore, accreting transition disks could be sustained by other
mechanisms. Gárate et al. (2021) showed that a disk with X-ray
photoevaporation in combination with dead zones can reproduce both the
accretion rates and gap sizes observed in transition disks. Their models
predict a long-lived compact inner disk, which is not seen in our ALMA images
under current resolution and sensitivity. MHD wind as another scenario, if the
large-scale magnetic field distribution is not significantly modified, could
generate a positive radial slope of the gas surface density (e.g., Suzuki et
al. 2016) which inhibits the inward drift of pebbles, and sustain the
accretion rates similar to classic T Tauri stars (e.g., Lesur 2021; Martel &
Lesur 2022).

##### Embedded planets

Giant planets can open gaps in the gas disks and form gas pressure bumps that
trap dust particles outside its orbit. Dust cavities will then form after the
inner dust disk materials are accreted onto the star. By applying the
analytical criterion for opening a gap in the gas (Crida et al., 2006), the
estimated masses of planets needed to open the small cavities in our disk
sample are: for viscous $\alpha\sim 10^{-4}$ and $10^{-3}$, 0.1 & 0.3 $M_{\rm
Jup}$ at $\sim 8$ au for 2M0434, 0.1 & 0.3 $M_{\rm Jup}$ at $\sim 10$ au for
2M0436, and 0.1 & 0.2 $M_{\rm Jup}$ at $\sim 6$ au for 2M0508. The locations
of planets are assumed to be five times the Hill radius
$r_{H}=r_{p}(M_{p}/3M_{*})^{1/3}$ away from the ring peak locations (e.g.,
Dodson-Robinson & Salyk 2011). Since only 2M0508 shows a tentative gas cavity,
these planet masses should be treated as upper limits. Moreover, if the disk
mass and angular momentum transport are dominated by MHD wind rather than
turbulent viscosity, planets can open wider and deeper gaps more easily (e.g.,
Elbakyan et al. 2022; Aoyama & Bai 2023; Wafflard-Fernandez & Lesur 2023). As
a result, the masses of embedded planets could be a factor of a few to ten
smaller than estimated above (Elbakyan et al. 2022).

Simulations of planet-disk interaction predict the segregation between
millimeter-sized dust continuum emission peak and gas emission peak (or
micron-sized dust, e.g., de Juan Ovelar et al. 2013; Facchini et al. 2018). We
find hints of such segregation between micron-sized and millimeter-sized dust
grains from SEDs (Section 4.2). While for gas emission, only 2M0508 shows
tentative inner gas depletion and the gas peak location overlap with the ring
location. Deeper sensitivity and higher angular resolution of molecular line
mapping are needed to resolve the gas depletion in the inner disk if present.

Taking disk dust masses in Section 3.4 and assuming a gas-to-dust mass ratio
of 100, the estimated current disk masses are 6.1 $M_{\rm Jup}$ for 2M0434,
0.5 $M_{\rm Jup}$ for 2M0436 and 0.8 $M_{\rm Jup}$ for 2M0508. If we assume a
constant accretion rate after disk formation, then the total mass accreted
onto the star from the disk is 1.5 $M_{\rm Jup}$ for 2M0434 and 4.5 $M_{\rm
Jup}$ for 2M0436 (using stellar ages with 50% spot coverage). The planet
masses needed to open cavities for 2M0434 and 2M0436 are within 10% of their
estimated total disk masses at initial stages. If 2M0508 has an accretion rate
$\sim 3\times 10^{-10}M_{\odot}\ {\rm yr}^{-1}$, similar to those in our
sample (Table 2.1), and that accretion rate has remained constant (most
likely, the accretion rate at earlier stages of formation was higher, e.g.,
Fischer et al. 2017; Fang et al. 2023a), its potential planet mass is also
within 10% of disk mass. Hence from the perspective of disk mass budget, the
three small cavity disks had enough materials to form those Saturn-mass
planets (e.g., Boss 2006; Lin et al. 2018).

If the observed cavities are carved out by giant planets, how does the current
detection rate of cavities around mid-to-late M stars compared to the
occurrence rate of exoplanets? The detection rate of cavities around mid-to-
late M stars in Taurus is around $24\%$, given the fraction of 7/13 for
cavities in disks (Section 4.2) and accounting for the disk fraction around
$45\%$ for M3-M6 stars in Taurus (Esplin & Luhman 2019; this latter correction
applies only if the age spread of the sample is much less than the disk
survival timescale). However, this is likely higher than the real cavity rate,
since the sample of 13 disks is slightly biased to the brighter end of disks
around M3-M6 stars. As for exoplanets around M-dwarfs, the CARMENES survey
obtained occurrence rates of giant planets with $M_{\rm p}\sin i$ >
$100M_{\oplus}$ and periods of 1-1000 days being $0.021^{+0.018}_{-0.011}$
planets per star and $0.045^{+0.021}_{-0.016}$ planets per star for stars less
and more massive than $0.337M_{\odot}$ respectively (Ribas et al. 2023). Wide-
orbit (>$10$ au) giant planets massive than Jupiter are found to be rare
around M-dwarfs (e.g., Bowler et al. 2015). Currently, the constraints on sub-
Jovian giants on orbits of a few au around M-dwarfs are limited. Measurement
around Sun-like stars show giant planet occurrence is enhanced by a factor of
four beyond 1 au compared to within 1 au (e.g., Fernandes et al. 2019; Fulton
et al. 2021), while it is unclear whether the result could be extrapolated to
M-dwarfs. Samples from both disk and exoplanet need to be developed to make a
robust comparison.

##### Binary companions

Binary stars could also carve out inner cavities in circumbinary disks, with
the cavity sizes expected to be 3-5 times the binary semimajor axis (e.g.,
Artymowicz & Lubow 1994; Miranda et al. 2017). Hence stellar companions within
4 au could be responsible for the 7-12 au cavities in the three disks. To our
knowledge, no evidence of stellar multiplicity have been found in our sample
(Kraus et al. 2011; Kraus & Hillenbrand 2012), although the scales of a few au
have yet to be explored. If the companions are sufficiently massive (mass
ratios $>$ 0.1-0.2), the binary would produce eccentric cavities (0.05-0.35
for circular binaries. Miranda et al. 2017; Ragusa et al. 2020) and cause a
shift between the stellar mass center and the geometric center of the cavity.
For the $\sim 10$ au cavities in our sample, resolutions better than the
binary separation are needed to detect the potential shift.

Among the mechanisms described above, photoevaporation and planets are more
preferred origins of small cavities (around 7 and 10 au) in our sample. To
distinguish between these possibilities, more deep and higher angular
resolution observations of molecular lines are needed. Furthermore, multi-
wavelength analysis of the spectral index at the cavity edge might help to
distinguish photoevaporation and planet scenarios (Picogna et al. 2023), where
a cavity created by photoevaporation is expected to have a higher spectral
index due to its lower dust filtering efficiency.

##  5 Summary

This paper presents high angular resolution ($\sim$50 mas, 8 au) ALMA Band 6
observations of six disks around mid-to-late type M stars in Taurus. We
characterize the disk continuum emission by fitting parametric models in the
visibility plane. We explore these disks’ global properties and the possible
origins of their substructures. The main findings are summarized as follows:

1\. We detect all six disks in millimeter continuum emission. 12CO is detected
in five disks with a marginal detection in 2M0450, and 13CO is detected in
2M0412 and 2M0434. Dust substructures are detected in 4 disks: 2M0412, 2M0434,
2M0436 and 2M0508, which are also the four brightest disks in our sample. We
perform fitting in the visibility plane, with rings modeled as radially
asymmetric Gaussian rings. With our current sensitivity ($\sim$ 40 mJy/beam)
and spatial resolution, 2M0412 show a large cavity surrounded by a ring at 60
au, followed by a gap and an outer dust ring at 110 au; the other 3 disks all
display a cavity surrounded by a dust ring at 10 au for 2M0434, 12 au for
2M0436 and 7 au for 2M0508.

2\. The Nuker profile fitting on the two compact smooth disks 2M0450 and CIDA
12 show that they may hold small cavities, with sizes being $\sim$ 1.7 au for
2M0450 and $\sim$ 5.7 au for CIDA 12. Current spatial resolution and
sensitivity may hide these small and shallow cavities. Nuker profiles can
reproduce the flux for CIDA 12 better than Gaussian profile, with the flux
from the Nuker profile being twice as much as Gaussian’s flux. Higher spatial
resolution is needed to resolve the potential cavities.

3\. The structured disks are brighter than the smooth disks in our sample, the
presence of dust trapping by substructures could account for this. 2M0436 and
2M0508 deviate from the flatter $L_{\rm mm}-M_{*}$ relation for structured
disks, which could be due to a combination of partially optical thickness and
boulder formation (Pinilla et al. 2018, 2020). Disks around 2M0436, 2M0450,
CIDA 12 and 2M0508 fall at the faint end of $L_{\rm mm}\propto R^{2}_{\rm
eff}$ relation, consistent with being radial drift dominated. After measuring
the ratios between gas disk size $R_{\rm gas,90\%}$ and the dust disk size
$R_{\rm dust,90\%}$, 2M0436 and 2M0450 show large ratios $R_{\rm
gas,90\%}/R_{\rm dust,90\%}\sim 5$, indicating very efficient dust radial
drift in these two disks.

4\. All structured disks in our sample show a central cavity, indicating the
clearing of millimeter-sized dust grains. Radiative transfer fits to SEDs
providing a rough constraint on the clearing of micron-sized dust particles.
Of the 4 structured disks, only 2M0436 does not show evidence of clearing of
micron dust, with the modeled cavity size being the dust sublimation radius.
2M0412 show a 0.9 au cavity size far smaller than its millimeter cavity.
2M0434 and 2M0508 show 4 au cavities for micron dust. For the 2 smooth disks,
CIDA 12 show a 0.3 au cavity consistent with its near-IR deficit. These SED-
derived cavity sizes are smaller than their counterpart in millimeter images.

5\. Various mechanisms could be responsible for the observed substructures.
For 2M0412, Long et al. (2023) shows that Saturn-mass planets are able to
carve out its cavity and gap. The large cavity alone can also be shaped by a
combination of dead zone and photoevaporation. 2M0412’s large rings are
unlikely to be related to icelines. For the other three disks, Saturn-mass
planets or photoevaporation could be responsible for their cavities around 10
au. Icelines could also play a role but their locations are highly uncertain.
High resolution mapping of gas emission and multi-wavelength analysis of the
spectral index at the edge of cavities, could help to distinguish different
mechanisms.

Our sample covers a wide range of millimeter luminosity, the four cavity disks
together with the other two harboring potential small cavities suggest that
substructures are likely ubiquitous in disks around mid-to-late M stars.
Current exoplanet statistics does not rule out that all these observed
cavities are created by giant planets.

High spatial resolution imaging of disks around mid-to-late M stars are still
very limited in number. The more frequent small cavities around 10 au around
those disks may originate from various mechanisms, which will make them great
laboratories for testing planet formation and disk physics.

##  6 Acknowledgments

G.J.H and Y.S. are supported by the National Key R&D Program of China No.
2019YFA0405100 and by grant 12173003 from the National Natural Science
Foundation of China. Support for F.L. was provided by NASA through the NASA
Hubble Fellowship grant #HST-HF2-51512.001-A awarded by the Space Telescope
Science Institute, which is operated by the Association of Universities for
Research in Astronomy, Incorporated, under NASA contract NAS5-26555. D. H. is
supported by Center for Informatics and Computation in Astronomy (CICA) grant
and grant number 110J0353I9 from the Ministry of Education of Taiwan. D. H.
also acknowledges support from the National Science and Technology Council of
Taiwan through grant number 111B3005191. P. P. acknowledges support from the
UK Research and Innovation (UKRI) under the UK government’s Horizon Europe
funding guarantee from ERC (under grant agreement No 101076489). E.R.
acknowledges financial support from the European Union’s Horizon Europe
research and innovation programme under the Marie Skłodowska-Curie grant
agreement No. 101102964 (ORBIT-D). E.R. also acknowledges financial support
from the European Research Council (ERC) under the European Union’s Horizon
2020 research and innovation programme (grant agreement No. 864965, PODCAST).
D.J. is supported by NRC Canada and by an NSERC Discovery Grant. C.F.M. is
Funded by the European Union (ERC, WANDA, 101039452). Views and opinions
expressed are however those of the author(s) only and do not necessarily
reflect those of the European Union or the European Research Council Executive
Agency. Neither the European Union nor the granting authority can be held
responsible for them. G.D.M. acknowledges support from FONDECYT project
11221206, from ANID — Millennium Science Initiative — ICN12_009, and the ANID
BASAL project FB210003. L.A.C. acknowledges support from FONDECYT project
1211656 from ANID and the Millennium Science Initiative, Center Code
NCN2021_080.

This paper makes use of the following ALMA data: ADS/JAO.ALMA#2019.1.00566.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. The National Radio
Astronomy Observatory is a facility of the National Science Foundation
operated under cooperative agreement by Associated Universities, Inc.

##  Appendix A CO Channel maps and moment maps

Moment maps of ${}^{12}\rm CO$ and ${}^{13}\rm CO$ (if detected) are shown in
Figure 7 and corresponding channel maps are shown in Figures 8-13.

Figure 7: Moment 0 and 8 maps for (marginally) detected 12CO and 13CO. The
$3\sigma$ white contours from dust continuum images are overlaid. For each
panel, beams are shown as white ellipses at lower left and 20 au white scale
bars are shown at lower right.

Figure 8: Channel maps of 12CO (top) and 13CO (bottom) for 2M0412. Central
velocities are given in the upper left corner of each channel in unit of $\rm
km\ s^{-1}$. The $3\sigma$ white contours from dust continuum images are
overlaid. Beams as white ellipses and 20 au white scale bars are shown at
lower left channel, at lower left and lower right corner correspondingly.

Figure 9: Same as Figure 8 for 2M0434. Figure 10: Same as Figure 8 for 2M0436.
Figure 11: Same as Figure 8 for 2M0450. Figure 12: Same as Figure 8 for CIDA
12. Figure 13: Same as Figure 8 for 2M0508.

##  Appendix B Central emission in 2M0436?

In this section, we test whether 2M0436 has unresolved central emission
(Section 3.1). A central point emission is added to the model of a radially
asymmetric Gaussian ring, with point source intensity formulated as
$F_{1}\delta(r)$. $F_{1}$ is set to distributed uniformly between 0 and 1 mJy.
Figure 14 shows the best-fit model images and intensity profiles for the two
models. The central point emission model reaches a convergence with the point
flux being $0.04^{+0.03}_{-0.03}$ mJy (Figure 15), indicating a non-detection
(detection at 1$\sigma$) of central emission. Section 3.1 mentioned that it
could be the unresolved central emission blends with the outer ring emission,
resulting in a ring with its inner width than outer width which is contrary to
dust trapping in pressure bump. However, in our attempts to add central
emission component the outer ring still shows wider inner width and the ratio
between the two widths is similar to that from ring-only model (Section
3.2.1).

Figure 14: Two different models for 2M0436. From left to right: (1) Best-fit
image from the model with only an asymmetric Gaussian ring. (2) Best-fit image
from the model with a central point and an asymmetric Gaussian ring. (3)
Azimuthall averaged radial profiles retrieved from beam convolved images.
Figure 15: Cornerplots of key parameters of the model includes a central point
and an asymmetric Gaussian ring. From left to right, the three parameters are
the flux of point emission, the inner width of the ring and the outer width of
the ring.

##  Appendix C Nuker fitting results on 2M0450 & CIDA 12

The modeling procedure follows Section 3.1, with Nuker parameters’ priors
following Tripathi et al. (2017): $p(\log_{10}\alpha)=\mathcal{U}(0,2),\
p(\beta)=\mathcal{U}(0,10)$ and $\gamma$ sampled approximated uniform in
$(-3,2)$:

$p(\gamma)\propto\frac{1}{1+e^{-5(\gamma+3)}}-\frac{1}{1+e^{-15(\gamma-2)}}.$
(C1)

Figure 16 shows the best-fit Nuker model images, visibilities and random
radial profiles retrieved from posterior distributions.

Figure 16: Nuker fitting results of disks that seem to lack structures in the
image plane. Columns from left to right: (1) unconvolved best-fit Nuker model
images, (2) radial profiles retrieved from best-fit (orange) and 1000 randomly
selected (gray) Nuker models, (3) real part of deprojected and binned
visibilities from observations and best-fit models.

##  Appendix D SED Fitting

The SED data points are collected from GALEX GR6+7 (Bianchi et al., 2017),
SDSS DR17 (Abdurro’uf et al., 2022), PAN-STARRS (Flewelling et al., 2020),
GAIA DR3 (Gaia Collaboration et al., 2022), 2MASS (Skrutskie et al., 2006),
WISE (Wright et al., 2010), Spitzer (Luhman et al., 2010), Herschel (Marton et
al., 2017), and our measurements at 1.3mm (Literature value at 1.3mm for CIDA
12 in Akeson et al. 2019).

We fit the photospheric emission of each host star using the BT-Settl models
(Allard et al. 2011), with the temperature and extinction in Table 2.1 and
assuming a gravity of $\log g=4.0$. With the fitted spectrum and stellar
luminosity as inputs, we then use the RADMC-3D code (Dullemond et al., 2012)
to perform self-consistent radiative transfer to obtain dust temperature
profiles and then simulate the SEDs. For a flared disk model with well-mixed
dust and gas (for small grains this is approximately correct), the disk mass
density profile is given as:

$\rho(r,z)=\frac{\Sigma(r)}{\sqrt{2\pi}H_{p}}\exp{\left(-\frac{z^{2}}{2H_{p}^{2}}\right)}.$
(D1)

where $r$ is the distance to the central star measured in the disk mid-plane,
$z$ is the distance to the disk mid-plane, $H_{p}$ is the pressure scale
height and $\Sigma$ is the surface density integrated over vertical direction.
We adopt simple power-law for $\Sigma$ and $H_{p}$:

$\displaystyle\Sigma(r)$ $\displaystyle=\Sigma_{0}\left(\frac{r}{r_{\rm
out}}\right)^{p}.$ (D2) $\displaystyle H_{p}(r)$
$\displaystyle=H_{100}\left(\frac{r}{100{\rm\ au}}\right)^{\beta}.$ (D3)

where $H_{100}$ is the scale height of small dust grains at 100 au and
$p,\beta$ are the slope.

Six parameters are needed to describe the simple flared disk model: $R_{\rm
in},R_{\rm out},p,\beta,H_{100},M_{\rm dust}$. We fix $p$ to -1 and $R_{\rm
out}$ to the $R_{90\%}$ in millimeter images (Section 3.4 and Table 3).
$M_{\rm dust}$ is first set to the estimated value in Table 3 and adjusted to
match the $F_{\rm 1.3mm}$. We note that $M_{\rm dust}$ has negligible effect
on the infrared SEDs (hence the parameter $R_{\rm in}$) since dust is
optically thick at that wavelength. We test each combination of the other
three parameters with $R_{\rm in}$ between dust sublimation radius and mm
cavity size, $\beta$ in [1.0, 1.2] and $H_{100}$ in [5,15]. As for the dust
prescription, we adopt the model in (Birnstiel et al., 2018) with grain sizes
ranging from 0.01 $\mu\rm m$ to 1000 $\mu\rm m$. The parameters of our best-
fit models are listed in Table 5, and Figure 5 shows the comparison between
the models and observations.

Table 5: Best-fit models for SED Parameters | 2M0412 | 2M0434 | 2M0436 | 2M0450 | CIDA 12 | 2M0508 | Unit
---|---|---|---|---|---|---|---
$R_{\rm in}$ | 0.90 | 4 | 0.03 | 0.02 | 0.3 | 4 | au
$R_{\rm out}$ | 126 | 82 | 14 | 6 | 20 | 21 | au
$p$ | -1 | -1 | -1 | -1 | -1 | -1 | –
$\beta$ | 1.15 | 1.20 | 1.15 | 1.15 | 1.20 | 1.15 | –
$H_{100}$ | 10.0 | 7.5 | 12.5 | 10.0 | 6.3 | 7.5 | au

## References

  * Abdurro’uf et al. (2022) Abdurro’uf, Accetta, K., Aerts, C., et al. 2022, ApJS, 259, 35, doi: 10.3847/1538-4365/ac4414
  * Akeson et al. (2019) Akeson, R. L., Jensen, E. L. N., Carpenter, J., et al. 2019, ApJ, 872, 158, doi: 10.3847/1538-4357/aaff6a
  * Allard et al. (2011) Allard, F., Homeier, D., & Freytag, B. 2011, in Astronomical Society of the Pacific Conference Series, Vol. 448, 16th Cambridge Workshop on Cool Stars, Stellar Systems, and the Sun, ed. C. Johns-Krull, M. K. Browning, & A. A. West, 91, doi: 10.48550/arXiv.1011.5405
  * ALMA Partnership et al. (2015) ALMA Partnership, Brogan, C. L., Pérez, L. M., et al. 2015, ApJ, 808, L3, doi: 10.1088/2041-8205/808/1/L3
  * Andrews (2020) Andrews, S. M. 2020, ARA&A, 58, 483, doi: 10.1146/annurev-astro-031220-010302
  * Andrews et al. (2013) Andrews, S. M., Rosenfeld, K. A., Kraus, A. L., & Wilner, D. J. 2013, ApJ, 771, 129, doi: 10.1088/0004-637X/771/2/129
  * Andrews et al. (2018a) Andrews, S. M., Terrell, M., Tripathi, A., et al. 2018a, ApJ, 865, 157, doi: 10.3847/1538-4357/aadd9f
  * Andrews et al. (2018b) Andrews, S. M., Huang, J., Pérez, L. M., et al. 2018b, ApJ, 869, L41, doi: 10.3847/2041-8213/aaf741
  * Ansdell et al. (2017) Ansdell, M., Williams, J. P., Manara, C. F., et al. 2017, AJ, 153, 240, doi: 10.3847/1538-3881/aa69c0
  * Ansdell et al. (2016) Ansdell, M., Williams, J. P., van der Marel, N., et al. 2016, ApJ, 828, 46, doi: 10.3847/0004-637X/828/1/46
  * Ansdell et al. (2018) Ansdell, M., Williams, J. P., Trapman, L., et al. 2018, ApJ, 859, 21, doi: 10.3847/1538-4357/aab890
  * Aoyama & Bai (2023) Aoyama, Y., & Bai, X.-N. 2023, ApJ, 946, 5, doi: 10.3847/1538-4357/acb81f
  * Artymowicz & Lubow (1994) Artymowicz, P., & Lubow, S. H. 1994, ApJ, 421, 651, doi: 10.1086/173679
  * Bae et al. (2023) Bae, J., Isella, A., Zhu, Z., et al. 2023, in Astronomical Society of the Pacific Conference Series, Vol. 534, Astronomical Society of the Pacific Conference Series, ed. S. Inutsuka, Y. Aikawa, T. Muto, K. Tomida, & M. Tamura, 423
  * Bai (2016) Bai, X.-N. 2016, ApJ, 821, 80, doi: 10.3847/0004-637X/821/2/80
  * Banzatti et al. (2019) Banzatti, A., Pascucci, I., Edwards, S., et al. 2019, ApJ, 870, 76, doi: 10.3847/1538-4357/aaf1aa
  * Barenfeld et al. (2016) Barenfeld, S. A., Carpenter, J. M., Ricci, L., & Isella, A. 2016, ApJ, 827, 142, doi: 10.3847/0004-637X/827/2/142
  * Barenfeld et al. (2017) Barenfeld, S. A., Carpenter, J. M., Sargent, A. I., Isella, A., & Ricci, L. 2017, ApJ, 851, 85, doi: 10.3847/1538-4357/aa989d
  * Bianchi et al. (2017) Bianchi, L., Shiao, B., & Thilker, D. 2017, ApJS, 230, 24, doi: 10.3847/1538-4365/aa7053
  * Birnstiel et al. (2018) Birnstiel, T., Dullemond, C. P., Zhu, Z., et al. 2018, ApJ, 869, L45, doi: 10.3847/2041-8213/aaf743
  * Boss (1997) Boss, A. P. 1997, Science, 276, 1836, doi: 10.1126/science.276.5320.1836
  * Boss (2006) —. 2006, ApJ, 643, 501, doi: 10.1086/501522
  * Boss & Kanodia (2023) Boss, A. P., & Kanodia, S. 2023, ApJ, 956, 4, doi: 10.3847/1538-4357/acf373
  * Bowler et al. (2015) Bowler, B. P., Liu, M. C., Shkolnik, E. L., & Tamura, M. 2015, ApJS, 216, 7, doi: 10.1088/0067-0049/216/1/7
  * Bryant et al. (2023) Bryant, E. M., Bayliss, D., & Van Eylen, V. 2023, MNRAS, 521, 3663, doi: 10.1093/mnras/stad626
  * Burn et al. (2021) Burn, R., Schlecker, M., Mordasini, C., et al. 2021, A&A, 656, A72, doi: 10.1051/0004-6361/202140390
  * CASA Team et al. (2022) CASA Team, Bean, B., Bhatnagar, S., et al. 2022, PASP, 134, 114501, doi: 10.1088/1538-3873/ac9642
  * Cieza et al. (2019) Cieza, L. A., Ruíz-Rodríguez, D., Hales, A., et al. 2019, MNRAS, 482, 698, doi: 10.1093/mnras/sty2653
  * Cieza et al. (2021) Cieza, L. A., González-Ruilova, C., Hales, A. S., et al. 2021, MNRAS, 501, 2934, doi: 10.1093/mnras/staa3787
  * Crida et al. (2006) Crida, A., Morbidelli, A., & Masset, F. 2006, Icarus, 181, 587, doi: 10.1016/j.icarus.2005.10.007
  * Curone et al. (2022) Curone, P., Izquierdo, A. F., Testi, L., et al. 2022, A&A, 665, A25, doi: 10.1051/0004-6361/202142748
  * Cutri et al. (2003) Cutri, R. M., Skrutskie, M. F., van Dyk, S., et al. 2003, 2MASS All Sky Catalog of point sources.
  * Daemgen et al. (2015) Daemgen, S., Bonavita, M., Jayawardhana, R., Lafrenière, D., & Janson, M. 2015, ApJ, 799, 155, doi: 10.1088/0004-637X/799/2/155
  * de Juan Ovelar et al. (2013) de Juan Ovelar, M., Min, M., Dominik, C., et al. 2013, A&A, 560, A111, doi: 10.1051/0004-6361/201322218
  * Diaz Trigo et al. (2019) Diaz Trigo, M., Carpenter, J., Maude, L., Miura, R., & Plunkett, A. 2019, ALMA Cycle 7 Proposer’s Guide, doi: 10.5281/zenodo.4511962
  * Dodson-Robinson & Salyk (2011) Dodson-Robinson, S. E., & Salyk, C. 2011, ApJ, 738, 131, doi: 10.1088/0004-637X/738/2/131
  * Dong et al. (2018) Dong, R., Liu, S.-y., Eisner, J., et al. 2018, ApJ, 860, 124, doi: 10.3847/1538-4357/aac6cb
  * Dullemond et al. (2001) Dullemond, C. P., Dominik, C., & Natta, A. 2001, ApJ, 560, 957, doi: 10.1086/323057
  * Dullemond et al. (2012) Dullemond, C. P., Juhasz, A., Pohl, A., et al. 2012, RADMC-3D: A multi-purpose radiative transfer tool, Astrophysics Source Code Library, record ascl:1202.015.  http://ascl.net/1202.015
  * Dullemond et al. (2018) Dullemond, C. P., Birnstiel, T., Huang, J., et al. 2018, ApJ, 869, L46, doi: 10.3847/2041-8213/aaf742
  * Elbakyan et al. (2022) Elbakyan, V., Wu, Y., Nayakshin, S., & Rosotti, G. 2022, MNRAS, 515, 3113, doi: 10.1093/mnras/stac1774
  * Esplin & Luhman (2019) Esplin, T. L., & Luhman, K. L. 2019, AJ, 158, 54, doi: 10.3847/1538-3881/ab2594
  * Facchini et al. (2018) Facchini, S., Pinilla, P., van Dishoeck, E. F., & de Juan Ovelar, M. 2018, A&A, 612, A104, doi: 10.1051/0004-6361/201731390
  * Facchini et al. (2020) Facchini, S., Benisty, M., Bae, J., et al. 2020, A&A, 639, A121, doi: 10.1051/0004-6361/202038027
  * Fang et al. (2023a) Fang, M., Pascucci, I., Edwards, S., et al. 2023a, ApJ, 945, 112, doi: 10.3847/1538-4357/acb2c9
  * Fang et al. (2023b) Fang, M., Wang, L., Herczeg, G. J., et al. 2023b, Nature Astronomy, 7, 905, doi: 10.1038/s41550-023-02004-x
  * Fernandes et al. (2019) Fernandes, R. B., Mulders, G. D., Pascucci, I., Mordasini, C., & Emsenhuber, A. 2019, ApJ, 874, 81, doi: 10.3847/1538-4357/ab0300
  * Fischer et al. (2017) Fischer, W. J., Megeath, S. T., Furlan, E., et al. 2017, ApJ, 840, 69, doi: 10.3847/1538-4357/aa6d69
  * Flewelling et al. (2020) Flewelling, H. A., Magnier, E. A., Chambers, K. C., et al. 2020, ApJS, 251, 7, doi: 10.3847/1538-4365/abb82d
  * Foreman-Mackey et al. (2013) Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306, doi: 10.1086/670067
  * Francis et al. (2020) Francis, L., Johnstone, D., Herczeg, G., Hunter, T. R., & Harsono, D. 2020, AJ, 160, 270, doi: 10.3847/1538-3881/abbe1a
  * Fulton et al. (2021) Fulton, B. J., Rosenthal, L. J., Hirsch, L. A., et al. 2021, ApJS, 255, 14, doi: 10.3847/1538-4365/abfcc1
  * Gaia Collaboration et al. (2022) Gaia Collaboration, Vallenari, A., Brown, A. G. A., et al. 2022, arXiv e-prints, arXiv:2208.00211, doi: 10.48550/arXiv.2208.00211
  * Gaia Collaboration et al. (2023) —. 2023, A&A, 674, A1, doi: 10.1051/0004-6361/202243940
  * Gárate et al. (2021) Gárate, M., Delage, T. N., Stadler, J., et al. 2021, A&A, 655, A18, doi: 10.1051/0004-6361/202141444
  * González-Ruilova et al. (2020) González-Ruilova, C., Cieza, L. A., Hales, A. S., et al. 2020, ApJ, 902, L33, doi: 10.3847/2041-8213/abbcce
  * Guo et al. (2018) Guo, Z., Herczeg, G. J., Jose, J., et al. 2018, ApJ, 852, 56, doi: 10.3847/1538-4357/aa9e52
  * Hardegree-Ullman et al. (2019) Hardegree-Ullman, K. K., Cushing, M. C., Muirhead, P. S., & Christiansen, J. L. 2019, AJ, 158, 75, doi: 10.3847/1538-3881/ab21d2
  * Hashimoto et al. (2021) Hashimoto, J., Dong, R., & Muto, T. 2021, AJ, 161, 264, doi: 10.3847/1538-3881/abf431
  * Hendler et al. (2020) Hendler, N., Pascucci, I., Pinilla, P., et al. 2020, ApJ, 895, 126, doi: 10.3847/1538-4357/ab70ba
  * Herczeg & Hillenbrand (2008) Herczeg, G. J., & Hillenbrand, L. A. 2008, ApJ, 681, 594, doi: 10.1086/586728
  * Herczeg & Hillenbrand (2014) —. 2014, ApJ, 786, 97, doi: 10.1088/0004-637X/786/2/97
  * Hildebrand (1983) Hildebrand, R. H. 1983, QJRAS, 24, 267
  * Huang et al. (2018a) Huang, J., Andrews, S. M., Pérez, L. M., et al. 2018a, ApJ, 869, L43, doi: 10.3847/2041-8213/aaf7a0
  * Huang et al. (2018b) Huang, J., Andrews, S. M., Dullemond, C. P., et al. 2018b, ApJ, 869, L42, doi: 10.3847/2041-8213/aaf740
  * Huang et al. (2020) —. 2020, ApJ, 891, 48, doi: 10.3847/1538-4357/ab711e
  * Jennings et al. (2022) Jennings, J., Tazzari, M., Clarke, C. J., Booth, R. A., & Rosotti, G. P. 2022, MNRAS, 514, 6053, doi: 10.1093/mnras/stac1770
  * Kraus & Hillenbrand (2012) Kraus, A. L., & Hillenbrand, L. A. 2012, ApJ, 757, 141, doi: 10.1088/0004-637X/757/2/141
  * Kraus et al. (2011) Kraus, A. L., Ireland, M. J., Martinache, F., & Hillenbrand, L. A. 2011, ApJ, 731, 8, doi: 10.1088/0004-637X/731/1/8
  * Kurtovic et al. (2021) Kurtovic, N. T., Pinilla, P., Long, F., et al. 2021, A&A, 645, A139, doi: 10.1051/0004-6361/202038983
  * Lantz et al. (2004) Lantz, B., Aldering, G., Antilogus, P., et al. 2004, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 5249, Optical Design and Engineering, ed. L. Mazuray, P. J. Rogers, & R. Wartmann, 146–155, doi: 10.1117/12.512493
  * Lauer et al. (1995) Lauer, T. R., Ajhar, E. A., Byun, Y. I., et al. 1995, AJ, 110, 2622, doi: 10.1086/117719
  * Lesur (2021) Lesur, G. R. J. 2021, A&A, 650, A35, doi: 10.1051/0004-6361/202040109
  * Lin et al. (2018) Lin, J. W., Lee, E. J., & Chiang, E. 2018, MNRAS, 480, 4338, doi: 10.1093/mnras/sty2159
  * Liu et al. (2019) Liu, B., Lambrechts, M., Johansen, A., & Liu, F. 2019, A&A, 632, A7, doi: 10.1051/0004-6361/201936309
  * Liu (2021) Liu, Y. 2021, Research in Astronomy and Astrophysics, 21, 164, doi: 10.1088/1674-4527/21/7/164
  * Long et al. (2018) Long, F., Pinilla, P., Herczeg, G. J., et al. 2018, ApJ, 869, 17, doi: 10.3847/1538-4357/aae8e1
  * Long et al. (2019) Long, F., Herczeg, G. J., Harsono, D., et al. 2019, ApJ, 882, 49, doi: 10.3847/1538-4357/ab2d2d
  * Long et al. (2022) Long, F., Andrews, S. M., Rosotti, G., et al. 2022, ApJ, 931, 6, doi: 10.3847/1538-4357/ac634e
  * Long et al. (2023) Long, F., Ren, B. B., Wallack, N. L., et al. 2023, ApJ, 949, 27, doi: 10.3847/1538-4357/acc843
  * Longarini et al. (2023) Longarini, C., Lodato, G., Bertin, G., & Armitage, P. J. 2023, MNRAS, 519, 2017, doi: 10.1093/mnras/stac3653
  * Luhman (2018) Luhman, K. L. 2018, AJ, 156, 271, doi: 10.3847/1538-3881/aae831
  * Luhman et al. (2010) Luhman, K. L., Allen, P. R., Espaillat, C., Hartmann, L., & Calvet, N. 2010, ApJS, 186, 111, doi: 10.1088/0067-0049/186/1/111
  * Luhman et al. (2017) Luhman, K. L., Mamajek, E. E., Shukla, S. J., & Loutrel, N. P. 2017, AJ, 153, 46, doi: 10.3847/1538-3881/153/1/46
  * MacGregor et al. (2018) MacGregor, M. A., Weinberger, A. J., Wilner, D. J., Kowalski, A. F., & Cranmer, S. R. 2018, ApJ, 855, L2, doi: 10.3847/2041-8213/aaad6b
  * Mairs et al. (2019) Mairs, S., Lalchand, B., Bower, G. C., et al. 2019, ApJ, 871, 72, doi: 10.3847/1538-4357/aaf3b1
  * Manara et al. (2023) Manara, C. F., Ansdell, M., Rosotti, G. P., et al. 2023, in Astronomical Society of the Pacific Conference Series, Vol. 534, Protostars and Planets VII, ed. S. Inutsuka, Y. Aikawa, T. Muto, K. Tomida, & M. Tamura, 539, doi: 10.48550/arXiv.2203.09930
  * Martel & Lesur (2022) Martel, É., & Lesur, G. 2022, A&A, 667, A17, doi: 10.1051/0004-6361/202142946
  * Marton et al. (2017) Marton, G., Calzoletti, L., Perez Garcia, A. M., et al. 2017, arXiv e-prints, arXiv:1705.05693, doi: 10.48550/arXiv.1705.05693
  * Mercer & Stamatellos (2020) Mercer, A., & Stamatellos, D. 2020, A&A, 633, A116, doi: 10.1051/0004-6361/201936954
  * Miguel et al. (2020) Miguel, Y., Cridland, A., Ormel, C. W., Fortney, J. J., & Ida, S. 2020, MNRAS, 491, 1998, doi: 10.1093/mnras/stz3007
  * Miotello et al. (2023) Miotello, A., Kamp, I., Birnstiel, T., Cleeves, L. C., & Kataoka, A. 2023, in Astronomical Society of the Pacific Conference Series, Vol. 534, Protostars and Planets VII, ed. S. Inutsuka, Y. Aikawa, T. Muto, K. Tomida, & M. Tamura, 501, doi: 10.48550/arXiv.2203.09818
  * Miranda et al. (2017) Miranda, R., Muñoz, D. J., & Lai, D. 2017, MNRAS, 466, 1170, doi: 10.1093/mnras/stw3189
  * Morales et al. (2019) Morales, J. C., Mustill, A. J., Ribas, I., et al. 2019, Science, 365, 1441, doi: 10.1126/science.aax3198
  * Mulders et al. (2021) Mulders, G. D., Drążkowska, J., van der Marel, N., Ciesla, F. J., & Pascucci, I. 2021, ApJ, 920, L1, doi: 10.3847/2041-8213/ac2947
  * Mulders et al. (2015) Mulders, G. D., Pascucci, I., & Apai, D. 2015, ApJ, 798, 112, doi: 10.1088/0004-637X/798/2/112
  * Oke & Gunn (1982) Oke, J. B., & Gunn, J. E. 1982, PASP, 94, 586, doi: 10.1086/131027
  * Okuzumi et al. (2016) Okuzumi, S., Momose, M., Sirono, S.-i., Kobayashi, H., & Tanaka, H. 2016, ApJ, 821, 82, doi: 10.3847/0004-637X/821/2/82
  * Owen (2020) Owen, J. E. 2020, MNRAS, 495, 3160, doi: 10.1093/mnras/staa1309
  * Owen et al. (2012) Owen, J. E., Clarke, C. J., & Ercolano, B. 2012, MNRAS, 422, 1880, doi: 10.1111/j.1365-2966.2011.20337.x
  * Paneque-Carreño et al. (2021) Paneque-Carreño, T., Pérez, L. M., Benisty, M., et al. 2021, ApJ, 914, 88, doi: 10.3847/1538-4357/abf243
  * Pascucci et al. (2011) Pascucci, I., Sterzik, M., Alexander, R. D., et al. 2011, ApJ, 736, 13, doi: 10.1088/0004-637X/736/1/13
  * Pascucci et al. (2016) Pascucci, I., Testi, L., Herczeg, G. J., et al. 2016, ApJ, 831, 125, doi: 10.3847/0004-637X/831/2/125
  * Pecaut & Mamajek (2013) Pecaut, M. J., & Mamajek, E. E. 2013, ApJS, 208, 9, doi: 10.1088/0067-0049/208/1/9
  * Pérez et al. (2016) Pérez, L. M., Carpenter, J. M., Andrews, S. M., et al. 2016, Science, 353, 1519, doi: 10.1126/science.aaf8296
  * Pérez et al. (2020) Pérez, S., Casassus, S., Hales, A., et al. 2020, ApJ, 889, L24, doi: 10.3847/2041-8213/ab6b2b
  * Picogna et al. (2021) Picogna, G., Ercolano, B., & Espaillat, C. C. 2021, MNRAS, 508, 3611, doi: 10.1093/mnras/stab2883
  * Picogna et al. (2019) Picogna, G., Ercolano, B., Owen, J. E., & Weber, M. L. 2019, MNRAS, 487, 691, doi: 10.1093/mnras/stz1166
  * Picogna et al. (2023) Picogna, G., Schäfer, C., Ercolano, B., et al. 2023, MNRAS, 523, 3318, doi: 10.1093/mnras/stad1504
  * Pinilla (2022) Pinilla, P. 2022, European Physical Journal Plus, 137, 1206, doi: 10.1140/epjp/s13360-022-03384-1
  * Pinilla et al. (2013) Pinilla, P., Birnstiel, T., Benisty, M., et al. 2013, A&A, 554, A95, doi: 10.1051/0004-6361/201220875
  * Pinilla et al. (2012) Pinilla, P., Birnstiel, T., Ricci, L., et al. 2012, A&A, 538, A114, doi: 10.1051/0004-6361/201118204
  * Pinilla et al. (2020) Pinilla, P., Pascucci, I., & Marino, S. 2020, A&A, 635, A105, doi: 10.1051/0004-6361/201937003
  * Pinilla et al. (2017a) Pinilla, P., Pohl, A., Stammler, S. M., & Birnstiel, T. 2017a, ApJ, 845, 68, doi: 10.3847/1538-4357/aa7edb
  * Pinilla et al. (2015) Pinilla, P., van der Marel, N., Pérez, L. M., et al. 2015, A&A, 584, A16, doi: 10.1051/0004-6361/201526655
  * Pinilla et al. (2017b) Pinilla, P., Pérez, L. M., Andrews, S., et al. 2017b, ApJ, 839, 99, doi: 10.3847/1538-4357/aa6973
  * Pinilla et al. (2018) Pinilla, P., Tazzari, M., Pascucci, I., et al. 2018, ApJ, 859, 32, doi: 10.3847/1538-4357/aabf94
  * Pinilla et al. (2021) Pinilla, P., Kurtovic, N. T., Benisty, M., et al. 2021, A&A, 649, A122, doi: 10.1051/0004-6361/202140371
  * Pollack et al. (1996) Pollack, J. B., Hubickyj, O., Bodenheimer, P., et al. 1996, Icarus, 124, 62, doi: 10.1006/icar.1996.0190
  * Rab et al. (2023) Rab, C., Weber, M. L., Picogna, G., Ercolano, B., & Owen, J. E. 2023, ApJ, 955, L11, doi: 10.3847/2041-8213/acf574
  * Ragusa et al. (2020) Ragusa, E., Alexander, R., Calcino, J., Hirsh, K., & Price, D. J. 2020, MNRAS, 499, 3362, doi: 10.1093/mnras/staa2954
  * Ribas et al. (2023) Ribas, I., Reiners, A., Zechmeister, M., et al. 2023, A&A, 670, A139, doi: 10.1051/0004-6361/202244879
  * Rosotti et al. (2019) Rosotti, G. P., Booth, R. A., Tazzari, M., et al. 2019, MNRAS, 486, L63, doi: 10.1093/mnrasl/slz064
  * Sanchis et al. (2020) Sanchis, E., Testi, L., Natta, A., et al. 2020, A&A, 633, A114, doi: 10.1051/0004-6361/201936913
  * Schlecker et al. (2022) Schlecker, M., Burn, R., Sabotta, S., et al. 2022, A&A, 664, A180, doi: 10.1051/0004-6361/202142543
  * Skrutskie et al. (2006) Skrutskie, M. F., Cutri, R. M., Stiening, R., et al. 2006, AJ, 131, 1163, doi: 10.1086/498708
  * Somers et al. (2020) Somers, G., Cao, L., & Pinsonneault, M. H. 2020, ApJ, 891, 29, doi: 10.3847/1538-4357/ab722e
  * Stefánsson et al. (2023) Stefánsson, G., Mahadevan, S., Miguel, Y., et al. 2023, Science, 382, 1031, doi: 10.1126/science.abo0233
  * Suzuki et al. (2016) Suzuki, T. K., Ogihara, M., Morbidelli, A., Crida, A., & Guillot, T. 2016, A&A, 596, A74, doi: 10.1051/0004-6361/201628955
  * Tazzari et al. (2018) Tazzari, M., Beaujean, F., & Testi, L. 2018, MNRAS, 476, 4527, doi: 10.1093/mnras/sty409
  * Tazzari et al. (2021) Tazzari, M., Clarke, C. J., Testi, L., et al. 2021, MNRAS, 506, 2804, doi: 10.1093/mnras/stab1808
  * Tobin et al. (2020) Tobin, J. J., Sheehan, P. D., Megeath, S. T., et al. 2020, ApJ, 890, 130, doi: 10.3847/1538-4357/ab6f64
  * Trapman et al. (2019) Trapman, L., Facchini, S., Hogerheijde, M. R., van Dishoeck, E. F., & Bruderer, S. 2019, A&A, 629, A79, doi: 10.1051/0004-6361/201834723
  * Tripathi et al. (2017) Tripathi, A., Andrews, S. M., Birnstiel, T., & Wilner, D. J. 2017, ApJ, 845, 44, doi: 10.3847/1538-4357/aa7c62
  * Tychoniec et al. (2020) Tychoniec, Ł., Manara, C. F., Rosotti, G. P., et al. 2020, A&A, 640, A19, doi: 10.1051/0004-6361/202037851
  * Valenti et al. (1993) Valenti, J. A., Basri, G., & Johns, C. M. 1993, AJ, 106, 2024, doi: 10.1086/116783
  * van der Marel & Mulders (2021) van der Marel, N., & Mulders, G. D. 2021, AJ, 162, 28, doi: 10.3847/1538-3881/ac0255
  * van der Marel et al. (2016) van der Marel, N., Verhaar, B. W., van Terwisga, S., et al. 2016, A&A, 592, A126, doi: 10.1051/0004-6361/201628075
  * van der Marel et al. (2013) van der Marel, N., van Dishoeck, E. F., Bruderer, S., et al. 2013, Science, 340, 1199, doi: 10.1126/science.1236770
  * van der Marel et al. (2018) van der Marel, N., Williams, J. P., Ansdell, M., et al. 2018, ApJ, 854, 177, doi: 10.3847/1538-4357/aaaa6b
  * van der Marel et al. (2022) van der Marel, N., Williams, J. P., Picogna, G., et al. 2022, arXiv e-prints, arXiv:2204.08225, doi: 10.48550/arXiv.2204.08225
  * van der Plas et al. (2016) van der Plas, G., Ménard, F., Ward-Duong, K., et al. 2016, ApJ, 819, 102, doi: 10.3847/0004-637X/819/2/102
  * Villenave et al. (2019) Villenave, M., Benisty, M., Dent, W. R. F., et al. 2019, A&A, 624, A7, doi: 10.1051/0004-6361/201834800
  * Wafflard-Fernandez & Lesur (2023) Wafflard-Fernandez, G., & Lesur, G. 2023, arXiv e-prints, arXiv:2305.11784, doi: 10.48550/arXiv.2305.11784
  * Wang et al. (2019) Wang, L., Bai, X.-N., & Goodman, J. 2019, ApJ, 874, 90, doi: 10.3847/1538-4357/ab06fd
  * Wang & Chen (2019) Wang, S., & Chen, X. 2019, ApJ, 877, 116, doi: 10.3847/1538-4357/ab1c61
  * Weber et al. (2023) Weber, P., Pérez, S., Zurlo, A., et al. 2023, ApJ, 952, L17, doi: 10.3847/2041-8213/ace186
  * Wright et al. (2010) Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868, doi: 10.1088/0004-6256/140/6/1868
  * Yamaguchi et al. (2021) Yamaguchi, M., Tsukagoshi, T., Muto, T., et al. 2021, ApJ, 923, 121, doi: 10.3847/1538-4357/ac2bfd
  * Zhang et al. (2015) Zhang, K., Blake, G. A., & Bergin, E. A. 2015, ApJ, 806, L7, doi: 10.1088/2041-8205/806/1/L7
  * Zhang et al. (2023) Zhang, S., Kalscheur, M., Long, F., et al. 2023, arXiv e-prints, arXiv:2305.03862, doi: 10.48550/arXiv.2305.03862
  * Zhu & Dong (2021) Zhu, W., & Dong, S. 2021, ARA&A, 59, 291, doi: 10.1146/annurev-astro-112420-020055
  * Zormpas et al. (2022) Zormpas, A., Birnstiel, T., Rosotti, G. P., & Andrews, S. M. 2022, A&A, 661, A66, doi: 10.1051/0004-6361/202142046