Improving Forecast of Severe Oceanic Mesoscale Convective Systems Using FY-4A Lightning Data Assimilation with WRF-FDDA

: The Fengyun-4A (FY-4A) geostationary satellite carries the Lightning Mapping Imager that measures total lightning rate of convective systems from space at high spatial and temporal resolutions. In this study, the performance of FY-4A lightning data assimilation (LDA) on the forecast of non-typhoon oceanic mesoscale convective systems (MCSs) is investigated by using an LDA method implemented in the Weather Research and Forecasting-Four Dimensional Data Assimilation (WRF-FDDA). With the LDA scheme, three-dimensional graupel mixing ratio ﬁelds retrieved from the FY-4A lightning data and the corresponding latent heating rates are assimilated into the Weather Research and Forecasting model via nudging terms. Two oceanic MCS cases over the South China Sea were selected to perform the study. The subjective evaluation results demonstrate that most of the oceanic convective cells missed by the control experiments are recovered in the analysis period by assimilating FY-4A lightning data, due to the promoted updrafts by latent-heat nudging, the more accurate and faster simulations of the cold pools, and the associated gust-fronts at the observed lightning locations. The cold pools and gust-fronts generated during the analysis period helped to maintain the development of the MCSs, and reduced the morphology and displacement errors of the simulations in the short-term forecast periods. The quantitative evaluation indicates that the most effective periods of the LDA for simulation enhancement were at the analysis time and the nowcasting (0–2 h forecast) periods.


