# Massive Star Formation in the Tarantula Nebula

Omnarayani Nayak, Alex Green, Alec Hirschauer, Rémy Indebetouw, Margaret Meixner, Tony Wong, Mélanie Chevance, Guido de Marchi, Vianney Lebouteiller, Min-Young Lee, et al.

► To cite this version:

Omnarayani Nayak, Alex Green, Alec Hirschauer, Rémy Indebetouw, Margaret Meixner, et al. Massive Star Formation in the Tarantula Nebula. The Astrophysical Journal, 2023, 944 (1), pp.26. 10.3847/1538-4357/acac8b. hal-04301596

HAL Id: hal-04301596

https://hal.science/hal-04301596v1

Submitted on 24 Nov 2023

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

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

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

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Massive Star Formation in the Tarantula Nebula

Omnarayani Nayak<fnref n="1" />, Alex Green<fnref n="2" />, Alec S. Hirschauer<fnref n="1" /> Rémy Indebetouw<fnref n="3" />,<fnref n="4" />, Margaret Meixner<fnref n="5" />, Tony Wong<fnref n="2" />, Mélanie Chevance<fnref n="6" />,<fnref n="7" />, Guido De Marchi<fnref n="8" />, Vianney Lebouteiller<fnref n="9" />, Min-Young Lee<fnref n="10" />, Leslie W. Looney<fnref n="2" />, Suzanne C. Madden<fnref n="11" />, Julia Roman-Duval<fnref n="1" />, Yasuo Fukui<fnref n="12" />, Alvaro Hacar<fnref n="13" />, K. E. Jameson<fnref n="14" />, Venu Kalari<fnref n="15" />, Luuk Oudshoorn<fnref n="16" />, Mónica Rubio<fnref n="17" />, and Elena Sabbi<fnref n="1" />

Received 2022 May 16; revised 2022 December 15; accepted 2022 December 15; published 2023 February 9

## Abstract

In this work, we present 299 candidate young stellar objects (YSOs) in 30 Doradus discovered using Spitzer and Herschel point-source catalogs, 276 of which are new. We study the parental giant molecular clouds in which these YSO candidates form using recently published Atacama Large Millimeter/submillimeter Array (ALMA) Cycle 7 observations of ¹²CO and ¹³CO. The threshold for star formation in 30 Doradus inferred by the LTE-based mass surface density is 178 M☉ pc⁻<fnref n="2" />, 40% higher than the threshold for star formation in the Milky Way. This increase in star formation threshold in comparison to the Milky Way and increase in line width seen in clumps 11 pc away in comparison to clumps 45 pc away from the R136 super star cluster could be due to injected turbulent energy, increase in interstellar medium pressure, and/or local magnetic field strength. Of the 299 YSO candidates in this work, 62% are not associated with ¹²CO molecular gas. This large fraction can be explained by the fact that 75%-97% of the H₂ gas is not traced by CO. We fit a Kroupa initial mass function to the YSO candidates and find that the total integrated stellar mass is 18,000 M☉ and that the region has a star formation rate (SFR) of 0.18 M☉ yr⁻<fnref n="1" />. The initial mass function determined here applies to the four 150'' × 150'' (37.5 pc × 37.5 pc) subfields and one 150'' × 75'' (37.5 pc × 18.8 pc) subfield observed with ALMA. The SFR in 30 Doradus has increased in the past few million years.

Unified Astronomy Thesaurus concepts: Star formation (1569); Interstellar medium (847); Young stellar objects (1834)

Supporting material: machine-readable tables

## 1. Introduction

Massive stars play an important role in their host galaxy's properties because stellar feedback regulates star formation (Shetty & Ostriker 2008; Fall et al. 2010; Krumholz 2014; Skinner & Ostriker 2015; Rathjen et al. 2021). The UV radiation from these stars creates H II regions that expand into the surrounding interstellar medium (ISM) owing to warm ionized gas pressure (Lopez et al. 2014; Barnes et al. 2020; Olivier et al. 2021). Stellar winds can sweep up surrounding molecular gas (Pabst et al. 2019, 2020; Lancaster et al. 2021). Supernova explosions are predicted to shred nearby molecular

clouds (Walch et al. 2015; Koo et al. 2020). Bipolar jets and outflows from massive young stellar objects (YSOs) inject energy into the ISM (van Dishoeck & Blake 1998; Matzner & McKee 2000; Bally 2016). Lopez et al. (2011) find several feedback mechanisms in 30 Doradus, a giant H II region located in the Large Magellanic Cloud (LMC), including pressure from warm ionized gas, shock-heated gas from stellar winds, radiation from stars, and dust-processed radiation. We aim to quantify the effects of feedback on the next generation of star formation through observations of YSOs and their nascent molecular clouds.

As the most luminous star-forming region in the Local Group, 30 Doradus provides a unique opportunity to study massive star formation and how it drives and responds to stellar feedback. At the heart of 30 Doradus lies R136, a super star cluster (SSC) with extraordinarily high stellar densities of 1.5 × 10⁴–10⁷ M☉ pc<fnref n="3" /> (Selman & Melnick 2013), containing the most massive stars known (Crowther et al. 2016).

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

<footnote n="1">¹Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA; inayak@stsci.edu</footnote>

<footnote n="2">²Astronomy Department, University of Illinois, Urbana, IL 61801, USA</footnote>

<footnote n="3">³Department of Astronomy, University of Virginia, P.O. Box 3818, Charlottesville, VA 22903, USA</footnote>

<footnote n="4">⁴National Radio Astronomy Observatory, 520 Edgemont Road, Charlottesville, VA 22903, USA</footnote>

<footnote n="5">⁵SOFIA-USRA, NASA Ames Research Center, MS 232-12, Moffett Field, CA 94035, USA</footnote>

<footnote n="6">⁶Zentrum für Astronomie der Universität Heidelberg, Institut für Theoretische Astrophysik, Albert-Ueberle-Str. 2, 69120 Heidelberg, Germany</footnote>

<footnote n="7">⁷Cosmic Origins Of Life (COOL) Research DAO¹⁸</footnote>

<footnote n="8">⁸European Space Research and Technology Centre, Keplerlaan 1, 2200 AG Noordwijk, The Netherlands</footnote>

<footnote n="9">⁹AIM, CEA, CNRS, Université Paris-Saclay, Université Paris Diderot, Sorbonne Paris Cité, F-91191 Gif-sur-Yvette, France</footnote>

<footnote n="10">¹⁰Korea Astronomy and Space Science Institute, 776 Daedeokdae-ro, 34055 Daejeon, Republic of Korea</footnote>

<footnote n="11">¹¹Département d'Astrophysique AIM/CEA Saclay, Orme des Merisiers, F-91191 Gif-sur-Yvette, France</footnote>

<footnote n="12">¹²Department of Physics, Nagoya University, Chikusa-ku, Nagoya 464-8602, Japan</footnote>

<footnote n="13">¹³Department of Astrophysics, University of Vienna, Türkenschanzstrasse 17, A-1180 Vienna, Austria</footnote>

<footnote n="14">¹⁴CSIRO, Space and Astronomy, PO Box 1130, Bentley, WA 6102, Australia</footnote>

<footnote n="15">¹⁵Gemini Observatory, NSF NOIRLab, Casilla 603, La Serena, Chile</footnote>

<footnote n="16">¹⁶Leiden Observatory, Leiden University, Niels Bohrweg 2, 2333CA Leiden, The Netherlands</footnote>

<footnote n="17">¹⁷Departamento de Astronomía, Universidad de Chile, Casilla 36-D, Santiago, Chile</footnote>

<footnote n="18">¹⁸ coolresearch.io</footnote>
H II regions and photon-dominated regions (PDRs) have been studied in detail, revealing that more than 75% of the molecular gas is CO-dark (Chevance et al. 2016, 2020), radiation pressure does not currently dominate (Pellegrini et al. 2011), and a hard radiation field exists with $T_{\text{rad}}$ between 30,000 and 85,000 K (Indebetouw et al. 2009). The metallicity of the LMC (0.5 $Z_{\odot}$) and the extreme mode of star formation in 30 Doradus provide a local laboratory to understand the universe's peak epoch of star formation at $z \sim 1.5-2.0$ (Madau & Dickinson 2014).

The extreme nature of the R136 SSC has led to multiple studies that try to understand what specifically has caused the extreme star formation episodes. Observational evidence shows four major star formation activities in 30 Doradus: an old starburst event that occurred 25 Myr ago and is currently known as Hodge 301 (Grebel & Chu 2000), an increase in star formation activity 12 Myr ago (De Marchi et al. 2011), a more recent event a few million years old as indicated from spectroscopy of O stars (Hunter et al. 1995), and another ongoing star formation episode occurring in the arc of molecular gas surrounding 30 Doradus (Rubio et al. 1992; Walborn et al. 1999). Powerful evidence of the ongoing star formation activity extending beyond the central cluster is provided by the Hubble Tarantula Treasury Program (HTTP) survey, which discovered more than 16,700 pre-main-sequence (PMS) stars (Sabbi et al. 2013; Walborn et al. 2013; Ksoll et al. 2018; Melnick et al. 2021).

Giant molecular clouds (GMCs), which are host to intense star formation activity, have been studied extensively in the LMC (Fukui et al. 2009; Minamidani et al. 2011; Wong et al. 2011, 2017). Theoretical predictions and observations show that GMCs are hierarchical in nature, with small and dense structures embedded in a larger, lower-density substrate. The process by which larger clumps fragment and form dense cores is dependent on the radiative and mechanical feedback of the environment. Radiative and mechanical energy from past star formation, strong stellar winds, and galaxy-wide interactions are all prevalent in 30 Doradus. Similar shock properties in 30 Doradus and regions 600 pc away suggest kiloparsec-scale energy being injected into the LMC owing to interactions from the Milky Way and the Small Magellanic Cloud (SMC; Lee et al. 2016, 2019). The intense radiation from R136 and the previous star formation episodes that took place provide feedback on local scales in the 30 Doradus region (De Marchi et al. 2011).

A typical GMC is $10^4-10^5 M_{\odot}$, and dense clouds within the GMC are $10^2-10^3 M_{\odot}$. Subparsec-size cores are even lower in mass, but they are extremely dense regions, $10^5-10^6 M_{\odot} \text{ pc}^{-2}$ (McKee & Ostriker 2007; Zinnecker & Yorke 2007). We aim to study the YSO candidates forming in these dense cores and the molecular gas that fuels the star formation process (Fukui et al. 2015; Nayak et al. 2016, 2018). YSOs forming within clouds and filaments in the LMC have been the subject of many studies, each using different strategies to identify the YSOs and separate likely YSO candidates from other dusty sources such as asymptotic giant branch (AGB) stars and background galaxies. Gruendl & Chu (2009) used a single cut in a mid-infrared color-magnitude diagram (CMD) and then visually inspected the source to classify it as a YSO candidate. Whitney et al. (2008) and Carlson et al. (2012) use a series of complex color cuts to maximize the number of YSO candidates and minimize background galaxy contaminants, as well as fit spectral energy distribution (SED) models to the photometry.

Seale et al. (2014) remove extended-looking sources from their catalog (most likely to be galaxies) and then use the 24 µm emission as a tracer of star formation to identify YSO candidates. Seventeen sources in this work have been identified as YSOs with Spitzer IRS spectra (Seale et al. 2009; Jones et al. 2017). However, these initial surveys missed the lower-mass and lower-luminosity YSO candidates owing to the stringent color cuts and different YSO candidate selection criteria.

In this work we use SED models to identify and create a list of 299 YSO candidates in 30 Doradus, in a manner similar to previous studies (Nayak et al. 2016, 2018). We compare these YSOs to the molecular gas properties derived from Atacama Large Millimeter/submillimeter Array (ALMA) observations (Wong et al. 2022). In Section 2, we describe the observations. Section 3 describes the YSO selection method and the dendrogram method in identifying clumps. Section 4 goes into further detail about the point sources in the region: those that are YSO candidates and those that do not fit our criteria to be considered a YSO. In Section 5, we discuss the properties of the population of molecular clouds and of the 299 YSO candidates, as well as their association with each other. We present our conclusions in Section 6.

## 2. Observations and Methods

### 2.1. Spitzer SAGE and Herschel HERITAGE Surveys

The Spitzer Survey of the Agents of Galaxy Evolution (SAGE; Meixner et al. 2006) and Herschel HERSchel Inventory of the Agents of Galaxy Evolution (HERITAGE; Meixner et al. 2013) imaged the LMC with multiwavelength photometry at 3.6, 4.5, 5.8, 8.0, 24, 70, 100, 160, 250, 350, and 500 µm (Fazio & Eisenhardt 2004; Rieke et al. 2004; Roellig et al. 2004). One of the main goals of SAGE was to detect infrared point sources in regions with dust, photodissociation, and ionized hydrogen (Meixner et al. 2006). The HERITAGE survey added 250, 350, and 500 µm wavelength photometry to better trace the dust that is inherent to the most embedded and young protostars (Meixner et al. 2013).

The J-, H-, and K-band photometry is from the Magellanic Clouds Point Source Catalog taken with the InfraRed Survey Facility (IRSF; Kato et al. 2007). The IRSF photometry was converted to the Two Micron All Sky Survey (2MASS) for compatibility with the SED fitter (Kato et al. 2007). Spitzer InfraRed Array Camera (IRAC) and IRSF point sources were cross-matched using a nearest neighbor criterion within 0.′′1. The nearest IRAC point source within 1.′′5 to a Multiband Imaging Photometer (MIPS) 24 µm point source is adopted to be the same source. For each of the wavelengths of the Herschel Photodetector Array Camera and Spectrometer (PACS) and Spectral and Photometric Imaging Receiver (SPIRE), point sources were extracted such that they are 5σ above the background, and then sources were matched between bands using a simple nearest match criterion. Seale et al. (2014) matched Herschel PACS and SPIRE points to the Spitzer 24 µm point-source catalog and then bootstrapped the IRAC point sources through the 24 µm matches. Herschel point sources do not directly match Spitzer point sources owing to the difference in beam sizes. The FWHM beam width is 1.′′−2.′′ for IRAC, 6.′′−40.′′ for MIPS, 5.′′−13.′′ for PACS, and 17.′′−35.′′ for SPIRE (Meixner et al. 2006, 2013).

