GUTS: Fast Calculation of the Likelihood of a Stochastic Survival...

GUTSR Documentation

Fast Calculation of the Likelihood of a Stochastic Survival Model

Description

GUTS (General Unified Threshold model of Survival) is a stochastic survival model for ecotoxicology. The package allows for the definition of exposure and survival time series as well as parameter values, and the fast calculation of the survival probabilities as well as the logarithm of the corresponding likelihood.

The package implements the GUTS-SIC (also called GUTS-RED) variants that assume a one-compartment model with first-order toxicokinetics.

Usage

guts_setup(C, Ct, y, yt, dist = "lognormal",
	model = "Proper",
	N = 1000L,
	MF = 100L,
	M = max(
		5000L,
		as.integer(ceiling(MF * length(union(Ct, yt)))),
		as.integer(ceiling(MF * max(union(Ct, yt))))
		),
	SVR = 1L,
	study = "", Clevel = ""
	)

guts_calc_loglikelihood(gobj, par, external_dist = NULL,
  use_multinomial_coefficient = FALSE)

guts_calc_survivalprobs(gobj, par, external_dist = NULL)

guts_report_damage(gobj)

guts_report_sppe(gobj)

guts_report_squares(gobj)

Arguments

C

Numeric vector of concentrations. Vector must contain at least 2 values and be of the same length as Ct.

Ct

Numeric vector of concentration time points. Vector must contain at least 2 values and be of the same length as C. Time points must start at 0, and contain unique values in ascending order.

y

Integer vector (counts) of survivors. Vector must contain at least 2 values and be of the same length as yt. y must not be ascending.

yt

Numeric vector of survivor time points. Vector must contain at least 2 values and be of the same length as y. Time points must start at 0, and contain unique values in ascending order. Survivor information at time points later than the latest concentration time point will be disregarded (with a warning).

dist

Distribution as character, either “lognormal” (default), “loglogistic”, “external” or “delta”.

model

Model as character, either “Proper” (for full model, the default), “IT” (for individual tolerance), or “SD” (for stochastic death).

MF

Integer. Multiplication factor for M. Must be greater than 1. MF is used only if “model = 'SD'” or “model = 'Proper'” and M is not specified. Setting MF automatically ensures that the number of points for time discretization M is at least the number of measurement time steps or the measurement time (which ever is larger) multiplied by MF. A minimum of M = 5000 is ensured.

M

Integer. Desired number of points for time discretization. Must be greater than 1. M is used only if “model = 'SD'” or “model = 'Proper'”

N

Integer. Sample length of individual tolerance thresholds. Must be greater than 2. N is used only, if “model = 'Proper'”

study

string with the name of the study

Clevel

character vector with names for each of the concentraton levels

SVR

Numeric surface-volume-ratio. A multiplication factor to kd.

gobj

GUTS object. The object to be updated (and used for the calculation).

par

Numeric vector of parameters. See details below.

external_dist

Numeric vector containing the distribution of individual thresholds. Only used if dist = 'external'. See details below.

use_multinomial_coefficient

If “TRUE” returns loglikelihood from the correct multinomial distribution. Defaults to ignoring the constant multinomial coefficient for performance reasons.

Details

Functions

Use guts_setup to define (or alter) a GUTS object. Various checks are applied to the data. On success, a GUTS object will be created.

Use guts_calc_loglikelihood to calculate the survival probabilities and the corresponding loglikelihood for a given set of parameters. The function is very fast and can be used in routines for parameter estimation. The function returns the loglikelihood, however it also updates the fields par, S, D, SPPE, squares, zt and LL of the GUTS-object.

guts_calc_survivalprobs is a convenience wrapper that can be used for predictions; it returns the survival probabilities, however it also updates the fields par, S, D, SPPE, squares, zt and LL of the GUTS-object.

guts_report_damage returns a data.frame with time grid points and the damage for each of these. The function reports the damage that was calculated in the previous call to guts_calc_loglikelihood or guts_calc_survivalprobs.

guts_report_squares returns the sum of squares. The function reports the sum of squares that was calculated in the previous call to guts_calc_loglikelihood or guts_calc_survivalprobs.

guts_report_sppe returns the survival-probability prediction error (SPPE). The function reports the SPPE that was calculated in the previous call to guts_calc_loglikelihood or guts_calc_survivalprobs.

Models, Parameters, and Distributions

The GUTS package provides three model types:

  • Proper: a GUTS-SIC-Proper (also called GUTS-RED-Proper) model using random individual tolerances and a stochastic death process, when individual tolerances are exceeded.

  • IT: a GUTS-SIC-IT (GUTS-RED-IT) individual tolerance model using random individual tolerances. If an individual's tolerance threshold is exceeded, the individual dies.

  • SD: a GUTS-SIC-SD (GUTS-RED-SD) stochastic death model using a stochastic death process above a population-wide tolerance threshold. The tolerance-threshold is the same for all individuals.

