The Kormendy relation of cluster galaxies in PPS regions

We study a sample of 936 early-type galaxies located in 48 low-z regular galaxy clusters with $M_{200}\geq 10^{14}~ M_\odot$ at $z<0.1$. We examine variations in the Kormendy relation (KR) according to their location in the projected phase space (PPS) of the clusters. We have used a combination of Bayesian statistical methods to identify possible differences between the fitted relations. Our results indicate that the overall KR is better fitted when we take into account the information about PPS regions. We also find that objects with time since infall $\geq 6.5$ Gyr have a significant statistical difference of the KR coefficients relative to objects that are more recent in the cluster environment. We show that giant central ellipticals are responsible for tilting the KR relation towards smaller slopes. These galaxies present a late growth probably due to cumulative preprocessing during infall, plus cannibalism and accretion of smaller stripped objects near the center of the clusters.


INTRODUCTION
The Kormendy relation (KR) is a relation between the effective radius   of Early-Type Galaxies (ETGs) and the surface brightness   at that radius (Kormendy 1977).It is a projection of the Fundamental Plane of galaxies, a three-dimensional parameter space correlating, in addition to   and   , the velocity dispersion  of galaxies (Dressler et al. 1987;Djorgovski & Davis 1987;Bender et al. 1992).The KR is a scaling relation indicating that more luminous objects are larger and have a lower characteristic surface brightness.
ETGs encompass objects that normally ceased their star formation, have red colors, small amounts of cold gas and dust, and which correspond morphologically to ellipticals and lenticulars (e.g.Kauffmann et al. 2004;Blanton & Moustakas 2009).Some studies indicate that the formation of massive ETGs has occurred in two phases (e.g.Oser et al. 2010).At an early stage, the gas collapses into dark matter halos and forms stars intensely for a short time interval (e.g.Thomas et al. 2005;Peng et al. 2010;Conroy et al. 2015).The second phase involves mass accumulation through a series of mergers (e.g.Naab et al. 2009;Feldmann et al. 2011;Johansson et al. 2012;Huang et al. 2016) that enrich galaxies with stars set ex-situ, thereby increasing their size and stellar mass continuously (e.g.Trujillo et al. 2007;Van Dokkum et al. 2010).
The structural analysis of ETGs through various scaling relations can be used to probe the evolutionary state of galaxies and constrain galaxy formation models and the assembly history of ETGs (e.g. ★ E-mail: albr@uesc.brTortorelli et al. 2018;Genel et al. 2018;Kuchner et al. 2022).In particular, the KR can be used to derive luminosity and size evolution of ETGs, but it requires homogeneous and complete data samples for clusters at different redshifts, wavelengths, magnitude range, taking into account differences between fitting methods, besides measurement and systematic errors (e.g.Ziegler et al. 1999;La Barbera et al. 2003;Kuchner et al. 2017;Tortorelli et al. 2018).All of this makes it difficult to compare results from different studies.
In this work, we address this problem by studying a homogenous sample of ETGs selected from regular galaxy clusters (Ribeiro et al. 2023).Our aim is to study the KR in subsamples defined in the projected phase space (PPS) of the clusters.The correlation between these regions and the time since infall in clusters is used to probe environmental and evolutionary effects on the coefficients of the KR.A brief description of our data is presented in Section 2. Our analysis and results are presented in Section 3. We summarize and list some conclusions in Section 5.

DATA
Our sample is selected from the extended version of FoF group catalog originally identified by Berlind et al. (2006) and contains 5352 groups with  > 5 and 0.03 ≤  ≤ 0.11, consisting of galaxies with absolute magnitudes   ≤ −20.5, and stellar masses in the range 10.4 < log ( * / ⊙ ) < 11.9, 1 with median ∼ 10 11  ⊙ .The list of members for each group was defined using the "shifting gapper" technique (Fadda et al. 1996;Lopes et al. 2009), extending to ∼ 4 Mpc around the group centers defined by La Barbera et al. (2010).The groups were then subject to the virial analysis, analogous to that described in Girardi et al. (1998); Popesso et al. (2005Popesso et al. ( , 2007)); Biviano et al. (2006); Lopes et al. (2009).This procedure yields estimates of velocity dispersion (  ), radii ( 500 ,  200 ) and masses ( 500 ,  200 ) for most of the groups from the FoF sample.Our final sample contains 107 massive clusters, with at least 20 galaxies within  200 , implying systems with  200 ≳ 10 14 M ⊙ extending to 2 200 .We reduce this sample even further, demanding that the systems do not have substructures according to the DS test (Dressler & Shectman 1988); and with Gaussian velocity distribution, according to the HD metric (see Ribeiro et al. 2013;de Carvalho et al. 2017).We also added morphology from the Domínguez Sánchez et al. (2018) catalog, and the Sérsic indices  in the -band from the catalog of Simard et al. (2011).Our final sample has 936 ETGs as objects with morphological parameter  < 0, probability of being a bulge-dominated object   > 0.5, and Sérsic index  > 2.5, located in 48 regular clusters.For more details, see Ribeiro et al. (2023).

