Add motifs() function to clustering module

This implements an efficient motif counting algorithm.

src/boost-workaround/boost/graph/copy.hpp

+//
+//=======================================================================
+// Copyright 1997-2001 University of Notre Dame.
+// Authors: Jeremy G. Siek, Lie-Quan Lee, Andrew Lumsdaine
+//
+// Distributed under the Boost Software License, Version 1.0. (See
+// accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+//=======================================================================
+//
+
+/*
+  This file implements the following functions:
+
+
+  template <typename VertexListGraph, typename MutableGraph>
+  void copy_graph(const VertexListGraph& g_in, MutableGraph& g_out)
+
+  template <typename VertexListGraph, typename MutableGraph,
+    class P, class T, class R>
+  void copy_graph(const VertexListGraph& g_in, MutableGraph& g_out,
+                  const bgl_named_params<P, T, R>& params)
+
+
+  template <typename IncidenceGraph, typename MutableGraph>
+  typename graph_traits<MutableGraph>::vertex_descriptor
+  copy_component(IncidenceGraph& g_in,
+                 typename graph_traits<IncidenceGraph>::vertex_descriptor src,
+                 MutableGraph& g_out)
+
+  template <typename IncidenceGraph, typename MutableGraph,
+           typename P, typename T, typename R>
+  typename graph_traits<MutableGraph>::vertex_descriptor
+  copy_component(IncidenceGraph& g_in,
+                 typename graph_traits<IncidenceGraph>::vertex_descriptor src,
+                 MutableGraph& g_out,
+                 const bgl_named_params<P, T, R>& params)
+ */
+
+
+#ifndef BOOST_GRAPH_COPY_HPP
+#define BOOST_GRAPH_COPY_HPP
+
+#include <boost/config.hpp>
+#include <vector>
+#include <boost/graph/graph_traits.hpp>
+#include <boost/property_map.hpp>
+#include <boost/graph/named_function_params.hpp>
+#include <boost/graph/breadth_first_search.hpp>
+#include <boost/type_traits/conversion_traits.hpp>
+
+namespace boost {
+
+  namespace detail {
+
+    // Default edge and vertex property copiers
+
+    template <typename Graph1, typename Graph2>
+    struct edge_copier {
+      edge_copier(const Graph1& g1, Graph2& g2)
+        : edge_all_map1(get(edge_all, g1)),
+          edge_all_map2(get(edge_all, g2)) { }
+
+      template <typename Edge1, typename Edge2>
+      void operator()(const Edge1& e1, Edge2& e2) const {
+        put(edge_all_map2, e2, get(edge_all_map1, e1));
+      }
+      typename property_map<Graph1, edge_all_t>::const_type edge_all_map1;
+      mutable typename property_map<Graph2, edge_all_t>::type edge_all_map2;
+    };
+    template <typename Graph1, typename Graph2>
+    inline edge_copier<Graph1,Graph2>
+    make_edge_copier(const Graph1& g1, Graph2& g2)
+    {
+      return edge_copier<Graph1,Graph2>(g1, g2);
+    }
+
+    template <typename Graph1, typename Graph2>
+    struct vertex_copier {
+      vertex_copier(const Graph1& g1, Graph2& g2)
+        : vertex_all_map1(get(vertex_all, g1)),
+          vertex_all_map2(get(vertex_all, g2)) { }
+
+      template <typename Vertex1, typename Vertex2>
+      void operator()(const Vertex1& v1, Vertex2& v2) const {
+        put(vertex_all_map2, v2, get(vertex_all_map1, v1));
+      }
+      typename property_map<Graph1, vertex_all_t>::const_type vertex_all_map1;
+      mutable typename property_map<Graph2, vertex_all_t>::type
+        vertex_all_map2;
+    };
+    template <typename Graph1, typename Graph2>
+    inline vertex_copier<Graph1,Graph2>
+    make_vertex_copier(const Graph1& g1, Graph2& g2)
+    {
+      return vertex_copier<Graph1,Graph2>(g1, g2);
+    }
+
+    // Copy all the vertices and edges of graph g_in into graph g_out.
+    // The copy_vertex and copy_edge function objects control how vertex
+    // and edge properties are copied.
+
+    template <int Version>
+    struct copy_graph_impl { };
+
+    template <> struct copy_graph_impl<0>
+    {
+      template <typename Graph, typename MutableGraph,
+        typename CopyVertex, typename CopyEdge, typename IndexMap,
+        typename Orig2CopyVertexIndexMap>
+      static void apply(const Graph& g_in, MutableGraph& g_out,
+                        CopyVertex copy_vertex, CopyEdge copy_edge,
+                        Orig2CopyVertexIndexMap orig2copy, IndexMap)
+      {
+        typename graph_traits<Graph>::vertex_iterator vi, vi_end;
+        for (tie(vi, vi_end) = vertices(g_in); vi != vi_end; ++vi) {
+          typename graph_traits<MutableGraph>::vertex_descriptor
+            new_v = add_vertex(g_out);
+          put(orig2copy, *vi, new_v);
+          copy_vertex(*vi, new_v);
+        }
+        typename graph_traits<Graph>::edge_iterator ei, ei_end;
+        for (tie(ei, ei_end) = edges(g_in); ei != ei_end; ++ei) {
+          typename graph_traits<MutableGraph>::edge_descriptor new_e;
+          bool inserted;
+          tie(new_e, inserted) = add_edge(get(orig2copy, source(*ei, g_in)),
+                                          get(orig2copy, target(*ei, g_in)),
+                                          g_out);
+          copy_edge(*ei, new_e);
+        }
+      }
+    };
+
+    // for directed graphs
+    template <> struct copy_graph_impl<1>
+    {
+      template <typename Graph, typename MutableGraph,
+        typename CopyVertex, typename CopyEdge, typename IndexMap,
+        typename Orig2CopyVertexIndexMap>
+      static void apply(const Graph& g_in, MutableGraph& g_out,
+                        CopyVertex copy_vertex, CopyEdge copy_edge,
+                        Orig2CopyVertexIndexMap orig2copy, IndexMap)
+      {
+        typename graph_traits<Graph>::vertex_iterator vi, vi_end;
+        for (tie(vi, vi_end) = vertices(g_in); vi != vi_end; ++vi) {
+          typename graph_traits<MutableGraph>::vertex_descriptor
+            new_v = add_vertex(g_out);
+          put(orig2copy, *vi, new_v);
+          copy_vertex(*vi, new_v);
+        }
+        for (tie(vi, vi_end) = vertices(g_in); vi != vi_end; ++vi) {
+          typename graph_traits<Graph>::out_edge_iterator ei, ei_end;
+          for (tie(ei, ei_end) = out_edges(*vi, g_in); ei != ei_end; ++ei) {
+            typename graph_traits<MutableGraph>::edge_descriptor new_e;
+            bool inserted;
+            tie(new_e, inserted) = add_edge(get(orig2copy, source(*ei, g_in)),
+                                            get(orig2copy, target(*ei, g_in)),
+                                            g_out);
+            copy_edge(*ei, new_e);
+          }
+        }
+      }
+    };
+
+    // for undirected graphs
+    template <> struct copy_graph_impl<2>
+    {
+      template <typename Graph, typename MutableGraph,
+        typename CopyVertex, typename CopyEdge, typename IndexMap,
+        typename Orig2CopyVertexIndexMap>
+      static void apply(const Graph& g_in, MutableGraph& g_out,
+                        CopyVertex copy_vertex, CopyEdge copy_edge,
+                        Orig2CopyVertexIndexMap orig2copy,
+                        IndexMap index_map)
+      {
+        typedef color_traits<default_color_type> Color;
+        std::vector<default_color_type>
+          color(num_vertices(g_in), Color::white());
+        typename graph_traits<Graph>::vertex_iterator vi, vi_end;
+        for (tie(vi, vi_end) = vertices(g_in); vi != vi_end; ++vi) {
+          typename graph_traits<MutableGraph>::vertex_descriptor
+            new_v = add_vertex(g_out);
+          put(orig2copy, *vi, new_v);
+          copy_vertex(*vi, new_v);
+        }
+        for (tie(vi, vi_end) = vertices(g_in); vi != vi_end; ++vi) {
+          typename graph_traits<Graph>::out_edge_iterator ei, ei_end;
+          for (tie(ei, ei_end) = out_edges(*vi, g_in); ei != ei_end; ++ei) {
+            typename graph_traits<MutableGraph>::edge_descriptor new_e;
+            bool inserted;
+            if (color[get(index_map, target(*ei, g_in))] == Color::white()) {
+              tie(new_e, inserted) = add_edge(get(orig2copy, source(*ei,g_in)),
+                                              get(orig2copy, target(*ei,g_in)),
+                                              g_out);
+              copy_edge(*ei, new_e);
+            }
+          }
+          color[get(index_map, *vi)] = Color::black();
+        }
+      }
+    };
+
+    template <class Graph>
+    struct choose_graph_copy {
+      typedef typename graph_traits<Graph>::traversal_category Trv;
+      typedef typename graph_traits<Graph>::directed_category Dr;
+      enum { algo =
+             (is_convertible<Trv, vertex_list_graph_tag>::value
+              && is_convertible<Trv, edge_list_graph_tag>::value)
+             ? 0 : is_convertible<Dr, directed_tag>::value ? 1 : 2 };
+      typedef copy_graph_impl<algo> type;
+    };
+
+    //-------------------------------------------------------------------------
+    struct choose_copier_parameter {
+      template <class P, class G1, class G2>
+      struct bind_ {
+        typedef const P& result_type;
+        static result_type apply(const P& p, const G1&, G2&)
+        { return p; }
+      };
+    };
+    struct choose_default_edge_copier {
+      template <class P, class G1, class G2>
+      struct bind_ {
+        typedef edge_copier<G1, G2> result_type;
+        static result_type apply(const P&, const G1& g1, G2& g2) {
+          return result_type(g1, g2);
+        }
+      };
+    };
+    template <class Param>
+    struct choose_edge_copy {
+      typedef choose_copier_parameter type;
+    };
+    template <>
+    struct choose_edge_copy<detail::error_property_not_found> {
+      typedef choose_default_edge_copier type;
+    };
+    template <class Param, class G1, class G2>
+    struct choose_edge_copier_helper {
+      typedef typename choose_edge_copy<Param>::type Selector;
+      typedef typename Selector:: template bind_<Param, G1, G2> Bind;
+      typedef Bind type;
+      typedef typename Bind::result_type result_type;
+    };
+    template <typename Param, typename G1, typename G2>
+    typename detail::choose_edge_copier_helper<Param,G1,G2>::result_type
+    choose_edge_copier(const Param& params, const G1& g_in, G2& g_out)
+    {
+      typedef typename
+        detail::choose_edge_copier_helper<Param,G1,G2>::type Choice;
+      return Choice::apply(params, g_in, g_out);
+    }
+
+
+    struct choose_default_vertex_copier {
+      template <class P, class G1, class G2>
+      struct bind_ {
+        typedef vertex_copier<G1, G2> result_type;
+        static result_type apply(const P&, const G1& g1, G2& g2) {
+          return result_type(g1, g2);
+        }
+      };
+    };
+    template <class Param>
+    struct choose_vertex_copy {
+      typedef choose_copier_parameter type;
+    };
+    template <>
+    struct choose_vertex_copy<detail::error_property_not_found> {
+      typedef choose_default_vertex_copier type;
+    };
+    template <class Param, class G1, class G2>
+    struct choose_vertex_copier_helper {
+      typedef typename choose_vertex_copy<Param>::type Selector;
+      typedef typename Selector:: template bind_<Param, G1, G2> Bind;
+      typedef Bind type;
+      typedef typename Bind::result_type result_type;
+    };
+    template <typename Param, typename G1, typename G2>
+    typename detail::choose_vertex_copier_helper<Param,G1,G2>::result_type
+    choose_vertex_copier(const Param& params, const G1& g_in, G2& g_out)
+    {
+      typedef typename
+        detail::choose_vertex_copier_helper<Param,G1,G2>::type Choice;
+      return Choice::apply(params, g_in, g_out);
+    }
+
+  } // namespace detail
+
+
+  template <typename VertexListGraph, typename MutableGraph>
+  void copy_graph(const VertexListGraph& g_in, MutableGraph& g_out)
+  {
+    if (num_vertices(g_in) == 0)
+      return;
+    typedef typename graph_traits<MutableGraph>::vertex_descriptor vertex_t;
+    std::vector<vertex_t> orig2copy(num_vertices(g_in));
+    typedef typename detail::choose_graph_copy<VertexListGraph>::type
+      copy_impl;
+    copy_impl::apply
+      (g_in, g_out,
+       detail::make_vertex_copier(g_in, g_out),
+       detail::make_edge_copier(g_in, g_out),
+       make_iterator_property_map(orig2copy.begin(),
+                                  get(vertex_index, g_in), orig2copy[0]),
+       get(vertex_index, g_in)
+       );
+  }
+
+  template <typename VertexListGraph, typename MutableGraph,
+    class P, class T, class R>
+  void copy_graph(const VertexListGraph& g_in, MutableGraph& g_out,
+                  const bgl_named_params<P, T, R>& params)
+  {
+    typename std::vector<T>::size_type n;
+      n = is_default_param(get_param(params, orig_to_copy_t()))
+        ? num_vertices(g_in) : 1;
+    if (n == 0)
+      return;
+    std::vector<BOOST_DEDUCED_TYPENAME graph_traits<MutableGraph>::vertex_descriptor>
+      orig2copy(n);
+
+    typedef typename detail::choose_graph_copy<VertexListGraph>::type
+      copy_impl;
+    copy_impl::apply
+      (g_in, g_out,
+       detail::choose_vertex_copier(get_param(params, vertex_copy_t()),
+                                    g_in, g_out),
+       detail::choose_edge_copier(get_param(params, edge_copy_t()),
+                                  g_in, g_out),
+       choose_param(get_param(params, orig_to_copy_t()),
+                    make_iterator_property_map
+                    (orig2copy.begin(),
+                     choose_const_pmap(get_param(params, vertex_index),
+                                 g_in, vertex_index), orig2copy[0])),
+       choose_const_pmap(get_param(params, vertex_index), g_in, vertex_index)
+       );
+  }
+
+  namespace detail {
+
+    template <class NewGraph, class Copy2OrigIndexMap,
+      class CopyVertex, class CopyEdge>
+    struct graph_copy_visitor : public bfs_visitor<>
+    {
+      graph_copy_visitor(NewGraph& graph, Copy2OrigIndexMap c,
+                         CopyVertex cv, CopyEdge ce)
+        : g_out(graph), orig2copy(c), copy_vertex(cv), copy_edge(ce) { }
+
+      template <class Vertex, class Graph>
+      void examine_vertex(Vertex u, const Graph& g_in) const {
+        typename graph_traits<NewGraph>::vertex_descriptor
+          new_u = add_vertex(g_out);
+        put(orig2copy, u, new_u);
+        copy_vertex(u, new_u);
+      }
+
+      template <class Edge, class Graph>
+      void examine_edge(Edge e, const Graph& g_in) const {
+        typename graph_traits<NewGraph>::edge_descriptor new_e;
+        bool inserted;
+        tie(new_e, inserted) = add_edge(get(orig2copy, source(e, g_in)),
+                                        get(orig2copy, target(e, g_in)),
+                                        g_out);
+        copy_edge(e, new_e);
+      }
+    private:
+      NewGraph& g_out;
+      Copy2OrigIndexMap orig2copy;
+      CopyVertex copy_vertex;
+      CopyEdge copy_edge;
+    };
+
+    template <typename Graph, typename MutableGraph,
+              typename CopyVertex, typename CopyEdge,
+              typename Orig2CopyVertexIndexMap, typename Params>
+    typename graph_traits<MutableGraph>::vertex_descriptor
+    copy_component_impl
+      (const Graph& g_in,
+       typename graph_traits<Graph>::vertex_descriptor src,
+       MutableGraph& g_out,
+       CopyVertex copy_vertex, CopyEdge copy_edge,
+       Orig2CopyVertexIndexMap orig2copy,
+       const Params& params)
+    {
+      graph_copy_visitor<MutableGraph, Orig2CopyVertexIndexMap,
+        CopyVertex, CopyEdge> vis(g_out, orig2copy, copy_vertex, copy_edge);
+      breadth_first_search(g_in, src, params.visitor(vis));
+      return get(orig2copy, src);
+    }
+
+  } // namespace detail
+
+
+  // Copy all the vertices and edges of graph g_in that are reachable
+  // from the source vertex into graph g_out. Return the vertex
+  // in g_out that matches the source vertex of g_in.
+  template <typename IncidenceGraph, typename MutableGraph,
+           typename P, typename T, typename R>
+  typename graph_traits<MutableGraph>::vertex_descriptor
+  copy_component(IncidenceGraph& g_in,
+                 typename graph_traits<IncidenceGraph>::vertex_descriptor src,
+                 MutableGraph& g_out,
+                 const bgl_named_params<P, T, R>& params)
+  {
+    typename std::vector<T>::size_type n;
+      n = is_default_param(get_param(params, orig_to_copy_t()))
+        ? num_vertices(g_in) : 1;
+    std::vector<typename graph_traits<IncidenceGraph>::vertex_descriptor>
+      orig2copy(n);
+
+    return detail::copy_component_impl
+      (g_in, src, g_out,
+       detail::choose_vertex_copier(get_param(params, vertex_copy_t()),
+                                    g_in, g_out),
+       detail::choose_edge_copier(get_param(params, edge_copy_t()),
+                                  g_in, g_out),
+       choose_param(get_param(params, orig_to_copy_t()),
+                    make_iterator_property_map
+                    (orig2copy.begin(),
+                     choose_pmap(get_param(params, vertex_index),
+                                 g_in, vertex_index), orig2copy[0])),
+       params
+       );
+  }
+
+  template <typename IncidenceGraph, typename MutableGraph>
+  typename graph_traits<MutableGraph>::vertex_descriptor
+  copy_component(IncidenceGraph& g_in,
+                 typename graph_traits<IncidenceGraph>::vertex_descriptor src,
+                 MutableGraph& g_out)
+  {
+    std::vector<typename graph_traits<IncidenceGraph>::vertex_descriptor>
+      orig2copy(num_vertices(g_in));
+
+    return detail::copy_component_impl
+      (g_in, src, g_out,
+       make_vertex_copier(g_in, g_out),
+       make_edge_copier(g_in, g_out),
+       make_iterator_property_map(orig2copy.begin(),
+                                  get(vertex_index, g_in), orig2copy[0]),
+       bgl_named_params<char,char>('x') // dummy param object
+       );
+  }
+
+} // namespace boost
+
+#endif // BOOST_GRAPH_COPY_HPP

