How to discover public receptors between different diseases?

[linqy-immune]

What kind of documentation would you like to request?

New tutorial for a specific use case

Please describe your idea or what you want to achieve

Hello! @vadimnazarov
I am very very excited about the wonderful newly released version and amazed by the powerful functions. Please receive my sincere gratitudes for your excellent work and efforts for updating and maintaining such helpful tools for immune repertoire researchers.
There are some questions and requests,

1. How should I load my own data (metadata and qecuence data)?
2. I have read the documentation carefully and find that the instructions are based on tumor immune repertoire analysis, so many analysis functions are integrated with single-cell RNA seq data. But for autoimmune disease immune repertoire research, not all researches include sc-RNA seq data, in such circumstance, how can we normalize the pure sc_TCR seq data?
3. For the sake of comparing sc_TCR seq data between AIDs and infectious diseases, whether can we analyze by immunarch 0.10.3? For data from differnt datasets, how can we normalize them?
4. In order to find cross-reactive CDR3 motif/sequence between different diseases, could I use Public receptor indices: quantifying repertoire overlap section? But I notice that public repertoire analysis - discover receptors which are suspiciously overabundant in one groups of repertoires in contrast to another groups of repertoires (e.g., control vs condition to search for TIL or antigen-specific receptors) is still in progress. Could you please make this function come out earlier?
   Looking forward to your reply.
   Thanks and best wishes,
   Linqy

[linqy-immune]

Hi,
I am sorry to bother you again as I encounter some other questions which I can't deal with.

1. Whether the single-cell processed data from MiXCR can be analyzed by immunarch 0.10? If not, how should I deal with them?
2. Is it feasible to load data from two sets of data pre-processed by two softwares in one time (maybe different make_default_preprocessing) ?

Looking forward to your reply!!!

[vadimnazarov]

Hey @linqy-immune

First of all, thank you very much for your kind words, I sincerely appreciate them. I moved this issue to immunarch as it is concerning data analysis rather than low-level details of ImmunData.

1. I have a tutorial on that: https://immunomind.github.io/docs/tutorials/single_cell/#4-reading-file-paths-from-metadata In short, you need a table with file paths to your files. What format do you use? Like, MiXCR, TRUST4, AIRR-C, etc.

2. immunarch / immundata works on immune repertoire data, meaning there is no difference whether it is from bulk or single-cell sequencing. The whole idea is to have a flexible tool that allows annotation of immune receptors, and it doesn't matter from what source those annotations came from.

3. What do you mean by "normalization" in that context?

4. Yes, sadly, it is still in progress – the statistics behind that are quite complex. I can try to accelerate things, could you please describe how your data looks like (like, how many repertoires, which groups, etc.) and what you want to do (e.g., find receptors which are significantly over-abundant in one repertoire / group in contrast to others).

5. Regarding the different datasets – this could be risky due to biological differences. Could you describe what datasets you have please? Currently, immundata does not support this, and I have no plans to support that, to be completely honest with you. But please describe your use case, and let me figure something.

[vadimnazarov]

If you used MiXCR, did you use the AIRR-C compliant output format? It is much better to do this

[linqy-immune]

Thank you for your urgent help! This means a lot to me!
There are two different formats of data which I wanna upload together, so I use the same code idata <- read_repertoires(path = inp_files, schema = schema, metadata = md_table, barcode_col = "barcode", locus_col = "locus", # Mind that we use locushere because of column renaming umi_col = "umis", preprocess = make_default_preprocessing("10x"), repertoire_schema = "Disease")
One is below, I download them from GEO website

The other is pre-processed by MiXCR from fastq data,

Afte I run the code above, one error comes out Error in agg_receptors(): ! Not all columns in the receptor schema present in the data: [cdr3, v_call, j_call, and locus]. Please double check and run again.
How should I do now?

[vadimnazarov]

Okay, so the first one is 10x format, the second one is MiXCR. Let's start by separately reading them to two different ImmunDatas.

10x:

 immundata::read_repertoires(
    path = YOUR_PATH or "<metadata>",
    schema = YOUR_SCHEMA,
    metadata = YOUR_METADATA,
    barcode_col = "barcode",
    locus_col = "locus",
    umi_col = "umis",
    preprocess = immundata::make_default_preprocessing("10x"),
    rename_columns = immundata::imd_rename_cols("10x")
  )

MiXCR:

immundata::read_repertoires(
    path = YOUR_PATH or "<metadata>",
    schema = YOUR_SCHEMA,
    metadata = YOUR_METADATA,
    barcode_code = "barcode",
    locus_col = "locus",
    umi_col = "umis",
    rename_columns = c(
      cdr3_aa = "aaSeqCDR3",
      cdr3_nt = "nSeqCDR3",
      v_call = "allVHitsWithScore",
      d_call = "allDHitsWithScore",
      j_call = "allJHitsWithScore",
      locus = "topChains",
      barcode = "cellGroup",
      umis = "uniqueMoleculeCount"
    )
  )
}

[vadimnazarov]

I would advise you to use this AIRR-C-compliant output for MiXCR https://mixcr.com/mixcr/reference/mixcr-exportAirr/

[linqy-immune]

Thank you so much!!! I would give them a try.

Hey @linqy-immune

First of all, thank you very much for your kind words, I sincerely appreciate them. I moved this issue to immunarch as it is concerning data analysis rather than low-level details of ImmunData.

1. I have a tutorial on that: https://immunomind.github.io/docs/tutorials/single_cell/#4-reading-file-paths-from-metadata In short, you need a table with file paths to your files. What format do you use? Like, MiXCR, TRUST4, AIRR-C, etc.
2. immunarch / immundata works on immune repertoire data, meaning there is no difference whether it is from bulk or single-cell sequencing. The whole idea is to have a flexible tool that allows annotation of immune receptors, and it doesn't matter from what source those annotations came from.
3. What do you mean by "normalization" in that context?
4. Yes, sadly, it is still in progress – the statistics behind that are quite complex. I can try to accelerate things, could you please describe how your data looks like (like, how many repertoires, which groups, etc.) and what you want to do (e.g., find receptors which are significantly over-abundant in one repertoire / group in contrast to others).
5. Regarding the different datasets – this could be risky due to biological differences. Could you describe what datasets you have please? Currently, immundata does not support this, and I have no plans to support that, to be completely honest with you. But please describe your use case, and let me figure something.

Aiming for finding receptors/clusters which are significantly cross-reactive between two diseases, it seems that I have to integrate different datastes? If convenient for you, could you give me some suggestions on achieving my goal? The datasets I have are just like above.
In the tutorial 1b. Number of receptors and barcodes for cell clusters, there is one sentence, "We can also normalize by the total number of barcodes per tissue for fairer cross-tissue comparison", whether is this applicable for data from different sets of data and even diseases? It is the first time I explore finding such cluster/receptor, could you please give me some suggestions and kind help.
Thank you again and best regards!!!

[linqy-immune]

Okay, so the first one is 10x format, the second one is MiXCR. Let's start by separately reading them to two different ImmunDatas.

10x:

immundata::read_repertoires(
path = YOUR_PATH or "",
schema = YOUR_SCHEMA,
metadata = YOUR_METADATA,
barcode_col = "barcode",
locus_col = "locus",
umi_col = "umis",
preprocess = immundata::make_default_preprocessing("10x"),
rename_columns = immundata::imd_rename_cols("10x")
)
MiXCR:

immundata::read_repertoires(
path = YOUR_PATH or "",
schema = YOUR_SCHEMA,
metadata = YOUR_METADATA,
barcode_code = "barcode",
locus_col = "locus",
umi_col = "umis",
rename_columns = c(
cdr3_aa = "aaSeqCDR3",
cdr3_nt = "nSeqCDR3",
v_call = "allVHitsWithScore",
d_call = "allDHitsWithScore",
j_call = "allJHitsWithScore",
locus = "topChains",
barcode = "cellGroup",
umis = "uniqueMoleculeCount"
)
)
}

