Bioinf
view release on metacpan or search on metacpan
my @poss_sub_dir_heads=('.'); ## <<<<------- This is critically important, when 'D' opt is not used!
if($each_seq_name=~/(\S+)_\d+\-\d+$/){
$each_seq_name_range=$each_seq_name;
$each_seq_name=$1;
@name_types=($each_seq_name, $each_seq_name_range);
}else{
@name_types=($each_seq_name);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Here I take chars from the sequ names, as dirs have fragments of chars
#_______________________________________________________________________________
for($s=1; $s <= $subdir_char_size ; $s++){ ## here, number 2 indicates, I check single or 2 char sub dir names
$sub_dir_head= substr($seq_names[$j], 0, $s);
push(@poss_sub_dir_heads, "\L$sub_dir_head") if (-d "\L$sub_dir_head" );
push(@poss_sub_dir_heads, "\U$sub_dir_head") if (-d "\U$sub_dir_head" );
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking all the possible subdirectories to crop all the sso files
#_______________________________________________________________________________
FOR1: for($p=0; $p <= @poss_sub_dir_heads; $p++){ ## Default has '.' will make things like '././fam_8_8.pbla.gz'
$num_of_seq_member=@list=split(/ +/, $clus{$keys[$i]}); # @list has (HIU001, HI002, HI333, MJ111, etc)
print "# $0: clu_to_sso_to_msp: No. of seq member=$num_of_seq_member after split \n" if $verbose;
FOR0: for($j=0; $j < @list; $j++){
my($sub_dir_head, $file_name_low, $file_name_up, $file_name_prot_low, @sub_dir_heads,
$file_name_prot_up, $file_name_low_gz, $file_name_up_gz,
$file_name_prot_low_gz, $file_name_prot_up_gz);
$each_seq_name=$list[$j];
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Here I take chars from the sequ names, as dirs have fragments of chars
#_______________________________________________________________________________
for($s=1; $s <=2 ; $s++){ ## here, number 2 indicates, I check single or 2 char sub dir names
$sub_dir_head= substr($list[$j], 0, $s);
push(@sub_dir_heads, "\L$sub_dir_head") if (-d "\L$sub_dir_head" );
push(@sub_dir_heads, "\U$sub_dir_head") if (-d "\U$sub_dir_head" );
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking all the possible subdirectories to crop all the sso files
#_______________________________________________________________________________
FOR1: for($p=0; $p < @sub_dir_heads; $p++){
$num_of_seq_member=@list=split(/ +/, $clus{$keys[$i]}); # @list has (HIU001, HI002, HI333, MJ111, etc)
print "# $0: clu_to_sso_to_msp: No. of seq member=$num_of_seq_member after split \n" if $verbose;
FOR0: for($j=0; $j < @list; $j++){
my($sub_dir_head, $file_name_low, $file_name_up, $file_name_prot_low, @sub_dir_heads,
$file_name_prot_up, $file_name_low_gz, $file_name_up_gz,
$file_name_prot_low_gz, $file_name_prot_up_gz);
$each_seq_name=$list[$j];
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Here I take chars from the sequ names, as dirs have fragments of chars
#_______________________________________________________________________________
for($s=1; $s <=2 ; $s++){ ## here, number 2 indicates, I check single or 2 char sub dir names
$sub_dir_head= substr($list[$j], 0, $s);
push(@sub_dir_heads, "\L$sub_dir_head") if (-d "\L$sub_dir_head" );
push(@sub_dir_heads, "\U$sub_dir_head") if (-d "\U$sub_dir_head" );
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking all the possible subdirectories to crop all the sso files
#_______________________________________________________________________________
FOR1: for($p=0; $p < @sub_dir_heads; $p++){
}else{ return(\$largest_file); }
}
#______________________________________________________________
# Title : get_sequence_complexity
# Usage : print "\n", ${&get_sequence_complexity(\$seq)};
# Function : caculates the single sequence's sequence complexity
# If the seq given is larger than 20, it divides it into
# frags of 20 aa and gets the average of it.
# Example : ${&get_sequence_complexity(\$seq)};
# while $seq='TTTTTACDEFGHIKLMNPQRSTVWYAAAAACCCADFADFA'
# Warning :
# Keywords : sequence_complexity, calc_sequence_complexity,
# calc_seq_complexity, get_seq_complexity, seg
# Options : _ for debugging.
# # for debugging.
# 'w=' for window size as the first arg
# Returns : Ref. for a scalar digit.
# Argument : ref. of string.
# Category :
# Version : 1.3
#--------------------------------------------------------------
sub get_sequence_complexity{
my ($complexity, @seq,$i, $j, @frag);
my $win=20;
if(ref($_[0]) eq 'ARRAY'){
@seq=@{$_[0]};
}else{
$seq=${$_[0]} || $_[0];
@seq=split(//, $seq);
}
if(defined($_[1])){ $win=${$_[1]} || $_[1]; }
if(@seq <= $win){
my (%seq, @keys);
for($i=0; $i< @seq; $i++){
$seq{$seq[$i]}++;
}
@keys= keys %seq;
$complexity=@keys/@seq;
}else{
my @frag=@{÷_array(\@seq, "s=$win")};
my @complexity=();
for($i=0; $i < @frag; $i++){
my (%seq, @keys);
my @arr=@{$frag[$i]};
for($j=0; $j< @arr; $j++){
$seq{$arr[$j]}++;
}
@keys=keys %seq;
push(@complexity, @keys/$win);
}
$complexity=${&array_average(\@complexity)};
}
return(\$complexity);
}
#"""""""""""""""""< handle_arguments{ head Ver 4.1 >"""""""""""""""""""
my(@A)=&handle_arguments(@_);my($num_opt)=${$A[7]};my($char_opt)=${$A[8]};
my(@hash)=@{$A[0]};my(@file)=@{$A[4]};my(@dir)=@{$A[3]};my(@array)=@{$A[1]};
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my ($connect_gap, @seq_frag, %digitized, $verbose, %hash, $best_block_opt);
my $margin=3;
my $threshold=0.8;
my $min_seqlet_size=25;
$connect_gap=5;
my @vars=keys %vars;
for($i=0; $i< @vars; $i++){
if($vars[$i] eq 'm'){
$margin=$margin_ori=$vars{$vars[$i]};
}elsif($vars[$i] eq 't'){
$threshold=$vars{$vars[$i]};
for($e=0; $e< @keys; $e++){
my $range="$keys[$e]\-$blocks{$keys[$e]}";
my $seq_let_leng=$blocks{$keys[$e]} - $keys[$e] + 1;
if($seq_let_leng < $min_seqlet_size){
next;
}else{
push(@RANGE, $range);
}
}
@seq_frag=&get_seq_fragments(\%hash, @RANGE,
"l=$min_seqlet_size", "$range_in_name");
return(\@seq_frag);
}
#______________________________________________________________
# Title : add_columns
# Usage :
# Function :
# Example :
# Warning : if the attached name is too long(over 12 char),
# it changes to 'Added_upX' while X is a numb.
for($j=$i; $j < ($factor+$i); $j++){
$temp_sum+=$ARRAY[$j]
}
push(@out, int($temp_sum/$factor) );
}
return(\@out);
}
#___________________________________________________________
# Title : get_seq_fragments
# Usage : @seq_frag=&get_seq_fragments(\%msf, @RANGE);
# Function : gets sequence(string) segments with defined
# ranges.
# Example :
# %test=('seq1', '1234AAAAAAAAAAAaaaaa', 'seq2', '1234BBBBBBB');
# @range = ('1-4', '5-8');
#
# %out = %{&get_seq_fragments(\%test, \@range)};
# %out => (seq1_5-8 AAAAA
# seq2_5-8 BBBBB
# seq1_1-4 1234
# seq2_1-4 1234 )
#
# Warning :
# Keywords : get_sequence_fragments,
# Options : _ for debugging.
# # for debugging.
# l= for min seqlet length
# r for adding ranges in the seq names
#
# Returns :
# Argument :
# Category :
# Version : 1.8
#-------------------------------------------------------
sub get_seq_fragments{
#"""""""""""""""""< handle_arguments{ head Ver 4.1 >"""""""""""""""""""
my(@A)=&handle_arguments(@_);my($num_opt)=${$A[7]};my($char_opt)=${$A[8]};
my(@hash)=@{$A[0]};my(@file)=@{$A[4]};my(@dir)=@{$A[3]};my(@array)=@{$A[1]};
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my $no_range_in_name=1;
for($i=0; $i< @vars; $i++){
if($vars[$i] eq 'l'){
$min_seqlet_size=$vars{$vars[$i]};
}
}
if($char_opt=~/v/){ print "\n \$char_opt is $char_opt @char_opt\n"; }
if($char_opt=~/n/){ $no_range_in_name = 1 }
if($char_opt=~/r/){ $no_range_in_name = 0 }
print "\nget_seq_fragments \$no_range_in_name is $no_range_in_name \n";
for($i=0; $i< @hash; $i++){
my (%out_frag, $frag_name, $range_start, $range_end, @out_hash);
my %seqs = %{$hash[$i]};
my @names = keys %seqs;
if(@names==1){
for($j=0; $j < @names; $j++){
my $seq_name = $names[$j];
my $seq = $seqs{$seq_name};
for($k=0; $k< @range; $k++){
my $range = $range[$k];
if($no_range_in_name==1){
$frag_name = "$seq_name";
}else{
$frag_name = "$seq_name\_$range";
}
#if(length($frag_name)>14 ){
# $frag_name ='x'."${j}_${range}";
#}
($range_start, $range_end)=$range=~/(\d+\.?\d*)\-(\d+\.?\d*)/;
my $frag_len = $range_end-$range_start+1;
if($frag_len < $min_seqlet_size){
next;
}
my $fragment = substr($seq, $range_start-1, $frag_len);
$out_frag{$frag_name}=$fragment;
}
}
push(@out_hash, \%out_frag);
}elsif(@names > 1){
for($k=0; $k< @range; $k++){
my %out_frag=();
my $range=$range[$k];
($range_start, $range_end)=$range=~/(\d+\.?\d*)\-(\d+\.?\d*)/;
my $frag_len = $range_end-$range_start+1;
if($frag_len < $min_seqlet_size){
next;
}
for($j=0; $j < @names; $j++){
my $seq_name=$names[$j];
my $seq = $seqs{$seq_name};
if($no_range_in_name==1){
$frag_name = "$seq_name";
}else{
$frag_name = "$seq_name\_$range";
}
#if(length($frag_name)>15 ){
# $frag_name ='x'."${j}_${range}";
#}
if($range_start==0){ $range_start++; } ## This is a bugfix
my $fragment = substr($seq, $range_start-1, $frag_len);
$out_frag{$frag_name}=$fragment;
}
push(@out_hash, \%out_frag);
}
}
}
if(@out_hash > 1){ return(@out_hash)
}elsif(@out_hash==1){ return($out_hash[0]) }
}
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
$ref_target_hash = shift @hash;
%target = %{$ref_target_hash};
@names = keys %target;
$name_target = $names[0];
if($name_target =~/\w+(\d+)/){ $ori_target_seq_len = length($1); }
@target_frag = split(/ +/, $target{$name_target} );
for($i =0; $i< @hash; $i ++){
%db=%{$hash[$i]};
@db_name = keys %db;
for($j=0; $j < @db_name; $j ++){
$name = $db_name[$j];
@db_frag = split( / +/, $db{$name} );
for($k=0; $k < @target_frag; $k ++){
if( ($target_frag[$k]=~/H(\d+)/i)&&($db_frag[$k]=~/H(\d+)/i) ){
$simple_match_output{$name}++;
$leng_diff = abs($1 - $2)/15;
$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}elsif( ($target_frag[$k]=~/E(\d+)/i)&&($db_frag[$k]=~/H(\d+)/i) ){
$simple_match_output{$name}--;
#$leng_diff = abs($1 - $2)/10;
#$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}elsif( ($target_frag[$k]=~/H(\d+)/i)&&($db_frag[$k]=~/E(\d+)/i) ){
$simple_match_output{$name}--;
#$leng_diff = abs($1 - $2)/10;
#$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}elsif( ($target_frag[$k]=~/E(\d+)/i)&&($db_frag[$k]=~/E(\d+)/i) ){
$simple_match_output{$name}++;
$leng_diff = abs($1 - $2)/15;
$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}
}
}
}
return(\%simple_match_output);
}
#_______________________________________________________________
# Title : divide_array
# Usage : &show_array(÷_array(\@input, 6));
# Function : divides any array to the denominator given.
# If you give array of 100 elem, with 5, you will
# get 5 arrays with 20 elem each.
# Example :
# Warning :
# Keywords : split_array_into_pieces, split_array, chop_array,
# fragment_array,
# Options : s= for dividing the array with sub array size
# eg) to get 20 elem length sub arrays from
# a big array
# @ar_ref=@{÷_array(\@array, 's=20')};
# Returns :
# Argument :
# Category :
# Version : 1.4
#-----------------------------------------------------------
sub divide_array{
if((@_ ==1)&&($denominator == 0)){
print "\n Denominator is 0, error, setting to 1\n";
$denominator = 1;
}
if($size_div==1){
while(@array){
push(@final_array_ref, [splice(@array, 0, $size)]);
}
}else{
my $frag_ar_size = int(@array/$denominator);
if($debug eq 1){ print "\n Frag arr size is : $frag_ar_size \n" }
$remaining = @array % $denominator;
if($debug eq 1){ print "\n Remnant elem size is : $remaining \n" }
for($i=0; $i < $denominator; $i++){
if($remaining > 0){
push(@final_array_ref, [splice(@array, 0, ($frag_ar_size+1),)] );
$remaining --;
}elsif(($remaining == 0)&&(@array>0)){
push(@final_array_ref, [splice(@array, 0, ($frag_ar_size),)] );
}
}
}
return(\@final_array_ref);
}
#__________________________________________________________________________
# Title : split_fasta_files
my($division_factor, $base_name);
$division_factor=2;
if($vars{'d'}=~/(\d+)/){ $division_factor=$1 };
for($i=0; $i< @file; $i++){
open(FILE, $file[$i]);
$base_name=${&get_base_names($file[$i])};
$extension=${&get_extension_names($file[$i])};
$total_line_num=@lines=<FILE>;
$splited_frag_size=int(@lines/$division_factor);
for($j=0; $j< $division_factor; $j++){
#$frag_file="$base_name\_s${j}\.split${j}";
$frag_file="$base_name\_s${j}\.$extension";
push(@splited_files, $frag_file);
open(FRAGMENT_FILE, ">$frag_file");
for($k=0; $k<= $splited_frag_size; $k++){
print FRAGMENT_FILE $lines[$k];
}
splice(@lines, 0, $splited_frag_size);
}
print FRAGMENT_FILE @lines; # Writing the very last remaining lines
close (FRAGMENT_FILE);
}
return(\@splited_files);
}
#_______________________________________________________________
# The output will be 'seq_1_half', '1234578'
# 'seq_2_half', '9ABCDEF'
# 'seq_3_half', 'HIJKLMN'
# Warning :
# Keywords : divide_string, split_string, chop_string, divide_sequence
# split_sequence(look at separate split_sequence sub),
# Options :
# $reverse_second_half=S by S -S
# $reverse_first_half =F by F -F
# $reverse_rest =R by R -R ## reversing all except the first
# $reverse_all =A by A -A # reverse all the fragments
# Returns :
# Argument :
# Category :
# Version : 1.3
#-----------------------------------------------------------
sub split_sequence{
my(@string, $frag_str_size, $remaining, @hash, $i, $j,
@seq_names, $denominator, %input_seq, @final_hash_ref,
$frag, $reverse_second_half, $reverse_first_half,
$reverse_rest, $reverse_all);
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking arguments
#________________________________________________
for($i=0; $i< @_; $i++){
if(ref($_[$i]) eq 'HASH'){
push(@hash, $_[$i]);
}elsif(ref($_[$i]) eq 'SCALAR'){
if(${$_[$i]}=~/^(\d+\.?\d*)$/){ # this can handle fraction number!!
for($i=0; $i< @hash; $i++){
my %input_seq=%{$hash[$i]};
@seq_names=keys %input_seq;
my (%out_seq);
for($k=0; $k < @seq_names; $k++){
$seq_name=$seq_names[$k];
@string = split(//, $input_seq{ $seq_name });
$frag_str_size = int(@string/$denominator);
if($debug eq 1){ print "\n Frag str size is : $frag_str_size \n" }
$remaining = @string % $denominator;
for($j=0; $j < $denominator; $j++){
$frag_number=$j+1;
if($remaining > 0){
$frag=join('', splice(@string, 0, ($frag_str_size+1) ) );
$remaining --;
$seq_name_split="$seq_name\_$frag_number";
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~```
# Before final, set the seq types
#_________________________________________________________________
if($frag_number==1 and $reverse_first_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number==2 and $reverse_second_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number >2 and $reverse_rest){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($reverse_all){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}
$out_seq{$seq_name_split}=$frag;
}elsif( $remaining == 0 and @string>0){
$seq_name_split="$seq_name\_$frag_number";
$frag=join('', splice(@string, 0, $frag_str_size,) );
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~```
# Before final, set the seq types
#_________________________________________________________________
if($frag_number==1 and $reverse_first_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number==2 and $reverse_second_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number >2 and $reverse_rest){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($reverse_all){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}
$out_seq{$seq_name_split}=$frag;
}
}
}
push(@final_hash_ref, \%out_seq);
}
wantarray? \@final_hash_ref : $final_hash_ref[0];
}
# Warning :
# Keywords : divide_string, split_string, chop_string, divide_sequence
# split_sequence(look at separate split_sequence sub),
# Options :
# Returns :
# Argument :
# Category :
# Version : 1.4
#-----------------------------------------------------------
sub divide_string{
my (@array, $i, $j, $denominator, @temp, @string, $frag_str_size,
$remaining, $frag);
for($i=0; $i< @_; $i++){
if(ref($_[$i]) eq 'ARRAY'){
push(@array, @{$_[$i]});
}elsif(ref($_[$i]) eq 'SCALAR'){
if(${$_[$i]} =~/^(\d+)$/){
$denominator = $1;
}else{
push(@array, ${$_[$i]});
}
push(@array, $_[$i]);
}
}
if($denominator == 0){
print "\n Denominator is 0, error, setting to 1\n";
$denominator = 1;
}
for($i=0; $i< @array; $i++){
my @temp;
my @string = split(//, $array[$i]);
$frag_str_size = int(@string/$denominator);
if($debug eq 1){ print "\n Frag str size is : $frag_str_size \n" }
$remaining = @string % $denominator;
for($j=0; $j < $denominator; $j++){
if($remaining > 0){
$frag=join('', splice(@string, 0, ($frag_str_size+1) ) );
push(@temp, $frag);
$remaining --;
}elsif( $remaining == 0 and @string>0){
$frag=join('', splice(@string, 0, $frag_str_size,) );
push(@temp, $frag);
}
}
push(@final_array_ref, \@temp);
}
wantarray? \@final_array_ref : $final_array_ref[0];
}
#____________________________________________________________
# Title : write_html_headbox
# Usage : &write_html_headbox($outfilename, \%entries);
# Example :
# Warning :
# Keywords :
# Options :
# Returns :
# Argument :
# Category :
# Version : 1.1
#-----------------------------------------------------------
sub make_template_from_sec_str{
my(%out, @name, @fragments, %in, $name, $leng, $frag_seq, $name2);
#""""""""""""""""""""""< handle_arguments{ head Ver 1.1 >""""""""""""""""""""""""""""""""
my(@A ) = &handle_arguments( @_ ); my( $num_opt )=${$A[7]};my( $char_opt )=${$A[8]};
my(@hash) =@{$A[0]};my(@file) =@{$A[4]};my(@dir ) =@{$A[3]};my(@array)=@{$A[1]};
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};my(@raw_string)=@{$A[9]};
my($i, $j, $c, $d, $e, $f, $g, $h, $k, $l, $p, $q, $r, $s, $t, $u, $v, $w, $x,$y,$z);
#"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
for($j = 0; $j < @hash ; $j ++){
%in = %{$hash[$j]};
@name = keys %in;
print "\@name is @name\n", if($debug eq 1);
for($t=0; $t < @name; $t++){
$name = $name[$t];
print "\$name is $name\n", if($debug eq 1);
$leng=length($in{$name});
print "$leng\n", if($debug eq 1);
$name2 ="$name"."$leng"; # to attach sequence length
@fragments = split(/_+/, $in{$name});
print "\@fragments is @fragments\n", if($debug eq 1);
for($i = 0; $i < @fragments; $i++){
if($fragments[$i] =~/(\w)\w+/){
$fraglength = length($fragments[$i]);
$frag_seq .= "$1"."$fraglength "; # space is delimiter 'H5 E3 E5 E4'
print "\$frag_seq is $frag_seq\n", if($debug eq 1);
}
$out{$name2}=$frag_seq;
}
}
}
return(\%out);
}
#_______________________________________________________________
# Title : calculate_protein_volume
$scop_class=$2;
$sequence{"${pdb_name}"} = $scop_class;
}
}
}
}
close FILE_1;
}
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
return(\%seq_fragments);
}elsif($char_opt=~/s/i){ # when SIZE return only option is set
@seq_Names=keys %sequence;
for($i=0; $i < @seq_Names; $i++){
$Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
}
return(\%Sizes);
}elsif(@hash >=1){
for($h=0; $h< @hash; $h++){
my %hash=%{$hash[$h]};
my @Keys=keys %hash;
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my (%sequence, %na_sequence, %HASH, @Keys, $seq_found1,
$S_start, $S_end, $seq_found,
$present_seq, @seq_Names, %Sizes, $bare_seq_name, $fasta_seq_idx_file,
%seq_fragments, $get_protein_seq, $reverse_seq, $translation_found);
$get_nucleic_acid_seq='n';
if($char_opt=~/r/){ $reverse_seq='r' }
if(@file<1){
print "\n \@file has less than 1 elem. There is no fileinput for open_embl_files\n";
exit
}
$seq=$na_seq='';
}
close FILE_1;
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`~~~~~~~~~~~~~~~~~~~~~`
# (3) When ranges information is given(via \@range), seq in those ranges are returned
#______________________________________________________________________________________
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
if($reverse_seq){
%seq_fragments=%{&reverse_sequences(\%seq_fragments)};
}
return(\%seq_fragments);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`
# (4) When only size is asked with -s option
#_____________________________________________________________________________
elsif($char_opt=~/s/){ # when SIZE(length of seq) return only option is set
#@seq_Names=keys %sequence;
#for($i=0; $i<@seq_Names; $i++){
# $Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
# }
return(\%Sizes);
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my (%sequence, %HASH, @Keys, $seq_found1, $S_start, $S_end, $seq_found,
$present_seq, @seq_Names, %Sizes, $bare_seq_name, $fasta_seq_idx_file,
%seq_fragments, $reverse_seq);
if($char_opt=~/r/){ $reverse_seq='r' }
if(@file<1){
print "\n \@file has less than 1 elem. There is no fileinput for open_fasta_files\n";
exit
}
print "\n# \@file are @file\n" if $verbose;
}# end of while
}
close FILE_1;
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`~~~~~~~~~~~~~~~~~~~~~`
# (3) When ranges information is given(via \@range), seq in those ranges are returned
#______________________________________________________________________________________
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
if($reverse_seq){
%seq_fragments=%{&reverse_sequences(\%seq_fragments)};
}
return(\%seq_fragments);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`
# (4) When only size is asked with -s option
#_____________________________________________________________________________
elsif($char_opt=~/s/){ # when SIZE(length of seq) return only option is set
#@seq_Names=keys %sequence;
#for($i=0; $i<@seq_Names; $i++){
# $Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
# }
return(\%Sizes);
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
#$| = 1;
my($sort_seq_names, $string, $string_leng, $na,$out_file_name_provided,
$write_pure_seq_only, $show_coil_region, %one_to_three, $each_reliability,
$each_seq_fragment, $each_str_fragment, $end_seq_point, $each_reliability_frag,
$protein_name, $omit_coil_region, $ext, $graphical_rep_of_str,
$show_on_screen_only, $seq_block_size);
my($output_file) ='default_out.pred'; ### when no output file name is given, this is used
$protein_name='pred_temp_name';
my $seq_block_size=80;
if($char_opt=~/r/){ $put_reliability_line='r' } ## putting additional reliability line
if($char_opt=~/c/){ $omit_coil_region='c' } ## putting additional reliability line
if($char_opt=~/g/){ $graphical_rep_of_str='g' } ## putting additional reliability line
if($char_opt=~/s/){ $show_on_screen_only='s' } ## putting additional reliability line
$sec_str=~s/[Cc]/_/g;
}
if($graphical_rep_of_str){
$sec_str=~s/[Hh]/>/g;
$sec_str=~s/[Ee]/=/g;
$sec_str=~s/[Cc]/~/g;
$sec_str=~s/_/ /g;
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~``
# writing down the 80 residue seq fragment
#______________________________________________________
$seq_ruler_line=' . . . . . . . .';
for($i=0; $i< @keys; $i+= $seq_block_size){
$each_seq_ruler_line=substr($seq_ruler_line, 0, @keys-$i-1);
$each_seq_fragment=substr($residue, $i, $seq_block_size);
$each_str_fragment=substr($sec_str, $i, $seq_block_size);
$end_seq_point=$i+length($each_seq_fragment);
unless($show_on_screen_only){
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_seq_ruler_line, ' ');
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n", $i+1, $each_seq_fragment, $end_seq_point);
}else{
printf ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_seq_ruler_line, ' ');
printf ("%5s %-${seq_block_size}s %-5s\n", $i+1, $each_seq_fragment, $end_seq_point);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# if you want to put reliability index line
#____________________________________________
if($put_reliability_line){
$each_reliability_frag=substr($reliability, $i, $seq_block_size);
if(! $show_on_screen_only){
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_str_fragment, ' ');
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_reliability_frag, ' ');
}else{
printf ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_str_fragment, ' ');
printf ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_reliability_frag, ' ');
}
}else{
if(! $show_on_screen_only){
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_str_fragment, ' ');
}else{
printf ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_str_fragment, ' ');
}
}
}
close PRED_WRITE unless $show_on_screen_only;
splice(@hash, $n, 1); $n--;
}
if( $out_file_name_provided != 1){
print "\n\n# You didnt give out file name, $output_file used\n" unless $show_on_screen_only;
}
if( -s $output_file ){
return(@out_hash_ref_list);
}
}
}
#__________________________________________________________________________
# Title : open_hmmls_files
# Usage : %out=%{&open_hmmls_files(\@file)};
# Function : hmmls matches the full length model to target seq. while, hmmfs
# does for fragments as well.
# Example :
# Keywords :
# Options :
# t=$thresh for bits score threshold
# e=$evalue_thresh for bits score threshold
# r for adding ranges
# m for making MSP file format output
# E=Enguiry_name for specifying enquiry seq name rather than 'HMM', the default
# Returns :
# Argument :
$sequence{$1} = $3;
}
}elsif(/^\>?(\S+)[\t ]+(\D\D+) *$/){ ## plain 'seqname ASDASDFAASSEQUECN'
$sequence{$1} = $2;
}else{ next }
}
}
close FILE_1;
}
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
return(\%seq_fragments);
}elsif($char_opt=~/s/i){ # when SIZE return only option is set
@seq_Names=keys %sequence;
for($i=0; $i < @seq_Names; $i++){
$Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
}
return(\%Sizes);
}elsif(@hash >=1){
for($h=0; $h< @hash; $h++){
my %hash=%{$hash[$h]};
my @Keys=keys %hash;
}
}
}
}
} #<<<<<-- for($i=0; $i< @file; $i++){
if($simple_default_out_flag >=1){ %matched_seq=%simple_default_output; }
if($simple_default_out_flag > 1){ print "\n# You have asked me to open more than one file. I merge output into one hash\n"; }
if(defined(@range) && ($char_opt !~ /n/)){ ## m opt is for matched seq names only
%seq_fragments=%{&get_seq_fragments(\%matched_seq, \@range)};
return(\%seq_fragments);
}elsif($char_opt=~/s/){ # when SIZE return only option is set
@seq_Names=keys %matched_seq;
for($i=0; $i<@seq_Names; $i++){
$Sizes{$seq_Names[$i]}=length($matched_seq{$seq_Names[$i]});
}
return(\%Sizes);
}else{
if($char_opt=~/n/){
return(\%matched_seq_names);
}else{
get_median get_median get_msp_enquiry_sequence get_msp_matched_sequence
get_msp_range get_multiple_array_entry get_occurances_of_char get_occurances_of_shift_type_hash
get_occurances_of_shift_type_hash_all get_overlapping_range get_overlapping_seq_match_size get_pair_homol_array
get_pair_homol_hash get_path_dirs_from_ENV get_pdb_file_start_number get_peptide_occurance
get_percent_homo_hash get_percent_homol_arr get_percentage get_perl_keywords
get_posi_diff get_posi_diff_abs get_posi_diff_and_rms_hash get_posi_diff_hash
get_posi_rates_hash_out get_posi_rates_hash_out get_posi_rates_hash_out_compact get_posi_rates_hash_out_jp
get_posi_rates_hash_out_msf get_posi_sans_gaps get_posi_shift_hash get_posi_shift_hash_rms
get_posi_shift_rate get_posi_shift_rms_hash get_posi_shift_rms_whole get_position_shift_rate
get_probable_half get_pwd_dir get_pwd_dir_name get_residue_error_rate
get_scop_correcting_pairs get_sd_of_length_diff get_segment_shift_rate get_seq_fragments
get_seq_hash_sans_gaps get_seq_identity get_seqblock get_sequence_complexity
get_sequence_number get_shortest_str_size get_smallest_file get_stat_FASTA_search_result_in_msp_0_files
get_sub_hash get_subdir_names get_subroutine_calls get_time
get_total_memory_size_in_linux get_unix_shell_name get_whole_pwd_path get_windows_cs_rate_array
get_windows_sc_rate_array get_wrong_segment_rate handle_arguments handle_arguments_old
hash_average hash_catenate hash_chk hash_common
hash_common2 hash_common_by_keys hash_no_common hash_output_chk
hash_stat_for_all hash_substract_by_keys hostname if_file_older_than_x_days
import_ENV_vars initialize_code insert_gaps_in_seq_hash insert_lines_anywhere
interm_lib_search is_html isroman key_ready
my @poss_sub_dir_heads=('.'); ## <<<<------- This is critically important, when 'D' opt is not used!
if($each_seq_name=~/(\S+)_\d+\-\d+$/){
$each_seq_name_range=$each_seq_name;
$each_seq_name=$1;
@name_types=($each_seq_name, $each_seq_name_range);
}else{
@name_types=($each_seq_name);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Here I take chars from the sequ names, as dirs have fragments of chars
#_______________________________________________________________________________
for($s=1; $s <= $subdir_char_size ; $s++){ ## here, number 2 indicates, I check single or 2 char sub dir names
$sub_dir_head= substr($seq_names[$j], 0, $s);
push(@poss_sub_dir_heads, "\L$sub_dir_head") if (-d "\L$sub_dir_head" );
push(@poss_sub_dir_heads, "\U$sub_dir_head") if (-d "\U$sub_dir_head" );
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking all the possible subdirectories to crop all the sso files
#_______________________________________________________________________________
FOR1: for($p=0; $p <= @poss_sub_dir_heads; $p++){ ## Default has '.' will make things like '././fam_8_8.pbla.gz'
$num_of_seq_member=@list=split(/ +/, $clus{$keys[$i]}); # @list has (HIU001, HI002, HI333, MJ111, etc)
print "# $0: clu_to_sso_to_msp: No. of seq member=$num_of_seq_member after split \n" if $verbose;
FOR0: for($j=0; $j < @list; $j++){
my($sub_dir_head, $file_name_low, $file_name_up, $file_name_prot_low, @sub_dir_heads,
$file_name_prot_up, $file_name_low_gz, $file_name_up_gz,
$file_name_prot_low_gz, $file_name_prot_up_gz);
$each_seq_name=$list[$j];
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Here I take chars from the sequ names, as dirs have fragments of chars
#_______________________________________________________________________________
for($s=1; $s <=2 ; $s++){ ## here, number 2 indicates, I check single or 2 char sub dir names
$sub_dir_head= substr($list[$j], 0, $s);
push(@sub_dir_heads, "\L$sub_dir_head") if (-d "\L$sub_dir_head" );
push(@sub_dir_heads, "\U$sub_dir_head") if (-d "\U$sub_dir_head" );
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking all the possible subdirectories to crop all the sso files
#_______________________________________________________________________________
FOR1: for($p=0; $p < @sub_dir_heads; $p++){
$num_of_seq_member=@list=split(/ +/, $clus{$keys[$i]}); # @list has (HIU001, HI002, HI333, MJ111, etc)
print "# $0: clu_to_sso_to_msp: No. of seq member=$num_of_seq_member after split \n" if $verbose;
FOR0: for($j=0; $j < @list; $j++){
my($sub_dir_head, $file_name_low, $file_name_up, $file_name_prot_low, @sub_dir_heads,
$file_name_prot_up, $file_name_low_gz, $file_name_up_gz,
$file_name_prot_low_gz, $file_name_prot_up_gz);
$each_seq_name=$list[$j];
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Here I take chars from the sequ names, as dirs have fragments of chars
#_______________________________________________________________________________
for($s=1; $s <=2 ; $s++){ ## here, number 2 indicates, I check single or 2 char sub dir names
$sub_dir_head= substr($list[$j], 0, $s);
push(@sub_dir_heads, "\L$sub_dir_head") if (-d "\L$sub_dir_head" );
push(@sub_dir_heads, "\U$sub_dir_head") if (-d "\U$sub_dir_head" );
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking all the possible subdirectories to crop all the sso files
#_______________________________________________________________________________
FOR1: for($p=0; $p < @sub_dir_heads; $p++){
}else{ return(\$largest_file); }
}
#______________________________________________________________
# Title : get_sequence_complexity
# Usage : print "\n", ${&get_sequence_complexity(\$seq)};
# Function : caculates the single sequence's sequence complexity
# If the seq given is larger than 20, it divides it into
# frags of 20 aa and gets the average of it.
# Example : ${&get_sequence_complexity(\$seq)};
# while $seq='TTTTTACDEFGHIKLMNPQRSTVWYAAAAACCCADFADFA'
# Warning :
# Keywords : sequence_complexity, calc_sequence_complexity,
# calc_seq_complexity, get_seq_complexity, seg
# Options : _ for debugging.
# # for debugging.
# 'w=' for window size as the first arg
# Returns : Ref. for a scalar digit.
# Argument : ref. of string.
# Category :
# Version : 1.3
#--------------------------------------------------------------
sub get_sequence_complexity{
my ($complexity, @seq,$i, $j, @frag);
my $win=20;
if(ref($_[0]) eq 'ARRAY'){
@seq=@{$_[0]};
}else{
$seq=${$_[0]} || $_[0];
@seq=split(//, $seq);
}
if(defined($_[1])){ $win=${$_[1]} || $_[1]; }
if(@seq <= $win){
my (%seq, @keys);
for($i=0; $i< @seq; $i++){
$seq{$seq[$i]}++;
}
@keys= keys %seq;
$complexity=@keys/@seq;
}else{
my @frag=@{÷_array(\@seq, "s=$win")};
my @complexity=();
for($i=0; $i < @frag; $i++){
my (%seq, @keys);
my @arr=@{$frag[$i]};
for($j=0; $j< @arr; $j++){
$seq{$arr[$j]}++;
}
@keys=keys %seq;
push(@complexity, @keys/$win);
}
$complexity=${&array_average(\@complexity)};
}
return(\$complexity);
}
#"""""""""""""""""< handle_arguments{ head Ver 4.1 >"""""""""""""""""""
my(@A)=&handle_arguments(@_);my($num_opt)=${$A[7]};my($char_opt)=${$A[8]};
my(@hash)=@{$A[0]};my(@file)=@{$A[4]};my(@dir)=@{$A[3]};my(@array)=@{$A[1]};
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my ($connect_gap, @seq_frag, %digitized, $verbose, %hash, $best_block_opt);
my $margin=3;
my $threshold=0.8;
my $min_seqlet_size=25;
$connect_gap=5;
my @vars=keys %vars;
for($i=0; $i< @vars; $i++){
if($vars[$i] eq 'm'){
$margin=$margin_ori=$vars{$vars[$i]};
}elsif($vars[$i] eq 't'){
$threshold=$vars{$vars[$i]};
for($e=0; $e< @keys; $e++){
my $range="$keys[$e]\-$blocks{$keys[$e]}";
my $seq_let_leng=$blocks{$keys[$e]} - $keys[$e] + 1;
if($seq_let_leng < $min_seqlet_size){
next;
}else{
push(@RANGE, $range);
}
}
@seq_frag=&get_seq_fragments(\%hash, @RANGE,
"l=$min_seqlet_size", "$range_in_name");
return(\@seq_frag);
}
#______________________________________________________________
# Title : add_columns
# Usage :
# Function :
# Example :
# Warning : if the attached name is too long(over 12 char),
# it changes to 'Added_upX' while X is a numb.
for($j=$i; $j < ($factor+$i); $j++){
$temp_sum+=$ARRAY[$j]
}
push(@out, int($temp_sum/$factor) );
}
return(\@out);
}
#___________________________________________________________
# Title : get_seq_fragments
# Usage : @seq_frag=&get_seq_fragments(\%msf, @RANGE);
# Function : gets sequence(string) segments with defined
# ranges.
# Example :
# %test=('seq1', '1234AAAAAAAAAAAaaaaa', 'seq2', '1234BBBBBBB');
# @range = ('1-4', '5-8');
#
# %out = %{&get_seq_fragments(\%test, \@range)};
# %out => (seq1_5-8 AAAAA
# seq2_5-8 BBBBB
# seq1_1-4 1234
# seq2_1-4 1234 )
#
# Warning :
# Keywords : get_sequence_fragments,
# Options : _ for debugging.
# # for debugging.
# l= for min seqlet length
# r for adding ranges in the seq names
#
# Returns :
# Argument :
# Category :
# Version : 1.8
#-------------------------------------------------------
sub get_seq_fragments{
#"""""""""""""""""< handle_arguments{ head Ver 4.1 >"""""""""""""""""""
my(@A)=&handle_arguments(@_);my($num_opt)=${$A[7]};my($char_opt)=${$A[8]};
my(@hash)=@{$A[0]};my(@file)=@{$A[4]};my(@dir)=@{$A[3]};my(@array)=@{$A[1]};
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my $no_range_in_name=1;
for($i=0; $i< @vars; $i++){
if($vars[$i] eq 'l'){
$min_seqlet_size=$vars{$vars[$i]};
}
}
if($char_opt=~/v/){ print "\n \$char_opt is $char_opt @char_opt\n"; }
if($char_opt=~/n/){ $no_range_in_name = 1 }
if($char_opt=~/r/){ $no_range_in_name = 0 }
print "\nget_seq_fragments \$no_range_in_name is $no_range_in_name \n";
for($i=0; $i< @hash; $i++){
my (%out_frag, $frag_name, $range_start, $range_end, @out_hash);
my %seqs = %{$hash[$i]};
my @names = keys %seqs;
if(@names==1){
for($j=0; $j < @names; $j++){
my $seq_name = $names[$j];
my $seq = $seqs{$seq_name};
for($k=0; $k< @range; $k++){
my $range = $range[$k];
if($no_range_in_name==1){
$frag_name = "$seq_name";
}else{
$frag_name = "$seq_name\_$range";
}
#if(length($frag_name)>14 ){
# $frag_name ='x'."${j}_${range}";
#}
($range_start, $range_end)=$range=~/(\d+\.?\d*)\-(\d+\.?\d*)/;
my $frag_len = $range_end-$range_start+1;
if($frag_len < $min_seqlet_size){
next;
}
my $fragment = substr($seq, $range_start-1, $frag_len);
$out_frag{$frag_name}=$fragment;
}
}
push(@out_hash, \%out_frag);
}elsif(@names > 1){
for($k=0; $k< @range; $k++){
my %out_frag=();
my $range=$range[$k];
($range_start, $range_end)=$range=~/(\d+\.?\d*)\-(\d+\.?\d*)/;
my $frag_len = $range_end-$range_start+1;
if($frag_len < $min_seqlet_size){
next;
}
for($j=0; $j < @names; $j++){
my $seq_name=$names[$j];
my $seq = $seqs{$seq_name};
if($no_range_in_name==1){
$frag_name = "$seq_name";
}else{
$frag_name = "$seq_name\_$range";
}
#if(length($frag_name)>15 ){
# $frag_name ='x'."${j}_${range}";
#}
if($range_start==0){ $range_start++; } ## This is a bugfix
my $fragment = substr($seq, $range_start-1, $frag_len);
$out_frag{$frag_name}=$fragment;
}
push(@out_hash, \%out_frag);
}
}
}
if(@out_hash > 1){ return(@out_hash)
}elsif(@out_hash==1){ return($out_hash[0]) }
}
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
$ref_target_hash = shift @hash;
%target = %{$ref_target_hash};
@names = keys %target;
$name_target = $names[0];
if($name_target =~/\w+(\d+)/){ $ori_target_seq_len = length($1); }
@target_frag = split(/ +/, $target{$name_target} );
for($i =0; $i< @hash; $i ++){
%db=%{$hash[$i]};
@db_name = keys %db;
for($j=0; $j < @db_name; $j ++){
$name = $db_name[$j];
@db_frag = split( / +/, $db{$name} );
for($k=0; $k < @target_frag; $k ++){
if( ($target_frag[$k]=~/H(\d+)/i)&&($db_frag[$k]=~/H(\d+)/i) ){
$simple_match_output{$name}++;
$leng_diff = abs($1 - $2)/15;
$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}elsif( ($target_frag[$k]=~/E(\d+)/i)&&($db_frag[$k]=~/H(\d+)/i) ){
$simple_match_output{$name}--;
#$leng_diff = abs($1 - $2)/10;
#$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}elsif( ($target_frag[$k]=~/H(\d+)/i)&&($db_frag[$k]=~/E(\d+)/i) ){
$simple_match_output{$name}--;
#$leng_diff = abs($1 - $2)/10;
#$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}elsif( ($target_frag[$k]=~/E(\d+)/i)&&($db_frag[$k]=~/E(\d+)/i) ){
$simple_match_output{$name}++;
$leng_diff = abs($1 - $2)/15;
$simple_match_output{$name} = $simple_match_output{$name}- $leng_diff;
}
}
}
}
return(\%simple_match_output);
}
#_______________________________________________________________
# Title : divide_array
# Usage : &show_array(÷_array(\@input, 6));
# Function : divides any array to the denominator given.
# If you give array of 100 elem, with 5, you will
# get 5 arrays with 20 elem each.
# Example :
# Warning :
# Keywords : split_array_into_pieces, split_array, chop_array,
# fragment_array,
# Options : s= for dividing the array with sub array size
# eg) to get 20 elem length sub arrays from
# a big array
# @ar_ref=@{÷_array(\@array, 's=20')};
# Returns :
# Argument :
# Category :
# Version : 1.4
#-----------------------------------------------------------
sub divide_array{
if((@_ ==1)&&($denominator == 0)){
print "\n Denominator is 0, error, setting to 1\n";
$denominator = 1;
}
if($size_div==1){
while(@array){
push(@final_array_ref, [splice(@array, 0, $size)]);
}
}else{
my $frag_ar_size = int(@array/$denominator);
if($debug eq 1){ print "\n Frag arr size is : $frag_ar_size \n" }
$remaining = @array % $denominator;
if($debug eq 1){ print "\n Remnant elem size is : $remaining \n" }
for($i=0; $i < $denominator; $i++){
if($remaining > 0){
push(@final_array_ref, [splice(@array, 0, ($frag_ar_size+1),)] );
$remaining --;
}elsif(($remaining == 0)&&(@array>0)){
push(@final_array_ref, [splice(@array, 0, ($frag_ar_size),)] );
}
}
}
return(\@final_array_ref);
}
#__________________________________________________________________________
# Title : split_fasta_files
my($division_factor, $base_name);
$division_factor=2;
if($vars{'d'}=~/(\d+)/){ $division_factor=$1 };
for($i=0; $i< @file; $i++){
open(FILE, $file[$i]);
$base_name=${&get_base_names($file[$i])};
$extension=${&get_extension_names($file[$i])};
$total_line_num=@lines=<FILE>;
$splited_frag_size=int(@lines/$division_factor);
for($j=0; $j< $division_factor; $j++){
#$frag_file="$base_name\_s${j}\.split${j}";
$frag_file="$base_name\_s${j}\.$extension";
push(@splited_files, $frag_file);
open(FRAGMENT_FILE, ">$frag_file");
for($k=0; $k<= $splited_frag_size; $k++){
print FRAGMENT_FILE $lines[$k];
}
splice(@lines, 0, $splited_frag_size);
}
print FRAGMENT_FILE @lines; # Writing the very last remaining lines
close (FRAGMENT_FILE);
}
return(\@splited_files);
}
#_______________________________________________________________
# The output will be 'seq_1_half', '1234578'
# 'seq_2_half', '9ABCDEF'
# 'seq_3_half', 'HIJKLMN'
# Warning :
# Keywords : divide_string, split_string, chop_string, divide_sequence
# split_sequence(look at separate split_sequence sub),
# Options :
# $reverse_second_half=S by S -S
# $reverse_first_half =F by F -F
# $reverse_rest =R by R -R ## reversing all except the first
# $reverse_all =A by A -A # reverse all the fragments
# Returns :
# Argument :
# Category :
# Version : 1.3
#-----------------------------------------------------------
sub split_sequence{
my(@string, $frag_str_size, $remaining, @hash, $i, $j,
@seq_names, $denominator, %input_seq, @final_hash_ref,
$frag, $reverse_second_half, $reverse_first_half,
$reverse_rest, $reverse_all);
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# Checking arguments
#________________________________________________
for($i=0; $i< @_; $i++){
if(ref($_[$i]) eq 'HASH'){
push(@hash, $_[$i]);
}elsif(ref($_[$i]) eq 'SCALAR'){
if(${$_[$i]}=~/^(\d+\.?\d*)$/){ # this can handle fraction number!!
for($i=0; $i< @hash; $i++){
my %input_seq=%{$hash[$i]};
@seq_names=keys %input_seq;
my (%out_seq);
for($k=0; $k < @seq_names; $k++){
$seq_name=$seq_names[$k];
@string = split(//, $input_seq{ $seq_name });
$frag_str_size = int(@string/$denominator);
if($debug eq 1){ print "\n Frag str size is : $frag_str_size \n" }
$remaining = @string % $denominator;
for($j=0; $j < $denominator; $j++){
$frag_number=$j+1;
if($remaining > 0){
$frag=join('', splice(@string, 0, ($frag_str_size+1) ) );
$remaining --;
$seq_name_split="$seq_name\_$frag_number";
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~```
# Before final, set the seq types
#_________________________________________________________________
if($frag_number==1 and $reverse_first_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number==2 and $reverse_second_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number >2 and $reverse_rest){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($reverse_all){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}
$out_seq{$seq_name_split}=$frag;
}elsif( $remaining == 0 and @string>0){
$seq_name_split="$seq_name\_$frag_number";
$frag=join('', splice(@string, 0, $frag_str_size,) );
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~```
# Before final, set the seq types
#_________________________________________________________________
if($frag_number==1 and $reverse_first_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number==2 and $reverse_second_half){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($frag_number >2 and $reverse_rest){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}elsif($reverse_all){
$frag=reverse($frag);
$seq_name_split="$seq_name_split\_rv";
}
$out_seq{$seq_name_split}=$frag;
}
}
}
push(@final_hash_ref, \%out_seq);
}
wantarray? \@final_hash_ref : $final_hash_ref[0];
}
# Warning :
# Keywords : divide_string, split_string, chop_string, divide_sequence
# split_sequence(look at separate split_sequence sub),
# Options :
# Returns :
# Argument :
# Category :
# Version : 1.4
#-----------------------------------------------------------
sub divide_string{
my (@array, $i, $j, $denominator, @temp, @string, $frag_str_size,
$remaining, $frag);
for($i=0; $i< @_; $i++){
if(ref($_[$i]) eq 'ARRAY'){
push(@array, @{$_[$i]});
}elsif(ref($_[$i]) eq 'SCALAR'){
if(${$_[$i]} =~/^(\d+)$/){
$denominator = $1;
}else{
push(@array, ${$_[$i]});
}
push(@array, $_[$i]);
}
}
if($denominator == 0){
print "\n Denominator is 0, error, setting to 1\n";
$denominator = 1;
}
for($i=0; $i< @array; $i++){
my @temp;
my @string = split(//, $array[$i]);
$frag_str_size = int(@string/$denominator);
if($debug eq 1){ print "\n Frag str size is : $frag_str_size \n" }
$remaining = @string % $denominator;
for($j=0; $j < $denominator; $j++){
if($remaining > 0){
$frag=join('', splice(@string, 0, ($frag_str_size+1) ) );
push(@temp, $frag);
$remaining --;
}elsif( $remaining == 0 and @string>0){
$frag=join('', splice(@string, 0, $frag_str_size,) );
push(@temp, $frag);
}
}
push(@final_array_ref, \@temp);
}
wantarray? \@final_array_ref : $final_array_ref[0];
}
#____________________________________________________________
# Title : write_html_headbox
# Usage : &write_html_headbox($outfilename, \%entries);
# Example :
# Warning :
# Keywords :
# Options :
# Returns :
# Argument :
# Category :
# Version : 1.1
#-----------------------------------------------------------
sub make_template_from_sec_str{
my(%out, @name, @fragments, %in, $name, $leng, $frag_seq, $name2);
#""""""""""""""""""""""< handle_arguments{ head Ver 1.1 >""""""""""""""""""""""""""""""""
my(@A ) = &handle_arguments( @_ ); my( $num_opt )=${$A[7]};my( $char_opt )=${$A[8]};
my(@hash) =@{$A[0]};my(@file) =@{$A[4]};my(@dir ) =@{$A[3]};my(@array)=@{$A[1]};
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};my(@raw_string)=@{$A[9]};
my($i, $j, $c, $d, $e, $f, $g, $h, $k, $l, $p, $q, $r, $s, $t, $u, $v, $w, $x,$y,$z);
#"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
for($j = 0; $j < @hash ; $j ++){
%in = %{$hash[$j]};
@name = keys %in;
print "\@name is @name\n", if($debug eq 1);
for($t=0; $t < @name; $t++){
$name = $name[$t];
print "\$name is $name\n", if($debug eq 1);
$leng=length($in{$name});
print "$leng\n", if($debug eq 1);
$name2 ="$name"."$leng"; # to attach sequence length
@fragments = split(/_+/, $in{$name});
print "\@fragments is @fragments\n", if($debug eq 1);
for($i = 0; $i < @fragments; $i++){
if($fragments[$i] =~/(\w)\w+/){
$fraglength = length($fragments[$i]);
$frag_seq .= "$1"."$fraglength "; # space is delimiter 'H5 E3 E5 E4'
print "\$frag_seq is $frag_seq\n", if($debug eq 1);
}
$out{$name2}=$frag_seq;
}
}
}
return(\%out);
}
#_______________________________________________________________
# Title : calculate_protein_volume
$scop_class=$2;
$sequence{"${pdb_name}"} = $scop_class;
}
}
}
}
close FILE_1;
}
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
return(\%seq_fragments);
}elsif($char_opt=~/s/i){ # when SIZE return only option is set
@seq_Names=keys %sequence;
for($i=0; $i < @seq_Names; $i++){
$Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
}
return(\%Sizes);
}elsif(@hash >=1){
for($h=0; $h< @hash; $h++){
my %hash=%{$hash[$h]};
my @Keys=keys %hash;
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my (%sequence, %na_sequence, %HASH, @Keys, $seq_found1,
$S_start, $S_end, $seq_found,
$present_seq, @seq_Names, %Sizes, $bare_seq_name, $fasta_seq_idx_file,
%seq_fragments, $get_protein_seq, $reverse_seq, $translation_found);
$get_nucleic_acid_seq='n';
if($char_opt=~/r/){ $reverse_seq='r' }
if(@file<1){
print "\n \@file has less than 1 elem. There is no fileinput for open_embl_files\n";
exit
}
$seq=$na_seq='';
}
close FILE_1;
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`~~~~~~~~~~~~~~~~~~~~~`
# (3) When ranges information is given(via \@range), seq in those ranges are returned
#______________________________________________________________________________________
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
if($reverse_seq){
%seq_fragments=%{&reverse_sequences(\%seq_fragments)};
}
return(\%seq_fragments);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`
# (4) When only size is asked with -s option
#_____________________________________________________________________________
elsif($char_opt=~/s/){ # when SIZE(length of seq) return only option is set
#@seq_Names=keys %sequence;
#for($i=0; $i<@seq_Names; $i++){
# $Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
# }
return(\%Sizes);
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
my (%sequence, %HASH, @Keys, $seq_found1, $S_start, $S_end, $seq_found,
$present_seq, @seq_Names, %Sizes, $bare_seq_name, $fasta_seq_idx_file,
%seq_fragments, $reverse_seq);
if($char_opt=~/r/){ $reverse_seq='r' }
if(@file<1){
print "\n \@file has less than 1 elem. There is no fileinput for open_fasta_files\n";
exit
}
print "\n# \@file are @file\n" if $verbose;
}# end of while
}
close FILE_1;
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`~~~~~~~~~~~~~~~~~~~~~`
# (3) When ranges information is given(via \@range), seq in those ranges are returned
#______________________________________________________________________________________
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
if($reverse_seq){
%seq_fragments=%{&reverse_sequences(\%seq_fragments)};
}
return(\%seq_fragments);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~`
# (4) When only size is asked with -s option
#_____________________________________________________________________________
elsif($char_opt=~/s/){ # when SIZE(length of seq) return only option is set
#@seq_Names=keys %sequence;
#for($i=0; $i<@seq_Names; $i++){
# $Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
# }
return(\%Sizes);
my(@string)=@{$A[2]};my(@num_opt)=@{$A[5]};my(@char_opt)=@{$A[6]};
my(@raw_string)=@{$A[9]};my(%vars)=%{$A[10]};my(@range)=@{$A[11]};
my($i,$j,$c,$d,$e,$f,$g,$h,$k,$l,$m,$n,$o,$p,$q,$r,$s,$t,$u,$v,$w,$x,$y,$z);
if($debug==1){print "\n\t\@hash=\"@hash\"
\@raw_string=\"@raw_string\"\n\t\@array=\"@array\"\n\t\@num_opt=\"@num_opt\"
\@char_opt=\"@char_opt\"\n\t\@file=\"@file\"\n\t\@string=\"@string\"\n" }
#""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
#$| = 1;
my($sort_seq_names, $string, $string_leng, $na,$out_file_name_provided,
$write_pure_seq_only, $show_coil_region, %one_to_three, $each_reliability,
$each_seq_fragment, $each_str_fragment, $end_seq_point, $each_reliability_frag,
$protein_name, $omit_coil_region, $ext, $graphical_rep_of_str,
$show_on_screen_only, $seq_block_size);
my($output_file) ='default_out.pred'; ### when no output file name is given, this is used
$protein_name='pred_temp_name';
my $seq_block_size=80;
if($char_opt=~/r/){ $put_reliability_line='r' } ## putting additional reliability line
if($char_opt=~/c/){ $omit_coil_region='c' } ## putting additional reliability line
if($char_opt=~/g/){ $graphical_rep_of_str='g' } ## putting additional reliability line
if($char_opt=~/s/){ $show_on_screen_only='s' } ## putting additional reliability line
$sec_str=~s/[Cc]/_/g;
}
if($graphical_rep_of_str){
$sec_str=~s/[Hh]/>/g;
$sec_str=~s/[Ee]/=/g;
$sec_str=~s/[Cc]/~/g;
$sec_str=~s/_/ /g;
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~``
# writing down the 80 residue seq fragment
#______________________________________________________
$seq_ruler_line=' . . . . . . . .';
for($i=0; $i< @keys; $i+= $seq_block_size){
$each_seq_ruler_line=substr($seq_ruler_line, 0, @keys-$i-1);
$each_seq_fragment=substr($residue, $i, $seq_block_size);
$each_str_fragment=substr($sec_str, $i, $seq_block_size);
$end_seq_point=$i+length($each_seq_fragment);
unless($show_on_screen_only){
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_seq_ruler_line, ' ');
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n", $i+1, $each_seq_fragment, $end_seq_point);
}else{
printf ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_seq_ruler_line, ' ');
printf ("%5s %-${seq_block_size}s %-5s\n", $i+1, $each_seq_fragment, $end_seq_point);
}
#~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
# if you want to put reliability index line
#____________________________________________
if($put_reliability_line){
$each_reliability_frag=substr($reliability, $i, $seq_block_size);
if(! $show_on_screen_only){
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_str_fragment, ' ');
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_reliability_frag, ' ');
}else{
printf ("%5s %-${seq_block_size}s %-5s\n", ' ', $each_str_fragment, ' ');
printf ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_reliability_frag, ' ');
}
}else{
if(! $show_on_screen_only){
printf PRED_WRITE ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_str_fragment, ' ');
}else{
printf ("%5s %-${seq_block_size}s %-5s\n\n", ' ', $each_str_fragment, ' ');
}
}
}
close PRED_WRITE unless $show_on_screen_only;
splice(@hash, $n, 1); $n--;
}
if( $out_file_name_provided != 1){
print "\n\n# You didnt give out file name, $output_file used\n" unless $show_on_screen_only;
}
if( -s $output_file ){
return(@out_hash_ref_list);
}
}
}
#__________________________________________________________________________
# Title : open_hmmls_files
# Usage : %out=%{&open_hmmls_files(\@file)};
# Function : hmmls matches the full length model to target seq. while, hmmfs
# does for fragments as well.
# Example :
# Keywords :
# Options :
# t=$thresh for bits score threshold
# e=$evalue_thresh for bits score threshold
# r for adding ranges
# m for making MSP file format output
# E=Enguiry_name for specifying enquiry seq name rather than 'HMM', the default
# Returns :
# Argument :
$sequence{$1} = $3;
}
}elsif(/^\>?(\S+)[\t ]+(\D\D+) *$/){ ## plain 'seqname ASDASDFAASSEQUECN'
$sequence{$1} = $2;
}else{ next }
}
}
close FILE_1;
}
if(defined(@range)){
%seq_fragments=%{&get_seq_fragments(\%sequence, \@range)};
return(\%seq_fragments);
}elsif($char_opt=~/s/i){ # when SIZE return only option is set
@seq_Names=keys %sequence;
for($i=0; $i < @seq_Names; $i++){
$Sizes{$seq_Names[$i]}=length($sequence{$seq_Names[$i]});
}
return(\%Sizes);
}elsif(@hash >=1){
for($h=0; $h< @hash; $h++){
my %hash=%{$hash[$h]};
my @Keys=keys %hash;
}
}
}
}
} #<<<<<-- for($i=0; $i< @file; $i++){
if($simple_default_out_flag >=1){ %matched_seq=%simple_default_output; }
if($simple_default_out_flag > 1){ print "\n# You have asked me to open more than one file. I merge output into one hash\n"; }
if(defined(@range) && ($char_opt !~ /n/)){ ## m opt is for matched seq names only
%seq_fragments=%{&get_seq_fragments(\%matched_seq, \@range)};
return(\%seq_fragments);
}elsif($char_opt=~/s/){ # when SIZE return only option is set
@seq_Names=keys %matched_seq;
for($i=0; $i<@seq_Names; $i++){
$Sizes{$seq_Names[$i]}=length($matched_seq{$seq_Names[$i]});
}
return(\%Sizes);
}else{
if($char_opt=~/n/){
return(\%matched_seq_names);
}else{
( run in 2.258 seconds using v1.01-cache-2.11-cpan-364913b4093 )