We use the published band-merged Spitzer SAGE point-source catalog to identify point sources within the ALMA
footprint of 30 Doradus. This catalog has been merged between IRSF, IRAC, PACS, and SPIRE bands as described above. There are 488 point sources from the Spitzer point-source catalog within the ALMA Cycle 7 observations of 30 Doradus. We aim to classify every source in the band-merged Spitzer point-source catalog as a YSO candidate or not. We perform new aperture photometry to extract the Spitzer and Herschel photometry for all the point sources in the ALMA footprint of the region to realize consistent flux calculations. The aperture size is circular and centers on the Spitzer/IRAC point-source coordinates. We set the aperture radius equal to the beam size in each band as described above. Flux uncertainties were determined from an annulus that is twice the radius of the aperture. If the flux of the point source is more than three times the rms noise in the annulus, we consider this to be a detection. This is the same aperture photometry algorithm used by Carlson et al. (2012). We check that our flux is within 10% of published photometric fluxes when there are common sources and bands (Whitney et al. 2008; Gruendl & Chu 2009; Carlson et al. 2012). Table 1 lists the *JHK* photometry from IRSF and the IRAC photometry from Spitzer for all 488 point sources in the 30 Doradus ALMA footprint taken during Cycle 7. Table 2 lists the Spitzer MIPS, Hershel PACS, and Hershel SPIRE photometry of the 488 point sources. The fluxes listed in Tables 1 and 2 are what we use in conjunction with the Robitaille (2017) YSO SED fitter and Castelli & Kurucz (2004) stellar photosphere SED fitter to classify the point sources.

## 2.2. ALMA $^{12}$CO and $^{13}$CO Observations

We use the ALMA Cycle 7 $^{12}$CO (2 − 1) and $^{13}$CO (2 − 1) observations (2019.1.00843.S) taken by Wong et al. (2022). The ALMA footprint of 30 Doradus covers the R136 SSC and its surrounding molecular gas. Four of these subfields map a region 150" × 150" (37.5 pc × 37.5 pc) in size, and the fifth subfield was 150" × 75" (37.5 pc × 18.8 pc) in size. The correlator was set up to cover the $^{12}$CO (2 − 1) and $^{13}$CO (2 − 1) lines with a spectral resolution of ~0.1 km s⁻¹ and spatial resolution of 1.75" (0.4 pc). ACA 7 m and total power (TP) cubes were combined using the feather task, and then this process was repeated to combine the 7 m + TP with the ALMA 12 m observations.

