<jupyter_start><jupyter_text>`seqp`: how to use sharded storage  IntroductionFor this example, we will be using DNA data to illustrate sharded storage with `seqp`.We will:1. Download a DNA data file in [FASTA format](https://en.wikipedia.org/wiki/FASTA_format).2. Parse the file with [biopython](http://biopython.org/DIST/docs/tutorial/Tutorial.htmlhtoc11).3. Store the sequences from the file in multiple HDF5 shards, annotating each DNA sequence with its identifier.4. Read the DNA data back. First, we will install some python libraries...<jupyter_code>pip install biopython tqdm numpy<jupyter_output>Note: you may need to restart the kernel to use updated packages.<jupyter_text>Download the data and take a look at it<jupyter_code>!wget 'ftp://ftp.ebi.ac.uk/pub/databases/ena/coding/release/con/fasta/rel_con_hum_r143.cds.fasta.gz' -O rel_con_hum_r143.cds.fasta.gz
!gunzip -f 'rel_con_hum_r143.cds.fasta.gz'
from Bio import SeqIO

file_name = 'rel_con_hum_r143.cds.fasta'

for k, seq_record in enumerate(SeqIO.parse(file_name, "fasta")):
    print("ID: {}".format(seq_record.id))
    print("   - Sequence: {}...{}".format(seq_record.seq[:20], seq_record.seq[-10:]))
    print("   - Lenght: {}".format(len(seq_record)))
    if k > 3:
        break<jupyter_output>ID: ENA|EAL24309|EAL24309.1
   - Sequence: atgaagcatgtgttgaacct...aagcatgtga
   - Lenght: 192
ID: ENA|EAL24310|EAL24310.1
   - Sequence: atggaggggccactcactcc...gctgtactga
   - Lenght: 1800
ID: ENA|EAL24311|EAL24311.1
   - Sequence: atggaggggccactcactcc...gctgtactga
   - Lenght: 1692
ID: ENA|EAL24312|EAL24312.1
   - Sequence: atggacccaaggacatccag...gacctcctga
   - Lenght: 276
ID: ENA|EAL24313|EAL24313.1
   - Sequence: atggccaggcatggctgtct...agacctgtga
   - Lenght: 2082<jupyter_text>Read the DNA data and store it in HDF5 with `seqp` as we goWe want to use `Hdf5RecordWriter` to write DNA sequences to files. We also want to write to multiple HDF5 files, each one containing up to a maximum amount of records, so we make use of a `ShardedWriter` decorator.Once we have our writer, we iterate over the FASTA file sequences and store them in the writer. Each nucleotide is saved as a byte-sized integer number obtained by subtracting the ASCII index of 'a' to the nucleotide letter.Once all the sequences are written to files, we write a piece of metadata with a dictionary from the protein name and the index within the files.<jupyter_code>import json
import numpy as np
from seqp.hdf5 import Hdf5RecordWriter
from seqp.record import ShardedWriter
from tqdm import tqdm

def nucleotide2num(letter: str) -> int:
    """ Converts a nucleoide letter to an integer"""
    return ord(letter.lower()) - ord('a')

protein2idx = dict()
output_file_template = "dna_example_{:02d}.hdf5"

with ShardedWriter(Hdf5RecordWriter,
                   output_file_template,
                   max_records_per_shard=5000) as writer:

    for seq_record in tqdm(SeqIO.parse(file_name, "fasta")):
        _, _, protein = seq_record.id.split('|')
        sequence = [nucleotide2num(letter) for letter in seq_record.seq]
        idx = writer.write(np.array(sequence, dtype=np.uint8))
        protein2idx[protein] = idx

    writer.add_metadata({'protein_idx': json.dumps(protein2idx)})<jupyter_output>65183it [01:13, 888.76it/s]<jupyter_text>Read the HDF5 records backWe open the HDF5 files with a `Hdf5RecordReader`. First, we read back the dictionary with the indexes of each protein sequence, and then we retrieve the sequences associated with some specific target proteins.<jupyter_code>from glob import glob
import json
from seqp.hdf5 import Hdf5RecordReader

target_proteins = ['EAL24309.1', 'EAL24312.1']

def num2nucleotide(num: int) -> str:
    """ Converts an integer to a nucleoide letter"""
    return chr(num + ord('a'))

with Hdf5RecordReader(glob('dna_example_*.hdf5')) as reader:
    loaded_protein2idx = json.loads(reader.metadata('protein_idx'))
    indexes = set(reader.indexes())
    for protein in target_proteins:
        sequence = reader.retrieve(protein2idx[protein])
        sequence = "".join([num2nucleotide(n) for n in sequence.tolist()])
        print("{} : {}...{}".format(protein, sequence[:20], sequence[-10:]))<jupyter_output>EAL24309.1 : atgaagcatgtgttgaacct...aagcatgtga
EAL24312.1 : atggacccaaggacatccag...gacctcctga<jupyter_text>The same but with records with fieldsWe just need to specify the `sequence_field` and `fields` parameters, and provide a dictionary with such keys to function `write`.<jupyter_code>output_file_template = "fields_dna_example_{:02d}.hdf5"

SEQUENCE_FIELD = 'seq'
OTHER_FIELD = 'other'
FIELDS = [SEQUENCE_FIELD, OTHER_FIELD]

with ShardedWriter(Hdf5RecordWriter,
                   output_file_template,
                   sequence_field=SEQUENCE_FIELD,
                   fields=FIELDS,
                   max_records_per_shard=5000) as writer:

    for seq_record in tqdm(SeqIO.parse(file_name, "fasta")):
        _, _, protein = seq_record.id.split('|')
        sequence = [nucleotide2num(letter) for letter in seq_record.seq]
        record = {'seq': np.array(sequence, dtype=np.uint8),
                  'other': np.array(sequence[:10], dtype=np.uint8)}
        idx = writer.write(record)
        protein2idx[protein] = idx

    writer.add_metadata({'protein_idx': json.dumps(protein2idx)})<jupyter_output>65183it [01:26, 757.52it/s]<jupyter_text>When reading back the records, we don't need to specify the fields, just `retrieve` would return a dictionary instead of a plain `np.ndarray`.<jupyter_code>with Hdf5RecordReader(glob('fields_dna_example_*.hdf5')) as reader:
    loaded_protein2idx = json.loads(reader.metadata('protein_idx'))
    indexes = set(reader.indexes())
    for protein in target_proteins:
        record = reader.retrieve(protein2idx[protein])
        sequence = record[SEQUENCE_FIELD]
        sequence = "".join([num2nucleotide(n) for n in sequence.tolist()])
        print("{} : {}...{}".format(protein, sequence[:20], sequence[-10:]))<jupyter_output>EAL24309.1 : atgaagcatgtgttgaacct...aagcatgtga
EAL24312.1 : atggacccaaggacatccag...gacctcctga
{'fraction_non_alphanumeric': 0.10981423639651487, 'fraction_numerical': 0.03600197271083347, 'mean_word_length': 4.253886010362694, 'pattern_counts': {'":': 0, '<': 19, '<?xml version=': 0, '>': 22, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '604971', 'license': 'MIT', 'n_tokens_mistral': 2220, 'n_tokens_neox': 2032, 'n_words': 532, 'path': 'examples/sharded_storage.ipynb', 'repo_name': 'noe/seqp', 'size': '6083'}