utils/fginfo.cpp

   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
   4
   5     This file is part of libDAI.
   6
   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.
  11
  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.
  16
  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
  20 */
  21
  22
  23 #include <iostream>
  24 #include <cstdlib>
  25 #include <dai/factorgraph.h>
  26 #include <dai/jtree.h>
  27
  28
  29 using namespace std;
  30 using namespace dai;
  31
  32
  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 );
  38         }
  39         else if( length > 1 )
  40             findLoopClusters( fg, allcl, newcl | *in, root, length - 1, ind / newcl );
  41     }
  42 }
  43
  44
  45 size_t countLoops( const FactorGraph & fg, size_t loopdepth ) {
  46     set<VarSet> loops;
  47     for( vector<Var>::const_iterator i0 = fg.vars().begin(); i0 != fg.vars().end(); i0++ ) {
  48         VarSet i0d = fg.delta(*i0);
  49         if( loopdepth > 1 )
  50             findLoopClusters( fg, loops, *i0, *i0, loopdepth - 1, i0d );
  51     }
  52     return loops.size();
  53 }
  54
  55
  56 bool hasShortLoops( const std::vector<Factor> &P ) {
  57     bool found = false;
  58     vector<Factor>::const_iterator I, J;
  59     for( I = P.begin(); I != P.end(); I++ ) {
  60         J = I;
  61         J++;
  62         for( ; J != P.end(); J++ )
  63             if( (I->vars() & J->vars()).size() >= 2 ) {
  64                 found = true;
  65                 break;
  66             }
  67         if( found )
  68             break;
  69     }
  70     return found;
  71 }
  72
  73
  74 bool hasNegatives( const std::vector<Factor> &P ) {
  75     bool found = false;
  76     for( size_t I = 0; I < P.size(); I++ )
  77         if( P[I].hasNegatives() ) {
  78             found = true;
  79             break;
  80         }
  81     return found;
  82 }
  83
  84
  85 int main( int argc, char *argv[] ) {
  86     if( argc != 3 ) {
  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;
  90         return 1;
  91     } else {
  92         // Read factorgraph
  93         FactorGraph fg;
  94         char *infile = argv[1];
  95         int calc_tw = atoi(argv[2]);
  96         fg.ReadFromFile( infile );
  97
  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;
 106         if( calc_tw ) {
 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;
 110         }
 111         double stsp = 1.0;
 112         for( size_t i = 0; i < fg.nrVars(); i++ )
 113             stsp *= fg.var(i).states();
 114         cout << "Total state space:   " << stsp << endl;
 115
 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()]++;
 124             else
 125                 cavsizes[di.size()] = 1;
 126             size_t Ds = fg.Delta(i).nrStates();
 127             if( Ds > max_Delta_size )
 128                 max_Delta_size = Ds;
 129             cavsum_lcbp += di.nrStates();
 130             for( VarSet::const_iterator j = di.begin(); j != di.end(); j++ )
 131                 cavsum_lcbp2 += j->states();
 132         }
 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 << ") ";
 139         cout << endl;
 140
 141         cout << "Type: " << (fg.isPairwise() ? "pairwise" : "higher order") << " interactions, " << (fg.isBinary() ? "binary" : "nonbinary") << " variables" << endl;
 142
 143         if( fg.isPairwise() ) {
 144             bool girth_reached = false;
 145             size_t loopdepth;
 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;
 149                 if( nr_loops > 0 )
 150                     girth_reached = true;
 151             }
 152             if( girth_reached )
 153                 cout << "Girth: " << loopdepth-1 << endl;
 154             else
 155                 cout << "Girth: infinity" << endl;
 156         }
 157
 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 ) )
 162                 facsizes[Isize]++;
 163             else
 164                 facsizes[Isize] = 1;
 165         }
 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 << ") ";
 169         cout << endl;
 170
 171         return 0;
 172     }
 173 }