3
$\begingroup$

I have a pyranges object that I am looking to iterate over. For each range, I would like to "resize" it to its centre (so start+end/2) and save that into a separate pyranges object. There will be additional processing associated based on what ranges already exist nearby, so I can't just use a resize function, I have to iterate over the object.

First, how would I go about manually creating a new pyranges object, then inserting one by one new ranges into it?

Second, when I try to iterate over my current pyranges object (using this line for myspot in pr.itergrs(hotspots, strand=False) :), it gives me this error: AttributeError: 'tuple' object has no attribute 'stranded'. This error occurs whether or not I have strands in my pyranges, and whether or not I specify strand=False. Why does this occur?


The data I am working with: a pyranges covering hotspots in the human genome (all chr1 for the moment)

+--------------+-----------+-----------+-----------+
| Chromosome   | Start     | End       | Score     |
| (category)   | (int32)   | (int32)   | (int64)   |
|--------------+-----------+-----------+-----------|
| chr1         | 154421820 | 154421821 | 26413     |
| chr1         | 30749469  | 30749470  | 14833     |
| chr1         | 55259503  | 55259504  | 9900      |
| chr1         | 14985097  | 14985098  | 7735      |
| ...          | ...       | ...       | ...       |
| chr1         | 247975751 | 247976497 | 1         |
| chr1         | 247977599 | 247981391 | 1         |
| chr1         | 248018012 | 248019454 | 1         |
| chr1         | 248177298 | 248178490 | 1         |
+--------------+-----------+-----------+-----------+

The hotspots are sorted by their score, and I'd like to go through and for each one, resize it to its centre, keeping the score, but if there is already a hotspot within a certain range (say +/-5kb), to discard it instead, so that I end up with clear chunks between hotspots. This is the "pseudocode" I was planning to do to implement this:

add centre of hotspot 1 to a new pyranges object
create ranges with score = -1 for 5kb flanking each side
for myhotspot in hotspots:
  resize to centre
  if not overlap with new pyranges object
    add centre
    create flanking regions
discard ranges with score<0

This may not be the most elegant solution, but I'll never have more than ~2000 ranges to process so it should work fine. Is this implementable with pyranges?

$\endgroup$

2 Answers 2

1
$\begingroup$

Quick and dirty (not to mention atrociously bad O(n)) solution:

import pyranges as pr
import numpy as np
np.random.seed(42 * 10)

# create large df to test on
gr = pr.random(int(1e5), length=10000, chromsizes={"chr1": 249250621})
gr.Score = np.random.randint(250, size=len(gr))

def remove_worst_scores_until_no_overlap(gr):
    df = gr.df
    old_length = -1
    new_length = len(df)
    while old_length != new_length:
        df = df.drop(["Cluster", "Count"], axis=1, errors="ignore")
        df = pr.PyRanges(df).cluster(count=True).df.sort_values("Score", ascending=False)
        df = df[df.duplicated('Cluster', keep='last') | (df.Count == 1)]
        old_length = new_length
        new_length = len(df)

    return pr.PyRanges(df.drop("Cluster Count".split(), axis=1)).sort()

hotspots_no_overlaps = remove_worst_scores_until_no_overlap(gr)
$\endgroup$
1
  • 1
    $\begingroup$ Quick and dirty, but it works! $\endgroup$
    – Whitehot
    Commented Mar 2, 2021 at 15:01
4
$\begingroup$

Answer to main question:

import pyranges as pr
assert int(pr.__version__.split(".")[2]) >= 93, "pip install pyranges==0.0.93"
import numpy as np
np.random.seed(42 * 10)

# create large df to test on
gr = pr.random(int(1e5), length=1, chromsizes={"chr1": 249250621})
gr.Score = np.random.randint(250, size=len(gr))

result = gr.max_disjoint(slack=5000)

# small example that is easy to validate
gr2 = pr.from_string("""Chromosome Start End Score
chr1 1 2 5
chr1 4999 5000 10
chr1 10001 10002 1
chr1 15000 15001 2
chr1 20100 20101 5
""")

gr2.max_disjoint(slack=5000)

# +--------------+-----------+-----------+-----------+
# | Chromosome   |     Start |       End |     Score |
# | (category)   |   (int32) |   (int32) |   (int64) |
# |--------------+-----------+-----------+-----------|
# | chr1         |         1 |         2 |         5 |
# | chr1         |     10001 |     10002 |         1 |
# | chr1         |     20100 |     20101 |         5 |
# +--------------+-----------+-----------+-----------+
# Unstranded PyRanges object has 3 rows and 4 columns from 1 chromosomes.
# For printing, the PyRanges was sorted on Chromosome.

Could you please post an example of what data you have and what you want to end up with?

  1. Inserting row after row is going to be slow, whether it is pandas or pyranges. It is best to collect all the data and then create a new dataframe or PyRanges.

  2. pr.itergrs is for iterating over multiple pyranges at the same time. It takes a list of pyranges while you seem to have a collection of tuples. The function iterates over every chromosome/strand pair in the list and for it extracts the data from every pyranges. See the code example at the end.

To iterate over a pyranges you can do for k, df in gr: # ....

Anyways, to do what you want:

import pyranges as pr
import numpy as np

f1 = pr.data.f1()
f1.Start = ((f1.Start + f1.End)/2).astype(np.int32) # int is 64 bit
f1.End = f1.Start + 1
  1. I understand I have not done the addititonal processing you speak of, but if you can post an example I'm sure I'd be able to help you :)

Code example pr.itergrs:

import pyranges as pr
f1 = pr.data.f1()
f2 = pr.data.f2()
for k, dfs in pr.itergrs([f1, f2], keys=True):
    print(k)
     for df in dfs:
         print("df:")
         print(df)

('chr1', '+')
df:
  Chromosome  Start  End       Name  Score Strand
0       chr1      3    6  interval1      0      +
2       chr1      8    9  interval3      0      +
df:
  Chromosome  Start  End Name  Score Strand
0       chr1      1    2    a      0      +
('chr1', '-')
df:
  Chromosome  Start  End       Name  Score Strand
1       chr1      5    7  interval2      0      -
df:
  Chromosome  Start  End Name  Score Strand
1       chr1      6    7    b      0      -
$\endgroup$
13
  • $\begingroup$ Thanks for the pointers, this already helps a lot. I'll post my code in the morning (it's 7PM for me now) $\endgroup$
    – Whitehot
    Commented Feb 24, 2021 at 18:03
  • $\begingroup$ Added the information you asked for :) $\endgroup$
    – Whitehot
    Commented Feb 25, 2021 at 14:37
  • 1
    $\begingroup$ I've now added a method max_disjoint to PyRanges in 0.0.93. Thanks for making me realise there was a fundamental interval operation I was missing. $\endgroup$ Commented Feb 25, 2021 at 20:53
  • 1
    $\begingroup$ Will fix. Could you post an issue on the pyranges github? $\endgroup$ Commented Feb 26, 2021 at 19:54
  • 1
    $\begingroup$ max_disjoint() is actually a well-known greedy algorithm: en.wikipedia.org/wiki/… It starts by sorting on the end column. I will continue thinking about the problem $\endgroup$ Commented Mar 1, 2021 at 18:02

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.