Humpback whales in Banderas Bay, Mexico: relative abundance and temporal patterns between 2004 and 2017

In 1986, after whaling nearly drove their populations to extinction, humpback whales (Megaptera novaeangliae) became globally protected by a ban on commercial whale hunting established by the International Whaling Commission (IWC). The abundance of whale populations began to recuperate (Clapham, 2016) as a result. Between 2004 to 2006, project “SPLASH: Structure of Populations, Levels of Abundance and Status of Humpback Whales” (Calambokidis et al., 2008) assessed several parameters of the North Pacific population of humpback whale. SPLASH results reported a yearly population increase of 6.8% and an estimate of 21,808 whales in the whole North Pacific, a number that surpassed pre-whaling population estimates (Barlow et al., 2011). The findings of Calambokidis et al. (2008) from the SPLASH project suggested that the Mexican Pacific was the second most important breeding region for humpback whales of the North Pacific (after the Hawaiian Archipelago), receiving about 40% of wintering humpback whale annually. Medrano et al. (2007) reported that humpback whale population abundance in Banderas Bay was 15 times greater than in the surrounding waters of the Eastern Tropical Pacific and the entrance to the Gulf of California. Later, Martinez-Aguilar (2011) assessed the population abundance of humpback whale in wintering areas of the North Pacific, particularly in Mexican waters (i.e. Baja California, mainland and Revillagigedo) and reported a higher recovery rate for humpback whales wintering in Mexico compared to those wintering in Hawaii. Humpback whale population increases are directly related to oceanographic conditions such as thermal fronts (impacted by climate change) and to El Niño-Southern Oscillations (ENSO) in high latitude summer feeding areas that either directly or indirectly affect prey and habitat availability and therefore, reproductive success (Salvadeo, et al., 2011; Guidino et al., 2014; Cartwright et al., 2019). North Pacific humpback whales winter in latitudes near 20°N, in shallow areas with warm waters ranging between Abstract


