Confusing Clone Definition

[yuyuleung]

Report

Hello,

I have tried running Scirpy with the following code to define BCR clonotypes using my BCR data:
(Rules: Same V gene, same J gene, and 85% sequence similarity at the nucleotide level of the junction region)

ir.pp.ir_dist(mdata, metric="normalized_hamming", cutoff=15, sequence="nt", histogram=False) ir.tl.define_clonotype_clusters( mdata, sequence="nt", metric="normalized_hamming", receptor_arms=sirpy_receptor_arms, dual_ir="all", same_v_gene=True, same_j_gene=True, partitions="fastgreedy", key_added="clone_id_85_similarity", )

However, I obtained a clonotype containing sequences with different junction lengths (both amino acid and nucleotide) and also different V genes.

Here is the AIRR file of the incorrectly assigned clonotype (clone 24):
test.airr.tsv

I have also tried only on this subset of airr.

It is confusing that contigs with different V genes were not compared (as I verified in the distance matrix—the junctions of index 0 and 1 were not compared to those of index 2 and 3), yet they were assigned the same clone_id.

Thank you for your help in debugging this issue.
Yuyu

Versions

[grst]

Hi,

I think the key is the parameter

receptor_arms=sirpy_receptor_arms

which I'd assume was set to "any" in your example. In that case, cells will be grouped into a clonotype with others if any of the two receptor arms matches below the distance cutoff. See also tl.clonotype_clusters docs.

Setting it to "all" (the default), produces the following table with your example

[yuyuleung]

Hi, thank you so much for your answer. The receptor_arms was set as "heavy" on my side. I think, in that case, cells will be grouped only considering the IGH receptor arms? However, these cells were grouped together.

[yuyuleung]

Sorry, The receptor_arms was set as "VDJ" on my side

[grst]

Ok, I think you are onto something. I'll need some time to figure out what's going on though!

[yuyuleung]

Ok, thank you for taking time to check it. :)

[yuyuleung]

By the way, when I used scirpy to analysis, I encountered an inconsistent behavior when trying to define clonotypes using only heavy chain ("VDJ") data.

When receptor_arms is set to "VDJ" (or heavy chain only), the clonotype defining should work using only the heavy chain information. Light chain information should not be required if it is not being considered for the analysis. However, even when receptor_arms is explicitly set to "VDJ", the function returns an error unless the paired light chain data is also provided. This seems contradictory because the light chain data is presumably ignored under this setting.

It seems to be a bug in the data validation or clonotype assignment logic, or not?

Thank you again for considering it.

[grst]

Hi @yuyuleung,

I took another look at your data, and it seems that the VDJ sequences are indeed the same between all cells in the test dataset? They are just swapped between VDJ_1 and VDJ_2 for different cells, but it's always the same pair:

In that case it would be expected that they all get assigned the same clonotype cluster, even with dual_ir="all"

When receptor_arms is set to "VDJ" (or heavy chain only), the clonotype defining should work using only the heavy chain information. Light chain information should not be required if it is not being considered for the analysis. However, even when receptor_arms is explicitly set to "VDJ", the function returns an error unless the paired light chain data is also provided. This seems contradictory because the light chain data is presumably ignored under this setting.

What do you mean with "VJ" data is not provided? Just remove all IGK/IGL chains before loading the data into scirpy? I think it should still work in that case, but it might be just a too restrictive data validation. Could you share the error message please?

[yuyuleung]

Hi,

Oh, I have understood! Thank you for testing it.

However, when there are two different lengths of junctions of each cell, it will return an error when I plot the CDR3-motif with following code.

ir.pl.logoplot_cdr3_motif(
mdata,
chains= "VDJ_1",
to_type="information",
ax=ax
)

How can I work on this? Alternatively, could we ​implement​ a feature in the function ​to sort​ the VDJ1 and VDJ2 sequences by their length?

Regarding the second question, I will check it again and give you feedback.

Thank you so much again for taking time on them.

[grst]

However, when there are two different lengths of junctions of each cell, it will return an error when I plot the CDR3-motif with following code.

There's an example here how you can filter by length:
https://scirpy.scverse.org/en/latest/tutorials/tutorial_5k_bcr.html#logoplots

But this is good feedback, I think we should make it easier to make logoplots of sequences with different lengths, e.g. by automatically generating separate plots for all lengths in the input dataset.

[yuyuleung]

Thank you for your guidance, I will ​give it another try. :)