Mapping TCR info to individual cells (scRNA-seq dataset)

[malonzm1]

Hi,

I would like to map the mixcr TCR output to individual cells (using scRNA-seq data). That is, determine which cells contribute to the identified TCR repertoire. Is this possible?

Thanks and good day.

[mizraelson]

Hi. Do you mean mapping bulk TCR repertoire sequencing (assembled clones without cell barcodes) to scRNA-seq?

[malonzm1]

Thanks. I generated TCR repertoire using mixcr with scRNA-seq data. I would like to know which cells actually contain the clonotypes. Does that make sense?

[mizraelson]

What MiXCR command did you use?

[malonzm1]

mixcr analyze rna-seq -s hsa --rna 1.fastq.gz 2.fastq.gz outdir and
mixcr analyze rna-seq -s hsa --rna {{CELL0ROW:a}}_1.fastq.gz {{CELL0ROW:a}}_2.fastq.gz outdir
but the data is for scRNA-seq (I run into errors when using the 10x and Smart-seq2 presets).

[mizraelson]

This preset is designed for a bulk RNA seq data.what protocol did you use for scRNA seq? Was is it 10x? What error did you get?

[malonzm1]

I used 10x and Smart-seq2. For 10x, there's no preset for 3'.

[mizraelson]

What chemistry did you use for 10x v2 or v3. Both are available for 3' and they differ in cell barcode length.
Also, did you encounter any errors using 'smart-seq2'? Can you share the text of the error?

[malonzm1]

What are the presets for 10x v2 or v3 for 3'? Can I map the TCR clonotypes to the individual cells?

[mizraelson]

We originally didn't include the 3' preset because the clonotype yield from this type of data is significantly low.

If you're using 3'10x V2, the barcodes are identical to 5', so you can use the following preset:

mixcr analyze 10x-sc-5gex \
--species hsa \
input_R1.fastq.gz \
input_R2.fastq.gz \
output

For 3'10x V3, use the command below:

mixcr analyze 10x-sc-5gex \
--species hsa \
-MrefineTagsAndSort.whitelists.CELL=builtin:3M-febrary-2018 \
--tag-pattern "^(CELL:N{16})(UMI:N{12})\^(R2:*)" \
input_R1.fastq.gz \
input_R2.fastq.gz \
output

In the resulting output tables, there will be a column showing the CELL barcoded sequence. You can utilize this to map the cells to their respective expression data.

Do you require assistance with smart-seq2 data?

[malonzm1]

Thanks. In the above examples, does input_R2.fastq.gz contain the barcode+UMI sequence?

[mizraelson]

Yes, everything is according to 10x protocol.

[malonzm1]

Thanks. What if there are two coding fastq files?

[mizraelson]

Sorry, I mislead you, in 10x CELL and UMI barcodes are in R1 and the "coding" sequence is in R2.

Do you mean that you have a longer R1 which also covers some sequence in addition to barcodes?

[malonzm1]

Thanks. Yes.

[mizraelson]

Then you can use --tag-pattern "^(CELL:N{16})(UMI:N{10})(R1:*)\^(R2:*)" and --tag-pattern "^(CELL:N{16})(UMI:N{12})(R1:*)\^(R2:*)" for V2 and V3 respectively.

[malonzm1]

I read here https://www.biostars.org/p/9529864/#9572200 that in cases where R1 is barcode+UMI, the rest of the read can be effectively ignored. What tag pattern should I use in that case?

[mizraelson]

--tag-pattern "^(CELL:N{16})(UMI:N{10})\^(R2:*)" and --tag-pattern "^(CELL:N{16})(UMI:N{12})\^(R2:*)" for V2 and V3 respectively.

[malonzm1]

Thanks!

[mizraelson]

Hi,

