A Uniformly Selected Sample of Low-mass Black Holes in Seyfert 1 Galaxies. III. Radio Sources from the SKA Pathfinders and Beyond

Occupying the intermediate-mass regime of the accretion-jet parameter space, radio continuum emission from active galactic nuclei (AGNs) with black hole mass M BH ≲ 106 M ⊙ (low-mass AGNs) is a valuable probe to the physics of relativistic jets. Yet the number of low-mass AGNs with radio detection is rather limited so far (≈40 in total). In this work, we make two efforts to search for radio counterparts for the largest sample of optically selected low-mass AGNs. First, we collect counterparts from the recent data releases of Square Kilometre Array (SKA) pathfinders such as LOFAR Two-meter Sky Survey (LoTSS). Additionally, we deeply mine in Faint Images of the Radio Sky at Twenty-Centimeters (FIRST), fitting the FIRST images of the optical AGNs with an elaborate procedure optimized to detect faint radio sources. We have obtained 151 radio sources (mainly from the SKA pathfinders), including 102 new reliable sources (signal-to-noise ratio, hereafter S/N, ≥ 5) and 23 new candidates (3.5 ≤ S/N < 5). The majority of these new sources (119 of 125) have flux densities lower than the threshold of the official FIRST catalog. The new sources have rest-frame 20 cm power (P 20 cm) from 1.98 × 1020 to 1.29 × 1023 W Hz−1. For low-z Seyfert galaxies, P 20 cm correlates with M BH intrinsically and positively, yet only marginally with Eddington ratio L/L Edd. In terms of the log N–log S relation for the expanding Universe, the limiting flux density for the completeness of our LoTSS sources turns out to be 0.45 mJy at 1.4 GHz; i.e., complete to such a flux-density level that is 4 times deeper than the official FIRST catalog.

1. INTRODUCTION The physics (e.g., launching, acceleration, collimation and propagation) of relativistic jets in active galactic nucleus (AGNs) is an outstanding hard problem for decades, and it also relates to various initial and boundary conditions of the supplied gas from the host galaxies (see Blandford et al. 2019 for a recent review).This complicated problem is beyond analytic approach, as can be clearly seen from the perspective of mathematics in which it is a large set of intricately nonlinear, time-evolving partial differential equations.It is also beyond numerical simulations and laboratory experiments in the near future, particularly considering the two-way connections between the accretion flow and jet (see Davis &Tchekhovskoy 2020 andBlackman &Lebedev 2022 for recent reviews).Yet progresses have been continually made over the decades; based on the inter-fertilization of physics, numerical simulations and particularly observations, theoretical models and new concepts have been kept improving and emerging (see the above three reviews).
A natural way to make further progresses, from the standpoint of astronomical observations, is to open up and explore new territory (or even new dimensions) in "parameter space".As far as AGN jets are concerned, a less explored space is radio AGNs in the low-mass end: for the time being the detections of radio emission in AGNs with BH mass M BH ≲ 10 6 M ⊙ (low-mass AGNs or called IMBH AGNs 1 ) are rather limited.
1 Following Greene & Ho (2007) and Dong et al. (2012), hereinafter we refer to BHs with M BH < 2 × 10 6 M ⊙ at the centers of galaxies as "lowmass" or "intermediate-mass" BHs (IMBHs); accordingly, for the ease of narration wherever it is not ambiguous, hereinafter we call AGNs hosting low-mass BHs as low-mass AGNs or IMBH AGNs.In the context concerning BHs yet without the term "AGN" occurring, normally we prefer "IMBHs" to "low-mass BHs" because of the possible confusion of the latter with the stellar-mass BHs in low-mass X-ray binaries (LMXBs), let alone the use in the literature of "low-mass BHs" in the so-called "low-mass gap" context concerning the smallest stellar BHs (see, e.g., Yang et al. 2018 and the mini-review therein).
So far, there are only 40 low-mass AGNs that have reliable radio detections in the literature.Previously, Greene & Ho (2007) listed 11 sources with FIRST detection in their optically selected sample of broad-line low-mass AGNs.Liu et al. (2018) listed 26 radio sources (including 4 sources in common with Greene & Ho 2007) by cross-matching the optically selected low-mass AGN samples of Dong et al. (2012) and Liu et al. (2018) with the official FIRST catalog.In addition, Qian et al. (2018) mined 35 more radio counterparts (with 3 < S/N < 5) in FIRST images for the Dong et al. (2012) sample; according to the strict criterion (S/N > 5) of the present paper for the definition of reliable radio sources, we regard Qian et al.'s radio sources as candidates only, and thus do not count them into the 40 detections.Besides the above small radio samples (from the FIRST survey), there are several more IMBH AGNs that are not in the FIRST catalog but have pointed observations (e.g., Gültekin et al. 2022).To sum up, by now the record of distinct IMBH AGNs with reliable radio detections numbers 40 in total,2 which is quite small (cf. the number of radio sources powered by supermassive BHs generally with M BH ≳ 10 7 M ⊙ ; let alone the radio galaxies and quasars having deep and high-resolution data conducted by say VLBI, cf.Blandford et al. 2019) disabling almost all statistically robust conclusions regarding radio properties in the IMBH regime.
It is not easy to expand radio sample of low-mass AGNs.Low-mass AGNs, just like the general population of low-z Seyfert galaxies, are of the so-called "radio quiet" (in terms of radio loudness) in general (Laor 2000;Greene et al. 2006;Greene & Ho 2007;Liu et al. 2018).After all, their radio power is expected to be low owing to their defining intermediate M BH (cf. Figure 7).It is thus difficult to detect their faint radio emission by existing radio facilities within a moderate integration time.
Fortunately, thanks to the international proposal of Square Kilometre Array (SKA; the next generation radio telescope now being constructed), several SKA pathfinders have been built and conducted continuum surveys (see, e.g., Norris et al. 2013).Some of the continuum surveys have released scientific data, which have better sensitivities than FIRST that is the deepest radio continuum survey to date among the finished "all-sky" surveys.Particularly, the pathfinder LOw Frequency ARray (LOFAR), recently releasing its DR2 data covering 5634 deg 2 for the LOFAR Two-metre Sky Survey (LoTSS), can reach very low frequencies from several tens to about two hundreds Mega Hertz (van Haarlem et al. 2013), which profits the detection of faint radio emission of synchrotron radiation nature.Thus, it is the very time to exploit the data of these SKA pathfinders now.
In the present work, as the first step for us to harvest from the aforementioned radio data, we collect radio counterparts as many as possible for the optical low-mass AGNs in two uniformly selected samples (Dong et al. 2012 andLiu et al. 2018;totaling 513 sources).By now the combined sample is the largest of low-mass AGNs (N.B. with M BH in the IMBH regime).Our goal is two-fold: (1) pushing the "parameter space" occupied by radio AGNs to the extremely small end along the M BH dimension, and (2) filling the IMBH regime of the radio-AGN parameter space with much more data points.
By searching in the recent data releases of the continuum surveys of SKA pathfinders, as well as deeply mining in FIRST images (as Qian et al. 2018 did), we find 102 brandnew radio counterparts (S/N > 5) and 23 new candidates (S/N > 3.5), about 3 times the number of previously known radio sources of low-mass AGNs (≈ 40).Besides increasing the number greatly, a second advantage is that these newly found sources are mostly below 1 mJy at 1.4 GHz (the fluxdensity limit of the official FIRST catalog).Our radio sample are complete down to 0.45 mJy at 1.4 GHz, deeper than the limit of completeness of the FIRST catalog (White et al. (1997)) by four times.
This paper is organized as follows.Section 2 describes the parent sample of optical-selected low-mass AGNs, the continuum surveys of SKA pathfinders and their recent data releases, and our deep mining of FIRST images.Section 3 presents the two catalogs of the radio sources from the SKA pathfinders and from our deep mining, respectively.Section 4 presents preliminary studies on radio properties and our investigation on the completeness and depth of our radio sample in terms of the log N − log S relation.Section 5 is the summary.Throughout the paper, we assume a cosmology with H 0 = 70 km s −1 Mpc −1 , Ω m = 0.3, and Ω Λ = 0.7.

DATA AND SAMPLE CONSTRUCTION
We search in the large-scale radio continuum surveys for radio counterparts of the optical low-mass AGNs in the uniformly selected samples of Dong et al. (2012) and Liu et al. (2018) (namely Papers I and II of this series), in order to compile a largest radio sample of low-mass AGNs.Specifically, we conduct our searches in the publicly released data of the continuum surveys of SKA pathfinders ( §2.2), as well as deeply mining the images of VLA FIRST ( §2.3).

