美文网首页全基因组/外显子组测序分析
7 比对到参考基因组输出bam文件

7 比对到参考基因组输出bam文件

作者: Y大宽 | 来源:发表于2019-06-11 23:41 被阅读33次

    总目录:三阴性乳腺癌全外显子分析(wes)


    接下来用 BWA mem把fastq map到参考基因组 hg38 版本。
    比对结果直接通过管道传给samtools处理,节省 I/O 时间。
    因为空间问题,比对好的文件放在
    /project/align/wes目录

    6.1设置好下面批量比对的数据文件

    kelly/wesproject/4_clean/wes目录下,也可以在align/wes目录下写完整路径

    ls *1.fq.gz> 1
    ls *2.fq.gz> 2
    paste 1 2 > config
    #vim config 写入第一列样本名,要以Tab分开
    cat 1|cut -d"_" -f 2,3 1>0
    paste 0 1 2 > config
    

    6.2 比对

    align/wes目录下
    根据前面的经验,一次并行比对10个文件

    (wes) pc@lab-pc:/home/kelly/wesproject/4_clean/wes$ cat config|head -10> config_10
    
    INDEX=/data/bigbiosoft/GATK/resources/bundle/hg38/bwa_index/gatk_hg38
    cat /home/kelly/wesproject/4_clean/wes/config10|while read id
    do 
    arr=($id)
    sample=${arr[0]}
    fq1=${arr[1]}
    fq2=${arr[2]}
    echo $sample $fq1 $fq2
    bwa mem -t 4 -R "@RG\tID:$sample/tSM:$sample\tLB:WGS\tPL:Illumina" $INDEX /home/kelly/wesproject/4_clean/wes/$fq1 /home/kelly/wesproject/4_clean/wes/$fq2 |samtools sort -@ 4 -o $sample.bam -  &
    done
    

    注意bam命令的-R参数,不加也可以运行,但是后面的gatk时会报错,但是也有解决办法,见后面。-t 和-@是线程数

    非常耗时的比对工作今天终于完成了。20190627

    相关文章

      网友评论

        本文标题:7 比对到参考基因组输出bam文件

        本文链接:https://www.haomeiwen.com/subject/gzjhfctx.html