Contents

1 Introduction

BERT (Batch-Effect Removal with Trees) offers flexible and efficient batch effect correction of omics data, while providing maximum tolerance to missing values. Tested on multiple datasets from proteomic analyses, BERT offered a typical 5-10x runtime improvement over existing methods, while retaining more numeric values and preserving batch effect reduction quality.

As such, BERT is a valuable preprocessing tool for data analysis workflows, in particular for proteomic data. By providing BERT via Bioconductor, we make this tool available to a wider research community. An accompanying research paper is currently under preparation and will be made public soon.

BERT addresses the same fundamental data integration challenges than the [HarmonizR][https://github.com/HSU-HPC/HarmonizR] package, which is released on Bioconductor in November 2023. However, various algorithmic modications and optimizations of BERT provide better execution time and better data coverage than HarmonizR. Moreover, BERT offers a more user-friendly design and a less error-prone input format.

Please note that our package BERT is neither affiliated with nor related to Bidirectional Encoder Representations from Transformers as published by Google.

Please report any questions and issues in the GitHub forum, the BioConductor forum or directly contact the authors,

2 Installation

Please download and install a current version of R (Windows binaries). You might want to consider installing a development environment as well, e.g. RStudio. Finally, BERT can be installed via Bioconductor using

if (!require("BiocManager", quietly = TRUE)){
    install.packages("BiocManager")
}
BiocManager::install("BERT")

which will install all required dependencies. To install the development version of BERT, you can use devtools as follows

devtools::install_github("HSU-HPC/BERT")

which may require the manual installation of the dependencies sva and limma.

if (!require("BiocManager", quietly = TRUE)){
    install.packages("BiocManager")
}
BiocManager::install("sva")
BiocManager::install("limma")

3 Data Preparation

As input, BERT requires a dataframe1 Matrices and SummarizedExperiments work as well, but will automatically be converted to dataframes. with samples in rows and features in columns. For each sample, the respective batch should be indicated by an integer or string in a corresponding column labelled Batch. Missing values should be labelled as NA. A valid example dataframe could look like this:

example = data.frame(feature_1 = stats::rnorm(5), feature_2 = stats::rnorm(5), Batch=c(1,1,2,2,2))
example
#>    feature_1  feature_2 Batch
#> 1 -0.7921692 -0.1067087     1
#> 2  1.4425722  0.4595086     1
#> 3  2.2972744  0.3883485     2
#> 4  1.8087603  0.0376165     2
#> 5 -0.3632649  0.3029304     2

Note that each batch should contain at least two samples. Optional columns that can be passed are

Note that BERT tries to find all metadata information for a SummarizedExperiment, including the mandatory batch information, using colData. For instance, a valid SummarizedExperiment might be defined as

nrows <- 200
ncols <- 8
expr_values <- matrix(runif(nrows * ncols, 1, 1e4), nrows)
# colData also takes all other metadata information, such as Label, Sample,
# Covariables etc.
colData <- data.frame(Batch=c(1,1,1,1,2,2,2,2), Reference=c(1,1,0,0,1,1,0,0))
dataset_raw = SummarizedExperiment::SummarizedExperiment(assays=list(expr=expr_values), colData=colData)

4 Basic Usage

BERT can be invoked by importing the BERT library and calling the BERT function. The batch effect corrected data is returned as a dataframe that mirrors the input dataframe2 In particular, the row and column names are in the same order and the optional columns are preserved..

library(BERT)
# generate test data with 10% missing values as provided by the BERT library
dataset_raw <- generate_dataset(features=60, batches=10, samplesperbatch=10, mvstmt=0.1, classes=2)
# apply BERT
dataset_adjusted <- BERT(dataset_raw)
#> 2024-02-28 16:06:17.842873 INFO::Formatting Data.
#> 2024-02-28 16:06:17.854089 INFO::Replacing NaNs with NAs.
#> 2024-02-28 16:06:17.864131 INFO::Removing potential empty rows and columns
#> 2024-02-28 16:06:18.274306 INFO::Found  600  missing values.
#> 2024-02-28 16:06:18.286837 INFO::Introduced 0 missing values due to singular proteins at batch/covariate level.
#> 2024-02-28 16:06:18.287525 INFO::Done
#> 2024-02-28 16:06:18.288051 INFO::Acquiring quality metrics before batch effect correction.
#> 2024-02-28 16:06:18.302702 INFO::Starting hierarchical adjustment
#> 2024-02-28 16:06:18.303832 INFO::Found  10  batches.
#> 2024-02-28 16:06:18.304423 INFO::Cores argument is not defined or BPPARAM has been specified. Argument corereduction will not be used.
#> 2024-02-28 16:06:20.686502 INFO::Using default BPPARAM
#> 2024-02-28 16:06:20.687153 INFO::Processing subtree level 1
#> 2024-02-28 16:06:22.792159 INFO::Processing subtree level 2
#> 2024-02-28 16:06:25.050168 INFO::Adjusting the last 1 batches sequentially
#> 2024-02-28 16:06:25.052703 INFO::Done
#> 2024-02-28 16:06:25.053671 INFO::Acquiring quality metrics after batch effect correction.
#> 2024-02-28 16:06:25.060941 INFO::ASW Batch was 0.441842244984297 prior to batch effect correction and is now -0.102003409929143 .
#> 2024-02-28 16:06:25.061859 INFO::ASW Label was 0.369041315939812 prior to batch effect correction and is now 0.843849474434861 .
#> 2024-02-28 16:06:25.063407 INFO::Total function execution time is  7.24793148040771  s and adjustment time is  6.74921202659607 s ( 93.12 )

BERT uses the logging library to convey live information to the user during the adjustment procedure. The algorithm first verifies the shape and suitability of the input dataframe (lines 1-6) before continuing with the actual batch effect correction (lines 8-14). BERT measure batch effects before and after the correction step by means of the average silhouette score (ASW) with respect to batch and labels (lines 7 and 15). The ASW Label should increase in a successful batch effect correction, whereas low values (\(\leq 0\)) are desireable for the ASW Batch3 The optimum of ASW Label is 1, which is typically however not achieved on real-world datasets. Also, the optimum of ASW Batch can vary, depending on the class distributions of the batches.. Finally, BERT prints the total function execution time (including the computation time for the quality metrics).

5 Advanced Options

5.1 Parameters

BERT offers a large number of parameters to customize the batch effect adjustment. The full function call, including all defaults is

BERT(data, cores = NULL, combatmode = 1, corereduction=2, stopParBatches=2, backend="default", method="ComBat", qualitycontrol=TRUE, verify=TRUE, labelname="Label", batchname="Batch", referencename="Reference", samplename="Sample", covariatename=NULL, BPPARAM=NULL, assayname=NULL)

In the following, we list the respective meaning of each parameter: - data: The input dataframe/matrix/SummarizedExperiment to adjust. See Data Preparation for detailed formatting instructions. - data The data for batch-effect correction. Must contain at least two samples per batch and 2 features.

  • cores: BERT uses BiocParallel for parallelization. If the user specifies a value cores, BERT internally creates and uses a new instance of BiocParallelParam, which is however not exhibited to the user. Setting this parameter can speed up the batch effect adjustment considerably, in particular for large datasets and on unix-based operating systems. A value between \(2\) and \(4\) is a reasonable choice for typical commodity hardware. Multi-node computations are not supported as of now. If, however, cores is not specified, BERT will default to BiocParallel::bpparam(), which may have been set by the user or the system. Additionally, the user can directly specify a specific instance of BiocParallelParam to be used via the BPPARAM argument.
  • combatmode An integer that encodes the parameters to use for ComBat.
Value par.prior mean.only
1 TRUE FALSE
2 TRUE TRUE
3 FALSE FALSE
4 FALSE TRUE

The value of this parameter will be ignored, if method!="ComBat".

  • corereduction Positive integer indicating the factor by which the number of processes should be reduced, once no further adjustment is possible for the current number of batches.4 E.g. consider a BERT call with 8 batches and 8 processes. Further adjustment is not possible with this number of processes, since batches are always processed in pairs. With corereduction=2, the number of processes for the following adjustment steps would be set to \(8/2=4\), which is the maximum number of usable processes for this example. This parameter is used only, if the user specified a custom value for parameter cores.

  • stopParBatches Positive integer indicating the minimum number of batches required at a hierarchy level to proceed with parallelized adjustment. If the number of batches is smaller, adjustment will be performed sequentially to avoid communication overheads.

  • backend: The backend to use for inter-process communication. Possible choices are default and file, where the former refers to the default communication backend of the requested parallelization mode and the latter will create temporary .rds files for data communication. ‘default’ is usually faster for small to medium sized datasets.

  • method: The method to use for the underlying batch effect correction steps. Should be either ComBat, limma for limma::removeBatchEffects or ref for adjustment using specified references (cf. Data Preparation). The underlying batch effect adjustment method for ref is a modified version of the limma method.

  • qualitycontrol: A boolean to (de)activate the ASW computation. Deactivating the ASW computations accelerates the computations.

  • verify: A boolean to (de)activate the initial format check of the input data. Deactivating this verification step accelerates the computations.

  • labelname: A string containing the name of the column to use as class labels. The default is “Label”.

  • batchname: A string containing the name of the column to use as batch labels. The default is “Batch”.

  • referencename: A string containing the name of the column to use as reference labels. The default is “Reference”.

  • covariatename: A vector containing the names of columns with categorical covariables.The default is NULL, in which case all column names are matched agains the pattern “Cov”.

  • BPPARAM: An instance of BiocParallelParam that will be used for parallelization. The default is null, in which case the value of cores determines the behaviour of BERT.

  • assayname: If the user chooses to pass a SummarizedExperiment object, they need to specify the name of the assay that they want to apply BERT to here. BERT then returns the input SummarizedExperiment with an additional assay labeled assayname_BERTcorrected.

5.2 Verbosity

BERT utilizes the logging package for output. The user can easily specify the verbosity of BERT by setting the global logging level in the script. For instance

logging::setLevel("WARN") # set level to warn and upwards
result <- BERT(data,cores = 1) # BERT executes silently

5.3 Choosing the Optimal Number of Cores

BERT exhibits a large number of parameters for parallelisation as to provide users with maximum flexibility. For typical scenarios, however, the default parameters are well suited. For very large experiments (\(>15\) batches), we recommend to increase the number of cores (a reasonable value is \(4\) but larger values may be possible on your hardware). Most users should leave all parameters to their respective default.

6 Examples

In the following, we present simple cookbook examples for BERT usage. Note that ASWs (and runtime) will most likely differ on your machine, since the data generating process involves multiple random choices.

6.1 Sequential Adjustment with limma

Here, BERT uses limma as underlying batch effect correction algorithm (method='limma') and performs all computations on a single process (cores parameter is left on default).

# import BERT
library(BERT)
# generate data with 30 batches, 60 features, 15 samples per batch, 15% missing values and 2 classes
dataset_raw <- generate_dataset(features=60, batches=20, samplesperbatch=15, mvstmt=0.15, classes=2)
# BERT
dataset_adjusted <- BERT(dataset_raw, method="limma")
#> 2024-02-28 16:06:25.126171 INFO::Formatting Data.
#> 2024-02-28 16:06:25.127303 INFO::Replacing NaNs with NAs.
#> 2024-02-28 16:06:25.129031 INFO::Removing potential empty rows and columns
#> 2024-02-28 16:06:25.132552 INFO::Found  2700  missing values.
#> 2024-02-28 16:06:25.179212 INFO::Introduced 0 missing values due to singular proteins at batch/covariate level.
#> 2024-02-28 16:06:25.180169 INFO::Done
#> 2024-02-28 16:06:25.180711 INFO::Acquiring quality metrics before batch effect correction.
#> 2024-02-28 16:06:25.195156 INFO::Starting hierarchical adjustment
#> 2024-02-28 16:06:25.196015 INFO::Found  20  batches.
#> 2024-02-28 16:06:25.196567 INFO::Cores argument is not defined or BPPARAM has been specified. Argument corereduction will not be used.
#> 2024-02-28 16:06:25.197161 INFO::Using default BPPARAM
#> 2024-02-28 16:06:25.197675 INFO::Processing subtree level 1
#> 2024-02-28 16:06:25.67324 INFO::Processing subtree level 2
#> 2024-02-28 16:06:26.071904 INFO::Processing subtree level 3
#> 2024-02-28 16:06:26.508667 INFO::Adjusting the last 1 batches sequentially
#> 2024-02-28 16:06:26.511161 INFO::Done
#> 2024-02-28 16:06:26.512017 INFO::Acquiring quality metrics after batch effect correction.
#> 2024-02-28 16:06:26.526957 INFO::ASW Batch was 0.450356828870141 prior to batch effect correction and is now -0.114965023930169 .
#> 2024-02-28 16:06:26.527959 INFO::ASW Label was 0.326079966512239 prior to batch effect correction and is now 0.82650368616144 .
#> 2024-02-28 16:06:26.52935 INFO::Total function execution time is  1.40303564071655  s and adjustment time is  1.31528520584106 s ( 93.75 )

6.2 Parallel Batch Effect Correction with ComBat

Here, BERT uses ComBat as underlying batch effect correction algorithm (method is left on default) and performs all computations on a 2 processes (cores=2).

# import BERT
library(BERT)
# generate data with 30 batches, 60 features, 15 samples per batch, 15% missing values and 2 classes
dataset_raw <- generate_dataset(features=60, batches=20, samplesperbatch=15, mvstmt=0.15, classes=2)
# BERT
dataset_adjusted <- BERT(dataset_raw, cores=2)
#> 2024-02-28 16:06:26.568378 INFO::Formatting Data.
#> 2024-02-28 16:06:26.569211 INFO::Replacing NaNs with NAs.
#> 2024-02-28 16:06:26.570369 INFO::Removing potential empty rows and columns
#> 2024-02-28 16:06:26.572804 INFO::Found  2700  missing values.
#> 2024-02-28 16:06:26.601013 INFO::Introduced 0 missing values due to singular proteins at batch/covariate level.
#> 2024-02-28 16:06:26.601885 INFO::Done
#> 2024-02-28 16:06:26.602515 INFO::Acquiring quality metrics before batch effect correction.
#> 2024-02-28 16:06:26.616975 INFO::Starting hierarchical adjustment
#> 2024-02-28 16:06:26.617969 INFO::Found  20  batches.
#> 2024-02-28 16:06:27.324286 INFO::Set up parallel execution backend with 2 workers
#> 2024-02-28 16:06:27.325238 INFO::Processing subtree level 1 with 20 batches using 2 cores.
#> 2024-02-28 16:06:30.370398 INFO::Adjusting the last 2 batches sequentially
#> 2024-02-28 16:06:30.371657 INFO::Adjusting sequential tree level 1 with 2 batches
#> 2024-02-28 16:06:32.31071 INFO::Done
#> 2024-02-28 16:06:32.311412 INFO::Acquiring quality metrics after batch effect correction.
#> 2024-02-28 16:06:32.321503 INFO::ASW Batch was 0.467878149378763 prior to batch effect correction and is now -0.124894505274773 .
#> 2024-02-28 16:06:32.322171 INFO::ASW Label was 0.31909634341049 prior to batch effect correction and is now 0.815830679676179 .
#> 2024-02-28 16:06:32.322949 INFO::Total function execution time is  5.75468635559082  s and adjustment time is  5.6924774646759 s ( 98.92 )

6.3 Batch Effect Correction Using SummarizedExperiment

Here, BERT takes the input data using a SummarizedExperiment instead. Batch effect correction is then performed using ComBat as underlying algorithm (method is left on default) and all computations are performed on a single process (cores parameter is left on default).

nrows <- 200
ncols <- 8
# SummarizedExperiments store samples in columns and features in rows (in contrast to BERT).
# BERT will automatically account for this.
expr_values <- matrix(runif(nrows * ncols, 1, 1e4), nrows)
# colData also takes further metadata information, such as Label, Sample,
# Reference or Covariables
colData <- data.frame("Batch"=c(1,1,1,1,2,2,2,2), "Label"=c(1,2,1,2,1,2,1,2), "Sample"=c(1,2,3,4,5,6,7,8))
dataset_raw = SummarizedExperiment::SummarizedExperiment(assays=list(expr=expr_values), colData=colData)
dataset_adjusted = BERT(dataset_raw, assayname = "expr")
#> 2024-02-28 16:06:32.376634 INFO::Formatting Data.
#> 2024-02-28 16:06:32.377298 INFO::Recognized SummarizedExperiment
#> 2024-02-28 16:06:32.377831 INFO::Typecasting input to dataframe.
#> 2024-02-28 16:06:32.415934 INFO::Replacing NaNs with NAs.
#> 2024-02-28 16:06:32.417102 INFO::Removing potential empty rows and columns
#> 2024-02-28 16:06:32.420653 INFO::Found  0  missing values.
#> 2024-02-28 16:06:32.427258 INFO::Introduced 0 missing values due to singular proteins at batch/covariate level.
#> 2024-02-28 16:06:32.42783 INFO::Done
#> 2024-02-28 16:06:32.428314 INFO::Acquiring quality metrics before batch effect correction.
#> 2024-02-28 16:06:32.43232 INFO::Starting hierarchical adjustment
#> 2024-02-28 16:06:32.432999 INFO::Found  2  batches.
#> 2024-02-28 16:06:32.433504 INFO::Cores argument is not defined or BPPARAM has been specified. Argument corereduction will not be used.
#> 2024-02-28 16:06:32.434074 INFO::Using default BPPARAM
#> 2024-02-28 16:06:32.43455 INFO::Adjusting the last 2 batches sequentially
#> 2024-02-28 16:06:32.435462 INFO::Adjusting sequential tree level 1 with 2 batches
#> 2024-02-28 16:06:32.485209 INFO::Done
#> 2024-02-28 16:06:32.485913 INFO::Acquiring quality metrics after batch effect correction.
#> 2024-02-28 16:06:32.49025 INFO::ASW Batch was 0.0021745411734109 prior to batch effect correction and is now -0.0989352634791682 .
#> 2024-02-28 16:06:32.490855 INFO::ASW Label was 0.0209102908725772 prior to batch effect correction and is now 0.0398206453788259 .
#> 2024-02-28 16:06:32.491648 INFO::Total function execution time is  0.114986658096313  s and adjustment time is  0.052311897277832 s ( 45.49 )

6.4 BERT with Covariables

BERT can utilize categorical covariables that are specified in columns Cov_1, Cov_2, .... These columns are automatically detected and integrated into the batch effect correction process.

# import BERT
library(BERT)
# set seed for reproducibility
set.seed(1)
# generate data with 5 batches, 60 features, 30 samples per batch, 15% missing values and 2 classes
dataset_raw <- generate_dataset(features=60, batches=5, samplesperbatch=30, mvstmt=0.15, classes=2)
# create covariable column with 2 possible values, e.g. male/female condition
dataset_raw["Cov_1"] = sample(c(1,2), size=dim(dataset_raw)[1], replace=TRUE)
# BERT
dataset_adjusted <- BERT(dataset_raw)
#> 2024-02-28 16:06:32.51918 INFO::Formatting Data.
#> 2024-02-28 16:06:32.519946 INFO::Replacing NaNs with NAs.
#> 2024-02-28 16:06:32.520841 INFO::Removing potential empty rows and columns
#> 2024-02-28 16:06:32.52293 INFO::Found  1350  missing values.
#> 2024-02-28 16:06:32.523873 INFO::BERT requires at least 2 numeric values per batch/covariate level. This may reduce the number of adjustable features considerably, depending on the quantification technique.
#> 2024-02-28 16:06:32.544372 INFO::Introduced 0 missing values due to singular proteins at batch/covariate level.
#> 2024-02-28 16:06:32.545112 INFO::Done
#> 2024-02-28 16:06:32.545665 INFO::Acquiring quality metrics before batch effect correction.
#> 2024-02-28 16:06:32.550964 INFO::Starting hierarchical adjustment
#> 2024-02-28 16:06:32.551697 INFO::Found  5  batches.
#> 2024-02-28 16:06:32.552231 INFO::Cores argument is not defined or BPPARAM has been specified. Argument corereduction will not be used.
#> 2024-02-28 16:06:32.552824 INFO::Using default BPPARAM
#> 2024-02-28 16:06:32.55334 INFO::Processing subtree level 1
#> 2024-02-28 16:06:32.732701 INFO::Adjusting the last 2 batches sequentially
#> 2024-02-28 16:06:32.734551 INFO::Adjusting sequential tree level 1 with 2 batches
#> 2024-02-28 16:06:32.788041 INFO::Done
#> 2024-02-28 16:06:32.788756 INFO::Acquiring quality metrics after batch effect correction.
#> 2024-02-28 16:06:32.794117 INFO::ASW Batch was 0.492773245691086 prior to batch effect correction and is now -0.0377157224767566 .
#> 2024-02-28 16:06:32.794734 INFO::ASW Label was 0.40854766060101 prior to batch effect correction and is now 0.895560693013661 .
#> 2024-02-28 16:06:32.795555 INFO::Total function execution time is  0.27644157409668  s and adjustment time is  0.236456155776978 s ( 85.54 )

6.5 BERT with references

In rare cases, class distributions across experiments may be severely skewed. In particular, a batch might contain classes that other batches don’t contain. In these cases, samples of common conditions may serve as references (bridges) between the batches (method="ref"). BERT utilizes those samples as references that have a condition specified in the “Reference” column of the input. All other samples are co-adjusted. Please note, that this strategy implicitly uses limma as underlying batch effect correction algorithm.

# import BERT
library(BERT)
# generate data with 4 batches, 6 features, 15 samples per batch, 15% missing values and 2 classes
dataset_raw <- generate_dataset(features=6, batches=4, samplesperbatch=15, mvstmt=0.15, classes=2)
# create reference column with default value 0.  The 0 indicates, that the respective sample should be co-adjusted only.
dataset_raw[, "Reference"] <- 0
# randomly select 2 references per batch and class - in practice, this choice will be determined by external requirements (e.g. class known for only these samples)
batches <- unique(dataset_raw$Batch) # all the batches
for(b in batches){ # iterate over all batches
    # references from class 1
    ref_idx = sample(which((dataset_raw$Batch==b)&(dataset_raw$Label==1)), size=2, replace=FALSE)
    dataset_raw[ref_idx, "Reference"] <- 1
    # references from class 2
    ref_idx = sample(which((dataset_raw$Batch==b)&(dataset_raw$Label==2)), size=2, replace=FALSE)
    dataset_raw[ref_idx, "Reference"] <- 2
}
# BERT
dataset_adjusted <- BERT(dataset_raw, method="ref")
#> 2024-02-28 16:06:32.837786 INFO::Formatting Data.
#> 2024-02-28 16:06:32.838547 INFO::Replacing NaNs with NAs.
#> 2024-02-28 16:06:32.839455 INFO::Removing potential empty rows and columns
#> 2024-02-28 16:06:32.84046 INFO::Found  60  missing values.
#> 2024-02-28 16:06:32.844309 INFO::Introduced 0 missing values due to singular proteins at batch/covariate level.
#> 2024-02-28 16:06:32.844835 INFO::Done
#> 2024-02-28 16:06:32.845345 INFO::Acquiring quality metrics before batch effect correction.
#> 2024-02-28 16:06:32.848343 INFO::Starting hierarchical adjustment
#> 2024-02-28 16:06:32.849016 INFO::Found  4  batches.
#> 2024-02-28 16:06:32.849537 INFO::Cores argument is not defined or BPPARAM has been specified. Argument corereduction will not be used.
#> 2024-02-28 16:06:32.850141 INFO::Using default BPPARAM
#> 2024-02-28 16:06:32.850672 INFO::Processing subtree level 1
#> 2024-02-28 16:06:32.944821 INFO::Adjusting the last 2 batches sequentially
#> 2024-02-28 16:06:32.946721 INFO::Adjusting sequential tree level 1 with 2 batches
#> 2024-02-28 16:06:32.976827 INFO::Done
#> 2024-02-28 16:06:32.977737 INFO::Acquiring quality metrics after batch effect correction.
#> 2024-02-28 16:06:32.981155 INFO::ASW Batch was 0.440355021914032 prior to batch effect correction and is now -0.087480278736629 .
#> 2024-02-28 16:06:32.981758 INFO::ASW Label was 0.373906827748893 prior to batch effect correction and is now 0.919791677398366 .
#> 2024-02-28 16:06:32.982553 INFO::Total function execution time is  0.144831657409668  s and adjustment time is  0.127890825271606 s ( 88.3 )

7 Issues

Issues can be reported in the GitHub forum, the BioConductor forum or directly to the authors.

8 License

This code is published under the GPLv3.0 License and is available for non-commercial academic purposes.

9 Reference

Please cite our manuscript, if you use BERT for your research: Yannis Schumann, Simon Schlumbohm et al., BERT - Batch Effect Reduction Trees with Missing Value Tolerance, 2023

10 Session Info

sessionInfo()
#> R Under development (unstable) (2024-01-16 r85808)
#> Platform: x86_64-pc-linux-gnu
#> Running under: Ubuntu 22.04.3 LTS
#> 
#> Matrix products: default
#> BLAS:   /home/biocbuild/bbs-3.19-bioc/R/lib/libRblas.so 
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_GB              LC_COLLATE=C              
#>  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
#>  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
#>  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
#> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
#> 
#> time zone: America/New_York
#> tzcode source: system (glibc)
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] BERT_0.99.16     BiocStyle_2.31.0
#> 
#> loaded via a namespace (and not attached):
#>  [1] tidyselect_1.2.0            blob_1.2.4                 
#>  [3] Biostrings_2.71.2           bitops_1.0-7               
#>  [5] fastmap_1.1.1               RCurl_1.98-1.14            
#>  [7] janitor_2.2.0               XML_3.99-0.16.1            
#>  [9] digest_0.6.34               timechange_0.3.0           
#> [11] lifecycle_1.0.4             cluster_2.1.6              
#> [13] statmod_1.5.0               survival_3.5-8             
#> [15] KEGGREST_1.43.0             invgamma_1.1               
#> [17] RSQLite_2.3.5               magrittr_2.0.3             
#> [19] genefilter_1.85.1           compiler_4.4.0             
#> [21] rlang_1.1.3                 sass_0.4.8                 
#> [23] tools_4.4.0                 yaml_2.3.8                 
#> [25] knitr_1.45                  S4Arrays_1.3.4             
#> [27] bit_4.0.5                   DelayedArray_0.29.8        
#> [29] abind_1.4-5                 BiocParallel_1.37.0        
#> [31] BiocGenerics_0.49.1         grid_4.4.0                 
#> [33] stats4_4.4.0                xtable_1.8-4               
#> [35] edgeR_4.1.17                iterators_1.0.14           
#> [37] logging_0.10-108            SummarizedExperiment_1.33.3
#> [39] cli_3.6.2                   rmarkdown_2.25             
#> [41] crayon_1.5.2                generics_0.1.3             
#> [43] httr_1.4.7                  DBI_1.2.2                  
#> [45] cachem_1.0.8                stringr_1.5.1              
#> [47] zlibbioc_1.49.0             splines_4.4.0              
#> [49] parallel_4.4.0              AnnotationDbi_1.65.2       
#> [51] BiocManager_1.30.22         XVector_0.43.1             
#> [53] matrixStats_1.2.0           vctrs_0.6.5                
#> [55] Matrix_1.6-5                jsonlite_1.8.8             
#> [57] sva_3.51.0                  bookdown_0.37              
#> [59] comprehenr_0.6.10           IRanges_2.37.1             
#> [61] S4Vectors_0.41.3            bit64_4.0.5                
#> [63] locfit_1.5-9.8              foreach_1.5.2              
#> [65] limma_3.59.3                jquerylib_0.1.4            
#> [67] annotate_1.81.2             glue_1.7.0                 
#> [69] codetools_0.2-19            lubridate_1.9.3            
#> [71] stringi_1.8.3               GenomeInfoDb_1.39.6        
#> [73] GenomicRanges_1.55.3        htmltools_0.5.7            
#> [75] GenomeInfoDbData_1.2.11     R6_2.5.1                   
#> [77] evaluate_0.23               lattice_0.22-5             
#> [79] Biobase_2.63.0              png_0.1-8                  
#> [81] memoise_2.0.1               snakecase_0.11.1           
#> [83] bslib_0.6.1                 Rcpp_1.0.12                
#> [85] SparseArray_1.3.4           nlme_3.1-164               
#> [87] mgcv_1.9-1                  xfun_0.42                  
#> [89] MatrixGenerics_1.15.0