Numerical forecast of the upper atmosphere and ionosphere using GAIA

Upper atmospheric conditions are crucial for the safe operation of spacecraft orbiting near Earth and for communication and positioning systems using radio signals. To understand and predict the upper atmospheric conditions, which include complex variations affected by both low altitude and upper surrounding environments, we have developed a quasi-real-time and forecast simulations using a physical global model, the Ground-to-topside model of Atmosphere and Ionosphere for Aeronomy (GAIA). The GAIA simulation system provides a global distribution of ionospheric total electron content (TEC) with background atmospheric and electric distributions including a few-days prediction. The prediction accuracy for the detection of significant ionospheric storms decreases with increasing lead time, i.e., the duration of the model simulation which is not constrained by realistic input parameters. Similar characteristic variations associated with sudden stratospheric warmings (SSWs) are reproduced with the full or limited input of meteorological data at least the prior 3 days. This is a first step toward the usage of GAIA for space weather forecasting.


Introduction
The upper atmosphere shows irregular variations in addition to regular daily, seasonal, and regional variations. These variations are results of the combination of (i) disturbances propagated from low atmospheres; (ii) interaction between the neutral and ionized atmosphere, i.e., ionosphere; (iii) energy input from the magnetosphere above, and (iv) radiation from the Sun. The effect of the upper atmospheric disturbance on today's society, which relies on highly advanced information techniques, has been increasing; for example, ionospheric disturbances affect the radio signals used for communication and positioning systems, and thermospheric density variation can affect operation of artificial satellites orbiting near Earth and space debris, the latter may damage satellites. "Space weather" of the upper atmosphere should cover not only the disturbances related to solar storms, but also the disturbances that are sometimes independent of solar storms, such as (i) and part of (ii) above.
Ionosphere and thermosphere models are expected to be essential tools for estimating and forecasting these upper atmospheric disturbances. There are both empirical and physical modeling approaches. The empirical models provide specific parameter/regions with high prediction accuracy. For example, neural networks are widely used for the prediction of the ionospheric variations (e.g., Nakamura et al. 2007;Sai Gowtam and Tulasi Ram 2017;Srivani et al. 2019).
Physical models are also applied to space weather forecasting. The Coupled Thermosphere Ionosphere Plasmasphere electrodynamics (CTIPe) model (Codrescu et al. 2012) provides results from real-time and 1-h forecasts run at their website (https ://www.swpc.noaa. gov/produ cts/ctipe -total -elect ron-conte nt-forec ast). Fernandez-Gomez et al. (2019) reported the difference between the real-time and research simulations using the CTIPe model. By necessity, they used different solar wind inputs to the real-time and research simulations. For the real-time simulation, they used real-time solar Open Access *Correspondence: chihiro.tao@nict.go.jp 1 National Institute of Information and Communications Technology (NICT), Tokyo, Japan Full list of author information is available at the end of the article wind observation data provided by the Advanced Composition Explorer (ACE) spacecraft. For the research simulation, they used OMNI post-processed data. Both neutral density and total electron content (TEC) from the operational simulation are underestimated compared with those from the research simulation and actual observation during a magnetospheric storm event on November 20, 2003. They found that the discrepancy in the ionospheric parameter is more significant relative to the thermospheric parameter. The Whole Atmosphere Community Climate Model with thermosphere and ionosphere extension (WACCM-X) (e.g., Liu et al. 2018) is another global whole atmosphere and ionosphere model. Real-time parameters from WACCM, the atmosphere part of the WACCM-X, are provided in their website (https ://www2.acom.ucar.edu/acres p/forec astsand-near-real-time-nrt-produ cts).
The empirical modeling approaches usually provide better prediction accuracy. On the other hand, the physical modeling approach is useful for its ability to predict various parameters including physical processes and sufficient spatial and temporal coverage/resolution beyond the observational limitation. In addition, approaches using different physical models would provide opportunities to test various model settings and behaviors.
The purpose of this study is to apply the other physical model, namely, Ground-to-topside model of Atmosphere and Ionosphere for Aeronomy (GAIA), for quasi-real-time and forecast uses. There are several target phenomena: (i) ionospheric storms; (ii) ionospheric disturbances during sudden stratospheric warmings (SSWs) (Jin et al. 2012); (iii) sporadic-E layers (Shinagawa et al. 2017), and (iv) plasma bubbles (Shinagawa et al. 2018). The performance of GAIA for phenomena (ii-iv) using GAIA research simulation was examined in those previous studies. The most challenging points in the application to real-time and forecast uses is achieving sufficient prediction accuracy under the limitation of input data. How much accuracy could be achieved in such a situation is an important question. In addition to the magnetospheric potential driven by the solar wind investigated by Fernandez-Gomez et al. (2019), the upper atmosphere is also affected by inputs from the low-altitude region. Here, we mainly focus on the effect of the low-altitude data limitation on the prediction accuracy of the upper atmosphere. In this paper, we target the (i) ionospheric storms in Japan and (ii) appearance of SSW events. Application to (iii) sporadic-E layers is investigated and reported by Shinagawa et al. (this issue). One of the objectives of these GAIA real-time and forecast simulations is to contribute to space weather forecasting carried out by the National Institute of Information and Communications Technology (NICT). This is the first system that provides a numerical forecast of the ionosphere/ thermosphere with a few days of lead time.
In "Model and real-time procedure" section, we introduce GAIA and model settings for the real-time and fewdays forecasting runs. Applications for the ionospheric storm events and SSW events are reported in "Target 1: ionospheric storms" and "Target 2: sudden stratospheric warming" sections, respectively, followed by the "Summary and future works" section.

Model and real-time procedure
Here, firstly we briefly introduce GAIA (e.g., Jin et al. 2011). Then we describe specific settings for the quasireal-time and forecast run.

GAIA
GAIA consists of three parts: a whole atmospheric general circulation model (GCM) part, an ionospheric part, and an electrodynamics part. These parts exchange physical parameters ( Fig. 1 of Jin et al. (2011)). We use a GCM Fig. 1 Schematics representing the timeline of GAIA quasi-realtime and forecast simulations. GAIA run on each day consists of 1-day simulation (blue thick arrow) including the meteorological data JRA and 4-day simulation without JRA input (light blue thick arrow). The WJRA calculation on the current day (middle ribbon) uses the final condition of the WJRA (with JRA) calculation on 1-day prior (top ribbon) as an initial condition. SWx denotes space weather, and UT denotes universal time. The simulation launch time following a data download is shown by filled stars. See details in the text version with a resolution of 2.8° in longitude and latitude and a vertical resolution of 0.2 scale height. The GCM provides neutral atmospheric parameters such as wind velocity and temperature information from the nonlinear primitive equations for continuity, momentum, and thermodynamics including photo-chemical processes (O 2 and O) in the thermosphere and important physical processes from the troposphere, stratosphere, mesosphere, and up to the thermosphere (e.g., Miyoshi and Fujiwara 2003). The ionospheric part provides the main ion compositions (N 2 + , O 2 + , O + , and NO + ), ion dynamics, and ion and electron temperatures, by solving equations of mass, momentum, and energy with 2.5° (longitude), 1° (latitude), and 10-100 km (altitude) resolution (e.g., Shinagawa 2009). By using the GCM neutral wind velocity and ionospheric conductivity information, the electrodynamics part estimates electric potential from a potential solver with 2.8° (longitude) and 0.2-0.6° (latitude) resolution at 75 km height (Jin et al. 2011).
The current version of GAIA includes variations of the solar extreme ultraviolet (EUV) emission and meteorological forcing. The variation of the solar EUV emission is treated by referring to a solar EUV flux model for aeronomic calculation (EUVAC) (Richards et al. 1994) as a function of the solar radio flux at 10.7 cm, i.e., the F10.7 index. Realistic meteorological forcing is included by assimilating the meteorological reanalysis data set, the Japanese 55-year Reanalysis (JRA-55), provided by the Japan Meteorological Agency (Kobayashi et al. 2015;Harada et al. 2016) at low altitude (< ~ 40 km) in the GAIA GCM part.
For external forcing from the magnetosphere, this GAIA version assumes a quiet convection electric field with a cross polar cap potential of 30 kV and weak auroral precipitation, as described in, e.g., Fujiwara and Miyoshi (2010), without variation. Therefore, this simulation does not contain the ionospheric and atmospheric responses to the solar wind and magnetospheric variations. The inclusion of these variations is one of the important future works. Figure 1 shows a timeline of GAIA quasi-real-time and forecast simulations. JRA data of the prior 2 days is uploaded by the Japan Meteorological Agency at universal time (UT) ~ 6:10 every day. We run GAIA for the prior 2 days assimilating JRA data for 24 h using outputs of the similar run of the previous day as the initial condition. This simulation setting is called "WJRA" (with JRA) hereafter. Then, we run GAIA without JRA data assimilation for 4 days to cover the forecast information; this setting is called the "WOJRA-1/2/3/4" setting for the 1st-/2nd-/3rd-/4th-day simulation without JRA, respectively. WOJRA simulations are necessitated by input latency and the desire for long lead-time prediction of a few days. GAIA itself includes meteorological processes in the troposphere. It takes ~ 3 h for this 5-day simulation sequence using the current supercomputer facility at NICT. The space weather forecast briefing at NICT is held daily at UT 5:30 (local time (LT) 14:30) for now (June 2020) to provide the space weather forecast within 24 h. WOJRA-3 and WOJRA-4 settings cover these periods. This GAIA run refers to the F10.7 index at UT 20:00 provided by Natural Resources Canada up to one day prior. The latest F10.7 value at one day prior is also input for WOJRA-2/3/4 simulations, i.e., constant solar EUV input over the WOJRA runs, for simplicity. This simulation sequence is iteratively operated from June 2019.

Procedure for real-time and forecast simulations
Both the meteorological and solar EUV information are limited for the simulation: the former is included until 2 days prior and the latter until one day prior. Usage of meteorological forecast data instead of the reanalysis data, which are not publicly open, is future work ("Summary and future works" section). Since model accuracy becomes worse with longer lead time for the prediction, how the accuracy behaves with lead time is an important information for both operational and scientific viewpoints. We investigate this for the ionospheric storm and SSW events in the following sections. Note that the ionospheric storms provided by this GAIA system do not originate from the magnetosphere. Some large ionospheric storms are associated with the magnetospheric storms, while a certain amount of them occur during magnetospheric low activity time. Figure 2 shows the F10.7 (input) and TEC (output) variations derived from the GAIA simulation run on January 6, 2020, example of an ionospheric storm. The F10.7 value was ~ 72 × 10 −22 Wm −2 Hz −1 during the interval (Fig. 2a). Figure 2b shows the variation of TEC averaged over 35-39° latitude and 125-147° longitude, for 5 days from January 4 to January 8. The TECs obtained under the WJRA setting (red line) and WOJRA settings (orange line) show a diurnal variation with enhancements during the local daytime. We refer to a 27-day median at the same local time as a reference profile. Compared with this reference profile (black line in Fig. 2b), several deviations are apparent, the largest one occurring in early 8 January. The grey hatches correspond to a scale indicating a significance of the deviation of the TEC profile from the reference profile. We measure the significance and detect ionospheric storms using the "I-scale" index (Nishioka et al. 2017) derived for GAIA (see Appendix). The shade of grey represents the deviation normalized by a standard deviation derived using long-term (1996-2016) simulation outputs. The significance levels are shown by the different darkness levels. The TEC variation (red or orange line) within the dark grey region is not significant, while that extending to the light grey or white region is a significant disturbance above the regular variation. From this index, two significant ionospheric storms are detected when the TEC values increased to the light grey region from the end of January 7 to UT 3:00 on January 8 and UT 8:00-18:00 on 8 January, as indicated by I P 2 labels in Fig. 2b. According to results based on observations using the Global Navigation Satellite System (GNSS) (Nishioka et al. 2017), the I P 2 positive storm was seen at UT 17:00-22:00 on January 7, which is slightly earlier than the first event.

Target 1: ionospheric storms
Figure 2c-e shows TEC global maps related to the large enhancement seen at LT midnight on January 7-8 at ~ 37° in Japan. The difference between the instantaneous profile (Fig. 2c) and the reference TEC, i.e., 27-day median at the same local time (Fig. 2d), provides a deviation TEC map (Fig. 2e). For this event, our simulation suggests that localized TEC enhancements around Japan a b c d e To understand the cause of the TEC enhancements, we also output electric field, neutral wind, and atmospheric composition. The east electric field lifts the ionospheric plasma by ExB drift. Southward (northward) neutral wind drags up the plasma along the slanted magnetic field in the middle latitudes in the northern (southern) hemisphere. Fewer recombinations take place at high altitude with fewer molecular ions. The O/N 2 ratio, a measure of ion production versus ion loss processes, also varies with season and magnetospheric activity caused by atmospheric dynamics.
The above parameters during this ionospheric event on January 8, 2020 are shown in Fig. 3. The time variation of the latitudinal profile of TEC at 135° longitude from the same GAIA run shows the extension of the equatorial anomaly with 30 TECU (1 TEC unit = 10 16 m −2 ) at 30° latitude and 20 TECU beyond 47° latitude in early January 8 (Fig. 3a). All parameters are dominated by diurnal variations, and some variations are added. The eastward electric field is dominant in the dayside at high latitudes owing to magnetospheric convection, as well as at mid-low latitudes (Fig. 3b). Further enhancement of the eastward electric field over a wide range of latitude is seen from late January 7 and continues LT daytime and evening time on January 8. The meridional neutral wind is northward at daytime and southward at night time at the northern mid-latitude of our interest (Fig. 3c). The O/ N 2 ratio at 300 km altitude is high in the northern winter hemisphere (Fig. 3d). The ratio rather decreases to < 5 at UT 0:00-8:00 on January 8 compared with O/N 2 of ~ 7 at the same time the previous day. Therefore, the TEC enhancement on January 8 would be mainly related to the enhancement of the east electric field, originated from the neutral wind dynamics and/or polarization field.
Next, we test the predictability of ionospheric storms using the WJRA simulation and WOJRA simulations from July 2019 to March 2020. Figure 4 shows the GAIA simulation results under these five settings with results based on observations using the GNSS (Nishioka et al. 2017). Here, flag = 1 is assigned when a strong positive (I P 2) or negative (I N 2) storm is detected in at least one of the TEC I-scales at the five latitudes in Japan (see Appendix), and flag = 2 is assigned when a severe positive (I P 3) or negative (I N 3) storm is detected in at least one of the TEC I-scales at the five latitudes. Other days are assigned flag = 0. Figure 4 also shows the daily maximum of the Kakioka K index for a reference of the magnetospheric activity (blue line). The observation shows flag enhancements sometimes over the entire year, and some of them are associated with the large K index, e.g., on day of year (DOY) ~ 216 and ~ 243. On the other hand, the current GAIA simulation assumes a quiet magnetospheric condition. The GAIA flag values appear slightly concentrated on days 270-300 and 340-400 from January 1, 2019, which would mainly reflect a seasonal variation. This more sensitive seasonal dependence of the model might be partly affected by the constant quiet magnetospheric setting. These concentrations of flag appearance are more significant under WOJRA-2/3/4 and seem to be improved under WOJRA-1 and WJRA settings.
To evaluate the estimation, we use accuracy, threat score (TS), and Heidke's skill score (HSS), which are calculated as (1) Accuracy = (TP + TN)/N , where TP (true positive) is the number of "hits", i.e., both the observation and model report events (flag ≥ 1). FP (false positive) is the number of "false alarms", i.e., the model forecasts events while the observation does not (flag = 0). FN (false negative) is the number of "misses", i.e., the model forecasts quiet while the observation shows events. TN (true negative) is the number of events that both model and observation forecasts as quiet. N = TP + FP + FN + TN. HSS represents the score excluding the random prediction. Sc is the hit rate by random prediction defined as follows.
These scores for the WJRA and WOJRA settings are shown in Fig. 5. All scores decrease with the lead time. The accuracy decreases from 0.80 to 0.72 (Fig. 5a). The TS decreases with time after the cessation of WJRA input: TS = 0.11 (0.096) for the WJRA setting, TS = 0.10 (0.090) for WOJRA-2, and TS = 0.061 (0.044) for WOJRA-4 of all case of K (limited daily maximum K < 4). Both the decreasing trend and values are similar to those in the analysis using magnetospheric quiet days with daily maximum K of less than 4 (blue dashed line in Fig. 5). Since the occurrence number of ionospheric storm events is small, the TS is considered to be a better evaluation than the accuracy for this case. The results of this analysis confirmed that the input of realistic low-altitude information improves the model prediction accuracy. As mentioned above, HSS represents the score excluding the random prediction, with positive value means better prediction compared to the random prediction. The derived score shows the positive value for the WJRA and WOJRA-1 ~ 3 settings, while it becomes negative for the WOJRA-4 setting (Fig. 6c). This indicates that the prediction ability is better than the random prediction up to the WOJRA-3 setting, while the accuracy is worse for the WOJRA-4 setting. Therefore, the maximum lead time is evaluated to be ~ 1-day.
On the other hand, if the prediction for tomorrow is assumed to be the same as today, the accuracy, TS, and HSS are 0.78, 0.19, and 0.20, respectively. The neural network method provides the prediction of ionospheric storms using the critical ionospheric frequency foF2 from the ionosonde at Tokyo over 1985-1996(Nakamura et al. 2007) with the accuracy of 0.63, TS of 0.22, and HSS of 0.11 based on their Table 15. Although prediction ability of our numerical model is less than their neural net approaches except for the accuracy value, we can obtain spatial and temporal variation with physical backgrounds, as described in Fig. 3. To contribute to space weather forecast efficiently, the improvements of our model and system improvements described in the last paragraph of "Summary and future works" section are required.

Target 2: sudden stratospheric warming
The changes in the low and middle atmospheres during SSWs cause significant disturbances in the ionosphere. Atmospheric tides are considered to be a primary driver for the ionospheric disturbances via the modulation of the E region dynamo in the equatorial and low-latitude regions and via direct propagation with changes in the thermosphere neutral composition in the middle-and high-latitude regions (e.g., Pedatella et al. 2016). The ionospheric storm flags (see details in "Target 1: ionospheric storms" section) often appear during SSW events; therefore, predicting SSW is also important for the ionospheric forecast. Karpechko et al. (2018) reported that two of the nine forecast models could predict the occurrence of a SSW in February 2018 at 12 days prior, and most of the models could predict the SSW at least 4 days prior (Karpechko et al. 2018). Here, we investigate the SSW signature in the GAIA with the WJRA and WOJRA settings, i.e., dependence on the low-altitude forcing. Figure 6 shows an SSW event that occurred in the southern hemisphere in August-September 2019. SSW in the southern hemisphere is rare compared with the northern hemisphere. As is typical of SSW, the local maximum of the zonal mean neutral temperature decreases its altitude (Fig. 6a) with an eastward neutral wind (Fig. 6b). A boundary where the zonal mean wind switches from easterly to westerly reaches ~ 1 hPa (up to ~ 5 hPa) on September 9 (September 17). The temperature (Fig. 6c) and velocity variations (Fig. 6d) under the WOJRA-4 setting show similar distributions to those in the WJRA simulation. This suggests that WOJRA-4 is sufficiently useful for forecasting including SSWs effects.
The time variations of the zonal mean temperature at 10 hPa and 88° latitude show a similar trend among the five settings, while there are modifications in the details. For example, the temperature on August 27 under the WOJRA-4 setting is 200 K, and this decreases to 196 K under the WJRA setting. On the other hand, the temperature on August 31 gradually increases from ~ 202 K under the WOJRA-4 setting to 205 K under the WJRA setting. This temperature difference is enough significant compare to a numerical error of ~ 0.4 K derived as a standard derivation of temperature distribution during this period. These modifications make the slope of the temperature steeper. The temperature increase after September 15 is greater than the threshold level of the SSW Fig. 5 a Accuracy, b threat score, and c Heidke's skill score obtained under the WJRA and WOJRA settings for the whole period from July 2019 to March 2020 (black solid line) and for the case of magnetospheric quiet condition with daily maximum K < 4 (blue dashed line). Error bars represent the sensitivity of the values derived from one hit event being moved from or into the true negative case, i.e., TP and TN in Eqs. (1,2,3,4) are replaced by TP ± 1 and TN ∓ 1, respectively definition, ~ 25 K in 1 week, as shown by the dotted line, under all settings. Figure 7 shows an SSW event that occurred in the northern hemisphere in February 2020. Similarities in the results under the WJRA and WOJRA-4 settings are again found. Note that white vertical lines in Fig. 7a-d show the results from the same simulation run which was carried out on DOY 38. The WJRA setting provides profiles up to DOY 36, as shown by the white lines in Fig. 7a the local maxima of the temperature and eastward velocity. The WOJRA-4 output indicates an increase in the altitude of the temperature peak for the following DOY 37-39 with the appearance of a new enhancement of eastward velocity at ~ 1 hPa on DOY 38-39. The cessation of SSW growth is foreseen for this event. Figure 8 shows the ionospheric storm flag over Japan estimated by observation and simulations for the same period as in the case of Fig. 7f. Observation shows disturbances on DOY 33 and 34 under weak magnetic activity with the daily maximum Kakioka K indexes of 3 and 1, respectively. The WOJRA-4 simulation predicts disturbances on DOY 30-31 and 35-36, which would be associated with the continuous increase in temperature over DOY 29-33 and a larger enhancement of temperature on DOY 34-35 (Fig. 7f ) under the WOJRA-4 setting. As the estimated temperature decreases, the disturbance flag shifts to DOY 31-32. Although no exact prediction, i.e., ionospheric disturbance on DOY 33-34, was achieved by a b c d e f Fig. 7 As in Fig. 6, except for parameters in the northern hemisphere during an SSW event in the northern hemisphere in February 2020 the system for this event, the occurrence of SSW is useful information for predicting ionospheric disturbances. The SSW is suggested to be more predictable for a few weeks prior (Karpechko et al. 2018). GAIA can also predict SSW growth and decay at least 2 days prior and it is expected to serve well in space weather forecasting. We briefly check the maximum lead time for a prediction of the SSW in 2019 using GAIA. Figure 9a shows the time variation of zonal mean neutral temperature at the southern pole at 10 hPa altitude from extended simulations carried out from August 26 (red line) to September 5 (light blue), 2019. Compared to the steep increase on September ~ 6 under the WJRA setting (black line), the temperature in the simulation run on August 26 (red) shows a smaller increase ~ 14 K (Fig. 9b) with less latitudinal shift of the temperature maxima (top left plot of Fig. 9c). A large temperature increase with a latitudinal shift indicating an occurrence of the SSW is firstly seen in the simulation run on August 28 and 29 (Fig. 9b). On the other hand, the amount of temperature increase decreases compared to the SSW occurrence criteria, 25 K within a week, in the simulation run on the following 4 days, from August 30 to September 2. The poleward extension of the temperature maxima also decreases  Fig. 9c). Then the large temperature increase is predicted again in the simulation runs on September 3 and later. Therefore, the lead time to forecast this SSW is 3 days with a precursor in 9 days prior. Note that the lead time varies with SSW events (e.g., Marshall and Scaife 2010), in addition to the various lead time among different models for a specific SSW (e.g., Karpechko et al. 2018). Application to more SSW events is also required to understand the behavior of the prediction system.

Summary and future works
We have applied real-time and forecast simulations using the physical global model GAIA. This version of GAIA includes variable forcing from low altitudes using the meteorological reanalysis data and the solar EUV flux referring to the F10.7 index. In this GAIA quasi-real-time and ~ 2-day forecast simulation, the latest information, i.e., the meteorological reanalysis data up to 2 days prior and the F10.7 index up to one day prior, is included.
The ionospheric storm events are detected by applying an ionospheric storm index called the I-scale to the GAIA simulation results. The GAIA simulation system provides global distributions of TEC with background atmospheric and electric distributions including a nearfuture prediction. The results of the detection of significant ionospheric storms under the full and limited (without meteorological and F10.7) input settings show that the model accuracy decreases with the duration of the simulation without the full inputs.
The characteristic temperature enhancement and neutral wind variation of SSWs are seen in the GAIA realtime and forecast simulations. The full and limited input settings provide similar predictions of SSW growth and decay at least 2 days prior.
Finally, we raise future work on the GAIA model and real-time system to contribute sufficiently to space weather forecast operations. One task is the inclusion of the variation in magnetospheric inputs such as highlatitude electric field and auroral precipitation. We are developing a new version of GAIA that includes these processes. Solar wind observation provides essential information to constrain these magnetospheric inputs. Both solar wind and magnetospheric inputs sometimes vary dramatically within a time scale of minutes. Solar wind observation is usually held at the L1 Lagrange point, which provides solar wind variation about one hour before the solar wind reaches Earth's magnetosphere. Consecutive inputs of the real-time solar wind and simulation run are required to reflect the highly variable information efficiently, while the current system runs the GAIA simulation once per day when the JRA information for the 2 days prior is available. The current version refers to the meteorological reanalysis data up to 2 days before. The utilization of the meteorological real-time and forecast data, which are not publicly open, in addition to the reanalysis data is also expected to improve the model prediction accuracy. Assimilating ionospheric observations in addition to the meteorological information is expected to increase the accuracy of the model. Ensemble methods with providing prediction accuracy are also expected. To extend the prediction lead time, the prediction information of the solar EUV flux (F10.7), solar wind, and magnetospheric inputs are also required. We also need to provide output that is useful for space weather operation and society rather than the model physical parameters. Work on predicting ionospheric phenomena such as the sporadic-E layer by GAIA simulation is ongoing (Shinagawa et al., submitted to this issue). The predictability of the ionospheric disturbances during SSWs should be tested further by applying the past SSW events. We are also testing the input of plasma density information into the radio propagation tool HF-START (High Frequency Simulator Targeting for All-users' Regional Telecommunications). Model evaluation with lead time information for these various data would be useful. We use TEC obtained from GAIA simulation results as the ionospheric parameter. The long-term  dataset is divided into 480 groups on the bases of combinations of (i) four seasons: February-April, May-July, August-October, and November-January; (ii) five latitudinal bands, and (iii) 24 LT bins (1, 2, …, 24 h). We use the five latitudinal bands over Japan used for the observation dataset: 27-31°N, 31-35°N, 35-39°N, 39-43°N, and 43-47°N at 125-147° longitude. We obtain the mean μ(λ, season, LT) and the standard deviation σ(λ, season, LT) values of each P-distribution group. The normalized P-distribution P is derived as The ionospheric storm level is defined by setting thresholds for the normalized numbers to seven levels: I P 3 (for P σ > 5σ), I P 2 (3 σ < P σ ≤ 5σ), I P 1 (1σ < P σ ≤ 3σ), I 0 (− 1σ ≤ P σ ≤ 1σ), I N 1 (− 2σ ≤ P σ < − 1σ), I N 2 (− 3σ ≤ P σ < − 2σ), and I N 3 ( P σ < − 3σ). We detect ionospheric storms I P 2, I P 3, I N 2, and I N 3 when these conditions continue for 2 h or longer. Figure 10 shows an example of P and P profiles derived using GAIA. The P profiles at nighttime (Fig. 10a) and at noon (Fig. 10b) show different shapes: small peak with more extension toward large values during nighttime compared with a sharp profile at noon. The ratios of |P|< 20% for the P at LT 20 h and 12 h are 56.0% and 83.3%, respectively. The difference in the ratio for the same criteria for the normalized P profile decreases, e.g., the ratios of | P|< 1 become 61.8% and 65.9%, respectively. Therefore, the normalized P provides similar significance information independent of LT, season, and latitude.
Variations in the standard deviation σ using GAIA at the 29° N and 41° N latitudinal bands are shown in Fig. 11. The deviation shows enhancements at pre-(LT = UT + 9 h, 19:00-23:00) and post-midnight (LT 3:00-5:00) compared with small values during daytime for all the seasons. The deviation during spring (black line) shows large values at early morning for 29° N and nighttime for 41° N bands, and that during winter (orange line) is largest at other LT. Compared with those derived from the observation (Nishioka et al. 2017), the seasonal and LT variations are similar. Since these variations would be related to the equatorial ionospheric anomaly, this similarity suggests a good ability of the model to reproduce the phenomena statistically. On the other hand, there are differences in detail, such as the larger amplitude of LT variation and less latitudinal difference in absolute values for GAIA than for the observation. This shows a statistical behavior of the model. The I-scale evaluation enables us to detect the significant events among the model variations.