Repository files navigation

The package xeredar is an R-package for analysis of the New ApproachMethodology (NAM) assays of XETA (Xenopus Eleutheroembryonic Thyroid),RADAR (Rapid Androgen Disruption Activity Reporter) and REACTIV (RapidEstrogen ACTivity In Vivo) for assessing endocrine effects of chemicalson the thyroid, androgen/steroid and estrogen axis. The functionality isbased on the SAS-script recommended in the Annex 13 of OECD testguideline No. 248 of the XETA assay(2019),written by John Green.

You can install xeredar using one of the followig commands:

devtools::install_github("basf/xeredar")pak::pkg_install("basf/xeredar")

Data frames that are supposed to be analyzed with xeredar needs tofulfill certain requirements. The data frame or tibble needs to containthe following column headers:

knitr::kable(head(xeredar::testDataSpiked))

The type of each column should be accordingly:

knitr::kable(purrr::map_df(xeredar::testDataSpiked, class))

Replicate (i.e. run),Treatment (i.e. a unique name for eachtreatment level in either spiked or unspiked mode) andRow(i.e. exposure vessel) can either befactor orcharacter columns,butFluor (i.e. measured fluorescence) must always benumeric andConc (i.e. concentration of test item) must always containorderedfactors. The order of the columns is not relevant.

Replicate,Treatment andRow can either befactor orcharacter columns, butFluor must always benumeric andConcmust always containordered factors. The order of the columns is notrelevant. It is important that the decimal separator is a period insteadof a comma.

When simply aiming to use thedata_prep() function, the data frameneeds to either contain spiked treatments or unspiked treatments. Whenstill having spiked and unspiked treatments in one data frame theyshould be separated. Imagine you have a XETA data frame (dat) whichcontains spiked and unspiked treatments as well as the T4 positivecontrol. The T3 or T4 additions are designated by “T3” and “T4” in theTreatment column. The spiked and unspiked datasets could quickly besubset using the following code:

datSpiked <- dat[grepl("T3",dat$Treatment),] # spiked data datUnspiked <- setdiff(dat,datSpiked)datUnspiked <- datUnspiked[which(datUnspiked$Treatment != "T4"),] # unspiked data

To demonstrate how to run XETA analysis, we will use one of the datasets by the French lab containing the spiked and unspiked measurementsfrom the XETA ring test as included in OECD test guideline No. 248 ofthe XETA assay(2019).

xeta_spiked <- xeredar::valid_data_xeta[["ptu_france_spiked"]]xeta_unspiked <- xeredar::valid_data_xeta[["ptu_france_unspiked"]]

The default XETA analysis can be run using thedata_prep() functionwith either spiked or unspiked data. This function automatically decideswhether trimming, outlier removal and/or transformations are conductedfollowing the manuscript (Spyridonov et al. 2025).The actual analysis is carried out by theana() function. Theana() function iscalled by thedata_prep() function and does not need to be calledseparately. For this dataset, the exposure well ID (Row of the 96 wellplate) is not recorded, therefore, we set therow argument toFALSE.In this case, we use the reduced mixed ANOVA model where the exposurewell ID is not included as a random effect. Please specifyrow=TRUE ifthe exposure well ID is recorded and you want to use the full mixedANOVA model.

xeta_spiked_result <- xeredar::data_prep(dataframe = xeta_spiked, row = FALSE)xeta_unspiked_result <- xeredar::data_prep(dataframe = xeta_unspiked, row= FALSE)

Here we use the spike data as an example to demonstrate the output ofthedata_prep() function.

The outputs of thedata_prep() function are lists containing thefollowing elements:

• A reasoning for the recommended transformation and trimming. The rawdata should be used for the analysis because the residuals of themixed ANOVA are normally distributed and show homogeneous variancesamong treatment groups.

xeta_spiked_result$Justify#> [1] "The raw data should be used for the analysis because\n                  the residuals of the mixed ANOVA are normally\n                  distributed and show homogeneous variances\n                  among treatment groups."

