A Deep Neural Network Based Reverse Radio Spectrogram Search Algorithm
Modern radio astronomy instruments generate vast amounts of data, and the
increasingly challenging radio frequency interference (RFI) environment
necessitates ever-more sophisticated RFI rejection algorithms. The "needle in a
haystack" nature of searches for transients and technosignatures requires us to
develop methods that can determine whether a signal of interest has unique
properties, or is a part of some larger set of pernicious RFI. In the past,
this vetting has required onerous manual inspection of very large numbers of
signals. In this paper we present a fast and modular deep learning algorithm to
search for lookalike signals of interest in radio spectrogram data. First, we
trained a B-Variational Autoencoder on signals returned by an energy detection
algorithm. We then adapted a positional embedding layer from classical
Transformer architecture to a embed additional metadata, which we demonstrate
using a frequency-based embedding. Next we used the encoder component of the
B-Variational Autoencoder to extract features from small (~ 715,Hz, with a
resolution of 2.79Hz per frequency bin) windows in the radio spectrogram. We
used our algorithm to conduct a search for a given query (encoded signal of
interest) on a set of signals (encoded features of searched items) to produce
the top candidates with similar features. We successfully demonstrate that the
algorithm retrieves signals with similar appearance, given only the original
radio spectrogram data. This algorithm can be used to improve the efficiency of
vetting signals of interest in technosignature searches, but could also be
applied to a wider variety of searches for "lookalike" signals in large
astronomical datasets.
Peter Xiangyuan Ma, Steve Croft, Chris Lintott, Andrew P. V. Siemion
# A Deep Neural Network Based Reverse Radio Spectrogram Search Algorithm

###### Abstract

We developed a fast and modular deep learning algorithm to search for lookalike signals of interest in radio spectrogram data. First, we trained an autoencoder on filtered data returned by an energy detection algorithm. We then adapted a positional embedding layer from classical Transformer architecture to a frequency-based embedding. Next we used the encoder component of the autoencoder to extract features from small (\(\sim 715\,\mathrm{Hz}\) with a resolution of \(2.79\,\mathrm{Hz}\) per frequency bin) windows in the radio spectrogram. We used our algorithm to conduct a search for a given query (encoded signal of interest) on a set of signals (encoded features of searched items) to produce the top candidates with similar features. We successfully demonstrate that the algorithm retrieves signals with similar appearance, given only the original radio spectrogram data.

keywords: Deep Learning - Machine Learning - Signal Processing - Radio Astronomy - Technosignatures

## 1 Introduction

### Interference Rejection in Radio Astronomy

The rejection of radio frequency interference (RFI) is a perennial challenge for radio astronomy, particularly given the increase in satellite constellations that transmit at a range of radio frequencies and are detectable even at remote observing sites. RFI rejection is traditionally performed using statistical or machine learning techniques that risk rejecting a potential signal of interest (SOI) in the process of flagging RFI signals (Pinchuk and Margot, 2022). This is a particular problem when searching for astrophysical transients or for technosignatures, since the signal morphology may not be known ahead of time and it may be rejected before a human has the opportunity to review it. Novel deep learning algorithms have been successfully used to address these kinds of problems (Ma et al., 2023).

An alternative approach, which is relatively unexplored, is to start with the SOI, search through some database of signals (including RFI), and if the signals match, to then decide whether to keep or discard the SOI. This is of particular interest in technosignatures searches where the vetting of candidate signals (e.g. BLC1; Sheikh et al., 2021) is rather onerous using existing methods. The ability to find additional examples of a given SOI at different times or different parts of the observed band can be invaluable in determining the nature of the signal. The focus of this paper is to develop an algorithm that allows us to reverse search a given SOI and return similar signals. We take inspiration from well known reverse image search algorithms.

### Classical Reverse Image Search Algorithms

