1 Introduction

King County’s Integrated Streams Program implemented a pilot study to evaluate whether periphyton can serve as a sensitive biological indicator of nutrient enrichment and stream health. Periphyton responds quickly to environmental change and may complement existing water quality and benthic invertebrate monitoring by providing earlier signals of nutrient- and development-related impacts.

Twelve ISP monitoring sites were selected to span a gradient of watershed development, from forested/low-development basins to highly urbanized streams. These sites represent first- through fifth-order streams across WRIAs 7, 8, 9, and 10 and were chosen to contrast conditions at the low and high ends of the County’s development spectrum.

The study’s objectives were to (1) compare traditional microscopy and eDNA metabarcoding approaches for characterizing periphyton communities, (2) establish baseline periphyton conditions across contrasting watershed types, and (3) assess the feasibility and value of incorporating periphyton monitoring into the County’s long-term stream assessment program.

To view the full SAP, visit: https://your.kingcounty.gov/dnrp/library/2025/kcr4156/kcr4156.pdf

2 Study Sites and Environmental Context

2.1 Site characteristics

SiteInfo %>%
dplyr::select(Stream, Grouping, BIBI_Locator, BIBI_Score_2024) %>%  # adjust column names as needed
arrange(Grouping, Stream) %>%
knitr::kable(digits = 1, caption = "Summary of site characteristics.")

Summary of site characteristics.
Stream	Grouping	BIBI_Locator	BIBI_Score_2024
Carey	Rural	08CED_bug2	96.1
Carey_FR	Rural	08CED_bug2	96.1
Cherry	Rural	07CHR_bug1	96.1
Griffin	Rural	07GRN077319	75.8
Judd	Rural	28A	59.0
Ravensdale	Rural	09COV1756	71.2
Rock	Rural	08CED4192	83.8
Juanita	Urban	08EAS_bug4	19.8
Longfellow	Urban	LF04	1.7
Longfellow_FR	Urban	LF04	1.7
McAleer	Urban	08WES0903	27.4
Miller	Urban	Miller_Club	31.9
Thornton	Urban	TMMA3168	10.7
Venema	Urban	PIVE1895	13.7

Rural sites had lower impervious cover and tended to have higher B-IBI scores.

2.2 Nutrients

Nutrient concentrations were summarized across the summer index period (Jul-Oct).

WQ_graph2

Nutrient concentrations were consistently higher in urban streams, often by an order of magnitude, relative to rural streams.

FinalNutrientTable <- FinalNutrientTest %>% 
  transmute(
    Parameter,
    Test = dplyr::recode(test_type,
                         "t" = "t-test",
                         "wilcox" = "Wilcoxon rank-sum"),
    `Mean (Rural)`  = mean_group1,
    `SD (Rural)`    = sd_group1,
    `n (Rural)`     = n_group1,
    `Mean (Urban)`  = mean_group2,
    `SD (Urban)`    = sd_group2,
    `n (Urban)`     = n_group2,
    `Test statistic` = statistic,
    `Effect size`    = effect_size,
    `p-value`        = p_value,
    Sig = significance
  )

FinalNutrientTable %>% 
  mutate(
    across(c(`Mean (Rural)`:`Effect size`), ~round(., 2)),
    `p-value` = signif(`p-value`, 3)
  ) %>% 
  kable(
    caption = "Comparison of nutrient concentrations between rural and urban streams.",
    align = "lcccccccccc",
    booktabs = TRUE
  ) %>% 
  kable_styling(full_width = FALSE, position = "center", font_size = 10)

Comparison of nutrient concentrations between rural and urban streams.
Parameter	Test	Mean (Rural)	SD (Rural)	n (Rural)	Mean (Urban)	SD (Urban)	n (Urban)	Test statistic	Effect size	p-value
Total Nitrogen	Wilcoxon rank-sum	2.57	0.29	25	2.99	0.15	27	47.00	1.82	1.00e-07
Nitrite + Nitrate Nitrogen	Wilcoxon rank-sum	2.32	0.53	25	2.89	0.24	27	96.00	1.68	1.02e-05
Ammonia Nitrogen	Wilcoxon rank-sum	0.71	0.16	25	1.17	0.30	27	71.00	1.75	1.10e-06
Total Phosphorus	t-test	1.22	0.38	25	1.87	0.14	27	-8.04	2.30	0.00e+00
Orthophosphate Phosphorus	t-test	0.93	0.41	25	1.75	0.13	27	-9.63	2.76	0.00e+00

