diff --git a/NAMESPACE b/NAMESPACE index e43ed5d..0823397 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -3,7 +3,6 @@ export(.extractAMRtable) export(.updateBVBRCdata) export(CDHIT2duckdb) -export(buildClusterFeatureMap) export(buildDyadFeatureMap) export(checkDataAvailability) export(cleanData) diff --git a/R/data_processing.R b/R/data_processing.R index 4825c16..9f58815 100644 --- a/R/data_processing.R +++ b/R/data_processing.R @@ -446,6 +446,7 @@ NULL tidyr::separate_rows(protein_ids, sep = ";") |> # dplyr::filter(!stringr::str_detect(protein_ids, "_pseudo")) |> dplyr::mutate(protein_ids = gsub("_pseudo", "", protein_ids)) |> + dplyr::mutate(protein_ids = gsub("_len", "", protein_ids)) |> DBI::dbWriteTable(conn = con, name = "genome_gene_protein", overwrite = TRUE) } @@ -1405,7 +1406,7 @@ CDHIT2duckdb <- function(duckdb_path, #' Default: `8`. #' @param n_workers Integer. Number of parallel HMMER jobs to run. Default: `8`. #' @param verbose Logical. Print progress messages. Default: `TRUE`.#' -#' @returns +#' @returns Invisibily returns completed HMMER Parquet files per requested database. #' #' @keywords internal .runHMMER <- function(duckdb_path, @@ -2235,7 +2236,7 @@ CDHIT2duckdb <- function(duckdb_path, #' @param ref_file_path Directory containing reference TSVs used by #' [cleanMetaData()] and [cleanData()] for metadata harmonization. #' Default: `"data_raw/"`. -#' +#' #' @export cleanMetaData <- function(duckdb_path, path, ref_file_path = "data_raw/") { duckdb_path <- normalizePath(duckdb_path) diff --git a/R/feature_to_cluster.R b/R/feature_to_cluster.R deleted file mode 100644 index 494f52e..0000000 --- a/R/feature_to_cluster.R +++ /dev/null @@ -1,171 +0,0 @@ -#' Build a protein cluster-to-feature mapping using DuckDB -#' -#' This function constructs a mapping between protein clusters and functional -#' features (e.g., gene families, domains, COGs, and ARGs) using a DuckDB-backed -#' workflow. The implementation is SQL-first and memory-efficient, leveraging -#' DuckDB views and Parquet output. -#' -#' @param duckdb_parquet_path Character. Path to the DuckDB database file with parquet file views -#' containing the input tables. -#' @param output_path Character or NULL. Directory where the output Parquet file -#' (\code{cluster_feature.parquet}) will be written. If NULL, the output is written -#' alongside the DuckDB database. -#' -#' @details -#' The function performs the following steps: -#' \itemize{ -#' \item Creates views for each feature type: -#' \itemize{ -#' \item Gene → protein features -#' \item Domain annotations -#' \item COG annotations -#' \item Antibiotic resistance gene (ARG) annotations -#' } -#' \item Combines all feature mappings into a unified protein–feature view -#' \item Joins protein–feature mappings to cluster membership -#' \item Writes the resulting cluster–feature mapping to a compressed Parquet file -#' } -#' -#' All joins and transformations are executed inside DuckDB, ensuring scalability -#' for large datasets without loading data into R memory. -#' -#' @return Invisibly returns the file path to the generated Parquet file. -#' -#' @examples -#' \dontrun{ -#' buildClusterFeatureMap( -#' duckdb_parquet_path = "data/Csp_parquet.duckdb", -#' output_path = NULL -#' ) -#' } -#' -#' @import DBI duckdb -#' @export -buildClusterFeatureMap <- function( - duckdb_parquet_path, - output_path = NULL -) { - con <- DBI::dbConnect( - duckdb::duckdb(), - normalizePath(duckdb_parquet_path) - ) - on.exit(DBI::dbDisconnect(con, shutdown = TRUE), add = TRUE) - - out_dir <- if (is.null(output_path)) { - dirname(duckdb_parquet_path) - } else { - normalizePath(output_path) - } - parquet_path <- file.path(out_dir, "cluster_feature.parquet") - - # ========================= - # Gene → protein features - # ========================= - DBI::dbExecute(con, " - CREATE OR REPLACE VIEW v_gene AS - SELECT DISTINCT - protein_ids AS protein_id, - REPLACE(Gene, '~', '.') AS feature - FROM genome_gene_protein - WHERE protein_ids IS NOT NULL - AND Gene IS NOT NULL - ") - - # ========================= - # Domain features - # ========================= - DBI::dbExecute(con, " - CREATE OR REPLACE VIEW v_domain AS - SELECT DISTINCT - AccNum AS protein_id, - \"DB.ID\" AS feature - FROM domain_names - WHERE AccNum IS NOT NULL - AND \"DB.ID\" IS NOT NULL - ") - - # ========================= - # Structural gene features - # (equivalent to separate_rows + inner_join(gp)) - # ========================= - # DBI::dbExecute(con, " - # CREATE OR REPLACE VIEW v_struct AS - # SELECT DISTINCT - # gp.protein_ids AS protein_id, - # s_gene AS feature - # FROM struct s - # JOIN genome_gene_protein gp - # ON gp.genome_ids = s.genome_id - # CROSS JOIN UNNEST(string_split(s.struct, '.')) AS t(s_gene) - # WHERE s.value = 1 - # ") - - # ========================= - # COG features - # ========================= - DBI::dbExecute(con, " - CREATE OR REPLACE VIEW v_cog AS - SELECT DISTINCT - query_name AS protein_id, - name AS feature - FROM protein_COG - WHERE query_name IS NOT NULL - AND name IS NOT NULL - ") - - # ========================= - # ARG features - # ========================= - DBI::dbExecute(con, " - CREATE OR REPLACE VIEW v_arg AS - SELECT DISTINCT - query_name AS protein_id, - REPLACE(REPLACE(name, '-NCBIFAM', ''), '-', '.') AS feature - FROM protein_ResFinder - WHERE query_name IS NOT NULL - AND name IS NOT NULL - ") - - # ========================= - # Union: protein → feature - # ========================= - DBI::dbExecute(con, " - CREATE OR REPLACE VIEW v_protein_feature AS - SELECT protein_id, feature FROM v_gene - UNION - SELECT protein_id, feature FROM v_domain - UNION - SELECT protein_id, feature FROM v_cog - UNION - SELECT protein_id, feature FROM v_arg - ") - - # ========================= - # Cluster → feature mapping - # ========================= - DBI::dbExecute(con, " - CREATE OR REPLACE VIEW cluster_feature AS - SELECT DISTINCT - cm.cluster, - pf.feature - FROM protein_members cm - JOIN v_protein_feature pf - ON pf.protein_id = cm.member - WHERE cm.cluster IS NOT NULL - AND pf.feature IS NOT NULL - ") - - # ========================= - # Write Parquet from DuckDB - # ========================= - DBI::dbExecute( - con, - sprintf( - "COPY cluster_feature TO '%s' - (FORMAT PARQUET, COMPRESSION ZSTD)", - parquet_path - ) - ) - - invisible(parquet_path) -} diff --git a/R/feature_to_head.R b/R/feature_to_head.R index 17bead1..9486096 100644 --- a/R/feature_to_head.R +++ b/R/feature_to_head.R @@ -411,7 +411,7 @@ buildDyadFeatureMap <- function( view_name = "v_pfam", parquet_dir = parquet_dir, dataset_name = "protein_Pfam", - feature_expr = "query_name" + feature_expr = "REPLACE(query_name, '-', '.')" ) } @@ -444,7 +444,7 @@ buildDyadFeatureMap <- function( view_name = "v_defensecas", parquet_dir = parquet_dir, dataset_name = "protein_DefenseCas", - feature_expr = "query_name" + feature_expr = "REPLACE(query_name, '-', '.')" ) } diff --git a/man/buildClusterFeatureMap.Rd b/man/buildClusterFeatureMap.Rd deleted file mode 100644 index 7b56363..0000000 --- a/man/buildClusterFeatureMap.Rd +++ /dev/null @@ -1,52 +0,0 @@ -% Generated by roxygen2: do not edit by hand -% Please edit documentation in R/feature_to_cluster.R -\name{buildClusterFeatureMap} -\alias{buildClusterFeatureMap} -\title{Build a protein cluster-to-feature mapping using DuckDB} -\usage{ -buildClusterFeatureMap(duckdb_parquet_path, output_path = NULL) -} -\arguments{ -\item{duckdb_parquet_path}{Character. Path to the DuckDB database file with parquet file views -containing the input tables.} - -\item{output_path}{Character or NULL. Directory where the output Parquet file -(\code{cluster_feature.parquet}) will be written. If NULL, the output is written -alongside the DuckDB database.} -} -\value{ -Invisibly returns the file path to the generated Parquet file. -} -\description{ -This function constructs a mapping between protein clusters and functional -features (e.g., gene families, domains, COGs, and ARGs) using a DuckDB-backed -workflow. The implementation is SQL-first and memory-efficient, leveraging -DuckDB views and Parquet output. -} -\details{ -The function performs the following steps: -\itemize{ -\item Creates views for each feature type: -\itemize{ -\item Gene → protein features -\item Domain annotations -\item COG annotations -\item Antibiotic resistance gene (ARG) annotations -} -\item Combines all feature mappings into a unified protein–feature view -\item Joins protein–feature mappings to cluster membership -\item Writes the resulting cluster–feature mapping to a compressed Parquet file -} - -All joins and transformations are executed inside DuckDB, ensuring scalability -for large datasets without loading data into R memory. -} -\examples{ -\dontrun{ -buildClusterFeatureMap( - duckdb_parquet_path = "data/Csp_parquet.duckdb", - output_path = NULL -) -} - -} diff --git a/man/class_abbr.Rd b/man/class_abbr.Rd index 7863770..dda73de 100644 --- a/man/class_abbr.Rd +++ b/man/class_abbr.Rd @@ -11,7 +11,7 @@ A data frame Internal reference data } \usage{ -class_abbr +data(class_abbr) } \description{ A dataset mapping drug classes to their abbreviations. diff --git a/man/clean_drug.Rd b/man/clean_drug.Rd index fc4b485..809676c 100644 --- a/man/clean_drug.Rd +++ b/man/clean_drug.Rd @@ -11,7 +11,7 @@ A data frame Internal reference data } \usage{ -clean_drug +data(clean_drug) } \description{ A dataset mapping original drug names to cleaned/standardized names. diff --git a/man/cleaned_bvbrc_countries.Rd b/man/cleaned_bvbrc_countries.Rd index 4edf9aa..dfda707 100644 --- a/man/cleaned_bvbrc_countries.Rd +++ b/man/cleaned_bvbrc_countries.Rd @@ -11,7 +11,7 @@ A data frame Internal reference data } \usage{ -cleaned_bvbrc_countries +data(cleaned_bvbrc_countries) } \description{ A dataset mapping raw country entries to cleaned and standardized names. diff --git a/man/dot-runHMMER.Rd b/man/dot-runHMMER.Rd index d693531..e950621 100644 --- a/man/dot-runHMMER.Rd +++ b/man/dot-runHMMER.Rd @@ -42,6 +42,9 @@ Default: \code{8}.} \item{verbose}{Logical. Print progress messages. Default: \code{TRUE}.#'} } +\value{ +Invisibily returns completed HMMER Parquet files per requested database. +} \description{ Wrapper for preparing HMM databases and running HMMER on protein sequences from duckdb and writing them. } diff --git a/man/drug_abbr.Rd b/man/drug_abbr.Rd index 39542f2..6912d4f 100644 --- a/man/drug_abbr.Rd +++ b/man/drug_abbr.Rd @@ -11,7 +11,7 @@ A data frame Internal reference data } \usage{ -drug_abbr +data(drug_abbr) } \description{ A dataset mapping drug names to their abbreviations. diff --git a/man/drug_class.Rd b/man/drug_class.Rd index 1cf412a..89e354a 100644 --- a/man/drug_class.Rd +++ b/man/drug_class.Rd @@ -11,7 +11,7 @@ A data frame Internal reference data } \usage{ -drug_class +data(drug_class) } \description{ A dataset mapping drugs to their drug classes. diff --git a/vignettes/BVBRC_stats.Rmd b/vignettes/BVBRC_stats.Rmd index cd6e380..b8c27ea 100644 --- a/vignettes/BVBRC_stats.Rmd +++ b/vignettes/BVBRC_stats.Rmd @@ -298,6 +298,27 @@ clean_all_char_cols <- function(df) { )) } +# Clean ALL List columns in the data frame +clean_list_col <- function(x) { + purrr::map( + x, + \(z) { + if (is.null(z) || length(z) == 0) + return(NA_character_) + + z <- as.character(unlist(z)) + z <- stringr::str_squish(z) + + z[is_placeholder_na(z)] <- NA_character_ + + if (all(is.na(z))) + return(NA_character_) + + z + } + ) +} + # --- Block-level summaries --------------------------------------------------- summarize_block <- function(df, cols, block_name) { if (length(cols) == 0) { @@ -337,6 +358,40 @@ bvbrc <- fetchCompleteBVBRCMetadataAPI() # Apply to your table bvbrc_clean <- clean_all_char_cols(bvbrc) +list_cols <- names(bvbrc_clean)[sapply(bvbrc_clean, is.list)] +bvbrc_clean <- bvbrc_clean |> + mutate( + across(all_of(list_cols), clean_list_col) + ) + +# Save the table to parquet + +flatten_list <- function(x) { + purrr::map_chr( + x, + \(z) { + if (is.null(z) || length(z) == 0) + return(NA_character_) + + paste(as.character(unlist(z)), collapse = ";") + } + ) +} + +bvbrc_clean <- bvbrc_clean |> + dplyr::mutate( + dplyr::across( + dplyr::all_of(list_cols), + flatten_list + ) + ) + +arrow::write_parquet( + bvbrc_clean, + "data/bvbrc_clean.parquet", + compression = "zstd" +) + # --- Per-column stats -------------------------------------------------------- col_stats <- tibble::tibble(column = colnames(bvbrc_clean)) |> dplyr::mutate(