CN103390119B - 一种转录因子结合位点识别方法 - Google Patents

Abstract

本发明属于分子生物信息检测领域，具体涉及一种基于条件随机场技术，融合ChIP-chip基因芯片数据和ChIP-seq？DNA测序数据的转录因子结合位点识别方法。本发明包括：建立条件随机场模型；获取ChIP-chip实验检测值，识别出对应状态值；获取ChIP-seq实验检测值识别出对应状态值；测试条件随机场模型的识别精度；加权融合识别结果的第n个DNA碱基片段识别概率；比较和识别转录因子结合位点。本发明利用条件随机场融合ChIP-chip和ChIP-seq实验的检测数据识别转录因子结合位点。通过实验验证，在识别准确率方面，本发明方法要高于采用单一技术的识别方法。

Description

一种转录因子结合位点识别方法

技术领域

本发明属于分子生物信息检测领域，具体涉及一种基于条件随机场技术，融合ChIP-chip基因芯片数据和ChIP-seqDNA测序数据的转录因子结合位点识别方法。

背景技术

转录因子作为调控基因转录的功能蛋白，其在DNA基因启动子区域上结合位点的识别对于基因转录调控机制的研究起着很大作用，它是生命科学研究中的重点和难点。随着信息科学实验技术的进步，用于获取分子生物信息的实验手段也日新月异。其中，在全基因组范围内检测DNA区域特定蛋白质结合状态的标准实验是染色质免疫沉淀反应(Chromatinimmunoprecipitaion,ChIP)，染色质免疫沉淀反应(ChIP)分别与基因芯片或者DNA测序两种检测技术相结合产生的ChIP-chip和ChIP-seq技术则是检测染色质免疫沉淀反应实验结果的两大主要技术。

由于检测原理的不同，ChIP-chip和ChIP-seq实验技术所产生的实验数据存在着各自的优缺点。基因芯片是对生物信息进行直接测量，它的敏感度较高，但是相比于高通量DNA测序，它的测试精度不高，测试特异性较低；高通量DNA测序则具有较高的检测分辨率，它的测试精度较高，测试特异性较好，但由于是对生物信息进行间接测量，它的敏感度较低。可见，分别通过ChIP-chip和ChIP-seq不同的实验途径获取的实验数据虽然具有各自不同的优缺点，但具有明显的互补性，因此，融合这两种技术获得的检测数据能够得到更为完整的检测信息。

作为针对一维序列的模式识别问题，分析当前DNA碱基片段是否是特定功能蛋白的结合位点，不仅要考虑当前碱基片段的生化特性是否适合特定功能蛋白结合，还需考虑其两侧多个碱基片段的生化特性对特定功能蛋白结合能力的影响。具体到转录因子结合位点的识别问题，就是在识别过程中应同时考虑结合位点待识别碱基片段及其附近多个碱基片段的特性。在众多针对一维序列的模式识别技术中，由于条件随机场(ConditionalRandomFields,CRF)技术具有特征选择无方向约束，特征位置支持长程关联等特性，使得它非常适合转录因子结合位点的识别。

发明内容

本发明的目的是提供一种更高精度的基于条件随机场的融合基因芯片数据和DNA测序数据的转录因子结合位点识别方法。

本发明的目的是这样实现的：

（1）建立条件随机场模型：

$p\left(y\right|x)=\frac{1}{Z\left(x\right)}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

$Z\left(x\right)=\underset{y}{\mathrm{\Σ}}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

其中，x={x₁,x₂,…,x_n}表示DNA碱基片段的实验检测值；y={y₁,y₂,…,y_n}是DNA碱基片段的对应状态值，1表示是转录因子结合位点，0表示不是转录因子结合位点；t_k(y_i-1,y_i,x,i)表示第i个碱基片段在当前实验检测值序列x下，状态为y_i的第i个碱基片段和状态为y_i-1的第i-1个碱基片段之间的转移特征函数；s_l(y_i,x,i)表示在当前实验检测值序列x下，第i个碱基片段的状态是y_i的状态特征函数；λ_k和μ_l分别是t_k(y_i-1,y_i,x,i)和s_l(y_i,x,i)对应的权值，表示各特征函数的重要性；Z(x)是规范化因子，使得p(y|x)位于[0,1]之间；

