TRUST4 on long read

[Chuang1118]

Hi mourisl,

I need help with your your personal point of view plz.

This discussion refer to my previous question #141 The custom 5' spatial BCR has well performance in term of pipeline, thanks to TRUST4. but less sensible, maybe dataset doesn't has enough sequencing depth.

Now, I have the new version dataset by PromethION ONT.

Long read library like below:

Approximative Library Size = 24 (P7) + 10 (i7) + 34 (TruSeq Read2) + 300 (VDJ H ) + 1500 (C H ) + 30 (polyA) + 12 (UMI) + 16 (Spatial BC) + 33 (TruSeq Read1) + 10 (i5) + 29 (P5) = 1998 nt.

So, from now I have short and long reads from the same tissue/slide. I have the truth Spatial BC(barcodeWhitelist).
There is my first time do long read analysis.
You mentioned the PacBio HiFi data in #39 (comment)

Here, I want to know which point can fit to TRUST4.
e.g. Can you tell me which tools you've used for PacBio HiFI?
Which format as input to TRUST4? Bam with tagged BC and UMI correction ?

I just know UMI and BC are after polyA, but those positions are dynamics/changed by reads captured, I don't have any idea for this data processing.

Splitting single-end read( cDNA + BC&UMI) to 2 files cDNA read + BC&UMI reads like short read as input TRUST4?

Thanks a lot,
Chuang

[mourisl]

I tested TRUST4 on the raw PacBio Hifi data. I would suggest directly run "fastq-extractor" and "annotator" because for long reads there is no need for assembly. This thread could be useful: #63 . Though this simple pipeline could not handle barcodes, I think renaming the read id, e.g. "@[barcode]_[number]" and add the option "--barcode" to "annotator", could get around this. Then, you can use the other downstream scripts like "trust-barcoderep.pl", "trust-airr.pl" normally.
Nevertheless, ONT has a much higher error rate, so you may need to run some error corrections first.

For the input to TRUST4, the bams with corrected BC and UMI would be nice. As you mentioned the BC, UMI positions can be varied and also have a lot of sequencing errors. If some pipelines have already cleaned the data and put them in the BAM file, I would recommend using the BAM input.

Hope this will help.

[Chuang1118]

Thank you very much.
I will test it tomorrow and get back to you later.

[Chuang1118]

Hi mourisl,

There is the bam file with corrected BC and UMI.
I have used Minimap2 mapping on all chromosome 14 then tagged BC & UMI.
Now reads in bam are reverse complement rather than library read initial.
With deleting p7 i7 & p5 i5 adapter.
Final like below:

TruSeq Read1 + 16 (Spatial BC)  + 12 (UMI) + 30 (polyT) + cDNA(VDJ H + C H) +  TruSeq Read2

e.g. 2 reads :
CB: corrected barcode    UB: corrected UMI
CR: uncorrected barcode   UR: uncorrected UMI
CY: barcode quality value   UY: umi quality
GN: gene name
TR: trimmed sequence ?

read 1: GN: RNASE1

read 2: GN: IGHA2

All the reads have tagged above,
but some read with e.g. GN:Z:chr14_104646000_104647000 TR:Z:- or GN:Z:chr14_105547000_105548000 TR:Z:ENST00000390547

run-trust4 can handle my case?

I waiting for pass arguments tag CB & UB and start from annotation.

pseudo CMD like:

run-trust4 -f /opt2/TRUST4/hg38_bcrtcr.fa -t 4 --ref /opt2/TRUST4/human_IMGT+C.fa -b test_1k_chr14.tagged.bam \
--barcode -b CB \
--UMI -b UB \
--barcodeWhitelist barcodes.txt \
-o 1k_longRead --od test_out --repseq --stage 2 

Can you correct me ? Maybe I need also BamExtractor.cpp

I have try some test CMD.
if my cmd like this without BC & UMI information.

run-trust4 -f /opt2/TRUST4/hg38_bcrtcr.fa -t 4 --ref /opt2/TRUST4/human_IMGT+C.fa -b test_1k_chr14.tagged.bam \
--barcode Whitelistbarcodes.txt \
-o 1k_longRead --od test_out --repseq --stage 2 

message error
Need to use -a to specify the assembly file.

if I just run simple cmd

run-trust4 -f /opt2/TRUST4/hg38_bcrtcr.fa -t 4 --ref /opt2/TRUST4/human_IMGT+C.fa -b test_1k_chr14.tagged.bam \
--barcodeWhitelist barcodes.txt \
-o 1k_longRead --od test_out --repseq 

message error

[Fri Dec  9 12:21:53 2022] TRUST4 begins.
[Fri Dec  9 12:21:53 2022] SYSTEM CALL: /opt2/TRUST4/bam-extractor -b test_1k_chr14.tagged.bam -t 4 -f /opt2/TRUST4/hg38_bcrtcr.fa -o test_out/1k_longRead_toassemble 
[Fri Dec  9 12:21:53 2022] Start to extract candidate reads from bam file.
Unknown genome name: chr1
system /opt2/TRUST4/bam-extractor -b test_1k_chr14.tagged.bam -t 4 -f /opt2/TRUST4/hg38_bcrtcr.fa -o test_out/1k_longRead_toassemble  failed: 256 at /opt2/TRUST4/run-trust4 line 49.

why chr1 I don't have chr1 in my bam
proof

samtools view test_1k_chr14.tagged.bam | grep -w "chr1" | wc -l
0

here is your example.bam

your demo bam is a particular bam?

I don't know if i need this paramters in Annotator.cpp
--needReverseComplement: reverse complement sequences on another strand (default: no)
I knew all the constant region in anti strand DNA of chr4.
But I have some confusion with my read in bam and initial fastq read which read I need work on.

Thank you for your help.
Chuang

[Chuang1118]

Thank you for your response.

Parameter '-k' make sense to me, but I don't know how TRUST4 takes my tagged BC & UMI as argument.
Here I don't understand. In other word, as my input is bam, I don't have barcode+umi.fastq with --barcodeRange/--umiRange.
I'm in situation of single-spot. I need TRUST consider that BC & UMI.

--barcode STRING: if -b, bam field for barcode; if -1 -2/-u, file containing barcodes (defaul: not used)
--UMI STRING: if -b, bam field for UMI; if -1 -2/-u, file containing UMIs (default: not used)

I imagine cmd like this.

run-trust4 -f /opt2/TRUST4/hg38_bcrtcr.fa -t 4 --ref /opt2/TRUST4/human_IMGT+C.fa -b test_1k_chr14.tagged.bam \
-k 53 \
--barcode -b CB \ # here doesn't work.
--UMI -b UB \ # here doesn't work.
--barcodeWhitelist barcodes.txt \
-o 1k_longRead --od test_out --repseq 

Another question:
It's possible building only 1 VDJ H constant region fasta file with 1 gtf annoationf file?
Maybe it's illegal, I see your hg38_bcrtcr.fa with annotation in header.
Now my step mapping with all chr14 using Minimap2, it is not efficient.

[mourisl]

1. You put an extra -b in the --barcode/--UMI options. The command should look like:
   run-trust4 -f /opt2/TRUST4/hg38_bcrtcr.fa -t 4 --ref /opt2/TRUST4/human_IMGT+C.fa -b test_1k_chr14.tagged.bam
   -k 53
   --barcode CB \ # here doesn't work.
   --UMI UB \ # here doesn't work.
   --barcodeWhitelist barcodes.txt
   -o 1k_longRead --od test_out --repseq

