跪求高手帮助解决perl问题

b114213903

回复 30# chenhao392


    因为我看到他“orthol-msu.txt”文件里的ID都没“.”，就直接用了/^（\w+）/，呵呵

chenhao392

用于多个isoform的code
那个ID的问题 fix了

1. #!/usr/bin/perl
   
2. 
   
3. use strict;
   
4. use warnings;
   
5. 
   
6. my %id;
   
7. my %seq;
   
8. 
   
9. #load orthoMCL result from MSU
   
10. open FILE, "<$ARGV[0]" or die "$!\n";
    
11. while(<FILE>){
    
12.     chomp;
    
13.     #skip header line
    
14.     if($. == 1){
    
15.         next;
    
16.     }
    
17. 
    
18.     #get ortholog groups
    
19.     my ($orthoID,undef,$speciesCount,@line)=split(/\s+/,$_);
    
20.     for(my $i=$speciesCount;$i<scalar(@line);$i++){
    
21.         $id{$orthoID}{$line[$i]}="";
    
22.         $seq{$line[$i]}="";
    
23.     }
    
24. 
    
25. }
    
26. close FILE;
    
27. 
    
28. 
    
29. #scanning through the seq
    
30. my $switch=0;
    
31. my $id="";
    
32. my $isoformID;
    
33. open FILE, "<$ARGV[1]" or die "$!\n";
    
34. while(<FILE>){
    
35.     chomp;
    
36.     chop;
    
37.     if($_ =~ />([\w\.]+)/){
    
38.         $switch =0;
    
39.         $id=$1;
    
40.         $isoformID=$id;
    
41.         $id=~s/(.{4,})[\._].*/$1/;
    
42.         if(defined $seq{$id}){
    
43.             $switch = 1;
    
44.             $seq{$id}.=">$isoformID\n";
    
45.         }
    
46.     }
    
47.     elsif($_ =~ /^[ATCGatcg]/ && $switch == 1){
    
48.         $seq{$id}.="$_\n";
    
49.     }
    
50. 
    
51. }
    
52. close FILE;
    
53. 
    
54. foreach my $orthoID(keys %id){
    
55.     open OUT, ">$orthoID.fa" or die "$!\n";
    
56.     foreach my $id(keys %{$id{$orthoID}}){
    
57.         print OUT "$seq{$id}";
    
58.     }
    
59.     close OUT;
    
60. }

复制代码

chenhao392

回复 31# b114213903


    囧..
    然后，你怎么知道Orthologs-msu.txt 可以先按 \t 切，再按 \s+ 切？
    我这里的测试数据，貌似找不到\t ， 我按\s+ 切了.