Use of multivariate analysis to identify phytoplankton bioindicators of stream water quality in the monomodal equatorial agroecological zone of Cameroon

The aquatic ecosystem is compromised by many contaminants that may reduce ecosystem functions and severely affect human health. This study aimed at determining suitable phytoplankton bioindicators of water quality for biomonitoring of freshwater streams in the monomodal agroecological zone of Cameroon. Water physicochemical and hydrological parameters, together with phytoplankton abundance and diversity, were measured from June 2016 to May 2017 along the Benoe Stream. Principal component analysis and redundancy analysis were used to determine phytoplankton spatial and temporal distribution and identify indicator species. The Shannon–Wiener diversity and Pielou’s evenness indices indicated a clean to mildly polluted stream with a diverse phytoplankton community consisting of 84 genera belonging to 51 families that was dominated by the Bacillariophyta (64%), followed by Chlorophyta (13%) and Cyanophyta (10%). The total dissolved solids, electrical conductivity, stream water velocity, and discharge were the most important stream characteristics affecting the abundance of the dominant phytoplankton genera. Seasonal variations in the stream characteristics as well as spatial community distribution along an urban—small-scale farming – large-scale farming gradient were unveiled and their influence on the phytoplankton relative abundances. Increased abundance of Synedra ulna was indicative of low TDS and EC, which was the contrary for Gyrosigma baltium dominance. High Pleurosira laevis abundance was associated with the urban zone while high Diatoma sp. and Oscillatoria sp. abundances were related to the large-scale farming zone of the stream. These phytoplankton species have good potential for use as bioindicators for stream water quality monitoring in the monomodal agroecological zone.


