Next generation ground-motion prediction equations for Indo-Gangetic Plains, India

Seismic hazard analysis (SHA) of Indo-Gangetic Plains (IGP), India, and its proximity to the Himalayas requires reliable ground motion prediction equations (GMPEs). This study attempts to derive the GMPEs for IGP using strong-motion accelerometer data recorded from 2005 to 2015 on IGP. For regression, peak ground acceleration (PGA) and pseudo-spectral acceleration (PSA) of 5% damped linear pseudo-absolute acceleration response spectra at 27 periods ranging from 0.01–10 s were used. Two-stage nonlinear regression trains the functional form of nonlinear magnitude scaling, distance scaling, and site conditions. The model includes a regionally independent geometric attenuation finite fault distance metric, style of faulting, shallow site response, basin response, hanging wall effect, hypocentre depth, regionally dependent anelastic attenuation, site conditions, and magnitude-dependent aleatory variability. Developed GMPEs are validated for active crustal continental earthquakes for epicentral distances REPI ranging from 10–1500 km, magnitude ranging from 3.3–7.9, and focal depths 1–70 km. The GMPEs developed are compared with the Campbell and Bozorgnia 2008, 2013 and 2014, and North Indian GMPEs, which are agreed upon consistently. The model can predict the horizontal ground motion for SHA on IGP, considering possible maximum magnitude earthquakes from the seismic gaps of the Himalayas.


