Missing CB from StarSolo ATACalignment

[akhst7]

HI,

I've been playing with 10x's PBMC from a Healthy Donor - Granulocytes Removed Through Cell Sorting (3k) and I used a following Solo script to generate alignments;

#!/opt/homebrew/bin/zsh

index=/Volumes/MySlateDrive/star_genome_index/index_human_Gencode_GRCh38_p13
whitelist=/Volumes/Bioinformatics/star_genome_index/3M-february-2018.txt
r1=/Volumes/MySlateDrive/pbmc_granulocyte_sorted_3k/gex/pbmc_granulocyte_sorted_3k_S1_L003_R1_001.fastq.gz
r2=/Volumes/MySlateDrive/pbmc_granulocyte_sorted_3k/gex/pbmc_granulocyte_sorted_3k_S1_L003_R2_001.fastq.gz
r3=/Volumes/MySlateDrive/pbmc_granulocyte_sorted_3k/gex/pbmc_granulocyte_sorted_3k_S1_L004_R1_001.fastq.gz
r4=/Volumes/MySlateDrive/pbmc_granulocyte_sorted_3k/gex/pbmc_granulocyte_sorted_3k_S1_L004_R2_001.fastq.gz

STAR --genomeDir $index \
--runThreadN 50 \
--readFilesCommand 'gzip -c -d' \
--readFilesIn $r2,$r4 $r1,$r3 \
--soloCBwhitelist $whitelist \
--soloType CB_UMI_Simple \
--soloCBstart 1 \
--soloCBlen 16 \
--soloUMIstart 17 \
--soloUMIlen 12 \
--soloStrand Forward \
--twopassMode Basic \
--clipAdapterType CellRanger4 \
--soloMultiMappers EM \
--soloFeatures Gene GeneFull \
--outFilterScoreMin 30 \
--soloCBmatchWLtype 1MM_multi_Nbase_pseudocounts \
--soloUMIfiltering MultiGeneUMI_CR \
--soloUMIdedup 1MM_CR \
--soloCellFilter  EmptyDrops_CR \
--outReadsUnmapped Fastx \
--outFileNamePrefix /Volumes/MySlateDrive/pbmc3Ksorted.Solo \
--outSAMtype None

A following is a final output;

Started job on |       Apr 07 17:28:08
                             Started mapping on |       Apr 07 17:28:46
                                    Finished on |       Apr 07 17:32:55
       Mapping speed, Million of reads per hour |       1116.37

                          Number of input reads |       77215568
                      Average input read length |       99
                                    UNIQUE READS:
                   Uniquely mapped reads number |       66533327
                        Uniquely mapped reads % |       86.17%
                          Average mapped length |       98.05
                       Number of splices: Total |       28112
            Number of splices: Annotated (sjdb) |       7168
                       Number of splices: GT/AG |       11499
                       Number of splices: GC/AG |       3037
                       Number of splices: AT/AC |       51
               Number of splices: Non-canonical |       13525
                      Mismatch rate per base, % |       0.22%
                         Deletion rate per base |       0.01%
                        Deletion average length |       1.54
                        Insertion rate per base |       0.00%
                       Insertion average length |       1.35
                             MULTI-MAPPING READS:
        Number of reads mapped to multiple loci |       3146380
             % of reads mapped to multiple loci |       4.07%
        Number of reads mapped to too many loci |       290961
             % of reads mapped to too many loci |       0.38%
                                  UNMAPPED READS:
  Number of reads unmapped: too many mismatches |       0
       % of reads unmapped: too many mismatches |       0.00%
            Number of reads unmapped: too short |       6951816
                 % of reads unmapped: too short |       9.00%
                Number of reads unmapped: other |       293084
                     % of reads unmapped: other |       0.38%
                                  CHIMERIC READS:
                       Number of chimeric reads |       0
                            % of chimeric reads |       0.00%

Also, flagstats;

samtools flagstats -@ 40 sorted.atac_SoloAligned.bam
158396246 + 0 in total (QC-passed reads + QC-failed reads)
139359414 + 0 primary
19036832 + 0 secondary
0 + 0 supplementary
0 + 0 duplicates
0 + 0 primary duplicates
158396246 + 0 mapped (100.00% : N/A)
139359414 + 0 primary mapped (100.00% : N/A)
139359414 + 0 paired in sequencing
69679707 + 0 read1
69679707 + 0 read2
139359414 + 0 properly paired (100.00% : N/A)
139359414 + 0 with itself and mate mapped
0 + 0 singletons (0.00% : N/A)
0 + 0 with mate mapped to a different chr
0 + 0 with mate mapped to a different chr (mapQ>=5)

When I looked at CR and CB, there are many CB missing as follows;

samtools view atac_SoloAligned.out.bam |cut -f14,15 |head -20
CR:Z:TTACCGGGTACTAGGA	CB:Z:-
CR:Z:TTACCGGGTACTAGGA	CB:Z:-
CR:Z:ATGTTGTTCCTAATCG	CB:Z:-
CR:Z:ATGTTGTTCCTAATCG	CB:Z:-
CR:Z:GCTTATGTCGCTCAGT	CB:Z:-
CR:Z:GCTTATGTCGCTCAGT	CB:Z:-
CR:Z:TGCGGTTTCATCAGCC	CB:Z:-
CR:Z:TGCGGTTTCATCAGCC	CB:Z:-
CR:Z:TAGCATGCAATGAGCC	CB:Z:-
CR:Z:TAGCATGCAATGAGCC	CB:Z:-
CR:Z:AATGGAACAACAAAGA	CB:Z:AATGGAACAACAAACA
CR:Z:AATGGAACAACAAAGA	CB:Z:AATGGAACAACAAACA
CR:Z:CCTGAATGTAGCAGGT	CB:Z:-
CR:Z:CCTGAATGTAGCAGGT	CB:Z:-
CR:Z:ATTGCGTGTACTTACA	CB:Z:-
CR:Z:ATTGCGTGTACTTACA	CB:Z:-
CR:Z:TGGCCAGCACAACATG	CB:Z:-
CR:Z:TGGCCAGCACAACATG	CB:Z:-
CR:Z:TAATCCGTCTTTAAGC	CB:Z:-
CR:Z:TAATCCGTCTTTAAGC	CB:Z:-

In fact, there are 3102544 reads with the barcodes, 0.002% of 1583962462 total reads. I appreciate any pointers.

[alexdobin]

Hi @akhst7

the CB:Z:- means that CB sequence was not matched to the passlist. Maybe this sample needs a different passlist?

[akhst7]

@alexdobin, I needed 737K-arc-v1.txt.gz

Thanks.