dnapars.c

一个神经网络工具箱

          like = bstlike2 = -tempsum->sumsteps;
          there = forknode;
          mulf = false;
        }
      }
    } else if (-tempsum->sumsteps > like) {
      like = -tempsum->sumsteps;
      if (-tempsum->sumsteps > bestyet) {
        there = forknode;
        mulf = false;
      }
    }
    trydescendants(item, forknode, forknode->back, tempf, false);
  }
  q = forknode->next;
  while (q != forknode) {
    if (q->back != item) {
      memcpy(temp2->base, q->base, endsite*sizeof(long));
      memcpy(temp2->numsteps, q->numsteps, endsite*sizeof(long));
      memcpy(temp2->numnuc, q->numnuc, endsite*sizeof(nucarray));
      temp2->numdesc = q->numdesc - 1;
      multifillin(temp2, temprm, -1);
      if (!q->back->tip) {
        trydescendants(item, forknode, q->back, temp2, true);
      } else {
        sumnsteps(temp1, q->back, tempadd, 0, endsite);
        sumnsteps2(tempsum, temp1, temp2, 0, endsite, threshwt);
        if (lastrearr) {
          if (-tempsum->sumsteps > bstlike2) {
            multf = false;
            bestever = true;
            savelocrearr(item, forknode, q->back, tmp, tmp1, tmp2, tmp3,
                         tmprm, tmpadd, &root, maxtrees, &nextree, multf,
                         bestever, &saved, place, bestrees, treenode,
                         &grbg, zeros);
            if (saved) {
              like = bstlike2 = -tempsum->sumsteps;
              there = q->back;
              mulf = false;
            }
          }
        } else if ((-tempsum->sumsteps) > like) {
          like = -tempsum->sumsteps;
          if (-tempsum->sumsteps > bestyet) {
            there = q->back;
            mulf = false;
          }
        }
      }
    }
    q = q->next;
  }
} /* trylocal */


void trylocal2(node *item, node *forknode, node *other)
{
  /* rearranges below forknode, below descendants of forknode when
     there are more than 2 descendants, then unroots the back of
     forknode and rearranges on its descendants.  Used if forknode
     has binary descendants */
  node *q;
  boolean bestever=0, multf, saved, belowbetter, trysave;

  memcpy(tempf->base, other->base, endsite*sizeof(long));
  memcpy(tempf->numsteps, other->numsteps, endsite*sizeof(long));
  memcpy(tempf->oldbase, forknode->base, endsite*sizeof(long));
  memcpy(tempf->oldnumsteps, forknode->numsteps, endsite*sizeof(long));
  tempf->numdesc = other->numdesc;
  if (forknode->back)
    trydescendants(item, forknode, forknode->back, tempf, false);
  if (!other->tip) {
    memcpy(temp->base, other->base, endsite*sizeof(long));
    memcpy(temp->numsteps, other->numsteps, endsite*sizeof(long));
    memcpy(temp->numnuc, other->numnuc, endsite*sizeof(nucarray));
    temp->numdesc = other->numdesc + 1;
    multifillin(temp, tempadd, 1);
    if (forknode->back)
      sumnsteps2(tempsum, forknode->back, temp, 0, endsite, threshwt);
    else
      sumnsteps2(tempsum, NULL, temp, 0, endsite, threshwt);
    belowbetter = true;
    if (lastrearr) {
      if (-tempsum->sumsteps >= bstlike2) {
        belowbetter = false;
        bestever = false;
        multf = true;
        if (-tempsum->sumsteps > bstlike2)
          bestever = true;
        savelocrearr(item, forknode, other, tmp, tmp1, tmp2, tmp3, tmprm,
                       tmpadd, &root, maxtrees, &nextree, multf, bestever,
                       &saved, place, bestrees, treenode, &grbg, zeros);
        if (saved) {
          like = bstlike2 = -tempsum->sumsteps;
          there = other;
          mulf = true;
        }
      }
    } else if (-tempsum->sumsteps >= like) {
      there = other;
      mulf = true;
      like = -tempsum->sumsteps;
    }
    if (forknode->back) {
      memcpy(temprm->base, forknode->back->base, endsite*sizeof(long));
      memcpy(temprm->numsteps, forknode->back->numsteps, endsite*sizeof(long));
    } else {
      memcpy(temprm->base, zeros, endsite*sizeof(long));
      memcpy(temprm->numsteps, zeros, endsite*sizeof(long));
    }
    memcpy(temprm->oldbase, other->back->base, endsite*sizeof(long));
    memcpy(temprm->oldnumsteps, other->back->numsteps, endsite*sizeof(long));
    q = other->next;
    while (q != other) {
      memcpy(temp2->base, q->base, endsite*sizeof(long));
      memcpy(temp2->numsteps, q->numsteps, endsite*sizeof(long));
      memcpy(temp2->numnuc, q->numnuc, endsite*sizeof(nucarray));
      if (forknode->back) {
        temp2->numdesc = q->numdesc;
        multifillin(temp2, temprm, 0);
      } else {
        temp2->numdesc = q->numdesc - 1;
        multifillin(temp2, temprm, -1);
      }
      if (!q->back->tip)
        trydescendants(item, forknode, q->back, temp2, true);
      else {
        sumnsteps(temp1, q->back, tempadd, 0, endsite);
        sumnsteps2(tempsum, temp1, temp2, 0, endsite, threshwt);
        if (lastrearr) {
          if (-tempsum->sumsteps >= bstlike2) {
            trysave = false;
            multf = false;
            if (belowbetter) {
              bestever = false;
              trysave = true;
            }
            if (-tempsum->sumsteps > bstlike2) {
              bestever = true;
              trysave = true;
            }
            if (trysave) {
              savelocrearr(item, forknode, q->back, tmp, tmp1, tmp2, tmp3,
                           tmprm, tmpadd, &root, maxtrees, &nextree, multf,
                           bestever, &saved, place, bestrees, treenode,
                           &grbg, zeros);
              if (saved) {
                like = bstlike2 = -tempsum->sumsteps;
                there = q->back;
                mulf = false;
              }
            }
          }
        } else if (-tempsum->sumsteps > like) {
          like = -tempsum->sumsteps;
          if (-tempsum->sumsteps > bestyet) {
            there = q->back;
            mulf = false;
          }
        }
      }
      q = q->next;
    }
  }
} /* trylocal2 */


