Skip to content

Export metadata and feature data tables - #34

Open
epbrenner wants to merge 18 commits into
mainfrom
export_tables
Open

Export metadata and feature data tables#34
epbrenner wants to merge 18 commits into
mainfrom
export_tables

Conversation

@epbrenner

Copy link
Copy Markdown
Contributor

Testing:

For this PR, I have included bundled test data in /inst/extdata/Campy_testdata.zip. Extract this into a single directory of your choosing. This bundled data is not intended to the ultimate example data for final package release.

For data_curation.R:
Run
prepareGenomes("Staphylococcus argenteus")
then
exportTables("data/Staphylococcus_argenteus/Sar.duckdb")

For data_processing.R:
Unzip the Campy test data first, then run
exportProcessedData("inst/extdata/Campy_testdata/Cco_parquet.duckdb")

There are several thrilling parameters to fiddle around with.

What this does:

Adds the ability to export basic stats into human readable CSVs from data_curation.R workflow.
Adds the ability to export all the processed feature data into a variety of formats from data_processing.R workflow.

What's the point?

Makes amRdata more independently useful to people even if they don't run the full amR suite. For example, you could export feature matrices for custom modeling workflows outside of ours, or just get an idea of what distributions of different QC and AST labels are across different taxa.

epbrenner and others added 16 commits July 22, 2026 11:17
Added ability to filter BV-BRC genomes based on genome stat metadata, including CheckM completeness and CheckM contamination as defaults. Also enabled optional checking of genome length, CDS count, and GC content based on deviations from median distribution per bug. Updated genome downloader functions to use ftpes_download_one instead of old ftp_download_one, and moved a dot function helper to inside the ftpes_download_one loop to avoid issues with BiocParallel workers not understanding what that helper is.
Uncommenting a helper that was still needed for several calls. Fixing stuff I broke myself!
Updating to remove BiocParallel implementation that caused bugs. Patching back to prior Future implementation instead. CLI version of genome downloading appears to have issues that may be separate and requires further troubleshooting.
data_processing.R runDataProcessing() missing(ref_file_path) check ignores if the default path is used. Updated to accept the default path if a user doesn't specify rather than erroring at the end of the processing steps.
Updated default workers to 8 for CLI and FTP downloads, made two-try FTP download option the default, added a helper to purge broken file sets (i.e., if .faa is missing for an AccID but .fna and .gff download successfully).
ftpes --> ftps
is_complete --> .is_complete
Cleaning up parallel implementation in data_curation.R, standardizing with furrr and dropping old duplicated within-loop helper function that was needed for earlier BiocParallel implementation. Standardized values for contamination and completeness across functions. Added a parameter in prepareGenomes to allow users to set number of workers. Reduced console clutter for final "out" print.
In data_curation.R, adding ability to export genome stats and QC metadata into human readable CSVs. `load_tables` allows you to run "example <- exportTables("data/Staphylococcus_aureus/Sau.duckdb", load_tables=TRUE)" and import the tables into R so you can play around with them directly. By default, it just saves the tables as exported CSVs. In data_processing.R, added a similar exporter that allows you to save out all the feature data and other DuckDB-connected tables, with the ability to output sequences too, and append AST phenotype data to the feature tables. Tables can be saved out as csv (default), tsv, Parquet, and xlsx. Added `writexl` to Suggests in DESCRIPTION. Fixed some reversions where I broke things in the PR this branch is based on. Will add those same fixes to that PR momentarily to keep things aligned.
Ran data processing workflow for Campylobacter coli. This can be unzipped and used as the test data for exportProcessedData() rather than having to run it all beforehand. Data curation step still requires that prepareGenomes() be run for a test bug, and Staphylococcus argenteus is a good use case for this.
@epbrenner epbrenner self-assigned this Jul 28, 2026
@epbrenner epbrenner added the enhancement New feature or request label Jul 28, 2026

@jananiravi jananiravi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Few quick comments. Will test it out tomorrow.

Comment thread R/data_curation.R
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R
)
return(all_complete)
if (isTRUE(verbose))
message("All genomes already complete.")

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not complete status but completed downloading? Is something missing here?

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We check the length of the diff of already (already downloaded) and ids (the full set needed), and if it's 0, then none require download. All genomes complete message is triggered.

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Makes sense. Download of all genomes complete.
Or something like that indicates what's complete?

Clarified a few names based on code review feedback.
jananiravi
jananiravi previously approved these changes Jul 29, 2026
@jananiravi

Copy link
Copy Markdown
Member

Few quick comments. Will test it out tomorrow.

@epbrenner prepareGenomes + exportTables --> they run locally with Sar.

@jananiravi

jananiravi commented Jul 29, 2026

Copy link
Copy Markdown
Member

After running the campy --> metadata.csv

image

what's TET doing in resistant_classes? genome_drug.resistance_phenotype marks them as susceptible. @epbrenner @AbhirupaGhosh

