df_pop 

mat_complete <- rbind(
  df_facs %>% 
    left_join(., df_meta %>% select(PatientID, `CITEseq`)) %>% 
    filter(`CITEseq`=="-") %>% 
    select(PatientID, Population, Prop=FACS),
  df_freq %>% 
    mutate(RNA=ifelse(is.nan(RNA), 0, RNA)) %>% 
    select(PatientID, Population, Prop=RNA)) %>% 
  filter(Population %in% df_pop$Population) %>% 
  left_join(., df_pop) %>% 
  select(-Population) %>% 
  pivot_wider(names_from = "IdentI", values_from = "Prop", values_fill = 0) %>% 
  column_to_rownames("PatientID") 

cl_order <- c(6,1,2,9,8,13,15,11,12,16,3,5,14,19)
names(labels_cl_parsed) <- as.character(cluster_order)

df_freqPlot <- 
  mat_complete %>% 
  rownames_to_column("PatientID") %>% 
  pivot_longer(cols=2:ncol(.), names_to = "IdentI", values_to = "Prop") %>% 
  add_entity() %>% 
  mutate(Entity=factor(Entity, levels = c("rLN", "DLBCL", "MCL", "FL", "MZL"))) %>% 
  mutate(IdentI=factor(IdentI, levels=cl_order)) %>% 
  mutate(Label=factor(IdentI, levels=cl_order, labels = labels_cl_parsed[as.character(cl_order)])) %>% 
  group_by(Entity, IdentI) %>% 
  mutate(outlier = (Prop > quantile(Prop, 0.75) + IQR(Prop) * 1.5) | (Prop < quantile(Prop, 0.25) - IQR(Prop) * 1.5))

df_medianLines <- df_freqPlot %>% 
  filter(Entity=="rLN") %>% 
  group_by(IdentI, Label) %>% 
  summarise(MedianProp=median(Prop))

df_freqPlot_pvalues <- 
  df_freqPlot %>% 
  group_by(IdentI) %>% 
  wilcox_test(data=., formula = Prop ~ Entity, detailed = T, ref.group = "rLN") %>% 
  adjust_pvalue(method = "BH") %>% 
  select(IdentI, Entity=group2, p.adj, estimate) %>% 
  mutate(Entity=factor(Entity, levels = c("rLN", "DLBCL", "MCL", "FL", "MZL"))) %>% 
  mutate(p.adj=ifelse(p.adj>0.05, NA, p.adj)) %>% 
  mutate(p.adj_s=format(p.adj, scientific = TRUE, digits=1)) %>% 
  mutate(p.adj_f=case_when(p.adj > 0.05 ~ "NA",
                           p.adj==0.05 ~ "0.05",
                           p.adj < 0.05 & p.adj > 0.001 ~ as.character(round(p.adj, 3)),
                           p.adj==0.001 ~ "0.001",
                           p.adj < 0.001 ~ p.adj_s)) %>% 
  filter(!is.na(p.adj)) %>% 
  left_join(., df_freqPlot %>% select(IdentI,  Label) %>% distinct) %>% 
  left_join(., data.frame(IdentI=factor(cl_order), height=c(52.5, 52.5, 21.5, 21.5, 22, 22, 20, 20, 26, 26, 80, 80, 21, 21)))

p <- list()

for(i in c(1:7)){
  
  y <- list(c(1,6),c(2,9),c(8,13),c(15,11),c(12,16),c(3,5),c(14,19))[[i]]
  ylim <- c(70,30,32,28,35,110,28)
  
  p[[i]] <- 
    ggplot(data=df_freqPlot %>% filter(IdentI %in% y) %>% 
             mutate(Label=factor(Label, levels = labels_cl_parsed[as.character(y)])), 
             aes(y=Prop, x=Entity, fill=IdentI))+
    geom_hline(data=df_medianLines %>%filter(IdentI %in% y), aes(yintercept=MedianProp),
               size=0.25, linetype="dashed", color="grey60")+
    geom_boxplot(width=0.4, outlier.shape = 21, outlier.size = 1, outlier.color = "white", 
                 outlier.alpha = 0, show.legend = F, size=0.25)+
    ggbeeswarm::geom_beeswarm(data = function(x) dplyr::filter_(x, ~ outlier), cex = 3, stroke=0.25, 
                              groupOnX = TRUE, shape = 21, size = 1, color = "white", alpha = 1)+
    geom_text(data=df_freqPlot_pvalues %>% filter(IdentI %in% y), 
              inherit.aes = F, aes(y=height, x=Entity, label=p.adj_f), hjust=0.1, size=2.3, angle=45)+
    scale_fill_manual(values = colors_umap_cl)+
    scale_y_continuous(name="% of total T-cells", limits=c(0,ylim[i]))+
    scale_x_discrete(expand = c(0.17,0.17))+
    facet_wrap(~Label, ncol = 2, labeller = label_parsed)+
    mytheme_1+
    theme(axis.title.x = element_blank(),
          strip.background = element_rect(color=NA),
          plot.margin = unit(c(0,0,0,0), units = "cm"),
          plot.title = element_text(margin = unit(c(0,0,0,0), units = "cm")),
          panel.border = element_rect(size = 0.5),
          axis.text.x = element_text(angle=45, hjust = 1))
  
  if(i!=7){
    p[[i]] <- p[[i]]+
      theme(axis.text.x = element_blank(),
            axis.ticks.x = element_blank())
  }
  
  if(i!=4){
    p[[i]] <- p[[i]]+
      theme(axis.title.y = element_blank())
    }
  
  if(i==1){
    p[[i]] <- p[[i]]+
      labs(tag = "A")+
      theme(plot.tag = element_text(vjust = -0.5))
}
  
  
}

