#!/usr/bin/perl

#use strict;
#use warnings;
use diagnostics;

# Program: GO_ancestors.pl
#
# Original Author: Mary Dolan
#
# Purpose: Given gene ontology in OBO format
#    generate mapping from GO ids to all ancestor ids
#
# Usage: perl GO_ancestors.pl
#

######### read in OBO terms
	&parse_obo(%nodes);

######### nodes (GO ids) to ancestors
	foreach $query_node (keys %nodes) {
	@ancestor_nodes = ();
	&find_ancestors($query_node);
	}

######### map OBO ids to OBO ancestors
open (GO_ANCS, ">GO_ancestors.txt");
#print header
    print GO_ANCS "GO_id\tGO_ancestor\n";

	foreach $query_node (keys %nodes) {
	@ancestor_nodes = @ { $nodes{$query_node}->{"ancestors"} };

	foreach $ancestor (@ancestor_nodes) {
		print GO_ANCS "$query_node\t$ancestor\n";
		}
	}


close (GO_ANCS) || die "GO_ANCS file was not closed\n";	
exit (0);


############
sub parse_obo()
{
    %nodes = ();

    open (OBO, "gene_ontology.obo") ||print "ontology open error";
    while (<OBO>) {

    if (/\[Term\]/) {
    #new term
	&get_node_info($node_id, $node_name, @parent_nodes);
	$term = {
		"name" => $node_name,
		"parents" => [ @parent_nodes ],
		"ancestors" => [  ],
		};

	$nodes{$node_id} = $term;
	} #end term test
     } #end read loop

     close (OBO);
}

############
sub get_node_info()
{
#clear info
    $node_id = "";
    $node_name = "";
    @parent_nodes = ();

    while (<OBO>) {
    chomp;
    $term_line = $_;

    #node_id
	if ($term_line =~ /^id: ([A-Z]+:\d+)/) {
	    $node_id = $1;
	    }

    #node name
	if ($term_line =~ /^name: (.*)/) {
	    $node_name = $1;
	    }

    #is_a
	if ($term_line =~ /^is_a: (\w+:\d+)/) {
	    $parent_id = $1;
	    $parent_type = "is_a";
	    $parent = {
		      "id" => $parent_id,
		      "type" => $parent_type,
		      };
	    push @parent_nodes, $parent;
	    }

    #other relationship
	if ($term_line =~ /^relationship: (\w+) (\w+:\d+)/) {
	    $parent_id = $2;
	    $parent_type = $1;
	    $parent = {
		      "id" => $parent_id,
		      "type" => $parent_type,
		      };
	    push @parent_nodes, $parent;
	    }

    #end of Term
	if ($term_line eq "") {
	    return();
	    }

    } #end read node
} #end get_node_info routine


#############
sub find_ancestors()
{
	#find parents/ancestors
	&find_parents($query_node);

		# select unique objects
		%seen = ();
    		$seen{$query_node}++;
		foreach $item (@ancestor_nodes) {
    			$seen{$item}++;
			}
		@ancestor_nodes = keys %seen;
	$nodes{$query_node}->{"ancestors"} = [ @ancestor_nodes ];
}

#############
sub find_parents()
{
    local $node_id = shift;
    local @parent_nodes = @ { $nodes{$node_id}->{"parents"} };

    #no more parents (stopping condition)
    if (!@parent_nodes) {
	return();
	}
    #
    foreach $parent (@parent_nodes) {
	$parent_id = $parent->{"id"};
	push @ancestor_nodes, $parent_id;
	&find_parents($parent_id);
	}
} #end find_parents routine


__END__

