上游分析

1. 安装linux

这就不多说了，自己搞一个虚拟机，我用的是Centos7。

ps：如果使用的是学校集群的话，注意在修改密码中改一下自己的密码，开启后账号为：root，密码自定义（注意是暗文，你敲进去是不会显示的）结束了enter即可

2. 预先安装

首先要安装anaconda，为了不污染环境

1
2
3
4
5
6
7
8
-- 安装linux安装包（如果报错自己去anaconda网站找地址）
wget https://repo.anaconda.com/archive/Anaconda3-2023.07-Linux-x86_64.sh

-- 解压anaconda，后面的就是刚刚安装好的名字
bash Anaconda3-2023.07-Linux-x86_64.sh

-- 更新环境变量
source ~/.bashrc

至此，anaconda安装完成。我们再创建一个环境，专门用于生信分析

tips：anaconda教程请查看TensorFlow1教程，需要包含的步骤为：

  1. create一个虚拟环境
  2. 设置镜像中除了最后一个，其他都跑
  3. 下载fastqc

3. 安装分析软件

这些软件因人而异，不一定需要全下载

1
2
3
4
5
6
7
8
9
10
11
12
-- 最主要的分析fast格式的软件
conda install fastqc

-- 有时候fastq不止一个文件，联合分析
conda install multigc

-- 高通量测序数据中去除接头序列和其他不需要的序列
conda install -c bioconda cutadapt

-- 用于将sra转化为fastq文件的脚本
conda install sra-tools

cellranger单独拎出来讲，这个软件应该使用到的频率比较大

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
-- 下载安装包（上官网找，复制粘贴即可，如果报错和anaconda一样）
wget -O cellranger-8.0.1.tar.gz "https://cf.10xgenomics.com/releases/cell-exp/cellranger-8.0.1.tar.gz?Expires=1724187641&Key-Pair-Id=APKAI7S6A5RYOXBWRPDA&Signature=Nym8cB6esJ3Zk7nueAU7BG3hhF2IZBp9FpD4OE5gW0a7C-m4ob7lChp0W-j7ydgnEBYafZ~igPGR~DzUq4CxsXS4XwkENuKhwn7Xr7RbdO~wGGh03fWzYyYt~Y7FK~V~73DzJEplvcls0p~KdbcQYvb7NflwtO9YMY~FnO4fB2VmFf5QBMdpXbPubMG~jNWE58ki7zr6ilsMGAfEnI6Po4cpZKKe1VCN7zJeipUKqS9qj~jGRpIGjHV9abluAPeodE-zCk0F-bsMpXSx6x~avzQfrN6ViJoRNEBemaB7anzMOTJ7L3XEP~QprOaJKKobFWZBGz4MBXokQBxByAZObQ__"

-- 解压，xyz为版本号
tar -xzvf cellranger-x.y.z.tar.gz

-- 配置环境变量
source ~/.bashrc

-- 查看是否成功
cellranger --version

-------如果失败，进入到vim中写，vim对于初学者来说比较复杂，建议认真读教程，不要多做也不要少做-------
-- 1. 打开文件
vim ~/.bashrc
-- 2. 进入编辑模式 ，按键盘i键
-- 3. 添加cellranger的路径，将光标移到最后一行，添加内容（-8.0.1是cellranger文件夹的名字，看自己下载的名字是什么）
export PATH=$PATH:/root/cellranger-8.0.1/bin
-- 4. 先按Esc退出编辑模式，然后输入:wq，最后按enter
-- 5. 重加载
source ~/.bashrc
-- 6.查看是否成功
cellranger --version

4. 测序数据质量评估

ok，开始正式的质量评估

4.1 cutadapt

去除接头序列：

在测序过程中，接头序列会附加到 DNA 或 RNA 片段的两端。cutadapt 可以从测序读段（reads）中识别并去除这些接头序列。

去除低质量碱基：

除了接头序列，还可以去除位于读段两端的低质量碱基，保证下游分析中使用的是高质量的数据。

处理指定长度的序列：

可以设置参数来过滤特定长度的序列，比如去除过短或过长的序列，这有助于控制下游数据分析中的数据质量。

双端测序数据处理：

cutadapt 支持双端测序数据，可以同时处理两端的序列，并保持双端读段的配对关系。