Introduction
Severe convective weather accounts for substantial property losses and casualties worldwide [1][2][3]. In modern meteorological services, numerical weather prediction (NWP) is one of the effective means for severe convection forecasting and warning. Previous studies have demonstrated the ability of convective-allowing NWP models to improve the forecasts of severe convective weather [4,5].
With the development of electromagnetic signal detection technology, many lightning detection networks have been made available around the world. The ground-based lightning detection platforms, such as the Lightning Mapping Arrays (LMA) [28][29][30], the Earth Network Total Lightning Network (ENTLN) [31,32], the World Wide Lightning Location Network (WWLLN) [33,34], the U.S. National Lightning Detection Network (NLDN) [35], and China's Advanced Time of Arrival and Direction System (ADTD) [36,37] are increasingly being used in NWP.
Fierro and Reisner [38] and Fierro et al. [31] developed a method to assimilate total lightning data at convection-allowing scales, wherein incremental increases in water vapor mass are applied in the mixed-phase region to locally enhance thermal buoyancy, which, ultimately, promoted the development of updrafts at observed lightning locations. The results show that the lightning data assimilation (LDA) improved the simulation of the convection and the associated precipitation. Additionally, the performance of the LDA scheme was compared with a radar data assimilation method on the short-term forecasts of a derecho event. The results indicated that the lightning assimilation was better able to capture the placement and intensity of the derecho, particularly with a relatively lower computational cost [39]. Marchand and Fuelberg [40] used a temperature nudging method to promote the development of convection where lightning flashes are observed. Qie et al. [41] developed an LDA method, which adjusts ice-phase hydrometeor mass where lightning was observed. Wang et al. [32] developed a method to retrieve graupel mixing ratio (q g ) and the corresponding latent heat releases using total lightning rates, and then nudge the simulated graupel contents and potential temperature toward the retrieved values. Their study shows that precipitation and cold pool were generated faster using the LDA method developed in their study compared to the LDA methods developed by Fierro et al. [31] and Marchand and Fuelberg [40], as the graupel, which was directly assimilated into the model, could affect the cold-cloud precipitation microphysics (i.e., graupel melts into raindrop).
In addition to the LDA method using nudging methods, the variational and ensemble techniques have also been applied to the LDA method in recent years [30,36,[42][43][44][45]. Such methods are able to promote the development of convection at lightning observed locations and suppress spurious convection developed in the model.
Since both radar and ground-based lightning detection systems are limited by detection range, satellite radiance data and satellite-based retrieval data are important data sources for the data assimilation of oceanic MCSs forecast [46][47][48][49][50]. In recent years, the lightning detection instruments mounted onboard the geostationary satellites have become available, e.g., the Geostationary Lightning Mapper (GLM) mounted onboard Geostationary Operational Environmental Satellite (GOES) -16 and -17 [51][52][53] and the Lightning Mapping Imager (LMI) mounted onboard Fengyun-4A (FY-4A) [54,55]. Some studies have demonstrated the benefits of assimilating geostationary satellite-based lightning data for the simulation and forecast of MCSs [44,[56][57][58]. However, the geostationary satellite-based lightning data have not been applied to the data assimilation for non-typhoon oceanic convection, despite its potential benefits.
The South China Sea is one of the most frequently affected regions in the world to be devastated by oceanic MCSs [59][60][61]. The accuracy of oceanic convection forecasts is essential for many sea-based activities, e.g., mariculture, navigation, and coastal tourism. An improved forecast can provide more time for preparation and better mitigate property losses and casualties. In this study, two non-typhoon oceanic convection cases which occurred over the South China Sea were selected for geostationary satellite-based LDA experiments to evaluate the potential benefits of assimilating FY-4A geostationary satellite lightning data on the simulation and forecast of severe non-typhoon oceanic convection.

Data for Assimilation and Evaluation
The data used for assimilation in this study are the total lightning data from LMI [62][63][64], which is mounted onboard the geostationary meteorological satellite FY-4A. LMI is composed of two 400 × 600 charge-coupled devices (CCDs) that locate lightning by receiving 2-ms energy integrations of the radiation generated by the sudden heating of the lightning channel. LMI uses ultra-narrow bandpass filters with a central wavelength of 777.4 nm and a bandwidth of 1 nm to ensure that the radiation energy is not affected by scattering in clouds [18]. The resolution of LMI is 7.8 km at nadir and 24.2 km at the edge of the full disk.
After being received, the lightning signals are preprocessed in a real-time event processor. Lightning signals are filtered by removing the background noise which is estimated from the radiance of adjacent frames. Similar to the Geostationary Lightning Mapper [65], LMI has a dynamic filtering threshold due to the large difference in the average background noise between the day and night periods. After filtering, possible lightning events are extracted, and Level 2 product datasets are generated after calibration and false signal elimination. In this study, the Level 2 flash data, which were generated from the Level 2 event data using the lightning clustering algorithm [65,66], were employed for data assimilation. The FY-4A LMI lightning data were interpolated onto the convectionallowing model grid (3 km × 3 km in this study) from the original latitude and longitude coordinates and accumulated over 15-min intervals (similar to the procedures in Wang et al. [32]).
The S-band radar composite reflectivity data from Haikou Meteorological Bureau (20.0 • N, 110.25 • E), with a detection range of 460 km, was utilized for model performance evaluation. Before the evaluation, a fuzzy logic approach was employed to identify and remove the ground clutters of composite reflectivity data by applying three features (the standard deviation of reflectivity, the vertical gradient of reflectivity, and the absolute value of radial velocity) to the membership functions. The detail of the approach is described in Cho et al. [67]. Afterwards, the radar composite reflectivity were interpolated into the respective domain grids of the model. In this study, only the region within the coverage (dotted circle in Figure 1) of the weather radar was analyzed.
Remote Sens. 2022, 14, x FOR PEER REVIEW 3 of 23 curred over the South China Sea were selected for geostationary satellite-based LDA experiments to evaluate the potential benefits of assimilating FY-4A geostationary satellite lightning data on the simulation and forecast of severe non-typhoon oceanic convection.

Data for Assimilation and Evaluation
The data used for assimilation in this study are the total lightning data from LMI [62][63][64], which is mounted onboard the geostationary meteorological satellite FY-4A. LMI is composed of two 400 × 600 charge-coupled devices (CCDs) that locate lightning by receiving 2-ms energy integrations of the radiation generated by the sudden heating of the lightning channel. LMI uses ultra-narrow bandpass filters with a central wavelength of 777.4 nm and a bandwidth of 1 nm to ensure that the radiation energy is not affected by scattering in clouds [18]. The resolution of LMI is 7.8 km at nadir and 24.2 km at the edge of the full disk.
After being received, the lightning signals are preprocessed in a real-time event processor. Lightning signals are filtered by removing the background noise which is estimated from the radiance of adjacent frames. Similar to the Geostationary Lightning Mapper [65], LMI has a dynamic filtering threshold due to the large difference in the average background noise between the day and night periods. After filtering, possible lightning events are extracted, and Level 2 product datasets are generated after calibration and false signal elimination. In this study, the Level 2 flash data, which were generated from the Level 2 event data using the lightning clustering algorithm [65,66], were employed for data assimilation. The FY-4A LMI lightning data were interpolated onto the convectionallowing model grid (3 km × 3 km in this study) from the original latitude and longitude coordinates and accumulated over 15-min intervals (similar to the procedures in Wang et al. [32]).
The S-band radar composite reflectivity data from Haikou Meteorological Bureau (20.0°N, 110.25°E), with a detection range of 460 km, was utilized for model performance evaluation. Before the evaluation, a fuzzy logic approach was employed to identify and remove the ground clutters of composite reflectivity data by applying three features (the standard deviation of reflectivity, the vertical gradient of reflectivity, and the absolute value of radial velocity) to the membership functions. The detail of the approach is described in Cho et al. [67]. Afterwards, the radar composite reflectivity were interpolated into the respective domain grids of the model. In this study, only the region within the coverage (dotted circle in Figure 1) of the weather radar was analyzed.

Data Assimilation Scheme
A nudging-based LDA scheme implemented in the Weather Research and Forecasting Four-Dimensional Data Assimilation (WRF-FDDA) system was utilized to assimilate the LMI lightning data [32]. In this scheme, the lightning rates derived graupel mixing ratio and latent heating rates are nudged toward the observations by adding tendency terms to the model prognostic equations [32,[68][69][70], as expressed by Equation (1): where f g is the model state variables, P(f g , t) is the original dynamical and physical terms of the WRF model prognostic equations, G f is the relaxation time scale that is used to prevent the occurrence of unrealistic noise wave due to the instantaneous change of model variable, W f is a spatial weight function, T f is a time weight representing a function of the time lags of the observation and the model state, and O g is the observation or analysis value on the model grids.
The empirical formula between total lightning flash rate and graupel mass derived by Carey et al. [71] was utilized to convert lightning flash rate into column-integrated graupel mass, as expressed by Equation (2): where FR is the total flash rate, and M g is the column-integrated graupel mass between the −10 • C and −40 • C isothermal layers.
Since the column-integrated graupel mass is determined by the graupel mixing ratio (q g ) of each grid within the model-column, the climatological profiles of q g for different bins of the column-integrated graupel masses were used to retrieve the three-dimensional q g fields. In this study, the climatological profiles of q g for different column-integrated graupel mass bins were obtained by the retrospective simulations of 10 summer severe convective events over the South China Sea. The physics and parameterization schemes used in the retrospective simulations were the same as those used in the experiments of this study (introduced in Section 2.3). The q g profiles of the retrospective simulation outputs were extracted hourly to construct the database of q g profiles. In the database, the column-integrated graupel mass of each q g profile was calculated, and the q g profiles were divided into different groups according to the column-integrated graupel mass bins. The climatological profile of q g for each bin of column-integrated graupel mass was calculated by averaging all the q g profiles within the corresponding group.
To consider the graupel in the adjacent regions of observed lightning locations [72], a distance-weighting function-based spread method [32] is employed to account for the q g in such regions, as expressed by Equation (3): where q g (x) is the q g value with zero lightning observation on the grid x, q g (x i ) is the q g value with nonzero lightning observation around the grid x within the influence radius R, N is the number of grids with nonzero lightning observations around grid-box x within the distance of R, and W i is the spatial-weight function proposed by Cressman et al. [73], which is expressed as: where r i is the horizontal distance between grid x and x i . Since the resolutions of FY-4A LMI lightning data are relatively lower than the typical ground-based lightning detection systems (e.g., the Earth Networks Total Lightning Network), the influence radius R is Remote Sens. 2022, 14, 1965 5 of 20 increased to 15 km (6 km was used in Wang et al. [32]). Figure 2 shows the examples of retrieved q g fields (Figure 2a,c) and the corresponding FY-4A lightning flash rates of two oceanic MCSs (Figure 2b,d). where ri is the horizontal distance between grid x and xi. Since the resolutions of FY-4A LMI lightning data are relatively lower than the typical ground-based lightning detection systems (e.g., the Earth Networks Total Lightning Network), the influence radius R is increased to 15 km (6 km was used in Wang et al. [32]). Figure   The latent heating rates associated with the formation of the graupel are computed by using an empirical relationship between qg increments and latent heating rates as in the following equation: The latent heating rates associated with the formation of the graupel are computed by using an empirical relationship between q g increments and latent heating rates as in the following equation: where ∆T is the temperature increments added to the tendency term, ∆q g is the increments of q g , L c is the latent heat of condensation at the preexisting temperature T 0 , L f is the latent heat of freezing at 0 • C, and C p is the specific heat of air at constant pressure.

Case Descriptions and Model Configuration
Two severe oceanic convective weather cases that occurred over the South China Sea were selected for FY-4A LDA experiments in this study. other two were merged into a squall line and moved southeast by 22:00 UTC 19 August. Low-level potential instability was observed according to the sounding data, and the convective available potential energy (CAPE) was higher than 3700 J kg −1 . Subsequently, the squall line dissipated at around 08:00 UTC 20 August.
Case 2. This case was composed by two mesoscale MCSs and one multi-cell convective system. They initiated at 17:00 UTC on 21 August 2019, and developed over several hours. The near-surface ambient temperature was approximately 29°C and the sea-level pressure was approximately 1015 hPa at the night of 21 August. By 04:00 UTC on 22 August, the three convective systems were gradually merged into a larger MCS and the CAPE reached 2000 J kg −1 . Afterward, the system moved northward and dissipated over a period of in 5 h.
Although both examples are typical of oceanic MCSs, the MCSs in Case 1 are relatively linear, while the MCSs in Case 2 are relatively scattered.
The numerical model used for this study is the three-dimensional compressible nonhydrostatic Weather Research and Forecasting model version 4.1.2 (WRF4.1.2) [74]. A double-nested scheme was adopted in this experiment. The horizontal grid spacings of the inner and outer grids were 3 and 9 km, respectively ( Figure 3). It is noted that the outer domains of the two oceanic convective cases are the same while the inner domains are different. This is because the convective cells are preferred to be covered in the inner domain in the whole simulation period, and the positions and development directions of convective cells are different in the two cases. Therefore, the coverage of the inner domains in the two cases is slightly different. The number of vertical layers of the inner and outer grids was 34, and the top layer was at 50 hPa. The National Centers for Environmental Protection Final Operational Global Analysis (NCEP-FNL) data with a resolution of 0.25 • × 0.25 • were used for the initial and lateral boundary conditions, and the Mercator projection was employed for the model area.  The physical parameterization schemes employed in this study include the Thompson scheme [75] for the microphysics, the RRTMG (Rapid Radiative Transfer Model for General Circulation Models) [76] for the shortwave and longwave radiation transfer, the Mellor-Yamada-Janjic scheme [77] for the planetary boundary layer, and the Noah scheme [78] for the land surface. The Grell-Freitas cumulus parameterization scheme [79] was utilized in the outer domain, while no cumulus parameterization scheme was activated in the inner domain with cloud-permitting resolution.
Lightning data were assimilated only into the inner domain. Two-way nesting between parent and inner nests was used, so the impact of data assimilation can feedback from the inner-most cloud-resolving domain to its parent domains. The flow diagrams for the experiments are shown in Figure 4. In the experiments with the LDA (ASML), the lightning data were assimilated during the first 3 h with a time interval of 15 min, thereafter, a 6-h forecast was made. For each case, a control experiment (CTRL, without the LDA) was conducted for comparison. The physical parameterization schemes employed in this study include the Thompson scheme [75] for the microphysics, the RRTMG (Rapid Radiative Transfer Model for General Circulation Models) [76] for the shortwave and longwave radiation transfer, the Mellor-Yamada-Janjic scheme [77] for the planetary boundary layer, and the Noah scheme [78] for the land surface. The Grell-Freitas cumulus parameterization scheme [79] was utilized in the outer domain, while no cumulus parameterization scheme was activated in the inner domain with cloud-permitting resolution.
Lightning data were assimilated only into the inner domain. Two-way nesting between parent and inner nests was used, so the impact of data assimilation can feedback from the inner-most cloud-resolving domain to its parent domains. The flow diagrams for the experiments are shown in Figure 4. In the experiments with the LDA (ASML), the lightning data were assimilated during the first 3 h with a time interval of 15 min, thereafter, a 6-h forecast was made. For each case, a control experiment (CTRL, without the LDA) was conducted for comparison.

Quantitative Evaluation Scheme
To evaluate the performances of the two sets of experiments, a neighborhood spatial verification scheme called the fractions skill score (FSS) is used in this study [80]. The FSS provides a reliable assessment of displacement errors and avoids the "double-penalty" problem caused by the small displacement both in space and time. Therefore, this scheme is more suitable for the evaluation of the simulation with fine scale grids, and it has been widely used to evaluate the performance of fine-scale forecasts of severe MCSs [32,36,43]. The FSS is calculated from fraction brier score divided by the sum of the mean squared observed and forecast fractions, as expressed by Equation (6): where N is the number of grid points in the verification domain, CRs and CRo is the fractional coverage of the studied elementary area by composite reflectivity that exceeds a given threshold value in forecast and observation, respectively. A perfect forecast has an FSS of 1, and a "no skill" forecast has a FSS of 0. For more details on the algorithm of the FSS, the reader is directed to Roberts and Lean [80]. FSSs for the different neighborhood radii (18 km, 24 km, and 36 km) at the different thresholds (25 dBZ, 35 dBZ, and 40 dBZ) were calculated in the study.

Performance in the Analysis Period
To evaluate the improvement of the LDA on the analysis and forecast fields of the oceanic convective weather case, the radar composite reflectivity simulated by ASML and CTRL were compared with the observations, respectively, as shown in Figure 5. At 21:00 UTC, the convective cells of Convection 1, 2, and 3, which formed over the sea were better captured by ASML than CTRL (Figure 5a,j,m). ASML was capable of simulating the morphology of the oceanic MCSs, especially for areas with radar composite reflectivity greater

Quantitative Evaluation Scheme
To evaluate the performances of the two sets of experiments, a neighborhood spatial verification scheme called the fractions skill score (FSS) is used in this study [80]. The FSS provides a reliable assessment of displacement errors and avoids the "double-penalty" problem caused by the small displacement both in space and time. Therefore, this scheme is more suitable for the evaluation of the simulation with fine scale grids, and it has been widely used to evaluate the performance of fine-scale forecasts of severe MCSs [32,36,43]. The FSS is calculated from fraction brier score divided by the sum of the mean squared observed and forecast fractions, as expressed by Equation (6): where N is the number of grid points in the verification domain, CR s and CR o is the fractional coverage of the studied elementary area by composite reflectivity that exceeds a given threshold value in forecast and observation, respectively. A perfect forecast has an FSS of 1, and a "no skill" forecast has a FSS of 0. For more details on the algorithm of the FSS, the reader is directed to Roberts and Lean [80]. FSSs for the different neighborhood radii (18 km, 24 km, and 36 km) at the different thresholds (25 dBZ, 35 dBZ, and 40 dBZ) were calculated in the study.

Performance in the Analysis Period
To evaluate the improvement of the LDA on the analysis and forecast fields of the oceanic convective weather case, the radar composite reflectivity simulated by ASML and CTRL were compared with the observations, respectively, as shown in Figure 5. At 21:00 UTC, the convective cells of Convection 1, 2, and 3, which formed over the sea were Remote Sens. 2022, 14,1965 8 of 20 better captured by ASML than CTRL (Figure 5a,j,m). ASML was capable of simulating the morphology of the oceanic MCSs, especially for areas with radar composite reflectivity greater than 35 dBZ (represented by black contour lines in Figure 5). The locations of MCSs simulated in ASML were consistent with the observations. However, these oceanic MCSs were not well simulated by CTRL (Figure 5j).
Remote Sens. 2022, 14, x FOR PEER REVIEW 9 of 23 than 35 dBZ (represented by black contour lines in Figure 5). The locations of MCSs simulated in ASML were consistent with the observations. However, these oceanic MCSs were not well simulated by CTRL (Figure 5j). The maximum updrafts, relative humidities, surface potential temperatures, and surface horizontal winds of ASML and CTRL were analyzed in Figure 6 to demonstrate the impact of the LDA on thermodynamic and dynamical environments. After the ingestion of latent heating rates, strong updrafts were promoted in the middle-upper levels. The averaged updrafts of Convection 1, 2, and 3 (black boxes in Figure 6b The positive qg increments in ASML could affect the cold-cloud precipitation microphysics (i.e., graupel melts into raindrops), enabling ASML to simulate the locations of the cold pools and the associated gust-fronts more accurately and faster (black boxes in Figure 6d) compared to CTRL. At 21:00 UTC, the cold pools of the Convection 2 and 3 have formed, and the cold pool of the Convection 1 has expanded further in ASML. The The positive q g increments in ASML could affect the cold-cloud precipitation microphysics (i.e., graupel melts into raindrops), enabling ASML to simulate the locations of the cold pools and the associated gust-fronts more accurately and faster (black boxes in Figure 6d) compared to CTRL. At 21:00 UTC, the cold pools of the Convection 2 and 3 have formed, and the cold pool of the Convection 1 has expanded further in ASML. The northeast cold-air outflow boundary of Convection 1 promoted the development of convection. Conversely, convergence regions of the observed MCSs were not accurately simulated by CTRL partly contributing to the poor simulations of these systems in this simulation.

Performance in the Forecast Period
During the forecast period, Convection 3 was gradually merged with Convection 1 in the following hours (Figure 5b). The merging process was well captured in ASML, and the simulation of merged MCS was consistent with the radar observation (Figure 5n). At 23:00 UTC, similar forecast improvements to those at 22:00 UTC were seen in ASML. On the contrary, only a few scattered, small convective cells near (18.  The impacts of the LDA scheme on the thermodynamic and dynamic environments were evaluated in Figure 8. In the locations of the observed MCSs, the temperature increments promoted the updrafts in the middle-upper levels, causing the averaged updrafts of Convection 1, 2, and 3 (black boxes in Figure 8b) in ASML were greater than 3 m/s (the maximum updraft to be around 7.9 m/s in Convection 2, not shown).. However, these updrafts failed to be simulated in CTRL (Figure 8a). Moreover, much stronger cold pools and the associated gust-fronts near (21.4 • N, 108.4 • E) and (21 • N, 112.3 • E) were simulated in ASML than CTRL, partly due to the q g increments-introduced cold-cloud precipitation processes (Figure 8c,d). Since the poor simulations of the initial convergence in the regions of observed MCSs, the corresponding MCSs and cold pools are not well simulated in CTRL.  maximum updraft to be around 7.9 m/s in Convection 2, not shown).. However, these updrafts failed to be simulated in CTRL (Figure 8a). Moreover, much stronger cold pools and the associated gust-fronts near (21.4°N, 108.4°E) and (21°N, 112.3°E) were simulated in ASML than CTRL, partly due to the qg increments-introduced cold-cloud precipitation processes (Figure 8c,d). Since the poor simulations of the initial convergence in the regions of observed MCSs, the corresponding MCSs and cold pools are not well simulated in CTRL.  (Figure 7b,k). Owing to the contribution of the cold pools and gust-fronts generated in the data assimilation period, the morphology and locations of the MCSs in ASML matched well with the radar observations (Figure 7b,n). At 23:00 UTC, ASML still outperformed CTRL Remote Sens. 2022, 14, 1965 13 of 20 in forecasting the morphology and locations of the MCSs (Figure 7c,l,o). Since the LDA mostly corrected the storm-scale background fields and had little impact on large-scale environments, simulation errors (e.g., location errors and spurious convection errors) were also presented in the forecast results of ASML.

Vertical Cross Sections of Radar Reflectivity Fields
The vertical cross sections of radar reflectivity through convective cores in Convection 3 of Case 1 and Convection 3 of Case 2 at the analysis time are shown in Figures 9 and 10. It is indicated that the simulated convective cells in ASML and the observed convective cells showed similar vertical reflectivity structure. For example, for the cross section shown in Figure 9, the heights of three convective cores were at~4.0 km in ASML, which were consistent with the observations and the cloud top heights simulated by ASML (~12.0 km, 14.3 km, and~13.5 km) were also close to the observation. However, the vertical ranges of the maximum reflectivity region showed differences between the observed and simulated results, especially for the cross section shown in Figure 10. It should be noted that the displacement errors in the simulation partly account for the differences.    Figures 11 and 12 show that ASML produced larger FSSs than CTRL at the analysis time and most of the forecast period, indicating an overall effective improvement of assimilating FY-4A lightning data. At the end of analysis time and the nowcasting period (1-2 h forecast), FSSs of CTRL were very small, due to the existence of the large intensity and displacement errors of the initially simulated MCSs. After the LDA, these errors were alleviated in ASML. In the late stage of forecasts, the differences between the two experiments gradually decreased, which was due to the LDA corrected only the convective-scale model states and had little impact on the large-scale environments, and the role of the large-scale environments became increasingly influential over time [36,81].  Figures 11 and 12 show that ASML produced larger FSSs than CTRL at the analysis time and most of the forecast period, indicating an overall effective improvement of assimilating FY-4A lightning data. At the end of analysis time and the nowcasting period (1-2 h forecast), FSSs of CTRL were very small, due to the existence of the large intensity and displacement errors of the initially simulated MCSs. After the LDA, these errors were alleviated in ASML. In the late stage of forecasts, the differences between the two experiments gradually decreased, which was due to the LDA corrected only the convectivescale model states and had little impact on the large-scale environments, and the role of the large-scale environments became increasingly influential over time [36,81].     For a more comprehensive quantitative evaluation of simulation results, we computed the FSSs for more neighborhood radii and radar reflectivity thresholds. The results are shown in Figures 11 and 12. In Roberts and Lean [80], a FSS greater than 0.5 is deemed "useful". For the threshold of 35 dBZ, the FSSs of CTRL hardly exceeded 0.5 when neighborhood radii were less than 36 km, while ASML increased the FSSs by 0.1-0.6 and exceeded 0.5 in forecast time when neighborhood radius was 36 km. The ASML showed improvement for all three radar reflectivity thresholds, indicating that the assimilation of FY-4A lightning data has a positive effect on the forecast of both convective and stratiform regions of MCSs.

Summary and Conclusions
In this study, the potential benefit of assimilating FY-4A lightning data on the forecast of non-typhoon oceanic convection was assessed with two oceanic convective case studies. To accomplish this work, an LDA scheme implemented in the WRF-FDDA was employed. With this LDA scheme, the process involves the calculation of the column-integrated graupel mass of each grid-column from the observed lightning rates from FY-4A LMI using an empirical formula between total flash rate and graupel mass [71]. Afterwards, three-dimensional q g fields are retrieved from the column-integrated graupel mass by using climatological profiles of q g for different integrated graupel mass bins. A distanceweighting spread method [73] is employed for the retrieved q g fields to consider the graupel contents in the adjacent regions of lightning observed points. The temperature increments are then derived from the latent heating rates which are proportional to the q g increments. The temperature and q g increments are then ingested into the model via the nudging terms.
Two oceanic convective weather cases over the South China Sea were selected to perform the study, and the simulations were compared subjectively and quantitatively with the radar observations to assess the performance of the LDA. The verification results demonstrated that the assimilation of FY-4A lightning data improves the performances on the analyses and forecasts of the oceanic MCSs. Most of the convective cells missed by the control experiments were recovered by using the LDA, due to the promoted updrafts by the temperature increments, and the more accurate and faster simulations of the cold pools and the associated gust-fronts.
The experiments with the LDA outperformed the control experiments in the forecast period because the cold pools and gust-fronts generated during the analysis period helped maintain the development of the MCSs, and reduced the morphology and displacement errors. The poor performance of the control experiments can be partly explained by the deviations in the initial convergence in the regions of observed MCSs, and the associated deviations in the locations and strength of cold pools and gust-fronts. The quantitative evaluation indicates that the most effective periods of the LDA for simulation enhancement were at the analysis time and the nowcasting (0-2 h forecast) periods.
The spatial resolution of LMI (7.8 km at nadir) was found to be lower than that of the cloud-allowing inner domain (3 km). Previous studies have found that the lower resolution of satellite-based lightning detection may lead to a larger spatial extent of MCSs compared to the observations [30,44]. However, benefits will still be obtained by assimilating such data, especially for marine areas not covered by high-resolution radar detection as indicated by the results of this study.
With the availability of lightning data with high spatial and temporal resolution from geostationary satellites, lightning data assimilation could provide a promising avenue to improve the short-term forecasts of oceanic MCSs. To make better use of the satellite-observed lightning data, future work should consider assimilating satellite-observed lightning data in combination with other observations and/or retrievals from geostationary satellites. Additionally, performance of assimilating FY-4A lightning data on the forecasts of oceanic MCSs with more advanced data assimilation methods (e.g., the ensemble Kalman filter or the hybrid ensemble variational method) is also subject to further study.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
This part provides the full names of the abbreviations used in this paper.