|
0
|
1 #!/usr/bin/perl
|
|
|
2
|
|
|
3 if(@ARGV<2) { print "usage: $0 iNIR_file iEUMA_file > output.iFVKM\n"; exit; }
|
|
|
4 my ($iNIR_file,$iEUMA_file) = @ARGV;
|
|
|
5
|
|
|
6 open IN, $iNIR_file or die "Can't open iNIR_file $sample_iNIR_file\n";
|
|
|
7 <IN>; # skip header
|
|
|
8 while(<IN>){
|
|
|
9 chomp;
|
|
|
10 split/\t/;
|
|
|
11 $gene = $_[0]."\t".$_[1];
|
|
|
12 $RNAseq{$gene}=$_[2];
|
|
|
13 }
|
|
|
14 close IN;
|
|
|
15
|
|
|
16 open IN, $iEUMA_file or die "Can't open iEUMA file $iEUMA_file\n";
|
|
|
17 while(<IN>){
|
|
|
18 chomp;
|
|
|
19 split/\t/;
|
|
|
20 $gene = $_[0]."\t".$_[1];
|
|
|
21 $iEUMA{$gene}=$_[2];
|
|
|
22 }
|
|
|
23 close IN;
|
|
|
24
|
|
|
25
|
|
|
26 print "gene\tisoform\tiFVK\tread.count\tiEUMA(bp)\n";
|
|
|
27 for my $gene (sort keys %iEUMA){
|
|
|
28 next if $iEUMA{$gene} ==0;
|
|
|
29 if(!exists $RNAseq{$gene}){ $RNAseq{$gene} = 0; }
|
|
|
30 my $iFVKM = $RNAseq{$gene} / $iEUMA{$gene} * 1000;
|
|
|
31 print "$gene\t$iFVKM\t$RNAseq{$gene}\t$iEUMA{$gene}\n";
|
|
|
32 }
|
|
|
33
|
|
|
34
|
|
|
35
|
