•     Land subsidence
•     Metro manila
•     InSAR time-series
•     Groundwater extraction

Extensive land subsidence in and around Metro Manila, largely from overuse of groundwater, is orders of magnitude more rapid than sea-level rise from global warming. It enhances the exposure of its residents to worsening floods and tidal incursions, and greatly exacerbates the storm-surge threat. Additionally, differential movements at pre-existing faults causes considerable damage to overlying properties. The subsidence has been recognized and quantified since the 1990s from the rates at which wells apparently rise as the ground around them sinks, the frequency with which roads have had to be raised, and social surveys of historical flooding and tide heights. To assess the impact of the subsidence and devise appropriate countermeasures requires that the subsiding areas be delineated and their subsidence rates be measured more precisely. Toward those ends, 2003-2010 Envisat and 2007-2011 ALOS PALSAR-1 imagery were processed using the PSInSAR method. During these periods, Manila, CAMANAVA, Rosario, San Pedro, Las Piñas and Dasmariñas subsided by at least 2-4.2 cm/yr, due mainly to over-extraction of groundwater. InSAR data also reveal subsidence in highly urbanized areas that also rely heavily on groundwater along southern segments of the West Marikina Valley Fault, but none on the northern segment. Land subsidence along the coast.

Eco RC, Rodolfo KS, Sulapas JJ, Morales Rivera AM, Lagmay AMF, et al. (2020) Disaster in Slow Motion: Widespread Land Subsidence in and Around Metro Manila, Philippines Quantified By Insar Time-Series Analysis. JSM Environ Sci Ecol 8(1): 1068.

Widespread land subsidence has been a persistent problem for many years in coastal Metro Manila and the adjacent coastal areas of Bulacan and Pampanga. This phenomenon has been recognized and quantified since 1997 from the rates at which wells apparently rise as the land surrounding them sinks, the frequency with which roads needed to be raised to avoid renewed inundation, social surveys of local communities, and emergence of ground fissures [1-4]. From 1991 to 2003 numerous areas around northern Manila Bay subsided at rates of 1.7-8.3 cm/yr [5]. More recently, Raucoules [6], used differential interferometric synthetic aperture radar (DInSAR), to measure subsidence rates of as much as 15 cm/yr from 1993-2010.

Areas away from the coast also suffer from land subsidence. Numerous ground fissures have emerged since the early 90s in Muntinlupa and San Pedro. Southern left-stepping en echelon segments of the West Marikina Valley Fault (WMVF), system are dominated by vertical movements; however, slip rates of about 2-20 cm/yr between 1996-1999 are higher than those that normally occur worldwide, suggesting that groundwater extraction is the major contributing factor [4,7]. Buildings atop these fissures have been extensively damaged and structurally compromised.

Anthropogenic land subsidence is also a concern in many other cities worldwide. In Shanghai, the added weight of new buildings has contributed 30 to 40 percent of more than a meter of subsidence from 2000 to 2010 [8]. Six major cities in Indonesia are subsiding at maximal rates of 22 cm/yr, mainly due to the extraction of groundwater and gas [9]. Bangkok, Thailand, which sank by as much as 12 cm/yr from 1978 to 1981, began to control usage by increasing the price of groundwater and piping in surface water. The city now subsides only up to 3 cm/yr, but flooding continues to worsen nonetheless [10,11]. In the Las Vegas Valley, combined InSAR, GPS analyses, and conventional leveling surveys have documented extensive land subsidence rates of at least 5-6 cm/yr since 1963. These rates have since declined from 2.5 cm/yr since 1991 to only a few millimeters per year by the 2000s [12,13]. Much of the deformation there is controlled by faulting [14]. Flooding from Hurricane Katrina was exacerbated when levees in New Orleans were overtopped because they had subsided too low to prevent floodwaters from rushing inland [15].

