knitr::opts_chunk$set( collapse = TRUE, comment = "#>", fig.align = "center", eval = TRUE )
library(zazou) library(tictoc)
data(chlamydiae)
pval_obs <- test_kruskalwallis(chlamydiae$X, chlamydiae$Y)$p.value zsco_obs <- p2z(pval_obs)
tree <- force_ultrametric(chlamydiae$tree) N_branch <- length(tree$edge.length) plot_shifts(tree, NA, obs_scores = zsco_obs)
tic() est_lbfgsb <- estimate_shifts(zscores = zsco_obs, lambda = 0.1, tree = tree, alphaOU = 1, method = "L-BFGS-B") toc() est_lbfgsb plot(est_lbfgsb)
set.seed(42) tic() est_shooting <- estimate_shifts(zscores = zsco_obs, lambda = 0.1, tree = tree, alphaOU = 1, method = "lasso") toc() est_shooting plot(est_shooting)
est_lbfgsb$objective_value est_shooting$objective_value est_lbfgsb$optim_info est_shooting$optim_info
library(ggplot2) theme_set(theme_minimal()) qplot(est_lbfgsb$zscores_obs, est_lbfgsb$zscores_est) + geom_abline() + geom_vline(xintercept = 0) qplot(est_shooting$zscores_obs, est_shooting$zscores_est) + geom_abline() + geom_vline(xintercept = 0) qplot(-est_shooting$zscores_est, -est_lbfgsb$zscores_est) + geom_abline() + scale_x_log10() + scale_y_log10()
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.