Bug report does not work anymore?

[richelbilderbeek]

Thanks for reporting a bug (or feature, we do not know yet) at nLTT issue 32.

I do not know what happened, but when I ran the (slightly adapted) code from your RPubs document today, I get the error message:

Error in rep(1, (length(wt) - 1)) : invalid 'times' argument

Can you reproduce the error? Just copy-paste the code below :-)

devtools::install_github("franciscorichter/dmea")
devtools::install_github("richelbilderbeek/nLTT", ref = "develop")

library(dmea)
library(ggplot2)
library(nLTT)

pp <- proc.time()
s <- sim_phyl(seed=2)
s2 <- phylo2p(drop.fossil(s$newick))
n_trees <- 5000
phylos <- vector('list',length=n_trees)

nLTT <- vector(length=n_trees)
gam <- vector(length=n_trees)
phylos2 <- vector('list',length=n_trees)
nLTT2 <- vector(length=n_trees)
gam2 <- vector(length=n_trees)
dt <- 0.01
for(i in 1:n_trees){
  r <- rec_tree(wt=s2$t)
  r1 <- p2phylo(r)
  phylos[[i]] <- read.tree(text=r1)
  nLTT[i] <- nLTTstat(s$newick, phylos[[i]], distance_method = "abs")
  gam[i] <- gammaStat(phylos[[i]])
  #
  r <- rec_tree(wt=s2$t,rec_method=2)
  r2 <- p2phylo(r)
  phylos2[[i]] <- read.tree(text=r2)
  nLTT2[i] <- nLTTstat(s$newick, phylos2[[i]], distance_method = "abs")
  gam2[i] <- gammaStat(phylos2[[i]])
}

[franciscorichter]

I will check it out

[richelbilderbeek]

Oops, typo:

devtools::install_github("richelbilderbeek/nLTT", dep = "develop")

should be

devtools::install_github("richelbilderbeek/nLTT", ref = "develop")

I updated this in the message above.

[franciscorichter]

this one should work (just change the for loop):

for(i in 1:n_trees){
  r <- rec_tree(wt=s2$wt)
  r1 <- p2phylo(r)
  phylos[[i]] <- r1
  nLTT[i] <- ltt_stat(s$newick, phylos[[i]])
  gam[i] <- gammaStat(phylos[[i]])
  r <- rec_tree(wt=s2$wt,rec_method=2)
  r2 <- p2phylo(r)
  phylos2[[i]] <- r2
  nLTT2[i] <- ltt_stat(s$newick, phylos2[[i]])
  gam2[i] <- gammaStat(phylos2[[i]])
}

I cannot use nLTTstat anymore though, that is because sometimes when you have non-ultrametric trees, the branching.times function gives negative branching times... that function does not work for non-ultrametric trees, so if you want nLTT working with non-ultrametric trees as well then you should use an alternative to branching.times... btw, ltt.plot.coords seems to work well.

> nLTT2[i] <- nLTTstat(s$newick, phylos2[[i]], distance_method = "abs")
Error in nLTT::nltt_diff(tree1, tree2, distance_method, ignore_stem) : 
  tree1 cannot have negative branching times

[richelbilderbeek]

Indeed this code does work. Going back to the issue.

devtools::install_github("franciscorichter/dmea")
devtools::install_github("richelbilderbeek/nLTT", ref = "develop")

library(dmea)
library(Hmisc)
library(ggplot2)
library(nLTT)

pp <- proc.time()
s <- sim_phyl(seed=2)
s2 <- phylo2p(drop.fossil(s$newick))
n_trees <- 5000
phylos <- vector('list',length=n_trees)

nLTT <- vector(length=n_trees)
gam <- vector(length=n_trees)
phylos2 <- vector('list',length=n_trees)
nLTT2 <- vector(length=n_trees)
gam2 <- vector(length=n_trees)
dt <- 0.01
for(i in 1:n_trees){
  r <- rec_tree(wt=s2$wt)
  r1 <- p2phylo(r)
  phylos[[i]] <- r1
  nLTT[i] <- ltt_stat(s$newick, phylos[[i]])
  gam[i] <- gammaStat(phylos[[i]])
  r <- rec_tree(wt=s2$wt,rec_method=2)
  r2 <- p2phylo(r)
  phylos2[[i]] <- r2
  nLTT2[i] <- ltt_stat(s$newick, phylos2[[i]])
  gam2[i] <- gammaStat(phylos2[[i]])
}

testit::assert(inherits(phylos[[1]], "phylo")) # Prerequisite for get_nltt_values

nltt_values <- get_nltt_values(phylos, dt = dt)
# Plot the phylognies, where the individual nLTT values are visible
qplot(t, nltt, data = nltt_values, geom = "point",
      ylim = c(0,1),
      main = "Average nLTT plot of phylogenies", size = I(0.1)
) + stat_summary(
  fun.data = "mean_cl_boot", color = "red", geom = "smooth"
)