# Foundations R Fundamentals

> Read/write R syntax for bioinformatics (vectors, data.frame, matrices, d/p/q/r distributions, DESeq2). Use when porting Python to R, debugging R from a paper/pipeline, or running DESeq2/edgeR/Seurat scripts.

- Skill: `pavel-kravchenko/foundations-r-fundamentals` (Agent Skill)
- Install (CLI): `npx skillmds@latest add pavel-kravchenko/foundations-r-fundamentals`
- Raw SKILL.md: https://api.skillmd.com/api/skills/pavel-kravchenko/foundations-r-fundamentals/raw
- Safety review: pending
- Works with: Claude Code, Claude.ai, OpenAI Codex
- Category: DevOps & Infra
- Author: pavel-kravchenko (https://skillmd.com/u/pavel-kravchenko)
- Updated: 2026-09-17
- Page: https://skillmd.com/skills/pavel-kravchenko/foundations-r-fundamentals

---


# R Fundamentals for Bioinformatics

## When to Use

- Reading or adapting R code from a published pipeline, DESeq2/edgeR/Seurat vignette, or supplementary script
- Translating a Python data-analysis snippet into R (or vice versa) and need the syntax mapping
- Debugging an off-by-one or "wrong subset" bug caused by R's 1-based, exclude-on-negative indexing
- Setting up a minimal Bioconductor workflow (install, DESeqDataSet, filtering results)
- Choosing/running a basic statistical test in R (t-test, Wilcoxon, chi-squared, correlation)

## Version Compatibility

R ≥ 4.2, Bioconductor ≥ 3.17 (BiocManager ≥ 1.30), DESeq2 ≥ 1.40, dplyr ≥ 1.1, ggplot2 ≥ 3.4. Base-R syntax here (`<-`, `data.frame`, `apply`) is stable across R versions back to R 3.x.

## Prerequisites

- R installed (`R --version`), or an R kernel (IRkernel) for Jupyter
- `install.packages("BiocManager")` for any Bioconductor package (DESeq2, GenomicRanges, ...)
- Comfort with Python's data types helps — this skill is written as a Python→R translation

## Python vs R Syntax Traps

**Goal:** avoid the bugs that come from assuming R behaves like Python.
**Approach:** memorize this table before writing/reading any R code — every row is a real, common bug source when porting between the two languages.

| Feature | Python | R |
|---------|--------|---|
| Indexing | Starts at 0 | Starts at **1** |
| Assignment | `=` | `<-` (preferred; `=` works but differs inside function calls) |
| Negative index | `x[-1]` = last element | `x[-1]` = **all except first** |
| Boolean values | `True` / `False` | `TRUE` / `FALSE` |
| Missing data | `None` | `NA` |
| Null | `None` | `NULL` |
| Auto-print | `print(x)` | Just type `x` (top level only) |
| AND/OR | `&`, `\|` (vectorized) | `&`, `\|` (vectorized); `&&`, `\|\|` (scalar, first element only) |

## Vectors and Data Frames

**Goal:** manipulate vectors and data frames the R way.
**Approach:** use vectorized operations and logical indexing instead of loops; use `data.frame` for tabular data and filter with boolean masks on columns.

```r
# --- Vectors ---
gene_expression <- c(2.5, 3.1, 4.2, 1.8, 5.6)
gene_names <- c("BRCA1", "TP53", "EGFR", "MYC", "KRAS")
positions <- 1:10
coverage <- seq(from = 0, to = 100, by = 10)
groups <- rep(c("control", "treatment"), each = 3)

# Named vector + logical indexing
gc_content <- c(BRCA1 = 0.42, TP53 = 0.38, EGFR = 0.55, MYC = 0.61)
gc_rich <- gc_content[gc_content > 0.5]

# Vectorized ops -- no explicit loops needed
log2fc <- log2(gene_expression + 1)          # pseudocount avoids log(0)
pvalue <- c(0.01, 0.15, 0.001, 0.02, 0.4)
sig_up <- log2fc > 1 & pvalue < 0.05          # vectorized filter
cat("Upregulated:", gene_names[sig_up], "\n")
which(sig_up)                                  # indices of TRUE values

# --- Data frames ---
gene_data <- data.frame(
    gene = c("BRCA1", "TP53", "EGFR"),
    log2fc = c(1.2, -0.5, 3.8),
    pvalue = c(0.01, 0.15, 0.001),
    stringsAsFactors = FALSE   # keep strings as character, not factor
)

gene_data$gene                  # column access
gene_data[1, ]                  # first row
gene_data[2, "log2fc"]          # specific cell
gene_data[gene_data$log2fc > 1 & gene_data$pvalue < 0.05, ]  # filter rows

gene_data$padj <- p.adjust(gene_data$pvalue, method = "BH")   # FDR (standard for genomics)
gene_data[order(gene_data$pvalue), ]                            # sort by p-value
```

## Matrices and Statistical Distributions

**Goal:** work with expression matrices and distribution functions.
**Approach:** use `matrix()` with `dimnames` for gene-by-sample data and `apply()` for row/column stats; use the `d/p/q/r` prefix family for any distribution.

```r
# --- Matrices (expression matrices) ---
expr_matrix <- matrix(
    c(5.2, 3.1, 8.5, 6.2, 4.8, 2.9, 7.1, 5.8),
    nrow = 2, byrow = TRUE,
    dimnames = list(c("Sample1", "Sample2"), c("BRCA1", "TP53", "EGFR", "MYC"))
)

expr_matrix[2, "EGFR"]          # specific cell
expr_matrix[1, ]                # full row
expr_matrix[, "BRCA1"]          # full column
apply(expr_matrix, 1, mean)     # row means (margin = 1)
apply(expr_matrix, 2, sd)       # column SDs (margin = 2)

# --- Statistical distributions: d(ensity) p(CDF) q(uantile) r(andom) ---
set.seed(42)                     # reproducibility
rnorm(n = 100, mean = 10, sd = 2)    # random samples
pnorm(q = 1.96, lower.tail = TRUE)   # CDF -> 0.975
qnorm(p = 0.975)                     # quantile -> 1.96
rnbinom(n = 100, mu = 50, size = 5)  # negative binomial -- what DESeq2 models RNA-seq counts with

# t-test / Wilcoxon / chi-squared -- the basic hypothesis-test toolkit
control <- rnorm(30, mean = 10, sd = 2)
treatment <- rnorm(30, mean = 12, sd = 2)
t.test(control, treatment)$p.value
wilcox.test(control, treatment)$p.value      # non-parametric alternative
chisq.test(table(c("Mut","WT","Mut","Mut"), c("Ctrl","Ctrl","Drug","Drug")))
```

## Minimal DESeq2 Workflow

**Goal:** run a minimal DESeq2 differential-expression analysis.
**Approach:** build a `DESeqDataSet` from a raw count matrix + sample metadata, run `DESeq()`, extract a contrast, then filter by FDR and effect size.

```r
# Install once
if (!requireNamespace("BiocManager", quietly = TRUE)) install.packages("BiocManager")
BiocManager::install("DESeq2")

library(DESeq2)

#' Run DESeq2 and return significant genes for a two-group contrast.
#' @param counts   integer matrix, genes x samples (raw counts, NOT normalized)
#' @param metadata data.frame, one row per sample, must include `condition`
#' @param alpha    FDR threshold (default 0.05)
#' @param lfc      minimum absolute log2 fold-change to call significant
#' @return data.frame of significant genes sorted by adjusted p-value
run_deseq2 <- function(counts, metadata, alpha = 0.05, lfc = 1) {
    dds <- DESeqDataSetFromMatrix(countData = counts, colData = metadata, design = ~condition)
    dds <- DESeq(dds)
    res <- results(dds, contrast = c("condition", "treatment", "control"))
    res_df <- as.data.frame(res)
    sig <- res_df[!is.na(res_df$padj) & res_df$padj < alpha & abs(res_df$log2FoldChange) > lfc, ]
    sig[order(sig$padj), ]
}

# Reading/writing tabular data
df <- read.csv("data.csv", stringsAsFactors = FALSE)
df <- read.delim("data.tsv", sep = "\t")
write.csv(df, "output.csv", row.names = FALSE)
```

```r
# Self-check: run_deseq2() end to end on synthetic count data
demo_deseq2 <- function() {
    set.seed(1)
    counts <- matrix(rnbinom(400, mu = 100, size = 5), nrow = 100, ncol = 4,
                      dimnames = list(paste0("gene", 1:100), paste0("s", 1:4)))
    metadata <- data.frame(condition = factor(c("control", "control", "treatment", "treatment")),
                            row.names = colnames(counts))
    sig <- run_deseq2(counts, metadata)
    stopifnot(is.data.frame(sig))
    stopifnot(all(c("log2FoldChange", "padj") %in% colnames(sig)))
    cat("OK:", nrow(sig), "significant genes out of", nrow(counts), "\n")
}
demo_deseq2()
```

## Pitfalls

- **R indexing starts at 1**: `x[1]` is first; `x[-1]` means "all except first" (NOT last element like Python)
- **Factors look like strings but are not**: a factor stores integer codes + labels; functions expecting strings will misbehave — use `stringsAsFactors = FALSE` in `data.frame()` or `as.character()` to convert
- **`library()` vs `require()`**: use `library()` in scripts — it errors loudly if a package is missing; `require()` returns `FALSE` silently and can hide broken pipelines
- **Data frames vs matrices**: Bioconductor functions often require numeric matrices; a data frame with character columns cannot be used directly — convert with `as.matrix(df[, numeric_cols])`
- **`<-` inside function arguments**: `f(x <- 1)` assigns `x` in the enclosing scope AND passes the value as the argument; use `=` for keyword arguments to avoid the side effect
- **Multiple testing**: `p.adjust(..., method = "BH")` (Benjamini-Hochberg FDR) is the genomics standard; Bonferroni is overly conservative for large gene sets
- **Scalar vs vectorized AND/OR**: `&&` and `||` only evaluate the first element — always use `&` and `|` when filtering vectors/data frames

## See Also

- `foundations-r-hypothesis-testing-and-nonparametrics` — Shapiro-Wilk, Wilcoxon, ANOVA, correlation tests in depth
- `foundations-r-regression-correlation-and-diagnostics` — linear models, `lm()`, diagnostic plots
- `foundations-biostatistics-fundamentals` — statistical concepts underlying these tests
- `bio-differential-expression-deseq2-basics` — full DESeq2 workflow beyond this quick-start

