© The Authors 2025

1 Introduction

Stars are the fundamental units of galaxies because they control the formation and evolution of their host galaxies. Most stars form in clustered environments (Lada & Lada 2003) within giant molecular clouds (GMCs) through fragmentation and hierarchical collapse (Larson 1981; Elmegreen & Scalo 2004; Mac Low & Klessen 2004; McKee & Ostriker 2007; Kruijssen 2012). Since stars in a given cluster form from the same molecular cloud and roughly at the same time, the stars of a cluster have similar ages and chemical compositions. Thus, star clusters hold unique information regarding the kinematics and evolution of star-forming regions (Karnath et al. 2019; Kuhn et al. 2019; Pang et al. 2020) as well as the origins and dynamics of host galaxies (Kroupa 2008; Ferraro et al. 2016).

Star clusters are ideal laboratories for understanding the origin of the initial mass function (IMF) (Elmegreen 2006; Bastian et al. 2010; Damian et al. 2021, 2023a), disc evolution, and planet formation scenarios in clustered environments (Piskunov et al. 2018; Concha-Ramírez et al. 2019; Reiter & Parker 2022; Damian et al. 2023b), as well as mass segregation (Meylan 2000; de Grijs et al. 2002; Olczak et al. 2011; Nony et al. 2021). Many clusters are also associated with massive O and B-type stars, which disrupt their surrounding environment, controlling the star formation ecosystem through the feedback effect of energetic stellar winds and H II regions (Samal et al. 2012, 2014; Jose et al. 2016; Das et al. 2017, 2021; Pandey et al. 2022). The dynamics of clusters and their formation processes regulate the global outcome of star formation efficiency and star formation rate of a molecular cloud (Dale et al. 2012, 2013; Walch et al. 2015; Patra et al. 2022).

Thus, clusters serve as the black boxes of the parental molecular cloud, encapsulating all of its star formation history. Analysing young stellar clusters represents the initial step towards understanding how individual stars and protoplanetary discs evolve in a clustered environment, and how the outcomes of star formation processes in a GMC depend on the stellar properties, kinematics, and feedback processes of associated clusters.

Therefore, the accurate detection of the member population of a star cluster becomes crucial to initiate such analyses. Membership analyses of clusters and star-forming complexes typically involve various techniques (Das et al. 2023, and references therein), each with its own mass limit determined by the sensitivity of the respective dataset. For instance, the Gaia data (Gaia Collaboration 2016) have revolutionised the identification of clusters and their members in the solar neighbourhood (Zucker et al. 2023). Such studies utilise the successful application of various machine learning (ML) algorithms (e.g. Gao 2018a,b; Das et al. 2023; Gupta et al. 2024). With Gaia, identifying high-confidence members is feasible using astrometric parameters such as parallax and proper motions. However, the sensitivity of Gaia does not extend to the detection of very low-mass objects, including brown dwarfs (<0.08 M_⊙) in relatively distant star-forming complexes (>500 pc, age ∼2 Myr, A_V ∼ 2 mag; Baraffe et al. 2015) (e.g. Damian et al. 2023a).

In this study, we focus on conducting a membership analysis of IC 1396, one of the most dominant Galactic bubbles powered by massive stars of the Trumpler 37 (Tr37) cluster (Trumpler 1930), using deep optical data from the Subaru-Hyper Supreme Cam (HSC) and applying ML techniques. This represents the deepest optical photometric data available for this complex to date. Following our previous membership analysis of this complex using Gaia-DR3 (Das et al. 2023), this work aims to identify more low-mass member stars within this star-forming complex.

This work is arranged as follows. We introduce the star-forming complex around Tr 37, the H II region IC 1396, in Section 2, along with the details of the membership analysis already carried out for this complex in the literature. The observational details and data reduction of Subaru-HSC data are given in Section 3, including the details of complete catalogue preparation and its properties, along with details regarding archival datasets used in this work. Section 4 presents the properties of images obtained from the HSC observations. We carry out the membership analysis in Section 5 with the details of the ML techniques implemented in this study. Results for the member population, their properties, and sub-cluster detection are discussed in Section 6. Using the detected member population, we discuss the star formation history of this complex, IC 1396, in Section 7. We summarise our work in Section 8.

2 Details of the star-forming complex IC 1396

2.1 IC 1396

IC 1396 (S131; Sharpless 1959) is a textbook example of a classical H II regions (Fig. 1), with a simple, circular morphology. This is one of the best-studied H II regions of the Milky Way Galaxy. IC 1396 is located in the Cepheus OB2 region (de Zeeuw et al. 1999) and has relatively low (A_V<5) foreground reddening (Sicilia-Aguilar et al. 2005; Getman et al. 2012; Nakano et al. 2012), associated with the cluster Tr 37 (Platais et al. 1998; Patel et al. 1995, 1998). This star-forming complex is ionised by massive multiple stellar systems (e.g. HD 206267; Stickland 1995) of spectral type O5–O9V, located at its near centre (Peter et al. 2012; Maíz Apellániz & Barbá 2020).

Several distance estimations have been made for the complex from several previous works. Based on main-sequence (MS) fitting, Contreras et al. (2002) obtained a distance of 870 pc, and later Sicilia-Aguilar et al. (2019) estimated a distance of ${945}_{-73}^{+90}$ pc based on Gaia-DR2. The two most recent distance estimations are 925 ± 73 pc (Pelayo-Baldárrago et al. 2023) based on Gaia-EDR3 and 917 ± 2.7 pc (Das et al. 2023) based on Gaia-DR3. Within the errors, all the distances fairly agree with each other. In this work, we use the distance of 917 pc obtained by Das et al. (2023).

This star-forming complex is relatively young, with an estimated age of 2–4 Myr (Sicilia-Aguilar et al. 2005) based on spectroscopically identified members and some pre-main sequence (PMS) tracks. The recent analysis based on Gaia-EDR3 and Gaia-DR3 also yields a similar average age for this complex (Pelayo-Baldárrago et al. 2023; Das et al. 2023).

IC 1396 is one of the dominant feedback-driven star formation complexes in the Milky Way. The feedback-driven activity of the central massive star(s) in the complex leads to the formation of several exciting structures such as bright-rimmed clouds (BRCs; Sugitani et al. 1991), fingertip structures, and elephant trunk structures. These structures are distributed in and around the complex (Schwartz et al. 1991; Froebrich et al. 2005; Saurin et al. 2012). The prominent BRCs of the complex IC 1396 A and IC 1396 N are often considered as the best examples of feedback-driven star formation (Sicilia-Aguilar et al. 2004, 2006a; Getman et al. 2007; Choudhury et al. 2010; Sicilia-Aguilar et al. 2013; Panwar et al. 2014; Sicilia-Aguilar et al. 2014, 2019).

In Fig. 1, we present the entire star-forming complex with a WISE 22 μm image. The energetic ultraviolet (UV) radiation and stellar wind of the central massive star(s) shape the complex, which is seen as a cavity of radius ∼1.5^∘ clearing up most of the gas and dust and leading to the formation of BRCs, fingertip structures, and elephant trunk structures around the boundary. By correctly identifying the stellar members, we can better understand the star formation history of the complex. Many studies have already been carried out to find the associated young stellar object (YSO) populations of the complex based on various datasets. In this work, we identify the low-mass member stars up to the brown dwarf using the much deeper Subaru-HSC optical catalogue.

2.2 Members from previous studies

This section briefly describes details of the member populations of those past studies. The population of IC 1396 consists of members identified using several datasets and methodologies, which is summarised in Table A.1. All the different studies have covered different regions of the star-forming complex. Along with the reference, we provide the method of the work, the area covered by the individual studies, and the number of stars retrieved by the individual works. A total of 2178 member stars were detected towards IC 1396 from all the non-Gaia based studies in the literature.

Gaia-based member populations have also been carried out for this complex. Cantat-Gaudin et al. (2018) conducted a membership analysis of a large sample of 1229 clusters of the Milky Way using Gaia-DR2 data and detected 460 candidate member stars towards IC 1396. Pelayo-Baldárrago et al. (2023) identified an additional 334 sources from Gaia-EDR3 catalogue. The work of Das et al. (2023) identified 1243 stars as members of the star-forming complex from Gaia-DR3 data. All the Gaia-based studies have revealed 1468 member stars of IC 1396.

The most recent study by Gupta et al. (2024) identifies brown dwarfs within the central 22^′ radius of the complex. In this work, we obtain the deepest membership catalogue of the complex within the circular region of radius 1.5^∘, as is shown in Fig. 1. The region of interest in this work is similar to the area covered in several previous studies (Barentsen et al. 2011; Nakano et al. 2012; Das et al. 2023). In Section 5.8, we compare the members identified in this work with those in the literature.

3 Observation, data reduction, and catalogue preparation

3.1 Details of observations

The Subaru Telescope (Iye et al. 2004) is an 8.2 m class telescope operated by the National Astronomical Observatory of Japan (NAOJ). The Hyper Supreme Cam¹ (HSC; Miyazaki et al. 2018; Komiyama et al. 2018; Kawanomoto et al. 2018; Furusawa et al. 2018; Huang et al. 2018; Coupon et al. 2018) is a 1.77 deg² mosaic CCD camera mounted at the prime focus of the telescope. The HSC consists of a total of 116 2k × 4 k CCDs, out of which 104 CCDs cover a 1.5^∘ field of view in diameter with a pixel scale of 0.17′′. A single exposure of the HSC is termed a ‘visit’, which is a pair of single exposures of the same pointing processed at the earliest stages. The datasets are identified by the combination of the visit ID and CCD ID.

The wide-field imaging of the region towards IC 1396 was conducted using Subaru-HSC on 17 September 2017 (PI: Jessy Jose, S17B0108N) in excellent imaging conditions (seeing ∼0.5–0.7^′′ and 1.0 ⩽ airmass ⩽ 2.0). This observation comprised two pointings, each with a radius of 0.75^∘ in three filters (r₂, i₂, and Y)². The two pointings of these observations are depicted in Fig. 1. Details regarding the number of frames and integration time in each filter are provided in Table 1.

3.2 Data reduction

The observed raw data was downloaded from STARS³ (Subaru Telescope archive system). We processed the observed raw data using the HSC Pipeline (HSCPIPE⁴) version 6.7. Further details on HSCPIPE V6.7 can be found in the description of data reduction in Gupta et al. (2021) for Cygnus OB2 complex. This data reduction followed the same procedure to obtain the final catalogue for the IC 1396 complex. In the following, we briefly explain the reduction process of HSCPIPE.

The complete data reduction process with the pipeline consists of multiple steps: (1) single-visit processing, (2) joint calibration, (3) co-addition, and (4) co-add processing and/or multi-band analysis. During the initial step of single-visit processing, the pipeline corrects each science frame by performing overscan subtraction, bias correction, dark subtraction, flat-fielding, and fringe subtraction. The pipeline also identifies and interpolates defects such as bad pixels, cross-talk, and saturated pixels through instrument signature removal (ISR). Following this, both the astrometry and photometry calibration processes commence. For the calibration process, HSCPIPE first selects point sources above a 50σ threshold and then performs calibration using the Pan-STARRS DR1 PV3 reference catalogue. By default, HSCPIPE employs the ‘Optimistic B’ matching algorithm for this process, which is suitable for low-density fields such as extragalactic regions. However, it performs poorly for high-density regions (Gupta et al. 2021). Given that our observed region is a highly dense Galactic star-forming region, we opt for the ‘Pessimistic B’ algorithm.

Then the process of co-addition is performed to create a deeper image from multiple images or visits. In this step, the pipeline warps the visits to the SkyMap co-ordinate system (warp) and then co-adds all visits together using the WCS and flux scale determined by mosaicking. Additionally, this step applies sky correction, writing a new background model. Following co-addition, HSCPIPE conducts multiband analysis to generate the final photometric catalogue for each band. During this stage, the pipeline extracts sources from the co-added images using a 5σ threshold value and retrieves photometry catalogues from the co-added images in each filter. To facilitate this process, footprints, which are regions above the threshold level and consist of many discrete sources, are generated. Following the methodology outlined in Gupta et al. (2021), we maintain the footprint size at 10¹⁰ pixels for the i₂ and Y filters, and 10¹¹ pixels for the r₂ filter to ensure maximum detection in the presence of nebulosities. Subsequently, we follow the same point source selection procedure by applying certain flags, as described by Gupta et al. (2021). Details regarding the catalogue flags are provided in Bosch et al. (2018), and Gupta et al. (2021) explains how applying certain flags removes spurious sources from the final catalogue. To prepare the final catalogue, we retain sources with magnitude errors <0.2, ensuring that more than 90% of sources lie within that error value.

To compile the final HSC source catalogue of IC 1396, we merged the individual photometric catalogues generated in the three bands. Our merged catalogue contains a total of 859 867 sources with photometry in all three bands. The final catalogue tends to saturate at approximately 14 mag in the r₂ and i₂ bands and around 12 mag in the Y band, similarly for the Cygnus OB2 complex (Gupta et al. 2021). As is depicted in Fig. 1, the HSC observations do not cover the entire 1.5^∘ circular region (as indicated by the dashed white circle) of the star-forming complex IC 1396. To account for stars in the region not covered by HSC observations and stars above the saturation limit of HSC, we incorporate stars from the Pan-STARRS catalogue, the details of which are provided in the following section. It is important to note that the inclusion of stars from the Pan-STARRS catalogue enhances the detection of member populations beyond the HSC survey area. We like to mention that preliminary results on the central part of the complex using this deep HSC photometry has been published in Gupta et al. (2024), which reports the presence of low-mass objects such as brown dwarfs. Further details regarding membership analysis of the entire complex are explained in subsequent sections.

3.3 Data quality and completeness

In Fig. 2, we present the magnitude distribution relative to error for the three bands. This plot illustrates two distributions primarily due to the photometry obtained from the long and short exposure frames. We also assess the quality of the retrieved photometric catalogue.

In Fig. 3, we depict the histogram distributions of the photometric magnitudes in individual bands. The limiting magnitudes are 27.6, 26.6, and 25.2 in the r₂, i₂, and Y bands, respectively. The HSC observations of the Cygnus OB2 complex also exhibit similar magnitude limits (Gupta et al. 2021). We determine the 90% completeness limits to be 24, 22, and 21 mag in the r₂, i₂, and Y bands, respectively.

3.4 Archival datasets

In this study, we have additionally used datasets from other archives within the circular region of radius 1.5^∘ centred at α = 21: 40: 39.28 and δ = +57: 49: 15.51 (as indicated by the dashed white circle in Fig. 1). Below, we present brief details of these datasets.

3.4.1 Optical data from Pan-STARRS

We use optical data in the r₂, i₂, z, and Y bands⁵ from the archive of the Pan-STARRS1 survey (PS1). For detailed information regarding the survey and data release, we refer to Chambers et al. (2016) and Flewelling et al. (2020). We downloaded the PS1 catalogue from the Vizier online data archive⁶ (Chambers & et al. 2017) within the circular region of radius 1.5^∘ (refer to Fig. 1).

For our analysis, we retained only the data of good quality based on specific criteria. Firstly, all stars must have photometric values in all four bands r₂, i₂, z, and Y bands. Additionally, the quality flag⁷ of all stars should fall within the range of 16 to 64, indicating good-quality stars. To ensure stars with a good photometric quality, we selected stars with an error in all bands of less than 0.2 and a photometric magnitude difference⁸ in PSF-mag and Kronmag of less than 0.05. Applying these constraints enabled us to obtain 800 125 stars of good photometric quality for further analysis.

To merge the PS1 catalogue with the HSC catalogue, we initially converted the PS1 stars into the HSC photometric system using the conversion equations provided in Gupta et al. (2021). During this merging process, we prioritised the photometry from HSC for stars that are common in both the HSC and PS1 catalogues. Additionally, we include additional stars from PS1 in regions not covered by the HSC observations. It is worth noting that the completeness of that region tends towards the brighter end due to Pan-STARRS photometry. As a result of this merging process, we incorporate 462 252 additional stars from PS1 into the HSC catalogue. Consequently, the combined catalogue encompasses 1 322 119 stars within the IC 1396 complex. For simplicity, henceforth, we refer to the complete HSC plus PS1 catalogue as the HSC catalogue.

3.4.2 NIR data from UKIDSS and 2MASS

In this study, we also use near-infrared (NIR) datasets in the JHK⁹ bands from the archives of the UKIRT Infrared Deep Sky Survey (UKIDSS; Lawrence et al. 2007) and the Two Micron All Sky Survey (2MASS; Skrutskie et al. 2006). We obtained stars from both catalogues within the same circular region of radius 1.5^∘.

We downloaded the UKIDSS¹⁰ 11PLUS Galactic Plane Survey (GPS) data. This data has a resolution of approximately 1^′′ and a typical limiting magnitude of K=18.4 mag (Hewett et al. 2006), with a saturation limit of K ∼ 13 mag (Lucas et al. 2008). To ensure that only stars of good quality are retained, we followed the criteria discussed in Lucas et al. (2008) and Damian et al. (2021). We filtered out sources classified as ‘noise’ in the catalogue and retained those with mergedClass values of −1 or −2, indicating a star or probable star, respectively. Additionally, to account for duplicate catalogue entries, we used the attribute PriOrSec. We retained stars with a magnitude in the K band greater than 13 mag.

To account primarily for the brighter end of the UKIDSS catalogue, we use NIR data in the JHK bands from the archive of 2MASS and merged it with the UKIDSS catalogue. Prior to merging, we converted the 2MASS photometry into the UKIDSS photometry system using the conversion equations provided by Carpenter (2001). For simplicity, we refer to this combined NIR catalogue as the UKIDSS catalogue.

There are 802 100 common stars between the HSC and UKIDSS catalogues¹¹ with photometry magnitudes in six bands. For simplicity in further analysis, we refer to this as the HSC+UKIDSS catalogue, which we utilise for membership analysis. Including UKIDSS data enhances the effectiveness and reliability of detecting member stars within the star-forming complex.

4 HSC image of IC 1396

In Fig. 4, we present the colour-composite image of the entire region observed by Subaru-HSC. The renowned BRCs of the complex, including IC 1396 A and IC 1396 N, are visible in the image. The image unveils the diffuse nebulosity and distinct, bright rims of the BRCs, predominantly in the r₂ band. Additionally, sharp nebular features are observable in the peripheral regions of the H II region.

Fig. 5 illustrates the colour-composite image depicting the area surrounding the BRC IC 1396 A. The signatures of feedback-driven activity within the star-forming complex are discernible across all the HSC bands, with the r₂ band primarily highlighting the head of IC 1396 A alongside the extended tail of the BRC. The morphology of the BRC in r₂ corresponds to the WISE 22 μm infrared image, revealing both the head and tails of the structure (refer to Fig. 1). Specifically, the r₂ band is ideal for examining the signatures of feedback and shock-driven activity, as this band covers the H_α line (656.3 nm), which is a tracer of ionised gas (Pandey et al. 2022).