2.3 Stream Width and Canopy Cover

StreamWidth_graph + Canopy_graph

FinalCanopyTable <- FinalCanopyTest %>% 
  transmute(
    Test = dplyr::recode(test_type,
                         "t" = "t-test",
                         "wilcox" = "Wilcoxon rank-sum"),
    `Mean (Rural)`  = mean_group1,
    `SD (Rural)`    = sd_group1,
    `n (Rural)`     = n_group1,
    `Mean (Urban)`  = mean_group2,
    `SD (Urban)`    = sd_group2,
    `n (Urban)`     = n_group2,
    `Test statistic` = statistic,
    `Effect size`    = effect_size,
    `p-value`        = p_value,
    Sig = significance
  )

FinalCanopyTable %>% 
  mutate(
    across(c(`Mean (Rural)`:`Effect size`), ~round(., 2)),
    `p-value` = signif(`p-value`, 3)
  ) %>% 
  kable(
    caption = "Comparison of stream width and canopy cover between rural and urban streams.",
    align = "lcccccccccc",
    booktabs = TRUE
  ) %>% 
  kable_styling(full_width = FALSE, position = "center", font_size = 10)

Comparison of stream width and canopy cover between rural and urban streams.
Parameter	Test	Mean (Rural)	SD (Rural)	n (Rural)	Mean (Urban)	SD (Urban)	n (Urban)	Test statistic	Effect size	p-value	Sig
NA	t-test	6.60	3.85	6	4.10	1.55	6	1.47	-0.85	0.1870
NA	t-test	82.38	12.42	6	95.22	3.07	6	-2.46	1.42	0.0521

Rural streams were wider and generally had less riparian cover than urban streams. However, the results were not statistically significant between groups.

3 Periphyton Biomass

Biomass was quantified by measuring ash free dry weight and chlorophyll a. Results were normalized by the area sampled.

BiomassPlot

FinalBiomassTable <- FinalBiomassTest %>% 
  transmute(
    Parameter,
    Test = dplyr::recode(test_type,
                         "t" = "t-test",
                         "wilcox" = "Wilcoxon rank-sum"),
    `Mean (Rural)`  = mean_group1,
    `SD (Rural)`    = sd_group1,
    `n (Rural)`     = n_group1,
    `Mean (Urban)`  = mean_group2,
    `SD (Urban)`    = sd_group2,
    `n (Urban)`     = n_group2,
    `Test statistic` = statistic,
    `Effect size`    = effect_size,
    `p-value`        = p_value,
    Sig = significance
  )

FinalBiomassTable %>% 
  mutate(
    across(c(`Mean (Rural)`:`Effect size`), ~round(., 2)),
    `p-value` = signif(`p-value`, 3)
  ) %>% 
  kable(
    caption = "Comparison of stream width and canopy cover between rural and urban streams.",
    align = "lcccccccccc",
    booktabs = TRUE
  ) %>% 
  kable_styling(full_width = FALSE, position = "center", font_size = 10)

Comparison of stream width and canopy cover between rural and urban streams.
Parameter	Test	Mean (Rural)	SD (Rural)	n (Rural)	Mean (Urban)	SD (Urban)	n (Urban)	Test statistic	Effect size	p-value	Sig
Ash Free Dry Weight	t-test	3.74	0.42	6	3.40	0.14	6	1.89	-1.09	0.10700
Chlorophyll a	t-test	1.54	0.48	6	0.55	0.44	6	3.76	-2.17	0.00378

4 Qualitative Assessment

Qual_graph

