Data challeng 1

library(tidyverse)
library(FRmatch)

##########
## data ##
##########

load("dat-challenge1.rda")

#################
## basic plots ##
#################

## cluster size
plot_clusterSize(sce.E1 = sce_query, sce.E2 = sce_ref, name.E1 = "query", name.E2 = "ref")

## original distribution
plot_exprDist(sce.E1 = sce_query, sce.E2 = sce_ref, name.E1 = "query", name.E2 = "ref")

###################
## barcode plots ##
###################

# ## create barcode plot folder
# folder_name <- "barcode_plots"
# if(!dir.exists(folder_name)){dir.create(folder_name)}
# 
# ## self barcode
# ref_clusters <- unique(colData(sce_ref)$cluster_membership)
# for(cl in ref_clusters){
#   plot_cluster_by_markers(sce_ref, cluster.name = cl, nsamp = 30, name.E1 = "",
#                           filename = paste0(folder_name,"/",cl,".pdf"))
# }
# query_clusters <- unique(colData(sce_query)$cluster_membership)
# for(cl in query_clusters){
#   plot_cluster_by_markers(sce_query, cluster.name = cl, nsamp = 30, name.E1 = "",
#                           filename = paste0(folder_name,"/",cl,".pdf"))
# }

###################
## normalization ## FINAL: decided not to do since query and ref data distribution look similar
###################

# ## min-max scaling only
# sce_query_norm <- normalization(sce_query)
# sce_ref_norm <- normalization(sce_ref)
# 
# 
# ## create barcode plot folder
# folder_name <- "barcode_plots_normalization"
# if(!dir.exists(folder_name)){dir.create(folder_name)}
# 
# ## self barcode
# ref_clusters <- unique(colData(sce_ref_norm)$cluster_membership)
# for(cl in ref_clusters){
#   plot_cluster_by_markers(sce_ref_norm, cluster.name = cl, nsamp = 30, name.E1 = "",
#                           filename = paste0(folder_name,"/",cl,".pdf"))
# }
# query_clusters <- unique(colData(sce_query_norm)$cluster_membership)
# for(cl in query_clusters){
#   plot_cluster_by_markers(sce_query_norm, cluster.name = cl, nsamp = 30, name.E1 = "",
#                           filename = paste0(folder_name,"/",cl,".pdf"))
# }

#############
## FRmatch ##
#############

rst_q2r <- FRmatch(sce_query, sce_ref, use.cosine=T) ### FINAL: subsamp.size=20, subsamp.iter=1000 (default)
rst_r2q <- FRmatch(sce_ref, sce_query, use.cosine=T) ### FINAL: subsamp.size=20, subsamp.iter=1000 (default)

plot_FRmatch(rst_q2r)

plot_FRmatch(rst_q2r, type="padj")

plot_FRmatch(rst_r2q)

plot_FRmatch(rst_r2q, type="padj")

out_FRmatch <- plot_bi_FRmatch(rst_q2r, rst_r2q, name.E1="", name.E2="", return.value = T) #sig.level = 0.05

plot_bi_FRmatch(rst_r2q, rst_q2r, name.E1="", name.E2="")

## output required lists
out_FRmatch <- out_FRmatch[-nrow(out_FRmatch),] #remove the unassigned row

## unassigned
(Q_unassigned <- names(which(colSums(out_FRmatch)==0)))

[1] "Q_cl_9"  "Q_cl_20" "Q_cl_50" "Q_cl_85" "Q_cl_22" "Q_cl_10" "Q_cl_8" 

(R_unassigned <- names(which(rowSums(out_FRmatch)==0)))

[1] "R_cl_21" "R_cl_40" "R_cl_15" "R_cl_18" "R_cl_6"  "R_cl_24" "R_cl_9" 

## one-to-many
(Q_multi <- names(which(colSums(out_FRmatch)>2)))

[1] "Q_cl_33"

(R_multi <- names(which(rowSums(out_FRmatch)>2)))

[1] "R_cl_11" "R_cl_32" "R_cl_84"

## one-to-one
matches <- get_values(out_FRmatch)
matches_1to1 <- matches %>% filter(!(rowname %in% R_multi)) %>% select(-value)
matches_1to1

## check if one-to-one
# table(matches_1to1$rowname)
# table(matches_1to1$colname)

## write output
write.csv(Q_unassigned, file="Q_unassigned.csv")
write.csv(R_unassigned, file="R_unassigned.csv")
write.csv(Q_multi, file="Q_multi.csv")
write.csv(R_multi, file="R_multi.csv")
write.csv(matches_1to1, file="matches_1to1.csv")