Contrasting lightning projection using the lightning potential index adapted in a convection-permitting regional climate model

Lightning climate change projections show large uncertainties caused by limited empirical knowledge and strong assumptions inherent to coarse-grid climate modeling. This study addresses the latter issue by implementing and applying the lightning potential index parameterization (LPI) into a fine-grid convection-permitting regional climate model (CPM). This setup takes advantage of the explicit representation of deep convection in CPMs and allows for process-oriented LPI inputs such as vertical velocity within convective cells and coexistence of microphysical hydrometeor types, which are known to contribute to charge separation mechanisms. The LPI output is compared to output from a simpler flash rate parameterization, namely the CAPE ×\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\times$$\end{document} PREC parameterization, applied in a non-CPM on a coarser grid. The LPI’s implementation into the regional climate model COSMO-CLM successfully reproduces the observed lightning climatology, including its latitudinal gradient, its daily and hourly probability distributions, and its diurnal and annual cycles. Besides, the simulated temperature dependence of lightning reflects the observed dependency. The LPI outperforms the CAPE ×\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\times$$\end{document} PREC parameterization in all applied diagnostics. Based on this satisfactory evaluation, we used the LPI to a climate change projection under the RCP8.5 scenario. For the domain under investigation centered over Germany, the LPI projects a decrease of 4.8%\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$4.8\%$$\end{document} in flash rate by the end of the century, in opposition to a projected increase of 17.4%\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$17.4\%$$\end{document} as projected using the CAPE ×\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\times$$\end{document} PREC parameterization. The future decrease of LPI occurs mostly during the summer afternoons and is related to (i) a change in convection occurrence and (ii) changes in the microphysical mixing. The two parameterizations differ because of different convection occurrences in the CPM and non-CPM and because of changes in the microphysical mixing, which is only represented in the LPI lightning parameterization.


