Present-day crustal deformation along the Philippine Fault in Luzon, Philippines

The Philippine Fault results from the oblique convergence between the Philippine Sea Plate and the Sunda Block/Eurasian Plate. The fault exhibits left-lateral slip and transects the Philippine archipelago from the northwest corner of Luzon to the southeast end of Mindanao for about 1200 km. To better understand fault slip behavior along the Philippine Fault, eight GPS surveys were conducted from 1996 to 2008 in the Luzon region. We combine the 12-yr survey-mode GPS data in the Luzon region and continuous GPS data in Taiwan, along with additional 15 International GNSS Service sites in the Asia-Pacific region, and use the GAMIT/GLOBK software to calculate site coordinates. We then estimate the site velocity from position time series by linear regression. Our results show that the horizontal velocities with respect to the Sunda Block gradually decrease from north to south along the western Luzon at rates of 85–49 mm/yr in the west–northwest direction. This feature also implies a southward decrease of convergence rate along the Manila Trench. Significant internal deformation is observed near the Philippine Fault. Using a two dimensional elastic dislocation model and GPS velocities, we invert for fault geometries and back-slip rates of the Philippine Fault. The results indicate that the back-slip rates on the Philippine Fault increase from north to south, with the rates of 22, 37 and 40 mm/yr, respectively, on the northern, central, and southern segments. The inferred long-term fault slip rates of 24–40 mm/yr are very close to back-slip rates on locked fault segments, suggesting the Philippine Fault is fully locked. The stress tensor inversions from earthquake focal mechanisms indicate a transpressional regime in the Luzon area. Directions of σ 1 axes and maximum horizontal compressive axes are between 90° and 110°, consistent with major tectonic features in the Philippines. The high angle between σ 1 axes and the Philippine Fault in central Luzon suggests a weak fault zone possibly associated with fluid pressure. Keywords Philippine Fault GPS velocity Crustal strain Interseismic deformation Dislocation model 1 Introduction The Philippine archipelago is a deformed orogenic belt resulting from the collage and collision of blocks of oceanic and continental affinities ( Karig, 1983 ). It is wedged between two converging plates: the oceanic northwest-moving Philippine Sea Plate in the east and the Sunda (Sundaland) Block/Eurasian Plate in the west ( Fig. 1 ). The east-dipping Manila Trench forms part of its western boundary and together with the Negros-Sulu-Cotabato Trench, absorbs the convergence along the western side. The northwestward motion of Philippine Sea Plate is absorbed in part, by the subduction of west-dipping Philippine Trench and the East Luzon Trough in the east, and in another by the Philippine Fault. Recent data derived from GPS ( Rangin et al., 1999; Simons et al., 1999; Kreemer et al., 2000; Bacolcol et al., 2005 ) and earthquake slip vectors ( Chamot-Rooke and Le Pichon, 1999 ) showed the Sunda Block to be a distinct entity and rotates clockwise with respect to the Eurasian plate ( Chamot-Rooke and Le Pichon, 1999; Michel et al., 2001 ). Simons et al. (2007) used a decade (1994–2004) of GPS data to characterize the Sunda Block boundaries and derived the rotation pole at 49.0°N–94.2°E, with a clockwise rotation rate of 0.34°/Myr. The convergence rate of about 80–90 mm/yr between the Philippine Sea Plate and the Eurasian Plate has been reported from the plate model ( Seno et al., 1993 ) and Global Positioning System (GPS) ( Yu et al., 1999 ). The oblique convergence between two plates is decomposed into a trench-parallel component of 20–25 mm/yr on the Philippine Fault ( Barrier et al., 1991 ) and a trench-perpendicular component of 40–90 mm/yr on the Philippine and Manila Trench ( Megawati et al., 2009 ). The Philippine Fault is a sinistral strike-slip fault which transects the Philippine archipelago from north to south for about 1200 km. In spite of its recognition as a major geological structure and sources of destructive earthquakes (M s 7.5 1973 Ragay Gulf earthquake; M s 7.9 1990 Luzon earthquake; M s 6.2 2002 Masbate earthquake), a number of characteristics, e.g., precise fault location, segmentation, fault slip rates, seismicity, and earthquake recurrence intervals, are poorly understood. The first quantitative measurement along the Philippine Fault was in Mindanao Island wherein a left-lateral displacement of about 28 km was found ( Gervacio, 1971 ). Since then, studies on the motion of the Philippine Fault have been proposed using various approaches ( Acharya, 1980; Karig, 1983; Hirano et al., 1986; Pinet, 1990; Barrier et al., 1991; Aurelio, 1992; Duquesnoy et al., 1994; Galgana et al., 2007 ). Based on GPS measurements, Duquesnoy et al. (1994) infer the slip rate of 26 ± 0.1 mm/yr on the creeping section of the Philippine Fault near the Leyte Island, consistent with slip rates of 23 and 36 mm/yr found in Masbate and Leyte ( Fig. 1 ), respectively ( Bacolcol, 2003; Bacolcol et al., 2005 ). On the other hand, the slip rate of the Philippine Fault is about 17–31 mm/yr in the Luzon area ( Yu et al., 1999 ). Geological and paleo-seismological investigations indicate that the slip rate on the Philippine Fault near central Luzon is generally between 9 and 17 mm/yr ( Daligdig, 1997 ), which is lower than the value computed from GPS data. This discrepancy will be discussed in Section 5 . The area of interest in this study is the segment of Philippine Fault in the Luzon Island, north of the Philippine archipelago. Luzon is part of the N–S trending Luzon arc, a 1200 km chain of mostly late Tertiary to Quaternary volcanics that extends from the Coastal Range of Taiwan (24°N) to Mindoro (13°N). In Luzon, the Philippine Fault acts as the tectonic boundary that separates the Northern Luzon volcanics and Eastern Luzon metamorphics from the Zambales-Angat ophiolites ( Karig, 1983; Karig et al., 1986 ). The Philippine Fault branches into several splays in central Luzon, including the San Jose Fault, the San Manuel Fault, the Gabaldon Fault and the Digdig Fault ( Nakata et al., 1977 ). Historic records indicate large events with M ∼ 7 on the Philippine Fault in the Luzon area occurred in 1901, 1937, and 1973 ( Acharya, 1980 ). The 1973 event occurred along the Guinyangan Fault near 13°N in southern Luzon ( Morante and Allen, 1973; Morante, 1974 ). The region has at least seven major events in the last two centuries with a recurrence interval of about 65 yrs ( Besana and Ando, 2005 ). The most recent large earthquake on the Philippine Fault, the M s 7.9 July 16, 1990 earthquake, occurred on the Digdig Fault segment in central Luzon and ruptured for about 120 km. The average left-lateral slip is 5.4 m from seismic inversion ( Yoshida and Abe, 1992 ) and 5.5–6.5 m from geodetic inversion ( Silcock and Beavan, 2001 ). Both studies infer the bottom depth of coseismic rupture is about 20 km. Additionally, Silcock and Beavan (2001) reported another fault segment of about 45 km-long to the north of the mapped rupture based on geodetic constraints. Daligdig (1997) examined geomorphic and paleoseismic data and found a recurrence interval of about 300–400 yrs along the Digdig Fault. Using satellite imagery, digital elevation models, and geophysical data, Galgana et al. (2007) delineated six tectonic blocks in Luzon and utilized a combination of earthquake slip vectors and GPS-derived horizontal surface velocities to invert for block rotations and elastic strain accumulated on the fault. They found that block rotations can explain the majority of regional deformation in Luzon; while fault-locking strain still makes a significant contribution to the observed GPS velocity field near the fault. The Philippine Fault is locked to partly-locked in the Luzon area and the locking depth is about 25 km. In this paper, we compute velocities of survey-mode GPS sites using data collected between 1996 and 2008 and discuss the implications to regional tectonics. We use an interseismic strain accumulation model ( Savage, 1983; Matsu’ura et al., 1986 ) to invert for the long-term slip rates (i.e., block motion) and back-slip rates (i.e., slip deficits, presumably to be repaid in some forthcoming earthquakes) on various segments along the Philippine Fault. Additionally, we use earthquake focal mechanisms to conduct stress tensor inversions and illustrate the stress status along the Philippine Fault. Our study aims to delineate a comprehensive deformation feature on the Philippine Fault in the Luzon area by combining available geological, geodetic, and seismic data. 2 GPS data acquisition and processing The survey-mode GPS network in Luzon was established in late 1995 and first measured in 1996 by the Institute of Earth Sciences, Academia Sinica, Taiwan in collaboration with the Philippine Institute of Volcanology and Seismology, Philippines. Initially, the network was composed of 15 stations, including 10 newly set up sites, two National Mapping and Resource Information Authority (NAMRIA) sites, and three Geodynamics of South and Southeast Asia (GEODYSSEA, Wilson et al., 1998; Michel et al., 2001 ) sites ( Yu et al., 1999 ). Since then 16 stations were added to densify the network in central and southern Luzon between 1998 and 1999. Succeeding GPS surveys were done annually until 2000. Three more campaigns were conducted in 2004, 2006, and 2008, respectively. In all of these surveys, most of the sites were occupied continuously for 2–3 days using dual-frequency, geodetic GPS receivers, with a 30-s sampling rate. The collected GPS data are processed by GAMIT/GLOBK software packages, version 10.3 ( Herring et al., 2009 ) using standard procedures based on double-difference phase observables, including tropospheric and ionospheric modeling. We fix the International GNSS Service (IGS) final precise ephemerides in the parameter estimation. The GPS data used in the processing includes 30 Luzon sites, 15 permanent IGS sites in the Asia-Pacific region, and 6 sites from Taiwan continuous GPS Array ( Yu et al., 1999 ) in southern Taiwan. These sites were integrated during the processing to obtain a more accurate and consistent regional deformation pattern of the Luzon Island. Fourteen IGS sites (COCO, DAEJ, DARW, GUAM, HOB2, IISC, IRKT, KUNM, PERT, SHAO, TOW2, TSKB, USUD, WUHN) with long observation history surrounding the studied region are adopted as the stabilization sites in GLOBK processing. The positions of these reference stations are constrained to their 2005 International Terrestrial Reference Frame (ITRF2005; Altamimi et al., 2007 ) coordinates. The station positions, variance–covariance matrices, and other parameters from the GAMIT solution are combined in GLOBK to produce estimates of station coordinates for other sites in the ITRF2005 reference frame. The position time series of stations are then extracted from the GLOBK output files. Fig. 2 shows an example of ITRF2005 position time series for east, north and up components of station LUZE. 3 Velocity field We perform a least squares linear fit individually to position time series of three components (north, east, and up) at each GPS station and estimate its average station velocity during the observation period. Outliers and anomalous data are removed prior to the final estimation of velocities. The best-fitting ITRF2005 velocities and their standard errors, station coordinates, and observation time periods for all sites are listed in Table 1 . We then transform the ITRF2005 velocities into a Sunda Block-fixed reference frame ( Simons et al., 2007 ). The horizontal velocities are shown in Table 2 and Fig. 3 . The station velocities with respect to Sunda Block are about 49–89 mm/yr, in the west–northwest (WNW) to northwest (NW) directions. At the Batanes Islands (BTS3) and the northwestern corner of Luzon Island (BRG1) velocities are 83–85 mm/yr in the directions of 286–301°. Along the western coast of Luzon, the station velocities reduce southward gradually, they are 78 mm/yr, in 287° at Santa (LUZD); 70 mm/yr, in 286° at Baguio City; 59 mm/yr, in 282° at Subic (LUZA) and 49 mm/yr, in 283° at Batangas (LUZP). This spatial variation of surface velocities also implies a southward decrease of convergence rate along the Manila Trench. Additionally, the station velocities on the east coast of Luzon are generally more northerly-directed than the corresponding sites on the west coast of Luzon and there are no significant changes on rates southward. The GPS velocities are 89 mm/yr, in the direction of 299° at the northeastern corner of Luzon (LUZI); 89 mm/yr, in 297° at Tuguegarao (LUZH); 85 mm/yr, in 300° at Santiago (LUZG); and 81 mm/yr, in 302° at Daet (CMN2). The relative velocities between four station pairs, LUZI-BRG1, LUZH-LUZD, LUZG-LUZC, and CMN2-LUZP are 19.8 mm/yr, in 10°; 17.6 mm/yr, in 344°; 23.3 mm/yr, in 343°; and 38.1 mm/yr, in 327°, respectively. These relative motions indicate a general picture of strain accumulation across the Philippine Fault. In order to investigate the internal deformation in the Luzon arc, the Luzon ITRF2005 velocities are also transformed to a reference frame fixed at Subic (LUZA), a GPS site located in western Luzon ( Fig. 4 ). For comparison and to be used in the following modeling studies, part of the velocity data published in Galgana et al. (2007) are also included in Fig. 4 . Their results are in general agreement with ours, while we have a denser spatial coverage and longer time history of GPS data across the Philippine Fault. We find that there are no significant velocity jumps across two sites of the Philippine Fault: along San Jose Fault (between PUGA and CRIS), the relative velocity rate is 5.9 ± 2.1 mm/yr; along Digdig Fault (between CRIS and LUZF), the velocity rate is 5.1 ± 1.5 mm/yr ( Fig. 4 ). This feature may imply near-surface aseismic creeping is absent on this portion of the Philippine Fault. Other known active faults found in the area also show insignificant velocity changes across the fault ( Fig. 4 ): – East Zambales Fault (between LUZA and LUZB): 3.1 ± 1.1 mm/yr. – Valley Fault System (between PHIV and LUZL): 0.8 ± 1.0 mm/yr. 4 Dislocation model Using GPS velocities derived in the previous section, we investigate the strain accumulation on the Philippine Fault. We use the interseismic crustal deformation model ( Savage, 1983; Matsu’ura et al., 1986 ) such that the interseismic velocity field can be represented as the sum of rigid block motion at long-term rates across the fault and back-slip rate (slip deficits) on the locked fault (negative dislocation). Because of an irregular spatial distribution of GPS sites in the Luzon area, we only choose three transects across the Philippine Fault wherein GPS data is sufficient to constrain fault slip rates. These transects from north to south are shown in Fig. 5 , Transect AA ′ is E–W directed and covered the northern Luzon region between 16°N and 17°N, Transect BB ′ is NE–SW directed with the width of 100 km and located in the middle portion of the Luzon Island wherein the 1990 M s 7.9 earthquake occurred, and Transect CC ′ is NE–SW directed with the width of 130 km and covered the southern portion of the Luzon Island. This region has experienced many large earthquakes in the past ( Besana and Ando, 2005 ). We project GPS data along these 2-D transects and estimate back-slip rate, fault width, dip, long-term slip rates using an elastic half-space model ( Okada, 1985 ). The positions of our model faults approximately follow the surface trace of the Philippine Fault. An initial value of fault dip is assumed to be 90°. We search for the optimal long-term slip rates and fault width using a grid search approach. In these 2-D models, we only consider the motion parallel to the Philippine Fault. The data are inverted using a weighted least-squares approach by minimizing the following functional: (1) F ( s , m ) = ∑ - 1 / 2 ( G ( m ) s - d ) 2 where ∑ - 1 / 2 is the inverse square root of the data covariance matrix; G ( m ) are Green’s functions, which depend on the fault geometry parameters m ; s is slip rate; d is the observed GPS velocities. The fit to the data is quantified from the mean of the normalized square residuals, χ r 2 , defined as the chi-square divide by the number of data points. A value of 1 for χ r 2 means that the model fits the data within uncertainties on average. Due to limitations of our knowledge of fault parameters, we invert for fault geometry, long-term fault slip rate, and back-slip slip rate using a gird search approach. We search for values of long-term rates from 15 to 45 mm/yr, fault widths from 10 km to 40 km, and fault dip from 60° to 90° based on apriori constraints from previous studies ( Acharya, 1980; Barrier et al., 1991 ). The inversions show that the modeling results are sensitive to the long-term slip rates but less sensitive to the fault width. However, our method does not directly provide correlation coefficients between all fault parameters. In order to estimate uncertainties of fault parameters, we apply a bootstrap method by re-sampling GPS data to generate 1000 synthetic data sets and compute fault parameters. The starting fault parameters are obtained from the results of grid search. Estimates of optimal fault parameters and their 95% confidence regions are given in Table 3 . 5 Results and discussion 5.1 Fault parameters The optimal models for three transects generally fit the data with the values of χ r 2 between 1 and 2.6. In Transect AA ′, the back-slip rate and fault width are 22 mm/yr and 15 km, respectively ( Table 3 and Fig. 6 a). However, a poor spatial coverage of GPS sites in the far field may limit the resolution on the fault parameters. In Transect BB ′, the back-slip rate and fault width are 37 mm/yr and 15 km, respectively ( Fig. 6 b). The spatial distribution of GPS sites in Transect BB ′ is the best among three profiles; therefore, the fault parameters are well determined. In Transect CC’ , the back-slip rate and fault width are 40 mm/yr and 28 km, respectively ( Fig. 6 c). The preferred fault dip is 70°, different from a nearly 90° dip obtained in Transects AA ′ and BB ′. The deformation near Transect CC ′ is more complex than other two sections and is possibly affected by the collision tectonics between Palawan and Mindoro ( Rangin et al., 1999 ). Fig. 5 shows the seismicity between 1977 and 2009 from US Geological Survey (USGS)/National Earthquake Information Center (NEIC) catalog. The seismicity near the Philippine Fault system in central Luzon is mainly related to the 1990 M s 7.9 Luzon earthquake and its aftershocks, while the seismicity activity is absent on the other segments of the Philippine Fault. The coupling ratio is defined as the back-slip rate relative to the long-term slip rate. If the ratio is large (close to 1.0), meaning fault is fully locked, it implies the earthquake potential is high in the area. In these three transects, the back-slip rates are close to the long-term fault slip rates ( Table 3 ), suggesting that the Philippine Fault is nearly fully locked in the Luzon area. The inferred fault slip rates increase from north to south, in agreement with the distribution of historic large earthquakes and the results of fault block models reported in Galgana et al. (2007) , which shows the slip rates near Transects AA ′- BB ′- CC ′ are 17–27, 25–29, and 29–40 mm/yr, respectively. Our estimate of the slip rate of about 37 mm/yr near the 1990 rupture zone is similar to 42 mm/yr derived from an elastic model ( Beavan et al., 2001 ). Additionally, the inferred fault slip rates are 20–25 mm/yr based on plate kinematics ( Barrier et al., 1991 ) and 35 mm/yr from GPS measurements ( Rangin et al., 1999; Thibault, 1999; Yu et al., 1999 ). However, fault slip rates inferred from geodetic data are larger than geologic slip rates. Geological and paleo-seismological investigations suggest that the slip rate on the Philippine Fault (near Transect BB ′) is generally between 9–17 mm/yr and the recurrence interval is 300–400 yrs ( Daligdig, 1997 ). The discrepancy between the long-term and short-terms slip rates presumably results from the contamination of postseismic deformation associated with the 1990 earthquake. The 1990 M s 7.9 earthquake produced a 110 km-long surface rupture and 5–6 m left-lateral slip near the Transect BB’ ( Yoshida and Abe, 1992; Silcock and Beavan, 2001 ). The mainshock is likely to induce viscoelastic relaxation of the lower lithosphere. Our preferred slip rate of 24–40 mm/yr on the Philippine Fault can be represented as an upper bound. Additionally, we compute the recurrence interval of about 162 yr near the rupture area of 1990 earthquake (Transect BB ′) based on the average coseismic slip of 6 m ( Yoshida and Abe, 1992; Silcock and Beavan, 2001 ) and the interseismic slip rate of 37 mm/yr. The return period of the 1990 type event derived from the ratio of peak slip to fault slip rate in this study can be represented as the lower bound of M ∼ 7 earthquakes in the seismic hazard analysis. 5.2 Stress tensor analysis To give a comprehensive understanding of deformation in the seismogenic zone, we investigate earthquake focal mechanisms between 1977 and 2009 from USGS/NEIC catalog to reveal stress status at depths. We choose events near the Luzon area (12–19°N, 118–125°E) with focal depths less than 40 km ( Fig. 7 a). The earthquake magnitude falls in the range of 4.7–7.7. Most events occur on the western branch of the Philippine Fault (16.5°N, 121°E) in northern Luzon and on the Verde Passage-Sibuyan Sea Fault (near 13–14°N, 120–122°E, VPSF in Fig. 4 ). We use the stress tensor inversion method proposed by Michael (1984, 1987) which assumes that the slip on the fault plane occurs in the direction parallel to the direction of resolved shear stress. Figs. 7b and 8 present inversion results based on a moving–window approach on the 0.5°-spacing grid. We include all events within a 1° × 1° rectangle centered at the node and estimate the stress tensor when there are at least 10 earthquakes within a given rectangular box ( Fig. 7 b). The present-day deformation shows predominately strike-slip faulting in Luzon with the stress ratio less than 0.5 ( Fig. 8 a), indicating a transpressional regime and corresponding with geological structures. The σ 1 axes are trending in the direction between WNW–ESE and E–W, while σ 2 axes are close to vertical oriented. Additionally, we compute the maximum horizontal compressive direction of about 90–110° in the Luzon area ( Fig. 8 b), consistent with major tectonic features in the Philippines ( Barrier et al., 1991; Aurelio et al., 1997; Rangin et al., 1999 ). Directions of σ 1 axes are sub-perpendicular to the Philippine Fault in central Luzon ( Fig. 7 b). A similar distribution of σ 1 axes has been found along the San Andreas Fault in California ( Zoback et al., 1987; Mount and Suppe, 1992 ). The σ 1 axes should be oriented about 30° to the fault plane for a strong fault ( Byerlee, 1978 ), thereby the high angle between σ 1 axes and fault plane in central Luzon suggests a weak fault zone. The mechanism for the weak fault might be related to high fluid pressure lowering the effective normal stress on a fault and decreasing shear stress ( Hubbert and Rubey, 1959; Sleep and Blanpied, 1992; Hardebeck and Hauksson, 1999 ). The fluid in the western central Luzon possibly release from the slab below the forearc region and facilitate partial melting and active volcanism ( Defant et al., 1988 ). However, the total number of earthquake focal mechanism is insufficient to resolve the spatial variation of stress axes near the fault zone. The analysis on small to moderate-sized earthquake focal mechanisms from the Luzon seismic network is required in the future. Based on modeling results in this study, we suggest to densify the present GPS network in the Luzon area. Additional sites are required to confidently resolve the long-term slip rates of the faults comprising the Philippine Fault system. Constructing new sites along known active structures like the Valley Fault System and East Zambales Fault ( Fig. 4 ) plays a crucial role in seismic hazard analysis. Measurements from various data sets would be gathered and eventually used towards a development of a comprehensive program in evaluating seismic risks in the Philippines. 6 Conclusions Using the survey-mode GPS data collected on eight campaigns from 1996 to 2008 and continuous GPS data from Taiwan Continuous GPS Array and IGS sites in the Asia-Pacific region, we derive a new velocity field in the Taiwan–Luzon region. The GPS velocities with respect to Sunda Block move at rates of 49–89 mm/yr in the west–northwest to northwest directions in Luzon. We observe a gradual decreasing of station velocity from 85 mm/yr at the northwest corner of Luzon to 49 mm/yr at Batangas of southwestern Luzon. This feature also suggests a southward decrease of converging rate along the Manila Trench. Using a two dimensional elastic dislocation model, we infer the fault slip rate falls in the range of 24–40 mm/yr and increases from northern to southern Luzon. The Philippine Fault is nearly fully locked in the Luzon area. We also perform stress tensor inversions using earthquake focal mechanisms and find a transpressional stress regime in Luzon. The σ 1 axes are in the direction of 90–110°, at high angle to a nearly N–S trending fault strike. This implies the segment of the Philippine Fault in central Luzon is possibly weak. Acknowledgements We thank many colleagues at the Philippine Institute of Volcanology and Seismology and the Institute of Earth Sciences, Academia Sinica who participated in Luzon GPS surveys. We also thank Ministry of the Interior (MOI) of Taiwan and International GNSS Service (IGS) community for providing the continuous GPS data in this study. Many figures in this paper were generated using the generic mapping tools (GMT) ( Wessel and Smith, 1998 ). We are indebted to two anonymous reviewers for their constructive suggestions and comments. This study was financially supported by Academia Sinica and National Science Council of Taiwan under Grant NSC 98-2119-M-001-030. This is a contribution of the Institute of Earth Sciences, Academia Sinica, IESAS1507. References Acharya, 1980 H.K. Acharya Seismic slip on the Philippine fault and its tectonic implications Geology 8 1980 40 42 Altamimi et al., 2007 Z. Altamimi X. Collilieux J. Legrand B. Garayt C. Boucher ITRF2005: a new release of the International Terrestrial Reference Frame based on time series of station positions and earth orientation parameters Journal of Geophysical Research 112 2007 19 21 Aurelio, 1992 Aurelio, M., 1992. Tectonique du segment central de la faille philippine, étude structurale, cinématique et évolution géodynamique, Ph.D., Univ. Paris 6, Paris. Aurelio et al., 1997 M. Aurelio E. Barrier R. Gaulon C. Rangin Deformation and stress states along the central segment of the Philippine Fault: implications to wrench fault tectonics Journal of Asian Earth Sciences 15 1997 107 119 Bacolcol, 2003 Bacolcol, T., 2003. Etudes geodesiques de la faille philippine a partir des donnees GPS, Thèse de Doctorat d’état, Université Pierre et Marie Curie. Bacolcol et al., 2005 T. Bacolcol A. Aguilar G. Jorgio R. de la Cruz M. Lasala GPS constraints on the Philippine fault slip rate in Masbate Island, Central Philippines Journal of the Geological Society of Philippines 60 2005 Barrier et al., 1991 E. Barrier P. Huchon M. Aurelio Philippine fault – a key for Philippine kinematics Geology 19 1991 32 35 Beavan et al., 2001 J. Beavan D. Silcock M. Hamburger E. Ramos C. Thibault R. Feir Geodetic constraints on postseismic deformation following the 1990 M-s 7.8 Luzon earthquake and implications for Luzon tectonics and Philippine Sea plate motion Geochemistry Geophysics Geosystems 2 2001 10.1029/2000GC000100 Besana and Ando, 2005 G.M. Besana M. Ando The central Philippine Fault Zone: location of great earthquakes, slow events, and creep activity Earth Planets and Space 57 2005 987 994 Byerlee, 1978 J. Byerlee Friction of rocks Pure and Applied Geophysics 116 1978 615 626 Chamot-Rooke and Le Pichon, 1999 N. Chamot-Rooke X. Le Pichon GPS determined eastward Sundaland motion with respect to Eurasia confirmed by earthquakes slip vectors at Sunda and Philippine trenches Earth and Planetary Science Letters 173 1999 439 455 Daligdig, 1997 Daligdig, J., 1997. Recent Faulting and Paleoseismicity Along the Philippine Fault Zone, North Central Luzon, Philippines, Ph.D. Thesis, Kyoto University. Defant et al., 1988 M.J. Defant J.Z. Deboer D. Oles The western central Luzon volcanic arc, the Philippines – 2 arcs divided by rifting Tectonophysics 145 1988 305 317 Duquesnoy et al., 1994 T. Duquesnoy E. Barrier M. Kasser M. Aurelio R. Gaulon R.S. Punongbayan C. Rangin B.C. Bautista E. Delacruz M. Isada S. Marc J. Puertollano A. Ramos M. Prevot A. Dupio I. Eto F.G. Sajona D. Rigor F.G. Delfin D. Layugan Detection of creep along the Philippine fault – 1st results of geodetic measurements on Leyte-island, central Philippine Geophysical Research Letters 21 1994 975 978 Galgana et al., 2007 G. Galgana M. Hamburger R. McCaffrey E. Corpuz Q.Z. Chen Analysis of crustal deformation in Luzon, Philippines using geodetic observations and earthquake focal mechanisms Tectonophysics 432 2007 63 87 Gervacio, 1971 F.C. Gervacio Geotectonic development of the Philippines Journal of the Geological Society of the Philippines 25 1971 18 38 Hardebeck and Hauksson, 1999 J.L. Hardebeck E. Hauksson Role of fluids in faulting inferred from stress field signatures Science 285 1999 236 239 Herring et al., 2009 Herring, T.A., King, R.W., McClusky, S.C., 2009. Documentation for the GAMIT Analysis Software , Release 10.3, Massachusetts Institute of Technology, Cambridge, MA, 183pp. Hirano et al., 1986 S. Hirano T. Nakata A. Sangawa Fault topography and quaternary faulting along the Philippine Fault zone, Central Luzon, the Philippines Journal of Geography 95 1986 1 23 Hubbert and Rubey, 1959 M.K. Hubbert W.W. Rubey Role of fluid pressure in mechanics of overthrust faulting. I. Mechanics of fluid-filled porous solids and its application to overthrust faulting Geological Society of America Bulletin 79 1959 115 166 Karig, 1983 D.E. Karig Accreted terranes in the northern part of the Philippine archipelago Tectonics 2 1983 211 236 Karig et al., 1986 D.E. Karig D.R. Sarewitz G.D. Haeck Role of strike-slip faulting in the evolution of allochthonous terranes in the Philippines Geology 14 1986 852 855 Kreemer et al., 2000 C. Kreemer W.E. Holt S. Goes R. Govers Active deformation in eastern Indonesia and the Philippines from GPS and seismicity data Journal of Geophysical Research 105 2000 663 680 Matsu’ura et al., 1986 M. Matsu’ura D.D. Jackson A. Cheng Dislocation model for aseismic crustal deformation at Hollister, California Journal of Geophysical Research 91 1986 2661 2674 Megawati et al., 2009 K. Megawati F. Shaw K. Sieh Z.H. Huang T.R. Wu Y.N. Lin S.K. Tan T.C. Pan Tsunami hazard from the subduction megathrust of the South China Sea: Part I. Source characterization and the resulting tsunami Journal of Asian Earth Sciences 36 2009 13 20 Michael, 1984 A.J. Michael Determination of stress from slip data - faults and folds Journal of Geophysical Research 89 1984 1517 1526 Michael, 1987 A.J. Michael Use of focal mechanisms to determine Stress - a control study Journal of Geophysical Research 92 1987 357 368 Michel et al., 2001 G.M. Michel Y.Q. Yu S.Y. Zhu C. Reigber M. Becker E. Reinhart W. Simons B. Ambrosius C. Vigny N. Chamot-Rooke X. Le Pichon P. Morgan S. Matheussen Crustal motion and block behaviour in SE-Asia from GPS measurements Earth and Planetary Science Letters 187 2001 239 244 Morante, 1974 E.M. Morante The Ragay Gulf earthquake of March 17, 1973, southern Luzon, Philippines Journal of the Geological Society of the Philippines 28 1974 69 93 Morante and Allen, 1973 E.M. Morante C.R. Allen Displacement of the Philippine fault during the Ragay Gulf earthquake of March 17, 1973 Geological Society of America Bulletin 5 1973 744 745 Mount and Suppe, 1992 V.S. Mount J. Suppe Present-day stress orientations adjacent to active strike-slip faults: California and Sumatra Journal of Geophysical Research 97 1992 11,995 91,2013 Nakata et al., 1977 T. Nakata A. Sangawa S. Hirano A report on tectonic landforms along the Philippine fault zone in Northern Luzon Philippines Science Reports 7th Series-Geography 27 1977 69 93 Okada, 1985 Y. Okada Surface deformation to shear and tensile faults in a half space Bulletin. Seismological Society of America 75 1985 1135 1154 Pinet, 1990 Pinet, N., 1990. Répartition des mouvements decrochants et chevauchements le long de frontière NW de la plaque Mer des Philippines, Paris. Rangin et al., 1999 C. Rangin X. Le Pichon S. Mazzotti M. Pubellier N. Chamotrooke M. Aurelio A. Walpersdorf R. Quebral Plate convergence measured by GPS across the Sundaland/Philippine sea plate deformed boundary: the Philippines and eastern Indonesia Geophysical Journal International 139 1999 296 316 Savage, 1983 J.C. Savage A dislocation model of strain accumulation and release at a subduction zone Journal of Geophysical Research 88 1983 4984 4996 Seno et al., 1993 T. Seno S. Stein A.E. Gripp A model for the motion of the Philippine Sea Plate consistent with NUVEL-1 and geological data Journal of Geophysical Research 98 1993 17941 17948 Silcock and Beavan, 2001 Silcock, D.M., Beavan, J., 2001. Geodetic constrains on coseismic rupture during the 1990 M-s 7.8 Luzon, Philippines, earthquake, Geochemistry Geophysics Geosystems 2, 2000GC000101. Simons et al., 1999 Simons, W.J.F., Ambrosius, B.A.C., Noomen, R., Angermann, D., Wilson, P., Becker, M., Reinhart, E., Walpersdorf, A., Vigny, C., 1999. Observing Plate Motions in S.E. Asia: Geodetic Results of the GEODYSSEA Project , edn, American Geophysical Union, Washington, DC, ETATS-UNIS. Simons et al., 2007 W.J.F. Simons A. Socquet C. Vigny B.A.C. Ambrosius S.H. Abu C. Promthong C. Subarya D.A. Sarsito S. Matheussen P. Morgan W. Spakman A decade of GPS in Southeast Asia: resolving Sundaland motion and boundaries Journal of Geophysical Research 112 2007 Sleep and Blanpied, 1992 N.H. Sleep M.L. Blanpied Creep, compaction and the weak rheology of major faults Nature 359 1992 687 692 Thibault, 1999 Thibault, C.A., 1999. GPS Measurements of Crustal Deformation in the Northern Philippine Island Arc, MS Thesis, Indiana University. Wessel and Smith, 1998 P. Wessel W.H.F. Smith New, improved version of generic mapping tools released EOS Transactions AGU 79 47 1998 579 Wilson et al., 1998 P. Wilson J. Rais C. Reigber E. Reinhart B.A.C. Ambrosius X. Le Pichon M. Kasser P. Suharto D. Majid D.A. Yaakub R. Almeda C. Boonphakdee Study provides data on active plate tectonics in Southeast Asia region Eos, Transactions AGU 79 545 1998 548 549 Yoshida and Abe, 1992 Y. Yoshida K. Abe Source mechanism of the Luzon, Philippines earthquake of July 16, 1990 Geophysical Research Letters 19 1992 545 548 Yu et al., 1999 S.B. Yu L.C. Kuo R.S. Punongbayan E.G. Ramos GPS observation of crustal deformation in the Taiwan–Luzon region Geophysical Research Letters 26 1999 923 926 Zoback et al., 1987 M.D. Zoback M.L. Zoback V.S. Mount J. Suppe J.P. Eaton J.H. Healy D. Oppenheimer P. Reasenberg L. Jones C.B. Raleigh I.G. Wong O. Scotti C. Wentworth New evidence on the state of stress of the San-Andreas fault system Science 238 1987 1105 1111