A search for kilonova radio flares in a sample of Swift/BAT short GRBs

The multi-messenger detection of GW170817 showed that binary neutron star (BNS) mergers are progenitors of (at least some) short gamma-ray bursts (GRBs), and that short GRB jets (and their afterglows) can have structures (and observational properties) more complex than predicted by the standard top-hat jet scenario. Indeed, the emission from the structured jet launched in GW170817 peaked in the radio band (cm wavelengths) at about 100 d since merger - a timescale much longer than the typical time span of radio follow-up observations of short GRBs. Moreover, radio searches for a potential late-time radio flare from the fast tail of the neutron-rich debris that powered the kilonova associated with GW170817 (AT2017gfo) have extended to even longer timescales (years after the merger). In light of this, here we present the results of an observational campaign targeting a sample of seven, years-old GRBs in the Swift/BAT sample with no redshift measurements and no promptly-identified X-ray counterpart. Our goal is to assess whether this sample of short GRBs could harbor nearby BNS mergers, searching for the late-time radio emission expected from their ejecta. We found one radio candidate counterpart for one of the GRBs in our sample, GRB111126A, though an origin related to emission from star formation or from an AGN in its host galaxy cannot be excluded without further observations.

Since the discovery of the binary neutron star (BNS) merger GW170817 in gravitational waves (GWs) and light at all wavelengths (Abbott et al. 2017a,b,c), the question of whether there may be a population of nearby, GW170817-like events in the known sample of short gamma-ray bursts (GRBs) has received particular attention (e.g., Horesh et al. 2016; Bartos et al. 2019;Dichiara et al. 2020;Schroeder et al. 2020;Matsumoto & Piran 2020;Ricci et al. 2021). Notably, GW170817 had a delayed electromagnetic afterglow (Abbott et al. 2017c), first detected in X-rays about 9 days after the merger and the prompt detection of γ-rays (Haggard et al. 2017;Margutti et al. 2017;Troja et al. 2017). A radio afterglow detection followed, about 15 days after the merger Hallinan et al. avery.eddins@ttu.edu 2017). This delayed, non-thermal, radio-to-X-ray emission was related to a structured jet launched after the merger and observed off-axis (e.g., Granot et al. 2018;Lazzati et al. 2017Lazzati et al. , 2018Mooley et al. 2018a,b;Nakar & Piran 2018;Ghirlanda et al. 2019;Makhathini et al. 2021). X-ray emission in excess to that expected from the structured jet was tentatively detected over 900 days after the GW170817 merger (Hajela et al. 2022;Troja et al. 2022). Possible explanations for this emission include a flat radio-to-X-ray spectrum afterglow from the slower kilonova ejecta that powered AT 2017gfo via the r-process (e.g., Arcavi et al. 2017;Chornock et al. 2017;Coulter et al. 2017;Cowperthwaite et al. 2017;Drout et al. 2017;Kasliwal et al. 2017;Kasen et al. 2017;Kilpatrick et al. 2017;Pian et al. 2017;Shappee et al. 2017;Smartt et al. 2017;Tanvir et al. 2017;Valenti et al. 2017;Metzger 2019;Radice et al. 2018Radice et al. , 2020, or radiation from accretion processes (Hajela et al. 2022;Metzger & Fernández 2021;Nedora et al. 2021;Troja et al. 2021). Balasubramanian et al. (2021) and Balasubramanian et al. (2022) found no evidence for radio emission in excess to that expected from a structured jet at late times, though efforts to detect a kilonova radio afterglow are still ongoing. O'Connor & Troja (2022) also did not find further evidence of an X-ray excess at 4.8 yrs since the GW170817 merger. Kilpatrick et al. (2021) observed the field of GW170817 in the optical band at late times, finding no remnant of the kilonova emission in their observations.
As evident from the above discussion, continued radio (and X-ray) observations are needed to confirm the presence of late-time emission from GW170817 and disentangle its origin. However, these results have spurred new interest in the search for so-called latetime radio flares from short GRBs (e.g. Nakar & Piran 2011;Hotokezaka & Piran 2015;Hotokezaka et al. 2018;Kathirgamaraju et al. 2019;Nedora et al. 2021), especially considering that the kilonova afterglow could be visible in radio for years after the initial event. Several recent efforts have targeted both well-localized short GRBs with known redshifts and short GRBs lacking accurate X-ray localizations and redshift measurements (e.g. Schroeder et al. 2020;Bruni et al. 2021;Grandorf et al. 2021;Ricci et al. 2021;Ghosh et al. 2022).
Searches targeting well-localized, cosmological short GRBs have been motivated by the expectation that radio flares can be brighter in the presence of a long-lived magnetar formed after the BNS merger and thus are potentially detectable also at cosmological distances (Ricci et al. 2021). For example, Schroeder et al. (2020) conducted 6 GHz observations of nine, low-redshift (z < 0.5) short GRBs with the Karl G. Jansky Very Large Array (VLA) on rest-frame timescales of ≈ 2 − 14 yr following the bursts. The lack of detections constrained the energy deposited into the ejecta to ≲ O(10 52 ) ergs (see also Ricci et al. 2021;Ghosh et al. 2022). A re-analysis of 27 short GRBs with GHz radio observations also enabled them to conclude that ≲ 50% of events could have formed a stable magnetar after the mergers. Grandorf et al. (2021) used the VLA to observe four short GRBs in the Swift/BAT sample without X-ray localizations. They found a previously uncatalogued radio source within the error region of GRB 130626 with a 3−6 GHz flux density consistent with a radio flare associated with a BNS at a distance of ≈ 100 Mpc. However, an origin related to a persistent radio source, not the GRB, could not be excluded.
Here, we present the results of continued efforts aimed at uncovering potentially nearby, GW170817-like short GRBs via a search for their late-time radio emission. Specifically, we carried out observations of 7 short GRBs in the Swift/BAT sample that lack accurate X-ray localizations and have unknown distances. Our short GRB selection criteria, as well as observations and analysis procedures, are very similar to those described in Grandorf et al. (2021). Our paper is organized as follows. In Section 2, we describe how we select our GRB sample. In Section 3, we describe our observations. In Section 4 we describe the counterpart radio candidates identification. In Sections 5 and 6, we present and discuss our results. In Section 7, we summarize and conclude.

GRB SAMPLE SELECTION
For this study, we consider only short GRBs in the Swift/BAT sample. We focus on short GRBs detected by the Swift/BAT because the VLA can cover the typical BAT localization error region in a single-field observation at 3−6 GHz. We consider events with T 90 ≤ 2 s as those have been traditionally linked to BNS merger progenitors (Nakar 2007), though see e.g. Rastinejad et al. (2022) for a recent discovery of a kilonova associated with the nearby, minute-long duration GRB 211211A.
We further select short GRBs with no (promptly) identified X-ray counterparts. Short GRBs in the Swift sample that have identified afterglows reside at larger distances at which the detection of long-term radio flares is unlikely with short VLA observations. Instead, we select GRBs with no redshift measurement because those could hide a population of nearby events. In fact, we note that without a GW detection, GW170817 would have remained one of the several short GRBs with no follow-up observations and no redshift measurement. This motivates our search.
We exclude from our short GRB sample events in regions of the sky inaccessible to the VLA, particularly those with declination below −40 • . In addition, we avoided the Clarke belt of satellites, ranging from a declination of −5 • to 15 • , because VLA observations of sources in the belt are subject to strong radio frequency interference (RFI) from satellites in C-band (4-8 GHz), X-band (8-12 GHz), and Ku-band (12-18 GHz).
Finally, we remove short GRBs for which existing observations would have already found evidence for kilonova emission had the GRB originated from a nearby BNS merger (e.g. Xu et al. 2014;Yurkov et al. 2012). See Bartos et al. (2019) for further details on the sample selection criteria. Table 1. VLA observations of the GRBs in our sample. All observations were carried out with the VLA in its C configuration (nominal FWHM of the synthesized beam of 3.5 ′′ at 6 GHz and 7.0 ′′ at 3 GHz). For each GRB field we report the observed central frequency, the image central RMS noise, the date and epoch of the VLA observation, the time spent on source with the VLA, the center of the BAT error region, the BAT position error radius, and the reference for the BAT position and position error. We carried out VLA observations of the 7 GRBs in our sample via project VLA/20A-239 (PI: Bartos). Each GRB was observed at nominal central frequencies of 6 GHz (C-band) and 3 GHz (S-band) with the VLA in its C configuration. The choice of frequency was motivated by the need to match the VLA field of view to the typical size of the error region of the Swift/BAT GRBs in our sample. Moreover, late-time radio flares from BNS mergers are expected to be powered by optically thin radio emission (Nakar & Piran 2011), thus lower frequencies offer better chances for discovery. On the other hand, the lower the frequency, the larger the expected number of unrelated sources in the field and the slower the time variability (hence, the harder the identification of false positives). Based on the above considerations, S-and/or C-band observations are overall advantageous compared to both L-band (nominal central frequency of 1.4 GHz) observations, or observations at frequencies above 6 GHz.
The VLA data were calibrated using the automated VLA calibration pipeline in CASA 1 . The calibrated data were then inspected manually for additional flagging, 1 https://casa.nrao.edu/index.shtml especially to mitigate RFI. The CASA task tclean was used in interactive mode to image the observed fields, with a robust parameter of 0.5 and Briggs weighting. Sources were identified in the cleaned images both by hand and by using Blobcat, a Python program designed to locate sources in an image, to ensure all sources were investigated (Hales et al. 2012a,b).
The central noise RMS of each image was determined using imstat by measuring the RMS of the residual image within a central circular region of radius 10 times the FWHM of the nominal synthesized beam. The efficiency-corrected RMS was then calculated by dividing this central RMS by the efficiency (primary beam correction) as measured at the location of each candidate. Our results are summarized in Table 1.

CANDIDATE RADIO COUNTERPART IDENTIFICATION
We visually inspected all of our calibrated images to search for candidate radio sources. Most sources identified via visual inspection were also confirmed by running Blobcat. We used the CASA task imfit to estimate the VLA position errors. Specifically, errors were calculated by dividing the clean beam semi-major axis, as measured using imfit, by the signal-to-noise ratio (SNR) Table 2. Candidate radio counterparts found within the BAT error regions of the GRBs in our sample. We list only radio detections that passed all criteria discussed in Section 4. For each candidate, we report the sky position, the class of the closest cataloged source, the epoch of the VLA detection, the central frequency of the observation, measured peak flux density, the offset between the measured VLA position and the position of the closest cataloged source, the estimated VLA position error, and the compactness parameter of the radio emission (see Section 4 for discussion).
101129A-Candidate-1 10h23m45.15s −17d35m51.20s IrS 58955.07 9.4 6.1 139 ± 11 .84 .25 1.12 ± .12 " " 58994.14 9.5 2.9 372 ± 29 " .48 1.25 ± .12 111126A-Candidate-1 18h24m20.56s +51d28m44.71s IrS 58962.53 8.4 6.0 105.2 ± 7.5 18 . The time between the BAT trigger and the mid-time of our VLA observation (see the "Epoch" column). The angular distance between the radio candidate and the closest known object in NED g VLA position error calculated as described in Section 4. h Compactness parameter calculated as described in Section 4. of a given source. Following Mooley et al. (2013), only radio sources found within the BAT error region of each GRB with SNR ≳ 7 were considered reliable detections. Hereafter, the SNR is defined as the ratio between the source peak flux density, and the peak flux density error calculated as described in what follows.
The peak flux density of each radio source was measured using imstat and a circular region centered on the source, with radius equal to the FWHM of the nominal synthesized beam. The peak flux density error was determined by adding the efficiency-corrected RMS in quadrature with the absolute flux calibration error. The last is estimated to be 5% of the peak flux density for observations using 3C 286 as absolute flux calibrator, and 10% of the peak flux density if the flux calibrator was 3C 48 (due to a recent flare).
For each of the radio sources with SNR≳ 7, we further used the CASA task imfit to calculate the integrated flux density. The integrated flux density error was determined by adding the error of the integrated flux density as reported by imfit and 5% (or 10%) of the imfit integrated flux density in quadrature, similarly to what done for the peak flux density.
Next, we derive the compactness parameter by dividing the integrated flux by the imstat peak flux density. We used the compactness parameter to determine whether each radio source is more likely to be a compact object, like a merger, or an extended object, like a galaxy (Itoh et al. 2020). We take 0.9 < C < 1.5 as the range of compactness values for a point-like source (Mooley et al. 2013). To ensure point-like morphology and avoid contamination from sidelobes, following Mooley et al. (2013), we also require that the size of the detected radio source, as reported by imfit, is smaller than 1.5× the size of the clean beam FWHM on corresponding axes and that the source is located at a distance ≳ 20× the geometric mean of the clean beam from sources with peak flux density ≳ 500 µJy and extended sources. We require all of the above conditions to be met at both 3 GHz and 6 GHz. This helps ensure that the 3 GHz and 6 GHz flux are both likely to be dominated by the same emission process, and that the spectral indices β (where β is defined so that f ν ∝ ν β ) can be reliably estimated.
Finally, we narrowed down our remaining candidate radio counterparts by checking for the existence of previously cataloged sources at their location (Grandorf et al. 2021;Mooley et al. 2013). We consulted NED (NASA/IPAC Extragalactic Database 2019), the VLA FIRST catalog (Becker et al. 1994), the VLASS quick look image repository (via CIRADA; Lacy et al. 2020), and the Chandra Source Catalog 2.0 (Evans et al. 2010). NED notably queries NVSS among other catalogs. If a source we observed with the VLA was found to be within 2 ′′ of a cataloged radio source detected in catalog images taken before the GRB trigger time, it was discarded from further analysis (as this indicates that the radio source is likely persistent, not a transient).
The results of this selection process are reported in Table 2. As evident from this Table, only 2 of the 7 GRBs in our sample are associated with candidate radio counterparts that pass all of the selection criteria described here. In the following Section, we discuss each GRB in more detail.

RESULTS
After observing the fields of the GRBs in our sample (Section 2) with the VLA, we calibrated our observations as described in Section 3 and selected candidate radio counterparts as described in Section 4. In what follows, we detail our analysis for each of the GRBs in our sample.

GRB 080121A
GRB 080121A triggered the Swift/BAT at 21:29:55 UT on 2008 January 21 (Cummings et al. 2008). It was detected at a refined, ground-based location of α = 09h09m01.8s and δ = +41d50m21.3s (J2000), with an error region of radius 2.48 ′ with 90% confidence 2 . We observed the field of this GRB at 6 GHz with the VLA for an hour (total time including calibration and overhead), starting at 23:04:52.0 UT on 2020 April 17. Our 3 GHz observations of the same field started at 18:00:54 UT on 2020 May 22, and lasted about an hour.
The central noise RMS we measure for this field is ≈ 39 µJy at 6 GHz, and ≈ 110 µJy at 3 GHz (Table 1). The relatively high RMS is due to a very bright source outside the BAT error region, which resulted in a limited dynamic range. None of the radio sources identified within the BAT error region passed the criteria described in Section 3. However, we note that this GRB was investigated further in Dichiara et al. (2020) and Ricci et al. (2021) due to its proximity to possible host galaxies, both of which are within 200 Mpc.

GRB 090417A
GRB 090417A triggered the Swift/BAT at 13:17:23 UT on 2009 April 17 (Mangano et al. 2009). It was detected at a refined, ground-based location of α = 02h19m58.3s δ = −07d08m28.9s (J2000) with an uncertainty of 2.8 ′ and 90% containment ). We observed the field of GRB 090417A with the VLA for one hour total time, starting at 20:21:37 UT on 2020 April 14 at 6 GHz and at 19:00:42 UT on 2022 May 23. The central noise RMS we measure for this field is ≈ 9.1 µJy at 6 GHz, and ≈ 16 µJy at 3 GHz (Table 1). No radio sources were detected within the BAT error region of this GRB.

GRB 101129A
GRB 101129A triggered the Swift/BAT at 15:39:32 UT on 2010 November 29 3 . It was detected at a refined, ground-based location of α = 10h23m41.0s, δ = −17d38m42.0s (J2000), with an uncertainty of 3.1 ′ and 90% containment. The VLA observed the GRB field for an hour, including calibrations, starting at 01:11:57.0 UTC on April 16th, 2020 for 6 GHz and starting at 02:49:57 UT on 2020 May 25 at 3 GHz. The central noise RMS we measure for this field is ≈ 5.7µJy at 6 GHz and ≈ 21 µJy at 3 GHz. Several radio sources were identified in the images collected for this field, but only one candidate radio counterpart passed the selection criteria described in Section 3 (  (Cummings et al. 2011). It was detected at a refined, ground-based location of α = 18h24m07.1s and δ = +51d28m06.1s (J2000) with an uncertainty of 2.5 ′4 . We observed the field of this GRB for one hour (total time), starting at 12:12:02 UT on 2020 April 23 at 6 GHz and at 08:32:43.0 UT on 2020 May 23 at 3 GHz. The central noise RMS we measure for this field is ≈ 4.8 µJy and ≈ 7.4 µJy at 6 GHz and 3 GHz, respectively. Several radio sources were identified within the BAT error circle, but only two candidate radio counterparts passed the selection criteria described in Section 3 ( Table 2).

GRB 120403A
GRB120403A triggered Swift BAT at 01:05:23 UT on April 3rd, 2012 (Beardmore et al. 2012). It was located at 02h49m49.8s +40d29m21.8s (J2000) with an uncertainty of 2.3 arcminutes and 90% containment (Sakamoto et al. 2012). The VLA observed GRB120403A for an hour, including calibrations, starting at 00:03:17.0 UT on April 16th, 2020 for 6 GHz and starting at 11:37:19.0 UT on May 23rd, 2020 for 3 GHz. The central noise RMS we measure for this field is ≈ 5.2 µJy and ≈ 15 µJy at 6 GHz and 3 GHz, respectively (Table 1). No radio sources were detected within the BAT error region of this GRB at either of the observed radio frequencies.

GRB 140606A
GRB 140606A triggered the Swift/BAT at 10:58:13 UT on 2014 June 6 (Stroh et al. 2014). It was located at α = 13h27m11.7s and δ = +37d35m56.5s (J2000), with an uncertainty of 2.4 ′ and 90% containment . We observed the field of this GRB with the VLA for one hour (total time) at both 6 GHz and 3 GHz, starting at 01:52:38 UT on 2020 April 25 and at 7:44:30 UT on 2020 May 29, respectively. The detected central RMS values are 5.5 µJy for 6 GHz and 9.7 µJy for 3 GHz (Table 1). Several radio sources were identified within the BAT error circle of this GRB, but none passed the selection criteria described in Section 3.

GRB 160726A
GRB 160726A triggered the Swift/BAT at 01:34:07.67 UT on 2016 July 26 5 . It was located at α = 06h35m14.3s and δ = −06d37m1.4s (J2000) with an uncertainty of 1.3 ′ (90% containment). We observed the field of this GRB with the VLA for one hour (total time) in each frequency, starting at 20:20:51.0 UT on 2020 April 18 at 6 GHz and at 17:55:04 UT on 2020 May 24 at 6 GHz. The central noise RMS we measure for this field is ≈ 7.9 µJy at 6 GHz, and ≈ 13.4 µJy at 3 GHz. No radio source passed the selection criteria described in Section 3 ( Table 2). 6. DISCUSSION

Spectral energy distribution
Using our observations at 3 GHz and 6 GHz, we compute the spectral indices of the three candidate radio counterparts listed in Table 2. We find β 1 = −1.32 ± 0.15, β 2 = −0.88 ± 0.14, and β 3 = −1.31 ± 0.17 for the radio candidates 101129A-Candidate-1, 111126A-Candidate-1, and 111126A-Candidate-2, respectively. While all three of these spectral indices are compatible with optically-thin synchrotron emission, as expected from radio from BNS mergers (Nakar & Piran 2011), they are also compatible (within errors) with radio emission from star formation (for which typical spectral indices are −1.1 ≲ β ≲ −0.4; Seymour et al. 2008). We note that the spectral index of GW170817 was measured to be β = −0.584 ± 0.002 (Makhathini et al. 2021). Therefore, candidate 101129A-Candidate-1 with β 2 ≈ −0.88 is the closest to GW170817, though our results all suggest steeper spectral indices. None of the spectral index values we derive for our candidate radio counterparts are suggestive of emission from flatspectrum AGN (β > −0.6; Itoh et al. 2020).

Contamination from unrelated (persistent or variable) radio sources
In this Section, we discuss the likelihood that the candidate radio counterparts we have identified in the GRB error regions considered in this work are false positives i.e., radio sources whose origin is unrelated to the GRB itself.
Using the spectral index of each candidate counterpart, we estimate 1.4 GHx flux densities in the range ≈ 0.4 − 1 mJy. According to Mooley et al. (2013) and Huynh et al. (2005), the average number of persistent radio sources (of any morphology) with a 1.4 GHz flux density above 0.1 mJy is ≈ 910 deg −2 . Hence, within the BAT error regions for GRB 101129A and GRB 111126A we would expect ≲ 5 − 8 unrelated persistent radio sources for observations conducted at ≳ 1.4 GHz (having also applied further cuts on their morphology). As discussed in Mooley et al. (2013), about ≈ 1% of unresolved radio sources above 40µJy at 1.4 GHz are variable at the 4σ level. Hence, we would expect an average of ≲ 0.05 − 0.08 variable unrelated radio sources in the error regions of our GRBs 6 . The Poisson probability of finding one or more unrelated variable sources would then be ≲ 5% − 8%, which is sufficiently low to motivate further follow-up studies in the radio aimed at establishing the level of time variability of our radio candidates.
To test whether our radio candidates can be related to star formation in unresolved galaxies, we first consider the constraints arising from the fact that these radio sources have point-like morphologies in our images. Short GRB host galaxies in the cosmological sample have effective radii in the range 0.2 ′′ − 1.2 ′′ with a median size of 0.36 ′′ . For the short GRBs with known redshifts, the median physical size is ≈ 3.6 kpc . Hence, if we considered a short-GRB-like host galaxy located at ≈ 200 Mpc, its angular radius would be of ≈ 3.7 ′′ . Our radio candidates have an angular size of ≲ 3.5 ′′ at 6 GHz (in C band) due to the requirement we imposed on their morphology and considering the nominal FWHM of the VLA synthesized beam in C configuration. Therefore, if our radio candidates are related to star formation in an unresolved host galaxy, most likely such a galaxy would be located at ≳ 200 Mpc.
Radio emission associated with star formation at GHz frequencies is dominated by synchrotron radiation from electrons accelerated by supernovae. The following relation can be used to estimate the star formation rate (SFR) in the galaxy given the measured luminosity at 1.4 GHz (Murphy et al. 2011;Perley & Perley 2013): Assuming the candidate radio counterparts listed in Table 2 are located at a distance of 200 Mpc, their estimated flux densities of ≲ 0.4 − 1 mJy at 1.4 GHz imply L 1.4GHz ≲ (1.8 − 4.6) × 10 28 erg s −1 Hz −1 and, in turn, a SFR rate of ≲ (1.2 − 2.9) M ⊙ yr −1 . This is compatible with normal galaxies and with cosmological short GRB host galaxies (0.2 − 6 M ⊙ yr −1 ; Berger 2009).
On the other hand, if we assume that our candidate radio counterparts are located at the median short GRB redshift of z = 0.72 (Berger 2014), then the measured fluxes would correspond to radio luminosities in the range (5.5 − 14) × 10 30 erg s −1 Hz −1 at 1.4 GHz, favoring an AGN origin.
The AGN scenario can be further tested using All-WISE (Wright et al. 2019) color information. Indeed, AGNs are expected to fall in the wedge defined by the conditions W 2 − W 3 > 2.517 and W 1 − W 2 > 0.315 × (W 2 − W 3) − 0.222 where W 1, W 2, and W 3 are the instrumental profile-fit photometry magnitudes 7 in bands one (3.4 µm), two (4.6 µm), and three (12 µm; Wright et al. 2010;Mateos et al. 2012;Gürkan et al. 2014). The 101129A-Candidate-1 candidate found in the error region of GRB 101129A has W 1 − W 2 = 0.10 ± .29 and W 2 − W 3 = 3.72 ± .30, fulfilling only the first condition of the AGN wedge. Therefore, a SFR origin is more We note the flatter behavior predicted for the light curve best fit at 200 Mpc compared to the one at 40 Mpc. Hence, we expect future (after 2022) radio observations to be able to probe more significant radio flux changes if the candidate is at 40 Mpc. See Section 6.4 for discussion.
likely. The candidate radio counterparts found in the error region of GRB 111126A did not have any color information available in the AllWISE Source Catalog.

BNS Merger Origin
To test the possibility that the radio candidates identified in our analysis are related to radio counterparts of BNS mergers powering the corresponding GRBs, similarly to what done in Grandorf et al. (2021), we fit 3−6 GHz model light curves to our flux density measurements. These light curves are derived assuming a twocomponent, GW170817-like ejecta model with masses of 0.04 M ⊙ and 0.01 M ⊙ and respective speeds of 0.1c and 0.3c. We set the fractions ϵ e and ϵ B of energy going into accelerated electrons and magnetic fields, respectively, equal to their fiducial values of 10% each. With these model assumptions, we have two free parameters left: the distance to the source and the interstellar medium (ISM) number density. Because we only have two observations per candidate radio counterpart, for each of these we fix the distance to a grid of values in between 40 Mpc and 200 Mpc and fit for the ISM density. Overall, our fits on the observed radio flux for GRB 101129A-Candidate-1 and for GRB 111126A-Candidate-2 return χ 2 > 10. Allowing for ϵ e and ϵ B to vary in the range 10 −4 -10 −1 returns best fits with ISM densities ≳ 10 cm −3 . For a more realistic range of 0.01-0.1 cm −3 (given that ∼ 80-95% of short GRBs have densities < 1 cm −3 ; Fong et al. 2015), we get χ 2 > 10. Therefore, we rule out these candidates as promising radio counterparts for late-time radio flares. We obtain better results for GRB 111126A-Candidate-1, which we show in Figures 1-2. All of the ISM density values we determined for GRB 111126A-Candidate-1 as a function of distance are plausible values within a BNS merger scenario (see Figure 2). Indeed, short GRBs have been associated with ISM media with densities ranging from 10 −4 cm −3 to 1 cm −3 , and GW170817 had a circum-burst density of ≲ 0.03 cm −3 (e.g. Halli- The 6 GHz discovery image of GRB111126A Candidate 1 is shown in Figure 3. If the hypothetical merger behind GRB111126A Candidate 1 is nearby, around 40 Mpc, our best fit model predicts rising light curves at 3-6 GHz (Figure 1, red and blue curves) and 7 − 10σ flux variations throughout 2023 and 2024 (compared to our observations in 2020). However, if the merger is farther away, say at 200 Mpc, the best fit light curves at 3-6 GHz flatten at late times and we would not expect to see significant flux variations in 2023-2024 compared to our 2020 observations (Figure 1, orange and light blue curves). Hence, re-observing this candidate in 2024 would be a promising strategy for further constraining its nature.

SUMMARY AND CONCLUSION
We used the VLA to observe 7 GRBs without accurate X-ray localizations to identify potential nearby, GW170817-like events. In the 7 observed fields, we find a total of 3 candidate radio counterparts (one in the error region of GRB 101129A and two in the error region of GRB 111126A) passing all our cuts. We have discussed these findings in the context of expectations for false positives, as well as in the context of BNS late-time radio flare models. Overall, one of the radio candidate counterparts found in the error region of GRB 111126A, GRB 111126A-Candidate-1, appears worthy of further follow-up in the radio. Indeed, detecting radio variability at the level of ≈ 7 − 10 σ appears possible over the next couple of years (at epochs of about 12 − 13 years since the GRB). The existence of time variability in the radio would significantly decrease the odds of a false positive origin.
We conclude by stressing that, if 10% of short GRBs in the known sample of events lacking a redshift measurement is located within 200 Mpc (Gupte & Bartos 2018), searches like the one described here and in Grandorf et al. (2021) could quickly become profitable: following up with the VLA a sample of ≳ 23 short GRBs should result in a ≳90% chance of finding at least one nearby, GW170817-like event. So far with this work and with the observations presented in Grandorf et al. (2021), we have collected about half of the necessary sample size using only filler-mode observations with the VLA and with a relatively small use of observing resources. With a higher observing priority, we could cover the whole sample at a faster pace. This would also enhance chances for well-timed follow-up observations aimed at testing expectations for radio flux variability (that could be critical to help discern the nature of any potential radio candidate). Overall, we envision this type of observational programs to continue in the future.
A.E. and A.C. gratefully acknowledge support from the National Science Foundation via grant #1907975. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. The authors thank the University of Florida and Columbia University in the City of New York for their generous support. I.B. acknowledges the support of the National Science Foundation under grants PHY-1911796 and PHY-2110060, and the Alfred P. Sloan Foundation. The Columbia Experimental Gravity group acknowledges the support of the support of the National Science Foundation under grant PHY-2012035.