Add CD-HIT parsing function and update DB logic

Added a function to parse CD-HIT .clstr output into a long-format mapping of clusters to member feature ids. Updated database writing logic to include the new protein members table.

Author: AbhirupaGhosh

R/data_processing.R

@@ -779,6 +779,54 @@ runPanaroo2Duckdb <- function(duckdb_path,
   cluster_map
 }
 
+#' Parse CD-HIT `.clstr` output into a long-format mapping
+#'
+#' Reads a CD-HIT `.clstr` file and constructs a mapping of clusters to member feature ids.
+#'
+#' @param clustered_faa Base path to CD-HIT output (without `.clstr` extension).
+#'
+#' @return A tibble with columns `cluster` and `member`.
+#'
+#' @keywords internal
+.extractMembersInClusters <- function(clustered_faa) {
+  clstr <- paste0(clustered_faa, ".clstr")
+  if (!file.exists(clstr)) {
+    stop("CD-HIT cluster file not found: ", clstr,
+         "\nEnsure .runCDHIT() completed successfully and produced the .clstr file.")
+  }
+
+  lines <- data.table::fread(clstr, sep = "\n", header = FALSE)$V1
+  cluster_ids <- grep("^>Cluster", lines)
+  cluster_member <- data.table::data.table()
+
+  for (i in seq_along(cluster_ids)) {
+    start <- cluster_ids[i] + 1
+    end   <- if (i < length(cluster_ids)) cluster_ids[i + 1] - 1 else length(lines)
+    cluster_lines <- lines[start:end]
+
+    # This finds the reference cluster ID and names the cluster with it
+    ref_line <- grep("\\*$", cluster_lines, value = TRUE)
+    ref_id <- if (length(ref_line) > 0) {
+      stringr::str_extract(ref_line, "fig\\|[0-9]+\\.[0-9]+\\.peg(?:sc)?\\.[0-9]+")
+    } else {
+      paste0("Cluster_", i - 1)
+    }
+
+    # Pull genome IDs
+    members <- stringr::str_match(cluster_lines,
+                                         "fig\\|([0-9]+\\.[0-9]+)\\.peg(?:sc)?\\.[0-9]+")[, 1]
+    members <- members[!is.na(members)]
+
+    if (length(members) > 0) {
+      cluster_member <- data.table::rbindlist(list(
+        cluster_member,
+        data.table::data.table(cluster = ref_id, member = members)
+      ), use.names = TRUE)
+    }
+  }
+
+  tibble::as_tibble(cluster_member)
+}  
 
 #' Build genome-by-protein-cluster count matrix
 #'
@@ -877,6 +925,10 @@ CDHIT2duckdb <- function(duckdb_path,
     ),
     overwrite = TRUE
   )
+
+  cluster_member <- .extractMembersInClusters(cdhit_outputs$clustered_faa)
+  DBI::dbWriteTable(con, "protein_members", cluster_member, overwrite = TRUE)
+  
   invisible(TRUE)
 }
 
@@ -1433,6 +1485,7 @@ cleanData <- function(duckdb_path, path) {
   protein_names_parquet <- file.path(path, "protein_names.parquet")
 
   protein_cluster_seq_parquet <- file.path(path, "protein_seqs.parquet")
+  protein_cluster_member_parquet <- file.path(path, "protein_members.parquet")
 
   writeCompressedParquet <- function(df, path) {
     arrow::write_parquet(
@@ -1502,6 +1555,9 @@ cleanData <- function(duckdb_path, path) {
   DBI::dbReadTable(con, "protein_cluster_seq") |> writeCompressedParquet(protein_cluster_seq_parquet)
   DBI::dbExecute(con_new, sprintf("CREATE OR REPLACE VIEW protein_seqs AS SELECT * FROM read_parquet('%s')", protein_cluster_seq_parquet))
 
+  DBI::dbReadTable(con, "protein_members") |> writeCompressedParquet(protein_cluster_member_parquet)
+  DBI::dbExecute(con_new, sprintf("CREATE OR REPLACE VIEW protein_members AS SELECT * FROM read_parquet('%s')", protein_cluster_member_parquet))
+  
   DBI::dbReadTable(con, "genome_gene_protein") |> writeCompressedParquet(genome_gene_protein_parquet)
   DBI::dbExecute(con_new, sprintf("CREATE OR REPLACE VIEW genome_gene_protein AS SELECT * FROM read_parquet('%s')", genome_gene_protein_parquet))