Astrophysics, Computer Science

M. J. Smith

J. E. Geach

# Astronomia ex machina: a history, primer, and outlook on neural networks in
astronomy

Michael J. Smith and James E. Geach  Department of Physics, Astronomy &
Mathematics, School of Physics, Engineering and Computer Science, University
of Hertfordshire, Hatfield, AL10 9AB<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

In this review, we explore the historical development and future prospects of
artificial intelligence (AI) and deep learning in astronomy. We trace the
evolution of connectionism in astronomy through its three waves, from the
early use of multilayer perceptrons, to the rise of convolutional and
recurrent neural networks, and finally to the current era of unsupervised and
generative deep learning methods. With the exponential growth of astronomical
data, deep learning techniques offer an unprecedented opportunity to uncover
valuable insights and tackle previously intractable problems. As we enter the
anticipated fourth wave of astronomical connectionism, we argue for the
adoption of GPT-like foundation models fine-tuned for astronomical
applications. Such models could harness the wealth of high-quality, multimodal
astronomical data to serve state-of-the-art downstream tasks. To keep pace
with advancements driven by Big Tech, we propose a collaborative, open-source
approach within the astronomy community to develop and maintain these
foundation models, fostering a symbiotic relationship between AI and astronomy
that capitalizes on the unique strengths of both fields.

###### keywords:

neural networks, astrophysics, machine learning

##  1 Introduction

The concept of artificial intelligence (AI) can be traced back at least 350
years to Leibniz’s Dissertation on the Art of Combinations [(1)]. Inspired by
Descartes and Llull, Leibniz posited that, through the development of a
‘universal language,’ all ideas could be represented by the combination of a
small set of fundamental concepts, and that new concepts could be generated in
a logical fashion, potentially by some computing machine. Leibniz’s ambitious
vision (‘let us calculate’) has not yet been realised, but the quest to
emulate human reasoning, or at least to build a machine to mimic the
computational and data processing capabilities of the human brain, has
persisted to this day.

###### Contents

  1. 1 Introduction
  2. 2 A primer on artificial neurons
    1. 2.1 The perceptron
    2. 2.2 The multilayer perceptron
  3. 3 Astronomy’s first wave of connectionism
    1. 3.1 Classification problems
    2. 3.2 Regression problems
  4. 4 Contemporary supervised deep learning
    1. 4.1 Convolutional neural networks
    2. 4.2 Recurrent neural networks
    3. 4.3 Sidestepping the vanishing gradient problem
      1. 4.3.1 Non-saturating activation functions
      2. 4.3.2 Graphics processing unit acceleration
      3. 4.3.3 Gated recurrent neural networks and residual networks
    4. 4.4 Translation, attention, and transformers
  5. 5 Astronomy’s second wave of connectionism
    1. 5.1 Convolutional neural network applications
    2. 5.2 Recurrent neural network applications
    3. 5.3 Transformer applications
    4. 5.4 A problem with supervised learning
  6. 6 Deep generative modelling
    1. 6.1 (Variational) autoencoders
    2. 6.2 Generative adversarial networks
    3. 6.3 Score-based generative modelling and diffusion models
      1. 6.3.1 Forward process
      2. 6.3.2 Reverse process
      3. 6.3.3 Denoising diffusion implicit models
  7. 7 Representation learning
    1. 7.1 Contrastive learning
    2. 7.2 Learning representations via a surrogate task
  8. 8 Astronomy’s third wave of connectionism
    1. 8.1 Deep astronomical generative modelling
    2. 8.2 Self-supervised astronomical representation learning
  9. 9 Foundation models: a fourth astroconnectionist wave?
    1. 9.1 Foundation models
    2. 9.2 Scaling laws and data moats
    3. 9.3 The practical implications and uses of an astrofoundation model
      1. 9.3.1 Democratising foundation models
      2. 9.3.2 Possible astronomical use-cases
      3. 9.3.3 A new class of simulation
  10. 10 Connectionism’s caveats
    1. 10.1 Possible practical pitfalls
    2. 10.2 Connectionism’s carbon crisis
  11. 11 Final comments, or how we learnt to stop worrying and love astronomy’s Big Data Era

It might be fair to say that the roots of AI stretch even as far back as
Llull’s medieval philosophy that inspired Leibniz [(97), (143)]. However, if
we now consider AI to be a bona fide scientific discipline, then that
discipline clearly emerged in the post-war years of the twentieth century,
following Turing’s simple enquiry ‘can machines think?’ [(4)]. Somewhat
philosophical in nature, Turing’s 1950 question succinctly articulates the
ambition of AI, but from a nuts and bolts standpoint it took a further five
years from Turing’s query for what one might call the first AI program—the so-
called ‘Logic Theorist’—to be developed by Allen Newell, Cliff Shaw, and
Herbert Simon. Funded by the Research and Development (RAND) Corporation, the
Logic Theorist was designed, in part, to emulate the role of a human
mathematician, in that it could automate the proof of mathematical theorems.
This was a breakthrough in computer science and the Logic Theorist was
presented at the seminal Dartmouth Summer Research Project on Artificial
Intelligence (DSRPAI) conference in 1956, now regarded as the true birth of AI
as a field. Indeed, it was DSRPAI organiser John McCarthy who is credited with
coining the term ‘artificial intelligence’ [(78)].

$\displaystyle{1995}$$\displaystyle{2000}$$\displaystyle{2005}$$\displaystyle{2010}$$\displaystyle{2015}$$\displaystyle{2020}$Paper
submission
date$\displaystyle{0}$$\displaystyle{20}$$\displaystyle{40}$$\displaystyle{60}$$\displaystyle{80}$Submissions
per month Figure 1: Here we see the number of arXiv:astro-ph submissions per
month that have abstracts or titles containing one or more of the strings:
‘machine learning’, ‘ML’, ‘artificial intelligence’, ‘AI’, ‘deep learning’, or
‘neural network’. The raw data is in the public domain and is available at
https://www.kaggle.com/Cornell-University/arxiv.

Natural intrigue—and clearly a good deal of fear—of the idea of AI has
inspired popular culture no end, from Dick’s Do Androids Dream Of Electric
Sheep? to Crichton’s Westworld, Terminator’s ‘Skynet’ and beyond. Iain M.
Banks’s Galactic civilisation known as ‘The Culture’ imagines a society run by
powerful ‘Minds’ whose intelligence and wisdom far exceeds that of humans, and
where biological beings and machines of equivalent sentience generally co-
exist peacefully, cooperatively, and equitably. Science fiction
notwithstanding, if these dreams are even possible, we are still years away
from a machine that can genuinely think for itself [(111), (197)].
Nevertheless, the question of how one mathematically (and algorithmically)
models the workings and inter-relationships of biological neurons—neural
networks—and the subsequent exploration of how they can find utility as tools
in the data analyst’s workshop is really what is being referred to when most
people use the term ‘AI’ today111And the term is regularly misused, not only
erroneously, but often cynically.. While we must always be wary of hype and
buzzwordism, it is the application of neural networks—and the possibility of
tackling hitherto intractable problems—that offers genuine reason for
excitement across many disparate fields of enquiry, including astronomy.

Astronomers have made use of artificial neural networks (ANNs) for over three
decades. In 1994, Ofer Lahav, an early trailblazer, wryly identified the
‘neuro-skeptics’—those resistant to the use of such techniques in serious
astrophysics research—and argued that ANNs ‘should be viewed as a general
statistical framework, rather than as an estoteric approach’ [(40)].
Unfortunately, this skepticism has persisted. This is despite the recent
upsurge in the use of neural networks (and machine learning in general) in the
field, as illustrated in Fig. 1. This skepticism also stands contrary to
achievements within astronomy that would not be possible without the use of
ANNs, such as photometric redshift estimation [(64), (63), e.g.], astronomical
object identification and clustering at scale [(243), e.g.], and entirely
data-driven simulation [(275), (301), e.g.]. Most of the criticism of machine
learning techniques, and deep learning222Deep learning referring to the use of
a network constructed of many layers of artificial neurons. in particular, is
levelled at the perceived ‘black box’ nature of the methodology. In this
review we provide a primer on how deep neural networks are constructed, and
the mathematical rules governing their learning, which we hope will serve as a
useful resource for neuro-skeptics. Nevertheless, we must recognise that a
unified theoretical picture of how deep neural networks work does not yet
exist. This remains a point of debate even within the deep learning community.
For example, Yann LeCun responding to Ali Rahimi’s ‘Test of Time’ award talk
at the 31st Conference on Neural Information Processing Systems (NIPS)
remarked:

> _Ali gave an entertaining and well-delivered talk. But I fundamentally
> disagree with the message. The main message was, in essence, that the
> current practice in machine learning is akin to ‘alchemy’ (his word). It’s
> insulting, yes. But never mind that: It’s wrong! Ali complained about the
> lack of (theoretical) understanding of many methods that are currently used
> in ML, particularly in deep learning … Sticking to a set of methods just
> because you can do theory about it, while ignoring a set of methods that
> empirically work better just because you don’t (yet) understand them
> theoretically is akin to looking for your lost car keys under the street
> light knowing you lost them someplace else. Yes, we need better
> understanding of our methods. But the correct attitude is to attempt to fix
> the situation, not to insult a whole community for not having succeeded in
> fixing it yet. This is like criticizing James Watt for not being Carnot or
> Helmholtz._ [(161)]

Philosophical concerns aside, LeCun’s fundamental point is that deep learning
‘works’ and therefore we should use it, even if we do not fully understand it.
If one were being uncharitable, we could make similar arguments about the
$\Lambda$CDM paradigm.

It is clear that in every field that deep learning has infiltrated we have
seen a reduction in the use of specialist knowledge, to be replaced with
knowledge automatically derived from data. We have already seen this process
play out in many ‘applied deep learning’ fields such as computer Go [(151)],
protein folding [(246)], natural language processing [(211)], and computer
vision [(216)]. We argue that astronomy’s data abundance corrals it onto a
path no different to that trodden by other applied deep learning fields. This
abundance is not a passing phase; the total astronomical data volume is
already large and will increase exponentially in the coming years. We
illustrate this in Fig. 2, where we present a selection of astronomical
surveys and their estimated data volume output over their lifetimes [(139)].
And this is not even considering data associated with ever larger and more
detailed numerical simulations [(180), (234), (271), e.g.]. The current scale
of the data volume already poses an issue for astronomy as many classical
methods rely on human supervision and specialist expertise, and the increasing
data volume will make exploring and exploiting these surveys through
traditional human supervised and semi-supervised means an intractable problem.
Of serious concern is the possibility that we will miss—or substantially
delay—interesting and important discoveries simply due to our inability to
accurately and consistently interrogate astronomical data at scale. Deep
learning has shown great promise in automating information extraction in
various data intensive fields, and so is ideally poised as a solution to the
challenge of processing ultra-large scale astronomical data. But we do not
need to stop there. This review’s outlook ventures a step further, and argues
that astronomy’s wealth of data should be considered a unique opportunity, and
not merely an albatross.

Since astronomical connectionism’s333 Since its inception, AI research can be
broadly categorised into two schools: ‘symbolic’ and ‘connectionist’.
Symbolists see the mind as a collection of fully-formed representations, and
attempt to mimic human reasoning through a logical rule-based processing of
these symbols. This approach contrasts with connectionist (or neural network-
based) AI, which takes a bottom-up approach and simulates cognition by
mimicking the way neurons in the human brain work.  humble beginnings in the
late 1980s, there have been numerous excellent reviews on the application of
artificial neural networks to astronomy [(32), (89), (318), e.g.]. We take an
alternative approach to previous literature reviews and survey the field
holistically, in an attempt to paint astronomical connectionism’s ‘Big
Picture’ with broad strokes. While we cannot possibly include all works within
astronomical connectionism444 We refer the reader to Fig. 1! , we hope that
this review serves as a historical background on astronomy’s ‘three waves’ of
increasingly automated connectionism, as well as presenting a general primer
on neural networks that may assist those seeking to explore this fascinating
topic for the first time.

$\displaystyle{1995}$$\displaystyle{2000}$$\displaystyle{2005}$$\displaystyle{2010}$$\displaystyle{2015}$$\displaystyle{2020}$$\displaystyle{2025}$$\displaystyle{2030}$Year
seeing first
light$\displaystyle{10^{2}}$$\displaystyle{10^{4}}$$\displaystyle{10^{6}}$$\displaystyle{10^{8}}$Data
volume over lifetime
[TB]DPOSS2MASSGBTGALEXSDSSSkymapperPanSTARRS*LSST*SKA**expected$\displaystyle
10\cdot 2^{\,(y-1997)/1.36}$ Figure 2: The data volume output of a selection
of astronomical surveys over their lifetimes. We can see the astronomical
survey data volume doubles every 16 months. Data is taken from [(139)].

In §2 and §3 we explore initial work on multilayer perceptrons within
astronomy, where models required manually selected emergent properties as
input. In §4 and §5 we explore the second wave, which coincided with the
dissemination of convolutional neural networks and recurrent neural
networks—models where the multilayer perceptron’s manually selected inputs are
replaced with raw data ingestion. In the third wave that is happening now we
are seeing the removal of human supervision altogether with deep learning
methods inferring labels and knowledge directly from the data, and we explore
this wave in §6–§8. Finally, in §9, we look to the future and predict that we
will soon enter a fourth wave of astronomical connectionism. We argue that if
astronomy follows the pattern of other applied deep learning fields we will
see the removal of expertly crafted deep learning models, to be replaced with
fine-tuned versions of an all-encompassing ‘foundation’ model. As part of this
fourth wave we argue for a symbiosis between astronomy and connectionism, a
symbiosis predicated on astronomy’s relative data wealth and deep learning’s
insatiable data appetite. Many ultra-large datasets in machine learning are
proprietary or of poor quality, and so there is an opportunity for astronomers
as a community to develop and provide a high quality multimodal public
dataset. In turn, this dataset could be used to train an astronomical
foundation model to serve state-of-the-art downstream tasks. Due to foundation
models’ hunger for data and compute, a single astronomical research group
could not bring about such a model alone. Therefore, we conclude that
astronomy as a discipline has slim chance of keeping up with a research pace
set by the Big Tech goliaths—that is, unless we follow the examples of
EleutherAI and HuggingFace and pool our resources in a grassroots open source
fashion.

Before moving on, we must first admit to our readers that we have not been
entirely honest with them. The abstract of this review has not been written by
us. It was generated by prompting OpenAI’s generative pretrained transformer 4
(‘GPT-4’) neural network-based foundation model with this paper’s introduction
[(322), (313)]. To be precise, we prompted the GPT-4 engine provided by
‘ChatGPT Plus’ with all the text in §1 up until this paragraph in raw LaTeX
format. We then appended the following prompt to the introduction text:

> _Write an abstract for the above text that will catch the reader’s eye, and
> make them interested in the paper. Make the abstract 160 words or less, and
> touch on the value of GPT-like models in astronomy._

We did not alter the GPT generated output whatsoever. We explore these
foundation models and their possible astronomical uses in more detail in §9.

##  2 A primer on artificial neurons

In 329 ref_mcculloch1943 proposed the first computational model of a
biological neuron (MP neuron; ref_mcculloch1943, 329). Their model consisted
of a set of binary inputs $x_{i}\in\\{0,1\\}$ and a single binary output
$y\in\\{0,1\\}$. Their model also defines a single ‘inhibitory’ input
$\mathcal{I}\in\\{0,1\\}$ that blocks output if $\mathcal{I}=1$. If the sum of
the inputs exceeds a threshold value $\Theta$, the MP neuron ‘fires’ and
outputs $y=1$. Mathematically we can write the MP neuron function as

$\text{MP}(\mathbf{x})=\begin{cases}\mbox{1}&\mbox{if
}\sum^{n}_{i=1}x_{i}>\Theta\mbox{ and }\mathcal{I}=0,\\\
\mbox{0}&\mbox{otherwise.}\end{cases}$

