# Quantifying Non-parametric Structure of High-redshift Galaxies with Deep
Learning

C. Tohill University of Nottingham, School of Physics & Astronomy, Nottingham,
NG7 2RD, UK  L. Ferreira  University of Nottingham, School of Physics &
Astronomy, Nottingham, NG7 2RD, UK  C. J. Conselice  University of Nottingham,
School of Physics & Astronomy, Nottingham, NG7 2RD, UK  Jodrell Bank Centre
for Astrophysics, University of Manchester, Oxford Road, Manchester UK  S. P.
Bamford  University of Nottingham, School of Physics & Astronomy, Nottingham,
NG7 2RD, UK  F. Ferrari  Instituto de Matemática Estatística e Física,
Universidade Federal do Rio Grande, 96203-900 Rio Grande, Brasil

(Received December 16, 2020; Revised March 31, 2021; Accepted May 18, 2021)

###### Abstract

At high redshift, due to both observational limitations and the variety of
galaxy morphologies in the early universe, measuring galaxy structure can be
challenging. Non-parametric measurements such as the CAS system have thus
become an important tool due to both their model-independent nature and their
utility as a straightforward computational process. Recently, convolutional
neural networks (CNNs) have been shown to be adept at image analysis, and are
beginning to supersede traditional measurements of visual morphology and
model-based structural parameters. In this work, we take a further step by
extending CNNs to measure well known non-parametric structural quantities:
concentration ($C$) and asymmetry ($A$). We train CNNs to predict $C$ and $A$
from individual images of $\sim 150,000$ galaxies at $0<z<7$ in the CANDELS
fields, using Bayesian hyperparameter optimisation to select suitable network
architectures. Our resulting networks accurately reproduce measurements
compared with standard algorithms. Furthermore, using simulated images, we
show that our networks are more stable than the standard algorithms at low
signal-to-noise. While both approaches suffer from similar systematic biases
with redshift, these remain small out to $z\sim 7$. Once trained, measurements
with our networks are $>10^{3}$ times faster than previous methods. Our
approach is thus able to reproduce standard measures of non-parametric
morphologies and shows the potential of employing neural networks to provide
superior results in substantially less time. This will be vital for making
best use of the large and complex datasets provided by upcoming galaxy
surveys, such as Euclid and Rubin-LSST.

Astronomy, Galaxies — Machine Learning — Deep Learning — High Redshift

††journal: AJ††software: Astropy (Astropy Collaboration et al., 2018),
Morfometryka (Ferrari et al., 2015), Tensorflow (Abadi et al., 2016), GPyOpt
(The GPyOpt authors, 2016)

##  1 Introduction

A galaxy’s morphology is a useful indicator of its assembly, interaction and
star-formation history. Morphological studies have therefore proven invaluable
for tracing the evolution of the galaxy population over cosmic time. However,
the faintness and small angular size of galaxies at high redshift ($z\gtrsim
1$) makes them difficult to classify in the same manner as those nearby.
Cosmological dimming causes more subtle features, such as spiral arms, to
rapidly disappear with increasing redshift, leaving only the brightest galaxy
components detectable (Barden et al., 2008).

Furthermore, the traditional Hubble sequence is of limited applicability at
high redshift. At early times, higher rates of star-formation and merging
increase the prevalence of more varied and irregular morphologies (Abraham et
al., 1996a; Elmegreen et al., 2005; Conselice & Arnold, 2009; Mortlock et al.,
2013). For studies of distant galaxies we need to consider more general and
robust approaches to characterising galaxy structure.

Galaxy structure can be studied using both parametric and non-parametric
methods. Parametric approaches fit analytic models, such as the Sérsic profile
(Sérsic, 1963), to a galaxy’s light distribution (e.g. Peng et al. 2002;
Buitrago et al. 2008; Simard et al. 2011; Häußler et al. 2013; Robotham et al.
2017). Such parametric methods are valuable for classifying symmetrical
Hubble-type galaxies. However, they break down for more irregular, peculiar-
type galaxies, as they assume a smooth light distribution. Non-parametric
methods make no such assumptions. They are therefore more applicable to the
variety of galaxies seen in the more distant universe, such as those with
‘clumpy’ morphologies (Abraham et al. 1994; Noguchi 1998; Bershady et al.
2000).

Motivated by these considerations, a number of authors, including Abraham et
al. (1994, 1996b), Schade et al. (1995), and Conselice (1997) focused on two
such non-parametric parameters, the concentration (C) and asymmetry (A) of a
galaxy’s light distribution. It has been shown that the concentration
parameter correlates with the bulge-to-disk ratio ($B/D$) of a galaxy, while
the asymmetry parameter is a good indicator of the merger history of the
galaxy (Conselice 2003; Lotz et al. 2008; Nevin et al. 2019). Using these
parameters they were able to separate galaxies into their morphological type
based on their position in this $C-A$ plane. Conselice (2003) expanded on this
by introducing a third parameter, the smoothness (S) of a galaxy’s light
distribution, creating the CAS system which has become one of the most common
non-parametric measures of galaxy structure. This system has since been used
in many investigations of galaxy structure across a wide range of redshifts
(e.g. Yagi et al. 2006; Hoyos et al. 2012). A variety of similar non-
parametric statistics are also in use (e.g. Lotz et al. 2004b; Freeman et al.
2013).

With the imminent arrival of large imaging surveys from new facilities, such
as the Euclid, Rubin and Roman telescopes, it is of paramount importance to
look into the efficacy of existing methods for measuring galaxy structure. For
example, parametric structural measurements are often very time-consuming to
apply to large surveys. Non-parametric measurements are generally faster, but
the algorithms are still typically applied to individual galaxies in series.
While the problem is ‘embarrassingly parallel’, significant computational
resources are required to measure large numbers of galaxies in a timely
fashion. With the future of extragalactic astronomy moving to extremely large
surveys, it is useful to explore more computationally efficient approaches.

One increasingly popular technique, which has already proved useful in a
number of areas of astronomy (Frontera-Pons et al. 2017; D’Isanto & Polsterer
2018; Pearson et al. 2019), is machine learning. In particular, deep learning,
utilizing neural networks, can apply sophisticated analyses to large datasets
at a much faster rate than conventional methods (e.g. Tuccillo et al. 2018).
Deep learning has been applied to the morphological classification of both
nearby (Dieleman et al. 2015; Cheng et al. 2020) and distant galaxies
(Huertas-Company et al., 2015; Ferreira et al., 2020). It has also been shown
to be very effective at reproducing parametric structural measurements
(Tuccillo et al., 2018). However, as yet, deep learning has not been applied
to determine the non-parametric CAS parameters. Given the arguments above,
this could be a highly valuable tool for studying the local and high-redshift
galaxy population in the next generation of surveys.

In this work we therefore create neural networks capable of predicting
concentration and asymmetry parameters from a galaxy’s image. (For now we
neglect the smoothness parameter as it is more difficult to measure at high
redshifts and needs a separate treatment.) We show that our networks are
consistent with conventional algorithms in their output, and demonstrate
reliable behaviour down to very low signal-to-noise ratios. Furthermore, we
find that our trained network is able to analyse $\sim 10,000$ galaxies in
under 1.5 seconds, much faster than convention methods, making it well-suited
to the large number of galaxies in future surveys.

This paper is organised as follows. In §2 we introduce the imaging data used
in this work and describe how the conventional CAS parameters are measured
using the Morfometryka software (Ferrari et al., 2015). The pre-processing of
the data and all data augmentation is detailed in §3.1 and 3.2. In §3.3 and
3.4 we describe the architecture and optimization of our neural networks. The
resulting performance of these networks is demonstrated through a number of
tests in §4, concluding with a brief summary in §5.

##  2 Data

###  2.1 CANDELS Fields

All of the images used in this project were taken with the Wide Field Camera 3
(WFC3) of the Hubble Space Telescope (HST) as part of the Cosmic Assembly
Near-infrared Deep Extragalactic Legacy Survey (CANDELS). We use data from all
5 CANDELS fields: the Great Observatories Origins Deep Survey (GOODS)-North
and GOODS-South fields, COSMOS, Extended Groth Strip (EGS) and Ultra-Deep
Survey (UDS).

The CANDELS/Deep survey ($5\sigma$ point-source limit $H=27.7$ mag) covers an
area of $\sim 125~{}\textrm{arcmin}^{2}$ with a resolution of $0.06\arcsec$
per pixel (Grogin et al., 2011; Koekemoer et al., 2011). In total we have
$\sim 150,000$ galaxy postage-stamp images. These galaxies have photometric
redshifts covering $z=0$–$7$, with many parameters already calculated,
including star formation rates (SFR) (Duncan et al., 2019) and CAS values. The
apparent magnitude–size distribution of our sample is shown in Fig. 1. In this
paper we use imaging from the $H$-band ($F160$W), as it provides the most
complete deep-coverage over all five CANDELS fields.

Figure 1: Distribution of the apparent magnitude in the H-band (F160W) vs half
light radius, with a resolution of $0.06\arcsec$ per pixel, for the CANDELS
galaxy sample used in this work.

###  2.2 Concentration and asymmetry

As mentioned in the introduction, non-parametric methods have been used for
many years to analyse the light distributions of distant galaxies, in order to
better understand their structure (Conselice 2003; Lotz et al. 2004a; Sazonova
et al. 2020). Such methods make very few assumptions, and so can be applied to
peculiar and irregular galaxies as well as to classic Hubble types.

In this paper we utilise a subset of the CAS (Concentration, Asymmetry and
Clumpiness) system as defined in Conselice (2003). This is a robust, non-
parametric method for classifying galaxy structure, in a manner that is
sensitive to their ongoing and past formation modes. In this paper, only
concentration and asymmetry are considered. The concentration ($C$) is based
on the measurement first established by Bershady et al. (2000), which was
found to correlate with both galaxy bulge-to-disk ratio (B/D) and the
effective radius of the bulge. This quantity is defined as

$C=5\log_{10}\left(\frac{r_{80}}{r_{20}}\right).$ (1)

where $r_{20}$ and $r_{80}$ are the radii containing 20% and 80% of the total
light of the galaxy, respectively. The value of $C$ is simply a measure of how
concentrated the light in the central region is relative to the galaxy’s
overall size. Galaxies with higher concentrations are typically ellipticals,
early-type disks, and edge-on disks. In this manner, it shares similarities
with the Sérsic index (Graham et al., 2005).

Figure 2: Distributions of asymmetry and concentration values for our selected
sample of CANDELS galaxies, measured with Morfometryka. The solid line shows
the selected training sample and the dotted line shows the training sample
after our augmentation process.

Galaxy asymmetry was first used in a basic form by Schade et al. (1995), when
trying to classify distant galaxies imaged with HST. Asymmetry ($A$) is
determined by rotating a galaxy $180^{\circ}$ about its center and then
subtracting from the original image. The centre of rotation is determined by
an iterative process that finds the minimum asymmetry. Further algorithmic
details are described in Conselice (2000) and Conselice (2003). The absolute
values of the residuals are summed and normalized by the original galaxy flux.
The resulting asymmetry contains a contribution from the background noise.
This is accounted for by subtracting a background term, determined by
computing the asymmetry for small areas of sky near the galaxy. The basic
calculation for the asymmetry is therefore given by

$A=\frac{\sum|I-I_{180}|}{I}-A_{\mathrm{bkg}}\,,$ (2)

where $I$ is the original galaxy image, $I_{180}$ is the rotated galaxy image,
and $A_{\mathrm{bkg}}$ is the background asymmetry (discussed further below).

Asymmetry can be used to identify a number of interesting galaxy classes, such
as mergers and starburst galaxies (Conselice 1997, Conselice 2000, Bluck et
al. 2012). These types of galaxies have a higher $A$ value than regular
ellipticals and disk galaxies, due to distributed areas of increased star
formation.

###  2.3 Morfometryka

CAS measurements were originally obtained using IRAF. However, a more modern
implementation, in Python, is provided by Morfometryka111The results in this
paper are based on Morfometryka version 8.2 (Ferrari et al., 2015).
Morfometryka extracts a number of features from astronomical images, such as
non-parametric morphology (including the CAS parameters) and Sérsic profiles.
Full details of the software can be found in Ferrari et al. (2015), however we
will briefly describe how the parameters used in this paper were calculated.

Figure 3: _Top row:_ Images of the galaxies with high concentration values.
These galaxies appear to be compact, spheroidal and have no close neighbours.
These galaxies also typically have low asymmetry values, although there are
some objects with high concentration and asymmetry. _Bottom row:_ Images of
galaxies with high asymmetry values. Many of these galaxies appear to be
undergoing mergers and have tidal effects present, although there are
occasional cases of line-of-sight projection. The concentration ($C$) and
asymmetry ($A$) are indicated above each galaxy stamp.

Morfometryka calculates the concentration, $C$, as explained in 2.2, with the
exception that the factor of 5 in Eq. 1 is omitted. However, in order to
remain consistent with previous studies, this factor was re-applied to our
concentration values.

The asymmetry, $A$, is also determined as described in Section 2.2, by
applying Eq. 2 within a $1.5$ $\times$ Petrosian-radius elliptical aperture
centred on the galaxy. However, the background term in Eq. 2 is computed in a
way that slightly deviates from the original CAS implementation. The standard
approach utilises a single background region. Originally Morfometryka did not
include the background asymmetry correction term. For our measurements, we
construct a $10\times 10$ pixel grid over the image area outside the galaxy
segmentation map. We then measure the asymmetry for each cell in the grid,
according to the first term of Eq. 2. Finally, we select the median asymmetry
across all the cells as our background term, $A_{\mathrm{bkg}}$. This ensures
a robust and accurate background correction, improving upon the original
background subtraction by eliminating the bias inherent in choosing only one
background area. This is now incorporated into Morfometryka.

The errors on the concentration values are derived from those of the
individual size measurements, which assume Poisson distributed fluxes. The
typical error on $C$ is $\pm 0.23$. The error on the asymmetry values were
calculated using the method described in Conselice (2003). We find that the
typical error on $A$ is $\pm 0.072$ for our sample.

We applied Morfometryka to all of the images in our dataset. We then select
suitable galaxies for our analysis based upon the the steps described in §3.1.
The distributions of concentration and asymmetry values for our selected
sample are shown in Fig. 2.

A subset of galaxy images were inspected to check that the measurements
correspond to visual expectations. As can be seen from the top row of Fig. 3,
galaxies with high $C$ values appear compact and spheroidal. Such galaxies
typically have low $A$ values, reflecting a broad anti-correlation between $C$
and $A$ for normal Hubble types. Galaxies with high asymmetries are shown in
the bottom row of Fig. 3. The contrast between the two sets of galaxies is
clear, with high $A$ galaxies appearing disrupted, or possessing features
associated with merging, such as tidal tails and multiple bright sources. Note
that high asymmetry galaxies span a range of concentrations. This is a
reassuring reconfirmation of how these parameters have been seen to behave in
past studies (Conselice et al., 2008, 2011).

##  3 Method

###  3.1 Pre-Processing

The initial images used in this analysis are $101\times 101$ pixel cutouts,
with the target galaxy in the center of each stamp. As we are only interested
in training the network to predict the $A$ and $C$ values for the target
galaxy, we need to remove any other sources. In order to remove neighbouring
sources from the cutouts, the galclean algorithm (Ferreira et al., 2018) was
utilised. This algorithm removes any non-central sources at a certain
threshold above the background level. These masked areas are replaced with
values sampled randomly from the background distribution to ensure they do not
leave shapes which could be picked up by the network.

The majority of our galaxies have a half light radius of $\sim 10$ pixels. For
computational efficiency, the individual galaxy images are therefore further
reduced in size to $60\times 60$ pixels, centered on the galaxy.

Since we are interested in measuring structure irrespective of overall galaxy
brightness, we individually normalize our images. The pixel values of each
image are rescaled so that the maximum pixel value of each image has a value
of 1. This is also a standard pre-processing procedure for deep learning. It
improves learning efficiency by ensuring that the inputs to the networks are
compatible with the domain of the activation functions used within the model.

We wish to consider only reliable galaxy detections, for which structural
parameters can reasonably be obtained. We therefore limit our sample to
galaxies above a minimum signal-to-noise. We define the average signal-to-
noise per pixel for each galaxy as

$\mathit{SNR}_{p}=\frac{L_{\mathrm{tot}}}{q\pi R_{p}^{2}\sigma_{sky}}$ (3)

where $L_{\mathrm{tot}}$ is the total integrated flux within the Petrosian
region (with semi-major axis $R_{p}$), $q$ is the axis ratio measured from the
intensity distribution using the image moments, and $\sigma_{sky}$ is the
standard deviation of the sky background. By visual inspection we define a
selection for our galaxy sample of $\mathit{SNR}_{p}>2$. We also limit our
sample to galaxies with $R_{p}>5$ to ensure they are properly resolved.

Once these steps have been completed, we are left with 94,192 galaxy images
with a median $\mathit{SNR}_{p}\sim 4.5$. These images were then split
randomly into training (80%), testing (10%) and validation (10%) datasets to
apply to our machine learning methods. With over 9,000 galaxies in each of our
testing and validation sets, our performance estimates will be both accurate
and precise.

###  3.2 Data Augmentation

Unbalanced datasets, whereby there are many more galaxies at one particular
value compared with others, can cause issues when dealing with both regression
and classification problems in machine learning. The relative frequency of
classes in the training set acts as a prior; the network may therefore be
biased against identifying rare cases. In extreme circumstances, the network
may fail to learn to identify rare cases at all. One way to combat this issue
is by data augmentation (Shorten & Khoshgoftaar, 2019).

Data augmentation is primarily used as a way of creating a larger training
sample, which more finely samples the space of possible inputs. It is a form
of regularisation and hence helps to prevent overfitting. By selectively
expanding the size of the potential training set, augmentation can also help
to balance the prevalence of different classes, while still using all of the
input data.

Looking at Fig.2, there is a large imbalance in the CAS values for our sample,
such that very high asymmetries are not common, nor are very low or high
concentrations. As we want our model to be accurate across all concentrations
and asymmetries, we selectively apply augmentation to create a more balanced
training set. That is, we need to supplement the images that occupy the
parameter space where there are few galaxies. For the range of $C$ or $A$
values where there are around half the number of images compared to the median
value, we rotated each image by $90^{\circ}$ once. Where there are relatively
fewer images, we apply a greater variety of augmentations: rotating by
$90^{\circ}$ 3 times and mirroring along both axes. These images were then
shuffled and added to the training set. After data augmentation, our training
sample increases in size from 75,353 to 141,453 images.

###  3.3 Convolutional Neural Networks

Figure 4: Architecture of the asymmetry network after optimization of its
hyperparameters. This network takes input images of size $60\times 60$ pixels,
followed by 3 blocks, each containing 2 convolutional layers with 32, 64 and
128 features, respectively. Each block is followed by an average pooling layer
of size 2 and a batch normalisation layer. Finally there are 2 fully-connected
layers of size 1024 each, which is followed by dropout of 0.5 before the
output value, i.e the network’s prediction of the asymmetry.

The purpose of this project is to efficiently and robustly predict CAS values
of a galaxy from an image. We chose to implement a Convolutional Neural
Network (CNN), as these are known to perform well when dealing with spatial
structured data. CNNs are made up of convolutional layers, which are able to
extract features from images by applying multiple filters (convolutional
kernels) to the image. Individually, these filters can detect simple features.
However, successive layers act hierarchically, identifying increasingly
complex patterns. One major advantage of CNNs for image classification
problems is the fact that they are able to exploit the spatial structure of
the data which in turn reduces the number of parameters and allows the
recognition of location invariant features.

CNNs were first popularised for image recognition/classification problems with
the creation of LeNet-5 (Lecun et al., 1998), a network trained to classify
handwritten digits. From this, CNNs have been applied in a range of fields,
addressing a number of different problems and are becoming increasingly
popular in astronomy.

CNNs were first utilised for galaxy classification by Dieleman et al. (2015)
using data from the Galaxy Zoo project (Willett et al., 2013). While many
others had applied different machine learning (ML) techniques to address this
problem (e.g., Storrie-Lombardi et al. 1992; Naim et al. 1995; Huertas-Company
et al. 2008; Banerji et al. 2010), these all required an earlier step of
extracting features (often including CAS parameters or similar) from the
images. The advent of CNNs provided a technique for efficiently extracting
high-quality information directly from images. CNNs have since seen wide usage
in extra-galactic astronomy, including morphological classification (e.g.,
Domínguez Sánchez et al. 2018; Cheng et al. 2020; Barchi et al. 2020),
performing photometry (Tuccillo et al. 2018; Boucaud et al. 2020), and
estimating merger rates (Ferreira et al. 2020).

There are many factors to consider when choosing the optimum architecture for
a network. Many early studies based their architecture on previous studies
(Huertas-Company et al. 2015; Domínguez Sánchez et al. 2018; Aniyan & Thorat
2017), trial-and-error (Dieleman et al. 2015; Feinstein et al. 2020), and
arbitrary choices. However, there are a number of optimisation techniques that
allow these choices to be optimised in a more satisfactory manner for the
problem at hand. The variety of network architectures we consider, and our
method for selecting from these, are described in the following section.

To evaluate how well our networks are performing we compute the mean absolute
error (MAE) and root mean squared error (RMSE) of the network’s predictions.
The RMSE metric also serves as our loss function. The MAE is simply a measure
of the average magnitude of error between the network’s prediction and the
expected result,

$\mathrm{MAE}=\frac{1}{n}\sum_{i=1}^{n}|y_{i}-\hat{y}_{i}|\,,$ (4)

where $n$ is the number of samples, $y_{i}$ is the expected value and
$\hat{y}_{i}$ is the network’s prediction. The RMSE is similar to the MAE, but
it is more sensitive to large errors and so can indicate if there are many
outliers present. It is calculated as

$\mathrm{RMSE}=\sqrt{\frac{1}{n}\sum_{i=1}^{n}{(y_{i}-\hat{y}_{i})}^{2}}\,.$
(5)

###  3.4 Bayesian Optimization

The various choices that must be made before training a network can be
considered as hyperparameters. These include aspects of the network
architecture, such as the number of convolutional layers and number of filters
in each layer, and of the training, such as the update algorithm, learning
rate and batch size. Varying these choices can significantly alter the
performance of the trained network. The problem of determining which
combination of hyperparameters will be best suited to a given problem
typically involves a trial and error process, which is often only partially
explored, or entirely neglected, resulting in a non-optimal solution.

To avoid this, many optimization techniques have been developed, from
simplistic random or grid-based searches (Bergstra et al., 2011), to more
advanced techniques such as random forests (Hutter et al., 2011). The aim of
these techniques is to find the optimum hyperparameters that will minimise the
average loss. Traditionally, these techniques can be computationally
expensive, as each variation in the hyperparameters results in a new version
of the network which must be trained and then evaluated. Bayesian Optimisation
(Snoek et al., 2015) provides a more efficient solution: a record of past
evaluation results are kept and used to form a probabilistic model, which the
method builds upon, reducing the time to converge on a optimal model.

Our networks comprise a number of convolutional blocks, between 1 and 3, with
each block having either 1 or 2 convolutional layers. Each convolutional layer
in a block has the same number of filters between 8 and 256 in powers of 2.
The kernel sizes are all fixed to $3\times 3$. Each convolutional block is
followed by a BatchNormalisation layer and an AveragePooling layer of fixed
(2$\times$ 2) size. Originally we started with MaxPooling layers, however we
found that the networks’ performance improved when using AveragePooling
layers. (Similar behaviour was found by Pasquet et al. 2019 when analysing
SDSS images.) Following the convolutional blocks we add some fully connected
layers, with their number and size as hyperparameters. The number of fully
connected layers ranges from 1 to 4, with each layer having the same number of
filters between 128 and 1024 in powers of 2. We include a dropout layer before
our output layer as a form of regularisation, allowing the dropout rate to
vary as another hyperparameter. The dropout rate is allowed to vary
continuously between 0.25 and 0.60. The activation function is fixed to the
common ReLu (Nair & Hinton, 2010) non-linearity.

When training a network, an optimization algorithm adjusts the weights to
minimise the cost function. With a plethora of optimizers now available, we
have included the choice as a hyperparameter, selecting from a pool of those
most commonly used, we include Adam, Adadelta, RMSprop, SGD and Adamax. We
also set the learning rate as a hyperparameter, where we evaluate 5 values,
0.001, 0.005, 0.01, 0.05 and 0.1.

The parameters we defined as hyperparameters and their optimised values are
displayed in Table (1).

Each network was trained for a maximum of 300 epochs, but we applied ”early
stopping” to halt the training when the validation loss had converged, which
was typically after $\sim 100$ epochs.

Our Bayesian Optimization was carried out using the GPyOpt python package (The
GPyOpt authors, 2016), with the aim to minimise the RMSE of the networks. Each
network created during the optimization was trained and validated using the
samples defined in §2. The MAE, RMSE and the Pearson coefficient were
monitored for each iteration in the optimization. The network that had the
lowest MAE and RMSE was selected as the optimum architecture for our network.

Hyperparameter | Optimum value
---|---
| Asymmetry | Concentration
batch size | 512 | 512
convolutional blocks | 3 | 3
conv. layers per block | 2 | 2
fully-connected layers | 2 | 2
fully-connected layer size | 1024 | 512
number of filters | 32 | 64
optimization | Adamax | Adam
learning rate | 0.001 | 0.001
dropout | 0.50 | 0.55
Table 1: Summary of the hyperparameters selected by the Bayesian Optimization
technique.

The architecture of the CNN selected for our asymmetry network is shown in
Fig.4. To ensure that the choice of optimum architecture is robust, we retrain
multiple times, and compare the variation in the loss to the variation
observed between different networks. The variation in MAE for the optimum
asymmetry network is quite stable and varies by $\sim 0.001$. Comparing the
top 10 network architectures, we find that the MAE varies by $\sim 0.002$. The
hyperparameters of these networks are quite similar, although the number of
fully connected layers varies between 1 and 2, the dropout rate between 0.47
and 0.56, and the optimizer varies between Adamax and Adadelta. These
parameters are not as significant in determining the optimum network.

The selected concentration network has a similar architecture, with some
slight variations. The MAE loss variation across different training runs is
$\sim 0.002$, and the MAE of the top 10 architectures vary by $\sim 0.003$,
very similar to above. Looking at the variation in the architectures which
give equivalent performance, we see that the number of fully connected layers
varies between 1 and 2 layers, the batch size between 256 and 512, the dropout
rate from 0.3 to 0.6, and the number of convolutional blocks varies between 2
and 3.

Following our use of Bayesian Optimisation and the above tests, we can be
confident that our final selected networks are well-optimised. However, it is
also reassuring that the performances we report below are robust to minor
variations in network architecture and training.

##  4 Results

Figure 5: Our final network’s predictions of concentration values for the test
set versus those measured by Morfometryka. The network does not see any of the
images used in this comparison during its training. The MAE of the network was
0.15, with a RMSE of 0.21 (see §4). The Pearson correlation coefficient of
0.75 indicates that there is a strong correlation between the two.  Figure 6:
Our final network’s predictions of the asymmetries for the test set versus
those measured by Morfometryka. The network did not use any of these images
during its training. The MAE of the network was 0.045, with a RMSE of 0.065
(see §4). The Pearson correlation coefficient of 0.76, indicates that there is
a strong correlation between the two.

As explained above, we train our networks on a subset (80%) of the images and
select our optimal model by its performance on a validation set (10%). To then
evaluate our selected, trained network, we use an additional independent test
set (10%; 9,420 galaxies). This ensures that the metric used to evaluate the
network’s performance is not biased by over-fitting the hyperparameters. We
find that the networks perform similarly on both the test and validation sets:
another indication that the selected network architecture is robust.

###  4.1 Model performance

The basic results of the concentration network after the hyperparameter
optimisation can be seen in Fig.5. The network’s predictions correlate
strongly with the Morfometryka measurements, with a MAE value of $0.15$ with a
RMSE of $0.21$. This error is lower than the average error on the
concentration measurements. This shows that our machine learning regression
can measure these parameters just as well as the direct measurement method.
Hence, the values from the network can be utilised with a similar level of
confidence as the original algorithm However, the scatter does get larger at
parameter values where there are fewer galaxies.

The results for the asymmetry network, again after the hyperparameter
optimisation, can be seen in Fig.6. The network’s predictions for the
asymmetries have a MAE of 0.045 with a RMSE of 0.065. As before, this error is
lower than the average error on the asymmetry measurements, showing that our
networks can be used to reliably measure both the concentration and asymmetry
values for a galaxy.

Overall, both networks perform well, achieving low residuals between the
measured values and the network’s predictions. Looking at the images of
galaxies where there was a large difference between our networks and the
Morfometryka-measured CAS values, we find that they are quite noisy, with
$\mathit{SNR}_{p}<4$. From this we decided to further investigate the impact
of noise on both our network predictions and the directly-measured CAS values.

###  4.2 Impact of noise

Our networks’ ability to accurately predict CAS values is potentially
dependent on the noise level in a given galaxy image. To investigate this, we
consider how the residuals (network prediction $-$ Morfometryka value) of each
network depend on the signal-to-noise per pixel, $\mathit{SNR}_{p}$. This is
shown in Fig.7. Confirming the results from Figs.5 & 6, we see that for
galaxies with moderate and high $\mathit{SNR}_{p}$ the residuals are close to
zero. The random scatter is also fairly constant with $\mathit{SNR}_{p}$,
indicating that our networks are reliable across a broad range of
$\mathit{SNR}_{p}$. There are a small number of galaxies with large deviations
at the higher $\mathit{SNR}_{p}$ however, when inspecting these images we find
that most contain another source in the image that was not removed by the
galclean algorithm. This could explain why the measurements for these galaxies
from Morfometryka and our networks varied. Within the low $\mathit{SNR}_{p}$
regime we find a slight bias where the networks, on average, under-predict the
values measured by standard algorithms.

Figure 7: _Top:_ The residuals of the asymmetry network and the measured
values from Morfometryka versus the $\mathit{SNR}_{p}$ of each galaxy image.
_Bottom:_ The fractional residuals between the concentration network and the
measured values from Morfometryka versus the $\mathit{SNR}_{p}$ of each galaxy
image. In both panels, black points represent the means plotted for each bin
with the darker shaded region representing $\pm 1$ standard deviation and the
lighter shaded region shows the area containing $2\sigma$. The trend indicates
that there is a bias at low $\mathit{SNR}_{p}$, where the networks will, on
average, slightly under-predict the values measured by the standard
algorithms.

The origin of this systematic trend at low signal-to-noise is interesting. Our
networks have been trained to reproduce the measured values, and are clearly
doing so in the majority of cases. So why the deviation at low
$\mathit{SNR}_{p}$? This could be seen as a failure of our model to capture
the details of the measurements. On the other hand, we apply regularisation
and optimise the hyperparameters to avoid over-fitting, with the aim of
producing a generally applicable model, capable of accurate measurements for a
wide variety of images. One optimistic possibility is that our networks are
able to learn a model which is more robust than the regular methods. This is
not inconceivable, since the regular methods must make a series of algorithmic
‘decisions’ (masking, fitting elliptical isophotes, recentering, etc.). The
networks, instead, consider all of these issues within a single ‘holistic’
calculation.

In order to determine if this low signal-to-noise trend is a bias in our
networks or in the standard algorithm (as implemented in Morfometryka), we
investigate how noise impacts these two approaches in an independent manner.
For this test, we select a sub-sample of 622 high $\mathit{SNR}_{p}$ galaxies
($\mathit{SNR}_{p}>10$), with low asymmetry residuals ($|\Delta A|<0.01$) from
the validation set. These galaxies also have low residuals in their
concentration values. These galaxies are chosen as both the network and
Morfometryka predicted these galaxies to have similar parameters, and hence we
can assume these to be the true values for the purpose of this test.

We then produce versions of each galaxy image with varying $\mathit{SNR}_{p}$
values. To do so, we first measure the mean and standard deviation in the
background of the original galaxy image, then create an image with
corresponding Gaussian noise. To this simulated background image, we add the
original image with the overall flux scaled, such that we achieve our desired
$\mathit{SNR}_{p}$. Finally, the image is normalized in the usual manner
(§3.1). Except for the variation in $\mathit{SNR}_{p}$, each galaxy image
remains identical to its original version. An example of these simulated noisy
images can be seen in Fig.8.

Figure 8: An example of our simulated noisy galaxy images. The top-left panel
shows the original image, while the remainder show the same galaxy at
different simulated $\mathit{SNR}_{p}$.

The asymmetry and concentration values for these galaxies are then re-measured
at each $\mathit{SNR}_{p}$, using both Morfometryka and our trained networks.
The variation from the values measured in the $\mathit{SNR}_{p}=10$ image is
plotted as a function of decreasing $\mathit{SNR}_{p}$ in Figs.9 & 10.

For asymmetry, it can be seen that at both high and moderate
$\mathit{SNR}_{p}$ the values recovered by our network are very similar to the
‘true’ values. For $\mathit{SNR}_{p}>5$ the recovered values vary with an
average standard deviation of $0.025$, reflecting the uncertainties due to
shot noise.

Furthermore, this scatter is significantly lower for our network than the
standard algorithm. The average scatter in the asymmetry measurements, at
$\mathit{SNR}_{p}\geq 5$, is $0.025$ compared to $0.037$ for Morfometryka.
Since we have already seen that our networks accurately recover Morfometryka
measurements, this suggests that our network is using information in these
moderately-noisy images that is not utilized by the Morfometryka algorithm.

At low $\mathit{SNR}_{p}$, we find a bias present in both the network and
Morfometryka, such that the $A$ values are, on average, overestimated.
However, it can be seen from Fig.9 that Morfometryka has a larger bias at
these low $\mathit{SNR}_{p}$, indicating that our network is slightly more
accurate in the low $\mathit{SNR}_{p}$ regime. At $\mathit{SNR}_{p}=3$ the
bias in the network’s asymmetry values is $0.016$ compared to $0.028$ for
Morfometryka.

While there is scatter in the individual measurements, on average our network
is able to accurately estimate asymmetry, with little bias from the ‘true’
value, at a lower $\mathit{SNR}_{p}$ than the original algorithm. This is
useful for merger fraction estimates, especially at high redshift, as we can
now include galaxy images down to a $\mathit{SNR}_{p}$ as low as 3 while still
retrieving unbiased measurements with our network. This means that we are able
to measure reliable CAS parameters for more galaxies using deep learning and
then we can with a direct measurement.

For concentration, plotted in Fig.10, we again see a a difference in the
variation of the measurements for moderately-noisy images, with our network
producing a significantly lower scatter than the standard algorithm. The
average scatter in the concentration measurements at $\mathit{SNR}_{p}\geq 5$
are $0.04$ for the network compared to $0.07$ for Morfometryka. Both the
network and Morfometryka slightly overestimate the concentration at
$\mathit{SNR}_{p}\leq 5$.

We also investigate the ‘catastrophic’ fraction ($f_{c}$) of both Morfometryka
and our network, that is the number of galaxies that fall outside of 2 sigma
deviation from the mean. Again the network performs marginally better, with
$f_{c}$ at $\mathit{SNR}_{p}\geq 5$ being 4.3% for Morfometryka compared to
3.3% for our network.

Based on these results, we conclude that our deep-learning approach is
performing at least as well as traditional measurements of non-parametric
structure.

Figure 9: _Top:_ Deviation in the asymmetry measurements from
$\mathit{SNR}_{p}=10$ for Morfometryka as the signal-to-noise is decreased in
bins of $\mathit{SNR}_{p}$. At $\mathit{SNR}_{p}\leq 5$ there is a consistent
bias, such that, on average, Morfometryka overestimates the values. Black
points indicate the mean residual at each bin, with error-bars showing $\pm 1$
standard deviation. The median deviations are shown by open points. _Bottom:_
As the top panel, but for our trained network. A similar trend can be seen at
low $\mathit{SNR}_{p}$, but the scatter and systematic deviation is lower.
This shows that the network is more stable and accurate than Morfometryka at
these low $\mathit{SNR}_{p}$ values.

Figure 10: _Top:_ Deviation in the concentration measurements from
$\mathit{SNR}_{p}=10$ for Morfometryka as the signal-to-noise is decreased in
bins of $\mathit{SNR}_{p}$. Black points indicate the mean residual at each
bin, with error-bars showing $\pm 1$ standard deviation. The median deviations
are shown by open points. _Bottom:_ As the top panel, but for our trained
network. At $\mathit{SNR}_{p}\leq 5$ there is a slight bias, such that, on
average, Morfometryka and the network overestimate the values of the sample.
It can be seen that the $C$ measurements are quite stable for
$\mathit{SNR}_{p}\geq 3$. It can also be seen that the network produces a
significantly lower scatter than Morfometryka

.

Figure 11: An example of our redshifted Frei galaxy images. In our redshifting
technique we include effects due to cosmological dimming, luminosity
evolution, size evolution and geometric scaling. Top left panel shows the
reference image used in our tests. We start at z = 0.2 so that the whole
galaxy fits within the $60\times 60$ crop which is the input for our networks.
The simulated redshift is indicated in each image.

Figure 12: Variations in the $C$ and $A$ measurements from our trained
networks and Morfometryka as a function of simulated redshift. Our simulations
include geometric scaling (with accompanying reduction in resolution), size
evolution, cosmological dimming and luminosity evolution. The original
asymmetry measurements range from $0.01<A<0.48$, while the concentration
measurements vary from $2.4<C<5.0$. There are a total of 100 galaxies from the
Frei sample in this analysis. The variations for a sample of individual
galaxies in the sample are shown by faint lines. These span a range of
morphological types and initial asymmetry and concentration values. Black
points indicate the mean variation and the error bars show the standard
deviation of the individual measurements.

###  4.3 Impact of redshift effects

While the previous test examined how our networks fare with respect to noise
alone, here we combine the effects of signal-to-noise and resolution to
determine our networks’ performance for galaxies at high redshifts.

It is known that at higher redshifts, cosmological dimming and decreasing
apparent size result in galaxies appearing more symmetric and less
concentrated than they otherwise would (Conselice, 2003). However, through the
use of simulations, we can model and correct for this variation. We thus
quantify the extent to which $C$ and $A$ values estimated by our networks are
biased by these issues, and hence the level of any correction which should be
applied when comparing galaxies at different redshifts.

To investigate how the performance of our networks varies with redshift, we
take a sample of nearby galaxies, with reliably measured concentration and
asymmetry values, and simulate how they would appear at higher redshift. For
this test, we selected objects from the Frei catalogue of nearby galaxies
(Frei et al., 1996). These are a well studied sample of regular, nearby
galaxies, containing all Hubble types and with previously measured $A$ and $C$
values (Conselice, 2003). We simulate the appearance of these galaxies as if
they were observed at a range of redshifts, from $z=0.1$–$7$.

There are a number of effects that need to be considered when artificially
redshifting galaxies. The first effect we address is geometric scaling,
whereby the apparent size of the galaxy will decrease when viewed at a higher
redshift. We follow the same procedure as described in Conselice (2003) and de
Albernaz Ferreira & Ferrari (2018) to reduce the sizes of the galaxies to how
they would appear at higher redshifts. Previous simulation work has kept the
physical size of the galaxies constant. However, it is well known that
galaxies of a given stellar mass are intrinsically smaller at higher
redshifts, reducing in size by a factor of $\sim 5$ between $z=0$ and $z=3$
(Trujillo et al., 2007; Buitrago et al., 2008). We therefore introduce size
evolution to better represent the properties of high-redshift galaxies.

We use the size evolution determined by Whitney et al. (2019), which is based
on the same $H$-band data from CANDELS GOODS North and South fields that we
use in our training sample. They measured how the average physical Petrosian
radius, $R_{p}$, changes with redshift for a mass selected sample, finding
that it varies according to $R_{p}(z)=\frac{R_{p}(z=0)}{(1+z)^{\beta}}$, with
$\beta=-0.97$. We therefore multiply the geometric scaling factor by this
value to correct for the size evolution in our simulation.

After the (flux-preserving) geometrical scaling, we apply cosmological dimming
Tolman (1930), according to

$I(z)=\frac{I(z=0)}{(1+z)^{4}}$ (6)

where $I(z)$ is the observed intensity. This is one of the major issues when
detecting high redshift galaxies, as it introduces a bias such that only the
brightest, most compact galaxies are detectable. The intrinsic brightness of
galaxies, with a given stellar mass, varies with redshift. We therefore
implement an evolution in the surface brightness of the galaxies as outlined
in Whitney et al. (2020). They found that the correction for the intrinsic
surface brightness follows

$\mu(z)=\mu(z=0)(1+z)^{\alpha}$ (7)

where $\mu=-2.5log(I)$ and $\alpha=-0.13$. The value of $\alpha$ was found to
vary from $-0.09$ to $-0.18$, but this will not result in much variation in
our results. The value of -0.13 was the value found for their size corrected
sample.

To complete our simulations, the galaxies are convolved with the HST PSF in
the H160-band filter and placed in an actual CANDELS background.

We do not account for morphological or magnitude k-corrections. Instead, we
test how a galaxy image would vary in restframe optical wavelengths, i.e.
choosing appropriate observed filters for different redshifts. While we are
currently only able to probe restframe optical up to $z=3$ with HST, future
surveys, such as JWST, will be able to probe up to $z=7$. Furthermore, it has
been found that the CAS parameters do not vary much between the UV and optical
for star forming galaxies (Conselice, 2003).

We select the brightest galaxies for this test, such that they are above
$\mathit{SNR}_{p}>2$ in all images out to $z=7$, as this was the cut off used
in our training sample. This leaves us with 100 out of the original 112
galaxies in the Frei sample. The original asymmetry values of the sample range
from $0.01<A<0.48$, while the concentration measurements vary from
$2.4<C<5.0$. The $A$ and $C$ parameters of each galaxy were remeasured at each
redshift. We only consider $z=0.2$ onward, since we need the whole galaxy to
fit within a $60\times 60$ pixel image, for input to our networks. An example
of one of the redshifted galaxies is shown in Fig.11.

The variations of the CAS values measured by our networks and Morfometryka are
plotted against redshift in Fig.12, with the full distributions at a sample of
redshifts shown in Fig.13 for comparison. As expected, at higher redshifts
both methods measure the galaxies to be more symmetric and less concentrated
than at $z=0.2$. As the outer regions of a galaxy fade below the background
noise, and their apparent size approaches the resolution limit of the PSF,
they appear more symmetric, as has been found previously (Conselice, 2003).

While the Frei sample of galaxies appear somewhat different to those the
networks were trained on, we see that the networks still perform well,
measuring values similar to Morfometryka. The average change in asymmetry at
$z=1.0$ compared to $z=0.2$ is $0.074$, which is similar to the average error
on the high-$z$ asymmetry measurements. This is important for merger
estimates, as this variation is small enough to avoid a merger appearing as a
non-merger and vice versa. At redshifts higher than $z=1$ the average
variation is around twice the average error on the measurements. While we
cannot be sure what the equivalent $z\sim 0$ $A$ value of an individual galaxy
would be, if investigating the galaxy population at high redshift, the average
$A$ value could be corrected.

Figure 13: Investigating the effect of redshift on the asymmetry (_top_) and
concentration (_bottom_) measurements from both our networks and Morfometryka.
Plots show the distributions of the deviations in the measurements for a
sample of redshifts. It can be seen that there are some extreme outliers in
the measurements made by Morfometryka, especially with the asymmetry
measurements. The range in the first asymmetry plot is reduced to better see
both distributions. The original asymmetry measurements range from
$0.01<A<0.48$, while the concentration measurements vary from $2.4<C<5.0$.
There are a total of 100 galaxies from the Frei sample in this analysis. We
also show the results from Conselice (2003) (the shaded region) demonstrating
that our trends agree well with those found previously.

We have included the equivalent results from Conselice (2003) in Fig.13. In
their redshift test they included 82 galaxies from the Frei sample and
investigated redshift effects up to $z=3$. They implement no size correction
in their redshifting procedure. Nevertheless, we see similar behaviour in both
the bias and scatter, illustrating that the qualitative behaviour is
insensitive to the details of our simulations.

These results indicate that our networks’ measurements at high redshift are
better behaved than Morfometryka, with Morfometryka having a broader range of
values. While both methods show systematic biases, the reduced scatter and
lower prevalence of outliers suggests one could more confidently correct
high-z $C$ and $A$ values based on the trained networks.

The differences seen in this section are greater than might be anticipated
purely due to $SNR$ from the results in Sec. 4.2. This indicates that there
are other factors affecting Morfometryka more than the networks. In the $SNR$
tests, simple uncorrelated Gaussian noise was added to the original images.
However, in the artificial-redshifting procedure, the galaxies are placed in
an apparently empty region of a real CANDELS image. The resulting background
is more realistic, containing low-level structure due to pixel covariances
introduced during the reduction and faint background galaxies. This could
affect the asymmetry and background calculations by Morfometryka, especially
at higher redshifts where the background subtraction becomes more significant.
Our networks, which we have shown to be less susceptible to noise, appear to
be less sensitive to these effects. The result is more stable measurements,
which may be applied at high-redshift with greater confidence.

Looking at the concentration measurements, we see that the network performs
significantly better than the standard algorithm at recovering the original
$C$ values. The network variations are, on average, a factor of two lower than
those measured by Morfometryka. The average change in the concentration
measurements between $z=1$–$3$ is 15% compared to 26% for Morfometryka.
However, as mentioned above, such systematic trends could be corrected. More
importantly, the scatter in the network’s measurements is somewhat lower and
more consistent than those measured by the standard algorithm.

Our networks have been trained on CANDELS data, but successfully applied to
data with a simple noise-degradation, and to artificially-redshifted, ground-
based data. The individual images cover a very wide variety of appearances,
and yet we recover reliable measurements across our test set. This indicates
that our networks are not particularly sensitive to the details of the
observations. They can be applied to roughly comparable datasets with similar
performance to standard methods, without the need for retraining.

The primary reason for this flexibility, is that both $C$ and $A$ measurements
are determined from an image alone, without requiring any other information,
such as the PSF, noise characteristics, etc. The network has learned to
calculate a statistic from the image pixel values, irrespective of the
observational details. It is therefore expected, but still pleasing to see,
that the networks remain accurate when applied to a wide variety of images.

###  4.4 Computational efficiency

We now briefly turn to the efficiency of our CNNs, compared to CAS
measurements using conventional measurements. Running both our trained
networks and Morfometryka on a single computational core, for comparison, our
CNNs are able to produce measurements $\sim 3,000$ times faster. However, for
a modern workstation, containing a single high-end consumer GPU (e.g. an
NVIDIA GeForce GTX 1080 Ti) and 16 CPU cores, the results are even more
striking. On such a system, our trained networks can analyse $\sim 10,000$
galaxies in under 1.5 seconds, while it would take 2 hours to perform these
measurements using the Morfometryka code. Thus, our networks could measure all
1.5 billion resolved galaxies in the Euclid survey (Laureijs et al., 2011) on
a single machine in a little over an hour. To do the same with Morfometryka
would take several weeks on a 1000-CPU cluster! Even with highly-optimized
software, using conventional algorithms would require significant time on a
computing cluster.

In the previous section we have argued that our networks may be applied to
other datasets without needing to retrain, providing the data characteristics
are reasonably similar (which will be the case for any intermediate- to high-
redshift galaxy surveys). However, should retraining be deemed necessary, this
need not be an onerous process. In Sec. 3.4 we show that the network
performance is consistent for moderate variations around the optimum. We
expect that the selected hyperparameters will be a suitable choice for a
variety of datasets. There should be no need to rerun the Bayesian
Optimisation process again.

Given the performance we see for our networks, only a few tens of thousands of
galaxies would be required to retrain the network. Using the optimal
architecture, fully training the network with 75,000 training examples takes
only around 30 minutes. Transfer learning is also a possibility, but the
training time should be no longer. In any case, the network training time is
short compared to that required to prepare the training set, which itself is
substantially faster than applying conventional methods to a large dataset.

A more general argument in favour of moving towards deep learning techniques
for these kinds of calculations, is that there is potential for many,
currently required, preparatory steps to be avoided. Pre-processing steps such
as creating segmentation maps and cleaning neighbouring objects could, in
principle, be performed by the network itself. In this paper we have not
explored this, and have instead applied our networks to the data prepared for
Morfometryka. However, an indication of the networks’ robustness is provided
by its stability when applied to artificially-redshifted galaxies. This could
significantly reduce the computational and human time spent preparing the data
to run these measurements.

It should be noted that Morfometryka performs a number of additional
measurements that are complementary to those discussed in this work. However,
one could train a network, in the same manner presented here, to predict these
parameters. Indeed, this has already been done for Sérsic profiles (Tuccillo
et al., 2018).

A further outstanding issue is that of uncertainties. We have not attempted to
produce uncertainties on individual measurements output by our networks,
beyond examining the scatter relative to Morfometryka. However, Pearson et al.
(2021) give a detailed explanation of how the estimation of uncertainties can
be incorporated into CNNs.

Extending our network to measure a wider variety of parameters, with
uncertainties, is beyond the scope of this present work. However, we hope that
this paper demonstrates the benefits for upcoming ‘Big Data’ surveys, such as
Euclid and Rubin-LSST. Deep-learning has the potential to improve over
conventional approaches in terms of the efficiency, accuracy and flexibility
with which the next generation of surveys can be analysed.

##  5 Summary

In this paper we trained two convolutional neural networks to perform
concentration ($C$) and asymmetry ($A$) measurements based on individual
galaxy input images. Our trained networks reproduce measurements by standard
algorithms with an average absolute error on the $C$ and $A$ values of $0.15$
and $0.045$, respectively. These are lower than the average uncertainties on
those measurements using conventional methods. Our networks can therefore be
used to measure these quantities with a similar level of confidence to
existing algorithms. Analysing these quantities for large samples of galaxies
can provide an estimate of the merger fraction, and help us understand the
transition from peculiar/irregular galaxies at high redshift to the well-
defined Hubble sequence we observe locally.

We have shown how both our networks’ and Morfometryka’s measurements are
impacted by noise, but find that our networks’ estimates are more stable in
the low signal-to-noise regime, in terms of both lower scatter and systematic
bias. By artificially-redshifting a sample of local galaxies from the Frei
catalogue, we investigate trends in the measurements due to redshift effects.
Again, we find that our networks produce measurements with a lower level of
random variation, compared to the conventional algorithms. While the measured
$A$ and $C$ values are slightly biased at high-redshift, our networks and
Morfometryka are both affected in similar manner, and consistent with
behaviour seen previously (Conselice, 2003). Furthermore, the systematic
offsets are comparable to the random uncertainty on individual galaxy
measurements, and so relatively minor.

Our trained networks are up to several thousand times faster than previous
non-parametric measurement algorithms, presenting a substantial advantage for
upcoming surveys. Our trained networks are made public with this work
222https://github.com/cbtohill/CASNET. The future of extragalactic astronomy
consists of ‘Big Data’ surveys, which will image billions of galaxies. Current
state of the art computational methods for analysing these surveys will become
impractical due to the computational resources and time they need. While
detailed analyses will be required for certain measurements, machine learning
techniques can replace many current algorithms. CNN-based approaches are more
efficient and, as we have shown for measuring CAS parameters, can be more
accurate and reliable than traditional measurements. Measuring non-parametric
morphologies in upcoming galaxy surveys, including those by the Euclid, Rubin,
and Roman observatories, will greatly benefit from the methods presented in
this paper. In addition, the high accuracy of our CNN-based measurements make
them equally suitable for use on smaller samples from deeper surveys, such as
those by JWST.

The authors would like to thank the Centre for Astronomy and Particle Theory
of the University of Nottingham for providing the computational infrastructure
needed to produce the networks used in this paper. CT acknowledges funding
from the Science and Technology Facilities Council (STFC). LF acknowledges
funding from Coordenação de Aperfeiçoamento de Pessoal de Nível Superior -
Brasil (CAPES) - Finance Code 001. We thank the anonymous referee for their
thorough review, which helped to significantly improve the presentation of
this work.

## References

  * Abadi et al. (2016) Abadi, M., Agarwal, A., Barham, P., et al. 2016, arXiv e-prints, arXiv:1603.04467.  https://arxiv.org/abs/1603.04467
  * Abraham et al. (1996a) Abraham, R. G., Tanvir, N. R., Santiago, B. X., et al. 1996a, MNRAS, 279, L47, doi: 10.1093/mnras/279.3.L47
  * Abraham et al. (1994) Abraham, R. G., Valdes, F., Yee, H. K. C., & van den Bergh, S. 1994, ApJ, 432, 75, doi: 10.1086/174550
  * Abraham et al. (1996b) Abraham, R. G., van den Bergh, S., Glazebrook, K., et al. 1996b, ApJS, 107, 1, doi: 10.1086/192352
  * Aniyan & Thorat (2017) Aniyan, A. K., & Thorat, K. 2017, ApJS, 230, 20, doi: 10.3847/1538-4365/aa7333
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123, doi: 10.3847/1538-3881/aabc4f
  * Banerji et al. (2010) Banerji, M., Lahav, O., Lintott, C. J., et al. 2010, MNRAS, 406, 342, doi: 10.1111/j.1365-2966.2010.16713.x
  * Barchi et al. (2020) Barchi, P. H., de Carvalho, R. R., Rosa, R. R., et al. 2020, Astronomy and Computing, 30, 100334, doi: 10.1016/j.ascom.2019.100334
  * Barden et al. (2008) Barden, M., Jahnke, K., & Häußler, B. 2008, ApJS, 175, 105, doi: 10.1086/524039
  * Bergstra et al. (2011) Bergstra, J., Bardenet, R., Bengio, Y., & Kégl, B. 2011, in Advances in Neural Information Processing Systems, Vol. 24, 25th Annual Conference on Neural Information Processing Systems (NIPS 2011), ed. J. Shawe-Taylor, R. Zemel, P. Bartlett, F. Pereira, & K. Weinberger (Granada, Spain: Neural Information Processing Systems Foundation).  https://hal.inria.fr/hal-00642998
  * Bershady et al. (2000) Bershady, M. A., Jangren, A., & Conselice, C. J. 2000, AJ, 119, 2645, doi: 10.1086/301386
  * Bluck et al. (2012) Bluck, A. F. L., Conselice, C. J., Buitrago, F. o., et al. 2012, ApJ, 747, 34, doi: 10.1088/0004-637X/747/1/34
  * Boucaud et al. (2020) Boucaud, A., Huertas-Company, M., Heneka, C., et al. 2020, MNRAS, 491, 2481, doi: 10.1093/mnras/stz3056
  * Buitrago et al. (2008) Buitrago, F., Trujillo, I., Conselice, C. J., et al. 2008, ApJ, 687, L61, doi: 10.1086/592836
  * Cheng et al. (2020) Cheng, T.-Y., Conselice, C. J., Aragón-Salamanca, A., et al. 2020, MNRAS, 493, 4209, doi: 10.1093/mnras/staa501
  * Conselice (1997) Conselice, C. J. 1997, PASP, 109, 1251, doi: 10.1086/134004
  * Conselice (2000) —. 2000, arXiv e-prints, astro.  https://arxiv.org/abs/astro-ph/0012454
  * Conselice (2003) —. 2003, ApJS, 147, 1, doi: 10.1086/375001
  * Conselice & Arnold (2009) Conselice, C. J., & Arnold, J. 2009, MNRAS, 397, 208, doi: 10.1111/j.1365-2966.2009.14959.x
  * Conselice et al. (2011) Conselice, C. J., Bluck, A. F. L., Ravindranath, S., et al. 2011, MNRAS, 417, 2770, doi: 10.1111/j.1365-2966.2011.19442.x
  * Conselice et al. (2008) Conselice, C. J., Rajgor, S., & Myers, R. 2008, MNRAS, 386, 909, doi: 10.1111/j.1365-2966.2008.13069.x
  * de Albernaz Ferreira & Ferrari (2018) de Albernaz Ferreira, L., & Ferrari, F. 2018, MNRAS, 473, 2701, doi: 10.1093/mnras/stx2266
  * Dieleman et al. (2015) Dieleman, S., Willett, K. W., & Dambre, J. 2015, MNRAS, 450, 1441, doi: 10.1093/mnras/stv632
  * D’Isanto & Polsterer (2018) D’Isanto, A., & Polsterer, K. L. 2018, A&A, 609, A111, doi: 10.1051/0004-6361/201731326
  * Domínguez Sánchez et al. (2018) Domínguez Sánchez, H., Huertas-Company, M., Bernardi, M., Tuccillo, D., & Fischer, J. L. 2018, MNRAS, 476, 3661, doi: 10.1093/mnras/sty338
  * Duncan et al. (2019) Duncan, K., Conselice, C. J., Mundy, C., et al. 2019, ApJ, 876, 110, doi: 10.3847/1538-4357/ab148a
  * Elmegreen et al. (2005) Elmegreen, D. M., Elmegreen, B. G., Rubin, D. S., & Schaffer, M. A. 2005, ApJ, 631, 85, doi: 10.1086/432502
  * Feinstein et al. (2020) Feinstein, A. D., Montet, B. T., Ansdell, M., et al. 2020, AJ, 160, 219, doi: 10.3847/1538-3881/abac0a
  * Ferrari et al. (2015) Ferrari, F., de Carvalho, R. R., & Trevisan, M. 2015, ApJ, 814, 55, doi: 10.1088/0004-637X/814/1/55
  * Ferreira et al. (2020) Ferreira, L., Conselice, C. J., Duncan, K., et al. 2020, ApJ, 895, 115, doi: 10.3847/1538-4357/ab8f9b
  * Ferreira et al. (2018) Ferreira, L., Ferrari, F., & Griffiths, A. 2018, galclean: v1.0.0, v1.0.0, Zenodo, doi: 10.5281/zenodo.4004571
  * Freeman et al. (2013) Freeman, P. E., Izbicki, R., Lee, A. B., et al. 2013, MNRAS, 434, 282, doi: 10.1093/mnras/stt1016
  * Frei et al. (1996) Frei, Z., Guhathakurta, P., Gunn, J. E., & Tyson, J. A. 1996, AJ, 111, 174, doi: 10.1086/117771
  * Frontera-Pons et al. (2017) Frontera-Pons, J., Sureau, F., Bobin, J., & Le Floc’h, E. 2017, A&A, 603, A60, doi: 10.1051/0004-6361/201630240
  * Graham et al. (2005) Graham, A. W., Driver, S. P., Petrosian, V., et al. 2005, AJ, 130, 1535, doi: 10.1086/444475
  * Grogin et al. (2011) Grogin, N. A., Kocevski, D. D., Faber, S. M., et al. 2011, ApJS, 197, 35, doi: 10.1088/0067-0049/197/2/35
  * Häußler et al. (2013) Häußler, B., Bamford, S. P., Vika, M., et al. 2013, MNRAS, 430, 330, doi: 10.1093/mnras/sts633
  * Hoyos et al. (2012) Hoyos, C., Aragón-Salamanca, A., Gray, M. E., et al. 2012, MNRAS, 419, 2703, doi: 10.1111/j.1365-2966.2011.19918.x
  * Huertas-Company et al. (2008) Huertas-Company, M., Rouan, D., Tasca, L., Soucail, G., & Le Fèvre, O. 2008, A&A, 478, 971, doi: 10.1051/0004-6361:20078625
  * Huertas-Company et al. (2015) Huertas-Company, M., Gravet, R., Cabrera-Vives, G., et al. 2015, ApJS, 221, 8, doi: 10.1088/0067-0049/221/1/8
  * Hutter et al. (2011) Hutter, F., Hoos, H. H., & Leyton-Brown, K. 2011, in Learning and Intelligent Optimization, ed. C. A. C. Coello (Berlin, Heidelberg: Springer Berlin Heidelberg), 507–523, doi: 10.1007/978-3-642-25566-3_40
  * Koekemoer et al. (2011) Koekemoer, A. M., Faber, S. M., Ferguson, H. C., et al. 2011, ApJS, 197, 36, doi: 10.1088/0067-0049/197/2/36
  * Laureijs et al. (2011) Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, arXiv e-prints, arXiv:1110.3193.  https://arxiv.org/abs/1110.3193
  * Lecun et al. (1998) Lecun, Y., Bottou, L., Bengio, Y., & Haffner, P. 1998, Proceedings of the IEEE, 86, 2278
  * Lotz et al. (2008) Lotz, J. M., Jonsson, P., Cox, T. J., & Primack, J. R. 2008, MNRAS, 391, 1137, doi: 10.1111/j.1365-2966.2008.14004.x
  * Lotz et al. (2004a) Lotz, J. M., Madau, P., Giavalisco, M., & Primack, J. 2004a, in American Astronomical Society Meeting Abstracts, Vol. 205, American Astronomical Society Meeting Abstracts, 163.07
  * Lotz et al. (2004b) Lotz, J. M., Primack, J., & Madau, P. 2004b, AJ, 128, 163, doi: 10.1086/421849
  * Mortlock et al. (2013) Mortlock, A., Conselice, C. J., Hartley, W. G., et al. 2013, MNRAS, 433, 1185, doi: 10.1093/mnras/stt793
  * Naim et al. (1995) Naim, A., Lahav, O., Sodre, L., J., & Storrie-Lombardi, M. C. 1995, MNRAS, 275, 567, doi: 10.1093/mnras/275.3.567
  * Nair & Hinton (2010) Nair, V., & Hinton, G. E. 2010, in Proceedings of the 27th International Conference on International Conference on Machine Learning, ICML’10 (Madison, WI, USA: Omnipress), 807–814
  * Nevin et al. (2019) Nevin, R., Blecha, L., Comerford, J., & Greene, J. 2019, ApJ, 872, 76, doi: 10.3847/1538-4357/aafd34
  * Noguchi (1998) Noguchi, M. 1998, Nature, 392, 253, doi: 10.1038/32596
  * Pasquet et al. (2019) Pasquet, J., Bertin, E., Treyer, M., Arnouts, S., & Fouchez, D. 2019, A&A, 621, A26, doi: 10.1051/0004-6361/201833617
  * Pearson et al. (2019) Pearson, J., Li, N., & Dye, S. 2019, MNRAS, 488, 991, doi: 10.1093/mnras/stz1750
  * Pearson et al. (2021) Pearson, J., Maresca, J., Li, N., & Dye, S. 2021, MNRAS, 505, 4362, doi: 10.1093/mnras/stab1547
  * Peng et al. (2002) Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2002, AJ, 124, 266, doi: 10.1086/340952
  * Robotham et al. (2017) Robotham, A. S. G., Taranu, D. S., Tobar, R., Moffett, A., & Driver, S. P. 2017, MNRAS, 466, 1513, doi: 10.1093/mnras/stw3039
  * Sazonova et al. (2020) Sazonova, E., Alatalo, K., Lotz, J., et al. 2020, ApJ, 899, 85, doi: 10.3847/1538-4357/aba42f
  * Schade et al. (1995) Schade, D., Lilly, S. J., Crampton, D., et al. 1995, ApJ, 451, L1, doi: 10.1086/309677
  * Sérsic (1963) Sérsic, J. L. 1963, Boletin de la Asociacion Argentina de Astronomia La Plata Argentina, 6, 41
  * Shorten & Khoshgoftaar (2019) Shorten, C., & Khoshgoftaar, T. 2019, Journal of Big Data, 6, 1, doi: 10.1186/s40537-019-0197-0
  * Simard et al. (2011) Simard, L., Mendel, J. T., Patton, D. R., Ellison, S. L., & McConnachie, A. W. 2011, ApJS, 196, 11, doi: 10.1088/0067-0049/196/1/11
  * Snoek et al. (2015) Snoek, J., Rippel, O., Swersky, K., et al. 2015, arXiv e-prints, arXiv:1502.05700.  https://arxiv.org/abs/1502.05700
  * Storrie-Lombardi et al. (1992) Storrie-Lombardi, M. C., Lahav, O., Sodre, L., J., & Storrie-Lombardi, L. J. 1992, MNRAS, 259, 8P, doi: 10.1093/mnras/259.1.8P
  * The GPyOpt authors (2016) The GPyOpt authors. 2016, GPyOpt: A Bayesian Optimization framework in Python, http://github.com/SheffieldML/GPyOpt
  * Tolman (1930) Tolman, R. C. 1930, Proceedings of the National Academy of Science, 16, 511, doi: 10.1073/pnas.16.7.511
  * Trujillo et al. (2007) Trujillo, I., Conselice, C. J., Bundy, K., et al. 2007, MNRAS, 382, 109, doi: 10.1111/j.1365-2966.2007.12388.x
  * Tuccillo et al. (2018) Tuccillo, D., Huertas-Company, M., Decencière, E., et al. 2018, MNRAS, 475, 894, doi: 10.1093/mnras/stx3186
  * Whitney et al. (2019) Whitney, A., Conselice, C. J., Bhatawdekar, R., & Duncan, K. 2019, ApJ, 887, 113, doi: 10.3847/1538-4357/ab53d4
  * Whitney et al. (2020) Whitney, A., Conselice, C. J., Duncan, K., & Spitler, L. R. 2020, ApJ, 903, 14, doi: 10.3847/1538-4357/abb824
  * Willett et al. (2013) Willett, K. W., Lintott, C. J., Bamford, S. P., et al. 2013, MNRAS, 435, 2835, doi: 10.1093/mnras/stt1458
  * Yagi et al. (2006) Yagi, M., Nakamura, Y., Doi, M., Shimasaku, K., & Okamura, S. 2006, MNRAS, 368, 211, doi: 10.1111/j.1365-2966.2006.10144.x


# Quasirandom groups enjoy interleaved mixing

Harm Derksen Partially supported by NSF grant DMS 2147769. Emanuele Viola
Supported by NSF grant CCF-2114116.

###### Abstract

Let $G$ be a group such that any non-trivial representation has dimension at
least $d$. Let $X=(X_{1},X_{2},\ldots,X_{t})$ and
$Y=(Y_{1},Y_{2},\ldots,Y_{t})$ be distributions over $G^{t}$. Suppose that $X$
is independent from $Y$. We show that for any $g\in G$ we have

$\left|\mathrm{\mathbb{P}}[X_{1}Y_{1}X_{2}Y_{2}\cdots
X_{t}Y_{t}=g]-1/|G|\right|\leq\frac{|G|^{2t-1}}{d^{t-1}}\sqrt{\mathrm{\mathbb{E}}_{h\in
G^{t}}X(h)^{2}}\sqrt{\mathrm{\mathbb{E}}_{h\in G^{t}}Y(h)^{2}}.$

Our results generalize, improve, and simplify previous works.

_Quasirandom groups_ , introduced by Gowers [Gow08], are groups whose non-
trivial representations have large dimension. Multiplication in such groups is
known to behave like a random function in several respects. The prime example
of this is that if $X$ and $Y$ are independent, high-entropy distributions
over a quasirandom group then $XY$ (i.e., sample from each and output the
product) becomes closer to uniform in $L_{2}$ norm. For a discussion of this
result and its many proofs we refer to Section 13 of [Gow17]. Other random-
like behaviors are known with respect to, for example, _progressions_[BHR22]
and _corners_[Aus16] (cf. [Vio19]).

In this work we are interested in a question posed by Miles and Viola [MV13].
Let $X=(X_{1},X_{2})$ and $Y=(Y_{1},Y_{2})$ be high-entropy distributions over
$G^{2}$ such that $X$ is independent from $Y$ (but $X_{1}$ needs not be
independent from $X_{2}$ and $Y_{1}$ needs not be independent from $Y_{2}$).
They asked if the _interleaved product_ $X_{1}Y_{1}X_{2}Y_{2}$ “mixes,” i.e.,
if it is close to uniform, for suitable groups $G$. Their question was
motivated by an application to cryptography (which follows from a positive
answer to a more general question they asked).

Gowers and Viola give a positive answer to this question for non-abelian
simple groups, which are known to be quasirandom. For the special case of
$G=SL(2,q)$ they prove a strong error bound. A simpler exposition of the
latter proof appears in [Vio19]. A follow-up paper by Shalev [Sha16] gives
stronger error bounds for non-abelian simple groups.

These proofs are somewhat complicated and use substantial machinery, and they
only apply to simple groups. Here we give a very short and elementary proof
that applies to any quasirandom group, as stated in the abstract. (But if one
is interested only in $G=SL(2,q)$ and is not willing to use representation
theory, the proof in [GV19] may be more accessible, especially with the
simplification presented in [Vio19].)

To illustrate the bound in the abstract, suppose that $X$ is uniform over a
set of density $\alpha$ and $Y$ is uniform over a set of density $\beta$. Then
the right-hand side is $|G|^{2t-1}\cdot
d^{-t+1}\cdot(\alpha\beta)^{-1/2}/|G|^{2t}=|G|^{-1}\cdot
d^{-t+1}\cdot(\alpha\beta)^{-1/2}.$ Our results also slightly improve the
parameters in the cases where interleaved mixing could be established. For
example for $t>2$ the bounds in [GV19] and [Sha16] have $(\alpha\beta)^{-1}$
instead of $(\alpha\beta)^{-1/2}$.

The paper [GV19] also shows that from interleaved mixing there follow a number
of other results (including the solution to the more general question in
[MV13], thus enabling the motivating application). Hence our results yield
these applications for any quasirandom group. Since this is an immediate
composition of proofs in [GV19] and this paper, we refer the reader to [GV19]
for precise statements.

## Proof of statement in the abstract

We follow standard notation for non-abelian Fourier analysis, see for example
Section 13 of [Gow17] or [GV22]. It suffices to prove the theorem for
$g=1_{G}$. Let $Z$ be a distribution over $G$. By Fourier inversion, and using
that $\rho(1)=I$ and $\widehat{Z}(1)=1/|G|$ we have

$|\mathrm{\mathbb{P}}[Z=1]-1/|G||=|\sum_{\rho}d_{\rho}\mathrm{{tr}}(\widehat{Z}(\rho)\rho(1)^{T})-1/|G||=|\sum_{\rho\neq
1}d_{\rho}\mathrm{{tr}}(\widehat{Z}(\rho))|\leq\sum_{\rho\neq
1}d_{\rho}|\mathrm{{tr}}(\widehat{Z}(\rho))|,$ (1)

where $\rho$ ranges over irreducible representations.

The main claim is that if $Z$ is the interleaved product
$X_{1}Y_{1}X_{2}Y_{2}\cdots X_{t}Y_{t}$ then for any $\rho$

$|\mathrm{{tr}}(\hat{Z}(\rho))|\leq|G|^{2t-1}|\hat{X}(\rho^{\otimes
t})|_{2}|\hat{Y}(\rho^{\otimes t})|_{2}.$ (2)

Assuming the claim the proof is completed as follows. Plugging Inequality (2)
into (1) and multiplying by $d_{\rho}^{t}/d^{t}$ which is $\geq 1$ for
$\rho\neq 1$, the error is at most

$\frac{|G|^{2t-1}}{d^{t-1}}\sum_{\rho\neq
1}\left(d_{\rho}^{t/2}\left|\hat{X}(\rho^{\otimes
t})\right|_{2}\right)\left(d_{\rho}^{t/2}\left|\hat{Y}(\rho^{\otimes
t})\right|_{2}\right).$

