Landscape Pattern Impacts on the Population Density and Distribution of Black Shama (Copsychus cebuensis Steere) in Argao Watershed Reserve, Argao, Cebu, Philippines

1 Cebu Technological University, Cebu Campus, Argao 6021, Cebu, Philippines 2 Institute of Renewable and Natural Resources, College of Forestry and Natural Resources, University of the Philippines Los Banos, Laguna 4031, Philippines 3 Institute of Biological Sciences, College of Arts and Sciences and School of Environmental Science & Management (SESAM), University of the Philippines Los Banos, Laguna 4031, Philippines 4 School of Environmental Science and Management, University of the Philippines Los Banos, Laguna 4031, Philippines


Introduction
Landscape ecology explains landscape pattern and interprets its possible ecological impacts particularly on energy and material flow on flora and fauna found at the landscape [1].Argao watershed is currently experiencing these various changes that threaten biodiversity.It is located at the southeastern part of Cebu, Philippines.Some of the limited natural forest patches of the watershed are sporadically distributed on upland barangays, which are pervasively threatened from sustained illegal cutting of trees, fuel wood gathering, and converting forest to other land uses such as agriculture [2].
In 2001, the International Union for the Conservation of Nature (IUCN) declared Mt.Lantoy, which is located within Argao watershed reserve (AWR), as 71st of the 117 Important Bird Areas (IBA) in the country [3], because of the presence of endangered Cebu Black Shama (Copsychus cebuensis) and Cebu Flowerpecker (Dicaeum quadricolor).
C. cebuensis is an endemic bird species to the island of Cebu in the Philippines.It has always been considered rare and endangered as it has very small range and population [3].Conducting detailed ecological research to determine threats and their level of impact, as a case study in Argao watershed (AWR), will facilitate conservation planning and

Description of the Study Site
2.1.1.Topography and Slope.The general topography of the area is steep to very steep, with rugged terrain ranging from 12% to 60% slope in any direction [4].The highest elevation is about 1,000 m above sea level which is located in Barangay Ablayan, Dalaguete, Cebu, in the southwestern part of the watershed where headwaters emanate.On the other hand, one of the highest mountain peaks within the watershed area is situated in Mt.Lantoy in Barangay Tabayag, Argao, with an elevation of about 593 m above sea level.Supervised using maximum likelihood (ML), the algorithm was used in image classification.Then ground truthing or ground validation activities between the two simulation processes were conducted.

Ground Truthing or Validation.
High-resolution digital map like Google Earth was used, and actual study site visitation based on randomly GIS generated geographic coordinates was done to validate the images.This was important to clearly identify the areas covered by heavy clouds or shadows produced on the SPOT-5 satellite image.
2.2.3.Supervised Classification.This is a process used for the extraction of quantitative information from remotely sensed image data with the analyst having prior knowledge of the available classes.This process is performed after the completion of the ground truthing activities, where the information obtained from the latter is used as the basis for the supervised classification technique.Finally, the number of pixels and area for each class and the landscape metrics or indices using ArcGIS/Patch analyst extensions [14] were estimated, and the statistics for each class were shown in graphs and/or tables.Figure 2 summarizes the foregoing processes.

Area Distribution Analysis of Different Land Use/Land
Cover (LULC) Classes.To determine the areal distribution of the different LULC classes, the newly generated LULC map was used in the analysis as input.The area of each LULC class was determined using GIS, specifically the spatial analyst tools.Subsequently, the results were tabulated and qualitatively analyzed.

Landscape Pattern Analysis (LPA)
2.3.1.Landscape Metrics (LM) Selection.Table 1 summarizes the landscape metrics that were used in quantifying the landscape pattern (LP) in AWR.Specifically, shape metrics intend to measure landscape configuration (LC), where LC refers to the nature of shape of patches in certain class or the entire watershed landscape.It is also an overall measure of how complex or irregular the shape of all the patches is in a class or landscape [5].

