LensExtractor: A Convolutional Neural Network in Search of Strong Gravitational Lenses

TL;DR

LensExtractor is a convolutional neural network designed to efficiently identify strong gravitational lenses in large astronomical surveys, significantly reducing manual inspection efforts.

Contribution

We developed a novel CNN-based algorithm for gravitational lens detection that outperforms classical methods in efficiency and scalability.

Findings

Successfully identified lensing candidates in COSMOS field data.

Demonstrated high accuracy and computational efficiency of LensExtractor.

Suitable for application to future large-scale surveys like LSST and WFIRST.

Abstract

In this work, we present our classification algorithm to identify strong gravitational lenses from wide-area surveys using machine learning convolutional neural network; LensExtractor. We train and test the algorithm using a wide variety of strong gravitational lens configurations from simulations of lensing events. Images are processed through multiple convolutional layers which extract feature maps necessary to assign a lens probability to each image. LensExtractor provides a ranking scheme for all sources which could be used to identify potential gravitational lens candidates significantly reducing the number of images that have to be visually inspected. We further apply our algorithm to the \textit{HST}/ACS i-band observations of the COSMOS field and present our sample of identified lensing candidates. The developed machine learning algorithm is much more computationally efficient than classical lens identification algorithms and is ideal for discovering such events across wide areas from current and future surveys such as LSST and WFIRST.

Figures

Tables

Table 1. Table 1 Catalog of displayed lenses in Figure 9 .

Lens	Object ID	Right Ascension	Declination	Einstein Radius	Magnitude†
		(deg)	(deg)	(arcsec)	(AB)
$1^{\ast }$	COSMOS0108+4029	150.2849	2.6749	1.51	19.08
$2$	COSMOS5844+3753	149.6860	1.6315	1.53	21.65
$3^{\ast }$	COSMOS5831+4334	149.6298	1.7263	1.17	20.76
$4$	COSMOS0142+5447	150.4285	1.9133	1.49	18.67
$5^{\ast }$	COSMOS0124+5121	150.3522	1.8559	0.79	22.31
$6^{\ast }$	COSMOS0047+5023	150.1986	1.8398	1.71	20.64
$7$	COSMOS0027+0513	150.1140	2.0871	1.80	18.87
$8^{\ast }$	COSMOS0208+1422	150.5355	2.2396	1.54	20.13
$9^{\ast }$	COSMOS0056+1226	150.2366	2.2072	2.04	18.78
$10$	COSMOS0237+2652	150.6560	2.4478	1.88	20.81
${11}^{\ast }$	COSMOS0012+2015	150.0526	2.3377	0.92	19.37
$12$	COSMOS5918+1911	149.8272	2.3198	1.89	19.86
$13$	COSMOS0211+2955	150.5486	2.4986	1.91	20.82
$14$	COSMOS0157+3510	150.4879	2.5863	1.20	20.68
${15}^{\ast }$	COSMOS0018+3845	150.0770	2.6460	1.26	23.62
$16$	COSMOS5844+3753	149.6860	1.6315	1.45	21.65
$17$	COSMOS0247+5601	150.6974	1.9337	0.82	20.24
$18$	COSMOS0229+1441	150.6260	2.2448	0.51	26.67
$19$	COSMOS5954+1128	149.9791	2.1913	0.97	21.19
${20}^{\ast }$	COSMOS0038+4133	150.1595	2.6927	0.55	20.51
${21}^{\ast }$	COSMOS5921+0638	149.8407	2.1106	0.72	20.43

Equations

$$a_i^l=f(\sum _j{a}_j^{l-1}{w}_{j\to i}^l+b_i^l).$$

$$\leavevmode\resizebox{399.2946pt}{}{$C(\boldsymbol{x},\boldsymbol{y})=\sum_{j=\text{unlensed},\quad\\ \text{lensed}}{y_{j}\ln a_{j}^{L}+(1-y_{j})\ln(1-a_{j}^{L})}$}.$$

$$p_{ij}\leftarrow p_{ij}-\frac{{\sum }_{ij}{p}_{ij}}{i_{max}{j}_{max}}$$

$$p_{ij}\leftarrow \frac{p_{ij}}{\sqrt{{\sum }_{ij}{p}_{ij}^2}}$$