FinalQualitativeTable <- FinalQualitativeTest %>% 
  transmute(
    Parameter,
    Test = dplyr::recode(test_type,
                         "t" = "t-test",
                         "wilcox" = "Wilcoxon rank-sum"),
    `Mean (Rural)`  = mean_group1,
    `SD (Rural)`    = sd_group1,
    `n (Rural)`     = n_group1,
    `Mean (Urban)`  = mean_group2,
    `SD (Urban)`    = sd_group2,
    `n (Urban)`     = n_group2,
    `Test statistic` = statistic,
    `Effect size`    = effect_size,
    `p-value`        = p_value,
    Sig = significance
  )

FinalQualitativeTable %>% 
  mutate(
    across(c(`Mean (Rural)`:`Effect size`), ~round(., 2)),
    `p-value` = signif(`p-value`, 3)
  ) %>% 
  kable(
    caption = "Comparison of rapid bioassessment scores between rural and urban streams.",
    align = "lcccccccccc",
    booktabs = TRUE
  ) %>% 
  kable_styling(full_width = FALSE, position = "center", font_size = 10)

Comparison of rapid bioassessment scores between rural and urban streams.
Parameter	Test	Mean (Rural)	SD (Rural)	n (Rural)	Mean (Urban)	SD (Urban)	n (Urban)	Test statistic	Effect size	p-value
NA	Wilcoxon rank-sum	1.33	2.16	6	0.25	0.61	6	22.0	0.94	0.462
NA	Wilcoxon rank-sum	2.58	1.11	6	2.17	1.13	6	23.5	0.86	0.408
NA	Wilcoxon rank-sum	2.50	1.64	6	1.50	0.55	6	24.0	0.83	0.342

5 eDNA Community Metrics

5.1 Summary of Overall Dataset

eDNA_richness_summary %>%
  mutate(Grouping = if_else(Stream == "All streams", "", Grouping)) %>%
  arrange(desc(Stream == "All streams"), Grouping, Stream) %>%
  kable(
    caption = "eDNA taxon richness (global and by stream and grouping).",
    booktabs = TRUE,
    align   = "llcccccc",
    digits  = 0
  ) %>%
  kable_styling(
    full_width  = FALSE,
    position    = "center",
    font_size   = 10,
    stripe_color = "#f7f7f7"
  ) %>%
  # draw a horizontal rule after the first data row
  row_spec(
    1,
    extra_css        = "border-bottom: 1px solid black;",   # for HTML
    extra_latex_after = "\\midrule"                         # for PDF/LaTeX
  )

eDNA taxon richness (global and by stream and grouping).
Stream	Grouping	Species	Genera	Families	Orders	Classes	Phyla
All streams		119	73	60	50	20	7
Carey	Rural	27	21	22	18	10	6
Carey_FR	Rural	29	21	16	10	4	2
Cherry	Rural	33	30	28	21	11	7
Griffin	Rural	34	22	22	18	9	4
Judd	Rural	19	19	20	18	8	4
Ravensdale	Rural	21	27	27	23	12	6
Rock	Rural	6	7	7	6	2	1
Juanita	Urban	58	36	30	21	13	6
Longfellow	Urban	26	23	19	14	8	5
Longfellow_FR	Urban	47	35	30	23	15	7
McAleer	Urban	25	26	28	24	12	5
Miller	Urban	13	15	18	15	7	4
Thornton	Urban	32	22	24	20	11	5
Venema	Urban	12	16	20	17	9	5

5.2 Site-Level Richness & Diversity metrics

site_metrics_pretty %>%
  kable(
    caption = "eDNA-based richness, diversity, and dominance metrics by stream.",
    booktabs = TRUE,
    digits = c(NA, NA, 0, 0, 0, 0, 2, 1, 1),
    align = "lcccccccccc"
  ) %>%
  kable_styling(
    full_width = FALSE,
    position   = "center",
    font_size  = 10
  ) %>%
  footnote(
    general = "Shannon diversity index, calculated from the proportional read abundance of collapsed eDNA taxa (taxon_final). Higher values indicate greater diversity and more even distribution of taxa.",
    general_title = ""
  )

