genealogy_t Class Reference

Encodes a genealogy. More...

#include <genealogy.h>

Inheritance diagram for genealogy_t:

Collaboration diagram for genealogy_t:

Public Member Functions
size_t	ndeme (void) const
	number of demes
size_t &	ndeme (void)
	number of demes
size_t	bytesize (void) const
	size of serialized binary form
	genealogy_t (double t0=R_NaReal, size_t ndeme=0)
	genealogy_t (raw_t *o)
	constructor from serialized binary form
	genealogy_t (SEXP o)
	constructor from RAW SEXP (containing binary serialization)
	genealogy_t (const genealogy_t &G)
	copy constructor
genealogy_t &	operator= (const genealogy_t &G)
	copy assignment operator
	genealogy_t (genealogy_t &&)=default
	move constructor
genealogy_t &	operator= (genealogy_t &&)=default
	move assignment operator
	~genealogy_t (void)
	destructor
slate_t &	time (void)
	view/set current time.
slate_t	time (void) const
	view current time.
slate_t &	timezero (void)
	view/set zero time.
slate_t	timezero (void) const
	get zero time.
void	lineage_count (double *tout, int *deme, int *ell, int *sat, int *etype) const
SEXP	lineage_count (void) const
	lineage count and saturation
void	gendat (double *tout, int *anc, int *lin, int *sat, int *type, int *deme, int *index, int *child) const
	genealogy information in list format
SEXP	gendat (void) const
	genealogy information in list format
size_t	nsample (void) const
	number of samples
size_t	nroot (void) const
	number of roots
string_t	yaml (string_t tab="") const
	human/machine-readable info
SEXP	structure (void) const
	R list description.
string_t	newick (bool extended=true) const
	put genealogy at current time into Newick format.
void	valid (void) const
	check the validity of the genealogy.
bool	check_genealogy_size (size_t grace=0) const
	check the size of the genealogy (to prevent memory exhaustion).
ball_t *	birth (ball_t *a, slate_t t, name_t d)
	birth into deme d
ball_t *	birth (node_t *p, name_t d)
	birth of second or subsequent sibling into deme d
void	death (ball_t *a, slate_t t)
	death
ball_t *	graft (slate_t t, name_t d)
	graft a new lineage into deme d
void	sample (ball_t *a, slate_t t)
	insert a sample node
void	sample_death (ball_t *a, slate_t t)
	insert a sample node and simultaneously terminate the lineage
void	migrate (ball_t *a, slate_t t, name_t d=0)
	movement into deme d
void	sample_migrate (ball_t *a, slate_t t, name_t d=0)
	insert a sample node and simultaneously migrate the lineage
genealogy_t &	prune (void)
	prune the tree (drop all black balls)
genealogy_t &	obscure (void)
	erase all deme information
void	curtail (slate_t tnew, slate_t troot)
genealogy_t &	operator+= (const genealogy_t &other)
	merges two genealogies, adjusting time, t0, and ndeme as needed
genealogy_t &	insert_zlb (void)
	insert zero-length branches for samples where needed
genealogy_t &	parse (const string_t &)
	Parse a Newick string and create the indicated genealogy.
void	time_rescale (slate_t scale, slate_t origin=0)
genealogy_t &	parse_cblv (const double *, const double *, int, double)
	parse a CBLV representation in the vectors x and y.
Public Member Functions inherited from nodeseq_t
	~nodeseq_t (void)
	destructor
size_t	bytesize (void) const
	size of serialized binary form
void	sort (void)
	order nodes in order of increasing time
pocket_t *	colored (color_t col) const
	Get all balls of a color.
size_t	ntime (slate_t t) const
	Number of distinct timepoints.
size_t	length (void) const
	Number of nodes in the sequence.
node_t *	position (int n)
	traverse to nth node, retrieve pointer
void	move (ball_t *b, node_t *p, node_t *q)
	move ball b from p to q
void	swap (ball_t *a, ball_t *b)
	swap balls a and b, wherever they lie
void	attach (node_t *p, node_t *q)
void	detach (node_t *p)
void	add (node_t *p, ball_t *a)
void	drop (ball_t *a)
void	destroy_node (node_t *p)
	remove a dead root node
void	weed (void)
	drop all dead roots
void	comb (void)
std::unordered_map< name_t, std::vector< node_t * > >	children_map (void) const
	map nodes onto vector of children
std::vector< node_t * >	sorted_roots (const std::unordered_map< name_t, slate_t > &height) const
	collect root nodes and sort them in order of decreasing subtree height
std::vector< node_t * >	ladderize (std::unordered_map< name_t, std::vector< node_t * > > &children) const
void	trace_lineages (void)
string_t	yaml (string_t tab="") const
	human/machine-readable info
SEXP	structure (void) const
	R list description.
string_t	newick (slate_t t, slate_t te, bool showdeme, bool extended) const
	put genealogy at time t into Newick format.

Private Member Functions
name_t	unique (void)
	get the next unique name
void	clean (void)
	clean up
node_t *	make_node (name_t d=undeme)
void	reuniqify (name_t shift)
	shifts names to avoid overlap
void	clip_zlb (void)
	clip out all zero-length branches