By Cauchy-Schwarz this is at most

$\frac{|G|^{2t-1}}{d^{t-1}}\sqrt{\sum_{\rho\neq
1}d_{\rho}^{t}\left|\hat{X}(\rho^{\otimes
t})\right|_{2}^{2}}\sqrt{\sum_{\rho\neq
1}d_{\rho}^{t}\left|\hat{Y}(\rho^{\otimes t})\right|_{2}^{2}}.$

Note that $d_{\rho}^{t}$ is the dimension of $\rho^{\otimes t}$. Each sum can
be bounded above by summing over all irreducible representations. Hence by
Parseval the sum with $X$ is at most $\mathrm{\mathbb{E}}_{h\in
G^{t}}X^{2}(h)$ and the same for $Y$, proving the theorem.

Next we verify Inequality (2). By definition we have

$\hat{Z}(\rho)=\mathrm{\mathbb{E}}_{g}Z(g)\overline{\rho(g)}=\mathrm{\mathbb{E}}_{g}\sum_{g_{1},g_{2},\ldots,g_{2t}:\prod
g_{i}=g}X(g_{1},g_{3},\ldots,g_{2t-1})Y(g_{2},g_{4},\ldots,g_{2t})\overline{\rho(g)}.$

This summation is the same as summing over all $g_{i}$ and setting $g$ to be
the product. Further, because $\rho$ is a representation one has
$\rho(\prod_{i}g_{i})=\prod_{i}\rho(g_{i})$. Hence we get

$\hat{Z}(\rho)=\frac{1}{|G|}\sum_{g_{1},g_{2},\ldots,g_{2t}}X(g_{1},g_{3},\ldots,g_{2t-1})Y(g_{2},g_{4},\ldots,g_{2t})\overline{\prod_{i\leq
2t}\rho(g_{i})}.$

And now the critical equation:

$\displaystyle\mathrm{{tr}}\hat{Z}(\rho)$
$\displaystyle=\sum_{i}\frac{1}{|G|}\sum_{g_{1},g_{2},\ldots,g_{2t}}X(g_{1},g_{3},\ldots,g_{2t-1})Y(g_{2},g_{4},\ldots,g_{2t})\sum_{i_{2},i_{3},\ldots,i_{2t}}\bar{\rho}(g_{1})_{i,i_{2}}\bar{\rho}(g_{2})_{i_{2},i_{3}}\cdots\bar{\rho}(g_{2t})_{i_{2t},i}$
$\displaystyle=\frac{1}{|G|}\sum_{i,i_{2},i_{3},\ldots,i_{2t}}\left(\sum_{g_{1},g_{3},\ldots,g_{2t-1}}X(g_{1},g_{3},\ldots,g_{2t-1})\bar{\rho}(g_{1})_{i,i_{2}}\cdot\bar{\rho}(g_{3})_{i_{3},i_{4}}\cdots\bar{\rho}(g_{2t-1})_{i_{2t-1},i_{2t}}\right)$
$\displaystyle\
\cdot\left(\sum_{g_{2},g_{4},\ldots,g_{2t}}Y(g_{2},g_{4},\ldots,g_{2t})\bar{\rho}(g_{2})_{i_{2},i_{3}}\cdot\bar{\rho}(g_{4})_{i_{4},i_{5}}\cdots\bar{\rho}(g_{2t})_{i_{2t},i}\right)$
$\displaystyle=|G|^{2t-1}\sum_{i,i_{2},i_{3},\ldots,i_{2t}}\left(\hat{X}(\bar{\rho}^{\otimes
t})_{i,i_{2},i_{3},\ldots,i_{2t}}\right)\left(\hat{Y}(\bar{\rho}^{\otimes
t})_{i,i_{2},i_{3},\ldots,i_{2t}}\right).$

Inequality (2) now follows by applying the Cauchy-Schwarz inequality.

## References

  * [Aus16] Tim Austin.  Ajtai-Szemerédi theorems over quasirandom groups.  In Recent trends in combinatorics, volume 159 of IMA Vol. Math. Appl., pages 453–484. Springer, [Cham], 2016.
  * [BHR22] Amey Bhangale, Prahladh Harsha, and Sourya Roy.  Mixing of 3-term progressions in quasirandom groups.  In Mark Braverman, editor, ACM Innovations in Theoretical Computer Science conf. (ITCS), volume 215 of LIPIcs, pages 20:1–20:9. Schloss Dagstuhl - Leibniz-Zentrum für Informatik, 2022.
  * [Gow08] W. T. Gowers.  Quasirandom groups.  Combinatorics, Probability & Computing, 17(3):363–387, 2008\.
  * [Gow17] W. T. Gowers.  Generalizations of Fourier analysis, and how to apply them.  Bull. Amer. Math. Soc. (N.S.), 54(1):1–44, 2017.
  * [GV] W. T. Gowers and Emanuele Viola.  The multiparty communication complexity of interleaved group products.  SIAM J. on Computing.
  * [GV15] W. T. Gowers and Emanuele Viola.  The communication complexity of interleaved group products.  In ACM Symp. on the Theory of Computing (STOC), 2015.
  * [GV19] W. T. Gowers and Emanuele Viola.  Interleaved group products.  SIAM J. on Computing, 48(3):554–580, 2019.  Special issue of FOCS 2016.
  * [GV22] W. T. Gowers and Emanuele Viola.  Mixing in non-quasirandom groups.  In ACM Innovations in Theoretical Computer Science conf. (ITCS), 2022.
  * [MV13] Eric Miles and Emanuele Viola.  Shielding circuits with groups.  In ACM Symp. on the Theory of Computing (STOC), 2013.
  * [Sha16] Aner Shalev.  Mixing, communication complexity and conjectures of Gowers and Viola.  Combinatorics, Probability and Computing, pages 1–13, 6 2016.  arXiv:1601.00795.
  * [Vio19] Emanuele Viola.  Non-abelian combinatorics and communication complexity.  SIGACT News, Complexity Theory Column, 50(3), 2019.


# Evolution of language driven by social dynamics

Moirangthem Shubhakanta Singh1 and R.K. Brojen Singh2<EMAIL_ADDRESS>(Corresponding author) 1 Department of Physics, Manipur University,
Canchipur-795003, Manipur, India.
2School of Computational & Integrative Sciences, Jawaharlal Nehru University,
New Delhi-110067, India.

###### Abstract

The survival of endangered languages in complex language competition depends
on socio-cultural status and honour endowed (by itself and by the other) among
them. The restriction in the endorsement of this honour leads to language
extinction of one language, and rise of the other. Endorsing proper mutual
honour each other trigger the co-existence of language speakers and can save
both languages from extinction. The lost of respect to each other drives the
death of both languages. We found a minimal or critical mutual honour
(a=0.9635) which protects the two languages from extinction. The increase in
mutual honour from this minimal value allows increase in the populations of
the two languages speakers. The state of co-existence of competing languages
abolishes the concept of minority and majority in language competition which
can be obtained by mutual honour. Further, excess biased honour to a
particular language (minority or majority) force the language to extinct. In
mean-field approximation of language competition, magnetization parameter can
be taken as an indicator of survival of a language.

##  I Introduction

Language carries socio-cultural status of the speakers fish , and therefore
trigger competition among the languages mine  which leads to the extinction of
some languages and dominance by some other languages abra . This language
competition could be affected by various factors, namely, decisions in social
consensus migu , interaction of heterogeneous agents in the society wang  and
many other factors. The two language competition modeled by Abrams and
Strogatz abra  proposed the diversity of languages caused by the competition,
which endangers one of the languages leading to extinction. The model has
limitations in describing dynamics of individual speakers considering only the
behavior of population as a whole. The model was extended to bilanguage
competition, which takes into account some individual speakers of both
languages in the dynamics mine . The incorporation of bilangual speakers as a
community in the two-language competition model was then proposed which could
lead to the multi-language modeling mira . On the other hand there have
language models where the two language competition was studied as predator-
prey model pina .

People carry language, and language spread out from the inhabitant region
through interaction, which is proportional to the population size of that
language speakers gong . Competition between any two languages is greatly
influenced by various factors, namely, ability to influence on a language to
the others, inheriting capacity of a language to the people, social and
cultural status of the language. Bilangualism model is one which incorporates
bilangual speakers in the competition of two mutually intelligible or
genetically similar languages allowing to maintain the endangered language
mira . It was proposed, within this model, that even though there is no
endangered monolingual speakers, the existence of bilangual speakers preserve
the endangered language preventing from extinction mine . However, if there is
no linguistic advantage of the two competing languages, the disappearance of
endangered language speakers will lead to extinction of the language due to
conversion of bilangual speakers to the other existing language mira ; mine .
In these models, existence of bilangual speakers in language competition is an
indicator of preserving endangered languages from extinction.

Language competition is greatly influence by social status, cultural
inheritance and hounor given among the speakers. These influencing factors are
so strong that minority languages, which are far different from majority
language, can able to survive if these factors are provided among them. We
address these issues in this work in order to show the probable preventive
measures to be taken up to save endangered languages from extinction. We
studied the bilangual model mine ; mira , which is an extended version of two
language competition due to Abrams and Strogatz abra , to understand the
fundamental roles of the factors in language dynamics and evolution.

Figure 1: The dynamics of two competing languages. Extinction and survival of
languages due to biased social status to the languages (left panel) for
$c_{XB}=0.035,c_{YB}=0.035,c_{BX}=1.0,c_{BY}=1.0$. Switching of language
survival due to change of social status (upper right panel). Low social status
drives both languages to extinct (lower right panel).

##  II Bilangual model of two language competition

The two language competition model due to Abrams and Strogatz abra  was
expended by Minett and Wang mine  by incorporating the emergence of bilangual
speakers from the two language speakers. If $X$ and $Y$ are two competing
languages which permits to emerge a group of people $B$ who speak both the
languages, the competition allows these basic transitions in due course:
$X\rightarrow Y$, $Y\rightarrow X$, $X\rightarrow B$, $B\rightarrow X$,
$Y\rightarrow B$ and $B\rightarrow Y$ mine . The first two transition can be
taken as rare transitions mine , and therefore it can be assumed that only the
last four transitions are trivially observed transitions during competition
wang . If $N$ is the population of the monolinguals and bilangual speakers,
then $X+Y+B=N$, so that one can define $x=X/N$, $y=Y/N$ and $bB/N$ such that
$x+y+b=1$. The transition probability of a transition from state $B$ to $X$ is
given by the power law, $P_{BX}=\mu c_{BX}s_{X}x^{u}$, where, $c_{BX}$ is
attractiveness of the bilangual speakers towards $X$ speakers, $s_{X}$ is
status of the language $X$ and $\mu$ is the mortality rate of $B$, i.e.
replacing rate of $B$ by $X$. Similarly, the probability that the transition
$X\rightarrow B$ occurs is given by, $P_{XB}=(1-\mu)c_{XB}s_{Y}y^{u}$, where,
$s_{X}+s_{Y}=1$. Then the time evolutions of $x$ and $y$ are given by the net
change in $x$ and $y$ at state $X$ and $Y$ respectively given by,

$\displaystyle\frac{dx}{dt}$ $\displaystyle=$ $\displaystyle zP_{BX}-xP_{XB}$
(1) $\displaystyle=$ $\displaystyle\mu
c_{BX}s_{X}(1-x-y)x^{u}-(1-\mu)c_{XB}s_{Y}xy^{u}$ $\displaystyle\frac{dy}{dt}$
$\displaystyle=$ $\displaystyle zP_{BY}-yP_{YB}$ (2) $\displaystyle=$
$\displaystyle\mu c_{BY}s_{Y}(1-x-y)y^{u}-(1-\mu)c_{YB}s_{X}yx^{u}$

The critical points of the above coupled differential equations (1) and (2)
from the conditions, $\frac{dx}{dt}=0,\frac{dy}{dt}=0$. The critical points
($x^{*},y^{*}$) are calculated from this conditions, and found to be (0,0) as
one critical point, and other can be obtained from the equations and the
following expression,