• A data frame of the processed data (e.g. raw data, trimmed,transformed, or outlier removed) used for actual statistical testingfollowing the reasoning. The box plots of the processed data perrun/replicate (i.e. each panel represents each run/replicate) arealso provided for visual inspection.

knitr::kable(head(xeta_spiked_result$ProcessedData))

xeta_spiked_result$BoxPlots

• Summary tables of the processed data (per replicate and overall)

knitr::kable(xeta_spiked_result$SummaryDF_Rep, caption="Summary statistics of fluorescence in different concentrations of test item per replicate")

knitr::kable(xeta_spiked_result$SummaryDF, caption="Summary statistics of fluorescence in different concentrations of test item of all replicates")

• Tables of results evaluated using increasing/decreasing Williamstest and/or Dunnett’s test, if applicable.

In the the Williams’ test result tables,Y.Tilde is the amalgamatedmean of the fluorescence in each treatment group,Y0 is the mean ofthe control fluorescence,DIFF is the estimated difference between thetreatment and the control,SE_DIFF is the standard error of theWilliams’ test,DF is the degrees of freedom for Williams’ test,WILLIncr orWill Decr is the Williams’ test statistic,crit Val is thecritical value of Williams distribution,Sign suggests if there issignificant difference between the treatment and the control, and%Incr is the percent increase of the fluorescence compared to thecontrol.

In the Dunnett’s test result table,Estimate is the estimateddifference between the treatment and the control,SE is the standarderror of the mixed ANOVA model,t value is the Dunnett’s teststatistic,adj p is the adjusted p value and%Incr is the percentincrease of the fluorescence compared to the control.

knitr::kable(xeta_spiked_result$WilliamsIncrease, caption="Increasing Williams' test")

knitr::kable(xeta_spiked_result$WilliamsDecrease, caption="Decreasing Williams' test")

knitr::kable(xeta_spiked_result$Dunnetts, caption="Dunnett's test")

• Further information about the normality test (Shapiro-Wilk), thehomogeneity of variance test (Levene’s test) of residuals of themixed ANOVA model, the monotonicity test and the model fit.

xeta_spiked_result$NormalityTest#> #>  Shapiro-Wilk normality test#> #> data:  stats::resid(mixedaov)#> W = 0.99186, p-value = 0.05343xeta_spiked_result$LeveneTest#> # A tibble: 1 × 4#>   statistic p.value    df df.residual#>       <dbl>   <dbl> <int>       <int>#> 1     0.759   0.580     5         342xeta_spiked_result$`Monotonicity Test`#>        Test t value Pr(>|t|) Significance#> 1    Linear    6.49  <0.0001          ***#> 2 Quadratic    2.14   0.0335            *xeta_spiked_result$MixedAnova#> Linear mixed model fit by REML ['lmerMod']#> Formula: #> Fluor ~ Conc + (1 | Replicate) + (1 | Replicate:Conc)#>    Data: dataframe#> REML criterion at convergence: 5608.294#> Random effects:#>  Groups         Name        Std.Dev.#>  Replicate:Conc (Intercept) 241.1   #>  Replicate      (Intercept) 350.5   #>  Residual                   827.7   #> Number of obs: 348, groups:  #> Replicate:Conc, 18; Replicate, 3#> Fixed Effects:#> (Intercept)       Conc.L       Conc.Q       Conc.C  #>      4824.2        758.5        264.5         52.6  #>      Conc^4       Conc^5  #>       149.8       -270.8

The list output from running the data_prep() function can besummarized with the data_summary() function.

xeredar::data_summary(xeta_spiked_result) |> knitr::kable()

To demonstrate how to run RADAR analysis, we will use one of the datasets by the Pos_mDHT_Fraunhofer_RADAR containing the spiked andunspiked measurements from the RADAR study validation in the labFraunhofer with an androgen axis active chemical.

