5: Biological File I/O & Bioconductor Formats

Chapter 5: Biological File I/O & Bioconductor Formats

Biological discovery relies on standardized file formats: FASTA (sequences), FASTQ (reads + quality), BED (intervals), GTF/GFF3 (gene models), and VCF (genomic variants). Reading and writing these efficiently prevents memory bottlenecks.

1. Parsing Multi-Record FASTA Files

Reading FASTA sequences in base R with header-sequence mapping:

read_fasta <- function(filepath) {
    lines <- readLines(filepath)
    headers <- c()
    sequences <- c()
    cur_seq <- ""
    for (line in lines) {
        line <- trimws(line)
        if (startsWith(line, ">")) {
            if (cur_seq != "") { sequences <- c(sequences, cur_seq); cur_seq <- "" }
            headers <- c(headers, substring(line, 2))
        } else {
            cur_seq <- paste0(cur_seq, line)
        }
    }
    if (cur_seq != "") sequences <- c(sequences, cur_seq)
    data.frame(Header = headers, Sequence = sequences, Length = nchar(sequences), stringsAsFactors = FALSE)
}

2. Parsing BED Interval Coordinates

Standard BED files use 0-based start and 1-based end coordinates:

bed_df <- read.table("exons.bed", sep = "\t", header = FALSE,
                     col.names = c("chrom", "chromStart", "chromEnd", "name", "score", "strand"))

# Calculate feature widths in bp
bed_df$width <- bed_df$chromEnd - bed_df$chromStart
print(head(bed_df))

3. Saving & Restoring Serialized R Objects (saveRDS)

Saving preprocessed genomic objects in compressed binary format preserves factors, S4 classes, and metadata exactly:

# Save compressed RDS
saveRDS(bundle, "genomic_bundle.rds", compress = "xz")

# Restore instantly in any future session
bundle <- readRDS("genomic_bundle.rds")