India plate angular velocity and contemporary deformation rates from continuous GPS measurements from 1996 to 2015

We estimate a new angular velocity for the India plate and contemporary deformation rates in the plate interior and along its seismically active margins from Global Positioning System (GPS) measurements from 1996 to 2015 at 70 continuous and 3 episodic stations. A new India-ITRF2008 angular velocity is estimated from 30 GPS sites, which include stations from western and eastern regions of the plate interior that were unrepresented or only sparsely sampled in previous studies. Our newly estimated India-ITRF2008 Euler pole is located significantly closer to the plate with ~3% higher angular velocity than all previous estimates and thus predicts more rapid variations in rates and directions along the plate boundaries. The 30 India plate GPS site velocities are well fit by the new angular velocity, with north and east RMS misfits of only 0.8 and 0.9 mm/yr, respectively. India fixed velocities suggest an approximate of 1–2 mm/yr intra-plate deformation that might be concentrated along regional dislocations, faults in Peninsular India, Kachchh and Indo-Gangetic plain. Relative to our newly-defined India plate frame of reference, the newly estimated velocities for 43 other GPS sites along the plate margins give insights into active deformation along India’s seismically active northern and eastern boundaries.

Global Positioning System (GPS) measurements at sites located in the Indian plate and along its seismically active and tectonically complex boundaries provide critical, present-day kinematic boundary conditions for a wide variety of geo-scientific studies, including studies of the tectonics and earthquake cycles of the numerous faults within and along the boundaries of the India plate. Plate tectonic estimates of India plate motion, which are based on reconstructions of the Carlsberg Ridge, seafloor spreading center magnetic anomaly sequences, strikes of transform faults and earthquake slip vectors 1-6 span millions of years and are thus more useful for studying the long-term evolution of Indian plate tectonics than its present-day deformation.
Previous estimates of India plate motion and deformation internal to the plate with space geodetic observations include 7 inversion of the velocities of 12 episodic GPS (eGPS) sites located south of 14°N (South India) and Delhi to estimate India plate angular velocities relative to the geodetic reference frame ITRF96 and the Eurasia plate. Their India-Eurasia angular velocity predicted convergence rates along the Himalayan frontal thrust are ~14% slower when compared to the rates predicted by the NUVEL-1A 3-Myr-average India-Eurasia angular velocity. This is consistent with a previous report 8 of a significant discrepancy between geodetic and geologic estimates of India plate motion. Subsequent estimates of India plate motion have relied increasingly on velocities from continuous GPS sites (cGPS) and eGPS which include: 3 cGPS,ITRF97 9 ; 2 cGPS + 4 eGPS located in South Nepal, ITRF2000 10 ; 3 cGPS + 2 eGPS located in South Nepal, ITRF2000 11 ; 7 cGPS, ITRF2000 12 ; 10 cGPS + 2 eGPS, ITRF2000 13 ; 4 cGPS, ITRF2005 14 ; and 13 cGPS, ITRF2008 15 . All of the above studies are affected to varying degrees by the sparse distribution of GPS sites in the plate interior and often-short observation time spans at many of their GPS sites. In particular, previous authors use few or no GPS sites from western and eastern region, thereby limiting the accuracy with which the India plate rotation pole could be determined.
In this study, we use data from 27 cGPS and 3 eGPS sites from stable regions of the India plate (more than any previous study), including for the first time cGPS sites from western and eastern regions of continental plate interior, to estimate the India plate angular velocity in ITRF2008. The data samples more of the plate interior sites and span more time than any previous study and thus offer the best opportunity to date to establish an upper limit on any deformation internal to the plate. For the first time velocities from seven sites in western India are used for estimation of angular velocity. cGPS observations from another 43 sites located in the Indian subcontinent are used to quantify deformation along the seismically active northern and eastern boundaries of Indian tectonic plate.

