gramod2.c

数据挖掘中的一算法 ines算法 c下实现的。适合初学习数据挖掘者借鉴

  csel.ptcnt     = (maxtest < n) ? n   : maxtest;
  csel.comps     = NULL;        /* clear for a proper error cleanup */
  csel.step      = 0;
  csel.verbose   = verbose;
  csel.edges = e = (EDGE*)  malloc(n *n *sizeof(EDGE));
  csel.etab      = (EDGE**) malloc(n *n *sizeof(EDGE*));
  csel.ptrees    = (PTREE**)calloc(csel.ptcnt, sizeof(PTREE*));
  if (!e || !csel.etab || !csel.ptrees) {
    _csel_clean(&csel); return GM_ERROR; }
  csel.seles     = csel.etab +n;
  for (e += n*n, i = n; --i >= 0; ) {
    for (k = n; --k >= 0; ) {   /* traverse the table columns */
      (--e)->y = k; e->x = i; e->mark = 0; e->eval = 0; }
    csel.etab[i] = e;           /* initialize the edge variables and */
  }                             /* organize the edges into a table */
  switch (search) {             /* evaluate the search method */
    case GM_OWST:   i = _owst  (&csel); break;
    case GM_EXTST:  i = _extst (&csel); break;
    case GM_TOPORD: i = _topord(&csel); break;
    case GM_NOLOOP: i = _noloop(&csel); break;
    default:        i = _topord(&csel); break;
  }                             /* induce a network structure */
  _csel_clean(&csel);           /* clean up the induction variables */
  if (i != 0) return GM_ERROR;  /* and check for an error */
  gm_clear(gm, 1);              /* clear the graphical model */
  if ((gm_estim(gm, table, 0) != 0)
  ||  (gm_eval(gm, measure, params) == GM_ERROR))
    return GM_ERROR;            /* rebuild the graphical model */
  if ((search != GM_OWST) && (search != GM_EXTST)
  &&  (lwise >= 0) && (gm_local(gm, lwise, mrgimp) == GM_ERROR))
    return GM_ERROR;            /* redo local structure learning */
  return gm->eval;              /* return the evaluation result */
}  /* gm_csel() */

/*----------------------------------------------------------------------
  Simulated Annealing
----------------------------------------------------------------------*/

static void _sian_clean (SIAN *sian)
{                               /* --- clean up induction variables */
  assert(sian);                 /* check the function argument */
  if (sian->table != sian->test)
    tab_delete(sian->table, 0); /* delete the estimation table */
  if (sian->best) ht_delete(sian->best);
  if (sian->curr) ht_delete(sian->curr);
  if (sian->tmp)  ht_delete(sian->tmp);
}  /* _sian_clean() */          /* delete the hypertrees */

/*--------------------------------------------------------------------*/

static int _ht2gm (HTREE *ht, GRAMOD *gm)
{                               /* --- build graphical model */
  int i, k, n;                  /* loop variables */
  int *nodes;                   /* to traverse the node identifiers */
  int *marks;                   /* node markers */

  assert(ht && gm               /* check the function argument */
      && (ht_nodecnt(ht) == gm->attcnt));
  if (ht_reduce(ht) != 0)       /* do Graham reduction to determine */
    return -1;                  /* a construction sequence */
  for (i = gm->attcnt; --i >= 0; )
    gm_conrem(gm, i, -1);       /* remove all conditions */
  for (marks = gm->buf +(i = gm->attcnt); --i >= 0; )
    *--marks = 0;               /* initialize the node markers */
  if (gm->type == GM_PROB) {    /* --- if in probabilistic mode */
    for (n = 0; n < ht_hecnt(ht); n++) {
      nodes = ht_henodes(ht, n);/* traverse the construction sequence */
      for (i = ht_hesize(ht, n); --i >= 0; ) {
        if (marks[nodes[i]])    /* traverse all unprocessed */
          continue;             /* nodes of the hyperedge */
        for (k = ht_hesize(ht, n); --k >= 0; ) {
          if (!marks[nodes[k]]) /* traverse all processed */
            continue;           /* nodes of the hyperedge */
          if (gm_conadd(gm, nodes[i], nodes[k]) < 0) return -1;
        }                       /* add all processed nodes */
        marks[nodes[i]] = 1;    /* as conditions/parents and */
      }                         /* mark the node as processed */
    } }
  else {                        /* --- if in possibilistic mode */
    for (n = 0; n < ht_hecnt(ht); n++) {
      nodes = ht_henodes(ht, n);/* traverse the construction sequence */
      for (i = ht_hesize(ht, n); --i >= 0; )
        if (!marks[nodes[i]])   /* traverse the nodes and */
          break;                /* find an unprocessed node */
      assert(i >= 0);           /* check whether a node was found */
      for (k = ht_hesize(ht, n); --k >= 0; ) {
        marks[nodes[k]] = 1;    /* mark the nodes as processed, */
        if (k == i) continue;   /* but then skip the chosen node */
        if (gm_conadd(gm, nodes[i], nodes[k]) < 0) return -1;
      }                         /* add all other nodes */
    }                           /* as conditions/parents */
  }
  return 0;                     /* return 'ok' */
}  /* _ht2gm() */

/*--------------------------------------------------------------------*/

double gm_sian (GRAMOD *gm, TABLE *table, int trials, int maxsz,
                double keep, double frac, double ppwgt, int distuv,
                double randfn (void), int verbose)
{                               /* --- simulated annealing */
  SIAN   sian;                  /* induction variables */
  int    step  = 0;             /* simulated annealing step */
  double worst =  DBL_MAX;      /* worth of worst hypertree found */
  double best  = -DBL_MAX;      /* worth of best hypertree found */
  double curr  = -DBL_MAX;      /* worth of current hypertree */
  double tmp;                   /* worth of test hypertree */
  double range;                 /* range of evaluation results */
  double fact = 1, mul;         /* factor for deviation */
  double res[3];                /* evaluation results */

  assert(gm && table            /* check the function arguments */
      && (trials > 0) && (maxsz >= 1) && (keep >= 0) && (keep <= 1)
      && (frac  >= 0) && (ppwgt >= 0));
  sian.table = sian.test = table;  /* note the data table and */
  sian.best  = sian.curr = NULL;   /* create an empty hypertree */
  sian.tmp   = ht_create(gm->attcnt);
  if (!sian.tmp) { _sian_clean(&sian); return GM_ERROR; }
  if (gm->type == GM_POSS) {    /* if to induce a poss. network */
    sian.table = tab_dup(table, 0);
    if (!sian.table || (tab_opc(sian.table, TAB_COND) != 0)) {
      _sian_clean(&sian); return GM_ERROR; }
  }                             /* compute one point coverages */
  gm_setnorm(gm, tab_getwgt(table, 0, INT_MAX));
  if (frac <= 0) frac = 0.001;  /* check and adapt the fraction */
  mul = pow(frac, 1.0/trials);  /* and compute the multiplier */
  while (1) {                   /* simulated annealing loop */
    if (((++step & 0x0f) == 0) && verbose)
      fprintf(stderr, "%8d\b\b\b\b\b\b\b\b", step);
    if ((_ht2gm(sian.tmp, gm) != 0)
    ||  (gm_estim(gm, sian.table, distuv) != 0)
    ||  (gm_evalx(gm, sian.test, res)     <  0)) {
      _sian_clean(&sian); return GM_ERROR; }
    tmp = ((gm->type & GM_POSS) ? -res[2] : res[0])
        - ppwgt *gm_parsum(gm); /* evaluate the graphical model */
    if (tmp < worst)            /* if it is worse than the worst, */
      worst = tmp;              /* note its evaluation */
    if (tmp > best) {           /* if it is better than the best */
      if (sian.best) ht_delete(sian.best);
      sian.best = ht_dup(sian.tmp);
      if (!sian.best) { _sian_clean(&sian); return GM_ERROR; }
      best = tmp;               /* copy the hypertree */
    }                           /* and its evaluation */
    range = ((step +1.0) /step) *(best -worst);
    if (range <= 0) range = 1;  /* compute estimate of eval. range */
    if ((tmp > curr)            /* alway select better result and */
    || ((tmp < curr)            /* select worse result randomly */
    &&  (randfn() *range *fact > curr -tmp))) {
      if (sian.curr) ht_delete(sian.curr);
      sian.curr = sian.tmp; sian.tmp = NULL;
      curr = tmp;               /* copy the hypertree */
    }                           /* and its evaluation */
    if (sian.tmp) { ht_delete(sian.tmp); sian.tmp = NULL; }
    if (step >= trials) break;  /* if all trials done, abort loop */
    sian.tmp = ht_modify(sian.curr, maxsz, keep, randfn);
    if (!sian.tmp) { _sian_clean(&sian); return GM_ERROR; }
    fact *= mul;                /* create modified hypertree and */
  }                             /* decrease the probability factor */
  if (sian.curr)  { ht_delete(sian.curr); sian.curr = NULL; }
  if (!sian.best) { _sian_clean(&sian); return GM_ERROR; }
  if ((_ht2gm(sian.best, gm) != 0)
  ||  (gm_estim(gm, sian.table, distuv) != 0)
  ||  (gm_evalx(gm, sian.test, res)     <  0)) {
    _sian_clean(&sian); return GM_ERROR; }
  best = ((gm->type & GM_POSS) ? -res[2] : res[0])
       - ppwgt *gm_parsum(gm);  /* reevaluate the best model */
  _sian_clean(&sian);           /* clean up the induction variables */
  return best;                  /* return the evaluation result */
}  /* gm_sian() */

