Spatiotemporal Analysis of Extreme Rainfall and Meteorological Drought Events over the Angat Watershed, Philippines

Abstract

Understanding the spatiotemporal distribution of extreme rainfall and meteorological drought on a watershed scale could be beneficial for local management of any water resources system that supports dam operation and river conservation. This study considered the watershed of Angat as a case, given its economic importance in the Philippines. A series of homogeneity tests were initially conducted on each rainfall dataset from monitoring stations in and near the watershed, followed by trend analysis to determine the rate and direction of change in the annual and seasonal rainfall extreme indices in terms of intensity, duration, and frequency. Three indices, using the rainfall deviation method (%DEV), percent of normal rainfall index (PNRI), and Standardized Precipitation Index (SPI), were also used to identify meteorological drought events. Generally, rainfall in the watershed has an increasing annual PCPTOT (4–32 mm/year), with increasing frequency and intensity in heavy rainfall and wet days. A significant increasing trend (α = 5%) in the seasonal PCPTOT (7–65 mm/year) and R10mm (1.7–10.0 days/decade) was particularly observed in all stations during the Amihan Monsoon Season (Dec–Feb). The observed increasing rainfall intensity and frequency, if it continues in the future, could have an implication both for the water resources operation to satisfy the multiple objectives of Angat Reservoir and for the flood operation that prevents damage in the downstream areas. The effect of each ENSO (El Niño- Southern Oscillation) phase on the rainfall is unique in magnitude, intensity, and duration. The seasonal reversal of the ENSO in the extreme rainfall and meteorological drought signals in Angat Watershed was also evident. The identified meteorological drought events in the watershed based on SPI-12 persisted up to 12–33 months, could reduce more than 60% (PNRI < 40%) of the normal rainfall. Insights from the study have implications for the hydrology of the watershed that should be considered for the water resources management of the Angat Reservoir.

1. Introduction

The Angat Watershed plays a pivotal role in the provision of water supply of the Philippine’s National Capital Region (NCR) and surrounding provinces. It is considered as one of the most important water resources systems in the country, not just for water supply but also for power generation and flood control. The multipurpose Angat Dam in the watershed allocates water for (i) domestic consumption (through the Ipo Dam) that accounts for more than 90% of the NCR’s water supply [1,2], (ii) irrigation use (through the Bustos Dam) for about 17,500 hectares (ha) and 24,000 ha of farmland during the wet and dry cropping seasons, respectively [3], (iii) hydropower generation with a total capacity of 218 megawatts (MW) [4], and (iv) flood operation, which prevents damage in the downstream areas along the Angat River and the afterbay re-regulation dams of the Bustos and Ipo Dams [5]. Flood control operation is extended to Barangay Tibag in Pulilan, Bulacan, which is the target watershed in the Angat Dam Flood Operation Rule.

The annual average natural inflow of the watershed is about 58.40 cubic meters per second (m³/s), but since the operation of the Umiray–Angat Transbasin Project (UATP) in June 2000, the reservoir receives an additional average annual flow of 11.60 m³/s, raising the average inflow to the dam to 70.0 m³/s. Rainfall over the Angat Watershed is a key factor influencing the inflow and water availability in the Angat Reservoir [6,7,8,9]. Pre-UATP studies [10,11] revealed that the watershed inflow to the reservoir has strong sensitivity to local rainfall variability. The rainfall in the watershed follows a distinct pattern of wet and dry seasons and is greatly influenced by the monsoons and various weather systems in the Philippines, such as tropical cyclones [12]. There is also the influence on the inter-annual climatic variability of El Niño Southern Oscillation (ENSO), a large-pattern atmospheric and ocean circulation process that triggers the fluctuation between above-normal (El Niño) and below normal (La Niña) sea surface temperature (SST) over the tropical Pacific [13,14,15]. These seasonal and inter-annual climatic variation adversely affected both the agricultural production [16,17,18], and water resources [19,20,21,22] in the country. In the case of Angat Watershed, the heavy rainfall during the wet season is crucial for replenishing the water supply in the Angat Reservoir [10], while during the dry season, rainfall variations over the watershed can have a significant impact on water availability and the risk of drought [19]. Given the economic worth of Angat Watershed in the country, there is a need for a thorough analysis of the spatiotemporal analysis of its local rainfall to properly assess the availability of water for its competing users, and effectively mitigate water-related hazards downstream.

In recent years, there has been a growing emphasis on studying extreme rainfall trends, particularly in the context of climate variability and climate change. Most peer-reviewed studies related to rainfall variability in the country were done on a nationwide [23,24,25,26,27,28,29] or regional scale [30,31,32,33]. According to the Intergovernmental Panel on Climate Change (IPCC) Fifth Assessment Report [34], an increase of 22 millimeters (mm)/decade on the total wet-day rainfall (PCPTOT) was reported over Southeast Asia, but trends vastly differ across the region and seasons. Endo et al. [35] who determined trends of extreme rainfall over Southeast Asia from the 1950s to 2000s based on the indices defined by the Expert Team on Climate Change Detection and Indices (ETCCDI), reported an increasing trend in the annual PCPTOT in southern Vietnam and Luzon Island in the Philippines. The same study observed an increase in heavy rainfall (R50mm) in most stations for the Luzon and Visayas islands. Overall, an increase in heavy rainfall was reported over Southeast Asia in the latest Assessment Report of the IPCC [36,37]. Understanding the dynamics of rainfall trends on this scale empowered national authorities and policy makers to make informed decisions related to disaster risk reduction, water resource management, and sustainable development [23,24,38]. While these findings are useful to get the bigger picture, a watershed-scale analysis could be most beneficial for the local management of any water resources system that supports dam operation and river conservation.

The continuous negative rainfall anomaly, on the other hand, is usually associated with meteorological drought. Unlike extreme rainfall events that have an associated water-related hazard (e.g., flood), and are abrupt and almost near real-time, the risk brought by meteorological drought could manifest over an extended period. Assessment of drought using multiple indexes and indicators is recommended by the World Meteorological Organization (WMO) [39]. In the Philippines, meteorological droughts were monitored and reported using the Percent of Normal Rainfall Index (PNRI) and the Standard Precipitation Index (SPI) [40]. Past drought events usually coincided with El Niño events [19,24]. Competition between the irrigation and domestic water allocation of Angat Reservoir caused farmers of AMRIS to cope with an “imposed water scarcity” especially during El Niño events [41]. It is evident with the past El Niño events where minimum releases were made for irrigation to give priority to domestic water use [6,8,42,43,44,45].

This study was carried out to assess the historical spatiotemporal variability of extreme rainfall and meteorological drought over the Angat Watershed, given its economic importance in the country. Specifically, this study aims to characterize rainfall in the watershed by conducting (1) trend analysis on the different extreme rainfall indices and (2) drought analysis using a rainfall deviation method, percent of normal method, and Standardized Precipitation Index (SPI). The coherency of the results with the national and regional context was included in the discussions.

2. Materials and Methods

2.1. Study Area and Preliminary Data Analysis

