1 /* Copyright (C) 2006-2008 Joris Mooij [joris dot mooij at tuebingen dot mpg dot de]
2 Radboud University Nijmegen, The Netherlands /
3 Max Planck Institute for Biological Cybernetics, Germany
5 This file is part of libDAI.
7 libDAI is free software; you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation; either version 2 of the License, or
10 (at your option) any later version.
12 libDAI is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
15 GNU General Public License for more details.
17 You should have received a copy of the GNU General Public License
18 along with libDAI; if not, write to the Free Software
19 Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 #include <dai/factorgraph.h>
26 #include <dai/jtree.h>
33 void findLoopClusters( const FactorGraph
& fg
, std::set
<VarSet
> &allcl
, VarSet newcl
, const Var
& root
, size_t length
, VarSet vars
) {
34 for( VarSet::const_iterator in
= vars
.begin(); in
!= vars
.end(); in
++ ) {
35 VarSet ind
= fg
.delta( *in
);
36 if( (newcl
.size()) >= 2 && ind
.contains( root
) ) {
37 allcl
.insert( newcl
| *in
);
40 findLoopClusters( fg
, allcl
, newcl
| *in
, root
, length
- 1, ind
/ newcl
);
45 size_t countLoops( const FactorGraph
& fg
, size_t loopdepth
) {
47 for( vector
<Var
>::const_iterator i0
= fg
.vars().begin(); i0
!= fg
.vars().end(); i0
++ ) {
48 VarSet i0d
= fg
.delta(*i0
);
50 findLoopClusters( fg
, loops
, *i0
, *i0
, loopdepth
- 1, i0d
);
56 bool hasShortLoops( const std::vector
<Factor
> &P
) {
58 vector
<Factor
>::const_iterator I
, J
;
59 for( I
= P
.begin(); I
!= P
.end(); I
++ ) {
62 for( ; J
!= P
.end(); J
++ )
63 if( (I
->vars() & J
->vars()).size() >= 2 ) {
74 bool hasNegatives( const std::vector
<Factor
> &P
) {
76 for( size_t I
= 0; I
< P
.size(); I
++ )
77 if( P
[I
].hasNegatives() ) {
85 int main( int argc
, char *argv
[] ) {
87 cout
<< "Usage: " << argv
[0] << " <in.fg> <tw>" << endl
<< endl
;
88 cout
<< "Reports some characteristics of the .fg network." << endl
;
89 cout
<< "Also calculates treewidth (which may take some time) unless <tw> == 0." << endl
;
94 char *infile
= argv
[1];
95 int calc_tw
= atoi(argv
[2]);
96 fg
.ReadFromFile( infile
);
98 cout
<< "Number of variables: " << fg
.nrVars() << endl
;
99 cout
<< "Number of factors: " << fg
.nrFactors() << endl
;
100 cout
<< "Connected: " << fg
.isConnected() << endl
;
101 cout
<< "Tree: " << fg
.isTree() << endl
;
102 cout
<< "Has short loops: " << hasShortLoops(fg
.factors()) << endl
;
103 cout
<< "Has negatives: " << hasNegatives(fg
.factors()) << endl
;
104 cout
<< "Binary variables? " << fg
.isBinary() << endl
;
105 cout
<< "Pairwise interactions? " << fg
.isPairwise() << endl
;
107 std::pair
<size_t,size_t> tw
= treewidth(fg
);
108 cout
<< "Treewidth: " << tw
.first
<< endl
;
109 cout
<< "Largest cluster for JTree has " << tw
.second
<< " states " << endl
;
112 for( size_t i
= 0; i
< fg
.nrVars(); i
++ )
113 stsp
*= fg
.var(i
).states();
114 cout
<< "Total state space: " << stsp
<< endl
;
116 double cavsum_lcbp
= 0.0;
117 double cavsum_lcbp2
= 0.0;
118 size_t max_Delta_size
= 0;
119 map
<size_t,size_t> cavsizes
;
120 for( size_t i
= 0; i
< fg
.nrVars(); i
++ ) {
121 VarSet di
= fg
.delta(i
);
122 if( cavsizes
.count(di
.size()) )
123 cavsizes
[di
.size()]++;
125 cavsizes
[di
.size()] = 1;
126 size_t Ds
= fg
.Delta(i
).nrStates();
127 if( Ds
> max_Delta_size
)
129 cavsum_lcbp
+= di
.nrStates();
130 for( VarSet::const_iterator j
= di
.begin(); j
!= di
.end(); j
++ )
131 cavsum_lcbp2
+= j
->states();
133 cout
<< "Maximum pancake has " << max_Delta_size
<< " states" << endl
;
134 cout
<< "LCBP with full cavities needs " << cavsum_lcbp
<< " BP runs" << endl
;
135 cout
<< "LCBP with only pairinteractions needs " << cavsum_lcbp2
<< " BP runs" << endl
;
136 cout
<< "Cavity sizes: ";
137 for( map
<size_t,size_t>::const_iterator it
= cavsizes
.begin(); it
!= cavsizes
.end(); it
++ )
138 cout
<< it
->first
<< "(" << it
->second
<< ") ";
141 cout
<< "Type: " << (fg
.isPairwise() ? "pairwise" : "higher order") << " interactions, " << (fg
.isBinary() ? "binary" : "nonbinary") << " variables" << endl
;
143 if( fg
.isPairwise() ) {
144 bool girth_reached
= false;
146 for( loopdepth
= 2; loopdepth
<= fg
.nrVars() && !girth_reached
; loopdepth
++ ) {
147 size_t nr_loops
= countLoops( fg
, loopdepth
);
148 cout
<< "Loops up to " << loopdepth
<< " variables: " << nr_loops
<< endl
;
150 girth_reached
= true;
153 cout
<< "Girth: " << loopdepth
-1 << endl
;
155 cout
<< "Girth: infinity" << endl
;
158 map
<size_t,size_t> facsizes
;
159 for( size_t I
= 0; I
< fg
.nrFactors(); I
++ ) {
160 size_t Isize
= fg
.factor(I
).vars().size();
161 if( facsizes
.count( Isize
) )
166 cout
<< "Factor sizes: ";
167 for( map
<size_t,size_t>::const_iterator it
= facsizes
.begin(); it
!= facsizes
.end(); it
++ )
168 cout
<< it
->first
<< "(" << it
->second
<< ") ";