2. I'm not sure I understand your question. You can use "grep --nogroupseparator -A 1 chr14 bcrtcr.fa > chr14_bartcr.fa" to get the chr14 only VDJH region.

[Chuang1118]

Thank you for your correction.
Now, I have used master TRUST4 version.
I have motified hg38_bcrtcr.fa to hg38_chr14_bcr.fa
and human_IMGT+C.fa to humanIGH_IMGT+C.fa [option]
to get only IGH

need delete extension .txt for get right format.
hg38_chr14_bcr.fa.txt
humanIGH_IMGT+C.fa.txt
test_5prime_5k_chr14.tagged.bam.txt

barcodes.txt

CMD

~/nice_env/TRUST4-master121222/run-trust4 -f ~/nice_env/TRUST4-master121222/hg38_chr14_bcr.fa -t 4 --ref ~/nice_env/TRUST4-master121222/humanIGH_IMGT+C.fa -b test_test_5prime_5k_chr14.tagged.bam \
-k 53 \
--barcode CB \
--UMI UB \
--barcodeWhitelist barcodes.txt \ #option
--repseq \
-o A_102_0023_longRead --od test_out3

without --barcode CB \ & --UMI UB \ it's work fine. but this isn't what I want.
I got error using cmd above:

[Mon Dec 12 12:45:25 2022] TRUST4 begins.
[Mon Dec 12 12:45:25 2022] SYSTEM CALL: /home/dong/nice_env/TRUST4-master121222/bam-extractor -b test_5prime_5k_chr14.tagged.bam -t 4 -f /home/dong/nice_env/TRUST4-master121222/hg38_chr14_bcr.fa -o TRUST_test_5prime_5k_chr14_toassemble  --barcode CB --UMI UB
[Mon Dec 12 12:45:25 2022] Start to extract candidate reads from bam file.
[Mon Dec 12 12:45:25 2022] Finish extracting reads.
[Mon Dec 12 12:45:25 2022] SYSTEM CALL: /home/dong/nice_env/TRUST4-master121222/trust4  -f /home/dong/nice_env/TRUST4-master121222/hg38_chr14_bcr.fa -t 4 -k 53 --trimLevel 2 --skipMateExtension -o TRUST_test_5prime_5k_chr14 -u TRUST_test_5prime_5k_chr14_toassemble.fq --barcode TRUST_test_5prime_5k_chr14_toassemble_bc.fa --UMI TRUST_test_5prime_5k_chr14_toassemble_umi.fa
[Mon Dec 12 12:45:25 2022] Found 596 reads.
[Mon Dec 12 12:45:25 2022] Finish sorting the reads.
[Mon Dec 12 12:45:25 2022] Finish rough annotations.
[Mon Dec 12 12:45:25 2022] Assembled 528 reads.
[Mon Dec 12 12:45:25 2022] Try to rescue 0 reads for assembly.
[Mon Dec 12 12:45:25 2022] Rescued 0 reads.
[Mon Dec 12 12:45:26 2022] SYSTEM CALL: /home/dong/nice_env/TRUST4-master121222/annotator -f /home/dong/nice_env/TRUST4-master121222/human_IMGT+C.fa -a TRUST_test_5prime_5k_chr14_final.out -t 4 -o TRUST_test_5prime_5k_chr14 --barcode --UMI -r TRUST_test_5prime_5k_chr14_assembled_reads.fa --airrAlignment > TRUST_test_5prime_5k_chr14_annot.fa
[Mon Dec 12 12:45:26 2022] Start to annotate assemblies.
[Mon Dec 12 12:45:27 2022] Start to realign reads for CDR3 analysis.
 has not barcode fieldsystem /home/dong/nice_env/TRUST4-master121222/annotator -f /home/dong/nice_env/TRUST4-master121222/human_IMGT+C.fa -a TRUST_test_5prime_5k_chr14_final.out -t 4 -o TRUST_test_5prime_5k_chr14 --barcode --UMI -r TRUST_test_5prime_5k_chr14_assembled_reads.fa --airrAlignment > TRUST_test_5prime_5k_chr14_annot.fa failed: 256 at /home/dong/nice_env/TRUST4-master121222/run-trust4 line 51.

for has not barcode fieldsystem I have checked TRUST_test_5prime_5k_chr14_assembled_reads.fa beacuse

TRUST4/Annotator.cpp

Line 707 in 8e5515d

if ( hasBarcode )

cat A_102_0023_longRead_assembled_reads.fa | grep "barcode" | head

>0f956c53-0298-4c3f-b892-b382aabd4e7a_0 1 3 39 barcode:GGACTAAGTCAGGAGT umi:353
>7ac0bb44-93c1-48cf-9cab-39e8b19813ef_0 1 2 8 barcode:TGATCCCAGCATTAGT umi:492
>268948dd-67c8-4b04-bafb-290e29958343_0 1 2 6 barcode:CTAAACGGGTGTAATC umi:398
>ffdfbd89-731a-4876-9cc9-850d7ca45664_0 1 2 2 barcode:TCTGGGTAGCGCTCAT umi:449
>be38587d-1a7a-4a4d-bbb0-9201bd01c26d_0 1 2 2 barcode:TCTGGGTAGCGCTCAT umi:449
>bd8d792d-fb29-4064-8a5d-27fc9396bb26_0 1 1 119 barcode:GAAACTCGTGCGATGC umi:20
>c5797b53-71b3-40d3-ab39-caaa667d52b0_0 1 1 118 barcode:GCTGGCAGGTGCCGTG umi:165
>0b3e452b-4729-47dd-b2aa-b0947b8f53bc_0 1 1 118 barcode:TCGGGATTCAAACATA umi:85
>58394321-9941-4975-be56-af38055ef3d7_0 1 1 117 barcode:AAACCGGGTAGGTACC umi:135
>f44b37d3-3fcd-43f8-9abd-5dbc001de758_0 1 1 117 barcode:TTACCATTGATTACCC umi:188
>fb596176-4dad-440d-bfcb-50de1c476446_0 1 1 117 barcode:CCCTGACTAACAAATT umi:172
>6b3333c4-8d63-4959-b5b0-3f404ea1ef99_0 1 1 117 barcode:GGCTTGGCTCTCACCT umi:196
>a8f67584-0867-4ed2-b4c8-8e33ef5e2eb1_0 1 1 117 barcode:GCACACGCCCATGGTC umi:50

All the header, I have barcode:
need delete extension .txt for get right format.
A_102_0023_longRead_assembled_reads.fa.txt
A_102_0023_longRead_final.out.txt

This error is permanence both local or singularity version.

if I used trust4 v1.0.7 only for this step
cmd below is your SYSTEM CALL: with trust4 v1.0.7

singularity exec -B /mnt/DOSI:/mnt/DOSI ccbr_trust4_1.0.7b.sif \
annotator -f /home/dong/nice_env/TRUST4-master121222/human_IMGT+C.fa -a test_out3/A_102_0023_longRead_final.out -t 4 -o test_out3/A_102_0023_longRead --barcode --UMI -r test_out3/A_102_0023_longRead_assembled_reads.fa --airrAlignment > A_102_0023_longRead_annot_old.fa