eDNA-based richness, diversity, and dominance metrics by stream.
Stream	Grouping	Total genera	Total species	Diatom genera	Diatom species	Shannon H*	% reads in top 5	% reads in top 10
Carey	Rural	18	25	16	24	2.96	72.3	92.2
Cherry	Rural	28	32	27	32	2.89	82.5	96.4
Griffin	Rural	20	34	20	34	3.09	76.9	94.6
Judd	Rural	17	18	14	17	2.68	80.0	94.4
Ravensdale	Rural	23	18	19	17	2.94	76.6	88.8
Rock	Rural	8	7	8	7	1.41	97.6	100.0
Juanita	Urban	34	57	30	55	3.29	73.3	90.2
Longfellow	Urban	21	25	20	24	2.17	71.0	90.2
McAleer	Urban	20	22	18	21	2.97	73.8	90.2
Miller	Urban	12	11	10	10	2.70	86.9	99.2
Thornton	Urban	19	30	16	29	2.79	80.2	96.1
Venema	Urban	12	10	11	10	2.78	76.8	98.8
Shannon diversity index, calculated from the proportional read abundance of collapsed eDNA taxa (taxon_final). Higher values indicate greater diversity and more even distribution of taxa.

Shannon diversity (H′) was calculated from collapsed eDNA taxa (taxon_final) derived by aggregating ASVs to the lowest confidently assigned taxonomic level. For each stream, relative read abundance (pᵢ) was computed as the proportion of total reads assigned to each taxon. Shannon diversity was then calculated as H′ = −Σ pᵢ ln pᵢ, which incorporates both richness and evenness components of community structure.

5.3 Site-Level Pollution Tolerance Metrics

Pollution tolerance metrics were calculated using genus- and species-level pollution tolerance index (PTI) scores from Lange-Bertalot, 1979; Bahls, 1993; and Barbour et al., 1999.

pti_results_pretty %>%
  kable(
    caption = "Periphyton tolerance metrics (PTI components) by stream.",
    booktabs = TRUE,
    digits = c(
      NA,  # Stream
      NA,  # Grouping
      2,   # WAT
      2,   # Mean PTV
      1,   # % sensitive
      1    # % tolerant
    ),
    align = "llcccc"
  ) %>%
  add_header_above(
    c(" " = 2,
      "Tolerance metrics" = 4)
  ) %>%
  kable_styling(
    full_width = FALSE,
    position   = "center",
    font_size  = 10
  ) %>%
  footnote(
    general = paste(
      "Weighted average tolerance (WAT) is calculated using abundance-weighted periphyton tolerance values.",
      "Mean periphyton tolerance value is the unweighted mean of taxon tolerance scores.",
      "% sensitive taxa and % tolerant taxa represent the proportion of taxa in the lowest and highest",
      "tolerance categories, respectively."
    ),
    general_title = ""
  )

Periphyton tolerance metrics (PTI components) by stream.
		Tolerance metrics
Stream	Grouping	Weighted average tolerance (WAT)	Mean periphyton tolerance value	% sensitive taxa	% tolerant taxa
Carey	Rural	2.75	2.70	19.9	14.2
Carey_FR	Rural	2.54	2.81	38.7	11.4
Cherry	Rural	2.58	2.75	28.5	14.8
Griffin	Rural	2.96	2.65	25.2	29.7
Judd	Rural	4.10	3.05	1.4	49.8
Ravensdale	Rural	2.82	2.74	11.3	14.7
Rock	Rural	3.01	2.89	8.4	12.3
Juanita	Urban	2.96	2.95	16.2	23.5
Longfellow	Urban	4.10	3.44	4.4	61.7
Longfellow_FR	Urban	3.44	3.36	9.5	34.6
McAleer	Urban	3.54	3.63	12.1	43.4
Miller	Urban	3.56	3.62	5.6	25.5
Thornton	Urban	3.73	3.27	2.4	39.5
Venema	Urban	3.41	3.00	7.1	23.8
Weighted average tolerance (WAT) is calculated using abundance-weighted periphyton tolerance values. Mean periphyton tolerance value is the unweighted mean of taxon tolerance scores. % sensitive taxa and % tolerant taxa represent the proportion of taxa in the lowest and highest tolerance categories, respectively.

6 eDNA-inferred Diatom Metrics

6.1 Global Diatom/Tolerance Metrics

DiatomMetrics_Graph