Angat watershed in Bulacan, Philippines, lies between 121°9′ and 121°21′ east longitudes, and 14°50′ and 15°12′ north latitudes (Figure 1). Angat Watershed elevation ranges between 183 and 1257 meters (m) above mean sea level (AMSL), with an average elevation of 503 m AMSL and hillslope of 23.4%. The watershed stores the Angat Reservoir and has a total drainage area of 546.5 square kilometers (sq.km) upstream of the Angat Dam—a 131 m high rock-fill dam with a 630 m long dam crest that curves at a radius of 620 m [5]. The watershed is as one of the most well-forested and managed forest reserves in the country, with 95.5% (506.62 sq.km) of the area covered with combined forest and brush, while the remaining 4.5% (24.20 sq.km) comprises the inland water of Angat Reservoir and upstream tributaries. It consists of three sub-watersheds, namely, the Upper Angat or Maputi watershed, Matulid Watershed, and Talaguio Watershed. The Sierra Madre Mountain Range on the eastern side of Angat watershed separates it from Umiray Watershed. Additional water from Umiray River is conveyed through a 13-kilometer (km) tunnel (4.30 m diameter) to Macua River across the Sierra Madre Mountains.

The study considered five (5) rainfall monitoring stations located within and outside the Angat Watershed (Figure 1). The stations of Talaguio, Maputi, Angat, and Matulid are maintained by the National Power Corporation (NPC), while Umiray station is monitored by the Metropolitan Waterworks and Sewerage System (MWSS). These stations are located in different climate zones: Type I with pronounced dry season from November to April, and wet during the rest of the year, which includes the stations in Talaguio and Angat; Type III, with a relatively dry season from November to April, and wet during the rest of the year, which includes the station in Matulid and Maputi; and lastly, Type IV is characterized with more or less even rainfall distribution throughout the year, which includes the stations of Umiray. No station belongs to Type II, which has no dry season, with a very pronounced maximum rain from December to February. The Philippine Atmospheric, Geophysical and Astronomical Services Administration (PAGASA) and the national meteorological and hydrological services (NHMS) agency of the Philippines classify these four climate zones following the Modified Corona Climate Classification based on average monthly rainfall totals [21,46,47].

The raw rainfall datasets were initially subjected to quality control through homogeneity assessment prior to any analysis [48]. Missing values were observed at all stations except Angat station (

Table 1. Descriptions of the monitoring rainfall stations in Angat Watershed.

Station	Climate Type *	Elevation (m AMSL)	Temporal Duration	No. of Years	% Missing Values
Maputi	III	557.99	1987–2021	35	22.67
Talaguio	I	644.38	1987–2021	35	7.07
Matulid	III	494.33	1987–2021	35	19.54
Angat	I	323.55	1987–2021	35	-
Umiray	IV	258.69	2022–2021	20	3.11

* Based on the Modified Corona Climate Classification.

), as there are no road networks to immediately access and regularly secure the other stations. To address the issue on data quality, two step data imputation were performed to treat missing data: (1) autoregressive integrated moving average (ARIMA) with Kalman smoothing for 10-day missing continuous gap; and (2) Inverse Distance Weighing (IDW) method using at least (2) nearby stations’ continuous recorded data for more than 10 days continuous data gap.

In this study, four homogeneity tests were employed, following the studies of Fernandez et al. [19], to identify the break periods of each rainfall series, namely: Pettitt’s Test, Standard Normal Homogeneity Test (SNHT), Buishand’s Range Test, and Von Neumann Ratio Test. The resultswere categorized into three classes based on the number of tests rejecting the null hypothesis at 5% significance level [49]

• Class 1 (Useful), which rejects one or zero null hypothesis of the four tests. Under this category, the series is considered homogenous and can be used for further analysis.
• Class 2 (Doubtful), which rejects two null hypotheses of the four tests with indication of inhomogeneity in the series. Results of trend analysis should be critically inspected due to possible presence of inhomogeneities.
• Class 3 (Suspect), when three or all null hypotheses are rejected. In this category, the series lacks credibility. Rainfall series that fall under this category were compared to the nearby stations for cross-validation of the inhomogeneity of the series. Necessary adjustments to the identified break point/s were made.

The continuous records of each station were used to estimate the total (PCPTOT) and coefficient of variation (CV) of rainfall across timescale (monthly, seasonal, and annual). The degree of variability based on CV is as follows: low for CV< 0.2, moderate for 0.2 < CV < 0.3, and high for CV > 0.3. To get a better perspective of the distribution of the PCPTOT and its variability over the watershed, spatial maps were generated using the IDW method. The seasonality index (SI) and the precipitation concentration index (PCI_annual) were also determined to understand the seasonality and shift in rainfall concentration patterns in the watershed:

$$\mathrm{S}\mathrm{I}=\frac{1}{{\mathrm{R}}_{\mathrm{j}}}\sum _{\mathrm{i}=1}^{12}\left({\mathrm{X}}_{\mathrm{i}}-\frac{{\mathrm{R}}_{\mathrm{j}}}{12}\right)$$

(1)

$${\mathrm{P}\mathrm{C}\mathrm{I}}_{\mathrm{a}\mathrm{n}\mathrm{n}\mathrm{u}\mathrm{a}\mathrm{l}}=\frac{\sum _{\mathrm{i}=12}^{12}{\mathrm{X}}_{\mathrm{i}}^2}{\sum _{\mathrm{i}=12}^{12}{\left({\mathrm{X}}_{\mathrm{i}}\right)}^2}$$

(2)

where ${\mathrm{R}}_{\mathrm{j}}$ is the total annual rainfall for a particular year j while ${\mathrm{X}}_{\mathrm{i}}$ is the total rainfall for ${\mathrm{i}}^{\mathrm{t}\mathrm{h}}$ month of year j. SI measures the spread of the monthly rainfall with respect to an ideally uniform monthly distribution in a year. SI values can range from 0 (all months have similar rainfall) to 1.83 (all rainfall incidences occur in a single month). PCI_annual indicates the degree of rainfall concentration during a given year. The values of PCI_annual can be lower than 10 (uniform rainfall distribution) to over 20 (strongly irregular rainfall distribution), with maximum value up to 100.

To describe the occurrence of extreme rainfall events, the study considered nine of the eleven rainfall extreme indices of ETCCDI [50] developed and recommended by the joint Commission for Climatology/World Climate Research Program/Joint Technical Commission for Oceanography and Marine Meteorology (CCI/WCRP/JCOMM) as shown in

Table 2. Description of rainfall extreme indices used in the study.

Index	Definition *	Units
PCPTOT	Total amount of rainfall in wet days	mm
RX1day	Maximum 1-day rainfall	mm
RX5day	Maximum 5-day rainfall	mm
CDD	Maximum length of dry days	day
CWD	Maximum length of wet days	day
R10mm	Count of days when rainfall ≥ 10 mm	day
R20mm	Count of days when rainfall ≥ 20 mm	day
R75mm	Count of days when rainfall ≥ 75 mm	day

* Wet days refer to days with rainfall ≥ 1 mm, while dry days refer to when rainfall < 1 mm.

. These indices consider the three major attributes of extreme events, namely: intensity (PCTTOT, RX1day, and RX5day), duration (CDD and CWD), and frequency (R10mm, R20mm, and R75mm). This study particularly added R75mm to the analysis because of the established daily threshold of 75 mm for any of NPC-monitored stations as an indicator for the beginning of the flood pre-caution period of downstream areas of Angat Dam [5].