Calculation and Analysis of Landscape Metrics (LM) or
Indices.The newly generated LULC map of AWR was used as base map to generate and calculate LM or shape metrics using PATCH ANALYST extensions for ArcMap Version 4.2 [14].The indices or metrics were calculated only at class level since fragmentation is best analyzed at this level rather than at the landscape [5].Furthermore, the calculation for the percentage occupied by each class with respect to the entire watershed landscape was undertaken with an excel format.Subsequently, qualitative analysis of LP in AWR was undertaken on the basis of calculated metrics or indices value and available literature.Equations used in the calculation of indices or metrics are as follows.
(1) Mean Shape Index (MSI).MSI is equal to 1 when all patches are circular (for polygons) or square (for rasters or grids), and it increases with increasing patch shape irregularity [5].
where   is the perimeter of the th patch and   is the area of the th patch.
(2) Area-Weighted Mean Shape Index (AWMSI).It is equal to 1 when all patches are circular (for polygons) or square (for rasters (grids)), and it increases with increasing patch shape irregularity [5].
where   is thae perimeter of the th patch and   is the area of the th patch. total = is the summed area of all patches (Pn  1/4 1   ).
(3) Mean Patch Fractal Dimension (MPFD).Index valueapproaches one for shapes with simple perimeters and approaches two when shapes are more complex [5].
where   is the perimeter (m) of patch  and   is the area (m 2 ) of patch .

Regression
Analysis.An enter method multiple regression analysis was applied to determine the relationship between the landscape variables (25 predictors were tested in the landscape level analysis and only 10 predictors were finally included in the sampling site) and population density of C. cebuensis using where  is the dependent variable (population density of C. cebuensis);  1 ,  2 , . . .,   are the independent variables or predictors (landscape variables);  0 ,  1 ,  2 , . . .,   are the partial regression coefficients of independent variables or predictors; and  is the Random error.All statistical analyses were undertaken with SPSS 11.5 license to UPLB College, Laguna, Philippines.4 there were three landscape elements or land uses generated, namely, (a) cultivated areas with a total land area of about 3,399 ha or 45%, (b) forest areas with an area of about 3,002 ha or 40%, and (c) built-up areas with 1,050 ha or 15%.It can be observed that AWR is dominated by cultivated areas followed by forest and builtup.This can be attributed to the fact that Argao watershed has higher human population living within and outside the perimeter of the watershed area.According to [16,17] about 24,096 people (from the municipalities of Argao and Dalaguete) reside within the perimeter of AWR and about 108,646 people reside outside the watershed area [18].In fact, Argao municipality was considered as one of the biggest municipalities in southern Cebu comprising 45 barangays, and these barangays are mostly located in the upland areas.
The difference between cultivated and forest areas is only about 375 ha or 5%.This figure may be considered as small fraction, yet this is, however, significant due to an increasing human population in the upland or in the watershed area as evident in the present population level.As population grows demand for watershed resources will likely grow especially areas devoted for cultivation or for agriculture, there might be an increasing encroachment into the remaining forest patches, there will be more land use conversion can be expected that would surely tantamount to further forest deforestation and degradation.Ultimately, this condition may seriously affect and would endanger native forest habitat for C. cebuensis and other natural resources in the watershed area.

Landscape Pattern Analysis (Using Shape Metrics)
(1) Area-Weighted Mean Shape Index (AWMSI).Figure 4 shows the calculated index values in AWMSI for three categories or classes.The results show that cultivated areas have the highest computed AWMSI of 65.45, followed by forest class with 15.49, and the lowest were built-up areas with 5.01.All calculated index values are greater than two which serve as benchmark values for patches having more complex and irregular shape [5].These findings suggest that all three classes or categories have very irregular and/or complex patch shape.These could be attributed to the intense human activities especially in cultivated areas encroaching ultimately into forest areas as observed in this study, thereby exposing the latter to pasture or agriculture, microclimatic alterations that usually occur in the nearby surrounding forest edges through increased penetration of sunlight and wind.
(2) Mean Shape Index (MSI).Figure 5 shows the calculated index values for MSI for all three classes.It can be observed that built-up areas got the highest value of 1.53, while forest and cultivated areas have the same value of 1.39.It can be observed that the values for all the three classes greater than one indicate that the patch shape had increasing irregularity and complexity.This indicates a more fragmented and more heterogeneous watershed landscape [5].Likewise, this can be attributed to the high level of human activities that usually occur along edges and within the forest interior (i.e., kaingin making, infrastructure development, and conversion of forest lands into other uses) which have been observed in the study.
(3) Mean Patch Fractal Dimension (MPFD).Figure 6 presents the computed index value for the three land use classes, of which forest areas had the highest value of 1.85, followed by cultivated lands of 1.84 and built-up areas with 1.73.It was observed that the value for three classes was greater than one and nearly approaches two, suggesting a more complex shape, high patch shape irregularity, and a highly fragmented watershed landscape.These findings further indicate that Argao watershed is more of an edged habitat which can bring impacts on the wildlife species especially forest-interior  species like most of insectivorous birds' especially C. cebuensis [19].
(4) Area-Weighted Mean Patch Fractal Dimension (AWMPFD).Figure 7 presents the results of the calculated index values for AWMPFD.It showed that cultivated areas got the highest index value of 1.65, while forest areas are on the second place of 1.57, and the lowest are built-up areas with 1.53.All the three values were greater than one and nearly close to two.These results suggest that the patch shapes of three classes were complex, further indicating that AWR had exposed edges, fragmented habitat, and patch shape that tend to be more irregular and complex and more heterogeneous.Barangay Tabayag with 81 individuals-ha −1 (Table 4).It can be observed that, though Barangay Canbantug and Barangay Cansuje have the same population density estimates, they had different %CV of 40% and 18%, respectively.It was observed further that the higher the number of contacts, the lower the %CV, whereas the lower the no. of contacts the higher %CV.According to [6], for a distance sampling method to be adapted especially in point count survey, it needs to have at least 80 to 100 individual contacts in order to achieve more reliable estimates of the density of population for certain organism under consideration.Based on these findings it can be inferred that the C. cebuensis are still widely distributed among the fragments of forest habitat patches remaining in Argao watershed.

