v1.0.0 (#221)

new features:
- countquartetsintrees (was in v0.16 but not v0.17)
- tablequartetCF (was writeTableCF in v0.16, now returns NamedTuple and has option to skip 4-taxon sets with no informative genes)
- tiplabels for vector of phylogenies: union

breaking changes: renaming with small case

---------

Co-authored-by: Joshua Justison (jjustison)

Author: cecileane

Project.toml

@@ -1,7 +1,7 @@
 name = "PhyloNetworks"
 uuid = "33ad39ac-ed31-50eb-9b15-43d0656eaa72"
 license = "MIT"
-version = "0.17.0"
+version = "1.0.0"
 
 [deps]
 BioSequences = "7e6ae17a-c86d-528c-b3b9-7f778a29fe59"
@@ -19,22 +19,24 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 NLopt = "76087f3c-5699-56af-9a33-bf431cd00edd"
 Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
 Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c"
+StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 
 [compat]
-BioSequences = "2.0, 3"
-BioSymbols = "4.0, 5"
+BioSequences = "2, 3"
+BioSymbols = "4, 5"
 CSV = "0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.10"
-Combinatorics = "0.7, 1.0"
+Combinatorics = "0.7, 1"
 DataFrames = "1.3"
 DataStructures = "0.9, 0.10, 0.11, 0.12, 0.13, 0.14, 0.15, 0.16, 0.17, 0.18"
 Distributions = "0.15, 0.16, 0.17, 0.18, 0.19, 0.20, 0.21, 0.22, 0.23, 0.24, 0.25"
 FASTX = "1.1, 2"
 Functors = "0.2, 0.3, 0.4, 0.5"
 NLopt = "0.5.1, 0.6, 1"
+StaticArrays = "0.8.3, 0.9, 0.10, 0.11, 0.12, 1"
 StatsBase = "0.26, 0.27, 0.28, 0.29, 0.30, 0.31, 0.32, 0.33, 0.34"
-julia = "1.5, 1.6, 1.7"
+julia = "1.5"
 
 [extras]
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"

src/PhyloNetworks.jl

@@ -22,6 +22,7 @@ module PhyloNetworks
     using DataFrames # innerjoin new in v0.21
     using DataStructures # for updateInCycle with priority queue
     using Distributions #for RateVariationAcrossSites
+    using StaticArrays # for quartets
     using FASTX
     using Functors: fmap
     using NLopt # for branch lengths optimization
@@ -122,7 +123,10 @@ module PhyloNetworks
         #maxParsimonyNet # broken after v0.17 refactoring: fix network search
         readfastatodna,
         # neighbor joining
-        nj
+        nj,
+        # quartets
+        tablequartetCF,
+        countquartetsintrees
 
     ##Constants
     const fAbsBL = 1e-10
@@ -148,5 +152,6 @@ module PhyloNetworks
     include("graph_components.jl")
     include("deprecated.jl")
     include("nj.jl")
+    include("quartets.jl")
 
 end #module

src/auxiliary.jl

@@ -572,7 +572,7 @@ level-1, with valid internal fields `e.inte1` (to track which cycle
 function deleteEdge!(net::HybridNetwork, e::Edge; part::Bool=true)
     if part
         if e.inte1 == -1 && !e.hybrid && !isempty(net.partition) && !isTree(net)
-            ind = whichPartition(net,e)
+            ind = whichpartition(net,e)
             indE = getIndex(e,net.partition[ind].edges)
             deleteat!(net.partition[ind].edges,indE)
         end
@@ -1460,50 +1460,50 @@ end
 # better to return the index than the partition itself, because we need the index
 # to use splice and delete it from net.partition later on
 # cycle: is the number to look for partition on that cycle only
-function whichPartition(net::HybridNetwork,edge::Edge,cycle::Integer)
+function whichpartition(net::HybridNetwork, edge::Edge, cycle::Integer)
     !edge.hybrid || error("edge $(edge.number) is hybrid so it cannot be in any partition")
     edge.inte1 == -1 || error("edge $(edge.number) is in cycle $(edge.inte1) so it cannot be in any partition")
     @debug "search partition for edge $(edge.number) in cycle $(cycle)"
     in(edge,net.edge) || error("edge $(edge.number) is not in net.edge")
     for i in 1:length(net.partition)
         @debug "looking for edge $(edge.number) in partition $(i): $([e.number for e in net.partition[i].edges])"
-        if(in(cycle,net.partition[i].cycle))
+        if in(cycle, net.partition[i].cycle)
             @debug "looking for edge $(edge.number) in partition $(i), with cycle $(cycle): $([e.number for e in net.partition[i].edges])"
             if in(edge,net.partition[i].edges)
                 @debug "partition for edge $(edge.number) is $([e.number for e in net.partition[i].edges])"
                 return i
             end
         end
     end
-    @debug begin; printPartitions(net); "" end
+    @debug begin; printpartitions(net); "" end
     error("edge $(edge.number) is not hybrid, nor part of any cycle, and it is not in any partition")
 end
 
 # function to get the partition where an edge is
 # returns the index of the partition, or error if not found
 # better to return the index than the partition itself, because we need the index
 # to use splice and delete it from net.partition later on
-function whichPartition(net::HybridNetwork,edge::Edge)
+function whichpartition(net::HybridNetwork, edge::Edge)
     !edge.hybrid || error("edge $(edge.number) is hybrid so it cannot be in any partition")
     edge.inte1 == -1 || error("edge $(edge.number) is in cycle $(edge.inte1) so it cannot be in any partition")
     @debug "search partition for edge $(edge.number) without knowing its cycle"
     in(edge,net.edge) || error("edge $(edge.number) is not in net.edge")
     for i in 1:length(net.partition)
         @debug "looking for edge $(edge.number) in partition $(i): $([e.number for e in net.partition[i].edges])"
-        if(in(edge,net.partition[i].edges))
+        if in(edge,net.partition[i].edges)
             @debug "partition for edge $(edge.number) is $([e.number for e in net.partition[i].edges])"
             return i
         end
     end
-    @debug begin printPartitions(net); "printed partitions" end
+    @debug begin printpartitions(net); "printed partitions" end
     error("edge $(edge.number) is not hybrid, nor part of any cycle, and it is not in any partition")
 end
 
-# function that will print the partition of net
-function printPartitions(net::HybridNetwork)
-    println("partition.cycle\t partition.edges")
+printpartitions(x) = printpartitions(stdout::IO, x)
+function printpartitions(io::IO, net::HybridNetwork)
+    println(io, "partition.cycle\t partition.edges")
     for p in net.partition
-        println("$(p.cycle)\t\t $([e.number for e in p.edges])")
+        println(io, "$(p.cycle)\t\t $([e.number for e in p.edges])")
     end
 end
 

src/deprecated.jl

@@ -6,6 +6,7 @@
 @deprecate readBootstrapTrees readmultinewick_files
 @deprecate writeTopology      writenewick
 @deprecate writeMultiTopology writemultinewick
+Base.@deprecate_moved snaq! "SNaQ"
 # minor name changes: camelCase to smallcase
 @deprecate printEdges printedges
 @deprecate printNodes printnodes

src/descriptive.jl

@@ -14,9 +14,30 @@ function tiplabels(net::HybridNetwork)
     return String[l.name for l in net.leaf] # AbstractString does not work for use by tree2Matrix
 end
 
-# function that we need to overwrite to avoid printing useless scary
-# output for HybridNetworks
-# PROBLEM: writenewick changes the network and thus show changes the network
+# extract & sort the union of taxa of list of gene trees
+function tiplabels(trees::Vector{HybridNetwork})
+    taxa = reduce(union, tiplabels(t) for t in trees)
+    return sort_stringasinteger!(taxa)
+end
+
+"""
+    sort_stringasinteger!(taxa)
+
+Sort a vector of strings `taxa`, numerically if
+elements can be parsed as an integer, alphabetically otherwise.
+"""
+function sort_stringasinteger!(taxa)
+    sortby = x->parse(Int,x)
+    try
+        parse.(Int,taxa)
+    catch
+        sortby = identity
+    end
+    sort!(taxa, by=sortby)
+    return taxa
+end
+
+# 'show' should not (and does not) modify the object
 function Base.show(io::IO, obj::HybridNetwork)
     disp = "$(typeof(obj)), "
     if obj.isrooted
@@ -39,7 +60,10 @@ function Base.show(io::IO, obj::HybridNetwork)
     end
     par = ""
     try
-        # par = writenewick(obj,round=true) # but writenewick changes the network, not good
+        #= *not* 'writenewick(obj,round=true)' because writenewick changes
+        the network (runs directedges! and if it fails, tries other rootings).
+        writesubtree! does *not* modify the network.
+        =#
         s = IOBuffer()
         writesubtree!(s, obj, false,true, true,3,true)
         par = String(take!(s))
@@ -53,8 +77,6 @@ function Base.show(io::IO, obj::HybridNetwork)
     println(io, disp)
 end
 
-
-
 function Base.show(io::IO, obj::Node)
     disp = "$(typeof(obj)):"
     disp = disp * "\n number:$(obj.number)"