Skip to main content

Effect of river flow, temperature, and water chemistry on proliferations of the benthic anatoxin-producing cyanobacterium Phormidium


Proliferations of the benthic anatoxin-producing cyanobacterium Phormidium are increasing in prevalence in cobble-bed rivers worldwide. This proliferation is of particular concern when rivers are used as sources of drinking water or for recreation. Little is known about the physicochemical variables promoting proliferations, and our existing knowledge is based on data from only a few rivers. We assessed Phormidium cover, physicochemical variables, and anatoxin concentrations at 10 sites in 7 New Zealand rivers every week for 2 y. Generalized additive mixed models (GAMMs) identified dissolved inorganic N (DIN) over the accrual period <0.8 mg/L, dissolved reactive P accrual <0.005 mg/L, water temperatures >15°C, and conductivity as having positive and statistically significant effects on % Phormidium cover. Flow intensity, expressed relative to the long-term median, had a positive effect up to 0.4× the median flow and a negative effect when >0.5× the median flow. Quantile regression models showed marked variability among sites in relation to the flow intensity required to reduce % Phormidium cover (90th percentile ranged 0.65–249× the long-term median flow). Anatoxins were detected in variable concentrations in samples from 7 of the 10 sites. GAMMs identified strong relationships between elevated toxin concentrations and low conductivity and increasing % Phormidium cover, and significantly lower toxin concentrations when DIN was <0.2 mg/L. These data demonstrate that multiple physicochemical variables influence Phormidium proliferations and toxin concentrations and indicate that the relative importance of these variables differs among rivers and sites.

Planktonic cyanobacterial blooms have become strongly associated with eutrophic lentic habitats where they affect aquatic life and are a health risk to animals and humans (Codd et al. 2005, Kouzminov et al. 2007, Havens 2008). However, benthic cyanobacteria in lentic and lotic environments pose a similar suite of issues (Quiblier et al. 2013). An increase in proliferations of mat-forming, toxin-producing Phormidium spp. in cobble-bed rivers worldwide is of particular concern (Gugger et al. 2005, Wood et al. 2007, Fetscher et al. 2015, McAllister et al. 2016, Echenique-Subiabre et al. 2016). In contrast to planktonic cyanobacterial blooms, Phormidium proliferations generally occur at sites with relatively good water quality (McAllister et al. 2016).

In New Zealand, periphyton communities in cobble-bed rivers are usually dominated by diatoms (Bacillariophyta) or green algae (Chlorophyta), with cyanobacteria reported in only low abundance in a nationwide survey in the late 1980s (Biggs 1990). However, proliferations of Phormidium-dominated (hereafter Phormidium) mats in human-influenced cobble-bed rivers have increased in the past decade (see McAllister et al. 2016 for a review of current distribution). These mats are dominated by Phormidium, but they contain other organisms including bacteria, other cyanobacteria, and algae, which are jointly adhered to the substrate by extracellular polymeric substances (Hart et al. 2013, Brasell et al. 2015). Investigators have demonstrated in polyphasic studies of Phormidium mats in New Zealand that the dominant species is Phormidium autumnale (Heath et al. 2010, Wood et al. 2012, Harland et al. 2014), but this taxonomic assignment is the subject of ongoing revision (Strunecký et al. 2013). Phormidium mats are very distinctive (several millimeters thick, cohesive, and black–brown), thereby enabling macroscopic identification and assessment of their abundance.

Phormidium proliferations can cause shifts in macroinvertebrate abundance and composition (Wood et al. 2014) and cause fish caught from affected rivers to be tainted with an ‘earthy/musty’ taste, and the proliferations are aesthetically undesirable. Phormidium proliferations have gained most notoriety because of the ability of some strains to produce powerful neuromuscular-blocking toxins, collectively known as anatoxins (Wood et al. 2007). These toxins pose a threat to humans and animals when consumed or when they come in contact with contaminated water (Wood et al. 2011). Anatoxins have killed ~100 dogs in New Zealand in the last 5 to 10 y and have led to issuance of health warnings advising against contact with rivers experiencing proliferations. Four chemical structural variants are found in Phormidium mats in New Zealand: anatoxin-a (ATX), homoanatoxin-a (HTX), and the reduced derivatives dihydro-anatoxin-a (dhATX) and dihydro-homoanatoxin-a (dhHTX). In Phormidium proliferations, ATX, HTX, dhATX, and dhHTX are almost always detected simultaneously, but their relative concentrations vary (McAllister et al. 2016). Anatoxin concentrations can change rapidly (Heath et al. 2011), a finding that has been attributed at least partly to the presence and relative abundance of toxic and nontoxic genotypes in the mats (Wood et al. 2012).

