Modeling total microcystin production by Microcystis aeruginosa using multiple regression

Cyanobacteria that produce toxic or taste and odor metabolites have become a widespread and critical problem in aquatic environments around the world (Saraf et al. 2018; Zhang et al. 2018). Intensive agricultural practices, the release of untreated domestic and industrial wastewater into water bodies, and climate change have boosted the growth rates of cyanobacteria biomass (Jacinavicius 2015). The annual global estimates of the socioeconomic impacts of cyanobacteria metabolites can reach billions of US dollars when considering water monitoring, treatment and analysis as well as adverse impacts on recreational use and fishing (Wiltsie et al. 2018).

One of the most studied species of cyanobacteria in lakes is Microcystis aeruginosa – MA (Chaffin et al. 2018). It can produce a group of toxic monocyclic heptapeptides called microcystins – MCs (Minasyan et al. 2018), the most common toxin found in natural waters (Walls et al. 2018). In addition, these cyanobacteria are problematic because it rapidly increases biomass on the water surface, limiting sunlight penetration, reducing dissolved oxygen concentration and making water treatment more laborious and costly (Ger et al. 2010). MCs can accumulate in the trophic network and cause harmful effects to mammalian health, especially to the liver, inhibiting essential enzymatic functions of the hepatic tissue (Yoshizawa et al. 1990; Massey et al. 2018; Mohamed et al. 2018). As a consequence, the World Health Organization proposed a maximum concentration of 1 μg L⁻¹ in finished water (WHO 2003, 2011), limit adopted by Brazil and Uruguay. Countries such as Australia and Canada have adopted MCs limit concentration of 1.3 and 1.5 μg L⁻¹, respectively (Canada 2002; NHMRC & NRMMC 2011).

The release of MCs occurs preferentially due to cellular senescence, when membrane lysis occurs naturally, rather than continuous excretion (Zi et al. 2018). Lysis may also occur during the application of chemicals for the treatment or control of cyanobacteria, especially algaecides (Chorus & Bartran 1999). After cyanobacteria lysis, the intracellular contents are dissolved into water, making them available to organisms (Schmidt et al. 2013) and accumulating in the aquatic systems (Pereira et al. 2018). Dissolved MCs are more difficult and costly to remove in conventional water treatment plants (Van Apeldoorn et al. 2007) since extra treatment steps must be incorporated into the treatment train (Hamilton et al. 2013). In addition, to minimize risks, frequent analysis of MCs and cyanobacteria cell density is required, further increasing the costs and time consumed. Frequently, by the time analysis confirms the presence of toxins, the finished water has already been delivered to consumers.

In order to overcome the high costs and time required for cell count and cyanotoxin analysis, some attempts to develop predictive methods were carried out to evaluate how environmental variables are connected to the growth of toxic cyanobacteria and toxin production (Bortoli et al. 2014). Harris & Graham (2017) observed that pigments, iron and dissolved oxygen presented a high correlation with MCs concentration. Marion et al. (2012) developed a multivariate logistic regression model using real-time measurements of phycocyanin in vivo and Secchi depth to evaluate the risks of MCs concentration exceeding 4 μg L⁻¹ in waters of a recreation pond. The final model reached high statistical significance as well as a good calibration, proving to be effective in identifying the risks involved with bathing. McQuaid et al. (2011) used in vivo analyses of phycocyanin in an attempt to evaluate the abundance of cyanobacteria and the associated MCs concentration in water sources. The results of this study demonstrated that there is a potential for the prediction of cyanobacteria biovolume and the concentration of MCs using phycocyanin probes.

Phycocyanin is a water-soluble and strongly fluorescent compound, making it viable to detect it with portable instruments such as probes (McQuaid et al. 2011). Probes based on phycocyanin fluorescence are applied for the rapid assessment of cyanobacteria in drinking water sources (Bastien et al. 2011; McQuaid et al. 2011), offering an alternative to monitoring cyanobacteria and hence cyanotoxins. However, most probes still have limitations when transforming pigment concentration signals into other cyanobacterial parameters, such as biomass and cell densities, and only perform accurately within very limited concentration ranges (Glibert et al. 2018). In this sense, the development of an efficient model that precisely correlates pigment concentration and cyanotoxins would be beneficial to the water industry.

