# Mid-Infrared Variability of Young Stellar Objects on Timescales of Days to
Years

Sieun Lee  School of Space Research, Kyung Hee University, 1732, Deogyeong-
daero, Giheung-gu, Yongin-si, Gyeonggi-do 17104, Korea  Jeong-Eun Lee
Department of Physics and Astronomy, Seoul National University, 1 Gwanak-ro,
Gwanak-gu, Seoul 08826, Korea  SNU Astronomy Research Center, Seoul National
University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826, Republic of Korea  Carlos
Contreras Peña Department of Physics and Astronomy, Seoul National University,
1 Gwanak-ro, Gwanak-gu, Seoul 08826, Korea  Doug Johnstone NRC Herzberg
Astronomy and Astrophysics, 5071 West Saanich Road, Victoria, BC, V9E 2E7,
Canada  Department of Physics and Astronomy, University of Victoria, 3800
Finnerty Road, Elliot Building, Victoria, BC, V8P 5C2, Canada  Gregory Herczeg
Kavli Institute for Astronomy and Astrophysics, Peking University, Yiheyuan 5,
Haidian Qu, 100871 Beijing, China  Seonjae Lee Department of Physics and
Astronomy, Seoul National University, 1 Gwanak-ro, Gwanak-gu, Seoul 08826,
Korea

###### Abstract

Variability in the brightness of Young Stellar Objects (YSOs) is a common
phenomenon that can be caused by changes in various factors, including
accretion, extinction, disk morphology, interactions between the disk and the
stellar photosphere, and the rotation of hot or cold magnetic spots on the
stellar photosphere. Analyzing the variability on different timescales gives
insight into the mechanisms driving the brightness changes in YSOs. We
investigate the variability of YSOs on both long and short timescales using
two mid-infrared datasets: the NEOWISE 7.5-year W2 (4.6$\mu$m) data and the
YSOVAR 40-day Spitzer/IRAC2 (4.5$\mu$m) data, respectively. We classify the
variability types in each timescale following Park et al. (2021). We find a
higher detection rate of variable sources in the short-term (77.6%) compared
to the long-term (43.0%) due to the higher sensitivity of the Spitzer
observations. In addition, the higher cadence of the YSOVAR data results in
the weeks-long short-term variability being mostly secular, while the years-
long long-term variability explored with the coarsely sampled NEOWISE data is
mostly stochastic. By cross-matching the two catalogs, we statistically
analyze the variability types exhibited by YSOs across both timescales. The
long-term variability amplitude is mostly three times (up to ten times)
greater than the short-term variability. Furthermore, we evaluate variability
on very short (1–2 day) timescales and recover a trend of the increasing
amplitude of variability as the timescales increase. By comprehensively
analyzing the variability of YSOs over various timescales, we contribute to a
deeper understanding of the underlying mechanisms driving their variability.

##  1 Introduction

Variability is a common characteristic observed in Young Stellar Objects
(YSOs). Extensive time-series observations have been conducted in regions of
star formation to investigate the underlying physical mechanisms driving the
variability of these systems. These studies have spanned various timescales
and have covered a wide range of wavelengths, including optical, near-
infrared, mid-infrared, and sub-mm observations (e.g., Herbst et al., 1994;
Carpenter et al., 2001; Morales-Calderón et al., 2011; Lee et al., 2021). The
observed changes in the brightness of YSOs can be attributed to several
mechanisms that affect the star, its disk, and its envelope. The physical
cause of these changes includes rotational modulation of cold or hot spots on
the stellar surface (Attridge & Herbst, 1992; Choi & Herbst, 1996), unstable
accretion processes (Romanova et al., 2013; Stauffer et al., 2014), variations
in extinction due to surrounding material (Bouvier et al., 2013), disk
instabilities (Bino et al., 2023), and the presence of eclipsing binary
systems (Stassun et al., 2004, 2006; Cargile et al., 2008).

Time-series photometry at mid-infrared (mid-IR) or longer wavelengths is a
valuable tool for studying the properties of the young protostars that are
still deeply embedded within their nascent envelopes. A recent study by Park
et al. (2021), hereafter referred to as P21, focused on analyzing 6000 known
Young Stellar Objects (YSOs) in nearby star-forming regions using 6.5 years of
mid-IR (3.4$\mu$m and 4.6$\mu$m) time-series photometry from the NEOWISE
mission provided at a 6 months cadence. Using the averaged single-epoch
exposures from NEOWISE, P21 discovered and characterized the long-term
lightcurves of 1700 variable YSOs. Their analysis revealed that changes in
extinction along the line of sight and variable accretion were the primary
mechanisms driving the mid-IR variability in these YSOs. Due to the cadence of
the NEOWISE epoch-averaged observations used in their study, P21 did not
investigate variability at timescales shorter than 180 days, although several
YSOs with significant variability on short (1-2 day) timescales were noted.

Short timescale (up to $\sim$ 40 days) mid-IR variability of YSOs has been
extensively studied through the Young Stellar Object VARiability (YSOVAR)
program (Morales-Calderón et al., 2011; Cody et al., 2014; Stauffer et al.,
2014; Wolk et al., 2018). For instance, Morales-Calderón et al. (2011)
conducted a comprehensive monitoring of YSOs in the Orion Nebula Cluster (ONC)
region using high-cadence mid-IR Spitzer IRAC1 (3.6$\mu$m) and IRAC2
(4.5$\mu$m) data over a 40-day period. Their work, along with the broader
YSOVAR program, demonstrated that a significant fraction of YSOs exhibits
variability on day and month timescales, attributed to rotational modulation
from stellar spots, variable extinction (analogous to AA Tau stars), and
stochastic bursts (see also Cody et al., 2014; Stauffer et al., 2014).

Our study aims to establish connections between YSOs’ short- and long-term
mid-IR variability by collating data from the Spitzer and WISE missions.
Specifically, we analyze the mid-IR data of YSOs in the ONC with photometric
observations available in both the YSOVAR (Morales-Calderón et al., 2011) and
NEOWISE (P21) datasets. Through this analysis, we aim to gain a deeper
understanding of the underlying physical mechanisms driving YSO variability
across timescales ranging from 1-2 days to up to 7.5 years ($\sim 2500$ days).

##  2 DATA

To analyze the variability of YSOs over many timescales, we analyze data from
two surveys: YSOVAR and NEOWISE. Here, we provide an overview of these surveys
and the methods employed for data analysis.

###  2.1 NEOWISE : Long-term (years) data

The NEOWISE (Near-Earth Object WISE) survey extends the WISE mission. Since
its reactivation in 2013, NEOWISE has been collecting all-sky photometric data
at 3.4$\mu$m (W1) and 4.6$\mu$m (W2) (Mainzer et al., 2014). During a single
epoch of observation, about 10 $\sim$ 20 exposures are taken over a span of a
few days, typically ranging over 1 to 2 days. These observations are repeated
every 6 months.

In P21, a catalog of 5398 nearby YSOs was compiled, encompassing nearby star-
forming regions such as Orion A, B (Megeath et al., 2012), Taurus (Esplin &
Luhman, 2019), and the rest of the Gould Belt (Dunham et al., 2015). The
NEOWISE W1 and W2 photometric data acquired between 2013 and 2020 were
extracted for this YSO sample. For each object, P21 averaged 1-2 days of
exposures within each epoch to create lightcurves with a single photometric
point every six months. This resulted in a maximum of 14 epochs per source,
enabling the characterization of the long-term 6.5-year mid-IR photometry. The
NEOWISE 2022 data release111Available at : NEOWISE Team (2020) includes
photometry from two additional epochs that were not available at the time of
the P21 analysis. This new data release expands the observational baseline of
NEOWISE to 7.5 years.

For our study, we used the single-epoch catalog of NEOWISE and performed a
query within a 3″ radius around the coordinates of each YSO, using the same
star-forming region catalogs as P21. For this set of measurements, we
calculated the mean and standard deviation of the distances ($sd_{d}$) between
the NEOWISE single exposure coordinates and the known coordinates of the YSO.
We included exposures that met the following criteria (Figure 1):

  * •