The abundance and composition of periphyton in rivers is constrained by flow, water chemistry (e.g., nutrients), light, temperature, and substrate and, hence, is affected by surrounding land use (Biggs and Smith 2002, Pan et al. 2004, O’Brien and Wehr 2010, Klose et al. 2012). The importance of the hydrological regime and nutrients in regulating periphyton development has been highlighted in New Zealand cobble-bed rivers (e.g., Biggs and Close 1989). Biggs (2000) concluded that when the effects of accrual time were taken into account, the best predictors of periphyton biomass were dissolved inorganic N (DIN) and dissolved reactive P (DRP) concentrations, with positive linear relationships observed. However, Phormidium mats do not follow this pattern and tend to proliferate at low to moderate nutrient concentrations (DIN > 0.1–0.2 mg/L and DRP < 0.01 mg/L; Heath et al. 2015, McAllister et al. 2016). Internal processes also may play a role in maintaining growth once mats are established (Wood et al. 2015). Phormidium mats can tolerate higher flows than other periphyton types (Hart et al. 2013) and are associated with stable substrates (Heath et al. 2015), but little is known about the magnitude and frequency of flushing events required to reduce Phormidium mats. Improved knowledge of the physicochemical variables regulating Phormidium accrual and dispersal cycles will enable new insights into Phormidium ecology, provide explanations regarding their recent increase in abundance and prevalence, and ultimately may lead to mitigation strategies.

We assessed % Phormidium cover and collected samples for water chemistry and toxin analysis weekly from 10 sites (7 rivers) experiencing a gradient of Phormidium proliferations over a 12- to 18-mo period. Our objectives were to: 1) assess the relative importance of physicochemical variables and river flow in regulating Phormidium proliferations, 2) investigate the intensity of flushing flow required to remove or reduce Phormidium mats and establish whether the same flow intensity would reduce mats in all rivers, and 3) test for relationships between anatoxin concentrations and environmental variables.


Site characteristics

This study was undertaken in the Manawatu River catchment, and all sites were on the mainstem or tributaries of the Manawatu River (Fig. 1). One or 2 sites at each of 7 rivers (Makakahi, Manawatu, Mangatainoka, Oroua, Oruakeretaki, Tiraumea, and Tokomaru) were selected for Phormidium monitoring (Table 1, Fig. 1). Sites were chosen to represent rivers with and without historical Phormidium proliferations, with similar geology, and a range of nutrient and flow regimes (Table 1). Land cover was primarily heavy pastoral or native vegetation, with urban accounting for <1% in catchments (Table 1). Riparian or overhead shading from vegetation was negligible at all sites. Sites 2/3 and 6/7 were paired sites within 1.9 and 0.62 km of each other but situated down- and upstream (respectively) of sewage treatment plant (STP) discharges (Fig. 1). Surveying and sampling were undertaken approximately weekly at 7 sites (sites 1, 4–7, 9, and 10) from 6 January 2012 to 26 or 28 June 2013, at 2 sites (sites 2 and 3) from 18 July 2012 to 28 June 2013, and at site 8 from 17 August 2012 to 28 June 2013.

Figure 1.
Figure 1.

Locations of sampling sites and flow gages in New Zealand.

Table 1.

Physical and hydrological conditions at study sites. Long-term median flow is from the gaging site closest to each sampling site (see Fig. 1). u/s = upstream, d/s = downstream, STP = sewage treatment plant. Catchment data were retrieved from the Freshwater Ecosystems of New Zealand database (Leathwick et al. 2011).

RiverSiteCatchment area (km2)Stream orderHeavy pastoral land use (%)Native vegetation land use (%)Other land use (%)Median flow (m3/s)
Makakahi River11635791833.18
Manawatu River d/s STP2402277321673.4
Manawatu River u/s STP3391577420673.4
Mangatainoka River411409822.13
Mangatainoka River54136762048.90
Oroua River d/s STP65856751697.10
Oroua River u/s STP75826761687.10
Oruakeretaki River8544673031.42
Tiraumea River97446831257.21
Tokomaru River10564177221.25

Site surveys and sampling

