R/collapsed.intervals.R

   1 ## collapsed.intervals.R (2002-09-12)
   2
   3 ##   Collapsed coalescent intervals (e.g. for the skyline plot)
   4
   5 ## Copyright 2002 Korbinian Strimmer
   6
   7 ## This file is part of the R-package `ape'.
   8 ## See the file ../COPYING for licensing issues.
   9
  10 # construct collapsed intervals from coalescent intervals
  11 collapsed.intervals <- function(ci, epsilon=0.0)
  12 {
  13   if (class(ci) != "coalescentIntervals")
  14     stop("object \"ci\" is not of class \"coalescentIntervals\"")
  15
  16   sz <- ci$interval.length
  17   lsz <- length(sz)
  18   idx <- c <- 1:lsz
  19
  20   p <- 1
  21   w <- 0
  22
  23   # starting from tips collapes intervals
  24   # until total size is >= epsilon
  25   for (i in 1:lsz)
  26   {
  27     idx[[i]] <- p
  28     w <- w + sz[[i]]
  29     if (w >= epsilon)
  30     {
  31       p <- p+1
  32       w <- 0
  33     }
  34   }
  35
  36   # if last interval is smaller than epsilon merge
  37   # with second last interval
  38   lastInterval <- idx==p
  39   if ( sum(sz[lastInterval]) < epsilon )
  40   {
  41     p <- p-1
  42     idx[lastInterval] <- p
  43   }
  44
  45   obj <- list(
  46      lineages=ci$lineages,
  47      interval.length=ci$interval.length,
  48      collapsed.interval=idx, # collapsed intervals (via reference)
  49      interval.count=ci$interval.count,
  50      collapsed.interval.count = idx[[ci$interval.count]],
  51      total.depth =ci$total.depth,
  52      epsilon = epsilon
  53     )
  54   class(obj) <- "collapsedIntervals"
  55
  56   return(obj)
  57 }