On the Origin of the Variety of Velocity Dispersion Profiles of Galaxies

Observed and simulated galaxies exhibit a significant variation in their velocity dispersion profiles. We examine the inner and outer slopes of stellar velocity dispersion profiles using integral field spectroscopy data from two surveys, SAMI (for z < 0.115) and CALIFA (for z < 0.03), comparing them with results from two cosmological hydrodynamic simulations: Horizon-AGN (for z = 0.017) and NewHorizon (for z ≲ 1). The simulated galaxies closely reproduce the variety of velocity dispersion slopes and stellar mass dependence of both inner and outer radii (0.5 r 50 and 3 r 50) as observed, where r 50 stands for half-light radius. The inner slopes are mainly influenced by the relative radial distribution of the young and old stars formed in situ: a younger center shows a flatter inner profile. The presence of accreted (ex situ) stars has two effects on the velocity dispersion profiles. First, because they are more dispersed in spatial and velocity distributions compared to in situ formed stars, it increases the outer slope of the velocity dispersion profile. It also causes the velocity anisotropy to be more radial. More massive galaxies have a higher fraction of stars formed ex situ and hence show a higher slope in outer velocity dispersion profile and a higher degree of radial anisotropy. The diversity in the outer velocity dispersion profiles reflects the diverse assembly histories among galaxies.


INTRODUCTION
Most galaxies in the Universe are expected to exist inside a dark matter halo that extends far beyond the size of the galaxy (e.g., Rubin et al. 1980;Bosma 1981).In the standard cold dark matter (CDM) paradigm, dark matter halos are assumed to be affected by various gravitational processes.Therefore, their structural properties are expected to reflect their assembly and evolutionary history (Navarro et al. 1996;Wechsler et al. 2002).For example, brightest cluster galaxies (BCGs) and brightest group galaxies (BGGs) located at the center of a large system are surrounded by exceptionally massive and extended dark matter halos (Moster et al. 2010;Pillepich et al. 2018a) presumably as a result of numer-ous mergers and accretion events during their formation (Gao et al. 2004;De Lucia & Blaizot 2007).However, satellites within such systems are prone to tidal stripping which preferentially removes loosely-bound outer parts of the dark halos, resulting in truncated dark matter halos (Smith et al. 2016).
It has been also claimed that the inner dark matter distributions of dwarf galaxies show flat profiles at the center, deviating from what is expected from cold dark matter properties (Flores & Primack 1994;Simon et al. 2005;Oh et al. 2015).These dark matter profiles with flat centers have revealed the necessity for additional physical processes beyond the simple CDM-based modeling of the Universe.The effect of baryonic feedback has been proposed as a possible solution.When the feedback is sufficiently strong, it can generate fluctuations in the gravitational potential of the halo, leading to the flattening of the inner dark matter profiles (e.g., Peirani et al. 2008;Pontzen & Governato 2012;Teyssier et al. 2013).An alternative scenario that has been suggested involves the self-interaction of dark matter.In this scenario, the collision of dark matter particles within the halo leads to the heating of the central cusp, transforming it into a shallower profile (Spergel & Steinhardt 2000;Rocha et al. 2013;Tulin & Yu 2018).
Therefore, inferring the dark matter profiles of galaxies is useful to understand the evolution of dark halos and their relationship with galaxies, as well as for probing the physical properties of dark matter.These processes serve as essential tests for the standard cosmological models.
Recent developments in the integral field unit (IFU) technique have enabled the detailed study of the kinematics of galaxies in large quantities and great detail (Sánchez et al. 2012;Ma et al. 2014;Bryant et al. 2015).Specifically, Santucci et al. (2022) used the Sydney-AAO Multi-object Integral-field unit system (SAMI, Croom et al. 2012;Bryant et al. 2015;Croom et al. 2021) data to derive dark matter fractions in passive galaxies.Recent studies suggest that galaxies display a variety of velocity dispersion profiles (Neumann et al. 2017;Falcón-Barroso et al. 2017;Pulsoni et al. 2018;Loubser et al. 2018;Veale et al. 2018;Mogotsi & Romeo 2019;Lu et al. 2020;Edwards et al. 2020).The velocity dispersion profile is found to be dependent on the morphology and bulge-tototal ratio of the galaxy (Neumann et al. 2017) and tends to vary more significantly at larger radii (Pulsoni et al. 2018).However, the degeneracy between the orbital velocity anisotropy and the total mass profile makes it difficult to reconstruct the dark matter distribution (Binney & Mamon 1982).It is known that this issue can be overcome by using the high-order Gauss-Hermite moment of the line-of-sight velocity distribution h 4 (Gerhard 1993;Read & Steger 2017) that is coupled with the velocity anisotropy of the system.However, there have been reports of certain massive early-type galaxies exhibiting rising velocity dispersion profiles (Edwards et al. 2020) and positive values of h 4 (Krajnović et al. 2008;Veale et al. 2017;Loubser et al. 2018Loubser et al. , 2020).This appears to be contradictory because rising profiles are usually associated with tangentially-biased anisotropy, whereas positive values of h 4 are considered hints of radially biased anisotropy.As a possible solution to this problem, Veale et al. (2018) proposed a variation in the total mass profile, and Loubser et al. (2020) suggested the contribution of intracluster light.
The degeneracy is mainly caused by projection effects, which makes it difficult to interpret the observed velocity dispersion profiles of the galaxies correctly.Numerical simulations are useful for investigating this is-sue because they provide detailed kinematic and spatial information on all the particles.Cosmological simulations provide tests for the ΛCDM cosmology based on the initial conditions of primordial density fluctuations that are observable in cosmic microwave background data.Advances in computational resources and simulation codes have enabled the emergence of large, high-resolution cosmological simulations with sophisticated subgrid prescriptions (e.g., Tremmel et al. 2019;Pillepich et al. 2019;Dubois et al. 2021).Recent numerical simulations have successfully reproduced various observed trends in velocity dispersion profiles (e.g., Pulsoni et al. 2020;Lu et al. 2020;Wang et al. 2022;Cannarozzo et al. 2023).However, the full understanding of the physical processes that shape galaxy velocity dispersion remains limited.Variations in the kinematic properties of galaxies are believed to result from the interplay between various processes in the hierarchical assembly paradigm.Understanding how these processes are encoded in the kinematic properties is crucial.Conversely, velocity profiles can be used to reconstruct the past assembly and evolution histories of galaxies.
Stars in galaxies typically assemble in two ways: insitu and ex-situ star formation.Two types of stellar components exhibit distinct kinematic characteristics.In-situ formed stars, or simply in-situ stars, exhibit lowvelocity dispersion as they retain dynamically cold properties inherited from the viscous nature of the cold gas component from which they originate.Ex-situ formed stars, or ex-situ stars, are defined as stars that originate from external galaxies and become incorporated into the main galaxy through mergers or accretion events.During the coalescence process, ex-situ stars lose information about the dynamical properties they had in their original galaxy.Understanding the distinction between these two components is crucial for comprehending how the velocity dispersion profiles of galaxies are shaped.
This study adopts two cosmological hydrodynamic simulations to measure the velocity dispersion profiles of galaxies using a methodology similar to that of IFU surveys.By tracing the history of star particles, we aimed to determine the main processes that shaped the current forms of velocity dispersion profiles.In Section 2, we describe the sample selection process used for simulations and observations.In Section 3, we compare the slopes of the velocity dispersion profiles obtained from simulations and observations.In Section 4, we use the Jeans' equation to perform a kinematic analysis and determine the key structural parameters that influence the velocity dispersion profiles of galaxies.In Section 5.1, we discuss the primary factors and processes that determine the velocity dispersion profiles of galaxies.This study assumes the standard ΛCDM cosmology based on Komatsu et al. (2011) (h=0.704,Ω m =0.272, Ω Λ =0.728) and the stellar initial mass function of Chabrier (2003).
Horizon-AGN (hereafter HAGN) involves a periodic cube with a comoving volume of (142 Mpc) 3 and has a "best" spatial resolution of ∼ 1 kpc.The vast scale of the simulation allows access to numerous galaxy samples from various environments.However, the limited gravitational and hydrodynamic spatial resolutions of HAGN, which are comparable to the typical value of the effective radii of dwarf galaxies, result in poor reproduction of sub-kpc scale structures such as thin disks.Additionally, the mass resolution of ∼ 3.5 × 10 6 M ⊙ limits the accurate representation of galaxy kinematics, leading to significant statistical noise in the analysis, particularly when measuring spatially resolved kinematics.Because these resolution issues become more severe for smaller galaxies, we limit the stellar mass range of our sample to log(M * /M ⊙ ) > 10.75.This corresponds to more than 15,000 star particles in each galaxy.We aim to obtain a statistically reliable number of star particles for measuring the velocity dispersion in 100 equally numbered radial (shell) bins around a galaxy, which corresponds to ≥ 150 particles in each bin.As a result, we managed to include 5,829 galaxies from the HAGN simulation at the final snapshot of z = 0.017, which includes 40 galaxies with log(M * /M ⊙ ) > 12 and 405 galaxies with log(M * /M ⊙ ) > 11.5.
NewHorizon (hereafter NH), which has the best spatial resolution of ∼ 34 pc, consists of a spherical zoom-in region with a comoving diameter of 20 Mpc.We utilize all the snapshots down to z = 0.017, while the introductory paper by Dubois et al. (2021) for NH used only down to z = 0.25.Because of the use of a zoom-in technique in NH, there is a possibility that some galaxies, usually near the boundary regions of the sphere, contain dark matter particles from low-resolution regions in the initial conditions.This can increase the shot noise in density, leading to a negative effect on the gravitational stability of the system.To minimize the contamination effect, we use galaxies without lowresolution dark matter particles for our sample.To include as many galaxies as possible in our analysis, we extract the galaxies from multiple snapshots in the NH simulation by selecting time intervals of approximately 0.5 Gyr.In total, we selected 12 snapshots in the redshift range of 0.17 < z < 0.96.We used a mass cut of log(M * /M ⊙ ) > 9.5 for the NH sample for easy comparison with the observed data that is described in Section 2.2.This means that, with the typical masses of star particles of ∼ 10 4 M ⊙ , each galaxy can be resolved into a minimum of 300,000 star particles.The sample includes 2,104 galaxies in total, with 278 galaxies having log(M * /M ⊙ ) > 10, and 86 galaxies having log(M * /M ⊙ ) > 10.5.While most of our galaxies reside in the field environments, 126 galaxies are located in two group-size halos of the masses of log(M vir /M ⊙ ) = 12.83 and log(M vir /M ⊙ ) = 12.96.In both simulations, we use the AdaptaHOP algorithm (Aubert et al. 2004) to detect galaxies based on the density distribution of star particles.The evolutionary track of each galaxy is traced by following the main progenitor branch of the merger tree.
NH and HAGN utilize the adaptive mesh refinement (AMR) technique implemented in the RAMSES code, where the hydro and gravitational solvers operate on octree grids.With each eightfold increase in density, the grid is subdivided into smaller grids.This approach enables RAMSES to assign smaller grids in higher-density regions to conduct more precise calculations.In NH, the supplementary refinement criterion is triggered when the Jeans length of gas becomes smaller than four times the size of the cell.The simulation stops refining grids when it reaches the smallest grid size, i.e., the best spatial resolution.
All stars in simulated galaxies are classified as either "ex-situ" or "in-situ" stars.Ex-situ stars are defined as stars that were older than 250 Myr at the time of accretion.The accretion is defined as the moment they were recognized as members of the hosting galaxy by the galaxy finder for the first time.The utilization of a time threshold helps exclude transient stellar systems within the galaxy, which may sometimes be identified as separate systems through the galaxy finder algorithm.The choice of a 250 Myr time threshold is comparable to the orbital period of the galaxy, ensuring sufficient time for the dissolution of clumps through differential rotation.The ex-situ fraction includes both minor and major mergers.Of the two, major mergers are considered to have a more dramatic effect on the velocity dispersion and even morphology of the remnant galaxy (e.g., Bournaud et al. 2005;Hilz et al. 2012).However, we do not distinguish them when identifying ex-situ stars because they have largely the same tendency to increase velocity dispersion, a key property of our investigation.
Finally, a mock line-of-sight velocity dispersion (σ los ) profile is generated for each galaxy using stars grouped

