Question about UMI information in the output files vdjca and clns

[silvia1234567890]

Hello,

for the alignment, I ran this code:

mixcr align --preset generic-amplicon \
> --library Scophthalmus_maximus-IGH.json \
> --species Scophthalmus_maximus \
> --rna \
> --rigid-left-alignment-boundary \
> --rigid-right-alignment-boundary C \
> -OsaveOriginalReads=true \
> seqs.fasta \
> alignments.vdjca

for the clonotype:
mixcr assemble alignments.vdjca clones.clns

and to export:

mixcr exportAlignments alignments.vdjca alignments.tsv
mixcr exportClones clones.clns clones.tsv

In my clonotype table there is no "UMI count" column, and when I want to check the barcode coverage statistics with:
mixcr exportQc tags clones.clns tags.pdf
it gives me this message in the terminal:
No tag refinement report for clones.clns; did you run refineTagsAndSort command?

My question is: is there any parameter i forgot to put in the mixcr align code to get this information?

Thank you in advance!

[mizraelson]

Hi,

The commands you used do not include UMI parsing and processing, which is why that information is not appearing in the output. Could you please share the library structure and the protocol used for your data so I can assist you with the appropriate set of commands.

Sincerely,
Mark

[silvia1234567890]

Hi,
I send you as attached file the library that I built with MiXCR.

I pre-processed the sequences with pRESTO (Immcantation tool) and to do the alignment and assembly I used MiXCR with the commands I shared above.
The structure of my sequences is similar to the ones used in this pRESTO workflow (which is what I followed to pre-process my sequences): https://presto.readthedocs.io/en/stable/workflows/VanderHeiden2017_Workflow.html

I also tried to do the complete analysis with MiXCR using the raw sequences, and I used this code:

mixcr analyze generic-amplicon-with-umi \
> --library Scophthalmus_maximus-IGH.json \
> --species Scophthalmus_maximus \
> --rna \
> --tag-pattern "^(R1:*) \ ^(UMI:N{15})N{10}(R2:*)" \
> --rigid-left-alignment-boundary \
> --floating-right-alignment-boundary C \
> input_R1.fastq \
> input_R2.fastq \
> results

For this case (after do the whole analysis with MiXCR), when I want to run this code:
mixcr exportQc tags clones.clns tags.pdf
I get this message in the terminal:

App version: 4.5.0; built=Fri Sep 22 14:39:05 CEST 2023; rev=cdb24b4fb7; lib=repseqio.v3.0.1
Libraries have different checksums.

Could you help me please?

Thank you in advance!

[mizraelson]

I do recommend using MiXCR on raw data. Did you run analyze and exportQc on the same MiXCR version? Can you share the reports from that analyze command.

[silvia1234567890]

Hello,

yes, I'm using the last version of MiXCR (4.5.0) through anaconda. I am sending you as an attachment the reports obtained with the mixcr analyze command.
analyze reports.txt

One of the questions is if it is possible to get UMI information in the alignments/clones output tables having preprocessed the raw sequences with another tool and then using the mixcr align and mixcr assemble commands, as well as the info obtained with mixcr exportQc tags.
And other problem I have is that when I used the raw sequences with mixcr analyze I got the UMI information in the output tables but when I run the command mixcr exportQc tags it doesn't work...

Could you help me please? Thank you!

[mizraelson]

Hi,
If your data has already been preprocessed and collapsed by UMI using other tools, MiXCR unfortunately cannot provide UMI statistics, as UMIs are no longer present in the data. To obtain these statistics, the data must be processed with MiXCR from the start. Additionally, MiXCR has more advanced barcode processing algorithms compared to other tools, including pRESTO.

Regarding the issue in the assembly step in your reports:

Reads used in clonotypes, percent of total: 377943 (34,32%)
Reads dropped due to the lack of a clone sequence, percent of total: 398086 (36,15%)

A significant portion of reads (36%) does not cover the CDR3 region. This could be an issue with the identification of one of the CDR3 boundaries. To investigate this further, I recommend manually examining the alignments using the following command:

mixcr exportAlignmentsPretty -s 1000 -n 15 input.vdjca

This command will skip 1000 records and export the next 15 alignments. You can then see why CDR3 is not clearly identified in some cases.

Regarding the error you encountered:

App version: 4.5.0; built=Fri Sep 22 14:39:05 CEST 2023; rev=cdb24b4fb7; lib=repseqio.v3.0.1
Libraries have different checksums.