void	cap_tips (void)
void	cap_roots (void)
	roots are added at zero time if needed
node_t *	scan_branch_label (string_t::const_iterator b, string_t::const_iterator e, node_t *p, slate_t bl)
std::pair< std::vector< slate_t >, std::vector< slate_t > >	cblv (void) const

Private Attributes
name_t	_unique
	The next unique name.
slate_t	_t0
	The initial time.
slate_t	_time
	The current time.
size_t	_ndeme
	The number of demes (excluding the undeme).

Static Private Attributes
static const name_t	magic = 1123581321

Friends
raw_t *	operator>> (const genealogy_t &G, raw_t *o)
	binary serialization
raw_t *	operator>> (raw_t *o, genealogy_t &G)
	binary deserialization
SEXP	cblv (genealogy_t &A)
	return the CBLV representation as an R matrix

Additional Inherited Members
Static Public Member Functions inherited from nodeseq_t
static bool	compare (node_t *p, node_t *q)

Detailed Description

Encodes a genealogy.

A genealogy consists of a sequence of nodes and the current time.

Definition at line 20 of file genealogy.h.

Constructor & Destructor Documentation

genealogy_t() [1/5]

genealogy_t::genealogy_t ( double t0 = R_NaReal, size_t ndeme = 0 )

inline

basic constructor for genealogy class

• t0 = initial time
• ndeme = number of demes (excluding the undeme)

