Closed shajoezhu closed 3 years ago
If I change n3 = readTopology("(((HS03:6.12,TMS01:6.12):0.71,(NBLZG:5.98,#H18:0.0::0.21):0.85):0.47,((AYEM:7.13,(XYBX:5.98)#H18:1.15::0.79):0.1,ZGQCS:7.23):0.07,((JYYJFS:0.0)#H19:6.2::0.89,(YNEM:5.64,(DLYQCS:0.0,#H19:0.0::0.11):5.64):0.56):1.1);")
to
n3 = readTopology("(((HS03:6.12,TMS01:6.12):0.71,(NBLZG:5.98,#H18:0.0::0.21):0.85):0.47,((AYEM:7.13,(XYBX:5.98)#H18:1.15::0.79):0.1,ZGQCS:7.23):0.07,((JYYJFS:0.0)#H19:6.2::0.89,(YNEM:5.64,(DL:0.0,#H19:0.0::0.11):5.64):0.56):1.1);")
it would work ...
DLYQCS -> DL
Problem with the condition to determine adding_to_parent 2.
Hi @cecileane
$ hybrid-Lambda -spcu '(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((JYYJFS:0.0)H19#0.89:6.2,(YNEM:5.64,(DLYQCS:0.0,H19#0.89:0.0)I6:5.64)I7:0.56)I8:1.1)I9;' -num 10 -seed 42 -o n3
Default Kingman coalescent on all branches.
Default population size of 10000 on all branches.
Random seed: 42
Produced gene tree files:
n3_coal_unit
Seems to work now for branch 38. Do you want to give it a try?
Thank you so much @shajoezhu ! I tried all the networks above and the smaller networks in #36, and a few more networks, with edges rotated around a few nodes. Everything looks good now, in the sense that after making sure that my input networks were ultrametric, I didn't get any warnings from hybrid-Lambda about non-ultrametric gene trees. 🎉 Thank you!
I noticed that your fixes branched off an older commit, and their version file has an older version v0.6.1-beta. You might want to bump this up to v0.6.3-beta before merging your fixes from this branch into your main branch.
Thanks for following up on this, Joe!
On Sun, Nov 29, 2020 at 8:33 AM Joe Zhu notifications@github.com wrote:
All the examples above run well with the version on your new branch. So I tried on 9-taxon network that had lead me to the 3-taxon example above. Unfortunately I hit the same bug again... and a new one. (sorry for being the messenger of bad news!). Here is a reproducible example.
I used 4 networks: all with the same topology: 9 taxa, 2 reticulations, level 1, ultrametric.
- n2: everything works fine
- n3: same as n2, but with daughter edges "rotated" around a particular node
- n2s: same as n2, but with species names modified to be "s1" through "s9"
- n3s: same as n3, but with species names modified to be "s1" through "s9"
[image: nets_n2n3n2sn3s] https://user-images.githubusercontent.com/9857294/100530415-7bc1b280-31b7-11eb-8122-1177e5360830.png
Simulations look good on n2 but the bug appears on n3: see warnings below. (this is using your bug fix on commit 781eb35 https://github.com/hybridLambda/hybrid-Lambda/commit/781eb352bc7eacaec0800fb2b97f8c01910af708 )
% hybrid-Lambda -spcu '(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((YNEM:5.64,(DLYQCS:0.0,H19#0.11:0.0)I6:5.64)I7:0.56,(JYYJFS:0.0)H19#0.11:6.2)I8:1.1)I9;' -num 10 -seed 42 -o n2
Default Kingman coalescent on all branches. Default population size of 10000 on all branches. Random seed: 42 Produced gene tree files: n2_coal_unit
[image: gt_n2] https://user-images.githubusercontent.com/9857294/100530560-3bfbca80-31b9-11eb-9135-bed5745bb4f7.png
% hybrid-Lambda -spcu '(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((JYYJFS:0.0)H19#0.89:6.2,(YNEM:5.64,(DLYQCS:0.0,H19#0.89:0.0)I6:5.64)I7:0.56)I8:1.1)I9;' -num 10 -seed 42 -o n3
Default Kingman coalescent on all branches. Default population size of 10000 on all branches. WARNING! NOT ULTRAMETRIC!!! Random seed: 42 WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric Produced gene tree files: n3_coal_unit
[image: gt_n3] https://user-images.githubusercontent.com/9857294/100530570-52098b00-31b9-11eb-9490-30da8db1b1f3.png
The taxon JYYJFS (just below the "second" reticulation) is not aligned with the other taxa: the external branch to this taxon is either too short or too long by 6.2 exactly. This looks like the same bug as before, as if the parents edges at that reticulations had their branch lengths switched (6.2 and 0).
I then wanted to modify the tip labels to write this bug report. And then surprise: the behavior changed after changing the tip labels. On n2s, no warnings (as expected). But on n3s, the warnings disappeared and the gene trees were ultrametric!
%hybrid-Lambda -spcu '(((s1:6.12,s2:6.12)I1:0.71,(s3:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((s4:7.13,(s5:5.98)H18#0.21:1.15)I4:0.1,s6:7.23)I5:0.07,((s7:5.64,(s8:0.0,H19#0.11:0.0)I6:5.64)I7:0.56,(s9:0.0)H19#0.11:6.2)I8:1.1)I9;' -num 10 -seed 42 -o n2s
Default Kingman coalescent on all branches. Default population size of 10000 on all branches. Random seed: 42 Produced gene tree files: n2s_coal_unit
[image: gt_n2s] https://user-images.githubusercontent.com/9857294/100530576-5a61c600-31b9-11eb-8778-0e6ea4d533a0.png
% hybrid-Lambda -spcu '(((s1:6.12,s2:6.12)I1:0.71,(s3:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((s4:7.13,(s5:5.98)H18#0.21:1.15)I4:0.1,s6:7.23)I5:0.07,((s9:0.0)H19#0.89:6.2,(s7:5.64,(s8:0.0,H19#0.89:0.0)I6:5.64)I7:0.56)I8:1.1)I9;' -num 10 -seed 42 -o n3s
Default Kingman coalescent on all branches. Default population size of 10000 on all branches. Random seed: 42 Produced gene tree files: n3s_coal_unit
[image: gt_n3s] https://user-images.githubusercontent.com/9857294/100530587-75343a80-31b9-11eb-8322-27ce9387d7b8.png
Why the taxon labels affect the simulations?? It looks like whether the branch lengths are switched, or not, is affected by the branch labels.
For completeness, here is the julia code that I used to make the plots, with networks in extended newick format. (I also used julia to modify taxon names, rotate edges around a node, and get the network string for hybrid-Lambda).
using PhyloNetworks # v0.12.0using RCall # v0.13.10using PhyloPlots # add PhyloPlots#master
n2 = readTopology("(((HS03:6.12,TMS01:6.12):0.71,(NBLZG:5.98,#H18:0.0::0.21):0.85):0.47,((AYEM:7.13,(XYBX:5.98)#H18:1.15::0.79):0.1,ZGQCS:7.23):0.07,((YNEM:5.64,(DLYQCS:0.0,#H19:0.0::0.11):5.64):0.56,(JYYJFS:0.0)#H19:6.2::0.89):1.1);")hybridlambdaformat(n2) # to get format required by hybrid-Lambda: "(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((YNEM:5.64,(DLYQCS:0.0,H19#0.11:0.0)I6:5.64)I7:0.56,(JYYJFS:0.0)H19#0.11:6.2)I8:1.1)I9;" n3 = readTopology("(((HS03:6.12,TMS01:6.12):0.71,(NBLZG:5.98,#H18:0.0::0.21):0.85):0.47,((AYEM:7.13,(XYBX:5.98)#H18:1.15::0.79):0.1,ZGQCS:7.23):0.07,((JYYJFS:0.0)#H19:6.2::0.89,(YNEM:5.64,(DLYQCS:0.0,#H19:0.0::0.11):5.64):0.56):1.1);") n2s = readTopology("(((s1:6.12,s2:6.12):0.71,(s3:5.98,#H18:0.0::0.21):0.85):0.47,((s4:7.13,(s5:5.98)#H18:1.15::0.79):0.1,s6:7.23):0.07,((s7:5.64,(s8:0.0,#H19:0.0::0.11):5.64):0.56,(s9:0.0)#H19:6.2::0.89):1.1);") n3s = readTopology("(((s1:6.12,s2:6.12):0.71,(s3:5.98,#H18:0.0::0.21):0.85):0.47,((s4:7.13,(s5:5.98)#H18:1.15::0.79):0.1,s6:7.23):0.07,((s9:0.0)#H19:6.2::0.89,(s7:5.64,(s8:0.0,#H19:0.0::0.11):5.64):0.56):1.1);") R"layout"([1 2; 3 4])plot(n2,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n2", line=-2, adj=0.1);plot(n3,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n3", line=-2, adj=0.1);plot(n2s,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n2s", line=-2, adj=0.1);plot(n3s,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n3s", line=-2, adj=0.1); R"layout"([1 2]); R"par"(mar=[0,0,0,1]);for gt in readMultiTopology("n2_coal_unit")[1:2] plot(gt,:R; showEdgeLength=true, useEdgeLength=true)endR"mtext"("first 2 gene trees from n2", line=-1, adj=0.05, outer=true);# and similar for gene trees from n3, n2s and n3s.
Originally posted by @cecileane https://github.com/cecileane in #36 (comment) https://github.com/hybridLambda/hybrid-Lambda/issues/36#issuecomment-735321005
— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/hybridLambda/hybrid-Lambda/issues/38, or unsubscribe https://github.com/notifications/unsubscribe-auth/ACHBVJ3J2FLXDAJVFT2NIQTSSJSUHANCNFSM4UGSQVPQ .
Thanks for following up on this, Joe! … On Sun, Nov 29, 2020 at 8:33 AM Joe Zhu @.**> wrote: All the examples above run well with the version on your new branch. So I tried on 9-taxon network that had lead me to the 3-taxon example above. Unfortunately I hit the same bug again... and a new one. (sorry for being the messenger of bad news!). Here is a reproducible example. I used 4 networks: all with the same topology: 9 taxa, 2 reticulations, level 1, ultrametric. - n2: everything works fine - n3: same as n2, but with daughter edges "rotated" around a particular node - n2s: same as n2, but with species names modified to be "s1" through "s9" - n3s: same as n3, but with species names modified to be "s1" through "s9" [image: nets_n2n3n2sn3s] https://user-images.githubusercontent.com/9857294/100530415-7bc1b280-31b7-11eb-8122-1177e5360830.png Simulations look good on n2 but the bug appears on n3: see warnings below. (this is using your bug fix on commit 781eb35 <781eb35> ) % hybrid-Lambda -spcu '(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((YNEM:5.64,(DLYQCS:0.0,H19#0.11:0.0)I6:5.64)I7:0.56,(JYYJFS:0.0)H19#0.11:6.2)I8:1.1)I9;' -num 10 -seed 42 -o n2 Default Kingman coalescent on all branches. Default population size of 10000 on all branches. Random seed: 42 Produced gene tree files: n2_coal_unit [image: gt_n2] https://user-images.githubusercontent.com/9857294/100530560-3bfbca80-31b9-11eb-9135-bed5745bb4f7.png % hybrid-Lambda -spcu '(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((JYYJFS:0.0)H19#0.89:6.2,(YNEM:5.64,(DLYQCS:0.0,H19#0.89:0.0)I6:5.64)I7:0.56)I8:1.1)I9;' -num 10 -seed 42 -o n3 Default Kingman coalescent on all branches. Default population size of 10000 on all branches. WARNING! NOT ULTRAMETRIC!!! Random seed: 42 WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric WARNING: Gene tree is not ultrametric Produced gene tree files: n3_coal_unit [image: gt_n3] https://user-images.githubusercontent.com/9857294/100530570-52098b00-31b9-11eb-9490-30da8db1b1f3.png The taxon JYYJFS (just below the "second" reticulation) is not aligned with the other taxa: the external branch to this taxon is either too short or too long by 6.2 exactly. This looks like the same bug as before, as if the parents edges at that reticulations had their branch lengths switched (6.2 and 0). I then wanted to modify the tip labels to write this bug report. And then surprise: the behavior changed after changing the tip labels. On n2s, no warnings (as expected). But on n3s, the warnings disappeared and the gene trees were ultrametric! %hybrid-Lambda -spcu '(((s1:6.12,s2:6.12)I1:0.71,(s3:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((s4:7.13,(s5:5.98)H18#0.21:1.15)I4:0.1,s6:7.23)I5:0.07,((s7:5.64,(s8:0.0,H19#0.11:0.0)I6:5.64)I7:0.56,(s9:0.0)H19#0.11:6.2)I8:1.1)I9;' -num 10 -seed 42 -o n2s Default Kingman coalescent on all branches. Default population size of 10000 on all branches. Random seed: 42 Produced gene tree files: n2s_coal_unit [image: gt_n2s] https://user-images.githubusercontent.com/9857294/100530576-5a61c600-31b9-11eb-8778-0e6ea4d533a0.png % hybrid-Lambda -spcu '(((s1:6.12,s2:6.12)I1:0.71,(s3:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((s4:7.13,(s5:5.98)H18#0.21:1.15)I4:0.1,s6:7.23)I5:0.07,((s9:0.0)H19#0.89:6.2,(s7:5.64,(s8:0.0,H19#0.89:0.0)I6:5.64)I7:0.56)I8:1.1)I9;' -num 10 -seed 42 -o n3s Default Kingman coalescent on all branches. Default population size of 10000 on all branches. Random seed: 42 Produced gene tree files: n3s_coal_unit [image: gt_n3s] https://user-images.githubusercontent.com/9857294/100530587-75343a80-31b9-11eb-8322-27ce9387d7b8.png Why the taxon labels affect the simulations?? It looks like whether the branch lengths are switched, or not, is affected by the branch labels. For completeness, here is the julia code that I used to make the plots, with networks in extended newick format. (I also used julia to modify taxon names, rotate edges around a node, and get the network string for hybrid-Lambda). using PhyloNetworks # v0.12.0using RCall # v0.13.10using PhyloPlots # add PhyloPlots#master n2 = readTopology("(((HS03:6.12,TMS01:6.12):0.71,(NBLZG:5.98,#H18:0.0::0.21):0.85):0.47,((AYEM:7.13,(XYBX:5.98)#H18:1.15::0.79):0.1,ZGQCS:7.23):0.07,((YNEM:5.64,(DLYQCS:0.0,#H19:0.0::0.11):5.64):0.56,(JYYJFS:0.0)#H19:6.2::0.89):1.1);")hybridlambdaformat(n2) # to get format required by hybrid-Lambda: "(((HS03:6.12,TMS01:6.12)I1:0.71,(NBLZG:5.98,H18#0.21:0.0)I2:0.85)I3:0.47,((AYEM:7.13,(XYBX:5.98)H18#0.21:1.15)I4:0.1,ZGQCS:7.23)I5:0.07,((YNEM:5.64,(DLYQCS:0.0,H19#0.11:0.0)I6:5.64)I7:0.56,(JYYJFS:0.0)H19#0.11:6.2)I8:1.1)I9;" n3 = readTopology("(((HS03:6.12,TMS01:6.12):0.71,(NBLZG:5.98,#H18:0.0::0.21):0.85):0.47,((AYEM:7.13,(XYBX:5.98)#H18:1.15::0.79):0.1,ZGQCS:7.23):0.07,((JYYJFS:0.0)#H19:6.2::0.89,(YNEM:5.64,(DLYQCS:0.0,#H19:0.0::0.11):5.64):0.56):1.1);") n2s = readTopology("(((s1:6.12,s2:6.12):0.71,(s3:5.98,#H18:0.0::0.21):0.85):0.47,((s4:7.13,(s5:5.98)#H18:1.15::0.79):0.1,s6:7.23):0.07,((s7:5.64,(s8:0.0,#H19:0.0::0.11):5.64):0.56,(s9:0.0)#H19:6.2::0.89):1.1);") n3s = readTopology("(((s1:6.12,s2:6.12):0.71,(s3:5.98,#H18:0.0::0.21):0.85):0.47,((s4:7.13,(s5:5.98)#H18:1.15::0.79):0.1,s6:7.23):0.07,((s9:0.0)#H19:6.2::0.89,(s7:5.64,(s8:0.0,#H19:0.0::0.11):5.64):0.56):1.1);") R"layout"([1 2; 3 4])plot(n2,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n2", line=-2, adj=0.1);plot(n3,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n3", line=-2, adj=0.1);plot(n2s,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n2s", line=-2, adj=0.1);plot(n3s,:R; showEdgeLength=true, useEdgeLength=true);R"mtext"("n3s", line=-2, adj=0.1); R"layout"([1 2]); R"par"(mar=[0,0,0,1]);for gt in readMultiTopology("n2_coal_unit")[1:2] plot(gt,:R; showEdgeLength=true, useEdgeLength=true)endR"mtext"("first 2 gene trees from n2", line=-1, adj=0.05, outer=true);# and similar for gene trees from n3, n2s and n3s. Originally posted by @cecileane https://github.com/cecileane in #36 (comment) <#36 (comment)>* — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#38>, or unsubscribe https://github.com/notifications/unsubscribe-auth/ACHBVJ3J2FLXDAJVFT2NIQTSSJSUHANCNFSM4UGSQVPQ .
No problem @jamdeg ! :) Thanks!
just fyi: hybrid-Lambda v0.6.3-beta was built under a number of platforms here and is being used in this julia package as the underlying simulator to do a goodness-of-fit test. Thanks for the bug fix! :smiley:
That's a really nice resource! Many thanks! @cecileane
All the examples above run well with the version on your new branch. So I tried on 9-taxon network that had lead me to the 3-taxon example above. Unfortunately I hit the same bug again... and a new one. (sorry for being the messenger of bad news!). Here is a reproducible example.
I used 4 networks: all with the same topology: 9 taxa, 2 reticulations, level 1, ultrametric.
Simulations look good on n2 but the bug appears on n3: see warnings below. (this is using your bug fix on commit 781eb35)
The taxon JYYJFS (just below the "second" reticulation) is not aligned with the other taxa: the external branch to this taxon is either too short or too long by 6.2 exactly. This looks like the same bug as before, as if the parents edges at that reticulations had their branch lengths switched (6.2 and 0).
I then wanted to modify the tip labels to write this bug report. And then surprise: the behavior changed after changing the tip labels. On n2s, no warnings (as expected). But on n3s, the warnings disappeared and the gene trees were ultrametric!
Why the taxon labels affect the simulations?? It looks like whether the branch lengths are switched, or not, is affected by the branch labels.
For completeness, here is the julia code that I used to make the plots, with networks in extended newick format. (I also used julia to modify taxon names, rotate edges around a node, and get the network string for hybrid-Lambda).
Originally posted by @cecileane in https://github.com/hybridLambda/hybrid-Lambda/issues/36#issuecomment-735321005