Observation data
This section describes the observational data collected from the SAMI and Calar Alto Legacy Integral Field Area (CALIFA) surveys that are used to compare gradients of the stellar velocity dispersion.We employed the third release data from the SAMI survey (Croom et al. 2021).The SAMI (Croom et al. 2012;Bryant et al. 2015) employs 13 fused optical fiber bundles (hexabundle), each containing 61 1.6 ′′ diameter fibers (Bland-Hawthorn et al. 2011;Bryant et al. 2014), feeding the AAOmega dual-arm spectrograph mounted on the Anglo-Australian Telescope (Sharp et al. 2006).The blue and red arms use 580V and 1000R gratings, covering, 3750-5750 Å at a resolution of R=1808, and 6300-7400 Å at a resolution of R=4304, respectively.The SAMI survey includes more than 3,000 galaxies with a stellar mass range of log(M * /M ⊙ ) = 8-12 and a redshift range of 0.004 < z < 0.115 (Bryant et al. 2015;Croom et al. 2021).More than two-thirds of the SAMI samples are from three equatorial fields (G09, G12 and G15) of the Galaxy And Mass Assembly survey (Driver et al. 2011).In addition, galaxies from eight clusters have been observed to complete the environmental matrix (Owers et al. 2017).We use spatially-resolved line-ofsight stellar velocity dispersions published by the SAMI team (Croom et al. 2021).van de Sande et al. (2017) described the measurement of SAMI stellar kinematics using Penalised piXel-Fitting software (pPXP; Cappellari & Emsellem 2004;Cappellari 2017).We select 2,146 galaxies with stellar mass log(M * /M ⊙ ) > 9.5 to reliably estimate the velocity dispersion gradient.We use the stellar masses provided by the SAMI catalog (Bryant et al. 2015), derived from i-band magnitudes and g-i colors.
In addition to the SAMI and CALIFA survey data, we use published data from two other sources.Veale et al. (2018) examined 90 early-type galaxies at the fixed radii of 2 kpc and 20 kpc using the MASSIVE IFU survey.We also use the line-of-sight velocity dispersion (σ los ) slopes and K-band magnitudes derived through long-slit spectroscopy of BCGs and BGGs from Loubser et al. (2018, Figure 5).For these catalogs, we employ the stellar mass estimated from K-band magnitude using the relation provided by Cappellari (2013).
We note that in our main results, we did not apply any morphological classifications on the simulated galaxies.For observations, SAMI and CALIFA encompass all morphological types of galaxies, whereas other sources, such as Veale et al. (2018) and Loubser et al. (2018), include only early-type galaxies.This bias is somewhat alleviated at the high-mass end, where both the SAMI and CALIFA samples are dominated by early-type galaxies.Although galaxy morphology has not been directly employed in our main results, discussions on the dependence on morphology are included in Appendix B. Table 1 summarizes the overall properties of simulation and observation data, presenting the total number of galaxies, redshift range, stellar mass range, spatial resolution, and type of galaxies.This section describes the post-processing procedures conducted on the simulation galaxies for producing 2D velocity dispersion maps for comparison with the IFUobserved data.In the simulations, a spherical boundary is defined for each galaxy with a radius where the projected surface brightness drops below 26.5 mag/arcsec 2 .This is analogous to the conventional cut for the boundary of galaxies in photometry.The luminosity of each star particle is evaluated based on the population synthesis model of Bruzual & Charlot (2003), in combination with the stellar initial mass function of Chabrier (2003).The radius is measured as the average of three different line-of-sight directions (x, y, and z).
For each galaxy, we measure the stellar line-of-sight velocity dispersion projected onto a 2D plane.First, 24 random directions were selected for each galaxy.In each projected view, the star particles are binned into segments of Voronoi tessellation.We design the Voronoi cells to contain an approximately equal number of star particles, following Cappellari & Copin (2003).The effective number of stars per cell is 5,000 for NH and 100 for HAGN.The former and latter correspond to the stellar masses of the ∼ 5 × 10 7 M ⊙ and ∼ 3.5 × 10 8 M ⊙ , respectively, which are comparable to the typical sizes of the SAMI spaxels.For each cell, the r-band fluxweighted standard deviation of the stellar line-of-sight velocities is measured.It is worth noting that the velocity dispersion we examine is consistent with that mea- sured in IFU observations, which differs from the traditional term of that derived from aperture spectroscopy.
The former measures the standard deviation of stellar velocity in a localized region, while the latter measures the integrated second velocity moment across the entire galaxy, considering spatial gradients in the mean velocity as well.

