graphtreepartition.cc.svn-base

Probabilistic graphical models in matlab.

// This file will be used to implement the tree Partition

#include <iostream>
#include <boost/graph/adjacency_list.hpp>
#include <boost/graph/graph_traits.hpp>
#include <boost/graph/subgraph.hpp>
#include <boost/graph/copy.hpp>
#include <queue>

#include <GraphTreePartition.h>
#include <Potential.h>
#include <RandomVariable.h>
#include <Node.h>
#include <Generic.h>
#include <OutputGraph.h>

using namespace std;
using namespace boost;

// Explicit template instantiations

typedef property< edge_index_t, unsigned int, DiscretePotential *> DiscretePairwiseEdgeProperty;
typedef adjacency_list<vecS, vecS, undirectedS, MessageNode <DiscreteRandomVariable, boost::adjacency_list_traits<vecS, vecS, undirectedS>::vertex_descriptor> *, DiscretePairwiseEdgeProperty > DiscretePairwiseGraph;
typedef subgraph <DiscretePairwiseGraph> SubDiscretePairwiseGraph;

typedef property< edge_index_t, unsigned int> EdgeProperty;
typedef adjacency_list<vecS, vecS, undirectedS, MessageNode <DiscreteRandomVariable, boost::adjacency_list_traits<vecS, vecS, undirectedS>::vertex_descriptor> *, EdgeProperty > DiscreteFactorGraph;
typedef subgraph <DiscreteFactorGraph> SubDiscreteFactorGraph;

typedef property< edge_index_t, unsigned int, ContinuousPotential *> NPEdgeProperty;
typedef adjacency_list<vecS, vecS, undirectedS, GraphicalModelNode *, NPEdgeProperty > NPGraph;
typedef subgraph <NPGraph> SubNPGraph;

// End of the typedefs

template void partition_graph_into_trees< SubDiscretePairwiseGraph > ( SubDiscretePairwiseGraph &);
template void partition_graph_into_trees< SubDiscreteFactorGraph > (SubDiscreteFactorGraph &);

template void partition_mrf_graph_into_chains < subgraph <NPGraph> > (subgraph <NPGraph> &);
template void partition_mrf_graph_into_chains < subgraph <DiscretePairwiseGraph> > (subgraph <DiscretePairwiseGraph> &);

template void partition_simple_chain < subgraph <NPGraph> > (subgraph <NPGraph> &);
template void partition_simple_chain < SubDiscretePairwiseGraph > (SubDiscretePairwiseGraph &);

// End of Explicit template instantiations

//********** Implementation of GraphPartitionQueueOrder ****************//

template <class Graph, class VertexDescriptor>
GraphPartitionQueueOrder<Graph, VertexDescriptor>::GraphPartitionQueueOrder(const Graph & a) : g(a)
{
}

// This is the comparison function. If it returns "true" it means t_1 will be LATER in the queue, ie t_1 > t_2
// If we want it to.

template <class Graph, class VertexDescriptor>
bool GraphPartitionQueueOrder<Graph, VertexDescriptor>::operator()(myTuple t_1, myTuple t_2) const
{
	// 1: choices
	// 2: overall degree
	// 3: step
	// 4: random
	
	if ( tuples::get<1>(t_1)  == tuples::get<1>(t_2) ) // If the number of choices is the same
	{
		
		if ( tuples::get<2>(t_1) == tuples::get<2>(t_2) ) // If the overall degree is the same	
		{
			
			if ( tuples::get<3>(t_1) == tuples::get<3>(t_2) ) // If the steps are the same
			{
				
				//return ( tuples::get<4>(t_1) < tuples::get<4>(t_2) ); // quantity representing the desirability of going into that vertex
				
				//at random
				
				return (uniform_distribution_0_1() < 0.5);
			}
			
			else // we decide by the steps
			{
				return ( tuples::get<3>(t_1) < tuples::get<3>(t_2) ); 
			}
			
		}
		
		else // we decide by the lowest overall degree
		{
			return ( tuples::get<2>(t_1) > tuples::get<2>(t_2) );
		}
		
	}
	
	else // we decide by the number of choices
	{
		return ( tuples::get<1>(t_1)  > tuples::get<1>(t_2) ); 
	}
}

/*
 template <class Graph, class VertexDescriptor>
 bool GraphPartitionQueueOrder<Graph, VertexDescriptor>::operator()(myTuple t_1, myTuple t_2) const
 {
	 // Totally random
	 
	 return (uniform_distribution_0_1() < 0.5);
 }
 */

//********** End of Implementation of GraphPartitionQueueOrder ****************//