In the close-up view of the BRC head (Fig. 5), sharp diffuse emissions are prominently visible around the head, accompanied by numerous small-scale filamentary structures behind the sharp, bright rim of the BRC head. The formation of the head in IC 1396 A arises from combined feedback-driven activity, primarily influenced by the massive star(s) HD 206267 and the winds from the Herbig Ae star V 390 Cep. A distinctive cavity surrounding the Herbig Ae star is observed in HSC images, attributed to local feedback-driven activity, wherein the energetic radiations from the Herbig Ae star disperse the gas, resulting in a void visible in the HSC images. Such a property of the Herbig Ae star has also been observed previously, where it was reported that, due to feedback-driven activity, this intermediate-mass star is creating an ionisation hole around its vicinity (Sicilia-Aguilar et al. 2019).

A Class 0 star was identified within the BRC head through Herschel PACS data, supporting the notion of multi-episodic star formation activity within the BRC. Sicilia-Aguilar et al. (2014) speculated that the Class 0 star might have been triggered via radiative-driven implosion (RDI) induced by the massive star(s) HD 206267. However, the Class 0 star remains undetected in HSC images because it is deeply embedded, as the HSC images primarily capture diffuse emissions from the warm dust surrounding the BRC head.

In Fig. B.1, we display the colour-composite images of several intriguing regions within the star-forming complex. We showcase the central part of the complex, focusing on the area surrounding the massive star(s) HD 206267. As anticipated, this region is densely populated, with the central massive star(s) exerting significant influence over the surrounding star formation activity. Notably, this central area hosts the prominent central cluster of IC 1396. Additionally, this figure highlights the BRC complexes IC 1396 N, SFO 35, and IC 1396 G.

In addition to point sources, the HSC images capture the bright diffuse emissions emanating from the BRCs. The HSC r₂ image also unveils intriguing morphologies, including diffuse nebulae, filaments, and proplyds. The filter coverage of r₂ band includes the H_α emission line at 656.3 nm, making it well-suited for tracing proplyds.

5 Membership analysis

This study utilises the ML algorithm to identify the low-mass membership population within the IC 1396 complex, using the HSC+UKIDSS catalogue. Das et al. (2023) provided a concise review of membership analysis methodologies found in the literature and determined the member population of this complex from Gaia-DR3 using RF classifier of the ML algorithm.

The random forest (RF; Breiman 2001; Pedregosa et al. 2012) classifier is part of the supervised ML algorithms, and due to its efficiency, RF classifier has been widely used in astronomy (Dubath et al. 2011; Brink et al. 2013; Liu et al. 2017; Lin et al. 2018; Plewa 2018; Gao 2018b,a; Mahmudunnobe et al. 2021; Gupta et al. 2024). A proper training set is essential for this ML algorithm to work efficiently. Additionally, the parameters used to filter the member population from the field population play a crucial role in such analyses.

In this study, we use the RF classifier to obtain the membership population of the whole IC 1396 complex. We describe the whole procedure in the following.

5.1 Dataset preparation for ML

5.1.1 Preparing the test dataset

Before proceeding with the entire ML procedure, it is sensible to filter out as many contaminated field stars as possible from the HSC+UKIDSS catalogue. Utilising Gaia-DR3 data can aid in the removal of these contaminated stars. Previous Gaia-based studies of this region have demonstrated that the member population of the IC 1396 complex typically falls within the parallax range of ∼0.8–1.6 mas (Sicilia-Aguilar et al. 2019; Pelayo-Baldárrago et al. 2023; Das et al. 2023). Therefore, stars lying outside of this parallax range are likely contaminants and should be excluded.

Initially, we obtained the Gaia-DR3 catalogue within the area of interest. To ensure the inclusion of only stars with good photometric quality, we set a criterion of G ⩽ 19 mag. Among these stars, 186485 exhibit parallax values outside the range of ∼0.8–1.6 mas. These stars are presumed to be contaminants, and their counterparts in the HSC+UKIDSS catalogue can be identified and removed. Subsequently, we identified 164 505 matches within the HSC+UKIDSS catalogue and eliminated them. After the removal of these contaminated stars, we are left with 637 595 stars in the HSC+UKIDSS catalogue, which constitutes the dataset used to derive the member population of IC 1396.

5.1.2 Preparing the training dataset

A suitable training set is crucial for effectively using the RF classifier algorithm. Das et al. (2023) derived the membership catalogue of 1243 stars for the entire complex from Gaia-DR3 data. In this study, we use the HSC+UKIDSS counterparts of the 1243 Gaia-based members. These counterparts constitute the member population of our training set. As the training set also requires a non-member population, we utilise the HSC+UKIDSS counterparts of the Gaia-based non-member population identified by Das et al. (2023). However, caution is necessary in selecting the non-member population.

Das et al. (2023) illustrated the distribution of the nonmember population on the colour-magnitude diagram (CMD). The distribution on the CMD (see their Fig. 15), indicates that stars with P_RF ⩽ 10 fall primarily to the left of the 10 Myr isochrone. In this analysis, we consider these stars as the more reliable non-member population for our training set. Consequently, our training set comprises 16 430 stars, including 1080 candidate members and 15 350 non-member stars. Fig. 6 displays the distribution of the training set stars on the r₂ vs r₂–i₂ and r₂ vs r₂–Y CMDs, illustrating the distinction between the two populations. As is visible from the plot, the magnitude of training set sources extends up to ∼20 mag in the r₂ band. This is primarily due to limitations in the photometry of the Gaia catalogue.

Many candidate member populations have been reported in the literature, identified by various methods, including Gaia data (see Section 2.2). However, for the training set, we use only the Gaia DR3-based member and non-member sources from Das et al. (2023). The primary reason for this choice is that, in membership analysis using an ML algorithm, a proper training set requires both member and non-member populations. Although member populations are available from other Gaia-based studies, there is no corresponding non-member population. Additionally, there is a significant overlap between all the Gaia-based candidate member catalogues, indicating that the additional Gaia sources are not crucial.

The non-Gaia-based studies report a large number of candidate member populations. Sources identified by different methods have varying properties. The literature-based sources include objects such as variable stars and evolved stars. However, in this analysis, we do not use them for our training set. This is because there is no adequate non-member population available, and several sources show a wide spread on the CMDs, which could indicate evolved older discless stars. In such cases, preparing a proper corresponding non-member catalogue is highly desirable. Due to these potential issues, we use the results of Das et al. (2023) to prepare the training set for this analysis.

5.2 RF classifier training procedure

The next phase of this analysis involves training the machine using the training set and then applying it to the complete HSC+UKIDSS catalogue to identify the candidate members of IC 1396. In this ML technique, we require parameters that possess distinctive features to successfully extract members from the extensive pool of field stars. For instance, in the Gaia-based membership analysis conducted by Das et al. (2023), certain colour parameters and proper motions played a significant role in membership identification. In our scenario, we need to identify parameters that can facilitate the more efficient retrieval of the HSC+UKIDSS-based member population.

In our case, we have spatial parameters (RA and Dec) and photometric magnitudes in the HSC and UKIDSS bands (r₂, i₂, and Y; J, H, and K). In the CMD plot (Fig. 6), we illustrate the distribution of the training set stars. However, due to the sensitivity limitations of Gaia, the training set is restricted to magnitudes up to ∼20 mag in the r₂ band. Consequently, in this analysis, we cannot heavily rely on photometric magnitudes alone since the final result might be biased towards the magnitude range of the training set.

The next set of parameters we consider are the colour terms. However, using only colour parameters for extracting the member population is insufficient. Therefore, in addition to the colour terms, we need to incorporate a few other parameters that can assist us in effectively identifying member stars. Hence, we aim to derive several additional parameters for all stars based on their photometric magnitudes and colours.

Below, we provide details of the new parameters. We opted not to utilise the spatial co-ordinate parameters (RA and Dec) as they may introduce biases towards certain parts of the complex.

5.2.1 New parameters used in RF classifier

After several attempts, we have derived the following parameters. The utilisation of these parameters yields a satisfactory outcome in our analysis. We briefly describe these parameters below.

First set of new parameters:

In this work, we derived a set of parameters using the relation AAB = A/(A-B), where A and B represent magnitudes in bands A and B, respectively. Therefore, in this relation, we derived the parameter by dividing a magnitude by the colour term of a band. For instance, RRY can be defined as r₂/(r₂–Y), and similarly for other terms. By employing this approach, we obtained twelve parameters from the HSC+UKIDSS dataset.

Second set of new parameters:

Similarly, we derived another set of parameters using a relation defined as AAB1 = A – (A – B)², where A and B represent magnitudes in bands A and B, respectively. In this case, the parameter is defined as the subtraction of the square of a colour from the magnitude of a band. For example, RRY1 = r₂–(r₂–Y)². By employing this approach, we obtained another twelve parameters through various combinations of magnitudes and colours from the HSC+UKIDSS catalogue.

Due to the combination of both colour and photometric magnitudes, these new parameters provide distinct information for each star. The values of these parameters vary significantly for a pair of stars with very different magnitudes and colours, and this unique property further aids in their effective segregation. The utility of these parameters in distinguishing between member and non-member populations can be observed in Fig. 7, in which we illustrate pair plots of some of these parameters. Additionally, pair plots of several more combinations of these parameters are presented in Fig. C.1.

The effectiveness of ML techniques arises from the characteristics embedded within parameters of the training set, and the inclusion of the derived parameters proved instrumental in identifying member populations within the complex. The resulting outcomes are statistically influenced by the composition of the training dataset; thus, the selection of an appropriate training set holds paramount importance in ML endeavours. Subsequent sections of the analysis explain how these parameters significantly contribute to achieving satisfactory results.

5.2.2 Total list of parameters