This issue could potentially be due to a file that has the same name as your library but contains slightly different content, which might be present in one of the MiXCR library directories. I suggest conducting the analysis again in a new, clean folder, starting from the reference preparation stage. Use a different name for the library to ensure there are no duplicates. If you can share one of the samples with us, I would be able to attempt reproducing the error and provide assistance in finding a solution.

Sincerely,
Mark

[silvia1234567890]

Thank you so much!!
Just other question: in some sequences, the segments V, D and J are detected but not CDR3 region.
For example, for the read_Id 1189,

When I ran the code:
mixcr exportAlignmentsPretty -s 1000 -n 15 input.vdjca

I get this information:

==============================================================================================

>>> Read ids: 1189

>>> Tags:
>>> UMI: GAAATGCAATTAGTG

>>> Target sequences (input sequences):

Sequence0:
Contains features: 

     0 ACAGCACGTCAGATTCAGCACAAACCATGTTTCCTGTAGCTCTGCTGCTGTTGGCAGCTGGATCATGTGTGAAGGGTGAA 79
       GGGGGGGGGGGGGGGGGGGGGGGGFGGFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGFGGGGGG

    80 CAGTTGACACAGCCGACCTCTGTGACTGTGCAGCCAGGTCAACGTCTGACCATCACCTGTCAGGTCTCTTATTCTCTTGG 159
       GGGGGGEGGGGGGGGGGGGGFFFGGGGGGFGGGGGGGGGGGGGGGGGGGGGGGFGGGFGGGGGFFGGGGGGGCGFG,@E9

   160 TAGCTACTA 168
       CGGG9999B

Sequence1:
Contains features: FR3, FR2_TO_FR4, CDR2, FR2, DRightPSegment, FR3_TO_FR4, DCDR3Part, 
VCDR3Part, CDR2_TO_FR4, VDJunction

     0 TACACAGCTTGGATCAGACAGCCTGCAGGAAAAGGACTGGAGTGGATTGGAATGAGATCTACTGGAGTTTCATACTACAA 79
       FFFFFFFEFEFFFFFE<;AFFA;FF@FFFGFFGFGGGGGGFGGGGCFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGG

    80 AGATTCATTAAAGAACAAGTTCAGTATCGACTTAGACACTTCCAGCAAAACTGTGACTCTAAATGGACAGAATGTGCAGC 159
       GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGDGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG

   160 CTGAAGACACTGCTGTGTATTACTGTGCCAGATACTGGGAGTACCTTGACTACTGGGGGAAAGGCACCATGGTCACCATC 239
       GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFF7GGGGGGGGGGGGGGGGGGGGGGGGG

   240 ACCTCGGCCACCCCAAAGGGACCAACTGTGTTTCCTCTGATGCCATGTGGTTCTGGAACT 299
       GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGCCCCC

>>> Gene features that can be extracted from this paired-read: 
FR3, FR2_TO_FR4, CDR2, FR2, DRightPSegment, FR3_TO_FR4, DCDR3Part, VCDR3Part, CDR2_TO_FR4, 
VDJunction

>>> Alignments with V gene:

IGHV3-18-F (total score = 1862.0)
Alignment of Sequence1 (score = 1862.0):
    101 TTCACAGCTTGGATCAGACAGCCTGCAGGAAAAGGACTGGAGTGGATTGGAATGAGATCTACTGGAGCTTCATACTACAA 180
        | ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||| ||||||||||||
      0 TACACAGCTTGGATCAGACAGCCTGCAGGAAAAGGACTGGAGTGGATTGGAATGAGATCTACTGGAGTTTCATACTACAA 79
        FFFFFFFEFEFFFFFE<;AFFA;FF@FFFGFFGFGGGGGGFGGGGCFGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGG

    181 AGATTCATTAAAGAACAAGTTCAGTATCGACTTAGACACTTCCAGCAAAACTGTGACTCTAAATGGACAGAATGTGCAGC 260
        ||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
     80 AGATTCATTAAAGAACAAGTTCAGTATCGACTTAGACACTTCCAGCAAAACTGTGACTCTAAATGGACAGAATGTGCAGC 159
        GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGFGGGGGGGGGGDGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG

    261 CTGAAGACACTGCTGTGTATTACTGTGCCAGA 292
        ||||||||||||||||||||||||||||||||
    160 CTGAAGACACTGCTGTGTATTACTGTGCCAGA 191
        GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG

>>> Alignments with D gene:

IGHDM1 (total score = 60.0)
Alignment of Sequence1 (score = 60.0):
     14 TACTGGGAGTAC 25
        ||||||||||||
    192 TACTGGGAGTAC 203
        GGGGGGGGGGGG

>>> Alignments with J gene:

IGHJ4 (total score = 354.0)
Alignment of Sequence1 (score = 354.0):
     23 TTGACTACTGGGGGAAAGGCACCATGGTCACCATCACCTCGGGTATG 69
        ||||||||||||||||||||||||||||||||||||||||||  |  
    205 TTGACTACTGGGGGAAAGGCACCATGGTCACCATCACCTCGGCCACC 251
        GGGGGGGFF7GGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGGG

>>> Alignments with C gene:

No hits.

==============================================================================================

Any thoughts about why CDR3 region is not covered here for example?

Thank you in advance

[mizraelson]

Can you please rerun that command with MiXCR v4.6 :

 mixcr exportAlignmentsPretty -s 1000 -n 15 input.vdjca

There was a bug in the output that we have fixed.

This should also answer the question about the clones with no CDR3. Most likely its due to the absence (not aligned) end of FR3 region.

[silvia1234567890]

Thank you, I just updated the MiXCR version.
I ran the command above with the last MiXCR version and I got this (for the same read_Id):

>>> Read ids: 1189

>>> Tags:
>>> UMI: GAAATGCAATTAGTG


                                                                                             
Quality   77777777777777777777777777777777777777777777777777777777777777777777777777777777   
Target0 0 ACAGCACGTCAGATTCAGCACAAACCATGTTTCCTGTAGCTCTGCTGCTGTTGGCAGCTGGATCATGTGTGAAGGGTGAA 79  Score

                                                                                               
Quality    77777777777777777777777777777777777777777777777777777777777777777777777767772674    
Target0 80 CAGTTGACACAGCCGACCTCTGTGACTGTGCAGCCAGGTCAACGTCTGACCATCACCTGTCAGGTCTCTTATTCTCTTGG 159  Score

                         
Quality     677744446    
Target0 160 TAGCTACTA 168  Score

               R1><FR2                                        FR2><CDR2           CDR2><FR3        
                Y  T  A  W  I  R  Q  P  A  G  K  G  L  E  W  I  G  M  R  S  T  G  V  S  Y  Y  K    
   Quality     77777777777777775567765776777777777777777777767777777777777777777777777777777777    
   Target1   0 TACACAGCTTGGATCAGACAGCCTGCAGGAAAAGGACTGGAGTGGATTGGAATGAGATCTACTGGAGTTTCATACTACAA 79   Score
IGHV3-18-F 101 tTcacagcttggatcagacagcctgcaggaaaaggactggagtggattggaatgagatctactggagCttcatactacaa 180  1862

                                                                                                   
                 D  S  L  K  N  K  F  S  I  D  L  D  T  S  S  K  T  V  T  L  N  G  Q  N  V  Q      
   Quality     77777777777777777777777777777777777777777777777777777777777777777777777777777777    
   Target1  80 AGATTCATTAAAGAACAAGTTCAGTATCGACTTAGACACTTCCAGCAAAACTGTGACTCTAAATGGACAGAATGTGCAGC 159  Score
IGHV3-18-F 181 agattcattaaagaacaagttcagtatcgacttagacacttccagcaaaactgtgactctaaatggacagaatgtgcagc 260  1862

                                                        DP>                                        
                                  FR3><CDR3  V><D      D><DP                                       
               P  E  D  T  A  V  Y  Y  C  A  R  Y  W  E  Y  L  D  Y  W  G  K  G  T  M  V  T  I     
   Quality     77777777777777777777777777777777777777777777777777777747777777777777777777777777    
   Target1 160 CTGAAGACACTGCTGTGTATTACTGTGCCAGATACTGGGAGTACCTTGACTACTGGGGGAAAGGCACCATGGTCACCATC 239  Score
IGHV3-18-F 261 ctgaagacactgctgtgtattactgtgccaga                                                 292  1862
    IGHDM1  14                                 tactgggagtac                                     25   60
     IGHJ4  23                                              ttgactactgggggaaaggcaccatggtcaccatc 57   354

                    FR4>                                                    
             T  S  A  T                                                     
Quality     777777777777777777777777777777777777777777777777777777766666    
Target1 240 ACCTCGGCCACCCCAAAGGGACCAACTGTGTTTCCTCTGATGCCATGTGGTTCTGGAACT 299  Score
  IGHJ4  58 acctcggGTaTG                                                 69   354