Surveys were conducted in runs and riffles. The length of the reach surveyed varied from 30 to 100 m among sites. At each site, 4 transects perpendicular to the water’s edge and extending to a maximum depth of 0.6 m were surveyed. Percent Phormidium cover was estimated visually at 5 points along each transect with the aid of an underwater viewer (giving a total of 20 views/site; Ministry for the Environment and Ministry of Health 2009). Substrate composition (macrophytes, moss, mud [<0.2 mm], sand [0.2–2 mm], fine gravel [2–8 mm], gravel [8–64 mm], cobble [64–264 mm], boulder [>264 mm], and bedrock) was assessed visually at the same 20 viewpoints.

Water temperature, conductivity, and pH were measured once at each site on each sampling occasion using a handheld multiparameter PCTESTR 35 sonde (Oakton Instruments, Vernon Hills, Illinois), and water samples were collected for dissolved nutrient analysis. When Phormidium mats were present, a single sample was collected from 10 rocks at each site by scraping mat material from a defined circular area (6.25-cm diameter). The 10 samples were pooled, chilled, and stored in darkness for transport to the laboratory where they were stored frozen (−20°C; within 5 h) until toxin analysis.

Continuous river flow data were obtained from permanent flow-gaging stations in close proximity to each sampling site (Fig. 1). Daily mean flow was calculated for each site. Long-term median flows at each site were obtained from Horizons Regional Council (Table 1).

Laboratory analysis

Nutrients in filtered water were analyzed on a Lachat QuickChem® Flow Injection Analyser (FIA+ 8000 Series, Zellweger Analytics, Loveland, Colorado) based on the methods provided by APHA (2005). The nutrients analyzed were: NH4-N, NO2-N, NO3-N, and dissolved reactive P (DRP). The accredited detection limits were 0.01 mg/L NH4-N, 0.002 mg/L NO2-N and NO3-N, and 0.005 mg/L DRP.

Phormidium mat samples were lyophilized and homogenized. Lyophilized material (100 mg) was suspended in 10 mL of Milli-Q water (MQ; Millipore, Billerica, Massachusetts) containing 0.1% formic acid and sonicated (30 min; Cole Parmer 8890; Biolab, Auckland, New Zealand) on ice. Samples were centrifuged (3000g, 10 min). The supernatants were transferred to a septum-capped vial and analyzed for ATX, HTX, dhATX and dhHTX based on liquid chromatography-mass spectrometry (LC-MS) as described by Wood et al. (2016). External standards were used for ATX calibration, using dilutions of a certified reference material for ATX (National Research Council, Canada) to prepare working standards that ranged from 0.5 to 20 ng/mL in 0.1% formic acid. A relative response factor of 1, with ATX as the calibration reference, was used to quantify HTX, dhATX, and dhHTX.

Data analyses

Generalized additive mixed models (GAMMs; Hastie and Tibshirani 1990) were used to model nonlinear trends in % Phormidium cover and total toxin concentration (sum of ATX, HTX, dhATX, and dhHTX in mg/kg of lyophilized material) in relation to time (wk) of the year and environmental variables. The % Phormidium cover data were means from 20 views obtained at each approximately weekly sampling occasion (n = 677). Total toxin models were fitted only with observations from sampling occasions when Phormidium samples were collected for toxin analysis (n = 176). Both % Phormidium cover and total toxin data were strictly positive, highly right-skewed, over-dispersed, and with a large proportion of 0 values, so models were constructed with log-normal errors. Seven noncollinear predictor variables were selected on the basis of initial data exploration and ecological knowledge and were included in the models. Collinearity among predictor variables was checked using the variance inflation factor (VIF < 3; Zuur et al. 2010). Week of the year was included as a continuous covariate (1–52) to account for seasonal trends based on a cyclic cubic spline to allow the continuity between weeks 1 and 52. Biggs and Close (1989) suggested that simple point-by-point correlations between nutrient concentrations and periphyton biomass do not provide a true indication of the water-column nutrients to which periphyton has been exposed over its accrual period. To address this problem, cumulative mean DIN (the sum of NH4-N, NO2-N, and NO3-N) and DRP concentrations over the accrual period were calculated for each time point (i.e., the mean concentrations since the previous sampling point at which the % Phormidium cover was 0). ‘Accrual’ DIN and DRP were √(x)- and log(x)-transformed, respectively. The influence of flow was considered by including log(x)-transformed ‘times median flow’, i.e., the maximum flow magnitude (obtained from daily mean flow data) between sampling periods divided by the long-term median flow for a given site. Water temperature, accrual DIN, accrual DRP, times median flow, and conductivity were included as continuous covariates, and site as a categorical variable with 10 levels. Percent Phormidium cover also was included as a continuous covariate in the total toxin model. To account for temporal autocorrelation in the weekly time-series data, models were fitted based on autoregressive moving average correlation structure (ARMA) of order 1. Models were selected with a stepwise procedure based on the generalized Akaike Information Criteria (GAIC) and were validated by inspecting the deviance residuals. Final models are presented as partial effects plots, which show the effect of each predictor variable conditional to others in the model. The partial effects of each predictor are displayed as cubic splines showing either negative or positive effects relative to the overall mean of the response variable centered on 0. Partial plots also show standard errors around the fitted spline and partial residuals for each observation. The nonlinear interaction between accrual DIN and DRP on % Phormidium cover was further investigated with a tensor product smoother in a GAMM model. The interactive effect of DIN and DRP on % Phormidium cover was visualized in a contour plot of the predicted values.