It's work

[Mon Dec 12 15:01:37 2022] Start to annotate assemblies.
[Mon Dec 12 15:01:41 2022] Start to realign reads for CDR3 analysis.
[Mon Dec 12 15:01:41 2022] Compute CDR3 abundance.
[Mon Dec 12 15:01:41 2022] Finish annotation.

Since you've motified this file Annotator.cpp [13 days ago], now this doesn't work for me.

My previous question 2:
I knew how to extract IGH in fasta file, but I don't how to build gtf file.
Since I just interested a small part of chr14 IGH constant region e.g.

TRUST4/human_IMGT+C.fa

Line 6416 in 8e5515d

>IGHG4 chr14 105620506 105626066 -

One giving IGHA2 have a multiple positions.
Here I just want to the acceleration of the process mapping by minimap2 with extract chr14 IGH constant region and a file gtf corrsponding. Any suggestion ?

Thanks a lot for your help.
Chuang

[mourisl]

1. Thank you for the detailed information. I think this is a bug I introduced a few days ago to remove a constraint on read length. I just fixed this issue and pushed it to the master branch. Could you please give it a try?
2. I'm not sure I get your question. Do you mean you want to extract the constant gene sequence based on a GTF file? You can run the BuildDatabaseFa.pl file. The command can be:

perl BuildDatabaseFa.pl reference_genome.fa genome_annotation.gtf file.txt > IGHC.fa

file.txt has the gene names of the constant gene: IGHM, IGHD, and so forth.

[Chuang1118]

I confirm TRUST4 master version work fine!
I really appreciate your help and your kind answer.

For my second question, since it works, just mapping on fully chr14 being less efficient.
In other words, IGH constant and variable region in a small part of chr14. I will do it later.
I have tried mapping only IGH constant region, that's why I want to gtf which has only IGH constant region.
Yesterday, I have checked the distribution of length in long read, it has one pic around 800 bases and we expected get the library around 2,000 bases. Maybe mapping not only constant region, but also variable region is a good strategy.

I will ask later if I had some issue for running TRUST4.

Thank you again,
Chuang

[Chuang1118]

Hi mourisl,

TRUST4 worked fine with my test dataset (head 4k), but I have some issues about my true dataset.

CMD

~/nice_env/TRUST4-master131222/run-trust4 -f ~/nice_env/TRUST4-master131222/hg38_chr14_bcr.fa \
-t 24 \
--ref ~/nice_env/TRUST4-master131222/humanIGH_IMGT+C.fa -b A_102_0023_chr14.tagged.bam \
-k 53 \
--barcode CB \
--UMI UB \
--barcodeWhitelist barcodes_A.txt \
--repseq \
-o A_102_0023_longRead --od trust4 

I have 4 libraries, named A , B, C, D.

Question 1:
Now, C has finished, according to log file.
nohup_C.out.txt

In my opinion, TRUST4 doesn't processe all the reads. Since Found 3137630 reads, Processed 1400000 reads.

[Mon Dec 19 12:58:55 2022] Read in and count kmers for 2900000 reads.
[Mon Dec 19 12:59:14 2022] Read in and count kmers for 3000000 reads.
[Mon Dec 19 12:59:30 2022] Read in and count kmers for 3100000 reads.
[Mon Dec 19 13:37:09 2022] Found 3137630 reads.
[Mon Dec 19 13:38:18 2022] Finish sorting the reads.
[Mon Dec 19 13:45:34 2022] Finish rough annotations.
[Mon Dec 19 13:48:57 2022] Processed 100000 reads (99656 are used for assembly).
[Mon Dec 19 13:57:14 2022] Processed 200000 reads (199143 are used for assembly).
[Mon Dec 19 14:15:53 2022] Processed 300000 reads (298632 are used for assembly).
[Mon Dec 19 14:49:18 2022] Processed 400000 reads (396027 are used for assembly).
[Mon Dec 19 15:40:29 2022] Processed 500000 reads (495769 are used for assembly).
[Mon Dec 19 16:45:48 2022] Processed 600000 reads (595611 are used for assembly).
[Mon Dec 19 18:05:36 2022] Processed 700000 reads (695479 are used for assembly).
[Mon Dec 19 19:11:43 2022] Processed 800000 reads (795326 are used for assembly).
[Mon Dec 19 19:52:53 2022] Processed 900000 reads (895024 are used for assembly).
[Mon Dec 19 20:12:37 2022] Processed 1000000 reads (993744 are used for assembly).
[Mon Dec 19 20:40:54 2022] Processed 1100000 reads (1093530 are used for assembly).
[Mon Dec 19 21:32:32 2022] Processed 1200000 reads (1192878 are used for assembly).
[Mon Dec 19 22:27:30 2022] Processed 1300000 reads (1292568 are used for assembly).
[Mon Dec 19 22:51:58 2022] Processed 1400000 reads (1392090 are used for assembly).
[Mon Dec 19 22:55:49 2022] Assembled 1434448 reads.
[Mon Dec 19 22:55:49 2022] Try to rescue 0 reads for assembly.
[Mon Dec 19 22:55:57 2022] Rescued 0 reads.

e.g. output C assembled_reads.fa

