Bienvenue sur IndexError.

Ici vous pouvez poser des questions sur Python et le Framework Django.

Mais aussi sur les technos front comme React, Angular, Typescript et Javascript en général.

Consultez la FAQ pour améliorer vos chances d'avoir des réponses à vos questions.

Utilisation de bitarrays

+1 vote

J'ai eu dans l'idée pour m'amuser de faire une bitarray qui coderait pour une séquence nucléique (4 bases possibles : ATGC). L'idée était de faire rentrer un caractère dans 2 bits tout en gardant une séquence mutable pour manipuler facilement la séquence. Pour cela j'ai pensé à prendre un bitarray et de l'encoder comme je le souhaitais, un peu comme ça :

from bitarray import bitarray
from random import choice

translation = {"A" : bitarray("00"), "T" : bitarray("01"), 
     "G" : bitarray("10"), "C" : bitarray("11")}


sequence = "".join(choice("ATGC") for _ in range(10))
print (sequence)

# coding
sequence_array = bitarray(sum([translation[base] for base in sequence], bitarray()))
print (sequence_array)

# decoding
sequence_01 = sequence_array.to01()
revert_translate = {v.to01() : k for k, v in translation.items()}
print ("".join(revert_translate[sequence_01[i:i+2]] for i in range(0, len(sequence_01), 2)))

Ca marche mais je suis persuadé que ce n'est pas du tout le bon moyen de faire les choses et j'ai plusieurs questions : Est-ce que c'est carrément une mauvaise idée ? Est-ce qu'il faudrait utiliser une autre lib ? Est-ce qu'il n'y aurait pas moyen de rendre l'encoding et le decoding bien plus rapide ? J'ai regardé du coté de struct mais j'imagine que ce n'est pas possible. Ou alors il faut passer sur du Cython voir même du C à ce niveau là ?

demandé 27-Avr par Jev (388 points)

2 Réponses

0 votes
 
Meilleure réponse

l'idée en soit

C'est pas mauvais, c'est même standard. Pour des raisons qui vont au-delà du simple gain de perfs, les formats binaires sont souvent un poil plus complexes que ça, mais l'idée de base est ici.

Concernant les struct, ces dernières sont un moyen plus bas niveau de manipuler la mémoire, surtout utilisées pour les dialogues C/Python quand on joue avec l'ABI.

Clairement, passer par un module tiers pour représenter des bits comme des booléens, c'est un peu tiré par les cheveux : en général on utilise plutôt des masques binaires, et l'API de bytearray (qui du coups donnerait 4 nucléotides / bytes).

Enfin, de mon point de vue, les compressions décompressions ne seront pas facilement simplifiables, surtout si c'est pour gagner en perf. Si c'est vraiment les perfs qui t'intéresse, fait d'abord ça en C, et dans un second temps apprend à interfacer C et Python.

Et pas la peine de passer à Cython, ce serais un bulldozer pour défourailler une mouche. types ou cffi font très bien le boulot pour des cas comme ceux-là.

standardisation

Plusieurs formats binaires existent en bioinformatique. L'un d'entre eux, formalisé par l'UCSC, gros acteur académique du domaine, se nomme le .2bit et un package python non-officiel nommé twobitreader permet de le manipuler.

D'autres existent. Les mot-clefs sont binary encoding of nucleotid sequence.

Si tu veux t'intéresser aux formats en général, et les manipuler de manière un peu «pro», Biopython est un package de référence en bioinformatique, présent dans tous les labos qui bossent en python sur des données bio, mais bizarrement il ne semble pas gérer les formats binaires pour les données génétiques. Peut-être est-ce juste un manque de la doc, qui est un bordel à peine masqué.

Dans tous les cas, les codes critiques sont généralement en C, avec un wrapper Python tout ce qu'il y a de plus standard.

répondu 27-Avr par lucas (2,292 points)
sélectionné 30-Avr par Jev

Merci pour la réponse. J'ai tenté de faire de quelque chose rapidement en utilisant la méthode encode et decode des bitarray, qui donne un résultat assez satisfaisant pour le moment. La prochaine étape sera de passer par du C quand j'aurais le temps. Je ne connaissais pas le format .2bit, c'est intéressant mais je reste un peu sceptique, c'est assez courant d'avoir des bases inconnues (N) dans les séquences. Sinon biopython, je ne connais que trop bien sa documentation ... Elle est correcte pour les modules les plus utilisées mais devient horrible à certains moments, particulièrement quand on doit chercher une info dans le code source.

Quitte à gérer les N, fait ça propre et gère tout IUPAC. Ça fait 16 signes à gérer, du coups tu peux compresser ça par groupe de 2 sur un byte. Je connais pas de format qui fait ça en bioinfo, mais je serais très surpris si ça n'existait pas (ça revient à de l'hexa).

+1 vote

Pour la performance (et la simplicité aussi), tu as loupé deux méthodes essentielles de bitarray : encode et decode. C'est beaucoup plus efficace que ton code.

translation = {
    "A": bitarray("00"),
    "T": bitarray("01"), 
    "G": bitarray("10"),
    "C": bitarray("11")
}

N = 100
sequence = [choice("ATGC") for _ in range(N)]

# encodage
arr = bitarray()
arr.encode(translation, sequence)

# decodage
seq = arr.decode(translation)

assert seq == sequence  # note que ce sont deux listes

Note: Il y a certainement moyen de faire plus rapide en C, vu que ton codage est de taille fixe, alors que ces méthodes prennent n'importe quel code préfixe.

répondu 9-Mai par yoch (2,312 points)
edité 9-Mai par yoch
...