2.2. Trend Analysis

Detecting significant monotonic trends present in each series of rainfall extremes was performed using the Mann–Kendall (MK) test [51,52] with Sen’s Slope (SS) estimator [53]. The application of these statistical tests had been demonstrated on previous works [24,26,29,33,35] for Philippine setting. In general, the null hypothesis that the series come from a sequence with independent realizations and are identically distributed (no trend) is evaluated against the alternative hypothesis, which implies a presence of significant trend. This is done at both the 5% and 10% significance levels.

2.3. Drought Analysis

Three drought indices using the percent of normal rainfall index (PNRI), rainfall deviation index (%DEV), and Standardized Precipitation Index (SPI) were employed in this study to describe the past meteorological drought events in the watershed.

2.3.1. Rainfall Deviation and Percent of Normal Rainfall Index

Meteorological drought can be understood from the rainfall deviation from the long-term mean of a rainfall series. The rainfall deviation index (%DEV) and percent of normal rainfall index (PNRI) in monthly scale was calculated using the formula:

$$\mathrm{\%}\mathrm{D}\mathrm{E}\mathrm{V}=\frac{X_i-X_{n\left(i\right) }}{X_{n\left(i\right) }}\times 100$$

(3)

$$\mathrm{P}\mathrm{N}\mathrm{R}\mathrm{I}=\frac{X_i}{X_{n\left(i\right) }}\times 100$$

(4)

where $X_i$ is the total rainfall for $i^{th}$ month, while $X_{n\left(i\right)}$ is the climate normal of total rainfall for the $i^{th}$ month covering the last 30 years (1991–2020). Station-based and watershed mean $X_{n\left(i\right)}$ was established following the guidelines of WMO [54]. PNRI has been used extensively to indicate past drought conditions in the Philippines. PNRI is used mainly by the DOST-PAGASA for monthly climate assessment and issuance of probabilistic rainfall forecasts [40]. The following classification suggested by the DOST-PAGASA was adopted in this study: rainfall way below normal for PNRI < 40% (red), rainfall below normal for 41% < PNRI < 80% (yellow), rainfall near normal for 80% < PNRI < 120% (green), and rainfall above the normal for PNRI > 120% (blue). The study of de los Reyes and David [20] use %DEV to describe the onset, duration, and magnitude of drought events in the Philippines for post-2000 El Niño events. Based on PNRI, drought event was defined as when three consecutive months are observed to have way below normal rainfall (red) or when five consecutive months have below normal rainfall (yellow) [40] (

Table 3. Criteria for drought, dry spell, and dry conditions based on PNRI and %DEV.

Classification	Criteria
Drought	3 months of consecutive records of way below normal rainfall (equivalent to less than −60%DEV); or 5 months of consecutive records of below normal rainfall (equivalent to −21 to −60%DEV)
Dry Spell	2 months of consecutive records of below normal rainfall (equivalent to less than −60%DEV); or 3 months of consecutive records of way below normal rainfall (equivalent to −21 to −60%DEV);
Dry Condition	2 months of consecutive records of below normal rainfall (equivalent to −20 to −60%DEV)

).

2.3.2. Standardized Precipitation Index (SPI)

The SPI developed by McKee et al. [55] was endorsed by the WMO in 2009 to be used by all NHMS around the world as the standard for determining the existence of meteorological drought [56,57] because of its applicability for various climate regime and computational simplicity. The total monthly rainfall series of each station were initially fitted with gamma distribution using coefficients by Thom [58] to compute for the cumulative distribution function (CDF) that is converted to Z-score distributions [55] to obtain the monthly SPI index. SPI was developed upon the relationship of drought frequency, duration, and scales [55]. A drought event is defined as the period in which SPI is continuously negative and reaches a value of −1.0 or less [57]. A drought event starts when the SPI first drops below zero and ends with a positive SPI value following an SPI of −1 or less. A drought event magnitude (positive sum of the SPI values for all months within a drought event) and intensity (ratio between drought magnitude and its duration) [59] are then calculated. This study calculated the SPI values on various timescales (3-, 6-, and 12-month), but focused on the analysis of SPI-12 since this period is normally tied to streamflow and reservoir levels [57] and is commonly used in the assessment of hydrologic impacts of meteorological droughts [39].

2.4. The Teleconnection of El Niño–Southern Oscillation (ENSO)