Signal-to-noise (SNR)$>7$;

  * •

Frame quality score (qual_frame) $>$0;

  * •

The source in the NEOWISE single exposure was located within 2*$sd_{d}$ of the
mean location overall NEOWISE measurements (where $sd_{d}$ is the standard
deviation).

To ensure adequate data for each investigated source, we selected targets with
more than 32 NEOWISE exposures, corresponding to an average of at least 2
exposures per single epoch. Subsequently, we averaged the single-epoch
exposures, typically obtained within 1 to 2 days, to generate regularly
sampled lightcurves with measurements available every 6 months.

Following the approach of P21, we applied additional criteria to construct a
high-quality dataset to analyze YSO variability. For our analysis, the targets
had to satisfy the following conditions:

  * •

Detected in more than 14 epochs in W2;

  * •

Standard deviation of distance from the known YSO coordinates smaller than
0.3″;

  * •

Mean W2 uncertainty smaller than 0.2 magnitudes.

After applying these criteria, we retained 5398 sources for further analysis.
We note that specific criteria used in analyzing the raw NEOWISE data for our
study differ from those employed in the P21 study. Our criteria are stricter
to select only the most robust measurements. Here, we note that the total
number of retained sources in this study and P21 is coincidentally the same,
but individual sources are not identical.

P21 employed less strict criteria than ours for individual exposures. Thus,
they chose to eliminate outlying data points, considering only 70% of the
data. During our analysis, however, we recognized that the excluded 30% of the
data demonstrated sufficient variability potential. As a result, we introduced
the aforementioned strict conditions and retain the entire dataset, as
detailed in Appendix A.

Figure 1: The excluded NEOWISE low-quality exposures. (Left) The plot
represents the original photometric data for W1 (red) and W2 (blue). Solid
circles indicate the exposures that satisfied the distance and signal-to-noise
ratio criteria, while open circles represent the data points that did not meet
the criteria. (Right) The x-axis corresponds to the right ascension (RA), and
the y-axis corresponds to the declination (Dec). The red dot at the center
represents the coordinates extracted from the original catalog, and the red
circle indicates the 3″search radius. The black dots represent all the
original detections that are subsequently removed due to our criteria, the
green circle represents the 2*$sd_{d}$ region, and the green dots represent
the exposures that satisfy this cut.

###  2.2 YSOVAR : Short-term (weeks) data

YSOVAR is a survey that used the Spitzer Space Telescope to investigate the
mid-IR photometric variability of YSOs. The pioneering work of Morales-
Calderón et al. (2011) focused on the YSOVAR survey conducted in the ONC
region, covering an area of 0.9 deg2. The observations were carried out in the
Fall of 2009, with time-series imaging in the mid-IR bands of IRAC1 (3.6
$\mu$m) and IRAC2 (4.5 $\mu$m). The observations were conducted over 40
consecutive days, with a cadence of approximately two observations per day,
resulting in 81 exposures, in total, in each band. This high-cadence
monitoring allowed for a detailed analysis of short-term variability of YSOs
within the ONC region. The resulting lightcurve catalog from Morales-Calderón
et al. (2011) contained 2067 sources, which served as valuable data for
investigating the short-term variability of YSOs. This dataset was used in our
study to analyze the short-term mid-IR variability of YSOs and establish
connections with the long-term variability observed in the NEOWISE data.

##  3 YSO VARIABILITY CLASSIFICATION

In this study, we adopted the YSO variability classification methodology
proposed by P21 for both the 7.5-year NEOWISE data and the 40-day YSOVAR data.
This classification scheme categorizes YSO variability into six distinct
types: secular variables as Linear, Curved, and Periodic, and stochastic
variables as Burst, Drop, and Irregular. The classification is based on
various parameters derived from the W2 lightcurves, including the standard
deviation of fluxes, the amplitude of variability, and the identification of
periodicity using the Lomb-Scargle Periodogram (Lomb, 1976; Scargle, 1989).
Linearly rising(+) or declining(–) trends in the lightcurves are examined.
Figure 2 depicts this variability classification methodology with a flowchart.

We made certain modifications to the definition of secular variability
(Linear, Curved, and Periodic sources) compared to P21 in order to account for
the different timescales studied using YSOVAR and the NEOWISE 2022 data
release. The basis of the timescale used for the classification in P21 was the
6.5 years covered by the NEOWISE lightcurves adopted in the study. Given the
different timescales of our datasets, we adjusted the period ranges to be more
appropriate for the respective data sets. The subsequent sections describe the
specific values and justifications for these new period ranges.

Figure 2: Variability Classification Methodology Scheme. $\Delta mag$ represents the difference between the maximum and minimum magnitudes, while $\sigma$(mag) signifies the mean value of magnitude uncertainty. $FAP_{LSP}$ denotes the False Alarm Probability from the Lomb-Scargle Periodogram, and $FAP_{Lin}$ indicates the False Alarm Probability from Linear Least Square fitting. SD stands for flux standard deviation, and $\sigma$(flux) represents the mean value of flux uncertainty. For periods, A is 200 days for NEOWISE data and 0 days for YSOVAR data. B is 1370 days for NEOWISE data and 20 days for YSOVAR data. C is 5480 days for NEOWISE data and 160 days for YSOVAR data.  Figure 3: Examples of NEOWISE lightcurves for each long-term variability type. The source name is shown at the top of each panel, and each variable type is indicated in the upper right corner. The first observation date (MJD) is marked on the very left x-axis of each panel. Curved and Periodic are depicted with LSP fittings shown as red lines.  Table 1: Long-term Variable Type | Linear | Curved | Periodic | Burst | Drop | Irregular | Non-var type | Total
---|---|---|---|---|---|---|---|---
NEOWISE | 86 [1.6]aaThe numbers in front are the count of variables for each type, while the numbers in the square bracket are the fractions (%) of variables relative to their total samples. The sum of the number in square bracket is 100. | 439 [8.1] | 154 [2.9] | 114 [2.1] | 21 [0.4] | 1505 [27.9] | 3079 [57.0] | 5398 [100]
+PeriodicbbCombined variable type | - | 4 | - | - | - | - | - | 4 (0.1) ccThe numbers in front are the total count for each combined variable type, while the numbers in parentheses are the fractions (%) of combined variables relative to total variable samples.
+Burst | - | 6 | 2 | - | - | - | - | 8 (0.1)
+Drop | 3 | - | - | - | - | - | - | 3 (0.1)
+Irregular | 20 | 138 | 5 | - | - | - | - | 163 (3.0)

###  3.1 Long-term data

In our classification of long-term variability using the NEOWISE W2
observations, we followed the criteria outlined in P21. The first criterion
was based on the magnitude domain, specifically looking for YSOs with a
magnitude difference $\Delta W2$ (${\rm=W2_{max}-W2_{min}}$) greater than
three times the mean uncertainty $\sigma(W2)$ in W2 magnitudes. We identified
4218 YSOs that satisfied this condition. For the classification of stochastic
variability, we followed the same criteria defined in P21. Among the targets
that satisfied the first condition, the criterion for classifying as Burst was
(median magnitude - minimum magnitude) $>$ $0.8\times$$\Delta$W2. For Drop,
the criterion was (maximum magnitude - median magnitude) $>$
$0.8\times$$\Delta$W2. For Irregular, the criterion was SD/$\sigma>3$, where
SD represents the standard deviation of flux, and $\sigma$ represents the mean
flux uncertainty.

Due to the longer timescales involved in our study, we adjusted the periods
used to define the different secular classes. YSOs with Lomb-Scargle
Periodogram (LSP) periods longer than 5480 days (15 yrs) have timescales
longer than twice the 7.5 yrs of NEOWISE monitoring and were therefore
classified as Linear variables. YSOs with periods ranging from 1370 days to
5480 days were classified as Curved variables. Finally, sources with periods
between 200 days and 1370 days were classified as Periodic variables.

