# The Evolution of AGN Activity in Brightest Cluster Galaxies

T. Somboonpanyakul, M. Mcdonald, A. Noble, M. Aguena, S. Allam, A. Amon, F. Andrade-Oliveira, D. Bacon, M. B. Bayliss, E. Bertin, et al.

► To cite this version:

T. Somboonpanyakul, M. Mcdonald, A. Noble, M. Aguena, S. Allam, et al. The Evolution of AGN Activity in Brightest Cluster Galaxies. The Astronomical Journal, 2022, 163, 10.3847/1538-3881/ac5030. insu-03748284

HAL Id: insu-03748284

https://insu.hal.science/insu-03748284

Submitted on 9 Aug 2022

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_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# The Evolution of AGN Activity in Brightest Cluster Galaxies

T. Somboonpanyakul<sup>1,2</sup> ID, M. McDonald<sup>1</sup> ID, A. Noble<sup>3</sup>, M. Aguena<sup>4</sup> ID, S. Allam<sup>5</sup>, A. Amon<sup>2</sup>, F. Andrade-Oliveira<sup>4,6</sup>, D. Bacon<sup>7</sup>, M. B. Bayliss<sup>8</sup> ID, E. Bertin<sup>9,10</sup> ID, S. Bhargava<sup>11</sup>, D. Brooks<sup>12</sup> ID, E. Buckley-Geer<sup>5,13</sup>, D. L. Burke<sup>2,14</sup> ID, M. Calzadilla<sup>1</sup> ID, R. Canning<sup>15</sup>, A. Carnero Rosell<sup>4,16,17</sup> ID, M. Carrasco Kind<sup>18,19</sup> ID, J. Carretero<sup>20</sup> ID, M. Costanzi<sup>21,22,23</sup> ID, L. N. da Costa<sup>4,24</sup> ID, M. E. S. Pereira<sup>25</sup>, J. De Vicente<sup>26</sup> ID, P. Doel<sup>12</sup>, P. Eisenhardt<sup>27</sup>, S. Everett<sup>28</sup>, A. E. Evrard<sup>25,29</sup> ID, I. Ferrero<sup>30</sup> ID, B. Flaugher<sup>5</sup>, B. Floyd<sup>31</sup> ID, J. García-Bellido<sup>32</sup> ID, E. Gaztanaga<sup>33,34</sup> ID, D. W. Gerdes<sup>25,29</sup> ID, A. Gonzalez<sup>35</sup> ID, D. Gruen<sup>36</sup> ID, R. A. Gruendl<sup>18,19</sup> ID, J. Gschwend<sup>4,24</sup> ID, N. Gupta<sup>37</sup> ID, G. Gutierrez<sup>5</sup>, S. R. Hinton<sup>38</sup> ID, D. L. Hollowood<sup>28</sup> ID, K. Honscheid<sup>39,40</sup>, B. Hoyle<sup>36</sup>, D. J. James<sup>41</sup>, T. Jeltema<sup>28</sup> ID, G. Khullar<sup>13,42</sup> ID, K. J. Kim<sup>8</sup> ID, M. Klein<sup>36</sup>, K. Kuehn<sup>43,44</sup> ID, M. Lima<sup>4,45</sup>, M. A. G. Maia<sup>4,24</sup>, J. L. Marshall<sup>46</sup> ID, P. Martini<sup>39,40,47</sup> ID, P. Melchior<sup>48</sup> ID, F. Menanteau<sup>18,19</sup> ID, R. Miquel<sup>20,49</sup> ID, J. J. Mohr<sup>36,50</sup>, R. Morgan<sup>51</sup> ID, R. L. C. Ogando<sup>4,24</sup> ID, A. Palmese<sup>5,42</sup> ID, F. Paz-Chinchón<sup>18,19</sup> ID, A. Pieres<sup>4,24</sup> ID, A. A. Plazas Malagón<sup>48</sup> ID, K. Reil<sup>14</sup>, A. K. Romer<sup>11</sup> ID, F. Ruppin<sup>1</sup> ID, E. Sanchez<sup>26</sup> ID, A. Saro<sup>21,22,23,52</sup>, V. Scarpine<sup>5</sup>, M. Schubnell<sup>25</sup>, S. Serrano<sup>33,34</sup> ID, I. Sevilla-Noarbe<sup>26</sup> ID, P. Singh<sup>22,23</sup>, M. Smith<sup>53</sup> ID, M. Soares-Santos<sup>25</sup> ID, V. Strazzullo<sup>22,23,54</sup> ID, E. Suchyta<sup>55</sup>, M. E. C. Swanson<sup>18</sup> ID, G. Tarle<sup>25</sup> ID, C. To<sup>2,14,56</sup> ID, D. L. Tucker<sup>5</sup> ID, and R. D. Wilkinson<sup>11</sup>

<sup>1</sup> Kavli Institute for Astrophysics and Space Research, Massachusetts Institute of Technology, 77 Massachusetts Avenue, Cambridge, MA 02139, USA; taweewat@stanford.edu

<sup>2</sup> Kavli Institute for Particle Astrophysics & Cosmology, P.O. Box 2450, Stanford University, Stanford, CA 94305, USA

<sup>3</sup> School of Earth and Space Exploration, Arizona State University, Tempe, AZ 85287, USA

<sup>4</sup> Laboratório Interinstitucional de e-Astronomia—LIneA, Rua Gal. José Cristino 77, Rio de Janeiro, RJ—20921-400, Brazil

<sup>5</sup> Fermi National Accelerator Laboratory, P. O. Box 500, Batavia, IL 60510, USA

<sup>6</sup> Instituto de Física Teórica, Universidade Estadual Paulista, São Paulo, Brazil

<sup>7</sup> Institute of Cosmology and Gravitation, University of Portsmouth, Portsmouth, PO1 3FX, UK

<sup>8</sup> Department of Physics, University of Cincinnati, Cincinnati, OH 45221, USA

<sup>9</sup> CNRS, UMR 7095, Institut d'Astrophysique de Paris, F-75014, Paris, France

<sup>10</sup> Sorbonne Universités, UPMC Univ Paris 06, UMR 7095, Institut d'Astrophysique de Paris, F-75014, Paris, France

<sup>11</sup> Department of Physics and Astronomy, Pavement Building, University of Sussex, Brighton, BN1 9QH, UK

<sup>12</sup> Department of Physics & Astronomy, University College London, Gower Street, London, WC1E 6BT, UK

<sup>13</sup> Department of Astronomy and Astrophysics, University of Chicago, 5640 South Ellis Avenue, Chicago, IL 60637, USA

<sup>14</sup> SLAC National Accelerator Laboratory, Menlo Park, CA 94025, USA

<sup>15</sup> Institute of Cosmology and Gravitation, University of Portsmouth, PO1 2UP, UK

<sup>16</sup> Instituto de Astrofísica de Canarias, E-38205 La Laguna, Tenerife, Spain

<sup>17</sup> Universidad de La Laguna, Dpto. Astrofísica, E-38206 La Laguna, Tenerife, Spain

<sup>18</sup> Center for Astrophysical Surveys, National Center for Supercomputing Applications, 1205 West Clark Street, Urbana, IL 61801, USA

<sup>19</sup> Department of Astronomy, University of Illinois at Urbana-Champaign, 1002 West Green Street, Urbana, IL 61801, USA

<sup>20</sup> Institut de Física d'Altes Energies (IFAE), The Barcelona Institute of Science and Technology, Campus UAB, E-08193 Bellaterra (Barcelona), Spain
<sup>21</sup> Astronomy Unit, Department of Physics, University of Trieste, via Tiepolo 11, I-34131 Trieste, Italy

<sup>22</sup> Osservatorio Astronomico di Trieste, via G. B. Tiepolo 11, I-34143 Trieste, Italy

<sup>23</sup> Institute for Fundamental Physics of the Universe, Via Beirut 2, I-34014 Trieste, Italy

<sup>24</sup> IObservatório Nacional, Rua Gal. José Cristino 77, Rio de Janeiro, RJ—20921-400, Brazil

<sup>25</sup> Department of Physics, University of Michigan, Ann Arbor, MI 48109, USA

<sup>26</sup> Centro de Investigaciones Energéticas, Medioambientales y Tecnológicas (CIEMAT), Madrid, Spain

<sup>27</sup> Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, M/S 169-327, Pasadena, CA 91109, USA

<sup>28</sup> Santa Cruz Institute for Particle Physics, Santa Cruz, CA 95064, USA

<sup>29</sup> Department of Astronomy, University of Michigan, Ann Arbor, MI 48109, USA

<sup>30</sup> Institute of Theoretical Astrophysics, University of Oslo. P.O. Box 1029 Blindern, NO-0315 Oslo, Norway

<sup>31</sup> Department of Physics and Astronomy, University of Missouri-Kansas City, 5110 Rockhill Road, Kansas City, MO 64110, USA