plot_freq <- wrap_plots(p, ncol = 1)
plot_freq

ggsave(width = 9, height = 21, units = "cm", filename = "Figure3_p1.pdf")

# Create data frame
data <- data.frame(
  level1="all",
  level2=c("rLN", 
           "MZL",
           "MCL",
           "FL",
           "DLBCL")
)

edges_level1_2 <- data %>% select(level1, level2) %>% unique %>% rename(from=level1, to=level2)
edge_list=rbind(edges_level1_2)
vert <- data.frame(
  name=unique(c(data$level1, data$level2))) %>% 
  mutate(label=c(NA, "rLN", "MZL", "MCL", "FL", "DLBCL"))

# Create graph object
mygraph_lasso <- graph_from_data_frame( edge_list ,vertices = vert)

# Plot dendrogramm
ggraph(mygraph_lasso, layout = 'tree', circular = FALSE)+ 
  geom_edge_diagonal(strength = 1.4, edge_width=0.25)+
  geom_node_point(shape=21, size=3.5, color="white", stroke=2, alpha=c(0,1,1,1,1,1),
                  fill=c(NA, brewer.pal(name = "Paired", 5)[c(5,4,3,2,1)]))+
  coord_flip(clip = "off")+
  scale_y_reverse()+
  theme_void()+
  theme(legend.position = "right",
        plot.margin = margin(0.25,0.25,0.25,0, unit = "cm"),  
        plot.title = element_text(hjust=0.4, size=7, face = "bold"),)

#ggsave(width = 3.5, height = 3.4, units = "cm", filename = "Figure3_p3.pdf")

total <- mat_complete %>% 
  rownames_to_column("PatientID") %>% 
  left_join(., df_meta %>% select(PatientID, Tcells)) %>% 
  add_entity() %>% 
  mutate_if(is.numeric, ~./100)

cell_types <- total %>% select(-Entity, -PatientID) %>% colnames()
gt <- my_glmnet(total)

entities <- c("DLBCL", "MCL", "FL", "MZL", "rLN")

tbl <- gt$confusion_table
class(tbl) = "matrix"
tbl = tbl / rowSums(tbl) # convert to probability estimates
tbl = tbl[entities, entities]

tbl %>% data.frame() %>% 
  rownames_to_column("truth") %>% 
  pivot_longer(cols = 2:ncol(.), names_to = "predicted", values_to = "Prop") %>% 
  ggplot(aes(x=truth, y=predicted, fill=Prop))+
  geom_tile()+
  scale_x_discrete(limits=rev(entities), expand = c(0,0), position = "top")+
  scale_y_discrete(limits=entities, expand = c(0,0), name="Truth")+
  geom_hline(yintercept = c(1.5,2.5,3.5,4.5),size=0.25, color="black")+
  geom_vline(xintercept = c(1.5,2.5,3.5,4.5),size=0.25, color="black")+
  scale_fill_gradientn(colours = colorRampPalette(RColorBrewer::brewer.pal(9, "BuPu"))(100), name="Prop", limits=c(0,0.8))+
  xlab("Reference")+
  mytheme_1+
  coord_fixed()+
  theme(panel.border = element_blank(),
        legend.position = "right",
        legend.key.height = unit(0.3, "cm"),
        legend.key.width = unit(0.3, "cm"),
        legend.box.spacing = unit(0.1, "cm"),
        legend.box.margin = unit(c(0,-0.25,0,0.05), units = "cm"),
        plot.margin = unit(c(0.1,0.1,0.1,0.1), units = "cm"),
        axis.text.x = element_text(angle=45, hjust = 0),
        axis.ticks = element_blank(),
        axis.title.y = element_blank())

#ggsave(width = 5.5, height = 5.5, units = "cm", filename = "Figure3_p4.pdf")