NOWCASTING EARTHQUAKES IN THE NORTHWEST HIMALAYA AND SURROUNDING REGIONS

With the rapid increase and availability of seismic data, an automatic, transparent and regular way of earthquake hazard estimation strategy is highly desirable in many seismically active large geographical regions. In this paper, we implement a novel method of nowcasting (Rundle et al., 2016) that can indirectly assess the current progression of a region through its earthquake cycle of large events. Nowcasting differs from the method of forecasting in which future earthquake probabilities are calculated. Using statistics of natural times, counts of small earthquakes between large earthquakes in a defined region, nowcasting provides an earthquake potential score (EPS) to enable scientists and city planners a snapshot of the current level of earthquake hazard in the region. Applied to a number of selected major cities in the northwest Himalaya and surrounding regions, we found that the EPS values corresponding to 6 M  events in New Delhi, Chandigarh, Dehradun and Shimla reach about 0.56, 0.87, 0.85 and 0.88, respectively. These estimated scores thus indicate that New Delhi is about half-way through its cycle for magnitude 6.0 or higher earthquakes, while Dehradun is about 85 percent of the way through its cycle. Towards the end, we discuss some implications and applications of these nowcast values to improve the present earthquake hazard assessment practice in the study region.


INTRODUCTION
In a seismically active large geographical region like northwest Himalaya, faults are almost everywhere (Scholz, 1990;Yin, 2006).Some faults occur at the Earth's surface, whereas some others may be buried enough to be escaped from standard geological mapping.Some faults based on their size and potential are marked as active, whereas a few of them are assumed to be relatively benign (Yin, 2006).Yet, often times, many damaging earthquakes occur in unexpected places, where neither records of previous large earthquakes exist, nor any major geological faults have been identified so far.In such situations where geodetic or geophysical observations are notoriously limited to yield much information on current and future hazard of a region to the driven complex system of earthquakes, statistical analysis of classifying seismicity and related empirical hazard estimation can prove to be effective (Jordan, 2006;Pasari and Dikshit, 2015).
In this study, we aim to implement a new method of nowcasting (Rundle et al., 2016(Rundle et al., , 2018) ) to statistically work out hazard of large earthquakes at current time in a few selected cities embedded in the larger northwest Himalaya and adjacent regions.The term "nowcasting" usually appears in meteorology, economics and political sciences where the goal is to assess the uncertain current state of a dynamic system in the present, the very near future and the very recent past (Rundle et al., 2016).The .* Corresponding author nowcasting idea in statistical seismology refers to the determination of the progression of a region through its regional earthquake cycles.In other words, the nowcasting method uses proxy earthquake data to provide a simple yet transparent estimation of stress and strain accumulation in a region since the last large event (Luginbuhl et al., 2018).In the past, this determination of tectonic stress regime has focused on direct or indirect, physical or remote observations at selected places, its relation to failure strength of active faults (through estimated slip rate) and finally measuring the state of earthquake hazards in a region (Scholz, 1990).Such a program of study, many a times, is not only expensive but also requires complex geodetic inversion (Rundle et al., 2016).The nowcasting approach thereby evolves to serve a hierarchy of roles to bridge the gap between the detailed physical modeling of a complex dynamical system and the pressing societal need of earthquake hazard estimation in large seismic prone areas.
The nowcasting idea in seismology is grounded on two innovative concepts.First, it considers the formulation of "natural" times, rather than the traditional calendar or clock times, to develop the underlying probability distribution of observed seismicity in a large geographic region.Second, it assumes "earthquake cycle" as the recurring events in a region comprising several mapped or unmapped faults, rather than the usual focus on repeating earthquakes on some selected active faults.The natural times, which are the counts of small earthquakes between successive large earthquakes in a defined region, are independent to the background seismicity level, as long as the b-value in the Gutenberg-Richter (G-R) frequency-magnitude relation turns out close to a constant (Rundle et al., 2016).As a result, the nowcasting technique is equally applicable to smaller as well as larger seismic regions provided we have sufficient data to span at least ~20 large earthquake cycles in the region (Rundle et al., 2016(Rundle et al., , 2018)).Neither declustering of catalog nor homogeneity in earthquake magnitude is necessary in nowcasting method (Luginbuhl et al., 2018).The crucial trick in nowcasting lies in that it accumulates natural time statistics from larger spatial region, such as the entire northwest Himalayan region, and then apply the same statistics to smaller areas, such as cities, embedded in that larger region.This assumption seems realistic in light of the ergodicity property in the statistical mechanics of earthquakes, in which the statistics of smaller regions over longer times is considered to be similar as the statistics of broader regions over large spatial domains and longer periods (Rundle et al., 2016(Rundle et al., , 2018)).Using the best-fit data-derived probability distribution of observed natural times in a region, the earthquake potential score (EPS) or the nowcast value is computed through the cumulative probability     P X n t  of the current number of small earthquake counts   nt at (calendar) time t.The numerical EPS scores derived in this way have several direct or indirect applications, such as the hazard based relative ranking of cities, scenario earthquake risk modeling, revision of insurance premiums and most importantly, nationwide earthquake risk reduction and resilience (Rundle et al., 2016(Rundle et al., , 2018)).