It's a little pity that I tried this just now, but found it still unsuccessful (the same error).
Does this mean that I need to place my two kinds of datasets into two pathes, and load them seperately?

[vadimnazarov]

Yes, two different folders, two different metadatas

[vadimnazarov]

And if you encounter some issue, please copy and paste the error message, it helps a lot in the identification of root causes of bugs

[linqy-immune]

Hi,
I have tried to load the data seperately, but it is a little pity that there still exists such error when I run the code you provide above for 10X data:
Error in filter(): ! This operation cannot be carried out by DuckDB, and the input is a stingy duckplyr frame. ℹ Use compute(prudence = "lavish") to materialize to temporary storage and continue with duckplyr. ℹ See vignette("prudence") for other options. Caused by error in rel_filter(): ! {"exception_type":"Binder","exception_message":"Could not choose a best candidate function for the function call \"r_base::==(BOOLEAN, STRING_LITERAL)\". In order to select one, please add explicit type casts.\n\tCandidate functions:\n\tr_base::==(BOOLEAN, INTEGER) -> BOOLEAN\n\tr_base::==(BOOLEAN, DOUBLE) -> BOOLEAN\n\tr_base::==(BOOLEAN, BOOLEAN) -> BOOLEAN\n"} Run rlang::last_trace() to see where the error occurred. even after I revise the data to make them identical to the test data Ib6 on the website in format.
And for MiXCR data, after I run the command you provide above, the error Error in agg_receptors(): ! Not all columns in the receptor schema present in the data: [cdr3, v_call, j_call, and locus]. Please double check and run again. Run rlang::last_trace() to see where the error occurred. occurs again.
I am confused that after I run the two commands, whether the two sets of data can be merged or integrated in one immundata/repertoire? Could I analyze them together>
How should I do now? Could you please give me some instructions?
The two files are attached here.

MiXCR.csv
10X.csv

[linqy-immune]

Maybe something wrong happens in the machine cause an error ── Saving the newly created ImmunData to disk ℹ Writing the receptor annotation data to [/public/software/R/lib64/R/library/immundata/extdata/single_cell/immundata-lb6.csv/annotations.parquet] DuckDB progress: 0%Error in compute_parquet(): ! {"exception_type":"IO","exception_message":"Failed to write '/public/software/R/lib64/R/library/immundata/extdata/single_cell/immundata-lb6.csv/annotations.parquet': Cannot open file \"/public/software/R/lib64/R/library/immundata/extdata/single_cell/immundata-lb6.csv/tmp_annotations.parquet\": Permission denied","errno":"13"} Run rlang::last_trace()to see where the error occurred. occurs when I run the test-data.

[vadimnazarov]

@linqy-immune I will look into the formats, but the error about the writing is pretty clear – you don't have permissions to write a file on the disk, it's not on the immunarch side. read_repertoires has an output_folder parameter to specify where to write the, well, output repertoires: https://immunomind.github.io/immundata/reference/read_repertoires.html

The file MiXCR you sent is completely corrupted, and immunarch isn't going to help you – it's out of the scope for it. Please use MiXCR exportAIRR function to export to AIRR-C format. I don't have plans to support anything except AIRR-C format and maybe 10x – the reason being is pretty much illustrated here, which could be easily solved by tools fully embracing the AIRR-C format.

MiXCR exportAIRR – https://mixcr.com/mixcr/reference/mixcr-exportAirr/

AIRR-C format – https://docs.airr-community.org/en/stable/datarep/rearrangements.html

The 10x file you sent has different values for productive sequences. "TRUE" in your file vs. "True" or "true" in the specification by 10x genomics: https://www.10xgenomics.com/support/software/cell-ranger/7.2/analysis/outputs/cr-5p-outputs-annotations-vdj

I attached your MiXCR file, but with the "corrupted" parts removed. There is no functionality to merge files together, as merging is an extremely rare case. I'm thinking about how to implement it at the moment.

I will ping you here once I have some news on

1. reading the files in the formats you send me
2. merging the data together

MiXCR.csv