diff --git a/src/cool_seq_tool/mappers/exon_genomic_coords.py b/src/cool_seq_tool/mappers/exon_genomic_coords.py index a03e604..82d4a72 100644 --- a/src/cool_seq_tool/mappers/exon_genomic_coords.py +++ b/src/cool_seq_tool/mappers/exon_genomic_coords.py @@ -1169,14 +1169,28 @@ def _get_adjacent_exon( :param end: Genomic coordinate of breakpoint :return: Exon number corresponding to adjacent exon. Will be 0-based """ - for i in range(len(tx_exons_genomic_coords) - 1): + # Check if a breakpoint occurs before/after the transcript boundaries + bp = start if start else end + exon_list_len = len(tx_exons_genomic_coords) - 1 + + if strand == Strand.POSITIVE: + if bp < tx_exons_genomic_coords[0].alt_start_i: + return 0 + if bp > tx_exons_genomic_coords[exon_list_len].alt_end_i: + return exon_list_len + if strand == Strand.NEGATIVE: + if bp > tx_exons_genomic_coords[0].alt_end_i: + return 0 + if bp < tx_exons_genomic_coords[exon_list_len].alt_start_i: + return exon_list_len + + for i in range(exon_list_len): exon = tx_exons_genomic_coords[i] if start == exon.alt_start_i: break if end == exon.alt_end_i: break next_exon = tx_exons_genomic_coords[i + 1] - bp = start if start else end if strand == Strand.POSITIVE: lte_exon = exon gte_exon = next_exon diff --git a/tests/mappers/test_exon_genomic_coords.py b/tests/mappers/test_exon_genomic_coords.py index 07ff3ff..c779953 100644 --- a/tests/mappers/test_exon_genomic_coords.py +++ b/tests/mappers/test_exon_genomic_coords.py @@ -840,6 +840,32 @@ async def test_get_adjacent_exon( ) assert resp == 4 + # Check cases where breakpoint occurs in before/after transcript boundaries + resp = test_egc_mapper._get_adjacent_exon( + tx_exons_genomic_coords=nm_001105539_exons_genomic_coords, + start=80486220, + strand=Strand.POSITIVE, + ) + assert resp == 0 + resp = test_egc_mapper._get_adjacent_exon( + tx_exons_genomic_coords=nm_001105539_exons_genomic_coords, + start=80526285, + strand=Strand.POSITIVE, + ) + assert resp == 5 + resp = test_egc_mapper._get_adjacent_exon( + tx_exons_genomic_coords=nm_152263_exons_genomic_coords, + end=154192110, + strand=Strand.NEGATIVE, + ) + assert resp == 0 + resp = test_egc_mapper._get_adjacent_exon( + tx_exons_genomic_coords=nm_152263_exons_genomic_coords, + end=154161809, + strand=Strand.NEGATIVE, + ) + assert resp == 9 + def test_is_exonic_breakpoint(test_egc_mapper, nm_001105539_exons_genomic_coords): """Test is breakpoint occurs on exon"""