hybridLambda / hybrid-Lambda

Hybrid-Lambda is a software package that can simulate gene trees within a rooted species network or a rooted species tree under the Kingman's coalescent or Lambda coalescent process.
http://hybridlambda.github.io/
GNU General Public License v3.0
6 stars 4 forks source link

Level-2 networks sometimes cause segfault or warnings about ultrametricity #42

Open gaballench opened 3 years ago

gaballench commented 3 years ago

Hello,

We have been simulating phylogenetic networks with SiPhyNetwork in R, and are trying to use hybrid-Lambda to simulate gene trees on these. However, when using some of our level-2 networks as input for hybrid-Lambda we get (1) segmentation faults, (2) a warning on lack of ultrametricity, or (3) both at the same time.

We have used the approach below for network generation:

### script.R

library(SiPhyNetwork)
library(ape)
set.seed(10)

# simulate networks under the birth-death hybridization process
numbsim = 500
n1 <- 15
bdNS <- SiPhyNetwork::sim.bdh.taxa.ssa(
    n = n1,
    numbsim = numbsim,
    lambda = 0.9,
    mu = 0,
    nu = 0.03,
    hybprops = c(0.5, 0.25, 0.25),
    hyb.inher.fxn = SiPhyNetwork::make.beta.draw(1, 1),
    frac = 1,
    mrca = FALSE,
    complete = TRUE,
    stochsampling = FALSE,
    hyb.rate.fxn = NULL,
    trait.model = NULL
    )

# here some code testing whether the networks are level-1 or not
out <- # indices to networks that are not level-1

# check that all non-level-1 networks are ultrametric using ape
if (sum(sapply(X = bdNS[out], FUN = ape::is.ultrametric) == length(out))) {
    cat("All networks to be written to nets_not_level1.extnex are ultrametric\n")
}

# write all the networks above level-1 to a file (all of them are)
SiPhyNetwork::write.net(net = bdNS[out], file = "nets_not_level1.extnex")

The content of net_not_level1.extnex is found attached to this issue.

Then, we use PhyloNetworks for conversion of format between extended newick and hybridlambda and write them to the file nets_not_level1.hybridlambda (also available in the attachment):

### script.jl

using PhyloNetworks
nets_not_level1 = readMultiTopology("nets_not_level1.extnex")

# bulk format conversion of networks read in nets_not_level1
open("nets_not_level1.hybridlambda", "w") do io    
    for i = 1:length(nets_not_level1)
        hyblambFormatted = hybridlambdaformat(nets_not_level1[i])
        write(io, "$hyblambFormatted\n")
        flush(io)
    end
end

And finally iterate over each network to simulate a set of three loci with some arbitrary parameters:

### script.sh
#!/usr/bin/bash

# run the r script
Rscript script.r

# run the julia script
julia script.jl

# iterate over the content of nets_not_level1.hybridlambda and
# run hybrid-Lambda on each line in order to check issues
# with ultrametricity
while read line
do
    let LNUM+=1
    echo "Processing network in line $LNUM:"
    hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
done < nets_not_level1.hybridlambda

The above is run redirecting std out and err to a file (also attached):

./script.sh > outerr.txt 2>&1

The whole procedure returns 63 networks stored in nets_not_level1.hybridlambda, some of which show the problems already mentioned:

We would like to know if we are missing something on the way hybrid-Lambda should be run, e.g., whether level-2 networks are not guaranteed to work.

Also, we would like to understand better the warning about lack of ultrametricity, and whether it is safe to ignore it due to the fact that the simulated network is by definition ultrametric, and tested positive with ape.

Finally, we would like to understand why the segfault is happening sometimes, and whether it can be solved somehow.

Thanks in advance for your time and feedback.

Best wishes,

Gustavo (@gaballench), Carlos (@acosta187), and Claudia (@crsl4). nets_not_level1_extnew.txt nets_not_level1_hybridlambda.txt outerr.txt

shajoezhu commented 3 years ago

Hi Gustovo

Thanks for the email. I am a bit busy this week. Will take a closer look at this this weekend

Best Joe

On Wed, 14 Jul 2021, 22:10 Gustavo A. Ballen, @.***> wrote:

Hello,

We have been simulating phylogenetic networks with SiPhyNetwork in R, and are trying to use hybrid-Lambda to simulate gene trees on these. However, when using some of our level-2 networks as input for hybrid-Lambda we get (1) segmentation faults, (2) a warning on lack of ultrametricity, or (3) both at the same time.

We have used the approach below for network generation:

script.R

library(SiPhyNetwork) library(ape) set.seed(10)

