The Emerging Stellar Complex in Mon R2: Membership and Optical Variability Classification

Monoceros R2 (Mon R2) is one of the closest large active star-forming regions. This extremely young and partially embedded region provides an excellent laboratory for studying star formation and the early evolution of young stellar objects (YSOs). In this paper, we conduct an optical study of the greater Mon R2 region. Beginning with 1690 previously identified candidate YSOs, we used 496 sources with good proper motions and parallaxes from Gaia Data Release (DR) 3 to determine the astrometric properties for likely members of Mon R2. We then used both astrometric and photometric (isochronal and variability) criteria to determine that 308 of these stars are highly probable members. Using the same criteria, we considered a broad area search around Mon R2 in Gaia DR3 and separated candidate members from field stars. In total, we selected 651 likely new cluster members that had been missed in the previous X-ray and infrared excess selection techniques used in the past to establish cluster membership. Revised astrometric properties of the cluster were found using the combined sample of ∼959 highly probable member stars. For the literature plus the new candidate member list, optical light curves were compiled from the Zwicky Transient Facility. For 470 identified variable sources, we attempted classification based on the flux asymmetry (M) and quasiperiodicity (Q) metrics. We find that Mon R2 is dominated by quasiperiodic symmetric variables, with aperiodic sources also a significant population. A few tens of large-amplitude variables are identified that may be of interest for further study.


INTRODUCTION
Stars form from the collapse and fragmentation of molecular clouds.These dust-rich gas clouds have enough material to produce anywhere from tens to tens of thousands of stars, often resulting in dense stellar systems known as star clusters.Detailed study of young stellar populations, and their clustering properties, in such active regions can help inform the processes and conditions during star formation from dense interstellar material.Before nuclear fusion begins within its core, a forming star goes through core collapse, a main accretion phase, and finally slow contraction, eventually reaching a quasi-equilibrium.During the early stages of a star's pre-main sequence evolution, the object is denoted as a Young Stellar Object (YSO).
One defining property of YSOs is their spectroscopic and photometric variability (Joy 1945).As a proto-star with a circumstellar disk, the object experiences many dynamic phenomena including accretion of material, ejection of material in outflows and jets, magnetic activity manifest as cool and hot spots on the surface, and occasional outbursts (Semkov 2011;Fischer et al. 2023).These physical phenomena each result in differing amplitudes and timescales of the star's variability and cause brightness fluctuations in different regions in the electromagnetic spectrum.Here, we study the optical variability of YSOs toward the Mon R2 region.
The Monoceros R2 Molecular Cloud, often referred to as just "Mon R2", is one of the closest large active starforming regions from the sun.The traditional distance estimate of 830 ± 50 pc (Racine 1968;Herbst & Racine 1976) was updated based on high-precision VLBA measurements to 893 ± 42 pc by Dzib et al. (2016).The Monoceros R2 region contains several early-type stars and reflection nebulae, and deeply embedded star clus-ters (see Carpenter & Hodapp 2008, for a review).It also contains the YSO outburst source V899 Mon, an EX Lup type object.
The main star cluster in Mon R2 was first identified by Beckwith et al. (1976) and further revealed by Hodapp (1994).Although embedded, it is associated with NGC 2170.A second dense embedded cluster associated with GGD 12-15 is just to the east (Gutermuth et al. 2005), and a third associated with GGD 16-17 is even further east1 .A network of filamentary gas and dust structures spiral outwards from the central cluster (Rayner et al. 2017;Treviño-Morales et al. 2019;Kumar et al. 2022).In this and other examples of hub-and-spoke type structure in the dense interstellar medium, it is thought that the hub maintains star formation through the accretion of material that flows from the filaments to the central hub (Rayner et al. 2017), and that these central regions are more likely to produce massive stars and clusters, relative to star formation in the fragmenting filaments.Other examples of hub-filament systems are presented in e.g.Mookerjea et al. (2023), Nagai et al. (1998), Myers (2009), Wang et al. (2020).
Whereas most of the previous literature on Mon R2 region has presented data and analysis on the gas (e.g.Loren 1977;Pilleri et al. 2012;Ginard et al. 2012;Treviño-Morales et al. 2016;Keown et al. 2019) or the dust (e.g.Pokhrel et al. 2016;Sokol et al. 2019;Hwang et al. 2022), our focus is on the young stellar population.The Mon R2 region has a very low inclination on the sky, allowing us to see the entire cluster face-on (Kumar et al. 2022).Young embedded star clusters provide an opportunity to study both star formation and early-stage stellar evolution.Using photometric measurements of the large sample of young, bright, and active stars in these clusters, we can better understand the surrounding environment and conditions for star formation.
Past work on the properties of stars includes that of (Carpenter et al. 1997;Andersen et al. 2006) in the main Mon R2 cluster and (Maaskant et al. 2011) in the GGD 12-15 cluster.However, little examination has been made of stellar properties of the YSOs in the larger region beyond the main clusters.Detailed work beyond census building (e.g.Gutermuth et al. 2009), specifically the establishment of membership of candidate young stars, has not been done.Nor has there been a large-scale analysis of optical variability using long-term light curves.Mid-infrared variability in the main cluster was characterized by Rebull et al. (2014), and a recent preprint by Orcajo et al. (2023) studies optical variability on relatively short timescales.
In this paper, we reconsider the census of the greater Mon R2 region, starting with a star list based on the existing literature.We then use Gaia DR3 astrometry and photometry (Gaia Collaboration et al. 2022), and light curves from the Zwicky Transient Facility (ZTF; Bellm et al. 2019), to firmly establish the Mon R2 census.We then study the optical photometric variability using the ZTF light curves, which sample timescales from days to years, allowing us to capture the short-term and longterm variability of the YSOs.Using thresholding metrics, we then test whether each star's brightness over time is changing and hence whether it is variable or not.After identifying variable stars in Mon R2 we can then classify the variability properties based on established variability metrics, which can suggest the physical phenomena that are causing the lightcurves to vary.Section 2 describes the assembled data sets and our data culling processes and quality control steps.Section 3 details how we assess cluster membership.Section 4 presents our analysis of the optical light curves.Section 5 contains our discussion and conclusions.

