CN110556184B - 基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法 - Google Patents

Abstract

本发明公开了一种基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法，属于系统生物学领域，主要包括三个步骤：步骤一，分别计算非编码RNA的高斯谱核相似性和疾病的高斯谱核相似性；步骤二，使用迭代求解算法计算非编码RNA‑疾病关联对预测得分；步骤三，根据算出的非编码RNA‑疾病关联对得分进行排序，给出最终预测结果。本发明通过Hessian正则化细致刻画数据内在流形结构，有效利用了阴性样本的信息；l2,1范数约束和近似正交约束保证了编码矩阵的群稀疏性，能减弱噪声数据的影响。本方法能获得比较可靠的预测结果，有效地解决了生物实验方法耗时长、成本高的问题。

Description

基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测 方法

技术领域

本发明涉及系统生物学领域，更具体地说，本发明涉及一种基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法。

背景技术

非编码RNA(Non-coding RNA，ncRNA)指转录组中不编码蛋白质的RNA分子，常见的包括microRNA、lncRNA和circRNA等。

微小RNA(miRNA)是长度约为22个核苷酸的内源性单链RNA，存在于多种物种中，包括植物，动物和某些病毒。作为一种重要的转录后调节因子，它们通过与靶RNA的3'非翻译区(untranslated regions,UTR)碱基配对而抑制基因表达并促进mRNA的降解。它们在多种生物学过程中发挥着关键作用，如细胞分裂，分化，发育，代谢，感染，衰老，细胞凋亡和信号转导。实验证据表明，miRNA的异常表达与诸多人类疾病有关。例如，miRNA 181a的上调表达可能引发人类1型糖尿病的进展。此外，高胆固醇血症与动脉粥样硬化小鼠肝脏miR-223水平的增加密切相关。此外，已有证明表明miR-21，miR-494和miR-1973是经典霍奇金淋巴瘤中的疾病应答生物标志物。

长链非编码RNA(long noncoding RNA，lncRNA)是一类长度＞200bp的RNA，参与多种生物学过程的调控，包括基因组表观遗传修饰、调控转录后翻译、增强子RNA作用等，从而对细胞的增殖、分化、迁移、凋亡、免疫等发挥调控作用。有实验表明，lncRNA AC006449.2在卵巢癌细胞中可能发挥抑癌因子的作用。此外，lncRNA H19高表达的肝癌细胞以外泌体方式,增强邻近肝癌细胞的增殖、迁移和侵袭能力,促进肝癌的发生与发展。通过大数据分析发现，lncRNA RP11-214F16.8在乳腺癌中高表达，且促进乳腺癌细胞增殖，进而推动乳腺癌进程。

环状RNA(circRNA)是一类由反向剪接形成的、没有5’端帽子和3’端多聚腺苷酸尾巴的环状闭合RNA分子，具有保守、稳定、组织特异性、时空特异性等特点。大量研究发现它可以通过多种机制参与动物生长发育调控，疾病等的发生和发展。研究发现，在ISO诱导的心肌肥大小鼠中强制表达circRNA HRCR可显著缓解心肌肥大。实验发现，circRNA Cdr1as会影响胰岛素分泌和胰岛B细胞的更新。结直肠癌相关研究表明，hsa_circ_001988在癌组织中减少，与肿瘤细胞分化的程度及预后相关。

由于非编码RNA影响多种人类复杂疾病的发生和进展，因此鉴定潜在的ncRNA-疾病关联可以更好地理解ncRNA水平的疾病发病机制，进而利于疾病的诊断和治疗。然而，由于通过实验方法揭示关联是昂贵且耗时的，因此需要用于关联预测的新颖且有效的计算方法。已开发出的方法常见的不足之处有：未能考虑全局相似性，涉及过渡组分的假阳性较高或使用随机化的未验证样本视作阴性导致近似替代的不确切等问题。

发明内容

为了克服现有技术的上述缺陷，本发明提出了一种基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法(Hessian Regularized Non-negative MatrixFactorization method for ncRNA-Disease Association prediction,简称HRNMF-RDA)，通过在NMF框架中引入Hessian正则化捕获数据内在流形结构，增加l2,1范数约束和近似正交约束以保证编码矩阵的群稀疏性，进而提高预测精度。本发明构建的模型属于半监督模型，不依赖于负样本，只需要阳性样本和未标记的样本，大大降低了构建模型的难度，模型有效提升了预测性能，从而获得比较可靠的预测结果。可以通过计算模型预测和排列潜在的ncRNA-疾病关联，将可能性较大的关联对优先用于进一步的实验验证可加速实验验证过程。

为实现上述目的，本发明提供如下技术方案：基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法，具体包括如下六个步骤：

步骤一，分别计算疾病间的高斯相互作用谱核相似性和ncRNA间的高斯相互作用谱核相似性：

若某一个疾病d(i)与非编码RNA之间存在关联，则对应位置记为1否则记为0，形成一个1×nm大小的0或1构成的行向量，记之为疾病d(i)的交互谱IP(d(i))。然后，计算疾病d(i)和d(j)之间的高斯相互作用谱核相似性：

KD(d(i),d(j))＝exp(-γ_d||IP(d(i))-IP(d(j))||²)

上式中，参数γ_d用于控制核带宽，通过归一化新带宽参数γ‘_d获得：

以类似的方式定义非编码RNA m(i)和m(j)之间的高斯相互作用谱核相似性：

KM(m(i),m(j))＝exp(-γ_m||IP(m(i))-IP(m(j))||²)

其中nd表示疾病的数量，nm表示非编码RNA的数量，取γ’_d＝γ’_m＝1；

步骤二，计算Hessian正则化：

通过既有的Hessian矩阵计算方法，输入矩阵KD，输出矩阵B_d，输入矩阵KM，输出矩阵B_m，参数设置为邻近元素的个数K＝25，PCA分解后的左奇异矩阵取其前d列中的d＝6；

步骤三，增加稀疏约束：

将“l_2,1范数正则化”添加到编码矩阵V，强制V中的某些行趋于零，对于新的表示矩阵V，将行稀疏正则化项引入目标函数，将V中的一些行向量收缩为零，保留重要功能并删除不相关的功能；

矩阵V的l_2,1范数定义为：

其中v_j.代表矩阵V的第j行；

步骤四，增加判别约束：

利用给定指示矩阵Y＝{0,1}_N×K，假定其第i个数据点属于第j类别，则Y_ij＝1，缩放的指示矩阵定义为

其中F的每列是：

其中n_j是第j组中的样本数，使用新矩阵表示V获取F中的判别信息，即

其中ε是任意小的常数；

其中，缩放指标矩阵是严格正交的；

其中I_k是k×k的单位矩阵，然后放松正交约束并让V近似正交，即

步骤五，构建目标函数：

将上述约束综合后，可得

式中，(1)

是非负矩阵分解项，刻画了重构误差的大小，其中X是非编码RNA-疾病关联的邻接矩阵，U和V是待求的基矩阵及系数矩阵；(2)tr(VB_mV^T)和tr(UB_dU^T)是Hessian正则化项，可捕获数据内在流形结构；(3)

为矩阵V的l2,1范数，v_j.代表矩阵V的第j行。行稀疏正则化目的是将V中的一些行向量收缩为零，保留重要功能并删除不相关的功能；(4)判别约束项

中I_k是k×k的单位矩阵，目的是让V近似正交进而获得判别信息，其中，λ，μ和γ是正则化参数，可通过网格搜索进行优化。

由拉格朗日乘数法和KKT条件(Karush–Kuhn–Tucker conditions)，可得下面的更新公式：

其中

其中R是一个对角矩阵，其第i个元素为：

类似地，

在一个优选地实施方式中，上述参数取λ₁＝λ₂＝0.7，μ＝0.5，γ₁＝γ₂＝0.3。

步骤六，得到矩阵U、V后，获取评分矩阵F_score＝UV，然后进行排序得到最终结果。

本发明的技术效果和优点：

1、本发明能够较准确的预测ncRNA-疾病关系，将可能性较大的关联对优先用于进一步的实验验证可加速实验验证过程，有效减少生物化学实验大量人力物力的消耗；

2、通过在NMF框架中引入Hessian正则化更好地捕获数据内在流形结构，有效利用阴性样本的信息提高了预测精度；

3、通过增加l_2,1范数约束和近似正交约束保证了编码矩阵的群稀疏性，能减弱噪声数据的影响；

4、本发明构建的模型属于半监督模型，不依赖于负样本，只需要少量阳性样本和未标记的样本，大大降低了构建模型的难度，有效提升了预测性能。

附图说明

图1为本发明总体流程图。