Parent sample
Our parent sample, of which we search for radio counterparts, is the combination of the two optically selected samples of Dong et al. (2012) and Liu et al. (2018).Dong et al. (2012) searched for type-1 AGNs having broad Hα line from the SDSS DR4 spectra, and obtained 309 low-mass AGNs with virial BH mass in the range 8 × 10 4 < M BH < 2 × 10 6 M ⊙ .Their virial mass estimator and threshold of low-mass AGNs (2 × 10 6 M ⊙ ) are the same as Greene & Ho (2007).Dong et al. (2012) aimed at so-called uniformly selected broad-line AGNs by elaborately designing a set of automated procedure and quantitative criteria, and attempted to probe accretion rates (namely Eddington ratio) as low as pos-sible (with the resulting Eddington ratios of their low-mass AGNs ranging from ≲ 0.01 to ≈ 1).Using this same procedure, Liu et al. (2018) extended the search to SDSS DR7, and found 204 additional low-mass AGNs.These additional sources have 1 × 10 5 < M BH < 2 × 10 6 M ⊙ , with a similar Eddington-ratio distribution to Dong et al. (2012).To sum up, our parent sample is to date the largest sample of optically selected low-mass AGNs, totaling 513 sources (at redshifts z < 0.35).They are selected from the spectroscopic data set of SDSS DR7 Legacy survey, with a footprint area of 8032 deg 2 (Abazajian et al. 2009).The M BH values range from 8×10 4 to 2×10 6 M ⊙ , and Eddington ratios from ≲ 0.01 to ≈ 1.

Continuum surveys of SKA pathfinders
So far, several radio continuum imaging sky surveys of SKA pathfinders have publicly released their data, such as the LOFAR Two-metre Sky Survey (LoTSS; Shimwell et al. 2017Shimwell et al. , 2019Shimwell et al. , 2022) ) and its deep fields version (hereafter LoTSS-deep;Tasse et al. 2021;Sabater et al. 2021), Apertif imaging survey wide-shallow tier (hereafter Apertif-shallow; Adams et al., in prep.), and the Karl G. Jansky Very Large Array Sky Survey (VLASS; Lacy et al. 2020).
LoTSS is for the entire northern sky conducted by LOFAR high-band antenna working in 120-168 MHz.Two data sets with catalogs have been released, namely the first and second data releases (called DR1 and DR2).The DR1 focuses on a small area of 424 square degrees in HETDEX Spring Field, with a median sensitivity of 71 µJy beam −1 at 144 MHz and an angular resolution of 6 ′′ (Shimwell et al. 2019).The DR2 covers 5634 square degrees for two regions (called RA-1 and RA-13; RA-13 encompassing the DR1 region), with a median sensitivity of 83 µJy beam −1 at 144 MHz and an angular resolution of 6 ′′ (Shimwell et al. 2022).
LoTSS-deep has released its first data set with the catalog recently.It covers three fields (Boötes, Lockman Hole and ELAIS-N1) with a total area about 30 square degrees, in which the rms noise is less than ≈30 µJy beam −1 and the resolution is 6 ′′ , with respect to the working frequency of about 150 MHz (Tasse et al. 2021;Sabater et al. 2021;Kondapally et al. 2021).
Apertif-shallow is conducted by the new Phased Array Feed (PAF) receiver system Apertif installed on Westerbork Synthesis Radio Telescope (WSRT), and plans to survey northern extragalactic sky in the frequency range 1130-1430 MHz.Recently it has released the images of its first data (Adams et al., in prep.) 3 , which covers about 1000 square degrees (part of its planned survey area).In this release, only data in the frequency range 1280-1430 MHz are processed, reaching a median rms depth of 40 µJy beam −1 and a resolution of about 15 ′′ × 15 ′′ / sin(δ).
VLASS is a survey in three epochs, conducted by Karl G. Jansky Very Large Array (JVLA; the upgraded version of VLA).Every epoch will cover the same sky (the sky north 3 http://hdl.handle.net/21.12136/B014022C-978B-40F6-96C6-1A3B1F4A3DB0 of declination −40 • , about 33,885 square degrees), with a frequency range of 2-4 GHz, angular resolution of 2.5 ′′ and sensitivity of about 120 µJy beam −1 .The sensitivity would be reduced to 70 µJy beam −1 when three epochs are combined (Lacy et al. 2020).Recently VLASS has released the Quick Look images of epoch 1, with a catalog published by Canadian Initiative for Radio Astronomy Data Analysis (CIRADA).The Quick Look images, not achieving the full quality of single-epoch images, possess mean sensitivities about 128-145 µJy beam −1 (Gordon et al. 2021).
LoTSS DR1 and DR2, LoTSS-deep DR1, VLASS epoch 1 Quick Look (VLASS1QL) have provided their catalogs, we thus search in them for the radio counterparts of our parent sample by simple positional matching with a match radius of 5 ′′ .In the case of VLASS1QL, we limit our matches to sources with flags Duplicate_flag < 2 and Quality_flag = 0, which indicate non-duplicates and high reliability (Gordon et al. 2021), respectively.It turns out that 116 distinct sources have radio counterparts, some of which have records in more than one catalogs; see Table 2.
Regarding Apertif-shallow, only the compound-beam images have been released without primary-beam correction.We first conduct the primary beam correction, and then fit the images corresponding to the optical AGNs in our parent sample.The correction and the selection of best compound-beam images for every sources follow the procedure described by Adams et al. (in preparation) 3 .We fit the Apertif-shallow images with the same code as to fit FIRST images, which is described in detail in §2.3.1 below.It turns out that there are six radio sources satisfying S/N ⩾ 5, five of which are already found in the above surveys.

Deep mining of VLA FIRST images
In addition to finding radio counterparts in the surveys of SKA pathfinders, we also mine deeply in the VLA FIRST images to extract radio sources fainter than the flux criteria of the official FIRST catalog (see Footnote 21 of Liu et al. 2017, andQian et al. 2018).The official catalog is constructed traditionally: firstly, detect and fit possible radio sources automatedly with software; secondly, pick up reliable sources in terms of a set of criteria.The official criteria are (1) the best-fit peak flux density being not less than 5 times rms; and (2) the peak flux density being not fainter than the threshold, 1 mJy (i.e., the best-fit peak flux density being not fainter than 0.75 mJy before accounting for the CLEAN bias; see §4 of White et al. 1997).
FIRST survey has homogeneous sensitivity of about 130 µJy beam −1 and resolution of 5.4 ′′ at 1.4 GHz over the vast of its coverage of 10,575 deg 2 ; the pixel size of the im-ages is 1.8 ′′ (Becker et al. 1995;White et al. 1997;Helfand et al. 2015).It covers nearly all the sources of our parent sample, while LoTSS-DR2, LoTSS-deep or Apertif-shallow cover half at most.Although VLASS (epoch 1) has a larger coverage area than FIRST, the sensitivity of the current data (namely VLASS1QL) distributes non-homogeneously (Gordon et al. 2021) and gets worse than FIRST in some sky regions.Thus, currently FIRST still has its advantage in terms of both survey area and sensitivity.
We take two measures to make full use of FIRST images: (1) an fairly elaborate procedure is applied during fitting the radio images of the optical AGNs, with greater thoroughness and more steps than Qian et al. (2018), to handle correlated errors and systematical biases or defects (such as residual structures and even sidelobe patterns of faint sources; Thompson et al. 2017;White et al. 2007); (2) an flux criterion (S/N ⩾ 3.5) is carefully selected for plausible radio sources, 4after being tested and assessed quantitatively.The two measures are described in detail below ( §2. 3.1 and §2.3.2).
We would like to point out that, to detect a possible radio source at a specific position in the sky, a criterion of S/N less than 5 is not unreasonable when it has been already known that there is a broad-line AGN detected in the optical (see also White et al. 2007).We can understand this point from two perspectives.Statistically, for a source, assuming random Gaussian noise5 the chance probability of false detection associated with the S/N ⩾ 5 criterion is 1 , which is of course reliable strictly.Here G(x) is the Gaussian probability density function, zero-centered and with a standard deviation σ.As for S/N ⩾ 3.5, that probability is 2E-4, which is actually still reliable statistically.More important is to think in a Bayesian way: the construction of the official FIRST catalog is essentially a blind-source detection procedure, without any prior information (existent fact) for a specific sky position; instead, in our case for every position in investigation there exists an optical broad-line AGN.It is well known that low-z AGNs emit radio (at centimeter wavelengths) more or less (cf.§5.2 of Ho 2008).

The fitting procedure
We develop a fairly elaborate, interactive procedure to fit the 513 FIRST cutout images centering at positions of the optical sources.As described in detail by White et al. (2007), faint sources are generally below the official 5σ threshold (0.75 mJy, without accounting for the CLEAN bias) and therefore susceptible to those systematical errors mentioned above, as well as being blurred by rms noise.
First of all, we visually inspect the radio images with the corresponding optical images overlaid to exclude those with potential contamination of background sources, characteristic sidelobe patterns, and spot-or ridge-like structures (see White et al. 1997White et al. , 2007;;Thompson et al. 2017 for the detailed descriptions of those features).The images of 124 sources are left, passing to the fitting in next step.
Our fitting code, based on Python MPFIT package6 , is used to fit every image interactively.First of all, we note that during the FIRST image processing the CLEAN procedure stops once the CLEANing limit of 0.5 mJy is reached (Becker et al. 1995;White et al. 2007), and therefore the positive and negative sidelobe patterns of faint sources remain in the images.Thus we must include a whole region for each source, rather than only an "island" of bright pixels as the customary strategy used in pipelined softwares to fit bright radio sources (e.g., PyBDSF).Besides, we must visually check the fitting regions to avoid serious contamination caused by the aforementioned systematical errors and defects, and determine a fitting region interactively on the image.The selection procedure of the fitting region is iterated until the histogram of flux densities of the enclosed pixels (excluding the central pixels in the case of bright sources) is close to a Gaussian distribution ("bell"-like shape).
We fit an elliptical Gaussian model to every fitting region.Most of our targets have a single fitting region and thus one Gaussian component, while some other targets have multiple fitting regions and thus multiple Gaussians.From the best-fit Gaussian models, we derive radio parameters: peak flux density F ′ p , integrated flux density F ′ i , right ascension and declination (RA and DEC; J2000), undeconvolved position angle of major axis (measured east of north), undeconvolved major and minor axes (FWHM).The peak and integrated flux densities are corrected for CLEAN bias or a similar bias (called "snapshot bias" by White et al. 2007), according to the empirical bias correction formula (Equation 1 of White et al. 2007; also cf.White et al. 1997).We calculate the S/N of each Gaussian component to be F ′ p /rms (see the FIRST catalog description7 ).For the multi-component sources, we take the S/N of their brightest components as their S/N.
Of the above 124 sources, 86 thus fitted turn out to have S/N ⩾ 3 (and 33 having S/N ⩾ 5), and are kept for our consideration.Their S/N histogram is illustrated in the top panel of Figure 1, where the cumulative number N(S/N > x) versus S/N (the x-axis) is plotted.Because there are not many sources with S/N ⩾ 4, in the figure we use bigger incremental intervals for the range 4 ⩽ S/N < 5 than for the range 3 ⩽ S/N < 4.

S/N criterion to detect deep FIRST sources
As stated in the beginning of §2.3, as the second measure we need to choose an appropriate S/N threshold to deeply mining reliable sources and candidates out of the FIRST images.The idea is this: if we know the numbers of true and false detections versus S/N thresholds, we can pick up an appropriate threshold value that maximizing the number of true detections at a small cost (i.e., with a small fraction of false detections).
To this purpose, we exploit the broad-line AGNs in the footprint of Deep VLA Observations of SDSS Stripe 82 (Hodge et al. 2011;hereafter VLA-S82).VLA-S82 has sensitivity and resolution three times better than FIRST.The broad-line AGNs are from the largest broad-line AGN sample (Liu et al. 2019) covering the Stripe 82 region.We apply our fitting procedure to the FIRST images of these Stripe-82 AGNs.If a thus-fitted FIRST source is recorded in the VLA-S82 catalog, we regard it as true detection; otherwise, a false detection.We can sense that this is a conservative definition of true detections.The detection rate is defined as the ratio of the "true detections" number confirmed by VLA-S82 to the FIRST source number above a specific S/N threshold.By setting the S/N threshold starting from 3 progressively to 5, we obtain the detection-rate curve as a function of S/N threshold, shown in the bottom panel of Figure 1 (red dash-dotted line corresponding to the right y-axis; all in red).
Then, we apply this detection-rate curve of Stripe-82 AGNs to the 86 (probable) radio sources of our parent sample (deeply mined from the FIRST images with S/N ⩾ 3; see §2.3.1), to estimate the numbers of true and false detections at different S/N thresholds.The estimated (i.e., predicted) numbers are shown in the bottom panel of Figure 1, in the same cumulative way as the top panel (squares and circles, respectively, corresponding to the left y-axis; all in blue).For the ease of devising a reasonable threshold, we also plot the curve (blue solid line) that is the number of true detections minus that of false detections.This blue line increases with increasing threshold (as expected), and gets flat around S/N = 3.5.Thus we choose S/N = 3.5 as the threshold to detect deep FIRST sources.
For safety, we make a further check of the 53 FIRST counterparts with 3 ⩽ S/N < 5 of our parent sample, by employing image stacking (cf.White et al. 2007).We divide those sources into five groups: 3 ⩽ S/N < 3.17 (11 sources), 3.17 ⩽ S/N < 3.50 (12), 3.50 ⩽ S/N < 3.80 ( 12 28, 9.93, 10.16, 12.21 and 10.38, respectively (from left to right).The significant central emissions in these stacked images (even of the 3 ⩽ S/N < 3.17 group) support the faint yet genuine existence of, at least a considerable fraction of, the radio sources.

Accuracy of the best-fit parameters
We check the accuracy of our fitting by comparing the bestfit parameters with certain references.To quantify these comparisons, we use two indices: absolute difference p our − p ref , and relative difference (p our − p ref )/p ref , where p stands for any parameters such as F p , F i and so on.
First, we conduct comparison with the official FIRST catalog.This comparison is based on the same data; on the other hand, it is limited to bright sources (only the 26 previously known sources, which is recorded in the official cata- log, can be used).The comparison results are presented in Table 1.The absolute differences of all parameters are small: the mean differences are close to zero, and the standard deviations are far less than their reference parameter values.
Relative differences are also small, with all the mean values less than 2% and standard deviations less than 4%.Second, we conduct comparison by exploiting the small survey VLA-S82 (see §2.3.2), to test the accuracy of our fitting to the FIRST images of faint sources (3 ⩽ S/N < 5).For this purpose, the 10 broad-line AGNs in §2.3.2, which have 3 ⩽ S/N < 5 fits based on the FIRST images and are recorded in the VLA-S82 catalog (noting than VLA-S82 is NOTE-This table lists the mean, standard deviation, maximum and minimum of absolute and relative differences of seven radio parameters, between our best-fit parameter values and those in the official FIRST catalog of the 26 known sources. deeper than FIRST), are also used here.Because the spatial resolution is different between VLA-S82 and FIRST and the sample size is small, we only consider F p (the directly fitted and important quantity).The mean of the F p difference is −0.041 mJy/beam; the standard deviation is 0.104 mJy/beam.Both are smaller than the rms noise of FIRST images (0.14 mJy/beam median).
3. NEW CATALOGS As described in the preceding section, we obtain two radio subsamples of low-mass AGNs.In the subsample from the surveys of SKA pathfinders (SKAPs subsample in short), there are 117 reliable sources (S/N ⩾ 5), among which 100 are brand-new.In the subsample deeply mined from the FIRST images (deep FIRST subsample in short), there are 63 sources (S/N ⩾ 3.5), among which 33 are deemed as reliable sources with S/N ⩾ 5 (including the 26 previously known sources recorded in the official FIRST catalog), and the remaining 30 are deemed as candidate sources with 3.5 ⩽ S/N < 5. We compile two respective catalogs, as presented in Table 2 and Table 3.There are 29 sources in common between the two subsamples.Combining the two subsamples, there are 151 distinct sources in total, including the 26 known sources, 102 reliable sources, and 23 candidate sources.
In addition, considering that LoTSS-DR2 is both deep and wide (now covering about half of the sky area of our parent sample), we present upper limits for the 104 low-mass AGNs that are covered by LoTSS-DR2 yet under detection.The upper-limit flux density is set to be 5 times the local rms of each source, where rms is from the official rms map of LoTSS-DR2.The data of these upper-limit sources are placed in Appendix C (see Table C1).Applications of these sources are deferred to future work; in this study, we only display the upper-limit distributions of their radio power and loudness (in Figures 5 and 6) as reference.We do not consider upper limits based on the FIRST data since it is relatively shallow now in this SKA-pathfinder era.
(3)  4), frequency range of the surveys.Column ( 5), the parent catalog; '1' standing for Dong et al. (2012) and '2' for Liu et al. (2018).Columns ( 6) to ( 9), right ascension and declination (J2000) measured for the radio sources, and the positional errors.Columns (10) to ( 13), peak and integrated flux densities, and their errors.Columns ( 14) to ( 19), undeconvolved position angle (east of north), major and minor axes (FWHM), and their errors.Column ( 20) is the additional note: the sources having detection in the deep FIRST subsample are marked with their identification numbers in that subsample; the sources listed in the official FIRST catalog, marked with 'FCat'; the source having multiple Gaussian components, marked with 'M'.(This table is available in its entirety in a machine-readable form in the online journal.A portion is shown here for guidance regarding its form and content.) Table 3.
Expanded radio low-mass AGNs deeply mined from FIRST images (3)  2) to ( 7) and ( 13) to ( 18), the same definitions as in the above Table 2; see the text for the methods of deriving the parameter values.Column ( 12), the total integrated flux density, which is the sum of the integrated flux densities of the individual Gaussian components for multi-component sources.Column (19), the additional note: the sources having detection in the SKAPs subsample are marked with its identification number in that subsample; the sources in the official FIRST catalog, marked with 'FCat'.(This table is available in its entirety in a machine-readable form in the online journal.A portion is shown here for guidance regarding its form and content.)

