gjuggler / ape Goto Github PK

Hi,

I am following a tutorial for phylogenetic reconstractions but when I load my fasta sequence alignment and I transform the object into phydat I cannot see the species names but it will return only numbers.

seq_alignment <- read.dna("ectos_alignment.txt", format="fasta")
alignment_phyDat <- phyDat(seq_alignment, type = "DNA", levels = NULL)

alignment_phyDat[[1]], it says [[1]], instead of the species names.

Do you maybe know how to fix this?
thanks a lot for your tutorial!

Hello there,

I am getting the following error when trying cophenetic(tree)

Error in double(nm * nm) : vector size cannot be NA
In addition: Warning message:
In nm * nm : NAs produced by integer overflow

My tree is unrooted, with 251141 tips and 251139 nodes.
My tip labels are like OTU_123, OTU_456, OTU_789 and my node labels like 0.123, 0.456, 0.789

My version of ape is 3.2

Any clue on that?

Thanks a lot,
Stelios

Code to reproduce:

library(ape)
D <- dist(1:500)
test <- ape::pcoa(D)

...and we get an error:

Error in array(STATS, dims[perm]) : 'dims' cannot be of length 0

What the problem is:

So this matrix only has one eigenvalue and its stupid to run pcoa. But still, theoretically it should work. We should get a single vector as output. I was trying to demonstrate this to a student.

pcoa() calculates k as the number of eigenvalues above epsilon, then subsets the eigenvectors to get columns 1:k. This fails in the case where k==1, because yall forgot to use drop=FALSE when you subset D.eig$vectors. Its necessary because by default in R when pulling one column out of a matrix, R decides that it isn't a matrix anymore and demotes it to a plain old vector. an example of this behavior:

M <- matrix(1:9, ncol=3)
M[ ,1]             # vector
M[ ,1, drop=FALSE] # a matrix

Easy Fix:

Current line in the function where error happens:

vectors <- sweep(D.eig$vectors[, 1:k], 2, sqrt(eig[1:k]), FUN = "*")

How it should look:

vectors <- sweep(D.eig$vectors[, 1:k, drop=FALSE], 2, sqrt(eig[1:k]), FUN = "*")

...at least that fixed it for me.

I wonder, how I can extract the newich format from a phylo object.

newicktree = '((((o1:1,o2:1)x1:1,(o3:1,o4:1)x2:1)x9:1,((o5:1,o6:1)x3:1,(o7:1,o8:1)x4:1)x10:1)x13:1,(((o9:1,o10:1)x5:1,(o11:1,o12:1)x6:1)x11:1,((o13:1,o14:1)x7:1,(o15:1,o16:1)x8:1)x12:1)x14:1)n1:1;'
cat(newicktree, file = "ex.tre", sep = "\n")
phy = read_tree(treefile = "ex.tre")

phy looks like

Phylogenetic tree with 16 tips and 15 internal nodes.

Tip labels:
o1, o2, o3, o4, o5, o6, ...
Node labels:
n1, x13, x9, x1, x2, x10, ...

Rooted; includes branch lengths.

and
str(phy)
returns

List of 6
$ edge : int [1:30, 1:2] 17 18 19 20 20 19 21 21 18 22 ...
$ Nnode : int 15
$ tip.label : chr [1:16] "o1" "o2" "o3" "o4" ...
$ edge.length: num [1:30] 1 1 1 1 1 1 1 1 1 1 ...
$ node.label : chr [1:15] "n1" "x13" "x9" "x1" ...
$ root.edge : num 1

• attr(*, "class")= chr "phylo"
• attr(*, "order")= chr "cladewise"

however, I would like to get back the actual newich format; if that is possible.