* ComponentSizeHistogram() becomes LabelComponents()

* fixed saving of undirected filtered graphs in graph_io.cc
* fixed correlated community detection in graph_community.cc


git-svn-id: https://svn.forked.de/graph-tool/trunk@53 d4600afd-f417-0410-95de-beed9576f240

src/graph-tool

@@ -102,7 +102,7 @@ statistics.add_option("--average-harmonic-distance", action="callback", callback
 statistics.add_option("--sampled-distance-histogram", action="callback", callback=push_option, type="string", metavar="[WEIGHT|]SAMPLES|SEED|FILE", help="get the sampled distance histogram")
 statistics.add_option("--average-sampled-distance", action="callback", callback=push_option, type="string", metavar="[WEIGHT|]SAMPLES|SEED|FILE", help="get the average sampled distance")
 statistics.add_option("--average-sampled-harmonic-distance", action="callback", callback=push_option, type="string", metavar="[WEIGHT|]SAMPLES|SEED|FILE", help="get the average sampled harmonic distance")
-statistics.add_option("--component-size-histogram", action="callback", callback=push_option, type="string", metavar="FILE", help="get the component size histogram")
+statistics.add_option("--label-components", action="callback", callback=push_option, type="string", metavar="PROPERTY", help="label components to PROPERTY")
 statistics.add_option("--average-vertex-property", action="callback", callback=push_option, type="string", metavar="PROPERTY|FILE", help="get the average of the vertex property")
 statistics.add_option("--average-edge-property", action="callback", callback=push_option, type="string", metavar="PROPERTY|FILE", help="get the average of the edge property")
 statistics.add_option("--reciprocity", action="callback", callback=push_option, type="string", metavar="FILE", help="get the edge reciprocity")
@@ -406,10 +406,10 @@ def parse_option(opt, just_file=False):
         hist = graph.GetSampledDistanceHistogram(weight, samples, seed)
         avg, err = hist_mean(dict((1.0/k,v) for k,v in hist.iteritems()))
         return ("%f\t%f" % (1.0/avg,1.0/err) , file_name)
-    elif opt.name == "component-size-histogram":
+    elif opt.name == "label-components":
         if just_file:
-            return opt.value
-        return (graph.GetComponentSizeHistogram(), opt.value)
+            return
+        graph.LabelComponents(opt.value)
     elif opt.name == "vertex-histogram":
         values = parse_values(opt.value)
         if len(values) != 2:

src/graph/graph.cc

@@ -278,60 +278,46 @@ GraphInterface::hist_t GraphInterface::GetEdgeHistogram(string property) const
 }
 
 //==============================================================================
-// GetComponentSizeHistogram()
+// LabelComponents(property)
 //==============================================================================
 
-struct strong_component
+struct label_components
 {
-    template <class Graph, class ComponentMap>
-    void operator()(const Graph &g, ComponentMap comp_map) const
+    template <class Graph, class CompMap>
+    void operator()(const Graph &g, CompMap comp_map) const
+    {   
+	get_components(g, comp_map, typename is_convertible<typename graph_traits<Graph>::directed_category, directed_tag>::type());
+    }
+
+    template <class Graph, class CompMap>
+    void get_components(const Graph& g, CompMap comp_map, boost::true_type is_directed) const
     {
 	strong_components(g, comp_map);
     }
-};
 
-struct connected_component
-{
-    template <class Graph, class ComponentMap>
-    void operator()(const Graph &g, ComponentMap comp_map) const
+    template <class Graph, class CompMap>
+    void get_components(const Graph& g, CompMap comp_map, boost::false_type is_directed) const
     {
 	connected_components(g, comp_map);
     }
+    
 };
 