head C*_longRead_assembled_reads.fa
>a90f8e84-2306-41b7-a3fb-91e7960e0cfa_0 1 648 376629 barcode:CATTCGTCGTAGCGGT umi:72753
CTGGGGGGTCCCTGAGACTCTCCTGTGAAGTCTCTGGACTTACCTTCAATATGTATACCATGCACTGGGTCCGGCAGGCTACAGGAAAGGGGCTGGAGTGGGTCGCGTCGATTACAAGTGGTAGTCATTATATCTTCTATGCAGACTCATTGAGGGGCCGATTCACCATCTCCAGAGACAACGCCAGGAATTTACTCTACCTTCAGATGACCAACCTGAGAGCCGAGGACACATCTGTGTATTACTGTGTGAGGGACTTTACAAATGTGGCCTCGCTTGACTACTGGGGCCAGGGAACCCCGGTCACCGTCTCCTCAGCATCCCCGACCAGCCCCAAGGTCTTCCCGCTGAGCCTCTGCAGCACCCAGCCAGATGGGAACGTGGTCATCGCCTGCCTGGTCCAGGGCTTCTTCCCCCAGGAGCCACTCAGTGTGACCTGGAGCGAAAGCGGACAGGGCGTGACCGCCAGAAACTTCCCACCCAGCCAGGATGCCTCCGGGGACCTGTACACCACGAGCAGCCAGCTGACCCTGCCGGCCACACAGTGCCTAGCCGGCAAGTCCGTGACATGCCACGTGAAGCACTACACGAATCCCAGCCAGGATGTGACTGTGCCCTGCCCAGTTCCCTCAACTCCACCTACCCCATCTCCCTCAACTCCACCTACCCCATCTCCCTCATGCTGCCACCCCCGACTGTCACTGCACCGACCGGCCCTCGAGGACCTGCTCTTAGGTTCAGAAGCGAACCTCACGTGCACACTGACCGGCCTGAGAGATGCCTCAGGTGTCACCTTCACCTGGACGCCCTCAAGTGGGAAGAGCGCTGTTCAAGGACCACCTGAGCGTGACCTCTGTGGCTGCTACAGCGTGTCCAGTGTCCTGCCGGGCTGTGCCGAGCCATGGAACCATGG
>93550f04-ce6b-478c-9b17-50fae669473a_0 1 648 373577 barcode:CATTCGTCGTAGCGGT umi:72753
CTGGGGGGTCCCTGAGACTCTCCTGTGAAGTCTCTGGACTTACCTTCAATATGTATACCATGCACTGGGTCCGGCAGGCTACAGGAAAGGGGCTGGAGTGGGTCGCGTCGATTACAAGTGGTAGTCATTATATCTTCTATGCAGACTCATTGAGGGGCCGATTCACCATCTCCAGAGACAACGCCAGGAATTTACTCTACCTTCAGATGACCAACCTGAGAGCCGAGGACACATCTGTGTATTACTGTGTGAGGGACTTTACAAATGTGGCCTCGCTTGACTACTGGGGCCAGGGAACCCCGGTCACCGTCTCCTCAGCATCCCCGACCAGCCCCAGGTCTTCCCGCTGAGCCTCTGCAGCACCCAGCCAGATGGGAACGTGGTCATCGCCTGCCTGGTCCAGGGCTTCTTCCCCCAGGAGCCACTCAGTGTGACCTGGAGCGAAAGCGGACAGGGCGTGACCGCCAGAAACTTCCCACCCAGCCAGGATGCCTCCGGGGACCTGTACACCACGAGCAGCCAGCTGACCCTGCCGGCCACACAGTGCCTAGCCGGCAAGTCCGTGACATGCCACGTGAAGCACTACACGAATCCCAGCCAGGATGTGACTGTGCCCTGCCCAGTTCCCTCAACTCCACCTACCCCATCTCCCTCAACTCCACCTACCCCATCTCCCTCATGCTGCCACCCCCGACTGTCACTGCACCGACCGGCCCTCGAGGACCTGCTCTTAGGTTCAGAAGCGAACCTCACGTGCACACTGACCGGCCTGAGAGATGCCTCAGGTGTCACCTTCACCTGGACGCCCTCAAGTGGGAAGAGCGCTGTTCAAGGACCACCTGAGCGTGACCTCTGTGGCTGCTACAGCGTGTCCAGTGTCCTGCCGGGCTGTGCCGAGCCATGGAACCATGG
>45e8287a-6789-416c-81a2-47e3fc768478_0 1 648 364817 barcode:CATTCGTCGTAGCGGT umi:72753
CTGGGGGGTCCCTGAGACTCTCCTGTGAAGTCTCTGGACTTACCTTCAATATGTATACCATGCACTGGGTCCGGCAGGCTACAGGAAAGGGGCTGGAGTGGGTCGCGTCGATTACAAGTGGTAGTCATTATATCTTCTATGCAGACTCATTGAGGGGCCGATTCACCATCTCCAGAGACAACGCCAGGAATTTACTCTACCTTCAGATGACCAACCTGAGAGCCGAGGACACATCTGTGTATTACTGTGTGAGGGACTTTACAAATGTGGCCTCGCTTGACTACTGGGGCCAGGGAACCCCGGTCACCGTCTCCTCAGCATCCCCGACCAGCCCCAAGGTCTTCCCGCTGAGCCTCTGCAGCACCCAGCCAGATGGAAACGTGGTCATCGCCTGCCTGGTCCAGGGCTTCTTCCCCCAGGAGCCACTCAGTGTGACCTGGAGCGAAAGCGGACAGGGCGTGACCGCCAGAAACTTCCCACCCAGCCAGGATGCCTCCGGGACCTGTACACCACGAGCAGCCAGCTGACCCTGCCGGCCACACAGTGCCTAGCCGGCAAGTCCGTGACATGCCACGTGAAGCACTACACGAATCCCAGCCAGGATGTGACTGTGCCCTGCCCAGTTCCCTCAACTCCACCTACCCCATCTCCCTCAACTCCACCTACCCCATCTCCCTCATGCTGCCACCCCCGACTGTCACTGCACCGACCGGCCCTCGAGGACCTGCTCTTAGGTTCAGAAGCGAACCTCACGTGCACACTGACCGGCCTGAGAGATGCCTCAGGTGTCACCTTCACCTGGACGCCCTCAAGTGGGAAGAGCGCTGTTCAAGGACCACCTGAGCGTGACCTCTGTGGCTGCTACAGCGTGTCCAGTGTCCTGCCGGGCTGTGCCGAGCCATGGAACCATGG
>10eb0e8b-6d06-4af1-b413-00c54ec71ac4_0 1 614 373057 barcode:CATTCGTCGTAGCGGT umi:72753
CTGGGGGGTCCCTGAGACTCTCCTGTGAAGTCTCTGGACTTACCTTCAATATGTATACCATGCACTGGGTCCGGCAGGCTACAGGAAAGGGGCTGGAGTGGGTCGCGTCGATTACAAGTGGTAGTCATTATATCTTCTATGCAGACTCATTGAGGGGCCGATTCACCATCTCCAGAGACAACGCCAGGAATTTACTCTACCTTCAGATGACCAACCTGAGAGCCGAGGACACATCTGTGTATTACTGTGTGAGGGACTTTACAAATGTGGCCTCGCTTGACTACTGGGGCCAGGGAACCCCGGTCACCGTCTCCTCAGCATCCCCGACCAGCCCCAAGGTCTTCCCGCTGAGCCTCTGCAGCACCCAGCCAGATGGGAACGTGGTCATCGCCTGCCTGGTCCAGGGCTTCTTCCCCCAGGAGCCACTCAGTGTGACCTGGAGCGAAAGCGGACAGGGCGTGACCGCCAGAAACTTCCCACCCAGCCAGATCCTCCGGGGACCTGTACACCACGAGCAGCCAGCTGACCCTGCCGGCCACACAGTGCCTAGCCGGCAAGTCCGTGACATGCCACGTGAAGCACTACACGAATCCCAGCCAGGATGTGACTGTGCCCTGCCCAGTTCCCTCAACTCCACCTACCCCATCTCCCTCAACTCCACCTACCCCATCTCCCTCATGCTGCCACCCCCGACTGTCACTGCACCGACCGGCCCTCGAGGACCTGCTCTTAGGTTCAGAAGCGAACCTCACGTGCACACTGACCGGCCTGAGAGATGCCTCAGGTGTCACCTTCACCTGGACGCCCTCAAGTGGGAAGAGCGCTGTTCAAGGACCACCTGAGCGTGACCTCTGTGGCTGCTACAGCGTGTCCAGTGTCCTGCCGGGCTGTGCCGAGCCATGGAACCATGG
>a706be34-6dac-4ffe-8302-7841d4908936_0 1 547 4250 barcode:CCATAGGTTGGCGTGG umi:79532
CGGGGAGTCTCTGAAGATCTCCTGTAAGGGTTTTGGATACATTTTTACCAGTTACTGGATCGGCTGGGTGCGCCAGATGCCCGGTCAAGGCCTGGAGTGGATGGGGATCATCTATCCTAGTGACTCTGATGCCCGATACAGCCCGTCCTTCCAAGGCCAGGTCACCATGTCAGCCGACAAGGCCATCAGCACCGCCTACCTGCACTGGACCAGCCTGAAGGCCTCGGACACCGCCATGTATTACTGTGCGAGACATCCCTTATCCCCGTGTAGTTCTAATCATTGTTTCAGAGACTACCAGCCGGGGAGATTCGACTCCTGGGGCCAGGGGACCCTGGTCACCGTCTCCTCAGCACCCACCAAGGCTCCGGATGTGTTCCCCATCATATCAGGGTGCAGACACCCAAAGGATAACAGCCCTGTGGTCCTGGCATGCTTGATAACTGGGTACCACCCAACGTCCGTGACTGTCACCTGGTACATGGGGACACAGAGCCAGCCCCAGAGAACCTTCCCTGAGATACAAAGACGGGACAGCTACTACATGACAAGCAGCCAGCTCTCCACCCCCCTCCAGCAGTGGCGCCAAGGCGAGTACAAATGCGTGGTCCAGCACACCGCCAGCAAGAGTAAGAAGGAGATCTTCCGCTGGCCAGAGTCTCCAAAGGCACAGGCCTCCTCAGTGCCCACTGCACAACCCCAAGCAGAGGGCAGCCTCGCCAAGGCAACCACAGCCCCAGCCACCACCCGTAACACAGGAAGAGGAGGAGAAGAGAAGAAGAAGGAGAAGGAGAAAGAGGAACAAGAAGAGAGAGAGACAAAGACACCAG