src/graph/clustering/Makefile.am

@@ -16,8 +16,12 @@ libgraph_tool_clustering_la_LDFLAGS = $(MOD_LDFLAGS)
 
 libgraph_tool_clustering_la_SOURCES = \
     graph_clustering.cc \
-    graph_extended_clustering.cc
+    graph_extended_clustering.cc \
+    graph_motifs.cc
+
 
 libgraph_tool_clustering_la_include_HEADERS = \
     graph_clustering.hh \
-    graph_extended_clustering.hh
+    graph_extended_clustering.hh \
+    graph_motifs.hh
+

src/graph/clustering/graph_clustering.cc

@@ -56,10 +56,14 @@ void local_clustering(GraphInterface& g, boost::any prop)
 using namespace boost::python;
 
 void extended_clustering(GraphInterface& g, python::list props);
+void get_motifs(GraphInterface& g, size_t k, python::list subgraph_list,
+                python::list hist, python::list p, bool comp_iso,
+                bool fill_list, size_t seed);
 
 BOOST_PYTHON_MODULE(libgraph_tool_clustering)
 {
     def("global_clustering", &global_clustering);
     def("local_clustering", &local_clustering);
     def("extended_clustering", &extended_clustering);
+    def("get_motifs", &get_motifs);
 }

src/graph/clustering/graph_motifs.cc

