{'fraction_non_alphanumeric': 0.03201116008049085, 'fraction_numerical': 0.02940631828450798, 'mean_word_length': 4.497387173396675, 'pattern_counts': {'":': 0, '<': 30, '<?xml version=': 0, '>': 8, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Investigating a population of infrared-bright gamma-ray burst host galaxies

We identify and explore the properties of an infrared-bright gamma-ray burst (GRB) host population. Candidate hosts are selected by coincidence with sources in WISE, with matching to random coordinates and a false alarm probability analysis showing that the contamination fraction is approx 0.5. This methodology has already identified the host galaxy of GRB 080517. We combine survey photometry from Pan-STARRS, SDSS, APASS, 2MASS, GALEX and WISE with our own WHT/ACAM and VLT/X- shooter observations to classify the candidates and identify interlopers. Galaxy SED fitting is performed using MAGPHYS, in addition to stellar template fitting, yielding 13 possible IR-bright hosts. A further 7 candidates are identified from previously published work. We report a candidate host for GRB 061002, previously unidentified as such. The remainder of the galaxies have already been noted as potential hosts. Comparing the IR-bright population properties including redshift z, stellar mass M*, star formation rate SFR and V-band attenuation Av to GRB host catalogues in the literature, we find that the infrared-bright population is biased toward low z, high M* and high Av. This naturally arises from their initial selection - local and dusty galaxies are more likely to have the required IR flux to be detected in WISE. We conclude that while IR-bright GRB hosts are not a physically distinct class, they are useful for constraining existing GRB host populations, particularly for long GRBs.


INTRODUCTION
Gamma-ray bursts (GRBs) occur when relativistic jets are launched by a newly formed neutron star or black hole, along our line of sight, in the immediate aftermath of a cataclysmic event such as the collapse of a massive star or merger of two compact objects. The distribution of GRB T 90 durations (the time over which 90 per cent of the gamma-ray radiation arrives) is indicative of two main populations. Long bursts (T 90 2s) are thought to arise from the collapse of particularly massive, rapidly rotating stars (Woosley & Bloom 2006;Levan et al. 2016;Schady 2017), while short GRBs (SGRBs) are associated with the merger of compact objects ( Troja et al. 2017;Goldstein et al. 2017;Tanvir et al. 2017b). Throughout, we denote short and long bursts in the text as SGRBs and LGRBs respectively. Because the progenitors of LGRBs are massive stars, and these only exist for a short while after formation, we expect to find that their host galaxies are actively star-forming. Multi-wavelength observations of LGRB hosts have found this to be the case. Within their hosts, LGRBs appear to trace the regions of highest UV luminosity, and hence star formation (Fruchter et al. 2006;Svensson et al. 2010;Starling et al. 2011;Blanchard et al. 2016;Lyman et al. 2017). Similarly, we would expect the comoving rate density of LGRBs to scale with the star formation rate density throughout cosmic history (Hopkins & Beacom 2006;Madau & Dickinson 2014). However, LGRBs have an apparent aversion to massive, luminous galaxies at fixed SFR, which has been interpreted as a host metallicity bias (Fruchter et al. 2006;Perley et al. 2013;Hjorth et al. 2012;Perley et al. 2016a). The bias has been confirmed spectroscopically (Krühler et al. 2015). For a massive star to produce a GRB, it is thought that the progenitor must have sufficient angular momentum (see Blandford & Znajek (1977), and Williams (1995)). This is difficult if the stellar metallicity is high, because metals increase the opacity of the stellar envelope, driving mass and angular momentum loss through line-driven winds (Vink & de Koter 2005;Vink & Harries 2017). Furthermore, if pre-burst mass loss produces a dense, high metallicity circumstellar medium, the opacity to gamma-rays is increased and may stifle the jet (Hjorth 2013).
Despite these theoretical considerations, the aversion to massive, dusty galaxies may also be due to observational biases (Perley et al. 2013). The picture of LGRBs favouring low mass, low metallicity hosts is complicated by the existence of dark bursts. These are GRBs which have a steeper X-ray to optical slope, β OX , than would be expected from extrapolation of their X-ray spectra (Rol et al. 2005). The implication is that the optical emission is suppressed. Studies have shown that ∼25 per cent of all bursts fall into this category (Jakobsson et al. 2004;Selsing et al. 2018). The proposed explanations for the high β OX values in this population include intrinsic faintness, the burst being at very high redshift, and dust obscuration within the host. The first of these is hard to explain on theoretical grounds, as the GRB and subsequent afterglow are though to arise from the same outflow of material interacting with the ISM; as such GRBs which are bright in X-rays ought to also be optically bright (e.g. Dainotti et al. 2015). If the burst is at very high redshift, observed optical bands correspond to rest-frame ultraviolet (UV), which is absorbed along the line of sight by neutral Hydrogen. However, such bursts are expected to be rare, and indeed only a small fraction of bursts have been confirmed at z > 5 ( Ruiz-Velasco et al. 2007;Tanvir et al. 2009;Salvaterra et al. 2009;Greiner et al. 2009;Cucchiara et al. 2011;McGuire et al. 2016;Tanvir et al. 2017a;Bolmer et al. 2017). The majority of dark bursts are thought to be optically faint due to dust extinction. Early studies of GRB host populations systematically missed dark burst hosts because localisation was performed with optical afterglows, particularly before the launch of Swift. More recent studies have attempted to account for this by selecting hosts of bursts which have an undetected optical afterglow, using X-rays for localisation instead. However, Krühler et al. (2011) and Perley et al. (2013) found that even when considering dark burst hosts, which are typically more massive and dusty, the overall GRB population still shows a bias towards fainter, less massive systems than the typical star forming galaxy population at the same epoch. Moving to host identification in the infrared (IR) may provide a route to further reducing this bias, since ultraviolet light from young stellar populations is preferentially absorbed by dust, and re-emitted in the IR. This only aids the identification of dark GRB hosts, however, if they are dark due to extinction from galaxy-wide dust. If the extinction is local to the burst site, or exclusively along the line of sight, then a burst might be optically suppressed in an otherwise IR-faint galaxy with little dust re-emission.
All-sky infrared surveys, such as the Wide-field Infrared Survey Explorer (WISE, Wright et al. 2010), can be used for the purpose of host identification. Such surveys are shal-low, favouring the identification of nearby (Kovács & Szapudi 2015) or luminous and dusty hosts. A small number of GRB hosts have been confirmed at z < 0.1, which we define as local. These include LGRBs 051109B, 060218, 100316D and 111005A (Perley et al. 2006;Liang et al. 2007;Starling et al. 2011;Micha lowski et al. 2016), while the host of LGRB 080517 was studied by Stanway et al. (2015a) following initial selection through coincidence with a notably bright source in the WISE bands. Subsequent follow-up resulted in characterisation of the stellar population and star formation rate in this galaxy through a number of indicators, including radio emission. It also secured the detection of molecular gas for only the third time in a GRB host, constraining the gas consumption timescale (Stanway et al. 2015b).
In general, the benefits of identifying local GRB hosts are threefold. First, as discussed, proximity makes observation at radio, submillimetre and infrared wavelengths more feasible (e.g., Micha lowski et al. 2015). This is exemplified by the recent identification of infrared molecular hydrogen emission lines in low-redshift GRB host galaxies (Wiersema et al, in prep). Secondly, local galaxies will tend to have greater angular extent and thus be easier to spatially resolve for GRB environment studies, increasingly using IFUs (e.g. Christensen et al. 2008;Starling et al. 2011). Finally, a rare class of low-luminosity long GRB (LLGRB) has emerged thanks to their low redshift identification (e.g. Galama et al. 1998;Stanway et al. 2015a). Because the supernovae (SNe) associated with LLGRBs appear typical of GRB-SNe across the full range of LGRB energies, it seems unlikely that the progenitors of LLGRBs are different to 'regular' LGRBs (Schady 2017). The question then is, what factors can produce the wide range of inferred LGRB isotropic energies, while influencing the range of SN energies much less? Suggestions have included the effect of viewing angle, differences between central engine activity duration versus the shock breakout time, and progenitor metallicity having an impact on burst efficiency (Hjorth 2013;Levan et al. 2016;Schady 2017). Studies of a large sample of low redshift LGRB hosts will be invaluable in determining the conditions capable of producing LGRBs, including low-luminosity bursts, as well as for studying the evolution of LGRB hosts over cosmic time.
In this paper we explore the properties of a population of IR-bright GRB host galaxies, detected in WISE. The hosts are localised in X-rays and are selected from all bursts detected between 2005-2016 inclusive. Such galaxies may be nearby, or extremely luminous and dusty. The paper is structured as follows. Section 2 describes the sample and selection criteria used. Section 3 details observations of a sub-sample of these candidate hosts, including VLT/Xshooter and WHT imaging and spectroscopy, in addition to ATCA radio observations. We compile archival and survey data in section 4. In section 5, previously studied hosts are identified. SED fitting is performed in section 6. Section 7 presents our results and discussion, with the broader implications considered in section 8. Our conclusions follow in section 9. Where required, the standard ΛCDM cosmology is used, with h = 0.7, Ω M = 0.3 and Ω Λ = 0.7. Magnitudes are quoted in the AB system (Oke & Gunn 1983).

Rationale
The Swift observatory (Gehrels et al. 2004), which has detected the bulk of GRBs since 2005, is mounted with an X-ray telescope (XRT, Burrows et al. 2004), a Gamma-ray Burst Alert Telescope (BAT, Barthelmy et al. 2005) and an Ultraviolet and Optical Telescope (UVOT, Roming et al. 2005), as on-board instruments. UVOT provides the best localisations, however only one third of bursts with an X-ray detection have a UVOT determined position. At the other extreme, all detected GRBs have a BAT detection by definition, however the localisation is no better than a few arcminutes. The best balance between the number of detections and the ability to locate a host is therefore provided by X-rays, for which ∼ 98 per cent of bursts have a localisation. We have identified a sample of infrared-bright gamma-ray burst host galaxies by cross-matching the GRB X-ray afterglow coordinates with the ALLWISE IR all-sky catalogue from WISE (the Widefield Infrared Survey Explorer, Wright et al. 2010). This provides aperture matched photometry in four wavebands, W1 − 4, at 3.4, 4.6, 12 and 22 µm. Any cross match procedure between these catalogues will identify both genuine matches and spurious matches to unassociated sources. Because the WISE dataset is relatively shallow, we expect to see nearby or very luminous extragalactic sources, in addition to Galactic stellar contaminants.

Initial Cross-Matching and Cuts
Data for GRBs (detected by Swift, INTEGRAL, Konus-Wind and the IPN) in the years 2005-2016 inclusive were downloaded from NASA's GRB catalogue 1 . Swift positions were checked against the Swift XRT-GRB catalogue 2 . The data include positions in the gamma-ray, X-ray and UV/optical bands, with their associated 90 per cent confidence error radii, in addition to the T 90 durations. We do not, at this stage, differentiate between long and short bursts. The total sample contains 1001 bursts, which are used for the following analysis. TOPCAT 3 (Taylor 2005) was used to cross-match the X-ray locations with sources in the ALL-WISE catalogue. Matching is primarily to the W1 band, i.e. all of our sources have at least a W1 detection. No significance cut in W1 was made at this stage, since the quality of the sources are determined through flags and visual inspection, as described later.
We perform an initial cross-matching analysis with a fixed radius for all bursts. In order to determine the expected contamination fraction, we also match to a catalogue of positions created by shifting the 1001 GRB positions by ± 1 and 2 deg in each of Right Ascension (RA) and Declination (Dec). Because the search radii used are of the order arcseconds, and the X-ray positional uncertainties are also on this scale, shifting by 1 or 2 deg removes all physical correlations and creates a random sample of coordinates. Crucially however, the broad distribution of points in galactic latitude and hence surface density is preserved. Various trial radii from 1 to 10 arcsec are tested. The difference between the number of matches to actual GRB coordinates, and to our 8008 new pseudo-random coordinates, is used to estimate a significance through the Poisson cumulative distribution, where the smallest P(≥N A |N R ) corresponds to the best matching radius. We find that r = 2.5 arcsec minimizes the contamination fraction. The significance is further improved by employing cuts. This includes the removal of sources contaminated by diffraction spikes, optical ghosts and similar data artefacts, using the ALLWISE contamination and confusion flag (CCF). If one of the W1 or W2 bands is dominated by contaminating flux, the match is rejected, or if both of these bands are contaminated (but not dominated), the match is rejected. Matched sources which do not satisfy these criteria cannot be considered robust or reliable. The random matches include both brighter and fainter objects than the actual matches, therefore we limit the random matches to the same range in apparent magnitude to evaluate the probability of selecting the same population by chance. The cuts effectively act as a signal-to-noise filter, with the lowest W1 SNR after cuts of ∼ 4.8. This gives a maximum centroiding error of ∼0.5 arcsec (given σ centroid = FWHM/2.35 * SNR). The X-ray positional uncertainties therefore dominate the cross-matching. A final cut was made by removing objects whose cross-matched counterparts were ambiguous or blended in W1 on visual inspection, but which otherwise satisfied the confusion flag cut. These are SGRB 060801 and LGRB 061007. Image cutouts of these (and all other) burst locations in the W1 band are given in the appendix. Given that visual inspection is both time consuming and subjective, this last cut was not applied to the randoms and as such all numerical comparisons between actual and random samples were made before this stage. With these cuts, and with a 2.5 arcsec matching radius, we find 45 matches to actual GRB positions and 23 to random coordinates. The corresponding Poisson cumulative probability P(≥N A |N R ) is ∼10 −5 . Using this methodology, we estimate a contaminant fraction f c of 0.51, with a Poisson 95 per cent confidence interval covering the region 0.36 ≤ f c ≤ 0.67.

Consideration of Burst Error Radii and Local Background Densities
In the previous section, we use the same matching radius for all bursts to give an estimate of the contamination fraction. However, this fails to consider two important factors. First, while 2.5 arcsec is the best matching radius when averaging over the GRB sample, individual burst error radii vary and can be larger than this, so we may be missing genuine matches which lie further out. Second, because GRBs occur in galaxies which tend to exist as members of groups and clusters, our previous analysis considering only the effect of galactic latitude on chance alignment probability is incomplete. The true chance of random alignment may be greater than suggested by averaging over degree-scales be- cause GRBs should preferentially occur in over-densities, which have angular scales much smaller than this.
To address the first issue concerning the tailoring of cross-matching radius to each GRB, we use a radius of 1.5×R 90 . This is approximately the 99 per cent error radius, assuming a Gaussian profile for the X-ray probability function. From 1001 GRBs, this yields 60 GRBs with one IR match and 4 with two or more. Some WISE sources are included in this count, and not in the 45 discussed previously, because their associated GRB has 1.5R 90 >2.5 arcsec and a match at r>2.5 arcsec. Others are not included because the matched radius from the previous analysis is greater than 1.5R 90 . These are GRBs 050716, 060428B, 070208, 120119A, 120612A and 161108A. We add these 6 bursts back into the sample -a small X-ray uncertainty is not used at this stage to reject an otherwise good match, because the source of the IR flux could plausibly be extended.
In order to estimate the chance that each association is genuine, we perform a false alarm probability (FAP) analysis. The surface density of sources in the entire ALLWISE catalogue is visualised in figure 1, as a function of latitude and W1 magnitude. The matched bursts are indicated. Clearly the apparent magnitude and galactic latitude both affect the probability that a match is spurious. However, sky object density also varies on small scales, in addition to the broad galactic latitude trend. To sample the local surface density Σ around each burst, we cross-match the X-ray coordinates for each burst with a 3 arcmin radius. Given that galaxy clusters have typical sizes of ∼10Mpc, a 6 arcmin diameter is a sufficiently small angular scale to sample density variation due to clustering and cosmic variance. This is demonstrated in figure 2, which relates redshift to angular extent θ for physical scales d of 2, 6 and 10Mpc, through θ = d D di a m , where D diam is the angular diameter distance. For the region around each burst, the local density of sources of magnitude equal to the WISE match or brighter is given by, where N(N m≤m g ) is the number of sources within 3 arcmin of the burst of W1 magnitude m g or brighter, and r 3 = 3 arcmin. The probability of a match at angular distance r being genuine and not a false alarm can be written as, which tends to 1 as r tends to 0, and tends to 0 as r tends to ∞, as required. Using this method, the false alarm probability (FAP) is given by 1 − P chance . The CCF flag cuts from section 2.2 are again used. After these are made, a cut of FAP<0.05 is chosen. This cut, when applied to the matching of random positions to ALLWISE, yields a theoretical maximum of 50 matches by chance. However, the distribution continues well below 0.05. The average FAP is therefore is much lower, and the number of false matches will also be lower. This is backed up by the addition of only 7 matches when going from a FAP cut of 0.025 to 0.05. Included in these 7 is LGRB 100316D, which has a previously noted z=0.059 host galaxy. Using a FAP cut of 0.05 allows us to catch hosts which have larger projected sizes, such as that of LGRB 100316D. In addition, it allows us to identify the hosts of bursts with large X-ray uncertainties, provided the field is not crowded and the WISE source is sufficiently bright. As in section 2.2, SGRB 060801 and LGRB 061007 are removed due to possible blending in WISE, in addition to LGRB 120811A. Three GRBs with more than one IR match (GRBs 060223B, 071007 and 071109) had all of their candidates rejected due to CCF or FAP cuts. The fourth example with more than one match, LGRB 050117, lies in the galactic plane (|b| = 3) and has two matches almost equidistant at ∼7 arcsec, with similar false alarm probabilities. We cannot distinguish which IR source is more likely to be associated, and the line of sight extinction meant there was no optical afterglow reported for this burst, precluding a improved localisation. Therefore, we reject LGRB 050117. This leaves us with a final sample of 55 bursts, each with one matched WISE source. This differs from the sample derived in section 2.2, in that ten extra bursts are Table 1. The 55 GRB X-ray positions for which a catalogued WISE source is identified within 1.5R 90 of the X-ray error circle centroid, and is not rejected based on FAP, CCF or blend cuts. Our classification of the candidates is given in the final column, and is discussed throughout the paper in the appropriate sections. SGRB 050724 is a disguised short burst, despite the long T 90 .

GRB
T LG -identified as an IR-bright host galaxy by comparison to the published literature. CA -rejected due to possible or confirmed chance alignment.

Sample Summary
The final sample includes candidate host galaxies for 55 GRBs. These are listed in table 1, with the candidates rejected for CCF flags, WISE blending and high false alarm probabilities given in table 2. The tables give the X-Ray coordinates, WISE coordinates, T 90 estimates, the 90 per cent confidence interval on the X-ray position, the X-ray-WISE separation and a false alarm probability for association with the WISE source. Given the analysis in section 2.2, we are confident that around a third to two-thirds of the associations are spurious. However, some will be Galactic stars, and others may be galaxies with properties inconsistent with being a GRB host. These contaminants can be identified as such through their photometric and spectroscopic properties, as well as through better burst localisation. Observations and archival searches for these observations are thus the next objective of this analysis.

OBSERVATIONS
In order to investigate the true hosts and determine which matches are spurious, we have observed subsets of the sample with WHT/ACAM, VLT/X-shooter and ATCA. 7 targets were observed with ACAM/WHT, 5 with VLT/Xshooter and 14 with ATCA.

WHT Imaging
Observations of 7 candidate hosts were taken over two nights (2015 January 19 and 20) with the auxiliaryport camera (ACAM) on the William Herschel Telescope (WHT). These were associated with programme WHT/2015A/34. Both nights were severely affected by Table 3. WHT/ACAM observations, taken on 2015 Jan 19/20. If a magnitude uncertainty is not given, the value corresponds to the 2 σ limit at the position of the WISE source. LGRB 070208, LGRB 061002 and LGRB 080307 were observed on 2015 January 20. Conditions were clearer but still windy, with poor seeing (∼2 arcsec). Sloan g, r, i and z filters were used. The images were reduced with standard IRAF procedures and aperture photometry performed on the candidate hosts. Aperture sizes were chosen to be ∼2 times the seeing FWHM (of the largest band) if the target was a point source, or else ∼2 times the FWHM of the object of interest. Aperture sizes were kept constant for each source. Details of the observations and measured quantities are listed in table 3, along with 2 σ depths and the seeing. The quoted magnitudes are in broad agreement with archival data where available.

WHT Spectroscopy
The candidate hosts of SGRB 111222A and LGRB 140331A were observed on 2015 January 19 and 20 respectively, using the V400 grism and a 1.5 arcsec slit on ACAM. The position of the slit with respect to the 90 per cent XRT error circle and WISE source are shown in figure 3, overlaid on r-band images. Given the poor seeing, slit losses were significant. The observations are listed in table 4. The LGRB 140331A candidate counterpart was not detected. The SGRB 111222A IR counterpart is identified as an M dwarf, as shown in figure 4.

VLT Spectroscopy
We observed five GRB host candidates using the echelle spectrograph on VLT/X-shooter (Vernet et al. 2011). Observations were associated with programme 096.D-0260(A) (PI: Stanway) and are detailed in table 4. Images in the rband and the position of the slit with respect to the XRT and WISE positions are also shown in figure 3. The XS images are used primarily for visualising the slit placement. The 2 σ depths of the XS images corresponding to GRBs 091102, 120224A, 120612A, 140331A and 141212A are 23.5, 22.6, 24.1, 22.9 and 21.7 respectively. In each case, the slit placement was chosen to overlap with the WISE source. The spectra were reduced using the standard ESO pipeline in Gasgano. Of the five targets, two were marginally detected (LGRBs 120224A and 140331A), one was detected with prominent emission lines (SGRB 141212A), and two were found to be foreground stars (LGRBs 091102 and 120612A). The two stellar spectra are shown in figures 5 and 6 respectively. We note that servicemode observations of the counterpart of LGRB 091102 has been misaligned, likely due to the misidentification of a selected offset and alignment star. Thus the spectroscopic identification of this as a star is irrelevant to the GRB and is presented here to avoid confusion in future studies of archival data for these observations. This is the only case where we are required to attempt photometry on the XS imaging, measuring an r-band magnitude of the faint WISE aligned source of ∼21.4 (anchored to the APOP magnitude for the nearby bright star, Qi et al. 2015). However, it is not detected in the other bands available, u and z. Combined with GALEX and 2MASS non-detections, we deem there to be insufficient data for fitting the SED of this object. We classify it as having no coverage in table 1.
The spectra of the LGRB 120224A and 140331A targets are featureless, with only marginally detected continuum flux and no readily identifiable absorption or emission lines. LGRB 120224A's candidate host has also been observed with X-shooter in a different programme, with similar results (Wiersema et al. 2012). If these are indeed the LGRB host galaxies, this requires them either to be mature stellar systems without nebular emission, or else heavily dust enshrouded. It should be noted that the WCS of the X-shooter images are misaligned with ALLWISE and the X-ray positions, leading to small offsets from their true positions. For example, the slit for LGRB 140331A has been deliberately placed over the fainter object south-east of the error circle, because the centroid of the IR flux aligns with it, suggesting that it corresponds to the source of the IR emission. While the other, brighter object might be the true host, it is likely not IR-bright, and would therefore be out of place in our sample.
The potential host of SGRB 141212A has a weak continuum with Hα, Hβ and O iii emission lines. The wavelength of these correspond to a redshift of 0.596±0.001. This is in agreement with Chornock et al. (2014), who observed an object within the enhanced XRT error circle one day post burst with the Gemini-N spectrograph. They found that, out of two objects near the error circle, the likely host has a redshift z=0.596. Portions of the 2D spectrum covering key emission lines are shown in figure 7. Emission line measurements are listed in table 5. Owing to the low signal-to-noise ratio, meaningful constraints on the Hα/Hβ ratio are not possible. However, we note that the presence of these lines is in qualitative agreement with the star-forming best-fitting SED as discussed later in section 6.

ATCA Radio Observations
Radio observations of 14 candidate hosts were made at central frequencies of 5.5 GHz and 9.0 GHz and a bandwidth of  Flux Hβ λ4862 [OIII] Doublet Mg λ5176.7 Hα λ6563 SII Doublet Mg λ517.7 Hα λ656.5 Ca Triplet Figure 5. The VLT/X-shooter spectrum of the object associated with LGRB 091102 target. Hβ, Mg, Hα and Ca absorption lines at negligible redshift confirm that this is an M-star. However, the slit was misaligned with the IR source, so we do not consider this identification any further.

Ca Triplet
Fe/Ca Blend Ca+ Doublet Figure 6. The VLT/X-shooter optical spectrum of the LGRB 120612A target. Hβ, Hα and various metal absorption lines at z≈ 0 indicate that this is a foreground star. None of the targets were detected. The observations are listed in table 6. We also list the synthesized beam size, which varied significantly from source to source given their wide range of declinations, and the final image RMS noise level. Where a redshift for the source is known, we use the 1.4 GHz flux to star formation rate (SFR) calibration of Kennicutt & Evans (2012) to estimate a 3 σ upper limit on the star formation rate (assuming a radio spectral slope of -1).

ARCHIVAL AND SURVEY DATA
Additional information for the candidate GRB host counterparts was gathered from archival surveys and the literature. The main source of optical photometry is the Pan-STARRS survey (DR1, Chambers et al. 2016), which covers the whole sky north of -30 • declination down to 3 σ depths of g, r, i, z, y < 23.3,23.2,23.1,22.3,21.3. Cross-matching between the Pan-STARRS 1 science archive and the 55 WISE counterparts was performed with a 2.5 arcsec matching radius, producing 27 matches. Of these, a small subset have more than one possible optical counterpart, and these are carefully considered in section 6.1. The candidate hosts for SGRB 141212A are visible in Pan-STARRS but below the cataloging threshold, so we measure the magnitudes from image cutouts, complementing our WHT photometry presented in section 3.1. The Sloan Digital Sky Survey (SDSS, Alam et al. 2015, we use DR12), VST/ATLAS (Shanks et al. 2015) and APASS (Henden & Munari 2014) were also searched by matching to within 2.5 arcsec of the WISE source. All of these surveys extend south of -30 • declination. Two GRB locations have matches in VST/ATLAS, SGRBs 070724A and 150101B, which are also covered by Pan-STARRS, and only one in APASS (LGRB 140927A). There are 12 matches in SDSS, which provides the only optical photometry for LGRB 161108A. In the remaining 11 matched cases, Pan-STARRS data also exists, and we use the best available combination of photometry. At the very least, u band limits are used from the SDSS matches. Overall, we have a total of 29 optical survey detections.
The remaining 26 positions may fail to obtain a match because the source lies outside the Pan-STARRS survey region (14 objects), or because the galaxy is optically faint (12 instances). Where a GRB falls in a field covered by one of the surveys used, but no object is detected at that position, 2 σ upper limits are used. Four host candidates lacking coverage have been studied in the literature (GRBs 050219A, 070429B, 100316D and 130515A), so 10 are classified 'NC' (no coverage) in table 1.
A total of 6 of the 12 optically undetected sources have been well studied in the literature. In these cases we use the results of those works. The total number of WISE sources for which we have survey coverage but are lacking information on an optical counterpart, is therefore 6. The 12 undetected sources are discussed further in sections 5 and 7.
The GALEX All Sky Imaging Survey (AIS, Martin & GALEX Science Team 2005) provides UV photometry or limits for all of our objects at near-constant depth (Martin et al. 2003). GALEX has two photometric bands, the Far and Near UV (FUV and NUV), with effective wavelengths of 1528Å and 2271Å respectively. Ten of the sample have a NUV source within 5.3 arcsec, the NUV PSF FWHM. Ex- panding the search radius to 10 arcsec yields only one more match (at 9 arcsec), suggesting that those matches identified are genuine. Of these, four (LGRB 080405, LGRB 080517, LGRB 100316D and SGRB 150101B) also have a FUV detection. Where we have no detection, we use the AIS mean 2 σ upper limits, m FUV = 20.89 and m NUV = 21.79. The 2 Micron All Sky Survey (2MASS) was used to provide NIR data or limits for our sample (Skrutskie et al. 2006). A catalogued 2MASS source was identified for 14 of the 55 WISE matches. Image cutouts for all the GRB positions were also inspected, and we measure JHK 2 σ upper limits for the remainder of the sample. We search the FIRST (Becker et al. 1994) and NRAO VLA Sky Survey (NVSS, Condon et al. 1998) radio surveys, cross-matching to the WISE coordinates, and find only one match. The host of SGRB 150101B is detected in NVSS at 1.4 GHz. The full table of FUV to W4 photometry derived from this compilation of observations and archival data is given in table A1 of the appendix.
Initial checks were performed with the available data to discern the physical nature of the sources. The first method employed uses the difference between PSF and Kron magnitudes (Farrow et al. 2014, and references therein), which is a recommended technique for star-galaxy separation in Pan-STARRS. Because the Kron radii vary with the light distribution of the object in question, and the PSF does not, extended sources such as galaxies (or saturated stars) show discrepancies between the two magnitudes. A plot showing star-galaxy separation in our sample using this method is given in figure 8. Beyond an i-band apparent magnitude of ∼21, the separation becomes unreliable. Additionally, PSF-Kron positions towards the lower right of the galaxy region might be contaminants, and we do not use positioning in this region as grounds for galaxy classification. For all sources, a visual check for extension was also made. It should be noted that while the PSF-Kron method can confirm an object as a galaxy, it cannot definitively classify stars. In particular, compact, dwarf or distant galaxies may be unresolved at Pan-STARRS resolution and appear 'starlike' by this classifier. Where the only information we have is an insecure galaxy or star classification using PSF-Kron, and the object is not obviously extended by eye, we use SED fitting to distinguish the possibilities (see section 6).
An alternative method, proposed by Kovács & Szapudi (2015), separates stars and galaxies with the aid of the Jband. This is demonstrated in figure 9. At fixed W1 magnitude, objects that are bluer in W1-J are more likely to be galaxies. Kovács & Szapudi (2015) found that a cut at W1-J = 0.09 (AB magnitudes) is an effective star-galaxy separator, with a stellar contamination on the galaxy side of only ∼1.8 per cent. However, we caution that this technique was applied to brighter sources than we are dealing with, with the galaxies lying at a median redshift of 0.14, much lower than our in our sample. Because it is unknown how the W1-J colours vary for lower mass, star forming galaxies across a range of redshift, the use of the cut here is a suggestive, but not decisive, diagnostic.
We also check for proper motion in WISE classifying those with notable proper motion (we require a total proper motion of at least 2 σ significance) as stars. The HSOY catalogue (Altmann et al. 2017) is also searched. This is a precursor to the full Gaia DR2. We identify several optical Table 7. Initial star galaxy separation results using archival photometry and catalogued data products, for those sources that had sufficient data for at least one of the separation tests. Note that sources with PSF-Kron positions towards the lower right of the galaxy region may actually be stellar contaminants. The 5 objects listed below the line lack coverage in Pan-STARRS and SDSS, but are classified in other ways.

GRB
PM WISE PM HSOY PSF-Kron W1-J Type sources, all within 1 arcsec of their respective WISE match, which have proper motions allowing us to rule them out as stars. A small number of sources have notable PM in HSOY but not WISEṪhese might be chance alignments between a foreground star and a background IR source. However, due to large uncertainties on the WISE proper motions, we are unable to tell. Where this scenario arises, we assume that the optical counterparts are associated with the IR sources, given that this is most likely. Table 7 gives the results for those sources that have data available for at least one of the star-galaxy checks discussed.

PREVIOUSLY REPORTED GRB HOSTS AND OBSERVATIONS
Many of the GRBs in our sample of 55 have previously been studied, yielding useful information for our analysis. In this section, we compile reported afterglow positions, detailed host studies and other noteworthy information. This allows us to rule out some IR sources as chance alignments, and add some hosts to our sample for which we lacked the required observations. We split these into 3 categories: matched IR sources with an optical detection, those without, and those lying outside the optical surveys searched in this paper.

Optical Survey Detections
LGRB 060428B: The candidate host represents a single IR source corresponding to what appears to be a single optical source. However, Perley et al. (2007) has suggested that a compact blue galaxy lies underneath the foreground elliptical's light, at an offset of 2.6 arcsec. This corresponds to the foreground object's Einstein radius, as such they claim that LGRB 060428B is likely a gravitationally lensed event originating from the higher redshift, bluer galaxy. While this work has not been fully published, it seems a plausible explanation. We therefore take a conservative approach and exclude this source from later analysis.
LGRB 160703A: Zheng et al. (2016) observed the afterglow with Keck-I in the g and r bands. The improved positional certainty over XRT suggests that the IR source is not the host. The object has proper motion (see table 7), confirming this interpretation. Therefore we remove the source from further analysis.
LGRB 161214A: The object associated with the matched WISE source was observed by Malesani et al. (2016), who obtained spectroscopy identifying it as a K or M star.
There are cases where an improved burst localisation (e.g. in the optical) strengthens the IR source association, rather than ruling it out. There are also hosts where we have sufficient photometry for SED fitting, in addition to previously reported host parameters. We use this information to compare to our SED fitting results. These studies and observations are referenced for each burst in table 8 of section 6.

Optical Survey Non-Detections
LGRB 050716: No object is detected in optical imaging. Rol et al. (2007) deduced a burst redshift > 2 based on the optical to X-ray afterglow SED, however the matched WISE source has proper motion -assuming that the flux is entirely from this star, we reject this association as a chance alignment.
LGRBs 070208 and 120119A: Blanchard et al. (2016) report HST imaging of the burst locations in these cases, with the optical afterglow positions indicating that these may be chance alignments with the IR-bright sources. This is backed up by our analysis in section 6.
LGRB 080207: We successfully identify this host, which has been extensively studied in the literature (Hunt et al. 2011;Svensson et al. 2012a;Hjorth et al. 2012;Krühler et al. 2012b;Arabsalmani et al. 2017) as an example of a red, dusty luminous infrared galaxy. There is a 1 arcsec separation between the WISE centroid and the Chandra position provided by Svensson et al. (2012a), which itself is clearly placed over the galaxy in question. Therefore, we include the reported parameters for this galaxy in our analysis. The optical faintness of this galaxy shows that such sources can have steep optical to IR spectral slopes, and demonstrates that other optical non-detections could be similar in nature.
LGRB 080307: This burst has an X-ray detected AGN a few arcseconds away (Page et al. 2009). Both the WISE IR emission and X-ray flux have levels consistent with expectations for local AGN (Eckart et al. 2010). Therefore, we suggest that the AGN is the most likely source of the IR flux. While it is possible that the AGN resides in the host galaxy (as is the case with SGRB 150101B, see Fong et al. (2016)), we use MAGPHYS for our SED fitting in the next section, which does not have a prescription for AGN. Because the IR flux is consistent with being AGN dominated, this would lead to incorrect parameters when the SED is fitted by MAG-PHYS. Therefore, we take a cautious approach and remove LGRB 080307 from further analysis.
LGRB 080605: In the imaging provided by Krühler et al. (2012a) and Blanchard et al. (2016), we can see that the WISE flux centroid is centered on a z = 1.64 galaxy undetected in Pan-STARRS, rather than either of the two bright sources. The burst redshift of z = 1.64 was determined from afterglow spectroscopy. Because this galaxy is therefore confirmed as the host, we reject the nearby bright objects and reclassify this host as a Pan-STARRS non-detection. Krühler et al. (2012a) provide estimated parameters for the host galaxy, which we employ.
GRBs 100816A and 110918A: For these sources we use the physical parameters reported by Pérez-Ramírez et al. (2013) and Elliott et al. (2013) respectively. In both cases, the IR bright source is aligned with the reported host galaxy.
GRB 120224A: This object is only marginally detected in our X-shooter spectroscopy, however Selsing et al. (2018) found a 2 σ emission line in similar X-shooter data. If the line is Hα, this corresponds to a redshift of 1.1. We tentatively assume this to be the case going forward, treating the source as an optically undetected galaxy at z = 1.1.
LGRBs 130527A, 150323C and 161007A: The positions of these bursts were observed in the optical. In each case, an extended object was seen, allowing galaxy classifications (Cano et al. 2013;Malesani et al. 2015;Heintz et al. 2016). There is insufficient photometry for SED fitting, so we continue to treat these as non-detections, but note the IR flux likely originates from faint galaxies (rather than stellar contaminants).
LGRB 130528A: While undetected in Pan-STARRS, this source is revealed to be in a crowded region in deeper imaging, with multiple objects in a 10 arcsec region (Jeong et al. 2014). Because we cannot assign the IR flux to a single object with any certainty, we classify this burst as a potential chance alignment and do not consider it any further.

Lacking Coverage
GRBs 050219A, 070429B, 100316D: We use the host galaxy physical parameters reported by Rossi et al. (2014), Cenko et al. (2008) and Starling et al. (2011) respectively. In these cases, we compare the reported host coordinates to the WISE positions. For all three, the IR flux is aligned with the galaxies identified as a potential hosts. Because these are good quality IR sources which satisfy our FAP cut, we include the reported host parameters in our analysis.
LGRB 130515A: Levan & Tanvir (2013) observed the brightest source in the X-ray error circle with VLT/FORS2, find-ing it to be an M star. The position of the star is consistent with the WISE source, so we discount this association.

Multiple Candidate Hosts
In our sample, there are cases where multiple optical sources lie inside the search radius used to match optical imaging to WISE. Visual inspection reveals 3 GRBs (141212A, 140331A and 150120A) where the candidate host IR flux is not uniquely associated with a single optical source. In these cases, the optical and IR images were aligned to check for astrometric offsets. For LGRB 140331A, the IR emission originates from the fainter of two optical sources within the XRT error circle. Inspection of Chandra X-ray imaging with sub-arcsecond afterglow localisation suggests that the burst may not, in fact, be associated with either optical galaxy, but we cannot rule this out (Chrimes et al., in prep). The optical source aligned with the IR flux matched to LGRB 140331A does not reach the threshold for cataloging in Pan-STARRS. Instead, and in addition to the WHT photometry (see section 3.1), we measure magnitudes from Pan-STARRS cutouts. It is this photometry which appears in the appendix, and is used for SED fitting.
In the case of SGRB 141212A, the IR source centroid lies closest to the host identified by Malesani et al. (2014), although there may be some blending with neighbouring objects. We download and perform aperture photometry on the Pan-STARRS images in addition to the WHT images previously discussed. The photometry from both is in good agreement, and the Pan-STARRS measurements are again given in the appendix.
The IR source associated with SGRB 150120A may have either of two optical counterparts, both catalogued in Pan-STARRS DR1. We fit their SEDs separately, in each case assuming that the entire IR flux is associated with the galaxy under consideration. The better of the two fits us used in the subsequent analysis.

Galaxy SED fitting with MAGPHYS
Spectral Energy Distribution (SED) fitting on this sample was performed by χ 2 minimisation using templates derived from MAGPHYS and its high redshift update (da Cunha et al. 2008(da Cunha et al. , 2015. MAGPHYS was chosen for fitting our sample due to its careful implementation of dust extinction and re-emission. Given that our sources are detected in the infrared, we expect dust to be an important influence on their spectra. The population synthesis models of Bruzual & Charlot (2003) are employed, with the dust absorption and re-emission model of Charlot & Fall (2000). The stellar populations are built up by assuming a range of ages distributed evenly from 0.1 Gyr to the maximum age permitted at a given redshift (i.e. the age of the Universe). The star formation rate is modelled as declining proportional to e −γt , where γ is the star formation timescale and t is the time elapsed since the onset of star formation. Random constant SFR bursts are overlaid, with durations evenly distribution between 3×10 7 yr and 3×10 8 yr. The probability of a burst is such that 50 per cent of the model galaxies have experienced a starburst phase in the last 2 Gyr. The amplitude of these bursts A is defined as the ratio of stellar mass formed in the burst and all stellar mass assembled since the galaxy was formed at time t form ; this parameter is distributed logarithmically from 0.03 to 4.00.
Many of the candidates already have confirmed or likely redshifts, either from afterglow absorption lines or from the host. These were used to fix the redshift where possible. We also derive a photometric redshift for each galaxy. At each step in a grid of trial redshifts, we make use of the internal Bayesian fitting code built into MAGPHYS to determine the redshift interval for which χ 2 ≤ χ 2 min + ∆ χ 2 (Avni 1976). We did this at redshift intervals of 0.05 over the range 0 ≤ z ≤ 1, or up to z = 3 if no acceptable solution is found at lower redshifts. The distribution of χ 2 over redshift is then minimized. In this way, MAGPHYS is effectively used a photometric redshift code (da Cunha et al. 2015), with the redshift treated as an additional free parameter. Our photometric redshifts are generally in agreement with spectroscopic redshifts where available, as shown in table 8. Uncertainties on SED fitting parameters include the effect of photometric redshift uncertainty where a spectroscopic (or SDSS) redshift is not used.
MAGPHYS fits for stellar mass, current star formation rate, star formation history (age, timescale and burst amplitude), metallicity, and dust extinction, amongst other parameters. We caution that the number of free parameters is greater than the number of data points available for a given galaxy and that the fits may be overconstrained. MAG-PHYS does not routinely report uncertainties on metallicity, as this is often a poorly constrained parameter, so we simply state the best-fitting metallicity given by MAGPHYS.

Galaxy SED fitting Summary
Of the 55 GRBs with candidate WISE counterparts, 29 have optical-NIR photometry (from Pan-STARRS, SDSS, APASS, 2MASS, WHT observations or some combination of the above) in addition to at least WISE band 1. We correct this observed photometry for the Galactic extinction of Schlafly & Finkbeiner (2011) and a Fitzpatrick reddening law with R V = 3.1, using the IRSA dust reddening and extinction service 4 and the York Extinction Solver (YES, McCall 2004). We fit 28 SEDs (excluding LGRB 060429B, as this has been rejected as probable chance alignment). For the host of SGRB 050120A, we fit the two optical components separately, assigning the entire WISE flux each time. 14 of the 28 objects are best-fitting, or otherwise confirmed, as galaxies. Objects identified as galaxies or stars solely through the quality of fitting to stellar or galaxy templates are classified as photometric galaxies or stars ('PG' or 'Ps') in table 1. The limited number of stellar templates available makes it difficult to decide whether an object is best-fitting as a star or galaxy based on the reduced χ 2 alone. Therefore, we note the reduced χ 2 values for galaxy fitting, but make the assumption that a visually good star fit indicates that the object is a star, when the corresponding galaxy χ 2 value is poor. 161108A log 10 (νL ν / L ) Figure 10. SEDs for the objects which were best-fitting (or else confirmed) as galaxies. Some upper limits are too high on this scale to be visible. The data for the objects best-fitting as galaxies, together with the best-fitting galaxy templates, are shown in figure 10 and table 8. There are some instances where the FUV/NUV flux appears inconsistent with the SED. This might be because we have matched to a different objectthe matching to WISE is independent for the UV, optical and NIR. However, this unlikely as we would expect a nearby UV source to also be seen in the optical. They might be UV upturns, a phenomena seen in otherwise red elliptical galaxies, or the fitting may simply be failing to properly account for complex stellar populations.
Some fits produced no clear minimum in χ 2 , and therefore lack a robust photometric redshift. In these cases, we use a spectroscopic or SDSS photometric redshift where available. For LRGB 061002, we obtain a photometric redshift of 0.10 +0.45 −0.10 . However, the SDSS photometric redshift, which is calculated using machine learning techniques, is much more precise at z = 0.564. Therefore we fix our redshift to the SDSS value, which has sufficiently small uncertainties that the fit parameters are unaffected by variation with this redshift range. Similarly, the photometric redshift for LGRB 131122A is unconstrained, likely because there is no clear Balmer or Lyman break evident from the photometric points. However, this object also has an SDSS photometric redshift, at z = 0.399. This produces a visually acceptable fit, so we employ this redshift.
Another object, the candidate host of LGRB 070208, is correctly identified as a galaxy in our analysis but our photometric redshift of z = 0.16 +0.37 −0.16 is inconsistent with a previous spectroscopic redshift of z = 1.165 (Cucchiara et al. 2007). The foreground object is not detected in the spectroscopy reported by Cucchiara et al. (2007), but there is a marginally detected object seen offset from the foreground galaxy, nearer to the afterglow centroid. The HST imaging of Blanchard et al. (2016) confirms this interpretation. Given this information, we do not include LGRB 070208 in the discussion of our sample properties.

Stellar Fitting
We also perform fitting to a library of 131 stellar spectra (Pickles 1998), with the expectation that a significant fraction of our objects will be best-fitting as stars. The spectra span the range 1150-25000Å, allowing fitting from the FUV to K bands. In most of the 14 cases where MAGPHYS does not provide a good fit, stellar fitting does. The best-fitting stellar templates for these sources are shown in figure 11. Here, we show the SEDs for those objects with proper motion or other star diagnostics, in addition to sources where a stellar SED is a better fit than a galaxy SED. We note that the sources which fit best to stellar templates include all objects for which a spectroscopic or proper motion confirmation as a star is available.

SED Fitting Results
In table 8 we present the results of our SED fitting for the 13 candidate or confirmed GRB hosts that match best to galaxy templates (the 14 shown in figure 10, excluding the likely interloper LGRB 070208), and also a compilation of information from the literature for seven of the GRB hosts discussed in section 5.
In comparing the SED fit results to independent measures of the host properties, potential discrepancies arises between them. The SED derived SFR, ∼ 0.15 M yr −1 , for the host of SGRB 150101B disagrees with the NVSS 1.4 GHz radio SFR, ∼ 300 M yr −1 and the detection of this source at 5.8µm in the W3 band. This is due to the presence of an AGN . Unlike LGRB 080307 which was rejected due to coincidence with a nearby AGN, the association of SGRB 150101B with this galaxy has been secured with spectroscopic observations. Therefore, we keep it in the analysis. Only two other sources in our sample are detected in the W3 band. While this is sometimes used as a SFR indicator or AGN discriminator, the two cannot be unambiguously differentiated without a reliable (and low) redshift (Davies et al. 2017). We do not consider this data any further at this point.

Short GRBs
Before considering the redshift distribution and other properties of this sample, it is important to consider the selection effects that will shape any comparison we make. The first issue to consider is whether we may be particularly biased towards long or short bursts. Since their progenitor mechanisms differ, we expect their host properties to also differ. As such it is essential to consider the short vs long divide. Of the 55 GRB locations identified as having a WISE counterpart, 11 of the associated bursts had an observed T 90 < 2 s. Of these, 8 are reported in table 8 as galaxies. Two are identified as a Galactic stars.
The lack of a clear divide between the two populations can lead to ambiguity (e.g. Berger 2014, and references therein). Two of our targets have a 2 < T 90 < 5 s and might be classified as short by some proposed criteria. Of these, one is identified in our galaxy sample, while one has insufficient photometry. In the following, we treat these "intermediate" sources as long bursts.
If we consider the 20 GRBs for which galaxy properties are assembled in table 8, the fraction of SGRBs is 40 per cent (8/20) and would be as high as 50 per cent if the intermediate bursts were included. This compares to the short burst fraction in the entire GRB catalog of only 6 per cent (Berger 2014), suggesting that we may be preferentially selecting short bursts. This may reflect the difference in the underlying redshift distribution of these sources, exemplified by the mean redshift of z = 0.45 for the SGRBs in our sample. Short GRBs are typically of lower isotropic-equivalent luminosity and their distribution is biased towards low redshifts ( z < 0.8) relative to long GRBs ( z ∼ 2), due to the differences in both their progenitors and detection probabilities (Berger 2014). Given the relatively shallow depth of the W1 band imaging, we might expect a low redshift, and therefore SGRB, excess in our sample.

Redshift Distribution
We now consider the detailed properties of the IR-bright LGRB host population. Figure 12 shows the cumulative redshift distribution for the 13 candidate host galaxies that form our IR-bright LGRB host population (the 12 in table 8 plus LGRB 120224A at z = 1.1). We compared these to the SHOALS sample of LGRB hosts (Perley et al. 2016a,b), all Swift LGRBs with a known redshift, and a sample of dark LGRB hosts (Perley et al. 2013). In each case we indicate uncertainty on the cumulative distribution by performing an analysis in which each value is permitted to vary by addition of its associated random error, drawn from a skewed Gaussian distribution, in order to account for asymmetric errorbars. The scale and alpha parameters are chosen in each case such that the asymmetric distributions are reproduced. We show the standard deviation of 1000 realisations of the perturbed cumulative distribution as a shaded region.
In our sample, 7 objects have survey coverage but are undetected, and lack extensive study in the literature. It is possible that some of these may be M, L or T dwarfs, which LGRBs with known redshift, the unbiased LGRB SHOALS host galaxy sample and the dark LGRB sample of Perley et al. (2013). We make the worst case assumption that the 5 optical-non detections that have not been ruled out as stars are higher redshift galaxies.
can be a few magnitudes brighter in the WISE bands than the optical (e.g., Best et al. 2013). Indeed, the optically undetected object associated with LGRB 050716 has significant WISE proper motion. Late L and T dwarfs have the reddest colours of these stars and would be most able to satisfy the criteria of W1 detection and optically non-detection, however these are also the rarest classes of these objects. Alternatively, the optical non-detections could be intrinsically faint or higher redshift galaxies. This possibility has been demonstrated by several such examples in section 5. Three of the optically undetected sources have been identified as extended sources in deeper imaging, as described in section 5. One of the non-detections, LGRB 120224A, may be at z = 1.1 (Selsing et al. 2018). We can see in figure 12 that our sample is biased towards low redshifts when we do not included the non-detections. Therefore, the most extreme scenario is that all 5 of the remaining nondetections are in fact galaxies at higher redshift than the highest confirmed LGRB in our sample. Setting these 5 nondetections to an arbitrarily high redshift, we are still unable to match the slower rise of the SHOALS and dark samples, demonstrating that there is a low z bias even in the 'worst case' scenario. In reality, their redshifts could well be lower, given that the optical undetected hosts of LGRBs 080605, 100816A, 110918A and 120224A all lie at z < 2.
We must also consider the possibility that chance alignments remain in our sample. LGRB 140331A is at high risk in this regard, since the XRT error circle position favours a different optical source, unaligned with the WISE flux. In general, foreground chance alignments preferentially select lower redshift or foreground objects, and this remains a possible explanation for the difference in cumulative distribution between our sample and others. However, all but three of these sources have spectroscopic redshifts, and it is unlikely that a spectroscopic redshift measurement would be unaffected by a foreground interloper. We do, however, note that chance alignment is the most likely scenario for LGRB 070208. Table 8. SED fitting results for the IR-bright sources with photometry best-fitting to galaxy templates. Additionally, where we have non-detections or no coverage, we list well studied host galaxies identified in the literature in the last 7 rows. References are given below the table and provide spectroscopic redshifts in most cases, or additionally, parameter values for galaxies lacking the photometry needed for SED fitting.

GRB
z ph ot z s p e c SFR (M yr −1 ) M (10 10 M ) sSFR (10 −10 yr −1 )  We further quantify the effect of the small number statistics in this sample by bootstrap resampling of the SHOALS sample. We extract a subsample of sources from one of the reference samples matching the IR bright sample in size, and calculate its redshift distribution. This is done 100,000 times to explore the frequency with which the subsample realisation matched the observed distribution. We define a match as a scenario in which an appropriate fraction of the sample lies at z = 1.24 (our highest optically detected galaxy redshift) or lower. We also consider the Kolmogorov-Smirnov (KS) statistic for the same distribution. As the results in table 9 demonstrate, the KS-test is unable to the reject the null hypothesis that IR-bright sample is drawn from the same population as SHOALS, while the bootstrapping estimate gives a ∼0.06 per cent chance of drawing this redshift distribution from SHOALS. Assuming instead that the undetected sources are dusty, lower redshift galaxies, or that they are stars, only increases the disparity between the samples (see later).
The differences between the bootstrap and KS-test results are pronounced. The bootstrapping method supports a much stronger identification of IR-bright LGRB hosts being biased in redshift, compared to the KS-test. Fundamentally the two tests are exploring different aspects of the data. The KS-test is primarily sensitive to the shape of the distribution, Table 9. Bootstrap and KS-test results for the LGRB redshift, stellar mass, V-band attenuation and star formation rate distributions, compared to unbiased samples over the same redshift range.

Property Boostrap Target
Bootstrap KS-test %age prob. and at z < 1.24 these are similar. However, the entire distribution for the IR-bright galaxies is shifted towards lower redshifts, producing the clear bootstrap results indicated in table 9. The IR-bright hosts possess a biased distribution in redshift, but since the distribution shape is similar, this does not necessarily imply a distinct underlying population.  Figure 13. The cumulative distribution of stellar masses in the IR-bright host population, the SHOALS sample, and the dark burst sample. The IR-bright hosts appear to be more massive than even the dark burst population.

Masses, Dust Extinction and SFR
We compare the stellar masses derived for our IR-bright LGRB host population (11 objects -6 from our SED fitting and 5 literature values), to existing samples over the same redshift range, in figure 13. The IR-bright and SHOALS distributions yield a KS-test p-value that passes the 2σ threshold for significance. This, and the corresponding bootstrap result, are given in table 9. In comparing the mass distributions, we assume that the optically-faint subsample is not significantly biased in mass. Given that we are selecting in the infrared, we expect the dust extinction (and hence re-emission at long wavelengths) to be a parameter of interest. A constraint on this is obtained from the SED fitting, parameterised by the Vband attenuation A V . The cumulative distribution of A V in our LGRB sample is shown in figure 14. This is compared to the LGRB host distribution of A V (where we have restricted the sample to z < 2), as determined by Perley et al. (2013) by correcting previous optically biased studies.
The bootstrap and KS-test results are again given in table 9. 50 and 75 per cent bootstrap targets are provided to demonstrate the significance of the divergence of the distributions around A V = 1. We make no a priori assumption about the extinction in optical non-detections. As in the case of the redshift distribution, it is appropriate to consider the possibility that the 5 optically undetected sources may, in fact, be biased and have very high dust extinction values. This would only strengthen the conclusion that IR bright sources are dustier than the typical host galaxy, while also making them more extreme outliers in redshift.
Finally, we compare our population distribution in terms of the star formation rate in LGRB hosts. We compare the 7 IR-bright LGRB sources at z < 1 for which we have SFR constraints (4 from SEDs, 3 from literature values) against the z < 1 distribution reported by Salvaterra et al. (2009) andJapelj et al. (2016) for the BAT6 LGRB subsample in figure 15. The KS-test p-value and boostrap result in table 9 fail to reject the null hypothesis that the IR

Host Luminosity
Given a large range of source redshifts, any single photometric band samples a range of rest frame wavelengths, complicating an analysis of their rest frame magnitudes. Perley et al. (2016b) obtained Spitzer 3.6 µm photometry for each of their targets in the unbiased SHOALS catalogue. We transform the 3.4µm W1 and Spitzer 3.6 µm apparent magnitudes to absolute magnitudes, without K-correction, while cautioning that these correspond to a rest-frame wavelength of ∼3.5 µm/(1 + z) in each case. The absolute W1/(1 + z) magnitudes of the LGRB sample are shown in figure 16, and are compared to the Perley et al. (2016b) distribution. Since this is a direct comparison with Spitzer data, the distribution in absolute magnitude at a given redshift is indepent of the K-correction uncertainties. At low redshift (z < ∼ 0.6), the  Figure 16. The absolute, rest frame magnitude versus redshift for the IR selected hosts (red) and SHOALS. No K-correction is made to either sample. The absolute magnitude corresponding to the WISE apparent magnitude 2 σ limit is given by the dashed curve.
W1 band is probing a fall off in the stellar flux and a rise in dust and PAH emission, rather than the stellar continuum. Above z ∼ 0.6 and below z ∼ 3, the W1 magnitude is probing a fairly flat region of stellar emission and so is a good mass indicator. Again, the comparison to literature work is valid because the SHOALS magnitudes also suffer from this effect. The faintest objects detected in WISE band 1 have an apparent magnitude m W 1 ∼ 20. Placing such a source at z = 2 corresponds to a 3.4µm/(1+z) absolute magnitude of ∼-24.5, which is at the very high end for GRB hosts as figure 16 shows. The dashed line in figure 16 indicates the 2 σ detection threshold for WISE as a function of redshift. It is clear that most SHOALS host galaxies, even at low redshift, fail to satisfy this threshold.
We note that we are only matching against very bright sources, satisfying the cut-off for inclusion in the ALLWISE catalog. Perley et al. (2016b) obtained their 3.6µm host galaxy fluxes by subtracting the flux from nearby bright sources to reveal an underlying host. This introduces the possibility that we may be overestimating the IR flux. However, in most cases we do not see another optical source that might be the true origin of the WISE flux. In other words, if the IR sources we match to are simply chance alignments, then many of the GRB hosts would have to be too optically faint for detection in Pan-STARRS. In addition, we know that IR-bright GRB hosts exist from previous work (Stanway et al. 2015a,b). Therefore, we have confidence that the differences between the 3.6µm magnitudes in this sample, and those of the SHOALS sample, are real.
The contrast between our sample and SHOALS is rather unsurprising, since our IR-bright sample is selected to be extreme in the W1 band. However, it does present the possibility that the bias towards low redshifts comes from sampling deeper into the host galaxy luminosity function. While our failure to identify higher redshift hosts may be a simple WISE data selection effect, our most luminous host galaxy would, theoretically, have been detectable out to z ∼ 3. At these redshifts, the W1 band probes the rest-frame near- LGRB 080517 is indicated by the lower limit. We have calculated E iso for LGRB 050219A and LGRB 061002, based on the E p and S bol values estimated by Butler et al. (2007). infrared and is unlikely to be strongly affected by dust, but the optical is probing the rest-frame ultraviolet, so dust extinction may account for observed optical non-detections.

Burst Luminosity
GRB studies have found a correlation between the restframe photon energy E p,i at peak prompt emission, and the isotropic equivalent energy E iso , known as the Amati relation (Amati et al. 2002;Amati 2006). One member of our IR-bright host population, LGRB 080517, has already been identified as a sub-luminous LLGRB on this relation (Stanway et al. 2015a) due to its low E iso and lower limit on E p,i . We thus consider whether other members of this population are low-luminosity bursts. A number of targets in our sample have literature constraints on E iso and E p,i . Six bursts had E iso and E p,i estimates in the literature (see the references in table 8 for details). For LGRBs 061002 and 050219A, we calculate E iso and the rest frame peak energy E p,i using the SDSS photometric redshift given in section 6, the spectroscopic redshift presented by Rossi et al. (2014) and the bolometric fluence and observed peak energy provided by Butler et al. (2007). We use a cut-off power law to extrapolate from the gammaray band and infer E iso .
These eight GRBs with IR-bright hosts are placed on the Amati relation in figure 17. Additional (non-IR bright) bursts are shown as grey circles in the background for reference, while the bursts from our sample are shown in red (long bursts) and blue (short) with error bars. While LGRB 061002 lies in a similar region to LGRB 080517, towards the low-luminosity region of parameter space, the uncertainties cannot rule out consistency with the Amati relation. Targets selected in the W1 band appear to cluster at the low E p end of the distribution, consistent with expectations given the low redshift bias in our sample.

INTERPRETATION
Our targets were selected to be bright in W1 imaging. At low redshifts, W1 lies above the bulk of stellar photospheric emission, while at intermediate redshift 0.6 < z < 3 it probes the bulk of the stellar mass. As figure 16 makes clear, we are strongly limited by the shallow depth of W1 imaging in WISE. This introduces a Malmquist bias (Malmquist 1922) which limits us to comparatively low redshift. Given that the star forming galaxy luminosity function favours less luminous systems, we are able to sample a larger fraction of the extant population at the lowest redshifts, and likely miss many of the galaxies in the more distant Universe.
However, many of the sources for which optical data is available are also extremely red (see appendix and figures 10 and 11). Red colours can arise from either old stellar populations or dusty populations. The Balmer break, which appears in a galaxy spectrum after the death of hot young O stars, is an age indicator and so may allow us to distinguish betwen these two cases. In 2 of the 14 galaxy SEDs presented in figure 10, there is no clear constraint on the Balmer break. In 5 cases there is a clear indication of a Balmer break. Even those galaxies with Balmer breaks have moderate star formation rates. Thus for the sources which are detected in the optical, we can rule out IR-bright, old galaxies as being dominant in this sample.
The population for which an IR-detection exists, but which remain undetected in optical surveys presents a slightly different challenge. It is possible that some of these may show a prominent 4000Å break, if the galaxy in question lies at z>0.8. However, given the attenuation distribution presented in figure 14, our failure to detect galaxies in the range 1 < z < 3 instead suggests that the rest-frame UV being probed at these redshifts has been attenuated by dust, resulting in observer-frame optical non-detections.
In figure 18 we present the multivariate distribution of properties for our LGRB hosts. While our W1 band magnitudes are bright by construction, the emission may arise for different reasons. The sources could be extremely dusty, or, if lying at 0.6 < z < 3, may instead be very massive. Figure 18 indicates that our sources are amongst the most massive LGRB hosts known at their redshifts (figure 18, left). In terms of attenuation and mass, mass appears to be the dominant factor (figure 18, centre). Around half of the IR-bright sample are dusty for their star formation rate (figure 18, right). In turn, this implies they have high mass for their SFR (Whitaker et al. 2017), or relatively low sSFRs for LGRB hosts (i.e. sit below what is known as the galaxy 'main sequence ', Noeske et al. 2007). This is also seen in the sSFRs given in table 8. While some hosts have high sSFRs of > 10 −9 yr −1 , as is expected from LGRB hosts generally, others have low values (e.g. 6×10 −12 yr −1 for LGRB 050219A and 3.8×10 −11 yr −1 for LGRB 080517).
This has implications for the supply and depletion of molecular gas in these systems. If our sample is preferentially selecting molecular gas (and therefore star formation) poor systems, this implies either that the hydrogen gas in these systems is in a different phase (i.e. heated by shocks, interactions or a large scale environment effect) or that there is very little gas, as might be seen in an old underlying system which accretes a small star forming satellite. This latter scenario has already been demonstrated in the most luminous example in this sample of IR-bright systems, the host of LGRB 080517 (Stanway et al. 2015b, and in prep). However a molecular gas detection has recently been reported for a second GRB host in our IR-bright selection, the dark burst LGRB 080207 (Arabsalmani et al. 2017), and two galaxies in our LGRB sample (080517 and 100316D) have detections of IR molecular Hydrogen emission lines (Wiersema et al., in prep). Arabsalmani et al. (2017) find that the host of LGRB 080207 appears to follow the normal scaling relations for star forming galaxies and is not noticeably poor in molecular gas. It is however a dark burst, in a host which has already been noted for being relatively massive and dustycriteria which bring it into our overall sample. This implies that we are not preferentially selecting molecular gas rich systems. Figure 18 shows that the IR-bright LGRB host sample may not have a single set of consistent properties. If we divide the sample into dusty (A V > 1) vs non-dusty, and massive (log 10 (M /M ) > 9) vs less massive, we find that of 11 LGRB host galaxies in our sample that have M and A V information, 6 would be classed as both massive and dusty, while 4 are simply massive, and one is on the border of being low in mass and not dusty.
A question thus arises: is this a distinct population, or is it sampling the tails of existing and known populations? The diversity of host galaxy properties in our sample favours the latter possibility. There is no clear separation between the IR luminous host galaxies and the bulk of the LGRB population in any projection of parameter space. While our sources are dustier than the bulk of the population over the same redshift range, they do not lie significantly outside of the range of typical LGRB host properties, with both A V and another extreme property (e.g. high stellar mass, low redshift) required for selection. This suggests that we are sampling the tails of previously established host distributions. As such, it is unlikely that implications for LGRB progenitors can be derived from this sample. None the less, selection in these relatively shallow infrared bands is useful for constraining the tails of existing distributions and therefore the full range of LGRB host properties. This analysis has identified 2 new candidate GRB hosts. These are LGRB 080517, which has already been extensively studied (Stanway et al. 2015a,b) following selection with this methodology, and LGRB 061002. This burst has a false alarm probability for association with the IR-bright galaxy of only 0.0003. We have also newly identified existing host candidates as being bright in the WISE bands. Generally, the low number of IR-bright GRB hosts we have found (only 20 from 1001 bursts), is consistent with expectations for the following reasons. First and foremost, a shallow survey such as WISE only allows us to sample low redshifts, covering a small cosmological volume. A secondary effect might be arising because the LGRB rate appears to be suppressed in the local Universe (e.g. Perley et al. 2016a), tracing the star formation rate density but also the increase in metallicity over cosmic time. Local and/or massive and dusty galaxies are typically detected in WISE, however massive and dusty galaxies are typically higher in metallicity and are less of- ten seen to host LGRBs. SGRBs could be hosted by these classes of galaxy, however short bursts are intrinsically less luminous and have a much lower observed redshift distribution. Finally, by the time we get to high enough redshifts that the rest frame bands LGRB hosts are bright in (e.g. UV) are redshifted to W1, the emission is too faint for detection by shallow surveys such as WISE. These factors conspire to produce only small sample of GRB hosts that can be detected in WISE.
In theory, using these data, it is possible to calculate the fraction of low redshift and massive, dusty GRB hosts that have been missed in previous follow up studies. This is primarily true for low redshift galaxies. As figure 16 shows, the WISE 2 σ limit corresponds to absolute magnitudes of ∼-18 at z=0.1. This is well below the knee of the galaxy luminosity function, probing galaxies of ≈0.1 L * and above (Babbedge et al. 2006). We would therefore naively expect to see the majority of galaxies at z < 0.1 in W ISE. However, our cross matching procedure fails to capture bursts at very large apparent offsets from their hosts, and this primarily affects short bursts and the lowest redshift galaxies in the GRB population. In fact, of the three known long bursts at z < 0.05, we identify zero as IR-bright. There are four long bursts known to be at 0.05 < z < 0.1, of these we detect two and one of those was selected as low z purely on the basis of its WISE detection. The remaining two GRBs constitute one source at large offset from a very extended galaxy, and one source in a highly sub-luminous host. This demonstrates that 3.6 micron selection is not a highly efficient method for identifying low redshift bursts in all cases, but nonetheless that for typical galaxies in its redshift sweet spot, we would expect to have identified any source present in the archival data. Thus, while we cannot rule out the continued presence of z < 0.05 galaxies in the archival data set, we can be reasonably confident that the fraction of GRB hosts without redshift identifications which lie at 0.05 < z < 0.10 is very low. With the exception of GRB 080517, which was selected as part of this sample, we identify no further candidates in this range and would have expected to detect at least half of those present, unless the GRB host luminosity function has a steep faint end slope. This suggests that the total number of 0.05 < z < 0.1 hosts still remaining unidentified in the archival sample is of order a few, at most.
We note that deeper IR data, for example that obtained with Spitzer or JWST, would detect galaxies further down the luminosity function, but would not overcome the issues of large projected offsets or galaxy sizes. It is therefore possible that future deeper infrared surveys will allow this analysis to be extended to slightly higher redshifts, and, for example, characterise the overlooked GRB hosts (if any) at 0.1 < z < 0.3.

CONCLUSIONS
A population of infrared-bright GRB host galaxies has been identified by cross-matching X-ray afterglow positions to the ALLWISE catalogue. Selective cuts in apparent magnitude and catalogue quality flags, in addition to a false alarm probability analysis, yield 55 IR-sources that are convincingly associated with a GRB X-ray position. Compiling photometry from surveys and our own observations, we perform SED fitting, finding that 14 sources fit well to galaxy templates and 14 are best-fittung to stellar SEDs. The remainder are either cut from the sample, optically undetected, lacking survey coverage or have been previously studied. Spectroscopy of 6 targets supplements the photometric data and allows us to rule out 3 stellar interlopers. Our methodology has identified the host of LGRB 080517, and potentially LGRB 061002. The former has already undergone extensive study, while the candidate host of LGRB 061002 is newly reported here. Focusing on LGRB hosts in particular, we find that the population is biased towards massive, dusty and low redshift galaxies, with respect to unbiased samples of LGRB hosts. The low redshift bias appears to be due to the depth of WISE, and within this low redshift population, galaxies with high stellar mass and dust content are selected. Dusty and local galaxies are the most frequently discovered hosts of dark and low-luminosity LGRBs respectively, classes of burst which are crucial to understand if we are to develop a full picture of the collapsar GRB phenomenon. We propose that the spatial association of IR-bright galaxies with LGRBs is therefore a useful technique for the identification of unusual host systems. Table A1. The full table of archival photometry for 30 GRB locations with reliable WISE detected counterparts, where optical data was available. No Galactic dust correction has been applied. If a magnitude is entered without an associated error, it is a 2 σ upper limit. If no value is entered, then photometry in that band was unavailable. GRB 150120A has two sets of photometry, one for each component (see section 6.1).   figure A1 we provide archival imaging postage stamp cutouts in the W1, J and r bands. Overlaid are the X-ray R 90 error circles (blue, varying size) and the WISE source centroids (red, fixed diameter). The images are 30 arcsec on each side and are centered on the GRB location. The figure extends over seven panels.