tail C*_longRead_assembled_reads.fa
>47f6c994-26a4-4a0d-ac58-8c56065c1469_0 1 1 1 barcode:CTGGCTGGTTGTCAGT umi:70557
CCTTCAGTAGCTTTTGAATGAGCCCAGTTCACCAGTCTGCAGGCAAGGGGCTGGAGTGAGTAATAGATATAAAAGATGATGGAAAAGTTAGATACACCATGCATACTTGTGAAGCAGATTTCTCATCTCTCAAAGACAATGCTAAAGACTCTCTGGATCCTAAATCAACAGTCAGAGAACTGAGGACATGGCTGTGTATGGAGCTGTACACAAGGTTCTGAAGTGAGGAAATATTGTGCCAATCCATACAATACTTCTGCAGAGAAAGAAGGAGTCTGGGCTAAGGTACTGAGACTGACAAACAGGTGCAGCCCCAAGGCAGGTGCAGATGGAGGGAGGGTAAGGATGCTTTTCCTTCAGGATCTGTGGTTTTCCTTGCTTCTAACAATTCCCCTTGAGCCTCTGACATTTATATTTTGTGCTGACTATGAAGGCTGCTATGGATCAGAAAACTAATTTGAAAGAGGAAAATATTCTCATATGTCCTAAAAACAGATGTAAGTATTGGAGGCAGAACCAGAGGAGCCGGGTGTGTCTGTAGACACTGCCACCCCCAATGCCAGACCCACAACTAGTGCTGGGAAGGTGGGAGTTCTATGGGAGCTTCCCTGACGACCCCGTGGTCCAAGCTAAGTCCAGCAAGGCCATTGGTGCCTCGCTGAGCACAGTTGTCCATCAGGGATCTCCCACGTGTCCCAGCAGCAGCCACGCCCAACCTCCCACTGTGCACAGCCATCGTCTGGGAGGAGCTCCCAGGATGGGTGTCTTTGGCACACACAGTTGATGGGTGTTAGAGTGCAGTGCAGCAGCTGGCTGCTGGGCCTCATCAGGCTCTGCACACTGGAGGGATTCGAGGTGATTCTCAGGGCCAGCATGCTTCTACGGTTTCCTTATATGAAAAACCATGTGATTTAATTCATTCTCTAGATGACATAGATAATTATAATACAACTGCAAACAAAAAAAAAAAAAAAAATACTGACAACCAGCCAGAGATCGGAAGAGTGTTGTGTA
>cf4ec41b-c7ff-4fe1-a531-46d657467846_0 1 1 1 barcode:GTACGTTTGCCCGTCA umi:2104
CTTCTGGATACAGCTTCACTGCCTACTTTATCCACTGGGGCTCTGCCACACCCACAATAAGTTTCTGTGGCATGGGTTGGGTCAACCCTCAAAGTGGTGACACAAACTATGCACAGAAGTTTCAGGGGCAGGGCCACCATTGGGACATGTCCATCACCACAGCCTTCATGAGCTGAGCAGAGAGCCTATCCCACACACACCCGGCCGTCTATTACTGTGCGAGAGAACCATCTCAAATCATAGATGTCCGCAATGATGCTTTTGTATGTCTGGGGTGCAAGGGACACTGGGTCTCGCCGTCTCTTTAGATCTCCCCGACCAGCCCCAAGGTCTTCCCCGCTGAGCCTCGACAGCACCCCCCAAGACTGGGGAATGTTTTTGTTGTTGCATGCCTGGTCCAGGGTCTTCTCCCGTGTTCCTGAGCCACTCAGTGTCTACCCCTGTGAAGCCAGCAGAACGTGACCGCCAGGAAATTTTCCCCTACCAGCCTAGGATGCCTCCGGGGACCCTAATATCACCAGCAGCCAGTGGACCCTGCCGGCCACAGAGTGCCCAGATTCTGAAGTCCGGTGACATGTCATGTGAAACTCATACTTGAATCCCAGCCAGGCATACTGTGCCCTGCCCATTTTTCCCCCACCTCCCCCACGCCTCACCCCCTGACTTGTCGCTGCACGTGGCCCTCCAATGACCTGCTCTTAGGTTCAGAAGCGAACCTCAGAGAGATGCCTTTCTGGTGGCTGGCTTCACGCCTGGTACTGCTCCAAACTGGAAGAGCGCTGTTCAAGGACCACCTGAGCGTGACCTCTGCGGCATTTACAGCGTGTCCAGTGTCCTGCCTGGCTG
>7d908e84-d809-4e17-bb81-15abd15ebf39_0 1 1 1 barcode:ATAGAGTACTGGGACA umi:83878
AGGCTTCTGCAGACACCTTCACCGACTACTTTATGAACTGGATGTGTGAGGCCCCTGGACAAAGGCTTGAGTGGATGGATGGATCAACGCTGGCAATGGTAACAGAAACAATATTCACAGAAGCTCCAGGGCAGAGCTTACTATTACCAGGGACACACTGCGAGACAGCCTACATGGAGCTGAGAGGCAGCTGACTGAGGAAGACGGCCACATGACTGTGCGAGACACAGAGTGAAAACCCACATCCTGAGAGTGTGAAGGAACCCCAGGGAGGAAGCAGCTGTACTGGCATGGAGGAAATGGACAAGATTATTATATTTATACTCTCCAACAAAAACTACATTAAGTCATTAAGGAAAAGAAACAATATAAATGTGTACTTGAGAAATTTTAATTACTCAGAGATTTTCATACAATATTTATTCTGCAAGCAAATTCAGGGATTGAATCAATAAAACTGATACAGTATTCTTTGTAGGTATCTGTGTAAACATCAATTTCTGAATCAGTGTTGTAATATTTTCTTAACACATACACAAATCACATCTTGATCTCATTTTTAAAAATGCGAAAAAACTCATTTTGTGCATGTAGCATTCTGAATTCGCACCATCAATGCATGATAGTTCTTGTGTTTTGACATTTATATTGCTATTTACAATTATGAAGAATCATGGTTGAATGATTCTAATAGGTGAGAATGGTACTGAACTTCAGTTAAATGCACATTCCTCCTTAATAAAAATCTACATTTAACATATTTTTCATATAGTTTTTGCTGAATATGCCTCTTCTCATATTTGGTTCATTTTTAACTGTATTGTTTTCTTTTGATTAGTTGTAAAGTTACTTGCATATCTATTATAAAAGTCATCTAACAAATTAAAATAATTCTATTTAACAAGATGTGACTTGGAGTATTTTCTCCAAGTCTGTGGCTGTCTTTTAATCCCCTTATCAGTGTGTATTGCAGCAAAAGTGTGTGTGTGTGTGTTTATACAAATTAGAATATTAAAAAAATAAAATTATTCATCTAGAAAAAAATGCAACAATCGATGTCCGAGTACTCATAGATCGGAAGAGTGTTGTGTA
>7697106b-e15a-4059-b7a2-6c7be5ad6a4c_0 1 1 1 barcode:ACTCTCTGACTTAGGT umi:86120
CCTCTGGATTCACCCTCAGTAGCTATGGCATGCACTGGGTCCGCCAGGCTCCAGGCAAGGGGCTGGAGTGGGTGGCAGTTATAATGATGGAAGTAATAAATACTATGCCAGAAACTGCATCTGGCCACATCACCATCTCCAGAGACAATTCCAAGAACACACTGTGCAAATGAACAGCCTGAGGGGGCACGGCTGTGTATTACTATAAAGACACAGTGAGGGGGAAGTGACTGTAGCCAGACAAACCTCCCCTGCAGGAACGCTGGCAGGGGAGTGAGGGGCGGGGGCACTAGAGGAGCCACTGATCAGAGTCAGCCCTGGAGGCAGGTGCAAGTGAAGCTGTTTCCTGTCAGGATGTATGGGACTTTGTCTTCTTCTAGCAAGTTGCCAGGAACCTCTAAATTTAGAAAACTGTCCTAACAATGTCTTCTCTATGCATATGAGGACCTTTTCTCCCTGGCACAAAATGCGGACTGACACTGACACAGATGAGATCCCTCAACCATGGTCACAAGGATCAGAGTCCTGAGTAACCTCAGGGCTTCCTGGTGAGTCTTCTCCCAATCAGACCCCAGGACAGGGACCTCCGTGAGATTCCCTGACTGGAACAGTCTATGGATCCTGGTCACAGACAGATAGAGAGGCTGAACCAGGGTTGAGCATTGTCATGTCAGAACCTCACAGATTTCACGTCTGATCCTTCTCCTGACACCAAAAATCTGTCAAATCAGTATCAGCACGATTGTGGTTGTTCTTCTGCTCCTAATCCATTTACTTTATTTTTTCATTCCATTTTCTTCCTTTTCCCATTTGTTTTTCCTGCTTTTCCAGAAAGGAAGATGTTTTCCCTGTGGGGGTAGGCAGGGGATGACAATTTGGGAGATGGCTGAAAATACAATATCCTCAGGGCCGGCCATCAGTAAGTGCAGGCTGGAAGTCTCAGAAAGAGCTGAAGCTGCTTAATCACCGTAGAGTACCTTCTCCAGTTCTGCTCTGATGGAATCAGGGCCAAGCAGGTTATCAATGATAATCTACCTAACATAGAGTCAACCAATTACAGGTTCAATAACATCTGTTTAAAAATTCACACCACCACCTGGATTAGTGTTTTGTAAAATAAATAACAGTGTTGTAGGCTGAGGTGAACCAAGAAACCAGCATTTACCATGGAGAAAAACATGACCTTAGATTATAGGTTCTTACAGTGTTAAAGGTGTAAAACTGATTATTAAAAATGAGGCTATTTTTCTTTTTGCTGTTGAGTTGTAAAAACTATTTTTGCATTTGGACATAAAACTTTTTGAGATATATAGCATATTGAAAAACGGTTAATTTACTAAGTCAGAGAGTAGATCAGAAAGCGTCGTGTGG
>5f6c999c-1351-4b64-a6d2-c6265042976c_1 1 1 1 barcode:TCGGAATGCGCTCTGA umi:88503
GTCACCATATCAGTAGACACGTCCAAGAACCAGTTCTCCCTGAAGCTGAGCTCTGGGCCCCACTGCGGACAGCGGCCATCTAGTACTCTGCGACACAGTGAGGGGAGGCAGAGAAGAGGCCATGTGGAGACGTAATGCTCTAGAAGGTGGGTTTACAAGACCTAGCCAACCAACCCCTTCCAGACATTAGATCTCAACATTCAGACTGCAGAATTTTGCAGAAAATCAATATTTGTTCTTTAAGCCACCCACTCCTGTCTCTTCTTATGAAGAATAACAACTAATACCACATAACTATCTTAGCGCTGTCCCCTGGATGCAAGTCAGCCCGCTGGGGCTGGGCACATCTCTCAGATTTCCACATAAAGTAAGCAAATTGTCTGATATAAAAATTTGTCATGTCCCTGTTGGCCAATTTCTGGGCAAGGTCTTTTAAAGAAGCCCTGGGATTTGTGCAGAGTTACCTTTATGATTTATTAGGACATAACTGAAAGAATAAGATGAGTACCAGTTGGATGGAGACTGACCCACTTGACCATCTTCTGCCTGTCTCCTAAGTATGCCACAGAAAACCAGCACCAACATTACTCTATGTCTTCAACTTTCTAAATTTGCACTGATTGGTATTTAAGGCAGGCCCAGCGTTAAATAACTCCTTAGTTTTGCTTCTCTGGGAAAGGTCTTATCTATCCTGGCCTTGGTCTTCAGAGTTTCAGCAATTCTGGGAAGCCAAGGACGCCTCTATCTCCTCCTCCATGCTCTGCAACTCACCTGAGAACAGCTTTTCTCCAATTCAAATGTCTCTGTATTAAAATGCAGGCTCCTGTTTCAGGCTTGGGTGCCTGAGTACACCTACTGGATGCAGCCAGGATTTGGAACTTTCCAGAACACATCACCTGAGAAATGACCAGTCACACTGCTACACTTTCACAATTTCCCCACTTCCCCTCATGGGGAGAAAATTAAAATTGCAGCTTTTTCAAAAAAAAAAAAATTAAGGAAGGCCTACAAGCGCATTCCGAAGATCGGAAGAGCGTGTGTAG

