createCommonString3dBlastDataset.pl

#!/usr/bin/perl -w
# Author: Abdullah Kahraman
# Date: 01.01.2000

###############################################################################
###############################################################################
### bla bla bla bla bla bla bla bla bla bla bla bla bla bla bla bla         ###
###############################################################################
###############################################################################

use strict;
use warnings;
use Getopt::Long;

my (
    # variable for parameters which are read in from commandline
    $help,
    $label,
   );

##############################################################################
### read all needed parameters from commandline ##############################

&GetOptions(
    "help!"      => \$help,     # print this help
    "label!"     => \$label,    # label entries that have been added by other dataset
) or die "\nTry \"$0 -h\" for a complete list of options\n\n";

##############################################################################

# help
if ($help) {printHelp(); exit}

##############################################################################
### SETTINGS #################################################################
##############################################################################

##############################################################################
### SUBROUTINES ##############################################################
############################################################################## 

###############################################################################
sub printHelp {
###############################################################################
    # prints a help about the using and parameters of this scripts 
    # (execute if user types commandline parameter -h)
    # param:  no paramaters
    # return: no return value

    my (
	$usage,
	$sourceCode,
	@rows,
	$row,
	$option,
	$scriptInfo,
	$example,
       );

    $usage = "$0\n";


    print "\nUsage: " .  $usage . "\n";

    print "Valid options are:\n\n";
    open(MYSELF, "$0") or
      die "Cannot read source code file $0: $!\n";
    $sourceCode .= join "", <MYSELF>;
    close MYSELF;
    $sourceCode =~ s/^.+?\&GetOptions\(\n//s;
    $sourceCode =~ s/\n\).+$//s;
    @rows = split /\n/, $sourceCode;
    foreach $row (@rows){
        $option = $row;
	$option =~ s/\s+\"//g;
	$option =~ s/\"\s.+\#/\t\#/g;
	$option =~ s/=./\t<value> [required]/;
	$option =~ s/:./\t<value> [optional]/;
	$option =~ s/!/\t<non value> [optional]/;

	$row =~ s/^.*//;
	print "\t";
	printf("%-1s%-30s%-30s\n", "-",$option,$row);

    } # end of foreach $row (@rows)
    print "\n";
    print "Options may be abreviated, e.g. -h for --help\n\n";

     $example  = "$0";
}

##############################################################################
### END OF SUBROUTINES########################################################
############################################################################## 


############
### MAIN ###
############

my $blastDsFile = "interactions_hensemble78_blast_result_predHumanFromNonHuman_testVsAllHumanPDB_sumMinMax_enspPfamBs_sel_uniq_nonHomolog_adjNegDs_sel.txt.gz";
my $stringDsFile = "string3d_human_all_consurf_predHumanFromNonHuman_testVsAllHumanPDB_full_sumMinMax_enspPfamBs_sel_uniq_nonHomolog_adjNegDs_sel_psiblastWater_string3d.txt.gz";
my $blastAllDsFile = "interactions_hensemble78_blast_result_predHumanFromNonHuman_testVsAllHumanPDB_sumMinMax_enspPfamBs_sel_uniq_nonHomolog.txt.gz";
my $stringAllDsFile = "string3d_human_all_consurf_predHumanFromNonHuman_testVsAllHumanPDB_full_sumMinMax_enspPfamBs_sel_uniq_nonHomolog.txt.gz";

my $blastNewDsFile = "interactions_hensemble78_blast_result_predHumanFromNonHuman_testVsAllHumanPDB_sumMinMax_enspPfamBs_sel_uniq_nonHomolog_adjNegDs_sel_filled.txt.gz";
my $stringNewDsFile = "string3d_human_all_consurf_predHumanFromNonHuman_testVsAllHumanPDB_full_sumMinMax_enspPfamBs_sel_uniq_nonHomolog_adjNegDs_sel_psiblastWater_string3d_filled.txt.gz";