simulate networks under the birth-death hybridization processnumbsim = 500n1 <- 15bdNS <- SiPhyNetwork::sim.bdh.taxa.ssa(

n = n1,
numbsim = numbsim,
lambda = 0.9,
mu = 0,
nu = 0.03,
hybprops = c(0.5, 0.25, 0.25),
hyb.inher.fxn = SiPhyNetwork::make.beta.draw(1, 1),
frac = 1,
mrca = FALSE,
complete = TRUE,
stochsampling = FALSE,
hyb.rate.fxn = NULL,
trait.model = NULL
)

here some code testing whether the networks are level-1 or notout <- # indices to networks that are not level-1

check that all non-level-1 networks are ultrametric using apeif (sum(sapply(X = bdNS[out], FUN = ape::is.ultrametric) == length(out))) {

cat("All networks to be written to nets_not_level1.extnex are ultrametric\n")

}

write all the networks above level-1 to a file (all of them are)SiPhyNetwork::write.net(net = bdNS[out], file = "nets_not_level1.extnex")

The content of net_not_level1.extnex is found attached to this issue.

Then, we use PhyloNetworks for conversion of format between extended newick and hybridlambda and write them to the file nets_not_level1.hybridlambda (also available in the attachment):

script.jl

using PhyloNetworks nets_not_level1 = readMultiTopology("nets_not_level1.extnex")

bulk format conversion of networks read in nets_not_level1open("nets_not_level1.hybridlambda", "w") do io

for i = 1:length(nets_not_level1)
    hyblambFormatted = hybridlambdaformat(nets_not_level1[i])
    write(io, "$hyblambFormatted\n")
    flush(io)
endend

And finally iterate over each network to simulate a set of three loci with some arbitrary parameters:

script.sh#!/usr/bin/bash

run the r script

Rscript script.r

run the julia script

julia script.jl

iterate over the content of nets_not_level1.hybridlambda and# run hybrid-Lambda on each line in order to check issues# with ultrametricitywhile read linedo

let LNUM+=1
echo "Processing network in line $LNUM:"
hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2done < nets_not_level1.hybridlambda

The above is run redirecting std out and err to a file (also attached):

./script.sh > outerr.txt 2>&1

The whole procedure returns 63 networks stored in nets_not_level1.hybridlambda, some of which show the problems already mentioned:

  • segfault networks: 5, 6, 14, 16, 20, 27, 28, 29, 32, 36, 40, 43, 44, 50, 56, 59, 63
  • ultrametric and segfault networks: 34, 52, 57, 58, 61

We would like to know if we are missing something on the way hybrid-Lambda should be run, e.g., whether level-2 networks are not guaranteed to work.

Also, we would like to understand better the warning about lack of ultrametricity, and whether it is safe to ignore it due to the fact that the simulated network is by definition ultrametric, and tested positive with ape.

Finally, we would like to understand why the segfault is happening sometimes, and whether it can be solved somehow.

Thanks in advance for your time and feedback.

Best wishes,

Gustavo @. https://github.com/gaballench), Carlos @. https://github.com/acosta187), and Claudia @.*** https://github.com/crsl4). nets_not_level1_extnew.txt https://github.com/hybridLambda/hybrid-Lambda/files/6818940/nets_not_level1_extnew.txt nets_not_level1_hybridlambda.txt https://github.com/hybridLambda/hybrid-Lambda/files/6818941/nets_not_level1_hybridlambda.txt outerr.txt https://github.com/hybridLambda/hybrid-Lambda/files/6818942/outerr.txt

— 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/42, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4FP7INMIUERZVWXAG5UKLTXX4KRANCNFSM5AMHOUIA .

shajoezhu commented 3 years ago

Dear @gaballench

Thanks for raising this issue. I think a bug was introduce in the latest release version, can you try version v0.6.2-beta and see how it goes

