Variations of phytoplankton chlorophyll in the Bay of Bengal: Impact of climate changes and nutrients from different sources

Phytoplankton biomass, quantified as the concentration of chlorophyll-a (CHL), is the base of the marine food web that supports fisheries production in the Bay of Bengal (BoB). Nutrients from river discharge, the ocean subsurface layer, and the atmosphere have been reported to determine CHL in the BoB. Which source of nutrients mainly determines CHL in different parts of the bay has not been determined. Furthermore, how climate variations influence nutrient inputs from different sources and their impacts on CHL have not been detailed. To address these questions, we used relationships between satellite-derived CHL and in situ river discharge data (a proxy for river-borne nutrients) from 1997 to 2016, physical variables, and modeled dust deposition (DD), a proxy for atmosphere-borne nutrients. Nutrients supplied from the ocean subsurface layer were assessed based on variations in physical parameters (i.e., wind stress curl, sea surface height anomaly, and sea surface temperature). We found that nutrients from the Ganges and Brahmaputra Rivers were important for CHL along the northern coast of the bay. By increasing rainfall and river discharge, La Niña extended high-CHL waters further southward. Nutrients from the ocean subsurface layer determine CHL variations mainly in the southwestern bay. We suggest that the variations in the supply of nutrients from the subsurface layer are related to the generation of mesoscale cyclonic eddies during La Niña, a negative Indian Ocean Dipole, or both. Climate-driven cyclonic eddies together with cyclones can intensify Ekman divergence and synergistically lead to a pronounced increase in CHL in the southwestern bay. Nutrients from the atmosphere mainly determine CHL in the central/eastern BoB. We further suggest that DD in the central/eastern BoB is influenced by ENSO with a 6–7-month time lag. CHL in the central/eastern bay responds to the ENSO 6–7 months after the ENSO peak because of the 6–7-month lag between ENSO and DD. This report provides valuable information needed to plan necessary actions for climate adaptation in local fisheries activities by elucidating how climate variations influence phytoplankton.


Introduction
The Bay of Bengal (BoB), among the largest of the Large Marine Ecosystems (LMEs) in the world, is bounded by the Asian continent in the north, the eastern coast of India in the west, and the Andaman Sea in the east. In the south, the bay is open to the Indian Ocean. The surrounding countries of India, Bangladesh, Sri Lanka, Myanmar, the Maldives, Thailand, Malaysia, and Indonesia have reached a consensus on a strategic action program to solve existing socioeconomic issues (e.g., Elayaperumal et al., 2019) because of the social and economic importance of the bay.
One of the key determinants of the services the BoB provides to nearby societies and economies is phytoplankton primary production, a common metric of which is phytoplankton biomass or chlorophyll-a concentration (CHL, mg m −3 ). Primary production is the production of organic matter at the base of the marine food web and supports the secondary production of biomass by marine organisms at higher trophic levels. Previous studies have reported a strong relationship between CHL, or primary production, and fisheries production in a wide variety of marine ecosystems in coastal waters and the open ocean (Nixon and Thomas, 2001;Nixon and Buckley, 2002). The production of the tropical Hilsa fishery, which contributes to the society and economy of India, Bangladesh, and Myanmar, is controlled by primary production in the BoB (Hossain et al., 2020).
Light, in addition to nutrients, is an important determinant of primary production in the BoB, especially in the northern part of the bay because of the high turbidity associated with high sediment inputs from the Ganges and Brahmaputra Rivers (GBR) . Because of the large influx of freshwater from the GBR, the water column in the northern BoB is strongly stratified. This stratification inhibits the influx of inorganic nutrients from the ocean subsurface layer to the euphotic zone, especially in the central part of the bay and during boreal summer (Kumar et al., 2002). In the part of the bay where light is not a limiting factor, the input of inorganic nutrients is thought to be the major determinant of primary production in the BoB.
There are at least three major sources of inorganic nutrients that support primary production in the BoB. Previous studies have mentioned how inputs of inorganic nutrients from rivers and ocean subsurface waters cause seasonal variations in CHL in the BoB (e.g., Gomes et al., 2000;Levy et al., 2007;Vinayachandran, 2009;Siswanto et al., 2022). The GBR, which globally ranks fourth with respect to freshwater discharge into the sea (an average of 1,032 km 3 year −1 ) (Dai and Trenberth, 2002), supplies nutrients that elevate CHL along the northern coast of the bay during the boreal summer (e.g., Vinayachandran, 2009;Kay et al., 2018). Inputs of nutrients associated with wind-driven upwelling and vertical mixing have been reported to elevate CHL and cause a strong seasonal cycle of CHL in mainly the western and/or southwestern bay (e.g., Shetye et al., 1991;Gomes et al., 2000;Vinayachandran and Mathew, 2003;Martin and Shaji, 2015;Xu et al., 2021). During the boreal summer, CHL variations along the western coast of the bay are influenced mainly by coastal upwelling due to local longshore wind stress (Shetye et al., 1991;Gomes et al., 2000). CHL enhancement during the fall is associated with the westward propagation of an upwelling Rossby wave (a reflection of an equatorial, upwelling Kelvin wave from the eastern boundary of the BoB) along the coast of the bay (e.g., Rao et al., 2010;Sreenivas et al., 2012;Gulakaram et al., 2018). During the winter, CHL concentrations are highest, and nutrients are entrained by strong, northeasterly, wind-driven Ekman pumping and/or vertical mixing (e.g., Xu et al., 2021).
In addition to nutrient inputs from rivers and ocean subsurface layers, atmospheric dust deposition (DD, mg m −2 s −1 ) is also a potential source of inorganic nutrients that can fertilize the bay. Grand et al. (2015) have pointed out that the concentration of nutrients in the surface waters of the BoB is highly consistent with the distribution of mineral atmospheric DD, except in the northern BoB, where nutrients are largely supplied by the GBR. Yadav et al. (2016) have recently confirmed that atmospheric DD increases CHL even in coastal areas where the input of river-borne nutrients is high.
While the relationships between CHL and nutrients have been widely studied on seasonal timescales, the relative contributions of different sources of nutrients to determining CHL, especially on interannual timescales, have not been delineated spatially. Deducing the impacts of climate changes on primary production in the BoB and hence the BoB LME requires understanding the major sources of nutrients that determine interannual variations in CHL in the bay. Climate variations such as the El Niño/Southern Oscillation (ENSO) and the Indian Ocean Dipole (IOD) do not directly drive primary production in the ocean but instead alter nutrient inputs by modifying river discharges, upwelling/downwelling, and atmospheric circulation.
Furthermore, because there is a strong relationship between phytoplankton primary production and fisheries production, understanding the impacts of climate changes on primary production in the BoB is important for fisheries management, i.e., for planning the actions of fishery activities/industries that are needed to adapt to climate change. In this study, we employed multi-year satellite-derived CHL observations (a proxy of phytoplankton primary production), in situ GBR discharge, and models of DD to assess the region-wise influence of nutrients from different sources (i.e., rivers, ocean subsurface layer, and the atmosphere) on interannual variations of CHL in the BoB, to identify the mechanisms that may introduce nutrients from different sources into the bay, and to link those mechanisms with climate variations. Figure 1A shows the area selected for this study. The satellitederived CHL in the offshore area is low (<0.5 mg m −3 ), but the CHL in the coastal area is high (>2 mg m −3 ) because of the input of terrigenous inorganic nutrients from the adjacent river systems. The wide area of high CHL in the northern coastal area is attributable to discharges from large rivers such as the GBR (e.g., Gomes et al., 2000). On average, GBR discharge peaks in late summer (August/ September, Figure 1B). Note that the summer season used in this manuscript refers to boreal summer (June-September). This high river discharge is associated with a large amount of rainfall that peaks 1-2 months earlier (June or July) than the peak of river discharge. CHL in the coastal region increases because of the summer peaks in rainfall and river discharge (hence in riverborne nutrient supply) and is then dispersed southward during the summer ( Figure 1C).

Study area
Because of the upwelled/entrained nutrients from the ocean subsurface layer, CHL is enhanced in the northern summer and peaks in the northern winter, especially over the western/ southwestern BoB. Increases in nutrient concentrations that enhance CHL are associated with local longshore wind stress (Shetye et al., 1991;Gomes et al., 2000), westward propagation of upwelling Rossby waves (Rao et al., 2010;Sreenivas et al., 2012;Gulakaram et al., 2018), and strong northeasterly winds (e.g., Xu et al., 2021) during the summer, fall, and winter, respectively.
In addition, the BoB receives large amounts of dust from fallout that also peaks around summer (July) ( Figure 1B). Aerosol dust over the BoB is transported from southwest Asia, northeast Asia, and the Middle East (Banerjee and PrasannaKumar, 2016;Banerjee et al., 2019). During summer, a high plume of DD intrudes into the BoB from the southwest ( Figure 1A). This plume is attributable to high emissions from dust sources, high rainfall over the BoB, and a slight deflection (southeastward to northeastward) of winds around the southern Indian Peninsula (Banerjee et al., 2019).

