R for Bioinformatics: Bioconductor, DESeq2, and What R Is Actually For

Quick answer: R is the statistics language of bioinformatics. Bioconductor — the field’s package repository with synchronized twice-yearly releases — standardizes data containers (a SummarizedExperiment holds the count matrix, the sample metadata and the gene annotation in one aligned object) and ships the analysis canon: DESeq2, edgeR and limma-voom for differential expression, Seurat for single-cell, DiffBind for ChIP-seq, clusterProfiler for enrichment. Nearly every transcriptomics methods section cites a Bioconductor tool, which is why R is non-negotiable for the statistical half of the field — while Python owns the file-wrangling and machine-learning half. Learn R through Bioconductor, not through generic data-science courses.

The companion to Python for bioinformatics and the Biopython tutorial on this blog: Python reads and transforms sequencing data; R tests and plots it. Most real pipelines use both, in that order.

Why bioinformatics standardized on R

R was built for statisticians — it descends from the S language, is vectorized by design, and treats a statistical model as a first-class object — and bioinformatics arrived at the moment genomics started needing statistics at population scale. Two decisions then locked the standard in place:

• Bioconductor, a domain repository with its own release discipline. Since 2001, Bioconductor has coordinated package releases twice a year, synchronized with R versions, with enforced standards for documentation, worked vignettes and — most importantly — shared data containers. Where CRAN is a bazaar, Bioconductor is a curated instrument collection: a Bioconductor package arrives with a tutorial run on real biological data.
• The methods-section economy. A statistical result is publishable when its method is citable. DESeq2, edgeR and limma have canonical papers with hundreds of thousands of citations between them; using them means your methods section writes itself and reviewers recognize every step. An equivalent script written from scratch in any other language has no citation and invites exactly the questions you do not want in review.

The practical consequence: whenever an analysis touches statistical testing on biological data, R is not a style choice — it is the field’s shared instrument, and the fluency is worth acquiring in that specific dialect.

The data model that makes R different: SummarizedExperiment

The single most important thing to understand before Bioconductor makes sense: data travel in typed containers, not loose data frames.

A SummarizedExperiment bundles three coordinated pieces:

• assay — the matrix, genes × samples: counts, intensities, methylation betas;
• colData — sample metadata, rows aligned to the assay’s columns: condition, batch, patient;
• rowData — feature annotation, rows aligned to the assay’s rows: gene symbols, biotypes, coordinates.

Subsetting is coherent by construction: pulling ten samples pulls their metadata with them; pulling fifty genes pulls their annotations. Batch, condition and pairing live with the data instead of in a parallel spreadsheet that silently desynchronizes — which is the most common way a real analysis quietly goes wrong. The family extends by assay type: SingleCellExperiment adds reduced-dimensional embeddings and cell-level metadata, VariantAnnotation wraps VCF content with transcript models, RaggedExperiment handles per-sample variant collections. Learn the container once and every downstream tool speaks it.

Installation has one entry point, and using the wrong one is the classic beginner failure:

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("DESeq2")
library(DESeq2)

BiocManager::install resolves the package version compatible with your Bioconductor release. Installing Bioconductor packages with plain install.packages produces version mismatches that fail deep inside dependency chains, usually on someone else’s machine.

The canon: DESeq2, edgeR, limma-voom

For differential expression, three tools cover virtually every design, and the choice is driven by scale and design complexity rather than fashion:

The workflow has the same shape in all three: build a count matrix (from Salmon via tximport or from STAR featureCounts — the counting step belongs to the bulk RNA-seq pipeline), attach sample metadata, fit a model with an explicit design formula, shrink fold changes, and report FDR-adjusted results. Every one of these packages returns Benjamini-Hochberg adjusted p-values by default because testing twenty thousand genes at α = 0.05 guarantees a thousand false positives by construction — the discipline that the multiple-comparisons guide covers in full, and the reason an unadjusted “significant” list is not a result.

Everywhere else R is the default

The same pattern — canonical container, canonical tool — repeats across the field:

• Single-cell RNA-seq: SingleCellExperiment plus Seurat for clustering, dimensionality reduction and marker testing — the computational subject of the scRNA-seq guide.
• ChIP-seq: DiffBind aggregates peak sets across samples and re-tests differential binding with DESeq2’s statistics under the hood, closing the loop on the peak-calling workflow.
• Enrichment analysis: clusterProfiler performs GO and KEGG over-representation and GSEA with correct background-universe handling and ontology-aware redundancy reduction — it is the standard execution path for the enrichment workflow.
• Methylation: minfi and bumphunter process arrays and bisulfite data into the beta matrices and region tests that the methylation workflow describes.
• Variant annotation: the VariantAnnotation package reads VCFs, filters on INFO/FORMAT fields, and overlaps variants with transcript models — the R-side companion to the VCF pipeline.

In each case R performs the same role: hold a domain data model and apply statistically principled, citable tests to it.

Figures: ggplot2 and the grammar of graphics

