Free Access
Issue
A&A
Volume 640, August 2020
Article Number A128
Number of page(s) 19
Section Catalogs and data
DOI https://doi.org/10.1051/0004-6361/201936984
Published online 28 August 2020

© ESO 2020

1. Introduction

Early studies about stars in the premain sequence (PMS) phase (e.g., Joy 1945) have recognized photometric variability as one of the main characteristics of young stars. Exploring PMS variability in different regions of the electromagnetic spectrum provides insights into the physical processes at work in these systems. Optical and infrared (IR) variability studies as the Young Stellar Object (YSO) Variability Spitzer Space Telescope program (YSOVAR Morales-Calderón et al. 2011; Rebull et al. 2014; Wolk et al. 2018) and the Coordinated Synoptic Investigation of NGC 2264 (e.g., Cody et al. 2014; Guarcello et al. 2017, 2019) have shown that most of the variability of a YSO is related to physical phenomena happening in the disk or due to the star-disk interaction. In disk-bearing PMS stars, the star-disk interaction and the accretion process causes variability due to the rotational modulation of hot spots at the stellar surface, which are produced at the base of the accretion column (e.g., Sousa et al. 2016; Venuti et al. 2017), instabilities in the accretion disk (e.g., Bouvier et al. 2007) or in the accretion shock (e.g., Koldoba et al. 2008), variable accretion and the evolution of hot spots (e.g., Vrba et al. 1988; Fernandez & Eiroa 1996; Scholz et al. 2009), and the obscuration of the PMS stellar photosphere by circumstellar material, as is the case of AA Tauri type stars (e.g., Bouvier et al. 2007; Alencar et al. 2010; Fonseca et al. 2014). Variability in PMS stars is not only restricted to those stars with disks. Rotational modulation by cool spots, caused by magnetic activity, can produce periodic optical and infrared photometric variability in both disk-bearing and diskless PMS stars (e.g., Artemenko et al. 2012; Grankin 2013).

In addition to probing the physical processes occurring during the PMS, variability can also be used to identify young stars. Knowing the typical variability characteristics of disk-bearing stars beforehand can lead to their identification in contrast to the low variability of field stars. This is especially powerful in the faint limits of optical and infrared surveys, where they are incomplete, and the young stars become blended with reddened background stars. In recent years, various variability surveys have used near-IR to unveil disk-bearing stars both in clusters and associations (e.g., Rice et al. 2012; Wolk et al. 2013; Meng et al. 2019) as well as in the field (e.g., Lucas et al. 2017). Recognizing the typical variability of disk-bearing stars requires determining the occurrence and prevalence of different variability mechanisms. Most of what we currently know about near-IR variability in the PMS is based on two studies by Carpenter et al. (2001, hereafter CHS01, JHKs survey with 2MASS) and by Rice et al. (2015, hereafter RRW15, JHK survey with the UKIRT), which altogether characterize the near-IR variability of ∼1800 young stars in the field of the Orion Nebula Cluster (ONC). In this paper, we contribute to the subject by characterizing the near-IR variability of ∼2500 PMS stars in the young and massive association Cygnus OB2.

Cygnus OB2 (hereafter CygOB2) is an OB association 1.33 pc away from the Sun (Kiminki et al. 2015) with estimated ages between 1 and 7 Myr for its stellar population (Hanson 2003; Wright et al. 2010). CygOB2 harbors a rich and well studied population of massive stars counting with more than 150 OB stars (e.g., Knödlseder 2000; Comerón et al. 2002; Rauw et al. 2015; Wright et al. 2015; Kiminki et al. 2015), along with a low mass population of several thousand stars, which has been uncovered in the last decade by studies with Chandra and Spitzer (Guarcello et al. 2013; Wright et al. 2014). Previous variability studies in CygOB2 have focused on its OB population (e.g., Pigulski & Kołaczkowski 1998) and before this survey, the only variability survey partially concerning the low mass population was performed by Henderson et al. (2011, hereafter HSP11) in two small fields where they observed in the R and I optical bands over two seasons of 19 and 18 observed nights and identified 121 variable stars. However their study only concerned the brightest low mass stars.

In this paper, we investigated the near-IR variability of 5083 previously known low mass members in CygOB2. This work is complementary to the study presented in Roquette et al. (2017, hereafter Paper I), where we focused only on the stars with periodic behavior. Here, we look at all types of variability, regardless of periodicity. We found that 2529 of the known CygOB2 low mass members are variable to a significant level in the JHK bands, and we present their variability characteristics and major trends. Not only the results we present are an important step toward understanding the low mass population of CygOB2, and they also represent the largest single body of near-IR variable stars studied in the literature by almost doubling the number of such variables studied in the ONC by CHS01 and RRW15.

This paper is organized as follows: We summarize the near-IR data used in Sect. 2 and present the methods applied for identifying and studying variable stars in Sect. 3. We investigated variability in terms of Stetson Variability Index (Sect. 3.1), amplitudes (Sect. 3.2), light curve morphology (Sect. 3.3), and variability in the color-space (Sects. 3.4.1 and 3.4.2). We discuss our results in Sect. 4. In Sect. 4.1 we show a comparison of the trends in color with the major physical mechanisms expected to produce variability in the near-IR. In Sect. 4.2 we discuss the newly discovered eclipsing binaries in the region. Finally, in Sect. 4.3, we discuss the major variability characteristics of the known disk-bearing stars in our survey and use those to identify new disk-bearing stars in our sample. Finally, Sect. 5 presents a summary of our results and concluding remarks.

2. Analyzed sample

This paper is based on observations obtained in 2007 with the 3.8 m United Kingdom Infra-Red Telescope (UKIRT), at Mauna Kea, Hawaii, equipped with the Wide Field Camera (WFCAM, Casali et al. 2007), programs U/07A/H16 and U/07B/H60 (P.I. Bo Reipurth). The complete dataset is composed of 115 observed nights in the J, H, and K near-IR filters (Hewett et al. 2006) and covers 217 days. The observations cover a field of 0.8 squared degrees, centered on α2000 = 20h33m, δ2000 = +41 ° 12′. We refer to Paper I for details on data reduction, calibration, and light curve production. Errors in the measurements used in this study are between 2% and 10% and their distributions are shown in Fig. 5 of Paper I.

The sample analyzed was composed of JHK light curves of previously known candidate members of CygOB2 that lay on the FOV of the survey (see Paper I Figs. 3 and 4 for details). The list of known candidate members was built based on previous membership evaluation available in the literature, which includes disk-bearing stars identified by Guarcello et al. (2013, hereafter GDW13) and candidate members identified based on their X-ray properties as part of the Chandra Cygnus OB2 Legacy Survey (Wright et al. 2014; Guarcello et al. 2015, hereafter GDW15). Out of 5083 objects analyzed, 7% (354) were in both X-ray sources, and disk-bearing candidate member lists, 18% (918) were only in the disk-bearing list, and 75% (3811) were only in the X-ray sources list. For this paper, we divided the analyzed sample in two: the sample of 1272 disk-bearing stars and the sample of 3811 members from X-ray not listed as disk-bearing stars (sample of diskless stars). Using GTC/OSIRIS riz photometry from Guarcello et al. (2012), we inferred in Paper I that the sample has masses mainly between 0.1 and 1.0 M.

Figure 1 shows a near-IR color-color diagram for the analyzed samples. The left and middle panels show median colors for each star in the sample. The right panel shows individual data-points for the star CygOB2-093173, the most variable star in our sample. Three regions are marked (following Rice et al. 2012; Itoh et al. 1996): region “P” is typically occupied by stars with normal photospheric colors affected only by the interstellar extinction (Class III stars) and by reddened disk-bearing stars that show only a weak near-IR excess (Class II stars with only weak near-IR-excess); region “D” is typically occupied by stars with near-IR excess characteristic of disk-bearing stars (Class II stars); and region “E” is typically occupied by stars with extreme near-IR-excess possibly due to the emission of a circumstellar envelope (Class I stars). The reddening vector was built using the interstellar extinction law for the near-IR bands estimated by Hanson (2003) for the massive stars in CygOB2, transformed from 2MASS to UKIRT photometric system using the transformations proposed by Hewett et al. (2006).

thumbnail Fig. 1.

Near-IR (JHK) color-color diagrams. Left and middle panels: median colors of the 3811 diskless stars (blue dots) and the 1272 disk-bearing stars (red dots) analyzed in the present study. Right panel: CygOB2-093173 variable star. Black dots show the color of the star in each observed night. The regions labeled as “P”, “D”, and “E” are the regions typically occupied by stars with normal photospheric colors, disk-bearing colors, and extreme near-IR excess respectively. The dark red arrow shows the effect of the interstellar reddening for the median visual extinction expected for stars in the CygOB2 association (AV = 4.1m). The gray solid curve shows an empirical dwarf sequence by Kraus & Hillenbrand (2007) and the orange line shows the Meyer et al. (1997) T Tauri locus. Dashed gray lines show extinction paths of AV = 10m for 0.6 M and 0.08 M stars, and the truncation point of the T Tauri locus.

The central panel in Fig. 1 shows a significant number of disk-bearing stars identified by GDW13, but located in the “P” region. The shorter wavelength in which a circumstellar disk produces IR-excess depends on the temperature distribution in the disk, and it may be the case that some disk-bearing stars do not present such excess in the near-IR but only at longer wavelengths. As GDW13 uses near- to far-IR to select disk-bearing stars, disk-bearing stars in the “P” region are likely stars that show an excess in the mid- and far-IR, but not in the near-IR.

The left panel in Fig. 1 also shows a small number of diskless stars falling in the regions “D” and “E”. These sources are identified as young by GDW15 and have near-IR excess typical of disk-bearing stars, but were missed by GDW13 presumably due to the different depths of the two surveys. Although GDW13 and GDW15 have similar completeness of about 0.9 M, the maximum depth of GDW13 is the one of the “Spitzer Legacy Survey of the Cygnus X region” (Beerer et al. 2010), ∼0.5 M, while GDW15 go as deep as ∼0.1 M. The possibility that these sources are indeed disk-bearing stars is discussed in Sect. 4.3 in the context of near-IR variability as a disk-diagnosis.

GDW13 classifies the variable star CygOB2-093173 shown in the right panel of Fig. 1 as Class II. However, because of the variability mechanism in action, the star describes a trajectory inside the color-color diagram, with a variability amplitude of ∼1 mag in the H − K color. It spends most of its observed time in the “D” region, with a near-IR-excess compatible with a Class II classification, and part of its observed time in the “E” region, with a near-IR-excess compatible with a Class I classification. Because of the ambiguities that the variability mechanism acting in disk-bearing stars can bring to the derivation of their evolutionary status based on single epoch observations, we do not adopt GDW13 evolutionary classes (Class I, Class II, or flat-spectrum) in our analysis. Instead, we considered that every star listed by them as a disk-bearing star can present variability caused by physical phenomena related to the disk or the star-disk interaction.

