pileup around centromeres

[ratheraarif]

Hi,

My apologies if this question sounds naive.

I actually want to perform the pile ups around the centromeres in cool file that I am working with. the coordinates look like this.

I tried to execute the example from https://cooltools.readthedocs.io/en/latest/notebooks/pileup_CTCF.html.

since it is mostly related to CTCF bindings, I am not able to get that..

I would appreciate if you could help with the code example.

[sergpolly]

The coordinates look rather small ... Are you working with human data or is it some other organism ? If it is a different organism - are centromeres mappable there i.e. is there data around the centromeres on the heatmap ?

Would you like to study:
(1) how your centromeres interact across different chromosomes
(2) or you're trying to generate an "average" heatmap around your centromeres within chromosomes ?

(1) - is a work in progress in cooltools, but is available in coolpuppy https://coolpuppy.readthedocs.io/en/latest/Examples/Walkthrough_CLI.html#trans-inter-chromosomal-pileups , or you could try using saddle plots (we can discuss later how)
(2) you should be able to do right now with cooltools alone - this should be conceptually the same as the "on-diagonal CTCF pileups" from the notebook you've mentioned above

[ratheraarif]

Thank you for the quick reply. Actually I am working with the yeast dataset. You guessed it right I am interested in both of the following analysis.

(1) How centromeres interact across different chromosomes.
(2) Generate an "average" heatmap around all the centromeres within chromosomes ?

when I run the pileup from cooltools stack = cooltools.pileup(clr, Centro_sites, view_df=yeast_arms, flank=300_000)

it just returns Nans

need help..

[Phlya]

If your view is arms and doesn't contain centromeres, it won't work... Remove the view argument. Are the centromeres repetitive in your species? Then also make sure that flank is sufficiently big to actually capture the actual sequence outside of the repetitive regions.

[ratheraarif]

Thank you for the reply. No centromeres are not repetitive, they have a fixed position but are present across all the chromosomes. I tried to set the flank size large enough but it didn't run. However, setting the flank size small works fine.
Should I then set a flank size to a small value. If yes how small it should be?

Best.

[Phlya]

I don't know what's the best size of flank regions in this specific case. If the centromeres are short, then that point is not very relevant.
Did removing the view argument solve the issue?

[ratheraarif]

Yes it solved the problem and I was able to get some plots, but plots didn't seem to be as expected. Any other suggestions please.

[sergpolly]

The number of centromeres/chromosomes in yeast shouldn't be terribly large ? If you suspect that the pileup does something strange, I would just try looking at them individually first - just to understand better what we can expect from the pileups : mat=clr.matrix().fetch((chrom, start, end)) and then plt.imshow(mat) - you can get more plotting ideas from the viz notebook from open2c_examples repo

You'd want to confirm things like - the coordinates are perfectly aligned with the centromeres , individually they look like they should , etc etc

Sorry for brevity - typing on my phone 📱

[ratheraarif]

Thank you. will give it a try and get back to you soon...

[Phlya]

Has this been solved @ratheraarif ?

[ratheraarif]

Not really @Phlya, @sergpolly

I want the output to look like this:

This plot was also generated using cooltools for the same coordinates as #491 (comment)

This is how many papers do pileups around centromeres using cooltools.

I would appreciate it if you could help with the code examples too.

Thank you!

[Phlya]

What code did you try and what output did you get?

[ratheraarif]

I use the following code snippet

expected = cooltools.expected_cis(clrDD,nproc=2,ignore_diags=2,  chunksize=1_000_000)
stack = cooltools.pileup(clrDD, centromere_sites, expected_df=expected, flank=10_000)
print(stack.shape)
mtx = np.nanmean(stack, axis=2)
print(mtx.shape)

flank = 10_000
resolution = 1000
plt.imshow(
    np.log10(mtx),
    vmax = 1,
    vmin = -1,
    cmap='coolwarm')

plt.colorbar(label = 'log2 mean obs/exp')
ticks_pixels = np.linspace(0, flank*2//resolution,5)
ticks_kbp = ((ticks_pixels-ticks_pixels[-1]/2)*resolution//1000).astype(int)
plt.xticks(ticks_pixels, ticks_kbp)
plt.yticks(ticks_pixels, ticks_kbp)
plt.xlabel('relative position, kbp')
plt.ylabel('relative position, kbp')

plt.show()

the output is as:

However, I Expect a dot at the centre as shown here .

Not really @Phlya, @sergpolly

I want the output to look like this:

This plot was also generated using cooltools for the same coordinates as [#491 (comment)](https://github.com//issues/491#issue-2104276553)

This is how many papers do pileups around centromeres using cooltools.

I would appreciate it if you could help with the code examples too.

Thank you!

[Phlya]

From the look if it, the plot that you want to achieve - is that interaction between centromeres? So trans pileups? This is not actually (yet) supported by cooltools pileup...
Are you sure you should see anything interesting and different from your current result when looking at local interactions around centromeres?
Did you try suggestions by @sergpolly above - i.e. looking at them manually, for example? Since there aren't so many of them...

[ratheraarif]

Please have a look on the individual plots of each chromosome
flank_size = 50000

plots.pdf

please note that authors claim to have done the above analysis (#491 (comment)) as:

Each interaction between centromeres was extracted from corresponding mcool files in JupiterNotebooks using cooltools, generating text files with a range of +/- 50 kbp and a bin size of 1 kbp. Within a window of +/- 10 kbp around the central pixel, the average interactions were calculated and normalized against the mean interactions across all corners in both Hi-C and Micro-C matrices. Following this normalization, the mean interactions of all centromere-centromere matrices were computed, and the disparity between TKO and WT experiments for each pixel was determined.

can you help?

[Phlya]

The explanation above is about trans-interactions between centromeres (note that it is not possible with cooltools, they must have just used cooler).
Where is this from?

The plots you show above are very low resolution, while the pileups you tried were very high resolution, and also very different flank size.
Regardless, if you have code that extracts them one by one, you can simply average these single snippets manually.