We also investigated the magnitude of river flow required to reduce % Phormidium cover at each site to <20%. Analysis was undertaken only for the 4 sites that had % Phormidium cover >20% on multiple occasions (sites 1, 3, 5, and 10). For each site, the maximum flow (obtained from daily means) divided by the long-term median flow on the 7 d prior to surveying was selected. A time frame of 7 d was selected because this is the length of time over which % Phormidium cover has been observed to increase from 0% to >20% cover (Heath et al. 2011, Wood and Young 2012). Greater than 20% Phormidium cover was selected because this level is the lower threshold of the amber alert level in the New Zealand Guidelines for Cyanobacteria in Recreational Freshwater (Ministry for the Environment and Ministry of Health 2009). If a site has a cover of 20%, cover probably will continue to increase in the absence of a flushing event. For each site, quantile regression models (80th, 85th, and 90th percentiles) were fitted to the % Phormidium cover data with log(x)-transformed river flow as the predictor. The river flow at which % Phormidium cover stayed at 0% was estimated from the regression lines.

All statistical analyses were performed with the software R (version 3.1.1; R Project for Statistical Computing, Vienna, Austria), GAMMs models were run with the gamlss package (Rigby and Stasinopoulos 2005), and quantile regressions were run with the quantreg package (Koenker 2012).


Substrate composition, spatial and temporal patterns in % Phormidium cover, physicochemical variables

The dominant substrate at all sampling sites was large (12–25 cm) and small (6–12 cm) cobbles (>42% at all sites; Fig. S1). Phormidium mats were present at all 10 river sites, but mean cover was spatially and temporally variable (Fig. 2). Five sites (1, 2, 5, 6, and 10) had mean Phormidium cover >40% on multiple occasions. Three sites (7, 8, and 9) had high % cover in single views (40–60%), but the mean cover was ≤15% (Fig. 2). At 2 sites (3 and 4) mean cover was <2% (Fig. 2).

Figure 2.
Figure 2.

Mean % Phormidium cover (± 95% confidence intervals) assessed weekly over the study period (6 January 2012 to 26 or 28 June 2013) at each site (see Table 1 for site names).

At the sites with high % Phormidium cover, the period with the greatest % cover was generally from early austral summer (November) to the end of autumn (May). Periods of high cover varied among sites over the study period. For example, high cover was recorded at site 5 and 10 from mid-April to mid-June 2013, but this pattern was not observed at other sites (Fig. 2). Although situated <2 km apart, % Phormidium cover differed notably between sites 2 and 3 and sites 6 and 7. In general, sites downstream of STPs had higher cover. This situation was particularly apparent at site 2 in March 2013 and site 6 in January 2013 (Fig. 2).

Median DIN concentrations varied markedly among sites. The lowest concentration was recorded at site 4 (0.12 mg/L) and the highest at site 6 (0.87 mg/L; Fig. 3A). Median DRP concentrations were generally low at all sites (0.01 mg/L). The 3 highest concentrations were recorded at sites 2, 6, and 7 (Fig. 3B). Median water temperatures ranged between 9.5 and 15.5°C and were comparable across sites except site 4, which was cooler (Fig. 3C). The lowest (4.6°C) and highest (26.2°C) water temperatures were recorded at site 4 in July 2012 and site 5 in February 2013, respectively. Median conductivity was variable among sites (range: 53–313 µS/cm at sites 4 and 9, respectively; Fig. 3D).

Figure 3.
Figure 3.

Box plots of √(dissolved inorganic N) (DIN) (A), log(dissolved reactive P) (DRP) (B), temperature (C), and conductivity (D) assessed at each site. Lines in boxes are medians, box ends are quartiles, and whiskers extend to the lowest or highest data point ≤1.5× the interquartile range. Black dots are outliers.