void tryrearr(node *p, boolean *success)
{
  /* evaluates one rearrangement of the tree.
     if the new tree has greater "likelihood" than the old
     one sets success = TRUE and keeps the new tree.
     otherwise, restores the old tree */
  node *forknode, *newfork, *other, *oldthere;
  double oldlike;
  boolean oldmulf;

  if (p->back == NULL)
    return;
  forknode = treenode[p->back->index - 1]; 
  if (!forknode->back && forknode->numdesc <= 2 && alltips(forknode, p))
    return;
  oldlike = bestyet;
  like = -10.0 * spp * chars;
  memcpy(tempadd->base, p->base, endsite*sizeof(long));
  memcpy(tempadd->numsteps, p->numsteps, endsite*sizeof(long));
  memcpy(tempadd->oldbase, zeros, endsite*sizeof(long));
  memcpy(tempadd->oldnumsteps, zeros, endsite*sizeof(long));
  if (forknode->numdesc > 2) {
    oldthere = there = forknode;
    oldmulf = mulf = true;
    trylocal(p, forknode);
  } else {
    findbelow(&other, p, forknode);
    oldthere = there = other;
    oldmulf = mulf = false;
    trylocal2(p, forknode, other);
  }
  if ((like <= oldlike) || (there == oldthere && mulf == oldmulf))
    return;
  recompute = true;
  re_move(p, &forknode, &root, recompute, treenode, &grbg, zeros);
  if (mulf)
    add(there, p, NULL, &root, recompute, treenode, &grbg, zeros);
  else {
    if (forknode->numdesc > 0)
      getnufork(&newfork, &grbg, treenode, zeros);
    else
      newfork = forknode;
    add(there, p, newfork, &root, recompute, treenode, &grbg, zeros);
  } 
  if (like > oldlike) {
    *success = true;
    bestyet = like;
  }
}  /* tryrearr */


void repreorder(node *p, boolean *success)
{
  /* traverses a binary tree, calling PROCEDURE tryrearr
     at a node before calling tryrearr at its descendants */
  node *q, *this;

  if (p == NULL)
    return;
  if (!p->visited) {
    tryrearr(p, success);
    p->visited = true;
  }
  if (!p->tip) {
    q = p;
    while (q->next != p) {
      this = q->next->back;
      repreorder(q->next->back,success);
      if (q->next->back == this)
        q = q->next;
    }
  }
}  /* repreorder */


void rearrange(node **r)
{
  /* traverses the tree (preorder), finding any local
     rearrangement which decreases the number of steps.
     if traversal succeeds in increasing the tree's
     "likelihood", PROCEDURE rearrange runs traversal again */
  boolean success=true;

  while (success) {
    success = false;
    clearvisited(treenode);
    repreorder(*r, &success);
  }
}  /* rearrange */


void describe()
{
  /* prints ancestors, steps and table of numbers of steps in
     each site */

  if (treeprint) {
    fprintf(outfile, "\nrequires a total of %10.3f\n", like / -10.0);
    fprintf(outfile, "\n  between      and       length\n");
    fprintf(outfile, "  -------      ---       ------\n");
    printbranchlengths(root);
  }
  if (stepbox)
    writesteps(chars, weights, oldweight, root);
  if (ancseq) {
    hypstates(chars, root, treenode, &garbage, basechar);
    putc('\n', outfile);
  }
  putc('\n', outfile);
  if (trout) {
    col = 0;
    treeout3(root, nextree, &col, root);
  }
}  /* describe */