Similar to the results in P21, some variable YSOs exhibit additional trends
besides their dominant long-term behavior. We conducted tests to determine
whether secular variables ($Linear$, $Curved$, $Periodic$) exhibit additional
variable types. Following the flow chart (Figure 2), we first removed the
dominant secular trend. We then tested for secondary secular or stochastic
variability by continuing the flow chart analysis, reinserting the residual
light curve at a lower decision box (diamond). Thus, for primary $Linear$ and
$Curved$ variables, secondary $Periodic$ and $Stochastic$ variations were
tested (by reinserting the residual light curve into the fourth diamond from
the top) while for primary $Periodic$ variables only $Stochastic$ secondary
variations were tested (by reinserting the residual into the fifth diamond
from the top). If a source shows both trends (Secular + additional trend), it
is classified as a combined variable type. This methodology is consistent with
that of P21. However, we found some variables with $Linear$ or $Curved$ types
that also exhibit a simultaneous periodic tendency. Therefore, we classified
these as $Secular+Periodic$, a category not present in P21.

Table 1 provides an overview of the long-term YSO variable type results based
on the NEOWISE data. Out of the 4218 sources that satisfy the first criterion
$\Delta\text{W2}/\sigma(\text{W2})>3$, 2319 were classified as having secular
or stochastic variability types. Among the variable YSOs, most were classified
as Irregular (64.9%), indicating unpredictable and stochastic variability.
Additionally, the combined variable type Secular+Irregular was dominant among
the YSOs exhibiting additional trends on top of the long-term behavior. To
visually illustrate the different long-term variable types, Figure 3 presents
example NEOWISE lightcurves for each variability type. These examples provide
a glimpse into the diverse patterns and behaviors exhibited by YSOs over many
years.

####  3.1.1 Comparison between P21 and this work

This section will examine the reasons behind the changes in the number and
types of variables between P21 and this study. In P21, out of 5398 YSOs, 1734
were identified as variables, while in this study, 2319 were identified as
variables, despite classifying variables by the same criteria. We attribute
this increase mainly to the following two factors:

  * •

Data coverage per epoch: P21 used only 70% of the data in each epoch, whereas
this study used 100% of the data points (see Appendix A).

  * •

Time coverage: P21 used 6.5 years of NEOWISE data, while this study used an
extended data set of 7.5 years.

One more difference between P21 and this work is the criterion for detected
epochs during initial sample determination, i.e., 5 epochs versus 14 epochs.
However, the stricter criterion in this work is expected to reduce the overall
number of samples and, consequently, the number of variables, which is the
opposite of what was found. To examine the difference, we investigate the
selection criterion of detection epochs under identical conditions (7.5 years
time coverage and 100% data coverage). According to our examination, the
number of targets detected in 5 or more epochs was 5696, with 2446 variables.
In contrast, the number of targets detected in 14 or more epochs was 5398,
with 2319 variables, indicating a 5% decrease in variables. As expected, the
detection epoch criterion does not cause the increase of variables in this
work.

Therefore, we analyze the two main factors of data coverage per epoch and time
coverage in more detail and compare the results with a heatmap, contrasting
them with the results of P21 (Figure 4). The left panel of Figure 4 shows the
differences in classification that arise from using 6.5 years of data but
analyzing 70% (as in P21) and 100% of the data in individual epochs (this
study). Each column represents the variability type identified when using 100%
of the data, and each row represents the variability type found when using 70%
of the data. If none of the sources changed their type between these analyses,
only the diagonal boxes would have non-zero values. For instance, in the box
located in the first row and second column, ‘10’ indicates the number of
sources that changed their variability type from $Linear$ when using 70% of
the data to $Curved$ when using 100% of the data.

On the other hand, the right panel shows the differences in using 70% of the
data but varying the time coverage from 6.5 to 7.5 years. In this panel, each
column shows the variability type found after using 7.5 years data, and each
row shows the variability after using 6.5 years data. As a result, while P21
identified 1734 variables out of 5398 YSOs, changing the coverage of data
points and time increased the number of variables to 2246 (left panel) and
2351 (right panel), respectively. Moreover, in both cases, there was a
significant increase in the number of variables changing from Non-variable to
Irregular types.

Therefore, in this study, which uses 100% data points in individual epochs and
extends the time coverage to 7.5 years, more variables are identified compared
to P21. Consequently, with the prospect of acquiring more data points with
time, the number of variables will continue to increase.

Figure 4: (Left) Comparing two analyses that used 6.5 years data but only
differed in data points per epoch (70% data (P21) and 100% data). (Right)
Comparing two analyses that used 70% data but differed in the time coverage
(6.5 years and 7.5 years). The value within each box indicates the sources
associated with that joint type.

###  3.2 Short-term data

To classify the short-term YSO variability, we followed a similar procedure as
described above and applied the methods of P21 to the YSOVAR IRAC2 catalog
presented by Morales-Calderón et al. (2011). Given the 40-day timescale of the
YSOVAR observations, we established LSP periods to differentiate the secular
variables in the short-term. YSOs with LSP periods longer than 160 days were
classified as Linear, those with periods between 20 days and 160 days were
classified as Curved, and those with periods shorter than 20 days were
classified as Periodic.

Table 2 provides an overview of the results obtained from the classification
of short-term YSO variability. Similar to the long-term classification, we
employed the first criterion for short-term classification, which requires
$\Delta\text{IRAC2}/\sigma(\text{IRAC2})>3$ in the magnitude domain. Here,
$\Delta\text{IRAC2}$ represents the difference between the maximum and minimum
magnitudes, and $\sigma(\text{IRAC2})$ represents the mean uncertainty of the
IRAC2 magnitudes. Out of the 2067 YSOVAR sources, 1807 satisfied this first
criterion, and among these, 1603 YSOs (corresponding to 77.6% of the total)
were classified into different variability types. Among the variable sources,
the majority were categorized as secular variables (85.5%), with the Curved
type (49.5%) being the most dominant among the YSO variability types.
Furthermore, a significant portion of the objects classified as secular
variables exhibited additional trends in their lightcurves, resulting in
classification as combined variables. Additional trends were observed in
67.6%, 70%, and 53.3% of the Linear, Curved, and Periodic YSOs, respectively.
Figure 5 displays example YSOVAR lightcurves for each short-term variable
type. These examples provide a visual representation of the different patterns
and behaviors exhibited by YSOs over many weeks.

Figure 5: Examples of YSOVAR lightcurves for each short-term variability type.
The source name is shown at the top of each panel, and each variable type is
indicated in the bottom left corner. Curved and Periodic are depicted with LSP
fittings shown as red lines.  Table 2: Short-term Variable Type

| Linear | Curved | Periodic | Burst | Drop | Irregular | Non-var type | Total
---|---|---|---|---|---|---|---|---
YSOVAR | 185 [9.0]aaVariable types. Y. refers to YSOVAR variable types | 793 [38.4] | 392 [24.5] | 5 [0.2] | 13 [0.6] | 215 [10.4] | 464 [22.4] | 2067 [100]
+PeriodicbbVariable types. N. refers to NEOWISE variable types | 64 | 357 | - | - | - | - | - | 421 (20.4)ccThe numbers in front are the count of total variables for each type, while the numbers in the square bracket are the fractions (%) of variables relative to their total samples.
+Burst | - | 1 | 1 | - | - | - | - | 2 (0.1)
+Drop | - | 2 | 1 | - | - | - | - | 3 (0.1)
+Irregular | 61 | 195 | 207 | - | - | - | - | 463 (22.4)

**footnotetext: The table note mark is same as Table 1.

####  3.2.1 Comparison between YSOVAR and this work

This section discusses the differences between the variability analysis
conducted in the YSOVAR study and the present study.