Reverse image searches are commonly employed in search engines and social media. These algorithms take an input image as a search term and return the locations of the same or similar images from across the web. These sites use algorithms like Scale-Invariant Feature Transforms (SIFT; Lowe, 1999), Maximally Stable Extremal Regions (MSER; Matas et al., 2004) or Bag of Words/Vocabulary Tree (Csurka et al., 2004) to power their search. Here we will briefly review some of these approaches and extend them to our work.

#### 1.2.1 Scale-Invariant Feature Transforms

A Scale-Invariant Feature Transform (SIFT) is an algorithm that takes an image, locates special "locally distinct points", and describes the features around such points using measurements such as the gradient in intensity or maximum intensity, etc. (Lowe, 1999). This process produces pairs of vectors: the coordinates of the locally distinct point, and the corresponding descriptor vector that contains the features of that point.

The locally distinct points are produced from a process of difference of Gaussians performed at varying resolutions of the images. To compute the descriptor vector we measure the gradients in pixel intensity over a region about the locally distinct points. Together this builds a feature extractor. Finally, to perform a search, we match the key points between images and check the similarity between the descriptor vectors. More concretely we can use a technique called a Vocabulary Tree, described in section 1.2.3.
#### 1.2.2 Maximally stable extremal regions

The maximally stable extremal regions (MSER) algorithm, like SIFT, attempts to find key points in the image. Specifically, it looks for objects called "blobs", defined as areas of an image that have connected elements, contrasting backgrounds, and close-to-uniform pixel intensities (Matas et al., 2004).

MSER works by taking various thresholds in the range (0, 255) and blacking or whitening out pixels depending on this threshold. The blacking and whitening of pixels create the "blobs". We perform a connected component analysis and track how the blobs evolve as we adjust the threshold. MSER helps pull out small regions of images that contain distinctive boundaries which we can use to label distinct points and compute descriptor vectors for each blob. These can then be used to compare and match.

#### 1.2.3 Bag of Words / Vocabulary Tree

The Vocabulary Tree approach in computer vision effectively implements a Bag of Words (BoW) model to compare images to each other (Csurka et al., 2004). We first run a local feature extractor like SIFT (section 1.2.1). We then construct a "codebook" containing codewords that describe several similar patches identified by SIFT. One simple means of determining codewords is to use K-Means clustering to produce centroids, which we denote as codewords. We can then compare and contrast images by comparing the corresponding codebooks using a variety of algorithms. However, a key limitation in this approach is that SIFT ignores the spatial relationships between local patches, which we know is incredibly important in describing an image accurately. This leaves room for improvement.

### Deep Learning for Computer Vision

Deep learning overcomes some of the shortcomings of the aforementioned algorithms and has been proven to effectively solve a wide range of image problems with the advent of Convolutional Neural Networks (CNN; LeCun et al., 1999). CNN's are simple in that they are a traditional neural network with the addition of convolutional layers. These layers operate by performing convolution operations between input data and a kernel. These operations have built-in inductive bias. For example, these operations are equivalent to translations that the algorithm exploits as fundamental properties in images. More concretely, this is because in images, small translations often do not change the prediction outcomes (Goodfellow et al., 2016). For example, if we take a picture of a dog and move the dog 5 pixels to the right, the image is still that of a dog. By using a convolutional layer we are baking this assumption into the model without getting the model to learn it from scratch. More formally this allows us to restrict the priors on the weights thus reducing the size of the model and allowing vast improvements in both performance and scale. More specifically CNN's have been proven to be efficient feature extractors of both local and global features which address the shortcomings described in Section 1.2.3. Due to the exceptional performance of CNN's on a wide range of image tasks, they have been widely adopted in industry for a variety of computer vision challenges (Krizhevsky et al., 2017; Koul et al., 2019) including for reverse image search algorithms (Singh and Gowdar, 2021).

#### 1.3.1 Autoencoders

Autoencoders are special cases of CNN's, employing a symmetrical CNN architecture that takes in an input and outputs a reconstruction of the original input. However, these models have a bottleneck that constricts the flow of data (Baldi, 2012). This bottleneck divides the network into two: an encoder and a decoder. The encoder takes an input image and compresses the data through the bottleneck and the decoder attempts to reverse that process. Intuitively this builds an efficient feature extractor since the encoder is attempting to reduce the input into the most important features such that the decoder can reconstruct the original input (Goodfellow et al., 2016). This kind of automatic feature extraction technique will serve as the core of our approach.