India plate tectonic setting
The India plate borders the Eurasia plate on its northern and eastern boundary; Arabian plate on its, western boundary; Somalia, Capricorn, and Australia plates to the south (Fig. 1). Relevant to this study, the ~2500-km-long Himalayan Arc at the northern limit of the India plate accommodates NNE-directed convergence between India and Eurasia, increasing from ~35 mm/yr at the arc's western end in Kashmir (~75°E) to ~50 mm/yr at its eastern end near the Eastern Himalayan Syntaxis 4 . The western, sinistral shear-dominated plate boundary of Indian tectonic plate consists of the Owen fracture zone and Murray Ridge between India and Arabia 16,17 and the Chaman fault zone between the Makran subduction zone and Himalayan Arc 18,19 . Motion between India and the Sunda Block 20 (which rotates clockwise with respect to Eurasia) along their dextral, shear-dominated boundary is variously partitioned between convergence along the Andaman and Arakan trenches and structures farther inboard, including the Andaman back-arc spreading center and the Sagaing fault of Myanmar 10 . The 2500 km Himalayan Arc is characterized by several thrust faults that sole into the basal detachment of the Himalayan wedge or the Main Himalayan Thrust which marks the upper boundary of under-thrusting Indian plate [21][22][23][24] . Main Frontal Thrust (MFT) marks the southern boundary of the Himalayas and the surface trace of MHT coincides with MFT. Further north (Figs 2 and 5) of it is Main Boundary Thrust (MBT) which marks the southern limit of Lesser Himalaya, Main Central Thrust (MCT) transition between the Lesser and Higher Himalaya, South Tibet detachment (STD) is the transition between Higher and Tethyan Himalaya and Indus Suture Zone (ISZ) is the northern boundary 25,26 . Major seismic activity (Fig. 2) in the Indian tectonic plate is mostly confined to northern and eastern convergence zones 27,28 . Seismic activity in peninsular India, except Kachchh region (Figs 2 and 3a), is low with very rare occurrence of earthquakes of magnitude 5.0 29 . Kachchh region recorded several felt earthquakes 30 throughout 19 th century and the notable earthquakes in this region are 1819 Allah bund earthquake 31 , 1956 Anjar 32 and 2001 Bhuj earthquake 33,34 . The seismic events in peninsular India are isolated from each other and are mainly due to the movement along the regional dislocations and faults [35][36][37] .