In this analysis, we have several parameters. These include three HSC photometric magnitudes (r₂, i₂, and Y), three UKIDSS photometric magnitudes (J, H, and K), twelve colour parameters (r₂ – Y, r₂ – I, i₂ – Y, r₂ – J, r₂ – H, r₂ – K, i₂ – J, i₂ – H, i₂ – K, Y – J, Y – H, and Y – K), the first set of twelve new parameters (RRY, RRI, IIY, RRJ, RRH, RRK, IIJ, IIH, IIK, YYJ, YYH, and YYK), and the second set of twelve new parameters (RRY1, RRI1, IIY1, RRJ1, RRH1, RRK1, IIJ1, IIH1, IIK1, YYJ1, YYH1, and YYK1). In total, we have forty-two parameters for use in this analysis.

5.3 Applying the RF algorithm

We have already discussed that our current training set obtained from previous Gaia analysis lacks coverage towards the fainter end of HSC bands. Therefore, our initial objective is to strengthen our training set, ensuring it encompasses both member and non-member populations across the entire photometric range of HSC bands. The main rationale behind this effort is that although we have derived a set of new parameters and can obtain member sources using them, the results are still distorted towards the fainter end of HSC bands. This bias occurs due to the prevalence of field stars compared to member populations and due to the absence of training set towards the fainter end.

To address this issue, we initially implemented the RF algorithm on the HSC+UKIDSS catalogue corresponding to a small circular region with a radius of 30^′, centred at the location of the massive star(s) HD 206267. As the central cluster is situated around the massive star(s), conducting the analysis in this limited region aids in more effectively distinguishing the member population from the field star population. The selection of a 30^′ radius was based on multiple trials, in which we determined it to be the optimal choice. Within this radius, the dominance of field stars is relatively reduced, facilitating a better separation of the cluster population from the field population.

5.4 Result from the 30^′ circle

We began by training the machine using a training set corresponding to this 30^′ circular region. Initially, there were 14 685 HSC+UKIDSS sources within this circular region, of which 3089 sources were eliminated following the same process outlined in Section 5.1.1. The remaining 11 596 sources were utilised for the analysis.

Given that the colour range of training set stars is limited compared to the entire dataset, at this stage, we utilised only the twenty-four newly derived parameters. Further details regarding the efficiency test of the machine are explained in Appendix C. The primary three leading parameters during the training procedure are RRK1, RRH1, and RRJ1. Subsequent to training the machine, we applied the RF classifier to the 11 596 HSC+UKIDSS stars and identified 399 stars with P_RF ⩾ 0.8, constituting the member populations within the 30^′ circular region. The reliability of the retrieved member population is evident from the CMD plots illustrated in Fig. 8. The distribution of these 399 sources in the CMDs conforms to a pre-main sequence locus. Of the 11 596 stars, there are 11 057 stars with P_RF ⩽ 0.5. These sources are most likely non-members or field stars within the circular region.

At this stage, this is not a definitive detection of cluster members due to the limitations of the training set and the exclusion of some parameters. However, applying the RF classifier to this small central region helps us obtain a member and non-member catalogue whose magnitudes and colours align with the complete HSC+UKIDSS catalogue. Now, we merged these derived member and non-member sources within the 30^′ circular region with the Gaia-based training set. Through this process, our final training set comprised 27 423 stars, including 1250 member stars and 26 173 non-member stars. We used this training set to retrain the machine and subsequently applied it across the entire HSC+UKIDSS catalogue of 637 595 sources (see Section 5.1.1), as is described in the next section.

5.5 Final result from RF classifier

After establishing a suitable training set, we initiated the ML procedure to identify member stars throughout the entire IC 1396 complex. The magnitude and colour range of the final training dataset aligns closely with the test dataset. Therefore, we trained the machine using all parameters, including the six photometric magnitude bands, twelve colour parameters, and twenty-four new parameters. This comprehensive training ensures the robustness of the machine with forty-two input parameters.

The efficiency test of the machine with the RF classifier is elaborated in Appendix C. The first five parameters with the highest importance are IIH1, RRY1, IIK1, IIY, and RRH. This highlights the critical role of these new parameters in segregating the member and non-member populations.

The machine demonstrates an efficiency¹² of 99.6% during its training. Following the successful training of the machine, we apply the RF classifier to the entire HSC+UKIDSS catalogue and identify 2425 stars with P_RF ⩾ 0.8. Details of the co-ordinates, photometric magnitudes in HSC+UKIDSS bands, and the P_RF value of these 2425 candidate members are provided in Table 2.

As our membership analysis relies on the photometry of sources, we cannot eliminate the possibility of a small fraction of contamination in the membership catalogue. Spectroscopic and proper-motion observations of the sources would be beneficial in further confirming their membership. Additionally, due to limitations on the input parameters, and the training set, the possibility of a small fraction of undetected of candidate members, particularly towards the fainter end, cannot be ruled out.

Furthermore, the HSC observations do not cover the entire 1.5^∘ circular region, as illustrated in Fig. 1. Consequently, the member population outside of the HSC observed area will have sensitivity corresponding to the Pan-STARRS catalogue.

Fig. 9 displays the spatial distribution of 2425 stars on the WISE 22 μm image. This distribution is consistent with earlier studies (Nakano et al. 2012; Das et al. 2023; Pelayo-Baldárrago et al. 2023). Most stars form a diagonal pattern from north to south, connecting the BRC IC 1396 N with the central cluster. The western part of the star-forming complex exhibits slightly higher stellar density. Additionally, a small fraction of stars appear randomly distributed throughout the complex with minor clumping. This spatial distribution suggests the presence of several sub-clusters associated with the complex. Previous studies have also identified sub-clusters in IC 1396 using member stars detected from various surveys (Nakano et al. 2012; Pelayo-Baldárrago et al. 2023; Das et al. 2023). In the next section, we explore into the analysis of derived properties of the member stars and identify the sub-clusters within the IC 1396 complex.

5.6 Impact of extinction in the analysis

Extinction has not yet been considered as an input parameter in this analysis, although variable extinction could potentially be present along the line of sight towards the star-forming complex. Therefore, it is important to test the impact of extinction in this membership analysis by incorporating visual extinction (A_V) alongside the other parameters. For this purpose, we obtained the A_V values (both mean and median) towards the location of each star in our HSC+UKIDSS catalogue from the Bayestar 19¹³ 3-D dust reddening map (Green et al. 2019). This dust reddening map, derived with the inclusion of Gaia parallaxes and stellar photometry from Pan-STARRS1 and 2MASS, created a 3-D dust map of the sky, out to a distance of a few kiloparsecs.

In Fig. 10, we show the histogram distribution of both mean and median values of A_V obtained for all stars. The mean and median values of A_V lie in the range of 0.96–3.69 and 0.86–4.53, respectively, consistent with the extinction values reported in the literature for the star-forming complex (Contreras et al. 2002; Sicilia-Aguilar et al. 2005; Nakano et al. 2012; Pelayo-Baldárrago et al. 2023).

To assess the impact of extinction on the analysis, we incorporated a mean A_V for each star and repeated the entire analysis. This resulted in 2310 sources with P_RF ⩾ 0.8, approximately 92% of which match the existing catalogue of 2425 member stars. The final results differ by only 5% compared to our previous findings. Fig. E.1 illustrates the distribution of the 2310 stars alongside the 2425 stars on CMDs, further demonstrating the minimal impact of including A_V in the analysis. We also tested the results using the median A_V values. However, this adjustment also did not lead to any significant changes, with the final result comprising 2241 stars, approximately 93% of which match the results obtained without considering A_V.

In this exercise, we observe that A_V does not have major impact on the training of the machine. The primary reason for this could be the use of the newly derived parameter sets. As several combinations of magnitude and colours are involved in their derivation, these parameters help in reducing the temperature-extinction degeneracy. Consequently, these parameters hold significantly greater importance in the analysis compared to extinction. Additionally, it is important to note that the Bayestar19 dust map has a resolution of 3.4^′–13.7^′ (Green et al. 2015, 2018), which is substantially larger than the resolutions of the datasets used in this analysis. As a result, stars in our catalogue that fall within the same pixel will share the same A_V values, potentially introducing bias for some sources. However, this analysis demonstrates that the inclusion of A_V does not significantly impact the final outcome. For further analysis in the paper, we used the 2425 candidate member stars initially obtained without considering A_V.

5.7 CMDs of member stars

In Fig. 11, we present the r₂ vs r₂–i₂ and r₂ vs r₂–Y CMD plots. These plots show the distribution of the 2425 candidate members of IC 1396 identified in this study. The plots demonstrate that the detected 2425 stars closely follow to the anticipated pre-main-sequence locus, with the majority of candidate members situated near 2–4 Myr isochrones. Comparison with the CMDs presented in Fig. 6 suggests that the distribution of this member population matches with the Gaia-DR3-based stars. This further confirms the authenticity of the member populations identified in this study. In the subsequent sections, we compared our identified membership catalogue with the literature to further assess its reliability and the effectiveness of the method adopted in this work.

5.8 Comparison with literature

Section 2.2 provides a brief overview of the literature-based membership analyses conducted for the star-forming complex IC 1396. There are 2178 stars (see Section 2.2) reported from various studies, excluding Gaia data. Analyses utilising Gaia data have reported a total of 1468 member stars in this complex. We compare these two literature catalogues with the member stars identified in our study.

5.8.1 Comparison with non-Gaia based literature sources