3. Variability analysis

3.1. Variable star selection

As in Paper I, the initial selection of variable stars was done using the Stetson variability index (Stetson 1996). The Stetson variability index (S) attributes S ∼ 0 to uncorrelated nonvariable stars and S ≥ 1 to significantly variable stars. Figure 2 shows S as a function of median magnitudes in the H filter. The left, middle, and right plots show the full sample of candidate members, the disk-bearing stars, and the diskless stars, respectively. The mean, standard deviation and median value of S for all the candidate members are S(μ, σ, ν) = 0.9, 1.1, 0.5. The bulk of the S distribution of disk-bearing stars is shifted toward larger values compared with the distribution of diskless stars. The former has S(μ, σ, ν) = 1.8, 1.7, 1.3, and the latter has S(μ, σ, ν) = 0.6, 0.6, 0.4. The difference between the diskless and disk-bearing samples is supported by a KS-test, resulting in a close to zero probability that the two distributions come from the same parent distribution, and corroborating the idea that disk-bearing stars are more variable in the near-IR than diskless stars.

thumbnail Fig. 2.

Distribution of Stetson variability index as a function of H median magnitude for the full sample of candidate members (left), disk-bearing stars (middle), and diskless stars (right). The red line shows S = 1, and the dashed black line shows S = 0.25.

Different authors adopt different Stetson variability index limits to select variable stars. For example, CHS01 adopts S ≥ 0.55. Because in Paper I we were interested in selecting low-amplitude periodic stars, we chose the limit S = 0.25 for initial variable selection, which selects 4079 stars (80% of the full sample). In the present study we identified as significant variables stars with S >  1 (sample S1, 1280 sources). We also identified variables with 0.25 ≤ S <  1.0 if they presented at least one of the following characteristics. First if their light curves presented significant periodicity and they were among the 1154 sources discussed in Sect. 3.3.1 and Paper I. Second if their light curves presented occultation features: five sources discussed in Sect. 3.3.2. Third if their light curves presented χ2 ≥ 4.5, 4.5, 5.71 for the J, H and K filters respectively: 90 sources discussed in Sect. 3.3.3.

This procedure accounts for the selection of 1280 variable stars in the sample S1, 1249 variable stars in the sample S2, and the rejection of 1547 stars with 0.25 ≤ S <  1.0. For future reference, we show in Table 1 the general characteristics (coordinates, magnitudes, and colors) of the 2554 stars not considered as variable stars within this study. Table 2 presents the general characteristics for the 2529 sources considered as variable stars in our study. Table 3 presents the variability flags assigned and described in the remaining of this section along with Stetson variability index and, for periodic stars, the measured period. In the rest of the paper, we only analyze the stars that we considered as variable stars.

Table 1.

General information for the nonvariable stars identified within this survey.

Table 2.

General information for the Variable stars identified and studied in the present paper.

Table 3.

Information about variability classification for the stars studied in the present paper.

3.2. Variability amplitude

Peak to peak (ptp) and root main square (rms) amplitudes were used to characterize the variability amplitudes in magnitude and color, and the values measured for each variable star are presented in Table 4. The differences between the photometric variability in disk-bearing and diskless stars are evident from a simple comparison between typical values for the two samples, and we verified that disk-bearing stars amplitudes are typically higher than diskless stars amplitudes. For example, the maximum, mean, and standard deviation of ptp amplitudes in the J band were 1.79, 0.27, and 0.22 mag for disk-bearing stars, and 0.94, 0.10, and 0.07 mag for diskless stars. For the H − K color, the maximum, mean, and standard deviation of ptp amplitudes were 0.97, 0.14, and 0.10 mag for disk-bearing stars and 0.41, 0.06, and 0.04 for diskless stars. Additionally, Fig. 3 shows that there is a steep increase in the disk-fraction for larger amplitudes and that nearly all the stars with an amplitude larger than about 0.5 mag are disk-bearing stars.

thumbnail Fig. 3.

Disk fraction as a function of peak to peak amplitude for each near-IR band (left) and color (right). In both plots, circles show the center of each ptp amplitude bin and horizontal bars indicate the size of each bin. Vertical bars show the standard errors of a Poisson counting in that bin.

Table 4.

Peak to peak and root mean square amplitudes in each observed band and color for each of the variable stars studied in the present paper.

3.3. Light curve morphology

We visually inspected the light and color curves and classified the variable stars according to their morphology as shown in Fig. 4. The example in Fig. 4 shows the disk-bearing star CygOB2-093173, whose color-color diagram presented in Fig. 1. We describe in the following sections the details of different types of variability morphology present in our sample.

thumbnail Fig. 4.

Star CygOB2-093173 is a disk-bearing (Class II) star and an example of non periodic variable star, with S = 15.64, and that presents reddening linear trajectories in the color-color and color-magnitude diagrams. Left, from top to bottom: J light curve, a color bar showing the color attributed for each observation date, K light curve, and H − K color curve. The dotted lines show maximum and minimum values in each light curves, and dashed lines show the median values. Right, from top to bottom: J vs. J − H color-magnitude diagram, K vs. H − K color-magnitude diagram, and color-color diagram with P, D, and E regions (discussed in Sect. 2). The dark red arrows show the effect of an interstellar reddening of AV = 2m. The vertical and horizontal magenta line segments show the median error in each axis. In each plot, error bars are shown in gray. When visible and close to the position of the star, a gray curve shows a Kraus & Hillenbrand (2007) empirical isochrone in the color-color diagram, and the orange line shows the T Tauri Locus. Each diagram has a range of 5σ around the median color or magnitude for that star.

3.3.1. Periodic variability

In Paper I, we identified stars with periodic variability produced by the rotational modulation of a spotted stellar surface. As an initial step, the first morphological class defined was of candidate periodic stars. We considered a total of 1679 stars as candidate periodic because their light curves were visibly oscillating between max and min values, or because they presented significant peaks in their Lomb-Scargle periodogram (Lomb 1976; Scargle 1982), in each photometric-band. In Paper I, we measured periods for 1256 of these stars (75% of the candidate periodic sample), but due to varied contamination issues introducing ambiguities in the period determination process – such as limitations due to the time sampling – only 894 had reliable measured periods. Finally, periodicity seems more common among low amplitude variable stars, and only 31% of the candidate periodic stars were variables of type S1.

Figures 5 and 6 present two examples of periodic stars with their light curves folded to the periods obtained in Paper I. The star CygOB2-000348 (Fig. 5) is diskless and periodic. The star CygOB2-113855 (Fig. 6) is disk-bearing and periodic in the light curves as well as in the color-curves. Only 39 stars (∼2% of the candidate periodic sample) presented periodicity in their phase folded color curves. This lower incidence of periodicity in color does not mean that the physical mechanisms causing periodic variability are colorless, but their amplitudes of variability in the near-IR colors are often below our limits of detection.

thumbnail Fig. 5.

Star CygOB2-000348 is an example of a diskless periodic variable, with S = 1.62, a ptp amplitude of 0.13 and 0.11 mag in the J and K bands respectively, and that does not show periodicity in color. The plots, colors, and symbols are the same as in Fig. 4, but with the light curves and color curve phased with the period of the star, P = 5.1 d.

thumbnail Fig. 6.

Star CygOB2-113855 is an example of a disk-bearing (Class II) periodic star. With S = 2.77 and ptp amplitude of 0.31 and 0.36 mag in the J and K bands, the star is also periodic in color. It is also an example of a compound variable with a variable IR-excess. The plots, colors, and symbols are the same as in Fig. 4, but with the light curves and color-curve phased with the period of the star, P = 10.3 d.

3.3.2. Occultation variability

We based our definition of occultation variability on the existence of asymmetry in the light curve: these were stars that spent most of their observed time around the median value of magnitude but presented drops in brightness that can be periodic or not. Based on the visual inspection of the light curves, we identified 125 stars with occultation variability: 15 presented only a single or a couple of dips in their light curves; 110 presented multiple dips and we previously examined their periodicity in the context of Paper I.

Figures 7 and 8 show two examples of occultation variables. The first is the disk-bearing star CygOB2-124656, which presented only two occultation features in the first season of observations. The second is the diskless star CygOB2-041510, which has periodic occultation features and is discussed in Sect. 4.2 as an eclipsing binary candidate.

thumbnail Fig. 7.

Disk-bearing (Class II) star CygOB2-124656 is an example of an occultation variable that presents reddening linear trajectories in the color-color and color-magnitude diagrams. The star has S = 3.59 and ptp amplitude of 0.43 and 0.17 in the J and K bands. The plots, colors, and symbols are the same as in Fig. 4.

thumbnail Fig. 8.

Diskless star CygOB2-041510 is an example of a periodic occultation variable considered as an eclipsing binary candidate. It is also an example of a compound variable with occultation data-points with distinct color behavior. The star has P = 4.68 days, S = 3.06, and ptp amplitude of 0.28 in both J and K bands. The plots, colors, and symbols are the same as in Fig. 4.

In Paper I, we used Monte Carlo simulations to investigate the limitations of our data-sampling on detecting periodic stars with different waveforms. We concluded that our dataset was not appropriated for measuring reliable periods in periodic occultation wave-forms and that most of the periods measured were harmonics of the correct period. We, therefore, decided not to try to distinguish between eclipsing binary variables and other types of variables with occultation features (like dippers) based on their periodicity and level of asymmetry. Section 4.2 presents a further discussion about the identification of eclipsing binary candidates in our sample.

3.3.3. Other types of variability

A total of 726 stars were identified as variables but presented neither periodicity nor occultation features. Based on the visual inspection of the light curves, we verified the existence of 90 long-term variables (timescales longer than a month) as, for example, the star CygOB2-093173 in Fig. 4. We also identified 128 mixed-variables (variability with different time-scales) and 508 variables that were neither long nor mixed and were therefore called short-term variables (time-scales from days to weeks), for example, the star CygOB2-113911 in Fig. 9.

thumbnail Fig. 9.

Disk-bearing (flat-spectrum) star CygOB2-113911 is an example of a long-term variable with S = 3.67, that presents a blueing linear trajectory in the color-magnitude diagrams. The plots, colors, and symbols are the same as in Fig. 4.

3.3.4. light curve asymmetry

We applied the asymmetry metric (M-index) as defined by Cody et al. (2014) to statistically asses the level of asymmetry in the light curves. We defined a threshold for selecting asymmetric light curves by estimating the values of M-index in a large sample of sinewaves with the same time-sampling and uncertainty of our survey and verifying that all sinusoidal light curves were in the interval −0.5 ≤ M ≤ +0.5.

