DNA图谱分析的一种SPDP问题求解模型及其算法

本发明涉及一种模型求解算法领域，主要涉及dna图谱分析的一种spdp问题求解模型及其算法。

背景技术：

1、在绘制dna限制性图谱时，知道dna分子长度十分必要。但目前存在的实验技术无法直接测量dna分子长度，只能通过其他方法：从dna片段入手，切开dna分子，对其分别进行局部测量，即测量dna的片段，然后通过其他手段间接地得出dna分子长度。但由于dna的特殊性，切开的dna片段无法辨别其原本在dna中的位置，dna分子的排列顺序会丢失，这也就引入了一个关键问题：dna的测序。想要知道dna分子的排序，构造一张限制性图谱，目前科学家仅能通过间接的图谱信息来实现。

2、现有一种方法：简化的部分消化(spdp)法，其思路是采用限制性酶消化dna分子，获得所需某些点位的片段，根据片段长度来重构。具体操作就是，首先将测序的dna分子复制成n+1份，其中一份采取所有限制性位点切开，获得n+1个细分的dna片段数据组，称n+1数据组；剩下的n份则每一份在不同的一个限制性位点切开，切成两段，获得2n个成对的dna片段数据组，称2n数据组。dna分子长度是一定的，dna片段的顺序不知道，现在的问题就是如何利用两组数据重构出原来的dna分子顺序，称spdp问题。如何解决spdp问题，关键在于建立数学模型，并采用合适的算法进行计算求解。

技术实现思路

1、针对dna分子重构的spdp问题，本发明提供一种求解模型及其算法，根据spdp实验方法所得的n+1数据组和2n数据组，建立数学模型，并提出了对应的两种算法。

2、建立数学模型如下：

3、模型依据的数据来自spdp实验方法，现有一个dna分子，采用spdp实验方法进行限制性酶消化，获得所需片段数据：n+1数据组和2n数据组，如图1所示。

4、作为优选，已知dna分子有n个限制性位点，以节点1，…，n，表示，令dna分子两个端点分别为0和n+1，其中0表示起点，n+1表示终点；设di,j是第i个点到第j个点之间的距离，di,j恒为正，此时i、j＝0,1,…,n,n+1。解释两个实验数据组：n+1数据组中，每一个dna片段长度都是相邻两个点之间的距离，即di,i+1；2n数据组中，每一个dna片段长度和另外某一个dna片段长度成对，是某个限制性位点到0和n+1的距离，即di,0和di,n+1，如图2所示。

5、作为优选，n+1数据组可以表示为集合an+1＝{di,i+1|i≤＝n且i∈n}，2n数据组可以表示为集合b2n＝{(di,0,di,n+1)|i＝1,2,,n}。

6、作为优选，从n+1数据组入手，重新排列an+1中的元素，考虑到b2n中元素的递增性，排列时的第一个元素直接选择b2n的最小数，假设排列之后的集合为an+1'＝{a1,a2,,an+1}。

7、作为优选，按spdp实验方法，消化集合an+1'重构所得的dna分子，得到每个限制性位点到两个端点的距离，公式如下：

8、

9、其中，(di,0',di,n+1')是计算出的成对数据，i＝1,2,…,n；当i＝0或i＝n+1时，表示两个端点到自身的距离，故初始化

10、作为优选，根据公式(1)计算出成对的dna片段数据(di,0',di,n+1')，将其和集合b2n中的任一元素(dj,0,dj,n+1)作比较，若有相同的就表明第i个元素的排序可能正确，继续进行下一个限制性位点的判断；反之，则说明排列不符合要求。

11、作为优选，对于集合an+1'＝{a1,a2,…,an+1}来说，因为b2n元素的唯一性，只要有一个元素ai对应计算出的(di,0',di,n+1')不满足集合b2n的内容，那么集合{ai|i＝1,2,…,n}就可以排除。建模思路如图3所示。

12、作为优选，最终，比较筛选出所有满足条件的排列集合an+1'，按照集合中元素的排列顺序，对应排列n+1数据组就是求解结果。

13、作为优选，在数学模型的基础上，提出了两种算法来进行计算求解。

14、算法1——dna序列穷举法：

15、作为优选，全排列集合an+1＝{di,i+1|i≤＝n且i∈n}中的所有元素，第一个元素确定为b2n＝{(di,0,di,n+1)|i＝1,2,…,n}中最小数，即minb2n。设ai,j表示第i次排列的第j个数据，全排列结果可以表示为如下矩阵：

