Analysis_meeting_11_19

1 Functions

library(stringi)
source_from_github(repositoy = "DEG_functions",version = "0.2.47")

ℹ SHA-1 hash of file is f5bb1cd741d13bded83fe3b6fd43169e29731216
Welcome to enrichR
Checking connection ... 
Enrichr ... Connection is Live!
FlyEnrichr ... Connection is available!
WormEnrichr ... Connection is available!
YeastEnrichr ... Connection is available!
FishEnrichr ... Connection is available!

source_from_github(repositoy = "cNMF_functions",version = "0.4.0",script_name = "cnmf_functions_V3.R")

ℹ SHA-1 hash of file is 7b39dfd215cf7e1d29ee976ecc6bfa0e3f1532f6
Loading required package: reticulate

source_from_github(repositoy = "sc_general_functions",version = "0.1.28",script_name = "functions.R")

ℹ SHA-1 hash of file is 737683e4e0a82ffa22769bbd4a0842a271fe63fc
Loading required package: RColorBrewer

2 Data

library("readxl")
acc_all = readRDS(file = "./Data/acc_tpm_nCount_mito_no146_15k_alldata.rds")
acc_cancer_pri = readRDS(file = "./Data/acc_cancer_no146_primaryonly15k_cancercells.rds")
acc_cancer = readRDS(file = "./Data/acc_tpm_nCount_mito_no146_15k_cancercells.rds")


neuronal_signatures <- read_excel("./Data/Neuronal Signatures.xlsx",col_names =F)
neuronal_signatures  %<>%  t() %>% as.data.frame() %>% janitor::row_to_names(1) %>%  filter(!row_number() == 1)
rownames(neuronal_signatures) <- NULL
colnames(neuronal_signatures)   %<>%  gsub(replacement = "", pattern = "_\\d.*") #remove any _numbers
colnames(neuronal_signatures)   %<>%  gsub(replacement = "", pattern = "\\(.*") #rename "(" to the end
colnames(neuronal_signatures)   %<>%  gsub(replacement = "_", pattern = " ") #rename all spaces

neuronal_pathways <- read_excel("./Data/Pathway analysis Gene sets.xlsx",col_names =F)
neuronal_pathways  %<>%  t() %>% as.data.frame() %>% janitor::row_to_names(1) %>%  filter(!row_number() == 1)  %>%  as.list() 
neuronal_pathways = lapply(neuronal_pathways, na.omit)
neuronal_pathways = lapply(neuronal_pathways, as.character)

3 Neuronal signatures

3.1 Violin plots

for (neural_name in colnames(neuronal_signatures)) {
  genes = neuronal_signatures[,neural_name,drop=T]
  # Assuming df is your data frame

  
  pathways_scores = FetchData(object = acc_all,vars = genes,slot = "data") %>% 
    mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_all  %<>% AddMetaData(metadata = pathways_scores$score,col.name = neural_name)
}

Warning in FetchData.Seurat(object = acc_all, vars = genes, slot = “data”) : The following requested variables were not found: SNHG29, NOP53 Warning in FetchData.Seurat(object = acc_all, vars = genes, slot = “data”) : The following requested variables were not found: NOP53 Warning in FetchData.Seurat(object = acc_all, vars = genes, slot = “data”) : The following requested variables were not found: NOP53

plt = VlnPlot(object = acc_all,features = colnames(neuronal_signatures))
plt[[1]] = plt[[1]]+ ylab("TPM")
plt[[4]] = plt[[4]]+ ylab("TPM")
plt[[7]] = plt[[7]]+ ylab("TPM")
print_tab(plt = plt,title = "expression with RP genes",subtitle_num = 3)

expression with RP genes

NA

rp_genes = lapply(neuronal_signatures %>% as.list(),  function(x) {x[startsWith(x = x,prefix = "RP")]}) %>% unlist() %>% unique() %>% as.data.frame()
print_tab(rp_genes,title = " RP genes names",subtitle_num = 3)

RP genes names

NA

neuronal_signatures_no_RP = lapply(neuronal_signatures %>% as.list(),  function(x) {x[!startsWith(x = x,prefix = "RP")]})
names(neuronal_signatures_no_RP) =  paste0(names(neuronal_signatures_no_RP),"_noRP")
for (neural_name in names(neuronal_signatures_no_RP)) {
  genes = neuronal_signatures_no_RP[[neural_name]]
  # Assuming df is your data frame

  pathways_scores = FetchData(object = acc_all,vars = genes,slot = "data") %>% 
    mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_all  %<>% AddMetaData(metadata = pathways_scores$score,col.name = neural_name)
}