We then used M ≥ +0.5 to identify and flag 115 variables with dimming asymmetries. Of those, 67 were visually identified as occultation variables, including 13 of the 14 stars discussed as eclipsing binaries in Sect. 4.2. Two examples are the stars CygOB2-124656 (Fig. 7) and CygOB2-041510 (Fig. 8). This criterion also flagged 29 periodic variables and 19 stars with other types of variability. We used M ≤ −0.5 to select 11 stars with burst asymmetries. The stars CygOB2-093173 (Fig. 4) and CygOB2-113911 (Fig. 9) are examples of it. All the stars with burst asymmetry were in the samples of other types of variability.

3.4. Variability in color-space

We mentioned in Sect. 3.3.1 that around 2% of the candidate periodic sample also showed periodicity in their folded color curves. In this section, we identify other behaviors in the color space that could be used to infer the physical processes responsible for the variability.

3.4.1. Correlations between color and magnitude variations

Variable stars may show strong correlations between their variability in different bands, which causes the star to move along linear trajectories in the color space. We investigated these correlations by measuring the slopes of the linear trajectories in the color-magnitude diagrams and and in the color-color diagram .

As weak correlations may be biased by uncertainties in the data, before estimating the slopes, it was mandatory to first investigate how strong was the correlation between the two variables in question. The Spearman nonparametric rank correlation test (Press et al. 1992) was used to evaluate the statistical significance of the correlations in the color-color and color-magnitude diagrams. Spearman correlation measures how well the relationship between two variables can be described by a monotonic function. It can assume values between −1 and 1 with the extremes occurring when the relationship between the variables is perfectly monotonic. For stars with the absolute value of Spearman correlation larger than 0.7, the orthogonal distance regression method2 was used to estimate the slopes of the trajectories in the color-color and color-magnitude diagrams and their uncertainties. This method takes the uncertainty of each data point into account while performing a linear fit. Finally, we overplotted the slopes measured along with the real data in each of the color spaces, and visually inspected them in order to remove stars that were not well described by a linear relation.

Figure 10 and Table 5 present the slopes for the 335 variable stars with trajectories in color space that could be well fitted by a linear relation. The majority of those stars showed higher ptp amplitudes than the median values for the whole sample of variable stars, or belonged to the S1 sample. This suggests that our analysis of the correlations between color and magnitude variability is biased toward high amplitude variable stars.

thumbnail Fig. 10.

Color-magnitude and color-color diagrams with the slopes for the stars presenting a linear variation. Red lines show disk-bearing stars, and blue lines show stars without discs. The length of each line segment reflects the amplitude of the variability. The dark red arrows show the median interstellar extinction toward CygOB2, and the gray lines show as reference of unreddened normal photosphere colors, an empiric isochrone from Kraus & Hillenbrand (2007).

Table 5.

Parameters for the stars with measured slopes in the color-color diagram or color-magnitude diagrams.

Only 50, out of 335 stars with linear trajectories were diskless. Among the 144 stars with measured slopes, 29 had positive slopes and 115 had negative slopes. Among the 196 with measured slopes, 193 had positive slopes, and 3 had negative slopes. Finally, 143 stars had measured slopes and all of them positive.

The stars CygOB2-124656 (Fig. 7), CygOB2-093173 (Fig. 4) and CygOB2-113911 (Fig. 9) are examples of stars presenting linear trajectories in the color-diagram. We further discuss in Sect. 4.1 how the slopes of their trajectories can help distinguishing between the different possible physical mechanisms causing variability.

In the discussion that follows, it is crucial to understand how the slopes are defined. In the case of color-magnitude diagrams (left and middle plots in Fig. 10), a positive slope is defined as a positive change in color, and a positive change in magnitude, which results in a vector pointed toward the bottom-right of the color-magnitude diagram due to the inverted scale of magnitude, and the slope grows clockwise. In the case of the color-color diagram (right plot in Fig. 10), a positive slope is defined as positive changes in both colors, the resultant vector points to the top right, and the slope grows counterclockwise.

3.4.2. Compound variables

Often the variability behavior was not confined into a single monotonic trajectory with a positive or negative slope in the color-space. For example, the disk-bearing periodic star CygOB2-113855 (Fig. 6) has an approximately fixed amplitude in its periodic variability, but presents a decreasing flux in K responsible for a varying near-IR excess that shifted the star across the K vs. H − K and color-color diagrams. Such a complex color behavior may indicate multiple variability mechanism in action in the same star. RRW15 describes similar behaviors as a “compound color behavior”. We verified three types of compound color behavior in our sample:

Asymmetry with a distinct color. The occultation variable CygOB2-041510 (Fig. 8) presents occultation datapoints with distinct color behavior than the rest of the stellar points. Similarly, 20% (25) of the occultation variables presents similar behavior, which is associated with eclipsing binary variability in Sect. 4.2.

We also verified this behavior in variable stars with burst asymmetry in the light curve. For example, the star CygOB2-041756 (Fig. 11) shows a mixture of variability modes. Overall, the star is dimming while describing a linear trajectory with a negative slope in the CMD. However, in the first 50 epochs of observation, the star has several burst-features with varying amplitudes up to one magnitude in both J and K, showing a distinct behavior in the color space.

thumbnail Fig. 11.

Disk-bearing (Class II) star CygOB2-041756 (S = 10.4) is an example of a mixed variable with a long term variability mixed with burst features in the first epochs observed. The star is a compound variable with asymmetries in the light curve having a distinct color behavior. The plots, colors, and symbols are the same as in Fig. 4.

Variable near-IR excess. Some stars present an increasing or decreasing near-IR excess, causing the data points observed in different seasons to be clumped in different regions of the color-color and CMDs. The periodic star CygOB2-113855 (Fig. 6) is an example of it: while the amplitudes of periodic variability are unchanged during the observed time, and no significant change is observed in the J band, the median magnitude in the K band diminishes shifting from the “P” to the “D” region of the color-color diagram. We note that whatever is the physical mechanism causing the decrease of brightness in the K band, it seems to be independent of the dominant physical mechanism producing the periodic variability. This type of variability can be associated with what Rice et al. (2012) called “transient near-IR excesses”, which is related to a variability mechanism that changes the amount of near-IR excess and some times makes the star transition between the “P”, “D”, and “E” regions of the color-color diagram.

Changing slope. Some stars show an evident change in the dominant variability mechanism, marked by a change in the slope of variability in the color-diagram. For example, the star CygOB2-021188 (Fig. 12) gets bluer when fainter during the first ∼100 observed nights, and it presents a linear trajectory with a negative slope in the K vs. H − K CMD during this time. However, during the final ∼100 nights of observations, the star presents a shorter term variability in which it gets redder when fainter, and it presents a linear trajectory with a positive slope in the J vs. J − H CMD during this time. Note that the star also shows a variable near-IR excess. As a consequence, it transitions between the “D” and “E” parts of the color-color diagram during the observed time.

thumbnail Fig. 12.

Disk-bearing (Class II) star CygOB2-021188 (S = 1.0) is an example of compound variability showing a variable near-IR excess and a changing slope during the observed time. The plots, colors, and symbols are the same as in Fig. 4.

We identified 149 stars with one or multiple modes of compound variability. Figure 13 shows more examples of each case. As in star CygOB2-021188 (Fig. 12), often variable near-IR excess and changing slope happened at the same time.

thumbnail Fig. 13.

Examples of stars with a compound behavior are shown inside color-color diagrams or color-magnitude diagrams. 1- Asymmetry with a distinct color; 2- Variable IR-excess; 3- Changing slope; 5- Mix of 2 and 3. Dark red arrows show the effect of interstellar reddening for AV = 1 mag. The color of each data point reflects its observational date following the color-bar in Fig. 4.

4. Results and discussion

The results obtained for the 2529 variable stars studied in this paper are available as electronic tables at CDS. Table 2 presented general information on the light curves used in this study, along with literature comparison for each variable star on our list. Table 3 presented the information related to the variability classification. Table 6 summarizes how we categorized disk-bearing and diskless stars as discussed in the previous sections.

Table 6.

Summary of the variability classification performed in Sect. 3 for both disk-bearing and diskless samples.

4.1. Variability trajectories in the color-space and the physics behind photometric variability of young stars

In Sect. 3.4.1 we estimated the slopes of the linear trajectories described by 355 variable stars inside the color diagrams in Fig. 10. In this section, we put such slopes in context by comparing them to reference values collected in the literature and discussing their relation with the physical mechanisms responsible for variability in the near-IR. Figure 14 shows the distributions of , , and slopes measured in Sect. 3.4.1. When examining Fig. 14, we stress that the distributions of slopes in each color space do not necessarily correspond to the same stars. Out of 355 stars that had at least one slope measured: 20 stars had all three slopes measured; 12 had both and , but no ; 88 had both and , but no ; nine had both and , but no ; 24 had only ; 79 had only ; 103 had only .

thumbnail Fig. 14.

Distribution of slopes measured in the color-color and color-magnitude diagrams. The figure shows slopes for J − H vs. H − K (top), J vs. J − H (middle), and K vs. H − K (bottom) diagrams. We plotted the expected slopes for each physical mechanism as line segments, following the ranges in Table 7. The X mark shows the exact position of the T Tauri Locus. The black circle shows the slope expected for the typical dust opacity index in the ISM. Black bars show the stars in S1 sample and with precision in the light curve better than 5%, empty bars show stars in the S2 sample or with precision worse than 5%.

4.1.1. Positive slopes in the CMD: spot and extinction variability

Stellar spots and extinction variability can cause positive slopes in the near-IR CMD. Stellar spots can be hotter or colder than the stellar photosphere. They produce variability because the fraction of the stellar surface covered by the spots and seen by the observer changes as the star rotates. Cold spots are similar to sunspots (e.g., Bouvier et al. 1993), they are colder than the photosphere and arise from magnetic activity. In young stars, they can cover up to ∼60% of the stellar surface (Venuti et al. 2015). They can be very stable and persist over thousands of rotation periods (Herbst et al. 1994). Cold spots are likely the most common source of photometric variability in low-mass stars at all ages. Indeed, in Paper I we found 1119 diskless stars with periodic variability likely associated to cold spots, which accounts for 65% of the diskless variables in our study.

Hot spots are only possible in disk-bearing accreting stars as they form in the regions where the accretion streams hit the stellar surface. Stable accretion streams are generated by a misalignment between the rotation and magnetic axis, and the hot spots are typically located at high latitude (Kulkarni & Romanova 2008; Romanova et al. 2013). Unlike cold spots, they only persist over a few rotation periods. In case of unstable accretion, hot spots may be located at low latitudes, and the accretion bursts occur stochastically and last only a few hours, typically less than a stellar rotational period (Romanova et al. 2012; Kurosawa & Romanova 2013; Stauffer et al. 2014; Colombo et al. 2016).