Results and Discussions
India-ITRF2008 angular velocity. Twenty-five continuous and three campaign sites included in our analysis are located on the continental peninsular Indian shield and Indo-gangetic plain and two IGS stations SGOC in Sri Lanka and MALD in the Maldives both located on islands in Indian Ocean (Figs 1, 2 and 3   previous studies demonstrated that little or no intraplate deformation occurs 7,12,13,15 . We thus estimated the India plate angular velocity in ITRF2008 by minimizing the ITRF2008 GPS velocities of these 30 GPS sites using GLORG. Table 2 gives the pole, angular rotation rate, and covariances for our newly estimated angular velocity, as well as angular velocities estimated by previous authors. Figure 3 shows the residual velocities of the 30 India plate sites with respect to the GPS site velocities that are predicted with our new angular velocity in map view (Fig. 3a) and in velocity space (Fig. 3b). Reduced chi-square, the least-squares misfit normalized by the degrees of freedom, is 1.24. The estimated site velocity uncertainties are thus, on average, nearly equal to the velocity misfits. The weighted root-mean-square (WRMS) misfits to the north and east velocity components are 0.79 and 0.90 mm/yr, with ninety percent of the residual velocities lower than 2 mm/yr (Fig. 3b). The uncertainties in the newly estimated angular velocity are small, partly because it is estimated from more than twice as many continuous station velocities as for any previous study and partly due to the superior geographic spread of those stations in the plate interior relative to previous works.
We next describe the motions of all the GPS sites used in the analysis in our newly estimated India plate frame of reference (India-ITRF08) which are given in Table 1. To facilitate the discussion, the deformation is subdivided by region, including continental portions of the plate interior (i.e. peninsular India and the Indo-Gangetic plain) (Figs 3 and 4), the Himalayan Arc (Fig. 5), Northeast India and the Indo-Burmese Arc (Fig. 7) and the Andaman Arc (Fig. 8).
Plate Interior (4-29°N; 68-88°E). Peninsular India and the Indo-Gangetic plains north of the peninsula constitute nearly all continental portions of the plate interior, excepting only the island of Sri Lanka (Fig. 3a). The tectonic features of peninsular India are complex with varied geology and numerous, possibly active, but poorly understood faults [35][36][37][38][39][40] . Since 1900, the level of seismic activity in peninsula (except Kachchh region) can be rated as below moderate to low with very rare occurrence of earthquakes of Magnitude 5.0 and above 29 . Two significant earthquakes have affected this nominally-stable portion of the plate interior in the recent past, the  (Table 1). Earthquakes with magnitude ≥5 from 1970 to 2016 are plotted as solid Red circles (source: http://ds.iris.edu/).  (a) Residual velocities tipped with 70% confidence error ellipse of the 30 GPS sites whose velocities were used to estimate the best-fitting India-ITRF2008 angular velocity in Table 2. ACSZ: Achankovil Shear zone, Blue lines are minor lineaments near PUNE GPS site. Focal mechanisms (Source: http://www.globalcmt. org/CMT) are given for two major earthquakes mentioned in the text, including the 1993 Latur earthquake (Mw 6.3) and the 2001 Bhuj earthquake (Mw 7.6). Earthquakes with magnitude ≥5 from 1970 to 2016 are plotted as solid red circles (source: http://ds.iris.edu/). Figure was created using GMT (generic mapping tool) software 72 . (b) North and east residual site velocity components with respect to velocities estimated with the India-ITRF2008 angular velocity of present study given in Table 2. Red and Black circle represent the velocity of ±2 and ±3 mm/yr respectively. Figure was created using qtiplot 73 .
Scientific RePORTS | 7: 11439 | DOI:10.1038/s41598-017-11697-w Mw = 6.3 Latur earthquake in 1993 38 and the Mw = 7.6 Kachchh/Bhuj/Gujarat earthquake in 2001 33,41 . Although the thrust-faulting focal mechanisms for both earthquakes (Fig. 3a) are consistent with slow ~N-S shortening of the plate interior, the velocities of the mostly continuous GPS sites that span the plate interior show no evidence for an organized deformation pattern at rates faster than ~1 to 2 mm/yr, the approximate resolution of our velocity field. India-fixed composite velocities of peninsular India continuous (except TVM_ and PUNE) and campaign sites are ~0.8 to 1.8 ± 1.5 mm/yr. The Kanpur, Lucknow, Varanasi and Delhi GPS sites on the Indo-Gangetic plain north of peninsular India also indicate composite velocities of ~0.6 to 1.4 mm/yr ± 1.2 mm/ yr. Baseline-length changes of the longest North-South baseline KODI-DELH is 1.1 ± 1.3 mm/yr and East-West baseline KHAV-DURG is 1.7 ± 2 mm/yr represent the deformation rates in the continental portion of plate interior.
Strains are computed at all the 28 continental plate interior sites using strain_zero program of Grid-Strain 42, 43 and plotted in Fig. 4. The maximum principal strain rates for plate interior sites (except TVM_, PUNE and BELP) vary between ~4e-10 to 6e-9 (extension) and the minimum principal strain rates vary from -6e-9 (compression) to 2e-9 (extension) which is consistent with the strain rates reported earlier i.e negligible regional dilatational and shear strain changes in the southernmost 530 km of India 44 ; statistically insignificant strain rate of 2.1 ± 6.1 and 2.6 ± 8.4 · 10-9 yr-1 in the north-south and east-west directions respectively 7 ; southern peninsula move as a rigid block with the velocity of Indian tectonic plate 45 . The evidence for slow or no deformation in these parts of the plate interior support the hypothesis 15 that the plate interior is not significantly affected by the Narmada Son failed rift region or any other fault. By corollary, Indian plate motion is well described by a single angular velocity. This study confirms that Indian plate interior region moves as a rigid block with the velocity of Indian tectonic plate with no significant strain accumulation and the occurrence of Intra-plate earthquakes are due to localized regional deformation specific to the active dislocations and faults in the region.  The two isolated cGPS sites TVM_ and PUNE (Fig. 3) indicate ~3 mm/yr India-fixed velocity which cannot be termed as significant in terms of GPS resolution but nevertheless cannot be ignored considering the GPS data span of 3-5 yrs and needs to be taken in to account for further studies in these regions as detailed below. TVM_ cGPS site located in the southern tip of Indian peninsula has India-fixed velocity of ~3.2 ± 1.2 mm/yr NE. This site is located in Achankovil Shear zone which extends 120 km laterally with a width of 20-30 km (Figs 3 and 4). The shear zone strikes NW-SE with distinct rock types to the north and the south 39 . Baseline length between TVM_ and the KODI, nearest cGPS site indicate shortening of 3.7 ± 1.3 mm/yr and the 3 campaign sites (CHEN, MANP, PLNI) north of KODI also indicate shortening of ~2.6-3.6 mm/yr which may be due to the active regional deformation in this shear zone. Strain rate at Trivandrum indicates compression of 1.2e-8 which is consistent with the shortening in this region.This needs further corroboration by dense network of GPS measurements in Achankovil Shear zone. PUNE GPS site located in the western margin of Deccan Plateau to the north of western ghats of peninsular India indicates India-fixed northward velocity of 3 ± 1 mm/yr in Indian reference frame. Though several faults and lineaments are located in this region (Figs 3 and 4), Panvel flexure is the most conspicuous lineament located in this region which strikes in the north-northwest direction parallel to the west coast and extends approximately between latitudes 16 to 21°N 29 . The origin of Panvel flexure is related to west-coast rifting, subsidence and uplift of Western Ghats. Recent micro-seismic studies indicate existence of active fault systems beneath the Panvel flexure possibly related to west-coast tectonics 40 . GPS motion of ~3 mm/yr in this region and strain rate (compression) of 1.1e-8 may be due to regional deformation related to this active fault system which needs further in-depth study.
The cGPS sites KHAV, BELP, RADP are located in Kachchh rift basin of Gujarat ~50-100 km north of the 2001 Mw = 7.6 Bhuj earthquake epicenter 34 . ISRR GPS site is located in Gandhinagar close to Ahmadabad (~160 km east of epicenter of Bhuj earthquake) and DHAR, MABU are located ~100 and 200 km further north of Ahmadabad. The transient post seismic deformation for Bhuj earthquake is found to be very low and attenuated rapidly within 3-4 years of the earthquake and is much low during 2007-2009 suggesting a weak mantle in this region 46,47 . India fixed velocities and baseline length changes of these continuous GPS sites located in KRB indicate N-S rate of −0.71 ± 2.1 mm/yr between ISRR -MABU and E-W rate of 2.3 ± 2.11 mm/yr between KHAV-DHAR. These measurements are during 2010-2012 and the residual velocities of 1.8 ± 1.5 mm/yr and strain rate (extension) of ~1e-08 at KHAV and BELP GPS sites located close to the epicenter of Bhuj earthquake do indicate the direction pointing towards the epicenter indicating that there may be visco-elastic effects of the 2001 earthquake which the current resolution of our GPS velocity field is not able to quantify.

Himalayan Arc. The Himalayan Arc is seismically active between Kashmir in the west and the Eastern
Syntaxis in the east (Fig. 5) due to active under-thrusting of India tectonic plate below Eurasian plate 48 . The Himalayan arc can be segmented from west to east into Kashmir, Ladakh, Gharwal, Kumaon, Nepal, Sikkim, Bhutan, Arunachal Himalaya and Eastern Syntaxis. It is characterized by several thrust faults (MFT, MBT, MCT, STD, ISZ) that sole in to the basal detachment i.e MHT which constitutes the surface over which Tibet together with Himalayan wedge moves southward over the Indian plate. From north to south each segment of Himalaya is classified as Trans-Himalaya, Tethyan Himalaya, Higher Himalaya, Lesser Himalaya and frontal Himalaya (Fig. 5) which are bounded by these thrust faults.  (Figs 5 and 6). Relative to our newly-estimated India plate frame of reference, the 35 site velocities vary from ~0.2 to 28 mm/yr. Figure 5 shows in map view the velocities of all 35 sites and other sites proximal to the arc and their approximate location from Plate boundary. Figure 6a,b document variations in the velocity components orthogonal and parallel to the Himalayan Arc from west to east. The GPS sites in the foreland immediately south of the arc move southward toward the plate interior at 0.2 to 3 mm/yr, consistent with the elastic effects associated with locking of the Main Frontal Thrust. The arc-normal rates for sites within the Himalayan Arc vary along the arc from west to east and variously move southward at rates of 3-15 mm/yr in the northwestern Himalaya and Nepal Himalaya, ~9 mm/yr in Sikkim; ~6 mm/yr in Bhutan Himalaya; ~16 to 28 mm/yr in Eastern Himalaya (Fig. 6a). The locking depth, width, and slip of Main Himalayan Thrust (MHT) have been estimated previously 11,24,[49][50][51][52] by modeling the GPS surface deformation rates in different segment of Himalayas. Arc-normal convergence rates reported in this study broadly support the hypothesis of the previous study 24 , that the slip and locking depth of MHT varies for different segments of Himalaya.
Average arc-parallel rate in Himalayan arc is ~3 mm/yr for different segments of Himalaya from Kashmir-Himachal to Arunachal. Ladakh Himalaya cGPS sites to the west and Lhasa IGS site in south Tibet to the east indicate average arc-parallel extension rate of ~6 to 9 mm/yr which is broadly consistent with the E-W GPS extension rates reported in Tibet 24, 53 . Spatially non-uniform E-W extension of ~9.7 mm/year that occurs between 78°E and 92°E in ~N110°E direction across the sub-structures in southern Tibet was reported earlier 49 . Their model also suggests that a additional ~3 mm/year deformation over the same distance and direction results from the curvature of locked MHT. Three-dimensional block model of the greater Tibetan Plateau region (20-45°N; 76-108°E) was given 54 by taking into account the inter-seismic GPS deformation, mapped fault geometry, earthquake cycle and moment release. The arc-parallel extension rates in 2500 km Himalayan arc (Fig. 6b) reported here from 35 cGPS sites support this hypothesis 49, 54 . Northeast India and the Indo-Burmese Arc. This seismically active and tectonically complex region is bounded by two convergence plate boundaries i.e Himalayan Arc and Indo-Burmese Arc (Fig. 7). The distinct tectonic domains of northeast India are elevated Shillong Plateau, Brahmaputra (Assam) valley, Indo-Burmese Arc, Eastern Himalaya (Sikkim, Bhutan, Arunachal) and Assam/Eastern Syntaxis. This region has several faults (Fig. 7), with some of the faults being reported to be active. India -fixed velocities of Shillong plateau and Assam valley cGPS sites during the study period indicate southward motion of ~7 mm/yr with respect to the stable Indian shield which is consistent with the rates reported earlier 12,55 . Campaign GPS measurements 55,56 indicate clockwise rotation of Shilling plateau and Assam valley which they attribute to their locations between the Indo-Eurasian and India-Sunda convergence zones. The cGPS site, TZPR in Assam Valley (Fig. 7) moves ~3 mm/ yr southward relative to cGPS sites CSOS/SHIL and GHTU in the Shillong Plateau and Assam Valley. The motion between these sites is attributable to the slip along Kopili fault which fragments the Assam Valley and is consistent with previously reported rates 15,55 .
The Indo-Burmese Arc (IBA) is the western extension of eastern boundary of Indian plate, which is defined by dextral Sagaing fault and Sunda-Arakan trench. This region is segmented by several transverse/oblique faults and has salient-recess topography consisting of foreland Tripura salient in the south and Naga salient in the north which are connected by Manipur recess 12 . India-fixed velocities (Fig. 7) of the three cGPS sites (LUMA, IMPH, AIZW) located in the salient-recess of Indo-Burmese Arc indicate SSW motion of ~16 mm/yr at the northern most site LUMA (Naga Salient); SW motion ~19 mm/yr at IMPH (Manipur Recess) site and ~9.5 mm/yr SW motion at the southernmost site AIWL (Tripura Salient). Indo-Burmese Arc is accommodating about half of the 36 mm/yr relative motion between the India and Sunda plates 10,[57][58][59] . Baseline lengths between these sites indicate extension of ~3 mm/yr between LUMA and IMPH and shortening of ~9 mm/yr between IMPH and Scientific RePORTS | 7: 11439 | DOI:10.1038/s41598-017-11697-w AIZW indicating segmented deformation related to the salient-recess topography and faults in this region 12,55 . An approximately E-W baseline (GHTU-LUMA) between Shillong Plateau and IBA indicates shortening of ~3.6 mm/yr in this region which is due to the subduction of IBA. The cGPS station ANIN northeast of the Lohit thrust fault in the Eastern Himalayan Syntaxis (Figs 6a and 7) indicates ~20 mm/yr oblique convergence between Indian and Eurasian plate and ~8 mm/yr convergence to the south of frontal Himalaya as indicated by the southward motion of CSOS, SHIL, GHTU, TZPR and LUMA in Shillong Plateau, Assam Valley and IBA indicating that faults in this region accommodate 30% of the ~28 mm/yr convergence rate in this region.
Andaman Arc. GPS data of ~four years (2012-2015) from the Port Blair IGS site (PBR2) in the Andaman Arc (Fig. 2, Table 1) to south of the IBA give an India-fixed velocity of 47.9 ± 1.1 mm/yr toward S62°W (i.e ITRF 2008 velocity of 15.5 ± 1.1 mm/yr N and 3.6 ± 1.1 mm/yr E). Figure 8 gives the position time series of PBR2 in north, east and up direction for ~4 years. For comparison time series of the north, east, up component from our earlier analysis were plotted for four epochs of episodic GPS measurements made at Port Blair before earthquake (1996)(1997)(1998)(1999) and at two cGPS sites in Port Blair soon after the earthquake (2005)(2006)(2007)(2008) in Fig. 8. Inter-seismic motion measured from 1996 to 1999 at a nearby episodic site in Port Blair (Fig. 8)

Conclusions
We estimate a new angular velocity for the India plate in ITRF2008 (Table 2) from the velocities of 27 continuous and 3 long-term episodic sites well distributed in the stable plate interior, including the offshore Srilanka and Maldives IGS sites. Our new India-ITRF2008 pole is located ~3 degrees east of the previous estimate 15 , with a slightly (~3%) higher angular velocity considering that both the poles were estimated using ITRF2008 velocities. This difference may be attributed to using 30 GPS sites, more than any previous study and also for the first time GPS sites from the west and east of the plate interior were used in the angular velocity estimate as well as long span of data. Velocity residuals are less than 2 mm/yr for most of the sites located in Indian Plate interior. Strain rates at the continental plate interior sites give maximum extension and compression rates of ~6e-09 indicating that there is no segmentation of the Indian plate along failed rifts or faults and the plate motion is best described with a single angular velocity. Intra-plate earthquakes are due to localized deformation along major faults in the plate interiors which does not contribute to the Indian plate motion.
In the Indian reference frame, TVM_ site in Achankovil Shear zone moves ~3 mm/yr N72°E and PUNE site in Panvel flexure zone moves ~3 mm/yr northward both the sites indicating strain rate (compression) of 1e-08 which may be related to localized active deformation along major faults in these region. This needs further investigation with dense network of GPS measurements to delineate the regional deformation in these regions. KHAV and BELP cGPS near the epicenter of Bhuj earthquake record maximum extension of ~1e-08 compatible with a post-seismic viscoelastic response to the 2001 earthquake which could not be quantified by the current resolution of GPS velocity field.
For the first time arc-normal and arc-parallel rates were estimated using 35 cGPS sites covering the entire 2500 Km Himalayan arc from west to east. Himalayan foreland sites move 0.2-3 mm/yr southward and the rest of the sites move southward at a rate of 3 to 28 mm/yr along the main thrust zones of the different segments of Himalayas i.e NW Himalaya (~15 mm/yr), Sikkim Himalaya (~9 mm/yr), Bhutan Himalaya (~6 mm/yr), Arunachal Himalaya (~20 mm/yr) to maximum rate (~28 mm/yr) in Assam syntaxis. In addition, this total convergence in each segment is accommodated differently in Lesser, Higher, Tethyan and Trans Himalaya. Average arc-parallel motion of ~3 mm/yr recorded in the Himalayan arc is a result of curvature of locked MHT. Arc-parallel extension rate of ~9 mm/yr recorded by Ladakh cGPS and Lhasa site located in South Tibet is due to the extension of Tibet Plateau.
Shillong Plateau and Assam Valley in Northeast India behaves as a distinct block moving southward at a rate of ~7 mm/yr with reference to the plate interior and consistent with the clock wise rotation of this block reported earlier 55 . Shortening and extension rates recorded by cGPS sites in Shillong Plateau, Assam Valley and IBA are consistent with the fragmentation of Assam Valley by Kopili fault and segmented deformation related to the salient-recess topography and faults in Indo Burmese Arc 12,55 . To the south of IBA, Port Blair site in the Andaman subduction zone to the east of Indian plate moves ~48 mm/yr SW with ~16 mm/yr uplift during ~4 year period (2012-2015), thus defining the upper bound of the ongoing active post-seismic relaxation in this region. Indian plate boundary in north and the east is seismically active 62-67 with cluster of seismic events (M < 5) recorded along the MCT, IBA and Andaman Arc which indicate that the stress accumulated due to the present-day active deformation in this region is being intermittently released by these events. Inverse modeling of these GPS derived deformation rates at the plate boundaries when combined with active seismicity and seismic gaps of the region give valuable information on earthquake hazard estimation.

Methodology. GPS Data and Analysis.
Continuous data from 44 Indian GPS sites established as part of Indian GPS network, episodic GPS data from 3 sites in southern India (spanning 17 years), 16 cGPS sites in Nepal and Bhutan Himalaya, and the Indian IGS station data were analysed for the period 1996 to 2015 using GAMIT/GLOBK and the IGS stations shown in Fig. 1. Detailed information about the sites and data used is given in Table 1. Quality control on all RINEX observation files was done with TEQC software 68 . Observation sessions shorter than 12 hours or with numerous cycle slips or poor multipath statistics were eliminated from the analysis. Phase data of all GPS stations were analyzed following the methodology 69 using GAMIT/GLOBK software to obtain loosely constrained daily station positions 69,70 . A data sampling interval of 30 s and elevation cut off angle of 15° was used for all GPS sites. Loosely constrained, daily site positions from GAMIT were refined by minimizing the errors contributed by satellite and receiver clock errors, phase ambiguities, atmosphere, receiver phase center variations and multipath. Zenith tropospheric delays were estimated every two hours using a piecewise linear (PWL) model with 0.5-m priori constraints in zenith delay. Ambiguity-free and ambiguity-fixed solutions were performed with the ionosphere-free linear combination to account for carrier phase ambiguities and signal delay due to the ionosphere. Daily loose GAMIT solutions with their covariance matrix were combined using GLORG for reference frame definition by stabilizing the IGS sites to their ITRF08 coordinates and velocities. Errors in modeling the orbits and atmosphere in the long time series and time-correlated sources of errors in position estimates including monument instability were accounted for by including random-walk noise of 2 mm/ yr in the station coordinate time series during the reference frame stabilization. Table 1 gives the ITRF2008 velocities and their uncertainties from our analysis. The same velocities and their associated error ellipses are shown in Fig. 2.
Strain Estimation at Plate interior sites. Crustal strain rates are estimated from India-fixed GPS velocities of plate interior sites using modified least-squares (MLS) approach 71 and coded in MATLAB as grid-strain program 43 . At any GPS site P surrounded by N GPS sites whose positions and displacements are x n and u n , n = 1, 2, … N, the functional model for Least Square is u = AI + e where matrix A contains the information related to positions of GPS points, I is the strain vector and e is the residual vector. In the MLS approach, computations are performed emphasizing the effect of nearest sites based on adjustment of least-squares covariance matrix by using a scale factor. This gives the choice to use different scales of analysis and evaluate the scale-dependent behaviour of the observed system and select a scale factor that best fits the observed system. We have used the strain_zero program which is a simplified version of grid-strain to calculate strain at GPS points with a scale factor of 150 which best fits the deformation of plate interior sites. This approach gives an option to reduce or exclude a particular GPS point for strain computation which facilitate in estimating the strains of unevenly distributed cluster of GPS points. Maximum and minimum principal strains and the directions thus estimated for the continental plate interior sites are plotted in Fig. 4.