STUDY AREA AND EARTHQUAKE DATA
The present study area covers the northwest portion of Himalayan orogen, part of alluvial Indo-Gangetic plains, east-northeast part of Pakistan and its surrounding regions bounded by the latitudes 26 0 N to 36 0 N and longitudes 72 0 E to 84 0 E (Figure 1).It is a region of active tectonics and high seismicity comprising a series of major thrust faults (e.g., main central thrust, main boundary thrust, Himalayan frontal thrust) and their imbricate branching faults (e.g., Jawalamukhi thrust, Palampur thrust, Barsar thrust, Panjal thrust, and Nalagarh thrust), ridges (e.g., Delhi-Haridwar ridge; Delhi-Lahore-Sargodha ridge), duns (e.g., Dehradun, Pinjore dun, Kota dun), uplifts, anticlines (e.g., Janauri and Sarkaghat), synclines (e.g., Mussoorie) and lineaments (Yin, 2006).From north to south, the study region covers a part of the snow-capped Tethys Himalaya, Higher Himalaya, Lesser Himalaya, Siwalik, and the foreland basin of Indo-Gangetic plains (Yin, 2006).The area has experienced a number of devastating earthquakes in the recent past, notably the 1803 Garhwal  (Yin, 2006;Pasari, 2018).The built-in infrastructure and the houses in these regions are hardly equipped with earthquake resistance, although the population has dramatically increased due to the rapid pace of urbanisation and growing dense city population (Pasari and Dikshit, 2015).A number of developing cities and towns including the national capital of India, New Delhi, national capital of Pakistan, Islamabad, and other Indian state capitals, such as Chandigarh, Dehradun, Jaipur, Jammu/Srinagar, Lucknow and Shimla lie in the study region.Economy in these cities is thriving, although most of the regions are highly vulnerable to earthquakes.occurred in the study region, among which 27 earthquakes had magnitude 6.0 or above.Therefore, we obtain 26 interevent small earthquake counts between large events.However, in 1975, two large earthquakes occurred on 19 January in the study region almost simultaneously without any small events between them.Therefore, we obtain 25 positive interevent counts of small earthquakes (natural time) for further analysis.