Question 2:
Others libraries aren't finish, This step assembly is very slow and doesn't have multi-threads supported

TRUST4/main.cpp

Line 1361 in 5cdff13

PrintLog( "Processed %d reads (%d are used for assembly).", i + 1, assembledReadCnt ) ;

Now every 100,000 reads over 3 hours.
I have e.g. Found 5,412,778 reads in library A.

dong@gauss01:~/nice_env$ tail ~/nice_env/*_OutputB/B*/bams/no*
[Mon Dec 19 20:23:26 2022] Processed 900000 reads (892585 are used for assembly).
[Mon Dec 19 21:39:39 2022] Processed 1000000 reads (992276 are used for assembly).
[Mon Dec 19 22:58:25 2022] Processed 1100000 reads (1092016 are used for assembly).
[Tue Dec 20 04:26:53 2022] Processed 1200000 reads (1191804 are used for assembly).
[Tue Dec 20 07:44:04 2022] Processed 1300000 reads (1291662 are used for assembly).
[Tue Dec 20 10:34:19 2022] Processed 1400000 reads (1391486 are used for assembly).
[Tue Dec 20 12:39:14 2022] Processed 1500000 reads (1491280 are used for assembly).
[Tue Dec 20 13:42:26 2022] Processed 1600000 reads (1590726 are used for assembly).
[Tue Dec 20 14:35:00 2022] Processed 1700000 reads (1687966 are used for assembly).
[Tue Dec 20 15:07:35 2022] Processed 1800000 reads (1787699 are used for assembly).
dong@gauss01:~/nice_env$ tail ~/nice_env/*_OutputA/A*/bams/no*
[Mon Dec 19 16:06:58 2022] Processed 1000000 reads (986307 are used for assembly).
[Mon Dec 19 16:35:11 2022] Processed 1100000 reads (1084467 are used for assembly).
[Mon Dec 19 16:58:31 2022] Processed 1200000 reads (1182196 are used for assembly).
[Mon Dec 19 17:33:26 2022] Processed 1300000 reads (1279995 are used for assembly).
[Mon Dec 19 17:47:20 2022] Processed 1400000 reads (1379355 are used for assembly).
[Mon Dec 19 18:50:39 2022] Processed 1500000 reads (1466675 are used for assembly).
[Mon Dec 19 23:46:44 2022] Processed 1600000 reads (1563329 are used for assembly).
[Tue Dec 20 06:33:45 2022] Processed 1700000 reads (1661531 are used for assembly).
[Tue Dec 20 11:02:22 2022] Processed 1800000 reads (1759975 are used for assembly).
[Tue Dec 20 14:32:38 2022] Processed 1900000 reads (1858779 are used for assembly).
dong@gauss01:~/nice_env$ tail ~/nice_env/*_OutputD/D*/bams/no*
[Mon Dec 19 15:28:44 2022] Processed 600000 reads (594800 are used for assembly).
[Mon Dec 19 15:49:14 2022] Processed 700000 reads (693959 are used for assembly).
[Mon Dec 19 16:55:15 2022] Processed 800000 reads (792579 are used for assembly).
[Mon Dec 19 17:27:15 2022] Processed 900000 reads (884203 are used for assembly).
[Mon Dec 19 18:56:59 2022] Processed 1000000 reads (982581 are used for assembly).
[Mon Dec 19 20:19:44 2022] Processed 1100000 reads (1081553 are used for assembly).
[Mon Dec 19 21:46:40 2022] Processed 1200000 reads (1180930 are used for assembly).
[Mon Dec 19 23:21:31 2022] Processed 1300000 reads (1280487 are used for assembly).
[Tue Dec 20 08:07:24 2022] Processed 1400000 reads (1380028 are used for assembly).
[Tue Dec 20 13:59:56 2022] Processed 1500000 reads (1479734 are used for assembly).