YSOC: Young Stellar Objects Previously
Associated with Mon R2 As described in Hillenbrand (2021), YSOC (the Young Stellar Object Corral) is an online database of nearby young stars containing astrometric, kinematic, photometric, spectroscopic, and other observed stellar properties, as well as associated membership information.The data is compiled from various papers studying young stellar objects from x-rays to millimeter wavelelengths, and complemented with data from survey archives, The site contains stars associated with various different regions, such as the North American-Pelican Nebula, Orion, Ophiuchus, etc.
For Mon R2, the historical literature that has established the stellar census is mainly based on infrared observations, namely the near-infrared surveys from the ground of Carpenter et al. (1997) and from HST of Andersen et al. (2006), and the mid-infrared surveys with Spitzer of Gutermuth et al. (2009) and Kryukova et al. (2012).In addition, X-ray surveys by Kohno et al. (2002), Nakajima et al. (2003) and Getman et al. (2017) also contributed unique candidate members.Our initial list of plausible members of the Mon R2 region is the unique list of 1690 stars compiled from the above sources.A region 2.5 deg × 2.5 deg in size (corresponding to about 36.3 x 36.3 pc 2 , for our adopted distance of 833 parsecs) was queried in Gaia DR3 around the position RA: 06h08m35.5s,Dec: -06d13m42.60s,resulting in a return of 107,557 stars.This provides the material to search for new members of the greater Mon R2 region that were potentially missed by previous work that relied on less diagnostic selection techniques.

Gaia DR3: Astrometry and Optical Photometry
For our Gaia query, we applied parallax quality cuts to eliminate nonphysical parallax (requiring plx > 0) and bad parallax measurements (requiring plxSN R > 2).This left 41,493 stars from our original 107,557 sources.Figure 1 shows the region of Mon R2 we studied, with the sources resulting from the initial Gaia (red) and the sources previously associated with Mon R2 as compiled in YSOC (blue).The right-hand side shows various distributions comparing all Gaia stars in the field with the previously associated stars.

Zwicky Transient Facility: Optical Time Series Photometry
The Zwicky Transient Facility (ZTF; Bellm et al. 2019) is a survey using the Samuel Oschin 48-inch telescope (P48) at Palomar Observatory.ZTF catalogs over a billion stars in the northern sky and measures photometric time series, primarily in the r-band (λ ef f = 6340 Å) and g-band (λ ef f = 4722 Å).Due to the red nature of our sources, we consider only r-band lightcurves due to their higher signal-to-noise.We utilize photometry from ZTF Data Release 12, spanning a 4-year period from March 2018 to May 2022.
We queried ZTF using the Caltech/IPAC InfraRed Science Archive (IRSA).For stars of interest, we used the RA and Dec coordinates from Gaia or, for sources not present in Gaia, from YSOC (which have their provenance in Gaia then 2MASS then WISE).We searched a 1.5 arcsec radius around each star's coordinate for ZTF counterparts.The search radius was chosen to minimize contamination in the lightcurve from nearby stars.We checked the quality of the matches by computing the difference of the RA and Dec coordinates of the star, and the median RA and Dec for all apparitions of the source in the ZTF light curve.After calculating the euclidean separation ( √ ∆RA 2 + ∆Dec 2 ), we examined these differences to validate that each lightcurve corresponds to the intended source.We also examined the stars with separation > 1 ′′ using CD-SPortal 2 to check for contamination by closely projected stars (e.g.cases where the optically visible ZTF source is not the same object as the embedded claimed Mon R2 member).
Separately, we identified cases of different sources in close proximity to one another having the same lightcurve, and assigned the lightcurve to one of the sources if there was no ambiguity, but eliminated objects where no unique attribution could be made.
The retrieved ZTF light curves were considered in our analysis below only if the source had a median magnitude r < 21, with > 25 observations in the time series.This restricts our analysis to the higher-quality light curves.We also note that ZTF saturates at approximately 13'th mag, so the brightest few sources in our sample do not have lightcurves.

CLUSTER MEMBERSHIP METHODS AND RESULTS
Ultimately, we are interested in studying the optical light curves of Mon R2 members with ZTF.In order to establish a membership list, we assembled data and applied cuts designed to leave us with high-quality members.In this section, we describe the steps taken to separate contaminating foreground and background field stars from likely young stellar members of the Mon R2 region.In the subsections below we design and apply criteria based on astrometry, photometry, and photometric variability.We then weight the different membership criteria and select a final list of candidate members from the stars assembled from YSOC and the wide-field Gaia query.

Astrometric and Photometric Properties of Mon R2 Established from the Literature Sample
Because star clusters consist of stellar groupings that are spatially and kinematically associated, their members should be located around the same distance from us, and be co-moving in the plane of the sky along similar directions.Using these assumptions, we can assess the parallax and proper motion values of all the stars in 2 http://cdsportal.u-strasbg.fr the greater Mon R2 field to find stars which are moving coherently, and identify those which have discrepant motions.To establish the astrometric parameters of the cluster we use the list of stars previously associated with Mon R2 as compiled from the YSOC database.
As stated above, we started with 1690 stars associated in the literature with the Mon R2 Region.We find that 994 of these stars lack proper motions or parallax values in Gaia DR3.This is due to the highly embedded nature of the core (or "hub") of Mon R2, where most of the members have been identified, as well as high extinction in some of the "spoke" regions.The stars without astrometric values in Gaia DR3 are dimmer in optical photometric bands compared to those having high quality astrometry.Specifically, based on matching to PanSTARRS, which is deeper than ZTF, we find that at least 90% of the stars with parallax and proper motion measurements have r < 21 mag.
We assert that it is reasonable to calculate the cluster parameters using the optically brighter sources in the cluster, but acknowledge potential bias if the more embedded sources have different kinematics.The stars without astrometry are not eliminated from our consideration, as we can still assess their membership through other criteria.However, we note that a brightness of r < 21 mag was also the limit imposed above on the ZTF data set, so the brighter population is in fact the ultimate target of our investigation.
For the remaining 696 stars from YSOC having parallax values in Gaia DR3, we applied parallax quality cuts to eliminate nonphysical parallax (requiring plx > 0) and bad parallax measurements (requiring plxSN R > 2).This left 496 stars with quality parallax measurements.As a first step in creating membership probability criteria to select likely members of MonR2, apart from field stars, we found the mean and modal astrometric values for this sample of likely member stars, established from the literature.
To calculate the cluster proper motions, we first plot the one-dimensional distributions of Proper Motion (RA) and Proper Motion (Dec), and fit Gaussian functions to extract the mean and standard deviation (Figure 2).It should be noted that we fit only data within the range -10 to 10 mas/yr, for both proper motions, which encloses the medians of the proper motions (µ α,med =-2.86, µ σ,med =0.70), and avoids stars in the initial catalog with extreme values.From our fits, we find that suggested members of the Mon R2 cluster have mean proper motions and standard deviations: • Proper Motion RA: −3.05 ± 0.69 mas/yr • Proper Motion Dec: 0.77 ± 0.61 mas/yr  with errors of ±0.03 mas/yr in RA and ±0.02 mas/yr in Dec.For parallax, we use the same Gaussian fitting process after eliminating stars with parallax error over parallax less than 0.38, or the 90th percentile parallax error over parallax value (Figure 3b).We thus ensure working with stars with both physical and accurate parallax measurements.Figure 3a shows the Gaussian fit, which has mean and standard deviation: • Parallax: 1.20 ± 0.17 mas with error of ±0.01 mas.The inferred distance of 833.33 +103.40  −137.54 parsecs matches well with previously determined distances to the cluster (see references in Sections 1 and 5).
We emphasize that the above kinematic properties are derived from candidate cluster members identified in previous literature.In the next subsection, we use the parallaxes and proper motions of individual stars as one factor in assessing their probability of membership in the cluster.
We also compare photometry for the literature sample with theoretical isochrones.The location of the bulk of the observed distribution relative to computed pre-main sequence isochrones allows us to determine a maximum age for likely members of the cluster.Stars born from the same molecular cloud are expected to lie along sequences in color-magnitude diagrams, though especially for young embedded systems like Mon R2, reddening plays a significant role in the colors and magnitudes of our stars.
In standard color-magnitude diagrams, pre-main sequence stars should lie above and to the right of main sequence stars at the same distance.We use photometry from PanSTARRS (Flewelling et al. 2020) and Gaia (Gaia Collaboration et al. 2016b) to assess this for our Mon R2 sample.Typical of star-forming regions, Mon R2 has an estimated age of around 10 6 to 10 7 years (Herbst & Racine 1976;Carpenter & Hodapp 2008).In Figure 4 we confirm this age range is consistent with our sample by comparing the photometry in different colormagnitude planes to 1, 10, and 100 Myr isochrones from Dotter (2016); Choi et al. (2016); Paxton et al. (2011).The majority of our candidate member stars are in the expected location.However, some are not, and in the next sub-section we use this information to help assess membership likelihood for individual stars.

