Loop구문 cluster별로 분석

 Idents(seurat_combined1) <- seurat_combined1$seurat_clusters

# 결과 저장용 리스트 초기화

comparison_markers <- list()

comparison_pathways <- list()

comparison_GO <- list()


# 0부터 19까지 모든 클러스터에 대해 반복 작업 수행

for (cluster_id in 0:19) {

  

  cat("Processing: Cluster", cluster_id, "vs others (Clusters excluding", cluster_id, ")\n")

  

  # 기준 클러스터 (예: Cluster 0)와 나머지 클러스터를 "others"로 합침

  reference_cluster <- as.character(cluster_id)

  

  # 나머지 클러스터들을 "others"로 설정 (현재 클러스터 제외)

  seurat_combined1$combined_cluster <- ifelse(seurat_combined1$seurat_clusters == reference_cluster, reference_cluster, "others")

  

  # Idents를 combined_cluster로 설정하여 "others"와 특정 클러스터를 기준으로 분석

  Idents(seurat_combined1) <- "combined_cluster"

  

  # ==========================

  # 마커 분석

  # ==========================

  markers <- FindMarkers(seurat_combined1, ident.1 = reference_cluster, ident.2 = "others", only.pos = TRUE)

  comparison_markers[[paste0("Cluster_", cluster_id, "_vs_others")]] <- markers

  

  # ==========================

  # GSEA pathway 분석 (Hallmark)

  # ==========================

  rank_vector <- markers$avg_log2FC

  names(rank_vector) <- rownames(markers)

  rank_vector <- sort(rank_vector, decreasing = TRUE)  # 내림차순 정렬

  

  fgsea_result <- fgsea(pathways = hallmark_pathways, stats = rank_vector, nperm = 1000)

  fgsea_result$leadingEdge <- sapply(fgsea_result$leadingEdge, function(x) paste(x, collapse = ", "))

  comparison_pathways[[paste0("Cluster_", cluster_id, "_vs_others")]] <- fgsea_result

  

  # ==========================

  # GO 분석 (Biological Process)

  # ==========================

  go_result <- gseGO(

    geneList = rank_vector,

    OrgDb = org.Hs.eg.db,

    ont = "BP",

    keyType = "SYMBOL",

    pvalueCutoff = 0.05,

    verbose = FALSE

  )

  comparison_GO[[paste0("Cluster_", cluster_id, "_vs_others")]] <- go_result@result

  

  # ==========================

  # 결과 파일로 저장

  # ==========================

  write.table(markers, 

              file = paste0("Cluster_", cluster_id, "_vs_others_Markers.txt"), 

              sep = "\t", quote = FALSE, row.names = TRUE)

  

  write.table(fgsea_result, 

              file = paste0("Cluster_", cluster_id, "_vs_others_Pathways.txt"), 

              sep = "\t", quote = FALSE, row.names = TRUE)

  

  write.table(go_result@result, 

              file = paste0("Cluster_", cluster_id, "_vs_others_GO_BP.txt"), 

              sep = "\t", quote = FALSE, row.names = TRUE)

}


cat("All clusters processed successfully!")


댓글

이 블로그의 인기 게시물

#single cell sequencing 기초 분석 - #1 R 설치 및 package 설치

리눅스 기초 #10 GATK calling을 사용하기 위하여, reference file indexing하는 방법

Single cell 분석을 위한 package 소개