（2）获取DNA碱基片段的ChIP-chip实验检测值x={x₁,x₂,…,x_n}，根据条件随机场模型，识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}；

（3）获取DNA碱基片段的ChIP-seq实验检测值x={x₁,x₂,…,x_n}，根据条件随机场模型，识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}；

（3）测试条件随机场模型的识别精度：

$S_n=\frac{\mathrm{TP}}{\mathrm{TP}+\mathrm{FN}},$

$S_P=\frac{\mathrm{TN}}{\mathrm{TN}+\mathrm{FP}},$

$A_c=\frac{S_n+S_p}{2},$

其中，敏感性S_n，特异性S_p，准确率A_c，TP代表转录因子结合位点被预测正确的个数；FN表示转录因子结合位点被预测错误的个数；TN表示非转录因子结合位点被预测正确的个数；FP表示非转录因子结合位点被预测错误的个数；

（4）对第n个DNA碱基片段，通过ChIP-chip实验检测值识别为转录因子结合位点的概率用表示，识别为非转录因子结合位点的概率用表示，通过ChIP-seq实验检测值识别为转录因子结合位点的概率用表示，识别为非转录因子结合位点的概率用表示，加权融合识别结果的第n个DNA碱基片段识别概率表示为：

$p_1^{\left(n\right)}=p_{\mathrm{lchip}}^{\left(n\right)}\×w_1+p_{1\mathrm{seq}}^{\left(n\right)}\×w_2,$

$p_0^{\left(n\right)}=p_{0\mathrm{chip}}^{\left(n\right)}\×w_1+p_{0\mathrm{seq}}^{\left(n\right)}\×w_2,$

融合权值w1和w2为与识别结果的识别准确率，w₁+w₂=1；

（5）比较和则此碱基片段识别为转录因子结合位点；则此碱基片段识别为转录因子非结合位点。

本发明的有益效果在于：

本发明利用条件随机场融合ChIP-chip和ChIP-seq实验的检测数据识别转录因子结合位点。通过实验验证，在识别准确率方面，本发明方法要高于采用单一技术的识别方法。

附图说明

图1为本发明的转录因子结合位点识别方法流程图；

图2为条件随机场模型的训练流程图。

具体实施方式

下面结合附图对本发明做进一步描述：

基于条件随机场的融合基因芯片数据和DNA测序数据的转录因子结合位点预测方法，包括下列步骤：

(1)采用条件随机场识别技术采用ChIP-chip基因芯片数据，对转录因子结合位点进行初步识别；

(2)采用条件随机场识别技术采用ChIP-seqDNA测序数据，对转录因子结合位点进行初步识别；

(3)采用决策级融合的方式，对转录因子结合位点进行最终识别。即对同一位置的碱基片段，加权融合分别来自ChIP-chip基因芯片数据和ChIP-seqDNA测序数据的转录因子结合位点初步识别结果，并得出最终的识别结果。其中，融合权值与初步识别结果的识别准确率成正比，且两个融合权值的和为1。

本方法所采用的条件随机场的表示形式如下：

$p\left(y\right|x)=\frac{1}{Z\left(x\right)}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

$Z\left(x\right)=\underset{y}{\mathrm{\Σ}}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

