© The Authors 2026

1 Introduction

Star clusters are typically groups of stars bound together by gravitational interactions (Bergond et al. 2001; Lada & Lada 2003; Piskunov et al. 2006), simply classified into two main types based on their morphology: globular clusters (GCs; Fujii et al. 2024) and open clusters (OCs; Perryman 2026). Open clusters that form from the same giant molecular cloud (GMC; e.g., Meynet et al. 1993; Alves & Bouy 2012; Megeath et al. 2016) are initially in a nesting phase, characterized by high-density gas accompanied by strong gas loss. After that, the number of members of OCs gradually changes over time, usually decreasing. For instance, two-body relaxation (e.g., Bergond et al. 2001; Mamon et al. 2019; Broggi et al. 2024) will lead to energy redistribution, causing some low-mass member stars to migrate to the tidal tail or even escape from the OCs (Portegies Zwart et al. 2010; Tang et al. 2019; Akhmetali 2026). Since the morphology of OCs is concerned with the number of their members, changes in the number of member stars will inevitably affect the OC morphology. This demonstrates that the morphology of OCs is closely related to their evolution. Therefore, studying the morphology of OCs can enhance intuitive understanding of their evolution (e.g., Meingast & Alves 2019; Dalessandro et al. 2021; Pang et al. 2021; Hu et al. 2021a,b; Tarricq et al. 2021, 2022; Pang et al. 2022; Della Croce et al. 2023).

The release of high-precision astrometry data from the Gaia satellite¹ (Gaia Collaboration 2016, 2023) has enabled the precise reconstruction of the 2D or 3D structure of OCs. Studies have revealed that some OCs exhibit morphological asymmetries, core–halo structures (Kadam et al. 2024), extended tidal tails (Sharma et al. 2025), and even signs of imminent disintegration (Almeida et al. 2025). Moreover, there is also an extensive body of literature on the statistical morphology of OCs. For example, Zhong et al. (2022) conducted a structural analysis of 256 OCs, revealing the radial density profile (Camargo et al. 2012) in their outer halo regions deviating significantly from the King model profile (King 1962). Subsequently, Hu et al. (2023) employed the rose diagram overlay method to study the morphological structures of 88 OCs, primarily analyzing their three 2D projections in the heliocentric Cartesian coordinate system to investigate layered structures and found that 74 of the 88 sample OCs present clear 3D layered structures. In 2024, Hu et al. (2024) defined the “morphological coherence” for the first time to quantify the 3D morphological difference between the inner and outer regions of 132 OCs located within 1 kpc of the Sun, and investigated, in a Cartesian reference frame, the full 3D manifestation of this morphological coherence. Hu et al. (2025b) examined the structural parameters of 105 OCs, suggesting that OCs with more member stars within their tidal radii possess stronger gravitational binding forces. In the same year, Xu et al. (2025) performed a 3D principal component analysis on eight OCs, revealing that the tidal structure orientation of certain OCs is perpendicular to their direction of motion. In addition, Lang et al. (2025) performed N-body simulations of 279 nearby OCs, correlating member star spatial distributions with layered structures, and suggested that OCs with fewer members typically lack such structures.

Although significant progress has been made in studying the morphological structure of OCs in recent years, most studies have been concerned with 3D spatial distributions, limited by sample size, especially the research on the morphological stability of OCs. As there are symmetrical errors in the parallax measurement of Gaia (Gaia Collaboration 2016, 2023), a direct propagation to distant distances through an inverting parallax (1/ω) could lead to a significant artificial elongation of the 3D spatial distribution of OCs along the line of sight (see, e.g., Smith & Eichhorn 1996; Bailer-Jones 2015; Luri et al. 2018; Carrera et al. 2019; Hu et al. 2023). Generally, to recover the true 3D spatial distribution of OCs, one restricted the OCs to d ≤ 500 pc precisely to keep their distance uncertainties below ~6 pc after distance correction (e.g., Hu et al. 2023). However, OC samples within 500 pc are limited. Therefore, to extend the OC sample and avoid the fake elongation, we adopt, in this work, a 2D projection in the Galactic spherical coordinate system to explore the morphological stability of a large number of OCs that are almost beyond 500 pc.

While this study did not directly investigate the morphological stability of OCs in 3D space, we can still obtain and study their morphological stability in specific projection planes in 3D space through the projection effect. In this paper, we systematically analyzed the morphological stability of 1490 OCs by the rose diagram construction, which was more than 10 times the sample number in previous studies (Hu et al. 2023, 2024, 2025b).

The paper is organized as follows. In Sect. 2, we introduce the data used in our work. We then describe the methodology adopted in this study in Sect. 3. Section 4 presents all results of our entire work. We discuss the sensitivities and potential systematic biases of our approaches in Sect. 5. In the end, Sect. 6 provides a summary and conclusion.

2 Data

Based on Gaia Data Release 3 (DR3) data, van Groeningen et al. (2023) provides a catalog of OC member stars complete to G > 18 mag. Hunt & Reffert (2024) recently determined 5647 reliably bound OCs, marked by an “o” symbol, from a sample of 6956 star cluster candidates. These two catalogs serve as the primary data sources for this work. We selected all OCs flagged with o from Hunt & Reffert (2024) and cross-matched them with the corresponding member star catalog from van Groeningen et al. (2023) by cluster names in topcat². In order to obtain reliable member stars, we retained only those members with membership Pmemb = 1 from van Groeningen et al. (2023). Previous studies and catalogs (e.g., Cantat-Gaudin et al. 2018; Cantat-Gaudin 2022; Hu et al. 2023) have typically considered OCs with at least 50 members to be robust for analysis. Thus, we limited the sample OCs to those with a number of member stars greater than or equal to 50 members, and a total of 1490 OCs remained. Finally, the member stars of 1490 OCs were cross-matched with the Gaia DR3 data to obtain their fundamental parameters.

3 Methods

To investigate the morphological stability of OCs, we first calculated their morphological stability parameters by the rose diagram construction. The specific approaches are given below.

3.1 Construction of the rose diagram

We adopted the rose diagram method developed by Hu et al. (2025b) to quantify the morphologies of all OCs in our sample. We first determined a plane distributed by its members for each OC. Due to the inclination of the equatorial plane of the celestial coordinate system relative to the Galactic midplane, taking an OC as an example, we used the position coordinates and parallax of its members and their uncertainties (RA_ICRS, DE_ICRS, Plx, e_RAdeg, e_DEdeg, e_Plx) to transform them into the Galactic coordinates. We employed a Monte Carlo sampling to properly propagate the observational uncertainties in this coordinate transformation. For every member of the OC, we performed 100 samples by drawing RA_ICRS, DE_ICRS, and Plx from their respective error distributions to obtain their Galactic spherical coordinates and uncertainties (l, b, e_l, e_b) via the astropy³ package (Astropy Collaboration 2013, 2018). For every sampling, we took the median of the Galactic coordinates (l, b) as the centric coordinate (l₀, b₀) of the OC. At the same time, the median of sampled parallax (Plx) was converted into distance (d).