Although it has been demonstrated that the concentration of MCs is related to the number and the growth phase of cells (McQuaid et al. 2011), studies comparing the productivity of metabolites at different stages of growth and with different MCs variants are still limited. Crettaz Minaglia et al. (2017), for example, mathematically modeled M. aeruginosa growth and microcystin-LR production in BG11 media at different temperatures. The authors established that there is a relationship between some environmental variables and one variant of MC. However, a mathematical model that estimates the concentration of total microcystin (McT) using easily measurable parameters, such as pigments, has not yet been found in the literature.

In this context, the main objective of this work was to evaluate the production of intra- and extracellular MCs throughout the life cycle of an MA culture under controlled conditions and to correlate McT concentration and cell age and density, concentrations of chlorophyll-a, phycocyanin and pheophytin in order to develop a mathematical model to estimate the concentration of McT using few parameters with relatively easy and rapid analysis.

The suspension using an MA culture and an ASM1 medium was prepared in one batch and, after that, it was subdivided into five individual flasks of 1 L each. Each flask was completely used for that week's analyses, so no changes in the initial volume occurred. The flasks were submitted to the same conditions. The M. aeruginosa strain (CIAR 03) was cultured in 1 L flasks in an ASM-1 medium (Gorham et al. 1964) using 12 : 12 h (L : D) photoperiod at 22 μ E m⁻² s⁻¹ (Digital Lux Tester YF-1065), pH 8, 25 °C (±1 °C) and under continuous aeration. In each of the subsequent five weeks (cell age – CA), one flask was used in full for analysis of the following parameters: MA cell density (Dens), chlorophyll-a (Chlo.a), phycocyanin (Phyc) and pheophytin (Pheo), called the explanatory variables, and Total Microscystin – McT, the response variable. The experiment was carried out in triplicate for 5 weeks, totaling 15 analysis. Aliquots of 1 mL were used for cell density analysis. For the determination of chlorophyll-a (Chlo.a) and pheophytin (Pheo), an aliquot of 100 mL was used; for phycocyanin (Phyc), an aliquot of 100 mL was used and for total microcystin – McT, an aliquot of 100 mL was used.

For microcystin analysis, the aliquot of 100 mL was filtered through 110 mm GF/C discs (Whatman, Maidstone, Kent, UK). The GF/C discs with the cell-bound or intracellular microcystin were freeze–thawed and placed in glass beakers containing 10 mL of methanol (100%) (J.T. Baker – Grade UV-HPLC) to extract the intracellular content for 1 h at room temperature (25 °C). The extract was removed and evaporated at 50 °C in a vacuum until dry. The residue was resuspended in 1 mL of methanol (80%) prior to analysis by high-performance liquid chromatography (HPLC). For the extracellular microcystin, the filtered samples (100 mL) were passed through C18 Cartridges (SupelClean LC-18 – SUPELCO 1 g and 3 mL) using a vacuum manifold system (Varian Analytical Instruments, San Fernando, CA, USA). The cartridges were conditioned before sample concentration using 10 mL of methanol (100%) (J.T. Baker – Grade UV-HPLC) followed by 10 mL of water. Samples were applied to the cartridges at a flow rate of 1 drop min⁻¹. After the sample was passed through, the cartridge was washed with 10 mL of 10, 20 and 30% v/v aqueous methanol in sequence. The air was drawn through the cartridge for about 10 min to minimize the amount of aqueous methanol eluted from the previous step. Then, the cartridge was eluted with 0.8 mL of methanol (80%). The analysis was performed using high-performance liquid chromatography (Agilent Technologies, model 1260 Infinity), with a reverse-phase column (Waters Symmetry C18, 5 μm, 2.1 × 150 mm) at 40 °C with a DAD detector set at wavelengths 238 and 222 nm. The mobile phases used were (i) 0.05% trifluoroacetic acid – TFA acid (Dynamic – Grade UV-HPLC) in water and (ii) 0.05% TFA in acetonitrile (J.T. Baker – Grade UV-HPLC) with a flow rate of 0.5 mL min⁻¹ and a total run time of 35 min. Gradient elution was applied using the following program: linear increase from 25 to 100% B within 22 min, hold at 100% B for 2 min and subsequently decrease to 25% B within 11 min. A calibration curve for each variant (Mc-LR, Mc-LA and Mc-LY) was produced using standards provided by the Cyanosol Laboratory (UK). The equipment's limit of detection (LOD) and quantification (LQD) was 0.25 and 0.5 mg L⁻¹, respectively, but the method's limit of detection (LOD) and quantification (LQD) was based on those limits and the degree of sample concentration detailed previously. This method was adapted from Lawton et al. (1994).