Measuring the σ slope
To measure the radial velocity dispersion profile for each projected view of observed and simulated galaxies, an effective ellipse is computed to enclose half of the total flux of the galaxy.The Kinemetry method (Krajnović et al. 2006) is employed to determine the effective ellipse that follows the isophotal lines.The log r-band flux is assumed to be even moments (n = 0, 2, 4) of the harmonic function.While maintaining the orientation and ellipticity of the effective ellipse, from 20 to 50 concentric ellipse bins with different sizes are selected.The radial velocity dispersion profiles are then measured as a function of the semi-major axes of the ellipse bins.
For the systematic measurement of the velocity dispersion profile, we adopt the broken power-law fit of Veale et al. (2018), where r represents the semi-major axis of the ellipse bin, γ 1 and γ 2 represent the inner and outer asymptotic power-law gradients of the fitted curve, respectively, and r b is the break radius, and σ 0 is the normalization of the profile (value of σ f at r = r b ).The fit is performed using optimize.curvefit function in scipy module, by minimizing the chi-square in the log-log plane using four parameters, γ 1 , γ 2 , r b , and σ 0 .To reduce potential issues arising from beam-smearing in observations and spatial resolution in simulations, we exclude the data points within a radius of 0.1r 50 , where r 50 denotes the semi-major axis of the effective ellipse1 The slope of the fitting function in the log-log plane, is measured at two radial points.The inner slope is defined as the slope at 0.5 r 50 , denoted as γ in = γ(0.5 r 50 ), and the outer slope is defined as the slope at the 3 r 50 , denoted as γ out = γ(3 r 50 ).
We use inner and outer slopes measured from Veale et al. (2018).Because they presented slopes only at 2 and 20 kpc, while the break radius, r b , was fixed at 5 kpc, the broken power-law curve can be recovered by using the following equations, where γ i and γ o denote slopes measured at 2 and 20 kpc, respectively.After recovering γ 1 and γ 2 , we re-measured slopes at 0.5 r 50 and 3 r 50 to ensure consistency with our definitions of the inner and outer radii.We derived the half-light radius of each galaxy from the published catalog of Greene et al. (2019).
We also utilize the slopes of velocity dispersion profiles presented in Loubser et al. (2018), Figure 5.Given that the slope has been measured from a single power-law fit, we use the same values for both the inner and outer radii in our analysis.Figure 3 demonstrates the variation in the σ los profiles among elliptical galaxies.The three sample galaxies extracted from the CALIFA survey exhibit flat, rising, and falling profiles.Even for similar morphologies and masses, the profiles are significantly different.
We select the final sample from the broken power-law fitted data based on the reduced chi-square and mean σ/error ratio, where "error" means the error in σ los .This secures the reliability of the data and avoids over or under-fitted samples due to poor measurement of the data or complex behaviors of the σ los profile caused by interlopers, mergers, and close companions.Detailed procedures and criteria are described in Appendix A. The final dataset comprises 2,424 galaxies, with 2,101 Figure 4 shows the inner slope, γ in , measured in simulations and observations as a function of stellar mass.The black and green lines represent galaxies from NH and HAGN, respectively.Observational data are also presented.The most notable feature of this diagram is the "V" shape.The inner slope of the σ los profile becomes steeper (more negative slope) with increasing stellar mass until log(M * /M ⊙ ) ∼ 11, after which it starts increasing.This fall-and-rise trend appears in both the simulated and observed samples.The most representative samples in this diagram are the SAMI and CALIFA data because they contain a large number of galaxies.They both show a "V" trend.
The "V" trend is visible in the simulations as well.NH does not show the upturn above log(M * /M ⊙ ) ∼ 11 because the field-environment simulation is almost exclusively composed of late-type galaxies.If only late-type galaxies are selected from the observation, the trend shows a monotonic decrease with stellar mass down to γ in ∼ −0.2, which is in good agreement with NH.The resulting figure is shown in Figure B1 turn trend in the high mass range.In contrast, the kinematic structure of its low-mass galaxies below our mass cut (log(M * /M ⊙ ) < 10.75) may not be suitable for this analysis, as we explained in Section 2.1; therefore, we excluded them from the diagram.If we accept the low-mass galaxies from the high-resolution NH simulation and massive galaxies from the HAGN simulation as a combination, we can reproduce the observed fall-andrise trend.The origins of this trend are discussed in the following section.
One may find the offset of the observed data from MASSIVE uncomfortably large.The vertical offset may have originated from a horizontal offset, i.e., in stellar mass estimates for observed galaxies.Alternatively, the HAGN simulation may be incorrect by as much as the vertical offset, which is not inconceivable considering its relatively poor (1 kpc) resolution.On the other hand, we use only the most massive galaxies from HAGN, exactly worrying about the issue, and the offset persists even for the most massive galaxies for which the signals (or calculations) should be more reliable than for less massive galaxies.
Figure 5 shows the outer slope, γ out , as a function of stellar mass.The black and green lines represent the projected outer slopes and the 1σ scatters for NH and HAGN galaxies, respectively.The observed data is also presented.We show the combined data of SAMI and CALIFA in this figure because SAMI and CALIFA by themselves show extremely large scatters in outer slopes.The number of galaxies in the drawing sample is reduced to 322 (190 from SAMI and 132 from CALIFA) after selecting observations that consist of data reaching up to 3 r 50 from the center.The large scatters may be indicative of a large variation in the physical properties in the outskirts of galaxies.In the mass range of log(M * /M ⊙ ) ≲ 11, NH exhibits no distinct trend.SAMI and CALIFA observations show a weak decreasing trend, but the errors are so large that the trend may not be statistically relevant.In addition, the "trend" might have originated from the continuation of the broken power-law fits that are dominated (in terms of chisquare evaluation) by the brighter inner profiles.More massive galaxies (log(M * /M ⊙ ) ≳ 11) show a hint of an increasing trend with stellar mass, and it is consistently visible for both observed and simulated galaxies.We will discuss the origin of this trend in the following sections.