The MP neuron is quite a powerful abstraction. Single MP neurons can calculate
simple boolean functions, and more complicated functions can be calculated
when many MP neurons are chained together. However, there is one show-stopping
issue: the MP neuron is missing the capacity to learn. ref_perceptron (332)
addressed this by combining the MP neuron with Hebb’s neuronal wiring
theory555Also known by the mantra ‘cells that fire together wire together’.
(ref_hebb1949, 330), and we will explore a related training formulation in the
next subsection.

###  2.1 The perceptron

This subsection aims to provide the reader a foundation and intuition for the
gradient-based learning that dominates contemporary neural network
architectures. Therefore, we diverge from ref_perceptron’s original learning
algorithm and instead describe a gradient-based training algorithm. The
interested reader will find an analysis of ref_perceptron’s original learning
algorithms in the ‘Mathematical analysis of learning in the perceptron’
section of ref_perceptron (332).

Like the MP neuron, the perceptron takes a number of numeric inputs ($x_{i}$).
However, unlike the MP neuron each one of these inputs is multiplied by a
corresponding weight ($w_{i}$) signifying the importance the perceptron
assigns to a given input. As shown in Fig. 3, we can then sum this list of
products and pass it into an ‘activation function’. Let us use the Heaviside
step function as our activation function:

$-5$$0$$5$$0.5$$1$$\mathbf{w}\cdot\mathbf{x}$

$\displaystyle\text{prediction}=H(\mathbf{w}\cdot\mathbf{x})=\begin{cases}\mbox{0}&\mbox{if
}\mathbf{w}\cdot\mathbf{x}<0,\\\ \mbox{1}&\mbox{if
}\mathbf{w}\cdot\mathbf{x}\geq 1,\end{cases}$ (1)

where $\mathbf{x}$ is a set of inputs, and $\mathbf{w}$ is a set of ‘weights’
that represent the importance of each input.

To concretise how we could train our perceptron we will use an example. Let us
say that we want to automatically label a set of galaxy images as either
‘spiral’ or ‘elliptical’. To do this we first need to compile a training
dataset of galaxy images. This training set would consist of spiral and
elliptical galaxies, and each image would have a ground truth label $y$—say
‘0’ for a spiral galaxy and ‘1’ for an elliptical. To train our perceptron we
randomly choose one image from the training set, and feed it to the
perceptron, with the numerical value of each pixel corresponding to an input
$\\{x_{1},\ldots,x_{N}\\}$. These inputs are multiplied by their corresponding
weight $\\{w_{1},\ldots,w_{N}\\}$. A bias term
$(b=w_{0}\,x_{0},\text{where}\,x_{0}=1)$ is also added to the inputs, which
allows the neuron to shift its activation function linearly. Since we do not
want our perceptron to have any prior knowledge of the task, we initialise the
weights at random. The resulting products are then summed. Finally, our
activation function $H$ transforms $\mathbf{w}\cdot\mathbf{x}$ and produces a
prediction $p$. We then compare $p$ to $y$ via a ‘loss function,’ which is a
function that measures the difference between $p$ and $y$. The loss can be any
differentiable function, so for illustration purposes we will define it here
as the L1 loss: $\mathcal{L}(y,p)=\lvert y-p\rvert$. Now that we can compare
to the ground truth, we need to work out how a change in one of our weights
affects the loss (that is, we want to find
$\partial\mathcal{L}/\partial\mathbf{w}$). We can calculate this change with
the chain rule

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}}=\frac{\partial\mathcal{L}}{\partial
p}\frac{\partial p}{\partial\mathbf{w}},$ (2)

and since $p=H(\mathbf{w}\cdot\mathbf{x})$ and $\partial
p/\partial\mathbf{w}=H^{\prime}\mathbf{x}^{T}$ we get

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}}=\frac{\partial\mathcal{L}}{\partial
p}\,\odot(H^{\prime}\mathbf{x}^{T})$

where $\odot$ is the distributive Hadamard product. Thus we can update the
weights to decrease the loss function:

$\begin{split}\textbf{w}_{\text{next}}&=\mathbf{w}-\eta\frac{\partial\mathcal{L}}{\partial\mathbf{w}}\\\
&=\mathbf{w}-\eta\frac{\partial\mathcal{L}}{\partial
p}\odot(H^{\prime}\mathbf{x}^{T}),\end{split}$

where $\eta$ is the learning rate666 The eagle-eyed reader may have noticed
that since the derivative of the Heaviside step function is the Dirac delta
function, we will only update the perceptron’s weights on an incorrect
prediction. If we want to also learn from positive examples, we need to use a
smoothly differentiable activation function. This is explored in the next
subsection. . If we repeat this process our perceptron will get better and
better at classifying our galaxies!

Figure 3: A single neuron (or perceptron) with a bias $w_{0}$, inputs
$x_{1},x_{2},\ldots,x_{N}$, and weights $w_{1},w_{2},\ldots,w_{N}$.

While we provide the above example for illustrative purposes, we will need a
more powerful algorithm to produce a useful classifier of galaxy morphology.
This need is perhaps most famously discussed in _Perceptrons: An Introduction
to Computational Geometry_ (e.g. §13.0; ref_minsky1969, 335). ref_minsky1969
show that the single layer perceptron is only able to calculate linearly
separable functions, among other limitations. Their book (alongside a
consensus that AI had failed to deliver on its early grandiose promises)
delivered a big blow to the connectionist school of artificial intelligence777
See ref_olazaran1996 (374) and ref_metz2021 (578) for a closer look at the
conflicts and personalities that shaped AI. . In the years following
ref_minsky1969 (335) governmental and industry funding was pulled from
connectionist research laboratories, ushering in the first ‘AI winter’888 At
least, in the Western world. Connectionism continued in earnest in the Soviet
Union (ref_ivakhnenko1965, 334, 337). .

Yet, as exemplified in ref_rosenblatt1962 (§5.2, theorem 1; 333) it was known
at the time that multilayer perceptrons could calculate non-linearly separable
functions (such as the ‘exclusive or’). We can prove intuitively that a set of
neurons can calculate _any_ function: a perceptron can perfectly emulate a
NAND gate (Fig. 4), and the singleton set $\\{\texttt{NAND}\\}$ is
functionally complete. Since we can combine a set of NAND gates to calculate
any function, _we must also be able to combine a set of neurons to calculate
any function_. This result is also explored in a more formal proof by both
ref_cybenko1989 (345) and ref_kurt1991 (352). They show that an infinitely
wide neural network can calculate any function. Similarly, ref_lu2017width
(489) show that an infinitely deep neural network is a universal approximator.
Such a group of neurons is known as the multilayer perceptron (MLP).
Unfortunately, we cannot simply stack perceptrons together as we are missing
one vital ingredient: a way to train the network! At the time of
ref_minsky1969’s treatise on perceptrons there was no widely known algorithm
(in the West; see ref_ivakhnenko1965, 334) that could train such a multilayer
network. In ref_minsky1969’s own words:

> _Nevertheless, we consider it to be an important research problem to
> elucidate (or reject) our intuitive judgment that the extension [from one
> layer to many] is sterile. Perhaps some powerful convergence theorem will be
> discovered, or some profound reason for the failure to produce an
> interesting ‘learning theorem’ for the multilayered machine will be found._
> (§13.2; ref_minsky1969 (335), on MLPs)

The field had to wait almost two decades for such an algorithm to become
widespread. In the next subsection we will explore backpropagation, the
algorithm that ultimately proved ref_minsky1969’s intuition wrong.

$x_{1}$ | $x_{2}$ | $\neg(x_{1}\land x_{2})$ |  $p=H(\mathbf{w}\cdot\mathbf{x}$)
---|---|---|---
0 | 0 | 1 | $H(1.5+(-1)\cdot 0+(-1)\cdot 0)=1$
0 | 1 | 1 | $H(1.5+(-1)\cdot 0+(-1)\cdot 1)=1$
1 | 0 | 1 | $H(1.5+(-1)\cdot 1+(-1)\cdot 0)=1$
1 | 1 | 0 | $H(1.5+(-1)\cdot 1+(-1)\cdot 1)=0$

Figure 4: If we define $H(\mathbf{w}\cdot\mathbf{x})$ as in Eq. 1 we can set a
perceptron’s weights so that it is equivalent to the NAND gate.

###  2.2 The multilayer perceptron

Grouping many artificial neurons together may result in something resembling
Fig. 5. This network consists of an input layer, two intermediate ‘hidden’
layers, and an output layer. As in the previous section, let us say that we
want a classifier that can classify a set of galaxy images into elliptical and
spiral types. In an MLP similar to Fig. 5 a neuron would be assigned to each
pixel in a galaxy image. Each neuron would take the numeric value of that
pixel, and propagate that signal forward into the network. The next layer of
neurons does the same, with the input being the previous layer’s output. This
process continues until we reach the output layer. In a binary classification
task like our galaxy classifier this layer outputs a value between zero and
one. Thus, if we define a spiral galaxy as zero, and an elliptical galaxy as
one, we would want the network output to be near zero for a spiral galaxy
input (and vice versa).

Figure 5: The multilayer perceptron, or artificial neural network. The
depicted network has two hidden layers. It takes $N$ inputs
$x_{1},x_{2},\ldots,x_{N}$, and outputs a prediction $p_{L}$. Note that here
we omit the explicit bias terms (i.e. $w_{0}$).
$-5$$5$$-1$$1$Sigmoid$-5$$5$$-1$$1$Tanh$-5$$5$$-1$$1$ReLU$-5$$5$$-1$$1$ELU$-5$$5$$-1$$1$Swish$-5$$5$$-1$$1$Mish
Figure 6: A curated selection of activation functions. In all plots, the x
axis is the input, and the y axis is the output. The rectified linear unit
(ReLU) activation function was first introduced in the context of neural
networks in ref_neocognitron (339) and later rediscovered, named, and
popularised in ref_nair2010relu (420). The exponential linear unit (ELU),
Swish and Mish activations were respectively introduced in ref_clevert2015
(468), ref_ramachandran2017swish (492), and ref_mish (523).

In §2.1 we found the change we needed to apply to a single neuron’s weights to
make it learn from a training example. We can train an MLP in a similar way by
employing the reverse mode of automatic differentiation (or backpropagation)
to learn from our galaxy training data set (ref_backpropenglish, 338, 340,
342)999 Some controversy surrounds backpropagation’s discovery. The Finnish
computer scientist ref_backpropenglish proposed the reverse mode of automatic
differentiation and adapted the algorithm to run on computers in their 1970
(Finnish language) thesis (ref_backpropfinnish, 336). They first published
their findings in English in 338. ref_werbos1981 then proposed applying an
adaptation of ref_backpropenglish’s method to artificial neural networks.
ref_rumelhart1986 (342) showed experimentally that backpropagation can
generate meaningful internal representations within a neural network, and
popularised the method. Here we will err on the side of caution and cite all
three manuscripts. For further reading we recommend ref_schmidhuber2014 (445)
and ref_baydin2018backprop (496). . We want our network to learn when it makes
both a correct and incorrect prediction, so we define our activation function
as a smoothed version of ref_perceptron’s perceptron activation. This ensures
that a signal is present in the derivative no matter which values are input.
This activation function is known as the ‘sigmoid’ function, and is shown in
Fig. 6. As in §2.1 we define a loss function $\mathcal{L}(y,p)$ that describes
the similarity between a ground truth ($y$) and a prediction ($p$). We also
define a neuron’s activation function as $\varphi(\mathbf{w}\cdot\mathbf{x})$
where $\mathbf{w}\cdot\mathbf{x}$ is the weighted sum of a neuron’s inputs.
Following from Eq. 2:

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}_{l}}=\frac{\partial\mathcal{L}}{\partial\mathbf{p}_{l}}\frac{\partial\mathbf{p}_{l}}{\partial\mathbf{w}_{l}}$

where $l$ is a layer in the MLP. In the same way as in §2.1 we can calculate
an MLP’s final layer’s ($l=L$) weight updates in terms of known values:

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}_{L}}=\frac{\partial\mathcal{L}}{\partial
p_{L}}\odot\left(\varphi_{L}^{\prime}\mathbf{p}^{T}_{L-1}\right),$ (3)

where $\mathbf{p}_{L-1}$ are the outputs from the previous layer. To calculate
the $(L-1)$th layer’s weight updates we use the chain rule:

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}_{L-1}}=\frac{\partial\mathcal{L}}{\partial
p_{L}}\frac{\partial
p_{L}}{\partial\mathbf{p}_{L-1}}\frac{\partial\mathbf{p}_{L-1}}{\partial\mathbf{w}_{L-1}}.$

Likewise for the $(L-n)$th layer:

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}_{L-n}}=\frac{\partial\mathcal{L}}{\partial
p_{L}}\left(\prod^{n}_{i=1}\frac{\partial\mathbf{p}_{L+1-i}}{\partial\mathbf{p}_{L-i}}\right)\frac{\partial\mathbf{p}_{L-n}}{\partial\mathbf{w}_{L-n}}.$

Now we can start plugging in some known values. Since
$\mathbf{p}_{l}=\varphi_{l}(\mathbf{w}_{l}\cdot\mathbf{p}_{l-1})$, it follows
that
$\partial\mathbf{p}_{l}/\partial\mathbf{p}_{l-1}=\varphi^{\prime}_{l}\mathbf{w}^{T}_{l}$,
and
$\partial\mathbf{p}_{l}/\partial\mathbf{w}_{l}=\varphi^{\prime}_{l}\mathbf{p}^{T}_{l-1}$.
So:

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}_{L-n}}=\frac{\partial\mathcal{L}}{\partial
p_{L}}\odot\left(\prod^{n}_{i=1}\varphi^{\prime}_{L-i}\mathbf{w}^{T}_{L-i}\right)\,\left(\varphi^{\prime}_{L-n}\mathbf{p}^{T}_{L-n-1}\right).$
(4)

Combining Eq. 3 with Eq. 4 we get the weight update algorithm for the
$(L-n)$th layer of the MLP:

$\mathbf{w}_{\text{next}}=\mathbf{w}-\eta\begin{cases}\frac{\partial\mathcal{L}}{\partial
p_{L}}\odot\left(\varphi_{L}^{\prime}\mathbf{p}^{T}_{L-1}\right),&\text{for
}n=0,\\\ \frac{\partial\mathcal{L}}{\partial
p_{L}}\odot\left(\prod^{n}_{i=1}\varphi^{\prime}_{L-i}\mathbf{w}^{T}_{L-i}\right)\,\left(\varphi^{\prime}_{L-n}\mathbf{p}^{T}_{L-n-1}\right),&\text{for
}n>0.\end{cases}$ (5)

With this equation101010 If we examine Eq. 5 carefully, we can see why we add
nonlinearities between the MLP layers; without activation functions Eq. 5
collapses to the equivilent of a single layer MLP!  in hand we can use the
same technique described earlier in this section and in §2.1 to update the
network’s weights with each galaxy image to decrease the loss function
$\mathcal{L}$. Again, as $\mathcal{L}$ is minimised, our MLP will classify our
elliptical and spiral galaxy images with increasing accuracy.

##  3 Astronomy’s first wave of connectionism

Connectionism was first discussed within astronomy in the late 1980s, after
the popularisation of backpropagation (see footnote 9) and the consequent
passing of the first ‘AI winter’. Two radical studies emerged in 1988 that
recognised areas where astronomy could benefit from the use of ANNs
(ref_rappaport1988, 344, 343). Together, they identified that astronomical
object classification111111 Specifically, galaxies were discussed in
ref_rappaport1988 (344) and point sources observed with the Infra-Red
Astronomical Satellite (IRAS) were discussed in ref_adorf1988 (343).  and
telescope scheduling could be solved through the use of an ANN. These studies
were followed by a rapid broadening of the field, and the application of
connectionism to many disparate astronomical use cases (ref_miller1993, 359,
and references therein). In this section, we will outline areas where MLPs
found an early use in astronomy.

###  3.1 Classification problems