Is this okay?
Any change in any of the parameters I used with the prior version (build library, full analysis with MiXCR, alignments or assembly)?

[mizraelson]

It seems that the automatic algorithm responsible for constructing the reference library failed to identify the FR4 Begin for this J gene, which led to the CDR3 sequence not being recognized. Additionally, there are missing reference points for the IGHV6-22-F and IGHV4-1-F genes. This issue can occur when these genes significantly differ from those in other species. We will manually update the library file that you have created to address this and I'll provide you with a revised version that you can use.

No changes in the parameters are required after the update to the latest MiXCR version.

In the meantime, were you able to successfully execute the command mixcr exportQc tags clones.clns tags.pdf?

[mizraelson]

Do you need the fasta files with the segments I built the library with? or is it enough with the library? Sorry, I haven't executed the mixcr exportQc tags command yet, I will do it tomorrow (here in a few hours) and let you know if it worked. Thanks again for your help

After reviewing the library, it would be beneficial if you could share the initial FASTA files. Additionally, could you clarify which portion of the V genes is included? I noticed the use of VRegion parameter in buildLibrary, but in the library file, there are sequences that contain L2. If all the sequences are similar, it is advisable to use --v-gene-feature VExon2.

[silvia1234567890]

last questions (I promise).

just to be sure I'm understanding well the clonotype output table.

For example, in the picture attached, for the cloneId 0:

• "target sequence" is the sequence of this clone
• "readCount" means that 1182 reads contributed to this clone <-- would this be the same as saying "the number of different sequences with the same UMI"?
• "UMI count" means that 1143 UMIs (different UMIs) could have the "target sequence" of this clone?

Am I right?

Thank you in advance!

[mizraelson]

Hi,

The target sequence refers to the sequence of the region used for clone assembly, typically the CDR3 region, unless otherwise specified using the --assemble-clonotypes-by option.

The readCount represents the total number of reads contributing to the clone, while the UMI count indicates the number of unique UMI groups that have contributed. The close numbers between reads and UMIs suggest that, on average, each UMI group is represented by approximately one read.

Feel free to reach out if you have more questions :)

Best regards,
Mark

[silvia1234567890]

oh, I understand.

I would like to know if it's possible with mixcr analyze or another command to get the information below after pre-processing and before alignment steps.
I mean, after getting the unique sequences in the preprocessing step, for each sequence I would like to know its abundance. And with pRESTO, for each unique sequence (with its UMI) it gives me two parameters: CONSCOUNT (the number of reads with the same UMI) & DUPCOUNT (count of UMIs that has the same sequence).

Example (attached picture):

Is it possible to get this information with MiXCR?

I ask you this because I did the comparison of pRESTO (preprocessing) + MiXCR (align and clonotype) vs MiXCR (whole process, with mixcr analyze) and I got better results with mixcr analyze (more sequences etc).

Thank you

[mizraelson]

The thing is, MiXCR works differently, and that's why it gives better results. Initially, MiXCR aligns all reads to the reference sequence, extracts a UMI sequence for each read, and assigns it to the alignments. Later, during the refineTagsAndSort step, the UMI barcode sequences with PCR and sequencing errors are corrected, and automatic filtering filters out the UMI groups that are covered by a low number of reads (in your case, MiXCR most likely didn't filter any UMIs because most of them are only covered by one read, and by default, MiXCR tries to preserve at least 80% of the data). Then, alignments are assembled first into consensuses within each UMI group and later into clones across different UMIs. Such alignment-based UMI consensus assembly is one of the important advantages of the MiXCR pipeline.

Thus, consensuses are only formed from already aligned reads (alignments).

You can follow this process by exporting alignments with mixcr exportAlignments from different intermediate files:

• The 'sample.vdca' file contains reads successfully aligned to the reference, with a UMI barcode assigned to each alignment.
• The 'sample.refined.vdca' holds the alignments with corrected and filtered UMIs.
• The 'sample.clna' file (which is only available if you have specified -Massemble.clnaOutput=true) can show you which alignments with what UMIs were used in the assembly process of each clone. Don't forget to add -cloneId to the mixcr exportAlignments command to backtrack final clonotypes to the alignments. The alignments, where cloneId equals -1 were dropped due to various reasons (e.g. low sequencing quality etc.)

Sincerely,
Mark

[silvia1234567890]

thank you so much for such detailed explanations!! I now have a much better understanding of how MiXCR works.

Best regards,
Silvia