Percent Phormidium cover in relation to environmental variables

The most parsimonious GAMM explained 56.2% of the total deviance in % Phormidium cover and included all predictor variables (Table 2; p < 0.05). Week of the year affected predicted % Phormidium cover positively for the first half of the year (January–May; p < 0.001; Fig. 4A). The model predicted greater cover with increasing accrual DIN concentrations up to 0.8 mg/L, after which the effect remained constant (p < 0.001; Fig. 4B). Accrual DRP had a weak, marginally significant negative effect on % Phormidium cover when DRP > 0.005 mg/L (p = 0.03; Fig. 4C). Water temperatures >15°C had a strongly positive effect, whereas colder temperatures did not affect % Phormidium cover (p < 0.001; Table 2, Fig. 4D). Times median flow had a strongly positive effect for flows <0.4 and a negative effect for flows >0.5 (p < 0.001; Fig. 4E). Percent Phormidium cover increased monotonically with conductivity (p < 0.001; Fig. 4F). Percent Phormidium cover differed strongly among sites (p < 0.05; Fig. 4G).

Table 2.

Results of the log-normal generalized additive mixed model for % Phormidium cover data (n = 677). The model was fitted using an autoregressive moving-average correlation structure of order 1. Fit degrees of freedom = 35. DIN = dissolved inorganic N, DRP = dissolved reactive P, SE = standard error. * = p < 0.05, ** = p < 0.01, *** = p < 0.001.

Week of the year−0.040.01−7.06***
DIN accrual2.220.514.36***
DRP accrual−0.360.17−2.11*
Water temperature0.120.034.54***
Times median flow−0.500.09−5.41***
Site 2−1.390.49−2.83**
Site 3−3.730.47−7.90***
Site 4−1.120.46−2.45*
Site 5−0.500.36−1.36 
Site 6−3.250.46−7.10***
Site 7−3.080.42−7.34***
Site 8−1.250.44−2.87**
Site 9−4.400.69−6.42***
Site 102.540.406.29***
Figure 4.
Figure 4.

Partial plots of seasonal (week of the year) (A), accrual dissolved inorganic N (DIN) (B), accrual dissolved reactive P (DRP) (C), water temperature (D), times median flow (E), conductivity (F), and site (G) effects on % Phormidium cover. Solid lines represent cubic splines (± SE, dashed lines) fitted based on log-normal generalized additive mixed model. See methods for description of x-axis partial effect scale.

Sites could be grouped loosely into 3 categories based on the interactive effects of accrual DRP and DIN on % Phormidium cover (Fig. 5A): 1) low DRP (<0.01 mg/L) and elevated DIN (>0.1 mg/L; sites 1 and 5), 2) downstream of STPs with elevated DIN and DRP (>0.1 mg/L and 0.01 mg/L; sites 2 and 6), and 3) low DRP (<0.01 mg/L) and DIN (<0.06 mg/L; site 10). The GAMMs predicted % Phormidium cover >20% when DIN was between 0.09 and 0.27 mg/L and DRP was between 0.004 and 0.005 mg/L (Fig. 5B).

Figure 5.
Figure 5.

Relationship between log(accrual dissolved reactive P) (DRP) and √(accrual dissolved inorganic N) (DIN) and observed % Phormidium cover at each site (A) and predicted by generalized additive modeling (B).

Predicting the flow required to reduce % Phormidium cover to <20%

Quantile regression based on the 90th percentile of the maximum daily mean flow (divided by the long-term median) between the sampling periods predicted that the flow required to reduce % Phormidium cover to <20% would vary markedly among the sites 1, 3, 5, and 10 (~2.5, ~0.7, ~5.4, and ~249× times the long-term median, respectively; Fig. 6A–D). At site 10, the magnitude of required flushing flow predicted based on the 90th percentile was large, so we also explored the 80th and 85th percentiles. These regressions predicted required flushing flows that were ~6.75 and ~16.5× the long-term median, respectively, at site 10 (Fig. 6D).

Figure 6.
Figure 6.

Quantile regressions for mean % Phormidium cover vs the maximum low (obtained from daily means) divided by the long-term median flow on the 7 d prior to surveying (times median flow) at sites 1(A), 2 (B), 5 (C), and 10 (D). Fitted lines show the predicted times median flow required to reduce % Phormidium cover to <20% as derived by 80th, 85th, and 90th percentile quantile regression.

Anatoxin concentrations and variant composition