The influence of the large-pattern atmospheric and ocean circulation of ENSO using the Oceanic Niño Index (ONI) 3.4 on the detrended series of extreme rainfall and drought indices over Angat watershed was examined through correlation analysis. This study uses the Spearman correlation at 5% level of significance as in other studies [60,61]. This study conducted thorough discussions on the characteristics of meteorological drought events on the Angat Watershed based on SPI-12 and PNRI after the onset or development of warm episodes in tropical Pacific (El Niño) from 1990 to 2020. The monthly teleconnection index of ONI 3.4 was collected from the Climate Prediction Center of the National Oceanic and Atmospheric Administration (NOAA, https://origin.cpc.ncep.noaa.gov/products/analysis_monitoring/ensostuff/ONI_v5.php, accessed on 1 September 2023).

3. Results and Discussion

3.1. Rainfall Characteristics

Table 4. Rainfall time series with significant break periods (α = 5%).

Station	Timeseries	Classification	Pettitt’s Test	SNHT	Buishand’s Test	Von Neumann’s
Angat	January	Doubtful	2006	2021	2006	0.2141
Maputi	January	Doubtful	2010	2017	2010	0.0492
Matulid	January	Doubtful	1999	1999	1999	0.2437
Matulid	March	Doubtful	1999	1999	1999	0.0127
Matulid	February	Suspect	2000	2000	2000	0.0801
Matulid	DJF	Suspect	1999	1999	1999	0.0582

shows the list of rainfall stations with identified break periods. The rainfall series of Talaguio and Umiray stations are found to be homogenous across timescale. The rainfall series of Angat, Maputi, and Matulid for the month of January falls under Class 2, with identified break periods in 2006, 2010, and 1999, respectively. However, these were found to be not associated with inhomogeneity after cross-validating with homogenous series of the same month from other stations. The February and December–January–February (DJF) rainfall series of Matulid station falls under Class 3 and is found to have a significant break period in 1999 and 2000, respectively. Figure 2 shows the identified breakpoint of the DJF rainfall series of Matulid station along with nearby homogenous series, which is disregarded for being high as compared to other stations and is corrected using the IDW method.

The mean PCI_annual values over the study period varies from 12.70 in Maputi to 15.56 in Angat station (

Table 5. Annual mean statistics and extreme indices covering the study period of each station.

Station	PCI	SI	CV	PCPTOT	RX1day	RX5day	CDD	CWD	R10mm	R20mm	R75mm
Angat	15.56	0.76	0.19	2619	177	334	37	20	66	40	5
Maputi	12.70	0.58	0.23	3919	224	409	16	34	95	55	8
Matulid	13.93	0.65	0.26	3734	238	444	21	27	88	53	9
Talaguio	14.78	0.72	0.21	3067	229	402	24	25	71	44	6
Umiray	12.71	0.54	0.16	5478	272	538	15	31	127	78	15

Note: PCPTOT, RX1day, and RX5day are in mm/year; CDD, CWD, R10mm, R20mm, and R75mm are in day/year.

), indicating a moderately distributed total annual rainfall taking place during at least six months of any year. This happened for more than 80% of the observation years considered for the Maputi and Umiray stations (

Table 6. Percent of the observation years per PCI classification of each station.

PCIannual	Classification	Angat	Maputi	Matulid	Talaguio	Umiray
<10	Uniform	-	-	2.86	-	4.76
10–15	Moderate	48.57	88.57	68.57	60.00	80.95
16–20	Irregular	42.86	11.43	28.57	37.14	14.29
>20	Strongly Irregular	8.57	-	-	2.86	-

). Furthermore, almost all the stations exhibit SI values less than 0.79, which confirms the presence of seasonality in rainfall over the watershed, with a short drier season particularly for the stations of Maputi (SI = 0.58) and Umiray (SI = 0.54).

The spatial mean annual and seasonal rainfall variations over the watershed are shown in Figure 3. The annual mean total rainfall (PCPTOT) over the Angat watershed is 3733 mm, with a CV of 0.24 implying a moderate interannual rainfall variability. The Amihan (DJF) season or the Philippines northeast monsoon contributes to the rainfall in the watershed’s eastern side, ranging from 209 to 1610 mm, or 18.66% of the annual total rainfall. The strong influence of the Philippines southwest monsoon or Habagat season can be seen on the western side of the watershed, which usually occurs in June, with its peak in August (JJA). This contributes a substantial amount of rainfall, especially in the western part of the watershed, ranging from 1125 to 1366 mm, or 34.48% of the annual rainfall. The summer season covers March–April–May (MAM), and the monsoon transition season occurs in September–October–November (SON); they have the lowest and highest contribution to the total annual rainfall of 11.70 and 35.16%, respectively. The variability of the seasonal rainfall in the watershed based on CV can be classified as high across seasons (0.34 to 0.59). The results showed that the seasons with the least amount of rainfall (DJF and MAM) showed higher interannual variability compared to seasons with usually heavy rainfall (JJA and SON).

The highest recorded one-day and five-day cumulative rainfall from the NPC-monitored stations over the study period are both observed in Matulid station: 510 mm during Typhoon Ruby in October 1988 (RX1day) and 807 mm brought by Typhoon Kai-Tak in July 2000 (RX5day). During these times, cold episodes were observed in the tropical Pacific, prompting a La Niña alert. On average, the side of the Matulid station receives most of the intense rainfall within the watershed, with an annual average RX1day and RX5day of 238 mm and 444 mm, respectively (Table 5). Maputi station in the northern part of the watershed is the wettest, with records on average of 34 consecutive wet days (CWD) and 16 consecutive dry days (CDD) per year. On average, the eastern side of the watershed, which falls under a Type III climate, receives more frequent heavy rainfall (R10mm, R20mm, R75mm): Maputi station (95, 55, 8 days/year) and Matulid station (88, 53, 9 days/year). Regular maintenance of these isolated monitoring stations is therefore necessary for flood dam operation, as well as for watershed management (e.g., the monitoring of rainfall erosivity for reservoir sediment management).

3.2. Trend Analysis of Seasonal and Annual Rainfall Series

Seasonal and annual trends of rainfall extremes are analyzed for each monitoring station of the Angat Watershed using the MK trend test with SS estimator at a confidence interval of 90% or higher. Autocorrelation for all rainfall series proved to be statistically insignificant.

Table 7. The rate of change based on the Sen’s Slope estimator of rainfall extreme indices per station. The sign of the estimates indicates if the annual/seasonal series shows increasing (positive values) or decreasing (negative values) trend over the study period.

Time Scale	Station	PCPTOT	RX1day	RX5day	CWD	CDD	R10mm	R20mm	R75mm
Annual	Angat	13.07	0.66	2.27	0.1	−0.27 **	0.33 **	0.13	0.08
	Maputi	31.83 **	0.58	2.32	0.2	0	0.86 *	0.58 *	0.09 **
	Matulid	7.44	−0.89	0.44	−0.47 *	−0.04	0.05	0.14	0
	Talaguio	3.87	0.44	0.32	−0.04	−0.11	0.25	0.11	−0.06
	Umiray	40.85	−4.61	−5.07	0.37	−0.33	1.44 *	1 *	0.08
Amihan	Angat	7.35 *	1.1 *	2.3 **	0.05	−0.5 *	0.17 *	0.1 *	0
Season	Maputi	16.48 *	2.69 *	4.09	0.15 *	−0.04	0.38 *	0.28 *	0
	Matulid	19.41 *	2.52 *	5.59*	−0.13	−0.13	0.42 *	0.29 *	0
	Talaguio	7.26 *	1.07	2.28	0	0	0.2 *	0.09 *	0
	Umiray	64.28 *	6.92 *	16.78 *	0.37	−0.07	1.09 *	0.79	0.33 *
Summer	Angat	2.29	0.43	−0.18	0	0.04	0.05	0	0
Season	Maputi	−2.11	0	−1.25	0.06	0	0	−0.04	0
	Matulid	3.96	0.72	0.6	0	−0.1	0.06	0	0
	Talaguio	0.05	−0.65	−0.79	0.03	0	0.05	0	0
	Umiray	3.56	1.09	−0.33	0	−0.25	0.36	0	0
Habagat	Angat	2.37	0.7	1.46	0.11	0	0.08	0.06	0
Season	Maputi	1.29	−0.45	−0.17	0	0	0.16	0.1	0
	Matulid	−10.75	−0.64	−1	0	0.04	−0.26	−0.13	0
	Talaguio	−0.14	0.22	−0.71	0.18	0	0.04	0	0
	Umiray	1.52	−0.41	−2.03	0.43 **	0	0.19	−0.18	0
Monsoon	Angat	−1.04	−0.33	0.86	−0.08	−0.05	−0.04	−0.04	0
Transition	Maputi	14.27	1.18	3.33	0.14	0	0.33 *	0.25 *	0.08 *
Season	Matulid	−0.76	−1	0.38	−0.17 **	0	0	0	0
	Talaguio	−6.31	0.57	1.16	-0.15 **	0	-0.07	-0.06	0
	Umiray	−14.03	−8.27	−7.2	0.14	0	0.29	0.17	0

* Significant at α = 5% based on MK trend test; ** Significant at α = 10% based on MK trend test. Note: PCPTOT, RX1day, and RX5day are in mm/year; CDD, CWD, R10mm, R20mm, and R75mm are in day/year.

shows the annual and seasonal rate of change of each rainfall extreme indices per station, supplemented with the following illustrations: Figure 4 shows the trend direction of the seasonal and annual extreme indices associated with the rainfall intensity (PCPTOT, RX1day, RX5day), Figure 5 for trend directions of indices associated with rainfall duration (CDD and CWD), and Figure 6 for trend directions of indices associated with frequency (R10mm, R20mm, R75mm). In each figure, an increasing rate is represented by a (+) sign, and a significant increasing rate with (▲) symbol; while, decreasing rate are shown using a (−) sign, significant decreasing rate with a (▼) symbol, and no trend with (●) symbol.

All NPC-monitored rainfall stations have an increasing trend of annual PCPTOT ranging from 3.87 to 31.83 mm/year over the 35-year study period (1987–2021). A significant trend at a 90% confidence interval was observed particularly in the Maputi annual series, which holds the highest rate of change among the stations (Figure 4). Umiray station, with a 20-year temporal record (2001–2021), also shows an increasing but not statistically significant trend at a rate of 40.85 mm/year. On average, the watershed area rainfall is increasing at a rate of 18.95 mm/year.

A significant increasing trend (α = 5%) in the seasonal PCPTOT was observed during the Amihan Monsoon season (DJF) for all stations, with an area watershed average rate of 15.87 mm/year. The rate of increase of Amihan total rainfall was observed to intensify at the eastward side of the watershed. Angat (7.35 mm/year) and Talaguio (7.26 mm/year) stations, located on the western side and with a Type I climate, have the lowest rate of increase, followed by the Maputi (16.48 mm/year) and Matulid (19.41 mm/year) stations, which are under a Type III climate. Umiray station, with a Type IV climate, has a significant total rainfall increase at a rate of 64.28 mm/year during the Amihan season. The finding Amihan Monsoon season is contrary to the findings for the monsoon transition period (SON), where four out of five (4/5) stations show a decreasing but insignificant trend (α = 5%) in the seasonal PRCTOT. Maputi station has the only increasing rate (14.27 mm/year) of SON seasonal PRCTOT, while the rest are decreasing (−0.76 to −14.03 mm/year). During the Habagat Monsoon season (JJA), insignificant trends were observed: decreasing for Matulid (−10.75 mm/year) and Talaguio (−0.14 mm/year), while increasing for Angat (2.37 mm/year), Maputi (1.29 mm/year), and Umiray (1.52 mm/year). Summer seasonal rainfall (MAM) has an increasing but insignificant trend in four out of five (4/5) monitoring stations. Overall, there is no spatial coherency in terms of the seasonal trend direction of PRCPTOT in the watershed except during Amihan season. Same findings were observed for other rainfall indices.

An increase in the rainfall intensity during the Amihan season is also observed in all stations in terms of RX1day and RX5day (Figure 4). A significant increasing trend (α = 5%) is specifically observed for seasonal RX1day (RX5day), with a rate of 2.52 mm/year (5.59 mm/year) for Matulid and 6.92 mm/year (16.78 mm/year) for Umiray. On average, the Angat watershed received an increasing RX1day and RX5day brought by the Amihan monsoon of 2.35 mm/year and 4.29 mm/year, respectively. The Habagat monsoon rainfall intensity, on the other hand, shows a declining but insignificant trend for RX1day (3 stations) and RX5day (4 stations). The watershed decreases of Habagat’s RX1day and RX5day are 0.37 mm/year and 0.44 mm/year, respectively. During the summer, watershed area average rainfall intensity has an increasing RX1day (0.19mm/year) but decreasing RX5day (0.50 mm/year). Annually, an increasing trend for RX1day (3 stations) and RX5day (4 stations) was observed over the watershed.

For the rainfall duration over the watershed, the results show an overall decreasing annual trend on consecutive dry days (CDD), which is significant (α = 5%) for Angat station (Figure 5). The decrease on dry days was observed mostly in Amihan season for four out of five (4/5) stations. No trend in CDD was observed from the JJA and SON seasons, except the statistically insignificant increasing trend in Matulid (JJA) and a decreasing trend in Angat station (SON). The annual consecutive wet days (CWD) are increasing for the stations of Angat, Maputi, and Umiray, but decreasing for the stations of Talaguio and Matulid. Maputi station has an increasing rate of CWD across seasons.

An increase in the annual frequency of heavy rainfall (R10mm and R20mm) is observed in all stations as shown in Figure 6. Particularly, a significant annual upward trend (α = 5%) is particularly detected for the R10mm (R20mm) of Maputi station, with a rate of 8.6 days/decade (5.8 days/decade), and for Umiray with 14.4 days/decade (10 days/decade). A significant (α = 5%) increasing trend in the seasonal R10mm is observed in all stations (1.7 to 10.09 days/decade) during the Amihan monsoon season and in the Maputi station (2.5 days/decade) during monsoon transition season. For R75mm, there is no trend observed for all stations from March to August, but a significant increasing trend (α = 5%) for Maputi station during the rainy season of SON, which reflects its annual R75mm series.

3.3. Meteorological Drought Analysis in Angat Watershed

Quantifying the occurrence of the monthly drought requires establishing the rainfall departure from the long-term mean over the Angat Watershed. Figure 7 shows the long-term monthly total rainfall per each NPC-monitored station, as well as for the watershed areal average computed using both the IDW method and Thiesen polygon method.

Rainfall starts to increase at the end of the summer season (May) throughout the watershed. Stations located in a Type I climate (Talaguio and Angat) peak only in the month of August, while stations located in Type III (Maputi and Matulid) have two maxima every August and November. In general, the watershed areal average rainfall follows the monthly behavior of Maputi (Thiessen weight of 0.478) and Matulid (Thiessen weight of 0.346), which peaks twice a year. The drought indices of PNRI (Figure 8) were established using the normal monthly rainfall over the Angat watershed.

There were eleven (11) El Niño events which lasted from 7 up to 18 months over the study period (Figure 9a). Two very strong El Niño events were recorded from May 1997 to May 1998, and November 2014 to April 2016; two strong events in July 1987 to February 1988 and May 1991 to June 1992; three moderately strong events in September 1994 to March 1995, June 2002 to February 2003, and July 2009 to March 2010; three weak events happened on July 2004 to February 2005, September 2006 to 2007, and September 2018 to June 2019. On average, there are 34 months intervals between El Niño events, but it still varies from 18 up to 50 months.

Droughts in the Philippines usually coincide with El Niño events [19,62], and Angat Watershed is not an exemption. During the manifestation of the very strong 1997–1998 El Niño, a drought event based on the PNRI classification was observed from October 1997 to February 1998 over the Angat Watershed, with PNRI values ranging from 11.18 to 40.40% (Figure 8). On the other hand, during the very strong 2014–2016 El Niño, the following events were detected: a dry spell from October to December 2014, dry conditions in May and June 2015, and a drought event from March to September 2016, with a lowest PNRI record of 35.69%. The longest drought events in the watershed based on PNRI were influenced by the strong 1991–1992 El Niño, namely: the 9-month drought from October 1991 to June 1992, with average (lowest) PNRI record of 40.59% (6.97%), and the 10-month drought event from December 1992 to September 1993, with an average (lowest) PNRI record of 38.73% (2.86%). The lowest record of PNRI was during the moderately strong 2009–2010 El Niño, with 1.18% on February 2010.

Above normal rainfall (PNRI > 120%) was also observed even during El Niño months (Figure 8, red highlight). For example, the presence of heavy rainfall separated each drought and dry event during the recent very strong 2014–2016 El Niño, since the months of April, October, and December of 2015 were detected to have PNRI records of 121.91%, 138.33%, and 149.22% respectively. During these months, the country was hit by tropical cyclones that contributed to the high rainfall records in the watershed, such as the Super Typhoon Koppu (“Lando”) with RX1day of 136.2 mm on October 17, 2015, and the Typhoon Melor (“Nona”) with RX1day of 154.9 mm on December 15, 2015. On the other hand, there were also neutral months (Figure 8, green phases) and La Niña months (Figure 8, blue phases) which have a rainfall record that was way below normal (PNRI < 40%), especially during the summer season. For instance are the 85.5% observed reduction in rainfall (PNRI = 14.5%) in March 2008, classified as a La Niña month, and the 89.3% reduction (PNRI = 10.7%) occurred in April 2004, a neutral month.

This study also identified the past meteorological drought over the Angat Watershed using SPI at various timescales as shown in Figure 9b–d. With an increasing timescale, the range and duration improved but the representation of the dryness magnitude deteriorated. Details of each drought event according to the computed SPI across timescales are summarized in

Table 8. Meteorological drought events in Angat Watershed were identified from 1987 to 2021 based on SPI values of different intervals (3-, 6-, 12-month) as shown in Figure 9.

Start		End	Duration	Magnitude	Intensity	Type	Interval
SPI-3
1988-Jul **	-	1989-Feb **	8	8.60	1.08	Extreme	-
1989-Nov	-	1990-Mar	5	6.33	1.27	Severe	9
1990-Nov	-	1991-May *	7	8.14	1.16	Severe	8
1991-Oct*	-	1992-Jul	10	12.71	1.27	Severe	5
1992-Dec	-	1993-Nov	12	18.16	1.51	Extreme	5
1997-Jan	-	1997-Apr	4	2.26	0.57	Moderate	38
1997-Aug *	-	1999-Jan **	18	18.86	1.05	Extreme	4
2000-Oct **	-	2000-Nov **	2	2.19	1.09	Severe	21
2002-May	-	2002-Jun *	2	1.50	0.75	Mild	18
2004-Apr	-	2004-Jul *	4	5.45	1.36	Extreme	22
2005-Feb *	-	2005-Nov **	10	9.59	0.96	Severe	7
2008-Feb **	-	2008-Apr **	3	2.37	0.79	Severe	27
2009-Dec *	-	2010-Nov **	12	8.67	0.72	Severe	20
2013-Dec	-	2015-Feb *	15	8.68	0.58	Extreme	37
2015-May *	-	2015-Sep *	5	2.53	0.51	Mild	3
2016-Mar *	-	2016-Oct **	8	7.95	0.99	Severe	6
2020-Aug **	-	2020-Sep **	2	2.05	1.03	Moderate	46
2021-May **	-	2021-Jul	3	3.32	1.11	Severe	8
SPI-6
1988-Jul **	-	1989-Feb **	8	10.25	1.28	Extreme	-
1989-Dec	-	1990-May	6	6.05	1.01	Severe	10
1990-Dec	-	1991-Aug *	9	6.96	0.77	Severe	7
1991-Oct*	-	1992-Jul	10	13.07	1.31	Extreme	2
1993-Feb	-	1994-Jan	12	17.88	1.49	Extreme	7
1997-Apr	-	1999-Mar **	24	27.68	1.15	Extreme	39
2004-May	-	2004-Oct *	6	4.90	0.82	Moderate	62
2005-May	-	2006-Jan **	9	10.69	1.19	Extreme	7
2009-Dec*	-	2011-Jan **	14	9.29	0.66	Severe	60
2014-Feb	-	2015-Nov *	22	14.01	0.64	Moderate	37
2016-Jun	-	2016-Dec **	7	9.13	1.30	Severe	7
SPI-12
1988-Aug **	-	1989-Jul	12	10.63	0.89	Moderate	-
1991-May*	-	1992-Oct	18	17.84	0.99	Extreme	23
1993-Apr	-	1994-Jul	16	14.56	0.91	Extreme	7
1997-Oct*	-	1999-Jul **	22	37.08	1.69	Extreme	40
2005-Jul	-	2006-Jul	13	9.39	0.72	Severe	73
2010-May	-	2011-Apr **	12	7.97	0.66	Severe	47
2014-Aug	-	2017-Apr	33	23.30	0.71	Moderate	41

* El Niño month, ** La Niña month; Note: (1) Duration is in months. (2) Type is based on the highest SPI value reach per event. (3) Interval refers to the time gap in months before the end of the last drought event.

. Based on SPI-3 (Figure 9b), 18 meteorological drought events in the Angat Watershed were identified, which lasted from 2 up to 18 months. There are 16-month intervals on average between drought events based on SPI-3, although they vary from 3 up to 46 months. Based on SPI-6 Figure 9c), 11 meteorological drought events in the Angat Watershed were identified which lasted from 6 up to 24 months. Intervals between SPI-6 drought events vary from 2 up to 62 months.