其中，x={x₁,x₂,…,x_n}表示DNA碱基片段的实验检测值；y={y₁,y₂,…,y_n}是DNA碱基片段的对应状态值，这里只有两种状态值，1表示是转录因子结合位点，0表示不是转录因子结合位点；t_k(y_i-1,y_i,x,i)表示第i个碱基片段在当前实验检测值序列x下，状态为y_i的第i个碱基片段和状态为y_i-1的第i-1个碱基片段之间的转移特征函数；s_l(y_i,x,i)表示在当前实验检测值序列x下，第i个碱基片段的状态是y_i的状态特征函数；λ_k和μ_l分别是t_k(y_i-1,y_i,x,i)和s_l(y_i,x,i)对应的权值，表示各特征函数的重要性；Z(x)是规范化因子，使得p(y|x)位于[0,1]之间。本发明所选取的条件随机场特征函数包含状态特征函数和转移特征函数。其中，状态特征选择窗口为1，即结合当前碱基片段和前后各一个碱基片段的检测值构成特征。转移特征简单取为相邻的片段特征为1，否则为0。特征函数则选为每个特征的自身、平方和开方，最后由权重来决定每个特征函数的重要性。采用上述模型，当获取DNA碱基片段的实验检测值x={x₁,x₂,…,x_n}后，即可识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}，用于判断每个碱基片段是否是转录因子结合位点。

本发明方法具体实施方式如图1所示，具体如下：

1.碱基片段ChIP-chip数据CRF模型转录因子结合位点初步识别结果

碱基片段ChIP-chip数据转录因子结合位点初步识别的条件随机场CRF模型设计如下：

$p\left(y\right|x)=\frac{1}{Z\left(x\right)}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

$Z\left(x\right)=\underset{y}{\mathrm{\Σ}}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

其中，x={x₁,x₂,…,x_n}表示DNA碱基片段的ChIP-chip实验检测值；y={y₁,y₂,…,y_n}是DNA碱基片段的对应状态值，这里只有两种状态值，1表示是转录因子结合位点，0表示不是转录因子结合位点；t_k(y_i-1,y_i,x,i)表示第i个碱基片段在当前实验检测值序列x下，状态为y_i的第i个碱基片段和状态为y_i-1的第i-1个碱基片段之间的转移特征函数；s_l(y_i,x,i)表示在当前实验检测值序列x下，第i个碱基片段的状态是y_i的状态特征函数；λ_k和μ_l分别是t_k(y_i-1,y_i,x,i)和s_l(y_i,x,i)对应的权值，表示各特征函数的重要性；Z(x)是规范化因子，使得p(y|x)位于[0,1]之间。本发明所选取的条件随机场特征函数包含状态特征函数和转移特征函数。其中，状态特征选择窗口为1，即结合当前碱基片段和前后各一个碱基片段的检测值构成特征。转移特征简单取为相邻的片段特征为1，否则为0。特征函数则选为每个特征数据的自身、平方和开方，最后由权重来决定每个特征函数的重要性。

采用上述模型，当获取DNA碱基片段的ChIP-chip实验检测值x={x₁,x₂,…,x_n}后，即可识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}，作为每个碱基片段转录因子结合位点初步识别结果。

2.碱基片段ChIP-seq数据CRF模型转录因子结合位点初步识别结果

碱基片段ChIP-seq数据转录因子结合位点初步识别的条件随机场CRF模型设计如下：

$p\left(y\right|x)=\frac{1}{Z\left(x\right)}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

$Z\left(x\right)=\underset{y}{\mathrm{\Σ}}\exp (\underset{i,k}{\mathrm{\Σ}}{\mathrm{\λ}}_k{t}_k(y_{i-1},y_i,x,i)+\underset{i,l}{\mathrm{\Σ}}{\mathrm{\μ}}_l{s}_l(y_i,x,i))$

其中，x={x₁,x₂,…,x_n}表示DNA碱基片段的ChIP-seq实验检测值；y={y₁,y₂,…,y_n}是DNA碱基片段的对应状态值，这里只有两种状态值，1表示是转录因子结合位点，0表示不是转录因子结合位点；t_k(y_i-1,y_i,x,i)表示第i个碱基片段在当前实验检测值序列x下，状态为y_i的第i个碱基片段和状态为y_i-1的第i-1个碱基片段之间的转移特征函数；s_l(y_i,x,i)表示在当前实验检测值序列x下，第i个碱基片段的状态是y_i的状态特征函数；λ_k和μ_l分别是t_k(y_i-1,y_i,x,i)和s_l(y_i,x,i)对应的权值，表示各特征函数的重要性；Z(x)是规范化因子，使得p(y|x)位于[0,1]之间。本发明所选取的条件随机场特征函数包含状态特征函数和转移特征函数。其中，状态特征选择窗口为1，即结合当前碱基片段和前后各一个碱基片段的检测值构成特征。转移特征简单取为相邻的片段特征为1，否则为0。特征函数则选为每个特征数据的自身、平方和开方，最后由权重来决定每个特征函数的重要性。