<sup>32</sup> Instituto de Física Teórica UAM/CSIC, Universidad Autonoma de Madrid, E-28049 Madrid, Spain

<sup>33</sup> Institut d'Estudis Espacials de Catalunya (IEEC), E-08034 Barcelona, Spain

<sup>34</sup> Institute of Space Sciences (ICE, CSIC), Campus UAB, Carrer de Can Magrans, s/n, E-08193 Barcelona, Spain

<sup>35</sup> Department of Astronomy, University of Florida, Gainesville, FL 32611, USA

<sup>36</sup> Faculty of Physics, Ludwig-Maximilians-Universität, Scheinerstrasse 1, D-81679 Munich, Germany

<sup>37</sup> CSIRO Astronomy & Space Science, PO Box 1130, Bentley, WA 6102, Australia

<sup>38</sup> School of Mathematics and Physics, University of Queensland, Brisbane, QLD 4072, Australia

<sup>39</sup> Department of Astronomy, The Ohio State University, Columbus, OH 43210, USA

<sup>40</sup> Center of Cosmology and Astro-Particle Physics, The Ohio State University, Columbus, OH 43210, USA

<sup>41</sup> Center for Astrophysics | Harvard & Smithsonian, 60 Garden Street, Cambridge, MA 02138, USA

<sup>42</sup> Kavli Institute for Cosmological Physics, University of Chicago, 5640 South Ellis Avenue, Chicago, IL 60637, USA

<sup>43</sup> Australian Astronomical Optics, Macquarie University, North Ryde, NSW 2113, Australia

<sup>44</sup> Lowell Observatory, 1400 Mars Hill Road, Flagstaff, AZ 86001, USA

<sup>45</sup> Departamento de Física Matemática, Instituto de Física, Universidade de São Paulo, CP 66318, São Paulo, SP, 05314-970, Brazil

<sup>46</sup> George P. and Cynthia Woods Mitchell Institute for Fundamental Physics and Astronomy, and Department of Physics and Astronomy, Texas A&M University, College Station, TX 77843, USA

<sup>47</sup> Radcliffe Institute for Advanced Study, Harvard University, Cambridge, MA 02138, USA

<sup>48</sup> Department of Astrophysical Sciences, Princeton University, Peyton Hall, Princeton, NJ 08544, USA

<sup>49</sup> Institució Catalana de Recerca i Estudis Avançats, E-08010 Barcelona, Spain

<sup>50</sup> Max Planck Institute for Extraterrestrial Physics, Giessenbachstrasse, D-85748 Garching, Germany

<sup>51</sup> Physics Department, 2320 Chamberlin Hall, University of Wisconsin-Madison, 1150 University Avenue Madison, WI 53706-1390, USA

<sup>52</sup> National Institute for Nuclear Physics, Via Valerio 2, I-34127 Trieste, Italy

<sup>53</sup> School of Physics and Astronomy, University of Southampton, Southampton, SO17 1BJ, UK
Received 2021 July 20; revised 2022 January 12; accepted 2022 January 13; published 2022 March 3

Abstract

We present the results of an analysis of Wide-field Infrared Survey Explorer (WISE) observations of the full 2500 deg² South Pole Telescope (SPT)-Sunyaev-Zel'dovich cluster sample. We describe a process for identifying active galactic nuclei (AGN) in brightest cluster galaxies (BCGs) based on WISE mid-IR color and redshift. Applying this technique to the BCGs of the SPT-SZ sample, we calculate the AGN-hosting BCG fraction, which is defined as the fraction of BCGs hosting bright central AGNs over all possible BCGs. Assuming an evolving single-burst stellar population model, we find statistically significant evidence (>99.9%) for a mid-IR excess at high redshift compared to low redshift, suggesting that the fraction of AGN-hosting BCGs increases with redshift over the range of $0 < z < 1.3$. The best-fit redshift trend of the AGN-hosting BCG fraction has the form $(1+z)^{4.1\pm1.0}$. These results are consistent with previous studies in galaxy clusters as well as in field galaxies. One way to explain this result is that member galaxies at high redshift tend to have more cold gas. While BCGs in nearby galaxy clusters grow mostly by dry mergers with cluster members, leading to no increase in AGN activity, BCGs at high redshift could primarily merge with gas-rich satellites, providing fuel for feeding AGNs. If this observed increase in AGN activity is linked to gas-rich mergers rather than ICM cooling, we would expect to see an increase in scatter in the $P_{cav}$ versus $L_{cool}$ relation at $z > 1$. Last, this work confirms that the runaway cooling phase, as predicted by the classical cooling-flow model, in the Phoenix cluster is extremely rare and most BCGs have low (relative to Eddington) black hole accretion rates.

Unified Astronomy Thesaurus concepts: Galaxy clusters (584); Galaxy evolution (594)

1. Introduction

