The Liuxi River is regarded as the mother river (of high cultural significance) of Guangzhou, the capital city of Guangdong Province in southern China. The study area has a typical subtropical monsoon climate. The mean annual precipitation of the watershed is 1800 mm, mostly occurring in April–September. The Liuxi River watershed (with an area of 2300 km²) is situated in the northeastern corner of the Pearl River Delta, which has experienced rapid development in the last two decades. The watershed spans four county-level districts (Huadu, Luogang, Baiyun, and Conghua), which occupy 70% of the watershed area ( Figure 1 ). The water of the Liuxi River is used for a wide range of purposes, such as drinking water, agriculture, industry, and recreation. However, over the past several decades, the Liuxi River (especially the lower reaches) has experienced rapid agricultural, industrial, and urban development. The water quality and habitat integrity of river ecosystems are disturbed by rapid population growth and associated anthropogenic activities, leading to drastic extraneous interference from humans.

Fish specimens and water samples were collected concurrently from the headwaters to the lower reaches at the 14 sites along the Liuxi River ( Figure 1 ). Site information (e.g., mainstream and tributary locations) and environmental data, including habitat characteristics and physicochemical parameters of water quality, were provided in Supplementary Table S1 . Forests are the dominant land use type in the watershed, accounting for 55% of the total area. Approximately 33% of the watershed is used for agriculture, including orchards (20%), paddy rice fields (9%) and vegetable growing areas (4%). The remaining watershed is occupied by built-up areas (9%) and water (3%). The development in the watershed is spatially unbalanced. The upper reaches (sites #1–#4) are covered by dense forests, while agricultural activities are concentrated in the middle reaches (sites #5–#10). The area close to the lower reaches (sites #11–#14) is highly urbanized; in particular, site 13 was located in industrialized area where a large amount of sewage was discharged into the river.

Fish samples were collected in 2022 from the headwaters to the lower reaches (sites #1–#14 in Figure 1 ) along the subtropical Liuxi River during the rainy (June to July) and dry (December) seasons. Each site was sampled two times over a season following the basic guidelines of Barbour et al. (1999) and Hauer and Lamberti (2007). Electrofishing equipment was used to effectively stun and collect fish (individual weight< 10 kg) in a 2 m wide × 2 m long × 3.5 m deep water column. Due to varying water levels, two electrofishing operations were conducted as follows: 1) At wadeable sites, single-pass backpack electrofishing was performed simultaneously by two operators moving in a zig-zag fashion. Electrofishing equipment was adjusted at low voltage and mixed frequency, and the walking speed was controlled to ensure a sampling effort of approximately 8 m² min^-1 over 30 min. 2) At nonwadeable sites, a 6-m-long welded hull boat was used for boat-electrofishing, and a bamboo quant was used to propel the boat to eliminate noise disturbance to fish. Electrofishing equipment was adjusted at a high voltage and main frequency, and the paddling speed was controlled to ensure a sampling effort of approximately 6 m² min^-1. Boat-electrofishing was conducted over a distance of 500 m, spanning both riverbanks at a depth of 1–3 m (Flotemersch et al., 2006). With the help of local fishermen, gill nets (mesh sizes of 10, 20, 30, 40 mm between adjacent nodes, high × long areas of 1.0 m× 60 m, 1.5 m × 100 m, 2.5 m × 150 m, 3.5 m × 200 m) and hoop nets (mesh sizes of 5, 10, 15 mm, volume of 0.25 m high × 0.35 m wide × 10 m long) were used as passive methods to supplement the requisite specimen in cases where the latter was precluded by high depth and large fish sizes.

All statistical analyses were conducted in R 4.0.5 with the primary packages cluster, factoextra, phyloseq, vegan, microeco, PerformanceAnalytics, and ggplot2. Cluster analysis for grouping spatial zones and scoring items relied on the Bray−Curtis dissimilarities of fish abundance (individual number) and indicators’ ecological properties, respectively. Statistically significant cluster groupings were identified using a bootstrap randomization technique in which the nonzero values were resampled and used to generate pseudovalues of Bray−Curtis dissimilarities under the null hypothesis. A frequency distribution of pseudovalues was generated from 1000 randomizations of the data matrix, and the 95th percentile was used as the critical value to determine between-group significance. Stepwise forward selection was performed to reduce the linearly correlated environmental factors, and the variables showed significant differences in their values.

Linear discriminant analysis effect size (LEfSe) is an algorithm for high-dimensional indicator discovery that identifies taxa by characterizing the differences between two or more biological conditions (Segata et al., 2011). LEfSe emphasizes both statistical significance and biological relevance, allowing researchers to identify discriminative features that are significantly different among biological classes. The nonparametric factorial Kruskal−Wallis sum-rank test is first used to detect features with significant differential abundance with respect to the class of interest. Second, LEfSe uses linear discriminant analysis to estimate the effect size of each differentially abundant feature and rank the feature accordingly (Liu et al., 2021). Spearman’s correlation analysis was used to measure the strength and direction of monotonic association between fish indicator scores and traditional biodiversity indices as well as environmental variables.

One hundred thirteen species belonging to 10 orders, 32 families, and 83 genera were found along the Liuxi River, including 97 native species and 11 alien species (see details in Supplementary Table S2 ). Seven spatial zones (i.e., zones I to VII shown in Figure 2A ) were grouped by clustering analysis based on the individual number of fish assemblages. Sites 1−2 in zone I were located in the headwaters ( Figure 2B ), where fish assemblages were dominated by Cyprinidae 62.7% (mainly composed of Danioninae 30.2% and Barbinae 24.8%), Gastromyzontidae (18.0%), and Gobiidae (11.6%). Sites 3−4 in zone II were located in the upper mainstream, where fish assemblages were dominated by Cyprinidae 68.6% (mainly composed of Cyprininae 14.2%, Rhodeinae 13.9%, Gobioninae 12.8%, and Culterinae 12.4%), Cichlidae 15.7%, and Cobitidae 9.1%. Sites 6−7 in zone III were located in the midstream tributaries, where fish assemblages were dominated by Cyprinidae 75.2% (mainly composed of Gobioninae 32.7% and Danioninae 31.4%), Cobitidae (8.1%), and Gobiidae (8.0%). Zone IV included mainstream site 5 and tributary sites 9−10, where fish assemblages were dominated by Cyprinidae 68.1% (mainly composed of Danioninae 28.7%, Gobioninae 11.0%, and Rhodeinae 8.9%), Gobiidae (13.7%), and Cobitidae (7.0%). Sites 8 and 11 in Zone V were located in the middle-lower mainstream, where fish assemblages were dominated by Cyprinidae 61.9% (mainly composed of Culterinae 23.0%, Labeoninae 8.2%, Cyprininae 6.5%, and Xenocyprininae 6.4%) and Cichlidae 18.4%. Sites 12 and 14 in Zone VI were located in the lower mainstream, where brackish fishes appeared and fish assemblages were dominated by Cyprinidae 54.5% (mainly composed of Culterinae 25.3% and Labeoninae 22.1%), Cichlidae 17.2%, Loricariidae 6.7%, and Engraulidae 5.8%. Site 13 in Zone VII was located in the downstream polluted tributary, where fish assemblages were dominated by Cichlidae 52.6%, Loricariidae 26.3%, and Clariidae 10.5%.

The differences in spatial distribution of fishes in the Liuxi River are mainly reflected at the genus and species levels ( Figure 3 ), suggesting that genus- and species-level indicators were effective in distinguishing fish assemblages in each zone. The difference at the family level was lower than that at the genus level because Cyprinidae contained more species than other families. LEfSe analysis results showed that the fish assemblages in Zone I were marked by Balitoridae, those in Zone III were marked by Cobitidae, those in Zone IV were marked by Cyprinidae, those in Zone V were marked by Mastacembelidae, those in Zone VI were marked by Engraulidae, and those in Zone VII were marked by Loricariidae. There were no family-level fish indicators in Zone II. At the genus and species levels, a total of 4 genera and 30 species were selected by LEfSe, including Pseudogastromyzon changtingensis, Vanmanenia (V. pingchowensis and V. gymnetrus), Acrossocheilus parallens, and Zacco platypus in Zone I; Rhodeus sinensis, Rasbora steineri, Hemiculter leucisculus, Acheilognathus barbatulus, Cobitis (C. sinensis and C. arenae), Channa maculata, and Pelteobagrus (P. fulvidraco, P. intermedius, and P. vachelli) in Zone II; Pseudogobio vaillanti, Microphysogobio (M. kiatingensis and M. fukiensi), Sarcocheilichthys nigripinnis, Opsariichthys bidens, Abbottina rivularis, Micronemacheilus pulcher, Rhinogobius giurinus, and Glyptothorax fokiensis in Zone III; Cirrhinus molitorella, Squalidus argentatus, and Hemibagrus guttatus in Zone IV; Culter alburnus, Parabramis pekinensis, Hemibarbus medius, Xenocypris davidi, Schistura fasciolatus, Mastacembelus armatus, and Silurus asotus in Zone V; Hypophthalmichthys molitrix, Ctenopharyngodon idella, and Coilia grayi in Zone VI; and Coptodon zillii and Hypostomus plecostomus in Zone VII.

The core design of our fish-based index is to use indicators’ accumulative positive feedback to evaluate river ecological status and food web attributes (e.g., trophic structure, transfer efficiency, and network connectance). To assess the positive feedback of each fish indicator to local environment (i.e., environmentally indicative function), 6 focuses were set up for the evaluation criteria (see details in Supplementary Table S3 ): 1) prioritize rare, endangered, and protected species; 2) rheophilic species with specific feeding and spawning/nursing requirements for habitats (e.g., pristine floodplains and riparian zones), substrate types (e.g., pebble, grave, and sand), and flow regimes (e.g., riffles, pools, and runs); 3) migratory species that can indicate both lateral and longitudinal river connectivity; 4) keystone species occupying different trophic levels, which reflect the integrity and diversity of local food web; 5) top predators at high trophic levels since most of them have a long life history and are exposed to multiple stressors (e.g., water pollution, habitat degradation, and overfishing); and (6) alien species with strong tolerance (e.g., Cichlidae and Loricariidae in zone VII) are assigned a score of 0 due to their negative impacts on native aquatic communities. On this basis, the final evaluation score of each fish indicator was calculated as positive feedback level × trophic level, where positive feedback level = feeding score + spawning/nursing score +migrating score + keystone store + endangered/protected score (see scoring diagram in Figure 4 ).

To clarify the environmentally indicative function, thirty-four fish indicators were clustered into 19 ecological patterns ( Figure 5 ) according to their scores for environmental evaluation. Fish indicators with evaluation scores of 5 – 8 have two characteristics: 1) they are widely distributed and able to swim across lotic and lentic water bodies (e.g., from river to reservoir), such as Culter recurviceps and Z. platypus, and 2) they have a certain tolerance to the degradation of water and habitat, such as Rhodeus sinensis and Rhinogobius giurinus ( Table 2 ). Fish indicators with scores > 10 were benthic top predators at high trophic levels (e.g., S. asotus, C. maculata, and Pelteobagrus), the appearance of which could indicate a long food chain and integrated food web. In particular, the native protected species (e.g., H. guttatus) and endangered species on the red list (e.g., Anguilla japonica) are not only species protected by the state but also important migratory species that indicate ‘stream – lake’ or ‘river − estuary’ connectivity; thus, they obtained the highest score of 15 – 19. Fish indicators with scores< 3 were those at low trophic levels and indicated short food chains, such as filter-feeding H. molitrix and scrape-feeding C. molitorella. C. zillii and H. plecostomus are dominant alien species in the tropics and subtropics of China (Wang et al., 2019), which compete with native species for limited resources. Given that fish species belonging to the same genus or subfamily commonly have analogous ecological characteristics and exhibit consistency in environmentally indicative functions, it is necessary to provide alternative species as candidates for the 34 target indicators, which not only enriches the evaluation system, but also facilitates field monitoring efficiency (see eight alternative species in Table 2 ).

