## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.dim = c(6, 4)
)

## ----setup--------------------------------------------------------------------
library(EZbakR)
library(dplyr)

## ----echo = F, results = 'hide'-----------------------------------------------

# Example metadf:
metadf <- tibble(
  sample = c("sampleA", "sampleB", "sampleC",
             "sampleD", "sampleE", "sampleF"),
  tl = c(0, 2, 2,
         0, 2, 2),
  genotype = c("WT", "WT", "WT",
               "KO", "KO", "KO")
)

metadf

## ----echo = FALSE, warning = FALSE--------------------------------------------

knitr::kable(metadf, "pipe")


## -----------------------------------------------------------------------------
metadf <- tibble(sample = c("sampleA", "sampleB", "sampleC",
                            "sampleD", "sampleE", "sampleF"),
                 tl = c(0, 2, 2,
                        0, 2, 2),
                 genotype = c("WT", "WT", "WT",
                               "KO", "KO", "KO"))


# Simulate kdeg and ksyn differences
set.seed(43)
simdata <- EZSimulate(500,
                      metadf = metadf %>% dplyr::rename(label_time = tl),
                      mean_formula = ~genotype -1,
                      pdiff_ks = 0.1)


ezbdo <- EZbakRData(simdata$cB,
                    metadf)


ezbdo <- EstimateFractions(ezbdo)

ezbdo <- EstimateKinetics(ezbdo)

### Linear modeling in EZbakR:

# Estimate average log(kdeg) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              formula_mean = ~ genotype)

# Estimate average log(ksyn) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              parameter = "log_ksyn",
                              formula_mean = ~ genotype)


# Calculate difference in log(kdeg) between the two genotypes
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           design_factor = "genotype",
                           reference = "WT",
                           experimental = "KO")

# Calculate differe in log(ksyn) between the two genotypes
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_ksyn",
                           design_factor = "genotype",
                           reference = "WT",
                           experimental = "KO")


## ----results = 'hold'---------------------------------------------------------
# Assess log(kdeg) differences
EZVolcanoPlot(ezbdo,
              parameter = "log_kdeg")
EZMAPlot(ezbdo,
              parameter = "log_kdeg")

## ----results = 'hold'---------------------------------------------------------
# Assess log(ksyn) differences
EZVolcanoPlot(ezbdo,
              parameter = "log_ksyn")
EZMAPlot(ezbdo,
              parameter = "log_ksyn")

## ----echo = F, results = 'hide'-----------------------------------------------

# Example metadf:
metadf <- tibble(
  sample = c("sampleA", "sampleB", "sampleC",
             "sampleD", "sampleE", "sampleF",
             "sampleG", "sampleH", "sampleI",
             "sampleJ", "sampleK", "sampleL"),
  tl = c(0, 2, 2,
         0, 2, 2,
         0, 2, 2,
         0, 2, 2),
  genotype = c("WT", "WT", "WT",
               "WT", "WT", "WT",
               "KO", "KO", "KO",
               "KO", "KO", "KO"),
  treatment = c("nodrug", "nodrug", "nodrug",
                "drug", "drug", "drug",
                "nodrug", "nodrug", "nodrug",
                "drug", "drug", "drug")
)

metadf

## ----echo = FALSE, warning = FALSE--------------------------------------------

knitr::kable(metadf, "pipe")


## -----------------------------------------------------------------------------
metadf <- tibble(
  sample = c("sampleA", "sampleB", "sampleC",
             "sampleD", "sampleE", "sampleF",
             "sampleG", "sampleH", "sampleI",
             "sampleJ", "sampleK", "sampleL"),
  tl = c(0, 2, 2,
         0, 2, 2,
         0, 2, 2,
         0, 2, 2),
  genotype = c("WT", "WT", "WT",
               "WT", "WT", "WT",
               "KO", "KO", "KO",
               "KO", "KO", "KO"),
  treatment = c("nodrug", "nodrug", "nodrug",
                "drug", "drug", "drug",
                "nodrug", "nodrug", "nodrug",
                "drug", "drug", "drug")
)


# Simulate kdeg and ksyn differences
set.seed(43)
simdata <- EZSimulate(500,
                      metadf = metadf %>% dplyr::rename(label_time = tl),
                      mean_formula = ~genotype:treatment -1,
                      pdiff_ks = 0.1)


ezbdo <- EZbakRData(simdata$cB,
                    metadf)


ezbdo <- EstimateFractions(ezbdo)

ezbdo <- EstimateKinetics(ezbdo)

### Linear modeling in EZbakR:

# Estimate average log(kdeg) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              formula_mean = ~ genotype:treatment)

# Estimate average log(ksyn) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              parameter = "log_ksyn",
                              formula_mean = ~ genotype:treatment)


# Calculate difference in log(kdeg) between WT drug and no drug
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           design_factor = "",
                           reference = "genotypeWT:treatmentnodrug",
                           experimental = "genotypeWT:treatmentdrug")

# Calculate difference in log(kdeg) between KO drug and no drug
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           design_factor = "",
                           reference = "genotypeKO:treatmentnodrug",
                           experimental = "genotypeKO:treatmentdrug")

# Calculate difference in log(kdeg) between KO and WT nodrug
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           design_factor = "",
                           reference = "genotypeWT:treatmentnodrug",
                           experimental = "genotypeKO:treatmentnodrug")


### You could do same for log_ksyn


## -----------------------------------------------------------------------------
ezbdo$metadata$averages$logkdeg_feature$fit_params

