Project: Wrangle a Real Data File

Bioinformaticians spend a surprising amount of time just interrogating files from the command line before any heavy analysis - how big is it, how many records, does it contain what I expect? That's exactly this project. You'll grab a real FASTA file of three human genes and answer real questions about it.

1Get the data

Make a workspace and download three real human gene sequences (insulin, beta-globin, and APP) straight from NCBI:

mkdir fasta-project
cd fasta-project

curl "https://eutils.ncbi.nlm.nih.gov/entrez/eutils/efetch.fcgi?db=nuccore&id=NM_000207,NM_000518,NM_000484&rettype=fasta&retmode=text" -o genes.fasta

(If curl isn't found, try wget with the same URL and -O genes.fasta.)

No internet, or NCBI blocked? Create a small real FASTA by hand instead - paste this whole block into your terminal and press Enter:

cat > genes.fasta <<'EOF'
>gene1 insulin
ATGGCCCTGTGGATGCGCCTCCTGCCCCTGCTGGCGCTGCTG
>gene2 hemoglobin_beta
ATGGTGCATCTGACTCCTGAGGAGAAGTCTGCCGTTACTGCC
>gene3 albumin
ATGAAGTGGGTAACCTTTATTTCCCTTCTTTTTCTCTTTAGC
EOF

Every command below works the same on this file.

2First look

head genes.fasta     # peek at the top - don't dump the whole thing
wc -l genes.fasta    # how many lines total?

3How many sequences?

Every FASTA record starts with a > header line. So counting headers counts sequences:

grep -c ">" genes.fasta

That's the single most common FASTA sanity-check in all of bioinformatics - and you can now do it instantly.

4What are they?

grep ">" genes.fasta            # list just the header lines
grep ">" genes.fasta > headers.txt  # save them to a file

5How long is the sequence, in total?

Here's a pipeline: strip out the header lines, remove the line breaks, and count the remaining characters - the total number of bases.

grep -v ">" genes.fasta | tr -d '\n' | wc -c

Reading the pipeline

grep -v ">" keeps every line that is not a header (the -v inverts the match). tr -d '\n' deletes the newlines so the bases form one continuous string. wc -c counts the characters. Three tiny tools, one real answer.

6GC content - from the command line

Count the G and C bases, then the total, then compute the percentage:

# Count just G and C characters in the sequence
grep -v ">" genes.fasta | tr -d '\n' | tr -cd 'GC' | wc -c

# Do the whole calculation in one go
gc=$(grep -v ">" genes.fasta | tr -d '\n' | tr -cd 'GC' | wc -c)
total=$(grep -v ">" genes.fasta | tr -d '\n' | wc -c)
awk "BEGIN{printf \"GC%% = %.1f\n\", $gc/$total*100}"

You just computed a real biological metric on real data with a few command-line tools - no programming language required. That's the power (and the surprise) of Bash for bioinformatics.

🚀 Make it your own

Download a different set of genes by editing the id= list in the URL (any NCBI nucleotide accessions).
Use grep -c "N" to count ambiguous bases (Ns) - a quick quality check.
Count how many headers mention a keyword: grep ">" genes.fasta | grep -c "insulin".
Save your whole session as a shell script (analyze.sh) and run it with bash analyze.sh - your first automation.
Portfolio move: put the script and a short README on GitHub.

What you just did

You downloaded real biological data and answered concrete questions about it - count, identity, length, composition - using only the command line. This is not a toy version of bioinformatics; it's the real daily reality of inspecting files before and after every analysis. The aligners and variant callers in Track 2 are just bigger commands that fit exactly this pattern.

Why tr -cd 'GC' and not tr -d?

tr -cd 'GC' means "delete every character that is not G or C" (the -c complements the set). What's left is only the G and C characters, so wc -c counts exactly those.

Check your understanding

How do you count the sequences in a FASTA file from the command line?

Tools like grep, wc and sort let you:

In the long-sequence pipeline, what does grep -v ">" do?

The -v flag inverts the match, so grep -v keeps every line that is not a header, leaving just the sequence lines.

Why does the GC-content step use tr -cd 'GC' rather than tr -d 'GC'?

The -c complements the set, so tr -cd 'GC' deletes every character that is not G or C; what remains is only the G and C bases for wc -c to count.

The project shows that, before heavy analysis, bioinformaticians spend a lot of time doing what?

Much real work is just inspecting files from the command line (size, record count, contents) before any heavy analysis, which is exactly what this project practices.

You've finished the foundations 🎉