Establishing Membership Probability
With the properties of probable cluster members established above, we now create a set of criteria to estimate membership probabilities that will be used to select likely cluster stars from among all stars projected on the sky in the vicinity of Mon R2.For each star, we evaluate three membership probability metrics, and calculate a total membership probability using weighted averaging.In addition to the astrometric and colormagnitude information discussed above, we also include photometric variability in the ZTF data as an additional factor in the membership probability weighting.We label our membership metrics as probabilities, though they are not strictly such in a rigorous statistical sense, but instead a simple way of numerically ranking stars as more likely to less likely cluster members.
For our astrometric metric, we use proper motion in RA, proper motion in DEC, and parallax, all relative to the mean cluster values quoted above.We evaluate the standard score (also called Z-score) for each, defined as Z = X−µ σ , where X is the variable of interest, µ is the variable mean, and σ is the variable standard deviation.We then assign a probability of cluster membership, adopting unity for values within 1 standard deviation and 1/Z for values more than 1 standard deviation away from the mean.Finally, we calculate the rootsum-square of the three individual proper motions and parallax Z-scores (Z pmra , Z pmdec , Z plx ) and produce an astrometric probability (P astrom ) by taking the inverse of this as Figure 5 shows the results of this criterion for the Gaia candidates.We isolate stars with P astrom > 0.25, which  indicates that they should be close to the cluster mean value in at least one of the astrometric parameters.We also use color-magnitude diagrams based on Gaia photometry to check whether our stellar candidates are consistent with the age of the cluster and pre-main sequence stars.A G-RP vs RP color-magnitude diagram (Figure 6) is used to assign a membership probability P iso based on the location relative to the isochrones discussed in Section 3.1.In making this assessment, we consider that many of the stars are young stellar objects surrounded by gas and dust, which will affect their colors on an individual basis, and we allow for reddening up to about A V = 5 since the main cluster is still embedded.We assume all stars below the main sequence at the distance of Mon R2 are field stars, in practice allowing the 10 8 yr isochrone as the oldest potential cluster member.Equation 2shows the assigned membership probabilities.
below 10 8 yr isochrone 0.25, between 10 7 yr and 10 8 yr isochrone 1, between 10 6 yr and 10 7 yr isochrone 0.75, above 10 6 yr isochrone (2) Finally, in addition to astrometry and photometry, we assign a photometric variability P var metric for cluster membership.Young pre-main sequence stars are ubiquitously variable (Fischer et al. 2023), and non-variable stars at the photometric precision of ZTF are highly likely to be older than 10-100 Myr, which would mean they are less likely to be members of the young Mon R2 cluster.We assign variability probability membership values based on the measured scatter in the optical light curve (σ mag ) relative to the median error (err mag ) in the data, according to Equation 3.
We thus have astrometric, isochronal, and variabilitybased membership probabilities for individual stars.Some objects have all three available, while others have only two or one of these values established.

Identifying New Candidates Members from Gaia DR3
In this section we apply the membership probabilities defined and assigned above to the initial candidate list from our broad 2.5×2.5 deg Gaia query of 41,493 stars around Mon R2.In order to only focus on the most likely candidates, we first require our sources to have astrometric probability P astrom > 0.25.This cut alone leaves us with 2551 stars, for which we then calculate a weighted average total membership probability based on all three criteria (Equation 4) to find the most highly probable members.We apply weights of w 1 = 3, w 2 = 1, w 3 = 1, to the astrometric, isochronal, and variability probabilities in deriving the final membership probability.These weights are used, since the parameter of astrometry is composed of three separate parameters, in contrast to the remaining probabilities which are composed of one parameter each.
As noted above, not all stars have all measurements.For the sample sourced from Gaia, while all stars have astrometric values, 42 of these do not have G or RP magnitudes, and 97 do not have ZTF lightcurves.As such, for this <5% of the sample, we set the weights to 0 when information is missing so as to not bias the final P total .Figure 7 shows our sample of 2551 Gaia stars with P astrom > 0.25, in proper motion and parallax space, with their final probability membership value, P total indicated. .The distribution of total membership probability P total , as given by the colorbar, across proper motion and parallax space for the 2551 Gaia stars with Pastro > 0.25.Given the weighting of the membership selection criteria, most of the highly probable members lie close to the estimated median kinematic parameters for the cluster.However, the additional criteria based on photometry relative to isochrones, and photometric variability, mean that some stars identified as candidate cluster members are further away from the kinematic centers.
We select those sources with total membership probabilities P total > 0.5 as likely members of Mon R2. Figure 8 shows the distribution of probabilities for each one of our assessed properties (astrometric, isochronal, and variability), and the combined probability used to separate members from non-members.The cutoff value of 0.5 is not well-justified, especially given the fairly flat distribution in P .However, experimentation with different cutoff choices shows that 0.5 represents a balance between being too conservative (e.g.eliminating many objects that are centrally located in clusters and have obvious YSO-like variability) and too liberal (e.g.including stars that are widely distributed across the region and not preferentially adding to the clusters).
The process results in 921 stars that are highly probable members.We illustrate these 921 sources in Figure 9 relative to those previously known from the candidate member list sourced from YSOC.In §3.5 we will assess how many of the 921 are newly identified vs. previously suspected members.

