Quantifying day-to-day variability of O/N2 and its correlation with geomagnetic activity using GOLD

We quantify the short-term (<30 day) variability of column O/N2 measured by GOLD from January 2019 to August 2022 for various geomagnetic activity conditions. We find enhanced variabilities at high latitudes during active (Kp ≥ 3.0) times and weak but statistically significant variabilities at low latitudes. For active times, the largest absolute variability of O/N2 ratio is 0.14 and the largest relative variability is 20.6% at ∼60.0°N in Fall, which are about twice those of quiet times. The variability at higher latitudes can be larger than that of lower latitudes by a factor of 5–8. We further quantify contributions of magnetospheric forcing to O/N2 variability in the Ionosphere-Thermosphere region by correlating O/N2 perturbations with Dst. During geomagnetic active times, positive correlations as large as +0.66 and negative correlations as large as −0.65 are found at high and low latitudes, respectively, indicative of storm-induced O and N2 upwelling at high latitudes and down welling at low latitudes. During quiet times, correlations between O/N2 perturbations and Dst become insignificant at all latitudes, implying a more substantial contribution from below. O/N2 variabilities maximize in Fall and decrease towards Summer, while correlations maximize in Spring/Summer and decrease in Winter/Spring, which may be related to seasonal variations of geomagnetic activity and mean circulation.