σ profiles and assembly history
In the previous section, we used 2D properties to compare the simulations with observations and demonstrated their similar behavior with respect to stellar mass.However, projection effects cause complications in analysis and interpretation.To measure the intrinsic kinematic properties of the galaxies in 3-D space, we set 100 spherical bins (shells) centered on each simulated galaxy with equally numbered star particles.The luminosity-weighted velocity dispersion in each shell is measured as where σ r , σ θ , and σ ϕ are the velocity dispersion of stars in three directions in spherical coordinates.
To figure out the origin of the different shapes in the velocity dispersion profile, it is important to trace the assembly and formation history of the stellar components of galaxies.For this purpose, we plot Figures 6 and  7, which show the velocity dispersion profiles and mass fractions of stellar components from different origins.In Figure 6, we plot the σ 3D profiles (panels (a) and (b)) and mass fractions (panels (c) and (d)) of the in-situ and ex-situ stars as a function of the radial distance in the NH galaxies.The in-situ stars are divided into young and old populations based on an arbitrary age cut of 3 Gyr.The left and right columns represent the two different ranges of stellar masses.Ex-situ stars generally exhibit higher velocity dispersion at all radii and thus increase the mean velocity dispersion in the region.Because ex-situ stars are more predominant in the outer regions of the galaxy (Rodriguez-Gomez et al. 2016), the strongest impact occurs in the outer regions.This  is consistent with the formation of dynamically hot stellar components via accretion (Abadi et al. 2006;Dubois et al. 2016;Park et al. 2019;Zhu et al. 2022).
In contrast, in-situ stars show relatively lower velocity dispersion overall (Dubois et al. 2016).When comparing samples of different ages, older stars exhibit higher values of σ 3D .This result can be interpreted in the context of dynamical heating (Quinn et al. 1993;Yi et al. 2024).Newly formed stars are likely to retain dissipative gas dynamics with low σ 3D ; however, their velocity dispersion is increased with time due to continuous gravitational perturbations from internal or external sources such as spiral arms, giant molecular clouds, mergers, and galaxy encounters.Consequently, older stars tend to exhibit dynamically hotter kinematics (Wielen 1977;Aumer et al. 2016;Park et al. 2021;Sharma et al. 2021).
Panels (c) and (d) show that old and young stars have different radial distributions for the two samples with different masses.In low-mass galaxies (panel (c)), the inner region is increasingly occupied by young stars, reducing σ 3D in the inner region, and the fraction of old stars increases with increasing radial distance, raising σ 3D gradually.As a result, the σ 3D slope of the inner region is flattened.
In the more massive galaxies (panel (d)), the trend is reversed: the inner region is dominated by old in-situ stars, which increases the central σ 3D , and the fraction of young stars remain relatively stable (or slightly increases) with the radial distance, leading to a steeply falling σ 3D profile.This may indicate a transition in the order of the star formation process from outsidein to inside-out with increasing stellar mass (Pan et al. 2015;Lin et al. 2019).Therefore, the inner σ 3D profiles are determined by the relative radial distributions of the young and old stellar populations.This can also be seen in mock IFU images.Figure 1 shows a massive galaxy containing an old core with high σ los , surrounded by a young disk with low σ los (panel (e)), exhibiting a steep inner gradient in σ los map (panel (c)).In contrast, Figure 2 shows a lower-mass galaxy containing young in-situ stars in the central region and age increases with radial distance (panel (e)), exhibiting a shallow inner gradient in σ los map (panel (c)).
Figure 7 is the counterpart of Figure 6 using HAGN galaxies and different mass ranges.The relative importance between the three components, as shown in panels (a) and (b), is similar to that of NH, with the ex-situ stellar component having the highest and young in-situ stars having the lowest.Similar to NH, the fraction of ex-situ stars increases with radial distance in panels (c) and (d).Flattened σ 3D profiles are observed in more massive galaxies (panel (b)), primarily following the dis- tribution of ex-situ stars that dominate the composition of the galaxy.For the in-situ components, the relative fraction of old stars over young stars increases with increasing radial distance, which is understandable in the context of dynamical heating and migration.This is in contrast to the high-mass sample of NH shown in panel (d) of Figure 6, which has a similar stellar mass range.This discrepancy has two possible explanations.First, the morphologies of high-mass end galaxies in NH are dominated by late-type galaxies which contain an old bulge with a star-forming disk.In contrast, HAGN includes a substantial number of early-type galaxies in this mass regime, which likely have central star formation rather than extended disk-mode star formation.Second, HAGN lacks sufficient resolution to accurately distinguish between bulge and disk components, particularly for lower-mass galaxies (log(M * /M ⊙ ) ≲ 11).
To quantify the effect of external accretion, we define ex-situ fraction as the mass fraction of ex-situ stars.Figure 8 shows the impact of ex-situ fraction on the outer slope in the same format as Figure 5, based on HAGN galaxies.Panels (a) and (b) show the slopes of the intrinsic (σ 3D ) and projected (σ los ) profiles, respectively.The same broken power-law fit and slope measurement described in Equations 1 and 2 are applied to σ 3D by using r as the radial distance from the center.Panel (b) also presents the observed data.As shown in Figure 5, a mild positive correlation exists between γ out and the stellar mass in Panel (b).The mass dependence is even more evident in the intrinsic 3D case, as shown in Panel (a).The mass trend is related to the impact of ex-situ fraction.In Panel (b), it is apparent that more massive galaxies have higher ex-situ fraction values.This is expected in the hierarchical paradigm, in which more massive galaxies are likely to have experienced a larger number of mergers and accretions (Lee & Yi 2013;Dubois et al. 2016;Rodriguez-Gomez et al. 2016;Davison et al. 2020;Remus & Forbes 2022).The apparent correlation between γ out and ex-situ fraction can also be explained by the fact that ex-situ stars are more spatially dispersed and kinematically hotter than in-situ stars (as shown in Figure 6).This dominance of ex-situ stars at the outskirts leads to an increase in γ out towards zero.The presence of a vertical gradient, although relatively less pronounced in Panel (b) due to the projection effect, suggests that the variation of exsitu fraction among galaxies causes the diversity of σ profiles of outer halos.