ref_odewahn1992 (355) classified astronomical objects into star and galaxy
types. These were taken from the Palomar Sky Survey Automated Plate Scanner
catalogue (ref_pennington1993, 361). To compile their dataset, they first
extracted a set of emergent image parameters from the scanned observations.
These parameters included the diameter, ellipticity, area, and plate
transmission. The parameters were then used to train both a linear perceptron
and a feedforward MLP to classify the objects into stars or galaxies.
ref_odewahn1992 (355) found that their best performing model could classify
galaxies with a completeness of $95\%$ for objects down to a magnitude
$<19.5$. This work was followed by many more studies on the star/galaxy
classification problem (e.g. ref_odewahn1993, 360, 379, 371, 382). Galaxy
morphological type classification was explored in the early 1990s.
ref_storrie1992 (356) describe an MLP that takes an input a selected set of
thirteen galaxy summary statistics, and uses this information to classify a
galaxy into one of five morphological types. ref_storrie1992 (356) report a
top one accuracy of 64%, and a top two accuracy of 90%. This pilot study was
followed by several studies from the same group that confirmed that MLPs are
effective automatic galaxy morphological classifiers (ref_lahav1995, 368, 372,
370, 369, 373, 392, see §5 for a continuation of this line of research).

MLPs were also used in other classification tasks; here we highlight a few
further areas where MLPs were applied. ref_von1994 (365) classified stellar
spectra into temperature types, and ref_klusch1993 (358) did the same for
Morgan-Keenan System types. ref_chon1998 (380) described the use of an MLP to
search for and classify muon events (and therefore neutrino observations) in
the Sudbury Neutrino Observatory. Quasar classification has been explored in
several studies (ref_carballo2004, 394, 403, 409). Seminally, ref_carballo2004
(394) used an MLP to select quasar candidates given their radio flux,
integrated-to-peak flux ratio, photometry and point spread function in the red
and blue bands, and their radio-optical position separation. They found good
agreement between their model and that of the decision tree described in
ref_white2000 (385), confirming MLPs as a competitive alternative to more
traditional machine learning. As part of the Supernova photometric
Classification Challenge (SPCC; ref_kessler2010, 419), ref_karpenka2013 (432)
proposed the use of a neural network to classify supernovae into Type-1a/non-
Type-1a classes. To classify their light curves, they first used a hand-
crafted fitting function, and then trained their MLP on the fitted
coefficients. They found that their model was competitive with other, more
complex models trained on the SPCC dataset. From the studies discussed in this
section we can safely conclude that MLPs are effective classifiers of
astronomical data, when given important parameters extracted by an expert
guide.

###  3.2 Regression problems

MLPs have also been used in regression problems. ref_angel1990 (347) applied
them first to adaptive telescope optics. They trained their MLP on $250\,000$
simulated in focus and out of focus observations of stars as seen by the
Multiple Mirror Telescope (MMT). From the flattened $13\times 13$ pixel
observations, their network predicted the piston position and tilt required
for each of the MMT’s mirrors to bring the stars into focus. After the
application of these corrections, the authors were able to recover the
original profile. In follow up studies, ref_sandler1991 (353) and
ref_lloydhart1992 (354) proved that ref_angel1990’s MLP worked on the real
MMT.

Photometric redshift estimation was explored in many concurrent studies (e.g.
ref_firth2003, 390, 391, 398, 395, 392). ref_firth2003 (390) train a neural
network to predict the redshift of galaxies contained in the Sloan Digital Sky
Survey (SDSS) early data release (ref_stoughton2002, 389). The galaxies were
input to the neural network as a set of summary parameters, and the output was
a single float representing the galaxy redshift. They found their network
attained a performance comparable to classical techniques. Extending and
confirming the work by ref_firth2003 (390), ref_ball2004 (392) used an MLP to
predict the redshift of galaxies contained in the SDSS’s first data release
(ref_sdss, 386). They also showed that MLPs were capable of predicting the
galaxies’ spectral types and morphological classifications.

Of course, MLPs have been used more widely in astronomical regression tasks.
Here we will cherry pick a few studies to show the MLP’s early breadth of use.
Sunspot maxima prediction was carried out by ref_koons1990 (349). They found
their MLP based method was capable of predicting the number of sunspots when
trained on previous cycles. ref_bailer1997 (376) predicted the effective
temperature of a star from its spectrum. ref_auld2007 (406) and ref_auld2008
(408) applied MLPs to cosmology, demonstrating that MLPs are capable of
predicting the cosmic microwave background power spectra and matter power
spectra when given a set of cosmological parameters. ref_norgaardnielsen2008
(412) used an MLP to remove the foreground from microwave temperature maps.
From the studies discussed in this section we can see that MLPs are effective
regressors of astronomical data, when given significant parameters extracted
by an expert guide.

##  4 Contemporary supervised deep learning

There are some issues with MLPs. Primarily they do not scale well to high
dimensional datasets. For example, if our dataset consists of images with a
$128\times 128$ pixels, we will need $16\,384$ neurons in the MLP’s input
layer alone! As we move into the hidden layers, this scaling issue only gets
worse. Also, since MLPs must take an unrolled image as an input, they
disregard any spatial properties of their training images, and so either need
a substantial amount of training data to classify or generate large
images121212 At the height of the convolutional neural network architecture’s
popularity in the mid 2010s these were real problems. However, with the growth
of computing power and data in recent years we are seeing a resurgence of the
more general MLP model (e.g. ref_tolstikhin2021, 590, 591, 575, 577). This
follows the prevailing trend in AI where the removal of human-crafted features
and biases ultimately results in more expressive models that learn such
features and biases dieectly from data (ref_sutton2019, 534, 601). , or an
expert to extract descriptive features from the data in a preprocessing step.
We can see this issue writ large in the previous section—most of the MLP
applications described in §3 require an expert to extract features from the
data for the network to then train on! This drawback is not ideal; what if
there are features within the raw data that are not present in these cherry
picked statistics? In that case, it would be preferable to let the neural
network take in the raw data as input, and then learn which features are the
most descriptive. We will discuss neural network architectures that solve both
the MLP scaling problem and the expert reliance problem in this section. After
we have explored these architectures in general, we will discuss their
application to astronomical problems in §5.

###  4.1 Convolutional neural networks

Unlike the MLP described in the previous section, convolutional neural
networks (CNNs; introduced in ref_neocognitron (339) and first combined with
backpropagation in ref_lecun1989 (346)) do not entirely consist of fully
connected layers, where each neuron is connected to every neuron in the
previous and subsequent layers. Instead, the CNN (such as the one depicted in
Fig. 13) uses convolutional layers in place of the majority (or all) of the
dense layers.

Figure 7: A convolutional neural network classifying a spiral galaxy
image131313All astronomical objects shown in the neural network diagrams
within this manuscript are generated via text prompts fed into a latent
diffusion neural network model (ref_rombach2021, 583)..

We can think of a convolutional layer as a set of learnt ‘feature filters’.
These feature filters perform a local transform on input imagery. In classical
computer vision, these filters are hand crafted, and perform a predetermined
function, such as edge detection or blurring. In contrast, a CNN learns the
optimal set of filters for its task (say, galaxy classification). Eq. 6 shows
two different convolution141414 We must note that in Eq. 6 we follow most deep
learning libraries and perform a cross-correlation and _not_ a convolution.
However, since the weights are learnt, this does not matter; the neural
network will simply learn a flipped representation of the cross-correlation.
operators being performed on an array.

$\begin{split}\begin{bmatrix}39&57&86&9&26\\\ 90&74&63&87&98\\\
79&34&26&16&46\\\ 67&61&96&1&79\\\
33&47&15&49&29\end{bmatrix}\star\begin{bmatrix}0&0&0\\\ 0&0&0\\\ 0&0&1\\\
\end{bmatrix}&=\begin{bmatrix}26&16&46\\\ 96&1&79\\\ 15&49&29\end{bmatrix}\\\
\\\ \begin{bmatrix}39&57&86&9&26\\\ 90&74&63&87&98\\\ 79&34&26&16&46\\\
67&61&96&1&79\\\ 33&47&15&49&29\end{bmatrix}\star\begin{bmatrix}1&0&0\\\
0&1&0\\\ 0&0&1\\\ \end{bmatrix}&=\begin{bmatrix}139&136&219\\\ 220&101&158\\\
155&179&56\end{bmatrix}\end{split}$ (6)

In the above equation the operation is represented as a matrix. In a CNN the
matrix is a set of neuronal weights. As seen in Fig. 13 there are multiple
feature maps in a convolutional layer, each containing a set of weights
independent to the other feature maps, and learning to extract a different
feature. Due to the convolution operator’s inbuilt translational equivarience,
these features can be detected by the convolutional layer no matter where they
are in the image. As in the MLP described in the previous section, the weights
are updated using backpropagation to minimise a loss function. We will discuss
astronomical applications of CNNs in §5, after we introduce modern CNN
architectures.

###  4.2 Recurrent neural networks

Standard feedforward neural networks like the MLP (§2.2) and CNN (§4.1)
generate a fixed size vector given a fixed size input151515 As with any rule
there are exceptions, such as CNNs containing a global average pooling layer
(ref_lin2013, 434). . But, what if we want to classify or generate a variably
sized vector? For example, we might want to classify a galaxy’s morphology
given its rotation curve. A rotation curve describes the velocity of a
galaxy’s visible stars versus their distance from the galaxy’s centre. Fig. 8
shows a possible rotation curve for Messier 81. A rotation curve’s length
depends on the size of its galaxy, and due to this variable length, and the
fact that MLPs take a fixed size input, we cannot easily use an MLP for
classification. Recurrent neural networks (RNNs), however, can take a variable
length input and produce a variable length output.

Distance from centreRotational velocity Figure 8: An example of a galaxy
rotation curve, plotted over an image of Messier 81 (ref_m81, 451).

An RNN differs from a feed forward MLP by having a hidden state that acts as a
‘memory’ store of previously seen information. As the RNN encounters new data,
its weights are altered through the backpropagation through time algorithm
(BPTT; ref_werbos1990, 350, and references therein. Also see footnote 9).

We can use an RNN similar to Fig. 9 to classify our rotation curves. We
express the rotation curve as a list $\\{x_{1},x_{2},\ldots,x_{N}\\}$, with
each $x$ being a measurement of the rotational velocity at a certain radius.
Then we feed this list into the RNN sequentially in the same way as shown in
Fig. 9. The RNN will produce an output for each $x$ fed to it, but we ignore
those until we feed in $x_{N}$, the rotational velocity furthest from the
galaxy’s centre. When we feed in $x_{N}$, the RNN produces a prediction
$p_{N}$, which we can then compare to a ground truth $y_{N}$ via a loss
function $\mathcal{L}_{N}$. In our case, $y$ is an integer label representing
the galaxy’s morphological class. The comparison
$\mathcal{L}_{N}(y_{N},p_{N})$ is a function that represents the distance
between the RNN prediction and the ground truth. We can then reduce
$\mathcal{L}_{N}(y_{N},p_{N})$ by updating the RNN’s weights through BPTT so
that the weights $\\{\mathbf{w}_{x},\mathbf{w}_{p},\mathbf{w}_{h}\\}$ follow
$\nabla\mathcal{L}_{N}$ downwards. As we do this, our RNN will improve its
galaxy classifications.

Figure 9: A recurrent neural network with weights
$\\{\mathbf{w}_{x},\mathbf{w}_{p},\mathbf{w}_{h}\\}$, a hidden state
$\mathbf{h}_{n}$, inputs $\mathbf{x}$, and a prediction $p_{n=N}$ is unrolled
into its constituent processes.

BPTT’s mathematical derivation is akin to the one we explored in §2.2, and we
will quickly derive it here for posterity. Let us first look at the forward
propagation equations:

$\displaystyle\mathcal{L}_{n}$ $\displaystyle=\lvert y_{n}-p_{n}\rvert,$
$\displaystyle p_{n}$
$\displaystyle=\varphi(\mathbf{w}_{p}\cdot\mathbf{h}_{n}),\,\text{and}$
$\displaystyle\mathbf{h}_{n}$
$\displaystyle=\phi(\mathbf{w}_{h}\cdot\mathbf{h}_{n-1}+\mathbf{w}_{x}\cdot\mathbf{x}_{n}).$

From these we see that we need to express
$\partial\mathcal{L}_{n}/\partial\mathbf{w}_{p}$,
$\partial\mathcal{L}_{n}/\partial\mathbf{w}_{h}$, and
$\partial\mathcal{L}_{n}/\partial\mathbf{w}_{x}$ as known values to train the
network. $\partial\mathcal{L}_{n}/\partial\mathbf{w}_{p}$ is relatively easy;
via the chain rule, and the fact that $\partial
p_{n}/\partial\mathbf{w}_{p}=\varphi^{\prime}\mathbf{h}_{n}^{T}$:

$\displaystyle\frac{\partial\mathcal{L}_{n}}{\partial\mathbf{w}_{p}}$
$\displaystyle=\frac{\partial\mathcal{L}_{n}}{\partial p_{n}}\frac{\partial
p_{n}}{\partial\mathbf{w}_{p}},$
$\displaystyle=\frac{\partial\mathcal{L}_{n}}{\partial
p_{n}}\odot\varphi^{\prime}\mathbf{h}_{n}^{T}.$ (7)

$\partial\mathcal{L}_{n}/\partial\mathbf{w}_{h}$ is more tricky, so we will go
step by step. We already know that

$\frac{\partial\mathcal{L}_{n}}{\partial\mathbf{w}_{h}}=\frac{\partial\mathcal{L}_{n}}{\partial
p_{n}}\frac{\partial
p_{n}}{\partial\mathbf{h}_{n}}\frac{\partial\mathbf{h}_{n}}{\partial\mathbf{w}_{h}}.$
(8)

However, we see in Fig. 9 that $\mathbf{h}_{n}$ depends on $\mathbf{h}_{n-1}$,
which depends on $\mathbf{h}_{n-2}$ (and so on). We also notice that all the
hidden states depend on $\mathbf{w}_{h}$. We therefore rewrite Eq. 8 to make
this explicit:

$\displaystyle\frac{\partial\mathcal{L}_{n}}{\partial\mathbf{w}_{h}}$
$\displaystyle=\frac{\partial\mathcal{L}_{n}}{\partial p_{n}}\frac{\partial
p_{n}}{\partial\mathbf{h}_{n}}\sum^{n}_{j=1}\frac{\partial\mathbf{h}_{n}}{\partial\mathbf{h}_{j}}\frac{\partial\mathbf{h}_{j}}{\partial\mathbf{w}_{h}},$
$\displaystyle=\frac{\partial\mathcal{L}_{n}}{\partial p_{n}}\frac{\partial
p_{n}}{\partial\mathbf{h}_{n}}\sum^{n}_{j=1}\left(\prod_{i=j+1}^{n}\frac{\partial\mathbf{h}_{i}}{\partial\mathbf{h}_{i-1}}\right)\frac{\partial\mathbf{h}_{j}}{\partial\mathbf{w}_{h}}.$

We can now substitute in some known values:

$\frac{\partial\mathcal{L}_{n}}{\partial\mathbf{w}_{h}}=\frac{\partial\mathcal{L}_{n}}{\partial
p_{n}}\odot\varphi^{\prime}\mathbf{h}_{n}^{T}\sum^{n}_{j=1}\left(\prod_{i=j+1}^{n}\phi^{\prime}\mathbf{w}_{h,i}^{T}\right){\phi^{\prime}\mathbf{h}_{j-1}^{T}}.$
(9)

Finally, $\partial\mathcal{L}_{n}/\partial\mathbf{w}_{x}$ is derived in the
same way as $\partial\mathcal{L}_{n}/\partial\mathbf{w}_{h}$:

$\displaystyle\frac{\partial\mathcal{L}_{n}}{\partial\mathbf{w}_{x}}$
$\displaystyle=\frac{\partial\mathcal{L}_{n}}{\partial p_{n}}\frac{\partial
p_{n}}{\partial\mathbf{h}_{n}}\sum^{n}_{j=1}\left(\prod_{i=j+1}^{n}\frac{\partial\mathbf{h}_{i}}{\partial\mathbf{h}_{i-1}}\right)\frac{\partial\mathbf{h}_{j}}{\partial\mathbf{w}_{x}},$
$\displaystyle=\frac{\partial\mathcal{L}_{n}}{\partial
p_{n}}\odot\varphi^{\prime}\mathbf{h}_{n}^{T}\sum^{n}_{j=1}\left(\prod_{i=j+1}^{n}\phi^{\prime}\mathbf{w}_{h,i}^{T}\right){\phi^{\prime}\mathbf{x}_{j}^{T}}.$
(10)

