<?xml-model href='http://www.tei-c.org/release/xml/tei/custom/schema/relaxng/tei_all.rng' schematypens='http://relaxng.org/ns/structure/1.0'?><TEI xmlns="http://www.tei-c.org/ns/1.0">
	<teiHeader>
		<fileDesc>
			<titleStmt><title level='a'>COLDz: Deep 34 GHz Continuum Observations and Free–Free Emission in High-redshift Star-forming Galaxies</title></titleStmt>
			<publicationStmt>
				<publisher></publisher>
				<date>05/01/2021</date>
			</publicationStmt>
			<sourceDesc>
				<bibl> 
					<idno type="par_id">10330881</idno>
					<idno type="doi">10.3847/1538-4357/abe6a5</idno>
					<title level='j'>The Astrophysical Journal</title>
<idno>0004-637X</idno>
<biblScope unit="volume">912</biblScope>
<biblScope unit="issue">1</biblScope>					

					<author>H. S. Algera</author><author>J. A. Hodge</author><author>D. Riechers</author><author>E. J. Murphy</author><author>R. Pavesi</author><author>M. Aravena</author><author>E. Daddi</author><author>R. Decarli</author><author>M. Dickinson</author><author>M. Sargent</author><author>C. E. Sharon</author><author>J. Wagg</author>
				</bibl>
			</sourceDesc>
		</fileDesc>
		<profileDesc>
			<abstract><ab><![CDATA[Abstract                          The high-frequency radio sky has historically remained largely unexplored due to the typical faintness of sources in this regime, and the modest survey speed compared to observations at lower frequencies. However, high-frequency radio surveys offer an invaluable tracer of high-redshift star formation, as they directly target the faint radio free–free emission. We present deep continuum observations at 34 GHz in the COSMOS and GOODS-North fields from the Karl G. Jansky Very Large Array (VLA), as part of the COLD              z              survey. The deep COSMOS mosaic spans                                                                                                            down to              σ              = 1.3              μ              Jy beam              −1              , while the wider GOODS-N observations cover                                                                                                            to              σ              = 5.3              μ              Jy beam              −1              . We detect a total of 18 galaxies at 34 GHz, of which nine show radio emission consistent with being powered by star formation; although for two sources, this is likely due to thermal emission from dust. Utilizing deep ancillary radio data at 1.4, 3, 5, and 10 GHz, we decompose the spectra of the remaining seven star-forming galaxies into their synchrotron and thermal free–free components, and find typical thermal fractions and synchrotron spectral indices comparable to those observed in local star-forming galaxies. We further determine free–free star formation rates (SFRs), and show that these are in agreement with SFRs from spectral energy distribution-fitting and the far-infrared/radio correlation. Our observations place strong constraints on the high-frequency radio emission in typical galaxies at high redshift, and provide some of the first insights into what is set to become a key area of study with future radio facilities, such as the Square Kilometer Array Phase 1 and next-generation VLA.]]></ab></abstract>
		</profileDesc>
	</teiHeader>
	<text><body xmlns="http://www.tei-c.org/ns/1.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xmlns:xlink="http://www.w3.org/1999/xlink">
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="1.">Introduction</head><p>Deep radio observations offer an invaluable view on star formation in the high-redshift universe. With current facilities, such as the upgraded National Science Foundation's (NSF's) Karl G. Jansky Very Large Array (VLA), both star-forming galaxies and faint active galactic nuclei (AGNs) can now be studied down to microjansky flux densities at gigahertz frequencies. However, a regime that remains substantially understudied is the faint radio population at high frequencies (&#957; &#61577; 10 GHz), which is in large part the result of the comparative inefficiency at which high-frequency radio surveys can be executed. First, for a fixed telescope size, the field of view of a single pointing decreases steeply with frequency as &#957; -2 . Second, radio sources are generally intrinsically fainter at high radio frequencies, S &#957; &#8733; &#957; &#945; , where &#945; &#8764; -0.7 <ref type="bibr">(Condon 1992)</ref>, and third, typically only a modest fraction of telescope observing time is suitable for highfrequency observations due to more stringent requirements on the observing conditions. As a result, the survey speed at 34 GHz is &#61577;5000 times smaller compared to observations at the more commonly utilized frequency of 1.4 GHz. Despite these observational difficulties, high-frequency radio observations provide complementary insight into both starforming galaxies and AGNs. Historically, low-frequency radio observations have been used as a tracer of star formation activity through the far-infrared/radio correlation (FIRRC; van der <ref type="bibr">Kruit 1971</ref><ref type="bibr">Kruit , 1973;;</ref><ref type="bibr">de Jong et al. 1985;</ref><ref type="bibr">Helou et al. 1985)</ref>. This correlation, which has been shown to hold over several orders of magnitude in terms of luminosity <ref type="bibr">(Yun et al. 2001;</ref><ref type="bibr">Bell 2003)</ref>, as well as to high redshift <ref type="bibr">(z &#8764; 5;</ref><ref type="bibr">Calistro Rivera et al. 2017;</ref><ref type="bibr">Delhaize et al. 2017;</ref><ref type="bibr">Algera et al. 2020a;</ref><ref type="bibr">Delvecchio et al. 2021)</ref>, relates the predominantly nonthermal synchrotron emission of a star-forming galaxy to its far-infrared (FIR) luminosity. The latter has been well-calibrated as a tracer of star formation, as at FIR-wavelengths, dust re-emits the light absorbed from young, massive stars (e.g., <ref type="bibr">Kennicutt 1998)</ref>. The synchrotron emission, instead, emanates from cosmic rays accelerated by supernova-induced shocks, and as such constitutes a tracer of the end product of massive-star formation <ref type="bibr">(Condon 1992;</ref><ref type="bibr">Bressan et al. 2002)</ref>. However, a second process is expected to dominate the radio spectral energy distribution (SED) at high frequencies (&#957; &#61577;30 GHz): radio free-free emission (FFE). Unlike radio synchrotron radiation, FFE is a much more direct star formation rate (SFR) tracer, as it originates from the H II regions in which massive stars have recently formed. In addition, unlike other commonly used probes of star formation such as UV continuum emission or the hydrogen Balmer lines, FFE constitutes a tracer of star formation that is largely unbiased by dust extinction. These characteristics establish FFE as one of the most reliable tracers of star formation, both in the local and the high-redshift universe.</p><p>Locally, the radio spectra of both individual star-forming regions and star-forming galaxies have been well characterized, and have established FFE as a means of calibrating other tracers of star formation <ref type="bibr">(Murphy et al. 2011</ref>). In addition, in nearby galaxies, FFE can typically be separated spatially from nonthermal synchrotron emission, as individual extragalactic star-forming regions can be resolved <ref type="bibr">(Tabatabaei et al. 2013;</ref><ref type="bibr">Querejeta et al. 2019;</ref><ref type="bibr">Linden et al. 2020)</ref>. However, in the high-redshift universe, FFE has remained elusive, despite the observational advantage that high-frequency continuum emission redshifts into radio bands that are more easily accessible from Earth, facilitating the sampling of the free-free-dominated regime of the radio spectrum. The comparative faintness of high-redshift galaxies, however, complicates the usage of FFE as a tracer of star formation at early cosmic epochs. Indeed, current detections of high-frequency continuum emission in distant galaxies remain limited to bright or gravitationally lensed starbursts <ref type="bibr">(Thomson et al. 2012;</ref><ref type="bibr">Aravena et al. 2013;</ref><ref type="bibr">Riechers et al. 2013;</ref><ref type="bibr">Wagg et al. 2014;</ref><ref type="bibr">Huynh et al. 2017;</ref><ref type="bibr">Penney et al. 2020</ref>). In addition, most of these studies lacked the ancillary low-frequency data required to robustly disentangle FFE from the overall radio continuum, which requires observations at a minimum of three frequencies (e.g., <ref type="bibr">Tabatabaei et al. 2017;</ref><ref type="bibr">Klein et al. 2018)</ref>, or probed rest-frame frequencies dominated by thermal emission from dust (&#957; &#61577; 200 GHz; <ref type="bibr">Condon 1992)</ref>. Despite this observational complexity, one of the key science goals for upcoming radio facilities such as the next-generation VLA is to systematically use FFE as a probe of star formation in the high-redshift galaxy population <ref type="bibr">(Barger et al. 2018)</ref>. As such, it is already of considerable interest to explore this high-frequency parameter space with current radio facilities.</p><p>While radio continuum observations are invaluable in characterizing high-redshift star formation in a dust-unbiased manner, radio surveys are additionally capable of detecting what fuels this process, namely molecular gas. For a clear census of the molecular gas reservoir of the universe, blind surveys are crucial, as they do not suffer from any biases arising from follow-up radio observations of known highredshift sources. The first such blind surveys have recently been completed, such as the ALMA Spectroscopic Survey in the Hubble Ultra-Deep Field (ASPECS; <ref type="bibr">Walter et al. 2016;</ref><ref type="bibr">Decarli et al. 2016)</ref>, and the CO Luminosity Density at High Redshift survey (COLDz; <ref type="bibr">Pavesi et al. 2018;</ref><ref type="bibr">Riechers et al. 2019</ref><ref type="bibr">Riechers et al. , 2020))</ref>, which targets low-J CO observations at 34 GHz using the VLA. Due to the large bandwidth of its highfrequency receivers, deep VLA surveys of molecular gas result in sensitive continuum images essentially "for free." In this work, we describe the deep continuum observations of the COLDz survey. Our main goal is to constrain the radio spectra of typical sources in a frequency range that has not been widely explored, and extend our analysis to a new parameter space of faint AGNs and star-forming galaxies, down to the microjansky level. The COLDz survey covers a region of two well-studied extragalactic fields, COSMOS <ref type="bibr">(Scoville et al. 2007</ref>) and GOODS-North <ref type="bibr">(Giavalisco et al. 2004)</ref>, and hence allows for a multiwavelength perspective on this faint population.</p><p>The outline of this paper is as follows. In Section 2 we describe our 34 GHz VLA observations and the creation of deep continuum images, as well as the existing ancillary multiwavelength data. We outline the detection and source extraction of sources at 34 GHz in Section 3, and assign the radio sample multiwavelength counterparts and redshifts in Section 4. We present our deep radio number counts at 34 GHz in Section 5, and separate the sample into radio AGNs and starforming galaxies in Section 6. In Section 6, we additionally decompose the spectra of the star-forming galaxies into radio synchrotron and FFE, and compare the latter as a tracer of star formation with more commonly adopted tracers at high redshift. Section 7 provides an outlook for the future, and discusses how upcoming radio facilities will revolutionize high-redshift studies of radio star formation. Finally, we summarize our findings in Section 8. Where necessary, we assume a standard &#923;-cold dark matter cosmology, with H 0 = 70 km s -1 Mpc -1 , &#937; m = 0.30, and &#937; &#923; = 0.70. Magnitudes are quoted in the AB-system, and a <ref type="bibr">Chabrier (2003)</ref> initial mass function is assumed.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="2.">Observations and Data Reduction</head></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="2.1.">COLDz</head><p>The COLDz survey <ref type="bibr">(Pavesi et al. 2018;</ref><ref type="bibr">Riechers et al. 2019</ref><ref type="bibr">Riechers et al. , 2020) )</ref> was designed to blindly probe the low-J CO-lines in high-redshift galaxies, which requires high-frequency observations spanning a large bandwidth. Both are provided by the Ka band of the upgraded VLA, which allows for a total of 8 GHz of contiguous frequency coverage, tuneable to be within the range of 26.5-40.0 GHz. The COLDz observing strategy is presented in <ref type="bibr">Pavesi et al. (2018)</ref>, to which we refer the reader for additional details.</p><p>For the optimal constraints on the CO luminosity function, COLDz combines a deep yet small mosaic in the COSMOS field with wider but shallower observations in the GOODS-N field. The COSMOS observations constitute a seven-pointing mosaic, and account for a total on-source time of 93 hr. The data were taken in the VLA D (82 hr) and DnC (11 hr) configurations, and span a total frequency range of 30.969-39.033 GHz, using dual polarization. The observations of the GOODS-N field comprise a total of 122 hr of on-source time, and make up a 57-pointing mosaic. The data were predominantly taken in the D-configuration (83.5 hr), with additional observations taking place in the D &#8594; DnC (4.4 hr), DnC (30.8 hr), and DnC &#8594; C (3 hr) configurations. The observations cover a frequency range between 29.981-38.001 GHz, and in what follows, we will refer to both the COSMOS and GOODS-N observations by their typical central frequency of 34 GHz.</p><p>Calibration of the observations was done in CASA version 4.1.0, with extensive use being made of a modified version of the VLA pipeline, the details of which can be found in <ref type="bibr">Pavesi et al. (2018)</ref>. As one of the main goals of COLDz is to detect spectral lines, both Hanning smoothing and RFLAG, used to remove radio frequency interference (RFI), were switched off. Instead, some persistent RFI at 31.5 GHz was flagged manually, and some occasions of narrow noise spikes were flagged via the methods detailed in <ref type="bibr">Pavesi et al. (2018)</ref>. The data were then concatenated, and, for the continuum observations presented here, subsequently averaged in time (9 s) and frequency (16 channels), in order to reduce the size of the data set prior to imaging. The imaging of the calibrated 34 GHz observations was carried out in CASA version 4.3.1, using the "mosaic" mode of CASA task CLEAN. For both the COSMOS and the GOODS-N mosaics, a multifrequency synthesis algorithm was employed to take into account the large bandwidth of the observations. A natural weighting was further adopted, in order to maximize the sensitivity of the data. The data were imaged iteratively, by cleaning all sources at &gt; 6&#963; down to the 2&#963; level.</p><p>We present the 34 GHz continuum maps across the COSMOS and GOODS-N fields, as well as the corresponding rms noise maps, in Figures <ref type="figure">1</ref> and<ref type="figure">2</ref>. The COSMOS mosaic covers a field of view of 9.6 arcmin 2 out to 20% of the peak primary beam sensitivity. The central rms noise in the image equals 1.3 &#956;Jy beam -1 , with the typical rms increasing to 1.5 and 1.9 &#956;Jy beam -1 , within 50% and 20% of the peak primary beam sensitivity, respectively. The synthesized beam of the COSMOS observations is well described by an elliptical Gaussian of 2 70 &#215; 2 41 with a position angle of -0&#176;.7.</p><p>The GOODS-N mosaic spans an area of 51 arcmin 2 out to 20% of the peak primary beam sensitivity, covering roughly  one-third of the "traditional" GOODS-N survey <ref type="bibr">(Giavalisco et al. 2004)</ref>. The typical noise level varies slightly across the mosaic, with a single, deep 34 GHz pointing in the field reaching a central noise level of 3.2 &#956;Jy beam -1 . Across the entire GOODS-N field, within 50% and 20% of the peak primary beam sensitivity, respectively, the typical rms noise equals 5.3 &#956;Jy beam -1 and 5.5 &#956;Jy beam -1 . The synthesized beam is well described by an elliptical Gaussian of 2 19 &#215; 1 84, with a position angle of 75&#176;.3, after smoothing the different pointings to a common beam. <ref type="bibr">Pavesi et al. (2018)</ref> also provide a version of the GOODS-N radio map where all pointings have been imaged at their native resolution. While the resulting mosaic cannot be described by a single beam, it allows for the search for unresolved sources at slightly higher signal-to-noise ratios (S/Ns), as no smoothing or tapering was required. We use this "unsmoothed" map for source detection in Section 3, in addition to the regular mosaic with a homogenized beam.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="2.2.">Ancillary Radio Observations</head><p>The COLDz observations in the COSMOS field overlap in their entirety with deep 3 and 10 GHz VLA pointings from the COSMOS-XS survey <ref type="bibr">(Algera et al. 2020b;</ref><ref type="bibr">van der Vlugt et al. 2021)</ref>. These images have a resolution similar to COLDz, with a synthesized beam of 2 14 &#215; 1 81 2 14 &#215; 1 81 and 2 33 &#215; 2 01 at 3 and 10 GHz, respectively. At 3 GHz, the observations reach a central rms of 0.53 &#956;Jy beam -1 , with a typical primary beam sensitivity of 90% of the maximum across the COLDz mosaic. The 10 GHz data reach a central rms sensitivity of 0.41 &#956;Jy beam -1 (typical primary beam sensitivity of 80%), and were centered on the COLDz observations by design. In total, 70 sources detected at 3 GHz at 5&#963; fall within the ~10 arcmin 2 COLDz field of view. A subset of 40 are additionally detected at 10 GHz. In addition, the central region of the COSMOS field, spanning &#8776;1 deg 2 , is covered by VLA observations at 1.4 GHz <ref type="bibr">(Schinnerer et al. 2007</ref><ref type="bibr">(Schinnerer et al. , 2010))</ref>. These observations reach a typical rms sensitivity of 12 &#956;Jy beam -1 across the COLDz field of view, at a resolution of 1 5 &#215; 1 4. In total, two sources within the COLDz field of view are detected at 1.4 GHz. Adopting a typical spectral index of &#945; = -0.70, the COLDz and 1.4 GHz COSMOS observations are approximately of equal relative depth. <ref type="foot">12</ref>The GOODS-N field is similarly well covered by radio observations from the VLA. At 1.4 GHz, the entire field has been imaged in a single pointing by <ref type="bibr">Owen (2018)</ref>, down to a central rms of 2.2 &#956;Jy beam -1 with little variation across the field of view of the COLDz mosaic. The resolution of the radio data equals 1 6, and a total of 186 sources detected at 1.4 GHz by <ref type="bibr">Owen (2018)</ref> fall within the 20% power point of the GOODS-N COLDz image. The field has additionally been mapped at 5 GHz by <ref type="bibr">Gim et al. (2019)</ref>. They detect 52 sources down to an rms of 3.5 &#956;Jy beam -1 across two VLA pointings covering a total area of 109 arcmin 2 , which fully overlaps with the COLDz footprint. Their angular resolution of 1 47 &#215; 1 42 is similar to that of the deeper 1.4 GHz observations, and <ref type="bibr">Gim et al. (2019)</ref> find that all 5 GHz detections have a counterpart in this lower-frequency map. Finally, <ref type="bibr">Murphy et al. (2017)</ref> present a single pointing at 10 GHz in the GOODS-N field (primary beam FWHM of &#162; 4.25) at a native resolution of 0 22. At this high angular resolution, their observations reach a sensitivity of 0.57 &#956;Jy beam -1 in the center of the pointing. <ref type="bibr">Murphy et al. (2017)</ref> further provide two tapered images with a resolution of 1&#8243; and 2&#8243;, similar to the ancillary radio data in the field. These images reach an rms noise of 1.1 and 1.5 &#956;Jy beam -1 , respectively. In total, <ref type="bibr">Murphy et al. (2017)</ref> recover 38 sources across the combined high-resolution and tapered images. Approximately 75% of the COLDz GOODS-N data overlap with the smaller pointing at 10 GHz, when imaged out to 5% of the primary beam FWHM <ref type="bibr">(Murphy et al. 2017)</ref>.</p><p>We summarize the detection limits of the various radio observations across both COSMOS and GOODS-N in Figure <ref type="figure">3</ref>, and compare these with the typical radio spectrum of a star-forming galaxy (SFR = 100 M e yr -1 at z = 1). This assumes the <ref type="bibr">Kennicutt (1998)</ref> conversion between SFR and infrared luminosity, adapted for a Chabrier IMF, as well as a simple optically thin dust SED with &#946; = 1.8 and T dust = 35 K. We further adopt the <ref type="bibr">Condon (1992)</ref> model for the radio spectrum of star-forming galaxies, and assume the FIRRC from <ref type="bibr">Delhaize et al. (2017)</ref>. Under these assumptions, a galaxy with SFR = 100 M e yr -1 can be directly detected at 34 GHz out to z = 2 (z = 1) in the COSMOS (GOODS-N) field.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="2.3.">Ultraviolet to Submillimeter Observations</head><p>Both the COSMOS and GOODS-North fields have been targeted by a wealth of multiwavelength observations, spanning the full X-ray to radio regime. For the COSMOS field, we adopt the multiwavelength matching procedure for the COSMOS-XS survey <ref type="bibr">(Algera et al. 2020b)</ref>, as all 34 GHz continuum detections have radio counterparts in this survey (Section 3). This cross-matching procedure invokes the recent "Super-deblended" catalog from <ref type="bibr">Jin et al. (2018)</ref>, who adopt a novel deblending technique to address confused mid-infrared to submillimeter observations. The Super-deblended catalog provides photometry from Spitzer/IRAC 3.6 &#956;m to MAMBO 1.2 mm (for a full list of references, see <ref type="bibr">Jin et al. 2018)</ref>, as well as radio data at 1.4 and 3 GHz from <ref type="bibr">Schinnerer et al. (2010)</ref> and <ref type="bibr">Smol&#269;i&#263; et al. (2017)</ref>, respectively. However, we directly Figure <ref type="figure">3</ref>. The detection limits of the various radio observations across COSMOS and GOODS-N, superimposed on the radio spectrum of a starforming galaxy (SFR = 100 M e yr -1 at z = 1; see the text for details). The blue (green) lines indicate the detection limit of the COLDz 34 GHz continuum data, scaled with a typical &#945; = -0.70. The different constituents of the radio spectrum are shown, including free-free emission, which is expected to become dominant at rest-frame &#957; &#61577;30 GHz. As such, the COLDz observations directly target the strongly free-free-dominated regime.</p><p>adopt the radio fluxes from the 1.4 GHz catalog from <ref type="bibr">Schinnerer et al. (2007</ref><ref type="bibr">Schinnerer et al. ( , 2010))</ref>, and use the deeper COSMOS-XS 3 GHz data in favor of the observations from <ref type="bibr">Smol&#269;i&#263; et al. (2017)</ref>. In addition, <ref type="bibr">Algera et al. (2020b)</ref> cross-match with the z ++ YJHK s -selected COSMOS2015 catalog from <ref type="bibr">Laigle et al. (2016)</ref> containing far-ultraviolet to near-infrared photometry, in order to complete the coverage of the SED. Finally, we search for Atacama Large Millimeter/submillimeter Array (ALMA) counterparts as part of the AS2COSMOS <ref type="bibr">(Simpson et al. 2020</ref>) and A3COSMOS surveys <ref type="bibr">(Liu et al. 2019)</ref>, which constitute a collection of individual pointings across the COSMOS field.</p><p>A wealth of multiwavelength data similarly exist across the GOODS-N field. In order to obtain optical/near-infrared photometry for the 34 GHz continuum detections, we adopt the 3D-HST photometric catalog from <ref type="bibr">Skelton et al. (2014)</ref>. Source detection was performed in a deep combined F125W+F140W +F160W image, and further photometry is carried out in the wavelength range of 0.3-8 &#956;m, including Spitzer/IRAC observations from <ref type="bibr">Dickinson et al. (2003)</ref> and <ref type="bibr">Ashby et al. (2013)</ref>. For further details, we refer the reader to <ref type="bibr">Skelton et al. (2014)</ref>. We obtain additional mid-and far-infrared photometry from the Super-deblended catalog across the GOODS-N field by <ref type="bibr">Liu et al. (2018)</ref>. They adopt combined priors from Spitzer/IRAC, Spitzer/ MIPS 24 &#956;m and VLA 1.4 GHz observations, and utilize these to deblend the photometry at more strongly confused wavelengths. The Super-deblended catalog provides additional photometry from the Spitzer/Infrared Spectropgraph (IRS), Herschel/Photodetector Array Camera and Spectrometer (PACS; <ref type="bibr">Magnelli et al. 2013)</ref> and Herschel/Spectral and Photometric Imaging REceiver (SPIRE; <ref type="bibr">Elbaz et al. 2011)</ref>, as well as data from JCMT/SCUBA-2 850 &#956;m <ref type="bibr">(Geach et al. 2017)</ref> and AzTEC+MAMBO 1.2 mm <ref type="bibr">(Penner et al. 2011</ref>).</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="2.4.">X-Ray Observations</head><p>While a radio-based selection renders one sensitive to AGN activity at radio wavelengths, &#8764; 20%-25% of the faint radio population (S 1.4 &#61576;1 mJy; equivalent to S 34 &#61576;100 &#956;Jy given &#945; = -0.70) is thought to consist of radio-quiet AGNs <ref type="bibr">(Bonzini et al. 2013;</ref><ref type="bibr">Smol&#269;i&#263; et al. 2017)</ref>. These sources show no substantial AGN-related emission at radio wavelengths, but are instead classified as AGNs based on signatures at other wavelengths. Strong X-ray emission, in particular, forms an unmistakable manifestation of AGN activity. As such, we make use of the deep X-ray coverage over both the COSMOS and GOODS-N fields to characterize the nature of the 34 GHz continuum detections.</p><p>The COSMOS field is covered in its entirety by the 4.6 Ms Chandra COSMOS Legacy survey <ref type="bibr">(Civano et al. 2016)</ref>, with the individual Chandra pointings accounting for a typical &#8776;160 ks of exposure time. The area covered by the COLDz survey contains three X-ray detections, at a typical detection limit of the survey of &#8764; 2 &#215; 10 -15 erg cm -2 s -1 in the full range of 2-10 keV. The catalog provided by <ref type="bibr">Marchesi et al. (2016)</ref> further includes X-ray luminosities for X-ray detections with robust optical and infrared counterparts.</p><p>The GOODS-North field is similarly covered by deep 2 Ms Chandra observations as part of the Chandra Deep Field North Survey <ref type="bibr">(Xue et al. 2016)</ref>. Across the COLDz field of view, the survey identifies 189 X-ray sources, and attains a flux limit of &#8764; 3.5 &#215; 10 -17 erg cm -2 s -1 in the 0.5-7 keV energy range -nearly a factor of 50 deeper than the COSMOS data, when converted to the same energy range adopting &#915; = 1.8. The catalog provided by <ref type="bibr">Xue et al. (2016)</ref> further includes X-ray luminosities for the entries with reliable redshift information.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="3.">Continuum Sources</head><p>We run source detection on our 34 GHz radio maps using PYBDSF <ref type="bibr">(Mohan &amp; Rafferty 2015)</ref>, prior to correcting the images for the primary beam. This has the benefit that the noise properties are uniform across the mosaics, which facilitates source detection, and additionally ensures that fewer spurious sources arise around the noisy edges of the maps. <ref type="foot">13</ref> In our source detection procedure, we can afford to set a liberal detection threshold, as both the COSMOS and GOODS-N fields contain additional low-frequency radio data of greater depth relative to the COLDz observations, and as such we expect any 34 GHz detections to have radio counterparts. In this work, we therefore adopt a 3&#963; peak detection threshold for both images.</p><p>In the COSMOS field, we compare to the deep S-and X-band observations from van der <ref type="bibr">Vlugt et al. (2021)</ref>. For a source to be undetected in this 10 GHz map, yet detected at 34 GHz, requires a highly inverted spectral index of a &#61577; 0.55 34 10 . In the GOODS-N field, we compare with the 1.4 GHz observations from Owen (2018), and find that 34 GHz sources require an inverted spectral index of a &#61577; 0.  to remain undetected in the lower-frequency map. We note that, as we perform source detection at low S/Ns, the source sizes calculated via PYBDSF may be affected by the local noise properties within the images. In particular, sources with fitted radio sizes smaller than the synthesized beam will have an integrated flux density smaller than the peak value, and will have an uncertainty on the latter that is smaller than the local rms noise in the map <ref type="bibr">(Condon 1997)</ref>. We show in Appendix A.2 that all but one source are likely to be unresolved at &#8764;2 5 resolution. For all unresolved sources, we adopt the peak brightness, while the integrated flux density is used otherwise. Following the discussion above, we redefine the uncertainty on the peak brightness and conservatively adopt the maximum of the calculated uncertainty from PYBDSF and the local rms at the source position.</p><p>In the COSMOS map, we detect 57 peaks at 3&#963; within 20% of the peak primary beam sensitivity (Figure <ref type="figure">4</ref>). We can match six to COSMOS-XS counterparts within 0 7, where we expect to have N false = 0.2 false matches based on randomly shifting the coordinates of all 3&#963; peaks in the mosaic, and repeating the matching a large number of times. This radius was chosen to include all close associations to COSMOS-XS sources, while minimizing the number of expected false matches (N false = 1). To further verify the robustness of the six close counterparts, we run our source detection procedure on the inverted radio map (i.e., multiplied by -1), and find a total of 60 negative "sources," all of which are by definition spurious. None of these can be matched to COSMOS-XS counterparts within 0 7, indicating the real matches are likely to be robust. However, out of the six associations to COSMOS-XS galaxies within 0 7, we find that one candidate source (S/N = 3.3 at 34 GHz) is detected solely at 3 GHz while a 10 GHz counterpart is also expected, implying it is likely to be spurious. As such, we discard it from our sample and retain five sources that form the robust COSMOS continuum sample.</p><p>Three of these have an S/N larger than the highest peak S/N in the inverted radio map of 3.7&#963;. The remaining two have a relatively low S/N of &#8764;3.5, but are deemed robust due to their multiwavelength associations.</p><p>We adopt an identical source detection procedure for the GOODS-N 34 GHz map. In total, PYBDSF identifies 236 peaks above 3&#963; in the map, of which we match 12 with lowfrequency radio counterparts in the catalog from Owen (2018) at 0 7, where we expect N false = 0.3 incorrect identifications. We find 263 sources in the inverted radio map, none of which are matched to lower-frequency radio counterparts within 0 7. However, upon cross-matching our 12 robust sources with the 5 GHz catalog from <ref type="bibr">Gim et al. (2019)</ref>, we find only 11 matches within 0 7. The single unmatched source (S/N = 3.3 at 34 GHz) also falls within the field of view of the single deep X-band pointing from <ref type="bibr">Murphy et al. (2017)</ref>, but is additionally undetected at 10 GHz. As such, this source is likely to be spurious, and we discard it from further analysis. The S/N distributions of all peaks identified by PYBDSF in the GOODS-N image and its inverted counterpart are shown in Figure <ref type="figure">5</ref>. The maximum S/N in the inverted map is S/N&#8776;5.0, which implies that 5/11 robust 34 GHz detections lie below the most significant spurious source. This emphasizes the added value of the deep, low-frequency data, which allow us to identify faint sources that would not have been recovered in a blind source detection procedure. We further perform source detection via PYBDSF on the unsmoothed mosaic (Section 2), using the same detection threshold as adopted for the regular GOODS-N mosaic, in order to maximize the number of recovered sources. We find two additional matches with low-frequency counterparts at both 1.4 and 5 GHz, resulting in a total of 13 sources identified in GOODS-N. We adopt the peak brightness for these two sources, given the lack of a common beam in the unsmoothed mosaic, and use the local rms in the map as the appropriate uncertainty.</p><p>Given that we require any identifications at 34 GHz to have robust low-frequency radio counterparts, we may be missing sources with highly unusual inverted spectra. To investigate this possibility, we median-stack in the low-frequency radio maps on the positions of the highest signal-to-noise peaks at 34 GHz for which no radio counterpart was found (3 and 10 GHz in COSMOS, 1.4 GHz in GOODS-N, using the publicly available radio map from <ref type="bibr">Morrison et al. 2010)</ref>. In neither the COSMOS nor the GOODS-N fields do we see any evidence for a positive signal in the stacks, indicating that most-and likely all-of these S/N&#8776;3-5 peaks at 34 GHz are spurious. In addition, we crossmatch the peaks in the 34 GHz map that do not have radio counterparts with optical/near-infrared-selected sources from the COSMOS2015 and 3D-HST catalogs in COSMOS and GOODS-N, respectively, within 0 7. We find two matches with the COSMOS2015 catalog, at &#8764;0 4. However, based on a visual inspection, there is no hint of emission at either 3 or 10 GHz for these sources, indicating that they are likely to be spurious. We find 37 matches within 0 7 between 3&#963; peaks in the 34 GHz GOODS-N radio map without radio counterparts and the 3D-HST catalog. In addition to a visual inspection, we stack these sources at  1.4 GHz, but find no detection in the stack, and place a 3&#963; upper limit on the average 1.4 GHz emission of S 1.4 2.7 &#956;Jy beam -1 , which corresponds to a highly inverted a &#61577; 0.6 34 1.4</p><p>. This further substantiates that the majority of the low-S/N peaks identified at 34 GHz are likely to be spurious.</p><p>We summarize the radio properties of the 18 sources detected at 34 GHz in Table <ref type="table">1</ref>, and show postage stamps on top of HST images in Figure <ref type="figure">6</ref>. The tabulated flux densities are further corrected for flux boosting, as outlined in Appendix A.3. For a source detected at 3&#963; (4&#963;), the typical correction factor is 20% (10%), whereas at S/N &#61577;5, the effects of flux boosting are found to be negligible.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="4.">Multiwavelength Properties</head></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="4.1.">Multiwavelength Counterparts and Redshifts</head><p>Across the combined COSMOS and GOODS-N 34 GHz mosaics, we identify a total of 18 robust high-frequency continuum detections. In this section, we detail the association of multiwavelength counterparts to this radio-selected sample. For the five COSMOS sources, we follow the cross-matching procedure from <ref type="bibr">Algera et al. (2020b)</ref>. We first match the 34 GHz continuum detections to the Super-deblended catalog containing FIR photometry, exploring matching radii up to 0 9. However, we find that all five sources can be matched to Super-deblended counterparts within 0 3, where we expect a negligible number of false associations. We further find that the five sources can be additionally cross-matched to galaxies in the COSMOS2015 catalog, similarly within 0 3. We subsequently cross-match with the AS2COSMOS and A3COSMOS ALMA catalogs, finding a single match at 0 1 that appears in both catalogs, and adopt the photometry from the former. We additionally match with the robust and tentative catalog of COLDz CO-emitters from <ref type="bibr">Pavesi et al. (2018)</ref> and recover two matches within 0 3 within the robust set of blind COdetections. <ref type="foot">14</ref>We then extract the optimal photometric or spectroscopic redshift from these catalogs. We prioritize redshifts in the following order: (1) a spectroscopic value based on a COLDz CO-line, (2) a spectroscopic redshift within the Superdeblended catalog, and (3) a photometric redshift within the COSMOS2015 catalog. In total, we find that four out of five COLDz COSMOS detections have a spectroscopically confirmed redshift, while the remaining source has a wellconstrained photometric redshift measurement (COLDz-cont-COS-3 at z = 0.98 &#177;0.01). We note that one of our 34 GHz continuum detections (COLDz-cont-COS-2) is the well-studied submillimeter galaxy AzTEC.3 at z = 5.3, which is additionally detected in CO-emission in the COLDz survey. For AzTEC.3, we further compile additional available ALMA continuum photometry at 230 and 300 GHz from <ref type="bibr">Pavesi et al. (2016)</ref>.</p><p>For the GOODS-N field, we follow a similar procedure. We first cross-match the COLDz continuum detections with the 3D-HST survey <ref type="bibr">(Brammer et al. 2012;</ref><ref type="bibr">Skelton et al. 2014;</ref><ref type="bibr">Momcheva et al. 2016)</ref>, adopting the radio positions at 1.4 GHz from <ref type="bibr">Owen (2018)</ref>, as these are of higher signal-to-noise than the 34 GHz data, and as such are less susceptible to the local noise properties. We find that all 13 radio sources have counterparts in the 3D-HST survey within 0 3, where no spurious matches are expected. We subsequently cross-match with the Super-deblended catalog, and find that all but one of the 34 GHz continuum sources have a counterpart at FIR wavelengths. All 12 cross-matches have a Super-deblended counterpart within 0 3, while there are no further matches within 3 0 of the single unmatched entry. We further find three cross-matches within 0 1 with the catalog of submillimeterselected galaxies from <ref type="bibr">Pope et al. (2005)</ref>, which we identify with GN12, GN16, and GN26, the latter of which has a spectroscopically measured redshift of z = 1.22 based on a CO(2-1) detection from <ref type="bibr">Frayer et al. (2008)</ref>. Additional matching with the COLDz catalog of line emitters in the GOODS-N field does not result in further matches. As such, the z = 5.3 submillimeter galaxy GN10 <ref type="bibr">(Daddi et al. 2009;</ref><ref type="bibr">Riechers et al. 2020)</ref>, detected in COLDz as the brightest CO-emitter <ref type="bibr">(Pavesi et al. 2018)</ref>, remains undetected in deep 34 GHz imaging. Upon performing photometry at the known position of GN10, we determine a peak brightness of S 34 = 11.0 &#177;5.6 &#956;Jy beam -1 (&#8776;2&#963;), placing it well below our survey detection limit. We compile the optimal redshifts for the 13 GOODS-N 34 GHz continuum sources in a similar manner as for the COSMOS field. In order of priority, we adopt a spectroscopic redshift from the Super-deblended catalog, or the best redshift from 3D-HST <ref type="bibr">(Momcheva et al. 2016)</ref>, the latter being either a spectroscopic redshift from HST grism, or a photometric redshift from <ref type="bibr">Skelton et al. (2014)</ref>. We further overwrite a single spectroscopic value for COLDz-cont-GN-10 at z = 4.424 with the updated value of z = 2.018, based on the discussion in <ref type="bibr">Murphy et al. (2017)</ref>. Overall, 10/13 sources have a spectroscopically confirmed redshift, while the remaining three sources have a well-constrained photometric redshift.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="4.2.">Spectral Energy Distributions</head><p>We use SED-fitting code MAGPHYS (da Cunha et al.  2008, 2015) to determine the physical properties of the 34 GHz continuum detections. MAGPHYS adopts an energy balance technique to couple the stellar emission at ultraviolet to near-infrared wavelengths to thermal dust emission at longer wavelengths, and as such, it models the ultraviolet to FIR SED in a self-consistent way. This is additionally useful for sources lacking FIR photometry; as in this case, shorter wavelengths will still provide constraints on the total dust emission. While MAGPHYS is also capable of modeling the radio spectrum of star-forming galaxies, we do not utilize any observations at wavelengths beyond 1.3 mm in our SED-fitting procedure, as any emission from an AGN at radio wavelengths is not incorporated in the fitting. We additionally add a 10% uncertainty in quadrature to the cataloged flux density uncertainties blueward of Spitzer/MIPS 24 &#956;m, following, e.g., <ref type="bibr">Battisti et al. (2019)</ref>. This accounts for uncertainties in the photometric zero-points, and further serves to guide the fitting into better constraining the FIR part of the SED. Via MAGPHYS, we obtain several physical properties for our galaxies, including FIR luminosities, SFRs, and stellar masses. We show the fitted SEDs for all 18 COLDz continuum detections in Figure <ref type="figure">17</ref> in Appendix B, and tabulate their physical parameters in Table <ref type="table">4</ref>.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="4.3.">X-Ray and Mid-infrared AGN Signatures</head><p>Given that the COLDz continuum detections in the COSMOS field are all identified in the COSMOS-XS survey, we adopt the results from <ref type="bibr">Algera et al. (2020b)</ref>, who match the COSMOS-XS 3 GHz continuum sources with the X-ray catalog from <ref type="bibr">Marchesi et al. (2016)</ref>. None of the five COLDz COSMOS continuum sources, however, have counterparts in X-ray emission within a separation of 1 4, and as such, we calculate upper limits on their X-ray luminosity adopting &#915; = 1.8. However, as the resulting limits are of modest depth, we cannot state definitively whether the COLDz COSMOS sources are X-ray AGNs. We note, however, that all of them fall below the typical X-ray emission seen in submillimeter galaxies <ref type="bibr">(Alexander et al. 2005)</ref>, which in turn are thought to predominantly be star formation dominated. We additionally adopt the criteria from <ref type="bibr">Donley et al. (2012)</ref> in order to identify AGNs at z 2.7 through the mid-infrared signature ascribed to a dusty torus surrounding the accreting black hole, and similarly find no signature of AGN-related emission at midinfrared wavelengths. We limit ourselves to this redshift range, as at higher redshifts, the <ref type="bibr">Donley et al. (2012)</ref> criteria are susceptible to false positives in the form of dusty star-forming galaxies (e.g., <ref type="bibr">Stach et al. 2019)</ref>.</p><p>For the GOODS-N field, we find 11/13 matches with the X-ray catalog from <ref type="bibr">Xue et al. (2016)</ref>, within a separation of 1 0. At this matching radius, no false identifications are expected. We compare the X-ray luminosities with the emission expected from star formation, adopting the relations from <ref type="bibr">Symeonidis et al. (2014)</ref>. This comparison classifies nine sources as AGNs based on their X-ray emission, despite two of these sources having [0.5-8] keV X-ray luminosities below the typical threshold for AGNs of L X = 10 42 erg s -1 . We additionally find that only a single X-ray detected source exhibits midinfrared colors that place it within the <ref type="bibr">Donley et al. (2012)</ref> wedge. Overall, we conclude that the majority of 34 GHz continuum detections in the GOODS-N field are X-ray AGNs. Due to the COSMOS X-ray data being comparatively shallower, we cannot assess definitely whether the five 34 GHz continuum sources in the COSMOS field are similarly AGNs based on their X-ray emission. However, a lower AGN fraction in COSMOS is expected, as its radio observations are significantly deeper than in GOODS-N, and the incidence rate of AGNs is a strong function of radio flux density <ref type="bibr">(Smol&#269;i&#263; et al. 2017;</ref><ref type="bibr">Algera et al. 2020b)</ref>. Regardless, if X-ray AGNs are present in the COLDz/COSMOS detections, they are unlikely to be dominating the SED.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="5.">34 GHz Source Counts</head><p>Radio number counts, while historically used as a probe for the cosmology of the universe, remain a useful tool for comparing surveys, in addition to visualizing the onset of different radio populations. At the bright end (S 1.4 ? 1 mJy), the Euclidean number counts decline smoothly toward lower flux densities, and are dominated by luminous radio AGNs (e.g., <ref type="bibr">Condon &amp; Mitchell 1984)</ref>. At 1.4 GHz, the number counts show a flattening at S 1.4 &#8776;1 mJy, believed to be the advent of star-forming galaxies and radio-quiet AGNs as the dominant radio populations (e.g., <ref type="bibr">Rowan-Robinson et al. 1993;</ref><ref type="bibr">Seymour et al. 2004;</ref><ref type="bibr">Padovani et al. 2009;</ref><ref type="bibr">Smol&#269;i&#263; et al. 2017)</ref>. For a typical spectral index of &#945; = -0.70, this flattening should arise around S 34 &#8776;100 &#956;Jy, and hence is covered in the range of flux densities probed in this work. We note that these flux densities may be "contaminated" by thermal emission from dust, or in the case of two COSMOS sources, by bright COemission, which are typically not an issue in low-frequency radio source counts. However, while we correct for this when examining the radio spectra of the star-forming COLDz continuum detections in detail (Section 6.2), here we compute the number counts based on the raw observed flux densities.</p><p>As we adopt deep radio observations as prior positions for possible 34 GHz continuum sources, we expect our sample to be fully reliable, i.e., not to contain any spurious detections. However, our sample may still be incomplete, in particular as a result of the rms of our radio maps increasing rapidly toward the image edges due to the enhanced primary beam attenuation (Figure <ref type="figure">2</ref>). In turn, in these regions we may miss faint 34 GHz continuum sources that would have been observed had the rms been constant across our field of view. As such, we adopt the fractional incompleteness as a function of radio flux density as determined from inserting mock sources into our radio maps (Appendix A.1). Denoting the completeness at flux density S &#957; for field i as f i (S &#957; ), the correction factor per field is simply</p><p>For multiple fields, we then adopt the full completeness to be</p><p>where &#937; i (S &#957; ) denotes the area in arcmin 2 across which a source of flux density S &#957; can be detected in field i at 3&#963; significance.</p><p>As such, the overall completeness is the area-weighted average of the completeness in the COSMOS and GOODS-N fields.</p><p>We present the completeness-corrected Euclidean-normalized number counts at 34 GHz in Figure <ref type="figure">7</ref>, and tabulate the results in Table <ref type="table">2</ref>. The number counts combine the continuum detections across both the COSMOS and GOODS-N fields, using the same area-weighting that is adopted for the completeness calculation (Equation ( <ref type="formula">1</ref>)). The uncertainties on the individual points constitute the combination of the error on the counting statistics from Gehrels (1986)-which are more appropriate than simple Poissonian errors for bins with few sources-and the error on the completeness. Cosmic variance is not included in the uncertainties, but its magnitude is discussed below.</p><p>While the high-frequency (&#957; &#61577;30 GHz) radio sky has been explored at the millijansky level (e.g., <ref type="bibr">Mason et al. 2009)</ref>, the COLDz survey provides the first constraints on the 34 GHz number counts in the regime where star-forming galaxies are expected to emerge as the dominant population, complicating any direct comparisons to the literature. At the bright end of our 34 GHz observations, we compare with <ref type="bibr">Murphy &amp; Chary (2018)</ref>, who perform a stacking analysis in Planck observations at 28.5 GHz, based on priors at 1.4 GHz from the NRAO VLA Sky Survey (NVSS; <ref type="bibr">Condon et al. 1998)</ref>. We rescale the number counts from 28.5 to 34 GHz adopting &#945; = -0.70, and find them to be in agreement with the COLDz observations at S 34 &#61577;200 &#956;Jy.</p><p>As no 34 GHz observations with a similar sensitivity to COLDz exist in the literature, the next best solution is to compare to the highest-frequency number counts available, and scale the counts to 34 GHz. For this, we adopt the recent 10 GHz number counts from van der Vlugt et al. (2021), as part of the COSMOS-XS survey, which we scale to 34 GHz via &#945; = -0.70. As the parent survey constitutes a very deep single pointing ( &#8764;0.40 &#956;Jy across ~30 arcmin 2 ), these data probe down to slightly fainter flux densities than those probed in this work. However, within the flux density range in common, we find that the number counts are in good agreement, despite the inherent uncertainties associated with the required frequency scaling.</p><p>While observationally little is known about the faint highfrequency radio sky, this is no less true for simulations. Modeling of the radio sky has predominantly been performed at low frequencies, with simulations by <ref type="bibr">Wilman et al. (2008)</ref> and more recently by <ref type="bibr">Bonaldi et al. (2019)</ref> only extending to 18 and 20 GHz, respectively. As such, we again invoke a frequency scaling to 34 GHz, adopting as before a spectral index of &#945; = -0.70. We focus on the recent 20 GHz simulations by <ref type="bibr">Bonaldi et al. (2019)</ref>, who model the radio population for two distinct classes of sources: star-forming galaxies that follow the FIRRC, and radio AGNs that show a strong excess in radio power compared to this correlation. We show both the individual and combined contributions of the two populations in Figure <ref type="figure">7</ref>, and find that the simulations predict that below S 34 &#61576;100 &#956;Jy, star-forming galaxies should make up the bulk of the radio population. Our measurements are in good agreement with the simulated number counts, indicating that a scaling from 20 to 34 GHz with a fixed spectral index is likely to be appropriate. We caution, however, that the uncertainties on our number counts are large as a result of the small number of sources detected in the 34 GHz continuum maps.</p><p>Furthermore, cosmic variance constitutes an additional uncertainty on our number counts, as our radio observations probe a relatively small field of view. We show the area probed  as a function of flux density in the bottom panel of Figure <ref type="figure">7</ref>.</p><p>Out to S 34 &#61576;15 &#956;Jy, the total area is dominated by the deeper COSMOS mosaic, and accounts for approximately 10 arcmin 2 . At flux densities above S 34 &#61577;30 &#956;Jy, our field of view increases to approximately 60 arcmin 2 , as now sources can be detected across the full GOODS-N mosaic as well. We quantify the magnitude of cosmic variance via the <ref type="bibr">Bonaldi et al. (2019)</ref> simulations, following <ref type="bibr">Algera et al. (2020b)</ref>. Briefly, for each of the flux density bins adopted in our computation of the 34 GHz number counts, we determine the effective area of the COLDz survey in which sources in the given bin can be detected. We then sample 200 independent circular regions of equivalent area from the 20 GHz <ref type="bibr">Bonaldi et al. (2019)</ref> simulations, accounting for the flux density scaling with &#945; = -0.70. We determine the number counts for all independent areas, and compute the 16-84th and 5-95th percentiles, which we adopt to be the 1&#963; and 2&#963; uncertainties due to cosmic variance, shown as the dark and light gray regions in Figure <ref type="figure">7</ref>, respectively. We note that such a calculation of the cosmic variance encapsulates two effects: at low flux densities, we have a relatively small field of view, which naturally increases the magnitude of cosmic variance. At large flux densities, our field of view constitutes the full ~60 arcmin 2 , but due to the relative paucity of bright radio sources, cosmic variance similarly constitutes an appreciable uncertainty. The typical magnitude of cosmic variance between 5&#61576;S 34 &#61576;60 &#956;Jy induces an additional uncertainty on the number counts of &#8764;0.1-0.2 dex-comparable to the Poissonian uncertainties-although this value rapidly increases for higher flux densities. Overall, we therefore conclude that our number counts are in good agreement with both observed and simulated counts from the literature at lower frequencies.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="6.">34 GHz Continuum Source Properties</head></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="6.1.">Radio AGNs</head><p>In the local universe, the existence of a linear correlation between the total FIR and radio emission of star-forming galaxies has been well established <ref type="bibr">(Yun et al. 2001;</ref><ref type="bibr">Bell 2003)</ref>. This FIRRC has been shown to hold over a wide range of luminosities, from dwarf galaxies to dust-obscured starbursts <ref type="bibr">(Bell 2003)</ref>. The correlation is commonly expressed via parameter q IR , first introduced by <ref type="bibr">Helou et al. (1985)</ref>, and defined as We show the FIRRC for the eighteen 34 GHz continuum detections as a function of redshift in Figure <ref type="figure">8</ref>. In total, half of the COLDz sample are identified as radio AGNs, comprising eight AGNs in GOODS-N, and one in COSMOS (Table <ref type="table">4</ref>). The fact that radio AGNs make up a large fraction of the bright radio population is evidenced by the relatively large contribution of AGNs in the wider but shallower GOODS-N observations. Among the galaxies classified as star-forming, Figure <ref type="figure">6</ref> suggests the existence of an extended tail at 34 GHz in COLDz-cont-COS-2, at face value indicative of an AGN jet. However, such extended features are absent in the deeper 3 and 10 GHz observations of this source, and at these frequencies, the galaxy is consistent with being a point source at 2&#8243; resolution. As the spectrum of jetted AGNs tends to steepen toward higher frequencies (e.g., <ref type="bibr">Mahatma et al. 2018)</ref>, and should therefore be detectable in the lower-frequency ancillary radio data, we interpret the extended tail at 34 GHz as simply being due to noise. As we adopt the peak brightness for COS-2, its flux density measurement is unlikely to be substantially affected by this noisy region in the radio map.</p><p>In total, six out of eight radio AGNs in the GOODS-N field are additionally classified as AGNs through their strong X-ray emission. While constituting only a small number of sources, this is a relatively large fraction, as the overlap between radio AGNs selected at lower frequencies and X-ray AGNs is typically found to be small (&#61576;30%; e.g., <ref type="bibr">Delvecchio et al. 2017;</ref><ref type="bibr">Smol&#269;i&#263; et al. 2017;</ref><ref type="bibr">Algera et al. 2020b</ref>). This difference, however, may in part be due to the deeper X-ray data available in GOODS-N, compared to the COSMOS field where these studies were undertaken.</p><p>We further show the long-wavelength spectra of the nine radio AGNs in Figure <ref type="figure">9</ref>. The median 1.4-34 GHz spectral index of the sample equals &#945; = -0.77 &#177;0.13, which is consistent with the commonly assumed synchrotron slope of &#945; = -0.70. However, the modest sample still spans a wide range of spectral slopes between 1.4 and 34 GHz, ranging from nearly flat (&#945; = -0.04 &#177;0.02) to steep (&#945; = -1.47 &#177;0.06). In addition, the AGNs exhibit relatively smooth spectra, with the median spectral index between 1.4-3 (or 5) GHz of a = -&#61617; 0.70 0.12 3 5 1.4</p><p>being consistent with the typical high- . Only two AGNs exhibit strong evidence for steepening of their radio spectra toward higher frequencies, although sources with strongly steepening spectra are more likely to be missed in a selection at high radio frequencies. Such spectral steepening is expected to occur due to synchrotron aging losses increasing toward higher frequencies, and in turn relates to the age of the AGNs (e.g., <ref type="bibr">Carilli et al. 1991</ref>).</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="6.2.">Radio Spectral Decomposition for Star-forming Galaxies</head><p>Detecting radio FFE in high-redshift star-forming galaxies is challenging due to its expected faintness, and the presence of a radio AGNs only further hinders the detection of this already elusive component in the radio spectrum. As such, we now turn our attention to the star-formation-powered sources detected in the COLDz survey. The radio spectrum of star-forming galaxies is frequently assumed to be the superposition of two power laws arising from nonthermal synchrotron and thermal FFE <ref type="bibr">(Condon 1992;</ref><ref type="bibr">Murphy et al. 2017;</ref><ref type="bibr">Tabatabaei et al. 2017)</ref>. Denoting their spectral indices as &#945; NT and &#945; FF , respectively, the radio flux density at a given frequency &#957; may be written as</p><p>given a reference frequency &#957; 0 , as well as the thermal and nonthermal flux densities n S FF 0 and n S NT 0 , respectively, evaluated at this frequency. We rewrite this equation by introducing the thermal fraction, defined through</p><p>, as is common in the literature (e.g., <ref type="bibr">Condon 1992;</ref><ref type="bibr">Tabatabaei et al. 2017</ref>). As such, we rewrite the radio spectrum as</p><p>This further assumes that the spectral index for thermal FFE is fixed at &#945; FF = -0.10 (Condon 1992; <ref type="bibr">Murphy et al. 2011)</ref>. We then determine the remaining free parameters, a n f , th NT 0 and n S 0 , using a Monte Carlo Markov Chain (MCMC) based fitting routine. In addition, we adopt an observer-frame frequency of 1.4 GHz as the reference frequency &#957; 0 , which defines the frequency where the thermal fraction is normalized. Where necessary, we convert the thermal fraction from observedframe frequency &#957; to a rest-frame frequency n&#162; via We adopt flat priors in our fitting routine for the normalization n S 0 and thermal fraction f th , while we adopt Gaussian priors for the nonthermal spectral index (following, e.g., <ref type="bibr">Linden et al. 2020)</ref>. For the former two parameters, we require that -&lt; &lt; n S 1 Jy 1 Jy 0 and -0.5 &lt; f th &lt; 1.5, that is, we allow both the normalization and the thermal fraction to take on unphysical, negative values, as artificially bounding both to be greater than zero will not fully capture the uncertainties in the MCMC sampling. In addition, allowing negative values in the thermal fraction will demonstrate the necessity for the thermal component, as well as limitations inherent to the simple model we adopt for the radio spectrum (see also <ref type="bibr">Tabatabaei et 2017)</ref>.</p><p>The Gaussian priors adopted on the nonthermal spectral index are motivated by a degeneracy that manifests between the synchrotron slope and thermal fraction at low signal-to-noise. This degeneracy occurs because it is impossible to accurately distinguish between an overall flat spectrum as being due to dominant FFE, or an intrinsically shallow synchrotron slope. To partially alleviate this degeneracy, we adopt prior knowledge from the local universe that the synchrotron spectral indices are typically distributed around &#945; NT &#8776; -0.85 <ref type="bibr">(Niklas et al. 1997;</ref><ref type="bibr">Murphy et al. 2011)</ref>. In large radio-selected samples, the scatter around the typical low-frequency spectral index equals 0.3-0.5 dex <ref type="bibr">(Calistro Rivera et al. 2017;</ref><ref type="bibr">Smol&#269;i&#263; et al. 2017;</ref><ref type="bibr">Gim et al. 2019)</ref>, and is approximately Gaussian. While these spectral index measurements include both the synchrotron and free-free components, the low-frequency nature of these data ensure the radio fluxes are likely dominated by nonthermal emission, and as such, the overall variation in the spectral index constitutes a proxy for the scatter in the typical synchrotron spectral index. To encompass the full observed scatter in the synchrotron slopes, we therefore adopt a Gaussian prior on &#945; NT centered on a mean value of -0.85, with a scatter of 0.50.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="6.3.">Line and Dust Continuum Subtraction</head><p>Prior to fitting the radio spectra of the star-forming COLDz continuum detections, we need to ensure the 34 GHz flux density is not contaminated by thermal emission from dust, or by strong line emission. In particular, the 34 GHz flux density of the two COLDz COSMOS sources detected in <ref type="bibr">a.k.a. AzTEC.3,</ref>, may be boosted by their respective CO-lines. We re-extract the peak brightness of these two sources after removing the channels contaminated by the CO-emission, and, as a sanity check, repeat this for the three COSMOS sources that do not show any evidence for strong line emission. While for the latter, the flux densities are . However, the AGNs exhibit relatively smooth radio spectra, with only two sources showing strong evidence for spectral curvature.</p><p>unaffected by this procedure, we find line-uncontaminated flux densities for COS-2 and COS-4 of S 34 = 5.2 &#177;1.3 &#956;Jy and S 34 = 3.0 &#177;1.4 &#956;Jy, respectively, which are lower than the original cataloged flux densities by &#8764;25% and &#8764;40% (Table <ref type="table">1</ref>). In turn, this correction brings the 34 GHz flux density from COS-4 below the formal detection limit (S/N&#8776;2).</p><p>The 34 GHz continuum flux densities may further contain a contribution from thermal emission by dust, which, at least for local normal star-forming galaxies, is thought to dominate the radio spectrum beyond rest-frame &#957; &#61577;100-200 GHz <ref type="bibr">(Condon 1992)</ref>. To determine the extent of this contribution, we fit-where available-the FIR observations of our galaxy sample with both optically thin and thick modified blackbody spectra, and extrapolate the resulting dust SED to observed-frame 34 GHz. We note that this methodology is rather sensitive to how well the global dust properties (e.g., temperature and emissivity) can be constrained, and as a result the predicted flux densities are quite uncertain. Nevertheless, we find that the 34 GHz emission of COS-2 and COS-4 is likely to be dominated by emission from dust, with predicted dust contributions of m -+ 4.9 Jy , respectively. As a result, the full 34 GHz flux densities of these two sources are consistent with being powered by the combination of CO-emission and dust. As we probe restframe frequencies of n&#162; &#187; 210 GHz and n&#162; &#187; 120 GHz for COS-2 and COS-4, respectively, this finding is consistent with the typical model for the long-wavelength SED of star-forming galaxies <ref type="bibr">(Condon 1992)</ref>.</p><p>In what follows, we will discard COS-2 and COS-4 from our sample, as any remaining contribution from free-free or synchrotron emission to the measured 34 GHz flux density is not statistically significant, such that no robust spectral decomposition for these two sources can be performed. One source in GOODS-N, GN-7 at z = 2.95 (n&#162; &#187; 130 GHz), may have &#8764;25% of its continuum flux density contaminated by dust emission. However, due to the aforementioned uncertainties in the fitting of the dust SED, we do not correct for this potential contribution. This analysis indicates that, even at low flux densities (S 34 &#61576;25 &#956;Jy), a 34 GHz-selected sample does not automatically yield a free-free-dominated population.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="6.4.">The Radio Spectra of High-redshift Star-forming Galaxies</head><p>We show the radio spectra of the seven remaining star-forming sources across the COSMOS and GOODS-N fields in Figure <ref type="figure">10</ref>. All sources can be well described by the combination of a synchrotron and free-free component, although for three sources (COS-5, GN-7, and GN-11)-while some contribution from FFE is preferred-the fitted thermal fractions are consistent with zero within 1&#963;. In turn, for these sources, a single power law representing synchrotron emission is sufficient to match the observed flux densities. We additionally emphasize that there is considerable covariance between the thermal fraction and synchrotron slope, and as such, any quoted 1D uncertainties are not fully representative of the multidimensional posterior distribution (Figure <ref type="figure">11</ref>). We therefore utilize these full posterior distributions in order to propagate the uncertainties into physical quantities such as free-free SFRs (Section 6.5).</p><p>The fitted spectral parameters are presented in Figure <ref type="figure">11</ref> and Table <ref type="table">3</ref>. Our bootstrapped median thermal fraction, scaled to 1.4 GHz rest frame as is common in the literature, equals f th = 0.06 &#177; 0.03, with a standard deviation of &#963; = 0.05. None of the seven star-forming galaxies exhibit thermal fractions of f th &#61577;0.20, indicating a fairly narrow distribution of f th at restframe 1.4 GHz, even among a sample showing substantial variation in SFRs. Our average thermal fraction is slightly lower than the average value observed by <ref type="bibr">Tabatabaei et al. (2017)</ref> of f th = 0.10 for star-forming galaxies in the local universe, though they report a large scatter of &#963; = 0.09. Additionally, the typical thermal fraction is similar to what was observed by <ref type="bibr">Niklas et al. (1997)</ref>, who determined f th = 0.08 &#177; 0.01 at 1 GHz, with a scatter of &#963; = 0.04 across 74 local galaxies.</p><p>We further determine an average thermal fraction at observed-frame 34 GHz (i.e., probing 34&#215; (1 + z) GHz rest frame) of f th = 0.78 &#177;0.07, with a range of f th = 0.45 -0.95, and a scatter of &#963; = 0.20 (Figure <ref type="figure">12</ref>). As such, we find that, even at rest-frame frequencies &#957; &#61577;60 GHz, the radio spectrum is not fully dominated by thermal FFE, though we caution that the uncertainties on the individual thermal fractions are large. We first compare these results with two local studies, both of which map FFE on subkiloparsec scales. At a typical resolution of &#8776;0.9 kpc, <ref type="bibr">Murphy et al. (2012)</ref> find an average thermal fraction across 103 star-forming regions of f th = 0.76 at restframe 33 GHz, with a scatter of &#963; = 0.24. Extrapolating this value via the simple model from <ref type="bibr">Condon (1992)</ref>, the typical thermal fraction at &#8764;60 GHz is expected to be &#8764;0.85-0.90, which is slightly higher than the thermal fraction we find for high-redshift star-forming galaxies. In addition, <ref type="bibr">Linden et al. (2020)</ref> recently measured a typical thermal fraction at 33 GHz of f th = 93 &#177;0.8% across 118 star-forming complexes in local galaxies, at a resolution of &#8776;0.2 kpc. Their typical thermal fraction is both higher than that determined by <ref type="bibr">Murphy et al. (2012)</ref> and the values measured in this work. This, however, is not surprising, given that the thermal fractions presented in this work are integrated over the entire galaxy. As FFE is predominantly produced in star-forming regions, spatial variations in the thermal fraction across a galaxy are naturally expected, with the thermal fraction peaking in star-forming complexes.</p><p>At high redshift, no previous studies have directly targeted the free-free-dominated regime (&#957; &#61577;30 GHz) in blindly selected galaxy samples, at a depth where star-forming galaxies are expected to dominate the radio population. However, particularly in bright dusty star-forming galaxies, some works have serendipitously detected high-frequency radio continuum emission, typically as a byproduct when targeting the CO(1-0) line. <ref type="bibr">Thomson et al. (2012)</ref> detect FFE in two z &#8764;2.9 lensed submillimeter galaxies, and determine thermal fractions of f th &#8764;0.3-0.4 at 34 GHz. Other studies of highly star-forming galaxies <ref type="bibr">(Aravena et al. 2013;</ref><ref type="bibr">Huynh et al. 2017</ref>) have additionally detected radio continuum emission at observedframe &#8764;30-35 GHz, but had to assume fixed synchrotron spectral indices due to a lack of ancillary data. Nevertheless, they estimate thermal fractions between f th &#8764;40%-70%. Overall, these studies find thermal fractions that are broadly consistent with, albeit typically slightly lower than, what we determine for the star-forming sample detected in our nontargeted 34 GHz observations. Finally, we compare our results with the 10 GHz pointing from <ref type="bibr">Murphy et al. (2017)</ref> in GOODS-N. They determine thermal fractions from 1.4-10 GHz spectral indices for &#8764;25 galaxies, under the assumption of a fixed synchrotron slope of &#945; NT = -0.85. At a typical rest-frame frequency of &#957; &#8764;20 GHz, they find a median thermal fraction of f th &#8776;50%. Extrapolating this value to rest-frame 60 GHz, this would imply a thermal fraction of f th &#8764;0.7, similar to what we observe among the COLDz star-forming sample.</p><p>We further determine a median nonthermal spectral index for the star-forming COLDz sample of a = - -+ 0.99 NT 0.37 0.19 , with a standard deviation of &#963; = 0.25. This typical value is consistent with the value observed by <ref type="bibr">Tabatabaei et al. (2017)</ref> for local starforming galaxies of &#945; NT = -0.97 &#177;0.16, but is slightly steeper than that of individual star-forming regions in NGC 6946, where <ref type="bibr">Murphy et al. (2011)</ref> find a typical value of &#945; NT = -0.81 &#177;0.02. This is not surprising, as these observations directly target the acceleration sites of cosmic rays, where the spectrum should be flatter. However, our median synchrotron slope is additionally slightly steeper than the value obtained by <ref type="bibr">Niklas et al. (1997)</ref>, who determine an average &#945; NT = -0.83 &#177;0.02 (&#963; = 0.13) across 74 local galaxies. The slightly steeper nonthermal spectral index we find for the COLDz sample may be the result of the higher rest-frame frequencies probed in this work, compared to the aforementioned local studies. Synchrotron cooling losses increase toward high rest-frame frequencies, and result in the steepening of the nonthermal spectral index (e.g., <ref type="bibr">Thomson et al. 2019)</ref>. While our data do not have the constraining power to determine whether Figure <ref type="figure">10</ref>. The radio spectra of the seven star-forming galaxies in the COSMOS (first two panels) and GOODS-North (last five) fields. We show the decomposition of the spectra into their synchrotron and free-free components, with the shaded regions indicating the 1&#963; confidence region on the fits. We find the radio emission for four out of seven galaxies (COS-1, GN-4, GN-9, and GN-12) to be dominated by free-free emission at observed-frame 34 GHz, whereas for the remainder, only a relatively minor thermal contribution is preferred. such spectral aging occurs, as this requires additional sampling of the radio spectrum, such increased high-frequency losses would be fitted by a relatively steep synchrotron slope in our twocomponent model, and may plausibly contribute to our moderately lower value for &#945; NT compared to local studies.</p><p>Based on the spectral parameters we determine for the COLDz sample, we calculate the rest-frame frequency n &#162; 50 where the thermal fraction reaches 50% via</p><p>While the uncertainties on this value are substantial for the individual star-forming galaxies due to its dependence on both the thermal fraction and synchrotron slope, we determine a median value of n &#162; = -+ 11.5 50 7.4 20.0 GHz, with a standard deviation of &#963; = 5.9 GHz. This is slightly lower than, albeit still consistent with, the canonically assumed value of n &#162; &#187; -25 30 50 GHz (e.g., <ref type="bibr">Condon 1992)</ref>, which is likely due to the relatively steep synchrotron slopes we are finding for the COLDz star-forming galaxies.</p><p>While the differences between the spectral parameters determined for COLDz and those observed in local starforming galaxies are of minor statistical significance, a highfrequency selection of star-forming galaxies is likely to bias the sample toward having overall shallow radio spectra. Combining the 10-34 GHz spectral indices for the star-forming galaxies across COSMOS and GOODS-N, we find an average slope of a = - &#61617; 0.47 0.10 34 10 (&#963; = 0.28), which is substantially shallower than the canonical &#945; = -0.70 assumed at lower frequencies. Shallow spectra are naturally expected in the freefree dominated regime, in particular for young starburst galaxies. Such sources should exhibit large thermal fractions, as synchrotron emission lags the onset of a starburst by &#61577;30 Myr (e.g., <ref type="bibr">Bressan et al. 2002)</ref>, and the galaxies should hence be dominated by FFE across the entire radio spectrum. However, our sample is fully comprised of sources with typical or low thermal fractions ( f th &#61576;0.20 at 1.4 GHz), and rather steep  </p><p>COLDz-cont-COS-1 64 Notes.</p><p>a The rest-frame frequency probed by the COLDz 34 GHz observations. b The thermal fraction at rest-frame 1.4 GHz. c The thermal fraction at rest-frame frequency n&#162;.</p><p>d The uncertainties on the free-free SFRs are propagated from those on the thermal fraction and radio luminosity. The sources where these SFRs are consistent with zero within 1&#963; are plotted as upper limits in Figure <ref type="figure">13</ref>.</p><p>Figure <ref type="figure">12</ref>. Violin diagrams of the thermal fraction at 34 GHz (observed-frame, hence probing rest-frame 34&#215; (1 + z) GHz) as a function of redshift. The width of the shaded regions represents the probability distribution of the thermal fraction, while the data points indicate the median and 16-84th percentiles. The probability distributions of the thermal fractions are characterized by a typically high probability of having a large thermal fraction f th &#61577;0.7 -0.8, with a long tail extending toward lower values. The model from <ref type="bibr">Condon (1992)</ref>, adopting &#945; NT = -0.85 and f th = 0.10 at 1.4 GHz, is shown as the solid black line, and predicts similar high-frequency thermal fractions as observed across the starforming galaxies.</p><p>synchrotron spectra (&#945; NT &#61576;-0.85). This, in turn, indicates we are likely detecting relatively mature starbursts, as opposed to young star-forming galaxies. In particular, galaxies with a declining star formation history are expected to exhibit only modest levels of FFE at low radio frequencies, and may show steepening of their synchrotron spectra toward higher frequencies. This scenario is qualitatively consistent with the spectral parameters we are finding for the COLDz star-forming galaxies, and may therefore be typical for a high-frequencyselected sample.</p><p>Alternatively, the relatively low thermal fractions observed for our star-forming sample at low frequencies may be the result of residual AGN contamination. An AGN will contribute additional synchrotron emission in excess of that arising from star formation, while the overall contribution from FFE is mostly unaffected. In this case, one may expect typical synchrotron spectra, in combination with low thermal fractions. However, we find no strong evidence for systematic residual radio AGN activity based on the values of q IR of our starforming galaxies, as we find that half of our sample falls onto or above the FIRRC for star-forming galaxies <ref type="bibr">(Delhaize et al. 2017;</ref><ref type="bibr">Algera et al. 2020a</ref>).</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="6.5.">Free-Free Star Formation Rates</head><p>Free-free emission is one of the most robust tracers of star formation, as it constitutes a direct probe of the ionizing photons emitted by recently formed (&#61576;10 Myr) massive stars. As such, it does not rely on reprocessed starlight, or emission produced by stellar remnants, such as, respectively, FIR and radio synchrotron emission. However, some caveats still apply, as carefully summarized in <ref type="bibr">Querejeta et al. (2019)</ref>. In particular, any ionizing photons absorbed by dust within the H II region will not contribute to the ionization of hydrogen atoms, and as such will reduce the free-free luminosity at a fixed SFR (e.g., <ref type="bibr">Inoue et al. 2001;</ref><ref type="bibr">Dopita et al. 2003)</ref>. Alternatively, if there is substantial leakage of ionizing photons, FFE will similarly be suppressed. Finally, since FFE is only sensitive to the massive end of the initial mass function, any variations in the IMF may substantially affect the calculated SFRs. We note that these caveats also apply to SFRs estimated using the Balmer lines (e.g., H&#945;, H&#946;), with the clear advantage of FFE being that it is fully dust-insensitive on scales beyond the H II region wherein the star formation occurs.</p><p>With these caveats in mind, we now set out to calculate freefree SFRs for the seven star-forming galaxies detected at 34 GHz. The calibration from <ref type="bibr">Murphy et al. (2012)</ref>, adapted to a Chabrier IMF, is given by Here T e is the electron temperature, which we assume to equal T e = 10 4 K. However, we note that our results are somewhat insensitive to the precise value adopted, given the modest exponent of - T e 0.45 . The SFR is further directly proportional to the product of the thermal fraction and the radio luminosity, which we evaluate at a rest-frame frequency of &#957; = 1.4 GHz.</p><p>We derive typical free-free SFRs between SFR&#8776;30-350 M e yr -1 for the four out of seven sources for which we have robustly constrained thermal fractions. The remaining sources instead have a thermal fraction that is consistent with zero within 1&#963;, such that we can only provide upper limits on their free-free SFRs. We compare the free-free SFRs with those derived from MAGPHYS, as well as from low-frequency radio synchrotron emission in Figure <ref type="figure">13</ref>. For the latter, we adopt the FIRRC from <ref type="bibr">Delhaize et al. (2017)</ref>, which is suitable for the radio-detected star-forming population. While we have derived individual q IR -values for the star-forming COLDz sample in Section 4, we adopt a fixed FIRRC from the literature to ensure that the synchrotron-derived SFRs are independent of those from FIR emission.</p><p>Overall, we observe a reasonable agreement between the SFRs from FFE, and those from the more commonly adopted tracers we use for comparison. This likely implies that the various aforementioned caveats do not greatly affect our calculated SFRs. However, the correlation between the two radio-based tracers appears tighter than the one between the SFRs from FFE and SED-fitting, which is surprising, as the timescale for FIR emission should more closely resemble that of FFE than synchrotron emission. Two sources in particular, Figure <ref type="figure">13</ref>. Left: a comparison of the star formation rates obtained from free-free emission, vs. those from MAGPHYS. The one-to-one relation is shown through the solid black line, and the vertical error bars represent the propagated uncertainties on the radio luminosity and thermal fraction. For three sources, we can only place upper limits on their free-free SFRs. The diamonds show the SFRs that we infer from SED-fitting when we instead convert the FIR luminosity into a star formation rate, for the two sources where these values are discrepant (see the text for details). Right: free-free star formation rates vs. those obtained from radio synchrotron emission, adopting the FIRRC from <ref type="bibr">Delhaize et al. (2017)</ref>. Despite the uncertainties on SFR FF being large as a result of the low signal-to-noise at 34 GHz, and the limited available sample size, we find that the free-free star formation rates are in reasonable agreement with the values derived from SED-fitting and low-frequency radio emission. ) and GN-9 (z = 0.87) appear to have wellconstrained free-free SFRs that exceed the ones from SEDfitting by a factor of 5.6 &#177;0.1 and 7.0 &#177;0.3, respectively. If dust extinction within the H II region, or leakage of ionizing photons were a concern, the free-free SFRs should be suppressed, in contrast to what we observe for these two sources. Instead, for GN-4, this offset is likely related to the SFR determined from SED-fitting. While the infrared SED of this galaxy is well constrained (Figure <ref type="figure">17</ref>), MAGPHYS predicts that a substantial fraction of the infrared emission from GN-4 originates from an older stellar population, with the ratio between its infrared-based SFR-assuming the conversion from <ref type="bibr">Kennicutt (1998)</ref> adjusted for a Chabrier IMF-and the fitted SFR equaling 3.2 &#177;0.1. Given that the SFRs from FFE and the 1.4 GHz luminosity are in good agreement for GN-4, it is likely that the contribution from old stars to the FIR luminosity is overestimated in the SED-fitting.</p><p>A similar discrepancy can be seen between the free-free and infrared SFRs derived for GN-9. While this source too has a modest contribution to its dust luminosity from older stars, upon accounting for this, its ratio between the free-free and SED-fitted SFRs remains a factor of - + 3.8 0.3 0.1 . Instead, the radio flux densities of this source are likely to be boosted by the emission from an AGN. We find an FIRRC parameter of q IR = 1.99 &#177;0.04 for GN-9 (Section 6.1), which implies that it lies &#8764;0.5 dex below the median value for radio-selected starforming galaxies at its redshift of z = 0.87 <ref type="bibr">(Delhaize et al. 2017)</ref>, though still above the threshold we adopt for identifying radio AGNs. Nevertheless, GN-9 is &#8764;3.5&#215; radio-bright with respect to the median FIRRC, which fully accounts for the difference between its radio and SED-fitted SFRs.</p><p>Interestingly, we can only place upper limits on the free-free SFRs for two of the brightest star-forming galaxies in GOODS-N (SFR SED &#8764;400 M e yr -1 ), GN-7 and GN-11. While naively these galaxies might be expected to exhibit strong FFE, recent studies of starburst galaxies have found that their radio spectra might steepen toward higher frequencies <ref type="bibr">(Thomson et al. 2019;</ref><ref type="bibr">Tisani&#263; et al. 2019</ref>). This may be related either to a deficit in FFE, or to the steepening of their synchrotron spectra. However, these studies were limited to rest-frame frequencies &#957;&#61576;20 GHz, and hence did not probe the regime where FFE is expected to dominate. While our study does directly target this high-frequency regime, with the current sampling of the radio spectrum, we are unable to distinguish between a deficit in FFE, or more complex behavior in the synchrotron emission. Nevertheless, our lack of a robust detection of FFE in these strongly star-forming sources provides support for the existence of more complex radio spectra, which may be constrained by probing their radio emission at intermediate frequencies, using, for example, the VLA K-band (22 GHz).</p><p>Taking the aforementioned caveats into account, our derived free-free SFRs are in good agreement with those from SEDfitting and the low-frequency FIRRC. While the uncertainties on the free-free SFRs remain large as a result of the typical faintness of star-forming galaxies at high rest-frame frequencies, our analysis indicates that a deep, 34 GHz-selected sample, supplemented by deep ancillary radio observations, can be used to accurately constrain star formation at high redshift. This, in turn, will be possible for significantly larger galaxy samples in the future, with the increased sensitivity of nextgeneration radio facilities.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="7.">Free-Free Emission with the SKA and ngVLA</head><p>The next large radio telescope to come online is the Square-Kilometer Array Phase 1 (SKA1), with SKA1-Mid set to cover a frequency range of 0.35-15 GHz. As such, observations with the highest-frequency band of the SKA1-Mid (at a central frequency of &#957; c = 12.5 GHz) will start probing the regime where FFE dominates in star-forming galaxies at z &#61577;1. A galaxy with SFR = 10 M e yr -1 at z = 1 will have a flux density of approximately S 12.5 &#8776;1.5 &#956;Jy in this band, assuming the FIRRC from <ref type="bibr">Delhaize et al. (2017)</ref>, and the calibration between star formation and FFE from <ref type="bibr">Murphy et al. (2012)</ref>. This, in turn, requires &#8764;15-20 hr of telescope time for a 5&#963; detection, based on the SKA1 sensitivity estimates from <ref type="bibr">Braun et al. (2019)</ref>. In order to robustly probe the FFE in such a modestly star-forming galaxy, additional sampling of its low-frequency radio spectrum is crucial. In particular, for a similar 5&#963; detection at 1.4 and 6.7 GHz, a further &#8764;10-15 hr of total telescope time is required. Given the &#8764;4&#215; larger field of view at 6.7 GHz compared to at 12.5 GHz, a possible observing strategy for the detection of FFE in faint star-forming galaxies is to combine two single SKA1-Mid pointings at 1.4 and 6.7 GHz with a five-pointing mosaic at 12.5 GHz, covering the entire 6.7 GHz field of view. With a total telescope time of &#8764;100 hr, this allows for the mapping of FFE in all z &#61577;1 star-forming galaxies at S 12.5 &#61577;1.5 &#956;Jy across an area of ~120 arcmin 2 . Adopting the <ref type="bibr">Bonaldi et al. (2019)</ref> simulations of the radio sky, developed specifically for the SKA1, typical &#8764;1100-1200 galaxies are expected at S 12.5 &#61577;1.5 &#956;Jy within this field of view. In particular, approximately &#8764;68 &#177;2% (&#8764;12 &#177;1%) of this sample is expected to lie at a redshift z 1 (z 3), allowing for the robust sampling of the free-free dominated regime. For comparison, the 200 hr VLA COSMOS-XS survey (van der Vlugt et al. 2021), reaches a similar depth to these template SKA1-Mid observations at 3 and 10 GHz, but covers a smaller area of ~30 arcmin 2 . The increased survey speed of SKA1-Mid, therefore, allows for a &#61577;8&#215; quicker mapping of FFE up to &#8764;15 GHz, compared to the VLA. However, as the frequency coverage of SKA1-Mid is not fully optimized to directly probe the high-frequency radio emission in star-forming galaxies, significant synergy with the VLA remains, as it allows for the extension of the spectral coverage from the SKA1 to higher frequencies.</p><p>The ngVLA <ref type="bibr">(Murphy et al. 2018;</ref><ref type="bibr">McKinnon et al. 2019</ref>), however, is set to truly transform our understanding of the high-frequency radio spectrum in distant galaxies, and will allow for the usage of FFE as a high-redshift SFR-tracer on an unprecedented scale. The current sensitivity estimates <ref type="bibr">(Butler et al. 2019)</ref> indicate that the ngVLA will attain a typical rms of &#963;&#8776;0.3 &#956;Jy beam -1 in one hour of band four (&#957; = 20.5-34.0 GHz) observations. This, in turn, translates to a 5&#963; detection of a galaxy forming stars at 100 M e yr -1 at z = 2.5, probing restframe 80 GHz, similar to the frequency range probed in this work for z&#8776;2 star-forming galaxies. However, with the expected improvement in sensitivity the ngVLA provides over the current VLA, high-redshift sources may more easily be targeted at relatively low observing frequencies, enabling wider surveys while still allowing for the free-free dominated regime to be probed. For example, <ref type="bibr">Barger et al. (2018)</ref> propose a survey at 8 GHz (ngVLA band two) to 0.2 &#956;Jy beam -1 across a large area of 1 deg 2 , which requires just an hour per pointing. By adopting a typical wedding cake strategy, deeper observations across a smaller area can further be used to target fainter star-forming galaxies. As an example, a star-forming galaxy of SFR = 25 M e yr -1 at z = 3 (z = 5) can in principle be detected at 8 GHz in only &#8764;3 hr ( &#8764;15 hr), modulo, of course, the large uncertainties on the typical thermal fraction in faint, star-forming sources, and the nature of the FIRRC in this population.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head n="8.">Conclusions</head><p>We have presented a deep continuum survey of the highfrequency radio sky with the VLA, which probes the microjansky galaxy population at 34 GHz. This regime has historically remained largely unexplored due to the relatively low survey speed of radio telescopes at these frequencies, as well as the expected faintness of sources. However, the high radio frequencies hold one of the most reliable tracers of star formation, radio FFE, and as such are set to become a key area of study with next-generation radio facilities.</p><p>We employ deep observations at 34 GHz from the COLDz project <ref type="bibr">(Pavesi et al. 2018;</ref><ref type="bibr">Riechers et al. 2019</ref><ref type="bibr">Riechers et al. , 2020))</ref>, which cover the well-studied COSMOS (10 arcmin 2 ) and GOODS-North (50 arcmin 2 ) fields to a typical depth of &#8764;1.5 &#956;Jy beam -1 and &#8764;5.3 &#956;Jy beam -1 , respectively. We perform source detection on the images down to a liberal 3&#963; detection threshold, aided by deep ancillary radio data across both fields, resulting in the detection of high-frequency continuum emission in 18 galaxies. We cross-match these detections with additional deep radio observations at 1.4, 3, and 10 GHz in the COSMOS field, as well as data at 1.4, 5, and 10 GHz in GOODS-N. In addition, we leverage the wealth of multiwavelength data across both fields to fully sample the SEDs of the galaxies from the X-ray to radio regime. The COLDz continuum sample spans a redshift range of z = 0.50-5.30, and lies at a median (mean) redshift of = -+ z 1.12 0.15 0.52 ( = -+ z 1.55 0.35 0.41</p><p>). The sample contains six sources at z 2, and includes the well-studied submillimeter galaxy AzTEC.3 at z = 5.3. Our main findings are as follows:</p><p>1. We present the first constraints on the radio number counts at 34 GHz in the regime where star-forming galaxies dominate the radio population (Figure <ref type="figure">7</ref>), and find that these are in good agreement with lowerfrequency number counts in the literature, both from observations (van der Vlugt et al. 2021) and simulations <ref type="bibr">(Bonaldi et al. 2019</ref>). 2. We use the FIRRC to divide the 34 GHz continuum sample into star-forming galaxies and AGNs (Figure <ref type="figure">8</ref>). In total, half of the sample (nine sources) shows AGN activity at radio wavelengths (Figure <ref type="figure">9</ref>), while the radio emission of the remainder is consistent with being powered predominantly through star formation. All but one of the faintest galaxies in our sample (S 34 &#61576;20 &#956;Jy) show radio emission of a star-forming origin, which is qualitatively consistent with the small fraction of radio AGNs found in deep observations at lower frequencies (e.g., <ref type="bibr">Algera et al. 2020b)</ref>. Two sources, including AzTEC.3 at z = 5.3, likely have their continuum emission at 34 GHz dominated by thermal emission from dust, leaving 7/18 sources ( &#8764;40%) of the sample with highfrequency radio continuum emission dominated by the combination of synchrotron and FFE. 3. We use the wealth of ancillary radio data across the COSMOS and GOODS-N fields to construct radio spectra of the star-forming galaxies, covering four frequencies in the range 1.4-34 GHz (Figure <ref type="figure">10</ref>). We fit the radio spectra with a combination of free-free and synchrotron emission, and determine thermal fractions and nonthermal spectral indices for our sample <ref type="bibr">(Figures 11 and 12)</ref>, which are consistent with the values observed in local galaxies. We further determine freefree SFRs for seven star-forming galaxies, and find good agreement with those obtained from SED-fitting and the FIRRC (Figure <ref type="figure">13</ref>).</p><p>With the 34 GHz continuum data from the COLDz survey, we have directly targeted FFE in faint star-forming sources at high redshift. While currently limited to a modest sample, nextgeneration radio facilities are set to significantly increase the number of galaxies for which the full radio spectrum is constrained, and will transform our understanding of highfrequency radio continuum emission in star-forming galaxies. However, combined with the wealth of ancillary data in the COSMOS and GOODS-N fields, the COLDz observations already allow for a census of FFE in the typical star-forming population via a multifrequency radio stacking analysis, which will be presented in a forthcoming publication (H. S. B. <ref type="bibr">Algera et al. 2021, in preparation)</ref>.</p><p>density S &#957; as A(r)&#215; S &#957; , where A(r) represents the primary beam sensitivity at position r within the mosaic. We randomize source positions within both mosaics, above A(r) 0.20, and draw flux densities from a power-law distribution to ensure low flux densities-where incompleteness will be the largest-are amply sampled. For COSMOS, we insert 50 mock sources in each run, for a total of 200 runs. As the GOODS-N mosaic is substantially larger, we instead insert 100 mock sources per run, for 100 runs total. In both cases, mock sources are required to be 2.5 beam sizes away from both real sources and other mock sources. We then repeat the source detection procedure described in Section 3, and cross-match the recovered sources to the inserted ones, using a matching radius of 0 7. We record their inserted flux density, as well as their recovered peak and integrated flux densities. We define the completeness in a given flux density bin i as C i = N rec,i /N ins,i , that is, as the ratio of the number of inserted and recovered mock sources with a flux density that falls within the ith bin. We determine the corresponding uncertainty via a bootstrap analysis, whereby we resample from the inserted flux densities, with replacement, and determine for each flux density bin the fraction of this sample that was recovered in our source detection procedure. The uncertainty then represents the 16th-84th percentile of the bootstrapped completeness analyses. We show the completeness in both the COSMOS and GOODS-N mosaics in Figure <ref type="figure">14</ref>.</p><p>In the COSMOS field, we reach 50% and 80% completeness at flux densities of S 34 = 7.3 &#956;Jy beam -1 (&#8776;6&#963;, where &#963; represents the typical rms in the map) and S 34 =13.3 &#956;Jy beam -1 (&#8776;10&#963;), respectively. In GOODS-North, we reach completeness fractions of 50% and 80% at S 34 = 19.6 &#956;Jy beam -1 (&#8776;4&#963;) and S 34 = 29.3 &#956;Jy beam -1 (&#8776;6&#963;), respectively. These differences may be explained by the nonuniform exposure map of the GOODS-N mosaic, allowing for faint sources to still be detected in a small portion of the mosaic.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>A.2. Peak versus Integrated Fluxes</head><p>In order to assign a flux density to the detected radio sources, we need to establish if they are resolved. Our observations have a typical beam size of &#8764;2 5, and as such, we expect most sources to be unresolved, based on the typical (sub-)arcsecond radio sizes of star-forming galaxies at low frequencies <ref type="bibr">(Cotton et al. 2018;</ref><ref type="bibr">Jim&#233;nez-Andrade et al. 2019)</ref>, and the finding that these sources are more compact at higher frequencies <ref type="bibr">(Murphy et al. 2017;</ref><ref type="bibr">Thomson et al. 2019</ref>). To verify this, we use our runs of inserted mock sources, which, by construction, are unresolved, and compare their peak and integrated flux densities as a function of signal-to-noise. We show the results in Figure <ref type="figure">15</ref>, for both the COSMOS and GOODS-N fields. We divide the results into logarithmically spaced bins in S/N, and determine the integrated/ peak ratio encompassing 95% of sources per bin (following, e.g.,  Both panels show a power-law fit to the upper 95th percentile of each bin in S/N (red points) via the red, dashed line. Sources below this line are taken to be unresolved, and for these, the peak brightness is adopted. The robust sources detected in both fields are shown in blue. A single source in GOODS-N is consistent with being resolved, whereas the remaining continuum detections are unresolved. Since the mock sources are inserted into the mosaic uncorrected for the primary beam, the rms is highly uniform, and hence the S/N can be mapped into a peak brightness (upper horizontal axis). <ref type="bibr">Bondi et al. 2008;</ref><ref type="bibr">van der Vlugt et al. 2021)</ref>. These percentiles are fitted with a power law, with the region below the best fit defining the limiting integrated/peak flux density where sources are taken to be unresolved. It is clear that at a modest signal-to-noise, S/ N&#61576;5, sources can have S int /S peak &#61577;2, despite being unresolved. This is a result of nearby noise peaks elongating the source, or throwing off the fitting, which mostly affects the integrated flux density. We find that all sources, with the exception of one bright (S/N &#8764;25) detection in GOODS-N, show an integrated/peak ratio that is consistent with the source being unresolved. As such, we adopt the integrated flux density for this single source, and the peak brightness for the rest.</p></div>
<div xmlns="http://www.tei-c.org/ns/1.0"><head>A.3. Flux Boosting</head><p>At low signal-to-noise, the peak brightness may be "boosted" as a result of noise properties in the image. To establish whether this is affecting the flux densities of our 34 GHz detections, we compare the recovered and inserted flux densities of our mock source analysis. The results are shown in Figure <ref type="figure">16</ref>, for both COSMOS and GOODS-N. At S/N &#61577;5, the median ratio of recovered-to-inserted flux density is consistent with unity, with a spread of less than&#61576;20%. At S/N&#61576;5, which is the typical signal-to-noise at which the faintest 34 GHz sources are detected, the level of flux boosting steadily increases, with a typical correction of &#8764;20% at S/N&#8776;3 in either field. In this low-S/N regime, the typical spread on the ratio of recovered and inserted flux densities similarly increases strongly. Following, e.g., <ref type="bibr">Stach et al. (2019)</ref>, we correct the flux densities of the sources detected at 34 GHz by the median level of flux boosting at their observed S/N. The uncertainty on the corrected flux density includes the propagated bootstrapped error on the median. </p></div><note xmlns="http://www.tei-c.org/ns/1.0" place="foot" xml:id="foot_0"><p>The Astrophysical Journal, 912:73 (23pp), 2021 May 1 Algera et al.</p></note>
			<note xmlns="http://www.tei-c.org/ns/1.0" place="foot" n="12" xml:id="foot_1"><p>However, we note that we adopt a 3&#963; detection threshold for the COLDz data in Section 3, whereas a 5&#963; threshold was adopted by<ref type="bibr">Schinnerer et al. (2007)</ref> for the 1.4 GHz observations. As such, sources detected at 3&#963; in the COLDz survey require a spectral index of &#945;&#61576;-0.85 to additionally be detected at 1.4 GHz.</p></note>
			<note xmlns="http://www.tei-c.org/ns/1.0" place="foot" n="13" xml:id="foot_2"><p>We have verified that, after applying the primary beam correction, the flux densities are consistent with those obtained from running source detection on the primary-beam-corrected map, with a typical ratio of S uncorr /S corr = 1.05 &#177; 0.07.</p></note>
			<note xmlns="http://www.tei-c.org/ns/1.0" place="foot" n="14" xml:id="foot_3"><p>As such, the observed 34 GHz continuum emission may be in part due to these bright CO-lines. We investigate this in Section 6.3.</p></note>
		</body>
		</text>
</TEI>
