Skip to content

Instantly share code, notes, and snippets.

This file has been truncated, but you can view the full file.
#NEXUS
begin DATA;
dimensions ntax=235 nchar=137324;
format datatype=nucleotide missing=? gap=-;
matrix
Acanthisitta_chloris_genome NNCTCTCAGGATCCAGACTTCACACAGGTT--AGAAACTACT-C-T-GCCTTCTGGATCC-T-TGCT----------------TT-ATGCCC---T-GTCCTG-TTTTATTGTGACAGTTGCCTG--------------TTTTT-CACAGAA--------------AC-TAAG-----------GCCTAGTGA---------CTCCTTATCCTTTGG--TTCGGGCA--ACACAGCTGTC-TTATCTTAAGGCCCAGTGGAATGA--------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------G-TCAGCTGGCTC-CTCAGGTGTGA--------AAAC------------TTGGG-CATAG------GTATAGATATGCTCT-TAATTGACCTC-TAGCTG------------GTG--
@roblanf
roblanf / code.R
Last active January 3, 2020 10:58
A look at gene- and site-concordance factors. See also https://www.biorxiv.org/content/early/2018/12/05/487801
library(viridis)
library(ggplot2)
library(dplyr)
library(ggrepel)
library(GGally)
library(entropy)
# read the data
d = read.delim("concord.cf.stat", header = T, skip = 15)
# infer the concatentated tree with 1000 ultrafast bootstraps and an edge-linked fully-partitioned model (-p is the same as -spp from version 1.7 onwards)
iqtree -s alignment.nex -p alignment.nex --prefix concat -bb 1000 -nt AUTO
# infer the 88 single-locus trees
iqtree -s alignment.nex -S alignment.nex --prefix loci -nt 50
# calculate concordance factors
iqtree -t concat.treefile --gcf loci.treefile -s alignment.nex --scf 100 --prefix concord
library(viridis)
library(ggplot2)
library(dplyr)
library(ggrepel)
library(GGally)
library(entropy)
# read the data
d = read.delim("concord.cf.stat", header = T, comment.char=‘#')
ID gCF gDF1 gDF2 gN sCF sDF1 sDF2 sN bootstrap branchlength
327 1.15 3.45 0 87 37.34 29.98 32.68 350.19 100 0.00304084
ID gCF gDF1 gDF2 gN sCF sDF1 sDF2 sN bootstrap branchlength gEF_p sEF_p
14 249 78.57 5.36 5.36 56 56.14 17.86 26.00 383.86 100 0.02134550 1.000000000 1.605600e-02
33 268 62.07 8.62 13.79 58 39.02 27.36 33.62 624.87 100 0.02370450 0.405380556 4.506223e-02
34 269 87.50 0.00 0.00 56 69.10 13.40 17.50 771.13 100 0.10990200 1.000000000 4.051930e-02
36 271 67.35 0.00 8.16 49 46.61 25.97 27.42 426.92 100 0.02632530 0.045500264 6.818018e-01
39 274 34.33 17.91 17.91 67 39.51 23.62 36.87 258.95 100 0.00694057 1.000000000 6.116524e-03
@roblanf
roblanf / mre_for_monophy.rmd
Created September 12, 2022 02:04
What is happening with is.monophyletic()
``` r
library(ape)
# make 1K trees to ensure we get some that seem odd...
trees = lapply(rep(4, 1000), rtree, rooted = FALSE)
# test monophyly of all six possible pairs
m1 = unlist(lapply(trees, is.monophyletic, tips = c("t1", "t2")))
m2 = unlist(lapply(trees, is.monophyletic, tips = c("t1", "t3")))
m3 = unlist(lapply(trees, is.monophyletic, tips = c("t1", "t4")))

Commit Message Guidelines

Short (72 chars or less) summary

More detailed explanatory text. Wrap it to 72 characters. The blank
line separating the summary from the body is critical (unless you omit
the body entirely).

Write your commit message in the imperative: "Fix bug" and not "Fixed
bug" or "Fixes bug." This convention matches up with commit messages