The sensitivity of Imaging Atmospheric Cherenkov Telescopes (IACTs) used to carry out Very-High-Energy (VHE; E>100 GeV) gamma-ray astrophysics strongly depends on the ability to reject cosmic-ray (hadron) background events in favor of gamma rays. Since cosmic-ray initiated Extensive Air Showers (EAS) dominate those initiated by gamma rays by several orders of magnitude, the ability to accurately distinguish between gamma-ray or hadron-initiated showers is a long-standing problem within the IACT community. Motivated by the physical differences in gamma-ray and hadron EAS, some existing work in this őeld has focused on implementing deep learning techniques to solve this classiőcation problem. The predominant deep learning approach has been to train models in a supervised fashion on simulated EAS data, which has encountered issues when transitioning from simulation training data to real EAS data. We take a novel deep learning approach focused on unsupervised learning with real data from the VERITAS IACT to learn spatial relations and temporal correlations of the EAS. We implemented a Two-Dimensional Convolutional Long-Short Term Memory Autoencoder network (2DConvLSTM-AE network) given its strong performance in both spatial-and time-relationed data. The autoencoder architecture enables us to encode a latent space mapping of the generalized features for a downstream classiőcation. We őnd that while the 2DConvLSTM-AE is capable of producing faithful reconstruction of EAS, the ability to differentiate EAS by their origin particle has not yet been demonstrated but provides a promising avenue for future research.

38th International Cosmic Ray Conference (ICRC2023) 26 July -3 August, 2023 Nagoya, Japan

Very-high-energy (VHE) gamma rays are photons with energies greater than 100 GeV, forming the highest energy range of the electromagnetic spectrum. VHE gamma rays are created by nonthermal processes, usually from violent events such as supernova and jets of material emanating from supermassive black holes [1]. Increasing the availability of VHE gamma-ray data will enable VHE astrophysics researchers to further understand these extreme events of our universe. The primary instrumentation used to capture VHE gamma ray events are Imaging Atmospheric Cherenkov Telescopes (IACTs). When a primary gamma ray or a hadron collides with atmospheric particles, it initiates a cascade of particles and Cherenkov photons, known as Extensive Air Showers (EAS). IACTs capture the images of these EAS, thereby indirectly observing the incoming particles. There are key physical differences in the EAS created by gamma rays and hadrons; EAS caused by gamma rays have a narrower particle spread through their development and primarily consist of positrons, electrons, and lower energy gamma rays, whereas EAS caused by hadronic particles have a wider particle spread and contain muons formed from the decay of charged pions and kaons [2]. The duration from an EAS can range from ∼ 10 -100 nanoseconds, and its developments are imaged on the IACT cameras as an ellipse. The result of these different particle compositions gives a gamma-ray induced event a visibly 'tighter' shower as compared to the wider shower produced by hadron induced events. These differences in EAS properties stemming from the different causal particles makes it possible to determine which particle caused a shower via analysis of Cherenkov radiation captured by IACTs. It is estimated that ∼ 99.9% of the EAS captured by IACTs are originated by hadrons as opposed to gamma rays [2]. The sensitivity of IACTs used to carry out VHE gamma-ray astrophysics strongly depends on the ability to reject hadron background events in favor of gamma rays.

The traditional approach to the gamma/hadron classiőcation problem is the implementation of the Hillas Parameters, established by [3]. This work analyzes images of EAS simulations and deőnes a set of key properties to describe the events, such as: length and width of the image ellipse, total charge integrated over the full shower development, nominal distance, azimuthal angle, and orientation angle. Figure 1 shows a mapping of the length (left) and width (right) parameters for gamma rays (signal) and hadrons (background) [4]. While the Hillas parameters have reasonable discriminating power, it is clear from Figure 1 that improvements in class separation could be made. One relevant limitation of these parameters

Grant Oie is that the calculation method is being performed on a time-integrated form of the EAS temporal evolution, losing information about the timing of EAS development and evolution which may be useful in increasing the classiőcation power.

