搜索
查看: 3862|回复: 9

[其他] 统计蛋白质互作hug gene degree(每个基因连接数)

[复制链接]

6

主题

36

帖子

465

积分

中级会员

Rank: 3Rank: 3

积分
465
发表于 2017-1-4 14:17:30 | 显示全部楼层 |阅读模式
用了stringDB这个包做完了蛋白质互作后得到了如下表格,
现在想统计每个目标基因链接的其他数目,像其他文章里找hug gene,尝试其他办法未解,无奈准备写一个小脚本解决。

于是就想用perl搞定。(就是统计该文本中第二列相同gene的数目)
这是我的初始代码:
[Perl] 纯文本查看 复制代码
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
        
}

我把第二列写到了同一数组@arry中,但无后续思路,太菜是这样的。。。
随后请教群里大神,大神用哈希很快就解决该问题了
附上大神代码:
[Perl] 纯文本查看 复制代码
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;

短短几行搞定。。。强。。。

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
回复

使用道具 举报

6

主题

36

帖子

465

积分

中级会员

Rank: 3Rank: 3

积分
465
 楼主| 发表于 2017-2-28 21:41:17 | 显示全部楼层
此贴已错成狗,千万勿抄
回复 支持 1 反对 0

使用道具 举报

6

主题

36

帖子

465

积分

中级会员

Rank: 3Rank: 3

积分
465
 楼主| 发表于 2017-1-4 14:30:41 | 显示全部楼层
盘龙 发表于 2017-1-4 14:25
很好 发帖的时候 可以点击这里插入代码 方便观看

OK,以后发会注意的
回复 支持 反对

使用道具 举报

634

主题

1182

帖子

4030

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
4030
发表于 2017-1-4 14:35:49 | 显示全部楼层
你们两个的代码量,是一样的呀!
你这个问题很复杂,需要打赏,请点击 http://www.bio-info-trainee.com/donate 进行打赏,谢谢
回复 支持 反对

使用道具 举报

6

主题

36

帖子

465

积分

中级会员

Rank: 3Rank: 3

积分
465
 楼主| 发表于 2017-1-4 14:42:35 | 显示全部楼层
Jimmy 发表于 2017-1-4 14:35
你们两个的代码量,是一样的呀!

我没有写完啊,我后面不会了。。。写到数组里面之后。。。
回复 支持 反对

使用道具 举报

0

主题

1

帖子

56

积分

注册会员

Rank: 2

积分
56
发表于 2017-1-4 15:58:34 | 显示全部楼层
本帖最后由 张小凡_2016 于 2017-1-4 16:01 编辑

[Perl] 纯文本查看 复制代码
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;

更加简洁的版本,大家可以参考一下
回复 支持 反对

使用道具 举报

0

主题

5

帖子

87

积分

注册会员

Rank: 2

积分
87
发表于 2017-1-24 18:16:06 | 显示全部楼层
你好,您在使用stringDB这个包的时候对输入文件有什么要求吗?为什么我这边总显示输入文件读取错误!
回复 支持 反对

使用道具 举报

6

主题

36

帖子

465

积分

中级会员

Rank: 3Rank: 3

积分
465
 楼主| 发表于 2017-2-28 22:10:29 | 显示全部楼层
因为我只统计了第一列重复的个数,结果我发现这是错的,今天导入cytoscape后发现不对,因为第二列和第一列是不重复的。。。
[Python] 纯文本查看 复制代码
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")  
回复 支持 反对

使用道具 举报

0

主题

8

帖子

65

积分

注册会员

Rank: 2

积分
65
发表于 2017-3-4 17:57:06 | 显示全部楼层
用R最方便~library(tidyverse) 里面的group_by,然后sammarize就好了
回复 支持 反对

使用道具 举报

5

主题

32

帖子

484

积分

管理员

Rank: 9Rank: 9Rank: 9

积分
484
发表于 2017-6-6 20:08:55 | 显示全部楼层
用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))
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|手机版|小黑屋|生信技能树 ( 粤ICP备15016384号  

GMT+8, 2019-9-16 09:09 , Processed in 0.039078 second(s), 31 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.