Skip to content

Commit

Permalink
update vignettes
Browse files Browse the repository at this point in the history
  • Loading branch information
sdgamboa committed Dec 4, 2023
1 parent 567f95a commit b09c522
Show file tree
Hide file tree
Showing 2 changed files with 51 additions and 55 deletions.
23 changes: 1 addition & 22 deletions vignettes/articles/uncertainty_binary.Rmd
Original file line number Diff line number Diff line change
@@ -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
Expand Down Expand Up @@ -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))
```

Expand Down
83 changes: 50 additions & 33 deletions vignettes/articles/uncertainty_multistate.Rmd
Original file line number Diff line number Diff line change
@@ -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}
Expand All @@ -11,6 +10,7 @@ knitr::opts_chunk$set(
```

```{r setup}
options(rstudio.viewer.autorefresh = FALSE)
library(ape)
library(phytools)
library(dplyr)
Expand All @@ -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) |>
Expand All @@ -50,6 +48,18 @@ myFun <- function(mat, per_uncertain = 0.1) {
}
```


## Some practical examples


```{r}
```





# No uncertainty

```{r}
Expand All @@ -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
Expand All @@ -112,4 +130,3 @@ title(main = '70% uncertain tips', line = -1)
sessioninfo::session_info()
```


0 comments on commit b09c522

Please sign in to comment.