With cold spots, the star gets redder and fainter when the spot is visible. With hot spots, the star gets bluer and brighter when the spots are visible. In both cases, the resultant trajectory in the CMD will have a positive slope. The amplitude and slope of the trajectories in the color-space depend on the difference of temperature between the spot and the star, and also on the size of the spot, but cold and hot spots are expected to produce a similar range of slopes in the CMDs. On the color-color diagrams, cold spots usually describe shallower trajectories than hot spots.

In Appendix A.1, we used a spot model to infer typical near-IR amplitudes and slopes for spot variability. We fed the model with values of spot temperatures in YSOs inferred from optical and UV observations in the similarly aged cluster NGC 2264 (Venuti et al. 2015), along with the effective temperature of the same stars, and obtained expected amplitudes in the JHK bands. We present the range of slopes found in Table 7, along with literature values.

Table 7.

Slopes of variability in the near-IR expected from models and other observational studies.

Extinction variability can arise from inhomogeneities in the absorbing material moving across the line of sight. These inhomogeneities can be in the inner circumstellar environment or in the ambient molecular cloud. For example, in AA Tauri variable stars (e.g., Bouvier et al. 2007), the magnetospheric star-disk interaction, via which accretion occurs, can lift some dust above the plane of the disk, creating a warp in the inner disk. As the dust is optically thick, this warp causes extinction occultation of the stellar photosphere.

As with spots, extinction variability also makes the star redder as it gets fainter. The slopes of the trajectories described by the cold spots, hot spots and extinction can be very similar, making it difficult to distinguish between the three phenomena in the color-space. Our estimations of typical extinction slopes are presented in Table 7 and detailed in Appendix A.2. In the remaining discussion, we call stars describing linear trajectories with position slopes “reddening variables”.

For example, the disk-bearing star CygOB2-093173 (Fig. 4) presents a reddening variability with a 42.9° slope in the J − H vs. H − K diagram, 74° in the J vs. J − H CMD and 55° in the K vs. H − K CMD. The star has high ptp amplitudes (1.79 in J and 0.77 in K), and has a borderline compound behavior in the K vs. H − K diagram, where the positive slope of its trajectory became shallower after the 80 first nights of observation, suggesting a mix of spot variability and extinction.

The disk-bearing occultation variable star CygOB2-124656 (Fig. 7) is another example of reddening variability. Over the course of the occultation feature (first ∼80 nights of observation), the star presents a reddening slope, which is more remarkable in the J vs. J − H diagram, with a 72° inclination. In the K vs. H − K diagram, the star does not present a definite linear trajectory in the second season of the observations when occultation features are no longer observed. The reddening variability correlated to the occultation features, added to the fact the star has a disk, suggest the variability is caused by dust and the star is likely an AA Tauri variable.

4.1.2. Negative slopes in the CMD: disk emission and accretion variability

Disks can contribute to the IR emission of YSOs by absorbing the optical and ultraviolet radiation from the star and reradiating it at longer wavelengths, and by releasing energy as the material in the disk is radially transported through a viscous optically thick accretion disk (Lynden-Bell & Pringle 1974; Meyer et al. 1997). The resultant IR-excess can occur in different wavelengths from near-IR to far-IR. The region in the star-disk system producing an excess in each wavelength range depends on the temperature distribution in the star-disk system. Typically, the near-IR excess originates in the hotter inner parts of the disk, and the shorter wavelength at which the IR-excess is observed depends both on the temperature of the inner disk material and on the inner hole size. An increasing accretion rate, for example, changes the release of energy and temperature distribution in the inner disk which alters the amount of absorbed and reprocessed radiation (Scholz et al. 2009). As the temperatures in the inner disk become larger than the stellar surface temperature, the flux from the disk overcomes the stellar one and the star gets brighter in the near-IR and redder. On the other hand, holes in the inner disk alters the amount of absorbed and reprocessed radiation, diminishing the disk contribution to the overall near-IR flux. As the disk emission becomes stronger at longer wavelengths, one can expect higher amplitude variations in the K band, than in the J band, while the variations make the star bluer (approaching the true stellar color) as it fades.

In both cases, variability related to the disk emission will produce negative (“blueing”) slopes inside the color-magnitude diagrams and positive slopes inside the color-color diagram. Table 7 shows literature values for the slopes produced by these type of physical mechanisms. CHS01 reports amplitudes as large as 1 mag in K band for this type of variability.

The disk-bearing star CygOB2-113911 (Fig. 9) is an example of blueing variable. The star described linear trajectories in all diagrams with a 48.72° slope in the J − H vs. H − K diagram, −74° in the K vs. H − K diagram, and −49° in the J vs. J − H diagram. The datapoints dominating the linear trajectory are related mainly to two burst-features in the light curve: a lower amplitude one between nights 120 and 160, and a very intense one after night 170. The behavior of these burst features suggests that its variability is caused by variable accretion episodes.

4.1.3. Distributions of slopes

Figure 14 showed the distribution of slopes of the linear trajectories described by 355 variable stars inside the color diagrams. The distributions can be compared to the ones for the ONC in CHS01 (see their Fig. 20) and RRW15 (see their Fig. 20). To facilitate the comparison, we adopted the same terms as these authors, and directly compare reddening variables (positive slopes in the CMDs related to spot or extinction variability) with blueing variables (negative slopes in the CMDs related to variability caused by changes in the inner disk or variable accretion).

In the distribution of slopes in the J − H vs. H − K diagram (top panel in Fig. 14), the distinction between types of variability is more complicated, as all physical mechanisms investigated produce positive slopes inside the color-color diagram. Even though blueing variability is expected to produce slopes shallower than about 35°, there is still a significant overlap between the slopes of spotted stars and the TTauri locus. With 28 stars with shallower slopes in the range 21°– 35° and 116 in the range 35°– 58°, this diagram indicates a mixture of disk blueing variability and spot variability.

In the J vs. J − H diagram (middle panel of Fig. 14), we found mostly reddening stars, with 154 stars presenting variability with slopes between 67° and 81°, and only three stars presenting blueing variability with slopes between −49° and −58°. This indicates that this diagram is composed of a mixture of extinction and spot variability. We found 39 stars with slopes shallower than 67°, which is a larger number of shallow positive slopes than CHS01 or RRW15. However, most of those are from the S2 sample or had large photometric uncertainty between 5 and 10%.

In the K vs. H − K diagram (bottom panel of Fig. 14), we found four times more blueing variables than reddening ones, with 115 stars presenting blueing variability with slopes between −51° and −80°, and 29 stars presenting reddening variability with slopes between 55° and 74°. This distribution differs a lot from RRW15, where they found only 1.2 times more blueing variables than reddening ones. As RRW15 used S ≥ 1 for selecting variable stars, if we considered in the comparison only S1 sources with photometric uncertainty lower than 5% (black bars in Fig. 14), we still find three times more blueing variables than reddening variables.

The fact that few stars presented blueing variability in the J vs. J − H CMD is explained by the nature of the near-IR excesses in disk-bearing stars. As the amount of light re-emitted by the disk in the near-IR is larger for longer wavelengths, one may expect that the amplitudes of variability due to variable accretion rate and changes in the inner disk will also be more significant for the K magnitude and H − K color, than for J magnitude and J − H color, thus also explaining why more blueing variability is observed in the K vs. H − K diagram. The flux and variability in the J band is dominated mainly by the stellar flux, therefore explaining why the distribution of slopes in the J vs. J − H CMD is dominated by reddening related to spot variability, which acts directly in the photosphere, and extinction, which dims the light emitted by the star itself.

Considering both J vs. J − H and K vs. H − K diagrams, CHS01 found six times more reddening type variables than blueing type variables, RRW15 found 1.5 times more reddening than blueing types. Considering only stars in the S1 sample with photometric uncertainty lower than 5%, our dataset presents more similar results to RRW15 than with CHS01, as we found 2.1 times more reddening than blueing variables. We caution the reader on interpreting the slope distributions. While the study of the trajectories in color is especially useful for selecting high amplitude variable stars with a single dominant physical mechanism producing variability in their light curves, the frequency in which each physical mechanism is observed in a dataset is most likely not representative of the recurrence of such physical mechanism in nature. RRW15 remarked that the larger fraction of blueing variables that they found compared to Carpenter et al. (2001) may be because the differences in the total time-span of each survey may favor the detection of certain types of variables. For example, surveys developed over a couple of months, like Carpenter et al. (2001) and Cody et al. (2014, who surveyed NGC 2264 for variability in the optical and mid-IR) might favor the detection of variable stars such as AA Tauri, which present variability with time-scales of days to weeks. On the other hand, as blueing variability has typically long time-scales, surveys as long as Rice et al. (2015, 2.5 years), Rice et al. (2012, 1.5 years) or the present study may favor the detection of this type of variable.

We identified four major differences between our data analysis and other studies that could introduce bias to our results and explain the differences with results from the literature. First, in our survey we considered stars with photometric uncertainties up to 10%. Second, even when considering only stars with uncertainty lower than 5%, our survey goes deeper than RRW15 by 0.2 mag in the J band, and by 0.3 mag in the H and K bands. Third, a collateral effect of the definition of our S2 sample is the introduction of variables with large amplitude, but with high uncertainties. Fourth, both RRW15 and CHS01 searched for variability among ∼18 000 good quality light curves in the FOV of ONC. In our case, we only analyzed light curves extracted for previously known YSO in CygOB2. As discussed in Sect. 4.1.2, the variability mechanisms producing blueing variability are intimately related to the infrared-excess in disk-bearing stars. This means that by analyzing the light curves for disk-bearing stars identified by GDW13 we are preferentially selecting stars that have augmented IR-excess and more likely to show blueing variability.

We address the first and third points by looking only at the stars in the S1 sample with uncertainty lower than 5%, which accounts for 26 reddening stars and 84 blueing stars in the K vs. H − K diagram, still yielding 3.2 times more blueing variables. Therefore confirming our previous results. On the other hand, the second and fourth points are likely the most important ones to keep in mind when comparing our results with the literature, specially when looking at the K vs H − K CMD, which is mostly dominated by disk-variability. Rice et al. (2015) found 1.2 times more blueing variables in this diagram than reddening ones, in a sample of ∼15k stars in the FOV of ONC. We found 3 more times blueing variables in a sample of ∼5k YSOs, 1/4 of which known to have strong IR-excesses.

4.2. Eclipsing binary candidates

Another source of variability not yet discussed is the occurrence of eclipses in a binary system. We verified the occurrence of eclipsing binaries in our dataset in two ways: First, we further analyzed our list of occultation variables to identify eclipsing binary candidates based on their variability characteristics. Second, we matched our dataset with the list of eclipsing binaries in the field of view of CygOB2 previously identified by HSP11.