Data sources and analysis
We used satellite-derived monthly CHL blended from multiple ocean color sensors with 4-km resolution downloaded from the ESA OC-CCI (https://esa-oceancolour-cci.org/), monthly sea surface height anomaly (SSHA, cm) with 25-km resolution downloaded from the AVISO+ Satellite Altimetry Data (https:// www.aviso.altimetry.fr/), monthly rainfall with 25-km resolution downloaded from the Asia-Pacific Data Research Center (http:// apdrc.soest.hawaii.edu/), and daily sea surface temperature (SST,°C ) with 4-km resolution retrieved by the Advanced Very High-Resolution Radiometer and downloaded from NOAA (https:// www.nodc.noaa.gov/). The daily SST data were then averaged to produce monthly mean SSTs. Variations in SST and SSHA are used widely to assess nutrient inputs from ocean subsurface layers associated with physical processes (e.g., Siswanto, 2015). We also calculated the wind stress curl (curl, N m −3 ) following Large and Pond (1981) from monthly wind fields acquired from the Cross-C a l i b r a t e d M u l t i -P la t f o r m ( C C M P ) p r o j e c t ( h t t p : / / podaac.jpl.nasa.gov) as an indicator of upwelling.
One of the objectives of this study was to determine whether climate variations influence CHL by modulating supplies of nutrients from different sources. Because climate variations have an interannual timescale, the seasonal cycles and long-term trends of CHL and other geophysical variables must be removed from their time series before the effects of changes in sources of nutrients can be detected. To assess whether the ocean subsurface layer or atmospheric deposition was a more important source of inorganic nutrients in different parts of the offshore regions of the BoB, we conducted a multiple linear regression analysis (MLRA) that related CHL to other environmental variables following the approach used by Siswanto (2015) and Nathans et al. (2012). The approach included removing seasonal variations and trends from variables and standardizing them so that they all had means of zero and standard deviations of one. Standardization facilitates the comparison of the importance of one independent variable to that of other variables in determining the variance of CHL. The MLRA required that no data be missing. An empirical, orthogonal functionbased data interpolation (e.g., Alvera-Azcarate et al., 2005) was used to fill in missing data in the CHL and SST satellite imageries. Figure   S1 in the Supplementary Information shows the numbers and percentages of valid CHL and SST pixels used for interpolation.

Spatial variability of the importance of nutrients from different sources
High CHL in the northern BoB was confined mainly to the area north of~21°N or within~200 km offshore of the mouth of the GBR ( Figures 1A, C). The time series of meridional CHL variation derived from a meridional transect south of the GBR mouth showed southward dispersion of high CHL during summer, when rainfall intensity and river discharge peak ( Figure 2A). Tan et al. (2006) and Siswanto et al. (2011) have reported that CHL in the coastal region is systematically overestimated by satellite-retrieved CHL because of the high concentration of total suspended matter, which is generally quantified by satellite-retrieved normalized water-leaving radiance at a wavelength of 555 nm (nLw 555 ). They have found that the variation of satellite-derived CHL in the coastal region can be considered to reflect real variations in phytoplankton biomass if nLw 555 <2.0 mW cm −2 µm −1 sr −1 . Figure 2C shows temporal variations of nLw 555 along the meridional transect line off the mouth of the GBR shown in Figure 1A. During the peak discharge of the GBR in summer, the area where the nLw 555 exceeded 2.0 mW cm −2 µm −1 sr −1 was confined mainly to the area north of 21°N. Variations of CHL south of 21°N, including the high CHL that is dispersed further southward, can be considered to reflect high phytoplankton biomass due to elevated primary productivity modulated by nutrients from rivers, especially from the GBR.
Along the northern coast of the BoB, nutrients from rivers are likely important for phytoplankton in coastal regions only as far as about 200 km from the river mouths. The distance from the river mouths would depend on the intensity of rainfall and interannual variations of river discharge; thus, CHL variation in the BoB offshore region is likely attributable to variations in other sources of nutrients.
To identify which source of nutrients was most important in determining the offshore variations of CHL in the BoB on an interannual basis, we used a MLRA with CHL as the dependent variable. The independent variables were DD and a physical variable There is much seasonal variability in CHL as well as in the independent variables SST, SSHA, curl, and DD. Prior to running the MLRAs, we removed seasonal cycles and long-term trends from the variables (see Section 2.2) to focus the analysis on interannual variations. Our results therefore concern interannual variability unrelated to seasonal cycles.
Negative values of b SST , which were especially common in the southwestern BoB indicated by box A in Figure 3A, indicated an increase in CHL associated with surface cooling (hence an increase in nutrient concentrations due to upwelling, vertical mixing, or both) (e.g., Shetye et al., 1991;Gomes et al., 2000;Vinayachandran and Mathew, 2003;Martin and Shaji, 2015;Xu et al., 2021). The hypothesis that nutrients introduced from the ocean subsurface layer by upwelling led to increases in CHL in the southwestern part of the bay was also supported by the negative b SSHA ( Figure 3D) and positive b Curl ( Figure 3G). A negative b SSHA implies that an increase in CHL was associated with a negative SSHA due to upwelling. A positive b Curl means that an increase in CHL was associated with a positive curl, which is also indicative of upwelling. As shown in Table 1, the values of the mean b SST (−0.25), mean b SSHA (−0.20), and mean b curl (0.26) averaged over box A were larger in magnitude than the mean b DD (0.14-0.15). The implication is that variations in the supply of nutrients from the ocean subsurface layer contributed more than variations in the supply of nutrients from atmospheric deposition to the interannual variations of CHL in the southwestern BoB. The mean R 2 (0.15-0.18) resulting from the three MLRAs (Figures 3C, F, I, Table 1) implies that variations in the supply of nutrients from the ocean subsurface layer explain 15%-18% of the interannual variations of CHL in the southwestern BoB.
In the central/eastern BoB, nutrients from upwelled/entrained water and atmospheric deposition could together explain about 15% (R 2 = 0.15-0.16, Table 1) of the interannual variance of CHL. In

Possible linkages with climate variations 3.2.1 BoB northern coastal region
There is interannual variation over the seasonal cycle of high CHL during summer/fall along the northern coast of the BoB (Figures 1C, 2A). The CHL interannual variation can be seen more obviously in terms of the variation of the meridional CHL anomaly along the meridional transect line (19. 0°N-21.5°N , Figure 2B). Expansions of high-CHL waters further southward (up to about 300 km south of the GBR mouth) were obvious during the summer/fall of 1998, 2005, 2007, and 2011. Those years were characterized by anomalously high rainfall over the catchment area and high discharges from the GBR ( Figure 2B). The years 1998The years , 2005The years , and 2011 were La Niña years and negative-IOD

BoB southwestern region
To assess whether nutrient inputs from the ocean subsurface layer that determine interannual variations of CHL in the southwestern BoB were correlated with climate variations, we plotted time series of CHL, SST, curl, and SSHA along with the Nino3.4 and DMI (Figures 4A-C). In general, CHL varied inversely with SST but positively with curl. Such an inverse relationship between CHL and SST and a positive relationship between CHL and curl can also be seen in the corresponding scatter plots (Figures 5A, B). The correlation coefficients between CHL and SST (−0.43, DF = 224) and between CHL and curl (0.43, DF = 224) were statistically significant (two-tailed t-test, p <0.0001) ( Table 2). The negative correlation between CHL and SST indicated that entrained cold (hence nutrient-rich) water was associated with high CHL. The entrainment could have been driven by upwelling and/or vertical mixing. The positive curl would be associated with upwelling or Ekman divergence processes. The positive correlation between CHL and curl therefore indicated that upwelling was the major physical process that entrained nutrients from the ocean's subsurface layer to enhance CHL.  Interestingly, in general, the negative correlation between CHL and SST ( Figures 4B, 5A) and the positive correlation between CHL and curl ( Figures 4B, C, 5B) were likely linked with ENSO, IOD, or both. For instance, negative CHL, positive SST, and negative curl variations were obvious during the concurrent El Niño and positive IOD (hereafter El Niño/pIOD) in 1997, El Niño/pIOD in 2009, and the 2017 El Niño. In contrast, positive CHL, negative SST, and positive curl anomalies were evident during the La Niña/nIOD of 1998, the La Niña of 1999, the La Niña/nIOD of 2005, and the La Niña/nIOD of 2013. This pattern indicates that atmosphere-ocean interactions associated with interannual climate variations of ENSO and IOD drive the supply of inorganic nutrients from the ocean subsurface layer, and an increase in the supply of nutrients (indicated by low SST) to the surface layer leads to increases in CHL.
Compared to SST and curl, SSHA was more weakly correlated with CHL (0.21, DF = 224), but the correlation was still statistically significant (p <0.005) ( Table 2). The strong seasonality of SSHA in the BoB is associated with westward-propagating Rossby waves (Sreenivas et al., 2012;Gulakaram et al., 2018) that likely weaken the relationship between CHL and SSHA on an interannual timeframe. Westward-propagating Rossby and coastal Kelvin waves (that traverse along the BoB coast) are generated at the eastern boundary of the BoB. They are feedback from eastwardpropagating equatorial Kelvin waves.

BoB central/eastern region
In the central/eastern BoB, the time series of CHL, SST, and DD anomalies show that CHL generally varies inversely with SST but positively with DD. This pattern is apparent in the corresponding scatter plots (Figures 4D, 5C, D). The strongest interannual association was between CHL and DD, for which the correlation coefficient was 0.30 (p <0.0001, DF = 244, Table 2, Figure 5D). CHL and SST were also significantly correlated (r = −0.21, p <0.005, DF = 224, Table 2, Figure 5C). CHL was significantly correlated with curl (r = 0.14, DF = 224) ( Figure 4E) but with a lower level of confidence (p <0.05). There was no significant correlation between CHL and SSHA (r = −0.001, p >0.05, DF = 224, Figures 4D, E, Table 2).
Despite the significant correlation between CHL and DD in the central/eastern BoB, neither CHL nor DD was clearly associated with either ENSO or IOD. Neither CHL nor DD were significantly correlated with those climate indices (Table 3). The lack of association between CHL, DD, and variations of climate is likely due to the time lag between the time when the climate event peaked and the time of DD entering the BoB (see Discussion). In contrast, as shown in Table 3, the significant correlations between SSHA and Nino3.4 (r = −0.18), between SSHA and DMI (r = −0.35), between curl and Nino3.4 (r = −0.32), and between curl and DMI (r = −0.20) indicate that curl and SSHA were influenced by both ENSO and IOD.
Climate variations had opposite effects on SSHA and curl ( Table 3). The negative correlations between the SSHA and climate indices (Table 3) indicate that El Niño and pIOD led to upwelling, but the negative correlations between the curl and climate indices (Table 3) indicate that they led to downwelling. Regardless of these apparently contradictory relationships, one of the reasons for the lack of association between CHL and climate variations in the central BoB was probably the formation of a freshwater lens associated with the influx of freshwater from river discharge that restricted climate-related atmosphere-ocean interactions (e.g., Kumar et al., 2002;Vinayachandran et al., 2002;Felton et al., 2014;Kay et al., 2018).

Factors responsible for anomalously high CHL
Here we will elaborate on the atmosphere-land-ocean interactions that were probably responsible for the remarkably high CHL anomalies in the coastal region of the northern BoB (Figure 2), box A in the southwestern BoB, and box B in the central/ eastern BoB (Figure 3). We first discuss the mechanisms responsible for the high CHL in the coastal region, followed by the mechanisms probably responsible for the high CHL in the southwestern and central/eastern BoB.
As reported by Perves and Henebry (2015), the amount of rainfall over the GBR basin increases by more than 100% (compared to the long-term mean rainfall) during years of La Niña alone or La Niña/ nIOD. The anomalously high CHL and/or southward dispersion of high-CHL waters along the northern coast of the bay during the La Niña/nIOD or La Niña years of 1998, 2005, 2007, and 2011 can be attributed to high river discharges and hence large influxes of river-  1997, 2002, 2009, and 2015, there was no apparent southward dispersion of high-CHL waters (Figures 2A, B) because of the low amounts of rainfall and low river discharges during those El Niño years (Perves and Henebry, 2015). The implication is that there is a teleconnection between basin-scale climate variations of ENSO and/or IOD and land-ocean interaction, i.e., variations of river-borne inorganic nutrient supply and primary production in coastal waters. Between September 1997 and May 2016, there were two obvious, anomalous peaks of CHL in the southwestern BoB that were accompanied by negative SST, negative SSHA, and positive curl anomalies (Figures 4B, C). The first CHL peak occurred during December 2005 ( Figure 6B), a year influenced by La Niña/nIOD ( Figure 6A) (see also Li et al., 2017;Chen and Li, 2018). That CHL peak resulted from a spatially large, positive CHL anomaly ( Figure 6C-2nd row). That spatially large, positive CHL anomaly was enabled by nutrients supplied from the ocean subsurface layer, as indicated by the positive curl, negative SST, and negative SSHA anomalies ( Figure 6C-3rd, 4th, and 5th rows). Previous work has suggested that the high CHL in December 2005 resulted from the passage of cyclones (Ali et al., 2007;Chen et al., 2013;Sridevi et al.,

2019
). An input of nutrients from the atmosphere is unlikely the explanation for the high CHL because the DD from October 2005 to February 2006 was normal ( Figure 6C-1st row), and DD is generally low in December ( Figure 1B) (see also Solmon et al., 2015).
A similar combination of anomalously high CHL accompanied by positive curl, negative SST, and negative SSHA was also apparent in December 2013 ( Figures 7B, C-2nd, 3rd, 4th, and 5th rows). These observations indicate that the upwelling of nutrients enhances CHL. Atmospheric inputs were not a cause of the enhancement of CHL because the DD anomalies were negative ( Figure 7C-1st row), and December is in the season of low DD ( Figure 1B) (Solmon et al., 2015). Oceanographic conditions during this period were also influenced by the 2013 La Niña/nIOD ( Figure 7A). Jayaram et al. (2019) have suggested that two cyclones that passed over adjacent areas caused the highest CHL anomaly observed in December 2013.
Regardless of the inorganic nutrient composition of DD, the fact that previous studies have shown that increased phytoplankton growth in various ocean basins can be associated with the SPRINTARS-modeled DD (Kitajima et al., 2009;Calil et al., 2011;Fukushima, 2014) indicates that the deposition of inorganic nutrients from the atmosphere is an important source of nutrients for phytoplankton growth. We therefore expected that the SPRINTARS-modeled DD would to some degree alter CHL in the BoB, particularly in the central/eastern bay, where nutrient supplies from river discharge are generally considered to be low. Indeed, we found a significant correlation between CHL and DD in the central/ eastern BoB.
However, the role of atmospheric deposition, which is suggested by the significant, positive correlation between CHL and DD (r = 0.30, Table 2), must be interpreted with caution because there was a significant negative (positive) correlation between CHL and SST (CHL and curl) ( Table 2); therefore, nutrient inputs from the ocean subsurface layer may also have played an important role. For instance, from May 2012 to May 2014, positive DDs were accompanied not only by positive CHL but also by negative SST (Figure 4D). During this study, only the DD in July-August 2006 likely indicated atmospheric deposition that enhanced CHL in August 2006 ( Figures 8B, C-1st and 2nd  rows). However, DD might contribute some nutrients because

Discussion
In general, CHL in the southwestern BoB was negatively and significantly correlated with Nino3.4 (r = −0.27, p <0.0001, Table 3). CHL therefore tended to be low (high) during El Niño (La Niña). Similarly, there was a negative (but insignificant) correlation between CHL and DMI, i.e., CHL was low (high) during pIOD (nIOD). A similar negative but insignificant CHL-DMI correlation in the southwestern part of the bay has also been reported by Xu et al. (2021). Our analysis indeed showed that ENSO rather than the DMI seemed to have a larger influence on the SST and curl (Table 3). We therefore hypothesize that elucidation of the impacts of climate variations on the biogeophysical variability of the BoB should include consideration of ENSO cycles.
Except for the southwestern BoB, much of the BoB during the La Niña/nIOD period from October 2005 to February 2006 was characterized by a positive SSHA ( Figure 6C-5th row). This positive SSHA is known to be associated with a strengthened downwelling Kelvin wave in the equatorial Indian Ocean that is driven by a strong westward wind stress anomaly during La Niña/nIOD (Sreenivas et al., 2012;Roman-Stork et al., 2021). This downwelling equatorial Kelvin wave is reflected westward from the eastern coast of the bay as a downwelling Rossby wave and is characterized by positive SSHA (Rao et al., 2010). However, the SSHA tends to be negative over the southwestern BoB because during La Niña/nIOD, a strong westward wind stress anomaly creates cyclonic wind stress over the southwestern BoB (see Sreenivas et al., 2012).
The 2005 La Niña/nIOD-driven cyclonic wind stress led to the generation of a mesoscale cyclonic eddy over the same area of high CHL in the southwestern BoB ( Figures 6C-2nd (Ali et al., 2007;Sridevi et al., 2019). However, during October-November 2005, before the arrival of cyclone Fanoos, the positive curl and negative SSHA indicated that a mesoscale cyclonic eddy already existed ( Figure 6C-3rd and 5th rows; see also Sridevi et al., 2019). An enhancement of CHL was already apparent in October 2005 ( Figure 6C-2nd row). An Ekman divergence associated with La Niña/nIOD was already occurring and was likely strengthened by cyclone Fanoos. The result was a pronounced positive CHL anomaly in December 2005. The increase in CHL could not be attributed to atmospheric input of nutrients because the DD was normal during the period from October 2005 to February 2006 ( Figure 6C-1st row).
Analysis of the CHL enhancement during the December 2013 La Niña/nIOD requires consideration of the fact that two cyclones passed over the western/southwestern BoB from late November to early December 2013: Cyclone Leher (24-28 November) and Cyclone Madi (5-12 December) (Jayaram et al., 2019). The high CHL in December 2013 (Figures 7B, C-2nd row) was the result of a CHL bloom from 25 November to 2 December 2013 that preceded Cyclone Madi's arrival (Jayaram et al., 2019). The implication is that Cyclone Madi was not the only driver of the high CHL in December 2013. Although patches of increased CHL along the track of cyclone Leher are apparent in the study of Jayaram et al. (2019), the locations and spatial patterns that show Cyclone Leher-driven Ekman pumping and that show increased CHL are not in the same places. It is therefore unlikely that the CHL bloom in the southwestern BoB was associated with the passage of Cyclone Leher. The positive curl and negative SSHA ( Figure 7C-3rd and 5th rows) suggest that a cyclonic eddy already existed in October-November 2013 (see also Sridevi et al., 2019). We therefore suggest that the cyclonic eddy that formed in association with the 2013 La Niña/nIOD and the cyclones may have synergistically enhanced Ekman divergence and produced a prominent, positive CHL anomaly in December 2013. The presence of a cyclonic eddy that preceded cyclone Madi has also B A been mentioned by Chen et al. (2013) and Chowdhury et al. (2022), but those authors did not discuss the connection between the cyclonic eddy that preceded cyclone Madi and the 2013 La Niña/nIOD. Mahala et al. (2015) and Roose et al. (2022) have recently mentioned that more cyclones are generated in the BoB during La Niña, nIOD, or both. This pattern is associated with anomalous equatorial westerly winds during La Niña, nIOD, or both that favor cyclone formation. Equatorial westerlies also favor cyclonic eddy formation in the southwestern BoB (see Sreenivas et al., 2012). The above-mentioned synergistic effect of mesoscale cyclonic eddies and tropical cyclones over the southwestern BoB is therefore likely to be more frequent during La Niña, nIOD, or both, and hence the CHL in the southwestern BoB will be greatly enhanced during those times.
The negative anomalies of DD from October 2013 to February 2014 ( Figure 7C-1st row) and the seasonal minimum of DD during this period ( Figure 1B) (see also Solmon et al., 2015) confirmed that nutrients upwelled from the ocean subsurface layer by climatedriven cyclonic eddies and cyclone passages were the only sources of nutrients that enhanced CHL in December 2013. The enhancements of CHL during the 1998 La Niña/nIOD, 1999 La Niña, and 2010 La Niña/nIOD were less remarkable than the enhancements during the 2005 and 2013 La Niña/nIOD. One reason was that the areas where phytoplankton bloomed during those climatic events were not entirely within the southwestern BoB (maps not shown). The entire CHL bloom was therefore not counted when the CHL within box A was averaged in the southwestern BoB.
Based on field experiments, Yadav et al. (2016) have also mentioned that inorganic nutrients from dust enhance CHL by factors of 1.5-4 in the coastal region of the western BoB. Patches of positive b DD around western bay coastal areas in our study that were approximately the same as the study area of Yadav et al. (2016) confirmed that phytoplankton production in coastal waters also seemed to be sensitive to atmospheric DD.
We observed no correlations between DD in the central/eastern BoB and climate variations in the same months (no time lag) (Table 3). There may, however, be a time lag between climate variations and DDs (e.g., Banerjee and PrasannaKumar, 2016). To determine whether there was a time lag between climate variations and DDs in the BoB, we investigated the relationship between Nino3.4 and DD with different time lags. Figure 10A shows the correlations and reveals that DDs were highly correlated with Nino3.4 variations that occurred 6-8 months earlier (r <−0.20, p <0.01). The highest correlations were observed with a 7-month lag (Nino3.4 leading, r = −0.24, DF = 217, p <0.005). Interestingly, high correlations between Nino3.4 and CHL (r <−0.20, p <0.01) were also observed with time lags of 6-7 months ( Figure 10B).
The probable explanation for the negative correlations between Nino3.4-DD and Nino3.4 versus CHL 6-7 months later is that La Niña causes anomalously low precipitation in the dust source regions (i.e., southwestern Asia, the Middle East, and Eastern Africa) that leads to more dust generation during the following summer (Banerjee and PrasannaKumar, 2016). High DD that partially enhanced CHL in August 2006 (Figures 8B, C-1st and 2nd rows) was apparent following the 2005 La Niña peak (January, Figure 8A). Similar La Niña-DD lagged co-variations were also apparent for other years: a La Niña peak in December 1998 preceded a DD peak in June 1999; a La Niña peak in January 2008 preceded high DD in May-June 2009; a La Niña peak in December 2011 preceded high DD in June-July 2012; and a La Niña in December 2012 preceded high DD in May-June 2013 ( Figures 4A, D). In contrast and consistent with Banerjee and PrasannaKumar (2016), El Niño peaks were followed by low DD in the following summers: an El Niño peak in December 1997 preceded low DD in July 1998, and an El Niño peak in December 2009 preceded a very low DD in July 2010.
Because the deposition of dust generated from the same sources mentioned above is important in determining primary production in the Arabian Sea (Tandule et al., 2022) and because we found a significant correlation between DD and CHL in the central/eastern BoB (r = 0.30, Table 2), we suggest that DD plays an important role in determining primary production in the BoB, especially in the central/eastern bay. In addition, the Nino3.4-DD and Nino3.4-CHL ( Figure 10) lagged correlations indicated that ENSO influences DD (and hence CHL) in the BoB 6-8 months after the peak of ENSO.
However, when assessing the impact of DD on CHL in the central/eastern bay, it is important to note that there may also have been a contribution of nutrients from other sources, even during the high DD in August 2006 ( Figure 8B, C-1st row). During August, nutrient inputs from the ocean subsurface layer are likely to be limited by the formation of a barrier layer associated with the flux of freshwater from the GBR that restricts climate-related atmosphereocean interactions (e.g., Vinayachandran et al., 2002;Felton et al., 2014;Kay et al., 2018). Moderate increases in CHL in the central/ eastern BoB (box B) may therefore be largely associated with DD. But episodic cyclonic eddies, such as the eddies that form off the western coast ( Figure 9B), may also supply nutrients from the ocean subsurface layer and enhance CHL in the coastal region. Coastal waters with high CHL may be transported into the central/eastern bay by surface circulation. A more quantitative analysis would be required to assess the contributions of nutrients from the atmosphere and ocean subsurface layer separately (e.g., Siswanto, 2015).
Understanding how climate impacts phytoplankton CHL in the BoB by altering nutrient inputs from different sources would be valuable for local fisheries management because fisheries production in the BoB is largely controlled by phytoplankton primary production, and a common metric of phytoplankton production is CHL (Hossain et al., 2020). Climate change may reduce or increase local fisheries production by similarly changing CHL. This report may therefore be used as a reference by fisheriesrelated stakeholders in planning necessary adaptations to fisheries activities in response to climate change.

Conclusions
This study revealed that variations of nutrient inputs (inferred from variations of physical parameters) from different sources (river discharge, ocean subsurface layer, and atmosphere) cause variations of phytoplankton CHL in different parts of the BoB from year to year. Nutrients from the GBR are important for CHL mainly in the coastal region at distances up to~200 km from the northern coast of the BoB. By increasing the amount of rainfall (hence the supply of river-borne nutrients from the GBR), La Niña can extend high-CHL waters further southward (~300 km from the coast). Variations in the supply of nutrients from the ocean subsurface layer determine interannual variations of phytoplankton CHL, mainly in the southwestern BoB. This variation in nutrient supply is linked with ENSO, IOD, or both through the generation of mesoscale cyclonic eddies during La Niña/nIOD. Such La Niña/nIOD-related mesoscale cyclonic eddies, together with the passage of cyclones, may synergistically intensify Ekman divergence, which results in a pronounced CHL enhancement in the southwestern bay. Nutrients deposited from the atmosphere make a nontrivial contribution to the interannual variations in CHL in the central/eastern BoB. This study also suggests that there is a lag of 6-7 months between the peaks of ENSO and the variations of DD and associated enhancements of CHL in the central/eastern BoB. Because phytoplankton are the base of the marine food web, local fishery production is also expected to be affected by climate variations that modulate CHL.

Author contributions
ES conceptualized and led the study. MS assisted with river discharge data acquisition and analysis. BP assisted with ocean surface circulation computation from global drifter and altimetry data. TT modeled and provided atmospheric dust deposition data. TH assisted with discussion concerning climate changes and teleconnection of physical oceanographic process. KM, FT, and MH assisted with discussion on the impacts of atmospheric deposition on ocean primary production. MH assisted with financial support and led the research project. All authors contributed to the article and approved the submitted version.