我所理解的cellranger軟件理想原始輸入數(shù)據(jù)就是SRA格式适瓦，然后利用sra-tools分為read谱仪、barcode+UMI、index三個fastq.gz文件嗦随。最后直接利用cellranger即可敬尺。但總會發(fā)現(xiàn)文獻提供的數(shù)據(jù)格式并非如此，就要花費一些心思了署恍。

• 文獻：https://www.nature.com/articles/s43018-020-00139-8
• 單細胞測序數(shù)據(jù)：https://trace.ncbi.nlm.nih.gov/Traces/sra/?study=SRP261877
  
  scRNA-seq data
  
  https://trace.ncbi.nlm.nih.gov/Traces/sra/?study=SRP261877
  
  scRNA-seq data

如圖盯质，基本了解到做了三組實驗袭蝗，每組兩個重復。但關鍵作者提供的是bam文件格式朵逝，就要想辦法轉(zhuǎn)換為cellranger所需要的文件格式乡范。

1配名、下載數(shù)據(jù)

• 還是使用ebi數(shù)據(jù)庫下載晋辆，詳見實戰(zhàn)1的介紹。
• https://www.ebi.ac.uk/ena/browser/view/PRJNA632854

cat > bam.txt
fasp.sra.ebi.ac.uk:/vol1/run/SRR117/SRR11798249/ICBtreated_Brca2null_rep1.bam
fasp.sra.ebi.ac.uk:/vol1/run/SRR117/SRR11798250/ICBtreated_Brca1null_rep2.bam
fasp.sra.ebi.ac.uk:/vol1/run/SRR117/SRR11798251/ICBtreated_Brca1null_rep1.bam
fasp.sra.ebi.ac.uk:/vol1/run/SRR117/SRR11798257/ICBtreated_Parental_rep2.bam
fasp.sra.ebi.ac.uk:/vol1/run/SRR117/SRR11798258/ICBtreated_Parental_rep1.bam
fasp.sra.ebi.ac.uk:/vol1/run/SRR117/SRR11798259/ICBtreated_Brca2null_rep2.bam

conda activate download
cat fq.txt |while read id
do ascp -QT -l 300m -P33001  \
-i ~/miniconda3/envs/download/etc/asperaweb_id_dsa.openssh \
era-fasp@$id  .
done
conda deactivate 

rawdata

conda環(huán)境配置可見實戰(zhàn)1芋膘，就是使用下ascp軟件

2、bam轉(zhuǎn)fastq

2.1 special bam of cellranger

• 這里的bam序列比對結果文件應該是作者使用cellranger后的產(chǎn)生的比對結果臂拓。

  
  
  cellranger
• 經(jīng)過一番查詢后习寸，知道了cellranger產(chǎn)生的bam文件里是帶有barcode與UMI的，儲存在tag標簽里霞溪。
  https://bioinformatics.stackexchange.com/questions/7096/bam-to-gene-expression-matrix-umi-counts-per-gene-per-cell-10x
  https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/output/bam

samtools view ICBtreated_Parental_rep1.bam | less -SN
samtools view ICBtreated_Parental_rep1.bam | head -3 | tr "\t" "\n" | cat -n

bam tag

相關標簽具體解釋見結尾鸯匹，這里知道 CB或者CR代表barcode、UB或者UR代表UMI即可幼东。

2.2 cellranger bamtofastq

• 知道上述知識點后科雳，我一開始想手動提取下bam文件里的barcode與UMI序列，也查了很多l(xiāng)inux字符處理方法简逮。但這樣之后還要自己組裝fastq尿赚，并且使header一致”辏總之很費事冰寻。
• 后來知道 cellranger也有bamtofastq功能，就試試看和其它軟件的轉(zhuǎn)換有什么不同轻腺。
• 如下圖是一個bam文件的處理結果划乖。其很強大的功能：自動根據(jù)bam文件，生成配套的三種文件误算。而且還修改為規(guī)范命名！I1,R1,R2的含義見實戰(zhàn)1
  
  ICBtreated_Brca1null_rep1.bam

three type fastq

I1代表的index序列筒占，一般是用于區(qū)分混合樣品的蜘犁，二代測序通配这橙。這里為空导披，表示bam文件里沒有。而且的確也不需要鹰晨，就提供一個空序列文件即可止毕。

• 批量處理

cat > bamtofastq.sh
bin=/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger
cat name.list |while read id
do
$bin bamtofastq $id ./fastq/${id}
done

