Confirmation of a Substantial Discrepancy between Radio and UV–IR Measures of the Star Formation Rate Density at 0.2 < z < 1.3

We present the initial sample of redshifts for 3839 galaxies in the MeerKAT DEEP2 field—the most sensitive ∼1.4 GHz radio field yet observed with σ n = 0.55 μ Jy beam−1, reaching the confusion limit. Using a spectrophotometric technique combining coarse optical spectra with broadband photometry, we obtain redshifts with σ z ≲ 0.01(1 + z), as determined from repeat observations. The resulting radio luminosity functions between 0.2 < z < 1.3 from our sample of 3839 individual galaxies are in remarkable agreement with those inferred from previous modeling of radio source counts, confirming a ≳50% excess in radio-based star formation rate density (SFRD) (z) measurements at 0.2 < z < 1.3 compared to those from the UV–IR. Several sources of systematic error are discussed—totalling ∼0.13 dex when added in quadrature. Even in the event that all systematic errors work to decrease the radio-based SFRD values, they are incapable of reconciling differences between the radio-based measurements with those from the UV–IR at 0.5 < z < 1.3. We conclude that significant work remains to have confidence in a full accounting of the star formation budget of the Universe.


INTRODUCTION
Since it was first discovered in the 1990s that the comoving star-formation rate density (SFRD, usually specified in M ⊙ yr −1 Mpc −3 ) around z ∼ 1 dwarfs that of today, immense progress has been made in understanding the general shape of SFRD evolution with cosmic time.Collections of SFRD measurements across wide redshift ranges place the peak of star formation activity near z ∼ 2 and show a decline by a factor of ∼10 to present day (see Hopkins & Beacom 2006;Madau & Dickinson 2014, for reviews of the topic).
Despite consistent agreement among the UV/optical and infrared (IR) communities on the SFRD through z ∼ 2, there has been established tension between the SFRD evolution derived at shorter wavelengths (e.g.UV and optical) with that derived at longer wavelengths (e.g.radio and sub-millimeter).Studies with ALMA have shown that UV measurements of the SFRD miss the obscured, dusty, and heavily star-forming galaxies more common in the distant universe (e.g.Casey et al. 2018;Bouwens et al. 2020).Even in the more "recent" past, Whitaker et al. (2017) found that > 80% of star formation is obscured at all redshifts z < 2.5 in galaxies having stellar masses log[M * /M ⊙ ] ≥ 10.FIR/sub-mm observations are essential to understand dusty galaxies, but single-dish telescopes are only sensitive to the massive, tip-of-the-iceberg sources.Sub-mm interferometers such as ALMA are limited to a very small field-of-view, making it expensive and impractical to amass a deep sample over a large enough area to minimize cosmic variance.
Radio observations are immune to dust obscuration, insensitive to contamination by older stellar populations, and available across wide sky areas.While active galactic nuclei (AGNs) are primarily responsible for powering strong 1.4 GHz radio sources, star-forming galaxies dominate the radio emission below 450 µJy (Algera et al. 2020;Matthews et al. 2021a).Star-forming galaxies produce radio emission through two processes: (1) thermal bremsstrahlung from H II regions ionized and heated by massive stars and (2) synchrotron radiation from cosmic-ray electrons accelerated in the shocks of supernova remnants (Condon 1992).
The tight relationship between the radio synchrotron and FIR luminosities of star-forming galaxies, described by the "q" parameter q ≡ log FIR/(3.75 × 10 12 Hz S 1.4 GHz (W m −2 Hz −1 ) , where 1.26 × 10 −14 [2.58 S 60 µm (Jy) + S 100 µm (Jy)] (3) (Helou et al. 1988) was first derived in the local universe (Helou et al. 1985).It has been shown to evolve to varying degrees (from not-at-all to moderately) through redshift z ∼ 4 (e.g.Ivison et al. 2010;Magnelli et al. 2015;Pannella et al. 2015;Delhaize et al. 2017).More recently, Delvecchio et al. (2021) found that while the FIR/radio correlation is largely invariant with redshift, there is a statistically significant dependence on stellar mass dq/d log(M * ) = −0.148± 0.013.Given the established relationship between stellar mass and star formation rate (e.g.Whitaker et al. 2012), this result is largely consistent with the nonlinear local FIR/radio correlation dq/d log(L ν ) = −0.147 of Matthews et al. (2021a), derived from the largest, clean (free from AGNdominated radio sources) sample of star-forming galaxies (∼ 4, 300) in the local universe.Further, previous studies have also found non-linearity in the FIR/radio correlation that may explain much of the perceived evolution in the q parameter as a function of redshift (e.g.Basu et al. 2015;Molnár et al. 2021).Unfortunately, galaxies whose radio emission is powered primarily by star-formation are very faint.It is necessary to detect sub-µJy sources to account for most of the star formation through z ∼ 1 − 3 when galaxies built up the majority of their stellar mass.Using the MeerKAT DEEP2 field-the deepest ∼1.4 GHz radio image yet taken- Matthews et al. (2021b) measured brightness-weighted radio source counts S 2 n(S) down to 0.25 µJy.As shown in Condon & Matthews (2018), the brightness-weighted source counts S 2 n(S) of either SFGs, AGNs, or both are proportional to the luminosity function of the respective population integrated over all redshift: (4) where D H0 ≡ c/H 0 is the Hubble distance, α ≡ +d ln S/d ln ν is the spectral index, L ν = 4πD 2 C (1 + z) 1−α S is the luminosity per unit frequency, D C is the comoving distance, E(z) = [Ω m (1 + z) 3 + Ω Λ + Ω r (1 + z) 4 ] 1/2 , and the energy density function is (5) where g(z) and f (z) represent density and luminosity evolution, respectively.If the source counts are computed for SFGs and AGNs separately using independently derived evolutionary functions, the total source counts-the quantity typically observed-are simply the sum of the source counts from the two populations.Matthews et al. (2021a) determined the evolutionary functions f (z) and g(z) such that the local luminosity function ρ dex (L ν |0) evolved backwards matches the observed source counts.At any redshift z, the SFRD is proportional to the product f (z)g(z) through the FIR/radio correlation.In this way, Matthews et al. (2021a) constrained the star formation history of the universe for a global population (i.e.there was no information on individual galaxies) by modeling the source counts down to 0.25 µJy.The resulting SFRD evolution measured from the deep MeerKAT radio observations is not only stronger than SFRD evolution based on UV/optical measurements, but also ≳50% stronger than combined UV-IR measurements across all redshiftsdeepening the divide between longer and shorter wavelength pictures of the star formation history of the universe.
Upgrades to the Very Large Array and the advent of science observations with the MeerKAT telescope have unlocked the potential to use radio continuum as a star formation tracer to cosmologically significant redshifts.Recent catalogs of radio continuum sources have been used to probe the SFRD, and there is an emerging scatter amongst the measurements and their implications for the star formation history of the universe.While some radio studies find agreement with UV-IR measurements of SFRD(z) (e.g.Novak et al. 2017;Ocran et al. 2020;van der Vlugt et al. 2022), others have found an increase over the UV-IR measurements (e.g.Leslie et al. 2020;Matthews et al. 2021a;Enia et al. 2022;Cochrane et al. 2023).These studies illustrate many of the difficulties in using radio continuum to trace star formation across redshift: stacking is usually necessary to reach the flux density sensitivity needed to detect star-forming galaxies (SFGs), relying on photometric redshifts to derive luminosity functions and galaxy characteristics, and observing at lower frequencies where the relationship between radio luminosity and SFR is less understood.
This work introduces a comprehensive multiwavelength follow-up campaign of the MeerKAT DEEP2 field.By surveying the same field whose source counts implied stronger SFRD evolution than UV-IR measurements we minimize systematics and robustly test the evolutionary models of Matthews et al. (2021a).Precise redshifts from fitting a combination of low-resolution spectra and optical/NIR photometry of 3,839 galaxies with corresponding MeerKAT detections confirm the cosmic radio SFRD modeling from the source counts and its discrepancy with the UV-IR picture.The results presented here establish a rich framework for dustunbiased tests of SFR evolution and its diagnostics down to faint, normal galaxies responsible for amassing most of the stellar mass in the universe.
This paper is organized as follows.Section 2 details the multiwavelength, spectrophotometric data of the MeerKAT DEEP2 field.Redshifts are derived from a novel SED fitting technique described in Section 3. We calculate radio luminosity functions from 0.2 < z < 1.3 and compare them with models in Section 4. We discuss the implications of confirming radio-based models of stronger SFRD evolution in Section 5.
Absolute quantities were calculated for the flat ΛCDM universe with H 0 = 70 km s −1 Mpc −1 and Ω m = 0.3.
Our spectral-index sign convention is α ≡ +d ln S/d ln ν.

MULTIWAVELENGTH OBSERVATIONS OF THE MEERKAT DEEP2 FIELD
The MeerKAT DEEP2 field was observed as part of commissioning the MeerKAT array (Mauch et al. 2020).The field covers a Θ = 69.2′ diameter halfpower circle centered on J2000 α = 04:13:26.4,δ = −80:00:00 with a θ 1/2 = 7. ′′ 6 circular Gaussian synthesized beam.The location of the field was optimized to avoid bright radio sources whose large flux densities-combined with systematic position and gain uncertainties of the telescope-limit the achievable dynamic range.About 160 hours of integration resulted in a final thermal noise of σ n = 0.56 ± 0.01 µJy beam −1 .The wideband DEEP2 image is the average of 14 narrow subband images weighted to maximize the signalto-noise (S/N) of sources with spectral index α = −0.7,resulting in an effective frequency of ν = 1.266GHz.The confusion distribution of faint sources has such a long tail that it is not well-described by its rms.Using the traditional definition of the confusion limit as the flux density at which there are ≳25 beam solid angles per source, the "rms" confusion noise is σ c ≈ 2.6 µJy beam −1 .Matthews et al. (2021b) describes the construction of the MeerKAT-DEEP2 radio source catalog of 17,350 components down to 10 µJy.In brief, Matthews et al. (2021b) applied the Obit (Cotton 2008) task FndSou, which decomposes islands of contiguous pixels into circular Gaussian components.Most sub-mJy radio sources have angular diameters ϕ ≪ 1 ′′ (Cotton et al. 2018;Murphy et al. 2011), so all sources were treated as unresolved point-sources.This approximation was supported through qualitative and quantitative comparisons of point-source-only simulated images with the observed MeerKAT-DEEP2 field.The radio component catalog extends down to S 1.266 GHz = 10 µJy and includes corrections to positions and flux-densities due to the presence of confusion.Matthews et al. (2021b) used the confusion amplitude distribution to constrain the counts of unresolved sources as faint as S = 0.25 µJy.We refer the reader to Matthews et al. (2021b) for more details.
Here we present initial results from a suite of photometric and spectroscopic observations of the MeerKAT DEEP2 field from the optical through near-infrared.A detailed account of the multiwavelenth observations and a corresponding catalog will be available in Matthews et al. 2024 (in prep).Below is a brief summary of the multiwavelength data contributing to the initial results presented here.

Optical Prism Spectroscopy
Obtaining redshifts for thousands of galaxies typically implies a photometric redshift approach; model or empirical template SEDs are fitted to fluxes in broadband filters to determine the physical properties of galaxies.With an effective spectral resolution of R < 10, current and near-future large imaging surveys achieve uncertainties of at best σ z ≳ 0.025(1 + z) (see Newman & Gruen 2022, for a review of the topic).
A spectrophotometric approach adds low-resolution (R ∼ 30) optical prism spectra-which can be obtained for thousands of galaxies in one observation-to the fluxes in broadband filters.The continuous wavelength coverage in the optical with broadband filters extending into the NIR enable redshift uncertainties σ z ≤ 0.01(1 + z), a method developed and implemented successfully in Patel et al. (2009); Coil et al. (2011); Kelson et al. (2014).We adopt this spectrophotometric approach needed to accurately characterize the ∼17,000 galaxies in our radio-selected sample.
Fourteen pointings of the IMACS f/2 camera with the UDP covered all of the MeerKAT DEEP2 half-power circle (the dashed line in Figure 1).Each of the fourteen masks contained between 1,100 and 1,500 objects, with several objects duplicated between masks of overlapping regions.We demanded that each radio target have a counterpart in the 3.6 µm Spitzer image (described in Section 2.3) to register the radio positions to optical/IR astrometry.In total, 11,671 unique radio sources were observed and are identified with yellow crosses in Figure 1.The median exposure time of the final sample is 12,560 seconds with a minimum exposure time of 4,000 seconds and 95% of the objects having exposure times greater than 7,200 seconds.
To guarantee that most of the slits placed on radio sources returned optical spectra, the initial round of objects (the first 7 of 14 masks) were chosen from the subset of the Matthews et al. (2021b) MeerKAT DEEP2 radio source catalog that had Spitzer cross-identifications within 1 ′′ .This criterion applied for all of the first 7 of 14 masks and 2/3 of the objects on the last 7 of 14 masks.For the remaining 1/3 of objects on the last 7 of 14 masks, we relaxed the criterion to include a random selection from a preliminary "deblended" radio catalog made using the XID+ algorithm by Hurley et al. (2017).Since these objects make up ≤ 15% of the targets (and an even smaller fraction of the reduced spectra) details of the deblending and cross-identifications will be described in a future work.
Observations were taken across 8.5 nights spanning February 2022 through January 2023.Typical exposure times for each mask ranged between 160 and 210 minutes.

Ground-based Optical and Near-Infrared Photometry
As part of programs 2022A-771331 and 2022B-132648, the MeerKAT DEEP2 field was observed by the Dark Energy Camera (DECam) in filters ugrizY .Observations spanned 3.5 nights from February 2022 through December 2022.We used the DECam Community Pipeline stacked image product for photometric measurements (Valdes et al. 2014).
Poor atmospheric conditions persisted through most of the observations and limited the resulting seeing FWHM to ≳ 1. ′′ 2 in all but i band (where it reached an average seeing of ∼ 0. ′′ 9).This significantly degraded the achieved depth of the images, resulting in ∼5σ limiting AB magnitudes of 24.5, 25, 24.75, 24.5, 24.0, and 23.0 for ugrizY , respectively.
Observations were taken in the J-band filter (1.1-1.4 µm) over 2.5 nights from February 2022 through December 2023.The average seeing over multiple nights of observation was 1. ′′ 1.Each pointing was dithered 9 times in a square pattern around the central position with dither offsets of 26 ′′ to cover the 19 ′′ gap between the detectors.Data reduction and image processing was done using the custom pipeline FourCLift developed and described in Kelson et al. (2014).
2.3.Warm-Spitzer 3.6 µm and 4.5 µm imaging The observations were carried out in twelve distinct epochs (Astronomical Observation Requests or AORs) from 8 July 2019 through 11 December 2019 and totaled 69.4 hours (Spitzer PID:14246).Each AOR consisted of a 14×14 mosaic of 100 s frames, using a single point from the cycling dither pattern to break up the mosaic grid pattern.We included small (∼ 10 ′′ ) offsets in the central position, and, because the field is close to the southern continuous viewing zone (CVZ), there was a spread in position angle of the individual mosaics.All together, this ensured full coverage of the circular field containing the FWHM of the MeerKAT primary beam sensitivity (although, sources can be detected beyond this).The reduction of the IRAC photometry was completed in the standard way using the MOPEX data reduction package (Makovoz & Khan 2005).

Photometric Measurements
All astrometry was registered to Gaia DR3.Our deepest optical data were taken in the DECam i filter, so we used this image to detect optical counterparts in the MeerKAT DEEP2 field.We used Source Extractor (Bertin & Arnouts 1996) to calculate and subtract the background and subsequently detect sources using the default parameters.We ran aperture photometry on the resulting ∼ 10 6 detections with a 4 ′′ diameter aperture to ensure all source flux of galaxies at intermediate redshifts and beyond is encompassed.The magnitude zeropoints in grizY J were simultaneously fit to match the synthesized stellar locus according to the algorithms in Kelly et al. (2014).The u-band zeropoint was then calibrated using the "blue tip" of halo stars as described in Liang & von der Linden (2023).

AN EIGENVECTOR APPROACH TO SED MODELING
Constraining galaxy evolution necessitates accurate physical parameters from SED fitting.The accuracy and complexity can come at the expense of longer computation times, particularly for large samples of thousands of galaxies.We developed a new algorithm for SED fitting based upon the fact that SEDs can be constructed from a linear combination of basis functions.The details of this method will be fully described in Kelson et al. (in prep).In summary, we choose to abandon the idea that each basis function has a physical meaning and rather utilize the basis functions as eigenvectors that can be combined in various magnitudes to match a wide variety of galaxy SEDs.

Construction of Eigenvectors
We first establish a grid of redshift, metallicity, and reddening values.There are 501 redshifts at intervals of ∆ log z = 0.005 between −2 ≤ log z ≤ 0.5, six metallicity values spaced at 0.3 dex intervals from -1.2 to 0.3, and eleven increments of log A V spaced at 0.25 dex intervals between −2 ≤ log A V ≤ 0.5.Each point in the redshift-metallicity-reddening grid has 100,000 star-formation histories generated stochastically via the method introduced in Kelson (2014).We constructed our model SEDs using the Flexible Stellar Population Synthesis (FSPS; Conroy et al. 2009;Conroy & Gunn 2010) models with the Chabrier IMF (Chabrier 2003).and adopted the reddening law of Cardelli et al. (1989), which gave the best results (on average).In future iterations of this novel fitting method, we will be incorporating finer increments in metallicity and reddening, and exploring variation in the reddening law.
Singular value decomposition (SVD) of the resulting galaxy SEDs produces a plethora of basis functions, but many fall below observational noise and are discarded (for example, the signal-to-noise of the present sample requires three to six).To capture the diversty of stellar populations to a part in a thousand requires 9 eigenvectors.The spectroscopy-and photometry-here, like most survey data, are not accurate to a part in a thousand with a full accounting of flux calibration and photometric uncertainties, allowing us to restrict the fits to fewer basis functions without any loss of information.

Data Fitting
Every observed galaxy SED is fit using a linear superposition of up to N ≤ 9 eigenvectors plus five nonnegative emission line vectors.The maximum number of eigenvectors N per galaxy is determined when the reduced χ 2 is no longer improved by adding more eigenvectors.The superposition of eigenvectors accurately recreates the full diversity of galaxy stellar populations and we include the following five (non-negative, Gaussian) emission lines in order to capture the full diversity of galaxy SEDs (see, e.g., Fig. 2): Hα, [N II], [O III] (a superposition of 5007 Å and 4959 Å in the ratio 3:1), Hβ, the [O II] doublet at 3727 Å (treated as two equal-magnitude Gaussians separated by 2.7 Å in the restframe), and Mg 2959 Å.
The final number of source SEDs fit using our novel eigenvector approach is 9,416, down from the original 11,671 due to criteria on deep photometry in multiple bands.Specifically, detection in [3.6], J-, and i-band were required as well as measurements in r, i, and z brighter than 25, 24.5, and 24 AB mag, respectively, so as to ensure that the individual flux-calibration polyno-mials could be well constrained.After SED-fitting, we imposed additional constraints to be conservative in this first analysis of the data.These constraints included: (a) the 25th and 75th percentile of the S/N spectrum be greater than 0.25 and 1.0, respectively-this demands that the spectrum is not overly corrupted by random data issues that can corrupt sky subtraction (e.g., slit edge artifacts or scattered light from neighboring alignment star); (b) the reduced χ 2 at a given S/N be less than 2 standard deviations from the median at that S/N ratio-to keep objects within the nominal data quality distribution; (c) objects do not lie on or near the stellar locus, as some stars were not excluded from the catalog prior to observation; (d) objects have narrow P (z) distributions (we consider objects with overly broad P (z) to have had too low S/N to be considered information of value; and (e) that the reduced χ 2 < 20, which not only removes bright stars, but also high-redshift quasars.
After these conservative cuts, our final sample is 3,839 galaxies.Four sources are shown for context in Figure 2. We are encouraged by the strong combination of efficiency and accuracy using this method.

Redshift Uncertainties
The location of the MeerKAT-DEEP2 field was chosen for its potential to achieve unparalleled sensitivityregardless of the availability of multiwavelength data.Unfortunately, being a field not previously covered by wide-field legacy surveys means that there is no extensive data set with high-resolution optical spectra or redshifts.Instead, we derive uncertainties using the large number of objects with repeat measurements.
We re-ran the SED fitting using UDP spectra of individual mask observations rather than using the combined spectra for objects with repeated measurements.There are 1,470 objects with repeated observations.We examined the biweight scatter in duplicate redshift measurements in redshift slices and find that the redshift errors are ≲ 2% in (1 + z) for all observations.The vast majority of combined spectra of sources used in the subsequent analysis have S/N > 10.If we limit the duplicate measurements to this regime, the median redshift error is σ z ∼ 0.006(1 + z).
The tails of the distribution of the redshift differences between two observations encode information on the catastrophic failures.The total fraction of objects where the duplicate redshift is more than 10% off in (1 + z) is 5%.For duplicate observations with S/N > 10representative of the deeper spectra when all the data have been combined-the catastrophic failure rate is < 4% for galaxies with 0.4 < z < 1.In our lowest and highest redshift slices it is < 8%-a result of the poor seeing, and resulting insufficient depth of the DECam imaging.

RADIO SOURCE EVOLUTION
Comparisons of the space densities of galaxies as functions of luminosity (i.e.luminosity functions) at different points in time trace cosmological evolution.Models correctly describing this evolution thus constrain the buildup of stellar mass over cosmic time.Below, we describe the completeness of our sample-in both radio (Section 4.1) and optical (Section 4.2)-and the corrections we made to overcome inevitable incompleteness.Finally, in Section 4.3 we recount our luminosity function derivation and discuss the effects (or lack thereof) due to cosmic variance (Section 4.4).

Radio Reliability and Completeness
The sensitivity of the MeerKAT DEEP2 image is limited by point-source confusion (rms σ c ∼ 2.6 µJy every-where), not by noise (rms σ n = 0.56 µJy at the pointing center increasing to 1.12 µJy at the primary beam halfpower circle).Catalogs of individual sources fainter than S ∼ 10 µJy ∼ 4σ c become increasingly incomplete and unreliable.In the case of MeerKAT DEEP2, the catalog was cut off at S 1.266 GHz = 10 µJy (∼18σ n ) to ensure high completeness and reliability as determined through injecting sources into mock images (Matthews et al. 2021b).Because the limiting flux density of the radio source catalog is independent of the thermal noise, the catalog sensitivity and completeness is uniform across the image.This translates to exceedingly simple radio completeness corrections that depend on flux density alone.
Ten mock images of the MeerKAT DEEP2 field were created by Matthews et al. (2021b) and cataloged with the same source finding algorithm employed on the real data.The completeness of the real source catalog was determined by calculating the fraction of input sources that were recovered and cataloged by the source-finding algorithm.The catalog at the faint limit of S 1.266 GHz = 10 µJy is 57% complete and quickly jumps to ≳ 96% complete by S 1.266 GHz = 20 µJy.
The radio sample reliability 0 ≤ R ≤ 1 is less than one because some of the sources with measured S 1.266 GHz ≳ 10 µJy are actually fainter than the sample limit and should not be in the sample.This is true for any radio catalog, but is of particular concern for confusion-limited images.Using the ten mock images presented in Matthews et al. (2021b), we estimate the reliability in 0.2 dex wide bins of log S centered on log S 1.266 GHz = −4.9,−4.7, −4.5, • • • .For cataloged sources in each log S bin, we calculate the fraction of the cataloged sources whose input source flux density is truly greater than 9.45 µJy (thereby excluding objects whose flux was enhanced by more than the rms thermal noise).We add an additional constraint that for each cataloged source in a log S bin, the input source flux density truly lies within the log S bin in question.In this way, we estimate not only the reliability that a source is truly above the sample limit, but also the reliability that sources truly belong in their resident flux density bin (i.e.accounting for situations when the true flux density would place the source in the next highest flux density bin, but the cataloged flux was underestimated (e.g.due to negative sidelobes).The radio source completeness and reliability is shown in Figure 3.
The radio catalog limit of S 1.266 GHz = 10 µJy sets the minimum luminosity-and therefore the minimum SFR-this survey is sensitive to as a function of redshift.Assuming our radio sources are drawn from the population of normal star-forming galaxies and that SFR is correlated with stellar mass, we can use the parameterization of Whitaker et al. (2012) to find the stellar mass that corresponds to our survey flux limit as a function of redshift (shown in Figure 4).We assume the characteristic mass of the galaxy stellar mass function is log[M * (M ⊙ )] ∼ 10.85, the middle of the range log[M * (M ⊙ )] ∼ 10.7 − 11 reported by Weaver et al. (2023) for the COSMOS2020 sample.At a z = 1.3, our radio flux density limit corresponds to a galaxy with stellar mass ∼ 0.63 M * , meaning our completeness of normal galaxies drawn from a representative sample approaches ∼ 50% at z = 1.3.Galaxies below this characteristic mass threshold do make it into the sample, but only if they have abnormally higher SFRs compared to the locus defined in Whitaker et al. (2012).By redshifts z ∼ 1.3, galaxies enter the sample in ways that are not representative of the full distribution of galaxies that fully contribute to the SFRD.

Optical Completeness
The poor seeing prevented our DECam observations from reaching the intended i = 25 magnitude limit.Nonetheless, because the seeing and depth of the iband data are superior to the other bands, we use it to define the positions of galaxies associated with the radio detections for the photometry behind the SEDs.Cross-identifying optical counterparts in the i band image yields 14,679 matches within 2. ′′ 5 (smaller than half the 7. ′′ 6 FWHM of the MeerKAT-DEEP2 beam) of the 17,150 radio galaxies.The Poisson probability that one or more unrelated radio sources lie within 2. ′′ 5 of an i band source is where ρ is the sky density of sources with flux density brighter than a given value (N (> S)) and r s is the crossmatching radius.For the sky density of sources brighter than S 1.266 = 10 µJy, P (≥ 1) ≈ 0.028.The rate of false matches is not expected to impact the results.The fraction of these 14,679 galaxies with photometric measurements and successful SED fitting represents the optical completeness of the radio sample.3,839 galaxies had both IMACS spectra free of data artifacts, high enough S/N, sufficient spectroscopic and photometric data quality and fidelity, and photometric measurements in a majority of the filters to fit an SED.We partitioned the 14,679 galaxies into log S bins of width 0.2 from log[S (Jy)] = −5 to log[S (Jy)] = −2.For each bin, the completeness is defined as the ratio of the number of galaxies with a successful SED fit to the total number within each bin.The resulting optical completeness curve is shown in Figure 3.We find a slight decreasing trend of optical completeness with increasing radio flux density, ranging from ∼40% complete at the lower limit log[S (Jy)] = −5 to ∼20% at log[S (Jy)] = −3.1.
Comparing the optical completeness with i-band magnitude, we are most complete (∼44%) at a i = 21 and the completeness falls near-linearly to ∼5% at i = 24.We limit our luminosity function calculations to i = 23.5, the point at which the completeness falls to 12.5%, just over a quarter of its peak value.

Radio luminosity functions
We calculate luminosity functions in several contiguous redshift slices using the 1/V max method (Schmidt 1968).We define redshift bins to be wide enough such that the Poisson counting errors in that redshift slice are at most ∼ 5% (N ≈ 350), but narrow enough to mitigate the impact of evolving stellar populations.If redshift slices are too wide, some galaxies that may have sufficient photometry at the front end of a slice would have had photometric properties-owing to both mass growth and stellar evolution-at the back end of the redshift slice that leave them excluded by the DECam imaging depth and complicate the estimation of V max .Our sample has a median redshift error σ z ≲ 0.01(1+z), which means the average uncertainties on z are at most ∼ 17.5% (for redshift bin 0.7 < z < 0.8) of the width of the redshift bins, meaning objects are unlikely to scatter into neighboring redshift slices at a level affecting the luminosity function or SFRD measurement.
The maximum comoving volume in which a galaxy can be observed V max = V (z max ) − V (z min ) is the volume of  the redshift bin over which the source is observable.It is bounded on one side by z min , the lower bound of the redshift bin containing the galaxy.On the high end, the maximum redshift z max is the minimum of the following: (1) the redshift at which the galaxy falls out of the radio sample (S 1.266 GHz < 10 µJy), (2) the redshift at which the galaxy falls out of the optical sample (i > 23.5-and assuming a galaxy's color does not evolve over the time span of its respective redshift slice), or (3) the maximum redshift of the redshift slice containing the galaxy.
The spectral luminosity L 1.4 GHz at redshift z depends on the observed flux density of the source and the spectral index of the source population where we assume the standard α = −0.7 for radio populations at ∼1.4 GHz.
Unlike in the radio regime, the SED of a galaxy in the optical does not follow a power law and depends sensitively on the K-correction.
where m i is the apparent magnitude of the source at redshift z, M i is the i band absolute magnitude, DM(z) is the distance modulus, and K i (z) is the K-correction at i band.We modeled the K-correction using the kcorrect software by Blanton & Roweis (2007).We assume the galaxy SED-the best-fit spectral template from kcorrect-does not change within its redshift bin and calculate the maximum redshift for which m i ≤ 23.5.
For each luminosity bin, the space density of radio sources is calculated by the following: where ∆ log L is the width of the luminosity bin, V max,i is the maximum comoving volume over which the ith galaxy could be observed, where Ω is the survey area in steradians, and C i is the completeness correction factor (equal to 1/completeness) for the ith galaxy.Luminosity bins are of width 0.2 dex centered on log[L 1.4 GHz (W Hz −1 )] = 21.1, 21.3, . . ., 25.3.The scenario in which our sample is 100% completewhere we have measured redshifts for every radio source at z < 1.3 and that all other sources are at z > 1.3is illustrated by the open circles in Figure 5. Assuming all unmeasured objects are at z < 1.3 overestimates the incompleteness and leads to measured luminosity functions far above the predicted luminosity functions calculated using the radio evolution models of Matthews et al. (2021a).Here, we presume completeness corrections assuming that the true fraction of radio sources at z < 1.3 is that implied by the modeling from Matthews et al. (2021a).As a function of flux density, the fraction of all radio sources with z < 1.3 ranges from 52% at log[S (Jy)] = −5, rises to a peak of 81% at log[S (Jy)] = −3.5, and plateaus around ∼70% for sources with flux densities approaching 100 mJy.The completeness correction factor c i for the ith galaxy includes a term for the reliability of the radio catalog as a function of flux density R r (S), for the completeness of the radio survey as a function of flux density C r (S), for the optical/redshift completeness of the radio sample as a function of flux density C opt (S), and for the fraction of objects at flux density S that lie within z < 1.3 by the evolutionary models of Matthews et al. (2021a).The final completeness correction factor for the ith galaxy is as follows: The errors are the quadrature sum of the expected fractional cosmic variance and rms Poisson counting errors for independent galaxies.If the number of galaxies in a luminosity bin is small (N < 5), the counting errors are taken from the 84% confidence limits tabulated in Gehrels (1986).The uncertainty caused by cosmic variance is described below in Section 4.4.

Cosmic Variance
We estimate cosmic variance using three independent methods: with the Cosmic Variance Calculator from Trenti & Stiavelli (2008), from the radio sky simulations of Bonaldi et al. (2019), and following the results of Driver & Robotham (2010).Trenti & Stiavelli (2008) model the cosmic variance using a combination of analytic estimates via the two-point correlation function of dark-matter halos in an extended Press-Schechter formalism (Press & Schechter 1974), and numerically using synthetic catalogs from N -body simulations of largescale structure formation.To connect the cosmic variance of the radio source population to the predicted clustering of dark matter halos from Trenti & Stiavelli (2008), we use the stellar-to-halo mass relation presented in Leauthaud et al. (2012) and estimate the average halo mass associated with the stellar-masses of the radio sources determined from the SED modeling, which range from ⟨log M * ⟩ = 10 to ⟨log M * ⟩ = 10.5 from the lowest to highest redshift slice.
With the estimated average halo masses corresponding to our radio sample, we use the Cosmic Variance Calculator to estimate the 1σ fractional cosmic variance in each redshift bin (Trenti & Stiavelli 2008).For the survey area of the MeerKAT-DEEP2 field (∼1 deg 2 ), the fractional count error due to cosmic variance is consistently less than 0.16 across all redshift slices in our sample.Tweaking the Cosmic Variance Calculator parameters (e.g.halo filling factor) has little effect on the output cosmic variance values and therefore inconsequential effects on the space densities of radio galaxies determined in this work.
We independently calculate the cosmic variance using the Tiered Radio Extragalactic Continuum Survey (T-RECS) presented in Bonaldi et al. (2019).T-RECS models two main populations of radio sources: SFGs and AGN over the 150 MHz to 20 GHz range.From the 25 deg 2 "medium" deep tier, we extract nonoverlapping circular sample areas of ∼1 deg 2 -the size of the MeerKAT-DEEP2 field.This amounts to 16 independent sky areas.In each simulated sky area, we combine both SFGs and AGNs-as we are currently unable to distinguish between the two populations in our observations of the MeerKAT DEEP2 field-and calculate number counts in the same 1.4 GHz luminosity bins used to construct the observed luminosity functions.We find that the cosmic variance is consistent with the values estimated by the Cosmic Variance Calculator through z = 0.7 and falls below the Cosmic Variance Calculator estimations for z > 0.7.For z > 0.7, the cosmic variance determined by comparing the number counts in the 16 sky areas of the Bonaldi et al. (2019) simu-lations is ≲ 10% and contributes only 0.04 dex to the error budget.
Finally, we refer to the cosmic variance estimates of Driver & Robotham (2010), who used the Sloan Digital Sky Survey (SDSS) to empirically develop formulae for estimating cosmic variance based on survey shape and volume.It is important to note that the estimations from Driver & Robotham (2010) were derived for a M * ± 1 mag population of galaxies.This assumption captures most-but not all-of the sample of radio sources we present here.For the median redshift in each slice, we calculate the comoving transverse distance of the MeerKAT DEEP2 field radius and the comoving radial depth for the width of that redshift slice.Using these values and Equation 4 of Driver & Robotham (2010), we calculate the percent cosmic variance in each redshift slice.For bins of ∆z = 0.1 dex, the cosmic variance ranges from 24-28%-much higher than that predicted by the other two methods.Wanting to be conservative with our systematic uncertainties, we adopt Equation 4 of Driver & Robotham (2010) to estimate the cosmic variance in our SFRD measurements and increase the minimum width of our redshift bins to ∆z = 0.2 dex for the SFRD calculation.

Estimating radio-based SFRD(z)
Converting radio luminosities to SFRs deserves extensive care and a thorough investigation of possibly uncertainties (as will be outlined in Section 5).However, to guide such discussion and quantify the possible disagreement in the SFRD(z) between radio and UV-IR studies, we present an initial conversion between the calculated radio luminosity functions and SFRD(z).
The SFRD at any redshift z is related to the total radio energy density produced by star-forming galaxies.In each redshift slice, we fit the projected luminosityweighted luminosity functions (e.g.energy density functions) of SFGs and AGNs combined to the observed data.As an initial estimate, we assume the form of the luminosity functions modeled by Matthews et al. (2021a) and only allow a renormalization of the SFG energy density function with respect to both axes (i.e. a scale and a shift).
To properly compare the observed SFRD in each redshift slice with that predicted by the Matthews et al. (2021a) evolutionary models, we adopt their prescription to convert radio luminosity to SFR.In short, we assume the relationship between SFR and integrated IR luminosity (8 < λ(µm) < 1000) of Murphy et al. (2011) with the average ratio between total IR luminosity and FIR luminosity (42.5 < λ(µm) < 122.5) measured by Bell et al. (2003).The conversion from the total radio spectral power at 1.4 GHz at any time t to SFRD ψ(t) for a Kroupa (2001) IMF is where ⟨q(L ν )⟩ is the ratio of the FIR and radio luminosities from Matthews et al. (2021a): We calculate U 1.4 GHz directly by summing L 1.4 GHz /(C V max )-where C is the same completeness correction factor defined in Section 4.3-over the unbinned sample of galaxies within each redshift slice (0.2 < z < 0.4, 0.4 < z < 0.6, 0.6 < z < 0.8, 0.8 < z < 1.0, and 1.0 < z < 1.3).We restrict the summation to galaxies with L 1.4 GHz ≤ 24.25 W Hz −1 , the luminosity corresponding to an SFR of ∼ 1000 M ⊙ yr −1 .This upper limit is chosen to minimize AGN contamination.
Our calculation of U 1.4 GHz is independent of any assumptions about the shape or evolution of the radio luminosity function, providing an essential test on the radio luminosity and density evolution models derived in Matthews et al. (2021a).However, the imposed flux limit corresponds to a different minimum luminosity in each redshift slice and probes different fractions of the total energy output by star-forming galaxies.To homogenize our SFRD measurements, we use the models derived in Matthews et al. (2021a) to estimate the fraction of the total energy-density our survey is sensitive to for each redshift slice.We divide by this fraction-ranging from 70% in the highest redshift slice to 92% in the lowest-to correct for this incompleteness.For transparency, we include three versions of our radio-based SFRD measurement in Figure 6 at the median redshift of each slice: (1) the raw measurements, uncorrected for incompleteness of any kind, (2) the measurements corrected for the completeness and reliability of the optical and radio surveys, and (3) the corrected measurements homogenized to probe an equivalent fraction of the total energy density.
There is excellent agreement between our independent SFRD measurements and both the SFRD predictions from the models in Matthews et al. (2021a) and SFRD measurements from other radio-based studies (Cochrane et al. 2023;Enia et al. 2022;Leslie et al. 2020) in the range 0.2 < z < 0.7.However, there is notable disagreement between the SFRD measurements presented here with other radio-based studies (Novak et al. 2017;Ocran et al. 2020;van der Vlugt et al. 2022).Further work to (1) push to fainter flux densities (and therefore further down the faint end of the luminosity function) and (2) carefully separate SFGs and AGNs in the MeerKAT DEEP2 sample will expand the integration limits and place stronger constraints on the radio-SFRD(z).

POTENTIAL SOURCES OF DISCREPANCY IN THE SFRD EVOLUTION AND IMPLICATIONS
At the end of the previous section, we showed star formation rate densities-inferred from our luminosity functions-in Figure 6.These SFRDs-shown by the filled triangles-are compared to those from other radio studies (e.g.Leslie et al. 2020;Enia et al. 2022;Cochrane et al. 2023), and to the Madau & Dickinson (2014) model based on a large sample of UV-IR estimates.While this study advances our understanding of the radio-SFRD by coupling the the unprecedented sensitivity of the MeerKAT DEEP2 radio image with a large spectrophotometric sample, the tension with UV-IR measurements has been seen across multiple radio studies (such as those listed above).Given that the luminosity functions agreed with those predicted by Matthews et al. (2021a), Figure 6 also shows that our SFRD measurements agree with the SFRD evolution of the Matthews et al. (2021a) modeling.
The radio values for SFRD between (0.2 ≤ z ≤ 1.3) are 50% to 100% higher than canonical estimates from not only UV/optical, but also, seemingly, over combined UV-IR measurements.This discrepancy suggests more than a mild tension, but rather a fundamental disagreement, in one of the most important metrics for understanding the evolution of our universe.Further, there is emerging scatter amongst the radio-based SFRD measurements-with some in agreement with UV-IR measurements (e.g.Novak et al. 2017;Ocran et al. 2020;van der Vlugt et al. 2022) and others significantly discrepant in all or selected redshift regimes (e.g.Cochrane et al. 2023;Matthews et al. 2021a).Here we discuss sources of systematic uncertainty on measurements of luminosity and space densities, quantify their impact on SFRD calculations, and explore possible reasons for disagreement between radio and UV-IR based SFRD(z).

Budget of systematic uncertainties
A full budget of systematic uncertainties must include many potential issues, some of which have been mentioned already.In Section 4.4 the uncertainty due to cosmic variance was estimated using three independent methods, with resulting fractional cosmic variance uncertainties ranging from 23% in the lowest-z bin to 14% in the highest-z one, according to the most conservative of the three methods (Driver & Robotham 2010).Other potential sources of uncertainty include confusion effects, completeness correction factors, contamination by AGN, and the FIR/radio correlation parameter q.
The volume surveyed here, and the modest 7. ′′ 6 PSF, are both large enough to produce more than zero incidents of coincidence between galaxies and background radio sources along the line of sight.We generated mock MeerKAT images by populating the image with sources at redshifts and luminosities drawn from the model radio luminosity functions of Matthews et al. (2021a).After a blind cataloging-without prior knowledge of injected source positions-of the image, we simulated the observed luminosity functions.Comparing the mock catalog with the positions and luminosities of input sources, we determined the occurance rate and effects of object coincidence per luminosity bin and as a function of redshift.At z < 0.7 and z > 1 the effect is negligible, while at 0.7 < z < 0.9 the observed luminosity functions may be contaminated by line of sight superposition at a level that can, on average, skew the luminosities of sources by +0.05 dex.
There is an additional systematic uncertainty that arises from the modeling of spectroscopic incompleteness, as we have only securely measured redshifts for a fraction of the radio catalog.We estimate this uncertainty by summing the galaxy weights in our final catalog of 3, 839 galaxies and comparing to the size of what ought to be the parent catalog-the number of radio galaxies in the MeerKAT DEEP2 catalog multiplied by the fraction of those expected to have z < 1.3.This latter term is model dependent and we adopt a flux-density dependent fraction from Matthews et al. (2021a).We find the sum of galaxy weights (i.e.inverse completeness) to be within ∼ 10% of the size of our parent catalog.As such, we estimate errors in the completeness correction factor to add ∼10% uncertainty (0.04 dex) to our SFRD measurements.
The fraction of AGN in our radio sample will decrease with decreasing flux density.Further, since AGN components are omitted in our SED fitting method, AGNs exhibiting strong spectral features in the optical will result in poor SED fits and are therefore cut from our sample following our conservative data-quality criteria.Nonetheless, radio emission powered by AGNs is present in our sample at some level.After identifying AGNs with a wide variety of multiwavelength diagnostics, Algera et al. (2020) found that at S 3 GHz = 30 µJy the fraction of radio sources powered primarily by starformation reaches 90%.At 1.266 GHz, S 3 GHz = 30 µJy equals ∼ 55 µJy (assuming α = −0.7).While 82% of our sample lie below this flux density value, we conservatively assume 10% of our radio emission is contributed by AGN.As such, AGN contamination add ∼0.04 dex to the total error budget on our SFRD measurements.

Sources of disagreement in SFR conversions
Disagreement between can arise for two reasons: problems in the measurements (UV, IR, and/or radio) themselves, or problems in the conversions from these measurements to SFRs.This work will focus on areas of improvement for radio-based determinations of the SFRD(z).
Particularly important to our work is the calibration of radio continuum as a tracer of SFR through the FIR/radio correlation.We have adopted the sub-linear relation calculated for SFGs in the local universe by Matthews et al. (2021a).Had a luminosity-independent ⟨q⟩ = 2.34 ± 0.01 been adopted from the seminal work of Yun et al. (2001)-using galaxies in the local universe- The blue points are the direct SFRD measurements from this work using the sub-linear FIR/radio correlation from Matthews et al. (2021a).As a comparison, the gold points show the resulting SFRD measurements when the stellarmass-dependent correlation found by Delvecchio et al. (2021) is applied to our sample.
Figure 7 demonstrates the range of resulting SFRD values that can occur when adopting various FIR/radio correlation prescriptions.
Of the more modernand robust-relations for q, the sub-linear relationship (L FIR ∝ L 0.85 1.4 GHz ) presented by Matthews et al. (2021a) and the stellar-mass dependent correlation found by Delvecchio et al. (2021) yield the largest SFRD values, but are consistent with the redshift-dependent q relationships presented by Magnelli et al. (2015) and Delhaize et al. (2017) in 0.2 < z < 1.3.
When added in quadrature the multiple sources of systematic error in our SFRD measurements yield ∼ 0.13 dex.Given the combined random and systematic uncertainties, our measurements substantially agree with the model by Matthews et al. (2021a) of SFRD evolution (the redshift slices at z ∼ 0.7 − 1 yield SFRD measure-ments ∼ 2 − 3σ high) and in significant disagreement with the SFRD characterized by "canonical" UV-IR estimators.We believe cosmic variance in the DEEP2 area is driving the jump in SFRD around z ∼ 0.8.Preliminary analysis of the spatial distribution of galaxies in each redshift slice suggest the presence of multiple galaxy groups or clusters around z ∼ 0.8.A detailed quantification of this overdensity and the implications will be addressed in a subsequent paper.
While some sources of systematic errors cannot be mediated (e.g.cosmic variance), others can be improved.Further analysis of the MeerKAT DEEP2 field will explore the connection between radio luminosity and SFR across a wide variety of galaxy stellar masses and SFRs.Comparison of SFRs derived from SED fitting (some of which are supplemented by emission line fluxes measured directly from the spectra) with those by radio luminosity will illuminate areas of galaxy parameter space where these diagnostics disagree.Particular attention to resolving the cause behind these SFR disagreements will inevitably improve our understanding of dust-reddening evolution, cosmic-ray diffusion, and stellar mass buildup in star-forming galaxies.

SUMMARY
We present initial results from the spectrophotometric multiwavelength follow-up of the MeerKAT DEEP2 field.Using a novel SED-fitting method on a combination of low-resolution spectra and optical/NIR photometry, we determined redshifts for 3,839 galaxies spanning 0.2 < z < 1.3.This sample provides the first tests of critical assumptions-such as the shape of the luminosity function remaining constant with timeand predictions-such as the redshift distributions of radio sources-underlying the modeling of counts by Matthews et al. (2021a).As can be seen from the luminosity functions in Figure 5, the assumed shape of the luminosity function and predicted number density of radio sources was not grossly in error.We conclude with the following: • The luminosity functions measured from the collection of 3,839 individual radio galaxies agree remarkably well with the luminosity and density evolutionary models derived from only global (i.e.lacking any individual galaxy-level information) radio source counts (Matthews et al. 2021a).
• This agreement confirms the findings of Matthews et al. (2021a)-there is strong evidence for ≳ 50% increased SFRD evolution at radio frequencies than what has been measured in the combined UV-IR.The source of this discrepancy remains unknown and is the focus of our continued work with these data.
• Radio measurements of the SFRD show increasing scatter, with some studies confirming enhanced SFRD evolution over UV-IR measurements and others in agreement with the UV-IR SFRD(z) (i.e.Madau & Dickinson 2014).We hope that our continued multi-wavelength follow up of the MeerKAT-DEEP2 radio sourcesspectroscopic and photometric-will provide further insights into the connections between the radio and UV-IR emission from star forming galaxies.

Figure 1 .
Figure 1.The 1.266 GHz MeerKAT DEEP2 field is shown in the left panel-uncorrected for primary beam attenuation-with the half-power circle of the primary beam marked by the white dashed circle.Positions of the radio galaxies targeted for optical prism spectra are shown as yellow crosses.The right panel shows an enlarged version of a central ∼ 3. ′ 4× ∼ 2. ′ 6 region.The initial data release of 3,839 redshifts represents a small fraction of the total radio sources.

Figure 2 .
Figure 2. Four examples of radio galaxy SEDs are shown with their multiwavelength photometry (as 14 ′′ × 14 ′′ cutouts) and redshift probability distribution.Observed spectra and photometric flux values are shown in black.Our best-fitting model spectra and predicted photometric flux values are shown as the magenta line and blue points, respectively.

Figure 3 .
Figure 3. Radio completeness Cr(S), radio reliability Rr(S), and optical completeness Copt(S)/fz<1.3(S)as a function of flux density.The optical completeness includes the correction fz<1.3(S)accounting for the fraction of objects at S that are expected to lie at z < 1.3, for the purpose of avoiding over-correcting for incompleteness and inflating the resulting luminosity functions.

Figure 4 .
Figure 4.The redshift distribution of the final 3,839 MeerKAT DEEP2 sources presented in this work (orange curve).The median stellar mass of a star-forming galaxy with S1.4 GHz = 10 µJy as determined by the star-forming "main sequence" (taken from Whitaker et al. 2012) is shown as a solid blue line.

Figure 5 .
Figure5.Observed number densities of radio sources (filled circles)-corrected for both radio and optical incompleteness-as a function of spectral radio luminosity.Open circles represent the number densities of the same sample without completeness corrections.Error bars reflect 1σ uncertainties.Solid colored line: radio luminosity functions predicted byMatthews et al. (2021a) for all galaxies (starforming and AGN) at the median redshift of the observed galaxies in each redshift slice.The predicted radio luminosity functions of only SFGs or only AGNs at the median redshift are shown as the dashed and dotted lines, respectively.Thick gray line: local radio luminosity function of star-forming galaxies fromCondon et al. (2019).
Figure6.The SFRD from present day to a lookback time of tL 9.5 Gyr.A sample of UV-IR(Madau & Dickinson 2014;Marchetti et al. 2016) and recent radio(Leslie et al. 2020;Enia et al. 2022;Cochrane et al. 2023) continuum SFRD measurements are shown for comparison.SFRD measurements from our sample of 3,839 galaxies in the MeerKAT DEEP2 field are shown as filled, thick triangles.Unfilled triangles show the SFRD measurements from the MeerKAT DEEP2 sample without any corrections for sample incompleteness.The uncertainties in each redshift slice due to cosmic variance-estimated from Driver & Robotham (2010)-are shown as the rectangular outlines around the dark blue triangles.All measurements have been converted to a Salpeter IMF(Salpeter 1955) for ease of comparison with theMadau & Dickinson (2014) relation.

Figure 7 .
Figure 7.The SFRD(z) is shown as a function of lookback time.Different prescriptions for the conversion from radio to FIR/TIR luminosity are explored as different colored curves.The blue points are the direct SFRD measurements from this work using the sub-linear FIR/radio correlation fromMatthews et al. (2021a).As a comparison, the gold points show the resulting SFRD measurements when the stellarmass-dependent correlation found byDelvecchio et al. (2021) is applied to our sample.