Publication figures in biology are disproportionately R-made because ggplot2 encodes the grammar of graphics: a plot is data, plus aesthetics mapped to variables, plus geometric layers, plus scales. A volcano plot is not a special plot type — it is a scatterplot with log₂ fold change on x, −log₁₀ adjusted p on y, and colour mapped to a significance threshold:

library(ggplot2)
res <- as.data.frame(results(dds, alpha = 0.05))
res$gene <- rownames(res)
ggplot(res, aes(log2FoldChange, -log10(padj))) +
  geom_point(aes(color = padj < 0.05), size = 0.6, alpha = 0.6) +
  scale_color_manual(values = c("grey70", "firebrick")) +
  geom_vline(xintercept = c(-1, 1), linetype = 2) +
  geom_hline(yintercept = -log10(0.05), linetype = 2) +
  theme_bw()

Because the plot is a declarative specification, it recomposes when the data change, and the same grammar drives heatmaps, PCA plots and annotation layers with no new API to learn. How to read what those axes actually mean — why fold change alone misleads and adjusted p alone misleads — is the subject of the effect-size and p-value guide; the volcano plot is that guide made visible.

The tidyverse middle ground

Between base R and Bioconductor sits the tidyverse — dplyr for verbs (filter, group_by, summarise), tidyr for reshaping, purrr for iteration — and it is where count tables become clean inputs:

library(dplyr)
top <- as.data.frame(results(dds)) %>%
  tibble::rownames_to_column("gene") %>%
  filter(!is.na(padj), padj < 0.01, abs(log2FoldChange) > 1) %>%
  arrange(padj) %>%
  head(20)

Pipe-native code reads like the analysis it describes, and that legibility is a reproducibility feature: a collaborator can audit the logic line by line. Bioconductor objects and tidy verbs interoperate cleanly — fit with DESeq2, convert with as.data.frame, filter with dplyr, present with ggplot2.

Where R stops and Python takes over

An honest boundary saves you from learning the wrong half of each language deeply:

• Raw sequencing files. R is not where you parse FASTQ at scale or align reads — that is Biopython, samtools and fastp territory (see the Biopython tutorial). Import counts or matrices into R; do not build them there.
• Machine learning beyond classical statistics. Deep learning on sequences, images or protein structures is Python’s ecosystem entirely — the machine-learning guide covers what actually applies in bioinformatics. R’s ML packages exist, but the frameworks, the frontier and the hiring market are Python.
• Pipeline orchestration. Production pipelines wrap R scripts as processes inside Nextflow or Snakemake. R is a stage, not the conductor.
• Long-lived software. R is a language for analyses; large engineering codebases are better served elsewhere.

The practical synthesis for a beginner: Python first for programming fundamentals and file handling (start here), then R for statistics and figures — or, if your entry point is expression analysis, R first with Python picked up alongside. Both, eventually; neither alone.

A first session: counts to results in ten lines

The shape of a complete DESeq2 session, from a counts matrix and a sample table:

library(DESeq2)
cts <- read.csv("counts.csv", row.names = 1)
col <- read.csv("samples.csv", row.names = 1)   # condition column: treat / control
dds <- DESeqDataSetFromMatrix(cts, col, ~ condition)
dds <- DESeq(dds)
res <- results(dds, contrast = c("condition", "treat", "control"))
res <- lfcShrink(dds, contrast = c("condition", "treat", "control"), type = "apeglm")
summary(res)          # how many genes at FDR < 0.05, up and down
write.csv(as.data.frame(res), "deseq2_results.csv")

Read that sequence once — construct, fit, test, shrink, report — and you have the skeleton of most expression papers. Everything beyond it is design-matrix care and interpretation discipline.

Key takeaways

• R is bioinformatics’ statistics layer; Bioconductor’s release discipline, shared containers and citable tools made it the field standard.
• Data travel in SummarizedExperiment-family containers, keeping counts, sample metadata and feature annotation aligned through every subset.
• DESeq2, edgeR and limma-voom cover essentially every differential-expression design — choose by scale and complexity, and always report FDR-adjusted results.
• Seurat, DiffBind, clusterProfiler, minfi and VariantAnnotation repeat the container-plus-canon pattern across single-cell, ChIP-seq, enrichment, methylation and variants.
• ggplot2’s grammar makes publication figures declarative; a volcano plot is a mapped scatterplot, not a special plot.
• Python owns raw-file wrangling and modern ML, R owns testing and figures — real pipelines use both, orchestrated by a workflow engine.

Challenge

Take a public count matrix from GEO’s supplementary files — the download guide shows exactly where they live — and run the ten-line DESeq2 session above end to end. Then answer three questions from your own output: how many genes pass FDR 0.05, which gene has the largest shrunken fold change, and how does the volcano plot change when you plot unshrunken instead of shrunken log₂ fold changes? That third question teaches more about shrinkage than any documentation page.

Want guided practice turning expression matrices into defensible statistical results? Ampersand Academy offers one-to-one mentoring in bioinformatics — part of Mahadhi, a development and digital marketing company.

FAQ