## ----setup, include = FALSE---------------------------------------------------
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>",
crop = NULL
)
## ----getPackage, eval=FALSE---------------------------------------------------
# ## install MIRit from Bioconductor
# if (!requireNamespace("BiocManager", quietly = TRUE)) {
# install.packages("BiocManager")
# }
# BiocManager::install("MIRit")
## ----getPackageDevel, eval=FALSE----------------------------------------------
# ## install the development version from GitHub
# if (!requireNamespace("BiocManager", quietly = TRUE)) {
# install.packages("BiocManager")
# }
# BiocManager::install("jacopo-ronchi/MIRit")
## ----loadPackage, message=FALSE-----------------------------------------------
## load MIRit
library(MIRit)
## ----example------------------------------------------------------------------
## load count matrix for genes
data(geneCounts, package = "MIRit")
## load count matrix for miRNAs
data(mirnaCounts, package = "MIRit")
## ----geneExpr, echo=FALSE-----------------------------------------------------
## print a table for gene expression matrix
knitr::kable(geneCounts[seq(5), seq(5)], digits = 2, caption = "Gene expression matrix. An expression matrix containing normalized values, with row names according to hgnc.")
## ----mirnaExpr, echo=FALSE----------------------------------------------------
## print a table for miRNA expression matrix
knitr::kable(mirnaCounts[seq(5), seq(5)], digits = 2, caption = "MiRNA expression matrix. An expression matrix containing normalized values, with row names according to miRBase.")
## ----colnames-----------------------------------------------------------------
## print sample names in geneCounts
colnames(geneCounts)
## print sample names in mirnaCounts
colnames(mirnaCounts)
## check if samples in geneCounts are equal to those in mirnaCounts
identical(colnames(geneCounts), colnames(mirnaCounts))
## ----metadata-----------------------------------------------------------------
## create a data.frame containing sample metadata
meta <- data.frame(primary = colnames(mirnaCounts),
mirnaCol = colnames(mirnaCounts),
geneCol = colnames(geneCounts),
disease = c(rep("PTC", 8), rep("NTH", 8)),
patient = c(rep(paste("Sample_", seq(8), sep = ""), 2)))
## ----mirnaObj-----------------------------------------------------------------
## create the MirnaExperiment object
experiment <- MirnaExperiment(mirnaExpr = mirnaCounts,
geneExpr = geneCounts,
samplesMetadata = meta,
pairedSamples = TRUE)
## ----mds, fig.wide=TRUE, fig.cap="MDS plots for miRNAs and genes. Both plots show that the main source of variability is given by the disease condition."----
geneMDS <- plotDimensions(experiment,
assay = "genes",
condition = "disease",
title = "MDS plot for genes")
mirnaMDS <- plotDimensions(experiment,
assay = "microRNA",
condition = "disease",
title = "MDS plot for miRNAs")
ggpubr::ggarrange(geneMDS, mirnaMDS,
nrow = 1, labels = "AUTO", common.legend = TRUE)
## ----model--------------------------------------------------------------------
## design the linear model for both genes and miRNAs
model <- ~ disease + patient
## ----diffexp------------------------------------------------------------------
## perform differential expression for genes
experiment <- performGeneDE(experiment,
group = "disease",
contrast = "PTC-NTH",
design = model,
pCutoff = 0.01)
## perform differential expression for miRNAs
experiment <- performMirnaDE(experiment,
group = "disease",
contrast = "PTC-NTH",
design = model,
pCutoff = 0.01)
## ----accessDe-----------------------------------------------------------------
## access DE results for genes
deGenes <- geneDE(experiment)
## access DE results for miRNAs
deMirnas <- mirnaDE(experiment)
## ----volcano, fig.wide=TRUE, fig.cap="Volcano plots of gene and miRNA differential expression. (A) shows the differentially expressed genes, while (B) displays differentially expressed miRNAs."----
## create a volcano plot for genes
geneVolcano <- plotVolcano(experiment,
assay = "genes",
title = "Gene differential expression")
## create a volcano plot for miRNAs
mirnaVolcano <- plotVolcano(experiment,
assay = "microRNA",
title = "miRNA differential expression")
## plot graphs side by side
ggpubr::ggarrange(geneVolcano, mirnaVolcano,
nrow = 1, labels = "AUTO", common.legend = TRUE)
## ----thyroidBars, fig.wide=FALSE, fig.cap="Differential expression bar plots for different thyroid genes. Differential expression analysis demonstrated how TG, TPO, DIO2 and PAX8 result downregulated in thyroid cancer."----
## create a bar plot for all thyroid features
plotDE(experiment,
features = c("TG", "TPO", "DIO2", "PAX8"),
graph = "barplot", linear = FALSE)
## ----species------------------------------------------------------------------
## list available species for Reactome database
supportedOrganisms("Reactome")
## ----ora----------------------------------------------------------------------
## perform over-representation analysis with GO
ora <- enrichGenes(experiment,
method = "ORA",
database = "GO",
organism = "Homo sapiens")
## ----gsea---------------------------------------------------------------------
## set seed for reproducible results
set.seed(1234)
## perform gene-set enrichment analysis with KEGG
gse <- enrichGenes(experiment,
method = "GSEA",
database = "KEGG",
organism = "Homo sapiens")
## ----gseaTab, echo=FALSE------------------------------------------------------
## display the results of GSEA
knitr::kable(enrichmentResults(gse), digits = 2, caption = "GSEA results. A table listing the affected KEGG pathways in thyroid cancer.")
## ----oraDot, fig.wide=FALSE, fig.cap="ORA results for downregulated genes. The enrichment of downregulated genes through the gene sets provided by GO database."----
## create a dot plot for ORA
enrichmentDotplot(ora$downregulated, title = "Depleted functions")
## ----gsePlot, fig.wide=FALSE, fig.cap="GSEA-style plot for Thyroid hormone synthesis. This type of plot shows the running sum that GSEA uses to determinate the enrichment score for each pathway."----
## create a GSEA plot
gseaPlot(gse, "Thyroid hormone synthesis", rankingMetric = TRUE)
## ----targets, results='hide', eval=FALSE--------------------------------------
# ## retrieve miRNA target genes
# experiment <- getTargets(experiment)
## ----targetsTMP, include=FALSE------------------------------------------------
experiment <- loadExamples()
## ----correlation--------------------------------------------------------------
## perform a correlation analysis
experiment <- mirnaIntegration(experiment, test = "correlation")
## ----correlationResults-------------------------------------------------------
## extract correlation results
integrationResults <- integration(experiment)
## take a look at correlation results
head(integrationResults)
## ----corPlot, fig.wide=TRUE, fig.cap="Correlation between miRNAs and key thyroid genes. These plots suggest that the upregulation of miR-146b-5p and miR-146b-3p may be responsible for the downregulation of DIO2 and PAX8, respectively."----
## plot the correlation between miR-146b-5p and DIO2
cor1 <- plotCorrelation(experiment,
mirna = "hsa-miR-146b-5p",
gene = "DIO2",
condition = "disease")
## plot the correlation between miR-146b-3p and PAX8
cor2 <- plotCorrelation(experiment,
mirna = "hsa-miR-146b-3p",
gene = "PAX8",
condition = "disease")
## plot graphs side by side
ggpubr::ggarrange(cor1, cor2, nrow = 1,
labels = "AUTO", common.legend = TRUE)
## ----association--------------------------------------------------------------
## perform a one-sided inverse association
exp.association <- mirnaIntegration(experiment,
test = "association",
associationMethod = "fisher-midp",
pCutoff = 0.2,
pAdjustment = "none")
## ----fry----------------------------------------------------------------------
## perform the integrative analysis through 'fry' method
exp.fry <- mirnaIntegration(experiment,
test = "fry",
pAdjustment = "none")
## ----sessionInfo, echo=FALSE--------------------------------------------------
sessionInfo()