-template <class ComponentSelector>
-struct get_component_size_histogram
+void GraphInterface::LabelComponents(string prop)
 {
-    template <class Graph, class IndexMap, class Hist>
-    void operator()(const Graph &g, IndexMap index_map, Hist &hist) const
-    {   
-	typedef HashedDescriptorMap<IndexMap, size_t> comp_map_t;
-	comp_map_t comp_map(index_map);
+    typedef HashedDescriptorMap<vertex_index_map_t, size_t> comp_map_t;
+    static comp_map_t comp_map(_vertex_index);
 
-	_components(g, comp_map);
+    check_filter(*this, bind<void>(label_components(), _1, comp_map), reverse_check(), directed_check());
 
-	tr1::unordered_map<size_t, size_t> comp_count;
-	typename graph_traits<Graph>::vertex_iterator v, v_begin, v_end;
-	tie(v_begin, v_end) = vertices(g);
-	for(v = v_begin; v != v_end; ++v)
-	    comp_count[comp_map[*v]]++;
-	for(typeof(comp_count.begin()) iter = comp_count.begin(); iter != comp_count.end(); ++iter)
-	    hist[iter->second]++;
+    try
+    {
+	find_property_map(_properties, prop, typeid(graph_traits<multigraph_t>::vertex_descriptor));
+	RemoveVertexProperty(prop);
     }
-    ComponentSelector _components;
-};
-
-GraphInterface::hist_t GraphInterface::GetComponentSizeHistogram() const
-{
+    catch (property_not_found) {}
 
-    hist_t hist;
-    if (_directed)
-	check_filter(*this, bind<void>(get_component_size_histogram<strong_component>(), _1, _vertex_index, var(hist)),
-		     reverse_check(),always_directed()); 
-    else
-	check_filter(*this, bind<void>(get_component_size_histogram<connected_component>(), _1, _vertex_index, var(hist)),
-		     reverse_check(),always_undirected());
-    return hist;
+    _properties.property(prop, comp_map);
 }
 
 //==============================================================================

src/graph/graph.hh

@@ -90,7 +90,7 @@ public:
     void SetExtendedClusteringToProperty(std::string property_prefix, size_t max_depth);
 
     // other
-    hist_t GetComponentSizeHistogram() const;
+    void   LabelComponents(std::string property);
     hist_t GetDistanceHistogram(std::string weight) const;
     hist_t GetSampledDistanceHistogram(std::string weight, size_t samples, size_t seed) const;
     double GetReciprocity() const;

src/graph/graph_bind.cc

@@ -222,7 +222,7 @@ BOOST_PYTHON_MODULE(libgraph_tool)
 	.def("GetNumberOfEdges", &GraphInterfaceWrap::GetNumberOfEdges)
 	.def("GetVertexHistogram", &GraphInterfaceWrap::GetVertexHistogram)
 	.def("GetEdgeHistogram", &GraphInterfaceWrap::GetEdgeHistogram)
-	.def("GetComponentSizeHistogram", &GraphInterfaceWrap::GetComponentSizeHistogram)
+	.def("LabelComponents", &GraphInterfaceWrap::LabelComponents)
 	.def("GetCombinedVertexHistogram", &GraphInterfaceWrap::GetCombinedVertexHistogram)
 	.def("GetAverageCombinedVertexCorrelation", &GraphInterfaceWrap::GetAverageCombinedVertexCorrelation)
 	.def("GetVertexCorrelationHistogram", &GraphInterfaceWrap::GetVertexCorrelationHistogram)

src/graph/graph_clustering.cc

@@ -171,7 +171,7 @@ struct set_clustering_to_property
 void GraphInterface::SetLocalClusteringToProperty(string property)
 {
     typedef HashedDescriptorMap<vertex_index_map_t,double> clust_map_t;
-    static clust_map_t clust_map(_vertex_index);
+    clust_map_t clust_map(_vertex_index);
 
     bool directed = _directed;
     _directed = false;

src/graph/graph_community.cc

@@ -410,7 +410,6 @@ public:
 	_Ks[s] += k; 
     }
 
-
     double operator()(size_t k, size_t s) const 
     {
 	size_t ks = 0;
@@ -432,66 +431,112 @@ class NNKSCorr
 public:
     NNKSCorr(Graph &g, CommunityMap s): _g(g)
     {
+	unordered_set<size_t> spins;
+
 	typename graph_traits<Graph>::vertex_iterator v,v_end;
 	for (tie(v,v_end) = vertices(_g); v != v_end; ++v)
 	{
 	    size_t k = out_degree_no_loops(*v, _g);
 	    _Nk[k]++;
-	    _kNs[k][s[*v]]++;
-	    unordered_set<typename graph_traits<Graph>::vertex_descriptor> neighbours;
-	    neighbours.insert(*v);
-	    typename graph_traits<Graph>::adjacency_iterator a,a_end;
-	    for (tie(a,a_end) = adjacent_vertices(*v, _g); a != a_end; ++a)
-		if (neighbours.find(*a) == neighbours.end())
-		{
-		    _Nkk[k][out_degree_no_loops(*a, _g)]++;
-		    neighbours.insert(*a);
-		}
+	    _Nks[k][s[*v]]++;
+	    spins.insert(s[*v]);
 	}
-
-	for (typeof(_Nk.begin()) k_iter = _Nk.begin(); k_iter != _Nk.end(); ++k_iter)
+	
+	size_t E = 0;
+	typename graph_traits<Graph>::edge_iterator e,e_end;
+	for (tie(e,e_end) = edges(_g); e != e_end; ++e)
 	{
-	    _degs.push_back(k_iter->first);
-	    for (tie(v,v_end) = vertices(_g); v != v_end; ++v)
+	    typename graph_traits<Graph>::vertex_descriptor s, t;
+	    s = source(*e,g);
+	    t = target(*e,g);
+	    if (s != t)
 	    {
-		size_t k = out_degree_no_loops(*v, _g);
-		_kNNs[k_iter->first][s[*v]] += _Nkk[k_iter->first][k] * _kNs[k][s[*v]]/double(_Nk[k_iter->first]*_Nk[k]);
+		size_t k1 = out_degree_no_loops(s, g);
+		size_t k2 = out_degree_no_loops(t, g);
+		_Pkk[k1][k2]++;
+		_Pkk[k2][k1]++;
+		E++;
 	    }
 	}
+
+	for (typeof(_Pkk.begin()) iter1 = _Pkk.begin(); iter1 != _Pkk.end(); ++iter1)
+	{
+	    double sum = 0;
+	    for (typeof(iter1->second.begin()) iter2 = iter1->second.begin(); iter2 != iter1->second.end(); ++iter2)
+		sum += iter2->second;
+	    for (typeof(iter1->second.begin()) iter2 = iter1->second.begin(); iter2 != iter1->second.end(); ++iter2)
+		iter2->second /= sum;
+	}
+
+	for (typeof(_Nk.begin()) k_iter = _Nk.begin(); k_iter != _Nk.end(); ++k_iter) 
+	{
+	    size_t k1 = k_iter->first;
+	    _degs.push_back(k1);
+	    for (typeof(spins.begin()) s_iter = spins.begin(); s_iter != spins.end(); ++s_iter)
+		for (typeof(_Nk.begin()) k_iter2 = _Nk.begin(); k_iter2 != _Nk.end(); ++k_iter2) 
+		{
+		    size_t k2 = k_iter2->first;
+		    if (_Nks[k2].find(*s_iter) != _Nks[k2].end())
+			_NNks[k1][*s_iter] +=  k1*_Pkk[k1][k2] * _Nks[k2][*s_iter]/double(_Nk[k2]);
+		}
+	}
     }
 
     void Update(size_t k, size_t old_s, size_t s)
     {
+
 	for (size_t i = 0; i < _degs.size(); ++i) 
 	{
-	    typeof(_Nkk[k].begin()) iter = _Nkk[k].find(_degs[i]);
-	    if (iter != _Nkk[k].end())
-	    {
-		size_t nkk = iter->second;
-		_kNNs[_degs[i]][old_s] -= nkk * _kNs[k][old_s]/double(_Nk[_degs[i]]*_Nk[k]);
-		if (_kNNs[_degs[i]][old_s] == 0.0)
-		    _kNNs[_degs[i]].erase(old_s);
-		_kNNs[_degs[i]][s] += nkk * _kNs[k][s]/double(_Nk[_degs[i]]*_Nk[k]);
-	    }
+	    size_t k1 = _degs[i], k2 = k;
+	    if (_Pkk.find(k1) == _Pkk.end())
+		continue;
+	    if (_Pkk.find(k1)->second.find(k2) == _Pkk.find(k1)->second.end())
+		continue;
+	    _NNks[k1][old_s] -=  k1*_Pkk[k1][k2] * _Nks[k2][old_s]/double(_Nk[k2]);
+	    if (_NNks[k1][old_s] == 0.0)
+		_NNks[k1].erase(old_s);
+	    if (_Nks[k2].find(s) != _Nks[k2].end())
+                _NNks[k1][s] -=  k1*_Pkk[k1][k2] * _Nks[k2][s]/double(_Nk[k2]);
+	    if (_NNks[k1][s] == 0.0)
+                _NNks[k1].erase(s);
 	}
+
+	_Nks[k][old_s]--;
+	if (_Nks[k][old_s] == 0)
+	    _Nks[k].erase(old_s);
+	_Nks[k][s]++;
+
+	for (size_t i = 0; i < _degs.size(); ++i) 
+	{
+	    size_t k1 = _degs[i], k2 = k;
+	    if (_Pkk.find(k1) == _Pkk.end())
+		continue;
+	    if (_Pkk.find(k1)->second.find(k2) == _Pkk.find(k1)->second.end())
+		continue;
+	    _NNks[k1][old_s] +=  k1*_Pkk[k1][k2] * _Nks[k2][old_s]/double(_Nk[k2]);
+	    if (_NNks[k1][old_s] == 0.0)
+		_NNks[k1].erase(old_s);
+	    _NNks[k1][s] +=  k1*_Pkk[k1][k2] * _Nks[k2][s]/double(_Nk[k2]);
+	}
+
     }
 
-    double operator()(size_t k, size_t s)
+    double operator()(size_t k, size_t s) const
     {
-	double knns = 0.0;
-	typeof(_kNNs[k].begin()) iter = _kNNs[k].find(s);
-	if (iter != _kNNs[k].end())
-	    knns = iter->second;
-	return knns;
+	const typeof(_NNks[k])& nnks = _NNks.find(k)->second;
+	const typeof(nnks.begin()) iter = nnks.find(s);
+	if (iter != nnks.end())
+	    return iter->second;
+	return 0.0;
     }
 
 private:
     Graph& _g;
     vector<size_t> _degs;
-    ManagedUnorderedMap<size_t, size_t> _Nk;
-    ManagedUnorderedMap<size_t, ManagedUnorderedMap<size_t, size_t> > _Nkk;
-    ManagedUnorderedMap<size_t, ManagedUnorderedMap<size_t,size_t> > _kNs;
-    ManagedUnorderedMap<size_t, ManagedUnorderedMap<size_t, double> > _kNNs;
+    unordered_map<size_t, size_t> _Nk;
+    unordered_map<size_t, unordered_map<size_t,double> > _Pkk;
+    unordered_map<size_t, ManagedUnorderedMap<size_t,size_t> > _Nks;
+    unordered_map<size_t, ManagedUnorderedMap<size_t,double> > _NNks;
 };
 
 
@@ -583,29 +628,25 @@ struct get_modularity
 	size_t E = 0;
 	double W = 0;
 
-	typename graph_traits<Graph>::edge_iterator e,e_end;
+	typename graph_traits<Graph>::edge_iterator e, e_end;
 	for (tie(e,e_end) = edges(g); e != e_end; ++e)
 	    if (target(*e,g) != source(*e,g))
 	    {
 		W += get(weights, *e);
 		E++;
+		if (get(s, target(*e,g)) == get(s, source(*e,g)))
+		    modularity += 2*get(weights, *e);
 	    }
 
-	typename graph_traits<Graph>::vertex_iterator v,v_end;
+	unordered_map<size_t, size_t> Ks;
+
+	typename graph_traits<Graph>::vertex_iterator v, v_end;
 	for (tie(v,v_end) = vertices(g); v != v_end; ++v)
-	{
-	    size_t k = out_degree_no_loops(*v,g);
-	    typename graph_traits<Graph>::out_edge_iterator e,e_end;
-	    for(tie(e, e_end) = out_edges(*v,g); e != e_end; ++e)
-	    {
-		typename graph_traits<Graph>::vertex_descriptor t = target(*e,g);
-		if (t != *v)
-		{
-		    if (get(s, t) == get(s, *v))
-			modularity += get(weights, *e) - k*out_degree_no_loops(t,g)/double(2*E);
-		}
-	    }
-	}
+	    Ks[get(s, *v)] += out_degree_no_loops(*v, g);
+	
+	for (typeof(Ks.begin()) iter = Ks.begin(); iter != Ks.end(); ++iter)
+	    modularity -= (iter->second*iter->second)/double(2*E);
+
 	modularity /= 2*W;
     }
 };