Relationship between Landscape Variables and Population Density of Copsychus cebuensis (1) Population Density Estimates for Habitat or Landscape Level (Mixed and Natural Forest).
(

3) Regression Analysis between Landscape Variables and Population Density of C. cebuensis for Habitat or Landscape
Level (Mixed and Natural Forest).Table 5 presents the results of regression analysis between the population density of C. cebuensis at habitat or landscape level.The final model shows that adjusted  2 = 0.345, (22, 107) = 4.093,  < 0.000 (using the enter method).By looking at the final model (particularly the adjusted  2 value) this shows that only about 34.5% of the variances of the population density of C. cebuensis were being accounted for by the model.Only three out 25 predictors have been able to explain the behavior of the population density of C. cebuensis, namely, (1) relative humidity (RH) and canopy cover (TCPY) having positive high significant correlation (beta) with population density, 0.400 and 0.261, respectively, with  values at 0.000 and 0.001 with, respectively, which is <0.01 significant level.These findings indicate that with high RH and TCPY the population density of C. cebuensis may increase.While tree basal area (TBA) it has a negative high significant correlation with C. cebuensis population density with coefficient (beta) value of −0.660 and a  value of 0.000 significant at 0.01 level.This finding indicates that if TBA increases, the population of density of C. cebuensis decreases.The study on insectivorous birds and environmental factors in Malaysia's tropical rainforest showed that terrestrial insectivores were sensitive to forest edge and could indicate the quality of forest interior habitats associated with high humidity, dense canopy cover, and deep litter depth [20].

(4) Regression Analysis between Pooled Population Density Estimates of C. cebuensis and Landscape
Variables for Sampling Sites.Based on the multiple linear regression analysis (using enter method) the results have shown a final model of adjusted  2 = 0.212, (10, 119) = 4.474,  < 0.000 with only 10 predictors included in the final analysis.This indicates that the model has accounted only for 21.2% of the variation of population density of C. cebuensis which is lower by 13.3% compared to what is accounted for in the regression done at the habitat level.Table 6 presents the results of the regression done at this level.
It has been observed further that there are four predictors that have explained the behavior of the population density of C. cebuensis, namely, (a) elevation, (b) slope, (c) canopy cover,  6).This means that with increasing elevation the density of population of C. cebuensis will likely increase.According to [21] variations in the number of bird occurrences along elevation gradient in the natural and urbanized habitats in western Himalaya, Nainital district of Uttarakhand, India, found out that bird species richness (BSR) varied considerably with elevation.This difference in BSR among study sites could be due to elevation and vegetation differences associated with elevation and not caused by the presence of a group of mid-altitude specialists.For slope and shrub cover it is negatively correlated with population density of C. cebuensis with coefficients of −0.178 and −0.165 and with values of 0.032 and 0.047, significant at 0.05 level, respectively.These findings indicate that with increasing slope and shrub cover the population density of C. cebuensis will likely decrease.Lastly, for canopy cover it had positive significant correlation with C. cebuensis population density with coefficient of 0.177 and a  value of 0.043 which is significant at 0.01 level.This means that with increasing canopy cover the population density of C. cebuensis will likely increase.

