Introduction
Bioconductor objects are useful for their shared structure within the project, and that they enable rich metadata.
How then should one save a Bioconductor object so that a collaborator can load and use it? Or so that it persists reliably across time?
The answer depends on several factors:
- Longevity: will the file need to be easily readable in 1 year, or 5-10 years from now?
- Language interoperability: does the recipient use R, Python, or both?
- Object size: is the object small enough to serialize entirely, or does it contain large on-disk arrays?
- Is it easy enough to write out in common formats (BED) or plaintext?
- Reproducibility: should the saved form capture provenance and metadata alongside data?
This vignette walks through the main options, their trade-offs, and some recommendations for common scenarios.
Acknowledgments
The content of this vignette has been informed by discussions on the Bioconductor community Zulip, including contributions from Kevin Rue-Albrecht, Johannes Rainer, Lori Shepherd, Jayaram Kancherla, Aaron Lun, Hervé Pagès, Laurent Gatto, Vince Carey, Sean Davis, Robert Castelo, Hugo Gruson, and Luke Zappia.
Version Info
R version: R version 4.6.1 (2026-06-24)
Bioconductor version: 3.23
Package version: 0.99.1
R serialization
saveRDS and readRDS
The simplest approach is R’s built-in binary serialization. saveRDS() saves a single object to an RDS file (the acronym is not defined in the man pages but likely stands for “R Data Serialization”); save() bundles one or more named objects into an .RData (or .rda) file.
library(SummarizedExperiment)
se <- SummarizedExperiment(
assays = list(counts = matrix(1:12, nrow = 3)),
colData = DataFrame(
condition = c("A", "A", "B", "B"),
row.names=1:4
),
rowData = DataFrame(
gene = c("gene1","gene2","gene3"),
row.names=1:3
)
)
# Single object
tmp_rds <- tempfile(fileext = ".rds")
saveRDS(se, file = tmp_rds)
se_from_rds <- readRDS(tmp_rds)
se_from_rds
class: SummarizedExperiment
dim: 3 4
metadata(0):
assays(1): counts
rownames(3): 1 2 3
rowData names(1): gene
colnames(4): 1 2 3 4
colData names(1): condition
# Multiple objects in one file
tmp_rda <- tempfile(fileext = ".RData")
save(se, file = tmp_rda)
load(tmp_rda) # restores 'se' by name into the current environment
Prefer saveRDS() for most use cases. save() / load() silently overwrites any object in the calling environment that shares a name, which is a common source of confusion. The main remaining use case for save() is .rda data files shipped inside R packages (under data/).
Advantages:
- Zero setup — works with any R object.
- Compression applied automatically (gzip by default).
- Round-trips perfectly: the loaded object is identical to the saved one.
Disadvantages:
- R-only: Python or other languages cannot read these files without a bridge library.
- Large objects are loaded entirely into RAM; there is no lazy / on-disk access.
When to use it: quick sharing between R users on the same project, saving intermediate objects in a pipeline, anything under ~1 GB.
If your workflow downloads RDS files from a remote URL, consider using BiocFileCache to cache them locally so they are only fetched once.
Reading RDS files in Python
BiocPy is a Python ecosystem that brings Bioconductor’s core data structures to Python, including BiocFrame, IRanges, GenomicRanges, SummarizedExperiment, SingleCellExperiment, and MultiAssayExperiment.
Python users can read RDS files with the rds2py package from the BiocPy ecosystem. Standard R types map to NumPy/SciPy equivalents (e.g. numeric vectors become numpy.ndarray), and Bioconductor classes such as SummarizedExperiment, SingleCellExperiment, GRanges, and MultiAssayExperiment are converted to their BiocPy counterparts. For unrecognised S4 classes the object falls back to a dictionary so no data is lost. Support for writing RDS files from Python is also in development.
Cross-release stability
There is no guarantee that an S4 object serialized today will work with older versions of Bioconductor. The most common reasons are that the S4 class is not defined at all in the earlier version, or that it exists but its definition has changed — new slots added, slots renamed or removed, or infrastructure moved between packages.
A concrete example: in Bioconductor 3.22, Seqinfo was moved out of GenomicRanges into its own package. A GRanges serialized under BioC ≥ 3.22 can still be loaded on a machine running BioC < 3.22, and many basic operations work — show(), seqnames(), ranges(), mcols(), and [ among them. But operations that require the Seqinfo package (such as shift() or reduce()) will fail with a confusing “package not available” error. The practical advice is to upgrade to BioC ≥ 3.22. The reverse was also true: older GRanges objects loaded into a newer session sometimes required updateObject() to migrate the internal representation:
# not evaluated — requires an object saved under an older Bioconductor release
gr <- readRDS("old_granges.rds")
gr <- updateObject(gr, verbose = TRUE)
The general lesson is: if you can save your data in a format that is not tied to a particular version of Bioconductor — or better, not tied to R at all — you should. For a GRanges, for instance, a simple TSV is often enough:
library(GenomicRanges)
gr <- GRanges(
seqnames = "chr1",
ranges = IRanges(start = c(100, 200, 300), width = 50),
seqinfo = Seqinfo(seqnames = "chr1", seqlengths = 248956422,
isCircular = FALSE, genome = "hg38")
)
names(gr) <- c("peak1", "peak2", "peak3")
gr$score <- c(500, 800, 300) # standard BED score column
gr$log2fc <- c(1.2, -0.5, 2.1) # extra metadata column
tmp_tsv <- tempfile(fileext = ".tsv")
write.table(as.data.frame(gr), tmp_tsv, sep = "\t", quote = FALSE)
gr_from_tsv <- makeGRangesFromDataFrame(
read.table(tmp_tsv, header = TRUE, sep = "\t"),
keep.extra.columns = TRUE
)
gr_from_tsv
GRanges object with 3 ranges and 2 metadata columns:
seqnames ranges strand | score log2fc
<Rle> <IRanges> <Rle> | <integer> <numeric>
peak1 chr1 100-149 * | 500 1.2
peak2 chr1 200-249 * | 800 -0.5
peak3 chr1 300-349 * | 300 2.1
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
This round-trip survives any Bioconductor version and is readable from Python or the command line. The alabaster ecosystem (described below) applies the same principle more systematically and with better support for complex objects, using HDF5 and JSON as the underlying storage formats.
HDF5-backed storage
Saving with HDF5Array
For large assay matrices (e.g., single-cell count matrices with millions of cells), it is impractical to hold the entire object in RAM. The HDF5Array package provides array classes backed by HDF5 files, enabling lazy loading and out-of-memory computation.
class: SummarizedExperiment
dim: 3 4
metadata(0):
assays(1): counts
rownames(3): 1 2 3
rowData names(1): gene
colnames(4): 1 2 3 4
colData names(1): condition
The saved directory contains an HDF5 file with the assay data and an RDS file for the non-assay metadata.
Advantages:
- Assay data is stored on disk; only the chunks you access are read into RAM.
- HDF5 is a widely used binary format with readers in Python (
h5py, anndata), Julia, C/C++, and more.
- Good for objects with tens of gigabytes of assay data.
Disadvantages:
- The output is a directory, not a single file, which complicates transfer (use
tar or zip before sharing).
- The RDS envelope for metadata is still R-specific.
- Write performance can be slower than
saveRDS() for small objects.
When to use it: large single-cell or spatial datasets where you want on-disk access; workflows shared between R users who need memory efficiency.
SummarizedExperiment and AnnData
For single-cell workflows that move between R and Python, the .h5ad format used by scanpy and related tools is often the most convenient path when collaborators are working in Python. Like HDF5Array, .h5ad is an HDF5-based format. The trade-off relative to alabaster is that .h5ad is AnnData-specific rather than a general Bioconductor serialization format.
The recommended starting point is the anndataR package, a more recent and complete implementation that is actively maintained as part of the scverse ecosystem:
class: SingleCellExperiment
dim: 3 4
metadata(0):
assays(1): counts
rownames(3): 1 2 3
rowData names(1): gene
colnames(4): 1 2 3 4
colData names(1): condition
reducedDimNames(0):
mainExpName: NULL
altExpNames(0):
The zellkonverter package is an alternative that also converts directly between SingleCellExperiment and .h5ad. It remains actively maintained and has advantages in some cases, though at the cost of managing a Python environment via basilisk.
# not evaluated — zellkonverter installs a full Python environment via basilisk
# on first use, which takes too long in CI
library(zellkonverter)
writeH5AD(sce, file = "sce.h5ad")
sce_from_h5ad <- readH5AD("sce.h5ad")
sce_from_h5ad
The alabaster ecosystem
saveObject and readObject
The alabaster family of packages is part of the broader ArtifactDB project, which provides a multi-language system for storing and retrieving analysis-ready Bioconductor objects. The core idea is to save objects as directories of standard files (HDF5, JSON, CSV) whose format is defined by explicit, versioned specifications — meaning the saved form is readable without R, and can evolve over time without breaking previously saved objects.
class: SummarizedExperiment
dim: 3 4
metadata(0):
assays(1): counts
rownames(3): 1 2 3
rowData names(1): gene
colnames(4): 1 2 3 4
colData names(1): condition
The alabaster umbrella package pulls in support for the most common Bioconductor classes. Individual sub-packages cover specific classes: alabaster.se for SummarizedExperiment, alabaster.sce for SingleCellExperiment, and so on.
Validation with takane
A key part of the ArtifactDB design is that saved directories can be independently validated against the format specification. This is handled by takane, a C++ library that maintains separate, versioned specifications for 30+ Bioconductor object types. Calling takane::validate() on a saved directory checks that all files conform to the expected layout and types, which means a collaborator or downstream tool can verify the integrity of a saved object without needing to load it into R. This makes alabaster directories suitable for deposition in data repositories where format conformance needs to be auditable.
The Python counterpart to alabaster is the dolomite family of packages, which reads and writes the same on-disk format. An object saved with alabaster in R can be read with dolomite in Python, and vice versa, with no conversion step.
Advantages:
- Truly language-agnostic: Python readers exist via the
dolomite family of packages, enabling seamless R ↔︎ Python interoperability.
- Built on open standards (HDF5, JSON); inspectable without R.
- Versioned format specifications with independent validation via takane.
- A good choice for archives or data portals.
Disadvantages:
- Newer ecosystem; not all Bioconductor classes have
alabaster support yet.
- Requires installing the relevant
alabaster.* sub-package for each class.
- Like HDF5Array, the output is a directory.
When to use it: archival storage, data portal submissions, cross-language workflows, or any situation where you want the saved format to be readable without R.
Writing BED files
When the object is a GRanges or similar ranges object and the goal is interoperability with other tools (genome browsers, Python, command-line utilities), exporting to BED format is often more useful than R-specific serialization. Our gr has range names, a standard BED score column, and an extra metadata column log2fc:
GRanges object with 3 ranges and 2 metadata columns:
seqnames ranges strand | score log2fc
<Rle> <IRanges> <Rle> | <numeric> <numeric>
peak1 chr1 100-149 * | 500 1.2
peak2 chr1 200-249 * | 800 -0.5
peak3 chr1 300-349 * | 300 2.1
-------
seqinfo: 1 sequence from hg38 genome
Both rtracklayer and plyranges can write BED files:
When reading back, the standard score column is preserved, but log2fc is silently dropped — BED has no mechanism to carry arbitrary metadata columns. Range names are stored in the BED name field but come back as a $name metadata column rather than as R names on the object; restore them manually:
gr_rtracklayer <- import(tmp_bed)
names(gr_rtracklayer) <- gr_rtracklayer$name
gr_rtracklayer$name <- NULL
gr_rtracklayer
GRanges object with 3 ranges and 1 metadata column:
seqnames ranges strand | score
<Rle> <IRanges> <Rle> | <numeric>
peak1 chr1 100-149 * | 500
peak2 chr1 200-249 * | 800
peak3 chr1 300-349 * | 300
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
gr_plyranges <- read_bed(tmp_bed2)
names(gr_plyranges) <- gr_plyranges$name
gr_plyranges$name <- NULL
gr_plyranges
GRanges object with 3 ranges and 1 metadata column:
seqnames ranges strand | score
<Rle> <IRanges> <Rle> | <numeric>
peak1 chr1 100-149 * | 500
peak2 chr1 200-249 * | 800
peak3 chr1 300-349 * | 300
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
If preserving all metadata columns is the priority and BED compatibility is not required, the simplest approach is the TSV round-trip shown earlier: write.table(as.data.frame(gr), ...) followed by makeGRangesFromDataFrame(..., keep.extra.columns = TRUE) restores all mcols in one step without a sidecar.
When you do need a BED file (e.g. for a genome browser or a tool that expects BED input), extra mcols can be preserved by writing them to a sidecar file. Here using plyranges to read the BED back, then reattaching from the sidecar:
tmp_meta <- tempfile(fileext = ".tsv")
write.table(
data.frame(name = names(gr), log2fc = gr$log2fc),
tmp_meta, sep = "\t", quote = FALSE, row.names = FALSE
)
gr_restored <- read_bed(tmp_bed2)
meta <- read.table(tmp_meta, header = TRUE, sep = "\t")
gr_restored$log2fc <- meta$log2fc
gr_restored
GRanges object with 3 ranges and 3 metadata columns:
seqnames ranges strand | name score log2fc
<Rle> <IRanges> <Rle> | <character> <numeric> <numeric>
[1] chr1 100-149 * | peak1 500 1.2
[2] chr1 200-249 * | peak2 800 -0.5
[3] chr1 300-349 * | peak3 300 2.1
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
Saving Seqinfo separately
BED files do not store chromosome lengths or genome build information, so Seqinfo is silently dropped on export. To preserve it, write it out alongside the BED file and restore it on load:
Seqinfo object with 1 sequence from hg38 genome:
seqnames seqlengths isCircular genome
chr1 248956422 FALSE hg38
Object-level metadata stored in metadata(object) — things like processing parameters, provenance notes, or experiment descriptors — is lost in any format that only encodes the ranges or assay data. This applies whether you are writing a BED file, an HDF5 matrix, or any other non-R format. Write it to a JSON sidecar file so it travels with the data:
$timestamp
[1] "2020-01-01 12:00:00"
$pipeline
[1] "v2.1"
$n_samples
[1] 4
toJSON handles simple R types (lists, vectors, data frames) well, but complex objects (S4 instances, environments) need to be simplified or omitted before serializing.
Summary and recommendations
In most new projects we recommend defaulting to saveRDS() for convenience and upgrading to alabaster when cross-language access or archival stability becomes a priority. Be aware that any R-serialized Bioconductor object may require updateObject() when loaded under a different Bioconductor release.
Session info
R version 4.6.1 (2026-06-24)
Platform: x86_64-pc-linux-gnu
Running under: Ubuntu 24.04.5 LTS
Matrix products: default
BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so; LAPACK version 3.12.0
locale:
[1] LC_CTYPE=C.UTF-8 LC_NUMERIC=C LC_TIME=C.UTF-8
[4] LC_COLLATE=C.UTF-8 LC_MONETARY=C.UTF-8 LC_MESSAGES=C.UTF-8
[7] LC_PAPER=C.UTF-8 LC_NAME=C LC_ADDRESS=C
[10] LC_TELEPHONE=C LC_MEASUREMENT=C.UTF-8 LC_IDENTIFICATION=C
time zone: UTC
tzcode source: system (glibc)
attached base packages:
[1] stats4 stats graphics grDevices utils datasets methods
[8] base
other attached packages:
[1] jsonlite_2.0.0 plyranges_1.32.0
[3] dplyr_1.2.1 rtracklayer_1.72.0
[5] alabaster.se_1.12.0 alabaster.base_1.12.1
[7] SingleCellExperiment_1.34.0 anndataR_1.2.1
[9] HDF5Array_1.40.0 h5mread_1.4.1
[11] rhdf5_2.56.1 DelayedArray_0.38.2
[13] SparseArray_1.12.2 S4Arrays_1.12.0
[15] abind_1.4-8 Matrix_1.7-5
[17] SummarizedExperiment_1.42.0 Biobase_2.72.0
[19] GenomicRanges_1.64.0 Seqinfo_1.2.0
[21] IRanges_2.46.0 S4Vectors_0.50.3
[23] BiocGenerics_0.58.1 generics_0.1.4
[25] MatrixGenerics_1.24.0 matrixStats_1.5.0
[27] BiocManager_1.30.27
loaded via a namespace (and not attached):
[1] rjson_0.2.23 xfun_0.61 lattice_0.22-9
[4] rhdf5filters_1.24.1 vctrs_0.7.3 tools_4.6.1
[7] bitops_1.1-0 curl_8.0.0 parallel_4.6.1
[10] tibble_3.3.1 pkgconfig_2.0.3 cigarillo_1.2.1
[13] lifecycle_1.0.5 compiler_4.6.1 Rsamtools_2.28.0
[16] Biostrings_2.80.2 codetools_0.2-20 htmltools_0.5.9
[19] RCurl_1.98-1.20 alabaster.matrix_1.12.0 yaml_2.3.12
[22] pillar_1.11.1 crayon_1.5.3 BiocParallel_1.46.0
[25] tidyselect_1.2.1 digest_0.6.39 purrr_1.2.2
[28] restfulr_0.0.17 fastmap_1.2.0 grid_4.6.1
[31] cli_3.6.6 magrittr_2.0.5 XML_3.99-0.24
[34] rmarkdown_2.32 XVector_0.52.0 httr_1.4.9
[37] otel_0.2.0 reticulate_1.47.0 png_0.1-9
[40] evaluate_1.0.5 knitr_1.52 BiocIO_1.22.0
[43] rlang_1.3.0 Rcpp_1.1.2 glue_1.8.1
[46] alabaster.ranges_1.12.0 alabaster.schemas_1.12.0 R6_2.6.1
[49] Rhdf5lib_2.0.0 GenomicAlignments_1.48.0