-
Notifications
You must be signed in to change notification settings - Fork 3
/
probe_score.pl
127 lines (112 loc) · 3.43 KB
/
probe_score.pl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
#!/usr/bin/perl -w
use strict;
use warnings;
use Getopt::Long;
use Data::Dumper;
use FindBin qw($Bin $Script);
use File::Basename qw(basename dirname);
require "$Bin/score.pm";
require "$Bin/common.pm";
require "$Bin/self_lib.pm";
require "$Bin/math.pm";
require "$Bin/io.pm";
my $BEGIN_TIME=time();
my $version="1.0.0";
#######################################################################################
# ------------------------------------------------------------------
# GetOptions
# ------------------------------------------------------------------
my ($foligo,$fkey,$outdir);
my $min_len=20;
my $max_len=36;
my $opt_tm=70;
my $NoFilter;
my $Methylation;
my $NoSpecificity;
GetOptions(
"help|?" =>\&USAGE,
"i:s"=>\$foligo,
"maxl:s"=>\$max_len,
"minl:s"=>\$min_len,
"opttm:s"=>\$opt_tm,
"k:s"=>\$fkey,
"Methylation:s"=>\$Methylation,
"NoSpecificity:s"=>\$NoSpecificity,
"NoFilter:s"=>\$NoFilter,
"od:s"=>\$outdir,
) or &USAGE;
&USAGE unless ($foligo and $fkey);
$outdir||="./";
`mkdir $outdir` unless (-d $outdir);
$outdir=AbsolutePath("dir",$outdir);
my $max_bound_num=100;
open(F,">$outdir/$fkey.probe.filter") or die $!;
open(O, ">$outdir/$fkey.probe.score") or die $!;
#$slen, $stm, $sself, $ssnp, $spoly, $sbound, $sCGd
print O "##ScoreInfo: ". &score_des("Probe", $Methylation), "\n";
print O &score_head($Methylation, $NoSpecificity), "\n";
open(P, $foligo) or die $!;
while(<P>){
chomp;
#($id, $seq, $len), ($tm, $gc, $hairpin, $dimert, $dimers, $snp, $poly, $cpgs, $cs), ($bnum, $btm, $binfo);
my ($abase, $afeature, $ameth, $aspec, $bnumtm)=&read_evaluation_info(0, $_, $Methylation, 1);
my ($id, $seq, $len)=@{$abase};
my $tm=$afeature->[0];
## filter
if(!defined $NoFilter && ($tm<$opt_tm-8 || $tm>$opt_tm+8)){
print F "TM\t$_\n";
next;
}
my ($is_G5, $CGd, $Gnum) = &G_content($seq);
if(!defined $NoFilter && $is_G5){
print F "5endG\t$_\n";
next;
}
if(!defined $NoFilter && $Gnum>=3){
print F "dupGnum\t$_\n";
next;
}
# if($bnum>$max_bound_num){
# print F "Bound\t$_\n";
# next;
# }
## score
my ($sadd, $score_info);
if(!defined $Methylation){
($sadd, $score_info)=&probe_oligo_score($opt_tm, $len, @{$afeature}, @{$ameth}, $bnumtm, $is_G5, $CGd, $Gnum);
}else{
($sadd, $score_info)=&probe_meth_score($opt_tm, $len, @{$afeature}, @{$ameth}, $bnumtm, $is_G5, $CGd, $Gnum);
}
my @out=(@{$abase}, $sadd, $score_info, @{$afeature}, @{$ameth}, @{$aspec});
print O join("\t", @out),"\n";
}
close(P);
close(O);
close(F);
#######################################################################################
print STDOUT "\nDone. Total elapsed time : ",time()-$BEGIN_TIME,"s\n";
#######################################################################################
# ------------------------------------------------------------------
# sub function
# ------------------------------------------------------------------
sub USAGE {#
my $usage=<<"USAGE";
Program:
Version: $version
Contact:zeng huaping<huaping.zeng\@genetalks.com>
Usage:
Options:
-i <file> Input oligo evaluation file, forced
-k <str> Key of output file, forced
--Methylation Design methylation oligos
--NoFilter Not filter any probes
--NoSpecificity Not evalue specificity of oligos
-minl <int> min len of probe, [$min_len]
-maxl <int> max len of probe, [$max_len]
-opttm <int> opt tm of probe, [$opt_tm]
-od <dir> Dir of output file, default ./
-h Help
USAGE
print $usage;
exit;
}