Worsening floods and tidal incursions in coastal Metro Manila and adjacent Bulacan are usually blamed on upland deforestation, urbanization, channel encroachment by informal settlers and fishponds, and garbage dumping in estuaries. Although not widely recognized, anthropogenic land subsidence from excessive groundwater extraction for coastal aquaculture and mining are significant contributors [5,9,16]. Around northern Manila Bay, for example, excessive groundwater extraction from deltaic muds and intercalated sand and gravel has caused subsidence of centimeters to more than a decimeter per year [1,2].

In this study, we use interferometric synthetic aperture radar (InSAR) time-series spanning 2003-2011 period to identify and characterize the magnitude, spatial extent, and temporal evolution of land subsidence in Metro Manila and surrounding regions. This highlights the advantages of using satellite-based remote sensing in detecting and monitoring land subsidence, especially in urban areas.

We used 35 Envisat ASAR scenes acquired from 2003 to 2010 and 19 ALOS PALSAR-1 scenes acquired from 2007 to 2011, respectively of descending and ascending tracks (Figure 1). The images were first focused to generate single-look complex images using ROI_PAC [17]. Pairs were then selected based on minimum perpendicular baselines relative to “master” images: 25 August 2004 for Envisat, and 28 August 2009 for ALOS PALSAR-1 (Figure 2). Using Doris software [18], the master images were then paired with the rest of images to generate interferograms. The interferograms were unwrapped using a statistical-cost approach with SNAPHU [19,20].

The displacement time-series and velocities were produced using the Stanford Method for Persistent Scatterers (StaMPS), which isolates places with radar-bright features such as houses and buildings to record slow-moving displacements [21]. This method minimizes temporal decorrelation and phase contributions caused by atmospheric variations, which are particularly problematic in tropical environments. Linear tropospheric correction was applied based on the phase and topographic data [22].

The PSInSAR processing derives displacements along line-of-sight (LOS), towards the satellite. In areas with both ascending and descending data, the LOS displacement to its vector components can be decomposed. However, the near-polar track of the Envisat and ALOS orbits render the displacement observations insensitive to north-south movements [23]. Inasmuch as the main objective of this study is to identify and evaluate subsiding areas, we disregard this vector component. We assume that subsidence occurs at a constant rate during the study period and compute the subsidence rate (vertical velocity component vv ) as [24]

where vA and vD are the LOS velocities from the ascending and descending data, respectively, and θA and θD are incidence angles. In our case, the average θ values are 20.1º and 39.2º for Envisat and ALOS-1, respectively. Mentions of displacements in this report will be in terms of LOS unless otherwise specified.

The maps in Figure 3 show the spatio-temporal deformation patterns yielded by the descending Envisat (2003-2010), and ascending ALOS PALSAR-1 (2007-2011), data. Positive LOS velocities (blue color) depict movement towards the satellite, such as uplift, whereas negative values (red color) represent movement away from the satellite such as subsidence. We obtain good coherence across most of the study area owing to the widespread and contiguous distribution of built-up areas in and around Metro Manila, especially in the north. To the south, urban centers clusters within a large agricultural and vegetated region with little coherence.

Our results show clear deformation (Figure 4). Most significant are LOS range increases in Caloocan, Malabon, Navotas, Valenzuela-collectively referred to as “CAMANAVA”- as well as Manila, Las Piñas, Rosario, San Pedro, Biñan, and Dasmariñas cities. All these locations are highly urbanized and densely populated. The mean velocities and the associated standard deviations at the reference point provides an estimate of the background noise: -0.008 ±0.01 cm/yr and -0.01 ± 0.005 cm/yr from the ascending and descending data, respectively. Equation 1 yields the vertical rates from the combined ascending and descending data, dominated by subsidence signals from the six areas (Table 1). Observations in these areas are summarized in the subsections below.

Manila and CAMANAVA