Applying Membership Probability to YSOC Catalog Stars
In addition to the new members we found using our Gaia query of the larger field, we can also investigate the membership probability of the stars in the YSOC catalog to see if there are any outlier candidates whose membership should be re-examined.However, many of the stars in the YSOC catalog do not have astrometric values, due to their optical faintness, but can still be likely members of the cluster.As such, we apply our membership criteria to the 496 stars with good parallax measurements and the 994 stars with no astrometry measurements.Similarly applying our process above for the stars sourced from Gaia, to the literature sample sourced from YSOC, we find 308 stars that are highly probable members with P total > 0.5.Figure 8 also shows the distribution of probabilities for the YSOC stars for each criterion and for the total membership probability, which can be compared to the Gaia-selected sample.

Comparing the New Gaia Members to Previously Known Mon R2 Members
In order to assess which of the 921 Gaia-selected Mon R2 members found in §3.3 are newly determined as such, versus already suspected Mon R2 members from the previous literature, we cross-matched to the list of 308 highly probability members from YSOC that we reconfirmed in §3.4.Stars with a Gaia versus YSOC coordinate separation of < 0.25 ′′ are taken as confirmed matches, and a small number (7) of matched sources with separations 0.25 − 1 ′′ were individually examined and are also determined to be matches.All objects in the Gaia-selected members list more than 1 ′′ from any YSOC source were taken as new members.We thus find that 270 of the 921 Gaia-selected members were previously known literature members.We have therefore newly identified 651 sources from Gaia that match the kinematics of the Mon R2 region.

Conclusions on the Kinematics of Mon R2
Thus far, we have used a catalog of previously associated stars from the YSOC database in order to estimate the parameters of the Mon R2 cluster and establish a set of astrometric, isochronal, and variability membership criteria.We then used these criteria both to reconfirm the highly probable members in our YSOC list, while also identifying likely nonmembers, and to search for new highly probable members in the field discovered using the Gaia DR3 catalog.Lastly, we crossmatched our new Gaia stars with the YSOC stars to eliminate re-   dundant sources.In the end, we are left with 959 highly probable optical cluster members, 308 from YSOC and 651 from Gaia.Table 1 shows the selection criteria along with the number of stars that remained at each step, leading to our final list.Table 4 shows our final members list, including the membership probabilities for all applied criteria.In this sub-section, we use our final list of highly probable members to re-estimate the kinematic parameters of Mon R2.By repeating our methods of §3.1, where we fit one-dimensional Gaussians to the proper motions and the parallax (Figure 10), we find the following means and standard deviations: • Proper Motion RA: −3.16 ± 0.78 mas/yr • Proper Motion Dec: 0.62 ± 0.53 mas/yr • Parallax: 1.17 ± 0.13 mas with errors of ±0.04 mas/yr in RA, ±0.02 mas/yr in Dec, and < 0.005 mas in parallax.The parallax corresponds to a distance of 854 +85 −107 parsecs.We note that our elimination of some previously claimed Mon R2 candidate members, and our addition of new, previously unassociated, stars has left the standard deviations on the proper motions and parallax the same or even smaller than reported above for the literature sample.

OPTICAL VARIABILITY METHODS
In this section, we classify the optical variability of the stars we have determined are members of Mon R2, using light curves from ZTF.We use both variability metrics and a periodicity search to select a subset of the cluster members list for further analysis.We then classify these stars based on the Q-M classification established by Cody et al. (2014).
We work only with 889 stars from the 959 star membership list, based on requiring the ZTF light curve to have >25 observations so that we can reliably calculate the variability metrics.Furthermore, while we used variability above as one of the criteria in assessing new member candidates, with a relative weighting of only 20%, there are some candidates with high membership probabilities that turn out to be not formally variable, and for which we do not pursue the Q-M classification.
Additionally, the median of the magnitude error over the time series (err mag ) is calculated to identify baseline noise for the variability assessment.Our metrics are those established by the ZTF Source Classification Project (Coughlin et al. 2021).
For some light curves, there are outlier points that are inconsistent with the rest of the lightcurve.In most cases, these are spurious noise that would bias the calculation of many of our metrics.To identify such contamination, we use the kurtosis metric or a measure of the "tailedness" of the data distribution, selecting all light curves with kurtosis > 4 for examination.We check whether the high kurtosis is caused by a few outlier points, or more systematic outliers that are likely astrophysical.For example, large amplitude variables with short timescale changes would have high kurtosis values.In total, there are 25 stars with high kurtosis due to outlier points; these sources have their variabil- ity metrics recalculated using only observations within three standard deviations of the median brightness.For the stars with high kurtosis that is likely due to astrophysical outliers, we keep all data points of the light curve.
We classify objects as variable stars if both of the following criteria are met.(1) the standard deviation (σ mag ) of the light curve is larger than twice the median of the magnitude error (err mag ), and (2) the inverse vonneumann statistic is less than unity.Figure 11 shows the distribution of median magnitude vs median magnitude error (black) and standard deviation (green) for all our stars.Stars with a red outline are not considered variable as they do not fit both our criteria.

Period Search
We used the astropy LombScargle periodogram package to find and analyze periodic signals in our ZTF light curves.We set our initial period search to find periodicity between 2 and 250 days.This range avoids the cadences from Earth's rotation (daily) and motion around the sun (annual), but special attention does need to be paid to possible confusion of astrophysical timescales with the lunar orbit cadence (≈30 days).The resulting periodograms are used to find the most likely period and its significance, as quantified by the standard false alarm probability (FAP).
Of the 889 stars for which we performed the period search, we found 466 significant (FAP < 0.01) periods.For these stars, we visually examined all observed and phased light curves, and determined whether the periods were justified, or if an expanded period search below 2 days or above 250 days would be informative.For stars where none of the plausible periods seemed appropriate, but the star met the variability criteria set in §4.1, we manually label the source as Aperiodic and assign it the timescale having the highest power.
We ended up eliminating 59 sources with approximately monthly or 250+ day timescales that seem unreliable.We also readjusted the periods of approximately 32 stars, some to ranges < 2 or > 250 days when such periods provided better phased light curves than those from our initial search.Figure 12 shows the distribution of significant/believable periods.The majority of stars in our sample have a period < 10 days, suggestive of typical stellar rotation timescales for a young population.All sources with significant periods are considered for our variable star sample, even if the object does not meet both variability criteria articulated in subsection 4.1.Ultimately, we retain those added based on their significant periodicity only if they also have a low value of the quasi-periodicity metric Q, as described below.