In Morales-Calderón et al. (2011), the authors classified the variability in
their targets using the Stetson index (Welch & Stetson, 1993) to search for
correlated variability in both IRAC channels. They also included YSOs with
determined periods, even if they did not exhibit high values of the Stetson
index. The classification into different types was supported by visual
inspection of the lightcurves and complemented by contemporaneous photometry
at shorter near-IR wavelengths. The variability in their sample was attributed
to physical mechanisms such as photospheric spots, changes in disk geometry,
variable extinction (AA Tau-like), and slow changes in the accretion rate.

In our study, we performed a statistical classification of the lightcurves
without relying on visual inspection. Our method aims to describe the shape of
the lightcurve without immediately associating the behavior with a specific
physical mechanism. Additionally, we did not use additional photometry at
shorter wavelengths. Figure 6 displays the fits and classification for six
sources from Figure 1 by Morales-Calderón et al. (2011), which were used to
illustrate the different types of variability identified in their work. Based
on our classification alone, it is impossible to discern between the different
mechanisms proposed by Morales-Calderón et al. (2011). However, all six YSOs
were classified as Periodic (or Periodic+Stochastic) in our analysis,
indicating variability at timescales shorter than 20 days, which is likely
associated with the phenomena in the stellar photosphere and inner disk, as
proposed by Morales-Calderón et al. (2011). Furthermore, we note that the
longer timescale trends observed in the YSOVAR study are revealed by the
Linear classification in our analysis.

Regarding the identification of PMS eclipsing binaries, Morales-Calderón et
al. (2011) reported nine cases, four of which were newly detected systems.
However, when applying our methods, we find it more challenging to distinguish
these four candidate PMS binaries compared to their results. Figure 7
illustrates that we are only able to classify one of their candidate YSOs as
Periodic, but when looking at the fitting, it is difficult to consider it as
an eclipsing binary.

When comparing all the periodic variables, however, the analysis of the
previous study and the present study yield very similar results. As shown in
Figure 8, the periods derived using our method align well with the periods
determined by Morales-Calderón et al. (2011), except for five out of 104
sources. Most periods lie along the one-to-one relation (blue line). We also
indicate objects that exhibit two periods in our analysis, classified as
Curved+Periodic. This represents one of the main differences from the work by
Morales-Calderón et al. (2011), as these additional trends were not identified
in their study.

Figure 6: Examples of YSOVAR lightcurves with their LSP fits. The gray line
represents the dominant variable trend while the red line includes a shorter
variable trend if such a period exists. These sources are the same as in
Figure 1 by Morales-Calderón et al. (2011). (a) and (b) Curved(gray line) +
Periodic sources (red line); (c) Periodic+Irregular source; (d)
Linear+Periodic source; (e) Periodic source; (f) Periodic + Irregular source.
Figure 7: Candidate PMS eclipsing binary lightcurves identified by Morales-
Calderón et al. (2011). These are sources listed in Table 5 of Morales-
Calderón et al. (2011) where they identified a stochastically repeating
variability. We are unable to reproduce their results by our analysis. Using
our criteria, these sources are classified as (a) Irregular, (b)
Curved+Periodic, (c) Periodic, and (d) Drop  Figure 8: Scatter plot comparing
the periods extracted from the YSOVAR data by the previous and current
studies. The x-axis represents the periods calculated in this study, while the
y-axis represents the periods calculated by Morales-Calderón et al. (2011).
The plotted data points represent sources classified as Periodic or +Periodic
as a combined variable. For sources classified as Curved+Periodic (15
sources), which exhibit two distinct periods corresponding to the Curved
classification and the +Periodic classification, they are depicted by pink
triangles and brown upside-down triangles, respectively. The two periods are
connected by gray dotted lines only for three sources as examples. Sources
classified solely as Periodic (84 sources) are denoted by blue dots, while the
+Periodic periods of sources classified as Linear+Periodic are represented by
olive dots (5 sources). The blue line represents the one-to-one relation
(y=x).

##  4 Relating SHORT- AND LONG-TERM VARIABILITY

To study the connection between short- and long-term, weeks and years,
variability in YSOs, we crossmatched the two catalogs described in the
previous sections.

###  4.1 Cross-matched sources

We conducted a cross-match between the NEOWISE (Section 3.1) and YSOVAR
(Section 3.2) catalogs using a matching radius of three arcseconds, resulting
in 703 positive matches. Table 3 provides the number of these YSOs classified
into different variability types based on both the short- and long-term
analyses. A significant portion of the sources (58.2%) do not have a NEOWISE
(long-term) variability classification. In contrast, the YSOVAR (short-term)
classifications reveal that the majority of sources (93.0%) are classified as
variable YSOs. This discrepancy can be attributed to the higher sensitivity of
the Spitzer telescope compared to the WISE telescope enabling the detection of
variability types with lower amplitudes. For the cross-matched sources, the
median magnitude uncertainty was 0.005 in YSOVAR and 0.04 in NEOWISE,
indicating that the sensitivity of YSOVAR is approximately 8 times better than
that of NEOWISE. Some YSOs likely exhibit high-frequency and low-amplitude
variability that remains relatively stable over time, without additional
mechanisms causing significant flux changes over longer observation periods.
As such, NEOWISE, with its poorer sensitivity, would not classify these YSOs
as variable stars in the long-term, whereas the YSOVAR data would still be
capable of capturing the variability in these systems.

The higher cadence of the YSOVAR observations provides an opportunity to study
the complexity of YSO variability in more detail. In the short-term analysis,
a higher percentage of YSOs are classified as combined variables compared to
the long-term analysis (see Table 4). Out of the 703 sources, 509 (72.4%) are
classified as combined variables, with the majority of them falling into the
Curved + Periodic class. On the other hand, in the long-term data, only 21
sources (3.0%) exhibit combined variable types, and they are predominantly
classified as Curved + Irregular.

We note that the Orion YSO catalog adopted by P21 lacks the class of PMS+E,
which denotes the evolutionary stage of Class III, possibly contaminated by
evolved stars. As a result, in the cross-matched sources (Table 3), the
‘periodic’ variables almost disappear; only 0.4% of cross-matched sources are
classified as periodic in long-term. The fraction of periodic sources
classified from all the NEOWISE sources is 2.9% as shown in Table 1 (refer to
P21 for additional discussion of the periodic classification).

Table 3: Variable Type : Cross-matched Sources | Y.LinearaaVariable types. Y. refers to YSOVAR variable types | Y.Curved | Y.Periodic | Y.Burst | Y.Drop | Y.Irregular | Y.Non-variable | NEOWISE Total
---|---|---|---|---|---|---|---|---
N.LinearbbVariable types. N. refers to NEOWISE variable types | 3ddThe number of sources having variable type linear in both YSOVAR(short-term) and NEOWISE(long-term). | 3 | 2 | - | - | 1 | - | 9 [1.3]ccThe numbers in front are the count of total variables for each type, while the numbers in the square bracket are the fractions (%) of variables relative to their total samples.
N.Curved | 3 | 36 | 10 | - | - | 5 | 1 | 55 [7.8]
N.Periodic | - | 2 | 1 | - | - | - | - | 3 [0.4]
N.Burst | 7 | 4 | - | - | - | 4 | 1 | 16 [2.3]
N.Drop | - | 3 | 1 | - | - | - | - | 4 [0.6]
N.Irregular | 35 | 112 | 53 | - | - | 5 | 2 | 207 [29.4]
N.Non-variable | 37 | 194 | 84 | - | 1 | 48 | 45 | 409 [58.2]
YSOVAR Total | 85 [12.1] | 354 [50.4] | 151 [21.5] | 0 [0.0] | 1 [0.1] | 63 [9.0] | 49 [7.0] | 703 [100]
Table 4: Combined Variable Type : Cross-matched Sources |  | YSOVAR |  |  |  | NEOWISE |  |
---|---|---|---|---|---|---|---|---
| Y.Linear | Y.Curved | Y.Periodic | Y.Total | N.Linear | N.Curved | N.Periodic | N.Total
$+$Periodic | 39 | 190 | - | 229 [45.0] | - | 1 | - | 1 [4.8]
$+$Burst | - | - | 1 | 1 [0.2] | - | 1 | - | 1 [4.8]
$+$Drop | - | - | 1 | 1 [0.2] | 1 | - | - | 1 [4.8]
$+$Irregular | 37 | 115 | 126 | 278 [54.5] | 4 | 14 | - | 18 [85.7]
Total | 76 [14.9] | 305 [59.9] | 128 [25.1] | 509 (72.4)aaThe numbers in front are the total count of combined variables, while the numbers in the parentheses are the fractions (%) of combined variables relative to their total samples. | 5 [23.8] | 16 [76.2] | 0 [0.0] | 21 (3.0)aaThe numbers in front are the total count of combined variables, while the numbers in the parentheses are the fractions (%) of combined variables relative to their total samples.

