# The VLA-COSMOS 3 GHz Large Project: The infrared-radio correlation of star-forming galaxies and AGN to z less than or similar to 6

J. Delhaize, V. Smolcic, I. Delvecchio, M. Novak, M. Sargent, N. Baran, B.

Magnelli, G. Zamorani, E. Schinnerer, E. J. Murphy, et al.

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

► To cite this version:

J. Delhaize, V. Smolcic, I. Delvecchio, M. Novak, M. Sargent, et al. The VLA-COSMOS 3 GHz Large Project: The infrared-radio correlation of star-forming galaxies and AGN to z less than or similar to 6. Astronomy & Astrophysics - A&A, 2017, 602, 10.1051/0004-6361/201629430. [hal-01678393](https://hal.science/hal-01678393v1)

HAL Id: hal-01678393

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

Submitted on 14 Aug 2025

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# The VLA-COSMOS 3 GHz Large Project:
The infrared-radio correlation of star-forming galaxies and AGN
to $z \le 6$

J. Delhaize<fnref n="1" />, V. Smolčić<fnref n="1" />, I. Delvecchio<fnref n="1" />, M. Novak<fnref n="1" />, M. Sargent<fnref n="2" />, N. Baran<fnref n="1" />, B. Magnelli<fnref n="3" />, G. Zamorani<fnref n="4" />,
E. Schinnerer<fnref n="5" />, E. J. Murphy<fnref n="6" />, M. Aravena<fnref n="7" />, S. Berta<fnref n="1" />,*, M. Bondi<fnref n="8" />, P. Capak<fnref n="9" />, C. Carilli<fnref n="10" />,<fnref n="11" />, P. Ciliegi<fnref n="4" />,
F. Civano<fnref n="12" />, O. Ilbert<fnref n="13" />, A. Karim<fnref n="3" />, C. Laigle<fnref n="14" />, O. Le Fèvre<fnref n="13" />, S. Marchesi<fnref n="15" />, H. J. McCracken<fnref n="14" />,
M. Salvato<fnref n="16" />, N. Seymour<fnref n="17" />, and L. Tasca<fnref n="13" />

Received 29 July 2016 / Accepted 17 February 2017

## ABSTRACT

We examine the behaviour of the infrared-radio correlation (IRRC) over the range $0 < z \le 6$ using new, highly sensitive 3 GHz observations with the *Karl G. Jansky Very Large Array (VLA)* and infrared data from the *Herschel Space Observatory* in the 2 deg<fnref n="2" /> COSMOS field. We distinguish between objects where emission is believed to arise solely from star-formation, and those where an active galactic nucleus (AGN) is thought to be present. We account for non-detections in the radio or in the infrared using a doubly-censored survival analysis. We find that the IRRC of star-forming galaxies, quantified by the infrared-to-1.4 GHz radio luminosity ratio ($q_{\text{TIR}}$), decreases with increasing redshift: $q_{\text{TIR}}(z) = (2.88 \pm 0.03)(1+z)^{-0.19\pm0.01}$. This is consistent with several previous results from the literature. Moderate-to-high radiative luminosity AGN do not follow the same $q_{\text{TIR}}(z)$ trend as star-forming galaxies, having a lower normalisation and steeper decrease with redshift. We cannot rule out the possibility that unidentified AGN contributions only to the radio regime may be steepening the observed $q_{\text{TIR}}(z)$ trend of the star-forming galaxy population. We demonstrate that the choice of the average radio spectral index directly affects the normalisation, as well as the derived trend with redshift of the IRRC. An increasing fractional contribution to the observed 3 GHz flux by free-free emission of star-forming galaxies may also affect the derived evolution. However, we find that the standard (M82-based) assumption of the typical radio spectral energy distribution (SED) for star-forming galaxies is inconsistent with our results. This suggests a more complex shape of the typical radio SED for star-forming galaxies, and that imperfect K corrections in the radio may govern the derived trend of decreasing $q_{\text{TIR}}$ with increasing redshift. A more detailed understanding of the radio spectrum is therefore required for robust K corrections in the radio and to fully understand the cosmic evolution of the IRRC. Lastly, we present a redshift-dependent relation between rest-frame 1.4 GHz radio luminosity and star formation rate taking the derived redshift trend into account.

**Key words.** galaxies: evolution – galaxies: star formation – radio continuum: galaxies – infrared: galaxies

## 1. Introduction

A tight correlation between the total infrared luminosity of a galaxy and its total 1.4 GHz radio luminosity, extending over at least three orders of magnitude, has been known to exist

for some time (e.g. van der Kruit 1971, 1973; de Jong et al. 1985; Helou et al. 1985; Condon 1992; Yun et al. 2001). This correlation exists for star-forming late-type galaxies, early-type galaxies with low levels of star formation and even for some merging systems (e.g. Dickey & Salpeter 1984; Helou et al. 1985; Wrobel & Heeschen 1988; Domingue et al. 2005).

<footnote n="1">¹ Department of Physics, Faculty of Science, University of Zagreb, Bijenička cesta 32, 10000 Zagreb, Croatia
e-mail: jacinta@phy.hr</footnote>

<footnote n="2">² Astronomy Centre, Department of Physics & Astronomy, University of Sussex, Brighton, BN1 9QH, UK</footnote>

<footnote n="3">³ Argelander Institute for Astronomy, University of Bonn, Auf dem Hügel 71, 53121 Bonn, Germany</footnote>

<footnote n="4">⁴ INAF-Osservatorio Astronomico di Bologna, via Piero Gobetti 93/3, 40129 Bologna, Italy</footnote>

<footnote n="5">⁵ Max Planck Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany</footnote>

<footnote n="6">⁶ National Radio Astronomy Observatory, 520 Edgemont Rd, Charlottesville, VA 22903, USA</footnote>

<footnote n="7">⁷ Núcleo de Astronomia, Facultad de Ingeniería y Ciencias, Universidad Diego Portales, Av. Ejército 441, Santiago, Chile</footnote>

<footnote n="8">⁸ INAF-Istituto di Radioastronomia, via Gobetti 101, 40129 Bologna, Italy</footnote>

<footnote n="9">⁹ Department of Astronomy, California Institute of Technology, MC 249-17, 1200 East California Blvd, Pasadena, CA 91125, USA</footnote>

<footnote n="10">¹⁰ National Radio Astronomy Observatory, Socorro, NM, USA</footnote>

<footnote n="11">¹¹ Cavendish Laboratory, Cambridge, UK</footnote>

<footnote n="12">¹² Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA 02138, USA</footnote>

<footnote n="13">¹³ Aix Marseille Université, CNRS, LAM (Laboratoire d'Astrophysique de Marseille) UMR 7326, 13388 Marseille, France</footnote>

<footnote n="14">¹⁴ Institut d'Astrophysique de Paris, Sorbonne Universités, UPMC Univ. Paris 06 et CNRS, UMR 7095, 98bis bd Arago, 75014 Paris, France</footnote>

<footnote n="15">¹⁵ Department of Physics and Astronomy, Clemson University, Kinard Lab. of Physics, Clemson, SC 29634-0978, USA</footnote>

<footnote n="16">¹⁶ Max-Planck-Institut für Extraterrestrische Physik (MPE), Postfach 1312, 85741 Garching, Germany</footnote>

<footnote n="17">¹⁷ International Centre for Radio Astronomy Research, Curtin University, Perth, WA 6102, Australia</footnote>

<footnote>* Visiting scientist.</footnote>
The so-called infrared-radio correlation (IRRC) has been used to identify and study radio-loud active galactic nuclei (AGN; e.g. Donley et al. 2005; Norris et al. 2006; Park et al. 2008; Del Moro et al. 2013) and to estimate the distances and temperatures of high-redshift submillimetre galaxies (e.g. Carilli & Yun 1999; Chapman et al. 2005). Another important application of the IRRC is to calibrate radio luminosities for use as indirect, dust-unbiased star formation rate (SFR) tracers (e.g. Condon 1992; Bell 2003; Murphy et al. 2011, 2012). This is particularly relevant considering the powerful new capabilities of the recently upgraded radio astronomy facilities (such as the Karl G. Jansky Very Large Array; VLA) and the next generation of radio telescopes coming online in the near future (such as MeerKAT, the Australian SKA Pathfinder and the Square Kilometer Array). Sensitive radio continuum surveys with these instruments will have simultaneously good sky coverage and excellent angular resolution and will thus have the potential to act as powerful SFR tracers at high redshifts. However, this relies on a proper understanding of whether, and how, the IRRC evolves with redshift.

Star-formation in galaxies is thought to be responsible for the existence of the IRRC, although the exact mechanisms and processes at play remain unclear. Young, massive stars emit ultraviolet (UV) photons, which are absorbed by dust grains and re-emitted in the infrared (IR), assuming the interstellar medium is optically-thick at UV wavelengths. After a few Myr, these massive stars die in supernovae explosions which produce the relativistic electrons that, diffusing in the galaxy, are responsible for synchrotron radiation traceable at radio wavelengths (e.g. Condon 1992). Several theoretical models attempt to explain the IRRC on global scales, such as the Calorimeter model proposed by Voelk (1989), the conspiracy model (e.g. Bell 2003; Lacki et al. 2010) and the optically-thin scenario (Helou & Bicay 1993). Models such as the small-scale dynamo effect (Schleicher & Beck 2013; Niklas & Beck 1997) attempt to explain the correlation on more local scales. However, none of these models successfully reproduce all observational constraints.

As to whether the IRRC evolves with redshift, several different theoretical predictions exist. Murphy (2009) predict a gradual increase in the infrared-to-radio luminosity ratio with increasing redshift due to inverse Compton scattering off the cosmic microwave background resulting in reduced synchrotron cooling, although this is dependent on the magnetic field properties of galaxy populations. Schober et al. (2016) model the evolving synchrotron emission of galaxies and also find a decreasing IRRC towards higher redshifts. On the other hand, Lacki & Thompson (2010) predict a slight decrease in the infrared-to-radio luminosity ratio with redshift (of the order of 0.3 dex) by $z \sim 2$ due to changing cosmic ray scale heights of galaxies.

Observationally, a lack of sensitive infrared and/or radio data has, until recently, restricted the redshift range of studies of the cosmic evolution of the IRRC. Several observation-based studies have concluded that the IRRC does not appear to vary over at least the past 10–12 Gyr of cosmic history, in that it is linear over luminosity (e.g. Sajina et al. 2008; Murphy et al. 2009). Sargent et al. (2010) found no significant evolution in the IRRC out to $z \sim 1.5$ using VLA imaging of the Cosmic Evolution Survey (COSMOS; Scoville et al. 2007) field at 1.4 GHz with rms ~15 µJy (Schinnerer et al. 2007, 2010). Using a careful survival analysis, Sargent et al. (2010) demonstrate that selecting sources only in the radio or in the infrared for flux-limited surveys can introduce a selection bias that can artificially indicate evolution.

Several other studies (e.g. Garrett 2002; Appleton et al. 2004; Garn et al. 2009; Jarvis et al. 2010; Mao et al. 2011; Smith et al. 2014) have similarly found no significant evidence for evolution of the IRRC out to $z \sim 2$, and out to $z \sim 3.5$ by Ibar et al. (2008).

More recently, studies of the IRRC evolution towards higher redshifts have been facilitated by the revolutionary data products provided by the Herschel Space Observatory (Pilbratt et al. 2010) at far-infrared wavelengths. For example, Magnelli et al. (2015) performed a stacking analysis of Herschel, VLA and Giant Metre-wave Radio Telescope radio continuum data to study the variation of the IRRC over $0 < z < 2.3$. They find a slight, but statistically-significant ($\sim 3\sigma$) evolution of the IRRC. Similarly, Ivison et al. (2010) find some evidence for moderate evolution of the IRRC to $z \sim 2$ using Herschel and VLA data, however their sample selection in the mid-infrared may introduce some bias.

In this paper, we conduct a careful analysis of thousands of galaxies to examine the IRRC out to $z \sim 6$ using deep Herschel observations of the COSMOS field in combination with the VLA-COSMOS 3 GHz Large Project (Smolčić et al. 2017b) – a new, highly sensitive, high-angular resolution radio continuum survey with the VLA. These are the most sensitive data currently available over a cosmologically-significant volume and are thus ideal for such studies. With the wealth of deep, multiwavelength data (from X-ray to radio) available in the COSMOS field, we can conduct a sophisticated separation of galaxy populations into AGN and non-active star-forming galaxies. This allows us to examine the behaviour of the IRRC for each population separately.

In Sect. 2 of this paper we describe our data, the construction of the jointly-selected source sample and the identification of AGN. In Sect. 3 we present our analysis of the IRRC as a function of redshift. In Sect. 4 we discuss our results with respect to the literature and examine the various biases involved. We present our conclusions in Sect. 5. We assume $H_0 = 70 \text{ km s}^{-1} \text{ Mpc}^{-1}$, $\Omega_M = 0.3$ and $\Omega_A = 0.7$ and a Chabrier (2003) initial mass function (IMF), unless otherwise stated. Magnitudes and colours are expressed in the AB system.

## 2. Data

### 2.1. Radio- and infrared-selected samples

It has been shown in Sargent et al. (2010) that studies using solely radio-selected or solely IR-selected samples are biased towards low and high average measurements of the IRRC, respectively, with the difference (in the ratio of infrared to radio luminosities) being roughly 0.3 dex. Therefore, an unbiased study of the IRRC requires a sample jointly selected in the radio and infrared. This section details the construction of the radio-selected and infrared-selected samples and the union of the two, constituting the jointly-selected sample.

#### 2.1.1. Radio-selected sample

The 3 GHz COSMOS Large Project survey was conducted over 384 h with the VLA between November 2012 and May 2014 in A and C configurations. The observations, data reduction and source catalogue are fully described in Smolčić et al. (2017b). The data cover the entire $2\deg^2$ COSMOS field to an average sensitivity of $2.3$ µJy beam$^{-1}$ and an average beamwidth of $0.75''$. In total, 10 830 individual radio components with $S/N \ge 5$ have been identified in the field.

We have searched for optical and/or near-IR (hereafter optical) counterparts to the 8696 radio sources in regions of
the COSMOS field containing good-quality photometric data (i.e. the unmasked regions presented in Laigle et al. 2016). The matching process is identical to that described in detail in Smolčić et al. (2017a) and is briefly summarised here. The best-matching optical counterpart was identified via a position cross-match with the multi-band COSMOS2015 photometry catalogue of Laigle et al. (2016)<fnref n="1" /> using a search radius of 1.2". After rejection of objects with false-match probabilities greater than 20%, the predicted fraction of spurious matches is <1% on average (Smolčić et al. 2017a). We find optical associations for 7729 (89%) of radio sources. These constitute our photometry-matched radio-selected sample.

### 2.1.2. Infrared-selected sample

We use a prior-based catalogue of Herschel-detected objects in the COSMOS field to construct our infrared-selected sample. The Herschel Photodetector Array Camera and Spectrometer (PACS) data at 100 and 160 µm are provided by the PACS Evolutionary Probe (PEP; Lutz et al. 2011) survey. The Herschel Spectral and Photometric Imaging Receiver (SPIRE) data at 250, 350 and 500 µm are available from the Herschel Multi-tier Extra-galactic Survey (HerMES; Oliver et al. 2012). The entire 2 deg<fnref n="2" /> COSMOS field is fully covered by both surveys.

The use of a prior-based, rather than a blind, Herschel source catalogue minimises blending issues. The priors come from the 24 µm Spitzer MIPS (Sanders et al. 2007; Le Floc'h et al. 2009) catalogue of >60 µJy detections, matched to the COSMOS2015 photometric catalogue within a search radius of 1". A source enters our infrared-selected sample if a ≥5σ detection is present in at least one Herschel band at the position of a prior. We have chosen to use a 5σ Herschel detection threshold in order to match the sensitivity level of the radio data. This will be discussed further in Sect. 2.3.4. See Laigle et al. (2016) for a detailed description of the MIPS/COSMOS2015 matching process and the extraction of fluxes from the PEP and HerMES maps. We find 8458 such infrared-detected objects with optical COSMOS2015 counterparts and these constitute our photometry-matched infrared-selected sample.

### 2.1.3. Jointly-selected sample

The jointly-selected sample consists of the union of the radio- and infrared-selected samples and contains 12 333 sources. As can be seen in the upper left panel of Fig. 1, 31% of objects are detected in both radio and infrared, 31% just in radio and 37% just in infrared.

Of the radio-selected sample, 50% are detected in the infrared. However, it is curious that the majority (54%) of objects in the infrared-selected sample, and hence star-forming, are not detected at ≥5σ in the radio. This cannot be explained by a difference in the sensitivities of the Herschel and VLA 3 GHz data, since the two are comparable, as will be shown in Sect. 2.3.4.

However, this can be partially explained by so-called resolution bias (see Smolčić et al. 2017b). Extended or diffuse objects may fall below the detection threshold of the 3 GHz mobile due to the high resolution of the data (0.75"). We have therefore convolved the 3 GHz map to several resolutions between 0.75" and 3.0" and searched for detections in each. This

will be discussed in further detail in Sect. 2.3.2. Of the 4604 objects present in the infrared-selected sample but not present in the radio-selected sample (i.e. undetected in the original, unsmoothed 3 GHz map), 455 are detected at ≥5σ in a mosaic of lower resolution. Hence, 51% of the infrared-selected sample are detected in the radio. The final distribution of objects detected in the infrared, radio or both can be seen in the right-hand panels of Fig. 1.

### 2.1.4. Spectroscopic and photometric redshifts

We require redshifts for all sources in our jointly-selected sample in order to conduct spectral energy distribution (SED) fitting and to compute luminosities. For 35% (4354) of optical counterparts, highly-reliable spectroscopic redshifts are available in the COSMOS spectroscopic redshift master catalogue (Salvato et al., in prep.), with redshifts coming mainly from the zCOSMOS survey (Lilly et al. 2007), DEIMOS runs (Capak et al., in prep.), and the VUDS survey (Le Fèvre et al. 2015; Tasca et al. 2017). Photometric redshifts were available for the remaining sources. For 7607 objects, these are taken from the COSMOS2015 photometric redshift catalogue of Laigle et al. (2016) and were generated using LEPHARE SED fitting (Ilbert et al. 2013). The remaining 372 objects have X-ray counterparts and for these it is more appropriate to use the photometric redshifts produced via LEPHARE SED fitting incorporating AGN templates (Salvato et al. 2009, 2011).

## 2.2. Identification and exclusion of AGN

We wish to consider the relationship between infrared and radio properties due solely to star-formation. Therefore, we identify galaxies likely to host AGN and exclude them from our sample. We exclude a source if it displays evidence of radiatively-efficient AGN emission based on the following criteria:

(i) it displays power-law like emission in the mid-infrared and the IRAC colours satisfy the criteria of Donley et al. (2012) to predict the presence of a dusty AGN torus (as in Smolčić et al. 2017a), and/or

(ii) it has an X-ray counterpart detected in the combined Chandra-COSMOS and COSMOS Legacy surveys (Elvis et al. 2009; Civano et al. 2012, 2016; Marchesi et al. 2016) with a full intrinsic ([0.5–8] keV) X-ray luminosity $L_X > 10^{42} \text{ erg s}^{-1}$ (as in Smolčić et al. 2017a); and/or

(iii) when fitting the SED of the object using both a purely star-forming template and a separate AGN template, the AGN component of the SED is found to be significant based on a Fisher test (Delvecchio et al. 2014). This multi-component SED fitting process is conducted using SED3FIT<fnref n="2" /> (Berta et al. 2013) and is discussed in detail in Delvecchio et al. (2017).

Using these three criteria, we identify 1967 objects from the jointly-selected sample as likely AGN. We refer to these objects as moderate-to-high radiative luminosity AGN (HLAGN). A discussion of this nomenclature can be found in Smolčić et al. (2017a) and Delvecchio et al. (2017), the latter of which also provides a discussion of the relative fraction of AGN identified by each criterion and the extent of overlap.

We further identify an object as an AGN and exclude it from our sample if it does not appear in the IR-selected sample (and

<footnote n="1">¹ An exhaustive list of all available COSMOS multiwavelength data and enhanced data products (such as photometric and redshift catalogues) can be found at http://cosmos.astro.caltech.edu/page/astronomers</footnote>

<footnote n="2">² The multi-component SED fitting code SED3FIT is described in Berta et al. (2013) and is publicly available from http://cosmos.astro.caltech.edu/page/other-tools</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Left: number and fraction of sources present in the radio-selected and/or infrared-selected samples for: all objects in the jointly-selected sample (top), only objects classified as star-forming galaxies (bottom). The grey boxes to the left (right) in each image show the fractions relevant to the infrared- (radio-) selected sample only. Right: same as for the left but including radio detections identified in convolved 3 GHz maps. These samples therefore show which objects are detected (as opposed to selected) in the infrared and/or radio (see Sect. 2.1.3).</figcaption></figure>

<table><caption>**Table 1.** Number of objects in the jointly-selected sample within each galaxy type classification.</caption><thead><tr><th></th><th></th><th>Optical colour selection</th><th>IR-detected only</th><th>Radio-detected only</th><th>IR & radio detected</th><th>Total</th><th>Total</th></tr></thead><tbody><tr><td rowspan="3">Star-forming</td><td>Blue</td><td>M<sub>NUV</sub> - M<sub>r</sub> &lt; 1.2</td><td>2490</td><td>1392</td><td>2331</td><td>6213</td><td rowspan="3">9575</td></tr><tr><td>Green</td><td>1.2 &lt; M<sub>NUV</sub> - M<sub>r</sub> &lt; 3.5</td><td>1228</td><td>783</td><td>1150</td><td>3161</td></tr><tr><td>Red</td><td>M<sub>NUV</sub> - M<sub>r</sub> &gt; 3.5</td><td>100</td><td>-</td><td>101</td><td>201</td></tr><tr><td rowspan="2">AGN</td><td>HLAGN</td><td>N/A</td><td>331</td><td>909</td><td>727</td><td>1967</td><td rowspan="2">2758</td></tr><tr><td>MLAGN</td><td>M<sub>NUV</sub> - M<sub>r</sub> &gt; 3.5</td><td>-</td><td>791</td><td>-</td><td>791</td></tr><tr><td colspan="2">Total</td><td></td><td>4149</td><td>3875</td><td>4309</td><td>12 333</td><td>12 333</td></tr></tbody></table>

**Notes.** The number of objects which are present only in the infrared-detected sample are shown in Col. 4, those present only in the radio-detected sample in Col. 5 and those present in both in Col. 6. The total in each class is also shown in Col. 8. Subsets in italics are considered AGN and are excluded from the star-forming sample.

thus displays no evidence of appreciable star-formation activity), displays red optical rest-frame colours ($M_{NUV} - M_r$) > 3.5 (and is hence considered “passive” in the classification scheme of Ilbert et al. 2009) and is radio-detected (i.e. present in the radio-selected sample). The colour-selection method is described in detail in Smolčić et al. (2017a) and ($M_{NUV} - M_r$) colours are defined in the COSMOS2015 catalogue (Laigle et al. 2016). Considering the lack of observed star formation, the majority of the radio synchrotron emission in such sources is expected to arise from AGN processes. These objects are likely to be low-to-moderate radiative luminosity AGN (MLAGN hereafter), sometimes referred to as low-excitation radio galaxies (LERGs; e.g. Sadler et al. 2002; Best et al. 2005). We note that these objects are referred to as quiescent MLAGN in Smolčić et al. (2017a). We find 791 such objects.

The remaining 9575 sources in the jointly-selected sample display no evidence of AGN presence and we therefore consider

their infrared and radio emission to arise predominantly from star-formation. The distribution of these between the infrared- and radio-selected samples can be seen in the lower panel of Fig. 1. A summary of the classification of all objects in the jointly-selected sample is presented in Table 1. Figure 2 shows the redshift distribution of the star-forming and AGN populations separately. The median redshifts of the star-forming and AGN samples are 1.02 and 1.14, respectively.

All further analysis will focus solely on the star-forming population, unless otherwise stated.

## 2.3. Radio and infrared luminosities

### 2.3.1. Radio spectral indices and 1.4 GHz luminosities

We calculate the spectral index ($\alpha$, where $S_v \propto v^\alpha$) of radio sources by comparing the 3 GHz fluxes to those in the 1.4 GHz VLA COSMOS data (Schinnerer et al. 2004, 2007, 2010). Of
<figure><img src="image_4.png" /><figcaption>Fig. 2. Redshift distribution of the star-forming population (blue solid line) and AGN population (red dashed line) in the jointly-selected sample.</figcaption></figure>

star-forming objects in the radio-selected sample, 1212 (23%) are detected in both the 3 GHz map and the shallower 1.4 GHz map. Figure 3 shows the individual measured spectral indices for these objects. The 5σ lower limit on the spectral index is also shown for all 3 GHz-detected objects without detections at 1.4 GHz. We use a single-censored survival analysis to calculate the median value of $\alpha_{1.4\,GHz}^3$ within several redshift bins. See Sect. 3.1 for details on the binning process. This uses the Kaplan-Meier estimator to incorporate the lower limits when computing the median (Kaplan & Meier, 1958). As seen in Fig. 3, no evolution of the spectral index with redshift is evident. The median in redshift bins at $z < 2.0$ are consistent with $\alpha_{1.4\,GHz}^3 = -0.7$, and is also consistent with that found for all objects in the full 3 GHz source catalogue (Smolčić et al. 2017b). In the two $z > 2$ bins, the median spectral index is more consistent with $\alpha = -0.8$ (see also Fig. 21 in Sect. 4.4). For simplicity, we assume $\alpha = -0.7$ for all objects undetected at 1.4 GHz, however we examine the impact of a particular choice of spectral index on the results in Sect. 4.4. We note that the use of $\alpha = -0.7$ predicts a 1.4 GHz flux that is inconsistent with the 1.4 GHz limit in only 3% of cases.

We convert the observer-frame 3 GHz fluxes ($S_{3\,GHz}$; W Hz$^{-1}$ m$^{-2}$) into 1.4 GHz luminosities ($L_{1.4\,GHz}$; W Hz$^{-1}$) via:

$$L_{1.4\,GHz} = \frac{4\pi D_L^2}{(1+z)^{\alpha+1}} \left(\frac{1.4}{3}\right)^{\alpha} S_{3\,GHz}, \quad (1)$$

where $D_L$ is the luminosity distance to the object in metres.

For any object with no $\ge 5\sigma$ detection in the original 3 GHz mosaic, $L_{1.4\,GHz}$ is calculated by replacing $S_{3\,GHz}$ in Eq. (1) by the flux measured from a lower resolution 3 GHz mosaic, or by the 5σ 3 GHz flux upper limit. The following section will describe how such fluxes and flux limits are determined.

### 2.3.2. 3 GHz detections and flux limits from convolved mosaics

As introduced in Sect. 2.1.3, the high resolution of the 3 GHz data (0.75") means that extended and/or diffuse emission may fall below the detection threshold of the mosaic, corresponding to a peak flux density of five times the local rms. It is therefore possible that some 3 GHz counterparts to objects in the

<figure><img src="image_5.png" /><figcaption>Fig. 3. 3 GHz to 1.4 GHz spectral indices ($\alpha_{1.4\,GHz}^3$) of the star-forming population as a function of redshift. Red points show direct measurements, while yellow triangles show 5σ lower limits for objects not detected at 1.4 GHz. The median within redshift bins are shown by black squares and have been calculated using a single-censored survival analysis, which incorporates the lower limits. The median $\alpha_{1.4\,GHz}^3$ of the star-forming population is consistent with $\alpha_{1.4\,GHz}^3 = -0.7$ (indicated by the horizontal dashed line), at least at $z \le 2$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 4. 3 GHz peak flux measured in each convolved 3 GHz map for three different objects with 1σ error bars. Fluxes associated with a given object are colour-coded and joined by a line. In each case, the indicated point shows where the flux ceases to change significantly with increasing convolution. The mosaic resolution at which this occurs is considered optimal. See text (Sect. 2.3.2) for further details.</figcaption></figure>

infrared-selected sample are missed, particularly at low redshift. To overcome this issue, we convolve the 3 GHz map to several resolutions between 0.75" and 3.0" (i.e. between one and four times the original beamwidth). The convolution increases the average rms of the map, but allows for the potential detection of sources with extended radio emission but missed in the 0.75" mosaic.

If an infrared-detected object is not detected at $\ge 5\sigma$ in the original 0.75" resolution radio mosaic, there are two possibilities:

(i) the object is detected at $\ge 5\sigma$ in one or more convolved radio mosaics; or

(ii) the object remains undetected in all convolved radio mosaics.
<table><caption>**Table 2.** Resolution (i.e. beamwidth) of each convolved 3 GHz mosaic, the average rms and the number of sources for which the 5σ flux measurement is taken from that particular mosaic.</caption><thead><tr><th>Mosaic resolution</th><th>⟨rms⟩ (μJy)</th><th>N (detected)</th></tr></thead><tbody><tr><td>0.75" (original)</td><td>2.34</td><td>7729</td></tr><tr><td>0.9"</td><td>2.5</td><td>199</td></tr><tr><td>1.0"</td><td>2.66</td><td>89</td></tr><tr><td>1.2"</td><td>3.08</td><td>80</td></tr><tr><td>1.4"</td><td>3.57</td><td>29</td></tr><tr><td>1.5"</td><td>3.84</td><td>14</td></tr><tr><td>1.6"</td><td>4.13</td><td>11</td></tr><tr><td>1.8"</td><td>4.77</td><td>11</td></tr><tr><td>2.0"</td><td>5.49</td><td>9</td></tr><tr><td>2.2"</td><td>6.04</td><td>7</td></tr><tr><td>2.4"</td><td>6.73</td><td>3</td></tr><tr><td>2.6"</td><td>7.64</td><td>1</td></tr><tr><td>2.8"</td><td>8.32</td><td>1</td></tr><tr><td>3.0"</td><td>9.13</td><td>1</td></tr></tbody></table>

We calculate the $S_{3\ GHz}$ measurement (or limit) differently for each of these two cases, as follows.

For case (i), we use the integrated flux density from the 3 GHz mosaic with the highest resolution (i.e. smallest beamwidth) where the object is detected at $\ge 5\sigma$. It is appropriate to use the integrated flux density since it is found to be stable with changing resolution, while the peak flux would be underestimated for such extended sources. Table 2 shows the number of sources per mosaic from which the flux measurement is taken. As mentioned in Sect. 2.1.3, 3 GHz counterparts to an additional 455 (5% of) infrared-detected objects are found in lower resolution mosaics.

While we are justified in using the measured 3 GHz flux for these 455 objects with prior positions in the infrared (Herschel and 24 µm), we do not allow the additional objects detected in convolved 3 GHz maps to enter our radio-selected sample. This would result in a highly incomplete sample due to the significantly changing rms with increasing level of convolution and would require additional complex completeness and false detection rate tests (see Smolčić et al. 2017b) which are beyond the scope of this paper.

For case (ii), the 5σ 3 GHz flux limit is taken as five times the value at the corresponding pixel position in the noise map associated with the most appropriate convolved mosaic. The most appropriate mosaic is chosen as follows. For all sources in a given redshift bin, which are detected at $\ge 5\sigma$ in at least one 3 GHz mosaic (i.e. any object in the radio-selected sample or satisfying case (i)), we track how the peak flux (surface brightness) changes with the level of convolution. Several examples are shown in Fig. 4. For each source, the optimal map resolution is that where the peak flux ceases to change significantly with increased convolution. i.e. the first data point which is inconsistent (considering the 1σ errors) with the native point (the measurement from the highest resolution map) but is consistent with all data points in lower resolution maps. This is considered to be the resolution at which all emission from the source is contained within a single map pixel.

For a given object undetected at 3 GHz, the mosaic from which to calculate the 3 GHz flux limit is chosen by sampling from the distribution of optimal resolutions in that redshift bin using a Monte-Carlo-like approach. Examples of the sampled distributions are shown in Fig. 5. The average rms of each convolved map is reported in Table 2. This technique for

<figure><img src="image_7.png" /><figcaption>**Fig. 5.** Normalised distribution of optimal convolved mosaic resolutions for 3 GHz detections in a given redshift bin. For clarity, this is shown only for three redshift bins. See text (Sect. 2.3.2) for explanation of how the optimal resolution is chosen.</figcaption></figure>

determining 3 GHz upper limits overcomes much of the resolution bias in our data.

### 2.3.3. Infrared luminosities

The total infrared luminosities ($L_{TIR}$) of each source in the jointly-selected sample are found by integrating the best-fitting galaxy template to the SED between 8–1000 µm in rest-frame. The data available over this range in the full COSMOS field include *Spitzer* MIPS 24 µm data and the five *Herschel* PACS and SPIRE bands. For 95 star-forming galaxies, sub-millimetre data was also available from various instruments including AzTEC and ALMA (Casey et al. 2013; Scott et al. 2008; Aretxaga et al. 2011; Bertoldi et al. 2007; Smolčić et al. 2012; Miettinen et al. 2015; Aravena et al., in prep.). The SED fitting to the COSMOS multiband photometry was conducted using MAGPHYS (da Cunha et al. 2008) and is presented in Delvecchio et al. (2017).

As discussed in Sect. 2.1.2, we require a $\ge 5\sigma$ detection for a source to enter the infrared-selected sample. This is for the purpose of sensitivity matching with the radio. Of the star-forming galaxies in the infrared-selected sample, 53% of objects are detected at $\ge 5\sigma$ in only one Herschel band, while 1% are detected in all bands. However, catalogued infrared photometry is also available for $3 \le S/N \le 5$ objects. We use this photometry for SED fitting where it is available as it provides better constraints compared with the use of a limit. We have confirmed that this does not result in any bias towards higher luminosities due to noise-induced flux boosting at the faint flux end.

If a source has $S/N < 3$ in a particular Herschel band, we constrain the SED fit using the corresponding 3σ upper limit to the flux. A single value for this limit is used for each band, and full details of this process are provided in Sect. 3 of Delvecchio et al. (2017).

In cases where the source is undetected at $\ge 5\sigma$ in all Herschel bands, integrating the resulting best-fit SED provides only an upper limit on the $L_{TIR}$. This is the case for the 2175 star-forming objects not in the infrared-selected sample. However, we note that the SED template fit, and therefore the $L_{TIR}$ limit, will still be somewhat constrained in the infrared regime since a 24 µm detection is available in 59% of such cases and also due
<figure><img src="image_8.png" /><figcaption>Fig. 6. Total infrared versus 1.4 GHz luminosity for star-forming objects in the jointly-selected sample. Black points show objects directly detected in both the radio and infrared data, red arrows indicate objects in the radio-detected sample with upper limits in the infrared and blue arrows indicate objects in the infrared-detected sample with upper limits in the radio.</figcaption></figure>

to the optical/infrared energy balance performed by MAGPHYS (da Cunha et al. 2008).

Figure 6 shows the $L_{TIR}$ versus $L_{1.4 GHz}$, including limits, for the star-forming sources in the jointly-selected sample. The $L_{TIR}$ and $L_{1.4 GHz}$ versus redshift are shown in Fig. 7.

We have verified that the particular choice of template suite used for SED fitting does not have a significant impact upon the derived infrared luminosities. For a random subsample of 100 objects, we have recomputed the $L_{TIR}$ by using the SED template library of Dale & Helou (2002). We find good agreement with the MAGPHYS-derived $L_{TIR}$, with a median difference of 0.027 dex and a scatter of 0.39 dex. Furthermore, we verified that the $L_{TIR}$ estimates derived from MAGPHYS are consistent with those calculated by using SED templates from Chary & Elbaz (2001), which rely on the 24 µm detection as a proxy for the $L_{TIR}$ at $z < 2$. We found no offset and a $1\sigma$ dispersion of ~0.3 dex between the two $L_{TIR}$ estimates. This agreement has also been determined in previous papers (e.g. Berta et al. 2013; Delvecchio et al. 2017, and references therein).

Along with the total infrared luminosity (which we have defined as 8–1000 µm), we also calculate the far-infrared (FIR) luminosity ($L_{FIR}$) by integrating the star-forming template over 42–122 µm in the rest-frame. The median difference between the total and far-infrared luminosities is 0.30 dex i.e. $\langle \log(L_{TIR}) - \log(L_{FIR}) \rangle = 0.30$. The direct calculation of $L_{FIR}$ allows for ease of comparison with previous studies of the IRRC in the literature which have limited their analyses to the FIR in order to avoid AGN contamination at the shorter wavelengths (e.g. Magnelli et al. 2015; Yun et al. 2001). We are not inhibited by this issue due to our extensive AGN identification process and our ability to decompose the origin of the emission using the multi-component SED3FIT fitting for such objects (see Sect. 2.2).

### 2.3.4. Survey sensitivity comparison

The luminosity limits of the infrared and radio surveys are compared in Fig. 8. The dashed, coloured lines show the 5σ detection limits in each Herschel band. These have been calculated assuming a “typical” $z = 0$ galaxy template found by averaging

the models of Béthermin et al. (2013)<fnref n="3" /> for normal star-forming objects on the galaxy main sequence. The solid black line traces the lowest dashed, coloured line at each redshift. For comparison, the equivalent line assuming $z = 5$ templates is also shown but does not differ significantly to the $z = 0$ case. This represents the lower limit for a galaxy to enter our infrared-selected sample as it must be detected at $\ge 5\sigma$ in at least one Herschel band. However, we stress that this serves only as a rough guide for comparison. In reality, different best-fitting galaxy templates apply to different sources, meaning that it is possible for the $L_{TIR}$ of a particular object to be significantly lower than the predicted limit, while still being present in our infrared-selected sample.

The dashed black line in Fig. 8 shows the sensitivity of the 3 GHz data, assuming a spectral index of $\alpha = -0.7$ and a local conversion factor of $q_{TIR} = 2.64$ (Bell 2003; see Sect. 3.1 below for the definition of $q_{TIR}$). We see that the sensitivities of the 3 GHz and Herschel surveys are well-matched out to high redshift. However, the 24 µm data, which have been used as a prior catalogue for the infrared-selected sample, are more sensitive than both the radio and Herschel data. In fact, 85% of star-forming galaxies in the radio-selected sample are detected in this 24 µm data. Thus, most radio-detected objects are in fact detected to some extent in the infrared, as expected.

## 3. Results

### 3.1. IR-radio correlation redshift trends

The IRRC can be quantified by the parameter $q_{TIR}$, defined as the logarithmic ratio of the total infrared (8–1000 µm) and 1.4 GHz luminosities:

$$q_{TIR} = \log\left(\frac{L_{TIR}}{3.75 \times 10^{12} \text{ Hz}}\right) - \log\left(\frac{L_{1.4 GHz}}{\text{W Hz}^{-1}}\right). \quad (2)$$

We note that the $L_{TIR}$ (in unit W) are divided by the central frequency of $3.75 \times 10^{12}$ Hz such that $q_{TIR}$ becomes dimensionless.

Figure 9 shows the $q_{TIR}$ of all 9575 star-forming galaxies in our jointly-selected sample, as a function of redshift. We have a well-populated sample out to $z \sim 3$, with direct detections in both the infrared and radio data. Upper and lower limits on $q_{TIR}$ are also indicated in the plot. We split the data into ten redshift bins such that they contain equal numbers of galaxies. To incorporate the lower and upper limits when calculating the median $q_{TIR}$ in each bin, we have employed a doubly-censored survival analysis, as presented in Sargent et al. (2010). The basic principle of this method is that the code (written in Perl/PDL by MTS) redistributes the limits, assuming they follow the underlying distribution of the directly-constrained values. This results in a doubly-censored distribution function, as described in Schmitt (1985). An example of the cumulative distribution function and associated 95% confidence interval determined for several redshift bins are shown in Fig. 10.

We use a bootstrap approach to estimate uncertainties on $q_{TIR}$ in each redshift bin by repeating the survival analysis 100 times. In each instance, the values of $L_{TIR}$ are randomly sampled from a Gaussian distribution with a mean equal to the directly-constrained nominal value and a dispersion equal to the measurement error on the nominal value. The $S_{3 GHz}$ measurements are also sampled in the same manner, and the flux limits are again sampled from the distribution of optimal mosaic resolutions (see Sect. 2.3.2 and Fig. 5). These values are then used

<footnote n="3">³ Galaxy templates by Béthermin et al. (2013) at $0 \le z \le 5$ are publicly available at ftp://cdsarc.u-strasbg.fr/pub/cats/J/A%2BA/557/A66/</footnote>
<figure><img src="image_9.png" /><figcaption>**Fig. 7.** Total infrared luminosity (**left**) and the 1.4 GHz luminosity (**right**) versus redshift for all star-forming objects in the jointly-selected sample. Limits are shown as arrows for non-detections.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Total infrared luminosity limit of various data sets. Dashed, coloured lines show the limit in various Herschel bands assuming a $z = 0$ galaxy template (see text, Sect. 2.3.4, for details). The black line traces the lowest coloured line at each redshift and represents the sensitivity limit of the infrared-selected sample. The magenta line is the equivalent using $z = 5$ templates. The $5\sigma$ sensitivity limit of the Spitzer 24 µm data is shown as the black dot-dashed line. The sensitivity limit of the VLA 3 GHz Large Project (dashed black line) is also shown, assuming $q_{\text{TIR}} = 2.64$ (Bell 2003) and a radio spectral index of $\alpha = -0.7$.</figcaption></figure>

for the calculation of the $q_{\text{TIR}}$ measurements or limits and the doubly-censored survival function is regenerated.

The median statistic in a given instance is the value of the 50th percentile of the survival distribution of $q_{\text{TIR}}$ (middle dotted line in Fig. 10). Figure 11 shows an example of the resultant distribution of the 100 median $q_{\text{TIR}}$ measurements in a particular redshift bin and a Gaussian fit to this distribution. The mean of this Gaussian fit provides the final average $q_{\text{TIR}}$ measurement within the redshift bin. The $1\sigma$ dispersion of the Gaussian (~0.01 on average) is combined in quadrature with the statistical error on the median output from the survival analysis (indicated by the shaded regions in Fig. 10; ~0.05 on average) to give the final uncertainty on the average $q_{\text{TIR}}$. These average values and uncertainties are reported in Table 3 and shown in Fig. 9. The 16th and 84th percentiles of the survival function are also quoted, as well as the spread of $q_{\text{TIR}}$ (i.e. $P_{84}-P_{16}$) in each bin. We note that the survival analysis does not constrain some of these parameters

in some redshift bins, due to the number and distribution of the limits in that bin.

We fit a power-law function to the average values of $q_{\text{TIR}}$, weighting by the uncertainty, and find a small but statistically-significant variation of $q_{\text{TIR}}$ with redshift: $q_{\text{TIR}}(z) = (2.88 \pm 0.03)(1+z)^{-0.19 \pm 0.01}$. The errors here are the $1\sigma$ uncertainty from the power-law fit.

# 4. Discussion

## 4.1. Comparison with previous studies

Here we compare our results to several other studies in the literature. To reduce systematics introduced by converting between measurements of $q_{\text{TIR}}$ and the FIR-radio correlation ($q_{\text{FIR}}$), we have compared our results separately to those quoted using TIR and those using FIR. As described in Sect. 2.3.3, we are able to directly measure the $L_{\text{TIR}}$ and $L_{\text{FIR}}$ as a result of the SED fitting process, and therefore can directly calculate both $q_{\text{TIR}}$ and $q_{\text{FIR}}$. For ease of comparison, we have also assumed a spectral index of $\alpha = -0.8$ for non-detections at 1.4 GHz when calculating $L_{1.4}$, as was assumed in Sargent et al. (2010) and Magnelli et al. (2015). Artificial discrepancies could be introduced if different studies assumed different spectral indices, as will be demonstrated in Sect. 4.4.

As shown in Fig. 12, our calculated median values of $q_{\text{TIR}}$ at $z < 1.4$ are consistent with those of Sargent et al. (2010), who also employ a doubly-censored survival analysis to incorporate non-detections into their measurements. At higher redshift, the increase of $q_{\text{TIR}}$ with redshift found by Sargent et al. (2010) is not consistent with our results; a possible reason for this discrepancy is the fact that, as noted by Sargent et al. (2010), high-quality photometric redshifts were not available to them over this range. Sargent et al. (2010) fit a linear relation with redshift to their data up to $z = 1.4$: $q_{\text{TIR}}(z) = (-0.268 \pm 0.115)z + (2.754 \pm 0.074)$. For ease of comparison with our adopted functional form of the fit, we also fit a power-law relation in $(1+z)$ to their data: $q_{\text{TIR}}(z) = (2.78 \pm 0.07)(1+z)^{-0.15 \pm 0.04}$. The slope of this best fit is slightly flatter than, but consistent within $2\sigma$, with our results based on a doubly-censored survival analysis using $\alpha = -0.8$: $q_{\text{TIR}}(z) = (2.85 \pm 0.03)(1+z)^{-0.22 \pm 0.01}$.

The redshift trend that we find is also in agreement with the recent results of Magnelli et al. (2015), as shown in Fig. 13. These authors use a stacking analysis to examine the evolution of
<figure><img src="image_11.png" /><figcaption>**Fig. 9.** IRRC ($q_{\text{TIR}}$) versus redshift for star-forming galaxies. Objects with detections in both the infrared and radio have directly-constrained values of $q_{\text{TIR}}$ and are shown as yellow points. Objects only detected in the radio are upper limits and shown as red triangles. Objects only detected in the infrared are lower limits and shown as blue triangles. A doubly-censored survival analysis has been used to calculate the median $q_{\text{TIR}}$ within redshift bins, indicated by the magenta points. Error bars (smaller than the magenta points here) represent the $1\sigma$ dispersion calculated via the bootstrap method. The magenta line shows the power-law fit to these. The black dotted line and grey shaded area are the local value of Bell (2003) ($q_{\text{TIR}}(z \approx 0) = 2.64 \pm 0.02$) and associated spread (0.26), respectively. In the right-hand panel, the $q_{\text{TIR}}$ distribution is shown separately for direct measurements, upper limits and lower limits.</figcaption></figure>

<table><caption>**Table 3.** Median value of $z$ and $q_{\text{TIR}}$ and number of star-forming galaxies in each redshift bin.</caption><thead><tr><th>z</th><th>Median(z)</th><th>N<sub>total</sub></th><th>N<sub>meas</sub></th><th>N<sub>upper</sub></th><th>N<sub>lower</sub></th><th>q<sub>TIR</sub></th><th>P<sub>16</sub></th><th>P<sub>84</sub></th><th>(P<sub>84</sub> - P<sub>16</sub>)</th></tr></thead><tbody><tr><td>0.005–0.346</td><td>0.23</td><td>958</td><td>425 (44%)</td><td>63 (7%)</td><td>470 (49%)</td><td>2.76<sup>+0.07</sup><sub>-0.06</sub></td><td>2.35</td><td>&gt;3.34</td><td>&gt;0.99</td></tr><tr><td>0.346–0.527</td><td>0.428</td><td>957</td><td>421 (44%)</td><td>117 (12%)</td><td>419 (44%)</td><td>2.66<sup>+0.03</sup><sub>-0.05</sub></td><td>2.28</td><td>3.27</td><td>0.99</td></tr><tr><td>0.527–0.704</td><td>0.626</td><td>957</td><td>402 (42%)</td><td>153 (16%)</td><td>402 (42%)</td><td>2.59<sup>+0.03</sup><sub>-0.06</sub></td><td>2.17</td><td>2.94</td><td>0.77</td></tr><tr><td>0.704–0.877</td><td>0.801</td><td>958</td><td>376 (39%)</td><td>163 (17%)</td><td>419 (44%)</td><td>2.59<sup>+0.03</sup><sub>-0.05</sub></td><td>2.09</td><td>2.92</td><td>0.83</td></tr><tr><td>0.877–1.020</td><td>0.944</td><td>957</td><td>361 (38%)</td><td>172 (18%)</td><td>424 (44%)</td><td>2.56<sup>+0.04</sup><sub>-0.04</sub></td><td>2.07</td><td>3.14</td><td>1.07</td></tr><tr><td>1.020–1.245</td><td>1.149</td><td>958</td><td>312 (33%)</td><td>239 (25%)</td><td>407 (42%)</td><td>2.52<sup>+0.03</sup><sub>-0.06</sub></td><td>1.45</td><td>2.93</td><td>1.48</td></tr><tr><td>1.245–1.509</td><td>1.381</td><td>957</td><td>299 (31%)</td><td>284 (30%)</td><td>374 (39%)</td><td>2.44<sup>+0.04</sup><sub>-0.06</sub></td><td>&lt;-0.54</td><td>2.84</td><td>&gt;3.38</td></tr><tr><td>1.509–1.835</td><td>1.657</td><td>958</td><td>295 (31%)</td><td>322 (34%)</td><td>341 (36%)</td><td>2.35<sup>+0.08</sup><sub>-0.05</sub></td><td>&lt;-1.83</td><td>2.81</td><td>&gt;4.64</td></tr><tr><td>1.835–2.247</td><td>1.995</td><td>957</td><td>355 (37%)</td><td>282 (29%)</td><td>320 (33%)</td><td>2.34<sup>+0.06</sup><sub>-0.05</sub></td><td>0.18</td><td>2.73</td><td>2.55</td></tr><tr><td>2.247–5.653</td><td>2.746</td><td>958</td><td>336 (35%)</td><td>380 (40%)</td><td>242 (25%)</td><td>2.19<sup>+0.10</sup><sub>-0.07</sub></td><td>&lt;-2.03</td><td>2.73</td><td>&gt;4.76</td></tr></tbody></table>

**Notes.** A break-down of the number of sources into those with directly measured $q_{\text{TIR}}$ values ($N_{\text{meas}}$), upper limits on $q_{\text{TIR}}$ ($N_{\text{upper}}$) or lower limits on $q_{\text{TIR}}$ ($N_{\text{lower}}$) is shown, with the fraction of the total shown in brackets. The $q_{\text{TIR}}$ is calculated using a doubly-censored survival analysis to incorporate lower and upper limits. Uncertainties on $q_{\text{TIR}}$ are calculated using a bootstrap approach and incorporate statistical, measurement and systematic errors. A radio spectral index of $\alpha = -0.7$ has been assumed where it is unknown. The 16th and 84th percentiles ($P_{16}$ and $P_{84}$) on the measurement of $q_{\text{TIR}}$ in each redshift bin are given, as determined via the cumulative distribution function output by the survival analysis. ($P_{84} - P_{16}$) is quoted to indicate the spread of the population. We note that limits are given when a value is not constrained by the survival analysis.

the FIR-radio correlation. They find $q_{\text{FIR}}(z) = (2.35 \pm 0.08)(1+z)^{-0.12 \pm 0.04}$. Although our measurements within each redshift bin for star-forming galaxies using a survival analysis are largely consistent with those of Magnelli et al. (2015), the fitted trend we derive has a slightly higher normalisation and steeper slope (although within $2\sigma$): $q_{\text{FIR}}(z) = (2.52 \pm 0.03)(1+z)^{-0.21 \pm 0.01}$. This trend is also in agreement with that found by Ivison et al. (2010) using Herschel data in the GOODS-North field: $q_{\text{FIR}}(z) \propto (1+z)^{-0.26 \pm 0.07}$. We note that Calistro Rivera et al. (2017) find a

similarly decreasing trend of $q_{\text{TIR}}(z)$ for a radio-selected sample of star-forming galaxies in the Boötes field.

As can be seen in Figs. 12 and 13, our measurements of $q_{\text{TIR}}$ and $q_{\text{FIR}}$ in the lowest redshift bin are slightly higher (by more than the $1\sigma$ uncertainty) than the local values of Bell (2003) and Yun et al. (2001), respectively. While we have attempted to account for resolution bias in the radio data, it is possible that we still miss emission from the most extended sources, which are likely to be present at the lowest redshifts. However, our low
<figure><img src="image_12.png" /></figure>

**Fig. 10.** Cumulative distribution functions produced via the doubly-censored survival analysis within the first, fifth and tenth redshift bins. The plots show the fraction of data with $q_{\text{TIR}}$ values less than the value indicated on the lower axis. Shaded regions indicate the 95% confidence interval. The 16th, 50th and 84th percentiles are indicated by the bottom, middle and top dotted lines, respectively.

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

**Fig. 11.** Distribution of the median statistic of the doubly-censored survival function generated by resampling $q_{\text{TIR}}$ 100 times. This particular distribution is for the $0.005 < z < 0.346$ redshift bin. A Gaussian function is fit to the distribution and used to determine the final average value of $q_{\text{TIR}}$ and its uncertainty.

redshift measurements are consistent with those of Sargent et al. (2010) who used radio data at a lower resolution (~1.5") and are therefore less affected by resolution bias. It is therefore unlikely that our results are significantly impacted by resolution bias. It is also possible that our results are affected by issues related to blending in the Herschel maps.

If we exclude the first redshift bin from the fitting procedure, we find that the $q_{\text{TIR}}(z)$ trend is not altered within $1\sigma$, as seen in Fig. 14. To examine the effect of including the local value in the analysis, we include a $q_{\text{TIR}}(z=0) = 2.64$ data point when performing the fit to $q_{\text{TIR}}(z)$. As shown in Fig. 14, the resulting $q_{\text{TIR}}(z)$ trend is slightly flatter: $q_{\text{TIR}}(z) = (2.78 \pm 0.04)(1+z)^{-0.15 \pm 0.02}$. To examine the extreme case, we “anchor” the trend to the local value by fitting the expression $q_{\text{TIR}}(z) = 2.64(1+z)^x$, where $x$ is the free parameter. We still find a decrease in $q_{\text{TIR}}$ with redshift to a $5\sigma$ significance level. This suggests that a decreasing trend of $q_{\text{TIR}}(z)$ is always observed, with the exponent

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

**Fig. 12.** Evolution of $q_{\text{TIR}}$ in comparison with the results of Sargent et al. (2010). The magenta points and fit show the results from this work using a full survival analysis, as in Fig. 9, however a spectral index of $\alpha = -0.8$ has now been assumed for objects not detected at 1.4 GHz. The measurements of Sargent et al. (2010) and their linear fit are shown by the green points and line. A power-law evolution to the individual measurements of Sargent et al. (2010) is shown by the blue line, for ease of comparison. The shaded magenta and blue regions show the $1\sigma$ uncertainty regions calculated by propagating the errors on the corresponding fitting parameters. The local measurement and spread (grey shading) of Bell (2003) are also shown.

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

**Fig. 13.** FIR-radio correlation ($q_{\text{FIR}}$) versus redshift for star-forming galaxies. The evolving fit generated via a survival analysis in this work, assuming $\alpha = -0.8$ for objects not detected at 1.4 GHz, is shown by the magenta line. The evolution found by Magnelli et al. (2015) using a stacking analysis is shown by the green points and curve. The shaded magenta and green regions show the $1\sigma$ uncertainty. The local value of Yun et al. (2001; $2.34 \pm 0.01$) and associated spread (0.26) are shown by the dashed line and grey shaded area, respectively.

of $(1+z)$ between -0.20 and -0.09, regardless of the treatment of the low-redshift measurement.

## 4.2. Impact of upper and lower limits

In Table 3 it can be seen that the fraction of upper and lower limits on $q_{\text{TIR}}$ in a given bin changes with redshift. It is possible that the apparent decrease in $q_{\text{TIR}}$ with increasing redshift could be somehow driven by the changing fraction of limits. To examine
<figure><img src="image_16.png" /></figure>

**Fig. 14.** $q_{\text{TIR}}$ versus redshift for star-forming galaxies and associated power-law fits derived using: all data points (magenta points and solid line; $1\sigma$ uncertainty region shaded), excluding the lowest redshift bin (black dashed line), including the local value of Bell (2003); i.e. fitting to all the magenta points as well as the green point; green dotted line with $1\sigma$ uncertainty region shaded, and anchoring to the local value of Bell (2003) by fitting to the function $2.64(1+z)^x$ where $x$ is the free parameter (blue dot-dashed line). Also shown is the median in each redshift bin calculated using only directly-measured values (i.e. without applying a survival analysis; cyan squares), and the associated fit (solid cyan line).

the extreme case, we ignore all limits and calculate the median of only directly-constrained values of $q_{\text{TIR}}$ in each redshift bin. These values are shown in Fig. 14 with error bars representing the standard error on the median. Using these measurements we find a trend of $q_{\text{TIR}}(z) = (2.59 \pm 0.02)(1+z)^{-0.09 \pm 0.01}$. This fit is flatter than that found when non-detections are correctly accounted for using a survival analysis, producing smaller $q_{\text{TIR}}$ values particularly at lower redshifts. This indicates that accounting for non-detections (limits) in such an analysis has a profound impact on the results.

It is interesting to note that the exponent of the $q_{\text{TIR}}$ trend found when excluding limits agrees with that found in Sect. 4.1 through anchoring to the local value while incorporating limits. It is perhaps worth noting that these studies at $z \sim 0$ also dealt only with direct detections and not with limits. Overall, our conclusion again is that a decrease in $q_{\text{TIR}}$ with redshift is always observed, with the value of the $(1+z)$ exponent varying between -0.20 and -0.09, depending on the particular treatment of non-detections and low-redshift data.

We also note that our survival analysis produces results consistent with those of Magnelli et al. (2015) who accounted for limits using the independent approach of stacking. Mao et al. (2011) also find that the use of a survival analysis and a stacking analysis to account for limits in studies of $q_{\text{TIR}}(z)$ give similar results. Of course, the optimal solution would be to have direct detections available for a complete sample. However, such data are not yet available. Thus, despite our attempts to account for the non-detections through a survival analysis, we acknowledge that our results could still be affected by the sensitivity limitations of the data.

Related to this, we also acknowledge the strong trend between redshift and luminosity of objects in our sample, resulting from the data sensitivity limits. We have performed a partial correlation analysis (see e.g. Macklin 1982) to determine whether a correlation between $q_{\text{TIR}}$ and redshift exists when the

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

**Fig. 15.** Star formation rate predicted from the infrared emission of the Herschel-detected star-forming galaxies in our sample, compared to that predicted via X-ray stacking. The grey region encloses a factor of two around the 1:1 relation, and corresponds to the observed scatter of the $L_X$-SFR relation presented by Symeonidis et al. (2014). No excess is seen in the X-rays, indicating no appreciable contribution from AGN.

dependence on radio or infrared luminosities are removed. However, our results are inconclusive due to biases introduced by the flux limit of our sample. Breaking this degeneracy would require a well-populated, complete sample spanning several orders of magnitude in both radio and infrared luminosity at each redshift. We therefore emphasise that the results we present in this paper are based upon the assumption of a luminosity-independence of $q_{\text{TIR}}$ at all redshifts.

### 4.3. AGN contributions

#### 4.3.1. Are many moderate-to-high radiative luminosity AGN misclassified as star-forming galaxies?

We wish to determine the extent to which AGN contamination could be influencing our results. Although we have used all information at hand to identify objects that are very likely to host AGN, it is still possible that some sources in our star-forming sample have been misclassified or contain low levels of AGN activity. We can investigate the extent to which our sample is contaminated by misclassified HLAGN via X-ray stacking. If misclassified AGN are present, the stacked X-ray flux of the full sample should exceed that expected purely from star formation processes. To test this, we used the publicly-available CSTACK<fnref n="4" /> tool to stack Chandra soft ([0.5–2] keV) and hard band ([2–8] keV) X-ray images of all objects within each redshift bin. The stacked count rate is converted into a stacked X-ray luminosity by assuming a power law spectrum with a slope of 1.4, consistent with the X-ray background (e.g. Gilli et al. 2007). We then apply the conversion between X-ray luminosity and SFR derived by Symeonidis et al. (2014). This conversion was calibrated on Herschel galaxies, both detected and undetected in X-ray, for a better characterisation of the average $L_X$-SFR correlation in inactive star-forming galaxies<fnref n="5" />. Figure 15 shows the SFR derived from X-ray stacking compared to the

<footnote n="4">⁴ CSTACK was created by Takamitsu Miyaji and is available at http://lambic.astrosen.unam.mx/cstack/</footnote>

<footnote n="5">⁵ We note that we have scaled the relation of Symeonidis et al. (2014) to match the X-ray bands and spectral slope chosen here.</footnote>
<figure><img src="image_18.png" /></figure>

**Fig. 16.** Evolution of the IRRC for different source populations. The magenta curve (and points) is the power-law relation found for star-forming galaxies only, while the green curve (and squares) is that found when AGN are included (i.e. star-forming galaxies plus all AGN). The red curve (and triangles) is found when only HLAGN are considered. The cyan curve (and points) is found for the star-forming population of galaxies, excluding those with radio excess. See text (Sect. 4.3.3) for the definition of radio excess. Shading shows the $1\sigma$ uncertainty regions.

SFR derived from infrared luminosities. The latter was found using the conversion of Kennicutt (1998) assuming a Chabrier (2003) IMF and is not expected to be significantly affected by AGN activity and therefore solely attributable to star formation. We find no excess in the X-ray-derived SFR with respect to the IR-derived SFR, indicating that there are very few misclassified HLAGN in our star-forming sample of galaxies.

### 4.3.2. Infrared-radio correlation of AGN

Despite the fact that we expect minimal numbers of misclassified HLAGN, we nonetheless investigate how the emission arising from AGN activity, rather than star-formation processes, could impact the results. Figure 16 shows the resulting $q_{\text{TIR}}$ as a function of redshift if we apply the survival analysis, described in Sect. 3.1, to all objects in the jointly-selected sample. That is, to all star-forming galaxies as well as all HLAGN and MLAGN (see Sect. 2.2). We find only a slight ($<2\sigma$) decrease in the normalisation of the power law fit and steepening of the slope when compared to star-forming galaxies only. This indicates that the inclusion or exclusion of known AGN (which only constitute 22% of the full sample) does not significantly impact the overall $q_{\text{TIR}}(z)$ trend found.

If we consider only objects in the HLAGN category, the inferred trend of $q_{\text{TIR}}$ with redshift for this population appears significantly steeper than that for star-forming galaxies only, although is affected by large uncertainties at higher redshifts. Overall, this suggests that the dependence with redshift of the IRRC of HLAGN is different to that of star-forming galaxies. We note that for this analysis, the $L_{\text{TIR}}$ of HLAGN has been calculated by integrating only the star-forming galaxy component of the multi-component SED template fit determined by SED3FIT. That is, we exclude the AGN component and its contribution to the $L_{\text{TIR}}$. See Sect. 2.2 and Delvecchio et al. (2017) for further details.

We note that, by definition, only upper limits on $q_{\text{TIR}}$ are available for the MLAGN (see Sect. 2.2) and therefore we cannot directly investigate the behaviour of this population alone.

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

**Fig. 17.** Distribution of direct $q_{\text{TIR}}$ measurements (solid green line), lower limits (blue dot-dashed line) upper limits (red dashed line) shown separately for the populations of star-forming galaxies, HLAGN and MLAGN, as indicated.

Figure 17 shows the distributions of direct $q_{\text{TIR}}$ measurements and limits separately for the star-forming galaxies and the two classes of AGN. Although the two classes of AGN comprise only 22% of the full sample, they are responsible for many of the extreme measurements (or limits) of $q_{\text{TIR}}$. In particular, the upper limits of the MLAGN largely sit towards lower $q_{\text{TIR}}$ values (i.e. have radio-excess) with respect to the $q_{\text{TIR}}$ distribution of star-forming galaxies. The lower median $q_{\text{TIR}}$, and large fraction of upper limits, of AGN may be explained by the presence of significant AGN contribution to the radio continuum, with a potentially lower fractional contribution in the infrared. In particular, the far-infrared Herschel bands should be relatively free of AGN contamination, as the thermal emission from the dusty torus peaks in the mid-IR (e.g. Dicken et al. 2009; Hardcastle et al. 2009). Furthermore, we find no obvious bias in the directly-detected $L_{\text{TIR}}$ distribution of the AGN compared to the star-forming population. We again note that any AGN contribution to the $L_{\text{TIR}}$ should have been excluded via the SED-fitting decomposition mentioned above. It is therefore possible that AGN contamination only in the radio regime could be contributing to the observed decrease of $q_{\text{TIR}}$ with redshift.
<figure><img src="image_20.png" /></figure>

**Fig. 18.** Probability distribution in a given redshift bin used to identify objects with radio-excess. The probability distribution function (red line) is generated by taking the derivative of the survival function (a cumulative distribution) in a given redshift bin and is fitted with a Gaussian function (black dashed line).

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

**Fig. 19.** Evolution of the IRRC found when using (i) the real spectral index, where it is known, otherwise using $\alpha = -0.7$ (magenta; $1\sigma$ uncertainty region shaded); (ii) a spectral index of $\alpha = -0.7$ for all sources (black); and (iii) a spectral index of $\alpha = -0.8$ for all sources (cyan). The green points and line show the result of sampling $\alpha$ (where it is unknown) from a Gaussian distribution with $\mu = -0.7$ and $\sigma = 0.35$. The red dashed line shows the use of $\alpha = -0.7$ (at $z < 2$) and $-0.8$ (at $z > 2$) where it is unknown.

### 4.3.3. Radio-excess objects

It is notoriously difficult to separate AGN and star-formation contributions to the radio when no AGN identifiers are available at other wavelengths. Although we have identified MLAGN based upon their red optical colours and lack of Herschel detections (see Sect. 2.2), it is still possible that some objects in our star-forming galaxy sample may also contain MLAGN which contribute only in the radio. Such objects may be expected to show radio excess in their $q_{\text{TIR}}$ values. We therefore again examine the trend of $q_{\text{TIR}}$ versus redshift for the star-forming population of galaxies, this time excluding objects displaying a radio excess. We define an appropriate cut to exclude such objects in each redshift bin as follows: we take the derivative of the survival function and then fit a Gaussian profile to the resulting probability distribution function. An example of this is shown in Fig. 18. The dispersion ($\sigma$) and the mean ($\mu$) of this Gaussian function are used to define radio-excess objects as those with $q_{\text{TIR}} < (\mu - 3\sigma)$. The median value of $\sigma$ across the redshift bins is 0.34. We then rerun the survival analysis excluding these 510 radio excess objects (5% of the star-forming sample). The result, as seen in Fig. 16, is inconsistent with the inclusion of these objects (i.e. the full star-forming sample), having a shallower slope: $q_{\text{TIR}}(z) = (2.83 \pm 0.02)(1+z)^{-0.15 \pm 0.01}$. Thus, sources with appreciable radio excess may play a role in the observed $q_{\text{TIR}}(z)$ trend of the star-forming sample. It is also possible that an appreciable fraction of objects in this star-forming sample are in fact composite systems containing (currently unidentified) MLAGN which contribute to the radio regime, perhaps impacting the observed $q_{\text{TIR}}(z)$ behaviour. Investigating this possibility further will be the subject of an upcoming paper.

### 4.4. Systematics in the computation of radio luminosity

In this section we investigate how the assumptions concerning the exact spectral shape of the emission in the radio regime may affect the derived IRRC.

#### 4.4.1. Influence of the radio spectral index

We firstly examine the impact of the choice of the spectral index ($\alpha$) on the IRRC. As the IRRC is defined via a rest-frame 1.4 GHz luminosity (see 2), which we here infer from the observed-frame 3 GHz flux density (see 1), the choice of spectral indices determines the K corrections<fnref n="6" />. As detailed in Sect. 2.3.1 we have made standard assumptions, i.e. that the radio spectrum is a simple power law ($S_v \propto v^\alpha$). This is supported by the inferred average spectral index of -0.7, approximately constant across redshift (see Fig. 3), and consistent with that typically found for star-forming galaxies, ($\alpha$ = -0.8 to -0.7; e.g. Condon 1992; Kimball & Vezić 2008; Murphy 2009). We have therefore assumed $\alpha = -0.7$ for our 3 GHz sources which are undetected in the shallower 1.4 GHz survey, while for the remainder of the sources we have computed their spectral indices using the flux densities at these two frequencies. From the expression for rest-frame 1.4 GHz luminosity (Eq. (1)) it follows that the change in $q_{\text{TIR}}$ (Eq. (2)), when assuming two different average spectral indices ($\alpha_1$ and $\alpha_2$, respectively), is $\Delta q_{\text{TIR}}(\alpha_1, \alpha_2) = -(\alpha_1 - \alpha_2)[\log(1+z) - \log(\frac{1.4}{3})]$. For $\alpha_1 = -0.7$, and $\alpha_2 = -0.8$, $\Delta q_{\text{TIR}} = -0.1 \log(1+z) + 0.033$. This is illustrated in Fig. 19 where we show $q_{\text{TIR}}$ as a function of redshift derived i) using the measured spectral index where it exists, otherwise setting $\alpha = -0.7$; ii) with an assumed $\alpha = -0.7$ for all sources; and iii) with an assumed $\alpha = -0.8$ for all sources. A change of 0.1 in the assumed spectral index (-0.7 → -0.8) systematically lowers $q_{\text{TIR}}$ and steepens the (1+z) redshift dependence. Thus, the choice of the average spectral index directly affects the normalisation, as well as the derived trend with redshift of the IRRC. As discussed in Sect. 2.3.1, the average spectral index of sources in the two redshift bins at $z > 2$ are consistent with $\alpha = -0.8$, rather than $\alpha = -0.7$. We therefore also show in Fig. 19 the $q_{\text{TIR}}$ trend found when assuming $\alpha = -0.7$ at $z < 2$ and $\alpha = -0.8$ at $z > 2$ for sources where the spectral index is unknown. Although

<footnote n="6">⁶ We note that this is the case for any observing frequency even if, for example, an observed 1.4 GHz flux density is used.</footnote>
<figure><img src="image_22.png" /></figure>

**Fig. 20.** Fractional contribution to 3 GHz flux from free-free emission (top) and synchrotron emission (middle) as a function of redshift, assuming 10, 20, 30, and 40% contributions of free-free emission at 1.4 GHz rest-frame frequency (see legend in bottom panel). The bottom panel shows the power-law evolution of $q_{\text{TIR}}$ determined in Sect. 3.1 (solid line), and the corrected evolution when the free-free emission contribution is properly taken into account.

slightly steeper, this is fully consistent with the use of $\alpha = -0.7$ at all redshifts.

Finally, Fig. 19 also shows the results of sampling the undefined spectral indices from a Gaussian distribution centred at $\mu = -0.7$ and with a dispersion $\sigma = 0.35$. This is the distribution reported in Smolčić et al. (2017b) for all objects detected in the VLA-COSMOS 3 GHz Large Project. This $\alpha$ sampling very slightly steepens the slope and increases the normalisation (albeit within the uncertainties) due to the non-linear dependence of $L_{1.4\,\text{GHz}}$ on $\alpha$ (see Eq. (1)).

### 4.4.2. Influence of free-free contributions

We next test whether the assumption of a simple power-law is a realistic description of the spectral energy distribution in the radio regime. Synchrotron emission is a major component of typical radio SEDs for star-forming galaxies at rest-frequencies of ~1–20 GHz. At higher frequencies, free-free (Bremsstrahlung) emission begins to contribute substantially (see e.g. Fig. 1 in Condon 1992). Both emission processes can be described as power-law radio spectra ($S_{\nu} \propto \nu^{\alpha}$), with a spectral index of –0.8 (synchrotron emission), and –0.1 (free-free emission). For low redshift galaxies, the observing frequencies probe the rest-frame part of the spectrum dominated by the synchrotron emission. However, towards higher redshifts the free-free contributions at rest-frame frequencies become increasingly significant.

In Fig. 20 we show the expected fractional contribution of free-free emission as a function of redshift, assuming various (10–40%) fractional contributions of free-free emission at

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

**Fig. 21.** Radio spectral index ($\alpha_{1.4\,\text{GHz}}^3$) of star-forming sources detected at 3 GHz as a function of redshift. As in Fig. 3, the median values within each redshift bin, derived from a single-censored survival analysis, are shown by the black squares. The predicted evolution in the spectral index due to the contamination of free-free emission based on the M82 model of Condon (1992) is shown by the thick, dashed blue line. Our assumed value of $\alpha = -0.7$ for non-detections at 1.4 GHz is shown by the black dashed line.

1.4 GHz rest-frame. The corresponding synchrotron fractions are also shown as a function of redshift.

The bottom panel of Fig. 20 shows $q_{\text{TIR}}(z)$ if we exclude the free-free contribution and calculate $q_{\text{TIR}}$ using only the synchrotron contribution to the total observed radio emission. The slope of $q_{\text{TIR}}(z)$ is flatter, however a declining trend with redshift is still observed when a 10% contribution of free-free emission at rest-frame 1.4 GHz frequency is assumed (consistent with Condon 1992; Murphy 2009). However, the local $q_{\text{TIR}}$ value is then at the high end of that locally derived by numerous studies (e.g. Bell 2003).

Examining the variation of the spectral index as a function of redshift may also provide information on the extent of the free-free contribution. If we again assume a simplistic radio SED with $\alpha = -0.8$ for synchrotron emission and $\alpha = -0.1$ for free-free emission, then we expect a flattening of the average observed radio spectral index towards higher redshifts. A higher rest-frame frequency is sampled at higher redshifts, given a fixed observing frequency. Since the fractional contribution of free-free emission is larger at higher frequencies, the measured total flux will be larger and hence the spectral index flatter.

Assuming a 10% contribution of free-free emission to the total radio flux density at rest-frame 1.4 GHz, we find that the change of the average spectral index amounts to $\Delta\alpha(z) = \alpha(z = 4.0) - \alpha(z = 0.2) = 0.11$ only. We note that the average spectral index is, under these assumptions, consistent with the local average, $\alpha(z = 0.2) = -0.7$ value inferred using the real data. If we assume free-free emission contributions to the total radio spectrum at rest-frame 1.4 GHz frequency of 20%, 30%, and 40%, we infer an increase (i.e. flattening) of the observed spectral index of only $\Delta\alpha = 0.17$ (albeit with a steeper local spectral index than inferred for the real data). However, the flattening of the average radio spectral index expected under the given assumptions is not supported by our data, as can be seen in Fig. 21.

The general conclusion is that the fractional contribution of free-free emission to the observed radio spectrum with standard, simple assumptions is inconsistent with the derived decreasing trend of $q_{\text{TIR}}$ with increasing redshift. This suggests a
<figure><img src="image_24.png" /><figcaption>**Fig. 22.** IRRC defined at 5.95 GHz versus $L_{\text{FIR}}$ for star-forming objects in our sample at $0.9 < z < 1.1$. The dashed line indicates the median value: $q_{\text{FIR}}^{5.95 \text{ GHz}} = 2.68 \pm 0.02$.</figcaption></figure>

more complex radio SED for star-forming galaxies, compared to the usual assumptions of a superposition of $\alpha = -0.8$ and $\alpha = -0.1$ power-law synchrotron and free-free spectra, respectively, such that at rest-frame 1.4 GHz the free-free contribution amounts to 10% of the total radio emission (e.g. Condon 1992; Yun & Carilli 2002; Bell 2003; Murphy 2009; Galvin et al. 2016).

### 4.4.3. Comparison with local (U)LIRGs

The radio SED for star-forming galaxies was studied by Leroy et al. (2011) who obtained VLA observations of local ($z \sim 0$) (ultra-) luminous infrared galaxies – (U)LIRGs – in C-band (5.95 GHz). They calculated the IRRC in this band and found a median value of $q_{\text{FIR}}^{5.95 \text{ GHz}} = 2.8$ with a scatter of 0.16 dex. At $z \sim 1$, rest-frame 5.95 GHz corresponds to observed-frame 3 GHz. This means that we can use our 3 GHz data to calculate $q_{\text{FIR}}^{5.95 \text{ GHz}}$ with no, or very little, K correction required for objects in our sample at $z \sim 1$. Figure 22 shows $q_{\text{FIR}}^{5.95 \text{ GHz}}$ versus $L_{\text{FIR}}$ for objects in our sample at $0.9 < z < 1.1$ and with $\log(L_{\text{FIR}}) > 11.5$ $L_{\odot}$ for the sake of completeness and a fair comparison (although we note that this restricts us to a luminosity range of $\sim 1$ dex). We find a median $q_{\text{FIR}}^{5.95 \text{ GHz}} = 2.68 \pm 0.02$ with a scatter of 0.24 dex. The $L_{\text{FIR}}$ range of these objects matches closely with the (U)LIRGs sample of Leroy et al. (2011). Therefore, we can directly compare the two samples. We find that minimising K corrections in the radio band flattens the observed trend of decreasing $q_{\text{TIR}}$ with increasing redshift. The inferred $q_{\text{FIR}}^{5.95 \text{ GHz}}$ value at $z=1$ is consistent with a trend $\propto (1+z)^{-0.06}$ (rather than with $q_{\text{FIR}}(z) \propto (1+z)^{-0.21}$ as derived in Sect. 4.1). This suggests that the observed redshift trend of $q_{\text{TIR}}$ may be at least partially attributable to uncertainties in the K corrections applied to the radio flux. Therefore, further investigations into the radio spectra of various star-forming galaxy populations are required for robust determinations of K corrections in the radio regime, having particular relevance for high-redshift star-forming galaxies.

### 4.5. Other physical factors

Along with uncertainty in the radio SED shape and the possible contribution from AGN, it is possible that other physical

mechanisms could be driving a decrease in $q_{\text{TIR}}$ towards higher redshifts. While a thorough investigation of these is beyond the scope of this paper, we nonetheless mention several mechanisms here. One possible driver of an evolving $q_{\text{TIR}}(z)$ is the changing magnetic field properties of galaxies. An increasing magnetic field strength would increase the flux of synchrotron radiation and thereby decrease the measured $q_{\text{TIR}}$. While galaxy-scale magnetic fields are thought to build up over time (e.g. Beck et al. 1996) perhaps from turbulent seed fields (Arshakian et al. 2009), the mean magnetic field strength in a galaxy undergoing a global starburst may be elevated. Tabatabaei et al. (2017) argue that the amplification of magnetic fields within star-forming regions in galaxies with high SFRs could result in a decrease of the infrared-radio correlation. Such a decrease may be stronger at higher redshifts due to the detection bias towards objects with higher SFRs.

It is also thought that major mergers of galaxies can enhance synchrotron emission through various processes and thus result in a decreased measurement of $q_{\text{TIR}}$. For example, Kotarba et al. (2010) performed a magnetohydrodynamical simulation of NGC 4038 and NGC 4039 (the Antennae galaxies) and found evidence for amplification of magnetic fields within merging systems due to compression and shear flows. As discussed above, an increased magnetic field strength would increase the synchrotron emission from pre-existing cosmic rays. Murphy (2013) studied a sample of nearby steep-spectrum infrared-bright starburst galaxies and argue that gas bridges between the interacting taffy-like systems could be the site of enhanced synchrotron radiation which is not related to star formation. Furthermore, Donevski & Prodanović (2015) argue that, in addition to the effects of enhanced magnetic fields, shocks generated by galactic interactions will accelerate electrons and thus further boost synchrotron emission.

The timescale for merger-enhanced infrared emission (due to shock-heating of gas and dust) is expected to be on the order of $\sim 10$ Myr, followed by the enhanced synchrotron emission phase which is expected to last from hundreds of Myrs up to a Gyr (Donevski & Prodanović 2015). If this is the case, then it is statistically more likely that a flux-limited sample contains more galaxies in the phase of synchrotron boosting (Prodanović; priv. comm.). Thus, it is possible that an increasing major merger fraction with redshift, such as that presented by Conselice et al. (2014) to $z \sim 3$, could partially explain a decreasing $q_{\text{TIR}}(z)$.

We note that this is not a comprehensive list of the many physical processes which could be driving an evolving $q_{\text{TIR}}(z)$ and that a number of competing mechanisms, such as inverse Compton energy losses towards higher redshifts (e.g. Murphy 2009), could also be at play.

### 4.6. Radio as a star-formation rate tracer

We have determined that $q_{\text{TIR}}$ decreases with increasing redshift, consistent with previous results in the literature (e.g. Ivison et al. 2010; Sargent et al. 2010; Magnelli et al. 2015). In Sect. 4.4, we have shown that this trend may partly be due to uncertainties in the K correction in the radio due to the overly simplistic assumptions that the radio spectrum can be well-described by a simple power-law. Nevertheless, regardless of the origin of the observed trend, we can make use of it to recalibrate radio luminosity as a SFR tracer as a function of redshift.

In the local Universe, 1.4 GHz rest-frame radio luminosity is anchored to the SFR via the $q_{\text{TIR}}$ parameter (e.g. Condon 1992; Yun et al. 2001; Bell 2003). Following Yun et al. (2001),
we make use of the Kennicutt (1998) calibration for total IR luminosity based SFR:

$$SFR [M_{\odot}/yr] = f_{\text{IMF}} 10^{-10} L_{\text{IR}} [L_{\odot}], \quad (3)$$

where SFR is the star formation rate in units of $M_{\odot}/yr$, $f_{\text{IMF}}$ is a factor accounting for the assumed initial mass function (IMF, $f_{\text{IMF}} = 1$ for a Chabrier IMF, $f_{\text{IMF}} = 1.7$ for a Salpeter IMF), and $L_{\text{IR}}$ is the total IR luminosity in units of Solar luminosities. Relating the SFR to the rest-frame 1.4 GHz luminosity through Eq. (2) and accounting for the redshift and radio spectral index dependences then yields:

$$SFR [M_{\odot}/yr] = f_{\text{IMF}} 10^{-24} 10^{q_{\text{TIR}}(z,\alpha)} L_{1.4 \text{ GHz}} [\text{W/Hz}], \quad (4)$$

where

$$q_{\text{TIR}}(z) = \begin{cases} (2.88 \pm 0.03)(1+z)^{-0.19 \pm 0.01} & \text{for } \langle\alpha\rangle = -0.7 \\ (2.85 \pm 0.03)(1+z)^{-0.22 \pm 0.01} & \text{for } \langle\alpha\rangle = -0.8, \end{cases}$$

where $\langle\alpha\rangle$ is the average assumed spectral index of the star-forming galaxy population, and

$$L_{1.4 \text{ GHz}} = L_{1.4 \text{ GHz}}(z, \alpha) = \frac{4\pi D_L^2}{(1+z)^{\alpha+1}} \left(\frac{1.4}{v_{\text{obs}}}\right)^{\alpha} S_{v_{\text{obs}}}, \quad (5)$$

where $v_{\text{obs}}$ is the observing frequency in units of GHz, here tested and verified for $v_{\text{obs}} = 1.4$ and 3 GHz, and $\alpha = -0.7$, and $-0.8$. It is important to note that the above is valid only for samples of star-forming galaxies selected similarly to those studied here and under the assumptions: (i) of a luminosity-independent IRRC; (ii) of simple K corrections of the radio spectrum ($S_v \propto v^\alpha$) as presented in Eq. (5); and (iii) that the infrared luminosity accurately traces the SFR with redshift.

## 5. Conclusions

We use the new, sensitive VLA-COSMOS 3 GHz Large Project and infrared data from Herschel and Spitzer to push studies of the infrared data correlation (IRRC) out to $z \sim 6$ over the 2 deg² COSMOS field. The excellent sensitivity of the 3 GHz data allows us to directly detect objects down to the micro-Jansky regime. We detect 7729 sources in the 3 GHz data with optical counterparts and redshifts available in the COSMOS database. We identify 8458 sources detected in the Herschel PEP and HERMES surveys with counterparts in Spitzer MIPS 24 µm data and in the optical. Our final sample, jointly-selected in both the radio and infrared, consists of 12 333 unique objects.

We take advantage of the plethora of high-quality multiwavelength data available in the COSMOS field, as well as our ability to perform a multi-component SED fitting process, to separate our sample into (non-active) star-forming galaxies, moderate-to-high radiative luminosity AGN (HLAGN) and low-to-moderate radiative luminosity AGN (MLAGN). We study the IRRC for each of these populations separately.

We examine the behaviour of the IRRC, characterised by the $q_{\text{TIR}}$ parameter, as a function of redshift using a doubly-censored survival analysis to account for non-detections in the radio or infrared along with a bootstrap approach to incorporate measurement errors. A slight, but statistically significant, trend of $q_{\text{TIR}}$ with redshift is found for the population of star-forming galaxies: $q_{\text{TIR}}(z) = (2.88 \pm 0.03)(1+z)^{-0.19 \pm 0.01}$. This is in good agreement with several other results from the literature, although is biased slightly high compared to studies of the local Universe. To examine biases introduced by the sensitivity limits of our

data, we perform various tests incorporating these local measurements, and/or ignoring non-detections. In all cases we find a statistically-significant decrease of $q_{\text{TIR}}$ with increasing redshift, with the slope (i.e. $(1+z)$ exponent) ranging between $-0.20$ and $-0.09$.

When examined separately, we find that AGN have $q_{\text{TIR}}$ measurements biased towards lower values, suggesting that radio wavelengths are more likely than the infrared to be influenced by emission from active processes. It is possible that AGN contributions only to the radio regime could be influencing (i.e. steepening) the observed $q_{\text{TIR}}(z)$ trend, particularly if this occurs in an appreciable fraction of star-forming host galaxies.

We find that the choice of radio spectral index used for the K correction of the 3 GHz flux can influence both the shape and normalisation of the $q_{\text{TIR}}(z)$. The increasing contribution of free-free emission towards higher radio frequencies may also influence the redshift trend, however our results are inconsistent with a typical (M82-based) model of the radio SED. We conclude that a better understanding of the SED of star-forming galaxies is needed for a comprehensive physical interpretation of the apparent redshift evolution of the IRRC. Other physical mechanisms which could potentially drive a decreasing $q_{\text{TIR}}(z)$ include changing galaxy magnetic field strengths and major merger fractions.

Finally, we present a redshift-dependent relation between rest-frame 1.4 GHz luminosity and star formation rate.

**Acknowledgements.** We thank the anonymous referee for useful comments which have helped to improve this paper. We also thank Tijana Prodanović for valuable discussions. This research was funded by the European Unions Seventh Frame-work programme under grant agreement 337 595 (ERC Starting Grant, “CoSMass”). This research was supported by the Munich Institute for Astro- and Particle Physics (MIAPP) of the DFG cluster of excellence “Origin and Structure of the Universe”. N.B. acknowledges the European Unions Seventh Framework programme under grant agreement 333 654 (CIG, AGN feedback). M.B. and P. Ciliegi acknowledge support from the PRIN-INAFF 2014. A.K. acknowledges support by the Collaborative Research Council 956, sub-project A1, funded by the Deutsche Forschungsgemeinschaft (DFG). M.T.S. acknowledges support from a Royal Society Leverhulme Trust Senior Research Fellowship. Support for B.M. was provided by the DFG priority programme 1573 “The physics of the interstellar medium”. M.A. acknowledges partial support from FONDECYT through grant 1140 099.

## References

* Appleton, P. N., Fadda, D. T., Marleau, F. R., et al. 2004, ApJS, 154, 147

* Aretxaga, I., Wilson, G. W., Aguilar, E., et al. 2011, MNRAS, 415, 3831

* Arshakian, T. G., Beck, R., Krause, M., & Sokoloff, D. 2009, A&A, 494, 21

* Beck, R., Brandenburg, A., Moss, D., Shukurov, A., & Sokoloff, D. 1996, ARA&A, 34, 155

* Bell, E. F. 2003, ApJ, 586, 794

* Berta, S., Lutz, D., Santini, P., et al. 2013, A&A, 551, A100

* Bertoldi, F., Carilli, C., Aravena, M., et al. 2007, ApJS, 172, 132

* Best, P. N., Kauffmann, G., Heckman, T. M., et al. 2005, MNRAS, 362, 25

* Béthermin, M., Wang, L., Doré, O., et al. 2013, A&A, 557, A66

* Calistro Rivera, G., Williams, W. L., Hardcastle, M. J., et al. (2017), MNRAS, submitted [arXiv:1704.06268]

* Carilli, C. L., & Yun, M. S. 1999, ApJ, 513, L13

* Casey, C. M., Chen, C.-C., Cowie, L. L., et al. 2013, MNRAS, 436, 1919

* Chabrier, G. 2003, PASP, 115, 763

* Chapman, S. C., Blain, A. W., Smail, I., & Ivison, R. J. 2005, ApJ, 622, 772

* Chary, R., & Elbaz, D. 2001, ApJ, 556, 562

* Civano, F., Elvis, M., Brusa, M., et al. 2012, ApJS, 201, 30

* Civano, F., Marchesi, S., Comastri, A., et al. 2016, ApJ, 819, 62

* Condon, J. J. 1992, ARA&A, 30, 575

* Conselice, C. J., Bluck, A. F. L., Mortlock, A., Palamara, D., & Benson, A. J. 2014, MNRAS, 444, 1125

* da Cunha, E., Charlot, S., & Elbaz, D. 2008, MNRAS, 388, 1595

* Dale, D. A., & Helou, G. 2002, ApJ, 576, 159

* de Jong, T., Klein, U., Wielebinski, R., & Wunderlich, E. 1985, A&A, 147, L6

* Del Moro, A., Alexander, D. M., Mullaney, J. R., et al. 2013, A&A, 549, A59

* Delvecchio, I., Gruppioni, C., Pozzi, F., et al. 2014, MNRAS, 439, 2736
Delvecchio, I., Smolčić, V., Zamorani, G., et al. 2017, A&A, 602, A3 (VLA-COSMOS SI)

Dicken, D., Tadhunter, C., Axon, D., et al. 2009, ApJ, 694, 268

Dickey, J. M., & Salpeter, E. E. 1984, ApJ, 284, 461

Domingue, D. L., Sulentic, J. W., & Durbala, A. 2005, AJ, 129, 2579

Donevski, D., & Prodanović, T. 2015, MNRAS, 453, 638

Donley, J. L., Rieke, G. H., Rigby, J. R., & Pérez-González, P. G. 2005, ApJ, 634, 169

Donley, J. L., Koekemoer, A. M., Brusa, M., et al. 2012, ApJ, 748, 142

Elvis, M., Civano, F., Vignali, C., et al. 2009, ApJS, 184, 158

Galvin, T. J., Seymour, N., Filipović, M. D., et al. 2016, MNRAS

Garn, T., Green, D. A., Riley, J. M., & Alexander, P. 2009, MNRAS, 397, 1101

Garrett, M. A. 2002, A&A, 384, L19

Gilli, R., Comastri, A., & Hasinger, G. 2007, A&A, 463, 79

Hardcastle, M. J., Evans, D. A., & Croston, J. H. 2009, MNRAS, 396, 1929

Helou, G., & Bicay, M. D. 1993, ApJ, 415, 93

Helou, G., Soifer, B. T., & Rowan-Robinson, M. 1985, ApJ, 298, L7

Ibar, E., Cirasuolo, M., Ivison, R., et al. 2008, MNRAS, 386, 953

Ilbert, O., Capak, P., Salvato, M., et al. 2009, ApJ, 690, 1236

Ilbert, O., McCracken, H. J., Le Fèvre, O., et al. 2013, A&A, 556, A55

Ivison, R. J., Magnelli, B., Ibar, E., et al. 2010, A&A, 518, L31

Jarvis, M. J., Smith, D. J. B., Bonfield, D. G., et al. 2010, MNRAS, 409, 92

Kennicutt, Jr., R. C. 1998, ApJ, 498, 541

Kimball, A. E., & Ivezić, Ž. 2008, AJ, 136, 684

Kotarba, H., Karl, S. J., Naab, T., et al. 2010, ApJ, 716, 1438

Lacki, B. C., & Thompson, T. A. 2010, ApJ, 717, 196

Lacki, B. C., Thompson, T. A., & Quataert, E. 2010, ApJ, 717, 1

Laigle, C., McCracken, H. J., Ilbert, O., et al. 2016, ApJS, 224, 24

Le Fèvre, O., Tasca, L. A. M., Cassata, P., et al. 2015, A&A, 576, A79

Le Floc'h, E., Aussel, H., Ilbert, O., et al. 2009, ApJ, 703, 222

Leroy, A. K., Evans, A. S., Momjian, E., et al. 2011, ApJ, 739, L25

Lilly, S. J., Le Fèvre, O., Renzini, A., et al. 2007, ApJS, 172, 70

Lutz, D., Poglitsch, A., Altieri, B., et al. 2011, A&A, 532, A90

Macklin, J. T. 1982, MNRAS, 199, 1119

Magnelli, B., Ivison, R. J., Lutz, D., et al. 2015, A&A, 573, A45

Mao, M. Y., Huynh, M. T., Norris, R. P., et al. 2011, ApJ, 731, 79

Marchesi, S., Civano, F., Elvis, M., et al. 2016, ApJ, 817, 34

Miettinen, O., Smolčić, V., Novak, M., et al. 2015, A&A, 577, A29

Murphy, E. J. 2009, ApJ, 706, 482

Murphy, E. J. 2013, ApJ, 777, 58

Murphy, E. J., Chary, R.-R., Alexander, D. M., et al. 2009, ApJ, 698, 1380

Murphy, E. J., Condon, J. J., Schinnerer, E., et al. 2011, ApJ, 737, 67

Murphy, E. J., Bremseth, J., Mason, B. S., et al. 2012, ApJ, 761, 97

Niklas, S., & Beck, R. 1997, A&A, 320, 54

Norris, R. P., Afonso, J., Appleton, P. N., et al. 2006, AJ, 132, 2409

Oliver, S. J., Bock, J., Altieri, B., et al. 2012, MNRAS, 424, 1614

Park, S. Q., Barmby, P., Fazio, G. G., et al. 2008, ApJ, 678, 744

Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1

Sadler, E. M., Jackson, C. A., Cannon, R. D., et al. 2002, MNRAS, 329, 227

Sajina, A., Yan, L., Lutz, D., et al. 2008, ApJ, 683, 659

Salvato, M., Hasinger, G., Ilbert, O., et al. 2009, ApJ, 690, 1250

Salvato, M., Ilbert, O., Hasinger, G., et al. 2011, ApJ, 742, 61

Sanders, D. B., Salvato, M., Aussel, H., et al. 2007, ApJS, 172, 86

Sargent, M. T., Schinnerer, E., Murphy, E., et al. 2010, ApJS, 186, 341

Schinnerer, E., Carilli, C. L., Scoville, N. Z., et al. 2004, AJ, 128, 1974

Schinnerer, E., Smolčić, V., Carilli, C. L., et al. 2007, ApJS, 172, 46

Schinnerer, E., Sargent, M. T., Bondi, M., et al. 2010, ApJS, 188, 384

Schleicher, D. R. G., & Beck, R. 2013, A&A, 556, A142

Schober, J., Schleicher, D. R. G., & Klessen, R. S. 2016, ApJ, 827, 109

Scott, K. S., Austermann, J. E., Perera, T. A., et al. 2008, MNRAS, 385, 2225

Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1

Smith, D. J. B., Jarvis, M. J., Hardcastle, M. J., et al. 2014, MNRAS, 445, 2232

Smolčić, V., Aravena, M., Navarrete, F., et al. 2012, A&A, 548, A4

Smolčić, V., Delvecchio, I., Zamorani, G., et al. 2017a, A&A, 602, A2 (VLA-COSMOS SI)

Smolčić, V., Novak, M., Bondi, M., et al. 2017b, A&A, 602, A1 (VLA-COSMOS SI)

Symeonidis, M., Georgakakis, A., Page, M. J., et al. 2014, MNRAS, 443, 3728

Tabatabaei, F. S., Schinnerer, E., Krause, M., et al. 2017, ApJ, 836, 185

Tasca, L. A. M., Le Fevre, O., Ribeiro, B., et al. 2017, A&A, 600, A110

van der Kruit, P. C. 1971, A&A, 15, 110

van der Kruit, P. C. 1973, A&A, 29, 263

Voelk, H. J. 1989, A&A, 218, 67

Wrobel, J. M., & Heeschen, D. S. 1988, ApJ, 335, 677

Yun, M. S., & Carilli, C. L. 2002, ApJ, 568, 88

Yun, M. S., Reddy, N. A., & Condon, J. J. 2001, ApJ, 554, 803