With $\partial\mathcal{L}_{n}/\partial\mathbf{w}_{p}$,
$\partial\mathcal{L}_{n}/\partial\mathbf{w}_{h}$, and
$\partial\mathcal{L}_{n}/\partial\mathbf{w}_{x}$ in hand we can apply the same
update rule shown in Eq. 5.

Aside from many-to-one encoding, RNNs can produce many predictions given many
inputs, or act similarly to an MLP and produce one or many outputs given a
single input. We will discuss the application of recurrent neural networks to
astronomical data in §5, after we introduce gated recurrent neural networks.

###  4.3 Sidestepping the vanishing gradient problem

In the early 1990s, researchers identified a major issue with the training of
deep neural networks through backpropagation. ref_hochreiter1991 first
formally examined the ‘vanishing gradient’ problem in their diploma thesis
(ref_hochreiter1991 (351), see also later work by ref_bengio1994 (364)). Due
to the vanishing gradient problem, it was widely believed that training very
deep artificial neural networks from scratch via backpropagation was
impossible. In this section we will explore what the vanishing gradient
problem is, and how contemporary end-to-end trained neural networks sidestep
this issue.

First let us remind ourselves of the sigmoid activation function introduced in
Fig. 6:

$-5$$0$$5$$0.5$$1$$\mathbf{x}$$\varphi(\mathbf{x})\quad$$\varphi^{\prime}(\mathbf{x})$

$\varphi(\mathbf{x})=1/(1-e^{-\mathbf{x}}).$ (11)

Eq. 11 and its accompanying plot shows the output of a sigmoid function
$\varphi$ and its derivative $\varphi^{\prime}$, when given an input
$\mathbf{x}$.

Now, let us revisit the weight update rule for the $(L-n)$th layer of a
feedforward MLP (Eq. 4):

$\frac{\partial\mathcal{L}}{\partial\mathbf{w}_{L-n}}=\frac{\partial\mathcal{L}}{\partial
p_{L}}\odot\underbrace{\left(\prod^{n}_{i=1}\varphi^{\prime}_{L-i}\mathbf{w}^{T}_{L-i}\right)}_{\lim\limits_{n\to\infty}\prod^{n}_{i=1}\varphi^{\prime}_{L-i}\mathbf{w}^{T}_{L-i}=\,0}\,\left(\varphi^{\prime}_{L-n}\mathbf{p}^{T}_{L-n-1}\right).$
(12)

If $\varphi^{\prime}$ is typically less than one (as in Eq. 11 and most other
saturating nonlinearities) the product term in the above equation becomes an
issue. In that case, we can see that the product rapidly goes to zero as $n$
(the number of layers) becomes large161616 Likewise, if $\varphi^{\prime}$ is
typically greater than one, the product term rapidly ‘explodes’ to infinity.
This is known as the ‘exploding gradient’ problem, also first identified in
ref_hochreiter1991 (351). . If we study Eq. 9, we can see the same problem
also plagues RNNs as we backpropagate through hidden states:

$\frac{\partial\mathcal{L}_{n}}{\partial\mathbf{w}_{h}}=\frac{\partial\mathcal{L}_{n}}{\partial
p_{n}}\odot\varphi^{\prime}\mathbf{h}_{n}^{T}\sum^{n}_{j=1}\underbrace{\left(\prod_{i=j+1}^{n}\phi^{\prime}\mathbf{w}_{h,i}^{T}\right)}_{\lim\limits_{n\to\infty}\prod^{n}_{i=j+1}\phi^{\prime}\mathbf{w}^{T}_{h,i}=\,0}{\phi^{\prime}\mathbf{h}_{j-1}^{T}}.$
(13)

Let us solidify this issue by reminding ourselves of Eq. 5—the weight update
rule for a network trained through backpropagation:

$\mathbf{w}_{\text{next}}=\mathbf{w}-\eta\frac{\partial\mathcal{L}}{\partial\mathbf{w}}.$
(14)

Combining Eq. 14 and the limits defined in Eq. 12 and Eq. 13 results in the
below weight update rule in the limit $n\to\infty$.

$\lim\limits_{n\to\infty}\mathbf{w}_{\text{next}}=\mathbf{w}.$ (15)

Eq. 15 shows that learning via backpropagation slows as we move deeper into
the network. This problem once again caused a loss of faith in the
connectionist model, ushering in the second AI winter. It took until 2012 for
a new boom to begin. In the following three subsections we will explore some
of the proposed partial solutions to the vanishing gradient problem and show
how they came together to contribute to the current deep learning boom.

####  4.3.1 Non-saturating activation functions

We can see in Eq. 13 and Eq. 12 that if $\varphi^{\prime}=1$ then the product
term does not automatically go to zero or infinity. If this is the case, why
not simply design our activation function around this property? The rectified
linear unit (ReLU; ref_neocognitron, 339, 420) is an activation function that
does precisely this171717 $\text{ReLU}^{\prime}$ is always zero if its inputs
are $<0$, removing any signal for further training. This is known as the
‘dying ReLU’ problem, but is not as big of an issue as it first seems. Since
contemporary deep neural networks are greatly overparameterised (see for
example ref_frankle2018 (499) and other work on the ‘lottery ticket
hypothesis’) backpropagation through the ReLU activation function can act as a
pruning mechanism, creating sparse representations within the neural network
and thus reducing training time even further (ref_glorot2011, 425). :

$-1$$0$$1$$0.5$$1$$\mathbf{x}$$\text{ReLU}(\mathbf{x})\quad$$\text{ReLU}^{\prime}(\mathbf{x})$

$\text{ReLU}(\mathbf{x})=\max(\mathbf{x},0).$ (16)

The gradient of ReLU is unity if the inputs are above zero, exactly the
property we needed to mitigate the vanishing gradient problem. Similar non-
saturating activation functions also share the ReLU gradient’s useful
property, see for example the Exponential Linear Unit, Swish, and Mish
functions in Fig. 6.

####  4.3.2 Graphics processing unit acceleration

If we can speed up training, we can run an inefficient algorithm (such as
backpropagation through saturating activations) to completion in less time.
One way to speed up training is by using hardware that is specifically suited
to the training of neural networks. Graphics processing units (GPUs) were
originally developed to render video games and other intensive graphical
processing tasks. These rendering tasks require a processor capable of massive
parallelism. We have seen in the previous sections that neural networks
trained through backpropagation also require many small weight update
calculations. With this in mind, it is natural to try to accelerate deep
neural networks using GPUs.

In 2004, ref_oh2004anngpu (397) were the first to use GPUs to accelerate an
MLP model, reporting a $20\times$ performance increase on inference with an
‘ATI RADEON 9700 PRO’ GPU accelerated neural network. Shortly after,
ref_steinkrau2005 (401) showed that backpropagation can also benefit from GPU
acceleration, reporting a three-fold performance increase in both training and
inference. These two breakthroughs were followed by a flurry of activity in
the area (e.g. ref_chellapilla2006, 402, 415, 418, 423), culminating in a
milestone victory for GPU accelerated neural networks at ImageNet 2012.
AlexNet (ref_alexnet, 430) won the ImageNet classification and localisation
challenges (ref_russakovsky2015, 462), scoring an unprecedented top-5
classification error of $16.4\%$, and a single object localisation error of
$34.2\%$. In both challenges AlexNet scored over $10\%$ better than the models
in second place. ref_alexnet’s winning network was a CNN (ref_neocognitron,
339) trained through backpropagation (ref_backpropenglish, 338, 346), with
ReLU activation (ref_nair2010relu, 420), and dropout (ref_dropout, 446) as a
regulariser181818 Dropout reduces the amount of neural network
overfitting—where a network performs well on the training set at the expense
of performance on data it has not yet seen. One performs dropout by randomly
removing a set of neurons at each training step, and using all neurons at test
time. This set up essentially trains a large ensemble of sub-models, whose
average prediction outperforms that inferred by a single model. . The
performance increase afforded by GPU accelerated training enabled the network
to be trained from scratch via backpropagation in a reasonable amount of time.
The discovery that it is possible to train a neural network from scratch by
using readily available hardware ultimately resulted in the end of
connectionism’s second winter, and ushered in the Cambrianesque deep learning
explosion of the mid-to-late 2010s and the 2020s (Fig. 10).

19501960197019801990200020102020Publication
date$\displaystyle{10^{3}}$$\displaystyle{10^{7}}$$\displaystyle{10^{11}}$$\displaystyle{10^{15}}$$\displaystyle{10^{19}}$$\displaystyle{10^{23}}$$\displaystyle{10^{27}}$$\displaystyle{10^{31}}$Training
compute [FLOPs]First era $\displaystyle\leftarrow$$\displaystyle\rightarrow$
DL eraDoubling every 15 monthsDoubling every 5 months Figure 10: If we plot
the total number of floating point operations (FLOPs) required to train a
neural network model, and compare it to the model’s publication date, we can
see a change in trend at around 430. This corresponds to the popularisation of
GPU accelerated training of very deep neural networks, with 430 demarcating
AI’s ‘Deep Learning Era’ and the beginning of astronomy’s second wave of
connectionism (§5). Data is taken from ref_sevilla2022 (624).

####  4.3.3 Gated recurrent neural networks and residual networks

The long short-term memory unit (LSTM; ref_lstm, 377, 383)191919 Compare also
the gated recurrent unit (GRU; ref_gru, 440).  mitigates the vanishing
gradient problem by introducing a new hidden state, the ‘cell state’
$(\mathbf{c}_{n})$, to the standard RNN architecture. This cell state allows
the network to learn long range dependencies, and we will show why this is the
case via a brief derivation202020 Here we loosely follow ref_bayer2015 (450,
§1.3.4). . First, as always, let us study Fig. 11 and write down the forward
pass equation for updating the cell state:

$\mathbf{c}_{n}=f(\mathbf{c}_{n-1},\mathbf{h}_{n-1},\mathbf{x}_{n})+g(\mathbf{h}_{n-1},\mathbf{x}_{n})$

where
$f(\mathbf{c}_{n-1},\mathbf{h}_{n-1},\mathbf{x}_{n})=\mathbf{c}_{n-1}\odot\varphi(\mathbf{h}_{n-1},\mathbf{x}_{n})$.
For brevity we define $\varphi_{n}=\varphi(\mathbf{h}_{n-1},\mathbf{x}_{n})$.

Figure 11: A set of sequential data $\mathbf{x}_{n}$ is input into an LSTM
network. Inside the cell  denotes elementwise operations and  denotes neuronal
layers. $\varphi$ is the sigmoid activation function, and Tanh is the
hyperbolic tangent activation function. $\oplus$ is an elementwise addition,
$\odot$ is the Hadamard product, and line mergers are concatenations.
$\mathbf{c}_{n}$ is the cell state, and $\mathbf{h}_{n}$ is the hidden state.

Like the RNN case (Eq. 9 and Eq. 10), we will need to find
$\partial\mathbf{c}_{n}/\partial\mathbf{c}_{n-1}$ to calculate
$\nabla\mathcal{L}$. Therefore,

$\displaystyle\frac{\partial\mathbf{c}_{n}}{\partial\mathbf{c}_{n-1}}$
$\displaystyle=\frac{\partial
f(\mathbf{c}_{n-1},\mathbf{h}_{n-1},\mathbf{x}_{n})}{\partial\mathbf{c}_{n-1}}+\cancelto{0}{\frac{\partial
g(\mathbf{h}_{n-1},\mathbf{x}_{n})}{\partial\mathbf{c}_{n-1}}},$
$\displaystyle=\frac{\partial\mathbf{c}_{n-1}\odot\varphi_{n}}{\partial\mathbf{c}_{n-1}},$
$\displaystyle=\mathbf{c}_{n-1}\,\cancelto{0}{\frac{\partial\varphi_{n}}{\partial\mathbf{c}_{n-1}}}+\cancelto{1}{\frac{\partial\mathbf{c}_{n-1}}{\partial\mathbf{c}_{n-1}}}\,\varphi_{n},$
$\displaystyle=\varphi_{n}.$

Thus, if we want to backpropagate to a cell state deep in the network we must
calculate

$\frac{\partial\mathbf{c}_{n}}{\partial\mathbf{c}_{N}}=\prod_{i=1}^{n-N}\varphi_{i},\quad
n>N.$ (17)

The product term above does not depend on the derivative of a saturating
activation function, and so does not automatically vanish as $N$ goes to
$\infty$. This means that a gradient signal can be carried through the LSTM
cell state without losing amplitude and vanishing212121 Which is great in
theory. In practice, LSTMs still have trouble learning very long range
dependencies due to their reliance on recurrent processing (ref_seq2seq, 447).
Transformer networks (ref_aiayn, 494) are an architecture that uses the
concept of attention to address this issue. We will discuss transformer
networks in §4.4. .

We can use a technique derived from the LSTM to solve our vanishing gradient
problem for deep feedforward neural networks (as studied in §2.2). ref_highway
(464) do this by applying the concept of the LSTM’s cell state to their deep
convolutional ‘highway network’. The highway network uses gated connections to
modulate the gradient flow back through neuronal layers. Later work by
ref_resnet (454) introduces the residual network (ResNet) by taking a highway
network and simplifying its connections. They apply an elementwise addition
(or ‘residual connection’) in place of the highway network’s gated connection
(Fig. 12(a)). One can go even further with residual connections, as ref_unet
(461) demonstrate with their U-Net model. The U-Net combines residual
connections with an autoencoder-like architecture (Fig. 12(b)). The U-Net has
gone on to become the _de facto_ network for many tasks that require an input
and output of the same size (such as segmentation, colourisation, and style
transfer).

(a) A single residual connection is applied within a neural network.

(b) The U-Net, a network that was originally developed to segment biological
imagery uses the residual connection.

Figure 12: The left subfigure shows the residual connection as originally
introduced in ref_resnet (454). The right subfigure shows an application of
the residual connection to an autoencoder like architecture (ref_unet, 461),
in this case colourising an astronomical object. Here, $\mathbf{z}$ is a
compressed shared representation of $\mathbf{x}$ and $\mathbf{\hat{x}}$.

###  4.4 Translation, attention, and transformers

Theoretically, gated RNNs (GRNNs) such as the LSTM can learn very long range
dependencies (see Eq. 17 and its accompanying text). In practice, GRNNs tend
to forget information about distant inputs. This is because the GRNN lacks
unmediated access to inputs beyond the immediate antecedent as a consequence
of its recurrent architecture. The problem is especially apparent in neural
machine translation tasks that require knowledge of an entire sequence to
produce an output, such as language to language translation. Fig. 13 shows
such a sequence to sequence (Seq2Seq; ref_seq2seq, 447) model. Seq2Seq
translates between two sets of sequential data by sharing a hidden state
between two GRNN units. In Fig. 13 we can see that the shared information is
bottlenecked by the hidden state. Therefore, to resolve the GRNN ‘forgetting
problem’ we must find a way to avoid any recursion, or serial processing of
input and output. We can do this by providing the neural network access to all
input while it is calculating an output. This was the primary motivation
behind the transformer architecture (ref_bahdanau2014, 437, 494).

Figure 13: A sequence to sequence (Seq2Seq; ref_seq2seq, 447) model. A
sequence $\mathbf{x}$ is input into a GRNN. The final hidden state
($\mathbf{h}$) of the input network is then passed into a second GRNN. The
second GRNN then unrolls to predict an output sequence $\mathbf{p}$. Due to
the hidden state acting as a intermediary, $\mathbf{x}$ and $\mathbf{p}$ need
not be of equal length.