First, we compare the member stars identified in this study with the non-Gaia-based catalogue of 2178 stars. Out of which 1224 stars have HSC+UKIDSS counter parts; however, of these, 514 stars are not considered for comparison because their Gaia parallaxes do not fall within the valid range of the star-forming complex, as explained in Section 5.1.1. Eventually, out of these 2178 stars, 710 have corresponding entries in the final HSC+UKIDSS catalogue, upon which we have applied the RF method. We find that, out of these 710 stars, we retrieve 403 (∼57%). Upon closer inspection (see Fig. 11), it becomes evident that the majority of the stars retrieved in this work belong to the brighter end of the HSC catalogue (around 20 mag in the r₂ band). Due to the higher sensitivity of Subaru-HSC observations, many fainter sources were detected in this work compared to previous studies. A few stars from the 403 stars are found towards the fainter end, with magnitudes reaching up to 26 mag in the r₂ band. Additionally, we can clearly see that the locus of the literature-based member population overlaps with the locus of the candidate member stars retrieved in this work, enhancing the reliability of this membership analysis.

In Table D.1, we provide a detailed comparison of our identified sources with member populations identified through various methods in the literature. The table shows that, except for the literature-based variability sources (Meng et al. 2019) and the YSOs detected through NIR, MIR and X-ray (Silverberg et al. 2021), sources detected through all other methods show a high retrieval percentage. Although the YSOs from Silverberg et al. (2021) also include X-ray detections, they represent a small fraction of the NIR population. Of their 139 sources (Table D.1) used in our RF analysis, only 18 sources (∼13%) have X-ray counterparts. Therefore, most of the YSOs not retrieved in our work from Silverberg et al. (2021) and Meng et al. (2019) are of NIR origin. A large fraction of sources from spectroscopy-based studies (Contreras et al. 2002; Sicilia-Aguilar et al. 2006b, 2013) and other surveys (Spitzer, H_α, and X-ray) are retrieved in this work. This higher retrieval suggest the effectiveness of the method adopted in this work.

Out of the 710 stars, 353 lie within the 90% completeness limit, which corresponds to 18.5 mag in the r₂ band. Of these, we successfully retrieve 273 stars (∼77%), indicating the robustness of the method adopted in this work. Among the remaining 80 stars within the 90% completeness range, many are located to the left of the 10 Myr isochrone, while those that fall to the right of the 10 Myr isochrone still have P_RF ⩾ 0.5 (see Fig. E.2). Beyond the 90% completeness magnitude, only a few stars from the literature are retrieved by our method, and these lie to the right of the 10 Myr isochrone.

In the r₂ band, 22 mag corresponds to a mass of approximately 0.1 M_⊙ (see Section 6.1), roughly marking the boundary between stars and brown dwarfs. Beyond 22 mag in the r₂ band, we identify 253 stars, indicating a significant improvement in the detection of very low-mass objects within the star-forming complex. However, from non-Gaia based literature sources, only 12 such very low-mass stars are retrieved by us.

The remaining 307 stars that were not recovered have a membership probability of less than 80%. The random distribution of these 307 stars on CMDs (refer to Fig. E.2) indicates their poor membership probability based on this membership analysis. Most of them are observed to fall on the left side of the 10 Myr isochrone on the CMDs. These 307 stars do not satisfy the properties of membership. This is because, in the current method, the machine is trained with an input training set that focuses on the population of the complex, with a potential age of less than 10 Myr. The majority of these 307 stars belong to the older population of the star-forming complex. This can be considered as a potential limitation of the method, as it is not able to identify the older populations of the complex, which we discuss briefly in Section 5.9. Most of the 307 stars have been detected through observations targeting variable stars (Meng et al. 2019) and YSOs from Silverberg et al. (2021). These 307 stars only include a small fraction of stars identified from spectroscopic studies and observations involving H_α and X-ray. Thus, most of the 307 stars are of NIR origin, and hence spectroscopic studies will be helpful in revealing their true nature.

5.8.2 Comparison with Gaia-based literature sources

Next, we compare our membership catalogue with all the Gaia-based member stars. There are 1468 Gaia-based member stars detected in IC 1396. Out of these, we have 1217 stars with HSC+UKIDSS counterparts over which the RF method has been applied. Our comparison reveals that out of these 1217 stars, 1020 (∼84%) stars are retrieved in this work. This high retrieval of the Gaia-based member stars further confirms the reliability of the candidate member stars obtained in this work. This is because the Gaia-based member population offers a higher reliability of candidate membership due to their astrometric parameters, such as parallax and proper motions. The mean, median, and standard deviation in parallax of the Gaia-retrieved sources are [1.089 ± 0.004, 1.078, 0.116] mas, respectively. Similar statistics are [−1.179 ± 0.002, −1.193, 0.284] mas yr⁻¹ and [−4.277 ± 0.004, −4.473, 0.684] mas yr⁻¹ for μ_α cos δ and μ_δ, respectively. These values agree with the results obtained by Das et al. (2023) from their Gaia-DR3 member stars. The statistics of parallax for the retrieved 1020 Gaia members also match other works (Sicilia-Aguilar et al. 2019; Pelayo-Baldárrago et al. 2023; Das et al. 2023).

The remaining 197 Gaia stars that were not recovered by us are primarily (∼52%) from the survey of Cantat-Gaudin et al. (2018), in which they included all stars with a membership probability above 50%. This current study could not recover their low-probability stars, which possibly have a higher chance of contamination. Of the 1080 Gaia-based sources used in the training set, 994 (∼92%) stars are retrieved with P_RF ⩾ 80%. For the remaining 86 stars, a large fraction of stars still have P_RF ⩾ 50%.

In classification tasks, the RF classifier computes the probability of a data point belonging to each class by averaging the probabilities assigned by individual decision trees. In the context of the membership analysis, the RF classifier assigns a probability to each star being a member of the IC 1396 complex based on its photometric properties and other derived parameters. The probability represents the confidence of the model in the classification of the star as a member or non-member, with higher probabilities indicating greater confidence in the classification. However, due to the statistical nature involved in assigning the probabilities, a small degree of mismatch between the final outcome and the training set is still expected, as observed in this study.

For easier understanding, we present the entire comparison of catalogues in the form of a Venn diagram in Fig. 12. Out of the 2425 candidate member stars, 1331 are newly identified in this study. In Table 2, the new member stars are marked with an asterisk symbol. Additionally, in one of the most recent analyses of the complex, Gupta et al. (2024) conducted a search for low-mass objects within a 22^′ region centred around the central massive star HD 206267. They employed ML techniques on the HSC datasets and identified a population of 458 stars within that specific region. In our work, we are able to retrieve approximately 86% of their sources, and the small discrepancy is primarily due to the different approaches employed by both analyses.

5.9 Caveats and points to ponder

In this section, we discuss the limitations or caveats associated with our membership analysis. As previously mentioned, conducting a membership analysis of a star-forming complex based solely on photometric information presents challenges. The primary difficulty arises due to the unavailability of a training set extending to the fainter end of the HSC+UKIDSS catalogue. To address this issue, we have derived several new parameters using the magnitudes and colours of the sources. These new parameters, when combined with magnitudes and colours, enable us to identify a reliable member population from a large pool of sources.

It is likely that some candidate members may remain undetected, especially at the fainter end. This is primarily due to the lack of sources in the training set at the fainter end and in addition to the consideration of a high probability of P_RF ⩾ 0.8. This lack of sources in the training set towards the fainter end thereby introduces a bias. Hence, statistically, the chance of non-detection of candidate members is higher at the fainter end compared to the brighter end. This limitation affects the detection of very low-mass objects such as brown dwarfs in the star-forming complex.

Furthermore, this analysis only identifies the young sources of the complex lying rightward of the 10 Myr isochrone. The older populations of the complex of >10 Myr (if any), or the YSOs that are scattered on the CMD and lying leftward of the 10 Myr isochrone would have potentially been removed in our analysis (See Appendix E). There are several reasons for the broad distribution of young sources on the CMD, such as differential reddening, variability, age spread, etc (e.g. Jose et al. 2017; Panwar et al. 2017). Such members are difficult to detect in this analysis due to a lack of proper knowledge about these sources and the absence of a reliable training set.

The non-detection of the older population can also be attributed to the training set used in this work (see Section 5.1.2), which includes Gaia-based stars with ages of less than 10 Myr, thus introducing a potential bias. The outcome of ML algorithms is highly dependent on the nature of the training set. Therefore, the detection of member populations younger than 10 Myr and possibly undetected stars towards the fainter end can be attributed to the bias present in the training set. However, we assume that the bias is less significant, as only a relatively small percentage of spectroscopic (13%) and X-ray (23%) stars lie to the left of the 10 Myr isochrone. We are keen on improving such membership analysis by addressing the known limitations in future studies.

It is worth noting that integrating UKIDSS with the HSC catalogue enhances the reliability of our membership detection. Therefore, we advocate for incorporating more bands with photometry to yield more dependable results from such analyses. The rationale for using additional wavelength bands is similar to spectral energy distribution (SED) analysis, where more data points contribute to improved accuracy. In a similar study, Rebull et al. (2023) identified YSO populations in IC 417 by analysing the SED obtained from datasets across several optical and infrared bands. Although our approach differs, we have retrieved the reliable member populations of IC 1396.

Additionally, it is true that even though the final membership detection is more reliable, it is less complete because the final HSC+UKIDSS catalogue is restricted to the common sources from both catalogues. Hence, if a source is missing in either one of the catalogues, for example, the UKIDSS catalogue, it will also be missing from the final HSC+UKIDSS catalogue and eventually from the final candidate member detection. This can also be considered as a limitation in such membership analysis.

