Awk, Sed, and Bash One-Liners Every Bioinformatician Should Know

Every analysis eventually stalls on the same unglamorous problem: the file is almost in the right shape, but not quite. The chromosome column has a chr prefix your reference lacks, the FASTA headers carry junk you need gone, or you just want the mean read length before committing to a two-hour alignment. Firing up Python for these tasks is overkill. The Unix text-processing trio, awk, sed, and grep, glued together with bash, will do each of them in a single line that runs in seconds on a multi-gigabyte file. This post teaches the patterns, not just the incantations, so you can adapt them instead of pasting them blindly.

Why These Tools Beat a Script Here

Genomic text formats, FASTA, FASTQ, VCF, BED, and SAM, are almost all line-oriented and tab- or space-delimited. That is exactly the data model awk was designed for in 1977, and nothing has dethroned it. awk reads input line by line, splits each line into fields ($1, $2, and so on) on whitespace by default, and runs your code on every line. sed is a stream editor built for find-and-replace and line surgery. grep filters lines by pattern. Because all three stream their input, they use near-constant memory regardless of file size, so a 10 GB FASTQ costs the same RAM as a 10 KB one.

The key mental model for awk is the pattern { action } rule: if the pattern is true for a line, the action runs. Omit the pattern and the action runs on every line; omit the action and matching lines are printed. Everything below is a variation on that one idea.

Reading Field-Delimited Files Correctly

The single most common mistake with genomic data is letting awk split on any whitespace when your file is strictly tab-delimited. VCF, BED, SAM, and GTF all use tabs, and some fields legitimately contain spaces (a sample description, a gene name). Always set the field separator explicitly:

# -F'\t' sets the INPUT field separator to a literal tab
awk -F'\t' '$1 == "chr1"' variants.vcf

To make the output tab-delimited too, set OFS (output field separator) in a BEGIN block, which runs once before any line is read:

awk 'BEGIN{FS=OFS="\t"} {print $1, $2, $3}' regions.bed

FASTA: Extracting and Reformatting Sequences

FASTA files alternate header lines (starting with >) and sequence lines. To count sequences, you do not need awk at all, just count headers with grep -c, where -c reports a count and ^> anchors the match to the start of the line:

grep -c "^>" genome.fa

To pull one record by name, sed can print a range, but a cleaner approach uses awk with a flag variable. The pattern below prints the matching header and every sequence line until the next header:

awk '/^>chr2$/{f=1; print; next} /^>/{f=0} f' genome.fa

Here /^>chr2$/ matches the exact header, sets flag f=1, prints it, and next skips to the next line. The rule /^>/{f=0} turns the flag off at the next header. The bare f at the end is a pattern with no action: when f is nonzero (true), the line prints. To linearize a multi-line FASTA (one sequence per line, which makes downstream awk trivial):

awk '/^>/{if(seq)print seq; print; seq=""; next} {seq=seq$0} END{if(seq)print seq}' genome.fa

FASTQ: Read Counts and Length Statistics

FASTQ is rigidly four lines per read: header, sequence, +, quality. That regularity makes NR, the built-in line-number variable, your best friend. The sequence sits on every line where NR % 4 == 2. To count reads, count lines and divide by four:

echo $(( $(wc -l < reads.fastq) / 4 ))

For read-length statistics, operate only on sequence lines and accumulate. length($0) gives the character count of the current line:

awk 'NR%4==2 {sum+=length($0); n++; if(length($0)>max)max=length($0)}
     END {printf "reads=%d mean=%.1f max=%d\n", n, sum/n, max}' reads.fastq

The END block runs once after the last line, which is where you emit summaries. For gzipped input, do not gunzip to disk; stream it with zcat (or gzcat on macOS) piped into awk. If any of this four-line logic feels unfamiliar, the FASTA/FASTQ formats explainer walks through why the structure is rigid in the first place.

VCF and BED: Filtering and Coordinate Surgery

VCF files carry ## metadata headers, a #CHROM column header, then one variant per line. To filter without destroying the header, let any line starting with # pass through, then apply your condition. The QUAL score lives in column 6:

awk -F'\t' '/^#/ || $6 >= 30' variants.vcf > filtered.vcf

Coordinate conversion is where these tools shine. BED is 0-based half-open; many other formats are 1-based inclusive. To convert a 1-based, 2-column position list into valid BED intervals, subtract one from the start:

awk 'BEGIN{OFS="\t"} {print $1, $2-1, $2}' positions.txt > sites.bed

The perennial reference-mismatch headache, chr1 versus 1, is a one-character sed job. Anchor with ^ so you only touch the start of the line, never a chr buried inside an INFO field:

# Strip the chr prefix
sed 's/^chr//' with_chr.vcf > no_chr.vcf
# Add it back
awk -F'\t' 'BEGIN{OFS="\t"} /^#/{print;next} {$1="chr"$1; print}' no_chr.vcf

For anything beyond these edits, reach for bedtools, which handles interval intersection and merging that pure awk should not attempt. To understand the columns you are slicing, keep the VCF format explainer open in another tab.

SAM/BAM: Filtering With samtools and awk

SAM is tab-delimited with @-prefixed headers. You almost never touch BAM directly with awk; instead you decode it with samtools view and pipe the SAM stream in. MAPQ (mapping quality) sits in column 5. To count primary alignments with MAPQ >= 20:

samtools view file.bam | awk -F'\t' '$5 >= 20' | wc -l

The FLAG field (column 2) is a bitmask, and this is the one place a naive awk comparison goes wrong: you cannot test bits with a plain numeric comparison. Let samtools view -f and -F do the bit filtering (for example -F 0x904 to drop secondary, supplementary, and unmapped reads), then use awk for the numeric columns it handles well. Mixing responsibilities like this, the right tool for each job in one pipeline, is the whole point.

bioawk: awk That Speaks Biology

Plain awk does not know what a FASTQ record is; you juggle NR % 4 yourself. bioawk, a fork of awk by Heng Li, adds format-aware parsing. With -c fastx it exposes named fields like $name, $seq, $qual, and $comment, and automatically handles both FASTA and FASTQ, gzipped or not, treating each record as one logical unit no matter how many physical lines it spans:

# GC content per sequence, no manual line counting
bioawk -c fastx '{gc=gsub(/[GgCc]/,"",$seq); print $name, gc/length($seq)}' genome.fa

Here is how the tools divide the labor:

bioawk ships on Bioconda, so conda install -c bioconda bioawk into a dedicated environment gets you a working binary without touching your system awk.

Composing One-Liners Into Pipelines

The real power comes from chaining with the pipe (|), where each tool's output feeds the next. To find the ten longest reads in a gzipped FASTQ, prepend the length, sort numerically descending, and cut it back off:

zcat reads.fastq.gz \
  | awk 'NR%4==2 {print length($0)"\t"$0}' \
  | sort -k1,1nr \
  | head

sort -k1,1nr sorts on field 1 (-k1,1), numerically (n), reversed (r). Two habits keep pipelines robust: put set -euo pipefail at the top of any bash script so a failure anywhere aborts the whole chain instead of silently producing a truncated file, and quote your variables ("$file") so filenames with spaces do not shatter into multiple arguments.

Practical Takeaways

• Always set -F'\t' for VCF, BED, SAM, and GTF; whitespace-splitting silently corrupts fields that contain spaces.
• Use NR % 4 == 2 to grab FASTQ sequence lines, and flag variables (f=1 ... f) to extract FASTA records by header.
• Anchor sed substitutions with ^ when editing chromosome prefixes so you never touch text mid-line, and remember awk needs a field reassignment before OFS takes effect.
• Let samtools and bedtools own bitmask FLAG filtering and interval math; awk handles the plain numeric columns.
• Reach for bioawk -c fastx whenever per-record sequence logic gets tangled in manual line counting.
• Chain tools with | and guard scripts with set -euo pipefail and quoted variables.

These patterns will carry most day-to-day wrangling, but they are a scalpel, not a Swiss Army knife: for heavy interval math and BAM manipulation, the dedicated tools above will always be more correct than a hand-rolled awk reimplementation. Learn the one-liners for the quick edits, and hand the structured operations to purpose-built software.