Bayesian change point quantile regression approach to enhance the understanding of shifting phytoplankton-dimethyl sulfide relationships in aquatic ecosystems

Highlights

• We proposed a novel Bayesian change point quantile regression model (BCPQR).

• BCPQR illustrates shifting CHL-DMS relationships at multiple quantiles.

• The inequality of variance was revealed for the relationship in seas.

• The pH is not likely to affect CHL-DMS relationships in aquatic ecosystems.

Abstract

Dimethyl sulfide (DMS) serves as an anti-greenhouse gas, plays multiple roles in aquatic ecosystems, and contributes to the global sulfur cycle. The chlorophyll a (CHL, an indicator of phytoplankton biomass)-DMS relationship is critical for estimating DMS emissions from aquatic ecosystems. Importantly, recent research has identified that the CHL-DMS relationship has a breakpoint, where the relationship is positive below a CHL threshold and negative at higher CHL concentrations. Conventionally, mean regression methods are employed to characterize the CHL-DMS relationship. However, these approaches focus on the response of mean conditions and cannot illustrate responses of other parts of the DMS distribution, which could be important in order to obtain a complete view of the CHL-DMS relationship. In this study, for the first time, we proposed a novel Bayesian change point quantile regression (BCPQR) model that integrates and inherits advantages of Bayesian change point models and Bayesian quantile regression models. Our objective was to examine whether or not the BCPQR approach could enhance the understanding of shifting CHL-DMS relationships in aquatic ecosystems. We fitted BCPQR models at five regression quantiles for freshwater lakes and for seas. We found that BCPQR models could provide a relatively complete view on the CHL-DMS relationship. In particular, it quantified the upper boundary of the relationship, representing the limiting effect of CHL on DMS. Based on the results of paired parameter comparisons, we revealed the inequality of regression slopes in BCPQR models for seas, indicating that applying the mean regression method to develop the CHL-DMS relationship in seas might not be appropriate. We also confirmed relationship differences between lakes and seas at multiple regression quantiles. Further, by introducing the concept of DMS emission potential, we found that pH was not likely a key factor leading to the change of the CHL-DMS relationship in lakes. These findings cannot be revealed using piecewise linear regression. We thereby concluded that the BCPQR model does indeed enhance the understanding of shifting CHL-DMS relationships in aquatic ecosystems and is expected to benefit efforts aimed at estimating DMS emissions. Considering that shifting (threshold) relationships are not rare and that the BCPQR model can easily be adapted to different systems, the BCPQR approach is expected to have great potential for generalization in other environmental and ecological studies.

Introduction

Dimethyl sulfide (DMS) was recognized as an anti-greenhouse gas because its oxidized products acted as cloud condensation nuclei, which reflected solar irradiation and thereby contributed to the reduction of earth temperature (i.e., the CLAW hypothesis, Charlson et al., 1987). Despite recent debates or rejection of the CLAW hypothesis (Cropp, Gabric, van Tran, Jones, Swan, Butler, 2018, Quinn, Bates, 2011), the hypothesis is likely relevant in some regions like the Southern and Arctic Oceans (Krüger, Graßl, 2011, Levasseur, 2013). DMS also plays multiple essential roles in aquatic ecosystems, such as serving as an antioxidant for phytoplankton (Sunda et al., 2002) and facilitating a tritrophic mutualism between primary producers and top predators (Savoca and Nevitt, 2014). In addition, DMS is important to the global sulfur cycle (Eyice et al., 2015), accounting for about 80$\%$ of global biogenic sulfur emissions to the atmosphere (Kettle and Andreae, 2000). Phytoplankton, indicated by chlorophyll a (CHL) (Bates, Kiene, Wolfe, Matrai, Chavez, Buck, Blomquist, Cuhel, 1994, Zhang, Yang, Zhu, 2008), is the major producer of DMS (Charlson, Lovelock, Andreae, Warren, 1987, Gondwe, Krol, Gieskes, Klaassen, de Baar, 2003). Therefore, understanding the CHL-DMS relationship is critical for estimating regional or global DMS emissions from aquatic ecosystems (Anderson, Spall, Yool, Cipollini, Challenor, Fasham, 2001, Galí, Devred, Levasseur, Royer, Babin, 2015, Simó, Dachs, 2002).

