--- a
+++ b/FigS5C.R
@@ -0,0 +1,63 @@
+GIT_HOME="/research/users/ppolonen/git_home/ImmunogenomicLandscape-BloodCancers/"
+source(file.path(GIT_HOME, "common_scripts/visualisation/plotting_functions.R"))
+source(file.path(GIT_HOME, "common_scripts/statistics/functions_statistics.R"))
+source(file.path(GIT_HOME, "common_scripts/statistics/statistics_wrappers.R"))
+source(file.path(GIT_HOME, "common_scripts/pathway_analysis/functions.GSEA.R"))
+source(file.path(GIT_HOME, "common_scripts/scRNA/functions.scRNA.analysis.R"))
+source(file.path(GIT_HOME, "common_scripts/statistics/useful_functions.R"))
+
+library(Matrix)
+library(Seurat)
+library(data.table)
+library(ComplexHeatmap)
+library(circlize)
+library(parallel)
+library(ggplot2)
+
+setwd("/research/groups/sysgen/PROJECTS/HEMAP_IMMUNOLOGY/petri_work/HEMAP_IMMUNOLOGY/Published_data_figures")
+
+
+# analyze FIMM and Galen AML and find costim associations to certain clusters.
+co.stim=data.table::fread("costim_ligands_final.txt", data.table = F)[,c(1,3,5)]
+co.stimR=data.table::fread("costim_ligands_final.txt", data.table = F)
+co.stimR.f=c(unlist(strsplit(co.stimR$`Receptor gene`, ", ")), "LAG3")
+co.stim=rbind(co.stim, c("FGL1", "LAG3", "Inhibitory"))
+
+load("AML_Galen_scRNA.Rdata")
+galen=scmat
+
+load("HCA_scRNA.Rdata")
+HCA=scmat
+
+load("FIMM_AML_scRNA.Rdata")
+
+# plot from each dataset all significant genes:
+subtype.order=c("MDS-like", "Progenitor-like", "Monocyte-like", "Monocyte-like-MLL", "CEBPA", "RUNX1-RUNX1T1", "CBFB-MYH11", "PML-RARA")
+
+fimm=scmat[,scmat[["SingleR.label"]][,1]%in%c("HSC", "MPP", "GMP", "CMP", "MEP", "Monocytes", "Erythrocytes")]
+galen.sub=galen[,galen[["SingleR.label"]][,1]%in%c("HSC", "MPP", "GMP", "CMP", "MEP", "Monocytes", "Erythrocytes")]
+hca=HCA[,HCA[["SingleR.label"]][,1]%in%c("HSC", "MPP", "GMP", "CMP", "MEP", "Monocytes", "Erythrocytes")]
+
+Idents(fimm)=fimm[["SingleR.label"]][,1]
+Idents(galen.sub)=factor(galen.sub[["SingleR.label"]][,1], levels=c("HSC", "MPP", "GMP", "CMP", "MEP", "Monocytes", "Erythrocytes"))
+Idents(hca)=hca[["SingleR.label"]][,1]
+
+# do DE gene tests:
+DE.FIMM=FindAllMarkers(object = fimm, features = genelist[genelist%in%rownames(scmat)], only.pos = T, logfc.threshold = 0.15)
+DE.galen=FindAllMarkers(object = galen.sub, features = genelist[genelist%in%rownames(galen)], only.pos = T, logfc.threshold = 0.15)
+DE.HCA=FindAllMarkers(object = hca, features = genelist[genelist%in%rownames(HCA)], only.pos = T, logfc.threshold = 0.15)
+
+all=rbind(DE.FIMM, DE.galen, DE.HCA)
+all=all[order(match(all$cluster, c("HSC", "MPP", "GMP", "CMP", "MEP", "Monocytes", "Erythrocytes"))),]
+
+a=table(all$cluster, all$gene)
+genelist.filt=c(colnames(a[,apply(a, 2, function(v)any(v>1))]), "C10orf54")
+
+# these were selected for the figure:
+genelist.filt=c("CD34", "CLEC2B", "TNFRSF14", "CD84", "VSIR", "C10orf54", "CD68", "CD48","CD86","ENTPD1") 
+
+pdf("FigS5C_VISTA_FIMM_dotplot.pdf", width = 4.5, height = 2.5)
+plot.DotPlot(fimm, features = genelist.filt[genelist.filt%in%rownames(scmat)], cols=c("grey75", "red"), scale.max = 50, scale.min = 10, dot.scale = 4)
+plot.DotPlot(galen.sub, features = genelist.filt[genelist.filt%in%rownames(galen)], cols=c("grey75", "red"), scale.max = 50, scale.min = 10, dot.scale = 4)
+plot.DotPlot(hca, features = genelist.filt[genelist.filt%in%rownames(HCA)], cols=c("grey75", "red"), scale.max = 50, scale.min = 10, dot.scale = 4)
+dev.off()