Next, the converted coordinates and distances were substituted into Eq. (1) to project the positions onto the Galactic projection that is perpendicular to the line of sight: $$\{\begin{array}{l}\mathrm{\Delta }\mathrm{x}=\mathrm{d}\cdot \sin (1-{\mathrm{l}}_0)\cdot \cos \mathrm{b}\\ \mathrm{\Delta }\mathrm{y}=\mathrm{d}\cdot (\cos {\mathrm{b}}_0\cdot \sin \mathrm{b}-\sin {\mathrm{b}}_0\cdot \cos \mathrm{b}\cdot \cos (1-{\mathrm{l}}_0)),\end{array}$$(1)

where l₀ and b₀ represent the coordinates of the center of the OC, with d being the distance from the OC. The mean and standard deviation of the 100 samples of (Δx, Δy) are taken as the final coordinates and their uncertainties, respectively. Figure 1 shows the Galactic distribution (left) and the flattened projection distribution of the members of NGC_6791 in physical space (right).

It is well known that only member stars within the tidal radius are bound by gravitation from their host OC. Therefore, we determined the tidal radius of each OC in the flattened spatial projection by the King model density profile fitting (King 1962): $$\rho (r)=k{[\frac{1}{\sqrt{1+{\left(r/r_{\mathrm{c}}\right)}^2}}-\frac{1}{\sqrt{1+{\left(r_{\mathrm{t}}/r_{\mathrm{c}}\right)}^2}}]}^2+c,$$(2)

where r_c is the core radius, r_t is the tidal radius, and k and c are model fitting parameters. The fitting was performed by the king_model function from the diagnostics subpackage of gaia_oc_amd⁴. This process was repeated for each set of 100 samples, gaining a distribution of tidal radii for each OC. As is shown in Fig. 1, the mean value (marked by a dashed black circle in the right panel) and standard deviation of tidal radii from the 100 sampling for NGC_6791 give a tidal radius of R_t = 26.95 ± 2.85 pc.

The rose diagram construction is a quantitative approach to analyze the morphological structure of OCs (Hu et al. 2023, 2025b). By partitioning and statistically analyzing the distribution of OCs’ members on a 2D projection plane according to azimuth angles within the spatial extent defined by the tidal radius (R_t) and combining it with polar coordinate visualization, one can effectively identify the morphological characteristics and asymmetries of OCs. The rose diagram method was first proposed by Hu et al. (2023). For the specific operation, readers can refer to this literature.

Prior to applying the rose diagram construction, we decomposed the flattened spatial projection planes of all sample OCs onto the planes perpendicular to the Galactic disk plane by multiplying the vertical coordinates of members on the flattened projections by the cosine of the median Galactic latitude of these members. In this study we obtained the rose diagram of sample OCs through the following steps. We first centered the coordinates in the plane corrected above of each OC with the kernel density estimation (Seleznev 2016) density peak, then divided the projected plane into 12 sectors of 30° each. For each sector, we calculated two normalized parameters: the number of member stars, n, and the median radial distance, d. The radius of each sector was then calculated as $\mathrm{R}=\sqrt{n^2+d^2}$, denoted by R_i, which is a dimensionless normalized quantity in the ith sector. The minimum value among these 12 radii was defined as the core radius (R_core), and the corresponding area $\pi R_{\text{core}}^2$ represents the core area (S_core). The total area of the irregular outer region (S_outer) was calculated as the sum of all areas in the sectors minus the core area. Their equations are as follows: $${\mathrm{R}}_{\text{core}}=min({\mathrm{R}}_{\mathrm{i}}),\mathrm{i}=[1,12],$$(3) $${\mathrm{S}}_{\text{core\hspace{0.17em}}}=\pi \cdot {\mathrm{R}}_{\text{core}}^2,$$(4)

and $${\mathrm{S}}_{\text{outer\hspace{0.17em}}}=\sum _{\mathrm{i}=1}^{n=12}\frac{\pi \cdot {\mathrm{R}}_{\mathrm{i}}^2}{12}-{\mathrm{S}}_{\text{core}}.$$(5)

Specifically, based on the mean tidal radius, R_t, of each sample OC and its error, σ_Rt, we constructed a normal distribution, $\mathcal{N}({\mathrm{R}}_{\mathrm{t}},{\sigma }_{{\mathrm{R}}_{\mathrm{t}}}^2)$. From this distribution, we independently drew ten tidal radius values. For each sampling iteration, j, we used ${\mathrm{R}}_{\mathrm{t}}^{(\mathrm{j})}$ as the boundary to reselect OC member stars, calculated their projected coordinates relative to the cluster center, and quantified their 2D spatial distribution through the rose diagram method. Repeating this process ten times, we obtained ten independent estimates of the core radius, R_core, core area, S_core, and outer area, S_outer. Finally, we reported the sample mean and sample standard deviation of these parameters as the best estimates of the sample OCs’ morphological parameters and their uncertainties, as compiled in Table 1.

By the processes mentioned above, we can obtain the rose diagrams of all sample OCs. Figure 2 shows the rose diagram of NGC_6791, with S_core marked in dark gray.

Fig. 1 Distribution of NGC_6791 in the Galactic coordinate system (left panel) and its flattened projection distribution (right panel). The scattered blue points in each panel represent the member stars of the OC, with the gray bars being the errors of the position coordinates of its members. The orange circles correspond to the tidal radii obtained from 100 Monte Carlo sampling runs. The dashed black circle represents the mean value of the 100 results.

Fig. 2 Rose diagram of the NGC_6791 on the 2D corrected projection, derived from one of ten samplings. The rose-colored sectors represent the normalized directional radial radius values obtained by the Rose diagram method. Each sector spans 30°, with the first sector starting at 0° corresponding to the horizontal line on the right of this diagram. The second sector begins at a position 30° counterclockwise from the first sector, and subsequent sectors follow this pattern. The semitransparent gray sectors are the average extension levels in each direction, with the scattered green points indicating the member stars.

3.2 The measurement of morphological stabilities

Hu et al. (2025b) first proposed the ratio of the core to the outer areas (S_core/S_outer, hereafter referred to as the morphological stability parameter I) to describe the morphological stability of the OCs. The essence of this parameter is consistent with the instability discussed by Hu et al. (2023), as both characterize the stability of an OC based on the structural ratio between the outer and the core region.

For an evolving OC, two-body relaxation (Portegies Zwart et al. 2001; Angelo et al. 2025), combined with the tidal stripping of weakly bound stars from the periphery, drives both the development of mass segregation and the global mass loss of the OC (Vázquez-Semadeni et al. 2017; Grudić et al. 2018). The macroscopic outcome is an overall contraction of the cluster scale, with the central region becoming increasingly dense and gravitationally dominant relative to the outer envelope. A sub-virial core (η < 1) in a contracting state deepens the central gravitational potential (Foster et al. 2015), accelerating the inward migration of massive stars and further enhancing the central density (Della Croce et al. 2025). Simultaneously, the outer envelope, heated by two-body relaxation and truncated by tidal forces, is progressively depleted as member stars approaching the escape velocity are preferentially stripped. Based on this theoretical expectation, we propose a direct observational diagnostic: the ratio of the number of member stars in the core to that in the outer region, N_core/N_outer (hereafter referred to as the morphological stability parameter II). This parameter is designed to quantify the degree of structural tightening a cluster has undergone due to internal relaxation and external stripping. A high N_core/N_outer value signifies a core that is gravitationally bound strongly enough to retain its stars, and an outer region that the combined effects of relaxation and tides have significantly depleted.

Morphological stability parameter I, S_core/S_outer, relies on the areas of the core and outer regions being determined. Morphological stability parameter II, N_core/N_outer, is calculated by directly counting the number of member stars located within these same respective regions.

4 Results

4.1 Distribution of sample clusters

Figure 3 shows the distribution of our sample OCs in Galactic coordinates. The majority of the sample OCs are concentrated at |b| ≤ 20°. Moreover, owing to the complex and dynamically active environment of the thin disk in which these OCs reside, most among them have only a few members and weak morphological stability. On the contrary, the sample OCs that are far from the disk usually have a larger number of member stars and exhibit stronger morphological stability. This means that in addition to the number of members, the morphological stability of OCs seems to be related to the spatial positions.

Figure 4 shows the histograms of the core and outer areas of the sample OCs. Here, we have only considered the histogram of areas because the sample OCs differ greatly in terms of the number of members in the core regions. We found from this figure that the core areas (S_core) and outer areas (S_outer) of our samples span a range from 0 to 4.5. The outer areas of most OCs lie between 1.0 and 2.5, while the majority of core areas are below 1.5. This suggests that for most clusters in our sample, the ratio S_core/S_outer is less than 1, which is consistent with Hu et al. (2025b), who report that more than 70 % of their sample OCs have a 2D morphological stability parameter below 1. This also indicates that most of the OCs have lower morphological stability. It can be supported by the sample distribution in Fig. 3, especially those with Log(S_core/S_outer) < 0. Moreover, it is apparent that the mean area of OC cores is smaller than that of OC outer regions, which suggests that a typical OC has a compact core and an extensive external structure.

We note that in this study the sample OCs with S_core = 0 or N_outer = 0 should be excluded from our samples, because the morphological stability parameters I and II derived from them are devoid of physical meaning. (Hu et al. 2025b) established that when the core area of an OC is zero, it signifies that the OC exhibits no discernible core structure in the projected plane. We found that our sample does not contain cases in which S_core = 0 but does contain cases in which N_outer = 0. Finally, a total of 16 OCs with N_outer = 0 were excluded. Their parameters can be found in Table 1.

Fig. 3 Map of the sample OCs in the Galactic coordinate system. The size of each color-filled circle is proportional to the number of members of each OC, with its color being coded by morphological stability parameter I.

Fig. 4 Histogram of the areas of the core (Score, blue) and the outer (Souter, red) of our sample OCs. Each set of data was fit with a normal distribution (solid colored lines), with the mean positions being marked with vertical dashed lines.

4.2 Relationships between the morphological stabilities and the number of members

Hu et al. (2025b) have established a linear relationship between morphological stability parameter I (S_core/S_outer) and the number of member stars (N) of OCs. In this section we verify this relationship and extend it further by exploring the correlation between the morphological stability parameter II (N_core/N_outer) and N.

To this end, we performed the linear fitting for them with the Bayesian linear fitting model (e.g., Anders et al. 2017; Hu & Soubiran 2025) and the Markov chain Monte Carlo (MCMC) method (e.g., Foreman-Mackey et al. 2013; Alfonso et al. 2024; Lang et al. 2025; Qiu et al. 2026). The Bayesian linear model is detailed in Appendix A. The fitting results are displayed in Fig. 5. We can see from this figure that there is indeed a strong linear correlation between S_core/S_outer and N for our sample OCs. Moreover, the correlation of N_core/N_outer with N suggests the number of core member stars (N_core) is a new metric for the gravitational binding in the morphological stability of OCs. As the number of member stars of OCs increases, particularly the number of member stars in the core, the gravitational binding capacity of the OCs strengthens, which conforms with the previous finding by Hu et al. (2025b). Furthermore, the slope derived from morphological stability parameter II is steeper than that from stability parameter I, with nearly twice the slope within the margin of error of their slopes. And the Pearson correlation coefficients are 0.5721 ± 0.0020 and 0.6052 ± 0.0006; these were derived from the median and standard deviation of Pearson regression of the fit line 1000 times for parameters I and II, respectively. Meanwhile, their P_values are separately 2.84 × 10⁻¹³⁰ ± 1.72 × 10⁻¹²⁶ and 5.35 × 10⁻¹⁴⁸ ± 2.61 × 10⁻¹⁴⁷, indicating that both correlations are strong. Thus, the ratio of member stars in the core to those in the outer seems to be a somewhat more sensitive probe of the morphological stability of OCs than the ratio of core area to outer area, although the difference between them within the uncertainties is modest.

Next, we further explored the relationship between the morphological stability parameters I and II. We fit linear parameter II against parameter I for the entire sample, as shown in Fig. 6. The relationship, Log(N_core/N_outer) = (1.244 ± 0.022) × Log(S_core/S_outer) + (1.275 ± 0.011), displays how the stellar density in the core increases more rapidly than the core area itself. This can support the hypothesis of the clusters with a more concentrated structure possessing stronger gravitational binding in their central regions.

Interestingly, from the top panel of Fig. 6, we can observe that the samples with the logarithm of the parameter I greater than 0 exhibit a noticeable deviation from the linear fitting line. This tells us that the correlation in this region may actually be stronger than that in other regions where the logarithm of the parameter I is less than 0. In addition, most of our samples have the logarithm of the parameter II greater than 0. Based on this, we divided the sample into four quadrants defined by the signs of Log(S_core/S_outer) and Log(N_core/N_outer), marked by two dashed lines in the top panel of Fig. 6. Due to the lack of samples in the quadrants with the logarithm of the parameter II less than 0, this work mainly focused on two quadrants with the logarithm of the parameter II greater than 0: Quadrants I (Q1, right upper) and II (Q2, left upper). It can be seen that the OC samples in Q1 have a large core area and numerous core members, while those in Q2 have a small core area and numerous core members. Compared to the samples in Q1, the OCs from Q2 probably possess a dense core surrounded by a loose envelope and may be in an advanced stage of disruption.

Furthermore, in order to figure out how different they are, we performed a quantitative analysis for the OCs from Q1 and Q2. It comes out that the slope of Q1 is approximately 2.8 times steeper than that of Q2, as shown in Fig. 7. This indicates that the structural concentration presented in Q1 samples is much more responsive to the density change in core members because of some processes such as dynamical heating (Camacho et al. 2025) or tidal stripping (Escala et al. 2025) that are altering their structure.

Moreover, the turnoff position along the morphological stability parameter I in Fig. 7 that we took to divide samples, Log(S_core/S_outer) = 0, illustrates that when the ratio of the core area to the outer area is greater than 1, the trend of the members’ number in the core greater than that in the outer would be obvious. Meanwhile, when morphological stability parameter I is greater than 1, the OC samples are typically with strong morphological stability (Hu et al. 2025b); the opposite holds true otherwise. However, regardless of whether parameter I is greater than 1 or less than 1, the morphological stability parameter II almost remains greater than 1. This indicates that the morphological stability parameter I possibly relies on the premise of the number of core members exceeding that of external members, since a larger number of core members within an OC enhances its potential capacity for self-constraint.

Theoretically, it is almost impossible for a typical OC with fewer members in its core than in its outer region to possess sufficient binding potential to maintain its structure, leading to eventual dissolution. This is why there is a lack of such sample OCs, as is illustrated by the field below the horizontal dashed line in Fig. 6. Even if such OCs exist, the estimated errors for their morphological stability parameters are large, most likely arising from the uncertainties of observational data and the extremely sparse distributions of their members.

After understanding the difference between Q1 and Q2, we further performed a linear fitting on the parameters I and II for the samples from these two regions against the number of members, as is shown in Figs. 8 and 9. All fitting parameters are listed in Table 2. Figure 8 displays the linear fitting of parameters I and II with respect to the number of members for the Q1 sample. It can be found that the slope of Log(N_core/N_outer) versus Log(N) is 1.106 ± 0.117, with the one of Log(S_core/S_outer) versus Log(N) being 0.135 ± 0.030. This illustrates that as the number of member stars increases, the proportion of core member stars increases sharply, supporting strong gravitational binding. However, unlike it, the core region may have approached its maximum physical density, limiting further contraction and resulting in slow area expansion.

The same method was applied to the Q2 sample, as is shown in Fig. 9. Due to the Q2 sample having Log(S_core/S_outer) less than 0, the slope of Log(N_core/N_outer) versus Log(N) is smaller than that of Q1 sample, with 0.627 ± 0.065, indicating that the proportion of member stars of the core continues to increase, but at a slower rate than in Group Q1. This can be explained by the fact that although the number of core members is still increasing, the spatial extent of the core grows faster, resulting in a decreasing stellar density. This behavior is commonly associated with internal dynamical interactions or external tidal perturbations, as Della Croce et al. (2024) demonstrated that a large fraction of young OCs experience pronounced global expansion during the first 30 Myr of their evolution, driven by the violent relaxation and expulsion of residual gas, and accompanied by a substantial loss of mass from their outer regions. In addition, the slope of Log(S_core/S_outer) versus Log(N) is 0.409 ± 0.045. The combination of a moderate increase in core stars coupled with a faster relative expansion of the core area is consistent with the scenario potentially driven by dynamical interactions or tidal forces (Krumholz et al. 2019).

In conclusion, our analysis confirms that morphological stability parameter II (N_core/N_outer) exhibits a stronger dependence on the number of members (N) than parameter I. The disparity between the Q1 and Q2 samples validates these parameters to diagnose the stability and evolutionary state of the OCs. The clear linear relationship between Log(N_core/N_outer) and Log(N) across the sample further underscores the critical role of the core stellar population in maintaining the structural integrity of OCs.

Fig. 5 Relationships between the morphological stability parameters (Score/Souter: the first two panels from left to right; Ncore/Nouter: the last two panels) and the number of members (N). For ease of visualization, a logarithmic scale was adopted. Blue scatters in linear fitting panels represent our sample OCs, with the gray bars denoting their parameter errors. The faint gray lines show 1000 MCMC samples, with the thick red lines (the median values) and the shaded yellow bands (1σ) corresponding to our best parameter estimates. The two corner plots display the marginal posterior distributions (diagonal) and the 2D joint posterior contours (below diagonal) of the fit parameters (m, b, σ). Contours enclose the 1σ and 2σ confidence regions, and the median together with the 16th and 84th percentiles are quoted in each panel.

Fig. 6 Relationship between Score/Souter and Ncore/Nouter for the whole sample through a logarithmic scale. The dotted black lines in the top panel correspond, respectively, to the zero positions on the horizontal and vertical axes, representing the scenario where both morphological stability parameters I and II are equal to 1. The color bar is coded by the number of member stars of each sample OC within its tidal radius. The corner legends in the bottom panel follow Fig. 5.

Fig. 7 Relationships between Log(Ncore/Nouter) and Log(Score/Souter) (Quadrant-wise Fitting). The panels are arranged from left to right as: Q1 corner, the linear fit plot, and Q2 corner. The corner panels’ legends follow Fig. 5.

Fig. 8 Relationships between the morphological stability parameters and the number of members (N) in the Q1 region: Score/Souter (the first two panels from left to right) and Ncore/Nouter (the last two panels). Other panels’ legends follow Fig. 5.

Fig. 9 Relationships between the morphological stability parameters and the number of members (N) in the Q2 region, in the same style as Fig. 8.

4.3 Potential variations in the morphological stabilities across different spatial locations

Our sample OCs are almost embedded within the Galactic disk, as is shown in Fig. 3, and undoubtedly, they are influenced by the gravitational potential at their positions within the disk, with the most significant effect being Galactic tidal forces (Balaguer-Núñez et al. 2020). Therefore, in this section we explore whether the morphological stability of the sample OCs changes in different spatial positions.

Hu et al. (2025b) report that morphological stability parameter I of OCs exhibits a correlation with the Y axis in the Y–Z projection of the heliocentric Cartesian coordinate system. This suggests possible changes in the morphological stability along the Y axis, which aligns with the radial direction under the Galactic Cartesian coordinate system. In this work, the morphology of sample OCs with a low height (|Z|) in the radial direction is, in fact, similar to that in the Y–Z projection, while the one in the tangential direction is almost the same as that in the X–Z projection. Thus, we first investigated potential variations in the morphological stability of sample OCs in the radial and tangential directions.

To this end, we grouped the sample OCs according to their locations in the Galactic Cartesian coordinate system that is centered on the Galactic center, where the negative X axis points toward the Sun, the positive Y axis points in the direction of Galactic rotation, and the positive Z axis points northward from the Galactic plane. The Sun’s position was taken from Xu et al. (2025) as [X, Y, Z] = (−8200.0, 0.0, 14.0) pc. The positional parameters of the sample OCs were taken directly from van Groeningen et al. (2023). To ensure that the radial and tangential sample projections correspond as closely as possible to the morphologies observed in the Y–Z and X–Z projection planes, respectively, we restricted all samples to the range |Z| ≤ 200 pc. We initially divided the solar neighborhood into two key regions: the radial region, including a stripe region with −250 pc ≤ Y ≤ 250 pc and the tangential region representing a range of −8450 pc ≤ X ≤ 7950 pc.

To verify whether the radial and tangential sample projections are almost equal to the morphologies observed in the Y–Z and X–Z projection planes, respectively, we compared the morphological stability parameter I in this study with those calculated by Hu et al. (2025b). By cross-matching, we obtained 36 and 45 common OCs between their samples and our samples in the radial and tangential directions, respectively, as is shown in Fig. 10. We can see from this figure a relatively good agreement between our samples and their samples. But there is still some dispersion around the dashed red lines of Fig. 10. This may be because our sample projection corrected still exhibits a certain tilt angle relative to the X–Z and Y–Z projection planes. Nevertheless, our approach remains effective overall.

Figure 11 shows the distributions of the two subsamples in the X–Y (the Galactic disk plane) and the X–Z planes. The size of the solid circles is proportional to the number of members, with their color being coded by the logarithm of the morphological stability parameter I. This distribution does not show any discernible trend. However, it is certain that the morphological stability of OCs is inevitably influenced by external forces at their spatial locations. We investigated how to measure this influence. We already know that an OC’s morphological stability is directly proportional to the number of its member stars, as demonstrated by the linear relationship presented earlier. If the morphological stability of OCs at certain spatial locations is affected, this will inevitably weaken this linear proportional relationship. Based on this speculation, we performed the linear fitting between the morphological stability parameters I and II and the number of members for the radial and tangential samples separately. The parameters of their linear relationships are compiled in Table 3.

As we expected, the linear relationships mentioned above change. For morphological stability parameters I and II, the radial samples show stronger correlations than the tangential samples. This demonstrates that the morphological stabilities of the sample OCs vary at different spatial positions, in which the influence of external forces on their morphology is variable. The radial samples in this study primarily show us their morphology in the Y–Z projection plane, which is predominantly influenced by shear forces from the disk’s differential rotation. In contrast, what we can see for the tangential samples is their morphology in the X–Z projection plane, which is mainly shaped by tidal forces directed roughly toward the Galactic center. Based on this, it can be concluded that the shear forces generated by the differential rotation possibly exert a weaker influence on the morphological stability of the sample OCs than the tidal forces directed toward the Galactic center.

In addition to the radial and tangential samples primarily influenced by shear forces and tidal forces, respectively, samples at other locations must be subjected to these two external forces simultaneously. To more precisely reveal the influence of other spatial positions on the morphological stability of the sample OCs, we thus divided the samples based on azimuthal angles (e.g., in 15° intervals) centered on the Sun and extracted the samples located at angular intervals of 0° (radial), 15°, 30°, 45°, 60°, 75°, 90° (tangential), 105°, 120°, 135°, 150°, and 165°. The number of sample OCs for most groups is around 100, as is shown in Table 3.

We then first studied the samples along the two directions (45° named by “tilt-right” and 135° named by “tilt-left”), as shown in Fig. 12, because they theoretically are subject to roughly equal tidal and shear forces. If the Galactic disk is symmetrical, the morphological stability of the sample OCs in both directions is equally affected by external forces. To test this hypothesis, we performed linear fits of two morphological stability parameters against the number of members for both samples. The fit parameters are listed in Table 3. We found that the slopes of the morphological stability parameters with the number of members for the tilt-right sample are different from those for the tilt-left sample. Therefore, in general, the morphological stability of the samples in these two directions may be affected to varying degrees, indicating that the external force environment influencing the morphology of sample OCs is asymmetric in these two directions.

We next systematically fit all the other subgroups. All results are also summarized in Table 3, and their trends are clearly displayed in Fig. 13. We found that the slopes of the radial sample in both morphological stability parameters I and II are relatively large within the error range. On the basis of previous speculation, the morphological stability of these sample OCs is relatively stable. However, the slope of the tangential samples is relatively small among all of the samples, suggesting that the morphological stability of the tangential sample is comparatively weaker than that of samples oriented in other directions. Moreover, from the overall trend in Fig. 13, the slopes of samples at different positions vary, indicating that the external force environment on the Galactic disk region around the Sun is actually asymmetric.

Furthermore, we can also see from Fig. 13 that the slopes of the sample OCs within tilt angles of less than 90° are almost smaller than those of other samples. We investigate why this is happening. To explore whether environmental factors are associated with this asymmetry, we plotted the map of our sample on the Galactic disk, as shown in Fig. 14, colored by the logarithm of their morphological stability II. Based on this figure, it can be easily explained that the OC samples within tilt angles of less than 90° are closer to the Galactic bar of the Milky Way than those with tilt angles greater than 90°, and suffer from more disruption of the external force environment. The OCs with tilt angles less than 90° are concentrated in a sector region with a Galactocentric distance of 4–6 kpc and an azimuthal angle of approximately 0°-90°. This region highly overlaps with the orientation of the Galactic bar (Joharle et al. 2026; Kalita et al. 2026) and the peak density area of the Gaia-Sausage-Enceladus merger remnant (Helmi et al. 2018). Recent lustrisTNG50 simulations further demonstrated that major mergers can leave azimuthal offsets in stellar warps that persist even after differential rotation correction, providing a dynamical memory of external perturbations lasting several gigayears (Thulasidharan et al. 2025). Therefore, the combined perturbations from bar torque, merger potential, and local arm shear in the region where low tilt angle OCs are located may continuously strip member stars from their outer regions, leading to a systematically smaller ratio of N_core/N_outer. Consequently, these OCs are subjected to more disruptive external environmental forces. We also found that most samples with a small parameter II are distributed along the local arm, suggesting that the external force environment around the local arm is likely complicated.

It is well known that young OCs are usually not far away from their birthplaces, which were likely close to the spiral arms. Therefore, we specifically studied the slope between morphological stability parameter II and the number of members in OC samples of different ages. Della Croce et al. (2024) point out that 80% of OCs under the age of 30 Myr are in an explosive phase, leading us to speculate that young OCs have poor morphological stability. Hao et al. (2021, see their Fig. 5) mention that 72% of the young OCs are located in the spiral arms. Qiu et al. (2026, see their Fig. 6) conclude that the further an OC is from the Galactic disk, the more compact its center. Therefore, we divided these sample OCs into six age groups (age ≤ 30 Myr, 30 Myr < age ≤ 100 Myr, 100 Myr < age ≤ 200 Myr, 200 Myr < age ≤ 400 Myr, 400 Myr < age ≤ 800 Myr, and age > 800 Myr), as is displayed in Fig. 15, ensuring a balanced number of OCs in each group. These linear fitting parameters of six groups are listed in Table 4, in which the weakest correlation is observed for young OCs with fewer than 30 Myr among them, with a fitting slope of 0.751 ± 0.166 for Log(N_core/N_outer) versus Log(N). The overall trend is illustrated in Fig. 16. This suggests that the gravitational disturbances from their birthplaces have a significant impact on the morphology of the youngest sample OCs. As the OCs age, they are gradually moving from their birthplaces, and the influence of the gravitational disturbances from the birthplaces on their morphology is also getting weaker. As expected, the slopes of the other five groups with ages greater than 30 Myr are, on the whole and within the range of error margins, larger than that of the youngest sample OCs group. In addition, the same holds true for the slopes of morphological stability parameter I against the number of members, as is displayed in Table 4.

The observed slopes of the morphological stability parameters increasing with the number of members as a function of age can be consistently explained from a statistical perspective by Coenda et al. (2025) and Liu et al. (2025). The majority of young OCs with an age of ≤ 30 Myr mostly exist in “grouped” configurations, with intra-group separations of only a few dozen parsecs (Coenda et al. 2025, see the left panel of their Fig. 3). Their overlapping potential wells generate a rapidly time-varying internal tidal field. At the same time, these young clusters are generally located in the local spiral arm region, which is still part of a growing arm where macroscopic turbulence and density waves continuously inject kinetic energy (Liu et al. 2025). Therefore, the combined internal and external disturbances constantly break the stellar field, likely resulting in the poorest morphological stability of our young sample OCs (Elmegreen 2008, 2018). It can also be inferred from the work of (Hu et al. 2025a) that a young pair (ASCC_19 and ASCC_21) is undergoing mutual tidal interaction processes.

The OCs in the 30 Myr < age ≤ 200 Myr range exist mostly as binary clusters, because Coenda et al. (2025) note that the median age for paired OCs is log(Age/yr) = 8.01 ± 0.03, which corresponds to approximately 100 Myr. Although tidal perturbations between the two clusters remain, their strength is significantly weaker than the multidirectional disturbances within young cluster groups. Meanwhile, as these OCs begin to move away from the most turbulent gas layers in the spiral arm (e.g., Binney & Tremaine 2008; Reid et al. 2019), and the external environmental noise also diminishes. Hence, their morphological stability is higher than that of young OCs. For OCs in the 200 Myr < age ≤ 400 Myr range, they have largely transitioned to more isolated states. The external tidal sources are progressively simplified, while internal dynamical processes begin to dominate the cluster evolution. The morphological stability parameters show a further increase, reflecting the ongoing stabilization process. Old sample OCs exhibit the strongest morphological stability, as they mostly exist as single clusters and are located farther from the galactic disk (Coenda et al. 2025). The external tidal sources are simplified to smooth the Galactic disk potential, while internal two-body relaxation rapidly erases any remaining substructure, ultimately forming a stabilized radial density profile. For the oldest OCs (age > 800 Myr), they maintain high morphological stability as single clusters located at significant distances from the galactic plane. However, the slopes show a slight decrease compared to the 400 Myr < age ≤ 800 Myr group, possibly due to the long-term effects of dynamical evolution and mass segregation processes that have been operating over extended timescales (Moreira et al. 2025). These OCs typically reside at high galactic latitudes far from the disk plane (Wu et al. 2024), where the non-axisymmetric gravitational potential generated by the Galactic warp introduces additional kinematic perturbations to their orbital motions (Jónsson & McMillan 2024; Thulasidharan et al. 2025; Reshetnikov et al. 2025). Therefore, as the age group increases, the morphological stability of the sample OCs generally exhibits a stronger trend.

Moreover, Hu et al. (2025b) suggest that there is no disruption for the morphological stability of OCs in the local environment around the Sun. So to investigate the influence of the local environment on our samples, we plotted the distribution of the morphological stability parameter I of left-tilt and right-tilt sample OCs against X, Y, and Z axes, as is shown in Fig. 17. The results show that sample OCs closer to the Sun (marked by a vertical dotted red line) exhibit stronger morphological stability, with the strongest stability observed within approximately 200 pc. This suggests that the solar vicinity may be a relatively “quiet” region, lacking complex additional external forces (e.g., strong gravitational perturbations from other massive structures and intense interactions with GMCs) that severely disrupt OC equilibrium (Rangwal et al. 2025). This result is consistent with the local bubbles mentioned in Gozha & Marsakov (2015). Besides, among the different directions, the morphological stability is generally higher in the Z direction (perpendicular to the Galactic plane) compared to the X and Y directions (within the plane). This is likely due to stronger gravitational restoring forces vertical to the disk, whereas in-plane directions are more susceptible to influences such as spiral arm crossings, tidal forces, and interactions with interstellar clouds.

Fig. 10 Left panel: radial Score/Souter comparison. Right panel: tangential Score/Souter comparison. The comparative samples are drawn from (Hu et al. 2025b), comprising 36 in the radial direction (Y–Z plane) and 45 in the tangential direction (X–Z plane). The dotted red line in each panel represents the y = x reference line, with the error bars indicating the dispersions of Score/Souter.

Fig. 11 Distributions of the radial and tangential sample OCs in the Milky Way: X–Y plane distribution (left panel) and X–Z plane distribution (right panel). The size of the solid circles indicates the number of member stars (N), with their colors being coded by the Log(Score/Souter). The arrow marks the direction of the Galactic center (GC).

Fig. 12 Distributions of the right-tilt and left-tilt sample OCs in the Milky Way: X–Y (left panel), X–Z (right panel), in the same style as Figure 11, except for color.

Fig. 13 Comparative graph of slopes (top panel: Log(Score/Souter) vs. Log(N), bottom panel: Log(Ncore/Nouter) vs. Log(N)) across different groups in Table 3. The inset is for a better understanding of the changing trend of the slopes.

Fig. 14 Morphological stability (Ncore/Nouter) map of sample OCs in our galaxy. The red cross marks the position of the Sun. We used the library available at mw_plot: https://milkyway-plot.readthedocs.io/en/stable/matplotlib_faceon.html#milkyway-bird-s-eye-view

Fig. 15 Distributions of sample OCs with varying ages in the Milky Way (from left to right and then from top to bottom: age ≤ 30 Myr, 30 Myr < age ≤ 100 Myr, 100 Myr < age ≤ 200 Myr, 200 Myr < age ≤ 400 Myr, 400 Myr < age ≤ 800 Myr, age > 800 Myr). The inset in each panel indicates the absolute distance from the Galactic plane. The age parameters used in this work were taken from Cavallo et al. (2024).

Fig. 16 Same as Fig. 13 but those from different age groups in Table 4. The first five age groups are labeled by their median ages, with the last group being marked by its lower age limit (800 Myr).

Fig. 17 Distributions of the morphological stability I Score/Souter in the X, Y, and Z directions for the sample OCs in the left-tilt and right-tilt regions (from top to bottom: X, Y, and Z). Every ten sample OCs in each panel form a bin. Blue scattered points represent the bootstrap (Kawano & Higuchi 1995) median estimate of each bin, with error bars indicating the 95% confidence intervals derived from bootstrap resampling (B = 1000). The vertical red line marks the position of the Sun. A bootstrap methodology was employed to provide robust uncertainty estimates for the median values, given the nonparametric nature of the data distribution.

5 Discussion

The morphological stability derived from the irregular morphology of quantitatively OCs in rose diagrams is a key parameter in this work, which depends on the calculation of the core radius in this approach. Therefore, it is necessary to discuss the systematic biases and sensitivities affecting the determination of the core radius (R_core).

In addition to the fluctuation of the tidal radius within the error range affecting the determination of the core radius, we also considered the influence of starting angle randomness in the rose diagram construction method on the calculation of the core radius. We conducted ten random starting angle experiments to verify the effect of the member star distribution on the core radius calculation. The specific method is as follows: based on the mean tidal radius, R_t, of each sample OC, we fixed R_t as the boundary to reselect its members, calculated their projected coordinates relative to the cluster center, and quantified their 2D spatial distribution through the rose diagram method. During the construction of the rose diagram, a random starting angle, θ_j ~ 𝒰(0°, 360°), between 0° and 360° was independently generated for each iteration to rotate the reference direction of the sector division. Repeating this process ten times, we obtained ten independent estimates of the core radius, R_core, core area, S_core, and outer area, S_outer. With this process, we can plot the sampling rose diagrams of all sample OCs, such as NGC_6791 presented in Fig. B.3. Finally, we compared the two morphological stability parameters obtained from the previous method and this way, as is shown in Fig. 18.

To further comprehensively evaluate the combined impact of tidal radius uncertainty and starting angle randomness on the estimation of morphological parameters, we designed a Monte Carlo experiment integrating the two aforementioned randomization processes. The specific procedure is as follows: based on the mean tidal radius, R_t, of each OC and its standard deviation, σ_Rt, a normal distribution, $\mathcal{N}({\mathrm{R}}_{\mathrm{t}},{\sigma }_{{\mathrm{R}}_{\mathrm{t}}}^2)$, was constructed. In each iteration, k, a tidal radius value, ${\mathrm{R}}_{\mathrm{t}}^{(\mathrm{k})}$, was first independently drawn from this distribution and used as the boundary to reselect cluster member stars, followed by their projected coordinates being calculated. Subsequently, when quantifying their 2D spatial distribution using the rose diagram method, an independent random starting angle, θ_k ~ 𝒰(0°, 360°), was generated to determine the reference direction for sector division. This combined sampling process was repeated ten times, yielding ten independent sets of estimates for the core radius, ${\mathrm{R}}_{\text{core}}^{(\mathrm{k})}$, core area, ${\mathrm{S}}_{\text{core}}^{(\mathrm{k})}$, and outer area, ${\mathrm{S}}_{\text{outer}}^{(\mathrm{k})}$, under the dual randomization of both tidal radius and starting angle. This method aims to simultaneously reflect the potential influence of uncertainties in both boundary delineation and directional definition on the calculation of morphological stability parameters. All results from the combined random sampling for our samples can be obtained, such as the rose diagrams of NGC_6791 shown in Fig. B.4. We also plotted the comparison of the two morphological stability parameters derived from this way and our original approach, as is displayed in Fig. 19. By comparing the distributions and dispersions of the morphological parameters obtained from all methods, as shown in Figs. 18 and 19, we found that the sampling analysis results of all random starting angles are consistent with the relevant conclusions based on the fixed starting angle at 0° within error ranges.

Fig. 18 Comparison of the morphological stability parameters (left panel: Ncore/Nouter, right panel: Score/Souter) derived from two methods. Method A denotes sampling the tidal radius ten times, with Method B sampling the starting angle ten times. Dashed red lines indicate 1:1 correspondence. Error bars represent standard deviations.

Fig. 19 Same as Fig. 18, but Method C samples both the tidal radius and the starting angle ten times.

6 Summary and conclusions

In this paper, we investigated the 2D morphological stability parameters of 1490 OCs in the Galactic spherical coordinate frame by the rose diagram construction and explored their relationships with the number of member stars and their potential variations at different spatial locations by employing the Bayesian fitting method. The data of sample OCs used in this work are sourced from the member stars cross-matched by van Groeningen et al. (2023) from Gaia DR3. Besides, instead of stacking multiple line-of-sight projections, we projected each sample OC only once along its individual line of sight, corrected for distortions in Galactic longitude b, and forced the representation onto the plane that is perpendicular to the Galactic disk. Finally, we discussed the potential influence of systematic biases inherent in the rose diagram method on the results. Our main findings and conclusions are as follows:

1. This work defined for the first time a new morphological stability parameter of N_core/N_outer. Within the 2D projection, both morphological stability parameters are positively correlated with the number of members. This all indicates that within a certain range, the more member stars an OC has, the stronger its binding ability will be, and the denser the structure, the stronger its stability will be.

2. The morphological stability parameter II N_core/N_outer has a higher slope in the Q1 region than in the Q2 region, while the morphological stability parameter I S_core/S_outer shows an inverse trend. This suggests that the morphological stability parameters I and II jointly characterize the structural configuration and evolutionary stage of OCs.

3. The two morphological stability parameters exhibit lower values in the tangential OCs compared to the radial region, indicating that the influence of the tidal forces towards the Galactic center on the tangential sample OCs is likely stronger than that of the shear force resulting from the Galactic differential rotation on the radial sample OCs.

4. The slope of the morphological stability parameters against the number of members in the tilt-right region is generally lower than that in the tilt-left region. This indicates that the influence of external environments on the sample OCs is asymmetric, which may be related to their proximity to the Galactic center. The closer they are to the Galactic center, the stronger the galactic tidal force they experience.

5. The slopes of the morphological stability of the sample OCs against the number of members show a significant increase in the trend from young to old. This implies that compared to older sample OCs, younger sample OCs experience more severe external environmental disturbances.

6. The highest morphological stability parameters observed in OCs near the Sun suggest that external perturbations in this region are relatively weak, implying the existence of a locally quiescent environment.

We established a connection between the internal evolution and external environment of the sample OCs by their morphological stability parameters. Specifically, our results reveal that a lower morphological stability parameter is indicative of a more complex and dynamically active evolutionary state for the OCs. To deepen our understanding, future studies could perform a stratified analysis considering OC age and mass to examine whether morphological stability trends of the OCs still hold across different subpopulations in the local environment. Extending the sample to greater distances, particularly toward spiral arms or the Galactic center, would help identify critical scales at which the morphological stability changes markedly, shedding light on the global influence of Galactic structure on cluster evolution.

Data availability

The full version of Table 1 is available at the CDS via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/708/A238

Acknowledgements

We thank an anonymous reviewer for valuable suggestions that have improved our methodology, results, and presentation. This work is supported by the National Natural Science Foundation of China (NSFC) under grant Nos. 12303037 and 12573035, the Natural Science Foundation of Sichuan Province (No. 2026NSFSC0742), and the Fundamental Research Funds of China West Normal University (CWNU, No. 23KE024). Qingshun Hu would like to acknowledge the financial support provided by the China Scholarship Council program (Grant No. 202308510136). This study has made an indirect use of Gaia data, through large catalogues of open clusters and stellar parameters data published in the third data release. Gaia is operated by the European Space Agency (ESA) (https://www.cosmos.esa.int/gaia). The preparation of this work has made extensive use of Topcat (Taylor 2005), of the Simbad and VizieR databases at CDS, Strasbourg, France, and of NASA’s Astrophysics Data System Bibliographic Services.

References

Appendix A Bayesian linear regression with MCMC

We implemented a uniform Bayesian linear-regression framework to quantify empirical relations between any two observables (X, Y) that are both subject to measurement errors and to an unknown intrinsic scatter. The method is written in Python and based on the affine-invariant ensemble sampler emcee (Foreman-Mackey et al. 2013), also employed by Anders et al. (2017) and Hu & Soubiran (2025).

Appendix A.1 Model parameterization

The linear model is $$\mathrm{Y}=m\mathrm{X}+b+\epsilon ,\epsilon \sim \mathcal{N}(0,{\sigma }^2),$$(A.1)

where m is the slope, b the intercept, and σ the intrinsic scatter around the linear relation.

Appendix A.2 Likelihood with heteroscedastic errors

Each datum carries Gaussian uncertainties (σ_X,i, σ_Y,i). Following Hogg & Foreman-Mackey (2018), the total variance for point i is $${\sigma }_{\mathrm{t}\mathrm{o}\mathrm{t},i}^2={\sigma }_{\mathrm{Y},i}^2+{\left(m{\sigma }_{\mathrm{X},i}\right)}^2+{\sigma }^2.$$(A.2)

The log-likelihood is $$\ln \mathcal{L}(m,b,\sigma )=-\frac{1}{2}\sum _{i=1}^N[\frac{{({\mathrm{Y}}_i-m{\mathrm{X}}_i-b)}^2}{{\sigma }_{\text{tot\hspace{0.17em}},i}^2}+\ln \left(2\pi {\sigma }_{\text{tot\hspace{0.17em}},i}^2\right)].$$(A.3)

Appendix A.3 Prior distributions

Uninformative uniform priors were adopted: $$m\sim \mathcal{U}(m_{min},m_{max}),$$(A.4) $$b\sim \mathcal{U}(b_{min},b_{max}),$$(A.5) $$\sigma \sim \mathcal{U}(0,{\sigma }_{max}).$$(A.6)

These ranges are chosen physically wide enough to contain any realistic solution, and the values used in this work are listed in Table A.1.

Appendix A.4 Posterior distribution

According to Bayes’ theorem, our goal is to maximize the posterior distribution. By combining the likelihood function and the prior distribution, we obtained the posterior distribution of the parameters m and b. The logarithmic posterior distribution is $$L_p=-\frac{1}{2}\sum _{i=1}^N[\ln {\sigma }_i^2+\frac{{(y_i-{\hat{y}}_i)}^2}{{\sigma }_i^2}-2\ln w_i],$$(A.7)

• σ_i the error of the i-th data point,

• w_i the weight of each data point,

• ${\hat{y}}_i$ the predicted value based on the regression model.

Appendix A.5 MCMC implementation

We sampled the posterior with emcee using 24 walkers evolved for 10 000 steps, with the first 500 discarded as burn-in. Convergence was checked by eye and via the Gelman–Rubin statistic ($\hat{R}<1.05$). The posterior chains yield the marginalized 16th, 50th, and 84th percentiles, which define the best-fit value and its 1σ uncertainty. Visualization is performed with the corner package (Foreman-Mackey 2016).

Appendix A.6 General applicability

The same pipeline is applied to any (X, Y) pair of interest. Here, we fit three linear relationships (log₁₀(N_core/N_outer) versus log₁₀(N), log₁₀(S_core/S_outer) versus log₁₀(N), and log₁₀(N_core/N_outer versus log₁₀(S_core/S_outer)) with this method without the modification of the likelihood or priors.

Appendix B Corner panels and Bayesian fit panels for each group

Fig. B.1 Relationships between the morphological stability parameters (Ncore/Nouter and Score/Souter) and the number of members (N) of the sample OCs in the radial and tangential region, in the same style as Fig. 5.

Fig. B.2 Same as Fig. B.1 but for the sample OCs in the “left_tilt” and “right_tilt” regions.

Fig. B.3 Rose diagrams of NGC_6791. These were obtained by randomly sampling the starting angle ten times, but with a fixed tidal radius Rt.

Fig. B.4 Same as Fig. B.3 but with ten random sampling iterations conducted simultaneously for the starting angle and tidal radius Rt.

All Tables

All Figures

Fig. 1 Distribution of NGC_6791 in the Galactic coordinate system (left panel) and its flattened projection distribution (right panel). The scattered blue points in each panel represent the member stars of the OC, with the gray bars being the errors of the position coordinates of its members. The orange circles correspond to the tidal radii obtained from 100 Monte Carlo sampling runs. The dashed black circle represents the mean value of the 100 results.

In the text

Fig. 2 Rose diagram of the NGC_6791 on the 2D corrected projection, derived from one of ten samplings. The rose-colored sectors represent the normalized directional radial radius values obtained by the Rose diagram method. Each sector spans 30°, with the first sector starting at 0° corresponding to the horizontal line on the right of this diagram. The second sector begins at a position 30° counterclockwise from the first sector, and subsequent sectors follow this pattern. The semitransparent gray sectors are the average extension levels in each direction, with the scattered green points indicating the member stars.

In the text

	Fig. 3 Map of the sample OCs in the Galactic coordinate system. The size of each color-filled circle is proportional to the number of members of each OC, with its color being coded by morphological stability parameter I.
In the text

	Fig. 4 Histogram of the areas of the core (Score, blue) and the outer (Souter, red) of our sample OCs. Each set of data was fit with a normal distribution (solid colored lines), with the mean positions being marked with vertical dashed lines.
In the text

Fig. 5 Relationships between the morphological stability parameters (Score/Souter: the first two panels from left to right; Ncore/Nouter: the last two panels) and the number of members (N). For ease of visualization, a logarithmic scale was adopted. Blue scatters in linear fitting panels represent our sample OCs, with the gray bars denoting their parameter errors. The faint gray lines show 1000 MCMC samples, with the thick red lines (the median values) and the shaded yellow bands (1σ) corresponding to our best parameter estimates. The two corner plots display the marginal posterior distributions (diagonal) and the 2D joint posterior contours (below diagonal) of the fit parameters (m, b, σ). Contours enclose the 1σ and 2σ confidence regions, and the median together with the 16th and 84th percentiles are quoted in each panel.

In the text

Fig. 6 Relationship between Score/Souter and Ncore/Nouter for the whole sample through a logarithmic scale. The dotted black lines in the top panel correspond, respectively, to the zero positions on the horizontal and vertical axes, representing the scenario where both morphological stability parameters I and II are equal to 1. The color bar is coded by the number of member stars of each sample OC within its tidal radius. The corner legends in the bottom panel follow Fig. 5.

In the text

	Fig. 7 Relationships between Log(Ncore/Nouter) and Log(Score/Souter) (Quadrant-wise Fitting). The panels are arranged from left to right as: Q1 corner, the linear fit plot, and Q2 corner. The corner panels’ legends follow Fig. 5.
In the text

	Fig. 8 Relationships between the morphological stability parameters and the number of members (N) in the Q1 region: Score/Souter (the first two panels from left to right) and Ncore/Nouter (the last two panels). Other panels’ legends follow Fig. 5.
In the text

	Fig. 9 Relationships between the morphological stability parameters and the number of members (N) in the Q2 region, in the same style as Fig. 8.
In the text

	Fig. 10 Left panel: radial Score/Souter comparison. Right panel: tangential Score/Souter comparison. The comparative samples are drawn from (Hu et al. 2025b), comprising 36 in the radial direction (Y–Z plane) and 45 in the tangential direction (X–Z plane). The dotted red line in each panel represents the y = x reference line, with the error bars indicating the dispersions of Score/Souter.
In the text

	Fig. 11 Distributions of the radial and tangential sample OCs in the Milky Way: X–Y plane distribution (left panel) and X–Z plane distribution (right panel). The size of the solid circles indicates the number of member stars (N), with their colors being coded by the Log(Score/Souter). The arrow marks the direction of the Galactic center (GC).
In the text

	Fig. 12 Distributions of the right-tilt and left-tilt sample OCs in the Milky Way: X–Y (left panel), X–Z (right panel), in the same style as Figure 11, except for color.
In the text

	Fig. 13 Comparative graph of slopes (top panel: Log(Score/Souter) vs. Log(N), bottom panel: Log(Ncore/Nouter) vs. Log(N)) across different groups in Table 3. The inset is for a better understanding of the changing trend of the slopes.
In the text

	Fig. 14 Morphological stability (Ncore/Nouter) map of sample OCs in our galaxy. The red cross marks the position of the Sun. We used the library available at mw_plot: https://milkyway-plot.readthedocs.io/en/stable/matplotlib_faceon.html#milkyway-bird-s-eye-view
In the text

Fig. 15 Distributions of sample OCs with varying ages in the Milky Way (from left to right and then from top to bottom: age ≤ 30 Myr, 30 Myr < age ≤ 100 Myr, 100 Myr < age ≤ 200 Myr, 200 Myr < age ≤ 400 Myr, 400 Myr < age ≤ 800 Myr, age > 800 Myr). The inset in each panel indicates the absolute distance from the Galactic plane. The age parameters used in this work were taken from Cavallo et al. (2024).

In the text

	Fig. 16 Same as Fig. 13 but those from different age groups in Table 4. The first five age groups are labeled by their median ages, with the last group being marked by its lower age limit (800 Myr).
In the text

Fig. 17 Distributions of the morphological stability I Score/Souter in the X, Y, and Z directions for the sample OCs in the left-tilt and right-tilt regions (from top to bottom: X, Y, and Z). Every ten sample OCs in each panel form a bin. Blue scattered points represent the bootstrap (Kawano & Higuchi 1995) median estimate of each bin, with error bars indicating the 95% confidence intervals derived from bootstrap resampling (B = 1000). The vertical red line marks the position of the Sun. A bootstrap methodology was employed to provide robust uncertainty estimates for the median values, given the nonparametric nature of the data distribution.

In the text

	Fig. 18 Comparison of the morphological stability parameters (left panel: Ncore/Nouter, right panel: Score/Souter) derived from two methods. Method A denotes sampling the tidal radius ten times, with Method B sampling the starting angle ten times. Dashed red lines indicate 1:1 correspondence. Error bars represent standard deviations.
In the text

	Fig. 19 Same as Fig. 18, but Method C samples both the tidal radius and the starting angle ten times.
In the text

	Fig. B.1 Relationships between the morphological stability parameters (Ncore/Nouter and Score/Souter) and the number of members (N) of the sample OCs in the radial and tangential region, in the same style as Fig. 5.
In the text

	Fig. B.2 Same as Fig. B.1 but for the sample OCs in the “left_tilt” and “right_tilt” regions.
In the text

	Fig. B.3 Rose diagrams of NGC_6791. These were obtained by randomly sampling the starting angle ten times, but with a fixed tidal radius Rt.
In the text

	Fig. B.4 Same as Fig. B.3 but with ten random sampling iterations conducted simultaneously for the starting angle and tidal radius Rt.
In the text