Commit bf3c8e36 authored by Tiago Peixoto's avatar Tiago Peixoto

uncertail_blockmodel: implement marginal_multigraph_sample and marginal_graph_sample

parent 36917e1f
......@@ -67,7 +67,13 @@ void collect_marginal_count_dispatch(GraphInterface& gi, GraphInterface& ui,
boost::any aecount, boost::any aexs,
boost::any aexc);
void marginal_count_entropy(GraphInterface& gi, boost::any aexc, boost::any aeh);
double marginal_count_entropy(GraphInterface& gi, boost::any aexc, boost::any aeh);
double marginal_multigraph_sample(GraphInterface& gi, boost::any axs,
boost::any axc, boost::any ax, rng_t& rng);
double marginal_graph_sample(GraphInterface& gi, boost::any ap,
boost::any ax, rng_t& rng);
void export_uncertain_state()
{
......@@ -123,4 +129,6 @@ void export_uncertain_state()
def("collect_xmarginal", &collect_xmarginal_dispatch);
def("collect_marginal_count", &collect_marginal_count_dispatch);
def("marginal_count_entropy", &marginal_count_entropy);
def("marginal_multigraph_sample", &marginal_multigraph_sample);
def("marginal_graph_sample", &marginal_graph_sample);
}
......@@ -81,30 +81,92 @@ void collect_marginal_count_dispatch(GraphInterface& gi, GraphInterface& ui,
ui.get_graph_view());
}
void marginal_count_entropy(GraphInterface& gi, boost::any aexc, boost::any aeh)
double marginal_count_entropy(GraphInterface& gi, boost::any aexc, boost::any aeh)
{
typedef eprop_map_t<double>::type ehmap_t;
auto eh = any_cast<ehmap_t>(aeh);
double S_tot = 0;
gt_dispatch<>()
([&](auto& g, auto exc)
{
parallel_edge_loop
(g,
[&](auto& e)
{
auto& S = eh[e];
S = 0;
size_t N = 0;
for (auto n : exc[e])
{
S -= xlogx_fast(n);
N += n;
}
if (N == 0)
return;
S /= N;
S += safelog_fast(N);
#pragma omp atomic
S_tot += S;
});
},
all_graph_views(), edge_scalar_vector_properties())
(gi.get_graph_view(), aexc);
return S_tot;
}
double marginal_multigraph_sample(GraphInterface& gi, boost::any axs, boost::any axc,
boost::any ax, rng_t& rng)
{
double L = 0;
gt_dispatch<>()
([&](auto& g, auto& xs, auto& xc, auto& x)
{
for (auto e : edges_range(g))
{
auto& S = eh[e];
S = 0;
size_t N = 0;
for (auto n : exc[e])
typedef std::remove_reference_t<decltype(xs[e][0])> val_t;
std::vector<double> probs(xc[e].begin(), xc[e].end());
Sampler<val_t> sample(xs[e], probs);
x[e] = sample.sample(rng);
size_t Z = 0;
size_t p = 0;
for (size_t i = 0; i < xs[e].size(); ++i)
{
S -= xlogx_fast(n);
N += n;
auto m = xs[e][i];
if (m == x[e])
p = xc[e][i];
Z += xc[e][i];
}
if (N == 0)
continue;
S /= N;
S += safelog_fast(N);
L += std::log(p) - std::log(Z);
}
},
all_graph_views(), edge_scalar_vector_properties())
(gi.get_graph_view(), aexc);
all_graph_views(), edge_scalar_vector_properties(),
edge_scalar_vector_properties(), writable_edge_scalar_properties())
(gi.get_graph_view(), axs, axc, ax);
return L;
}
double marginal_graph_sample(GraphInterface& gi, boost::any ap,
boost::any ax, rng_t& rng)
{
double L = 0;
gt_dispatch<>()
([&](auto& g, auto& p, auto& x)
{
for (auto e : edges_range(g))
{
std::bernoulli_distribution sample(p[e]);
x[e] = sample(rng);
if (x[e] == 1)
L += std::log(p[e]);
else
L += std::log1p(-p[e]);
}
},
all_graph_views(), edge_scalar_properties(),
writable_edge_scalar_properties())
(gi.get_graph_view(), ap, ax);
return L;
}
......@@ -187,6 +187,8 @@ __all__ = ["minimize_blockmodel_dl",
"bethe_entropy",
"microstate_entropy",
"marginal_multigraph_entropy",
"marginal_multigraph_sample",
"marginal_graph_sample",
"PartitionHist",
"BlockPairHist",
"half_edge_graph",
......
......@@ -1518,3 +1518,21 @@ def marginal_multigraph_entropy(g, ecount):
_prop("e", g, ecount),
_prop("e", g, eh))
return eh
def marginal_multigraph_sample(g, ew, ecount):
w = g.new_ep("int")
L = libinference.marginal_multigraph_sample(g._Graph__graph,
_prop("e", g, ew),
_prop("e", g, ecount),
_prop("e", g, w),
_get_rng())
return w, L
def marginal_graph_sample(g, ep):
w = g.new_ep("int")
L = libinference.marginal_graph_sample(g._Graph__graph,
_prop("e", g, ep),
_prop("e", g, w),
_get_rng())
return w, L
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