It is important to note that the stars in the training set are available up to ∼20 mag in the r₂ band. Additionally, other non-Gaia-based studies in the literature primarily correspond to the brighter end (see Fig. 11). Of the 2425 candidate stars identified in this work, 1461 are brighter than 20 mag in the r₂ band. Moreover, there is a high retrieval rate with Gaia-based literature sources, as well as sources detected through methods based on optical spectroscopy, Spitzer, Hα/X-ray emissions, optical photometry, and 2MASS photometry (see Section 5.8). Therefore, these 1461 stars have a higher confidence level for being members of the star-forming complex. The distribution of the 2425 candidate populations in the CMD (see Fig. 13) suggests the presence of very low-mass stars (e.g. brown dwarfs).

It is also important to consider that the identified candidate population, particularly at the low-mass end, may include contamination from different origins. For example, very faint foreground stars and distant background objects may be contaminating the retrieved sample. The nebular emission, particularly the bright rims near the BRCs, is associated with line emissions such as Hα and [S II] (Beltrán et al. 2009; Pelayo-Baldárrago et al. 2023). These line emissions may blend with the stars and alter their apparent properties. The nebular emissions may also appear as compact, star-like sources. Additionally, shock emissions near the BRCs (Nisini et al. 2001; Beltrán et al. 2002; Caratti o Garatti et al. 2006; Pelayo-Baldárrago et al. 2023) might further affect the observed characteristics of the stars. Such contaminants are difficult to filter out using a photometry-based membership analysis.

Therefore, additional studies, such as spectroscopic and proper motion analyses, will be essential to further confirm the membership of the identified candidate population in this work, including brown dwarfs. A detailed analysis of these brown dwarfs will provide crucial insights into the formation of low-mass objects in clustered environments within the star-forming complex.

As was explained in previous sections, we initially applied the RF technique to a 30^′ circular area in our methodology, which facilitated the enhancement of our training set. While deriving the member stars, we utilised twenty-four and forty-two parameters in the 30^′ circle and the entire complex, respectively. However, we experimented with reducing the parameters by leveraging the correlation matrix and retaining only the highly uncorrelated parameters that would provide maximal information. It is worth noting that, through this approach, our results varied by approximately ten per cent. We opted to use more parameters to ensure the retrieval of the most reliable outcome.

The applicability of this method employed in this work may vary from one star-forming region to another and for different datasets. Nevertheless, this method serves as a helpful or guiding tool in such membership analyses. Our intention is to apply this methodology to other star-forming regions using diverse datasets.

6 Results

6.1 Age and mass completeness

Previous studies have provided an age estimate of ∼2–4 Myr for IC 1396 (Sicilia-Aguilar et al. 2005; Pelayo-Baldárrago et al. 2023; Das et al. 2023). From the identified candidate members, we also attempt to derive an age estimate of the complex. Additionally, to gain insight into the mass limit of the candidate population, we use all the candidate members identified in Section 5.5, and the r₂ vs r₂–Y CMD is shown in Fig. 13. After correcting for extinction and distance using reddening laws of Wang & Chen (2019), we plot isochrones along with the mass evolutionary tracks from Baraffe et al. (2015). Earlier works have reported sparse visual extinction towards IC 1396. Using NIR and optical data, Sicilia-Aguilar et al. (2005) estimated an average visual extinction of the entire complex to be A_V = 1.5 ± 0.5 mag, which also matches estimations from a few other studies (Contreras et al. 2002; Nakano et al. 2012). A more recent analysis of the complex derived an average extinction of A_V = 1.4 ± 0.52 mag from Gaia-EDR3 colours of the member stars (Pelayo-Baldárrago et al. 2023).

In the star-forming complex, the feedback effect of central massive stars clears the dust and gas mass and creates a cavity, which can be seen in the WISE 22 μm image (see Fig. 1). Thus, the inner part of the complex is expected to be associated with lesser visual extinction compared to the surrounding regions. Due to these reasons, Das et al. (2023) have used the minimum extinction value of A_V = 1 mag obtained from Nakano et al. (2012). In this analysis, we have also used the same extinction value and the distance of 917 pc derived from Gaia-DR3-based member stars to correct the isochrones.

The Fig. 13 presents the density distribution of the 2425 candidate members on the CMD. The plot shows that the maximum number of sources lie within an age range of ∼2–4 Myr, and the mean age would be ∼3 Myr. This matches with all the previous studies (Sicilia-Aguilar et al. 2005; Pelayo-Baldárrago et al. 2023; Das et al. 2023). However, a detailed SED analysis will help deriving a more accurate age estimation for individual stars and the whole complex.

In Fig. 13, we have overplotted evolutionary tracks of various masses starting from 0.02 M_⊙. In this section, we derive a quantitative estimate of the limiting mass and the 90% mass completeness limit of the member populations. We use the HSC-r₂ band to estimate them since its depth is maximum compared to the other two bands. HSC-r₂ band magnitude limit of the membership catalogue is ∼27 mag. Using an isochrone of age 2 Myr, this magnitude limit in the r₂ band corresponds to a mass ∼0.02 M_⊙. Detection of low-mass objects of IC 1396 becomes possible due to the higher sensitivity of the Subaru-HSC observations. Next, we derive the 90% mass completeness of the member population. In Section 3.3, we have derived the 90% completeness of the HSC catalogue to be ∼24 mag in r₂ band. Considering the isochrone of 2 Myr, the 24 mag in r₂ corresponds to a mass of ∼0.04–0.05 M_⊙. This analysis shows possible detection of very low-mass objects (e.g. brown dwarfs) in IC 1396. It is to be noted that this membership analysis includes the UKIDSS catalogue with the HSC catalogue, so the final HSC+UKIDSS has 90% completeness of ∼22 mag in r₂ band, and so, the 90% mass completeness is ∼0.1 M_⊙.

We note that the isochrones plotted in Fig. 13 do not correctly fit the low-mass stars and deviate for mass below 0.1 M_⊙. The origin of this discrepancy is still unclear and probably an inherent limitation of isochrone models for low-mass stars (Khalaj & Baumgardt 2013; Li et al. 2020). This discrepancy between observations and multiple isochrone models has also been reported in several observational analysis such as in Taurus, Ophiuchus, Chamaeleon, Serpens-south, NGC 2244 etc. (e.g. Mužić et al. 2013; Ribas et al. 2017; Zhang et al. 2018; Jose et al. 2020; Almendros-Abad et al. 2023). Such caveats in the isochrone models can affect determining the individual properties of stars. However, the overall cluster properties, such as age and mass, will not be affected because those properties are mainly constrained by the location of the majority of sources in the CMD.

6.2 Sub-clustering in IC 1396

Previous studies have identified several sub-clusters within the IC 1396 complex (Nakano et al. 2012; Pelayo-Baldárrago et al. 2023; Das et al. 2023). In this section, we aim to identify the subclusters associated with IC 1396, using the more reliable set of 1461 member stars (see Section 5.9), whose spatial distribution (see Fig. 14) reveals the presence of multiple sub-clusters within the star-forming complex. To accomplish this, we first generate a surface density map using the member stars and applying the nearest neighbour (NN) method (Casertano & Hut 1985; Schmeja 2011). This method enables us to derive the j-th nearest neighbour density of a star using the following expression: $${\rho }_{\mathrm{j}}=\frac{\mathrm{j}-1}{\mathrm{S}({\mathrm{r}}_{\mathrm{j}})},$$(1)

where r_j is the distance to the j-th nearest neighbour, and S(r_j) is the surface area with radius r_j. Many studies have indicated that j = 20 is an optimal value for cluster identification (Schmeja et al. 2008; Ramachandran et al. 2017; Damian et al. 2021). In our analysis, we also use j = 20 and obtain the nearest neighbour density distribution of member stars. From the derived density values of all member stars, we generate a stellar density map with a pixel size of 0.1 pc (20.5′′). In Fig. 14, we present the distribution of stellar density in the form of contours overlaid on the WISE 22 μm map. The lowest contour is set at 1.5 stars pc⁻², representing the threshold above which the maximum number of sources is observed. These stellar density contours also reveal the presence of sub-clusters within the IC 1396 complex.

Using the astrodendro algorithm (Robitaille et al. 2019) of Python, we identify sub-clusters within IC 1396. This algorithm constructs tree structures starting from the brightest pixels of the dataset and progressively includes fainter ones. The input parameters for the algorithm are the threshold flux value (minimum value), a contour separation value (minimum delta), and the minimum number of pixels. In this analysis, we set the threshold and minimum delta to 1.5 and 0.5 stars pc⁻², respectively. For potential cluster detection, we specify the minimum number of pixels as 200. These parameter settings, determined through multiple trials, enable optimal cluster detection. Applying the astrodendro algorithm with these inputs, we identify 9 leaf structures in the dataset, which correspond to the sub-clusters of the IC 1396 complex.

In Fig. 14, we illustrate the detected sub-clusters along with their respective names. All identified sub-clusters are found to be associated with the BRCs and PDRs of IC 1396. Notably, we detect cluster C-1 around the central massive star(s) HD 206267, which exhibits two distinct sub-structures labelled as C–1A and C-1B. Several of these sub-clusters, including C-1A, C-1B, C-2, C-4, and C-6 have been reported in previous studies (Nakano et al. 2012; Das et al. 2023). For instance, sub-cluster C-4 is associated with BRC IC 1396 N, while C-6 is linked with BRCs SFO 39 and SFO 41. Moreover, Das et al. (2023) identified a sub-cluster near BRC SFO 35 based on Gaia-DR3 member stars, which again detected by us.

In our study, we detect a few new sub-clusters within the IC 1396 complex. We identify sub-cluster C-3 along the diagonal line connecting the central cluster and IC 1396 N. Furthermore, a new sub-cluster, C-5 is detected towards the eastern edge of IC 1396, linked with BRC IC 1396 G (Pottasch 1956).