Q and M Metrics and Variability Classification
For Q-M classification, we require the star to either be a statistically significant photometric variable, defined above based on the lightcurve having a standard deviation more than twice the median error, and an inverse von Neumann statistic ν < 1, or the lightcurve having a significant period, with false alarm probability of less than 1%.
If a star has a period that is deemed not significant, and it does not meet the variable star criteria, it is removed from our list for Q-M classification.In total, among 889 stars, we end up with 470 stars that are variable, significantly periodic, or both.
For each of these stars, we used the Q and M metrics first codified in Cody et al. (2014), and also applied in Cody & Hillenbrand (2018) and Cody et al. (2022) to young star lightcurves that were obtained from spacecraft (CoRoT and K2, respectively).Application of this same technique to ZTF data by Hillenbrand et al. (2022) was taken as proof of concept for the present study.Q − M classification of young star light curves results in a parameter space of periodicity to aperiodicity (Q), and fading to bursting lightcurves (M ).
The Flux Asymmetry, or M metric, measures whether a light curve has brightened or dimmed over time, and is defined as: where ⟨m 10% ⟩ is the is the mean magnitude of the top 90% and bottom 10% percentile m med is the median magnitude, and σ mag is the standard deviation of the magnitude.
The Quasi-Periodicity, or Q metric, is a measure of whether a light curve contains purely periodic signals, quasi-periodic signals, or aperiodic signals.To measure Q, we apply the formula: by first determining the dominant timescale identified in the Lomb-Scargle period finding process described above, and subtracting this purported periodic signal from the observed lightcurve to produce the σ resid term.The σ phot term is the mean photometric error in the lightcurve and σ mag is the standard deviation of the magnitude measurements.
The initial run produced a small number of stars with Q values outside the range of 0 to 1. Upon visual inspection, it was determined that these sources had either a small number of measurements, or did not meet the formal variability criteria and had been included based only on their low-FAP periods.As these periods were close to the lunar cadence, all such sources were eliminated from the Q − M analysis.
Stars with M values > 0.25, have predominantly dimming light curves and are considered dippers, while stars with M values < -0.25 are predominantly brightening and considered bursters, while M values between -0.25 to 0.25 have relatively symmetric lightcurves and are considered symmetric.Stars with low Q values have a more strictly periodic signal, with most of the lightcurve subtracting out well from the smoothed periodic signal while a star with a Q value of 1 would be very stochastic or aperiodic, having a large signal even after subtraction of the periodic signal.We considered stars with Q < 0.45 as periodic, 0.45 < Q < 0.82 as quasi-periodic, and Q > 0.82 as aperiodic.These boundaries were validated by visual inspection.We note that very few significant periods (FAP < 1%) are found at Q > 0.48, and stars with such high Q should be taken to have "timescales" rather than true periods.
With M and Q calculated for each star, we assigned the Q − M categories as in Cody et al. (2014): Periodic (P), Quasi-Periodic Symmetric (QPS), Quasi-Periodic Dipper (QPD), Aperiodic Dipper (APD), Stochastic (S) -which corresponds to an aperiodic-symmetric signaland Burster (B), a class that includes both aperiodic and quasi-periodic bursters.We then manually inspected all of the phased lightcurves to validate the dominant timescale reported from the Lomb-Scargle analysis.If we determined in the manual inspection that the chosen period was inaccurate, we looked for other (usually shorter) plausible periodogram peaks that phased well, and reassigned the Q value, and possibly the classification between aperiodic, quasi-periodic, and periodic.If none were found, we removed the Q value calculated for the star and assigned it a classification in one of the aperiodic categories.
Figure 13 shows an example star for each of the classification groups, with both the entire light curve and the phase-folded version illustrated.The distribution of our final Q − M classifications for Mon R2 is shown in Figure 14 and the complete list of sources and their Q − M classes is given in Table 5.Table 2 provides a summary of the different variable star classes.

Stars with Additional Timescales
In our periodograms, the most significant (tallest) peak is used to calculate the Q metric.However, some lightcurves show additional timescales that we describe in this section.We flagged stars that we considered to have long-timescale or truly multi-periodic behavior, and added a secondary class to the primary Q−M classification indicating multi-periodic (MP) or long timescale (L) behavior.
We identify potential multi-periodic signals by flagging stars with more than one peak above the 1% significance power threshold.To validate such candidates, we extract all significant periods from the periodogram and manually inspect each phase-folded light curve.If there are multiple periods that phase well and are not aliases of one another, we give the light curve a secondary classification of Multi-Periodic (MP).In total, we have only two stars with the secondary classification of MP.
To identify long timescale stars, we searched for light curves that have systematic changes in brightness over a significant portion or even the entirety of the fouryear lightcurve time span.The stars with a secondary classification of long timescale (L) are all named as such through manual inspection.
Figure 15 shows one example of each type of secondary classification.

Correlations with the Q-M Classes
In this section, we consider how well the Q − M variability classes align with standard statistical metrics of variability.We also investigate the variability behavior of the population with respect to physical location in the Mon R2 region.
For each source, we first simplified its variability classification to just the M metric designation as a burster, symmetric, or dipper.Then, we examine the variability metric distributions across the simplified classes.
We find that dippers have on average a lower inverse von-Neumann statistic than symmetrics and bursters, making them more readily identified as variables.This is further evident from examining the variability level as measured by standard deviation over median magnitude error, where we find the most variable stars are also most likely to be dippers, while symmetrics make up the majority of the least variable stars, likely due to their relatively low amplitudes.The burster sources tend to have dimmer light curves, on average, than the symmetric and dipper sources.It may be that these stars are picked out as variables based on the burst behavior, whereas comparably faint stars with symmetric or dipping variability would be harder to identify as variables based on small-amplitude changes relative to their median brightness especially given that the noise levels are increasingly higher for fainter stars.We also note that there are no dippers or bursters with Q < 0.4; in other words, all of the truly periodic sources at low Q have relatively symmetric This is consistent with other studies, and indicative that the mechanisms causing the aperiodicity in the lightcurves also cause the skew, measured here as elevated values of the |M | metric.In these cases, the presumed culprit is circumstellar material.
Another trend in the variable star population is that dipper sources appear to be concentrated toward the central cluster of Mon R2, and toward the clusters found to the east and to the north.Under the assumption that dipper behavior is likely associated with circumstellar material, which is known to be correlated with stellar youth, the clustering of the dipper population supports another common hypothesis that clustering is related to stellar youth.

