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

rubyish

ref

NM_000721        chr1        100        110        +        246
NM_020394        chr1        110        300        -        2470

lnc

lnc-FAM20C-5:2        chr1        100        150        +

输出文件1
OR
输出文件2

yanglc2013

回复 43# rubyish

lnc和ref比较，ref作为数据库。输出四个文件，输出的内容是lnc中的序列：

文件一：contain：就是ref图中的红色区间包含lnc图中红色区间。
文件二：overlap：就是lnc图中的红色区间和ref图中红色区间相交。
文件三：gapped：就是lnc图中的红色区间和ref图中红色区间包含，也不相交，是独立分开的。
文件四：reverse：就是lnc图中的红色区间包含ref图中的红色区间，和文件一刚好相反。

   

yanglc2013

回复 42# picbhan

大哥，其实是这样的，lnc中每个区间都要单独和ref中的区间比较，每一个染色体（比如lnc-HES5-1:2        chr3        14 100        -），它和NM_016505        chr3        50 128        +        827是overlap的，但是它和NR_026818        chr3        110 150        -        585却是non_include，所以输出的non_include文件中应该也包含lnc-HES5-1:2        chr3        14 100        - 啊。（就是lnc中同一个区间，他可能有很多种情况，都要分别输出。这也没有考虑正负链啊）

yanglc2013

回复 42# picbhan
我意思是对每种情况都输出啊。就是一个region，它可能是contain，可能是gapped，可能是overlap也可能lnc包含ref，如果这四种情况都存在，那么输出的四个文件中，都要包含这个region。大哥，真抱歉啊，我没表达清楚。你那个程序，我还没研究透，不是很会改你的程序。还是会出错。


   

rubyish

比如
lnc-HES5-1:2        chr3        14 100        -

NM_016505        chr3        50 128        +        827
NR_026818        chr3        110 150        -        585
NR_026835        chr3        12 150        -        586
NR_026822        chr3        16 50        -        534
NR_0268200       chr3        6 12        -        531

输出文件?

yanglc2013

回复 47# rubyish
不是啊，输出的是lnc中的数据。（比如lnc-HES5-1:2        chr3        14 100        -）中的 [14,100]这个区间和ref中的（NM_016505        chr3        50 128        +        827
）中的[50,128]这个区间比较，它们有一个相交的区间[50,100]，所以在overlap文件中输出。如果和（NR_026818        chr3        110 150        -        585）比较，即区间[110,150]，没有公共部分，就输出到文件gapped。如果和（NR_026835        chr3        12 150        -        586）比较，lnc包含在ref区间中，输出到contain中，如果反过来，ref包含在lnc中，就输出到reserve文件中，输出的文件都是（lnc-HES5-1:2        chr3        14 100        -），如果四种情况都存在，那么，输出的四个文件中，都应该有（lnc-HES5-1:2        chr3        14 100        -）
   

rubyish

{:3_188:}

1. 
   
2. Delete
   
3. 有错误
   

复制代码

rubyish

{:3_188:}

1. 
   
2. Delete
   
3. 有错误
   

复制代码

picbhan

回复 45# yanglc2013


    你觉得这样子输出有具体意义吗？inc中任何一个region都会被输出到non_include中去，只要你的refGene足够多。对于后面三种有交集的情况，如果要输出到多个文件中，你把下面的代码替换一下就行。

1. my %result; # one of ref_include_inc, inc_include_ref, overlap or undef
   
2.             while (@ref_data) {
   
3.                 my @rd = @{+shift @ref_data};
   
4.                 # igore genes on the left side of inc region
   
5.                 next if $rd[3] < $id[2];
   
6.                 push @store, [@rd];
   
7.                 # stop when this gene is on the right side of inc region
   
8.                 last if $rd[2] > $id[3];
   
9.                 # this gene must overlap with this inc region
   
10.                 if ($id[2] >= $rd[2] && $id[3] <= $rd[3]) {
    
11.                     $result{ref_include_inc} = 1;
    
12.                 }
    
13.                 elsif ($rd[2] >= $id[2] && $rd[3] <= $id[3]) {
    
14.                     $result{inc_include_ref} = 1;
    
15.                 }
    
16.                 else {
    
17.                     $result{overlap} = 1;
    
18.                 }
    
19.             }
    
20.             # output result
    
21.             if (%result) {
    
22.                 if ($result{ref_include_inc}) {
    
23.                     print {$ref_include_inc} join("\t", @id), "\n";
    
24.                 }
    
25.                 if ($result{inc_include_ref}) {
    
26.                     print {$inc_include_ref} join("\t", @id), "\n";
    
27.                 }
    
28.                 if ($result{overlap}) {
    
29.                     print {$overlap} join("\t", @id), "\n";
    
30.                 }
    
31.             }
    
32.             else { # if $result is undef, means no overlapped genes
    
33.                 print {$non_include} join("\t", @id), "\n";
    
34.             }
    
35.             
    
36.             # put this genes back into @ref_data, use to check next inc region
    
37.             unshift @ref_data, @store;
    
38.             @store = ();
    
39.         }

复制代码

yanglc2013

回复 51# picbhan

不会出现这样情况的。因为ref中每个染色体，比如很多chr1都是来源于有一个DNA，这样发生的剪切，如果存在contain，那么另一种情况只能是gapped。