my %blastDs = ();
my %blastDsY = ();
open(F1, "zcat $blastDsFile |") or die "\nERROR: Failed to unzip and read blast dataset file \"$blastDsFile\": $!\n\n";
while(my $l = <F1>) {
    chomp($l);
    next if($l =~ /^#/);
    my ($id1, $id2, @tmp) = split(/\t/, $l);

    if($id1 ge $id2) {
	my @a = split(/\t/, $l);
	$id1 = $a[1];
	$id2 = $a[0];
	$l = "$id1\t$id2";
	for(my $i=2; $i<=$#a; $i++) {
	    $l .= "\t$a[$i]";
	}
    }

    $blastDs{$id1}->{$id2} = $l;
    $blastDsY{$id1}->{$id2} = $tmp[-1];
    if(exists $blastDs{$id2}->{$id1}) {
	print STDERR "WARNING: reverse interaction found in blast dataset between $id2 and $id1\n";
    }
}
close(F1);

my %stringDs = ();
my %stringDsY = ();
open(F2, "zcat $stringDsFile |") or die "\nERROR: Failed to unzip and read string dataset file \"$stringDsFile\": $!\n\n";
while(my $l = <F2>) {
    chomp($l);
    next if($l =~ /^#/);
    my ($id1, $id2, @tmp) = split(/\t/, $l);

    if($id1 ge $id2) {
	my @a = split(/\t/, $l);
	$id1 = $a[1];
	$id2 = $a[0];
	$l = "$id1\t$id2";
	for(my $i=2; $i<=$#a; $i++) {
	    $l .= "\t$a[$i]";
	}
    }

    $stringDs{$id1}->{$id2} = $l;
    $stringDsY{$id1}->{$id2} = $tmp[-1];
    if(exists $stringDs{$id2}->{$id1}) {
	print STDERR "WARNING: reverse interaction found in string dataset between $id2 and $id1\n";
    }
}
close(F2);

my %blastAllDs = ();
open(F3, "zcat $blastAllDsFile |") or die "\nERROR: Failed to unzip and read all blast dataset file \"$blastAllDsFile\": $!\n\n";
while(my $l = <F3>) {
    chomp($l);
    my @a = split(/\t/, $l);

    my $id1 = $a[8];
    my $id2 = $a[9];

    if($id1 ge $id2) {
	$id1 = $a[9];
	$id2 = $a[8];
    }

    $l = "$id1\t$id2";
    for(my $i=2; $i<=7; $i++) {
	$l .= "\t$a[$i]";
    }
    $l .= "\t$a[-3]\n";

    $blastAllDs{$id1}->{$id2} = $l;
    if(exists $blastAllDs{$id2}->{$id1}) {
	print STDERR "WARNING: reverse interaction found in all blast dataset between $id2 and $id1\n";
    }
    $blastAllDs{$id2}->{$id1} = $l;
}
close(F3);

my %stringAllDs = ();
open(F4, "zcat $stringAllDsFile |") or die "\nERROR: Failed to unzip and read all string dataset file \"$stringAllDsFile\": $!\n\n";
while(my $l = <F4>) {
    chomp($l);
    my @a = split(/\t/, $l);

    my $id1 = $a[-5];
    my $id2 = $a[-4];

    if($id1 ge $id2) {
	$id1 = $a[-4];
	$id2 = $a[-5];
    }

    $l = "$id1\t$id2";
    for(my $i=2; $i<=45; $i++) {
	$l .= "\t$a[$i]";
    }
    $l .= "\t$a[-3]\n";

    $stringAllDs{$id1}->{$id2} = $l;
    if(exists $stringAllDs{$id2}->{$id1}) {
	print STDERR "WARNING: reverse interaction found in all string dataset between $id2 and $id1\n";
    }
    $stringAllDs{$id2}->{$id1} = $l;
}
close(F4);

##########################################
# now create comparable list

my %string = ();
foreach my $id1 (sort keys %stringDs) {
    foreach my $id2 (sort keys %{$stringDs{$id1}}) {
	$string{$id1}->{$id2} = $stringDs{$id1}->{$id2};
    }
}
foreach my $id1 (sort keys %blastDs) {
    foreach my $id2 (sort keys %{$blastDs{$id1}}) {
	if(!exists $stringDs{$id1}->{$id2} and !exists $stringDs{$id2}->{$id1}) {
	    if(exists $stringAllDs{$id1}->{$id2}) {
		$string{$id1}->{$id2} = ($label ? "-" : "").$stringAllDs{$id1}->{$id2};
	    } else {
		$string{$id1}->{$id2} = ($label ? "-" : "")."$id1\t$id2\t0\t0\t99\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t999\t999\t0\t0\t0\t999\t0\t0\t99\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t0\t999\t999\t0\t0\t0\t999\t".$blastDsY{$id1}->{$id2};
	    }
	}
    }
}

my %blast = ();
foreach my $id1 (sort keys %blastDs) {
    foreach my $id2 (sort keys %{$blastDs{$id1}}) {
	$blast{$id1}->{$id2} = $blastDs{$id1}->{$id2};
    }
}
foreach my $id1 (sort keys %stringDs) {
    foreach my $id2 (sort keys %{$stringDs{$id1}}) {
	if(!exists $blastDs{$id1}->{$id2} and !exists $blastDs{$id2}->{$id1}) {
	    if(exists $blastAllDs{$id1}->{$id2}) {
		$blast{$id1}->{$id2} = ($label ? "-" : "").$blastAllDs{$id1}->{$id2};
		print O2 
	    } else {
		$blast{$id1}->{$id2} = ($label ? "-" : "")."$id1\t$id2\t0\t0\t0\t0\t0\t0\t".$stringDsY{$id1}->{$id2}."\n";
	    }
	}
    }
}

#####################################
# now output comparable list to files

open(O1, "| gzip > $stringNewDsFile") or die "\nERROR: Failed to zip and new string dataset file \"$stringNewDsFile\": $!\n\n";
open(O2, "| gzip > $blastNewDsFile") or die "\nERROR: Failed to zip and new blast dataset file \"$blastNewDsFile\": $!\n\n";

print O1 "#id1\tid2\tbitScoreSum\tprobisZscoreSum\tprobisRmsdSum\tprobisInsideAlignSum\tprobisInOutRatioSum\thhsearchSimSum\thhsearchIdentSum\thhsearchProbSum\thhsearchInsideAlignSum\thhsearchInOutRatioSum\tconsurfTemplateAvgSum\tconsurfTargetAvgSum\tconsurfTemplateSumSum\tconsurfTargetSumSum\tprobisBlosumSum\tprobisBsSizeSum\tprobisOutsideAlignSum\thhsearchEvalueSum\thhsearchScoreSum\thhsearchAlignColsSum\thhsearchSumProbsSum\thhsearchOutsideAlignSum\tbitScoreMin\tprobisZscoreMin\tprobisRmsdMax\tprobisInsideAlignMin\tprobisInOutRatioMin\thhsearchSimMin\thhsearchIdentMin\thhsearchProbMin\thhsearchInsideAlignMin\thhsearchInOutRatioMin\tconsurfTemplateAvgMin\tconsurfTargetAvgMin\tconsurfTemplateSumMin\tconsurfTargetSumMin\tprobisBlosumMin\tprobisBsSizeMin\tprobisOutsideAlignMax\thhsearchEvalueMax\thhsearchScoreMin\thhsearchAlignColsMin\thhsearchSumProbsMin\thhsearchOutsideAlignMax\ty\n";
print O2 "#id1\tid2\talignLengthSum\tidentSum\tbitScoreSum\talignLenthMin\tidentMin\tbitScoreMin\ty\n";
foreach my $id1 (sort keys %string) {
    foreach my $id2 (sort keys %{$string{$id1}}) {
	print O1 $string{$id1}->{$id2}."\n";
	print O2 $blast{$id1}->{$id2}."\n";
    }
}

close(O1);
close(O2);