xavierdidelot
commented
3 years ago Hi Feng,
From the output res=bactdate(...) you can calculate the local rate of evolution of each branch. This is equal to the number of substitutions on a branch (stored in res$tree$subs) divided by the estimated duration of that branch (stored in res$tree$edge.length). You could see if this rate seems to be consistently different in your lineage A-E as opposed to the remainder of the tree. These local rates are computed and displayed when you run the command plot(res,'scatter') and you could look at the code for this function to see how it is done. Cf from line 64 at https://github.com/xavierdidelot/BactDating/blob/master/R/methods.R
Best wishes, Xavier
Sorry it's me again, and I'm wondering if there is a way or if it is appropriate to estimate (actually just to work out) mu, sigma and alpha of a specific branch/lineage/subtree/clade from the results obtained for the whole tree?
e.g. we have 10 strains (A to J) , and the whole analysis worked fine (effective size all greater than 500). Also it showed pstrict=0 indicating strains evolved in different rates. We knew that strains A, B, C, D and E are clustered together on the tree (pairwised SNP<20) sharing a more recent common ancestor, which was distant (SNP>500) from the rest strains. Bactdating results analyzed from 10 strains showed mu=15[7;22], sigma=17[6;25], alpha=6[4;9].
I'm wondering if it is feasible and appropriate to work out the mu, sigma for the subtree representing strains A, B, C, D, E and their LCA from the existing result?
Starting a new analysis on these five strains only did not help (effective size remained low even ran for 1e8 in length), since they were isolated within a short period (3 months), also I'm worried that the negative dependency between substitution rate estimates and sampling time mentioned in the paper would give an incorrect result.
Sorry to bother you again and many thanks.