Modern transformer architectures consist of a series of self-attention layers
interspersed with other layer types222222 In the original transformer
formulation described in ref_aiayn (494), the network consisted of a connected
‘encoder’ and ‘decoder’ section much like a Seq2Seq model (Fig. 13). Later
work has found this to be an unnecessary complication. For example, the
generative pretrained transformer (GPT) 2 and 3 models (ref_radford2019gpt2,
527, 538) consist of only decoder layers, and the bidirectional encoder
representations from transformers (BERT) model consists of only encoder layers
(ref_bert, 516). . Self-attention as described in ref_aiayn (494) is shown in
Fig. 14. Intuitively, it captures the relationships between quanta within a
data input. To perform self-attention we first take an input sequence

$\mathbf{x}=\begin{bmatrix}x_{1}&x_{2}&\cdots&x_{n}\\\ \end{bmatrix},$

where $\mathbf{x}$ can be any sequence, such as a sentence, a variable star’s
time series, or an unravelled galaxy image232323 One can go very general with
this, as DeepMind demonstrated with their ‘Gato’ transformer model
(ref_reed2022gato, 622). Gato can predict sequences for myriad tasks, from
operating a physical robotic arm, to completing natural language sentences, to
playing Atari games. . This sequence has a maximum length ($n$) that must be
defined at train time, but we can process shorter sequences by masking out any
surplus values so that they do not affect the loss. Here we will follow the
literature and refer to $[x_{1},\ldots,x_{n}]$ as tokens. As we can see in
Fig. 14 the input is passed through a trainable pair of weight matrices
$\mathbf{Q}$ (or ‘query’) and $\mathbf{K}$ (or ‘key’). The output matrices
$\mathbf{q}$ and $\mathbf{k}^{\dagger}$ are then multiplied together to yield

Figure 14: An input $(\mathbf{x})$ is fed into a self-attention mechanism. The
weights used to produce the query $(\mathbf{q})$, key $(\mathbf{k})$, and
value $(\mathbf{v})$ matrices are learnt via backpropagation. Here the learnt
weights are denoted as the capitalised versions of their child matrices.
$\mathbf{q}$ and $\mathbf{k}$ are normalised and multiplied together, and a
softmax nonlinearity ($\varsigma$) is applied. Finally, $\mathbf{v}$ is
multiplied with output of the upper path and the final output is fed forward
to the next neuronal layer. $\bigotimes$ denotes a matrix multiplication.

$(\mathbf{Q}\cdot\mathbf{x})(\mathbf{K}\cdot\mathbf{x})^{\dagger}=\mathbf{q}\mathbf{k}^{\dagger}=\begin{bmatrix}Q_{1}x_{1}K_{1}x_{1}&Q_{1}x_{1}K_{2}x_{2}&\cdots&Q_{1}x_{1}K_{n}x_{n}\\\
Q_{2}x_{2}K_{1}x_{1}&Q_{2}x_{2}K_{2}x_{2}&\cdots&Q_{2}x_{2}K_{n}x_{n}\\\
\vdots&\vdots&\ddots&\vdots\\\
Q_{n}x_{n}K_{1}x_{1}&Q_{n}x_{n}K_{2}x_{2}&\cdots&Q_{n}x_{n}K_{n}x_{n}\end{bmatrix}.$
(18)

We can see that Eq. 18 describes the relationships between tokens within
$\mathbf{x}$. For example, if $x_{1}$ is similar semantically to $x_{2}$, we
would expect $Q_{1}x_{1}K_{2}x_{2}$ and $Q_{2}x_{2}K_{1}x_{1}$ to have a high
value. We then normalise $\mathbf{q}\mathbf{k}^{\dagger}$ to mitigate
vanishing gradients242424See Footnote 16., and apply a softmax non-linearity
so that the maximum weighting (or similarity) is one, and the similarity
values sum to unity.

Meanwhile, the input sequence $\mathbf{x}$ is passed through the neuronal
layer $\mathbf{V}$, resulting in a weighted representation $\mathbf{v}$:

$\mathbf{V}\cdot\mathbf{x}=\mathbf{v}=\begin{bmatrix}V_{1}x_{1}&V_{2}x_{2}&\cdots&V_{n}x_{n}\\\
\end{bmatrix}.$

$\mathbf{v}$ is multiplied with the similarity matrix
$\varsigma(\mathbf{q}\mathbf{k}^{\dagger}/\sqrt{n})$. This process weighs
similar tokens within the sequence higher, increasing their relative
importance in later neuronal layers.

We will use an astronomical example to solidify our understanding of the self-
attention mechanism. Let us assume that our self-attention mechanism is
attending to a natural language caption describing a galaxy’s morphology that
has been provided by a citizen scientist. The caption could be something like:

> $\mathbf{x}=$ ‘ _A barred galaxy with five spiral arms_ ’,

with each word acting as a separate token.

Figure 15: We can think of $\mathbf{q}\mathbf{k}^{\dagger}$ within self-
attention as a graph of relationships between a prompt and itself. Each of the
edges in this graph represents the weight shared between a pair of tokens in
the input sequence.

Let us imagine that we put this prompt into our self-attention mechanism:

$(\mathbf{Q}\cdot\mathbf{x})(\mathbf{K}\cdot\mathbf{x})^{\dagger}=\mathbf{q}\mathbf{k}^{\dagger}=\bordermatrix{&\text{A}&\text{barred}&\text{galaxy}&\text{with}&\text{five}&\text{spiral}&\text{arms}\cr\text{A}&0.7&0.2&0.1&0.0&0.0&0.0&0.0\cr\text{barred}&0.2&0.5&0.3&0.0&0.0&0.0&0.0\cr\text{galaxy}&0.1&0.3&0.5&0.1&0.0&0.0&0.0\cr\text{with}&0.0&0.0&0.1&0.8&0.1&0.0&0.0\cr\text{five}&0.0&0.0&0.0&0.1&0.6&0.3&0.0\cr\text{spiral}&0.0&0.0&0.0&0.0&0.3&0.5&0.2\cr\text{arms}&0.0&0.0&0.0&0.0&0.0&0.2&0.8}.$

We can see that in the above matrix higher values have been assigned to pairs
of words that are more closely related within the sentence. For example, the
weight between ‘ _barred_ ’ and ‘ _galaxy_ ’ is relatively high (0.3), as the
term ‘ _barred_ ’ describes a feature of galaxy. Similarly, the weight between
‘ _five_ ’ and ‘ _spiral_ ’ is also high (0.3), as these words together define
the number of spiral arms in the galaxy. Conversely, lower weights have been
assigned to word pairs that are less related, such as ‘ _A_ ’ and ‘ _with_ ’
(0.0). As shown in Fig. 15, one can think of these relationships between
tokens within our sequence as a learnt mathematical graph252525 This view
demonstrates that transformers can be thought of as a class of graph neural
network—a network that is tasked with learning the relationships between nodes
in a graph. One can also approach this task with a feed-forward neural network
(§2.2; ref_gori2005, 399), convolutional architecture (§4.1; ref_bruna2013,
431, 487), or with a recurrent architecture (§4.3.3; ref_li2015, 457). . Now
that we have calculated $\mathbf{q}\mathbf{k}^{\dagger}$, we can use this
matrix to weigh our example sentence as shown in Fig. 14. This weighting gives
the subsequent layers in our neural network an awareness of the relationships
between the tokens in our sequence.

##  5 Astronomy’s second wave of connectionism

Compared to classical connectionist approaches262626This includes most MLP
applications in astronomy, see §3. deep learning as outlined in §4 does not
require an extraction of emergent parameters to train its models. CNNs in
particular are well suited to observing raw information within image-based
data. Likewise, RNNs are well suited to observing the full raw information
within a time series. Astronomy is rich with both types of data, and in this
section we will review the history of the application of CNN, RNN, and
transformer models to astronomical data.

###  5.1 Convolutional neural network applications

It did not take long after ref_alexnet (430) established CNNs as the _de
facto_ image classification network for astronomers to take notice: in 449
they were applied in the search for pulsars (ref_zhu2014, 449) as part of an
ensemble of methods. ref_zhu2014 (449) found that their ensemble was highly
effective, with 100% of their test set pulsar candidates being ranked within
the top 961 of the 90 008 test candidates. Shortly after, ref_hala2014 (443)
described the use of one dimensional CNNs for a ternary classification
problem. They find that their model is capable of classifying 1D spectra into
quasars, galaxies, and stars to an impressive accuracy. CNNs have been also
been extensively used in galaxy morphological classification. First on the
scene was ref_dieleman2015 (452). They used CNNs to classify galaxy morphology
parameters as defined in the Galaxy Zoo dataset (ref_galzoo, 421) from galaxy
imagery. They observed their galaxies via the SDSS, and found a 99% consensus
between the Galaxy Zoo labels, and the CNN classifications. ref_huertas2015
(455) showed that ref_dieleman2015’s CNN is equally applicable to
morphological classification of galaxies in the CANDELS fields
(ref_koekemoer2011, 426). Likewise, ref_aniyan2017 (483) showed that CNNs are
capable of classifying radio galaxies. The combined work of ref_dieleman2015
(452), ref_huertas2015 (455), and ref_aniyan2017 (483) confirms that CNNs are
equally applicable to visually dissimilar surveys, with little-to-no
modification. Looking a little further afield, ref_wilde2022 (636) used a deep
CNN model to classify simulated lensing events. They also applied some
interpretability techniques to their data, using occlusion mapping
(ref_zeiler2014, 448), gradient class activation mapping (ref_selvaraju2016,
477), and Google’s DeepDream to prove that the CNN was indeed classifying via
observing the gravitational lenses. Alternative CNN models have also been
used, such as the U-Net (Fig. 12(b)). The U-Net was initially developed to
segment biological imagery (ref_unet, 461). Its first use in astronomy was
related: ref_akeret2017 (482) use a U-Net (ref_unet, 461) CNN to isolate via
segmentation, and ultimately remove, radio frequency interference from radio
telescope data. Likewise, ref_berger2019 (514) used a three dimensional U-Net
(V-Net; ref_milletari2016, 475) to predict and segment out galaxy dark matter
haloes in simulations, and ref_aragon2019 (513) used a V-Net to segment out
the cosmological filaments and walls that make up the large scale structure of
the Universe. ref_hausen2020 (547) demonstrate that a U-Net is capable of
performing pixelwise semantic classification of objects in HST/CANDELS
imagery, thus proving that U-Nets are capable of useful work directly within
large imaging surveys, particularly in the deblending of overlapping objects,
which is a perennial challenge in deep imaging. The U-Net in ref_lauritsen2021
(574) is used to superresolve simulated submillimetre observations. They found
that the U-Net could successfully do this when using a loss comprising of the
L1 loss, and a custom loss that measures the distance between predicted and
ground truth point sources. ref_choma2018 (498) were the first to demonstrate
that graph convolutional neural networks (GCNNs) are useful within
astronomical context. They showed that their 3D GCNN could classify signals
from the IceCube neutrino observatory, and found that it outperformed both a
classical method, and a standard 3D CNN. ref_villanueva2021 (593, 631)
demonstrated that EdgeNet—a class of GCNN—can estimate halo masses when given
the positions, velocities, stellar masses, and radii of the host galaxies
(ref_wang2018graph, 510). The authors also demonstrated that EdgeNet can
estimate the halo masses of both Andromeda and the Milky Way. We must conclude
from the studies described in this subsection that CNNs are effective
classifiers and regressors of image-based astronomical data.

###  5.2 Recurrent neural network applications

RNNs were first applied in astronomy very close to home; ref_aussem1994 (363)
predicted atmospheric seeing for observations from ESO’s Very Large Telescope,
and the prediction of geomagnetic storms given data on the solar wind was also
explored in the mid-to-late 1990s and early 2000s (ref_wu1996 (375),
ref_lundstedt2002 (388), and other work from the same group;
ref_vassiliadis2000 (384)).

The first use of RNNs for classification in astronomy was carried out in a
prescient study by ref_brodrick2004 (393). They describe the use of an RNN-
like Elman network (ref_elman1990, 348). Their RNN was tasked with the search
for artificially generated narrowband radio signals that resemble those that
may be produced by an extraterrestrial civilisation. They found that their
model had a test set accuracy of 92%, suggesting that RNNs could be a useful
tool in the search for extraterrestrial intelligence. More than a decade after
ref_brodrick2004 (393), ref_charnock2017 (484) used an LSTM (Fig. 11) to
classify simulated supernovae. They describe two classification problems. One,
a binary classification between type-Ia and non type-Ia supernovae, and the
other a classification between supernovae types I, II, and III. For their best
performing model they report an accuracy of more than 95% for their binary
classification problem, and an accuracy of over 90% for their trinary
classification. This study cemented the usefulness of RNNs for classification
problems in astronomy. ref_charnock2017 (484) was followed by numerous
projects studying the use of RNNs for classification of time series
astronomical data. A non-exhaustive list of modern RNN use in astronomy
includes: stochastically sampled variable star classification (ref_naul2018,
503); exoplanet instance segmentation (ref_gonzalez2018, 500); variable
star/galaxy sequential imagery classification (ref_carrascodavis2019, 515);
and gamma ray source classification (ref_finke2021, 569). We must conclude
from these studies that RNNs are effective classifiers of astronomical time
series, provided that sufficient data is available.

Of course, recurrent networks are not limited to classification; they can also
be used for regression problems. First, ref_weddell2008 (413) successfully
used an echo state network (ref_jaeger2004, 396) to predict the point spread
function of a target object in a wide field of view. ref_capizzi2012 (429)
used an RNN to inpaint missing NASA Kepler time series data for stellar
objects. They found that their model could recreate the missing time series to
an excellent accuracy, suggesting that the RNN could internalise information
about the star it was trained on. As in the classification case, research into
the use of RNNs for regression problems picked up massively in the late 2010s,
and here we will highlight a selection of these studies that represent the
range of RNN use cases. ref_shen2017 (529) used both an LSTM and an
autoencoder based RNN to denoise gravitational wave data, and
ref_morningstar2019 (525) used a recurrent inference machine to reconstruct
gravitationally lensed galaxies. ref_liu2019 (522) used an LSTM to predict
solar flare activity. From these studies, similarly to the classification case
above, we can once again conclude that RNNs are effective regressors of
astronomical time series.

RNNs have also been used in cases that are a little more unconventional. For
example, ref_kugler2016 (473) used an autoencoding RNN (specifically an echo
state network) to extract representation embeddings of variable main sequence
stars. They find that these embeddings capture some emergent properties of
these variable stars, such as temperature, and surface gravity, suggesting
that clustering within the embedding space could result in semantically
meaningful variable star classification. We will revisit this line of research
when we explore representation learning within astronomy in detail in §8. An
example of more drastic cross-pollination between ideas within deep learning
and those within astronomy is ref_smith2021 (587). They use an encoder-decoder
network comprising of a CNN encoder and RNN decoder to predict surface
brightness profiles of galaxies. This class of neural network was previously
used extensively within natural language image captioning, and by treating
surface brightness profiles as ‘captions’ their model was capable of
prediction over $100\times$ faster than the previous classical, human-agent
based method.

###  5.3 Transformer applications

Although initially used for natural language, transformers have also been
adapted for use in imagery, first by ref_parmar2018imtrans (504), and also in
ref_dosovitskiy2020vit (543). To our knowledge, transformers have not yet been
applied to astronomical imagery, but they have started to find use in time-
series astronomy. ref_donosoolivia2022 (642) used BERT (ref_bert, 516) to
generate a representation space for light curves in a self-supervised manner.
ref_morvan2022 (616) use an encoding transformer to denoise light curves from
the Transiting Exoplanet Survey Satellite (TESS; ref_ricker2015, 460), and
show that the denoising surrogate task results in an expressive embedding
space. ref_pan2022 (618) also use a transformer model to analyse light curves
for exoplanets. Transformers have taken the fields of natural language
processing and computer vision by storm (§9), and so if we extrapolate from
trends in other fields we expect to see many more examples of transformers
applied to astronomical use cases in the near future. We will revisit the
transformer architecture in the context of foundation models
(ref_bommasani2021, 565, and references therein), and their possible future
astronomical applications in §9.

###  5.4 A problem with supervised learning

