Application of model fitting technique to enhance bacterial regrowth potential (BRP) measurement for drinking water supply monitoring

In drinking water distribution systems, microbial regrowth can lead to various problems, such as low aesthetic water quality (unwanted odour and colour) and issues related to infrastructure and operation (corroding pipe materials which can lead to difficulty in maintaining disinfectant residual (Kaplan et al. 1992; Camper 2014). One of the most important nutrients for microbial regrowth is the biodegradable component of organic carbon in the water distribution system. The readily biodegradable or assimilable organic carbon (AOC) was identified by van der Kooij et al. (1985), comprising only a very small proportion (0.1–0.9%) of total organic carbon in drinking water (Volk & LeChevallier 2000; Chowdhury 2012; van der Kooij et al. 2017; Huang et al. 2020). However, these concentrations of organic carbon are not high enough to be detected accurately by chemical analysis (Withers & Drikas 1998a; Ross et al. 2013). Hence, several microbial assay methods have been developed to measure this substrate as a parameter to monitor the microbial regrowth potential in drinking water distribution systems (van der Kooij 1990; Hammes & Egli 2005; Elhadidy et al. 2016), especially flow cytometry-based assay was widely applied in recent years (Elhadidy et al. 2017; Miller et al. 2020).

The bacterial regrowth potential (BRP) method was developed by Werner & Hambsch (1986) and was commonly used to determine the risk of microbial regrowth in drinking water distribution systems. This method utilises indigenous bacteria retrieved from water samples as an inoculum and microbial regrowth is assessed indirectly by the change in turbidity over time (Withers & Drikas 1998b; Ross et al. 2013; Pick et al. 2019). The correlation between turbidity and total cell number is high under different conditions in the bacterial growth phase and has been demonstrated by Werner & Hambsch (1986). The advantage of this method is that it incorporates indigenous bacteria which can broaden the substrate spectrum range to the nutrients present in the water sample (Withers & Drikas 1998a; Ross et al. 2013). Ross et al. (2013) have shown that BRP measurement is highly sensitive, and using the BRP method would result in the detection of a higher AOC concentration than using the ‘van der Kooij’ method.

Normally, a microbiological growth curve has a lag phase, maximum growth rate and asymptotic phases to describe its behaviour (Zwietering et al. 1990).

Hambsch & Werner (1993) determined the growth rate in the BRP method using Equation (2). However, the exponential phase is estimated with uncertainty, due to the difficulty of not being able to consistently determine the initial and final turbidity values on the curve and it can lead to relatively inaccurate determination of the growth rate. To develop a more consistent approach and improve accuracy, both the lag phase and stationary phase are required to be considered as part of the maximum growth rate calculation. Therefore, using the statistical modelling approach to develop a mathematical curve to describe the behaviour of BRP growth curve should improve the consistency of measurements of the maximum growth rate, lag phase and AOC values and enable more reliable comparisons between water samples (van der Kooij et al. 2017).

The Weibull regression has been used to model indeterminate growth for plants, fish and molluscs. López et al. (2004) and Swintek et al. (2019) indicated the Weibull model achieved a better correlation coefficient for microbial growth experimental data than the commonly used Gomperts model/sigmoid curves (Tjørve & Tjørve 2017).

The aim of this study was to develop a mathematical model using BRP data to improve understanding of the behaviour of bacterial regrowth in drinking water distribution systems. The objectives of this study are (1) to use experimental data to determine best-fitted models using the curve-fitting software, (2) to understand bacterial growth curve by developing/extracting microbial growth parameters from the determined model and (3) to validate and evaluate the developed model and compare it with the Weibull model.

A series of BRP experiments (82 sets) were conducted using water samples obtained from six different distribution systems to generate a database and use this to develop a mathematical model to further study the BRP technique.

Six drinking water distribution systems, Anstey Hill, Barossa, Myponga, Happy Valley, Hope Valley and Little Para, in the South Australian metropolitan area were selected to conduct this study. The drinking water in these systems was supplied by the corresponding treatment plant using conventional treatment processes (coagulation/flocculation, sedimentation or floatation and then sand filtration). Water samples were collected at the water treatment plant from the treated water tank (after filtration and prior to chlorination) and from the end of distribution system (consumers’ taps).

A Monitek Turbidimeter (Monitec mAOC Analyser Model 251-4 mAOC, Germany) was used for the acquisition of the BRP data. There were four units of this instrument, each one had a four-channel turbidimeter (total of 16 channels) connected with an external computer and automatically registered readings and analysing growth curve through an accompanying computer program. The accompanying program exported the growth curve data (raw data) to a text format file, as time and turbidity (two columns of data).

Water samples were collected in glass bottles with a quenching agent (4% of sodium thiosulphate) added and transported back to the laboratory with cold pack. The samples were then filtered through pre-rinsed glass fibre membrane filters. The sterilised (autoclaved) inorganic nutrient salt solution prepared using analytical grade chemicals, ammonium chloride, calcium nitrate, calcium chloride, magnesium sulphate, ferrous sulphate, sodium silicate nanohydrate, aluminium sulphate and potassium dihydrogen orthophosphate, with Milli-Q water was added to the water samples to ensure AOC was the limiting nutrient (Withers & Drikas 1998b). Subsequently, pre-rinsed polycarbonate membranes were used to filter the sample, and then these membranes were placed in isotonic saline solution to obtain natural consortia as inoculum (Withers & Drikas 1998b).

To inoculate the sample with the indigenous bacteria, the inoculum was washed off the polycarbonate membrane and was added into a cuvette, then placed in a channel of the instrument until the turbidity in the range of 0.02–0.04 ppm SiO₂. After the experiments were started, turbidity was automatically measured and recorded by the computer program every 30 min. The experiments were completed when growth curves reached a stationary phase, and it is determined by the experimental operator. The instrument was located in a constant temperature room with temperature maintained at 22 °C (Withers & Drikas 1998a, 1998b).

The curve-fitting software TableCurve 2D version 5.01 (Systat Software, Inc., Richmond, CA, USA) was selected to analyse data and determine ideal empirical equations representing the bacterial growth curve for different types of water. TableCurve 2D was applied to find the best models by fitting the experimental data to its database of a wide range of equations in its equation pool. To simplify the process, the automated fitting feature was selected, and 3,665 built-in equations were used to provide the ranks of R², degrees of freedom-adjusted R², maximum error and F-statistic.

Two-column data namely time and turbidity were processed using TableCurve 2D. Every set of experimental data was processed by TableCurve 2D and fitted to the selected Equation (4) (named as BRP-Hill model hereafter) to obtain the coefficient of the equation for the experimental data. For a few sets of data, a rational coefficient for parameter b from TableCurve 2D was not able to be determined, because the software output for the highest R² value leads to coefficients of parameter b being extremely high. Microsoft Excel Solver was applied to adjust coefficients for these sets of data. Although the R² values after adjustment were slightly lower than the TableCurve 2D fitted, the coefficients of parameters provided for accurate fitting.

The goodness of fit of the BRP-Hill model (Equation (4)) was statistically assessed by the coefficient of determination (R²) values.

When modelling some BRP data sets by the BRP-Hill model (Equation (4)), an issue was noted where differences between the model fitted AOC value and the instrument built-in determined value were relatively high. This is because the lag phases of some BRP experiments were relatively long and they did not reach the stationary phase during test. In other words, the experimental operator stopped the experiments too soon at the predetermined time, which followed the routine procedure. Hence, a mathematical approach (second derivative) was applied.

To validate this model, this study compared the instrument built-in method for AOC values with model-determined AOC values using regression by plotting instrument built-in vs. model-determined and both the slope and R² were used, and errors were enumerated for two types of water, treated water and tap water.

To describe a bacterial growth curve, the characteristics of the curve need to be studied. Four phases of the growth curve need to be evaluated: (1) the lag time (T_Lag) (the intercept of lower asymptote line and tangent line that passes through the inflection point), (2) the maximum growth rate (μm) (the slope of the tangent line through the inflection point), (3) lower asymptote (start value of the curve) and (4) upper asymptote (maximum values that the curve is able to reach) (Zwietering et al. 1990). In this model, those four phases are able to be calculated (as shown in Figure 1(c)), and this algorithm is according to Zwietering et al.’s (1990) method.

Examples of growth curves of two water samples showing details of the experimental data and the 95% confidence interval of the fitted growth curves are shown in Figure 2. The model-determined AOC values, maximum growth rates and lag times are shown in Figure 2. These water samples were collected from the Myponga drinking water distribution system (treated water from the plant and tap water from customer taps, respectively) in winter. This model provided good fitting to the data of these waters, and the R² values are all greater than 0.99.

The BRP-Hill model (Equation (4)) provided overall good fitting to the data for plant treated water (after filtration and prior to chlorination) and tap water samples, shown in Table 2. Although tap water's minimum R² value was smaller than 0.900, only one set of tap water's R² values are smaller than 0.900. These data indicate the BRP-Hill model (Equation (4)) fits suitably the experimental data for the BRP of the waters.

With pre-set analysis times, some cases of BRP analyses with relatively long lag phases were stopped before the plateau phase was reached. Hence, the original method for AOC value determination and model-determined values have unacceptable differences. For instance, Figure 3(a) shows an example of a large difference between the method built-in AOC value and the model-determined AOC value of 233.8 and 293.1 μg/L ACE for a data set acquired using Equations (5) and (3), respectively. These two AOC values were calculated from the model-determined turbidity difference (b) and experimental turbidity difference. It can be seen that the end point of predicted curve (Pred-Tur – solid curve)'s vertical axis value is 0.76, and the upper asymptotic value is 0.93 (coefficient of parameter a plus coefficient of parameter b). Hence, the model predicted curves have not reached the plateau phase, and that was the reason for the relatively large difference between the AOC values. It can be inferred that for such experiments the full monitoring period to plateau stage was not completed, and that made experimental AOC values lower than model-determined values, which was based on a plateau being attained.

When parameter d is smaller than −1, the BRP-Hill model (Equation (4)) is a sigmoidal curve with an upper asymptote (Panik 2014), and it has a bell-shaped first derivative (Basterretxea et al. 2004). As seen in Figure 4, for Example 2, that when the BRP experiment is run over a longer time to reach the stationary phase, the circular dot on x-axis (the cross point of tangent line 2 and x-axis) is a distinct point where the curve did not reach the plateau (when growth curve (top solid line)'s x-value was smaller than circular dot's x-value. Hence, if the experiment stops before this point (cross point of second tangent line and x-axis), the full monitoring period is not completed.

The second inflection point of first derivative curve was selected as the first point to calculate the average AOC value. Some sets of data were investigated, and the first derivative values of this point are not very small, which means there was modest gradient of this point on the model curve. Subsequently, calculating the average turbidity value from x″ to the end point, and obtain the difference between this average value and parameter a, then dividing by 3, gives a result that equals 198.0 μg/L ACE. The instrument built-in AOC value equals 212.0 μg/L ACE, the model-determined AOC value equals 234.7 μg/L ACE, and hence the result from this approach was the closest value to the instrument built-in AOC value.

The BRP-Hill model (Equation (4)) was investigated to determine the accuracy of how it represented the experimental data. The 82 sets of data were analysed by comparing the instrument built-in values and model-determined values in terms of AOC. There were seven sets of data's model-determined AOC values being calculated by applying the alternative algorithm. The model's parameter values can be utilised for determining the AOC value on the condition where experiments were stopped after arriving at the plateau phase, and the alternative algorithm was used to determine AOC values for experiments that were stopped before reaching the plateau phase. The linear relationships between the AOC value obtained from instrument built-in and the AOC value determined from the model's parameter values for the two types of water are presented in Figure 5. The slope values of the two regression lines are close to 1, and the R² values of the plant treated water and tap water were 0.97 and 0.98, respectively. These values indicate that the results determined by the instrument built-in method and model-determined results are in good agreement.

In this study, the Weibull model frequently used in predictive microbiology (Cunha et al. 1998) was selected for comparison for the proposed BRP-Hill model (Equation (4)). In this study, Gomperts model and another frequently used logistic model (Fujikawa et al. 2004) have been trialled to fit the experimental data, but every R² value obtained for those experimental data was smaller than 0.9. Hence, those two models were not been considered further in this study.

Table 3 provides details of six sets of data of plant treated waters and tap waters selected from six drinking water distribution systems (Anstey Hill, Barossa, Myponga, Happy Valley, Hope Valley and Little Para). It can be seen that both the BRP-Hill model (Equation (4)) and the Weibull model described well these growth curves. Half of the BRP-Hill model (Equation (4))'s determined maximum growth rates were closer to the observed growth rates, and half of the Weibull model's determined maximum growth rates were closer to observed growth rates (Table 3). For the lag time, all of the BRP-Hill model (Equation (4))'s determined lag times were closer to observed lag time than the Weibull model. There are four sets of data where the determined AOC values by the BRP-Hill model (Equation (4)) are closer to the observed AOC values than the Weibull model (Table 3).

There were 24 sets of data for comparing AICc values between the BRP-Hill model (Equation (4)) and the Weibull model selected randomly from all drinking water distribution system data for two types of water. Seventeen out of the 24 sets of data's BRP-Hill model (Equation (4)) had smaller AICc values, and the mean evidence ratio of BRP-Hill model and Weibull model for these 17 sets of data was 10.79 (data not shown). When some sets of data were fitted by Weibull model, their coefficient d of equation was smaller than 1. The curve is saturation shape when d ≤ 1, and the curve is sigmoid curve when d > 1 (Swintek et al. 2019). However, the microbiological growth curve generally is a sigmoidal curve (Fujikawa et al. 2004).

In summary, the BRP-Hill model (Equation (4)) was found to be better for describing bacterial regrowth curve in water than the Weibull model. Although both the Weibull model and the BRP-Hill model (Equation (4)) have good (high) R² values, the Weibull model cannot provide a good fit when experimental data plotted curve has a more saturation shape, and as well, the BRP-Hill model had better AICc scores.

A BRP-Hill model (Equation (4)) that provides microbial growth parameters to describe the bacterial growth curve has been developed with some alternative algorithms to adapt the large difference between the method built-in AOC value and the model-determined AOC value occurred for some experimental data. The coefficient of determination (R²) values obtained for the BRP-Hill model (Equation (4)) fitted data sets (n = 82) were close to 1. The average R² values were 0.987 for plant treated water and 0.983 for tap water (as shown in Table 2). Furthermore, the BRP method built-in AOC values and model-determined AOC values for the 82 sets of data were compared, and the result show this model fitted these data well. The BRP-Hill model (Equation (4)) parameters can help improving understanding of BRP and also of the AOC in waters. This suggests that the developed model can be a helpful tool to understand the behaviour of bacterial regrowth of drinking water and enable prediction where some data is not available. However, the BRP-Hill model (Equation (4)) was developed based on experiments carried out at a constant temperature (22 °C); thus, all the BRP experiments were conducted at this temperature. The BRP-Hill model (Equation (4)) may not be able to fit data from the BRP experiment conducted in a different environment such as at another constant temperature (30 °C).

The authors acknowledge the contribution of SA Water to support this study by providing the BRP instrument, result database and technical advice from their staff, Mary Drikas and Miriam Nedic.

All relevant data are included in the paper or its Supplementary Information.