Introduction
Water bodies harbour a very rich biodiversity (Biggs et al., 2017), which is unfortunately threatened by contamination risks from agricultural, industrial, and domestic activities (Egessa et al., 2020;Hanif et al., 2020;Kasonga et al., 2021;Shah et al., 2020). Aquatic pollutants from point and non-point sources such as surface run-off from urban areas and farmlands cause water quality deterioration globally (Shah et al., 2020).In developing countries, waste management is particularly challenging (Gallego Schmid & Tarpani, 2019) and as a consequence, water bodies often serve as dumping grounds for domestic wastes. In addition, agriculture is the backbone of the economy of Cameroon, like many other developing nations (ZEF FARA IRAD, 2017). The crops are grown either by small-scale farmers or agro-industrial complexes such as the Cameroon Development Corporation (CDC). CDC was created in 1947 with the objective to acquire, develop and operate intensive cultivation of tropical crops. This has led to the extensive use of huge amounts of chemical inputs over the years in the hot and humid coastal belt of Cameroon, which makes up the mono-modal equatorial agro-ecological zone. Several streams, including the Benoe Stream, flow through urban areas as well as the CDC banana plantations and are used for various purposes by the CDC and the local population. The Benoe Stream takes its rise from Mount Cameroon and flows through several villages and towns including Buea, Tiko and Mutengene with rapidly growing populations who rely on it for domestic, agricultural and other uses such as car washing, and swimming. Lambi and Shende (2009) revealed a direct relationship between rural development and water availability and its exploitation among villagers on the eastern slope of Mount Cameroon, showing the importance of the Benoe stream for the development of the villages and towns through which it flows. On the other hand, aquatic organisms in such streams are exposed to a large number of chemical contaminants in the water column, sediments, and in their food preys (Amiard & Amiard-Triquet, 2015). However, no appropriate environmental health or monitoring scheme has been put in place to ensure good water quality and prevention of adverse health effects on the aquatic community. Simple, reliable, and cost-effective stream bioindicators are needed for inclusion in biomonitoring schemes to ensure sustained aquatic ecosystem functioning. Any organism characteristic of a particular set of environmental conditions qualifies for the term bioindicator (Suthers et al., 2019). Algae are abundant in freshwater, forming the base of the aquatic food chain and playing a tremendous role in freshwater productivity, and have been used as ecological bio-indicators worldwide for a long time (Ahmed Mohammed & Nodhy Mahran, 2022;França et al., 2017;Inyang & Wang, 2020;Yusuf, 2020). The plankton community responds both to natural changes and anthropogenic pressures (Ganai & Parveen, 2014) and are highly sensitive to pollutants (Onana et al., 2014) as well as to nutrient changes as they accumulate over days the effects of hourly changes in water quality (Suthers et al., 2019). According to Yusuf (2020), phytoplankton communities give more evidence concerning alterations in water quality than nutrient or chlorophyll-a concentration.
Phytoplankton communities, like other biological communities, are characterized by significant temporal and spatial variability resulting from a complex interplay between various processes occurring at multiple scales (Dray et al., 2012;Shitikov & Zinchenko, 2019). Several biodiversity indices have been used to identify relationships between the taxonomic structure of biological communities with environmental factors which aim to quantify the richness, evenness, variety of species present and the trophic state of aquatic ecosystems. The most commonly used indices in phytoplankton studies include, species richness, Simpson's index, Shannon-Wiener index, Pielou index and Margaleff index (Motwani et al., 2014;Spatharis et al., 2011;Tas, 2021;Wei et al., 2022). In addition, bioassessment methods have led to the development of a plethora of biotic indices, many of which use diatoms as bioindicators to assess environmental conditions of aquatic ecosystems (Celekli et al., 2021). However, these biotic indices have been considered by some as extreme simplification that reduces multiple species information to a single synthetic value (Dray et al., 2012). In addition, Spatharis et al. (2011) reported high variability and non-linearity of the Shannon and Simpson indices along a eutrophication gradient and suggested that these measures of diversity may be inappropriate for use as water quality monitoring assessment tools. According to Buckley et al. (2021) multi-taxon, multivariate datasets can be more sensitive to changes in environmental variables than univariate summary statistics, such as taxon richness, diversity, or biomass. Thus, multivariate analysis represents a much richer alternative that takes better advantage of the multi-dimensional nature of ecological data (Dray et al., 2012). These methods can be used to decipher the occurrences and/or successions of diatom assemblages in different environmental conditions to determine their trophic weight and indicator value (Celekli et al., 2021). Biodiversity indices, however, have their merits due to their quantitative and dimensionless nature, independence of sample size, and the fact that they lend themselves to statistical analysis. Therefore, ecological indices currently used to evaluate stream health include biotic, multivariate, and multimetric indices (Mwedzi et al., 2022).
Few studies from Cameroon have examined the phytoplankton communities in running water (Anyinkeng et al., 2016;Fonge et al., 2015Fonge et al., , 2016Menye et al., 2012;Noumssi et al., 2021), lakes (Awo et al., 2020;Kemka et al., 2004Kemka et al., , 2009; the Douala estuary (Fonge et al., 2013) and floodplains (Fonge et al., 2012). Several of these studies used diversity indices like species richness, Shannon diversity, evenness and Sorenson similarity indices (Cibils-Martina et al., 2017;Fonge et al., 2016;Tas, 2021). Some other phytoplankton indices that have been previously used to assess the ecological status of aquatic ecosystems include: the Palmer's pollution index, Phytoplankton trophic index, phytoplankton assemblage index, Specific Pollution-Sensitive Index, Trophic Index, Eutrophication and/or pollution index-diatom, Diatom Species Index Australian Rivers, Trophic Diatom Index, Duero Diatom Index, Trophic Water Quality Index, Richmond River Diatom Index, Trophic Index Turkey, saprobic index and Diatom Ecological Quality Index (Anyinkeng et al., 2016;Celekli & Ozturk, 2014;Celekli & Yildiz, 2022;Celekli et al., 2021;Fonge et al., 2013). Previous studies on rivers in the monomodal agroecological zone (Fonge et al., 2015 and2016) focussed on the trophic structure of the streams and the studies were cross-sectional giving a snapshot of the phytoplankton community present at the time with no indication on the temporal variations within the stream. In addition, while it was possible to compare the phytoplankton communities at specific points between streams, the sampling strategy did not permit the study of the phytoplankton spatial distribution along a stream flowing through different zones influenced by various anthropogenic activities. Multivariate analysis methods such as principal component analysis (PCA), nonmetric multidimensional scaling, redundancy analysis (RDA), and canonical correspondence analysis (CCA), have been used in previous studies to assess differences in the composition and the structure of phytoplankton species in water bodies as well as to determine the relationship between phytoplankton species and environmental variables (Celekli & Ozturk, 2014;Celekli & Yildiz, 2022;Kumar et al., 2020;Zongo et al., 2019). The main objective of the present study was to assess the potential of using phytoplankton for biomonitoring of streams in the monomodal agroecological zone of Cameroon. Specifically, this research set out to 1) use phytoplankton diversity indices to determine the ecological status of the Benoe Stream; 2) use multivariate analysis to determine the effects of seasonality, land use patterns and stream water hydrological and physicochemical characteristics on the phytoplankton population dynamics in the Benoe Stream; and 3) identify potential phytoplankton indicator species in the Benoe Stream.

Material and methods
The study area and sampling stations The study was carried out in the Benoe Stream, a typical tropical freshwater ecosystem of the mono-modal agroecological zone of Cameroon, located in the Fako Division, South West Region of Cameroon. The Humid Forest Monomodal Agro-Ecological Zone (AEZ) is zone four (AEZ IV) of five AEZs in Cameron having equatorial forest with rainfall throughout the year peaking in July and August and a temperature range of 21 to 37 0 C (Yengoh et al., 2010). The dry season runs from November to February and the rainy season from March to October (WorldData.info, 2023;Yengoh et al., 2010). On average, rainfall is much more than the rest of the country; annual rainfall in Fako varies between 3000 and 3500 mm (Akoachere et al., 2018). Located at the base of Mount Cameron and close to the Atlantic coast of Cameroon, the Benoe Stream has its sources around Small-Soppo (Buea). It passes through Mutengene and Tiko where it enters the forest to join the sea around Bwenda after crossing CDC Ndongo (2275 ha) and CDC Holforth (more than 3000 ha) estates (Fig. 1). Based on land use, point sources, tributaries and the location of farms, eight stations were chosen and named Benoe 1-8 (B1-8). A GPS device (Garmin ETREX-H™) was used to get the coordinates and elevation of each sampling station. Station B1 was close to the source. Station B2 was located in the urban zone after the Mutengene abattoir just before the CDC mini dam where swimming and car washing activities take place. Located near small-scale farms, station B3 was after the CDC mini dam. Station B4 was located at the entrance of the CDC Ndongo farm i.e., before the large-scale (CDC) farms. Located between the CDC Ndongo and Holforth banana farm, station B5 receives wastes water from the Ndongo parking house, where bananas are processed after harvesting, and receives pollutant residues from the Ndongo farm. Toward the middle of the Holforth farm, station B6 is subject to wastewaters from Holforth parking house and a tributary from the Holforth farm. Situated towards the end of the Holforth farm, station B7 assessed the contribution of a tributary from the said plantation. Station B8 was after the "Cabbe Bridge" and was the exit point at which the stream enters the forests to join the sea (Fig. 1).

Analysis of the hydrological and physicochemical quality of water
From June 2016 to May 2017, triplicate samples were collected once a month at each of the eight sampling stations along the stream. The triplicate samples for each station were combined into a composite sample and hydrological and physicochemical parameters were measured in order to determine their relationships with the phytoplankton diversity and abundance. Temperature, total dissolved solids (TDS), and conductivity were measured in situ using a COM-100™ Multimeter. The amount of dissolved oxygen (DO) was measured in situ using a Milwaukee MW600™ Dissolved Oxygen Meter. Water pH was measured in situ using a pH meter (Dr. Meter™). Salinity was measured in situ with an EC70™ salinity meter. Flow velocity and flow rate or discharge (i.e., the volume of water passing a section of a stream per unit time) were determined in situ using the float method (Gordon et al., 2004;Hauet et al., 2008). For parameters measured in the laboratory, water samples were collected just below the surface using 2-l plastic bottles and transported to the laboratory in an icebox in order to maintain the samples around 4 °C. Turbidity (turbidity meter), colour, orthophosphates (molybdenum blue method), suspended solids (differential weighing), nitrogenous compounds (nitrates and nitrites) were determined following standard methods (APHA, 2005). Water color (in Pt-Co) was determined spectrophotometrically following the method of Pauwels et al. (1992) by establishing and using a calibration absorbance curve from demineralized water.

Phytoplankton collection and identification
From June 2016 to May 2017, phytoplankton samples were collected monthly using a LaMotte™ plankton net, mesh size 30 µm (Mahar, 2004) at the eight sampling stations, at the same time as water samples for hydrological and physicochemical parameters. Samples were preserved in a 3% formalin solution (Suthers et al., 2019), taken to the Zoology Laboratory (University of Buea) for identification and counting of the phytoplankton cells, using a Sedgwick-Rafter counting chamber under an Olympus BH-2 light microscope equipped with a micrometre and the abundance expressed in cells/m 3 . The identification was done using phytoplankton identification keys (Bellinger & Sigee, 2010;Mahar, 2004;Suthers et al., 2019;van Vuuren et al., 2006;Verlecar & Desai, 2004).

Data processing and analysis
Data was compiled using Microsoft Excel 2019 and summary statistics were done in SPSS 14.0 software. The list of identified plankton taxa was established and two main diversity indices were computed. The Shannon-Wiener diversity index, H′ was calculated according to Eq. 1 (Shannon, 1948). where: (1) -H′ = Shannon-Wiener Diversity Index -n i = Number of individuals of species i -S = total number of species or taxonomic richness -N = the total number of individuals of all species found (N) H′ takes into account species abundance and it is sensitive to species with low frequencies (Martensson, 2016). Typical values are generally between 1.5 and 3.5 in most ecological studies (Türkmen & Kazanci, 2010) and the index is rarely greater than 4. It increases as both the richness and the evenness of the community increase (Inyang & Wang, 2020) H′ can be interpreted as follows (Shekhar et al., 2008): -H′ > 4: clean water -3 ≤ H′ ≤ 4: mildly polluted water -2 ≤ H′ < 3: moderately polluted water -H′ < 2: heavily polluted water The Pielou Evenness Index (J) was computed according to Eq. 2 (Pielou, 1969). The value of J varies between 0 (one species dominates) and 1 (all the species tend to have the same abundance). where: -J is the Pielou Index -H′ is the Shannon-Wiener diversity index -S is the total number of species or taxonomic richness A correlation plot between biotic and abiotic variables was performed using the Corrplot package function in the R software (R Core Team, 2021).
Multivariate ordination methods such as principal component analysis (PCA) are useful for analysing temporal community dynamics as they allow optimal projection of data with a large number of variables into low-dimensional spaces to detect important ecological gradients in communities (Shitikov and Zinchenko, 2019). Multivariate data were analysed in the Canoco 5 software (ter Braak & Šmilauer, 2018). PCA was used to examine the relationships among all hydrological and physicochemical variables describing water samples from different sites but only variables that loaded above 15% were considered in further analyses. These variables included water turbidity, colour, velocity, discharge and conductivity, as well as the total suspended solids (TSS) and total dissolved solids (TDS), along with concentrations of orthophosphates (PO 4 3− ), nitrite (NO 2 − ) and nitrate (NO 3 − ). The data unit was represented by sampling events at each sampling site (n = 94). This analysis allowed us to reduce the dimensionality of variables while minimizing information loss by finding new variables that are linear functions of those in the original dataset.
Linear models can be used to regress ordination axis scores (as a response variable) generated by PCA against time and other explanatory variables of interest where the first axis represents the gradient explaining the greatest amount of variation in composition, followed by the second axis that explains the next greatest amount (Buckley et al., 2021). The redundancy analysis (RDA) analysis has been recommended as the best linear method for investigating community-environment relationships (Legendre & Fortin, 2010). A second-order partial RDA was performed to calculate the effect of environmental variables (i.e., independent variables) on the phytoplankton community structure (i.e., response variables). The environmental variables included scores from the first ordination axis of the previous PCA analysis (representing a gradient between water velocity, discharge, conductivity and TDS), and scores from the second ordination axis of the previous PCA analysis (representing a gradient of water colour, turbidity, amount of suspended solid particles, salinity, temperature, pH, dissolved oxygen as well as phosphate, nitrite and nitrate concentrations. In addition, the dry/ rainy season and zone of sampling site (urban, smallscale farming, large-scale farming) were included. The phytoplankton community was sorted according to the total abundance of each genus in the stream and the 7 most abundant genera (having a mean density > 900 cells/m 3 ) were selected for further analysis. Each phytoplankton genus was represented by its percentage abundance relative to the other genera present at a site (i.e., response variables). The percentage data was further log-transformed prior to the analysis.
All visible trends among percentages of algae genera and environmental variables were further tested by specific analyses using Statistica 13 software (TIBCO Software Inc., 2017). The non-parametric statistical methods, Mann-Whitney U test and Kruskal-Wallis test with post-hoc tests were used for comparisons. Linear or Non-linear estimation using loss function (y = a/(b + x)) was used to describe the effect of transformed scores from the first PCA axis on the percentage of genus for genera that were correlated with the scores from the first PCA. The transformation was done by shifting all the values to be only above zero allowing us to use the function correctly.

Results
Phytoplankton diversity indices and the ecological status of the Benoe Stream A total of 84 genera from 8 major phytoplankton phyla were identified and enumerated. These were distributed among 14 classes, 38 orders and 51 families (Table S1). Bacillariophyta was the most abundant phylum with 63% of all individuals (Fig. 2a) as well as the most diverse having 47% of the phytoplankton orders, 51% of the families and 50.6% of the species. Bacillariophyceae was the largest of its 4 classes, containing 72% of the 18 identified orders (Table A1), six of which made up 94% of the Bacillariophyte abundance with two of the orders (Triceratiales and Naviculales) making almost 50% of the Bacillariophyte abundance (Fig. 2b). In general, six of the ten most abundant phytoplankton species belonged to the Bacillariophyta. These were Pleurosira laevis Compère, 1982 (Triceratioceae), Synedra ulna (Nitzsch) Ehrenberg, 1832 (Fragilariaceae), Navicula rhynchocephala Kützing, 1844 (Naviculaceae), Gyrosigma baltium (Ehrenberg) Rabenhorst, 1853 (Coscinoiscoideae), Diatoma vulgaris Bory de Saint-Vincent, 1824 (Tabellariaceae) and Melosiria varians (Melosiraceae). The Chlorophyta, Cyanophyta and Dinophyta taken together only constituted a third of the Bacillariophyta abundance (Fig. 2a) and 38% of the number of species (Suppl . Table 1). However, two Chlorophyte genera, Desmidium sp. and Bacteriastrum sp. and one Cyanophyte genera, Oscillatoria sp. Gomont, 1893 (Oscillatoriaceae) were among the top ten most abundant taxa. Despite contributing only 0.2% of the total abundance, the phylum Charophyta contained one species, Spirogyra porticalis (O. F. Müller) Dumortier 1822 (Zygnemataceae), that featured among the top ten most abundant genera. It is also worth noting that some genera were only found in one sampling station, such as Bacteriastrum sp. (Chlorophyta) which was only found at station B2; Dinobryon sp. (Chrysophyta) only found at B3 and Dinophysis sp. (Dinophyta) at B5 only with 15 cells/m 3 (Table S1). The total number of genera per site ranged from 34 to 45 except for B4, which only had 29 genera. The overall phytoplankton diversity of the stream indicated by the Shannon-Wiener diversity index (H′) was high at all sampling stations (ranging from 3.5 to 4) indicating that the stream contained mildly polluted to clean water. This was further corroborated by the Pielou evenness index (J) which was equally high (ranging from 0.6 to 0.75) (Table 1) indicating a diverse phytoplankton community.
Effects of seasonality, land use patterns, and stream water hydrological and physicochemical characteristics on the phytoplankton population dynamics With respect to the physicochemical variables, only the mean annual temperature showed a significant spatial variation, increasing from 24.7 °C at B1 to 27 °C at B5 and then slowly decreasing to 25.8 °C at B8 (Table S2). With the exception of turbidity and orthophosphate concentration, physicochemical variables were within the WHO acceptable limit for drinking water (WHO, 2017) and did not show a significant spatial variation. However, there were significant temporal variations in the discharge, velocity, TDS, EC, and temperature ( Fig S2). The velocity ( Fig  S2a) and discharge (Table S2) increased from March to September which corresponded to the rainy season and reduced from October to February, which is the dry season having fewer rain days per month and lower precipitation. The opposite trend was observed for the EC, TDS and temperature (Fig S2d-f) which decreased from March to September and increased from October to February. The orthophosphate concentration increased from November to January, and reduced to May after which it was not detectable for the rest of the year (Fig S2c).
A correlation plot of water abiotic variables and dominant phytoplankton taxa (Fig. 3) shows strong negative relationships between velocity and discharge on the one hand and electrical conductivity and stream water discharge on the other hand. We also found significant positive correlations between suspended solids, colour, nitrite, nitrate and turbidity. These correlations were also uncovered by PCA  analysis, as two independent gradients among measured water hydrology and physicochemical characteristics (Fig. 4a). The first gradient along the first ordination axis was represented by negative correlations between TDS (correlation coefficient with the first ordination axis -1.06) and electrical conductivity-EC (-1.05) on one side and discharge (0.73) and stream velocity (0.61) on the other side. The second independent gradient (Fig. 4a) was found along the second ordination axis, which scores were positively correlated with concentrations of suspended solids (1.44), water colour (1.44), stream turbidity (1.64) and concentrations of nitrite NO 2 − (1.05) and nitrates NO 3 − (1.05). The first gradient explained 79.8% of the variation in water quality while the second gradient explained a further 9.7% giving a cumulative explained variation of 89.5%. These two independent gradients (i.e., scores from first and second ordination axis) were used for subsequent analysis.
The seven most abundant species in decreasing order were (mean cells/m 3 ± SD): Oscillatoria sp. (5164 ± 3018); Pleurosira laevis (3493 ± 4479); Synedra ulna (2224 ± 905); Spirogyra porticalis (1529 ± 1285); Diatoma vulgaris (1343 ± 1991); Navicula rhynchocephala (979 ± 855); Gyrosigma baltium (936 ± 931). The large standard deviations reflect the large differences in the phytoplankton cell densities at the different study sites. These species constituted 57.6% of the total phytoplankton abundance, and were used for further analysis to explore the influence of the stream hydrological and physicochemical characteristics, land use pattern and season on the phytoplankton dynamics. Synedra ulna and Gyrosigma baltium were negatively correlated with Synedra ulna dominating during the rainy season and Gyrosigma baltium dominating in the dry season (January to March). Oscillatoria sp. was found during the dry season, but the relationship was weak because this genus only makes a brief appearance in the peak months (January and February) and almost disappears for the rest of the year, barely maintaining a very low population.
RDA analysis showed that the phytoplankton community structure based on these dominant genera, was significantly affected by the gradient along the first ordination axis, the zone (urban, small-scale and large-scale farms) as well as the dry/rainy season (Table 2). In particular, we found marginally significant differences along the first ordination axis between dry and rainy season and statistically significant positive correlation with the PCA scores from the first ordination axis. The effect of the PCA scores from the second ordination axis were not statistically significant (P = 0.17). A second independent gradient was found along the second ordination axis from urban areas (B1 and B2), through small-scale farms (B3 and B4) to large-scale farms, i.e., CDC (B5 to B8) (Fig. 4b). These gradients were significantly linked with the main phytoplankton community structure. With the PCA scores from the first ordination axis, we found a positive correlation of the percentages of Synedra ulna (0.50) and Diatoma vulgaris (0.55) and a negative correlation of percentages of Spirogyra porticalis (-0.45) and Navicula rhynchocephala (-0.40). Percentages of genus Oscillatoria were negatively correlated with the scores from the second ordination axis (-0.64). Percentages of Gyrosigma baltium and Pleurosira laevis were negatively correlated with the PCA scores on the first ordination axis (-0.49 and -0.45 respectively) and positively (P. laevis 0.61) or negatively (G. baltium -0.82) correlated with the PCA scores from the second ordination axis (Fig. 4b).
Further detailed analyses on relationships between the three above-mentioned independent variables ( Table 2) and percentages of phytoplankton genera showed several clear patterns. Firstly, the percentages of Gyrosigma baltium and Synedra ulna clearly differed between the dry and rainy season. Percentages of genus Gyrosigma baltium were higher in the dry season compared to the rainy season (Mann-Whitney U test, U = 671.0, P = 0.014, Fig. 5a) and the opposite was true for Synedra ulna (Mann-Whitney U test, U = 423.0, P < 0.001, Fig. 5b). Second, for the genus Diatoma (Kruskal-Wallis test, H = 5.8, P = 0.056) and Pleurosira laevis (Kruskal-Wallis test, H = 16.4, P < 0.001), we found some patterns along the gradient from urban areas, through small-scale farms to largescale farms. Based on post-hoc tests, we found that percentages of genus Diatoma were marginally significantly higher (P = 0.095) in large-scale farms compared to urban areas (Fig. 5c). Other differences were not statistically significant (P at least 0.674). For the genus Pleurosira, post hoc tests uncovered that in urban areas the percentages were higher (Fig. 5d) compared to both small-scale (P = 0.019) and large-scale farms (P < 0.001). Comparison of small-scale and largescale farms was not statistically significant (P = 0.992). Finally, we found clear relationships between percentages of Gyrosigma baltium and Synedra ulna and the PCA scores from the first ordination axis. Percentages of Gyrosigma baltium were negatively correlated with the PCA scores from the first ordination axis (regression, R 2 = 0.14, F = 15.3, ß = -0.38, P < 0.001, Fig. 6a) and the opposite was true for Synedra ulna (regression, R 2 = 0.21, F = 24.7, ß = 0.46, P < 0.001, Fig. 6b). These results corroborate the correlation plot (Fig. 3) which shows the abundance of Gyrosigma baltium being negatively correlated to velocity and discharge while being positively correlated with electrical conductivity and total dissolved solids. In addition, the abundance of Gyrosigma baltium was negatively correlated with that of Synedra ulna. Table 2 Results of the RDA analysis on the effect of independent variables (PCA scores from the first ordination axis, zone and season) on the community structure of the dominant phytoplankton genera in the samples Independent variables % of explained variability pseudo-F P PCA scores from the first ordination axis 38.4 7.1 0.002 Zone (urban areas /small-scale farms/ large-scale farms) 35.4 6.9 0.002 Season (dry/rainy) 9.2 1.8 0.098 The temporal distribution of the major phytoplankton general (Fig S1) indicates that the phytoplankton community is largely dominated by Synedra ulna from April to September, which is the rainy season also corresponding to low values for TDS and EC ( Fig S2). In October and November (transition period from the rainy to the dry season), Pleurosira laevis, Diatoma vulgaris and Spirogyra porticalis became prevalent almost completely replacing Synedra ulna. However, some spatial distribution was noticed among these three genera during this period. While the Spirogyra was almost evenly distributed throughout the stream, Pleurosira was dominant upstream in the urban zone and Diatoma dominated downstream in the large-scale farming zone of the stream (Table S1). The abundance of Diatoma vulgaris drastically reduced in December and by January it gave way to Oscillatoria sp. which makes an ephemeral appearance for two months, as well as Gyrosigma baltium which dominated from January to March ( Fig  S1). Pleurosira laevis is abundant in the stream for most of the year but it seems to avoid the peak period of the rainy season (i.e., from June to September). Navicula made sporadic appearances throughout the year with no discernible pattern.

Variations in the Benoe Stream hydrological and physicochemical characteristics
This study set out to assess the effects of season, land use patterns as well as the hydrological and physicochemical characteristics of the Benoe Stream on its phytoplankton community in order to ascertain the stream water quality. Hydrological characteristics, physical and chemical properties, habitat conditions, are among the basic measures of ecological monitoring that are widely used in the assessment of water health (Zhu et al., 2021). Except for the orthophosphate concentration, all physicochemical characteristics measured were within the WHO limit for drinking water (WHO, 2017). A strong negative relationship was observed between velocity and discharge on one hand and EC and TDS on the other hand. TDS and EC give an indication of the salinity level of the stream (Rusydi, 2018). EC increases with increase in TDS and is known to be an approximation of water pollution level (Das et al., 2006). In addition, TDS and EC have been shown to be among the most important physical-chemical characteristics that explain the variability of water quality (Bertossi et al., 2013;Fetene & Teshager, 2020). We showed that TDS, EC, stream velocity and discharge explained up to 79.8% of the variation in physicochemical characteristics in the Benoe Stream. Stream velocity and discharge were mainly determined by the amount of precipitation. Therefore, the negative correlation between these characteristics and TDS and EC is most likely a result of dilution due to increased precipitation in the rainy season from mid-March to September with a peak in August in the monomodal agroecological zone in which the Benoe stream is located (Yengoh et al., 2010). This reduces the concentration of salts in the stream and consequently TDS and EC. On the other hand, TDS and EC are higher in the dry season (October to mid-March) due to increased evaporation resulting from increased temperature (WorldData.info, 2023;Yengoh et al., 2010). The suspended solids, water color, turbidity, phosphate, nitrite and nitrate concentrations were positively correlated with each other, and their combined effects explained only 9.7% of the variation in the Benoe Stream.
Ecological status of the Benoe Stream and the influence of stream characteristics, land use pattern and season on the phytoplankton dynamics Based on the high values obtained for Shannon-Wiener (H′) and Pielou's evenness (J), the Benoe Stream status ranged from healthy to mildly polluted and it contained a diverse phytoplankton community. Attempting to use all phytoplankton taxa present in a stream simultaneously to describe its ecological status is a herculean task because of the wide range of spatial scales and life-histories encompassed within these phototropic organisms (Kelly et al., 2008). Similar to Fonge et al. (2015), the phylum Bacillariophyta, commonly known as diatoms, were found to be the most abundant and diverse taxa making up 63.4% of the phytoplankton abundance and 50.6% of the identified species. The class Bacillariophyceae was the most abundant and diverse family. In contrast, Fonge et al., (2016) reported that Euglenaceae was the most abundant family, although their study was done in the mangrove. Diatoms have been described as the most species-rich group of autotrophic algae, found in fresh, brackish, and marine waters worldwide and may account for about 20% of global net primary production (Ali et al., 2018;França et al., 2017;Mann et al., 2017). This explains why diatoms are used as proxies of photosynthetic organisms in assessing the ecological status of aquatic ecosystems (Kelly et al., 2008) as well as the plethora of biotic indices that have been developed based on diatoms as bioindicators of aquatic ecological status (Celekli et al., 2021). Tropical freshwater streams are strongly influenced by seasonality as there is a well-defined wetdry climate, which has been found to affect the phytoplankton biomass and community structure (Aboim et al., 2019). Seasonal changes in the distribution and abundance of phytoplankton species have been reported by several authors, attributing them to environmental changes including water physicochemical parameters, land use and flow regime (Chauhan et al., 2015;Inyang & Wang, 2020;Qu et al., 2022). Given the relationships already established between the four stream characteristics (TDS, EC, velocity and discharge) and the rainy/dry seasons, it could be deduced that high TDS and EC and low stream velocity and discharge favoured the proliferation of Gyrosigma baltium, while Synedra ulna preferred the reverse. Variation in the abundance of Synedra sp. was not related to temperature changes. It can grow under a wide temperature range with high cell densities previously recorded between 7.2-23.9 °C (Youn et al., 2020), which encompass the temperature range found in the Benoe Stream. In addition, Synedra sp. grows well even in the absence of phosphorus but is limited by the absence of Si, which they need to synthesize their silicified cell wall (Chauhan et al., 2015). The Si concentration was not measured in the present study, but it would be interesting in future studies to find out the role of Si concentration in the abundance and distribution of the diatoms in the Benoe Stream.
Another important gradient uncovered that influenced the phytoplankton community structure was the Page 13 of 16 788 Vol.: (0123456789) urban to small-scale farming to large-scale farming (CDC) with high Pleurosira sp. abundance associated more with the urban zone, which was upstream, while Diatoma sp. was more abundant in the largescale farming zone downstream. Sometimes, land use can be more important in explaining differences in phytoplankton community than nutrients. Qu et al. (2022) reported a rise in phytoplankton abundance in a medium level residential urban area. Pleurosira sp. and Spirogyra sp. were abundant for most of the year, except during the peak period of the rainy season (June to September) when Synedra ulna was dominant.
Phytoplankton bioindicators in the Benoe Stream TDS and EC give an indication of pollution since the salts they measure come both from nature and anthropogenic activities (Rusydi, 2018). Therefore, in terms of bioindication, Synedra ulna was pollutant-sensitive since it dominated during periods of low TDS and EC while Pleurosira laevis, Spirogyra sp. and Gyrosigma baltium, were pollution-tolerant since they were negatively correlated with the PCA scores on the first ordination axis, with their abundance showing a preference for higher TDS and EC. With respect to the anthropogenic activities, percentages of Pleurosira laevis, which were positively correlated with scores from the second ordination axis prevailed more in the urban zone of the stream, while Diatoma sp. and Oscillatoria sp. percentages increased at large-scale farming zone downstream. It should be noted that Diatoma sp. and Oscillatoria sp. abundances do not peak at the same time. Diatoma sp. dominates in October and November, reducing in December, while Oscillatoria sp. takes over in January and February.

Conclusion
The phytoplankton community of the Benoe Stream was dominated by Bacillariophytes, also known as diatoms. The stream also showed high values for the Shannon-Wiener diversity (H′) and Pielou's evenness (J) indices indicating the presence of a diverse phytoplankton community, which was made up of eight phyla: Bacillariophyta, Chlorophyta, Cyanophyta, Dinophyta, Euglenophyta, Rhodophyta, Charophyta and Chrysophyta. H′ values also indicated good to mildly polluted water quality status for the Benoe Stream. Most of the physicochemical characteristics measured, including EC, TDS, DO, SS, and pH, were within WHO acceptable limits for drinking water. Multivariate statistical analysis, particularly PCA and RDA, clearly showed that the most important physicochemical parameters that influenced the phytoplankton community abundance were TDS, EC, stream water velocity and discharge. The multivariate analyses also brought out temporal variations in these stream parameters between the rainy and dry season as well as spatial variations from the urban through small-scale farming to large-scale farming zones of the stream and revealed their impact on the dominant phytoplankton genera. From these, bioindicator species were identified which can be used for the biomonitoring of the Benoe Stream and other similar streams in the monomodal agroecological zone of Cameroon. These included: Synedra ulna, Pleurosira laevis, Spirogyra sp., Gyrosigma baltium, Diatoma sp., and Oscillatoria sp.