We considered that an occultation feature was due to an eclipsing binary if:

  • 1.

    The star was classified as occultation variability (Sect. 3.3.2);

  • 2.

    The drops in brightness were sharp, deep and homogeneous:

    • (a)

      Given the combination of the cadence of about one observation per night in each filter, and the short periods expected for eclipsing binaries3, each drop in brightness is expected to be composed typically of one or two sets of JHK observations;

    • (b)

      The eclipse happened similarly in the three bands;

    • (c)

      Their typical depth was ≳0.1 mag;

  • 3.

    The star was classified as a compound variable of the first type (Sect. 3.4.2), with all the data points that are part of the eclipse being localized outside the major bulk or pattern of observations inside the CMDs, and not following a correlated slope as in Sect. 3.4.1.

Following these three criteria, we identified 14 eclipsing binary candidates. Their properties are presented in Table 8. To sum up: the identified eclipsing binary candidates presented magnitudes between 10.45 and 16.29 mag in the J band, and ptp amplitudes between 0.06 and 0.33 mag in the same band. None of them were disk-bearing stars or had measured slopes of variability in color. Their periods ranged from 0.83 to 7.30 days.

Table 8.

Information about the eclipsing binary candidates identified in the present study.

Among the eclipsing binary candidates identified, the four stars in Fig. 15 presented recurrent eclipses over the observed time, but these eclipses were not periodic inside the observational window of our data.

thumbnail Fig. 15.

Four eclipsing binary candidates that were not periodic in our study. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are written in the plot. Right: J vs. J − H CMDs. In both left and right plots, colors reflect the date of the observation as in Fig. 4.

As mentioned in the Sect. 3.3.2, even though we could detect periodicity in some of the occultation variables, we often could not recover their correct period. This is partially due to the limitations of the classical implementation of the Lomb-Scargle periodogram, which searches for pure sinusoidal signals, and partially due to the observational cadence of our observations, which was not ideal for detecting periods below two days, as further discussed in Paper I. To complement the study of the eclipsing binary candidates by taking a closer look into their periodicity, we applied two additional period search techniques alternative to the Lomb-Scargle periodogram: the fast-χ2 method (Palmer 2009)4, which is a higher order Fourier method capable of dealing with waveforms that are more complex than a pure sinewave; and the String-Length method (Clarke 2002), which is a nonparametric period search method that does not use any assumption regarding the shape of the light curve. When these two methods lead to different periods, or when we found distinct periods in the different observed bands, we visually inspected the folded light curves and their phase dispersion and based the choice of the adopted period on it.

For the stars CygOB2-041491 and CygOB2-124487, we verified the existence of two periods in the light curve. Figure 16 shows their light curves folded to the period found for the occultation features. Figure 17 shows the folded light curves for the secondary period, and in both cases, the secondary period seems to be reflecting the rotational modulation of spots in the primary star in the eclipsing binary system. We remark that in Table 8, we presented the star CygOB2-041491 with two possible primary periods because different period-search techniques gave different periods for the occultation features. The period adopted for building the plot in Fig. 16 was the one that minimized the phase dispersion of the occultation features.

thumbnail Fig. 16.

Two multiperiodic eclipsing binary candidates. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band, folded using the detected period attributed to the occultation feature. J, H, and K bands are shown as blue, green, and red respectively. Figure 17 shows the folded light curves for the secondary period. Right: J vs. J − H CMD. Colors reflect the date of the observation, as in Fig. 4.

thumbnail Fig. 17.

Phased light curve for each band, folded using the secondary detected period found for the eclipsing binary candidates presented in Fig. 16. J, H, and K bands are shown as blue, green, and red respectively.

Figure 18 shows complementary plots for the eclipsing binary candidate CygOB2-041510, already presented in Fig. 8, and which has a period of 4.69 days. Seven other eclipsing binary candidates also had a single period measured and five of them are presented in Fig. 19 and two of them in Fig. 20.

thumbnail Fig. 18.

Eclipsing binary candidate CygOB2-041510 (Fig. 8). Left: J, H, and K light curves for the star. The median magnitude and ptp amplitude are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Right: phased light curve for each band, folded using the detected period attributed to the occultation feature. J, H, and K bands are shown as blue, green, and red respectively.

thumbnail Fig. 19.

Five periodic eclipsing binary candidates. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band, folded using the detected period attributed to the occultation feature. J, H, and K bands are shown as blue, green, and red respectively. Right: J vs. J − H CMD. Colors reflect the date of the observation, as in Fig. 4.

thumbnail Fig. 20.

Stars CygOB2-000644 and CygOB2-072104 are both occultation and compound variables and were considered as eclipsing binary candidates in the present study. They are also reported as eclipsing binaries in HSP11. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are also written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band. J, H, and K bands are shown as blue, green, and red respectively. Right: J vs. J − H CMD for the same star. Colors also reflect the date of the observation, as in Fig. 4.

Finally, we also took into consideration the list of eclipsing binaries reported by HSP11. Among HSP11 sources, we recovered in our sample four out of ten eclipsing binaries and 13 out of 17 candidates. One of their eclipsing binaries, the star CygOB2-031292, was not considered a variable within the limits of our data. Two of their eclipsing binaries were also considered as eclipsing binaries by us, and they are shown in Fig. 20. We found a 3.21 day period for CygOB-000644 based on the classical Lomb-Scargle periodogram and on the string-length method, but a 1.447 day-period based on the fast-χ2 method. The shorter period is a harmonic of the 2.89 days period found by HSP11 for the same star, which we considered as the correct period since it looks better in terms of the phase dispersion in the folded light curve. In the case of the star CygOB2-072104, we found the same 0.83 days period as in HSP11. The forth HSP11 eclipsing binary in our data was the star CygOB2-114152, shown in Fig. 21. Even though this star presented a compound behavior in the CMDs, similar to other identified eclipsing binaries, it was identified as a periodic star, but it was eliminated from the list of stars with a measured rotation period in Paper I because it presented a different period from the one in HSP11. We found a period of 3.63 days based on the Lomb-Scargle periodogram and a 0.78 days period based on the fast-χ2 method. Both values are aliases of the 0.42330 days period reported by HSP11 for the star, but the correct period is below the limit of detection of our observations. Out of the 13 eclipsing binary candidates from HSP11, in our survey two were not considered as variable within the limits of our dataset, and the 11 remaining were reported as periodic variables, but none of them presented other characteristics indicating they could be eclipsing binaries. The differences in classification are probably due to the different observational windows in the two surveys, and to account for the possibility those are indeed eclipsing binaries, in Paper I we removed those stars from our sample of spotted stars.

thumbnail Fig. 21.

Star CygOB2-114152 is a periodic candidate and compound in this survey but it is reported as an eclipsing binary in HSP11. Left: J, H, and K light curve. The median magnitude and ptp amplitudes are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band. J, H, and K bands are shown as blue, green, and red respectively. The period reported by HSP11 (below our detection limit) was adopted in the phased light curve. Right: J vs. J − H CMD. Colors reflect the date of the observation, as in Fig. 4.

4.3. Near-infrared variability as a disk indicator

Finally, we address the reliability of using near-IR variability as a disk diagnosis. We start by looking at what are the typical variability behavior of the disk-bearing stars in our survey, and if we could use any of their typical behavior as a disk diagnosis.

Variability amplitude. Figure 3 shows that the disk fraction sharply increases with variability amplitude, suggesting that large amplitudes might indicate the presence of a disk. In our sample, 196 stars have ptp amplitude larger than 0.5 mag in at least one of the near-IR bands, 90% of which are disk-bearing. Out of the 12 diskless stars with such large amplitudes, 8 have slopes measured in Sect. 3.4.1. Their slopes are larger than ≳70° in the J vs. J − H CMD, or shallower than ≲35° in the color-color diagram, implying that their variability is either disk-variability or extinction by circumstellar material, and hinting they are disk-bearing stars.

Slopes in the color-space. The measurement of the slopes of the trajectories in the color diagrams is a promising disk indicator as 84.8% of all slopes we measured are in known disk-bearing stars. As discussed in Sect. 4.1, only disk variability can produce negative slopes, thus measuring those can lead to identifying disk-bearing stars. We found 115 stars with negative slopes: 110 known disk-bearing stars, and 10 newly discovered disk-bearing sources. Positive slopes can also be an indirect disk-indicator, as among their causes there are hots pots produced in the accretion process and circumstellar extinction. However, that variability caused by cold spots is independent of the presence of a disk and can also explain trajectories with positive slopes having a range of slopes that is degenerate with the ones related to disk variability. Therefore, positive slopes do not immediately imply the existence of a disk. The overlapping of the expected range of slopes for cold spots and other disk related variability is almost complete in the J vs. J − H CMD, but only partial in the K vs. H − K CMD. Therefore, it is still be possible to differentiate between disk-variability and cool spots given the positive slopes they produce in the K vs. H − K CMD, as the former causes shallower slopes than the latter (see Table 7). We found 196 stars with positive slopes in the CMDs, 157 of which are known disk-bearing stars. Unfortunately, all the 39 diskless stars with a positive slope in our sample had only measured, and further information is required to confirm if they could indeed be disk-bearing candidates.

Compound variability. In Sect. 3.4.2 we associated asymmetries with a distinct color in occultation variables with eclipsing binary variability, and those do not suggest the presence of a disk. Variable near-IR excess indicates the presence of a disk and accordingly the 24 compound variables of this type in our sample were known disk-bearing stars. Once eclipsing binary variability is ruled out, the only variability mechanism expected in the diskless stars is due to cold spots. As we associated compound stars with a changing slope with two or more alternating variability mechanisms, the 24 compound variable of this type must host a disk. Indeed, 21 of them were known disk-bearing stars and thus the remaining 3 changing slope variables are newly discovered disk-bearing (for example, the star CyOB2-000111 shown in Fig. 13). We also found 59 compound stars showing both variable near-IR excess and changing slopes (see, for example, the star CygOB2-021188 in Fig. 12), 3 of which are newly discovered disk-bearing stars.

Light curve morphology. Periodic variability is not a good indicator of the presence of a disk as it is the most common type of variability among the diskless stars (only 10.2% of the periodic stars had a known disk). For occultation variables, once eclipsing binary candidates were identified and removed, ∼85% of the remaining stars (83 stars) are known disk-bearing. Those are occultation variables that can be associated with circumstellar extinction (“dippers”), and therefore also related to a disk-presence. On the other hand, 15 are diskless: four are considered disk-bearing due to their changing slope compound variability, two are considered candidates due to their positive slopes of variability in the J vs. J − H CMD. The remaining nine sources are also considered disk-bearing candidates because of their light curve morphology.

Finally, 560 (77% of them) of the other types of variables are known disk-bearing stars. Among the 166 that are diskless: 10 are identified as new disk-bearing stars due to their negative slopes in the CMD and two due to their compound variability; 28 are identified as disk-bearing candidates due to their positive slopes in the CMDs. We thus consider that ruling out both eclipsing binary candidates and periodic variables from a list of variable stars seem to be a valuable first step toward selecting disk-bearing stars.

Considering the variability characteristics described in this section, we verified that 681 previously known disk-bearing stars in CygOB2 can be identified as disk-bearing stars based solely on their near-IR variability. We found 18 newly discovered disk-bearing stars among the list of members of Cygnus OB2, based on their negative slopes in the CMD or in their compound variability. Additionally, we found 187 stars with ambiguous variability characteristics that could be either caused by cold spots in diskless stars, or by hot spots or extinction in disk-bearing stars, but further information is needed to confirm those as disk-bearing systems.

5. Summary and conclusions

In this paper, we summarized and discussed the results for the first near-IR photometric variability survey focusing on the low mass population of the Cygnus OB2 association. We observed the region in 115 nights spanning 217 days (∼7 months) using WFCAM/UKIRT JHK photometry. We focused our study on a list of 5083 Cygnus OB2 low mass candidate members identified by previous studies in the literature, 25% of which were listed as disk-bearing stars. While characterizing them according to their JHK variability, we verified that:

  • 1.

    2529 (∼50%) candidate members are significantly variable in the JHK-bands.

    • (a)

      1280 of them are highly variable and present large values of the Stetson variability index (Sample S1:S ≥ 1);

    • (b)

      1249 are low amplitude variables and present moderate values of the Stetson variability index (Sample S2:1 <  S ≤ 0.25).

  • 2.

    Given the morphology of the light curves:

    • (a)

      1679 of variable stars show some periodicity. 10% of those are disk-bearing stars, 39 (2.3%) are also periodic in color, 29 (1.7%) have assymetric light curves and 1256 have a period measured in Paper I.

    • (b)

      124 variable stars show occultation variations; 71% of which are disk-bearing stars. Among these, 83 are disk-bearing stars likely associated with AA Tau variability.

    • (c)

      726 variable stars show other types of variability, including long-term and short-term nonperiodic variability, and mixed light curve morphologies. Of those, 560 (77%) are disk-bearing stars. Also, 11 are burst-variables.

  • 3.

    Disk-bearing stars are more likely to be high amplitude variables. The main results supporting this statement are:

    • (a)

      While ∼57% of the sample S1 is composed of disk-bearing stars, only ∼8% of the sample S2 have disks.

    • (b)

      There is a steep increase in the disk fraction for higher variability amplitude in all bands.

  • 4.

    Periodic variability is more common among diskless stars: while 88% of all variable diskless stars present periodic behavior in their light curves, only 21% of the variable disk-bearing stars show the same behavior.

We also explored their variability in the color-space and verified that:

  • 5.

    Some variable stars present correlated variability in the color space resulting in linear trajectories inside the color-color and color-magnitude diagrams. We measured the slopes of these trajectories for 335 of such stars. We compared the slopes measured with the slopes expected from a theoretical description of the physical mechanisms typically attributed as the causes of near-IR variability in PMS stars. Considering the slopes measured inside the color-magnitude diagrams, we found that overall, our sample of stars with measured slopes is composed of a mixture of disk-bearing variability and circumstellar extinction variability, with 1.7 times more stars with variability caused by extinction/spot than caused by changes in the inner disk or accretion.

  • 6.

    149 variable stars show a compound behavior in color. Those are variables in which the color-behavior is not restricted to a linear trajectory. We observed 3 different types of compound behavior:

    • (a)

      Occultation with distinct color, which are variables that present occultation features with a distinct color from the rest of the light curve (42 sources).

    • (b)

      Variable near-IR excess (Increasing or decreasing) that causes the data points observed in different seasons to be clumped in different regions of the color-color and color magnitude diagrams (24 sources).

    • (c)

      Changing slope, which are variables showing a clear change in the dominant variability mechanism, often marked by a change in the slope of variability in the color-magnitude diagrams (24 sources).

    • (d)

      A mix of (b) and (c) (59 sources).

  • 7.

    We identified 15 eclipsing binary candidates, four of them previously identified by HSP11 and two of them presenting two periods: a period due to the eclipse, and a second period probably due to the rotational modulation of spots in the primary star. We report periods for 11 of the eclipsing binary candidates.

Characterizing the variability of young stars in the near-IR is an efficient way of identifying disk-bearing stars and that 64.5% of the previously known disk-bearing stars in Cygnus OB2 were variable in the near-IR. We also verified that 53.4% (681 sources) of the known disk-bearing stars in the association can be identified as disk-bearing stars based solely on their near-IR variability. Additionally, among the member candidates of CygOB2 previously identified based on their X-ray properties, we found 18 newly identified disk-bearing stars.


1

The χ2-test is defined as , where mi is the ith observation in a certain band, is the weighted mean value, σi is the photometric uncertainty for the ith observation, and N is the total number of observations. We calculated χ2 for each near-IR band and color. See Cody & Hillenbrand (2010), Rebull et al. (2014), Sokolovsky et al. (2017) for examples of application of the χ2-test.

2

Both Spearman nonparametric rank correlation and Orthogonal Distance Regression are included in the SciPy package for Python 2.7.

3

The detectability of a binary system as an eclipsing binary depends on the inclination of the system and its orbital period and separation. Thus the chances of detecting eclipses are higher for shorter orbital periods, because the wider the separation of components is, the longer the orbital period gets, and the narrower is the range of inclinations in which the eclipse is visible. This way, we expect most of the eclipsing binaries detected to have short periods.

4

Available in Python as part of the astropy.timeseries: https://docs.astropy.org/en/stable/timeseries/index.html

Acknowledgments

This study was part of JR Ph.D. thesis, which was granted by CNPq (Conselho Nacional de Desenvolvimento Científico e Tecnológico) and CAPES (Coordenação de Aperfeiçoamento de Pessoal de Nível Superior). JR has also received funding from the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation programme (grant AWESoMeStars, agreement No. 682393). S.H.P.A. acknowledges financial support from CNPq, CAPES, and FAPEMIG. J.B. acknowledges the support of ANR grant 2011 Blanc SIMI5-6 020 01 Toupies: Toward understanding the spin evolution of stars (http://ipag.osug.fr/Anr_Toupies/). We also thank Laura Venuti, Francisco Maia, and Alana Souza for useful comments over the development of the work. This study has made use of NASA’s Astrophysics Data System Bibliographic Services. This research made extensive use of TOPCAT software Taylor (2005). This research made use of Astropy (http://www.astropy.org), a community-developed core Python package for Astronomy (Astropy Collaboration 2013, 2018). The authors also thank the referee, Scott Wolk, for his detailed comments which helped to improve this paper. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Mauna Kea has always had within the indigenous Hawaiian community. We are most fortunate to have the opportunity to use observations from this mountain.

References

  1. Alencar, S. H. P., Teixeira, P. S., Guimarães, M. M., et al. 2010, A&A, 519, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Artemenko, S. A., Grankin, K. N., & Petrov, P. P. 2012, Astron. Lett., 38, 783 [NASA ADS] [CrossRef] [Google Scholar]
  3. Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Astropy Collaboration (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123 [Google Scholar]
  5. Beckwith, S. V. W., & Sargent, A. I. 1991, ApJ, 381, 250 [NASA ADS] [CrossRef] [Google Scholar]
  6. Beerer, I. M., Koenig, X. P., Hora, J. L., et al. 2010, ApJ, 720, 679 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bouvier, J., Cabrit, S., Fernandez, M., Martin, E. L., & Matthews, J. M. 1993, A&A, 272, 176 [NASA ADS] [Google Scholar]
  8. Bouvier, J., Alencar, S. H. P., Boutelier, T., et al. 2007, A&A, 463, 1017 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  9. Carpenter, J. M., Hillenbrand, L. A., & Skrutskie, M. F. 2001, AJ, 121, 3160 [NASA ADS] [CrossRef] [Google Scholar]
  10. Casali, M., Adamson, A., Alves de Oliveira, C., et al. 2007, A&A, 467, 777 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Claret, A., & Bloemen, S. 2011, A&A, 529, A75 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Clarke, D. 2002, A&A, 386, 763 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  13. Cody, A. M., & Hillenbrand, L. A. 2010, ApJS, 191, 389 [NASA ADS] [CrossRef] [Google Scholar]
  14. Cody, A. M., Stauffer, J., Baglin, A., et al. 2014, AJ, 147, 82 [NASA ADS] [CrossRef] [Google Scholar]
  15. Colombo, S., Orlando, S., Peres, G., Argiroffi, C., & Reale, F. 2016, A&A, 594, A93 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  16. Comerón, F., Pasquali, A., Rodighiero, G., et al. 2002, A&A, 389, 874 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  17. Fernandez, M., & Eiroa, C. 1996, A&A, 310, 143 [NASA ADS] [Google Scholar]
  18. Fonseca, N. N. J., Alencar, S. H. P., Bouvier, J., Favata, F., & Flaccomio, E. 2014, A&A, 567, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  19. Froebrich, D., Ray, T. P., Murphy, G. C., & Scholz, A. 2005, A&A, 432, L67 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Garufi, A., Meeus, G., Benisty, M., et al. 2017, A&A, 603, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Grankin, K. N. 2013, Astron. Lett., 39, 251 [CrossRef] [Google Scholar]
  22. Guarcello, M. G., Wright, N. J., Drake, J. J., et al. 2012, ApJS, 202, 19 [NASA ADS] [CrossRef] [Google Scholar]
  23. Guarcello, M. G., Drake, J. J., Wright, N. J., et al. 2013, ApJ, 773, 135 [NASA ADS] [CrossRef] [Google Scholar]
  24. Guarcello, M. G., Drake, J. J., Wright, N. J., et al. 2015, ArXiv e-prints [arXiv:1501.03761] [Google Scholar]
  25. Guarcello, M. G., Flaccomio, E., Micela, G., et al. 2017, A&A, 602, A10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  26. Guarcello, M. G., Flaccomio, E., Micela, G., et al. 2019, A&A, 628, A74 [CrossRef] [EDP Sciences] [Google Scholar]
  27. Hanson, M. M. 2003, ApJ, 597, 957 [NASA ADS] [CrossRef] [Google Scholar]
  28. Henderson, C. B., Stanek, K. Z., Pejcha, O., & Prieto, J. L. 2011, ApJS, 194, 27 [NASA ADS] [CrossRef] [Google Scholar]
  29. Herbst, W., Herbst, D. K., Grossman, E. J., & Weinstein, D. 1994, AJ, 108, 1906 [NASA ADS] [CrossRef] [Google Scholar]
  30. Hewett, P. C., Warren, S. J., Leggett, S. K. T. D., & Hodgkin, S. T. 2006, MNRAS, 367, 454 [NASA ADS] [CrossRef] [Google Scholar]
  31. Itoh, Y., Tamura, M., & Gatley, I. 1996, ApJ, 465, L129 [NASA ADS] [CrossRef] [Google Scholar]
  32. Joy, A. H. 1945, ApJ, 102, 168 [NASA ADS] [CrossRef] [Google Scholar]
  33. Kiminki, D. C., Kobulnicky, H. A., Vargas Álvarez, C. A., Alexander, M. J., & Lundquist, M. J. 2015, ApJ, 811, 85 [NASA ADS] [CrossRef] [Google Scholar]
  34. Knödlseder, J. 2000, A&A, 360, 539 [NASA ADS] [Google Scholar]
  35. Koldoba, A. V., Ustyugova, G. V., Romanova, M. M., & Lovelace, R. V. E. 2008, MNRAS, 388, 357 [NASA ADS] [CrossRef] [Google Scholar]
  36. Kraus, A. L., & Hillenbrand, L. A. 2007, AJ, 134, 2340 [NASA ADS] [CrossRef] [Google Scholar]
  37. Kulkarni, A. K., & Romanova, M. M. 2008, MNRAS, 386, 673 [NASA ADS] [CrossRef] [Google Scholar]
  38. Kurosawa, R., & Romanova, M. M. 2013, MNRAS, 431, 2673 [NASA ADS] [CrossRef] [Google Scholar]
  39. Lomb, N. R. 1976, Ap&SS, 39, 447 [Google Scholar]
  40. Lucas, P. W., Smith, L. C., Contreras Peña, C., et al. 2017, MNRAS, 472, 2990 [CrossRef] [Google Scholar]
  41. Lynden-Bell, D., & Pringle, J. E. 1974, MNRAS, 168, 603 [NASA ADS] [CrossRef] [Google Scholar]
  42. Mathis, J. S. 1990, ARA&A, 28, 37 [NASA ADS] [CrossRef] [Google Scholar]
  43. Meng, H. Y. A., Rieke, G. H., Kim, J. S., et al. 2019, ApJ, 878, 7 [CrossRef] [Google Scholar]
  44. Meyer, M. R., Calvet, N., & Hillenbrand, L. A. 1997, AJ, 114, 288 [NASA ADS] [CrossRef] [Google Scholar]
  45. Miyake, K., & Nakagawa, Y. 1993, Icarus, 106, 20 [NASA ADS] [CrossRef] [Google Scholar]
  46. Morales-Calderón, M., Stauffer, J. R., Hillenbrand, L. A., et al. 2011, ApJ, 733, 50 [NASA ADS] [CrossRef] [Google Scholar]
  47. Palmer, D. M. 2009, ApJ, 695, 496 [NASA ADS] [CrossRef] [Google Scholar]
  48. Pigulski, A., & Kołaczkowski, Z. 1998, MNRAS, 298, 753 [NASA ADS] [CrossRef] [Google Scholar]
  49. Pinte, C., Padgett, D. L., Ménard, F., et al. 2008, A&A, 489, 633 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  50. Press, W. H., Flannery, B. P., Teukolsky, S. A., & Veterling, W. T. 1992, Numerical Recipes in C: The Art of Scientific Computing, 2nd edn. (Cambridge: Cambridge University Press) [Google Scholar]
  51. Rauw, G., Nazé, Y., Wright, N. J., et al. 2015, ApJS, 221, 1 [NASA ADS] [CrossRef] [Google Scholar]
  52. Rebull, L. M., Cody, A. M., Covey, K. R., et al. 2014, AJ, 148, 92 [NASA ADS] [CrossRef] [Google Scholar]
  53. Rice, T. S., Wolk, S. J., & Aspin, C. 2012, ApJ, 755, 65 [NASA ADS] [CrossRef] [Google Scholar]
  54. Rice, T. S., Reipurth, B., Wolk, S. J., Vaz, L. P., & Cross, N. J. G. 2015, AJ, 150, 132 [NASA ADS] [CrossRef] [Google Scholar]
  55. Rodmann, J., Henning, T., Chandler, C. J., Mundy, L. G., & Wilner, D. J. 2006, A&A, 446, 211 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  56. Romanova, M. M., Ustyugova, G. V., Koldoba, A. V., & Lovelace, R. V. E. 2012, MNRAS, 421, 63 [NASA ADS] [Google Scholar]
  57. Romanova, M. M., Ustyugova, G. V., Koldoba, A. V., & Lovelace, R. V. E. 2013, MNRAS, 430, 699 [NASA ADS] [CrossRef] [Google Scholar]
  58. Roquette, J., Bouvier, J., Alencar, S. H. P., Vaz, L. P. R., & Guarcello, M. G. 2017, A&A, 603, A106 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  59. Sadavoy, S. I., Stutz, A. M., Schnee, S., et al. 2016, A&A, 588, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Scargle, J. D. 1982, ApJ, 263, 835 [Google Scholar]
  61. Scholz, A., Xu, X., Jayawardhana, R., et al. 2009, MNRAS, 398, 873 [NASA ADS] [CrossRef] [Google Scholar]
  62. Sokolovsky, K. V., Gavras, P., Karampelas, A., et al. 2017, MNRAS, 464, 274 [NASA ADS] [CrossRef] [Google Scholar]
  63. Sousa, A. P., Alencar, S. H. P., Bouvier, J., et al. 2016, A&A, 586, A47 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  64. Stauffer, J., Cody, A. M., Baglin, A., et al. 2014, AJ, 147, 83 [NASA ADS] [CrossRef] [Google Scholar]
  65. Stetson, P. B. 1996, PASP, 108, 851 [NASA ADS] [CrossRef] [Google Scholar]
  66. Taylor, M. B. 2005, in Astronomical Data Analysis Software and Systems XIV, eds. P. Shopbell, M. Britton, & R. Ebert, ASP Conf. Ser., 347, 29 [Google Scholar]
  67. Venuti, L., Bouvier, J., Irwin, J., et al. 2015, A&A, 581, A66 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. Venuti, L., Bouvier, J., Cody, A. M., et al. 2017, A&A, 599, A23 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  69. Vrba, F. J., Herbst, W., & Booth, J. F. 1988, AJ, 96, 1032 [NASA ADS] [CrossRef] [Google Scholar]
  70. Williams, J. P., & Cieza, L. A. 2011, ARA&A, 49, 67 [NASA ADS] [CrossRef] [Google Scholar]
  71. Wolk, S. J., Rice, T. S., & Aspin, C. 2013, ApJ, 773, 145 [NASA ADS] [CrossRef] [Google Scholar]
  72. Wolk, S. J., Günther, H. M., Poppenhaeger, K., et al. 2018, AJ, 155, 99 [NASA ADS] [CrossRef] [Google Scholar]
  73. Wright, N. J., Drake, J. J., Drew, J. E., & Vink, J. S. 2010, ApJ, 713, 871 [NASA ADS] [CrossRef] [Google Scholar]
  74. Wright, N. J., Drake, J. J., Guarcello, M. G., et al. 2014, ArXiv e-prints [arXiv:1408.6579] [Google Scholar]
  75. Wright, N. J., Drew, J. E., & Mohr-Smith, M. 2015, MNRAS, 449, 741 [NASA ADS] [CrossRef] [Google Scholar]

Appendix A: Estimation of typical amplitudes and slopes of spot and extinction variability in the near-IR colors

A.1. Spot variability

To investigate the near-IR variability amplitudes produced by spots in the stellar surface we used a simple model in which two blackbodies with different temperatures describe the stellar surface and the spot, and in which the limb-darkening effects are considered. In such model, the variation in magnitude caused by a spot of temperature Tspot at the surface of a star of photospheric temperature T* can be estimated as:

(A.1)

where f is the fraction of the stellar surface covered by spots (filling factor), Bλ(T) is Plank Law for the radiation of a blackbody of temperature T at wavelength λ, and μ is the limb-darkening coefficient. This model does not take into account the shape, position, and number of spots.

Previous studies (e.g., Carpenter et al. 2001; Scholz et al. 2009) apply similar models for estimating amplitudes and slopes of spot variability in the near-infrared color space using a small set of parameters. Here, instead, we used the parameters retrieved from observations of spotted stars. We estimated Δm for each of the near-infrared bands by feeding the model with the parameters Teff, Tspot and f estimated by Venuti et al. (2015) for 360 stars in NGC 2264 using optical and UV observations. We only used the parameters of stars with spectral types between F2 and M2, which covers approximately the mass interval 0.3 M to 1.5 M. This same mass range contains most of the stars in our member candidate list. Venuti et al. (2015) sample includes spotted stars with ΔT (difference of temperature between the stellar photosphere and the spot) up to −6820 K for hot spots, and 7030 K for cold spots. The limb-darkening coefficients for the WFCAM JHK photometric system were provided by Dr. Antônio Claret, and were estimated as described in Claret & Bloemen (2011).

For cold spots, we found typical (median) amplitudes of 0.10, 0.08, and 0.07 mag for the J, H, and K bands. These values are in good agreement with the median amplitudes for the sample of diskless periodic stars analyzed in Paper I, which were 0.10, 0.08, and 0.08 mag, supporting the idea that the periodic sample was composed mainly of stars with cold spots. On the other hand, we found that the stars with hot spots from Venuti et al. (2015) have median near-IR amplitudes of 0.07, 0.05, 0.04 mag for J, H, and K bands respectively. These values are startling low compared to the amplitudes of periods in disk-bearing sample, whose median values were 0.31, 0.30, 0.20 mag.

We showed the range of slopes of variability we found in Table 7 and discussed it in Sect. 4.1. Overall we found slopes for cold and hot spot variability similar to previous studies. However, we found typical higher amplitudes for cold spots, than for hot spots, with median and max values of 0.09 and 0.61 mags in the J band for cold spots, and 0.05 and 0.5 for hot spots. The amplitudes we found were systematically lower than in previous studies for both cold and hot spots. The explanation for that is that among the spotted stars from Venuti et al. (2015) there is a lack of hot spotted stars with both a very hot spot and a large spot coverage (We encourage the interested reader to refer to Venuti et al. 2015 Sect. 3.4, and in particular to their Fig. 10). Thereby, the hypothetical star producing the max values of typical amplitudes for spotted stars in Carpenter et al. (2001) – a spotted star with a photospheric temperature of 4000 K, a spot temperature of 8000 K, and spot coverage of 30% – is not possible within the observational scenario evaluated by Venuti et al. (2015).

It may be the case that hot spotted stars with large ΔT and spot coverage are more common than reflected by Venuti et al. (2015) sample of spotted stars, and in this case, Carpenter et al. (2001) and Scholz et al. (2009) models and parameters are reflecting more realistic amplitude ranges for this spot variability than our results. If this is not the case, then the sample of disk-bearing periodic stars presented in Paper I cannot be composed only of hot spotted stars, and the observed variations may be caused instead by other types of periodic variability mechanisms such as variable circumstellar extinction. Either way, the spot model used and the set of parameters estimated by Venuti et al. (2015) show that there are sets of observed spot temperature and spot coverage that result in cold spot variability amplitude as high as those produced by hot spots. Hence, it is not possible to differentiate hot spots from cold spots in the near-IR based only in their variability amplitudes. It is also not possible to distinguish between hot and cold spots using their slopes, and Table 7 shows a significant overlap for the typical slopes in the two cases.

A.2. Extinction variability

The wavelength dependence of extinction in the near-IR and submillimeter can be described as a power law: Aλ ∼ λβ (Scholz et al. 2009), where β is the dust opacity index and represents the efficiency at which dust grains radiate at long wavelengths (e.g., Sadavoy et al. 2016). A typical value for β in the interstellar medium in the IR is β ∼ 1.7 (Mathis 1990). β can vary even inside the same cloud (Froebrich et al. 2005), as it is a function of the size of the grain producing the opacity, and it can evolve with both density and temperature. In circumstellar disks, a typical value is β ∼ 1 (Beckwith & Sargent 1991; Miyake & Nakagawa 1993; Williams & Cieza 2011), but it can present values as low as 0.4 and as high as 1.6 (e.g., Rodmann et al. 2006; Pinte et al. 2008; Garufi et al. 2017). We used values in this range to estimate the slopes described by variable extinction in the CMD (Table 7 and Fig. 14). Due to the nature of this type of variability, it can assume arbitrary amplitudes.

All Tables

Table 1.

General information for the nonvariable stars identified within this survey.

Table 2.

General information for the Variable stars identified and studied in the present paper.

Table 3.

Information about variability classification for the stars studied in the present paper.

Table 4.

Peak to peak and root mean square amplitudes in each observed band and color for each of the variable stars studied in the present paper.

Table 5.

Parameters for the stars with measured slopes in the color-color diagram or color-magnitude diagrams.

Table 6.

Summary of the variability classification performed in Sect. 3 for both disk-bearing and diskless samples.

Table 7.

Slopes of variability in the near-IR expected from models and other observational studies.

Table 8.

Information about the eclipsing binary candidates identified in the present study.

All Figures

thumbnail Fig. 1.

Near-IR (JHK) color-color diagrams. Left and middle panels: median colors of the 3811 diskless stars (blue dots) and the 1272 disk-bearing stars (red dots) analyzed in the present study. Right panel: CygOB2-093173 variable star. Black dots show the color of the star in each observed night. The regions labeled as “P”, “D”, and “E” are the regions typically occupied by stars with normal photospheric colors, disk-bearing colors, and extreme near-IR excess respectively. The dark red arrow shows the effect of the interstellar reddening for the median visual extinction expected for stars in the CygOB2 association (AV = 4.1m). The gray solid curve shows an empirical dwarf sequence by Kraus & Hillenbrand (2007) and the orange line shows the Meyer et al. (1997) T Tauri locus. Dashed gray lines show extinction paths of AV = 10m for 0.6 M and 0.08 M stars, and the truncation point of the T Tauri locus.

In the text
thumbnail Fig. 2.

Distribution of Stetson variability index as a function of H median magnitude for the full sample of candidate members (left), disk-bearing stars (middle), and diskless stars (right). The red line shows S = 1, and the dashed black line shows S = 0.25.

In the text
thumbnail Fig. 3.

Disk fraction as a function of peak to peak amplitude for each near-IR band (left) and color (right). In both plots, circles show the center of each ptp amplitude bin and horizontal bars indicate the size of each bin. Vertical bars show the standard errors of a Poisson counting in that bin.

In the text
thumbnail Fig. 4.

Star CygOB2-093173 is a disk-bearing (Class II) star and an example of non periodic variable star, with S = 15.64, and that presents reddening linear trajectories in the color-color and color-magnitude diagrams. Left, from top to bottom: J light curve, a color bar showing the color attributed for each observation date, K light curve, and H − K color curve. The dotted lines show maximum and minimum values in each light curves, and dashed lines show the median values. Right, from top to bottom: J vs. J − H color-magnitude diagram, K vs. H − K color-magnitude diagram, and color-color diagram with P, D, and E regions (discussed in Sect. 2). The dark red arrows show the effect of an interstellar reddening of AV = 2m. The vertical and horizontal magenta line segments show the median error in each axis. In each plot, error bars are shown in gray. When visible and close to the position of the star, a gray curve shows a Kraus & Hillenbrand (2007) empirical isochrone in the color-color diagram, and the orange line shows the T Tauri Locus. Each diagram has a range of 5σ around the median color or magnitude for that star.

In the text
thumbnail Fig. 5.

Star CygOB2-000348 is an example of a diskless periodic variable, with S = 1.62, a ptp amplitude of 0.13 and 0.11 mag in the J and K bands respectively, and that does not show periodicity in color. The plots, colors, and symbols are the same as in Fig. 4, but with the light curves and color curve phased with the period of the star, P = 5.1 d.

In the text
thumbnail Fig. 6.

Star CygOB2-113855 is an example of a disk-bearing (Class II) periodic star. With S = 2.77 and ptp amplitude of 0.31 and 0.36 mag in the J and K bands, the star is also periodic in color. It is also an example of a compound variable with a variable IR-excess. The plots, colors, and symbols are the same as in Fig. 4, but with the light curves and color-curve phased with the period of the star, P = 10.3 d.

In the text
thumbnail Fig. 7.

Disk-bearing (Class II) star CygOB2-124656 is an example of an occultation variable that presents reddening linear trajectories in the color-color and color-magnitude diagrams. The star has S = 3.59 and ptp amplitude of 0.43 and 0.17 in the J and K bands. The plots, colors, and symbols are the same as in Fig. 4.

In the text
thumbnail Fig. 8.

Diskless star CygOB2-041510 is an example of a periodic occultation variable considered as an eclipsing binary candidate. It is also an example of a compound variable with occultation data-points with distinct color behavior. The star has P = 4.68 days, S = 3.06, and ptp amplitude of 0.28 in both J and K bands. The plots, colors, and symbols are the same as in Fig. 4.

In the text
thumbnail Fig. 9.

Disk-bearing (flat-spectrum) star CygOB2-113911 is an example of a long-term variable with S = 3.67, that presents a blueing linear trajectory in the color-magnitude diagrams. The plots, colors, and symbols are the same as in Fig. 4.

In the text
thumbnail Fig. 10.

Color-magnitude and color-color diagrams with the slopes for the stars presenting a linear variation. Red lines show disk-bearing stars, and blue lines show stars without discs. The length of each line segment reflects the amplitude of the variability. The dark red arrows show the median interstellar extinction toward CygOB2, and the gray lines show as reference of unreddened normal photosphere colors, an empiric isochrone from Kraus & Hillenbrand (2007).

In the text
thumbnail Fig. 11.

Disk-bearing (Class II) star CygOB2-041756 (S = 10.4) is an example of a mixed variable with a long term variability mixed with burst features in the first epochs observed. The star is a compound variable with asymmetries in the light curve having a distinct color behavior. The plots, colors, and symbols are the same as in Fig. 4.

In the text
thumbnail Fig. 12.

Disk-bearing (Class II) star CygOB2-021188 (S = 1.0) is an example of compound variability showing a variable near-IR excess and a changing slope during the observed time. The plots, colors, and symbols are the same as in Fig. 4.

In the text
thumbnail Fig. 13.

Examples of stars with a compound behavior are shown inside color-color diagrams or color-magnitude diagrams. 1- Asymmetry with a distinct color; 2- Variable IR-excess; 3- Changing slope; 5- Mix of 2 and 3. Dark red arrows show the effect of interstellar reddening for AV = 1 mag. The color of each data point reflects its observational date following the color-bar in Fig. 4.

In the text
thumbnail Fig. 14.

Distribution of slopes measured in the color-color and color-magnitude diagrams. The figure shows slopes for J − H vs. H − K (top), J vs. J − H (middle), and K vs. H − K (bottom) diagrams. We plotted the expected slopes for each physical mechanism as line segments, following the ranges in Table 7. The X mark shows the exact position of the T Tauri Locus. The black circle shows the slope expected for the typical dust opacity index in the ISM. Black bars show the stars in S1 sample and with precision in the light curve better than 5%, empty bars show stars in the S2 sample or with precision worse than 5%.

In the text
thumbnail Fig. 15.

Four eclipsing binary candidates that were not periodic in our study. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are written in the plot. Right: J vs. J − H CMDs. In both left and right plots, colors reflect the date of the observation as in Fig. 4.

In the text
thumbnail Fig. 16.

Two multiperiodic eclipsing binary candidates. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band, folded using the detected period attributed to the occultation feature. J, H, and K bands are shown as blue, green, and red respectively. Figure 17 shows the folded light curves for the secondary period. Right: J vs. J − H CMD. Colors reflect the date of the observation, as in Fig. 4.

In the text
thumbnail Fig. 17.

Phased light curve for each band, folded using the secondary detected period found for the eclipsing binary candidates presented in Fig. 16. J, H, and K bands are shown as blue, green, and red respectively.

In the text
thumbnail Fig. 18.

Eclipsing binary candidate CygOB2-041510 (Fig. 8). Left: J, H, and K light curves for the star. The median magnitude and ptp amplitude are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Right: phased light curve for each band, folded using the detected period attributed to the occultation feature. J, H, and K bands are shown as blue, green, and red respectively.

In the text
thumbnail Fig. 19.

Five periodic eclipsing binary candidates. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band, folded using the detected period attributed to the occultation feature. J, H, and K bands are shown as blue, green, and red respectively. Right: J vs. J − H CMD. Colors reflect the date of the observation, as in Fig. 4.

In the text
thumbnail Fig. 20.

Stars CygOB2-000644 and CygOB2-072104 are both occultation and compound variables and were considered as eclipsing binary candidates in the present study. They are also reported as eclipsing binaries in HSP11. Left: J, H, and K light curves for each star. The median magnitude and ptp amplitude are also written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band. J, H, and K bands are shown as blue, green, and red respectively. Right: J vs. J − H CMD for the same star. Colors also reflect the date of the observation, as in Fig. 4.

In the text
thumbnail Fig. 21.

Star CygOB2-114152 is a periodic candidate and compound in this survey but it is reported as an eclipsing binary in HSP11. Left: J, H, and K light curve. The median magnitude and ptp amplitudes are written in the plot, and the colors reflect the date of the observation, as in Fig. 4. Middle: phased light curve for each band. J, H, and K bands are shown as blue, green, and red respectively. The period reported by HSP11 (below our detection limit) was adopted in the phased light curve. Right: J vs. J − H CMD. Colors reflect the date of the observation, as in Fig. 4.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.