Anatoxins were detected in variable concentrations in samples from 6 of the 10 sites (Fig. 7A, B, Table S1). The highest median toxin concentration was measured at site 10 (0.20 mg/kg). The highest toxin concentration measured in a single sample was at site 1 (82.4 mg/kg, 17 February 2012). When % Phormidium cover was >~10%, the relationship between cover and toxin concentrations was close to linear (Fig. 7B).

Figure 7.
Figure 7.

A.—Box plot of log(total anatoxin concentration) at each site. See Fig. 3 for interpretation of box plots. B.—Relationship between log(% Phormidium cover) and log(total anatoxin concentration) (sum of homoanatoxin-a [HTX], dihydro-anatoxin-a [dhATX], and dihydro-homoanatoxin-a [dhHTX]). The black line is a local polynomial regression (LOESS)-smoothed local polynomial regression (±SE, gray shading). Colors identify sites and match those in panel A. C.—Percentage of each anatoxin variant present at each site over all samples.

No ATX variant was detected in any sample (Fig. 7C). Over the entire sampling period, HTX was the dominant variant at 4 sites (2, 8, 9, and 10), where it accounted for >50% of the total toxin concentration (Fig. 7C). dhHTX was dominant at 2 sites (1 and 5; Fig. 7C). The dominant variant showed some temporal variability among sites, but no obvious patterns were observed (Table S1).

Total anatoxin concentration and relationship with physicochemical variables and % Phormidium cover

The most parsimonious GAMM for total toxins explained 50.8% of the total deviance and included 5 predictor variables: week of the year, accrual DIN, times median flow, conductivity, and % Phormidium cover (Fig. 8A–E). The effect of week of the year was weak and marginally significant (p = 0.06; Fig. 8A), and followed a similar pattern to that observed for the % Phormidium cover GAMMs. The effects of accrual DIN and flow were less pronounced than those for % Phormidium cover, but were significant (p < 0.001), with the highest toxin concentrations predicted at mid-range DIN and during low flows (Fig. 8B, C). Low conductivity ~50 µS/cm) resulted in high toxin concentrations. Higher conductivity produced a similar pattern but had a smaller effect on toxin concentrations (Fig. 8D). Percent Phormidium cover had a strongly positive and almost linear effect on toxin concentration (p < 0.001; Fig. 8E).

Figure 8.
Figure 8.

Partial plots of the effects of week of the year (A), √(accrual dissolved inorganic N) (DIN) (B), times median flow (C), conductivity (D), and % Phormidium cover (E) on log(x)-transformed total anatoxin (sum of anatoxin-a, homoanatoxin-a, dihydro-anatoxin-a, and dihydro-homoanatoxin-a). Solid lines represent cubic splines (±SE, dashed lines) fitted based on log-normal generalized additive mixed model. See methods for description of x-axis partial effect scale.


Drivers of % Phormidium proliferations

The P. autumnale-complex (Strunecký et al. 2010) is one of the most widespread cyanobacterial taxa in the world and has been identified in a variety of habitats and nutrient regimes including eutrophic and oligotrophic ponds, lakes, and rivers, and on damp soils (Komárek and Anagnostidis 2005, Richter et al. 2009, Wood et al. 2016). Nevertheless, it and other Phormidium species have received international scientific attention only recently because of their ability to form toxin-containing benthic proliferations (Gugger et al. 2005, Bouma-Gregson and Higgins 2015, McAllister et al. 2016). Research on environmental drivers of Phormidium proliferations in cobble-bed rivers has been spatially and temporally limited but has highlighted the importance of variables including river flow, fine sediment, and water-column nutrients (Heath et al. 2011, 2015, Hart et al. 2013, Wood et al. 2015, McAllister et al. 2016). Our study was limited in its geographic extent, but spanned the range of river water-chemistry conditions observed in New Zealand rivers that experience Phormidium proliferations. The GAMMs indicated that week of the year, water chemistry, temperature, and flow are important regulators of Phormidium proliferations.

The basic model for control of periphyton biomass in cobble-bed rivers identifies hydrologic disturbance as the primary regulator, within which nutrients operate by influencing the rate of accrual during periods of stable flow (Biggs 1995, Biggs et al. 2005). Our findings differ from those of other investigators who showed that Phormidium was present in higher abundance than other algal species when DRP concentrations were elevated (Horner et al. 1990), was prevalent in rivers with high DIN and DRP concentrations (Loza et al. 2013), and was out-competed by other cyanobacteria when grown in competitive assays under low-nutrient regimes (Loza et al. 2014). However, these investigators did not identify Phormidium to species level, and different species and strains within genera may vary in their optimal nutrient requirements (Heath et al. 2016). Moreover, the mats studied by Horner et al. (1990) and Loza et al. (2013, 2014) contained many cyanobacterial or algal species. In mats similar to those we observed, Phormidium was volumetrically abundant (Brasell et al. 2015) and formed thick and cohesive mats. This growth structure probably is critical to understanding the prevalence of Phormidium mats in rivers with low-nutrient concentrations.