Correlation analysis and ordinary linear regression have been the most widely used methods to explore the CHL-DMS relationship. Most studies revealed a positive effect of CHL on DMS (e.g., a significantly positive correlation coefficient or regression slope) (Gao, Yang, Zhang, Liu, 2017, Iverson, Nearhoof, Andreae, 1989, Lana, Simó, Vallina, Dachs, 2011, Law, Smith, Harvey, Bell, Cravigan, Elliott, Lawson, Lizotte, Marriner, McGregor, Ristovski, Safi, Saltzman, Vaattovaara, Walker, 2017, Lizotte, Levasseur, Galindo, Gourdal, Gosselin, Tremblay, Blais, Charette, Hussherr, 2020, Tan, Wu, Liu, Yang, 2017, Tortell, Guguen, Long, Payne, Lee, Ditullio, 2011, Walker, Harvey, Bury, Chang, 2000, Yang, 1999, Yang, 2000, Yang, Tsunogai, 2005, Yang, Zhang, Su, Zhou, 2009, Yang, Zhang, Zhou, Yang, 2011, Zhang, Yang, Zhang, Yang, 2014), while several studies reported a negative relationship (Froelichd et al., 1985) or no relationship at all (Nemcek, Ianson, Tortell, 2008, Watanabe, Yamamoto, Tsunogai, 1995). We note that maximum CHL concentrations in studies deducing positive CHL-DMS relationships were always much lower than those deducing negative or no relationships. For example, CHL concentrations in a series of studies on Chinese seas (Yang, 1999, Yang, 2000, Yang, Tsunogai, 2005, Yang, Zhang, Su, Zhou, 2009, Yang, Zhang, Zhou, Yang, 2011) were all lower than 4 $\mu$g/L. In contrast, the CHL concentration can reach approximate 60 $\mu$g/L in Froelichd et al. (1985).

A recent study examining the CHL-DMS relationship across a broad range of CHL concentrations implemented a change point model to capture the ascending and descending limbs of this relationship (Deng et al., 2020). The change point model aims to determine one or more unknown change points at which the stressor-response relationship changes. In Deng et al. (2020), the authors used 246 paired observations of CHL and DMS from 100 Chinese lakes and collected 426 paired observations from global oceans. They applied a piecewise linear regression model (Muggeo, 2003) to detect thresholds of CHL concentration, at which CHL-DMS relationships significantly changed. Benefiting from the novel application of piecewise regression, the authors revealed hump-shaped CHL-DMS relationships in both lakes and seas, which were expected to increase the estimation accuracy of global DMS emissions from aquatic ecosystems (Deng et al., 2020). The hump-shaped relationship also seemed to resolve the contradiction of the sign of the CHL-DMS relationship in previous studies, whose deductions might have been constrained by a relatively smaller sample size, a narrow range of sampled CHL concentration, or an application of a overly simplified linear regression model.

Although many informative studies have investigated the CHL-DMS relationship, we note that those studies mainly used mean regression methods (e.g., the ordinary linear regression or piecewise regression), by which the relationship between CHL and the mean of DMS distribution was estimated. A practically important alternative to classical mean regression methods is quantile regression (QR) (Koenker and Bassett, 1978). To the best of our knowledge, QR has not been used to examine CHL-DMS relationships.

QR explores the effect of one or more predictors on any quantile of the response variable distribution (Das, Krzywinski, Altman, 2019, Koenker, Bassett, 1978). Compared with mean regression methods, QR can provide a more complete view of possible causal relationships and can reveal useful predictive relationships at some parts of the response variable distribution, even when there is a weak or no predictive relationship between the predictor(s) and the mean of the response variable distribution (Cade and Noon, 2003). In addition, QR appears more robust to outliers of the response variable (Scharf et al., 1998) and is not constrained by the equal variance assumption (Cade, Noon, 2003, Das, Krzywinski, Altman, 2019). QR has been successfully applied to environmental and ecological studies. QR has been used to 1) illustrate a relatively complete view of stressor-response relationships at multiple regression quantiles (Cade, Terrell, Porath, 2008, Liang, Xu, Qiu, Liu, Lu, Wagner, 2021, Muller, Cade, Schwarzkopf, 2018, Niinemets, Valladares, 2006, Simkin, Allen, Bowman, Clark, Belnap, Brooks, Cade, Collins, Geiser, Gilliam, Jovan, Pardo, Schulz, Stevens, Suding, Throop, Waller, 2016, Xu, Schroth, Isles, Rizzo, 2015), 2) obtain reliable prediction intervals of the response variable (Heiskary, Bouchard, 2015, Kampichler, Sierdsema, 2018), and 3) reveal the limiting effect of the stressor on the response variable via the upper boundary of the stressor-response relationship (Fornaroli, Cabrini, Zaupa, Bettinetti, Ciampittiello, Boggero, 2016, Keeley, Macleod, Forrest, 2012, Youngflesh, Jenouvrier, Li, Ji, Ainley, Ballard, Barbraud, Delord, Dugger, Emmerson, Fraser, Hinke, Lyver, Olmastroni, Southwell, Trivelpiece, Trivelpiece, Lynch, 2017). The upper boundary of a stressor-response relationship illustrates the behavior of response variable when the stressor is the limiting factor (Cade, Terrell, Schroeder, 1999, Sankaran, Hanan, Scholes, Ratnam, Augustine, Cade, Gignoux, Higgins, Roux, Ludwig, Ardo, Banyikwa, Bronn, Bucini, Caylor, Coughenour, Diouf, Ekaya, Feral, February, Frost, Hiernaux, Hrabar, Metzger, Prins, Ringrose, Sea, Tews, Worden, Zambatis, 2005).