It's possible speed up ? I've add argument --repseq & -k 53

Question 3:
if impossible to Question 2. how to use *_longRead_toassemble.fq, *_longRead_toassemble_bc.fa & *_longRead_toassemble_umi.fa
to get *_longRead_assembled_reads.fa
I don't understand the nomenclature of output assembled_reads.fa in header.
e.g. current read 0 1 1 1 to next read 1 1 1 1 and umi with the number.

Thank you,
Chuang

[mourisl]

Q1: the candidate the reads are the reads with some amount of k-mer match against the reference database. So they may not really from the VDJ region and may be filtered before assembly.

Q2: I recently added another parameter "--minHitLen", which controls the size of a valid read-contig overlap. Maybe you can try some large values for this, e.g. 1000 .

Q3: the numbers in the header of assembled_reads are the strand and some k-mer count information for the read (min, medium, and max).

[Chuang1118]

Hi,
Thank you for your reply.

Q1:
In your paper figure 1 a1

Candidate reads link to fastq-extractor or bam-extractor in my first impression.
if not, I don't understand why you add argument e.g. TRUST4/mouse/GRCm38_bcrtcr.fa

According to your explanation, this step candidate reads is also in step Contigs(kmer).
In case of short read, I always observe a a reasonable log file.
e.g.

[Wed Dec 21 16:53:33 2022] Read in and count kmers for 5500000 reads.
[Wed Dec 21 16:53:35 2022] Read in and count kmers for 5600000 reads.
[Wed Dec 21 16:53:38 2022] Read in and count kmers for 5700000 reads.
[Wed Dec 21 16:55:42 2022] Found 5769294 reads.
[Wed Dec 21 16:55:55 2022] Finish sorting the reads.
[Wed Dec 21 16:57:38 2022] Finish rough annotations.
[Wed Dec 21 16:57:40 2022] Processed 100000 reads (100000 are used for assembly).
[Wed Dec 21 16:57:41 2022] Processed 200000 reads (200000 are used for assembly).
[Wed Dec 21 16:57:41 2022] Processed 300000 reads (300000 are used for assembly).
[Wed Dec 21 16:57:41 2022] Processed 400000 reads (400000 are used for assembly).
....
[Wed Dec 21 16:58:46 2022] Processed 5500000 reads (5181489 are used for assembly).
[Wed Dec 21 16:58:53 2022] Processed 5600000 reads (5275710 are used for assembly).
[Wed Dec 21 16:59:00 2022] Processed 5700000 reads (5363664 are used for assembly).
[Wed Dec 21 16:59:03 2022] Assembled 5381065 reads.
[Wed Dec 21 16:59:03 2022] Try to rescue 6458 reads for assembly.

However, I've lost a half reads in case long read. It's abnormal for me?
It's there have the limitation on maximum characters ?

Q2:
parameter "--minHitLen" can be speed up contig step.
But I don't know why called minHitLen, why not maxHitLen?
Does less length read need to contig building not?

Q3:
I want to just avoid write own script to adapter format *_assembled_reads.fa

Just one suggestion. like I've a lot of project dependent TRUST4.
One paper in preparing, others in processing.
Could you please check TRUST4 work correctly then push in master.

I've pull a master verison 1.0.7 6 July 2022.
All the v_identity and j_identity have same values in barcode_airr.tsv.

v_identity	 j_identity
95.5	95.5
97.92	97.92
97.28	97.28
100	100
94.52	94.52
97.62	97.62
NA	NA
98.99	98.99
100	100
96.26	96.26
97.59	97.59
NA	NA
98.61	98.61
98.98	98.98
99.66	99.66
100	100
98.63	98.63

master verison 1.0.8 12 Dec 2022.
with error

has not barcode fieldsystem

master verison 1.0.8 21 Dec 2022 works fine. mostly v_identity and j_identity have the distinct values

Thank you
Chuang

[mourisl]

Q1: the candidate reads are extracted by the fastq-extractor or bam-extractor. The process is to check the alignment coordinate or map the reads against the reference genes, e.g. GRCm38_bcrtcr.fa. The mapping is done by chaining concordant k-mer matches. The k-mer size is quite small for this stage, and the long reads have a higher chance to get a valid chain by chance. As a result, many candidate reads in long read might be really from VDJ region.

There is no limit on the maximum characters in the recent updates (in master branch but not in the release)

Q2: --minHitLen is the minimal overlap size to add the read to the contig, or the minimum hit length to consider a read can be aligned to a contig.

