Overlay_projection.R

library(Seurat)
library(dplyr)
library(ggplot2)


# Run SCTransform on the reference to match the query's assay type
# This creates an 'SCT' assay in your obj_10x
obj_10x <- SCTransform(obj_10x, verbose = FALSE)

# Re-run PCA on the new SCT assay
obj_10x <- RunPCA(obj_10x, verbose = FALSE)

# Ensure UMAP is still available (using the author's coordinates you injected previously)
# If you haven't run this yet, do it now:
obj_10x <- RunUMAP(obj_10x, dims = 1:30, reduction = "pca", return.model = TRUE)

# --- PART A: THE MATH (SNAP) ---
# 1. Find Anchors
# We change normalization.method to "SCT" to match your obj_smart
anchors <- FindTransferAnchors(
  reference = obj_10x, 
  query = obj_smart, 
  dims = 1:30, 
  reference.reduction = "pca",
  normalization.method = "SCT"
)

# 2. Map Query
# This 'snaps' your query cells into the reference UMAP space
obj_smart <- MapQuery(
  anchorset = anchors, 
  reference = obj_10x, 
  query = obj_smart,
  refdata = list(celltype = "celltype"), 
  reference.reduction = "pca", 
  reduction.model = "umap"
)

# --- PART B: THE VISUALIZATION ---
# 1. Grab coordinates
ref_emb   <- Embeddings(obj_10x, "umap")
query_emb <- Embeddings(obj_smart, "ref.umap") # This was created by MapQuery

colnames(ref_emb)[1:2]   <- c("UMAP_1","UMAP_2")
colnames(query_emb)[1:2] <- c("UMAP_1","UMAP_2")

# 2. Build plotting data frames
ref_df <- data.frame(ref_emb[, 1:2], celltype = obj_10x$celltype, source = "Reference")
query_df <- data.frame(query_emb[, 1:2], celltype = obj_smart$predicted.celltype, source = "Query")

all_coords <- rbind(ref_df, query_df)

# 3. Final Plot with Large Legend
ggplot(all_coords, aes(UMAP_1, UMAP_2)) +
  geom_point(data = subset(all_coords, source == "Reference"),
             aes(color = celltype), size = 0.5, alpha = 0.3) +
  geom_point(data = subset(all_coords, source == "Query"),
             color = "red", size = 2, alpha = 0.9) +
  theme_classic() +
  ggtitle("Query Projection onto Liu et al. (2021) Reference Map") +
  # TWEAK: Override legend size so you can see the colors clearly
  guides(color = guide_legend(override.aes = list(size = 5, alpha = 1))) +
  theme(legend.text = element_text(size = 12),
        legend.title = element_text(face = "bold", size = 14))

#===output===

# 1. Extract the relevant data from your query object (obj_smart)
# We pull the cell barcode (row names), the predicted lineage, and the UMAP coordinates
output_data <- data.frame(
  cell_barcode = rownames(obj_smart@meta.data),
  predicted_lineage = obj_smart$predicted.celltype,
  # We also include the prediction score (how confident the model is)
  prediction_score = obj_smart$predicted.celltype.score,
  # The UMAP coordinates on the Liu et al. reference map
  UMAP_1 = Embeddings(obj_smart, "ref.umap")[,1],
  UMAP_2 = Embeddings(obj_smart, "ref.umap")[,2],
  # Include your original condition (TE, PE, IM) for comparison
  original_condition = obj_smart$cluster_factor,
  stringsAsFactors = FALSE
)

# 2. Save the results to a CSV file
write.csv(output_data, "query_cells_mapping_to_Liu_reference.csv", row.names = FALSE)

# 3. Print a summary to the console to verify
print("Mapping Summary:")
table(obj_smart$predicted.celltype)