The subsiding area in Manila and CAMANAVA is contiguous and round-shaped that overlaps the administrative boundaries of the said cities. Both Envisat and ALOS 1 datasets showed a linear temporal evolution of the subsidence; maximum subsidence rates computed for Manila and CAMANAVA were 2.5 and 4.2 cm/ yr, respectively. The total affected area was about 67 km2 .

Las Piñas

The deformation zone in Las Piñas was irregularly shaped and affected about 21 km2 of the city. Subsidence rates from 2003 to 2011 were about 2 cm/yr; rates were highest from 2003 to 2006 then slightly decreased from 2007 to 2011, as reflected in both Envisat and ALOS 1 time-series graphs.

San Pedro and Biñan

The southern segment of the West Marikina Valley Fault constrained most of the northern and western extent of subsidence in San Pedro and Biñan. The southern extent exhibited a linear, NE-SW trending deformation pattern; however, there is no known fault that constrains that section of the deformation area. The area subsided about 3 cm/yr during the study period and affected about 27 km2 of the region. Here, the Envisat and ALOS PALSAR-1 showed slightly different trends. Envisat data documented a steady increase in LOS range until around 2007, whereby subsidence slowed drastically. Such slowing down was not detected in the ALOS PALSAR-1 time-series.

Dasmariñas

Subsidence in Dasmariñas was mainly concentrated in the urban center of the city, extending to about 40 km2 in a bowlshaped pattern. Both the Envisat and ALOS PALSAR-1 data documented a steady increase in LOS range during the entire observation period. Combined ascending and descending data yielded cumulative subsidence rates of 3 cm/yr.

Rosario

In Rosario, both the Envisat and ALOS PALSAR-1 time-series show that LOS range increased rapidly until around 2006, then began to decrease. The maximum cumulative subsidence rate was 1.1 cm/yr affecting about 13 km2 . Similar to Dasmariñas, the deformation pattern was bowl-shaped. By the later part of the time-series, it appears that subsidence may have been significantly reduced or even stopped altogether.

Table 1: Subsidence rates (in LOS and vertical displacement rates) in key areas in and around Metro Manila.

Groundwater extraction and land subsidence

During deposition, fluvial sediments undergo natural “autocompaction”, compressed because the weight of each successive sediment deposit squeezes water out of the layers beneath them [2]. This natural process causes the surface to subside no more than a few millimeters per year, such as in the Pampanga River Delta northwest of Manila, where natural autocompaction rates are about 2.5-5 mm/yr, roughly the same magnitude as global sea-level rise [25]. Without human impact, the surface is built back up again by the accumulation of new river deposits, similarly a few millimeters thick.

Pumping water out of the ground faster than it is replenished lowers the water table of an unconfined aquifer, and lowers the piezometric surface around the well in a confined aquifer, creating a radial-shaped cone of depression; the larger the extraction, the larger the cone [26,27]. The causal relationship between land subsidence and over-extraction of groundwater has long been recognized [28-30], and clearly described by Galloway [31].

The National Water Resources Board (NWRB), identified several cones of depression in and around Metro Manila [32]. Five of these are within our study area: one each in CAMANAVA, Las Piñas, and Imus, and two in central Metro Manila (Figure 4). Among these areas, only central Metro Manila did not subside.

Figure 4 shows the areal extent of subsiding areas with the NWRB groundwater surface map and cones of depression. In CAMANAVA, the subsiding areas do not exactly coincide with the center of the cone of depression but are nevertheless still within its overall extent. The Las Piñas cone of depression coincides with the location of land subsidence. In Imus, the cone of depression is displaced a few kilometers away from the subsiding area in Dasmariñas. Additionally, concentric piezometric contour lines starting from -20 m above sea level inland going down to -100 masl surrounding Rosario has not been considered as a cone of depression, but do document land subsidence.

Surface geology and land subsidence

Three geological units have been mapped within the study area: Guadalupe Formation, Taal Tuff, and recent alluvial and fluvial deposits (Figure 4). The Guadalupe Formation underlies most of inland Metro Manila. It has two members: Diliman Tuff and Alat Conglomerate. The presence of Stegodon fossils and other vertebrate remains indicate a Pliestocene age for this formation [33].

