取特征序列代码改进

第一部分：整理逻辑细化判断条件。

整理了取特征序列的代码逻辑，主体逻辑是：同时看两条序列，每次迭代向下一位。同时看本条和下一条序列。如果两条同名且同时read1或者read2，使用同一种判断方式；在read不同或者序列易名时，前者一定是read1，后者一定是read2，且作为read1或者read2交汇处、不同序列交汇处只有一条且未最后一条，一定是含有H的，不能是S。总的来说相当于分块考虑，分为交汇处和非交汇处。见图如下：

总逻辑

总结：分块考虑交汇处（r2 和r2交汇、两个不同的测序数据的交汇）<特殊>和只有read1或者read2的部分<一般>

            方法是考虑特殊与一般，特点是使用两条序列进行迭代，但是只考虑第一条，第二条是提供位置信息的（看在块内还是在交汇处）

整理完逻辑修改了代码的两个部分，细化了在不同部分的判断条件：从上可知，主要的判断部位是交汇处和非交汇处，交汇处应只有H，非交汇处应是((S or H) and M)

以下第二部分：查看文件寻找问题

接着查看了r1 r2比对条数都是2时候的情况，可见蓝圈以外的是对的，蓝圈内是有问题的，从代码中能找到问题：

问题 问题代码

我是要求有r2时才做存储，如果只有r1，他会继续等待直到符合条件的r2出现才存储。是个重大错误。可以预测的是：这导致许多本应是one_pair的序列被放到two_pair，损失大量单条多比对的序列。这在最后的文件大小处也得到验证。

改进

改进后查看one_pair.sam 发现仍有错，问题在于最后没把不取出的清零。

one_pair问题

清零之后抽样调查结果暂时无误

清零 大小

大小而言，单条的文件大小都增大，两条的文件大小在减小，也是比较合理的情况。

至此，对取特征序列的代码改进暂时结束。还有再细化的空间，需待再多观察取出的序列。

还有一个小的优化，当序列next到最后一条时，会因为没有下一条而导致提前退出，少判断了文件的最后一条。就逻辑而言修补它比较麻烦，简单的办法是对j.next（）  try except ，except里给j赋值，因为j的值不用来取只是用来判断的，所以只要使得其qname与i.qnme不相等即可。可以保证到i的最后一条也进入判断。之后对i同理try except exit（）。程序就完整了。当然现在不优化也能使用，只是少一条。

序列next到最后一条

改进结束，但今天的工作还没有结束，接下来需要筛查每个文件，并转格式放igv观察。