16、

17、矩阵a中，每一行表示为一次排列结果，总共有n！。

18、作为优选，an+1＝{di,i+1|i≤＝n且i∈n}中可能会出现重复元素，这使得全排列出的结果也会出现重复现象。筛除重复的全排列行，假设筛除了k行。

19、作为优选，针对每一次an+1元素排列，实际上就是一次dna分子重构，根据公式(1)计算，得到：b2n'＝{(di,0',di,n+1')|i＝1,2,…,n}。然后，拆开成对数据，得到集合：b＝{d1,0′,d1,n+1′,d2,0′，d2,n+1′,…，di,0′,di,n+1′,…,dn,0′,dn,n+1′}，再进行从小到大排序，又得到新一个数据集。

20、作为优选，设第一次排列(矩阵第一行)并计算，得到从小到大的数据集为：{l1,1,…,l1,j,…,l1,2n}；以此类推，第i次排列得到集合{li,1,li,2,…,li,j,…,li,2n}，li,j指的是第i次排列之后，计算并排序得到的第j个数据；最终，得到每次排列之后，每个限制性位点到两个端点的距离，可以表示为如下矩阵：

21、

22、其中，li,j≤li,j+1；i＝1,2,…,n！-k，j＝1,2,…,2n；每一行对应一次排列。

23、作为优选，为了便于比较，将2n数据组，即b2n＝{(di,0,di,n+1)|i＝1,2,,n}也进行转换，得到数据集：l2n＝{l1,l2,…,li,…,l2n}，li表示排序之后的第i个数据并且li≤li+1，i＝1,2,…,2n-1。

24、作为优选，对l2n和矩阵l中每一行元素做差并求绝对值，计算公式如下：

25、δli,j＝|lj-li,j|  (4)

26、式中，δli,j表示的第i次排列计算并排序所得的第j个片段长度数据与2n数据组的绝对差是，最终得到绝对差值矩阵如下：

27、

28、作为优选，判断矩阵δl中每一行，由唯一性，若某一行求和为零，则该行对应的排列就是原dna分子序列；反之则不是。直到判断完所有行，满足条件的即为spdp实验方法所求结果，可能存在多解。

29、算法2——dna序列树搜索法：

30、作为优选，根据集合b2n＝{(di,0,di,n+1)|i＝1,2,,n}递增性和唯一性，知道minb2n就是dna分子序列的第一个片段对应的长度数据，从而可以确定an+1＝{di,i+1|i≤＝n且i∈n}中d0,1元素是哪个。

31、作为优选，按照集合an+1＝{di,i+1|i≤＝n且i∈n}，改进公式(1)得到：

32、

33、作为优选，知道d0,1＝minb2n，由公式(6)计算可得(d1,0',d1,n+1')，因为d0,1位置是确定的，此时直接将an+1中的d0,1和b2n中与(d1,0',d1,n+1')一致的数据剔除。

34、作为优选，按照顺序，接下来确定d1,2。搜索an+1＝{di,i+1|i≤＝n且i∈n}中剩余元素，依次代入到d1,2位置，按照公式(6)计算(d2,0',d2,n+1')，然后判断该成对的dna片段数据是否在b2n＝{(di,0,di,n+1)|i＝1,2,,n}里面。如果(d2,0',d2,n+1')在b2n当中，则说明该元素可能是d1,2；反之，则该元素不是d1,2，判断下一个元素，直到所有剩余元素搜索结束。

35、作为优选，因为判断出的d1,2可能有多个，那么在搜索d2,3时就需要分别进行讨论。搜索d1,2所有可能，当d1,2代入一个元素时，此时便将该d1,2对应的数据分别从an+1和b2n中剔除，然后判断d2,3；判断完再回溯到d1,2继续搜索，代入下一个元素继续找d2,3的其他可能元素。直到搜索完d1,2，搜索并判断d2,3结束。

36、作为优选，以此类推，最终找出符合条件的所有{d0,1,d1,2,…,di,i+1,…dn-1,n}。

37、作为优选，集合{d0,1,d1,2,…,di,i+1,…dn-1,n}中元素排序，就是spdp实验方法下原dna分子的序列，结果可能存在多解。

38、作为优选，算法的实现实际上是树结构的体现，如图4所示。

39、作为优选，两种算法都具备可靠性，准确度一致，算法2的效率更高；两种算法结果是顺序形式的dna分子序列，只需将其逆序再呈现出来即可；结果存在多解，是由于spdp实验方法本身具有一定的缺陷。