$$p_{ij}\leftarrow (1-e^{-((i-100{)}^2+(j-100{)}^2)/1000})p_{ij}$$

$$p_{ij}\leftarrow (p_{ij}+1.01\min (\mathbf{p}){)}^{0.15}$$

$$p_{ij}\leftarrow \min (p_{ij},0)$$

$${\sigma }_c(\mathbf{Z})=\frac{e^{Z_c}}{e^{Z_{\text{unlensed}}}+e^{Z_{\text{lensed}}}}.$$

$$\text{accuracy}\equiv \frac{N(\text{lensed}|p_{\text{lens}}>0.5)}{N(\text{lensed})+N(\text{unlensed})}$$

Full text

††journal: Submitted to the Astrophysical Journal

1. Introduction

Gravitational lensing, a prediction of Einstein’s general theory of relativity, is a very powerful tool in cosmological studies. It has been used extensively to understand various aspects of galaxy formation and evolution (e.g. Refsdal & Bondi (1964); Blandford & Narayan (1992); Nayyeri et al. (2016); Postman et al. (2012); Atek et al. (2015)). This involves accurate cosmological parameter estimation (Treu 2010), studies of dark matter distribution from weak gravitational lensing events (Kaiser & Squires 1993; Velander et al. 2014), black-hole physics (Peng et al. 2006) and searches for the most distant galaxies (Coe et al. 2012; Oesch et al. 2015), among others.

