美文网首页
【Read质量】Read质量值计算

【Read质量】Read质量值计算

作者: adell898 | 来源:发表于2023-07-11 11:51 被阅读0次

对测序fastq数据,碱基质量值(base quality)的本质是体现测序错误率,在fastq文件中,碱基质量值以字符形式存储在fastq文件中。字符的碱基质量值体系通常是Phred33,即碱基质量值Q = 字符的ASCII码– 33。而Q和碱基的错误率之间又具有如下对应关系:

Q = -10logP

Q:碱基质量值;

P:碱基测序错误率;

同理,Read的质量值表示的是Read中碱基平均测序错误率,基于此,Reads的质量值和错误率应具有如下关系,

Q_r=-10logP_r

Qr:Read碱基质量值;

Pr:Reads碱基平均测序错误率;

由此,当需要计算Read的质量值时,其计算步骤如下,

[if !supportLists]1) [endif]获得Read各碱基的测序错误率Pi

P_i={10}^{-\ \frac{Qi}{\ 10}}

Pi:各碱基错误率;

Qi:各碱基质量值;

i = 1 .. len(Read)

[if !supportLists]1) [endif]计算Read碱基平均测序错误率

P_r=\frac{\sum_{1}^{n}P_i}{n}

n:read长度;

Pi:各碱基错误率

[if !supportLists]1)[endif]计算Read质量值

Q_r=-10logP_r

注:

当计算Read质量值时,如果直接用base的质量值(即Q值)求均值,获得的结果是不对的。

```python

import os,sys

import numpy as np

def cal_p(qual_str):

  Q = ord(qual_str) - 33

  p = 10 ** (-Q/10)

  reutrn p

qual = '####???'

p_list = np.array([])

for q  in qual:

  p = cal_p(q)

  np.append(p_list,p)

mean_p = p_list.mean()

read_Q = -10np.log10(mean_p)

print('Reads Qual: ')

print(read_Q)

```

参考资料:

https://labs.epi2me.io/quality-scores/

相关文章

网友评论

      本文标题:【Read质量】Read质量值计算

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