radar_spiked <- xeredar::RADAR_valid_data_table_spiked_unspiked[["mDHTFRAUNH_Spiked"]]radar_unspiked <-  xeredar::RADAR_valid_data_table_spiked_unspiked[["mDHTFRAUNH_Unspiked"]]

The default radar analysis can be run using thedata_prep() functionwith either spiked or unspiked data. This function automatically decideswhether trimming, outlier removal and/or transformations are conductedfollowing the manuscript (Spyridonov et al. unpublished). The actualanalysis is carried out by theana() function. Theana() function iscalled by thedata_prep() function and does not need to be calledseparately. the analysis the RADAR assay follows the description of theMethod 2 (the mixed ANOVA approach) in the Annex 8: methods for thestatistical analysis of RADAR assay data of the OECD TG 251(2022).For this dataset, the exposure well ID (Row of the 96 well plate) is notrecorded, therefore, we set therow argument toFALSE. In this case,we use the reduced mixed ANOVA model where the exposure well ID is notincluded as a random effect. Please specifyrow=TRUE if the exposurewell ID is recorded and you want to use the full mixed ANOVA model.Trimming is not required.

radar_spiked_result <- xeredar::data_prep(dataframe = radar_spiked, row = FALSE, trimming=FALSE)radar_unspiked_result <- xeredar::data_prep(dataframe = radar_unspiked, row= FALSE, trimming=FALSE)

Here we use the spike data as an example to demonstrate the output ofthedata_prep() function.

The outputs of thedata_prep() function are lists containing thefollowing elements:

• A reasoning for the recommended transformation and trimming. The rawdata (without trimming or outlier removal) where the fluorescencevalues are log transformed should be used for the analysis becauseonly after log transformation, the residuals of the mixed ANOVA arenormally distributed and have homogeneous variances among treatmentgroups.

radar_spiked_result$Justify#> [1] "The raw data (without trimming or outlier removal) where\n                    the fluorescence values are log transformed should be used\n                    for the analysis because only after log transformation,\n                    the residuals of the mixed ANOVA are normally distributed\n                    and have homogeneous variances among treatment groups."

• A data frame of the processed data (e.g. raw data, transformed, oroutlier removed) used for actual statistical testing following thereasoning. The box plots of the processed data per run/replicate(i.e. each panel represents each run/replicate) are also providedfor visual inspection.

knitr::kable(head(radar_spiked_result$ProcessedData))

radar_spiked_result$BoxPlots

• Summary tables of the processed data (per replicate and overall)

knitr::kable(radar_spiked_result$SummaryDF_Rep, caption="Summary statistics of fluorescence in different concentrations of test item per replicate")

knitr::kable(radar_spiked_result$SummaryDF, caption="Summary statistics of fluorescence in different concentrations of test item of all replicates")

• Tables of results evaluated using increasing/decreasing Williamstest and/or Dunnett’s test, if applicable.

In the the Williams’ test result tables,Y.Tilde is the amalgamatedmean of the fluorescence in each treatment group,Y0 is the mean ofthe control fluorescence,DIFF is the estimated difference between thetreatment and the control,SE_DIFF is the standard error of theWilliams’ test,DF is the degrees of freedom for Williams’ test,WILLIncr orWill Decr is the Williams’ test statistic,crit Val is thecritical value of Williams distribution,Sign suggests if there issignificant difference between the treatment and the control, and%Incr is the percent increase of the fluorescence compared to thecontrol.

In the Dunnett’s test result table,Estimate is the estimateddifference between the treatment and the control,SE is the standarderror of the mixed ANOVA model,t value is the Dunnett’s teststatistic,adj p is the adjusted p value and%Incr is the percentincrease of the fluorescence compared to the control.

knitr::kable(radar_spiked_result$WilliamsIncrease, caption="Increasing Williams' test")

knitr::kable(radar_spiked_result$WilliamsDecrease, caption="Decreasing Williams' test")

knitr::kable(radar_spiked_result$Dunnetts, caption="Dunnett's test")

• Further information about the normality test (Shapiro-Wilk), thehomogeneity of variance test (Levene’s test) of residuals of themixed ANOVA model, the monotonicity test and the model fit.

