Commit 412abc9a authored by Tiago Peixoto's avatar Tiago Peixoto
Browse files

BlockState.sample_graph(): Use correct probabilities for undirected graphs

parent d305ec12
......@@ -28,7 +28,7 @@ void generate_sbm(GraphInterface& gi, boost::any ab, boost::python::object ors,
{
auto rs = get_array<int64_t, 1>(ors);
auto ss = get_array<int64_t, 1>(oss);
auto probs = get_array<double, 2>(oprobs);
auto probs = get_array<double, 1>(oprobs);
typedef vprop_map_t<int32_t>::type bmap_t;
auto b = any_cast<bmap_t>(ab).get_unchecked();
......@@ -37,6 +37,6 @@ void generate_sbm(GraphInterface& gi, boost::any ab, boost::python::object ors,
auto in_deg = any_cast<dmap_t>(ain_deg).get_unchecked();
auto out_deg = any_cast<dmap_t>(aout_deg).get_unchecked();
run_action<graph_tool::detail::always_directed_never_reversed>()
run_action<>()
(gi, [&](auto& g) { gen_sbm(g, b, rs, ss, probs, in_deg, out_deg, rng); })();
}
......@@ -64,12 +64,14 @@ void gen_sbm(Graph& g, VProp b, IVec& rs, IVec& ss, FVec probs, VDProp in_deg,
v_out_sampler.emplace_back(rvs[r], v_out_probs[r]);
}
for (size_t i = 0; i < rs.shape()[0]; ++i)
{
size_t r = rs[i];
size_t s = ss[i];
double p = probs[0][i];
double p = probs[i];
if (!is_directed::apply<Graph>::type::value && r == s)
p /= 2;
if (r >= v_out_sampler.size() || v_out_sampler[r].prob_sum() == 0 ||
s >= v_in_sampler.size() || v_in_sampler[s].prob_sum() == 0)
......
......@@ -851,13 +851,16 @@ def generate_sbm(b, probs, out_degs=None, in_degs=None, directed=False):
probs : two-dimensional :class:`numpy.ndarray` or :class:`scipy.sparse.spmatrix`
Matrix with edge propensities between groups. The value ``probs[r,s]``
corresponds to the average number of edges between groups ``r`` and
``s``.
``s`` (or twice the average number if ``r == s`` and the graph is
undirected).
out_degs : iterable or :class:`numpy.ndarray` (optional, default: ``None``)
Expected out-degree for each node. If not provided, a constant value
will be used.
Out-degree propensity for each node. If not provided, a constant value
will be used. Note that the values will be normalized inside each group,
if they are not already so.
in_degs : iterable or :class:`numpy.ndarray` (optional, default: ``None``)
Expected in-degree for each node. If not provided, a constant value
will be used. This parameter is ignored if ``directed == False``.
In-degree propensity for each node. If not provided, a constant value
will be used. Note that the values will be normalized inside each group,
if they are not already so.
directed : ``bool`` (optional, default: ``False``)
Whether the graph is directed.
......@@ -879,12 +882,19 @@ def generate_sbm(b, probs, out_degs=None, in_degs=None, directed=False):
P({\boldsymbol A}|{\boldsymbol \lambda},{\boldsymbol \theta},{\boldsymbol b})
= \prod_{i<j}\frac{e^{-\lambda_{b_ib_j}\theta_i\theta_j}(\lambda_{b_ib_j}\theta_i\theta_j)^{A_{ij}}}{A_{ij}!}
\times\prod_i\frac{e^{-\lambda_{b_ib_i}\theta_i^2}(\lambda_{b_ib_i}\theta_i^2)^{A_{ij}/2}}{(A_{ij}/2)!},
\times\prod_i\frac{e^{-\lambda_{b_ib_i}\theta_i^2/2}(\lambda_{b_ib_i}\theta_i^2/2)^{A_{ij}/2}}{(A_{ij}/2)!},
where :math:`\lambda_{rs}` is the edge propensity between groups :math:`r`
and :math:`s`, and :math:`\theta_i` is the propensity of node i to receive
edges, which is proportional to its expected degree. For directed graphs,
likelihood is analogous:
edges, which is proportional to its expected degree. Note that in the
algorithm it is assumed that the node propensities are normalized for each
group,
.. math::
\sum_i\theta_i\delta_{b_i,r} = 1.
For directed graphs, the probability is analogous:
.. math::
......@@ -931,7 +941,6 @@ def generate_sbm(b, probs, out_degs=None, in_degs=None, directed=False):
if not directed:
if out_degs is None:
out_degs = in_degs = g.new_vp("double", 1)
else:
out_degs = in_degs = g.new_vp("double", out_degs)
else:
......@@ -949,7 +958,9 @@ def generate_sbm(b, probs, out_degs=None, in_degs=None, directed=False):
idx = r <= s
r = r[idx]
s = s[idx]
p = probs[r, s]
p = numpy.squeeze(probs[r, s])
g.set_directed(directed)
libgraph_tool_generation.gen_sbm(g._Graph__graph,
_prop("v", g, b),
......@@ -959,8 +970,6 @@ def generate_sbm(b, probs, out_degs=None, in_degs=None, directed=False):
_prop("v", g, in_degs),
_prop("v", g, out_degs),
_get_rng())
g.set_directed(directed)
return g
def predecessor_tree(g, pred_map):
......
......@@ -2054,8 +2054,8 @@ class BlockState(object):
in_degs = self.g.degree_property_map("in").fa
else:
in_degs = None
g = generate_sbm(b=self.b.fa,
probs=adjacency(self.bg, weight=self.mrs),
probs = adjacency(self.bg, weight=self.mrs)
g = generate_sbm(b=self.b.fa, probs=probs,
in_degs=in_degs, out_degs=out_degs,
directed=self.g.is_directed())
if not multigraph:
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment