转录组入门学习(三)

学习——数据下载及质控

数据预处理

1. 准备工作
1.1 构建项目目录
  • “基因组文件”与“注释信息文件”放在同一个目录中
  • 原始数据目录
  • 结果文件目录:包括“比对结果”与“定量结果”


    阿里云终端截图,不能安装tree包
#解决办法
sudo apt-get update #进行升级
#之后安装tree
apt install tree
#但是当时没考虑用conda直接安装tree
#安装完毕后,运行tree时发现有乱码,解决方法:
alias tree='tree --charset ASCII'
解决之后的tree命令状态
1.2 参考序列下载
  • 参考基因组 fasta
  • 注释信息 gtf/gff
1.3 原始数据上传
  • 检查数据完整性:md5值
  • 给自己的文件生成md5值:
  • 比对已有的md5值:
md5sum *gz > md5.txt    #生成md5.txt文件
md5sum -c md5.txt       #检查md5文件的一致性
2. 质量控制:两个软件
2.1 FastQC
  • 安装:conda
#FastQC为java软件,需要预先配置Linux java 环境
getconf LONG_BIT  #查看系统是32位还是64位
wget -c http://download.oracle.com/otn-pub/java/jdk/8u151-b12/e758a0de34e24606bca991d704f6dcbf/jdk-8u151-linux-x64.tar.gz
tar zxf jdk-8u151-linux-x64.tar.gz  #解压缩
mkdir java
mv ./jdk1.8.0_151/ ./java  #将解压缩的文件移入新文件夹

#在profile中加入下面的内容
vi .profile   
JAVA_HOME=/usr/java/jdk1.8.0_151
export JRE_HOME=/usr/java/jdk1.8.0_151/jre
export CLASSPATH=.:$JAVA_HOME/lib:$JRE_HOME/lib:$CLASSPATH
export PATH=$JAVA_HOME/bin:$JRE_HOME/bin:$PATH

reboot   #重启计算机
  • 如何批量处理多个样品
#方法1:
fastqc sample*gz     #fastq文件不需要解压
rm -f *html    #删除所有html文件(-f (--force): 强制删除文件或目录,没有提示)
rm -f *zip     #删除所有zip文件

#方法2:
for i in 'ls *gz'; do fastqc $i; done #依次遍历以gz结尾的文件,并进行fastqc命令

#方法3:
ls *gz | xargs -I [] echo 'nohup fastqc [] &' > fastqc.sh
#并行处理多个样本
#解释:
#ls *gz |: 列出文件列表,并通过管道操作符传入后面
#xargs -I: xargs命令是给其他命令传递参数的一个过滤器,也是组合多个命令的一个工具。把前一个命令的输出结果当做后一个命令的输入文件。
#它擅长将标准输入数据转换成命令行参数,xargs能够处理管道或者stdin并将其转换成特定命令的命令参数。
#[]: 替代符
2.2 MultiQC
  • 安装:conda
  • 特点:
    1. 自动检测已有输出结果
    2. 支持多种结果的整合
  • 同时展示多个结果文件
3. 质量过滤
  • 软件:Trimmomatic
  • 特点:
    1. 支持多线程,处理数据速度快
    2. 主要用来去除 Illumina 平台接头
    3. 根据碱基质量值对 fastq 进行筛选
    4. 支持 SE 和 PE 测序数据,支持 gzip 和 bzip2 压缩文件
  • 过滤依据
    1. Illumianclip: 过滤reads 中的illumina 接头
    2. LEADING: 从 reads 开头切除质量值低于阈值的碱基
    3. TRAILING: 从 reads 末尾切除质量值低于阈值的碱基
    4. SLIDINGWINDOW: 从 reads 的 5' 端开始,进行滑窗过滤,切掉碱基质量平均值低于阈值的滑窗
    5. MINLEN: 丢弃经过剪切后长度低于阈值的这条 reads
    6. TOPHRED33: 将 reads 的碱基质量体系转为 phred-33
    7. TOPHRED64: 将 reads 的碱基质量体系转为 phred-64
  • 使用
    1. 接头序列的选择:
      "Illumian Single End" / "Illumian Paired End": "TruSeq2-SE.fa" and "TruSeq2-PE.fa"
      "TruSeq Universal Adapter" / "TruSeq Adapter, Index ...": "TruSeq3-SE.fa" and "TruSeq3-PE.fa"
    2. 去接头参数的选择: true; false
#参数:
trimmomatic PE -threads 4 \
sample1_R1.fastq.gz sample1_R2.fastq.gz \
../02clean_data/sample1_paired_clean_R1.fastq.gz \
../02clean_data/sample1_unpair_clean_R1.fastq.gz \
../02clean_data/sample1_paired_clean_R2.fastq.gz \
../02clean_data/sample1_unpair_clean_R2.fastq.gz \
ILLUMINACLIP:/home/leon/miniconda2/share/trimmomatic-0.36-5/adapters/TruSeq3-PE-2.fa:2:30:10:1:true \
LEADING:3 TRAILING:3 \
SLIDINGWINDOW:4:20 MINLEN:50 TOPHRED33

trimmomatic PE -threads 4 \
sample2_R1.fastq.gz sample2_R2.fastq.gz \
../02clean_data/sample2_paired_clean_R1.fastq.gz \
../02clean_data/sample2_unpair_clean_R1.fastq.gz \
../02clean_data/sample2_paired_clean_R2.fastq.gz \
../02clean_data/sample2_unpair_clean_R2.fastq.gz \
ILLUMINACLIP:/home/leon/miniconda2/share/trimmomatic-0.36-5/adapters/TruSeq3-PE-2.fa:2:30:10:1:true \
LEADING:3 TRAILING:3 \
SLIDINGWINDOW:4:20 MINLEN:50 TOPHRED33

#注意:冒号 ':' 之后没有空格,否则会报错
sample1运行结果
sample2运行结果