#endif
/*----------------------------------------------------------------------
  Other Main Functions
----------------------------------------------------------------------*/

int gm_estim (GRAMOD *gm, TABLE *table, int distuv)
{                               /* --- estimate model parameters */
  int   i, k;                   /* loop variables */
  GMATT *att;                   /* to traverse the attributes */
  float wgt;                    /* weight of the instantiation */

  assert(gm && table);          /* check the function arguments */
  for (i = tab_tplcnt(table); --i >= 0; ) {
    tpl_toas(tab_tpl(table,i)); /* traverse the tuples */
    wgt = as_getwgt(gm->attset);/* get the instantiation weight */
    for (att = gm->atts +(k = gm->attcnt); --k >= 0; )
      if ((--att)->ptree && (pt_aggr(att->ptree, wgt, distuv) < 0))
        return -1;              /* aggregate the instantiation weight */
  }                             /* for all attributes */
  return 0;                     /* return 'ok' */
}  /* gm_estim() */

/*--------------------------------------------------------------------*/

double gm_eval (GRAMOD *gm, int measure, double *params)
{                               /* --- evaluate a graphical model */
  int    i;                     /* loop variable */
  GMATT  *att;                  /* to traverse the attributes */

  assert(gm);                   /* check the function arguments */
  gm->eval = 0;                 /* init. the evaluation result */
  for (att = gm->atts +(i = gm->attcnt); --i >= 0; ) {
    if (!(--att)->ptree) continue;
    gm->eval += att->eval = pt_eval(att->ptree, measure, params);
    if (att->eval <= PT_ERROR) return GM_ERROR;
  }                             /* evaluate all prob./poss. trees */
  return gm->eval;              /* return sum of evaluation results */
}  /* gm_eval() */

/*--------------------------------------------------------------------*/

double gm_evalx (GRAMOD *gm, TABLE *table, double res[])
{                               /* --- evaluate a graphical model */
  int    i;                     /* loop variable */
  double tsr[3];                /* tuple specific results */
  double imposs;                /* number of impossible tuples */
  float  wgt;                   /* instantiation weight */

  assert(gm && table);          /* check the function arguments */
  if (gm_check(gm, -1, 0) != 0) /* if there are loops, */
    return -1;                  /* abort the function */
  res[0] = res[1] = res[2] = imposs = 0;
  for (i = tab_tplcnt(table); --i >= 0; ) {
    tpl_toas(tab_tpl(table,i)); /* traverse the table tuples */
    gm_execx(gm, tsr);          /* execute the graphical model */
    wgt = as_getwgt(gm->attset);/* get the instantiation weight */
    if (gm->type & GM_POSS) {   /* aggregate possibilities */
      res[0] += wgt *tsr[0];    /* (simply sum the weighted */
      res[1] += wgt *tsr[1];    /* evaluation results) */
      res[2] += wgt *tsr[2]; }
    else {                      /* aggregate probabilities */
      if (tsr[0] > 0) res[0] += wgt *log(tsr[0]);
      if (tsr[1] > 0) res[1] += wgt *log(tsr[1]);
      if (tsr[2] > 0) res[2] += wgt *log(tsr[2]);
    }                           /* (compute log-likelihood) */
    if (tsr[2] <= 0)            /* if the maximum prob./poss. is 0, */
      imposs += wgt;            /* count the tuple as impossible */
  }
  return imposs;                /* return the number of */
}  /* gm_evalx() */             /* impossible tuples */

/*--------------------------------------------------------------------*/

double gm_local (GRAMOD *gm, int lwise, double minimp)
{                               /* --- do local structure learning */
  int    i;                     /* loop variable */
  GMATT  *att;                  /* to traverse the attributes */

  assert(gm);                   /* check the function arguments */
  gm->eval = 0;                 /* init. the evaluation result */
  for (att = gm->atts +(i = gm->attcnt); --i >= 0; ) {
    if (!(--att)->ptree) continue;
    gm->eval += att->eval = pt_local(att->ptree, lwise, minimp);
    if (att->eval <= PT_ERROR) return GM_ERROR;
  }                             /* evaluate all prob./poss. trees */
  return gm->eval;              /* return sum of evaluation results */
}  /* gm_local() */

/*--------------------------------------------------------------------*/

const char* gm_sname (int search)
{                               /* --- get name of search method */
  if ((search < GM_NONE) || (search > GM_UNKNOWN))
    search = GM_UNKNOWN;        /* check the search method index */
  return snames[search];        /* return the search method name */
}  /* gm_sname() */