Because the CHL-DMS relationship represents a stressor-response relationship (McDowell et al., 2018), QR seems applicable and helpful to enhance the understanding of the CHL-DMS relationship. Considering the recent finding on the shifting nature of CHL-DMS relationships (Deng et al., 2020), a simple linear QR might not be adequate. A QR method with the ability to detect a change point is required but has rarely been explored (an exploration of this approach could be found in Zhou et al. (2015) who proposed a sequential change point detection method for linear QR).

In this study, we propose a novel Bayesian change point quantile regression (BCPQR) approach to investigate the CHL-DMS relationship in aquatic ecosystems. The BCPQR model integrates two well-developed Bayesian models: a Bayesian change point (BCP) model (Barry, Hartigan, 1993, Erdman, Emerson, 2007) and a Bayesian quantile regression (BQR) model (Benoit, Van den Poel, 2017, Yu, Moyeed, 2001). Both the BCP model (Beckage, Joseph, Belisle, Wolfson, Platt, 2007, Liang, Qian, Wu, Chen, Liu, Yu, Yi, 2019, Thomson, Kimmerer, Brown, Newman, Nally, Bennett, Feyrer, Fleishman, 2010) and the BQR model (Barneche, Kulbicki, Floeter, Friedlander, Allen, 2016, Uranchimeg, Kim, Kim, Kwon, Lee, 2018, Yu, Zou, Wang, 2019, Zou, Shi, 2020) have been recently introduced and applied to develop a stressor-response relationship in environmental and ecological fields. However, to our knowledge, this is the first proposal of BCPQR model in environmental and ecological studies.

There are several features of the BCPQR model that makes it desirable for ecological investigations. First, the BCPQR model inherits advantages of the BQR model and the BCP model. It is expected to be able to provide a complete view on the stressor-response relationship (Muller, Cade, Schwarzkopf, 2018, Xu, Schroth, Isles, Rizzo, 2015). The detection of any change point in the regression intercept, slopes, and/or variance of residuals is possible (Beckage, Joseph, Belisle, Wolfson, Platt, 2007, Liang, Qian, Wu, Chen, Liu, Yu, Yi, 2019). Second, the Bayesian framework would provide the convenience for parameters estimation. We can straightforwardly incorporate the change point into the BQR model structure. Parameter estimation of BCPQR model could then be achieved using Markov-chain Monte Carlo (MCMC) methods (Qian et al., 2003). Moreover, the parameter estimation framework would allow for the calculation of probability densities representing the uncertainty of parameters (including the change point and the other model parameters) (Ellison, 2004, Gende, Hendrix, Harris, Eichenlaub, Nielsen, Pyare, 2011, Underwood, Rizzo, Schroth, Dewoolkar, 2017). In addition, based on posterior distributions of parameters, comparing parameters is straightforward (Alameddine, Qian, Reckhow, 2011, Qian, Craig, Baustian, Rabalais, 2009). Finally, prior information – if available – could be used during model development (Ellison, 1996, Ellison, 2004).

Our objective was to examine whether or not the BCPQR approach can enhance the understanding of CHL-DMS relationships in lakes and seas. We applied the proposed BCPQR model to reevaluate shifting CHL-DMS relationships revealed in Deng et al. (2020). We separately fitted BCPQR models at five regression quantiles. To avoid over confidence in the change point model, a common practice is comparing the change point model with a model without any change point (Cahill, Rahmstorf, Parnell, 2015, Liang, Qian, Wu, Chen, Liu, Yu, Yi, 2019). Therefore, we also fitted a BQR model at each regression quantile and compared performances of the two models as a means to select the best model for characterizing the CHL-DMS relationship.