plt = VlnPlot(object = acc_all,features = names(neuronal_signatures_no_RP))
for (i in 1:(length(plt$patches$plots)+1)) {

  plt[[i]] = plt[[i]] + theme(plot.title = element_text(size=12))
  if (i %in% c(1,4,7)) {
    plt[[i]] = plt[[i]]+ylab("TPM")
  }
}
print_tab(plt = plt,title = "expression without RP genes",subtitle_num = 3)

expression without RP genes

NA

3.2 UMAPS

plt1 = FeaturePlot(object = acc_all,features = colnames(neuronal_signatures))& theme(plot.title = element_text(size=13))
plt2 = FeaturePlot(object = acc_all,features = names(neuronal_signatures_no_RP)) & theme(plot.title = element_text(size=12))
print_tab(plt = plt1,title = "UMAP with RP",subtitle_num = 3)

UMAP with RP

print_tab(plt = plt2,title = "UMAP without RP",subtitle_num = 3)

UMAP without RP

NA

3.3 High Sympathetic_cholinergic_neuron

high_cholinergic_neuron = list( high_cholinergic_neuron = colnames(acc_all)[acc_all$Sympathetic_cholinergic_neuron_noRP>1000])
DimPlot(object = acc_all, cells.highlight = high_cholinergic_neuron, cols.highlight = "red", cols = "gray", order = TRUE)

4 Pathway comparison

for (pathway_name in names(neuronal_pathways)) {
  genes = neuronal_pathways[[pathway_name]]
  pathways_scores = FetchData(object = acc_cancer,vars = genes,slot = "data") %>% 
    mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = paste0(pathway_name,"_TPM"))
}

for (pathway_name in names(neuronal_pathways)) {
  genes = neuronal_pathways[[pathway_name]]
  pathways_scores = FetchData(object = acc_cancer,vars = genes,slot = "data") %>% 
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = pathway_name)
}

library(reshape2) 
data = FetchData(object = acc_cancer,vars = names(neuronal_pathways)[c(1:2,7,8)])
my.labels <- c("GOBP_NOREPINEPHRINE\nSECRETION",
               "GOBP_NOREPINEPHRINE\nTRANSPORT", 
               "GOBP_ADRENERGIC_RECEPTOR\nSIGNALING_PATHWAY",
               "GOMF_BETA_2_ADRENERGIC\nRECEPTOR_BINDING") 

data =  reshape2::melt(data)

No id variables; using all as measure variables

names (data) = c("pathway","log(TPM)")
ggboxplot(data, x = "pathway", y ="log(TPM)",
              palette = "jco",
              add = "jitter")+ theme( axis.text.x  = element_text(size=10))+
    stat_summary(fun.data = function(x) data.frame(y=5, label = paste("Mean=",round(mean(x),digits = 
                                    2))), geom="text")  + scale_x_discrete(labels = my.labels)+ggtitle("ACC cells pathways")

4.1 pathways UMAP

plt = FeaturePlot(object = acc_cancer,features = paste0(names(neuronal_pathways),"_TPM"),pt.size = 0.5) & theme(plot.title = element_text(size=10.5)) & labs(color='TPM') & scale_color_gradientn(colours = rainbow(5)) 

print_tab(plt,title = "TPM",subtitle_num = 3)

TPM

NA

plt1 = FeaturePlot(object = acc_cancer,features = names(neuronal_pathways),pt.size = 0.5)&
  theme(plot.title = element_text(size=10.5)) & labs(color='log (TPM+1)') & scale_color_gradientn(colours = rainbow(5))

Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale. Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale.

plt2 = FeaturePlot(object = acc_cancer,features = names(neuronal_pathways),pt.size = 0.5)&
  theme(plot.title = element_text(size=10.5)) & labs(color='log (TPM+1)')

print_tab(plt2,title = "log(TPM+1) ",subtitle_num = 3)

log(TPM+1)

print_tab(plt1,title = "log(TPM+1) rainbow",subtitle_num = 3)

log(TPM+1) rainbow

NA

