pycisTopic
pycisTopic copied to clipboard
aertslab
SCENIC+: pycistopic QC errors
Hi, I'm trying to run the new SCENIC+ using this tutorial: https://pycistopic.readthedocs.io/en/latest/notebooks/human_cerebellum.html
When I get to the QC steps after inferring consensus peaks, pycistopic qc does not run. The error looks like this:
/.....site-packages/pycisTopic/topic_binarization.py:319: RuntimeWarning: invalid value encountered in divide mean1 = np.cumsum(hist * bin_centers) / weight1
..../site-packages/pycisTopic/topic_binarization.py:320: RuntimeWarning: invalid value encountered in divide mean2 = (np.cumsum((hist * bin_centers)[::-1]) / weight2[::-1])[::-1]
The strange thing is, it ran once already and I got some plots for sample level QC. But the TSS enrichment look very odd (attached image), and different from previous runs on the old SCENIC+ version....
When I ran the following code it seemed to work earlier:
pycistopic qc -f /rds/general/............../atac_fragments.tsv.gz -r /rds/general/................../consensus_peak_calling/consensus_regions.bed -t /rds/general/.................../consensus_peak_calling/tss_hg38.bed -o /rds/general/............../outs/qc/MCSF_16h_batch1 --min_fragments_per_cb 100 --tss_flank_window 1000 --tss_smoothing_rolling_window 10 --tss_minimum_signal_window 100 --tss_window 50
I added the the following patraemters to match the analysis I did in the old SCENIC+: --min_fragments_per_cb 100 --tss_flank_window 1000 --tss_smoothing_rolling_window 10 --tss_minimum_signal_window 100 --tss_window 50`
I'm uncertain if they make sense or if it's appropriate for the dataset.
Any suggestions or ideas on what could have gone wrong? The data is the single-nuclear ATAC modality from 10x multiome.
Thanks in advance!
What are the values for weight1 and weight2? The error seems to indicate that there are zeros there (or infinity or some other weird values).
Hi! Thanks for your message.
This is what the values are:
Lots of zeros..
I don't understand why..also not sure if things went wrong at the peak calling step or at consensus_peaks.to_bed.
The bed file looks like this:
Any thoughts on what's going wrong?
