20积分悬赏，求个perl程序，新手，写了好几个都无果。

yanglc2013

回复 9# blackfur

意思就是ref中的红色区间包含lnc那个红色区间啊。其实就是两个相同的染色体随机切割后。有的会没有交集，有的有，有的一条包含在另一条内。

   

yanglc2013

回复 9# blackfur
楼下正解啊


   

picbhan

昨晚写的那个程序忘了三行，结果漏掉了一点没有交集的行。修改如下：

1. #!/usr/bin/perl
   
2. use strict;
   
3. use warnings;
   
4. 
   
5. my (%ref, %data);
   
6. 
   
7. open my $in , '<', 'refGene.txt'
   
8.     or die "Can't open 'refGene.txt' for reading. $!\n";
   
9. while (<$in>) {
   
10.     chomp;
    
11.     my @t = split;
    
12.     push @{$ref{$t[1]}}, [@t];
    
13. }
    
14. close $in;
    
15. 
    
16. open $in , '<', 'lncipedia_location.txt'
    
17.     or die "Can't open 'lncipedia_location.txt' for reading. $!\n";
    
18. while (<$in>) {
    
19.     chomp;
    
20.     my @t = split;
    
21.     push @{$data{$t[1]}}, [@t];
    
22. }
    
23. close $in;
    
24. 
    
25. open my $include , '>', 'include.txt'
    
26.     or die "Can't create 'include.txt' for writing. $!\n";
    
27. open my $non_include , '>', 'non_include.txt'
    
28.     or die "Can't create 'non_include.txt' for writing. $!\n";
    
29. open my $overlap , '>', 'overlap.txt'
    
30.     or die "Can't create 'overlap.txt' for writing. $!\n";
    
31. 
    
32. 
    
33. while (my ($chr, $v) = each %ref) {
    
34.     my @ref_data = sort {$a->[2] <=> $b->[2]} @$v;
    
35.     next unless exists $data{$chr};
    
36.     my @inc_data = sort {$a->[2] <=> $b->[2]} @{$data{$chr}};
    
37.     for (@ref_data) {
    
38.         my @rd = @$_;
    
39.         last unless @inc_data;
    
40.         while (1) {
    
41.             last unless @inc_data;
    
42.             my @id = @{+shift @inc_data};
    
43.             if ($id[3] < $rd[2]) {
    
44.                 print {$non_include} join("\t", @id), "\n";
    
45.                 next;
    
46.             }
    
47.             if ($id[2] > $rd[3]) {
    
48.                 unshift @inc_data, [@id];
    
49.                 last;
    
50.             }
    
51.             if ($id[2] >= $rd[2] && $id[3] <= $rd[3]) {
    
52.                 print {$include} join("\t", @id), "\n";
    
53.             }
    
54.             else {
    
55.                 print {$overlap} join("\t", @id), "\n";
    
56.             }
    
57.         }
    
58.     }
    
59.     for (@inc_data) {
    
60.         print {$non_include} join("\t", @$_), "\n";
    
61.     }
    
62.     delete $data{$chr};
    
63. }
    
64. 
    
65. while (my ($chr, $v) = each %data) {
    
66.     my @inc_data = @$v;
    
67.     for (@inc_data) {
    
68.         print {$non_include} join("\t", @$_), "\n";
    
69.     }
    
70. }
    
71. 
    
72. close $include;
    
73. close $non_include;
    
74. close $overlap;

复制代码

picbhan

回复 8# dahe_1984


    为毛我跑的时候貌似一秒钟就完了，数据不大啊。

yanglc2013

回复 4# picbhan
这个程序好像无论ref还是lnc中重复的chr没有比对啊。那不能忽略，其起始和终止位置不一样，还是得单独比对。输出的文件中，可能三个文件中有相同的name，但是其区间不一样。


   

picbhan

回复 15# yanglc2013


    我不明白你说的什么意思。我是按照你在帖子开头说的意思来写代码的，即只在inc和ref之间比较。所以你如果觉得不是这个意思，麻烦说清楚你到底要什么结果。从生物的角度说就行，你这每个区间是个基因？

picbhan

回复 15# yanglc2013


    你说的name是指第一列吗？同一个文件里面相同染色体上的区间也要比较？

climby

如果让我迅速的解决这个问题，鉴于数据的工整性，我会将数据倒入数据库，利用数据库的查询功能，会非常迅速的得到想要的结果。
如果非要用Perl 来解决这个问题，我看了一下数据， ref文件有4万多条记录， lnc文件有3万多条记录，一次读进内存然后比较的话，对机器的要求还是比较高的。
我的处理方法是
1. 先按照第二列对文件进行分组处理。
   读取ref 文件，按照第二列进行分组，生成一些小文件，放在ref目录下。
   分组后的文件如下：
   对lnc文件进行类似处理
   

2.然后对分组后的文件一对一进行对比分析。

这样处理起来，应该对机器的压力会小很多。由于我上午抽不出时间来完成所有代码，仅列出方法供大家参考

rubyish

是否:

1. 
   
2. Delete
   
3. 有错误
   

复制代码

yanglc2013

回复 18# climby
这个方法我也试过，但是比较麻烦，因为我后面还有很多这样的数据。如果每个都得细分，那很耗时间的。