1. Preset for Smartseq2 is smart-seq2-vdj. You can read about it here.
2. Use the original commands provided above (where the pattern doesn't have R1) to skip the rest of R1 read.
3. I couldn't find any info regarding 3' v1 on their website. If you can share the barcode structure and the whitelist I can help you with that.

[malonzm1]

3. I couldn't find any info regarding 3' v1 on their website. If you can share the barcode structure and the whitelist I can help you with that.

Thanks. https://www.biostars.org/p/462568/

10x v1
Whitelist, 737K-april-2014_rc.txt
CB length, 14
UMI start, 15
UMI length, 10 (courtesy ATpoint)

[malonzm1]

Hi,

How can I run mixcr with SureCell sequences?

SureCell (18 bp barcode, 8 bp UMI): surecell, ddseq, biorad

Thanks and good day.

[mizraelson]

Hi, are you certain you're using 3' v1? To my knowledge, that kit version was discontinued a long time ago.

From the link you provided, the tag pattern appears as follows:
--tag-pattern "^(CELL:N{14})(UMI:N{10})\^(R2:*)"

The corresponding command would be:

mixcr analyze 10x-sc-5gex \
--species hsa \
-MrefineTagsAndSort.whitelists.CELL=builtin:737K-april-2014_rc \
--tag-pattern "^(CELL:N{14})(UMI:N{10})\^(R2:*)" \
input_R1.fastq.gz \
input_R2.fastq.gz \
output

However, based on this link I found, the cell barcode is actually located in the Index read.

Given this, the command should be:

mixcr analyze 10x-sc-5gex \
--species hsa \
-MrefineTagsAndSort.whitelists.CELL=builtin:737K-april-2014_rc \
--tag-pattern "^(UMI:N{10})\^(R2:*)\^(CELL:N{14})" \
input_R1.fastq.gz \
input_R2.fastq.gz \
input_I1.fastq.gz \
output

Illumina Bio-Rad SureCell 3' WTA for ddSEQ

We haven't tried this protocol before. Since it's a 3' transcriptome, we don't anticipate a substantial number of TCR/BCR clones from it. But, based on the link shared above, the preset should look like this:

mixcr analyze generic-single-cell-gex-with-umi \
--species hsa \
--tag-pattern "^(CELL1:N{6})tagccatcgcattgc(CELL2:N{6})tacctctgagctgaa(CELL3:N{6})acg(UMI:N{8})\^(R2:*)" \
input_R1.fastq.gz \
input_R2.fastq.gz \
output

Sincerely,
Mark

[malonzm1]

Thanks!

Do you also have a preset for Fluidigm C1?

[mizraelson]

Can you share the protocol that you used fro cDNA library construction? To my knowledge Fluidigm C1 platform can be used with different chemistry.

If you can provide the library structure I will help you with the command. UMIs, CELL barcodes any technical sequences, are all cells in the single pair of FASTQ files or each cells' sequences are in a distinct pair of FSTQ files, etc.

[mizraelson]

Is it related to this #1105

[malonzm1]

Yes. Many thanks!

[malonzm1]

Hi again,

Re #498, can I use the following command:

mixcr analyze 10x-sc-5gex -s hsa --rna --threads 4 --force-overwrite <(cat SRR7159836/fastq/bamtofastq_R1.fastq.gz SRR7159837/fastq/bamtofastq_R1.fastq.gz) <(cat SRR7159836/fastq/bamtofastq_R2.fastq.gz SRR7159837/fastq/bamtofastq_R2.fastq.gz) GSM3140593/GSM3140593

[malonzm1]

Using the command above raises the following error:

sh: -c: line 0: syntax error near unexpected token `('
sh: -c: line 0: `mixcr analyze 10x-sc-5gex -s hsa --rna --threads 4 --force-overwrite <(cat SRR7159836/fastq/bamtofastq_R1.fastq.gz SRR7159837/fastq/bamtofastq_R1.fastq.gz) <(cat /SRR7159836/fastq/bamtofastq_R2.fastq.gz SRR7159837/fastq/bamtofastq_R2.fastq.gz) GSM3140593/GSM3140593'

Also when using zcat instead of cat. Please advise.

[mizraelson]

Is SRR7159837 and SRR7159836 the same biological sample?

[malonzm1]

Yes

[malonzm1]

I made a work-around by using softlinks and file name expansion.

[mizraelson]

Yes, the correct way would be to use wildcards like {{a}} to aggregate multiple samples. There is no need to 'cat' the file.
#1598

[malonzm1]

Hi,

I get the following Exception during refineTagsAndSort.

java.lang.IllegalStateException: Unexpected buffer behaviour. This might be a sign that you ran out of storage in the output or tmp folder.
        at com.milaboratory.o.Bt.completed(SourceFile:1582)
        at com.milaboratory.o.CO$a.completed(SourceFile:1163)
        at java.base/sun.nio.ch.Invoker.invokeUnchecked(Invoker.java:127)
        at java.base/sun.nio.ch.SimpleAsynchronousFileChannelImpl$3.run(SimpleAsynchronousFileChannelImpl.java:389)
        at java.base/java.util.concurrent.ForkJoinTask$RunnableExecuteAction.exec(ForkJoinTask.java:1426)
        at java.base/java.util.concurrent.ForkJoinTask.doExec(ForkJoinTask.java:290)
        at java.base/java.util.concurrent.ForkJoinPool$WorkQueue.topLevelExec(ForkJoinPool.java:1020)
        at java.base/java.util.concurrent.ForkJoinPool.scan(ForkJoinPool.java:1656)
        at java.base/java.util.concurrent.ForkJoinPool.runWorker(ForkJoinPool.java:1594)
        at java.base/java.util.concurrent.ForkJoinWorkerThread.run(ForkJoinWorkerThread.java:183)

I still have plenty of disk space. Please advise.

[mizraelson]

Hi, do you run in on a cluster?

[malonzm1]

Yes. Br, Get Outlook for Android<https://aka.ms/AAb9ysg>

…

________________________________ From: mizraelson ***@***.***> Sent: Wednesday, March 27, 2024 8:49:35 PM To: milaboratory/mixcr ***@***.***> Cc: Malonzo Maia ***@***.***>; Author ***@***.***> Subject: Re: [milaboratory/mixcr] Mapping TCR info to individual cells (scRNA-seq dataset) (Discussion #1412) Hi, do you run in on a cluster? — Reply to this email directly, view it on GitHub<#1412 (reply in thread)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/AHXYMJOAJNGUNHUYRYTMSE3Y2K557AVCNFSM6AAAAABFAUALNCVHI2DSMVQWIX3LMV43SRDJONRXK43TNFXW4Q3PNVWWK3TUHM4DSMRXHA4DQ>. You are receiving this because you authored the thread.Message ID: ***@***.***>

[malonzm1]

But there are also samples that completed with no problems.

[mizraelson]

Hi, so it seems like the issue is with the way analysis is run.
/scratch/work/malonzm1/.conda_envs/mixcr_env/bin/mixcr: line 196: 31502 Killed $java -Dmi.license.file="$dir/mi.license" -D${app}.path="$dir" -D${app}.command=${app} "${javaArgs[@]}" -jar "$jar" "${appArgs[@]}" This is a memory error which can happen when you try to allocate more memory than available in your system.

The other errors looks like a result of an analysis when you run the same sample several times in parallel writing the outputs to the same directory. Can you try running the problematic sample, the one that gave an error before, separately? And make sure not to allocate more memory than available.

[malonzm1]

Thanks. You're right. I mistakenly ran the same samples in parallel. I've corrected this now.

[malonzm1]

Hi, re getting "Killed", I allocated 480Gb (4 cpu with 120Gb per cpu) and use -Xmx300g but the process still gets killed. I'm now trying -Xmx200g.

[mizraelson]

Is it possible that you allocate more memory than available ?

[malonzm1]

I'm using SLURM, so I believe the SLURM job wouldn't run unless the memory requested is allocated. I think. Br, Get Outlook for Android<https://aka.ms/AAb9ysg>

…

________________________________ From: mizraelson ***@***.***> Sent: Monday, April 22, 2024 10:39:12 AM To: milaboratory/mixcr ***@***.***> Cc: Malonzo Maia ***@***.***>; Author ***@***.***> Subject: Re: [milaboratory/mixcr] Mapping TCR info to individual cells (scRNA-seq dataset) (Discussion #1412) Is it possible that you allocate more memory than available ? — Reply to this email directly, view it on GitHub<#1412 (reply in thread)>, or unsubscribe<https://github.com/notifications/unsubscribe-auth/AHXYMJI3AGFYSVZPD4J6WBTY6RZ5BAVCNFSM6AAAAABFAUALNCVHI2DSMVQWIX3LMV43SRDJONRXK43TNFXW4Q3PNVWWK3TUHM4TCOBTGA2DK>. You are receiving this because you authored the thread.Message ID: ***@***.***>

[mizraelson]

Cell grouping tells us about a reliable set of clones in a number of cells. 'Undefined' means there were not enough information across the dataset to assign a certain clone to a group. This might still be a valid cell, or maybe two cells in one well, but we can't really tell for sure.

[malonzm1]

Thanks! Why are there rows with the same barcodes with cellGroup in both "contamination" and a valid cellGroup?

[mizraelson]

If a certain clone has been found in more then 80% of cells, which already contain their specific clones, such clone is marked as a contamination. Thus, it is a characteristic of a clone, not a cell barcode.

[malonzm1]

Thanks. So "undefined" is also a characteristic of a clone? Is that why a single barcode can be associated with both "undefined" and a valid cellGroup?

Is this a clone "TRBV200 TRBJ2-600 TRBC2*00"?

[mizraelson]

Yes, undefined means we can't safely assign a clone to any cell.

[malonzm1]

It's hard to say without looking at the actual data and the expression levels. MiXCR only identifies the TCR/BCR regardless of the expression of other genes. For example, if you had some plasma cells (full of BCR RNA) burst, you might end up with other cells being contaminated with these sequences.

Could this be caused by, e.g., doublets?

[mizraelson]

That is possible too, yes. Also depends on the software you use for cell type annotation and the way it deals with doublets.

[malonzm1]

Thanks. From file *.clone.groups.tsv, there are columns IGH.primary.allVHitsWithScore IGH.second.allVHitsWithScore. What is the difference between primary and second?

Also, why can there be multiple entries, e.g. TRBV4-3*00(1260),TRBV4-2*00(1231)? So TRBV4-3*00(1260) is the highest hit, but why is there also TRBV4-2*00(1231)?

[mizraelson]

1. Both primary and secondary chains can be found in this cell group. Typically, the primary chain has a higher count, but the secondary chain is also present in substantial amounts.

2. The sequences of the V genes are very similar. Consequently, MiXCR often reports several hits. The number in the brackets represents the alignment score; a higher score indicates greater certainty regarding the identification of the specific gene.

[malonzm1]

Thanks. What is the difference between the information in clone.groups.tsv, clone.groups_TRAB\TRAG\IG.tsv and clone.groups_mixed.tsv?

[malonzm1]

Are the information in clone.groups_TRAB\TRAG\IG\mixed.tsv just a subset of clone.groups.tsv?

[mizraelson]

yes, those files are the subsets of the clone.groups.tsv

[malonzm1]

Is it ok to use different versions of MiXCR on different samples from the same dataset?

[mizraelson]

It is always better to use the same version for the same dataset.

[malonzm1]

Thanks. Can I assume that the clones found in a cellGroup are present in all cells associated with the cellGroup (via barcode)?

[mizraelson]

Yes