bug fix, for v0.14.2 (#165)

bug fixes for inline strings from CSV 0.9 better plots in docs
JuliaPhylo · Oct 23, 2021 · a8155c8 · a8155c8 · cecileane · Oct 23, 2021
1 parent e0838bc
commit a8155c8
Show file tree

Hide file tree

Showing 7 changed files with 57 additions and 50 deletions.
diff --git a/Project.toml b/Project.toml
@@ -1,7 +1,7 @@
 name = "PhyloNetworks"
 uuid = "33ad39ac-ed31-50eb-9b15-43d0656eaa72"
 license = "MIT"
-version = "0.14.1"
+version = "0.14.2"
 
 [deps]
 BioSequences = "7e6ae17a-c86d-528c-b3b9-7f778a29fe59"

diff --git a/docs/Project.toml b/docs/Project.toml
@@ -13,8 +13,4 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 StatsModels = "3eaba693-59b7-5ba5-a881-562e759f1c8d"
 
 [compat]
-BioSymbols = "4.0"
-CSV = "0.4, 0.5, 0.6, 0.7, 0.8"
-DataFrames = "1.0"
-Documenter = "~0.26"
-RCall = "0.13"
+Documenter = "~0.27"
diff --git a/docs/src/man/dist_reroot.md b/docs/src/man/dist_reroot.md
@@ -81,7 +81,7 @@ network below. We plotted the edge numbers, because we will want to use them
 later to place the root.
 
 ```@example dist_reroot
-net7taxa = readTopology("(C,D,((O,(E,#H7:::0.196):0.314):0.664,(B,((A1,A2))#H7:::0.804):10.0):10.0);")
+net7taxa = readTopology("(C,D,((O,(E,#H7:::0.196):0.314):0.664,(((A1,A2))#H7:::0.804,B):10.0):10.0);")
 R"svg(name('reroot_net7taxa_1.svg'), width=4, height=4)" # hide
 R"par"(mar=[0,0,0,0]) # hide
 plot(net7taxa, :R, showGamma=true, showEdgeNumber=true, tipOffset=0.2);
@@ -98,7 +98,7 @@ or with the A clade (on edge 11), will fail with a `RootMismatch` error:
 ```julia
 rootatnode!(net7taxa, "A1"); # ERROR: RootMismatch: non-leaf node 5 had 0 children. ...
 rootatnode!(net7taxa, "A2"); # ERROR: RootMismatch (again)
-rootonedge!(net7taxa, 11);   # ERROR: RootMismatch (again)
+rootonedge!(net7taxa, 10);   # ERROR: RootMismatch (again)
 ```
 
 In this case, however, it is possible to root the network on either parent edge
@@ -109,9 +109,10 @@ We get these 2 rooted versions of the network:
 R"svg(name('reroot_net7taxa_2.svg'), width=7, height=4)"; # hide
 R"layout(matrix(1:2,1,2))";
 R"par"(mar=[0,0,0.5,0]); # hide
-rootonedge!(net7taxa, 12);
-plot(net7taxa, :R, showGamma=true, tipOffset=0.2);
-R"mtext"("rooted on hybrid edge 12 (major)", line=-1)
+rootonedge!(net7taxa, 11);
+rotate!(net7taxa, -5)
+plot(net7taxa, :R, showGamma=true, tipOffset=0.2, showNodeNumber=true);
+R"mtext"("rooted on hybrid edge 11 (major)", line=-1)
 rootonedge!(net7taxa, 5);
 plot(net7taxa, :R, showGamma=true, tipOffset=0.2);
 R"mtext"("rooted on hybrid edge 5 (minor)", line=-1);
@@ -201,13 +202,13 @@ truenet = readTopology("((((D:0.4,C:0.4):4.8,((A:0.8,B:0.8):2.2)#H1:2.2::0.7):4.
 hardwiredClusterDistance(net1, truenet, true)
 ```
 ```@example dist_reroot
-R"svg(name('truenet.svg'), width=4, height=4)" # hide
+R"svg(name('truenet_sim.svg'), width=4, height=4)" # hide
 R"par"(mar=[0,0,0,0]) # hide
 plot(truenet, :R, useEdgeLength=true, showGamma=true);
 R"dev.off()" # hide
 nothing # hide
 ```
-![truenet](../assets/figures/truenet.svg)
+![truenet](../assets/figures/truenet_sim.svg)
 
 Our estimated network is not the same as the true network:
 - the underlying tree is correctly estimated

diff --git a/docs/src/man/fitDiscrete.md b/docs/src/man/fitDiscrete.md
@@ -22,7 +22,7 @@ The simplest way is to use a vector of species names with a data frame of traits
 
 ```@repl fitdiscrete_trait
 # read in network
-net = readTopology("(O:4,(A:3,((B:0.4)#H1:1.6::0.92,((C:0.4,#H1:0::0.08):0.6,(D:.2,E:.2):0.8):1):1):1);");
+net = readTopology("(O:4,(A:3,((B:0.4)#H1:1.6::0.92,((#H1:0::0.08,C:0.4):0.6,(D:.2,E:.2):0.8):1):1):1);");
 # read in trait data
 species = ["C","A","D","B","O","E"];
 dat = DataFrame(trait=["hi","lo","lo","hi","lo","lo"])
@@ -49,13 +49,11 @@ isequal(species[o], tipLabels(net)) # true :)
 traitcolor = map(x -> (x=="lo" ? "grey" : "red"), dat.trait[o])
 leaves = res[13][!,:lea]
 R"points"(x=res[13][leaves,:x] .+0.1, y=res[13][leaves,:y], pch=16, col=traitcolor, cex=1.5); # adds grey & red points
-R"legend"(x=1, y=2, legend=["hi","lo"], pch=16, col=["red","grey"],
+R"legend"(x=1, y=7, legend=["hi","lo"], pch=16, col=["red","grey"],
           title="my trait", bty="n",var"title.adj"=0);
-# next: add arrow to show gene flow edge, and proportion γ of genes affected
-hi = findfirst([!e.isMajor for e in net.edge]) # 6 : "h"ybrid "i"ndex: index of gene flow edge (minor hybrid) in net
-(hx1, hx2, hy1, hy2) = (res[i][hi] for i in 9:12); # coordinates for minor hybrid edge; 1=start, 2=end
-R"arrows"(hx1, hy1, hx2, hy2, col="deepskyblue", length=0.08, angle=20); # adds the arrow
-R"text"(res[14][hi,:x]-0.2, res[14][hi,:y]+0.1, res[14][hi,:gam], col="deepskyblue", cex=0.75); # add the γ value
+# next: add to gene flow edge the proportion γ of genes affected
+hi = findfirst([!e.isMajor for e in net.edge]) # 6 : "h"ybrid "i"ndex: index of gene flow edge (minor hybrid) in net: horizontal segment
+R"text"(res[14][hi,:x]-0.3, res[14][hi,:y]-0.1, res[14][hi,:gam], col="deepskyblue", cex=0.75); # add the γ value
 R"dev.off"(); # hide
 nothing # hide
 ```

diff --git a/docs/src/man/fixednetworkoptim.md b/docs/src/man/fixednetworkoptim.md
@@ -147,7 +147,7 @@ Now imagine that our outgroup is taxon A.
   (rooted network on the left below).
 - For the second best network in our list, there are 2 ways to root it
   with A: on the external edge 8 to A (top right), or on its parent edge 10
-  (bottom right). These 2 options give quite different rooted versions
+  These 2 options give quite different rooted versions
   of the network, one of which requires the existence of an unsampled taxon,
   sister to BOECD, that would have contributed to introgression into
   an ancestor of E. The second rooted version says that an ancestor of
@@ -162,14 +162,18 @@ R"par"(mar=[0,0,0.5,0]) # hide
 rootonedge!(netlist[1], 10); # root best net to make A outgroup
 rotate!(netlist[1], -4); # to 'un-cross' edges
 rotate!(netlist[1], -6);
+rotate!(netlist[1], -5);
 plot(netlist[1], :R, showGamma=true, tipOffset=0.1);
 R"mtext"("best net, score=28.3", line=-1);
 global_logger(NullLogger()); # hide
 rootatnode!(netlist[2], "A"); # net with modified direction: first way to make A outgroup
 global_logger(baselogger);   # hide
+rotate!(netlist[2], -4) # to 'un-cross' edges
+rotate!(netlist[2], -6)
 plot(netlist[2], :R, showGamma=true, tipOffset=0.1);
 R"mtext"("second best in list, score=31.5\nrequires unsampled population", line=-2);
 rootonedge!(netlist[2], 10) # net with modified direction: second way to make A outgroup
+for i in [9,-7] rotate!(netlist[2], i); end; # to 'un-cross' edges
 plot(netlist[2], :R, showGamma=true, tipOffset=0.1);
 R"mtext"("second best in list, score=31.5\ndifferent root position", line=-2);
 R"dev.off()"; # hide

diff --git a/docs/src/man/snaq_plot.md b/docs/src/man/snaq_plot.md
@@ -7,6 +7,7 @@ astralfile = joinpath(dirname(pathof(PhyloNetworks)), "..","examples","astral.tr
 astraltree = readMultiTopology(astralfile)[102] # 102th tree = last tree here
 net0 = readTopology(joinpath(dirname(pathof(PhyloNetworks)), "..","examples","net0.out"))
 net1 = readTopology(joinpath(dirname(pathof(PhyloNetworks)), "..","examples","net1.out"))
+rotate!(net1, -6)
 net2 = readTopology(joinpath(dirname(pathof(PhyloNetworks)), "..","examples","net2.out"))
 net3 = readTopology(joinpath(dirname(pathof(PhyloNetworks)), "..","examples","net3.out"))
 net0.loglik = 53.53150526187732
@@ -66,8 +67,8 @@ We can visualize the estimated network and its inheritance values γ, which
 measure the proportion of genes inherited via each parent at a reticulation event
 (e.g. proportion of genes inherited via gene flow).
 ```@example snaqplot
-R"svg(name('snaqplot_net1_1.svg'), width=4, height=3)" # hide
-R"par"(mar=[0,0,0,0]) # hide
+R"svg(name('snaqplot_net1_1.svg'), width=4, height=3)"; # hide
+R"par"(mar=[0,0,0,0]); # hide
 plot(net1, :R, showGamma=true);
 R"dev.off()"; # hide
 nothing # hide
@@ -306,8 +307,8 @@ then we tell R to wrap up and save its image file.
 ```@example snaqplot
 using PhyloPlots # to visualize networks
 using RCall      # to send additional commands to R like this: R"..."
-R"name = function(x) file.path('..', 'assets', 'figures', x)" # function to create file name in appropriate folder
-R"svg(name('snaqplot_net1_2.svg'), width=4, height=3)" # starts image file
+imagefilename = "../assets/figures/snaqplot_net1_2.svg"
+R"svg"(imagefilename, width=4, height=3) # starts image file
 R"par"(mar=[0,0,0,0]) # to reduce margins (no margins at all here)
 plot(net1, :R, showGamma=true, showEdgeNumber=true); # network is plotted & sent to file
 R"dev.off()"; # wrap up and save image file
@@ -322,11 +323,17 @@ major edge with γ>0.5), and edges were annotated with their internal numbers.
 
 Type `?` to switch to the help mode
 of Julia, then type the name of the function, here `plot`.
-Edge colors can be modified, for instance.
+Below are two visualizations.
+The first uses the default style (`:fulltree`) and modified edge colors.
+The second uses the `:majortree` style.
+That style doesn't have an arrow by default for minor hybrid edges,
+but we can ask for one by specifying a positive arrow length.
 ```@example snaqplot
-R"svg(name('snaqplot_net1_3.svg'), width=4, height=3)" # hide
+R"svg(name('snaqplot_net1_3.svg'), width=7, height=3)" # hide
 R"par"(mar=[0,0,0,0]) # hide
-plot(net1, :R, showEdgeLength=true, minorHybridEdgeColor="tan")
+R"layout(matrix(1:2,1,2))";
+plot(net1, :R, showEdgeLength=true, minorHybridEdgeColor="tan");
+plot(net1, :R, style=:majortree, arrowlen=0.07);
 R"dev.off()"; # hide
 nothing # hide
 ```

diff --git a/src/readData.jl b/src/readData.jl
@@ -47,7 +47,7 @@ end
 writeTableCF(d::DataCF) = writeTableCF(d.quartet)
 
 function writeTableCF(quartets::Vector{QuartetT{T}},
-                      taxa=Vector{String}()::AbstractVector{String}) where T<:AbstractVector
+                      taxa::AbstractVector{<:AbstractString}=Vector{String}()) where T<:AbstractVector
     V = eltype(T) # would expect Float64, but Int would be reasonable if counts are not normalized
     V <: Real || error("CFs need to take real values")
     df = DataFrames.DataFrame(t1=String[],t2=String[],t3=String[],t4=String[],
@@ -247,7 +247,6 @@ end
 
 # function to list all quartets for a set of taxa names
 # return a vector of quartet objects, and if writeFile=true, writes a file
-# warning: taxon has to be vector of String, vector of AbstractString do not work
 function allQuartets(taxon::Union{Vector{<:AbstractString},Vector{Int}}, writeFile::Bool)
     quartets = combinations(taxon,4)
     vquartet = Quartet[];
@@ -305,7 +304,7 @@ end
 
 # function that will not use randQuartets(list of quartets,...)
 # this function uses whichQuartet to avoid making the list of all quartets
-# fixit: i think we should write always the file, but not sure
+# fixit: writeFile is not used. remove the option?
 function randQuartets(taxon::Union{Vector{<:AbstractString},Vector{Int}},num::Integer, writeFile::Bool)
     randquartets = Quartet[]
     n = length(taxon)
@@ -516,10 +515,10 @@ function calculateObsCFAll_noDataCF!(quartets::Vector{Quartet}, trees::Vector{Hy
 end
 
 """
-    countquartetsintrees(trees [, taxonmap=Dict{String,String}]; which=:all, weight_byallele=true)
+    countquartetsintrees(trees [, taxonmap]; which=:all, weight_byallele=true)
 
 Calculate the quartet concordance factors (CF) observed in the `trees` vector.
-If present, `taxonmap` should map each allele name to it's species name.
+If present, `taxonmap` should be a dictionary that maps each allele name to it's species name.
 To save to a file, first convert to a data frame using [`writeTableCF`](@ref).
 When `which=:all`, quartet CFs are calculated for all 4-taxon sets.
 (Other options are not implemented yet.)
@@ -628,9 +627,9 @@ julia> show(writeTableCF(q,t), allcols=true)
 ```
 """
 function countquartetsintrees(tree::Vector{HybridNetwork},
-                           taxonmap=Dict{String,String}()::Dict{String,String};
-                           whichQ=:all::Symbol, weight_byallele=false::Bool,
-                           showprogressbar=true::Bool)
+                           taxonmap::Dict=Dict{String,String}();
+                           whichQ::Symbol=:all, weight_byallele::Bool=false,
+                           showprogressbar::Bool=true)
     whichQ in [:all, :intrees] || error("whichQ must be either :all or :intrees, but got $whichQ")
     if isempty(taxonmap)
         taxa = unionTaxa(tree)
@@ -691,7 +690,7 @@ function countquartetsintrees(tree::Vector{HybridNetwork},
 end
 function countquartetsintrees!(quartet::Vector{QuartetT{MVector{4,Float64}}},
             tree::HybridNetwork, whichQ::Symbol, weight_byallele::Bool, nCk::Matrix,
-            taxonnumber::Dict{String,Int}, taxonmap::Dict{String,String})
+            taxonnumber::Dict{<:AbstractString,Int}, taxonmap::Dict{<:AbstractString,<:AbstractString})
     tree.numHybrids == 0 || error("input phylogenies must be trees")
     # next: reset node & edge numbers so that they can be used as indices: 1,2,3,...
     resetNodeNumbers!(tree; checkPreorder=true, type=:postorder) # leaves first & post-order
@@ -875,7 +874,10 @@ function readInputData(trees::Vector{HybridNetwork}, quartetfile::AbstractString
 end
 
 
-function readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{<:AbstractString}, Vector{Int}}, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
+function readInputData(treefile::AbstractString, whichQ::Symbol=:all, numQ::Integer=0,
+        taxa::Union{Vector{<:AbstractString}, Vector{Int}}=unionTaxaTree(treefile),
+        writetab::Bool=true, filename::AbstractString="none",
+        writeFile::Bool=false, writeSummary::Bool=true)
     if writetab
         if(filename == "none")
             filename = "tableCF.txt" # "tableCF$(string(integer(time()/1000))).txt"
@@ -891,18 +893,13 @@ function readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer,
     readInputData(trees, whichQ, numQ, taxa, writetab, filename, writeFile, writeSummary)
 end
 
-readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}, writetab::Bool) = readInputData(treefile, whichQ, numQ, taxa, writetab, "none", false, true)
-readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}) = readInputData(treefile, whichQ, numQ, taxa, true, "none",false, true)
-readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool, filename::AbstractString) = readInputData(treefile, whichQ, numQ, unionTaxaTree(treefile), writetab, filename,false, true)
 readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer, writetab::Bool) = readInputData(treefile, whichQ, numQ, unionTaxaTree(treefile), writetab, "none",false, true)