FinalDiatomMetricsTable <- pti_tests_final %>% 
  transmute(
    `Metric` = metric_name,
    Test = dplyr::recode(test_type,
                         "t" = "t-test",
                         "wilcox" = "Wilcoxon rank-sum"),
    `Mean (Rural)`  = mean_group1,
    `SD (Rural)`    = sd_group1,
    `n (Rural)`     = n_group1,
    `Mean (Urban)`  = mean_group2,
    `SD (Urban)`    = sd_group2,
    `n (Urban)`     = n_group2,
    `Test statistic` = statistic,
    `Effect size`    = effect_size,
    `p-value`        = p_value,
    Sig = significance
  )%>%
  mutate(
    # Match the same metric naming used in your plot
    Metric = dplyr::recode(
      Metric,
      mean_PTV         = "Mean PTV",
      pct_sensitive    = "% sensitive taxa",
      pct_tolerant     = "% tolerant taxa",
      richness_genus   = "Genus richness",
      richness_species = "Species richness",
      WAT              = "Weighted average tolerance (WAT)"
    ),
    # Apply correct facet/table ordering
    Metric = factor(
      Metric,
      levels = c(
        "Genus richness",
        "Species richness",
        "% sensitive taxa",
        "% tolerant taxa",
        "Mean PTV",
        "Weighted average tolerance (WAT)"
      )
    )
  ) %>%
  arrange(Metric)

FinalDiatomMetricsTable %>% 
  mutate(
    across(c(`Mean (Rural)`:`Effect size`), ~round(., 2)),
    `p-value` = signif(`p-value`, 3)
  ) %>% 
  kable(
    caption = "Comparison of diatom-inferred and pollution tolerance metrics between rural and urban streams.",
    align = "lcccccccccc",
    booktabs = TRUE
  ) %>% 
  kable_styling(full_width = FALSE, position = "center", font_size = 10)

Comparison of diatom-inferred and pollution tolerance metrics between rural and urban streams.
Metric	Test	Mean (Rural)	SD (Rural)	n (Rural)	Mean (Urban)	SD (Urban)	n (Urban)	Test statistic	Effect size	p-value
Genus richness	t-test	22.00	7.28	7	25.71	8.32	7	-0.89	0.48	0.39200
Species richness	t-test	25.14	9.77	7	31.43	16.98	7	-0.85	0.45	0.41700
% sensitive taxa	t-test	19.05	12.89	7	8.21	4.77	7	2.09	-1.12	0.07210
% tolerant taxa	Wilcoxon rank-sum	21.00	14.11	7	36.01	13.81	7	9.00	1.78	0.05530
Mean PTV	t-test	2.80	0.13	7	3.32	0.27	7	-4.60	2.46	0.00138
Weighted average tolerance (WAT)	Wilcoxon rank-sum	2.96	0.53	7	3.53	0.35	7	9.00	1.78	0.05530

7 Community Composition and Multivariate Analyses (eDNA-Inferred Diatoms)

This section evaluates whether diatom communities differ between rural and urban stream types and identifies the biological and environmental factors that best explain those differences. Analyses include tests of community dissimilarity (PERMANOVA), assessments of community dispersion, constrained ordinations (RDA and db-RDA), and variation partitioning.

7.1 PERMANOVA (Bray–Curtis): Do Urban and Rural Streams Support Different Diatom Communities?

PERMANOVA asks whether overall community composition differs between two or more groups—in this case, rural vs. urban streams. It does this by comparing multivariate distances (Bray–Curtis) between and within groups. A significant result means the groups have different community “centroids” (i.e., different typical communities).

permanova_diatom

## Permutation test for adonis under reduced model
## Permutation: free
## Number of permutations: 999
## 
## adonis2(formula = diatom_dist ~ Grouping, data = diatom_meta, permutations = 999)
##          Df SumOfSqs      R2      F Pr(>F)  
## Model     1   0.5649 0.15488 1.8326  0.039 *
## Residual 10   3.0823 0.84512                
## Total    11   3.6472 1.00000                
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Diatom community composition differed significantly between rural and urban streams (Bray–Curtis PERMANOVA, p = 0.01). Approximately 15% of community variation was explained by the urban–rural grouping, which is a relatively strong effect size for ecological community data. This indicates that land-use differences are associated with distinct diatom assemblages.

