Forecasting Equatorial Ionospheric Convective Instability With ICON Satellite Measurements

Measurements from the Ionospheric Connections Explorer satellite (ICON) form the basis of direct numerical forecast simulations of plasma convective instability in the postsunset equatorial F region ionosphere. ICON data are selected and used to initialize and force the simulations and then to test the results one orbit later when the satellite revisits the same longitude. Data from the IVM plasma density and drifts instrument and the MIGHTI red‐line thermospheric winds instrument are used to force the simulation. Data from IVM are also used to test for irregularities (electrically polarized plasma depletions). Fourteen datasets from late March 2022, were examined. The simulations correctly predicted the occurrence or non‐occurrence of irregularities 12 times while producing one false positive and one false negative. This demonstrates that the important telltales of instability are present in the ICON state variables and that the important mechanisms for irregularity formation are captured by the simulation code. Possible refinements to the forecast strategy are discussed.

Day-to-day physics-based ESF forecasting has not been widely attempted (see Retterer (2005Retterer ( , 2010aRetterer ( , 2010b e.g.,), possibly due to logistical difficulties in assembling suitable forecast datasets. Recently, Hysell et al. (2022) presented a study based on a campaign at the Jicamarca Radio Observatory in which both ionospheric state parameters (causes) and ionospheric irregularities (effects) were measured in the postsunset sector. A numerical simulation initialized and forced by the former was used to reproduce the latter. The effort and other similar efforts preceding it were largely successful (see Hysell et al. (2022) and references therein). The crucial aspect of the simulation was the accurate reproduction of background flows and currents priming irregularity development. Accurate forecasts depended on capturing the background electric field and neutral winds which drive destabilizing currents. Rapid and sustained irregularity formation furthermore depends on the sequencing of the forcing, something not captured by linear, local growth rate analysis (Hysell et al., 2022).
There being few incoherent scatter radars near the magnetic equator, the aforementioned forecast efforts have been limited in scope mainly to the Peruvian sector. In this paper, we attempt to extrapolate the forecast method to other longitudes. The basis for the simulations is data from the Ionospheric Connections Explorer satellite (ICON) and the Ion Velocity Meter (IVM-A), Michelson Interferometer for Global High-resolution Thermospheric Imaging (MIGHTI), and Far Ultraviolet (FUV) instruments in particular (see Harding et al. (2017); Heelis et al. (2017); Mende et al. (2017) for descriptions). The IVM-A, MIGHTI, and FUV data considered here are revisions v06, v05, and v05, respectively, filtered using the internal data quality flags. ICON has a circular orbit with an altitude of 590 km and an inclination of 27°. For this study, ICON data from March 2022, were selected, parameterized, and incorporated in regional simulations of plasma convective instability at different longitudes. Below, we describe the methods, results, and implications of the work.

IVM Data Validation
In the aforementioned forecast studies using incoherent scatter radar data, as in numerous other studies, the postsunset vertical plasma drifts and the prereversal enhancement in particular were salient precursors of instability (e.g., Huang & Hairston, 2015). It is therefore important to verify the accuracy of the vertical drifts measured by the IVM-A instrument on ICON. Toward this end, we compare the averaged drifts with the established climatology. The global quiet-time (Kp ≤ 3) climatology of the vertical plasma drifts in the equatorial F region was modeled by Fejer et al. (2008) on the basis of measurements from the Ionospheric Plasma and Electrodynamics Probe Instrument (IPEI) on the ROCSAT-1 satellite. ROCSAT-1 had a circular orbit with an altitude of 600 km and an inclination of 35°. The Fejer et al. model incorporates IPEA data from July 1999, through June 2004, spanning a range of solar flux conditions above 100 sfu. It is noteworthy that the ROCSAT-1 drifts climatology in the American sector closely resembles average vertical drifts derived from ISR measurements at Jicamarca despite the fact that the averaging kernels involved are very different and that the ISR measurements generally reflect lower altitudes.
For comparison with the climatology and for use in this study, we consider ICON IVM-A data acquired in 2022. This choice obviates possible problems associated with low solar flux conditions that prevailed before this time and an excessive hydrogen ion population at ICON's altitude which could compromise the IVM-A drift measurements. The climatological comparison is further restricted to solar local times between 1700 and 2,200, electron densities greater than 1.0 × 10 5 cm −3 , O + ion concentrations greater than 90% of the total, and magnetic latitudes between ±7.5°. Large vertical drifts (>80 m/s) presumed to be associated with plasma density depletions are also omitted from the analysis.
The IVM-A drifts measurements exhibit a gradually-varying persistent bias which causes them to appear to violate the irrotational electric field condition. To remove the bias, a sinusoidal curve with a 24-hr period and a constant offset is fit to each daily vertical drift measurement as a function of local time. The offset is then subtracted from the daily measurement. The fitting excludes the local time interval between 0500 and 1130 SLT when spurious signals associated with photo emissions from the IVM grid may be present (Heelis et al., 2022). Figure 1 shows a comparison between the ROCSAT-1 vertical drifts climatology modeled by Fejer et al. (2008) (blue symbols) and ICON IVM-A vertical drifts averages (red symbols) sorted by longitude sector and seasons. Results for December solstice and September equinox can be discounted due to the relatively small volume of data available for comparison at the time of writing. The agreement between the ICON-IVM data and the established climatology during March equinox and June solstice is rather close, inspiring confidence in the new drifts data. The prereversal enhancement in the postsunset vertical drifts is clearly reproduced in the March equinox data in particular in all longitude sectors. This would seem to be a promising season for investigating and forecasting ionospheric instability.

Data Selection
The strategy for this forecast study hinges on careful data selection. We sought orbital segments where (a) the ICON satellite was near the magnetic equator in a local-time window suitable for making measurements needed to initialize a forecast simulation and (b) the satellite was again near the magnetic equator and able to detect irregularities when it revisited the same longitude approximately 104 min later. These conditions are met hundreds of times per year, most often but not always on the orbital ascending nodes. For the reasons described in the previous paragraph, we will concentrate here on data from March equinox.
Care must be taken in combining in situ and remotely sensed data from ICON. The winds remotely sensed by MIGHTI are attributed to the tangent point of the line of sight on the limb. For a tangent altitude of 250 km, this point lies a horizontal distance of ∼2,100 km away from the spacecraft. MIGHTI A and B are oriented ∼45 and ∼135 deg from the spacecraft velocity vector, so the vector wind derived from their combination traces a path that lies ∼1,500 km away from the spacecraft path. In nominal operating mode, the displacement is to the north. The geometry is clearly illustrated in Harding et al. (2017), Immel et al. (2018).
We refer to the forecast and validation segments as the (a) and (b) segments, respectively. The condition for selecting the segments are:

Forecast segment (a):
• ICON is at a postsunset local time between 1800 and 1900 when meaningful forecast simulations of postsunset irregularities can be initialized. • ICON is outside the South Atlantic Anomaly (SAA) region between 30 and 90° west longitude where remote sensing is impossible. • The midpoint between ICON's in situ and remote sensing (northward) observations is close to the dip equator. This represents a "happy medium" where all the ICON measurements can be regarded as equatorial.

Validation segment (b):
• ICON is very close to the geomagnetic equator approximately 104 min later when it returns to the longitude of (a) and therefore well positioned to detect plasma density and drift variations with IVM-A indicative of topside irregularities with magnetic apex heights reaching the satellite altitude.
A tool for applying the aforementioned conditions was built and run using ICON two-line elements taken from Space-Track.Org. The tool found 86 suitable segment pairs in the first 6 months of 2022. Segments tend to cluster around a few favorable days each month but are reasonably uniformly distributed across months. Segments occur in all longitude sectors. Not all segments are suitable for analysis, however, because of intermittent data quality issues such as dim airglow emissions for MIGHTI retrievals and low O + densities for IVM measurements, for example, Moreover, some of the datasets identified this way are ambiguous, for example, in cases where irregularities are already present in segment (a), rendering the forecast analysis moot. In other cases, dropouts or interference in one or another data set inhibited measurement interpretation and parameterization. Ultimately, about half the segments identified by the tool were suitable for forecast analysis.  The blue and red lines denote segments (a) and (b), respectively. The longitude for both segments was close to zero degrees in this particular case. During segment (a), IVM-A recorded a strong, distinct prereversal enhancement which is evident in panel II. During segment (b), meanwhile, deep density depletions containing rapidly-ascending plasma were detected by IVM-A as seen in panels I and II. This signifies plasma convective instability leading to irregularities ascending to ICON's altitude by the time of segment (b).

Data Parameterization
In the forecast studies described by Hysell et al. (2022), a three-dimensional simulation volume was populated with plasma state parameters inferred from radar observations and extrapolated from several physics-based and empirical models including SAMI2 for initial plasma number density (Huba et al., 2000), HWM-14 for thermospheric winds (Drob et al., 2015), the IRI-2016 model for ionospheric composition (Bilitza et al., 2016), and a simple model parameterization for the background vertical drifts.
For the plasma density initialization, parameters in SAMI2 were tuned so as to reproduce the electron density profile observed over Jicamarca at the initiation time of the simulation. SAMI2 is a 2D model for a single longitude, and extrapolation to other longitudes was achieved by equating local-time and longitude variations. Likewise, local time and longitudinal variations in the background vertical plasma drifts were also taken to be equivalent.
The initial intent of this study was to parametrize SAMI2 so as to reproduce as closely as possible the electron density profiles measured by the FUV instrument on ICON at the earliest local time the (nighttime) data became available Mende et al. (2017)). However, this first step revealed several instrumental complications associated with deriving the bottomside scale height that will need to be resolved. For the present study, we exclude the FUV data and refrain from tuning SAMI2 which has demonstrated the ability to reproduce climatological background electron density profiles reasonably accurately even without tuning (e.g., Hysell et al., 2015). The FUV data will be revisited in future forecast studies whereas this one will focus on nudging the model using ICON background electric fields and thermospheric winds.
A parameterization of the time history of the vertical plasma drifts measured by IVM-A is crucial for the forecast simulations. Note that we neglect vertical variations in the vertical plasma drifts in this parametrization. The orange curve in panel II represents a fit to the vertical plasma drifts versus local time during segment (a). The parameterization model here is a simple sine function with a specific amplitude, period, phase, and offset, viz.
where t is the local time. We equate local time with longitude when using this model to populate the simulation volume and neglect variations with height.
The behavior of the zonal winds measured by the MIGHTI red-line instrument around the time of segment (a) must also be captured (parameterized) to drive the forecast simulations. Here, we make use of the HWM-14 model as a parameterization tool. HWM-14 provides time-varying profiles, and our parameterization identifies the scaling factor that maximizes the congruity of the zonal HWM-14 winds with the MIGHTI zonal wind measurements. At local times in the vicinity of segment (a), a the best fit multiplicative scaling factor f(t) is determined, taking measurement errors into account. A linear trend with local time/longitude is then fit to the factor, that is, where t is once again local time. This combination of parameters affords some control not only over the strength of the winds but also the timing and rapidity of the evening wind reversal. The trend is passed to the numerical simulation which is forced by accordingly scaled HWM-14 winds.

Numerical Simulations
The numerical simulation used for this study was described in detail by Hysell et al. (2022). The simulation has two components. One solves the quasineutrality condition for the electrostatic potential fully in three dimensions. This is a linear elliptic partial differential equation which is solved using a preconditioned stabilized conjugate gradient algorithm in a sparse framework (Saad, 1990). Ion inertia is neglected. The transport coefficients and background forcing required are derived from a collection of physics-based and empirical models and on other relevant ionospheric datasets as are available. These include NRLMSIS 2.0 for neutral atmospheric parameter specification, IRI-2016 for initial ion composition, SAMI2 for initial background electron densities, and HWM14 for neutral winds (Bilitza et al., 2016;Drob et al., 2015;Emmert et al., 2021;Huba et al., 2000).
The other component advances the continuity equation forward in time for four ion species (NO + , O + 2 , O + , and H + ) using a second-order Runge Kutta scheme. The initial densities are seeded with white noise (20% RMS). The time advance includes limited chemistry (charge-exchange and dissociative recombination reactions). The simulation is cast in magnetic dipole coordinates with spatial gridding with voxels a few km on a side. The time advance is 7.5 s. Convection is handled using a flux assignment scheme based on the total variation diminishing (TVD) condition. It incorporates monotone upwind schemes for conversation laws (MUSCLs) with flux limiting and a second-order TVD scheme to minimize diffusion and dispersion. Dimensional splitting is used to extend the technique to three dimensions. See (Trac & Pen, 2003) and references therein for additional discussion.

Results
Fourteen numerical simulations were performed based on 14 ICON orbital segment pairs selected in the last week of March 2022. The corresponding ICON data parameterizations are listed in Table 1. These include the vertical plasma drift parameterizations, the thermospheric wind parameterizations, and electron density profile parameterizations derived from the FUV instrument (but not used further). The parametrizations describe Chapman profiles, that is, They are included here for reference but not incorporated in the simulations.
Also shown in the last two columns of the table are simple yes-or-no assessments regarding the detection and prediction, respectively, of ionospheric irregularities at 590 km close to the time of segment (b). Table 1. The forecast simulation panels show plasma density and composition, with red, green and blue hues representing molecular ion density, atomic oxygen density, and hydrogen ion density, according to the legends shown. The IVM vertical drifts are presented in the same format as in Figure 3. The occurrence of irregularities is signified by compact regions with large (several hundred m/s) vertical drifts in the vicinity of the vertical red lines. These drifts are indicative of large polarization electric fields in narrow depletion plumes or channels characteristic of convective instability and ESF conditions. Note that the IVM density data for these segments (not shown) exhibit deep density depletions wherever there are large vertical velocity enhancements. Figures 4 and 5   Note. Each row represents an orbital segment (a) used for a forecast simulation run. The first set of figures are the UT date and time, the second set the vertical drift parameterizations, the third set the wind parameterizations, the fourth set the electron density profile parameterizations, and the final set the observed and predicted ESF occurrences. The IVM datasets considered here share some remarkable characteristics. For one, the prereversal enhancement of the zonal electric field is clearly evident in each vertical drift measurement around twilight. This is true both in ESF events and non-ESF events. This echoes the fact that the prereversal enhancement is prevalent in all longitude sectors during March equinox during moderate solar flux conditions and inspires confidence in the IVM vertical drift measurements.

Figures 4-7 show IVM measurements and forecast simulation snapshots for each of the 14 orbital pairs in
For another, irregularities associated with ESF are also prevalent in the data. Evidently, depletions plumes penetrated into the topside to an altitude of 590 km or more in half the datasets considered. This fraction exceeded expectations; at Jicamarca, plumes penetrate to the topside about one quarter of the time during equinox moderate

Space Weather
HYSELL ET AL.
10.1029/2023SW003427 11 of 15 solar flux conditions. The ICON data were not selected on the basis of irregularity occurrence, and the high irregularity occurrence seems somewhat anomalous. This is possibly a byproduct of the relatively small database in question.
The thermospheric wind parameterizations α and β point to a systematic discrepancy between the MIGHTI wind measurements and HWM-14 around twilight. Namely, values of β consistently less than unity indicate that the measured winds are smaller than the model winds around 18 LT. Uniformly positive values of α, meanwhile, suggest that the measured winds largely "recover" to HWM-like values within a few hours of 18 LT. The overall pattern suggests that the postsunset reversal of the zonal winds occurs somewhat later in the MIGHTI data than the HWM model would predict. This is important, as the time of the evening reversal has been seen to be a key influence on overall postsunset stability in our simulations (Hysell et al., 2022).
Alongside each IVM data set in the figures is a snapshot of the forecast simulation output. The snapshots show plasma number density and composition using the scales shown. Here, red, green, and blue hues represent molecular ions, atomic oxygen ions, and hydrogen ions, respectively. The times shown above each snapshot are the elapsed simulations times since initialization which was undertaken at the local time and longitude of segment (a), under the associated conditions, and with the associated background vertical drifts and thermospheric winds. The times chosen for the snapshots are the lesser of 1) the time when plumes first exceeded 590 km altitude, the upper boundary of the simulation or 2) 104 min, whichever comes first. Recall that the ICON satellite was orbiting at 590 km, establishing the ceiling for the simulation.
In six of seven cases shown in Figures 4 and 5, plumes were predicted to have reached ICON's altitude in time for segment (b). The marginal case is the one at the bottom of Figure 4 when plumes were nearly but not quite to 590 km after 104 min. This was a false negative.
In six of the seven cases in Figures 6 and 7 when no plumes were detected by ICON in segment (b), no plumes were predicted to have reached ICON's altitude. In the case at the top of Figure 6, irregularities ascended nearly to ICON's altitude by 104 min, making this a marginal but correct prediction, strictly speaking. In the case immediately below that one, the numerical simulation produced a distinct false positive, predicting irregularities at 590 km within just 72 min. This false alarm is consistent with the robust prereversal enhancements in the zonal electric field evident in the (a) segment of the corresponding IVM data set.