### Problem Statement

Despite these algorithms being industry-standard, out-of-the-box solutions fail to satisfy our needs. Firstly, these traditional (not deep learning) techniques disregard large-scale spatial features, as described in section 1.2.3. Traditional algorithms fail to capture complex structures within the data. This is partially the reason why traditional computer vision techniques are being superseded by deep learning in tackling computer vision problems. Secondly, under certain circumstances, we want to bias our model towards particular features, for example, frequency. Currently, there is no encoding of information on the radio frequency of the signal in any of these methods.

During our search process, sometimes we wish to match signals that not only appear similar to each other, but also are located in similar regions of the frequency band. On the surface, this appears to be a trivial problem: simply filter out the candidates after the search is completed. However, there exists a trade-off between visually similar signals and signals that are just close by in frequency. Thus if one wants to search for both similar appearing signals and signals close by in frequency we need to "bake in" the information about the frequency into the search process rather than apply it as a secondary filter. Classical search algorithms do not have an out-of-the-box solution to embedding complex features such as signal frequency.

The central questions we want to answer are: _Given a signal in a radio spectrogram, can we find all lookalike signals? Can we also find similar signals in similar regions of the frequency band? Can we make these approaches modular?_

## 2 Methods

We begin by outlining the structure of our approach. First, we build a feature extractor using an autoencoder. To do so we describe our data source in section 2.1, the cleaning of the data in section 2.2, and data preprocessing in section 2.3. We then move to build the model in 2.4. We train and test the model in section 2.5. We extend the model capabilities using a _modular_ frequency embedding strategy discussed in 2.6, and finally assemble the algorithm in section 2.7.

### Data Source

The training, validation, and testing dataset are derived from real observations from Breakthrough Listen's Green Bank Telescope \(1.8-2.9\,\mathrm{GHz}\) dataset, sourced from several different observational campaigns (Enriquez et al., 2017; Price et al., 2018; Price et al., 2020; Price et al., 2021). We used the high-frequency resolution data product where each frequency bin is \(2.79\,\mathrm{Hz}\) wide (Lebofsky et al., 2019), since our goal is to tackle data that contains RFI, and RFI can have
fine resolution in the frequency domain. The data is open source and are available from the Breakthrough Listen Open Data Archive1.