Radio subsample from SKA pathfinders
The data of the 117 radio sources in the SKAPs subsample are listed in Table 2.The columns are as follows: Column (1), identification number (ID) assigned for the radio sources in this paper.Different detections of a source have the same ID.Prefix 'S' represents the surveys of SKA pathfinders.
Column (3), data releases of the surveys of SKA pathfinders.
Column ( 4), frequency range of the surveys (in units of MHz).
Columns ( 6) to ( 9), right ascension (RA) and declination (DEC) in decimal degrees (J2000.0)measured for the radio sources, and their measurement errors in arcseconds; the values are from the SKA pathfinders' catalogs (except the Apertif sources; see below).
Columns ( 10) to ( 13), peak flux density F p and its error in mJy beam −1 ; integrated flux density F i and its error in mJy.
Column ( 20) is the additional note.If a source also having a detection in the deep FIRST subsample, we mark its identification number in that subsample; if a source is already listed in the official FIRST catalog, marked with 'FCat'; if a source has multiple Gaussian components, marked with 'M'.
Regarding the parameters in Columns ( 6) to ( 19): for the sources selected from surveys of LoTSS, LoTSS-deep and VLASS1QL, the parameter values come from their official radio catalogs; for the sources selected from Apertif-shallow continuum images, the values are given by MPFIT in the 2D Gaussian fitting (except that the error of peak flux density is the local rms noise of the images, and the error of integrated flux density is calculated through error propagation).

Radio subsample deeply mined from FIRST
The data of the 63 radio sources in the deep FIRST subsample are listed in Table 3.For the sources fitted with multiple Gaussian components, following White et al. (1997) we list every individual components in the catalog.The columns are as follows: Column (1), identification number assigned in this paper, with decimals in digit representing the individual Gaussian components of those multi-component sources.Prefix 'F' represents FIRST survey.
Columns ( 2) to ( 11) and ( 13) to ( 18), the same definitions as in §3.1.Note that the peak and integrated flux densities of the deep FIRST sources have been corrected for CLEAN bias or "snapshot bias" (cf.§2.3.1).
Column (12), the total integrated flux density F i_tot in mJy, which is the sum of the integrated flux densities of the individual Gaussian components for multi-component sources, following Best et al. (2005).
Column ( 19), the note.The sources also having detection in the SKAPs subsample are marked with their identification numbers in that subsample; the sources in the official FIRST catalog, marked with 'FCat'.
The parameter values in columns ( 4) to ( 18) are given by MPFIT in the 2D Gaussian fitting (except the errors to the peak and integrated flux densities; see §3.1).
In the deep FIRST subsample, the 29 sources in common with the SKAPs subsample are as follows: 17 previously known sources (in the official FIRST catalog), 5 newly found reliable sources (S/N ⩾ 5), and 7 candidate sources (3.5 ⩽ S/N < 5).Four of the 12 new sources (5 + 7) have data of pointed observations in the literature, which are marked in Table A1.We noticed that among the 37 new FIRST sources (i.e., excluding the 26 ones recorded in the official FIRST catalog), 16 are covered by the footprints of LoTSS-DR2, giving a confirmation fraction by LOFAR of 75% (namely 12/16).
Except the above 63 radio candidates out of the 86 (probable) radio sources deeply mined from the FIRST images for our parent low-mass AGN sample ( §2.3.1),there are 23 remaining sources with 3 ⩽ S/N < 3.5.They are optimal targets of radio low-mass AGNs for deep and high-resolution observations in the future, and thus we list their information (including the best-fit parameters based on their FIRST images) in Appendix B.

SAMPLE PROPERTIES
As stated above, this work has obtained 151 radio lowmass AGNs, of which 128 are reliable detections (with S/N ⩾ 5).Either of the two numbers is more than 3 times the number of known radio low-mass AGNs (totaling 40).Of the 151 sources, 125 are newly found (from the surveys of SKA pathfinders and our deep mining of FIRST images).Among our newly found sources, 119 sources are fainter than the flux threshold of the official FIRST catalog (1 mJy at 1.4 GHz; regarding sources from SKA pathfinders, for comparison their flux densities are converted to those at 1.4 GHz assuming a radio spectral index of 0.46).Thus our radio sample extends deeper into the faint-flux regime than those previous small radio samples of low-mass AGNs (see §1).This encourages us to investigate the statistical properties of the radio sources, and assess the completeness and depth of the present radio sample.
We note that, due to the insufficient spatial resolution of the radio data to identify radio emissions from the Seyfert nuclei, strictly speaking, these radio sources are of plausible AGN origin; thus the results of our investigations on jet emission are instructive but not conclusive.Considering the hybrid parentage and not decent spatial resolutions of these radio sources, here we take a conservative way to identify "extended" radio sources: we only classify these sources to be either single-component or multicomponent, in terms of the number of best-fit Gaussian components.It turns out that 132 of the 151 sources have only one Gaussian component, and the rest 19 have two or more Gaussian components.Ten of the multi-component sources are in the SKAPs subsample (we mark them with 'M' in Table 2), and nine in the deep FIRST subsample (denoted by decimals in their identification numbers in Table 3).The (probably) most extended sources, in terms of the maximum separation among the Gaussian components, have 10.9 kpc in the SKAPs subsample, and 20.9 kpc in the deep FIRST subsample.
Brief notes are in order for two special multi-component sources.J122342.81+581446.1 in the SKAPs subsample has two Gaussian components, but the optical image indicates one of the two radio components sits on a tiny background galaxy (J122343.72+581454.3 at z = 0.588 as identified from its SDSS image and spectrum)8 ; it is not clear at this point if there is any connection between that radio component and the background galaxy.Thus we mark this source with 'M?' in Table 2.The second is J134738.24+474301.9 in the deep FIRST subsample, for which we need two Gaussian components to fit.Interestingly, it has been actually recorded also in the LoTSS-DR2 catalog yet without a multicomponent tag.We inspect its LoTSS image, and see multiple radio sources within a radius 28 ′′ (namely, a projected physical radius of 33.6 kpc) without any other optical counterparts.Thus the spatially distributed radio emission is of the same source plausibly, with the brightest component J134738.24+474301.9 located in the center; all these components appear to constitute a jet, roughly in a linear configuration.The LoTSS-DR2 and FIRST images, with its optical image from the SDSS (shown as contours) overlaid, are shown in Figure 3.The low-mass AGNs with extended jets are very few so far, and valuable for jet physics (see, e.g., SDSS J095418.15+471725.1 by Gültekin et al. 2022).Thus J134738.24+474301.9 deserves deep and highresolution pointed radio observations in the future.

Spectral index
There are 25 sources in our sample that have detections in both LoTSS-DR2 and deep FIRST subsamples.There is a large frequency difference between LoTSS and FIRST, centering at 144 MHz and 1.4 GHz, respectively, which profits the investigation of radio spectral slope.Moreover, the angular resolutions of the two surveys are similar (6 ′′ and 5.4 ′′ , respectively), minimizing the potential bias in flux due to resolution mismatch.Assuming there is no considerable variability between different epochs (at least statistically), we investigate the two-point radio spectral indices of the 25 sources.The indices are calculated by using their integrated flux densities and adopting the spectral form f ν ∝ ν −αr , where α r is the radio spectral index.
As shown in Figure 4, the derived α r values span from 0.03 to 0.98, with a standard deviation of 0.26, a mean of 0.47 and a median of 0.45, consistent with Ho & Ulvestad (2001).Adopting the generally used division α r = 0.5 between flat and steep radio spectra, 14 of the 25 sources are flat (α r < 0.5).In the context of AGNs, flat radio spectra are from compact radio cores, and steep spectra from extended jets or lobes (see §1.3.1 of Peterson 1997).More importantly, while it is also possible for steep radio spectra to originate from star formation activity, it is generally believed that flat-spectrum radio sources cannot be associated with star formation.Interestingly, 5 of the 11 steep-spectrum sources are multi-component in morphology as described in §4.1.1,while only 3 of the 14 flat-spectrum sources are fitted with multiple Gaussians.With α r being 0.45 and 0.37 two of them are near the boundary to be steep-spectrum, and the third one is J134738.24+474301.9 with an extended jet as shown in Figure 3, of which the low-frequency flux (LoTSS) might be underestimated owing to the extended structure (see §4.1.1)and thus its true spectral index could be steeper.