Total microcystin concentration (McT) comprised of the sum of the Mc-LR, Mc-LA and Mc-LY concentrations (intra- and extracellular). For the cell quota calculation, the total intracellular microcystin concentration (sum of the intracellular Mc-LR, Mc-LA, and Mc-LY) was divided by the number of cells or biomass, in the same base volume.

After preserving MA samples with Lugol, cell density (Dens) at different growth stages was measured using a Sedgewick-rafter chamber and an inverted optical microscope with a magnification of 40× (Zeiss, Model Vert. A1; APHA 2012). For the determination of cell density, each unicellular organism was considered as an individual (Yoshida et al. 2006).

The spectrophotometric method (10200H), proposed in APHA (2012), was used for the determination of chlorophyll-a (Chlo.a) and pheophytin (Pheo). Samples were analyzed at wavelengths of 664, 665 and 750 nm, respectively. Phycocyanin (Phyc) was analyzed using the spectrophotometric method at 665 and 750 nm wavelengths proposed by Bennett & Bogorad (2011), with adaptations (Genesys 10-S, Thermo-Scientific, USA).

The hypothesis of data normality and homogeneity of variance was verified using the Shapiro–Wilk (P_S-W) test and the Levene test. For comparisons between the analyzed values at each week after using the analysis of variance (ANOVA) test, the Tukey test was performed for samples with normal distribution, and pairwise comparisons using the Wilcoxon signed-rank test with Bonferroni adjustment for samples without normal distribution. The Spearman correlation coefficient (R_s), as well as its intensity (weak, moderate and strong) and significance, was determined to verify the degree of association between all variables. For this study, it was considered a strong association when |R_s| ≥ 0.7 and moderate when 0.7 > |R_s| ≥ 0.5 (Zar 1996).

The explanatory variables (CA, Dens, Chlo.a, Phyc and Pheo) that correlated significantly with the response variable (McT) were incorporated into the model using the Stepwise method. The regression model that contained easily determined explanatory variables was preferred among the models with the highest Adjusted R². The selected regression model was then evaluated by both normality and homoscedasticity of the residues (evaluated graphically – Figure 4) and coefficient significance. All the statistical tests were performed using a significance level of 5% (α = 5%) and the RStudio program (R Core Team 2018). Details of the data analyses mentioned are available in the Supplementary material.

Until the fourth week of the experiment (28 days), there was a significant (p < 0.05) increase in cell density (Figure 1). However, no significant difference (p > 0.05) was observed between the fourth and fifth weeks indicating that the culture reached the stationary phase. There was no need for additional points after that since they would not improve the model precision. Bortoli et al. (2014) studied the growth of a MA strain isolated from the Billings Reservoir (Brazil) under different nutritional conditions and observed that the stationary phase was reached after 21 days, differing from our experiment probably because of the different experimental conditions and MA strain. According to Bittencourt-Oliveira et al. (2005), MA growth and microcystin production vary according to a biochemical activity regulator, which is different for each strain.

