iSEEde 1.5.0
In this example, we use the ?airway
data set.
library("airway")
data("airway")
This section demonstrates one of many possible workflows for adding annotations to the data set. Those annotations are meant to make it easier for users to identify genes of interest, e.g. by displaying both gene symbols and ENSEMBL gene identifiers as tooltips in the interactive browser.
First, we make a copy of the Ensembl identifiers – currently stored in the
rownames()
– to a column in the rowData()
component.
rowData(airway)[["ENSEMBL"]] <- rownames(airway)
Then, we use the org.Hs.eg.db package to map the
Ensembl identifiers to gene symbols. We store those gene symbols as an
additional column of the rowData()
component.
library("org.Hs.eg.db")
rowData(airway)[["SYMBOL"]] <- mapIds(
org.Hs.eg.db, rownames(airway),
"SYMBOL", "ENSEMBL"
)
Next, we use the uniquifyFeatureNames()
function of the
scuttle package to replace the rownames()
by a
unique identifier that is generated as follows:
library("scuttle")
rownames(airway) <- uniquifyFeatureNames(
ID = rowData(airway)[["ENSEMBL"]],
names = rowData(airway)[["SYMBOL"]]
)
airway
#> class: RangedSummarizedExperiment
#> dim: 63677 8
#> metadata(1): ''
#> assays(1): counts
#> rownames(63677): TSPAN6 TNMD ... APP-DT ENSG00000273493
#> rowData names(12): gene_id gene_name ... ENSEMBL SYMBOL
#> colnames(8): SRR1039508 SRR1039509 ... SRR1039520 SRR1039521
#> colData names(9): SampleName cell ... Sample BioSample
To generate some example results, we first use edgeR::filterByExpr()
to remove
genes whose counts are too low to support a rigorous differential expression
analysis. Then we run a standard Limma-Voom analysis using edgeR::voomLmFit()
,
limma::makeContrasts()
, and limma::eBayes()
; alternatively, we could have
used limma::treat()
instead of limma::eBayes()
.
The linear model includes the dex
and cell
covariates, indicating the
treatment conditions and cell types, respectively. Here, we are interested in
differences between treatments, adjusted by cell type, and define this
comparison as the dextrt - dexuntrt
contrast.
The final differential expression results are fetched using limma::topTable()
.
library("edgeR")
design <- model.matrix(~ 0 + dex + cell, data = colData(airway))
keep <- filterByExpr(airway, design)
fit <- voomLmFit(airway[keep, ], design, plot = FALSE)
contr <- makeContrasts("dextrt - dexuntrt", levels = design)
fit <- contrasts.fit(fit, contr)
fit <- eBayes(fit)
res_limma <- topTable(fit, sort.by = "P", n = Inf)
head(res_limma)
#> gene_id gene_name entrezid gene_biotype gene_seq_start gene_seq_end seq_name
#> CACNB2 ENSG00000165995 CACNB2 NA protein_coding 18429606 18830798 10
#> DUSP1 ENSG00000120129 DUSP1 NA protein_coding 172195093 172198198 5
#> PRSS35 ENSG00000146250 PRSS35 NA protein_coding 84222194 84235423 6
#> MAOA ENSG00000189221 MAOA NA protein_coding 43515467 43606068 X
#> STEAP2 ENSG00000157214 STEAP2 NA protein_coding 89796904 89867451 7
#> SPARCL1 ENSG00000152583 SPARCL1 NA protein_coding 88394487 88452213 4
#> seq_strand seq_coord_system symbol ENSEMBL SYMBOL logFC AveExpr t
#> CACNB2 1 NA CACNB2 ENSG00000165995 CACNB2 3.205606 3.682244 37.68303
#> DUSP1 -1 NA DUSP1 ENSG00000120129 DUSP1 2.864778 6.644455 28.50569
#> PRSS35 1 NA PRSS35 ENSG00000146250 PRSS35 -2.828184 3.224885 -28.10830
#> MAOA 1 NA MAOA ENSG00000189221 MAOA 3.256085 5.950559 27.66135
#> STEAP2 1 NA STEAP2 ENSG00000157214 STEAP2 1.894559 6.790009 27.40834
#> SPARCL1 -1 NA SPARCL1 ENSG00000152583 SPARCL1 4.489009 4.166904 27.34820
#> P.Value adj.P.Val B
#> CACNB2 1.115938e-10 1.881472e-06 14.45518
#> DUSP1 1.148539e-09 4.309046e-06 13.01074
#> PRSS35 1.291089e-09 4.309046e-06 12.44376
#> MAOA 1.475547e-09 4.309046e-06 12.75540
#> STEAP2 1.592916e-09 4.309046e-06 12.71210
#> SPARCL1 1.622327e-09 4.309046e-06 12.33402
Then, we embed this set of differential expression results in the airway
object using the embedContrastResults()
method and we use the function contrastResults()
to display the contrast results stored in the airway
object.
library(iSEEde)
#> Loading required package: iSEE
airway <- embedContrastResults(res_limma, airway,
name = "dextrt - dexuntrt",
class = "limma"
)
contrastResults(airway)
#> DataFrame with 63677 rows and 1 column
#> dextrt - dexuntrt
#> <iSEELimmaResults>
#> TSPAN6 <iSEELimmaResults>
#> TNMD <iSEELimmaResults>
#> DPM1 <iSEELimmaResults>
#> SCYL3 <iSEELimmaResults>
#> FIRRM <iSEELimmaResults>
#> ... ...
#> ENSG00000273489 <iSEELimmaResults>
#> ENSG00000273490 <iSEELimmaResults>
#> ENSG00000273491 <iSEELimmaResults>
#> APP-DT <iSEELimmaResults>
#> ENSG00000273493 <iSEELimmaResults>
contrastResults(airway, "dextrt - dexuntrt")
#> iSEELimmaResults with 63677 rows and 18 columns
#> gene_id gene_name entrezid gene_biotype gene_seq_start gene_seq_end
#> <character> <character> <integer> <character> <integer> <integer>
#> TSPAN6 ENSG00000000003 TSPAN6 NA protein_coding 99883667 99894988
#> TNMD NA NA NA NA NA NA
#> DPM1 ENSG00000000419 DPM1 NA protein_coding 49551404 49575092
#> SCYL3 ENSG00000000457 SCYL3 NA protein_coding 169818772 169863408
#> FIRRM ENSG00000000460 C1orf112 NA protein_coding 169631245 169823221
#> ... ... ... ... ... ... ...
#> ENSG00000273489 NA NA NA NA NA NA
#> ENSG00000273490 NA NA NA NA NA NA
#> ENSG00000273491 NA NA NA NA NA NA
#> APP-DT NA NA NA NA NA NA
#> ENSG00000273493 NA NA NA NA NA NA
#> seq_name seq_strand seq_coord_system symbol ENSEMBL SYMBOL
#> <character> <integer> <integer> <character> <character> <character>
#> TSPAN6 X -1 NA TSPAN6 ENSG00000000003 TSPAN6
#> TNMD NA NA NA NA NA NA
#> DPM1 20 -1 NA DPM1 ENSG00000000419 DPM1
#> SCYL3 1 -1 NA SCYL3 ENSG00000000457 SCYL3
#> FIRRM 1 1 NA C1orf112 ENSG00000000460 FIRRM
#> ... ... ... ... ... ... ...
#> ENSG00000273489 NA NA NA NA NA NA
#> ENSG00000273490 NA NA NA NA NA NA
#> ENSG00000273491 NA NA NA NA NA NA
#> APP-DT NA NA NA NA NA NA
#> ENSG00000273493 NA NA NA NA NA NA
#> logFC AveExpr t P.Value adj.P.Val B
#> <numeric> <numeric> <numeric> <numeric> <numeric> <numeric>
#> TSPAN6 -0.464216 5.02559 -6.589673 0.000136603 0.00171857 0.899306
#> TNMD NA NA NA NA NA NA
#> DPM1 0.125077 4.60191 1.632802 0.139297605 0.25279527 -6.274792
#> SCYL3 -0.042107 3.47269 -0.438929 0.671768770 0.77787936 -7.230811
#> FIRRM -0.228517 1.40857 -0.977961 0.355374953 0.49644090 -6.208897
#> ... ... ... ... ... ... ...
#> ENSG00000273489 NA NA NA NA NA NA
#> ENSG00000273490 NA NA NA NA NA NA
#> ENSG00000273491 NA NA NA NA NA NA
#> APP-DT NA NA NA NA NA NA
#> ENSG00000273493 NA NA NA NA NA NA
In this example, we use iSEE::panelDefaults()
to specify rowData()
fields to
show in the tooltip that is displayed when hovering a data point.
The application is then configured to display the volcano plot and MA plot for the same contrast.
Finally, the configured app is launched.
library(iSEE)
panelDefaults(
TooltipRowData = c("SYMBOL", "ENSEMBL")
)
app <- iSEE(airway, initial = list(
VolcanoPlot(ContrastName = "dextrt - dexuntrt", PanelWidth = 6L),
MAPlot(ContrastName = "dextrt - dexuntrt", PanelWidth = 6L)
))
if (interactive()) {
shiny::runApp(app)
}