Radio power
We make quick statistics to the powers at rest-frame 20 cm P 20cm of the radio sources.To calculate the powers, we use the integrated flux densities specified as follows.For the sources having multiple detections in the SKAPs subsample, their flux densities from LoTSS-DR2 are used; for the 29 sources having detections in the both subsamples, the flux densities from the SKAPs subsample are used for the 12 newly-found ones (because their SKAPs fluxes are above 5σ), and their FIRST flux densities are used for the 17 previously known (measured at the rest-frame 20 cm wavelength and with S/N ⩾ 5); the total integrated flux densities are used for the multi-component sources.For flux densities at frequencies different from 1.4 GHz, we convert them to 1.4 GHz (namely at the 20 cm wavelength) by assuming a radio spectral index α r = 0.46 (Ho & Ulvestad 2001;Ulvestad & Ho 2001).The log P 20cm values of our radio sources are recorded in Table 4.The P 20cm values of the 125 newly found sources range from 1.98 × 10 20 to 1.29 × 10 23 W Hz −1 with a median of 1.12 × 10 22 W Hz −1 , while the P 20cm of the 26 known sources (i.e., listed in the official FIRST catalog) from 2.40 × 10 21 to 3.72 × 10 23 W Hz −1 with a median of 3.34 × 10 22 W Hz −1 .
For comparison, we also calculate the rest-frame P 20cm of the sample of broad-line AGNs in SDSS DR7 (Liu et al. 2019; dubbed BLAGNs hereinafter, excluding the 26 known radio low-mass AGNs of Liu et al. 2018 since they have been included in our sample).For those BLAGNs we use their integrated flux densities from the FIRST catalog.
The log P 20cm histograms of our radio sample and the socalled BLAGNs are plotted in Figure 5.The histograms are normalized to have unity probability (namely the area under every histogram).The distribution of our 151 low-mass AGNs is in the lower end of the BLAGN distribution, which has a larger median of 1.58 × 10 23 W Hz −1 .

Radio loudness
Radio loudness R, an indicator of the relative intensity of jet to accretion-disk emission, is defined as the ratio in the rest frame, R ≡ f ν (6 cm)/ f ν (4400Å).We calculate the rest-frame f ν (6 cm), by using the same integrated flux densities, radio spectral index α r and rest-frame correction as used above to calculate the rest-frame P 20cm .The restframe f ν (4400Å) is converted from monochromatic luminosity L λ (5100Å), by assuming an optical spectral index α o = 0.44 ( f ν = ν −αo ; Vanden Berk et al. 2001).The L λ (5100Å) data are derived from the luminosity of broad Hα emission line (L Hα B ).The calculated log R data are recorded in Table 4.The R values of the 125 newly found sources range from 0.9 to 181.3 with a median of 12.0, while the R of the 26 known sources from 6.3 to 314.1 with a median of 29.6.
For comparison, similar to the P 20cm case, we also calculate the radio loudness of the BLAGNs.The normalized histograms of log R of our sample and the BLAGNs are plotted in Figure 6.Both distributions are similar.But here a caveat is in order: the distributions of both samples are only for the detected radio sources from the surveys such as FIRST (and LoTSS for our sample); if the upper-limit sources are taken into account, the distributions certainly would be moved and skewed to small-R end (cf.Figure 6).
If the conventional radio-loud criterion R > 10 is taken, 23 of the 26 known sources (Liu et al. 2018) and 70 of our 125 brand-new sources are radio-loud.Such a high radio-loud fraction may be indicative of considerable radio contributions from their host galaxies; this is particularly possible considering the low optical luminosity of these low-mass AGNs.

Dependence or not of radio power and loudness on M BH and L/L Edd in Seyferts
A linear correlation in logarithmic space has been found among BH mass M BH , and radio and X-ray powers, which forms a so-called fundamental plane (FP) of BH activity (Merloni et al. 2003;Falcke et al. 2004;Körding et al. 2006;Gültekin et al. 2009Gültekin et al. , 2014)).Similarly, the FP can be translated into the plane of radio loudness with M BH and the bolometric Eddington ratio L/L Edd (see, e.g., Xie & Yuan 2017;Qian et al. 2018).Safely and strictly speaking, the BH FP is restricted to AGNs with radiatively inefficient accretion flows The normalized histograms of the so-called BLAGNs (black dashed line) and the upper-limit sources (blue dotted line) are also plotted in the same way as in Figure 5.
In light of the data quality of the present sample (e.g., concerning the accurate radio emission from the jets), to avoid any possible confusion or even misleading here we do not attempt to model/fit the FP relationship in Seyferts or IMBH AGNs, but conduct correlation and partial correlation tests instead, which are less demanding on data quality than fitting a FP relation.
For this purpose, we use all the reliable radio data for both our low-mass AGN sample and the BLAGNs described in §4.1.3and §4.1.4.Note that the BLAGNs denote the other broad-line AGNs (namely, more massive Seyfert 1s) in the sample of Liu et al. (2019), and that their radio data are retrieved from the official FIRST catalog.The data of M BH and L/L Edd are obtained from the respective catalogs of (Dong et al. 2012;Liu et al. 2018Liu et al. , 2019)).
Viewing from the left panel of Figure 7, there are apparent correlations between any two quantities of log P 20cm , log L/L Edd and log M BH .The one between log L/L Edd and log M BH , is well known in any optically selected samples and caused by selection effects and cosmic downsizing (see Dong et al. 2012 for the detail).The other two concern P 20cm correlating positively with M BH and L/L Edd , respectively.
In order to minimize the effects caused by inter-dependence among the three quantities and to pin down the intrinsically significant correlations, we perform various Spearman partial correlation tests to the above two correlations concerning radio power.First, we conduct the following tests: for the correlation between P 20cm and L/L Edd controlling for M BH , r s (Spearman correlation coefficient) is 0.48 and P null (probability for null hypothesis that there is no correlation) is 1.2 × 10 −104 ; for that between P 20cm and M BH controlling for L/L Edd , r s = 0.67 and P null = 2.4 × 10 −245 .Besides, these two correlations may be affected by the dependence of the three variables on redshift in common.Thus we further control the redshift in the above two partial correlation tests.As a result, P null becomes 1.1 × 10 −3 for that between P 20cm and L/L Edd (controlling for both M BH and redshift), whereas r S = 0.26 and P null = 9.8 × 10 −29 for that between P 20cm and M BH (controlling for both L/L Edd and redshift).Thus we conclude that for Seyfert galaxies radio power correlates intrinsically and positively with M BH , yet marginally at most with L/L Edd .
We also have test correlations among log R, log L/L Edd and log M BH .There are no strong correlations of log R with either M BH or L/L Edd .This can be seen easily from the right panel of Figure 7, in which log R nearly distributes uniformly in the plane of log L/L Edd and log M BH .This may be due to large scatter in the log R data compared with its small dynamic range (log R mainly in the range between 0 and 2.5).Thus we cannot say anything about the correlations concerning radio loudness.
Note that in Figure 7 the gap between low-mass AGNs and other Broad-line AGNs are due to the difference in the viral mass formulae used in the three catalog papers.But the gap has little influences on the above Spearman tests because Spearman method is rank ordered.We are still working on an expanded sample of broad-line AGNs (W.-J.Liu et al. 2023, in preparation), in the present work we thus do not investigate further the radio properties of the whole Seyfert-1 population.

The log N − log S test and predictions
Motivated by the log N − log S exploration in the FIRSTcatalog paper (White et al. 1997), with great interest we assess the completeness and depth of our radio sample by carrying out the similar exploration but down to a fainter limit (S).In addition, with a certain intellectual curiosity we pre-   2 and 3. Columns ( 2) to ( 5), the redshift, virial mass of black hole, Eddington ratio and luminosity of broad Hα emission line, retrieved from the two parent optical catalogs.Column (6), the continuum luminosity at rest-frame 5100 Å (L5100 ≡ λL λ at λ = 5100 Å), derived from L Hα B .Column ( 7), the calculated radio power at rest-frame 20 cm.Column (8), radio loudness.(This table is available in its entirety in a machine-readable form in the online journal.A portion is shown here for guidance regarding its form and content.) dict the detectability of the entire 513 optical selected lowmass AGNs in the parent sample in terms of the log N − log S relation.
In a non-expanding spatially flat universe, the so-called Euclidean universe, this relation is N ∝ S −1.5 if the observed population has a constant space density, where N is the number of all the sources having flux above the limiting flux S (see §10.1 of Peterson 1997 for the detail).Our Universe is spatially flat yet expanding, and thus the log N − log S relation is not so simple (actually without an analytic expression).But once the cosmological parameters are given, the relation is certain (see Appendix D).Thus, if we know a priori that the space density of the sample of interest is constant, we can use the log N − log S relation to judge to what a limiting flux the sample is complete.
In this respect, we regard the radio emission of the whole galaxies (Seyfert nuclei plus their host galaxies) to be from the same population.The practical reason is that in light of the spatial resolution we cannot separate the nuclei emission from the host-galaxy one.The second and positive reason is that the radio properties of general Seyfert galaxies (unlike quasars that can either have large-scale jets or not; see §4; Ho 2008), as well as their host-galaxy properties (see, e.g., Dong et al. 2012), appear indeed belong to a same population.

