if(!requireNamespace('BiocManager', quietly = TRUE))
install.packages('BiocManager')
BiocManager::install("BioNERO")
# Load package after installation
library(BioNERO)
set.seed(12) # for reproducibility
Comparing different coexpression networks can reveal relevant biological
patterns. For instance, seeking for consensus modules can identify
coexpression modules that occur in all data sets regardless of natural
variation and, hence, are core players of the studied phenotype.
Additionally, module preservation within and across species can reveal
patterns of conservation and divergence between transcriptomes.
In this vignette, we will explore consensus modules and module preservation
analyses with BioNERO
. Although they seem similar, their goals are opposite:
while consensus modules identification focuses on the commonalities, module
preservation focuses on the divergences.
We will use RNA-seq data of maize (Zea mays) and rice (Oryza sativa) obtained from Shin et al. (2020).
data(zma.se)
zma.se
## class: SummarizedExperiment
## dim: 10802 28
## metadata(0):
## assays(1): ''
## rownames(10802): ZeamMp030 ZeamMp044 ... Zm00001d054106 Zm00001d054107
## rowData names(0):
## colnames(28): SRX339756 SRX339757 ... SRX2792103 SRX2792104
## colData names(1): Tissue
data(osa.se)
osa.se
## class: SummarizedExperiment
## dim: 7647 27
## metadata(0):
## assays(1): ''
## rownames(7647): Os01g0100700 Os01g0100900 ... Os12g0641400 Os12g0641500
## rowData names(0):
## colnames(27): SRX831140 SRX831141 ... SRX263041 SRX1544234
## colData names(1): Tissue
All BioNERO
’s functions for consensus modules and module preservation
analyses require the expression data to be in a list. Each element of the
list can be a SummarizedExperiment object (recommended) or an expression data
frame with genes in row names and samples in column names.
The most common objective in consensus modules identification is to find core modules across different tissues or treatments for the same species. For instance, one can infer GCNs for different types of cancer in human tissues (say prostate and liver) and identify modules that occur in all sets, which are likely core components of cancer biology. Likewise, one can also identify consensus modules across samples from different geographical origins to find modules that are not affected by population structure or kinship.
Here, we will subset 22 random samples from the maize data twice and find consensus modules between the two sets.
# Preprocess data and keep top 2000 genes with highest variances
filt_zma <- exp_preprocess(zma.se, variance_filter = TRUE, n = 2000)
## Number of removed samples: 1
# Create different subsets by resampling data
zma_set1 <- filt_zma[, sample(colnames(filt_zma), size=22, replace=FALSE)]
zma_set2 <- filt_zma[, sample(colnames(filt_zma), size=22, replace=FALSE)]
colnames(zma_set1)
## [1] "SRX3804716" "SRX2792102" "SRX2527287" "SRX2792108" "SRX3804723"
## [6] "SRX2792107" "SRX2792103" "SRX2792104" "SRX3804715" "SRX339808"
## [11] "SRX339758" "SRX339756" "SRX339809" "SRX2792105" "SRX339757"
## [16] "SRX2641029" "SRX339762" "SRX2792110" "SRX339807" "SRX3804718"
## [21] "ERX2154032" "SRX339764"
colnames(zma_set2)
## [1] "SRX2792111" "SRX339756" "SRX3804715" "SRX2792108" "SRX2792103"
## [6] "SRX339762" "SRX339807" "ERX2154030" "SRX2792105" "SRX2792107"
## [11] "SRX2527287" "ERX2154032" "SRX2792104" "SRX2527288" "SRX339809"
## [16] "SRX3804718" "SRX3804716" "SRX2792109" "SRX3804723" "SRX2792102"
## [21] "SRX339808" "SRX339758"
# Create list
zma_list <- list(set1 = zma_set1, set2 = zma_set2)
length(zma_list)
## [1] 2
As described in the first vignette, before inferring the GCNs, we need to
identify the optimal \(\beta\) power that makes the network closer to a scale-free
topology. We can do that with consensus_SFT_fit()
.
cons_sft <- consensus_SFT_fit(zma_list, setLabels = c("Maize 1", "Maize 2"),
cor_method = "pearson")
## Power SFT.R.sq slope truncated.R.sq mean.k. median.k. max.k.
## 1 5 0.658 -0.650 0.656 107.00 90.50 304.0
## 2 6 0.761 -0.706 0.774 79.70 62.90 244.0
## 3 7 0.833 -0.757 0.863 61.30 45.10 200.0
## 4 8 0.825 -0.830 0.868 48.30 32.80 169.0
## 5 9 0.822 -0.916 0.887 38.70 24.50 145.0
## 6 10 0.838 -0.975 0.910 31.50 18.90 126.0
## 7 11 0.824 -1.040 0.914 26.00 15.00 111.0
## 8 12 0.837 -1.090 0.932 21.70 12.20 97.6
## 9 13 0.851 -1.130 0.944 18.30 10.10 86.6
## 10 14 0.842 -1.190 0.941 15.50 8.42 77.3
## 11 15 0.844 -1.230 0.947 13.30 7.17 69.4
## 12 16 0.842 -1.270 0.951 11.50 6.13 62.5
## 13 17 0.855 -1.300 0.962 9.99 5.24 56.5
## 14 18 0.864 -1.320 0.967 8.73 4.54 51.3
## 15 19 0.869 -1.330 0.970 7.67 3.92 46.7
## 16 20 0.875 -1.340 0.977 6.78 3.38 42.7
## Power SFT.R.sq slope truncated.R.sq mean.k. median.k. max.k.
## 1 5 0.592 -0.621 0.574 94.70 82.40 265.0
## 2 6 0.702 -0.666 0.710 69.90 56.80 209.0
## 3 7 0.757 -0.713 0.771 53.20 40.30 167.0
## 4 8 0.829 -0.729 0.853 41.50 29.50 136.0
## 5 9 0.855 -0.782 0.884 32.90 22.10 113.0
## 6 10 0.868 -0.826 0.911 26.60 17.20 96.6
## 7 11 0.873 -0.868 0.931 21.80 13.70 83.2
## 8 12 0.837 -0.943 0.914 18.10 11.30 73.0
## 9 13 0.826 -1.000 0.917 15.20 9.30 64.7
## 10 14 0.815 -1.070 0.917 12.80 7.73 57.6
## 11 15 0.810 -1.120 0.926 11.00 6.53 51.5
## 12 16 0.823 -1.150 0.943 9.42 5.45 46.3
## 13 17 0.815 -1.210 0.942 8.16 4.57 41.8
## 14 18 0.827 -1.240 0.954 7.11 3.92 37.8
## 15 19 0.833 -1.260 0.964 6.23 3.42 34.3
## 16 20 0.831 -1.280 0.965 5.49 2.99 31.2
This function returns a list with the optimal powers and a summary plot,
exactly as SFT_fit()
does.
powers <- cons_sft$power
powers
## set1 set2
## 7 8
cons_sft$plot
Now, we can infer GCNs and identify consensus modules across data sets.
consensus <- consensus_modules(zma_list, power = powers, cor_method = "pearson")
## ..connectivity..
## ..matrix multiplication (system BLAS)..
## ..normalization..
## ..done.
## ..connectivity..
## ..matrix multiplication (system BLAS)..
## ..normalization..
## ..done.
## ..done.
## multiSetMEs: Calculating module MEs.
## Working on set 1 ...
## Working on set 2 ...
names(consensus)
## [1] "consModules" "consMEs" "exprSize"
## [4] "sampleInfo" "genes_cmodules" "dendro_plot_objects"
head(consensus$genes_cmodules)
## Genes Cons_modules
## 1 ZeamMp030 royalblue
## 2 ZeamMp044 darkred
## 3 ZeamMp092 darkred
## 4 ZeamMp108 saddlebrown
## 5 ZeamMp116 darkred
## 6 ZeamMp158 saddlebrown
Finally, we can correlate consensus module eigengenes to sample metadata (here, plant tissues).1 NOTE: Blank grey cells in the heatmap represent correlation values that have opposite sign in the expression sets. For each correlation pair, consensus correlations are calculated by selecting the minimum value across matrices of consensus module-trait correlations.
consensus_trait <- consensus_trait_cor(consensus, cor_method = "pearson")
head(consensus_trait)
## ME trait cor pvalue
## 1 MEblack endosperm -0.297224435 0.1815879
## 2 MEblack pollen 0.284525339 0.2021641
## 3 MEblack whole_seedling 0.007301901 0.9746086
## 4 MEblue endosperm -0.164695363 0.4687656
## 5 MEblue pollen 0.295979347 0.1835410
## 6 MEblue whole_seedling -0.016437184 0.9428769
Users can also transpose the heatmap as in module_trait_cor()
.
Module preservation is often used to study patterns of evolutionary conservation and divergence across transcriptomes, an approach named phylotranscriptomics. This way, one can investigate how evolution shaped the expression profiles for particular gene families across taxa.
To calculate module preservation statistics, gene IDs must be shared by
the expression sets. For intraspecies comparisons, this is an easy task,
as gene IDs are the same. However, for interspecies comparisons, users need
to identify orthogroups between the different species and collapse the
gene-level expression values to orthogroup-level expression values.
This way, all expression sets will have common row names. We recommend
identifying orthogroups with OrthoFinder (Emms and Kelly 2015), as it is simple
to use and widely used.2 PRO TIP: If you identify orthogroups with OrthoFinder,
BioNERO
has a helper function named parse_orthofinder()
that parses
the Orthogroups.tsv file generated by OrthoFinder into a suitable data frame
for module preservation analysis. See ?parse_orthofinder
for more details. Here, we will compare maize and rice expression
profiles. The orthogroups between these species were downloaded from the
PLAZA 4.0 Monocots database (Van Bel et al. 2018).
data(og.zma.osa)
head(og.zma.osa)
## Family Species Gene
## 1548 ORTHO04M000001 osa Os01g0100700
## 1549 ORTHO04M000001 osa Os01g0100900
## 4824 ORTHO04M000001 zma Zm00001d009743
## 4854 ORTHO04M000001 zma Zm00001d020834
## 4874 ORTHO04M000001 zma Zm00001d026672
## 4921 ORTHO04M000001 zma Zm00001d039873
As you can see, the orthogroup object for BioNERO
must be a data frame with
orthogroups, species IDs and gene IDs, respectively. Let’s collapse gene-level
expression to orthogroup-level with exp_genes2orthogroups()
. By default, if
there is more than one gene in the same orthogroup for a given species, their
expression levels are summarized to the median. Users can also summarize
to the mean.
# Store SummarizedExperiment objects in a list
zma_osa_list <- list(osa = osa.se, zma = zma.se)
# Collapse gene-level expression to orthogroup-level
ortho_exp <- exp_genes2orthogroups(zma_osa_list, og.zma.osa, summarize = "mean")
# Inspect new expression data
ortho_exp$osa[1:5, 1:5]
## SRX831140 SRX831141 SRX831137 SRX831138 SRX831134
## ORTHO04M000001 6.909420 7.258330 94.20870 92.85195 123.01060
## ORTHO04M000002 9.203498 8.709974 66.45512 44.97913 33.86936
## ORTHO04M000003 9.417930 9.444861 42.57513 66.02237 55.37741
## ORTHO04M000004 9.019436 8.920091 96.22074 62.56506 109.32262
## ORTHO04M000005 40.845040 41.844234 52.33474 31.31474 22.42236
ortho_exp$zma[1:5, 1:5]
## SRX339756 SRX339757 SRX339758 SRX339762 SRX339763
## ORTHO04M000001 26.02510 15.07917 14.91571 13.82989 8.080476
## ORTHO04M000002 19.28281 13.94254 13.57854 12.96032 13.522525
## ORTHO04M000003 45.17294 48.63796 54.22404 42.12135 10.779117
## ORTHO04M000004 28.05475 38.53734 39.48070 27.13272 2.978207
## ORTHO04M000005 67.58868 34.87009 21.46280 12.79565 7.452068
Now, we will preprocess both expression sets and keep only the top 1000 orthogroups with the highest variances for demonstration purposes.
# Preprocess data and keep top 1000 genes with highest variances
ortho_exp <- lapply(ortho_exp, exp_preprocess, variance_filter=TRUE, n=1000)
## Number of removed samples: 2
## Number of removed samples: 2
# Check orthogroup number
sapply(ortho_exp, nrow)
## osa zma
## 1000 1000
Now that row names are comparable, we can infer GCNs for each set. We will do that iteratively with lapply.
# Calculate SFT power
power_ortho <- lapply(ortho_exp, SFT_fit, cor_method="pearson")
## Power SFT.R.sq slope truncated.R.sq mean.k. median.k. max.k.
## 1 3 0.00503 -0.174 0.792 196.00 196.00 277.0
## 2 4 0.10300 -0.767 0.874 133.00 131.00 210.0
## 3 5 0.15600 -0.861 0.903 94.40 92.90 165.0
## 4 6 0.18900 -0.873 0.929 69.90 68.30 132.0
## 5 7 0.22000 -0.936 0.915 53.40 51.70 109.0
## 6 8 0.32100 -1.070 0.918 41.90 40.10 90.9
## 7 9 0.35400 -1.030 0.903 33.60 31.70 77.0
## 8 10 0.39800 -0.991 0.909 27.50 25.50 66.0
## 9 11 0.46900 -0.997 0.914 22.80 20.90 57.2
## 10 12 0.51800 -0.969 0.917 19.20 17.30 50.0
## 11 13 0.58400 -0.944 0.935 16.40 14.50 44.1
## 12 14 0.66000 -0.913 0.969 14.20 12.50 39.3
## 13 15 0.70900 -0.906 0.977 12.30 10.60 35.5
## 14 16 0.76000 -0.920 0.983 10.80 9.14 32.2
## 15 17 0.77900 -0.945 0.929 9.57 7.79 29.3
## 16 18 0.81000 -0.954 0.938 8.52 6.81 26.9
## 17 19 0.82600 -0.972 0.954 7.64 5.97 25.1
## 18 20 0.84700 -1.000 0.954 6.88 5.26 23.6
## Power SFT.R.sq slope truncated.R.sq mean.k. median.k. max.k.
## 1 3 0.6910 0.7290 0.6260 269.0 284.0 407.0
## 2 4 0.4540 0.3460 0.3270 199.0 209.0 329.0
## 3 5 0.0802 0.1070 -0.0486 153.0 158.0 273.0
## 4 6 0.0459 -0.0746 -0.0302 121.0 123.0 233.0
## 5 7 0.2620 -0.2140 0.3480 98.1 97.1 201.0
## 6 8 0.4070 -0.3270 0.4820 80.7 78.1 176.0
## 7 9 0.5030 -0.4180 0.5960 67.4 63.4 155.0
## 8 10 0.4930 -0.4950 0.5950 57.0 52.5 138.0
## 9 11 0.5730 -0.5520 0.6790 48.8 43.8 124.0
## 10 12 0.5770 -0.6090 0.6860 42.1 37.2 112.0
## 11 13 0.6560 -0.6500 0.7490 36.7 31.6 101.0
## 12 14 0.6740 -0.7010 0.7630 32.2 27.5 92.1
## 13 15 0.7240 -0.7330 0.7930 28.4 23.8 84.2
## 14 16 0.7800 -0.7460 0.8470 25.3 20.6 77.3
## 15 17 0.8330 -0.7540 0.9060 22.6 18.0 71.2
## 16 18 0.8560 -0.7740 0.9200 20.3 15.8 65.7
## 17 19 0.8690 -0.7870 0.9320 18.3 13.9 60.9
## 18 20 0.8810 -0.8120 0.9240 16.6 12.3 56.6
# Infer GCNs
gcns <- lapply(seq_along(power_ortho), function(n)
exp2gcn(ortho_exp[[n]], SFTpower = power_ortho[[n]]$power,
cor_method = "pearson")
)
## ..connectivity..
## ..matrix multiplication (system BLAS)..
## ..normalization..
## ..done.
## ..connectivity..
## ..matrix multiplication (system BLAS)..
## ..normalization..
## ..done.
length(gcns)
## [1] 2
Initially, module preservation analyses were performed with WGCNA’s
algorithm (Langfelder and Horvath 2008). However, the summary preservation statistics
used by WGCNA rely on parametric assumptions that are often not met.
For this purpose, the NetRep algorithm (Ritchie et al. 2016) is more accurate
than WGCNA, as it uses non-parametric permutation analyses.
Both algorithms are implemented in BioNERO
for comparison purposes,
but we strongly recommend using the NetRep algorithm.
Module preservation analysis can be performed with a
single function: module_preservation()
.
# Using rice as reference and maize as test
pres <- module_preservation(ortho_exp,
ref_net = gcns[[1]],
test_net = gcns[[2]],
algorithm = "netrep")
## [2021-08-29 06:50:43 EDT] Validating user input...
## [2021-08-29 06:50:43 EDT] Checking matrices for problems...
## [2021-08-29 06:50:44 EDT] Input ok!
## [2021-08-29 06:50:44 EDT] Calculating preservation of network subsets from
## dataset "osa" in dataset "zma".
## [2021-08-29 06:50:44 EDT] Pre-computing network properties in dataset
## "osa"...
## [2021-08-29 06:50:46 EDT] Calculating observed test statistics...
## [2021-08-29 06:50:46 EDT] Generating null distributions from 1000
## permutations using 1 thread...
##
##
0% completed.
98% completed.
100% completed.
##
## [2021-08-29 06:50:48 EDT] Calculating P-values...
## [2021-08-29 06:50:48 EDT] Collating results...
## [2021-08-29 06:50:50 EDT] Done!
## None of the modules in osa were preserved in zma.
None of the modules in rice were preserved in maize. This can be either due to the small number of orthogroups we have chosen or to the natural biological variation between species and sampled tissues. You can (and should) include more orthogroups in your analyses for a better view of transcriptional conservation between species.
Finally, BioNERO
can identify singletons and duplicated genes
with is_singleton()
. This function returns logical vectors indicating if
each of the input genes is a singleton or not.
# Sample 50 random genes
genes <- sample(rownames(zma.se), size = 50)
is_singleton(genes, og.zma.osa)
## Zm00001d000240 Zm00001d006084 Zm00001d007075 Zm00001d007888 Zm00001d008982
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d009683 Zm00001d013642 Zm00001d013678 Zm00001d013730 Zm00001d014286
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d014770 Zm00001d014993 Zm00001d015585 Zm00001d016894 Zm00001d016988
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d018142 Zm00001d018636 Zm00001d019283 Zm00001d021464 Zm00001d022182
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d023452 Zm00001d023655 Zm00001d025946 Zm00001d027337 Zm00001d027505
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d027976 Zm00001d028661 Zm00001d030955 Zm00001d031310 Zm00001d033532
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d034152 Zm00001d035039 Zm00001d035933 Zm00001d037263 Zm00001d038003
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d038361 Zm00001d038685 Zm00001d039684 Zm00001d039816 Zm00001d042108
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d042506 Zm00001d043483 Zm00001d044525 Zm00001d045195 Zm00001d045763
## TRUE TRUE TRUE TRUE TRUE
## Zm00001d046937 Zm00001d047187 Zm00001d048627 Zm00001d049325 Zm00001d050347
## TRUE TRUE TRUE TRUE TRUE
This vignette was created under the following conditions:
## R version 4.1.1 (2021-08-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.2 LTS
##
## Matrix products: default
## BLAS: /home/biocbuild/bbs-3.13-bioc/R/lib/libRblas.so
## LAPACK: /home/biocbuild/bbs-3.13-bioc/R/lib/libRlapack.so
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_GB LC_COLLATE=C
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] BioNERO_1.0.4 BiocStyle_2.20.2
##
## loaded via a namespace (and not attached):
## [1] readxl_1.3.1 backports_1.2.1
## [3] circlize_0.4.13 Hmisc_4.5-0
## [5] plyr_1.8.6 igraph_1.2.6
## [7] splines_4.1.1 GENIE3_1.14.0
## [9] BiocParallel_1.26.2 GenomeInfoDb_1.28.2
## [11] ggnetwork_0.5.10 ggplot2_3.3.5
## [13] sva_3.40.0 digest_0.6.27
## [15] foreach_1.5.1 htmltools_0.5.2
## [17] magick_2.7.3 GO.db_3.13.0
## [19] fansi_0.5.0 magrittr_2.0.1
## [21] checkmate_2.0.0 memoise_2.0.0
## [23] cluster_2.1.2 doParallel_1.0.16
## [25] limma_3.48.3 openxlsx_4.2.4
## [27] ComplexHeatmap_2.8.0 fastcluster_1.2.3
## [29] Biostrings_2.60.2 annotate_1.70.0
## [31] matrixStats_0.60.1 jpeg_0.1-9
## [33] colorspace_2.0-2 ggrepel_0.9.1
## [35] blob_1.2.2 haven_2.4.3
## [37] xfun_0.25 dplyr_1.0.7
## [39] crayon_1.4.1 RCurl_1.98-1.4
## [41] jsonlite_1.7.2 genefilter_1.74.0
## [43] impute_1.66.0 survival_3.2-13
## [45] iterators_1.0.13 glue_1.4.2
## [47] gtable_0.3.0 zlibbioc_1.38.0
## [49] XVector_0.32.0 GetoptLong_1.0.5
## [51] DelayedArray_0.18.0 car_3.0-11
## [53] shape_1.4.6 BiocGenerics_0.38.0
## [55] abind_1.4-5 scales_1.1.1
## [57] edgeR_3.34.0 DBI_1.1.1
## [59] rstatix_0.7.0 Rcpp_1.0.7
## [61] xtable_1.8-4 htmlTable_2.2.1
## [63] clue_0.3-59 foreign_0.8-81
## [65] bit_4.0.4 preprocessCore_1.54.0
## [67] Formula_1.2-4 stats4_4.1.1
## [69] htmlwidgets_1.5.3 httr_1.4.2
## [71] RColorBrewer_1.1-2 ellipsis_0.3.2
## [73] farver_2.1.0 pkgconfig_2.0.3
## [75] XML_3.99-0.7 nnet_7.3-16
## [77] sass_0.4.0 locfit_1.5-9.4
## [79] utf8_1.2.2 dynamicTreeCut_1.63-1
## [81] labeling_0.4.2 reshape2_1.4.4
## [83] tidyselect_1.1.1 rlang_0.4.11
## [85] AnnotationDbi_1.54.1 cellranger_1.1.0
## [87] munsell_0.5.0 tools_4.1.1
## [89] cachem_1.0.6 generics_0.1.0
## [91] RSQLite_2.2.8 statnet.common_4.5.0
## [93] broom_0.7.9 evaluate_0.14
## [95] stringr_1.4.0 fastmap_1.1.0
## [97] yaml_2.2.1 RhpcBLASctl_0.20-137
## [99] knitr_1.33 bit64_4.0.5
## [101] zip_2.2.0 purrr_0.3.4
## [103] KEGGREST_1.32.0 nlme_3.1-152
## [105] compiler_4.1.1 rstudioapi_0.13
## [107] curl_4.3.2 png_0.1-7
## [109] ggsignif_0.6.2 minet_3.50.0
## [111] tibble_3.1.4 statmod_1.4.36
## [113] geneplotter_1.70.0 bslib_0.2.5.1
## [115] stringi_1.7.4 highr_0.9
## [117] forcats_0.5.1 lattice_0.20-44
## [119] Matrix_1.3-4 vctrs_0.3.8
## [121] networkD3_0.4 pillar_1.6.2
## [123] lifecycle_1.0.0 BiocManager_1.30.16
## [125] jquerylib_0.1.4 GlobalOptions_0.1.2
## [127] cowplot_1.1.1 data.table_1.14.0
## [129] bitops_1.0-7 GenomicRanges_1.44.0
## [131] R6_2.5.1 latticeExtra_0.6-29
## [133] bookdown_0.23 network_1.17.1
## [135] rio_0.5.27 gridExtra_2.3
## [137] IRanges_2.26.0 codetools_0.2-18
## [139] assertthat_0.2.1 SummarizedExperiment_1.22.0
## [141] DESeq2_1.32.0 rjson_0.2.20
## [143] S4Vectors_0.30.0 GenomeInfoDbData_1.2.6
## [145] intergraph_2.0-2 mgcv_1.8-36
## [147] hms_1.1.0 parallel_4.1.1
## [149] grid_4.1.1 rpart_4.1-15
## [151] tidyr_1.1.3 NetRep_1.2.4
## [153] coda_0.19-4 rmarkdown_2.10
## [155] carData_3.0-4 MatrixGenerics_1.4.3
## [157] Cairo_1.5-12.2 ggpubr_0.4.0
## [159] ggnewscale_0.4.5 Biobase_2.52.0
## [161] WGCNA_1.70-3 base64enc_0.1-3
Emms, David M., and Steven Kelly. 2015. “OrthoFinder: solving fundamental biases in whole genome comparisons dramatically improves orthogroup inference accuracy.” Genome Biology 16 (1): 1–14. https://doi.org/10.1186/s13059-015-0721-2.
Langfelder, Peter, and Steve Horvath. 2008. “WGCNA: an R package for weighted correlation network analysis.” BMC Bioinformatics 9 (1): 559. https://doi.org/10.1186/1471-2105-9-559.
Ritchie, Scott C., Stephen Watts, Liam G. Fearnley, Kathryn E. Holt, Gad Abraham, and Michael Inouye. 2016. “A Scalable Permutation Approach Reveals Replication and Preservation Patterns of Network Modules in Large Datasets.” Cell Systems 3 (1): 71–82. https://doi.org/10.1016/j.cels.2016.06.012.
Shin, Junha, Harald Marx, Alicia Richards, Dries Vaneechoutte, Dhileepkumar Jayaraman, Junko Maeda, Sanhita Chakraborty, et al. 2020. “A network-based comparative framework to study conservation and divergence of proteomes in plant phylogenies.” Nucleic Acids Research, 1–23. https://doi.org/10.1093/nar/gkaa1041.
Van Bel, Michiel, Tim Diels, Emmelien Vancaester, Lukasz Kreft, Alexander Botzki, Yves Van De Peer, Frederik Coppens, and Klaas Vandepoele. 2018. “PLAZA 4.0: An integrative resource for functional, evolutionary and comparative plant genomics.” Nucleic Acids Research 46 (D1): D1190–D1196. https://doi.org/10.1093/nar/gkx1002.