Conclusion and Recommendations.
AWR is dominated by cultivated land, and forest area at AWR is fragmented and heterogeneous.On the basis of the landscape pattern analysis conducted in AWR using four landscape metrics/indices (shape metric) shows that the landscape pattern in AWR could be characterized as more fragmented and more heterogeneous with patches at classes (especially forest areas) that are increasingly complex and irregular.This condition could affect C. cebuensis subpopulations to become more fragmented or isolated.However, results show that C. cebuensis subpopulations are still widespread or widely distributed within these remaining forest patches especially those being sampled.
In regression analysis, there are only three out of 25 and four out of 10 have been able to predict or explain the behavior of the population density of C. cebuensis at habitat (landscape) and sampling site levels.Predictors which are positively highly significant at the habitat level were the following: (a) relative humidity (0.400, 0.01 <  value = 0.000) and (b) canopy cover (0.261, 0.01 <  value = 0.001) while tree basal area is negatively highly correlated (−0.660, 0.01 <  value = 0.000).At sampling site, elevation has high positive significant correlation (0.398, 0.01 <  value = 0.000) with C. cebuensis population density, and canopy cover has positive significant correlation (0.177, 0.05 <  value = 0.043), whereas slope and shrub cover have negative significant correlation (−0.178, 0.05 <  value = 0.032 and −0.165, 0.05 <  value = 0.047).Based on the final model and given the two levels of analysis, at landscape, results have shown that 34.5% of the variances of the population density of C. cebuensis have been accounted for, while at the sampling site it accounts only for 21.2%.
Study suggested that regression analysis at landscape or habitat level is better than same analysis at sampling site.Study concludes with a series of sound policies and management recommendations to preserve AWR.

Policy Measures to Protect C. cebuensis Native Habitats.
Integrating the use of indigenous tree species in rehabilitation and reforestation programs whether it is undertaken by public agencies or private entities within study site to preserve native habitat for C. cebuensis taking into account the remaining forest within the watershed area.Promoting and adopting landscape approach in rehabilitating degraded habitats.The Department of Environment and Natural Resources (DENR) may have been practicing this approach in rehabilitating degraded areas.This approach, however, focused on the relationship at various scales from species up to landscape.Each and every element or patch is taken into account as a basic unit of the entire system.
3.4.2.Management Aspect.Management and rehabilitation of remaining forest habitat patches within Argao watershed should be based the findings of landscape pattern analysis conducted in this study.Results show that the higher landscape metric value indicates that the forest patches in Argao watershed were highly fragmented.Institute connections on isolated patches, mainly through corridor establishment using riparian forest, is an excellent option of patches connection is imperative.There is a need to rehabilitate riparian vegetation especially those that are highly degraded and sparsely stocked with vegetation.Results show that most of the contacts with C. cebuensis were in the valley bottom where relative humidity is high due to the presence of vegetation.

Further Research on Copsychus cebuensis.
Follow-up study on the foraging behavior of the target species to ascertain more on the type of insects being preyed about.Follow-up study on the effect of patch area on the population density and/or abundance of the target species.This study was able to survey selected forest habitat patches within AWR.However, the effect of the patch area with the population density of C. cebuensis has not been explored.
Continue the monitoring populations in all forest patches in Argao watershed to quantify all existing subpopulations.Since this study had covered only limited number of the remaining forest patches within Argao watershed, it is advised to continue monitor all the forest patches within the area to be able to quantify all existing subpopulations.

Figure 2 :
Figure 2: GIS simulation processes for land use/land cover classes generation in AWR.
Figure 3 presents the newly generated land use/land cover of Argao watershed.As shown in Table