src/graph/graph_io.cc

@@ -341,7 +341,7 @@ void GraphInterface::WriteToFile(string file, string format)
 	    typedef tr1::unordered_map<graph_traits<multigraph_t>::vertex_descriptor, size_t>  map_t;
 	    map_t vertex_to_index;
 	    associative_property_map<map_t> index_map(vertex_to_index);
-	    check_filter(*this,bind<void>(generate_index(),_1,index_map),reverse_check(),always_directed());
+	    check_filter(*this, bind<void>(generate_index(), _1, index_map), reverse_check(), directed_check());
 	    if (graphviz)
 	    {
 		try
@@ -355,12 +355,12 @@ void GraphInterface::WriteToFile(string file, string format)
 	    }
 	    if (GetDirected())
 	    {
-		check_filter(*this,bind<void>(write_to_file(),var(stream),_1,index_map,var(dp),graphviz),
-			     reverse_check(),always_directed());
+		check_filter(*this,bind<void>(write_to_file(),var(stream), _1, index_map, var(dp), graphviz),
+			     reverse_check(), always_directed());
 	    }
 	    else
 	    {
-		check_filter(*this,bind<void>(write_to_file_fake_undir(),var(stream),_1,index_map,var(dp),graphviz),
+		check_filter(*this,bind<void>(write_to_file_fake_undir(), var(stream), _1, index_map, var(dp), graphviz),
                              never_reversed(), always_undirected());
 	    }
 	}