library(survival)
library(reactable)
select <- c("ENSG00000135365", "ENSG00000162747", "ENSG00000180871", "ENSG00000186529", "ENSG00000105835",
"ENSG00000081237", "ENSG00000090104", "ENSG00000211772", "ENSG00000117091","ENSG00000145649",
"ENSG00000106066", "ENSG00000140968", "ENSG00000100079", "ENSG00000127951", "ENSG00000128815",
"ENSG00000112149", "ENSG00000167604", "ENSG00000166501", "ENSG00000184588", "ENSG00000173193",
"ENSG00000158050", "ENSG00000111796", "ENSG00000188042", "ENSG00000101109", "ENSG00000277734")
Adeno_bulk <- readRDS("~/new_data/Adeno_bulk.rds")
bulk <- Adeno_bulk
bulk <- bulk[select,]
bulk <- t(bulk)
bulk <- as.data.frame(bulk)
LUAD_BULK_CLINIC <- readRDS("~/new_data/LUAD_BULK_CLINIC.rds")
select_col <- c("overall_survival", "vital_status")
clinic <- LUAD_BULK_CLINIC[,select_col]
all_clinic <- clinic
data <- cbind(all_clinic,bulk)
data$vital_status <- ifelse(data$vital_status == "Alive", 1, 0)
data_clean <- na.omit(data)
univ_cox <- function(variable) {
cox_model <- coxph(Surv(overall_survival, vital_status) ~ variable, data = data_clean)
cox_summary <- summary(cox_model)
c(
coef = cox_summary$coefficients[1],
HR = exp(cox_summary$coefficients[1]),
p_value = cox_summary$coefficients[5],
CI_lower = exp(cox_summary$conf.int[ ,3]),
CI_upper = exp(cox_summary$conf.int[ ,4])
)
}
univ_results <- lapply(data_clean[, 3:ncol(data_clean)], univ_cox)
univ_results_df <- do.call(rbind, univ_results)
univ_results_df <- as.data.frame(univ_results_df)
rownames(univ_results_df) <- colnames(data_clean[, 3:ncol(data_clean)])
reactable(univ_results_df)
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000135365"] <- "neutrophils"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000162747"] <- "neutrophils"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000180871"] <- "neutrophils"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000186529"] <- "neutrophils"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000105835"] <- "neutrophils"
#"ENSG00000105835" "ENSG00000135365" "ENSG00000162747" "ENSG00000180871" "ENSG00000186529"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000081237"] <- "t_cell_CD8_effector_memory"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000117091"] <- "t_cell_CD8_effector_memory"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000090104"] <- "t_cell_CD8_effector_memory"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000211772"] <- "t_cell_CD8_effector_memory"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000145649"] <- "t_cell_CD8_effector_memory"
#"ENSG00000081237" "ENSG00000090104" "ENSG00000211772" "ENSG00000145649", ENSG00000117091
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000158050"] <- "t_cell_CD8_activated"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000111796"] <- "t_cell_CD8_activated"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000188042"] <- "t_cell_CD8_activated"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000101109"] <- "t_cell_CD8_activated"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000277734"] <- "t_cell_CD8_activated"
#"ENSG00000158050","ENSG00000111796", "ENSG00000188042", ENSG00000101109, ENSG00000277734
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000112149"] <- "B_cell"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000167604"] <- "B_cell"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000166501"] <- "B_cell"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000184588"] <- "B_cell"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000173193"] <- "B_cell"
# "ENSG00000112149", "ENSG00000081237", "ENSG00000166501", "ENSG00000184588", "ENSG00000173193"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000106066"] <- "cDC1"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000140968"] <- "cDC1"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000100079"] <- "cDC1"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000127951"] <- "cDC1"
univ_results_df$cell_type[rownames(univ_results_df) == "ENSG00000128815"] <- "cDC1"
# "ENSG00000106066", "ENSG00000140968", "ENSG00000100079", "ENSG00000127951", "ENSG00000128815"
univ_results_df$ln_B <- abs(log(univ_results_df$coef))
#"ENSG00000105835" "ENSG00000135365" "ENSG00000162747" "ENSG00000180871" "ENSG00000186529":neu
univ_results_df$si[rownames(univ_results_df) == "ENSG00000105835"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000162747"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000180871"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000186529"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000135365"] <- 1
#"ENSG00000081237" "ENSG00000117091" "ENSG00000090104" "ENSG00000211772" "ENSG00000145649": eff
univ_results_df$si[rownames(univ_results_df) == "ENSG00000081237"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000090104"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000145649"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000211772"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000117091"] <- -1
#"ENSG00000158050", "ENSG00000111796", "ENSG00000188042", "ENSG00000101109", "ENSG00000277734": act
univ_results_df$si[rownames(univ_results_df) == "ENSG00000158050"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000111796"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000188042"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000101109"] <- 1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000277734"] <- 1
# "ENSG00000112149", "ENSG00000167604", "ENSG00000166501", "ENSG00000184588", "ENSG00000173193": B
univ_results_df$si[rownames(univ_results_df) == "ENSG00000112149"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000167604"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000166501"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000184588"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000173193"] <- -1
# "ENSG00000106066", "ENSG00000140968", "ENSG00000100079", "ENSG00000127951", "ENSG00000128815": cDC1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000106066"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000140968"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000100079"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000127951"] <- -1
univ_results_df$si[rownames(univ_results_df) == "ENSG00000128815"] <- -1
univ_results_df$beta <- univ_results_df$si * univ_results_df$ln_B
univ_results_df$ens <- rownames(univ_results_df)
library(biomaRt)
mart <- useMart("ensembl", dataset = "hsapiens_gene_ensembl")
ensembl_ids <- rownames(univ_results_df)
gene_conversion <- getBM(attributes = c("ensembl_gene_id", "hgnc_symbol"),
filters = "ensembl_gene_id",
values = ensembl_ids,
mart = mart)
univ_results_df$Gene_Name <- gene_conversion$hgnc_symbol[match(rownames(univ_results_df), gene_conversion$ensembl_gene_id)]
univ_results_df <- as.data.frame(univ_results_df)
univ_results_df <- univ_results_df[!is.na(univ_results_df$Gene_Name), ]
rownames(univ_results_df) <- univ_results_df$Gene_Name
univ_results_df$Gene_Name <- NULL
reactable(univ_results_df)
final <- data_clean[, c("ENSG00000135365", "ENSG00000162747", "ENSG00000180871", "ENSG00000186529", "ENSG00000105835",
"ENSG00000081237", "ENSG00000090104", "ENSG00000211772", "ENSG00000117091","ENSG00000145649",
"ENSG00000106066", "ENSG00000140968", "ENSG00000100079", "ENSG00000127951", "ENSG00000128815",
"ENSG00000112149", "ENSG00000167604", "ENSG00000166501", "ENSG00000184588", "ENSG00000173193",
"ENSG00000158050", "ENSG00000111796", "ENSG00000188042", "ENSG00000101109", "ENSG00000277734")]
si <- univ_results_df$si
ln_B <- univ_results_df$ln_B
final$risk_scores <- apply(final, 1, function(x) sum(x * si * ln_B))
head(final$risk_scores)
## [1] -221725.164 -56358.127 -92280.171 187195.936 -57500.873 8218.054
library(survival)
library(survminer)
## Loading required package: ggplot2
## Loading required package: ggpubr
##
## Attaching package: 'survminer'
## The following object is masked from 'package:survival':
##
## myeloma
final$overall_survival <- data_clean$overall_survival
final$vital_status <- data_clean$vital_status
final$deceased <- ifelse(final$vital_status == "1", FALSE, TRUE)
final$overall_survival <- as.numeric(final$overall_survival)
plot_list <- list()
res.cut <- surv_cutpoint(final, time = "overall_survival", event = "deceased",
variables = c("risk_scores"))
cutpoint <- res.cut[["cutpoint"]][["cutpoint"]]
final$risk_scores <- ifelse(final$risk_scores > cutpoint, "HIGH", "LOW")
fit_risk_scores <- survfit(Surv(overall_survival, deceased) ~ risk_scores, data = final)
ggsurvplot(fit_risk_scores,
final,
pval = TRUE,
)