The MA strain cultivated produced Mc-LR, Mc-LA and Mc-LY which were detected in the chromatograms at 4.7, 16.1 and 19.7 min, respectively. The intracellular Mc-LR variant was predominant during the whole duration of the experiment (Figure 2). In the fourth week, intracellular Mc-LR concentration was 70 μg L⁻¹, while extracellular concentration was 1.25 μg L⁻¹, close to the concentrations observed by Turner et al. (2018). On the other hand, Almuhtaram et al. (2018), studying a reservoir used for water supply, observed that intra- and extracellular MC concentrations were approximately 4 and 0.5 μg L⁻¹, respectively. This difference between intra- and extracellular proportions may be explained by the fact that extracellular MC modifies the intracellular MC synthesis (Schatz et al. 2007). The authors suggested that Microcystis cell lysis enhances MC concentration in the medium, inducing accumulation of mcyB (genetic subunits involved in the synthesis of MCs) and increasing the intracellular production of MC in the remaining intact cells. With the release of MCs and other oligopeptides, microgynin and micropeptin, the intact population responded by increasing the intracellular production of MC (Schatz et al. 2007), a feature more easily identified in axenic cultures.

In natural waters, MCs are preferentially found inside the MA cells (Westrick et al. 2010; Merel et al. 2013) and released in greater amounts during the cell senescence or death phase (Merel et al. 2013; Ye et al. 2017). This behavior may also be observed in MA cells submitted to controlled conditions (Neilan et al. 2013). The concentration of MC may be correlated with the cell growth rate, with the highest concentrations of MCs observed at the end of the exponential phase of growth. In the stationary phase, concentrations of MCs tend to stabilize and then to decrease in the cell senescence phase (Repka et al. 2004). This reduction of MC concentration can be explained by the reduced activity of enzymes complex NRPS/PKS, responsible for the biosynthesis of secondary metabolites (Nikolouli & Mossialos 2012).

According to Crettaz Minaglia et al. (2017), Mc-LR cell quota is an important parameter in determining the cellular physiology of MC production. The authors also observed average cell quota values (20 ± 10 fg cell⁻¹) during MA growth similar to the ones found in our study (20 ± 10 fg cell⁻¹).

The cell quota, a relation between total intracellular MC concentration and cell density or biomass, showed no significant difference (p > 0.05) along all growth periods, suggesting that the cell quota was constant throughout this experiment (Table 1). Long et al. (2001) found that to be a common behavior during the growth phase of MA. This linear correlation between MC production and the cyanobacterial growth rate has also been reported by Orr & Jones (1998) and Bortoli et al. (2014). Puddick et al. (2019) evaluated 187 MC quotas in New Zealand and found that the highest value detected was 13 pg cell⁻¹ (equivalent to 13,000 fg cell⁻¹), much higher than the quota found in our study (32 fg cell⁻¹). This may be due to the fact that Microcystis strains found in natural environments have shown the capacity to produce higher cell quota than strains cultivated under controlled conditions (Wood et al. 2012).

Therefore, it is feasible to consider that the number of Microcystis cells is an indication of the concentration of MCs in both cultured and natural waters. It is important to remember, however, that the cellular quota can vary between different Microcystis taxa, environmental conditions and in different natural reservoirs (Salmaso et al. 2014).

The first line of Figure 3 shows that all explanatory variables were strongly correlated (p < 0.001) to the response variable (McT). In the left column, the scatter plots of the McT values versus the covariates (explanatory variables) are presented with their respective linear regression lines. The visual verification of the scatter plots can be used to estimate the McT values from each explanatory variable and to evaluate how closely the experimental results are to the modeled values. Looking at the linear regression adjustment lines (graphs below the main diagonal of Figure 3), it is noted that the explanatory variables alone do not explain the McT production adequately because the vertical distances between the points and the line are large. This justifies the search for a more elaborate model to estimate McT concentration. Thus, according to Figure 3, the most relevant variables are CA, Dens, Chlo.a, Phyc and Pheo with R_s > 0.7, indicating that these variables could be used to estimate McT.

