#!/usr/bin/perl

#########################################################
#  Written by Justin Elser 4/29/2009                    #
#                					#
#	This program takes the genes known to interact	#
#	in the long day flowering process in 		#
#	Arabidopsis and predicts the similar short day	#
#	interaction in rice				#
#							#
#	Uses the orthologs first, then paralogs		#
#							#
#########################################################


use strict;
use warnings;

# print usage
if ($#ARGV !=4) {
	print "usage: flowering_interaction.pl input_file interaction_network_file ortholog_file paralog_file output_file_prefix\n";
	exit;
}

my $in_file = $ARGV[0];
my $inter_file = $ARGV[1];
my $ortho_file = $ARGV[2];
my $para_file = $ARGV[3];
my $out_file = $ARGV[4] . ".sif";
my $out_file_inter = $ARGV[4] . "_inter.sif";
my $out_file_ortho = $ARGV[4] . "_ortho.sif";
my $out_file_para = $ARGV[4] . "_para.sif";

my %in_hash;


open (input_file, "$in_file") || die "Error: file '$in_file' can not be opened \n";
while(<input_file>) {
	my $entry = $_;
	# strip off newline characters
	$entry =~ s/\r//g;
	$entry =~ s/\n//g;
	# split the columns into 2 separate variables (using either multiple congruent spaces or a tab)
	my ($gene_name, $gene_id) = split /\s{2,}|\t/, $entry;
	$in_hash{$gene_id} = $gene_name;
}
close (input_file);

open (interaction_file, "$inter_file") || die "Error: file '$inter_file' can not be opened\n";
my %inter_hash;
while(<interaction_file>) {
	my $entry = $_;
	# strip off newline characters
	$entry =~ s/\r//g;
	$entry =~ s/\n//g;
	# split the columns into 3 separate variables
	my ($gene_1, $inter_type, $gene_2) = split ("\t", $entry);
	# change gene ids to all caps (for matching)
	$gene_1 =~ tr/a-z/A-Z/;
	$gene_2 =~ tr/a-z/A-Z/;
	
	# because it doesn't matter which way the interaction goes, need to do the interaction in both directions
	#      ie. as gene_1 as the key then as gene_2 as the key
	if(defined($inter_hash{$gene_1})){
		# if gene already in hash, just add to the list
		my $oldhash = $inter_hash{$gene_1};
		$oldhash->{'type'} = "$oldhash->{'type'}\t$inter_type";
		$oldhash->{'inter_gene'} = "$oldhash->{'inter_gene'}\t$gene_2";
		$inter_hash{$gene_1} = $oldhash;
		# if new gene, create hash element
	}else{
		my %hash;
		$hash{'map_gene'} = $gene_1;
		$hash{'type'} = $inter_type;
		$hash{'inter_gene'} = $gene_2;
		$inter_hash{$gene_1} = \%hash;
	}
	
	if(defined($inter_hash{$gene_2})){
		# if new gene in the hash, just add to the list
		my $oldhash = $inter_hash{$gene_2};
		$oldhash->{'type'} = "$oldhash->{'type'}\t$inter_type";
		$oldhash->{'inter_gene'} = "$oldhash->{'inter_gene'}\t$gene_1";
		$inter_hash{$gene_2} = $oldhash;
		# if new gene, create hash element
	}else{
		my %hash;
		$hash{'map_gene'} = $gene_2;
		$hash{'type'} = $inter_type;
		$hash{'inter_gene'} = $gene_1;
		$inter_hash{$gene_2} = \%hash;
	}
}
close(interaction_file);



my $min_ident=0;

my %ortho_hash;
#my @ortho_array;

open (ortholog_file, "$ortho_file") || die "Error: file '$ortho_file' can not be opened\n";
while(<ortholog_file>) {
	my $entry = $_;
	# strip off newline characters
	$entry =~ s/\r//g;
	$entry =~ s/\n//g;
	# split the columns into 3 separate variables
	my ($gene_id, $ortho_id, $ident_id) = split("\t", $entry);
	# change gene ids to all caps
	$gene_id =~ tr/a-z/A-Z/;
	$ortho_id =~ tr/a-z/A-Z/;
	# strip off the TAIR-G suffix on ids
	$gene_id =~ s/\-TAIR-G//g;
	
	# check to make sure the %identity is above the minimum defined above
	if($ident_id =~ /\d/ && $ident_id>=$min_ident){
		# if gene_id already found, add the new ortholog to hash
		if(defined($ortho_hash{$gene_id})){
			my $oldhash = $ortho_hash{$gene_id};
			$oldhash->{'ortho'} = "$oldhash->{'ortho'}\t$ortho_id";;
			$ortho_hash{$gene_id} = $oldhash;
		# if new gene_id, create hash element
		}else{
			my %hash;
			$hash{'ortho'} = $ortho_id;
			$hash{'gene'} = $gene_id;
			$ortho_hash{$gene_id} = \%hash;
			#push @ortho_array, $gene_id;
		}
	}
}
close (ortholog_file);