A thorough investigation of SPI-12 was conducted, since this timescale is normally tied to streamflow and reservoir levels [57] and commonly used in the assessment of the hydrologic impacts of meteorological droughts [39]. Based on SPI-12 (Figure 9d), seven meteorological drought events in the Angat Watershed were identified which lasted from 12 up to 33 months. There are 39-month intervals on average between drought events, although they vary from 7 up to 73 months. Based on SPI-12, three extreme drought events were identified, which all happened before the year 2000: from May 1991 to October 1992, with a drought intensity of 0.99; from April 1993 to July 1994, with a drought intensity of 0.91; and from October 1997 to July 1999, with a drought intensity of 1.69.

The effect of each El Niño event on the meteorological drought based on SPI-12 over the Angat Watershed is unique in magnitude, intensity, and duration. The onset of the strong 1991–1992 El Niño, and the drought over the watershed, were detected simultaneously in the same month (May 1991). The drought event ended four months later (October 1992) than the end of the warm episode in the Pacific (June 1992). It is followed by the 16-month extreme drought event which happened during the neutral months of April 1993 to July 1994. This event developed ten months after the termination of the strong 1991–1992 El Niño. In October 1997, a 22-month drought event developed, which is five months later than the onset of the very strong 1997–1998 El Niño. Extreme dry months with SPI less than −2 were observed for five consecutive months from May to October 1998, despite the fact that there was already a shift from warm to cold episodes in the Pacific during that time.