Sand-Jensen (1983) noted that the relationships between water-column nutrients and algal mats are most important and most evident during early growth stages before mats mature and internal processes, such as nutrient recycling, begin. Several mechanisms might underlie the ability of Phormidium to form expansive mats under low-nutrient conditions; e.g., formation of geochemical conditions (e.g., pH and dissolved O2) within mats that enable release of DRP bound to sediments entrapped in the mat (Wood et al. 2015), internal microbial processes (i.e., N fixation by bacteria in the mats), and possible advection of DIN and DRP from turbulent near-bed flow that could compensate for low-nutrient concentrations. Future investigators should explore nutrient concentrations in Phormidium mats directly to provide insights to this possible disconnection between nutrients in the water column and actual nutrient availability.

Previous investigators have suggested that water-column DIN concentrations >0.10 to 0.20 mg/L promote Phormidium proliferations (Wood and Young 2012, Heath et al. 2015, McAllister et al. 2016). This threshold held true for 4 of the 5 sites in our study that experience Phormidium proliferation. However, DIN concentrations at site 10 are much lower than the 0.10 mg/L threshold, but the site experiences prolonged proliferations (Fig. 2). Reasons for this phenomenon are unknown but may be related to the higher intensity of flow required to remove Phormidium mats at this site than at other sites (see Discussion below), which would provide longer periods for mats to become established or would result in larger starting inoculums after flushing events. These relic patches might already contain functioning microbial communities, e.g., N-fixing organisms, that would enable rapid mat expansion.

Conductivity had a strongly positive effect on % Phormidium cover. Conductivity is often correlated with periphyton biomass and is strongly related to catchment geology and water source in New Zealand rivers (Biggs 1990, 1995). It also is linked integrally to river flow, and generally increases during periods of low flow because of lack of dilution from water with lower ionic content, a trend also evident in our data set (data not shown). Many of the major ions contributing to conductivity, e.g., Ca2+, Cl, Na+, K+, and Mg2+, influence important metabolic processes or enhance/suppress growth in cyanobacteria (Seale et al. 1987, Parker et al. 1997, Carneiro et al. 2011). Individual ions/compounds that contribute to conductivity may play a role in regulating Phormidium accrual, but further research is required to test this suggestion.

Our analysis highlighted the importance of week of the year, with higher % Phormidium cover predicted in summer and early autumn than at other times. At most sites, this period coincides with periods of stable flow and elevated temperature, variables that also were related to elevated cover. Our model predicted greater % Phormidium cover when water temperature was >15°C. A recent analysis of a 20-y data set of benthic algae in 2 pristine Norwegian streams also forecast greater Phormidium cover with increasing temperatures (Schneider 2015). River flow equal to approximately half the long-term median flow was projected to favor % Phormidium cover, with a ‘hump-back curve’ apparent (Fig. 4E). A habitat suitability criteria model developed for Phormidium based on velocity, depth, and substrate type showed a similar nearly parabolic curve, with highest cover observed at moderate velocities (1 m/s; Heath et al. 2015). Gas bubbles caused by photosynthesis are often visible in Phormidium mats. Bubble formation is more likely to occur under low flow because diffusion of O2 will be slowed by the thick boundary layer (Hawes et al. 2014). Therefore, Phormidium growth in slower flow may be coupled with more frequent autogenic detachment, which would reduce biomass accrual in low flows (Boulêtreau et al. 2006). Thus, flows that are sufficient to enhance nutrient and gas flux but insufficient to generate shear stress that results in biomass loss are postulated to be optimal for Phormidium accrual (sensu Biggs and Thomsen 1995).

Intensity of flushing flow required to reduce Phormidium proliferations

High-energy flow events ‘reset’ periphyton communities (Clausen and Biggs 1997). Some management authorities involved in issuing human health warnings related to Phormidium have used ‘length of time since a flushing flow’ as an early warning indicator of elevated risk. For example, Milne and Watts (2007, p. 24) advocated for the criteria “…that if there has not been a flushing flow (defined in their study as 3× the median flow [sic]) for two weeks and low river flows (set at the lowest 10th percentile flow for each river).” Exceeding these criteria triggered weekly surveys of sites prone to Phormidium proliferations.

