Dataset DOI: 10.5061/dryad.1ns1rn96w

Description of the data and file structure

This dataset contains the spatial data, processed ecological data, statistical outputs, and R scripts used to generate the main and supplementary figures associated with the study Riparian land use reduces multitrophic community stability in rivers.

The study examined how riparian land use, particularly impervious surface and cropland cover, was associated with the taxonomic richness, community composition, temporal dynamics, and stability of six aquatic organism groups in the Dongjiang River, southern China. The organism groups included fish, aquatic insects, protists, algae, cyanobacteria, and bacteria. Community data were obtained using environmental DNA metabarcoding.

The dataset is organized according to figure number. The archive for Figure 1 contains ArcGIS-compatible spatial files used to generate the study-area map. The archives for Figures 2–6 contain the R scripts, processed data tables, and statistical results required to reproduce the corresponding analyses and figures. Supplementary figures are included with the most closely related main figure.

Raw sequencing data are not included in this Dryad dataset. They are publicly available from the NCBI Sequence Read Archive under BioProject accession PRJNA1282341.

Files and variables

File: Figure_1.zip

Description: This archive contains the ArcGIS-compatible spatial data used to generate Figure 1A, including the Dongjiang River network, sampling-site locations, elevation information, geographic regions, and land-use classifications.

Figure 1A presents the location of the Dongjiang River in southern mainland China and the distribution of the sampling sites. Digital elevation model data were used to extract river-network structure and elevation. Land-use types were interpreted from Landsat-8 Operational Land Imager imagery at a spatial resolution of 30 m. The land-use layer represents average land-cover conditions from 2021 to 2023.

The map illustrates the longitudinal transition from relatively forested and minimally disturbed upstream reaches, through regulated or meandering tributaries such as the Lijiang River, to wider and increasingly urbanized downstream reaches.

Figure 1B presents the conceptual workflow and analytical framework of the study. It summarizes the links evaluated among riparian land use, community composition, taxonomic richness, temporal stability, and multitrophic community stability.

File: Figure_2.zip

Description: This archive contains the taxonomic records used to construct the phylogenetic trees and identify dominant families for Figure 2.

Figure 2 presents phylogenetic trees for the six riverine organism groups detected using environmental DNA metabarcoding. Taxonomic records were projected onto an a priori taxonomic framework. Taxa that did not have valid identification numbers in the NCBI Taxonomy Database or Genome Taxonomy Database were not included in the phylogenetic trees.

Text and shading colors distinguish the six organism groups. Symbols identify the two or three families with the highest relative sequence abundance in each group. Temporal changes in relative sequence abundance are shown relative to the initial survey conducted in November 2021. Red arrows indicate increases and blue arrows indicate decreases relative to this initial survey.

Phylogenetic relationships were inferred using PhyloT v2. Trees were exported in Newick format and visualized and annotated using Interactive Tree Of Life, iTOL v6.

File: Figure_3.zip

Description: This archive contains the statistical outputs and R scripts used to generate Figure 3 and the associated supplementary figures. The workbook contains six worksheets representing different taxonomic groups. Each row represents a sample record within the corresponding group, and each column contains either sample metadata or a taxon-specific numerical value.

The sp column names are indexed taxon identifiers rather than scientific names. Interpret each identifier together with its worksheet name; matching column labels across worksheets should not be assumed to represent the same taxon.

Figure 3 describes spatiotemporal variation in the community structure of fish, aquatic insects, protists, algae, cyanobacteria, and bacteria.

Community dissimilarity was calculated from presence–absence environmental DNA data using the Jaccard dissimilarity index. Principal coordinate analysis was used to reduce the dissimilarity matrices to site-level ordination scores. Linear discriminant analysis was then applied to the principal-coordinate scores to display differences among sampling times and geographic regions.

Colors represent sampling times, whereas plotting symbols represent geographic regions.

A significant PERMANOVA result indicates that community composition differs among years, regions, or year-by-region combinations. Because PERMANOVA may also be affected by unequal within-group dispersion, the corresponding betadisper() results are provided to help distinguish differences in group centroids from differences in multivariate variability.

File: Figure_4_Figure_S7_S9-11.zip

Description: This archive contains the generalized additive model outputs used to generate Figure 4 and the associated supplementary figures.

Figure 4A shows temporal trajectories of taxonomic richness for the six organism groups. Richness values were standardized within organism groups using z-score transformation and modeled separately using generalized additive models. Sampling time was fitted as a smooth predictor.

The lower x-axis reports the number of months elapsed since the initial survey in November 2021.

Figure 4B shows longitudinal variation in standardized community stability along the river system. Distance to the estuary was used as the smooth spatial predictor. Colored curves represent organism-group-specific relationships, whereas the black curve represents the multitrophic or multi-taxa relationship.

Adjusted R² and P values are reported for each fitted model.

File: Figure_5_Figure_S12-15.zip