###  4.2 Comparison of the variabilities

To compare the degree of variability between the short-term and long-term
observations, we analyzed two parameters: the Normalized Lomb-Scargle
Periodogram amplitude ($N_{amp}$) and the normalized flux root-mean-square
value ($Nf_{rms}$). To calculate $N_{amp}$, we selected the 49 YSO sources
that were classified as either Curved or Periodic in both the short-term and
long-term analyses. For each of these sources, we obtained the long-term and
short-term amplitudes from the Lomb-Scargle periodogram fitting and normalized
them by dividing them by the mean flux determined during the fitting (see
Figure 9).

The comparison of $N_{amp}$ between the short-term and long-term is presented
in Figure 10. The scatter plot illustrates that the sources tend to lie to the
right of the $y=x$ line, indicating that the amplitude of long-term
variability is generally larger than that of short-term variability. The
histogram (right of Figure 10) shows the distribution of the amplitude ratio
between NEOWISE and YSOVAR. All sources have ratios above 1, clearly
signifying that long-term amplitude exceeds short-term amplitude, with a peak
at $\sim$3\. The largest ratio goes up to 10. This difference in amplitude can
be attributed to the different underlying mechanisms driving the short-term
and long-term changes. Short-term variability is primarily influenced by
physical mechanisms related to the stellar photosphere or inner disk, which
typically result in lower amplitudes. On the other hand, long-term changes are
often driven by mechanisms involving variations in geometry and accretion
rate, which can produce larger amplitude variations.

Figure 9: An example of calculating $N_{amp}$ for the YSO source
ISOY_J053522.36-050739.2 (also known as ’M2341’ in the NEOWISE catalog). The
top panel shows the YSOVAR lightcurve, while the bottom panel displays the
NEOWISE lightcurve. The red curved lines represent the Lomb-Scargle
periodogram fitting for each dataset, with the corresponding amplitudes
indicated by red arrows. The green horizontal lines represent the mean flux of
the fitting in each dataset. The $N_{amp}$ values were calculated as the ratio
of the LSP amplitude to the mean flux of the fitting.  Figure 10: (Left)
Scatter plot of Normalized Amplitude ($N_{amp}$) for 49 sources classified as
Curved or Periodic in both long-term and short-term analyses. In the legend,
”Y.” denotes YSOVAR, and ”N.” denotes NEOWISE. For instance, ”Y.curved &
N.periodic” signifies that a source is classified as Curved in YSOVAR data and
Periodic in NEOWISE data. (Right) Histogram illustrating the ratio of NEOWISE
normalized amplitude to YSOVAR normalized amplitude. All values are greater
than 1, indicating that long-term amplitude exceeds short-term amplitude. The
peak is centered around 3, and the distribution extends up to 10.

In a second method, we analyzed the “normalized flux root-mean-square value
($Nf_{rms}$)” for all 703 sources. For each source, we calculated the flux rms
(root-mean-square) value from the long- and short-term lightcurves and for
normalization, we divided it by the median flux. Figure 11 presents the
scatter plot and distribution histograms of $Nf_{rms}$ for the cross-matched
sources.

The leftmost panel of Figure 11 shows the scatter plot of the Normalized Flux
rms ($Nf_{rms}$) for all 703 cross-matched sources. The top two histogram
(panels (b) and (d)) illustrates the $Nf_{rms}$ values for YSOVAR variability
types, revealing that the distribution peak of secular variables (red) is at
larger values compared to stochastic (blue) and non-variable (black) YSOs.
This suggests that (quasi-)periodic variability related to hot spots or
obscuration events dominates among the YSO variable population during the
40-day observations of the YSOVAR program. Interestingly, as these same
sources have the highest $Nf_{rms}$ values in the NEOWISE data (panel (b)),
they are also associated with long-term (NEOWISE) secular (yellow) and
stochastic (green) variability (compare with panel (c)).

In the bottom panel (panels (c) and (e)), focusing on the long-term
classification, secular variables (yellow) exhibit larger $Nf_{rms}$ values in
the NEOWISE data ($N.Nf_{rms}$) compared to the YSOVAR data ($Y.Nf_{rms}$).
This indicates that YSOs dominated by long-term, high-amplitude variability,
potentially associated with significant changes in accretion or geometrical
variations in the disk, tend to have relatively constant or minimally varying
brightness over a few weeks of observations. As further support, the long-term
non-variables (black) peak at roughly the same $Nf_{rms}$ value for both data
sets. Notably, the distributions of long-term secular and long-term non-
varying sources are almost indistinguishable in panel (e). This is not the
case for long-term stochastic sources (green). These sources exhibit less of
decrease in $Nf_{rms}$ values between the NEOWISE and YSOVAR datasets and show
the largest $Nf_{rms}$ values in panel (e).

As noted above, the distribution peaks for NEOWISE stochastic variables
(green) align with the $Nf_{rms}$ values of secular variables (red) classified
from YSOVAR data (upper panel). This suggests that many of the objects
displaying secular variability over 40 days of observations exhibit a
consistent level of variability over years-long observations. Due to the
cadence of NEOWISE, these YSOs would be observed at random epochs in their
lightcurves, and their short-term variability would likely appear irregular
over the long term, leading to their classification as stochastic.

Figure 11 also highlights the higher sensitivity of Spitzer observations. In
the region where $N.Nf_{rms}<$ 0.03 (panels (b) and (c)), many YSOs classified
as secular or stochastic based on short-term data are categorized as non-
variable in the long-term data analysis.

Figure 11: (a) Scatter plot depicting the Normalized Flux rms ($Nf_{rms}$) for
all 703 cross-matched sources. Here, ’Y.’ denotes YSOVAR, and ’N.’ represents
NEOWISE. The ’sec’ category encompasses secular variables, while ’sto’
corresponds to stochastic variables. ’non’ signifies the Non-variable type.
(b) and (c) present histograms of NEOWISE normalized flux rms. (b) illustrates
the classification results based on YSOVAR, while (c) showcases the
classification results based on NEOWISE. Within the region of $N.Nf_{rms}<$
0.03, some sources that were classified as secular or stochastic in the short-
term data are classified as non-variable in the long-term data. (d) and (e)
display histograms of YSOVAR normalized flux rms. (d) demonstrates the
classification outcomes based on YSOVAR, whereas (e) illustrates the
classification outcomes based on NEOWISE. The histogram’s peak locations
reveal that sources classified as secular in the short-term exhibit the
highest variability at both of these timescales.

###  4.3 Degrading the sensitivity of YSOVAR observations

To assess the potential impact of lower NEOWISE sensitivity on comparing
short- and long-term variability in YSOs, we degraded the sensitivity of the
YSOVAR data to examine the effect on detecting variability. In the previous
section, we established that most YSOs exhibit short-term variability and
suggested that this could be attributed to the superior sensitivity of the
Spitzer observations compared to WISE. To verify this hypothesis, we applied
the mean uncertainty of NEOWISE W2 flux to the YSOVAR flux for each source
when determining YSOs as a ‘variable’.