radar_spiked_result$NormalityTest#> #>  Shapiro-Wilk normality test#> #> data:  stats::resid(mixedaov)#> W = 0.99071, p-value = 0.05007radar_spiked_result$LeveneTest#> # A tibble: 1 × 4#>   statistic p.value    df df.residual#>       <dbl>   <dbl> <int>       <int>#> 1     0.474   0.796     5         300radar_spiked_result$`Monotonicity Test`#>        Test t value Pr(>|t|) Significance#> 1    Linear   -1.57   0.1167            .#> 2 Quadratic   -1.48   0.1400            .radar_spiked_result$MixedAnova#> Linear mixed model fit by REML ['lmerMod']#> Formula: #> log(Fluor) ~ Conc + (1 | Replicate) + (1 | Replicate:Conc)#>    Data: dataframe#> REML criterion at convergence: 1218.932#> Random effects:#>  Groups         Name        Std.Dev.#>  Replicate:Conc (Intercept) 0.000   #>  Replicate      (Intercept) 0.466   #>  Residual                   1.756   #> Number of obs: 306, groups:  #> Replicate:Conc, 18; Replicate, 3#> Fixed Effects:#> (Intercept)       Conc.L       Conc.Q       Conc.C  #>     8.35799     -0.37806     -0.37347      0.01863  #>      Conc^4       Conc^5  #>     0.68924     -0.25055  #> optimizer (nloptwrap) convergence code: 0 (OK) ; 0 optimizer warnings; 1 lme4 warnings

The list output from running the data_prep() function can besummarized with the data_summary() function.

xeredar::data_summary(radar_spiked_result) |>  knitr::kable()

To demonstrate how to run REACTIV analysis, we will use one artificialdata set containing the spiked and unspiked measurements.

reactiv_spiked <- xeredar::REACTIV_valid_data_table_spiked_unspiked[["Anastrozole_UK"]]$Spiked reactiv_unspiked <- xeredar::REACTIV_valid_data_table_spiked_unspiked[["Anastrozole_UK"]]$Unspiked

The default REACTIV analysis can be run using thedata_prep() functionwith either spiked or unspiked data. This function automatically decideswhether trimming, outlier removal and/or transformations are conductedfollowing the manuscript (Spyridonov et al. unpublished). The actualanalysis is carried out by theana() function. Theana() function iscalled by thedata_prep() function and does not need to be calledseparately. The analysis the REACTIV assay follows the description ofthe Method 2 (the mixed ANOVA approach) in the Annex 8: methods for thestatistical analysis of REACTIV assay data of the Amended Draft new TestGuideline for the REACTIV assay for second WNT-review(30.01.2024).For this assay, therow argument should be set toFALSE. Trimming isnot required. In case there are residuals deviate from normality andvariance homogeneity, outlier removal (e.g. by applying the Tukey rule(Green et al., 2018) and data transformation (for example log- orsquare-root) can be conducted.

reactiv_spiked_result <- xeredar::data_prep(dataframe = reactiv_spiked, row = FALSE, trimming=FALSE, boxcox = FALSE)reactiv_unspiked_result <- xeredar::data_prep(dataframe = reactiv_unspiked, row= FALSE, trimming=FALSE, boxcox = FALSE)

Here we use the spiked data as an example to demonstrate the output ofthedata_prep() function.

The outputs of thedata_prep() function are lists containing thefollowing elements:

• A reasoning for the recommended transformation and trimming. Thedata from which outliers were removed with the Tukey-rule where thefluorescence values are square-root transformed, should be used forthe analysis because only after outlier removal and sqrttransformation, the residuals of the mixed ANOVA are normallydistributed and have homogeneous variances among treatment groups

reactiv_spiked_result$Justify#> [1] "The data from which outliers were removed with the\n                    Tukey-rule where the fluorescence values are square-root\n                    transformed, should be used for the analysis because only\n                    after outlier removal and sqrt transformation, the\n                    residuals of the mixed ANOVA are normally distributed\n                    and have homogeneous variances among treatment groups"

• A data frame of the processed data (e.g. raw data, transformed, oroutlier removed) used for actual statistical testing following thereasoning. The box plots of the processed data per run/replicate(i.e. each panel represents each run/replicate) are also providedfor visual inspection.

knitr::kable(head(reactiv_spiked_result$ProcessedData))

reactiv_spiked_result$BoxPlots

• Summary tables of the processed data (per replicate and overall)

knitr::kable(reactiv_spiked_result$SummaryDF_Rep, caption="Summary statistics of fluorescence in different concentrations of test item per replicate")

knitr::kable(reactiv_spiked_result$SummaryDF, caption="Summary statistics of fluorescence in different concentrations of test item of all replicates")

• Tables of results evaluated using increasing/decreasing Williamstest and/or Dunnett’s test, if applicable.

In the Williams’ test result tables,Y.Tilde is the amalgamated meanof the fluorescence in each treatment group,Y0 is the mean of thecontrol fluorescence,DIFF is the estimated difference between thetreatment and the control,SE_DIFF is the standard error of theWilliams’ test,DF is the degrees of freedom for Williams’ test,WILLIncr orWill Decr is the Williams’ test statistic,crit Val is thecritical value of Williams distribution,Sign suggests if there issignificant difference between the treatment and the control, and%Incr is the percent increase of the fluorescence compared to thecontrol.

In the Dunnett’s test result table,Estimate is the estimateddifference between the treatment and the control,SE is the standarderror of the mixed ANOVA model,t value is the Dunnett’s teststatistic,adj p is the adjusted p value and%Incr is the percentincrease of the fluorescence compared to the control.

knitr::kable(reactiv_spiked_result$WilliamsIncrease, caption="Increasing Williams' test")

knitr::kable(reactiv_spiked_result$WilliamsDecrease, caption="Decreasing Williams' test")

knitr::kable(reactiv_spiked_result$Dunnetts, caption="Dunnett's test")

• Further information about the normality test (Shapiro-Wilk), thehomogeneity of variance test (Levene’s test) of residuals of themixed ANOVA model, the monotonicity test and the model fit.

reactiv_spiked_result$NormalityTest#> #>  Shapiro-Wilk normality test#> #> data:  stats::resid(mixedaov)#> W = 0.97621, p-value = 0.0227reactiv_spiked_result$LeveneTest#> # A tibble: 1 × 4#>   statistic p.value    df df.residual#>       <dbl>   <dbl> <int>       <int>#> 1     0.829   0.532     5         123reactiv_spiked_result$`Monotonicity Test`#>        Test t value Pr(>|t|) Significance#> 1    Linear   -6.34  <0.0001          ***#> 2 Quadratic    6.05  <0.0001          ***reactiv_spiked_result$MixedAnova#> Linear mixed model fit by REML ['lmerMod']#> Formula: #> sqrt(Fluor) ~ Conc + (1 | Replicate) + (1 | Replicate:Conc)#>    Data: wt_outlier#> REML criterion at convergence: 1607.975#> Random effects:#>  Groups         Name        Std.Dev.#>  Replicate:Conc (Intercept)   0.00  #>  Replicate      (Intercept)  51.74  #>  Residual                   151.56  #> Number of obs: 129, groups:  #> Replicate:Conc, 18; Replicate, 3#> Fixed Effects:#> (Intercept)       Conc.L       Conc.Q       Conc.C  #>      551.13      -471.74       437.64      -271.05  #>      Conc^4       Conc^5  #>      175.76       -30.96  #> optimizer (nloptwrap) convergence code: 0 (OK) ; 0 optimizer warnings; 1 lme4 warnings

The list output from running the data_prep() function can besummarized with the data_summary() function.

xeredar::data_summary(reactiv_spiked_result) |>  knitr::kable()

xeredar contains an integrated shiny app that is available to the usersby writingrun_app() into the console. The only requirement to use theapp is the successful installation of xeredar. When the app started, theuser has a couple of options to analyse the data. These options can beadjusted by clicking on the little gear sign next toInputs. Thefollowing options are available:

• Use d’Agostino test?
  • The default here is that this is not selected meaning theShapiro-Wilk test is utilized to check for residual normality.However, the RADAR TG mentions the d’Agostino test which is whyit is also available to the users.
• Apply 10% Trimming
  • The default here is that this is selected, meaning that 10%Trimming is conducted when the ANOVA assumptions are notfulfilled by the raw data. Please be aware that this does notmean that 10% Trimming is always conducted. It is onlyconducted, when the raw data is violating the residual normalityand variance homogeneity assumptions, tested with the respectivetest with an adjustable alpha level.
• Remove outliers
  • The default here is that this is selected, meaning that outlierremoval is conducted when the ANOVA assumptions are notfulfilled by the raw data or by 10% Trimming, if it is selectedabove.
• Try Box-Cox transformation
  • The default here is that this is selected, meaning that box-coxtransformation is carried out when the ANOVA assumptions are notfulfilled by the raw or processed data (10% Trimming or outlierremoval) nor by log- or square-root transformation. Box-Coxtransformation is not mentioned in any of the TGs which is whyit is left to the choice of the user.
• alpha for residual normality and variance homogeneity tests
  • The default here is 0.05. In the TGs the alpha level isdiscussed so the user is adviced to inspect the requirements ofthe respective study to analyze.

The next option in the little box on the top-left of the app is to set ahook or remove the hook to decide whether the exposure well ID isregarded as random term in the underlying mixed ANOVA model. In case aREACTIV study is investigated the hook should be removed. For RADAR andXETA, the hook should be placed as long as information about theexposure well ID was documented. Of course, reducing the complexity ofthe random term might also make sense for XETA and RADAR studies in caseno variance is explained by the exposure well-ID. However, this choiceis left to the user.

When own data is supposed to be analyzed xlsx files with either spikedor unspiked data can be uploaded. Please make sure that the uploadeddata fulfills the Data requirements explained above. In case ofinsecurity, inspect the uploaded data structure of the pre-loadedexample data.

The app contains several output boxes about the required dataprocessing, a conclusion table, the uploaded data, residual diagnosticsplots, boxplots, a data summary table, the output of the residualnormality and homogeneity tests, the monotonicity test result, theDunnett’s and increasing and decreasing Williams test result tablesalong with the summary table of the underlying mixed ANOVA model.

The depicted information can be downloaded in a simple report byclicking on the button ** Word report**. It takes a couple ofseconds until the final docx file is produced and ready to download.Avoid clicking the button several times as this can lead to long waitingtimes and the production of several reports.

OECD. 2019a. Validation Report of the Xenopus Eleutheroembryonic ThyroidSignaling Assay (XETA) for the Detection of Thyroid Active Substances.OECD.

OECD. 2019b. TG 248: Xenopus Eleutheroembryonic Thyroid Assay (XETA).OECD.https://doi.org/10.1787/a13f80ee-en.

OECD. 2022b. Test No. 251: Rapid Androgen Disruption Activity Reporter(RADAR) Assay. OECD.https://doi.org/10.1787/da264d82-en.

OECD. 2022a. Rapid Estrogen Activity in Vivo (REACTIV) Assay (OECD DraftTG): Guideline for the Testing of Chemicals, Section 2: Effects onBiotic System. OECD.

Inka Marie Spyridonov, Lijuan Yan, Eduard Szöcs, Ana Filipa Pereira Miranda, Carsten Lange,Andrew Tindall, David Du Pasquier, Gregory Lemkine, Lennart Weltje, Maike Habekost,Pernille Thorbek. 2025. Xeredar: An open-source R-package for the statistical analysisof endocrine new approach methods (NAMs) using fish or amphibian eleutheroembryos,Environmental Toxicology and Chemistry.https://doi.org/10.1093/etojnl/vgaf056