Introduction
The earthquake mechanism is very complex. As we know, the prediction of occurrences of earthquakes is uncertain. The best possible way to save life and property is by assessing the earthquake hazard and preparing for the risk of future Earthquakes. Ground motion prediction equations (GMPEs) play a crucial role in the seismic hazard assessment in a region. GMPEs are a mathematical expression that relates a strong-motion parameter like peak ground acceleration (PGA) and pseudo-spectral acceleration (PSA) of ground shaking to seismotectonic parameters of the earthquake source, the propagation path between the source and the site, also the surface and subsurface geological layers beneath the site. GMPEs are broadly used to predict the extent of ground shaking at different frequencies. In most regions, site-speciBc GMPEs are unavailable; therefore, the ground motions are predicted from GMPEs developed in regions where ground-motion data is available. The purpose of this study is to derive coefBcients of GMPEs for the Indo-Gangetic Plains (IGP) in India and reduce variability in ground-motion prediction models. In the IGP region, there is an insufBciency of recorded strong-motion data, so it has been a core reason for the difBculty in developing GMPEs for IGP India.
With the development of instrumentation and installation of a strong motion network, reliable past GMPEs are developed by Singh et al. (1996), Sharma (1998), Nath et al. (2005), Sharma et al. (2009), Iyengar (NDMA report 2010), Anbazhagan et al. (2013), Harbindu et al. (2014), and Raghukanth and Kavitha (2014). Most developed GMPEs are based on simulated data and very few on recorded data. Recently, a program for excellence in strong motion studies (PESMOS), Indian Institute of Technology Roorkee installed strongmotion accelerometers over IGP and proximity to the Himalayas , the data provided via the National Centre of Seismology Delhi. Many authors have studied attenuation and seismotectonic parameters based on recorded strong ground motion data and developed GMPEs for the northwest, northeast, and northern India region at bedrock levels (Singh et al. 2017(Singh et al. , 2021Harinarayan and Kumar 2020;Kumar et al. 2021;Ramkrishnan et al. 2021). Harbindu et al. (2014) developed GMPEs for a magnitude range of 3.5-6.8 and a distance range of 250 rock sites and it is based on synthetic ground motion. Raghukanth and Kavitha (2014) have used a total of 236 records from 62 earthquakes with magnitudes ranging from M 3.4-7.8 are used to develop the seismological model. Singh et al. (2017) used data from the 2015 Gorkha, Nepal earthquake (M 7.9) and some of its more significant aftershocks to derive GMPEs for earthquakes along the Western Himalayan arc. Kumar et al. (2021), Harinarayan and Kumar (2020) GMPEs are applicable for the magnitude (M) range of M 3.5-7.8, up to an epicentral distance of 250 km and for NEHRP site classes A/B, C, D, and E separately. Ramkrishnan et al. (2021) GMPEs applicable to a more extensive magnitude range of M 4.1-7.8 and distance range up to 1560 km. Singh et al. (2021) used recordings from the 2015 Gorkha, Nepal earthquake and some of its more significant aftershocks to derive GMPEs for earthquakes along the Western Himalayan arc.
The previous models have not considered all the seismotectonic parameters included in other worldwide models, such as Campbell andBozorgnia (2013 and. In the present study, all these seismotectonic parameters have been accounted. We consider the GMPEs developed in this study to be a reliable model for predicting PGA and PSA from Himalayan earthquakes on IGP. The GMPEs can be used for regional and local-level seismic hazard analysis for sites located on alluvial soil, soft rock, or hard rock conditions in the IGP and adjacent Shiwalik hills. This paper provides a brief description of the database, functional forms, and analyses used to develop GMPEs for IGP, followed by a comparison of the GMPEs model with Campbell and Bozorgnia (2008;hereafter, CB08) and Campbell and Bozorgnia (2014;hereafter, CB14). We also present practical guidance on how to use the model in seismic hazard analysis and engineering applications.

Indo-Gangetic Plains seismic hazard scenario
The Himalayas is a convergent plate margin with a high rate of strain accumulation. The nearby IGP consists of four significant sedimentary basins plains, such as the Ganga, Sindhu, Indus, and Brahmaputra (Valdiya 2016 (Ambraseys and Douglas 2004;Bilham 2004Bilham , 2008Bilham , 2019Grandin et al. 2015;Kubo et al. 2016). According to Bilham et al. (2001) (Khattri and Wyss 1978;Seeber and Armbruster 1981;Khattri 1987Khattri , 1993. The region between the 1934 Bihar-Nepal earthquake and the 1905 Kangra earthquake is known as the Central seismic gap (Singh et al. 2002) (Bgure 1); the length of this gap is about 750 km. In 1505, a catastrophic earthquake occurred in the Central seismic gap, which caused significant damage over the Ganga basin. Other events that have occurred in the Central seismic gap are the large (M * 7.0-7.9) 1803 and 1833 earthquakes Mukhopadhyay 2014, 2015) and the major (M * 6.0-6.9) 1991 Uttarkashi (M 6.8) and 1999 Chamoli (M 6.6) (Thakur and Sushil 1994;Thakur and Kumar 2001). However, because of their relatively small size, these are not gap-Blling earthquakes (Bilham 1995;Khattri 1999). In 2015, a major earthquake occurred in Nepal (M 7.9) close to the Central seismic gap. It was not a Central gap-Blling earthquake because the rupture propagated east-southeast from the hypocentre away from the seismic gap for about 160 km (Fan and Shearer 2015). Great earthquakes of magnitude M [ 8.0 in the Himalayan region are possible in the three identiBed seismic gaps. Therefore, there is an urgent need to develop new GMPEs for the IGP, which will be helpful in the assessment of seismic hazards and risks in this highly populated region.

Database
The Indian Institute of Technology Roorkee (IITR) has a Program for Excellence in Strong Motion Studies (PESMOS), comprising the Himalayan mountain range from Jammu and Kashmir in the west to Meghalaya in the east. PESMOS has installed about 290 strong-motion accelerometer (SMA) stations on the IGP. Recordings from these stations have been processed and provided by PESMOS (2005PESMOS ( -2014, and the SMA data for Nepal 2015 earthquake is recorded by PESMOS instrumentation and retrieved and provided by the National Centre for Seismology (NCS) of the Ministry of Earth Sciences, Government of India. The database contains accelerograms recorded over 694 stations on IGP, but for the preparation of GMPEs, the earthquake data of M C 3.0 and focal depth \70 km recorded from 2005-2015 are used and the remaining were removed from the database. The geometric mean of the peaks of the two horizontal component values of PGA and PSA are utilised for regression. The amplitude of the seismic waves at the recording station depends on the site geology, as well as the path considering the far-Beld propagation eAect coming from longer distances.
The selected strong-motion database consists of 515 recordings from 140 earthquakes of magnitude ranging from M 3.0-7.9 with epicentral distance (R EPI ) ranging from 1-1555 km and calculated rupture distance (R RUP ) ranging from 1-1495 km from India, Nepal, Bhutan, Myanmar, and Iran-Pakistan border. There are recordings at 259 alluvium sites, 59 soft rock sites, and 197 hard rock sites based on the site categorization done by Mittal et al. (2012) and Kumar et al. (2012). Figure 2 shows the earthquake locations and recording stations along the Himalayas and IGP, which have been used to develop the GMPEs. The earthquake database used in the analysis shows magnitude, focal mechanism solution, depth, and location. All the parameters have been provided in Supplementary tables 1 and 2. Figure 3(a) shows a plot between PGA and the closest distance to rupture (R RUP ). In Bgure 3(a), there are only four points where the PGA exceeds 0.1 g, which might be considered strong shaking. Similarly, a plot between the magnitude and closest distance to rupture is shown in Bgure 3(b); it explains the distribution of data with the magnitudes and rupture distances for all recording stations. It also shows that smaller numbers of strong motion data were recorded in the study region for the magnitude M [ 6 and rupture distance R RUP \ 100 km. With this limitation in the data, modelling of the hazard in empty space is troubled by uncertainty. We have to develop a model with this dataset because the strong motion recording in India is very low, and we need an updated ground motion prediction equation for this region. A plot of magnitude and the focal depth is shown in Bgure 3(c). A plot between magnitude and the timeaveraged shear-wave velocity in the top 30 m of the site (V S30 ) is shown in Bgure 3(d). It can infer that the different magnitudes have been recorded on different site classes. A plot of basement depth and time-averaged shear-wave velocity in the top 30 m of the site is shown in Bgure 3(e). The accelerographs used for data recordings are installed in a room on the ground Coor of a one-or two-story building or in the free Beld inside a small fabricated housing . The thickness of sediment deposits above basement rock (Z 2.5 ) is taken from GSI (2000).

Ground motion model
In this study, the Anderson and Hough (1984) method is used to Blter strong ground motion data. The referenced GMPEs approach in the functional form of Campbell and Bozorgnia (2014) has been used for ground motion model development. A Bnal functional form has been selected based on the availability of the data and the parameters necessary for engineering applications such as probabilistic seismic hazard analysis (PSHA). Many regression methods have been proposed to develop GMPEs; two commonly used regression methods are the two-stage method, Brst proposed by Joyner and Boore (1981) and reBned by Joyner and Boore (1993), and the random-effects method, Brst proposed by Preisler (1984, 1985) and popularized by Abrahamson and Youngs (1992). The two methods have been shown to give virtually identical results when appropriate weights are used in the second stage (Joyner and Boore 1993). In this paper, the twostage nonlinear regression method is employed to obtain the Bnal set of coefBcients and standard deviations. The two-stage analysis allowed us to partition the variability into within-event and between-event residuals and standard deviations. In the Brst stage, nonlinear regression is used to Bt the model for the recording-speciBc terms f dis , f atn , f hng , f site , and f sed . The result is a set of amplitude factors (event terms) and their standard deviations. In the second stage, a regression analysis is conducted on the weighted event terms to Bt the earthquake-speciBc terms f mag , f Ct , and f hyp . These terms are deBned in the following sections of the article.

Strong motion intensity measures
Strong-motion parameters represent a speciBc attribute of an earthquake's time history or its frequency-domain equivalent. To describe the active ground motion in engineering seismology, the most common peak time-domain parameters used are PGA, to some extent, peak ground velocity (PGV), and peak ground displacement (PGD). Peak time-domain parameters represent the maximum absolute amplitude of ground motion from a recorded or synthetic accelerogram, velocity, or displacement time series.
For the peak frequency-domain parameters, the most common are the response-spectral parameters PSA and pseudo-spectral velocity (PSV). As seismic design procedures have become more advanced, engineers have begun to consider the natural period of artiBcial structures in their design using a response spectrum. The component of the strong motion intensity measures (IMs) used in developed GMPEs is the geometric mean (GM) of the two recorded horizontal components. The IMs parameters for this study are PGA and PSA for 27 oscillator periods (T) ranging from 0.01-10 s. The plot of the geometric mean of response spectra at 5% damping with the time period for all recorded accelerograms from all earthquake events data has been shown in Supplementary Bgures 1 and 2.

Ground motion model formulation
The natural logarithm of PGA (g) and PSA (g) (Y) in the Brst-stage regression is Bt by the equation The event terms amplitude factors from the Brststage regression is Bt in the second stage by the equation The initial functional forms of the terms in equations (1a and 1b) are those in CB14, as deBned below.

Magnitude term
where M is the moment magnitude.

Geometric attenuation term
where R RUP (km) is the closest distance to the coseismic rupture plane.

Style-of-faulting term
where F RV is an indicator variable representing reverse and reverse-oblique fault types.
F NM is an indicator variable representing normal and normal-oblique fault types.
and k (°) is the rake angle deBned as the average angle of slip on the fault.

4.2.4
Hanging-wall term ; 5:5\M 6:5 where a 2 is period-dependent, theoretically constrained model coefBcients given in the result table 3, same as (Campbell andBozorgnia 2013, 2014), and R JB km ð Þ is the closest distance to the surface projection of the co-seismic rupture plane (Joyner-Boore distance).

Shallow site response term
where V S30 (m/sec) is the time-averaged shearwave velocity in the top 30 m of the site, and A 1100 (g) is the median predicted value of PGA on the rock with V S30 = 1100 m/sec (rock PGA).

Basin response term
where Z 2:5 (km) has been taken as the basement depth of 2.5 km/sec shear wave velocity.
4.2.7 Hypocentral depth term where Z hyp (km) is the hypocentral depth of the earthquake.

Treatment of missing data
In case of missing predictor variables for recordings in the database, they are either estimated using proxies. Missing values of V S30 are taken from the proxies deBned in Allen and Wald (2009) Thompson and Worden (2018).

Aleatory variability model
Aleatory variability is deBned in terms of both the geometric mean (GM) of the two horizontal components and the arbitrary horizontal component of motion.

Geometric mean horizontal component
Consistent with the two-stage regression analysis used to derive the median value of Y, the aleatory variability model for the GM horizontal component is deBned by the equation: where g i is the between-event (inter-event) residual for an event i and y ij ,Y ij , and e ij are the observed value, predicted value, and within-event (intraevent) residual for recording j of event i, respectively. The independent normally distributed variables g i and e ij have zero means and estimated between-event and within-event standard deviations on reference rock or on soil represented by linear site response, s and /, given by the magnitude-dependent equations where s i and / i are empirically derived standard deviations.
The total aleatory standard deviation is given by combining the between-event and within-event standard deviations by the square root of the sum of squares (SRSS) by the equation

Arbitrary horizontal component
The aleatory variance of the arbitrary horizontal components (Boore and Joyner 1997;Boore 2005;Campbell and Bozorgnia 2014) can be calculated from the equation where subscripts 1 and 2 refer to the two orthogonal horizontal components, j is an index representing the recording, and N is the total number of recordings. This equation is used to calculate the values of / c associated with the database that is used to develop the model. The total standard deviation corresponding to this component is given by the equation: where r is the total GM standard deviation shown in equation (21).

Results and discussion
The coefBcients determined in the regression are c 1 to c 17 . In the shallow site response term, the periodindependent theoretically constrained coefBcients c = 1.88 and n = 1.18, and the theoretically constrained period-dependent coefBcients k 1 and k 2 are taken from CB14. The theoretically constrained period-dependent coefBcient k 3 in the basin response term is also from CB14, although technically, it is not needed since it is modiBed by the regression coefBcient c 14 . It is included to be consistent with CB14. The coefBcients of the twostage regression analysis are summarized from c 1 -c 7 in table 1, from c 8 -c 14 in table 2 and from c 15 -c 17 and period-dependent parameters k 1 , k 2 , k 3 , and a 2 in table 3. The aleatory standard deviations are summarized in table 4.

Evaluation of residuals
To validate the derived GMPEs, we analyse the between-event residuals (g i ) and the within-event residuals (e ij ) for PGA and PSA at spectral periods of 0.05, 0.2, 1, 3, and 5 sec against selected predictor variables. It has been tested for the range of over and under prediction of PGA and PSA at different periods. The plot of residuals for PGA with different predictor variables is shown in Bgures 4 and 5. In these plots, a positive residual indicates underestimation by the model and a negative residual indicates overestimation by the model. The plot for between-event and withinevent residuals for the PGA and PSA at different response periods like 0.05, 0.2, 1.0, 3.0, and 5.0 s with different model parameters utilized for GMPEs preparation has been given in Supplementary Bgures 3-10. Figure 4 shows between-event residuals for PGA as a function of magnitude and hypocentral depth. Figure 5 shows within-event residuals for PGA as a function of magnitude, rupture distance, site velocity, sediment depth, and rock PGA. The plots show no significant trends or biases in the residuals that would indicate that the model is inconsistent with the data.

Model evaluation
Figures 6-12 show how the predicted ground motion scales with magnitude, rupture distance, spectral period, and site eAects. The values of the predictor variables used to compute the predicted ground motions are given at the top of each plot. Figure 6(a) shows the scaling of PGA with distance (attenuation) for magnitudes of 3.5, 4.5, 5.5, 6.5, and 7.5 for a strike-slip fault. Figure 6(b) shows a similar plot for PSA at T = 1.0 sec. The scaling of PGA with distance for strike-slip faults for developed GMPEs compared with CB08 and CB14 are shown in Bgure 7. For the smaller magnitudes and smaller distances, the predicted PGA values from our model are lower, and for larger distances, the PGA values match well with the CB08 and CB14. For a magnitude of M 6.5, our GMPEs match well, but for larger magnitudes at smaller distances, our predicted values are higher than CB08 and CB14 for reverse faults (Bgures 8 and 9). In the present study, the closest epicentral distance data is very less in number, so the scaling at short distances is not well constrained.
Scaling of PGA and PSA at T = 1.0 sec with magnitude for strike-slip faults for different rupture distances (5, 15, 25, 50, 100, 300, and 500 km) for our GMPEs. Figure 10(a and b) show a decreasing trend in PGA/PSA with longer distances and an increasing trend with higher magnitudes. Figure 11 shows the scaling of PSA with periods for different magnitude ranges for the GMPEs, at rupture distances of 10, 20, 50, 80, 150, and 300 km. There is a shift in the peak of the spectra from spectral periods of about 0.1-1 sec at short distances as magnitude increases from around 3.5-5.5, where the peak becomes relatively constant. There is also a noticeable shift at larger distances where the spectral peak shifts to longer periods at small magnitudes and broaden considerably at large magnitudes.  Figure 12 shows the variation of PSA with respect to the median estimate of bedrock PGA (A 1100 ) from 0.1-0.7g. It also clearly shows the nonlinear soil eAects and the associated reduction in spectral amplitude for the softer site conditions (NEHRP D and E) as rock PGA increases. The between-event, within-event, and total aleatory standard deviations for M B 4.5 and M C 5.5 for IGP GMPEs are shown in Bgure 13. Within-event standard deviations (SD) are higher than between-events SD for M C 5.5 and viceversa for M B 4.5. There is a relatively constant increase in between-event and within-event standard deviations from 0.01-0.4 sec, decreasing for higher periods. For both magnitude ranges, the significant increase in standard deviations due to a smaller number of digital recordings and some outliers present in residual for between-event and within-event impact the value for standard deviation. Comparing the developed equations for recorded site conditions V S30 = 900 m/sec hard Table 3. Derived coefBcients for the GMPEs from c 15 -c 17 and period-dependent parameters k 1 , k 2 , k 3 , and a 2 .

User guidance
In this section, guidelines are provided for potential users on evaluating the model for seismic hazard analysis.

Limits for applicability
Indo-Gangetic Plains, ground motion model is considered valid for shallow focus continental  earthquakes with focal depths from 0-70 km and for epicentre distance (R EPI ) ranges from 10-1500 km. The developed GMPEs can be used for seismic hazard analysis for the source rupture distances greater than 10 km in active tectonic regimes with high accuracy. Due to scarcity (small number of recordings) in the strong motion data within a 10 km epicentral distance, PGA does not match well with a small distance and magnitude (Bgures 7-9). Other limitations in terms of magnitude, faulting types, site velocity, sediment depth, and hypocentral depths are listed below: • The minimum magnitude that can be used in seismic hazard analysis with high accuracy for these GMPEs is M 3.3. • The maximum magnitude that can be used in seismic hazard analysis with high accuracy using these GMPEs is different for different faulting types of earthquakes. For strike-slip faults, these GMPEs are applicable for M B 7.0. For reverse faults, it is applicable for magnitude M B 7.9, and for normal faults, it is applicable for M B 7.8. Further, these GMPEs can be used cautiously after testing for certain upper limits like (M [ 8.0) of magnitude ranges. • Shear-wave velocities of V S30 = 150-1620 m/sec as the data used are recorded in this range, but for regression, V S30 velocity was used from the United States Geological Survey (USGS) V S30 model. • Sediment depths of Z 2.5 = 0-12 km.
The IGP model is uniformly valid over the entire range of predictor variables listed above. When the model is extrapolated beyond the data limits of the predictor variable, the errors can become large and should be used with caution under such conditions. The applicable range of some predictor variables has not been extended beyond the limits of the data. Additional details and explanations are given in the following sections.

Magnitude
Magnitude scales deBne the size of an earthquake and it is of many different types. The most common magnitude used in attenuation relations is moment magnitude (M w ), body-wave magnitude (M b ), surface-wave magnitude (M s ), and local magnitude (M L ). There is a propensity to adopt M as the standard for measuring magnitude because it is derived from seismological properties like seismic moment, which measures the earthquake energy radiation (Hanks and Kanamori 1979). In this study, magnitude M provided by the PESMOS database has been used for the derivation of the IGP model. The upper magnitude limit for strikeslip is 0.2 units smaller than the limit given, and the upper limits for reverse and normal faults are the same as that used in the database. For the application of derived IGP model the lower limit of earthquake magnitude can be taken as M C 3.3.

Geometric attenuation and anelastic attenuation
The ground motion amplitudes decrease as the seismic wave propagates away from the earthquake source to the site. The 'Source-to-site distance' is used to know the decrease in ground motion amplitudes. One of the basis of the assumption of the earthquake as point sources or as Bnite rupture sources is grouped into two categories: epicentral distance (R EPI ) and hypocentral distance (R hyp ), which are the distances measured from point source assumption to the site. Hypocentral distance is measured from the point within the Earth where the earthquake rupture initiated to the site. Epicentre distance is the distance between the Earth's surface points from the hypocentre just above the rupture surface to the site. Finite source distances are generally used in attenuation relations Brst deBned by Joyner and Boore (1981), the horizontal distance to the vertical projection of the rupture plane from the site (R JB ), second the distance between the seismogenic parts of the rupture plane to the site (R SEIS ). In this study, the rupture distance (R RUP ) is used as a Bnite fault source approximation of the closest horizontal distance to the vertical projection of the rupture plane to the site, and R JB is used to constrain the hanging wall term. Figure 17 shows the distance metrics used in the development of worldwide GMPEs. Geometric attenuation is constrained for 1 to 80 km for up to the magnitude M 6.8 because higher magnitude data are also recorded beyond 80 km in the database, and their applicability for lower limits might be beyond 15 km but considered minimum up to  10 km. The attenuation rate is anelastic for the earthquake recordings at large distances ([80 km) which are constrained by higher magnitude earthquakes in the database. For the IGP model development the recorded accelerograms are used for the rupture distances (R RUP ) up to 1495 km. To constrain the anelastic attenuation term, we consider the GMPEs to be most valid at R RUP B 1000 km, beyond which the number of recordings decreases rapidly.
For IGP, R RUP in the 1000 km buAer range from site to earthquake sources should be considered for PSHA studies. For source zone characterization for IGP, it needs to consider the far-Beld eAect due to deeper sediment deposits in the proximity of the Himalayas. Many large-magnitude earthquakes from the Himalayas aAect the IGP due to far-Beld eAects, and presently, no proper GMPEs model is there to consider these eAects. The lack of proper GMPEs for the IGP region limits us to characterize source zones in the range of 300 km for PSHA (Keshri et al. 2020). When the source zones for higher magnitude at a distance greater than 300 km, the unrealistic PSHA results from Himalaya's earthquake zones using GMPEs based on simulated strong motion data.

Style-of-faulting
The orientation of slip on the fault plane is known as the rake and dip of the fault plane, and the type of faulting or focal mechanism characterizes it. The type of faulting is typically classiBed into slip (horizontal slip), reverse (dip-slip with the hanging-wall side up), thrust (same as reverse but with shallow dip), and normal (dip-slip with the hanging-wall side down). Campbell (1981) empirically demonstrated that reverse and thrust-faulting earthquakes have relatively higher ground motions than strike-slip or normal-faulting earthquakes. It has been common practice that strike-slip and normal faulting events are in the same category, and reverse and thrust faulting events are in another category. For the study of faulting type, we have used the GCMT catalog focal mechanism solution for magnitude M C 5, and for magnitude M \ 5 focal mechanism solution of local source region data used from the seismotectonic atlas of India (GSI 2000). Due to the unavailability of a focal mechanism and fault plain solution for   smaller events, we have kept the style of the faulting term in equation (6) as zero for smaller events like magnitude M B 5. So, the style of the faulting term has been modelled for the magnitude M [ 5 by taking care of the scarcity of data. The indicator variables F RV and F NM represents reverse, thrust, normal, strike-slip, and oblique faulting by using dummy variable depending on the average slip of the focal mechanism solution for M [ 5.

Hanging-wall
For the study, those sites on the hanging wall exhibit more significant acceleration than those on the footwall. The hanging wall eAect increased the inCuence of the seismic dynamic loading on structures. For the IGP model, we have modiBed the NGA-West, 14 hanging wall term. Only the eAect due to rupture distance and magnitude were used for regression, and the other term was removed because of insufBciency of the faults data. For constraining the hanging wall, some of the terms related to R X (km), which is the closest distance to the surface projection of the top edge of the coseismic rupture plane measured perpendicular to its average strike, dip-related terms were dropped from regression analysis because these data were not present in the database. For this study, hanging wall side stations are decided based on using the event station latitude-longitude and seismotectonic map of India, and geological information on the Arc-GIS platform also along with the visualization of the GCMT focal mechanism solution. The hanging wall side and foot wall side stations have been marked as dummy variables '1' and '0', respectively, in the database. For the earthquake of low magnitude M \ 5.5, the f hng;M term is zero in equation (11), so the hanging wall term f hng is zero in equation (9). As the slip term is already removed for modelling the hanging wall eAect, so the GMPEs model, the eAects are only based on magnitude and distance.

Shallow site response
Local site conditions are deBned in terms of nearsurface geology, shear-wave velocity, and sediment depth; it describes the earthy material lying directly beneath the site. Usually, local site conditions are categorized as soil and rock (soft or hard) types. Various GMPEs use this soil classiBcation, which would be further divided into shallow, soft, Brm, and very Brm soil. Moreover, in recent developments, the rock could be further divided into soft and hard rock (Campbell and Bozorgnia 2008; shear wave velocity V S30 is used for site parameters. In this study, we have used the V S30 for local site conditions as alluvium soil site class C, 180-375 m/sec for soils (alluvium soil (slope washed)), soft rock (sandstones/slates/limestones/dolomites), class B, 375-700 m/sec for soft to Brm rocks and hard rock (quartzite, dolomites, schist, granodiorites, and gneiss), class A, 700-1620 m/sec for Brm/hard rocks. In our selected IGP database, the V S30 data were calculated using the USGS model, which characterizes the recording sites in 180 B V S B 900. We caution our GMPEs for site category for V S30 B 180 because no single data is used for this model. The upper limit for V S30 is 900 m/sec of site class A, but data may be recorded in the range of 700-1620 m/sec. The shallow site term in the IGP model approximates an empirical estimate of site response for general site classiBcation, and it is recommended. It may be extrapolated with caution, and it is recommended to go for site response analysis.

Basin response
For the IGP model, the basement rock depth of sediments (Z 2.5 ) is used for the basin response to see how the basins aAect the ground motions. The IGP model may be valid for the basement depth of 10-12 km, shallow basement depth, and small V S30 . The estimate of basin ampliBcation is not predicted well enough because of the jelly-like behaviour of sediment and basement resonance. In this case, site response analysis is essential.

Hypocentral depth
There is a high strain accumulation rate in the Himalayas convergent plate margins, which releases from time to time in the region. The main fault types for the Himalayan region are mainly reverse/ thrust faults, and it comes in the category of shallow-crustal earthquakes in active tectonic areas. Sources in the Indo-Gangetic regions' deformation rates are generally low and localized in intraplate regions. This region will be highly inCuenced by the far-Beld eAect of the Hindukush earthquakes and the shallow thrust earthquakes in the central Himalayas and the subduction zone in the northeast region. In this study, the shallow focus earthquake is considered (Z hyp \ 70 km).
Based on an analysis of residuals, hypocentre depth shows its dependence on between-event residuals. Figure 4 (b) shows that the between-event residuals were generally unbiased down to shallow focus earthquake depth. We recommend that our IGP model not to be used for hypocentral depths or depths to the top of rupture greater than 70 km. This depth limit constrains the modelled earthquakes to occur in the shallow lithosphere up to 70 km. This model is not applicable for the Hindukush earthquakes, intermediate focal depth earthquakes (70-200 km).

Conclusions
In the absence of proper GMPEs in the IGP region, Earthquake Strong Motion recordings from the seismically active Himalayan belt recorded on 259 alluvium sites, 59 soft rock sites, and 197 hard rock sites of IGP is used to establish new GMPEs for the IGP region. The data set consists of 515 strong-motion seismic recordings from 140 earthquakes of magnitude ranging from M 3.0-7.9 with rupture distance (R RUP ) ranging from 1-1495 km. An attempt has been made in the present paper to develop a new reference GMPEs (IGP) same as CB13 and CB14. In previous models for northern India GMPEs, minimal model parameters were used. In the present study, period-dependent magnitude saturation, magnitudedependent style of faulting eAects, scaling with hypocentral depth, fault dip, geometric attenuation, regionally dependent anelastic attenuation, hangingwall eAects, shallow linear and nonlinear site response, shallow sediment response, and magnitudedependent nonlinear between-event and withinevent aleatory variability have been taken as model parameters to develop GMPEs based on two-stage nonlinear regression analysis.
Though strong motion data for the Indian region is less to develop GMPEs, the result shows that the explanatory variables explain the variation in PGA/ PSA. PGA/PSA values were calibrated with previous GMPEs for northern India and worldwide for the Btted model. It is also calibrated for the PSA at different periods for the past 1999 M 6.6 Chamoli earthquake for two different sites. For the alluvium and hard rock sites, the predicted and observed values are highly correlated. Further, the future IGP model can upgrade for recorded and simulated ground motion parameters from different Himalayan regions. GMP21 will help to predict ground motion parameters for the computation of seismic hazard risk on the IGP and the adjacent Siwalik region.