ANALYSIS
We study our stacked sample of cluster ETGs according to their loci in the PPS, which is built by normalizing the projected clustercentric distances and velocities by the virial radius,  200 , and velocity dispersion,   , respectively.To explore the PPS, we divide cluster galaxies into four subsamples, following the classification introduced by Rhee et al. (2017), see also Ribeiro et al. (2023).It is important to consider possible differences between our definitions and those used by Rhee et al. (2017), specifically our use of  200 and their use of   , which are generally not equivalent (see Rhee et al. 2020).We derive a first estimate of the virial mass from equation 5 of Girardi et al. (1998).After applying the surface pressure term correction to the mass estimate, we assume a NFW profile to obtain an estimate of M 200 .This is found after the interpolation (most cases) or extrapolation of the virial mass M  from R  to R 200 (R  is the radial offset of the most distant cluster member).Then, we use the definition of M 200 and a final estimate of R 200 is derived.Further details can be found in Lopes et al. (2009).With this definition,   ≈ 1.36 200 in the ΛCDM cosmology (see Gill et al. 2005;Reiprich et al. 2013), which may introduce a bias to our results.In order to assess for any systematic effects, we performed a series of tests and found that the present analysis remains robust and is not significantly impacted by this difference in the radial normalization of the clusters.Consequently, we have chosen to maintain PPS normalization using  200 .
In Figure 1 we see the distribution of ETG galaxies in the PPS regions of Rhee et al. (2017), following the approximation of Song et al. (2018).All the ancient (171 objects), intermediate (157 objects), recent (307 objects), and first infallers (301 objects) are shown in this figure .A schematic curve illustrates the expected trajectory of a galaxy through the PPS, and therefore indicates that an object in the "ancient" region requires a long time to delve into the potential of the cluster (tipically more than 6.5 Gyr -see e.g.Rhee et al. 2017).

The Kormendy relation
The Kormendy relation (KR) defines an observational correlation between effective radius and surface brightness of ETGs, usually expressed as   =  +  log   (Kormendy 1977).This relation provides information on the distribution of the light profiles and the sizes of ETGs, thus it can be used to study their size-evolution as a function of the environment, in particular the PPS regions.
The effective radii are calculated by where  deV and / are the semimajor axis length and the axis ratio from the de Vaucouleurs fit, respectively.We only use ETGs with / > 0.3 to avoid edge-on objects, and 100 <  < 420 km s −1 (see Ribeiro et al. 2023).These restrictions on  and / reduce the sample to 912 galaxies.The parameters  deV and / in -band are taken from the catalogs PhotObjAll and SpecObjAll of DR15 (Aguado et al. 2019).
We transformed the surface brightness   at   to ⟨  ⟩ using the relation recalling that we now have (e.g.Graham & Colless 1997;Longhetti et al. 2007).This is assumed for the rest of the analysis, but for simplicity we keep the nomenclature .

The KR fitting -all data
We begin our analysis by asking whether classifying the data according to PPS regions influences the result of the linear regression.We first ran the robust linear regression with M-estimator for all ETGs, regardless of the PPS region.Then, we ran the regression once more, now considering the region as a categorical predictor.
Hence, "region" is a categorical variable including the interaction log R e (kpc) effect of PPS loci on the coefficients of the KR relation.The results are presented in Figure 2. We note that both the slopes and intercepts are slightly different between the models:  = 19.62 ± 0.19 and  = 2.94 ± 0.26 (including the PPS regions); or  = 19.24± 0.11 and  = 3.47 ± 0.15 (not including the PPS regions).We compared the models using the Akaike Information Criterion (AIC) difference for the respective likelihoods, finding Δ AIC = 2.59, which provides strong support in favor of the model with influence from the PPS regions.We can assign a probability to the alternative model, without influence, as which provides a relative probability that the alternative model minimizes the AIC.We find p AIC = 0.27, which is low but not entirely convincing.Then we resample the data 1000 times and check how many times the alternative model minimizes the AIC.This procedure results in only 3.2% probability (p boot ) favoring the alternative model, which gives us confidence about the influence of the PPS regions to explain the KR.

