Published at 2020-05-19 10:48
Author：zhixy
View：1817
基因组de novo测序，没有参考序列，如要计算测序深度，可以将reads回贴到拼接好的contigs/scaffolds序列上，近而统计每个核酸位点上被reads覆盖的次数，即测序深度。
回贴mapping短序列reads到基因组上，成熟的软件有bwa、bowtie2等，下面以bwa为例介绍一种测序深度估算方法。
(base) [user@server ~]# bwa index contigs.fa
(base) [user@server ~]# bwa mem -t 60 -M contigs.fa strain_trim1.fastq strain_trim2.fastq > bwa.sam
-t 60，计算核心/线程数
(base) [user@server ~]# samtools view -@ 60 -bS bwa.sam > bwa.bam
-@ 60，计算核心/线程数
(base) [user@server ~]# samtools sort -@ 60 bwa.bam > bwa.sort.bam
(base) [user@server ~]# samtools depth bwa.sort.bam > bwa.depth
(base) [user@server ~]# cat bwa.depth | awk '{ sum ＋= ＄3; } END { print "coverage = " sum/NR }'
对第三列（$3）求和，并除以行数（NR），行数即基因组序列长度。
