Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

add error handling in tree.R #70

Merged
merged 5 commits into from
Oct 30, 2024
Merged
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
61 changes: 57 additions & 4 deletions R/tree.R
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,8 @@
#' be saved. Default is the path to "data/alns/pspa_snf7.tre".
#' @param fasttree_path Path to the FastTree executable, which is used to
#' generate the phylogenetic tree. Default is "src/FastTree".
#'
#' @importFrom rlang abort
#'
#' @return No return value. The function generates a tree file (.tre) from the
#' input FASTA file.
Expand All @@ -60,7 +62,31 @@ convertFA2Tree <- function(fa_path = here("data/alns/pspa_snf7.fa"),
# fa_path=here("data/alns/pspa_snf7.fa")
# tre_path=here("data/alns/pspa_snf7.tre")
# fasttree_path=here("src/FastTree")
print(fa_path)

# Check if the FASTA file exists
if (!file.exists(fa_path)) {
abort(paste("Error: The FASTA file does not exist at:", fa_path))
}

# Check if the FastTree executable exists
if (!file.exists(fasttree_path)) {
abort(paste("Error: The FastTree executable does not exist at:",
fasttree_path))
}

# Check if the output directory exists
tre_dir <- dirname(tre_path)
if (!dir.exists(tre_dir)) {
abort(paste("Error: The output directory does not exist:", tre_dir))
}

# Check if the output file already exists
if (file.exists(tre_path)) {
cat("Warning: The output file already exists and will be overwritten:",
tre_path, "\n")
}

message(fa_path)
system2(
command = fasttree_path,
args = paste(c(fa_path, ">", tre_path),
Expand Down Expand Up @@ -98,8 +124,18 @@ convertFA2Tree <- function(fa_path = here("data/alns/pspa_snf7.fa"),
#' generate_trees(here("data/alns/"))
#' }
convertAlignment2Trees <- function(aln_path = here("data/alns/")) {

# Check if the alignment directory exists
if (!dir.exists(aln_path)) {
abort(paste("Error: The alignment directory does not exist:", aln_path))
}
# finding all fasta alignment files
fa_filenames <- list.files(path = aln_path, pattern = "*.fa")
# Check if any FASTA files were found
if (length(fa_filenames) == 0) {
abort("Error: No FASTA files found in the specified directory.")
}

fa_paths <- paste0(aln_path, fa_filenames)
variable <- str_replace_all(basename(fa_filenames),
pattern = ".fa", replacement = ""
Expand Down Expand Up @@ -157,6 +193,23 @@ createFA2Tree <- function(fa_file = "data/alns/pspa_snf7.fa",
## SAMPLE ARGS
# fa_file="data/alns/pspa_snf7.fa"
# out_file="data/alns/pspa_snf7.tre"

# Check if the FASTA file exists
if (!file.exists(fa_file)) {
abort(paste("Error: The FASTA file does not exist at:", fa_file))
}

# Check if the output directory exists
out_dir <- dirname(out_file)
if (!dir.exists(out_dir)) {
abort(paste("Error: The output directory does not exist:", out_dir))
}

# Check if the output file already exists
if (file.exists(out_file)) {
cat("Warning: The output file already exists and will be overwritten:",
out_file, "\n")
}

###########################
## Approach 1
Expand All @@ -182,7 +235,7 @@ createFA2Tree <- function(fa_file = "data/alns/pspa_snf7.fa",
## Model Testing & Distance Matrices
## Comparison of different nucleotide or amino acid substitution models
mt <- modelTest(prot10, model = "all")
print(mt)
message(mt)

# estimate a distance matrix using a Jules-Cantor Model
dna_dist <- dist.ml(prot10, model = "JC69")
Expand All @@ -203,7 +256,7 @@ createFA2Tree <- function(fa_file = "data/alns/pspa_snf7.fa",
## Maximum likelihood and Bootstrapping
# ml estimation w/ distance matrix
fit <- pml(prot_NJ, prot10)
print(fit)
message(fit)
fitJC <- optim.pml(fit, model = "JC", rearrangement = "stochastic")
logLik(fitJC)
bs <- bootstrap.pml(fitJC,
Expand All @@ -216,7 +269,7 @@ createFA2Tree <- function(fa_file = "data/alns/pspa_snf7.fa",
prot10_dm <- dist.ml(prot10)
prot10_NJ <- NJ(prot10_dm)
fit2 <- pml(prot10_NJ, data = prot10)
print(fit2)
message(fit2)
fitJC2 <- optim.pml(fit2, model = "JC", rearrangement = "stochastic")
logLik(fitJC2)
bs_subset <- bootstrap.pml(fitJC2,
Expand Down
Loading