From b09c522945ddab3898a2fa99a86ff9ebffcaa0e9 Mon Sep 17 00:00:00 2001 From: sdgamboa Date: Mon, 4 Dec 2023 18:36:43 -0500 Subject: [PATCH] update vignettes --- vignettes/articles/uncertainty_binary.Rmd | 23 +---- vignettes/articles/uncertainty_multistate.Rmd | 83 +++++++++++-------- 2 files changed, 51 insertions(+), 55 deletions(-) diff --git a/vignettes/articles/uncertainty_binary.Rmd b/vignettes/articles/uncertainty_binary.Rmd index c36717c..b6968fb 100644 --- a/vignettes/articles/uncertainty_binary.Rmd +++ b/vignettes/articles/uncertainty_binary.Rmd @@ -1,5 +1,5 @@ --- -title: "Uncertainty in tips for binary attributes with 'TRUE' and 'FALSE' annotations" +title: "Uncertainty in tips with binary attributes" output: html_document: toc: true @@ -52,27 +52,6 @@ original <- data |> ) |> as.matrix() colnames(original) <- c('A--TRUE', 'A--FALSE') - - -# myFun <- function(mat, uncertainty = 0.7, input_value = 0) { -# m <- mat -# if (input_value == 0.5) { -# m[] <- input_value -# } else if (input_value == 0) { -# m[,1] <- 0 -# m[,2] <- 1 -# } -# n <- round(nrow(mat) * (1 - uncertainty)) -# a <- round(n / 3) -# b <- n - a -# nm1 <- sample(rownames(mat), b) -# nm2 <- sample(rownames(mat)[!rownames(mat) %in% nm1], a) -# m[nm1, 1] <- 1 -# m[nm1, 2] <- 0 -# m[nm2, 1] <- 0 -# m[nm2, 2] <- 1 -# return(m) -# } head(sort(table(sub('^(\\w+)_.*$', '\\1', rownames(original))), decreasing = TRUE)) ``` diff --git a/vignettes/articles/uncertainty_multistate.Rmd b/vignettes/articles/uncertainty_multistate.Rmd index 0f6bf88..892774b 100644 --- a/vignettes/articles/uncertainty_multistate.Rmd +++ b/vignettes/articles/uncertainty_multistate.Rmd @@ -1,6 +1,5 @@ --- -title: "Effect of uncertainty in tip annotations" -subtitle: 'Multistate-intersection' +title: "Uncertainty in tips with multistate(-intersection) annotations" --- ```{r, include = FALSE} @@ -11,6 +10,7 @@ knitr::opts_chunk$set( ``` ```{r setup} +options(rstudio.viewer.autorefresh = FALSE) library(ape) library(phytools) library(dplyr) @@ -25,8 +25,6 @@ data('primate.data') tree <- primate.tree data <- primate.data data <- data[tree$tip.label,] -rownames(data) <- paste0('taxon', 1:nrow(data)) -tree$tip.label <- paste0('taxon', 1:Ntip(tree)) original <- data |> tibble::rownames_to_column(var = 'Taxa') |> select(Taxa, Activity_pattern) |> @@ -50,6 +48,18 @@ myFun <- function(mat, per_uncertain = 0.1) { } ``` + +## Some practical examples + + +```{r} + +``` + + + + + # No uncertainty ```{r} @@ -61,49 +71,57 @@ plot(ace, args.plotTree = list(direction="upwards")) title(main = '0% uncertain tips', line = -1) ``` -# Uncertaintity - -## 99% uncertainty +## Uncertainty about 12% ```{r} -set.seed(123) -m99 <- myFun(original, 0.99) -m99 <- m99[tree$tip.label,] -fit99 <- fitMk(tree = tree, x = m99, - model = "ARD", pi = "fitzjohn", - lik.func = "pruning", logscale = TRUE) -ace99 <- ancr(fit99, tips=TRUE) -plot(ace99, args.plotTree = list(direction="upwards")) -title(main = '99% uncertain tips', line = -1) +macaca <- grep('^Macaca_', rownames(original)) # A +galago <- grep('^Galago_', rownames(original)) # B +eulemur <- grep('^Eulemur', rownames(original)) # C +paste0(round((length(macaca) + length(galago) + length(eulemur)) / 90 * 100), '%') ``` -## 95% uncertainty + ```{r} -set.seed(123) -m95 <- myFun(original, 0.95) -m95 <- m95[tree$tip.label,] -fit95 <- fitMk(tree = tree, x = m95, +m1 <- original +m1[] <- 1 / 3 +m1[macaca, ] <- c(rep(1,3), rep(0, 6)) +m1[galago, ] <- c(rep(0,4), rep(1, 4), rep(0, 4)) +m1[eulemur, ] <- c(rep(0,8), rep(1, 4)) + +fit1 <- fitMk(tree = tree, x = m1, model = "ARD", pi = "fitzjohn", lik.func = "pruning", logscale = TRUE) -ace95 <- ancr(fit95, tips=TRUE) -plot(ace95, args.plotTree = list(direction="upwards")) -title(main = '95% uncertain tips', line = -1) +ace1 <- ancr(fit1, tips=TRUE) +plot(ace1, args.plotTree = list(direction = "upwards")) +title(main = '12% uncertain tips', line = -1) +``` + + +Let's add some more taxa +```{r} +cebus <- grep('^Cebus_', rownames(original)) # A +microcebus <- grep('^Microcebus', rownames(original)) # B +paste0(round((length(macaca) + length(galago) + length(eulemur) + length(cebus) + length(microcebus)) / 90 * 100), '%') ``` -## 70% uncertainty ```{r} -set.seed(123) -m75 <- myFun(original, 0.7) -m75 <- m75[tree$tip.label,] -fit75 <- fitMk(tree = tree, x = m75, +m2 <- original +m2[] <- 1 / 3 +m2[macaca, ] <- c(rep(1,3), rep(0, 6)) +m2[galago, ] <- c(rep(0,4), rep(1, 4), rep(0, 4)) +m2[eulemur, ] <- c(rep(0,8), rep(1, 4)) +m2[cebus, ] <- c(rep(1,2), rep(0, 4)) +m2[microcebus, ] <- c(rep(0,2), rep(1, 2), rep(0, 2)) + +fit2 <- fitMk(tree = tree, x = m2, model = "ARD", pi = "fitzjohn", lik.func = "pruning", logscale = TRUE) -ace75 <- ancr(fit75, tips=TRUE) -plot(ace75, args.plotTree = list(direction="upwards")) -title(main = '70% uncertain tips', line = -1) +ace2 <- ancr(fit2, tips=TRUE) +plot(ace2, args.plotTree = list(direction = "upwards")) +title(main = '17% uncertain tips', line = -1) ``` # Session information @@ -112,4 +130,3 @@ title(main = '70% uncertain tips', line = -1) sessioninfo::session_info() ``` -