@AbhirupaGhosh

Copy link
Copy Markdown
Contributor

After running the campy --> metadata.csv

image what's TET doing in `resistant_classes`? genome_drug.resistance_phenotype marks them as susceptible. @epbrenner @AbhirupaGhosh

resistant_classes is a column that adds all the drug classes the genome is resistant to. If there are more than 1, say TET and FLQ, it will show TET_FLQ and it shows susceptible to the genomes that are resistant to none.

@jananiravi jananiravi left a comment

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

is this going to be merged to qc_filtering and then qc_filtering --> main? that might be helpful to avoid duplicate reviews.

Comment thread R/data_processing.R
gene_count = list(source = "gene_count", stem = "gene_count", appendable = TRUE),
protein_count = list(source = "protein_count", stem = "protein_count", appendable = TRUE),
domain_count = list(source = "domain_count", stem = "domain_count", appendable = TRUE),
struct = list(source = "struct", stem = "struct", appendable = TRUE),

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

struct or gene_struct?

Comment thread R/data_curation.R
Comment on lines +2486 to +2493
count_if_present("metadata"),
count_if_present("genome_data"),
count_if_present("amr_phenotype"),
count_if_present("metadata_qc"),
count_if_present("metadata_qc_rejections"),
count_if_present("filtered"),
count_if_present("files"),
count_if_present("bac_data")

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

how do these checks tie back to basic_tables? is the latter not exported? so should "files" be added to basic_tables?

Comment thread R/data_curation.R
"genome_data",
"amr_phenotype",
"metadata_qc",
"metadata_qc_rejections"

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

see below. add files?

@eboyer221 eboyer221 left a comment

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

When I tested this, it looks like genome IDs such as 263.1, 263.10, and 263.100 are being written to the CSV as plain numbers rather than text, since they happen to look numeric. When I read that export back (with readr::read_csv(), pandas, or Excel's defaults) all three collapsed to the same value, since trailing zeros after a decimal don't count (263.1 and 263.10 parse as numerically identical). So the three distinct genomes merged under one ID, and it'll happen for any taxon with more than ~9 genomes. I also tried forcing quotes around every value (quote = "all") to see if that would prevent it, but the collision still happened on read. Quoting only protects against commas breaking the file format, it doesn't tell the reader "treat this as text."

This affects every genome_id-keyed table both export functions produce (genome_data, metadata, amr_phenotype_wide, gene_count, etc.).
Two options I see: (1) add a clear warning telling users to force col_types/dtype=str on ID columns when reading CSV/TSV back, or (2) recommend Parquet as the default/safe export format for programmatic reuse, since arrow preserves the real column type - I verified this did not occur with Parquet.

Comment thread R/data_curation.R
#'
#' @return Invisibly returns a list with exported file paths and, if requested,
#' in-memory tables.
#' @keywords internal

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
#' @keywords internal
#' @export

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Without @export this stays internal-only after a normal install. library(amRdata); exportTables(...) fails even though that's exactly what the PR description tells testers to run.

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

After applying this change, run devtools::document() and commit the regenerated NAMESPACE/man/*.Rd. exportProcessedData() has the same problem (has @export in source but isn't in NAMESPACE since docs were never regenerated on this branch).

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

And that's why it didn't run the first time for me. Good catch, Emily!!

@epbrenner epbrenner left a comment

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh I have to submit my own review to my own PR to actually post these comments, huh? Oops.

Comment thread R/data_curation.R
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R
Comment thread R/data_curation.R Outdated
Comment thread R/data_curation.R
)
return(all_complete)
if (isTRUE(verbose))
message("All genomes already complete.")

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We check the length of the diff of already (already downloaded) and ids (the full set needed), and if it's 0, then none require download. All genomes complete message is triggered.

@jananiravi

Copy link
Copy Markdown
Member

When I tested this, it looks like genome IDs such as 263.1, 263.10, and 263.100 are being written to the CSV as plain numbers rather than text, since they happen to look numeric. When I read that export back (with readr::read_csv(), pandas, or Excel's defaults) all three collapsed to the same value, since trailing zeros after a decimal don't count (263.1 and 263.10 parse as numerically identical). So the three distinct genomes merged under one ID, and it'll happen for any taxon with more than ~9 genomes. I also tried forcing quotes around every value (quote = "all") to see if that would prevent it, but the collision still happened on read. Quoting only protects against commas breaking the file format, it doesn't tell the reader "treat this as text."

This affects every genome_id-keyed table both export functions produce (genome_data, metadata, amr_phenotype_wide, gene_count, etc.).
Two options I see: (1) add a clear warning telling users to force col_types/dtype=str on ID columns when reading CSV/TSV back, or (2) recommend Parquet as the default/safe export format for programmatic reuse, since arrow preserves the real column type - I verified this did not occur with Parquet.

@epbrenner @AbhirupaGhosh didnt we have this fix in already?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

enhancement New feature or request

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants