<figure><img src="image_1.png" /></figure>

HAL
open science

# Measuring the Intracluster Light Fraction with Machine Learning

Louisa Canepa, Sarah Brough, Francois Lanusse, Mireia Montes, Nina Hatch

► To cite this version:

Louisa Canepa, Sarah Brough, Francois Lanusse, Mireia Montes, Nina Hatch. Measuring the Intracluster Light Fraction with Machine Learning. The Astrophysical Journal, 2025, 980 (2), pp.245. 10.3847/1538-4357/adabc7. hal-05412127

HAL Id: hal-05412127

https://hal.science/hal-05412127v1

Submitted on 11 Dec 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Measuring the intracluster light fraction with machine learning

LOUISA CANEPA,¹,² SARAH BROUGH,¹,² FRANCOIS LANUSSE,³,⁴ MIREIA MONTES,⁵,⁶,⁷ AND NINA HATCH⁸

¹School of Physics, University of New South Wales, NSW 2052, Australia

²ARC Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO 3D)

³The Flatiron Institute, 162 5th Ave, New York, NY 10010, USA

⁴Université Paris-Saclay, Université Paris Cité, Paris 91190, France

⁵Institute of Space Sciences (ICE, CSIC), Campus UAB, Carrer de Can Magrans, s/n, 08193 08193, Cerdanyola del Valles, Spain

⁶Instituto de Astrofísica de Canarias, Calle Vía Láctea s/n, 38204, San Cristóbal de La Laguna, Tenerife, Spain

⁷Departamento de Astrofísica, Universidad de La Laguna, 38206, La Laguna, Tenerife, Spain

⁸School of Physics and Astronomy, University of Nottingham, University Park, Nottingham NG7 2RD, UK

ABSTRACT

The intracluster light (ICL) is an important tracer of a galaxy cluster's history and past interactions. However, only small samples have been studied to date due to its very low surface brightness and the heavy manual involvement required for the majority of measurement algorithms. Upcoming large imaging surveys such as the Vera C. Rubin Observatory's Legacy Survey of Space and Time are expected to vastly expand available samples of deep cluster images. However, to process this increased amount of data, we need faster, fully automated methods to streamline the measurement process. This paper presents a machine learning model designed to automatically measure the ICL fraction in large samples of images, with no manual preprocessing required. We train the fully supervised model on a training dataset of 50,000 images with injected artificial ICL profiles. We then transfer its learning onto real data by fine-tuning with a sample of 101 real clusters with their ICL fraction measured manually using the surface brightness threshold method. With this process, the model is able to effectively learn the task and then adapt its learning to real cluster images. Our model can be directly applied to Hyper Suprime-Cam images, processing up to 500 images in a matter of seconds on a single GPU, or fine-tuned for other imaging surveys such as LSST, with the fine-tuning process taking just 3 minutes. The model could also be retrained to match other ICL measurement methods. Our model and the code for training it is made available on GitHub.

1. INTRODUCTION

Galaxy groups and clusters contain a diffuse, low surface brightness component that extends throughout the cluster (e.g. Mihos et al. 2005). This component is known as the intracluster light (ICL; see Mihos 2019; Contini 2021; Montes 2022 for reviews), and is thought to be mainly made up of stars that have been ripped out of their original galaxies during the accretion history of the cluster, through mergers and interactions with other galaxies (e.g. Rudick et al. 2006; Contini et al. 2014). In-situ star formation might also be a contributor, albeit less significant, to the ICL (e.g. Barfety et al. 2022; Montenegro-Taborda et al. 2023). The ICL forms an important record of the formation history of the cluster and the evolution of the galaxies within it, particularly the brightest cluster galaxy (BCG), around which the ICL is concentrated (e.g. Rudick et al. 2011; Cañas et al. 2020).

One metric that is widely used to quantify the amount of ICL in a cluster is the ICL fraction, that is, the fraction of the total stellar light in the cluster that belongs to the ICL. The relationship of this quantity with host cluster properties such as redshift, mass, or dynamical state can provide key insights into the efficiency of interactions in different environments.

Simulations suggest an increasing amount of ICL over time (with decreasing redshift), as more stars are added to the ICL as the cluster evolves (e.g. Murante et al. 2007). However, this relationship could be non-monotonic, with accretion of other groups and clusters with lower ICL fractions causing the ICL fraction of the combined system to decrease, while disruptions of the individual galaxies already within the cluster cause the fraction to increase (Cañas et al. 2020). The overall relationship is then the result of the interplay between the timescales of these processes. Simulations disagree about the expected strength of this relationship how-
ever, with Montes (2022) showing that, for example, Rudick et al. (2011) found a much more gradual increase across redshift than Tang et al. (2018).

Some simulations suggest a weakly increasing correlation between ICL fraction and host cluster halo mass in the group-cluster mass range (e.g. Purcell et al. 2007; Cañas et al. 2020; Proctor et al. 2024), albeit with significant scatter, while others do not find any trend (e.g. Rudick et al. 2011; Contini et al. 2014; Contreras-Santos et al. 2024).

Simulations also suggest a relationship between the amount of ICL in a cluster and the cluster's dynamical state, where relaxed and hence evolved clusters have higher ICL fractions than their counterparts, although noting that simulations define the ICL fraction in mass rather than stellar light(e.g. Contini et al. 2023; Brough et al. 2024; Contreras-Santos et al. 2024).

Analyzing these trends with ICL fraction and cluster properties using observational data presents difficulties. The observational definition of the ICL is ambiguous, and it is unclear how to best separate the light in the ICL from the rest of the light in the cluster, in particular the BCG. There are a number of different methods currently used. These range from using a surface brightness threshold to define and separate the ICL (e.g. Feldmeier et al. 2004; Burke et al. 2015; Montes & Trujillo 2018; Furnell et al. 2021; Martínez-Lombilla et al. 2023), to modeling and subtracting cluster galaxies through either analytical models (e.g. Gonzalez et al. 2005; Spavone et al. 2017; Morishita et al. 2017; Ragusa et al. 2021) or orthonormal mathematical bases (e.g. Jiménez-Teja & Dupke 2016; Jiménez-Teja et al. 2023), to using wavelet-based algorithms (e.g. Da Rocha & De Oliveira 2005; Guennou et al. 2012; Adami et al. 2013; Ellien et al. 2021). These each come with their own assumptions and biases, which makes it particularly challenging to compare between different existing samples and study trends on a large scale.

The review by Montes (2022) compared measurements taken with the surface brightness threshold method, composite model method, and wavelet method. They found that the observed relationship between ICL fraction and redshift was strongly dependent on the measurement method chosen, with the surface brightness threshold method showing an increase in ICL fraction with decreasing redshift, the composite model method showing slight evolution in redshift, and the final method showing no correlation. There was also significant scatter between different samples due to the different depths, bands, and surface brightness thresholds applied in the images. These issues make it difficult to compare with expectations from simulations and to com-

pare different observations with one another. However, with the heterogeneous samples available, observations also seem to indicate that there is no trend of ICL fraction with host halo mass (e.g. Krick & Bernstein 2007; Montes 2022; Ragusa et al. 2023). Finally, there are also contradictory results from observations regarding the correlation with cluster relaxation – some studies (e.g. Da Rocha et al. 2008; Iodice et al. 2020; Ragusa et al. 2023) found that the ICL fraction appears to be higher in more dynamically evolved systems, whereas Jiménez-Teja et al. (2018) find the opposite trend. These results show the difficulty posed by drawing conclusions from heterogeneous samples.

Brough et al. (2024) compared eight different observational measurement methods applied to a homogeneous sample of simulations. They found agreement between the mean ICL fraction measured by the different methods, but found significant scatter present on an individual measurement basis. The largest difference in ICL fraction for a single cluster was 23% (ranging from 0.11 to 0.34). These results indicate that using a homogeneous sample does alleviate some of the difficulties in comparison, and that the measurement methods do show broad agreement, however the scatter involved still makes direct comparison between samples analyzed with different methods challenging.

These issues are compounded by the dependence of the ICL fraction on the photometric band and survey limiting depth used for the observation (e.g. Burke et al. 2015; Montes & Trujillo 2018). Constructing a large, homogeneous sample of ICL measurements using a single method would be the most effective way to remove these systematic issues and draw more robust conclusions about the ICL fraction's dependence on host cluster properties.

To generate such a sample, we need many homogeneous cluster images with the depth required to study this low surface brightness component. With the next generation of optical imaging surveys expected to reach unprecedented depths, we will soon have enough data to solve this issue. The Vera C. Rubin Observatory's Legacy Survey of Space and Time (LSST; Ivezić et al. 2019) will image the entire Southern sky to a depth sufficient to study the ICL (Montes 2019; Brough et al. 2020). Thus, to take full advantage of these data and avoid the problems presented by different measurement methods, we need to be able to efficiently measure the ICL fraction in hundreds of thousands of group and cluster images, as LSST is likely to find ~100,000 clusters and ~1 million groups up to z ~ 1.2 over the southern hemisphere (Brough et al. 2020). However, the measurement methods currently available are not readily
scalable to sample sizes of this magnitude, as many of them rely on some manual involvement and tuning on a per-cluster basis (e.g. Brough et al. 2024). Alongside investigation into making these methods applicable to significantly larger samples, we need to explore other possible options.

Machine learning has proven to be particularly useful for these types of data-intensive problems, as it requires little to no human involvement, making it an attractive solution for effectively processing large numbers of images. Machine learning models have been developed for a wide array of problems in astronomy, such as image denoising (e.g. Vojtekova et al. 2021) and outlier identification (e.g. Han et al. 2022). In the ICL field, Marini et al. (2022) used a machine learning model to classify star particles in simulated galaxy clusters as belonging to either the ICL or the BCG based on the given host cluster mass, normalized particle cluster-centric distance, and rest-frame velocity. Machine learning models have seen particular use in image classification problems, for example identification of galaxy mergers (Pearson et al. 2019), morphological classification of galaxies (Hayat et al. 2021), and classification of tidal features (Desmons et al. 2024). They have also seen use in image regression problems, for example in predicting photometric redshifts from input images (Hayat et al. 2021).

A difficulty presented in astronomy applications for machine learning, in particular ICL measurement, is the lack of labeled training samples. Supervised machine learning models often require tens of thousands of training samples in order to learn a task sufficiently well, and such datasets do not exist even for current ICL measurement methods. Some works (e.g. Hayat et al. 2021; Stein et al. 2022; Desmons et al. 2024) solve this problem by using self-supervised representation learning, which uses a combination of a large unlabeled dataset and a small labeled dataset to learn a task. However, this technique is not currently suitable for ICL measurement, owing to the limited numbers of deep images of clusters available to construct a large enough unlabeled dataset. Instead, here we explore a technique called transfer learning, which involves adapting learning from a different task or dataset to the target domain (Pan & Yang 2010). It is a well known issue that machine learning models struggle with different datasets. For example, Pearson et al. (2019) found that the accuracy of a model trained on real images fell from 91.5% to 53.0% when applied to simulated mock images. However, if the two datasets are similar enough, this issue can be mitigated with fine-tuning (e.g. Oquab et al. 2014; Walmsley et al. 2022). This involves taking a

model pretrained on a very large dataset related to the target task, and training it further on a small sample of images from the target dataset. Generally, either only the final few layers are trained or a low learning rate is used, to avoid destroying all the learned weights from the previous task.

In this work, we develop a machine learning model to efficiently measure the ICL fraction in input cluster images. Section 2 describes our data sources, outlines the methods used to construct our training and fine-tuning datasets, and describes the architecture of the machine learning model. Section 3 presents the results of our model, and Section 4 discusses the reliability of the results and potential areas of improvement for the model. Throughout this work we assume a standard cosmological model with parameters $H_0 = 70 \text{ km s}^{-1} \text{ Mpc}^{-1}$, $\Omega_m = 0.3$, and $\Omega_\Lambda = 0.7$.

## 2. METHOD

### 2.1. Data sources

We use images from the Hyper Suprime-Cam Subaru Strategic Program (HSC-SSP; Miyazaki et al. 2018) Public Data Release 2 (PDR2; Aihara et al. 2019). Although Public Data Release 3 (PDR3; Aihara et al. 2022) is now available, we choose to use PDR2 as it fulfils the requirements of our work and has been extensively tested for low surface-brightness studies (e.g. Huang et al. 2018; Huang et al. 2020; Li et al. 2022) whereas PDR3 introduces several differences in sky conditions, observing strategy, and data treatment pipeline that would need to be accounted for. HSC-SSP is carried out with the Hyper Suprime-Cam instrument on the 8.2m Subaru telescope, and is a three layered multi-band (grizy) imaging survey. The three layers have varying limiting depths: Wide ($m_r \sim 26.4$ mag), Deep ($m_r \sim 27.4$ mag), and Ultradeep ($m_r \sim 28.0$ mag). LSST is expected to image the entire southern sky to a surface brightness depth<fnref n="1" /> of $\mu_r \sim 30.3 \text{ mag/arcsec}^2$ ($3\sigma$, $10'' \times 10''$), which is comparable to the Ultradeep layer of HSC, which has a depth of $\mu_r \sim 29.8 \text{ mag/arcsec}^2$ (Martínez-Lombilla et al. 2023). The Deep and Ultradeep layers are imaged to a depth sufficient to study low surface brightness features such as ICL, and also has the advantage that it uses the same pipeline as LSST to reduce its data (Bosch et al. 2018), meaning that a model trained on HSC-SSP data should be easily adapted to LSST data as well. This is important as the need for an efficient measuring algorithm such as the one we present

<footnote n="1">¹ https://smtn-016.lsst.io/</footnote>
here is driven by the amount of data that LSST will provide.

In order to make manual measurements of the ICL, we require a catalogue of clusters. We use CAMIRA (Oguri et al. 2018), an optically-selected catalogue of clusters run on all HSC-SSP images. The algorithm, fully described in Oguri (2014), leverages the expectation that all massive clusters are expected to show a “red sequence” of galaxies (e.g. Gladders & Yee 2000) to identify likely red sequence galaxies using a stellar population synthesis model and group them into clusters.

The final product of the CAMIRA algorithm is a catalogue of clusters, where the central coordinates of each cluster are the coordinates of the BCG identified by the algorithm. We only use the Deep/Ultradeep clusters for our manual ICL measurements, as the depths of these layers are best suited for low surface brightness studies. The Wide clusters are used to define a BCG stellar mass distribution, which is used in the construction of our training sample as described below. In the Deep/Ultradeep layers 248 clusters are identified by the CAMIRA algorithm, and 7939 clusters are identified in the Wide layer. Photometric cluster redshifts are calculated by the CAMIRA algorithm based on the most likely redshift for a particular likely red sequence of galaxies, as described in detail in Oguri (2014). We restrict our sample to clusters that have a photometric redshift of less than 0.5, as low surface brightness features become extremely challenging to detect past this limit due to cosmological dimming. 130 clusters in the Deep/Ultradeep layers and 2709 clusters in the Wide layer fit this criterion. These clusters were identified in the HSC-SSP PDR3 survey, which covers more area than the PDR2 that we use for our image data (1298 square degrees compared to 1022 square degrees in the Wide layer, and 36 square degrees compared to 35 square degrees in the Deep/Ultradeep layers). Using PDR2 images reduces our sample size to 129 clusters in the Deep layer and 2379 clusters in the Wide layer. We also remove clusters from these layers that have bad photometry or have significant contamination from bright stars (identified by manual inspection of these images), leaving us with 101 clusters in the Deep/Ultradeep layers.

The redshifts of these clusters have a range of $0.1 < z < 0.5$, where 0.5 is our chosen upper limit as described previously and 0.1 is the lower limit of the CAMIRA catalogue due to the large angular sizes and brightness of member galaxies in HSC-SSP below this redshift complicating the cluster finding (Oguri et al. 2018). The parameter richness in CAMIRA is defined as the number of red member galaxies with stellar masses $M_* \gtrsim$

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Distribution of cluster richness (top) and photometric cluster redshift (bottom) of the final Wide sample (left) and Deep/Ultradeep sample (right). Values are as provided by the CAMIRA algorithm.</figcaption></figure>

$10^{10.2} M_{\odot}$, within an aperture of radius $R \sim 1h^{-1}$Mpc, where the stellar mass of all likely red sequence galaxies is predicted by the stellar population synthesis model (Oguri 2014). In our final sample, the cluster richness ranges from 15 to 93 in the Deep/Ultradeep layers and 15 to 170 in the Wide layer. This lower limit roughly corresponds to $M_{200m} \gtrsim 10^{14}h^{-1}M_{\odot}$ (Oguri et al. 2018). The richness and redshift distribution of the samples over all three layers is shown in Figure 1.

The CAMIRA algorithm also produces a catalogue of 2,191,333 galaxies from the Wide layer that fit the red sequence stellar population synthesis model well, called Luminous Red Galaxies (LRGs), which we also make use of in our work. As with the clusters, we exclude galaxies with redshifts higher than 0.5, leaving a sample of 479,435 LRGs.

Finally, in this work we make use of the photometric redshifts provided by HSC-SSP (Nishizawa et al. 2020), calculated using the Mizuki template fitting code (Tanaka et al. 2018), which gives photometric redshifts for all objects in the Deep/Ultradeep layers of HSC-SSP.
Machine learning models require significant amounts of data (ideally tens of thousands of samples) to effectively learn a task. Given that we only have a sample of 2480 cluster images over the entire HSC-SSP survey, and only 101 of these in the Deep/Ultradeep layer, it is currently not possible to assemble such a sample for training. We use transfer learning to overcome this issue by generating a large sample of artificial images to use as our main training dataset, and then transfer the model’s learning into the real domain using the smaller sample of real images.

Transfer learning is most effective when the training domain is as similar as possible to the target domain. For this reason, rather than using simulated data, which would give us realistic ICL morphologies based on real physical processes and the interaction history of the cluster, we choose to prioritize the realistic appearance of the images by using existing HSC-SSP images as a base and injecting simple artificial ICL profiles.

We use the LRGs identified by the CAMIRA algorithm in the Wide layer as our base images of galaxies. We use these instead of injecting many different artificial ICL profiles into repeat versions of the real BCG images to have more diversity in the training sample for more robust learning. Using the limited number of real BCGs in the Deep layers as the base images could lead to the model overfitting to those particular clusters, which would ultimately hurt its ability to generalize to unseen data. BCGs are typically a special case of the LRG population (e.g. Loh & Strauss 2006; Brough et al. 2007; Dalal et al. 2021), meaning that the generated images, centered on LRGs, will look more similar to the real images, centered on BCGs, compared to a sample of randomly drawn galaxies.

The stellar mass distribution of the LRG sample is significantly different to the combined sample of BCGs from all layers, as shown in Figure 2. The reason that the number of BCGs in high mass bins exceeds the number of LRGs is that the BCG sample is drawn from all three layers, whereas the LRG sample is only from the Wide layer. This difference in the two distributions means that a random sample drawn from images of LRGs will contain many more low mass galaxies than the BCG sample. In order to keep our datasets as similar as possible, we re-sample the LRGs to match the mass distribution of our BCG sample. To do this, we estimate the probability density function (PDF) of the BCG mass distribution and the LRG mass distribution using Gaussian kernel density estimation. We then weight the probability of sampling each galaxy by the ratio of the BCG mass PDF to the LRG mass PDF at that galaxy’s mass. In this way, the likelihood of drawing an LRG of

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Stellar mass distribution of the Wide layer LRG sample with redshift < 0.5, and the combined BCG sample in Wide, Deep, and Ultradeep layers.</figcaption></figure>

the same mass as a BCG increases. We sample 50,000 LRGs using this method, as this should be a sufficient number for the supervised training. After this weighted sampling of the LRGs, the stellar mass distribution of our drawn sample of 50,000 LRGs matches the mass distribution of the BCGs.

We use 600 × 600 kpc cutouts in r-band from the Wide layer of the HSC-SSP survey, converting this into an angular size for downloading each cluster. Given the HSC pixel scale of 0.168 arcsec, this covers a range of 1936 pixels at $z = 0.1$ to 585 pixels at $z = 0.5$. Current observations generally detect ICL out to radii of between 100 kpc and 300 kpc depending on the survey, before it becomes too faint compared to the background level to be detected (e.g. Montes & Trujillo 2014; DeMaio et al. 2015; Jiménez-Teja et al. 2018; Montes et al. 2021). This physical size should therefore be large enough to capture the ICL in the majority of cases. We then mask out bright stars in the image using the bright star masks provided by the HSC-SSP survey.

The method of data generation is illustrated in Figure 3. Firstly, we remove all light in the image fainter than a surface brightness of $\mu_r = 26$ mag/arcsec<fnref n="2" /> (e.g. Rudick et al. 2011; Montes & Trujillo 2018; Martínez-Lombilla et al. 2023), corrected for cosmological dimming and k-corrected following Equation 1, where the k-correction, $k(z)$, is estimated as a function of redshift using the “k-corrections calculator”<fnref n="2" /> (Chilingarian et al. 2010) and assuming an ICL color $g - r = 0.7$ (Martínez-

<footnote n="2">² http://kcor.sai.msu.ru/</footnote>
Lombilla et al. 2023).

$$ \mu_{obs}(z) = 26 + 2.5 \log_{10}(1+z)^4 + k(z) \quad (1) $$

This removes any existing low surface brightness features around the central galaxy, including any potential ICL (see panel 2 in Figure 3).

We then randomly generate an artificial ICL profile, using an exponential profile with randomized ellipticity and position angle (c.f. Watkins et al. 2024). We apply a surface brightness cut to find the 26 mag/arcsec² isophote of the LRG, and set the effective radius of the profile to a randomly selected pixel size between the semi-major and semi-minor axes of the isophote. Several studies have found that the ICL can be modeled using an exponential profile (e.g. Montes & Trujillo 2018; Kluge et al. 2021; Ragusa et al. 2021). We visually determined that our parameters automatically produce a reasonable ICL profile for a number of different LRG base images. We also tested the impact of modeling ICL with an exponential profile by training a model on images injected with Sérsic profiles with random indices between 1 and 8. We find that this does have a minor impact on the model’s performance, with a Mean Absolute Error (MAE) in its predictions of 0.0186 compared to our final model’s MAE of 0.0159. This means that for an ICL fraction of 0.2, the error on the fraction increases from 7.9% to 9.3%. This indicates that even using wildly unrealistic profiles during initial model training does not have a significant impact on the model’s ability to recognize the ICL profiles in real images. However, we note that the use of this profile could create the possibility that our model may not be as adept at identifying ICL with more asymmetric morphologies (e.g. Rudick et al. 2006; Presotto et al. 2014; Montes & Trujillo 2022). We then inject this profile into the thresholded image (see panel 3 in Figure 3).

In the first step, we removed all low surface brightness light in the image. This is intended to remove any low surface brightness light initially present around the central galaxy, which is replaced by our injected ICL profile. However, this will also remove any other low surface brightness light in the image associated with other sources. Removing this light gives these other sources an unnatural appearance that could hinder the model’s ability to process real images. To test the impact of this step, we trained a model on images that instead have Gaussian noise injected as the background. We find that this causes the model to fail to make predictions on real data. This shows that maintaining the realistic appearance and noise properties of the images is important to enabling the model to transfer its learning to the target dataset. We therefore want to preserve this

light in our generated training images. To do this, we create a version of the ICL profile which is normalized between 0 and 1 (where 1 corresponds to the maximum flux value, equivalent to 26 mag/arcsec²). This is then inverted, and used to scale the original subtracted background before adding it to the image. This means that where the ICL profile is close to 0, away from the central galaxy, the image has its original background, but as the profile gets brighter, the original background is suppressed, leaving our injected ICL profile around the central galaxy. Finally, so that the ICL does not appear unnaturally noiseless as in panel 3 of Figure 3, we calculate the 3σ clipped statistics of the original background, and add Gaussian noise with the same standard deviation to the image, scaled by the normalized ICL profile so that the noise across the new image remains consistent with the original. This new image is the final generated image (see panel 4 in Figure 3).

Running this process on all 50,000 images of LRGs gives us our final set of images for the initial stage of the model training. We refer to this dataset as the “artificial dataset”, in contrast to the “fine-tuning dataset” which consists of the sample of 101 real Deep/Ultradeep cluster images that we introduce below.

## 2.3. Measurement method

In order to train the model to predict the ICL fraction, we first need to measure the ICL fractions in both our artificial dataset and the fine-tuning dataset. There are multiple different possible methods to choose from to measure the ICL as discussed in Montes (2022) and Brough et al. (2024). In this case we use the surface brightness cut method to make these measurements. This method assumes that all flux fainter than a chosen surface brightness threshold is contributed by the ICL component. Brough et al. (2024) found that this method gave results consistent with other measurement methods. We choose this method for its relative simplicity and the fact that it makes no assumptions about the morphology of either the BCG or the ICL.

We measure the ICL using the surface brightness cut method for both our artificial and fine-tuning datasets, and describe our procedures in more detail below. The exact procedure differs slightly between the two. For the fine-tuning dataset, we manually tune the background subtraction and bright star masking, and also apply masks for non-cluster member galaxies on an individual cluster basis when conducting the measurement. The fine-tuning dataset is small enough that manual inspection and tuning of each cluster is feasible, and it is important to have accurate measurements for this dataset to achieve a good transfer of the model’s learning
<figure><img src="image_5.png" /><figcaption>**Figure 3.** From left to right the images show the steps of the training data generation procedure. (1) The original cutout of an LRG from HSC-PDR2 wide data, resized to 224x224 pixels. (2) The image with light fainter than 26 mag/arcsec² removed. (3) Randomly generated exponential profile added to the thresholded image. (4) Final generated image with extra noise and background added. (5) The difference between images (1) and (4) to better highlight the added low surface brightness light. All images have been arcsinh scaled.</figcaption></figure>

into the final image domain. However, for the artificial dataset, we prioritize making the measurements completely automatically, as it is not feasible to manually tweak the measurement for 50,000 clusters. The relative simplicity of the method makes it easy to automate the process to some extent, but it is quite reliant on effective masking and background subtraction, which often requires human evaluation and involvement and makes it problematic to apply automatically. Also, since the artificial dataset does not contain actual clusters, we cannot apply cluster membership information here. This means that the automatic measurements taken of the artificial dataset are likely to be slightly inaccurate for some clusters due to these factors. However, given that the purpose of the artificial dataset is to give the model only a general understanding of the task, and that inaccuracies in its understanding should be corrected in the fine-tuning stage, this level of noise in the measurements should be acceptable for our artificial dataset. We explore the impact of this assumption in the next section.

### 2.3.1. Measurement method for artificial dataset

The model accepts images of size 224 × 224 pixels. This is a standard size for computer vision models, to lower the computational cost of training and inference, so we first resize the cutouts to this size. This is done using `SCIKIT-IMAGE` (van der Walt et al. 2014), which interpolates the image to downsize the image. This resizing does mean that lower redshift clusters will be more heavily binned than higher redshift clusters due to their larger initial size, which could affect the measured ICL fraction. We return to this in Section 4.

We then perform an automatic background estimate and subtraction. The HSC-SSP data has already been coadded and sky subtracted using a global sky subtraction algorithm (Aihara et al. 2019), however this often requires some additional correction, particularly to correct for gradients across the cutouts which would be problematic for ICL measurement. To estimate the background, we divide the cutout into a 14 × 14 grid,

and use `PHOTUTILS` (Bradley et al. 2023) to create a coarse 2D background map using the $3\sigma$ clipped statistics in each box, interpolating over bright stars. This is likely to include any ICL and the BCG in the center of the image as part of the background, as they extend over large areas in the image. To avoid this, we keep only the background values from boxes around the edge of the image, removing the central region of the coarse background map. We first calculate the $3\sigma$ clipped median of the counts in these boxes and subtract this constant sky value. We then fit a bivariate B-spline surface to these edge background values and use this as an estimate of the background. This avoids including any large features such as ICL in the center of the image as part of the background, and will smooth over any remaining sharp peaks in the background values due to extended galaxies on the edge of the image, while still allowing us to capture any large-scale 2D gradients across the image. We then subtract this estimate of the background from the image. As mentioned above, it is possible that this automatic procedure could lead to some over or under-subtraction in some images, leading to an increase in the uncertainties. However our procedure aims to minimize these issues as much as possible.

To test the sensitivity of the model's learning to the choice of background subtraction procedure, we trained a model on measurements made with only a constant sky value subtracted, and another model with only the B-spline curve background estimate. After fine-tuning, the performance of these models was not significantly different from that of our final model, with an MAE of 0.0164 and 0.0158 compared to the final MAE of 0.0159. This indicates that the choice of background subtraction procedure in the artificial dataset does not have a strong effect on the model's performance.

Next, we locate and mask very small background galaxies that might otherwise fall below the surface brightness threshold and artificially inflate the ICL fraction. We do this by unsharp masking the original image to increase the contrast before segmenting the image to
create the “hot mask”, as done in Montes et al. (2021) and Martínez-Lombilla et al. (2023). First, we mask out the parts of the image that are brighter than the $\mu_r = 26$ mag/arcsec² threshold, as our goal here is to capture very small scale peaks or asymmetries fainter than this threshold that could be caused by background galaxies. The unsharp masking is then done by convolving the masked image with a Gaussian filter of $\sigma = 2$ pixels using ASTROPY (Astropy Collaboration et al. 2022), subtracting it from the original, and segmenting the resulting image using PHOTUTILS methods. This will identify faint, small objects in the image that we then mask for the ICL calculation.

We then estimate the $3\sigma$ surface brightness limits of the image. We calculate the standard deviation of the sky by calculating the $3\sigma$ clipped standard deviation in the same region around the edges of the image from which we derived the background estimate, again to avoid including any extended light in the center of the image in our calculation. We then use this to find the surface brightness limit of the image as in Equation 2 following Román et al. (2020), where $\mu_{lim}(3\sigma_{10'' \times 10''})$ is the $3\sigma$ surface brightness limit for an angular scale of $10 \times 10$ arcsecond², $\sigma$ is our calculated standard deviation of the sky, pix is the pixel scale, and ZP is the zero point of the image given by the HSC-SSP catalogue.

$$ \mu_{\lim}(3\sigma_{10'' \times 10''}) = -2.5 \times \log \left( \frac{3\sigma}{\text{pix} \times 10} \right) + ZP \quad (2) $$

Any pixels that fall below this calculated surface brightness limit are masked out to minimize background contamination in our calculation.

The ICL fraction is the flux from the ICL divided by the total flux from all galaxies and ICL in the cluster. However, the images from our artificial dataset are centered on LRGs rather than BCGs. So, for the artificial dataset measurements, we simplistically assume that all unmasked objects within the field of view belong to the cluster (which does not include bright stars or the small faint background galaxies as these have already been automatically masked). This is likely to lead to an overestimate of the cluster flux and a consequent underestimate of the ICL fraction, however, we leave the correction of this potential bias to the fine-tuning process. To investigate the effect of this assumption, we made a set of test measurements masking all galaxies other than the central LRG in the image, that is, assuming that the central LRG is the only cluster member, and a second set of test measurements where we classify galaxies with a photometric redshift within $3\sigma$ of the central LRG as member galaxies (the same method that we use for our real BCGs, detailed below). The first set of test measurements

results in worse performance after fine-tuning (an MAE of 0.0242 compared to 0.0159), because the model is trained to ignore satellite galaxies entirely during the initial training phase and this is too large of a difference for it to overcome during fine-tuning. On the other hand, the second test results in similar performance after fine-tuning (MAE of 0.0185 compared to our final MAE of 0.0159).

We choose a surface brightness threshold of $\mu_r = 26$ mag/arcsec², k-corrected and corrected for cosmological surface brightness dimming as for the data generation.

We apply the HSC-SSP provided bright star mask and the hot mask described above to the image. The ICL flux is calculated by summing the unmasked pixels below the chosen surface brightness threshold within a circular 300kpc radius aperture, covering the full radius of the image, which as previously described, should be large enough to capture the measurable extent of the ICL above the background in most cases. The total flux of the cluster is the sum of all unmasked pixels within the aperture. The ratio of these is then the final ICL fraction, which we use to label each of our artificial dataset images. For the artificial dataset, the measured ICL fractions range from 0.001 to 0.261, with a mean of $0.075 \pm 0.034$. This distribution arises naturally from our data generation process which prioritizes the realistic appearance of the generated ICL profiles.

### 2.3.2. Measurement method for fine-tuning dataset

The measurement method for the fine-tuning dataset is largely similar to the artificial dataset, but at each step we visually inspect the image and make manual adjustments as necessary. We note that although these manual adjustments are necessary for the manual measurement of the fine-tuning dataset to ensure accurate labeling of the images, the images provided to the machine learning model during training and fine-tuning are not adjusted manually, i.e. they are not sky subtracted and only apply the provided HSC-SSP bright star masks. The model thus naturally learns to ignore the background of the images when making its predictions, meaning that when the model is later applied to samples of images, background subtraction and manual masking of bright objects is not necessary. The images we use for the manual measurements are not resized to $224 \times 224$ pixels before measurement.

Firstly, we visually inspect the image to find any obvious problems with the bright object masks provided by HSC-SSP. At times, these masks are not large enough to fully cover bright stars, or there are stars that have not been masked out, so we manually create and apply new masks to the images in these cases.
We perform the background estimate as for the artificial dataset. After visual inspection of the result and manually checking the radial surface brightness profile of the image after subtraction, we then make manual adjustments to the background subtraction process as needed, such as adjusting the estimated gradient over the image or subtracting additional constant sky.

To calculate the total cluster light for this sample, we use the photometric redshifts for other objects in the cutout to construct a mask that includes only the galaxies belonging to that cluster. We classify a galaxy as a cluster member if its photometric redshift is within $3\sigma$ of the cluster's, where the $1\sigma$ accuracy of the photometric redshifts is $0.05(1+z_{phot})$ (Tanaka et al. 2018). After identifying non-cluster members, we use the PHOTUTILS segmentation map to create an initial mask covering these galaxies. We enlarge these non-cluster member masks by convolving them with a Gaussian kernel of size 0.3 pixels after finding this size kernel worked well for covering the faint outskirts of these sources at our image resolution. We apply this mask to the image, so that the non-cluster member flux is not included in either the estimation of the ICL or total cluster light.

We then construct the hot mask and mask below the surface brightness limits of the image as described previously, and finally calculate the ICL fraction as for the artificial dataset. The observational uncertainty on the ICL fraction is often estimated from the sky variance. However, the impact of background subtraction and masking on ICL fraction are significantly greater than the impact of sky variance in the images (mean = 0.001). We therefore use these steps to estimate the uncertainties on the fractions by measuring the images multiple times, spanning the range of parameters used in the manual measurements. We then take the standard deviation in the resulting ICL fraction as the uncertainty on the manually measured fraction.

## 2.4. Model architecture

We build a model that accepts a single band $224 \times 224$ pixel image as input, and produces a prediction for the ICL fraction for that image as output. A primary goal of this method is to remove the need for manual involvement. For this reason, the images given to the model are exactly as downloaded from HSC-SSP, with the only processing done being to downsize the images to the required resolution and applying the bright star masks as provided by HSC-SSP. No background subtraction or additional masking is performed.

The images are first standardized by arcsinh stretching the images, dividing them by the mean absolute deviation of the images, and then normalizing them between 0 and 1.

This suppresses the highest and lowest intensities of the image, and compresses the range to make it easier for the neural network to handle the inputs.

The architecture of the network consists of an augmenter, an encoder, and output layers to generate the final answer.

The augmenter introduces small random perturbations to input images that have no effect on the images' labels, and is key to showing the model what features of the image are unimportant to the task and preventing the model from overfitting to the training dataset. We randomly apply three image augmentations:

* Randomly flipping the image horizontally and/or vertically, each with 50% probability.

* Adding random Gaussian noise to the image with $\sigma$ equal to the mean absolute deviation of the unstandardized images.

* Randomly rotating the image by 0 (no rotation), 90, 180, or 270 degrees, each with equal probability.

We base our encoder on a ResNet-50 architecture (He et al. 2016) followed by a global average pooling layer to produce a 256-dimensional vector representation of the input image.

The output layers then use this vector to produce the final ICL prediction. We represent this answer as a probability density distribution over the possible ICL fractions. This is done by first passing the 256-dimensional representation through a dense hidden layer with 2048 neurons, and then another dense layer which generates 48 outputs. These outputs are used as parameters to the output distribution, which is a mixture of beta distributions. We choose beta distributions as they are bounded between 0 and 1, which is a useful constraint for this problem. The mode of this distribution will then be the most probable value of the ICL fraction for the corresponding input image, and we retain information about the model's confidence in this prediction through the full probability distribution.

## 2.5. Model training and fine-tuning

The model is built and trained using Tensorflow (Tensortflow Developers 2023). It aims to maximize the log likelihood of the correct ICL fraction under the output probability distribution. This is implemented as a negative log likelihood loss function using the Tensorflow Probability library (Dillon et al. 2017).

The model is compiled with the Adam optimizer (Kingma & Ba 2014) and is first trained on the 50,000
images in the artificial dataset using 90% (45,000) for training and 10% (5,000) for validation. The learning rate starts at $10^{-4}$ and is decreased using a learning rate scheduler to $10^{-5}$ when the validation loss has not decreased for 10 epochs. The training is then stopped once the validation loss has not decreased for 25 epochs. Our final model was trained for 180 epochs with these conditions, taking 21 hours and 45 minutes on a single Tesla V100-SXM2 GPU.

We then need to transfer the model's learning onto real data using the fine-tuning data. A conventional 80-20 split for training and testing here will result in a testing sample too small to be able to accurately judge the final performance of the model, given that our full sample only has 101 cluster images. We therefore perform 5-fold cross-validation by repeating the fine-tuning process five times, using different partitions of the full dataset for testing each time. The model is trained from the same initial checkpoint each time. This firstly allows us to verify that the performance of our fine-tuned model is stable across different sets of training and testing data, and secondly, by collating the results of all five fine-tuned models, allows us to evaluate the effect of fine-tuning and performance of the final models over a statistically significant sample of all 101 real clusters.

For the fine-tuning stage, we freeze all but the final two dense layers in the network, meaning that only the parameters in these final layers are allowed to change during fine-tuning, while the more fundamental weights of the encoder are kept fixed. This is because the model should only need small adjustments to transfer its learning, and allowing fewer parameters to vary will reduce the risk of the model overfitting to the small fine-tuning dataset. The model is then trained for 100 epochs at a learning rate of $10^{-6}$, which takes 3 minutes.

## 3. RESULTS

In this section we first show the model's performance on the artificial dataset after the initial training stage, showing that it has successfully learned the ICL fraction prediction task for this first dataset. We then compare the model's performance on the real data before and after the fine-tuning step to explore the effect of fine-tuning.

### 3.1. Artificial dataset

Figure 4 shows a plot of the model's final performance on the validation split (5,000 images) of the artificial dataset after the initial training was completed. The top plot shows the distribution of ICL fractions over the validation split. The bottom plot shows the ICL fraction value predicted by the model as a function of

<figure><img src="image_6.png" /><figcaption>Figure 4. Top plot: Distribution of the measured ICL fractions in the validation data labels. Bottom plot: Performance of the model on the artificial dataset, plotting predicted values against actual values. Points have been binned into 20 equally spaced bins, with the solid line showing the median, and the solid colors marking regions within which 75% (dark blue) and 90% (lighter blue) of the points fall. The black dashed line shows the ideal one to one relationship.</figcaption></figure>

the actual ICL fraction value, with the values binned into 20 equally spaced bins for clarity, and the scatter in each bin represented by the filled regions. A perfect model would predict every value correctly, and would be represented as a perfect one to one relationship on the plot, which is marked as a black dashed line. We see that in general, the model predicts very close to the one to one relationship with only a small amount of scatter. The mean absolute error (MAE) between the model's predictions and the actual fraction is 0.00509.

The mean predictions of the model as shown in Figure 4 appear to worsen at the higher range of ICL fractions, with the MAE in the second highest and highest bins being 0.0116 and 0.0447 respectively. However, given the small number of points in these bins (2 in each bin), it is not clear how representative these values are of the true ability of the model in this range. There are very few clusters in our sample that have an ICL fraction greater than 0.2 - only 129 clusters from our full artificial dataset are in this range (0.258% of the dataset). Although, as seen in Figure 4, the model does qualitatively appear to have been able to learn to distinguish these from lower ICL fraction images, it is clear from the MAE in the final bins that this range of fractions is
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Performance of the model on the real dataset before fine-tuning, plotting predicted values against actual values. Grey points show all data points, where vertical error bars are the 1σ confidence intervals calculated from the output probability distributions, and horizontal error bars show the uncertainty in the manual measurements. Red points are median values in five bins with equal numbers of points, with horizontal error bars showing the extent of each bin.</figcaption></figure>

less well characterized by the model than clusters with ICL fractions below 0.2 due to the lack of training samples. This could cause it to predict with more scatter, as well as preferentially predict slightly closer to the mean to minimize the loss penalty of mistakenly predicting a high ICL fraction.

Despite these small inaccuracies, we can see that the model has certainly learned how to distinguish high ICL images from low ICL images for this dataset, and generally achieves a high accuracy in its predictions.

## 3.2. Real data

Figure 5 shows the model's performance on the real data, before any fine-tuning has taken place. We also bin the data into 5 bins with equal numbers of points to better see the general trend in the model's predictions. There is a clear linear upwards trend in the model's predictions, showing that the model has been able to apply its learning to this new data domain to an extent, but the trend is flatter than the ideal one to one relationship shown as the black dashed line. The model also appears to be overconfident in its predictions, shown by the extremely small vertical error bars, which represent the 1σ confidence intervals calculated from the output probability distributions. The MAE in the model's predictions is 0.0295. Fine-tuning on some amount of real data should help the model adjust to the differences in domain, and reevaluate its confidence in its predictions.

Figure 6 shows the model's performance after fine-tuning on the real dataset. We use an 80-20 split for

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Performance of the model on the real dataset after fine-tuning, plotting predicted values against actual values. Top plot: Results from each fine-tuning run indicated by different colors, with different samples used for testing each time. Vertical error bars are the 1σ confidence intervals calculated from the output probability distributions, and horizontal error bars show the uncertainty in the manual measurements. Bottom plot: All data points are shown in grey. Red points showing median values in five bins with equal numbers of points, with horizontal error bars showing the extent of each bin.</figcaption></figure>

training and testing, however, due to the small size of the real dataset, this equates to only 20 points used for testing. We therefore run the fine-tuning process on 5 different splits of the data, with different points used for testing each time, in order to collate the results and evaluate the effect of fine-tuning on the model's performance in a statistically significant way. Each model is trained from the same initial point, using weights saved from the initial training on the artificial dataset. The top panel in Figure 6 shows which predictions come from which split in the data, whereas the bottom panel shows all the data points, binned into 5 bins, each with equal numbers of points.

The top panel shows that the splits are sufficiently random to give the model a good range of data points each time, and also shows that the fine-tuning of the model is stable, since each split shows a similar testing
performance despite the different partitions of training and testing data used.

We can see from the binned trend in the bottom panel that fine-tuning has successfully corrected the bias at higher ICL fractions that we observed in the pre-fine-tuning performance in Figure 5, with the binned points now falling very close to the one to one relationship. The MAE for the model now is 0.0159. The model has also reevaluated its confidence intervals after being trained on this new data. These now better reflect the actual scatter in the model's predictions. On average, the model now has excellent agreement with the fractions measured manually with the surface brightness method, however, we do observe some scatter in this relationship. In Figure 7 we show the extent of this scatter for our dataset. Generally it is of the order of a few percent, with a mean of -0.005 and rms of 0.025, but it does reach a maximum of 0.162 in the worst case.

This maximum is due to a significant outlier in the measured fraction, evident in Figure 6, which has a measured ICL fraction of 0.338±0.003. This is 0.134 higher than the next highest fraction. This outlier poses an issue for the model, which underpredicts the measured fraction significantly (prediction $0.176^{+0.031}_{-0.024}$). One explanation for this is that the model characterizes this cluster poorly due to it being a very large outlier. When this cluster is in the test split during fine-tuning (where the prediction shown in Figure 6 comes from), the model will never see any cluster with a fraction this high, and is therefore not able to recognize it correctly during testing. In actual use of the model, this would be mitigated since this cluster would form part of the fine-tuning set. The other possibility is that some element of this specific image poses a challenge for the model. We investigated whether the presence of a significant foreground galaxy is causing the model to mistakenly include it as part of the cluster, resulting in the underestimated ICL fraction. However, masking the galaxy in question did not result in a significant difference in the model's prediction (now predicting $0.196^{+0.032}_{-0.024}$, within the uncertainties of the original prediction). Without other examples of very high ICL fraction clusters, it is unclear whether this is an issue with this particular cluster, or a systematic issue with all very high ICL fraction clusters. In any case, we must always be cautious about trusting the particular prediction for one single cluster by the model, particularly for higher fractions, until the fine-tuning sample is sufficiently expanded to provide more training with high ICL fraction clusters.

# 4. DISCUSSION

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Scatter in ICL fraction predictions of the model with respect to actual fractions in the fine-tuning dataset. The distribution has a mean of -0.005 and rms of 0.025.</figcaption></figure>

We have shown that our machine learning model is capable of accurately predicting ICL fractions in cluster images, where the provided images require no manual preprocessing (i.e. the images given to the model are exactly as downloaded from HSC-SSP, with automatic application of the bright star masks provided) and the model is not provided any other external information. Even with only a small amount of available data and measurements (just 101 real cluster images), through the use of transfer learning the model was able to effectively learn the task and adapt its learning from an artificial domain to the real domain. However, there are concerns that come with any machine learning model that can influence how useful these models are in real science applications. In this section, we discuss how our model may be affected by these issues.

## 4.1. *Manual measurements*

A machine learning model is only as good as the data it is trained on, and without any real “ground truth” ICL dataset available for training, and with the wide variety of different measurement methods available, it is important to consider the systematics involved in the method we chose to use for our data labeling.

In this work we used the surface brightness cut method. This is one of the simplest methods to implement, making it relatively easy to automate for our large artificial dataset. It also has the advantage that it makes no assumptions about the morphology of either the BCG or the ICL. The clear drawback to this method is the relatively arbitrary choice of surface brightness threshold, and the fact that this method cannot measure the part of the ICL that falls over the BCG or any other galaxies.
Brough et al. (2024) made a comparison of the performance of eight different ICL measurement methods applied to 61 mock images of galaxy clusters from simulations. They found that, within the uncertainties, the mean ICL fraction measured by the eight different methods on the same cluster images were consistent, ranging from $0.09 \pm 0.02$ to $0.17 \pm 0.08$, with a mean over all methods of $0.13 \pm 0.05$. The two surface brightness threshold methods found a mean ICL fraction of $0.14 \pm 0.04$ (‘Martinez-Lombilla’) and $0.14 \pm 0.03$ (‘Montes’). Our measurement method, applied to our fine-tuning dataset, finds a mean ICL fraction of $0.075 \pm 0.055$, and for our artificial dataset the mean is $0.075 \pm 0.034$. These values are lower than the mean found by Brough et al. (2024). However, this is not surprising when we consider the differences between the datasets. Brough et al. (2024) are measuring the ICL fraction in mock images from simulations rather than real observations, and make their measurements within a 1 Mpc aperture, significantly larger than ours. They also simulated a limiting surface brightness of $\mu_r = 30.3 \text{ mag/arcsec}^2$ to match the estimated LSST 10 year limiting surface brightness. This limit is deeper than our images, which reach a maximum of $\mu_r = 29.8 \text{ mag/arcsec}^2$ for the Ultradeep layer of HSC (Martínez-Lombilla et al. 2023).

The only other work measuring the ICL fraction in HSC-SSP images is Furnell et al. (2021). They measured the ICL in 18 clusters detected from the XMM Cluster Survey data that fall in the deep HSC-SSP field. Five of these clusters also appear in our final fine-tuning sample. They also use the surface brightness cut method, however, their method deviates in other ways. Most significantly, they use a different surface brightness threshold to perform their measurement ($\mu_B = 25 \text{ mag/arcsec}^2$). Assuming an age of 2 Gyr and [Fe/H] of -0.4 for the ICL (e.g. Montes & Trujillo 2018) giving a color of $B-r = 0.74$ (Vazdekis et al. 2016), this corresponds to a surface brightness threshold of $24.26 \text{ mag/arcsec}^2$ in the $r$ band, significantly brighter than our threshold of $26 \text{ mag/arcsec}^2$. A brighter threshold will lead to more light from the brighter parts of the cluster being included in the ICL, naturally leading to higher fractions. They also make their measurements within varying radii, depending on the $R_{500}$ of their clusters as estimated from X-ray temperatures, with a mean over all clusters of $\sim 600 \text{ kpc}$, on average significantly larger than our chosen radius of 300 kpc. As expected, Furnell et al. (2021) find a mean fraction of $0.243 \pm 0.002$, much higher than ours. The differences in our measurement methods make it difficult to draw conclusions from these comparisons, despite having some of the same clusters in our samples.

<figure><img src="image_10.png" /><figcaption>**Figure 8.** ICL fraction as a function of redshift. Blue and purple circles represent our manually measured and model predicted fractions, respectively, and the black dashed line shows the line of best fit for our manual measurements. Teal squares are the fractions from Montes & Trujillo (2018), and yellow diamonds are the fractions from Furnell et al. (2021). Grey lines connect fractions that correspond to the same cluster.</figcaption></figure>

Martínez-Lombilla et al. (2023) also use HSC-SSP images, although they measured the diffuse light within a galaxy group (intragroup light, or IGL), rather than a cluster. They performed their measurements in $g$, $r$, and $i$ bands, and used both a composite model method and a surface brightness threshold method at two different thresholds. The method most comparable to ours (surface brightness threshold of $26 \text{ mag/arcsec}^2$ in the $r$ band, within a radius of $\sim 275 \text{ kpc}$), found an IGL fraction of $0.086 \pm 0.024$, which agrees well with the mean ICL fraction that we find in our sample. However, the IGL fraction that they find using other methods and other bands varies considerably in their studied group, from $0.016 \pm 0.043$ to $0.365 \pm 0.022$. This again demonstrates the importance of building a homogeneous sample of ICL measurements using a consistent measurement method, as even when using the same data, differences in measurement method can make it very difficult to compare between studies.

We can also consider any dependence on redshift that we observe in our sample. In Figure 8 we present the fractions we measure for our clusters, as well as the fractions predicted by the model. We also plot measurements from Furnell et al. (2021), and Montes & Trujillo (2018), as their clusters are within a similar redshift range to ours. Grey lines connect different measurements for the same cluster.