Introduction
In 1986, after whaling nearly drove their populations to extinction, humpback whales (Megaptera novaeangliae) became globally protected by a ban on commercial whale hunting established by the International Whaling Commission (IWC). The abundance of whale populations began to recuperate (Clapham, 2016) as a result. Between 2004 to 2006, project "SPLASH: Structure of Populations, Levels of Abundance and Status of Humpback Whales" (Calambokidis et al., 2008) assessed several parameters of the North Pacific population of humpback whale. SPLASH results reported a yearly population increase of 6.8% and an estimate of 21,808 whales in the whole North Pacific, a number that surpassed pre-whaling population estimates (Barlow et al., 2011). The findings of Calambokidis et al. (2008) from the SPLASH project suggested that the Mexican Pacific was the second most important breeding region for humpback whales of the North Pacific (after the Hawaiian Archipelago), receiving about 40% of wintering humpback whale annually. Medrano et al. (2007) reported that humpback whale population abundance in Banderas Bay was 15 times greater than in the surrounding waters of the Eastern Tropical Pacific and the entrance to the Gulf of California. Later, Martinez-Aguilar (2011) assessed the population abundance of humpback whale in wintering areas of the North Pacific, particularly in Mexican waters (i.e. Baja California, mainland and Revillagigedo) and reported a higher recovery rate for humpback whales wintering in Mexico compared to those wintering in Hawaii. Humpback whale population increases are directly related to oceanographic conditions such as thermal fronts (impacted by climate change) and to El Niño-Southern Oscillations (ENSO) in high latitude summer feeding areas that either directly or indirectly affect prey and habitat availability and therefore, reproductive success (Salvadeo, et al., 2011;Guidino et al., 2014;Cartwright et al., 2019).
North Pacific humpback whales winter in latitudes near 20°N, in shallow areas with warm waters ranging between Abstract Banderas Bay, Mexico is an important breeding and transit area for the North Pacific humpback whale (Megaptera novaeangliae) population. In this paper we estimated relative abundance (RA = number of whales/hours of navigation) as a proxy to assess population temporary patterns in the area. We analyzed data from 14 breeding seasons (2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013)(2014)(2015)(2016)(2017), collected between December and March each winter. A total of 8,013 whales were observed in 1,394.6 navigation hours. Average seasonal RA was 5.7 whales per hour with a maximum of 7.5 (2013) and a minimum of 4.0 (2016). Sea surface temperature (SST) averaged 25.1°C and remained within the range considered optimal for humpback whale reproduction areas. SST showed no significant correlation with RA (r = 0.183). Inter-seasonal RA values suggested an increase throughout the study period, although the increase was not statistically significant (R = 0.32; R 2 = 0.10; t = 1.15, p > 0.05). Intra-seasonal analyses showed that RA in December and January were significantly higher (U = 150, p < 0.05) than in February and March; this pattern was consistent throughout the seasons of study. These results represent a shift in the intra-seasonal abundance peak relative to previous studies when most whales were observed between January and February. It is important to recognize changes in 21.1 and 28.3°C (Rasmussen et al., 2007;Guidino et al., 2014). Banderas Bay (Jalisco-Nayarit) is considered to be part of the mainland coast of Mexico wintering aggregation (Rice, 1978;Urbán and Aguayo, 1987) and because of its bathymetric and oceanographic characteristics, this region constitutes an ideal spot for humpback whale reproductive activities 1 . Humpback whale photoidentification studies show that migrations between winter regions and summer areas do not follow a simple pattern. However, there is clear evidence that most of the whales from mainland Mexico prefer California, Oregon, Washington and southern British Columbia as their feeding grounds, whilst a small percentage travels to southeast Alaska (Urbán and Aguayo, 1987;Calambokidis et al., 2000;Urbán et al., 2000).
Humpback whales are affected by anthropogenic activities because of their coastal distribution. Boat collisions, entanglements, acoustic pollution, oil spills and prey reduction are their main identified threats globally (Fisheries and Oceans Canada, 2013;Clapham, 2016;SEMARNAT, 2018). Estimates of cetacean abundance, biomass, and population density are key to assessing the potential effects of anthropogenic perturbations on cetacean populations (Carretta et al., 2006). Moreover, identifying and understanding spatial and temporal patterns of behavior contributes to the prediction and mitigation of emerging threats (Ingman et al., 2021). Therefore, it is a priority to continue to evaluate population trends in sites that constitute important aggregation areas, such as Banderas Bay. In this study, new estimates of relative regional abundance were calculated in order to determine if there have been perceivable changes within Banderas Bay.
We hypothesized that the number of humpback whales sighted per time unit has increased between seasons, since it is the overall trend in the North Pacific population.
To assess the relative abundance of whales observed per hour, we used data collected from whale watching boat tours. The data was coupled with inter-and intra-seasonal analyses of variation during 2004-2017 humpback whale breeding seasons in Banderas Bay, Mexico.
It is an open bay which includes the municipalities of Cabo Corrientes, Puerto Vallarta and Banderas Bay, in the states of Nayarit and Jalisco. The bay has its limits to the north at Punta Mita, Nayarit and to the south at Cabo Corrientes, Jalisco ( Fig. 1). It extends 33 km from north to south and 40 km from season is hereafter referred to by the year that includes most of its duration (i.e., the 2004 season begins in December 2003 and continues through January, February and March 2004). For data management purposes, each season was divided into four months following the Gregorian calendar. All sightings were performed in accordance with Mexican whale watching regulation, NOM-131-SEMARNAT-2010. A single Sea Surface Temperature (SST) value (in Celsius degrees) for Banderas Bay per month was obtained from the NOAA NCEI Extended Reconstructed Sea Surface Temperature (NCEI, 2019).

Data analysis
The relative abundance (RA, whales/h) was calculated dividing the total number of whales observed by the search effort, which was defined as the time spent transiting between observations, excluding the time spent observing whales. Search effort was calculated by subtracting the sighting duration, determined from the beginning and end of the encounter, from the total duration of the trip.
The data did not show normal distribution, nor variance homogeneity. Therefore, to explore the population tendency, we ran a logarithmic regression where RA was the dependent variable and time (season) the independent variable. Regression results were tested for significance with F-Fisher statistics (n=14). Inter-and intra-seasonal variations in RA were analyzed using the Kruskal-Wallis test (H) with p < 0.05, n = 14; followed by a Mann-Whitney U test (p < 0.05) when significant differences were found. To evaluate the relationship between RA and SST, a Spearman Correlation analysis was conducted, with monthly values. All statistical analyses were done using the software Minitab ver. 19 (Minitab, Inc, 2019).

