一个Nphard问题，我都快把头想破了也没有想出来。(200分)

  • 主题发起人 主题发起人 qiaohj
  • 开始时间 开始时间
下载数据的网址
http://www.ncbi.nlm.nih.gov/Traces/trace.cgi?
在search框中添进去TRACE_TYPE_CODE='SHOTGUN'就可以得到我们现在要的实验数据，每一段就是一个片断，
returned 27217811 Items，使所有的片断数目。
 
修改一下查询语句
SPECIES_CODE='MUS MUSCULUS' and TRACE_TYPE_CODE='SHOTGUN'
 
先进入 http://www.tigr.org/，浏览一下内容，进入Genome Database下的Bacillus anthracis，再到NCBI Trace Archive 中下载shotgun raw data。
 
>gnl|ti|7378889 TUWAA1U0269
NTTTTAGAATTGGNTAGGTCGGNGAAGNCTTAAGGGNNANNGNAGCNNGNNNGAANGTNNTNGGGGGGTG
晕——怎么还有"N"？？
还有，我发现实际数据的基因分布很不符合随机分布的特性——看来我的索引匹配要更加
严格了。
 
最新进展：
产生一个完全随机的、长度为3000K的DNA序列。六条这样的序列共断裂成长度在1..1200
之间的共30668个片断，随后丢失信息（随机获得前500..1000个序列信息），共丢失序列信息
1441437个，占总信息量的约8%。对这些总长度为16990563的片断进行索引化处理——每9个
相邻序列为一组（共4^9即262144组索引）——耗时30.09秒（PIV 1.6A），内存204M。
正在写核心拼接算法，难啊...
还有，我发现网上的数据几乎没看到长度有小于300的，Why?是不是长度太小的片断无法被
检测出来？您总要对数据的特点加以说明吧。
 
昨天有些事情，没有来得及跟贴，序列中的N代表任意碱基都可以
>>您总要对数据的特点加以说明吧。
我现在正在和有关方面联系着方面的东西，找华大的人，要他们说明一下数据的特点。现在大家办事都比较拖沓，郁闷 的说
 
想找个时间跟您当面讨教，您在北京吗？
 
可惜呀，我在上海。（要不然就可以搓一顿了 :) ）实际上，像这种纯技术问题，我贴个
算法就可以了，没有什么东西需要“讨教”的。不过，如果你打长途不用自己掏钱的话，
可以考虑电话联系（19:00-0:00），你可以把Mail给我，我把号码给你。
 
qiaohj@panda.ioz.ac.cn
没准什么时候就去上海了，近期有可能去的。
 
最后再呆1天，然后就结题了
 
不用这么早吧... 请您再给我两周时间。
构造索引消耗时间：8.88 秒
[索引信息]
索引数=4194304
非零长度索引数=980703
非零比率=23.38%
索引长度和=5556671
平均非零索引长度=5.666
索引最大长度=30
[片断索引信息]
索引数=10060
总长度=5556304
平均长度=552.317
最大长度=990
片断索引匹配消耗时间：5.53 秒
[片断匹配信息]
片断数=10121
索引匹配数下限=200
拥有匹配的片断数=8327
拥有匹配片断比率=82.27%
总合格匹配数=49138
平均每片匹配数=5.901
最大匹配数=13
拼接片断消耗时间：10.92 秒 （——注：正在处理这个部分）

算法思路：
1.对片断进行索引化处理（我发现索引越大，后期处理就越快——现在采用的是4^11）
2.根据索引信息，找到重合索引数量达到指定值（比如200）的片断对
3.从这些拥有匹配对象的片断中任取一个，找到和它相匹配的片断们的相对偏移量，进而
获得所有匹配片段的相对覆盖范围
4.根据覆盖范围，选择位于最前和最后的片断，接着匹配，直到找不到匹配片断为止（我
估计此时已经超不多完成整个DNA序列中30%-70%的部分了）
5.将已经成功匹配的（也就是位于拼合后的DNA之内的）片断排除，回到过程3，直到片断
耗尽。
6.对获得的大片DNA进行再次拼合，构成完整的DNA。
 
由于无法获得完整的测试用实际数据，我目前只能用随机数模拟基因的碱基序列——这样
做的必然结果是：基因序列的分布是均匀的、完全随机的。我发现网上的基因数据分部非常
的不均匀，虽然总的说来，四种碱基的比例接近1:1:1:1，但是从局部看来，就完全不是这样
了，因此我需要基因中碱基对分布的不均匀参数。还有，您提到了“模糊匹配”，我现在写
的算法可以满足非精确匹配的要求，但我需要知道两条原始基因的一致比率。以上这两点很
重要，因为它们直接涉及到索引的可靠性。
我目前的初步拼接算法，已经可以将8000多个碎片中的99%以上拼接成不到160条大片断。
由于索引的精度关系，如果将拼接精度从6提高到10，最后形成的片断数就会急剧增长到400
多。
 
>>因此我需要基因中碱基对分布的不均匀参数
由于谁也没有得到实际情况的序列，这些参数是得不到的。
但是不均匀参数您暂时可以不考虑在内。如您所说在DNA中有很多重复序列，
有一些我们在生物学术语上面叫做“保守序列”，这些保守序列将是我们下一
步讨论的问题，即如何通过现有的DNA片断找到这些保守序列，如果能够找到这
些保守序列，我们就可以将这些保守序列进行替换，然后将DNA片断中所有的保
守序列用一定的符号替换出来，这样得到的新的片断就是现在您做实验数据用
的完全随即分布的DNA片断。
 