Diliman Tuff is of uncertain origin, and although some workers have attributed this unit to eruptions from the Laguna de Bay Caldera [34], the upper sections of this deposit is chemically distinct from other known deposits from Laguna de Bay [35]. This geologic unit is thin-to medium-bedded with regularly stratified layers that are almost flat-lying, although in some parts they gently dip to the west. It is composed mainly of fine-grained vitric tuffs and welded volcanic breccias with subordinate amount of tuffaceous fine- to medium-grained sandstones [36,37].

Alat Conglomerate is massive, poorly-sorted with wellrounded pebbles and small boulders cemented by coarsegrained, calcareous and sandy matrix. There are also interbedded sandstones that are loosely cemented and friable, and consists of tuffaceous, fine- to medium-grain sized particles. Also present in this unit are mudstones that are medium- to thin-bedded, soft, sticky, silty, and tuffaceous [33].

Taal Tuff is composed of fine- to medium-grained basaltic tuffs and minor flows, including volcanic breccias; as the name implies, these deposits are ascribed to calderagenic eruptions from the Taal Caldera. Recent alluvial and fluvial deposits along the western and eastern coastlines are mainly unconsolidated deposits of silt, sand, and gravel along valleys and coastal plains [33].

InSAR data show that land subsidence occurred regardless of the surface geology. Subsidence in CAMANAVA occurred mostly where there are recent alluvial and fluvial deposits, although parts in the northeast that are subsiding extend towards those underlain by the Diliman Tuff. Rosario, San Pedro and Biñan are also underlain by alluvial and fluvial deposits; Las Piñas is underlain by the Diliman Tuff, Dasmariñas by the Taal Tuff. Furthermore, the proximities of cones of depression and deformation patterns clearly indicate that the subsidence is anthropogenic, caused by groundwater extraction.

Fault-constrained subsidence

The West Marikina Valley Fault (WMVF), is the only officially recognized fault cutting through Metro Manila (Figure 6A). It is an oblique dextral fault trending almost N-S along the eastern edge of the metropolis [38]. The combined ascending and descending data show subsidence in the San Pedro and Biñan areas adjacent to the southern section of the WMVF (Figure 6B). The said fault constrains the areal extent of the deformation; the western portion of the deformation pattern trends N-S along the fault. Additionally, the south is bounded by a linear deformation pattern, which trends NE-SW and continues down to the coast of Laguna de Bay. A splay of the fault bounds this deformation, although at only a small section. The rest is not bounded by any known trace of the WMVF (Figure 6B).

The northeast portion of Figure 6A displays a fault splay that trends southwestward from the East Marikina Valley Fault (EMVF), to the WMVF. If extended, that trend coincides nicely with the sharply-defined southern edge of the CAMANAVA subsidence zone. The hypothesis that differential motion along that edge is also localized at a preexisting fault is worthy of future geophysical testing.

Ground fissures have emerged in parts of Muntinlupa, San Pedro, and Biñan where the southern segment of the WMVF is [4]. PHIVOLCS reported vertical slip rates of about 2-20 cm/ yr between 1996-1999, and concluded from the locations, extents, and geometries of the fissures that they are tectonically controlled, and are likely associated with the WMVF. Their slip rates are higher than known fault-related slip rates worldwide (e.g. Archaya 1980; Thatcher 1990; Rigenbach 1993), however, and deep wells in the area are drying up, suggesting that the exceptionally high differential movements are due predominantly to massive extraction of groundwater from industrial and domestic use [7].

