查看: 13661|回复: 5

[Perl语言] read count 转 FPKM 和cuflinks计算的值不一样

[复制链接]

中华鲟

Rank: 5Rank: 5

主题
13
注册时间
2018.1.3
在线时间
338 小时

热心会员活跃会员


发表于 2019.7.16 20:43:59 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

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

x
请各位高手看下我哪里是不是出错了,先谢过。
read count 文件如下形式
geneid    exon_Len    sampleA    sampleB    sampleC
1    124    66    55    34
2    866    65    43    23
3    649    34    76    88
我的perl脚本如下:
[Perl] 纯文本查看 复制代码
use strict;

open A,"$ARGV[0]";

<A>;
my @colsum;
while(<A>){
        chomp;
        my @line=split;
        shift @line;
        shift @line;
        for(0..$#line){
                $colsum[$_]+=$line[$_];
        }
}

open B,"$ARGV[0]";
my $head=<B>;
chomp($head);
print "$head\n";
while(<B>){
        chomp;
        my @line=split;
        my $gene=shift @line;
        my $length=shift @line;
        print "$gene\t";
        for(my $i=0;$i<@line;$i++){
                my $fpkm=$line[$i]*1000000*1000/($colsum[$i]*$length);
                print "$fpkm\t";
        }
        print "\n";
}

新的一天加油!
回复

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
13
注册时间
2018.1.3
在线时间
338 小时

热心会员活跃会员


 楼主| 发表于 2019.7.16 20:49:47 | 显示全部楼层
说明一下,我的read count值是用featureCounts 跑出来的,count值取的是结果文件最后一列的值,exon长度取的是结果文件倒数第二列的值
回复 支持 反对

使用道具 举报

中华鲟

Rank: 5Rank: 5

主题
13
注册时间
2018.1.3
在线时间
338 小时

热心会员活跃会员


 楼主| 发表于 2019.7.18 09:16:07 | 显示全部楼层
没有人啊
新的一天加油!
回复

使用道具 举报

管理员

Rank: 15Rank: 15Rank: 15Rank: 15Rank: 15

主题
402
注册时间
2018.4.19
在线时间
896 小时

推广达人宣传达人


发表于 2019.7.18 10:01:17 | 显示全部楼层
帮顶
回复

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
3
注册时间
2017.4.12
在线时间
52 小时

发表于 2020.2.15 21:57:26 | 显示全部楼层
大裤衩 发表于 2019.7.16 20:49
说明一下,我的read count值是用featureCounts 跑出来的,count值取的是结果文件最后一列的值,exon长度取 ...

使用counts计算FPKM值,取exon 长度直接使用featureCounts中的length应该是不对的,我看了featureCount软件中的length是全长,而不是exon的
新的一天加油!
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
3
注册时间
2017.4.12
在线时间
52 小时

发表于 2020.2.16 20:40:21 | 显示全部楼层
杜康华 发表于 2020.2.15 21:57
使用counts计算FPKM值,取exon 长度直接使用featureCounts中的length应该是不对的,我看了featureCount软 ...

你可以使用perl写一个脚本提取外显子长度,再进行计算。如果可以分享Perl提取外显子长度的脚本的话,那就非常好了,对我们这些不会写脚本的童鞋就非常友好了,谢谢
新的一天加油!
回复 支持 反对

使用道具 举报

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

本版积分规则

快速回复 返回顶部 返回列表