LLMs-Universal-Life-Science-and-Clinical-Skills- fastq-quality

<!--

install
source · Clone the upstream repo
git clone https://github.com/mdbabumiamssm/LLMs-Universal-Life-Science-and-Clinical-Skills-
Claude Code · Install into ~/.claude/skills/
T=$(mktemp -d) && git clone --depth=1 https://github.com/mdbabumiamssm/LLMs-Universal-Life-Science-and-Clinical-Skills- "$T" && mkdir -p ~/.claude/skills && cp -r "$T/Skills/Sequence_Analysis/sequence-io/fastq-quality" ~/.claude/skills/mdbabumiamssm-llms-universal-life-science-and-clinical-skills-fastq-quality && rm -rf "$T"
manifest: Skills/Sequence_Analysis/sequence-io/fastq-quality/SKILL.md
source content
<!-- # COPYRIGHT NOTICE # This file is part of the "Universal Biomedical Skills" project. # Copyright (c) 2026 MD BABU MIA, PhD <md.babu.mia@mssm.edu> # All Rights Reserved. # # This code is proprietary and confidential. # Unauthorized copying of this file, via any medium is strictly prohibited. # # Provenance: Authenticated by MD BABU MIA -->

name: bio-fastq-quality description: Work with FASTQ quality scores using Biopython. Use when analyzing read quality, filtering by quality, trimming low-quality bases, or generating quality reports. tool_type: python primary_tool: Bio.SeqIO measurable_outcome: Execute skill workflow successfully with valid output within 15 minutes. allowed-tools:

  • read_file
  • run_shell_command

FASTQ Quality Scores

Analyze and manipulate FASTQ quality scores using Biopython.

Required Imports

from Bio import SeqIO
from Bio.Seq import Seq

Accessing Quality Scores

Quality scores are stored in

letter_annotations['phred_quality']
as a list of integers.

for record in SeqIO.parse('reads.fastq', 'fastq'):
    qualities = record.letter_annotations['phred_quality']
    print(record.id, qualities[:10])  # First 10 quality scores

Quality Score Basics

Phred ScoreError ProbabilityAccuracy
101 in 1090%
201 in 10099%
301 in 100099.9%
401 in 1000099.99%

Code Patterns

Calculate Average Quality per Read

for record in SeqIO.parse('reads.fastq', 'fastq'):
    quals = record.letter_annotations['phred_quality']
    avg_qual = sum(quals) / len(quals)
    print(f'{record.id}: {avg_qual:.1f}')

Filter Reads by Mean Quality

def high_quality_reads(records, min_avg_qual=20):
    for record in records:
        quals = record.letter_annotations['phred_quality']
        if sum(quals) / len(quals) >= min_avg_qual:
            yield record

records = SeqIO.parse('reads.fastq', 'fastq')
good_reads = high_quality_reads(records, min_avg_qual=25)
SeqIO.write(good_reads, 'filtered.fastq', 'fastq')

Filter by Minimum Quality at Any Position

def all_bases_above(records, min_qual=20):
    for record in records:
        if min(record.letter_annotations['phred_quality']) >= min_qual:
            yield record