KINEMATIC ANALYSIS OF THE VELOCITY DISPERSION PROFILE
This section presents the investigations of the structural properties of galaxies that govern the shape of the σ profile.To address this in a quantitative way, we employ the Jeans' equation in spherical coordinates (Binney & Tremaine 2008, 4.215), where ρ * is the stellar density of the system, v 2 r is the second moment of velocity in the radial direction, Φ is the gravitational potential of the system, and β is the anisotropy parameter defined as Here, v 2 r , v 2 θ , and v 2 ϕ represent the second moment of velocities of stars in spherical coordinates.It is worth noting that our definition of β considers not only the anisotropy of the velocity dispersion but also the bulk tangential velocity of the system (i.e., rotation).This implies that β determines whether the system is supported tangentially or radially.
We applied an assumption that the system has an equilibrium of inflow and outflow in the radial direction, v r = 0, so that v 2 r = σ 2 r , and the coordinate axis is aligned to the direction of the total angular momentum v θ = 0, so that v 2 θ = σ 2 θ .However, this assumption cannot be applied to the meridional direction if the system has non-zero angular momentum, leading to v ϕ ̸ = 0.By converting the differential operators to log scale, Equation 5 can be rewritten as Equation 7 indicates that the radial velocity dispersion profile σ r (r) is primarily correlated with three parameters: the circular velocity profile v c (r) = GM (< r)/r, the stellar density gradient α * (r) = −d(ln ρ * )/d(ln r) and the anisotropy profile β(r).The σ r gradient γ r (r) = d(ln σ r )/d(ln r) depends on the σ r profile itself, and its change along the radial distance is typically small (±0.3 at maximum).Thus, we do not consider it as one of the important parameters that describe σ r (r).Accordingly, the role of each parameter affecting σ profiles can be summarized as follows.
1.The circular velocity, v c (r), directly correlates with σ r (r) and is a function of the enclosed mass profile, which can also be represented by the power-law slope of the total density profile α(r) = −d(ln ρ)/d(ln r).
2. The stellar density power-law slope, α * (r), has a negative correlation with the σ r (r) profile.
3. The velocity anisotropy, β(r), quantifies the dominance of either radial or tangential stellar motion and is positively correlated to the σ r (r) profile.When the anisotropy increases (i.e., the velocity becomes radially biased), the velocity dispersion increases.
In Section 3.3, we observed that the slope of velocity dispersion correlates with the stellar mass and ex-situ fraction.Here, we demonstrate how the 3 parameters mentioned above change with stellar mass and ex-situ fraction and present their significance to the velocity dispersion profile.
Figure 9 shows the radial profiles of the radialcomponent velocity dispersion (σ r ) in the first row, and the profiles of three parameters (α, α * , and β) in the second-fourth rows, as a function of radial distance.We present the total density profile slope, α, instead of the circular velocity, v c , for the reasons mentioned above.Each line with a different color represents the median profile of the galaxies subsampled based on their stellar mass.The left and right columns correspond to NH and HAGN, respectively.In the first row, the σ r profile of the NH galaxies exhibits significant mass dependence on the inner slope.In low-mass galaxies, the inner velocity dispersion exhibits a flattened profile.However, as the stellar mass increases, the slope of the velocity dispersion profile transitions to a steeper, falling profile, indicating a more rapid increase in the velocity dispersion toward the inner regions of the galaxy.
In contrast, the HAGN galaxies, which cover a higher mass range (log(M * /M ⊙ ) > 10.75), show a clear transition of the slope with stellar mass from a negative, falling slope to a positive, rising slope.The mass dependence of the radial σ r profile (first row: panels a and e) is consistent with the trend of the projected velocity dispersion slope, as shown in Figures 4 and 5.
The second and third rows present the power-law density slopes for the total and stellar masses of the galaxies, respectively, categorized by different stellar masses.In NH, the density profile exhibits a significant difference in the inner region.More massive galaxies maintain a steeper central density slope.In contrast to the NH galaxies in Panel (b), the HAGN galaxies in Panel (f) exhibit more distinct mass dependence on the outskirts.
More massive galaxies exhibit a shallower slope, indicating a more extended mass distribution on the outskirts.
The fourth row demonstrates the radial distribution of velocity anisotropy.There is a weak dependence of the velocity anisotropy on the stellar mass in NH and HAGN.However, there is a distinct difference in the mass dependence between NH and HAGN.In NH, more massive galaxies have more negative values of β and thus tangential anisotropy, indicating the presence of rotational motion.On the contrary, the HAGN galaxies typically have positive values of β indicative of radial anisotropy.Overall, the change in the velocity dispersion profile with respect to the stellar mass in panels (a) and (e) appears to be primarily driven by the change in the density profile.
There appears to be some degree of consistency between NH and HAGN in the outer velocity dispersion and anisotropy profiles for common mass bins (orange lines) which correspond to the most massive galaxies of our NH sample and the least massive galaxies of our HAGN sample.However, it is worth noting that the HAGN and NH samples exhibit different trends in the central velocity dispersion profile and density slope, even for similar mass ranges.This can be partly attributed to the difference in the environments they represent.HAGN, being a large-volume simulation, covers a wide variety of environments and thus includes all types of galaxies.In contrast, NH represents a field region; thus, its massive galaxies are mostly late-type, exhibiting a dynamically hot bulge and a rotating cold star-forming disk, leading to a steeper profile and tangential anisotropy in the inner region compared to its counterpart in HAGN.Additionally, the differences in spatial and mass resolutions between the two simulations may have led to variations in the stellar dynamics of the inner regions, which suggests that the inner kinematic information of HAGN galaxies, particularly for lower mass galaxies, may not be reliable.
In Figure 10, the median profiles of the parameters are plotted in each row in the same manner as in Figure 9, except that the galaxies are now subsampled by ex-situ fraction.To address the bias introduced by the correlation between the stellar mass and ex-situ fraction (as shown in Figure 8), we first narrow the stellar mass ranges to log(M * /M ⊙ ) = 9.5-10.25 for NH and log(M * /M ⊙ ) = 10.75-11.25 for HAGN.These mass ranges correspond to the subsamples with the lowest masses, as shown in Figure 9.In the first row, the σ r profiles at large radii clearly exhibit dependence on exsitu fraction in both simulations.In the second-third rows, no strong dependence on ex-situ fraction is observable for the total density or stellar density profiles. .Median radial velocity dispersion profiles and three parameters for simulated galaxies (same as Figure 9) with different ex-situ fractions.The variation of anisotropy profiles primarily drives the changes in the outer slopes of velocity dispersion profiles.
In contrast, we found a strong correlation between the anisotropy profile β(r) and ex-situ fraction (panels (d) and (h)).In both NH and HAGN, the galaxies show a nearly isotropic velocity distribution at the inner radii, regardless of ex-situ fraction.At large radii, high ex-situ fraction galaxies remain isotropic (β(r) ∼ 0), whereas low ex-situ fraction galaxies exhibit a tangentially biased motion (β(r) < 0).This negative slope of the velocity anisotropy profile leads to a reduction in the radial velocity dispersion (σ r (r)) at large radii, as shown in panels (a) and (e).This reduction in σ r (r) results in a deviation between the galaxies with different ex-situ fractions.
Because these three parameters contribute to the σ r profile in combination, it is not trivial to determine the direct effect of each parameter.Therefore, to analytically remove the effect of the radial variation of the velocity anisotropy on the σ r profile, we assume a flat velocity anisotropy.More explicitly, we multiply σ r (r) by ((α * − 2β − 2γ r )/(α * − 2γ r )) 1/2 .This corresponds to how the σ r profile would look for the case of isotropic velocity distribution (β(r) = 0).The results are shown in Figure 11 binned by the stellar mass (top row) and ex-situ fraction (bottom row), respectively.In the upper panels, a systematic variation is visible between subsamples, indicating that the stellar mass dependence still exists.This means that the mass dependence of the σ r profile cannot be explained based on the correlation between mass and velocity anisotropy alone.In contrast, the bottom panels show that the dependence on the exsitu fraction no longer exists after the anisotropy effect is removed.In summary, when galaxies with different stellar masses are compared, the mass profiles (both total and stellar) are the key factors driving the change in the velocity dispersion profile.However, when the stellar mass is fixed, a radial variation in the velocity anisotropy, which is closely related to ex-situ fraction, is the key driver that causes a change in the outer σ r profile gradient.