Advancements in the technological and theoretical foundations of machine learning have motivated the applications of deep learning techniques to the gamma/hadron classiőcation problem. Given that the Hillas parameters of the EAS leave room for improvement in class separation, deep learning approaches have been implemented in hopes of discovering additional parameters by which classiőcation becomes more effective. Prior machine learning work has included tools such as Boosted Decision Trees (BDTs) [2,4], Convolutional Neural Networks (CNNs) [2,5] and Long-Short-Term Memory (LSTM) networks [6]. The predominant deep learning approach has been to train models in a supervised fashion on simulated EAS data, which has encountered issues when transitioning from simulation training data to real EAS data. We believe that pursuing an unsupervised learning approach has the potential to solve this issue and provide the deep learning models with more information from which to learn relevant parameters for categorization. To accomplish this, we are implementing a Two-Dimensional Convolutional Long-Short Term Memory Autoencoder network (2DConvLSTM-AE network) given LSTM networks' proven strong performance in both spatial-and time-relationed data, with the autoencoder structure enabling unsupervised learning [6].

The dataset includes 23,994 cubic-interpolated images of EAS captured by VERITAS, primarily from the Crab Nebula source. The Crab Nebula is being used because it is the strongest VHE gamma-ray emitter, frequently used as a standard calibration source. The original data captured by each VERITAS telescope is in the form of a 499-pixel hexagonal image, with each pixel corresponding to a photo-multiplier-tube. We then take the 499-pixel hexagonal images and apply cubic interpolation to convert each frame to a 96×96 image, one for each of the four telescopes. The data is otherwise un-normalized, where each pixel value corresponse to the acquired raw charge value. When the VERITAS telescope array is triggered, it captures a set of 16 images in a time series. Each time frame captures 2 nanoseconds of sensor data, so when collated, each image cube covers 32 nanoseconds (16 frames x 2 nanoseconds). Preserving the temporal aspect of the EAS introduces new information to the deep learning model that has the potential to improve its performance. The resulting EAS data is formatted in an image 'hypercube' with dimensions of 4x16x96x96: 4 telescopes, 16 time frames, 96×96 pixel images.

Our selection of the 2DConvLSTM-AE framework was inspired by prior work with a supervised 2DConvLSTM architecture by [6], which evaluated the effectiveness of a 2DConvLSTM on parameters extracted from simulated hadron and gamma-ray events. Similar to other supervised learning approaches dealing with gamma/hadron separation, the 2DConvLSTM performed very well on simulated data but faltered when introduced to real data from VERITAS [6]. Given this observed issue common to many gamma/hadron approaches, we decided to take an unsupervised, data-driven approach that emphasized using real events from VERITAS while making as few modiőcations to the data as possible. To accommodate an unsupervised approach, the autoencoder architecture enables us to encode a latent space mapping of the generalized features for downstream classiőcation. The 2DConvLSTM-AE is trained on the time series images of shape 4×16×96×96 described above, applying convolutions to compress the sample down to a latent vector with 512 dimensions. The latent vector is extracted for downstream latent space analysis which is described fully in later sections. The decoder then takes the latent vector and applies transposed convolutions to reconstruct a set of time series images of shape 4×16×96×96. The 2DConvLSTM-AE is trained to minimize the mean square error loss between the original and reconstructed cube; once the model has been trained sufficiently to convergence, the individual MSE loss values were computed for each data sample and it is considered as an anomaly score. Conceptually, owing to a relative low incidence of gamma rays with respect to hadronic signals in the data, it is to be expected that the events that have a high anomaly score will correspond to poorly generalized ones by the machine model, i.e., the gamma ray events. The ability to encode a latent space is at the crux of our unsupervised deep learning approach, given that the latent space can potentially provide more information about the individual image, rather than just a binary separation into two classes. By compressing each sample down to a lower dimensional representation, there are a number of incremental inference pathways that become available by plotting this representation with a Principal Component Analysis (PCA), Uniform Manifold Approximation and Projection (UMAP) [7], or a combination of PCA and UMAP. With a visualized low dimensional representation of the feature space, sub-classiőcations of event types both between and within groups such as hadrons, gamma rays, and muons may become evident. Additionally, observing the spatial relationships between the classiőcation groups and overlaying the latent space with a number of relevant parameters such as anomaly score or Hillas parameters may provide additional insight relevant to classiőcation.