According to the principle of ‘simple, quick, precise, and practical’, we proposed an index of fish indicator (IFI) to evaluate local environment based on the abovementioned scoring criteria. The evaluation on local environment was based on the accumulation of each fish indicator’s positive feedback to local environment. Thus, the higher the final indicator-based assessment score, the better the local environmental quality. Equations (1–3) were developed and compared to optimize the evaluation effects, including a qualitative method using indicators’ species abundance and two quantitative methods weighted by indicators’ individual number and biomass:

At a sampling site, IFI is the fish indicator-based environmental evaluation score. E_i and n_i are the assigned ecological score (see Table 2 and Supplementary Table S3 ) and individual number of a single fish indicator i (i = 1 ~ m), respectively. N is the total individual number of sampled fish, including indicator and non-indicator species. IFI_Qis a qualitative index that equals the sum of the E_i value of each fish indicator species. IFI_N and IFI_B were two quantitative indexes, which were weighted by the relative contribution of fish indicators’ individual number (i.e., n i N ) and by the relative contribution of fish indicators’ individual biomass (i.e., b i B ), respectively. For example, if the number of fish indicators appearing at the sampling site is m, IFI_N = E 1 × n 1 N + E 2 × n 2 N +… + E m × n m N . After calculating the IFI, the criteria for evaluating local environment from a perspective of food web properties, e.g., predator-prey linkages, food chain length, and energy flows in trophic network, were explored in Table 1 .