Wong et al. (2022) used dendrograms to categorize the filaments and cores seen in the ALMA $^{12}$CO and $^{13}$CO observations in order to understand the correlation between the molecular gas and the formation of the YSO candidates. The dendrogram tree structure is a novel approach of identifying isosurfaces at various emission levels nested within other isosurfaces in a position–position–velocity data cube (Rosolowsky et al. 2008). The trunk on the dendrogram structure is the largest isosurface and represents the lowest emission level, which is determined from the user input. Trunks are connected to the branches of a dendrogram, which represent the nested isosurfaces with increasing emission levels. The leaves are the local emission maxima that meet the minimum pixel criteria input from the user. Leaves represent high-density cores that are traced by molecular gas such as HCN, HCO$^{+}$, and SO (Nayak et al. 2019). The assumptions for calculating the clump parameters, as well as the full catalog of clumps, are given in Wong et al. (2022). We use the 250 m s⁻¹ feathered mosaic data cubes from Wong et al. (2022) in this work, which has an rms noise per channel of ~0.26 K (35 mJy per 1.75" beam). The dendrogram code outputs properties such as radius,

area, and line width of individual structures, which are then used to calculate column density, virial mass, LTE mass, etc. (Rosolowsky et al. 2008).

Dendrograms have previously been used to study $^{12}$CO and $^{13}$CO clumps in very active star-forming regions like the central molecular zone (CMZ), N159 (in the LMC), and 30 Doradus. These studies show that high-mass star-forming regions have different properties than regions of less intense star formation (e.g., Kauffmann et al. 2010; Shetty et al. 2012; Nayak et al. 2016, 2018). There is an excess in line width measurements for a given clump size in very active star-forming regions, shifting the size–line width relation upward almost 1 dex (Shetty et al. 2012; Nayak et al. 2016). This shift in excess line widths is further discussed in Section 5.3.

The clump structures output from the dendrogram code are nested within each other and therefore are not an independent set. We use the Spectral Clustering for Interstellar Molecular Emission Segmentation (SCIMES) dendrogram method to obtain a set of independent structures (Colombo et al. 2015). The main idea of the dendrogram + SCIMES algorithm is to map the data as a graph onto a vector space, and if cluster properties exist, then these properties will be represented by different vectors (Colombo et al. 2015). Leaves that do not form isolated structures are removed in an iterative process until a complex dendrogram tree structure is only a few branches. Therefore, the SCIMES clusters are usually slightly smaller than the dendrogram trunks, but still large-scale isoemission structures that are independent from one another. In this work, we use the dendrogram method and the dendrogram + SCIMES method to analyze the $^{12}$CO and $^{13}$CO ALMA data.

## 3. YSO Selection Criteria

We investigate every point source from the published SAGE point-source catalog within the ALMA footprint of 30 Doradus to determine whether these 488 sources are likely YSO candidates. Previous galaxy-wide studies by Whitney et al. (2008), Gruendl & Chu (2009), and Carlson et al. (2012) implemented stringent color–color cuts, and it has been shown by Nayak et al. (2016, 2018) that massive YSOs were missed by such studies because of these cuts or the requirement of detection with the MIPS 24 µm band. We use the fluxes of the 488 point sources listed in Tables 1 and 2 and fit Robitaille (2017) models named “spbhmi”—i.e., models with ambient media, a power-law envelope, and bipolar cavities. Table 3 lists the point sources, the R.A. and decl. of each source, the reduced χ² fits to the Robitaille (2017) “spbhmi” models and Castelli & Kurucz (2004) photosphere models, and whether the source was previously studied or classified. The dust properties in the “spbhmi” models use a grain size distribution from Weingartner & Draine (2001) with the dust model by Draine (2003a, 2003b). Robitaille (2017) does not take into account PAH emission in the dust models.

The requirements to qualify as a likely YSO candidate in this work are as follows: (i) it has to be a good fit to Robitaille (2017) models (reduced χ²$_{robitaille}$ < 10), (ii) the source has to look point-like and not extended in the Hubble Space Telescope (HST) F555W band, and (iii) the source has previously not been identified as an OB star with Very Large Telescope (VLT) observations. The YSO models from Robitaille (2017) include a wide range of parameters (i.e., disk, envelope, ambient media, viewing angle) in order to
<table><caption>Table 1
2MASS and IRAC Point-source Photometry</caption><thead><tr><th>Name</th><th>R.A.</th><th>Decl.</th><th>J</th><th>eJ</th><th>H</th><th>eH</th><th>K</th><th>eK</th><th>I1</th><th>eI1</th><th>I2</th><th>eI2</th><th>I3</th><th>eI3</th><th>I4</th><th>eI4</th></tr></thead><tbody><tr><td>Y1</td><td>84.713666</td><td>-69.042588</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0505</td><td>0.0081</td><td>0.0000</td><td>0.0000</td><td>0.3870</td><td>0.0774</td><td>0.4660</td><td>0.2005</td></tr><tr><td>Y2</td><td>84.700704</td><td>-69.069917</td><td>1.3900</td><td>0.1390</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>12.5000</td><td>1.2610</td><td>17.2000</td><td>1.7320</td><td>39.9000</td><td>4.1120</td><td>98.1000</td><td>10.1670</td></tr><tr><td>Y3</td><td>84.641899</td><td>-69.081180</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>13.1000</td><td>1.3270</td><td>30.1000</td><td>3.0230</td><td>56.6000</td><td>5.7180</td><td>114.0000</td><td>11.6360</td></tr><tr><td>Y4</td><td>84.657806</td><td>-69.083230</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.8660</td><td>0.1369</td><td>0.4740</td><td>0.1004</td><td>4.8300</td><td>0.6271</td><td>0.0000</td><td>0.0000</td></tr><tr><td>Y5</td><td>84.631746</td><td>-69.037392</td><td>2.3900</td><td>0.2390</td><td>3.6200</td><td>0.3620</td><td>6.3700</td><td>0.6370</td><td>10.6000</td><td>1.0670</td><td>15.8000</td><td>1.5840</td><td>44.4000</td><td>4.4940</td><td>136.0000</td><td>13.7320</td></tr><tr><td>Y6</td><td>84.719300</td><td>-69.077300</td><td>1.4999</td><td>0.1490</td><td>1.3889</td><td>0.1380</td><td>2.8410</td><td>0.2840</td><td>5.1630</td><td>0.2350</td><td>5.2275</td><td>0.2452</td><td>17.4400</td><td>1.4786</td><td>38.1410</td><td>3.8483</td></tr><tr><td>Y7</td><td>84.684938</td><td>-69.053818</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.1190</td><td>0.0198</td><td>0.0848</td><td>0.0132</td><td>0.6270</td><td>0.1266</td><td>0.0000</td><td>0.0000</td></tr><tr><td>Y8</td><td>84.628583</td><td>-69.115335</td><td>1.8800</td><td>0.1110</td><td>0.0000</td><td>0.0000</td><td>2.8700</td><td>0.1860</td><td>1.1200</td><td>0.4220</td><td>1.9600</td><td>0.5960</td><td>1.9400</td><td>0.1890</td><td>0.0000</td><td>0.0000</td></tr><tr><td>Y9</td><td>84.702173</td><td>-69.067105</td><td>0.0767</td><td>0.0283</td><td>0.5730</td><td>0.0430</td><td>1.0900</td><td>0.0891</td><td>1.1500</td><td>0.1428</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td></tr><tr><td>Y10</td><td>84.799534</td><td>-69.092910</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.0000</td><td>0.1060</td><td>0.0183</td><td>0.1830</td><td>0.0405</td><td>0.0933</td><td>0.0464</td><td>0.3800</td><td>0.1259</td></tr></tbody></table>

**Note.** Column (1): name of source. Column (2): R.A. of source. Column (3): decl. of source. Column (4): *J* flux. Column (5): error in *J* flux. Column (6): *H* flux. Column (7): error in *H* flux. Column (8): *K* flux. Column (9): error in *K* flux. Column (10): IRAC 3.6 µm flux. Column (11): error in IRAC 3.6 µm flux. Column (12): IRAC 4.5 µm flux. Column (13): error in IRAC 4.5 µm flux. Column (14): IRAC 5.8 µm flux. Column (15): error in IRAC 5.8 µm flux. Column (16): IRAC 8.0 µm flux. Column (17): error in IRAC 8.0 µm flux. All fluxes and errors are in units of mJy. A value of 0 indicates a nondetection. The full table is available in machine-readable format.

(This table is available in its entirety in machine-readable form.)
<table><caption>**Table 2**

Mid-to-far-infrared Point-source Photometry</caption><thead><tr><th>Name</th><th>M1</th><th>eM1</th><th>M2</th><th>eM2</th><th>PACS2</th><th>ePACS2</th><th>PACS3</th><th>ePACS3</th><th>SPIRE1</th><th>eSPIRE1</th><th>SPIRE 2</th><th>eSPIRE2</th><th>SPIRE3</th><th>eSPIRE3</th></tr></thead><tbody><tr><td>Y1</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr><tr><td>Y2</td><td>4530.00</td><td>192.20</td><td>81000.00</td><td>5922.20</td><td>30000.00</td><td>3006.90</td><td>21000.00</td><td>2156.60</td><td>23900.00</td><td>2013.70</td><td>18900.00</td><td>1542.50</td><td>16800.00</td><td>1351.80</td></tr><tr><td>Y3</td><td>3790.00</td><td>183.28</td><td>18500.00</td><td>1510.70</td><td>12200.00</td><td>1228.40</td><td>9720.00</td><td>982.76</td><td>7190.00</td><td>637.99</td><td>3270.00</td><td>654.55</td><td>3280.00</td><td>425.53</td></tr><tr><td>Y4</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>5500.00</td><td>530.86</td></tr><tr><td>Y5</td><td>5460.00</td><td>218.99</td><td>59600.00</td><td>4189.20</td><td>29700.00</td><td>2971.10</td><td>20700.00</td><td>2073.20</td><td>18200.00</td><td>1467.50</td><td>10000.00</td><td>817.04</td><td>7180.00</td><td>589.15</td></tr><tr><td>Y6</td><td>6419.00</td><td>108.50</td><td>77519.00</td><td>1471.50</td><td>13521.00</td><td>537.22</td><td>21343.00</td><td>589.32</td><td>28562.00</td><td>510.39</td><td>17834.00</td><td>331.58</td><td>15762.00</td><td>152.22</td></tr><tr><td>Y7</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>1790.00</td><td>351.49</td></tr><tr><td>Y8</td><td>3510.00</td><td>73.20</td><td>81100.00</td><td>1130.00</td><td>4460.00</td><td>174.00</td><td>6650.00</td><td>274.00</td><td>13600.00</td><td>507.00</td><td>8260.00</td><td>326.00</td><td>8270.00</td><td>146.00</td></tr><tr><td>Y9</td><td>3350.00</td><td>345.25</td><td>50900.00</td><td>5389.00</td><td>18600.00</td><td>1875.70</td><td>19000.00</td><td>1917.30</td><td>21900.00</td><td>1904.20</td><td>15100.00</td><td>1287.00</td><td>14900.00</td><td>1201.90</td></tr><tr><td>Y10</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td><td>0.00</td></tr></tbody></table>

**Note.** Column (1): name of source. Column (2): MIPS 24 µm flux. Column (3): MIPS 24 µm flux error. Column (4): MIPS 70 µm flux. Column (5): MIPS 24 µm flux error. Column (6): PACS 100 µm flux. Column (7): PACS 100 µm flux error. Column (8): PACS 160 µm flux. Column (9): PACS 160 µm flux error. Column (10): SPIRE 250 µm flux. Column (11): SPIRE 250 µm flux error. Column (12): SPIRE 350 µm flux. Column (13): SPIRE 350 µm flux error. Column (14): SPIRE 500 µm flux. Column (15): SPIRE 500 µm flux error. All fluxes and errors are in units of mJy. A value of 0 indicates a nondetection. The full table is available in machine-readable format.

(This table is available in its entirety in machine-readable form.)
<table><caption>**Table 3**
Point Sources and Reduced $\chi^2$ Fits</caption><thead><tr><th>Name</th><th>R.A.</th><th>Decl.</th><th>Log10 (Reduced χ<sup>2</sup> Robitaille SED)</th><th>Log10 (Reduced χ<sup>2</sup> Kurucz SED)</th><th>Current Classification</th><th>Reason</th><th>References</th><th>Previous Classification</th></tr></thead><tbody><tr><td>Y1</td><td>84.713666</td><td>-69.042588</td><td>-1.65</td><td>1.39</td><td>YSO</td><td>SED fit</td><td></td><td></td></tr><tr><td>Y2</td><td>84.700704</td><td>-69.069917</td><td>1.75</td><td>3.29</td><td>YSO</td><td>Spitzer IRS</td><td>W08, G09, S14, J17</td><td>YSO/H II</td></tr><tr><td>Y3</td><td>84.641899</td><td>-69.081180</td><td>1.18</td><td>3.27</td><td>YSO</td><td>Spitzer IRS</td><td>W08, G09, S14, J17</td><td>YSO/H II</td></tr><tr><td>Y4</td><td>84.657806</td><td>-69.083230</td><td>0.47</td><td>1.60</td><td>YSO</td><td>SED fit</td><td></td><td></td></tr><tr><td>Y5</td><td>84.631746</td><td>-69.037392</td><td>1.96</td><td>1.28</td><td>YSO</td><td>Spitzer IRS</td><td>G09, S14, J17</td><td>YSO/H II</td></tr><tr><td>Y6</td><td>84.719300</td><td>-69.077300</td><td>2.05</td><td>3.98</td><td>YSO</td><td>Spitzer IRS</td><td>G09, S09, J17</td><td>YSO/H II</td></tr><tr><td>Y7</td><td>84.684938</td><td>-69.053818</td><td>0.22</td><td>1.34</td><td>YSO</td><td>SED fit</td><td></td><td></td></tr><tr><td>Y8</td><td>84.628583</td><td>-69.115335</td><td>1.58</td><td>3.55</td><td>other</td><td>SED fit</td><td></td><td></td></tr><tr><td>Y9</td><td>84.702173</td><td>-69.067105</td><td>0.02</td><td>3.03</td><td>YSO</td><td>SED fit</td><td></td><td></td></tr><tr><td>Y10</td><td>84.799534</td><td>-69.092910</td><td>-0.21</td><td>0.78</td><td>YSO</td><td>SED fit</td><td></td><td></td></tr></tbody></table>

**Note.** Column (1): name of the point source we use in this manuscript. Column (2): R.A. of the point source. Column (3): decl. of the point source. Column (4): log of the reduced $\chi^2$ fit to the “sphbmi” Robitaille (2017) model. Column (5): log of the reduced $\chi^2$ fit to the Castelli & Kurucz (2004) stellar photosphere model. Column (6): how we classify these sources in this work. Column (7): reason for classification in this work. Column (8): references if this source was previous identified as a YSO candidate: W08 is Whitney et al. (2008), G09 is Gruendl & Chu (2009), S09 is Seale et al. (2009), S14 is Seale et al. (2014), N16 is Nayak et al. (2016), J17 is Jones et al. (2017), S18 is Schneider et al. (2018). Column (9): previous classification. The full table is available in machine-readable format. Note: reduced $\chi^2$ values less than 1 are possible because occasionally a model fits the observed data well. For example, in the top right and middle left panels of Figure 1, the $\chi^2$ values are 0.343 and 2.415, respectively. There are three fitted points in the top right panel and five fitted points in the middle left panel, making the reduced $\chi^2$ values 0.114 and 0.483, respectively.

(This table is available in its entirety in machine-readable form.)

<table><caption>**Table 4**
Properties of Likely YSOs Ranked from Most Massive to Least Massive</caption><thead><tr><th>Name</th><th>Rad (R<sub>☉</sub>)</th><th>Rad Error (R<sub>☉</sub>)</th><th>Temp (K)</th><th>Temp Error (K)</th><th>Mass (M<sub>☉</sub>)</th><th>Mass Error (M<sub>☉</sub>)</th><th>Lum (L<sub>☉</sub>)</th><th>Lum Error (L<sub>☉</sub>)</th><th><sup>12</sup>CO Association</th><th><sup>13</sup>CO Association</th></tr></thead><tbody><tr><td>Y1</td><td>99.6</td><td>18.9</td><td>1.42e+04</td><td>7.27e+03</td><td>38.5</td><td>5.3</td><td>3.53e+05</td><td>6.69e+04</td><td></td><td></td></tr><tr><td>Y2</td><td>36.1</td><td>0.0</td><td>2.05e+04</td><td>0.00e+00</td><td>32.8</td><td>0.0</td><td>2.02e+05</td><td>0.00e+00</td><td>y</td><td></td></tr><tr><td>Y3</td><td>30.8</td><td>0.0</td><td>1.99e+04</td><td>0.00e+00</td><td>29.0</td><td>0.0</td><td>1.32e+05</td><td>0.00e+00</td><td>y</td><td>y</td></tr><tr><td>Y4</td><td>38.9</td><td>14.0</td><td>1.76e+04</td><td>4.76e+03</td><td>28.9</td><td>9.0</td><td>1.30e+05</td><td>2.32e+05</td><td>y</td><td>y</td></tr><tr><td>Y5</td><td>64.4</td><td>0.0</td><td>1.30e+04</td><td>0.00e+00</td><td>27.3</td><td>0.0</td><td>1.07e+05</td><td>0.00e+00</td><td>y</td><td></td></tr><tr><td>Y6</td><td>81.6</td><td>0.0</td><td>1.12e+04</td><td>0.00e+00</td><td>26.2</td><td>0.0</td><td>9.21e+04</td><td>0.00e+00</td><td>y</td><td>y</td></tr><tr><td>Y7</td><td>18.7</td><td>9.7</td><td>2.28e+04</td><td>5.06e+03</td><td>25.5</td><td>7.5</td><td>8.33e+04</td><td>5.37e+04</td><td></td><td></td></tr><tr><td>Y9</td><td>67.0</td><td>14.3</td><td>1.16e+04</td><td>1.90e+03</td><td>24.4</td><td>1.1</td><td>7.14e+04</td><td>1.31e+04</td><td>y</td><td>y</td></tr><tr><td>Y10</td><td>28.4</td><td>20.3</td><td>1.77e+04</td><td>7.24e+03</td><td>24.3</td><td>8.7</td><td>7.05e+04</td><td>9.39e+04</td><td></td><td></td></tr><tr><td>Y13</td><td>31.3</td><td>16.8</td><td>1.54e+04</td><td>3.02e+03</td><td>21.9</td><td>1.3</td><td>4.90e+04</td><td>6.39e+03</td><td>y</td><td></td></tr></tbody></table>

**Note.** Column (1): YSO candidate name we use in this work. Column (2): best-fit radius of YSO candidate. Column (3): best-fit temperature of YSO candidate. Column (4): mass of the YSO candidate calculated using $L \propto M^{3.5}$. Column (5): luminosity of YSO candidate calculated using the best-fit radius and best-fit temperature with the following formula: $L = 4 \pi R^2 \sigma T^4$. An error of 0 indicates that there was only one YSO model that fit the observed fluxes. Column (6): if the YSO candidate is associated with <sup>12</sup>CO. Column (7): if the YSO candidate is associated with <sup>13</sup>CO. The full table is available in machine-readable format.

(This table is available in its entirety in machine-readable form.)

accurately classify the newly forming stars. There are 299 YSO candidates within the ALMA footprint of 30 Doradus, 5 of which do not have HST coverage, and therefore we cannot confirm them to be point-like in the HST F555W band. As YSOs evolve, they are no longer enshrouded by dust, and more of the optical light from the star itself is observed. Therefore, a more evolved YSO and a star whose optical radiation is dominating can have very similar SEDs. We require an additional criterion for a point source to be a YSO candidate: (iv) the fit to the Robitaille (2017) YSO models is better than the fit to the Castelli & Kurucz (2004) stellar photosphere models. Previously identified YSOs that have spectra but do not meet the above four criteria are included in our final list of YSO candidates because the emission- and absorption-line features, such as ices, $H_2$, and polycyclic aromatic hydrocarbons (PAHs) in the spectra, are characteristic of YSOs. The

physical properties of the 299 likely YSO candidates are listed in Table 4. Six of the YSO candidates with good SED fits are shown in Figure 1, and 18 of the SED fits for point sources that are not YSO candidates are shown in Figures 2–4. For the purpose of our analysis, we use the 299 likely YSO candidates identified in this work. There are 33 point sources that were previously identified as YSO candidates, but only 23 fit the criteria in this work, and therefore these 23 are part of the 299 likely YSO candidates.

As evidenced by the CMDs in Figure 5, the 276 new YSO candidates identified in this work are generally less luminous than those that have been previously identified. Other dusty sources that could lead to false positives include AGB stars, planetary nebulae, and background galaxies. Sewiło et al. (2013) look at the same CMD filters to identify YSOs in the SMC, and they separate them from AGB stars and background
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Example SEDs of likely YSO candidates (i.e., those with reduced $\chi^2_{\text{robitaille}} < 10$). The best-fit SED is shown with the black line, and fits with $\chi^2 < 3$ relative to the best fit are shown with the gray lines. The best-fit line model is used to determine YSO properties such as mass and luminosity.</figcaption></figure>

galaxies. AGB stars are located in the same region as YSOs, evident from Figure 2 of Sewiło et al. (2013), which shows the same four CMDs we plot in Figure 5 in this work. However, the SED of a typical AGB star is different from that of a YSO. Jones et al. (2017) show the spectra of AGB stars in Figures 9–11 on their work, and none of the spectra have the characteristic rise toward the mid- and far-infrared wavelengths seen in the SEDs we show in Figure 1. Even though AGB stars

can indeed occupy a similar space in the CMD, we require the YSO candidates in this work to fit the Robitaille (2017) SED models. Therefore, we do not expect AGB stars masquerading as YSO candidates in our robust photometric list. Similar to AGB stars, planetary nebulae are not expected to be abundant in sites of young star formation. We look at the list of planetary nebulae by Hora et al. (2008) and find that there are none within the ALMA Cycle 7 footprint.
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Example SEDs of point sources unlikely to be YSO candidates with reduced $\chi^2_{\text{robitaille}}$ > 10 and <100.</figcaption></figure>

Background galaxies can be well fit by Robitaille SED models; however, they inhabit a different region of the CMD space than YSOs in this paper. The wide separation in color is optimal for distinguishing between background galaxies and YSOs. The new YSO candidates in this work are in this CMD range, where the [4.5]–[24] color ranges between 8 and 10 and [4.5] mag ranges between 14 and 16, as can be seen in the bottom right panel of Figure 5. Figure 2 from Sewiło et al. (2013) shows the [4.5]–[24] versus [24] CMD. There are no

background galaxies in the region where the [4.5]–[24] color ranges between 8 and 10 and [4.5] mag ranges between 14 and 16. Background galaxies compose up to 12% of the point sources in the LMC (Meixner et al. 2006). No background galaxies are seen by Sewiło et al. (2013) in the specific CMD region where the new YSO candidates in this work are located. All 276 new YSO candidates look point-like in the HST F555W band. Additionally, the extinction at the HST wavelengths is far too high to be able to detect background
<figure><img src="image_5.png" /><figcaption>Figure 3. Example SEDs of point sources unlikely to be YSO candidates with reduced $\chi^2_{\text{robitaille}}$ > 100 and <1000.</figcaption></figure>

galaxies that could masquerade as a YSO a few solar masses in size (i.e., we are not concerned with the lower-luminosity and lower-mass YSO candidates in this work being false-positive YSO candidates). Therefore, it is unlikely that the new YSO candidates we identify in this work are background galaxies.

We compare the CMDs in Figure 5 to the CMD from Carlson et al. (2012) to determine whether the YSO candidates in this work indeed fall in a similar parameter space to other YSO candidates in the LMC. Carlson et al. (2012) look for

YSO candidates in nine different GMCs in the LMC: N11, N44, N51, N105, N113, N120, N144, N160, and N206. They used a series of $\alpha$ and $\beta$ color cuts for their selection criteria. The $\alpha$ cut criteria can be applied to galaxy-wide surveys where low contamination is expected. This cut selects the more massive and more luminous YSO candidates in the region. The $\beta$ cut is applied only to star-forming regions where there is a higher probability that a color-selected source is a YSO. Because of the higher chance of contamination with the $\beta$ cuts,
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Example SEDs of point sources unlikely to be YSO candidates with reduced $\chi^2_{\text{robitaille}}$ > 1000.</figcaption></figure>

Carlson et al. (2012) required that a YSO candidate meet two different $\beta$ cuts if it did not meet the criteria for the $\alpha$ cut. Figure 5 shows the Carlson et al. (2012) $\alpha$ cut with the red line and the $\beta$ cut with the blue line.

Previous YSO candidates identified in the 30 Doradus region by Whitney et al. (2008), Gruendl & Chu (2009), Seale et al. (2009, 2014), and Jones et al. (2017) that have 4.5 and 8.0 µm photometry are shown with green circles in Figure 5. New YSO candidates we identify in this work are shown with orange

squares. The new YSO candidates we identify in this work are fainter (i.e., lower luminosity and less massive) than the previously identified YSO candidates owing to the stringent color–color cuts implemented. For example, Whitney et al. (2008) required that YSO candidates have IRAC [8.0] mag brighter than 8.0. However, when individually inspecting each of the sources, we find that several YSO candidates have IRAC [8.0] mag much fainter than what Whitney et al. (2008) required, with some even as faint as 14.0, as visible in the top
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /><figcaption>**Figure 5.** CMDs showing the α and β cuts from Carlson et al. (2012). The α is shown in red, and the β cut is shown in blue. Previous YSOs are shown with green circles. Newly identified YSOs from this work are shown with orange squares.</figcaption></figure>

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

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

right panel of Figure 5. The previous and new YSO candidates we plot in this panel have a color range between 0.5 and 5.5 and a magnitude range between 6 and 14. The new YSO candidates we identify here fall into a space of the CMD similar to that of Sewiło et al. (2013) and Carlson et al. (2012). Comparing where the new YSO candidates in this work fall in the CMD space to previous YSO candidates confirms that the 299 point sources we identify are likely YSO candidates.

There are 189 point sources in this work that are not YSO candidates based on the criteria we set. However, it is likely that some of them are indeed YSO candidates given that many of the YSOs with spectra do not meet the reduced $\chi^2$ criteria. We could relax our criteria to minimize the false-negative sources. However, that will allow for an increase in false positives in the catalog. The only true method in identifying sources is with spectra. In this work we have made an effort to compile a list of real YSO candidates with minimal contamination from other dusty sources.

# 4. Description of YSO Candidates

## 4.1. Previous YSO Candidates That Fit the Criteria to Still Be YSO Candidates

We have created an inventory of the most massive YSOs in the 30 Doradus region to better understand them in relation to their parental CO cloud. There are 23 YSO candidates from the

literature that satisfy the conditions from Section 3 in our final YSO list (Whitney et al. 2008; Gruendl & Chu 2009; Seale et al. 2014; Nayak et al. 2016), 17 of which have spectra (Seale et al. 2009; Jones et al. 2017). Table 3 lists all 488 point sources and has columns indicating the classification and whether the source was determined to be a YSO candidate using SED modeling or via spectra. We show the location of the 23 likely YSO candidates based on our criteria in Figure 6.

Out of the 17 YSO candidates with spectra, 9 of these sources have reduced $\chi^2_{\text{robitaille}} > 10$. The bad fits to the SED modeling for these nine sources can be explained if some of these sources are in a cluster, where the majority of stars are thought to form (Lada & Lada 2003). Nakajima et al. (2005) studied the N159/N160 region south of the 30 Doradus complex and found that several of the Spitzer YSO candidates are actually clusters. Additionally, Fleener et al. (2010) determined via spectral analysis that 5 of the 10 YSO candidates identified by Gruendl & Chu (2009) in NGC 2074, a star-forming region near 30 Doradus, are actually multiple sources. We can still use the Robitaille (2017) models even if these massive YSOs are indeed a small cluster of YSOs because the most massive and most luminous source in the cluster will dominate the SED (Whitney et al. 2008). The only method to definitively identify YSO candidates is with spectra. Modeling the point sources with SEDs and individually
<figure><img src="image_11.png" /><figcaption>**Figure 6.** Left: the 23 previously studied YSO candidates in the ALMA footprint. Right: the 276 new YSO candidates we identify in this work. The YSO locations are overplotted on the ALMA $^{12}$CO moment 0 map (integrated intensity in units of K km s$^{-1}$). The white dotted line shows the ALMA footprint.</figcaption></figure>

inspecting each source is the next best method in identifying
likely YSO candidates.

The youngest and most embedded YSOs, merely a few tens of thousands of years old, show ice emission and silicate absorption features in their spectra. As YSOs evolve and the intense radiation from the hot, burgeoning protostars creates an H II region, their spectra show fine-structure lines and PAH emission features (Seale et al. 2009). Seale et al. (2009) used Spitzer IRS spectra and categorized the YSOs into six different groups: The S and SE groups have silicate absorption features at 10 and 18 µm and fine-structure lines, respectively; the P and PE groups have PAH emission features and fine-structure lines, respectively; the E group contains candidates with very strong fine-structure lines; and the F group contains sources whose spectra look similar to YSOs but do not fit the criteria of the above groups. The six YSO candidates in this work with spectra from Seale et al. (2009) are categorized as PE, with strong PAH features as well as visible fine-structure line emission. These six YSO candidates are also studied by Jones et al. (2017), who determined that these are late-stage YSOs with H II regions, similar to the 11 additional YSOs from this study. Early-stage YSOs are deeply embedded in dust and gas and have strong silicate and ice absorption features. The presence of fine-structure lines and H II regions is consistent with classifications done by Seale et al. (2009) and Jones et al. (2017).

## 4.2. Previous YSO Candidates That Are Not YSO Candidates in This Work

The Spitzer and Herschel surveys identified millions of point sources in the LMC (Meixner et al. 2006, 2013). Implementing color cuts to the sources was the only feasible method in processing and categorizing them. In this work we look at the region of 30 Doradus observed with ALMA Cycle 7, where there are 488 point sources. We improve on previous methods of identifying YSOs by inspecting each source individually with our criteria described in Section 3. Some of the known

YSO candidates do not fit our criteria because of the assumptions made or color cuts implemented in previous works. There are 10 point sources that were previously classified as YSO candidates but do not meet our criteria in this work. We show the SED for 3 of the 10 sources: candidates “y217” in Figure 2, “y126” in Figure 3, and “y404” in Figure 4. The 7 other previously identified YSO candidates that do not meet our criteria in Section 3 and do not have spectra are “y31,” “y95,” “y134,” “y154,” “y159,” “y296,” and “y307.”

Six of the previously identified YSO candidates that do not meet the criteria in this work were identified by Gruendl & Chu (2009): candidates “y126,” “y134,” “y154,” “y217,” “y2961,” and “y307”. These six sources do not have a reduced $\chi^2_{\text{robitaille}} < 10$. Instead, the reduced $\chi^2_{\text{robitaille}}$ ranges from 14 for “y296” to 441 for “y307.”

Schneider et al. (2018) used the VLT Fiber Large Array Multi-Element Spectrograph (FLAMES) survey to classify 452 sources in 30 Doradus. We exclude three sources previously classified to be YSOs because they were confirmed to be OB stars by Schneider et al. (2018): candidates “y31,” “y95,” and “y159.”

Finally, point source “y404” was identified as a YSO candidate by Whitney et al. (2008) and Seale et al. (2014). Whitney et al. (2008) applied stringent color–magnitude cuts to the Spitzer SAGE point sources and then fit Robitaille et al. (2006, 2007) YSO SED models to their high-probability YSO candidates. Seale et al. (2014) use far-IR Herschel observations and Spitzer 24 µm observations to separate background galaxies from YSO candidates. This object, however, has reduced $\chi^2_{\text{robitaille}} = 12100$ with the Robitaille (2017) SED models, making it unlikely that this is a YSO candidate. The differences in the SED models in Robitaille et al. (2006, 2007) and the updated models from Robitaille (2017) that we use in this work include more evenly sampled parameter space, inclusion of cold dust cooler than 30 K, and a range of model complexity in the updated YSO SED models. Therefore, “oldy35” did fit the Robitaille et al. (2006, 2007) YSO SED
<figure><img src="image_12.png" /><figcaption>**Figure 7.** Radius, mass, luminosity, and temperature histograms of the 299 YSO candidates (23 previously known YSO candidates and 276 newly identified YSO candidates).</figcaption></figure>

models but does not have a good fit to the Robitaille (2017) SED models. In this work, point sources are labeled as YSO candidates only if they meet our criteria, regardless of whether they were previously identified as a YSO candidate or not.

### 4.3. New YSO Candidates

There are 276 point sources that we identify as new YSO candidates in this work, which are shown in the right panel of Figure 6. The stringent YSO color–color or color–magnitude cuts implemented by Whitney et al. (2008), Gruendl & Chu (2009), and Seale et al. (2014) led to these sources being removed from the final catalogs from the respective papers. Whitney et al. (2008) used a series of color–magnitude cuts in their query of the SAGE database using IRAC bands and the MIPS 24 µm band. Gruendl & Chu (2009), on the other hand, used two simple cuts: [4.5] – [8.0] > 2.0 to exclude normal and evolved stars, and [8.0] > 14 – ([4.5] – [8.0]) to exclude background galaxies. Seale et al. (2014) require the presence of 24 µm emission as a proxy for star formation to identify possible YSO candidates. The different methods led to different YSO candidates being identified in each of the respective catalogs, oftentimes with some overlap. For example, YSO candidate “y45” was identified by Gruendl & Chu (2009) and Seale et al. (2014). However, YSO candidate “y49” was identified by only Gruendl & Chu (2009).

We show the radius, mass, luminosity, and temperature distribution of the 276 new YSO candidates and 23 previously identified YSO candidates in Figure 7. These physical properties are from the best-fit Robitaille (2017) models to the SEDs. The radius, mass, and temperature distributions peak at 20 $R_⊙$, 10 $M_⊙$, and 7080 K, respectively. The luminosity distribution of the YSO candidates in this work is double peaked, with one peak at 1780 $L_⊙$ and another at 10,000 $L_⊙$, which can be seen in the bottom left panel of Figure 7. The average mass and luminosity of the 23 previously identified YSO candidates are 17.3 $M_⊙$ and 37,300 $L_⊙$, respectively, while for the 276 new YSO candidates the average mass and luminosity are 10.3 $M_⊙$ and 7900 $L_⊙$, respectively (20% of the average luminosity of the previously known massive YSO candidates). The method of selecting YSO candidates as described in Section 3 identifies lower-mass and lower-luminosity candidates missed by previous studies because we are able to sidestep the stringent and automated color–color cuts. ALMA Cycle 7 observations of CO molecular gas in 30 Doradus cover a region 40 times larger than the Cycle 0 observations of 30 Doradus, within which there are 10 YSO candidates (Gruendl & Chu 2009; Seale et al. 2009; Walborn et al. 2013; Seale et al. 2014; Nayak et al. 2016). Therefore, it is expected to find approximately 400 YSO candidates in a map of the 30 Doradus region that is 40 times larger. We find a total of 299 YSO candidates in the region.
### 4.4. Point Sources That Are Not YSO Candidates

There are 189 point sources in the ALMA footprint that are not YSO candidates, 59 of which are OB stars, based on VLT observations (Schneider et al. 2018). An additional 68 are good fits to the reddened stellar photosphere models (reduced $\chi^2_{kurucz} < 10$), which indicates that these 68 sources are likely stars. However, we label these 68 sources as “other” in Table 3 owing to the lack of additional information.

There are 62 sources that fit with neither YSO SED models nor stellar photosphere SED models, 20 of which have a characteristic steep decline in the near-IR wavelengths in their SEDs and then a jump in 24 µm flux that is comparable to the far-IR fluxes. A few examples of these SEDs can be seen in Figures 2, 3, and 4, all of which show unlikely YSO candidates. These 20 sources could be reflection nebulae. 30 Doradus is one of the most extreme star-forming regions in the LMC, with hundreds of massive protostars. The light from embedded YSOs or more evolved protostars could be reflecting off compact clouds, shifting the original flux to bluer wavelengths. There are another 42 point sources that we are unable to identify as a YSO candidate, a star, or a reflection nebula. Without further evidence, we cannot verify these sources and therefore exclude them from our final list of YSO candidates.

## 5. Results and Discussion

### 5.1. CO Clumps with and without YSOs

We aim to understand the physical processes and feedback mechanisms that lead some molecular gas clumps to gravitationally collapse and form stars and others to not form any stars. Figure 8 shows the $^{12}$CO virial mass, $^{12}$CO virial mass density, $^{12}$CO luminosity-based mass, $^{12}$CO luminosity-based mass density, $^{13}$CO LTE mass, and $^{13}$CO LTE mass density distributions for the dendrogram clumps. Figure 9 shows the $^{12}$CO virial mass, $^{12}$CO virial mass density, $^{12}$CO luminosity-based mass, $^{12}$CO luminosity-based mass density, $^{13}$CO LTE mass, and $^{13}$CO LTE mass density distributions for the SCIMES clusters. More massive stars form in more massive clumps, as can be seen in both Figures 8 and 9. The distribution of YSO candidates greater than 20 $M_\odot$ (seen in green) is similar to the distribution of all YSO candidates (seen in orange) in Figure 8. Furthermore, the distribution of the YSO candidates in Figure 8 is different from the distribution of the CO clumps. This implies that our list of YSO candidates in this work are indeed real YSOs.

We assume that a YSO candidate is associated with a molecular gas clump if it falls within twice the radius of the clump. This assumption is due to the dendrogram algorithm not interpolating the clump radius down to the zero-intensity level; rather, the algorithm clips emission above a certain level determined by the user. We look at Figure 9 to determine the threshold for star formation as inferred by the lowest surface density clump with a YSO candidate. We find that YSO candidates greater than 4 $M_\odot$ form in clumps with virial mass greater than 178 $M_\odot$, CO luminosity mass greater than 6 $M_\odot$, and LTE mass greater than 316 $M_\odot$. YSO candidates greater than 4 $M_\odot$ form in virial mass surface density greater than 32 $M_\odot$ pc$^{-2}$, CO luminosity mass surface density greater than 6 $M_\odot$ pc$^{-2}$, and LTE mass surface density greater than 178 $M_\odot$ pc$^{-2}$. When looking at YSO candidates that are greater than 20 $M_\odot$, these protostars are more likely to form in higher-

mass and higher mass density clumps. There are 22 YSO candidates greater than 20 $M_\odot$ in this work, 16 of which are associated with $^{12}$CO molecular gas. The top right panel of Figure 9 shows that the virial mass surface density needed to form protostars greater than 20 $M_\odot$ is 100 $M_\odot$ pc$^{-2}$. In this work we find that the threshold for forming protostars greater than 20 $M_\odot$ is three times the threshold for forming protostars greater than 4 $M_\odot$ as inferred by the virial mass surface density.

The mass surface density threshold should decrease as metallicity decreases. CO is known to trace the bulk of $H_2$ in solar-metallicity environments (Bolatto et al. 2013). At the lower-metallicity environment, where there are lower dust column densities, $H_2$ remains efficiently self-shielded, whereas CO is more photodissociated. The transition of $C^0/C^+/CO$ is shifted deeper into the cloud, leading to a larger fraction of the CO not tracing the $H_2$ (Wolfire et al. 2010; Bolatto et al. 2011, 2013). Therefore, in low-metallicity environments there exists a substantial amount of $H_2$ outside the CO-emitting regions, which is referred to as “CO-dark” molecular gas. The threshold in the Milky Way is 129 $M_\odot$ pc$^{-2}$ for protostars 0.5 $M_\odot$ in size, as calculated by Heiderman et al. (2010), who use Spitzer surveys and the Gould Belt Survey. Kalari et al. (2020) use ALMA $^{12}$CO (1–0) observations to study the threshold for star formation in the Magellanic Bridge, which has a metallicity of 0.2 $Z_\odot$. They find that the threshold of star formation for 0.4–5 $M_\odot$ protostars in the Magellanic Bridge is 17 $M_\odot$ pc$^{-2}$. The YSOs in this work range from 4 to 38 $M_\odot$, higher protostellar masses than previous studies. We find that the YSO candidates in 30 Doradus form in clumps with virial mass surface density greater than 32 $M_\odot$ pc$^{-2}$, almost twice the threshold in the Magellanic Bridge. It is assumed that all the clumps are in virial equilibrium when using the virial-based mass surface density, which is not true. Using the LTE mass surface density is a more accurate representation of the true local threshold. We find that the LTE-based mass surface density threshold for star formation is 178 $M_\odot$ pc$^{-2}$ for YSO candidates greater than 4 $M_\odot$. It is possible for YSOs smaller than 4 $M_\odot$ to form in gas clumps less dense than what we measure using $^{13}$CO observations. Alternatively, added turbulence from R136 or an increase in ISM pressure in the region is causing the local threshold, as determined by $^{13}$CO LTE analysis, to be almost 40% higher than the star formation threshold in the Milky Way. Further observations of protostars down to 0.5 $M_\odot$ in 30 Doradus and $^{13}$CO molecular gas in the Magellanic Bridge are necessary to determine whether CO is indeed getting increasingly photodissociated with decreasing metallicity, therefore leading to a lower star formation threshold.

Mouschovias & Spitzer (1976) predicted that there is a certain surface density above which magnetic fields can no longer oppose gravitational collapse:

$$ \Sigma_{\text{crit}} > \left( \frac{80}{M_{\odot} \text{ pc}^{-2}} \right) \times \left( \frac{B}{30 \mu\text{G}} \right). $$

For a star formation threshold of 32 $M_\odot$ pc$^{-2}$ inferred from the virial mass surface density, the magnetic field strength needed to balance out gravitational collapse is 12 $\mu$G. This is consistent with the 4–1878 $\mu$G magnetic field strength needed to explain the observed $\nu - \Sigma$ relationship (as calculated in Section 5.1). Observations of the magnetic field strength in 30 Doradus are necessary to determine a critical virial mass surface density above which star formation occurs.
<figure><img src="image_13.png" /><figcaption>**Figure 8.** The ¹²CO virial mass, ¹²CO virial mass density, ¹²CO luminosity-based mass, ¹²CO luminosity-based mass density, ¹³CO LTE mass, and ¹³CO LTE mass density distributions for the SCIMES clusters. Blue represents all the clumps, orange represents clumps associated with the YSOs in this work, and green represents YSOs greater than 20 M⊙.</figcaption></figure>

## 5.2. YSOs Associated with CO Clumps

The YSO candidates in this work are young objects that formed in the past $10^5$ yr and should be embedded within their natal gas cloud. However, not every YSO candidate is associated with CO molecular gas. Of the 23 YSOs previously

studied, 19 are associated with ¹²CO molecular gas, while 16 are associated with ¹³CO molecular gas. Out of the 276 new YSO candidates in this work, 94 YSO candidates are associated with ¹²CO molecular gas, and 18 of the 94 candidates are also associated with ¹³CO molecular gas. In this work 34% of the YSO candidates are associated with ¹²CO and 7% are
<figure><img src="image_14.png" /><figcaption>**Figure 9.** The log of ¹²CO virial mass, ¹²CO virial mass surface density, ¹²CO luminosity-based mass, ¹²CO luminosity-based mass surface density, ¹³CO LTE mass, and ¹³CO LTE mass surface density distributions for all the dendrogram leaves associated with a YSO candidate. Orange represents all leaves associated with the YSOs in this work, and green represents all leaves associated YSOs greater than 20 M⊙.</figcaption></figure>

associated with ¹³CO. Overall, only 38% of the protostars are associated with ¹²CO, leaving the remaining 62% of them seemingly not associated with molecular clouds. Chevance et al. (2020) use Stratospheric Observatory for Infrared Astronomy (SOFIA) Far Infrared Field-Imaging Line Spectrometer (FIFI-LS) observations of [O III] 52 µm, [O III] 88 µm,

[O I] 145 µm, and [C II] 158 µm of 30 Doradus to determine the effects of a strong radiation field from a nearby SSC in a low-metallicity environment. They find that 75%-97% of the H₂ is undetected in the 30 Doradus region (CO-dark molecular gas), similar to our findings that 62% of the protostars greater than 4 M⊙ are not associated with CO. It is indeed possible that
<figure><img src="image_15.png" /><figcaption>**Figure 10.** Correlation between the YSO candidate mass and ¹²CO and ¹³CO clump properties. A YSO candidate is often associated with multiple dendrogram clumps. We chose the smallest dendrogram clump associated with the YSO since we are assuming that these massive protostars form in the dense cores. We indicate the slope and Pearson correlation coefficient.</figcaption></figure>

these protostars are associated with H₂ that is not traced by CO. Alternatively, these massive YSOs could have formed in isolation and quickly dissipated their parental molecular cloud.

In Figure 10, we show how the YSO candidate masses relate to the ¹²CO virial mass, ¹²CO virial mass surface density, ¹²CO luminosity mass, ¹²CO luminosity mass surface density, ¹³CO LTE mass, and ¹³CO LTE mass surface density. In Table 5 we show the slope, intercept, and correlation values. Tables 6 and

7 indicate ¹²CO and ¹³CO clump and YSO candidate associations. The ¹²CO luminosity mass surface density and YSO mass relation have a slope of 9.75 ± 1.14 and a Pearson's r-value of 0.31. The ¹³CO LTE mass surface density–YSO mass relation has a slope of 9.48 ± 2.36 and a Pearson's r-value of 0.34. The stronger correlation between YSO candidate mass and mass surface densities is indicative that the mass densities of the molecular gas clouds play a more important
<table><caption>**Table 5**
Power-law Fit Parameters: $\log Y = a_1 \log X + a_0$</caption><thead><tr><th>X-axis</th><th>Y-axis</th><th>Data Set</th><th>a<sub>1</sub></th><th>a<sub>0</sub></th><th>r-value</th></tr></thead><tbody><tr><td>R</td><td>σ</td><td>All <sup>12</sup>CO</td><td>0.50 ± 0.01</td><td>0.06 ± 0.01</td><td>0.76</td></tr><tr><td>R</td><td>σ</td><td><sup>12</sup>CO with YSOs</td><td>0.49 ± 0.02</td><td>0.06 ± 0.01</td><td>0.73</td></tr><tr><td>R</td><td>σ</td><td>All <sup>13</sup>CO</td><td>0.76 ± 0.05</td><td>0.05 ± 0.01</td><td>0.57</td></tr><tr><td>R</td><td>σ</td><td><sup>13</sup>CO with YSOs</td><td>0.59 ± 0.09</td><td>0.07 ± 0.02</td><td>0.52</td></tr><tr><td>Σ</td><td>ν</td><td>All <sup>12</sup>CO</td><td>0.38 ± 0.01</td><td>-0.76 ± 0.02</td><td>0.42</td></tr><tr><td>Σ</td><td>ν</td><td><sup>12</sup>CO with YSOs</td><td>0.59 ± 0.02</td><td>-1.25 ± 0.05</td><td>0.46</td></tr><tr><td>Σ</td><td>ν</td><td>All <sup>13</sup>CO</td><td>0.24 ± 0.02</td><td>-0.61 ± 0.06</td><td>0.44</td></tr><tr><td>Σ</td><td>ν</td><td><sup>13</sup>CO with YSOs</td><td>0.32 ± 0.08</td><td>-0.85 ± 0.22</td><td>0.49</td></tr><tr><td><sup>12</sup>CO virial mass</td><td>YSO mass</td><td>All <sup>12</sup>CO</td><td>5.49 ± 0.88</td><td>-3.14 ± 2.94</td><td>0.10</td></tr><tr><td><sup>12</sup>CO virial mass density</td><td>YSO mass</td><td>All <sup>12</sup>CO</td><td>15.3 ± 2.17</td><td>-26.0 ± 5.58</td><td>0.21</td></tr><tr><td><sup>12</sup>CO luminosity mass</td><td>YSO mass</td><td>All <sup>12</sup>CO</td><td>4.50 ± 0.60</td><td>4.20 ± 1.33</td><td>0.17</td></tr><tr><td><sup>12</sup>CO luminosity mass density</td><td>YSO mass</td><td>All <sup>12</sup>CO</td><td>9.75 ± 1.14</td><td>-3.67 ± 1.96</td><td>0.31</td></tr><tr><td><sup>13</sup>CO LTE mass</td><td>YSO mass</td><td>All <sup>13</sup>CO</td><td>4.25 ± 1.09</td><td>8.32 ± 2.60</td><td>0.21</td></tr><tr><td><sup>13</sup>CO LTE mass density</td><td>YSO mass</td><td>All <sup>12</sup>CO</td><td>9.48 ± 2.36</td><td>-3.78 ± 5.64</td><td>0.34</td></tr></tbody></table>

**Note.** Column (1): X-axis. Column (2): Y-axis. Column (3): data used for the fit. Column (4): slope of the best-fit line. Column (5): intercept of the best-fit line. Column (6): Pearson's correlation coefficient r-value between X and Y.

<table><caption>**Table 6**
YSO Candidates Associated with the Smallest <sup>12</sup>CO Virial Mass Clump</caption><thead><tr><th>YSO Name</th><th>YSO Mass</th><th><sup>12</sup>CO Clump ID</th><th><sup>12</sup>CO Virial Mass (M<sub>☉</sub>)</th><th><sup>12</sup>CO Virial Mass Err (M<sub>☉</sub>)</th><th><sup>12</sup>CO Virial Mass Density (M<sub>☉</sub> pc<sup>-2</sup>)</th><th><sup>12</sup>CO Virial Mass Density Err (M<sub>☉</sub> pc<sup>-2</sup>)</th><th><sup>12</sup>CO Lum Mass (M<sub>☉</sub>)</th><th><sup>12</sup>CO Lum Mass Err (M<sub>☉</sub>)</th><th><sup>12</sup>CO Lum Mass Density (M<sub>☉</sub> pc<sup>-2</sup>)</th><th><sup>12</sup>CO Lum Mass Density Err (M<sub>☉</sub> pc<sup>-2</sup>)</th></tr></thead><tbody><tr><td>Y2</td><td>32.8</td><td>987</td><td>1.44e+05</td><td>2.27e-02</td><td>1.51e+03</td><td>2.27e-02</td><td>2.88e+04</td><td>1.00e-01</td><td>3.03e+02</td><td>1.00e-01</td></tr><tr><td>Y3</td><td>29.0</td><td>1719</td><td>6.54e+01</td><td>7.73e-01</td><td>2.34e+02</td><td>7.73e-01</td><td>7.15e+01</td><td>1.00e-01</td><td>2.56e+02</td><td>1.00e-01</td></tr><tr><td>Y4</td><td>28.9</td><td>744</td><td>5.25e+02</td><td>2.29e-01</td><td>4.26e+02</td><td>2.29e-01</td><td>2.58e+02</td><td>1.00e-01</td><td>2.09e+02</td><td>1.00e-01</td></tr><tr><td>Y5</td><td>27.3</td><td>1906</td><td>1.27e+04</td><td>5.25e-02</td><td>3.34e+02</td><td>5.25e-02</td><td>6.23e+03</td><td>1.00e-01</td><td>1.64e+02</td><td>1.00e-01</td></tr><tr><td>Y6</td><td>26.2</td><td>910</td><td>1.45e+04</td><td>7.71e-02</td><td>1.84e+03</td><td>7.71e-02</td><td>4.86e+03</td><td>1.00e-01</td><td>6.18e+02</td><td>1.00e-01</td></tr><tr><td>Y9</td><td>24.4</td><td>2290</td><td>1.91e+03</td><td>1.15e-01</td><td>3.23e+02</td><td>1.15e-01</td><td>1.74e+03</td><td>1.00e-01</td><td>2.93e+02</td><td>1.00e-01</td></tr><tr><td>Y13</td><td>21.9</td><td>358</td><td>5.61e+01</td><td>5.41e-01</td><td>3.94e+01</td><td>5.41e-01</td><td>1.63e+01</td><td>1.00e-01</td><td>1.14e+01</td><td>1.00e-01</td></tr><tr><td>Y16</td><td>23.3</td><td>1706</td><td>1.84e+04</td><td>3.66e-02</td><td>2.56e+02</td><td>3.66e-02</td><td>1.15e+04</td><td>1.00e-01</td><td>1.60e+02</td><td>1.00e-01</td></tr><tr><td>Y17</td><td>21.4</td><td>1869</td><td>1.30e+02</td><td>4.82e-01</td><td>3.69e+02</td><td>4.82e-01</td><td>5.91e+00</td><td>1.00e-01</td><td>1.68e+01</td><td>1.00e-01</td></tr><tr><td>Y19</td><td>21.2</td><td>118</td><td>2.40e+03</td><td>2.05e-01</td><td>1.51e+03</td><td>2.05e-01</td><td>1.37e+02</td><td>1.00e-01</td><td>8.62e+01</td><td>1.00e-01</td></tr></tbody></table>

**Note.** Column (1): name of YSO candidate. Column (2): mass of YSO candidate. Column (3): <sup>12</sup>CO clump ID from Wong et al. (2022). This is the 250 m s<sup>-1</sup> feathered mosaic dendrogram output. Column (4): virial mass of the <sup>12</sup>CO clumps. Column (5): error of the virial mass. Column (6): surface density calculated by dividing the virial mass by the size of the dendrogram leaf. Column (7): error in the surface density calculated from the virial mass. Column (8): CO-based mass calculated assuming alpha-scale = 3. Column (9): error in the CO-based mass. Column (10): surface density calculated by dividing the CO-based mass by the size of the dendrogram leaf. Column (11): error in the surface density calculated from the CO-based mass. The full table is available in machine-readable format.
(This table is available in its entirety in machine-readable form.)

role than the masses of the clouds for star formation (i.e., a clump that is 50 M<sub>☉</sub> but has a mass surface density of 1000 M<sub>☉</sub> pc<sup>-2</sup> is more likely to form a protostar than a clump that is 1000 M<sub>☉</sub> but has a density of 50 M<sub>☉</sub> pc<sup>-2</sup>). However, the correlation is not very strong (indicated by the Pearson's r-value), and it could be possible that the individual mass of each protostar is not directly correlated with the clump it forms in. Protostars that do not have properties that correlate with their natal molecular cloud properties may be a consequence of the merger between two different star clusters, or gravitationally unbound stars eventually becoming bound to R136. All 299 YSO candidates in this work are assumed to be 100 Myr old or younger, and we find that they are scattered throughout the region (Figure 6). We examine the clustering of YSOs with respect to the [4.5]–[8.0] color in Figure 11; however, we find no evidence of more embedded or less embedded (i.e., younger or older) YSOs to cluster together. Our findings are consistent with Schneider et al. (2018), who used VLT-FLAMES survey data to search for stars of all ages in the 30 Doradus region,

even at the core of R136 and find no correlation between age and location.

### 5.3. Size, Line Width, and Mass Surface Density

The size-line width relation (or Larson's equation) is an empirical relationship for quiescent molecular clouds with little to no star formation activity, whereby the line width (σ) is proportional to the radius of the cloud raised to the one-half power (r<sup>0.5</sup>; Larson 1981; Lada et al. 1985; Bolatto et al. 2008; Heyer et al. 2009). This relation is thought to be due to the cascade of turbulent energy from large-scale to small-scale clouds (McKee & Ostriker 2007). We investigate the size-line width relation of <sup>12</sup>CO and <sup>13</sup>CO molecular gas clouds in 30 Doradus and compare this relation to others from Galactic and extragalactic studies in Figure 12. For this plot, the ALMA beam size of 30 Doradus observations taken in Cycle 7 is 1"75, and the ALMA beam size of Cycle 0 observations is 2"3 × 1"5. In this work, we use the dendrogram clump and
<table><caption>**Table 7**

YSO Candidates Associated with the Smallest <sup>13</sup>CO Clump Calculated Using the LTE Mass</caption><thead><tr><th>YSO Name</th><th>YSO Mass</th><th><sup>13</sup>CO Clump ID</th><th><sup>13</sup>CO LTE Mass (M<sub>☉</sub>)</th><th><sup>13</sup>CO LTE Mass Err (M<sub>☉</sub>)</th><th><sup>13</sup>CO LTE Mass Density (M<sub>☉</sub> pc<sup>-2</sup>)</th><th><sup>13</sup>CO LTE Mass Density Err (M<sub>☉</sub> pc<sup>-2</sup>)</th></tr></thead><tbody><tr><td>Y3</td><td>29.0</td><td>98</td><td>2.055e+03</td><td>5.540e-03</td><td>5.573e+02</td><td>5.540e-03</td></tr><tr><td>Y4</td><td>28.9</td><td>10</td><td>4.006e+02</td><td>1.390e-02</td><td>3.246e+02</td><td>1.390e-02</td></tr><tr><td>Y6</td><td>26.2</td><td>30</td><td>8.613e+03</td><td>3.090e-03</td><td>1.086e+03</td><td>3.090e-03</td></tr><tr><td>Y9</td><td>24.4</td><td>370</td><td>3.360e+03</td><td>2.560e-03</td><td>7.915e+02</td><td>2.560e-03</td></tr><tr><td>Y22</td><td>20.6</td><td>342</td><td>3.237e+02</td><td>1.800e-02</td><td>2.181e+02</td><td>1.800e-02</td></tr><tr><td>Y23</td><td>20.2</td><td>288</td><td>9.969e+02</td><td>5.350e-03</td><td>7.457e+02</td><td>5.350e-03</td></tr><tr><td>Y26</td><td>20.2</td><td>446</td><td>4.143e+01</td><td>6.070e-02</td><td>8.295e+01</td><td>6.070e-02</td></tr><tr><td>Y27</td><td>20.2</td><td>629</td><td>2.030e+02</td><td>2.420e-02</td><td>2.709e+02</td><td>2.420e-02</td></tr><tr><td>Y34</td><td>18.5</td><td>338</td><td>1.151e+03</td><td>6.120e-03</td><td>1.206e+03</td><td>6.120e-03</td></tr><tr><td>Y36</td><td>18.5</td><td>364</td><td>3.112e+02</td><td>8.980e-03</td><td>6.053e+02</td><td>8.980e-03</td></tr></tbody></table>

Note. Column (1): name of YSO candidate. Column (2): mass of YSO candidate. Column (3): <sup>13</sup>CO clump ID from Wong et al. (2022). This is the 250 m s<sup>-1</sup> feathered mosaic dendrogram output. Column (4): LTE mass of the <sup>13</sup>CO clumps. Column (5): error of the LTE mass. Column (6): surface density calculated by dividing the LTE mass by the size of the dendrogram leaf. Column (7): error in the surface density calculated from the LTE mass. The full table is available in machine-readable format.

(This table is available in its entirety in machine-readable form.)

SCIMES cluster properties from Wong et al. (2022). The slope and intercept in the size-line width relation for all <sup>12</sup>CO clumps are 0.50 ± 0.01 and 0.06 ± 0.01, respectively. For <sup>12</sup>CO clumps associated with YSOs the slope and intercept are similar, 0.49 ± 0.02 and 0.06 ± 0.01, respectively. This is given in Table 5. Nayak et al. (2016) looked at ALMA Cycle 0 observations of <sup>12</sup>CO clumps within 11 pc of the R136 SSC and found the slope and intercept of 0.65 ± 0.04 and 0.22 ± 0.01, respectively. Clumps that are within 11 pc of the R136 SSC are therefore offset to higher line widths (intercept of 0.22 on a log scale) in comparison to all the <sup>12</sup>CO clumps in the 30 Doradus region (intercept of 0.06), which can be seen in Figure 12. We observe that the clumps closer to the SSC have higher line widths for a given size in comparison to clumps farther away (up to 45 pc away). The R136 SSC has hundreds of massive O stars and intense radiation and stellar winds. One possible scenario is that energy from R136 is being injected into clumps, such that those more nearby to the SSC are injected with proportionally more energy than those lying a greater distance away.

This trend of higher line widths for a given clump size for clumps that are closer to the SSC seen in <sup>12</sup>CO is also true for <sup>13</sup>CO (Figure 12). Nayak et al. (2016) find the slope and intercept for <sup>13</sup>CO clumps that are within 11 pc of R136 to be 0.97 ± 0.12 and 0.19 ± 0.05, respectively. We find the slope and intercept for all <sup>13</sup>CO clumps within 45 pc to be 0.76 ± 0.05 and 0.05 ± 0.01, respectively (Table 5). The <sup>13</sup>CO clumps that are within 11 pc of R136 are offset to 60% larger line widths in comparison to <sup>13</sup>CO clumps in the entire 30 Doradus region mapped in Cycle 7 (intercept of 0.19 versus intercept of 0.05 on a log scale). Even when looking at just clumps with star formation activity, the <sup>13</sup>CO clumps in the immediate vicinity of R136 are offset to line widths 30% higher (intercept of 0.18 ± 0.07 versus 0.07 ± 0.02). The slope of the best-fit line for <sup>13</sup>CO clumps within 11 pc of R136 is also steeper than the slope when looking at clumps in a much larger region of 30 Doradus (slope of 1.02 ± 0.17 versus slope of 0.76 ± 0.05). Clumps that are farther away from the R136 SSC are not as affected by the intense star formation occurring in the region in comparison to clumps that are closer. The intercepts and the

slopes of the best-fit line in the size-line width relation indicate that more energy might be getting injected into the clumps closer to R136. We assume that the line widths are dominated by nonthermal broadening (i.e., turbulence). Observations by Kalari et al. (2018) of three clouds in 30 Doradus within 20 pc of R136 indicate that feedback from massive stars also contributes to the measured line widths. Kalari et al. (2018) present <sup>12</sup>CO observations using data from the Swedish-ESO Submillimetre Telescope (SEST). The line widths they observed are large for a given cloud size in comparison to Larson's empirical relation, and this increase in observed line width is likely because the molecular gas is being pushed away by the ionizing radiation produced by the massive stars in R136 (Kalari et al. 2018).

Comparisons of high-mass star-forming regions similar to R136 are in agreement with our findings that line widths are higher for a given clump size. NGC 604 is located in M33 and has the highest star formation rate (SFR) in the entire galaxy, with over 200 O stars surrounded by filaments and H II shells (Relaño & Kennicutt 2009; Miura et al. 2012). Though not compact enough to be considered an SSC, the population of O stars in NGC 604 is similar to that of R136. Phiri et al. (2021) study the <sup>13</sup>CO molecular gas in NGC 604 and use dendrograms to analyze the clump structure and properties. While the slope and the intercept of their size-line width relation are not reported, using the provided data, we have determined the slope and intercept to be 2.89 ± 0.59 and -2.44 ± 0.48, respectively. The slope of the size-line width relation in NGC 604 is 2.8 times higher than the slope of the <sup>13</sup>CO clumps within 11 pc of R136, indicating high amounts of turbulence in NGC 604. Even in Milky Way environments, sites of star formation show an increased line width for a given clump size relative to quiescent regions. Ma et al. (2021) study molecular clouds in the second quadrant of the Milky Way, a relatively quiescent region, as well as molecular clouds in Perseus. They find $\sigma \propto r^{0.27}$ for <sup>12</sup>CO clumps and $\sigma \propto r^{0.29}$ for <sup>13</sup>CO clumps in the Milky Way second quadrant. The slope of the size-line width relation for clumps in the star-forming region Perseus is 50% larger, with a slope of 0.44 for <sup>12</sup>CO clumps and a slope of 0.43 for <sup>13</sup>CO clumps. The increase in
<figure><img src="image_16.png" /><figcaption>**Figure 11.** We plot the location of the YSO candidates from the image in Figure 5 based on their [4.5]–[8.0] color. The color scale is in units of K km s⁻¹.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure 12.** Top left: the size–line width relation of the $^{12}$CO dendrogram clumps. We assume that the best-fit linear equation is given by $log Y = a_1 log X + a_0$. The best-fit line is given by $(a_1, a_0) = (0.50, 0.06)$ for all $^{12}$CO clumps from Wong et al. (2022), while the best-fit line is given by (0.49, 0.06) for $^{12}$CO clumps associated with YSOs. Top right: The size–line width relation of the $^{13}$CO dendrogram clumps. The best-fit line is given by (0.76, 0.05) for all $^{13}$CO clumps from Wong et al. (2022), and the best-fit line is given by (0.59, 0.07) for $^{13}$CO clumps associated with YSOs. We apply spatial and spectral resolution limits of 0.344 pc and 0.106 km s$^{-1}$ to the $^{12}$CO and $^{13}$CO clumps in this figure. Bottom: size–line width relation from previous studies and this work.</figcaption></figure>

the size–line width slope in Perseus in comparison to the second quadrant of Milky Way is similar to the increase in slope when observing CO clumps in the immediate vicinity of R136 in comparison to the CO clumps in the entire map of 30 Doradus.

Simulations indicate that stellar winds play an important role in the early formation stages for stellar clusters. Rathjen et al. (2021) simulate the ISM with both thermal and nonthermal processes such as stellar winds, hydrogen ionizing radiation, supernovae, and cosmic rays. They find that radiation and winds from stellar clusters play a major role in regulating star formation on timescales of 100 Myr. High-mass star-forming regions with hundreds of O stars such as R136 and NGC 604 have higher line widths for a given clump size than nearby quiescent clouds like those studied by Solomon et al. (1987). Turbulence from the strong radiation from high-mass star-forming regions is likely getting injected into CO clumps, and we see this effect in the size–line width relation.

A different approach to illustrating Larson's equation is to consider the quantity $\sigma/r^{0.5}$ (we will define as $\nu$) versus the mass surface density $\Sigma$. Assuming that clumps are self-gravitating, the virial mass can be given by $M_{vir} = 5\sigma^2r/G$. Because $\Sigma = M_{vir}/\pi r^2$, we can substitute $\Sigma$ for $M_{vir}$. Solving

for line width $\sigma$, we find

$$ \sigma = \left(\frac{\pi G}{5}\right)^{0.5} \Sigma^{0.5} r^{0.5}. $$

Subsequently, $\sigma/r^{0.5} = \nu = (\pi G\Sigma/5)^{0.5}$. By Larson's scaling relations, $\Sigma$ is constant for all clouds, and therefore the relationship of $\nu$ against $\Sigma$ should yield a flat slope. In the left panel of Figure 13, we show the $\nu-\Sigma$ relationship for $^{12}$CO, which has a slope of $0.38 \pm 0.01$, and in the right panel we show the relationship for $^{13}$CO, which has a slope of $0.24 \pm 0.02$.

The presence of turbulence and radiation from R136 may explain the observed deviance from a flat slope as predicted by Larson's $\nu-\Sigma$ relation. Another possibility is the propagation of large-amplitude, long-wavelength Alfvén waves through the clouds (Arons & Max 1975). The difference in the slope in the $\nu-\Sigma$ relation in different environments could additionally be due to the difference in the local magnetic field strength (Heyer et al. 2009). Myers & Goodman (1988) assume equipartition between magnetic, kinetic, and gravitational energy and determine the magnetic field strength to be $B \approx (45/G)^{0.5}(\sigma^2/r)$. Using the sizes and line widths of all the $^{12}$CO clumps in 30 Doradus as identified by Wong et al. (2022) and the previous equation, the magnetic field in 30 Doradus
<figure><img src="image_18.png" /><figcaption>**Figure 13.** Left: the $v-\Sigma$ relation of the $^{12}\text{CO}$ dendrogram clumps. According to Larson's equations, the slope is supposed to be flat. However, we find in this work that regions of high-mass star formation are different from the quiescent CO cloud studies by Larson. The best-fit line is given by $(0.38, -0.76)$, assuming an $(a_1, a_0)$ format, for all $^{12}\text{CO}$ clumps from Wong et al. (2022). And the best-fit line is given by $(0.59, -1.25)$ for $^{12}\text{CO}$ clumps associated with YSOs. We use the surface density calculated using the $^{12}\text{CO}$ virial mass divided by the area of the clump. Right: the $v-\Sigma$ relation of the $^{13}\text{CO}$ dendrogram clumps. The best-fit line is given by $(0.24, -0.61)$ for all $^{13}\text{CO}$ clumps from Wong et al. (2022), and the best-fit line is given by $(0.32, -0.85)$ for $^{13}\text{CO}$ clumps associated with YSOs. We use the surface density calculated using the $^{13}\text{CO}$ LTE mass divided by the area of the clump.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Figure 14.** Left: relation between the virial parameter $\alpha$ and the mass derived from the $^{12}\text{CO}$ luminosity. We define $\alpha$ to be equal to $\frac{5\sigma^2r}{GM}$. Black squares are all the $^{12}\text{CO}$ clumps, and the cyan circles are the $^{12}\text{CO}$ clumps with a YSO candidate. Right: relation between the virial parameter $\alpha$ and the mass derived from the SCIMES $^{12}\text{CO}$ cluster luminosity.</figcaption></figure>

ranges from 4 to 1878 µG, with a median of 95 µG. Con- firmation of the magnetic field strength in 30 Doradus requires polarimetric observations. As the measurements of Crutcher (1999) found values between 10 and 3100 µG in star-forming molecular clouds, it is therefore possible that the 4–1878 µG required to explain the $v-\Sigma$ relationship in 30 Doradus is present.

### 5.4. Energy Balance

The $^{12}\text{CO}$ virial mass density is indicative of the gravita- tional potential energy, while the luminosity mass density is indicative of the kinetic energy of the clumps. In this work, we find that the luminosity mass density is on average larger than the virial mass density, $\Sigma_{\text{vir}} \propto \Sigma_{\text{CO}}^{0.91}$, consistent with Wong et al. (2022).

We analyze the dynamics of the clumps and define the virial parameter as $\alpha = 5\sigma^2r/GM$, which is the ratio of the kinetic energy and the gravitational potential energy. In Figure 14, we plot the virial parameter against the virial mass of all $^{12}\text{CO}$

clumps on the top panels and that of only the $^{12}\text{CO}$ SCIMES clusters in the bottom panels. When clumps have virial parameters $\alpha < 2$, they are considered to be gravitationally bound and are more likely to collapse and form stars. When clumps have virial parameters $\alpha > 2$, they are considered to be unbound and can either expand or be mixed in with the diffuse ISM.

### 5.5. Virial Analysis

The virial parameter of $^{12}\text{CO}$ clumps in 30 Doradus ranges from 0.4 to 197. About 80% of the $^{12}\text{CO}$ clumps in 30 Doradus are not gravitationally bound, and 70% of the star-forming $^{12}\text{CO}$ clumps are also not bound. The virial parameter for the $^{13}\text{CO}$ clumps ranges from 0.3 to 13. Approximately 60% of the $^{13}\text{CO}$ clumps in 30 Doradus are gravitationally bound, whereas the majority of $^{12}\text{CO}$ are gravitationally unbound. This difference in the virial parameter between the two CO isotopologues is because $^{12}\text{CO}$ has a critical density of $2.3 \times 10^3 \text{ cm}^{-3}$, whereas $^{13}\text{CO}$ has a critical density of $2.0 \times 10^4 \text{ cm}^{-3}$. $^{13}\text{CO}$ molecular
<figure><img src="image_20.png" /><figcaption>**Figure 15.** Top: the histogram of the YSOs (blue) and the best-fit Kroupa IMF (orange; for YSOs with $M > 8 M_\odot$). We assume that we are complete for YSO candidates above $8 M_\odot$ because that is the peak of the distribution. The total YSO mass in the ALMA footprint region of 30 Doradus is $1.8 \times 10^4 M_\odot$. Assuming a typical YSO formation timescale of $10^5$ yr, the SFR is then $0.18 M_\odot$ yr$^{-1}$. Bottom: the histogram of the 23 previous YSO candidates and the best-fit Kroupa IMF. The total YSO mass assumed from the IMF fit is $1.2 \times 10^3 M_\odot$, and the SFR is $0.012 M_\odot$ yr$^{-1}$.</figcaption></figure>

gas is a better indicator of the high-density knots in which star formation takes place, and therefore the majority of the $^{13}$CO clumps are gravitationally bound.

## 5.6. Initial Mass Function and Efficiency of Star Formation

We fit a three-part Kroupa initial mass function (IMF; Kroupa 2001) to the YSO candidates in this work. The sample of YSO candidates is assumed to be complete above $8 M_\odot$, which is the peak of the distribution in the top panel of Figure 15. We fit this observed peak as the expected completeness value with the following IMF: $\xi \propto M^{-0.3}$ for $M < 0.5 M_\odot$, $\xi \propto M^{-1.2}$ for $0.5 M_\odot < M < 1.0 M_\odot$, and $\xi \propto M^{-1.7}$ for $M > 1.0 M_\odot$ (Kroupa 2001; Romano et al. 2005). Other common IMFs used are Salpeter (1955) and Chabrier (2003). A comparison between Kroupa, Salpeter, Chabrier, and other IMFs can be found in Figure 1 of Romano et al. (2005). The Salpeter and Chabrier IMFs lie above the Kroupa IMF for protostellar masses ranging between 10 and 100 $M_\odot$; therefore, the Kroupa IMF gives the best fit to the YSO candidates in this work.

The IMF is incomplete at the low-mass end because we cannot observe protostars smaller than a few solar masses at the distance of the LMC, and it is incomplete at the high-mass end

owing to observational bias effects such as saturation and crowding. The total YSO mass derived from integrating the best-fit IMF to the YSO candidates in the ALMA footprint of 30 Doradus is 18,000 $M_\odot$. Assuming a typical YSO formation timescale of $10^5$ yr (McKee & Tan 2002), the SFR is then approximately $0.18 M_\odot$ yr$^{-1}$. Carlson et al. (2012) calculate the total YSO mass and SFR via integration over a Kroupa IMF and assuming a formation timescale of 1 Myr. If the SFR is scaled from Table 8 in Carlson et al. (2012) for a formation timescale of $10^5$ yr, then the SFR in their work ranges from $0.008 M_\odot$ yr$^{-1}$ in N144 to $0.068 M_\odot$ yr$^{-1}$ in N160. The four known episodes of star formation and the current radiation from the massive stars have created conditions for the 30 Doradus region to have one of the highest SFRs in the LMC, 2.6 times higher than the SFR in N160 (another high-mass star-forming GMC adjacent to and toward the south of 30 Doradus).

We fit a Kroupa IMF to the 23 previously known YSO candidates in the ALMA footprint of 30 Doradus (bottom panel of Figure 15). Integrating over the Kroupa IMF while adopting the same formation timescale and completeness limit as the previous paragraph, we find that the total mass of YSOs would be $1300 M_\odot$ and the SFR would be $0.013 M_\odot$ yr$^{-1}$, five times smaller than the SFR in N160. The 30 Doradus region is one of
the most luminous giant H II regions in the Local Group (Kennicutt 1984) and is known to host numerous O stars (Melnick 1985; Pehlemann et al. 1992). By identifying an additional 276 YSO candidates, we improve on the inferred SFR.

Ochsendorf et al. (2017) calculated the SFR using $H\alpha$, which they assumed to be a proxy for the average SFR over the past 4 Myr. They find $SFR_{H\alpha} = 0.03 M_{\odot} yr^{-1}$, 17% of the $0.18 M_{\odot} yr^{-1}$ rate we calculated in this work. An SFR calculated using YSOs is indicative of the rate averaged over the past 0.5 Myr. Comparing the SFR from YSOs in this work to that calculated using $H\alpha$ in Ochsendorf et al. (2017), we find that the SFR has increased over time in 30 Doradus. This is consistent with evidence of ongoing star formation activity that shows 16,700 PMS stars in the 30 Doradus region (Sabbi et al. 2013; Walborn et al. 2013).

## 5.7. GMC and Stellar Cluster Origin

An alternative to energy being injected into the clumps is that the ISM pressure can confine and compress the observed clumps, increasing their surface density and observed line widths (Hughes et al. 2013; Utomo et al. 2015; Sun et al. 2020). This pressure can be from intercloud thermal pressure or the pressure of the entire large-scale clump acting on the clumps within it (Indebetouw et al. 2013). The average pressure of the $^{13}$CO clumps we measure in this work can be given by

$$ \frac{P}{k_B} = 1.5 \left(\frac{M}{M_{\odot}}\right)^2 \left(\frac{r}{\text{pc}}\right)^{-4} \text{cm}^{-3} \text{K.} $$

We use the $^{13}$CO clump masses and radii from Wong et al. (2022) to find an average pressure in 30 Doradus of $6.5 \times 10^6 \text{cm}^{-3} \text{K}$. The average pressure we calculate using $^{13}$CO molecular gas observations is consistent with pressure inferred from modeling optical diagnostic emission lines $H\alpha$, $H\beta$, $H\gamma$, [O III], [N II], [S II], and [S III] (Pellegrini et al. 2011). Pellegrini et al. (2011) use observations from HST and the Magellanic Cloud Emission Line Survey (MCELS), model the emission lines using the CLOUDY photoionization tool, and find that the pressure in 30 Doradus ranges from $9.4 \times 10^5 \text{cm}^{-3} \text{K}$ to $9.7 \times 10^6 \text{cm}^{-3} \text{K}$. For comparison, Faesi et al. (2018) find an ISM pressure in the Milky Way in the range of $(8.5-19.8) \times 10^3 \text{cm}^{-3} \text{K}$ and an ISM pressure in M51 in the range of $(0.5-10) \times 10^5 \text{cm}^{-3} \text{K}$ on scales of 10–100 pc. The pressure in 30 Doradus is comparable to the pressure in the CMZ, which is estimated to be between $6 \times 10^6 \text{cm}^{-3} \text{K}$ and $1 \times 10^8 \text{cm}^{-3} \text{K}$ (Spergel & Blitz 1992; Miyazaki & Tsuboi 2000).

The elevated line widths for a given clump size we calculate in this work can be due to the molecular gas clumps being situated in a high-pressure environment. Figure 12 shows that there is not difference in the size–line width relation of all $^{12}$CO and $^{13}$CO clumps versus the star-forming clumps. The presence of YSOs does not affect the clump itself, and it is not what is creating this turbulence inferred from the elevated line widths for a given clump size. Turbulence is not driven by pressure from the star formation activity in 30 Doradus.

We find the 299 YSO candidates to be scattered throughout the complex and clumpy CO cloud structure in 30 Doradus, a region that has undergone four major star formation episodes. A simple model of H II shells expanding into an ambient medium is not an accurate scenario that can be used to explain how star formation is getting triggered or quenched from the

feedback. Walch et al. (2013) use a high-resolution numerical simulation to determine the nature of the molecular clouds in which star formation takes place. They find that for clouds in which the density structure is dominated by large-scale fluctuations, the resulting H II regions have shell-like structures along the borders. These shell-like structures break up into a small number of massive high-density clumps, which then spawn star clusters. Star formation then occurs quickly and results in higher stellar masses. Walch et al. (2013) find that for clouds where small-scale fluctuations become important the border of the H II region is dominated by a large number of pillars. These pillars contain compact dense clumps and tend to spawn single stars or individual multiple systems. Star formation then occurs later, the stellar masses are lower, and the stars are more widely distributed. A collision between large-scale filaments or H I shells can explain the turbulence we infer from the elevated line width measurements. The most recent star formation episode took place 2 Myr ago when two H I clouds with an estimated mass of $1 \times 10^6 M_{\odot}$ started colliding, and at the intersection of this collision lies the SSC R136 (Fukui et al. 2017). Of the two H I components colliding in 30 Doradus, one is shaped like a shell/arc and the other is shaped like a pillar/ridge (Fukui et al. 2017). The merging of these two H I shells suggests that a mixture of the two modes of star formation predicted by Walch et al. (2013) is occurring in 30 Doradus, resulting in the formation of the 299 YSO candidates we observe with ALMA Cycle 7.

## 6. Conclusions

We present a catalog of 488 point sources within the ALMA Cycle 7 footprint of 30 Doradus. From fitting SEDs we determine that 299 of these point sources are likely YSO candidates, 276 of which are newly identified protostars. Our method of selecting YSO candidates sidesteps the issue of removing potential good candidates with stringent color cuts. There are 189 point sources in the ALMA footprint that are not YSO candidates, 20 of which may be reflection nebulae. These 20 objects have a steep decline in near-IR wavelengths with a jump at 24 µm flux. We use dendrograms to analyze the $^{12}$CO and $^{13}$CO molecular gas (Wong et al. 2022) and determine the following:

1. We find that 113 of 299 YSO candidates are associated with $^{12}$CO molecular gas, leaving the remaining 62% seemingly not associated with the observed clouds. However, these YSO candidates may potentially be in “CO-dark” clouds, those $H_2$ reservoirs not traced by CO. Indeed, up to 97% of the $H_2$ is undetected by CO in 30 Doradus (Chevance et al. 2020). Another potential scenario is that these YSOs could have formed in isolation and quickly dissipated their parental molecular cloud.

2. The threshold for star formation in 30 Doradus is $32 M_{\odot} \text{pc}^{-2}$ for $4 M_{\odot}$ stars and $100 M_{\odot} \text{pc}^{-2}$ for stellar masses greater than $20 M_{\odot}$ as inferred by the virial mass surface density, which is less than the threshold for low-mass star formation in the Milky Way. This may be due to added turbulence and pressure. By comparison, the threshold for star formation inferred using LTE surface mass density is $178 M_{\odot} \text{pc}^{-2}$.

3. If the magnetic field strength in 30 Doradus is 12 µG, then clouds with mass surface density threshold of
32 $M_{\odot}$ pc$^{-2}$ or higher can undergo gravitational collapse. Observations of the magnetic field, as well as CO observations of other low-metallicity star-forming regions, are needed for direct comparisons.

4. The $^{12}$CO size–line width relation for all clumps has a slope of $0.50 \pm 0.01$ and intercept of $0.06 \pm 0.01$. For $^{13}$CO, the slope and intercept are $0.76 \pm 0.05$ and $0.05 \pm 0.01$, respectively, in the size–line width relation. Nayak et al. (2016) analyzed the ALMA Cycle 0 observations of molecular gas within 11 pc of R136 and find that the intercept of the best fit is $0.22 \pm 0.01$. The clumps that are within 11 pc have elevated line widths by a factor of 1.4 in comparison to all clumps within 45 pc of R136. It is possible that turbulent energy from R136 is being injected into the molecular gas clumps.

5. Larson's relations imply that the slope should be flat when plotting $v$ versus $\Sigma$. Large-amplitude, long-wavelength Alvén waves propagating through the cloud can explain why we do not find the slope to be flat. Assuming equipartition between magnetic, kinetic, and gravitational energy, we calculate that a magnetic field strength of 4–1878 $\mu$G is enough to explain the observed $v-\Sigma$ relation.

6. We fit a Kroupa IMF to the YSO candidates in this work and find that the total stellar mass within the ALMA Cycle 7 footprint of 30 Doradus is $1.8 \times 10^4 M_{\odot}$. If we assume a formation timescale of $10^5$ yr, the SFR is therefore $0.18 M_{\odot} yr^{-1}$. The SFR calculated using $H\alpha$, proxy for the average SFR over the past 4 Myr, is $0.03 M_{\odot} yr^{-1}$. Comparing the SFR from $H\alpha$ to that from YSOs, we find that the star formation activity has increased over the past few million years.

7. Pressure in the ISM from thermal sources or as a consequence of large-scale clumps affecting smaller clumps within can compress the molecular gas. The pressure needed to produce the observed properties for $^{13}$CO is $6.5 \times 10^6$ cm$^{-3}$ K, similar to the pressure inferred from molecular gas observations of the central molecular zone.

8. R136 is thought to have formed from the collision of two H I flows, one that is shaped like an arc and another that is shaped like a pillar (Fukui et al. 2017). Comparing numerical simulations by Walch et al. (2013) to observational results by Fukui et al. (2017), we theorize that a combination of high-mass stars, which form in shell-like structures, and low-mass stars, which form in pillar-like structures, exists in 30 Doradus because of the collision of the two H I flows.

This paper makes use of the following ALMA data: ADS/JAO.ALMA #2019.1.00843.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSC 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.

This work makes use of data collected by the Herschel Space Observatory. Herschel is an ESA space observatory with science instruments provided by European-led Principal

Investigator consortia and with important participation from NASA. This work is also based in part on observations made with the Spitzer Space Telescope, which is operated by the Jet Propulsion Laboratory, California Institute of Technology under a contract with NASA.

This research made use of SCIMES, a Python package to find relevant structures into dendrograms of molecular gas emission using the spectral clustering approach.

O.N. was supported by STScI Director's Discretionary Fund. T.W., M.M., and R.I. acknowledge support from collaborative NSF AAG awards 2009849, 2009544, and 2009624. M.M. acknowledges support from NSF, grant No. 2054178. M.C. gratefully acknowledges funding from the Deutsche Forschungsgemeinschaft (DFG) in the form of an Emmy Noether Research Group (grant Nos. KR4801/1-1 and CH 2137/1-1) and the DFG Sachbeihilfe (grant No. KR4801/2-1), and from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program via the ERC Starting Grant MUSTANG (grant agreement No. 714907). COOL Research DAQ is a Decentralized Autonomous Organization supporting research in astrophysics aimed at uncovering our cosmic origins. M.R. wishes to acknowledge support from ANID(CHILE) through FONDECYT grant No. 1190684.

## ORCID iDs

Omnarayani Nayak https://orcid.org/0000-0001-6576-6339

Alex Green https://orcid.org/0000-0002-8432-3362

Alec S. Hirschauer https://orcid.org/0000-0002-2954-8622

Rémy Indebetouw https://orcid.org/0000-0002-4663-6827

Tony Wong https://orcid.org/0000-0002-7759-0585

Mélanie Chevance https://orcid.org/0000-0002-5635-5180

Guido De Marchi https://orcid.org/0000-0001-7906-3829

Vianney Lebouteiller https://orcid.org/0000-0002-7716-6223

Min-Young Lee https://orcid.org/0000-0002-9888-0784

Leslie W. Looney https://orcid.org/0000-0002-4540-6587

Suzanne C. Madden https://orcid.org/0000-0003-3229-2899

Julia Roman-Duval https://orcid.org/0000-0001-6326-7069

Yasuo Fukui https://orcid.org/0000-0002-8966-9856

Alvaro Hacar https://orcid.org/0000-0001-5397-6961

K. E. Jameson https://orcid.org/0000-0001-7105-0994

Venu Kalari https://orcid.org/0000-0002-4641-2532

Mónica Rubio https://orcid.org/0000-0002-5307-5941

Elena Sabbi https://orcid.org/0000-0003-2954-7643

## References

Arons, J., & Max, C. E. 1975, ApJ, 196, L77

Bally, J. 2016, ARA&A, 54, 491

Barnes, A. T., Longmore, S. N., Dale, J. E., et al. 2020, MNRAS, 498, 4906

Bolatto, A. D., Leroy, A. K., Jameson, K., et al. 2011, ApJ, 741, 12

Bolatto, A. D., Leroy, A. K., Rosolowsky, E., Walter, F., & Blitz, L. 2008, ApJ, 686, 948

Bolatto, A. D., Wolfire, M., & Leroy, A. K. 2013, ARA&A, 51, 207

Carlson, L. R., Sewilo, M., Meixner, M., Romita, K. A., & Lawton, B. 2012, A&A, 542, A66

Castelli, F., & Kurucz, R. L. 2004, A&A, 419, 725

Chabrier, G. 2003, PASP, 115, 763

Chevance, M., Madden, S. C., Fischer, C., et al. 2020, MNRAS, 494, 5279

Chevance, M., Madden, S. C., Lebouteiller, V., et al. 2016, A&A, 590, A36

Colombo, D., Rosolowsky, E., Ginsburg, A., Duarte-Cabral, A., & Hughes, A. 2015, MNRAS, 454, 2067
Crowther, P. A., Caballero-Nieves, S. M., Bostroem, K. A., et al. 2016, MNRAS, 458, 624

Crutcher, R. M. 1999, ApJ, 520, 706

De Marchi, G., Paresce, F., Panagia, N., et al. 2011, ApJ, 739, 27

Draine, B. T. 2003a, ApJ, 598, 1017

Draine, B. T. 2003b, ARA&A, 41, 241

Faesi, C. M., Lada, C. J., & Forbrich, J. 2018, ApJ, 857, 19

Fall, S. M., Krumholz, M. R., & Matzner, C. D. 2010, ApJL, 710, L142

Fazio, G., & Eisenhardt, P. 2004, Spitzer spectroscopy of 8 micron Bright, Optically Faint Objects Identified in the IRAC Shallow Survey, Spitzer Proposal, 111

Fleener, C. E., Payne, J. T., Chu, Y.-H., Chen, C.-H.-R., & Gruendl, R. A. 2010, AJ, 139, 158

Fukui, Y., Harada, R., Tokuda, K., et al. 2015, ApJL, 807, L4

Fukui, Y., Kawamura, A., Wong, T., et al. 2009, ApJ, 705, 144

Fukui, Y., Tsuge, K., Sano, H., et al. 2017, PASJ, 69, L5

Grebel, E. K., & Chu, Y.-H. 2000, AJ, 119, 787

Gruendl, R. A., & Chu, Y.-H. 2009, ApJS, 184, 172

Heiderman, A., Evans, Neal, J., et al. 2010, ApJ, 723, 1019

Heyer, M., Krawczyk, C., Duval, J., & Jackson, J. M. 2009, ApJ, 699, 1092

Hora, J. L., Cohen, M., Ellis, R. G., et al. 2008, AJ, 135, 726

Hughes, A., Meidt, S. E., Colombo, D., et al. 2013, ApJ, 779, 46

Hunter, D. A., Shaya, E. J., Scowen, P., et al. 1995, ApJ, 444, 758

Indebetouw, R., Brogan, C., Chen, C.-H.-R., et al. 2013, ApJ, 774, 73

Indebetouw, R., de Messières, G. E., Madden, S., et al. 2009, ApJ, 694, 84

Jones, O. C., Woods, P. M., Kemper, F., et al. 2017, MNRAS, 470, 3250

Kalari, V. M., Rubio, M., Elmegreen, B. G., et al. 2018, ApJ, 852, 71

Kalari, V. M., Rubio, M., Saldaño, H. P., & Bolatto, A. D. 2020, MNRAS, 499, 2534

Kato, D., Nagashima, C., Nagayama, T., et al. 2007, PASJ, 59, 615

Kauffmann, J., Pillai, T., Shetty, R., Myers, P. C., & Goodman, A. A. 2010, ApJ, 716, 433

Kennicutt, R. C., Jr. 1984, ApJ, 287, 116

Koo, B.-C., Kim, C.-G., Park, S., & Ostriker, E. C. 2020, ApJ, 905, 35

Kroupa, P. 2001, MNRAS, 322, 231

Krumholz, M. R. 2014, PhR, 539, 49

Ksoll, V. F., Gouliermis, D. A., Klessen, R. S., et al. 2018, MNRAS, 479, 2389

Lada, C. J., & Lada, E. A. 2003, ARA&A, 41, 57

Lada, E. A., Magnani, L., & Blitz, L. 1985, BAAS, 17, 616

Lancaster, L., Ostriker, E. C., Kim, J.-G., & Kim, C.-G. 2021, ApJ, 914, 89

Larson, R. B. 1981, MNRAS, 194, 809

Lee, M.-Y., Madden, S. C., Le Petit, F., et al. 2019, A&A, 628, A113

Lee, M.-Y., Madden, S. C., Lebouteiller, V., et al. 2016, A&A, 596, A85

Lopez, L. A., Krumholz, M. R., Bolatto, A. D., et al. 2014, ApJ, 795, 121

Lopez, L. A., Krumholz, M. R., Bolatto, A. D., Prochaska, J. X., & Ramirez-Ruiz, E. 2011, ApJ, 731, 91

Ma, Y., Wang, H., Li, C., et al. 2021, ApJS, 254, 3

Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415

Matzner, C. D., & McKee, C. F. 2000, AAS Meeting, 195, 135.07

McKee, C. F., & Ostriker, E. C. 2007, ARA&A, 45, 565

McKee, C. F., & Tan, J. C. 2002, Natur, 416, 59

Meixner, M., Gordon, K., Indebetouw, R., et al. 2006, AJ, 132, 2268

Meixner, M., Panuzzo, P., Roman-Duval, J., et al. 2013, AJ, 146, 62

Melnick, J. 1985, A&A, 153, 235

Melnick, J., Tenorio-Tagle, G., & Telles, E. 2021, A&A, 649, A175

Minamidani, T., Tanaka, T., Mizuno, Y., et al. 2011, AJ, 141, 73

Miura, R. E., Kohno, K., Tosaki, T., et al. 2012, ApJ, 761, 37

Miyazaki, A., & Tsuboi, M. 2000, ApJ, 536, 357

Mouschovias, T. C., & Spitzer, L., Jr. 1976, ApJ, 210, 326

Myers, P. C., & Goodman, A. A. 1988, ApJL, 326, L27

Nakajima, Y., Kato, D., Nagata, T., et al. 2005, AJ, 129, 776

Nayak, O., Meixner, M., Fukui, Y., et al. 2018, ApJ, 854, 154

Nayak, O., Meixner, M., Indebetouw, R., et al. 2016, ApJ, 831, 32

Nayak, O., Meixner, M., Sewilo, M., et al. 2019, ApJ, 877, 135

Ochsendorf, B. B., Zinnecker, H., Nayak, O., et al. 2017, NatAs, 1, 784

Olivier, G. M., Lopez, L. A., Rosen, A. L., et al. 2021, ApJ, 908, 68

Pabst, C., Higgins, R., Goicoechea, J. R., et al. 2019, Natur, 565, 618

Pabst, C. H. M., Goicoechea, J. R., Teyssier, D., et al. 2020, A&A, 639, A2

Pehlemann, E., Hofmann, K. H., & Weigelt, G. 1992, A&A, 256, 701

Pellegrini, E. W., Baldwin, J. A., & Ferland, G. J. 2011, ApJ, 738, 34

Phiri, S. P., Kirk, J. M., Ward-Thompson, D., Sansom, A. E., & Bendo, G. J. 2021, MNRAS, 504, 4511

Rathjen, T.-E., Naab, T., Girichidis, P., et al. 2021, MNRAS, 504, 1039

Relaño, M., & Kennicutt, R. C., Jr. 2009, ApJ, 699, 1125

Rieke, G. H., Mather, J. C., Young, E. T., et al. 2004, Proc. SPIE, 5487, 50

Robitaille, T. P. 2017, A&A, 600, A11

Robitaille, T. P., Whitney, B. A., Indebetouw, R., & Wood, K. 2007, ApJS, 169, 328

Robitaille, T. P., Whitney, B. A., Indebetouw, R., Wood, K., & Denzmore, P. 2006, ApJS, 167, 256

Roellig, T. L., Mather, J. C., Werner, M. W., et al. 2004, Proc. SPIE, 5487, 38

Romano, D., Chiappini, C., Matteucci, F., & Tosi, M. 2005, A&A, 430, 491

Rosolowsky, E. W., Pineda, J. E., Kauffmann, J., & Goodman, A. A. 2008, ApJ, 679, 1338

Rubio, M., Roth, M., & Garcia, J. 1992, A&A, 261, L29

Sabbi, E., Anderson, J., Lennon, D. J., et al. 2013, AJ, 146, 53

Salpeter, E. E. 1955, ApJ, 121, 161

Schneider, F. R. N., Ramírez-Agudelo, O. H., Tramper, F., et al. 2018, A&A, 618, A73

Seale, J. P., Looney, L. W., Chu, Y.-H., et al. 2009, ApJ, 699, 150

Seale, J. P., Meixner, M., Sewilo, M., et al. 2014, AJ, 148, 124

Selman, F. J., & Melnick, J. 2013, A&A, 552, A94

Sewilo, M., Carlson, L. R., Seale, J. P., et al. 2013, ApJ, 778, 15

Shetty, R., Beaumont, C. N., Burton, M. G., Kelly, B. C., & Klessen, R. S. 2012, MNRAS, 425, 720

Shetty, R., & Ostriker, E. C. 2008, ApJ, 684, 978

Skinner, M. A., & Ostriker, E. C. 2015, ApJ, 809, 187

Solomon, P. M., Rivolo, A. R., Barrett, J., & Yahil, A. 1987, ApJ, 319, 730

Spergel, D. N., & Blitz, L. 1992, Natur, 357, 665

Sun, J., Leroy, A. K., Schinnerer, E., et al. 2020, ApJL, 901, L8

Utomo, D., Blitz, L., Davis, T., et al. 2015, ApJ, 803, 16

van Dishoeck, E. F., & Blake, G. A. 1998, ARA&A, 36, 317

Walborn, N. R., Barbá, R. H., Brandner, W., et al. 1999, AJ, 117, 225

Walborn, N. R., Barbá, R. H., & Sewilo, M. M. 2013, AJ, 145, 98

Walch, S., Girichidis, P., Naab, T., et al. 2015, MNRAS, 454, 238

Walch, S., Whitworth, A. P., Bisbas, T. G., Wünsch, R., & Hubber, D. A. 2013, MNRAS, 435, 917

Weingartner, J., & Draine, B. T. 2001, ApJ, 548, 296

Whitney, B. A., Sewilo, M., Indebetouw, R., et al. 2008, AJ, 136, 18

Wolfire, M. G., Hollenbach, D., & McKee, C. F. 2010, ApJ, 716, 1191

Wong, T., Hughes, A., Ott, J., et al. 2011, ApJS, 197, 16

Wong, T., Hughes, A., Tokuda, K., et al. 2017, ApJ, 850, 139

Wong, T., Oudshoorn, L., Sofovich, E., et al. 2022, ApJ, 932, 47

Zinnecker, H., & Yorke, H. W. 2007, ARA&A, 45, 481