**1. Introduction**

The bimodal duration distribution of Gamma-ray Bursts (GRBs) suggests the separation of GRBs at T<sup>90</sup> ≈ 2 s into short/hard and long/soft classes [1]. The association of long GRBs with star forming galaxies [2] and Type Ic supernovae (Galama et al. [3], Woosley and Bloom [4]; for a review, see Cano et al. [5]) provides an observational link between long GRBs and the deaths of massive stars, supporting the collapsar scenario [6]. There is substantial evidence to support compact object mergers (neutron star–neutron star or neutron star–black hole) as the progenitors of short GRBs [7,8]. The location offsets of short GRBs from their host galaxies [9,10], their proximity to elliptical galaxies [11], and the association of GRB 170817A, an unusual short GRB, with the neutron star merger event GW 170817 detected by aLIGO [12–14], all support the merger hypothesis for the origin of short GRBs.

Other formation scenarios for short GRBs include the accretion-induced collapse of a white dwarf, double white dwarf mergers, or neutron star–white dwarf mergers [15–17], possibly leading to an unstable magnetar remnant. There are notable exceptions to the shortmerger/long-collapsar paradigm, such as the short-collapsar event GRB 200826A [18–20], and GRB 060614, a long GRB without a supernova [21]. It has been suggested that many of the short duration GRBs of high redshift arise from collapsars [22]. Consideration of additional GRB characteristics, such as late X-ray flares in some short GRBs, and the non-detection of a supernova associated with some long GRBs [23], led to the suggestion of a new classification scheme [21], with Type I (massive star origin) and Type II (compact object merger origin) GRBs defined by many multiple observational criteria beyond the traditional duration and hardness [22,24]. Lü et al. [25] suggested a new parameter *e*,

**Citation:** Salmon, L.; Hanlon, L.; Martin-Carrillo, A. Two Dimensional Clustering of *Swift*/BAT and *Fermi*/GBM Gamma-ray Bursts. *Galaxies* **2022**, *10*, 77. https:// doi.org/10.3390/galaxies10040077

Academic Editors: Elena Moretti and Francesco Longo

Received: 27 April 2022 Accepted: 22 June 2022 Published: 25 June 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

based on the isotropic equivalent energy and peak energy, to classify bursts. Additionally, Donaghy et al. [26] considered 10 observational criteria for HETE-2 bursts, concluding that the best criteria to classify GRBs as 'short population' or 'long population' bursts are host galaxy properties, spectral lag, and the presence of a long-soft bump or gravitational waves.

In view of the diversity in GRB phenomenology, a definitive classification of GRBs based on duration alone is challenging. Several studies have found evidence for an additional 'intermediate' duration class of GRBs, first identified through Gaussian fits to the duration distribution of GRBs in the Third BATSE catalogue [27] and, subsequently, in fits to the GRB duration distributions of BeppoSAX [28], RHESSI [29], and *Swift* [30–35]. This class appears as an additional Gaussian 'component' required for the best-fit solution. However, the observed duration distribution can be recovered by modelling it as two skewed distributions [31,36,37], without requiring a third component.

GRB catalogues provide a set of standard parameters measured for each GRB, including duration (T90), hardness ratio (HR), fluence (S), peak flux (PF), peak energy (Epeak), and spectral fit parameters, including the low and high energy spectral indices of the Band function [38], which fits the keV-MeV GRB spectrum, typically denoted in the literature as *α* and *β*, respectively. In the case of *Fermi*/GBM, the catalogue contains over 300 parameters for each GRB [39,40]. The availability of such large GRB catalogues allows the application of bivariate and multidimensional analyses to the data.

Table 1 summarises the previous studies, along with the resulting number of components identified for different GRB datasets. Between two and five classes of GRBs are found, depending on the sample, parameters, and methods used. Clustering of the durationhardness plane of the final BATSE GRB catalogue identified three [41–43] or five [44–46] classes of GRBs separated by their duration, fluence, and hardness. Unsupervised neural network analysis also revealed an intermediate class [47] or two classes [48,49]. However, only two classes were found in the BATSE sample using self-organising maps [50] and fits to the duration-hardness plane with skewed bivariate distributions [51,52].

The clustering of the duration and hardness of *Swift*/BAT GRBs [53,54] and the clustering of light curve shape indicators [55] identified three classes of bursts. Gaussian Mixture Modelbased (GMM) clustering applied to the *Fermi*/GBM sample revealed that GRB 170817A fit within the intermediate class in the duration-hardness plane [56], and that five classes could be identified by clustering spectral fit parameters, fluences, and durations [57]. Principal Component Analysis (PCA) also identified three classes in *Fermi*/GBM [58] and BATSE [59] samples.

**Table 1.** Methods and resulting components identified in clustering, fitting, and dimensionality reduction techniques applied to GRB populations. HR denotes Hardness Ratio, PF denotes Peak Flux, and S represents fluence. Studies which consider intrinsic properties such as redshift-corrected duration and hardness are marked with a \*.



**Table 1.** *Cont.*

Observational bias has been suggested as a possible origin of the putative intermediate class. Bias caused by short temporal trigger windows favours short low-fluence bursts (fluence-duration bias; Hakkila et al. [49]), while the low signal-to-noise ratios of long faint bursts can cause them to be mistaken for short bursts ('tip-of-the-iceberg' effect; Lü et al. [71]). However, neither of these effects have been able to reproduce the third class in simulations. It has been shown that the third class can arise as a consequence of fitting symmetrical models to the GRB duration distribution, which may be skewed rather than symmetrical [31,36,37,51], possibly as a result of the GRB pulse shapes [72].

The significant number of GRBs with measured redshift in the *Swift* and *Fermi* samples has allowed studies of intrinsic properties, which have pointed to the existence of two classes in the *Fermi*/GBM sample [32]. For the *Swift*/BAT sample of bursts, one [30,34], two [35,61,63,65], or three [32,54] classes have been identified. However, cosmological time dilation applied to GRB durations has not been found to transform a rest-frame twocomponent Gaussian duration distribution to the observed skewed one [73]. While there are now more than 400 *Swift* GRBs with measured redshift, there are only 25 short duration bursts with T90,obs < 2 s. The rest-frame studies outlined in Table 1 note that the short duration sample is not statistically significant, and a larger sample is required [54,65].

This paper reports on an updated two-dimensional clustering analysis in the durationhardness plane of the large *Fermi*/GBM and *Swift*/BAT GRB samples. Advancing previous studies, the analysis presented here makes use of an entropy criterion to identify 'excess' components that may be identified in the standard GMM clustering of data but which arise from the application of Gaussian models to non-Gaussian underlying data [74]. This method has been applied in other astrophysical contexts, for example in the clustering of stars [75]. As the number of short GRBs with redshift has not grown significantly since previous studies, this paper focuses on GMM clustering using observed, rather than intrinsic, properties.

Section 2 outlines the sample construction, while Section 3 provides details of the methods applied to perform clustering. The results and discussion are presented in Sections 4 and 5 respectively, while the conclusions are outlined in Section 6.
