aggr代码 cellranger_cellranger使用的初步探索（1）

本篇笔记是按照单细胞天地公众号教程里的代码练习的，因为之前了解到cell ranger运行需要大量的内存，无奈自己的电脑配置不行，现在有了服务器，那就要把这一节补起来~

参考文章：

1.单细胞实战(二) cell ranger使用前注意事项

2.单细胞实战(三) Cell Ranger使用初探

3.单细胞实战(四) Cell Ranger流程概览

练习数据：GSE117988。首先在GEO网站上看一下作者对数据的处理方法：

对于这几个原始数据来说，其中包括：Read1-26个碱基(UMI)，read2-8个碱基(Index),Read2-98个碱基(RNA read)。对原始的BCL文件使用cell ranger进行mkfastq，然后用cell ranger的count进行分析，使用STAR进行比对(基因组GRCh38和HM011556.1)。

(一)下载sra数据

#!/bin/bash

#module load sratoolkit/2.9.1(这是我调用服务器里的软件的命令，如果是自己的电脑，请忽略这一行)

cat SRR_Acc_List.txt | while read i

do

prefetch $i -O `pwd` && echo "**${i}.sra done**"

done

(二)提取fastq文件

#!/bin/bash

# module load sratoolkit/2.9.1

for i in SRR77229*

do

fastq-dump --gzip --split-files ./$i

done

这一步每一个SRA文件都会生成3个fastq文件。分别标注_1, _2, _3。现在来看看每一个fastq文件都是什么样子的：

$ zless -SN SRR7722942_1.fastq.gz | head #这是每一个_1的read，这里显示每一条read都是8个碱基长度，那么就是上面我们提到的Index。

@SRR7722942.1 SN367:911:HKMNCBCXY:2:1101:1192:1900 length=8

CGCTATGT

+SRR7722942.1 SN367:911:HKMNCBCXY:2:1101:1192:1900 length=8

GGGGGIII

@SRR7722942.2 SN367:911:HKMNCBCXY:2:1101:1112:1988 length=8

CGCTATGT

+SRR7722942.2 SN367:911:HKMNCBCXY:2:1101:1112:1988 length=8

GGGGGIII

@SRR7722942.3 SN367:911:HKMNCBCXY:2:1101:1404:1952 length=8

CGCTATGT

$ zless -SN SRR7722942_2.fastq.gz | head #标注_2的fastq文件长度是26个碱基，也就是上面说的细胞barcode(UMI)

@SRR7722942.1 SN367:911:HKMNCBCXY:2:1101:1192:1900 length=26

CGTTGGGGTCTGCGGTAAATAGGCCA

+SRR7722942.1 SN367:911:HKMNCBCXY:2:1101:1192:1900 length=26

GAGGGIGGIIGGIIIIIGIGGGGGGI

@SRR7722942.2 SN367:911:HKMNCBCXY:2:1101:1112:1988 length=26

CTTAACTTCTCGATGAAGGGGTCTCG

+SRR7722942.2 SN367:911:HKMNCBCXY:2:1101:1112:1988 length=26

GGGGGIIIIIIIIIIGIIIIIIGIII

@SRR7722942.3 SN367:911:HKMNCBCXY:2:1101:1404:1952 length=26

AAGACCTTCGCCCTTATACCGGTCCC

$ zless -SN SRR7722942_3.fastq.gz | head #标注_3的fastq文件每一个Read是98个碱基，也就是RNA的read

@SRR7722942.1 SN367:911:HKMNCBCXY:2:1101:1192:1900 length=98

NNNNNGTGGTATCAACGCAGAGTACATGGGGCNCTTACCGCCATCTTGGCTCCTGTGGNTGNCTGCTGGGAACGGGACTTCTAAAAGNNNNTATGTCT

+SRR7722942.1 SN367:911:HKMNCBCXY:2:1101:1192:1900 length=98

#####<.>

@SRR7722942.2 SN367:911:HKMNCBCXY:2:1101:1112:1988 length=98

CTCAGCCTTAGCCCTGTGCCCCCTGAAACAGCTGCCACCATCACTCGCAAGAGAATCCCCTCCATCTTTGGGAGGGGTTGATGCCAGACATCACCAGG

+SRR7722942.2 SN367:911:HKMNCBCXY:2:1101:1112:1988 length=98

@SRR7722942.3 SN367:911:HKMNCBCXY:2:1101:1404:1952 length=98

GGAGGAGATGGGCAGCATTGTTAAGAGATTGGTACCACTGAGCAAATATGTTGAGAATGATGATGGCAAGGTTTCTCCCTGTTAGAGAAGGTATTTGT

(三)修改fastq文件名

根据10x官网上的说明，在后续处理数据之前，最好把fastq文件的名字改一下：here

改成如下格式：

这里就可以对应上面文章里描述的，index文件改成I1，Read1(UMI)改成R1，Read2(RNA read)改成R2。所以这里fastq_1文件应该改成I1, fastq_2文件改成R1，fastq_3文件改成R2：

$ cat SRR_Acc_List.txt

SRR7722937

SRR7722941

SRR7722938

SRR7722939

SRR7722940

SRR7722942

$ cat SRR_Acc_List.txt | while read i ;do (mv ${i}_1*.gz ${i}_S1_L001_I1_001.fastq.gz;mv ${i}_2*.gz ${i}_S1_L001_R1_001.fastq.gz;mv ${i}_3*.gz ${i}_S1_L001_R2_001.fastq.gz);done

改后文件名就都变成了下面这样：

(四)fas