In the search for a mathematical model capable of estimating McT concentration based on an explanatory variable or a set of them, a multiple regression of the data was performed. Multiple regression models were generated, with their respective correlation coefficient (Adjusted R²), involving all possible combinations (Table 2).

Models containing only pigments presented Adjusted R² values lower than those containing Dens and CA (Table 2). Beversdorf et al. (2015, 2018) showed that there is a strong correlation between MA cell density and MCs production in isolated MA cultures, indicating that it is possible to estimate MCs concentration using the number cells. Pereira et al. (2018) demonstrated that increased cell density of MA induced an increase in chlorophyll and MCs.

On the other hand, although the model selected to describe the response variable (McT) should have high Adjusted R², for simplicity's sake, it is preferable that it requires the smallest number of variables and that these variables are easily analyzed. According to Table 2, the presence of a single variable in the mathematical model does not explain the production of McT satisfactorily. This can also be observed by the small Adjusted R² of single-variable models (Models 12, 19, 28, 30 and 31 – Table 2). Determining cell age (CA) and density (Dens) in environmental samples is not an easy task, and therefore, models based only on pigments were selected since they can be analyzed by simple multiparametric probes. On this basis, the models selected for determining McT were models 25 and 26 with an Adjusted R² = 0.824 and 0.814, respectively. The adjusted coefficient of determination, also called Adjusted R², is a criterion used for the selection of an optimal model with the largest Adjusted R² with a reduced set of explanatory variables. This choice is corroborated by studies that have shown a strong correlation between pigments and McT concentrations in environmental samples (Lee et al. 2000; Downing et al. 2005; Francy et al. 2016).

Table 3 shows the coefficients of the selected multiple regression models, as well as p-values and upper and lower limits for each estimated coefficient. Models 25 and 26 were then evaluated using the linear regression method to choose the most appropriate one. From Table 3, it is possible to observe that model 26 is not the most appropriate because Chlo.a p-value is higher (p-value = 0.546) than the adopted significance level (α = 5%) and the confidence interval for each coefficient goes through zero. On the other hand, model 25, which includes the concentration of Phyc and Pheo pigments, is a reasonable choice since its p-values are below the adopted significance level (α = 5%) and the confidence interval for each coefficient does not go through zero. Arranging the variables and their respective coefficients, Equation (1) is obtained.

Figure 4 shows the analysis of the residues obtained from the model selected (25) for the prediction of McT. According to Figure 4(a), the residues were found to be normally distributed, a required condition for the regression model validity (Mann 2010; Zar 2010). Figure 4(b) shows the distribution of the residues, indicating a good fit since the residues are randomly scattered around zero. This graph is important for verifying the assumption that residues are randomly distributed and have constant variance (Mann 2010; Zar 2010). Therefore, there are no indications of a violation of assumptions to certify the quality of the model according to Cook & Weisberg (1982).

A comparison between the actual values obtained in the experiment and those predicted by the model selected for the McT variable (Equation (1)) is presented in Figure 4(c). Based on this figure, it is possible to state that the evidence of normality and homoscedasticity of the residues was confirmed by adopting an Adjusted R² of 0.824. Thus, the multiple regression model proved to be valid and was able to explain 82.4% of the McT variable.

The multiple regression model was able to predict McT concentration in cultured MA using easily measured variables such as pheophytin and phycocyanin. Although the model containing only pigments presented an Adjusted R² value lower than the one with all the studied variables, it was considered a more practical solution and showed good precision.

Although differences between calculated and measured McT concentrations were less than 5%, demonstrating that the model can be applied for our cultured MA, its use for environmental samples still depends on further investigations. The combination of sensor technology and accurate mathematical modeling can be an important tool to alert for toxins in water and to the need for complementary treatments such as activated carbon application or, in extreme cases, the complete stop of the water distribution.

The Supplementary Material for this paper is available online at https://dx.doi.org/10.2166/aqua.2020.128.