Distribution on the log N − log S plane
For our radio detected low-mass AGNs, we cannot treat them as a single sample to make the log N − log S test, because they come from different radio surveys with different sensitivities and different survey areas.Below we use two groups of our radio sources, one from the FIRST survey (including our deeply mined and the officially cataloged; hereafter called FIRST sources for short) and the other from LoTSS-DR2, because both surveys have large areas and their sources dominate our radio sample.The integrated flux densities (physically, nucleus + host emission) are used.
We plot the log N − log S distributions of both FIRST and LoTSS-DR2 sources in the same diagram (right panel of Figure 8), for the ease of checking the completeness of the two groups and the overall behavior of our radio sample.For putting them in the same figure, the two groups need to have flux densities at the same frequency, and have source numbers corresponding to an equal survey area.First, we convert flux densities of FIRST sources at 1.4 GHz to that at 0.15 GHz (the central frequency of LoTSS) adopting the spectral index of 0.46 used above.Second, we scale the number of LoTSS-DR2 sources to that corresponding to the survey area of FIRST, by multiplying the ratio of their survey areas (7813/4111).Here, 7813 (deg 2 ) is the overlapped survey area of FIRST and the Spectroscopic Legacy of SDSS-DR7 (where our parent sample come); 4111 (deg 2 ) is the overlapped survey area of LoTSS-DR2 and SDSS-DR7 Spectroscopic Legacy.
We also plot the log N − log S distributions of LoTSS-DR1 along with FIRST sources (left panel of Figure 8); this is used to verify the accuracy of predicting source number with the log N − log S relation (see §4.3.4).

Fitting the log N − log S relation
On the theoretical side, however, the situation is a little complicated because our real Universe is a spatially flat yet expanding universe, not an Euclidean one.In this Universe, the theoretical log N − log S relation deviates from the simple power law (namely the canonical N ∝ S −1.5 relation of the Euclidean universe).Here we give a brief description, and refer the reader to Appendix D for the detail.First of all, because our parent sample (the optically selected low-mass AGNs) is limited to low redshifts z < 0.35, we can simply assume that our radio sources have a constant comoving density n 0 .Then the N − S relation is expressed as follows, where ϕ(L)dL is the luminosity function; z max and d co are described in Appendix D. Note that z max , or equivalently d co (z max ), is a function of S/L (refer to Equation D3), not of L explicitly; i.e., the functional forms of z max and d co (z max ) can be constructed in a way that they are not sensitive to the concrete values of luminosity (see Appendix D).Moreover, the distribution of the radio luminosities of our sample is somehow narrow and well peaked (see Figure 5).Taking the above two points into account, in order to derive a theoretical, expansion-modified log N − log S relation we thus simplify the luminosity distribution of our sample with an effective luminosity (e.g., around the median luminosity).Then Equation 1 is fully represented by the system of Equations D3 and D5.As discussed in Appendix D, this equation system can be solved numerically only, but we can model the numerical N(S) relation with a carefully constructed function of the modified double power-law form (Equation D7).This model N(S) curve is devised thus that the curve shape is controlled by three parameters (β, γ and d t ) that are independent of the limiting flux S and the luminosity L (when S and L are normalized; see the text around Equation D7 for the detail).
We would like to note that the steep radio luminosity function of AGNs in the related luminosity range (Heckman & Best 2014;Nagar et al. 2005;Mauch & Sadler 2007;Yuan et al. 2017;Saikia et al. 2018) further supports the above simplification of a single effective luminosity.The intrinsic radio luminosity function of low-redshift AGNs in our luminosity range has a fairly steep slope, α LF ≈ 0.7-0.78,if a single power-law (n ∝ L −αLF ) is fitted to the data (Nagar et al. 2005;Saikia et al. 2018).By Estimating in terms of Equation D6, we can sense that with decreasing L the power-law increase of the n 0 term and the decrease of L term in that equation counterbalance each other.We notice also that those studies reported an intrinsic turnover begins at around P 20cm ≈ 10 19.5 W Hz −1 toward low-luminosity end of the radio luminosity function of low-redshift AGNs.But this turnover luminosity is six times smaller than the smallest P 20cm of our sample (see Figure 5).As the above studies stated, the radio luminosity function in the neighboring (higher) luminosity range continues with the same slope down to at least P 20cm ≈ 10 20 W Hz −1 .Thus, our present sample is immune to this lowluminosity turnover in the radio luminosity function.
We now fit the model N(S) function to the distributions of our source counts as plotted in Figure 8.For the ease of narration, we put the model function here (see also Equation D7), Here the shape parameters β, γ and d t are fixed to the theoretical values obtained in Appendix D; β = −8.339× 10 −3 , γ = 0.392, and d t = 12.36.The normalization factors S 0 and L 0 are fixed to 1 mJy and 10 22.57 W•Hz −1 at 0.15 GHz (the median luminosity of our data plotted in Figure 8), the same as adopted in Appendix D (cf. Figure D1).The effective luminosity, L ′ , is fixed to the median luminosity (being equal to L 0 ); this brings little error considering that the distribution of the radio luminosities of our sample is narrow and fairly peaked, and that model N(S) is not sensitive to the value of L ′ (note that L ′ just shifts the N(S) curve horizontally).C, mainly related to the unknown quantity n 0 and responsible for shifting the model N(S) vertically, is a free parameter and obtained by fitting the model N(S) to the log N − log S distributions of our sources.
A technical detail should be described, concerning the fitting range of the data.Obviously the low-flux ends of both FIRST and LoTSS data are not complete, and thus should not be included in the fitting.We determine the fitting range in this way: first we get the highest redshift of an observed data set (either FIRST or LoTSS), treat it as the maximum redshift z max (corresponding to the limit flux S) of the theoretical log N − log S curve (see Figure D1), and calculate the slope of log N − log S at this z max from the theoretical curve, which is regarded to be the flattest of the possible pointwise slopes of the complete source counts range (called 'the possible flattest slope' for short); then we fit our log N − log S data with a modified double power law model and calculate the slopes of this model at fluxes starting from the low-flux end upwards, and find at which flux level the calculated slope gets steeper than 'the possible flattest slope'; we take this flux level as the lower bound of the fitting range; at the high flux end, we only discard the points with so small source counts that large fluctuations appear (these points distribute unnaturally horizontally).As for the LoTSS-DR2 data set, the highest redshift is 0.33, 'the possible flattest slope' is −1.262, and thus the fitting range is 2.49 mJy ⩽ S ⩽ 12.48 mJy; likewise, the fitting range of the FIRST data set is 3.06 mJy ⩽ S ⩽ 12.20 mJy (corresponding to the filled symbols in the right panel of Figure 8, while the open symbols are not used for the fitting).In the fitting, the two data sets are used jointly.
The best-fit log N − log S curve to our data set is plotted in right panel of Figure 8, with best-fit C = 3.83 × 10 6 .The canonical Euclidean log N ∝ −1.5 log S relation is also plotted as comparison.

Assessing the completeness and depth
With the best-fit log N − log S curve in hand, we can assess the completeness of our radio sample, following the approach of White et al. (1997, see their §6).Just as we noted in §4.3.2, the log N − log S distribution of the present sample, including its turnover (getting flatter) toward lower-flux end, is not affected by the low-luminosity turnover of radio luminosity function.From the right panel of Figure 8, the distribution of the sources from the FIRST deep mining agrees well with the best-fit log N − log S relation downwards till flux density of 2.79 mJy at 0.15 GHz (corresponding to 1 mJy at 1.4 GHz),9 and gets flatter than the relation below this flux density, which indicates that the limiting flux density for the source completeness of our FIRST deep mining is 1 mJy at 1.4 GHz.The thus-determined limiting flux density is half of that of the FIRST official catalog (≈ 2 mJy at 1.4 GHz; see §6 of White et al. 1997).Note that the limiting flux density of completeness is not the same concept as the detection limit defined in the FIRST official catalog (1 mJy) that is 5 times the rms noise plus the CLEAN bias (0.25 mJy).Besides, the limiting flux density of completeness concerns the integrated flux density of a source, whereas the detection limit is generally applied to the peak flux density (see, e.g., the definition of the official FIRST catalog; §4 of White et al. 1997); the integrated and peak flux densities of a source are not the same, when the source is resolved.Likewise, the limiting flux density for the completeness of our LoTSS-DR2 sources turns out to be 1.26 mJy at 0.15 GHz (namely, 0.45 mJy at 1.4 GHz), albeit conservatively determined which is broadly consistent with the LoTSS-DR2 official estimation, 1.1 mJy for the 95% completeness (see §3.6 of Shimwell et al. 2022).
Again, note that this completeness limit should not be confused with the 5σ detection limit (0.415 mJy).The nominal detection limit concerns peak flux density, and is just a median because the rms noise varies across the LoTSS-DR2 sky region (see §3.6 of Shimwell et al. 2022 for the detail).