Galaxy clusters are the most massive gravitationally bound and collapsed objects in the universe (Voit 2005). Because of their extremely deep potential wells, the temperature of the intracluster medium (ICM) is high enough to radiate X-rays. The central parts of clusters, which have the densest X-ray emitting gas, often have shorter cooling times than the Hubble time, implying that the hot X-ray gas should have had enough time to cool and form large inward flows of cooling material, known as cooling flows (Sarazin 1986; Fabian 1994). However, multiwavelength observations have only seen a fraction of the massive cooling flows that are expected from standard cooling models (e.g., O'Dea et al. 2008; Donahue et al. 2015; McDonald et al. 2018). This is referred to as “the cooling-flow problem,” and active galactic nucleus (AGN) feedback is thought to be responsible for preventing the hot gas from cooling by propagating energy from the supermassive black hole (SMBH) to the ICM. The two primary modes of AGN feedback are the kinetic mode, in which relativistic jets push the hot gas aside and create cavities, and the quasar mode, in which the radiation comes from the accretion disk (see reviews by Fabian 2012; McNamara & Nulsen 2012).

With the recent development of galaxy cluster surveys that use the Sunyaev-Zel'dovich (SZ) effect (Sunyaev & Zeldovich 1972), such as the South Pole Telescope (SPT; Carlstrom et al. 2011; Bleem et al. 2015, 2020; Huang et al. 2020) and the Atacama Cosmology Telescope (ACT; Hilton et al. 2018, 2021), the number of known high-$z$ ($z > 1$) clusters with good mass estimates has increased dramatically. This has enabled many studies of the evolution of AGN feedback in clusters over cosmic time (McDonald et al. 2013, 2017; Gupta et al. 2020). However, the evolution of AGN feedback in galaxy clusters

with redshift remains poorly understood. In particular, Hlavacek-Larrondo et al. (2015) found no evidence for evolution in jetted power generated by AGN feedback from X-ray cavities over the past 7 Gyr ($z = 0.8$). An earlier study by Hlavacek-Larrondo et al. (2013) suggested that the fraction of brightest cluster galaxies (BCGs) with X-ray bright nuclei is decreasing with time (or increasing with redshift), suggesting a strong evolution in radiative mode feedback. In contrast, a recent study looking for nuclear BCG X-ray emission in Chandra archival data instead found no evidence for an evolution between two redshift bins ($\langle z \rangle \sim 0.25$ and $\langle z \rangle \sim 0.65$; Yang et al. 2018). The disagreement between various studies about the evolution of AGN feedback restricts our ability to fully understand this issue.

In this work, we calculate the AGN-hosting BCG fraction by identifying BCGs in the 2500 deg² SPT-SZ cluster samples (Bleem et al. 2015) and by classifying whether they are AGNs based on mid-IR data. The SZ cluster catalogs allow for an effectively mass-selected sample of clusters, making it possible to study the evolution of galaxy clusters over time. In addition, the SZ catalogs typically have less contamination compared to optical/IR catalogs. The fraction of BCGs hosting luminous AGNs is an important indicator for AGN fueling processes, availability of cold clumpy gas in the centers of clusters, the duration and duty cycle of the AGNs, and how BCGs and the host clusters grow and evolve together. This is because additional physical mechanisms are often required to explain the transport of the cold gas, which serves as the primary fuel source for the central black holes. The fact that we find a relative absence of AGNs in the centers of clusters has led us to study many physical processes, including ram pressure stripping (Gunn & Gott 1972), tidal effects from the cluster gravitational potential (Merritt 1983), and the lack of new infall of cold gas (Larson et al. 1980). Similarly, understanding the evolution of AGN activities in BCGs will help us understand the evolution mechanism of galaxy clusters, and how the feedback might play a role in it.

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

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

<footnote n="54">⁵⁴ Osservatorio Astronomico di Brera, Via Brera 28, I-20121 Milano, Italy</footnote>

<footnote n="55">⁵⁵ Computer Science and Mathematics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA</footnote>

<footnote n="56">⁵⁶ Department of Physics, Stanford University, 382 Via Pueblo Mall, Stanford, CA 94305, USA</footnote>
Our goal for this paper is to study the redshift evolution of the AGN-hosting BCG fraction up to $z = 1.3$ to understand the fueling processes in the centers of clusters, determine when AGN feedback is fully established, and identify whether there are any more extreme AGN-hosting BCGs in the sample, similar to the Phoenix cluster. The paper is organized as follows. In Sections 2 and 3 we summarize the data and additional information used in this paper. The results and their implications are presented in Sections 4 and 5, respectively. We conclude our work in Section 6. We assume $H_0 = 70 \text{ km s}^{-1} \text{ Mpc}^{-1}$, $\Omega_m = 0.3$ and $\Omega_\lambda = 0.7$. All errors are $1\sigma$ unless noted otherwise.

## 2. Data

### 2.1. The SPT-SZ 2500 deg² Cluster Sample

We use the full 2500 deg² SPT-SZ cluster sample from Bleem et al. (2015) with the improvement in the cluster redshift estimates from Bocquet et al. (2019) by incorporating new spectroscopic and improved photometric measurements (Bayliss et al. 2016; Khullar et al. 2019). The survey spans a contiguous 2500 deg² area within a boundary of R.A. = 20 hr–7 hr and decl. = −65° to −40°. Once we limit the redshift range to $0 < z < 1.3$, the total number of clusters in our sample is 475.

### 2.2. Position of BCGs

Given the diversity of BCG colors, morphologies, and assembly state as a function of redshift, typical identification algorithms may be biased when they select BCGs based on single-band fluxes. We have instead developed a novel BCG identification pipeline that uses the full probability distribution of redshift and stellar mass for every object within 500 projected kiloparsec of the SZ cluster center to assign BCG likelihoods. Photometry is provided by the Dark Energy Survey (Year 3) catalogs (Jarvis et al. 2021; Sevilla-Noarbe et al. 2021), cross-correlated with unWISE (Lang 2014), which is a combination of WISE and NEOWISE images. Various cuts and flags are used to avoid stars and objects with poor photometric measurements. WISE is an IR-wavelength space satellite with four IR filters, including W1 ($\lambda_{cen} = 3.6$ µm), W2 (4.3 µm), W3 (12 µm), and W4 (22 µm; Wright et al. 2010). The satellite operated for two years with cryogen until 2011, before it was reactivated and resumed operations as NEOWISE in 2013. It has continued to observe ever since (Mainzer et al. 2014).

Probability distributions of photometric redshift and stellar mass for each source are estimated with EAZY (Brammer et al. 2008) and FAST (Kriek et al. 2009), respectively. We then randomly sample from each distribution to find the most massive cluster galaxy at the cluster redshift within each field, iterating this process $10^5$ times to build up a BCG likelihood for each galaxy. In this way, all galaxies are assigned a value between 0% and 100% probability of being the BCG within each cluster. Full details about the pipeline, along with the BCG catalog, will be provided in Noble et al. (in preparation). The top three panels of Figure 1 show optical images of example SPT galaxy clusters with identified BCGs in white squares. This demonstrates that the algorithm selects likely BCGs that match the galaxies that typical or traditional visual BCG identification methods would select over a wide range of redshift.

### 2.3. Data for AGN Selection

Most photometric techniques for identifying AGNs are severely biased toward unobscured (type 1) AGNs since their nuclear emissions dominate the emission of the host galaxies, making these AGNs easily identifiable. This implies that most obscured (type 2) AGNs are often underrepresented in most studies. The most promising techniques for identifying both obscured and unobscured AGNs include radio, hard X-ray, and mid-IR selections. However, not all AGNs are radio loud (e.g., Stern et al. 2000), and the current hard X-ray satellites remain limited in their sensitivity and field of view. This leaves mid-IR selection as a popular technique to quickly identify large AGN populations (obscured and unobscured). The idea of mid-IR selection is to separate between the power-law AGN spectrum and the blackbody stellar spectrum of galaxies, which peaks at rest frame 1.6 µm. The power-law spectra of the AGNs are due to the thermal emission from the warm-hot dust in the torus, which is heated by absorbing shorter-wavelength photons from the accretion disk (Stern et al. 2012; Hickox & Alexander 2018). This implies that the emission is not strongly suppressed by the dusty torus, unlike at the UV to near-IR wavelength, allowing this technique to detect more obscured AGNs. Additionally, with the first all-sky data release of the Wide-field Infrared Survey Explorer (WISE; Wright et al. 2010) in 2012, mid-IR selection became one of the top methods of probing the AGN population over the entire sky without additional observations.

One drawback of the mid-IR selection technique is that the host galaxy is still bright at these wavelengths, limiting the detection of low-luminosity AGN, which have to compete with a bright stellar continuum (Stern et al. 2012). This means that AGNs that are selected by their mid-IR color tend to be brighter relative to the host galaxies than those selected by other techniques. For example, Assef et al. (2013) found that in the sample of relatively luminous AGNs ($L_{AGN}/L_{host} > 0.5$), the luminosity of mid-IR AGNs tends to be greater than ~5 × $10^{44}$ erg s⁻¹, taking into account the bolometric correction from Singal et al. (2016). Assuming the efficiency of turning accreting matter into energy $\epsilon_{acc} = 0.1$ and a typical mass of the SMBH $M_{BH} \sim 10^8–10^9 M_\odot$ (Russell et al. 2013), the black hole mean accretion rate ($\dot{M}_{BH} = L_{bol,nuc}/\epsilon_{acc}c^2$) of mid-IR selected AGNs should be greater than $4 \times 10^{-3}–4 \times 10^{-2} \dot{M}_{Edd}$, where $\dot{M}_{Edd}$ is the limiting Eddington accretion rate. This level of accretion is relatively high compared to typical optical/radio AGNs, which have an accretion rate of about $10^{-6}–10^{-2} \dot{M}_{Edd}$ (McDonald et al. 2021). This implies that the mid-IR technique will mainly identify a brighter and more massive AGN population. From now on, AGNs mentioned in this paper mean the mid-IR selected AGN population.

#### 2.3.1. Mid-infrared Data from WISE

Instead of using the main source catalog from WISE (AllWISE), which only includes the data obtained from 2010 to 2011, we make use of CatWISE to obtain the best photometry with available data. CatWISE is an updated all-sky IR source catalog that combines the 2010 and 2011 data from WISE with the 2013 through 2016 NEOWISE data (Eisenhardt et al. 2020). The caveat is that the CatWISE catalog only includes 3.6 (W1) and 4.3 (W2) µm data. In this work, we use the Preliminary CatWISE catalog to obtain the
<figure><img src="image_3.png" /><figcaption>**Figure 1.** The top three panels show DES gri optical images of three SPT galaxy clusters, including SPT-CLJ0106-5943, SPT-CLJ0106-5355, and SPT-CLJ2341-5119. The white boxes show the location of the two highest-probability BCG candidates for each cluster, while the white numbers show the probability (in percentage) for each object to be a true BCG. The three examples are ranging from z = 0.35 to z = 1.00, demonstrating the ability of this method to find BCG candidates up to high redshift. The bottom three panels show the corresponding WISE images from the W1 channel for the three SPT clusters. The green boxes show the same location of the two highest-probability BCG candidates, suggesting that the BCG candidates are detectable in mid-IR.</figcaption></figure>

WISE color for each BCG, which has the advantage that more
than four times as many exposures are included as were used
for the AllWISE catalog while using the same AllWISE
software. This makes a comparison with previous works more
straightforward.

For every identified BCG (with a probability >5%) from
Section 2.2, we search for mid-IR counterparts in the CatWISE
catalog within a radius of 3" from the identified BCG because
the typical FWHMs for W1 and W2 are 6"'08 and 6"'84,
respectively. Both W1 and W2 are converted from Vega into
AB magnitudes using the correction from the Explanatory
Supplement to the WISE Products.<fnref n="57" /> The bottom three panels
of Figure 1 show W1 images of the three SPT clusters, showing
that their BCG candidates can be detected with WISE data.

2.3.2. *CatWISE Color Correction*

We perform a comparison test between the AllWISE and
CatWISE catalogs. The test is carried out by comparing the
(W1 − W2)AB color of bright objects (16.5 < W1AB < 18)
between the two catalogs within 10 arcmin of one field. A
histogram of the color differences is plotted in Figure 2,
showing the offset between the color from AllWISE and
CatWISE to be around 0.042. Further investigation reveals that
this is due to the gradual diminishing of the W2 throughput
with time, leading to a bluer (W1 − W2)AB color compared to
AllWISE. We apply this correction of 0.042 mag to CatWISE
W1−W2 colors.

3. Method

3.1. AGN Selection

With the WISE satellite, Stern et al. (2012) developed a
well-known formula for quickly identifying AGN candidates
with a simple color criterion, (W1 − W2)vega ≥ 0.8 or

<figure><img src="image_4.png" /><figcaption>**Figure 2.** This figure shows a comparison of W1–W2 between AllWISE and CatWISE. The red line is the median of the difference at 0.042 mag.</figcaption></figure>

(W1 − W2)AB ≥ 0.16. One benefit of mid-IR selected samples
is that both unobscured (type 1) and obscured (type 2) AGN
can be identified (Lacy et al. 2004; Stern et al. 2005, 2012).
However, because the colors of galaxies drastically change
over a large redshift range, this simple criterion is not accurate
enough to characterize a large population of AGNs. To increase
the number of AGNs we can identify, we develop a new color
criterion that depends on the redshifts of galaxies.

3.1.1. *EzGal Galaxy Color Model*

To determine whether each BCG harbors an AGN, we
calculate the expected color for typical elliptical galaxies using
EzGal.<fnref n="58" /> EzGal calculates the magnitude evolution as a
function of redshift from evolving the spectral energy
distribution (SED) models of a stellar population with time

<footnote n="57">⁵⁷ wise2.ipac.caltech.edu/docs/release/allsky/expsup/sec4_4h.html</footnote>

<footnote n="58">⁵⁸ www.baryons.org/ezgal/</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Top: W1–W2 color for each BCG candidate as a function of cluster redshift. The solid orange line shows the expected color as a function of redshift for our elliptical galaxy model using EzGal, as described in Section 3.1.1. The dashed green line shows our criterion for selecting AGNs, which is derived from the orange line, and the dotted pink line shows the cut from the previous work by Stern et al. (2012). Bottom: the W1–W2 color difference between each BCG candidate and the expected color. The dashed green line shows our selection with the residual >0.2. Every object with a residual greater than 0.2 is likely to be an AGN. The vertical dashed gray lines show the binning for the results in Figure 6. The redshift bins are defined such that each bin contains roughly the same number of systems, making for uniform counting statistics across all redshifts. The two colored stars are known BCGs in galaxy clusters with a high SFR, showing that our selection criterion does not select these starburst BCGs, while the two clusters with luminous AGNs (H 1821+643 and IRAS 09104+4109) are clearly above our criterion.</figcaption></figure>

and projecting them through filters (Mancone & Gonzalez 2012). This calculation takes into account both the stellar evolution of a galaxy as young stars evolve and the wavelength shift due to the distance of a galaxy. To find the model that best describes our overall sample, we perform a grid search between three stellar population model sets (i.e., Bruzual & Charlot 2003; Maraston 2005; Conroy et al. 2009), various formation redshifts, two different initial mass functions (IMF) (i.e., Salpeter 1955; Chabrier 2003), star formation history as a single exponential decaying burst of star formation with an e-folding time parameter ($\tau$) between 0.1 and 10 Gyr, and the representative metallicity (Z) for our galaxy sample from 0.001 to 0.03. Ultimately, the best-fit model (based on the chi-square test) is a Bruzual & Charlot (2003) stellar model with a formation redshift of ($z_{\tau}$) 3.5, a Salpeter (1955) IMF, $\tau$ = 0.1 Gyr for the star formation history, and a metallicity

of 0.016. The solid orange line in the top panel of Figure 3 shows the expected W1–W2 color evolution, generated from the EzGal model with this particular set of parameters. The bottom panel shows the residual from the expected value of W1–W2 for each BCG. It demonstrates that the scatter is distributed around zero with a relatively weak redshift dependence, implying that we have successfully removed the continuum contribution. In this work, galaxies that are redder (a residual greater than 0.2) than typical elliptical galaxies based on the EzGal model are considered AGN candidates. For a range of threshold values from 0.2 to 0.4, the highest redshift bin has more AGN-hosting BCG than the lowest redshift bin, implying that an increase in the fraction of AGN-hosting BCGs is independent of this choice. To further test this notion, we perform a heteroscedasticity test, specifically, the Breusch-Pagan test, in the bottom panel of Figure 3, which shows whether the scatter of the IR residual depends on redshift, regardless of the choice of a threshold value. The test results in a P-value = 0.0047, meaning that the heteroscedasticity is present and the scatter of W1–W2 residuals depends on the redshift, implying that an increase in the fraction is feasible regardless of the threshold choice.

One assumption that we apply in this section is that we only consider a single-burst stellar population model with a single formation redshift, star formation history, and metallicity. In Figure 4 we consider both a single-burst stellar population model and a more complicated two-age stellar population (old and young) model with a wide range of parameters for both models. The model with two stellar populations keeps the same parameter sets from a single-burst model for the “old” population, while a “young” population is represented by a 50 Myr old stellar population at all redshifts. Even though these two models are likely not sufficient to describe our data, more sophisticated models would be unconstrained by the data that we have available. Based on these single-age and two-age models, we find no combination of formation time, metallicity, and IMF that can fully account for the observed evolution in the mid-IR excess, as shown in the right panel of Figure 4. We propose that this mid-IR excess comes from a dusty torus, which is a signature of an actively accreting SMBH. It is difficult to imagine other astronomical sources for this emission because star formation typically yields significantly cooler dust temperatures with a peak brightness of ~100 µm instead of ~1–10 µm. On the other hand, it could also be that our current population models are not adequate to describe the data. We should keep this caveat in mind when we discuss the implication of our results.

As a test to see how a starburst can affect the mid-IR color, we consider Abell 1835 (Ehlert et al. 2011) and RX J1532.9+3021 (Hlavacek-Larrondo et al. 2013), which are the most star-forming BCGs known (SFR ~100 $M_{\odot}$ yr$^{-1}$; McDonald et al. 2018) that also lack evidence of a strong AGN. The two colored stars in Figure 3 demonstrate that even though a star-forming BCG would have boosted mid-IR emission due to dust, polycyclic aromatic hydrocarbon molecules (PAH), and molecular gas, the emission is not as strong as the power-law spectra of AGNs, and our selection does not include these two BCGs. On the other hand, the two clusters with the most luminous AGNs (H 1821+643; Russell et al. 2010, and IRAS 09104+4109; O’Sullivan et al. 2012) are easily detected with our criterion.
<figure><img src="image_6.png" /></figure>

**Figure 4.** Left: W1–W2 colors (black dots) for each BCG candidate in this sample. The solid lines show the expected color as a function of redshift using EzGal models for both a single-burst stellar population (blue) and the two-age stellar populations (orange). The single-burst models are plotted from a wide range of parameters, including two IMF (Salpeter 1955; Chabrier 2003), three formation redshifts ($z_f = 1.5, 2.5, \text{and } 3.5$), and a span of metallicity from 0.75 to 1.0 $Z_\odot$. To obtain the best fit, we perform a grid search between the IMF and the three formation redshifts before fitting the metallicity to minimize the chi-square. The old and young stellar models use the same parameter sets, but include an additional “young” population, which is represented by a 50 Myr old stellar population at all observed redshifts. Right: The fraction of AGN-hosting BCGs as a function of redshift for all models. The figure demonstrates that the fraction of AGN-hosting BCGs increases with redshift regardless of our choice of stellar population model.

### 3.1.2. Spitzer Color Verification

Because the point-spread functions (PSFs) of the two WISE bands are not small (PSF $w_1 = 6''08$ and PSF $w_2 = 6''84$), we compare the results from WISE mid-IR color with those from the Spitzer Space Telescope. Spitzer is an IR telescope with the Infrared Array Camera (IRAC; Fazio et al. 2004) as one of its main science instruments. IRAC is a four-channel imaging camera capable of taking simultaneous images at wavelengths of 3.6, 4.5, 5.8, and 8.0 $\mu$m. Thus, channels 1 and 2 ([3.6] and [4.5]) on IRAC are roughly equivalent with W1 and W2 from WISE, but with the benefit of having a much better PSF at 1"95 and 2"02, respectively (Fazio et al. 2004).

A certain fraction (~35%, predominantly at $z>0.8$) of the SPT cluster sample has been observed with IRAC. For verification, we compare the [3.6]–[4.5] colors of our AGN-hosting BCG candidates with their W1–W2 colors. If the Spitzer color, which has a higher angular resolution, is bluer (smaller) than the WISE color, it shows that there is a contamination from nearby galaxies within the WISE aperture. On the other hand, if the Spitzer color is redder (larger), it implies that the object is even more likely to be an AGN. Figure 5 shows the comparison between WISE's W1–W2 (gray squares) and Spitzer's [3.6]–[4.5] color (circles) for our AGN candidates. We find that most AGN candidates have a difference of Spitzer and WISE color that is either compatible (60% of the sources have a difference of within ±0.07 mag, which is roughly the mean of the WISE color uncertainty) or that Spitzer is slightly redder (33% of the sources have a Spitzer color larger than the WISE color by ~0.2 mag). This suggests that most of our AGN candidates are likely to be real quasars. One clear exception is SPT-CL J2146-4633, which has a WISE color much redder than Spitzer. Further investigation shows that there is a point-like source near the location of the object, but not at the BCG location, meaning that the WISE

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

**Figure 5.** Residual plot similar to the bottom panel of Figure 3. The blue points are the W1–W2 color difference between each BCG candidate and the expected color. The gray squares emphasize the W1–W2 color from WISE for the sample that has been observed by Spitzer, while the green circles show the color from Spitzer, demonstrating that most of the objects we classify as AGNs have a difference of Spitzer and WISE color that is either compatible (60% of the sources have a difference of within ±0.07 mag) or Spitzer is slightly redder (33% of the sources have a Spitzer color larger than the WISE color by ~0.2 mag). This means that the color and our results are not strongly impacted by WISE's larger PSF.

color is probably contaminated by a nearby AGN, while the Spitzer color is not. This object has been removed from the further analysis.
<table><caption>**Table 1**
Total Number of BCGs and Those Hosting AGN in Each Redshift Bin</caption><thead><tr><td>Redshift Bin</td><td>AGN-hosting BCG</td><td>All BCG</td></tr></thead><tbody><tr><td>0.00–0.39</td><td>0</td><td>141</td></tr><tr><td>0.39–0.55</td><td>3</td><td>142</td></tr><tr><td>0.55–0.71</td><td>7</td><td>143</td></tr><tr><td>0.71–1.30</td><td>12</td><td>137</td></tr></tbody></table>

### 3.2. *Calculating the AGN fraction*

We compute the fraction of galaxy clusters with AGN-hosting BCGs based on the number of BCGs whose mid-IR colors are redder than expected, as described in Section 3.1.1, in four redshift bins to study the redshift evolution. With the probability estimated in Section 2.2, we first include in Figure 3 all the BCGs that are identified with a probability that is higher than 20%, meaning that some clusters will have more than one BCG candidate. The fraction of AGN-hosting BCGs is calculated over the total number of BCG candidates, instead of over the total number of clusters in each bin. We discuss these particular choices of calculating AGN-hosting BCG fraction in Section 4. The bins are defined from $z=0–1.3$ in such a way that each bin contains roughly the same number of BCG candidates (~140 BCGs; see Table 1 for the exact number). This choice of binning yields uniform counting statistics across all redshifts. The uncertainties associated with the AGN fractions are estimated from the Wilson interval, which remains accurate for fractions near 0 and 1 (Brown et al. 2001).

## 4. Results and Verification

Table 1 shows the number of AGN-hosting BCG in the four redshift bins, whereas the blue points in Figure 6 show the fraction of AGN-hosting BCGs in the four bins with their corresponding 68% confidence intervals. We observe that the fraction is increasing with redshift in the SPT sample. Gray crosses show the fraction of points that have residuals lower than -0.2. If the scatter is truly random, there should be a similar number of points below -0.2 compared to those above +0.2, which are classified as AGNs. This result implies that the increasing trend of the AGN-hosting BCGs fraction is not a result of the data quality. Such a trend has been suggested and shown in previous works (e.g., Hlavacek-Larrondo et al. 2013; McDonald et al. 2016; Bîrzan et al. 2017; Gupta et al. 2020; Mo et al. 2020). In particular, we show that the fraction is ~2% at $z \approx 0.5$, which is consistent with what Somboonpanyakul et al. (2021a) found from studying extreme central BCGs in clusters. We note that because some AGN-dominated galaxies will have poor photometric redshift constraints as they are estimated from the stellar spectrum and not from the power-law spectrum of the AGN, we might misidentify these galaxies in our BCG-finding algorithm. This implies that the number of BCGs with central AGNs found in this work gives a lower limit for the AGN-hosting BCG fraction, and the actual evolution could be even stronger.

One way to evaluate the observed result with the chosen bins is to perform Fisher's Exact Test in order to see whether there are any associations between two categorical variables. The result shows that we can reject the null hypothesis of independence with P-value = 0.00045, meaning that there is

<figure><img src="image_8.png" /><figcaption>**Figure 6.** The fractions of AGN-hosting BCGs as a function of redshift. Blue points show the fractions from this work, which come from the W1–W2 color residual in the SPT sample. The size of the error bar takes the binomial uncertainty into account. Gray crosses show the fraction of points that has residuals lower than -0.2, which we consider a “noise level.” This figure demonstrates that the fraction of AGN-hosting BCG increases with redshift.</figcaption></figure>

a statistically significant association (>99.9%) between redshift and whether or not BCGs host AGNs.

The approaches taken in this work are as follows: (i) we include all BCG candidates with a probability higher than 20% in our sample instead of picking only one BCG per cluster, and (ii) we calculate the fraction of BCGs with AGNs over the total number of BCG candidates, and not the fraction of clusters with AGN-hosting BCGs over the total number of clusters. The reason for these two assumptions is that we want to include AGN-hosting BCGs from systems with more than one obvious BCG, which are typical for merging systems such as the Coma cluster (Zwicky 1933), and the Bullet cluster (Markevitch et al. 2004). We perform consistency checks to address both of these assumptions. Figure 7 shows the fraction of AGN-hosting BCGs when we consider the most likely BCG candidates, every BCG candidate with a probability higher than 20%, and every BCG candidate with a probability higher than 10%. This figure shows that the increasing trend of the fraction of AGN-hosting BCGs over redshift remains consistent in all three scenarios, regardless of how we select BCGs.

On the other hand, Figure 8 shows the results when we use different definitions of AGN fractions. The gray and blue points in Figure 8 are calculated with the total number of clusters as a denominator instead of the number of BCG candidates. For the blue points, we consider one BCG per cluster and include both the probability of being BCGs, as calculated in Section 2.2, and the uncertainty of the mid-IR color for each BCG to emphasize the fact that the uncertainties of identifying BCGs and BCG colors are higher at high redshift. The empty gray dots are the largest possible fractions, which are calculated from clusters for which any of their potential BCGs are considered as AGNs, while the empty gray squares are the smallest possible fractions by counting only clusters for which all of their BCG candidates are classified as AGNs. This figure illustrates that all of these definitions qualitatively give the same conclusion as our initial results.
<figure><img src="image_9.png" /><figcaption>**Figure 7.** All three panels show the fractions of AGN-hosting BCGs, similar to Figure 6. The left panel (triangles) shows the fraction when we consider only the most likely BCG for each cluster. The middle panel (dots) includes all possible BCGs with a probability higher than 20%, while the right panel (squares) includes those with a probability higher than 10%. In all scenarios, the fraction increases as a function of redshift.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 8.** The fractions of clusters with AGN-hosting BCGs, similar to Figure 6. The blue points assume one BCG per cluster and incorporate the probability of being a BCG, estimated from Section 2.2. The empty gray dots are for clusters for which any of their potential BCGs are considered as AGNs, which is equivalent to the maximum fraction. The empty gray squares only include clusters whose BCG candidates are all considered AGNs, which is the minimum fraction. This figure shows that all of these definitions qualitatively give the same conclusion as our initial results.</figcaption></figure>

We compare our results with the AGN fraction in field galaxies to determine whether there is a difference in the fractions between the two environments. The green and pink squares in the left panel of Figure 9 show the field X-ray AGN fractions from the zCOSMOS survey (Silverman et al. 2009) and the Chandra Multiwavelength Project results (ChaMP; Haggard et al. 2010), respectively. The results from our work are consistent with these two results, suggesting that the source of fuel required for AGN accretion in field galaxies could be similar to that in the brightest cluster galaxies. Additional evidence for the AGN fraction evolution in field galaxies has been seen in other works. For example, Lehmer et al. (2007) reported an evolution in early-type galaxies ($z \sim 0.7$) consistent

with the $(1+z)^3$ pure luminosity evolution model. The dotted gray line in Figure 9 shows the curve for $(1+z)^3$, although it is only intended to be illustrative because it is arbitrarily normalized. The dashed line instead shows the curve for $(1+z)^{5.3}$. This line was first suggested by Martini et al. (2009), who showed that the AGN fraction of cluster members increases as $(1+z)^{5.3}$ for AGN with an X-ray luminosity higher than $L_x > 10^{43}$ erg s$^{-1}$, hosted by luminous galaxies. We also fit the power-law model $(\propto(1+z)^\alpha)$ to the blue points in Figure 9 and find a power-law exponent $\alpha = 4.1 \pm 1.0$, as shown by the dash-dotted brwon line, which is consistent with the results from both Lehmer et al. (2007) and Martini et al. (2009). Nevertheless, there are caveats regarding the relation between cluster BCGs and field galaxies. One concern is that the AGN selection criteria for both BCGs and field galaxies are different, making it difficult to make a direct comparison between the two. In addition, according to the work about the evolution of AGN luminosity, which shows that AGNs in galaxies tend to be brighter at high redshift (Hasinger et al. 2005; Silverman et al. 2008), we would naturally expect to find a higher AGN fraction at high redshift because we usually selected AGNs based on a certain luminosity threshold.

## 5. Discussion

The results obtained in Section 4 demonstrate that given a single-burst stellar population model, there is an increase in the fraction of AGN-hosting BCGs with redshift. This finding is consistent with previously published studies (e.g., Hlavacek-Larrondo et al. 2013; Bîrzan et al. 2017), which focused on different samples with distinct selection effects. In particular, the results from Hlavacek-Larrondo et al. (2013) and ours, as shown in the middle panel of Figure 9, show the same trend of an increasing fraction of AGN-hosting BCGs with redshift. However, the normalizations are vastly different. Hlavacek-Larrondo et al. (2013) claimed that the fraction of active BCGs is 30% at $z \approx 0.1$ and 60% at $z \approx 0.5$. On the other hand, we show that the fraction is smaller than 20% at all redshift bins. A possible explanation stems from the very different samples that these works consider. While this study is based on an effectively mass-selected sample of clusters, without any consideration of X-ray properties, the sample used by Hlavacek-Larrondo et al. (2013) focused solely on highly
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Comparison of the fraction of clusters with AGN-hosting BCGs to other published works. Left: Green and pink squares are the field AGN fraction from Silverman et al. (2009) and Haggard et al. (2010), respectively. The dotted gray line indicates pure luminosity evolution $\propto (1+z)^3$ as suggested by Lehmer et al. (2007), while the dashed gray line shows the $(1+z)^{5.3}$ line suggested by Martini et al. (2009). The dash-dotted brown line displays the best-fit model in the form $(1+z)^{4.1\pm1.0}$. Middle: Orange diamonds show the estimated fractions of BCGs with bright AGNs ($\log_{10}(L_{X,\text{nuc}} [\text{erg s}^{-1}]) > 42.2$) from Hlavacek-Larrondo et al. (2013), and open orange diamond showing the scaled version by changing the total number of clusters to the larger parent sample that Hlavacek-Larrondo et al. (2013) had drawn from because this work only focuses on highly X-ray luminous clusters ($L_X > 3 \times 10^{44} \text{ erg s}^{-1}$). Right: Purple points show the fraction of BCGs with an SFR > 100 $M_{\odot} \text{ yr}^{-1}$ from McDonald et al. (2016).</figcaption></figure>

X-ray luminous clusters ($L_{X,\text{cluster}} > 3 \times 10^{44} \text{ erg s}^{-1}$), which show clear X-ray cavities. To place these two studies on the same scale, we modify the denominator used by Hlavacek-Larrondo et al. (2013) in calculating the AGN fraction to account for the full parent population of clusters from which their sample of 32 clusters was drawn. Given that this previous work included a subsample of clusters drawn from the REFLEX (Böhringer et al. 2004), eBCS (Ebeling et al. 2000), MACS (Ebeling et al. 2001), and SPT-XVP (McDonald et al. 2013) surveys, we consider these surveys in their entirety as the total population (the denominator) when calculating the AGN fraction. With this rescaling, we find that the fraction of BCGs hosting powerful AGN is consistent between our work and Hlavacek-Larrondo et al. (2013), and that the observed evolution is consistent in the two studies. We note that AGN are identified in very different ways in these two samples, but they both select at the extreme end—the few most mechanically powerful and the most IR-bright outbursts. The fact that these two works agree after the aforementioned rescaling is applied is reassuring, although a much more thorough analysis (and proper bias correction) is needed before conclusions about the coevolution of jet power and mid-IR emission can be made.

An increase in the fraction of BCGs hosting central AGNs with redshift suggests that the accretion rates of the SMBHs in the BCGs are higher at high redshift because the AGN luminosity is proportional to the accretion rate. Several works about the relation between the mean black hole accretion rates and the cavity (kinetic)/quasar (radiative) power of the central AGN (Churazov et al. 2005; Russell et al. 2013) have shown that as black hole accretion increases in the BCGs, the cavity power of the AGN also increases to counteract the cooling from the accretion in the form of a negative feedback cycle. However, as the black hole accretion rates rise to near the Eddington limit, the cavity/jet power seems to be saturated, and the radiative power tends to dominate at this level of accretion. The fact that the radiative power from AGNs usually promotes more cooling in the ICM instead of preventing it suggests that a well-regulated feedback system between the

central black hole and its host cluster is no longer possible at high accretion, implying that some galaxies might not have a fully established AGN feedback loop at this redshift range.

A similar conclusion has been reached in the work related to star-forming galaxies (Webb et al. 2015; McDonald et al. 2016; Bonaventura et al. 2017), which showed that the fraction of starburst BCGs is higher at high redshift ($z > 1$). Specifically, McDonald et al. (2016) found the fraction of BCGs with SFR over $10 M_{\odot} \text{ yr}^{-1}$ to be $34 \pm 5\%$ at $0.25 < z < 1.25$, compared to $\sim 1\% - 5\%$ at $z \sim 0$. The right panel of Figure 9 compares the fraction of AGN-hosting BCGs in this work with the fraction of starburst BCGs ($\text{SFR}_{\text{BCG}} > 100 M_{\odot} \text{ yr}^{-1}$) from McDonald et al. (2016), demonstrating that in the center of cluster environments, both massive starburst galaxies and bright AGNs behave similarly. These two results strongly indicate that AGN feedback might not be as effective to prevent overcooling at high redshift as we have previously thought.

All of these results lead us to suspect that the reason for the observed redshift trend and the breakdown of AGN feedback at high redshift comes from the fact that there is an abundance of cold gas at that redshift. In the local universe, BCGs typically grow by merging with gas-poor satellites without triggering any AGN activity. However, BCGs at high redshift could grow by merging with gas-rich members instead. Cold gas from the mergers could be a source of fuel for increasing AGN activities in the center of clusters. This is consistent with the picture we obtain from the studies of starburst BCGs (Webb et al. 2015; McDonald et al. 2016) because cold dense clouds from gas-rich mergers could provide enough matter required for creating new stars. Further evidence supporting a gas-rich merger explanation includes the prevalence of cluster galaxies with massive CO or cold gas reservoirs at high redshift (Noble et al. 2017, 2019; Hayashi et al. 2018; Markov et al. 2020) and the detections of molecular gas in many BCGs (Dunne et al. 2021). This scenario can also explain recent studies about the cool-core (CC) fraction, which shows no sign of evolution over the same redshift range (McDonald et al. 2017; Ruppin et al. 2021). If AGN feedback breaks down at high redshift, one
would expect that the CC fractions of clusters would be higher because more gas should have been cooled near the center. However, if black hole accretion and star formation in high-redshift BCGs are fueled by something other than cooling of the hot gas, such as gas-rich mergers (Barnes & Hernquist 1991; Hopkins et al. 2006), it would be reasonable to think that the trends of AGN- and starburst-hosting BCG fractions would be different from the trend of the CC fraction. If this observed increase in AGN activity is linked to gas-rich mergers rather than ICM cooling, we would expect to see an increase in scatter in the $P_{cav}$ versus $L_{cool}$ relation (Rafferty et al. 2006) at $z > 1$.

Another possible scenario to explain the trend of high AGN-hosting BCG fraction at high redshift has to do with cluster mergers. It has been shown both in simulations (Fakhouri et al. 2010) and observations (McDonald et al. 2017) that the cluster merger rate is significantly higher at high redshift. Major mergers between two clusters have the potential of disrupting a tightly regulated AGN feedback loop and promote black hole accretion and star formation by potentially increasing the local turbulence of the system. This is consistent with the turbulent picture in the precipitation model for AGN feedback, called “chaotic cold accretion (CCA),” which states that turbulence is a key component in driving nonlinear thermal instability and extended condensation (Voit et al. 2015; Gaspari et al. 2020). Turbulent forcing can help stimulate precipitation and condensation by increasing the velocity dispersion of the ambient medium, resulting in more black hole accretion and star formation (Voit 2018). With a more energetic environment in the early universe, it is reasonable to assume that the turbulence will be higher at high redshift, resulting in higher black hole accretion rates. The recent discovery of CHIPS 1911 + 4455, a merging galaxy cluster with a massive starburst at the center, provides strong evidence that mergers can indeed increase star formation (Somboonpanyakul et al. 2021a, 2021b). With the development of next-generation X-ray observatories such as Athena and Lynx, we will be able to directly measure motions in the hot gas and determine whether mergers of groups/clusters can boost cooling via an increase in turbulence.

Last, Figure 5 shows that the BCG of the Phoenix cluster remains the most extreme AGN in the entire SPT-SZ sample, which covers a 2500 deg² area and spans all redshifts. In combination with the recent work from the CHIPS survey (Somboonpanyakul et al. 2021a), which confirmed that the Phoenix cluster hosts the most extreme BCG with the strongest CC at $z < 0.7$, the runaway cooling phase, as we have seen in the Phoenix cluster, is indeed extremely rare.

## 6. Conclusion

In this work, we present results for the mid-IR colors of BCGs in SPT-selected galaxy clusters at $0 < z < 1.3$. This study allows us to track the evolution of BCG properties over ~9 Gyr of cluster growth. In particular, we focus our work on black hole accretion in BCGs, which turns these central galaxies into bright AGNs. Our findings are summarized as follows:

1. Assuming a single-burst stellar population model, we find statistically significant evidence ($> 99.9\%$) for a mid-IR excess in high-redshift BCGs compared to low-redshift BCGs, suggesting an increase with redshift in the fraction of AGN-hosting BCGs in galaxy clusters over $0 < z < 1.3$. For the lower-redshift bins ($z < 0.6$), an

increase is not statistically significant, and the results are compatible with the noise level. On the other hand, we see an increase in the fraction of BCGs with AGNs in high-redshift bins ($z > 0.6$), similar to what others have found in previous works (Hlavacek-Larrondo et al. 2013; McDonald et al. 2016; Bîrzan et al. 2017).

2. We show that our results are consistent with both the evolution of the fraction of AGNs in field galaxies (Silverman et al. 2009; Haggard et al. 2010) and the fraction of starburst BCGs (Webb et al. 2015; McDonald et al. 2016), suggesting that the reason for the evolution of both AGN and starburst fraction could come from the fact that more cold gas is available in the early universe. This should lead to a higher level of gas-rich mergers in BCGs, which could fuel both AGN activity and star formation at the center of clusters. Some caveats about the direct comparison between cluster and field galaxies remain and range from selection criteria to the evolution of AGN luminosity.

3. Another possible explanation for the increase in the fraction of AGN-hosting BCGs with redshift could be a higher level of local turbulence from dynamically active galaxy clusters at high redshift, leading to elevated cooling and subsequent black hole accretion. However, for this scenario, it is difficult to explain the similarity to the trends in field galaxies.

4. We do not see any additional cluster with a BCG that is as extreme in the mid-IR color as the Phoenix cluster. In other words, the Phoenix cluster likely hosts the most extreme central AGN in the SPT sample.

An enhancement of AGN activity in BCGs at high redshift compared to low redshift, which is similar in magnitude to the increase observed in field galaxies and cluster members, suggests that this increased AGN activity is not related to cooling flows, but rather to accretion of gas-rich satellites at early times. These accretion events ought to throw off the precise cooling/feedback balance in the centers of clusters that responsible for preventing runaway cooling flows, which leads to a less tightly regulated feedback loop at early times. Further studies with deeper and higher angular resolution mid-IR imaging, such as the upcoming James Webb Space Telescope (JWST; Gardner et al. 2006), will be required to better understand the evolution of AGN feedback and its impact on galaxy clusters.

T.S. and M.M. acknowledge support from the Kavli Research Investment Fund at MIT, and from NASA through Chandra grant GO5-16143.

F.R. acknowledges financial support provided by NASA through SAO Award Number SV2-82023 issued by the Chandra X-Ray Observatory Center, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of NASA under contract NAS8-03060.

This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration.

The South Pole Telescope is supported by the National Science Foundation through grant PLR-1248097. Partial support is also provided by the NSF Physics Frontier Center grant PHY-1125897 to the Kavli Institute of Cosmological
Physics at the University of Chicago, the Kavli Foundation, and the Gordon and Betty Moore Foundation grant GBMF 947.

Funding for the DES Projects has been provided by the U.S. Department of Energy, the U.S. National Science Foundation, the Ministry of Science and Education of Spain, the Science and Technology Facilities Council of the United Kingdom, the Higher Education Funding Council for England, the National Center for Supercomputing Applications at the University of Illinois at Urbana-Champaign, the Kavli Institute of Cosmological Physics at the University of Chicago, the Center for Cosmology and Astro-Particle Physics at the Ohio State University, the Mitchell Institute for Fundamental Physics and Astronomy at Texas A&M University, Financiadora de Estudos e Projetos, Fundação Carlos Chagas Filho de Amparo à Pesquisa do Estado do Rio de Janeiro, Conselho Nacional de Desenvolvimento Científico e Tecnológico and the Ministério da Ciência, Tecnologia e Inovação, the Deutsche Forschungsgemeinschaft and the Collaborating Institutions in the Dark Energy Survey.

The Collaborating Institutions are Argonne National Laboratory, the University of California at Santa Cruz, the University of Cambridge, Centro de Investigaciones Enrgeticas, Medioambientales y Tecnológicas-Madrid, the University of Chicago, University College London, the DES-Brazil Consortium, the University of Edinburgh, the Eidgenössische Technische Hochschule (ETH) Zürich, Fermi National Accelerator Laboratory, the University of Illinois at Urbana-Champaign, the Institut de Ciències de l'Espai (IEEC/CSIC), the Institut de Física d'Altes Energies, Lawrence Berkeley National Laboratory, the Ludwig-Maximilians-Universität München and the associated Excellence Cluster Universe, the University of Michigan, the National Optical Astronomy Observatory, the University of Nottingham, The Ohio State University, the OzDES Membership Consortium, the University of Pennsylvania, the University of Portsmouth, SLAC National Accelerator Laboratory, Stanford University, the University of Sussex, and Texas A&M University.

**Facilities:** Pan-STARRS, Chandra X-ray Observatory (ACIS), Hubble Space Telescope (ACS/WFC), the Nordic Optical Telescope (ALFOSC).

**Software:** astropy (Astropy Collaboration et al. 2018), CIAO (Fruscione et al. 2006), pandas (McKinney 2010), seaborn (Waskom et al. 2016).

## ORCID iDs

T. Somboonpanyakul https://orcid.org/0000-0003-3521-3631

M. McDonald https://orcid.org/0000-0001-5226-8349

M. Aguena https://orcid.org/0000-0001-5679-6747

M. B. Bayliss https://orcid.org/0000-0003-1074-4807

E. Bertin https://orcid.org/0000-0002-3602-3664

D. Brooks https://orcid.org/0000-0002-8458-5047

D. L. Burke https://orcid.org/0000-0003-1866-1950

M. Calzadilla https://orcid.org/0000-0002-2238-2105

A. Carnero Rosell https://orcid.org/0000-0003-3044-5150

M. Carrasco Kind https://orcid.org/0000-0002-4802-3194

J. Carretero https://orcid.org/0000-0002-3130-0204

M. Costanzi https://orcid.org/0000-0001-8158-1449

L. N. da Costa https://orcid.org/0000-0002-7731-277X

J. De Vicente https://orcid.org/0000-0001-8318-6813

A. E. Evrard https://orcid.org/0000-0002-4876-956X

I. Ferrero https://orcid.org/0000-0002-1295-1132

B. Floyd https://orcid.org/0000-0003-4175-571X

J. García-Bellido https://orcid.org/0000-0002-9370-8360

E. Gaztanaga https://orcid.org/0000-0001-9632-0815

D. W. Gerdes https://orcid.org/0000-0001-6942-2736

A. Gonzalez https://orcid.org/0000-0002-0933-8601

D. Gruen https://orcid.org/0000-0003-3270-7644

R. A. Gruendl https://orcid.org/0000-0002-4588-6517

J. Gschwend https://orcid.org/0000-0003-3023-8362

N. Gupta https://orcid.org/0000-0001-7547-4241

S. R. Hinton https://orcid.org/0000-0003-2071-9349

D. L. Hollowood https://orcid.org/0000-0002-9369-4157

T. Jeltema https://orcid.org/0000-0001-6089-0365

G. Khullar https://orcid.org/0000-0002-3475-7648

K. J. Kim https://orcid.org/0000-0001-6505-0293

K. Kuehn https://orcid.org/0000-0003-0120-0808

J. L. Marshall https://orcid.org/0000-0003-0710-9474

P. Martini https://orcid.org/0000-0002-4279-4182

P. Melchior https://orcid.org/0000-0002-8873-5065

F. Menanteau https://orcid.org/0000-0002-1372-2534

R. Miquel https://orcid.org/0000-0002-6610-4836

R. Morgan https://orcid.org/0000-0002-7016-5471

R. L. C. Ogando https://orcid.org/0000-0003-2120-1154

A. Palmese https://orcid.org/0000-0002-6011-0530

F. Paz-Chinchón https://orcid.org/0000-0003-1339-2683

A. Pieres https://orcid.org/0000-0001-9186-6042

A. A. Plazas Malagón https://orcid.org/0000-0002-2598-0514

A. K. Romer https://orcid.org/0000-0002-9328-879X

F. Ruppin https://orcid.org/0000-0002-0955-8954

E. Sanchez https://orcid.org/0000-0002-9646-8198

S. Serrano https://orcid.org/0000-0002-0211-2861

I. Sevilla-Noarbe https://orcid.org/0000-0002-1831-1953

M. Smith https://orcid.org/0000-0002-3321-1432

M. Soares-Santos https://orcid.org/0000-0001-6082-8529

V. Strazzullo https://orcid.org/0000-0001-7975-2894

M. E. C. Swanson https://orcid.org/0000-0002-1488-8552

G. Tarle https://orcid.org/0000-0003-1704-0781

C. To https://orcid.org/0000-0001-7836-2261

D. L. Tucker https://orcid.org/0000-0001-7211-5729

## References

Assef, R. J., Stern, D., Kochanek, C. S., et al. 2013, ApJ, 772, 26

Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123

Barnes, J. E., & Hernquist, L. E. 1991, ApJL, 370, L65

Bayliss, M. B., Ruel, J., Stubbs, C. W., et al. 2016, ApJS, 227, 3

Bírzan, L., Rafferty, D. A., Brüggen, M., et al. 2017, MNRAS, 471, 1766

Bleem, L. E., Bocquet, S., Stalder, B., et al. 2020, ApJS, 247, 25

Bleem, L. E., Stalder, B., de Haan, T., et al. 2015, ApJS, 216, 27

Bocquet, S., Dietrich, J. P., Schrabback, T., et al. 2019, ApJ, 878, 55

Böhringer, H., Schuecker, P., Guzzo, L., et al. 2004, A&A, 425, 367

Bonaventura, N. R., Webb, T. M. A., Muzzin, A., et al. 2017, MNRAS, 469, 1259

Brammer, G. B., van Dokkum, P. G., & Coppi, P. 2008, ApJ, 686, 1503

Brown, L. D., Cai, T. T., & DasGupta, A. 2001, StaSc, 16, 101

Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000

Carlstrom, J. E., Ade, P. A. R., Aird, K. A., et al. 2011, PASP, 123, 568

Chabrier, G. 2003, PASP, 115, 763

Churazov, E., Sazonov, S., Sunyaev, R., et al. 2005, MNRAS, 363, L91

Conroy, C., Gunn, J. E., & White, M. 2009, ApJ, 699, 486

Donahue, M., Connor, T., Fogarty, K., et al. 2015, ApJ, 805, 177

Dunne, D. A., Webb, T. M. A., Noble, A., et al. 2021, ApJL, 909, L29

Ebeling, H., Edge, A. C., Allen, S. W., et al. 2000, MNRAS, 318, 333

Ebeling, H., Edge, A. C., & Henry, J. P. 2001, ApJ, 553, 668

Ehlert, S., Allen, S. W., von der Linden, A., et al. 2011, MNRAS, 411, 1641

Eisenhardt, P. R. M., Marocco, F., Fowler, J. W., et al. 2020, ApJS, 247, 69
Fabian, A. C. 1994, ARA&A, 32, 277

Fabian, A. C. 2012, ARA&A, 50, 455

Fakhouri, O., Ma, C.-P., & Boylan-Kolchin, M. 2010, MNRAS, 406, 2267

Fazio, G. G., Hora, J. L., Allen, L. E., et al. 2004, ApJS, 154, 10

Fruscione, A., McDowell, J. C., Allen, G. E., et al. 2006, Proc. SPIE, 6270, 62701V

Gardner, J. P., Mather, J. C., Clampin, M., et al. 2006, SSRv, 123, 485

Gaspari, M., Tombesi, F., & Cappi, M. 2020, Nature Astronomy, 4, 10

Gunn, J. E., & Gott, J. R. I. 1972, ApJ, 176, 1

Gupta, N., Pannella, M., Mohr, J. J., et al. 2020, MNRAS, 494, 1705

Haggard, D., Green, P. J., Anderson, S. F., et al. 2010, ApJ, 723, 1447

Hasinger, G., Miyaji, T., & Schmidt, M. 2005, A&A, 441, 417

Hayashi, M., Tadaki, K., Kodama, T., et al. 2018, ApJ, 856, 118

Hickox, R. C., & Alexander, D. M. 2018, ARA&A, 56, 625

Hilton, M., Hasselfield, M., Sifón, C., et al. 2018, ApJS, 235, 20

Hilton, M., Sifón, C., Naess, S., et al. 2021, ApJS, 253, 3

Hlavacek-Larrondo, J., Allen, S. W., Taylor, G. B., et al. 2013, ApJ, 777, 163

Hlavacek-Larrondo, J., Fabian, A. C., Edge, A. C., et al. 2013, MNRAS, 431, 1638

Hlavacek-Larrondo, J., McDonald, M., Benson, B. A., et al. 2015, ApJ, 805, 35

Hopkins, P. F., Hernquist, L., Cox, T. J., et al. 2006, ApJS, 163, 1

Huang, N., Bleem, L. E., Stalder, B., et al. 2020, AJ, 159, 110

Jarvis, M., Bernstein, G. M., Amon, A., et al. 2021, MNRAS, 501, 1282

Khullar, G., Bleem, L. E., Bayliss, M. B., et al. 2019, ApJ, 870, 7

Kriek, M., van Dokkum, P. G., Labbé, I., et al. 2009, ApJ, 700, 221

Lacy, M., Storrie-Lombardi, L. J., Sajina, A., et al. 2004, ApJS, 154, 166

Lang, D. 2014, AJ, 147, 108

Larson, R. B., Tinsley, B. M., & Caldwell, C. N. 1980, ApJ, 237, 692

Lehmer, B. D., Brandt, W. N., Alexander, D. M., et al. 2007, ApJ, 657, 681

Mainzer, A., Bauer, J., Cutri, R. M., et al. 2014, ApJ, 792, 30

Mancone, C. L., & Gonzalez, A. H. 2012, PASP, 124, 606

Maraston, C. 2005, MNRAS, 362, 799

Markevitch, M., Gonzalez, A. H., Clowe, D., et al. 2004, ApJ, 606, 819

Markov, V., Mei, S., Salomé, P., et al. 2020, A&A, 641, A22

Martini, P., Sivakoff, G. R., & Mulchaey, J. S. 2009, ApJ, 701, 66

McDonald, M., Allen, S. W., Bayliss, M., et al. 2017, ApJ, 843, 28

McDonald, M., Benson, B. A., Vikhlinin, A., et al. 2013, ApJ, 774, 23

McDonald, M., Gaspari, M., McNamara, B. R., et al. 2018, ApJ, 858, 45

McDonald, M., McNamara, B. R., Calzadilla, M. S., et al. 2021, ApJ, 908, 85

McDonald, M., Stalder, B., Bayliss, M., et al. 2016, ApJ, 817, 86

McKinney, W. 2010, in Proc. IX Python in Science Conf., ed. S. van der Walt & J. Millman (Austin, TX: SciPy), 51

McNamara, B. R., & Nulsen, P. E. J. 2012, NJPh, 14, 055023

Merritt, D. 1983, ApJ, 264, 24

Mo, W., Gonzalez, A., Brodwin, M., et al. 2020, ApJ, 901, 131

Noble, A. G., McDonald, M., Muzzin, A., et al. 2017, ApJL, 842, L21

Noble, A. G., Muzzin, A., McDonald, M., et al. 2019, ApJ, 870, 56

O'Dea, C. P., Baum, S. A., Privon, G., et al. 2008, ApJ, 681, 1035

O'Sullivan, E., Giacintucci, S., Babul, A., et al. 2012, MNRAS, 424, 2971

Rafferty, D. A., McNamara, B. R., Nulsen, P. E. J., et al. 2006, ApJ, 652, 216

Ruppin, F., McDonald, M., Bleem, L. E., et al. 2021, ApJ, 918, 43

Russell, H. R., Fabian, A. C., Sanders, J. S., et al. 2010, MNRAS, 402, 1561

Russell, H. R., McNamara, B. R., Edge, A. C., et al. 2013, MNRAS, 432, 530

Salpeter, E. E. 1955, ApJ, 121, 161

Sarazin, C. L. 1986, RvMP, 58, 1

Sevilla-Noarbe, I., Bechtol, K., Carrasco Kind, M., et al. 2021, ApJS, 254, 24

Silverman, J. D., Green, P. J., Barkhouse, W. A., et al. 2008, ApJ, 679, 118

Silverman, J. D., Kovač, K., Knobel, C., et al. 2009, ApJ, 695, 171

Somboonpanyakul, T., McDonald, M., Bayliss, M., et al. 2021b, ApJL, 907, L12

Somboonpanyakul, T., McDonald, M., Gaspari, M., et al. 2021a, ApJ, 910, 60

Stern, D., Assef, R. J., Benford, D. J., et al. 2012, ApJ, 753, 30

Stern, D., Djorgovski, S. G., Perley, R. A., et al. 2000, AJ, 119, 1526

Stern, D., Eisenhardt, P., Gorjian, V., et al. 2005, ApJ, 631, 163

Sunyaev, R. A., & Zeldovich, Y. B. 1972, A&A, 20, 189

Voit, G. M. 2005, RvMP, 77, 207

Voit, G. M. 2018, ApJ, 868, 102

Voit, G. M., Donahue, M., Bryan, G. L., et al. 2015, Natur, 519, 203

Waskom, M., Botvinnik, O., Drewokane, et al. 2016, Seaborn: v0.7.1, Zenodo, doi:10.5281/zenodo.54844

Webb, T. M. A., Muzzin, A., Noble, A., et al. 2015, ApJ, 814, 96

Wright, E. L., Eisenhardt, P. R. M., Mainzer, A. K., et al. 2010, AJ, 140, 1868

Yang, L., Tozzi, P., Yu, H., et al. 2018, ApJ, 859, 65

Zwicky, F. 1933, AcHPh, 6, 110