2021.01.07丨使用fastp统计样品质量结果

• 各位小伙伴在对测序样品进行质控的时候，首选基本上都是fastQC，他能能够生成许多图片直观地展示质控结果。
  • 
• 然而，当我们有多个样品，希望对其结果以表格形式进行展示的时候，fastQC能提供的信息就比较少了，比如GC含量精确到小数点，或者Q30等等
  • fastQC能统计到的基本信息
  • 
  • 我们希望得到的统计结果
    • 
• 那么如何能够批量统计到更详细的质控信息呢？fastp工具和这篇文章脚本的必要性就产生了，它可以统计测序数据较多的信息并以.json形式进行展示。
  • 我们用Editplus打开fastp.json文件，可以看到fastp对测序数据有详细的统计情况
    • 
• 数据整理流程
  • step.1下载并安装fastp

    • conda install fastp

       

  • step.2使用fastp获取统计文件

    • for i in *_val_1.fq.gz; #遍历所有样品R1测序数据
      do
      i=${i%_val_1.fq.gz*};#获取样品名称
      fastp -w 8 -i ${i}_val_1.fq.gz -o ${i}_1 -I ${i}_val_2.fq.gz -O ${i}_2 -j ${i}.json #
      rm -rf ${i}_1 ${i}_2 #删除指控后结果（之前没做过trim可以保留）
      done

       

  • step.3获取统计结果

    • import re
      import os
      newfile_name = 'fastp_stat'
      desktop_path = '/home/yangxin/project/rna_seq/refer_rna_seq/6_8_CGRSH201026_20201215/01.data/trim_QC/'
      file_path = desktop_path+newfile_name+'.txt'
      newfile = open(file_path,'w')
      for root,dirs,files in os.walk(r"/home/yangxin/project/rna_seq/refer_rna_seq/6_8_CGRSH201026_20201215/01.data/trim_QC"):
      for file in files:
      if ".json" in file:
      #获取文件路径
      print(os.path.join(root,file))
      fastp_file = open(os.path.join(root,file),'r')
      genome_line = fastp_file.readlines()
      total_reads = re.sub("\D","",genome_line[14])
      total_bases = re.sub("\D","",genome_line[15])
      Q30_tmp_rates = re.sub("\D","",genome_line[19])
      Q30_rates = Q30_tmp_rates[2:]
      Q30_rates = int(Q30_rates) *0.0001
      GC_connects = re.sub("\D","",genome_line[22])
      GC_connects = int(GC_connects) *0.0001
      file_name = file.replace('.json','')
      result_stat = file_name+"\t"+total_reads+"\t"+total_bases+"\t"+'%.2f'%Q30_rates+"\t"+'%.2f'%GC_connects+"\n"
      print(result_stat)
      newfile.write(result_stat)
      newfile.close()
      fastp_file.close()

       

  • 文件结果展示
    • 
• 最后直接复制粘贴到RNA-seq的表格模板上面即可。