Table 3 provides statistics including the radius, number of stars, and average stellar density of the sub-clusters. To calculate the physical radius (R_cluster = (A_cluster/π)^0.5; Das et al. 2017), we utilised the apertures obtained from astrodendro. Area of each sub-cluster (A_cluster) is determined as the product of the number of pixels (N) and the area of each pixel (A_pixel). The central cluster C–1 is the largest, with a radius of 3.06 pc, hosting the maximum number of stars (312), and with an average stellar density of 10.62 stars pc⁻². Sub-clusters C-1A and C-1B, nestled within the central cluster, exhibit the highest stellar densities among all sub-clusters, with values of 15.32 and 13.95 stars pc⁻², respectively. Considering the stellar density and distribution of sub-clusters, we explore the star formation history of the complex in subsequent sections.

7 Star formation history of IC 1396

IC 1396, nearby, is one of the best-studied star-forming complexes, where the feedback-driven star formation activity dominates the entire complex. The strong stellar wind from the central massive star(s) creates a giant cavity of radius ∼1.5^∘ (∼25 pc) by clearing most of the gas and dust from the inner region. At the same time, the energetic UV radiations shape the surrounding region by creating several structures such as BRCs, filaments, and fingertip structures within and around the complex (Schwartz et al. 1991; Froebrich et al. 2005; Saurin et al. 2012). Star formation activity of the individual BRCs and also the entire complex have been carried out in many past studies (Ikeda et al. 2008; Sicilia-Aguilar et al. 2014; Panwar et al. 2014; Pelayo-Baldárrago et al. 2023; Das et al. 2023). In this section, we attempt to briefly discuss the over-all star formation history of the entire complex.

7.1 IC 1396: An expanding bubble

In Fig. 15, we show the distribution of 1020 candidate members that have Gaia counterparts. The direction of the proper motions relative to the central sub-cluster C–1 is illustrated, indicating a random motion of the stars. Previous studies have shown that the complex IC 1396 is gravitationally unbound and is expanding (Patel et al. 1995; Pelayo-Baldárrago et al. 2023). In this section, we explore the kinematics of IC 1396 using the astrometric details of these 1020 stars.

The variation in RA with respect to μ_α cos δ(μ_{α *}) and Dec with respect to μ_δ shows a positive correlation (Fig. 16). This correlation suggests the presence of velocity gradient in IC 1396 and indicates the expanding nature of the complex, as observed in other regions (e.g. Wright et al. 2019). To further investigate the expanding feature of IC 1396, we conduct a simple 2-D analysis, as a 3-D analysis is not feasible due to the lack of radial velocities for all 1020 stars.