From July 2005 to July 2006, there was a severe meteorological drought based on SPI-12 over the watershed, despite the occurrence of a cold episode in tropic Pacific from November 2005 to March 2006 (ONI < −0.6). This drought event developed five months after the end of the weak 2004–2005 El Niño. A severe meteorological drought event over the watershed developed a month after the termination of the moderately strong 2009–2010 El Niño. This drought lasted over a year until April 2011, coinciding with the strong 2010–2011 La Niña. The longest meteorological drought in the watershed based on SPI-12 was recorded from August 2014 to April 2017. This drought event developed earlier (by three months) and terminated much later (by twelve months) relative to the start and end of the very strong 2014–2015 El Niño, respectively. No drought events, according to the classification of SPI-12, were identified during the moderately strong 1994–1995 and 2002–2003 El Niño. It is the same with the PNRI classification, where only a dry spell was recorded from January to April 1995, with lowest PNRI values of 27.24%, while a dry condition from May to June 2002 provided lowest PNRI records of 36.6%

3.4. Correlation Analysis of Extreme Rainfall Indices and Drought Indices to ENSO

Spearman correlation analysis was conducted to assess the relationship of the large-scale atmospheric and ocean circulation patterns of ENSO to the local characteristics of seasonal extreme rainfall and meteorological drought over the Angat Watershed from 1987 to 2021. The correlation coefficients between the ONI 3.4 and the extreme and drought indices are shown in

