# Using ssd4mosaicâ€™s functions in R

library(ssd4mosaic)

When using the MOSAIC SSD web application, a code is provided after each analysis to reproduce the same results directly in R. Here is an example of censored data species sensitivity distribution analysis using {ssd4mosaic} functions.

## Defining the inputs

# Data creation
# Most often, you would archive the same result by reading a table file with a
data <- ssd4mosaic::fluazinam

# Which distribution to fit to the data.
# See get_fits function documentation for possible options
distributions <- list("lnorm")
# Whether to display the results plots with a logscale x-axis
logscale <- TRUE
# Concentration unit for plots labels
unit <- "\u03bcg/L"

## Fitting to the data

## model fitting
fits <- ssd4mosaic::get_fits(data, distributions, TRUE)

## bootstrapping
bts <- ssd4mosaic::get_bootstrap(fits)[[1]]

## Extracting information from the fit

## Model parameters
lapply(fits, summary)
#> [[1]]
#> Fitting of the distribution ' lnorm ' By maximum likelihood on censored data
#> Parameters
#>         estimate Std. Error
#> meanlog 4.976920  0.7422075
#> sdlog   2.687785  0.6056713
#> Loglikelihood:  -72.81266   AIC:  149.6253   BIC:  150.9034
#> Correlation matrix:
#>           meanlog     sdlog
#> meanlog 1.0000000 0.1350239
#> sdlog   0.1350239 1.0000000

## HCx values
lapply(bts, quantile, probs = c(0.05, 0.1, 0.2, 0.5))
#> [[1]]
#> (original) estimated quantiles for each specified probability (censored data)
#>            p=0.05    p=0.1    p=0.2    p=0.5
#> estimate 1.743522 4.629205 15.10194 145.0271
#> Median of bootstrap estimates
#>            p=0.05    p=0.1    p=0.2    p=0.5
#> estimate 1.941508 5.110783 16.45573 148.3947
#>
#> two-sided 95 % CI of each quantile
#>            p=0.05     p=0.1     p=0.2     p=0.5
#> 2.5 %   0.3247236  1.084853  4.172938  36.89388
#> 97.5 % 18.6651384 35.903303 83.634981 787.67532
## CDF plot with confidence intervals
p <- ssd4mosaic::base_cdf(fits, unit = unit, logscale = logscale)
ssd4mosaic::options_plot(fits, unit, logscale, data, use_groups = TRUE)