Results
Analyses were based on data collected during 887 whale watching tours. Each trip had a 3-hour total duration, however, navigation effort varied based on the number of sightings. On average, navigation effort was 1.6 navigation hour per trip (s.d. east to west covering an area of approximately 1,000 km 2 (Plata and Filonov, 2007). Oceanographically, it is a transition area which also has unique geological and ecological characteristics. The northern area of the bay has a wide continental shelf; the gradual slope has an average value of 0.012 with a depth of 100 m that is reached 8 km from shore (Plata et al., 2006). At the southern coast, the average slope value is 0.080 and the 100 m isobath is 1.5 km from shore (Plata et al., 2006). There is a submarine fossa in the southern part of the bay that reaches a maximum depth of 1,436 m (INEGI-SPP, 1983) located towards the center of the bay, 8 km from the south coast (Plata and Filonov, 2007).

Data collection
In order to assess humpback whale relative abundance, we used data collected by trained biologists working for a local tour company (Ecotours de México) during every whale watching season from 2004 to 2017. Navigations were conducted onboard an inflatable 7 m zodiac-type boat with two 90 HP outboard engines. Launching points for the boats were Punta Mita at the north end of the bay and Puerto Vallarta in the center of the bay (Fig. 1). The direction for each navigation was random, using the method described by Urbán (1983) as "cetacean observations at sea". In the aforementioned method two observers search for whales scanning the horizon 360° with their naked eyes. Effort was mostly directed to the northern part of the bay since previous studies have shown that this is the region with the highest humpback whale concentration 2 (Ladrón de Guevara-Porras, 1995, Espinoza-Rodríguez, 2009). The collected data included the trip schedule: date, port of departure, times of departure and return to port, and time of sightings, including initial and final time of  ±0.6) and the remaining time was dedicated to whale observation. A total of 8,013 whales (7,324 adults and 689 calves) were observed in 3,289 sightings over 1,394.6 navigation hours (Table 1). Overall, an average seasonal RA of 5.7 (s.d. ± 1.0), whales per hour was calculated for the whole study period, with a seasonal maximum of 7.5 and a minimum of 4.0 whales per hour (Table 1). The logarithmic regression showed no significant tendency for the seasonal RA (R 2 = 0.10; F = 1.36, p > 0.05; Fig. 2), nor the seasonal comparison (H = 10.5, p > 0.05). Monthly relative abundance varied throughout the study period with a maximum of 9.5 and a minimum of 1.2 whales per hour (Table 2). Significant differences were found in the intra-seasonal RA analysis (H = 28.31, p < 0.05). RA values were significantly higher for the first two months of the season across all seasons of study. Intra-seasonal comparisons showed that December and January had similar RA values (U = 150, p < 0.05) which differed from RA values during February and March (Fig. 3).
December and January consistently had the highest SST during the study period averaging 26.4°C and 25.1°C, respectively, while February and March had the lowest with 24.5°C and 24.3°C ( Table  2). The Spearman Correlation Analysis did not show a significant relation between the RA and the SST (r = 0.183; Table 1, Fig. 4).