We plot the vertical displacement map computed from the Envisat and ALOS 1 data with the trace of the WMVF mapped by PHIVOLCS [39] (Figure 6B). It appears that segments adjacent to San Pedro and Biñan still experienced differential movement during the observation period, while those in Muntinlupa did not. This may indicate either a temporary quiescence or possibly the depletion of groundwater supply. We note that the spatial extent of deformation during this time period suggests that the NE-SW creeping segments may be more extensive than what had been mapped by PHIVOLCS [7], and Rimando [4].

Land subsidence and implications to hazards

In 2011, Typhoon Nesat raised water levels up to about 1.28m above mean sea level along the coast of Metro Manila, with a peak storm surge of 0.78m [40]. Its nearest point to Manila was about 200 km north; by that time, it had already weakened to a Category 2 typhoon (Saffir-Simpson scale). Despite these, it still caused considerable damages because much of coastal Manila is low-lying. Should a Haiyan-strength tropical cyclone traverse Manila Bay, it could generate storm surges 2-4 meters high, which would cause much more damage to coastal communities [41].

Based on historical records, at least two tsunami events have affected Manila: 9 November 1828 and 3 June 1863 [42]. Water heights were estimated to have reached 1-2 m. Considering the water levels from Typhoon Nesat, tsunami waves this high could cause significant damage, notwithstanding other earthquakerelated effects such as liquefaction.

Land subsidence increases exposure to storm surges, tsunami, and flooding, tidal incursion, and potential decline in groundwater quality and availability in coastal communities. In areas with faults, ground fissures that emerged from fault reactivation threaten the structural integrity of buildings and structures built on top of them. Infrastructure projects intended to protect residents from these hazards will be useless without detailed knowledge of where and how fast the land subsides. Climate change is expected to increase the frequency of stronger typhoons [43], enhancing the impacts of floods and storm surges on coastal communities.

The combined ascending and descending InSAR data have established the rates and extents of land subsidence affecting Metro Manila and surrounding areas during 2003-2011. Large areas have subsided at rates of about 2-4.2 cm/yr, most notably the coastal cities of Manila and CAMANAVA. While the InSAR time-series data showed continued subsidence in most areas, subsidence rates in Rosario appears to have slowed down.

In general, subsidence is mostly linear, although in some areas the rates decreased towards the end of the study period; however, critical areas such as those in coastal Manila Bay, faulttraversed San Pedro and Biñan, and the rapidly urbanizing Las Piñas and Dasmariñas continued to subside. We attribute the land subsidence to over-extraction of groundwater extraction in these areas.

We thank the Philippines’ Department of Science and Technology and the USAID-funded PEER Science Program for supporting this research. The ALOS PALSAR 1 data are copyright of the Japan Aerospace Exploration Agency and the Japanese Ministry of Economy, Trade, and Industry, which were made available through the Alaska Satellite Facility.

1. Siringan F, Rodolfo KS. Relative sea level changes and worsening floods in the western Pampanga Delta: Causes and some possible mitigation measures. Science Diliman. 2013; 15: 1-12.

2. Rodolfo KS, Siringan FP. Global sea-level rise is recognised, but flooding from anthropogenic land subsidence is ignored around northern Manila Bay, Philippines. Disasters. 2006; 30: 118-139.

3. Lagmay AMF, Eco RN, Agdeppa J. PSInSAR detection of ground subsidence and fault movement in Muntinlupa City, MM and Binan, Laguna, In: J.C. Thouret, S.C. Liew and A. Gupta (eds) Remote sensing, natural hazards and environmental change. Center for Remote Sensing and Processing, NUS and Laboratoire Magmas et Volcans, CNRS, CLERVOLC, Universite Blaise Pascal, France. 2011; 63-66.

4. Rimando RE. Neotectonic and paleoseismic study of the Marikina Valley Fault System, Philippines. Dissertation, State University of New York. 2002.

5. Rodolfo KS, Siringan FP, Remotigue CT, Lamug CB. Worsening floods around northern Manila Bay, Philippines: Research-based analysis from physical and social science perspectives. Philippine Sociological Review. 2003; 51: 17-40.

