diff --git a/intact/intact.py b/intact/intact.py index 5ed3130..485586e 100644 --- a/intact/intact.py +++ b/intact/intact.py @@ -710,6 +710,19 @@ def analyse_single_sequence(holistic, sequence, blast_rows): reference = subtype_choices[reference_name] aligned_subtype = AlignedSequence(this=reference, reference=st.HXB2()) + forward_aligned_sequence = AlignedSequence(this=sequence, reference=aligned_subtype.this) + reverse_aligned_sequence = forward_aligned_sequence.reverse() + + forward_score = forward_aligned_sequence.alignment_score() + reverse_score = reverse_aligned_sequence.alignment_score() + if forward_score >= reverse_score: + aligned_sequence = forward_aligned_sequence + else: + print("Reversing sequence " + sequence.id + "; forward score " + + str(forward_score) + "; reverse score " + str(reverse_score)) + aligned_sequence = reverse_aligned_sequence + sequence = aligned_sequence.this + # convert ORF positions to appropriate subtype forward_orfs, reverse_orfs, small_orfs = [ [ @@ -719,10 +732,6 @@ def analyse_single_sequence(holistic, sequence, blast_rows): for orfs in [hxb2_forward_orfs, hxb2_reverse_orfs, hxb2_small_orfs] ] - # convert PSI locus and RRE locus to appropriate subtype - psi_locus = [ReferenceIndex(x).mapto(aligned_subtype) for x in hxb2_psi_locus] - rre_locus = [ReferenceIndex(x).mapto(aligned_subtype) for x in hxb2_rre_locus] - holistic.orfs_start = min(forward_orfs, key=lambda e: e.start).start holistic.orfs_end = max(forward_orfs, key=lambda e: e.end).end clamp = lambda p: max(min(p, holistic.orfs_end), holistic.orfs_start) @@ -730,18 +739,9 @@ def analyse_single_sequence(holistic, sequence, blast_rows): blast_matched_orfs_slen = holistic.orfs_end - holistic.orfs_start holistic.blast_sseq_orfs_coverage = aligned_reference_orfs_length / blast_matched_orfs_slen - forward_aligned_sequence = AlignedSequence(this=sequence, reference=aligned_subtype.this) - reverse_aligned_sequence = forward_aligned_sequence.reverse() - - forward_score = forward_aligned_sequence.alignment_score() - reverse_score = reverse_aligned_sequence.alignment_score() - if forward_score >= reverse_score: - aligned_sequence = forward_aligned_sequence - else: - print("Reversing sequence " + sequence.id + "; forward score " - + str(forward_score) + "; reverse score " + str(reverse_score)) - aligned_sequence = reverse_aligned_sequence - sequence = aligned_sequence.this + # convert PSI locus and RRE locus to appropriate subtype + psi_locus = [ReferenceIndex(x).mapto(aligned_subtype) for x in hxb2_psi_locus] + rre_locus = [ReferenceIndex(x).mapto(aligned_subtype) for x in hxb2_rre_locus] alignment = aligned_sequence.get_alignment()