The predictions
Based on the best-fit log N − log S curve, we can also make some predictions about the respective numbers of radio detections of our parent sample corresponding to various detection limits (as well as specific survey areas).Certainly for this purpose there is a premise that as far as the total radio emission is concerned these Seyfert 1 galaxies can be regarded as a same population; as stated in the beginning of this Subsection ( §4.3), this is a reasonable assumption.
First, we predict how many the detected LOFAR counterparts will be once LoTSS is completed.This can be made by keeping the current limiting flux density of completeness of the LoTSS-DR2 sources (1.26 mJy; see the above) and simply expanding the survey area.That is, multiplying the source number on the best-fit log N − log S curve at the limiting flux density of 1.26 mJy (see the vertical axis of the right panel of Figure 8) by a factor of the area ratio (8032/7813).Here 8032 deg 2 is the total area of Spectroscopic Legacy of SDSS-DR7, which our parent sample covers; the number 7813 is explained in §4.3.1.As a result, we predict that once LoTSS is completed about 149 radio counterparts of our parent sample will be detected.
We have tested and verified this prediction in the following way: using the sources from LoTSS-DR1 (cf.Table 2) to "predict" the number of the sources from LoTSS-DR2.For this purpose, we fit the log N − log S model to the data set of LoTSS-DR1 and FIRST sources (see the last paragraph of §4.3.1, and the left panel of Figure 8), and find that the limiting flux density of completeness for LoTSS-DR1 is 0.79 mJy at 0.15 GHz (cf. §4.3.3).We now scale vertically this best-fit curve by the area ratio 4111/7813 to correspond to the actual sky area of LoTSS-DR2 (see §4.3.1 for the two numbers and the detail of Figure 8).Reading from this scaled curve (at the limiting flux density of completeness 0.79 mJy), we can "predict" (i.e., based on LOTSS-DR1) that we can find 130 counterparts in LoTSS-DR2.This "prediction" agrees reasonably with the actual number of LoTSS-DR2 sources (112), considering a worse (and not homogeneous) complete-ness of LoTSS-DR2 than LoTSS-DR1 (see Shimwell et al. 2022Shimwell et al. , 2019)).If we use the nominal limiting flux density of completeness for LoTSS-DR2 determined in §4.3.3, i.e. "predicting" based on LoTSS-DR1 a posteriori, then the predicted number of LoTSS-DR2 counterparts would be 68, almost the same as the actual number of LoTSS-DR2 counterparts with flux density ⩾ 1.26 mJy (64 sources).
The second prediction is about the question: to what a limiting flux density of completeness can we find (in the statistical sense) all radio counterparts of the parent sample?First, we scale vertically the log N − log S curve presented in §4.3.1 (this curve corresponds to a sky area of 7813 deg 2 ), to match the total sky area of 8032 deg 2 .Then, along the curve, a limiting flux density is found corresponding to the source number of 513, which is 0.5 mJy at 0.15 GHz.Here we should note that the actual limiting flux density may turn out to be somehow lower; this is because our parent sample being optically selected Seyfert 1s, just like normal optical surveys being magnitude-limited, is not complete at high-redshift end.In addition, just as stated in §4.3.2, the low-luminosity radio sources with P 20cm ≲ 10 19.5 W Hz −1 may have a smaller space density than the power-law extrapolation of the luminosity function of their neighbors not too luminous away.If so, and if the radio powers of some yetto-detect sources in our parent sample are below the turnover luminosity, then the actual limit of completeness should be further lower.
5. SUMMARY In this work, we have built a radio sample of low-mass AGNs, by searching in large-scale radio continuum surveys for the radio counterparts of 513 optically selected low-mass AGNs (Dong et al. 2012;Liu et al. 2018;the parent sample).By exploiting the recent data releases of the surveys of SKA pathfinders, particularly the second data release of LO-FAR Two-metre Sky Survey (LoTSS DR2), we have compiled a catalog of 117 radio low-mass AGNs with 5σ detection (called the SKAPs subsample; Table 2), which constitute the main body of our radio sample.We have also deeply mined the FIRST images with an elaborate procedure specialized for fitting faint radio sources, and compiled a catalog of 63 radio low-mass AGNs with S/N ⩾ 3.5 (called deep FIRST subsample; Table 3).
Combining the two catalogs, our radio sample comprises 151 distinct radio low-mass AGNs, including 102 new reliable sources (S/N ⩾ 5) and 23 new candidates (3.5 ⩽ S/N < 5), as well as 26 known sources (see Liu et al. 2018).Thus our sample is four times the record of radio low-mass AGNs (totaling about 40 sources previously).Moreover, the majority of our newly found radio sources (119 of 125) are fainter than the flux threshold of the official FIRST catalog (1 mJy at 1.4 GHz), extending deeper into the low-flux regime than those previous small samples of radio low-mass AGNs.
We have investigated the radio properties.Only a small fraction (19/151) of the radio sources exhibit "extended" morphology; in particular, J134738.24+474301.9 appears to have a jet on tens of kpc scale.Regarding the 125 newly found sources, their P 20cm ranges from 1.98 × 10 20 to 1.29 × 10 23 W Hz −1 with a median of 1.12 × 10 22 W Hz −1 , and their R ranges from 0.9 to 181.3 with a median of 12.0.By combining the radio low-mass AGNs and the radio counterparts in the FIRST catalog of a larger sample of broad-line AGNs (i.e., radio Seyfert 1s at z < 0.35), we have carried out partial correlation tests among P 20cm (or R), M BH and L/L Edd , and found that for low-z Seyfert 1s P 20cm correlates intrinsically and positively with M BH , yet only marginally with L/L Edd .
We have investigated the depth and completeness of our radio sample ( §4. 3.3), in terms of the log N − log S relation for the expanding, spatially flat Universe.For this purpose, we carefully designed a modified double power-law function as the log N − log S model, and fit it to the distributions of our LoTSS and FIRST sources.The limiting flux densities of completeness of the deeply mined FIRST sources is pushed down to half the value of the FIRST official catalog (2 mJy; note the difference between the limit of completeness and the detection limit, and see §6 of White et al. 1997).The limiting flux density for the completeness of our LoTSS-DR2 sources turns out to be 1.26 mJy at 0.15 GHz (namely, 0.45 mJy at 1.4 GHz).Thus our radio sample is complete to such a fluxdensity level that is deeper than White et al. (1997) by a factor of 4.
We also predict in terms of the log N − log S relation that in total 149 radio counterparts of the parent sample will be detected once LoTSS is completed.Further, when the limiting flux density of completeness is achieved to 0.5 mJy at 0.15 GHz, all radio counterparts of the 513 optically selected low-mass AGNs in our parent sample would be detected.
The present work kicks off our massive search for radio low-mass AGNs in the surveys of SKA pathfinders.More and deeper data are to be released in near future, probably even from SKA-1.We believe that it will be feasible to con-struct a sample of low-mass AGNs, well defined and uniformly selected both in the optical and in the radio.Besides, it will be valuable to search for new (and fainter) low-mass AGNs directly in the radio.
As introduced in §1, there are 40 reliable radio sources of low-mass AGNs prior to the present work, of which 26 are included in our radio sample (see Footnote 2).For the other 14 sources, we summarized their information from the literature in Table A2 below.B. PROBABLE RADIO SOURCES WITH 3 ⩽ S/N < 3.5 FROM VLA-FIRST IMAGES Below we list the 23 probable radio counterparts with 3 ⩽ S/N < 3.5 of the optically selected broad-line AGNs, which are deeply mined from the FIRST images as described in §2.3.Because there exist optical AGNs at their well-specified sky positions, according to Bayesian argument the genuineness of the radio emission is more plausible than traditional blind-source detections with the same S/N (see also White et al. 2007).Anyway, as optimal, probable candidates of radio low-mass AGNs, these sources deserve deep and high-resolution observations in the future.3) and ( 5) to ( 13), see Table 3 for the definitions.Column ( 10) is the same as the Fp_e column in Table 3, both being the local rms noise of the images.Column (4), redshift; retrieved from the parent optical catalogs.Note that J152637.37+065941.7 has been detected by deep VLA observation (Gültekin et al. 2022).(This table is available in its entirety in a machine-readable form in the online journal.A portion is shown here for guidance regarding its form and content.)C. LOW-MASS AGNS WITH RADIO UPPER LIMITS FROM LOTSS-DR2 LoTSS-DR2 is both deep and wide, covering now about half of the sky area of our parent sample.There are 104 low-mass AGNs in the parent sample that are covered by LoTSS-DR2 yet under detection.Here we provide the information of those upper-limit sources.The upper-limit flux density is set to be 5 times the local rms of each source, where rms is retrieved from the official rms map of LoTSS-DR2.1), identification number assigned in this paper for the upper-limit sources (Prefix 'U' means 'Upper-limit').Column ( 2) to ( 4), the same as in the above Table B1.Column ( 5), upper limit of the flux density at 0.15 GHz, set to be five times the rms value at the respective location of every source.Columns ( 6) and ( 7), rest-frame 20 cm power and radio loudness derived with the values in Column (5).Columns ( 8) to (10), virial mass of black hole, Eddington ratio and luminosity of broad Hα emission line, retrieved from the parent catalogs.Column (11), continuum luminosity at rest-frame 5100 Å (L5100 ≡ λL λ at λ = 5100 Å), derived from L Hα B .(This table is available in its entirety in a machine-readable form in the online journal.A portion is shown here for guidance regarding its form and content.)D. log N − log S RELATION FOR THE EXPANDING, SPATIALLY FLAT UNIVERSE Although the real Universe is likely flat for its 3-dimensional space with measured k = 0, usually called spatially flat or even spatially Euclidean in the literature, it is yet expanding or called time evolving (i.e., its 4-dimensional spacetime is curved).Thus we must assess the impact of the cosmological expansion (namely redshift effects) onto the log N − log S test; it concerns not only the measured frequency and fluxes related to the so-called K-correction, but also the geometrical relation between comoving volume (and thus density) and luminosity distance (see §11.1.1 of Peterson 1997 for the derivations). 10Here we present the necessary formulation for the analyses in §4.3.
Generally for an expanding universe, the flux measured at some observer-frame frequency ν 0 has the following relation with luminosity distance (d L ) or redshift (z): where we assume the radio continuum is modeled as a power law, L ν ∝ ν −α .On the other hand, for an isotropic universe, the total number of sources per unit solid angle (dΩ) detected above a specific limiting flux S (at the observed frequency ν 0 ) can be expressed as the following volume integral: where z max is the maximum redshift, at which the sources with luminosity L ν (ν 0 ) can be detected at the limit S. That is, z max (L ν , S) satisfying Equation D1, which reads where d L (z max ) is the luminosity distance corresponding to z max , the two being definitely linked as prescribed by the cosmological parameters.Because of the complex formula of d L (z max ), the relation between S and z max (Equation D3) can be solved numerically only, and generally the dependence of S on distance deviates from the simple inverse-square relation of the non-expanding Euclidean universe.For an spatially flat universe (k = 0), Equation D2 can be expressed easily using comoving density and comoving volume.Particularly, when a uniform comoving source density n(r) = n 0 is further assumed, it becomes in a simple form Again, as discussed above, because z max (S, L ν ) is a complicated function of S (Equation D3), generally Equation (D6) deviates from N ∝ S −1.5 , more or less.We numerically solve the equation system (D3 and D5) for the ΛCDM cosmology adopted in this work with the AGN continuum slope we used (α = 0.46), and obtain the relation between log N (or log d 3 co equivalently) and log S; see the data points (solid circles) plotted in Figure D1.
In fact, we use log N − log S/S0 L/L0 relation as displayed in the figure, after careful considerations.Our primary consideration, as discussed in §4.3, is to keep the overall shape of the N(S) curve independent of the concrete values of S and L as well as n 0 .In other words, we do our best to achieve the following goal: the effects of S, L and n 0 to the curve shape is just to shift the N(S) curve vertically or horizontally in log-log scale.Thus we make the following two manipulations.(1) To use S/L instead of S, which is in fact implied by Equation D3.By doing so, in log-log scale the N(S) curve only shifts horizontally with changing S or L, keeping its overall shape unchanged.(2) To normalize S and L. The concrete values of the normalization factors S 0 and L 0 , just like the above point (1), do not impact the curve shape, and are somehow arbitrary in principle.In practice, a good choice is to fix them to be the mean values of the data set under consideration, because this choice is a natural rule of thumb of numerical calculation with computers, and makes easy for us to incarnate a simple (and elegant) analytical form to approximate N(S) (see below).To make Figure D1, we set S 0 and L 0 to be 1 mJy and 10 22.57 W•Hz −1 , respectively, adopt a single luminosity L to be L 0 simply (see §4.3), and set n 0 arbitrarily to be 3 Mpc −3 ; the data points are calculated by scanning z max (i.e., log S/S0 L/L0 equivalently) from 0.005 to 50.The first power law with the slope fixed to −1.5 is the only dominant component in the low-z end (i.e., high-S), which asymptotes to the Euclidean-universe case.The second, modified power law component of this model gets more and more significant as S approaches toward the low end, and works together with the first component to make the log N − log S slope get flatter and flatter and asymptotic to 0 as S tends to 0. The characteristic parameter d t accounts for the turnover of the curve, and the exponent γ controls the abruptness degree of the turnover.We fit the modified double power-law model (D7) to the numerically calculated data points of the expansion version of the log N − log S relation, and the fit is excellent indicating the aptness of this carefully devised model (see the solid line overplotted in Figure D1).The two normalization factors, S 0 and L 0 , are fixed to the mean values of the data set we used (see §4.3), as stated in the above.In principle, S 0 and L 0 can be arbitrary; but if so, the parameter d t would not be regarded as invariant, yet would have to be shifted by a quantity like (S 0 /L 0 ) γ when applying the above model to any data set.The parameters C is a free parameter to be fitted to match N. The best-fitted values of the three shape parameters, β, γ and d t , are listed in §4.3 already.Here we explain once again that we can directly apply the best-fitted d t value to §4.3, because we use the same S 0 and L 0 here and there.
We also plot the canonical log N − log S/S0 L/L0 relation (namely N ∝ S −1.5 , the static version) in Figure D1 as comparison, with the same S 0 , L 0 , L and n 0 values.Compared with the canonical relation, the expansion-modified version differs little in the bright S end (i.e., for the cases with very small z max ; see Equation D3) as expected, yet gets flatter (namely with less and less sources) toward the faint S end (i.e., for higher z max ).In the case of z max = 0.35 that is the maximum redshift of our parent sample (the optical low-mass AGNs), for the same flux limit (S) the expansion-modified version predicts 0.52 times the number (N) of the canonical relation.This deviation is too significant to be neglected.
), 3.80 ⩽ S/N < 4.38 (12), and 4.38 ⩽ S/N < 5 (6).The five stacked images based on the FIRST images are shown in Figure 2. The S/N of their central region are calculated, which are 8.