6. Raucoules D, Cozannet GL, Woppelman G, de Michele M, Gravelle M, Daag A, et al. High nonlinear urban ground motion in Manila (Philippines) from 1993 to 2010 observed by DInSAR: Implications for sea-level measurement. Remote Sensing of Environment. 2003; 139: 386-397.

7. PHIVOLCS. Ground Fissures in Muntinlupa City and Vicinity: Tectonic or Man-Made? Unpublished report. 1999: 30.

8. Dong S, Samsonov S, Yin H, Ye S, Cao Y. Time-series analysis of subsidence associated with rapid urbanization in Shanghai, China measured with SBAS InSAR method. Environmental Earth Sciences. 2014; 72: 677-691.

9. Chaussard E, Amelung F, Abidin H, Hong SH. Sinking cities in Indonesia: ALOS PALSAR detects rapid subsidence due to groundwater and gas extraction. Remote Sensing of Environment. 2013; 128: 150-161.

10. Aobpaet A, Cuenca MC, Hooper A, Trisirisatayawong I. InSAR timeseries analysis of land subsidence in Bangkok, Thailand. Int J Remote Sensing. 2013; 34: 2969-2982.

11. Phien-Wej N, Giao PH, Nutalaya P. Land subsidence in Bangkok, Thailand. Engineering Geology. 2006; 82: 187-201.

12. Bell JW, Amelung F, Ramelli AR, and Blewitt G. Land subsidence in Las Vegas, Nevada, 1935-2000: New geodetic data show evolution, revised spatial patterns, and reduced rates. Environmental and Engineering Science. 2002; 8: 155-174.

13. Bell JW, Amelung F, Ferretti A, Bianchi M, Novali F. Permanent scatterer InSAR reveals seasonal and long-term aquifer-system response to groundwater pumping and artificial recharge. Water Resour Res. 2008; 44: W02407.

14. Amelung F, Galloway DL, Bell JW, Zebker HA, Laczniak RJ. Sensing the ups and downs of Las Vegas: InSAR reveals structural control of land subsidence and aquifer-system deformation. Geology. 1999; 27: 483- 486.

15. Dixon TH, Amelung F, Ferretti A, Novali F, Rocca F, Dokka R, et al. Space geodesy: Subsidence and flooding in New Orleans. Nature. 2006; 441: 587-588.

16. Chaussard E, Kerosky S. Characterization of Black Sand Mining Activities and Their Environmental Impacts in the Philippines Using Remote Sensing. Remote Sensing. 2016; 8: 100.

17. Rosen PA, Hensley S, Peltzer G, Simons M. Updated repeat orbit interferometry package released. Eos Transactions, American Geophysical Union. 2004; 85: 47-47.

18. Kampes B, Hanssen R, Perski Z. Radar Interferometry with Public Domain Tools, In: Proceedings of FRINGE 2003, December 1-5, Frascati, Italy. 2003.

19. Hooper A. A Statistical-Cost Approach to Unwrapping the Phase of InSAR Time Series. In: The FRINGE Workshop 09, ESA-ESRIN, Frascati (Rome), Italy, November 30–December 4, 2009. 2010.

20. Chen CW, Zebker HA. Two-dimensional phase unwrapping with use of statistical models for cost functions in nonlinear optimization. Journal of the Optical Society of America. A, Optics, image science, and vision. 2001; 18: 338-351.

21. Hooper A, Bekaert D, Spaans K, Arikan M. Recent advances in SAR interferometry time series analysis for measuring crustal deformation. Tectonophysics. 2012; 514-517.

22. Bekaert DPS, Walters RJ, Wright TJ, Hooper AJ, Parker DJ. Statistical comparison of InSAR tropospheric correction techniques. Remote Sensing of Environment. 2015; 170: 40-47.

23. Wright TJ, Parsons B, England PC, Fielding EJ. InSAR observations of low slip rates on the major faults of western Tibet. Science. 2004; 305: 236-239.