Q3: The VJ gene similarity was indeed an issue in 1.0.7 and was fixed later. I think the current master branch works well, and I will draft a new release (v1.0.9) later.

[Chuang1118]

There is no limit on the maximum characters in the recent updates (in master branch but not in the release)
There is an interesting point.
I work with TRUST4 master branch version. I think my run is in the case of the limitation on maximum characters.
How about the mechanism if my run pass the maximum characters ?
TRUST4 crash ? or write a message in log file ? or like below?
This is other sample with parameter --minHitLen 1000

[Wed Dec 21 11:22:48 2022] Read in and count kmers for 3400000 reads.
[Wed Dec 21 11:23:16 2022] Read in and count kmers for 3500000 reads.
[Wed Dec 21 12:10:21 2022] Found 3588880 reads.
[Wed Dec 21 12:11:31 2022] Finish sorting the reads.
[Wed Dec 21 12:19:25 2022] Finish rough annotations.
[Wed Dec 21 12:19:59 2022] Processed 100000 reads (99387 are used for assembly).
...
[Thu Dec 22 00:45:27 2022] Processed 3100000 reads (3061306 are used for assembly).
[Thu Dec 22 01:03:00 2022] Processed 3200000 reads (3159777 are used for assembly).
[Thu Dec 22 01:10:18 2022] Assembled 3213987 reads.

Giving kmers for 3500000 reads, As expected processed 3500000 reads( * are used for assembly).
However, last line write processed 3200000 reads (3159777 are used for assembly) in log file.

I have checked these 2 cases by file size, they are not same.
I just want to make sure that I correctly use TRUST4.

Could you tell me the limitation on the maximum characters in which script and which line please ?

Chuang

[mourisl]

Just want to confirm, does the maximum character limitation mean the length of a read? I'm also not sure I understand your question.

[Chuang1118]

I mean the length of a read in the buffer c++.
I don't think that I have passed the maximum limitation.
Instead of I have the problem understand log file.
e.g.

[Mon Dec 19 12:58:55 2022] Read in and count kmers for 2900000 reads.
[Mon Dec 19 12:59:14 2022] Read in and count kmers for 3000000 reads.
[Mon Dec 19 12:59:30 2022] Read in and count kmers for 3100000 reads.
[Mon Dec 19 13:37:09 2022] Found 3137630 reads.
[Mon Dec 19 13:38:18 2022] Finish sorting the reads.
[Mon Dec 19 13:45:34 2022] Finish rough annotations.
[Mon Dec 19 13:48:57 2022] Processed 100000 reads (99656 are used for assembly).
[Mon Dec 19 13:57:14 2022] Processed 200000 reads (199143 are used for assembly).
[Mon Dec 19 14:15:53 2022] Processed 300000 reads (298632 are used for assembly).
[Mon Dec 19 14:49:18 2022] Processed 400000 reads (396027 are used for assembly).
[Mon Dec 19 15:40:29 2022] Processed 500000 reads (495769 are used for assembly).
[Mon Dec 19 16:45:48 2022] Processed 600000 reads (595611 are used for assembly).
[Mon Dec 19 18:05:36 2022] Processed 700000 reads (695479 are used for assembly).
[Mon Dec 19 19:11:43 2022] Processed 800000 reads (795326 are used for assembly).
[Mon Dec 19 19:52:53 2022] Processed 900000 reads (895024 are used for assembly).
[Mon Dec 19 20:12:37 2022] Processed 1000000 reads (993744 are used for assembly).
[Mon Dec 19 20:40:54 2022] Processed 1100000 reads (1093530 are used for assembly).
[Mon Dec 19 21:32:32 2022] Processed 1200000 reads (1192878 are used for assembly).
[Mon Dec 19 22:27:30 2022] Processed 1300000 reads (1292568 are used for assembly).
[Mon Dec 19 22:51:58 2022] Processed 1400000 reads (1392090 are used for assembly).
[Mon Dec 19 22:55:49 2022] Assembled 1434448 reads.

Thank you again for all your answers before.

[Chuang1118]

One of my run has segmentation fault is similar to #29, #116 or #22

[Fri Dec 23 05:31:07 2022] SYSTEM CALL: /opt2/TRUST4/annotator -f /home/cdong/TRUST4-master131222/humanIGH_IMGT+C.fa -a trust4/B_102_0023_longRead_final.out -t 24 -o trust4/B_102_0023_longRead --barcode --UMI -r trust4/B_102_0023_longRead_assembled_reads.fa --airrAlignment > trust4/B_102_0023_longRead_annot.fa
[Fri Dec 23 05:31:07 2022] Start to annotate assemblies.
Segmentation fault
system /opt2/TRUST4/annotator -f /home/cdong/TRUST4-master131222/humanIGH_IMGT+C.fa -a trust4/B_102_0023_longRead_final.out -t 24 -o trust4/B_102_0023_longRead --barcode --UMI -r trust4/B_102_0023_longRead_assembled_reads.fa --airrAlignment > trust4/B_102_0023_longRead_annot.fa failed: 35584 at /opt2/TRUST4/run-trust4 line 51.

1,4G -rw-r--r-- 1 cdong cdong 1,4G 23 déc.  06:06 B_102_0023_longRead_airr_align.tsv
3,4G -rw-r--r-- 1 cdong cdong 3,4G 23 déc.  06:06 B_102_0023_longRead_annot.fa
2,9G -rw-r--r-- 1 cdong cdong 2,9G 23 déc.  05:17 B_102_0023_longRead_assembled_reads.fa
 24G -rw-r--r-- 1 cdong cdong  24G 23 déc.  05:30 B_102_0023_longRead_final.out
 24G -rw-r--r-- 1 cdong cdong  24G 23 déc.  05:17 B_102_0023_longRead_raw.out
253M -rw-r--r-- 1 cdong cdong 253M 22 déc.  18:56 B_102_0023_longRead_toassemble_bc.fa
 11G -rw-r--r-- 1 cdong cdong  11G 22 déc.  18:56 B_102_0023_longRead_toassemble.fq
235M -rw-r--r-- 1 cdong cdong 235M 22 déc.  18:56 B_102_0023_longRead_toassemble_umi.fa

I found out one strange read in last line of B_102_0023_longRead_annot.fa

This read has a super length and has cdr3 abnormal in header. Without carriage return in last line.

The context of this read in B_102_0023_longRead_assembled_reads.fa [before/After 10 lines at this read]

Cmd

run-trust4 -f ~/TRUST4-master131222/hg38_chr14_bcr.fa \
-t 24 \
--ref ~/TRUST4-master131222/humanIGH_IMGT+C.fa -b B_102_0023_chr14.tagged.bam \
-k 53 \
--minHitLen 1000 \
--barcode CB \
--UMI UB \
--barcodeWhitelist barcodes_B.txt \
--repseq \
-o B_102_0023_longRead --od trust4

Thank you very much for your help.
Happy holidays,
Chuang

[mourisl]

This could be a bug in TRUST4. Could you please share the sequence with me? I'll take a look.

[Chuang1118]

Hello,

I try get 5%, 10% and 20% from dataset, TRUST4 completed perfect.
Without the full dataset, I can't reproduce the issue.
I will send you a private message with the original dataset.
I found your mail by https://mourisl.github.io

Chuang