diff --git a/ChangeLog.md b/ChangeLog.md
--- a/ChangeLog.md
+++ b/ChangeLog.md
@@ -5,6 +5,16 @@
## Unreleased changes
+## Version 0.6.0.0
+
+- **elynx-tree:** remove parallel folds with layers (`parBranchFoldMapWithLayer`
+ too special and slow).
+- **elynx-tree:** fix various tree instances; add zip trees with appropriate
+ instances.
+- Remove `monad-logger` dependency and implement lighter alternative.
+- Significant changes to the tool chain.
+
+
## Version 0.5.1.0
- **elynx-tree:** new functions `isValidPath`, `isLeaf`, `depth`; add conversion
diff --git a/README.md b/README.md
--- a/README.md
+++ b/README.md
@@ -2,7 +2,7 @@
# The ELynx Suite
-Version: 0.5.1.0.
+Version: 0.6.0.0.
Reproducible evolution made easy.
@@ -39,7 +39,7 @@
# Installation
-ELynx is written in [Haskell](https://www.haskell.org/) and can be installed with [Stack](https://docs.haskellstack.org/en/stable/README/).
+ELynx is written in [Haskell](https://www.haskell.org/) and can be installed with [cabal-install](https://cabal.readthedocs.io/en/3.4/cabal-commands.html) or [Stack](https://docs.haskellstack.org/en/stable/README/).
1. Install Stack with your package manager, or directly from the web
page.
@@ -67,15 +67,15 @@
[PATH](https://en.wikipedia.org/wiki/PATH_(variable)) environment variable. Then, they can be used directly.
-# SLynx
+# Get help
-Handle evolutionary sequences.
+For example:
- stack exec slynx -- --help | head -n -16
+ slynx --help
- ELynx Suite version 0.5.1.0.
+ ELynx Suite version 0.6.0.0.
Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
+ Compiled on September 3, 2021, at 20:56 pm, UTC.
Usage: slynx [-v|--verbosity VALUE] [-o|--output-file-basename NAME]
[-f|--force] [--no-elynx-file] COMMAND
@@ -84,7 +84,7 @@
Available options:
-h,--help Show this help text
-V,--version Show version
- -v,--verbosity VALUE Be verbose; one of: Quiet Warning Info Debug
+ -v,--verbosity VALUE Be verbose; one of: Quiet Warn Info Debug
(default: Info)
-o,--output-file-basename NAME
Specify base name of output file
@@ -112,107 +112,33 @@
- Protein (amino acids)
- ProteinX (amino acids; including gaps)
- ProteinS (amino acids; including gaps, and translation stops)
-
-
-## Concatenate
-
-Concatenate multi sequence alignments.
-
- stack exec slynx -- concatenate --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: slynx concatenate (-a|--alphabet NAME) INPUT-FILE
- Concatenate sequences found in input files.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -a,--alphabet NAME Specify alphabet type NAME
- INPUT-FILE Read sequences from INPUT-FILE
- -h,--help Show this help text
-
-
-## Examine
-
-Examine sequence with `slynx examine`.
-
- stack exec slynx -- examine --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: slynx examine (-a|--alphabet NAME) INPUT-FILE [--per-site]
- Examine sequences. If data is a multi sequence alignment, additionally analyze columns.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -a,--alphabet NAME Specify alphabet type NAME
- INPUT-FILE Read sequences from INPUT-FILE
- --per-site Report per site summary statistics
- -h,--help Show this help text
-
-
-## Filter
-
-Filter sequences with `filer-rows`.
-
- stack exec slynx -- filter-rows --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: slynx filter-rows (-a|--alphabet NAME) INPUT-FILE [--longer-than LENGTH]
- [--shorter-than LENGTH] [--standard-characters]
- Filter rows (or sequences) found in input files.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -a,--alphabet NAME Specify alphabet type NAME
- INPUT-FILE Read sequences from INPUT-FILE
- --longer-than LENGTH Only keep sequences longer than LENGTH
- --shorter-than LENGTH Only keep sequences shorter than LENGTH
- --standard-characters Only keep sequences containing at least one standard
- (i.e., non-IUPAC) character
- -h,--help Show this help text
-
-Filter columns of multi sequence alignments with `filter-columns`.
-
- stack exec slynx -- filter-columns --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
+ - ProteinI (amino acids; including gaps, translation stops, and IUPAC codes)
+ ELynx
+ -----
+ A Haskell library and tool set for computational biology. The goal of ELynx is
+ reproducible research. Evolutionary sequences and phylogenetic trees can be
+ read, viewed, modified and simulated. The command line with all arguments is
+ logged consistently, and automatically. Data integrity is verified using SHA256
+ sums so that validation of past analyses is possible without the need to
+ recompute the result.
- Usage: slynx filter-columns (-a|--alphabet NAME) INPUT-FILE
- [--standard-chars DOUBLE]
- Filter columns of multi sequence alignments.
+ slynx Analyze, modify, and simulate evolutionary sequences.
+ tlynx Analyze, modify, and simulate phylogenetic trees.
+ elynx Validate and redo past analyses.
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -a,--alphabet NAME Specify alphabet type NAME
- INPUT-FILE Read sequences from INPUT-FILE
- --standard-chars DOUBLE Keep columns with a proportion standard (non-IUPAC)
- characters larger than DOUBLE in [0,1]
- -h,--help Show this help text
+ Get help for sub commands:
+ slynx examine --help
-## Simulate
+## Sub command
-Simulate sequences with `slynx simulate`.
+The documentation of sub commands can be accessed separately:
- stack exec slynx -- simulate --help
+ slynx simulate --help
- ELynx Suite version 0.5.1.0.
+ ELynx Suite version 0.6.0.0.
Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
+ Compiled on September 3, 2021, at 20:56 pm, UTC.
Usage: slynx simulate (-t|--tree-file Name) [-s|--substitution-model MODEL]
[-m|--mixture-model MODEL] [-e|--edm-file NAME]
@@ -281,315 +207,4 @@
LG exchangeabilities with site profiles (Phylobayes) given in FILES.
-m "EDM(LG-Custom)" -p FILES
For special mixture models, mixture weights are optional.
-
-
-## Sub-sample
-
-Sub-sample columns from multi sequence alignments.
-
- stack exec slynx -- sub-sample --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: slynx sub-sample (-a|--alphabet NAME) INPUT-FILE
- (-n|--number-of-sites INT)
- (-m|--number-of-alignments INT) [-S|--seed [INT]]
- Sub-sample columns from multi sequence alignments.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -a,--alphabet NAME Specify alphabet type NAME
- INPUT-FILE Read sequences from INPUT-FILE
- -n,--number-of-sites INT Number of sites randomly drawn with replacement
- -m,--number-of-alignments INT
- Number of multi sequence alignments to be created
- -S,--seed [INT] Seed for random number generator; list of 32 bit
- integers with up to 256 elements (default: random)
- -h,--help Show this help text
-
- Create a given number of multi sequence alignments, each of which contains a given number of random sites drawn from the original multi sequence alignment.
-
-
-## Translate
-
-Translate sequences.
-
- stack exec slynx -- translate --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: slynx translate (-a|--alphabet NAME) INPUT-FILE (-r|--reading-frame INT)
- (-u|--universal-code CODE)
- Translate from DNA to Protein or DNAX to ProteinX.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -a,--alphabet NAME Specify alphabet type NAME
- INPUT-FILE Read sequences from INPUT-FILE
- -r,--reading-frame INT Reading frame [0|1|2].
- -u,--universal-code CODE universal code; one of: Standard,
- VertebrateMitochondrial.
- -h,--help Show this help text
-
-
-# TLynx
-
-Handle phylogenetic trees in Newick format.
-
- stack exec tlynx -- --help | head -n -16
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx [-v|--verbosity VALUE] [-o|--output-file-basename NAME]
- [-f|--force] [--no-elynx-file] COMMAND
- Compare, examine, and simulate phylogenetic trees.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -v,--verbosity VALUE Be verbose; one of: Quiet Warning Info Debug
- (default: Info)
- -o,--output-file-basename NAME
- Specify base name of output file
- -f,--force Ignore previous analysis and overwrite existing
- output files.
- --no-elynx-file Do not write data required to reproduce an analysis.
-
- Available commands:
- compare Compare two phylogenetic trees (compute distances and branch-wise differences).
- connect Connect two phylogenetic trees in all ways (possibly honoring constraints).
- distance Compute distances between many phylogenetic trees.
- examine Compute summary statistics of phylogenetic trees.
- shuffle Shuffle a phylogenetic tree (keep coalescent times, but shuffle topology and leaves).
- simulate Simulate phylogenetic trees using a birth and death or coalescent process.
-
-
- Available tree file formats:
- - Newick Standard: Branch support values are stored in square brackets after branch lengths.
- - Newick IqTree: Branch support values are stored as node names after the closing bracket of forests.
- - Newick RevBayes: Key-value pairs is provided in square brackets after node names as well as branch lengths. XXX: Key value pairs are ignored at the moment.
-
-
-## Compare
-
-Compute distances between phylogenetic trees.
-
- stack exec tlynx -- compare --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx compare [-n|--normalize] [-b|--bipartitions] [-t|--intersect]
- [-f|--newick-format FORMAT] NAMES
- Compare two phylogenetic trees (compute distances and branch-wise differences).
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -n,--normalize Normalize trees before comparison
- -b,--bipartitions Print and plot common and missing bipartitions
- -t,--intersect Compare intersections; i.e., before comparison, drop
- leaves that are not present in the other tree
- -f,--newick-format FORMAT
- Newick tree format: Standard, IqTree, or RevBayes;
- default: Standard; for detailed help, see 'tlynx
- --help'
- NAMES Tree files
- -h,--help Show this help text
-
-
-## Connect
-
-Connect two phylogenetic tree in all ways (possibly honoring constraints).
-
- stack exec tlynx -- connect --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx connect [-f|--newick-format FORMAT] [-c|--contraints CONSTRAINTS]
- TREE-FILE-A TREE-FILE-B
- Connect two phylogenetic trees in all ways (possibly honoring constraints).
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -f,--newick-format FORMAT
- Newick tree format: Standard, IqTree, or RevBayes;
- default: Standard; for detailed help, see 'tlynx
- --help'
- -c,--contraints CONSTRAINTS
- File containing one or more Newick trees to be used
- as constraints
- TREE-FILE-A File containing the first Newick tree
- TREE-FILE-B File containing the second Newick tree
- -h,--help Show this help text
-
-
-## Distancce
-
-Compute distances between many phylogenetic trees.
-
- stack exec tlynx -- distance --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx distance (-d|--distance MEASURE) [-n|--normalize] [-t|--intersect]
- [-s|--summary-statistics]
- [-m|--master-tree-file MASTER-TREE-File]
- [-f|--newick-format FORMAT] [INPUT-FILES]
- Compute distances between many phylogenetic trees.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -d,--distance MEASURE Type of distance to calculate (available distance
- measures are listed below)
- -n,--normalize Normalize trees before distance calculation; only
- affect distances depending on branch lengths
- -t,--intersect Compare intersections; i.e., before comparison, drop
- leaves that are not present in the other tree
- -s,--summary-statistics Report summary statistics only
- -m,--master-tree-file MASTER-TREE-File
- Compare all trees to the tree in the master tree
- file.
- -f,--newick-format FORMAT
- Newick tree format: Standard, IqTree, or RevBayes;
- default: Standard; for detailed help, see 'tlynx
- --help'
- INPUT-FILES Read tree(s) from INPUT-FILES; if more files are
- given, one tree is expected per file
- -h,--help Show this help text
-
- Distance measures:
- symmetric Symmetric distance (Robinson-Foulds distance).
- incompatible-split[VAL] Incompatible split distance. Collapse branches with (normalized)
- support less than 0.0<=VAL<=1.0 before distance calculation;
- if, let's say, VAL>0.7, only well supported differences contribute
- to the total distance.
- branch-score Branch score distance.
-
-
-## Examine
-
-Compute summary statistics of phylogenetic trees.
-
- stack exec tlynx -- examine --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx examine INPUT-FILE [-f|--newick-format FORMAT]
- Compute summary statistics of phylogenetic trees.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- INPUT-FILE Read trees from INPUT-FILE
- -f,--newick-format FORMAT
- Newick tree format: Standard, IqTree, or RevBayes;
- default: Standard; for detailed help, see 'tlynx
- --help'
- -h,--help Show this help text
-
-
-## Shuffle
-
-Shuffle a phylogenetic tree (keep coalescent times, but shuffle topology and
-leaves).
-
- stack exec tlynx -- shuffle --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx shuffle [-f|--newick-format FORMAT] [-n|--replicates N] TREE-FILE
- [-S|--seed [INT]]
- Shuffle a phylogenetic tree (keep coalescent times, but shuffle topology and leaves).
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -f,--newick-format FORMAT
- Newick tree format: Standard, IqTree, or RevBayes;
- default: Standard; for detailed help, see 'tlynx
- --help'
- -n,--replicates N Number of trees to generate
- TREE-FILE File containing a Newick tree
- -S,--seed [INT] Seed for random number generator; list of 32 bit
- integers with up to 256 elements (default: random)
- -h,--help Show this help text
-
-
-## Simulate
-
-Simulate phylogenetic trees using birth and death processes.
-
- stack exec tlynx -- simulate --help
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: tlynx simulate (-t|--nTrees INT) (-n|--nLeaves INT) PROCESS
- [-u|--sub-sample DOUBLE] [-s|--summary-statistics]
- [-S|--seed [INT]]
- Simulate phylogenetic trees using a birth and death or coalescent process.
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
- -t,--nTrees INT Number of trees
- -n,--nLeaves INT Number of leaves per tree
- -u,--sub-sample DOUBLE Perform sub-sampling; see below.
- -s,--summary-statistics For each branch, print length and number of children
- -S,--seed [INT] Seed for random number generator; list of 32 bit
- integers with up to 256 elements (default: random)
- -h,--help Show this help text
-
- Available processes:
- birthdeath Birth and death process
- coalescent Coalescent process
-
- See, for example, 'tlynx simulate birthdeath --help'.
- Sub-sample with probability p:
- 1. Simulate one big tree with n'=round(n/p), n'>=n, leaves;
- 2. Randomly sample sub-trees with n leaves.
- (With p=1.0, the same tree is reported over and over again.)
-
-
-# ELynx
-
-Validate and (optionally) redo past ELynx analyses.
-
- stack exec elynx -- --help | head -n -16
-
- ELynx Suite version 0.5.1.0.
- Developed by Dominik Schrempf.
- Compiled on June 12, 2021, at 14:54 pm, UTC.
-
- Usage: elynx COMMAND
- Validate and redo past ELynx analyses
-
- Available options:
- -h,--help Show this help text
- -V,--version Show version
-
- Available commands:
- validate Validate an ELynx analysis
- redo Redo an ELynx analysis
diff --git a/bench/Bench.hs b/bench/Bench.hs
--- a/bench/Bench.hs
+++ b/bench/Bench.hs
@@ -23,6 +23,7 @@
import Length
import Lens
import System.Random.MWC
+import Tree
treeFileMany :: FilePath
treeFileMany = "data/Many.trees"
@@ -45,8 +46,10 @@
main :: IO ()
main = do
!ts <- getManyTrees
- !ht <- first getLen <$> hugeTree
- let mr1 = hugeTreeCalcPar 0 ht
+ !ht <- first getLength <$> hugeTree
+ let !pht = lengthToPhyloTree ht
+ !ht' = flipLabels ht
+ mr1 = hugeTreeCalcPar 0 ht
mr2 = hugeTreeCalcPar 1 ht
if mr1 == mr2
then putStrLn "Map OK."
@@ -55,7 +58,7 @@
print mr2
error "Map wrong."
let fr1 = (foldl' (+) 0 . branches) ht
- fr2 = parBranchFoldMap 1 id (+) ht
+ fr2 = parBranchFoldMap 3 id (+) ht
if 1e-8 > abs (fr1 - fr2)
then putStrLn "Fold OK."
else do
@@ -90,6 +93,16 @@
bench "sum with accessor function" $ nf sumWithAccessorFunction [0 .. 1000000 :: Length],
bench "sum with setter and getter" $ nf sumWithSetter [0 .. 1000000 :: Length],
bench "sum with modify and accessor functions" $ nf sumWithModifyFunction [0 .. 1000000 :: Length]
+ ],
+ bgroup
+ "[mono|bi][functor|foldable|traversable]"
+ [ bench "normal fmap tree" $ nf fmapNormalFunctor ht',
+ bench "mono fmap tree" $ nf fmapFunctor ht,
+ bench "bi fmap tree" $ nf fmapBifunctor ht,
+ bench "mono fold tree" $ nf totalBranchLengthFoldable ht,
+ bench "bi fold tree" $ nf totalBranchLengthBifoldable ht,
+ bench "mono traverse tree" $ nf toLengthTreeTraversable pht,
+ bench "bi traverse tree" $ nf toLengthTreeBitraversable pht
]
]
@@ -128,3 +141,120 @@
-- mean 368.3 ms (345.6 ms .. 394.6 ms)
-- std dev 25.80 ms (7.999 ms .. 35.08 ms)
-- variance introduced by outliers: 20% (moderately inflated)
+
+-- benchmarking length/length sum foldl' safe
+-- time 123.4 ms (117.8 ms .. 136.4 ms)
+-- 0.993 R² (0.986 R² .. 1.000 R²)
+-- mean 120.3 ms (118.1 ms .. 124.8 ms)
+-- std dev 4.572 ms (2.556 ms .. 6.511 ms)
+-- variance introduced by outliers: 11% (moderately inflated)
+
+-- benchmarking length/length sum foldl' unsafe
+-- time 82.04 ms (77.92 ms .. 86.95 ms)
+-- 0.994 R² (0.988 R² .. 0.999 R²)
+-- mean 88.56 ms (85.87 ms .. 92.70 ms)
+-- std dev 6.079 ms (4.635 ms .. 7.548 ms)
+-- variance introduced by outliers: 19% (moderately inflated)
+
+-- benchmarking length/length sum foldl' num instance
+-- time 102.4 ms (89.16 ms .. 124.8 ms)
+-- 0.964 R² (0.896 R² .. 0.994 R²)
+-- mean 98.94 ms (91.15 ms .. 106.4 ms)
+-- std dev 11.90 ms (8.779 ms .. 15.29 ms)
+-- variance introduced by outliers: 32% (moderately inflated)
+
+-- benchmarking length/double sum foldl'
+-- time 83.87 ms (71.82 ms .. 104.0 ms)
+-- 0.954 R² (0.878 R² .. 0.996 R²)
+-- mean 103.1 ms (95.66 ms .. 112.2 ms)
+-- std dev 14.21 ms (10.77 ms .. 17.84 ms)
+-- variance introduced by outliers: 43% (moderately inflated)
+
+-- benchmarking length/double sum
+-- time 102.6 ms (88.13 ms .. 116.0 ms)
+-- 0.970 R² (0.938 R² .. 0.997 R²)
+-- mean 94.08 ms (89.94 ms .. 100.9 ms)
+-- std dev 8.071 ms (5.233 ms .. 10.58 ms)
+-- variance introduced by outliers: 21% (moderately inflated)
+
+-- benchmarking lens/sum with getter
+-- time 4.320 ms (4.303 ms .. 4.352 ms)
+-- 1.000 R² (0.999 R² .. 1.000 R²)
+-- mean 4.328 ms (4.320 ms .. 4.346 ms)
+-- std dev 38.06 μs (18.94 μs .. 58.05 μs)
+
+-- benchmarking lens/sum with accessor function
+-- time 4.314 ms (4.073 ms .. 4.509 ms)
+-- 0.912 R² (0.795 R² .. 0.986 R²)
+-- mean 4.795 ms (4.546 ms .. 5.264 ms)
+-- std dev 1.056 ms (562.7 μs .. 1.691 ms)
+-- variance introduced by outliers: 90% (severely inflated)
+
+-- benchmarking lens/sum with setter and getter
+-- time 16.34 ms (14.93 ms .. 18.00 ms)
+-- 0.955 R² (0.913 R² .. 0.983 R²)
+-- mean 17.06 ms (16.30 ms .. 18.25 ms)
+-- std dev 2.137 ms (1.656 ms .. 2.794 ms)
+-- variance introduced by outliers: 59% (severely inflated)
+
+-- benchmarking lens/sum with modify and accessor functions
+-- time 14.87 ms (14.11 ms .. 15.68 ms)
+-- 0.977 R² (0.946 R² .. 0.993 R²)
+-- mean 15.54 ms (15.15 ms .. 16.23 ms)
+-- std dev 1.319 ms (882.4 μs .. 2.030 ms)
+-- variance introduced by outliers: 40% (moderately inflated)
+
+-- Mono vs bi.
+--
+-- The Bifoldable instance is much slower.
+--
+-- The Bifunctor instance is as fast as the normal Functor instance.
+--
+-- The Bitraversable instance is fine, although slightly slower.
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/normal fmap tree
+-- time 30.57 ms (29.87 ms .. 31.16 ms)
+-- 0.999 R² (0.998 R² .. 1.000 R²)
+-- mean 29.84 ms (29.66 ms .. 30.23 ms)
+-- std dev 528.5 μs (358.8 μs .. 725.0 μs)
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/mono fmap tree
+-- time 28.98 ms (28.54 ms .. 29.36 ms)
+-- 0.999 R² (0.999 R² .. 1.000 R²)
+-- mean 28.73 ms (28.44 ms .. 29.40 ms)
+-- std dev 901.2 μs (411.2 μs .. 1.583 ms)
+-- variance introduced by outliers: 10% (moderately inflated)
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/bi fmap tree
+-- time 29.53 ms (28.27 ms .. 31.23 ms)
+-- 0.993 R² (0.984 R² .. 1.000 R²)
+-- mean 28.93 ms (28.57 ms .. 29.55 ms)
+-- std dev 1.103 ms (436.8 μs .. 1.577 ms)
+-- variance introduced by outliers: 10% (moderately inflated)
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/mono fold tree
+-- time 81.88 ms (80.24 ms .. 83.03 ms)
+-- 0.999 R² (0.996 R² .. 1.000 R²)
+-- mean 83.68 ms (82.70 ms .. 85.33 ms)
+-- std dev 2.207 ms (581.3 μs .. 2.915 ms)
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/bi fold tree
+-- time 124.1 ms (120.7 ms .. 129.4 ms)
+-- 0.999 R² (0.996 R² .. 1.000 R²)
+-- mean 123.9 ms (121.8 ms .. 127.5 ms)
+-- std dev 4.128 ms (2.139 ms .. 6.374 ms)
+-- variance introduced by outliers: 11% (moderately inflated)
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/mono traverse tree
+-- time 118.5 ms (115.2 ms .. 122.6 ms)
+-- 0.999 R² (0.995 R² .. 1.000 R²)
+-- mean 122.2 ms (120.4 ms .. 126.8 ms)
+-- std dev 4.086 ms (1.704 ms .. 5.961 ms)
+-- variance introduced by outliers: 11% (moderately inflated)
+
+-- benchmarking [mono|bi][functor|foldable|traversable]/bi traverse tree
+-- time 122.6 ms (116.8 ms .. 127.2 ms)
+-- 0.999 R² (0.997 R² .. 1.000 R²)
+-- mean 121.3 ms (120.2 ms .. 123.0 ms)
+-- std dev 2.377 ms (1.043 ms .. 3.289 ms)
+-- variance introduced by outliers: 11% (moderately inflated)
diff --git a/bench/Length.hs b/bench/Length.hs
--- a/bench/Length.hs
+++ b/bench/Length.hs
@@ -15,11 +15,19 @@
lengthSumFoldl'NumInstance,
doubleSumFoldl',
doubleSum,
+ fmapNormalFunctor,
+ fmapFunctor,
+ fmapBifunctor,
+ totalBranchLengthFoldable,
+ totalBranchLengthBifoldable,
)
where
+import Data.Bifoldable
+import Data.Bifunctor
import Data.Foldable
import ELynx.Tree.Length
+import ELynx.Tree.Rooted
lengthSumFoldl' :: [Length] -> Length
lengthSumFoldl' =
@@ -42,3 +50,18 @@
doubleSum :: [Double] -> Double
doubleSum = foldl' (+) 0
+
+fmapNormalFunctor :: HasLength a => Tree e a -> Tree e a
+fmapNormalFunctor = fmap (modifyLength cos)
+
+fmapFunctor :: HasLength e => Tree e a -> Tree e a
+fmapFunctor = getZipBranchTree . fmap (modifyLength cos) . ZipBranchTree
+
+fmapBifunctor :: HasLength e => Tree e a -> Tree e a
+fmapBifunctor = first (modifyLength cos)
+
+totalBranchLengthFoldable :: HasLength e => Tree e a -> Length
+totalBranchLengthFoldable = totalBranchLength
+
+totalBranchLengthBifoldable :: HasLength e => Tree e a -> Length
+totalBranchLengthBifoldable = bifoldl' (+) const 0 . first getLength
diff --git a/bench/Lens.hs b/bench/Lens.hs
--- a/bench/Lens.hs
+++ b/bench/Lens.hs
@@ -22,7 +22,7 @@
import Lens.Micro
len :: HasLength a => Lens' a Length
-len = lens getLen (flip setLen)
+len = lens getLength (flip setLength)
sumWithGetter :: HasLength a => [a] -> Length
sumWithGetter = foldl' (\x y -> x ^. len + y ^. len) 0
@@ -31,7 +31,7 @@
sumWithSetter = sumWithGetter . map (\x -> x & len %~ (+ 10))
sumWithAccessorFunction :: HasLength a => [a] -> Length
-sumWithAccessorFunction = foldl' (\x y -> getLen x + getLen y) 0
+sumWithAccessorFunction = foldl' (\x y -> getLength x + getLength y) 0
sumWithModifyFunction :: HasLength a => [a] -> Length
-sumWithModifyFunction = sumWithAccessorFunction . map (modLen (+ 10))
+sumWithModifyFunction = sumWithAccessorFunction . map (modifyLength (+ 10))
diff --git a/bench/Tree.hs b/bench/Tree.hs
new file mode 100644
--- /dev/null
+++ b/bench/Tree.hs
@@ -0,0 +1,38 @@
+-- |
+-- Module : Tree
+-- Description : Benchmark some tree algorithms
+-- Copyright : (c) 2021 Dominik Schrempf
+-- License : GPL-3.0-or-later
+--
+-- Maintainer : dominik.schrempf@gmail.com
+-- Stability : experimental
+-- Portability : portable
+--
+-- Creation date: Fri Jul 2 13:10:48 2021.
+module Tree
+ ( toLengthTreeTraversable,
+ toLengthTreeBitraversable,
+ )
+where
+
+import Data.Bifunctor
+import Data.Bitraversable
+import ELynx.Tree.Length
+import ELynx.Tree.Phylogeny
+import ELynx.Tree.Rooted
+
+-- Compare speed of bitraversal and special traversal instances.
+
+toLengthTreeTraversable :: Tree Phylo a -> Tree Length a
+toLengthTreeTraversable = either error id . toLengthTree
+
+fromMaybeWithError :: String -> Maybe a -> Either String a
+fromMaybeWithError s = maybe (Left s) Right
+
+toLengthTreeBitraversable :: HasMaybeLength e => Tree e a -> Either String (Tree Length a)
+toLengthTreeBitraversable t =
+ fromMaybeWithError "toLengthTree: Length unavailable for some branches." $ bisequenceA t'
+ where
+ t' = modifyStem cleanLength $ bimap getMaybeLength pure t
+ cleanLength Nothing = pure 0
+ cleanLength x = x
diff --git a/elynx-tree.cabal b/elynx-tree.cabal
--- a/elynx-tree.cabal
+++ b/elynx-tree.cabal
@@ -1,19 +1,18 @@
-cabal-version: 2.2
-name: elynx-tree
-version: 0.5.1.1
-license: GPL-3.0-or-later
-license-file: LICENSE
-copyright: Dominik Schrempf (2021)
-maintainer: dominik.schrempf@gmail.com
-author: Dominik Schrempf
-homepage: https://github.com/dschrempf/elynx#readme
-bug-reports: https://github.com/dschrempf/elynx/issues
-synopsis: Handle phylogenetic trees
-description:
- Examine, compare, and simulate phylogenetic trees in a reproducible way. Please see the README on GitHub at .
+cabal-version: 2.2
+name: elynx-tree
+version: 0.6.0.0
+synopsis: Handle phylogenetic trees
+description: Examine, compare, and simulate phylogenetic trees in a reproducible way. Please see the README on GitHub at .
+category: Bioinformatics
+homepage: https://github.com/dschrempf/elynx#readme
+bug-reports: https://github.com/dschrempf/elynx/issues
+author: Dominik Schrempf
+maintainer: dominik.schrempf@gmail.com
+copyright: Dominik Schrempf (2021)
+license: GPL-3.0-or-later
+license-file: LICENSE
+build-type: Simple
-category: Bioinformatics
-build-type: Simple
extra-source-files:
README.md
ChangeLog.md
@@ -35,111 +34,118 @@
data/UltraMetric.tree
source-repository head
- type: git
- location: https://github.com/dschrempf/elynx
+ type: git
+ location: https://github.com/dschrempf/elynx
library
- exposed-modules:
- ELynx.Topology
- ELynx.Topology.Phylogeny
- ELynx.Topology.Rooted
- ELynx.Tree
- ELynx.Tree.Bipartition
- ELynx.Tree.Distance
- ELynx.Tree.Length
- ELynx.Tree.Name
- ELynx.Tree.Parallel
- ELynx.Tree.Partition
- ELynx.Tree.Phylogeny
- ELynx.Tree.Rooted
- ELynx.Tree.Splittable
- ELynx.Tree.Support
- ELynx.Tree.Zipper
- ELynx.Tree.Distribution.BirthDeath
- ELynx.Tree.Distribution.BirthDeathCritical
- ELynx.Tree.Distribution.BirthDeathCriticalNoTime
- ELynx.Tree.Distribution.BirthDeathNearlyCritical
- ELynx.Tree.Distribution.CoalescentContinuous
- ELynx.Tree.Distribution.TimeOfOrigin
- ELynx.Tree.Distribution.TimeOfOriginNearCritical
- ELynx.Tree.Distribution.Types
- ELynx.Tree.Export.Newick
- ELynx.Tree.Export.Nexus
- ELynx.Tree.Import.Newick
- ELynx.Tree.Import.Nexus
- ELynx.Tree.Simulate.Coalescent
- ELynx.Tree.Simulate.PointProcess
-
- hs-source-dirs: src
- other-modules: Paths_elynx_tree
- autogen-modules: Paths_elynx_tree
- default-language: Haskell2010
- ghc-options: -Wall -Wunused-packages
- build-depends:
- aeson >=1.5.6.0,
- attoparsec >=0.13.2.5,
- base >=4.14.1.0,
- bytestring >=0.10.12.0,
- comonad >=5.0.8,
- containers >=0.6.2.1,
- deepseq >=1.4.4.0,
- double-conversion >=2.0.2.0,
- elynx-nexus >=0.5.1.1,
- math-functions >=0.3.4.2,
- mwc-random >=0.15.0.1,
- parallel >=3.2.2.0,
- primitive >=0.7.1.0,
- statistics >=0.15.2.0
+ exposed-modules:
+ ELynx.Topology
+ ELynx.Topology.Phylogeny
+ ELynx.Topology.Rooted
+ ELynx.Tree
+ ELynx.Tree.Bipartition
+ ELynx.Tree.Distance
+ ELynx.Tree.Length
+ ELynx.Tree.Name
+ ELynx.Tree.Parallel
+ ELynx.Tree.Partition
+ ELynx.Tree.Phylogeny
+ ELynx.Tree.Rooted
+ ELynx.Tree.Splittable
+ ELynx.Tree.Support
+ ELynx.Tree.Zipper
+ ELynx.Tree.Distribution.BirthDeath
+ ELynx.Tree.Distribution.BirthDeathCritical
+ ELynx.Tree.Distribution.BirthDeathCriticalNoTime
+ ELynx.Tree.Distribution.BirthDeathNearlyCritical
+ ELynx.Tree.Distribution.CoalescentContinuous
+ ELynx.Tree.Distribution.TimeOfOrigin
+ ELynx.Tree.Distribution.TimeOfOriginNearCritical
+ ELynx.Tree.Distribution.Types
+ ELynx.Tree.Export.Newick
+ ELynx.Tree.Export.Nexus
+ ELynx.Tree.Import.Newick
+ ELynx.Tree.Import.Nexus
+ ELynx.Tree.Simulate.Coalescent
+ ELynx.Tree.Simulate.PointProcess
+ other-modules:
+ Paths_elynx_tree
+ autogen-modules:
+ Paths_elynx_tree
+ hs-source-dirs: src
+ ghc-options: -Wall -Wunused-packages
+ build-depends:
+ aeson
+ , attoparsec
+ , base
+ , bytestring
+ , comonad
+ , containers
+ , data-default-class
+ , deepseq
+ , double-conversion
+ , elynx-nexus
+ , math-functions
+ , mwc-random
+ , parallel
+ , primitive
+ , statistics
+ default-language: Haskell2010
test-suite tree-test
- type: exitcode-stdio-1.0
- main-is: Spec.hs
- hs-source-dirs: test
- other-modules:
- ELynx.Topology.RootedSpec
- ELynx.Tree.Arbitrary
- ELynx.Tree.BipartitionSpec
- ELynx.Tree.DistanceSpec
- ELynx.Tree.PartitionSpec
- ELynx.Tree.PhylogenySpec
- ELynx.Tree.RootedSpec
- ELynx.Tree.SupportSpec
- ELynx.Tree.Export.NewickSpec
- ELynx.Tree.Export.NexusSpec
- ELynx.Tree.Import.NewickSpec
- ELynx.Tree.Import.NexusSpec
- Paths_elynx_tree
-
- default-language: Haskell2010
- ghc-options: -Wall -Wunused-packages
- build-depends:
- QuickCheck >=2.14.2,
- attoparsec >=0.13.2.5,
- base >=4.14.1.0,
- bytestring >=0.10.12.0,
- containers >=0.6.2.1,
- elynx-tools >=0.5.1.1,
- elynx-tree -any,
- hspec >=2.7.10
+ type: exitcode-stdio-1.0
+ main-is: Spec.hs
+ other-modules:
+ ELynx.ClassLaws
+ ELynx.Topology.Arbitrary
+ ELynx.Topology.RootedSpec
+ ELynx.Tree.Arbitrary
+ ELynx.Tree.BipartitionSpec
+ ELynx.Tree.DistanceSpec
+ ELynx.Tree.PartitionSpec
+ ELynx.Tree.PhylogenySpec
+ ELynx.Tree.RootedSpec
+ ELynx.Tree.SupportSpec
+ ELynx.Tree.Export.NewickSpec
+ ELynx.Tree.Export.NexusSpec
+ ELynx.Tree.Import.NewickSpec
+ ELynx.Tree.Import.NexusSpec
+ Paths_elynx_tree
+ autogen-modules:
+ Paths_elynx_tree
+ hs-source-dirs: test
+ ghc-options: -Wall -Wunused-packages -threaded -rtsopts -with-rtsopts=-N
+ build-depends:
+ attoparsec
+ , base
+ , bytestring
+ , containers
+ , data-default
+ , elynx-tools
+ , elynx-tree
+ , hspec
+ , QuickCheck
+ , quickcheck-classes
+ default-language: Haskell2010
benchmark tree-bench
- type: exitcode-stdio-1.0
- main-is: Bench.hs
- hs-source-dirs: bench
- other-modules:
- Length
- Lens
- Paths_elynx_tree
-
- default-language: Haskell2010
- ghc-options:
- -Wall -Wunused-packages -threaded -rtsopts -with-rtsopts=-N
-
- build-depends:
- base >=4.7 && <5,
- criterion >=1.5.9.0,
- elynx-tools >=0.5.1.1,
- elynx-tree -any,
- microlens >=0.4.12.0,
- mwc-random >=0.15.0.1,
- parallel >=3.2.2.0
+ type: exitcode-stdio-1.0
+ main-is: Bench.hs
+ other-modules:
+ Length
+ Lens
+ Tree
+ Paths_elynx_tree
+ autogen-modules:
+ Paths_elynx_tree
+ hs-source-dirs: bench
+ ghc-options: -Wall -Wunused-packages -threaded -rtsopts -with-rtsopts=-N
+ build-depends:
+ base >=4.7 && <5
+ , criterion
+ , elynx-tools
+ , elynx-tree
+ , microlens
+ , mwc-random
+ , parallel
+ default-language: Haskell2010
diff --git a/src/ELynx/Topology/Phylogeny.hs b/src/ELynx/Topology/Phylogeny.hs
--- a/src/ELynx/Topology/Phylogeny.hs
+++ b/src/ELynx/Topology/Phylogeny.hs
@@ -10,46 +10,36 @@
--
-- Creation date: Sat Jul 18 13:15:49 2020.
--
--- THIS MODULE IS INCOMPLETE.
---
-- A topology, as it is used in phylogenetics is a 'Topology' with unique leaf
-- labels, and the order of the topologies in the sub-forest is considered to be
-- meaningless.
--
--- Internally, however, the underlying 'Topology' data structure stores the
--- sub-forest as a (non-empty) list, which has a specific order. Hence, we have
--- to do some tricks when comparing topologies, and topology comparison is slow.
---
--- Also, the uniqueness of the leaves is not ensured by the data type, but has
--- to be checked at runtime. Functions relying on the tree to have unique leaves
--- do perform this check, and return 'Left' with an error message, if the tree
--- has duplicate leaves.
---
--- Note: 'Topology's are rooted.
---
--- Note: 'Topology's encoded in Newick format correspond to rooted topologies.
--- By convention only, a topology parsed from Newick format is usually thought
--- to be unrooted, when the root node is multifurcating and has three or more
--- children. This convention is not enforced here. Newick topologies are just
--- parsed as they are, and a rooted topology is returned.
---
--- THIS MODULE IS INCOMPLETE.
+-- NOTE: The functions in this module are defined using the functions in
+-- "ELynx.Tree.Phylogeny". This induces a runtime overhead, but greatly reduces
+-- the probability of additional bugs.
module ELynx.Topology.Phylogeny
( equal,
equal',
+ intersect,
+ bifurcating,
+ outgroup,
+ midpoint,
+ roots,
)
where
-import Data.List
+import Data.Default.Class
+import Data.List hiding (intersect)
import Data.Maybe
+import qualified Data.Set as S
import ELynx.Topology.Rooted
+import ELynx.Tree.Length
+import qualified ELynx.Tree.Phylogeny as T
+import qualified ELynx.Tree.Rooted as T
-- | The equality check is slow because the order of children is considered to
-- be arbitrary.
--
--- NOTE: The equality check is only meaningful if the topologies have unique
--- leaves.
---
-- Return 'Left' if a topology does not have unique leaves.
equal :: (Eq a, Ord a) => Topology a -> Topology a -> Either String Bool
equal tL tR
@@ -59,7 +49,8 @@
-- | Same as 'equal', but assume that leaves are unique.
equal' :: Eq a => Topology a -> Topology a -> Bool
-equal' (Leaf lbL) (Leaf lbR) = lbL == lbR
+equal' (Leaf lbL) (Leaf lbR) =
+ lbL == lbR
equal' (Node tsL) (Node tsR) =
(length tsL == length tsR)
&& all (`elem'` tsR) tsL
@@ -67,12 +58,70 @@
elem' t ts = isJust $ find (equal' t) ts
equal' _ _ = False
--- TODO.
+-- | Intersection of topologies.
+--
+-- | See 'ELynx.Tree.Phylogeny.intersect'.
+intersect ::
+ (Ord a) => [Topology a] -> Either String [Topology a]
+intersect ts
+ | S.null lvsCommon = Left "intersect: Intersection of leaves is empty."
+ | otherwise = case sequence [dropLeavesWith (predicate ls) t | (ls, t) <- zip leavesToDrop ts] of
+ Nothing -> Left "intersect: A topology is empty."
+ Just ts' -> Right ts'
+ where
+ -- Leaf sets.
+ lvss = map (S.fromList . leaves) ts
+ -- Common leaf set.
+ lvsCommon = foldl1' S.intersection lvss
+ -- Leaves to drop for each topology in the forest.
+ leavesToDrop = map (S.\\ lvsCommon) lvss
+ -- Predicate.
+ predicate lvsToDr l = l `S.member` lvsToDr
--- A multifurcating root node can be resolved to a bifurcating root node with
--- 'outgroup'.
+-- | Check if topology is bifurcating.
--
--- The bifurcating root node can be changed with 'outgroup' or 'midpoint'.
+-- | See 'ELynx.Tree.Phylogeny.intersect'.
+bifurcating :: Topology a -> Bool
+bifurcating (Leaf _) = True
+bifurcating (Node ts) = (length ts == 2) && all bifurcating ts
+
+-- Perform a computation over the 'Tree' data type.
+overTree ::
+ (Default a, Functor f) =>
+ (T.Tree Length a -> f (T.Tree Length a)) ->
+ Topology a ->
+ f (Topology a)
+overTree f = goBack . f . goThere
+ where
+ goThere = toBranchLabelTreeWith (toLengthUnsafe 1.0) def
+ goBack = fmap fromBranchLabelTree
+
+-- | Root topology using an outgroup.
--
--- For a given topology with bifurcating root node, a list of all rooted
--- topologies is returned by 'roots'.
+-- See 'ELynx.Tree.Phylogeny.outgroup'.
+outgroup :: (Default a, Ord a) => S.Set a -> Topology a -> Either String (Topology a)
+outgroup xs = overTree (T.outgroup xs)
+
+-- | Root topology at the midpoint.
+--
+-- See 'ELynx.Tree.Phylogeny.midpoint'.
+--
+-- Use 'depth' to measure topology height.
+--
+-- If the midpoint is ambiguous because the sum of the left and right depths is
+-- odd, the depth of the left sub-topology will be set to be one node greater
+-- than the one of the right sub-topology.
+midpoint :: Default a => Topology a -> Either String (Topology a)
+midpoint = overTree T.midpoint
+
+-- | For a rooted tree with a bifurcating root node, get all possible rooted
+-- trees.
+--
+-- See 'ELynx.Tree.Phylogeny.roots'.
+roots :: Default a => Topology a -> Either String [Topology a]
+roots = goBack . T.roots . goThere
+ where
+ -- We have to use a special 'overTree' function here, since a list of
+ -- topologies is returned.
+ goThere = toBranchLabelTreeWith (toLengthUnsafe 1.0) def
+ goBack = (fmap . fmap) fromBranchLabelTree
diff --git a/src/ELynx/Topology/Rooted.hs b/src/ELynx/Topology/Rooted.hs
--- a/src/ELynx/Topology/Rooted.hs
+++ b/src/ELynx/Topology/Rooted.hs
@@ -13,34 +13,31 @@
--
-- Creation date: Sat Jul 11 10:28:28 2020.
--
--- THIS MODULE IS INCOMPLETE.
---
-- A rooted 'Topology' differs from a classical rooted rose 'Data.Tree.Tree' in
-- that it does not have internal node labels. The leaves have labels.
--
--- For rooted trees with branch labels, please see "ELynx.Tree.Rooted". Please
--- also see the note about tree traversals therein.
---
--- THIS MODULE IS INCOMPLETE.
+-- For rooted trees with branch labels, see "ELynx.Tree.Rooted".
module ELynx.Topology.Rooted
( -- * Data type
Topology (..),
Forest,
- fromTree,
- fromLabeledTree,
- toLabeledTreeWith,
+ fromRoseTree,
+ fromBranchLabelTree,
+ toBranchLabelTreeWith,
-- * Access leaves, branches and labels
leaves,
duplicateLeaves,
+ setLeaves,
identify,
-- * Structure
degree,
+ depth,
prune,
dropLeavesWith,
- zipTreesWith,
- zipTrees,
+ zipTopologiesWith,
+ zipTopologies,
)
where
@@ -50,6 +47,7 @@
import Data.Aeson
import Data.Data
import Data.Foldable
+import Data.Functor
import Data.List.NonEmpty (NonEmpty)
import qualified Data.List.NonEmpty as N
import Data.Maybe
@@ -68,7 +66,7 @@
| Leaf {label :: a}
deriving (Eq, Read, Show, Data, Generic)
--- | A shortcut.
+-- | Shorthand.
type Forest a = NonEmpty (Topology a)
instance Functor Topology where
@@ -89,30 +87,26 @@
traverse g (Node ts) = Node <$> traverse (traverse g) ts
traverse g (Leaf lb) = Leaf <$> g lb
--- TODO: This type checks, but I doubt the implementation is bug-free.
instance Applicative Topology where
pure = Leaf
- (Node tsF) <*> tx = Node $ fmap (<*> tx) tsF
+ (Node tsF) <*> tx = Node $ tsF <&> (<*> tx)
(Leaf lbF) <*> tx = lbF <$> tx
liftA2 f (Node tsX) ty = Node $ fmap (\tx -> liftA2 f tx ty) tsX
liftA2 f (Leaf lbX) (Node tsY) = Node $ fmap (f lbX <$>) tsY
liftA2 f (Leaf lbX) (Leaf lbY) = Leaf $ f lbX lbY
- (Node tsX) *> ty@(Node tsY) = Node $ tsY <> fmap (*> ty) tsX
+ (Node tsX) *> ty = Node $ tsX <&> (*> ty)
(Leaf _) *> (Node tsY) = Node tsY
- _ *> (Leaf lbY) = Leaf lbY
+ (Leaf _) *> (Leaf y) = Leaf y
- (Node tsX) <* ty = Node $ fmap (<* ty) tsX
- (Leaf lbX) <* _ = Leaf lbX
+ (Node tsX) <* ty = Node $ tsX <&> (<* ty)
+ (Leaf x) <* ty = x <$ ty
--- TODO: This type checks, but I doubt the implementation is bug-free.
instance Monad Topology where
(Node ts) >>= f = Node $ fmap (>>= f) ts
- (Leaf lb) >>= f = case f lb of
- Node ts' -> Node ts'
- Leaf lb' -> Leaf lb'
+ (Leaf lb) >>= f = f lb
instance NFData a => NFData (Topology a) where
rnf (Node ts) = rnf ts
@@ -122,32 +116,29 @@
instance FromJSON a => FromJSON (Topology a)
--- TODO: Provide and fix tests, provide arbitrary instances.
-
-- | Convert a rooted rose tree to a rooted topology. Internal node labels are lost.
-fromTree :: T.Tree a -> Topology a
-fromTree (T.Node lb []) = Leaf lb
-fromTree (T.Node _ xs) = Node $ fromTree <$> N.fromList xs
+fromRoseTree :: T.Tree a -> Topology a
+fromRoseTree (T.Node lb []) = Leaf lb
+fromRoseTree (T.Node _ xs) = Node $ fromRoseTree <$> N.fromList xs
--- | Convert a rooted, branch-labeled rose tree to a rooted topology. Branch
--- labels and internal node labels are lost.
-fromLabeledTree :: R.Tree e a -> Topology a
-fromLabeledTree (R.Node _ lb []) = Leaf lb
-fromLabeledTree (R.Node _ _ xs) = Node $ fromLabeledTree <$> N.fromList xs
+-- | Convert a rooted, branch-label tree to a rooted topology. Branch labels and
+-- internal node labels are lost.
+fromBranchLabelTree :: R.Tree e a -> Topology a
+fromBranchLabelTree (R.Node _ lb []) = Leaf lb
+fromBranchLabelTree (R.Node _ _ xs) = Node $ fromBranchLabelTree <$> N.fromList xs
--- | Convert a rooted topology to a rooted, branch-labeled rose tree. Use the
--- given node label at internal nodes.
-toLabeledTreeWith :: a -> Topology a -> R.Tree () a
-toLabeledTreeWith _ (Leaf lb) = R.Node () lb []
-toLabeledTreeWith x (Node ts) = R.Node () x $ map (toLabeledTreeWith x) $ N.toList ts
+-- | Convert a rooted topology to a rooted, branch-label tree. Use the given
+-- node label at internal nodes.
+toBranchLabelTreeWith :: e -> a -> Topology a -> R.Tree e a
+toBranchLabelTreeWith b _ (Leaf lb) = R.Node b lb []
+toBranchLabelTreeWith b l (Node ts) = R.Node b l $ map (toBranchLabelTreeWith b l) $ N.toList ts
--- TODO: Maybe use foldr similar to 'flatten'.
--- | Set of leaves.
+-- | List of leaves.
leaves :: Topology a -> [a]
leaves (Leaf lb) = [lb]
leaves (Node ts) = concatMap leaves ts
--- -- TODO: Check if this implementation of 'leaves' is faster.
+-- -- NOTE: This implementation of 'leaves' may be faster.
-- -- | Return leaf labels in pre-order.
-- flatten :: Topology a -> [a]
-- flatten t = squish t []
@@ -165,8 +156,16 @@
duplicateLeaves :: Ord a => Topology a -> Bool
duplicateLeaves = duplicates . leaves
--- TODO: This is the same as in ELynx.Tree.Rooted.
--- | Label the leaves with unique integers starting at 0.
+-- | Set leaf labels in pre-order.
+--
+-- Return 'Nothing' if the provided list of leaf labels is too short.
+setLeaves :: Traversable t => [b] -> t a -> Maybe (t b)
+setLeaves xs = sequenceA . snd . mapAccumL setLeafM xs
+ where
+ setLeafM [] _ = ([], Nothing)
+ setLeafM (y : ys) _ = (ys, Just y)
+
+-- | Label the leaves in pre-order with unique indices starting at 0.
identify :: Traversable t => t a -> t Int
identify = snd . mapAccumL (\i _ -> (i + 1, i)) (0 :: Int)
@@ -175,7 +174,18 @@
degree (Node ts) = (+ 1) $ length ts
degree (Leaf _) = 1
+-- | Depth of a topology.
+--
+-- See 'ELynx.Tree.Rooted.depth'.
+depth :: Topology a -> Int
+depth = maximum . go 1
+ where
+ go n (Leaf _) = [n]
+ go n (Node xs) = concatMap (go (n + 1)) xs
+
-- | Prune degree two nodes.
+--
+-- See 'ELynx.Tree.Rooted.prune'.
prune :: Topology a -> Topology a
prune (Node ts)
| singleton ts = Node $ fmap prune $ forest $ N.head ts
@@ -184,9 +194,7 @@
-- | Drop leaves satisfying predicate.
--
--- Degree two nodes may arise.
---
--- Return 'Nothing' if all leaves satisfy the predicate.
+-- See 'ELynx.Tree.Rooted.dropNodesWith'.
dropLeavesWith :: (a -> Bool) -> Topology a -> Maybe (Topology a)
dropLeavesWith p (Leaf lb)
| p lb = Nothing
@@ -194,25 +202,23 @@
dropLeavesWith p (Node ts) =
if null ts'
then Nothing
- else -- XXX: May be slow, unnecessary conversion to and from list.
+ else -- NOTE: Unnecessary conversion to and from list?
Just $ Node $ N.fromList ts'
where
ts' = catMaybes $ N.toList $ fmap (dropLeavesWith p) ts
-- | Zip leaves of two equal topologies.
--
--- Return 'Nothing' if the topologies are different.
-zipTreesWith :: (a1 -> a2 -> a) -> Topology a1 -> Topology a2 -> Maybe (Topology a)
-zipTreesWith f (Node tsL) (Node tsR) =
+-- See 'ELynx.Tree.Rooted.zipTreesWith'.
+zipTopologiesWith :: (a1 -> a2 -> a) -> Topology a1 -> Topology a2 -> Maybe (Topology a)
+zipTopologiesWith f (Node tsL) (Node tsR) =
if N.length tsL == N.length tsR
- then -- XXX: May be slow, unnecessary conversion to and from list.
- zipWithM (zipTreesWith f) (N.toList tsL) (N.toList tsR) >>= Just . Node . N.fromList
+ then -- NOTE: Unnecessary conversion to and from list?
+ zipWithM (zipTopologiesWith f) (N.toList tsL) (N.toList tsR) >>= Just . Node . N.fromList
else Nothing
-zipTreesWith f (Leaf lbL) (Leaf lbR) = Just $ Leaf $ f lbL lbR
-zipTreesWith _ _ _ = Nothing
+zipTopologiesWith f (Leaf lbL) (Leaf lbR) = Just $ Leaf $ f lbL lbR
+zipTopologiesWith _ _ _ = Nothing
--- | Zip leaves of two equal topologies.
---
--- Return 'Nothing' if the topologies are different.
-zipTrees :: Topology a1 -> Topology a2 -> Maybe (Topology (a1, a2))
-zipTrees = zipTreesWith (,)
+-- | See 'zipTopologiesWith'.
+zipTopologies :: Topology a1 -> Topology a2 -> Maybe (Topology (a1, a2))
+zipTopologies = zipTopologiesWith (,)
diff --git a/src/ELynx/Tree/Bipartition.hs b/src/ELynx/Tree/Bipartition.hs
--- a/src/ELynx/Tree/Bipartition.hs
+++ b/src/ELynx/Tree/Bipartition.hs
@@ -11,18 +11,6 @@
-- Portability : portable
--
-- Creation date: Fri Aug 30 15:28:17 2019.
---
--- 'Bipartition's are weird in that
--- > Bipartition x y == Bipartition y x
--- is True.
---
--- Also,
--- > Bipartition x y > Bipartition y x
--- is False, even when @x > y@.
---
--- That's why we have to make sure that for
--- > Bipartition x y
--- we always have @x >= y@.
module ELynx.Tree.Bipartition
( groups,
@@ -59,13 +47,34 @@
-- | A bipartition of a tree is a grouping of the leaves of the tree into two
-- non-overlapping, non-empty sub sets.
--
--- For example, each branch of a tree partitions the leaves of the tree into two
--- subsets, or a bipartition. Also, a bifurcating root induces a bipartition;
--- see 'bipartition'.
+-- For unrooted trees:
--
--- The order of the two subsets of a 'Bipartition' is meaningless. We ensure by
--- construction that the smaller subset comes first, and hence, that equality
--- checks are meaningful.
+-- - Each branch partitions the leaves of the tree into two subsets, or a
+-- bipartition.
+--
+-- For rooted trees:
+--
+-- - A bifurcating root node induces a bipartition; see 'bipartition'.
+--
+-- - Each inner node induces a bipartition by taking the leaves of the sub tree
+-- and the complement leaf set of the full tree.
+--
+-- The order of the two subsets of a 'Bipartition' is meaningless. That is,
+-- 'Bipartition's are weird in that
+--
+-- > Bipartition x y == Bipartition y x
+--
+-- is 'True'. Also,
+--
+-- > Bipartition x y > Bipartition y x
+--
+-- is False, even when @x > y@. That's why we have to make sure that for
+--
+-- > Bipartition x y
+--
+-- we always have @x >= y@. We ensure by construction that the larger subset
+-- comes first, and so that equality checks are meaningful; see 'bp' and
+-- 'bpUnsafe'.
newtype Bipartition a = Bipartition
{ fromBipartition :: (Set a, Set a)
}
@@ -73,7 +82,7 @@
-- | Create a bipartition from two sets.
--
--- Ensure that the smaller set comes first.
+-- Ensure that the larger set comes first.
--
-- Return 'Left' if one set is empty.
bp :: Ord a => Set a -> Set a -> Either String (Bipartition a)
@@ -84,7 +93,7 @@
-- | Create a bipartition from two sets.
--
--- Ensure that the smaller set comes first.
+-- Ensure that the larger set comes first.
bpUnsafe :: Ord a => Set a -> Set a -> Bipartition a
bpUnsafe xs ys = if xs >= ys then Bipartition (xs, ys) else Bipartition (ys, xs)
@@ -98,7 +107,8 @@
-- > read . show = id
--
-- This identity is met by the derived instance anyways. A more human readable
--- instance would most likely violate the identity.
+-- instance would most likely violate the identity. However, I provide separate
+-- functions to convert bipartitions into human readable strings.
-- | Show a bipartition in a human readable format. Use a provided function to
-- extract information of interest.
@@ -109,10 +119,6 @@
setShow :: Show a => Set a -> String
setShow = intercalate "," . map show . S.toList
--- -- | Map a function over all elements in the 'Bipartition'.
--- bpMap :: Ord b => (a -> b) -> Bipartition a -> Bipartition b
--- bpMap f (Bipartition (x, y)) = bp (S.map f x) (S.map f y)
-
-- | For a bifurcating root, get the bipartition induced by the root node.
--
-- Return 'Left' if
@@ -163,7 +169,7 @@
-- Since the induced bipartitions of the daughter branches of a bifurcating root
-- node are equal, the branches leading to the root have to be combined in this
-- case. See http://evolution.genetics.washington.edu/phylip/doc/treedist.html
--- and how unrooted trees should be handled.
+-- and how unrooted trees are handled.
--
-- Further, branches connected to degree two nodes also induce the same
-- bipartitions and have to be combined.
diff --git a/src/ELynx/Tree/Distance.hs b/src/ELynx/Tree/Distance.hs
--- a/src/ELynx/Tree/Distance.hs
+++ b/src/ELynx/Tree/Distance.hs
@@ -54,7 +54,8 @@
-- XXX: Comparing a list of trees recomputes bipartitions.
symmetric :: Ord a => Tree e1 a -> Tree e2 a -> Either String Int
symmetric t1 t2
- | S.fromList (leaves t1) /= S.fromList (leaves t2) = Left "symmetric: Trees contain different leaves."
+ | S.fromList (leaves t1) /= S.fromList (leaves t2) =
+ Left "symmetric: Trees contain different leaves."
| otherwise = do
bps1 <- bipartitions t1
bps2 <- bipartitions t2
@@ -101,7 +102,7 @@
incompatibleSplits :: (Show a, Ord a) => Tree e1 a -> Tree e2 a -> Either String Int
incompatibleSplits t1 t2
| S.fromList (leaves t1) /= S.fromList (leaves t2) =
- Left "incompatibleSplits: Trees do not have equal leaf sets."
+ Left "incompatibleSplits: Trees contain different leaves."
| otherwise = do
-- Bipartitions.
bs1 <- bipartitions t1
@@ -126,10 +127,11 @@
-- XXX: Comparing a list of trees recomputes bipartitions.
branchScore :: (HasLength e1, HasLength e2, Ord a) => Tree e1 a -> Tree e2 a -> Either String Double
branchScore t1 t2
- | S.fromList (leaves t1) /= S.fromList (leaves t2) = Left "branchScoreWith: Trees do not have equal leaf sets."
+ | S.fromList (leaves t1) /= S.fromList (leaves t2) =
+ Left "branchScoreWith: Trees contain different leaves."
| otherwise = do
- bpToBr1 <- bipartitionToBranch $ first (Sum . getLen) t1
- bpToBr2 <- bipartitionToBranch $ first (Sum . getLen) t2
+ bpToBr1 <- bipartitionToBranch $ first (Sum . getLength) t1
+ bpToBr2 <- bipartitionToBranch $ first (Sum . getLength) t2
let dBs = M.unionWith (-) bpToBr1 bpToBr2
dsSquared = foldl' (\acc e -> acc + e * e) 0 dBs
return $ sqrt $ fromLength $ getSum dsSquared
diff --git a/src/ELynx/Tree/Export/Newick.hs b/src/ELynx/Tree/Export/Newick.hs
--- a/src/ELynx/Tree/Export/Newick.hs
+++ b/src/ELynx/Tree/Export/Newick.hs
@@ -1,3 +1,5 @@
+{-# OPTIONS_GHC -Wno-orphans #-}
+
-- |
-- Module : ELynx.Tree.Export.Newick
-- Description : Export tree objects to Newick format
@@ -29,6 +31,14 @@
import ELynx.Tree.Rooted
import ELynx.Tree.Support
+-- Allow export of trees having branches with lengths only.
+instance HasMaybeSupport Length where
+ getMaybeSupport = const Nothing
+
+-- Allow export of trees having branches with support values only.
+instance HasMaybeLength Support where
+ getMaybeLength = const Nothing
+
buildBrLen :: Length -> BB.Builder
buildBrLen bl = BB.char8 ':' <> BB.doubleDec (fromLength bl)
@@ -36,7 +46,7 @@
buildBrSup bs = BB.char8 '[' <> BB.doubleDec (fromSupport bs) <> BB.char8 ']'
-- | See 'toNewick'.
-toNewickBuilder :: HasName a => Tree Phylo a -> BB.Builder
+toNewickBuilder :: (HasMaybeLength e, HasMaybeSupport e, HasName a) => Tree e a -> BB.Builder
toNewickBuilder t = go t <> BB.char8 ';'
where
go (Node b l []) = lbl b l
@@ -45,25 +55,32 @@
<> mconcat (intersperse (BB.char8 ',') $ map go ts)
<> BB.char8 ')'
<> lbl b l
- mBrSupBuilder x = maybe mempty buildBrSup (brSup x)
- mBrLenBuilder x = maybe mempty buildBrLen (brLen x)
+ mBrSupBuilder x = maybe mempty buildBrSup (getMaybeSupport x)
+ mBrLenBuilder x = maybe mempty buildBrLen (getMaybeLength x)
lbl x y =
BB.lazyByteString (fromName $ getName y)
<> mBrLenBuilder x
-- After reading several discussions, I go for the "more semantical
-- form" with branch support values in square brackets.
<> mBrSupBuilder x
+{-# SPECIALIZE toNewickBuilder :: Tree Length Name -> BB.Builder #-}
+{-# SPECIALIZE toNewickBuilder :: Tree Length Int -> BB.Builder #-}
+{-# SPECIALIZE toNewickBuilder :: Tree Phylo Name -> BB.Builder #-}
+{-# SPECIALIZE toNewickBuilder :: Tree Phylo Int -> BB.Builder #-}
--- | General conversion of a tree into a Newick 'BL.Bytestring'. Use provided
--- functions to extract node labels and branch lengths builder objects. See also
--- Biobase.Newick.Export.
+-- | General conversion of a tree into a Newick 'BL.ByteString'.
--
-- Functions to write key value pairs for nodes are not provided. Those can just
--- be set as node names. For example, the posterior density and the confidence
--- interval of a node can be encoded by setting the node name to:
+-- be set as node labels. For example, the posterior density and the confidence
+-- interval of a node can be encoded by setting the node label to a
+-- 'BL.ByteString':
--
-- @
-- "ACTUALNAME[posterior=-2839.2,age_95%_HPD={4.80804,31.6041}]"
-- @
-toNewick :: HasName a => Tree Phylo a -> BL.ByteString
+toNewick :: (HasMaybeLength e, HasMaybeSupport e, HasName a) => Tree e a -> BL.ByteString
toNewick = BB.toLazyByteString . toNewickBuilder
+{-# SPECIALIZE toNewick :: Tree Length Name -> BL.ByteString #-}
+{-# SPECIALIZE toNewick :: Tree Length Int -> BL.ByteString #-}
+{-# SPECIALIZE toNewick :: Tree Phylo Name -> BL.ByteString #-}
+{-# SPECIALIZE toNewick :: Tree Phylo Int -> BL.ByteString #-}
diff --git a/src/ELynx/Tree/Export/Nexus.hs b/src/ELynx/Tree/Export/Nexus.hs
--- a/src/ELynx/Tree/Export/Nexus.hs
+++ b/src/ELynx/Tree/Export/Nexus.hs
@@ -19,13 +19,17 @@
import qualified Data.ByteString.Lazy.Char8 as BL
import ELynx.Export.Nexus
import ELynx.Tree.Export.Newick
+import ELynx.Tree.Length
import ELynx.Tree.Name
-import ELynx.Tree.Phylogeny
import ELynx.Tree.Rooted
+import ELynx.Tree.Support
-- | Export a list of (NAME, TREE) to a Nexus file.
-toNexusTrees :: HasName a => [(BL.ByteString, Tree Phylo a)] -> BL.ByteString
+toNexusTrees ::
+ (HasMaybeLength e, HasMaybeSupport e, HasName a) =>
+ [(BL.ByteString, Tree e a)] ->
+ BL.ByteString
toNexusTrees ts = toNexus "TREES" (map tree ts)
-tree :: HasName a => (BL.ByteString, Tree Phylo a) -> BL.ByteString
+tree :: (HasMaybeLength e, HasMaybeSupport e, HasName a) => (BL.ByteString, Tree e a) -> BL.ByteString
tree (n, t) = " TREE " <> n <> " = " <> toNewick t
diff --git a/src/ELynx/Tree/Import/Newick.hs b/src/ELynx/Tree/Import/Newick.hs
--- a/src/ELynx/Tree/Import/Newick.hs
+++ b/src/ELynx/Tree/Import/Newick.hs
@@ -45,6 +45,30 @@
import GHC.Generics
import Prelude hiding (takeWhile)
+-- IDEA: Key-value pairs in Newick files.
+--
+-- After some thinking I believe the best way to go is RevBayes-like key-value
+-- pairs after branch lengths and node labels by default and the option to
+-- import IqTree-like trees.
+--
+-- I can not really provide a general parser for key-value pairs, but I can
+-- provide appropriate export functions for reasonably general key-value pairs
+-- such as branch support values. This could look like so:
+--
+-- @
+-- fromNewickG :: Parser (Tree (Maybe Length, BL.ByteString) (Name, BL.ByteString))
+--
+-- fromNewick :: Parser (Tree (Maybe Length) Name)
+--
+-- toNewickG :: Tree BL.ByteString BL.ByteString -> BL.ByteString
+--
+-- toNewick :: (HasMaybeLength e, HasMaybeSupport e, HasName a) => Tree e a -> BL.ByteString
+-- @
+--
+-- In this case, I would also rename RevBayes to KeyVal (or provide a separate
+-- function for IqTree-like trees). I would not ignore the key values but just
+-- provide the whole string to be parsed by the user.
+
-- | Newick tree format.
--
-- - Standard: Branch support values are stored in square brackets after branch
@@ -54,7 +78,7 @@
-- bracket of forests.
--
-- - RevBayes: Key-value pairs are provided in square brackets after node names
--- as well as branch lengths. XXX: Key value pairs are ignored at the moment.
+-- as well as branch lengths. NOTE: Key value pairs are ignored.
data NewickFormat = Standard | IqTree | RevBayes
deriving (Eq, Show, Read, Bounded, Enum, Generic)
@@ -187,10 +211,14 @@
-- RevBayes.
-- RevBayes uses square brackets and key-value pairs to define information
--- about nodes and branches. Parse a single Newick tree. Also succeeds when more
--- trees follow.
+-- about nodes and branches.
--
--- TODO: Key value pairs are ignored at the moment.
+-- Parse a single Newick tree. Also succeeds when more trees follow.
+--
+-- NOTE: Key value pairs are ignored. In my opinion, it is just not a good
+-- option to import key values pairs in this form. Key value pairs can still be
+-- exported by first converting them to a ByteString, and then performing a
+-- normal export.
newickRevBayes :: Parser (Tree Phylo Name)
newickRevBayes =
skipWhile isSpace
@@ -228,7 +256,7 @@
b <- optional branchLengthRevBayes
return $ Node (Phylo b Nothing) n []
--- Drop anything between brackets.
+-- NOTE: Drop anything between brackets.
brackets :: Parser ()
brackets = ( "brackets") $ do
_ <- char '['
diff --git a/src/ELynx/Tree/Length.hs b/src/ELynx/Tree/Length.hs
--- a/src/ELynx/Tree/Length.hs
+++ b/src/ELynx/Tree/Length.hs
@@ -4,7 +4,7 @@
-- |
-- Module : ELynx.Tree.Length
--- Description : Measurable labels
+-- Description : Labels having a length
-- Copyright : (c) Dominik Schrempf 2021
-- License : GPL-3.0-or-later
--
@@ -21,6 +21,7 @@
Length (fromLength),
toLength,
toLengthUnsafe,
+ HasMaybeLength (..),
HasLength (..),
height,
rootHeight,
@@ -37,8 +38,9 @@
import Control.DeepSeq
import Data.Aeson
-import Data.Bifoldable
import Data.Bifunctor
+import Data.Default.Class
+import Data.Foldable
import Data.Semigroup
import ELynx.Tree.Rooted
import ELynx.Tree.Splittable
@@ -72,7 +74,19 @@
-- @
newtype Length = Length {fromLength :: Double}
deriving (Read, Show, Generic, NFData)
- deriving (Enum, Eq, Floating, Fractional, Num, Ord, Real, RealFloat, RealFrac) via Double
+ deriving
+ ( Default,
+ Enum,
+ Eq,
+ Floating,
+ Fractional,
+ Num,
+ Ord,
+ Real,
+ RealFloat,
+ RealFrac
+ )
+ via Double
deriving (Semigroup, Monoid) via Sum Double
instance Splittable Length where
@@ -82,10 +96,13 @@
instance FromJSON Length
+instance HasMaybeLength Length where
+ getMaybeLength = Just
+
instance HasLength Length where
- getLen = id
- setLen = const
- modLen f = f
+ getLength = id
+ setLength = const
+ modifyLength f = f
-- | Return 'Left' if negative.
toLength :: Double -> Either String Length
@@ -97,18 +114,18 @@
toLengthUnsafe :: Double -> Length
toLengthUnsafe = Length
--- | A data type with measurable and modifiable values.
-class HasLength e where
- -- | Get length.
- getLen :: e -> Length
-
- -- | Set length.
- setLen :: Length -> e -> e
+-- | Class of data types that may have a length.
+class HasMaybeLength e where
+ getMaybeLength :: e -> Maybe Length
- -- For computational efficiency.
+instance HasMaybeLength () where
+ getMaybeLength = const Nothing
- -- | Modify length.
- modLen :: (Length -> Length) -> e -> e
+-- | Class of data types with measurable and modifiable length.
+class HasMaybeLength e => HasLength e where
+ getLength :: e -> Length
+ setLength :: Length -> e -> e
+ modifyLength :: (Length -> Length) -> e -> e
-- | The maximum distance between origin and leaves.
--
@@ -125,22 +142,22 @@
--
-- The distances include the branch length of the stem.
distancesOriginLeaves :: HasLength e => Tree e a -> [Length]
-distancesOriginLeaves (Node br _ []) = [getLen br]
-distancesOriginLeaves (Node br _ ts) = map (getLen br +) (concatMap distancesOriginLeaves ts)
+distancesOriginLeaves (Node br _ []) = [getLength br]
+distancesOriginLeaves (Node br _ ts) = map (getLength br +) (concatMap distancesOriginLeaves ts)
-- | Total branch length of a tree.
totalBranchLength :: HasLength e => Tree e a -> Length
-totalBranchLength = bifoldl' (+) const 0 . first getLen
+totalBranchLength = foldl' (+) 0 . fmap getLength . ZipBranchTree
-- | Normalize branch lengths so that the sum is 1.0.
normalizeBranchLengths :: HasLength e => Tree e a -> Tree e a
-normalizeBranchLengths t = first (modLen (/ s)) t
+normalizeBranchLengths t = first (modifyLength (/ s)) t
where
s = totalBranchLength t
-- | Normalize height of tree to 1.0.
normalizeHeight :: HasLength e => Tree e a -> Tree e a
-normalizeHeight t = first (modLen (/ h)) t
+normalizeHeight t = first (modifyLength (/ h)) t
where
h = height t
@@ -163,5 +180,5 @@
where
h = height t
go :: HasLength e => Length -> Tree e a -> Tree e a
- go h' (Node br lb []) = let dh = h - h' - getLen br in Node (modLen (+ dh) br) lb []
- go h' (Node br lb ts) = let h'' = h' + getLen br in Node br lb $ map (go h'') ts
+ go h' (Node br lb []) = let dh = h - h' - getLength br in Node (modifyLength (+ dh) br) lb []
+ go h' (Node br lb ts) = let h'' = h' + getLength br in Node br lb $ map (go h'') ts
diff --git a/src/ELynx/Tree/Name.hs b/src/ELynx/Tree/Name.hs
--- a/src/ELynx/Tree/Name.hs
+++ b/src/ELynx/Tree/Name.hs
@@ -1,4 +1,5 @@
{-# LANGUAGE DerivingVia #-}
+{-# LANGUAGE OverloadedStrings #-}
{-# LANGUAGE TypeSynonymInstances #-}
-- |
@@ -23,7 +24,10 @@
import qualified Data.ByteString.Builder as BB
import qualified Data.ByteString.Char8 as BS
import qualified Data.ByteString.Lazy.Char8 as BL
-import Data.Double.Conversion.ByteString as BC
+import Data.Default.Class
+-- TODO: 2021-09-02: Native conversion is being implemented at the moment.
+-- Remove external library when this is available.
+import qualified Data.Double.Conversion.ByteString as BC
import Data.String
-- | Node name.
@@ -33,6 +37,9 @@
deriving (Show, Read, Eq)
deriving (Ord, Monoid, Semigroup, IsString, NFData) via BL.ByteString
+instance Default Name where
+ def = Name ""
+
-- XXX: This is pretty lame, but I need those instances. At the moment, I just
-- go via 'String', but this is certainly not the best solution.
@@ -57,7 +64,7 @@
getName = Name . BB.toLazyByteString . BB.intDec
instance HasName Double where
- getName = Name . BL.fromStrict . toShortest
+ getName = Name . BL.fromStrict . BC.toShortest
instance HasName Char where
getName = Name . BB.toLazyByteString . BB.char8
diff --git a/src/ELynx/Tree/Parallel.hs b/src/ELynx/Tree/Parallel.hs
--- a/src/ELynx/Tree/Parallel.hs
+++ b/src/ELynx/Tree/Parallel.hs
@@ -19,8 +19,7 @@
module ELynx.Tree.Parallel
( parTree,
parBranchFoldMap,
- parBranchFoldMapWithLayer,
- parNodeFoldMap,
+ parLabelFoldMap,
)
where
@@ -51,35 +50,25 @@
branchFoldMap :: (e -> f) -> (f -> f -> f) -> Tree e a -> f
branchFoldMap f op (Node br _ ts) = foldl' op (f br) $ map (branchFoldMap f op) ts
+-- IDEA: Use and benchmark branch and node specific instances with parFoldMaps.
+--
+-- @
+-- parFoldMap' = blabla
+-- parBranchFoldMap' = parFoldMap' . ZipBranchTree
+-- parNodeFoldMap' = parFoldMap' . ZipNodeTree
+-- @
+
-- | Map and fold over branches. Evaluate the sub trees up to given layer in parallel.
parBranchFoldMap :: NFData f => Int -> (e -> f) -> (f -> f -> f) -> Tree e a -> f
parBranchFoldMap n f op t@(Node br _ ts)
| n >= 1 = foldl' op (f br) (map (parBranchFoldMap (n - 1) f op) ts `using` myParList rdeepseq)
| otherwise = branchFoldMap f op t
-branchFoldMapWithLayer :: Int -> (Int -> e -> f) -> (f -> f -> f) -> Tree e a -> f
-branchFoldMapWithLayer d f op (Node br _ ts) =
- foldl' op (f d br) (map (branchFoldMapWithLayer (d + 1) f op) ts)
-
--- | Map and fold over branches.
---
--- The used function has access to the layer of the node to which the handled
--- branch is attached to.
---
--- Evaluate the sub trees up to given layer in parallel.
-parBranchFoldMapWithLayer :: NFData f => Int -> (Int -> e -> f) -> (f -> f -> f) -> Tree e a -> f
-parBranchFoldMapWithLayer = go 0
- where
- go d n f op t@(Node br _ ts)
- | n >= 1 =
- foldl' op (f d br) (map (go (d + 1) (n - 1) f op) ts `using` myParList rdeepseq)
- | otherwise = branchFoldMapWithLayer d f op t
-
nodeFoldMap :: (a -> b) -> (b -> b -> b) -> Tree e a -> b
nodeFoldMap f op (Node _ lb ts) = foldl' op (f lb) $ map (nodeFoldMap f op) ts
--- | Map and fold over nodes. Evaluate the sub trees up to given layer in parallel.
-parNodeFoldMap :: NFData b => Int -> (a -> b) -> (b -> b -> b) -> Tree e a -> b
-parNodeFoldMap n f op t@(Node _ lb ts)
- | n >= 1 = foldl' op (f lb) (map (parNodeFoldMap (n - 1) f op) ts `using` myParList rdeepseq)
+-- | Map and fold over labels. Evaluate the sub trees up to given layer in parallel.
+parLabelFoldMap :: NFData b => Int -> (a -> b) -> (b -> b -> b) -> Tree e a -> b
+parLabelFoldMap n f op t@(Node _ lb ts)
+ | n >= 1 = foldl' op (f lb) (map (parLabelFoldMap (n - 1) f op) ts `using` myParList rdeepseq)
| otherwise = nodeFoldMap f op t
diff --git a/src/ELynx/Tree/Partition.hs b/src/ELynx/Tree/Partition.hs
--- a/src/ELynx/Tree/Partition.hs
+++ b/src/ELynx/Tree/Partition.hs
@@ -9,9 +9,6 @@
-- Portability : portable
--
-- Creation date: Thu Dec 12 12:58:49 2019.
---
--- A multifurcation induces a 'Partition', similar to branches inducing
--- 'ELynx.Tree.Bipartition's.
module ELynx.Tree.Partition
( -- * Data type
Partition (fromPartition),
@@ -36,14 +33,20 @@
-- | A partition of a tree is a grouping of the leaves of the tree into
-- non-overlapping, non-empty sub sets.
--
--- For example, each branch of a tree partitions the leaves of the tree into two
--- subsets (see 'ELynx.Tree.Bipartition'). In a similar way, each bifurcating
--- internal node (excluding the root node) partitions the leaves into three
--- subsets called a 'Partition'. If the tree is multifurcating, and a specific
--- node has more than two children, the number of subsets induced by this node
--- is larger than three. Partitions are interesting in that we can use them for
--- calculating incompatible splits, see 'ELynx.Tree.Distance'.
+-- For unrooted trees:
--
+-- - For example, each branch of an unrooted tree partitions the leaves of the
+-- tree into two subsets (see 'ELynx.Tree.Bipartition').
+--
+-- For rooted trees:
+--
+-- - In a similar way, each bifurcating internal node (excluding the root node)
+-- partitions the leaves into three subsets called a 'Partition'. If the tree
+-- is multifurcating, and a specific node has more than two children, the
+-- number of subsets induced by this node is larger than three. Partitions are
+-- interesting in that we can use them for calculating incompatible splits,
+-- see 'ELynx.Tree.Distance'.
+--
-- The order of the subsets of a 'Partition' is meaningless. We ensure by
-- construction that the subsets are ordered, and hence, that equality checks
-- are meaningful.
@@ -52,12 +55,10 @@
}
deriving (Eq, Ord, Show, Read)
--- TODO: Check that list is not empty after filtering.
-
-- | Create a partition.
pt :: Ord a => [Set a] -> Either String (Partition a)
pt xs = case filter (not . S.null) xs of
- [] -> Left "mp: Empty list."
+ [] -> Left "pt: Empty list."
xs' -> Right $ ptUnsafe xs'
-- | Create a partition.
diff --git a/src/ELynx/Tree/Phylogeny.hs b/src/ELynx/Tree/Phylogeny.hs
--- a/src/ELynx/Tree/Phylogeny.hs
+++ b/src/ELynx/Tree/Phylogeny.hs
@@ -23,8 +23,8 @@
-- 1. All trees are rooted. Unrooted trees can be treated with a rooted data
-- structure, as it is used here. However, some functions may be meaningless.
--
--- 2. Changing branch labels, node labels, or the topology of the tree are slow
--- operations, especially, when the changes are close to the leaves of the tree.
+-- 2. Changing branch labels, node labels, or the topology of the tree is slow,
+-- especially when the changes are close to the leaves of the tree.
--
-- 3. Internally, the underlying 'Tree' data structure stores the sub-forest as
-- an ordered list. Hence, we have to do some tricks when comparing phylogenies
@@ -35,21 +35,13 @@
-- perform this check, and return 'Left' with a message, if the tree has
-- duplicate leaves.
--
--- Note: 'Tree's are rooted.
+-- NOTE: 'Tree's are rooted.
--
--- Note: 'Tree's encoded in Newick format correspond to rooted trees. By
+-- NOTE: 'Tree's encoded in Newick format correspond to rooted trees. By
-- convention only, a tree parsed from Newick format is usually thought to be
-- unrooted, when the root node is multifurcating and has three or more
-- children. This convention is not used here. Newick trees are just parsed as
-- they are, and a rooted tree is returned.
---
--- A multifurcating root node can be resolved to a bifurcating root node with
--- 'outgroup'.
---
--- The bifurcating root node can be changed with 'outgroup' or 'midpoint'.
---
--- For a given tree with bifurcating root node, a list of all rooted trees is
--- returned by 'roots'.
module ELynx.Tree.Phylogeny
( -- * Functions
equal,
@@ -62,12 +54,16 @@
-- * Branch labels
Phylo (..),
+ toPhyloLabel,
toPhyloTree,
- measurableToPhyloTree,
- supportedToPhyloTree,
- phyloToLengthTree,
- phyloToSupportTree,
- phyloToSupportTreeUnsafe,
+ lengthToPhyloLabel,
+ lengthToPhyloTree,
+ supportToPhyloLabel,
+ supportToPhyloTree,
+ toLengthTree,
+ toSupportTree,
+
+ -- * Explicit branch labels
PhyloExplicit (..),
toExplicitTree,
)
@@ -75,9 +71,8 @@
import Control.DeepSeq
import Data.Aeson
-import Data.Bifoldable
import Data.Bifunctor
-import Data.Bitraversable
+import Data.Default.Class
import Data.List hiding (intersect)
import Data.Maybe
import Data.Monoid
@@ -107,11 +102,11 @@
(brL == brR)
&& (lbL == lbR)
&& (length tsL == length tsR)
- && all (elem' tsR) tsL
+ && all (`elem'` tsR) tsL
where
- elem' ts t = isJust $ find (equal' t) ts
+ elem' t ts = isJust $ find (equal' t) ts
--- | Compute the intersection of trees.
+-- | Intersection of trees.
--
-- The intersections are the largest subtrees sharing the same leaf set.
--
@@ -137,7 +132,7 @@
-- Predicate.
predicate lvsToDr l = l `S.member` lvsToDr
--- | Check if a tree is bifurcating.
+-- | Check if tree is bifurcating.
--
-- A Bifurcating tree only contains degree one (leaves) and degree three nodes
-- (internal bifurcating nodes).
@@ -146,79 +141,88 @@
bifurcating (Node _ _ [x, y]) = bifurcating x && bifurcating y
bifurcating _ = False
--- I believe that manual treatment with 'outgroup' is preferable.
-
--- -- | Remove multifurcations.
--- --
--- -- A caterpillar like bifurcating structure is used to resolve all
--- -- multifurcations on a tree.
--- --
--- -- Multifurcating nodes are copied and branches are 'split'.
--- resolve :: Splittable e => Tree e a -> Tree e a
--- resolve t@(Node _ _ []) = t
--- resolve (Node br lb [x]) = Node br lb [resolve x]
--- resolve (Node br lb [x, y]) = Node br lb $ map resolve [x, y]
--- resolve (Node br lb (Node brL lbL xsL : xs)) = Node br lb [Node brL' lbL (map resolve xsL), Node brL' lb (map resolve xs)]
--- where brL' = split brL
-
--- | Root the tree using an outgroup.
---
--- If the current root node is multifurcating, a bifurcating root node with the
--- empty label is introduced by 'split'ting the leftmost branch. The 'Monoid'
--- instance of the node label and the 'Splittable' instance of the branch length
--- are used.
---
--- NOTE: In this case, the degree of the former root node is decreased by one!
+-- | Root tree using an outgroup.
--
--- Given that the root note is bifurcating, the root node is moved to the
--- required position specified by the outgroup.
+-- If the root note is bifurcating, the root node is moved to the position
+-- specified by the outgroup.
--
--- Branches are connected according to the provided 'Semigroup' instance.
+-- If the root node is multifurcating, a new root node is introduced using the
+-- 'Default' instance of the node labels. Thereby, the degree of the original
+-- root node is reduced by one.
--
--- Upon insertion of the root node at the required position, the affected branch
--- is 'split' according to the provided 'Splittable' instance.
+-- Branches are connected and split according to the provided 'Semigroup' and
+-- 'Splittable' instances.
--
-- Return 'Left' if
--
--- - the root node is not multifurcating;
+-- - the root node is a leaf;
--
+-- - the root node has degree two;
+--
-- - the tree has duplicate leaves;
--
--- - the provided outgroup is not found on the tree or is polyphyletic.
-outgroup :: (Semigroup e, Splittable e, Monoid a, Ord a) => Set a -> Tree e a -> Either String (Tree e a)
+-- - the provided outgroup is polyphyletic or not found on the tree.
+outgroup ::
+ (Semigroup e, Splittable e, Default a, Ord a) =>
+ Set a ->
+ Tree e a ->
+ Either String (Tree e a)
outgroup _ (Node _ _ []) = Left "outgroup: Root node is a leaf."
outgroup _ (Node _ _ [_]) = Left "outgroup: Root node has degree two."
-outgroup o t@(Node _ _ [_, _]) = do
+outgroup o t = do
bip <- bp o (S.fromList (leaves t) S.\\ o)
rootAt bip t
-outgroup o (Node b l ts) = outgroup o t'
+
+-- Root the tree at the branch defined by the given bipartition. The original
+-- root node is moved to the new position.
+rootAt ::
+ (Semigroup e, Splittable e, Eq a, Default a, Ord a) =>
+ Bipartition a ->
+ Tree e a ->
+ Either String (Tree e a)
+rootAt b t
+ -- Do not use 'duplicateLeaves' here, because we also need to compare the leaf
+ -- set with the bipartition.
+ | length lvLst /= S.size lvSet = Left "rootAt: Tree has duplicate leaves."
+ | toSet b /= lvSet = Left "rootAt: Bipartition does not match leaves of tree."
+ | otherwise = do
+ ts <- roots t
+ case find (\x -> bipartition x == Right b) ts of
+ Nothing -> Left "rootAt': Bipartition not found on tree."
+ Just t' -> Right t'
where
- (Node brO lbO tsO) = head ts
- -- Introduce a bifurcating root node.
- t' = Node b mempty [Node (split brO) lbO tsO, Node (split brO) l (tail ts)]
+ lvLst = leaves t
+ lvSet = S.fromList $ leaves t
--- The 'midpoint' algorithm is pretty stupid because it calculates all rooted
--- trees and then finds the one minimizing the difference between the heights of
--- the left and right sub tree. Actually, one just needs to move left or right,
--- with the aim to minimize the height difference between the left and right sub
--- tree.
+-- NOTE: The 'midpoint' algorithm has not been optimized. All rooted trees are
+-- calculated and then the one minimizing the difference between the heights of
+-- the left and right sub tree is chosen. Better: Move left or right minimizing
+-- the height difference between the left and right sub tree.
--- | Root tree at the midpoint.
+-- | Root tree at midpoint.
--
+-- Branches are connected and split according to the provided 'Semigroup' and
+-- 'Splittable' instances.
+--
-- Return 'Left' if
--
--- - the root node is not bifurcating.
-midpoint :: (Semigroup e, Splittable e, HasLength e) => Tree e a -> Either String (Tree e a)
+-- - the root node is a leaf;
+--
+-- - the root node has degree two.
+midpoint ::
+ (Semigroup e, Splittable e, HasLength e, Default a) =>
+ Tree e a ->
+ Either String (Tree e a)
midpoint (Node _ _ []) = Left "midpoint: Root node is a leaf."
midpoint (Node _ _ [_]) = Left "midpoint: Root node has degree two."
-midpoint t@(Node _ _ [_, _]) = roots t >>= getMidpoint
-midpoint _ = Left "midpoint: Root node is multifurcating."
+midpoint t = roots t >>= getMidpoint
-- Find the index of the smallest element.
findMinIndex :: Ord a => [a] -> Either String Int
findMinIndex (x : xs) = go (0, x) 1 xs
where
go (i, _) _ [] = Right i
+ -- Indices with respect to original list: i is index of z, j is index of y.
go (i, z) j (y : ys) = if z < y then go (i, z) (j + 1) ys else go (j, y) (j + 1) ys
findMinIndex [] = Left "findMinIndex: Empty list."
@@ -232,60 +236,68 @@
Node
br
lb
- [ applyStem (modLen (subtract dh)) l,
- applyStem (modLen (+ dh)) r
+ [ modifyStem (modifyLength (subtract' dh)) l,
+ modifyStem (modifyLength (+ dh)) r
]
- -- Explicitly use 'error' here, because roots is supposed to return trees with
- -- bifurcating root nodes.
- Right _ -> error "getMidpoint: Root node is not bifurcating; please contact maintainer."
+ Right _ -> error "getMidpoint: Root node is not bifurcating?"
Left e -> Left e
where
dhs = map getDeltaHeight ts
+ -- Find index of minimum. Take this tree and move root to the midpoint of
+ -- the branch.
t = (ts !!) <$> findMinIndex dhs
-
--- find index of minimum; take this tree and move root to the midpoint of the branch
+ -- Subtract, and check that larger equal 0 with a precision close to the
+ -- machine precision of roughly 1e-16.
+ subtract' dx x =
+ let x' = subtract dx x
+ in case compare x' 0 of
+ LT -> if x' < 1e-14 then error "getMidpoint: Length less than zero." else 0
+ _ -> x'
-- Get delta height of left and right sub tree.
getDeltaHeight :: HasLength e => Tree e a -> Length
getDeltaHeight (Node _ _ [l, r]) = abs $ height l - height r
--- Explicitly use 'error' here, because roots is supposed to return trees with
--- bifurcating root nodes.
-getDeltaHeight _ = error "getDeltaHeight: Root node is not bifurcating; please contact maintainer."
+getDeltaHeight _ = error "getDeltaHeight: Root node is not bifurcating?"
--- | For a rooted tree with a bifurcating root node, get all possible rooted
--- trees.
+-- | Get all rooted trees with bifurcating root nodes.
--
--- The root node (label and branch) is moved.
+-- If the root node of the original tree is bifurcating, the root node (label
+-- and branch) is moved, and the original tree is part of the result.
--
--- For a tree with @l=2@ leaves, there is one rooted tree. For a bifurcating
--- tree with @l>2@ leaves, there are @(2l-3)@ rooted trees. For a general tree
--- with a bifurcating root node, and a total number of @n>2@ nodes, there are
--- (n-2) rooted trees.
+-- If the root node of the original tree is multifurcating, a new root node is
+-- introduced using the 'Default' instance of the node labels. Thereby, the
+-- degree of the original root node is reduced by one. The original,
+-- multifurcating tree is not part of the result.
--
--- A bifurcating root is required because moving a multifurcating root node to
--- another branch would change the degree of the root node. To resolve a
--- multifurcating root, please use 'outgroup'.
+-- Branches are connected and split according to the provided 'Semigroup' and
+-- 'Splittable' instances.
--
--- Connect branches according to the provided 'Semigroup' instance.
+-- For a tree with @n@ nodes we have:
--
--- Split the affected branch into one out of two equal entities according the
--- provided 'Splittable' instance.
+-- - @n-2@ rooted trees if the root node is bifurcating;
--
--- Return 'Left' if the root node is not 'bifurcating'.
-roots :: (Semigroup e, Splittable e) => Tree e a -> Either String (Forest e a)
+-- - (n-1) rooted trees if the root node is multifurcating.
+roots :: (Semigroup e, Splittable e, Default a) => Tree e a -> Either String (Forest e a)
roots (Node _ _ []) = Left "roots: Root node is a leaf."
roots (Node _ _ [_]) = Left "roots: Root node has degree two."
roots t@(Node b c [tL, tR]) = Right $ t : descend b c tR tL ++ descend b c tL tR
-roots _ = Left "roots: Root node is multifurcating."
+roots (Node b c ts) = roots $ Node b def [tL, tR]
+ where
+ (Node bL lL tsL) = head ts
+ bL' = split bL
+ tL = Node bL' lL tsL
+ tR = Node bL' c $ tail ts
complementaryForests :: Tree e a -> Forest e a -> [Forest e a]
complementaryForests t ts = [t : take i ts ++ drop (i + 1) ts | i <- [0 .. (n -1)]]
where
n = length ts
--- From the bifurcating root, descend into one of the two pits.
+-- Descend into the downward tree.
--
--- descend splitFunction rootBranch rootLabel complementaryTree downwardsTree
+-- @
+-- descend rootBranch rootLabel complementaryTree downwardsTree
+-- @
descend :: (Semigroup e, Splittable e) => e -> a -> Tree e a -> Tree e a -> Forest e a
descend _ _ _ (Node _ _ []) = []
descend brR lbR tC (Node brD lbD tsD) =
@@ -301,60 +313,15 @@
tC' = tC {branch = brC'}
cfs = complementaryForests tC' tsD
--- Root a tree at a specific position.
---
--- Root the tree at the branch defined by the given bipartition. The original
--- root node is moved to the new position.
---
--- The root node must be bifurcating (see 'roots' and 'outgroup').
---
--- Connect branches according to the provided 'Semigroup' instance.
---
--- Upon insertion of the root, split the affected branch according to the
--- provided 'Splittable' instance.
---
--- Return 'Left', if:
---
--- - the root node is not bifurcating;
---
--- - the tree has duplicate leaves;
---
--- - the bipartition does not match the leaves of the tree.
-rootAt ::
- (Semigroup e, Splittable e, Eq a, Ord a) =>
- Bipartition a ->
- Tree e a ->
- Either String (Tree e a)
-rootAt b t
- -- Tree is checked for being bifurcating in 'roots'.
- --
- -- Do not use 'duplicateLeaves' here, because we also need to compare the leaf
- -- set with the bipartition.
- | length lvLst /= S.size lvSet = Left "rootAt: Tree has duplicate leaves."
- | toSet b /= lvSet = Left "rootAt: Bipartition does not match leaves of tree."
- | otherwise = rootAt' b t
- where
- lvLst = leaves t
- lvSet = S.fromList $ leaves t
-
--- Assume the leaves of the tree are unique.
-rootAt' ::
- (Semigroup e, Splittable e, Ord a) =>
- Bipartition a ->
- Tree e a ->
- Either String (Tree e a)
-rootAt' b t = do
- ts <- roots t
- case find (\x -> bipartition x == Right b) ts of
- Nothing -> Left "rootAt': Bipartition not found on tree."
- Just t' -> Right t'
-
-- | Branch label for phylogenetic trees.
--
-- Branches may have a length and a support value.
+--
+-- Especially useful to export trees to Newick format; see
+-- 'ELynx.Tree.Export.Newick.toNewick'.
data Phylo = Phylo
- { brLen :: Maybe Length,
- brSup :: Maybe Support
+ { pBranchLength :: Maybe Length,
+ pBranchSupport :: Maybe Support
}
deriving (Read, Show, Eq, Ord, Generic, NFData)
@@ -364,114 +331,108 @@
(getSum <$> (Sum <$> mBL) <> (Sum <$> mBR))
(getMin <$> (Min <$> mSL) <> (Min <$> mSR))
+instance HasMaybeLength Phylo where
+ getMaybeLength = pBranchLength
+
+instance HasMaybeSupport Phylo where
+ getMaybeSupport = pBranchSupport
+
instance ToJSON Phylo
instance FromJSON Phylo
--- | Set all branch lengths and support values to 'Just' the value.
---
--- Useful to export a tree with branch lengths in Newick format.
-toPhyloTree :: (HasLength e, HasSupport e) => Tree e a -> Tree Phylo a
+-- | Set branch length and support value.
+toPhyloLabel :: (HasMaybeLength e, HasMaybeSupport e) => e -> Phylo
+toPhyloLabel x = Phylo (getMaybeLength x) (getMaybeSupport x)
+
+-- | See 'toPhyloLabel'.
+toPhyloTree :: (HasMaybeLength e, HasMaybeSupport e) => Tree e a -> Tree Phylo a
toPhyloTree = first toPhyloLabel
-toPhyloLabel :: (HasLength e, HasSupport e) => e -> Phylo
-toPhyloLabel x = Phylo (Just $ getLen x) (Just $ getSup x)
+-- | Set branch length. Do not set support value.
+lengthToPhyloLabel :: HasMaybeLength e => e -> Phylo
+lengthToPhyloLabel x = Phylo (getMaybeLength x) Nothing
--- | Set all branch lengths to 'Just' the values, and all support values to
--- 'Nothing'.
---
--- Useful to export a tree with branch lengths but without branch support values
--- to Newick format.
-measurableToPhyloTree :: HasLength e => Tree e a -> Tree Phylo a
-measurableToPhyloTree = first measurableToPhyloLabel
+-- | See 'lengthToPhyloLabel'.
+lengthToPhyloTree :: HasMaybeLength e => Tree e a -> Tree Phylo a
+lengthToPhyloTree = first lengthToPhyloLabel
-measurableToPhyloLabel :: HasLength e => e -> Phylo
-measurableToPhyloLabel x = Phylo (Just $ getLen x) Nothing
+-- | Set support value. Do not set branch length.
+supportToPhyloLabel :: HasMaybeSupport e => e -> Phylo
+supportToPhyloLabel x = Phylo Nothing (getMaybeSupport x)
--- | Set all branch lengths to 'Nothing', and all support values to 'Just' the
--- values.
---
--- Useful to export a tree with branch support values but without branch lengths
--- to Newick format.
-supportedToPhyloTree :: HasSupport e => Tree e a -> Tree Phylo a
-supportedToPhyloTree = first supportedToPhyloLabel
+-- | See 'supportToPhyloLabel'.
+supportToPhyloTree :: HasMaybeSupport e => Tree e a -> Tree Phylo a
+supportToPhyloTree = first supportToPhyloLabel
-supportedToPhyloLabel :: HasSupport e => e -> Phylo
-supportedToPhyloLabel x = Phylo Nothing (Just $ getSup x)
+fromMaybeWithError :: String -> Maybe a -> Either String a
+fromMaybeWithError s = maybe (Left s) Right
-- | If root branch length is not available, set it to 0.
--
-- Return 'Left' if any other branch length is unavailable.
-phyloToLengthTree :: Tree Phylo a -> Either String (Tree Length a)
-phyloToLengthTree =
- maybe (Left "phyloToLengthTree: Length unavailable for some branches.") Right
- . bitraverse brLen pure
- . cleanStemLength
-
-cleanStemLength :: Tree Phylo a -> Tree Phylo a
-cleanStemLength (Node (Phylo Nothing s) l f) = Node (Phylo (Just 0) s) l f
-cleanStemLength t = t
+toLengthTree :: HasMaybeLength e => Tree e a -> Either String (Tree Length a)
+toLengthTree (Node br lb ts) =
+ case traverse go ts of
+ Nothing -> Left "toLengthTree: Length unavailable for some branches."
+ Just ts' -> Right $ Node br' lb ts'
+ where
+ br' = fromMaybe 0 $ getMaybeLength br
+ go t = getBranchTree <$> traverse getMaybeLength (BranchTree t)
-- | Set branch support values of branches leading to the leaves and of the root
-- branch to maximum support.
--
-- Return 'Left' if any other branch has no available support value.
-phyloToSupportTree :: Tree Phylo a -> Either String (Tree Support a)
-phyloToSupportTree t =
- maybe
- (Left "phyloToSupportTree: Support value unavailable for some branches.")
- Right
- $ bitraverse brSup pure $
- cleanLeafSupport m $
- cleanRootSupport m t
- where
- m = getMaxSupport t
-
--- | Set all unavailable branch support values to maximum support.
-phyloToSupportTreeUnsafe :: Tree Phylo a -> Tree Support a
-phyloToSupportTreeUnsafe t = cleanSupport m t
+toSupportTree :: HasMaybeSupport e => Tree e a -> Either String (Tree Support a)
+toSupportTree t@(Node br lb ts) =
+ fromMaybeWithError "toSupportTree: Support value unavailable for some branches." $
+ getBranchTree <$> sequenceA (BranchTree (Node br' lb $ map go ts))
where
m = getMaxSupport t
+ br' = cleanSupportWith m br
+ go (Node b l []) = Node (cleanSupportWith m b) l []
+ go (Node b l xs) = Node (getMaybeSupport b) l (map go xs)
-- If all branch support values are below 1.0, set the max support to 1.0.
-getMaxSupport :: Tree Phylo a -> Support
-getMaxSupport = fromJust . max (Just 1.0) . bimaximum . bimap brSup (const Nothing)
-
-cleanRootSupport :: Support -> Tree Phylo a -> Tree Phylo a
-cleanRootSupport maxSup (Node (Phylo b Nothing) l xs) = Node (Phylo b (Just maxSup)) l xs
-cleanRootSupport _ t = t
-
-cleanLeafSupport :: Support -> Tree Phylo a -> Tree Phylo a
-cleanLeafSupport s (Node (Phylo b Nothing) l []) = Node (Phylo b (Just s)) l []
-cleanLeafSupport s (Node b l xs) = Node b l $ map (cleanLeafSupport s) xs
+getMaxSupport :: HasMaybeSupport e => Tree e a -> Support
+getMaxSupport = fromJust . max (Just 1.0) . maximum . fmap getMaybeSupport . ZipBranchTree
-cleanSupport :: Support -> Tree Phylo a -> Tree Support a
-cleanSupport maxSup (Node (Phylo _ s) l xs) = Node (fromMaybe maxSup s) l $ map (cleanSupport maxSup) xs
+cleanSupportWith :: HasMaybeSupport e => Support -> e -> Maybe Support
+cleanSupportWith m x = case getMaybeSupport x of
+ Nothing -> Just m
+ Just y -> Just y
-- | Explicit branch label with branch length and branch support value.
data PhyloExplicit = PhyloExplicit
- { sBrLen :: Length,
- sBrSup :: Support
+ { eBranchLength :: Length,
+ eBranchSupport :: Support
}
deriving (Read, Show, Eq, Ord, Generic)
instance Semigroup PhyloExplicit where
PhyloExplicit bL sL <> PhyloExplicit bR sR = PhyloExplicit (bL + bR) (min sL sR)
+instance HasMaybeLength PhyloExplicit where
+ getMaybeLength = Just . eBranchLength
+
instance HasLength PhyloExplicit where
- getLen = sBrLen
- setLen b pl = pl {sBrLen = b}
- modLen f (PhyloExplicit l s) = PhyloExplicit (f l) s
+ getLength = eBranchLength
+ setLength b pl = pl {eBranchLength = b}
+ modifyLength f (PhyloExplicit l s) = PhyloExplicit (f l) s
instance Splittable PhyloExplicit where
- split l = l {sBrLen = b'}
+ split l = l {eBranchLength = b'}
where
- b' = sBrLen l / 2.0
+ b' = eBranchLength l / 2.0
+instance HasMaybeSupport PhyloExplicit where
+ getMaybeSupport = Just . eBranchSupport
+
instance HasSupport PhyloExplicit where
- getSup = sBrSup
- setSup s pl = pl {sBrSup = s}
- modSup f (PhyloExplicit l s) = PhyloExplicit l (f s)
+ getSupport = eBranchSupport
+ setSupport s pl = pl {eBranchSupport = s}
+ modifySupport f (PhyloExplicit l s) = PhyloExplicit l (f s)
instance ToJSON PhyloExplicit
@@ -479,12 +440,15 @@
-- | Conversion to a 'PhyloExplicit' tree.
--
--- See 'phyloToLengthTree' and 'phyloToSupportTree'.
-toExplicitTree :: Tree Phylo a -> Either String (Tree PhyloExplicit a)
+-- See 'toLengthTree' and 'toSupportTree'.
+toExplicitTree ::
+ (HasMaybeLength e, HasMaybeSupport e) =>
+ Tree e a ->
+ Either String (Tree PhyloExplicit a)
toExplicitTree t = do
- lt <- phyloToLengthTree t
- st <- phyloToSupportTree t
+ lt <- toLengthTree t
+ st <- toSupportTree t
case zipTreesWith PhyloExplicit const lt st of
- -- Explicit use of error, since this case should not happen.
- Nothing -> error "toExplicitTree: Can not zip two trees with different topologies; please contact maintainer."
+ -- Explicit use of error, since this case should never happen.
+ Nothing -> error "toExplicitTree: Can not zip two trees with different topologies."
Just zt -> return zt
diff --git a/src/ELynx/Tree/Rooted.hs b/src/ELynx/Tree/Rooted.hs
--- a/src/ELynx/Tree/Rooted.hs
+++ b/src/ELynx/Tree/Rooted.hs
@@ -1,5 +1,7 @@
{-# LANGUAGE DeriveDataTypeable #-}
{-# LANGUAGE DeriveGeneric #-}
+{-# LANGUAGE GeneralizedNewtypeDeriving #-}
+{-# LANGUAGE StandaloneDeriving #-}
-- |
-- Module : ELynx.Tree.Rooted
@@ -13,8 +15,8 @@
--
-- Creation date: Thu Jan 17 09:57:29 2019.
--
--- Rooted 'Tree's differes from a classical rose 'Data.Tree.Tree' in that it has
--- labeled branches.
+-- Rooted 'Tree's differs from a classical rose 'Data.Tree.Tree's in that they
+-- have labeled branches.
--
-- For rooted topologies, please see 'ELynx.Topology.Rooted'.
--
@@ -36,8 +38,9 @@
--
-- This means, that the word 'Node' is reserved for the constructor of a tree,
-- and that a 'Node' has an attached 'branch', a 'label', and a sub-'forest'.
--- The value constructor /Node/ and the record function /label/ are not to be
--- confused. The elements of the sub-forest are often called /children/.
+-- The terms /Node/ and /label/ referring to the value constructor 'Node' and
+-- the record function 'label', respectively, are not to be confused. The
+-- elements of the sub-forest are often called /children/.
--
-- In mathematical terms: A 'Tree' is a directed acyclic graph without loops,
-- with vertex labels, and with edge labels.
@@ -57,9 +60,10 @@
-- functions such as 'branches', or 'labels' which is the same as 'toList'.
-- Please let me know, if post-order algorithms are required.
module ELynx.Tree.Rooted
- ( -- * Data type
+ ( -- * Tree with branch labels
Tree (..),
Forest,
+ fromRoseTree,
toTreeBranchLabels,
toTreeNodeLabels,
@@ -67,14 +71,13 @@
leaves,
duplicateLeaves,
setStem,
- applyStem,
+ modifyStem,
branches,
setBranches,
setLabel,
- applyLabel,
+ modifyLabel,
labels,
setLabels,
- applyRoot,
identify,
-- * Structure
@@ -85,6 +88,12 @@
dropLeavesWith,
zipTreesWith,
zipTrees,
+ flipLabels,
+
+ -- * Newtypes with specific instances
+ ZipTree (..),
+ BranchTree (..),
+ ZipBranchTree (..),
)
where
@@ -92,13 +101,13 @@
import Control.Comonad
import Control.DeepSeq
import Control.Monad
-import Control.Monad.Fix
import Data.Aeson
import Data.Bifoldable
import Data.Bifunctor
import Data.Bitraversable
import Data.Data
import Data.Foldable
+import Data.Functor
import Data.List
import Data.Maybe
import qualified Data.Set as S
@@ -110,8 +119,6 @@
-- Unary instances such as 'Functor' act on node labels, and not on branch
-- labels. Binary instances such as 'Bifunctor' act on both labels (`first` acts
-- on branches, `second` on node labels).
---
--- Lifted instances are not provided.
data Tree e a = Node
{ branch :: e,
label :: a,
@@ -119,13 +126,13 @@
}
deriving (Eq, Read, Show, Data, Generic)
--- | A shorthand.
+-- | Shorthand.
type Forest e a = [Tree e a]
-- | Map over node labels.
instance Functor (Tree e) where
fmap f ~(Node br lb ts) = Node br (f lb) $ map (fmap f) ts
- x <$ ~(Node br _ ts) = Node br x (map (x <$) ts)
+ lb <$ ~(Node br _ ts) = Node br lb (map (lb <$) ts)
-- | The function 'first' acts on branch labels, 'second' on node labels.
instance Bifunctor Tree where
@@ -150,50 +157,25 @@
instance Bitraversable Tree where
bitraverse f g ~(Node br lb ts) = Node <$> f br <*> g lb <*> traverse (bitraverse f g) ts
--- The following code provides a zip-like applicative instance. However,
--- the zip-like instance makes the Monad instance meaningless. So, either we
--- provide only 'Applicative' in zip-like form, or we use the classic instance
--- for 'Applicative' and 'Monad'.
-
--- -- | Note: The 'Applicative' instance of 'Tree' is similar to the one of
--- -- 'Control.Applicative.ZipList', and differs from the instance of
--- -- 'Data.Tree.Tree'!
--- --
--- -- >>> let t = Node "" 0 [Node "" 1 [], Node "" 2 []] :: Tree String Int
--- -- >>> let f = Node "+3" (+3) [Node "*5" (*5) [], Node "+10" (+10) []] :: Tree String (Int -> Int)
--- -- >>> f <*> t
--- -- Node {branch = "+3", label = 3, forest = [Node {branch = "*5", label = 5, forest = []},Node {branch = "+10", label = 12, forest = []}]}
--- --
--- -- Note: The 'Monoid' instance of the branch labels determines how the branches
--- -- are combined. For example, distances can be summed using the
--- -- 'Data.Monoid.Sum' monoid.
--- instance Monoid e => Applicative (Tree e) where
--- pure lb = Node mempty lb []
--- ~(Node brF lbF tsF) <*> ~(Node brX lbX tsX) =
--- Node (brF <> brX) (lbF lbX) (zipWith (<*>) tsF tsX)
--- liftA2 f ~(Node brX lbX tsX) ~(Node brY lbY tsY) =
--- Node (brX <> brY) (f lbX lbY) (zipWith (liftA2 f) tsX tsY)
--- ~(Node brX _ tsX) *> ~(Node brY lbY tsY) =
--- Node (brX <> brY) lbY (zipWith (*>) tsX tsY)
--- ~(Node brX lbX tsX) <* ~(Node brY _ tsY) =
--- Node (brX <> brY) lbX (zipWith (<*) tsX tsY)
-
-- | The 'Semigroup' instance of the branch labels determines how the
-- branches are combined. For example, distances can be summed using
-- 'Data.Semigroup.Sum'.
--
-- The 'Monoid' instance of the branch labels determines the default branch
-- label when using 'pure'.
+--
+-- This instance is similar to the one provided by 'Data.Tree.Tree'. For an
+-- alternative, see 'ZipTree'.
instance (Semigroup e, Monoid e) => Applicative (Tree e) where
pure lb = Node mempty lb []
~(Node brF lbF tsF) <*> ~tx@(Node brX lbX tsX) =
- Node (brF <> brX) (lbF lbX) (map (lbF <$>) tsX ++ map (<*> tx) tsF)
+ Node (brF <> brX) (lbF lbX) (map (bimap (brF <>) lbF) tsX ++ map (<*> tx) tsF)
liftA2 f ~(Node brX lbX tsX) ~ty@(Node brY lbY tsY) =
- Node (brX <> brY) (f lbX lbY) (map (f lbX <$>) tsY ++ map (\tx -> liftA2 f tx ty) tsX)
+ Node (brX <> brY) (f lbX lbY) (map (bimap (brX <>) (f lbX)) tsY ++ map (\tx -> liftA2 f tx ty) tsX)
~(Node brX _ tsX) *> ~ty@(Node brY lbY tsY) =
- Node (brX <> brY) lbY (tsY ++ map (*> ty) tsX)
+ Node (brX <> brY) lbY (map (first (brX <>)) tsY ++ (tsX <&> (*> ty)))
~(Node brX lbX tsX) <* ~ty@(Node brY _ tsY) =
- Node (brX <> brY) lbX (map (lbX <$) tsY ++ map (<* ty) tsX)
+ Node (brX <> brY) lbX (map (bimap (brX <>) (const lbX)) tsY ++ (tsX <&> (<* ty)))
-- | The 'Semigroup' instance of the branch labels determines how the branches
-- are combined. For example, distances can be summed using
@@ -203,10 +185,11 @@
-- label when using 'return'.
instance (Semigroup e, Monoid e) => Monad (Tree e) where
~(Node br lb ts) >>= f = case f lb of
- Node br' lb' ts' -> Node (br <> br') lb' (ts' ++ map (>>= f) ts)
+ Node br' lb' ts' -> Node (br <> br') lb' (map (first (br <>)) ts' ++ map (>>= f) ts)
--- -- Cannot provide MonadZip instance because branch labels cannot be
--- -- recovered from combined label.
+-- -- NOTE: We cannot provide a MonadZip instance because branch labels cannot
+-- -- be recovered from the combined label.
+--
-- instance Monoid e => MonadZip (Tree e) where
-- mzipWith f (Node brL lbL tsL) (Node brR lbR tsR) =
-- Node (brL <> brR) (f lbL lbR) (mzipWith (mzipWith f) tsL tsR)
@@ -215,24 +198,27 @@
-- where
-- (tsL, tsR) = munzip (map munzip ts)
-instance Monoid e => MonadFix (Tree e) where
- mfix = mfixTree
+-- -- NOTE: I don't really know much about 'MonadFix', and so do not provide the
+-- -- instance.
+--
+-- instance Monoid e => MonadFix (Tree e) where
+-- mfix = mfixTree
-mfixTree :: (a -> Tree e a) -> Tree e a
-mfixTree f
- | Node br lb ts <- fix (f . label) =
- Node
- br
- lb
- ( zipWith
- (\i _ -> mfixTree ((!! i) . forest . f))
- [0 ..]
- ts
- )
+-- mfixTree :: (a -> Tree e a) -> Tree e a
+-- mfixTree f
+-- | Node br lb ts <- fix (f . label) =
+-- Node
+-- br
+-- lb
+-- ( zipWith
+-- (\i _ -> mfixTree ((!! i) . forest . f))
+-- [0 ..]
+-- ts
+-- )
instance Comonad (Tree e) where
duplicate t@(Node br _ ts) = Node br t (map duplicate ts)
- extract (Node _ lb _) = lb
+ extract = label
{-# INLINE extract #-}
instance (NFData e, NFData a) => NFData (Tree e a) where
@@ -242,6 +228,10 @@
instance (FromJSON e, FromJSON a) => FromJSON (Tree e a)
+-- | Conversion from 'T.Tree'.
+fromRoseTree :: T.Tree a -> Tree () a
+fromRoseTree (T.Node l ts) = Node () l $ map fromRoseTree ts
+
-- | Conversion to 'T.Tree' using branch labels.
toTreeBranchLabels :: Tree e a -> T.Tree e
toTreeBranchLabels (Node br _ ts) = T.Node br (map toTreeBranchLabels ts)
@@ -250,11 +240,12 @@
toTreeNodeLabels :: Tree e a -> T.Tree a
toTreeNodeLabels (Node _ lb ts) = T.Node lb (map toTreeNodeLabels ts)
--- TODO: Maybe use foldr similar to 'labels'.
--- | Get leaves.
+-- | List of leaves.
leaves :: Tree e a -> [a]
-leaves (Node _ lb []) = [lb]
-leaves (Node _ _ ts) = concatMap leaves ts
+leaves t = squish t []
+ where
+ squish (Node _ lb []) xs = lb : xs
+ squish (Node _ _ ts) xs = foldr squish xs ts
duplicates :: Ord a => [a] -> Bool
duplicates = go S.empty
@@ -270,9 +261,9 @@
setStem :: e -> Tree e a -> Tree e a
setStem br (Node _ lb ts) = Node br lb ts
--- | Change the root branch of a tree.
-applyStem :: (e -> e) -> Tree e a -> Tree e a
-applyStem f t = t {branch = f $ branch t}
+-- | Modify the stem of a tree.
+modifyStem :: (e -> e) -> Tree e a -> Tree e a
+modifyStem f t = t {branch = f $ branch t}
-- | Get branch labels in pre-order.
branches :: Tree e a -> [e]
@@ -290,13 +281,13 @@
setBranch (y : ys) _ = (ys, Just y)
noChange ys z = (ys, Just z)
--- | Set the label to a given value.
+-- | Set label.
setLabel :: a -> Tree e a -> Tree e a
setLabel lb (Node br _ ts) = Node br lb ts
--- | Change the root branch of a tree.
-applyLabel :: (a -> a) -> Tree e a -> Tree e a
-applyLabel f t = t {label = f $ label t}
+-- | Modify the root label of a tree.
+modifyLabel :: (a -> a) -> Tree e a -> Tree e a
+modifyLabel f t = t {label = f $ label t}
-- | Return node labels in pre-order.
labels :: Tree e a -> [a]
@@ -313,11 +304,7 @@
setLabelM [] _ = ([], Nothing)
setLabelM (y : ys) _ = (ys, Just y)
--- | Change the root label of a tree.
-applyRoot :: (a -> a) -> Tree e a -> Tree e a
-applyRoot f t = t {label = f $ label t}
-
--- | Label the nodes with unique integers starting at the root with 0.
+-- | Label the nodes in pre-order with unique indices starting at 0.
identify :: Traversable t => t a -> t Int
identify = snd . mapAccumL (\i _ -> (i + 1, i)) (0 :: Int)
@@ -336,14 +323,16 @@
-- is 1.
depth :: Tree e a -> Int
depth = maximum . go 1
- where go n (Node _ _ []) = [n]
- go n (Node _ _ xs) = concatMap (go (n+1)) xs
+ where
+ go n (Node _ _ []) = [n]
+ go n (Node _ _ xs) = concatMap (go (n + 1)) xs
-- | Prune degree two nodes.
--
--- The information stored in a pruned node is lost. The branches are combined
--- according to their 'Semigroup' instance of the form @\daughterBranch
--- parentBranch -> combinedBranch@.
+-- The label of a pruned node is lost. The branches are combined according to
+-- their 'Semigroup' instance of the form
+--
+-- @\daughterBranch parentBranch -> combinedBranch@.
prune :: Semigroup e => Tree e a -> Tree e a
prune t@(Node _ _ []) = t
prune (Node paBr _ [Node daBr daLb daTs]) = Node (daBr <> paBr) daLb daTs
@@ -353,9 +342,13 @@
--
-- Degree two nodes may arise.
--
--- Also drop parent nodes of which all daughter nodes are dropped.
+-- Also drop nodes of which all daughter nodes are dropped.
--
--- Return 'Nothing' if the root node satisfies the predicate.
+-- Return 'Nothing' if
+--
+-- - The root node satisfies the predicate.
+--
+-- - All daughter nodes of the root are dropped.
dropNodesWith :: (a -> Bool) -> Tree e a -> Maybe (Tree e a)
dropNodesWith p (Node br lb ts)
| p lb = Nothing
@@ -370,7 +363,7 @@
--
-- Degree two nodes may arise.
--
--- Also drop parent nodes of which all leaves are dropped.
+-- Also drop nodes of which all daughter nodes are dropped.
--
-- Return 'Nothing' if all leaves satisfy the predicate.
dropLeavesWith :: (a -> Bool) -> Tree e a -> Maybe (Tree e a)
@@ -386,6 +379,10 @@
-- | Zip two trees with the same topology.
--
+-- This function differs from the 'Applicative' instance of 'ZipTree' in that it
+-- fails when the topologies don't match. Further, it allows specification of a
+-- zipping function for the branches.
+--
-- Return 'Nothing' if the topologies are different.
zipTreesWith ::
(e1 -> e2 -> e) ->
@@ -399,8 +396,195 @@
zipWithM (zipTreesWith f g) tsL tsR >>= Just . Node (f brL brR) (g lbL lbR)
else Nothing
--- | Zip two trees with the same topology.
---
--- Return 'Nothing' if the topologies are different.
+-- | See 'zipTreesWith'.
zipTrees :: Tree e1 a1 -> Tree e2 a2 -> Maybe (Tree (e1, e2) (a1, a2))
zipTrees = zipTreesWith (,) (,)
+
+-- | Flip the branch and node lables.
+flipLabels :: Tree e a -> Tree a e
+flipLabels (Node x y zs) = Node y x $ map flipLabels zs
+
+-- | This newtype provides instances acting on the branch labels, and not on the
+-- node labels as it is the case in 'Tree'.
+newtype BranchTree a e = BranchTree {getBranchTree :: Tree e a}
+ deriving (Eq, Read, Show, Data, Generic)
+
+-- | Map over branch labels.
+instance Functor (BranchTree a) where
+ fmap f ~(BranchTree (Node br lb ts)) =
+ BranchTree $ Node (f br) lb $ map (getBranchTree . fmap f . BranchTree) ts
+ br <$ ~(BranchTree (Node _ lb ts)) =
+ BranchTree $ Node br lb (map (getBranchTree . (br <$) . BranchTree) ts)
+
+-- | Combine branch labels in pre-order.
+instance Foldable (BranchTree a) where
+ foldMap f ~(BranchTree (Node br _ ts)) =
+ f br <> foldMap (foldMap f . BranchTree) ts
+ null _ = False
+ {-# INLINE null #-}
+ toList = branches . getBranchTree
+ {-# INLINE toList #-}
+
+instance Traversable (BranchTree a) where
+ traverse g ~(BranchTree (Node br lb ts)) =
+ assemble lb <$> fbr' <*> fts'
+ where
+ assemble lb' br' ts' = BranchTree $ Node br' lb' ts'
+ fbr' = g br
+ fts' = map getBranchTree <$> traverse (traverse g . BranchTree) ts
+
+instance Comonad (BranchTree a) where
+ duplicate (BranchTree t@(Node _ lb ts)) =
+ BranchTree $
+ Node (BranchTree t) lb $
+ map (getBranchTree . duplicate . BranchTree) ts
+ extract = branch . getBranchTree
+
+instance Monoid a => Applicative (BranchTree a) where
+ -- Infinite layers with infinite subtrees.
+ pure br = BranchTree $ Node br mempty []
+ (BranchTree ~(Node brF lbF tsF)) <*> tx@(BranchTree ~(Node brX lbX tsX)) =
+ BranchTree $
+ Node
+ (brF brX)
+ (lbF <> lbX)
+ ( map (bimap brF (lbF <>)) tsX
+ ++ map (getBranchTree . (<*> tx) . BranchTree) tsF
+ )
+ liftA2 f (BranchTree ~(Node brX lbX tsX)) ty@(BranchTree ~(Node brY lbY tsY)) =
+ BranchTree $
+ Node
+ (f brX brY)
+ (lbX <> lbY)
+ ( map (bimap (f brX) (lbX <>)) tsY
+ ++ map (\tx -> getBranchTree $ liftA2 f (BranchTree tx) ty) tsX
+ )
+ (BranchTree ~(Node _ lbX tsX)) *> ty@(BranchTree ~(Node brY lbY tsY)) =
+ BranchTree $
+ Node
+ brY
+ (lbX <> lbY)
+ ( getBranchTree
+ <$> ( map (BranchTree . second (lbX <>)) tsY
+ ++ map ((*> ty) . BranchTree) tsX
+ )
+ )
+ (BranchTree ~(Node brX lbX tsX)) <* ty@(BranchTree ~(Node _ lbY tsY)) =
+ BranchTree $
+ Node
+ brX
+ (lbX <> lbY)
+ ( map (bimap (const brX) (lbX <>)) tsY
+ ++ map (getBranchTree . (<* ty) . BranchTree) tsX
+ )
+
+-- | This newtype provides a zip-like applicative instance, similar to
+-- 'Control.Applicative.ZipList'.
+--
+-- The default applicative instance of 'Tree' is not zip-like, because the
+-- zip-like instance makes the Monad instance meaningless (similar to the
+-- behavior observed with lists).
+newtype ZipTree e a = ZipTree {getZipTree :: Tree e a}
+ deriving (Eq, Read, Show, Data, Generic)
+
+deriving instance Functor (ZipTree e)
+
+deriving instance Foldable (ZipTree e)
+
+instance Traversable (ZipTree e) where
+ traverse f (ZipTree t) = ZipTree <$> traverse f t
+
+instance Comonad (ZipTree e) where
+ duplicate (ZipTree t) = ZipTree $ second ZipTree $ duplicate t
+ extract = label . getZipTree
+
+-- | The 'Monoid' instance of the branch labels determines the default branch
+-- label, and how the branches are combined. For example, distances can be
+-- summed using the 'Data.Monoid.Sum' monoid.
+--
+-- >>> let t = ZipTree $ Node "" 0 [Node "" 1 [], Node "" 2 []] :: ZipTree String Int
+-- >>> let f = ZipTree $ Node "+3" (+3) [Node "*5" (*5) [], Node "+10" (+10) []] :: ZipTree String (Int -> Int)
+-- >>> f <*> t
+--
+-- ZipTree {getZipTree = Node {branch = "+3", label = 3, forest = [Node {branch = "*5", label = 5, forest = []},Node {branch = "+10", label = 12, forest = []}]}}
+instance Monoid e => Applicative (ZipTree e) where
+ -- Infinite layers with infinite subtrees.
+ pure lb = ZipTree $ Node mempty lb $ repeat (getZipTree $ pure lb)
+ (ZipTree ~(Node brF lbF tsF)) <*> (ZipTree ~(Node brX lbX tsX)) =
+ ZipTree $ Node (brF <> brX) (lbF lbX) (zipWith f tsF tsX)
+ where
+ f x y = getZipTree $ ZipTree x <*> ZipTree y
+ liftA2 f (ZipTree ~(Node brX lbX tsX)) (ZipTree ~(Node brY lbY tsY)) =
+ ZipTree $ Node (brX <> brY) (f lbX lbY) (zipWith g tsX tsY)
+ where
+ g x y = getZipTree $ liftA2 f (ZipTree x) (ZipTree y)
+ (ZipTree ~(Node brX _ tsX)) *> (ZipTree ~(Node brY lbY tsY)) =
+ ZipTree $ Node (brX <> brY) lbY (zipWith f tsX tsY)
+ where
+ f x y = getZipTree $ ZipTree x *> ZipTree y
+ (ZipTree ~(Node brX lbX tsX)) <* (ZipTree ~(Node brY _ tsY)) =
+ ZipTree $ Node (brX <> brY) lbX (zipWith f tsX tsY)
+ where
+ f x y = getZipTree $ ZipTree x <* ZipTree y
+
+-- | Like 'ZipTree' but act on branch labels; see 'BranchTree'.
+newtype ZipBranchTree a e = ZipBranchTree {getZipBranchTree :: Tree e a}
+ deriving (Eq, Read, Show, Data, Generic)
+
+-- | Map over branch labels.
+instance Functor (ZipBranchTree a) where
+ fmap f ~(ZipBranchTree (Node br lb ts)) =
+ ZipBranchTree $ Node (f br) lb $ map g ts
+ where
+ g = getZipBranchTree . fmap f . ZipBranchTree
+ br <$ ~(ZipBranchTree (Node _ lb ts)) =
+ ZipBranchTree $ Node br lb (map f ts)
+ where
+ f = getZipBranchTree . (br <$) . ZipBranchTree
+
+-- | Combine branch labels in pre-order.
+instance Foldable (ZipBranchTree a) where
+ foldMap f ~(ZipBranchTree (Node br _ ts)) =
+ f br <> foldMap g ts
+ where
+ g = foldMap f . ZipBranchTree
+ null _ = False
+ {-# INLINE null #-}
+ toList = branches . getZipBranchTree
+ {-# INLINE toList #-}
+
+instance Traversable (ZipBranchTree a) where
+ traverse g ~(ZipBranchTree (Node br lb ts)) =
+ assemble lb <$> fbr' <*> fts'
+ where
+ assemble lb' br' ts' = ZipBranchTree $ Node br' lb' ts'
+ fbr' = g br
+ fts' = map getZipBranchTree <$> traverse (traverse g . ZipBranchTree) ts
+
+instance Comonad (ZipBranchTree a) where
+ duplicate (ZipBranchTree t@(Node _ lb ts)) =
+ ZipBranchTree $
+ Node (ZipBranchTree t) lb $
+ map (getZipBranchTree . duplicate . ZipBranchTree) ts
+ extract = branch . getZipBranchTree
+
+-- | See the 'Applicative' instance of 'ZipTree'.
+instance Monoid a => Applicative (ZipBranchTree a) where
+ -- Infinite layers with infinite subtrees.
+ pure br = ZipBranchTree $ Node br mempty $ repeat (getZipBranchTree $ pure br)
+ (ZipBranchTree ~(Node brF lbF tsF)) <*> (ZipBranchTree ~(Node brX lbX tsX)) =
+ ZipBranchTree $ Node (brF brX) (lbF <> lbX) (zipWith f tsF tsX)
+ where
+ f x y = getZipBranchTree $ ZipBranchTree x <*> ZipBranchTree y
+ liftA2 f (ZipBranchTree ~(Node brX lbX tsX)) (ZipBranchTree ~(Node brY lbY tsY)) =
+ ZipBranchTree $ Node (f brX brY) (lbX <> lbY) (zipWith g tsX tsY)
+ where
+ g x y = getZipBranchTree $ liftA2 f (ZipBranchTree x) (ZipBranchTree y)
+ (ZipBranchTree ~(Node _ lbX tsX)) *> (ZipBranchTree ~(Node brY lbY tsY)) =
+ ZipBranchTree $ Node brY (lbX <> lbY) (zipWith f tsX tsY)
+ where
+ f x y = getZipBranchTree $ ZipBranchTree x *> ZipBranchTree y
+ (ZipBranchTree ~(Node brX lbX tsX)) <* (ZipBranchTree ~(Node _ lbY tsY)) =
+ ZipBranchTree $ Node brX (lbX <> lbY) (zipWith f tsX tsY)
+ where
+ f x y = getZipBranchTree $ ZipBranchTree x <* ZipBranchTree y
diff --git a/src/ELynx/Tree/Simulate/Coalescent.hs b/src/ELynx/Tree/Simulate/Coalescent.hs
--- a/src/ELynx/Tree/Simulate/Coalescent.hs
+++ b/src/ELynx/Tree/Simulate/Coalescent.hs
@@ -49,7 +49,7 @@
i <- uniformR (1, n - 1) g
-- The time of the coalescent event.
t <- toLengthUnsafe <$> genContVar (coalescentDistributionCont n) g
- let trs' = map (applyStem (+ t)) trs -- Move time 't' up on the tree.
+ let trs' = map (modifyStem (+ t)) trs -- Move time 't' up on the tree.
tl = trs' !! (i - 1)
tr = trs' !! i
-- Join the two chosen trees.
diff --git a/src/ELynx/Tree/Simulate/PointProcess.hs b/src/ELynx/Tree/Simulate/PointProcess.hs
--- a/src/ELynx/Tree/Simulate/PointProcess.hs
+++ b/src/ELynx/Tree/Simulate/PointProcess.hs
@@ -52,8 +52,6 @@
epsNearCriticalPointProcess = 1e-5
-- Also the distribution of origins needs a Tailor expansion for near critical values.
---
--- TODO: Check why the two epsilons are chosen differently.
epsNearCriticalTimeOfOrigin :: Double
epsNearCriticalTimeOfOrigin = 1e-8
@@ -130,7 +128,7 @@
Gen (PrimState m) ->
m (PointProcess Int Double)
simulateRandom n l m g
- | -- XXX. There is no formula for the over-critical process.
+ | -- There is no formula for the over-critical process.
m > l =
error
"simulateRandom: Please specify height if mu > lambda."
@@ -139,7 +137,7 @@
m =~= l =
do
!vs <- replicateM (n - 1) (D.genContVar (BDCNTD l) g)
- -- XXX: The length of the root branch will be 0.
+ -- The length of the root branch will be 0.
let t = maximum vs
return $ PointProcess [0 .. (n - 1)] vs t
| -- For the near critical process, we use a special distribution.
@@ -219,7 +217,7 @@
flattenIndices :: [Int] -> [Int]
flattenIndices is = snd $ mapAccumL fAcc [] is
--- TODO: This is the bottleneck for simulating large trees.
+-- NOTE: fAcc is the speed bottleneck for simulating large trees.
--
-- The accumulating function. Count the number of indices which are before the
-- current index and lower than the current index.
@@ -293,7 +291,7 @@
!heights = S.replicate (length ps) 0
!treeRoot = toReconstructedTree' isSorted vsSorted l lvs heights
!h = last vsSorted
- !treeOrigin = applyStem (+ (o - h)) treeRoot
+ !treeOrigin = modifyStem (+ (o - h)) treeRoot
-- Move up the tree, connect nodes when they join according to the point process.
toReconstructedTree' ::
@@ -317,8 +315,8 @@
!hr = hs `S.index` (i + 1)
!dvl = v - hl
!dvr = v - hr
- !tl = applyStem (+ dvl) $ trs `S.index` i
- !tr = applyStem (+ dvr) $ trs `S.index` (i + 1)
+ !tl = modifyStem (+ dvl) $ trs `S.index` i
+ !tr = modifyStem (+ dvr) $ trs `S.index` (i + 1)
!h' = hl + dvl -- Should be the same as 'hr + dvr'.
!tm = Node 0 l [tl, tr]
!trs'' = (S.take i trs S.|> tm) S.>< S.drop (i + 2) trs
diff --git a/src/ELynx/Tree/Splittable.hs b/src/ELynx/Tree/Splittable.hs
--- a/src/ELynx/Tree/Splittable.hs
+++ b/src/ELynx/Tree/Splittable.hs
@@ -24,5 +24,5 @@
class Splittable e where
split :: e -> e
-instance Splittable Double where
- split = (/ 2)
+instance Splittable () where
+ split = id
diff --git a/src/ELynx/Tree/Support.hs b/src/ELynx/Tree/Support.hs
--- a/src/ELynx/Tree/Support.hs
+++ b/src/ELynx/Tree/Support.hs
@@ -18,6 +18,7 @@
Support (fromSupport),
toSupport,
toSupportUnsafe,
+ HasMaybeSupport (..),
HasSupport (..),
-- * Functions on trees
@@ -28,7 +29,6 @@
import Control.DeepSeq
import Data.Aeson
-import Data.Bifoldable
import Data.Bifunctor
import Data.List
import Data.Semigroup
@@ -54,34 +54,43 @@
instance FromJSON Support
+instance HasMaybeSupport Support where
+ getMaybeSupport = Just
+
+instance HasMaybeSupport () where
+ getMaybeSupport = const Nothing
+
instance HasSupport Support where
- getSup = id
- setSup = const
- modSup f = f
+ getSupport = id
+ setSupport = const
+ modifySupport f = f
-- | Return 'Left' if negative.
toSupport :: Double -> Either String Support
-toSupport x | x < 0 = Left $ "Support is negative: " ++ show x ++ "."
- | otherwise = Right $ Support x
+toSupport x
+ | x < 0 = Left $ "Support is negative: " ++ show x ++ "."
+ | otherwise = Right $ Support x
-- | Do not check if value is negative.
toSupportUnsafe :: Double -> Support
toSupportUnsafe = Support
--- | A data type with measurable and modifiable values.
-class HasSupport e where
- getSup :: e -> Support
- setSup :: Support -> e -> e
+-- | Class of data types that may have a support value.
+class HasMaybeSupport e where
+ getMaybeSupport :: e -> Maybe Support
- -- For computational efficiency.
- modSup :: (Support -> Support) -> e -> e
+-- | Class of data types with measurable and modifiable support values.
+class HasMaybeSupport e => HasSupport e where
+ getSupport :: e -> Support
+ setSupport :: Support -> e -> e
+ modifySupport :: (Support -> Support) -> e -> e
-- | Normalize branch support values. The maximum branch support value will be
-- set to 1.0.
normalizeBranchSupport :: HasSupport e => Tree e a -> Tree e a
-normalizeBranchSupport t = first (modSup (/ m)) t
+normalizeBranchSupport t = first (modifySupport (/ m)) t
where
- m = bimaximum $ bimap getSup (const 0) t
+ m = maximum $ getSupport <$> ZipBranchTree t
-- | Collapse branches with support lower than given value.
--
@@ -94,7 +103,7 @@
-- A leaf has full support.
highP :: HasSupport e => Support -> Tree e a -> Bool
highP _ (Node _ _ []) = True
-highP th (Node br _ _) = getSup br >= th
+highP th (Node br _ _) = getSupport br >= th
-- See 'collapse'.
collapse' :: HasSupport e => Support -> Tree e a -> Tree e a
diff --git a/test/ELynx/ClassLaws.hs b/test/ELynx/ClassLaws.hs
new file mode 100644
--- /dev/null
+++ b/test/ELynx/ClassLaws.hs
@@ -0,0 +1,72 @@
+-- |
+-- Module : ELynx.ClassLaws
+-- Description : Unit tests for ELynx.ClassLaws
+-- Copyright : (c) 2021 Dominik Schrempf
+-- License : GPL-3.0-or-later
+--
+-- Maintainer : dominik.schrempf@gmail.com
+-- Stability : experimental
+-- Portability : portable
+--
+-- Creation date: Thu Jul 22 20:39:58 2021.
+module ELynx.ClassLaws
+ ( prop_appl_right,
+ prop_appl_left,
+ prop_appl,
+ prop_appl_func,
+ filterLaws,
+ lawsCheckResult,
+ lawsCheckSpec,
+ )
+where
+
+-- import Control.Comonad
+import Control.Applicative
+import Data.Traversable
+import ELynx.Tree.Arbitrary ()
+import Test.Hspec
+import Test.QuickCheck
+import Test.QuickCheck.Classes
+
+prop_appl_right :: (Applicative f, Eq (f a)) => f a -> f a -> Bool
+prop_appl_right u v = (u *> v) == ((id <$ u) <*> v)
+
+prop_appl_left :: (Applicative f, Eq (f a)) => f a -> f a -> Bool
+prop_appl_left u v = (u <* v) == liftA2 const u v
+
+prop_appl :: (Applicative f, Eq (f a)) => (a -> a -> a) -> f a -> Bool
+prop_appl f x = liftA2 f x x == (f <$> x <*> x)
+
+prop_appl_func :: (Applicative f, Eq (f b)) => (a -> b) -> f a -> Bool
+prop_appl_func f x = fmap f x == (f <$> x)
+
+filterLaws :: [String] -> Laws -> Laws
+filterLaws xs (Laws tn ps) = Laws tn [(n, p) | (n, p) <- ps, n `notElem` xs]
+
+lawsCheckResult :: Laws -> IO Bool
+lawsCheckResult (Laws className properties) =
+ and <$> do
+ for properties $ \(name, p) -> do
+ putStr (className ++ ": " ++ name ++ " ")
+ isSuccess <$> quickCheckResult p
+
+lawsCheckSpec :: Laws -> Spec
+lawsCheckSpec (Laws className properties) =
+ parallel $
+ describe className $
+ mapM_ (\(name, p) -> it name (property p)) properties
+
+-- -- TODO: Comonad laws.
+-- --
+-- -- See https://hackage.haskell.org/package/comonad/docs/Control-Comonad.html.
+--
+-- -- Requires: {-# LANGUAGE QuantifiedConstraints #-}
+-- comonadLaw ::
+-- forall proxy f.
+-- ( Comonad f,
+-- Functor f,
+-- forall a. Eq a => Eq (f a),
+-- forall a. Show a => Show (f a),
+-- forall a. Arbitrary a => Arbitrary (f a)
+-- )
+-- comonadLaw = undefined
diff --git a/test/ELynx/Topology/Arbitrary.hs b/test/ELynx/Topology/Arbitrary.hs
new file mode 100644
--- /dev/null
+++ b/test/ELynx/Topology/Arbitrary.hs
@@ -0,0 +1,48 @@
+{-# OPTIONS_GHC -Wno-orphans #-}
+
+-- |
+-- Module : ELynx.Topology.Arbitrary
+-- Description : Arbitrary instances for topologies
+-- Copyright : (c) 2021 Dominik Schrempf
+-- License : GPL-3.0-or-later
+--
+-- Maintainer : dominik.schrempf@gmail.com
+-- Stability : experimental
+-- Portability : portable
+--
+-- Creation date: Thu Jul 22 21:02:46 2021.
+module ELynx.Topology.Arbitrary
+ (
+ )
+where
+
+import qualified Data.List.NonEmpty as N
+import Data.Traversable
+import ELynx.Topology
+import Test.QuickCheck
+
+instance Arbitrary1 Topology where
+ liftArbitrary arbL = go
+ where
+ go = sized $ \n -> do
+ pars <- frequency [(1, pure [1, 1]), (4, arbPartition (n - 1))]
+ case pars of
+ [] -> Leaf <$> arbL
+ xs -> do
+ frst <- for xs $ \i -> resize i go
+ return $ Node $ N.fromList frst
+ arbPartition :: Int -> Gen [Int]
+ arbPartition k = case compare k 1 of
+ LT -> pure []
+ EQ -> pure [1]
+ GT -> do
+ first <- elements [1 .. k]
+ rest <- arbPartition $ k - first
+ return $ first : rest
+
+instance Arbitrary a => Arbitrary (Topology a) where
+ arbitrary = arbitrary1
+
+instance (CoArbitrary a) => CoArbitrary (Topology a) where
+ coarbitrary (Node frst) = coarbitrary $ N.toList frst
+ coarbitrary (Leaf val) = coarbitrary val
diff --git a/test/ELynx/Topology/RootedSpec.hs b/test/ELynx/Topology/RootedSpec.hs
--- a/test/ELynx/Topology/RootedSpec.hs
+++ b/test/ELynx/Topology/RootedSpec.hs
@@ -14,10 +14,29 @@
)
where
+import Data.Proxy
+import ELynx.ClassLaws
+import ELynx.Topology
+import ELynx.Topology.Arbitrary ()
import Test.Hspec
+import Test.QuickCheck
+import Test.QuickCheck.Classes
+type T = Topology Double
+
spec :: Spec
-spec =
- describe "TODO" $
- it "returns the original number when given a positive input" $
- (1 :: Int) `shouldBe` 1
+spec = do
+ describe "Topology" $ do
+ it "has reasonable applicative take right instance" $
+ property (prop_appl_right :: T -> T -> Bool)
+ it "has reasonable applicative take left instance" $
+ property (prop_appl_left :: T -> T -> Bool)
+ it "has reasonable applicative liftA2 instance" $
+ property (prop_appl (*) :: T -> Bool)
+ it "has reasonable applicative and functor instances" $
+ property (prop_appl_func (+ 3) :: T -> Bool)
+ lawsCheckSpec (functorLaws (Proxy :: Proxy Topology))
+ lawsCheckSpec (foldableLaws (Proxy :: Proxy Topology))
+ lawsCheckSpec (traversableLaws (Proxy :: Proxy Topology))
+ lawsCheckSpec (applicativeLaws (Proxy :: Proxy Topology))
+ lawsCheckSpec (monadLaws (Proxy :: Proxy Topology))
diff --git a/test/ELynx/Tree/Arbitrary.hs b/test/ELynx/Tree/Arbitrary.hs
--- a/test/ELynx/Tree/Arbitrary.hs
+++ b/test/ELynx/Tree/Arbitrary.hs
@@ -2,7 +2,7 @@
-- |
-- Module : ELynx.Tree.Arbitrary
--- Description : Arbitrary instance, needed for QuickCheck
+-- Description : Arbitrary instances for trees
-- Copyright : (c) Dominik Schrempf 2021
-- License : GPL-3.0-or-later
--
@@ -23,7 +23,7 @@
-- Of course, the boundaries for branch support and length have been chosen
-- pretty arbitrarily :).
--
--- XXX: This instance does not produce values without branch lengths nor branch
+-- NOTE: This instance does not produce values without branch lengths nor branch
-- supports.
instance Arbitrary Phylo where
arbitrary =
@@ -38,7 +38,7 @@
-- Sized is the size of the trees.
br <- arbB
val <- arbN
- pars <- frequency [(1, pure [1, 1]), (3, arbPartition (n - 1))] -- can go negative!
+ pars <- frequency [(1, pure [1, 1]), (1, pure [1]), (4, arbPartition (n - 1))]
frst <- for pars $ \i -> resize i go
return $ Node br val frst
arbPartition :: Int -> Gen [Int]
@@ -50,14 +50,6 @@
rest <- arbPartition $ k - first
return $ first : rest
- liftShrink2 _ shrN = go
- where
- go (Node br val frst) =
- frst
- ++ [ Node br e fs
- | (e, fs) <- liftShrink2 shrN (liftShrink go) (val, frst)
- ]
-
instance (Arbitrary e, Arbitrary a) => Arbitrary (Tree e a) where
arbitrary = arbitrary2
@@ -67,3 +59,33 @@
instance Arbitrary Length where
arbitrary = toLengthUnsafe . getPositive <$> arbitrary
+
+instance Arbitrary2 BranchTree where
+ liftArbitrary2 ga gb = BranchTree <$> liftArbitrary2 gb ga
+
+instance (Arbitrary a, Arbitrary e) => Arbitrary (BranchTree a e) where
+ arbitrary = arbitrary2
+
+instance (CoArbitrary a, CoArbitrary e) => CoArbitrary (BranchTree a e) where
+ coarbitrary (BranchTree (Node br val frst)) =
+ coarbitrary br . coarbitrary val . coarbitrary frst
+
+instance Arbitrary2 ZipTree where
+ liftArbitrary2 ga gb = ZipTree <$> liftArbitrary2 ga gb
+
+instance (Arbitrary e, Arbitrary a) => Arbitrary (ZipTree e a) where
+ arbitrary = arbitrary2
+
+instance (CoArbitrary e, CoArbitrary a) => CoArbitrary (ZipTree e a) where
+ coarbitrary (ZipTree (Node br val frst)) =
+ coarbitrary br . coarbitrary val . coarbitrary frst
+
+instance Arbitrary2 ZipBranchTree where
+ liftArbitrary2 ga gb = ZipBranchTree <$> liftArbitrary2 gb ga
+
+instance (Arbitrary a, Arbitrary e) => Arbitrary (ZipBranchTree a e) where
+ arbitrary = arbitrary2
+
+instance (CoArbitrary a, CoArbitrary e) => CoArbitrary (ZipBranchTree a e) where
+ coarbitrary (ZipBranchTree (Node br val frst)) =
+ coarbitrary br . coarbitrary val . coarbitrary frst
diff --git a/test/ELynx/Tree/BipartitionSpec.hs b/test/ELynx/Tree/BipartitionSpec.hs
--- a/test/ELynx/Tree/BipartitionSpec.hs
+++ b/test/ELynx/Tree/BipartitionSpec.hs
@@ -82,5 +82,5 @@
it "creates a map from bipartitions to branch lengths" $
do
simpleTrees <- getSimpleTrees
- (phyloToLengthTree (simpleTrees !! 2) >>= bipartitionToBranch)
+ (toLengthTree (simpleTrees !! 2) >>= bipartitionToBranch)
`shouldBe` Right bipartitionToBranchAnswer
diff --git a/test/ELynx/Tree/DistanceSpec.hs b/test/ELynx/Tree/DistanceSpec.hs
--- a/test/ELynx/Tree/DistanceSpec.hs
+++ b/test/ELynx/Tree/DistanceSpec.hs
@@ -180,19 +180,19 @@
]
incSplitTree1a :: Tree Phylo Name
-incSplitTree1a = parseByteStringWith (oneNewick IqTree) "((a,b)0.7,(c,d));"
+incSplitTree1a = parseByteStringWith (oneNewick IqTree) "((a,b)0.7,(c,d)1.0);"
incSplitTree1b :: Tree Phylo Name
incSplitTree1b = parseByteStringWith (oneNewick IqTree) "((a,b)0.7,c,d);"
incSplitTree2 :: Tree Phylo Name
-incSplitTree2 = parseByteStringWith (oneNewick IqTree) "((a,c),(b,d));"
+incSplitTree2 = parseByteStringWith (oneNewick IqTree) "((a,c)1.0,(b,d)1.0);"
incSplitTree3 :: Tree Phylo Name
-incSplitTree3 = parseByteStringWith (oneNewick IqTree) "(((a,b)0.7,c),(d,e));"
+incSplitTree3 = parseByteStringWith (oneNewick IqTree) "(((a,b)0.7,c)1.0,(d,e)1.0);"
incSplitTree4 :: Tree Phylo Name
-incSplitTree4 = parseByteStringWith (oneNewick IqTree) "(((a,c),b),(d,e));"
+incSplitTree4 = parseByteStringWith (oneNewick IqTree) "(((a,c)1.0,b)1.0,(d,e)1.0);"
-- Compute distances between adjacent pairs of a list of input trees. Use given
-- distance measure.
@@ -233,11 +233,11 @@
it "calculates correct distances for sample trees with branch support" $ do
incompatibleSplits incSplitTree1a incSplitTree2 `shouldBe` Right 2
incompatibleSplits incSplitTree1b incSplitTree2 `shouldBe` Right 2
- let t1a = phyloToSupportTreeUnsafe incSplitTree1a
- t1b = phyloToSupportTreeUnsafe incSplitTree1b
- tr2 = phyloToSupportTreeUnsafe incSplitTree2
- tr3 = phyloToSupportTreeUnsafe incSplitTree3
- tr4 = phyloToSupportTreeUnsafe incSplitTree4
+ let t1a = either error id $ toSupportTree incSplitTree1a
+ t1b = either error id $ toSupportTree incSplitTree1b
+ tr2 = either error id $ toSupportTree incSplitTree2
+ tr3 = either error id $ toSupportTree incSplitTree3
+ tr4 = either error id $ toSupportTree incSplitTree4
incompatibleSplits (collapse 0.7 t1a) tr2 `shouldBe` Right 2
incompatibleSplits (collapse 0.71 t1b) tr2 `shouldBe` Right 0
incompatibleSplits (collapse 0.71 tr3) tr4 `shouldBe` Right 0
@@ -248,7 +248,7 @@
describe "branchScore" $ do
it "calculates correct distances for sample trees" $ do
manyTrees <- getManyTrees
- let ts = map (either error id . phyloToLengthTree) manyTrees
+ let ts = map (either error id . toLengthTree) manyTrees
let ds = map (either error id) $ each 2 $ adjacent branchScore ts
ds `shouldSatisfy` nearlyEqListWith 1e-5 branchScoreAnswers
it "is zero for a collection of random trees" $
diff --git a/test/ELynx/Tree/PhylogenySpec.hs b/test/ELynx/Tree/PhylogenySpec.hs
--- a/test/ELynx/Tree/PhylogenySpec.hs
+++ b/test/ELynx/Tree/PhylogenySpec.hs
@@ -17,6 +17,7 @@
)
where
+import Data.Default
import Data.Either
import qualified Data.Set as S
import ELynx.Tree
@@ -50,66 +51,13 @@
Node () "i" [Node () "j" [Node () "z" [], Node () "x" []], Node () "y" []]
]
-instance Splittable () where
- split = id
-
--- Skip leaves and trees with multifurcating root nodes.
-prop_roots :: Tree () a -> Bool
-prop_roots t@(Node _ _ [_, _])
- | length (leaves t) == 2 = (length <$> roots t) == Right 1
- | otherwise = (length <$> roots t) == Right (length (labels t) - 2)
-prop_roots _ = True
-
--- -- Skip leaves and trees with multifurcating root nodes.
--- prop_connect :: a -> Tree () a -> Tree () a -> Bool
--- prop_connect n l@(Node _ _ [_, _]) r@(Node _ _ [_, _])
--- | length (leaves l) < 3 = (length <$> connect n l r) == Right (length (flatten r) - 2)
--- | length (leaves r) < 3 = (length <$> connect n l r) == Right (length (flatten l) - 2)
--- | otherwise =
--- (length <$> connect n l r)
--- == (Right $ (length (flatten l) - 2) * (length (flatten r) - 2))
--- prop_connect _ _ _ = True
-
--- -- | Determine compatibility between a bipartition and a set.
--- --
--- -- If both subsets of the bipartition share elements with the given set, the
--- -- bipartition is incompatible with this subset. If all elements of the subset
--- -- are either not in the bipartition or mapping to one of the two subsets of the
--- -- bipartition, the bipartition and the subset are compatible.
--- --
--- -- See also 'ELynx.Tree.Partition.compatible'.
--- bipartitionCompatible :: (Show a, Ord a) => Either String (Bipartition a) -> Set a -> Bool
--- -- compatible (Bipartition (l, r)) ss = sintersection l ss `sdisjoint` sintersection r ss
--- bipartitionCompatible (Left _) _ = False
--- bipartitionCompatible (Right p) s = S.null lOverlap || S.null rOverlap
--- where
--- (l, r) = fromBipartition p
--- lOverlap = S.intersection l s
--- rOverlap = S.intersection r s
-
--- compatibleAll :: (Show a, Ord a) => Tree e a -> [Set a] -> Bool
--- compatibleAll (Node _ _ [l, r]) cs =
--- all (bipartitionCompatible (bipartition l)) cs && all (bipartitionCompatible (bipartition r)) cs
--- compatibleAll _ _ = error "Tree is not bifurcating."
-
--- compatibleWith ::
--- (Show b, Ord b) => (a -> b) -> [Set a] -> Tree e a -> Bool
--- compatibleWith f cs t = compatibleAll (fmap f t) (map (S.map f) cs)
-
--- -- Get groups induced by multifurcations. Collect the leaves of all trees
--- -- induced by multifurcations.
--- multifurcatingGroups :: Tree e a -> [[a]]
--- multifurcatingGroups (Node _ _ []) = []
--- multifurcatingGroups (Node _ _ [x]) = multifurcatingGroups x
--- multifurcatingGroups (Node _ _ [x, y]) = multifurcatingGroups x ++ multifurcatingGroups y
--- multifurcatingGroups t = leaves t : concatMap multifurcatingGroups (forest t)
-
--- -- TODO.
--- prop_bifurcating_tree
--- :: (Ord a, HasLength a, HasName a, HasSupport a) => Tree a -> Bool
--- prop_bifurcating_tree t = partitions (resolve t) == empty
+prop_roots :: Default a => Tree () a -> Bool
+prop_roots t@(Node _ _ []) = isLeft $ roots t
+prop_roots t@(Node _ _ [_]) = isLeft $ roots t
+prop_roots t@(Node _ _ [_, _]) = (length <$> roots t) == Right (length (labels t) - 2)
+prop_roots t = (length <$> roots t) == Right (length (labels t) - 1)
-prop_roots_total_length :: Tree Length a -> Bool
+prop_roots_total_length :: Default a => Tree Length a -> Bool
prop_roots_total_length t@(Node _ _ [_, _]) =
all (\x -> abs (totalBranchLength x - l) < 1e-8) $
either error id $
@@ -135,39 +83,15 @@
roots tcherry `shouldBe` Right [tcherry]
it "correctly handles simple trees" $
either error id (roots simpleTree1) `shouldBe` simpleSol
- modifyMaxSize (* 100) $
+ modifyMaxSize (* 100) $ do
it "returns the correct number of rooted trees for arbitrary trees" $
property (prop_roots :: (Tree () Int -> Bool))
- describe "outgroup" $
- modifyMaxSize (* 100) $
- it "correctly handles simple trees" $
- do
- let p = fst $ fromBipartition $ either error id $ bipartition simpleTree1
- outgroup p simpleTree1 `shouldBe` Right simpleTree1
- let l = S.singleton "x"
- either error id (outgroup l simpleTree1) `equal` (simpleSol !! 1)
- `shouldBe` Right True
- describe "rootsWithBranch" $
- modifyMaxSize (* 100) $
- it "does not change the tree height" $
+ it "does not change the total tree length" $
property (prop_roots_total_length :: Tree Length Int -> Bool)
-
--- -- TODO: Move this test to the executable.
--- describe "connect" $
--- modifyMaxSize (* 100) $ do
--- it "returns the correct number of rooted trees for arbitrary trees" $
--- property (prop_connect :: Int -> Tree () Int -> Tree () Int -> Bool)
--- it "correctly connects sample trees without and with constraints" $ do
--- a <- parseFileWith (oneNewick Standard) "data/ConnectA.tree"
--- b <- parseFileWith (oneNewick Standard) "data/ConnectB.tree"
--- c <- parseFileWith (someNewick Standard) "data/ConnectConstraints.tree"
--- let ts =
--- either error id $
--- connect "ROOT" (first (const ()) a) (first (const ()) b)
--- cs =
--- map S.fromList $
--- concatMap (multifurcatingGroups . first (const ())) c ::
--- [Set ByteString]
--- ts' = filter (compatibleWith getName cs) ts
--- length ts `shouldBe` 63
--- length ts' `shouldBe` 15
+ describe "outgroup" $
+ it "correctly handles simple trees" $ do
+ let p = fst $ fromBipartition $ either error id $ bipartition simpleTree1
+ outgroup p simpleTree1 `shouldBe` Right simpleTree1
+ let l = S.singleton "x"
+ either error id (outgroup l simpleTree1) `equal` (simpleSol !! 1)
+ `shouldBe` Right True
diff --git a/test/ELynx/Tree/RootedSpec.hs b/test/ELynx/Tree/RootedSpec.hs
--- a/test/ELynx/Tree/RootedSpec.hs
+++ b/test/ELynx/Tree/RootedSpec.hs
@@ -16,11 +16,19 @@
)
where
+import Data.Bifunctor
import qualified Data.ByteString.Lazy.Char8 as BL
import Data.Maybe
+import Data.Monoid
+import Data.Proxy
+import ELynx.ClassLaws
import ELynx.Tools
import ELynx.Tree
+import ELynx.Tree.Arbitrary ()
import Test.Hspec
+-- import Test.Hspec.QuickCheck
+import Test.QuickCheck hiding (labels)
+import Test.QuickCheck.Classes
node :: Int -> Tree () Int
node n = Node () n []
@@ -43,8 +51,120 @@
subSampleLargeTree :: Tree Phylo Name
subSampleLargeTree = fromJust $ dropLeavesWith ((/= 'P') . BL.head . fromName) largeTree
+-- Branch trees treat branches.
+prop_BranchTree_fmap :: (Eq e, Eq f) => (e -> f) -> Tree e e -> Bool
+prop_BranchTree_fmap f t = first f t == getBranchTree (f <$> BranchTree t)
+
+-- Check that the Traversable instances of Tree and BranchTree work the same. I
+-- am pretty confident that the Traversable instance of Tree is correct, so this
+-- should be enough.
+prop_BranchTree_traversable :: Eq e => Tree e a -> Bool
+prop_BranchTree_traversable t = identify t == bt
+ where
+ bt = flipLabels $ getBranchTree $ identify $ BranchTree $ flipLabels t
+
+-- Check that zipping works the same for both instances ZipTree and
+-- ZipBranchTree. However, this check does not verify that either works
+-- correctly :).
+prop_zip :: (Monoid e, Eq e, Eq a) => Tree e a -> Bool
+prop_zip t = flipLabels (getZipBranchTree zbt') == getZipTree znt'
+ where
+ zbt = ZipBranchTree $ flipLabels t
+ zbt' = (,) <$> zbt <*> zbt
+ znt = ZipTree t
+ znt' = (,) <$> znt <*> znt
+
+-- Same as above but for zip trees.
+prop_ZipTrees_traversable :: Eq e => Tree e a -> Bool
+prop_ZipTrees_traversable t = (t' == zbt) && (t' == znt)
+ where
+ t' = identify t
+ zbt = flipLabels $ getZipBranchTree $ identify $ ZipBranchTree $ flipLabels t
+ znt = getZipTree $ identify $ ZipTree t
+
+type T = Tree String Double
+
+type BT = BranchTree String Double
+
+type ZT = ZipTree String Double
+
+type ZBT = ZipBranchTree String Double
+
spec :: Spec
spec = do
+ -- Data types.
+ describe "Tree" $ do
+ it "[Applicative] Reasonable take right instance" $
+ property (prop_appl_right :: T -> T -> Bool)
+ it "[Applicative] Reasonable take left instance" $
+ property (prop_appl_left :: T -> T -> Bool)
+ it "[Applicative] Reasonable liftA2 instance" $
+ property (prop_appl (*) :: T -> Bool)
+ it "[Functor/Applicative] Reasonable fmap/pure functions" $
+ property (prop_appl_func (+ 3) :: T -> Bool)
+ lawsCheckSpec (functorLaws (Proxy :: Proxy (Tree Int)))
+ lawsCheckSpec (bifunctorLaws (Proxy :: Proxy Tree))
+ lawsCheckSpec (foldableLaws (Proxy :: Proxy (Tree Int)))
+ lawsCheckSpec (bifoldableLaws (Proxy :: Proxy Tree))
+ lawsCheckSpec (traversableLaws (Proxy :: Proxy (Tree Int)))
+ lawsCheckSpec (bitraversableLaws (Proxy :: Proxy Tree))
+ lawsCheckSpec (applicativeLaws (Proxy :: Proxy (Tree String)))
+ lawsCheckSpec (monadLaws (Proxy :: Proxy (Tree String)))
+ describe "BranchTree" $ do
+ it "[Functor] Treats branches correctly with fmap" $
+ property (prop_BranchTree_fmap (* 2) :: Tree Double Double -> Bool)
+ it "[Traversable] Equal traversable instance as Tree" $
+ property (prop_BranchTree_traversable :: Tree Int Int -> Bool)
+ it "[Applicative] Reasonable take right instance" $
+ property (prop_appl_right :: BT -> BT -> Bool)
+ it "[Applicative] Reasonable take left instance" $
+ property (prop_appl_left :: BT -> BT -> Bool)
+ it "[Applicative] Reasonable liftA2 instance" $
+ property (prop_appl (*) :: BT -> Bool)
+ it "[Functor/Applicative] Reasonable fmap/pure functions" $
+ property (prop_appl_func (+ 3) :: BT -> Bool)
+ lawsCheckSpec (functorLaws (Proxy :: Proxy (BranchTree Int)))
+ lawsCheckSpec (foldableLaws (Proxy :: Proxy (BranchTree Int)))
+ lawsCheckSpec (traversableLaws (Proxy :: Proxy (BranchTree Int)))
+ lawsCheckSpec (applicativeLaws (Proxy :: Proxy (BranchTree String)))
+ describe "ZipTree" $ do
+ it "[Applicative] Reasonable take right instance" $
+ property (prop_appl_right :: ZT -> ZT -> Bool)
+ it "[Applicative] Reasonable take left instance" $
+ property (prop_appl_left :: ZT -> ZT -> Bool)
+ it "[Applicative] Reasonable liftA2 instance" $
+ property (prop_appl (*) :: ZT -> Bool)
+ it "[Functor/Applicative] Reasonable fmap/pure functions" $
+ property (prop_appl_func (+ 3) :: ZT -> Bool)
+ lawsCheckSpec (functorLaws (Proxy :: Proxy (ZipTree Int)))
+ lawsCheckSpec (foldableLaws (Proxy :: Proxy (ZipTree Int)))
+ lawsCheckSpec (traversableLaws (Proxy :: Proxy (ZipTree Int)))
+ lawsCheckSpec
+ ( filterLaws ["Homomorphism"] $
+ applicativeLaws (Proxy :: Proxy (ZipTree String))
+ )
+ describe "ZipBranchTree" $ do
+ it "[Applicative] Reasonable take right instance" $
+ property (prop_appl_right :: ZBT -> ZBT -> Bool)
+ it "[Applicative] Reasonable take left instance" $
+ property (prop_appl_left :: ZBT -> ZBT -> Bool)
+ it "[Applicative] Reasonable liftA2 instance" $
+ property (prop_appl (*) :: ZBT -> Bool)
+ it "[Functor/Applicative] Reasonable fmap/pure functions" $
+ property (prop_appl_func (+ 3) :: ZBT -> Bool)
+ lawsCheckSpec (functorLaws (Proxy :: Proxy (ZipBranchTree Int)))
+ lawsCheckSpec (foldableLaws (Proxy :: Proxy (ZipBranchTree Int)))
+ lawsCheckSpec (traversableLaws (Proxy :: Proxy (ZipBranchTree Int)))
+ lawsCheckSpec
+ ( filterLaws ["Homomorphism"] $
+ applicativeLaws (Proxy :: Proxy (ZipBranchTree String))
+ )
+ describe "ZipTree and ZipBranchTree" $ do
+ it "[Applicative] Somewhat corresponding instances of <*>" $
+ property (prop_zip :: Tree (Sum Int) Int -> Bool)
+ it "[Traversable] Somewhat correspnding instances of traverse" $
+ property (prop_ZipTrees_traversable :: Tree Int Int -> Bool)
+ -- Functions.
describe "prune" $ do
it "leaves a normal tree untouched" $
prune largeTree `shouldBe` largeTree
@@ -53,7 +173,7 @@
it "leaves height constant for trees with branch lengths" $ do
let t' =
either error id $
- phyloToLengthTree subSampleLargeTree
+ toLengthTree subSampleLargeTree
height (prune t') `shouldBe` height t'
describe "dropLeavesWith" $ do
it "returns the same tree if no leaves satisfy predicate" $
diff --git a/test/ELynx/Tree/SupportSpec.hs b/test/ELynx/Tree/SupportSpec.hs
--- a/test/ELynx/Tree/SupportSpec.hs
+++ b/test/ELynx/Tree/SupportSpec.hs
@@ -21,7 +21,7 @@
import Test.Hspec
collapseTree :: Tree Phylo Name
-collapseTree = parseByteStringWith (oneNewick IqTree) "((a,b),(c,d));"
+collapseTree = parseByteStringWith (oneNewick IqTree) "((a,b)1.0,(c,d)1.0);"
collapseStarTree :: Tree Phylo Name
collapseStarTree = parseByteStringWith (oneNewick Standard) "(a[1.0],b[1.0],c[1.0],d[1.0])[1.0];"
@@ -30,8 +30,8 @@
spec = do
describe "collapse" $ do
it "creates a star tree for 1.0" $ do
- let t = phyloToSupportTreeUnsafe collapseTree
- s = phyloToSupportTreeUnsafe collapseStarTree
+ let t = either error id $ toSupportTree collapseTree
+ s = either error id $ toSupportTree collapseStarTree
collapse 0 t `shouldBe` t
collapse 0.01 t `shouldBe` t
collapse 0.99 t `shouldBe` t