采用上述模型，当获取DNA碱基片段的ChIP-seq实验检测值x={x₁,x₂,…,x_n}后，即可识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}，作为每个碱基片段转录因子结合位点初步识别结果。

3.转录因子结合位点的融合识别结果

本发明采用决策级融合的方式获取最终的转录因子结合位点识别结果。即对同一位置的碱基片段，采用加权的方式融合分别来自ChIP-chip基因芯片数据和ChIP-seqDNA测序数据的碱基片段CRF模型转录因子结合位点初步识别结果，并得出最终的转录因子结合位点识别结果。具体为：

首先，采用交叉验证的方式，对条件随机场训练模型的识别精度进行测试。其中，测试指标采用敏感性S_n(Sensitivity)，特异性S_p(Specificity)，准确率A_c，它们定义如下：

$S_n=\frac{\mathrm{TP}}{\mathrm{TP}+\mathrm{FN}}---\left(1\right)$

$S_P=\frac{\mathrm{TN}}{\mathrm{TN}+\mathrm{FP}}---\left(2\right)$

$A_c=\frac{S_n+S_p}{2}---\left(3\right)$

式中，TP代表转录因子结合位点被预测正确的个数；FN表示转录因子结合位点被预测错误的个数；TN表示非转录因子结合位点被预测正确的个数；FP表示非转录因子结合位点被预测错误的个数。条件随机场训练模型的准确率A_c将用于融合权值的求取。

其次，针对第n个DNA碱基片段，基于ChIP-chip数据初步识别为转录因子结合位点的概率用表示，初步识别为非转录因子结合位点的概率用表示；基于ChIP-seq数据初步识别为转录因子结合位点的概率用表示，初步识别为非转录因子结合位点的概率用表示。则融合上述初步识别结果的该碱基片段识别概率表示为：

$p_1^{\left(n\right)}=p_{\mathrm{lchip}}^{\left(n\right)}\×w_1+p_{1\mathrm{seq}}^{\left(n\right)}\×w_2---\left(4\right)$

$p_0^{\left(n\right)}=p_{0\mathrm{chip}}^{\left(n\right)}\×w_1+p_{0\mathrm{seq}}^{\left(n\right)}\×w_2,$

式中，融合权值w₁和w₂设计为与初步识别结果的识别准确率，即条件随机场训练模型的准确率A_c成正比，且两个融合权值的和为1。

最后，比较和假如则此碱基片段识别为转录因子结合位点；假如则此碱基片段识别为转录因子非结合位点。

4.实验验证

4.1数据预处理

4.1.1ChIP-chip数据预处理

本实验所采用的ChIP-chip基因芯片数据来源于国际大型公共生物信息数据库GEO数据库(GSE6892)，该数据检测了人体44个基因上转录因子STAT1的结合状态。每一个测得数据为50bp长，但由于采用了瓦片检测形式，相邻两个测得数据间存在部分重叠。我们对数据进行重新分割和组合，最终获得38bp的相连的碱基片段，并计算了每个片段上的ChIP-chip检测值。

4.1.2ChIP-seq数据预处理

本实验所采用的ChIP-seq测序数据来源于BCGSC数据库，也是针对人体转录因子STAT1结合状态的检测。在对ChIP-seq的数据进行预处理的时候，首先把原始数据转换成fastq格式，然后通过比对，将测得数据映射到人体基因组上，并计算出DNA上每一个碱基位点的测得数据覆盖次数，最后，对应ChIP-chip数据所获得的DNA碱基片段，计算出各碱基片段的ChIP-seq测得数据覆盖次数值。