The KR fitting in PPS regions
The KR fitting is sensitive to sample selection effects as well as to fitting procedures (e.g.Ziegler et al. 1999;La Barbera et al. 2003;Tortorelli et al. 2018).In particular, ETGs in clusters define the KR with a high intrinsic dispersion, probably due to measurement and systematic errors.(e.g.La Barbera et al. 2003;Nigoche-Netro et al. 2008) Then, we chose to perform a Bayesian linear regression to recover the whole range of inferential solutions, i.e. the posterior distributions of the KR coefficients, rather than a point estimate and a confidence interval as in classical regression.The model is defined as That is, we assign to the parameters a bivariate normal with noise being independent, normally distributed random variables,   ∼ N (0,  2 ).The normal priors are taken from the best fit we found in Section 3.2.From the normal-normal conjugacy we know that the posterior distributions for each regression is also normal distributed (e.g.Lock et al. 2020).The posteriors are presented in Figure 3, with the corresponding bivariate distributions in Figure 4 for 5000 Markov-Chain-Monte-Carlo (MCMC) samples.The differences observed in these figures need to be tested and quantified.For this purpose, we use Bayesian posterior distributions for hypothesis testing.Methods include the Bayes factor two sample t-test (e.g.Morey & Rouder 2011), the Bayesian p-value based on the density at the Maximum A Posteriori (MAP) (e.g.Mills 2018), the region of practical equivalence (ROPE) (e.g.Kelter 2022) and the probability of direction (PD) (e.g.Makowski et al. 2019).These methods are briefly described in Appendix A. In all cases the quantity to be tested is the difference of slopes and intercepts () between two posterior distributions from the Bayesian regressions with respect to different PPS regions.We define the PPS regions as 1 -ancient; 2 -interme- diate; 3 -recent; and 4 -first infall.Accordingly, the differences of the posteriors are written as    , where  and  refer to two different regions of the PPS.
Our results are summarized in Tables 1 and 2, where we notice an almost complete agreement between the methods.Red cells indicate cases in which the tests reject the null hypothesis that slopes and/or intercepts are equivalent.It is clear from the p-values in these tables that the most conservative test is the MAP-based one, and that all tests agree that there is no significant difference between galaxies taken from the intermediate and first infall regions.This result is a little surprising since, considering the time since infall, these are not regions in direct sequence.In fact, the steepest Kormendy relations are precisely those of intermediate and first infall galaxies.Objects in the recent region are in the middle, while the lowest slope comes from ancient galaxies, as we can also see in Figure 5. Being conservative, we assume the results of the MAP test and interpret that the strongest and most relevant difference between the regressions is the one that distinguishes ancient objects from the other subsamples.BCGs could be affecting this result (Ribeiro et al. 2023), but removing them does not change the result [( BCG = 19.25 ± 0.19); ( BCG = 2.95 ± 0.27)].Likewise, removing the lenticulars we found no significant difference in the fittings at the 95% confidence level.This suggests an environmental-evolutionary explanation to understand the particular behavior of the KR for ancient galaxies.

Mass-radius relation
One can study changes in the KR by considering the evolution of luminosity at a certain   (e.g.La Barbera et al. 2003), or through the relative growth of the size of objects in relation to their stellar mass (e.g.Van Dokkum et al. 2010;van der Wel et al. 2014).In this work, we explore this second approach.Galaxy size is a property that is found to significantly vary with galaxy mass, star-formation activity, and redshift (e.g.Shen et al. 2003;Genel et al. 2018).In this context, the mass-radius relationship is of fundamental importance for probing changes in the KR.Usually, it is assumed to be a power law (e.g.Shen et al. 2003;van der Wel et al. 2014), describing the growth of   with the stellar (or total) mass.According to Van Dokkum et al. (2010) the variation of the effective radius   with  * is given by We call this variation  (n), and we used it to investigate possible differences between our subsamples.The idea to be tested is whether the time since infall promotes different behaviors for the ETG in clusters.Using the interpolation of lagged and iterated differences of   and  * we compute the smoothed behaviour of  (n), presented in Figure 6.In this figure, we highlight three points: (i) ancient objects present a significant growth in size with respect the stellar mass; (ii) all the remaining objects present a nearly constant growth, a significant deviation from the power law expectation; and (iii) the Van Dokkum et al. (2010) growth line is steeper than the ancient behaviour.This result clearly indicates a stronger environmental effect on galaxies longer within clusters.According to Hilz et al. (2013), a key factor to understand the inside-out growth scenario is the massratio of merging galaxies, such that for minor mergers, the sizes grow significantly faster and the profile shapes change more rapidly.Ribeiro et al. (2023) found a significant change of both size and shape (given by the Sérsic index) for ancient satellites, suggesting the relevance of possible mergers.
It is important to note that major mergers are supposed not to occur once galaxies have established themselves in the ancient region, an environment unsuitable for such mergings.The process takes place over time across the other regions, possibly through preprocessing in infall groups (see e.g.Lopes et al. 2024).On the other hand, it is known that central galaxies gradually swallows their smaller companions as dynamical friction brings them to the cluster center, growing both in size and mass, while losing velocity dispersion (e.g.Goto 2005).The growth of   with the stellar mass for ancient ETGs seen in Figure 6 aligns with the likelihood of minor mergers and cannibalism, as supported by the aforementioned dynamical friction in the central regions of clusters.In fact, at least 30% of the mass accreted by central galaxies since  ≈ 1 is expected to come from mergers and accretion from stripped satellites (e.g Laporte et al. 2013;Nipoti et al. 2018;Ragone-Figueroa et al. 2018;Chu et al. 2021).This sets a lookback time of ≈7.6 Gyr, which is close to the lower limit of the time since the infall for the ancient objects (Rhee et al. 2017).

DISCUSSION
The Kormendy relation links the effective surface brightness of ETGs to their effective radius.Our findings indicate a relevant role for the environment in establishing this relationship.Basically, we see that objects in the ancient region of the clusters present a tilt towards smaller slopes of the linear fit.To better understand this result, we split the sample of ancient ETGs into BCGs, satellites with low mass ( * < 10 11  ⊙ -LM ANC SAT), and satellites with high mass ( * ≥ 10 11  ⊙ -HM ANC SAT).Ribeiro et al. (2023) verified a change in the behavior of the structural properties of galaxies close to this value, which is also approximately the median stellar mass of the entire sample.The subsample of BCGs is composed of the 48 most luminous ancient cluster ETGs (each taken from a cluster) having   < −21.4 and  * ≥ 10 11  ⊙ (see Ribeiro et al. 2023).The LM and HM ancient satellites samples contain 56 and 67 objects, respectively.We also consider a sample of isolated earlytype galaxies taken from the list of 1-member groups of the catalog of Yang et al. (2007).This field sample, composed of 6670 objects, is defined with the same criteria (except membership) used to select the cluster ETGs (see also Ribeiro et al. 2023).
We compared the linear coefficients of the KR relation for each subsample, following the procedure described in Section 3.3.The confidence ellipses of the parameters from the Bayesian linear linear fit are presented in Figure 7, where the bivariate distributions correspond to 5000 Markov-Chain-Monte-Carlo (MCMC) samples.Using the methods described in Section 3.3 and Appendix A, we find no significant difference between the intercepts of the subsamples.The result of the slope comparison is presented in Table 3. Numbers 1 to 4 indicate the datasets LM SAT, HM SAT, BCGs, and FIELD galaxies, respectively.
The strongest result is the statistical evidence of a smaller slope for HM ancient satellites with respect to the other subsamples, see Figure 7.For HM ANC SAT, ∼90% of objects have  < −2, so they are dominated by massive ellipticals.Our result indicates that this population is likely responsible for tilting the KR relation towards smaller slopes.On the other hand, the comparison for BCGs does not suggest a significant difference in relation to the other subsamples, except the HM ANC SAT.This result does not agree with Samir et al. ( 2020) who find smaller slopes for BCGs in relation to the field.This disagreement may be related to differences in the definition of BCGs, as well as to uncertainties in the brightness profiles of these objects, which are difficult to map due to the presence of many other galaxies close to them.In fact, for mergers occurring at late times, BCGs mainly accrete mass into their extended outskirts, beyond the observational photometry apertures (e.g.Inagaki et al. 2015).The multiple cores in BCGs at  ≲ 0.1 also suggests a complex stellar mass growth of these objects over the past 4.5 Gyr (Hsu et al. 2022) .This growth could be affecting the slope of the KR relation, but our results do not indicate this, either due to problems in determining their structural properties, or due to some mechanism regulating the growth in brightness and size of these galaxies.Differently, the slope variation for the HM ANC SAT is significant, which may be associated with both late accretion in the central region and the history of mergers in infall groups over the cluster evolution (e.g.Lopes et al. 2024).It is also important to mention that field galaxies do not present significantly different coefficients from the intermediate, recent and first infall samples, reinforcing the idea of evolution in ancient galaxies.Therefore, our work is consistent with a scenario where giant ellipticals show a late growth,  ≲ 0.1, either through infall group interactions or through interactions and accretion promoted by dynamic friction in the central region of the clusters.The combination of these evolutionary channels produce a significant effect on the Kormendy relation, indicating that it can be used as an important indicator of the late evolution of cluster galaxies.