Definition at line 103 of file genealogy.h.

                                                       {
    clean();
    _time = _t0 = slate_t(t0);
    _ndeme = ndeme;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

genealogy_t() [2/5]

genealogy_t::genealogy_t ( raw_t * o )

inline

constructor from serialized binary form

Definition at line 109 of file genealogy.h.

                         {
    o >> *this;
  };

genealogy_t() [3/5]

genealogy_t::genealogy_t ( SEXP o )

inline

constructor from RAW SEXP (containing binary serialization)

Definition at line 113 of file genealogy.h.

                       {
    if (LENGTH(o)==0)
      err("in %s: cannot deserialize a NULL.",__func__);
    PROTECT(o = AS_RAW(o));
    RAW(o) >> *this;
    UNPROTECT(1);
  };

genealogy_t() [4/5]

genealogy_t::genealogy_t ( const genealogy_t & G )

inline

copy constructor

Definition at line 121 of file genealogy.h.

                                     {
    raw_t *o = new raw_t[G.bytesize()];
    G >> o;
    o >> *this;
    delete[] o;
  };

Here is the call graph for this function:

genealogy_t() [5/5]

genealogy_t::genealogy_t ( genealogy_t && )

default

move constructor

Here is the call graph for this function:

~genealogy_t()

genealogy_t::~genealogy_t ( void )

inline

destructor

Definition at line 141 of file genealogy.h.

                      {
    clean();
  };

Here is the call graph for this function:

Member Function Documentation

birth() [1/2]

ball_t * genealogy_t::birth ( ball_t * a, slate_t t, name_t d )

inline

birth into deme d

Definition at line 243 of file genealogy.h.

                                                 {
    time() = t;
    node_t *p = make_node(a->deme());
    ball_t *b = new ball_t (p,p->uniq,black,d);
    p->insert(b);
    p->slate = time();
    add(p,a);
    return b;
  };

Here is the call graph for this function:

birth() [2/2]

ball_t * genealogy_t::birth ( node_t * p, name_t d )

inline

birth of second or subsequent sibling into deme d

Definition at line 253 of file genealogy.h.

                                      {
    ball_t *b = new ball_t(p,unique(),black,d);
    p->insert(b);
    return b;
  };

Here is the call graph for this function:

bytesize()

size_t genealogy_t::bytesize ( void ) const

inline

size of serialized binary form

Definition at line 72 of file genealogy.h.

                               {
    return 3*sizeof(name_t) +
      2*sizeof(slate_t) + nodeseq_t::bytesize();
  };

Here is the call graph for this function:

Here is the caller graph for this function:

cap_roots()

void genealogy_t::cap_roots ( void )

inlineprivate

roots are added at zero time if needed

Definition at line 382 of file genealogy.h.

                        {
    node_nit j = begin();
    while (j != end()) {
      if ((*j)->is_root() && (*j)->slate > timezero()) {
        node_t *q = make_node();
        q->slate = timezero();
        attach(q,*j);
        push_front(q);
      }
      j++;
    }
    sort();
  };

Here is the call graph for this function:

cap_tips()

void genealogy_t::cap_tips ( void )

private

• tips without descendants are reclassified as samples.
• tips with black balls are swapped out.

Definition at line 8 of file parse.cc.

{
  for (node_t *p : *this) {
    if (p->empty()) {
      ball_t *b = new ball_t(p,p->uniq,blue,p->deme());
      p->insert(b);
    } else if (p->holds(black))
      swap(p->last_ball(),p->green_ball());
  }
}

Here is the call graph for this function:

Here is the caller graph for this function:

cblv()

std::pair< std::vector< slate_t >, std::vector< slate_t > > genealogy_t::cblv ( void ) const

private

return the CBLV representation in the vectors x and y:

• x[i] = branch length added by the i-th sample leaf
• y[i] = height above t0 of the i-th internal branching event

Definition at line 72 of file cblv.cc.

{
  auto children = children_map();
  auto roots  = ladderize(children);
  std::vector<slate_t> x, y;
  x.reserve(nsample());
  y.reserve(nsample());
  std::unordered_map<name_t, bool> memo;
  memo.reserve(size());
  for (node_t* p : *this) memo[p->uniq] = false;
  slate_t t0 = timezero();
  for (node_t* p : roots) {
    p->cblv(x, y, memo, children, t0);
    y.push_back(slate_t(0));
    memo[p->uniq] = true;
  }
  return {x, y};
}

Here is the call graph for this function:

check_genealogy_size()

bool genealogy_t::check_genealogy_size ( size_t grace = 0 ) const

inline

check the size of the genealogy (to prevent memory exhaustion).

Definition at line 215 of file genealogy.h.

                                                     {
    static size_t maxq = MEMORY_MAX/(sizeof(node_t)+2*sizeof(ball_t));
    bool ok = true;
    if (size() > maxq+grace) {
      err("maximum genealogy size exceeded!"); // #nocov
    } else if (size() > maxq) {
      ok = false;               // #nocov
    }
    return ok;
  };

Here is the caller graph for this function:

clean()

void genealogy_t::clean ( void )

inlineprivate

clean up

Definition at line 51 of file genealogy.h.

                    {
    _unique = 0;
    _ndeme = 0;
    _t0 = _time = R_NaReal;
  };

Here is the caller graph for this function:

clip_zlb()

void genealogy_t::clip_zlb ( void )

private

clip out all zero-length branches

Definition at line 21 of file parse.cc.

{
  for (node_t *p : *this) {
    if (!p->is_root() &&
        p->slate == p->parent()->slate &&
        p->deme() == p->parent()->deme()) {
      while (!p->empty()) {
        ball_t *b = p->last_ball();
        p->erase(b); p->parent()->insert(b);
      }
      detach(p);
    }
  }
}

Here is the call graph for this function:

Here is the caller graph for this function:

curtail()

void genealogy_t::curtail	(	slate_t	tnew,
		slate_t	troot )

curtail the genealogy by removing nodes with times later than tnew and/or earlier than troot

Definition at line 11 of file curtail.cc.

{
  if (tnew < troot) troot = tnew;
  if (!empty() && tnew < time()) {
    node_t *p = back();
    while (!empty() && p->slate > tnew) {
      ball_t *b;
      while (p->size() > 1) {
        b = p->last_ball();
        switch (b->color) {
        case black:
          p->erase(b); delete b;
          break;
        case green: case blue: // #nocov
          assert(0);           // #nocov
          break;               // #nocov
        }
      }
      b = p->last_ball();
      switch (b->color) {
      case blue:
        b->color = black;
      case black:
        b->deme() = p->deme();
        swap(b,p->green_ball());
      case green:
        destroy_node(p);
        break;
      }
      if (!empty()) p = back();
    }
  }
  time() = tnew;
  if (!empty() && troot > timezero()) {
    node_t *p = front();
    node_t *q;
    while (!empty() && p->slate < troot) {
      ball_t *b;
      assert(p->is_root());
      while (p->size() > 1) {
        b = p->last_ball();
        switch (b->color) {
        case blue:
          p->erase(b); delete b;
          break;
        case black:
          q = make_node(b->deme());
          q->slate = troot;
          move(b,p,q); push_back(q);
          break;
        case green:
          q = b->child();
          if (q == p) {
            b = p->first_ball();
            q = b->child();
          }
          if (q->slate < troot) {
            move(b,p,q);
          } else {
            node_t *pp = make_node(b->deme());
            pp->slate = troot;
            move(b,p,pp); push_back(pp);
          }
          break;
        }
      }
      destroy_node(p);
      if (!empty()) p = front();
    }
    sort();
  }
  if (troot > timezero()) timezero() = troot;
}

Here is the call graph for this function:

Here is the caller graph for this function:

death()

void genealogy_t::death ( ball_t * a, slate_t t )

inline

death

Definition at line 259 of file genealogy.h.

                                    {
    time() = t;
    drop(a);
  };

Here is the call graph for this function:

gendat() [1/2]

void genealogy_t::gendat	(	double *	tout,
		int *	anc,
		int *	lin,
		int *	sat,
		int *	type,
		int *	deme,
		int *	index,
		int *	child ) const

genealogy information in list format

Definition at line 8 of file gendat.cc.

                               {
  int m, n, k;
  node_it i, j;
  for (k = 0, n = 0, i = begin(); i != end(); i++, n++) {
    node_t *p = *i;
    assert(!p->holds(black)); // tree should be pruned first
    tout[n] = p->slate;
    deme[n] = p->deme();
    if (p->is_root()) {
      type[n] = 0;            // root node
    } else if (p->holds(blue)) {
      type[n] = 1;            // sample node
      deme[n] = p->ball(blue)->deme();
    } else {
      type[n] = 2;            // internal node
    }
    lin[n] = p->lineage();    // 0-based indexing
    sat[n] = p->nchildren();
    index[n] = k;
    k += sat[n];
    child[n] = NA_INTEGER;
    if (p->is_root()) {
      anc[n] = n;             // 0-based indexing
    } else {
      for (m = 0, j = begin(); j != i; j++, m++) {
        node_t *q = *j;
        if (p->parent()->uniq == q->uniq) {
          anc[n] = m;
          break;
        }
      }
    }
  }
  tout[n] = time();
  for (k = 0, n = 0, i = begin(); i != end(); i++, n++) {
    node_t *p = *i;
    j = i; j++;
    for (m = n+1; j != end(); m++, j++) {
      node_t *q = *j;
      if (p->uniq == q->parent()->uniq) {
        child[k++] = m;
      }
    }
  }
}

Here is the call graph for this function:

Here is the caller graph for this function:

gendat() [2/2]

SEXP genealogy_t::gendat

(

void

)

const

genealogy information in list format

Definition at line 59 of file gendat.cc.

             {
  SEXP t0, tout, anc, lin, sat, type, deme, index, child, ns, nr, nn;
  SEXP out, outn;
  size_t n = length();
  PROTECT(t0 = NEW_NUMERIC(1));
  PROTECT(tout = NEW_NUMERIC(n+1));
  PROTECT(type = NEW_INTEGER(n));
  PROTECT(deme = NEW_INTEGER(n));
  PROTECT(lin = NEW_INTEGER(n));
  PROTECT(sat = NEW_INTEGER(n));
  PROTECT(index = NEW_INTEGER(n));
  PROTECT(child = NEW_INTEGER(n));
  PROTECT(anc = NEW_INTEGER(n));
  PROTECT(ns = NEW_INTEGER(1));
  PROTECT(nr = NEW_INTEGER(1));
  PROTECT(nn = NEW_INTEGER(1));
  PROTECT(out = NEW_LIST(12));
  PROTECT(outn = NEW_CHARACTER(12));
  set_list_elem(out,outn,t0,"t0",0);
  set_list_elem(out,outn,tout,"nodetime",1);
  set_list_elem(out,outn,type,"nodetype",2);
  set_list_elem(out,outn,deme,"deme",3);
  set_list_elem(out,outn,lin,"lineage",4);
  set_list_elem(out,outn,sat,"saturation",5);
  set_list_elem(out,outn,index,"index",6);
  set_list_elem(out,outn,child,"child",7);
  set_list_elem(out,outn,anc,"ancestor",8);
  set_list_elem(out,outn,ns,"nsample",9);
  set_list_elem(out,outn,nr,"nroot",10);
  set_list_elem(out,outn,nn,"nnode",11);
  SET_NAMES(out,outn);
  gendat(REAL(tout),INTEGER(anc),INTEGER(lin),INTEGER(sat),
         INTEGER(type),INTEGER(deme),INTEGER(index),INTEGER(child));
  *REAL(t0) = double(timezero()); // zero-time
  *INTEGER(ns) = nsample();       // number of samples
  *INTEGER(nr) = nroot();         // number of roots
  *INTEGER(nn) = length();        // number of nodes
  UNPROTECT(14);
  return out;
}

Here is the call graph for this function:

Here is the caller graph for this function:

graft()

ball_t * genealogy_t::graft ( slate_t t, name_t d )

inline

graft a new lineage into deme d

Definition at line 264 of file genealogy.h.

                                      {
    time() = t;
    node_t *p = make_node();
    ball_t *b = new ball_t (p,p->uniq,black,d);
    p->insert(b);
    p->slate = timezero();
    push_front(p);
    return b;
  };

Here is the call graph for this function:

insert_zlb()

genealogy_t & genealogy_t::insert_zlb ( void )

inline

insert zero-length branches for samples where needed

Definition at line 359 of file genealogy.h.

                                 {
    for (node_t *p : *this) {
      if (p->holds(green) && p->holds(blue)) {
        assert(!p->holds(black)); // genealogy should have already been pruned
        ball_t *b = p->last_ball();
        assert(b->is(blue));
        node_t *q = make_node(p->deme());
        q->slate = p->slate;
        swap(q->green_ball(),b);
        push_back(q);
      }
    }
    sort();
    return *this;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

lineage_count() [1/2]

void genealogy_t::lineage_count	(	double *	tout,
		int *	deme,
		int *	ell,
		int *	sat,
		int *	etype ) const

lineage count, saturation, and event-type. types are:

• 0 = non-event
• -1 = root
• 1 = sample
• 2 = non-sample node
• 3 = end of interval

Definition at line 39 of file lineages.cc.

{
  size_t nd = ndeme()+1;
  slate_t tcur = timezero();
  for (size_t j = 0; j < nd; j++) {
    tout[j] = tcur;
    deme[j] = j;
    sat[j] = ell[j] = 0;
    etype[j] = 0;
  }
  for (const node_t *p : *this) {
    if (tcur < p->slate) {
      tout += nd; ell += nd; sat += nd;
      deme += nd; etype += nd;
      tcur = p->slate;
      for (size_t j = 0; j < nd; j++) {
        tout[j] = tcur;
        deme[j] = j;
        ell[j] = (ell-nd)[j];
        sat[j] = 0;
        etype[j] = 0;
      }
    }
    p->lineage_incr(ell,sat,etype);
  }
  tout += nd; ell += nd; sat += nd;
  deme += nd; etype += nd;
  tcur = time();
  for (size_t j = 0; j < nd; j++) {
    tout[j] = tcur;
    sat[j] = ell[j] = 0;
    deme[j] = j;
    etype[j] = 3;
  }
}

Here is the call graph for this function:

lineage_count() [2/2]

SEXP genealogy_t::lineage_count

(

void

)

const

lineage count and saturation

Definition at line 79 of file lineages.cc.

{
  SEXP tout, deme, ell, sat, etype, out, outn;
  int nt = ntime(timezero())+1;
  int nl = (ndeme()+1)*nt;
  PROTECT(tout = NEW_NUMERIC(nl));
  PROTECT(deme = NEW_INTEGER(nl));
  PROTECT(ell = NEW_INTEGER(nl));
  PROTECT(sat = NEW_INTEGER(nl));
  PROTECT(etype = NEW_INTEGER(nl));
  PROTECT(out = NEW_LIST(5));
  PROTECT(outn = NEW_CHARACTER(5));
  set_list_elem(out,outn,tout,"time",0);
  set_list_elem(out,outn,deme,"deme",1);
  set_list_elem(out,outn,ell,"lineages",2);
  set_list_elem(out,outn,sat,"saturation",3);
  set_list_elem(out,outn,etype,"event_type",4);
  SET_NAMES(out,outn);
  lineage_count(REAL(tout),INTEGER(deme),INTEGER(ell),
                INTEGER(sat),INTEGER(etype));
  UNPROTECT(7);
  return out;
}

Here is the call graph for this function:

Here is the caller graph for this function:

make_node()

node_t * genealogy_t::make_node ( name_t d = undeme )

inlineprivate

create a node holding its own green ball. this does not insert into the genealogy.

Definition at line 230 of file genealogy.h.

                                        {
    check_genealogy_size(0);
    name_t u = unique();
    node_t *p = new node_t(u,_time);
    ball_t *g = new ball_t(p,u,green,d);
    p->green_ball() = g;
    p->insert(g);
    return p;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

migrate()

void genealogy_t::migrate ( ball_t * a, slate_t t, name_t d = 0 )

inline

movement into deme d

Definition at line 293 of file genealogy.h.

                                                    {
    time() = t;
    node_t *p = make_node(a->deme());
    p->slate = time();
    add(p,a);
    a->deme() = d;
  };

Here is the call graph for this function:

ndeme() [1/2]

size_t & genealogy_t::ndeme ( void )

inline

number of demes

Definition at line 64 of file genealogy.h.

                       {
    return _ndeme;
  };

ndeme() [2/2]

size_t genealogy_t::ndeme ( void ) const

inline

number of demes

Definition at line 60 of file genealogy.h.

                            {
    return _ndeme;
  };

Here is the caller graph for this function:

newick()

string_t genealogy_t::newick

(

bool

extended = true

)

const

put genealogy at current time into Newick format.

Definition at line 106 of file newick.cc.

{
  return nodeseq_t::newick(time(),timezero(),(ndeme() > 0),extended);
}

Here is the call graph for this function:

nroot()

size_t genealogy_t::nroot ( void ) const

inline

number of roots

Definition at line 193 of file genealogy.h.

                            {
    size_t n = 0;
    for (const node_t *p : *this) {
      if (p->is_root()) n++;
    }
    return n;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

nsample()

size_t genealogy_t::nsample ( void ) const

inline

number of samples

Definition at line 184 of file genealogy.h.

                              {
    size_t n = 0;
    for (const node_t *p : *this) {
      if (p->holds(blue)) n++;
    }
    return n;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

obscure()

genealogy_t & genealogy_t::obscure ( void )

inline

erase all deme information

Definition at line 322 of file genealogy.h.

                              {
    // erase deme information from black balls.
    pocket_t *blacks = colored(black);
    while (!blacks->empty()) {
      ball_t *a = *(blacks->begin());
      a->deme() = undeme;
      blacks->erase(a);
    }
    delete blacks;
    // erase deme information from nodes.
    for (node_t *p : *this) {
      p->deme() = undeme;
    }
    // drop superfluous nodes (holding just one ball).
    comb();
    ndeme() = 0;
    return *this;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

operator+=()

genealogy_t & genealogy_t::operator+=

(

const genealogy_t &

other

)

merges two genealogies, adjusting time, t0, and ndeme as needed

merge two genealogies:

1. the node-sequences are merged;
2. the root time advances as necessary;
3. the current time retreats as necessary;
4. the unique-name stack advances as necessary.

Definition at line 29 of file sum.cc.

{
  genealogy_t G = other;
  slate_t t0 = timezero();
  slate_t t = time();
  t0 = (t0 < G.timezero()) ? t0 : G.timezero();
  t = (t > G.time()) ? t : G.time();
  G.reuniqify(_unique);
  merge(G,compare);
  timezero() = t0;
  time() = t;
  ndeme() = (ndeme() > G.ndeme()) ? ndeme() : G.ndeme();
  _unique = G._unique;
  return *this;
}

Here is the call graph for this function:

operator=() [1/2]

genealogy_t & genealogy_t::operator= ( const genealogy_t & G )

inline

copy assignment operator

Definition at line 128 of file genealogy.h.

                                                {
    clean();
    raw_t *o = new raw_t[G.bytesize()];
    G >> o;
    o >> *this;
    delete[] o;
    return *this;
  };

Here is the call graph for this function:

operator=() [2/2]

genealogy_t & genealogy_t::operator= ( genealogy_t && )

default

move assignment operator

Here is the call graph for this function:

parse()

genealogy_t & genealogy_t::parse

(

const string_t &

s

)

Parse a Newick string and create the indicated genealogy.

Definition at line 153 of file parse.cc.

{
  node_t *p = 0, *q;
  slate_t tf = timezero();
  slate_t bl = 0.0;
  string_t::const_reverse_iterator f = s.crend(), e = s.crbegin(), b = e;
  bool open = false;            // branch-string reading-frame open?
  int stack = 0, sqstack = 0;
  if (!s.empty() && *b != ';')
    err("in '%s': invalid Newick: no final semicolon.",__func__);
  while (b != f) {
    switch (*b) {
    case ';':                   // root
      if (stack != 0)
        err("in '%s': invalid Newick: unbalanced parentheses.",__func__);
      if (open) {
        q = scan_branch_label(b.base(),e.base(),p,bl);
        tf = (q->slate > tf) ? q->slate : tf;
        bl = 0.0;
      }
      p = make_node();
      p->slate = timezero();
      push_front(p);
      b++; e = b;
      open = true;
      bl = 0.0;
      break;
    case ')':                   // internal node
      if (open) {
        q = scan_branch_label(b.base(),e.base(),p,bl);
        if (q->holds(black))
          err("in '%s': 'type=extant' on internal node.",__func__);
        tf = (q->slate > tf) ? q->slate : tf;
        bl = 0.0;
        p = q;
      } else {
        err("in '%s': invalid Newick: missing comma or semicolon.",__func__);
      }
      b++; e = b;
      stack++;
      open = true;
      break;
    case '(':                   // tip node, eldest sister
      if (open) {
        q = scan_branch_label(b.base(),e.base(),p,bl);
        tf = (q->slate > tf) ? q->slate : tf;
        bl = 0.0;
      }
      p = p->parent();
      b++;
      e = b;
      stack--;
      open = false;
      break;
    case ',':                   // tip node, younger sister
      if (stack <= 0)
        err("in '%s': invalid Newick: misplaced comma or unbalanced parentheses.",__func__);
      if (open) {
        q = scan_branch_label(b.base(),e.base(),p,bl);
        tf = (q->slate > tf) ? q->slate : tf;
        bl = 0.0;
      }
      b++; e = b;
      open = true;
      bl = 0.0;
      break;
    case ']':                   // skip metadata
      sqstack++;
      while (b != f && sqstack > 0) {
        b++;
        if (*b == ']') sqstack++;
        if (*b == '[') sqstack--;
      }
      if (sqstack != 0)
        err("in '%s': invalid Newick: unbalanced square brackets.",__func__);
      else
        b++;
      break;
    case '[':
      err("in '%s': invalid Newick: unbalanced square brackets.",__func__);
      break;
    case ':':
      if (open) {
        bl = scan_branch_length(b.base(),e.base());
        b++; e = b;
      } else {
        err("in '%s': invalid Newick: misplaced colon.",__func__);
      }
      break;
    default:
      b++;
      break;
    }
  }
  if (stack != 0)
    err("in '%s': invalid Newick: unbalanced parentheses.",__func__);
  if (open) {
    q = scan_branch_label(b.base(),e.base(),p,bl);
    tf = (q->slate > tf) ? q->slate : tf;
    bl = 0.0;
  }
  time() = tf;
  sort(); cap_tips(); clip_zlb(); weed();
  return *this;
}

Here is the call graph for this function:

Here is the caller graph for this function:

parse_cblv()

genealogy_t & genealogy_t::parse_cblv	(	const double *	x,
		const double *	y,
		int	nin,
		double	tin )

parse a CBLV representation in the vectors x and y.

Definition at line 93 of file cblv.cc.

{
  if (nin <= 0) err("invalid CBLV");
  size_t n = size_t(nin);
  slate_t t0 = timezero();
  slate_t t = slate_t(tin);
  if (t < t0+slate_t(x[0]))
    err("invalid CBLV: x[0] = %lg > %lg = time-t0", x[0], t-t0);
  time() = t;
  node_t* p = 0;
  for (size_t k = 0; k < n; k++) {
    if (p == 0) {               // new root node at time t0
      p = make_node();
      p->slate = t0;
      push_back(p);
    }
    if (x[k] < 0)
      err("invalid CBLV: negative tip-edge length in position %zu", k+1);
    node_t* q = make_node();    // new tip node
    q->slate = p->slate + slate_t(x[k]);
    attach(p, q);
    push_back(q);
    t = t0 + slate_t(y[k]);     // new internal node time
    if (y[k] < 0)
      err("invalid CBLV: negative internal branch-time in position %zu", k+1);
    if (t > t0) {
      node_t* i = q->parent();  // points to p
      node_t* j = q;
      if (j->slate < t) err("invalid CBLV: node %zu cannot attach.", k+1);
      while (j != i && i->slate > t) {
        j = i;
        i = i->parent();
      }
      assert(j != i);
      // Create new internal node at time t
      node_t* node = make_node();
      node->slate = t;
      attach(i, node);
      move(j->green_ball(), i, node);
      push_back(node);
      p = node;
    } else {
      p = 0;
    }
  }
  if (p != 0) err("invalid CBLV: last value of y is nonzero.");
  sort(); cap_tips(); clip_zlb(); weed();
  return *this;
}

Here is the call graph for this function:

Here is the caller graph for this function:

prune()

genealogy_t & genealogy_t::prune ( void )

inline

prune the tree (drop all black balls)

Definition at line 311 of file genealogy.h.

                            {
    pocket_t *blacks = colored(black);
    while (!blacks->empty()) {
      ball_t *b = *(blacks->begin());
      blacks->erase(b);
      drop(b);
    }
    delete blacks;
    return *this;
  };

Here is the call graph for this function:

Here is the caller graph for this function:

reuniqify()

void genealogy_t::reuniqify ( name_t shift )

private

shifts names to avoid overlap

Definition at line 20 of file sum.cc.

{
  this->_unique += shift;
  for (node_t *p : *this) p->reuniqify(shift);
}

Here is the call graph for this function:

Here is the caller graph for this function:

sample()

void genealogy_t::sample ( ball_t * a, slate_t t )

inline

insert a sample node

Definition at line 274 of file genealogy.h.

                                     {
    time() = t;
    node_t *p = make_node(a->deme());
    ball_t *b = new ball_t (p,p->uniq,blue,a->deme());
    p->insert(b);
    p->slate = time();
    add(p,a);
  };

Here is the call graph for this function:

sample_death()

void genealogy_t::sample_death ( ball_t * a, slate_t t )

inline

insert a sample node and simultaneously terminate the lineage

Definition at line 283 of file genealogy.h.

                                           {
    time() = t;
    node_t *p = make_node(a->deme());
    ball_t *b = new ball_t (p,p->uniq,blue,a->deme());
    p->insert(b);
    p->slate = time();
    add(p,a);
    drop(a);
  };

Here is the call graph for this function:

sample_migrate()

void genealogy_t::sample_migrate ( ball_t * a, slate_t t, name_t d = 0 )

inline

insert a sample node and simultaneously migrate the lineage

Definition at line 301 of file genealogy.h.

                                                           {
    time() = t;
    node_t *p = make_node(a->deme());
    ball_t *b = new ball_t (p,p->uniq,blue,a->deme());
    p->insert(b);
    p->slate = time();
    add(p,a);
    a->deme() = d;
  };

Here is the call graph for this function:

scan_branch_label()

node_t * genealogy_t::scan_branch_label ( string_t::const_iterator b, string_t::const_iterator e, node_t * parent, slate_t bl )

private

Scan the branch-label string. This has format s[&&PhyloPOMP deme=d type=s]s

Scan the branch string. This has format s[&&PhyloPOMP deme=d type=s]s:f

Definition at line 124 of file parse.cc.

{
  name_t deme = 0;
  color_t col = green;
  if (b != e) {
    std::smatch m;
    if (std::regex_match(b,e,m,std::regex("^.*?\\[&&PhyloPOMP.+?deme=(\\w+).*?\\].*$")))
      deme = scan_name(m[1].str());
    if (std::regex_match(b,e,m,std::regex("^.*?\\[&&PhyloPOMP.+?type=(\\w+).*?\\].*$")))
      col = scan_color(m[1].str());
  }
  node_t *q = make_node(deme);
  if (col != green) {
    ball_t *b = new ball_t(q,q->uniq,col,deme);
    q->insert(b);
  }
  q->slate = bl+parent->slate;
  attach(parent,q);
  push_back(q);
  ndeme() = (ndeme() > q->deme()) ? ndeme() : q->deme();
  return q;
}

Here is the call graph for this function:

Here is the caller graph for this function:

structure()

SEXP genealogy_t::structure

(

void

)

const

R list description.

Definition at line 84 of file structure.cc.

{
  SEXP O, On, T0, Time, Nodes, Ndeme;
  PROTECT(O = NEW_LIST(4));
  PROTECT(On = NEW_CHARACTER(4));
  PROTECT(Time = NEW_NUMERIC(1));
  *REAL(Time) = double(time());
  PROTECT(T0 = NEW_NUMERIC(1));
  *REAL(T0) = double(timezero());
  PROTECT(Ndeme = NEW_INTEGER(1));
  *INTEGER(Ndeme) = int(ndeme());
  PROTECT(Nodes = nodeseq_t::structure());
  set_list_elem(O,On,Time,"time",0);
  set_list_elem(O,On,T0,"t0",1);
  set_list_elem(O,On,Ndeme,"ndeme",2);
  set_list_elem(O,On,Nodes,"nodes",3);
  SET_NAMES(O,On);
  UNPROTECT(6);
  return O;
}

Here is the call graph for this function:

time() [1/2]

slate_t & genealogy_t::time ( void )

inline

view/set current time.

Definition at line 146 of file genealogy.h.

                       {
    return _time;
  };

Here is the caller graph for this function:

time() [2/2]

slate_t genealogy_t::time ( void ) const

inline

view current time.

Definition at line 150 of file genealogy.h.

                            {
    return _time;
  };

time_rescale()

void genealogy_t::time_rescale ( slate_t scale, slate_t origin = 0 )

inline

Definition at line 405 of file genealogy.h.

                                                        {
    timezero() = scale*(timezero()-origin);
    for (node_t *p : *this ) {
      p->slate = scale*(p->slate-origin);
    }
    time() = scale*(time()-origin);
  };

Here is the call graph for this function:

Here is the caller graph for this function:

timezero() [1/2]

slate_t & genealogy_t::timezero ( void )

inline

view/set zero time.

Definition at line 154 of file genealogy.h.

                           {
    return _t0;
  };

Here is the caller graph for this function:

timezero() [2/2]

slate_t genealogy_t::timezero ( void ) const

inline

get zero time.

Definition at line 158 of file genealogy.h.

                                {
    return _t0;
  };

unique()

name_t genealogy_t::unique ( void )

inlineprivate

get the next unique name

Definition at line 44 of file genealogy.h.

                       {
    name_t u = _unique;
    _unique++;
    return u;
  };

Here is the caller graph for this function:

valid()

void genealogy_t::valid ( void ) const

inline

check the validity of the genealogy.

Definition at line 213 of file genealogy.h.

{};

yaml()

string_t genealogy_t::yaml

(

string_t

tab = ""

)

const

human/machine-readable info

Definition at line 63 of file yaml.cc.

{
  string_t o;
  string_t t = tab + "  ";
  o = tab + "t0: " + double2string(timezero()) + "\n"
    + tab + "time: " + double2string(time()) + "\n"
    + tab + "ndeme: " + std::to_string(ndeme()) + "\n"
    + tab + "nodes:\n" + nodeseq_t::yaml(tab);
  return o;
}

Here is the call graph for this function:

Here is the caller graph for this function:

Friends And Related Symbol Documentation

cblv

SEXP cblv ( genealogy_t & A )

friend

return the CBLV representation as an R matrix

Definition at line 150 of file cblv.cc.

{
  const char *colnames[] = {"tip","node"};
  double *x, *y;
  size_t i, n;
  SEXP S;
  std::pair<std::vector<slate_t>, std::vector<slate_t>> rep;
  rep = A.prune().obscure().insert_zlb().cblv();
  n = rep.first.size();
  PROTECT(S = make_matrix(n,2,colnames));
  x = REAL(S);
  y = REAL(S)+n;
  for (i = 0; i < n; i++) {
    x[i] = rep.first[i];
    y[i] = rep.second[i];
  }
  UNPROTECT(1);
  return S;
}

operator>> [1/2]

raw_t * operator>> ( const genealogy_t & G, raw_t * o )

friend

binary serialization

Definition at line 77 of file genealogy.h.

                                                            {
    name_t A[3]; A[0] = magic; A[1] = G._unique; A[2] = name_t(G.ndeme());
    slate_t B[2]; B[0] = G.timezero(); B[1] = G.time();
    memcpy(o,A,sizeof(A)); o += sizeof(A);
    memcpy(o,B,sizeof(B)); o += sizeof(B);
    return reinterpret_cast<const nodeseq_t&>(G) >> o;
  };

operator>> [2/2]

raw_t * operator>> ( raw_t * o, genealogy_t & G )

friend

binary deserialization

Definition at line 85 of file genealogy.h.

                                                      {
    G.clean();
    name_t A[3];
    slate_t B[2];
    memcpy(A,o,sizeof(A)); o += sizeof(A);
    memcpy(B,o,sizeof(B)); o += sizeof(B);
    if (A[0] != magic)
      err("in %s: corrupted genealogy serialization.",__func__);
    G._unique = A[1]; G.ndeme() = size_t(A[2]);
    G.timezero() = B[0]; G.time() = B[1];
    return o >> reinterpret_cast<nodeseq_t&>(G);
  };

Field Documentation

_ndeme

size_t genealogy_t::_ndeme

private

The number of demes (excluding the undeme).

Definition at line 37 of file genealogy.h.

_t0

slate_t genealogy_t::_t0

private

The initial time.

Definition at line 33 of file genealogy.h.

_time

slate_t genealogy_t::_time

private

The current time.

Definition at line 35 of file genealogy.h.

_unique

name_t genealogy_t::_unique

private

The next unique name.

Definition at line 31 of file genealogy.h.

magic

const name_t genealogy_t::magic = 1123581321

staticprivate

Definition at line 39 of file genealogy.h.

---

The documentation for this class was generated from the following files:

• src/genealogy.h
• src/cblv.cc
• src/curtail.cc
• src/gendat.cc
• src/lineages.cc
• src/newick.cc
• src/parse.cc
• src/structure.cc
• src/sum.cc
• src/yaml.cc