Trim Low-Quality Ends (3' Trimming)

def trim_low_quality(record, min_qual=20):
    quals = record.letter_annotations['phred_quality']
    trim_pos = len(quals)
    for i in range(len(quals) - 1, -1, -1):
        if quals[i] >= min_qual:
            trim_pos = i + 1
            break
    return record[:trim_pos]

records = SeqIO.parse('reads.fastq', 'fastq')
trimmed = (trim_low_quality(rec) for rec in records)
SeqIO.write(trimmed, 'trimmed.fastq', 'fastq')

Sliding Window Quality Trim

def sliding_window_trim(record, window_size=5, min_avg_qual=20):
    quals = record.letter_annotations['phred_quality']
    for i in range(len(quals) - window_size + 1):
        window = quals[i:i + window_size]
        if sum(window) / window_size < min_avg_qual:
            return record[:i] if i > 0 else None
    return record

Quality Statistics Summary

import statistics

all_quals = []
for record in SeqIO.parse('reads.fastq', 'fastq'):
    all_quals.extend(record.letter_annotations['phred_quality'])

print(f'Mean quality: {statistics.mean(all_quals):.1f}')
print(f'Median quality: {statistics.median(all_quals):.1f}')
print(f'Min quality: {min(all_quals)}')
print(f'Max quality: {max(all_quals)}')

Per-Position Quality Profile

from collections import defaultdict

position_quals = defaultdict(list)
for record in SeqIO.parse('reads.fastq', 'fastq'):
    for i, q in enumerate(record.letter_annotations['phred_quality']):
        position_quals[i].append(q)

for pos in sorted(position_quals.keys())[:20]:
    quals = position_quals[pos]
    print(f'Position {pos}: mean={sum(quals)/len(quals):.1f}')

Count Reads by Quality Threshold

thresholds = [20, 25, 30, 35]
counts = {t: 0 for t in thresholds}

for record in SeqIO.parse('reads.fastq', 'fastq'):
    avg = sum(record.letter_annotations['phred_quality']) / len(record.seq)
    for t in thresholds:
        if avg >= t:
            counts[t] += 1

for t, c in counts.items():
    print(f'Q>={t}: {c} reads')

Remove N Bases and Low Quality Together

def clean_read(record, min_qual=20):
    seq = str(record.seq)
    quals = record.letter_annotations['phred_quality']
    keep = [(s, q) for s, q in zip(seq, quals) if s != 'N' and q >= min_qual]
    if not keep:
        return None
    new_seq, new_quals = zip(*keep)
    new_record = record[:0]  # Empty copy with same metadata
    new_record.seq = Seq(''.join(new_seq))
    new_record.letter_annotations['phred_quality'] = list(new_quals)
    return new_record

FASTQ Format Variants

VariantFormat StringQuality EncodingASCII Range
Sanger/Illumina 1.8+
'fastq'
Phred+33 (standard)33-126
Solexa
'fastq-solexa'
Solexa+6459-126
Illumina 1.3-1.7
'fastq-illumina'
Phred+6464-126

Most modern data uses standard

'fastq'
(Sanger/Illumina 1.8+).

Quality Score Conversion

For legacy data using different quality encodings:

from Bio.SeqIO.QualityIO import phred_quality_from_solexa, solexa_quality_from_phred

Convert Solexa to Phred

from Bio.SeqIO.QualityIO import phred_quality_from_solexa

# Convert single score
solexa_score = 10
phred_score = phred_quality_from_solexa(solexa_score)

# Convert list of scores
solexa_scores = [10, 20, 30, 40]
phred_scores = [phred_quality_from_solexa(s) for s in solexa_scores]

Convert Phred to Solexa

from Bio.SeqIO.QualityIO import solexa_quality_from_phred

phred_score = 30
solexa_score = solexa_quality_from_phred(phred_score)

Convert Between FASTQ Variants

from Bio import SeqIO

# Read old Illumina format, write standard format
records = SeqIO.parse('old_reads.fastq', 'fastq-illumina')
SeqIO.write(records, 'standard_reads.fastq', 'fastq')

# Read Solexa format, write standard format
records = SeqIO.parse('solexa_reads.fastq', 'fastq-solexa')
SeqIO.write(records, 'standard_reads.fastq', 'fastq')

Auto-Detect Quality Encoding

def detect_quality_encoding(filepath, sample_size=1000):
    '''Guess FASTQ quality encoding from ASCII values'''
    min_qual = 126
    max_qual = 0
    count = 0

    with open(filepath) as f:
        for i, line in enumerate(f):
            if i % 4 == 3:  # Quality line
                for char in line.strip():
                    min_qual = min(min_qual, ord(char))
                    max_qual = max(max_qual, ord(char))
                count += 1
                if count >= sample_size:
                    break

    if min_qual < 59:
        return 'fastq'  # Sanger/Illumina 1.8+ (Phred+33)
    elif min_qual < 64:
        return 'fastq-solexa'  # Solexa+64
    else:
        return 'fastq-illumina'  # Illumina 1.3+ (Phred+64)

Common Errors

ErrorCauseSolution
KeyError: 'phred_quality'
Not FASTQ or wrong variantCheck format, try 'fastq-illumina'
Quality scores all 0Wrong encoding assumedTry different FASTQ variant
Trimmed reads emptyToo aggressive trimmingLower quality threshold

Related Skills

  • read-sequences - Parse FASTQ files
  • filter-sequences - Filter reads by other criteria (length, content)
  • paired-end-fastq - Handle R1/R2 paired quality filtering
  • sequence-statistics - Generate summary statistics including quality
  • alignment-files - After filtering, align reads with bwa/bowtie2; quality scores in BAM
<!-- AUTHOR_SIGNATURE: 9a7f3c2e-MD-BABU-MIA-2026-MSSM-SECURE -->