Description: This archive contains the Mantel-test results and R scripts used to generate Figure 5 and the associated supplementary figures. 

The accompanying file, Land.csv, contains paired sampling-site identifiers and land-use variables, including cropland, impervious surfaces, and human land use at multiple buffer distances. Each row identifies a pair of sampling sites through the Site1 and Site2 columns. Land-use percentages are stored as decimal fractions ranging from 0 to 1; for example, a value of 0.30 represents 30%. Column names are case-sensitive and correspond exactly to those in Land.csv. In buffer-specific column names, Crop denotes cropland, Imper denotes impervious surfaces, and HLU denotes human land use. The prefix specifies the buffer distance: 250 and 500 indicate metres, whereas 1km, 5km, 10km, and 25km indicate kilometres.

The file SiRich.csv contains pairwise differences in richness for six aquatic biological groups: fish, insects, protists, algae, cyanobacteria, and bacteria. Each row represents a pair of sampling sites, identified by Site1 and Site2. The symbol Δ denotes the difference in richness between the two sampling sites. The letters F, I, P, A, C, and B identify the six biological groups.

Figure 5 evaluates relationships between riparian land use and changes in taxonomic richness or community stability.

Figures 5A and 5B show relationships between impervious-surface or cropland cover and taxonomic richness. Figures 5C and 5D show relationships between impervious-surface or cropland cover and community-stability metrics.

Land-use gradients were calculated as differences in land-cover values between pairs of sampling sites. Changes in richness and stability were also represented as pairwise differences between sites.

Mantel tests with 9,999 permutations were used to evaluate associations between the land-use difference matrices and the corresponding ecological difference matrices. The fitted lines and grey shaded areas are visual summaries of the relationships and their 95% confidence intervals. The Mantel correlation coefficient and permutation-based P value are reported in each panel.

The sign of a pairwise difference depends on the order in which the two sites are entered. A positive value indicates that the value at Site_1 is greater than the value at Site_2, whereas a negative value indicates the opposite. Users should retain the original site order when reproducing the analyses.

File: Figure_6_Figure_S16.zip

Description: This archive contains the processed richness, stability, land-use, and structural equation model data used to generate Figure 6. Data_for_Figure_6B.xlsx contains the processed piecewise structural equation model results supporting Figure 6B. The worksheet "Data for SEMs" includes model information, component-model path coefficients, coefficients of determination (R-squared), and a global model-fit test. Blank cells indicate that the corresponding field does not apply to that particular result type or that the statistic was not provided in this summary table. They should be treated as not available (NA), never as zero. Asterisks (*, **, ***) in the Significance column are significance annotations carried over from the source table, not missing or hidden numerical values.

Figure 6A evaluates relationships between taxonomic richness and community stability. Pairwise differences in richness were compared with pairwise differences in the corresponding community-stability metric. Mantel tests with 9,999 permutations were used to quantify these associations.

Figure 6B presents a structural equation model describing the statistical pathways linking impervious-surface cover and cropland cover to community stability. The model includes direct associations between land-use variables and stability and indirect associations mediated through ecological variables such as richness and the component processes contributing to stability.

Numbers adjacent to arrows are standardized path coefficients. R² values indicate the proportion of variance in each endogenous response variable explained by its predictors.

Because the study is observational, the structural equation model represents hypothesized statistical pathways and should not be interpreted as definitive experimental evidence of causation.

Code/software

The R scripts can be run using R and RStudio. Required R packages are listed at the beginning of each script or loaded within the relevant code sections.

The analyses include functions for:

community-data processing;

Jaccard dissimilarity calculation;

principal coordinate analysis;

linear discriminant analysis;

PERMANOVA;

multivariate-dispersion testing;

generalized additive models;

Mantel tests;

structural equation modeling;

statistical visualization.

Major R packages may include:

vegan, for community dissimilarities, PERMANOVA, dispersion tests, and Mantel tests;

mgcv, for generalized additive models;

ggplot2, for figure production;

packages used for data manipulation and spreadsheet import;

the structural equation modeling package specified in the Figure 6 script.

Package requirements and object names should be retained as provided in the scripts.

Phylogenetic relationships for Figure 2 were inferred using phyloT v2. Trees were exported in Newick format and visualized using Interactive Tree Of Life, iTOL v6.

ArcGIS-compatible files for Figure 1 can be opened in ArcGIS Pro, ArcMap, QGIS, or equivalent geographic information system software.

Access information

Other publicly accessible locations of the data:

• The raw data that support the findings of this study are openly available in the Sequence Read Archive of the NCBI database under the accession project number PRJNA1282341.

Contact information

Questions regarding the dataset should be directed to the corresponding author identified in the associated publication.

Readers who require additional data, clarification of variable definitions, assistance with reuse, or other information related to the dataset are encouraged to contact the corresponding author by email (Feilong Li, fl_li@gdut.edu.cn).