$\displaystyle x^{*}=\Lambda
y^{*};~{}~{}\Lambda=\left(\frac{c_{XB}c_{BY}s_{Y}^{2}}{c_{BX}c_{YB}s_{X}^{2}}\right)^{\frac{1}{2u-1}}$
(3)

This expression and the two equations from $\frac{dx}{dt}=0,\frac{dy}{dt}=0$
lead to the expressions for $x^{*}$ and $y^{*}$ given by,

$\displaystyle x^{*}$ $\displaystyle=$
$\displaystyle\frac{\Lambda}{1+\Lambda+\Gamma}$ (4) $\displaystyle y^{*}$
$\displaystyle=$ $\displaystyle\frac{1}{1+\Lambda+\Gamma}$
$\displaystyle\Gamma=\Lambda\left(\frac{1}{\mu}-1\right)\frac{c_{XB}s_{Y}}{\mu
c_{BX}s_{X}}$

The critical point given by (4) and (5) depends on language bias, social
status and hounor provided.

Figure 2: Mutual honour endorsed to two competing languages allows to co-
existence of both the languages. The dynamics co-existence of both the
languages as a function of $u$ is shown. Figure 3: Plots of two competing
language dynamics as a function of social status of one language $X$ for a
fixed value of $s_{X}$.

##  III Mean-field approximation of language competition

The infinite population limit allows the society fully connected due to
interaction of various language speakers cast . In this limit, one can define
$m=x-y$, which is closely similar to magnetization parameter. If $\nu$ is a
bias parameter, defined by $\nu=s_{X}-s_{Y}$, then putting $s_{X}=s$, we get,
$\nu=2s-1$. Then, the time evolution of magnetization parameter, $m$ can be
obtained from equations (1) and (2). Similarly, from $x+y+b=1$ we have,
$\frac{db}{dt}=-\frac{d}{dt}(x+y)$. The time evolution of $m$ and $b$ are
given by,

$\displaystyle\frac{dm}{dt}$ $\displaystyle=$
$\displaystyle\frac{\mu}{2^{u}}\left[c_{BX}s(1-b+m)^{u}\right.$ (6)
$\displaystyle\left.-c_{BY}(s-\nu)(1-b-m)^{u}\right]$
$\displaystyle+\frac{1-\mu}{2^{u+1}}\left[(1-b)^{2}-m^{2}\right]$
$\displaystyle\times\left[c_{XB}(s-\nu)(1-b-m)^{u-1}\right.$
$\displaystyle\left.-c_{YB}s(1-b+m)^{u-1}\right]$ $\displaystyle\frac{db}{dt}$
$\displaystyle=$
$\displaystyle\frac{\mu}{2^{u}}\left[c_{BX}s(1-b+m)^{u}\right.$ (7)
$\displaystyle\left.+c_{BY}(s-\nu)(1-b-m)^{u}\right]$
$\displaystyle+\frac{\mu-1}{2^{u+1}}\left[(1-b)^{2}-m^{2}\right]$
$\displaystyle\times\left[c_{XB}(s-\nu)(1-b-m)^{u-1}\right.$
$\displaystyle\left.+c_{YB}s(1-b+m)^{u-1}\right]$

The critical points $(m^{*},b^{*})$ of the above differential equations (6)
and (7) can be obtained using the conditions
$\frac{dm}{dt}=0,~{}\frac{db}{dt}=0$. The two equations obtained from these
two conditions give two critical points, (0,1) and the other critical point is
given by,

$\displaystyle m^{*}$ $\displaystyle=$
$\displaystyle\frac{\mu}{\mu-1}\left[\frac{c_{BY}(s-\nu)}{c_{YB}s}-\frac{c_{BX}s}{c_{XB}(s-\nu)}\right]$
(8) $\displaystyle b^{*}$ $\displaystyle=$ $\displaystyle
1-\frac{\mu}{1-\mu}\left[\frac{c_{BX}s}{c_{XB}(s-\nu)}+\frac{c_{BY}(s-\nu)}{c_{YB}s}\right]$
(9)

The above critical point is greatly affected by various parameters, such as,
status, bias and mortality rate of the languages. The population difference of
the two language speakers remains constant which means that there is no
competition between the two languages providing harmony by the bilangual
speakers.

##  IV Results and discussions

The language dynamics driven by various parameters are presented here by
solving the coupled differential equations provided in equations (1), (2), (8)
and (8) using standard fourth order Runge-kutta method of numerical
integration pres . The complicated language dynamics, affected by various
social and cultural factors, is studied and provide some of the possible
criteria to save endangered languages from extinction.

###  IV.1 Language dynamics driven by social status

The social status of the two competing languages, indicated by $s_{X}$ and
$s_{Y}$, influence greatly the attractiveness of each language leading to the
conversion of speakers of one language to the other, and vice versa. If $X$
and $Y$ are minority and majority languages respectively, and the social
status endorsed to both the languages has large difference
($s_{X}=4,~{}s_{Y}=8$), then the speakers of minority language will decrease
monotonically with time (Fig. 1 left panel) and will be extinct after a
particular time. The speakers of majority language will be increased for
certain range of time due to conversion of minority to majority language
speakers, and then remain stationary (Fig. 1 left panel). The extinction time
of the minority language slow if the social status difference, $\Delta
s=(s_{Y}-s_{X})\rightarrow 0$ with large $s_{X},s_{Y}$, and minority language
comes into existence if $\Delta s=0$.

The social status given in the majority language is smaller as compared to the
minority language, $s_{X}\rangle s_{Y}$, then the majority language will
become extinct and minority language can able to survive (Fig. 1 upper right
panel). Further, if $s_{X},s_{Y}\rightarrow small$, then even if $\Delta s=0$,
both minority and majority languages will become extinct (Fig. 1 right lower
panel). Hence the social status to be given to the competing languages should
be given large and same, so that both the language could able to survive.

###  IV.2 Language driven by mutual honour

The mutual honour (indicated by parameter $u$) endorsed to each other in
language competition plays an important role which can able to save both
minority and majority languages from extinction (Fig. 2). Low social status
($s_{X}=0.4,~{}s_{Y}=0.4$) with low mutual honour (large value of $u=1$) given
to each other drive both languages to extinct at different times (Fig. 1 and
2). Keeping this social status fixed, decreasing the value of $u$ allows to
increase the number of the language speakers. The maximum critical value of
$u$ is found to be $u=0.9635$ which can save the both the languages from
extinction (non-zero population). The populations of the two languages
speakers first decreases upto minimum non-zero populations, and then the two
populations become equal and constant as a function of time. After this time
of co-existence there is no distinction between majority and minority
languages because of the same populations of the two languages. This indicates
the co-existence of the two language speakers at this social status, needs
minimum honour between the two language speakers which can be used to increase
the attractiveness of each individual language, and can able to save the
languages from extinction (Fig. 2).

Now, if the value of $u$ decreases the co-existence of the two languages
occurs much earlier, and the numbers of the language speakers are increased.
This increase in the populations of two languages at co-existence is due to
proper honour to each other.

Figure 4: The dynamics of magnetization for different values of mutual honour
parameter $u$ for parameter values
$c_{XB}=0.008,c_{YB}=0.008,c_{BX}=0.035,c_{BY}=0.035,\mu=0.05,\nu=0.3,s=0.35$.

###  IV.3 Social status can save minority language

The population of majority and minority language speakers are calculated as a
function of $s_{X}$ for a fixed value of $s_{Y}$ (Fig. 3). The populations of
the speakers of the two languages fluctuate as a function of social status of
the minority language ($s_{X}$), but the population of minority language
speakers increases as $s_{X}$ increases. At the same time the population of
the majority language decreases as a function of $s_{X}$ with large
fluctuations due to social fluctuations with time.

The co-existence of the two languages happens only at few points, and then
they are parted away from the co-existence. These co-existing points are the
points at which both the languages can be save from extinction. Beyond
this/these points one of the languages will be extinct.

###  IV.4 Magnetization as indicator of language extinction

The magnetization $m$ is calculated as a function of time for different values
of $u$ for fixed values of other parameter values (Fig. 4). Magnetization for
large values of $u$ ($u>1$) shows slow variation with time and the behavior
become steady state. The dynamics of $m$ for small values of $u$ ($u<0.0382$)
become monotonically increased as $u$ decreases. Since $m=x-y$, the monotonic
increase in $m$ with time indicates the dominance of one language over the
other leading to the extinction of one language.

The increase in $m$ is due to language competition, and reaching steady state
indicates the co-existence of the two languages where conversion of one
language speaker to the other and vice versa does not takes place. The
decrease in $m$ for $u=2$ as a function of time indicates either both the
languages extinct or co-existence of the two languages (same populations of
the language speakers).

##  V Conclusion

Language dynamics is a complicated process, which involves competition,
change, evolution and struggle for survival. The dynamics of competing
languages is affected by various factors, namely, mutual honour, prestige,
social status and carrying culture which can drive the languages at different
states, i.e. extinct or survival. In competing two languages, the language
which has relatively low social status will quickly extinct and the language
having relatively higher social status will able to survive in long run.
However, if the social status of both the languages are smaller than a
critical value, both the competing languages will become extinct.

The survival of bilangual speakers in two language competition could be a
signature of maintaining minority language to survive by them along with
majority language. It could be due to the fact that the minority language is
preserved by these speakers along with the other language. However, if the
status of the minority language is very low as compared to the other, then the
bilingual speakers will be attracted and convert quickly towards the higher
status language, and the other language will be extinct. The competition
between the two languages is done via transition through bilingual speakers.

The co-existence between two languages is one peculiar state of competing
language dynamics which allows to survive both the languages with same
population. This is possible only when the mutual honour between each other is
increased above a critical value. Hence the survival of both the languages
could be possible if the two competing languages endorse mutual honour and
respect.

Magnetization in mean-field approximation can be used as an indicator of
language survival in language competition. The monotonic increase in
magnetization shows the survival of one language and the other language will
extinct. The steady state of dynamics of magnetization could be a state of co-
existence of the two language. The decrease in the value of this parameter may
indicate extinction of both languages or co-existence state. However, the
language competition is a complicated process, and there are lots of issues in
understanding this process so that one can save endangered languages. The
breakdown of a language into a number of dialects, which may become
independent languages, and understanding their complicated dynamics need to be
investigated to highlight importance of language competition.

Acknowledgments
The author thanks the Special Centre for Study of North-East India, Jawaharlal
Nehru University for inviting me to present this work and exciting discussion.

## References

  * (1) Fishman, J.A. Reversing Language Shift: Theoretical and Empirical Foundations of Assistance to Threatened Languages. Multilingual Matters, Clevedon, Philadelphia (1991).
  * (2) Minett J.W., William S-Y. Wang W.S-Y. Lingua 118, 19–45 (2008).
  * (3) Abrams, D. M. and Strogatz, S. H. (2003). Modelling the dynamics of language death. Nature, 424, 900.
  * (4) Miguel M.S., Eguiluz V., Toral R. and Klemm K. Binary and multi-variate stochastic models of consensus formation. Computer in Science and Engineering, 7, 67–73 (2005).
  * (5) Wang W. S.-Y. and Minett J. W.. The invasion of language: emergence, change and death. Trends in Ecology and Evolution, 20, 263-269 (2005).
  * (6) Mira J. and Paredes A. Europhys. Lett. 69, 1031-1035 (2005).
  * (7) Pinasco J. and Romanelli L. Physica A 361, 355 (2006).
  * (8) Gong T., Shuai L. and Zhang M. Physics of Life Reviews 11, 280–302 (2014).
  * (9) Castelló X., Loureiro-Porto L., Eguiluz V.M. and Miguel M.S. Advancing Social Simulation: The First World Congress, 83-94 (2007).
  * (10) Press WH ,Teukolsky SA, Vetterling WT, Flannery BP. Numerical Recipe in Fortran. Cambridge University Press. 1992;1. ISBN: 9780521430647