Large Amplitude Variables
While examining our dataset, we identified several tens of stars with light curves exhibiting substantial changes in magnitude.Some of these stars are well known and have been studied in detail based on their previously identified erratic behavior.A few such examples include: V899 Mon (Ninan et al. 2015;Park et al. 2021), V1818 Ori (Chiang et al. 2015), and [CMD97]-1031 (Jiang & Hillenbrand 2022).We labeled these and other sources as Large-Amplitude Variables (LAVs) by manually inspecting all lightcurves exhibiting peak-topeak magnitude changes of at least one magnitude, and having standard deviation greater than three times the median magnitude error.
In total, we identify 52 stars that are LAVs.Their full light curves are illustrated in Appendix A, segregated by the simplified burster, symmetric, and dipper classes according to their M values.There are 8 bursters, 12 symmetrics, and 32 dippers.

DISCUSSION
In this paper, we have presented new information on the optically visible stellar population of Mon R2.First, we have used Gaia DR3 to assign membership probabilities to stars known in the literature as potential members.Second, we have used ZTF to study the variability properties of these members.Our discussion centers around the topics of clustering in the Mon R2 region, and comparison to other studies of the Mon R2 population, including both kinematics and variability.

Clustering in the Mon R2 Region
Using our final member list, we examined the kinematic properties of three previously known sub-clusters in Mon R2.Cluster 1 is Mon R2 main, Cluster 2 comprises GGD 12-15, and Cluster 3 is the GGD 16-17 region.We defined three adjacent 0.5 x 0.5 deg rectangular regions and examined the proper motions and parallax distributions of the stars within each group.Figure 16 and Table 3 shows locations of the sub-clusters and the kinematic results.
The parallax distribution between the three groups is the same.However, the proper motions of stars in Cluster 1 are more positive in each direction than the motions measured for stars in Cluster 2 and Cluster 3.These values overlap within one standard deviation, which is expected given that we selected our members to have similar astrometric values.However, we do see that the majority of stars in subcluster 1 are moving on average in a different direction from stars in subclusters 2 and 3.This indicates both shared and distinct submotions of our final member stars, whose origins and futures could inform star formation histories and evolution in this region.We also looked at Gaia photometry for the subcluster stars, examining their relative distributions in the RP vs G-RP magnitude-color diagram.All three subclusters overlap significantly, meaning that they have similar ages.The more distributed population, that is, stars that are not part of any of the three subclusters, have somewhat lower luminosity.According to the isochrones, these non-clustered stars in the region are an older population than the more clustered stars.We explored the variability class population in each subcluster and found slightly more dippers proportionally in the central cluster (cluster 1, with 27.5%) than the other cluster regions (cluster 2 with 17.3%, and cluster 3 with 20%).This could be due to the more embedded nature of the central cluster region, as explained before, which causes more dust and material to block the light from the young star.

Comparison of Kinematic Results
In this section, we compare our kinematic results to two other works that have conducted kinematic studies on populations in Mon R2.Kuhn et al. (2019) used a sample of 97 candidate member stars identified in previous surveys and located towards the core of Mon R2 (the center of our Cluster 1).These authors calculated a weighted median value for proper motion RA = −2.91 ± 0.11 mas/yr, proper motion Dec. = 1.05 ± 0.18 mas,/yr, and parallax = 1.06 ± 0.04 mas, corresponding to distance = 948 +42 −38 pc.More recently, Orcajo et al. (2023) examined stars over the central 26 ′ × 26 ′ region in the Mon R2 region using the Las Cumbres Observatory Global Telescope Network (LCOGT).They matched stars they found in LCOGT to Gaia EDR3 stars, and selected a large rectangular region in PMRA-PMDec space defined by (-2.75, 1.15) mas/yr ± 1.25 mas/yr to study photometric variability.For parallax, they performed a Gaussian fit to the distribution of stars within their proper motion box, and found values within 0.888 to 1.42 mas as likely members, corresponding to distance = 825 ± 51 pc for the cluster.
The region examined in our work is much larger (2.5 deg × 2.5 deg) than the works mentioned above, both of which encompassed smaller regions centered on our Cluster 1.Although our membership methods differ substantially, and our parent sample size and is much larger, we find that our measurements of proper motion and parallax are consistent and contained within one standard deviation of the values found in these previous works.Comparing only our Cluster 1 values, we find even better agreement between our kinematic values and those found previously.

Comparison of Variability Results
Regarding variability properties of stars in Mon R2, we can compare our results to two other studies that have also used ZTF data to study the variability of stars in young star-forming regions using the Cody et al. (2014) Q − M statistics.Caveats to any such comparison include the different membership selection criteria used among the studies, as well as the different distances and extinctions of the clusters, which could mean that somewhat different populations are being probed with ZTF photometry.
We find that our population of stars in Mon R2 has similar proportional distribution among the Q − M classes as found in the ZTF study by Hillenbrand et al. (2022) of the North America and Pelican (NAP) Nebulae region.However, there are some differences in results between our study and that of the Perseus (Per) Molecular Cloud by Wang et al. (2023), also using optical photometry from ZTF.
Specifically, our population of bursters (12.5%) and aperiodic dippers (11.7%) is comparable to those found in the NAP region (B: 13.9% / APD: 9.6%).For Quasi-Periodic Dippers, Mon R2 (9.7%) has a more similar population to Per (7.5%), than the larger population of this variable type found in the NAP region (19.2 %).We find that there is a much lower population in Mon R2 of purely periodic variability, and a much greater population of stochastic variability compared to NAP and Per.Among all three regions and studies, Quasi-Periodic Symmetric dominates the variability classification.This is likely due to the optical photometry picking up days-to weeks-long variability that is modulated on the stellar rotation timescale, whether originating from the stellar surface itself, or in circumstellar regions near the co-rotation radius in the disk.
Next, we can compare our variability results to a recent search for variable stars also conducted in Mon R2 by Orcajo et al. (2023).These authors present data spanning 23 days over which 1500 photometric observations were aquired.Our ZTF data, by contrast, covers a much longer baseline (4 years) but with fewer samples (33-570 photometric observations per star).We crossmatched our 470 member variables to the 136 sources Orcajo et al. (2023) designated as having "rotationmodulated variability", and compare the classifications and periods.We find some apparent differences in astrometry, necessitating use of a 5 arcsec radius to match sources.We found 59 stars in common between the two studies, all of which were visually examined to ensure the same astrophysical source was being compared.Among the 59 matches, only 39 of the reported periods are within 1 day of each other.A handful of stars in the "rotation-modulated variability" list have reported periods of one day, which we had disregarded among our sample in order to avoid the diurnal signal.Further investigation showed that all of the sources with mismatched periods have very high Q values, making them aperiodic sources with "timescales", rather than truly rotationally-modulated behavior which is expected to be sinusoidal.For the period-mismatched stars, we find that using the period reported by Orcajo et al. (2023) results in a more aperiodic (higher Q) signal from our ZTF data, and as such we believe that our periods (or timescales) are more correct for our dataset.Out of the 59 matched stars, only 6 are classified by us as purely periodic, or the equivalent of a rotation-modulated variable, while another 29 are defined in our analysis as "Quasi-Periodic Symmetric", which is related to stellar rotation but likely measuring phenomena in the corotating circumstellar environment.

