## ----style, echo = FALSE, results = 'asis'------------------------------------
BiocStyle::markdown()

## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
    echo = TRUE,
    tidy.opts = list(width.cutoff = 80),
    tidy = TRUE
)
# long links
adj_link <- paste0(
    "https://thomazbastiaanssen.github.io/anansi/",
    "articles/adjacency_matrices.html"
)
library(S7)

## ----formula-web-full---------------------------------------------------------
# Setup AnansiWeb with some dummy data
library(anansi)
web <- krebsDemoWeb()

## ----y-x-formula,echo=FALSE---------------------------------------------------
y ~ x

## ----formulate, echo = FALSE--------------------------------------------------
f <- y ~ x * (group_ab + score_a)
f

## ----reformulate, echo=FALSE--------------------------------------------------
update.formula(f, . ~ .)

## ----formula-anansi-full------------------------------------------------------
out <- anansi(web,
    formula = ~ group_ab + score_a,
    groups = "group_ab"
)

## ----plot_full----------------------------------------------------------------
#|
library(patchwork)
library(ggplot2)

full_df <- data.frame(
    citrate = tableX(web)[, "citrate"],
    `cis-aconitate` = tableX(web)[, "cis-aconitate"],
    aconitase = tableY(web)[, "aconitase"]
)

ggplot(full_df) +
    aes(x = aconitase, y = citrate) +
    facet_wrap(~"citrate") +
    ggplot(full_df) +
    aes(x = aconitase, y = cis.aconitate) +
    labs(y = "cis-aconitate") +
    facet_wrap(~"cis-aconitate") +
    plot_layout(axes = "collect", axis_titles = "collect") +
    plot_annotation(
        subtitle =
        ) &
    geom_point(shape = 21, fill = "#ffffbf") &

    coord_equal(ratio = 1, xlim = c(-2.5, 2.5), ylim = c(-2.5, 2.5)) &
    theme_bw()

## ----disj_cat-----------------------------------------------------------------
#|
diff_df <- data.frame(
    group_ab = paste0("Treatment ", metadata(web)$group_ab),
    score_a = metadata(web)$score_a,
    score_a_cat = cut(
        metadata(web)$score_a,
        labels = c("0-20%", "20-40%", "40-60%", "60-80%", "80-100%"),
        quantile(metadata(web)$score_a, probs = seq(0, 1, 1 / 5)),
        include.lowest = TRUE
    ),
    repeated = metadata(web)$repeated,
    sample_id = metadata(web)$sample_id,
    isocitrate = tableX(web)[, "isocitrate"],
    ketoglutarate = tableX(web)[, "ketoglutarate"],
    isocitrate.dehydrogenase = tableY(web)[, "isocitrate dehydrogenase"],
    ketoglutarate.dehydrogenase = tableY(web)[, "ketoglutarate dehydrogenase"],
    succinyl.CoA = tableX(web)[, "succinyl-CoA"],
    succinate = tableX(web)[, "succinate"],
    succinyl.CoA.synthetase = tableY(web)[, "succinyl-CoA synthetase"],
    succinate.dehydrogenase = tableY(web)[, "succinate dehydrogenase"],
    fumarate = tableX(web)[, "fumarate"],
    fumarase = tableY(web)[, "fumarase"]
)

ggplot(diff_df) +
    aes(x = isocitrate.dehydrogenase, y = isocitrate, fill = group_ab) +
    facet_wrap(~group_ab, ) +
    labs(x = "isocitrate dehydrogenase") +
    geom_point(shape = 21, show.legend = FALSE) +
    scale_fill_manual(values = c(
        "Treatment a" = "#d73027",
        "Treatment b" = "#2166ac"
    )) +
    theme_bw() +
    coord_equal(ratio = 1, xlim = c(-2.5, 2.5), ylim = c(-2.5, 2.5))

## ----disj_cont----------------------------------------------------------------
#|
ggplot(diff_df) +
    aes(
        x = ketoglutarate.dehydrogenase,
        y = ketoglutarate,
        fill = pnorm(score_a)
    ) +
    facet_grid(~score_a_cat) +
    geom_point(shape = 21, show.legend = FALSE) +
    scale_fill_gradientn(
        name = "Example health index",
        colors = c("#a50026", "#f46d43", "#ffffbf", "#74add1", "#313695")
    ) +
    labs(x = "ketoglutatate dehydrogenase") +
    theme_bw() +
    guides(fill = guide_colourbar(
        position = "bottom",
        legend.direction = "horizontal"
    )) +
    coord_equal(ratio = 1, xlim = c(-2.5, 2.5), ylim = c(-2.5, 2.5))

## ----disj-lm------------------------------------------------------------------

score_lm <- pairwiseApply(web, FUN = function(x, y) {
    anova(
        lm(y ~ x * (group_ab + score_a), data = metadata(web))
    )$`F value`[5L]
})

group_lm <- pairwiseApply(web, FUN = function(x, y) {
    anova(
        lm(y ~ x * (score_a + group_ab), data = metadata(web))
    )$`F value`[5L]
})

all.equal(score_lm, out$disjointed_score_a_f.values)
all.equal(group_lm, out$disjointed_group_ab_f.values)

## ----emerg_cat----------------------------------------------------------------
#|
ggplot(diff_df) +
    aes(x = succinyl.CoA.synthetase, y = succinyl.CoA, fill = group_ab) +
    facet_wrap(~group_ab) +
    geom_point(shape = 21, show.legend = FALSE) +
    scale_fill_manual(values = c(
        "Treatment a" = "#d73027",
        "Treatment b" = "#2166ac"
    )) +
    labs(x = "succinyl-CoA synthetase", y = "succinyl-CoA") +
    theme_bw() +
    coord_equal(ratio = 1, xlim = c(-2.5, 2.5), ylim = c(-2.5, 2.5))

## ----emerg_cont---------------------------------------------------------------
#|
ggplot(diff_df) +
    aes(x = succinate.dehydrogenase, y = succinate, fill = pnorm(score_a)) +
    facet_grid(~score_a_cat) +
    geom_point(shape = 21, show.legend = FALSE) +
    scale_fill_gradientn(
        name = "Example health index",
        colors = c("#a50026", "#f46d43", "#ffffbf", "#74add1", "#313695")
    ) +
    labs(x = "succinate dehydrogenase") +
    theme_bw() +
    coord_equal(ratio = 1, xlim = c(-2.5, 2.5), ylim = c(-2.5, 2.5))

## ----intro-error--------------------------------------------------------------
outErr <- anansi(web,
    formula = ~ group_ab + score_a + Error(sample_id),
    groups = "group_ab"
)

## ----aov-error----------------------------------------------------------------
summary(
    aov(fumarate ~ fumarase * repeated + Error(sample_id), data = diff_df)
)

## ----aov-regular--------------------------------------------------------------
summary(
    aov(fumarate ~ sample_id + fumarase * repeated, data = diff_df)
)

# See ?stats::aov for more.

## ----session-info-------------------------------------------------------------
sessionInfo()