-readInputData(treefile::AbstractString, whichQ::Symbol, numQ::Integer) = readInputData(treefile, whichQ, numQ, unionTaxaTree(treefile), true, "none",false, true)
-readInputData(treefile::AbstractString) = readInputData(treefile, :all, 0, unionTaxaTree(treefile), true, "none",false, true)
-readInputData(treefile::AbstractString,taxa::Union{Vector{String}, Vector{Int}}) = readInputData(treefile, :all, 0, taxa, true, "none",false, true)
+readInputData(treefile::AbstractString, taxa::Union{Vector{<:AbstractString}, Vector{Int}}) = readInputData(treefile, :all, 0, taxa, true, "none",false, true)
 # above: the use of unionTaxaTree to set the taxon set
 #        is not good: need to read the tree file twice: get the taxa, then get the trees
 #        this inefficiency was fixed in readTrees2CF
 
-function readInputData(trees::Vector{HybridNetwork}, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{String}, Vector{Int}}, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
+function readInputData(trees::Vector{HybridNetwork}, whichQ::Symbol, numQ::Integer, taxa::Union{Vector{<:AbstractString}, Vector{Int}}, writetab::Bool, filename::AbstractString, writeFile::Bool, writeSummary::Bool)
     if(whichQ == :all)
         numQ == 0 || @warn "set numQ=$(numQ) but whichQ=all, so all quartets will be used and numQ will be ignored. If you want a specific number of 4-taxon subsets not random, you can input with the quartetfile option"
         quartets = allQuartets(taxa,writeFile)