One of the main goals of observational cosmology is to constrain the main cosmological parameters that dictate the evolution of the Universe (Tegmark et al. 2004; Komatsu et al. 2009; Weinberg et al. 2013). Strong gravitational lensing has been utilized over the past few years to estimate and contain these cosmological parameters (Broadhurst et al. 2005; Suyu et al. 2013, 2014; Goobar 2016; Agnello et al. 2017; More et al. 2017). This is achieved through accurate lens modeling of such events and comparing the model predictions with observations (such as with observations of lensing induced time delays (Eigenbrod et al. 2005; Treu 2010; Suyu et al. 2014; Rodney et al. 2016; Treu & Marshall 2016). In a recent study, for example, Suyu et al. (2013) used combined WMAP, Keck and HST data on gravitational time delays in two lensed sources to constrain the Hubble constant within $4\%$ in a $\rm\Lambda CDM$ cosmological framework.

One of the key aspects of gravitational lensing is its use as natural telescopes through boosting the observed signal and increasing the spatial resolution (Treu 2010). This is quite advantageous in searches for distant and/or faint objects at moderate observing costs and has been utilized extensively in various surveys in searches for such objects, the identification of which would not have been possible without it (Bolton et al. 2006; Heymans et al. 2012). Given that the number of identified lenses for different classes of galaxies rises sufficiently due to better lens finding algorithms, by modeling the lenses and using spectra stacking techniques, we can better understand the physical and chemical composition of farther and fainter galaxies which in turn would excel our understanding of galaxy evolution (Wilson et al. 2017; Timmons et al. 2016). In the past few years, deep diffraction limited observations have also taken advantage of gravitational lensing to extend the faint end of the luminosity function of galaxies by a few orders of magnitude (Atek et al. 2015) to produce the deepest images of the sky ever taken across multiple bands. Strong gravitational lensing events have been observed extensively in such surveys as galaxy-galaxy lensing in field surveys such as the Cosmic Assembly Near-infrared Deep Extragalactic Legacy Survey (CANDELS) (Grogin et al. 2011; Koekemoer et al. 2011) and the Cosmological Evolution Survey (COSMOS) (Scoville et al. 2007; Capak et al. 2007) or as cluster lensing from observations of nearby massive clusters (Postman et al. 2012; Treu et al. 2015; Lotz et al. 2017) with Hubble Space Telescope. These yield identification of the first generations of galaxies at $z>3$ (and out to $z\sim 11$; Oesch et al. 2015) to study galaxy formation and evolution at the epoch of re-ionization. This was, in fact, one of the main motivations behind Hubble cluster lensing studies such as CLASH and Frontier Fields (Postman et al. 2012; Lotz et al. 2017). These magnification provided by the strong lensing could potentially raise the observed flux by as much as a factor of hundred (depending on configuration). The power of gravitational lensing could also be used in the detection of low surface brightness emission from extended objects such as far-infrared and radio emissions from dust and molecular gas at $z\sim 2-3$. This is indeed one of the main techniques in observing these systems even in the era of powerful mm/radio telescopes such as ALMA.

Strongly lensed galaxies are normally targeted and identified from dedicated surveys (Bolton et al. 2006). Traditionally these lens identifications are either catalog-based, in which lensing events are identified by looking for objects in a lensing configuration, or pixel based, with the search starting from a set of pixels. These lensing searches are normally computationally challenging in that individual pixels are constantly compared with adjacent ones and they could be biased towards a given population and/or brighter objects. Recent far-infrared wide area observations (such as with Herschel) significantly advanced these searches for lensed galaxies by adopting a simple efficient selection technique of lensed candidates through observations of excessive flux in the far infrared (as an indication of strong lensing events supported by number count distributions; Nayyeri et al. 2016; Wardlow et al. 2012). However such surveys are also biased towards populations of red dusty star-forming galaxies (missing any blue lenses) and are not always available across the full sky (the Herschel surveys that were targeted had $\rm\sim 0.2-0.4\,deg^{-2}$ lensing events, much lower than expected from optical surveys). Given that tests of cosmological models require simple unbiased selection function, it is important to have a complete unbiased catalog of lensing events.

We have entered the era of big data astronomy. Sky surveys such as the LSST, Euclid, and WFIRST will produce more imaging data than humans can ever analyze by eye. The challenges of designing such surveys are no longer merely instrumentational, but they also demand powerful data analysis and classification tools that can identify astronomical objects autonomously. Fortunately, computer vision has drastically improved in the last couple of years to make autonomous astronomy possible. The past couple of years has been the most exciting era in the field of machine learning (ML). Researchers from both the public and the private sectors have achieved landmarks in developing image recognition/classification techniques. One of the most exciting recent events in the ML community was the release of TensorFlow by Google, a parallel processing platform designed for development of fast deep learning algorithms (Abadi et al. 2016). Packages like TensorFlow, Caffe, and others have enabled researchers to develop very complex and fast classification algorithms. Among these deep learning programs, ConvNets have deservingly received a lot of attention in many fields of science and industry in the past few years (Krizhevsky et al. 2012). Complex ConvNets such as GoogleNet and AlexNet, which are publicly available, have achieved superhuman performance on the task of image classification. Google’s TensorFlow has made it possible to easily develop parallelized deep learning algorithms which if integrated with Google’s Tensor Processing Units (TPUs), could address the data mining challenges in the field of astronomy. The field of astronomy and observational astrophysics should take advantage of these new image classification algorithms. For instance, our team has been working on developing LensFlow, a ConvNet that can be used to search for strong gravitational lenses. Our work will be publicly available on Github after publication.

Non-machine learning computer algorithms have been previously used for finding gravitational lenses (Alard 2006; More et al. 2012). For instance, More et al. (2012) use two algorithms called RingFinder and ArcFinder. The former uses color information and the latter detects arc-like pixels. Initially, ArcFinder polishes the images by convolving a smoothing kernel. For each pixel, an estimator of elongation is calculated by taking the ratio between the sum of the flux of a few pixels along the horizontal line and the maximum value of a few nearby pixels along the vertical line which pass through the pixel in hand. This process is repeated for all pixels and those with smaller than an specified elongation threshold are set to zero to create a sharp arc map. An arc map that satisfies thresholds on the arc properties such as the size and surface brightness will be selected as an arc candidate for further visual inspection. Unlike deep learning, such classical algorithms are not easily parallelized and they can be computationally more expensive depending on the deepness of the ConvNets used. They also require threshold tuning which may cause insensitivity to smaller arcs. However, they do not require a massive training dataset. Even though more challenging, creating a large dataset could eliminate biases toward certain arc-lens morphologies. As we will discuss in Section 4, these algorithms do suffer from the same lens contaminants as ours. The hope is that machine has the capacity to improve its performance with better training datasets and improved architecture while remaining computationally efficient but such improvements are not trivial regarding classical algorithms.

Other researchers (Petrillo et al. 2017; Jacobs et al. 2017; Lanusse et al. 2017) also find deep learning a suitable solution for finding gravitational lenses. Lanusse et al. (2017) use residual ConvNets with 46 layers. Residual ConvNets are modified ConvNet that do not suffer from layer saturation as ordinary ConvNets do. After adding more than 50 layers, the accuracy of ordinary ConvNets no longer improves and the training becomes more challenging. He et al. (2016) were able to overcome this issue by providing residual maps in between layers, which has been employed by Lanusse et al. (2017). They have simulated LSST mock observations in a single band and have trained and tested their network on these images. Jacobs et al. (2017) have trained their ConvNet using multiple color bands and have applied it to Canada-France-Hawaii Telescope Legacy Survey. Petrillo et al. (2017) have searched for lenses in Kilo Degree Survey by training their ConvNet on cataloged luminous red galaxies.

In our independently developed work, we focus on the morphology of the lenses and only rely on one color band, similar to Petrillo et al. (2017) and Lanusse et al. (2017). Our lens simulation method is very similar to Petrillo et al. (2017) where we both merge simulated arcs with real images of galaxies to preserve the complexity of the physical data. In contrast to others, we do not discriminate against different sources found in the COSMOS field. Artifacts, starts, and other sources have been included in our training dataset so LensFlow can be directly applied to fields without a need for a catalog with galaxy type information. The deepness of our ConvNet is comparable to Petrillo et al. (2017) and Jacobs et al. (2017) but it is shallower than Lanusse et al. (2017). As mentioned in Jacobs et al. (2017), the morphology of lenses are much simpler than the morphology of daily objects and human faces which extremely deep ConvNets are developed for. However, the cost to performance ratio of ConvNets with varying deepness has not been studied yet. The effectiveness of deeper ConvNets cannot be compared between ours (and Petrillo et al. 2017 ) and Lanusse et al. (2017) since they have not applied their algorithm to physical data. However, they have studied the change in the performance of their ConvNet by varying the Einstein radii and signal-to-noise ratio of their lenses.

This paper is organized as follow. In Section 2, we will explain the principal concepts underlying neural networks, supervised learning, and ConvNets. A supervised learning algorithm requires a large sample of labeled images, known as the training dataset. In Section 3.1, we will discuss the procedure we have taken to create our training and testing datasets. Before feeding the images to a ConvNet, they must be normalized and should be enhanced. The details of these methods are discussed in Section 3.2. Section 3.3 will lay out the architecture of LensFlow and Section 3.4 will illustrate LensFlow’s performance on the testing dataset. We will conclude by sharing and analyzing the scan results of the COSMOS field in Section 4. Throughout this paper, we assume a standard cosmology with $H_{0}=70\,\text{kms}^{-1}\text{Mpc}^{-1}$, $\Omega_{m}=0.3$ and $\Omega_{\Lambda}=0.7$. Magnitudes are in the AB system where $\text{m}_{\rm AB}=23.9-2.5\times\text{log}(f_{\nu}/1\mu\text{Jy})$ (Oke & Gunn 1983).

2. Deep Learning Algorithms

Artificial neural networks are inspired by biological neurons. Just like biological neurons, artificial neurons receive input signals and send out an output signal to other neurons (see Figure 1). The synaptic connections between neurons are known as weights and the output of a neuron is know as its activation. To reduce the computational time and simplify neural network models, neurons are placed in consecutive layers rather than having a connection with every other neuron. Neurons from one layer cannot talk to each other or to the neurons in arbitrary layers; they may only send their signal to the neurons in the succeeding layer. A neuron receives the weighted sum of the activation of all the neurons in the previous layer, adds an internal parameter known as the bias and maps this sum to a value computed by an activation function (e.g. sigmoid, hyperbolic tangent, rectilinear, softmax). This model can be stated mathematically by the following equation:

[TABLE]

Here, $a_{i}^{l}$ is the activation of the neuron in hand (i.e. the $i$’th neuron in the $l$’th layer), $f$ is the activation function of this neuron, $a^{l-1}_{j}$ is the activation of the neuron $j$ in layer $l-1$ (the previous layer), $w_{j\rightarrow i}^{l}$ is the synaptic weight connecting $i$’th neuron in layer $l$ to the $j$’th neuron in layer $l-1$, and $b_{i}^{l}$ is the bias of the neuron to adjust its activation sensitivity. The first layer, i.e. the input layer, in a deep learning neural net acts as a sensory layer, analogous to the retina. As it gets analyzed, the information from the input layer travels through multiple layers until it reaches the final layer called the classification layer. Each class of images corresponds to a classifying neuron. In our case, we have a neuron corresponding to unlensed and another to lensed images. The neuron with the highest output determines which class an input image is placed in.

A neural net learns how to classify images by adjusting the weights between its neurons and the biases within them, having one goal in mind: minimizing the loss function $C(\boldsymbol{x},\boldsymbol{y})$. The loss function, sometimes called the cost function, can take many forms but it has to captures the misfiring of the classification neurons, i.e. the deviation between the target class versus the predicted class. This is why such algorithms are known as supervised learning algorithms, in contrast to unsupervised techniques. A common choice for the loss function is the cross-entropy loss function with the following form (Nielsen 2016):

[TABLE]

$a_{j}^{L}$ is the activation of neurons in the final (classifying) layer. $\boldsymbol{x}$ is the input data in the vector form and $\boldsymbol{y}$ represents the desired activations of the two classifying neurons. Of course, this function depends on the architecture of the neural net, weights, and biases, but they have not been expressed explicitly. As an example, if an image is a lens, its target output has to be $(0.0,1.0)$, meaning the activation of the unlensed neuron should be zero and the activation of the lensed neuron should be unity. During the training process of a neural net, images from a training dataset are presented to the network and the weights and the biases are adjusted to minimize the loss function for those images. The parameter space is massive and a change in one of the parameters of a neuron will affect the activation of a series of neurons in other layers. The first challenge is solved by minimizing algorithms such as the stochastic gradient descent (SGD) and the second one is solved via back-propagation. We won’t go in the details of these two techniques, but it worths emphasizing on the stochasticity of SGD. Stochasticity refers to randomly selecting images from the training set and bundling them in one batch. The loss function for the batch is the average of the loss function for individual images. It is very important to use a batch rather than training the neural net with one image at a time. Loosely speaking, using a batch would drastically improve classification accuracy of the network since neurons learn the features in images rather than memorizing examples.

ConvNets are a class of neural networks with multiple convolutional layers. A convolutional layer consists of a set of convolving neurons (on the order of 10 neurons) which can be connected to a small rectangular region of an image. The set of weights of a convolving neuron is known as a filter and are subject to change as the ConvNet learns. A filter scans an entire image by striding (convolving with specified steps) over the image and assembling its output into an image knows as a feature map. Feature maps contain information such as texture and edges. See Figure 2 as an example of a set of filters in a LensFlow convolutional layer. Two extracted feature maps of a physical lens (not simulated) have been shown in Figure 3.

3. Methodology

3.1. Training and Testing Datasets

Neural networks learn to classify data by learning from examples, referred to as the training dataset. A training dataset usually contains a few thousands of pre-classified images. Our training dataset contains $15,000$ unlensed sources and $15,000$ lensed sources. To make this dataset, we used SExtractor to obtain a catalog of sources from a few high-quality central tiles as well as a few low-quality edge tiles to include noisy images and artifacts. A cutout of $200\times 200$ pixels was made around the identified sources and stored individually, labeled as unlensed. After training and finalizing the architecture of the ConvNet, using a source extraction software is not necessary since an entire tile can be divided into smaller tiles and scanned to identify locations with high lensing probability.

Creating lenses are more challenging and more involved. To create these lenses, we had two options; we could either artificially boost up the number of know lenses or simulate them. We tried both methods. For the first method, we used 17 out of 18 lenses discussed in Faure et al. (2011) paper.111One of the lens candidates from Faure et al. (2011) is more similar to a tidal interaction. Since the number of known lenses is very limited, we formed $200\times 200$-pixel cutouts by applying all combination of the following transformations: (a) Rotating from 0 to 360 degrees in steps of 30 degrees, (b) Shifting the center by [math], $\pm 10$, and $\pm 25$ pixels in horizontal and vertical directions, (c) Scaling by 0.8, 1, 1.1 and (d) Adding a small Gaussian noise after transformation. This combination results in 900 transformed images per lens, boosting up our training dataset to fifteen thousand.

However, the training dataset using the transformation method does not contain a wide range of lenses and lenses with different structures might be missed. For this reason, we have also created a simulated training dataset using LensTool (Kneib et al. 1996; Jullo et al. 2007; Jullo & Kneib 2009). LensTool receives input parameters of the lensing model and generates a lensed image of the background galaxy without a foreground lensing source. We will refer to these as arcs. Lens model parameters such as redshift of the background and foreground galaxies, their ellipticity, orientation, and their relative positions were randomly changed to create a wider range of arcs, from complete Einstein rings with different radii to arcs with different orientations and sizes. The code for generating these arcs which includes the range of model parameters will be included in our online distribution. After generating a wide range of arcs, they were merged with raw images. Both images, the raw sources, and the simulated arcs, were normalized and the pixels of the arc images were multiplied by a random number between $0.3$ and $0.8$ so a range of relative arc to foreground luminosity would be captured in the training set. For a limited number of generated lenses, other ranges were used to create extremely faint and extremely bright lenses. $85\%$ of the lenses were generated by selecting elliptical sources while the remaining were sources randomly selected from tile 55. A bias like this is necessary since the foreground of the know COSMOS lenses are all elliptical. After examining the simulated lenses by eye and eliminating lenses with unnatural geometry, an eightfold transformation was performed to increase the number of lenses up to fifteen thousand. These transformations come from 8 elements of the symmetry group of the square, namely: $0\,^{\circ}$, $90^{\circ}$, $180^{\circ}$ and $270^{\circ}$ rotations, and horizontal, vertical, diagonal and anti-diagonal mirroring.

At the end, 500 lenses and 500 raw images were removed from the training set and placed in a separate set called the testing dataset to measure the performance of the ConvNet as the network gets trained on the training dataset.

3.2. Image Normalization and Enhancement

Before inputting the images to our neural net, we normalize and enhance them. For deep learning purposes, there are different methods of image normalization to choose from, but not normalizing is not a choice. This is due to the fact that raw images come in a wide range of values. However, a limited number of neurons are not capable of handling such variations. On top of that, the goal of a ConvNet is to learn geometric features rather than the statistical properties of the pixels, which could be unwillingly introduced while constructing the training dataset. To address these two issues, all input images must be normalized. Our method of normalization consists of two steps. In the first step, pixel values are shifted by their average so their new mean would be zero (Equation 3). In the second step, pixel values in an image are divided by their standard deviation (Equation 4).

[TABLE]

[TABLE]

$p_{ij}$ stands for the value of the pixel at row $i$ and column $j$. $i_{max}$ and $j_{max}$ indicate the total number of rows and columns, respectively. In order to improve our classification accuracy, we made the arcs stand out by apply an anti-Gaussian function (Equation 5) on each image. This function will attenuate the central pixels which belong to the foreground galaxy (see Figure 4). The width of this function is about a quarter of the cutout width.

[TABLE]

At this step, the images are ready for contrast adjustments which is done using the gamma correction. The pixels must be shifted up in order to eliminate zero or negative pixel values.

[TABLE]

The normalization steps are applied again and at the end, the negative pixels are dropped.

[TABLE]

This improves the image quality by removing low-intensity pixels which mostly include noise. Now images may be inputted to the neural net for training or for classification of new data. See Figure 4 for an example of image normalization and enhancement.

3.3. Architecture of LensFlow ConvNet

The architecture of the data determines the dimensionality of the ConvNet layers. We use $200\times 200\times 1$ images where $1$ indicate the number of color channels222In this paper and in our code, we have adapted the $N\times H\times W\times C$ format from TensorFlow where $N$, $H$, $W$, and $C$ stand for number of input images in a batch, height, width, and number of color (or feature) channels.. In our code, we have localized these parameters as input variables to ease the transition from one to multiple color bands. Classifying lenses with multiple color bands will be easier and more accurate since foreground and background sources have a color contrast. However, we have chosen to use one color channel so our algorithm can be sensitive to geometry rather than color contrast in order to expand its applicability to a wider range of bands as well as eliminating its need for multi-band images when unavailable. As it can be seen in Figure 5, after normalizing and enhancing these single-channeled images, LensFlow applies an average-pooling of kernel $4\times 4$ and stride of $4$. This means that the image is divided into $4\times 4$ segments and the average of each segment will become a pixel of the down-sampled output image. These down-sampled images are then fed to the first convolutional layer which consists of sixteen $5\times 5$ filters. The hyperbolic tangent function is chosen as the activation function for the neurons in this layer as well as the upcoming convolutional layers. This layer outputs sixteen feature maps which should be interpreted as one image with 16 feature channels, i.e. a $50\times 50\times 16$ image using our notation above. These feature maps are down-sampled using a max-pooling layer of size $2\times 2$ and stride of $2$ resulting in a $25\times 25\times 16$ reduced feature map. This means the input of the max-pooling layer is divided into $2\times 2$ tiles whose maximum values are selected to form down-sampled feature maps. This $25\times 25\times 16$ feature map is inputted to the second convolutional layer with twenty-five $5\times 5$ filters where each filter will scan all channels by combining them linearly, hence outputting a $25\times 25\times 25$ feature map to be down-sampled to $13\times 13\times 25$. Another convolutional layer with thirty-six $4\times 4$ filters paired with a max-pooling layer is followed, reducing the data size to a $7\times 7\times 36$ feature map. This information is flattened into a 1-D array where it is fully connected to 128 sigmoid neurons. The output of these neurons is inputted to the classifying layer with 2 neurons with the softmax activation function:

[TABLE]

Here, each component of $\mathbf{Z}$ is the total weighted input plus the bias of each of the classifying neurons, $z_{c}=\sum_{j}{a_{j}^{L-1}w_{j\rightarrow c}^{L}}+b_{c}^{L}$ . $c$ specifies whether we are talking about the unlensed neuron or the lensed neuron. This function enhances the training process and it will assign a pseudo-probabilistic number to each class (since the probabilities across all classes add up to 1). The higher $\sigma_{c=\text{lensed}}$ for an image, the higher its probability to be a lens which can be used to rank the images based on their lensing probability.

3.4. Measuring Performance

To optimize our ConvNet, we have chosen a cross-entropy function as our loss function which we minimize using Adam’s Optimizer. During the training phase, 64 unlensed and 64 lensed images were placed in a batch of 128 images. This combination technique will prevent under- or over-representation of classes even if the training size for different classes contain a different number of examples. This process was iterated one thousand times. An example of accuracy as a function of the number of iterations has been graphed in Figure 6 where the accuracy plateaus at $91.5\%$. To obtain this accuracy, any image with a lensing probability lower/higher than $0.5$ was labeled lensed/unlensed. These predicted labels were compared against the true label of each image. Thus, accuracy is defined as the number of correct classifications over the number of images in the testing set.

[TABLE]

This result was achieved with a ”ReLu”333A rectified linear unit is a neuron with the following common activation function:$f(x)=\max(0,x)$. neurons in the convolutional layers. After switching to hyperbolic tangent, our results improved to $96\%$ meaning that for every 100 test images, 4 were misclassified. If the testing dataset contains simulated lenses, the reader should not translate these accuracies as a measure of physical lens classification performance, neither for our paper nor for others. Unlike physical lenses, simulated lenses may contain statistical artifacts and are limited in shape. Hence, the true performance of a classifier should be measured with its ability to locate physical lenses in a sufficiently large field.

To optimize computational time while addressing this issue, we have mixed 15 known lenses from the COSMOS field (Faure et al. 2011) with 3.5 thousand images from the same field and have assigned a lensing probability to each image that was obtained from the output of the lens class neuron. These images were ranked based on their lensing probabilities. The relative ranking for these lenses, i.e. their rank divided by the total number of scanned images, has been plotted in Figure 7. We have also tried eightfold scanning where the probabilities of eight transformations discussed above (Section 3.1) have been summed to improve the accuracy. This technique does not show an improved ranking of the lenses. The majority of the lenses fall under top $8\%$ where they can be further examined by eye. In the next section, we will discuss the images that contaminate the high-rank candidates. We will also discuss the remedies to further separate lenses from other sources.

4. Results & Discussion

A catalog of the strong gravitational lenses in the COSMOS field has been generated previously (Faure et al. 2008) by looking at early type bright galaxies in the redshift range of $0.2\leq z\leq 1.0$ and visually inspecting and cataloging sixty high and low-quality lens candidates. In contrast, we have examined all sources in HST/ACS i-band of the COSMOS field that are more extended than 200 pixels and are $1.5\sigma$ brighter than the background, i.e. 230,000 images. After scanning these images with LensFlow, we inspected the top outputs and selected 21 as good quality lens candidates which are presented in Figure 9. Their coordinates and other physical parameters are also listed in Table LABEL:tab:lens_coords. The lens candidates previously identified by Faure et al. (2008) are marked in Table 1. Among sixty lens candidates presented in Faure et al. (2008), 25 were more extended than our image cutout size of $200\times 200$ pixels which translates to larger than $3\arcsec\times 3\arcsec$ (e.g. COSMOS0208+1422, COSMOS0009+2455) while some others were considered as lens contaminates or low quality lenses rather than as good lens candidates during our visual selection process (e.g. COSMOS0055+3821). To resolve the former issue, a secondary scan can be performed by down-sampling the field to half its size to include larger lenses.

The main contaminants of our high lensing probability images fall into the main categories of, spiral galaxies with contamination arising from spiral arms and/or ring-like structures, tidally interacting galaxies and satellite/nearby galaxies in a lens-like configuration and image artifacts. Examples of these contaminants are shown in Figure 8. Given that these contaminants are similar to many of our training examples and that it is also time-consuming to separate them by human eye we created another training dataset which would eliminate such contaminants by only looking for perfect Einstein rings. After another scan of the field, we were able to identify the Einstein ring shown in the last two panels of Figure 9.

These indicate that perhaps, rather than increasing the number of layers in series, we must aggregate parallel ConvNets, each one trained to extract lenses from one class of contaminants only. Furthermore, reducing the down-sampling factor would, in turn, increase the resolution of the input images which might increase the gap between spirals and lenses at a greater computational cost. And perhaps, it would be more efficient to have separated training datasets for visually distinct lenses, as tried above, rather than having one large training dataset with a variety of lenses. We will study these algorithms in our future publications.

5. Summary

In this paper, we have emphasized the importance of gravitational lenses in the field of cosmology and have presented an introduction to neural networks including ConvNets. Furthermore, we have laid out the procedure for constructing simulated images for training and testing LensFlow. The importance and details of our normalization and enhancement methods have been discussed. Then, we have discussed the architecture of LensFlow, its accuracy on test data, and its performance on real data. The latter was done by scanning HST/ACS i-band images of the COSMOS field and listing the lens candidates that were visually separated from images with large assigned lens probability by LensFlow.

Acknowledgement

Financial support for this paper was provided by NSF grant AST-1213319 and GAANN P200A150121. Figure 1 and 5 were generated using www.draw.io. The backbone of our algorithm was initially inspired by Hvass Laboratories on Github (Pedersen 2016). Figure 2 was plotted using his code as well. We would like to thank Pierre Baldi for valuable feedback. We would also like to thank Noah Ghotbi and Aroosa Ansari for their assistance.

Appendix A Remarks on Present and Future of LensFlow Code

LensFlow has been written in Python 3.5.2 and to enhance the user interface, it was developed in a Jupyter Notebook environment which enables the user to easily modify the code, plot, and read the documentation alongside the code. LensFlow relies on Keras, a high-level neural networks API, written in Python and capable of running on top of TensorFlow(Chollet 2015). Our product consists of three main notebooks: Development notebook, Arc Maker notebook, LensFlow notebook.

The Development notebook contains the necessary code for file management, to automate source extraction, to merge raw images with simulated arcs, to create training and testing datasets, for automated handling of 81 COSMOS tiles, and includes other helpful functions. Currently, LensFlow is in its development stage and requires a list of $200\times 200$-pixel FITS cutouts. In future versions, we will eliminate this need and users could provide a large FITS tile with/without a catalog of identified sources which would require less/more time to scan the data. This can be done by convolving the ConvNet itself over a FITS tile and identify spots with high lensing probabilities. To make this fully functional, the training dataset must include off centered and cropped sources. However, for the purposes of training and optimizing, isolated cutouts are more useful since sources from different tiles can be mixed without facing memory shortage and since sources have to be accessed again and again as one would change the architecture of the ConvNet.

Arc Maker notebook automates the arc creation process and uses LensTool. LensFlow notebook contains the ConvNet and can be used to specify the architecture of the ConvNet, to train or to test the ConvNet, and to search through new data for lens candidates. This notebook also contains the functions that we have used to generate some of the plots in this paper as well as the feature to view multiple top lens candidates in one screen using DS9. This feature will help the user to quickly examine the output of the LensFlow. Due to the restrictions that TensorFlow and LensTool apply, LensFlow would only function in a Linux operating system.