In Sage-Bionetworks/sageseqr: Identify differentially expressed genes counts data

```{css, echo = FALSE} / Whole document: / body{ font-family: Lato; font-size: 12pt; } / Headers / h1,h2,h3,h4,h5,h6{ font-size: 24pt; }

```r
project <- config::get("report")

---

title: r project

knitr::opts_chunk$set(echo = FALSE, 
                      fig.align = 'left',
                      fig.width = 10, 
                      fig.height = 10)
library(sageseqr)

Data

Explore Metadata

targets::tar_load(clean_md)
var <- config::get("x_var")
summary <- dplyr::group_by_at(clean_md, var)
summary <- dplyr::summarise(summary, count = dplyr::n())
knitr::kable(summary)

Visualize the distribution of data.

targets::tar_read(boxplots)

if (!is.null(config::get("sex check"))) {
  cat("Visualize gene expression of sex chromosomes using XIST 
      as a X-chromosome marker and UTY as a Y-chromosome marker.")
  targets::tar_read(sex_plot)$plot
}

if (!is.null(config::get("sex check"))) {
   cat("Visualize gene expression across X and Y chromosomes by principal component analysis (PCA).")
   targets::tar_read(sex_plot)$plot
}

Sample swaps identified by discordant sex label, are r glue::glue_collapse(targets::tar_read(sex_plot_pca)$discordant, ", ", last = " and ").

if (!is.null(targets::tar_read(sex_plot_pca)$outliers)) {
  list <- glue::glue_collapse(targets::tar_read(sex_plot_pca)$outliers, ", ", last = " and ")
  glue::glue("Possible outliers, identified as 2 standard deviations from the mean of samples n split by {config::get('sex check')} are {list}.")
} else {
  cat("No outliers identified.")
}

Visualize the relationships between covariates.

correlation_input <- targets::tar_read(
  correlation_plot
  )$plot
col2 <- grDevices::colorRampPalette(rev(c("#67001F", "#B2182B", "#D6604D", "#F4A582",
                                 "#FDDBC7", "#FFFFFF", "#D1E5F0", "#92C5DE",
                                 "#4393C3", "#2166AC", "#053061")))
corrplot::corrplot(correlation_input, col = col2(200), tl.col = "#000000")

Filter Genes

Remove genes that have less than r config::get("cpm threshold") counts per million (CPM) in at least r config::get("percent threshold") of samples per specified condition.

r dim(targets::tar_read(filtered_counts))[1] genes are used in this analysis.

knitr::kable(targets::tar_read(biotypes))

Check distribution of correlation between genes.

targets::tar_read(gene_coexpression)

Identify Outliers

targets::tar_read(outliers)$plot

Outliers, based on logCPM expression, are r glue::glue_collapse(targets::tar_read(outliers)$outliers, ", ", last = " and ").

Identify Dropped Gene Features

Gene features dropped, because they are not annotated in BioMart, are r glue::glue_collapse(targets::tar_read(dropped), ", ", last = " and ").

Significant Covariates

Significant covariates are identified by the pearson correlation (p-value of 1%) between principal component analysis (PCA) of normalized transcripts and variables that meet a 0.1 false discovery rate (FDR) threshold. Significant covariates to adjust for are r glue::glue_collapse(targets::tar_read(significant_covariates_plot)$significant_covariates, ", ", last = " and ").

targets::tar_read(
  significant_covariates_plot
  )$pc_results

Model Identification

Covariates are added as fixed and random effects iteratively if model improvement by Bayesian Information Criteria (BIC) was observed.

if (!isTRUE(config::get("skip model"))) {
  summary <- targets::tar_read(model)
  glue::glue('The best model is: {glue::glue_collapse(as.character(summary$formula), " ")}')
} else {
  glue::glue(

    '\n\nStepwise regression model not quantified.\n\n'

    )
}

\n

if (!isTRUE(config::get("skip model"))) {
  knitr::kable(summary$to_visualize)
}

Differential Expression

plots <- targets::tar_read(plot_de_volcano)
for(p in plots) {
  print(p)
}

Sage-Bionetworks/sageseqr documentation built on June 13, 2024, 2:11 p.m.