Compared with the qualitative IFI_Q , IFI_N and IFI_B weighted by individual number and biomass could more notably distinguished the ecological status between upstream and downstream as well as between tributary and mainstream ( Figure 6 ). Generally, the environmental evaluation results showed higher scores in upper stream than in lower river and in tributaries than in mainstream. The highest evaluation scores of IFI_N and IFI_B were recorded at tributary sites #1–#2 and #6, whereas the highest evaluation scores of IFI_Q were recorded at mainstream sites #3–#5. Evaluation score at site #13 was 0 since there was no indicator species with positive feedback. Spearman correlation coefficients between IFI and traditional ecological indices (e.g., Shannon‒Wiener diversity of fish community and benthic index of biotic integrity) as well as environmental factors were calculated to evaluate the practical application effect of IFI ( Figure 7 ). The results showed that IFI_Q , IFI_N and IFI_Bexhibited significantly (p< 0.05) positive correlations with elevation (m), dissolved oxygen (mg/L), flow velocity (cm/s), Shannon‒Wiener diversity, and benthic index of biotic integrity and significantly negative correlations with total nitrogen (mg/L), NH₃-N (mg/L), and electronic conductivity (μS/cm). Generally, IFI calculated by a concise scoring formula could comprehensively reflect the relationships between biological properties and environmental conditions. Among the three indices, IFI_B exhibited the best performance — significantly correlated with most biotic and environmental factors. Compared with IFI_Q and IFI_N , the greatest advantage of IFI_B is its reasonable accounts for the composition of fish communities, especially for carnivorous and herbivorous fishes with large body sizes.

In the longitudinal gradient, the 34 fish indicators selected along the Liuxi River could reflect not only the aquatic fauna diversity and community structure but also the basic structure of river food web (e.g., predator-prey linkages and trophic energy flows). Moreover, the indicator-based evaluation scores could notably distinguished the differences in ecological status along the upper, middle, and lower reaches as well as the differences between the mainstream and tributaries (Petry and Schulz, 2006). For example, fish in zones I and III (e.g., Danioninae and Barbinae) indicated wadable mountain headwaters that had higher elevation and lower water temperature than zones located in downstream plains (Wang et al., 2018b). In a subtropical monsoon climate, perennial rainfall, fast flowing velocity, high dissolved oxygen, and scattered riffle and pebble substrates provide favorable conditions for benthic algae and invertebrates (Dudgeon, 2008). These food resources are the main energy supply for periphytivorous and insectivorous food chains, which could be indicated by rheophilic Vanmanenia, O. bidens, Z. platypus, and M. pulcher.

Rheophilic species in small streams are vulnerable to human disturbance, especially for those in rural areas with developed tourism (e.g., agritainment) and intensive livestock breeding (e.g., duck industry) (Wang et al., 2020a). For example, although the water quality at sites #7 and #9 were as good as headwaters #1 and #2, the original pebble-strewn substrates with coarse sand in midstream tributaries are illegally extracted or buried during agricultural land expansions (Lasne et al., 2007). The reduction in pebble coverage area led to the homogenization of the substrate type, leading to the loss of Nemacheilidae (Micronemacheilus and Schistura) and Botiidae feeding on aquatic insects and Gobioninae (Hemibarbus and Pseudogobio) hidden under sand (Wang et al., 2018a, b). Vanmanenia and Rhinogobius, which crawl on pebble surfaces, were also threatened. Moreover, because Rhinogobius are important food resources for carnivorous fish, their population decrease had negative impacts on top predators (e.g., S. asotus and C. maculata), thus affecting the whole food web through trophic cascading (Wang et al., 2019). This is the reason why evaluation scores for sites #7 and #9 were much lower than for sites #1 and #2.