void dnapars_coordinates(node *p, double lengthsum, long *tipy,
                        double *tipmax)
{
  /* establishes coordinates of nodes */
  node *q, *first, *last;
  double xx;

  if (p == NULL)
    return;
  if (p->tip) {
    p->xcoord = (long)(over * lengthsum + 0.5);
    p->ycoord = (*tipy);
    p->ymin = (*tipy);
    p->ymax = (*tipy);
    (*tipy) += down;
    if (lengthsum > (*tipmax))
      (*tipmax) = lengthsum;
    return;
  }
  q = p->next;
  do {
    xx = q->v;
    if (xx > 100.0)
      xx = 100.0;
    dnapars_coordinates(q->back, lengthsum + xx, tipy,tipmax);
    q = q->next;
  } while (p != q);
  first = p->next->back;
  q = p;
  while (q->next != p)
    q = q->next;
  last = q->back;
  p->xcoord = (long)(over * lengthsum + 0.5);
  if ((p == root) || count_sibs(p) > 2)
    p->ycoord = p->next->next->back->ycoord;
  else
    p->ycoord = (first->ycoord + last->ycoord) / 2;
  p->ymin = first->ymin;
  p->ymax = last->ymax;
}  /* dnapars_coordinates */


void dnapars_printree()
{
  /* prints out diagram of the tree2 */
  long tipy;
  double scale, tipmax;
  long i;
 
  if (!treeprint)
    return;
  putc('\n', outfile);
  tipy = 1;
  tipmax = 0.0;
  dnapars_coordinates(root, 0.0, &tipy, &tipmax);
  scale = 1.0 / (long)(tipmax + 1.000);
  for (i = 1; i <= (tipy - down); i++)
    drawline3(i, scale, root);
  putc('\n', outfile);
}  /* dnapars_printree */


void globrearrange()
{
  /* does global rearrangements */
  long j;
  double gotlike;
  boolean frommulti;
  node *item, *nufork;

  recompute = true;
  do {
    printf("   ");
    gotlike = bestlike;
    for (j = 0; j < nonodes; j++) {
      bestyet = -10.0 * spp * chars;
      if (j < spp)
        item = treenode[enterorder[j] -1];
      else 
        item = treenode[j];

      if ((item != root) && 
           ((j < spp) || ((j >= spp) && (item->numdesc > 0))) &&
           !((item->back->index == root->index) && (root->numdesc == 2)
              && alltips(root, item))) {
        re_move(item, &nufork, &root, recompute, treenode, &grbg, zeros);
        frommulti = (nufork->numdesc > 0);
        clearcollapse(treenode);
        there = root;
        memcpy(tempadd->base, item->base, endsite*sizeof(long));
        memcpy(tempadd->numsteps, item->numsteps, endsite*sizeof(long));
        memcpy(tempadd->oldbase, zeros, endsite*sizeof(long));
        memcpy(tempadd->oldnumsteps, zeros, endsite*sizeof(long));
        if (frommulti){
          oldnufork = nufork;
          getnufork(&nufork, &grbg, treenode, zeros);
        }
        addpreorder(root, item, nufork);
        if (frommulti)
          oldnufork = NULL;
        if (!mulf)
          add(there, item, nufork, &root, recompute, treenode, &grbg, zeros);
        else
          add(there, item, NULL, &root, recompute, treenode, &grbg, zeros);
      }
      if (progress) {
        putchar('.');
        fflush(stdout);
      }
    }
    if (progress) {
      putchar('\n');
#ifdef WIN32
      phyFillScreenColor();
#endif
    }
  } while (bestlike > gotlike);
} /* globrearrange */

void load_tree(long treei)
{ /* restores a tree from bestrees */
  long j, nextnode;
  boolean recompute = false;
  node *dummy;

  for (j = spp - 1; j >= 1; j--)
    re_move(treenode[j], &dummy, &root, recompute, treenode, &grbg,
              zeros);
  root = treenode[0];
  recompute = true;
  add(treenode[0], treenode[1], treenode[spp], &root, recompute,
    treenode, &grbg, zeros);
  nextnode = spp + 2;
  for (j = 3; j <= spp; j++) {
    if (bestrees[treei].btree[j - 1] > 0)
      add(treenode[bestrees[treei].btree[j - 1] - 1], treenode[j - 1],
            treenode[nextnode++ - 1], &root, recompute, treenode, &grbg,
            zeros);
    else
      add(treenode[treenode[-bestrees[treei].btree[j-1]-1]->back->index-1],
            treenode[j - 1], NULL, &root, recompute, treenode, &grbg,
            zeros);
  }
}

void grandrearr()
{
  /* calls global rearrangement on best trees */
  long treei;