Footnote 1: [http://seti.berkeley.edu/opendata](http://seti.berkeley.edu/opendata)

### Data Filtering Energy Detection

Since we are working with real data, we need to balance the training features. Large regions of the spectrograms consist of Gaussian noise, so if we did not filter the data to extract signals, then the resulting model would be biased towards generating noise rather than real signals. To perform this filtering we apply a simple energy detection process2.

Footnote 2: [https://github.com/FXi96/SEII-Energy-Detection](https://github.com/FXi96/SEII-Energy-Detection)

First we perform a simple bandpass filtering on the entire observation. This removes the polyphase filterbank shape imposed by the first of a two-stage channelization procedure (Lebofsky et al., 2019). This is done by collapsing the spectrogram in the time dimension and fitting a piece-wise polynomial spline to each \(\sim 3\,\mathrm{MHz}\)-wide coarse channel. We subtract the fitted polynomial from the data. Finally, we iterate through windows of size \(715\,\mathrm{Hz}\) and search for excess energy above the expectations of Gaussian random noise. We chose an S-score threshold of \(512\)(D'Agostino and Pearson, 1973).

We also perform a complementary energy detection process where we select regions consistent with Gaussian noise by inverting the threshold condition. When constructing the training set we use an equal number of "signal" and "noise" regions, thus balancing the dataset.

### Data and Preprocessing

To create the training set we randomly draw 30 observations from a total of \(12{,}000\), and draw an additional six observations for the test set (excluding those already selected as part of the training set). We chose \(30\) and \(6\) to keep computing time reasonable. Running energy detection required 30 minutes to process each file and with 30 examples it already provided more than 2 million training samples which is more than enough to demonstrate our algorithm's capabilities. The spectrograms have a time resolution of \(18.25\,\mathrm{s}\) and a frequency resolution of \(2.79\,\mathrm{Hz}\), giving the dataset \(3\times 10^{8}\) frequency channels.

We then split the band into \(\sim 715\,\mathrm{Hz}\) "snippets". This resulted in a training set of approximately 2 million snippets and a test set of 300,000 snippets. We then log normalize the data, add a constant to make all the values positive, and scale the data to have a final range within 0 and 1. Examples of the resulting snippets are shown in Figure 1. Note that the normalization is done per snippet independently of each other sample. These snippets are used as inputs to the autoencoder, and the targets are the same snippets.

### Autoencoder

As described in section 1.3.1, our autoencoder consists of an encoder and a decoder (Baldi, 2012). The encoder consists of five convolutional layers (LeCun et al., 1999), with filter sizes of \(16\), \(32\), \(32\), \(32\), and \(32\) respectively. In between each layer is a 2-D Maxpool layer (Krizhevsky et al., 2017) of size (1,2). A batch normalization layer (Ioffe and Szegedy, 2015) is included between the final maxpool layer and the convolutional layer. This is followed by three dense layers of size \(32\), \(16\), and \(5\) respectively. All the activations used are ReLu (Fukushima, 1975) activations. The model was built using Tensorflow (Abadi et al., 2015) and Keras (Chollet et al., 2015).

The decoder is similar but in reverse order. We once again have dense layers of size \(16\), and \(32\) where the output is then reshaped and fed into a convolutional transpose (Baldi, 2012) layer which upscales the image back to the same dimension of the input. Between the five convolutional layers, we have batch normalization and maxpool size of (2,1). The model architecture is shown in Figure 2.

### Training and Validation

The training scheme is standard. We fitted the model using the ADAM (Kingma and Ba, 2014) optimizer with a learning rate of \(1\times 10^{-4}\), utilizing an early stopping routine with the patience of 10 epochs trained in batches of 16 samples, for 100 epochs. We then visually evaluated how well the model reconstructed the original spectrogram. Figure 1 shows a few randomly drawn examples. The reconstruction achieves our desired level of accuracy in that visually the reconstructions appear similar to the inputs. We proceed with this model (Figure 3).

### Frequency Embedding

Signal morphology varies with frequency, given that transmitters of different types occupy different regions of the band. For example, WiFi and Bluetooth signals are common in our spectra around \(2.4\,\mathrm{GHz}\). However, during our training phase and in the construction of our model this frequency information is not preserved in the feature extractor. Should a user choose to search for signals that are also similar in frequency this information would be lost. However, we can add frequency information to the feature vector.

Initially, this appears trivial -- one can simply extend the feature vector by a dimension and add the frequency information here. However, this would mean that frequency is an orthogonal feature to all the other extracted features. We know, however, that signal morphology is correlated to some extent with frequency. Here we borrow techniques from Natural Language Processing (specifically transformer architectures) to encode frequency information in feature vectors called Positional Embeddings (Vaswani et al., 2017).

Positional embeddings work by taking encoded feature vectors and perturbing these vectors by small offsets based on their position to obtain an ideal balance between correct adjustment and over-adjustment which could lead to confusion with other signals. One can imagine this as trying to build miniature clusters in some high-dimensional feature space where each element of the cluster is unique but the elements of the cluster don't intersect with other clusters.

The adjustment vector is found using Equation 1, adapted from the paper "Attention Is All You Need" (Vaswani et al., 2017),

\[P(k,i)=\begin{cases}\sin(\frac{k}{n^{2}}),&i\,\,\mathrm{even}\\ \cos(\frac{k-k}{n^{2}}),&i\,\,\mathrm{odd}\end{cases} \tag{1}\]

The variable \(k\) is the index in position, so if we have a sequence of length \(L\), we pick an integer \(k\in[0,L]\). The index in the feature vector is denoted \(i\), and \(i\in[0,d]\) where \(d\) is the dimension of the embedding space. The variable \(n\) is a tunable hyperparameter, where Vaswani et al. choose \(n=10,000\). These equations ought to satisfy two conditions:

1. The adjustments are unique for a given frequency index \(k\).
2. The adjustments are bounded, in order to prevent the adjustment vector from over-adjusting.

The functions are all bounded between \([0,1]\) satisfying condition (ii). The sinusoidal functions are unique to each position, satisfying
Figure 1: A random sample of training examples. A variety of RFI signals are seen in these snippets.

Figure 3: Eight randomly drawn real observations (top row) and their corresponding autoencoder reconstructions (bottom row). The autoencoder can reconstruct the signals in the input data. We see that the reconstruction of the signal appears good, whereas the reconstruction of the noise is slightly poorer. This isn’t a concern since the main focus should be on the signal.

Figure 2: A representation of our autoencoder model. The encoder shows a progressive compression of data, forcing the network to decide the most important features of the original spectrogram.
condition (i). We can measure similarity in positions using metrics such as cosine similarity or Euclidean distance. The embeddings are visualized in Figure 4.

The choice of sequence length determines how many frequency chunks the band is broken into. Our choice of 1000 bins (giving \(\sim 1\,\mathrm{MHz}\) chunks) is a compromise between small chunks (which increase computational resource intensity) and large chunks (which reduce accuracy).

### The Search Algorithm

We now combine the components of our algorithm and use them to perform a search. The search has three elements. First, the feature extractor, using the encoder of the Autoencoder, extracts features from the spectrogram. Second, we index and construct the frequency embedding for each of the extracted features and add those to the encoded vectors. We repeat this process with the SOI. Finally, we compute similarity scores between the SOI and the images in the search list. This similarity score is the cosine similarity metric and is computed efficiently by performing a matrix multiplication and renormalizing by the respective vector norms. The procedure is visualized in Figure 5.

## 3 Results

We perform three sets of tests. First, we apply the classical SIFT + BoW algorithm to the problem discussed in Section 1.2. Second, we apply our algorithm to the problem but without the frequency embedding step. Lastly, we apply our algorithm with frequency embeddings to the problem and assess the results, in order to demonstrate the modularity of our design.

To assess the ability of the algorithm to return similar signals to those used as input, we select 10 filtered Energy Detection snippets at random and perform a search using the same observation set as used for the SOI observations. Images of the spectrograms used can be found in our GitHub repository3.

Footnote 3: [https://github.com/PetchMa/Reverse_Radio_Search](https://github.com/PetchMa/Reverse_Radio_Search)

In the first method, we use the classical SIFT + Bag of Words algorithm described in section 1.2. We run SIFT through each window of the spectrogram and generate a set of descriptor vectors for each image. Then we use K-Means clustering with 800 centroids to build a codebook. The codebook is built by taking the fitted centroids as the codewords and then looping back in the descriptor vectors and assigning them a codeword. Then for each image, we created a histogram of codewords, hence treating it as a BoW model. Finally, we performed a match between the histogram of the SOI and the histograms of each candidate signal. The matches use the K-Nearest Neighbours approach (Cover and Hart, 1967). Results are shown in Figure 6. One downside to this approach is that the algorithm is very slow. For each observation, we need to generate millions of descriptor vectors and repeat the SIFT operation millions of times. Unlike deep learning-based approaches, these feature vectors can be generated in parallel on the GPU using much more efficient methods.

Next, we applied our algorithm _without_ using frequency embedding techniques (Figure 7). We see that this approach was able to retrieve more convincing candidates than the classical SIFT+ BoW approach. When we investigate _where_ the candidates lie in frequency space, we see that many are close to the frequency of the true signal (Figure 7).

Finally, we applied our algorithm _with_ frequency embedding (Figure 8). The resulting signals show even higher visual similarity to the input signal. Additionally, the frequency distribution of the matches is closer to the frequency of the input signal, showing that the frequency embedding was successful.

## 4 Discussion

Our algorithm is more successful than classical approaches in returning similar candidates. The algorithm also runs approximately 20 times faster than our classical implementations at our data center because our algorithm can take advantage of GPUs with deep learning hardware acceleration such as tensor cores. Traditional matching techniques suffer an inability to generalize large-scale features as discussed in Section 1.2. Furthermore, the addition of frequency embedding results in candidates that have a better match in frequency to the input signal, which is desirable under certain situations. In addition, it appears that SIFT-based algorithms tend to match local features and thus have a tendency to match noise (see Figure 6), unlike the Deep Learning based approaches where models were able to learn more global features resulting in better matches overall.

Looking forward we believe this algorithm has a wide range of practical use in transient radio astronomy. One particularly exciting use case is to help automate signal verification steps for candidate technosignature signals such as BLC1 (Sheikh et al., 2021). By automatically searching for potential lookalikes based on signal morphology, rather than relying on more low-level parameters of a signal, such as drift rate, signal width, or signal-to-noise ratio, we can do a better job of distinguishing between true anomalies, and RFI.

More generally, this technique can be used to build an RFI database. Many RFI databases consist simply of frequency ranges that are excluded (for example, the frequency ranges corresponding

Figure 4: A visualization of the patterns in the embedding confirms that they are unique for each position we encode. We used a dimension of 512 and a sequence length of 100 for demonstration purposes; in the actual algorithm, we used a dimension of 4 and a sequence of 1000.

Figure 5: Visualization of the search process and the flow of data. We first extract features from the SOI and the set of possible candidates. Then we apply frequency embedding, and finally, we produce the similarity scores by matrix multiplication.
Figure 6: (a) The top 10 candidates with the closest similarity to the target [shown as the top left box], using the SIFT + BoW algorithm. The best candidate is shown in the top right box. Note: that this method did not find itself because this uses a KNN matching algorithm that excludes the trivial case (b) Frequency distribution of the top 10,000 most similar hits.

Figure 7: (a) The top 10 candidates with the closest similarity to the target [shown as the top left box] for our deep learning algorithm with no frequency embedding. The best match is unsurprisingly, the input signal [top right box]. (b) Frequency distribution of the top 10,000 most similar hits.
to GPS satellites). Our algorithm enables a finer-grained approach, comparing signal morphologies without necessitating the exclusion of large ranges in frequency, thus helping to make more efficient use of the spectrum and expanding the power of the search.

Another use case is to perform template searches. For example, one can build theoretical models that simulate some desired signal, and search the recorded observations for signals that match.

The next step in the development of our algorithm is to modify it to handle a wider range of data products. For example, developing a model to handle arbitrary _collections_ of spectrogram data (data from different receivers or even telescopes). Or perhaps building a model that can deal with cases where the spectral resolution varies (e.g. combining data products with 3 Hz resolution and 1 Hz resolution).

We also plan to extend our approach to series (or cadences) of multiple observations, which intersperse scans of the target star with comparison scans of neighboring targets. Signals such as BLC1, which appear only in scans of the primary target, are consistent with being spatially localized on the sky. However, the dataset in which BLC1 was discovered also contained "lookalike" signals at other frequencies, indicating that they were likely due to particularly pernicious RFI. By applying our new methodology to cadences of data, we can much more easily locate lookalike signals at other frequencies, or in scans of other targets, providing an additional powerful means of screening technosignature candidates.

## 5 Data Release and Code Availability

All source code is released here and the data is publicly released here.

## Acknowledgements

Breakthrough Listen is managed by the Breakthrough Initiatives, sponsored by the Breakthrough Prize Foundation4. We are grateful to the staff of the Green Bank Observatory for their help with installation and commissioning of the Breakthrough Listen backend instrument and extensive support during Breakthrough Listen observations. We thank Yuhong Chen for his helpful discussion on the Energy Detection Algorithm.

Footnote 4: [http://www.breakthroughinitiatives.org](http://www.breakthroughinitiatives.org)
