differential_outcomes
Brian M. Schilder
2023-03-23
library(data.table)
library(ggplot2)
library(ggridges)Can we say, that this phenotype, when it acts via this celltype, has a different clinical course?
Import data
d <- MultiEWCE::load_example_results("Descartes_All_Results_extras.symptoms.full_join.rds")## Registered S3 method overwritten by 'ggtree':
## method from
## fortify.igraph ggnetwork{
d <- HPOExplorer::add_disease(d)
d$Phenotype <-
HPOExplorer::harmonise_phenotypes(d$HPO_ID, keep_order = TRUE)
d[,DiseaseDB:=sapply(DatabaseID,function(x){strsplit(x,":")[[1]][1]})]
d <- HPOExplorer::add_death(d, agg_by = "DatabaseID")
d <- HPOExplorer::add_ndisease(d)
d <- HPOExplorer::add_ancestor(d)
d <- HPOExplorer::add_pheno_frequency(d)
d <- HPOExplorer::add_ont_lvl(d)
d <- HPOExplorer::add_onset(d, agg_by = c("DatabaseID","HPO_ID"))
d <- HPOExplorer::add_severity(d)
d <- HPOExplorer::add_tier(d)
} ## Annotating phenos with Disease## Importing existing file: ... phenotype.hpoa## Translating all phenotypes to names.## + Returning a vector of phenotypes (same order as input).## Annotating phenos with AgeOfDeath## Annotating phenos with n_diseases## Importing existing file: ... phenotype_to_genes.txt## Importing existing file: ... genes_to_phenotype.txt## Importing existing file: ... phenotype.hpoa## Adding level-3 ancestor to each HPO ID.## Annotating phenotype frequencies.## Getting absolute ontology level for 6,170 HPO IDs.## Annotating phenos with Onset.## Annotating phenos with Modifiers## Annotating phenos with Tiers.Find phenotypes associated w/ 20+ diseases & have metadata.
Then test whether the EWCE enrichment p-values (conditional on celltype) significantly affect whether the phenotype affects clinical outcomes (e.g. age of death).
The formula for the linear model is as follows: >
outcome ~ EWC_pvalue * celltype
Run tests
save_path <- here::here("reports","differential_outcomes.csv.gz")
if(file.exists(save_path)){
res_dt <- data.table::fread(save_path)
} else {
vars <- c("AgeOfDeath_score_min","AgeOfDeath_score_max","AgeOfDeath_score_mean",
"Onset_score_min","Onset_score_max","Onset_score_mean",
"pheno_freq_min","pheno_freq_max","pheno_freq_mean",
"Severity_score","tier_merge"
# "n_diseases","symptoms.pval"
)
set.seed(2023)
options(future.globals.maxSize = 8000*1024^2)
#### Iterate over each variable ####
res_dt <- lapply(stats::setNames(vars,vars),
function(v){
message("Testing: ",v)
d1 <- d[!is.na(get(v))]
d1[,n_diseases:=length(unique(DatabaseID)), by="HPO_ID"]
d1 <- d1[n_diseases>20]
ids <- unique(d1$HPO_ID)
#### Iterate over each HPO ID ####
furrr::future_map(.x = ids,
.progress = TRUE,
.f = function(hpo_id){
d2 <- d1[HPO_ID==hpo_id]
if(length(unique(d2[[v]]))<2 |
length(unique(d2[["CellType"]]))<2) {
return(NULL)
}
res <- lm(data = d2,
formula = paste(v,"~ p*CellType"))
data.table::data.table(broom::tidy(res))[,HPO_ID:=hpo_id]
}) |> data.table::rbindlist()
}) |> data.table::rbindlist(fill = TRUE, use.names = TRUE, idcol = "variable")
res_dt[,test_id:=paste(variable,HPO_ID,sep=".")]
data.table::fwrite(res_dt,save_path)
} Plot
Intercepts
Plotting the q-value of the intercept of each test basically tells us the number of phenotypes in which celltype was a strong modifier of each clinical course.
plot_variables <- function(res_dt,
type="barplot",
scales="free"){
dat <- res_dt[term=="(Intercept)"][,q.value:=stats::p.adjust(p = p.value,method = "fdr")]
if(type=="violin"){
ggplot(dat, aes(x=-log10(q.value), y=variable, fill=variable)) +
geom_violin(show.legend = FALSE, color=alpha("black",.3), na.rm = TRUE) +
geom_point(size=1, alpha=.1, show.legend = FALSE, na.rm = TRUE) +
geom_vline(xintercept = -log10(0.05), linetype="dashed", alpha=0.5,color="blue") +
annotate(x=-log10(0.05), y=-1.1,
label=expression("| q-value<0.05" %->% ""),
geom = "text", hjust=0,
size=3,color=alpha("black",.5)) +
coord_cartesian(ylim = c(0, 8), clip = "off") +
theme_bw()
} else {
ggplot(dat, aes(x=-log10(q.value),
fill=variable)) +
geom_histogram(bins = 75, show.legend = FALSE) +
geom_vline(xintercept = -log10(0.05), linetype="dashed", alpha=0.5, color="blue") +
facet_wrap(facets = variable~., scales = scales) +
scale_fill_viridis_d(option = "mako", begin = .25, end = .8) +
theme_bw() +
theme(strip.background = element_rect(fill="white"))
}
}Violin plot
Summarised as a violin plot:
vp <- plot_variables(res_dt = res_dt,
type = "violin")
vp## Warning: Groups with fewer than two data points have been dropped.## Warning in is.na(x): is.na() applied to non-(list or vector) of type
## 'expression'
Let’s cut off super significant values so that we can see the distributions better:
vp + xlim(0,110)## Warning: Groups with fewer than two data points have been dropped.## Warning in is.na(x): is.na() applied to non-(list or vector) of type
## 'expression'
Histograms
We can see the distributions a bit better with histograms.
The vast majority of tests came back as significant. Meaning celltypes very have an impact on all clinical course variables.
plot_variables(res_dt = res_dt,
type = "histogram")
Celltypes
What proportion of celltypes are each phenotype enriched for?
Finally, let’s plot the p-values of each celltype within the model. This shows that celltypes tend not to be indiscriminately associated with most phenotype.
In other words, celltypes are relatively specific to certain symptoms/phenotypes.
plot_celltypes <- function(res_dt,
v=NULL,
x_var="-log10(p.value)",
scales="free_y"){
dat <- res_dt[grepl("^CellType",term)][,term:=gsub("^CellType","",term)]
if(!is.null(v)){
dat <- dat[variable %in% v]
}
gp <- ggplot(dat,
aes_string(x=x_var,
fill="term",
color="term")) +
geom_density(adjust=.9, position = "stack", color=alpha("cyan",.75), linewidth=.01, na.rm = TRUE) +
scale_fill_manual(values = pals::ocean.phase(n = length(unique(dat$term)))) +
theme_bw() +
theme(legend.key.size = unit(1, 'lines'), #change legend key size
legend.key.height = unit(.5, 'lines'), #change legend key height
legend.key.width = unit(.5, 'cm'), #change legend key width
legend.title = element_text(size=6), #change legend title font size
legend.text = element_text(size=5),
legend.background = element_blank(),
strip.background = element_rect(fill="white"),
legend.direction = "horizontal",
legend.position="bottom"
) +
guides(fill=guide_legend(ncol=5)) +
facet_wrap(facets = variable~.,
scales = scales)
if(x_var=="-log10(p.value)"){
gp <- gp +
geom_vline(xintercept = -log10(0.05), linetype="dashed", alpha=0.5, color="blue")
}
return(gp)
}plot_celltypes(res_dt = res_dt)
What is the magnitude of effect of celltype identity on each clinical course?
plot_celltypes(res_dt = res_dt,
x_var = "log10(abs(estimate))",
scales="free")## Warning: Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.
## Groups with fewer than two data points have been dropped.## Warning: Removed 10 rows containing missing values (`position_stack()`).
Distributions
Now let’s highlight some specific examples where celltype can cause a differential clinical course.
sig_ids <- res_dt[term=="(Intercept)"][,q.value:=stats::p.adjust(p = p.value,method = "fdr")][q.value<0.05]$test_id
sig_dt <- res_dt[test_id %in% sig_ids][grepl("^CellType",term)][,CellType:=gsub("^CellType","",term)]
message(formatC(length(unique(sig_dt$HPO_ID)),big.mark = ","),
" phenotypes have a differential clinical course dependent on the celltype:")## 745 phenotypes have a differential clinical course dependent on the celltype:MultiEWCE::create_dt(head(sig_dt))## Loading required namespace: DTid.vars <- c("HPO_ID","Phenotype","CellType")
d_melt <-
(d[,unique(c(id.vars,unique(sig_dt$variable)) ),with=FALSE] |>
unique() |>
data.table::melt.data.table(
id.vars = id.vars,
measure.vars = unique(sig_dt$variable))
)[!is.na(value)]
dat <- data.table::merge.data.table(sig_dt,
d_melt,
all.x = TRUE,
all.y = TRUE,
by=c("HPO_ID","variable","CellType")) |> unique()
dat[is.na(p.value)]$p.value <- 1
dat[,modifying_celltype:=factor(p.value<.05, levels = c(TRUE, FALSE),
ordered = TRUE)]
dat[,value_mean:=mean(value),by="variable"]
dat[,direction:=ifelse(estimate>0,"+","-"),]
dat[,directional_diff:=ifelse(estimate>0,
value_mean+value,
value_mean-value)]
dat[,directional_diff:=ifelse(directional_diff<0,0,directional_diff)]
dat[,fdr:=stats::p.adjust(p.value,method = "fdr")]ggplot(dat, aes(x=modifying_celltype,
y= directional_diff ,
fill=modifying_celltype)) +
geom_violin(na.rm = TRUE) +
geom_boxplot(fill="transparent", color="turquoise",
alpha=.2,
outlier.alpha = .25, na.rm = TRUE) +
facet_wrap(facets = .~variable, scales = "free") +
scale_fill_viridis_d(option = "mako", begin = .25, direction = -1) +
labs(x="Directional difference", y=NULL) +
theme_bw() +
theme(axis.text.x = element_blank(),
strip.background = element_blank()) 
Let’s plot the distribution each metadata attribute overall.
ggplot(dat,
aes(x=value,
fill=variable,
y = -0.5)) +
geom_boxplot(position = position_dodge2(preserve = "single"),
show.legend = FALSE, na.rm = TRUE) +
geom_density(aes(x = value, fill = variable),
show.legend = TRUE, na.rm = TRUE, inherit.aes = FALSE,
alpha=.75, adjust=.5) +
facet_wrap(facets = variable~.,
scales="free_x",
ncol = 3) +
stat_boxplot(geom = "vline", aes(xintercept = ..xlower..),
width=position_dodge2(width = 1)) +
stat_boxplot(geom = "vline", aes(xintercept = ..xmiddle..),
width=position_dodge2(width = 1)) +
stat_boxplot(geom = "vline", aes(xintercept = ..xupper..),
width=position_dodge2(width = 1)) +
labs(y="Density") +
scale_y_discrete(drop=FALSE) +
scale_fill_viridis_d(option = "mako", alpha=.7,
begin = .25, end=.75,
drop=FALSE) +
theme_bw() +
theme(strip.background = element_blank())## Warning: The dot-dot notation (`..xlower..`) was deprecated in ggplot2 3.4.0.
## ℹ Please use `after_stat(xlower)` instead.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.## Warning: Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Caused by error in `$<-.data.frame`:
## ! replacement has 5 rows, data has 1
For each test result, what is the absolute difference in the quantitative clinical course between the modifying celltypes and the non-modifiying celltypes?
## Find instances where we have a pvalue but no estimate
p_not_e <- nrow(dat[(!is.na(p.value)) & is.na(estimate)])
# test_id = clinical course + phenotype
dat_diff <- dat[,list(
n_modifying=length(unique(CellType[modifying_celltype==TRUE])),
n_nonmodifying=sum(modifying_celltype==FALSE),
value_modifying=mean(value[modifying_celltype==TRUE]),
value_nonmodifying=mean(value[modifying_celltype==FALSE]),
estimate_modifying=mean(estimate[modifying_celltype==TRUE]),
estimate_nonmodifying=mean(estimate[modifying_celltype==FALSE])
),
by=c("variable","HPO_ID","test_id")][n_modifying>0 & n_nonmodifying>0][,value_diff:=abs(value_modifying-value_nonmodifying)][,estimate_diff:=abs(estimate_modifying-estimate_nonmodifying)] ggplot(dat_diff,
aes(x=value_diff,
fill=variable,
y = "Boxplot")) +
geom_density(aes(x = value_diff, fill = variable),
show.legend = TRUE, na.rm = TRUE, inherit.aes = FALSE,
alpha=.75, adjust=.5, color="white") +
geom_boxplot(position = position_dodge2(preserve = "single"),
show.legend = FALSE) +
facet_wrap(facets = variable~.,
scales="free_x",
ncol = 3) +
stat_boxplot(geom = "vline", aes(xintercept = ..xlower..),
width=position_dodge2(width = 1)) +
stat_boxplot(geom = "vline", aes(xintercept = ..xmiddle..),
width=position_dodge2(width = 1)) +
stat_boxplot(geom = "vline", aes(xintercept = ..xupper..),
width=position_dodge2(width = 1)) +
labs(y="Density") +
scale_y_discrete(drop=FALSE) +
scale_fill_viridis_d(option = "mako", alpha=.7,
begin = .25, end=.75,
drop=FALSE) +
theme_bw() +
theme(strip.background = element_blank())## Warning: Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Computation failed in `stat_boxplot()`
## Caused by error in `$<-.data.frame`:
## ! replacement has 5 rows, data has 1
Now let’s try to project these difference back onto the full data distribution so we can see how different that looks on the scale.
We will want to differentiate between +/- values to get a sense of how much higher or how much lower celltypes can make a given clinical course.
If the estimate is negative, subtract the difference between the population mean and the values from the population mean:
When the estimates are positive > lower = mean + abs(mean-values)
When the estimates are negative > lower = mean - abs(mean-values)
lvls <- c("1. Population distribution",
"2. All results",
"3. Non-modifying cell types\n(p-value>0.05)",
"4. Modifying cell types\n(p-value<0.05)",
"5. Modifying cell types\n(p-value<0.005)",
"6. Modifying cell types\n(FDR<0.05)",
"7. Directional difference")
ggplot(data = dat[modifying_celltype==TRUE,],
aes(x=directional_diff,
y = tail(lvls,1))) +
ggridges::stat_density_ridges(
aes(x = value,
y = lvls[1],
color = lvls[1],
fill = 0.5 - abs(0.5 - stat(ecdf))) ,
data=d_melt,
# color="purple",
geom = "density_ridges_gradient", calc_ecdf = TRUE,
quantiles = 4, quantile_lines = TRUE,
quantile_fun=function(x,...)mean(x),
alpha = 0.5) +
ggridges::stat_density_ridges(
aes(x = directional_diff,
y = lvls[2],
color=lvls[2],
fill = 0.5 - abs(0.5 - stat(ecdf))),
data = dat,
# color="red",
geom = "density_ridges_gradient", calc_ecdf = TRUE,
quantiles = 4, quantile_lines = TRUE,
quantile_fun=function(x,...)mean(x),
alpha = 0.5) +
ggridges::stat_density_ridges(
aes(x = directional_diff,
y = lvls[3],
color=lvls[3],
fill = 0.5 - abs(0.5 - stat(ecdf))),
data = dat[modifying_celltype==FALSE,],
# color="red",
geom = "density_ridges_gradient", calc_ecdf = TRUE,
quantiles = 4, quantile_lines = TRUE,
quantile_fun=function(x,...)mean(x),
alpha = 0.5) +
####
ggridges::stat_density_ridges(
aes(x = directional_diff,
y = lvls[4],
color = lvls[4],
fill = 0.5 - abs(0.5 - stat(ecdf))),
data = dat[modifying_celltype==TRUE & p.value<0.05,],
geom = "density_ridges_gradient", calc_ecdf = TRUE,
quantiles = 4, quantile_lines = TRUE,
quantile_fun=function(x,...)mean(x),
jittered_points = FALSE,
position = position_points_jitter(width = 0.05,
yoffset = -.1,
height = 0),
point_shape = '|', point_size = 3,alpha = 0.5) +
ggridges::stat_density_ridges(
aes(x = directional_diff,
y = lvls[5],
color=lvls[5],
fill = 0.5 - abs(0.5 - stat(ecdf))),
data = dat[modifying_celltype==TRUE & p.value<0.005,],
# color="red",
geom = "density_ridges_gradient", calc_ecdf = TRUE,
quantiles = 4, quantile_lines = TRUE,
quantile_fun=function(x,...)mean(x),
alpha = 0.5) +
ggridges::stat_density_ridges(
aes(x = directional_diff,
y = lvls[6],
color=lvls[6],
fill = 0.5 - abs(0.5 - stat(ecdf))),
data = dat[modifying_celltype==TRUE & fdr<0.05,],
# color="red",
geom = "density_ridges_gradient", calc_ecdf = TRUE,
quantiles = 4, quantile_lines = TRUE,
quantile_fun=function(x,...)mean(x),
alpha = 0.5) +
scale_fill_viridis_c(option = "mako", alpha=.8) +
scale_color_viridis_d(option = "plasma", alpha=.8, end=.8) +
facet_wrap(facets = variable~.,
scales="free",
ncol = 3) +
labs(y="Density", x="Directional difference") +
theme_bw() +
theme(strip.background = element_blank()) +
xlim(limits = c(0, NA))
# geom_boxplot(aes(color=direction),
# data = dat[modifying_celltype==TRUE,],
# position = position_identity(),
# show.legend = TRUE) +
# stat_boxplot(geom = "vline", aes(xintercept = ..xlower..),
# width=position_dodge2(width = 1),
# show.legend = FALSE) +
# stat_boxplot(geom = "vline", aes(xintercept = ..xmiddle..),
# width=position_dodge2(width = 1),
# show.legend = FALSE) +
# stat_boxplot(geom = "vline", aes(xintercept = ..xupper..),
# width=position_dodge2(width = 1),
# show.legend = FALSE) Session info
sessionInfo()## R version 4.2.1 (2022-06-23)
## Platform: x86_64-apple-darwin17.0 (64-bit)
## Running under: macOS Big Sur ... 10.16
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/4.2/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] ggridges_0.5.4 ggplot2_3.4.1 data.table_1.14.8
##
## loaded via a namespace (and not attached):
## [1] backports_1.4.1 AnnotationHub_3.6.0
## [3] BiocFileCache_2.6.1 plyr_1.8.8
## [5] lazyeval_0.2.2 orthogene_1.4.1
## [7] ewceData_1.7.1 crosstalk_1.2.0
## [9] GenomeInfoDb_1.34.9 ggnetwork_0.5.12
## [11] digest_0.6.31 yulab.utils_0.0.6
## [13] htmltools_0.5.4 RNOmni_1.0.1
## [15] fansi_1.0.4 magrittr_2.0.3
## [17] memoise_2.0.1 ontologyPlot_1.6
## [19] limma_3.54.2 Biostrings_2.66.0
## [21] matrixStats_0.63.0 GeneOverlap_1.34.0
## [23] R.utils_2.12.2 colorspace_2.1-0
## [25] blob_1.2.4 rappdirs_0.3.3
## [27] xfun_0.37 dplyr_1.1.1
## [29] crayon_1.5.2 RCurl_1.98-1.10
## [31] jsonlite_1.8.4 graph_1.76.0
## [33] ape_5.7-1 glue_1.6.2
## [35] pals_1.7 gtable_0.3.3
## [37] zlibbioc_1.44.0 XVector_0.38.0
## [39] HGNChelper_0.8.1 DelayedArray_0.24.0
## [41] car_3.1-1 Rgraphviz_2.42.0
## [43] SingleCellExperiment_1.20.1 maps_3.4.1
## [45] BiocGenerics_0.44.0 abind_1.4-5
## [47] scales_1.2.1 DBI_1.1.3
## [49] rstatix_0.7.2 Rcpp_1.0.10
## [51] viridisLite_0.4.1 xtable_1.8-4
## [53] gridGraphics_0.5-1 tidytree_0.4.2
## [55] mapproj_1.2.11 bit_4.0.5
## [57] DT_0.27 stats4_4.2.1
## [59] htmlwidgets_1.6.2 httr_1.4.5
## [61] ontologyIndex_2.10 gplots_3.1.3
## [63] HPOExplorer_0.99.7 RColorBrewer_1.1-3
## [65] ellipsis_0.3.2 farver_2.1.1
## [67] R.methodsS3_1.8.2 pkgconfig_2.0.3
## [69] sass_0.4.5 dbplyr_2.3.2
## [71] here_1.0.1 utf8_1.2.3
## [73] labeling_0.4.2 reshape2_1.4.4
## [75] ggplotify_0.1.0 tidyselect_1.2.0
## [77] rlang_1.1.0 later_1.3.0
## [79] AnnotationDbi_1.60.2 munsell_0.5.0
## [81] BiocVersion_3.16.0 tools_4.2.1
## [83] cachem_1.0.7 MultiEWCE_0.1.4
## [85] cli_3.6.0 generics_0.1.3
## [87] RSQLite_2.3.0 ExperimentHub_2.6.0
## [89] statnet.common_4.8.0 broom_1.0.4
## [91] evaluate_0.20 stringr_1.5.0
## [93] fastmap_1.1.1 yaml_2.3.7
## [95] ggtree_3.6.2 babelgene_22.9
## [97] knitr_1.42 bit64_4.0.5
## [99] caTools_1.18.2 purrr_1.0.1
## [101] KEGGREST_1.38.0 gprofiler2_0.2.1
## [103] nlme_3.1-162 mime_0.12
## [105] R.oo_1.25.0 grr_0.9.5
## [107] aplot_0.1.10 compiler_4.2.1
## [109] rstudioapi_0.14 plotly_4.10.1
## [111] filelock_1.0.2 curl_5.0.0
## [113] png_0.1-8 interactiveDisplayBase_1.36.0
## [115] ggsignif_0.6.4 paintmap_1.0
## [117] treeio_1.23.1 tibble_3.2.1
## [119] EWCE_1.7.2 bslib_0.4.2
## [121] homologene_1.4.68.19.3.27 stringi_1.7.12
## [123] highr_0.10 lattice_0.20-45
## [125] Matrix_1.5-3 vctrs_0.6.1
## [127] pillar_1.9.0 lifecycle_1.0.3
## [129] BiocManager_1.30.20 jquerylib_0.1.4
## [131] bitops_1.0-7 httpuv_1.6.9
## [133] patchwork_1.1.2 GenomicRanges_1.50.2
## [135] R6_2.5.1 promises_1.2.0.1
## [137] network_1.18.1 KernSmooth_2.23-20
## [139] IRanges_2.32.0 dichromat_2.0-0.1
## [141] gtools_3.9.4 SummarizedExperiment_1.28.0
## [143] rprojroot_2.0.3 withr_2.5.0
## [145] S4Vectors_0.36.2 GenomeInfoDbData_1.2.9
## [147] parallel_4.2.1 grid_4.2.1
## [149] ggfun_0.0.9 tidyr_1.3.0
## [151] coda_0.19-4 rmarkdown_2.20.1
## [153] MatrixGenerics_1.10.0 carData_3.0-5
## [155] ggpubr_0.6.0 piggyback_0.1.4
## [157] Biobase_2.58.0 shiny_1.7.4