An Atacama Large Millimeter/submillimeter Array Survey of Chemistry in Disks around M4–M5 Stars

Table 1 presents the stellar characteristics of the survey sample, which consists of five disks around M4–M5 stars. The disks were chosen from multiple star-forming regions, avoiding bias toward any one molecular cloud. Two of the disks are from the Taurus region, two from the Lupus region, and one from the Chamaeleon I region. All disks are within 128–192 pc (Gaia Collaboration et al. 2016, 2018). Millimeter dust disk sizes are ∼30–63 au (where the size denotes the radius that contains ∼90% of the emission, similar to the methods of Ansdell et al. 2018). We refer to the millimeter dust disk as the "pebble disk" throughout the remainder of the paper. The host stars of the disks are ∼1–3 Myr in age and range in spectral type, stellar luminosity, estimated stellar mass, and stellar effective temperature within M4–M5.5, ∼0.05–0.32 L_⊙, ∼0.14–0.23 M_⊙, and ∼3000–3300 K, respectively. All stellar characteristics were derived in the literature from continuum photometry, spectral energy distribution (SED) fits, scaling relations, and/or stellar evolutionary models ¹⁰ (Luhman et al. 2010; Andrews et al. 2013; Ward-Duong et al. 2018; Manara et al. 2017; Alcalá et al. 2017).

Figure 1 compares the dust continuum fluxes of our M4–M5 disk sample to the dust continuum fluxes of solar-type and Herbig Ae disks (compiled from the chemistry surveys of Huang et al. 2017, Bergner et al. 2019, 2020, and Pegues et al. 2020 and discussed further in Section 4.3). Figure 1 also compares our sample to the dust continuum flux detections from Andrews et al. (2013) and from Pascucci et al. (2016), which were extracted from the Taurus and Chamaeleon I star-forming regions, respectively. The M4–M5 disks in our sample allow us to probe a lower stellar mass regime than the previous chemistry surveys dominated by solar-type disks. The majority of our M4–M5 disks fall within the brighter regime of the dust continuum flux versus stellar mass relationship. Our sample thus provides an initial rather than representative look into resolved disk chemistry around M4–M5 stars and in the low-mass M-star regime.

Zoom In Zoom Out Reset image size

We note that the five disks in the M4–M5 disk sample were selected specifically for their bright ¹²CO (J = 3–2) and ¹³CO (J = 3–2) detections in existing ALMA disk surveys (e.g., Ansdell et al. 2016; Long et al. 2017; G. van der Plas et al. 2021, in preparation). These selection criteria are similar to those used by early exploratory chemistry surveys of solar-type and Herbig Ae disks (e.g., the DISCS survey, Öberg et al. 2010, 2011b), and so our M4–M5 disk sample can reasonably be compared to disks around more massive stars from these previous surveys. The solar-type and Herbig Ae disks form the core of our literature sample, which is overplotted in Figure 1. We compare our M4–M5 disk sample to this literature sample in Section 4.3.

We survey the overall chemistry of these disks by targeting 1.1 mm and 1.3 mm lines of CO isotopologs and precursor organic molecules. Our primary target molecular lines are ¹²CO (J = 2–1), ¹³CO (J = 2–1), C¹⁸O (J = 2–1), C₂H (N = 3–2, J = 7/2–5/2), DCN (J = 3–2), HCN (J = 3–2), and H₂CO 3₀₃–2₀₂. ¹¹ Together, these lines provide constraints on the C/N/O chemistry and organic inventories in the emitting layers of the disks (e.g., Bergin et al. 2016; Cleeves et al. 2018; Miotello et al. 2019). These lines have also been previously observed toward disks around solar-type stars, allowing direct comparison of the disk chemistry between these two types of stars. We additionally use observations of the hyperfine structure of the C₂H and HCN lines to estimate C₂H and HCN column densities, excitation temperatures, and optical depths. The molecular characteristics of the target lines are summarized in Table 2.

Table 3 presents the observational characteristics of the sample. All five disks were observed with ALMA during Project 2017.1.01107.S from 2017 December through 2018 September. Each individual observation (i.e., each execution) used a total of 41–46 antennas, providing a minimum baseline of 15 m and maximum baselines from 1241 to 2517 m. Individual on-source observation times spanned 16–40 minutes. The angular resolution and maximum angular scales ranged from 012–029 and 235–423, respectively.

Initial calibration (including flux, phase, and bandpass calibration) was performed by ALMA/NAASC using standard procedures. We then used the Common Astronomy Software Applications package (CASA) version 4.7.2 to self-calibrate each observation and image each molecular line. We self-calibrated each observation using the line-free continuum combined from all spectral windows. We uniformly used Briggs weighting and a robust value of 0.5. We performed two rounds of phase calibration and one round of amplitude calibration per observation when able. Solution intervals for the disks J1545-3417 and Sz 69 ranged from 10 to 100 s for each observation. These intervals were chosen as the minimum values in this range, which still maximized the number of solutions with a signal-to-noise ratio ≥2. Due to low continuum fluxes, we used solution intervals from 70 to 210 s for the five observations toward FP Tau. We were unable to self-calibrate on the continuum for J0432+1827 and J1100-7619.

We used CASA's uvcontsub and clean functions to (1) subtract the continuum from each spectral window and (2) image each molecular line. We created cleaning masks by hand to cover the molecular emission and little else in each channel. For fainter lines, we recycled cleaning masks from brighter lines and cleaned down to 3σ rather than to the deeper 2σ, where σ is the average rms across non-emission channels, to avoid creating image artifacts. In order to balance sensitivity and resolution, we imaged each line using Briggs weighting and a robust value of 0.5. The resulting synthesized beams range from 018 to 047 in size. We imaged the dust continuum and molecular line emission toward the disk observed with the highest resolution (J1100-7619) with a pixel size of 002, while we imaged the emission toward the other four disks in the sample (FP Tau, J0432+1827, J1545-3417, and Sz 69) with a pixel size of 004.

We imaged the following line+disk combinations at a channel resolution of 0.2 km s⁻¹: all line emission toward the bright disk J1100-7619; the ¹²CO 2–1, ¹³CO 2–1, and HCN 3–2 lines toward FP Tau and J0432+1827; and the ¹²CO 2–1 and ¹³CO 2–1 lines toward J1545-3417 and Sz 69. All other line+disk combinations were imaged at a channel resolution of 0.4 km s⁻¹ to increase the signal-to-noise ratio. We analyzed velocity channels for each disk+line only within a specific velocity range. We selected these ranges based on visual inspection of the Keplerian-masked spectra (Section 3.1) to encompass visible emission. These ranges also included additional "buffer" channels to ensure that we had incorporated all emission. Note that the velocity ranges for the C₂H 3–2 and HCN 3–2 lines were extended to include all hyperfine components.

We estimated the disk centers using two-dimensional Gaussian fits to the dust continuum emission. We used Keplerian masks (e.g., Rosenfeld et al. 2013; Yen et al. 2016; Salinas et al. 2017) to extract emission for each disk (software: jpegues 2019) across all emission channels for each molecular line. Keplerian masks decrease the amount of noise that is incorporated into the integrated emission. The mask parameters used for this survey are given in Appendix A. We used the Keplerian-masked channels to then generate spectra, velocity-integrated emission maps, radial profiles, and integrated fluxes.

We took the average rms across 1000 2'' × 2'' random samples of non-emission channels to be the channel rms. We estimated the noise for the integrated fluxes via bootstrapping over 1000 Keplerian-masked random samples of non-emission channels. For the velocity-integrated emission maps, we used the median of 1000 random rms maps (σ_map; Bergner et al. 2018) as a representation of the noise. For the radial profiles, we adapted the approach of Bergner et al. (2018) and estimated the noise per ring as ${\sigma }_{\mathrm{map}}/\sqrt{N}$. N is the number of independent measurements in the ring, assumed to be N = (# of pixels within each ring's width)/(# of pixels in the beam area). Once the count per ring was less than the circumference of the beam, N was fixed to the value of N using the circumference.

The C₂H (N = 3–2, J = 7/2–5/2), C₂H (N = 3–2, J=5/2–3/2), and HCN 3–2 lines exhibit hyperfine structure, and so the emission spectra can be decomposed into individual hyperfine emission components (e.g., the 262.00426 GHz and 262.00648 GHz lines given in Rows 4 and 5 of Table 2 are hyperfine components for the C₂H (N = 3–2, J = 7/2–5/2) line). Since the disk J1100-7619 is bright and well resolved, we were able to fit models to the resolved hyperfine structure (e.g., Hily-Blant et al. 2013; Estalella 2017; Pety 2018). We could then extract excitation temperatures, column densities, and optical depths from the parameters of the fitted models, through each pixel in the image. By fitting the spectrum per pixel rather than across a larger disk area, we reduced the line width, and therefore the line blending, of the hyperfine components in the spectrum. This in turn permitted tighter constraints on the model fit. We used a hyperfine pixel-by-pixel procedure based on the work of Bergner et al. (2019), which in turn was adapted from the fitting procedure of Estalella (2017) and the calculations of Mangum & Shirley (2015). See Appendix B for the model equations, the fitting procedure, examples of model fits, and an important discussion of intrinsic uncertainties.

For lines where we could not use the hyperfine fitting method to measure excitation temperatures, optical depths, and column densities directly, we could still estimate the column densities over a given disk region, using assumptions of the excitation temperature and total optical depth (e.g., from pixel-by-pixel hyperfine fits to the brightest disk), as well as a measurement of the total emission across all hyperfine components within that given disk region. This methodology is described in detail in Appendix C. Uncertainties are discussed in Appendix B.

We used the following criteria to determine whether or not a molecular line is detected toward a disk:

• 1.  
  Emission is ≥3σ in the velocity-integrated emission map.
• 2.  
  Peak emission within the Keplerian masks is ≥3σ in at least three velocity channels.

Lines that satisfy both criteria are detected, while lines that satisfy at least one criterion are tentatively detected. Lines that fail both criteria are nondetected.

Based on these criteria, we detected ¹²CO (J = 2–1), ¹³CO (J = 2–1), C¹⁸O (J = 2–1), HCN (J = 3–2), C₂H (N = 3–2, F = 7/2–5/2), and C₂H (N = 3–2, F = 5/2–3/2) toward all 5/5 M4–M5 disks. We detected DCN (J = 3–2) toward one disk (J1100-7619), and we detected H₂CO 3₀₃–2₀₂ toward two disks (J0432+1827 and J1100-7619). For H₂CO 3–2 toward two other disks (FP Tau and Sz 69), and for DCN 3–2 toward four other disks (FP Tau, J0432+1827, J1545-3417, and Sz 69), the emission was tentatively detected and would be worth a follow-up observation with deeper integration.

Channel maps for the detections and tentative detections are presented in Appendix D. Table 4 presents the emission measurements for the dust continuum, and Table 5 presents the flux measurements for the target molecular lines. All detections, tentative detections, and nondetections are marked within Table 5.

Table 5. Fluxes and Detection Upper Limits for Target Molecular Lines

Disk	90% Em.	Integrated	Peak Flux	Integrated	Kep. Mask	Channel	Channel	Beam
	Radius	Flux	(mJy beam−1	Velocity Range	Extent	Width	rms	Size
	(au)	(mJy km s−1)	× km s−1)	(km s−1)	('')	(km s−1)	(mJy beam−1)	(P.A.)
(1)	(2)	(3)	(4)	(5)	(6)	(7)	(8)	(9)
12CO (J = 2–1)
FP Tau	55	781 ± 11	247 ± 4.3	2.0–14.6	1.24	0.20	3.4	030 × 023 (−63)
J0432+1827	126	1847 ± 13	129 ± 2.9	0.6–10.6	2.16	0.20	3.3	024 × 020 (−148)
J1100-7619	186	1306 ± 6.1	69 ± 1.4	2.1–7.3	2.04	0.20	2.5	040 × 025 (−11)
J1545-3417	98	738 ± 10	141 ± 3.3	−1.0–10.0	2.00	0.20	3.7	027 × 024 (710)
Sz 69	138	2001 ± 15	192 ± 3.0	1.2–9.4	2.68	0.20	3.8	027 × 024 (730)
13CO (J = 2–1)
FP Tau	46	232 ± 8.3	100 ± 3.5	2.0–14.6	0.68	0.20	3.2	030 × 024 (−72)
J0432+1827	81	537 ± 9.4	81 ± 2.9	0.6–10.6	1.32	0.20	3.0	025 × 021 (−134)
J1100-7619	136	417 ± 6.0	38 ± 1.4	2.2–7.4	1.92	0.20	2.7	041 × 026 (−48)
J1545-3417	80	218 ± 7.7	57 ± 3.1	−1.0–10.0	1.32	0.20	3.3	028 × 025 (714)
Sz 69	136	345 ± 13	60 ± 2.6	1.2–9.4	2.68	0.20	3.5	028 × 025 (727)
C18O (J = 2–1)
FP Tau	44	91 ± 6.7	54 ± 3.3	1.9–14.7	0.56	0.40	2.1	030 × 024 (−71)
J0432+1827	58	148 ± 9.0	36 ± 2.4	0.4–10.8	1.32	0.40	2.0	025 × 021 (−135)
J1100-7619	116	71 ± 3.6	12 ± 1.2	2.2–7.4	1.36	0.20	2.1	042 × 027 (−23)
J1545-3417	123	89 ± 9.5	21 ± 2.9	−1.1–10.1	1.96	0.40	2.3	028× 025 (691)
Sz 69	37	38 ± 6.6	22 ± 2.5	1.4–9.4	0.96	0.40	2.2	029 × 025 (710)
C2H (N = 3–2, J = $\tfrac{7}{2}$–$\tfrac{5}{2}$)
FP Tau	63	140 ± 18	81 ± 4.7	1.1–14.7	1.44	0.40	3.2	044 × 034 (−1765)
J0432+1827	104	288 ± 19	39 ± 4.3	−1.6–10.0	1.32	0.40	3.0	034 × 024 (144)
J1100-7619	167	1180 ± 17	22 ± 1.9	−0.4–7.4	1.48	0.20	2.9	018 × 014 (−95)
J1545-3417	56	61 ± 12	47 ± 3.0	−1.9–9.7	1.16	0.40	2.6	037 × 033 (724)
Sz 69	39	28 ± 8.1	25 ± 3.7	−2.6–9.0	0.40	0.40	2.8	037 × 032 (729)
HCN (J = 3–2)
FP Tau	48	272 ± 17	164 ± 4.4	2.0–14.6	1.20	0.20	4.2	038 × 025 (−1738)
J0432+1827	85	600 ± 18	125 ± 2.9	−0.2–11.4	1.36	0.20	3.8	034 × 024 (143)
J1100-7619	137	1655 ± 24	67 ± 1.6	−0.2–9.8	1.64	0.20	2.6	018 × 014 (−91)
J1545-3417	88	120 ± 14	51 ± 2.6	−1.1–10.1	1.36	0.40	2.4	037 × 033 (701)
Sz 69	39	42 ± 8.6	40 ± 3.4	−0.2–11.0	0.44	0.40	2.6	037 × 033 (710)
DCN (J = 3–2)
FP Tau	⋯	≲5.7	≲12	1.9–14.7	0.12	0.40	2.0	031 × 024 (−97)
J0432+1827	⋯	≲18	≲9.1	0.4–10.8	0.48	0.40	2.0	026 × 021 (−165)
J1100-7619	62	15 ± 2.7	8.3 ± 1.0	2.1–7.3	0.72	0.20	1.6	047 × 034 (−33)
J1545-3417	⋯	≲14	≲8.7	−1.1–10.1	0.36	0.40	2.2	029 × 026 (733)
Sz 69	⋯	≲8.7	≲8.7	1.4–9.4	0.16	0.40	2.1	029 × 026 (750)
H2CO 303–202
FP Tau	⋯	≲15	≲15	1.9–14.7	0.36	0.40	1.8	030 × 024 (−100)
J0432+1827	136	80 ± 8.0	6.5 ± 2.1	0.4–10.8	1.48	0.40	1.7	026 × 021 (−177)
J1100-7619	219	100 ± 3.7	7.7 ± 1.0	2.2–7.4	2.12	0.20	1.6	041 × 027 (−28)
J1545-3417	⋯	<14	<7.2	−1.1–10.1	0.40	0.40	1.9	029 × 026 (733)
Sz 69	⋯	≲11	≲6.9	1.4–9.4	0.24	0.40	1.9	029 × 026 (747)

Notes. Em. and Kep. are abbreviations of Emission and Keplerian, respectively. The emitting radii (Column 2) are the radii containing 90% of the molecular line emission and are presented only for detections. The integrated and peak fluxes were measured within the Keplerian masks (Column 6; Appendix A). Note the difference in units between the two quantities. The errors for the integrated and peak fluxes and the channel rms were estimated via bootstrapping of 1000 random samples of non-emission channels. Samples for the integrated and peak fluxes were extracted from within the Keplerian masks, while the channel rms was extracted from within 2'' × 2'' regions. Uncertainties do not include ∼15% systematic flux calibration uncertainties. The 3σ upper limits are given for tentative detections and nondetections and are marked with ≲ and <, respectively.

Download table as:  ASCIITypeset image

Figure 2 displays the millimeter dust continuum emission and velocity-integrated molecular line emission maps. Figures 3 and 4 display the corresponding radial profiles and spectra, respectively. All dust continua appear smooth at this resolution (i.e., there are no cavities in the millimeter dust continuum emission). We use the 90% emitting radii of the dust continuum to represent the pebble disk edges. In terms of pebble disk size, J0432+1827 and J1100-7619 are the largest disks, while Sz 69 is the smallest disk. The 90% emitting radii for all CO isotopolog, C₂H 3–2, and HCN 3–2 emission are either comparable to or exceed the pebble disk sizes, indicating that the gas of the disks is typically extended relative to the millimeter dust.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We use simple two-dimensional Gaussian fits to the dust continuum emission to estimate the center of each disk, and we qualitatively characterize peaks of the molecular line emission with respect to that disk center. All five disks in our sample have centrally peaked CO isotopolog emission distributions. The ¹³CO and C¹⁸O emission distributions extend beyond the pebble disk edge (i.e., beyond ∼30–60 au) for 4/5 and 3/5 disks, respectively. The disks FP Tau, J1545-3417, and Sz 69 are significantly affected by cloud contamination, which is evident in the CO spectra and the CO channel maps (Appendix D). The CO emission substructure seen beyond the pebble disk edges may be a byproduct of this contamination but could also be a contribution from an extended disk structure or disk–envelope interaction. The disk J1100-7619 shows asymmetries in the ¹²CO spectrum. It is not clear what is causing these asymmetries, as the channel maps do not suggest cloud contamination.

Two of the five disks have dips or holes in the C₂H 3–2 emission. The inner edge of the hole toward J1100-7619 is aligned with the edge of the pebble disk, between 46 and 63 au based on the 1.1 and 1.3 mm dust continuum. This means that the majority of the C₂H emission is coming from beyond the pebble disk (e.g., Bergin et al. 2016). The hole in C₂H emission toward the second disk, J0432+1827, is off-center by ∼01, appearing in Figure 2 but not in Figure 3. J0432+1827 additionally shows an emission plateau beyond the pebble disk. C₂H 3–2 emission toward the remaining three disks in our sample appears centrally peaked but is likely barely resolved. We clearly distinguish each hyperfine emission component (Rows 4 and 5 of Table 2) in the C₂H 3–2 spectrum toward J1100-7619, while the C₂H and HCN hyperfine components for all other disks are blended (Figure 4).

Zoom In Zoom Out Reset image size

All disks have centrally peaked HCN 3–2 emission morphologies. The two well-resolved disks, J0432+1827 and J1100-7619, have more compact HCN emission compared to the C₂H emission. Both disks also show plateaus or shelves in the HCN emission beyond the pebble disk edge. The remaining three disks may also show substructures exterior to the pebble disk at higher spatial resolution. DCN 3–2 is only detected toward J1100-7619, where the emission is centrally peaked and extends just beyond the edge of the pebble disk.

When detected, H₂CO 3–2 emission shows a hole or depression near the disk center. J0432+1827 has a central hole in H₂CO emission and suggestive rings that extend past the dust continuum. J1100-7619 has a hole in H₂CO emission that is off-center by <01 (which appears in Figure 2 but not in Figure 3), with two rings that peak at and beyond the pebble disk edge.

We now compare relative line fluxes (i.e., the flux for one molecular line with respect to the flux of a different molecular line) for the M4–M5 disks in this work to relative line fluxes measured for solar-type and Herbig Ae disks in the literature. The literature disk sample was compiled from the chemistry surveys of Huang et al. (2017), Bergner et al. (2019, 2020), and Pegues et al. (2020). We use only the molecular line+disk pairs in these surveys that were detected with ALMA, leading to 13 unique disks in all with at least two detections of C¹⁸O, C₂H, DCN, HCN, and H₂CO line emission. C¹⁸O 2–1 was detected with ALMA toward 11/13 disks, C₂H 3–2 toward 10/13 disks, DCN 3–2 toward 10/13 disks, HCN 3–2 toward 7/13 disks, and either H₂CO 3₀₃–2₀₂ or H₂CO 4₀₄–3₀₃ toward 13/13 disks. We consider both H₂CO 3₀₃–2₀₂ and H₂CO 4₀₄–3₀₃ line fluxes from the literature, because both lines have similar flux behaviors (Pegues et al. 2020). Stellar masses and spectral types for the T Tauri disks in this combined literature sample range within ∼0.4–1.8 M_⊙ and M1–G7, respectively. Two Herbig Ae disks (HD 163296 and MWC 480, A-star disks with stellar masses of ∼1.8–2.0 M_⊙) are included in the literature disk sample. We include the two Herbig Ae disks in relevant figures and tables for completeness, but since there are only two of them, we focus mainly on comparisons with the solar-type disk sample in later discussion.

Figure 5 plots distance-normalized C¹⁸O fluxes as a function of stellar mass across the combined sample of disks. The combined sample appears to follow the same trend in C¹⁸O versus stellar mass. The addition of the M4–M5 disks suggests a positive dependence of C¹⁸O flux on stellar mass across the sample, which was not distinguishable from the solar-type and Herbig Ae disk samples alone.

Zoom In Zoom Out Reset image size

Figure 6 compares distance-normalized line fluxes between the M4–M5 and literature disk samples. The Spearman correlation coefficients (r), which show how well the data can be described by a monotonic function, are shown when the associated p-value is ≤0.01. All correlation coefficients and measures of statistical significance are presented and explained in Appendix E. The line fluxes are often lower for the M4–M5 disks than for the solar-type and Herbig Ae disks. J1100-7619, which is particularly bright in C₂H and HCN emission, is a notable exception. Despite differences in individual fluxes, trends found between fluxes for solar-type and Herbig Ae disks seem to apply to M4–M5 disks. The strongest correlations across the M4–M5, solar-type, and Herbig Ae disk samples are found between C₂H 3–2 and HCN 3–2 (r = 0.92, p-value < 0.001) and between H₂CO and HCN 3–2 (r = 0.92, p-value < 0.001). Notably C₂H 3–2 and HCN 3–2 are perfectly correlated for the M4–M5 disks alone (r_M4−M5 = 1.00, p-value < 0.001). For the combined sample, the trends of the C₂H 3–2 versus C¹⁸O 2–1 line fluxes and the HCN 3–2 versus C¹⁸O 2–1 line fluxes are similar to each other (r = 0.81, p-value < 0.001 and r = 0.88, p-value < 0.001, respectively). We also find a weak correlation between H₂CO line fluxes and C¹⁸O 2–1 (r = 0.69, p-value = 0.01). All other pairs of line fluxes do not appear to be significantly correlated.

Zoom In Zoom Out Reset image size

We note that while we expect H₂CO and DCN line emission to be optically thin, studies have shown that C₂H and HCN line emission is often optically thick (e.g., Bergner et al. 2019). It is thus possible that the C₂H 3–2 and HCN 3–2 molecular line fluxes trace the distribution size rather than the underlying molecular abundance. We investigate this possibility in Appendix F, and we find evidence that the strong correlation between the C₂H 3–2 and HCN 3–2 fluxes likely cannot be explained by optical depth effects alone.

We plot a subset of the disk flux ratios for the M4–M5, solar-type, and Herbig Ae disks against stellar mass in Figure 7. The flux ratio subset includes molecular line emission relative to C¹⁸O 2–1 emission, which tells us how the molecular line emission relates to the disk gas. We also include the strongly correlated C₂H 3–2 versus HCN 3–2 and the deuterated fraction DCN 3–2 versus HCN 3–2. We note that the C₂H 3–2 / HCN 3–2 flux ratios for the M4–M5 disks are similar in value to each other. The C₂H 3–2 / HCN 3–2 flux ratios for the M4–M5 disks are high relative to the typical C₂H 3–2 / HCN 3–2 flux ratios for the solar-type disks. Otherwise, we find no clear trend for any of the disk molecular line flux ratios in either the individual or the combined disk samples. The ratios across all samples appear flat relative to stellar mass and also to stellar luminosity (not shown). The lack of any trends between molecular line flux ratios and stellar mass, despite the clear trends between the millimeter dust continuum fluxes, C¹⁸O fluxes, and stellar masses in Figures 1 and 5, suggests that the underlying disk chemistry is similar for the M4–M5 and solar-type stars in the outer disk regions probed by our disk-integrated fluxes.

Zoom In Zoom Out Reset image size

J1100-7619 presents the brightest molecular emission in our sample and was observed with the highest spatial resolution. We therefore used J1100-7619 to explore how the chemistry changes with disk radius around an M4–M5 star.

4.4.1. Radial Fluxes

Figure 8 shows the disk-averaged and azimuthally averaged (radial) molecular flux ratios toward J1100-7619 for the same molecule pairs in Figure 7. The (C₂H 3–2 / C¹⁸O 2–1) and (H₂CO 3–2 / C¹⁸O 2–1) ratios both increase across the disk. There is a notable bump in the (H₂CO 3–2 / C¹⁸O 2–1) ratio that appears at ∼80 au, just beyond the edge of the pebble disk. The (DCN 3–2/HCN 3–2) and (HCN 3–2/C¹⁸O 2–1) ratios are both roughly constant across the disk where a significant signal-to-noise ratio exists, with values of ∼(1.5–3.0) × 10⁻² and ∼20–30, respectively. The (C₂H 3–2/HCN 3–2) ratio increases across the pebble disk and then flattens out beyond the pebble disk edge. Past ∼100 au, the (C₂H 3–2/HCN 3–2) ratio is notably constant, with values of ∼0.8–0.9.

Zoom In Zoom Out Reset image size

4.4.2. Excitation Temperatures, Column Densities, and Optical Depths

Where a sufficient signal-to-noise ratio exists for J1100-7619, we used the pixel-by-pixel hyperfine fitting procedure described in Section 3.2 and Appendix B. This procedure allowed us to estimate excitation temperatures, column densities, and optical depths for C₂H and HCN within each pixel. Although we proceeded with this procedure, we note that these hyperfine fits are subject to significant intrinsic uncertainties (see Appendix B).

To increase the velocity resolution of the hyperfine fits, we reimaged the C₂H 3–2, C₂H (N = 3–2, J = 5/2–3/2), and HCN 3–2 line emission at channel widths of 0.14 km s⁻¹ prior to performing the fits. After the fitting procedure, we obtained the weighted averages of the pixel values within each annulus around the host star. The weighting process is described in Appendix B. We used this procedure to create weighted radial profiles of the excitation temperatures, column densities, and optical depths, which are plotted in Figure 9.

Zoom In Zoom Out Reset image size

The C₂H excitation temperatures are ∼10–20 K and appear roughly constant across the disk. These low temperatures suggest either thermal emission from the disk midplane or subthermal emission from an upper disk layer. The C₂H optical depths range within ∼1–4 and are also roughly constant with the radius. The C₂H column densities are ∼(10–1) × 10¹⁴ cm⁻² over ∼20–100 au and appear roughly flat with the radius.

The signal-to-noise ratio for HCN is sufficient only for probing interior to the pebble disk edge (at ∼20–60 au). In this region, HCN excitation temperatures are roughly constant. These temperatures are slightly warmer but still consistent with the C₂H excitation temperatures. HCN is optically thick in this region, with optical depths decreasing from ∼15 to 7. The HCN column densities are ∼(30–3) × 10¹³ cm⁻². The HCN column density profile decreases steeply in comparison to the C₂H column density profile over the same portion of the disk.

Bergner et al. (2019) used their pixel-by-pixel hyperfine procedure to derive C₂H and HCN excitation temperature and column density constraints for their sample of solar-type and Herbig Ae disks. When we compare the J1100-7619 results with their solar-type and Herbig Ae disk sample, we find that our estimated excitation temperatures and optical depths across the pebble disk for J1100-7619 are consistent with estimates for the disks around more massive stars.

We now use the methodology described in Section 3.3 and Appendix C to estimate disk-averaged C₂H, DCN, and HCN column densities for the entire M4–M5 disk sample. For each molecular line and disk, we took the extent of the emission to be where the emission radial profile decreases to 1σ. We extracted emission only within that radius, and we used that radius to define the angular size Ω of the disk (see Appendix C).

We assumed that the disk-averaged excitation temperatures range within 10–30 K for both molecules, with a median value of 20 K. For disk-averaged optical depths, we assumed 1–5 with a median of 2 for C₂H and 4–12 with a median of 8 for HCN. We also estimated disk-averaged column densities for DCN, assuming that DCN is optically thin and has the same excitation temperature and conditions as C₂H and HCN.

Figure 10 displays the estimated disk-averaged C₂H and HCN column densities for the M4–M5 disk sample using the disk-averaged median, minimum, and maximum excitation temperature and optical depth constraints. The median C₂H column density estimates are (5–10) × 10¹³ cm⁻² across the sample. Changing the excitation temperature and opacity results in factor of 2–4 higher or lower estimates of the column density. These values, which were estimated using the C₂H (N = 3–2, J = 7/2–5/2) line, are consistent with those estimated using the C₂H (N = 3–2, J = 5/2–3/2) line (not shown). The median HCN column density estimates span a factor of 3, within (4–12) × 10¹² cm⁻². Changes in optical depth from 4 to 12 lead to factors ≲3 in difference between estimates, while changes in excitation temperature from 10 to 30 K lead to factors of ∼3–10 (excluding Sz 69). For the bright and well-resolved disk J1100-7619, DCN is detected, and the estimated column densities (not shown) span ∼(4–10) × 10¹¹ cm⁻² at 10–30 K with a median of ∼7 × 10¹¹ cm⁻².

Zoom In Zoom Out Reset image size

Figure 10 also displays the estimated disk-averaged C₂H/HCN and DCN/HCN column density ratios. These are plotted in comparison to the 16th, median, and 84th percentiles derived from a sample of disk-averaged ratios from the literature for solar-type and Herbig Ae disks (Huang et al. 2017; Bergner et al. 2020). We used only the disk-averaged ratios from these references that were observed with ALMA and had comparable excitation temperature and optical depth assumptions. Detailed descriptions of the excitation temperatures and optical depths assumed in the literature are given in Appendix C.

For all M4–M5 disks, estimates of the disk-averaged C₂H/HCN ratio do not change significantly with optical depth. Median estimates of the ratio range from ∼7 for FP Tau to ∼11 for J1545-3417. Notably all C₂H/HCN estimates for the M4–M5 disks exceed the 84th percentile of the solar-type and Herbig Ae disk C₂H/HCN estimates. The DCN/HCN estimates are relatively well constrained for J1100-7619. This disk has a median DCN/HCN value of ∼0.072 and ranges within ∼0.04–0.1 across the extremes. This single value exceeds the 84th percentile disk-averaged DCN/HCN ratio found for the solar-type and Herbig Ae disks, although we would need more DCN detections toward M4–M5 disks to determine if this finding is representative of a larger sample.

The molecular emission morphologies we observe are products of the local chemistry and environment. Here we discuss the observed emission morphologies across our M4–M5 disk sample in the context of chemical pathways and environmental conditions that may form them. For recent, in-depth discussions of formation pathways and conditions for each molecule, refer to Huang et al. (2017; HCN and DCN), Bergner et al. (2019; C₂H and HCN), and Pegues et al. (2020; H₂CO).

5.1.1. C₂H and HCN

Observational studies have found significant positive correlations between HCN and C₂H line fluxes for solar-type disks (stellar masses of ∼0.5–1.8 M_⊙; Bergner et al. 2019), along with significant positive correlations between CN and C₂H line fluxes for both low-mass M-star (stellar masses <0.5 M_⊙) and solar-type disks (Miotello et al. 2019). These studies have concluded that the production of these molecules can be traced to the same or a similar chemical and/or physical driver(s). Theoretical studies of T Tauri disks (Du et al. 2015; Bergin et al. 2016) have found that depletions of oxygen relative to carbon and heightened UV irradiation are likely major drivers, as they can cause increased production of C₂H, other hydrocarbons, and cyanides in the outer T Tauri disk regions. Simulated observations of C₂H have shown that this increased production can be observed as rings in the hydrocarbon emission (Bergin et al. 2016).

It is difficult to isolate physical/chemical drivers and identify how they are tied to emission morphologies without detailed modeling. That being said, we note two intruiging similarities between C₂H and HCN beyond the pebble disk edges of our sample: morphologies and relative fluxes and column densities, which together point toward a related chemistry beyond the pebble disk. First, the two well-resolved disks in the sample (J0432+1827 and J1100-7619; Figure 11) show slope changes and suggestive plateaus in both C₂H and HCN emission beyond the pebble disk. Second, the disk (C₂H 3–2/HCN 3–2) flux ratios (Figures 6 and 7) appear perfectly correlated across the M4–M5 disk sample. The radially resolved (C₂H 3–2/HCN 3–2) flux ratio toward J1100-7619 helps us understand this correlation, because the ratio is virtually constant beyond the pebble disk edge (Figure 8). The estimated disk-averaged C₂H/HCN column density ratios (Figure 10), which are influenced most heavily by contributions from the outer disk, are similar across the M4–M5 disk sample.

Zoom In Zoom Out Reset image size

It is important to note how the C₂H and HCN column density profiles measured across the pebble disk for J1100-7619 (Figure 9) are different in shape, with the HCN column densities decreasing more steeply across the pebble disk than the C₂H column densities. This difference may indicate that the C₂H and HCN chemistry is decoupled and/or affected by different physical/chemical processes interior to the pebble disk edge. This further suggests that the similarities and correlations between C₂H and HCN morphologies, fluxes, and disk-averaged column densities discussed here are driven by emission beyond the pebble disk.

5.1.2. H₂CO, CO, and C₂H

5.2.1. Hydrocarbons and Hydrocyanides

5.2.2. Formaldehyde

We estimated excitation temperatures, column densities, and optical depths for C₂H and HCN in the disk J1100-7619 using a pixel-by-pixel hyperfine fitting procedure (Sections 3.2, 4.4.2). This procedure is based on the work of Bergner et al. (2019), which in turn was adapted from the fitting procedure of Estalella (2017) and the calculations of Mangum & Shirley (2015). By fitting to the spectrum through each pixel, we reduced line blending of the hyperfine components of the spectra, which allowed for tighter constraints on the fit. We also avoided systematic errors due to changes in line width, line center, and line peak as a function of disk location.

The model assumes that the molecule is in local thermodynamic equilibrium, that the line width is the same for each hyperfine component, and that the excitation conditions are homogeneous along the line of sight within the given region. We fit the same model across both C₂H 3–2 and C₂H (N = 3–2, J = 5/2–3/2), because the additional emission from C₂H (N = 3–2, J = 5/2–3/2) further improves our model constraints for C₂H. Note that more hyperfine components exist for C₂H 3–2 beyond those listed in Table 2; however, those unlisted hyperfine components have very weak line intensities and likely would not significantly affect our spectrum fits if we included them.

The hyperfine model contains four parameters: an amplitude ($\hat{A}*$), a transformed optical depth for the main hyperfine component (${\hat{\tau * }}_{{\rm{m}}}$), a line width (${\hat{{\rm{\Delta }}V}}_{\mathrm{line}}$), and a central velocity reference point (${\hat{V}}_{\mathrm{LSR}}$). Note that we use the ^∧ symbol to differentiate these model parameters from constants and other quantities in the equations. The amplitude and transformed optical depth parameters are directly related to the excitation temperature, column density, and total optical depth of the hyperfine line at that pixel.

We now write down the final equations used for this procedure. First, the spectrum model is in brightness temperature units and is a function of velocity. We denote this model as $\hat{H}\{{V}_{k}\}$, where we again use the ^∧ symbol to differentiate between the model spectrum $\hat{H}$ and the actual spectrum H later on:

$$\begin{eqnarray}&&\hat{H}\{{V}_{k}\}=\frac{{\hat{A}}^{* }}{{\hat{\tau }}_{{\rm{m}}}^{* }}\left(1-\exp (-{\tau }_{k})\right),\end{eqnarray} \tag{ B1 }$$

where V_k is the central velocity of channel k. ${\hat{A}}^{* }$ is the amplitude parameter (in brightness temperature units) and ${\hat{\tau }}_{{\rm{m}}}^{* }$ is the transformed optical depth parameter (bounded such that $0\lt {\hat{\tau }}_{{\rm{m}}}^{* }\lt 1$). τ _k is the optical depth of the kth channel, summed over all L hyperfine components:

$$\begin{eqnarray}&&{\tau }_{k}=\displaystyle \sum _{i}^{L}\frac{{\tau }_{{\rm{m}}}({S}_{i}{\mu }^{2}/{S}_{{\rm{m}}}{\mu }^{2}){G}_{i}\{{V}_{k}\}}{{\rm{\Delta }}{V}_{\mathrm{chan}}},\end{eqnarray} \tag{ B2 }$$

where τ_m is the optical depth of the main hyperfine component, calculated as

$$\begin{eqnarray}{\tau }_{{\rm{m}}}=\left\{\begin{array}{lc}-\mathrm{ln}(1-{\hat{\tau }}_{{\rm{m}}}^{* }) & \mathrm{if}\,{\hat{\tau }}_{{\rm{m}}}^{* }\gt 0.01\\ {\hat{\tau }}_{{\rm{m}}}^{* }+\displaystyle \frac{{\left({\hat{\tau }}_{{\rm{m}}* }\right)}^{2}}{2}+\displaystyle \frac{{\left({\hat{\tau }}_{{\rm{m}}* }\right)}^{3}}{3} & \mathrm{if}\,{\hat{\tau }}_{{\rm{m}}* }\leqslant 0.01,\end{array}\right.\end{eqnarray} \tag{ B3 }$$

and (S _i μ²/S_m μ²) is the strength of the ith hyperfine component relative to the main hyperfine component (Table 2). ΔV_chan is the width of each velocity channel (identical across all channels). G _i {V _k } is the value of the line profile for the ith hyperfine component at velocity channel k, assumed to be

$$\begin{eqnarray}&&{G}_{i}\{{V}_{k}\}={\hat{{\rm{\Delta }}V}}_{\mathrm{line}}\sqrt{\frac{\pi }{16\mathrm{ln}(2)}}\left(\mathrm{erf}({x}_{i,k}^{+})-\mathrm{erf}({x}_{i,k}^{-})\right),\end{eqnarray} \tag{ B4 }$$

where ${\hat{{\rm{\Delta }}V}}_{\mathrm{line}}$ is the line width parameter of the model, assumed to be the same for each hyperfine component, and erf(z) represents the error function for input z. If $| {x}_{i,k}^{+}-{x}_{i,k}^{-}| \lt {10}^{-4}$, then G _i {V _k } is instead approximated as

$$\begin{eqnarray}\begin{array}{rcl}{G}_{i}\{{V}_{{k}}\} & \approx & \frac{{\hat{{\rm{\Delta }}V}}_{\mathrm{line}}}{2\sqrt{\mathrm{ln}(2)}}\left({x}_{i,k}^{+}-{x}_{i,k}^{-}\right)\\ & & \times \exp \left(-{\left(\frac{{x}_{i,k}^{+}+{x}_{i,k}^{-}}{2}\right)}^{2}\right).\end{array}\end{eqnarray} \tag{ B5 }$$

${x}_{i,k}^{+}$ and ${x}_{i,k}^{-}$ are calculated as

$$\begin{eqnarray}\begin{array}{rcl}{x}_{i,k}^{\pm } & = & 2\sqrt{\mathrm{ln}(2)}\left(\frac{{V}_{k}\pm ({\rm{\Delta }}{V}_{\mathrm{chan}}/2)-{\hat{V}}_{\mathrm{LSR}}-{\rm{\Delta }}{V}_{i}}{{\hat{{\rm{\Delta }}V}}_{\mathrm{line}}}\right),\end{array}\end{eqnarray} \tag{ B6 }$$

where ${\hat{V}}_{\mathrm{LSR}}$ is the reference point parameter of the model. ΔV _i is the shift in velocity between the ith and main hyperfine components and can be calculated from the frequencies (Table 2).

We estimated the four parameters of the model (the amplitude parameter ${\hat{A}}^{* }$, the transformed optical depth parameter ${\hat{\tau }}_{{\rm{m}}}^{* }$, the line width parameter ${\hat{{\rm{\Delta }}V}}_{\mathrm{line}}$, and the velocity reference parameter ${\hat{V}}_{\mathrm{LSR}}$) per pixel in the image using the EnsembleSampler function from the emcee package (Foreman-Mackey et al. 2013). The log-likelihood function used for emcee included error in the resolved peaks (i.e., three peaks for HCN 3–2 and four peaks across C₂H 3–2 and C₂H 3–2) and error across all model values:

$$\begin{eqnarray}&&\mathrm{ln}{ \mathcal L }=-\mathop{}\limits_{k}{\left(\frac{H\{{V}_{k}\}-\hat{H}\{{V}_{k}\}}{{\epsilon }_{\mathrm{chan}}}\right)}^{2}-({\epsilon }_{\mathrm{peak}}\times 100),\end{eqnarray} \tag{ B7 }$$

where H{V _k } and $\hat{H}\{{V}_{k}\}$ are the actual and model spectra, respectively, at velocity channel k. _chan is the channel rms converted to brightness temperature. _peak is the summed relative error in the peaks of the resolved hyperfine components (where "relative error" is the absolute difference between the actual value and model value, divided by the actual value) and applies an additional penalty when the model and actual spectrum peaks are misaligned. The actual spectrum is converted from flux density units (i.e., units of power per square distance per frequency) to brightness temperature units using the following (e.g., Condon & Ransom 2016):

$$\begin{eqnarray}&&T=\left(\displaystyle \frac{{ \mathcal F }{c}_{0}^{2}}{2{k}_{{\rm{B}}}{\nu }_{0}^{2}{{\rm{\Omega }}}_{\mathrm{area}}}\right),\end{eqnarray} \tag{ B8 }$$

where ${ \mathcal F }$ is the flux density, c₀ is the speed of light, k_B is the Boltzmann constant, and ν₀ is the frequency of the main hyperfine component. Ω_area is the solid angle of the emitting area (in this case, the solid angle of the pixel). Note that in practice, we calculate ${ \mathcal F }$ by converting our observed fluxes from units of janskys per beam to janskys per pixel and finally to (W m⁻² Hz⁻¹/pixel), before applying Equation (B8) in SI units.

We used 100 emcee walkers per pixel and walked them for 1000 steps. Initial estimates of each walker for ${\hat{A}}^{* }$, ${\hat{\tau }}_{{\rm{m}}}^{* }$, ${\hat{{\rm{\Delta }}V}}_{\mathrm{line}}$, and ${\hat{V}}_{\mathrm{LSR}}$ were selected from uniform distributions in ranges of [peak of the emission in brightness temperature units ±10%], (0, 1), [the disk's systemic velocity ±10%], and [300 m s⁻¹, 500 m s⁻¹], respectively. We disregarded the first 950 steps for each walker to account for the burn-in phase, creating a sampling distribution of 100 × 50 for each parameter. We took the median parameter of each sampling distribution, ${\hat{A}}^{* (50)}$, ${\hat{\tau }}^{* (50)}$, ${\hat{{\rm{\Delta }}V}}_{\mathrm{line}}^{(50)}$, and ${\hat{V}}_{\mathrm{LSR}}^{(50)}$, to be the final estimate of each pixel's four parameters. We used the difference in the 16th and 84th percentile values relative to the median values as estimates of the lower and upper errors, respectively.

Figures 12 and 13 show example hyperfine emcee chains and fits for the C₂H and HCN molecular lines, respectively. The bottom panel of Figure 13 illustrates the importance of the HCN 3–2 satellite components discussed previously in anchoring the fit. This same panel also illustrates the intrinsic limitations of the underlying model at high optical depths. We discuss these limitations later in this section.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Adapting from Mangum & Shirley (2015), Estalella (2017), and Bergner et al. (2019), we estimated the excitation temperature T_ex from the following series of calculations:

$$\begin{eqnarray}&&{{ \mathcal T }}_{\mathrm{ex}}=\displaystyle \frac{{\hat{A}}^{* (50)}}{{\hat{\tau }}^{(50)}f}+{{ \mathcal B }}_{\mathrm{cont}}\end{eqnarray} \tag{ B9 }$$

$$\begin{eqnarray}&&{T}_{\mathrm{ex}}=\displaystyle \frac{h{\nu }_{0}}{{k}_{{\rm{B}}}\mathrm{ln}(h{\nu }_{0}/({k}_{{\rm{B}}}{{ \mathcal T }}_{\mathrm{ex}})+1)},\end{eqnarray} \tag{ B10 }$$

where Equation (B9) determines the Planck-corrected excitation temperature from the model parameters and Equation (B10) calculates the excitation temperature in brightness temperature units from the Planck-corrected excitation temperature. We calculated ${{ \mathcal B }}_{\mathrm{cont}}$ from Equation (B8), using the dust continuum flux at that pixel for ${ \mathcal F }$ (set to a very small number if the continuum at that pixel is 0) and the solid angle of the pixel for Ω_area. h is the Planck constant, k_B is the Boltzmann constant, and f is the filling factor (assumed to be f = 1 for J1100-7619). Note that unlike the previous studies, we did not calculate the Planck-corrected temperature of the background continuum prior to Equation (B9), as this would be an overcorrection of the background brightness temperature. The effect of this change is negligible at these wavelengths, continuum fluxes, and levels of uncertainty.

We then estimated the total optical depth τ as

$$\begin{eqnarray}&&\tau =-\displaystyle \frac{\mathrm{ln}(1-{\hat{\tau }}^{* (50)})}{{R}_{{\rm{m}}}},\end{eqnarray} \tag{ B11 }$$

where R_m is the strength of the main hyperfine transition relative to other same-level J transitions. R_i is generally calculated for any hyperfine transition as ${R}_{i}=\tfrac{({S}_{i}{\mu }^{2})}{{\sum }_{Jn=Ji}({S}_{n}{\mu }^{2})}$. As an example for a main hyperfine line, R_m = 0.527 for C₂H (N = 3–2, J = 7/2–5/2) is calculated from $\tfrac{2.2870}{2.2870+1.7110}$, which is the line intensity of the main line (N = 3–2, J = 7/2–5/2) divided by the sum of the line intensities for all other J = 7/2–5/2 lines.

Finally, we calculated N_tot as (Mangum & Shirley 2015)

$$\begin{eqnarray}\begin{array}{rcl}{\bar{N}}_{\mathrm{tot}} & = & \left(\frac{3{hQ}\{{T}_{\mathrm{ex}}\}}{8{\pi }^{3}({S}_{{\rm{m}}}{\mu }^{2})}\right)\left(\frac{\exp ({E}_{{\rm{u}}}/{T}_{\mathrm{ex}})}{\exp (h{\nu }_{0}/({k}_{{\rm{B}}}{T}_{\mathrm{ex}}))-1}\right)\\ & & \times \left(\frac{{\displaystyle \int }_{k}{\hat{H}}_{{\rm{m}}}{dV}}{{{ \mathcal T }}_{\mathrm{ex}}-{{ \mathcal B }}_{\mathrm{cont}}}\right)\left(\frac{\tau {R}_{{\rm{m}}}}{(1-\exp (-\tau {R}_{{\rm{m}}}))}\right),\end{array}\end{eqnarray} \tag{ B12 }$$

where Q{T_ex} is the value of the partition function at temperature T_ex (taken from the CDMS database; Endres et al. 2016). S_m μ² and E_u are the absolute line intensity and upper-level energy, respectively, of the main hyperfine component (Table 2). Note that the column density equation in Mangum & Shirley (2015) explicitly includes the upper degeneracy level g_u. For our calculations, that factor of degeneracy is already implicitly included in the CDMS values for S μ² (Table 2). Finally, ${\int }_{k}{\hat{H}}_{{\rm{m}}}{dV}$ is the integral of only the main hyperfine component of the model across all velocity channels. ${\hat{H}}_{{\rm{m}}}$ can be calculated by using Equation (B2) for only i = m and then plugging that result into Equation (B1). Note that any ith hyperfine component can be used in this way to determine N_tot; we use the main hyperfine component here because it has the strongest emission.

To generate the weighted radial profiles of Figure 9, we first defined ${\sigma }_{\mathrm{pix},Q}^{+}=\sqrt{{({Q}_{84}-{Q}_{50})}^{2}+{({Q}_{50}\times {f}_{50})}^{2}}$ and ${\sigma }_{\mathrm{pix},Q}^{-}\,=\sqrt{{({Q}_{16}-{Q}_{50})}^{2}+{({Q}_{50}\times {f}_{50})}^{2}}$ as the upper and lower errors, respectively, for each pixel and each measurement. Here Q refers to the measured quantity (e.g., the excitation temperature), and the subscripts 16, 50, and 84 refer to the 16th, median, and 84th percentiles, respectively. f₅₀ is the difference between the integrated fluxes of the median and actual spectrum, relative to the integrated flux of the actual spectrum.

We next extracted all valid pixels within each deprojected annulus around the disk center. "Valid" pixels are those that fulfill all of the following three criteria: (1) The signal-to-noise ratio of the faintest resolved hyperfine component is at least 3 (there are three total resolved components for HCN and four for C₂H). (2) The relative errors for all four parameters for this pixel (${\hat{A}}^{* (50)}$, ${\hat{\tau }}^{* (50)}$, ${\hat{{\rm{\Delta }}V}}_{\mathrm{line}}^{(50)}$, and ${\hat{V}}_{\mathrm{LSR}}^{(50)}$) are ≤0.20. Here the relative error is the maximum over [∣P₁₆ − P₅₀∣/P₅₀, ∣P₈₄ − P₅₀∣/P₅₀], where P refers to each parameter. (3) f₅₀ ≤ 0.20.

We then obtained a weighted average on the valid pixels within each annulus to produce the profiles of Figure 9. The weight per pixel and per measurement is (1/σ²), where $\sigma =\sqrt{{({\sigma }_{\mathrm{pix},Q}^{+})}^{2}+{({\sigma }_{\mathrm{pix},Q}^{-})}^{2}}$. The upper and lower shaded error ranges for the measurements in Figure 9 are the rms across all ${\sigma }_{\mathrm{pix},Q}^{+}$ and ${\sigma }_{\mathrm{pix},Q}^{-}$, respectively, for each annulus.

We note that the pixel-by-pixel hyperfine fits and thus the weighted radial profiles have significant intrinsic uncertainties. First, the C₂H and HCN pixel-by-pixel spectra toward J1100-7619 have relatively small line widths (∼0.2–0.5 km s⁻¹) because of the disk's low inclination angle. This fact, coupled with the relatively coarse velocity resolution of the data (0.14 km s⁻¹), means that the hyperfine procedure has a limited number of data points to fit for each hyperfine component. Second, the pixel-by-pixel fitting procedure does not account for uncertainties due to correlations between neighboring pixels. Third, the error in the fits is intrinsically large. This is particularly true for HCN, where the procedure is dependent on the fainter satellite hyperfine components as an "anchor" for the overall fit. Finally, the underlying column density equation (Equations (B12) and (C1)) makes Gaussian assumptions about the molecular line shape that break down as the optical depth increases. At the optical depths we encounter for C₂H and HCN, uncertainties intrinsic to the model may become as large as ∼40%. We thus stress that Figure 9 should be interpreted as rough constraints, rather than precise measurements, of the values, as should the column density estimates determined in Section 4.5.

Assuming an excitation temperature T_ex and total optical depth τ for a hyperfine line, we can estimate the column density ${\overline{N}}_{\mathrm{tot}}$, which is averaged over the disk region, using the total emission across all hyperfine components within that disk region. We use a column density equation adapted from Mangum & Shirley (2015), which presents ${\overline{N}}_{\mathrm{tot}}$ as measured from a specific ith hyperfine emission component:

$$\begin{eqnarray}\begin{array}{rcl}{\bar{N}}_{\mathrm{tot}} & = & \left(\frac{3{hQ}\{{T}_{\mathrm{ex}}\}}{8{\pi }^{3}({S}_{i}{\mu }^{2})}\right)\left(\frac{\exp ({E}_{{\rm{u}}}/{T}_{\mathrm{ex}})}{\exp (h{\nu }_{i}/({k}_{{\rm{B}}}{T}_{\mathrm{ex}}))-1}\right)\\ & & \times \left(\frac{{\displaystyle \int }_{k}{H}_{i}{dV}}{{{ \mathcal T }}_{\mathrm{ex}}-{{ \mathcal B }}_{\mathrm{cont}}}\right)\left(\frac{\tau {R}_{i}}{(1-\exp (-\tau {R}_{i}))}\right),\end{array}\end{eqnarray} \tag{ C1 }$$

where h is the Planck constant, k_B is the Boltzmann constant, and Q{T_ex} is the value of the partition function at temperature T_ex (taken from the CDMS database; Endres et al. 2016). ν _i , S _i μ², and E_u are the frequency, absolute line intensity, and upper-level energy, respectively, of the ith hyperfine component. R_i is the line intensity of the ith component relative to the summed line intensities of all other J-level transitions (the equation and an example calculation are in Appendix B). ${{ \mathcal T }}_{\mathrm{ex}}$ and ${{ \mathcal B }}_{\mathrm{cont}}$ are the Planck excitation temperature and continuum brightness temperature, respectively, calculated using the conversion procedures in Appendix B. Finally, ∫_k H _i dV is the integral across k velocity channels of only the ith hyperfine component of emission in brightness temperature units. This quantity is converted from flux density units using the conversion procedure in Appendix B. Note that the column density equation in Mangum & Shirley (2015) explicitly includes a factor of the upper degeneracy level g_u. For our calculations, that factor of degeneracy is already included implicitly within the CDMS values for S _i μ².

Since we are assuming an excitation temperature T_ex and total optical depth τ, Equation (C1) reduces to the following:

$$\begin{eqnarray}&&\longrightarrow {\bar{N}}_{\mathrm{tot}}={D}_{i}\{{T}_{\mathrm{ex}},\tau \}\times {\int }_{k}{H}_{i}{dV},\end{eqnarray} \tag{ C2 }$$

where D _i {T_ex, τ} is a term dependent on T_ex, τ, and the molecular line characteristics of only the ith hyperfine component.

When the hyperfine components are blended, the quantity (∫_k H _i dV) is not known. We only know the total flux density integrated across all hyperfine components, which we convert to brightness temperature units (∫_k H_tot dV) using the brightness temperature conversion in Appendix B. Even so, we have enough information to set up a system of (L+1) equations with (L+1) unknowns:

$$\begin{eqnarray}\begin{array}{rcl}{\int }_{k}{H}_{\mathrm{tot}}{dV} & = & \sum _{i=1}^{L}{\int }_{k}{H}_{i}{dV},\\ {\bar{N}}_{\mathrm{tot}} & = & {D}_{1}\{{T}_{\mathrm{ex}},\tau \}\times {\int }_{k}{H}_{1}{dV},\\ {\bar{N}}_{\mathrm{tot}} & = & {D}_{2}\{{T}_{\mathrm{ex}},\tau \}\times {\int }_{k}{H}_{2}{dV},\\ & & \vdots \\ {\bar{N}}_{\mathrm{tot}} & = & {D}_{L}\{{T}_{\mathrm{ex}},\tau \}\times {\int }_{k}{H}_{L}{dV},\end{array}\end{eqnarray} \tag{ C3 }$$

where L is the total number of hyperfine components for this hyperfine line (i.e., L = 2 for C₂H (N = 3–2, J = 7/2–5/2), L = 2 for C₂H (N = 3–2, J = 5/2–3/2), and L = 6 for HCN 3–2). Note that while we combine the C₂H (N = 3–2, J = 7/2–5/2) and C₂H (N = 3–2, J = 5/2–3/2) hyperfine spectra in Section 3.2 and Appendix B, we treat the lines separately for this general analysis and use the C₂H (N = 3–2, J = 5/2–3/2) line to check consistency with the C₂H (N = 3–2, J = 7/2–5/2) results.

The unknowns of this system of equations are (∫_k H _i dV) for each ith hyperfine component and ${\overline{N}}_{\mathrm{tot}}$. The latter equations can be subtracted from each other to eliminate ${\overline{N}}_{\mathrm{tot}};$ e.g., (D _i {T_ex, τ} × ∫_k H _i dV) − (D _j {T_ex, τ} × ∫_k H _j dV) = 0 for hyperfine components i and j. These equations can then be written in the form of Mx = b, where M is the matrix (dimensions L × L) of known coefficients:

$$\begin{eqnarray*}{\boldsymbol{M}}=\left[\begin{array}{ccccccc}1 & 1 & 1 & 1 & ... & 1 & 1\\ {D}_{1} & -{D}_{2} & 0 & 0 & ... & 0 & 0\\ 0 & {D}_{2} & -{D}_{3} & 0 & ... & 0 & 0\\ 0 & 0 & {D}_{3} & -{D}_{4} & ... & 0 & 0\\ ... & ... & ... & ... & \ddots & \vdots & \vdots \\ 0 & 0 & 0 & ... & 0 & {D}_{L-1} & -{D}_{L}\end{array}\right],\end{eqnarray*}$$

(where we have temporarily dropped the explicit {T_ex, τ} notation simply to reduce the width of the matrix text). x is the vector (dimensions L × 1) of unknown flux components:

$$\begin{eqnarray*}&&{\boldsymbol{x}}=\left[\begin{array}{c}{\int }_{k}{H}_{1}{dV}\\ {\int }_{k}{H}_{2}{dV}\\ \vdots \\ {\int }_{k}{H}_{L}{dV}\end{array}\right],\end{eqnarray*}$$

and b is the vector (dimensions L × 1) of equation solutions (which are mostly 0):

$$\begin{eqnarray*}&&{\boldsymbol{b}}=\left[\begin{array}{c}{\int }_{k}{H}_{\mathrm{tot}}{dV}\\ 0\\ \vdots \\ 0\end{array}\right].\end{eqnarray*}$$

We calculate M and b, and then we use standard linear algebra practices to solve for x (e.g., numpy's linalg.solve function in Python). We then plug one of the solved hyperfine flux components from x (e.g., the main hyperfine flux component) back into Equation (C2) to recover ${\overline{N}}_{\mathrm{tot}}$.

We now describe the literature samples used to determine the 16th, median, and 84th percentiles of C₂H/HCN and DCN/HCN (Huang et al. 2017; Bergner et al. 2020) in Figure 10.

There were a total of seven disks with ALMA detections of both C₂H 3–2 and HCN 3–2 emission presented in Bergner et al. (2020). When estimating the disk-averaged column densities, Bergner et al. (2020) assumed optical depths for C₂H and HCN of ∼1.1 and ∼5.8, respectively, and an excitation temperature of 30 K. The 16th, median, and 84th percentiles of these estimates (2.4, 2.6, 4.7) are used in the C₂H/HCN panel of Figure 10.

There were a total of five disks with ALMA observations of the optically thin lines DCN 3–2 and H¹³CN (J = 3–2) and one disk with ALMA observations of DCN 3–2 and HCN 3–2 (assumed to be optically thin within the emitting area) presented in Huang et al. (2017). These observations were used to estimate the DCN/HCN abundance ratios for all six disks as proxies for the overall D/H ratios. There were also DCN/HCN column density ratios for three additional disks presented in Bergner et al. (2020). Huang et al. (2017) considered excitation temperatures of 15 and 75 K, while Bergner et al. (2020) assumed 30 K. Here we use the 16th, median, and 84th percentiles across the combined 15 and 30 K estimates (0.011, 0.018, 0.042) in the DCN/HCN panel of Figure 10. The column density ratios are insensitive to temperature changes at these excitation temperatures (Huang et al. 2017; Bergner et al. 2020).