The Proper GUTS model requires the following parameters par, while variants IT and SD are based on a reduced subset (as indicated in brackets). Parameter values in par must be ordered as listed here:

  • hb: background mortality rate (Proper, IT, SD)

  • ke: dominant rate constant (Proper, IT, SD)

  • kk: killing rate (Proper, SD)

  • further parameters for the tolerance threshold (in SD) or the threshold distribution dist (in Proper and IT)

For model type “SD” (stochastic death), required parameters par[1:4] are hb, ke, kk and mn, which is the population-wide tolerance threshold. For backwards compatibility this model type can be initiated setting dist = "Delta" and model = "Proper".

For model type “IT” (individual tolerance), required parameters par[1:2] are hb, ke, as well as respective distribution parameters (from par[3] onwards). Parameter (kk) is set internally to infinity and does not need to be provided.

For model type “Proper”, all parameters are needed. par[1:3] take hb, ke, kk, distribution parameters follow (from par[4] onwards).

For model types “Proper” and “IT” individual tolerance thresholds are created internally. Individual tolerances are drawn from the specified distribution dist. The parameter values required depend on the specified dist:

  • "lognormal": requires the parameters mn and sd which are the mean and standard deviation of the lognormal random distribution. In contrast to parameters meanlog and sdlog of function dlnorm, these parameters are not on the logscale. They relate in the following way:

    sdlog = \sqrt{\ln{\frac{1 + sd^2}{mn^2}}}

    meanlog = \ln{mn} - \frac{1}{2} * sdlog^2

  • "loglogistic": requires the parameters mn = scale = median and beta = shape.

  • "external": uses random variates provided to external_dist. With this option GUTS can be run with arbitrarily distributed individual tolerance thresholds. With the option “external” only parameters hb, ke and kk (Proper only) are required. Further, the thresholds sample length N is internally adjusted to the length of the external vector of random variates external_dist. The adjustment of N is notified by a warning.

For performance reasons the implemented distributions “lognormal” and “loglogistic” are approximated using importance sampling. The option “external” generally performs well, but might require a larger thresholds sample (i.e. length(external_dist) should be large).

The number of parameters is checked according to dist and model. Wrong number of parameters invokes an error, wrong parameter values (e.g., negative values) invoke a warning, and the loglikelihood is set to -Inf.

Field and Attribute Access

Fields and attributes of an object of class “GUTS” are read-only. It is not possible to directly modify single elements of the GUTS object. Instead use function guts_setup to create GUTS objects or modify fields on existing GUTS objects. Functions guts_calc_loglikelihood and guts_calc_survivalprobs update an object's fields par (parameters), D (damage), squares (sum of squares), SPPE (survival-probability prediction error), S (survival probabilities) and LL (the loglikelihood).

Value

guts_setup returns a list of class “GUTS” with the following fields:

C

Concentrations.

Ct

Concentration time points.

y

Survivors.

yt

Survivor time points.

dist

Distribution.

model

Model.

N

Sample length.

M

Time grid points.

par

Parameters.

S

Vector of survivor probabilities.

D

Vector of internal damage for each of the M time grid points.

squares

Sum of squares

SPPE

Survival-probability prediction error.

LL

The loglikelihood.

guts_calc_loglikelihood returns the loglikelihood.

guts_calc_survivalprobs returns the survival probabilities.

guts_report_damage returns the damage.

guts_report_squares returns the sum of squares.

guts_report_sppe returns the survival-probability prediction error (SPPE).

Note

The GUTS project web site can be found here: http://guts.r-forge.r-project.org. For questions and discussion, please subscribe to the mailing list there.

Author(s)

Carlo Albert carlo.albert@eawag.ch, Sören Vogel soeren.vogel@posteo.ch, Oliver Jakoby oliver.jakoby@rifcon.de, Alexander Singer alexander.singer@rifcon.de and Dirk Nickisch dirk.nickisch@rifcon.de

Maintainer: Oliver Jakoby oliver.jakoby@rifcon.de

References

Albert, C., Vogel, S., and Ashauer, R. (2016). Computationally efficient implementation of a novel algorithm for the General Unified Threshold Model of Survival (GUTS). PLOS Computational Biology, 12(6), e1004978. \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1371/journal.pcbi.1004978")}.

Jager, T., Albert, C., Preuss, T., and Ashauer, R. (2011). General Unified Threshold Model of Survival – a toxicokinetic toxicodynamic framework for ecotoxicology. Environmental Science & Technology, 45(7), 2529–2540, \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1021/es103092a")}.

Ashauer, R., Albert, C., Augustine, S., Cedergreen, N., Charles, S., Ducrot, V., Focks, A., Gabsi, F., Gergs, A., Goussen, B., Jager, T., Kramer, N.I., Nyman, A.-M., Poulsen, V., Reichenberger, S., Schäfer, R.B., Van den Brink, P.J., Veltman, K., Vogel, S., Zimmer, E.I., Preuss, T.G. (2016) Modelling survival: exposure pattern, species sensitivity and uncertainty. Scientific Reports, 6, 1, \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1038/srep29178")}.