Introduction
The Ionosphere-Thermosphere (IT) environment is strongly coupled to the magnetosphere above and the regions of the lower atmosphere below (Rishbeth et al., 2000). During active times, the short-term variability of the IT region is dominated by forcing from above in the form of enhanced magnetospheric forcing related to geomagnetic storms and substorms. During quiet times, the short-term variability of the IT region is modulated primarily by in-situ generated tides and/or upward-propagating waves originating from the lower atmosphere (Richmond, 1979;Rishbeth, 1997;Immel et al., 2006;Oberheide et al., 2011;Oberheide et al., 2013). Even during quiet times, however, solar geomagnetic activity can still have a significant impact on IT preconditioning (Cai et al., 2020), and it can take as long as 2 days for the IT to return to pre-storm conditions after a strong or moderate geomagnetic disturbance (Richmond et al., 2003;Lu et al., 2015;Yao et al., 2016). This makes the accurate characterization of the day-to-day variability difficult because it requires that we disentangle the effects of the different forcings in order to understand the relative contributions of each.
Strong agreement between the observed climatology of the IT region and the modeled climatology [from models like the Thermosphere-Ionosphere-Electrodynamics General Circulation Model (TIEGCM), its predecessor the Thermosphere General Circulation Model (TGCM), and the Whole Atmosphere Community Climate Model with thermosphere and ionosphere eXtension (WACCM-X)] has been demonstrated since as early as 1981 (Dickinson et al., 1981;Richmond et al., 1992;Fesen et al., 2002;Qian et al., 2014;Liu et al., 2018), but the day-to-day variations are comparatively difficult to model/study and therefore less understood (Zhou et al., 2020, Zhou et al., 2021Andoh et al., 2022). Many details of the state of the IT system can be captured by parameters such as the ratio of atomic Oxygen (O) to molecular Nitrogen (N 2 ) (O/N 2 ). The seasonal O/N 2 composition of the IT region depends, in particular, on global, vertical circulation in the ionosphere (Burns et al., 1995;Rishbeth et al., 2000;Qian et al., 2010;Luan et al., 2017;Li et al., 2022), the absorption of solar EUV during geomagnetic active periods, which changes the ratio via a combination of photodissociation and ionization processes (Kelley, 2009;Schmölter et al., 2020;Yu et al., 2020), and auroral heating-induced vertical winds which bring up atomic oxygen-poor air from the lower atmosphere (Hays et al., 1973).
The recently launched Global Observations of the Limb and Disk (GOLD) satellite provides full-disk ionospheric column O/N 2 as a Level 2 data product from a geostationary orbit centered over the mouth of the Amazon River basin (Eastes et al., 2017). GOLD is the first satellite to observe the IT region from a fixed, geostationary orbit, which allows for the unambiguous separation of spatial and temporal IT variability and therefore makes possible unprecedented analysis of the day-to-day variability of O/N 2 . Recent studies have utilized GOLD O/N 2 and temperature measurements to study seasonal variations of O/N 2 as well as the response of O/N 2 to both strong and weak geomagnetic activity (Cai et al., 2020;Laskar et al., 2021;Qian et al., 2022), but the seasonal, short-term response to varying geomagnetic activity levels has not yet been studied systematically across the full, presently available GOLD O/N 2 dataset. Efforts to assimilate GOLD disk temperatures into whole atmosphere models have already led to significant improvements in model effective temperature root-mean square error (RMSE) (Laskar et al., 2021), and it should be possible to use GOLD O/N 2 in a similar fashion to improve the abilities of models to capture day-to-day IT variability.
2 Data and methodology 2.1 GOLD data and data selection GOLD measures the ratio of the vertical column density of O relative N 2 defined at a standard reference N 2 depth of 10 17 cm −2 . This reference depth is chosen to minimize uncertainty in O/N 2 as a function of 135.6 nm/LBH and is equivalent to~4 nanobar, or~160 km (Correira et al., 2021). The O/N 2 ratio is retrieved directly from dayside disk measurements of the OI-135.6 nm and N2-LBH band emission lines in the~134-164 nm wavelength range. We use version 3 of the O/N 2 data product. The GOLD O/N 2 data product is not optimized for auroral latitudes, so we limit our analysis to ±60°lat. The observations have longitudinal coverage ranging from 113.5°W to 19.5°E.
The O/N 2 data from January 2019 to August 2022 and universal time (UTC) =~14:00-15:00 are used for the current study. This UTC time range is chosen because the full footprint of GOLD is illuminated at this time and because GOLD scanned at a cadence of one full disk every half hour from 2019-2021, and every hour in 2022. We omit several days in which the data are unreliable for instrumental reasons. We omit the data from the days surrounding a Gyrating Yaw Mechanism (GYM) actuation. GYM actuations, which lead to abrupt changes in 135.6 nm radiance and therefore make standard deviations unreliable, were applied twice in 2019 (4/ 26/2019 and 10/10/2019), twice in 2020 (3/20/2020 and 8/11/2020), and three times in 2021 (2/8/2021, 6/12/2021, and 8/27/2021) to mitigate the effects of detector burn-in. To quantify the uncertainties associated with each data point, we add the random, system, and model uncertainties reported for the column O/N 2 data product in quadrature. We then omit all data points in which the uncertainties exceed 0.3.
For a particular UTC, there is a 2-D snapshot image of the O/N 2 provided by GOLD for each individual day. In order to spatially suppress the noise, we take the average O/N 2 ratio in a 10°× 10°( longitude × latitude) grid box and attribute this value to the center grid, and we collect all the data from January 2019 to August 2022 to each grid box. In order to have enough data points for calculating the standard deviation with significance, we omit all grid points where less than half of the days had reliable data. This results in the removal of much of the equatorial region, where uncertainties are generally higher because of flat-field correction errors (see the GOLD instrument documentation, available at https://gold.cs.ucf.edu/).

Calculation of short-term day-to-day variability
After binning the data for each grid box each day, we have obtained a time series with a temporal resolution of 1 day. We apply a 30-day high-pass filter to the time series at each longitude-latitude grid point to obtain the perturbations with periods shorter than 30 days, which are defined as O/N 2 residuals hereinafter. After these perturbations are obtained, they are grouped together according to different seasons and their standard deviations (stds) are calculated. The seasonal and year-to-year variations are suppressed by the highpass filter and therefore do not contribute to the std. The std is treated as a proxy for the day-to-day variability shorter than 30 days and longer than 2 days within a given season as a result of the temporal sampling (1 day) and filtering. Winter is defined as December, January, and February, Spring is defined as March, April, and May, Summer is defined as June, July, and August, and Fall/Autumn is defined as September, October, and November.
We consider all pre-screened data that fall into each season (~270-360 points) first, then we divide the data into two activity classes categorized by Kp indices. We consider two separate activity thresholds based on Kp: days in which the maximum Kp ≥ 3.0 or 2.0 are considered to be active days, and days in which the maximum Kp < 3.0 or 2.0 were considered to be quiet days. To account for relatively long (~24-48 h) IT recovery phases, we omit 1 day after each active day from the quiet day data. We use Kp to establish our quiet-time threshold to be consistent with, e.g., Cai et al. (2020). The stds are then calculated for these different groups of data (quiet, Frontiers in Astronomy and Space Sciences frontiersin.org active, and all) and used to represent the variability under each condition. Finally, we calculate a relative std by dividing the absolute std by the quiet-time mean O/N 2 . Geomagnetic conditions with the statistics of quiet and active days are summarized in Figure 1. In total, there are 440 active days (Kp ≥ 3.0) and 641 quiet days (Kp < 3.0). For the threshold of Kp = 2.0, there are 846 active days (Kp ≥ 2.0) and 220 quiet days (Kp < 2.0). Plotted Kp values are the daily maxima of 3-h Kp, and the reported Dst and F10.7 values are daily means. We choose daily maximum Kp as our quiet/active cutoff threshold in order to ensure that all days with any geomagnetic activity are considered active (using daily means could lead to misclassifications of days with short spikes of geomagnetic activity but low average Kp). We plot the solar F10.7 flux to show that our daily maximum Kp threshold is sufficient to account for periods of relatively high F10.7, and we show daily mean Dst because that value is used to track the O/N 2 response to geomagnetic activity by correlating it with O/N 2 residuals (see Section 3). Figure 2 shows the full-disk absolute stds of O/N 2 perturbations as a proxy for the day-to-day variability with periods < 30 day in the IT region. Comparing different rows, variabilities are the strongest for Kp ≥ 3.0 (third row), with a maximum of 0.14 identified in Fall ( Figure 2D3). Variabilities are the weakest during Kp < 2.0 (fourth row), with a maximum std of 0.07 also in Fall ( Figure 2D4), which is about half of that in active time. Figure 3 is the same as Figure 2 except for showing the relative stds of O/N 2 perturbations. Relative stds are also the strongest when strong activities are considered (Kp ≥ 3.0) and decrease gradually towards moderate strong (Kp ≥ 2.0), moderate weak (Kp < 3.0) and the weakest (Kp < 2.0) conditions. We focus on the relative stds for the discussion of variabilities. For the Kp ≥ 3.0 case (third row of Figure 3), the averaged relative std is 4.9% across all seasons and grid points, and this value is 3.7% for Kp < 2.0, which gives a 28% difference between variabilities during most active and quiet times.

Observations and results
Variability is enhanced at high latitudes which may be related to auroral heating. For Kp ≥ 3.0 (third row), the highest variability observed is 20.6% at (70°W, 60.0°N) in Fall ( Figure 2D3) and the weakest variability is 2.4% at (12°E, 6°N) in Summer which differs by a factor of 8-9 between high and low latitude variability for this activity category. During quiet conditions (Kp < 3.0, second row), the high latitude variability enhancement is weakened but still observed. The maximum variability is found to be~11.1% at a high latitude (70°W, 60°N) in Fall ( Figure 3D2), which is about 4-5 times the minimum variability, that is, 2.4% at (8.0°E, 37°S) in Summer ( Figure 3C2). The highest relative variability (20.6%) in active time is also about twice of that in quiet time (11.1%). The maximum and averaged variabilities for the various geomagnetic conditions and four seasons are plotted in Figures 6A, B.  Frontiers in Astronomy and Space Sciences frontiersin.org 04 also decreases from 21% to 11% from Fall to Summer for Kp ≥ 3.0. To quantify the extent of the high variability region which may be affected by aurora, we count all grid points where relative O/N 2 stds exceed 7.7% (red contours in Figure 3) and convert that number to a percentage of the full GOLD field of view (fov). For each set of geomagnetic conditions, the extent of the region of influence is at maximum during Fall and Spring and at minimum during Winter and Summer. As an example, 5.0%, 13%, 2.5%, and 11% of the Frontiers in Astronomy and Space Sciences frontiersin.org 05 GOLD fov had relative O/N 2 stds > 7.7% for Winter, Spring, Summer, and Fall, respectively considering the full dataset. The extents of the high variability region for different geomagnetic conditions and four seasons are illustrated in Figure 6C. It should be noted that 3-4 years of GOLD O/N 2 data are used for this study: four Winters, Springs, and Summers and only three Falls. More data may help to better confirm the seasonal dependence.
From Figures 2, 3, variability increases with Kp index, which implies that magnetospheric forcing plays a role in the day-to-day variability of O/N 2 , especially at high latitudes. In order to better quantify this contribution, we calculate the correlation between the time series of O/N 2 residuals and the Dst index. Figure 4 shows the correlation between these two parameters at (45°W, 59°N, a relatively high latitude region) during Summer and (2°W, 20°S, a relatively low latitude region) during Spring. The active-time data points (red in Figures 4A1, A2) are more tightly clustered than their quiet time (light blue) counterparts, and a linear fit shows a steep, positive (+0.66) correlation at (45°W, 59°N) and a steep negative (−0.65) correlation at (2°W, 20°S). In contrast, correlations are small when only quiet days are considered: +0.21 for (45°W, 59°N) and +0.07 for (2°W, 20°S). Magnitudes of the correlations decrease compared to active times when both quiet and active days are considered: they become +0.46 for (45°W, 59°N) and −0.38 for (2°W, 20°S). Figures 4B1, B2 show the time series of O/N 2 residuals overplotted with Dst indices (both quiet and storm times are considered). The correlation at high latitudes and the anticorrelation at low latitudes are apparent.
Shown in Figure 5 are the calculated correlation coefficients for each grid point in the GOLD fov for each season. Data gaps are present where correlations are determined to be insignificant based on a Student's t-test and its associated p-value. We keep correlations where the p-value is less than 0.05, which means that the significance level of the correlation is greater than 95%. Figures 5A1-D1 show correlation maps where all reliable data for both quiet and storm conditions are considered. Red regions indicate positive correlations, and blue regions indicate negative correlations. Since Dst decreases in response to enhanced geomagnetic activity, a positive correlation means that the O/N 2 ratio decreases as Dst decreases which corresponds to a storm-induced depletion in O/N 2 . A negative correlation means O/N 2 ratio increases as Dst decreases, corresponding to a storm-induced enhancement in O/N 2 .
There are very few significant correlations for quiet conditions (Kp < 3.0 and Kp < 2.0, Figures 5A2-D2, A4-D4) except for the quiet Summer. However, the region of significant correlation during Kp < 3.0 ( Figure 5C2) shrinks significantly when only Kp < 2.0 is considered ( Figure 5C4), which indicates that the threshold of Kp < correlations are robust at both high and low latitudes during active times, variabilities are much larger at high latitudes, implying a more direct influence locally from aurora.   1 The maximum correlation and largest anti-correlation for each season and activity class which includes active time data. The position of each correlation is given in parentheses.

Winter
Spring Summer Fall Frontiers in Astronomy and Space Sciences frontiersin.org of largest anti-correlations for the three geomagnetic conditions (all data, Kp ≥ 3.0 and Kp ≥ 2.0) and four seasons are plotted in Figure 6F.

Discussion and summary
The dependence of variability and correlation on Kp index, and the statistically significant correlation with Dst index indicate that the day-to-day O/N 2 variabilities are strongly coupled to storm activities. Both auroral particle precipitation and Joule heating are enhanced during active times, with the auroral oval itself expanding and contracting with the activity level in a similar manner to the expansion and contraction of the high variability region (Figures 2,  3) and positive correlation region ( Figure 5) observed at high latitudes. The fact that the high latitude enhancement in O/N 2 variability is present (albeit weakened) during quiet times indicates that forcing from above can still dominate at high latitudes, regardless of geomagnetic activity levels.
The correlation maps shown in Figure 5 reflect the storminduced, global, thermospheric circulation patterns which influence the distribution of O and N 2 . Joule heating causes the upwelling of N 2rich air and O-depleted air from the lower atmosphere, which is then distributed over high and mid latitudes by neutral winds (Strickland et al., 2004;Zhang et al., 2004;Zhang et al., 2015). The high latitude, positive correlation regions in Figure 5 are indicative of this upwelling. The corresponding negative correlations at low latitudes can be explained by the storm induced change in global circulation which drives low latitude downwelling (Fuller-Rowell et al., 2002). The low latitude response required to maintain continuity of mass flow is both delayed and weakened relative to the high latitude response (Forbes, 2007), which explains the smaller variabilities at low latitudes seen in Figures 2, 3. Based on this interpretation, the correlation map for the storm-time matches very closely the accepted circulation patterns for O and N 2 . Luan et al. (2017) used data from the Global Ultraviolet Imager (GUVI) on board the TIMED (Thermosphere, Ionosphere, and Mesosphere, Energetics and Dynamics) satellite to obtain global maps of O/N 2 . The group reported larger values of O/N 2 at high latitudes, which increased with increasing F10.7. They also reported significant longitudinal variations in O/N 2 at a particular LT. These results compliment our own work, which instead focuses on the dayto-day variability of O/N 2 . Similar to the high latitude enhancement in mean O/N 2 reported by Luan et al. (2017), we report high latitude enhancements in the day-to-day variability of O/N 2 (Figures 2, 3). Luan et al. (2017) attributed the enhancement in mean O/N 2 to auroral heating, which agrees with our own conclusions about auroral heating as the source of enhanced variability in O/N 2 . Luan et al. (2017) also made use of TIE-GCM to further study the contribution of Joule heating to the observed O/N 2 latitudelongitudinal patterns. This study showed that stronger auroral precipitation causes a decrease in O/N 2 at high latitudes, which is consistent with the positive correlations between O/N 2 and Dst that we report at high latitudes ( Figure 5). Our work demonstrates that geomagnetic activity not only changes the absolute values of O/N 2 , it also contributes to the absolute and relative short-term temporal variability of O/N 2 . This implies that O/N 2 changes more dramatically during storms and substorms.
The seasonal dependence of O/N 2 variability and O/N 2 residuals/ Dst correlations may be partly explained by the seasonal dependence of auroral activity (Newell et al., 2010;Bower et al., 2022). O/N 2 variabilities are strongest during Fall and Winter and weakest during Summer, which is consistent with observations of stronger auroral hemispheric power (HP) during Winter months (Zheng et al., 2013). Although auroral activity is at a minimum during Summer months, geomagnetic active-time upper thermospheric winds are at a maximum during these months (Dhadly and Conde, 2017), which may cause N 2 rich and O poor air that has been brought up by upwelling to spread to lower latitudes. This may have the effect of increasing the auroral region of influence during Summer and enhancing the correlation, which is reflected in Figure 5.
In summary, the day-to-day (2-30 day) variability of GOLD O/N 2 column density has been quantified and its correlation with Dst has been evaluated for various geomagnetic conditions and seasons for the first time. The main findings include: 1. Variability is significantly enhanced at high latitudes and during active (Kp ≥ 3.0) times, during which the largest absolute and relative variabilities of O/N 2 ratio are found to be 0.14% and 21%, respectively, at~60.0°N in Fall. Variabilities are the weakest during Summer. Small but statistically significant variabilities are also observed at low latitudes and during quiet (Kp < 3.0) times. The maximum variability at high latitude during quiet times can decrease to half of the that during active times. The variability at low latitudes is smaller than that at high latitudes by a factor of 5-8 in general. This work has provided baseline/references for numerical models to capture the short-term variability of an important parameter (O/N 2 ) of the IT system. To further diagnose the contributions from each source and find quantitative explanations for the ways that chemistry, dynamics, and electrodynamics contribute to the observed variability, for instance geomagnetic origin versus lower-atmosphere forcing, we will involve model sensitivity runs in a future work.

Data availability statement
Publicly available datasets were analyzed in this study. This data can be found here: https://gold.cs.ucf.edu/data/search/.
Frontiers in Astronomy and Space Sciences frontiersin.org