CONCLUSION AND SUMMARY
To summarize, we gathered 1690 candidate Mon R2 members that had been identified in previous literature and examined the proper motion and parallax properties of the young stellar population.We derived their mean kinematic values.We then created a weighted membership scheme for testing new candidate members that included astrometry, location in color-magnitude diagrams, and optical variability.
From 107557 field stars queried in Gaia DR3 over a 2.5 deg × 2.5 deg region centered at Mon R2, we thus selected 921 highly probable member stars having the expected kinematic and photometric behavior of the young stellar population of Mon R2.We then reassessed the membership for our 1690 YSOC literature candidates using the same techniques, and crossmatched the highly probable members with the stars sourced from Gaia.Including both the previously proposed members sourced from the literature, and the newly identified members from the broader field, we confirmed a total of 959 Mon R2 members in this study.Due to the optical selection of our sources, our membership list is incomplete and does not include more embedded members, mainly those fainter than r > 21.
Our compiled Mon R2 member sample has the following kinematic properties: • Proper Motion RA: −3.16 ± 0.78 mas/yr • Proper Motion Dec: 0.62 ± 0.53 mas/yr • Parallax: 1.17 ± 0.13 mas Furthermore, our study has found for the first time that the previously known sub-clusters in the region have distinct kinematic properties.We studied the proper motions among three sub-clusters in the region, the central (hub) cluster of Mon R2, and two additional clusters located to the west of it.The central cluster kinematics we measure here agree well with values reported in two previous works that study this region.However, we additionally find different proper motion directions for the two western clusters (Cluster 2 and Cluster 3) than for the central cluster (Cluster 1), which impacts the kinematic values we measure for the Mon R2 region overall.Figure 16.Left: RA-Dec and PMRA-PMDEC spaces for Mon R2 member stars that are not clustered (black) and associated with the three sub-cluster regions (red, green, blue).Sub-clusters are chosen by grouping within equal area 0.5 x 0.5 deg regions, centered at coordinates provided in Table 3 Right: Density distribution for proper motions and parallax of the three clusters (red, green, blue curves) compared with those for all member stars (black dotted curves).
For our compiled member list for Mon R2, we then collected optical r-band lightcurves from ZTF.We identified 470 sources with sufficiently high photometric sampling in ZTF that could be classified as variable based on excess scatter (σ mag > 2 × err mag ) or significant periodicity (Lomb-Scargle False Alarm Probability < 1%) sources.We then classified the optical variability using Q − M classification which divides stars by their lightcurve repeatability and flux asymmetry.
We find that the majority of the variability occurs at timescales of less than 10 days, and is dominated by Quasi-Periodic Symmetric (32.5%) and Stochastic (25.1%) variability phenomena.We also find that our population contains 11.7% (55) Aperiodic Dippers, 12.5% (59) Bursters (Aperiodic and Quasi-Periodic), 8.2% (39) Periodics, and 9.7% (46) Quasi-Periodic Dippers.In addition, we find by visual inspection that 7.8% (37) of our stars have a longtime scale variability, while just 0.4% (2) exhibit possible multi-periodic variability.We also find that dipper sources (M > 0.25) are on average more variable and brighter than their burster counterparts.Finally, we identified a list of 52 largeamplitude variable stars whose fluxes go through greater than 1 magnitude of change over the span of the entire light curve.Through this work, we have expanded and updated the stellar census in Mon R2 through the use of kinematics from Gaia DR3 and optical photometry from Gaia and ZTF.We also examined the optical variability in Mon R2 in detail.The results of our study can be used to further understanding of the evolution of young stellar populations.Note-m med -median magnitude of ZTF light curve mavg -average magnitude of ZTF light curve errmag -median of the magnitude error σmag -standard deviation of the magnitude m skew -skew or measure of symmetry of magnitude χ 2 -Chi-Square of the magnitude η -inverse von-neumann statistic m kurt -kurtosis or "tailedness" of light curve, high kurtosis sources have light curve data points cut to only within 3 standard deviations of median (eliminate errors) Period -highest peak on Lomb-Scargle Periodogram, verified by visual inspection FAP -False Alarm Probability calculated by Lomb-Scargle Peridogram M -M-metric, or measure of Flux Asymmetry, calculated using variability metrics, Equation 5Q -Q-metric, or measure of Quad-Periodicity, calculated using variability and Period metrics, Equation 6Primary -primary QM classification based on Q and M value, See Table 2 for all names Secondary -secondary classification assigned through visual inspection (Longtime/LT or Multiperiodic/MP)

Figure 1 .
Figure 1.Left: Spatial distribution of sources considered in this study.Blue points show stars previously associated with the Mon R2 region from infrared and x-ray surveys, that are mainly concentrated in the central cluster and the two groupings to the east: GGD 12-15 and GGD 16-17.Red points result from our initial Gaia query and show the broad distribution of the optical sources that we tested for potential cluster membership.Right: histograms in proper motion RA, proper motion Dec, and parallax with red indicating the Gaia stars and blue bars showing the ±3σ range of the previously associated stars.Lower right panel is a color magnitude diagram in Gaia G-RP versus RP filters for the same two populations.Gaia is a European Space Agency Mission (Gaia Collaboration et al. 2016a) measuring high-precision astrometric, photometric, and spectroscopic data for over 1.5 billion sources.Gaia Data Release 3 (DR3) (Gaia Collaboration et al. 2022) occurred in June 2022 and provided updated position, proper motion, and parallax information for magnitudes as faint as G ≈ 21 mag, as well as photometry.A region 2.5 deg × 2.5 deg in size (corresponding to about 36.3 x 36.3 pc 2 , for our adopted distance of 833 parsecs) was queried in Gaia DR3 around the position RA: 06h08m35.5s,Dec: -06d13m42.60s,resulting in a return of 107,557 stars.This provides the material to search for new members of the greater Mon R2 region that were potentially missed by previous work that relied on less diagnostic selection techniques.For our Gaia query, we applied parallax quality cuts to eliminate nonphysical parallax (requiring plx > 0) and bad parallax measurements (requiring plxSN R > 2).This left 41,493 stars from our original 107,557 sources.Figure1shows the region of Mon R2 we studied, with the sources resulting from the initial Gaia (red) and the sources previously associated with Mon R2 as compiled