Table 9. Correlation between ONI 3.4 and the seasonal extreme rainfall and drought indices over Angat Watershed for the period 1987–2021.

Index	Season
	Amihan Monsoon	Summer	Habagat Monsoon	Monsoon Transition
Extreme Rainfall Index
PCPTOT	−0.42 *	−0.52 *	0.11	−0.46 *
RX1day	0.02	−0.13	0.17	−0.25
RX5day	−0.07	−0.23	0.09	−0.27
CWD	−0.36 *	−0.51 *	−0.08	−0.23
CDD	0.13	0.36 *	−0.14	0.51 *
R10mm	−0.44 *	−0.4 *	0.25	−0.59 *
R20mm	−0.56 *	−0.4 *	0.14	−0.46 *
R75mm	−0.43 *	−0.16	0.12	−0.56 *
Meteorological Drought Index
PNRI	−0.45 *	−0.55 *	0.09	−0.45 *
%DEV	−0.45 *	−0.55 *	0.09	−0.45 *
SPI-3	−0.55 *	−0.62 *	0.04	−0.24
SPI-6	−0.56 *	−0.69 *	−0.05	−0.08
SPI-12	−0.40 *	−0.53 *	−0.06	−0.10

* Significant at α = 5%.

.

ONI 3.4 had a negative correlation with the rainfall intensity (PCPTOT, RX1day. RX5day) and frequency (R10mm, R20mm, R75mm) over the Angat Watershed for all seasons except the Habagat season, where a positive correlation was observed. These negative correlations were statistically significant (α = 5%) for PRCTOT, R10mm, R20mm, and R75mm. The CWD which has negative correlation with ONI 3.4 across seasons is found to be statistically significant (α = 5%) during the Amihan season (−0.36) and summer season (−0.51). The CDD which has a positive correlation with ONI 3.4 across seasons except for the Habagat season are found to be significant during the summer (0.36) and monsoon transition season (0.51). A significant inverse relationship between the ENSO and the PNRI and SPI-3 drought index was detected to all seasons except during Habagat season. A negative correlation was observed for both SPI-6 and SPI-12 across seasons, which is found to be statistically significant (α = 5%) for the Amihan season (−0.56 and −0.4) and summer season (−0.69 and −0.53).

4. Discussion

This study conducted a watershed-scale rainfall variability analysis, considering the water resource system of Angat as a case study, given its economic worth in the Philippines for water supply, power generation, and flood control. Most of the past studies in the country related to spatiotemporal rainfall variability were done on a nationwide [23,24,25] or regional scale [28,29], utilizing the homogenous ground-based rainfall records obtained from the DOST-PAGASA [63] or rainfall derived from reanalysis or satellite products [64,65]. An homogeneity test to detect break/s in a climate record is essential in climatological research, particularly when data are to be used for validating climate models or satellite estimates, for long-term climate analysis such as quantifying trends and establishing climate normal, as discussed by the WMO [48,54,57]. Results of homogeneity tests of this study for rainfall records of NPC- and MWSS-monitored stations in the Angat Watershed showed that most of the data series proved to have a high degree of homogeneity. However, the February and Amihan seasonal rainfall series of Matulid station were found to have a significant break period (α = 5%) in 1999 and 2000, respectively. The results of the homogeneity test suggest that, for trend analysis (and for climate-related studies), one should consider the detected breakpoint to ensure that the only variations to be detected are caused by climate and not by non-climate factors [24,48]. Also, the use of long-term and reliable rainfall data is a prerequisite for the accurate identification of the drought tendency of any region [66].

The Angat Watershed inflow to the reservoir has a strong sensitivity to rainfall variability [10,19]. The local rainfall in the watershed shows the presence of seasonality with short drier season and moderately distributed total annual rainfall. The variations of rainfall likewise have implications for the water availability of the Angat Reservoir. The observed increasing rainfall intensity (PCPTOT, RX1day, RX5day) and frequency (R10mm, R20mm, R75mm), if it continues in the future, could have an implication both for the water resources operation to satisfy the multiple objectives of the Angat Reservoir and for the flood operation that prevents damage in the downstream areas. Generally, rainfall within the watershed has an increasing annual PCPTOT, with increasing frequency in terms of R10mm and R20mm. All NPC-monitored rainfall stations have an increasing trend of PCPTOT, ranging from 3.87 to 31.83 mm/year over the 35-year study period (1987–2021), while Umiray station, with a 20-year temporal record (2001–2021), has an increasing rate of 40.85 mm/year. On average, the watershed rainfall increases at a rate of 18.95 mm/year. An increase in the annual frequency of heavy rainfall is detected at all stations, with a significant annual upward trend (α = 5%) particularly for the R10mm (R20mm) of Maputi station with a rate of 8.6 days/decade (5.8 days/decade) and Umiray station with 14.4 days/decade (10 days/decade). Specifically, rainfall during the Amihan season shows an increasing trend both in terms of intensity and frequency. A significant (α = 5%) increasing trend in the Amihan rainfall extremes was observed in all stations, i.e., PCPTOT at a rate ranging from 7 to 65 mm/year and R10mm at a rate ranging from 1.7 to 10.09 days/decade. For R75mm, there was no trend observed for all NPC-monitored stations, except for the significant increasing trend (α = 5%) in Maputi station during the monsoon transition season (SON). Overall, the study suggests that spatial incoherency of the trend direction of extreme rainfall across seasons may also be observed even for a watershed-scale analysis.

The significantly increasing trend of PCPTOT in the eastern side of the Angat watershed, across the Sierra Madre Mountains Range, means additional river flow that could be diverted for Angat Reservoir through UATP, which has a total inflow capacity of 21 CMS and current annual average diversion of 11.6 CMS. There is also an ongoing Sumag River Diversion project tapping this side of the watershed [67]. Studies to determine the relationship of the watershed rainfall to the reservoir inflow by simulating available physically based hydrologic models or data-driven statistical models could be implemented to further quantify how much rainfall would have a substantial benefit and risk for reservoir storage.

These findings on the annual rainfall trend are consistent with the regional studies conducted in the past. According to the Fifth Assessment Report of the IPCC [34], an increase of 22 mm/decade of annual PCPTOT was reported over Southeast Asia (SEA), but trends vastly differ across the region and seasons. Endo et al. [27] who also used the MK trend test on the extreme values collected from the dense network of rainfall stations in SEA from the 1950s to 2000s, reported an increasing trend in the annual PCPTOT and R50mm for most stations considered in Luzon Island in the Philippines. A recent study of Cheong et al. [42], who analyzed the rainfall datasets of 146 stations in SEA from 1972 to 2010, revealed a statistically significant increase (α = 5%) in the annual PCPTOT of 59.6 mm/decade in the region. Cheong et al. [30] also observed an upward trend in the regional average annual RX1day (1.6 mm/decade), RX5day (5mm/day), and R20mm (1 day/decade). Another study by Salvacion et al. [26] who also applied the MK trend test with the SS estimator on the rainfall gridded data of the Climate Research Unit (CRU) from 1951 to 2015, reported an average monthly increase of 0.34 mm/year in the country. Cinco et al. [23] who analyzed the long-term trend of extreme rainfall from 34 DOST-PAGASA stations in the country from 1951 to 2010, found an increasing trend in the annual rainfall intensity and frequency (events greater than or equal to the 99th percentile each year) in almost all stations (31/34), which is statistically significant (α = 5%) for 6 stations including the Laoag, Baguio, and Infanta stations in Luzon Island.

Few studies considered the seasonal rainfall trend in the country. According to Cruz et al. [18], who conducted a thorough analysis of Habagat monsoon rainfall from 1961 to 2000, a decreasing trend in PCPTOT (0.026% to 0.075% per decade) was observed for six out of nine (6/9) considered stations in the country. Likewise, the intensity of Habagat rainfall over the Angat Watershed has a decreasing trend for 2/5 stations for PCPTOT, 3/5 stations for RX1day, and 4/5 stations for RX5day. Villafuerte et al. [10], who analyzed the trends of the seasonal rainfall of 35 stations from 1951 to 2010 using the seven ETCCDI extreme indices, suggested a drier condition from January to March (JFM) in the country as indicated by 4/35 stations showing a significant (α = 10%) decrease in PCPTOT. The same study, however, showed that 9/14 stations in Luzon Island have an increasing but statistically insignificant trend (α = 10%) in PCPTOT in the JFM PCPTOT. Villafuerte et al. [24] also revealed that the rainfall trends obtained in the shorter term could either consistently represent the continuous long-term trends or denote interdecadal variability.

Three drought indices using the PNRI, %DEV, and SPI were employed in this study to describe the past meteorological drought events in the watershed. PNRI and %DEV both utilize the climatological monthly normal over the watershed. PNRI has been used extensively by DOST-PAGASA to indicate drought and dry conditions in the country through the monthly climate assessment and issuance of rainfall forecasts [68]. Recently, the DOST-PAGASA launched the Southeast Asia Climate Monitoring (SEACM) Project [69] that incorporates the real-time monthly monitoring of SPI-1 utilizing satellite rainfall product. Meteorological drought over the Angat Watershed was observed to be unique in magnitude, intensity, and duration, and highly dependent on the timescale of the index considered. From 1987 to 2021, the number of meteorological drought events detected over the watershed based on SPI-3, SPI-6, and SPI-12 were 18, 11, and 7, respectively. It indicates that, with an increasing timescale, the range and duration improved but the representation of the dryness magnitude deteriorated. The same findings were observed in the study of Valete et al. [70] who examined the historical SPI in the country using the Tropical Rainfall Measuring Mission (TRMM) dataset and found that, from 1998 to 2019, there was a total of six and three drought events corresponding to SPI-3 and SPI-12, respectively.

The drought events over the Angat Watershed were investigated in detail from the results of SPI-12. Several studies [71,72,73] found that the reservoir system has a higher response and coherence with SPI at a higher timescale. SPI-12 is a commonly used timescale in the assessment of the hydrologic impacts of meteorological droughts [39] and is mostly associated with streamflow and reservoir level variation [57]. Based on SPI-12, there were seven meteorological drought events in the Angat Watershed from 1987 to 2021 which lasted from 12 up to 33 months. On average, there are 39-month intervals between these drought events, although they vary from 7 up to 73 months. The onset of the SPI-12 drought event relative to the timing of El Niño showed no distinct pattern. The drought events in the watershed could begin either earlier or later than the start of warm episodes in the Pacific, like what happened during the 2014–2015 and 1997–1998 El Niño events, respectively. A study by de los Reyes and David [29], who use %DEV to describe the influence of El Niño events in the rainfall anomaly from 1971 to 2000 in 73 rainfall stations in the country, revealed that the reduction in rainfall (%DEV < −50%) in the country was initially detected 3 to 5 months after the development of a warm episode in the tropical Pacific, with the recovery usually starting in Mindanao when the tropical Pacific returns to normal. The study of Jaranilla-Sanchez et al. [19] showed that there was a 2 to 7 months lag time between the rainfall anomaly over the Pampanga River Basin before this deficit occurred in the groundwater level.

The results of correlation analysis to assess the relationship of ENSO to the extreme rainfall and meteorological drought over the Angat Watershed were in agreement with the earlier findings of Lyon et al. [13,14] and Villafuerte et al. [24] on the seasonal reversal of ENSO in rainfall variability in the country. An El Niño (La Niña) event is normally associated with drier (wetter) than normal conditions which could cause severe events of drought (floods). However, Lyon et al. [13,14] found that the rainfall response in the country to ENSO during Habagat season tends to be above (below) the long-term normal rainfall during an El Niño (La Niña) event, while an exactly inverse observation occurs monsoon transition and Amihan Monsoon season. As the seasons approach the ENSO mature phase which is usually during the Amihan Monsoon season, Villafuerte et al. [24] suggested a statistically significant drier (wetter) condition over the Philippines was expected during El Niño (La Niña). A recent study by Liao et al. [74] concluded that, during the mature phase of El Niño throughout the Amihan Season, particularly in January and February, the rainfall in the country started to significantly decrease. Likewise, ONI 3.4 had a negative correlation with the rainfall intensity (PCPTOT, RX1day, RX5day) and frequency (R10mm, R20mm, R75mm) over the Angat Watershed for all seasons except the Habagat season, where a positive correlation was observed. A significant inverse relationship between the ONI 3.4 drought indices (PNRI and SPI-3) was also found across seasons except during the Habagat season. The results of the correlation analysis also showed that ONI 3.4 is best correlated with the local SPI-6 of the watershed during Amihan and summer seasons, indicating a 6-month lag between the two indices.

5. Conclusions

Efforts to analyze rainfall in the Angat Watershed are essential for effective water resource management and flood control of the Angat Dam. This study utilized the historical rainfall data from rainfall monitoring stations for the Angat watershed. Trend analysis was performed on the homogenous rainfall series to determine the rate and direction of change in extreme rainfall characteristics in the study area. These indices describe the three major attributes of extreme events, namely: intensity (PCPTOT, RX1day, and RX5day), duration (CDD and CWD), and frequency (R10mm, R20mm, and R75mm). Rainfall in the Angat Watershed generally has an increasing annual PCPTOT: all NPC-monitored stations at a rate ranging from 4 to 32 mm/year, and the MWSS-monitored Umiray station at a rate of 41 mm/year. Specifically, rainfall during the Amihan season shows an increasing trend both in terms of intensity and frequency. A significant (α = 5%) increasing trend in the Amihan rainfall extremes was observed in all stations, i.e., PCPTOT at a rate ranging from 7 to 65 mm/year and R10mm at a rate ranging from 1.7 to 10.09 days/decade. The effect of each El Niño event on the rainfall is unique in magnitude, intensity, and duration. The seasonal reversal of the ENSO on the extreme rainfall and meteorological drought signals in Angat Watershed was also evident. The identified drought events in the watershed based on SPI-12 persisted up to 12–33 months, and could deviate more than 60% from the normal rainfall (PNRI < 40%) during its peak. These insights on the spatio-temporal distribution of local rainfall over Angat Watershed could be handy for the operators and key-actors of Angat Dam for effective flood supervision and water allocation. Understanding these rainfall patterns and their association with El Niño events is pivotal in preparing for and mitigating the impact of future climate variations on water resource availability and flood risk management in the watershed.