所以，我的想法是：假设我们现在已经能够找到序列中的所有的保守序列，而
且已经将保守序列进行了替换。这样我们得到的序列将是一个完全随机的无序
的片断组合，通过这样的数据我们将能够更好的，更顺利地进行拼接过程，在
拼接结束以后，再将原来的保守序列替换回来，这样应该能够恢复DNA的原来的
面目。不知道我的表述是不是清楚。
 
>>我需要知道两条原始基因的一致比率
这个我需要您进一步描述一下原始基因的一致性比率是个什么概念，是说有多
少单位长度的一致就表示这两条DNA片断是接在一块儿的吗？这个没有一个明
确的概念，也就是说，在不同的长度的限制条件下，也许会得到不同的拼接结
果，对程序的要求就是能够对不同的标准作比较，得出最好的方案来：换句话
说就是多大长度的一致能够的最好的拼接结果。通俗一点就是如果定为一个碱
基匹配就认为是接在一块的，这样能得出无数种情况，而制定的长度标准过长，
就有可能没有结果，就是要找出一个最好的匹配结果。
 
>>只能用随机数模拟基因的碱基序列——这样做的必然结果是：基因序列的分布是均匀的、完全随机的。
这和我的想法完全符合
 
关于“一致性比率”：
我是针对您在前面提出的“并不要求精确匹配，字符串模糊匹配就行了，保证匹配率在95%
以上就行”——不知道您是不是这个意思：任取两条原始DNA，将它们从头至尾进行逐位比较，
必然有95%以上的位是相同的（如果是的话我的任务就有可能完成了），如果不是这样，也就
是说两个DNA中可能出现“错位”——那就麻烦了。
>>在不同的长度的限制条件下，也许会得到不同的拼接结果
我现在的拼接方法也是需要用户指定拼接精度的。原理如下：
1.对数据进行索引化处理
所谓“索引”，是指长度固定的连续基因序列。我在程序中采用的索引长度为11，而每个
元素有4种可能，总索引数量为4^11=2^22=4*10^6。考虑到单个DNA的长度约在2*10^6-3*10^6
左右，如果保证随机性的话，索引可以保证99.99%以上的唯一性。
索引化包含一个索引对应于哪些片断以及一个片段对应于哪些索引。通过实验我发现对于
有片断对应的索引，它们的平均对应片段数为5.7左右（注：原始DNA有6条），而平均每个片
段则对应了约635个索引，接近它们的平均长度。
2.根据索引化的结果进行片段匹配
在这一步，用户可以设定“索引匹配数下限”参数，默认值为100——这意味着，只有当两
个片段对应的索引集合中，有100个以上一致的时候，程序才认为两者可能有重合部分，考虑
到前面已经论证过的索引的唯一性，想从两个不相关的片断中找到5个以上重合的索引都几乎
是完全不可能的，因此，采用100作为匹配下限是绝对安全的。通过实验，在原始数据丢失率
为8%左右时，我可以实现高达99.98%的成功匹配率。
在这一步完成之后，我们可以得到每个片段对应于哪些和它相匹配的片断的列表。
3.初步片断拼接
在这一步中，我没有采用根据位于片断头部或尾部的一段较长序列信息进行定位的方法，
因为这要求几条原始DNA必须完全一致（99.999%以上一致也可以）。为了减小计算量，我先
针对被匹配的片断生成稀疏索引（这样生成的索引数量为第一步的1/5，并且可以保证采样点
在片断中的绝对均匀分布）；然后再计算每个匹配片断的索引和它的交集（我通过跟踪发现，
平均说来交集运算的结果索引有约10-30个）；接着，程序对这些结果索引进行片断定位，得
到每个索引在片断中的偏移量，这样形成两个数组；最后，程序对这些偏移量信息进行对比
分析，发现如果有足够数量的索引的相对位置一致，就认为片段匹配成功，并返回两个片断
的相对偏移量。用户可以控制最后一步的索引对小匹配对数，默认值为6，也就是说两个片段
中，至少有六组相对位置一样的索引才能作为匹配偏移量的可信证据。
说了大半天，不知说明白了没有...
 
>>任取两条原始DNA，将它们从头至尾进行逐位比较，必然有95%以上的位是相同的（如果是的话我的任务就有可能完成了）
(我说的后面这句话是个普遍规律，和这个算法无关，可以不看)原始DNA理论上讲是完全相同的，但是在测序过程中可能会出现误差。
我说的95%是指拼接算法结束以后，得到的拼接出来的字符串，和拼接以前的原始字符串相比，允许有95%的不同。
 
我找个明白人，我现在的浅薄的算法基础已经跟不上您的思路了
 
哈哈！乔兄，你跑不掉了！——我专程到北京吃你来了！ :)
 
您必须登录或注册才可回复。

Similar threads

不
B编码以及BT种子文件分析 - 不得闲
回复
0
查看
1K
不得闲
不
S
[原]详解Variant 的相关函数
回复
0
查看
1K
SUNSTONE的Delphi笔记
S
S
[原]详解Variant 的相关函数
回复
0
查看
900
SUNSTONE的Delphi笔记
S
D
[原]2015-07-06：发布OrangeUI For FMX 1.50
回复
0
查看
909
DelphiTeacher的专栏
D
D
[原]2015-07-06：发布OrangeUI For FMX 1.50
回复
0
查看
704
DelphiTeacher的专栏
D
后退
顶部