Table 5 compares the use of two different sensitivities. ‘Y.eflux’ denotes the
results obtained using the original YSOVAR flux errors, while ‘N.eflux’
represents the outcomes obtained by degrading the sensitivity of the YSOVAR
data to match the NEOWISE sensitivity. When determining short-term variability
using the original YSOVAR flux errors, 81.9% of YSOs were classified as
variables. However, when employing the WISE flux errors, 15.1% were identified
as variables. As expected, the increased value of $\sigma$ in the YSOVAR data
fails to capture the low-amplitude variability. Therefore, the telescope’s
sensitivity is critical in detecting YSO variability, suggesting that the
NEOWISE lightcurves miss a majority of variable YSOs with low amplitudes, even
in the long-term.

Table 5: Variable detection by different telescope  Short-term data | Y.efluxbb$\sigma$ = $\sigma$(IRAC2 flux unc) ; the mean uncertainty of IRAC2 flux. | N.efluxcc$\sigma$ = $\sigma$(W2 flux unc) ; the mean uncertainty of W2 flux.
---|---|---
VariableaaVariable = SD/$\sigma\textgreater 3$ | 576 [81.9] | 106 [15.1]
Non variable | 127 [18.1] | 597 [84.9]
Total | 703 [100] | 703 [100]ddThe numbers in the parentheses are the fractions (%) of variables relative to their total 703 samples.

##  5 Discussion

###  5.1 Variability Mechanisms

YSO variability can arise from various processes, including accretion,
extinction, changes in disk properties, interactions between the disk and
stellar photosphere, and hot$/$cold magnetic spots in the stellar photosphere
(see e.g., Herbst et al., 1994; Carpenter et al., 2001; Stauffer et al., 2014;
Contreras Peña et al., 2017; Cody et al., 2017; Park et al., 2021). Comparing
the short-term and long-term variability analyses from YSOVAR and NEOWISE
provides insights into the mechanisms driving brightness changes in YSOs.

Analyzing weeks-long YSOVAR data alone leads to similar conclusions as those
derived by Morales-Calderón et al. (2011). The majority of YSOs in the ONC
sample (93.0% among cross-matched YSOs in our analysis) exhibit some degree of
variability. The timescales of variability are related to processes affecting
the stellar photosphere-inner disk regions, such as short-term bursts caused
by magnetic instabilities and reconnection between the stellar magnetosphere
and the inner disk edge or periodic obscuration events similar to AA Tau-like
extinction ($P<15$ days, Bouvier et al., 2013; Cody et al., 2014). Within the
ONC sample, a significant fraction of YSOs (28.4%) are classified as secular
variables in the short-term analysis, which display stochastic ($Irregular$)
variability in the long-term NEOWISE data. The amplitude of short-term and
long-term variability in these objects shows similar values. This suggests
that the mechanisms driving variability in the short term also persist over
the long-term, with the classification as Irregular arising from the cadence
of NEOWISE observations. 51.8% of the sample exhibit short-term variability
but are classified as non-variable in the long-term. These objects tend to
have lower amplitudes, indicating that we are likely observing the same
variability mechanisms over the long-term but cannot classify them due to the
lower sensitivity of WISE observations.

YSOs classified as secular variables in the years-long NEOWISE data generally
exhibit larger amplitudes in the long term compared to the short-term
observations. In these cases (approximately 9.7% of the sample), the
variability may be attributed to structural changes in the disk such as
changes in the scale height of the inner disk (Covey et al., 2021), long-term
obscuration caused by a dusty wind (Davies et al., 2018), or variations in the
accretion rate (Fischer et al., 2022). Evidence for accretion variations
driving mid-IR variability comes from a direct comparison of sub-mm
variability and mid-IR variability in the same source, as shown by Contreras
Peña et al. (2020). Figure 12 illustrates the example of V418 Ori, a YSO
displaying Curved+Periodic short-term variability and long-term changes with
higher amplitudes classified as Linear+Irregular.

Figure 12: An example lightcurve of Disk source (V418 Ori). (a) The overall
lightcurves, including IRAC1 (green), IRAC2 (yellow), W1 (red), and W2 (blue)
bands, are represented on a magnitude scale. W1 and W2 data at $\sim$ MJD
55200 are AllWISE data. (b) The lightcurve specifically for the IRAC2 band,
showing short-term variability. This source is classified as Curved+Periodic
in the short-term analysis, with periods of 37.39 days and 6.92 days. Gray and
red lines represent the corresponding Lomb-Scargle Periodogram (LSP) fitting
for each period, with the y-axis representing the flux scale. (c) The W2
NEOWISE lightcurve, revealing long-term variability. This target is classified
as Linear+Irregular. The orange line represents the Linear Least Square
fitting for this classification, with the y-axis representing the normalized
flux, calculated as (flux - median flux) $/$ median flux.

###  5.2 Very short-term variability

The NEOWISE single-epoch catalog offers observations with durations mostly
ranging from 1 to 2 days, repeated every 6 months (see Section 2.1). In this
section, we focus on analyzing the variability at these timescales, which we
refer to as the very short-term.

There are limitations in analyzing very short-term data due to its brief
duration and a 6-month time gap between observations. Therefore, it becomes
challenging to determine whether the observed variability in the very short-
term data will persist beyond that specific time period. Additionally, due to
the limited number of data points within each NEOWISE epoch, it is crucial to
ensure the inclusion of only high-quality data. Consequently, we applied
additional criteria to select NEOWISE high-quality data, extending beyond the
considerations made in the analysis by P21 (see Section 2.1). With these
carefully chosen high-quality exposures, we used the individual exposures,
which were averaged for the long-term analysis, to examine the very short-term
variability.

To investigate the behavior of very short-term variability, we visually
examined the lightcurves to compare the observed variability with the expected
trends derived from fitting the short- and long-term variations. In cases
where a YSO was classified as Curved or Periodic based on the YSOVAR data
(representing short-term variability) and defined as a secular variable in the
NEOWISE data (representing long-term variability), we combined the YSOVAR
Lomb-Scargle Periodogram (LSP) fitting with the NEOWISE long-term fitting.
Subsequently, we examined whether the 1-2 day NEOWISE observations
(representing very short-term variability) followed this combined fitting
trend.

For instance, Figure 13 compares the observed and expected variability for the
YSO V418 Ori. In the NEOWISE 7.5-year data, this source is classified as
Linear+Irregular, whereas in the YSOVAR 40-day timescale data, it is
categorized as Curved+Periodic. The bottom panel of the figure illustrates
that the very short-term variability aligns well with the combined fitting
trend. However, certain epochs exhibit discrepancies, which could arise from
stochastic variability processes occurring at very short timescales, such as
bursts resulting from unstable accretion (see e.g., Cody et al., 2017).
Similar stochastic behavior can also be observed during the YSOVAR
observations (refer to Figure 12).

To assess the level of YSO variability across different timescales, we
examined the variability amplitude ($\Delta$) at IRAC2 ($4.5\mu$m) and W2
($4.6\mu$m) for the very short-term, short-term, and long-term periods. In the
case of the very short-term analysis, we calculated the flux standard
deviation (SD) within a specific epoch lasting 1-2 days. In this study, we
chose one epoch, as an illustrative example, which occurred at MJD $\sim$
57810, since data for all 703 sources that were cross-matched between the
NEOWISE and YSOVAR catalogs (refer to Section 4.1) were available during this
epoch. We conducted similar analyses for other very short-term epochs, and all
yielded comparable results to those observed in this epoch. The short-term SD
was determined based on the 40-day YSOVAR observations, while the long-term SD
was derived from the 7.5-year NEOWISE data. To ensure the same sensitivity to
variable identification, we adopt the mean flux uncertainty calculated from
the NEOWISE data for all three timescales.

