realign_with_best_hits: Realign top blast hit of multi-fasta file with that fasta...

In joelnitta/baitfindR: Find Baits for Sequence Capture

Description

Top hits and original fasta files are matched based on the first part of the filename separated by periods (i.e., the filename without any extension).

Usage

realign_with_best_hits(
  best_hits_dir,
  best_hits_pattern = "bestmatch",
  fasta_dir,
  fasta_pattern = "\\.fa$",
  ...
)

Arguments

best_hits_dir	Path to directory containing top blast hits.
best_hits_pattern	Pattern used for matching with grep. Only files with names matching the pattern will be included as the top blast hit.
fasta_dir	Path to directory containing fasta files for realignment.
fasta_pattern	Pattern used for matching with grep. Only files with names matching the pattern will be included for realignment.
...	Additional other arguments. Not used by this function, but meant to be used by drake_plan for tracking during workflows.

Value

List of lists, each of which is of class 'DNAbin'.

Examples

library(ape)

# Make temp dir for storing files
temp_dir <- fs::dir_create(fs::path(tempdir(), "baitfindR_example"))

# Write out ape::woodmouse dataset as DNA
data(woodmouse)
ape::write.FASTA(woodmouse, fs::path(temp_dir, "woodmouse.fasta"))
ape::write.FASTA(woodmouse, fs::path(temp_dir, "woodmouse2.fasta"))

# Make blast database
build_blast_db(
  fs::path(temp_dir, "woodmouse.fasta"),
  db_type = "nucl",
  out_name = "wood",
  parse_seqids = TRUE,
  wd = temp_dir)

# Blast the original sequences against the database
blast_n_list(
  fasta_folder = temp_dir,
  fasta_pattern = "fasta",
  database_path = fs::path(temp_dir, "wood")
)

# Extract the top BLAST hit for each fasta file.
extract_blast_hits(
  blast_results_dir = temp_dir,
  blast_results_pattern = "\\.tsv$",
  database_path = fs::path(temp_dir, "wood"),
  out_dir = temp_dir,
  out_ext = "bestmatch"
)

realign_with_best_hits(
  best_hits_dir = temp_dir,
  best_hits_pattern = "bestmatch",
  fasta_dir = temp_dir,
  fasta_pattern = "fasta"
)

# Cleanup.
fs::file_delete(temp_dir)

joelnitta/baitfindR documentation built on May 7, 2020, 6:21 p.m.