CONCLUSIONS
The study of ETGs is important to improve our understanding of the structure formation and the galaxy-environment connection.In this work, we investigated the KR of cluster ETGs according to their loci in the PPS.We have used a combination of statistical methods to identify possible differences between the fitted KR.From our analysis we can conclude that: • The KR is is better fitted when we take into account information about the PPS regions.
• There is a significant statistical difference between the KR coefficients of the ancient ETGs in relation to the others.
• The longer time within host clusters and the cumulative effect of mergers (in infall groups and/or in the central regions due to dynamical friction) plus accretion of stripped galaxies are the likely causes for the difference between the KR through the PPS regions.The relative importance of these mechanisms deserves a further study, to be presented in a forthcoming paper.
• HM ancient satellites are the component that most contributes to tilting the KR relation.
• A better determination of the brightness profile of BCGs and a larger sample of these objects can help to better understand the slope of the KR relation in the ancient region.
A possible caveat to this work is that we have not considered the effects of the velocity dispersion on the intrinsic scattering of the KR, a point to be examined in the future, in a complete study of the fundamental plane.

A3 Region of Practical Equivalence (ROPE)
The general idea of this test is to establish a region of practical equivalence around the null value  0 of the null hypothesis  0 :  =  0 , which express the range of parameter values that are equivalent to the null value for practical purposes (Kruschke & Liddell 2018).To test hypotheses via the ROPE we need do find the highest posterior density (HPD) interval, which is basically the shortest interval on a posterior density for some given confidence level.Kruschke (2018) proposed the following decision rule: • Reject the null value  0 specified by  0 :  =  0 ( = 0), if the 95% HPD falls entirely outside de ROPE.
• Accept the null value, if the 95% HPD falls enterily inside the ROPE.
Hence, if the 95% HPD falls entirely inside the ROPE, the parameter value is located inside the ROPE with at least 95% posterior probability.As a consequence, it is practically equivalent to the null value  0 and it is reasonable to accept  0 (see Makowski et al. 2019).

A4 Probability of Direction (PD)
The probability of direction (PD) is dvan2010growthefined as the proportion of the posterior distribution that is of the median's sign,is a measure of effect existence representing the certainty with which an effect is positive or negative.Therefore, the PD is simply the proportion of the posterior probability density which is of the median's sign.Based on the PD, one can test the null hypothesis  0 :  = 0 by requiring a specified amount of posterior probability density to be strictly positive or negative.Usually, more than 95% of the posterior can be used as threshold for deciding between  0 and  1 .The two-side p-value is obtained as (A4) (see Makowski et al. 2019).
This paper has been typeset from a T E X/L A T E X file prepared by the author.For the likelihood ratio, points in red are showed, and for the computation of the PD area, the arrow indicates the lower limit of the density distribution integral, which extends to −∞.

Figure 1 .
Figure 1.Distributions of ETGs in different PPS regions.Ancients objects are in red, intermediate in orange, recent in green, and the first infallers are in blue.A schematic curve illustrates the expected trajectory of a galaxy through the PPS.

Figure 2 .
Figure 2. Linear regressions including (magenta) or not (cyan) the region as an interaction term on the slope.

Figure 3 .Figure 4 .
Figure 3. Posterior regression parameters from the Bayesian fits.Ancients objects are in red, intermediate in orange, recent in green, and the first infallers are in blue.

Figure A1 .
Figure A1. Figure illustrating the most important aspects of the Bayesian methods used in this work.The gray shade area shows the ROPE, while the blue shaded area corresponds to the posterior probability inside the 95% HPD.For the likelihood ratio, points in red are showed, and for the computation of the PD area, the arrow indicates the lower limit of the density distribution integral, which extends to −∞.

Table 1 .
Results of the Bayesian tests for intercepts .Significant differences are highlighted in red.

Table 2 .
Results of the Bayesian tests for slopes .Significant differences are highlighted in red.

Table 3 .
Results of the Bayesian tests for slopes  of subsamples of LM SAT, HM SAT, BCGs, and FIELD galaxies.