Jager, T., Ashauer, R. (2018). Modelling survival under chemical stress. A comprehensive guide to the GUTS framework. Leanpub: https://leanpub.com/guts_book, http://www.debtox.info/book_guts.html

EFSA PPR Panel (EFSA Panel on Plant Protection Products and their Residues), Ockleford, C., Adriaanse, P., Berny, P., Brock, T., Duquesne, S., Grilli, S., Hernandez-Jerez, A.F., Bennekou, S.H., Klein, M., Kuhl, T., Laskowski, R., Machera, K., Pelkonen, O., Pieper, S., Smith, R.H., Stemmer, M., Sundh, I., Tiktak, A., Topping, C.J., Wolterink, G., Cedergreen, N., Charles, S., Focks, A., Reed, M., Arena, M., Ippolito, A., Byers, H. and Teodorovic, I. (2018). Scientific Opinion on the state of the art of Toxicokinetic/Toxicodynamic (TKTD) effect models for regulatory risk assessment of pesticides for aquatic organisms. EFSA Journal, 16(8):5377, 188 pp., \Sexpr[results=rd]{tools:::Rd_expr_doi("10.2903/j.efsa.2018.5377")}.

Gergs, A., Gabsi, F., Zenker, A., Preuss, T. (2011). Demographic toxicokinetic–toxicodynamic modeling of lethal effects. Environmental Science & Technology, 50(11), 6017–6024, \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1021/acs.est.6b01113")}.

Nickisch, D., Rall, B., Singer, A., Ashauer, R. (2022). Fish Species Sensitivity Ranking Depends on Pesticide Exposure Profiles. Environmental Toxicology & Chemistry, 41, 1732-1741, \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1002/etc.5348")}.

Singer, A., Nickisch, D., Gergs, A. (2023). Joint survival modelling for multiple species exposed to toxicants. Science of The Total Environment, 857, 159266, \Sexpr[results=rd]{tools:::Rd_expr_doi("10.1016/j.scitotenv.2022.159266")}.

See Also

diazinon, GUTS-package and the package vignettes for examples on how to calibrate and project GUTS-models.

Examples

data(diazinon)

# create GUTS object to calculate the Proper model
#  using a log-normal distribution of tolerance thresholds
gts.lognormal <- guts_setup(
  C = diazinon$C1, Ct = diazinon$Ct1,
  y = diazinon$y1, yt = diazinon$yt1,
  dist = "lognormal", model = "Proper")

# calculate likelihood of Proper model using log-normal distribution
guts_calc_loglikelihood(
  gts.lognormal,
  c(0.051, 0.126, 1.618, 19.099, 6.495))
gts.lognormal  # show GUTS object

# repeating calculation above
#  with threshold values from an external log-normal distribution.
#  Note, we need to account for the different parametrisations
#  used in the GUTS-package and in rlnorm
sigma2 <- log( 1 + 6.495^2 / 19.099^2)
mu <- log(19.099) - 0.5 * sigma2
lognormal.thresholds <- rlnorm(1000, meanlog = mu, sdlog = sqrt(sigma2))
gts.external <- guts_setup(
  C = diazinon$C1, Ct = diazinon$Ct1,
  y = diazinon$y1, yt = diazinon$yt1,
  dist = "external", model = "Proper")
guts_calc_loglikelihood(
  gts.external,
  c(0.051, 0.126, 1.618), external_dist = lognormal.thresholds)
# -> Results using external and internal distributions are comparable

# create GUTS object to calculate the Proper model
#  using a log-logistic distribution of tolerance thresholds
gts.loglogistic <- guts_setup(
  C = diazinon$C1, Ct = diazinon$Ct1,
  y = diazinon$y1, yt = diazinon$yt1,
  dist = "loglogistic", model = "Proper")
guts_calc_survivalprobs( # returning survival probabilities
  gts.loglogistic,
  c(0.01, 0.2, 0.3, 3, 2))

str(guts_report_damage(gts.loglogistic)) # returning damage

# calculate survival probabilities with IT model
#  using a log-logistic distribution of tolerance thresholds
guts_calc_survivalprobs(
  guts_setup(
  C = diazinon$C1, Ct = diazinon$Ct1,
  y = diazinon$y1, yt = diazinon$yt1,
  dist = "loglogistic", model = "IT"),
  c(0.01, 0.2, 3, 2))

# calculate survival probabilities with an SD model with a fixed tolerance threshold
guts_calc_survivalprobs(
  guts_setup(
  C = diazinon$C1, Ct = diazinon$Ct1,
  y = diazinon$y1, yt = diazinon$yt1,
  dist = "loglogistic", model = "SD"),
  c(0.01, 0.2, 0.3, 3))


## Not run: guts_calc_survivalprobs(gts.external, rep(.5, 3))
# Warning and no result, because no external distribution was specified

## Not run: guts_calc_survivalprobs(gts.loglogistic , 1:4 ) # Error.

## Not run: gts.loglogistic[["C"]] <- 1:3 # Error.

GUTS documentation built on Oct. 20, 2023, 9:08 a.m.