  Several works in the past decade have used the ratio between total (rest
8-1000$\mu$m) infrared and radio (rest 1.4~GHz) luminosity in star-forming
galaxies (q$_{IR}$), often referred to as the "infrared-radio correlation"
(IRRC), to calibrate radio emission as a star formation rate (SFR) indicator.
Previous studies constrained the evolution of q$_{IR}$ with redshift, finding a
mild but significant decline, that is yet to be understood. For the first time,
we calibrate q$_{IR}$ as a function of \textit{both} stellar mass (M$_{\star}$)
and redshift, starting from an M$_{\star}$-selected sample of $>$400,000
star-forming galaxies in the COSMOS field, identified via (NUV-r)/(r-J)
colours, at redshifts 0.1$<$z$<$4.5. Within each (M$_{\star}$,z) bin, we stack
the deepest available infrared/sub-mm and radio images. We fit the stacked IR
spectral energy distributions with typical star-forming galaxy and IR-AGN
templates, and carefully remove radio AGN candidates via a recursive approach.
We find that the IRRC evolves primarily with M$_{\star}$, with more massive
galaxies displaying systematically lower q$_{IR}$. A secondary, weaker
dependence on redshift is also observed. The best-fit analytical expression is
the following:
q$_{IR}$(M$_{\star}$,z)=(2.646$\pm$0.024)$\times$(1+z)$^{(-0.023\pm0.008)}$-(0.148$\pm$0.013)$\times$($\log~M_{\star}$/M$_{\odot}$-10).
The lower IR/radio ratios seen in more massive galaxies are well described by
their higher observed SFR surface densities. Our findings highlight that using
radio-synchrotron emission as a proxy for SFR requires novel
M$_{\star}$-dependent recipes, that will enable us to convert detections from
future ultra deep radio surveys into accurate SFR measurements down to low-SFR,
low-M$_{\star}$ galaxies.