2DConvLSTM-AE on Very-High-Energy Gamma-Ray Classiőcation Grant Oie

We assessed the reconstruction performance of the 2DConvLSTM-AE and achieved a satisfactory level of visual ődelity between the original data and the reconstructions, demonstrated in Figure 4, with the reconstructions closely resembling the original images when there is a strong shower image present in the frame. One signiőcant issue that persists in the reconstructions is a 'phantom' shower image pattern which appears to be remnants of the strongest shower image in the frame spilling into the other telescopes, as demonstrated in Figure 5.

We tested two dimension reduction methods when analyzing the 512-dimension latent space: a traditional UMAP with a latent space input, and a UMAP with a PCA-reduced latent space input (Figures 6, 7, 11). Given the high dimensionality and variance in the latent vector, applying two levels of dimension reduction with the PCA followed by a UMAP enabled a greater degree of organization among anomalous and normal events. In Figure 6 and 7 the events in the UMAP are colorized by the log of their anomaly score to better visualize where the anomalies exist in the latent space, with darker points corresponding to higher anomaly scores.

The highly anomalous events cluster around the central region of the 3-dimensional latent space mapping that form a loop of dark points in Figure 6. In order to analyze the regions composed of high anomaly score events, a 2-dimensional side perspective was generated in Figure 7. The following analysis examines the reconstructions from some the clusters of interest as highlighted in Figure 7. The clusters analyzed below consist of 50-200 events per cluster. The events shown below were randomly selected from the clusters. More extensive analysis of each clusters' event composition is required to draw concrete conclusions about their signiőcance.

Cluster 1 & Cluster 2 (Figure 8): Groupings of very high energy hadron-initiated events that triggered the low-gain mode on the VERITAS telescopes. Potential presence of muons.

2DConvLSTM-AE on Very-High-Energy Gamma-Ray Classiőcation PoS(ICRC2023)692 2DConvLSTM-AE on Very-High-Energy Gamma-Ray Classiőcation Grant Oie Cluster 4 (Figure 10): Potential clustering of very strong hadron-initiated EAS, with many of the events having a brighter background. As such, we highlight that the 2DConvLSTM-AE has learned meaningful representations in the data and the different clusters correspond to physical parameters in the dataset.

Prior work utilising the Boosted Decision Tree (BDT) machine learning algorithm on Hillas parameter data to perform gamma/hadron classiőcation has shown the ability to produce reliable values that indicate the likelihood of an event being a gamma-ray or hadron [8]. BDT scores close to 1 indicate a higher likelihood of a gamma-ray event, while values closer to -1 indicate hadron events. We expected that the high BDT score events would cluster in one portion of the latent space, but instead they are more homogeneously distributed. This is a key avenue we will explore in future work.

The 2DConvLSTM-AE architecture demonstrated the ability to generate relatively strong reconstructions of EAS images, with the exception of the 'phantom' shower image pattern examined above. The reduced-dimension mapping of the latent space partially indicates the ability to organize the events by high-and low-anomaly, as well as BDT scores. However, the high anomaly scores and high BDT scores composed inverse regions of the UMAP which is contrary to our expectation that the 2DConvLSTM-AE will assign higher anomaly score values to gamma ray events. By examining speciőc anomalous clusters in the latent space, visual inspection suggests that the model is likely producing higher anomaly scores for events with higher energies. The resulting sentiment is that a 2DConvLSTM-AE deep learning architecture coupled with unsupervised learning on real data from the VERITAS IACT offers a promising avenue for future research. Further experimentation could include deeper exploration of the latent space spatial composition, overlaying the latent space with a number of relevant metrics, including Hillas parameter values, applying energy range cuts to the data, as well as implementing a Variational Autoencoder and further őne-tuning of the 2DConvLSTM-AE parameters.