Figure 14 compares distributions of $\Delta$magnitude of variables in the
three different timescales. Applying the same mean flux uncertainty of NEOWISE
($\sigma(W2)$) to both YSOVAR and NEOWISE data of the 703 cross-matched YSOs,
variables were identified with the criteria of $SD/\sigma(W2)>3$ in three
different timescales: 52 very short-term variables from the NEOWISE
measurements within the epoch of MJD $\sim$ 57810, 106 short-term variables
from the YSOVAR data, and 264 long-term variables from the epoch-averaged
NEOWISE data.

A pair-wise Kolmogorov-Smirnov (KS) test, utilizing the scipy.stats.kstest
function from the Python scipy package (Virtanen et al., 2020) was performed
between the very-short/short and short/long results to determine if the
samples were drawn from a single distribution. In both cases this null
hypothesis was definitely rejected ($p<<0.05$), arguing that the three
observed distributions are independent. The histogram illustrates an increase
in the value of $\Delta$ as the observational baseline extends, aligning with
the findings of the previous section and previous YSO variability studies
(e.g., Scholz, 2012; Contreras Peña et al., 2019). Furthermore, a modest
increase in $\Delta$ between very short-term and short-term variability
exists. The change, however, is very small with respect to the breadths of
each distribution, suggesting that mechanisms driving 1-2 day variability are
prevalent among YSOs.

Figure 13: (Top) An example of lightcurve with YSOVAR IRAC2 (yellow), NEOWISE
W2 data (blue) and the combined fitting of long-term with short-term (red
line), corresponding to the same source as shown in Figure 12. (Bottom)
NEOWISE very short-term lightcurve, where the data points with an interval of
approximately 6 months are removed. The gray boxes indicate the interval time,
marked at the top of each box, with ’186d’ indicating 186 days. The red line
represents the combined fitting, which is the same as indicated in the top
panel.  Figure 14: Histogram of $\Delta$ magnitude (Max - Min) at different
timescales with Kernel Density Estimation (KDE) fitting. Each color
corresponds to a specific timescale, highlighting the observed variability
within each period. The histogram includes data for 264 long-term variables
(blue), 106 short-term variables (yellow), and 52 very short-term variables
(pink). The KDE fitting displays the peak value for each timescale, which are
0.3, 0.4, and 0.6, respectively. This indicates that the magnitude difference
($\Delta$ magnitude) tends to increase with longer timescales.

###  5.3 Interesting Individual sources

This section focuses on specific sources that exhibit intriguing lightcurve
behavior across the various timescales analyzed in this study.

One such source is V935 Ori (Figure 15). Across the 40-day YSOVAR lightcurve
it demonstrates a Linear + Periodic trend with a combined period of 13.72
days. A burst-like event is observed around MJD $\sim$ 55150 and lasting for a
few days, although it is not classified as a Burst source due to the
relatively small amplitude of the burst. In the NEOWISE long-term lightcurve,
a more pronounced burst is evident, leading to its classification as a Burst
source. Additionally, the AllWISE data at approximately MJD 55200 also
captures a burst event. These burst events are likely associated with unstable
accretion resulting from the interaction between the stellar magnetosphere and
the inner edge of the disk (Stauffer et al., 2014; Cody et al., 2017).

Figure 15: V935 Ori lightcurves in various timescales. It shows
Linear+Periodic with little burst in 40-day analysis (IRAC2, orange data) and
Burst in the long-term (W2, blue data).

Another notable source, V2528 Ori (Figure 16), exhibits an Irregular type in
the long-term analysis, indicating irregular fluctuations over an extended
period. However, in the 40-day analysis, it shows a Linear+Irregular type,
suggesting that the irregular fluctuations occur on top of a linear trend. The
very short-term variability appears similar to the irregular fluctuations
observed in the YSOVAR data. The complexity of the lightcurve in V2528 Ori is
likely attributed to a combination of physical mechanisms driving its
variability. The short-term irregularity, which appears nearly periodic on
timescales exceeding a few days based on visual inspection, may be influenced
by the presence of (hot) stellar spots. On the other hand, the dominant linear
trend observed in the long-term analysis, as evidenced by the long period of
approximately 160 days derived from the YSOVAR data, is most likely driven by
structural changes occurring in the accretion disk of the YSO.

Figure 16: V2528 Ori lightcurves in various timescales. It shows
Linear+Irregular type in 40-day analysis (IRAC2, orange data) and Irregu lar
type in the long-term (W2, blue data).

A last example is V2012 Ori (JW 61) (Figure 17). In the long-term analysis,
this source does not exhibit any variability. However, in the YSOVAR 40-day
analysis, it is classified as Periodic, displaying a periodic behavior with a
period of 6.25 days. The very short-term data further confirms the persistence
of this periodicity throughout the NEOWISE coverage of the source. In this
case, the short-term periodic variability is most likely attributed to
magnetic spots on the stellar photosphere. These spots, possibly associated
with magnetic activity, give rise to the observed periodic variations in the
lightcurve. However, it is worth noting that V2012 Ori is a binary system
(Reipurth et al., 2007), which could potentially contribute to its
variability.

Figure 17: V2013 Ori lightcurves in various timescales. It shows
Periodic+Irregular type with a period of 6.25 days in 40-day analysis (IRAC2,
orange data), but it is non-variable in the long-term (W2, blue data).

By examining YSO variability across these different timescales, we can gain a
more comprehensive understanding of the physical processes driving the
observed variability. This multi-timescale approach helps us to distinguish
between different types of variability, identify correlations between short-
term and long-term behaviors, and unravel the complex interplay of various
mechanisms influencing YSOs’ brightness variations.

##  6 SUMMARY

Variability in YSOs is a common phenomenon that can be caused by various
physical mechanisms affecting the stellar photosphere, accretion disk and
envelope in these systems. To gain a better understanding of these mechanisms
it is important to examine YSO variability on different timescales. In
particular, mid-IR data sets are powerful for evaluating protostellar
accretion as well as other physical mechanisms driving YSO variability.

In this study, we utilize mid-IR data arising from the NEOWISE and YSOVAR
surveys, which cover the wavelength range of 3 to 5 $\mu$m over 7.5 years and
40-day timescales respectively. From analysis of these data we found.

  1. 1.

The majority of YSOs exhibit variability on the 40-day timescales provided by
YSOVAR (77.6%). This variability is primarily driven by YSOs classified as
secular variables, i.e objects that display variability over timescales larger
than 20 days. On the other hand, less than half of the YSOs (43.0%) display
variability on the 7.5 year timescales from NEOWISE, primarily attributed to
objects classified as stochastic variables (Section 3).

  2. 2.

We found 703 YSOs in Orion which contain data sets from both surveys. Among
these YSOs, a higher percentage (93.0%) exhibited variability in the short-
term, while a lower percentage (41.8%) showed variability in the long-term
(Section 4.1). Additionally, the amplitude of variability in the short-term
was smaller than in the long-term, suggesting that the stellar photosphere-
inner disk interface causes the short-term variability, while accretion or
geometrical variations in the disk result in the long-term variability
(Section 4.2).

  3. 3.

We studied the possible effect that the lower sensitivity of NEOWISE data,
compared with YSOVAR observations, could have on our results. When we degraded
the YSOVAR data using the NEOWISE magnitude error, the fraction of variable
sources significantly decreased from 81.9% to 15.1%. This suggests that the
sensitivity of the telescope is crucial in detecting YSO variability (Section
4.3).

  4. 4.

Finally, we studied the YSO variability in a very short-term timescale (1-2
days) provided by the single-epoch NEOWISE observations. We find that the
variability amplitude ($\Delta$) increases as the observational baseline
extends (Section 5.2), in agreement with previous studies of the variability
of YSOs.

##  Appendix A NEOWISE Outlier data cut

In our analysis of the 7.5-year NEOWISE data, we have adopted a different
approach compared to the 6.5-year NEOWISE data analysis conducted in P21. In
P21, they utilized only 70% of the exposures in each epoch to mitigate
potential errors in the photometry. However, considering that most YSOs
exhibit variability on short timescales (Morales-Calderón et al., 2011),
removing 30% of the NEOWISE data as done in P21 could potentially hinder our
ability to detect very short-term variability within a single NEOWISE epoch.