4.2条件随机场模型训练

条件随机场模型的特征函数包含状态特征函数和转移特征函数。其中，状态特征选择窗口为1，即结合当前碱基片段和前后各一个碱基片段的检测值构成特征。转移特征简单取为相邻的片段特征为1，否则为0。特征函数则选为每个特征的自身、平方和开方。利用已知的DNA碱基片段上STAT1结合状态数据，对条件随机场模型进行训练，训练流程图如图2所示。

4.3条件随机场模型识别结果的融合

首先，采用5次交叉验证的方式，对条件随机场训练模型的识别精度进行测试。测试结果如表1和表2所示：

表1ChIP-chip数据交叉验证结果

表2ChIP-seq数据交叉验证结果

表1中，用五次交叉验证法测得采用ChIP-chip数据训练的CRF模型的平均准确率为57.68%。，表2中，用五次交叉验证法测得采用ChIP-seq数据训练的CRF模型的平均准确率为63.76%。

然后，根据分别基于ChIP-chip数据和ChIP-seq数据的CRF模型训练准确率计算融合权值为：w1=0.475，w2=0.525。并据此计算第n个DNA碱基片段为转录因子结合位点的概率和为非转录因子结合位点的概率

最后，比较和假如则第n个DNA碱基片段识别为转录因子STAT1的结合位点；假如则第n个DNA碱基片段识别为转录因子STAT1的非结合位点。采用五次交叉验证法对融合结果的精度测试如表3所示：

表3融合的交叉验证结果

可见，融合ChIP-chip数据和ChIP-seq数据得到的识别结果的平均准确率是0.7697，大幅高于单独使用ChIP-chip数据和ChIP-seq数据得到的识别结果的平均准确率。

Claims (1)

1.一种转录因子结合位点识别方法，其特征在于：

（1）建立条件随机场模型：

其中，x={x₁,x₂,…,x_n}表示DNA碱基片段的实验检测值；y={y₁,y₂,…,y_n}是DNA碱基片段的对应状态值，1表示是转录因子结合位点，0表示不是转录因子结合位点；t_k(y_i-1,y_i,x,i)表示第i个碱基片段在当前实验检测值序列x下，状态为y_i的第i个碱基片段和状态为y_i-1的第i-1个碱基片段之间的转移特征函数；s_l(y_i,x,i)表示在当前实验检测值序列x下，第i个碱基片段的状态是yi的状态特征函数；λ_k和μ_l分别是t_k(y_i-1,y_i,x,i)和s_l(y_i,x,i)对应的权值，表示各特征函数的重要性；Z(x)是规范化因子，使得p(y|x)位于[0,1]之间；

（2）获取DNA碱基片段的ChIP-chip实验检测值x={x₁,x₂,…,x_n}，根据条件随机场模型，识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}；

（3）获取DNA碱基片段的ChIP-seq实验检测值x={x₁,x₂,…,x_n}，根据条件随机场模型，识别出对应的DNA碱基片段的状态值y={y₁,y₂,…,y_n}；

（4）测试条件随机场模型的识别精度：

其中，敏感性S_n，特异性S_p，准确率A_c，TP代表转录因子结合位点被预测正确的个数；FN表示转录因子结合位点被预测错误的个数；TN表示非转录因子结合位点被预测正确的个数；FP表示非转录因子结合位点被预测错误的个数；

（5）对第n个DNA碱基片段，通过ChIP-chip实验检测值识别为转录因子结合位点的概率用表示，识别为非转录因子结合位点的概率用表示，通过ChIP-seq实验检测值识别为转录因子结合位点的概率用表示，识别为非转录因子结合位点的概率用表示，加权融合识别结果的第n个DNA碱基片段识别概率表示为：

融合权值w₁和w₂为与识别结果的识别准确率，w₁+w₂=1；

（6）比较和则此碱基片段识别为转录因子结合位点；则此碱基片段识别为转录因子非结合位点。