METHODOLOGY
In the driven complex system of earthquakes or any other power-law systems, such as typhoons, landslides and market economy crashes, the occurrence of events are associated with a frequency-magnitude distribution, indicating that several smaller magnitude events will occur between the powerful, large magnitude events in a defined geographic region (Scholz, 1990;Rundle et al., 2016).However, such interevent counts (natural times) of small events are random.Therefore, developing a statistical distribution of natural times is the key idea in nowcasting analysis.In other words, for a given set of observed random sample 12 , , , n X X X of natural times, we would like to develop a probability density function (pdf) or cumulative distribution function (cdf) to assign a unique numerical score to the selected geographical region at a given time to describe its current progression in earthquake cycle on the basis of elapsed small earthquake counts since the last large event in the region.
The methodology for the main body of nowcasting technique therefore comprises three main steps: description of candidate probability distributions, parameter estimation and selection of best-fit distribution.Of course, the descriptive measures of the natural times, such as histogram, mean, standard deviation, coefficient of variation, skewness and kurtosis of the sample dataset will guide for inferential statistics in many ways (Pasari and Dikshit, 2015).
In this study, we employ ten candidate probability distributions from an exhaustive set of timeindependent, time-dependent, heavy-tailed and exponentiated family of distributions to fit the observed natural times for nowcast value calculation.These distributions are exponential, gamma, lognormal, Weibull, inverse Gaussian, inverse Weibull, Maxwell, Pareto, exponentiated exponential and exponentiated Weibull distributions.Each of these distributions offer significant modeling benefits in statistical seismology (Pasari andDikshit, 2014, 2015;Pasari, 2018).For instance, time-independent exponential distribution provides constant hazard function, whereas gamma, Weibull and exponentiated exponential offer increasing, decreasing and constant hazard functions; lognormal, inverse Gaussian, inverse Weibull and exponentiated Weibull distributions provide monotone as well as nonmonotone hazard shapes, whereas the heavy-tailed Pareto distribution offers only decreasing hazard shapes (Pasari and Dikshit, 2014, 2015, 2018).
For estimating parameters of the studied distributions, we use the method of maximum likelihood estimation (MLE), whereas to select the best-fit distribution, we utilize two goodness-of-fit tests, namely the Akaike information criterion (AIC) and the Kolmogorov-Smirnov (K-S) minimum distance criterion.Once the best-fit distribution is identified, the EPS score will be established using the cumulative distribution evaluated for the current number of small events in selected "local" circular city regions within a radius of 300 km around their respective city centers.These "local" regions are of course embedded in the larger northwest Himalaya and adjacent regions.

RESULTS AND DISCUSSION
Before we turn to natural time statistics and nowcast evaluation, it is imperative to develop the G-R frequency-magnitude plot   10 log N a bM  of earthquakes to estimate the b-value for the study region.
It is also important to assess the magnitude dependency between the cumulative number of small and large earthquakes in the study region (Rundle et al., 2016(Rundle et al., , 2018a)).For this, let us denote the small and large earthquakes by    1, it is clear that the Weibull and exponentiated distributions provide the best fit for the observed natural times, although the fit from exponential, gamma and exponentiated Weibull models cannot be ignored.A deeper insight of the overall fitting analysis from the K-S graph further highlights the closeness of these five competitive models with the dataderived empirical distribution function.
After developing natural time statistics and selecting the best-fit distribution(s) for the study region, we now focus on nowcast calculation for a few selected cities, namely New Delhi, Chandigarh, Dehradun and Shimla.
For this, we tabulate earthquakes of magnitude greater or equal to 4 M   within a circle of 300 km around their geographic city centers (Figure 1).We assume that the natural time statistics of these "regional" circular regions are indifferent from the statistics of the entire study region (Rundle et al., 2016(Rundle et al., , 2018)).It is observed that for the New Delhi circular city region, a total of 145 events of 4 M  have occurred since 1963, among which there are two large events.The current number of small events in this region reaches 70, since the last large event on 28 March 1999 at Chamoli, about 290 km away from Delhi.Using the best-fit Weibull distribution, we found that the nowcast score as on 04 August 2018 for the New Delhi region turns out to be 0.56, indicating that the Delhi city region has progressed about 56% of the way to the next large   6 M  earthquake, since the last such event in 1999.
We repeat the above process to compute earthquake potential score (EPS) for other city regions.The results are summarized in The EPS scores in Table 2 show that the nowcast values for the Himalayan cities are consistent and they are significantly higher than that of the New Delhi region.Apart from the determination of current progression of the region, the emanated nowcast values can be used for ranking of cities as to their current exposure to seismic risk of large earthquakes and other seismological applications such as in insurance sectors, engineering and safety divisions, mineral exploration, cultural resource management or geospatial mapping in the study regions.Nowcasting coupled with earthquake forecasting thus may lead to an enhanced earthquake hazard management system in the study region.
The EPS score     P X n t  , which provides the nowcast value for a region, is nothing but a continuous and monotonically non-decreasing function that resets to 0, immediately after every large earthquakes in the selected region.The nowcast values or the current state of a regional fault system may remain in the same standard for many years, whereas for some active regions, the EPS values may change rapidly over a short span of time.Therefore, by and large, the nowcasting method does not provide any direct implication to the future seismic activity, although the concept of natural time may provide some latent information of the earthquake dynamics in a selected region (Rundle et al., 2016;Luginbuhl et al., 2018).