Impact of ex-situ fraction and radial anisotropy
Our results indicate that there is a wide variation in the outer velocity dispersion slopes among both simulated and observed galaxies.This diversity is primarily caused by the different fractions of ex-situ stars present Figure 12.Velocity anisotropy (top panels) and total mass density slope (bottom panels) at 3 r50 as a function of exsitu fraction (left panels) and outer velocity dispersion slope (right panels), drawn for HAGN galaxies (green lines) and NH galaxies (black markers).Shaded region and error bars indicate 1σ scatter.Velocity anisotropy shows a good correlation with ex-situ fraction and outer σ slope.
in the galaxies (see Figure 8).Furthermore, galaxies with higher ex-situ fractions exhibit less tangentially biased velocity distribution (see Figure 10-d & h).In Figure 12 we further investigate the roles of α and β.We take the HAGN galaxies as our main sample here because they cover a wide range of galaxy types and environments.In order to minimize the mass effect while securing a sufficiently large sample for statistics, we limit the HAGN galaxy sample to log M * /M ⊙ = 10.75-11.25.
For reference, we show the NH sample, too.We divide the NH sample by mass: the open circles show the same mass range as that used for the HAGN galaxies in this diagram, and the closed circles show the lowermass sample which represents the majority of NH.The figure shows anisotropy and density slope measured at 3 r 50 of simulated galaxies with respect to ex-situ fraction and γ out .We find a strong correlation between the anisotropy and ex-situ fraction (panel (a)), while the correlation between the density slope and ex-situ fraction is weak at best (panel (c)).Similar trends are observed when compared with the outer slope of the velocity dispersion profile (panels (b) and (d)).Panel (c) shows that ex-situ fraction does not affect the outer density slope, which has already been demonstrated in Figure 10 (second-third rows).A change in γ out occurs pri-marily through ex-situ stars via a change in anisotropy rather than a change in density slope.
In the case of most massive galaxies, it becomes difficult to investigate the effects of accretion because they exhibit little variation in ex-situ fraction among themselves.
For example, the HAGN galaxies of log(M * /M ⊙ ) ≳ 11.6 have an ex-situ fraction of 80-100% (see Figure 8).Their outer velocity dispersion profiles are affected both by density slopes (i.e., mass) and anisotropy (see Figure 9).As a result, they exhibit 1) higher slopes in outer velocity dispersion profiles compared to lower mass counterparts, accompanied by 2) radial anisotropy.This finding is consistent with observations of massive early-type galaxies, which also display higher slopes of velocity dispersion and positive values of the h 4 parameter (Veale et al. 2018;Loubser et al. 2018).
However, it is widely accepted that the radial anisotropy causes the projected velocity dispersion to have a falling profile (Gerhard 1993), which appears to contradict our finding of higher slopes in outer velocity dispersion profiles mentioned above.For example, Loubser et al. ( 2020) demonstrated an expected anticorrelation between the anisotropy parameter β and the velocity dispersion slope.However, they assumed a constant anisotropy, which differs from the radially varying anisotropy observed in our simulated galaxies.
In fact, given that all our simulated galaxies have a more isotropic velocity distribution (β = 0) toward the center, the effect of radial anisotropy gradient is the opposite of that expected from the projection effect.According to Equation 7, with the total and stellar mass profile fixed, the velocity dispersion should move towards a rising profile if radial anisotropy increases with increasing radii (Figure 9).This can be interpreted in a way that galaxies with high ex-situ fractions show an intrinsically (i.e., in 3-D) rising velocity dispersion profile that is significant enough to overwhelm the projection effect.It is natural to expect that the accreted stars will maintain momentum along the infalling path, making them follow a radially biased orbit, and also mainly reside in the outskirt due to the high kinetic energy gained from the potential of the accreting galaxy.
Figure 13 displays the same plot as Figure 8 but with the background color replaced with the anisotropy parameter at 3 r 50 .The intrinsic slope generally shows a vertical gradient of increasing β with the outer slope.This can be understood as a competitive contribution between the in-situ stars and the ex-situ stars, with the former having tangential anisotropy due to disk-like kinematics and the latter having radial anisotropy due to the infalling velocity during the accretion or merger.Outer anisotropy, (3r 50 ) Figure 13.Same as Figure 8, except for the background color pixels indicate anisotropy parameter measured at 3 r50.More massive galaxies exhibit more radial velocity anisotropy and higher velocity dispersion slopes, which are consistent with observations.
When the projection effect is taken into account, the outer slope of the most massive galaxies shifts toward a falling profile.The upper part of the figure (rising and flat slope) exhibits a weak gradient with different directions to that of ex-situ fraction.This is primarily driven by a collective, but differential shift of galaxies toward lower outer slopes depending on their degree of radial anisotropy.