for network 5, it seems to work (((((t15:1.184583142)H24#0.7410187677:0.8005621307,H24#0.7410187677:0.8005621307)I1:0.01518882639,H23#0.4876111529:1.196276338)I4:0.07817391523,(((t1:0.08702091214,(t4:0.08702091214)H44#0.9538987672:0.0)I5:0.1581596425,((t2:0.1427137787,t13:0.1427137787)I6:0.102466776)H41#0.7859932417:0.0)I7:0.9739095645,(((t6:0.3484408972,(H38#0.204212385:0.06880064649,H35#0.31131891:0.0)I8:0.003584901564)I9:0.3234586189,(t3:0.3002319667,((t9:0.2760553491,(t8:0.2760553491)H38#0.204212385:0.0)I10:0.01557614853,t10:0.2916314976)I11:0.008600469108)I12:0.3716675493)I13:0.1321582449,(((t7:0.3448559956)H35#0.31131891:0.01451935831,((t12:0.08702091214,H44#0.9538987672:0.0)I2:0.2723544418)H34#0.9136001719:0.0)I3:0.444682407)H23#0.4876111529:0.0)I14:0.4150323583)I15:0.8594178952)I16:0.2879563441,((t14:0.3593753539,H34#0.9136001719:0.0)I17:0.5942466018,(t11:0.6738791754,(t16:0.2451805547,H41#0.7859932417:0.0)I18:0.4286986207)I19:0.2797427803)I20:1.412842403)I21;

network5-1

For network 34, (((H32#0.2870094643:0.08858339566,t7:0.9057121006)I2:0.3863821978,((((t17:0.3577009352,(t8:0.1127078449,t11:0.1127078449)I3:0.2449930903)I4:0.4663350034,t1:0.8240359386)I5:0.3833084735,t9:1.207344412)I6:0.08474988634)H25#0.6319662712:0.0)I7:1.44765081,((t5:1.601079997,(t3:1.497177316,(((t12:0.6411242122,t14:0.6411242122)I8:0.2647376045,((t6:0.005253654298,t4:0.005253654298)I9:0.569506896,(t2:0.4245830973)H22#0.9227559222:0.150177453)I10:0.3311012664)I11:0.3862324817,H25#0.6319662712:0.0)I12:0.2050830174)I13:0.1039026816)I14:0.6035924175,(H22#0.9227559222:1.484776982,((t13:0.05635283412,t10:0.05635283412)I1:0.7607758709)H32#0.2870094643:1.092231374)I15:0.2953123361)I16:0.5350726934)I17;

network34-1

Could you check the network please,

I will try and find some time and look into this bug. Thanks!

gaballench commented 3 years ago

Hello Joe,

Thanks for taking a look into this. I checked out into release v0.6.2-beta (whose version returns v0.6.1-beta after building, though):

hybrid-Lambda/src$ git checkout v0.6.2-beta
HEAD is now at de99ef3 update
hybrid-Lambda/src$ git status
HEAD detached at v0.6.2-beta
hybrid-Lambda/src$ make
hybrid-Lambda/src$ ./hybrid-Lambda
hybrid-Lambda v0.6.1-beta
Usage:
...

I ran again the script and all previous instances of segfaults are now gone. However, some of these are instead reported as not ultrametric now, or a couple instances which were initially reported as ultrametric are no longer so:

diff --git a/outerr.txt b/outerr.txt
index e076b82..5e5f413 100644
--- a/BashScript_Network_NotLevelOne/outerr.txt
+++ b/BashScript_Network_NotLevelOne/outerr.txt
@@ -42,11 +42,19 @@ level2_2_coal_unit
 Processing network in line 5:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8781 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 6:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8784 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+WARNING! NOT ULTRAMETRIC!!!
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 7:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -90,7 +98,11 @@ Random seed: 2
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 14:
-./script.sh: line 23:  8802 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 15:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -98,11 +110,19 @@ Random seed: 2
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 16:
-./script.sh: line 23:  8812 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 17:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
 Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 18:
@@ -118,7 +138,15 @@ Random seed: 2
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 20:
-./script.sh: line 23:  8824 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 21:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -140,7 +168,11 @@ level2_2_coal_unit
 Processing network in line 24:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
 Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 25:
@@ -156,13 +188,31 @@ Random seed: 2
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 27:
-./script.sh: line 23:  8838 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
Processing network in line 28:
-./script.sh: line 23:  8844 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 29:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8852 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+WARNING! NOT ULTRAMETRIC!!!
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 30:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -176,7 +226,11 @@ Random seed: 2
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 32:
-./script.sh: line 23:  8865 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 33:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -187,7 +241,12 @@ Processing network in line 34:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
 WARNING! NOT ULTRAMETRIC!!!
-./script.sh: line 23:  8874 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 35:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -197,7 +256,9 @@ level2_2_coal_unit
 Processing network in line 36:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8883 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 37:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -219,7 +280,9 @@ level2_2_coal_unit
 Processing network in line 40:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8894 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 41:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -235,9 +298,15 @@ level2_2_coal_unit
 Processing network in line 43:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8904 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 44:
-./script.sh: line 23:  8909 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 45:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -247,7 +316,11 @@ level2_2_coal_unit
 Processing network in line 46:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
 Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 47:
@@ -269,7 +342,15 @@ Random seed: 2
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 50:
-./script.sh: line 23:  8925 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 51:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -279,8 +360,9 @@ level2_2_coal_unit
 Processing network in line 52:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-WARNING! NOT ULTRAMETRIC!!!
-./script.sh: line 23:  8934 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 53:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -290,7 +372,11 @@ level2_2_coal_unit
 Processing network in line 54:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
 Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 55:
@@ -302,19 +388,35 @@ level2_2_coal_unit
 Processing network in line 56:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
-./script.sh: line 23:  8947 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 57:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
 WARNING! NOT ULTRAMETRIC!!!
-./script.sh: line 23:  8950 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 58:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
 WARNING! NOT ULTRAMETRIC!!!
-./script.sh: line 23:  8958 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 59:
-./script.sh: line 23:  8966 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 60:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
@@ -325,15 +427,28 @@ Processing network in line 61:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
 WARNING! NOT ULTRAMETRIC!!!
-./script.sh: line 23:  8978 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 62:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.
+WARNING! NOT ULTRAMETRIC!!!
 Random seed: 2
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
+WARNING: Gene tree is not ultrametric
 Produced gene tree files:
 level2_2_coal_unit
 Processing network in line 63:
-./script.sh: line 23:  8988 Segmentation fault      (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2
+Default Kingman coalescent on all branches.
+Default population size of 10000 on all branches.
+Random seed: 2
+Produced gene tree files:
+level2_2_coal_unit
 Processing network in line 64:
 Default Kingman coalescent on all branches.
 Default population size of 10000 on all branches.

Hope these results help to solve the bug. We are still trying to understand what triggers the warning, and whether it is safe to use hybrid-Lambda in the presence of such lack of ultrametricity. We may discard these not-ultrametric networks if it is better to do so, but given that these were simulated and in principle should all be ultrametric, we are still puzzled by this warning.

By looking at the code of net.cpp (L323) we find that the function check_isUltrametric uses 0.000001 in square difference of two paths as a threshold for considering them different, and hence the network as not ultrametric. That makes absolute differences in path length greater than 0.001 to fail the test. Is there any reason for choosing such particular threshold and not perhaps a bigger value? That could help us understand whether the simulation of the networks might produce misleading results.

Thanks again,

Gustavo

shajoezhu commented 3 years ago

Hi Gustavo

I think there are two issues here.

  1. There is a known bug for parsing network, hence there was the fix and newer version, but it seems causing more problems for you in the new version. I will look into it and get a patch version for you. I am quite busy lately with work, perhaps if you can use the ones are working first, and we can workout a timeline see how I can best help.

  2. The network parsed into the program is not recognized as an ultrametic network. These are easily identified, and if you use the plotting function from the program, it will help you to compare if this network is what you think it should be, if not. We can dig deeper and see if this is a bug from the program, or the network string needs improvement as well.

Thanks Joe

On Tue, 20 Jul 2021, 23:40 Gustavo A. Ballen, @.***> wrote:

Hello Joe,

Thanks for taking a look into this. I checked out into release v0.6.2-beta (whose version returns v0.6.1-beta after building, though):

hybrid-Lambda/src$ git checkout v0.6.2-beta HEAD is now at de99ef3 update hybrid-Lambda/src$ git status HEAD detached at v0.6.2-beta hybrid-Lambda/src$ make hybrid-Lambda/src$ ./hybrid-Lambda hybrid-Lambda v0.6.1-beta Usage: ...

I ran again the script and all previous instances of segfaults are now gone. However, some of these are instead reported as not ultrametric now, or a couple instances which were initially reported as ultrametric are no longer so:

diff --git a/outerr.txt b/outerr.txt index e076b82..5e5f413 100644--- a/BashScript_Network_NotLevelOne/outerr.txt+++ b/BashScript_Network_NotLevelOne/outerr.txt@@ -42,11 +42,19 @@ level2_2_coal_unit Processing network in line 5: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8781 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 6: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8784 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+WARNING! NOT ULTRAMETRIC!!!+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 7: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -90,7 +98,11 @@ Random seed: 2 Produced gene tree files: level2_2_coal_unit Processing network in line 14:-./script.sh: line 23: 8802 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 15: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -98,11 +110,19 @@ Random seed: 2 Produced gene tree files: level2_2_coal_unit Processing network in line 16:-./script.sh: line 23: 8812 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 17: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!! Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric Produced gene tree files: level2_2_coal_unit Processing network in line 18:@@ -118,7 +138,15 @@ Random seed: 2 Produced gene tree files: level2_2_coal_unit Processing network in line 20:-./script.sh: line 23: 8824 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!!+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 21: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -140,7 +168,11 @@ level2_2_coal_unit Processing network in line 24: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!! Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric Produced gene tree files: level2_2_coal_unit Processing network in line 25:@@ -156,13 +188,31 @@ Random seed: 2 Produced gene tree files: level2_2_coal_unit Processing network in line 27:-./script.sh: line 23: 8838 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 28:-./script.sh: line 23: 8844 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!!+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 29: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8852 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+WARNING! NOT ULTRAMETRIC!!!+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 30: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -176,7 +226,11 @@ Random seed: 2 Produced gene tree files: level2_2_coal_unit Processing network in line 32:-./script.sh: line 23: 8865 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 33: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -187,7 +241,12 @@ Processing network in line 34: Default Kingman coalescent on all branches. Default population size of 10000 on all branches. WARNING! NOT ULTRAMETRIC!!!-./script.sh: line 23: 8874 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 35: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -197,7 +256,9 @@ level2_2_coal_unit Processing network in line 36: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8883 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 37: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -219,7 +280,9 @@ level2_2_coal_unit Processing network in line 40: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8894 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 41: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -235,9 +298,15 @@ level2_2_coal_unit Processing network in line 43: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8904 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 44:-./script.sh: line 23: 8909 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 45: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -247,7 +316,11 @@ level2_2_coal_unit Processing network in line 46: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!! Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric Produced gene tree files: level2_2_coal_unit Processing network in line 47:@@ -269,7 +342,15 @@ Random seed: 2 Produced gene tree files: level2_2_coal_unit Processing network in line 50:-./script.sh: line 23: 8925 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!!+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 51: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -279,8 +360,9 @@ level2_2_coal_unit Processing network in line 52: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-WARNING! NOT ULTRAMETRIC!!!-./script.sh: line 23: 8934 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 53: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -290,7 +372,11 @@ level2_2_coal_unit Processing network in line 54: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!! Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric Produced gene tree files: level2_2_coal_unit Processing network in line 55:@@ -302,19 +388,35 @@ level2_2_coal_unit Processing network in line 56: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.-./script.sh: line 23: 8947 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 57: Default Kingman coalescent on all branches. Default population size of 10000 on all branches. WARNING! NOT ULTRAMETRIC!!!-./script.sh: line 23: 8950 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 58: Default Kingman coalescent on all branches. Default population size of 10000 on all branches. WARNING! NOT ULTRAMETRIC!!!-./script.sh: line 23: 8958 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 59:-./script.sh: line 23: 8966 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 60: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.@@ -325,15 +427,28 @@ Processing network in line 61: Default Kingman coalescent on all branches. Default population size of 10000 on all branches. WARNING! NOT ULTRAMETRIC!!!-./script.sh: line 23: 8978 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+Produced gene tree files:+level2_2_coal_unit Processing network in line 62: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.+WARNING! NOT ULTRAMETRIC!!! Random seed: 2+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric+WARNING: Gene tree is not ultrametric Produced gene tree files: level2_2_coal_unit Processing network in line 63:-./script.sh: line 23: 8988 Segmentation fault (core dumped) hybrid-Lambda -spcu $line -num 3 -seed 2 -o level2_2+Default Kingman coalescent on all branches.+Default population size of 10000 on all branches.+Random seed: 2+Produced gene tree files:+level2_2_coal_unit Processing network in line 64: Default Kingman coalescent on all branches. Default population size of 10000 on all branches.

Hope these results help to solve the bug. We are still trying to understand what triggers the warning, and whether it is safe to use hybrid-Lambda in the presence of such lack of ultrametricity. We may discard these not-ultrametric networks if it is better to do so, but given that these were simulated and in principle should all be ultrametric, we are still puzzled by this warning.

By looking at the code of net.cpp (L323 https://github.com/hybridLambda/hybrid-Lambda/blob/e782a53947ba0aab3967edcc2e53bb9466a3304e/src/net.cpp#L323) we find that the function check_isUltrametric uses 0.000001 in square difference of two paths as a threshold for considering them different, and hence the network as not ultrametric. That makes absolute differences in path length greater than 0.001 to fail the test. Is there any reason for choosing such particular threshold and not perhaps a bigger value? That could help us understand whether the simulation of the networks might produce misleading results.

Thanks again,

Gustavo

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/hybridLambda/hybrid-Lambda/issues/42#issuecomment-883752128, or unsubscribe https://github.com/notifications/unsubscribe-auth/AA4FP7IEGD2NU6GXTYCHBGLTYX3NVANCNFSM5AMHOUIA .