如何计算Perl中两个不同结构的文件中的匹配数
如何计算Perl中两个不同结构的文件中的匹配数
提问于 2013-09-04 16:04:59
如果问题太简单，很抱歉。但是，当涉及到下面代码中的一个问题时，我的头脑就陷入了困境。我们会跟进这个问题。
1. open(my $go_file, "<", "gene_associations_go_human.txt") or die "Can't open the file!";
2. open(my $selected_genes, "<", "my_selected_genes.txt") or die "Can't open the file!";
3. open(my $output, ">", "output_go_file.txt") or die "Can't open the file!";
4. my %go_hash;
5. chomp(my @sel_genes=<$selected_genes>);
6. while(<$go_file>){
7. chomp($_);
8. my @go_line=split("\t", $_);
9. $go_hash{$go_line[4]}=[] unless exists $go_hash{$go_line[4]};
10. push @{$go_hash{$go_line[4]}}, $go_line[2];
11. }
12. foreach my $go_term (sort keys %go_hash){
13. my @genes=@{$go_hash{$go_term}};
14. @genes= uniq(@genes);
15. my $count=0;
16. foreach my $element(@genes){
17. my $score=grep{$element eq $_} @sel_genes;
18. $count++ if($score>0);
19 }
21. @genes=sort(@genes);
22. push(@genes, ($#genes+1, $count));
23. print $output($go_term."\t".join("\t",@genes)), "\n";
24. }
25. close($go_file);
26. close($selected_genes);
27. close($output);编辑:输入和输出文件示例
**$go_file:**
UniProtKB A0A183 LCE6A NA GO:0031424
UniProtKB A0A5B9 TRBC2 NA GO:0016021
UniProtKB A0AUZ9 KANSL1L NA GO:0000123
UniProtKB A0AV02 SLC12A8 NA GO:0006813
UniProtKB A0AV02 SLC12A8 NA GO:0015293
UniProtKB A0AV02 SLC12A8 NA GO:0016021
**$selected_genes:**
DOLPP1
SPIC1
KANSL1L
SLC12A8
TRAF1
CDF7
**$output should be like this:**
GO:0000123 KANSL1L 1 1
GO:0006813 SLC12A8 1 1
GO:0031424 LCE6A 1 0
GO:0015293 SLC12A8 1 1
GO:0016021 SLC12A8 TRBC2 2 1我正在制作%go_hash (基于$go_file)，它基于关联的$go_terms (文件的5th column)来保持基因数组(文件的3rd column)，因此%go_hash中的数组长度可以不同。我还有另一个文件，$selected_genes，它只有一列有5000多个独特的基因。我应该计算%go_hash的每个数组中的基因数量，并找到$selected_genes列表中存在的每个数组的基因数量(如果没有重叠，0应该在那里)。然后将这两个数字添加到哈希中相应数组的末尾，并创建新的$output文件。当打印到这个输出文件时，除了一个thinkg之外，最终结果中的所有内容都很好。反变量$count，即与$selected_genes重叠的每个阵列的基因数目，一直是导致0的原因。(实际上有许多重叠，所以不应该总是0 )。我尝试了很多方法，但没有改变，特别是在15到19之间的代码行中。也许问题在代码的其他部分。
我在哪里搞错了？有人能纠正我吗？在此之前，谢谢您的评论/帮助。
回答 1
Stack Overflow用户
回答已采纳
发布于 2013-09-05 22:56:45
open my $go_file, "gene_associations_go_human.txt" or die "Can't open the file!";
open my $selected_genes, "my_selected_genes.txt" or die "Can't open the file!";
my %sel_genes = map { chomp; $_ => 1 } <$selected_genes>;
my %result;
while(<$go_file>){
chomp;
my @r = split /\t/;
$result{$r[4]} = { count1=>0, count2=>0, data=>{}} if not defined $result{$r[4]};
$result{$r[4]}->{count1}++;
$result{$r[4]}->{data}->{$r[2]} = defined $sel_genes{$r[2]} ? 1 : 0;
$result{$r[4]}->{count2}++ if defined $sel_genes{$r[2]};
}
for my $r (keys %result) {
print $r . "\t" . join("\t", keys %{$result{$r}->{data}}) . "\t" . $result{$r}->{count1} . "\t" . $result{$r}->{count2} . "\n";
}页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接：
https://stackoverflow.com/questions/18618736
复制相关文章
相似问题
腾讯云开发者