5 Beta-2 receptor

5.1 UMAP

b2r_genes = c("RAF1", "PAMPK1", "PRKACA", "PIK3CA", "AKT1", "AKT2", "AKT3", "CREB1", "CREBBP", "BCL2", "BAD", "JUN", "FOS",  "MYC")
tpm_1 = FeaturePlot(object = acc_tpm,features = b2r_genes[1:9],pt.size = 0.5,slot = "counts") & scale_color_gradientn(colours = rainbow(5)) & labs(color='TPM')
tpm_2 = FeaturePlot(object = acc_tpm,features = b2r_genes[10:length(b2r_genes)],pt.size = 0.5,slot = "counts") & scale_color_gradientn(colours = rainbow(5)) & labs(color='TPM')

logTPM1 = FeaturePlot(object = acc_all,features = b2r_genes[1:9],pt.size = 0.5) & labs(color='log(TPM+1)')
logTPM2 = FeaturePlot(object = acc_all,features = b2r_genes[10:length(b2r_genes)]) & labs(color='log(TPM+1)')

print_tab(plt = tpm_1,title = "TPM part 1",subtitle_num = 3)

TPM part 1

print_tab(plt = tpm_2,title = "TPM part 2",subtitle_num = 3)

TPM part 2

print_tab(plt = logTPM1,title = "log TPM part 1",subtitle_num = 3)

log TPM part 1

print_tab(plt = logTPM2,title = "log TPM part 2",subtitle_num = 3)

log TPM part 2

NA

5.2 Signature

pathways_scores = FetchData(object = acc_all,vars = b2r_genes,slot = "data") %>% 
  mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_all  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature_TPM")

tpm = FeaturePlot(object = acc_all,features ="b2r_signature_TPM",pt.size = 0.5,slot = "counts") & scale_color_gradientn(colours = rainbow(5)) & labs(color='TPM')

Scale for ‘colour’ is already present. Adding another scale for ‘colour’, which will replace the existing scale.

pathways_scores = FetchData(object = acc_all,vars = b2r_genes,slot = "data") %>% 
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_all  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature")

logTPM = FeaturePlot(object = acc_all,features ="b2r_signature",pt.size = 0.5,slot = "counts")  & labs(color='log TPM')

print_tab(plt = tpm,title = "TPM",subtitle_num = 3)

TPM

print_tab(plt = logTPM,title = "log TPM",subtitle_num = 3)

log TPM

NA

5.3 b2r with pathways TPM

pathways_scores = FetchData(object = acc_cancer,vars = b2r_genes,slot = "data") %>% 
  mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature_TPM")

FeaturePlot(object = acc_cancer,features = c("b2r_signature_TPM",paste0(names(neuronal_pathways)[c(1:2,7,8)],"_TPM")),pt.size = 0.5,slot = "counts") & scale_color_gradientn(colours = rainbow(5)) & labs(color='TPM')&
  theme(plot.title = element_text(size=10.5))

Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.

5.4 b2r with pathways log(TPM+1)

pathways_scores = FetchData(object = acc_cancer,vars = b2r_genes,slot = "data") %>% 
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature")

FeaturePlot(object = acc_cancer,features = c("b2r_signature",names(neuronal_pathways)[c(1:2,7,8)]),pt.size = 0.5,slot = "counts") & scale_color_gradientn(colours = rainbow(5)) & labs(color='log (TPM+1)')&
  theme(plot.title = element_text(size=10.5))

Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.
Scale for 'colour' is already present. Adding another scale for 'colour', which will replace the existing scale.

5.5 Correlation TPM

pathways_scores = FetchData(object = acc_cancer,vars = b2r_genes,slot = "data") %>% 
  mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean

Warning in FetchData.Seurat(object = acc_cancer, vars = b2r_genes, slot = "data") :
  The following requested variables were not found: PAMPK1

acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature_TPM")


data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM",paste0(names(neuronal_pathways)[c(1:2,7,8)],"_TPM")))
cor_res = cor(data)
breaks <- c(seq(-1,1,by=0.01))
colors_for_plot <- colorRampPalette(colors = c("blue", "white", "red"))(n = length(breaks))

pheatmap(cor_res,color = colors_for_plot,breaks = breaks,display_numbers = T,main = "Pearson correlation of pathways") 

5.6 Correlation llog TPM