Supervised learning requires a high quality labelled dataset to train a neural
network. In turn, these datasets require labourious human intervention to
create, and so supervised data is in short supply. One can avoid this issue by
prompting the deep learning model to gather semantic information from entirely
unlabelled data. This learnt semantic information can then be accessed through
a hidden descriptive ‘latent space’, and then used for downstream tasks like
data generation, classification, and regression. Indeed, all of the networks
described previously in this review can be repurposed for non-supervised
tasks, and in §6 and §7 we will explore some deep learning frameworks that do
not require supervision.

##  6 Deep generative modelling

In this section we discuss generative modelling within the context of
astronomy. Unlike discriminative models, generative models explicitly learn
the distribution of classes in a dataset (Fig. 16). Once we learn the
distribution of data, we can use that knowledge to generate new synthetic data
that resembles that found in the training dataset. In the following
subsections we will explore in detail three popular forms of deep generative
model: the variational autoencoder (§6.1); the generative adversarial network
(§6.2); and the family of score-based (or diffusion) models (§6.3). Finally,
in §8 we discuss applications of deep generative modelling in astronomy.

Figure 16: Here we show a possible latent space representation of a set of
galaxies and a set of stars. A latent (or embedding) space is a compressed
representation of a set of objects where similar objects are clustered closer
together than dissimilar objects. While this space is often highly
dimensional, here we project our latent space onto two dimensions for
visualisation purposes. On the left we see a generative model attempting to
learn the probability distributions of the latent representation of a dataset
that contains a set of galaxies and a set of stars. On the right is a
discriminative model, which is attempting to learn the boundary that separates
the star and galaxy types.

###  6.1 (Variational) autoencoders

Autoencoders have long been a neural network architectural staple. In a sister
paper to backpropagation’s populariser, ref_rumelhart1986autoencoder (341)
demonstrate backpropagation within an autoencoder. Fig. 17 demonstrates the
basic neural network autoencoder architecture. An autoencoder is tasked with
recreating some input data, squeezing the input information ($\mathbf{x}$)
into a bottleneck latent vector ($\mathbf{z}$) via a neural network
$q(\mathbf{\mathbf{z}}|\mathbf{x})$. $\mathbf{z}$ is then expanded to an
imitation of the input data ($\mathbf{\hat{x}}$) by a second neural network
$p(\mathbf{\hat{x}}|\mathbf{z})$. The standard autoencoder is trained via a
reconstruction loss; $\mathcal{L}_{R}(\mathbf{x},\mathbf{\hat{x}})$, where
$\mathcal{L}_{R}(\mathbf{x},\mathbf{\hat{x}})$ measures the difference in
pixelspace between $\mathbf{x}$ and $\mathbf{\hat{x}}$.

Figure 17: An autoencoder (ref_rumelhart1986autoencoder, 341) attends to an
image of a black hole. $\mathbf{z}$ is a latent vector and $\mathbf{x}$ is a
sample from a training set. The encoder, $q$ learns to encode the incoming
data into a latent vector while the decoder $p$ takes as input $\mathbf{z}$
and attempts to recreate $\mathbf{x}$.

Naïvely, one would think that once trained, one could ‘just’ sample a new
latent vector, and produce novel imagery via the decoding neural network
$p(\mathbf{\hat{x}}|\mathbf{z})$. We cannot do this, as autoencoders trained
purely via a reconstruction loss have no incentive to produce a smoothly
interpolatable latent space. This means we can use a standard autoencoder to
embed and retrieve data contained in the training set, but cannot use one to
generate new data. To generate new data we require a smooth latent space,
which variational autoencoders (VAEs; Fig 18) produce by design (ref_vae,
433).

Figure 18: A variational autoencoder (ref_vae, 433) operates on a spiral
galaxy. $\mathbf{z}$ is a latent vector and $\mathbf{x}$ is a sample from the
training set. The encoder, $q$ learns to compress the incoming data into a
latent vector that encodes the normal distribution. The decoder $p$ takes as
input $\mathbf{z}$ and attempts to recreate $\mathbf{x}$.

A VAE differs from the standard autoencoder by enforcing a spread in each
training set samples’ latent vector. We can see in Fig. 18 how this is done;
instead of directly predicting $\mathbf{z}$ the encoder $q$ predicts two
vectors, $\bm{\mu}$ and $\bm{\sigma}$. $\mathbf{z}$ is then sampled
stochastically via the equation

$\mathbf{z}=\bm{\mu}+\bm{\sigma}\odot\bm{\epsilon},$ (19)

where $\odot$ is the Hadamard product, and $\bm{\epsilon}$ is noise generated
externally to the neural network graph272727 To avoid breaking the
backpropagation chain the VAE injects noise via an external parameter,
$\bm{\epsilon}$. This is described in ref_vae (433) as the ‘reparameterisation
trick’. . This spread results in similar samples overlapping within the latent
space, and therefore we end up with a smooth latent space that we can
interpolate through. However, currently there is no incentive for the neural
network to provide a coherent, compact global structure in the latent space.
For that we require a regularisation term in the loss. This regularisation is
provided via the Kullback-Leibler (KL) divergence, which is a measure of the
difference between two probability distributions. A standard VAE uses the KL
divergence to push the latent distribution towards the standard normal
distribution, incentivising a compact, continuous latent space. Hence, the
final VAE loss is a combination of the reconstruction loss and KL divergence:

$\mathcal{L}_{\text{VAE}}=\mathcal{L}_{R}(\mathbf{x},\mathbf{\hat{x}})+\text{KL}(q(\mathbf{x}|\mathbf{z})\|\rho),$
(20)

where $\rho$ is some prior. In a standard VAE
$\rho=\mathcal{N}(\mathbf{0},\mathds{1})$.

In practice, VAEs are able to generate smooth and coherent samples, as they
model the data distribution explicitly, which also means that we can perform
latent space arithmetic on the latent vector—such as interpolation,
reconstruction, and anonomly detection (ref_vae, 433). Their explicit learning
of the latent vector ($\mathbf{z}$) means that they can trivially be
repurposed for semi-supervised, self-supervised, and supervised downstream
tasks by manipulating $\mathbf{z}$ (ref_regier2015deep, 458, 559). However,
the quality of samples generated by VAEs is lower than that of generative
adversarial networks or score-based generative models (ref_dosovitskiy2016,
469). This reduction in quality is due to the VAE’s simple posterior
$q(\mathbf{z}|\mathbf{x})$, but one can mitigate this shortcoming by
iteratively approaching a more complex posterior282828 Interestingly, this
iterative approximation is similar to the approach used in the training of
score-based generative models and diffusion models (ref_zhao2017, 495), and
the similarities between the training methods of state-of-the-art in VAE
models and SBGMs are striking. For example, the Vector-Quantised VAE, Very
Deep VAE, and the Nouveau VAE all use a heirarchical architecture that
iteratively injects latent codes that are used to produce finer and finer
detail in the generated image (ref_oord2017, 490, 560, 542). . To regularise
the latent space, VAEs require an assumption of the prior distribution which
requires some knowledge of the dataset, although often this is can be set as
‘just’ a normal distribution as shown in Eq. 20.

###  6.2 Generative adversarial networks

Generative adversarial networks (GAN; ref_gan, 441) can be thought of as a
minimax game between two competing neural networks. If we anthropomorphise we
can gain an intuition for how a GAN learns: let us imagine an art forger, and
an art critic. The forger wants to paint paintings that are similar to famous
expensive works, and needs to fool the critic when selling these paintings.
Meanwhile, the critic wants to ensure that no reproductions are sold and so
they need to accurately determine whether any painting is an original or a
reproduction. At first, our forger is a poor painter, and so the critic can
easily identify our forger’s works. However, the forger learns from the
critic’s choices and produces more realistic paintings. As the forger’s
paintings improve, the critic also learns better methods for detecting
forgeries. This minimax game incentivises the critic to keep improving their
classifications, and the forger to keep improving their painting. If this
continues, we get to a point where the forger’s works are indiscernible from
the real thing—the forger has learnt to perfectly mimic the dataset! In a GAN,
we name the critic the discriminator ($D$), and we name the forger the
generator ($G$).

(a) A typical GAN according to ref_gan (441). $\mathbf{z}$ is a noise vector,
$\mathbf{x}$ is a sample from the training set. The discriminator learns to
classify the incoming images as either fake or real, and the generator learns
to fool the discriminator by producing realistic fakes.

(b) A Pix2Pix-like model with a U-Net generator (ref_pix2pix, 471, 461). The
discriminator learns to classify the incoming image tuples as either fake or
real. Meanwhile, the generator learns to fool the discriminator by
approximating the colourisation function mapping $\mathbf{x}\to\mathbf{y}$.
Line mergers denote channel-wise concatenations.

Figure 19: The GAN and Pix2Pix models.

In ref_gan’s original GAN formulation (Fig. 19(a)), $G$ and $D$ are neural
networks (typically CNNs, although other architectures can be used) that
compete during training in a minimax game where $G$ aims to maximise the
probability of $D$ mispredicting that a generated datapoint is sampled from
the real dataset. $G$ takes as input a randomly sampled latent vector z, and
outputs a synthetic datapoint $G(\textbf{z})$. $D$ takes either this synthetic
datapoint, or a real datapoint x, and outputs $D(G(\textbf{z}))$ or
$D(\textbf{x})$. This output is the probability that the datapoint is drawn
from the real dataset. To train the network we can write the GAN adversarial
loss like so:

$\displaystyle\mathcal{L}_{D}$
$\displaystyle=-(\mathbb{E}_{\textbf{x}}[\log(D(\textbf{x}))]+\mathbb{E}_{\textbf{z}}[\log(1-D(G(\textbf{z})))]),$
$\displaystyle\mathcal{L}_{G}$
$\displaystyle=\mathbb{E}_{\textbf{z}}[\log(1-D(G(\textbf{z})))],$

where here we attempt to minimise both $\mathcal{L}_{D}$ and
$\mathcal{L}_{G}$. In practice we train the networks by alternating freezing
the weights of $G$ and backpropagating $\mathcal{L}_{D}$, and then freezing
the weights of $D$ and backpropagating $\mathcal{L}_{G}$ for each training
batch. In this way the networks’ weights are updated to follow
$\nabla_{\mathbf{w}}\mathcal{L}_{G}$ and $\nabla_{\mathbf{w}}\mathcal{L}_{D}$
downwards until the distribution of $G(\textbf{z})$ closely resembles that of
the real dataset. Once trained, $G$ can be used to generate entirely novel
synthetic data that closely resembles (but is not identical to) the training
set data.

One can condition a GAN to guide the network towards a desired output image
(ref_mirza2014, 444). To do this we alter the adversarial loss so that it is
conditioned on a label y:

$\displaystyle\mathcal{L}_{D}$
$\displaystyle=-(\mathbb{E}_{\mathbf{x}}[\log(D(\textbf{x}|\textbf{y}))]+\mathbb{E}_{\mathbf{z}}[\log(1-D(G(\textbf{z}|\textbf{y})))]),$
$\displaystyle\mathcal{L}_{G}$
$\displaystyle=\mathbb{E}_{\mathbf{z}}[\log(1-D(G(\textbf{z}|\textbf{y})))].$

As an example, if we set $\mathbf{y}$ as the redshift of the galaxies in the
training set, we could use a conditional GAN to guide the network to generate
galaxies of a certain redshift. Furthermore, we are not restricted to
conditioning single values; GANs can also be conditioned on entire images. In
Fig. 19(b) we see that the GAN adversarial loss can be used to translate
between image domains (ref_pix2pix, 471). In ref_pix2pix’s Pix2Pix model, the
generator takes as input an image $\mathbf{x}$, and attempts to produce a
related image $\mathbf{y}$. Meanwhile, the discriminator attempts to discern
whether the $(\mathbf{x},\mathbf{y})$ pair that it is given is sampled from
the training set, or the generator. Otherwise, Pix2Pix is trained in the same
way as the standard GAN.

GANs are capable of generating high-quality, sharp, and realistic samples
(ref_brock2018, 497, 646). They have long been a sweetheart of the deep
generative learning community, having been used for various state-of-the-art
applications, such as data embedding (e.g. ref_cheng2016, 467), style transfer
(e.g. ref_karras2018, 501), superresolution (e.g. ref_ledig2016, 474), and
image inpaining and object removal (e.g. ref_yu2018, 512). Unfortunately
however, GANs have some downsides. They are quite difficult to train;
maintaining the balance between the generator and discriminator networks is
challenging and requires careful finetuning (ref_weng2019, 536). $G$ and $D$
must work in tandem and one cannot overpower the other or learning will cease.
One of the most famous symptoms of this imbalance is mode collapse, where $G$
only generates a limited variety of samples that reliably fool $D$. This
instability during training makes it quite a time-consuming task to find a
stable network architecture if one is designing a GAN themselves. Finally, the
GAN adversarial losses are relative and so are not representative of the image
quality. This is not the case for the VAE and SBGM families of models.

###  6.3 Score-based generative modelling and diffusion models

Diffusion models were introduced by ref_sohldickstein2015 (463) and were first
shown to be capable of producing high quality synthetic samples by ref_ho2020
(548). Diffusion models are part of a family of generative deep learning
models that employ denoising score matching via annealed Langevin dynamic
sampling (first explored by ref_hyvarinen2005 (400, 427). More recent work can
be found in ref_song2020 (558, 548, 550, 572, 588)). This family of score-
based generative models (SBGMs) can generate imagery of a quality and
diversity surpassing state of the art GAN models (ref_gan, 441), a startling
result considering the historic disparity in interest and development between
the two techniques (ref_song2021, 588, 579, 568, 621). SBGMs can super-resolve
images (ref_kadkhodaie2020, 551, 584), translate between image domains
(ref_sasaki2021, 586), separate superimposed images (ref_jayaram2020, 549),
and in-paint information (ref_kadkhodaie2020, 551, 588).

Diffusion models define a diffusion process that projects a complex image
domain space onto a simple domain space. In the original formulation, this
diffusion process is fixed to a predefined Markov chain
$q(\mathbf{x}_{t}\mid\mathbf{x}_{t-1})$ that adds a small amount of Gaussian
noise with each step. As Fig. 20 shows, this ‘simple domain space’ can be
noise sampled from a Gaussian distribution
$\mathbf{x}_{T}\sim\mathcal{N}(\mathbf{0},\mathds{1})$.

Figure 20: It is easy (and achievable without learnt parameters) to add noise
to an image, but more difficult to remove it. Diffusion models attempt to
learn an iterative removal process via training an appropriate neural network,
$p_{\theta}(\mathbf{x}_{t-1}\mid\mathbf{x}_{t})$.

####  6.3.1 Forward process

To slowly add Gaussian noise to our data we define a Markov chain

$q(\mathbf{x}_{0\ldots
T})=q(\mathbf{x}_{0})\prod^{T}_{t=1}q(\mathbf{x}_{t}\mid\mathbf{x}_{t-1}),$

where $\mathbf{x}_{0}$ is an image sampled from the training set. The amount
of noise added per step is controlled with a variance schedule
$\\{\beta_{t}\in(0,1)\\}^{T}_{t=1}$:

$q(\mathbf{x}_{t}\mid\mathbf{x}_{t-1})=\mathcal{N}(\mathbf{x}_{t};\sqrt{1-\beta_{t}}\,\mathbf{x}_{t-1},\,\beta_{t}\mathds{1}).$
(21)

This process is applied incrementally to the input image. Since we can define
the above equation such that it only depends on $\mathbf{x}_{0}$ we can
immediately calculate an image representation $\mathbf{x}_{t}$ for any $t$
(ref_ho2020, 548). If we define $\alpha_{t}=1-\beta_{t}$ and
$\bar{\alpha}_{t}=\prod^{t}_{i=1}\alpha_{i}$:

$\displaystyle\mathbf{x}_{t}$
$\displaystyle=\sqrt{\alpha_{t}}\,\mathbf{x}_{t-1}+\sqrt{1-\alpha_{t}}\,\mathbf{z}_{t-1}$
$\displaystyle=\sqrt{\alpha_{t}\alpha_{t-1}}\,\mathbf{x}_{t-2}+\sqrt{(1-\alpha_{t})+\alpha_{t}(1-\alpha_{t-1})}\,\bar{\mathbf{z}}_{t-2}$
$\displaystyle=\sqrt{\alpha_{t}\alpha_{t-1}\alpha_{t-2}}\,\mathbf{x}_{t-3}+\sqrt{(1-\alpha_{t}\alpha_{t-1})+\alpha_{t}\alpha_{t-1}(1-\alpha_{t-2})}\,\bar{\mathbf{z}}_{t-3}$
$\displaystyle=\ldots$
$\displaystyle=\sqrt{\bar{\alpha_{t}}}\mathbf{x}_{0}+\sqrt{1-\bar{\alpha}_{t}}\mathbf{z},$
(22)

where $\mathbf{z}_{t}\sim\mathcal{N}(\mathbf{0},\mathds{1})$ and
$\bar{\mathbf{z}}$ is a combination of Gaussians. Plugging the above
expression into Eq. 21 removes the $\mathbf{x}_{t-1}$ dependency and yields

$q(\mathbf{x}_{t}\mid\mathbf{x}_{0})=\mathcal{N}(\mathbf{x}_{t};\sqrt{\bar{\alpha}_{t}}\,\mathbf{x}_{0},\,(1-\bar{\alpha}_{t})\mathds{1}).$
(23)

####  6.3.2 Reverse process

Diffusion models attempt to reverse the forward process by applying a Markov
chain with learnt Gaussian transitions. These transitions can be learnt via an
appropriate neural network, $p_{\theta}$:

$\displaystyle p_{\theta}(\mathbf{x}_{0\ldots T})$
$\displaystyle=p(\mathbf{x}_{T})\prod^{T}_{t=1}p_{\theta}(\mathbf{x}_{t-1}\mid\mathbf{x}_{t}),$
$\displaystyle p_{\theta}(\mathbf{x}_{t-1}\mid\mathbf{x}_{t})$
$\displaystyle=\mathcal{N}(\mathbf{x}_{t-1};\bm{\mu}_{\theta}(\mathbf{x}_{t},t),\bm{\Sigma}_{\theta}(\mathbf{x}_{t},t)).$

While $\mathbf{\Sigma}_{\theta}(\mathbf{x}_{t},t)$ can be learnt292929See for
example ref_nichol2021 (579)., the ref_ho2020 (548) formulation fixes
$\mathbf{\Sigma}_{\theta}$ to an iteration-dependent constant
$\sigma_{t}^{2}\mathds{1}$, where $\sigma^{2}_{t}=1-\alpha_{t}$.

By recognising that diffusion models are a restricted class of hierarchical
VAE303030 Denoising autoencoders (§6.1) have an interesting relationship with
score-based generative (or diffusion) models. As a taster, ref_turner2021
(592) reframe diffusion models as a class of hierarchical denoising VAE, and
ref_dieleman2022 (607) show through a brief derivation that diffusion models
optimise the same loss as a denoising autoencoder. , we see that we can train
$p_{\theta}$ by optimising the evidence lower bound (ELBO, introduced in
ref_vae, 433) that can be written as a summation over the Kullback-Leibler
divergences at each iteration step313131 See Appendix B in
ref_sohldickstein2015 (463) and Appendix A in ref_ho2020 (548) for the full
derivation. :

$\displaystyle\mathcal{L}_{\text{ELBO}}$
$\displaystyle=\mathbb{E}_{q}\Bigl{[}D_{\text{KL}}(q(\mathbf{x}_{T}\mid\mathbf{x}_{0})\|p(\mathbf{x}_{T}))+$
$\displaystyle\hskip
50.00008pt\sum_{t>1}D_{\text{KL}}(q(\mathbf{x}_{t-1}\mid\mathbf{x}_{t},\mathbf{x}_{0})\|p_{\theta}(\mathbf{x}_{t-1}\mid\mathbf{x}_{t}))+\log
p_{\theta}(\mathbf{x}_{0}\mid\mathbf{x}_{1})\Bigr{]}.$ (24)

In the ref_ho2020 (548) formulation, the first term in Eq. 24 is a constant
during training and the final term is modelled as an independent discrete
decoder. This leaves the middle summation. Each summand can be written as

$\mathcal{L}(\bm{\mu}_{t},\bm{\mu}_{\theta})=\frac{1}{2\sigma_{t}^{2}}\|\bm{\mu}_{t}(\mathbf{x}_{t},\mathbf{x}_{0})-\bm{\mu}_{\theta}(\mathbf{x}_{t},t)\|^{2},$
(25)

where $\bm{\mu}_{\theta}$ is the neural network’s estimation of the forward
process posterior mean $\bm{\mu}_{t}$. In practice it would be preferable to
predict the noise addition in each iteration step ($\mathbf{z}_{t}$), as
$\mathbf{z}_{t}$ has a distribution that by definition is centred about zero,
with a well defined variance. To this end we can define $\bm{\mu}_{\theta}$ as

$\bm{\mu}_{\theta}(\mathbf{x}_{t},t)=\frac{1}{\sqrt{\alpha_{t}}}\left(\mathbf{x}_{t}-\frac{1-\alpha_{t}}{\sqrt{1-\bar{\alpha}_{t}}}\mathbf{z}_{\theta}(\mathbf{x}_{t},t)\right),$
(26)

and by combining Eqs. 25 and 26 we get

$\displaystyle\mathcal{L}(\mathbf{z}_{t},\mathbf{z}_{\theta})$
$\displaystyle=\frac{1}{2\sigma_{t}^{2}}\Biggl{\|}\frac{1}{\sqrt{\alpha_{t}}}\left(\mathbf{x}_{t}-\frac{1-\alpha_{t}}{\sqrt{1-\bar{\alpha}_{t}}}\mathbf{z}_{t}\right)-\frac{1}{\sqrt{\alpha_{t}}}\left(\mathbf{x}_{t}-\frac{1-\alpha_{t}}{\sqrt{1-\bar{\alpha}_{t}}}\mathbf{z}_{\theta}(\mathbf{x}_{t},t)\right)\Biggr{\|}^{2}$
$\displaystyle=\frac{(1-\alpha_{t})^{2}}{2\sigma_{t}^{2}\alpha_{t}(1-\bar{\alpha}_{t})}\|\mathbf{z}_{t}-\mathbf{z}_{\theta}(\mathbf{x}_{t},t)\|^{2}.$
(27)

ref_ho2020 (548) empirically found that a simplified version of the loss
described in Eq. 27 results in better sample quality. They use a simplified
version of Eq. 27 as their loss, and optimise to predict the noise required to
reverse a forward process iteration step:

$\mathcal{L}(\mathbf{z}_{t},\mathbf{z}_{\theta})=\|\mathbf{z}_{t}-\mathbf{z}_{\theta}(\textbf{x}_{t},\,t)\|^{2},\quad\text{where}~{}\mathbf{x}_{t}=\sqrt{\bar{\alpha}_{t}}\mathbf{x}_{0}+\sqrt{1-\bar{\alpha}_{t}}\mathbf{z}_{t}.$
(28)

By recognising that $\mathbf{z}_{t}=\sigma_{t}^{2}\nabla_{\mathbf{x}_{t}}\log
q(\mathbf{x}_{t}\mid\mathbf{x}_{t-1})$, we see that Eq. 28 is equivalent to
denoising score matching over $t$ noise levels (ref_vincent2011, 427). This
connection establishes a link between diffusion models and other SBGMs (such
as ref_song2019, 531, 558, 550).

To run inference for the reverse process, one progressively removes the
predicted noise $\mathbf{z}_{\theta}$ from an image. The predicted noise is
weighted according to a variance schedule:

$\mathbf{x}_{t-1}=\frac{1}{\sqrt{\alpha_{t}}}\left(\mathbf{x}_{t}-\frac{1-\alpha_{t}}{\sqrt{1-\bar{\alpha}_{t}}}\,\mathbf{z}_{\theta}(\mathbf{x}_{t},t)\right)+\bm{\sigma}_{t}\mathbf{z}.$

If we take
$p(\mathbf{x}_{T})\sim\mathcal{N}(\mathbf{x}_{T};\mathbf{0},\mathds{1})$, we
can use $p_{\theta}$ to generate entirely novel data that are similar—but not
identical to—those found in the training set.

In practice, diffusion models are trained by sampling an integer value of
$t\sim\mathcal{U}(1,T)$, where $T$ is a large value typically in the 1000s. We
then use Eq. 23 to sample an image $\mathbf{x}_{t}$ that has had noise added
to it $t$ times. The model then attempts to predict the exact noise required
to reverse a forward iteration time step—that is, the output of a neural
network323232 Typically a U-Net; see §4.3.3 for more detail.  of the form
$z_{\theta}(\mathbf{z}_{t}|\mathbf{x}_{t-1})$. As shown in Fig. 20, we can
estimate $\mathbf{x}_{t}$ by removing the predicted noise from
$\mathbf{x}_{t-1}$. To optimise the model $\mathbf{z}_{t}$ is compared via Eq.
28 to the actual noise required to reverse the forward iteration, and this is
the loss that is reduced during training. For a detailed astronomical example
with code we direct the reader to ref_smith2022 (628).

####  6.3.3 Denoising diffusion implicit models

ref_ho2020’s diffusion model performs inference at a rate orders of magnitude
slower than single shot generative models like the VAE (§6.1) or the GAN
(§6.2). This is because diffusion models need to sequentially reverse every
step in the forward process Markov Chain. Reducing the inference time for
diffusion models is an active area of research (ref_ajm2021, 572, 576, 633),
and here we will review one proposed solution to the problem; the denoising
diffusion implicit model (DDIM; ref_song2020ddim, 557).

ref_song2020ddim (557, §§3-4) propose the following reparameterisation of Eq.
22:

$\displaystyle\mathbf{x}_{t-1}$
$\displaystyle=\sqrt{\bar{\alpha}_{t-1}}\,\mathbf{x}_{0}+\sqrt{1-\bar{\alpha}_{t-1}-\sigma_{t}^{2}}\,\mathbf{z}_{\theta}^{(t)}+\sigma_{t}\mathbf{z}_{t}$
$\displaystyle=\sqrt{\bar{\alpha}_{t-1}}\,\underbrace{\left(\frac{\mathbf{x}_{t}-\sqrt{1-\bar{\alpha}_{t}}\,\mathbf{z}^{(t)}_{\theta}}{\sqrt{\bar{\alpha}_{t}}}\right)}_{\mathbf{x}_{0}\text{
prediction}}+\underbrace{\sqrt{1-\bar{\alpha}_{t-1}-\sigma_{t}^{2}}\,\mathbf{z}^{(t)}_{\theta}}_{\text{vector
towards
}\mathbf{x}_{t}}+\underbrace{\sigma_{t}\mathbf{z}_{t}}_{\text{noise}},$

where $(t)$ is noted as a superscript to denote the output of the neural
network $z_{\theta}$ at timestep $t$. Intuitively, the first term can be
thought of as the prediction of the input image $\mathbf{x}_{0}$, given an
iteration step $t$. The second term can be thought of as a vector from
$\mathbf{x}_{t-1}$ towards the current iteration step image $\mathbf{x}_{t}$.
The third term is random noise. If we substitute in $\mathbf{x}_{t}$ from Eq.
28 we make this intuition explicit:

$\mathbf{x}_{t-1}=\sqrt{\bar{\alpha}_{t-1}}\,\mathbf{x}_{0}+\sqrt{1-\bar{\alpha}_{t-1}-\sigma_{t}^{2}}\,\frac{\mathbf{x}_{t}-\sqrt{\bar{\alpha}_{t}}\,\mathbf{x}_{0}}{\sqrt{1-\bar{\alpha}_{t}}}+\sigma_{t}\mathbf{z}.$

If we then set $\sigma_{t}=0$, we remove the noise dependency and the forward
process becomes deterministic:

$q_{\text{DDIM}}(\mathbf{x}_{t-1}\mid\mathbf{x}_{t},\mathbf{x}_{0})=\sqrt{\bar{\alpha}_{t-1}}\,\mathbf{x}_{0}+\sqrt{1-\bar{\alpha}_{t-1}}\,\frac{\mathbf{x}_{t}-\sqrt{\bar{\alpha}_{t}}\,\mathbf{x}_{0}}{\sqrt{1-\bar{\alpha}_{t}}}.$
(29)

This means that DDIMs can deterministically map to and from the latent space,
and so inherit all the benefits of this property. For example, two objects
sampled from similar latent vectors share high level properties, latent space
arithmetic is possible, and we can perform meaningful interpolation within
this space. We demonstrate DDIM latent space interpolation in Fig. 21.

Figure 21: Meaningful latent space interpolation via a DDIM model
(ref_song2020ddim, 557, 628). This property comes ‘for free’ with most other
generative models, however the denoising diffusion probabilistic model
(ref_ho2020, 548) requires a tweak to its sampling scheme (Eq. 29).

We can also subsample every $\tau$ number of steps at inference time, where
$\tau$ is a set of evenly spaced steps between $0$ and $T$, the maximum number
of steps in the forward process:

$q_{\text{DDIM}}(\mathbf{x}_{\tau_{i-1}}\mid\mathbf{x}_{\tau_{i}},\mathbf{x}_{0})=\sqrt{\bar{\alpha}_{t-1}}\,\mathbf{x}_{0}+\sqrt{1-\bar{\alpha}_{t-1}}\,\frac{\mathbf{x}_{\tau_{i}}-\sqrt{\bar{\alpha}_{t}}\,\mathbf{x}_{0}}{\sqrt{1-\bar{\alpha}_{t}}}.$
(30)

As shown in ref_song2020ddim (557) this results in acceptable generations with
a $T/\tau$ inference speed up.

SBGMs have emerged as a promising alternative to GANs, VAEs, and other
generative models, showcasing their ability to generate high-quality samples
with a level of detail comparable to that of the previous state-of-the-art
(ref_dhariwal2021, 568, 579, 621). One of the key advantages of SBGMs is how
easy they are to train; they do not inherit any of the instability issues that
plague GANs. However, SBGMs do have their share of weaknesses. For instance,
the SBGM sampling process is computationally expensive and slow. This is
because generating a single sample requires a pass through a learnt Markov
chain (Fig. 20), which can limit their practicality in certain applications.
Finally, diffusion models and other SBGMs have not been as extensively
explored in the deep learning literature as VAEs and GANs (although this is
changing fast!). This leaves their applicability across various domains still
under investigation.

##  7 Representation learning

Self-supervised333333 A model that employs self-supervised learning is one
that obtains a supervisory signal from the data itself. ‘Self-supervised
learning’ as a descriptor has largely superseded the older term ‘unsupervised
learning’. This is because the older term suggests that there is no
supervisory signal at all—but the signal is there, just not explicitly defined
by a human expert!  representation learning has recently exploded in
popularity, with a slew of models being developed in rapid succession (e.g.
ref_simclr, 539, 540, 546, 519, 541, 544). At its core, representation
learning attempts to produce semantically meaningful compressed
representations (or embeddings) of complex highly dimensional data. Aside from
simply being a compression device, these embeddings can also be taken and used
in downstream tasks, like clustering, anomaly detection, or classification.

In this section we will describe two approaches to representation learning
that are popular within astronomy. The first approach uses contrastive
learning as defined by the SimCLR model. The second approach defines and uses
a ‘surrogate task’ (such as autoencoding or next value prediction) to train a
deep learning model, and extracts semantically meaningful representations from
the subsequent trained network.

###  7.1 Contrastive learning

Fig. 22 describes a simple contrastive learning model similar to SimCLR
(ref_simclr, 539). This model takes as input a sample $\mathbf{x}$ from the
training set, and augments it to produce $\mathcal{A}(\mathbf{x})$. This
augmentation is performed in such a way that $\mathcal{A}(\mathbf{x})$ shares
enough semantically meaningful data with $\mathbf{x}$ to belong to the same
class. In the contrastive learning literature
$(\mathbf{x},\mathcal{A}(\mathbf{x}))$ is known as a positive pair. This
positive pair is passed to a Siamese neural network $\Phi$, which projects the
high dimensional input data onto a lower dimensional ‘embedding space’. All
other training set samples are assumed to belong to a different class to
$\mathbf{x}$, and so can be combined with $\mathbf{x}$ to produce ‘negative
pairs’.