Comparison with other literature
Our result is qualitatively consistent with Wang et al. (2022), which employed a similar technique to ours in a cosmological simulation IllustrisTNG (Pillepich et al. 2018b) based on moving mesh hydrodynamics.They focused on ETGs for comparison with counterparts of the MASSIVE dataset (Ma et al. 2014).They found that the outer velocity dispersion slope, despite the difference in the definition of the radius where the slope is measured (1.5r eff ) compared to ours (3 r 50 ), is positively correlated with stellar mass, velocity anisotropy, and major merger fraction, which is in line with the ex-situ fraction.However, on the contradicting observation of positive σ gradient with positive h 4 parameter, they suggest h 4 is uncorrelated to the stellar orbital anisotropy.This differs from our explanations based on the variation of orbital anisotropy.This can be a subject of a future study.
The study of Cannarozzo et al. (2023) also utilizes the IllustrisTNG simulation and compares its velocity dispersion profile with the MaNGA IFU survey (Bundy et al. 2015;Yan et al. 2016), providing velocity dispersion profiles of in-situ and ex-situ stars in galaxies.Their result indicates that the overall velocity dispersion of exsitu stars shows a slightly higher value than in-situ stars (Figure 3, 4th row in their result), which is consistent with our result (Figure 6 and 7.However, for larger radii and more massive galaxies, the difference appears not to be as significant as ours.Further investigation is needed to explore this aspect in more detail.

Caveats
There are several caveats in this study.Beamsmearing is one of the major obstacles in measuring velocity dispersion profiles in observations.In our main results, we did not account for the beam-smearing effect on 2-D velocity dispersion profiles for simulated galaxies.This could potentially introduce bias to the measurement of slope in the velocity dispersion profile when compared with observations.Another related issue is the variation in beam size applied to different observational datasets, leading to potential inconsistencies in comparisons.We discuss in Appendix C the impact of the artificial beam-smearing effect applied to simulated galaxies.The result suggests that the impact of beamsmearing is not large enough to alter the main trends found in Figures 4 and 5. From the result, we also anticipate that the influence of differences in beam sizes between observations will be small.
Dust absorption is not considered in estimating the luminosity-weighted velocity dispersion profile from simulated galaxies.This could result in an unusually bright central region, reducing the effective radius of galaxies, as well as affecting the central velocity dispersion profile.This problem is particularly prominent for the highmass NH galaxies, which have both high resolution and bright central bulges.
The use of Voronoi cells may result in an overestimation of velocity dispersion in sparse regions, particularly where cell size increases and the spatial gradient of the stellar group velocity becomes significant.This issue is prevalent in the measurement of 2-D velocity dispersion profile in both observation and simulation.While it may be particularly influential for the outskirts of low-mass galaxies flattening the velocity dispersion profiles, the actual bias may not be very significant considering the large scatter in their outer σ slope, as shown in Figure 5.
Throughout the paper, we consistently employed a fixed shape and position angle for the ellipses of each galaxy in the 2-D velocity dispersion measurement, both for observed and simulated galaxies.The purpose of employing ellipses is to account for the effect of inclination on the binning of stellar information, ensuring that the binning optimally reconstructs the 3-D shells surrounding the galaxy.This is based on the assumption that visually elongated galaxies possess an oblate shape or are composed of ideal thin disks.It should be acknowledged that some galaxies have triaxial shapes (King 1978;Kimm & Yi 2007;Santucci et al. 2022) or warped disks (Sancisi 1976;Briggs 1990) that likely have a radial variation in shape and position angle when fitted by ellipses.We do not account for this effect, and as a result, our velocity dispersion measurements are subject to uncertainties related to it.
Throughout the kinematic analysis of the 3-D velocity dispersion profile, we assumed each simulated galaxy to have a spherically symmetrical mass distribution and to be in an equilibrium state where inflow and outflow are equal.This simplification may not fully represent the actual structure of galaxies.However, we observed that Jeans' equation reasonably well reproduces the expected velocity dispersion profile of simulated galaxies from a given density and velocity anisotropy profile, indicating that our assumptions are valid to some extent.

CONCLUSION
In this study, we have used two cosmological simulations, Horizon-AGN (for log(M * /M ⊙ ) > 10.75 at z = 0.017) and NewHorizon (for log(M * /M ⊙ ) > 9.5 at z ≲ 1), to investigate the velocity dispersion profiles of galaxies.We have measured the slopes of the velocity dispersion in the 2-D plane at the inner (0.5 r 50 ) and outer (3 r 50 ) radius and compared them with integral field spectroscopy data from two surveys, SAMI (for z < 0.115) and CALIFA (for z < 0.03) , for galaxies within the range of log(M * /M ⊙ ) > 9.5.Using the versatility provided by numerical simulations, we have investigated the variations in the velocity dispersion profiles in connection with the spatial age distribution, mass profile, and velocity anisotropy, and found their close connection to the stellar mass and ex-situ fraction of galaxies.Our main results can be summarized as follows.
1.The inner slope of the velocity dispersion profile shows a "V-shaped" trend with stellar mass in both simulated and observed samples.The slope becomes more negative (moving towards a falling profile) with increasing stellar mass until log(M * /M ⊙ ) ∼ 11, and appears to bounce back thereafter (Figure 4).The outer slopes of velocity dispersion profiles show a large variation.Lower-mass galaxies (log(M * /M ⊙ ) ≲ 11) exhibit no strong dependence on mass.There is a hint of an increasing slope with stellar mass in more massive galaxies (log(M * /M ⊙ ) ≳ 11), where the observational and simulation data are in agreement with each other (Figure 5).

Han et al.
2. The shape of the inner velocity dispersion profile of a galaxy is determined by the relative radial distribution of young and old stellar populations.The change in the inner slope of the velocity dispersion profile with increasing stellar mass is driven by a transition in the order of the star formation process, from outside-in to inside-out (Figure 6).
3. Ex-situ stars exhibit higher velocity dispersion than in-situ stars and are more predominant in the outer region of the galaxy, hence increasing the mean velocity dispersion at the galactic outskirts (Figure 6).As a result, galaxies with higher ex-situ fractions have higher outer velocity dispersion slopes (Figure 8).
4. The stellar mass dependence of the velocity dispersion profile is primarily driven by changes in the mass profile (Figure 9).However, the diverse outer velocity dispersion slope among galaxies with similar masses is primarily driven by the radial variation of anisotropy (Figure 10), which tightly correlates with ex-situ fraction (Figure 12).On the other hand, ex-situ fraction does not affect the outer density profile (Figure 12).
5. We found that most massive (log(M * /M ⊙ ) ≳ 11.6) galaxies exhibit a very high ex-situ fraction (80-100%, Figure 8), resulting in a higher outer velocity dispersion slope and radial anisotropy, which is consistent with observations.This investigation suggests that the recent star formation and past assembly history of a galaxy contribute to shaping its velocity dispersion profile.However, this perspective may be viewed as an oversimplification, given the various complications involved, such as environmental effects, the presence of supermassive black holes, and different modes of star formation.Pinning down the exact processes observationally proves challenging, but exploring them in numerical simulations could offer valuable insights.This topic holds promise for future studies.
The result presented in this paper is largely independent of redshift.This is indicated in Figures 4 and 5, by comparing the final snapshot (z = 0.17) and the full sample (z < 1) of NH galaxies.This is because, while redshift evolution may alter the absolute values of the velocity dispersion of galaxies (Cannarozzo et al. 2020), our analysis focuses on measuring the log-log slope of σ(r) at inner and outer radii.The trend we find is more dependent on the specific history of individual galaxies rather than on the global evolution of galaxies with cosmic time.In addition, introducing an artificial beamsmearing effect into the simulation by applying a Gaussian kernel-based measurement of the velocity dispersion profile does not significantly change our conclusion, as indicated in Appendix C.
We propose two parameters that can potentially serve as indicators for measuring the ex-situ fraction of galaxies: the outer slope of the velocity dispersion profile as a weak indicator, and the radial anisotropy of the velocity distribution as a strong one.The former can be measured relatively easily from stellar absorption line profiles, although it is still challenging to get it done on outer regions of galaxies.Measuring the latter presents an even greater challenge.Measuring the higher order Gauss-Hermite moment h 4 requires even more accurate spectroscopy.Future observations should aim to perform these measurements for a greater sample of galaxies, covering wider ranges of mass and environments, than what is available now.
Analyzing the kinematic properties of today's galaxies and inferring past events are pivotal to understanding the role of mergers and accretions in the evolution of galaxies.The measurement of ex-situ fraction for example links the kinematic properties of galaxies at the galactic scale and their hierarchical assembly histories in a large-scale environment.
One interesting feature of this work is the presence of an upturn in the trend of inner velocity dispersion slopes with stellar mass (Figure 4).Both observations and simulations consistently suggest negative correlations for lower masses and positive correlations for higher masses.This strongly indicates the presence of a turning point around log(M * /M ⊙ ) ∼ 11.The "V-shaped" trend is likely the result of two effects impacting the inner slope in different directions.Lower-mass galaxies, with low ex-situ fractions, are more influenced by the properties of in-situ stars.The differential distribution of young and old in-situ stars results in negative correlations of the inner slope with mass.However, above a certain stellar mass, ex-situ stars that are more abundant in more massive galaxies may infiltrate the inner regions of galaxies, influencing the inner velocity dispersion profiles, producing the upturn in the trend.
It is also intriguing that the stellar mass of the upturn in the trend corresponds to the point where the galaxy transitions from a star-formation-dominated phase to a merger-dominated phase in the size evolution (van Dokkum et al. 2015) This is consistent with observations showing that the upturn is not present for late-type galaxies, which are not star-formation quenched yet and have not undergone dry mergers (Figure B1).For earlytype galaxies, the turning point roughly corresponds to the stellar mass where the distinction between fast and slow rotators occurs (Emsellem et al. 2007(Emsellem et al. , 2011)).This is consistent with our scenario of in-situ to ex-situ transition, as slow rotators are expected to be a result of dry mergers (Choi & Yi 2017).
The use of two simulations (NH and HAGN) in combination in this study reproduced the observed upturn in the inner slope.However, the significant difference in astrophysical prescriptions and (spatial and mass) resolutions between the two simulations poses a major caveat, making it challenging to interpret the results and determine.This emphasizes the need for high-resolution simulations that are capable of reproducing a diverse galaxy population covering a wide range of stellar masses and environments under a single astrophysical prescription.If the current computing resources are limited to perform such a simulation in one piece, it would be a good strategy to run separate zoom-in simulations using identical physical ingredients and resolutions.
We thank the anonymous referee for the comments that clarified the manuscript greatly.This work was granted access to the HPC resources of CINES under the allocations c2016047637, A0020407637 and A0070402192 by Genci, KSC-2017-G2-0003 by KISTI, and as a "Grand Challenge" project granted by GENCI on the AMD Rome extension of the Joliot Curie supercomputer at TGCC.The large data transfer was supported by KRE-ONET which is managed and operated by KISTI. .The inner velocity dispersion slope of ETGs and LTGs measured at 0.5 r50 as a function of galaxies' stellar mass in the same format as Figure 4 in the main text.Lines in each panel show simulated galaxies in HAGN and NH classified as ETG and LTGs based on the criterion (v/σ)3D = 1 (see the text for more details), with shades indicating the 1σ scatter.The NH galaxies are in good agreement with the trend of late-type galaxies in observation.To ensure the consistency of kinematical classification, we excluded two galaxies with counter-rotating features from the NH sample.

B. INNER VELOCITY DISPERSION SLOPE OF GALAXIES WITH DIFFERENT MORPHOLOGY
Figure B1 shows for reference the inner profiles of velocity dispersion of the simulated galaxies in the same format as Figure 4 in the main text but for observed early-and late-type galaxies (ETGs and LTGs) separated.For SAMI and CALIFA, galaxy morphologies are determined through visual classification.We employed the morphology information from Cortese et al. (2016); Croom et al. (2021) for SAMI and Pak et al. (2019) for CALIFA.We considered ellipticals (E) and leticulars (S0) as ETGs, and the rest as LTGs.We use the criterion (v/σ) 3D = 1 to distinguish ETGs and LTGs (or dispersion-and rotation-dominated galaxies) following Dubois et al. (2016, MNRAS, 463, 3948), where (v/σ) 3D indicates the ratio between maximum stellar rotation speed and 3-dimensional stellar velocity dispersion.We note that the criterion is for the rough classification of simulated galaxies and may not be indicative of visual morphology in some cases.For example, a large fraction of HAGN late-type galaxies in Panel (b) may actually be fast-rotating elliptical or lenticular galaxies.The massive NH galaxies are mostly late-type in many aspects because of being in the field environment.The most massive galaxies in the HAGN simulations (green lines and shades) show similar trends regardless of the morphological classification used.

C. SLOPE OF VELOCITY DISPERSION WITH DIFFERENT PSF SIZE
To investigate the effect of beam-smearing, we present Figure C1, which shows the inner and outer slope of the velocity dispersion (same as Figure 4 and 5).The solid lines represent the original data from previous figures.The dashed and dotted lines represent simulated data with Gaussian smoothing applied, to consider the effect of two sizes of the point spread function (PSF) which correspond to the half-width at half-maximum (HWHM) of CALIFA (0.38 kpc) and SAMI (1 kpc), respectively.The PSF-applied data shows differences only at the high-mass end of NH galaxies.The inner and outer slopes are affected by ∼ 25% and ∼ 50%, respectively, at log(M * /M ⊙ ) ∼ 11 of NH, where the effect is most pronounced.This is likely to be a consequence of the changes to the velocity dispersion profile being particularly influential for the galaxies with a strong central peak in the velocity dispersion profile.Because the luminosity of the bright central parts of massive NH galaxies is dispersed to surrounding regions, the effect may also have been amplified, thereby influencing the luminosity-weighted chi-square fitting of the broken power-law profile.

FluxFigure 1 .
Figure 1.Mock IFU image of a galaxy in NewHorizon simulation.The face-on (upper row) and edge-on (lower row) projections are shown.Each column represents the projected distribution of r-band flux (a), line-of-sight mean velocity (b) and velocity dispersion (c), ex-situ fraction (d), age of in-situ stars (e).The side length of one panel is 40 kpc.The white and black lines indicate the effective ellipse measured in isophote fitting.The difference in velocity dispersion between the old bulge and the young disk leads to a steep increase in velocity dispersion in the inner region.

Figure 2 .
Figure 2. Mock IFU images simiar to Figure 1, with a less massive galaxy from NewHorizon.The side length of one panel is 40 kpc.Compared to its more massive counterpart in Figure 1, the galaxy contains younger stars in the inner region, resulting in a relatively flat velocity dispersion profile at the center.

Figure 3 .
Figure3.Examples of velocity dispersion profiles of CAL-IFA galaxies, with their photometric image.The left column shows the radial trend of velocity dispersion (black), with a fitted broken power-law curve (red).The right column shows the SDSS image of target galaxies with stellar mass, name, and morphology.Galaxies with similar visual morphology exhibit different trends in their velocity dispersion profiles.Galaxy morphologies are derived from the same dataset used inPak et al. (2019), where visual classification was performed.

Figure 4 .
Figure 4.The inner velocity dispersion slope measured at 0.5 r50 as a function of galaxies' stellar mass.Solid lines show full sample of our simulated galaxies, with shades indicating the 1σ scatter.The data for NH galaxies from the final snapshot (z = 0.17) is shown as a dotted line.Observed data are also presented: 1σ scatters are shown by arrows and the standard errors of the median are shown by bars.A fall-and-rise trend can be seen with a turning point at log(M * /M⊙) ∼ 11.

Figure 5 .
Figure5.The outer slope measured at 3 r50 as a function of galaxy's stellar mass.Similar schemes with Figure4are employed, except that the SAMI and CALIFA samples are combined (blue symbols).Overall, the outer slope exhibits a large variation.Massive galaxies (log(M * /M⊙) ≳ 11) exhibit a hint of an increasing slope with stellar mass.

Figure 6 .
Figure 6.Radial velocity dispersion (σ3D) profiles (panels (a) and (b)) and mass fraction profiles (panels (c) and (d)) of the NH galaxies.The left and right panels are for two different mass ranges.The whole stellar data (black dashed line) is divided into young in-situ (blue solid line), old in-situ (orange solid line), and ex-situ stars (red dash-dotted line).An age cut of 3 Gyr separates the young and old stars.The shades indicate the 1σ scatters.Low-mass and high-mass galaxies exhibit an inverted radial composition of young and old stars.

Figure 7 .
Figure 7. Similar to Figure 6, but based on the HAGN sample with different stellar mass ranges, log(M * /M⊙) = 10.75-11.5 and log(M * /M⊙) > 11.5.Most massive galaxies in the simulation are dominated by ex-situ stars, and the galaxy's overall velocity dispersion profile is constant.

Figure 8 .
Figure 8. Outer slopes of velocity dispersion profile as a function of stellar mass in the range of log(M * /M⊙) = 10.75-12.25,which corresponds to the massive end of Figure 5, for the HAGN galaxies.Outer slopes before (a) and after (b) applying projection effects are shown.Color pixels in the background represent the median ex-situ fraction.More massive galaxies exhibit higher ex-situ fractions.Galaxies with higher ex-situ fractions have increased outer velocity dispersion slopes.

Figure 9 .
Figure 9. Median radial velocity dispersion profile, normalized by its value at r50 (1st row) and three parameters (2nd-4th row) of the Jeans' equation, namely α(r) (total density slope), α * (r) (stellar density slope) and β(r) (velocity anisotropy) are shown as functions of the radial distance, normalized by r50.Each column represents NH (left) and HAGN (right) galaxies.Each line represents a distinct range of stellar mass, indicated by different colors.Shades indicate 1σ scatter.The total and stellar density profiles are dependent on the stellar mass, which primarily drives the changes in the velocity dispersion profile.
Figure10.Median radial velocity dispersion profiles and three parameters for simulated galaxies (same as Figure9) with different ex-situ fractions.The variation of anisotropy profiles primarily drives the changes in the outer slopes of velocity dispersion profiles.

Figure 11 .
Figure11.Same σ3D profiles as Figure9and Figure10, except for assuming β = 0 for all radii to eliminate effects from velocity anisotropy.The variation of outer velocity dispersion profiles in the function of ex-situ fraction is primarily caused by velocity anisotropy.
Figure B1.The inner velocity dispersion slope of ETGs and LTGs measured at 0.5 r50 as a function of galaxies' stellar mass in the same format as Figure4in the main text.Lines in each panel show simulated galaxies in HAGN and NH classified as ETG and LTGs based on the criterion (v/σ)3D = 1 (see the text for more details), with shades indicating the 1σ scatter.The NH galaxies are in good agreement with the trend of late-type galaxies in observation.To ensure the consistency of kinematical classification, we excluded two galaxies with counter-rotating features from the NH sample.

Figure C1 .
Figure C1.Similar to Figures4 and 5, the simulated data includes the results of Gaussian convolution with two different kernel sizes on the 2-D plane, applied to the line-of-sight velocity distribution.HWHM of 0.38 (dashed lines) and 1 kpc (dotted lines) correspond to the physical PSF sizes for CALIFA and SAMI galaxies, respectively, measured based on the mean redshift and angular resolution of the data.The solid lines represent the original slopes without the beam-smearing effect.

Table 1 .
Summary of sample data Observation data † Corresponds to the best spatial resolution for the case of simulations and rough estimates of HWHM at the mean distance for the case of observations.Note-For simulated 2-D datasets, sample sizes are multiplied by 24 after employing different line-ofsight projections.