Figure 1 .
Figure 1.Determination of the S/N threshold for the radio sources deeply mined from FIRST images.Top panel: Cumulative histogram of the S/N values of the 86 FIRST counterparts for the parent sample with S/N ⩾ 3. Bottom panel: Corresponding to the right yaxis (red) is the radio detection-rate curve of Stripe-82 AGNs (red dash-dotted line); the S/N values (related to x-axis) are based on the fitting of their FIRST images, and the "true detections" and detection rate is defined in terms of the deep VLA-S82 catalog (see §2.3.2 for the detail).Corresponding to the left y-axis (blue) are the estimated numbers of true (blue squares) and false detections (blue circles), and the number difference (the true minus false; blue solid line), at various S/N thresholds for the 86 FIRST counterparts displayed in the top panel.

Figure 3 .Figure 4 .
Figure 3.The LoTSS-DR2 (left) and FIRST (right) images of J134738.24+474301.9, with the contours of its optical image overlaid.The flux density in the images are color-coded according to the respective colorbars.

Figure 5 .Figure 6 .
Figure 5. Normalized histograms of the rest-frame power at 20 cm for our radio sample (151 sources; black solid line) and all the other radio sources in the SDSS-DR7 broad-line AGN sample (dubbed BLAGNs in the text; black dashed line).Also plotted is the distribution of the upper-limit radio sources for our parent sample of optical low-mass AGNs (104 sources; blue dotted line).The histograms are normalized to make the probability (namely the area under every histogram) being unity.

Figure 7 .
Figure 7. Distributions of the rest-frame power at 20cm (P20cm; left panel) and radio loudness (R; right panel) on the diagram of black hole mass versus Eddington ratio for the radio sample of 151 low-mass AGNs, together with all the other radio sources in the SDSS-DR7 broad-line AGN sample (dubbed BLAGNs in the text).The radio power and loudness are color-coded according to the respective colorbars.

Figure 8 .
Figure8.The log N − log S distributions of the LoTSS sources (DR1, red upward triangles on left panel; DR2, red downward triangles on right panel) and FIRST sources (black circles on both panels).The flux densities are at 0.15 GHz; the source numbers, scaled in the LoTSS cases, correspond to the sky area of our FIRST sources (see §4.3.1).The canonical log N − log S relation for the static Euclidean universe (N ∝ S −1.5 ; dashed line in every panel) is plotted for comparative purpose.We also plot the log N − log S model (solid curve; designed for the expanding, spatially flat Universe) that is fitted to the aforementioned distributions of our radio sample.Note that only the data points denoted with filled symbols are used for the fittings.The dotted lines indicate the limiting flux densities of completeness for LoTSS and FIRST sources.
n 0 and n(r) are quantities with respect to the comoving frame (see Page 157 of Peterson 1997 for the explanation), and d co (z max ) is the comoving distance (namely proper distance at z = 0) corresponding to z max , and d co (z max ) = d L (z max )/(1 + z max ).In order to demonstrate the complicated N − S relation, by combining Equations (D3) and (D5) we obtain z max (S, L ν )) −1.5 (1+α) S −1.5 .(D6)

Figure D1 .
Figure D1.The numerically computed log N − log S/S 0 L/L 0 relation for the expanding spatially flat Universe.This relation (black dots) is in fact the log d 3 co − log S/S 0 L/L 0 relation with an arbitrary space density n0.The canonical log N − log S/S 0 L/L 0 relation for the Euclidean universe (dashed line) is plotted just for comparison.The modified double power-law model, which is used to approximate the numerically computed relation, is also plotted (solid curve).The vertical dotted line indicates where zmax = 0.35 .

Table 1 .
Statistics of absolute and relative differences of radio parameters

Table 2 .
Radio low-mass AGNs from the surveys of SKA pathfinders

Table 4 .
Physical properties of the radio low-mass AGNs

Table A2 .
Other radio sources in the literature

Table B1 .
Low-mass AGNs deeply mined from the FIRST images with 3 ⩽ S/N < 3.5

Table C1 .
Low-mass AGNs in the parent sample with 5σ upper limits from LoTSS-DR2