We find a negative relationship between our measured ICL fractions and redshift. The line of best fit for the manual measurements is shown in Figure 8 by the black dashed line and is given by $f_{\text{ICL}} = (-0.33 \pm 0.04)z +$
(0.18 ± 0.01). The relationship given by the model predictions is $f_{\text{ICL}} = (-0.30 ± 0.04)z + (0.17 ± 0.01)$. This relationship is consistent within uncertainties between the manual measurements and the model predictions, showing that the model does not show a particular redshift-dependent bias in its performance. The strength of this negative relationship is similar to that found by Montes & Trujillo (2018), although our fractions are on average slightly higher than theirs. They also use a surface brightness threshold method, using a threshold of $\mu_V = 26 \text{ mag/arcsec}^2$ and images that are 600 kpc on each side. Furnell et al. (2021) measure fractions much higher than ours on average, and also find a stronger correlation with redshift. As discussed previously, this is likely due to several important differences between our method and theirs.

There are a number of assumptions that we made in our measurements, particularly of the artificial dataset, that could potentially affect our model's accuracy. We noted that our use of smooth ICL profiles in the dataset generation could impact our model's ability to identify ICL with more asymmetric morphologies, however we do not find this to be a particular issue after the fine-tuning of the model on real data that does include examples of more asymmetric ICL profiles. In our artificial dataset measurements, we were not able to use cluster membership information, and instead assumed that all galaxies in the cutout belonged to the cluster. This would generally lead to an overestimate of cluster flux and therefore could lead to underestimates in ICL fraction, and we did in fact observe a tendency of the model to underestimate ICL fraction before fine-tuning, as shown in Figure 5. However, fine-tuning appears to have been able to correct this issue as shown in Figure 6. Our automatic background subtraction method as used for the artificial dataset was not always perfect, and we did observe some over or under-subtraction that required adjustment in the manual measurements. However, the tests we present in Section 2 and the model's ability to easily learn to match the manual measurements with more accurate background subtraction shows that it is not overly sensitive to the exact background subtraction method. We do not find these potential issues with the data generation and automatic measurements to have lasting effects on the model's accuracy as shown by its higher performance after the fine-tuning stage.

There are other assumptions made while measuring the fine-tuning dataset that could affect the uncertainties in our manual measurements. Although we do include cluster membership information in our manual measurements, these are based on cuts in photometric redshift, a common way to estimate cluster members

in the absence of spectroscopic information (e.g. Burke et al. 2015; Montes & Trujillo 2018; Furnell et al. 2021), but which could lead to misidentifications.

Our background subtraction is also a source of measurement uncertainty, even though in the case of the manual measurements we visually inspect the result for each cluster. We calculate our measurement uncertainty based on the sensitivity of the ICL fraction to background subtraction parameter choice, and find that this aspect should not influence our measurements significantly.

Downsizing the images to 224×224 pixels for the model input may also be impacting its predictions as compared to our measurements, as images that are lower redshift are more heavily binned due to their larger original size than images at higher redshifts. However, Brough et al. (2024) found that binning images by different amounts before measurement did not have a significant impact on the ICL fractions found, and we do not observe sensitivity in the model to this as shown by the lack of dependence of the model's performance on redshift in Figure 8.

Finally, our choice of radius is another important assumption in our method. We chose a radius of 300 kpc, owing to the fact that this generally corresponds to the ICL detection limit (e.g. DeMaio et al. 2015; Jiménez-Teja et al. 2018; Montes et al. 2021) and the fact that we cannot estimate physical cluster sizes for all of our clusters with our available data. However, it is possible that using a fixed physical radius will have an impact on our measured values, particularly at low redshifts where the ICL is more extended on the sky. To investigate the impact of our chosen radius limit on our measurements, we manually measure ICL fractions for a sample of clusters out to the $R_{500}$ scale radius (the radius containing a volume with a mean density 500 times the critical density at that redshift) and compare the fractions to our measurements within 300 kpc. To estimate $R_{500}$ we crossmatch our Deep/Ultraleep CAMIRA clusters with the Galaxy and Mass Assembly (GAMA) Galaxy Group Catalogue (Robotham et al. 2011) and find 13 matches with $\log(M_{200})$ ranging from 13.7 $M_{\odot}$ to 15.0 $M_{\odot}$. We use a scaling relation to calculate $R_{500}$ from the measured velocity dispersions of the GAMA groups (Zhang et al. 2011) Over this subsample, the mean $R_{500}$ is 640 ± 220 kpc.

We find that the ICL fractions measured within $R_{500}$ increase by a mean of 0.026 and up to a maximum of 0.073. Clusters with larger $R_{500}$ have larger differences in the measured fraction compared to the measurements taken within a radius of 300 kpc, however the differences appear to be independent of cluster redshift. This indi-
cates that there is some ICL that is missed when measuring within a radius of 300 kpc as compared to $R_{500}$. However, we note that there is currently no consensus on the best radius to use for measuring the ICL fraction, and there are a number of approaches taken in the literature ranging from a consistent physical radius, to cluster-dependent radius measurements. We have chosen to use a physical radius of 300 kpc as it is easy to maintain consistency amongst a photometrically selected cluster sample where scale radius estimates are not available and balances with the machine learning requirements, where significantly larger radii (e.g. 1 Mpc) pose issues with the spatial resolution of pixels after they have been resized to the 224 x 224 pixels necessary for the model input. However, we note that our choice of a 300 kpc radius is a free parameter that could complicate comparison with measurements that are made with a different radius. Our focus is to enable the construction of a large homogeneous sample, so it is more important to ensure that the choice of radius can be easily applied to all of our cluster measurements.

There are significant difficulties with comparing manual ICL fraction measurements between different samples and different methods, as anything from differing surface brightness limits to choice in observed wavelength, surface brightness threshold, or radius can have a significant impact on the measured fraction. Our fractions are consistent with previous comparable studies, and although some assumptions we have made introduce some uncertainties that may be adopted by the model, we do not find these to contribute significantly to uncertainties in the model's predictions.

## 4.2. Interpretability of machine learning outputs

An issue in the application of machine learning models is the trustworthiness of the model's outputs. Models are often treated as “black boxes”, meaning that the user sometimes has no idea why the model has produced a particular answer for a particular input.

In particular for regression models, many models have no way of indicating confidence in their answer, meaning that there is no way to distinguish between trustworthy and untrustworthy outputs. We address this issue by using a probabilistic output. In the model's training, it aims to maximize the likelihood of the correct answer, and is penalized based on how unlikely the correct answer was under the predicted distribution. This means that if the model is unsure, it is incentivized to predict a more flat probability distribution to minimize the penalty. Since we retain this probability distribution, we therefore have information about how confident the

model is in its prediction. This is demonstrated in the error bars given in Figures 5 and 6.

Another concern about machine learning models is their ability to find and exploit shortcuts in data that they have been given, so although it can appear from the test results that the model has learned the task, in fact it may be exploiting peculiarities in the dataset to “cheat” the answer and not performing the task as intended. This is a well-known problem across many domains of machine learning, sometimes known as shortcut learning (see Geirhos et al. 2020 for a review). Given that our model is easily able to transfer its learning from one dataset to another (the artificial dataset to the real images), this is unlikely to be the case, and indicates that the model has in fact learned the real task of identifying the ICL and not overfit to some underlying pattern in the artificial dataset. However, to check this assumption, we can use saliency maps, which highlight pixels in input images that the model deems to be important when making its decision. These can be very useful to check that a model is not exploiting shortcuts in the data (e.g. Zech et al. 2018).

The method we use to generate these saliency maps is GradCAM (Selvaraju et al. 2020). GradCAM uses the magnitude of gradients flowing into the last convolutional layer in the network (the last layer to retain spatial information from the image) to assign importance to different features in the network. An activation map can then be generated over the original image that combines all feature map activations weighted by importance. This can be used to check that the model is focusing on areas of the image that make sense.

Figure 9 shows three GradCAM outputs for images from the fine-tuning dataset, where we ask it to highlight areas of the image that are important for predicting a high ICL fraction. Note that this is not the same as asking the model to highlight areas where ICL exists, rather, the model shows us what areas in the image it is using to distinguish between it having a low ICL fraction and it having a high ICL fraction. For the sake of comparison, we also show in Figure 9 a thresholded version of each input image to show where the low surface brightness light is actually located in the image. A shortcoming of GradCAM is that it does not tell us exactly how the model is using these areas in making its decision, only the relative importance. This means that as a tool, it is most suited for checking that the model is using areas of the image that indicate it is performing the task as intended.

In Figure 9 we see that the model is focusing primarily on the low surface brightness areas of the image while largely ignoring the bright, central parts of the galaxy,
<figure><img src="image_11.png" /><figcaption>**Figure 9.** GradCAM outputs on three example cluster images from the fine-tuning dataset. Top row: an example of a typical activation map. Middle row: an example of the model activating on a tidal tail and shell. Bottom row: an example of the model ignoring a background gradient due to a poorly masked star. Left column: Original images from the fine-tuning dataset. Middle column: Input images, masked and thresholded to bring out the low surface brightness light. Right column: Activation maps overlaid on the original images, where red indicates stronger activation.</figcaption></figure>

showing that it is considering relevant parts of the image to make an informed conclusion. Its attention to these areas correctly applied to each image makes it unlikely to be exploiting peculiarities in the dataset. There are some points of interest in these activation maps. Firstly, we see that there is a tendency to activate along the borders of the images. It is possible that these indicate the model checking for background gradients to ignore – in fact, we see in the example in the bottom row that the model is capable of disregarding the background gradient that is present in that image due to a large star that is not fully masked in the top corner. There is also the possibility that this is an artifact due to the model's architecture. Alsallakh et al. (2020) found similar artifacts in their images to be connected to padding in the model's convolutional layers, which has also been connected to a model's ability to learn spatial information (Semih Kayhan & van Gemert 2020). With only these saliency maps, the exact reason for the activation along the borders remains unclear. We also see that, interestingly, the model shows interest in tidal features – in

the middle example, the central galaxy has a clear tidal tail and shells, which correspond to bright spots in the saliency maps.

This type of visual interpretability is still quite new in machine learning, and as the tools improve we will be able to have greater insight into how models make their decisions. GradCAM allows us to see that the model is not exploiting shortcuts in the dataset to “cheat” the task, and really has learned some underlying rules in the images to using the low surface brightness light to predict the ICL fraction. This indicates that its high performance on the test data can be trusted.

## 5. CONCLUSION