图2为本发明与几种已报道方法在同一数据集上五倍交叉验证的结果图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一种关于miRNA的实施例，而不是全部的实施例(ncRNA还包括其他种类，如lncRNA、circRNA等)。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明的实施例中所使用的已知人类miRNA-疾病关联的数据是从数据库HMDDV2.0(网址http://www.cuilab.cn/hmdd)中检索然后下载的，在对下载的数据进行清洗、分类和标准化之后，可以获得5430个实验验证的人类miRNA-疾病关联，包括383种疾病和495种miRNA。

然后执行如图1所示的基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法，具体包括如下步骤：

步骤一，输入已知的miRNA-疾病关联对，构建邻接矩阵X：

得到一个大小为383×495元素为0或1的矩阵X；

步骤二，分别计算疾病间的高斯相互作用谱核相似性和miRNA间的高斯相互作用谱核相似性：

若某一个疾病d(i)与miRNA之间存在关联，则对应位置记为1否则记为0，形成一个1×495大小的0或1构成的行向量，记之为疾病d(i)的交互谱IP(d(i))。然后，计算疾病d(i)和d(j)之间的高斯相互作用谱核相似性：

KD(d(i),d(j))＝exp(-γ_d||IP(d(i))-IP(d(j))||²)

上式中，参数γ_d用于控制核带宽，通过归一化新带宽参数γ‘_d获得：

以类似的方式定义miRNA m(i)和m(j)之间的高斯相互作用谱核相似性：

KM(m(i),m(j))＝exp(-γ_m||IP(m(i))-IP(m(j))||²)

取γ’_d＝γ’_m＝1。

其中nd表示疾病的数量，此处为383，nm表示miRNA的数量，此处为495，此步计算后得到大小为383×383的对称矩阵KD和大小为495×495的对称矩阵KM，且KD和KM的每个元素都大于0且小于1；

步骤二，计算Hessian正则化：

通过既有的Hessian矩阵计算方法(详见PLoS ONE 9(9):e108474.doi:10.1371/journal.pone.0108474的第7页step1中的k取25step3中的d取6，或者Y.Ma et al.,Methods(2016),http://dx.doi.org/10.1016/j.ymeth.2016.06.017的第2页)，输入矩阵KD，输出矩阵B_d(大小为383×383的对称矩阵)，输入矩阵KM，输出矩阵B_m(大小为495×495的对称矩阵)，参数设置为邻近元素的个数K＝25，PCA分解后的左奇异矩阵取其前d列中的d＝6；

步骤三，增加稀疏约束：

将“l2,1范数正则化”添加到编码矩阵V，强制V中的某些行趋于零，对于新的表示矩阵V，将行稀疏正则化项引入目标函数，将V中的一些行向量收缩为零，保留重要功能并删除不相关的功能；

矩阵V的l_2,1范数定义为：

其中v_j.代表矩阵V的第j行；

步骤四，增加判别约束：

利用给定指示矩阵Y＝{0,1}_N×K，假定其第i个数据点属于第j类别，则Y_ij＝1，缩放的指示矩阵定义为

其中F的每列是：

其中n_j是第j组中的样本数，使用新矩阵表示V获取F中的判别信息，即

其中ε是任意小的常数；

其中，缩放指标矩阵是严格正交的；

其中I_k是k×k的单位矩阵，然后放松正交约束并让V近似正交，即

步骤五，构建目标函数：

将上述约束综合后，可得

式中，(1)

是非负矩阵分解项，刻画了重构误差的大小，其中X是miRNA-疾病关联的邻接矩阵，U和V是待求的基矩阵及系数矩阵；(2)tr(VB_mV^T)和tr(UB_dU^T)是Hessian正则化项，可捕获数据内在流形结构；(3)

为矩阵V的l_2,1范数，v_j.代表矩阵V的第j行。行稀疏正则化目的是将V中的一些行向量收缩为零，保留重要功能并删除不相关的功能；(4)判别约束项

中I_k是k×k的单位矩阵，目的是让V近似正交进而获得判别信息，其中，λ，μ和γ是非负的正则化参数，可以通过网格搜索法进行优化。

由拉格朗日乘数法和KKT条件(Karush–Kuhn–Tucker conditions)，可得下面的更新公式：

其中

其中R是一个对角矩阵，其第i个元素为：

类似地，

在上述的实施例中，经网格搜索优化后，选取参数值λ₁＝λ₂＝0.7，μ＝0.5，γ₁＝γ₂＝0.3；

使用matlab编程实现上述算法时，将矩阵U初始化为383行×100列的随机矩阵，V初始化为100行×495列的随机矩阵(对角矩阵R大小同矩阵V的行数，此处为100行×100列)，设定迭代循环次数为1000次或者满足

退出迭代循环，运行结束后得到矩阵U、V。

计算评分矩阵F_score＝UV，然后进行排序得到最终结果。

本发明的有效性验证：

如图1所示的基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法，采用五重交叉验证进行预测评估，以这种方式进行：将所有已知的miRNA-疾病关联随机平均分成5组，然后将5组中的每一组分别设为测试样本，其他组作为训练样本。

因此，使用训练样本作为本方法的输入得到预测结果，最后将该组中每个测试样品的预测分数与候选miRNA的分数进行比较。

为了减少在获得测试样品的过程中随机划分可能造成的影响，进行了五次交叉验证100次。

结果如图2所示，HRNMFMDA与现有最先进的几种疾病-miRNA关联预测模型之间的性能比较。

本方法在5折交叉验证中取得了0.9451±0.0175的ROC曲线下面积(AUROC)，表现出了比以往所有经典模型更加出色的预测性能。

另外一方面，对于具体的某一种疾病，如淋巴瘤(Lymphoma)，基于上述HMDD V2.0中的已知关联，使用HRNMFMDA做miRNA-淋巴瘤关联预测，所得结果中前50个miRNA有49个能够得到外部数据库的支持。

注：上表的第一列记录了前1-25个相关的miRNA，第三列记录了前26-50个相关的miRNA。表中Ⅰ,Ⅱ,Ⅲ分别代表dbDEMC，miR2Disease和HMDD v3.0三个外部数据库。

最后：以上所述仅为本发明的优选实施例而已，并不用于限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (2)

1.基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法，具体包括如下步骤：

步骤一，输入已知的疾病和非编码RNA关联对并构建邻接矩阵，分别计算疾病间的高斯相互作用谱核相似性和非编码RNA间的高斯相互作用谱核相似性：

若某一个疾病d(i)与非编码RNA之间存在关联，则对应位置记为1否则记为0，形成一个1×nm大小的0或1构成的行向量，记之为疾病d(i)的交互谱IP(d(i))；然后，计算疾病d(i)和d(j)之间的高斯相互作用谱核相似性：

KD(d(i),d(j))＝exp(-γ_d||IP(d(i))-IP(d(j))||²)

上式中，参数γ_d用于控制核带宽，通过归一化新带宽参数γ‘_d获得：

以类似的方式定义非编码RNAm(i)和m(j)之间的高斯相互作用谱核相似性：

KM(m(i),m(j))＝exp(-γ_m||IP(m(i))-IP(m(j))||²)

其中nd表示疾病的数量，nm表示非编码RNA的数量，取γ′_d＝γ′_m＝1；

步骤二，计算Hessian正则化：

通过既有的Hessian矩阵计算方法，输入疾病间的高斯相互作用谱核相似性矩阵KD，输出矩阵B_d，输入非编码RNA间的高斯相互作用谱核相似性矩阵KM，输出矩阵B_m，其中参数设置为邻近元素的个数K＝25，PCA分解后的左奇异矩阵取其前d列中的d＝6；

步骤三，构建目标函数：

目标函数是非负矩阵分解框架整合了Hessian正则化、行稀疏正则化和判别约束项，具体为：

式中，(1)

是非负矩阵分解项，刻画了重构误差的大小，其中X是疾病和非编码RNA构成的邻接矩阵，U和V是待求的基矩阵及系数矩阵；(2)tr(VB_mV^T)和tr(UB_dU^T)是Hessian正则化项，可捕获数据内在流形结构信息；(3)

为矩阵V的l_2,1范数，v_j.代表矩阵V的第j行，行稀疏正则化目的是将V中的一些行向量收缩为零，保留重要功能并删除不相关的功能；(4)判别约束项

中I_k是k×k的单位矩阵，目的是让V近似正交进而获得判别信息；其中‖ ‖_F指矩阵的F范数，tr()表示矩阵的迹，λ，μ和γ是非负的正则化参数，可通过网格搜索进行优化选取，参数组合是λ₁＝λ₂＝0.7，μ＝0.5，γ₁＝γ₂＝0.3；

步骤四，使用迭代算法计算矩阵U和V；

步骤五，得到矩阵U、V后，计算评分矩阵F_score＝UV，然后进行排序得到最终预测结果。

2.根据权利要求1所述的基于Hessian正则非负矩阵分解的非编码RNA与疾病关系预测方法，其特征在于，上述步骤四中求解U和V的迭代公式具体为：

其中

其中

R是一个对角矩阵，其第i个元素为：

Images