Introduction
Current lightning climate simulations mainly rely on parameterizations, which relate climate model output to observed lightning (Clark et al. 2017), but rarely closely reflect the physical mechanisms leading to lightning. Indeed the complexity, the spatio-temporal scales, and the diversity of mechanisms related to lightning flashes do not allow global climate models for its explicit representation.
As an example, the mechanism dominating thunderstorm electrification is known as the non-inductive charging mechanism (Reynolds et al. 1957;Takahashi 1978;Saunders 1993;Saunders and Peck 1998;Latham et al. 2007). It implies electric charge separation through rebounding collision between small ice crystals growing by water vapor diffusion and graupel pellets growing by accretion of supercooled water droplets (Yair 2008). The sign and the magnitude of the charge exchanged between ice crystals and graupel pellets depend on many parameters (e.g., temperature, cloud water content, rime accretion rate, etc. -Jayaratne et al. 1983;Saunders 2008). Once the charge is exchanged, the ice particles are dispersed in the clouds through gravitational processes, with heavy graupel pellets moving downward relatively to the ice crystals (Williams 1988;Houze 2014). On the macro scale, the cloud will, therefore, be characterized by multiple large cloud layers with homogeneous polarities (Stolzenburg et al. 1998a, b, c).
The non-inductive process is mainly occurring in areas with convective activity. Therefore, many parametrizations relate the occurrence, and sometimes the intensity, of the convective activity to derive flash densities (Finney et al. 2016). An example of such a parameterization is the convective available potential energy times the precipitation rates ( CAPE × PREC ) used to derive the lightning flash rate. With these simple and easily accessible variables, it is possible to obtain both the location where flashes are likely to occur, but also the intensity of the flash activity as high CAPE and high precipitation rate indicate high convective activity. Over the contiguous United States, Romps et al. (2014) found that this parameterization explains 77% of the temporal variance in the lightning flash rate. Another approach commonly used for parameterizing lightning in climate models is based on the observed relationship between the cloud top height (CTH) and the flash rates (Price and Rind 1992), the cloud height being a good indicator for the occurrence and the intensity of convective activity. While such a relationship compares well to current observations, it is unknown whether it will hold in a different climate. For example, while CTH is expected to increase in a warmer climate, so does the cloud base height (Chepfer et al. 2014), which may result in constant cloud depth. Yoshida et al. (2009) used cold cloud height to tackle this issue, but all these parameterizations are still highly dependent on the empirical relationship and may be far away from representing the electrification process, notably the non-inductive one.
Closer to the theoretical mechanism, some more advanced lightning parameterizations have been implemented by introducing variables accounting for the microscale processes. For example, the ice mass flux scheme (Finney et al. 2014(Finney et al. , 2018 uses the upward cloud ice flux at 440 hPa. Although the variables selected in more complex parameterizations are uncommon GCM output, they allow for clear improvements both spatially and temporally, over parameterization essentially based on the convective activity (Finney et al. 2014). Moreover, while most of the convective activitybased parameterizations shows an increase in lightning flashes in a warmer climate (e.g., 12 ± 5% per • C increase for CAPE × PREC (Romps et al. 2014) or from 7 to 13% per • C increase for CTH (Krause et al. 2014;Clark et al. 2017)), the first flash rate climate projection based on the ice mass flux scheme parameterization (Finney et al. 2018) show decreasing lightning toward the end of the century following the RCP 8.5. These results are promising for the climate community, but the ice mass flux scheme still suffers from significant inadequacies (Finney et al. 2018). First, accounting for the ice flux on a single level is still far from representing the non-inductive charging mechanism. Second, the variables from which the ice flux scheme is based on, are still strongly parameterized in current GCM.
For example, the variability of hydrometeors or that of convective up/down-drafts have scales that are currently far from being explicitly represented in global climate models (Prein et al. 2015). The updraft is particularly critical as it conditions hydrometeors' transport, the growth of graupel through riming, and the chances for collisions between ice crystals and graupel pellets (Houze 1993). Microphysical processes will remain parameterized in the next generation of global climate projections. Recently, Lopez (2016) developed a new lightning parameterization for the global convection-parameterizing NWP model IFS of ECMWF. It employs diagnosed vertical profiles of graupel and supercooled liquid from within the convection parameterization and CAPE to form an estimator for the flash rate. While physically more sound, its possible performance in current coarse-grid GCMs is not clear.
However, the first convection-permitting projections produced by regional climate models are now available (Prein et al. 2015). These models show realistic representations of convective processes in general (Leutwyler et al. 2017). On the hourly and sub-hourly scales, convection-permitting models (CPM) improved the representation of spatial precipitation patterns and extremes values (Prein et al. 2013;Brisson et al. , 2018. The representation of convective clouds also shows improved features Keller et al. 2016). Interestingly, the sensitivity of convective related variables derived from CPM to global warming does not always match that derived from coarser climate models, as observed by, e.g., Kendon et al. (2014), Vanden Broucke et al. (2019, Helsen et al. (2020).
Besides, the improved representation of convective processes allows for more complex parameterizations (Fierro et al. 2013). In McCaul et al. (2009, the flash rate is estimated based on the resolved upward flux of graupel in the mixed-phased region at −15 • C and the integration of solid hydrometeors over the full storm depth. This parameterization is related to that established latter in Finney et al. (2014) but fits closely to the non-inductive process as it focuses on graupel and not anymore on the full ice content. Furthermore, it is applied at a specific temperature level instead of a pressure one. One step further into the complexity is the lightning potential index (LPI). The LPI accounts not only for the presence of solid hydrometeors but also for the presence of both solid (with a stress on graupel) and supercooled liquid hydrometeors. Besides, it uses all atmospheric levels ranging from −20 to 0 • C. Proposed by Yair et al. (2010), this parameterization was applied to several case studies with CPMs. The LPI shows a high correlation with observed lightning Lynn and Yair 2010;Lagasio et al. 2017). To the authors' knowledge, the output of such parameterization has only been applied to short periods for numerical weather prediction purposes. With the rise of multi-decadal CPM simulations, such parameterization could now be applied on a climate time scale. According to Kendon et al. (2016), the climate community expects an improved representation of lightning with CPMs and that further developments/research are needed in that direction.
This study aims in this direction. While many paths are open, we choose to start with the following questions: Can more complex parameterizations such as the LPI provide a realistic representation of lightning? What is the added value of these parameterizations applied to CPMs when compared to simpler parameterization used in non-CPMs? Do the projections of flash rates as represented by such complex parameterizations differ from that of simpler non-CPM parameterizations? To provide elements of answers to these questions, the lightning potential index (LPI) will be implemented and applied in a CPM and compared to the CAPE × PREC parameterization applied to a non-CPM. The models and their parameterizations are detailed in the methods section, together with an observational dataset. Such an indirect comparison of a CPM with one lightning parameterization and a non-CPM using another might appear to be not ideal at first glance. The fact is that LPI is not easy to implement in a non-CPM (lack of explicit convective updrafts, supercooled liquid, and graupel), while CAPE × PREC has severe problems in a CPM (CAPE-removal by explicit convection leads to a null CAPE where there is PREC and vice versa). The validation of the parameterizations against observations and an investigation on flash rate projections are described in the result section. The added value of a parameterization like the LPI compared to current GCM parameterizations is discussed but is not the main focus of this paper. Finally, the conclusion provides elements of answers to the research questions listed above.

Model setup and lightning parameterizations
All simulations investigated in this study were performed using the Consortium for Small-scale Modelling in climate mode (COSMO-CLM) model. The COSMO-CLM model is a non-hydrostatic limited-area climate model, based on the COSMO model (Steppeler et al. 2003;Doms et al. 2011), a model designed by the Deutsche Wetterdienst (DWD) for operational weather predictions. The climate limited-area modeling (CLM) community adapted this model to perform climate projections (Böhm et al. 2006;Rockel et al. 2008). We used the version COSMO5.0 clm7.
In all COSMO-CLM simulations performed in this study, we use the 5 th order Runga-Kutta split-explicit time-stepping scheme (Wicker and Skamarock 2002), the lower boundary fluxes provided by the TERRA model (Doms et al. 2011), the radiative scheme after Ritter and Geleyn (1992), and the one-moment microphysical scheme (Steppeler et al. 2003;Doms et al. 2011). Besides, as recommended by , in the CPM, the one-moment microphysical scheme predicts the mass evolution of graupel, in addition to the four standard hydrometeor types (i.e., cloud droplets, raindrops, cloud ice, snow). Such configuration is close to the operational convection-allowing limited area configuration COSMO-DE of DWD (Baldauf et al. 2011). Finally, in the CPM parent nest, the wind's horizontal components, together with the temperature, are spectrally nudged to reduce boundary effect.
The finest nest is centered over central Germany (Fig. 1) at a resolution of 0.025 • ( ∼ 3 km). A double-step one-way nesting strategy is applied to reach this resolution. The global fields, i.e., ERA-Interim for the hindcast simulation , and the EC-Earth ensemble member r12i1p1 for both the historical period  and the RCP8.5 projection (2071-2100), are used to drive a 0.22 • ( ∼ 25 km) grid mesh COSMO-CLM simulation over a limited-area centered over Europe. In a second step, the resulting hourly outputs of this simulation are used to drive the 0.025 • COSMO-CLM simulation. While the 0.22 • simulations are performed with a parameterization of deep and shallow convection following Tiedtke (1989), in the 0.025 • simulations, only the shallow (non-precipitating) part of that convection parameterization is activated, allowing for the explicit simulation of convective cells with their associated strong up-and downdrafts and microphysical mixed-phase precipitation formation processes. In the innermost convection-allowing COSMO domain, we adopt the LPI parameterization of Lynn and Yair (2010) and Yair et al. (2010) as described in Sect. 2.1.1. This parameterization uses the local microphysical conditions for charge separation in updrafts as a necessary condition for lightning, in a relatively simple and computationally cheap way.
There are more accurate schemes in the literature, which explicitly predict charged hydrometeors, the dynamics of atmospheric charge distribution, and associated lightning discharge when the electric potential between any two points in space rises above a threshold, e.g., Barthe and Pinty (2007), Barthe et al. (2012). However, such schemes would have been much too computationally expensive for the present study.

The lightning potential index (LPI) parameterization
This parameterization is based on the LPI parameterization proposed by Lynn and Yair (2010) and Yair et al. (2010) was adapted and implemented into the COSMO model by one of the authors (U. Blahak) in 2014. Subsequently, it has been evaluated in an NWP context by Lopez (2016) and found to be a useful COSMO diagnostic for lightning forecasting. The LPI is defined as and represents the kinetic energy of the updraft -derived from the vertical wind speed -scaled by the potential for charge separation inherent to the microphysical mixing. The latter is given by the function with where q c , q r , q i , q s and q g are the mixing ratios of cloud water, rain water, cloud ice, snow and graupel, respectively. q L is the sum of the liquid species, while q F represents (weighted) contributions of the frozen hydrometeors. In essence, is the ratio of the geometric to the arithmetic mean of two quantities. It also appears twice in q F . takes a value of 1 when the liquid water and the solid hydrometeors masses are equal to 0 when all the water is in the same phase (e.g., solid or liquid). These values are vertically averaged over the primary charging zone, ranging from 0 • C to −20 • C . Moreover, non-zero values can only be attained if -besides the supercooled liquid water -graupel, and a second non-rimed species (in our case q i or q s ), are present at the same time, representing the non-inductive charge separation mechanism. In addition, the boolean function restricts the vertical averaging to grid boxes having updrafts above a threshold of 0.5 m/s .
To prevent spurious and unrealistic LPI signals in areas of weak, isolated updraft columns as well as in deep orographic gravity wave clouds (lenticularis), two additional filter functions f 1 and f 2 are applied. These functions are not directly related to the parameterized physical mechanisms of charge separation in the LPI and are briefly described in Appendices A and B.
The LPI has to be calibrated by observations to be an estimator for the flash rate.

The CAPE × PREC parameterization
The convection-permitting model's performance using LPI in representing lightning is compared to that of the coarser 0.22 • model. Applying the LPI parameterization is impossible in this non-CPM due to the lack of strong updrafts and graupel hydrometeors. Another parameterization, namely the CAPE × PREC product as a flash rate estimator, is applied here. This parameterization performed well in representing flash rates over the contiguous United States (Romps et al. 2014) and is the only parameterization which could be derived from the available output of our 0.22 • simulation.
As mentioned before, please note that the CAPE × PREC parameterization can hardly be directly applied to a CPM, as CAPE is rapidly depleted at locations of convective cells by latent heat release in moist adiabatic updrafts and evaporative cooling in precipitation below the cloud base. Therefore, a direct comparison of both parameterizations in the CPM framework is rarely applied in this study.
To a lesser degree, this depletion also happens in the Tiedtke parameterization of deep convection. Precipitation is produced simultaneously as the instability is gradually removed from the atmosphere leading to a lower CAPE in the grid column where convective precipitation is produced. Thus, multiplying the hourly CAPE with the hourly accumulated precipitation as output by the model results in relatively low CAPE × PREC values. To tackle this issue, the CAPE at time t was multiplied by the precipitation accumulated from time t to t + 1 h . This method was selected as it shows the closest fit to observed flash rates when compared to other methods (e.g., using daily average, using precipitation accumulated from time t − 30 min to t + 30 min).
Similarly, The CAPE × PREC values have to be calibrated by observations to be an estimator for the flash rate.

Validation dataset
The lightning observations were extracted from the BLIDS dataset (Blitz Information Dienst Siemens (2019)). The raw spatial resolution of the BLIDS data ranges from 200 m to 700 m, and the temporal resolution is 1 ms. The BLIDS network measures lightning since 1992 and improved several times since then (Schulz and Diendorfer 2004). Besides, the observation sensors were upgraded in 1999, resulting in an improved detection of strokes, especially intra-clouds ones (Schulz and Diendorfer 2004). This upgrade led to temporal inhomogeneities when integrating over the full period (i.e., before and after the upgrade of the sensors). To avoid these inhomogeneities, the model simulations are evaluated for the period 2000 to 2013 only. Besides, the evaluation domain was restricted to areas for which the German regional agencies TLUG (Thüringen), LAU (Sachsen-Anhalt), HLNUG (Hessen), and FAWF (Rheinland-Pfalz) provided us with the BLIDS dataset over each corresponding region.

Adjustment of parameterizations output
The observations and the output of lightning parameterizations do not match directly in terms of units, spatial and temporal scales. This mismatch prevents the comparison of the parameterization output directly with the observations. This issue is tackled by remapping the observations to the CPM grid using a conservative remapping. Besides, the flashes are grouped into 15-min bins for the LPI parameterization and into hourly bins for the CAPE × PREC parameterization by deriving the number of flashes occurring during a given time interval (i.e., ±7.5 min for the LPI, which was output every 15 minutes and ±30 min for the CAPE × PREC, which was output every hour).
Finally, we converted the parameterization's output to flash rates in a two-step procedure. First, the lowest values of the parameterization's output are set to 0, so that the occurrence of values above zero is adjusted to that of the observation. Second, a simple linear model is built that relates the observed flash rates to the parameterization output. The use of a linear adjustment is motivated by the high correlation found in Romps et al. (2014) suggesting a linear relation for the CAPExPREC parameterization. Besides, such linear models are not likely to result in overfitting with only two parameters to calibrate, having a large sample size. The adjusted values are therefore derived based on the following equations: where X is the parameterization output, and a, b, and c are three coefficients, which were calibrated using the accumulated observed flashes rates for each grid-point and the parameterization output from the ERA-Interim driven simulations. For the period 2000-2013, a = 0.405 , b = −0.073 and, c = 0.33 for the LPI parameterization and a = 0.000172 , b = 0.0118 and, c = 15.85 for the CAPE× PREC parameterization provide the best adjustments (i.e., the lowest root mean square error between the observed and simulated flash rate probability density). In addition, note that the adjustment does not alter the climate change signal described in Sect. 3.2. The latter is of similar amplitude for both the LPI and the LPI adj .
As defined above, the X adj provides a number of flashes per grid-cell. To avoid a dependency of the flash rate to cell size, two strategies were adopted. First, the unit flash rate per km 2 per time period is adopted through all this manuscript. Second, the evaluation was mostly performed using spatial averaging over the full domain. While this may limit the analysis of high local values (e.g., flashes from local storms), it strengthens the analysis' robustness. Besides, it allows for a fair comparison between the two parameterizations, which are applied to different spatial resolutions.
The X adj is written as X unless specified otherwise in this manuscript to ease the reading.

Skill scores
The performance of the model is evaluated using the mean square skill score (MSSS, Murphy (1988)), which is derived as where s and o respectively represent the simulation and observational datasets, n the sample size, and w the weighting function. In this manuscript, the weighting function is set to one unless specified otherwise. The MSSS is equal to 1 when the simulation is identical to the observation and is greater than 0 when the simulation is more realistic than the observed climatology ( ō).
In addition, the probability density functions (PDFs) are evaluated using the Perkins skill score (PSS, Perkins et al. 2007), which measures the common area between two PDFs by calculating their cumulative minimum for each group of binned values, where i is the bin index, n is the total number of bins, and Z s (i) and Z o (i) are the simulation-and observation-probability masses, respectively. Two PDFs are identical when the PSS is equal to 1. Besides, Perkins et al. (2007) suggests that two PDFs differ significantly from each other when their PSS is lower than ∼ 0.7.

Uncertainty and significance
The flash rate uncertainty is estimated through bootstrapping with replacement. This bootstrapping is performed by randomly selecting blocks from the original time series to create a new time series of similar length, hereafter referred to as bootstrap. Blocks of one month are selected to keep the temporal dependency of the data. Also, an equal number of each month is selected in each bootstrap (e.g., n January, n February, n March, etc.), to account for the seasonal variability. The uncertainty is derived using 100 bootstraps by taking the range between the 2.5th and 97.5th percentiles of the studied quantity. These bootstraps are also used to assess whether two simulations significantly differ from each other. For each of the two sets of 100 bootstraps, the studied quantity is derived (e.g., the average value for a given month). The differences between all pairs from these two sets are then computed, resulting in 10,000 values. A significant change at the l level is detected when the l/2 and the 100 − l∕2 percentiles of these 10,000 values are respectively greater or smaller than 0 such as Figure 2 compares the probability densities of the modeled daily lightning flash rate computed from the two parameterizations with the observations. Both parameterizations fairly reproduce the observed hourly and daily flash rate probability densities over the domain where observations are available. Still, the LPI outperforms the CAPE × PREC resulting in PSSs of respectively 0.93 and 0.77 for the hourly timescale and 0.84 and 0.73 for the daily timescale.

Evaluation
The observed flash rate is spatially heterogeneous (Fig. 3a). High flash rates are observed along the Rhine valley ( ∼ 9 • E , ∼ 49 to 52 • N ), as well as an area over near the southeastern corner of the evaluation domain while lower flash rates are found in the north and between 10 and 12 • E . A negative north to south gradient is also observed. The small scale spatial variability is poorly reproduced by both parameterizations (MSSS of −0.14 for the LPI and -0.25 for the CAPE × PREC). However, the north to south gradient is well captured by LPI with an MSSS of 0.61 (Fig. 4). The CAPE × PREC parameterization produces a north to south gradient that is too steep, resulting in a bias larger than the observed climatology, as indicated by the MSSS of −1.66.
One of the main challenges in non-CPM simulations is the representation of the convective activity's diurnal cycle. CPMs generally improve the representation of the diurnal cycle of convective precipitation (Prein et al. 2015). Here, the diurnal cycle of the LPI follows closely that of the observation with a minimum hourly flash rates in the morning and a maximum in the mid-afternoon (Fig. 5 -  parameterizations. The shaded areas represent the uncertainty derived through the bootstrapping technique with replacement introduced in the method section peak with too early occurrence of flash rates. The LPI largely corrects these biases with an improved timing of the midafternoon peak that falls within the sampling uncertainty. The difference in performance between the two parameterizations is similar to that found for precipitation when comparing CPM and non-CPM (Ban et al. 2014;). The investigations described above focus on the flash rate parameterizations performance in modeling the current climatology of lightning. However, limited information is provided on the model's robustness to provide reliable flash rate climatologies for different climates. As an indication of robustness, we studied the scaling of flash rates with 2 m temperature. As shown in Fig. 6 the temperature scaling of the observed flash rates can roughly be decomposed into three parts, each of them characterized by different exponential growth. The two related breakpoints occur at about 273 K and 294 K. The LPI reproduces well the observed temperature scaling with similar exponential growth rates and breakpoints. However, the model overestimates the mean daily flash rate for the lowest temperatures and simulated a drop of flash rates for the highest temperature (higher than 301 K) as opposed to observed rates resulting in an MSSS of 0.18. However, for these highest temperatures mean rates, the sampling size is small, as shown by the rug in Fig. 6. When weighting the MSSS with the sampling size of each bin, the MSSS increases to 0.83. A robust estimation of this possible misbehavior of the model would require a longer sampling period.
The CAPE × PREC is worse at reproducing the observed temperature scaling than the LPI resulting in a weighted MSSS of 0.57. Still, CAPE × PREC reproduces the tipping points and the general behavior of the temperature scaling. It should be noted that for the lowest and highest temperature bins, no flashes are produced by the CAPE × PREC as opposed to the LPI. Generally, the observed frequency of days with flash occurrence ( 61.8% ) is underestimated by both parameterization, with a larger underestimation for the CAPE × PREC ( 43.5% ) compared to that of the LPI 49.7%.

Climate projections
In this subsection, the evolution of the LPI and CAPE × PREC parameterizations with time under the RCP8.5 scenario is investigated. For this, the parameterizations are applied to COSMO-CLM driven by a historical (1975 to 2005) and a future (RCP8.5 from 2071 to 2100) simulation. These simulations are used to compare the climate change signal derived from the two parameterizations. While the wording projection is used in the following section, these simulations do not robustly indicate the evolution of lightning under a warmer climate for Germany, which would require a larger simulation ensemble. When comparing the future with the historical simulations, the two parameterizations show significantly different changes at the 5% level with a decrease of −4.8% for the LPI and an increase of +17.4% for the CAPE × PREC. Note that the changes for both parameterizations are independent of the latitude (Fig. S1 in supplemental material). Figure 7 shows that the climate change signal is not equally distributed on the full range of flash rates. For most of the flash rate range, no climate change signal is found. In contrast, the probability for high flash rates (i.e., flash rates above the historical 95th percentiles) increases for both parameterizations on both the daily and the hourly timescale. This increase is more substantial for the CAPE × PREC than for the LPI (e.g., the significant change at the 5% level of the 99.9th percentiles is 64% for the CAPE × PREC against 16% for the LPI).
As shown by the representation of the annual cycle (Fig. 8), the highest flash rates are mostly occurring in summer. For both parameterizations, the highest increase in flash rate is occurring in July. However, significant changes between the historical and the RCP8.5 simulations are found for other months (May for the LPI and January, February, October, and December for CAPE×PREC). The changes are of opposite signs with a decrease for the LPI and an increase for the CAPE×PREC.
While the annual cycles for the two historical and RCP8.5 simulations are relatively similar to each other (Fig. 8), the representation of the diurnal cycle by these two sets of simulations largely differ ( Fig. 9) with differences throughout the day except in the early morning with both parameterizations showing an increase in flash rate. The increasing occurrence of large-scale driven convection, which tends to occur homogeneously throughout the day in central Europe, or the occurrence of convective events that favors flash rate occurrence could explain this increase. The development of resampling techniques (e.g., resampling based on weather types), which are outside the scope of this study, may help  Fig. 6 The dependency of domain-averaged daily mean flash rates with daily temperature for the period 2000-2013. The observation (green squares) and the parameterizations (orange circle for the LPI and blue triangle for the CAPE × PREC) are binned to temperature classes (i.e., one value per 1 K class). In addition, the corresponding bars provide the sampling uncertainty derived from 100 bootstraps with replacement when the sample size is larger than one. The rug at the bottom of the plot shows the simulated temperature values for the investigated period test these hypotheses. During the afternoon, the two parameterizations disagree on the projection with a significant decrease for the LPI and a significant increase for the CAPE × PREC. This disagreement is observed for most months featuring a mid-afternoon peak in flash rate (Fig. 10).
The increase of temperature between the historical and the RCP8.5 simulations is 3.3 K (in both the CPS and the 0.22 • simulations) for the full year and of similar amplitude for the convective season (here taken from April to September included − 3.3 K). Assuming that the temperature scaling described in Sect. 3.1 stays unchanged in the RCP8.5 scenario, the temperature increase found in the RCP8.5 scenario should result in increased flash rates. This assumption mostly holds for the CAPE × PREC resulting in the total increase described in previous paragraphs. In contrast, the LPI scaling with temperature significantly differs in the RCP8.5 simulations (Fig. 11) compared to that in the historical simulation with lower flash rates for a given temperature. Most of the significant changes occur for temperatures ranging from 278 K to 294 K. This decrease in the scaling occurs mainly in the afternoon, while the decrease is limited to a few bins in the morning (Fig. S2 in supplemental material).     . 11 The dependency of the mean daily flash rate with daily mean temperature for the COSMO EC-Earth historical driven simulation (green) and the EC-Earth future driven simulation (purple). The rugs at the bottom and the top of the plot show the occurrence of a given temperature for the corresponding simulation. Flash rates are binned to temperature values (i.e., one point per 1 K). These values are plotted for both the LPI (left column) and the CAPE×PREC (right column) parameterizations. Stars indicate significant differences between the two simulations at the 5% level. This significance is only indicated for points for which the sample size exceeds 100 days Similarly, the CAPE×PREC parameterizations can be decomposed into different factors. It notably shows a high degree of similarity with the LPI no ; it incorporates a filter, through precipitation occurrence, and an indicator of the updraft intensity, through both precipitation intensity and CAPE. For the two selections of three years, the CAPE and precipitation variables were generated by both the 0.22 • RCM and the CPM. The occurrence of precipitation's accumulation greater than 1.5 mm/hour , a threshold that results in an equivalent modeled number of precipitation events as events with flashes for the evaluation period, is decreasing by 3.8% (Table 2) in the CPM projection. For these events, the mean CAPE increases by 7.8% in the CPM projection, resulting in a CAPE × PREC increase of 8.7% . In opposition to the CPM, in the 0.22 • RCM projection, the precipitation occurrence above 1.5 mm/hour shows an increase of 4.0% . Besides, the CAPE, for the events above 1.5 mm/hour , shows a more substantial increase in the 0.22 • RCM projection compared to that in the CPM, with +37.9% change resulting in a change of CAPE × PREC of +31.3% . This result suggests that the CPM, as opposed to the 0.22 • RCM, is partly responsible for the disagreement observed in Section 3.2 between the LPI and the CAPE × PREC.
Therefore, our results support the use of a CPM, as well as accounting for the microphysical mixing for projecting flash rates. Accounting for the microphysical mixing makes the LPI one of the most physically consistent flash rate parameterizations applied on climate scale at the convection-permitting scale. Interestingly, it provides a different climate change signal than the CAPE × PREC in this study. Generally, the present study agrees with previous findings. Indeed, while the literature considering convective activity as a proxy for lightning suggested a global increase in lightning activity under a warmer climate ranging from 5 to 10% increase per degree of warming (Price and Rind 1994;Price 2009;Krause et al. 2014;Romps et al. 2014;Rädler et al. 2019), the few studies using a scheme, which include microphysics related parameters showed little increase or even a global decrease in lightning activity (Jacobson and Streets 2009;Finney et al. 2018Finney et al. , 2020. However, the literature also highlights heterogeneous lightning changes around the globe. For example, for the area under investigation in the present study, Krause et al. (2014) find a small decrease in lightning despite considering convective activity as a proxy for lightning. Besides, the use of a single simulation and, therefore, single CPMs/GCMs significantly weakens the robustness of the projections computed in this study. Large ensembles of simulations that encompass multiple CPM and GCM combinations would be needed to validate the divergence of the two types of lightning schemes (i.e., scheme with or without microphysical parameters) under a warmer climate.
While the LPI may sound more advanced than a more simple parameterization, and therefore, more trustworthy, it is still suffering from some weaknesses. First, although being more physically consistent than other parameterizations, the LPI is still not fully representing the cloud electrification. Thus, more complex and more realistic parameterizations are needed. Second, the microphysics parameterizations still suffer from fundamental limitations in climate models even at convection-permitting scales. The use of different microphysics schemes, in addition to different lightning parameterizations, may be required to provide more robust conclusions as well as realistic projection ensembles of flash rates for the future.

Conclusion
Representing lightning in climate models is still a challenge. While many parameterizations exist, they strongly rely on the representation of convective cloud properties. Since about a decade, CPMs are applicable to climate time scales. These CPMs allow lightning parameterizations closer to process understanding than coarser-grid climate models that rely on a deep convection parameterization. In this study, such a lightning parameterization, namely the LPI, was implemented in the regional climate model COSMO-CLM. This model configuration was applied to perform three simulations on a convection-permitting scale. The first simulation uses ERA-Interim as driving data to adjust the LPI to observed flash rates and an evaluation at different spatiotemporal scales. The additional two additional simulations use EC-Earth global climate simulations as driving data allowing for the investigation of the LPI evolution under the RCP8.5 scenario. The output of the LPI parameterization for these simulations was compared to another flash rate parametrization, namely the CAPE × PREC, that was derived for non-CPM corresponding simulations.
Our results show that the LPI parameterizations reproduce the present-day daily and hourly probability densities, the latitudinal dependency, and the diurnal cycle of the observed flash rates. While the CAPE × PREC parameterization produces correct representations of the temporal probability densities and the latitudinal dependency, these representations are not as skillful as those derived with the LPI. Furthermore, the diurnal cycle is poorly represented  torical 1975torical -2005.5 2071-2100) indicates a climate change signal in mean flash rate by −4.8% as projected by the LPI and by +17.4% as projected by the CAPE × PREC. Both parameterizations generate more frequent extremely high flash rates (99.9th percentile), but this increase is more substantial for the CAPE × PREC. Both historical and future simulations show similar annual cycles with higher flash rates in summer than in winter. However, the LPI projects a decrease in flash rates while the CAPE × PREC shows a change of opposite sign. This decrease is related to a decrease in flash rates produced during midafternoon convection. Finally, the increase of flash rate found for the CAPE × PREC can be related to the temperature increase as the CAPE × PREC temperature-scaling of flash rates remains identical for both simulations. For the LPI, the temperature scaling is altered in the RCP8.5 simulation showing possible limitations of temperature based rescaling techniques for downscaling flash rates.
The investigation of potential sources of flash rate reduction revealed that both the use of a CPM and accounting for microphysical mixing changes are likely causes for the projected disagreements between the two flash rates parameterizations. Changes in the updraft velocity partly compensate for the changes induced by the microphysical mixing. These findings fit with the projected change in flash rate as simulated by the CAPE × PREC for which the increase in the vertical velocity is carried through the increase in CAPE and precipitation accumulation, but does not have any factors that can carry the signal of the microphysics mixing. Besides, applying the CAPE×PREC at the convection-permitting scale does not result in a convergence of the projected flash rate with the LPI. This absence of convergence points at the influence of the CPM on the differences observed between the two flash rates parameterizations. Quantifying this influence is an important question that still has to be answered. The use of a broader range of parameterizations applied in both multiple CPM and non-CPM model may help to provide elements of answers to this question.
Besides, applying the CAPE×PREC at the convectionpermitting scale does not result in a convergence of the projected flash rate with the LPI. This absence of convergence points at the influence of the CPM on the differences observed between the two flash rate parameterizations. Quantifying this influence is an important question that still has to be answered. The use of a broader range of parameterizations applied in both multiple CPMs and non-CPMs may help to provide elements of answers to this question. Furthermore, using a single CPM/RCM driven, by only one global climate model and one reanalysis, of one parameterization, and for a specific region does not allow for deriving robust lightning climate projections. Still, this study shows that explicitly accounting for the microphysics in flash rate parameterizations as well as representing explicitly and more robustly convective processes using convection-permitting models potentially improve current climate projections of lightning.
A filtering of noisy weak LPI signals ( f 1 ) As already noted in Yair et al. (2010), weak and noisy LPI signals caused by isolated single-grid-column updrafts may occur in km-scale models, which do not fully resolve smallscale convective updrafts. Because such situations do not really represent physically coherent convective updrafts, valid positive LPI values are restricted to grid-points for which the majority of grid columns in a certain horizontal neighborhood exhibits a maximum updraft speed above a certain threshold ). This criterion is adapted in the boolean filter function f 1 in Eq. (1) with the area fraction a in a 10 × 10 km 2 neighbourhood with column maximum updrafts larger than a threshold w (max,0) , w (max,0) is set to 1.1 m/s in our application, a value which was found to produce reasonable spatial LPI distributions in comparison to observed flash rates.

B Filtering of false LPI signals in strong orographic gravity wave clouds ( f 2 )
The filter f 1 is not enough to prevent spurious LPI signals in deep orographic wave clouds, which have been noted during the implementation work, for example, during Föhn events in the Alps in winter time. Large-amplitude gravity waves embedded in moist flow may lead to small but spurious graupel formation in cloud microphysics schemes, which in turn leads to false LPI signals. To prevent this, a second boolean function f 2 has been added, which is based on the fact that gravity waves require a rather stable stratification throughout the troposphere. f 2 is defined as the column-integrated buoyancy in a 20 × 20 km 2 neighbourhood: with where p is the pressure, p s the surface pressure, R d the gas constant of dry air, T v,parcel the virtual temperature of a moist adiabatic parcel ascent starting from 50 hPa above ground with average properties of the lowermost 100 hPa and T v,s local virtual temperature. B ML is formally similar to mixed layer CAPE, but with fixed integration over a 500 hPa layer starting at 50 hPa above ground. B ML is approximately 0 or slightly negative at locations of explicitly simulated convective cells, but attains large negative values for the stable conditions associated with orographic mountain waves. The threshold value and B 0 = − 1500 J/kg 2 was found by experimentation to reasonably separate these two regimes for our CPM domain over Central Europe.
Acknowledgements This research has been partially funded by TLUG, Thüringen, by LAU, Sachsen-Anhalt, by HLNUG, Hessen, by FAWF, Rheinland-Pfalz, by the French National Research Agency under the future investment program ANR-18-MPGA-0005, and by the EUCP project which is funded by the European Commission through the Horizon 2020 Programme for Research and Innovation (Grant Agreement 776613). The computational resources have been provided by the LOEWE-CSC.
Funding Open Access funding enabled and organized by Projekt DEAL.

Declarations
Data availability The data and materials used in this study are available upon request to the corresponding author except from the observed dataset and the COSMO-CLM model source code. The latter can be obtained by becoming a member of the CLM community.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.