# There’s no difference: Convolutional Neural Networks for transient detection
without template subtraction

Tatiana Acero-Cuellar1,2, Federica B. Bianco2,3,4,5, Gregory Dobler2,3,4,5,
Masao Sako6, Helen Qu6 and The LSST Dark Energy Science Collaboration
1Observatorio Astronómio Nacional, Universidad Nacional de Colombia, Bogotá,
Colombia
2University of Delaware Department of Physics and Astronomy 217 Sharp Lab
Newark, DE 19716 USA
3University of Delaware Joseph R. Biden, Jr. School of Public Policy and
Administration, 184 Academy St, Newark, DE 19716 USA
4University of Delaware Data Science Institute
5Center for Urban Science and Progress, New York University, 370 Jay St,
Brooklyn, NY 11201, USA
6Department of Physics and Astronomy, University of Pennsylvania,
Philadelphia, PA 19104, USA

(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

We present a Convolutional Neural Network (CNN) model for the separation of
astrophysical transients from image artifacts, a task known as “real-bogus”
classification, that does not rely on Difference Image Analysis (DIA) which is
a computationally expensive process involving image matching on small spatial
scales in large volumes of data. We explore the use of CNNs to (1) automate
the “real-bogus” classification, (2) reduce the computational costs of
transient discovery. We compare the efficiency of two CNNs with similar
architectures, one that uses “image triplets” (templates, search, and the
corresponding difference image) and one that adopts a similar architecture but
takes as input the template and search only. Without substantially changing
the model architecture or retuning the hyperparameters to the new input, we
observe only a small decrease in model efficiency (97% to 92% accuracy). We
further investigate how the model that does not receive the difference image
learns the required information from the template and search by exploring the
saliency maps. Our work demonstrates that (1) CNNs are excellent models for
“real-bogus” classification that rely exclusively on the imaging data and
require no feature engineering task; (2) high-accuracy models can be built
without the need to construct difference images. Since once trained, neural
networks can generate predictions at minimal computational costs, we argue
that future implementations of this methodology could dramatically reduce the
computational costs in the detection of genuine transients in synoptic surveys
like Rubin Observatory’s Legacy Survey of Space and Time by bypassing the DIA
step entirely.

††pubyear: 2021††pagerange: There’s no difference: Convolutional Neural
Networks for transient detection without template subtraction–LABEL:lastpage

##  1 Introduction

Modern observational astronomy has shown us that the Universe is not static
and immutable: to the contrary, it is a lively and dynamic system. We now know
and understand a variety of different phenomena that can give rise to
variations in the brightness and color of astrophysical objects, including
explosive stellar death (supernovae, kilonovae, Gamma Ray Bursts, etc), less
dramatic and powerful stellar variability (flares, pulsations), and
variability arising from geometric effects (such as planetary transients and
microlensing). The time scales for these phenomena range from seconds to
years. When the variations are terminal –like supernovae– or stochastic –like
stellar flares– they are referred to as “transients”. Detecting optical
astrophysical transients characteristically requires sequences of images
across a significant temporal baseline. Surveys like the Dark Energy Survey
(The Dark Energy Survey Collaboration, 2005, DES) or the future Rubin
Observatory Legacy Survey of Space and Time (Ivezić et al., 2019a, LSST)
search for spatially localized changes in brightness in patches of sky
previously observed.

Due to the rarity of astrophysical transients, a long baseline in conjunction
with the observation of a large area of the sky is typically required to
detect a statistically significant sample of transients like supernovae and
individual examples of rare events, such as kilonovae. The process itself is
arduous and requires considerable human intervention at multiple stages. The
first step is typically the creation of high quality templates of each region
of the sky that is to be searched; the templates are then subtracted from
nightly images, a process known in astrophysics as Difference Image Analysis
(DIA) that was initially pioneered by Tomaney & Crotts (1996) and then
formalized by Alard & Lupton (1998) and there is a rich history of subsequent
improvements in the efficiency and accuracy of DIA models. Templates (tmpl)
are typically constructed as stacks of high quality (favorable observing
conditions) sky images. These high quality images must then be aligned with
the “today” image, typically called the “search image” (srch), and degraded to
match its Point Spread Function (PSF, which we note may vary across a single
large field of view image) and scaled to match its brightness. The product
generated by the subtraction of the template from the search image is the so-
called “difference image” (diff; see for example a description of the DIA
processing pipeline for DES in Kessler et al. 2015. Transients can then be
detected as clusters of adjacent pixels deviating significantly from the sky
noise. However, even the best existing DIA algorithms produce difference
images with large pixel value deviation from an ideal 0-average that would be
expected if no changes occurred in a patch of sky. Transients, variable stars,
and moving objects will result in detections, but also a typically large
number of artifacts will be detected by these thresholding schemes.

Machine learning models offer an excellent opportunity to improve the
efficiency of transient detection at this stage, automating the classification
between “real” astrophysical transients and “bogus” artifacts: these models
are often referred to as “Real-Bogus” (hereafter, RB). Generally, the
applications of machine learning to this problem are based on the extraction
of features from the (difference) images that are then fed to models like
Random Forests, k-Nearest Neighbors, or Support Vector Machines (Goldstein et
al., 2015; Sánchez et al., 2019; Mong et al., 2020). These models have
achieved high accuracy and enabled the discovery of transients at scale in
larger synoptic surveys, for example in the Palomar Transient Factory (PTF
Bloom et al., 2012) and Zwicky Transient Facility (ZTF Mahabal et al., 2019).

The process of engineering features for machine learning models allows the
experts to embed domain knowledge in models. In the case of RB classification,
the features engineered to be used by models such as the ones described above
usually rely on visual inspection of a subset of the data. However, these
features may overlook abstract associations between image properties that can
be effective for classification, and in fact, may be biased towards human
perception and theoretical expectations. An alternative approach involves the
use of models that can learn features directly from the data, such as
Convolutional Neural Networks (LeCun et al., 1989, CNNs). Here, it is possible
to train a model using the images themselves, skipping the step of feature-
design entirely.

RB CNN-based models appeared in the literature as early as 2017 (Gieseke et
al., 2017), and they are typically based on the analysis of the difference
images arising from the DIA. While neural networks may be computationally
demanding in the training phase, the "feed-forward" classification that arises
from a pre-trained model is typically rapid and computationally
light,111although the model may require a large amount of memory to cache.
leaving DIA as the computational bottle-neck in the process of astrophysical
transients’ detection. Yet in principle, the entirety of the information
content embedded in the diff-tmpl-srch image triplet is also contained in the
tmpl-srch image pair. In this paper, we explore the use of CNNs as RB models
concentrating on the potential for building high-accuracy models that do not
require the construction of difference images.

This paper represents a first step in the process of conceptualizing a DIA-
free model for astrophysical transients. Here we demonstrate the potential for
a model that discriminates between astrophysical transients and artifacts,
while still relying on DIA for detection. In future work we intend to develop
models that will not rely on DIA for detection, and, finally, that will detect
and characterize transients (measure magnitude) from the tmpl-srch image pair
only.

This paper is organized as follows: in section 1, the introduction, we briefly
describe DIA process, examples of architectures, and input images for CNNs
models on the literature and describe the importance of our work. In section 2
we discuss the history, approaches, and details of DIA in the literature. In
section 3 we present the DES data that we will use to build our RB
classification models and the pre-processing steps. In section 4 we presented
some general aspects of CNN and illustrate the CNNs architectures used, in
subsection 4.3 we present, share and discuss in some insights into the model
obtained enabled by saliency maps, as well as the method used to interpret the
data decisions through saliency maps. In section 5 we show the results of
building a model that does not use the difference image as input, named noDIA-
based model and compare its performance to one that does, named DIA-based
model; we outline future work in section 6 and we conclude with a discussion
of the broader implications for this result in light of upcoming surveys.

This study is reproducible and all the code that supports the analysis
presented here is available on a dedicated GitHub repository.
222https://github.com/taceroc/DIA_noDIA.

Figure 1: Example of four DIA-sets from the first season of the Dark Energy
Survey (DES): from left to right the images correspond to the template (tmpl),
search (srch), and difference (diff) image; the diff is generated as the
subtraction of tmpl and srch. Each pair of tmpl and srch images is mapped to
the same color range. We refer to these 3-images sets as “image triplets” or
DIA-sets. A and B show artifacts, human-labeled as “bogus” (label = 1). C and
D show transients labeled as “real” (label = 0). Above each triplet is the
unique ID of the transient (see Goldstein et al. 2015)  Figure 2: Histogram of
pixel values for the data shown in Figure 1 (before scaling and
normalization). While all difference images (right) show bell-shaped
distribution, the template and search image in both “bogus” and “real” sets
have different behaviours. For example, B and D have similar pixel values
distribution, however, B is “bogus” and D is “real”. The vertical lines show
the $\mu\pm 3\sigma$ interval for the srch and tmpl images. The diff images
(last right column) are standardized individually to a mean of $0$ and a
standard deviation of $1$. The srch and tmpl are instead scaled, setting the
pixel contained inside the 3$\sigma$ interval (vertical lines on the
histograms) to the range 0-1. This allows to retain negative values or as well
values above $1$ while keeping the core of the distributions to within a
homogeneous range.  Figure 3: Horizontally stacked data for the examples shown
in Figure 1. On the left, the composite images used as input to the DIA-based
CNN model: the composite follows the order, from left to right: diff, srch,
tmpl. On the right, the composite images used for the noDIA model, composed of
srch and tmpl. Each images element was scaled or normalized following the
description given in section 3 and Figure 2 before combining them into a
single image. Above each composite are the unique transient “ID”, the original
label, and prediction made by our model. The four transients were classify
correctly by both DIA-based and noDIA models. Purple shades indicate negative
and green positive pixel values.

##  2 state of the art and traditional solutions to detection of transients

###  2.1 Difference Imaging

Difference images produced subtracting a template generated coadding multiple
images (e.g. Kessler et al., 2015) from a sky image are currently the basis
for most astrophysical transient search algorithms. The difference image
allows brightness changes to be detected even if embedded in Galaxy light, for
example, in the case of extra-galactic explosive transients. Great efforts
have been made to improve the quality and effectiveness of the difference
images. Although the name may suggest the process simply entails subtracting
images from each other, the procedure is in fact riddled with complications
because of the following reasons. First, the images used to build the template
and the search images are taken principally within different atmospheric
conditions (Zackay & Ofek, 2017) generating variations in the quality of the
images. The construction of a proper template is also a delicate task;
typically templates are built by stacking tens of images taken under
favourable sky conditions at different times. This improves the image quality
but also mitigates issues related to variability in the astrophysical objects
captured in the image (Hambleton et al., 2020): one wants to capture each
variable source at its representative brightness. Typically then, the template
image is of higher quality than the search image, and it is degraded to match
the search image PSF and scaled to match its brightness. Yet, the scaling and
PSF may vary locally in the image plane, especially for images from large
field of view synoptic surveys such as the 2.2 sq. degrees DES or $\sim 10$
sq. degrees Rubin images. Finally it is important that the images are
perfectly aligned, both in the creation of the template and in the subtraction
process to create the difference image. This implies accounting for rotation
as well as potentially different warping effects on the images and template.
Once the PSF match and the alignment is done, it is possible to subtract the
degraded template from the search images to obtain the difference image. To
degrade the image quality of the template to match the search image, a
convolution kernel that must be applied to the template needs to be determined
(Alard & Lupton, 1998),

$\displaystyle\text{{\it tmpl}}(x,y)\otimes\text{Kernel}(u,v)=\text{{\it
srch}}(x,y),$ (1)

where tmpl is the high-quality image, the template image, srchis the one night
image or search image and $\otimes$ is the convolutional operation. The
arguments $x$ and $y$ represent the coordinates of matrix that compose the
images; $u$ and $v$ the coordinates of the kernel matrix.

To solve the computationally expensive problem of matching PSFs, the kernel
can be decomposed in terms of simple functions, for instance Gaussian
functions, and the method of least squares can be used to determine the best
values for the kernel. The fitted solution of one search image can be
determined in a short computational time. However, the computational cost
scales with the image size and resolution. Surveys and telescopes constructed
with the goal of discovering new transients are generally designed to collect
tremendous amount of data to maximize event rate (detection of astrophysical
transients). For the DES, the computational cost per 2.2 sq. degree image
$\sim 15.5$ CPU hours (with roughly 2/3 of that time spent on PSF matching).
The upcoming Rubin LSST will collect more than 500 images every night each
with 3.2 Gigapixels. This process thus is bound to turn out to be very
expensive. In this paper, we trained our models on postage stamps where
transients were detected or simulated (see section 3), thus a direct
comparison of the computational cost is not trivial. A detailed discussion of
the computational cost of our models is included in subsection 5.2 and the CPU
Node hours required to train and generate predictions from our models are
reported in Table 1. We note here that, with a Deep Neural Network approach to
this problem, the computational cost is high in training, but the predictions
are essentially instantaneous.

A bad subtraction can occur either because of poor PSF matching, poor
alignment, or poor correction of image warping. In all of these cases, the
subtraction would lead to artifacts or “bogus” alerts, like the one show in
the difference image in Figure 1A and B. In particular: in Figure 1A, the
difference image shows a so-called “dipole”, where one side of a suspected
transient is dark and the other side bright: this typically arises in case of
mis-alignments, but it might also be caused by moving objects in the field, or
differential chromatic refraction (Carrasco-Davis et al., 2021). Conversely,
Figure 1B shows a “bogus” alert caused by an image artifact: a column of bad
pixels in the search image. At this location there is no astrophysical object
in the image thumbnail: no host galaxy or star that could give rise to
variations.

Panels Figure 1C and Figure 1D show genuine transients in our DES training
data: in both of these two examples, there is a clear transient in the diff
images (high pixel values at the center of the image).

###  2.2 Autoscan and other feature-based Real Bogus models

We developed our models on data collected in the first year of DES. Thus, a
direct precursor of our work is Goldstein et al. (2015), in which the authors
created an automated RB based on a Random Forest (RF) supervised learning
model (Ho, 1995) to detect transients, and particularly supernovae, in the DES
data, hereafter refered to as autoscan. For these kinds of models, the process
of selecting and engineering features is pivotal. Autoscan is based on 38
features derived from the diff, srch and tmpl images. The selection and
computation of these features was done attempting to represent quantitatively
what humans would leverage in visual inspections. For instance, r_aper_psf
distinguishes a bad subtraction of srch and tmpl that would lead to a diff
qualitative similar to Figure 1A; the feature diffsum measures the
significance of the detection by summing the pixel values in the center of the
diff image; the feature colmeds, indicating the CCD used for the detection, is
designed to identify artifacts specific of a CCD, like bad rows/columns of
pixels.

In other RB models, like in Sánchez et al. (2019), the feature selection is
performed purely statistically: features were initially selected based on
variance thresholds. In the same work, different techniques are explored to
reduce the number of features, and thus the complexity of the classification
problem. For example, a RF model was trained using all features. Then, feature
values were randomly modified and the performance on the “new” data set was
analyzed. Decreased accuracy associated to a change of a feature measured the
importance of the feature that was modified. In this way the authors reduced
the dimensionality of the problem by removing possible redundant or irrelevant
information. Examples of models based on features closer to the data include
Mong et al. (2020), where the features are simply the flux values of the
pixels around the center of the image.

###  2.3 AI approach

CNNs have demonstrated enormous potential in image analysis including object
detection, recognition, and classification across domains (Deng et al., 2009).
Examples of astrophysics applications of CNNs include Dieleman et al. (2015)
for galaxy morphology prediction, Kim & Brunner (2016) for star-galaxy
classification, Gabbard et al. (2018) for signal/background separation for
Gravitational Waves (GW) searches, where the GW time series are purposefully
encoded as images to be analyzed by a CNN, and many more.

CNNs are particularly well-suited to learning discriminating features from
image input data. CNNs can work on high dimensional spaces (here the
dimensionality of the input is as large as the number of pixels in the image)
due to the generalizability of the convolution operation to $n$ dimensions
while preserving relative position information. Vectors of raw pixel values
can theoretically be used to train traditionally feature-based models, such as
RFs, but pixel-to-pixel position data in higher dimensions is unequivocally
lost. Previous studies already compared feature-based supervised models, like
RF, and supervised CNNs for RB, demonstrating that CNN generally leads to
increased accuracy. In Gieseke et al. (2017) an accuracy of $\sim 0.984$ is
achieved with a RF model in the RB task, and it increases to $\sim 0.990$ when
applying a CNNs to the same data. In Cabrera-Vives et al. (2016) accuracy of
$\sim 0.9889$ achieved with a RF is increased to to $\sim 0.9932$ with a CNN.
In Cabrera-Vives et al. (2017), an RF gave 0.9896 and a CNN 0.9945. In Liu et
al. (2019) the accuracy improves from $\sim 0.9623$ with a RF to $\sim 0.9948$
with a CNN.

Beyond RB, a model for image-based transient classification through CNNs has
been prototyped by the Automatic Learning for the Rapid Classification of
Events team (ALeRCE Carrasco-Davis et al., 2021). The model classifies between
AGNs (Active Galactic Nuclei), SNe (SuperNovae), variable stars, asteroids,
and artifacts in The Zwicky Transient Factory (Bellm et al., 2018, ZTF) survey
data with a reported accuracy exceeding 95% for all types, except SNe (87%).
This CNN model was trained using a combination of the srch, tmpl, and diff.

The problem we propose to address here is exclusively the RB classification,
but we explore the potential of leveraging AI to bypass the DIA step. Existing
CNN RB models differ not only in their architecture (for example, a single or
multiple sequences of convolutional, pooling, dropout, and dense layers) but
also in the choice of input. For instance, Gieseke et al. (2017) used the
tmpl, srch, and diff images in combination for training their CNN; Cabrera-
Vives et al. (2016) augmented this image set with an image generated as the
diff divided by an estimate of the local noise; Cabrera-Vives et al. (2017)
trained an ensemble of CNNs on different rotations of this four-fold image
set. Conversely, Liu et al. (2019) used the diff as the sole input. Although,
all these attempts have shown good results with accuracy higher than $90\%$,
these models all rely on DIA to construct the diff. Taking into account that
diff are built from the tmpl and srch, and in principle, should carry no
additional information content than the search-template pair alone, a logical
step to follow would be to only consider the two latter images.

A first attempt in this direction is presented in Wardęga et al. (2020).
Wardęga et al. (2020) trained a model that could distinguish between optical
transient and artifacts using a search image from Dr. Cristina V. Torres
Memorial Astronomical Observatory (CTMO, a facility of the University of Texas
Rio Grande
Valley333https://www.utrgv.edu/physics/outreach/observatory/index.htm) and a
template image from the Sloan Digital Sky Survey (Gunn et al., 2006, SDSS).
They trained two Artificial Neural Network models, a CNN and a Dense Layer
Network on simulated data and tested the models using data from CTMO and SDSS.
The data used for training and testing had specific characteristics:
transients were a combination of a source in the CTMO images (or srch image)
and background in the SDSS image (or tmpl image); artifacts were a combination
of a source in both the CTMO and SDSS images. Within this dataset, both models
yield high accuracy ($>95\%$) and establish a precedent for the possibility of
avoiding the construction of the diff to reliably detect optical transients.
However, studies based on more diverse and realistic data, e.g., sources near
galaxies, or embedded in clusters, are needed to demonstrate the feasibility
of this approach. The data used in this paper fulfills this condition and is
described in section 3.

##  3 Data

This study is designed as a detailed comparison of CNN-based models, with and
without diff in input, with the results from one RB model, autoscan (Goldstein
et al., 2015). The choice of autoscan as our point of reference is motivated
by its application to the discovery of transients in a state of the art
facility, the Dark Energy Survey (The Dark Energy Survey Collaboration, 2005),
which can be considered a precursor of upcoming surveys like the Rubin Legacy
Survey of Space and Time (Ivezić et al., 2019b). The latter, expected to start
in 2024, will deliver $\sim 20~{}\mathrm{Tb}$ of high resolution sky image
data each night, covering a footprint of $\sim 20,000$ sq. deg. every $\sim 3$
nights, with expected millions of transients per night, demanding rapid
methodological and technical advances in accuracy and efficacy of transient
detection models.

The data used in this work consists of postage stamps of images collected by
the DES during its first observational season (Y1), August 2013 through
February 2014444The data can be found in
https://portal.nersc.gov/project/dessn/autoscan/# (Abbott et al., 2018). The
data corresponds to 898,963 DIA-sets, a template (tmpl) image, search (srch)
image, and their difference (diff). The construction of the templates images
for the DES Y1 leveraged the data collected in season two (Y2) as well as the
Science Verification images (observations collected prior to survey start in
order to evaluate the performance of the instrument). More information on the
DES DIA pipeline can be found in (Kessler et al., 2015). Of these DIA sets,
454,092 contain simulated SNe Ia, which constitute the “real” astrophysical
transients set (label = 0) and 444,871 are the actual data from DES, i.e., the
“bogus” set (label = 1). Each image is offered in two formats: “.gif”, and
“.fits”. The former is an 8-bit compressed format (convenient for visual
inspection as it can be opened with commonly available software); the latter,
the “Flexible Image Transport System”, is a common data format for
astronomical data sets which enables high precision with a large dynamic rage.
The data in “.fits” format was used in this work.555More information related
to how to manipulate this astronomical data format is available in
https://docs.astropy.org/en/stable/io/fits/ Each image is $51\times 51$ pixels
- corresponding to approx 0.26 arcseconds square of sky.

Some examples of the data are in Figure 1. Each transient is identified by a
unique “ID”. The metadata, includes the labels associated with each image as
well as the the $38$ features used for classification in Goldstein et al.
(2015). Because we only analyze postage stamps with detections, we implicitly
still rely on the DIA to enable the detection step at this stage of our work.
The tmpl and diff image in our postage stamps, however, are not PSF matched.

###  3.1 Scaling and normalization

A word about data preparation and normalization is in order as astrophysical
images are inherently very different from the images upon which CNNs have been
built. When training CNN models for image analysis, each image is typically
simply scaled to a common range ($0-1$). However, the dynamic range of an
astrophysical image is typically large and the distribution of pixel values is
generally very different from Gaussian, with the majority of pixels sitting at
low values (the sky) and a few pixels at or near saturation (which in some
cases may carry the majority of the information content). Furthermore, in the
DIA-set case, the pixel-value distribution of the diff differs qualitatively
from the tmpl and srch ones. While the tmpl and srch are typically naturally
positive valued with a long tail at the bright end (right-skewed because of
the presence of bright astrophysical sources such as galaxies that host
transients or stars that vary), the diff image is, in absence of variable or
transient sources, symmetric around 0.

The diff images were standardized to have a mean $\mu=0$ and a standard
deviation $\sigma=1$. The srch and tmpl images, instead, were scaled to map
the $\mu\pm 3\sigma$ interval of the original image to $0-1$. This scheme
allows us to retain resolution in the shape of the core of the distribution
while also retaining extreme pixel values. Figure 2 shows the distribution of
pixel values for four DIA sets, for the same data as in Figure 1, which
include two “real” and two “bogus” labels. Appendix A shows the pixel
distributions before and after scaling for the same data.

Figure 4: Architecture of the Neural Networks used in this project to classify
“real” and “bogus” transients. Left: the DIA-based model that uses image
triplets as input (diff, tmpl, srch). The input layer is $51\times 153$ (see
left Figure 3); a convolution layer $(5\times 5)$ learns $16$ filters; max
pooling $(2\times 2)$ and dropout; convolution $(5\times 5)$ learns $32$
filters; maximum pooling $(2\times 2)$ and dropout; convolution $(5\times 5)$
learns $64$ filters; maximum pooling $(2\times 2)$ and dropout; flatten layer,
Dense $(32)$ and the output is a Dense $(2)$-class layer. Right: the noDIA
model that uses the tmpl and srch images only. The input layer is $51\times
102$ (see right Figure 3); a convolution layer $(7\times 7)$ learns $1$
filter; maximum pooling $(2\times 2)$; convolution $(3\times 3)$ learns $16$
filters; maximum pooling $(2\times 2)$ and dropout; convolution $(3\times 3)$
learns $32$ filters; maximum pooling $(2\times 2)$ and dropout; flatten layer,
Dense $(32)$ and the output was a Dense $(2)$-class layer.The illustrations
were made using NN-SVG tool by LeNail (2019).

One further decision has to be made in combining the three images in the DIA
set to feed them to the CNN. We stacked the scaled diff, srch and tmpl
horizontally. Thus the size of the data in input to our CNN is ($N_{tr}\times
51\times 153$), where $N_{tr}$ is the number of transients to be considered.
Four examples of the data “triplets” in input to our DIA-based CNN are in the
left panel of Figure 3. Following this horizontal structure, we mimic closely
the way that human scanned this type of data for classification, and we will
take advantage of this scheme when examining the models’ decisions in
subsection 4.3, although stacking the images along a third depth axis has also
been tested, with no differences in the model accuracy.

But our goal is to test the performance of models that do not use the DIA.
Thus for our noDIA models the training and testing data set were constructed
in the same way as the previous triplets, but without the diff image. Some
examples are in the right panel of Figure 3.

##  4 Methodology

###  4.1 Basics of Neural Networks

Neural Networks are models that learn important features and feature
associations directly from data. These are typically supervised learning
models used for classification or regression. They consist of a series of
layers of linear combinations of the input data each with real value
parameters known as weights and biases, combined with activation functions
that enable learning non-linear, potentially very complex, relationships in
the data. The weights tell us the relevance of the input feature with respect
to the output, and the biases are the threshold value that determine the
output. Fitting these quantities to the data minimizes the loss, meaning the
prediction is as close as possible to the original target/label (Nielsen,
2015). For Deep Neural Networks (DNNs), “deep” refers to the fact that there
are multiple hidden layers, besides the input and output layer. With this
architecture the features learned by each layer do not follow a human
selection, rather, the features arise in the analysis of the data (LeCun et
al., 2015). Among DNNs, CNN models have layers that represent convolutional
filters. More information on DNNs and CNNs can be found in Agarap (2018);
Krizhevsky et al. (2012); Dieleman et al. (2015), as well as Gieseke et al.
(2017), among many others.

We used the
Keras666https://keras.io/api/layers/convolution_layers/convolution2d/.
implementation of the CNN (Chollet et al., 2015).

###  4.2 DIA and noDIA based models architecture

The input of our Neural Networks are the horizontally stacked images of size
$51\times 153$ for the DIA-based model (diff, srch, tmpl), and $51\times 102$
for noDIA model (srch, tmpl). For both the DIA-based model and the noDIA model
$100,000$ images were used, $80,000$ images for training and $20,000$ for
validation. The images were selected randomly from the $898,963$ and while
certainly training with a larger set can lead to higher accuracy, the amount
of images was sufficient for this comparison of DIA-based and noDIA models
while being conservative with limited computational resource. The data is
composed of $50,183$ images labeled as “bogus” and $49,817$ labeled as “real”.

The network architectures used for this work are shown in Figure 4, left
panel, for the DIA-based model, and right panel for the noDIA model. More
details about the architecture can be found in Appendix B. Both architectures
follow a similar structure. In designing the neural networks, we started with
the DIA-based model and developed an architecture that would match the
performance of Goldstein et al. (2015). With this model in hand, we created a
noDIA with a similar structure in order to enable direct comparison and
measure the effect of the change in the input data. We emphasize that in the
development of our models we followed two general guidelines:

  1. 1.

We pushed the accuracy of the DIA-based CNN model to match the accuracy of
autoscan. While a more exhaustive architectural exploration or a
hyperparameter grid search may well lead to increased efficacy, matching the
accuracy of autoscan (at 97%) is sufficiently for our demonstration. We
further note that this historical data set cannot be further validated (by,
for example, assessing the recurrence of transients at a sky position to
validate the bogus nature of a non-simulated detection) and it is subject to
potential label inaccuracies that may compromise the ability to achieve a
better performance. The final architecture used is one that reached the same
accuracy as Goldstein et al. (2015), and False Positive and False Negative
Rates most similar to the ones obtain in Goldstein et al. (2015) (see Figure 5
and section 5). Once we matched the autoscan performance, we focused on the
potential for removing the diff image from the input.

  2. 2.

While the architecture of the DIA-based CNN was designed with the attempt to
achieve a specific target performance, the architecture of noDIA is
deliberately kept as close as possible to that of the DIA-based model. This
enables a direct comparison of the effects of the removal of the diff image.
In other words, the hyperparameters and architecture of noDIA, shown on the
right in Figure 4 were not optimized explicitly for RB classification, rather
we chose a standard architecture, only modifying the original design to adapt
for the different dimensionality of the input data. One further deliberate
modification is implemented in the choice of a single-filter first layer for
noDIA. The CNN that is not offered the diff image needs to learn the image
PSF, which is constant across the postage stamp-size image, and can be modeled
with a single filter.

###  4.3 Saliency maps

Saliency maps quantify the importance of each pixel of an image in input to a
CNN in the training process. They provide some level of interpretability
through a process akin to feature importance analysis by enabling an
assessment of which are the pixels the model relies on the most for the final
classification. If the task were, for example, to identify cats and dogs in
images, the expectation will be to find that the most important pixels are
located within the the dog or cat bodies, and not in the surrounding, while if
the task were to identify activities performed by cats and dogs, one may find
the important pixels both within the dogs and in the surrounding, particularly
in objects associated with the performed tasks. Furthermore, some portion of
the subject’s body may be most distinctive (ears, nose) and we would expect
more importance being given to those pixels. The "importance": of a pixel
would be simply measured by the weight the trained NN assigns to that pixel in
the case of a single-layer perceptron, but in the case of Deep Neural
Networks, a highly non-linear operation is performed on the input data, and
the importance, or saliency, of an input feature, or pixel, is harder to
assess. We follow (Simonyan et al., 2014)’s definition of saliency: denoting
the class score with $S_{c}$, such that the NN output on image $\mathcal{I}$
is represented by $S_{c}(\mathcal{I})$, the saliency map is given by :

$\displaystyle S_{c}(\mathcal{I})\approx w^{T}\mathcal{I}+b,$ (2)
$\displaystyle w=\frac{\partial
S_{c}}{\partial\mathcal{I}}\biggr{\rvert}_{\mathcal{I}_{0}}.$ (3)

That is: $w^{T}$ is the weight and $b$ the bias of the model. $w$ is the
derivative of $S_{c}$ with respect to $\mathcal{I}$ calculated specifically in
the local neighborhood of pixel $\mathcal{I}_{0}$, and the approximation sign
indicates a first order Taylor expansion has been used to approximate the
solution.

Figure 5: Confusion matrix (CM) for our testing data, a set composed of
$10,078$ objects labeled as “real” and $9,922$ as “bogus”. Squares of the
matrix, from the top left in the clockwise direction, indicate: True Positive
(TP, label = 0, prediction = 0 ), False Negative (FN, label = 0, prediction =
1), True Negatives (TN, label = 1, prediction = 1), and False Positive (FP,
label = 1, prediction = 0). We note that here, somewhat unusually, 0
corresponds to “real” and “positive”, and 1 to “bogus” and “negative”, as we
chose to remain consistent with the original labeling of the data presented in
(Goldstein et al., 2015). Left: The CM of our DIA-based model shows that from
the $10,078$ transients labeled as “real”, $9564$ (i.e., 95%) were correctly
classified and the accuracy for the TN is even higher, at $97\%$. Right: The
CM for the noDIA model shows that of the $10,078$ transients labeled as
“real”, $9163$ (i.e., 91%) were predicted correctly and the TN accuracy is
$92\%$.  Figure 6: ROC curve for the $20,000$ images used for testing, Blue
line for the DIA-based model, the area under the curve is $0.992$. Orange
line, for the noDIA model, the area under the curve is $0.973$. This figure is
discussed in section 5. Figure 7: Loss (top) and accuracy (bottom) for both
models presented in this work. Purple lines correspond to the validation data
(20% or 20,000 images) and green to the training data (80% or 80,000 images).
Left: results for our DIA-based model (Figure 4). Right: results for the noDIA
model (Figure 4).

Each pixel in an image in the training set is associated with the
corresponding pixel in the saliency map, and the saliency score (the
importance) of that pixel measures the change in model output as a function of
changes in the value of that input pixel by back-propagation. The higher the
value of a pixel in a saliency map, the more influence that pixel has in the
final classification. In this work we also refer to these maps as maps of
pixel importance.

In our case, given the side-by-side organization of the three elements of the
input image set, the saliency maps can help us assess how much the DIA-based
model relied on the diff to enable correct classification, and, thus, provide
some intuition in the difficulty of the challenge offered to the noDIA model.

For the DIA-based model we have an expectation guided by intuition: that a
greater concentration of important pixels should be found in the diff image.
In section 5, we will consider the veracity of this hypothesis both
qualitatively by visually inspecting the saliency maps, and by designing a
saliency-based metric that enables a quantitative approach. We calculated the
normalized sum of the saliency pixel values for each third of the image-
triplet, corresponding to diff, srch, tmpl:

$\displaystyle\begin{split}I_{\mathrm{diff}}=\frac{\sum_{p_{d}}{s_{p_{d}}}}{\sum_{p}{s_{p}}};\\\
I_{\mathrm{srch}}=\frac{\sum_{p_{s}}{s_{p_{s}}}}{\sum_{p}{s_{p}}};\\\
I_{\mathrm{tmpl}}=\frac{\sum_{p_{t}}{s_{p_{t}}}}{\sum_{p}{s_{p}}};\end{split}$
(4)

where $I$ indicated the importance of a segment of the input image, $p$ is an
index running over all pixels and $s_{p}$ indicates the value of pixel $p$ in
the saliency map. The denominator normalizes each metric by the total sum of
the saliency pixel values. In the numerators, $p_{d}$ denotes the pixels in
the leftmost third of the image, corresponding to diff; $p_{s}$ the middle
third corresponding to the srch; and $p_{t}$ the right most third,
corresponding to the tmpl. This metric allows us to assess the relative
importance of the diff (srch or tmpl) component of the image in performing RB
classification. Results from these metrics are discussed in detail in
subsection 5.1.

##  5 Results

The accuracies of our models are presented in Table 1. The DIA-based model
reached, by design, the accuracy of our benchmark model autoscan: 97% on True
Negative (TN) rate and 95% on True Positive (TP) rate (see section 4). We
remind the reader that we use the autoscan convention for the definition of TN
and TP: Positive is a “real” transient (label = 0), negative is a “bogus”
(label = 1). There is a drop of $\sim 4\%$ between the accuracy of the DIA-
based and the noDIA models. Along with the accuracy, in Table 1 we present
computational costs of training on 20,000 images measured in CPU Node
hours,777The computational cost was calculated by using the formula given in
https://docs.nersc.gov/jobs/policy/ and the allocations information given by
the system Iris https://iris.nersc.gov. Training time is reported in CPU
hours. Prediction time is reported per “postage stamp” (1ms) and the
classification operation is trivially parallelizable (can be run independently
on each postage stamp). and the clock time for the prediction for one single
image.

Model | Accuracy | Train time7 | Prediction7
---|---|---|---
Train | Test | (CPU hours) | (Clock-time, ms)
DIA-based | 0.962 | $0.959\pm$ | $\sim$35 | $1.00\pm 0.03$
$0.004$
noDIA | 0.915 | $0.916\pm$ | $\sim$56 | $0.30\pm 0.01$
$0.007$
Table 1: Comparative table of the accuracy and computational cost for both
training and testing data for the DIA-based (Figure 4, left) and noDIA (Figure
4, right) models.

Confusion matrices for the testing set are shown for the DIA-based model in
left panel of Figure 5 and for the noDIA model in the right panel of the same
figure. In this figure, as in the following confusion matrices and histograms
that we will present in subsection 4.3, correct predictions are indicated in
shades of blue, incorrect in shades of orange, and darker shades are
associated with the true labels. The percent accuracy for each class, True
Positives (TP), True Negatives (TN), False Positives (FP), and False Negatives
(FN), as well as the number of images in each class are reported within the
figure.

The Receiver Operating Characteristic (ROC) curve shows the relation between
the True Positive Rate (TP / (TP + FN)) also know as recall, and the False
Positive Rate (FP / (FP + TN)) when changing the threshold value (e.g., a
threshold of 0.5 will indicated that values greater than 0.5 would be
classified as “bogus”). The ROC for the testing data for the DIA-based and
noDIA models are presented in Figure 6. The Area Under the Curve (AUC), which
can be used as a comprehensive metric of the aggregated classification
performance of a model (Hanley & McNeil, 1983; Hernández-Orallo et al., 2012),
is 0.992 and 0.973 for the DIA-based and noDIA models respectively,
demonstrating good performance for both models.

The loss and the accuracy curves in the left panel of Figure 7 for the DIA-
based model shows some evidence of overfitting (the validation curve flattens
compared to the training curve) starting just a few epochs before training
came to end of the 400-epochs; meanwhile, for the noDIA model in right Figure
7, after $700$ epochs there was no visual evidence of overfitting indicating
that the model is still learning generalized information from the data; yet
the accuracy improvements from epoch $400$ to $700$ were small.

The nature of the noDIA model leads to hypothesize that because the input data
contains less information, this model takes longer to learn features from the
data to be able to classify them. The noDIA model took in fact longer (more
epochs) to come to stable accuracy and loss values. The loss and accuracy
curves are also noisier for the validation of the noDIA model in right Figure
7 compared to the DIA-based model. This also can be explained with the same
argument: the noDIA model had a harder problem to solve and its learning is
less linear. We cannot however rule out, and in fact we suspect, that the
hyperparameters or the architecture selected for the noDIA CNN model could be
improved — as this is a proof of concept paper and the architecture was chosen
to enable a direct comparison with the DIA-based model, future work will
include experimenting with different architectures and an extensive grid-
search to optimize hyperparameters.

Figure 8: Saliency map for the four transients in Figure 1 for the DIA-based
model. On the left, in grey color, the original diff, srch, and tmpl images
are plotted in their natural flux scale (before normalization). On the right,
the saliency map for the combined image. The intensity of a pixel color in the
white-to-green scale indicates the pixel relative importance: the maps are
normalized to 1 individually, such that dark green corresponds to high
saliency score, with 1 corresponding to the most important pixel in the image
triplet. With the side-by-side organization of the input data, these maps
enable a visual understanding of the importance of each element of the
combined image in the real-bogus classification. We note how in some cases
(panel A) the decision is largely based on the tmpl, rather than the diff
image, and in some cases all three image elements contribute similarly to the
decision (panel B). This figure is discussed in more detail in subsection 4.3.

Figure 9: As Figure 8 but for the noDIA model: saliency map for the four
transients in Figure 1. This figure is discussed in more detail in subsection
4.3.

.

###  5.1 A peek into the model decisions through saliency maps

Figure 10: Confusion matrix reporting the proportion of transients for which
the highest concentration of important pixels is found in the diff, srch, or
tmpl portion of the input image for the DIA-based model results (left) and
noDIA model (right), Left: DIA-based. For 90% ($8630$) of the $9564$
transients classified correctly as “real”, the classification principally
relied on the the diff image ($I_{\it diff}>1/3$); for 9% ($847$) on the the
srch; and for 1% ($87$) on the the tmpl. For incorrect “real” classifications
86% of the images relied principally on diff 7% on srch, and 7% on tmpl. For
incorrect “bogus” classifications 95% of the images relied principally on diff
3% on srch, and 2% on tmpl. For correct “bogus” classifications 61% of the
images relied principally on diff 12% on srch, and 27% on tmpl. Right: noDIA.
Of the $9163$ transients classified correctly as “real”, for 68% of them, the
classification relied principally on the tmpl image. For the incorrect “real”
and correct “bogus” classifications 76% of the cases were principally based on
the tmpl. For the incorrect “bogus” classifications 66% of the cases were
principally based on the tmpl.

Figure 11: The distribution of $I_{\mathrm{diff}}$ (top), $I_{\mathrm{srch}}$
(middle), and $I_{\mathrm{tmpl}}$ (bottom) value for the 20,000 transients in
the test data (see Equation 4). The colors corresponds to the quadrants of the
confusion matrix to which the transient belong according to the model
prediction, DIA-based predictions in the first and second column from the left
to the right and noDIA in the third and fourth. Blue shades correspond to
correct predictions (TP and TN) and orange to FP and FN. Note that the $y$
axis values are different in each plot and the that FP and FN histograms
contain far fewer observations. On the $I_{\mathrm{tmpl}}$Left, for the DIA-
based model, the TP classification shows a preference for the diff images and
the distribution peaks at $I_{\mathrm{diff}}\sim 0.5$. For TN, however, small
$I_{\mathrm{diff}}$ values are more common, with an significant fraction of
observations in the $0<$$I_{\mathrm{diff}}$$<0.33$ region. This behaviour is
complemented by a long $I_{\mathrm{tmpl}}$tail in the TN distribution
($0.4<I_{\mathrm{tmpl}}<0.8$). The FP and FN distributions are qualatively
similar to the TP and TN, but noisier, as they contain fewer than 10% of the
objects. Right: noDIA model. All four classes have qualitative similar
distributions in $I_{\mathrm{tmpl}}$and $I_{\mathrm{srch}}$. Classifications
rely mostly on the tmpl in all cases ($I_{\mathrm{tmpl}}>0.5$)

In subsection 3.1, we described how the importance of individual image pixels
in the RB prediction performed by our models can be measured, and the design
of a saliency-based metric to assess which component of the image is most
important to perform the RB classification. Here we inspect the saliency maps,
both visually and quantitatively through the measured values of
$I_{\mathrm{diff}}$, $I_{\mathrm{srch}}$, and $I_{\mathrm{tmpl}}$.

In Figure 8 and Figure 9 we show the four transients we considered as examples
throughout this work, the same images used in Figure 1 and Figure 3, and the
corresponding saliency maps for the DIA-based model and noDIA model
respectively. Figure 10 and Figure 11 report the results of Equation 4 for the
objects in our training set.

Lets start with some considerations about the saliency maps for the four image
examples for the DIA-based model (Figure 8), and specifically from panels C
and D were the transients were correctly predicted as “real”. We observe that
the greatest concentration of important pixels for both these images is found
in the left-most third of the image: $I_{\mathrm{diff}}\sim 0.5$ for both. We
speculate, from our experience in labeling real/bogus by visual inspection and
consulting with some of the human scanners that labeled the original autoscan
images, that this behaviour is similar to what a human scanner would do: if in
the diff image there is a clearly real transient, the scanner would not need
to study in detail the srch and tmpl images.

Figure 8A and B show correctly classified “bogus” transients. In Figure 8A, a
“bogus” produced likely by a moving object displaying a classical dipole, the
majority of the important pixels are located in the tmpl
($I_{\mathrm{tmpl}}\sim 0.54$), concentrated around the location of the
central source and the location where its “ghost” image is (the coordinates
corresponding to the location the bright patch of pixels in the diff, but in
the tmpl portion of the composite image). In Figure 8B, there is no central
source and the detection is triggered by an image artifact. The important
pixels are found in all three image segments and are spread around a large
area of each image: the model has inspected the image in its entirety to
decide the classification. Following our consideration about the similarity
between the CNN and human decision process, we speculate that here too the CNN
mimics closely what a human scanner would do: because there is not a clear
central source (a “real” object) in the diff, the scanner would analyze the
srch and tmpl images to extract more information from the context to enable a
classification. However, it should be noted that no quantitative studies of
the features the human scanners use to classify transients has been done, thus
this remains simply an intriguing suggestion.

For the case of the noDIA model, the expectation was less clear: both the srch
and the tmpl images are necessary to “reconstruct” the information contained
in the diff, and while the pixels overlapping with the central transient are
obviously expected to be important, the pixels that surround it are necessary
to essentially reproduce the scaling and PSF-matching operations between tmpl
and srch that the DIA performs. Accordingly, the saliency maps presented in
Figure 9 are more difficult to interpret: in all four cases, important pixels
are found all over the composite images.

To explore how the choice of important pixels may depend on the image label
and on the correct classification we report the fraction of images for which
the diff (srch, tmpl) is the dominant source of important pixels within the
confusion matrix in Figure 10. To do this, we use a rough but intuitive
cutoff: if the normalized sum of the saliency pixels in a third of the image
is larger than $\frac{1}{3}$, then we deduce that the model principally used
that component for its decision. For example, where $I_{\mathrm{diff}}>0.33$,
we conclude that the model principally relied on diff to make the RB
classification. With this cut-off we can assess if there are differences in
the model behavior when classifying objects as a function of their labels or
their classification. In all four cases (all combinations of “real” and
“bogus” label and prediction) the concentration of important pixels is largest
in the portion of the image corresponding to the diff in the DIA-based model.
It is however interesting to note that, in order to correctly classify the
“bogus”, the DIA-based model uses the template and search image more heavily
than in all other cases (diff, srch, tmpl, = 61%, 12%, 27% for TN, while
$I_{\it diff}>85\%$ for TP, FP, and FN).

The cut-off method described above does not allow us to distinguish between
cases where multiple sections of the images were used jointly, perhaps with
similar importance, from cases where the model truly only relied on one
section of the image. For that, we take a closer look at the distribution of
saliency values. In Figure 11 we show the distribution of values of the three
metrics defined in Equation 4 for each of the four cases: TP and TN, in shades
of blue, FP and FN, in shades of orange, following the color-scheme adopted in
Figure 5 and Figure 10. For the DIA-based model (the two left-most columns),
for the majority of the 20,000 images in the training set
$I_{\mathrm{diff}}>0.33$, but there is a secondary pick in the
$I_{\mathrm{diff}}$ distribution near $I_{\mathrm{diff}}\sim 0.1$ populated
entirely by TN cases, complementary to a long right tail in the
$I_{\mathrm{tmpl}}$distribution ($I_{\mathrm{tmpl}}>0.4$). This confirms that
the correct classification in the presence of real transients relies on diff,
but tmpl and srch become important to correctly classify bogus transients,
just like we had seen in the exemplary cases in Figure 8. We note also that
the general shape of each distribution ($I_{\mathrm{diff}}$,
$I_{\mathrm{srch}}$, $I_{\mathrm{tmpl}}$) is similar for the TP and TN case
(blue) and for the FP and FN cases (orange).

For the noDIA model the important pixels are concentrated in the tmpl for most
images ($I_{\mathrm{tmpl}}>0.5$) for both correct and incorrect
classifications. This is somewhat counter-intuitive since the tmpl does not
contain the transient itself! However, one may speculate that this is because
the tmpl, a higer quality image, contains more accurate information about the
context in which the transient arises: e.g., if it is located near a galaxy or
not. This information is important to the classification. It is also
interesting to note that for the transients predicted as “real” in both TP and
FP, the fraction of images that leveraged primarily the tmpl is approximately
$2/3$, and for images predicted as “bogus”, in both TN and FN, it is
approximately $3/4$.

To help guide the interpretation of the saliency maps, a few more maps are
plotted in Appendix C. Where we provide 6 examples per each class of the
confusion matrix, for both the DIA-based and the noDIA models.

###  5.2 Computational cost of our models

The computational cost of our models, reported in CPU Node hours in Table 1,
confirms that while training a CNN model for RB can be computationally
expensive, and significantly more so if the diff is not used in input (noDIA),
the model prediction is instantaneous. Using a NN-based platform, the
computational costs are front-loaded. In transient detection this could mean
that the observation-to-transient discovery process is instantaneous, while
computation time can be spent during off-sky hours (principally to build
templates). Furthermore, while the training time is longer for the noDIA model
than for the model that uses the diff in input (DIA-based) the computational
cost of the forward-pass (prediction) scales superlinearly with the size of
the feature set (pixels), so that our noDIA model takes less than half the
time than the DIA-based one to perform the RB task. With a clock time of 0.3
ms per 51$\times$51 pixels postage stamp, predicting over the full DES focal
plane would take $\sim 1$ minute. However, we note that at this stage of our
work we still rely on the DIA in several ways: while the tmpl and srch in
input to noDIA are not PSF matched, this proof of concept is performed on
transients that were detected in diff images, and we leverage the alignment of
tmpl and srch and centering of the postage stamp that arose from the DIA (see
section 6). Conversely, predictions would only need to be done in
correspondence of the sources detected in the tmpl and srch images, and not
for the entire CCD plane.

##  6 Future work and limitation of this work

This work is targeted to the investigation of CNN RB model performance, with
and without diff in input, on a single data set, the same dataset upon which
the development of the random-forest-based autoscan was based (see section 3).
This approach enables a straightforward comparison, but it comports some
limitations.

The labels in our data set come from simulations of SNe (label = 0) and visual
inspection that classifies artifacts and moving objects (label = 1). We
reserve to future work the investigation of the efficacy of the model on
transients of different nature, including quasars (QSOs), strong lensed
systems, Tidal Disruption Events (TDEs), and supernovae of different types.
These transients may have characteristically different associations with the
host galaxy, including preferences for different galaxy types and locations
with respect to the galaxy center, compared to the simulated SNe events in our
training set.

Specifically thinking of Rubin LSST data, an additional source of variability
may be introduced by Differential Chromatic Refraction effects (Abbott et al.,
2018; Richards et al., 2018), or stars with significant proper motion which,
due to exquisite image quality of the Rubin images, would be detectable
effects.

While we demonstrated CNN model’s potential in the detection of transients
without DIA, we did not address the question of completeness as a function of
srch or tmpl depth or the potential for performing accurate photometry without
DIA.

Improvements in the noDIA model performance could be achieved by considering
alternative network architectures, more training epochs, tuned
hyperparameters, etc., and we leave these tasks as future work beyond the
proof-of-concept presented here.

Finally, we note that our models did in fact implicitly leverage some of the
information generated by the DIA even when they did not use the diff image
itself as input. First, the tmpl and srch images are de-warped. Since the
transient alerts are generated from aligned DIA images, the transient source
is always located at the image center in our data. We use postage stamps that
were, however, not PSF matched or scaled to match the template brightness. To
move beyond a proof-of-concept, in future work we will re-train and apply our
models to images whose alignment does not depend on the existence of a
transient.

##  7 Conclusions

In this work, we have demonstrated that high-accuracy models for classifying
true astrophysical transients from artifacts and moving objects (a task
generally known as “real-bogus”) can be built without leveraging the results
of a difference image analysis (DIA) pipeline that constructs a “template-
subtracted” image.

Starting from the Dark Energy Survey dataset that supported the creation of
the well known real-bogus autoscan model (Goldstein et al., 2015), we first
built a CNN-based model, dubbed DIA-based, that uses a sky template (tmpl), a
nightly image (srch), and a template-subtracted version of the nightly image
(diff) that performs real-bogus classification at the level of state-of-the
art models. Our DIA-based model reaches a state-of-the art $97\%$ accuracy in
the bogus classification with an Area-Under-the-Curve of 0.992 and does not
require human decision in the feature engineering or extraction phase.

We then created noDIA, a model which uses only the tmpl and srch images and
can extract information that enables the identification of bogus transients
with $92\%$ accuracy even with a naive model architecture that, to enable a
direct comparison, is purposefully designed to be minimally different from our
high accuracy DIA-based model.

We further investigated what information enables the real-bogus classification
in both the DIA-based and noDIA models and demonstrated that a CNN trained
with the DIA output primarily uses the information in the diff image to make
the final classification, and that the model examines a diff-srch-tmpl image-
set fundamentally differently in the cases where there is a transient, than in
the cases where there is not one. The noDIA model, conversely, takes a more
comprehensive look at both tmpl and srch images, but relies primarily on tmpl
to enable the reconstruction of the information found in the diff.

Implementation of this methodology in future surveys could reduce the time and
computational cost required for classifying transients by entirely omitting
the construction of the difference images.

##  8 Acknowledgment

TAC and FBB acknowledge the support of the LSST Corporation and the Enabling
Science Grant program that partially supported TAC through Grant No. 2021-040,
Universidad Nacional de Colombia and University of Delaware for the Beyond
Research Program (2020) and Summer Research Exchange (2021). MS and HQ were
supported by DOE grant DE-FOA-0002424 and NSF grant AST-2108094.

TAC: Prepared the data, conducted the analysis, created and maintains the
machine learning models, wrote the manuscript. FBB: Advised on data selection,
curation, and preparation, selection and design, revised the manuscript. GD:
Advised on data preparation and model selection and design, revised the
manuscript. MS: Advised on data selection, shared knowledge about the
compilation and curation of the original dataset, revised the manuscript. HQ:
Advised on data preparation and model design, revised the manuscript.

This paper has undergone internal review in the LSST Dark Energy Science
Collaboration. The internal reviewers were: Viviana Acquaviva, Suhail Dhawan,
and Michael Wood-Vasey.

The DESC acknowledges ongoing support from the Institut National de Physique
Nucléaire et de Physique des Particules in France; the Science & Technology
Facilities Council in the United Kingdom; and the Department of Energy, the
National Science Foundation, and the LSST Corporation in the United States.
DESC uses resources of the IN2P3 Computing Center (CC-IN2P3–Lyon/Villeurbanne
- France) funded by the Centre National de la Recherche Scientifique; the
National Energy Research Scientific Computing Center, a DOE Office of Science
User Facility supported by the Office of Science of the U.S. Department of
Energy under Contract No. DE-AC02-05CH11231; STFC DiRAC HPC Facilities, funded
by UK BIS National E-infrastructure capital grants; and the UK particle
physics grid, supported by the GridPP Collaboration. This work was performed
in part under DOE Contract DE-AC02-76SF00515.

This project used public archival data from the Dark Energy Survey (DES).
Funding for the DES Projects has been provided by the U.S. Department of
Energy, the U.S. National Science Foundation, the Ministry of Science and
Education of Spain, the Science and Technology Facilities Council of the
United Kingdom, the Higher Education Funding Council for England, the National
Center for Supercomputing Applications at the University of Illinois at
Urbana–Champaign, the Kavli Institute of Cosmological Physics at the
University of Chicago, the Center for Cosmology and Astro-Particle Physics at
the Ohio State University, the Mitchell Institute for Fundamental Physics and
Astronomy at Texas A&M University, Financiadora de Estudos e Projetos,
Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro,
Conselho Nacional de Desenvolvimento Científico e Tecnológico and the
Ministério da Ciência, Tecnologia e Inovação, the Deutsche
Forschungsgemeinschaft and the Collaborating Institutions in the Dark Energy
Survey. The Collaborating Institutions are Argonne National Laboratory, the
University of California at Santa Cruz, the University of Cambridge, Centro de
Investigaciones Enérgeticas, Medioambientales y Tecnológicas–Madrid, the
University of Chicago, University College London, the DES-Brazil Consortium,
the University of Edinburgh, the Eidgenössische Technische Hochschule (ETH)
Zürich, Fermi National Accelerator Laboratory, the University of Illinois at
Urbana-Champaign, the Institut de Ciències de l’Espai (IEEC/CSIC), the
Institut de Física d’Altes Energies, Lawrence Berkeley National Laboratory,
the Ludwig-Maximilians Universität München and the associated Excellence
Cluster Universe, the University of Michigan, the National Optical Astronomy
Observatory, the University of Nottingham, The Ohio State University, the
OzDES Membership Consortium, the University of Pennsylvania, the University of
Portsmouth, SLAC National Accelerator Laboratory, Stanford University, the
University of Sussex, and Texas A&M University.

Based in part on observations at Cerro Tololo Inter-American Observatory,
National Optical Astronomy Observatory, which is operated by the Association
of Universities for Research in Astronomy (AURA) under a cooperative agreement
with the National Science Foundation.

This research used resources of the National Energy Research Scientific
Computing Center (NERSC), a U.S. Department of Energy Office of Science User
Facility located at Lawrence Berkeley National Laboratory, operated under
Contract No. DE-AC02-05CH1123.

All of our code is written in Python using the following packages:

  * •

TensorFlow (Abadi et al., 2015)

  * •

Keras (Chollet et al., 2015)

  * •

numpy (Harris et al., 2020)

  * •

astropy (Astropy Collaboration et al., 2013, 2018)

  * •

pandas (Wes McKinney, 2010; pandas development team, 2020)

  * •

matplotlib (Hunter, 2007)

  * •

seaborn (Waskom, 2021)

##  9 Data availability

The data underlying this article are available in
https://portal.nersc.gov/project/dessn/autoscan/#tdata, and explain in more
detail in Goldstein et al. (2015).

## References

  * Abadi et al. (2015) Abadi M., et al., 2015, TensorFlow: Large-Scale Machine Learning on Heterogeneous Systems, https://www.tensorflow.org/
  * Abbott et al. (2018) Abbott T. M. C., et al., 2018, The Astrophysical Journal Supplement Series, 239, 18
  * Agarap (2018) Agarap A. F., 2018, CoRR, abs/1803.08375
  * Alard & Lupton (1998) Alard C., Lupton R. H., 1998, The Astrophysical Journal, 503, 325
  * Astropy Collaboration et al. (2013) Astropy Collaboration et al., 2013, A&A, 558, A33
  * Astropy Collaboration et al. (2018) Astropy Collaboration et al., 2018, AJ, 156, 123
  * Bellm et al. (2018) Bellm E. C., et al., 2018, Publications of the Astronomical Society of the Pacific, 131, 018002
  * Bloom et al. (2012) Bloom J., et al., 2012, Publications of the Astronomical Society of the Pacific, 124, 1175
  * Cabrera-Vives et al. (2016) Cabrera-Vives G., Reyes I., Förster F., Estévez P. A., Maureira J.-C., 2016, in 2016 International Joint Conference on Neural Networks (IJCNN). pp 251–258, doi:10.1109/IJCNN.2016.7727206
  * Cabrera-Vives et al. (2017) Cabrera-Vives G., Reyes I., Förster F., Estévez P. A., Maureira J.-C., 2017, The Astrophysical Journal, 836, 97
  * Carrasco-Davis et al. (2021) Carrasco-Davis R., et al., 2021, Alert Classification for the ALeRCE Broker System: The Real-time Stamp Classifier (arXiv:2008.03309)
  * Chollet et al. (2015) Chollet F., et al., 2015, Keras, https://github.com/fchollet/keras
  * Deng et al. (2009) Deng J., Dong W., Socher R., Li L.-J., Li K., Fei-Fei L., 2009, in 2009 IEEE Conference on Computer Vision and Pattern Recognition. pp 248–255, doi:10.1109/CVPR.2009.5206848
  * Dieleman et al. (2015) Dieleman S., Willett K. W., Dambre J., 2015, Monthly Notices of the Royal Astronomical Society, 450, 1441–1459
  * Gabbard et al. (2018) Gabbard H., Williams M., Hayes F., Messenger C., 2018, Phys. Rev. Lett., 120, 141103
  * Gieseke et al. (2017) Gieseke F., et al., 2017, Monthly Notices of the Royal Astronomical Society, 472, 3101–3114
  * Goldstein et al. (2015) Goldstein D. A., et al., 2015, The Astronomical Journal, 150, 82
  * Gunn et al. (2006) Gunn J. E., et al., 2006, The Astronomical Journal, 131, 2332–2359
  * Hambleton et al. (2020) Hambleton K., et al., 2020, Research Notes of the AAS, 4, 40
  * Hanley & McNeil (1983) Hanley J. A., McNeil B. J., 1983, Radiology, 148, 839
  * Harris et al. (2020) Harris C. R., et al., 2020, Nature, 585, 357
  * Hernández-Orallo et al. (2012) Hernández-Orallo J., Flach P., Ferri Ramírez C., 2012, Journal of Machine Learning Research, 13, 2813
  * Ho (1995) Ho T. K., 1995, in Proceedings of 3rd international conference on document analysis and recognition. pp 278–282
  * Hunter (2007) Hunter J. D., 2007, Computing in Science & Engineering, 9, 90
  * Ivezić et al. (2019a) Ivezić Ž., et al., 2019a, ApJ, 873, 111
  * Ivezić et al. (2019b) Ivezić Ž., et al., 2019b, The Astrophysical Journal, 873, 111
  * Kessler et al. (2015) Kessler R., et al., 2015, The Astronomical Journal, 150, 172
  * Kim & Brunner (2016) Kim E. J., Brunner R. J., 2016, Monthly Notices of the Royal Astronomical Society, 464, 4463–4475
  * Krizhevsky et al. (2012) Krizhevsky A., Sutskever I., Hinton G. E., 2012, in Proceedings of the 25th International Conference on Neural Information Processing Systems - Volume 1. NIPS’12. Curran Associates Inc., Red Hook, NY, USA, p. 1097–1105
  * LeCun et al. (1989) LeCun Y., et al., 1989, Connectionism in perspective, 19, 143
  * LeCun et al. (2015) LeCun Y., Bengio Y., Hinton G., 2015, Nature, 521, 436–444
  * LeNail (2019) LeNail A., 2019, Journal of Open Source Software, 4, 747
  * Liu et al. (2019) Liu W., et al., 2019, Research in Astronomy and Astrophysics, 19, 042
  * Mahabal et al. (2019) Mahabal A., et al., 2019, PASP, 131, 038002
  * Mong et al. (2020) Mong Y.-L., et al., 2020, Monthly Notices of the Royal Astronomical Society, 499, 6009–6017
  * Nielsen (2015) Nielsen M. A., 2015, Neural Network and Deep Learning. Determination Press
  * Richards et al. (2018) Richards G., Peters C., Martin B., UW J. B. K., UW J. P., UW I. S., Bauer F. E., 2018
  * Simonyan et al. (2014) Simonyan K., Vedaldi A., Zisserman A., 2014, Deep Inside Convolutional Networks: Visualising Image Classification Models and Saliency Maps (arXiv:1312.6034)
  * Sánchez et al. (2019) Sánchez B., et al., 2019, Astronomy and Computing, 28, 100284
  * The Dark Energy Survey Collaboration (2005) The Dark Energy Survey Collaboration 2005, arXiv e-prints, pp astro–ph/0510346
  * Tomaney & Crotts (1996) Tomaney A. B., Crotts A. P. S., 1996, The Astronomical Journal, 112, 2872
  * Wardęga et al. (2020) Wardęga K., Zadrożny A., Beroiz M., Camuccio R., Díaz M. C., 2020, Detecting optical transients using artificial neural networks and reference images from different surveys (arXiv:2009.14614)
  * Waskom (2021) Waskom M. L., 2021, Journal of Open Source Software, 6, 3021
  * Wes McKinney (2010) Wes McKinney 2010, in Stéfan van der Walt Jarrod Millman eds, Proceedings of the 9th Python in Science Conference. pp 56 – 61, doi:10.25080/Majora-92bf1922-00a
  * Zackay & Ofek (2017) Zackay B., Ofek E. O., 2017, The Astrophysical Journal, 836, 188
  * pandas development team (2020) pandas development team T., 2020, pandas-dev/pandas: Pandas, doi:10.5281/zenodo.3509134, https://doi.org/10.5281/zenodo.3509134

##  Appendix A Scaling astrophysical images for input to Neural Networks

To visualize and compare the behaviour of the flux of the srch and tmpl
images, the value distribution for the four transients presented in Figure 1
is shown as a violin plot in Figure 12. The first transient on the left is
labeled as “bogus”, the srch distribution pixel values are in general greater,
positive and non-zero center than the values for tmpl, the same behaviour is
observed to the last transient on the right, but this one is labeled as
“real”. This same comparison can be applied to the transients plotted on the
middle of Figure 12, both show similar distributions, however one is “bogus”
and the other is “real”. For the four transients presented the pixel values
have long tails, outside the $\mu\pm 3\sigma$ values. The scaling of the srch
and tmpl images for the four transients, according to the description given in
subsection 3.1 is visualized in Figure 13.

Figure 12: Violin plot of pixel values inside the $\mu\pm 3\sigma$ interval
for srch and temp image of the data shown in Figure 1 (before normalization).
These values were scaled between $0$ and $1$. Left and right curve correspond
to srch and temp distribution of values respectively. The distribution of
values correspond to the same shown in Figure 2 for srch and temp, the green
and purple colors has its correspondence to Figure 2. Figure 13: Violin plot
of pixel values for srch and temp image of the data shown in Figure 1 (after
normalization). Values between the vertical lines in $0$ and $1$ were mapped
to the values in Figure 12. Negative values and above $1$ correspond to the
scaled values outside the $3\sigma$ clip.

##  Appendix B Detailed architecture of DIA-based and noDIA models

We designed a network of 12 layers using `tensorflow.keras`in python, for the
DIA-based case as follows:

⬇

layer1 = keras.layers.Conv2D(16, kernel_size=(5, 5), padding="valid",
activation="relu",

input_shape=(51,153,1))

layer2 = keras.layers.MaxPooling2D((2, 2), strides=2)

layer3 = keras.layers.Dropout(0.4)

layer4 = keras.layers.Conv2D(32, kernel_size=(5, 5), padding="valid",
activation="relu")

layer5 = keras.layers.MaxPooling2D((2, 2), strides=2)

layer6 = keras.layers.Dropout(0.4)

layer7 = keras.layers.Conv2D(64, kernel_size=(5, 5), padding="valid",
activation="relu")

layer8 = keras.layers.MaxPooling2D((2, 2), strides=2)

layer9 = keras.layers.Dropout(0.4)

layer10 = keras.layers.Flatten()

layer11 = keras.layers.Dense(32, activation="relu")

layer12 = keras.layers.Dense(2, activation="softmax")

and a network of 11 layers using `tensorflow.keras`in python, for the noDIA
case as follows:

⬇

layer1 = keras.layers.Conv2D(1, kernel_size=(7, 7), padding="valid",
activation="relu",

input_shape=(51,102,1))

layer2 = keras.layers.MaxPooling2D((2, 2), strides=2)

layer3 = keras.layers.Conv2D(16, kernel_size=(3, 3), padding="valid",
activation="relu")

layer4 = keras.layers.MaxPooling2D((2, 2), strides=2)

layer5 = keras.layers.Dropout(0.4)

layer6 = keras.layers.Conv2D(32, kernel_size=(3, 3), padding="valid",
activation="relu")

layer7 = keras.layers.MaxPooling2D((2, 2), strides=2)

layer8 = keras.layers.Dropout(0.4)

layer9 = keras.layers.Flatten()

layer10 = keras.layers.Dense(32, activation="relu")

layer11 = keras.layers.Dense(2, activation="softmax")

To compile the models we used the optimizer `SGD` (stochastic gradient
descent) with a `learning_rate=0.01` the `"sparse_categorical_crossentropy"`
loss function and the `"accuracy"` metric. Our models are available in a
dedicated GitHub repository 888https://github.com/taceroc/DIA_noDIA.

##  Appendix C Saliency maps for various transients

We include a series of saliency maps in the following 8 figures. Several
interesting behavioral patterns can be observed. The figures are organized by
model and by classification as follows: TN —Figure 15 and Figure 15—, TP
—Figure 16 and Figure 17—, FN Figure 18 and Figure 19–––, and FP —Figure 20
and Figure 21.

In Figure 15 and Figure 15 we include the transient’s contours overplotted
onto the saliency maps for the DIA-based and noDIA model respectively to guide
the reader’s eye. Notice the offset of the DIA-based model’s focus in Figure
15A with respect to the transient: the model is principally inspecting the
diff and tmpl at the location corresponding to the bright patch in the diff.
In Figure 15C a correctly predicted “bogus”, characterized by a very tight
dipole probably arising from a poorly centered DIA, the saliency map inspect
the tmpl and even seem to reproduce traditional aperture photometry, with
pixel values measured at the core of the transients, and in the sky
surrounding the transients, but away from the tails of the bright pixels
cluster. The behavior is principally different for the same transients when
they are inspected, and also correctly predicted, by the noDIA model (Figure
15): the focus of the model is in all three cases away from the transients,
and shifted to the surrounding: the model is learning the transients’ context
and extracting information to enable the comparison of tmpl and diff
(essentially, to enable the image differencing). This behavior is generally
seen throughout all examples in the following figures. In addition, for each
figure we highlight potential reasons for failed predictions, and potential
inaccuracies in the labeling that may lead to an artificial lowering of our
measured accuracy.

Figure 14: Transients (diff-srch-tmpl) and their respective saliency map for
DIA-based model True Negative predictions (correctly identified “bogus”).
Contour plot of light intensity from the original images are overplotted,
delineating the bright sources. In panel A: $I_{\it diff}\sim 0.4$; B: $I_{\it
diff}\sim 0.69$ and C $I_{\it tmpl}\sim 0.42$. In panel C, the saliency for
the tmpl portion of the image looks strikingly similar to the selection of
pixels that is made in aperture photometry, with pixel values considered in
the core of the source, ignored in a region immediately around the source, and
again considered farther out to calculate the source’s background. This figure
is further discussed in Appendix C.

Figure 15: The same transients (srch-tmpl) shown in Figure 15 and their
respective saliency maps for noDIA model True Negatives (correctly identified
“bogus”). Important pixels are found at nearly all locations in the image,
rather than in a small region around the center. The model needs to learn
properties of the image at large to enable a comparison of the tmpl and diff.
This figure is further discussed in Appendix C.

Figure 16: Transients (diff-srch-tmpl) and their respective saliency map for
DIA-based model True Positives (correctly identified “real” astrophysical
transients). The important pixels are generally found in the diff portion of
the image for the DIA-based model, as discussed in section 5, but there are
exceptions: here we show several cases of True Positives (real transients)
classifications where the component of the image that was principally
leveraged by the model was the diff (the leftmost third) but two cases where
the noDIA relied principally on the srch and/or tmpl (bottom two panels).

Figure 17: Transients (srch-tmpl) and their respective saliency map for the
noDIA model True Positives (correctly identified “real” astrophysical
transients). Important pixels are found everywhere in the image, as the CNN
learns how to compare the diff and tmpl taking a synoptic look at the
properties of each image component.

Figure 18: Transients (diff-srch-tmpl) and their respective saliency map for
DIA-based model False Negatives (real transients identified as “bogus”). We
remind the reader that the labels are inherited from Goldstein et al. (2015)
and cannot be verified. Some level of label inaccuracy is expected. The “real”
transients in this dataset are implanted supernovae onto real DES images.
However, in this collection, several transient display DIA inaccuracies (row
1, 2, 4, 6 show “dipoles”, see subsection 2.1) that likely lead to the
incorrect classification. Two very low signal-to-noise detections are missed
(row 3 and 5) by our model. Important pixels are more commonly found in diff
portion of the image. In the srch saliency maps we see again that the core of
the central source is used in the classification, as well as pixels that
surround the source, but these two sets of important pixels are separate by it
by a gap, again reminiscent of the typical aperture photometry technique (top
two panels).

Figure 19:  Transients (srch-tmpl) and their respective saliency map for the
noDIA model False Negatives (astrophysical sources classified as “bogus”). In
all cases but row 5 it is not clear why the classification fails. In row 5,
another source dominates the image scaling (and pre-processing) reducing the
visibility of the transient (that is completely missed by human inspection).
We remind the reader that the “real” transients in this dataset are implanted
supernovae onto real DES images. Important pixels are found everywhere in the
image, as the CNN learns how to compare the diff and tmpl taking a synoptic
look at the properties of each image component.

Figure 20: Transients (diff-srch-tmpl) and their respective saliency map for
DIA-based model False Positives (“bogus” predicted as “real”). We remind the
reader again that the labels are inherited from Goldstein et al. (2015) and
cannot be verified. Some level of label inaccuracy is expected. Bogus
transients were labeled by human scanners among astrophysical images with
detection. However, in this collection, we cannot verify the nature of the
transient and we argue that in the cases presented here there is no obvious
evidence of its “bogus” nature. Important pixels are most commonly found in
the diff, but in tmpl and srch we see again the CNN analyzes the central
source and its surrounding, but avoiding the tail of the central source, in a
way similar to traditional aperture photometry techniques.

Figure 21: Transients (srch-tmpl) and their respective saliency map for the
noDIA model False Positives (“bogus” classified as “real”). Important pixels
are found everywhere in the image, as the CNN learns how to compare the diff
and tmpl taking a synoptic look at the properties of each image component.