Discussion and Summary
This preliminary ESF forecast effort was posed with a blunt objective-predicting whether ionospheric irregularities caused by plasma convective instability would reach an altitude of 590 km within 104 min based on ICON state-parameter observations from 1800 to 1900 LT onward. It did not consider the morphology of the irregularities, their longevity, or their disruptiveness. It did not look far into the future. All of the analyses were retroactive. It did not parse the underlying mechanisms underlying instability. Crucially, the study did not trace the sensitivity of irregularity development to specific aspects of the forcing, especially the neutral winds. (Such a sensitivity analysis was undertaken in the Hysell et al. (2022) forecast study preceding this one.).
Nonetheless, the forecast accuracy was encouraging; in 14 attempts, there were 12 successful forecasts, one false positive, and one false negative. The false alarm seems to have been triggered by a robust prereversal enhancement which failed to provoke strong instability for some reason. The reason is far from clear based on the analyses performed here, prompting the need for a more incisive analysis of the ICON data. However, the study demonstrates that the most important telltales of instability are present in the ICON state variable measurements and that the most important mechanisms for irregularity formation are captured by the simulation code.
The study produced marginal cases where simulated irregularities ascended just barely to ICON's altitude.
Marginal cases occur in studies using ground-based measurements as well, where simulated irregularities emerge close to but not within the ground-based instrument field of view. A strategy combining ground-and space-based irregularity measurements could yield more meaningful forecast assessments.
Several options for building on this foundation and improving the forecasts present themselves. For example, notice that the Y/N predictions found through simulations could have been estimated by asking whether the c parameter in the IVM drifts fits was positive or not. This points to the essential importance of bias removal in the IVM data set for accurate forecasting. The bias removal method employed in this study is rudimentary and requires further development. Obtaining unbiased vertical drifts measurements from spacecraft is challenging but a priority in space weather research going forward.
Another area for improvement is the parameterization of the thermospheric winds. Earlier simulation studies pointed to the timing of the evening reversal of the winds being a key element in determining postsunset ionospheric stability (Hysell et al., 2022). In the present study, there is no obvious relationship between irregularity occurrence, actual or simulated, and the wind parametrization. This may point to shortcomings in the parametrization method. The linear regression HWM weighting method applied here was expedient but also rudimentary. It would be straightforward to interpolate the actual MIGHTI zonal and meridional winds in terms of an appropriate expansion basis, rooted for example, in the vertical extensions of the Hough functions. Doing so is an obvious first step, but this would not capture information about meridional gradients in the thermospheric winds. Huba and Krall (2013) showed that meridional wind gradients could influence ionospheric instability. A method for inferring meridional wind gradients using successive satellite orbits should therefore be explored.
Finally, the fidelity of the numerical simulations should improve if the initial conditions could be rooted in observed local postsunset electron density profiles. Rather than using inverted FUV profiles, it could be expedient to solve a forward problem instead, estimating brightness profiles from candidate SAMI2 electron density distributions and tuning the model for agreement with the FUV low-level data. This strategy replaces the direct numerical inversion problem with a model-based inversion which is constrained by well-understood ionospheric physics with just a few control parameters and suggests a new modality for space-based remote sensing.

Data Availability Statement
ICON data are available through https://icon.ssl.berkeley.edu/Data.