DRomics (version 2.1-3)

bmdcalc: Computation of benchmark doses for responsive items

Description

Computes x-fold and z-SD benchmark doses for each responsive item using the best fit dose-reponse model.

Usage

bmdcalc(f, z = 1, x = 10)
# S3 method for bmdcalc
print(x, …)
# S3 method for bmdcalc
plot(x, BMDtype = c("zSD", "xfold"), 
            plottype = c("ecdf", "hist", "density"), 
            by = c("none", "trend", "model", "typology"),
            hist.bins = 30, …)

Arguments

f

An object of class "drcfit" returned by the function drcfit.

z

Value of z defining the BMD-zSD as the dose at which the response is reaching y0 +/- z * SD, with y0 the level at the control given by the dose-response fitted model and SD the residual standard deviation of the dose-response fitted model.

x

Value of x given as a percentage and defining the BMD-xfold as the dose at which the response is reaching y0 +/- (x/100) * y0, with y0 the level at the control given by the dose-response fitted model.

For print and plot functions, an object of class "bmdcalc".

BMDtype

The type of BMD to plot, "zSD" (default choice) or "xfold".

plottype

The type plot, "ecdf" for an empirical cumulative distribution plot (default choice), "hist" for a histogram or "density" for a density plot.

by

If different from "none" the plot is split by trend (if "trend"), by model (if "model") or by typology (if "typology").

hist.bins

The number of bins, only used for histogram(s).

…

further arguments passed to graphical or print functions.

Value

bmdcalc returns an object of class "bmdcalc", a list with 4 components:

res

a data frame reporting the results of the fit and BMD computation on each selected item sorted in the ascending order of the adjusted p-values returned by function itemselect. The different columns correspond to the identifier of each item (id), the row number of this item in the initial data set (irow), the adjusted p-value of the selection step (adjpvalue), the name of the best fit model (model), the number of fitted parameters (nbpar), the values of the parameters b, c, d, e and f, (NA for non used parameters), the residual standard deviation (SDres), the typology of the curve (typology, (twelve class typology described in the help of the drcfit function)), the rough trend of the curve (trend) defined with four classes (U, bell, increasing or decreasing shape), the theoretical value at the control (y0), the theoretical y range for x within the range of tested doses (yrange) and for biphasic curves the x value at which their extremum is reached (xextrem) and the corresponding y value (yextrem), the BMD-zSD value (BMD.zSD) and the BMD-xfold value (BMD.xfold).

z

Value of z given in input to define the BMD-zSD.

x

Value of x given in input as a percentage to define the BMD-xfold.

microarraydata

The corresponding object of class "microarraydata" given in input (component of itemselect).

Details

Two types of benchmark doses (BMD) were computed for each responsive item using the best fit dose-reponse model previously obtained using the drcfit function :

• the BMD-zSD defined as the dose at which the response is reaching y0 +/- z * SD, with y0 the level at the control given by the dose-response model, SD the residual standard deviation of the dose response model fit and z given as an input (z fixed to 1 by default),

• the BMD-xfold defined as the dose at which the response is reaching y0 +/- (x/100) * y0, with y0 the level at the control given by the dose-response fitted model and x the percentage given as an input (x fixed at 10 by default.)

When there is no analytical solution for the BMD, it is numerically searched along the fitted curve using the uniroot function.

In cases where the BMD cannot be reached due to the asymptote at high doses, NaN is returned. In cases where the BMD is not reached at the highest tested dose, NA is returned.

References

Larras F, Billoir E, Baillard V, Siberchicot A, Scholz S, Wubet T, Tarkka M, Schmitt-Jansen M and Delignette-Muller ML (2018). DRomics: a turnkey tool to support the use of the dose-response framework for omics data in ecological risk assessment. Environmental science & technology.https://doi.org/10.1021/acs.est.8b04752

See Also

See uniroot for details about the function used for the numerical search of the benchmark dose for cases where there is no analytical solution.

Examples

# NOT RUN {
# (1) a toy example (a very small subsample of a microarray data set) 
#
datafilename <- system.file("extdata", "transcripto_very_small_sample.txt", package="DRomics")

# to test the package on a small (for a quick calculation) but not very small data set
# use the following commented line
# datafilename <- system.file("extdata", "transcripto_sample.txt", package="DRomics")

(o <- microarraydata(datafilename, check = TRUE, norm.method = "cyclicloess"))
(s_quad <- itemselect(o, select.method = "quadratic", FDR = 0.01))
(f <- drcfit(s_quad, progressbar = TRUE))
(r <- bmdcalc(f))
plot(r) 

# changing the values of z and x for BMD calculation
# }
# NOT RUN {
(rb <- bmdcalc(f, z = 2, x = 50))
plot(rb)
# }
# NOT RUN {
# (2) an example on a microarray data set (a subsample of a greater data set) 
#
# }
# NOT RUN {
datafilename <- system.file("extdata", "transcripto_sample.txt", package="DRomics")

# to test the package on a small (for a quick calculation) but not very small data set
# use the following commented line
# datafilename <- system.file("extdata", "transcripto_sample.txt", package="DRomics")

(o <- microarraydata(datafilename, check = TRUE, norm.method = "cyclicloess"))
(s_quad <- itemselect(o, select.method = "quadratic", FDR = 0.01))
(f <- drcfit(s_quad, progressbar = TRUE))
(r <- bmdcalc(f))
plot(r) 
if (require(ggplot2))
  plot(r, plottype = "ecdf") + scale_x_log10() # with log10 dose scale 

# different plots of BMD-zSD

plot(r, plottype = "hist") 
plot(r, plottype = "density") 
plot(r, plottype = "density", by = "trend") 
plot(r, plottype = "ecdf", by = "trend") 
plot(r, plottype = "ecdf", by = "model") 
plot(r, plottype = "ecdf", by = "typology") 

# a plot of BMD-xfold (by default BMD-zSD is plotted)
plot(r, BMDtype = "xfold", plottype = "hist", by = "typology", hist.bins = 10) 
# }
# NOT RUN {
# }