+// graph-tool -- a general graph modification and manipulation thingy
+//
+// Copyright (C) 2007  Tiago de Paula Peixoto <tiago@forked.de>
+//
+// This program is free software; you can redistribute it and/or
+// modify it under the terms of the GNU General Public License
+// as published by the Free Software Foundation; either version 3
+// of the License, or (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+#include "graph_filtering.hh"
+
+#include "graph.hh"
+#include "graph_selectors.hh"
+#include "graph_properties.hh"
+
+#include "graph_motifs.hh"
+
+#include <boost/python.hpp>
+#include <boost/ref.hpp>
+
+using namespace std;
+using namespace boost;
+using namespace boost::lambda;
+using namespace graph_tool;
+
+struct null_copy
+{
+    template <class T1, class T2>
+        void operator()(const T1& t1, const T2& t2) const {}
+};
+
+struct append_to_list
+{
+
+    template <class Graph>
+    void operator()(Graph* gp, boost::any& list) const
+    {
+        typedef typename mpl::if_<typename is_directed::apply<Graph>::type,
+                                  d_graph_t,
+                                  u_graph_t>::type graph_t;
+        vector<graph_t>& glist = any_cast<vector<graph_t>&>(list);
+        glist.push_back(graph_t());
+        copy_graph(*gp, glist.back(), vertex_copy(null_copy()).
+                   edge_copy(null_copy()));
+    }
+};
+
+struct retrieve_from_list
+{
+    template <class Graph>
+    void operator()(Graph* gp, boost::any& list, bool& done) const
+    {
+        typedef typename mpl::if_<typename is_directed::apply<Graph>::type,
+                                  d_graph_t,
+                                  u_graph_t>::type graph_t;
+        vector<graph_t>& glist = any_cast<vector<graph_t>&>(list);
+        if (glist.empty())
+        {
+            done = true;
+            return;
+        }
+        copy_graph(glist.back(), *gp, edge_copy(null_copy()));
+        glist.pop_back();
+    }
+};
+
+void get_motifs(GraphInterface& g, size_t k, python::list subgraph_list,
+                python::list hist, python::list p, bool comp_iso,
+                bool fill_list, size_t seed)
+{
+    rng_t rng(static_cast<rng_t::result_type>(seed));
+
+    boost::any list;
+    if (g.GetDirected())
+        list = vector<d_graph_t>();
+    else
+        list = vector<u_graph_t>();
+    try
+    {
+        for (size_t i = 0; i < python::len(subgraph_list); ++i)
+        {
+            GraphInterface& sub =
+                python::extract<GraphInterface&>(subgraph_list[i]);
+            run_action<>()(sub, lambda::bind<void>(append_to_list(), lambda::_1,
+                                                   lambda::var(list)))();
+        }
+    }
+    catch (bad_any_cast&)
+    {
+        throw GraphException("All motif graphs must be either directed or "
+                             "undirected!");
+    }
+
+    vector<size_t> phist;
+    vector<double> plist;
+    double total = 1;
+    for (size_t i = 0; i < python::len(p); ++i)
+    {
+        plist.push_back(python::extract<double>(p[i]));
+        total *= plist[i];
+    }
+
+    boost::any sampler;
+    if (total == 1.0)