Discussion
Seasonal relative abundance Large, migratory predators are considered sentinel species for ecosystem processes and climate-related changes; therefore, their utility as indicators is dependent upon an understanding of their response to environmental variability (Fleming et al., 2016) underlying the importance of generating comprehensive and ongoing knowledge about the environmental factors influencing abundance in breeding grounds. During all 14 sampling seasons, oceanographic conditions in Banderas Bay were relatively stable; SST remained within the 21.1-28.3°C interval which is considered optimal for humpback whale reproduction areas ( Fig. 4; Rasmussen et al., 2007). Therefore, based on our correlation analysis results, we suggest that the presence of humpback whales in Banderas Bay during the study period was not influenced by SST, as reported by Urbán et al. 1 for the 2004-2007 period in the same breeding ground.
Overall, there was a gradual increase in RA across the 14 seasons of this study (Fig. 2). Even though this increase was not statistically significant it does seem to reflect the general humpback whale increase in the North Pacific population as reported by several authors 2 (Calambokidis et al., 2008;Barlow et al., 2011;Martínez-Aguilar, 2011). Furthermore, from 1989to 1999Ladrón de Guevara-Porras (2001 reported a maximum seasonal RA of 2.7 whales per hour for Banderas Bay (with a similar methodology -searching randomly for whales with outboard motor vessels) which is below all seasonal RA estimated in the present study, including the lowest (4.0 whales per hour).
Fluctuations in breeding grounds have been reported as a response to changes of environmental conditions in feeding grounds (Morete et al., 2008;Ávila et al. 2020), i.e., whales are affected directly or indirectly by oceanographic alterations and prey availability (Cartwright et al., 2019;Seyboth et al., 2021). While RA was not significantly correlated with SST in the study area, our results suggest a relation with SST changes in the feeding grounds. In 2013, we recorded the highest RA in Banderas Bay which coincides with an increase in abundance documented (with a different methodology) for the Hawaiian subpopulation of humpback whale (Cartwright et al., 2019), which also belongs to the North Pacific population. This could be a reflection of the 2008 to 2018 period when various large-scale climatic phenomena produced changes in the conditions of the northern Pacific region (Fleming et al., 2016;Cartwright et al., 2019). From 2009 to 2010, the Pacific Decadal Oscillation (PDO) began a weak negative phase which intensified from 2011 to the beginning of 2014. At the same time, the ENSO entered a weak negative phase. Both phenomena were associated with colder relative average water temperatures for the North Pacific (Fleming et al., 2016). This caused changes in upwellings and, in turn, increased food availability for whales at high latitudes (Fleming et al., 2016;Cartwright et al., 2019;Seyboth et al., 2021).
Beginning in the summer of 2014, environmental conditions changed in the North Pacific feeding grounds. The PDO entered a high positive phase (from 2014 to 2016) with the highest positive values in the 115-year record (Peterson et al., 2016). This phase had unprecedented impacts on trophic levels of the marine ecosystem, including low primary productivity (Whitney, 2015). Lower food availability, as a result of changing environmental conditions in high latitudes, could have translated into a period  of increased nutritional stress during the following reproductive seasons (Cartwright et al., 2019). Because proper energy stores are a requirement for female whales to maintain reproductive status (Lockyer, 1981;Christiansen et al., 2016), the preceding food scarcity during the summer may have caused ovulating females to become rare towards the end of the following reproductive seasons (2015 and 2016). This would have reduced mating opportunities and caused males to return to their feeding grounds earlier than usual (Cartwright et al., 2019). Essentially, years with low food availability during the feeding seasons were followed by shorter than usual reproductive seasons and lower levels of RA. Potentially this accounts for the drops in RA that we saw during the 2015 and 2016 reproductive seasons in Banderas Bay (Fig. 2).

Temporal patterns
In breeding areas, whale abundance may fluctuate throughout winter because of variable migration patterns among groups, local movements and changes in social behavior (Dawbin, 1966;Herman and Antinoja, 1977;Darling et al., 1983;Mobley and Herman, 1985;Mattila et al., 1994). According to previous reports, humpback whales are present in Banderas Bay from November to April (Ladrón de Guevara-Porras, 1995;Medrano-González et al., 2009). We analyzed monthly variation in the RA during the official whale watching season and found that the RA was consistently high and statistically similar during December and January across the study period with an average of 6.0 and 7.7 whales sighted respectively per hour. Afterwards, a gradual decline in RA was observed in February and March with a respective average of 4.9 and 3.6 whales per hour (Table 2, Fig. 3). Moreover, as the RA in December and January increases, the RA in March decreases (Fig. 3). These results show the peak of abundance in December and January, in contrast to peaks previously reported by Ladrón de Guevara-Porras (1995) for the 1989-1999 period, where the peak was the last two weeks of January and the first two weeks of February. This could be a reflection of changing temporal migration patterns to the feeding grounds. The timing of migration is also influenced largely by basin-scale environmental variables in the northeastern Pacific Ocean such as the ENSO, the PDO and the North Pacific Gyre Oscillation (Ingman et al., 2021). In California, Ingam et al. (2021) found that humpback whales arrived to their feeding grounds on average 120 days earlier in 2016 than they had in 1993, while departures to the breeding grounds showed little to no change. Furthermore, Ávila et al. (2020) found that, as the whale population increases, individuals tend to arrive earlier to their reproduction areas in warm regions. Our findings suggest that North Pacific humpback whale abundance has continued to increase, and as a result, their arrival to the wintering grounds in Banderas Bay occurs earlier. This idea is supported by the results of the present study and can be observed in the new peak in RA which occurs in the earlier months of the season (December and January).
If these changes in temporal patterns prevail, they could indicate notable plasticity for the species in response to global climate change (Ramp et al., 2015) which directly impacts oceanographic conditions and food availability for humpback whales.
In conclusion, this study demonstrates that the population tendencies of the humpback whale which breed and mate in Banderas Bay are changing. Two significant variations in RA were found between 2004 and 2017. The first was a gradual increase throughout the study period, most evident from 2004 to 2014; the observed drop in RA during the 2015 and 2016 years could have been caused by environmental changes during the summer in the feeding grounds. The second variation was a peak in humpback whale RA during the months of December and January which represents a shift compared to results reported for the previous decade.
It is important to consider fluctuations in abundance and temporal patterns of humpback whales for management purposes, since they have direct effects in regions which have developed whale watching and other associated activities for tourism. Changes in relative abundance support the need for long-term studies of absolute abundance using consistent methods (Mobley et al., 1999). Therefore, we recommend continuing to monitor changes in temporal patterns of humpback whales' occurrence as well as to carry out studies that allow estimating the closest to the absolute abundance, through population models in order to better understand if the increasing  trend in relative abundance and the shift in the seasonal abundance peak are maintained over time.
This study highlights the usefulness of data collected by the whale watching industry to conduct scientific research, provided it is collected and managed adequately. This is known as citizen science and it is a scientific and educational alternative which also promotes environmental awareness (Dickinson et al., 2010;Earp and Liconti, 2020) and leads to the conservation of marine ecosystems (Kelly et al., 2019). Generally, the whale-search effort of the whale-watching industry surpasses the effort which any research team can accomplish, and tour boats are ideal platforms to collect data. The information gained from this study will be useful for refining our knowledge on breeding humpback whales in Banderas Bay and will contribute to inform conservation and management strategies and policies in the region.