Conversation
There was a problem hiding this comment.
Pull request overview
This PR extends DisjointIntervalSequence (DIS) with additional coordinate/interval manipulation and mapping utilities (including genomic ↔ DIS coordinate conversion and DNA extraction), updates documentation accordingly, and adds tests for the new behaviors. It also removes the prior DIS “index direction” toggle and introduces a Genome.dna(...) convenience path intended to accept a DIS directly.
Changes:
- Added DIS methods:
expand_coord,_lower_coord,lower,genomic_span,_lift_position,lift_interval,intersect, anddna; removed the global index-direction toggle. - Updated
Genome.dnabehavior to dispatch when passed aDisjointIntervalSequence. - Added/expanded unit tests and documentation for the new DIS coordinate mapping and DNA behaviors.
Reviewed changes
Copilot reviewed 6 out of 6 changed files in this pull request and generated 1 comment.
Show a summary per file
| File | Description |
|---|---|
| tests/test_genome.py | Adds coverage for Genome.dna() when called with an Interval vs a DIS. |
| tests/test_diseq.py | Adds extensive tests for DIS coordinate expansion, lowering/lifting, intersection, span, and DNA. |
| genome_kit/genome.py | Changes Genome.dna implementation to dispatch for DIS inputs. |
| genome_kit/diseq.py | Implements new DIS methods (expand/lower/lift/intersect/dna) and removes index-direction configuration. |
| genome_kit/init.py | Exposes DisjointIntervalSequence in __all__. |
| docs-src/diseq.rst | Documents expand_coord, coordinate mapping (lower/lift_interval), and DIS DNA behavior (including Genome.dna(dis) convenience). |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
This reverts commit 572bc63.
There was a problem hiding this comment.
Pull request overview
Copilot reviewed 4 out of 4 changed files in this pull request and generated 3 comments.
Comments suppressed due to low confidence (1)
docs-src/diseq.rst:819
- The docs claim
Genome.dnaaccepts a DIS directly (genome.dna(dis)), butGenomeDNA.__call__only acceptsInterval-typed arguments viaPyAsInterval(seesrc/py_interval.h). Unless there’s an unshown overload, this example will raiseTypeError. Either remove this claim/example or implement DIS dispatch in the DNA API.
For convenience, :py:meth:`Genome.dna <genome_kit.Genome.dna>` accepts a
DIS directly and dispatches to ``DisjointIntervalSequence.dna``, so
either of the following is equivalent::
>>> dis.dna()
>>> genome.dna(dis)
| ``dnstream`` bases. The segment is expanded an equal amount in the | ||
| upstream/downstream direction as the coordinate intervals. |
There was a problem hiding this comment.
What happens when the segment isn't adjacent to the coord system edge?
There was a problem hiding this comment.
it still get's expanded. So if the segment is inset 10bp on both 5' and 3' ends and then you expand the coordinate, it will remain so after the expansion. This also means if you have a 0-length interval and expand the coordinate space, the segment will no longer be 0-length after expansion.
There was a problem hiding this comment.
@AliceDG is this the desired behavior?
(this comment is not blocking)
There was a problem hiding this comment.
I don't think this is the desired behaviour, unless start, end represents the interval spanning the full-length coord.
One use case where we need to expand the coord is when we're constructing a DIS coord from an annotated transcript with e.g. inaccurate 3' annotation. In this case we'll create the DIS coord from the transcript first, and then expand it by say (0, 500). Note that I have not introduced any intervals on DIS coord at this point.
Under the proposed API, I imagine user code would look like this for the above use case:
tx = genome.transcripts["foo"]
dis_old_coord = DisjointIntervalSequence.from_transcript(tx) # important: this will default start-end to the full length
dis_new_coord = dis_old.expand_coord(0, 500) # now start-end also got expanded
# now I'm ready to define my interval: (10, 20) on the coord space
# approach 1: making use of `coord_end5, coord_end3`
dis_1 = dis_new_coord.coord_end5.expand(-10, 20)
# approach 2: alternatively, making use of `end5, end3`
dis_2= dis_new_coord.end5.expand(-10, 20)
- In the above use case, it's conceptually easier to think in terms of
coord_end5, coord_end3(and ignoreend5, end3), because we think about manipulating the coordinate when callingexpand_coord. - I don't really see a use case where:
start, endis not full-length & user needs to expandstart, endaccording to the coord expansion.
There was a problem hiding this comment.
@AliceDG How about expand_coord() then only expands start and end when the DIS has a full-length interval, and otherwise leaves the segment as-is?
Also, another thought related to the examples you gave, how would you feel about a function cut() that would let you create/set the DIS segment based on absolute positions within the DIS coord-space? So your example instead could be:
tx = genome.transcripts["foo"]
dis_old_coord = DisjointIntervalSequence.from_transcript(tx)
dis_new_coord = dis_old.expand_coord(0, 500)
dis_1 = dis_new_coord.cut(10, 20) # Equivalent to dis_new_coord.coord_end5.expand(-10, 20)
There was a problem hiding this comment.
Implemented the change to expand_coord() so that segment would only expand if it is equal to the coord intervals. If the segment is not equal it remains as is, while the coord-space grows. We can discuss this behaviour and change in a future PR if needed.
There was a problem hiding this comment.
How about expand_coord() then only expands start and end when the DIS has a full-length interval, and otherwise leaves the segment as-is?
Yes I think that makes more sense. As long as we document the behaviour to avoid confusion it's good!
how would you feel about a function cut() that would let you create/set the DIS segment based on absolute positions within the DIS coord-space
I think this is super useful! In fact a lot of the use cases involving disjoint interval coordinate are related to CDS / codons, so it'll be very convenient to be able to do something like this:
dis_cds = DisjointIntervalSequence(genome.transcripts["foo"].cdss)
# I'm interested in the 3rd codon
n = 3
itv_dis_codon = dis_cds.cut(3*(n-1), 3*n)
|
|
||
| ``lower`` projects the segment back to genomic space. ``lift_interval`` | ||
| projects a genomic interval into the DIS's index space and clips it | ||
| against the segment. They are conceptual inverses, but each has to deal |
There was a problem hiding this comment.
are they mechanically inverses as well? I remember the old code merged adjacent intervals. this one returns in sorted order?
There was a problem hiding this comment.
I've since updated the constructor to merge adjacent intervals, so that behaviour should be the same as with the old API, but yes the new one also returns in sorted order.
As for mechanical inverses, sort of. More specifically, if you were to lift an interval that contained introns of the transcript, then lowered the resulting DIS, you'd get a List[Interval] returned which would be comprised of only exons. If the initial lifted Interval did only consist of exons, then the functions would behave as inverses. For what it's worth, there is discussion about whether the lift function should raise an error if the Interval doesn't cleanly map onto the DIS segment, so this may change.
| on_coordinate_strand: bool | ||
|
|
||
|
|
||
| class DisjointIntervalSequence: |
There was a problem hiding this comment.
is a DIS supposed to handle the same protocol as a regular interval? EG, you don't need to if/else on the type. If so, I think maybe there's some missing methods like spanning/midpoint/distance etc
There was a problem hiding this comment.
The current understanding is no, DIS and Interval are separate and not interchangeable, however they do have significant overlap. As of right now, methods like spanning/distance are poorly defined for a DIS (is the returned value on the genomic coord, some combined shared coord, coord-space of just 1 of the DISes?), and other methods like midpoint are lower priority.
| cumulative += len(iv) | ||
| assert False, "Position not found in any interval" | ||
|
|
||
| def lift_interval( |
There was a problem hiding this comment.
My understanding is that this method finds the maximal-compatible interval on DIS by intersecting with DIS coord. I wonder if it's better to use a different method name, especially in our old API lift_interval raises an error if it's not compatible (since it doesn't do the intersection).
There was a problem hiding this comment.
How about this function doesn't do the intersection then, and instead we raise an error if the lifted segment isn't within the reference DIS segment?
There was a problem hiding this comment.
Going to leave as-is (intersect Interval with DIS when lifting) for now. We can discuss again later if this is needed.
There was a problem hiding this comment.
I think it's quite useful to have a method to do "intersect Interval with DIS when lifting". Maybe lift_interval can default to raising an error if the lifted segment isn't within the reference DIS segment (i.e. consistent with old API), and we add an arg to allow "intersection"?
also: