统计蛋白质互作hug gene degree（每个基因连接数）

017中大普外

用了stringDB这个包做完了蛋白质互作后得到了如下表格，
现在想统计每个目标基因链接的其他数目，像其他文章里找hug gene，尝试其他办法未解，无奈准备写一个小脚本解决。

于是就想用perl搞定。（就是统计该文本中第二列相同gene的数目）
这是我的初始代码：[mw_shl_code=perl,true]use warnings;
use strict;

my $file=$ARGV[0];

open(RF,"$file") or die $!;

my @arry=();
my @spli=();
while(my $line=<RF>){
         chomp($line);
         @spli=split(/\t/,$line);
         my $lim=$spli[1];
         my @arry=$lim;
         foreach $_(@arry){
                 print $_."\n";
         }
}

foreach $_(@arry){
        my $count=0
        
}[/mw_shl_code]
我把第二列写到了同一数组@arry中，但无后续思路，太菜是这样的。。。
随后请教群里大神，大神用哈希很快就解决该问题了
附上大神代码：[mw_shl_code=perl,true]use strict;
use warnings;

open my $gene,"DEGppi.txt" or die $!;
my %hash_gene;
while (<$gene>) {
        chomp;
        next if(m/^#/);
        my @tmp = split /\t/,$_;
        $hash_gene{$tmp[1]}++;
}
close $gene;

open my $out,"> outfile.txt" or die $!;
foreach (sort {$a cmp $b} keys %hash_gene) {
        # body...
        print $out "$_\t$hash_gene{$_}\n"
}
close $out;[/mw_shl_code]
短短几行搞定。。。强。。。

017中大普外

此贴已错成狗，千万勿抄

017中大普外

盘龙 发表于 2017-1-4 14:25
很好 发帖的时候 可以点击这里插入代码 方便观看

OK,以后发会注意的

Jimmy

你们两个的代码量，是一样的呀！

017中大普外

Jimmy 发表于 2017-1-4 14:35
你们两个的代码量，是一样的呀！

我没有写完啊，我后面不会了。。。写到数组里面之后。。。

张小凡_2016

[mw_shl_code=perl,true]
use strict;
use warnings;

open my $gene,"DEGppi.txt" or die $!;
open my $out,"> outfile.txt" or die $!;
my %hash_gene;
map{$hash_gene{$_}++} (split /\t/, $_)[1] while(<$gene>);

foreach (sort {$a cmp $b} keys %hash_gene) {
        print $out "$_\t$hash_gene{$_}\n"
}
close $gene;
close $out;[/mw_shl_code]
更加简洁的版本，大家可以参考一下

173318328

你好，您在使用stringDB这个包的时候对输入文件有什么要求吗？为什么我这边总显示输入文件读取错误！

017中大普外

因为我只统计了第一列重复的个数，结果我发现这是错的，今天导入cytoscape后发现不对，因为第二列和第一列是不重复的。。。
[mw_shl_code=python,true]import os
import re
os.chdir('F:/RNA/ppi/HUG GENE')
dic={}
with open("1233.txt",'r') as f:
    for line in f:
        line=line.rstrip().split('\t')
        gene1=line[0]
        gene2=line[1]
        #print(gene2)
        if not gene1 in dic:
            dic[gene1]=1
        else:
            dic[gene1]+=1
        if not gene2 in dic:
            dic[gene2]=1
        else:
            dic[gene2]+=1
print(dic)
            
            
with open('1111.txt','w') as f3:
    for i,t in dic.items():
        f3.write(i+'\t'+str(t)+"\n")  [/mw_shl_code]

321wangke321

用R最方便~library(tidyverse) 里面的group_by，然后sammarize就好了

tss940504

用R 两句话
my_degree<-sort(table(c(all_ppi$from,all_ppi$to)),decreasing = T)
tmp=data.frame(STRING_id=names(my_degree),degrees=as.numeric(my_degree))