template <class Graph, class ColorMap>
unsigned int compute_possible_choices(const Graph & g, ColorMap & color, const typename graph_traits<Graph>::vertex_descriptor v)
{
	unsigned int a(0);
	
	typename graph_traits<Graph>::adjacency_iterator u, u_end;
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	for (tie (u, u_end) = adjacent_vertices(v, g); u != u_end; ++u)
	{
		if (get(color, *u) == Color::white())
		{
			++a;
		}
	}
	
	return a;
}

template <class Graph, class ColorMap>
bool prune_out_degree_0 (Graph & g, ColorMap color)
{
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	typename graph_traits <Graph>::vertex_iterator u, u_end;
	
	for (tie(u, u_end) = vertices(g); u != u_end; ++u)
	{
		if ( out_degree(*u, g) == 0)
		{
			put (color, *u, Color::red() );
			return true;
		}
	}
	
	return false;
}


template <class Graph, class ColorMap, class ParentMap>
void prune_out_degree_1 (Graph & g, ColorMap color, ParentMap parent)
{
	bool prune_again = true;
	
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	typename graph_traits <Graph>::vertex_iterator u, u_end,v, v_end;
	typename graph_traits <Graph>::vertex_descriptor w;
	
	while (prune_again)
	{
		prune_again = false;
		
		for (tie(u, u_end) = vertices(g); u != u_end; ++u)
		{
			
			if ( out_degree(*u, g) == 1)
			{	
				w = *(adjacent_vertices(*u,g).first);
				clear_vertex(*u,g);
				put(color, *u, Color::green());
				prune_again = true;
				
				// Need to record parent, that should do it
				
				for (tie(v, v_end) = vertices(g); v != v_end; ++v)
				{
					if (get(parent,*v) == *u)
					{
						put (parent, *v, w);
					}
				}
			}
		}
	}
}

template <class Graph, class ColorMap, class ParentMap>
void prune_out_degree_2 (Graph & g, ColorMap color, ParentMap parent)
{	
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	typedef typename graph_traits <Graph>::vertex_descriptor VertexDescriptor;
	
	typename graph_traits <Graph>::vertex_iterator u, u_end, x, x_end;
	typename graph_traits <Graph>::adjacency_iterator v, v_end, w;
	
	for (tie(u, u_end) = vertices(g); u != u_end; ++u)
	{
		
		if ( out_degree(*u, g) == 2 ) // we also have to check if there isn't an edge connecting the two neighbours
		{
			tie (v,v_end) = adjacent_vertices (*u,g);
			w = v;
			++w;
			
			VertexDescriptor a = *v;
			VertexDescriptor b = *w;
			
			if (edge(*v, *w, g).second) // this means there is an edge between v and w
			{
				continue;
			}
			
			//cout << "pruned " << *u+1 << endl;
			
			remove_edge(*u, a, g);
			
			remove_edge(*u, b, g);
			
			add_edge(a,b,g);
			
			
			put(color, *u, Color::green());
			
			
			// Need to record parent
			
			for (tie(x, x_end) = vertices(g); x != x_end; ++x)
			{
				if (get(parent,*x) == *u)
				{
					put (parent, *x, *w);
				}
			}
			
		}
		
	}
}



template <class Graph, class ColorMap>
bool check_black_neighbor(Graph & g, ColorMap & color, const typename graph_traits<Graph>::vertex_descriptor v)
{
	typename graph_traits<Graph>::adjacency_iterator u, u_end;
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	for (tie (u, u_end) = adjacent_vertices(v, g); u != u_end; ++u)
	{
		if (get(color, *u) == Color::black())
		{
			return true;
		}
	}
	
	return false;
}

template <class Graph, class ColorMap>
bool check_gray_neighbor(Graph & g, ColorMap & color, const typename graph_traits<Graph>::vertex_descriptor v)
{
	typename graph_traits<Graph>::adjacency_iterator u, u_end;
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	for (tie (u, u_end) = adjacent_vertices(v, g); u != u_end; ++u)
	{
		if (get(color, *u) == Color::gray())
		{
			return true;
		}
	}
	
	return false;
}