For this analysis, we first translate the celestial co-ordinates to the 2-D Cartesian co-ordinate system. Both positions and velocities in the Cartesian co-ordinates can be derived using the following equations (Gaia Collaboration 2018). $$\{\begin{array}{l}x=\cos \delta \sin (\alpha -{\alpha }_C)\\ y=\sin \delta \cos {\delta }_C-\cos \delta \sin {\delta }_C\cos (\alpha -{\alpha }_C)\\ {\mu }_x={\mu }_{\alpha \ast }\cos (\alpha -{\alpha }_C)-{\mu }_{\delta }\sin \delta \sin (\alpha -{\alpha }_C)\\ {\mu }_y={\mu }_{\alpha \ast }\sin {\delta }_C\sin (\alpha -{\alpha }_C)\\ +{\mu }_{\delta }(\cos \delta \cos {\delta }_C+\sin \delta \sin {\delta }_C\cos (\alpha -{\alpha }_C))\end{array}$$(2)

where, α_C and δ_C are the mean RA and Dec (324.668; 57.486) of stars within sub-cluster C–1, used as the reference position. The linear velocities (μ_x, μ_y) are derived as μ_x = (–k/ω) μ_{α *} and μ_y = (k/ω) μ_δ, where ω is the parallax of each star and k = 4.74 is the conversion factor from mas yr⁻¹ to km s⁻¹ at a distance of 1 kpc. We derived the positions and velocities for all stars in the X and Y directions using these equations.

After deriving the individual velocities in the X and Y directions, the relative velocities v(v_x, v_y) of each star with respect to the mean values of C–1 are obtained. The directions of the relative velocities indicate whether a star is moving towards or away from C-1. The final expansion velocity can be obtained as v_exp = v · r̂, where r̂ is the radial unit vector. If a star moves outward from C–1, then v_exp will be positive; otherwise, it will be negative.

In Fig. 17, we show the distribution of the 1020 stars in the Cartesian co-ordinate system along with their relative velocities. As seen from the plot, the majority of stars are moving away from C-1, thereby providing evidence of the expansion of the star-forming complex. To obtain a quantitative estimate of the expansion velocity, we calculate the mean velocity within radial distance bins of 1 pc from the reference position. The variation in v_exp with radial distance is shown in Fig. 18. There is very little expansion up to a distance of ∼5 pc, which encompasses the central sub-cluster C-1. Previous work by Das et al. (2023) shows very little expansion of ∼1 km s⁻¹ for sub-cluster located within the central region, based on their 3-D analysis of 68 cluster members with radial velocity information. However, Fig. 18 does not reveal expansion within the central cluster. This is because, in this analysis, we consider the relative velocities with respect to the mean values of C–1 itself. Thus, it is expected that stars within C-1 will exhibit both positive and negative v_exp, not reflecting the true expanding structure within C–1 itself. In this analysis, a clearer picture of expansion is obtained by moving away from cluster C–1. Beyond 5 pc, the expansion shows a gradual rise up to 12 pc, reaching an expansion velocity of 4.2 ± 0.5 km s⁻¹. A dip is observed up to 13 pc, followed by another rise in v_exp.. The expansion velocity increases to its maximum value of 5.7 ± 0.5 km s⁻¹ at a distance of 21 pc, which is close to the periphery of the star-forming complex.

This variation in v_exp with distance reveals the outward movement of distant stars from the central cluster. This provides evidence of a gravitationally unbound system, demonstrating a clear expanding feature. The gravitationally unbound and expanding nature of the complex has been reported in earlier works. Pelayo-Baldárrago et al. (2023) identified sub-clusters of the complex based on astrometry and observed the expanding nature by calculating pairwise projected velocity differences among the stars. From CO maps, Patel et al. (1995) suggested an expansion velocity of ∼5 km s⁻¹, which is also supported by Getman et al. (2012) and Sicilia-Aguilar et al. (2019). Although we used a different approach to analyse the kinematics, our findings agree with the conclusions reported in previous studies of the star-forming complex.

7.2 Hierarchical, sequential, and distributed mode of star formation

The spatial distribution of the member sources (Fig. 14) and the association of sub-clusters with the BRCs and PDRs suggest ongoing feedback-driven star formation activity within the complex. Table 3 provides the mean stellar densities of all subclusters. Higher stellar densities are observed around the central massive star(s) HD 206267, reaching a maximum within subcluster C-1A and then peaking within sub-cluster C-1B, near the head of the BRC IC 1396 A. The stellar density gradually decreases and shows localised peaks in surrounding sub-clusters. This distribution indicates that the majority of star formation is occurring within the central part around the massive star(s), with density decreasing gradually towards the surrounding area. The localised peaks in stellar densities can be attributed to the local molecular conditions, where the star formation efficiency may vary. The presence of the massive star(s) HD 206267 drastically influencing the star formation within the central cluster, resulting in its higher stellar density. However, all the sub-clusters of IC 1396 have nearly similar average ages. This suggests a scenario where the molecular cloud undergoes a hierarchical fragmentation process, leading to star formation occurring within a similar or nearly simultaneous timeframe (Bonnell et al. 2003; Grudić et al. 2018; Torniamenti et al. 2022). Such possibilities of the star formation process in IC 1396 have also been reported earlier (Das et al. 2023). However, it is also possible that the fragmentation process may continue, forming much younger stars and leading to sequential hierarchical star formation within the molecular cloud.

The distribution of the identified stellar population in this work shows the presence of many sources distributed randomly in the whole star-forming complex of IC 1396 (see Fig. 14). These stars are lying in regions of stellar density less than 1.5 stars pc⁻². This low density and the scattered distribution of stars are characteristic of a distributed mode of star formation, where stars form in a more isolated and spread-out manner rather than in dense, clustered regions (Lada & Lada 2003). This observation is similar to other regions known for distributed star formation, such as the Taurus Molecular Cloud (Luhman et al. 2006; Luhman 2018), and supports the notion that star formation can occur efficiently even in environments with relatively low stellar densities (Evans et al. 2009). Apart from distributed star formation, these low-mass objects might have been ejected from their parent sub-clusters through dynamic processes (Bate et al. 2002, 2003; Goodwin & Whitworth 2007).

In IC 1396, stars are expected to move at a velocity of ∼3 km s⁻¹ when they are ejected from their parent sub-clusters. This velocity is an average of values reported in previous studies of this complex. Das et al. (2023) estimated the expansion velocity of the central cluster to be 1.11 km s⁻¹. From CO maps, Patel et al. (1995) obtained the expansion velocity of the entire complex to be around 5 km s⁻¹, a result supported by later studies (Pelayo-Baldárrago et al. 2023; Getman et al. 2012; Sicilia-Aguilar et al. 2019). Our current kinematic analysis also reveals an expanding feature of IC 1396. An earlier study by Sterzik & Durisen (1995) estimated a velocity of ∼3–4 km s⁻¹ for ejected stars from clusters. If we assume an ejected velocity of ∼3 km s⁻¹ perpendicular to the line of sight, then over the average age of 4 Myr, a star might travel a distance of ∼12 pc. The radius of IC 1396 is 1.5^∘, corresponding to ∼25 pc. Thus, randomly distributed stars are expected to lie within IC 1396, even if they are ejected from the central cluster. Although we speculate about the ejection of stars, we are unable to determine the ejection rate, which requires a clear understanding of the total cluster mass and its distribution, necessitating more detailed analysis.

In this work, we provide a brief discussion of the possible modes of star formation processes in IC 1396. However, a much more dedicated analysis is essential to gain a deeper understanding of the ongoing star formation activities in the star-forming complex.

8 Summary

This work presents the deepest photometric observations of the star-forming complex IC 1396 with the Subaru-HSC. By combining the optical and NIR data from Pan-STARRS and UKIDSS catalogues with the HSC catalogue, we conduct a comprehensive membership analysis of this star-forming complex. In the following, we outline the significant findings of this study.

1. Our Subaru-HSC observations of IC 1396 were conducted with two pointings in the r₂, i₂, and Y bands. We processed the optical data using HSCPIPE V6.7 to generate the photometric catalogue and images. After applying various quality flags, we obtained a clean HSC catalogue ensuring good photometric quality, where the magnitude error is less than 0.2 for all three bands. The complete catalogue comprises 859 867 stars with photometric magnitudes in all three bands. The limiting magnitude of the catalogue is 27.6, 26.6, and 25.2 mag in the r₂, i₂, and Y bands, respectively. The 90% completeness of the catalogue is 24, 22, and 21 mag in the r₂, i₂, and Y bands, respectively;

2. We expanded our dataset by including additional sources from the Pan-STARRS catalogue to cover the complete circular region of radius 1.5^∘ along with to compensate the brighter end of the HSC catalogue. The merged catalogue consists of 1 322 119 stars, of which 802 100 have counterparts in the UKIDSS catalogue. The membership analysis was conducted using these stars;

3. We employed the RF classifier method of the ML algorithm and identify 2425 high-probability candidate members with P_RF ⩾ 0.8. The training set was largely limited to sources younger than approximately 10 Myr; consequently, the retrieved candidate members are also within this age range;

4. The CMD distribution of these 2425 candidate member stars indicates an average age of 2–4 Myr, consistent with previous studies. We determined the limiting and 90% mass completeness limit using the r₂ band and isochrones from Baraffe et al. (2015). The limiting mass of the members is approximately ∼0.01–0.02 M_⊙, with a 90% mass completeness limit of ∼0.04–0.05 M_⊙. Notably, this study marks the detection of many low-mass member objects in IC 1396;

5. The spatial distribution of the member stars reveals the presence of several sub-clusters. Using the 1461 stars brighter than 20 mag in the r₂ band, we identify sub-clusters in the complex. By applying the NN method, we detected eight sub-clusters. In addition to the central cluster, one subcluster (C-3) is observed on the diagonal connecting the BRC IC 1396 N and the central cluster. Moreover, the central cluster C-1 exhibits the presence of two sub-clusters, C–1 A and C–1 B, situated on two sides of the central massive star(s).;

6. Our kinematic analysis, based on the astrometry of Gaia counterparts of candidate member sources, reveals the expanding nature of the star-forming complex. Additionally, we briefly discuss the possibility of hierarchical, sequential, and distributed modes of star formation activities in IC 1396.

We wish to emphasise that the primary aim of this work is to identify the member populations of the star-forming complex IC 1396 using the deep, sensitive Subaru-HSC data. Therefore this work provides an opportunity to evaluate the strengths and limitations of ML techniques when applied to such astronomical datasets.

It is true that all membership detection methods have inherent limitations, including our own. In this work, we have derived several new parameters from various combinations of magnitudes and colours, which help mitigate some of these limitations, although not entirely. Due to these constraints, we are unable to detect candidate members of IC 1396 older than 10 Myr, primarily because of the absence of older candidate member populations in the training set. Additionally, the lack of sources in the training set towards the fainter end introduces a bias, and hence, statistically, the chance of non-detection of candidate member stars is higher at the fainter end compared to the brighter end. This limitation affects the detection of very low-mass objects such as brown dwarfs in the star-forming complex. Furthermore, the presence of nebular emission lines, along with foreground and background stars, can lead to potential contamination. Spectroscopic analyses will be crucial in further filtering out such contaminants identified through photometric studies.

We are keen to improve our detection method in future works. It would be interesting to test the results by combining physical parameters such as effective temperature, luminosity, age, and mass, with the existing parameters used in this study. Since ML methods are effective in handling large parameter sets, we expect that the use of more robust parameters could potentially help to better constrain the ML model. However, we must also focus on improving the training set, ensuring that it includes both member and non-member populations without age restrictions, as failing to do so may introduce bias into the final outcome. We are eager to contribute to such membership analyses and make the method more robust by addressing known limitations, enabling it to work without age restrictions across various datasets and star-forming complexes in diverse environments.

Data availability

Full Table 2 is available at the CDS via anonymous ftp to cdsarc.cds.unistra.fr (130.79.128.5) or via https://cdsarc.cds.unistra.fr/viz-bin/cat/J/A+A/698/A297

Acknowledgements

We thank the anonymous referee for his/her detailed comments, which significantly improved the paper. This research is based on data collected at Subaru Telescope with Hyper Suprime-Cam, which is operated by the National Astronomical Observatory of Japan. We are honoured and grateful for the opportunity of observing the Universe from Mauna Kea, which has the cultural, historical and natural significance in Hawaii. We are grateful to The East Asian Observatory which is supported by The National Astronomical Observatory of Japan; Academia Sinica Institute of Astronomy and Astrophysics; the Korea Astronomy and Space Science Institute; the Operation, Maintenance and Upgrading Fund for Astronomical Telescopes and Facility Instruments, budgeted from the Ministry of Finance (MOF) of China and administrated by the Chinese Academy of Sciences (CAS), as well as the National Key R&D Program of China (no. 2017YFA0402700). The authors thank the entire HSC staff and HSC helpdesk for their help. We would like to thank S. Mineo, H. Furusawa, Y. Yamada, and M. Kubo in HSC helpdesk team for useful discussions regarding the data reduction. We thank NAOJ for providing access to hanaco account which was used to perform some initial stages of data reduction. This work presents results from the European Space Agency (ESA) space mission Gaia. Gaia data are being processed by the Gaia Data Processing and Analysis Consortium (DPAC). Funding for the DPAC is provided by national institutions, in particular the institutions participating in the Gaia MultiLateral Agreement (MLA). The Gaia mission website is https://www.cosmos.esa.int/gaia. The Gaia archive website is https://archives.esac.esa.int/gaia. This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This work made use of various packages of Python programming language. SRD acknowledges support from the Fondecyt Postdoctoral fellowship (project code 3220162) and ANID BASAL project FB210003. SRD acknowledges membership in the BRICS-SO. JJ acknowledges the financial support received through the DSTSERB grant SPG 2021/003850. We gratefully acknowledge the use of high performance computing facilities at IUCAA, Pune for the HSC data reduction. We acknowledge the use of the AI tool ChatGPT in this research work.

Appendix A Candidate populations in the literature

Appendix B HSC images of IC 1396

Appendix C Efficiency test of RF classifier

C.1 For 30^′ region

For the 30^′ region, the training set consists of 435 stars, of which 236 are member stars, and 199 are non-member stars. These are the HSC+UKIDSS counterparts of the Gaia-DR3-based sources obtained by Das et al. (2023). We use only the twenty-four newly derived parameters to train the machine in this step. The reason for using only these new parameters is discussed in Section 5. Of the total 435 stars, 60% are used for training, and 40% are used to test its accuracy. With this setup, we train the machine and see a high accuracy of 0.99. In Fig. C.2, we show the confusion matrix, which demonstrates the high accuracy level of the machine during training. Out of the 174 stars used for testing, the machine successfully identified most of the member and non-member stars. The machine is confused to identify the true nature of only two sources. This indicates the higher accuracy of the machine in training. In Table C.1, we provide the relative importance of the parameters used in this analysis. Here, we only provide the top ten parameters. As discussed in Section 5, we can see the effectiveness of these new parameters in segregating the member and non-member sources.

C.2 For the whole complex

As explained in Section 5, for the entire complex, our training set comprises 27423 stars, with 1250 being member stars and 26173 being non-member or field stars. The RF classifier algorithm efficiently handles large dimensions; therefore, in this analysis, we utilised the complete set of forty-two parameters, as detailed in Section 5. Consequently, we construct the RF classifier using the forty-two-dimensional reliable training set and assess its accuracy. For this purpose, out of the 27423 stars in the input training set, 60% were utilised to train the RF classifier, with the remaining 40% used for testing accuracy. Thus, in this analysis, out of the 27423 stars, 16453 stars were employed for training, and the remaining 10970 stars were utilised to test the efficiency of the machine in distinguishing member stars from field stars. The machine autonomously selects the training and test sets randomly. Our training procedure yields a high accuracy of 0.996 with the RF method. We present the confusion matrix in Fig. C.3, demonstrating the high accuracy of the RF method. This confusion matrix illustrates how the machine identifies sources based on training. From the confusion matrix, we observe that out of the 10970 stars used for testing the accuracy of the machine, and it successfully identifies 10466 field stars and 465 cluster member stars. During classification, the machine encounters confusion with only a few field and cluster member stars. This exercise highlights the effectiveness of the RF method. Table C.2 provides the relative importance of the input parameters as determined by RF while providing the membership probability. Here, we present the top ten parameters. Similarly, the new parameters exhibit high relative importance, indicating their effectiveness in significant segregation between member stars and non-member stars.

Appendix D Details of the comparison with literature candidate stars

Appendix E CMD plots

References

---

All Tables

All Figures