7.2 Multivariate Dispersion: Are Differences Due to Spread or True Shifts in Composition?

The betadisper test checks whether groups differ in their within-group variability (dispersion). If dispersion differs, PERMANOVA significance can be driven by unequal spread rather than genuine compositional shifts.

anova(betadisper_res)

## Analysis of Variance Table
## 
## Response: Distances
##           Df   Sum Sq  Mean Sq F value Pr(>F)
## Groups     1 0.022072 0.022072  1.8643 0.2021
## Residuals 10 0.118390 0.011839

We did not detect differences in multivariate dispersion between rural and urban streams (F = 2.60, p = 0.13). Because dispersion is similar, the PERMANOVA result reflects a shift in community centroid, not an artifact of different within-group variability.

The NMDS biplot below visualizes these centroid differences and highlights taxa or gradients associated with each group.

gg_NMDS_vecs

Urban streams cluster toward vectors representing tolerant diatom taxa, whereas rural streams cluster toward vectors representing sensitive taxa and higher diversity.

7.3 Unconstrained Ordination of Diatom Communities (PCA-like RDA)

We used an unconstrained redundancy analysis (RDA) on Hellinger-transformed eDNA diatom data to summarize the major axes of community variation across sites. Because no environmental variables were included, this ordination is equivalent to a principal components analysis (PCA) of the transformed community matrix. The input matrix (diatom_hell) contained proportional abundances of eDNA-inferred diatom taxa (columns) for each stream (rows).

rda_site_plot

7.3.1 Site-Level Patterns

PC1 and PC2 explained 23% and 17% of total community variation, respectively (40% cumulative), based on eigenvalues from the unconstrained RDA..

The site-score ordination shows clear structure in diatom assemblages across the development gradient:

Urban sites tend to cluster more tightly, indicating greater similarity in their diatom communities.
Rural sites are more widely dispersed, reflecting greater ecological heterogeneity among low-development watersheds.
Axis 1 (23% of variance) captures the primary gradient separating these assemblages, whereas Axis 2 (17%) represents secondary differences that are not strongly aligned with the urban–rural contrast.

These interpretations are based on the distribution of site scores in the ordination diagram; the eigenvalues quantify how much variance each axis explains but do not by themselves indicate clustering or group separation.

7.3.2 Taxon-Level Patterns (PCA-like Biplot)

A PCA-like biplot was used to overlay diatom taxa vectors on the same ordination space, illustrating which taxa contribute most strongly to variation along PC1 and PC2. Direction of arrow shows which sites have relatively higher abundance of that taxon. Length of arrow shows how strongly that taxon contributes to variation

pca_biplot

In the biplot, arrows represent diatom taxa and points represent streams:

Rural sites are more strongly associated with taxa such as Reimeria sinuata, Melosira varians, several Gomphonema species, and certain Achnanthidium and Nitzschia ASVs. Although some Achnanthidium and Nitzschia species are often considered tolerant in classical indices, the specific ASVs detected here may represent generalist or moderately sensitive lineages, highlighting the finer taxonomic resolution of eDNA.
Urban sites cluster toward taxa such as Navicula spp., Cocconeis spp., and Amphora pediculus, which are commonly associated with nutrient enrichment and altered hydrology.
PC1 therefore represents a gradient from assemblages dominated by Navicula / Cocconeis (Urban) to assemblages characterized by Gomphonema, Reimeria, and Melosira (Rural), while PC2 captures additional within-group variation likely related to microhabitat or hydrologic differences among sites.

Together, the unconstrained ordination and biplot show that diatom communities naturally separate along a gradient consistent with watershed development intensity, even before environmental variables or tolerance metrics are explicitly included in the models.

7.4 dbRDA using Tolerance Metrics

Here, we ask whether diatom traits themselves (e.g., tolerance values, % tolerant taxa, richness) can explain observed community differences. db-RDA is a constrained ordination that uses dissimilarity measures (Bray–Curtis) rather than Euclidean distance.