We used quantile regression to predict the magnitude of flushing flow required to reduce % Phormidium cover and showed that this value differed among sites. Identifying the reasons for these differences is challenging. Phormidium mats are likely to be removed during periods of elevated flow through high shear stresses, abrasion by mobilized sediments, and grinding action of tumbling gravel/cobble substrata (Grimm and Fisher 1989, Horner et al. 1990, Biggs et al. 1999, Francoeur and Biggs 2006). Heath et al. (2015) showed in an intensive field survey that the best predictor for Phormidium presence was stable substrates unlikely to be disturbed by flood flows. Substrate heterogeneity also might be relevant if cobbles or boulders provide more refuges (i.e. cracks and crevices) during flushing flows and enable faster recolonization (Bergey 2005, Murdock and Dodds 2007). In addition, Phormidium mats undergo a series of developmental stages as they grow. These stages include initial attachment, maturation, and dispersal (McAllister et al. 2016). Dispersal can occur when bubbles of O2 become trapped in the mats and cause them to detach or slough off the substrate. Early-stage mats are more firmly adhered to substrates and, therefore, are likely to require more intense flushing to remove them.


Total anatoxin concentrations and the relative composition of structural variants showed high spatial and temporal variability among and within sites, an observation consistent with previous studies (Heath et al. 2011, Wood et al. 2010). Toxin concentrations were generally low compared to those in a nationwide data set (McAllister et al. 2016), and the complete absence of toxins in mats from sites 3, 6, and 7 is noteworthy.

An important caveat when assessing relationships between toxin concentrations in Phormidium mats and physicochemical variables is that mats contain a mixture of toxic and nontoxic genotypes (Heath et al. 2010, Wood et al. 2010), other organisms (Brasell et al. 2015), and inorganic debris (Wood et al. 2015). These variables can vary among sites. Moreover, even within toxic genotypes, the relative amount of toxin produced can vary 100-fold (Wood et al. 2012). This variability limits our ability to infer relationships among physicochemical variables and actual toxin production or shifts in genotypes. However, our analysis provides some indication of periods or conditions associated with greatest toxin content of mat samples and, therefore, periods of highest risk to river users.

We identified 3 prominent relationships in our analysis. First, total toxin concentrations were higher when conductivity was low than when it was high, a pattern opposite that observed for % Phormidium cover. A similar result was noted in a culture-based study where the anatoxin quota was lowest in high N and P conditions (Heath et al. 2016). The authors suggested that this result supports the growth–differentiation balance hypothesis, which states that actively dividing cells are less likely to produce secondary metabolites than are nondividing cells. Second, the relationship between % Phormidium cover and toxin concentrations was nearly linear, particularly when % Phormidium cover was >10% (Figs 7A–C, 8A–E). This relationship differs from findings of other studies that suggested no relationship between these variables (Wood et al. 2010, Heath et al. 2011) and culture-based studies that showed peak toxin quotas during the initial growth phase (Harland et al. 2013, Heath et al. 2016). Third, toxin concentrations were lower when DIN was <0.2 mg/L, a pattern similar to that observed for % cover and consistent with the observations by Heath et al. (2016) but not with the growth–differentiation balance hypothesis. Development of molecular techniques (e.g., quantitative polymerase chain reaction) that enable enumeration of Phormidium cells in environmental samples and distinguishing toxic and nontoxic genotypes is essential to allow insights into how these and other variables are related to toxin production and genotype succession. This ability may ultimately assist in understanding the ecological role of anatoxins in these mat communities.

Author contributions: SAW helped with study design and led the analysis of the data and manuscript preparation. JA, AW, KD, RGY, and IH assisted with analysis of the data and preparation of the manuscript. LB assisted with study design, coordinated and undertook the field work, and helped with preparation of the manuscript.

We thank Jon Roygard, Carol Nicholson, Michael Patterson, Josh Markham, and Manas Chakraborty (Horizons Regional Council) for sample and data collection, advice, and support throughout this project. We thank Michael Boundy and Bryan Stokes (Cawthron) for assistance with LC-MS analysis, and Eric Goodwin (Cawthron) for statistical advice. SAW, AW, and RGY thank the New Zealand Ministry for Business, Innovation and Employment Cumulative Effects program (CO1X0803) for funding. We acknowledge the New Zealand Ministry of Environment for funding and support.

Literature Cited