Figure 2 .
Figure 2. Top: Distribution of proper motion in RA (left) and proper motion in Dec (middle) for YSOC Catalog Stars.Gaussian fits to the distributions are overplotted in red, with ±3σ values marked with orange lines.The vector point diagram (right) shows the proper motion distribution, with an ellipse indicating the 3σ contour.Bottom: Same plot as above but zoomed in around the main distribution.

Figure 3 .
Figure 3. Left: Gaussian fit to parallax values of YSOC catalog stars with parallax error over parallax larger than the 90th percentile value (0.375).Right: Parallax versus inverse signal-to-noise in the parallax measurements; The boxed region shows data within 3σ of Gaussian mean.Bottom: Same plot as above but zoomed in around the main distribution.

Figure 4 .
Figure 4. Sloan filter i vs i − z, i vs r − i, and i vs g − i absolute magnitude-color diagrams for YSOC stars.Isochrones are from MIST (Choi et al. 2016; Dotter 2016; Paxton et al. 2011) with parameters [Fe/H] = 0, v/v critical = 0.4, and Av = 0. Red dashed line on g − i diagram shows cutoff for stars we deem to be likely background due to high reddening.

Figure 5 .
Figure 5. Left: Vector point diagram of proper motion distributions among all Gaia DR3 star over the 2.5 × 2.5 deg 2 query area centered on Mon R2.Colorbar indcates astrometric membership probability (Pastrom) as described in the text.Middle: Stars with astrometric membership probability Pastrom > 0.25.Right: Field stars excluded from the candidate member sample as a result of low astrometric membership probability.

Figure 6 .
Figure6.Isochrones compared to RP absolute magnitude and G-RP color for all stars in the initial Gaia query.Stars are colored by their location relative to the isochrones, with green points representing stars above the 10 6 year isochrone, cyan between 10 6 and 10 7 year, red between 10 7 and 10 8 year, and blue points stars below the 10 8 year isochrone.

Figure 7
Figure7.The distribution of total membership probability P total , as given by the colorbar, across proper motion and parallax space for the 2551 Gaia stars with Pastro > 0.25.Given the weighting of the membership selection criteria, most of the highly probable members lie close to the estimated median kinematic parameters for the cluster.However, the additional criteria based on photometry relative to isochrones, and photometric variability, mean that some stars identified as candidate cluster members are further away from the kinematic centers.

Figure 8 .
Figure8.Membership probability distribution for each membership assessment technique shown for stars from the YSOC catalog (left) and the Gaia query (right).The top three panels show the contributions from the astrometric, photometric isochronal, and lightcurve variability assessments, with the combined final result in the bottom panels.We consider all stars with P > 50% in the combined histogram (red vertical line) to be likely members of Mon R2.

Figure 9 .
Figure9.Same as Figure1, but now showing only the 921 Gaia-selected sources that have membership probability P total > 0.5.Blue points are members already suspected as such from the literature, while red points are members of Mon R2 identified from Gaia.

Figure 10 .
Figure 10.Gaussian fits (red lines) and 3σ limits (orange lines) for proper motion in RA (Left), proper motion in Dec (Middle), and parallax (Right) for the final catalog of highly probable literature plus newly identified member stars.

Figure 11 .
Figure 11.Left: Black points indicate median magnitude, m med versus median error, errmag, in the ZTF data.Red and green points indicate median magnitude versus lightcurve standard deviation, σ with red outline indicating stars that do not meet a threshold of standard deviation, σmag, divided by median magnitude error, 2 × errmag, larger than two, and green those that do.Right: excess variability, plotted as standard deviation (σmag) minus twice the median magnitude error (2 × errmag) versus inverse von-Neumann statistic.Variable stars are selected as those having inverse von-Neumann <1 and σmag > 2 × errmag

Figure 12 .
Figure 12.Distribution of most significant Lomb-Scargle peaks for all stars.Left panel shows the full distribution; center panel a zoom-in on the lower range of values that are most likely to be true rotation periods rather than a mix of rotation and other timescales associated with young star variability; right panel shows a logarithmic format of the full period distribution.

Figure 13 .
Figure 13.Example sources with different Q − M classification, as labeled.Top: Observed light curves for the entire four seasons of ZTF data; error bars are shown.Bottom: Phased light curves using the highest peak in a Lomb-Scargle periodogram (red points) with binned median magnitude (blue points); error bars are shown.Q and M values are indicated.

Figure 14 .
Figure 14.Flux Asymmetry (M ) versus Quasi-Periodicity (Q) plotted for our sample of 470 variable member sources.Primary variability classes are shown by the marker indicated in the legend.Sources with a secondary classification (multi-periodic or long timescale) are marked in addition to their primary class.Dotted lines show the limits in Q and M for the differing classes, with Q = 0.45 and Q = 0.82 dividing the periodic from quasi-periodic and quasi-periodic from stochastic sources, and M = −0.25 and M = +0.25 dividing the burster from symmetric and symmetric from dipper sources.

Figure 15 .
Figure 15.Same as Figure 13, but example sources with secondary classification of long timescale (L) or multiperiodic (MP).

Figure 17 .
Figure 17.Light Curves with Error for all Large Amplitude Variables that have Burster Classification (M < -0.25) ZTF 92.33017624 -6.69885234 or V899 Mon has been studied in more detail in Ninan et al. (2015); Park et al. (2021) ZTF 88.42731129 -10.40019459 or V1818 Ori has been studied in more detail inChiang et al. (2015)

Table 1 .
Sample size as a function of membership cuts for the literature sources (YSOC) and newly identified sources (Gaia)

Table 2 .
Summary of QM classification for Variable and Periodic Stars in Mon R2

Table 3 .
Table of Positions and Astrometric values of stars in each subcluster and for all stars in our sample

Table 4 .
Table of all final members from Literature (YSOC) and new query (Gaia) DESIGNATION RA Dec PMRA PMDEC PLX PISO PAST RO PV AR PT OT AL Note-DESIGNATION -object identifier from either previous literature obtained from YSOC or Gaia DR3 number PISO -probability membership value based on location on Gaia photometry isochrones PAST RO -probability membership value based on astrometry value in comparison to cluster's median values PV AR -probability membership value based on variability of object PT OT AL -total membership probability calculated with weighted average For more detail see Section 3

Table 5 .
Table of all variable stars and their variability metrics