suppressPackageStartupMessages({
if (!require("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install(version = "3.18", force = TRUE)
if (!require("SummarizedExperiment", quietly = TRUE))
BiocManager::install("SummarizedExperiment", quiet = TRUE)
library(SummarizedExperiment)
if (!require("edgeR", quietly = TRUE))
BiocManager::install("edgeR", quiet = TRUE)
library(edgeR)
if (!require("stringr", quietly = TRUE))
BiocManager::install("stringr", quiet = TRUE)
library(stringr)
if (!require("ggplot2", quietly = TRUE))
BiocManager::install("ggplot2", quiet = TRUE)
library(ggplot2)
if (!require("biomaRt", quietly = TRUE))
BiocManager::install("biomaRt", quiet = TRUE)
library(biomaRt)
if (!require("org.Mm.eg.db", quietly = TRUE))
BiocManager::install("org.Mm.eg.db", quiet = TRUE)
library(org.Mm.eg.db)
if (!require("org.Hs.eg.db", quietly = TRUE))
BiocManager::install("org.Hs.eg.db", quiet = TRUE)
library(org.Hs.eg.db)
if (!require("msigdbr", quietly = TRUE))
BiocManager::install("msigdbr", quiet = TRUE)
library(msigdbr)
if (!require("survival", quietly = TRUE))
BiocManager::install("survival", quiet = TRUE)
library(survival)
if (!require("ggpubr", quietly = TRUE))
BiocManager::install("ggpubr", quiet = TRUE)
library(ggpubr)
if (!require("ggfortify", quietly = TRUE))
BiocManager::install("ggfortify", quiet = TRUE)
library(ggfortify)
if (!require("ggrepel", quietly = TRUE))
BiocManager::install("ggrepel", quiet = TRUE)
library(ggrepel)
if (!require("survminer", quietly = TRUE))
BiocManager::install("survminer", quiet = TRUE)
library(survminer)
if (!require("sva", quietly = TRUE))
BiocManager::install("sva", quiet = TRUE)
library(sva)
if (!require("ComplexHeatmap", quietly = TRUE))
BiocManager::install("ComplexHeatmap", quiet = TRUE)
library(ComplexHeatmap)
if (!require("viridis", quietly = TRUE))
BiocManager::install("viridis", quiet = TRUE)
library(viridis)
if (!require("knitr", quietly = TRUE))
BiocManager::install("knitr", quiet = TRUE)
library(knitr)
})RNAseq analysis of the mouse mammary tumors derived from PROCR+ mammary stem cells, bulk tumors
II. Survival analysis
Introduction
In this post we will assess the PROCR metagene (signature) expression across triple-negative breast cancer (TNBC) samples and its prognostic valie for patients. The PROCR metagene was obtained via non-negative matrix factorization (NMF) from the RNA count matrix. We integrated TCGA and Metabric public data sets, selected breast tumors that were clinically classified as triple negative, and split them according to the exposure to the PROCR metagene. Then, we compared the survival of those with high and low exposure using Cox regression.
Reading the patients’ data
# when downloaded from Zenodo, data stored in input_data folder
z.Sig.Claudin <- read.csv("../RNAseq_PROCRpos-derived_tumor_pieces/input_data/Genes_Sig4.csv")
colnames(z.Sig.Claudin) <- c("gene", "logCPM")
getwd()[1] "/Users/michal/Documents/Projects/Marion/SalvadorM_PROCR_project/Presentation"Now, we loaded the Metabric RNAseq data from the input_data folder.
mb.rna <- read.table("../RNAseq_PROCRpos-derived_tumor_pieces/input_data/data_expression_median.txt",
header = T, sep = "\t", check.names = F)When exploring the data, we can see that in the first two columns correspond to the gene symbol and an Entrezid, respectively. The remaining columns represent patients in the cohort (1904). The data are also log2-transformed and median-normalized.
kable(mb.rna[1:5, 1:8])| Hugo_Symbol | Entrez_Gene_Id | MB-0362 | MB-0346 | MB-0386 | MB-0574 | MB-0503 | MB-0641 |
|---|---|---|---|---|---|---|---|
| RERE | 473 | 8.676978 | 9.653590 | 9.033589 | 8.814856 | 9.274265 | 9.286585 |
| RNF165 | 494470 | 6.075331 | 6.687887 | 5.910885 | 5.628740 | 5.908698 | 6.206729 |
| CD049690 | NA | 5.453928 | 5.454185 | 5.501577 | 5.471941 | 5.531743 | 5.372668 |
| BC033982 | NA | 4.994525 | 5.346010 | 5.247467 | 5.316523 | 5.244094 | 5.167365 |
| PHF7 | 51533 | 5.838270 | 5.600876 | 6.030718 | 5.849428 | 5.964661 | 5.783374 |
dim(mb.rna)[1] 24368 1906In the next steps, we cleaned the data set. First, we removed genes without an Entrezid or those with duplicated ids. Afterwards, we also removed genes with NA values. In the last steps, we set Entrezids as row identifiers and split the gene count matrix (mb.rna) and the gene annotation (mb.genes)
mb.rna <- mb.rna[!is.na(mb.rna$Entrez_Gene_Id), ]
mb.rna <- mb.rna[!duplicated(mb.rna$Entrez_Gene_Id), ]
mb.rna <- mb.rna[!is.na(rowSums(mb.rna[, 3:ncol(mb.rna)])), ]
rownames(mb.rna) <- mb.rna$Entrez_Gene_Id
mb.genes <- mb.rna[, c(1, 2)]
mb.rna <- mb.rna[, -c(1, 2)]Next, we loaded the clinical data of the patients. They contain information on the survival status and biomarker expression (ER: estrogen, PR: progesterone, HER2: human epidermal growth factor receptor 2 ). We completed the annotation by adding a status column to the data dividing the patients into clinical subgroups (triple-negative - TNBC, Her2 over-expressing, luminal A and B). We also removed the patients’ samples when the necessary information was missing.
metb <- read.csv("../RNAseq_PROCRpos-derived_tumor_pieces/input_data/brca_metabric_clinical_data_Aug2019.tsv",
sep = "\t", stringsAsFactors = F)
selInd <- !(is.na(metb$Overall.Survival..Months.) | is.na(metb$Overall.Survival.Status) |
is.na(metb$ER.Status) | is.na(metb$PR.Status) | is.na(metb$HER2.Status))
metb <- metb[selInd, ]
status <- rep("", nrow(metb))
status[metb$ER.Status == "Negative" & metb$PR.Status == "Negative" & metb$HER2.Status ==
"Negative"] <- "TNBC"
status[metb$ER.Status == "Negative" & metb$PR.Status == "Negative" & metb$HER2.Status ==
"Positive"] <- "Her2 over-expressing"
status[metb$ER.Status == "Positive" & metb$PR.Status == "Positive"] <- "Luminal A"
status[metb$ER.Status == "Positive" & metb$PR.Status == "Negative"] <- "Luminal B"
metb$status <- status
metb <- metb[metb$status != "", ]
rownames(metb) <- metb$Patient.ID
metb <- metb[metb$Overall.Survival.Status != "", ]
metb$Overall.Survival.Status[metb$Overall.Survival.Status == "DECEASED"] = 1
metb$Overall.Survival.Status[metb$Overall.Survival.Status == "LIVING"] = 0
metb$Overall.Survival.Status = as.numeric(metb$Overall.Survival.Status)Similarly, we loaded the TCGA expression and clinical data and performed similar steps.
# Expresssion Data
tc.rna <- read.csv("../RNAseq_PROCRpos-derived_tumor_pieces/input_data/brcarsemfpkmtcgat.txt.gz",
sep = "\t", stringsAsFactors = FALSE, check.names = FALSE)dim(tc.rna)[1] 19738 984kable(tc.rna[1:4, 1:8])| Hugo_Symbol | Entrez_Gene_Id | TCGA-AC-A2BM-01A-11R-A21T-07 | TCGA-A8-A075-01A-11R-A084-07 | TCGA-A2-A0EM-01A-11R-A034-07 | TCGA-A7-A3J1-01A-11R-A213-07 | TCGA-BH-A0BO-01A-23R-A12D-07 | TCGA-A8-A08L-01A-11R-A00Z-07 |
|---|---|---|---|---|---|---|---|
| KIF1C | 10749 | 974.50 | 1225.22 | 1304.15 | 1794.29 | 1561.89 | 480.04 |
| SWAP70 | 23075 | 387.02 | 575.03 | 733.19 | 633.73 | 1674.06 | 426.57 |
| PAK1 | 5058 | 748.61 | 693.58 | 1583.71 | 1059.11 | 1015.93 | 173.85 |
| KRTAP9-1 | 728318 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 | 0.00 |
colnames(tc.rna) <- substr(colnames(tc.rna), 1, 15)
tc.rna <- tc.rna[!duplicated(tc.rna$Entrez_Gene_Id), ]
rownames(tc.rna) <- tc.rna$Entrez_Gene_IdWe separated the first two columns with gene annotations and the expression matrix. Also, this time, the data were not log2 transformed, so we transformed them.
genes.tc <- tc.rna[, 1:2]
tc.rna <- tc.rna[, -c(1, 2)]
tc.rna <- log2(tc.rna + 1)Next, we downloaded the annotation file and assigned clinical status to the patients.
# TCGA Clinical Data
tcga <- read.csv("../RNAseq_PROCRpos-derived_tumor_pieces/input_data/brca_tcga_clinical_data_Aug2019.tsv",
sep = "\t", stringsAsFactors = F)
selInd <- !(is.na(tcga$Disease.Free..Months.) | is.na(tcga$Disease.Free.Status) |
is.na(tcga$ER.Status.By.IHC) | is.na(tcga$PR.status.by.ihc) | is.na(tcga$IHC.HER2))
tcga <- tcga[selInd, ]
status <- rep("", nrow(tcga))
status[tcga$ER.Status.By.IHC == "Negative" & tcga$PR.status.by.ihc == "Negative" &
tcga$IHC.HER2 == "Negative"] <- "TNBC"
status[tcga$ER.Status.By.IHC == "Negative" & tcga$PR.status.by.ihc == "Negative" &
tcga$IHC.HER2 == "Positive"] <- "Her2 over-expressing"
status[tcga$ER.Status.By.IHC == "Positive" & tcga$PR.status.by.ihc == "Positive"] <- "Luminal A"
status[tcga$ER.Status.By.IHC == "Positive" & tcga$PR.status.by.ihc == "Negative"] <- "Luminal B"
tcga$status <- status
tcga <- tcga[tcga$status != "", ]
rownames(tcga) <- tcga$Sample.ID
tcga$Overall.Survival.Status[tcga$Overall.Survival.Status == "DECEASED"] = 1
tcga$Overall.Survival.Status[tcga$Overall.Survival.Status == "LIVING"] = 0
tcga$Overall.Survival.Status = as.numeric(tcga$Overall.Survival.Status)Signature matrix from NMF as a projector, batch correction
The goal of this step was to search for the “fingerprints” of the identified PROCR-metagene in TNBC samples. Recalling the NMF decomposition of the RNA count matrix C \approx WH, \dim{C}\{n,m\}, the signature matrix W, \dim{W} =\{n,k\}, encodes the structure of all k metagenes (signatures) as a linear combination of the expression of all n genes with non-negative coefficients. Given an RNA count matrix C_u, \dim{C_u} = \{n',m'\} of m' unknown samples, we can extract how much are these samples exposed to individual metagenes using transposed W^T as a projector.
W^T C_u = H_u, \ \ \ H_u: \ \ {\rm exposure \ matrix \ of \ the \ unknown \ samples.} By this means, we can say how much each an unknown sample expresses the metagene of interest. Before performing the projection, we only need to make sure that the “genomic dimensions” are consistently aligned. It may happen, that the number of genes n detected in the in the original experiment for which NMF was performed differs from those detected in the unknown data set under investigation with n' detected genes. So we need to filter both matrices W and C_u in such a way that their rows correspond to the same shared genes in both data sets.
Next, we loaded the projection matrix W from the input_data folder again. It was extracted filtered to 1-to-1 human orthologs (the original decomposition was derived from mouse mammary tumors, here we tested the relevance in human).
procr_Wspecific <- read.csv(paste0("../RNAseq_PROCRpos-derived_tumor_pieces/input_data/W_projector.csv"),
row.names = 1)Next, we filtered the TNBC patients in both reference cohorts and sorted both expression matrices from TCGA and Metabric together with the W matrix to contain only shared genes in the same order.
# select only TNBC
selTumor <- "TNBC"
tumors <- intersect(rownames(metb), colnames(mb.rna))
tumors <- tumors[metb[tumors, "status"] == selTumor]
metb.1 <- metb[tumors, ]
mb.rna.sel <- mb.rna[mb.genes$Hugo_Symbol %in% toupper(rownames(procr_Wspecific)),
tumors]
rownames(mb.rna.sel) <- mb.genes$Hugo_Symbol[mb.genes$Hugo_Symbol %in% toupper(rownames(procr_Wspecific))]
procr_Wspecific <- procr_Wspecific[rownames(mb.rna.sel), ] #sort genes
dim(mb.rna.sel) # 3575 299[1] 3575 299dim(procr_Wspecific) #[1] 3575 4[1] 3575 4# ________________
tumors.tc <- intersect(rownames(tcga), colnames(tc.rna))
tumors.tc <- tumors.tc[tcga[tumors.tc, "status"] == selTumor]
tcga.1 <- tcga[tumors.tc, ]
mb.1.sel <- metb.1[, c("Overall.Survival..Months.", "Overall.Survival.Status")]
tc.1.sel <- tcga.1[, c("Overall.Survival..Months.", "Overall.Survival.Status")]
comb.1.sel <- rbind(mb.1.sel, tc.1.sel)
tc.rna.sel <- tc.rna[genes.tc$Hugo_Symbol %in% toupper(rownames(procr_Wspecific)),
tumors.tc]
rownames(tc.rna.sel) <- genes.tc$Hugo_Symbol[genes.tc$Hugo_Symbol %in% toupper(rownames(procr_Wspecific))]
dim(tc.rna.sel) #3554 92[1] 3554 92dim(procr_Wspecific) #[1] 3575 4[1] 3575 4Next, we organized the data according to the shared genes.
# intersect genes!
shared.genes <- intersect(rownames(mb.rna.sel), rownames(tc.rna.sel))
length(shared.genes) # 33554[1] 3554tc.rna.sel <- tc.rna.sel[shared.genes, ]
mb.rna.sel <- mb.rna.sel[shared.genes, ]
procr_Wspecific <- procr_Wspecific[shared.genes, ]In the next step, we regressed out the batch effect between the TCGA and Metabric data using ComBat
comb.rna.sel <- cbind(mb.rna.sel, tc.rna.sel)
batch <- c(rep("mb", dim(mb.rna.sel)[2]), rep("tc", dim(tc.rna.sel)[2]))
# parametric adjustment: ComBat
combat_comb.rna.sel = ComBat(dat = comb.rna.sel, batch = batch, mod = NULL, par.prior = TRUE,
prior.plots = FALSE)Lastly, we needed to transform back the log2 expression data by exponentiation (non-negative count matrix is needed to perform the projection via NMF)
countmatrix.comb.rna.sel <- 2^combat_comb.rna.selPROCR-metagene expression in the TNBC patients
We performed the projection and plotted the exposure matrix
ExposureTNBC <- t(as.matrix(procr_Wspecific)) %*% as.matrix(countmatrix.comb.rna.sel)
h_TNBC_4 <- Heatmap(ExposureTNBC, col = viridis(100), clustering_distance_columns = "pearson",
show_column_dend = TRUE, name = "TNBC, exposure", show_column_names = FALSE,
show_row_names = TRUE, cluster_rows = FALSE)
# exposure matrix of the unknown samples
print(h_TNBC_4)
The fourth signature is the identified PROCR-metagene which we further focused on.
ExposureTNBC.SigPROCR <- ExposureTNBC[4, ]
ExposureTNBC.SigPROCR <- sort(ExposureTNBC.SigPROCR, decreasing = TRUE)
summary(ExposureTNBC.SigPROCR) Min. 1st Qu. Median Mean 3rd Qu. Max.
217.6 556.0 690.5 770.2 916.6 2393.7 Next, we investigated the distribution of exposures among the patients (Figure S10 B).
max.point.eff <- unname(sort(ExposureTNBC.SigPROCR, decreasing = TRUE)[2])
min.point <- min(ExposureTNBC.SigPROCR)
dif <- max.point.eff - min.point
midpoint <- min(ExposureTNBC.SigPROCR) + dif/2
hist(ExposureTNBC.SigPROCR, breaks = 30)
abline(v = midpoint, col = "blue")
The histogram had a distinct tail towards higher values indicating the existence of a group of patients with high exposure. We then split the patients into two groups: PROCR-metagene low and high. To do that, we identified the compact support of the histogram (meaning without the outlier with 2393.7 exposure) and its midpoint served as a limiting value for deciding to which group a given patient would be assigned. In detail, max.point.eff set the highest value of the compact support (note that here we neglected the outlier), min.point was the minimal value. The midpoint is marked by the vertical blue line in the plot.
ExposureStatus <- ifelse(ExposureTNBC.SigPROCR >= midpoint, "High", "Low")
top.TNBC <- names(ExposureTNBC.SigPROCR[ExposureStatus == "High"])
low.TNBC <- names(ExposureTNBC.SigPROCR[ExposureStatus == "Low"])As a next step, we investigated differences in the molecular compositions between these groups of patients. Interestingly, The PROCR expression was elevated by circa 36% in the group of high PROCR-metagene exposure (Figure 3E).
df.procr.combat <- data.frame(patients = colnames(combat_comb.rna.sel), lognorm.procr = combat_comb.rna.sel["PROCR",
], status = ifelse(colnames(combat_comb.rna.sel) %in% low.TNBC, "low exposure",
"high exposure"))
table(df.procr.combat$status)
high exposure low exposure
71 320 p <- ggboxplot(df.procr.combat, y = "lognorm.procr", x = "status", color = "status",
palette = c("cadetblue3", "chocolate1"), add = "jitter")
p + stat_compare_means(method = "t.test") + ggtitle("PROCR expression in survival groups \n (TCGA and Metabric)") +
ylim(5.2, 9.2)
Within the Metabric cohort, the clinical data also contained information on PAM50 molecular subtyping (Parker et al. 2009) including the Claudin-low subtype. So next, we checked what is the composition of these subtypes in the patients’ tumors with high exposure PROCR metagene exposure (Figure 3D).
data <- data.frame(table(metb[top.TNBC, c("Pam50...Claudin.low.subtype")]))
names(data) <- c("class", "high_exp")
df <- data %>%
mutate(perc = high_exp/sum(high_exp)) %>%
mutate(labels = scales::percent(perc))
bp <- ggplot(df, aes(x = "", y = high_exp, fill = class)) + geom_bar(width = 1, stat = "identity")
pie.high <- bp + coord_polar("y", start = 0)
my_colors <- c("coral2", "cornflowerblue", "darkolivegreen4", "darkorchid4", "darkorange")
pie.high + scale_fill_manual(values = my_colors, labels = c("Basal", "claudin-low",
"Her2", "LumA", "Normal")) + theme_void() + theme(axis.text.x = element_blank()) +
geom_text_repel(aes(label = high_exp), position = position_stack(vjust = 0.6),
size = 5) + ggtitle("Patients with high exposure (total metabric patients 56)")
Conversly, we also checked the composition in the lowly-exposed group.
data.low <- data.frame(table(metb[low.TNBC, c("Pam50...Claudin.low.subtype")]))
names(data.low) <- c("class", "low_exp")
df <- data.low %>%
mutate(perc = low_exp/sum(low_exp)) %>%
mutate(labels = scales::percent(perc))
bp <- ggplot(df, aes(x = "", y = low_exp, fill = class)) + geom_bar(width = 1, stat = "identity")
pie.low <- bp + coord_polar("y", start = 0)
pie.low + scale_fill_manual(values = my_colors[c(1, 2, 3, 5)], labels = c("Basal",
"claudin-low", "Her2", "Normal")) + theme_void() + theme(axis.text.x = element_blank()) +
geom_text_repel(aes(label = low_exp), position = position_stack(vjust = 0.6),
size = 5) + ggtitle("Patients with low exposure (total metabric patients 243)")
In the “high exposure” group, the dominant class is Claudin-low which, again, corresponds with our observation from the mouse model that the PROCR+ MaSCs, when transformed give dominantly rise to Claudin-low tumors and the PROCR-metagene captures this relation.
Survival analysis
Finally, we performed a survival analysis and plotted the Kaplan-Meier curve (Figure 3C).
comb.1.sel$signature <- ifelse(names(ExposureTNBC.SigPROCR) %in% top.TNBC, "High", ifelse(names(ExposureTNBC.SigPROCR) %in% low.TNBC,"Low","Mid"))
fit <- survfit(Surv(Overall.Survival..Months., Overall.Survival.Status) ~ signature, data = comb.1.sel[comb.1.sel$signature %in% c("Low","High") , ])
ggsurvplot(fit, data = comb.1.sel[comb.1.sel$signature %in% c("Low","High"),], risk.table = TRUE, censor.size = 5,
pval = TRUE,
title = "Metabric + TCGA, TNBC",
legend.labs = c("High", "Low"), # Change legend labels
risk.table.height = 0.3, # Useful to change when you have multiple groups
xlab = "Time (months)")
Session Info
sessionInfo()R version 4.3.2 (2023-10-31)
Platform: aarch64-apple-darwin20 (64-bit)
Running under: macOS 26.2
Matrix products: default
BLAS: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRblas.0.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/4.3-arm64/Resources/lib/libRlapack.dylib; LAPACK version 3.11.0
locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
time zone: Europe/Zurich
tzcode source: internal
attached base packages:
[1] stats4 grid stats graphics grDevices utils datasets
[8] methods base
other attached packages:
[1] viridis_0.6.5 viridisLite_0.4.3
[3] ComplexHeatmap_2.18.0 sva_3.50.0
[5] BiocParallel_1.36.0 genefilter_1.84.0
[7] mgcv_1.9-0 nlme_3.1-168
[9] survminer_0.5.0 ggrepel_0.9.6
[11] ggfortify_0.4.19 ggpubr_0.6.3
[13] survival_3.8-3 msigdbr_25.1.1
[15] org.Hs.eg.db_3.18.0 org.Mm.eg.db_3.18.0
[17] AnnotationDbi_1.64.1 biomaRt_2.58.2
[19] ggplot2_4.0.2 stringr_1.6.0
[21] edgeR_4.0.16 limma_3.58.1
[23] SummarizedExperiment_1.32.0 Biobase_2.62.0
[25] GenomicRanges_1.54.1 GenomeInfoDb_1.38.8
[27] IRanges_2.36.0 S4Vectors_0.40.2
[29] BiocGenerics_0.48.1 MatrixGenerics_1.14.0
[31] matrixStats_1.5.0 BiocManager_1.30.27
[33] png_0.1-8 knitr_1.51
loaded via a namespace (and not attached):
[1] RColorBrewer_1.1-3 shape_1.4.6.1 rstudioapi_0.18.0
[4] jsonlite_2.0.0 magrittr_2.0.4 farver_2.1.2
[7] rmarkdown_2.30 GlobalOptions_0.1.3 zlibbioc_1.48.2
[10] vctrs_0.6.5 Cairo_1.7-0 memoise_2.0.1
[13] RCurl_1.98-1.17 rstatix_0.7.3 htmltools_0.5.9
[16] S4Arrays_1.2.1 progress_1.2.3 curl_7.0.0
[19] broom_1.0.12 SparseArray_1.2.4 Formula_1.2-5
[22] htmlwidgets_1.6.4 zoo_1.8-15 cachem_1.1.0
[25] lifecycle_1.0.5 iterators_1.0.14 pkgconfig_2.0.3
[28] Matrix_1.6-4 R6_2.6.1 fastmap_1.2.0
[31] clue_0.3-66 GenomeInfoDbData_1.2.11 digest_0.6.39
[34] colorspace_2.1-2 RSQLite_2.4.1 labeling_0.4.3
[37] filelock_1.0.3 km.ci_0.5-6 httr_1.4.8
[40] abind_1.4-8 compiler_4.3.2 bit64_4.6.0-1
[43] withr_3.0.2 doParallel_1.0.17 S7_0.2.1
[46] backports_1.5.0 carData_3.0-6 DBI_1.3.0
[49] ggsignif_0.6.4 rappdirs_0.3.3 DelayedArray_0.28.0
[52] rjson_0.2.23 tools_4.3.2 otel_0.2.0
[55] glue_1.8.0 gridtext_0.1.5 cluster_2.1.8.1
[58] generics_0.1.4 gtable_0.3.6 KMsurv_0.1-6
[61] tidyr_1.3.2 data.table_1.17.8 hms_1.1.4
[64] xml2_1.3.8 car_3.1-5 XVector_0.42.0
[67] foreach_1.5.2 pillar_1.11.1 babelgene_22.9
[70] circlize_0.4.17 splines_4.3.2 ggtext_0.1.2
[73] dplyr_1.1.4 BiocFileCache_2.10.2 lattice_0.22-7
[76] bit_4.6.0 annotate_1.80.0 tidyselect_1.2.1
[79] locfit_1.5-9.12 Biostrings_2.70.3 gridExtra_2.3
[82] xfun_0.52 statmod_1.5.1 stringi_1.8.7
[85] yaml_2.3.12 evaluate_1.0.5 codetools_0.2-20
[88] tibble_3.3.0 cli_3.6.5 xtable_1.8-8
[91] survMisc_0.5.1 Rcpp_1.1.0 dbplyr_2.5.2
[94] XML_3.99-0.18 parallel_4.3.2 assertthat_0.2.1
[97] blob_1.3.0 prettyunits_1.2.0 bitops_1.0-9
[100] scales_1.4.0 purrr_1.0.4 crayon_1.5.3
[103] GetoptLong_1.1.0 rlang_1.1.6 formatR_1.14
[106] KEGGREST_1.42.0 References
Parker, Joel S., Michael Mullins, Maggie C. U. Cheang, Samuel Leung, David Voduc, Tammi Vickery, Sherri Davies, et al. 2009. “Supervised Risk Predictor of Breast Cancer Based on Intrinsic Subtypes.” Journal of Clinical Oncology 27 (8): 1160–67. https://doi.org/10.1200/JCO.2008.18.1370.