(a) Possible application to imagery.

(b) Possible application to sequential data.

Figure 22:  A simple contrastive learning model is applied to both imagery and
sequential data. $\mathcal{A}$ is an augmentation pipeline. For imagery,
$\mathcal{A}$ could consist of random crops, noise addition, and colour
jitter. For sequential data, $\mathcal{A}$ could consist of noise addition,
stochastic temporal shifting, and random data deletion. $\Phi$ is a function
approximator that projects inputs onto an embedding space. $\Phi$ is typically
a neural network: when processing imagery, $\Phi$ could take the form of a
CNN, and when processing sequential data $\Phi$ could be an RNN. The loss
$\mathcal{L}$ measures the distance between the embeddings
$\Phi(\bm{\mathbf{x}})=\mathbf{z}_{i}$ and
$\Phi(\mathcal{A}(\bm{\mathbf{x}}))=\mathbf{z}_{j}$, and we train by
attempting to minimise this distance while maximising the distance between
dissimilar samples.

Once we produce some embeddings we need to define a loss that clusters similar
samples together, while simultaneously pushing away dissimilar samples.
ref_hadsell2006 (404) propose such a loss—the maximum margin contrastive loss:

$\mathcal{L}(\mathbf{z}_{i},\mathbf{z}_{j})=\delta_{y_{i}y_{j}}\,d(\mathbf{z}_{i},\mathbf{z}_{j})+(1-\delta_{y_{i}y_{j}})\max(0,m-d(\mathbf{z}_{i},\mathbf{z}_{j})),$

where $\delta$ is the Kronecker delta, $\mathbf{z}_{i}$ and $\mathbf{z}_{j}$
are embedding vectors343434 All embeddings in this subsection are normalised.
, $y_{i}$ and $y_{j}$ are the class labels for the embedding vectors, and $m$
is the margin. $d$ is a ‘distance metric’ (such as for example the L1 loss)
that reduces to zero in the case where its inputs are identical. If
$\mathbf{z}_{i}$ and $\mathbf{z}_{j}$ are a positive pair, the loss pulls the
embeddings closer, and if they are a negative pair the loss pushes the
embeddings away from each other. The margin imposes an upper distance bound on
dissimilar embeddings.

While useful, the maximum margin contrastive loss does not take into account
the embedding space beyond the pair it is attending to in each training step.
This limitation ultimately results in a less expressive embedding space. The
triplet loss (ref_chechik2010, 417) solves this issue by taking into account
the broader embedding space and simultaneously attracting a positive pair
while repulsing a negative pair with each training step:

$\mathcal{L}(\mathbf{z}_{i},\mathbf{z}_{j},\mathbf{z}_{k})=\max(0,\,d(\mathbf{z}_{i},\mathbf{z}_{j})-d(\mathbf{z}_{i},\mathbf{z}_{k})+m)$
(31)

where $\mathbf{z}_{k}$ is a sampled from a different class to
$\mathbf{z}_{i}$, and $\mathbf{z}_{j}$ is sampled from the same class as
$\mathbf{z}_{i}$.

If we study Eq. 31 we see that it is possible to generalise our loss even
further, taking into account an arbitrary number of negative samples. The
normalized temperature-scaled cross entropy loss (NT-Xent; ref_sohn2016, 479)
does precisely this:

$\mathcal{L}(\mathbf{z}_{i},\mathbf{z}_{j})=-\text{log}\left(\frac{\exp(d(\mathbf{z}_{i},\mathbf{z}_{j})/\mathcal{T})}{\sum^{2N}_{k=1}(1-\delta_{ik})\exp(d(\mathbf{z}_{i},\mathbf{z}_{k})/\mathcal{T})}\right),$
(32)

where $\mathbf{z}_{i}$ and $\mathbf{z}_{j}$ are a positive embedding pair, and
$\mathbf{z}_{i}$ and $\mathbf{z}_{k}$ are a negative pair. $\mathcal{T}$ is a
‘temperature’ hyperparameter introduced in ref_simclr (539) to help the model
learn from hard negatives (negatives closer to the anchor than the comparison
positive, see Fig. 23(b)).

(a) The triplet (Eq. 31) and NT-Xent (Eq. 32) losses simultaneously
incentivise attraction between embeddings sampled from the same class
($\mathbf{z}_{i}$ and $\mathbf{z}_{j}$), and repulsion between embeddings
sampled from different classes ($\mathbf{z}_{i}$ and $\mathbf{z}_{k}$).

(b) Types of negative embeddings. $\mathbf{z}_{i}$ and $\mathbf{z}_{j}$ form a
positive embedding pair. If a negative is closer than the current positive it
is considered a hard negative, if it lies within the margin it is considered a
semi-hard negative, and if it is beyond the margin it is considered an easy
negative.

Figure 23: More information on self-supervised embeddings. Fig. 23(a) depicts
the inner workings of the triplet and NT-Xent losses, and Fig. 23(b) shows the
three possible negative embedding types as described in the literature.

###  7.2 Learning representations via a surrogate task

One can also learn representations via a surrogate task. A surrogate task is
any task that is unrelated to the network’s final use. However, in the process
of learning to perform the surrogate task, the network learns what is
important, and what is unimportant about data within the training set. This
information can then be extracted in the form of learnt representations. If
the surrogate task is general enough, these representations will contain
useful semantic information about the items in the dataset, and can then be
used for downstream applications.

Let us concretise this process by revisiting an example that we previously
discussed in §4.2. Let us imagine we have a large set of galaxy rotation
curves that we want to extract embeddings from. We could train an LSTM model
(Fig. 24) on the task of predicting the next item in the rotation curve, with
the model only having access to the previous items in the profile. Once the
LSTM model is trained on this task, we can feed in a full, new rotation curve,
and repurpose the final hidden state as a representative embedding. Note that
this set up does not rely on any external labels, only on the rotation curve
itself353535 This self-supervised training set up is similar to that used to
train autoregressive foundation models. These models will be explored in
detail in §9.1. .

(a) While training we feed in the galaxy rotation curve, and predict the next
observation in its sequence.

(b) While inferring we feed in the full galaxy rotation curve, and extract the
LSTM hidden state as a compressed representation embedding of the curve.
Otherwise, we ignore whatever output (i.e. $\\{p_{1},\ldots,p_{N}\\}$) the
LSTM generates.

Figure 24: A hypothetical surrogate task for extracting rotation curve
representations is shown. $\\{x_{0},\ldots,x_{N}\\}$ is a set of observations
from a galaxy rotation curve, in order of radial distance from the galactic
centre. $\\{p_{1},\ldots,p_{N}\\}$ is the LSTM’s corresponding set of
predictions for $\\{x_{1},\ldots,x_{N}\\}$. $\mathbf{h}$ is the LSTM hidden
state vector. See Fig. 11 for more about the internal workings of the LSTM.

We can generate embeddings via an autoencoding task. Again, let us use an
astronomical example to specify this, and say that we want to extract
embeddings from a set of galaxy observations. We could repurpose a variational
autoencoder for this, training it as normal as described in §6.1. However,
once the model is trained we would discard the decoding part of the network,
and only consider the encoder. To generate embeddings we would then simply
pass in our galaxy images to the trained encoder. The same process can be
carried out by a GAN (§6.2). In the GAN case, we would discard the generator
after training, and use the discriminator’s penultimate layer outputs as our
embeddings.

Supervised networks can also be used to generate embeddings. If a network has
been trained in a supervised manner to classify or regress data, it will have
learnt some properties about that data that helps it to carry out its task. We
can access these learnt representations by taking the outputs from a trained
network’s penultimate layer as an embedding363636 Interestingly, this process
is used in the calculation for the Fréchet Inception Distance (FID)
(ref_heusel2017, 485, 556). The FID acts as a measurement of the visual
similarity between two datasets. The FID works by taking the penultimate layer
representations from a trained Inception-v3 model (ref_szegedy2016, 480) for
each dataset and calculating the distance between them. .

##  8 Astronomy’s third wave of connectionism

Since its astronomical debut in the mid-2010s (ref_regier2015deep, 458)373737
Also compare its companion paper (ref_regier2015celeste, 459). , deep
generative modelling has become a popular subfield within astronomical
connectionism. This popularity is driven by its inherent scalability; the lack
of a need for labelled data allows the methods to be repurposed for any
dataset that might be at hand. Self-supervised connectionism has been around
for longer (i.e. ref_serra1993, 362), but again has recently exploded in
popularity due to its usefulness in wrangling enormous unlabelled datasets.
This section is split into two major parts. We will first outline the history
of deep astronomical generative modelling in §8.1, and the history of
astronomical representation learning will be discussed in §8.2. Although
representation learning is the explicit goal for only the studies described in
§8.2, it must be stressed that representations _can also be extracted from all
the deep generative models described in § 8.1._

###  8.1 Deep astronomical generative modelling

Capturing genuine astronomical data demands accurate knowledge of telescope
behaviour, equipment features, environmental factors during observations, and
data reduction techniques. These complex steps are often tailored to
individual observation sets. However, there’s an alternative to classical
simulation: leveraging examples from a specific survey allows for the
development of a data-driven method to simulate not only the astronomical
signal but also the data’s inherent characteristics. In addition to this, deep
learning models trained to replicate astronomical observations are much
cheaper to run than classical simulation, and so can rapidly generate massive
amounts of data. Data that can then be used for astronomical pipeline
prototyping at scale, aiding the development of new analysis methods, and for
dataset augmentation. Data-driven simulation is made possible via the power of
deep generative models, and this section describes the history of their use
within astronomy.

Seminally, ref_regier2015deep (458) proposed the use of a VAE to model galaxy
observations. They trained their network on downscaled $69\times 69$ crops of
galaxies from a SDSS-sampled dataset containing 43 444 galaxies. They trained
their network in the same way as described in §6.1, and find that the network
is capable of generating galaxies similar to those found in the training set.
They also find that their network produces semantically meaningful embeddings,
noting that their galaxies are clustered by orientation and morphological
type. This same line of enquiry was followed by ref_ravanbakhsh2016 (476), who
showed that VAEs could be used to generate galaxies conditionally.
ref_ravanbakhsh2016 (476) also pioneered the use of GANs to generate galaxy
imagery. ref_spindler2020 (559) used a VAE combined with a Gaussian mixture
model prior (see Eq. 20 and accompanying text) to generate and cluster galaxy
images into morphological types. While the previous studies in this paragraph
used images with relatively small pixel dimensions in their training set,
ref_fussell2019 (518) and ref_holzschuh2022 (612) demonstrated that GANs are
capable of generating large high fidelity galaxy observations. ref_fussell2019
(518) achieved this with a stacked GAN architecture (ref_zhou2016cvpr, 481),
and ref_holzschuh2022 (612) use the related StyleGAN architecture
(ref_karras2018, 501) to the same end. ref_bretonniere2021 (602) use a flow-
based model383838 Flow-based models have not been discussed in detail in this
review, but see ref_weng2018flow (511) for a magisterial introduction to the
subject.  (ref_germain2015, 453, 491) to conditionally simulate galaxy
observations. They found that their approach could produce more accurate
simulations than the previous analytical approach, at the cost of inference
time. Relatedly, ref_smith2022 (628) use a diffusion model to generate large
high fidelity galaxies. They trained their network on two datasets comprised
of galaxies as observed by the Dark Energy Spectroscopic Instrument (DESI;
ref_desi, 517). One, a set of $306\,006$ galaxies catalogued in the SDSS Data
Release 7 (ref_sdss, 386, 414, 422), and the other a set of $1962$ late-type
galaxies, as catalogued in the Photometry and Rotation curve OBservations from
Extragalatic Surveys (PROBES; ref_stone2019, 532) dataset. PROBES contains
well resolved galaxies that exhibit spiral arms, bars, and other features
characteristic of late-type galaxies. They found that their model produces
galaxies that are both qualitatively and statistically indistinguishable from
those in the training set, proving that diffusion models are a competitive
alternative to the more established GAN and VAE models for astronomical
simulation. From all of these studies we can conclude that deep generative
models can internalise a model capable of physically and morphologically
describing galaxies.

Generative models have also been used to simulate astronomical data on larger
scales. In a use-case tangential to galaxy generation, ref_sagan (530) show
that a Spatial-GAN (ref_jetchev2016, 472) can simulate arbitrarily wide field
surveys. They train on the Hubble eXtreme Deep Field, and find that galaxies
‘detected’ within their model’s synthetic deep fields are statistically
indistinguishable from the real thing. Cosmological simulations have also been
explored, with ref_rodriguez2018 (506) using a GAN to generate cosmic web
simulations at pace, and ref_mustafa2019 (526) generating weak lensing
convergence maps at a pace faster than classic simulations. Beyond GANs,
ref_remy2020 (555)393939 This preliminary work has been subsequently extended
in ref_remy2022 (651).  trained a SBGM on simulated maps from MassiveNus
(ref_liu2018, 502), and found that their model was capable of replicating
these maps. They also demonstrated that their model was capable of producing a
likely spread in the posterior predictions. Finally, they demonstrate that a
SBGM is capable of predicting the mass map of the real Hubble Cosmic Evolution
Survey (COSMOS) field (ref_scoville2007, 407).

The image domain translation abilities of GANs in a Pix2Pix-like formulation
(ref_pix2pix, 471, also see Fig. 19(b).) is particularly useful in astronomy.
ref_schawinski2017 (493) demonstrated this use first by training a Pix2Pix-
like model to denoise astronomical data. They trained their network on 4550
galaxies sampled from SDSS. The galaxies were convolved to increase the
seeing, and speckle noise was added. The GAN was tasked with reversing this
process. They found that their method outperformed both blind deconvolution,
and Lucy-Richardson deconvolution. Generative models are also capable of
separating sources, as ref_stark2018 (508) demonstrate by using a Pix2Pix
model to deblend a quasar’s point source emission from the extended light of
its host galaxy. ref_reiman2019 (528) use a similar model to ref_stark2018
(508) to deblend overlapping galaxies.

At the time of writing there are only three examples of score-based (or
diffusion) modelling in the astronomy literature (ref_remy2020, 555, 628,
651)404040 Since the first posting of this review there have been several
workshop papers presented at the 36th Conference on Neural Information
Processing Systems (NeurIPS 2022) on the application of SBGMs to astronomical
problems (e.g. ref_adam2022, 594, 613, 617). Here we will highlight a
particularly neat example of diffusion model application: ref_karchev2022
(613) tackle the inverse problem of strong-lensing source reconstruction and
prove that a denoising diffusion restoration model (DDRM; ref_kawar2022, 614)
inference scheme alongside an off-the-shelf ‘AstroDDPM’ model (ref_smith2022,
628) can restore galaxies that have been through a lensing process.
Remarkably, they achieved this _without any retraining or finetuning of the
original AstroDDPM model_ , demonstrating that generalist pretrained score-
based models like that described in ref_smith2022 (628) can easily be
repurposed for seemingly out-of-distribution downstream tasks. We will revisit
the idea of pretrained models that can be repurposed for downstream tasks when
we discuss ‘foundation’ models in §9. . It is surprising that these studies
are the only examples of score-based modelling in astronomy, as SBGMs produce
generations that rival that of state-of-the-art GAN models, without drawbacks
present in other models (like blurring in the case of VAEs, or mode collapse
and training instability in the case of GANs). SBGMs also have some natural
uses in astronomical data pipelines. For example, an implementation similar to
ref_sasaki2021 (586) could be used for survey-to-survey photometry translation
similarly to ref_buncher2021 (567). The source image separation model
described in ref_jayaram2020 (549) has the obvious application as an
astronomical object deblender (i.e. ref_stark2018, 508, 528, 563). To
summarise, SBGMs are ripe for exploitation by the astronomical community and
we hope to see much interest in this area in the coming years.

###  8.2 Self-supervised astronomical representation learning

In 362, ref_serra1993 (362) proposed using an autoencoder to learn embeddings
for stars as observed by the Two Micron Galactic Survey (ref_calbet1993, 357).
They first proved that their autoencoder model worked better than principle
component analysis (PCA) on the toy problem of separating Gaussian