24. Salzer JT, Milillo P, Varley N, Perissin D, Pantaleo M, Walter TR. Evaluating links between deformation, topography and surface temperature at volcanic domes: Results from a multi-sensor study at Volcan de Colima, Mexico. Earth and Planetary Science Letters. 2017; 479: 354-365.

25. Soria JL, Siringan FP, Parreño PE. Compaction rates and paleo-sea levels along the delta complex north of Manila Bay, Luzon Island, Philippines. Science Diliman. 2005; 17.

26. Pielou EC. Fresh Water. University of Chicago Press. 1998.

27. Theis CV. The significance and nature of the cone of depression in ground-water bodies. Economic Geology. 1938; 33: 889-902.

28. Terzaghi K. Erdbaumechanik auf bodenphysikalischer Grundlage. Leipzig u. Wien: F. Deuticke. 1925.

29. Tolman CF, Poland JF. Ground-water, salt-water infiltration, and ground-surface recession in Santa Clara Valley, Santa Clara County, California. Eos, Transactions American Geophysical Union. 1940; 21: 23-35.

30. Poland JF. Guidebook to studies of land subsidence due to groundwater withdrawal. International Hydrological Programme, Working Group 8.4, UNESCO. 1984.

31. Galloway WE. Cenozoic evolution of sediment accumulation in deltaic and shore-zone depositional systems, northern Gulf of Mexico Basin. Marine and Petroleum Geology. 2001; 18: 1031-1040.

32. Piquero E. Identification of Groundwater Critical Areas in Metro Manila and Vicinity and Formulating an Integrated Approach to Aquifer Recharge and Protection. 2005.

33. Mines and Geosciences Bureau. Geology of the Philippines. Mines and Geosciences Bureau, Quezon City, Philippines, 2nd Edition. 2010.

34. Catane S, Arpa M. Large-scale eruptions of Laguna caldera: contributions to the accretion and other geomorphic developments of Metro Manila and provinces. Philippine Institute of Volcanology and Seismology. Internal report. 1998.

35. Arpa MCB, Patino LC, Vogel TA. The basaltic to trachydacitic upper Diliman Tuff in Manila: Petrogenesis and comparison with deposits from Taal and Laguna Calderas. J Volcanology and Geothermal Res. 2008; 177: 1020-1034.

36. Gonzales B, Ocampo V, Espiritu E. Geology of Southeastern Nueva Ecija and Eastern Bulacan Provinces, Luzon Central Valley. J Geol Soc Philipp. 1971; 25: 3-41.

37. Bureau of Mines and Geo-sciences. Geological Map of Silang Quadrangle, Sheet 3162 I. 1987.

38. Rimando RE, Knuepfer PL. Neotectonics of the Marikina Valley fault system (MVFS) and tectonic framework of structures in northern and central Luzon, Philippines. Tectonophysics. 2006; 415: 17-38.

39. PHIVOLCS. The Valley Fault System in Greater Metro Manila Area Atlas. 2014.

40. Morin VM, Warnitchai P, Weesakul S. Storm surge hazard in Manila Bay: Typhoon Nesat (Pedring) and the SW monsoon. Natural Hazards. 2016; 81: 1569-1588.

41. Lapidez JP, Tablazon J, Dasallas L, Gonzalo LA, Cabacaba KM, Ramos MMA, et al. Identification of storm surge vulnerable areas in the Philippines through the simulation of Typhoon Haiyan-induced storm surge levels over historical storm tracks. Nat Hazards Earth Syst Sci. 2015; 15: 1473-1481.

42. Bautista MLP, Baustista BC, Salcedo JC, Narag IC. Philippine tsunamis and seiches (1589-2012). Philippine Institute of Volcanology and Seismology. 2012.

43. Mei W, Xie SP. Intensification of landfalling typhoons over the northwest Pacific since the late 1970s. Nature Geoscience. 2016; 9: 753-757