ps：在做cutadapt之前，首先是需要通过fastqc的质控报告观测是否需要去除的，具体在Adapter Content这一栏中，如果是打钩就不需要，打叉的话分情况讨论：

  1. 数据数据集上down下来的数据一般会显示测序的平台是啥，这时候使用平台默认的就ok了
  2. 自己的数据也知道用的测序平台是啥，同上

如果用的是illumina的话，read1为AGATCGGAAGAGC，read2为AATGATACGGCGACC。如果用的是其他测序平台的话，建议翻看手册查询（一般都会给出的）

1
2
3
4
5
6
7
8
9
10
-- 检查接头序列
cutadapt --detect-adapters -o /path/to/trimmed_output.fq /path/to/input.fq

-- 本数据为read2（illumina），read1同理，反正都要处理
cutadapt -a AGATCGGAAGAGC -q 20 --minimum-length 50 \
         -j 4 \
         -o /root/rowdata/Rhesus-Liver-1_L1_2_trimmed.fq.gz \
         /root/rowdata/Rhesus-Liver-1_L1_2.fq.gz


在预处理完数据之后，就可以通过fastqc可视化数据质量了，可以百度网页结果

1
fastqc -t 4 -o ~/fastqc_output /root/rowdata/Rhesus-Liver-1_L1_2.fq.gz

4.2 cellranger

  1. 创建自定义参考基因组
1
2
3
4
5
6
7
8
9
10
11
12
13
-- 这一步是因为我的gtf中，Curated Genomic包含空格，因此需要替换为Curated_Genomic
sed 's/Curated Genomic/Curated_Genomic/g' GCF_003339765.1_Mmul_10_genomic.gtf > fixed_GCF_003339765.1_Mmul_10_genomic.gtf

-- 设置核心数，按照自己的配置来定，一般max-1就ok了
export CELLRANGER_NUM_THREADS=15

-- 如果没有处理过gtf的话，原名就ok，如果sed了，那就改成自己处理后的gtf文件名，output-dir为文件输出的位置，可自定义
cellranger mkref --genome=mmul_ref \
--fasta=/root/compare/GCF_003339765.1_Mmul_10_genomic.fna \
--genes=/root/compare/fixed_GCF_003339765.1_Mmul_10_genomic.gtf \
--output-dir=/root/compare/mmul_ref


  1. 使用 cellranger count 进行比对和分析
1
2
3
4
5
6
7
8
9
10
cellranger count \
  --id=liver_analysis \  # 输出文件夹的名称
  --transcriptome=/root/compare/mmul_ref \  # 参考基因组路径
  --fastqs=/root/liver_PE \  # 存放FASTQ文件的目录
  --sample=sample_name \  # 样本名称，一般为liver_PE-1或liver_PE-2这样的
  --expect-cells=1000 \  # 预期细胞数，可根据实验设定调整（也可不用，默认1k）
  --localcores=16 \  # 使用的核心数
  --localmem=180  # 使用的内存（GB）
  --create-bam=true  # 通常默认，但是当发生该参数缺少的报错时可以显式添加

需要注意的是，如果我有liver_PE-1和liver_PE-2，那需要分别跑两次，然后对其进行合并

  1. 使用cellranger arr对结果进行合并

首先船舰csv文件，列出每个样本的molecule_info.h5 文件路径，内容如下：

1
2
3
library_id,molecule_h5
liver_PE-1,/root/liver_analysis1/outs/molecule_info.h5
liver_PE-2,/root/liver_analysis2/outs/molecule_info.h5

随后对其进行合并

1
2
3
4
5
6
cellranger aggr \
  --id=liver_combined_analysis \
  --csv=aggregation.csv \
  --normalize=mapped \
  --localcores=16 \
  --localmem=180

生成一个新的文件夹 liver_combined_analysis，其中包含合并后的 filtered_feature_bc_matrix 文件夹。

ps:在此推荐两个小工具了，下载和安装就百度吧，基本上都有教程的

  1. winscp用于本地到服务器上的文件传输（可视化，很方便，并且也可以在这里增删改查文件）
  2. finalshell可以复制粘贴代码（学校集群cv不了），并且旁边有当前内存，储存空间，cpu占用等有用的信息
许可协议

本文采用 署名-非商业性使用-相同方式共享 4.0 国际 许可协议，转载请注明出处。

站内搜索

没有找到内容！