SUMMARY
In any power-law system of complex dynamics like earthquakes, a fixed (average) number of frequent small event counts can be associated with one large event to define a strategy to determine the current progression of a region towards the infrequent next powerful event (Rundle et al., 2016).Of course, this criterion will only be fulfilled in a statistical sense, given the natural variation in small event counts (Rundle et al., 2016(Rundle et al., , 2018)).Using statistical concepts, these small event counts can then be summarized in a simple yet effective manner to assess the current state of system in terms of their nowcast values.
In this paper, we have successfully implemented the nowcasting process for earthquakes to indirectly determine the current state of earthquake hazards in a few seismically active city regions in the northwest Himalaya and adjoining regions.The underlying statistics of natural times are developed using a series of ten probability distributions: exponential, gamma, lognormal, Weibull, inverse Gaussian, inverse Weibull, Maxwell, Pareto, exponentiated exponential and exponentiated Weibull.Using goodness-of-fit tests, it was found that the best fit comes from Weibull and exponentiated exponential distributions.Applying the best fit distributions, we found that the nowcast values corresponding to 6 M  events in New Delhi, Chandigarh, Dehradun and Shimla regions reach about 0.56, 0.87, 0.85 and 0.88, respectively.These nowcast values, like thermometer readings, enable scientists and city planners a snapshot of their current exposure to earthquake hazards, which otherwise remain elusive from limited, sparse geodetic or geophysical observations on "selected" faults.
As a summary, the nowcasting idea in the present analysis offers an automatic, transparent and regular way of earthquake hazard estimation strategy in the seismically active northwest Himalayan region.Involving multiple directions from statistics, earthquake physics and seismology, this study not only contributes to the development of innovative methods in statistical seismology but also enhances industry-academia collaboration in large-scale utility mapping, smart city planning, disaster risk reduction and seismic resilience.

Figure 1 .
Figure 1.Earthquake data in the northwest Himalaya and adjacent regions; selected cities are highlighted along with their circular city regions The earthquake data for the present nowcasting analysis is acquired from the global public seismic catalogs, such as the Advanced National Seismic System (ANSS) Composite Catalog and the International Seismological Centre (ISC) catalog.Due to an insufficient instrumental coverage and weak signals from deep seated earthquakes in the mountainous regions, we confine ourselves to 4 M  events having a maximum focal depth of 200 km (Rundle et al., 2016, 2018).In addition, as our concern in this study is on small interevent counts, we concentrate only on the instrumental catalog starting from 1963.It is observed that during 06 March 1963 to 04 August 2018, a total of 2270 events   4 M  occurred in the study

Figure 2 .
Figure 2. (a) G-R magnitude-frequency plot for the study region, (b) Dependence plot between the cumulative number of small earthquake counts and large earthquake counts in the study region It is observed that the sample natural times of 25

Figure 3 .
Figure 3. K-S graph for the most competitive probability distributions, namely exponential, gamma, Weibull, exponentiated exponential and exponentiated Weibull distributions Table 2 below.