二代测序均一度(Uniformity)的计算方法

均一度是指测序的reads在目标区域的测序深度均匀程度的度量。

目前用的较多的是大于平均深度的20% 的碱基位点占目标区域碱基位点总数的比例。

可以使用IGV去看测序reads在 hg19 基因组上的覆盖情况， 这样能更直观的看出均一度是怎么计算的。

上图是用IGV定位到 8号染色体的 145743076bp 附近

可以看出在Coverage 这个显示窗口，大部分位点的测序深度是在100 - 422 之间， 其中一个位点chr8：145742550 的测序深度是422， 有420个G， 2个A

该窗口内的reads 长度为绝大多数为150bp， 在RECQL4 的这几个外显子上覆盖的比较均匀。

计算均一度的方法，首先需要求出每个位点测序深度。

以下针对两种类型的测序数据 提供计算方法。

如果是杂交捕获的方法用picard 的CalculateHsMetrics 工具

java -jar /home/software/picard.jar CalculateHsMetrics \
R=/data/reference/hg19.2bit/hg19.fasta \
I=/data/testSample.dup.bam \
O=out.target.stat.xls \
BI=Target_gene.bed \
TI=Target_gene.bed \
PER_TARGET_COVERAGE=out.target.coverage \
PER_BASE_COVERAGE=testSample.base.coverage \
MQ=20 \
Q=20 \
COVMAX=100000

求出平均测序深度

sed '1d' testSample.base.coverage | awk '{sum+=$4} END {print "Average = ", sum/NR}'

计算大于平均深度的20%的 bases number （注意：ave_depth 需要自己填写）

grep -v 'chrom' testSample.base.coverage | awk '{if($4>0.2* ${ave_depth}){print}}' | wc -l

这个bases 个数除以 target区域的总bases 数，然后再乘以100%，就是Uniformity的值

如果是amplicon获得的测序数据，用samtools或者bedtools都可以计算

samtools depth -d 100000 -b Target_gene.bed testSample.dup.bam -o test_Sample_Depth.txt
coverageBed -b testSample.dup.bam -a Target_gene.bed -d > res_Sample_cov.txt

结果如下：

chr1    198866926 198867539 MIR181A1HG        1       3569
chr1    198866926 198867539 MIR181A1HG        2       3572
chr1    198866926 198867539 MIR181A1HG        3       3574
chr1    198866926 198867539 MIR181A1HG        4       3584
chr1    198866926 198867539 MIR181A1HG        5       3586
......
chr1    198866926 198867539 MIR181A1HG        15      3573
chr1    198866926 198867539 MIR181A1HG        16      3576
chr1    198866926 198867539 MIR181A1HG        17      3581

第一列是染色体名，

第二列和第三列是bed文件中各个区域的起点和终点，

第四列是基因名，

第五列是相应区域内的坐标索引，

最后一列是测序深度。

用这个结果可以进一步计算平均测序深度，以及Uniformity了。