The problem of systematics introduced by choice of measurement method is still an unsolved one in the field of ICL measurements, which is partly what motivates the construction of a large, homogeneous sample with consistent assumptions and systematic uncertainties. Here we present a model trained using one reliable and commonly used measurement method, showing that machine learning models are capable of learning to perform a task on par with manual measurement. As understanding of different ICL measurement methods evolves, the data that the model is trained on can be developed further. In the future, it will be possible to fine-tune the model on larger samples made with different methods or even combinations of methods, in order to compare the quality of the predictions made.

The main goal of developing a machine learning method such as this is to enable fast processing of vast amounts of data that more traditional measurement methods are unable to cope with, which will be available in new surveys such as LSST as well as Euclid (e.g. Kluge et al. 2024). From our results, we can see that our model is capable of replicating measurements made with traditional methods without any of the manual inspection and preprocessing that would normally be necessary.

A distinct advantage to the throughput of the model is its ability to process a batch of input images at once, in parallel. The size of the batch that can be processed at once is limited by the available memory. On the single GPU that we use for training the model, which has 31GB of memory, our model is able to process a maximum of 500 samples at once in a matter of seconds.

Of course, the model is not always perfect, and is not always able to perfectly replicate the manual measurements. This scatter is relatively small, but given that the model in itself is not physically motivated, for more detailed study and measurement of individual clusters, we will still need traditional measurement meth-
ods to get a complete picture of the ICL. More detailed
measurements using traditional methods will also pro-
vide the better understanding of the systematics and
assumptions involved in the different methods needed
to advance this field. With this in place we will also be
able to continue to train this model to make its predic-
tions more robust and accurate. Currently, the model
is trained on just 101 real cluster images and measure-
ments, and already achieves high performance, so even
a small amount of extra data could help to improve its
predictions. For now, the model's best use is in quickly
and automatically processing a very large number of im-
age samples in order to study large-scale trends in the
ICL fraction, and identify any particularly interesting,
outlying systems for detailed follow up studies.

In this work, we have shown that a machine learning
model can accurately predict the ICL fraction in clus-
ter images given only a small number of real training
images by leveraging transfer learning. Its ability to
quickly process many images at once makes it a good
candidate for dealing with the massive amounts of data
that large imaging surveys such as LSST and Euclid will
provide us. More data and development of the various
manual measurement methods used in the ICL field will
be used to refine the fine-tuning of the model and con-
tinue to improve its performance. All of the code used
for creating, training, and testing the model is made
available to anyone interested in using it on GitHub at
https://github.com/lpcan/MICL (Canepa et al. 2024).

ACKNOWLEDGEMENTS

We thank the anonymous referee for their thoughtful
report that has improved the paper.

SB acknowledges funding support from the Aus-
tralian Research Council through a Discovery Project
DP190101943. MM acknowledges support from
the project RYC2022-036949-I financed by the MI-
CIU/AEI/10.13039/501100011033 and by FSE+. Parts
of this research were supported by the Australian Re-

search Council Centre of Excellence for All Sky Astro-
physics in 3 Dimensions (ASTRO 3D), through project
number CE170100013. Parts of this research were sup-
ported by the Astronomical Society of Australia (ASA),
through the Student Travel Assistance Scheme (STAS).

This research includes computations using the com-
putational cluster Katana supported by Research Tech-
nology Services at UNSW Sydney.

The Hyper Suprime-Cam (HSC) collaboration in-
cludes the astronomical communities of Japan and Tai-
wan, and Princeton University. The HSC instrumenta-
tion and software were developed by the National As-
tronomical Observatory of Japan (NAOJ), the Kavli
Institute for the Physics and Mathematics of the Uni-
verse (Kavli IPMU), the University of Tokyo, the High
Energy Accelerator Research Organization (KEK), the
Academia Sinica Institute for Astronomy and Astro-
physics in Taiwan (ASIAA), and Princeton University.
Funding was contributed by the FIRST program from
the Japanese Cabinet Office, the Ministry of Educa-
tion, Culture, Sports, Science and Technology (MEXT),
the Japan Society for the Promotion of Science (JSPS),
Japan Science and Technology Agency (JST), the Toray
Science Foundation, NAOJ, Kavli IPMU, KEK, ASIAA,
and Princeton University.

This paper makes use of software developed for Vera
C. Rubin Observatory. We thank the Rubin Observa-
tory for making their code available as free software at
http://pipelines.lsst.io/.

This paper is based on data collected at the Subaru
Telescope and retrieved from the HSC data archive sys-
tem, which is operated by the Subaru Telescope and
Astronomy Data Center (ADC) at NAOJ. Data analy-
sis was in part carried out with the cooperation of Cen-
ter for Computational Astrophysics (CfCA), NAOJ. We
are honored and grateful for the opportunity of observ-
ing the Universe from Maunakea, which has the cultural,
historical and natural significance in Hawaii.

REFERENCES

Adami, C., Durret, F., Guennou, L., & Rocha, C. D. 2013, Astronomy & Astrophysics, 551, A20, doi: 10.1051/0004-6361/201220282

Aihara, H., AlSayyad, Y., Ando, M., et al. 2019, Publications of the Astronomical Society of Japan, 71, 114, doi: 10.1093/pasj/psz103

—. 2022, Publications of the Astronomical Society of
Japan, 74, 247, doi: 10.1093/pasj/psab122

Alsallakh, B., Kokhlikyan, N., Miglani, V., Yuan, J., &
Reblitz-Richardson, O. 2020, arXiv preprint
arXiv:2010.02178

Astropy Collaboration, Price-Whelan, A. M., Lim, P. L., et al. 2022, The Astrophysical Journal, 935, 167, doi: 10.3847/1538-4357/ac7c74

Barfety, C., Valin, F.-A., Webb, T. M. A., et al. 2022, The Astrophysical Journal, 930, 25, doi: 10.3847/1538-4357/ac61dd
Bosch, J., Armstrong, R., Bickerton, S., et al. 2018, Publications of the Astronomical Society of Japan, 70, S5, doi: 10.1093/pasj/psx080

Bradley, L., Sipőcz, B., Robitaille, T., et al. 2023, astropy/photutils: 1.10.0, Zenodo, doi: 10.5281/zenodo.1035865

Brough, S., Proctor, R., Forbes, D. A., et al. 2007, Monthly Notices of the Royal Astronomical Society, 378, 1507, doi: 10.1111/j.1365-2966.2007.11900.x

Brough, S., Collins, C., Demarco, R., et al. 2020, The Vera Rubin Observatory Legacy Survey of Space and Time and the Low Surface Brightness Universe, arXiv. http://arxiv.org/abs/2001.11067

Brough, S., Ahad, S. L., Bahé, Y. M., et al. 2024, Monthly Notices of the Royal Astronomical Society, 528, 771, doi: 10.1093/mnras/stad3810

Burke, C., Hilton, M., & Collins, C. 2015, Monthly Notices of the Royal Astronomical Society, 449, 2353, doi: 10.1093/mnras/stv450

Canepa, L., Brough, S., Lanusse, F., Montes, M., & Hatch, N. 2024, Ipcan/MICL: Paper version, Zenodo, doi: 10.5281/ZENODO.14376238

Cañas, R., Lagos, C. d. P., Elahi, P. J., et al. 2020, Monthly Notices of the Royal Astronomical Society, 494, 4314, doi: 10.1093/mnras/staa1027

Chilingarian, I. V., Melchior, A.-L., & Zolotukhin, I. Y. 2010, Monthly Notices of the Royal Astronomical Society, 405, 1409, doi: 10.1111/j.1365-2966.2010.16506.x

Contini, E. 2021, Galaxies, 9, 60, doi: 10.3390/galaxies9030060

Contini, E., De Lucia, G., Villalobos, Á., & Borgani, S. 2014, Monthly Notices of the Royal Astronomical Society, 437, 3787, doi: 10.1093/mnras/stt2174

Contini, E., Jeon, S., Rhee, J., Han, S., & Yi, S. K. 2023, The Astrophysical Journal, 958, 72, doi: 10.3847/1538-4357/acfd25

Contreras-Santos, A., Knebe, A., Cui, W., et al. 2024, Astronomy and Astrophysics, 683, A59, doi: 10.1051/0004-6361/202348474

Da Rocha, C., & De Oliveira, C. M. 2005, Monthly Notices of the Royal Astronomical Society, 364, 1069, doi: 10.1111/j.1365-2966.2005.09641.x

Da Rocha, C., Ziegler, B. L., & Mendes de Oliveira, C. 2008, Monthly Notices of the Royal Astronomical Society, 388, 1433, doi: 10.1111/j.1365-2966.2008.13500.x

Dalal, R., Strauss, M. A., Sunayama, T., et al. 2021, Monthly Notices of the Royal Astronomical Society, 507, 4016, doi: 10.1093/mnras/stab2363

DeMaio, T., Gonzalez, A. H., Zabludoff, A., Zaritsky, D., & Bradač, M. 2015, Monthly Notices of the Royal Astronomical Society, 448, 1162, doi: 10.1093/mnras/stv033

Desmons, A., Brough, S., & Lanusse, F. 2024, Monthly Notices of the Royal Astronomical Society, 531, 4070, doi: 10.1093/mnras/stae1402

Dillon, J. V., Langmore, I., Tran, D., et al. 2017, TensorFlow Distributions, arXiv. http://arxiv.org/abs/1711.10604

Ellien, A., Slezak, E., Martinet, N., et al. 2021, Astronomy & Astrophysics, 649, A38, doi: 10.1051/0004-6361/202038419

Feldmeier, J. J., Mihos, J. C., Morrison, H. L., et al. 2004, The Astrophysical Journal, 609, 617, doi: 10.1086/421313

Furnell, K. E., Collins, C. A., Kelvin, L. S., et al. 2021, Monthly Notices of the Royal Astronomical Society, 502, 2419, doi: 10.1093/mnras/stab065

Geirhos, R., Jacobsen, J.-H., Michaelis, C., et al. 2020, Nature Machine Intelligence, 2, 665, doi: 10.1038/s42256-020-00257-z

Gladders, M. D., & Yee, H. K. C. 2000, The Astronomical Journal, 120, 2148, doi: 10.1086/301557

Gonzalez, A. H., Zabludoff, A. I., & Zaritsky, D. 2005, The Astrophysical Journal, 618, 195, doi: 10.1086/425896

Guennou, L., Adami, C., Rocha, C. D., et al. 2012, Astronomy & Astrophysics, 537, A64, doi: 10.1051/0004-6361/201117482