bash bamtofastq.sh

3、cellranger count

• 經(jīng)過上一步的處理忍疾，接下來就比較簡單了谨朝。

find /home/shensuo/test/fastq/ | grep bam/out | grep -v XX/bam > fq.txt

bin=/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger
db=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/
target=bamtofastq

cat fq.txt |while read id
do
echo $bin count --id=${id:0-39:14} \
--localcores=4 \
--transcriptome=$db \
--fastqs=$id \
--sample=$target \
--expect-cells=3000 \
--nosecondary  
done

cat > cellranger.sh
/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger count --id=ICBtreated_Brca1null_rep1 --localcores=10 --transcriptome=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/ --fastqs=/home/shensuo/test/fastq/ICBtreated_Brca1null_rep1.bam/output_0_1_HG3HHDRXX --sample=bamtofastq --expect-cells=13009 --nosecondary
/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger count --id=ICBtreated_Brca1null_rep2 --localcores=10 --transcriptome=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/ --fastqs=/home/shensuo/test/fastq/ICBtreated_Brca1null_rep2.bam/output_0_1_HG3HHDRXX --sample=bamtofastq --expect-cells=21789 --nosecondary
/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger count --id=ICBtreated_Brca2null_rep1 --localcores=10 --transcriptome=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/ --fastqs=/home/shensuo/test/fastq/ICBtreated_Brca2null_rep1.bam/output_0_1_HG3HHDRXX --sample=bamtofastq --expect-cells=18684 --nosecondary
/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger count --id=ICBtreated_Brca2null_rep2 --localcores=10 --transcriptome=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/ --fastqs=/home/shensuo/test/fastq/ICBtreated_Brca2null_rep2.bam/output_0_1_HG3HHDRXX --sample=bamtofastq --expect-cells=16731 --nosecondary
/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger count --id=ICBtreated_Parental_rep1 --localcores=10 --transcriptome=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/ --fastqs=/home/shensuo/test/fastq/ICBtreated_Parental_rep1.bam/output_0_1_HG3HHDRXX --sample=bamtofastq --expect-cells=13292 --nosecondary
/home/shensuo/biosoft/cellranger/cellranger-4.0.0/bin/cellranger count --id=ICBtreated_Parental_rep2 --localcores=10 --transcriptome=/home/shensuo/biosoft/cellranger/test/refdata-gex-mm10-2020-A/ --fastqs=/home/shensuo/test/fastq/ICBtreated_Parental_rep2.bam/output_0_1_HG3HHDRXX --sample=bamtofastq --expect-cells=20718 --nosecondary

nohup  bash /home/shensuo/test/cr_out/fq.txt &
#需要使用全路徑

耐心等待結果即可字币。估計一夜肯定是需要的。
此外其中有幾個注意點士复，具體如下

• 一開始想只用while語句運行共苛，發(fā)現(xiàn)還是不能盡善盡美。就echo下澄峰，再根據(jù)實際情況修改辟犀，保存為cellranger.sh
• 關于--localcores=設置绸硕，可根據(jù)實際情況玻佩。我是設置10席楚，之后也單獨嘗試了32，24烦秩，20。結果還是24還比較適合目前得到服務器環(huán)境的最大承受兜蠕，如果想盡快跑完的話抛寝。
• 關于--expect-cells設置，主要參考原文獻的測序結果細胞數(shù)盗舰。
  
  paper result

4岭皂、導出最終結果

• 如前所述，cellranger count結果很多书劝，主要是需要其中如下圖的三個文件（每個樣品）

  
  
  three type result

cat name.list | while read id
do 
mkdir ./out/$id
cp $id/outs/filtered_feature_bc_matrix/* ./out/$id
done

接下來就可以购对，將數(shù)據(jù)導入到R中陶因，用seurat等包進行下游分析了。

附：cellranger indexed bam

https://support.10xgenomics.com/single-cell-gene-expression/software/pipelines/latest/output/bam

image.png

• CB、CR切诀、CY表示barcode搔弄，一般是16個堿基顾犹；
• UB褒墨、UR、UY表示UMI浑玛，一般是10個堿基噩咪。
• R一般代表原始測序數(shù)據(jù)，Y代表質(zhì)量分數(shù)涨享，而B代表校正后的R仆百，可能對應堿基質(zhì)量分數(shù)太低等因素。一般來說R與B都是相同的吁讨。