Skip to content

Commit

Permalink
update Rmd file
Browse files Browse the repository at this point in the history
  • Loading branch information
sdgamboa committed Sep 12, 2023
1 parent 08ea1ee commit 8edb291
Showing 1 changed file with 50 additions and 1 deletion.
51 changes: 50 additions & 1 deletion vignettes/articles/phytools.Rmd
Original file line number Diff line number Diff line change
Expand Up @@ -137,13 +137,62 @@ phys_data_ready <- left_join(set_with_ids, set_without_ids)

```{r}
tree = mpa()
tree_data = mpa(data = TRUE)
```


```{r}
tree_data <- as_tibble(mpa(data = TRUE)) |>
modify(as.character)
head(tree_data)
```


```{r}
tip_data <- left_join(
tree_data, phys_data_ready, by = c('taxid' = 'NCBI_ID'),
relationship = 'many-to-many'
)
m1 <- tip_data |>
filter(!is.na(Attribute)) |>
select(tip_label, Attribute, Score) |>
pivot_wider(
names_from = Attribute, values_from = 'Score', values_fill = 0
) |>
tibble::column_to_rownames(var = 'tip_label') |>
as.matrix()
no_annotated_tips <- tip_data |>
filter(!tip_label %in% rownames(m1)) |>
pull(tip_label)
m2 <- matrix(
data = rep(rep(1/ncol(m1), ncol(m1)), length(no_annotated_tips)),
nrow = length((no_annotated_tips)),
byrow = TRUE,
dimnames = list(rownames = no_annotated_tips, colnames = colnames(m1))
)
m3 <- rbind(m1, m2)
m3 <- m3[tree$tip.label,]
```


```{r}
system.time(
fit <- fitMk(
tree = tree,
x = m3,
model = 'ARD',
pi = "fitzjohn",
logscale = TRUE,
lik.func = "pruning"
)
)
```



Expand Down

0 comments on commit 8edb291

Please sign in to comment.