Han, Y., Zou, Z., Li, N., & Chen, Y. 2022, Research in Astronomy and Astrophysics, 22, 085006, doi: 10.1088/1674-4527/ac7386

Hayat, M. A., Stein, G., Harrington, P., Lukić, Z., & Mustafa, M. 2021, The Astrophysical Journal Letters, 911, L33, doi: 10.3847/2041-8213/abf2c7

He, K., Zhang, X., Ren, S., & Sun, J. 2016, in 2016 IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 770–778, doi: 10.1109/CVPR.2016.90

Huang, S., Leauthaud, A., Greene, J. E., et al. 2018, Monthly Notices of the Royal Astronomical Society, 475, 3348, doi: 10.1093/mnras/stx3200

Huang, S., Leauthaud, A., Hearin, A., et al. 2020, Monthly Notices of the Royal Astronomical Society, 492, 3685, doi: 10.1093/mnras/stz3314

Iodice, E., Spavone, M., Cattapan, A., et al. 2020, Astronomy & Astrophysics, 635, A3, doi: 10.1051/0004-6361/201936435

Ivezić, Ž., Kahn, S. M., Tyson, J. A., et al. 2019, The Astrophysical Journal, 873, 111, doi: 10.3847/1538-4357/ab042c
Jiménez-Teja, Y., & Dupke, R. 2016, The Astrophysical Journal, 820, 49, doi: 10.3847/0004-637X/820/1/49

Jiménez-Teja, Y., Dupke, R. A., Lopes, P. A. A., & Vílchez, J. M. 2023, Astronomy & Astrophysics, 676, A39, doi: 10.1051/0004-6361/202346580

Jiménez-Teja, Y., Dupke, R., Benítez, N., et al. 2018, The Astrophysical Journal, 857, 79, doi: 10.3847/1538-4357/aab70f

Kingma, D. P., & Ba, J. 2014, Adam: A Method for Stochastic Optimization, doi: 10.48550/arXiv.1412.6980

Kluge, M., Bender, R., Riffeser, A., et al. 2021, The Astrophysical Journal Supplement Series, 252, 27, doi: 10.3847/1538-4365/abcda6

Kluge, M., Hatch, N. A., Montes, M., et al. 2024, Euclid: Early Release Observations – The intracluster light and intracluster globular clusters of the Perseus cluster, arXiv. http://arxiv.org/abs/2405.13503

Krick, J. E., & Bernstein, R. A. 2007, The Astronomical Journal, 134, 466, doi: 10.1086/518787

Li, J., Huang, S., Leauthaud, A., et al. 2022, Monthly Notices of the Royal Astronomical Society, 515, 5335, doi: 10.1093/mnras/stac2121

Loh, Y.-S., & Strauss, M. A. 2006, Monthly Notices of the Royal Astronomical Society, 366, 373, doi: 10.1111/j.1365-2966.2005.09714.x

Marini, I., Borgani, S., Saro, A., et al. 2022, Monthly Notices of the Royal Astronomical Society, 514, 3082, doi: 10.1093/mnras/stac1558

Martínez-Lombilla, C., Brough, S., Montes, M., et al. 2023, Monthly Notices of the Royal Astronomical Society, 518, 1195, doi: 10.1093/mnras/stac3119

Mihos, J. C. 2019, Deep Imaging of Diffuse Light Around Galaxies and Clusters: Progress and Challenges, doi: 10.48550/arXiv.1909.09456

Mihos, J. C., Harding, P., Feldmeier, J., & Morrison, H. 2005, The Astrophysical Journal, 631, L41, doi: 10.1086/497030

Miyazaki, S., Komiyama, Y., Kawanomoto, S., et al. 2018, Publications of the Astronomical Society of Japan, 70, S1, doi: 10.1093/pasj/psx063

Montenegro-Taborda, D., Rodriguez-Gomez, V., Pillepich, A., et al. 2023, Monthly Notices of the Royal Astronomical Society, 521, 800, doi: 10.1093/mnras/stad586

Montes, M. 2019, The intracluster light and its role in galaxy evolution in clusters, arXiv. http://arxiv.org/abs/1912.01616

—. 2022, Nature Astronomy, 6, 308, doi: 10.1038/s41550-022-01616-z

Montes, M., Brough, S., Owers, M. S., & Santucci, G. 2021, The Astrophysical Journal, 910, 45, doi: 10.3847/1538-4357/abddb6

Montes, M., & Trujillo, I. 2014, ApJ, 794, 137, doi: 10.1088/0004-637X/794/2/137

Montes, M., & Trujillo, I. 2018, Monthly Notices of the Royal Astronomical Society, 474, 917, doi: 10.1093/mnras/stx2847

—. 2022, The Astrophysical Journal Letters, 940, L51, doi: 10.3847/2041-8213/ac98c5

Morishita, T., Abramson, L. E., Treu, T., et al. 2017, The Astrophysical Journal, 846, 139, doi: 10.3847/1538-4357/aa8403

Murante, G., Giovalli, M., Gerhard, O., et al. 2007, Monthly Notices of the Royal Astronomical Society, 377, 2, doi: 10.1111/j.1365-2966.2007.11568.x

Nishizawa, A. J., Hsieh, B.-C., Tanaka, M., & Takata, T. 2020, Photometric Redshifts for the Hyper Suprime-Cam Subaru Strategic Program Data Release 2, arXiv. http://arxiv.org/abs/2003.01511

Oguri, M. 2014, Monthly Notices of the Royal Astronomical Society, 444, 147, doi: 10.1093/mnras/stu1446

Oguri, M., Lin, Y.-T., Lin, S.-C., et al. 2018, Publications of the Astronomical Society of Japan, 70, S20, doi: 10.1093/pasj/psx042

Oquab, M., Bottou, L., Laptev, I., & Sivic, J. 2014, in 2014 IEEE Conference on Computer Vision and Pattern Recognition, 1717–1724, doi: 10.1109/CVPR.2014.222

Pan, S. J., & Yang, Q. 2010, IEEE Transactions on Knowledge and Data Engineering, 22, 1345, doi: 10.1109/TKDE.2009.191

Pearson, W. J., Wang, L., Trayford, J. W., Petrillo, C. E., & Van Der Tak, F. F. S. 2019, Astronomy & Astrophysics, 626, A49, doi: 10.1051/0004-6361/201935355

Presotto, V., Girardi, M., Nonino, M., et al. 2014, Astronomy & Astrophysics, 565, A126, doi: 10.1051/0004-6361/201323251

Proctor, K. L., Lagos, C. d. P., Ludlow, A. D., & Robotham, A. S. G. 2024, Monthly Notices of the Royal Astronomical Society, 527, 2624, doi: 10.1093/mnras/stad3341

Purcell, C. W., Bullock, J. S., & Zentner, A. R. 2007, The Astrophysical Journal, 666, 20, doi: 10.1086/519787

Ragusa, R., Spavone, M., Iodice, E., et al. 2021, Astronomy & Astrophysics, 651, A39, doi: 10.1051/0004-6361/202039921

Ragusa, R., Iodice, E., Spavone, M., et al. 2023, Astronomy & Astrophysics, 670, L20, doi: 10.1051/0004-6361/202245530
Robotham, A. S. G., Norberg, P., Driver, S. P., et al. 2011, Monthly Notices of the Royal Astronomical Society, 416, 2640, doi: 10.1111/j.1365-2966.2011.19217.x

Román, J., Trujillo, I., & Montes, M. 2020, Astronomy & Astrophysics, 644, A42, doi: 10.1051/0004-6361/201936111

Rudick, C. S., Mihos, J. C., & McBride, C. 2006, The Astrophysical Journal, 648, 936, doi: 10.1086/506176

Rudick, C. S., Mihos, J. C., & McBride, C. K. 2011, The Astrophysical Journal, 732, 48, doi: 10.1088/0004-637X/732/1/48

Selvaraju, R. R., Cogswell, M., Das, A., et al. 2020, International Journal of Computer Vision, 128, 336, doi: 10.1007/s11263-019-01228-7

Semih Kayhan, O., & van Gemert, J. C. 2020, in 2020 IEEE/CVF Conference on Computer Vision and Pattern Recognition (CVPR), 14262-14273, doi: 10.1109/CVPR42600.2020.01428

Spavone, M., Capaccioli, M., Napolitano, N. R., et al. 2017, Astronomy & Astrophysics, 603, A38, doi: 10.1051/0004-6361/201629111

Stein, G., Blaum, J., Harrington, P., Medan, T., & Lukić, Z. 2022, The Astrophysical Journal, 932, 107, doi: 10.3847/1538-4357/ac6d63

Tanaka, M., Coupon, J., Hsieh, B.-C., et al. 2018, Publications of the Astronomical Society of Japan, 70, S9, doi: 10.1093/pasj/psx077

Tang, L., Lin, W., Cui, W., et al. 2018, The Astrophysical Journal, 859, 85, doi: 10.3847/1538-4357/aabd78

Tensorflow Developers. 2023, TensorFlow, Zenodo, doi: 10.5281/zenodo.10126399

van der Walt, S., Schönberger, J. L., Nunez-Iglesias, J., et al. 2014, PeerJ, 2, e453, doi: 10.7717/peerj.453

Vazdekis, A., Koleva, M., Ricciardelli, E., Röck, B., & Falcón-Barroso, J. 2016, Monthly Notices of the Royal Astronomical Society, 463, 3409, doi: 10.1093/mnras/stw2231

Vojtekova, A., Lieu, M., Valtchanov, I., et al. 2021, Monthly Notices of the Royal Astronomical Society, 503, 3204, doi: 10.1093/mnras/staa3567

Walmsley, M., Scaife, A. M. M., Lintott, C., et al. 2022, Monthly Notices of the Royal Astronomical Society, 513, 1581, doi: 10.1093/mnras/stac525

Watkins, A. E., Kaviraj, S., Collins, C. C., et al. 2024, Monthly Notices of the Royal Astronomical Society, 528, 4289, doi: 10.1093/mnras/stae236

Zech, J. R., Badgeley, M. A., Liu, M., et al. 2018, PLOS Medicine, 15, e1002683, doi: 10.1371/journal.pmed.1002683

Zhang, Y.-Y., Andernach, H., Caretta, C. A., et al. 2011, Astronomy & Astrophysics, 526, A105, doi: 10.1051/0004-6361/201015830