@@ -958,7 +955,7 @@ See also:
 """
 function readTrees2CF(treefile::AbstractString; quartetfile="none"::AbstractString, whichQ="all"::AbstractString, numQ=0::Integer,
                       writeTab=true::Bool, CFfile="none"::AbstractString,
-                      taxa=Vector{String}()::Union{Vector{String},Vector{Int}},
+                      taxa::AbstractVector=Vector{String}(),
                       writeQ=false::Bool, writeSummary=true::Bool, nexus=false::Bool)
     trees = (nexus ? readNexusTrees(treefile, readTopologyUpdate, false, false) : readInputTrees(treefile))
     if length(taxa)==0        # unionTaxa(trees) NOT default argument:
@@ -969,7 +966,11 @@ function readTrees2CF(treefile::AbstractString; quartetfile="none"::AbstractStri
 end
 
 # same as before, but with input vector of HybridNetworks
-function readTrees2CF(trees::Vector{HybridNetwork}; quartetfile="none"::AbstractString, whichQ="all"::AbstractString, numQ=0::Integer, writeTab=true::Bool, CFfile="none"::AbstractString, taxa=unionTaxa(trees)::Union{Vector{String},Vector{Int}}, writeQ=false::Bool, writeSummary=true::Bool)
+function readTrees2CF(trees::Vector{HybridNetwork};
+        quartetfile="none"::AbstractString, whichQ="all"::AbstractString, numQ=0::Integer,
+        writeTab=true::Bool, CFfile="none"::AbstractString,
+        taxa::AbstractVector=unionTaxa(trees),
+        writeQ=false::Bool, writeSummary=true::Bool)
     whichQ == "all" || whichQ == "rand" ||
         error("whichQ should be all or rand, not $(whichQ)")
     if(quartetfile == "none")
@@ -1185,7 +1186,7 @@ end
 # function to determine the resolution of taxa picked from part1,2,3,4 and DataCF
 # names: taxa from part1,2,3,4
 # rownames: taxa from table of obsCF
-function resolution(names::Vector{String},rownames::Vector{String})
+function resolution(names::Vector{<:AbstractString},rownames::Vector{<:AbstractString})
     length(names) == length(rownames) || error("names and rownames should have the same length")
     length(names) == 4 || error("names should have 4 entries, not $(length(names))")
     bin = [n == names[1] || n == names[2] ? 1 : 0 for n in rownames]
@@ -1206,7 +1207,7 @@ end
 # this function is meant to replace extractQuartet! in calculateObsCFAll
 # input: Quartet, Matrix, vector of taxa names
 # returns 1 if quartet found is 12|34, 2 if 13|24, 3 if 14|23, and 0 if not found
-function extractQuartetTree(q::Quartet, M::Matrix{Int},S::Union{Vector{String},Vector{Int}})
+function extractQuartetTree(q::Quartet, M::Matrix{Int},S::Union{Vector{<:AbstractString},Vector{Int}})
     @debug "extractQuartet: $(q.taxon)"
     @debug "matrix: $(M)"
     inds = indexin(q.taxon, S)
@@ -1260,7 +1261,7 @@ end
 # input: list of taxa names, vector of integers (each integer corresponds to a taxon name),
 # num is the number of the Quartet
 # assumes taxa is sorted already
-function createQuartet(taxa::Union{Vector{String},Vector{Int}},qvec::Vector{Int}, num::Int)
+function createQuartet(taxa::Union{Vector{<:AbstractString},Vector{Int}}, qvec::Vector{Int}, num::Int)
     length(qvec) == 4 || error("a quartet should have only 4 taxa, not $(length(qvec))")
     names = String[]
     for i in qvec