View source: R/convert_animalcules.R
convert_animalcules | R Documentation |
Upon completion of the MetaScope pipeline, users can analyze and visualize
abundances in their samples using the animalcules package. This function
allows interoperability of metascope_id
output with both animalcules
and QIIME. After running this function, the user should save the returned MAE
to an RDS file using a function like saveRDS
to upload the output into
the animalcules
package.
convert_animalcules(
meta_counts,
annot_path,
which_annot_col,
end_string = ".metascope_id.csv",
qiime_biom_out = FALSE,
path_to_write = ".",
NCBI_key = NULL,
num_tries = 3
)
meta_counts |
A vector of filepaths to the counts ID CSVs output by
|
annot_path |
The filepath to the CSV annotation file for the samples.
This CSV metadata/annotation file should contain at least two columns,
one with names of all samples WITHOUT the extension listed in |
which_annot_col |
The name of the column of the annotation file
containing the sample IDs. These should be the same as the
|
end_string |
The end string used at the end of the metascope_id files. Default is ".metascope_id.csv". |
qiime_biom_out |
Would you also like a qiime-compatible biom file
output? If yes, two files will be saved: one is a biom file of the counts
table, and the other is a specifically formatted mapping file of metadata
information. Default is |
path_to_write |
If |
NCBI_key |
(character) NCBI Entrez API key. optional. See taxize::use_entrez(). Due to the high number of requests made to NCBI, this function will be less prone to errors if you obtain an NCBI key. |
num_tries |
(numeric) Number of attempts to get UID. |
Returns a MultiAssay Experiment file of combined sample counts data and/or biom file and mapping file for analysis with QIIME. The MultiAssay Experiment will have a counts assay ("MGX").
tempfolder <- tempfile()
dir.create(tempfolder)
# Create three different samples
samp_names <- c("X123", "X456", "X789")
all_files <- file.path(tempfolder,
paste0(samp_names, ".csv"))
create_IDcsv <- function (out_file) {
final_taxids <- c("273036", "418127", "11234")
final_genomes <- c(
"Staphylococcus aureus RF122, complete sequence",
"Staphylococcus aureus subsp. aureus Mu3, complete sequence",
"Measles virus, complete genome")
best_hit <- sample(seq(100, 1050), 3)
proportion <- best_hit/sum(best_hit) |> round(2)
EMreads <- best_hit + round(runif(3), 1)
EMprop <- proportion + 0.003
dplyr::tibble(TaxonomyID = final_taxids,
Genome = final_genomes,
read_count = best_hit, Proportion = proportion,
EMreads = EMreads, EMProportion = EMprop) |>
dplyr::arrange(dplyr::desc(.data$read_count)) |>
utils::write.csv(file = out_file, row.names = FALSE)
message("Done!")
return(out_file)
}
out_files <- vapply(all_files, create_IDcsv, FUN.VALUE = character(1))
# Create annotation data for samples
annot_dat <- file.path(tempfolder, "annot.csv")
dplyr::tibble(Sample = samp_names, RSV = c("pos", "neg", "pos"),
month = c("March", "July", "Aug"),
yrsold = c(0.5, 0.6, 0.2)) |>
utils::write.csv(file = annot_dat,
row.names = FALSE)
# Convert samples to MAE
outMAE <- convert_animalcules(meta_counts = out_files,
annot_path = annot_dat,
which_annot_col = "Sample",
end_string = ".metascope_id.csv",
qiime_biom_out = FALSE,
NCBI_key = NULL)
unlink(tempfolder, recursive = TRUE)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.