for (pathway_name in names(neuronal_pathways)) {
  genes = neuronal_pathways[[pathway_name]]
  pathways_scores = FetchData(object = acc_cancer,vars = genes,slot = "data") %>% 
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = pathway_name)
}

pathways_scores = FetchData(object = acc_cancer,vars = b2r_genes,slot = "data") %>% 
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_cancer  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature")


data = FetchData(object = acc_cancer,vars =  c("b2r_signature",names(neuronal_pathways)[c(1:2,7,8)]))
cor_res = cor(data)
breaks <- c(seq(-1,1,by=0.01))
colors_for_plot <- colorRampPalette(colors = c("blue", "white", "red"))(n = length(breaks))

pheatmap(cor_res,color = colors_for_plot,breaks = breaks,display_numbers = T,main = "Pearson correlation of pathways") 

6 Correlation all cancer cells

library(ggpubr)
data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM"))
p1 = ggplot(data, aes(x=b2r_signature_TPM, y=GOBP_NOREPINEPHRINE_TRANSPORT_TPM)) +  geom_smooth(method = lm)+
  geom_point()+ stat_cor(method = "pearson")

data = FetchData(object = acc_cancer,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT"))
p2 = ggplot(data, aes(x=b2r_signature, y=GOBP_NOREPINEPHRINE_TRANSPORT)) +  geom_smooth(method = lm)+
  geom_point()+ stat_cor(method = "pearson")

print_tab(p1,title = "TPM")

TPM

geom_smooth() using formula ‘y ~ x’

print_tab(p2,title = "log TPM")

log TPM

geom_smooth() using formula ‘y ~ x’

NA

7 Intra-patient correlation

library(ggpubr)

data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM","patient.ident"))
p1 = ggplot(data,
           aes(x = b2r_signature_TPM, y = GOBP_NOREPINEPHRINE_TRANSPORT_TPM)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson") + facet_wrap(~patient.ident,nrow = 2)

data = FetchData(object = acc_cancer,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT","patient.ident"))
p2 = ggplot(data,
           aes(x = b2r_signature, y = GOBP_NOREPINEPHRINE_TRANSPORT)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson") + facet_wrap(~patient.ident,nrow = 2)
print_tab(p1,title = "TPM")

TPM

geom_smooth() using formula ‘y ~ x’

print_tab(p2,title = "log TPM")

log TPM

geom_smooth() using formula ‘y ~ x’

print_tab(DimPlot(acc_cancer,group.by = "patient.ident"),title = "UMAP")

UMAP

NA

8 Intra-cluster correlation

data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM","seurat_clusters"))

data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM","seurat_clusters"))
p1 = ggplot(data,
           aes(x = b2r_signature_TPM, y = GOBP_NOREPINEPHRINE_TRANSPORT_TPM)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson") + facet_wrap(~seurat_clusters,nrow = 2)

data = FetchData(object = acc_cancer,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT","seurat_clusters"))
p2 = ggplot(data,
           aes(x = b2r_signature, y = GOBP_NOREPINEPHRINE_TRANSPORT)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson") + facet_wrap(~seurat_clusters,nrow = 2)
print_tab(p1,title = "TPM")

TPM

geom_smooth() using formula ‘y ~ x’

print_tab(p2,title = "log TPM")

log TPM

geom_smooth() using formula ‘y ~ x’

print_tab(DimPlot(acc_cancer,group.by = "seurat_clusters"),title = "UMAP")

UMAP

NA

9 Classification

data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM"))
b2r_thresh1 = 500
b2r_thresh2 = 1000

data %<>% mutate(b2r_level = case_when(
  b2r_signature_TPM <b2r_thresh1 ~ "low",
  b2r_signature_TPM >=b2r_thresh1 & b2r_signature_TPM <=b2r_thresh2 ~ "medium",
  b2r_signature_TPM > b2r_thresh2 ~ "high",
)) 

epi_thresh1 = 250
epi_thresh2 = 500
data %<>% mutate(epi_trans_level = case_when(
  GOBP_NOREPINEPHRINE_TRANSPORT_TPM <epi_thresh1 ~ "low",
  GOBP_NOREPINEPHRINE_TRANSPORT_TPM >=epi_thresh1 & GOBP_NOREPINEPHRINE_TRANSPORT_TPM <=epi_thresh2 ~ "medium",
  GOBP_NOREPINEPHRINE_TRANSPORT_TPM > epi_thresh2 ~ "high",
))

b2r low = b2r < 500

b2r medium= b2r >= 500 & <= 1000

b2r high= b2r < 1000

GOBP_NOREPINEPHRINE_TRANSPORT low = GOBP_NOREPINEPHRINE_TRANSPORT< 250

GOBP_NOREPINEPHRINE_TRANSPORTmedium= GOBP_NOREPINEPHRINE_TRANSPORT>= 500 & <= 500

GOBP_NOREPINEPHRINE_TRANSPORT low = GOBP_NOREPINEPHRINE_TRANSPORT< 500

data_to_test = data %>% select(c("b2r_level","epi_trans_level")) %>% filter(b2r_level %in% c("low","high") & epi_trans_level %in%  c("low","high"))
test = fisher.test(table(data_to_test))

library(ggstatsplot)

plt = ggbarstats(
  data_to_test,
  b2r_level,
  epi_trans_level,
  results.subtitle = FALSE,
  subtitle = paste0("Fisher's exact test", ", p-value = ",
                    round(test$p.value, 13))
)
plt

data_for_plot = table(data %>% select(c("b2r_level","epi_trans_level"))) %>% as.data.frame.matrix()
data_for_plot = data_for_plot[c("low","medium","high"),c("low","medium","high")]
colnames(data_for_plot) = c("epi_trans_low","epi_trans_medium","epi_trans_high")
rownames(data_for_plot) = c("b2r_low","b2r_medium","b2r_high")

pheatmap(data_for_plot,cluster_rows = F,cluster_cols = F)

10 Between-clusters correlation

data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM","seurat_clusters"))
average_data1  = data %>% group_by(seurat_clusters) %>%
    dplyr::summarize(b2r_mean = mean(b2r_signature_TPM, na.rm=TRUE))

average_data2  = data %>% group_by(seurat_clusters) %>%
    dplyr::summarize(epi_trans_mean = mean(GOBP_NOREPINEPHRINE_TRANSPORT_TPM, na.rm=TRUE))
average_all = cbind(average_data1,average_data2)
average_all = average_all[,c(-1)]
average_all$seurat_clusters = paste("cluster",average_all$seurat_clusters)
library(ggrepel)
p = ggplot(average_all,
           aes(x = b2r_mean, y = epi_trans_mean, label=seurat_clusters)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson")+geom_text_repel()
print_tab(plt = p,title = "TPM")

TPM

geom_smooth() using formula ‘y ~ x’

NA

data = FetchData(object = acc_cancer,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT","seurat_clusters"))
average_data1  = data %>% group_by(seurat_clusters) %>%
    dplyr::summarize(b2r_mean = mean(b2r_signature, na.rm=TRUE))

average_data2  = data %>% group_by(seurat_clusters) %>%
    dplyr::summarize(epi_trans_mean = mean(GOBP_NOREPINEPHRINE_TRANSPORT, na.rm=TRUE))
average_all = cbind(average_data1,average_data2)
average_all = average_all[,c(-1)]
average_all$seurat_clusters = paste("cluster",average_all$seurat_clusters)
library(ggrepel)
p = ggplot(average_all,
           aes(x = b2r_mean, y = epi_trans_mean, label=seurat_clusters)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson")+geom_text_repel()
print_tab(plt = p,title = "log TPM")

log TPM

geom_smooth() using formula ‘y ~ x’

NA

11 Between-patients correlation

data = FetchData(object = acc_cancer,vars =  c("b2r_signature_TPM","GOBP_NOREPINEPHRINE_TRANSPORT_TPM","patient.ident"))
average_data1  = data %>% group_by(patient.ident) %>%
    dplyr::summarize(b2r_mean = mean(b2r_signature_TPM, na.rm=TRUE))

average_data2  = data %>% group_by(patient.ident) %>%
    dplyr::summarize(epi_trans_mean = mean(GOBP_NOREPINEPHRINE_TRANSPORT_TPM, na.rm=TRUE))
average_all = cbind(average_data1,average_data2)
average_all = average_all[,c(-1)]
library(ggrepel)
p = ggplot(average_all,
           aes(x = b2r_mean, y = epi_trans_mean, label=patient.ident)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson")+geom_text_repel()
print_tab(plt = p,title = "TPM")

TPM

geom_smooth() using formula ‘y ~ x’

data = FetchData(object = acc_cancer,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT","patient.ident"))
average_data1  = data %>% group_by(patient.ident) %>%
    dplyr::summarize(b2r_mean = mean(b2r_signature, na.rm=TRUE))

average_data2  = data %>% group_by(patient.ident) %>%
    dplyr::summarize(epi_trans_mean = mean(GOBP_NOREPINEPHRINE_TRANSPORT, na.rm=TRUE))
average_all = cbind(average_data1,average_data2)
average_all = average_all[,c(-1)]
library(ggrepel)
p = ggplot(average_all,
           aes(x = b2r_mean, y = epi_trans_mean, label=patient.ident)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson")+geom_text_repel()
print_tab(plt = p,title = "log TPM")

log TPM

geom_smooth() using formula ‘y ~ x’

NA

12 Primary cancer cells

# add scores
pathways_scores = FetchData(object = acc_cancer_pri,vars = b2r_genes,slot = "data") %>% 
  mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_cancer_pri  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature_TPM")


pathways_scores = FetchData(object = acc_cancer_pri,vars = b2r_genes,slot = "data") %>% 
  rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
acc_cancer_pri  %<>% AddMetaData(metadata = pathways_scores$score,col.name = "b2r_signature")

for (pathway_name in names(neuronal_pathways)) {
  genes = neuronal_pathways[[pathway_name]]
  pathways_scores = FetchData(object = acc_cancer_pri,vars = genes,slot = "data") %>% 
    mutate_all(~ 2^.)%>% mutate_all(~ .-1) %>%  #covert log(TPM+1) to TPM
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_cancer_pri  %<>% AddMetaData(metadata = pathways_scores$score,col.name = paste0(pathway_name,"_TPM"))
}

for (pathway_name in names(neuronal_pathways)) {
  genes = neuronal_pathways[[pathway_name]]
  pathways_scores = FetchData(object = acc_cancer_pri,vars = genes,slot = "data") %>% 
    rowwise() %>% mutate(score = mean(c_across(everything()))) #mean
  acc_cancer_pri  %<>% AddMetaData(metadata = pathways_scores$score,col.name = pathway_name)
}

data = FetchData(object = acc_cancer_pri,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT","patient.ident"))
average_data1  = data %>% group_by(patient.ident) %>%
    dplyr::summarize(b2r_mean = mean(b2r_signature, na.rm=TRUE))

average_data2  = data %>% group_by(patient.ident) %>%
    dplyr::summarize(epi_trans_mean = mean(GOBP_NOREPINEPHRINE_TRANSPORT, na.rm=TRUE))
average_all = cbind(average_data1,average_data2)
average_all = average_all[,c(-1)]
# average_all$patient.ident = paste("cluster",average_all$patient.ident)
library(ggrepel)
ggplot(average_all,
           aes(x = b2r_mean, y = epi_trans_mean, label=patient.ident)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson")+geom_text_repel()

`geom_smooth()` using formula 'y ~ x'

acc_cancer_pri <- FindNeighbors(acc_cancer_pri, dims = 1:10,verbose = F) %>%  FindClusters(resolution = 0.7) 

library(ggrepel)
data = FetchData(object = acc_cancer_pri,vars =  c("b2r_signature","GOBP_NOREPINEPHRINE_TRANSPORT","seurat_clusters"))
average_data1  = data %>% group_by(seurat_clusters) %>%
    dplyr::summarize(b2r_mean = mean(b2r_signature, na.rm=TRUE))

average_data2  = data %>% group_by(seurat_clusters) %>%
    dplyr::summarize(epi_trans_mean = mean(GOBP_NOREPINEPHRINE_TRANSPORT, na.rm=TRUE))
average_all = cbind(average_data1,average_data2)
average_all = average_all[,c(-1)]
average_all$seurat_clusters = paste("cluster",average_all$seurat_clusters)
ggplot(average_all,
           aes(x = b2r_mean, y = epi_trans_mean, label=seurat_clusters)) +  geom_smooth(method = lm) +
  geom_point() + stat_cor(method = "pearson")+geom_text_repel()

DimPlot(acc_cancer_pri,group.by = "patient.ident")
DimPlot(acc_cancer_pri)