dbRDA_ptv_plot

Tolerance-based metrics explained 47% of variation in diatom community composition (R² = 0.47, p = 0.022).

Urban sites grouped strongly on the positive side of dbRDA1 and were associated with higher Weighted Average Tolerance (WAT), higher mean PTV, and greater % tolerant taxa.

Rural sites grouped on the negative side of dbRDA1 and were associated with higher diatom richness and greater % sensitive taxa.

dbRDA2 represented additional within-group variability not directly linked to tolerance structure. Overall, tolerance traits appear to be major drivers of community distinction along the urbanization gradient.

7.5 db-RDA with Environmental Predictors (forward selection)

Forward-selection db-RDA evaluates which environmental variables (nutrients, canopy cover, stream width, etc.) significantly explain community differences after controlling for model complexity. Variables are added only if they improve adjusted R².

gg_dbRDA   # the db-RDA (capscale) biplot with TN_z (and any other selected vars)

Forward selection identified total nitrogen (TN) as the only predictor that significantly improved model fit. Although multiple variables were considered, TN alone explained the largest share of constrained variation.

Urban streams with higher TN clustered on the right side of dbRDA1.

Rural streams with lower TN clustered on the left.

dbRDA2 represented unexplained residual structure unrelated to TN.

These results suggest that nutrient enrichment—specifically TN—is the dominant environmental driver of diatom community differences among the variables tested.

7.6 Variation Partitioning (Nutrients vs Canopy)

Variation partitioning quantifies how much of the explained community variation is uniquely attributable to different sets of predictors (e.g., nutrients vs. canopy) and how much is shared.

varpart_res

## 
## Partition of variance in RDA 
## 
## Call: varpart(Y = diatom_mat2, X = X_nutr, X_can)
## 
## Explanatory tables:
## X1:  X_nutr
## X2:  X_can 
## 
## No. of explanatory tables: 2 
## Total variation (SS): 743754180 
##             Variance: 67614016 
## No. of observations: 12 
## 
## Partition table:
##                      Df R.squared Adj.R.squared Testable
## [a+c] = X1            4   0.36628       0.00415     TRUE
## [b+c] = X2            1   0.15702       0.07272     TRUE
## [a+b+c] = X1+X2       5   0.45458       0.00006     TRUE
## Individual fractions                                    
## [a] = X1|X2           4                -0.07266     TRUE
## [b] = X2|X1           1                -0.00409     TRUE
## [c]                   0                 0.07681    FALSE
## [d] = Residuals                         0.99994    FALSE
## ---
## Use function 'rda' to test significance of fractions of interest

plot(varpart_res, bg = c("lightblue","lightgreen"))

The combined nutrient + canopy model had a raw explanatory power of 45.5% (R² = 0.4546), but after adjusting for small sample size and model complexity, the adjusted R² dropped to essentially zero (0.00006). As a result, the residual fraction was 99.99%, indicating that almost all variation in the Hellinger-transformed community matrix remained unexplained.

A small shared component (~7.7%) was attributable jointly to nutrient and canopy gradients, but neither variable explained unique variation once the other was accounted for (adjusted unique fractions were negative and therefore interpreted as zero).

These results are consistent with (1) the small sample size (n = 12), which strongly penalizes adjusted R² in constrained ordination, and (2) high site-to-site variability commonly observed in stream periphyton communities. Despite the large residuals, nutrient and canopy gradients still appeared ecologically meaningful in earlier RDA/db-RDA analyses, suggesting underlying structure that cannot be cleanly partitioned in this small dataset.

8 Taxonomy via Microscopy (coming soon…)

9 Discussion

Urban streams: high nutrients, high canopy, low biomass, distinct diatom assemblages with higher tolerance values. Rural streams: lower nutrients, somewhat more open canopy, higher biomass, richer and less tolerant diatom communities. Evidence that light limitation and hydromorphology override nutrient enrichment in controlling biomass. Management implications (e.g., riparian thinning vs nutrient control, depending on objectives).

Periphyton Pilot Study: Biomass, eDNA, and Environmental Drivers

Jeremy Walls

December 2025