## ----results = 'hide'---------------------------------------------------------
metadf$condition <- factor(paste0(metadf$genotype, metadf$treatment))
metadf

## ----echo = FALSE, warning = FALSE--------------------------------------------

knitr::kable(metadf, "pipe")


## -----------------------------------------------------------------------------

# factor() ensures that WT + nodrug is treated as the reference
# which probably makes more intuitive sense
metadf <- tibble(
  sample = c("sampleA", "sampleB", "sampleC",
             "sampleD", "sampleE", "sampleF",
             "sampleG", "sampleH", "sampleI",
             "sampleJ", "sampleK", "sampleL"),
  tl = c(0, 2, 2,
         0, 2, 2,
         0, 2, 2,
         0, 2, 2),
  genotype = factor(c("WT", "WT", "WT",
               "WT", "WT", "WT",
               "KO", "KO", "KO",
               "KO", "KO", "KO"),
               levels = c("WT", "KO")),
  treatment = factor(c("nodrug", "nodrug", "nodrug",
                "drug", "drug", "drug",
                "nodrug", "nodrug", "nodrug",
                "drug", "drug", "drug"),
                levels = c("nodrug", "drug"))
)


# Simulate kdeg and ksyn differences
set.seed(43)
simdata <- EZSimulate(500,
                      metadf = metadf %>% dplyr::rename(label_time = tl),
                      mean_formula = ~genotype*treatment -1,
                      pdiff_ks = 0.1)


ezbdo <- EZbakRData(simdata$cB,
                    metadf)


ezbdo <- EstimateFractions(ezbdo)

ezbdo <- EstimateKinetics(ezbdo)

### Linear modeling in EZbakR:

# Estimate average log(kdeg) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              formula_mean = ~ genotype*treatment)

# Estimate average log(ksyn) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              parameter = "log_ksyn",
                              formula_mean = ~ genotype*treatment)


# Viewing parameters will be easier in the future
ezbdo$metadata$averages$logkdeg_feature$fit_params


## -----------------------------------------------------------------------------
# Assess significance of interaction term
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           param_name = "genotypeKO:treatmentdrug")

# Assess significance of drug effect on WT
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           param_name = "treatmentdrug")

# Assess difference of WT and KO
ezbdo <- CompareParameters(ezbdo,
                           parameter = "log_kdeg",
                           design_factor = "genotype",
                           reference = "WT",
                           experimental = "KO")


## ----echo = F, hide = T-------------------------------------------------------
metadf <- tibble(sample = c("sampleA", "sampleB", "sampleC", "sampleD", "sampleE",
                            "sampleF", "sampleG", "sampleH", "sampleI", "sampleJ"),
                 tl = c(0, 2, 2, 2, 2,
                        0, 2, 2, 2, 2),
                 genotype = c("WT","WT", "WT", "WT", "WT",
                              "KO","KO", "KO", "KO", "KO"),
                 batch = c("A", "A", "B", "A", "B",
                           "A", "A", "B", "A", "B"))

## ----echo = FALSE, warning = FALSE--------------------------------------------

knitr::kable(metadf, "pipe")


## -----------------------------------------------------------------------------
# factor() ensures that WT is treated as the reference
# which probably makes more intuitive sense
metadf <- tibble(sample = c("sampleA", 
                            "sampleB", "sampleC", "sampleD", "sampleE",
                            "sampleF", 
                            "sampleG", "sampleH", "sampleI", "sampleJ"),
                 tl = c(0, 
                        2, 2, 2, 2,
                        0, 
                        2, 2, 2, 2),
                 genotype = factor(c("WT",
                              "WT", "WT", "WT", "WT",
                              "KO",
                              "KO", "KO", "KO", "KO"),
                              levels = c("WT", "KO")),
                 batch = c("A", 
                           "A", "B", "A", "B",
                           "A", 
                           "A", "B", "A", "B"))

# Simulate kdeg and ksyn differences
set.seed(43)
simdata <- EZSimulate(500,
                      metadf = metadf %>% dplyr::rename(label_time = tl),
                      mean_formula = ~genotype + batch -1,
                      pdiff_ks = 0.1)


ezbdo <- EZbakRData(simdata$cB,
                    metadf)


ezbdo <- EstimateFractions(ezbdo)

ezbdo <- EstimateKinetics(ezbdo)

### Linear modeling in EZbakR:

# Estimate average log(kdeg) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              formula_mean = ~ genotype + batch)

# Estimate average log(ksyn) for each genotype
ezbdo <- AverageAndRegularize(ezbdo,
                              parameter = "log_ksyn",
                              formula_mean = ~ genotype + batch)


# Viewing parameters will be easier in the future
ezbdo$metadata$averages$logkdeg_feature$fit_params


## ----echo = F, hide = T-------------------------------------------------------
metadf <- tibble(sample = c("sampleA", "sampleB", "sampleC", "sampleD", "sampleE",
                            "sampleF", "sampleG", "sampleH", "sampleI", "sampleJ"),
                 tl = c(0, 2, 2, 2, 2,
                        0, 2, 2, 2, 2),
                 genotype = c("WT","WT", "WT", "WT", "WT",
                              "KO","KO", "KO", "KO", "KO"),
                 batch = c("A", "A", "A", "A", "A",
                           "B", "B", "B", "B", "B"))

## ----echo = FALSE, warning = FALSE--------------------------------------------

knitr::kable(metadf, "pipe")