my %para_hash;
open(paralog_file, "$para_file") || die "Error: file '$para_file' can not be opened\n";
	
	my $min_ident_para = 0; # set the threshold for paralog identities
	while(<paralog_file>) {
		my $entry = $_;
		# strip off newline characters
		$entry =~ s/\r//g;
		$entry =~ s/\n//g;
		# split the columns into 3 separate variables
		my ($gene_1, $gene_2, $ident) = split("\t", $entry);
		# skip if $ident is empty or Nan
		next if(!defined($ident));
		next if($ident =~ /D/);
		next if($ident eq "");
		next if($ident =~ /\%/);
		# change gene ids to all caps
		$gene_1 =~ tr/a-z/A-Z/;
		$gene_2 =~ tr/a-z/A-Z/;
		
		# put the genes in ASCII order to help remove dupes
		if($gene_1 gt $gene_2) {
			my $temp = $gene_1;
			$gene_1 = $gene_2;
			$gene_2 = $temp;
		}
		
		# only use paralogs with high confidence values
		if($ident>=$min_ident_para) {
			# most genes will show up many times, so push values for hash key
			if(defined($para_hash{$gene_1})) {
				my $oldhash = $para_hash{$gene_1};
				$oldhash->{'gene2'} = "$oldhash->{'gene2'}\t$gene_2";;
				$para_hash{$gene_1} = $oldhash;
				#print "$oldhash\n";
			}else{
				my %hash;
				$hash{'gene1'} = $gene_1;
				$hash{'gene2'} = $gene_2;
				$para_hash{$gene_1} = \%hash;
			}
		}
		
	}

	close(paralog_file);

# open output files for writing
open(output_file, ">$out_file");
open(output_file_inter, ">$out_file_inter");
open(output_file_ortho, ">$out_file_ortho");
open(output_file_para, ">$out_file_para");



# initialize the counters
my $Ath_inter_count = 0;
my $rice_ortho_count = 0;
my $rice_para_count = 0;

foreach my $key (keys %in_hash) {
	my $name = $in_hash{$key};
	
	# find all interactions with the genes in %in_hash
	if(defined($inter_hash{$key})) {
		
		my @inter_array = split ("\t", $inter_hash{$key}->{'inter_gene'});
		my @inter_type_array = split ("\t", $inter_hash{$key}->{'type'});
		my $inter_index = @inter_type_array;
		for (my $i=0; $i<$inter_index; $i++) {
			my $in_gene_1 = $key;
			my $in_gene_2 = $inter_array[$i];
			# add one to the number of Ath interactions counted
			$Ath_inter_count++;
			# output the interactions to file
			print output_file_inter "$in_gene_1\t$inter_type_array[$i]\t$in_gene_2\n";
			
			# make the genes be listed in ASCII order to make sorting out dupes easier
			if ($in_gene_1 gt $in_gene_2) { # note that this should not be used if order matters for interaction type
				my $temp = $in_gene_1;
				$in_gene_1 = $in_gene_2;
				$in_gene_2 = $temp;
			}
			
			# next step is to find all the orthologs for both genes in the interaction
			if(exists($ortho_hash{$in_gene_2})) {
				my @ortho_gene_2_array = split ("\t", $ortho_hash{$in_gene_2}->{'ortho'});
				foreach my $ortho_gene_2 (@ortho_gene_2_array) {
					# output orthologs to file
					print output_file_ortho "$in_gene_2\tortho\t$ortho_gene_2\n";
					
					if(exists($ortho_hash{$in_gene_1})) {
						my @ortho_gene_1_array = split ("\t", $ortho_hash{$in_gene_1}->{'ortho'});
						foreach my $ortho_gene_1 (@ortho_gene_1_array) {
							#print output_file "$ortho_gene_1\t$inter_type_array[$i]\t$ortho_gene_2\n";
							print output_file "$ortho_gene_1\tortho\t$ortho_gene_2\n";
							# add one to the number of orthologs counted
							$rice_ortho_count++;
							# output the orthologs to file
							print output_file_ortho "$in_gene_1\tortho\t$ortho_gene_1\n";
							
							# now find all paralogs to the orthologs found
							if(exists($para_hash{$ortho_gene_1})) {
								my @para_gene_1_array = split ("\t", $para_hash{$ortho_gene_1}->{'gene2'});
								foreach my $para_gene_1 (@para_gene_1_array) {
									# output the paralogs to file
									print output_file_para "$ortho_gene_1\tpara\t$para_gene_1\n";
									
									if(exists($para_hash{$ortho_gene_2})) {
										my @para_gene_2_array = split ("\t", $para_hash{$ortho_gene_2}->{'gene2'});
										foreach my $para_gene_2 (@para_gene_2_array) {
											# add one to the number of paralogs counted
											$rice_para_count++;
											#print output_file "$para_gene_1\t$inter_type_array[$i]\t$para_gene_2\n";
											print output_file "$para_gene_1\tpara\t$para_gene_2\n";
											# output the paralogs to file
											print output_file_para "$ortho_gene_2\tpara\t$para_gene_2\n";

										}
									}
								}
							}
						}
					}
				}
			}
			
		}
	}
}

close(output_file);
print "Ath_inter_count = $Ath_inter_count\nrice_ortho_count = $rice_ortho_count\nrice_para_count = $rice_para_count\n";
print "Now removing duplicates\n";
# sort the file and get rid of duplicates
#   Note: this requires replacing the file, hence the mv command
system "sort $out_file | uniq > $out_file.tmp; mv $out_file.tmp $out_file";


close(output_file_inter);
close(output_file_ortho);
close(output_file_para);




