Hi,
I have to analyze 8 single-NUCLEI RNA-seq samples from patient-derived biopsies. 4 patients, each with pre-/post-treatment. The samples are from difficult-to-dissociate tissues stored a long time ago. On the wet lab everything looked fine and it seemed we obtained a decent amount of single nuclei under the microscope. We sequenced them with 10x, targeting 5k cells/sample. We used 2x pools of 4-sample On-Chip Multiplexing.
After sequencing, the results from cellranger multi looked disappointing. 10x support suggested to give a try to Cellbender, to recover more "usable" barcodes/cells by removing ambient RNA, however, the results are still pretty bad:
| raw barcodes | cellranger-multi | cellbender (whatever #cells) | cellbender (--expected-cells=5k) | |
|---|---|---|---|---|
| pool1A_p1_t1 | 240748 | 530 | 794 | 4234 |
| pool1B_p1_t2 | 224959 | 898 | 1406 | 3176 |
| pool1C_p2_t1 | 414563 | 3157 | 4244 | 4631 |
| pool1D_p2_t2 | 556365 | 5684 | 41782 | 40905 |
| pool2A_p3_t1 | 323891 | 972 | 1245 | 4832 |
| pool2B_p3_t2 | 800648 | 9519 | 43903 | 42810 |
| pool2C_p4_t1 | 89333 | 44656 | 26286 | 18037 |
| pool2D_p4_t2 | 98760 | 283 | 920 | 4675 |
Something went obviously wrong with these samples. Both samples from patient4 are unusable. And the samples with high number of barcodes (p2_t2 and p3_t2) have suspicious-looking barcode x UMI plots:
pool2A_p3_t1 (decent, 900~1200 cells):
pool2B_p3_t2 (unreliable, 9.5k~43k cells):
All (non-garbage) samples present this peak, at a different number of UMIs:
And these are the ranked barcode x number of UMIs:
I don't have much experience with Cellbender, but I assume this kind of peak is common, probably caused by ambient RNA, and non-empty cells lie to the left of it. I suspect that most of my samples have large amounts of ambient RNAs from the cytoplasm.
I am already re-running cellbender following the suggestions from the warning messages (halving the Learning Rate does nothing) and from the guidelines (manually selected --expected-cells and --total-droplets-included per sample.
I understand that cellbender is just a first step in QC, and it might be better to be lenient on its barcode selection, to later further filter them downstream with metrics like %mtRNA.
I was wondering what to do with metrics like the fraction of intronic reads. Droplets containing only nuclei should have high % of unspliced transcripts. This, combined with %mtRNA is key to selecting droplets with actual nuclei inside, instead of large amounts of (cytoplasmic) ambient RNA. However, I don't see any option in Cellbender to take these parameters into account, in order to better select empty droplets to create the ambient RNA profile.
Are there alternative tools to Cellbender designed specially for single nuclei sequencing? If not, in which order should I run the different tools? Should I run velocyto on the raw data, or only after filtering? Can I rely on Cellbender's ambient correction? or should I use it first only for non-empty barcode selection, then further filter the data by %mt and %intronic, and finally ambient removal + doublet detection?
Alternatively, is it possible that those "peaks" seen in the UMI counts x number of droplet plots are in fact the properly-isolated nuclei, and droplets to their left are droplets containing a "decent" number of nuclei reads + large number of ambient cytoplasmic reads, or even whole cells?
PS: yes, I am aware that the samples look bad, but I have to make the most out of them.
In reply to @Bastien Hervé's comment, as the forum seems to not like tables in comments.
From the BAM file, how is p4_t1 (90k raw barcodes -> 45k cells) different from p3_t2 (800k raw barcodes -> 9.5k cell)
The sample_alignments.bam (multi's equivalent to count's possorted_genome_bam.bam) have wildly different sizes, matching the differences in number of barcodes.
| sample_alignments.bam | |
|---|---|
| pool1A_p1_t1 | 2.0G |
| pool1B_p1_t2 | 1.9G |
| pool1C_p2_t1 | 5.0G |
| pool1D_p2_t2 | 9.0G |
| pool2A_p3_t1 | 3.3G |
| pool2B_p3_t2 | 22G |
| pool2C_p4_t1 | 227M |
| pool2D_p4_t2 | 389M |
The huge differences between p3_t2 and both p4_tx samples I assume are driven by those cells being predominant during sequencing. And similarly in pool1 with p2_t2.
I also ran cellranger count on the merged pools and obtained similar sizes for the two possorted_genome_bam.bam (18G for pool1, 27G for pool2).
These are the metrics from multi's report for the pooled samples (pre-demultiplexing) and my run of count. The only differences between methods are derived from calling more/less cells. The mapping metrics are identical.
| multi - pool1 | count - pool1 | multi - pool2 | count - pool2 | ||
|---|---|---|---|---|---|
| Cell Calling Quality | Cells | 10,269 | 15,029 | 55,430 | 8,097 |
| Conf. map. reads in cells | 43.40% | 53.80% | 67.10% | 64.00% | |
| Median UMI Counts per Cell | 2,279 | 9,070 | |||
| Median Genes per Cell | 1,538 | 3,591 | |||
| Total Genes Detected | 32,267 | 32,455 | |||
| Mapping Quality | Conf. map. to transcriptome | 67% | 67.00% | 70.60% | 70.60% |
| Map. to genome | 95.80% | 95.80% | 95.90% | 95.90% | |
| Conf. map. to genome | 83.80% | 83.80% | 86% | 86.00% | |
| Conf. map. to exonic regions | 24.50% | 24.50% | 26.40% | 26.40% | |
| Conf. map. to intronic regions | 53% | 53.00% | 53.30% | 53.30% | |
| Conf. map. to intergenic regions | 6.20% | 6.20% | 6.20% | 6.20% | |
| Conf. map. antisense | 9.60% | 9.60% | 8.20% | 8.20% | |
| Sequencing Quality | Number of reads | 299,137,975 | 299,137,975 | 445,813,475 | 445,813,475 |
| Mean reads per cell | 29,130 | 19,904 | 8,043 | 55,059 | |
| Sequencing saturation | 53.60% | 53.60% | 57.20% | 57.20% | |
| Valid barcodes | 95.30% | 95.30% | 96.20% | 96.20% | |
| Valid UMIs | 100% | 100.00% | 99.90% | 99.90% | |
| Q30 barcodes | 93.40% | 93.40% | 93.60% | 93.60% | |
| Q30 UMI | 91.90% | 91.90% | 91.80% | 91.80% | |
| Q30 RNA read | 89% | 89.00% | 88.90% | 88.90% |
The gross metrics per OCM of multi. There are huge differences in the number of UMIs per sample.
| Metrics per OCM Barcode | ||||
|---|---|---|---|---|
| pool | OCM Barcode ID | Sample ID | UMIs per OCM barcode | Cells per OCM barcode |
| pool1 | OB1 | pool1A_p1_t1 | 7,078,433 (7.9%) | 530 (5.2%) |
| pool1 | OB2 | pool1B_p1_t2 | 7,207,461 (8%) | 898 (8.7%) |
| pool1 | OB3 | pool1C_p2_t1 | 25,712,718 (28.6%) | 3,157 (30.7%) |
| pool1 | OB4 | pool1D_p2_t2 | 49,781,294 (55.4%) | 5,684 (55.4%) |
| pool2 | OB1 | pool2A_p3_t1 | 15,318,842 (11.7%) | 972 (1.8%) |
| pool2 | OB2 | pool2B_p3_t2 | 113,041,243 (86.5%) | 9,519 (17.2%) |
| pool2 | OB3 | pool2C_p4_t1 | 799,817 (0.6%) | 44,656 (80.6%) |
| pool2 | OB4 | pool2D_p4_t2 | 1,483,499 (1.1%) | 283 (0.5%) |
And the detailed metrics per sample from multi, with huge differences in the number of reads per sample:
| pool1A_p1_t1 | pool1B_p1_t2 | pool1C_p2_t1 | pool1D_p2_t2 | pool2A_p3_t1 | pool2B_p3_t2 | pool2C_p4_t1 | pool2D_p4_t2 | ||
|---|---|---|---|---|---|---|---|---|---|
| Cell Calling Quality | Cells | 530 | 898 | 3,157 | 5,684 | 972 | 9,519 | 44,656 | 283 |
| Conf. map. reads in cells | 32.10% | 31% | 45.80% | 44.90% | 37% | 70.70% | 96.80% | 7.40% | |
| Median genes per cell | 1,687 | 1,228 | 1,615 | 1,881 | 2,032 | 3,219 | 16 | 210 | |
| Median UMI counts per cell | 2,670 | 1,796 | 2,441 | 2,913 | 3,439 | 7,256 | 16 | 226 | |
| Total genes detected | 23,077 | 22,647 | 28,003 | 30,056 | 25,459 | 32,398 | 21,368 | 13,828 | |
| Number of reads from cells called from this sample | 5,902,300 | 6,254,289 | 38,276,995 | 69,053,447 | 14,913,734 | 261,457,737 | 3,041,221 | 308,013 | |
| Mapping Quality (Amongst Reads From Cells Assigned To Sample) | Conf. map. to transcriptome | 68% | 68.40% | 69.20% | 71.80% | 73.40% | 74.20% | 33% | 50.60% |
| Map. to genome | 95.10% | 95.30% | 96.60% | 96.70% | 96.20% | 96.70% | 82.30% | 93.10% | |
| Conf. map. to genome | 82.30% | 82.40% | 92.50% | 92.40% | 89.30% | 91.80% | 41.70% | 61.40% | |
| Conf. map. to exonic regions | 25.80% | 24% | 22.20% | 22.10% | 25% | 25.60% | 15.70% | 18.70% | |
| Conf. map. to intronic regions | 49.40% | 51.30% | 63% | 63% | 57.50% | 59.30% | 22% | 37% | |
| Conf. map. to intergenic regions | 7.10% | 7.10% | 7.30% | 7.40% | 6.80% | 7% | 4% | 5.80% | |
| Conf. map. antisense | 6.30% | 6% | 15% | 12.10% | 8.10% | 9.60% | 4.30% | 4.40% |
•
written
by
txema.heredia •
0 answers
No answers yet.
Log in to answer this question.
I think you are going a bit fast on graphical interpretation of your outputs.
I would check first :
Any reason of using cellranger multi over cellranger count ?
I don't know your wet lab protocol but in my knowledge nuclei is less prone to ambient RNA and doublets than cells, due to nuclei isolation.
Thanks for your reply
These are multiplexed samples using 10x Flex. 4 samples were pooled together using unique barcodes for each sample. I double-checked and that seemed to have worked fine. Running
cellranger countinstead just gives me all 4 samples pooled together.AFAIK, these samples come from metastases. All of the patients have the same cancer type. Some of the samples I think are liver metastases, but others might come from other organs, with different "contexts" affecting tissue dissociation. Including more or less surrounding healthy tissue.
biostars is not letting me paste tables in the comment. I'll try to add them to the post.
All in all, yes several of these samples are messed up, most probably each sample contained wildly different number and quality of nuclei from origin. That caused huge differences in the number of reads per sample.
However, my question still stands. How can I integrate Quality metrics/characteristics intrinsic from single nuclei, into cellbender or other alternatives, so it improves the tool's ability to call reliable non-empty, non-ambient droplets?
Did you load the same number of nuclei for each sample ?
Supposedly, yes. They checked under the microscope and could see single nuclei for each sample before FACS sorting them, which yielded ~10k nuclei/sample. We expected sequencing to get ~5k nuclei/sample.
Something clearly went wrong with the overall protocol across all samples. And something was very wrong with p4's samples. Yet still cellbender has trouble with the best sample in this data (p3_t2) whose UMI per droplet x num droplets plot has a very sharp "cut", which look very different from what I've seen in other data (10x single CELL in this case):
As you have received fastq files of 22Go and 220Mo while loading the same amount of nuclei and multiplexing, something went wrong in the library preparation, the nuclei sanity, the barcoding or the demultiplexing step, nothing to do with ambient RNA.
How are you running
cellranger multi? I've had similar issues withcellbenderwhen usingcellbenderwith theflexchemistry.cellbenderneeds to be run on a per-sample basis (see this Github issue https://github.com/broadinstitute/CellBender/issues/373cellranger multiwas run by the core facility that did the sequencing.By the
_cmdlinefile, they used:cellranger multi --id=pool1 --csv=[REDACTED]/ibraries_Pool1.csv --localcores=20I double-checked. The barcodes from each demultiplexed sample are unique within the pool, and match between pool1A-pool2A and such.
I ran
cellbenderin a per-sample basis, with the demultiplexed data.I also ran beforehand
cellranger counton the entire pools and used a couple of demultiplexing-by-genotyping tools (vireoandsouporcell). Both tools found 2 donors in the combined pool1 (patients #1 & #2), but could only find 1 patient (and/or hallucinations) in pool2 because the almost nonexistent number of cells from patient #4.