Figure 3 :
Figure 3: Newly generated land use/land cover map of AWR from SPOT-5.(Data sources: CENRO, Argao, 2003; CBRMP-LGU, Argao, Cebu, 1999; NAMRIA.) The climate at AWR is classified under Climate Type III.Rainy season is from May up to November, while dry season is from December to April.Figure1shows the climogram of Mactan, Cebu, the nearest weather station to AWR.The figure shows the average monthly temperature and rainfall data from 1973 to 2007.The dotted area represents the dry period.The wet season is represented by the grayshaded area, with rainfall less than 100 mm, and black area, with rainfall more than 100 mm.
[13]2.Climatic Type.Some remnants of natural forest fragments are found at the peak of Mt.Lantoy and at the slopes of Argao River.Isolated forest patches are still present in barangays Canbantug, Panadtaran, Usmad, Tabayag, Conalum, and Cansuje.Mother trees of Ipil (Intsia bijuga (Colebr.)Kuntze)andMolave(VitexparvifloraJuss)arenative to the area.Natural vegetation and plantation forests comprise 29%(1,119.54ha) of the watershed area[12].Mahogany (Swietenia macrophylla) and Yemane (Gmelina arborea) are the two exotic species commonly planted in the plantation forest.2.1.4.Collection/Procurement of Geographic Information System (GIS) Thematic Layers and Remotely Sensed (RS) Image.GIS thematic layers and RS image were obtained to be able to generate the land use and land cover of the study site as well as other watershed maps.These data were obtained from the National Mapping and Resource Information Authority or NAMRIA.NAMRIA provided three digital images like AVNIR, SPOT, and LANDSAT.Among these, the image from SPOT-5 was utilized for this study.SPOT-5 image was with 2.5 m resolution and taken 2006/2007 available during conduct of study.Other base maps were acquired from different agencies including DENR Regional Office-VII, CENRO-Argao, and LGU-Argao particularly CBRM Office.2.1.5.Collecting/Gathering of Other SupplementaryInformation.Related information on the study topic was collected from the internet, especially in the directory of open online journal in the UP system to acquire published articles/materials, books, monograph, and journals, among others which were related to the present work.softwareprogramversion 10[13], particularly in generating the land use/land cover classes (LULCC) of the study site.

Table 2 :
[15] of landscape variables.speciesusingdistancepointcount sampling method[6].A total of 130 circular point count stations, 87 of which were mixed and 43 were natural forest, respectively, with 20 meter radius have been established among the four sampling areas: twenty-six (26) plots in Mt.Lantoy (Brgy.Tabayag), 23 plots in Brgy.Canbantug, 38 in Brgy.Usmad, and 43 plots in Brgy.Cansuje.At each point survey station, the target bird species detected (i.e., either seen or heard) within 20 m radius were recorded for about 10 minutes[15].The distance of every C. cebuensis encountered was estimated when first detected and recorded.Other information also noted such as no. of clusters, time start, time observed, and type of contact, among others.The survey was conducted early in the morning and late in the afternoon (between 0600 H to 1000 H and 1600 H to 1800 H, and no survey was conducted during bad weather especially during rainfall).This period of time is suitable as most of forest birds are active in the early morning and late in the afternoon.All point survey stations were visited just only one time during the whole study period from April 7 to May 28, 2012 (for almost 2-month period).
(b) Dist.from sampling site to a nearest village/community Spatial analyst (GIS) (c) Dist. to agriculture/cultivation Spatial analyst (GIS) (3) Vegetation structures (a) Total height (TH) Bibby et al. [6] H a g a (b) Diameter at breast height (DBH) Bibby et al. [6]/Rosli et al. [7] D i a m e t e r t a p e (c) Canopy density Bibby et al. [6]2.4.Relationship between Landscape Variables and Population Density of C. cebuensis 2.4.1.Landscape Variables (LV).To determine the relationship between landscape variables and population density of ered simultaneously with the bird data, and some LV were gathered indirectly using field data as inputs and equations.Table 2 presents the variables that are considered in this study.2.4.2.Bird Sampling.Data on the population density of C. cebuensis was collected during the breeding season of the target

Table 3
, however, presents the estimated population density of C. cebuensis (using DISTANCE 6 release 2 software program) for two habitat

Table 3 :
Summary table for estimating population density on C. cebuensis.

Table 4 :
Pooled density estimates of population density of C. cebuensis in four sampling locations.
(2) Population Density Estimates at Sampling Sites.Among the sampling locations, Barangay Cabantug and Barangay Cansuje had the highest pooled population density estimate of C. cebuensis with 118 individuals ha −1 followed by Barangay Usmad with 114 individuals-ha −1 , and the lowest is

Table 5 :
Results of regression analysis between population density estimates of C. cebuensis at habitat level (with 25 predictors or landscape variables included).

Table 6 :
Results of the regression analysis between pooled population density estimates in four sampling sites (with 10 predictors included in the analysis).
and (d) shrub cover.Elevation has positive high significant correlation (beta) with population density of C. cebuensis with standardized coefficient of 0.398 and with  value of 0.000, significant at 0.01 level (Table