pyfastx

repository·master·Indexed 43 days ago

https://github.com/lmdu/pyfastx

A lightweight Python C extension for fast, random access to sequences in plain or gzipped FASTA and FASTQ files. It utilizes a sqlite3-based indexing system to provide efficient access with low memory consumption. Key features include the pyfastx.Fasta and pyfastx.Fastq classes for parsing and indexing, the pyfastx.Fastx class for high-performance iteration, and support for parallel random access via Python's multiprocessing module.

Tokens
12.7K
Snippets
59
Records
63
Agent score
52%

What's inside pyfastx

  1. Overview of pyfastx features and capabilities

    master

    pyfastx is a lightweight Python C extension designed for efficient, random access to sequences in both plain and gzipped FASTA/Q files. It uses a sqlite3 database to store indexes, allowing for fast random access without high memory consumption.

    Key Features:

    • Random Access: Quickly extract sequences from FASTA or reads from FASTQ files by identifier or index number.
    • Gzip Support: Fast random access even within compressed gzipped files.
    • Memory Efficiency: Lightweight parsing suitable for large genomic datasets.
    • Flexible Parsing: Supports both standard (fixed-length lines) and nonstandard (variable-length lines) FASTA formats.
    • Sequence Analysis:
      • Calculate assembly N50 and L50.
      • Calculate GC content and nucleotide composition.
      • Extract reverse, complement, and antisense sequences.
    • Iteration: Read FASTA sequences line by line.
  2. Manage FASTA keys with FastaKeys

    master

    The fa.keys() method returns a FastaKeys object, which is a list-like object containing all sequence names in a FASTA file. It supports advanced filtering and sorting.

    Basic Operations

    • len(ids): Get count of sequences.
    • ids[index]: Get key by index.
    • 'name' in ids: Check if a name exists.
    • list(ids): Convert to a standard Python list.

    Sorting

    You can sort the keys using .sort(by='...', reverse=True). Valid by values are:

    • 'id'
    • 'name'
    • 'length'

    Filtering

    You can filter keys using the .filter() method with various expressions:

    • By length: ids.filter(ids > 600) or ids.filter(500 < ids < 700)
    • By name (substring): ids.filter(ids % 'substring')
    • Combined: ids.filter(ids % 'JZ8226', ids > 550)

    Resetting

    To clear any applied filters or sort orders, use ids.reset().

    import pyfastx
    fa = pyfastx.Fasta('test.fa')
    ids = fa.keys()
    
    # Sort by length descending
    for name in ids.sort(by='length', reverse=True):
        print(name)
    
    # Filter by name and length, then sort
    filtered = ids.filter(ids % 'JZ8226', ids > 730).sort('length', reverse=True)
    print(list(filtered))
    
    # Reset
    ids.reset()
  3. Optimize sequence name lookups in Python

    master

    When using pyfastx.Fasta in Python, checking if a name exists using the in operator (e.g., if seqname in fa:) can be slow because it performs a database query for every check.

    To significantly speed up lookups, load all sequence names into a Python set first:

    import pyfastx
    fa = pyfastx.Fasta('tests/data/test.fa.gz')
    # Faster alternative: load all names into a set
    all_names = set(fa.keys())
    for seqname in seqnames:
        if seqname in all_names:
            do_something()
    fa = pyfastx.Fasta('tests/data/test.fa.gz')
    all_names = set(fa.keys())
    for seqname in seqnames:
        if seqname in all_names:
            do_something()
  4. Prepare test data for benchmarking

    master

    To perform benchmarks, you need FASTA and FASTQ files.

    1. FASTA files: Download genome FASTA files (e.g., from NCBI FTP). Place them in a directory (e.g., data/fastas), uncompress them, and ensure you retain the .gz compressed versions.
    2. FASTQ files: Download FASTQ files (e.g., from NGDC database). Place them in a directory (e.g., data/fastqs), uncompress them, and ensure you retain the .gz compressed versions.
  5. Fast iteration of FASTA files without indexing

    master

    If you need to iterate through a FASTA file quickly without the overhead of building an index, use build_index=False.

    • To get (name, seq) tuples: use build_index=False.
    • To use the full header line as the identifier: use build_index=False and full_name=True.
    • To iterate over Sequence objects: use the default build_index=True (which allows access to attributes like .name, .seq, and .description).
    import pyfastx
    
    # Fastest iteration (returns name, seq)
    for name, seq in pyfastx.Fasta('test.fa.gz', build_index=False):
        print(name, seq)
    
    # Full header as name
    for name, seq in pyfastx.Fasta('test.fa', build_index=False, full_name=True):
        print(name, seq)
    
    # Iteration with Sequence objects (default)
    for seq in pyfastx.Fasta('test.fa.gz'):
        print(seq.name, seq.seq)
  6. Read and iterate FASTQ files

    master

    You can read plain or gzipped FASTQ files. There are two primary ways to iterate:

    1. Fastest (No Index): Set build_index=False. This returns a tuple of (name, seq, qual) for each read. This is ideal for one-pass processing.
    2. Random Access (With Index): The default behavior. Iteration returns a Read object, allowing access to attributes like .name, .seq, .qual (ASCII), and .quali (integer list).
    import pyfastx
    
    # Method 1: Fastest iteration (returns tuples)
    for name, seq, qual in pyfastx.Fastq('test.fq.gz', build_index=False):
        print(name, seq, qual)
    
    # Method 2: Iteration with Read objects (allows random access)
    for read in pyfastx.Fastq('test.fq.gz'):
        print(read.name, read.seq, read.qual, read.quali)
  7. Benchmark random access for FASTA and FASTQ

    master

    Benchmarking random access involves two steps: generating random access targets (sequence names or intervals) and then performing the access.

    FASTA Random Access (Sequence Names):

    1. Get 30% of sequence names with seed 123: python3 get_random_sequence_names_from_fasta.py 123 0.3 data/fastas/*.fa
    2. Run benchmark: ./benchmark_fasta_random_access.sh 3 data/fastas/*.fa > benchmark_result_fasta_random_access.tsv
    3. Process matrix: python3 make_fasta_ggplot2_matrix.py benchmark_result_fasta_random_access.tsv > matrix_fasta_random_access.tsv

    FASTA Random Access (Subsequences):

    1. Generate 1,000 random 1Kb intervals with seed 123: python3 generate_random_interval.py 123 data/fastas/*.fa
    2. Run benchmark: ./benchmark_fasta_extract_subsequences.sh 3 data/fastas/*.fa > benchmark_result_fasta_extract_subsequences.tsv
    3. Process matrix: python3 make_fasta_ggplot2_matrix.py benchmark_result_fasta_extract_subsequences.tsv > matrix_fasta_extract_subsequences.tsv

    FASTQ Random Access:

    1. Get 10,000 read names with seed 123: python3 get_random_read_names_from_fastq.py 123 data/fastqs/*.fq
    2. Run benchmark: ./benchmark_fastq_random_access.sh 3 data/fastqs/*.fq > benchmark_result_fastq_random_access.tsv
    3. Process matrix: python3 make_fastq_ggplot2_matrix.py benchmark_result_fastq_random_access.tsv > matrix_fastq_random_access.tsv
  8. Benchmark index building for FASTA and FASTQ

    master

    Run the following commands to benchmark the time and memory required to build indices for FASTA and FASTQ files. The results are processed into matrices for plotting.

    For FASTA:

    ./benchmark_fasta_build_index.sh 3 data/fastas/*.fa > benchmark_result_fasta_build_index.tsv
    python3 make_fasta_ggplot2_matrix.py benchmark_result_fasta_build_index.tsv > matrix_fasta_build_index.tsv

    For FASTQ:

    ./benchmark_fastq_build_index.sh 3 data/fastqs/*.fq > benchmark_result_fastq_build_index.tsv
    python3 make_fastq_ggplot2_matrix.py benchmark_result_fastq_build_index.tsv > matrix_fastq_build_index.tsv
  9. Use pyfastx with multiprocessing for parallel random access

    master

    To speed up random access to sequences, you can use pyfastx with Python's multiprocessing module. Because pyfastx uses a SQLite3 database for its index files, it supports concurrent access from multiple processes.

    Critical Requirement: You must ensure that the index file is created in the main process before spawning subprocesses. This ensures the SQLite3 database is initialized and available for the workers to access.

    import random
    import pyfastx
    import multiprocessing as mp
    
    # process worker
    # randomly fetch five sequences and print to stdout
    def worker(woker_num, seq_counts):
        # recreate the Fasta object in subprocess
        fa = pyfastx.Fasta('test.fa')
    
        for i in random.sample(range(seq_counts), 5):
            print("worker {} print:\n{}".format(worker_num, fa[i].raw))
    
    if __name__ == '__main__':
        # ensure index file has been created in main process
        fa = pyfastx.Fasta('test.fa')
    
        # get sequence counts
        c = len(fa)
    
        # start the process pool
        pool = mp.Pool()
    
        # add five task workers to run
        for n in range(5):
            pool.apply_async(worker, args=(n, c))
    
        # wait for tasks to finish
        pool.close()
        pool.join()
  10. Install pyfastx from source

    master

    To install from source, clone the repository and run the setup script. Note that pyfastx depends on zlib, sqlite3, and indexed_gzip. In the latest versions, pyfastx will automatically download these libraries during the build process.

    1. Clone the repository:
      git clone https://github.com/lmdu/pyfastx.git
    2. Navigate to the directory and install:
      cd pyfastx
      python setup.py install
    git clone https://github.com/lmdu/pyfastx.git
    cd pyfastx
    python setup.py install