Best practices for long-read data

[olliecheng]

Hi,

Thanks for your work on umi-tools.

I’m trying to deduplicate some bulk long read data which has been genome-mapped. Since the start/end positions of each read are variable I am using featureCounts, which I have seen suggested in the documentation and the GitHub Issues:

featureCounts \
    -a /path/to/gencode.v49.annotation.gtf \
    -o path/to/results_counts \
    -L --primary --largestOverlap \
    -R CORE \
    /path/to/data.sam

Unfortunately, featureCounts is causing a segfault (died with <Signals.SIGSEGV: 11>). I’ve also been running into other issues with featureCounts on single-cell long-read data. Because -R BAM causes errors with long reads, I’m currently using -R CORE along with a custom script to write tags back to the reads, but some of my datasets are still causing an error despite this.

Do you have any recommendations on how best to run umi-tools for long-read data, both single-cell and bulk?

Thanks for your help,
Ollie

[IanSudbery]

Hi Olli,

Whether you use the start/end position, or the gene assignment not about how variable the ends, but whether PCR takes place before or after fragmentation. Many scRNA-seq procedures use 3'-end tagging protocols that add UMIs, PCR amplify, and then fragment - this means a transcript can end up represented by sequence reads with different mapping coordinates. If fragmentation happens before PCR, then you should be using mapping position, irrespective of how variable the ends are.

If you just need to assign long reads to genes (not transcripts), you could try extracting the "gene" entries from the gtf file, and using bedtools instead of featureCounts, and postprocessing the results. It seems you problem is memory, and if you sort both the bam and the gtf into chromosome/position order first, then bedtools should use very little memory, irrespective of the size of your files.