To address this concern, we have employed an alternative method to select
high-quality exposures based on the photometry parameters provided in the
NEOWISE database. Specifically, we have considered parameters such as the
signal-to-noise ratio (SNR) of each exposure and the distance to the
coordinates of the YSO (see Section 2.1). By applying these selection
criteria, we ensure that the data used for analyzing very short-term
variability are of sufficient quality.

Figure 18 illustrates two lightcurves, comparing the results obtained by
removing 30% of the NEOWISE data versus our approach of selecting high-quality
exposures. We have observed that in most sources, removing the 30% of data
would result in the loss of data points that exhibit genuine variability at
very short timescales.

Therefore, instead of discarding a fixed percentage of data, we opted for a
data selection strategy based on the photometry parameters in the NEOWISE
database. This approach allows us to retain valuable information regarding
very short-term variability while ensuring the reliability and quality of the
data used in our analysis.

Figure 18: Examples of very short-term lightcurve with marked 30% Outlier
(red). (Top) Lightcurve of a source showing secular trends in the very short-
term variable. Most of the variables have this trend. (Bottom) Lightcurve of a
source showing stochastic trends in the very short-term variable.

## Acknowledgement

This work was supported by the National Research Foundation of Korea (NRF)
grant funded by the Korea government (MSIT) (grant number 2021R1A2C1011718).
G.J.H. is supported by general grants 12173003 awarded by the National Natural
Science Foundation of China. D.J. is supported by NRC Canada and by an NSERC
Discovery Grant.

This publication makes use of data products from NEOWISE, which is a project
of the Jet Propulsion Laboratory/California Institute of Technology, funded by
the Planetary Science Division of the National Aeronautics and Space
Administration.

## References

  * Attridge & Herbst (1992) Attridge, J. M., & Herbst, W. 1992, ApJ, 398, L61, doi: 10.1086/186577
  * Bino et al. (2023) Bino, G., Basu, S., Dey, R., et al. 2023, arXiv e-prints, arXiv:2302.03742, doi: 10.48550/arXiv.2302.03742
  * Bouvier et al. (2013) Bouvier, J., Grankin, K., Ellerbroek, L. E., Bouy, H., & Barrado, D. 2013, A&A, 557, A77, doi: 10.1051/0004-6361/201321389
  * Cargile et al. (2008) Cargile, P. A., Stassun, K. G., & Mathieu, R. D. 2008, ApJ, 674, 329, doi: 10.1086/524346
  * Carpenter et al. (2001) Carpenter, J. M., Hillenbrand, L. A., & Skrutskie, M. F. 2001, AJ, 121, 3160, doi: 10.1086/321086
  * Choi & Herbst (1996) Choi, P. I., & Herbst, W. 1996, AJ, 111, 283, doi: 10.1086/117780
  * Cody et al. (2017) Cody, A. M., Hillenbrand, L. A., David, T. J., et al. 2017, ApJ, 836, 41, doi: 10.3847/1538-4357/836/1/41
  * Cody et al. (2014) Cody, A. M., Stauffer, J., Baglin, A., et al. 2014, AJ, 147, 82, doi: 10.1088/0004-6256/147/4/82
  * Contreras Peña et al. (2020) Contreras Peña, C., Johnstone, D., Baek, G., et al. 2020, MNRAS, 495, 3614, doi: 10.1093/mnras/staa1254
  * Contreras Peña et al. (2019) Contreras Peña, C., Naylor, T., & Morrell, S. 2019, MNRAS, 486, 4590, doi: 10.1093/mnras/stz1019
  * Contreras Peña et al. (2017) Contreras Peña, C., Lucas, P. W., Minniti, D., et al. 2017, MNRAS, 465, 3011, doi: 10.1093/mnras/stw2801
  * Covey et al. (2021) Covey, K. R., Larson, K. A., Herczeg, G. J., & Manara, C. F. 2021, AJ, 161, 61, doi: 10.3847/1538-3881/abcc73
  * Davies et al. (2018) Davies, C. L., Kreplin, A., Kluska, J., Hone, E., & Kraus, S. 2018, MNRAS, 474, 5406, doi: 10.1093/mnras/stx3150
  * Dunham et al. (2015) Dunham, M. M., Allen, L. E., Evans, Neal J., I., et al. 2015, ApJS, 220, 11, doi: 10.1088/0067-0049/220/1/11
  * Esplin & Luhman (2019) Esplin, T. L., & Luhman, K. L. 2019, AJ, 158, 54, doi: 10.3847/1538-3881/ab2594
  * Fischer et al. (2022) Fischer, W. J., Hillenbrand, L. A., Herczeg, G. J., et al. 2022, arXiv e-prints, arXiv:2203.11257, doi: 10.48550/arXiv.2203.11257
  * Herbst et al. (1994) Herbst, W., Herbst, D. K., Grossman, E. J., & Weinstein, D. 1994, AJ, 108, 1906, doi: 10.1086/117204
  * Lee et al. (2021) Lee, Y.-H., Johnstone, D., Lee, J.-E., et al. 2021, ApJ, 920, 119, doi: 10.3847/1538-4357/ac1679
  * Lomb (1976) Lomb, N. R. 1976, Ap&SS, 39, 447, doi: 10.1007/BF00648343
  * Mainzer et al. (2014) Mainzer, A., Bauer, J., Cutri, R. M., et al. 2014, ApJ, 792, 30, doi: 10.1088/0004-637X/792/1/30
  * Megeath et al. (2012) Megeath, S. T., Gutermuth, R., Muzerolle, J., et al. 2012, AJ, 144, 192, doi: 10.1088/0004-6256/144/6/192
  * Morales-Calderón et al. (2011) Morales-Calderón, M., Stauffer, J. R., Hillenbrand, L. A., et al. 2011, ApJ, 733, 50, doi: 10.1088/0004-637X/733/1/50
  * NEOWISE Team (2020) NEOWISE Team. 2020, NEOWISE-R Single Exposure (L1b) Frame Metadata Table, IPAC, doi: 10.26131/IRSA143
  * Park et al. (2021) Park, W., Lee, J.-E., Contreras Peña, C., et al. 2021, ApJ, 920, 132, doi: 10.3847/1538-4357/ac1745
  * Reipurth et al. (2007) Reipurth, B., Guimarães, M. M., Connelley, M. S., & Bally, J. 2007, AJ, 134, 2272, doi: 10.1086/523596
  * Romanova et al. (2013) Romanova, M. M., Ustyugova, G. V., Koldoba, A. V., & Lovelace, R. V. E. 2013, MNRAS, 430, 699, doi: 10.1093/mnras/sts670
  * Scargle (1989) Scargle, J. D. 1989, ApJ, 343, 874, doi: 10.1086/167757
  * Scholz (2012) Scholz, A. 2012, MNRAS, 420, 1495, doi: 10.1111/j.1365-2966.2011.20136.x
  * Stassun et al. (2006) Stassun, K. G., Mathieu, R. D., & Valenti, J. A. 2006, Nature, 440, 311, doi: 10.1038/nature04570
  * Stassun et al. (2004) Stassun, K. G., Mathieu, R. D., Vaz, L. P. R., Stroud, N., & Vrba, F. J. 2004, ApJS, 151, 357, doi: 10.1086/382353
  * Stauffer et al. (2014) Stauffer, J., Cody, A. M., Baglin, A., et al. 2014, AJ, 147, 83, doi: 10.1088/0004-6256/147/4/83
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Welch & Stetson (1993) Welch, D. L., & Stetson, P. B. 1993, AJ, 105, 1813, doi: 10.1086/116556
  * Wolk et al. (2018) Wolk, S. J., Günther, H. M., Poppenhaeger, K., et al. 2018, AJ, 155, 99, doi: 10.3847/1538-3881/aaa6c4