At mainstream sites #11 and #12 with low evaluation scores, the subfluvial pebbles were removed or buried due to channel reconstruction, and some were piled up on the bank for the manufactured landscape. In addition, tributaries in the lower urban reach have been subjected to severe anthropogenic modification (Wang et al., 2021). The reinforced bank led to a separation between the riparian zone and river channel, which destroyed the aquatic vegetation that was necessary for fish reproduction (Sosa-López et al., 2005). In particular, due to the discharge of domestic sewage, extremely high COD_Mn and low DO were recorded at tributary site #13 in zone VII, where only invasive alien species (e.g., C. zillii, H. plecostomus, and Pangasius sutchi) resistant to pollution survived, with no native fish indicators observed. The spawning process of H. guttatus and A. japonica requires a certain migration distance, which is totally dependent on hydrologic and hydrodynamic connectivity (Petry and Schulz, 2006). However, water conservancy projects destroyed the original fluvial morphology and obstructed the migratory route, leading to a decline in migratory fish stocks (Southerland et al., 2007).

The purpose of constructing this scoring system is to overcome some deficiencies of traditional ecological indices, such as Shannon‒Wiener diversity and the fish index of biotic integrity (Souza and Vianna, 2020). In terms of calculation, traditional indices often underestimate the environmental conditions in upstream tributaries, where the diversity of the fish community is low but water quality and habitat integrity are quite high (Carignan and Villard, 2002). For example, in headwater zone I, fish species richness was limited by harsh physical conditions (e.g., high velocity and shallow water); however, local environment remained pristine since there was little external interference (Wang et al., 2018b). Thus, the environmental evaluation results for zone I should be high, which used to be underestimated by traditional ecological indices (Bal et al., 2018). The advantage of our scoring and evaluation system is to imbue ecological significance (e.g., fishes’ living requirements) to mathematical statistics, through which abiotic conditions can be indicated by biotic feedback (Wang et al., 2023). Such a scoring method can be modified to meet the evaluation needs of various rivers and in different scales (e.g., regional, provincial, or national).

The successive presence of upstream fish indicators feeding on aquatic insects, periphyton, and gastropods; midstream indicators feeding on hydrophytes, bivalves, odonate larvae, shrimp, and fish; and downstream indicators feeding on plankton, polychaetes, crabs, and seaweeds are three key processes accounting for the spatial distribution of fish assemblages (Wang et al., 2018a, b). It could be expected that the cumulative addition of indicators that had already appeared in the tributaries of upper reaches would lead to peak richness of indicator species in the mainstreams of mid-lower reaches (Welcomme et al., 2006; Eick and Thiel, 2014). This is the reason why IFI_Q -based evaluation scores reached the highest at mainstream sites #3–#5. However, IFI_Q failed to reflect the relative importance of fish indicators in community composition and could not avoid the bias caused by sampling occasionality. For example, since invasive cichlids and loricariids were dominant in most tropical and subtropical rivers of China, the evaluation score would be overestimated if only the species richness rather than the community composition was considered. Our scoring method intended to provide a realistic judgement on ecological status by using IFI_N and IFI_B . These two quantitative indices weighted by fish indicators’ relative contribution performed better than qualitative IFI_Q in revealing the practical conditions of the local environment, which provides evidence for managers to make decisions on river conservation (Mujiyanto et al., 2021).

The spatial distribution and population dynamics of fish with specific feeding and reproductive requirements carry the comprehensive information about water quality (e.g., flow velocity, dissolved oxygen, and nutrient contents), substrate types (e.g., sand, gravel, and pebble), and habitat characteristics (e.g., riffle, pool, and run), which are critical environmental factors determining the fluvial biodiversity (Lasne et al., 2007; Souza and Vianna, 2020). In addition, the fish indicators selected in this study occupied different trophic levels and represented important trophic nodes within predator‒prey relationships (Wang et al., 2018b), which can indicate the integrity of food chains occurring in heterogeneous habitats. Thus, the contribution (e.g., individual number and biomass) of these fish indicators to community composition has great impacts on the structure and energy flows of local food web (Bal et al., 2018), which is important for maintaining health and functioning of aquatic ecosystem.

Notably, although endangered or national protected fishes are can be indicators of multiple functions in river environments, it is difficult for filed investigators to gather their information due to their limited distribution, low population density, and legal protection (Southerland et al., 2007). Therefore, we recommend that the monitoring of protected species be achieved in two ways: 1) substitution by using other species with similar environmentally indicative functions, e.g., using Hemibagrus macropterus to indicate habitat areas where H. guttatus (national second-class protected animals) is likely to occur, and 2) using environmental DNA technology to conduct long-term monitoring in areas where protected species may occur (Wang et al., 2023). The scoring and evaluation system constructed in this study are suitable not only for traditional ecological data (e.g., individual number and biomass) but also for operational taxonomic unit abundance obtained by high-throughput sequencing.