template <class Graph, class ColorMap>
void transmit_escape(typename graph_traits<Graph>::vertex_descriptor v, ColorMap & color, vector <unsigned int> & free_edges, vector<list <typename graph_traits<Graph>::vertex_descriptor> > & available_neighbors, Graph & g)
{
	//typename graph_traits<Graph>::adjacency_iterator u, u_end;
	
	typedef typename property_traits<ColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	typename graph_traits<Graph>::vertex_descriptor u = available_neighbors[v].front();
	
	//cout << "We decrease vertex " << g[u]->get_random_variable()->get_index() << " which now has " <<  free_edges [u] -1 << endl;
	
	free_edges [u] = free_edges [u] -1;
	available_neighbors[u].remove(v);
	
	if (free_edges [u] == 1)
	{
		transmit_escape(u, color, free_edges, available_neighbors, g);
	}
	
	
	/*
	 for (tie (u, u_end) = adjacent_vertices(v, g); u != u_end; ++u)
	 {
		 // We are not interested in propagating if neighbors are red or black.
		 
		 if ( (get(color, *u) != Color::red()) && (get(color, *u) != Color::black() ))
		 {
			 cout << "We decrease vertex " << *u+1 << " which now has " <<  free_edges [*u] -1 << endl;
			 
			 free_edges [*u] = free_edges [*u] -1;
			 
			 if (free_edges [*u] == 1)
			 {
				 transmit_escape(*u, color, free_edges, g);
			 }
		 }
	 }
	 */
}

//template <class Graph, class ColorMap, class ParentMap>



typedef iterator_property_map < vector < default_color_type >::iterator, property_map < SubDiscretePairwiseGraph, vertex_index_t>::type > DiscretePairwiseGraphColorMap;
typedef iterator_property_map < vector < graph_traits <SubDiscreteFactorGraph>::vertex_descriptor >::iterator, property_map < SubDiscretePairwiseGraph, vertex_index_t>::type > DiscretePairwiseGraphParentMap;

bool get_tree(DiscretePairwiseGraph & g, DiscretePairwiseGraphColorMap color, DiscretePairwiseGraphParentMap parent)
{
	typedef DiscretePairwiseGraph Graph;
	
	typedef graph_traits <Graph>::vertex_descriptor VertexDescriptor;
	typedef property_traits<DiscretePairwiseGraphColorMap>::value_type ColorValue;
	typedef color_traits<ColorValue> Color;
	
	graph_traits <Graph>::vertex_iterator v, v_end;
	graph_traits <Graph>::adjacency_iterator u, u_end;
	
	priority_queue< GraphPartitionQueueOrder<Graph>::myTuple , vector < GraphPartitionQueueOrder<Graph>::myTuple >, GraphPartitionQueueOrder <Graph> > gray_vertices ( (GraphPartitionQueueOrder <Graph> (g)) );
	
	list <VertexDescriptor> backtrack_list_black;
	list <VertexDescriptor> backtrack_list_gray;
	
	vector <unsigned int> free_edges;

	bool trapped_vertex;
	unsigned int a;
	unsigned int b(0);
	
	VertexDescriptor vertex_with_lowest_degree = VertexDescriptor();

	//cout << "******************************************************" << endl;
	//cout << "We are starting the growing of a new tree." << endl;

	
	tie(v, v_end) = vertices(g);
	
	// If there are no more vertices, we are finished, return false
	
	if (v == v_end)
	{
		return false;
	}
	
	// Initialize all vertices to White, and they should be their own parents
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		put(color,*v, Color::white());
		put(parent,*v, *v);
	}
	
	//cout << "end init_tree()" << endl;
	
	if (prune_out_degree_0(g, color))
	{
		return true;
	}
	
	prune_out_degree_1(g, color, parent);
	
	prune_out_degree_2(g, color, parent);
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		if (get(color,*v) == Color::green())
		{
			continue;
		}
		
		vertex_with_lowest_degree = *v;
		b = out_degree(vertex_with_lowest_degree, g);
		break;
	}
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		free_edges.push_back(out_degree(*v, g));
		
		if (get(color,*v) == Color::green())
		{
			continue;
		}
		
		a = out_degree(*v, g);
		
		if (a < b)
		{
			b = a;
			vertex_with_lowest_degree = *v;
		}
	}
	
	
	vector<list <VertexDescriptor> > available_neighbors (num_vertices(g), list <VertexDescriptor> () );
	
	for (tie(v, v_end) = vertices(g); v != v_end; ++v)
	{
		for (tie (u, u_end) = adjacent_vertices(*v, g); u != u_end; ++u)
		{
			available_neighbors[*v].push_back(*u);
		}
	}
	
	// Obtain vertex with lowest degree
	
	gray_vertices.push( make_tuple (vertex_with_lowest_degree, compute_possible_choices(g, color, vertex_with_lowest_degree), out_degree(vertex_with_lowest_degree, g), 0, 0) );
	
	put(color, vertex_with_lowest_degree, Color::gray());
	
	unsigned step = 0;