-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathtasks.py
66 lines (50 loc) · 2.59 KB
/
tasks.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
from Bio import SeqIO
from Bio.Blast import NCBIXML
from Bio.Blast import NCBIWWW
def run_blast(sequence, program, database, evalue, id): #TODO
# Prendre la première séquence pour l'exemple
sequences = list(SeqIO.parse(sequence, "fasta"))
sequence_ = sequences[0].seq
# Exécuter BLAST avec les paramètres spécifiés (Pour la première séquence)
result_handle = NCBIWWW.qblast(program, database, sequence_, expect=evalue)
# Sauvegarder les résultats dans un fichier XML
with open("blast_results.xml", "w") as out_handle:
out_handle.write(result_handle.read())
print("Les résultats BLAST ont été enregistrés dans 'blast_results.xml'.")
return result_handle
def parse_blast_results(blast_result_handle):
# Lire et parser le fichier XML
blast_records = NCBIXML.parse(blast_result_handle)
# Parcourir les résultats de BLAST
for blast_record in blast_records:
print(f"Query: {blast_record.query}")
# Vérifier si des hits ont été trouvés
if blast_record.alignments:
for alignment in blast_record.alignments:
print(f" Hit: {alignment.hit_id}")
print(f" Description: {alignment.hit_def}")
# Vérifier si 'hit_len' est disponible
if hasattr(alignment, 'hit_len'):
print(f" Length: {alignment.hit_len}")
else:
print(" Length: non disponible")
# Afficher les top alignments
for hsp in alignment.hsps:
print(f" E-value: {hsp.expect}")
print(f" Score: {hsp.score}")
print(f" Identity: {hsp.identities}/{hsp.align_length} ({(hsp.identities / hsp.align_length) * 100:.2f}%)")
print(f" Query start: {hsp.query_start}, Query end: {hsp.query_end}")
print(f" Hit start: {hsp.sbjct_start}, Hit end: {hsp.sbjct_end}")
print("-" * 50)
else:
print(" Aucun hit trouvé.")
print("-" * 50)
# Exemple d'utilisation
if __name__ == "__main__":
sequences = get_fasta_file()
if sequences:
program, database, evalue = get_blast_parameters()
blast_results = run_blast(sequences, program, database, evalue)
# Lire les résultats BLAST enregistrés pour les analyser
with open("blast_results.xml", "r") as result_handle:
parse_blast_results(result_handle)