Mycelia Documentation

AssemblyAction

Action representation for reinforcement learning decisions during assembly.

Fields

• decision::Symbol: Primary decision (:continuek, :nextk, :terminate)
• viterbi_params::Dict{Symbol, Float64}: Viterbi algorithm parameters
• correction_threshold::Float64: Quality threshold for error correction
• batch_size::Int: Batch size for processing (memory management)
• max_iterations::Int: Maximum iterations at current k before forced progression

Assembly configuration structure.

AssemblyEnvironment

Reinforcement learning environment for training assembly decision policies.

Fields

• current_state::AssemblyState: Current environment state
• training_datasets::Vector{String}: Paths to training FASTQ files
• validation_datasets::Vector{String}: Paths to validation FASTQ files
• episode_length::Int: Maximum steps per training episode
• step_count::Int: Current step in episode
• reward_history::Vector{Float64}: Reward history for current episode
• action_history::Vector{AssemblyAction}: Action history for experience replay
• assembly_cache::Dict{String, Any}: Cache for intermediate assembly results

Assembly method enumeration for unified interface.

Assembly result structure containing contigs and metadata.

AssemblyState

State representation for reinforcement learning environment containing all information needed to make assembly decisions.

Fields

• current_k::Int: Current k-mer size being processed
• assembly_quality::Float64: Current assembly quality score (QV-based)
• correction_rate::Float64: Rate of successful error corrections in recent iterations
• memory_usage::Float64: Current memory utilization (fraction of limit)
• graph_connectivity::Float64: Graph connectivity metric (proportion of strongly connected components)
• coverage_uniformity::Float64: Uniformity of k-mer coverage distribution
• error_signal_clarity::Float64: Clarity of error signal detection (sparsity-based)
• iteration_history::Vector{Float64}: Recent reward history for trend analysis
• k_progression::Vector{Int}: Sequence of k-mer sizes processed so far
• corrections_made::Int: Total corrections made at current k
• time_elapsed::Float64: Time spent on current k-mer size (seconds)

Benchmark configuration for consistent testing.

Edge data for variable-length BioSequence graphs.

Vertex data for variable-length BioSequence graphs.

BubbleStructure

Represents a bubble (alternative path) in the assembly graph.

ContigPath

Represents a linear path through the graph forming a contig.

DQNPolicy

Deep Q-Network policy for high-level assembly decisions.

This is a placeholder structure for the neural network architecture that will be implemented with a machine learning framework like Flux.jl or MLJ.jl.

Fields

• state_dim::Int: Dimension of state representation
• action_dim::Int: Number of possible actions
• hidden_dims::Vector{Int}: Hidden layer dimensions
• learning_rate::Float64: Learning rate for training
• epsilon::Float64: Exploration rate for epsilon-greedy policy
• experience_buffer::Vector{Any}: Experience replay buffer

Stores comprehensive FASTQ quality analysis results.

Graph mode for handling strand information.

• SingleStrand: Sequences are single-stranded (RNA, amino acids, or directional DNA)
• DoubleStrand: Sequences are double-stranded DNA/RNA with canonical representation

Represents a complete path through the k-mer graph.

Type-stable metadata for k-mer graph edges.

Edges represent valid strand-aware transitions between canonical k-mers. The transition is valid only if the strand orientations allow for proper overlap.

Fields:

• coverage: Vector of edge traversal observations with strand information
• weight: Edge weight/confidence score based on coverage
• src_strand: Required strand orientation of source k-mer for this transition
• dst_strand: Required strand orientation of destination k-mer for this transition

Type-stable metadata for k-mer graph vertices.

Vertices always represent canonical k-mers for memory efficiency and cleaner graphs. Strand information is tracked in the coverage data and edge transitions.

Fields:

• coverage: Vector of observation coverage data as (observationid, position, strandorientation) tuples
• canonical_kmer: The canonical k-mer (BioSequence type - NO string conversion)

PangenomeAnalysisResult

Results of k-mer based pangenome analysis.

Quality-aware BioSequence edge data.

Quality-aware BioSequence vertex data.

Stores quality distribution statistics for FASTQ analysis.

Edge data for quality-aware k-mer graphs.

Qualmer observation for tracking k-mer occurrences in sequences.

Vertex data for quality-aware k-mer graphs.

RepeatRegion

Represents a repetitive region in the assembly graph.

RewardComponents

Structured representation of reward signal components for training the RL agent.

Fields

• accuracy_reward::Float64: Primary reward based on assembly accuracy (weighted 1000x)
• efficiency_reward::Float64: Secondary reward for computational efficiency (weighted 10x)
• error_penalty::Float64: Penalty for false positives/negatives (weighted -500x)
• progress_bonus::Float64: Bonus for making meaningful progress
• termination_reward::Float64: Reward for appropriate termination timing
• total_reward::Float64: Weighted sum of all components

ScaffoldResult

Results from scaffolding analysis.

Strand orientation for k-mer observations and transitions.

• Forward: k-mer as observed (5' to 3')
• Reverse: reverse complement of k-mer (3' to 5')

ViterbiConfig

Configuration for Viterbi algorithm parameters.

ViterbiPath

Complete Viterbi path through the k-mer graph.

ViterbiState

State information for Viterbi algorithm on k-mer graphs.

Represents a step in a probabilistic walk through the graph.

JLD2_read_table(filename::String) -> Any

Read a DataFrame from a JLD2 file without needing to know the internal name. If the file contains multiple DataFrames, returns the first one found.

JLD2_write_table(; df, filename)

Write a DataFrame to a JLD2 file using a standardized internal name.

Add edges to BioSequence graph based on sequence overlaps.

Add edges between BioSequences based on k-mer path relationships.

_add_observation_to_graph!(
    graph,
    observation,
    obs_idx,
    canonical_kmers,
    graph_mode
)

Add a sequence observation to an existing k-mer graph with strand-aware edge creation.

Arguments

• graph: MetaGraphsNext k-mer graph with canonical vertices
• observation: FASTA/FASTQ record
• obs_idx: Observation index
• canonical_kmers: Vector of canonical k-mers in the graph
• graph_mode: SingleStrand or DoubleStrand mode

Add quality-weighted edges to BioSequence graph.

Add edges between quality-aware BioSequences based on qualmer path relationships.

Add edges to the qualmer graph based on k-mer adjacency in sequences.

Add simulated sequencing errors to a sequence string.

Helper function to add strand-aware coverage data to an edge.

This function creates edges that respect strand orientation constraints. Each edge represents a biologically valid transition between k-mers.

Helper function to add coverage data to a vertex.

BioSequence graph assembly implementation (variable-length simplified from k-mer graphs).

Hybrid OLC assembly (placeholder for future implementation).

K-mer graph assembly implementation (fixed-length k-mer foundation).

Multi-k assembly (placeholder for future implementation).

N-gram graph assembly implementation (fixed-length unicode character analysis).

Quality-aware BioSequence graph assembly implementation (variable-length simplified from qualmer graphs).

Quality-aware k-mer graph assembly implementation (fixed-length qualmer foundation).

String graph assembly implementation (variable-length simplified from N-gram graphs).

Calculate consensus quality from multiple observations at a specific position. Uses the joint probability and multiple observations to compute a consensus quality score.

_calculate_l_statistic(sorted_lengths, threshold)

Calculate L-statistic (number of contigs needed to reach a given percentage of total assembly length). For example, L50 is the number of contigs needed to reach 50% of the total assembly length.

Arguments

• sorted_lengths: Vector of contig lengths sorted in descending order
• threshold: Fraction of total length (e.g., 0.5 for L50, 0.9 for L90)

Returns

• Int: Number of contigs needed to reach the threshold

Calculate N-statistic (N50, N90, etc.) for contig lengths.

Calculate coverage for a path from constituent k-mers.

Calculate coverage for qualmer path.

Helper function to calculate normalized transition probabilities.

_check_conda_env_exists(env_name::AbstractString) -> Bool

Check if a conda environment exists.

Convert contigs to FASTA records for graph construction.

Helper function to create an empty k-mer graph.

_detect_sequence_extension(sequence_type::Symbol) -> String

Internal helper function to convert sequence type to file extension.

Arguments

• sequence_type: Symbol representing sequence type (:DNA, :RNA, or :AA)

Returns

• String: Appropriate file extensions

Determine appropriate k-mer type from observations.

Determine strand orientation by comparing original and canonical qualmers.

Find the heaviest (highest confidence) Eulerian path starting from a given vertex.

Find linear paths in a k-mer graph.

Find linear paths in a qualmer graph.

Find quality-weighted overlap between two sequences.

Find quality-aware paths in a qualmer graph.

Find overlap between two sequences.

Generate contigs from qualmer graph using probabilistic walks.

Generate contigs using probabilistic walks when Eulerian paths are not available.

Generate FASTQ contigs from qualmer graph using probabilistic walks when no paths are found.

_get_output_files(output_dir::AbstractString) -> Vector{String}

Get a list of all files in the output directory.

Helper function to get valid transitions from a vertex with given strand orientation.

Validate that the transition between two k-mers is biologically valid.

_is_valid_transition(
    src_kmer,
    dst_kmer,
    src_strand,
    dst_strand,
    k
) -> Any

Validate that a transition between two k-mers with given strand orientations is biologically valid.

For a transition to be valid, the suffix of the source k-mer must match the prefix of the destination k-mer when accounting for strand orientations.

Arguments

• src_kmer: Source canonical k-mer
• dst_kmer: Destination canonical k-mer
• src_strand: Strand orientation of source k-mer
• dst_strand: Strand orientation of destination k-mer
• k: K-mer size

Returns

• Bool: true if transition is biologically valid

Iterative Viterbi algorithm for finding optimal paths through qualmer graph. Uses dynamic programming with quality scores as emission/transition probabilities.

Convert a path of vertices to a DNA sequence.

Polish a single contig using Viterbi error correction.

Convert observations to FASTQ records for quality processing.

Prepare observations from various input formats (FASTA/FASTQ records or file paths).

Perform quality-weighted walk through qualmer graph.

Enhanced qualmer path to sequence with consensus quality calculation. Uses joint probability from multiple observations to compute consensus quality.

Enhanced qualmer path to FASTQ record conversion with quality propagation. This is the core function that enables quality-aware assembly output.

Convert qualmer path to DNA sequence.

Reconstruct sequence and quality from qualmer path.

Reconstruct a BioSequence from a path of k-mers.

Helper function to reconstruct sequence from a graph path.

Helper function to reconstruct shortest path from Dijkstra's algorithm.

Helper function to sample a transition based on probabilities.

_sequence_to_canonical_path(
    canonical_kmers,
    sequence,
    graph_mode
) -> Vector{<:Tuple{Any, Mycelia.StrandOrientation}}

Convert a sequence to a path through canonical k-mer space with strand awareness.

This is the key function that handles the distinction between single-strand and double-strand modes while maintaining canonical k-mer vertices.

Arguments

• canonical_kmers: Vector of canonical k-mers available in the graph
• sequence: DNA/RNA sequence to convert
• graph_mode: SingleStrand or DoubleStrand mode

Returns

• Vector of (canonicalkmer, strandorientation) pairs representing the path

Details

In DoubleStrand mode:

• Each observed k-mer is converted to its canonical form
• Strand orientation tracks whether the canonical form matches the observed k-mer
• Edges respect the biological constraint that transitions must maintain proper overlap

In SingleStrand mode:

• K-mers are used as-is (no reverse complement consideration)
• All strand orientations are Forward
• Suitable for RNA, amino acids, or directional DNA analysis

_setup_phageboost_environment(force_reinstall::Bool=false)

Set up the PhageBoost conda environment if it doesn't exist or if force_reinstall is true.

Simplify N-gram graph by removing unnecessary complexity.

Convert N-gram graph to variable-length string graph.

Simplify string graph by removing unnecessary complexity.

_simulate_fastq_reads_from_sequence(
    sequence,
    identifier::String;
    coverage,
    read_length,
    error_rate,
    sequence_type
) -> Vector{FASTX.FASTQ.Record}

Generate simulated FASTQ reads from a reference sequence.

Arguments

• sequence: Reference sequence (String or BioSequences type)
• identifier: Sequence identifier
• coverage::Int: Desired coverage depth
• read_length::Int: Length of individual reads
• error_rate::Float64: Simulated error rate
• sequence_type::String: Type of sequence ("DNA", "RNA", "AA")

Returns

• Vector{FASTX.FASTQ.Record}: Vector of simulated FASTQ reads with quality scores

Validate assembly against reference sequence (placeholder).

Find optimal path using Viterbi-like dynamic programming on quality scores.

accuracy(true_labels, pred_labels)

Returns the overall accuracy.

add_bioconda_env(pkg; force) -> Union{Nothing, Base.Process}

Create a new Conda environment with a specified Bioconda package.

Arguments

• pkg::String: Package name to install. Can include channel specification using

the format "channel::package"

Keywords

• force::Bool=false: If true, recreates the environment even if it already exists

Details

The function creates a new Conda environment named after the package and installs the package into it. It uses channel priority: conda-forge > bioconda > defaults. If CONDA_RUNNER is set to 'mamba', it will ensure mamba is installed first.

Examples

# Install basic package
add_bioconda_env("blast")

# Install from specific channel
add_bioconda_env("bioconda::blast")

# Force reinstallation
add_bioconda_env("blast", force=true)

Notes

• Requires Conda.jl to be installed and configured
• Uses CONDA_RUNNER global variable to determine whether to use conda or mamba
• Cleans conda cache after installation

add_edgemer_to_graph!(
    graph,
    record_identifier,
    index,
    observed_edgemer
) -> Any

Add an observed edgemer to a graph with its associated metadata.

Arguments

• graph::MetaGraph: The graph to modify
• record_identifier: Identifier for the source record
• index: Position where edgemer was observed
• observed_edgemer: The biological sequence representing the edgemer

Details

Processes the edgemer by:

1. Splitting it into source and destination kmers
2. Converting kmers to their canonical forms
3. Creating or updating an edge with orientation metadata
4. Storing observation details (record, position, orientation)

Returns

Modified graph with the new edge and metadata

Note

If the edge already exists, the observation is added to the existing metadata.

add_fastx_records_to_graph!(graph, fastxs) -> Any

Add FASTX records from multiple files as a graph property.

Arguments

• graph: A MetaGraph that will store the FASTX records
• fastxs: Collection of FASTA/FASTQ file paths to process

Details

Creates a dictionary mapping sequence descriptions to their corresponding FASTX records, then stores this dictionary as a graph property under the key :records. Multiple input files are merged, with later files overwriting records with duplicate descriptions.

Returns

The modified graph with added records property.

add_record_edgemers_to_graph!(graph) -> Any

Processes DNA sequence records stored in the graph and adds their edgemers (k+1 length subsequences) to build the graph structure.

Arguments

• graph: A Mycelia graph object containing DNA sequence records and graph properties

Details

• Uses the k-mer size specified in graph.gprops[:k] to generate k+1 length edgemers
• Iterates through each record in graph.gprops[:records]
• For each record, generates all possible overlapping edgemers
• Adds each edgemer to the graph with its position and record information

Returns

• The modified graph object with added edgemer information

Adjust quality scores based on likelihood improvement.

alphabet_to_biosequence_type(
    alphabet::Symbol
) -> Union{Type{BioSequences.LongAA}, Type{BioSequences.LongSequence{BioSequences.DNAAlphabet{4}}}, Type{BioSequences.LongSequence{BioSequences.RNAAlphabet{4}}}}

Determine the BioSequence type from an alphabet symbol.

Maps alphabet symbols to the corresponding BioSequences.jl type for type-safe sequence operations throughout the codebase.

Arguments

• alphabet::Symbol: The alphabet symbol (:DNA, :RNA, or :AA)

Returns

• Type{<:BioSequences.BioSequence}: The corresponding BioSequence type

Examples

alphabet_to_biosequence_type(:DNA)  # Returns BioSequences.LongDNA{4}
alphabet_to_biosequence_type(:RNA)  # Returns BioSequences.LongRNA{4}
alphabet_to_biosequence_type(:AA)   # Returns BioSequences.LongAA

amino_acids_to_codons(

) -> Dict{BioSymbols.AminoAcid, DataType}

Creates a mapping from amino acids to representative DNA codons using the standard genetic code.

Returns

• Dictionary mapping each amino acid (including stop codon AA_Term) to a valid DNA codon that encodes it

analyze_fastq_quality(fastq_file::String)

Analyzes quality metrics for a FASTQ file.

Calculates comprehensive quality statistics including read count, quality scores, length distribution, GC content, and quality threshold percentages.

Arguments

• fastq_file: Path to FASTQ file (can be gzipped)

Returns

FastqQualityResults with the following fields:

• n_reads: Total number of reads
• mean_quality: Average Phred quality score across all reads
• mean_length: Average read length
• gc_content: GC content percentage
• quality_distribution: QualityDistribution with Q20+, Q30+, Q40+ percentages

Example

quality_stats = Mycelia.analyze_fastq_quality("reads.fastq")
println("Total reads: $(quality_stats.n_reads)")
println("Mean quality: $(quality_stats.mean_quality)")
println("Q30+ reads: $(quality_stats.quality_distribution.q30_percent)%")

analyze_pangenome_kmers(genome_files::Vector{String}; kmer_type=Kmers.DNAKmer{21}, distance_metric=:jaccard)

Perform comprehensive k-mer based pangenome analysis using existing Mycelia infrastructure.

Leverages existing count_canonical_kmers and distance metric functions to analyze genomic content across multiple genomes, identifying core, accessory, and unique regions.

Arguments

• genome_files: Vector of FASTA file paths containing genome sequences
• kmer_type: K-mer type from Kmers.jl (default: Kmers.DNAKmer{21})
• distance_metric: Distance metric (:jaccard, :bray_curtis, :cosine, :js_divergence)

Returns

• PangenomeAnalysisResult with comprehensive pangenome statistics

Example

genome_files = ["genome1.fasta", "genome2.fasta", "genome3.fasta"]
result = Mycelia.analyze_pangenome_kmers(genome_files, kmer_type=Kmers.DNAKmer{31})
println("Core k-mers: $(length(result.core_kmers))")
println("Total pangenome size: $(size(result.presence_absence_matrix, 1)) k-mers")

Analyze a potential repeat region starting from a vertex.

annotate_aa_fasta(
;
    fasta,
    identifier,
    basedir,
    mmseqsdb,
    threads
)

Annotate amino acid sequences in a FASTA file using MMseqs2 search against UniRef50 database.

Arguments

• fasta: Path to input FASTA file containing amino acid sequences
• identifier: Name for the output directory (defaults to FASTA filename without extension)
• basedir: Base directory for output (defaults to current directory)
• mmseqsdb: Path to MMseqs2 formatted UniRef50 database (defaults to ~/workspace/mmseqs/UniRef50)
• threads: Number of CPU threads to use (defaults to system thread count)

Returns

• Path to the output directory containing MMseqs2 search results

The function creates a new directory named by identifier under basedir, copies the input FASTA file, and runs MMseqs2 easy-search against the specified database. If the output directory already exists, the function skips processing and returns the directory path.

annotate_fasta(
;
    fasta,
    identifier,
    basedir,
    mmseqsdb,
    threads
)

Perform comprehensive annotation of a FASTA file including gene prediction, protein homology search, and terminator prediction.

Arguments

• fasta::String: Path to input FASTA file
• identifier::String: Unique identifier for output directory (default: FASTA filename without extension)
• basedir::String: Base directory for output (default: current working directory)
• mmseqsdb::String: Path to MMseqs2 UniRef50 database (default: joinpath(homedir(), "workspace/mmseqs/UniRef50"))
• threads::Int: Number of CPU threads to use (default: all available). Note: This argument is not explicitly used by Pyrodigal or MMseqs2 in this version of the function, they might use their own defaults or require modifications to run_pyrodigal or run_mmseqs_easy_search to respect it.

Processing Steps

1. Creates output directory and copies input FASTA.
2. Runs Pyrodigal for gene prediction (nucleotide, amino acid, and GFF output).
3. Performs MMseqs2 homology search against UniRef50.
4. Predicts terminators using TransTerm.
5. Combines annotations into a unified GFF file.
6. Generates GenBank format output.

Returns

• String: Path to the output directory containing all generated files.

Files Generated (within the output directory specified by identifier)

• (basename(fasta)).pyrodigal.fna: Predicted genes (nucleotide) from Pyrodigal.
• (basename(fasta)).pyrodigal.faa: Predicted proteins from Pyrodigal.
• (basename(fasta)).pyrodigal.gff: Pyrodigal GFF annotations.
• (basename(fasta)).gff: Combined GFF annotations (MMseqs2 and TransTerm).
• (basename(fasta)).gff.genbank: Final GenBank format from the first combined GFF.
• (basename(fasta)).transterm_raw.gff: Combined GFF (MMseqs2 and a second TransTerm run).
• (basename(fasta)).transterm_raw.gff.genbank: Final GenBank format from the second combined GFF.

apply_learned_policy(policy::DQNPolicy, input_fastq::String; output_dir="rl_assembly")

Apply a trained RL policy to perform genome assembly.

This function uses a trained policy to make autonomous assembly decisions.

Arguments

• policy::DQNPolicy: Trained assembly policy
• input_fastq::String: Path to input FASTQ file
• output_dir::String: Output directory for assembly results

Returns

• Dict{String, Any}: Assembly results and metadata

Example

results = apply_learned_policy(trained_policy, "genome.fastq")

Check if two bubbles are equivalent.

assemble_genome(reads; method=StringGraph, config=AssemblyConfig()) -> AssemblyResult

Unified genome assembly interface using Phase 2 next-generation algorithms.

Arguments

• reads: Vector of FASTA/FASTQ records or file paths
• method: Assembly strategy (StringGraph, KmerGraph, HybridOLC, MultiK)
• config: Assembly configuration parameters

Returns

• AssemblyResult: Structure containing contigs, names, and assembly metadata

Details

This is the main entry point for the unified assembly pipeline, leveraging:

• Phase 1: MetaGraphsNext strand-aware graph construction
• Phase 2: Probabilistic algorithms, enhanced Viterbi, and graph algorithms
• Phase 3: Integrated workflow with polishing and validation

Examples

# Basic assembly with default parameters
reads = load_fastq_records("reads.fastq")
result = assemble_genome(reads)

# Custom assembly with specific k-mer size and error rate
config = AssemblyConfig(k=25, error_rate=0.005, polish_iterations=5)
result = assemble_genome(reads; method=KmerGraph, config=config)

# Access results
contigs = result.contigs
stats = result.assembly_stats

Assemble strings by performing graph walks from each connected component.

Generate a summary report of the assembly process.

assess_alignment(
    a,
    b
) -> @NamedTuple{total_matches::Int64, total_edits::Int64}

Aligns two sequences using the Levenshtein distance and returns the total number of matches and edits.

Arguments

• a::AbstractString: The first sequence to be aligned.
• b::AbstractString: The second sequence to be aligned.

Returns

• NamedTuple{(:total_matches, :total_edits), Tuple{Int, Int}}: A named tuple containing:
  • total_matches::Int: The total number of matching bases in the alignment.
  • total_edits::Int: The total number of edits (insertions, deletions, substitutions) in the alignment.

assess_alignment_accuracy(alignment_result) -> Any

Return proportion of matched bases in alignment to total matches + edits.

Calculate the accuracy of a sequence alignment by computing the ratio of matched bases to total alignment operations (matches + edits).

Arguments

• alignment_result: Alignment result object containing total_matches and total_edits fields

Returns

Float64 between 0.0 and 1.0 representing alignment accuracy, where:

• 1.0 indicates perfect alignment (all matches)
• 0.0 indicates no matches

assess_assembly_kmer_quality(; assembly, observations, ks)

Evaluate genome assembly quality by comparing k-mer distributions between assembled sequences and raw observations.

Arguments

• assembly: Input assembled sequences to evaluate
• observations: Raw sequencing data for comparison
• ks::Vector{Int}: Vector of k-mer sizes to analyze (default: k=17 to 23)

Returns

DataFrame containing quality metrics for each k-mer size:

• k: K-mer length used
• cosine_distance: Cosine similarity between k-mer distributions
• js_divergence: Jensen-Shannon divergence between distributions
• qv: MerQury-style quality value score

assess_assembly_quality(contigs_file)

Assess basic assembly quality metrics from a FASTA file.

Calculates standard assembly quality metrics including contig count, total length, and N50 statistic for assembly evaluation.

Arguments

• contigs_file: Path to FASTA file containing assembly contigs

Returns

• Tuple of (ncontigs, totallength, n50, l50)
  • n_contigs: Number of contigs in the assembly
  • total_length: Total length of all contigs in base pairs
  • n50: N50 statistic (length of shortest contig in the set covering 50% of assembly)
  • l50: L50 statistic (number of contigs needed to reach 50% of assembly length)

Example

n_contigs, total_length, n50, l50 = assess_assembly_quality("assembly.fasta")
println("Assembly has $n_contigs contigs, $total_length bp total, N50=$n50, L50=$l50")

See Also

• assess_assembly_kmer_quality: For k-mer based assembly quality assessment

assess_dnamer_saturation(
    fastx::AbstractString;
    power,
    outdir,
    min_k,
    max_k,
    threshold,
    kmers_to_assess
)

Analyzes k-mer saturation in a FASTA/FASTQ file to determine optimal k-mer size.

Arguments

• fastx::AbstractString: Path to input FASTA/FASTQ file
• power::Int=10: Exponent for downsampling k-mers (2^power)
• outdir::String="": Output directory for results. Uses current directory if empty
• min_k::Int=3: Minimum k-mer size to evaluate
• max_k::Int=17: Maximum k-mer size to evaluate
• threshold::Float64=0.1: Saturation threshold for k-mer assessment
• kmers_to_assess::Int=10_000_000: Maximum number of k-mers to sample

Returns

Dict{Int,Float64}: Dictionary mapping k-mer sizes to their saturation scores

assess_dnamer_saturation(
    fastxs::AbstractVector{<:AbstractString},
    kmer_type;
    kmers_to_assess,
    power,
    min_count
) -> Union{@NamedTuple{sampling_points::Vector{Int64}, unique_kmer_counts::Vector{Int64}}, NamedTuple{(:sampling_points, :unique_kmer_counts, :eof), <:Tuple{Vector, Vector{Int64}, Bool}}}

Assess k-mer saturation in DNA sequences from FASTX files.

Arguments

• fastxs::AbstractVector{<:AbstractString}: Vector of paths to FASTA/FASTQ files
• kmer_type: Type of k-mer to analyze (e.g., DNAKmer{21})
• kmers_to_assess=Inf: Maximum number of k-mers to process
• power=10: Base for exponential sampling intervals
• min_count=1: Minimum count threshold for considering a k-mer

Returns

Named tuple containing:

• sampling_points::Vector{Int}: K-mer counts at which samples were taken
• unique_kmer_counts::Vector{Int}: Number of unique canonical k-mers at each sampling point
• eof::Bool: Whether the entire input was processed

Details

Analyzes k-mer saturation by counting unique canonical k-mers at exponentially spaced intervals (powers of power). Useful for assessing sequence complexity and coverage. Returns early if all possible k-mers are observed.

assess_dnamer_saturation(
    fastxs::AbstractVector{<:AbstractString};
    power,
    outdir,
    min_k,
    max_k,
    threshold,
    kmers_to_assess,
    plot
)

Analyze k-mer saturation in DNA sequences to determine optimal k value.

Arguments

• fastxs: Vector of paths to FASTA/FASTQ files to analyze
• power: Base of logarithmic sampling points (default: 10)
• outdir: Optional output directory for plots and results
• min_k: Minimum k-mer size to test (default: 7)
• max_k: Maximum k-mer size to test (default: 17)
• threshold: Saturation threshold to determine optimal k (default: 0.1)
• kmers_to_assess: Maximum number of k-mers to sample (default: 10M)
• plot: Whether to generate saturation curves (default: true)

Returns

Integer representing the first k value that achieves saturation below threshold. If no k value meets the threshold, returns the k with minimum saturation.

Details

• Tests only prime k values between mink and maxk
• Generates saturation curves using logarithmic sampling
• Fits curves to estimate maximum unique k-mers
• If outdir is provided, saves plots as SVG and chosen k value to text file

assess_duplication_rates(fastq; results_table) -> Any

Analyze sequence duplication rates in a FASTQ file.

This function processes a FASTQ file to quantify both exact sequence duplications and canonical duplications (considering sequences and their reverse complements as equivalent). The function makes two passes through the file: first to count total records, then to analyze unique sequences.

Arguments

• fastq::String: Path to the input FASTQ file to analyze
• results_table::String: Optional. Path where the results will be saved as a tab-separated file. Defaults to the same path as the input file but with extension changed to ".duplication_rates.tsv"

Returns

• String: Path to the results table file

Output

Generates a tab-separated file containing the following metrics:

• total_records: Total number of sequence records in the file
• total_unique_observations: Count of unique sequence strings
• total_unique_canonical_observations: Count of unique canonical sequences (after normalizing for reverse complements)
• percent_unique_observations: Percentage of sequences that are unique
• percent_unique_canonical_observations: Percentage of sequences that are unique after canonicalization
• percent_duplication_rate: Percentage of sequences that are duplicates (100 - percentuniqueobservations)
• percent_canonical_duplication_rate: Percentage of sequences that are duplicates after canonicalization

Notes

• If the specified results file already exists and is not empty, the function will return early without recomputing.
• Progress is displayed during processing with a progress bar showing speed.

Example

# Analyze a FASTQ file and save results to default location
result_path = assess_duplication_rates("data/sample.fastq")

# Specify custom output path
result_path = assess_duplication_rates("data/sample.fastq", results_table="results/duplication_analysis.tsv")

assess_optimal_kmer_alignment(
    kmer,
    observed_kmer
) -> Tuple{@NamedTuple{total_matches::Int64, total_edits::Int64}, Union{Missing, Bool}}

Used to determine which orientation provides an optimal alignment for initiating path likelihood analyses in viterbi analysis

Compare alignment scores between a query k-mer and an observed k-mer in both forward and reverse complement orientations to determine optimal alignment.

Arguments

• kmer: Query k-mer sequence to align
• observed_kmer: Target k-mer sequence to align against

Returns

A tuple containing:

• alignment_result: The alignment result object for the optimal orientation
• orientation: Boolean indicating orientation (true = forward, false = reverse complement, missing = tied scores)

Details

• Performs pairwise alignment in both orientations using assess_alignment()
• Calculates accuracy scores using assess_alignment_accuracy()
• For tied alignment scores, randomly selects one orientation
• Uses BioSequences.reverse_complement for reverse orientation comparison

Attempt to correct a specific k-mer using probabilistic path finding. Returns true if correction was applied.

bam_to_fastq(; bam, fastq)

Convert a BAM file to FASTQ format with gzip compression.

Arguments

• bam: Path to input BAM file
• fastq: Optional output path. Defaults to input path with ".fq.gz" extension

Returns

• Path to the generated FASTQ file

Details

• Uses samtools through conda environment
• Automatically skips if output file exists
• Output is gzip compressed
• Requires samtools to be available via conda

bandage_visualize(; gfa, img)

Generate a visualization of a genome assembly graph using Bandage.

Arguments

• gfa: Path to input GFA (Graphical Fragment Assembly) file
• img: Optional output image path. Defaults to GFA filename with .png extension

Returns

• Path to the generated image file

batch_download_viroid_references(
    species_list::Vector{String};
    base_outdir,
    download_genome,
    download_cds,
    download_protein,
    max_per_species
) -> Dict{String, Any}

Download reference data for multiple viroid species in batch.

Arguments

• species_list::Vector{String}: List of viroid species to download
• base_outdir::String: Base output directory (subdirectories created per species)
• download_genome::Bool: Download genomic sequences (default: true)
• download_cds::Bool: Download CDS sequences (default: true)
• download_protein::Bool: Download protein sequences (default: true)
• max_per_species::Int: Maximum sequences per species per type (default: 5)

Returns

• Dict{String, NamedTuple}: Dictionary mapping species names to their downloaded file info

Examples

# Download data for all known viroid species
species_list = get_viroid_species_list()
results = batch_download_viroid_references(species_list, "viroid_database/")

# Download only genomes for a subset
pstv_like = ["Potato spindle tuber viroid", "Tomato planta macho viroid"]  
results = batch_download_viroid_references(pstv_like, "pstv_group/";
                                         download_cds=false, download_protein=false)

benchmark_graph_construction()
benchmark_graph_construction(
    config::Mycelia.BenchmarkConfig
)

Benchmark graph construction performance: Legacy vs Next-generation.

Compares:

• MetaGraphs.jl (legacy) vs MetaGraphsNext.jl (next-gen)
• Memory allocation patterns
• Construction time
• Type stability

Arguments

• config: BenchmarkConfig for test parameters

Returns

• NamedTuple with benchmark results

benchmark_memory_patterns()
benchmark_memory_patterns(config::Mycelia.BenchmarkConfig)

Benchmark memory usage patterns for different graph representations.

Compares memory usage of:

• Stranded vertices (legacy) vs canonical vertices (next-gen)
• Edge metadata structures
• Coverage tracking efficiency

Arguments

• config: BenchmarkConfig for test parameters

Returns

• NamedTuple with memory analysis results

benchmark_type_stability()
benchmark_type_stability(config::Mycelia.BenchmarkConfig)

Benchmark type stability and allocation patterns.

Measures:

• Type inference success
• Runtime allocations
• Performance predictability

Arguments

• config: BenchmarkConfig for test parameters

Returns

• NamedTuple with type stability metrics

bernoullipcaepca(M::AbstractMatrix{Bool}; k::Int=0)

Perform Bernoulli (logistic) EPCA on a 0/1 matrix M (features × samples).

When to use

Use for binary (0/1) data, such as presence/absence or yes/no features.

Returns

A NamedTuple with

• model : the fitted ExpFamilyPCA.BernoulliEPCA object
• scores : k×n_samples matrix of low‐dimensional sample scores
• loadings : k×n_features matrix of feature loadings

best_label_mapping(true_labels, pred_labels)

Finds the optimal mapping from predicted labels to true labels using the Hungarian algorithm, so that the total overlap (confusion matrix diagonal) is maximized. Returns the remapped predicted labels and the mapping as a Dict.

binary_matrix_to_jaccard_distance_matrix(binary_matrix::Union{BitMatrix, Matrix{Bool}})

Pairwise Jaccard distance between columns of a binary matrix (BitMatrix or Matrix{Bool}). Throws an error if the input is not strictly a binary matrix.

binomial_pca_epca(M::AbstractMatrix{<:Integer}; k::Int=0, ntrials::Int=1)

Perform Binomial EPCA on a count matrix M (features × samples).

When to use

Use for integer count data representing the number of successes out of a fixed number of trials (e.g., number of mutated alleles out of total alleles).

Keyword arguments

• k : desired number of latent dimensions; if k<1 defaults to min(n_samples-1, n_features, 10)
• ntrials : number of trials for the Binomial distribution (default=1)

Returns

NamedTuple with fields

• model : the fitted ExpFamilyPCA.BinomialEPCA object
• scores : k×n_samples matrix of sample scores
• loadings : k×n_features matrix of feature loadings

biosequences_to_counts_table(; biosequences, k)

Convert a collection of biological sequences into a k-mer count matrix.

Arguments

• biosequences: Vector of biological sequences (DNA, RNA, or Amino Acids)
• k: Length of k-mers to count

Returns

Named tuple with:

• sorted_kmers: Vector of all unique k-mers found, lexicographically sorted
• kmer_counts_matrix: Sparse matrix where rows are k-mers and columns are sequences

Details

• For DNA sequences, counts canonical k-mers (both strands)
• Uses parallel processing with Thread-safe progress tracking
• Memory efficient sparse matrix representation
• Supports DNA, RNA and Amino Acid sequences

biosequences_to_dense_counts_table(; biosequences, k)

Convert a collection of biological sequences into a dense k-mer count matrix.

Arguments

• biosequences: Collection of DNA, RNA, or amino acid sequences (BioSequence types)
• k::Integer: Length of k-mers to count (must be ≤ 13)

Returns

Named tuple containing:

• sorted_kmers: Vector of all possible k-mers in sorted order
• kmer_counts_matrix: Dense matrix where rows are k-mers and columns are sequences

Details

• For DNA sequences, counts canonical k-mers (both strands)
• For RNA and protein sequences, counts exact k-mers
• Uses parallel processing with threads

blastdb2table(
;
    blastdb,
    ALL_FIELDS,
    sequence_sha256,
    sequence_hash,
    sequence_id,
    accession,
    gi,
    sequence_title,
    blast_name,
    taxid,
    taxonomic_super_kingdom,
    scientific_name,
    scientific_names_leaf_nodes,
    common_taxonomic_name,
    common_names_leaf_nodes,
    leaf_node_taxids,
    membership_integer,
    ordinal_id,
    pig,
    sequence_length,
    sequence
)

Convert a BLAST database to an in-memory table with sequence and taxonomy information.

Arguments

• blastdb::String: Path to the BLAST database
• outfile::String="": Optional output file path. If provided, results will be saved to this file
• force::Bool=false: Whether to overwrite existing output file
• ALL_FIELDS::Bool=true: If true, include all fields regardless of other flag settings
• Field selection flags (default to false unless ALL_FIELDS is true):
  • sequence_sha256::Bool: Include SHA256 hash of the sequence
  • sequence_hash::Bool: Include sequence hash
  • sequence_id::Bool: Include sequence ID
  • accession::Bool: Include accession number
  • gi::Bool: Include GI number
  • sequence_title::Bool: Include sequence title
  • blast_name::Bool: Include BLAST name
  • taxid::Bool: Include taxid
  • taxonomic_super_kingdom::Bool: Include taxonomic super kingdom
  • scientific_name::Bool: Include scientific name
  • scientific_names_leaf_nodes::Bool: Include scientific names for leaf-node taxids
  • common_taxonomic_name::Bool: Include common taxonomic name
  • common_names_leaf_nodes::Bool: Include common taxonomic names for leaf-node taxids
  • leaf_node_taxids::Bool: Include leaf-node taxids
  • membership_integer::Bool: Include membership integer
  • ordinal_id::Bool: Include ordinal ID
  • pig::Bool: Include PIG
  • sequence_length::Bool: Include sequence length
  • sequence::Bool: Include the full sequence

Returns

• DataFrame: DataFrame containing the requested columns from the BLAST database

blastdb_to_fasta(
;
    blastdb,
    entries,
    taxids,
    outfile,
    force,
    max_cores
)

Convert a BLAST database to FASTA format.

Arguments

• blastdb::String: Name of the BLAST database to convert (e.g. "nr", "nt")
• dbdir::String: Directory containing the BLAST database files
• outfile::String: Path for the output FASTA file

Returns

• Path to the generated FASTA file as String

Compute the Bray-Curtis distance between columns of a count matrix.

Arguments

• M::AbstractMatrix{<:Integer}: Count matrix where rows are features and columns are samples

Returns

• Matrix{Float64}: Symmetric distance matrix with Bray-Curtis distances

Breadth-first sampling: sample at least one from each group, then sample remaining proportionally to group frequencies

Apply breadth-first sampling to a DataFrame

build_biosequence_graph(
    fasta_records::Vector{FASTX.FASTA.Record};
    graph_mode,
    min_overlap
) -> MetaGraphsNext.MetaGraph{Int64, Graphs.SimpleGraphs.SimpleDiGraph{Int64}, _A, _B, _C, Nothing, Mycelia.var"#267#268", Float64} where {_A, _B, _C}

Build a BioSequence graph directly from FASTA records.

This creates a variable-length BioSequence graph where vertices are the input sequences and edges represent relationships between sequences (e.g., overlaps, containment).

Arguments

• fasta_records: Vector of FASTA records
• graph_mode: SingleStrand or DoubleStrand mode (default: DoubleStrand)
• min_overlap: Minimum overlap length for creating edges (default: 100)

Returns

• MetaGraphsNext.MetaGraph with BioSequence vertices and overlap edges

Example

records = [FASTX.FASTA.Record("seq1", "ATCGATCGATCG"), 
           FASTX.FASTA.Record("seq2", "CGATCGATCGAA")]
graph = build_biosequence_graph(records)

build_directed_kmer_graph(; fastq, k, plot)

Constructs a directed graph representation of k-mer transitions from FASTQ sequencing data.

Arguments

• fastq: Path to input FASTQ file
• k: K-mer size (default: 1). Must be odd and prime. If k=1, optimal size is auto-determined
• plot: Boolean to display quality distribution plot (default: false)

Returns

MetaDiGraph with properties:

• assembly_k: k-mer size used
• kmer_counts: frequency of each k-mer
• transition_likelihoods: edge weights between k-mers
• kmermeanquality, kmertotalquality: quality metrics
• branchingnodes, unbranchingnodes: topological classification
• likelyvalidkmer_indices: k-mers above mean quality threshold
• likelysequencingartifact_indices: potential erroneous k-mers

Note

For DNA assembly, quality scores are normalized across both strands.

build_genome_distance_matrix(genome_files::Vector{String}; kmer_type=Kmers.DNAKmer{21}, metric=:js_divergence)

Build a distance matrix between all genome pairs using existing distance metrics.

Creates a comprehensive pairwise distance matrix using established k-mer distance functions, suitable for phylogenetic analysis and clustering.

Arguments

• genome_files: Vector of genome FASTA file paths
• kmer_type: K-mer type from Kmers.jl (default: Kmers.DNAKmer{21})
• metric: Distance metric (:js_divergence, :cosine, :jaccard)

Returns

• Named tuple with distance matrix and genome names

Example

genomes = ["genome1.fasta", "genome2.fasta", "genome3.fasta"]
result = Mycelia.build_genome_distance_matrix(genomes, kmer_type=Kmers.DNAKmer{31})
println("Distance matrix: $(result.distance_matrix)")

build_kmer_graph_next(
    kmer_type,
    observations::AbstractVector{<:Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}};
    graph_mode
) -> MetaGraphsNext.MetaGraph{Int64, Graphs.SimpleGraphs.SimpleDiGraph{Int64}, Label, VertexData, EdgeData, Nothing, WeightFunction, Float64} where {Label, VertexData, EdgeData, WeightFunction}

Create a next-generation, type-stable k-mer graph using MetaGraphsNext.

This implementation uses canonical k-mers as vertices with strand-aware edges that respect biological transition constraints for both single-strand and double-strand sequences.

Arguments

• kmer_type: Type of k-mer (e.g., DNAKmer{K})
• observations: Vector of FASTA/FASTQ records
• graph_mode: SingleStrand for directional sequences, DoubleStrand for DNA (default)

Returns

• MetaGraphsNext.MetaGraph with canonical vertices and strand-aware edges

Details

• Vertices: Always canonical k-mers (lexicographically smaller of kmer/reverse_complement)
• Edges: Strand-aware transitions that respect biological constraints
• SingleStrand mode: Only forward-strand transitions allowed
• DoubleStrand mode: Both forward and reverse-complement transitions allowed

build_quality_biosequence_graph(
    fastq_records::Vector{FASTX.FASTQ.Record};
    graph_mode,
    min_overlap,
    min_quality
) -> MetaGraphsNext.MetaGraph{Int64, Graphs.SimpleGraphs.SimpleDiGraph{Int64}, _A, _B, _C, Nothing, Mycelia.var"#277#278", Float64} where {_A, _B, _C}

Build a quality-aware BioSequence graph directly from FASTQ records.

This creates a variable-length quality-aware BioSequence graph where vertices are the input sequences with their quality scores and edges represent quality-weighted relationships between sequences.

Arguments

• fastq_records: Vector of FASTQ records with quality scores
• graph_mode: SingleStrand or DoubleStrand mode (default: DoubleStrand)
• min_overlap: Minimum overlap length for creating edges (default: 100)
• min_quality: Minimum mean quality to include sequence (default: 20)

Returns

• MetaGraphsNext.MetaGraph with quality-aware BioSequence vertices

Example

records = [FASTX.FASTQ.Record("seq1", "ATCGATCGATCG", "IIIIIIIIIIII")]
graph = build_quality_biosequence_graph(records)

build_qualmer_graph(
    fastq_records::Vector{FASTX.FASTQ.Record};
    k,
    graph_mode,
    min_quality,
    min_coverage
) -> MetaGraphsNext.MetaGraph{Int64, Graphs.SimpleGraphs.SimpleDiGraph{Int64}, Label, VertexData, EdgeData, Nothing, WeightFunction, Float64} where {Label, VertexData, EdgeData, WeightFunction}

Build a quality-aware k-mer graph from FASTQ records using existing Qualmer functionality. This function leverages the existing qualmer extraction functions and adds graph construction.

Arguments

• fastq_records: Vector of FASTQ records with quality scores
• k: K-mer size
• graph_mode: SingleStrand or DoubleStrand mode (default: DoubleStrand)
• min_quality: Minimum average PHRED quality to include k-mer (default: 10)
• min_coverage: Minimum coverage (observations) to include k-mer (default: 1)

Returns

• MetaGraphsNext.MetaGraph with QualmerVertexData and QualmerEdgeData

build_stranded_kmer_graph(
    kmer_type,
    observations::AbstractVector{<:Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}}
) -> MetaGraphs.MetaDiGraph{T, Float64} where T<:Integer

Create a weighted, strand-specific kmer (de bruijn) graph from a set of kmers and a series of sequence observations in FASTA format.

Calculate assembly accuracy metrics for reward function. Returns a comprehensive score based on multiple quality indicators.

Calculate assembly mapping metrics for validation reads.

calculate_assembly_quality_metrics(
    qualmer_graph;
    low_confidence_threshold
) -> NamedTuple{(:mean_coverage, :mean_quality, :mean_confidence, :low_confidence_fraction, :total_kmers), <:NTuple{5, Any}}

Calculate comprehensive assembly quality metrics for a qualmer graph.

Arguments

• graph: Qualmer graph (MetaGraphsNext with QualmerVertexData)
• low_confidence_threshold::Float64=0.95: Threshold for identifying low-confidence k-mers

Returns

• NamedTuple: Assembly quality metrics including coverage, quality, and confidence statistics

Details

Calculates mean values for coverage, quality scores, and joint probabilities. Identifies fraction of k-mers below confidence threshold as potential error indicators.

Calculate reward for current k-mer processing iteration. Higher rewards indicate better assembly quality progress.

Calculate complexity score for a bubble.

Calculate in-degrees and out-degrees for all vertices.

calculate_emission_probability(state::ViterbiState, observation::String, config::ViterbiConfig) -> Float64

Calculate emission probability for a state given an observation.

calculate_gc_content(sequence::AbstractString) -> Float64

Calculate GC content from a string sequence.

Convenience function that accepts string input and converts to appropriate BioSequence. Automatically detects DNA/RNA based on presence of T/U.

Arguments

• sequence::AbstractString: Input DNA or RNA sequence as string

Returns

• Float64: GC content as a percentage (0.0-100.0)

Examples

# Calculate GC content from string
gc_percent = calculate_gc_content("ATCGATCGATCG")

calculate_gc_content(
    sequence::BioSequences.LongSequence
) -> Float64

Calculate GC content (percentage of G and C bases) from a BioSequence.

This function calculates the percentage of guanine (G) and cytosine (C) bases in a nucleotide sequence. Works with both DNA and RNA sequences.

Arguments

• sequence::BioSequences.LongSequence: Input DNA or RNA sequence

Returns

• Float64: GC content as a percentage (0.0-100.0)

Examples

# Calculate GC content for DNA
dna_seq = BioSequences.LongDNA{4}("ATCGATCGATCG")
gc_percent = calculate_gc_content(dna_seq)

# Calculate GC content for RNA
rna_seq = BioSequences.LongRNA{4}("AUCGAUCGAUCG") 
gc_percent = calculate_gc_content(rna_seq)

calculate_gc_content(
    records::AbstractArray{T<:Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}, 1}
) -> Float64

Calculate GC content from FASTA/FASTQ records.

Processes multiple records and calculates overall GC content across all sequences.

Arguments

• records::AbstractVector{T} where T <: Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}`: Input records

Returns

• Float64: Overall GC content as a percentage (0.0-100.0)

Examples

# Calculate GC content from FASTA records
records = collect(FASTX.FASTA.Reader(open("sequences.fasta")))
gc_percent = calculate_gc_content(records)

Calculate support for a path (simplified version).

Calculate likelihood of a k-mer given its observed quality scores and qualmer vertex data. Leverages existing qualmer graph quality calculations.

Calculate likelihood of a FASTQ read given the current graph.

Calculate confidence in recommendation based on variance across folds.

Calculate confidence in repeat identification.

Calculate quality-aware likelihood of a sequence given the qualmer graph. Uses both k-mer presence and quality-based confidence from qualmer observations.

Calculate k-mer sparsity for a given k-mer size. Returns fraction of possible k-mers that are NOT observed. Higher sparsity indicates errors become singletons.

Calculate penalty for unseen k-mers based on their quality scores. High quality unseen k-mers get less penalty than low quality ones.

canonical(
    qmer::Mycelia.Qualmer{<:Kmers.Kmer{BioSequences.AminoAcidAlphabet}}
) -> Mycelia.Qualmer{<:Kmers.Kmer{BioSequences.AminoAcidAlphabet}}

canonical(
    qmer::Mycelia.Qualmer{KmerT<:Union{Kmers.DNAKmer, Kmers.RNAKmer}}
) -> Mycelia.Qualmer{KmerT} where KmerT<:Kmers.Kmer

Get the canonical representation of a DNA qualmer, considering both sequence and quality. For DNA/RNA, this involves potentially reverse-complementing the k-mer and reversing the quality scores accordingly.

canonicalize_kmer_counts!(kmer_counts) -> Any

Canonicalizes the k-mer counts in the given dictionary.

This function iterates over the provided dictionary kmer_counts, which maps k-mers to their respective counts. For each k-mer that is not in its canonical form, it converts the k-mer to its canonical form and updates the count in the dictionary accordingly. If the canonical form of the k-mer already exists in the dictionary, their counts are summed. The original non-canonical k-mer is then removed from the dictionary.

Arguments

• kmer_counts::Dict{BioSequences.Kmer, Int}: A dictionary where keys are k-mers and values are their counts.

Returns

• The input dictionary kmer_counts with all k-mers in their canonical form, sorted by k-mers.

canonicalize_kmer_counts(kmer_counts) -> Any

Normalize k-mer counts into a canonical form by creating a non-mutating copy.

Arguments

• kmer_counts: Dictionary or collection of k-mer count data

Returns

• A new normalized k-mer count collection

check_bioconda_env_is_installed(pkg) -> Bool

Check whether a named Bioconda environment already exists.

Arguments

• pkg::String: Name of the environment.

Returns

Bool indicating if the environment is present.

Check conditions for Eulerian path existence.

check_matrix_fits_in_memory(bytes_needed::Integer; severity::Symbol=:warn)

Checks whether the specified number of bytes can fit in the computer's memory.

• bytes_needed: The number of bytes required (output from estimate_dense_matrix_memory or estimate_sparse_matrix_memory).
• severity: What to do if there is not enough available memory. Can be :warn (default) or :error.

Returns a named tuple: (willfittotal, willfitavailable, totalmemory, freememory, bytes_needed) Where:

• will_fit_total: true if the matrix fits in total system memory.
• will_fit_available: true if the matrix fits in currently available (free) system memory.
• total_memory: Total system RAM in bytes.
• free_memory: Currently available system RAM in bytes.
• bytes_needed: Bytes requested for the matrix.

If will_fit_available is false, either warns or errors depending on severity.

Check if memory usage is within acceptable limits.

choose_top_n_markers(N)

Return a vector of the top N most visually distinct marker symbols for plotting.

Arguments

• N::Int: Number of distinct markers to return (max 17 for best differentiation).

Returns

• Vector{Symbol}: Vector of marker symbol names.

Example

markers = choose_top_n_markers(7)

chromosome_coverage_table_to_plot(cdf) -> Plots.Plot

Creates a visualization of chromosome coverage data with statistical thresholds.

Arguments

• cdf::DataFrame: Coverage data frame containing columns:
  • index: Chromosome position indices
  • depth: Coverage depth values
  • chromosome: Chromosome identifier
  • mean_coverage: Mean coverage value
  • std_coverage: Standard deviation of coverage
  • 3σ: Boolean vector indicating +3 sigma regions
  • -3σ: Boolean vector indicating -3 sigma regions

Returns

• A StatsPlots plot object showing:
  • Raw coverage data (black line)
  • Mean coverage and ±1,2,3σ thresholds (rainbow colors)
  • Highlighted regions exceeding ±3σ thresholds (red vertical lines)

Classify reads based on taxonomic alignments using individual alignment scoring.

This function takes taxonomy-aware alignment data and performs classification by:

1. Loading the taxonomy-aware alignment data
2. Analyzing alignment score distributions per read
3. Identifying dominant taxonomic assignments
4. Applying conservative taxonomy classification

Arguments

• taxonomy_aware_file: Path to taxonomy-aware alignment file (.tsv.gz or .arrow)
• min_relative_proportion::Float64=60.0: Minimum relative proportion threshold for accepting a taxonomic assignment
• verbose::Bool=true: Whether to print progress information

Returns

A DataFrame with taxonomic classification results including individual alignment metrics

Classify the type of repeat.

cleanup_directory(
    directory::AbstractString;
    verbose,
    force
) -> @NamedTuple{existed::Bool, files_deleted::Int64, bytes_freed::Int64, human_readable_size::String}

Clean up a directory by calculating its size and file count, then removing it.

Arguments

• directory::AbstractString: Path to the directory to clean up
• verbose::Bool=true: Whether to report cleanup results (default: true)
• force::Bool=false: Whether to proceed without confirmation for large directories

Returns

• Named tuple with fields:
  • existed: Whether the directory existed before cleanup
  • files_deleted: Number of files that were deleted
  • bytes_freed: Total bytes freed up
  • human_readable_size: Human-readable representation of bytes freed

Details

This function will:

1. Check if the directory exists and is non-empty
2. Calculate the total size and number of files recursively
3. Remove the directory and all its contents
4. Report the cleanup results unless verbose=false

For directories larger than 1GB or containing more than 10,000 files, confirmation is required unless force=true.

Examples

# Clean up a temporary directory with reporting
result = cleanup_directory("/tmp/myapp_temp")

# Silent cleanup
cleanup_directory("/tmp/cache", verbose=false)

# Force cleanup of large directory
cleanup_directory("/tmp/large_data", force=true)

codon_optimize(
;
    normalized_codon_frequencies,
    protein_sequence,
    n_iterations
)

Optimizes the DNA sequence encoding for a given protein sequence using codon usage frequencies.

Arguments

• normalized_codon_frequencies: Dictionary mapping amino acids to their codon frequencies
• protein_sequence::BioSequences.LongAA: Target protein sequence to optimize
• n_iterations::Integer: Number of optimization iterations to perform

Algorithm

1. Creates initial DNA sequence through reverse translation
2. Iteratively generates new sequences by sampling codons based on their frequencies
3. Keeps track of the sequence with highest codon usage likelihood

Returns

• BioSequences.LongDNA{2}: Optimized DNA sequence encoding the input protein

codons_to_amino_acids() -> Dict

Creates a mapping between DNA codons and their corresponding amino acids using the standard genetic code.

Returns a dictionary where:

• Keys are 3-letter DNA codons (e.g., "ATG")
• Values are the corresponding amino acids from BioSequences.jl

Collapse linear paths (vertices with one incoming and one outgoing edge) into a simpler graph where sequences are concatenated.

Compare statistical properties of intelligent vs iterative assemblies.

compare_genome_kmer_similarity(genome1_file::String, genome2_file::String; kmer_type=Kmers.DNAKmer{21}, metric=:js_divergence)

Compare two genomes using existing k-mer distance metrics.

Leverages existing distance metric functions to compare genomic similarity between pairs of genomes using various distance measures.

Arguments

• genome1_file: Path to first genome FASTA file
• genome2_file: Path to second genome FASTA file
• kmer_type: K-mer type from Kmers.jl (default: Kmers.DNAKmer{21})
• metric: Distance metric (:js_divergence, :cosine, :jaccard)

Returns

• Named tuple with similarity/distance metrics and k-mer statistics

Example

similarity = Mycelia.compare_genome_kmer_similarity(
    "genome1.fasta", "genome2.fasta", 
    kmer_type=Kmers.DNAKmer{31}, 
    metric=:js_divergence
)
println("JS divergence: $(similarity.distance)")
println("Shared k-mers: $(similarity.shared_kmers)")

concatenate_files(; files, file)

Join fasta files without any regard to record uniqueness.

A cross-platform version of cat *.fasta > joint.fasta

See mergefastafiles

Concatenate multiple FASTA files into a single output file by simple appending.

Arguments

• files: Vector of paths to input FASTA files
• file: Path where the concatenated output will be written

Returns

• Path to the output concatenated file

Details

Platform-independent implementation of cat *.fasta > combined.fasta. Files are processed sequentially with a progress indicator.

confusion_matrix(true_labels, pred_labels)

Returns the confusion matrix as a Matrix{Int}, row = true, col = predicted. Also returns the list of unique labels in sorted order and a heatmap plot.

contbernoulli_pca_epca(M::AbstractMatrix{<:Real}; k::Int=0)

Perform Continuous Bernoulli EPCA on a matrix M (features × samples).

When to use

Use for continuous data in the open interval (0, 1), such as probabilities or normalized intensities.

Keyword arguments

• k : desired number of latent dimensions; if k<1 defaults to min(n_samples-1, n_features, 10)

Returns

NamedTuple with fields

• model : the fitted ExpFamilyPCA.ContinuousBernoulliEPCA object
• scores : k×n_samples matrix of sample scores
• loadings : k×n_features matrix of feature loadings

contig_is_circular(
    graph_file::String,
    contig_name::String
) -> Any

Returns bool indicating whether the contig is a circle

graphfile = path to assembly graph.gfa file contigname = name of the contig

Determine if a contig represents a circular structure in the assembly graph.

A circular contig is one where the sequence forms a complete loop in the assembly graph, typically representing structures like plasmids, circular chromosomes, or other circular DNA elements.

Arguments

• graph_file::String: Path to the assembly graph in GFA format
• contig_name::String: Name/identifier of the contig to check

Returns

• Bool: true if the contig forms a circular structure, false otherwise

contig_is_cleanly_assembled(
    graph_file::String,
    contig_name::String
) -> Bool

Returns bool indicating whether the contig is cleanly assembled.

By cleanly assembled we mean that the contig does not have other contigs attached in the same connected component.

graphfile = path to assembly graph.gfa file contigname = name of the contig

Check if a contig exists in isolation within its connected component in an assembly graph.

Arguments

• graph_file::String: Path to the assembly graph file in GFA format
• contig_name::String: Name/identifier of the contig to check

Returns

• Bool: true if the contig exists alone in its connected component, false otherwise

Details

A contig is considered "cleanly assembled" if it appears as a single entry in the assembly graph's connected components. This function parses the GFA file and checks the contig's isolation status using the graph structure.

convert_legacy_gfa_to_next(
    gfa_file::AbstractString
) -> MetaGraphsNext.MetaGraph{Int64, Graphs.SimpleGraphs.SimpleDiGraph{Int64}, _A, _B, _C, Nothing, Mycelia.var"#147#148", Float64} where {_A, _B, _C}
convert_legacy_gfa_to_next(
    gfa_file::AbstractString,
    graph_mode::Mycelia.GraphMode
) -> MetaGraphsNext.MetaGraph{Int64, Graphs.SimpleGraphs.SimpleDiGraph{Int64}, _A, _B, _C, Nothing, Mycelia.var"#147#148", Float64} where {_A, _B, _C}

Convert a legacy MetaGraphs GFA to next-generation MetaGraphsNext format.

This convenience function reads a GFA file using the legacy parser and converts it to the new strand-aware format.

Arguments

• gfa_file: Path to GFA file
• graph_mode: GraphMode for the output graph

Returns

• MetaGraphsNext.MetaGraph in strand-aware format

convert_sequence(seq::AbstractString) -> Any

Converts the given sequence (output from FASTX.sequence) into the appropriate BioSequence type:

• DNA sequences are converted using BioSequences.LongDNA
• RNA sequences are converted using BioSequences.LongRNA
• AA sequences are converted using BioSequences.LongAA

copy_to_tempdir(file_path::String) -> String

Create a copy of a file in a temporary directory while preserving the original filename.

Arguments

• file_path::String: Path to the source file to be copied

Returns

• String: Path to the newly created temporary file

copy_with_unique_identifier(
    infile,
    out_directory,
    unique_identifier;
    force
) -> Any

Copy a file to a new location with a unique identifier prepended to the filename.

Arguments

• infile::AbstractString: Path to the source file to copy
• out_directory::AbstractString: Destination directory for the copied file
• unique_identifier::AbstractString: String to prepend to the filename
• force::Bool=true: If true, overwrite existing files

Returns

• String: Path to the newly created file

Perform error correction at the current k-mer size. Returns the number of corrections made.

correct_errors_next(graph::MetaGraph, sequences::Vector, config::ViterbiConfig) -> Vector{FASTX.FASTA.Record}

Correct errors in sequences using Viterbi algorithm and return corrected FASTA records.

count_canonical_kmers(_::Type{KMER_TYPE}, sequences) -> Any

Count canonical k-mers in biological sequences. A canonical k-mer is the lexicographically smaller of a DNA sequence and its reverse complement, ensuring strand-independent counting.

Arguments

• KMER_TYPE: Type parameter specifying the k-mer size and structure
• sequences: Iterator of biological sequences to analyze

Returns

• Dict{KMER_TYPE,Int}: Dictionary mapping canonical k-mers to their counts

count_kmers(
    _::Type{KMER_TYPE},
    fastx_file::AbstractString
) -> Any

Count k-mers in a FASTA/FASTQ file and return their frequencies.

Arguments

• KMER_TYPE: Type parameter specifying the k-mer type (e.g., DNAKmer{K})
• fastx_file: Path to input FASTA/FASTQ file

Returns

• Dict{KMER_TYPE, Int}: Dictionary mapping each k-mer to its frequency

count_kmers(
    _::Type{Kmers.Kmer{A<:BioSequences.AminoAcidAlphabet, K}},
    sequence::BioSequences.LongSequence
) -> OrderedCollections.OrderedDict{K, Int64} where K<:(Kmers.Kmer{BioSequences.AminoAcidAlphabet, _A, _B} where {_B, _A})

Count the frequency of amino acid k-mers in a biological sequence.

Arguments

• Kmers.Kmer{A,K}: Type parameter specifying amino acid alphabet (A) and k-mer length (K)
• sequence: Input biological sequence to analyze

Returns

A sorted dictionary mapping each k-mer to its frequency count in the sequence.

count_kmers(
    _::Type{Kmers.Kmer{A<:BioSequences.DNAAlphabet, K}},
    sequence::BioSequences.LongSequence
) -> Any

Count the frequency of each k-mer in a DNA sequence.

Arguments

• ::Type{Kmers.Kmer{A,K}}: K-mer type with alphabet A and length K
• sequence::BioSequences.LongSequence: Input DNA sequence to analyze

Returns

A sorted dictionary mapping each k-mer to its frequency count in the sequence.

Type Parameters

• A <: BioSequences.DNAAlphabet: DNA alphabet type
• K: Length of k-mers

count_kmers(
    _::Type{Kmers.Kmer{A<:BioSequences.RNAAlphabet, K}},
    sequence::BioSequences.LongSequence
) -> Any

Count the frequency of each k-mer in an RNA sequence.

Arguments

• Kmer: Type parameter specifying the k-mer length K and RNA alphabet
• sequence: Input RNA sequence to analyze

Returns

• Dict{Kmers.Kmer, Int}: Sorted dictionary mapping each k-mer to its frequency count

count_kmers(
    _::Type{KMER_TYPE},
    sequences::Union{FASTX.FASTA.Reader, FASTX.FASTQ.Reader}
) -> Any

Counts k-mer occurrences in biological sequences from a FASTA/FASTQ reader.

Arguments

• KMER_TYPE: Type parameter specifying the k-mer length and encoding (e.g., DNAKmer{4} for 4-mers)
• sequences: A FASTA or FASTQ reader containing the biological sequences to analyze

Returns

A dictionary mapping k-mers to their counts in the input sequences

count_kmers(
    _::Type{KMER_TYPE},
    record::Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}
) -> Any

Count the frequency of amino acid k-mers in a biological sequence.

Arguments

• Kmers.Kmer{A,K}: Type parameter specifying amino acid alphabet (A) and k-mer length (K)
• sequence: Input biological sequence to analyze

Returns

A sorted dictionary mapping each k-mer to its frequency count in the sequence.

count_kmers(
    _::Type{KMER_TYPE},
    fastx_files::AbstractArray{T<:AbstractString, 1}
) -> Any

Count k-mers across multiple FASTA/FASTQ files and merge the results.

Arguments

• KMER_TYPE: Type parameter specifying the k-mer length (e.g., DNAKmer{4} for 4-mers)
• fastx_files: Vector of paths to FASTA/FASTQ files

Returns

• Dict{KMER_TYPE, Int}: Dictionary mapping k-mers to their total counts across all files

count_kmers(
    _::Type{KMER_TYPE},
    records::AbstractArray{T<:Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}, 1}
) -> Any

Count k-mers across multiple sequence records and return a sorted frequency table.

Arguments

• KMER_TYPE: Type parameter specifying the k-mer length (e.g., DNAKmer{3} for 3-mers)
• records: Vector of FASTA/FASTQ records to analyze

Returns

• Dict{KMER_TYPE, Int}: Sorted dictionary mapping k-mers to their frequencies

count_matrix_to_probability_matrix(counts_matrix) -> Any

Convert a matrix of counts into a probability matrix by normalizing each column to sum to 1.0.

Arguments

• counts_matrix::Matrix{<:Number}: Input matrix where each column represents counts/frequencies

Returns

• Matrix{Float64}: Probability matrix where each column sums to 1.0

count_predicted_genes(gff_file)

Count the number of predicted genes from a GFF file.

Parses a GFF/GTF file and counts the number of CDS (coding sequence) features, which correspond to predicted genes.

Arguments

• gff_file: Path to GFF/GTF file

Returns

• Integer count of predicted genes (CDS features)

See Also

• run_pyrodigal: For gene prediction that generates GFF files
• parse_transterm_output: For parsing other annotation tool outputs

count_records(fastx) -> Int64

Counts the total number of records in a FASTA/FASTQ file.

Arguments

• fastx: Path to a FASTA or FASTQ file (can be gzipped)

Returns

• Number of records (sequences) in the file

countmap_columns(table)

Generate and display frequency counts for all columns in a DataFrame.

Arguments

• table::DataFrame: Input DataFrame to analyze

Details

Iterates through each column in the DataFrame and displays:

1. The column name
2. A Dict mapping unique values to their frequencies using StatsBase.countmap

create_assembly_environment(training_data, validation_data; episode_length=100)

Create a new reinforcement learning environment for assembly training.

Arguments

• training_data::Vector{String}: Paths to training FASTQ datasets
• validation_data::Vector{String}: Paths to validation FASTQ datasets
• episode_length::Int: Maximum steps per training episode (default: 100)

Returns

• AssemblyEnvironment: Initialized RL environment ready for training

Example

training_files = ["genome1.fastq", "genome2.fastq", "genome3.fastq"]
validation_files = ["validation1.fastq", "validation2.fastq"]
env = create_assembly_environment(training_files, validation_files)

create_curriculum_schedule(; stages=4, datasets_per_stage=10)

Create a curriculum learning schedule that progressively increases difficulty.

Arguments

• stages::Int: Number of curriculum stages (default: 4)
• datasets_per_stage::Int: Datasets per stage (default: 10)

Returns

• Vector{Dict}: Curriculum schedule with parameters for each stage

Example

curriculum = create_curriculum_schedule(stages=5, datasets_per_stage=15)

create_database(; database, address, username, password)

Creates a new Neo4j database instance if it doesn't already exist.

Arguments

• database::String: Name of the database to create
• address::String: Neo4j server address (e.g. "neo4j://localhost:7687")
• username::String: Neo4j authentication username (defaults to "neo4j")
• password::String: Neo4j authentication password

Notes

• Requires system database privileges to execute
• Silently returns if database already exists
• Temporarily switches to system database to perform creation

create_dqn_policy(; state_dim=11, action_dim=3, hidden_dims=[128, 64], learning_rate=0.001, epsilon=0.1)

Create a Deep Q-Network policy for assembly decisions.

Arguments

• state_dim::Int: Dimension of state representation (default: 11)
• action_dim::Int: Number of discrete actions (default: 3 for continue/next/terminate)
• hidden_dims::Vector{Int}: Hidden layer sizes (default: [128, 64])
• learning_rate::Float64: Learning rate (default: 0.001)
• epsilon::Float64: Exploration probability (default: 0.1)

Returns

• DQNPolicy: Initialized policy network

Example

policy = create_dqn_policy(hidden_dims=[256, 128, 64])

create_hmm_from_graph(graph::MetaGraph, config::ViterbiConfig) -> (states, transitions, emissions)

Create Hidden Markov Model parameters from a k-mer graph structure.

Create k-fold partitions with validation holdout for cross-validation.

create_node_constraints(
    graph;
    address,
    username,
    password,
    database
)

Creates unique identifier constraints for each node type in a Neo4j database.

Arguments

• graph: A MetaGraph containing nodes with TYPE properties
• address: Neo4j server address
• username: Neo4j username (default: "neo4j")
• password: Neo4j password
• database: Neo4j database name (default: "neo4j")

Details

Extracts unique node types from the graph and creates Neo4j constraints ensuring each node of a given type has a unique identifier property.

Failed constraint creation attempts are silently skipped.

create_tarchive(; directory, tarchive)

Creates a gzipped tar archive of the specified directory along with verification files.

Arguments

• directory: Source directory path to archive
• tarchive: Optional output archive path (defaults to directory name with .tar.gz extension)

Generated Files

• {tarchive}: The compressed tar archive
• {tarchive}.log: Contents listing of the archive
• {tarchive}.hashdeep.dfxml: Cryptographic hashes (MD5, SHA1, SHA256) of the archive

Returns

• Path to the created tar archive file

Generate cross-validation summary report.

current_unix_datetime() -> Int64

Get the current time as a Unix timestamp (seconds since epoch).

Returns

• Int: Current time as an integer Unix timestamp (seconds since January 1, 1970 UTC)

Examples

unix_time = current_unix_datetime()
# => 1709071368 (example value, will differ based on current time)

cypher(
    cmd;
    address,
    username,
    password,
    format,
    database
) -> Cmd

Constructs a command to execute Neo4j Cypher queries via cypher-shell.

Arguments

• cmd: The Cypher query command to execute
• address::String="neo4j://localhost:7687": Neo4j server address
• username::String="neo4j": Neo4j authentication username
• password::String="password": Neo4j authentication password
• format::String="auto": Output format (auto, verbose, or plain)
• database::String="neo4j": Target Neo4j database name

Returns

• Cmd: A command object ready for execution

dataframe_convert_dicts_to_json(df) -> Any

dataframe_replace_nothing_with_missing(df::DataFrames.DataFrame) -> DataFrames.DataFrame

Return the DataFrame with all nothing values replaced by missing.

dataframe_to_ndjson(df::DataFrame; outfile::Union{String,Nothing}=nothing)

Converts a DataFrame df into a newline-delimited JSON (NDJSON) string. Each line in the returned string represents one DataFrame row in JSON format, suitable for upload to Google BigQuery.

Keyword Arguments

• outfile::Union{String,Nothing}: If provided, writes the resulting NDJSON to the file path given.

Examples

```julia using DataFrames, Dates

Sample DataFrame

df = DataFrame( id = [1, 2, 3], name = ["Alice", "Bob", "Carol"], created = [DateTime(2025, 4, 8, 14, 30), DateTime(2025, 4, 8, 15, 0), missing] )

ndjsonstr = dataframetondjson(df) println(ndjsonstr)

Optionally, write to a file

dataframetondjson(df; outfile="output.ndjson")

deduplicate_fasta_file(in_fasta, out_fasta) -> Any

Remove duplicate sequences from a FASTA file while preserving headers.

Arguments

• in_fasta: Path to input FASTA file
• out_fasta: Path where deduplicated FASTA will be written

Returns

Path to the output FASTA file (same as out_fasta parameter)

Details

• Sequences are considered identical if they match exactly (case-sensitive)
• For duplicate sequences, keeps the first header encountered
• Input sequences are sorted by identifier before deduplication
• Preserves the original sequence formatting

detect_alphabet(seq::AbstractString) -> Symbol

Determines the alphabet of a sequence. The function scans through seq only once:

• If a 'T' or 't' is found (and no 'U/u'), the sequence is classified as DNA.
• If a 'U' or 'u' is found (and no 'T/t'), it is classified as RNA.
• If both T and U occur, an error is thrown.
• If a character outside the canonical nucleotide and ambiguity codes is encountered, the sequence is assumed to be protein.
• If neither T nor U are found, the sequence is assumed to be DNA.

detect_alphabet(sequence::BioSequences.LongAA) -> Symbol

Detect the alphabet of a LongAA sequence.

Always returns :AA.

detect_alphabet(sequence::BioSequences.LongDNA) -> Symbol

Detect the alphabet of a LongDNA sequence.

Always returns :DNA.

detect_alphabet(sequence::BioSequences.LongRNA) -> Symbol

Detect the alphabet of a LongRNA sequence.

Always returns :RNA.

detect_and_extract_sequence(
    record::Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}
) -> Tuple{Symbol, Union{BioSequences.LongAA, BioSequences.LongDNA{4}, BioSequences.LongRNA{4}}}

Detect alphabet and extract typed sequence from FASTX record in one step.

Convenience function that combines alphabet detection with type-safe sequence extraction, ideal for workflows that need to determine sequence type once at the beginning.

Arguments

• record::Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}: Input sequence record

Returns

• Tuple{Symbol, BioSequences.BioSequence}: (alphabetsymbol, typedsequence)

Examples

record = FASTX.FASTQ.Record("read1", "ATCG", "IIII")
alphabet, sequence = detect_and_extract_sequence(record)
# alphabet = :DNA, sequence = LongDNA{4} object

detect_bubbles_next(graph::MetaGraphsNext.MetaGraph, min_bubble_length::Int=2, max_bubble_length::Int=100) -> Vector{BubbleStructure}

Detect bubble structures (alternative paths) in the assembly graph.

detect_sequence_extension(
    record::Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}
) -> String

Detect sequence type from input and suggest appropriate file extension.

Arguments

• record: A FASTA/FASTQ record
• sequence: A string or BioSequence containing sequence data

Returns

• String: Suggested file extension:
  • ".fna" for DNA
  • ".frn" for RNA
  • ".faa" for protein
  • ".fa" for unrecognized sequences

determine_fasta_coverage_from_bam(bam) -> Any

Calculate per-base genomic coverage from a BAM file using bedtools.

Arguments

• bam::String: Path to input BAM file

Returns

• String: Path to the generated coverage file (.coverage.txt)

Details

Uses bedtools genomecov to compute per-base coverage. Creates a coverage file with the format: <chromosome> <position> <coverage_depth>. If the coverage file already exists, returns the existing file path.

Dependencies

Requires bedtools (automatically installed in conda environment)

determine_max_canonical_kmers(k, ALPHABET) -> Any

Calculate the maximum number of possible canonical k-mers for a given alphabet.

Arguments

• k::Integer: Length of k-mer
• ALPHABET::Vector{Char}: Character set (nucleotides or amino acids)

Returns

• Int: Maximum number of possible canonical k-mers

Details

• For amino acids (AA_ALPHABET): returns total possible k-mers
• For nucleotides: returns half of total possible k-mers (canonical form)
• Requires odd k-mer length for nucleotide alphabets

determine_max_possible_kmers(k, ALPHABET) -> Any

Calculate the total number of possible unique k-mers that can be generated from a given alphabet.

Arguments

• k: Length of k-mers to consider
• ALPHABET: Vector containing the allowed characters/symbols

Returns

• Integer representing the maximum number of possible unique k-mers (|Σ|ᵏ)

determine_primary_contig(qualimap_results) -> Any

Determines the contig with the greatest number of total bases mapping to it

Identify the primary contig based on mapping coverage from Qualimap results.

Arguments

• qualimap_results::DataFrame: DataFrame containing Qualimap alignment statistics with columns "Contig" and "Mapped bases"

Returns

• String: Name of the contig with the highest number of mapped bases

Description

Takes Qualimap alignment results and determines which contig has the most total bases mapped to it, which often indicates the main chromosomal assembly.

determine_read_lengths(
    fastq_file;
    total_reads
) -> Vector{Int64}

Calculate sequence lengths for reads in a FASTQ file.

Arguments

• fastq_file::String: Path to input FASTQ file
• total_reads::Integer=Inf: Number of reads to process (defaults to all reads)

Returns

• Vector{Int}: Array containing the length of each sequence read

dictvec_to_dataframe(dictvec::Vector{<:AbstractDict}; symbol_columns::Bool = true)

Convert a vector of dictionaries (with possibly non-uniform keys and any key type) into a DataFrame. Missing keys in a row are filled with missing.

Arguments

• dictvec: Vector of dictionaries.
• symbol_columns: If true (default), columns are named as Symbols (when possible), else as raw keys.

Returns

• DataFrames.DataFrame with columns as the union of all keys.

distance_matrix_to_newick(
;
    distance_matrix,
    labels,
    outfile
)

Create distance matrix from a column-major counts matrix (features as rows and entities as columns) where distance is a proportional to total feature count magnitude (size) and cosine similarity (relative frequency)

Convert a distance matrix into a Newick tree format using UPGMA hierarchical clustering.

Arguments

• distance_matrix: Square matrix of pairwise distances between entities
• labels: Vector of labels corresponding to the entities in the distance matrix
• outfile: Path where the Newick tree file will be written

Returns

• Path to the generated Newick tree file

Details

Performs hierarchical clustering using the UPGMA (average linkage) method and converts the resulting dendrogram into Newick tree format. The branch lengths in the tree represent the heights from the clustering.

document_frequency(
    documents
) -> Dict{T, Int64} where T<:(SubString{_A} where _A)

Calculate the document frequency of tokens across a collection of documents.

Arguments

• documents: Collection of text documents where each document is a string

Returns

• Dictionary mapping each unique token to the number of documents it appears in

Description

Computes how many documents contain each unique token. Each document is tokenized by splitting on whitespace. Tokens are counted only once per document, regardless of how many times they appear within that document.

download_bandage() -> String
download_bandage(outdir) -> Any

Downloads and installs Bandage, a bioinformatics visualization tool for genome assembly graphs.

Arguments

• outdir="/usr/local/bin": Target installation directory for the Bandage executable

Returns

• Path to the installed Bandage executable

Details

• Downloads Bandage v0.8.1 for Ubuntu
• Installs required system dependencies (libxcb-glx0, libx11-xcb-dev, libfontconfig, libgl1-mesa-glx)
• Attempts installation with sudo, falls back to root if sudo fails
• Skips download if Bandage is already installed at target location

Dependencies

Requires system commands: wget, unzip, apt

download_blast_db(; db, dbdir, source, wait)

Smart downloading of blast dbs depending on interactive, non interactive context

For a list of all available databases, run: Mycelia.list_blastdbs()

Downloads and sets up BLAST databases from various sources.

Arguments

• db: Name of the BLAST database to download
• dbdir: Directory to store the downloaded database (default: "~/workspace/blastdb")
• source: Download source - one of ["", "aws", "gcp", "ncbi"]. Empty string auto-detects fastest source
• wait: Whether to wait for download completion (default: true)

Returns

• String path to the downloaded database directory

download_genome_by_accession(
;
    accession,
    outdir,
    compressed
)

Downloads a genomic sequence from NCBI's nucleotide database by its accession number.

Arguments

• accession::String: NCBI nucleotide accession number (e.g. "NC_045512")
• outdir::String: Output directory path. Defaults to current directory
• compressed::Bool: If true, compresses output file with gzip. Defaults to true

Returns

• String: Path to the downloaded file (.fna or .fna.gz)

download_genome_by_ftp(; ftp, outdir)

Downloads a genome file from NCBI FTP server to the specified directory.

Arguments

• ftp::String: NCBI FTP path for the genome (e.g. "ftp://ftp.ncbi.nlm.nih.gov/.../")
• outdir::String: Output directory path. Defaults to current working directory.

Returns

• String: Path to the downloaded file

Notes

• If the target file already exists, returns the existing file path without re-downloading
• Downloads the genomic.fna.gz version of the genome

download_mmseqs_db(; db, dbdir, force, wait)

Downloads and sets up MMseqs2 reference databases for sequence searching and analysis.

Arguments

• db::String: Name of database to download (see table below)
• dbdir::String: Directory to store the downloaded database (default: "~/workspace/mmseqs")
• force::Bool: If true, force re-download even if database exists (default: false)
• wait::Bool: If true, wait for download to complete (default: true)

Returns

• Path to the downloaded database as a String

Available Databases

  Name                  Type            Taxonomy        Url                                                           
- UniRef100             Aminoacid            yes        https://www.uniprot.org/help/uniref
- UniRef90              Aminoacid            yes        https://www.uniprot.org/help/uniref
- UniRef50              Aminoacid            yes        https://www.uniprot.org/help/uniref
- UniProtKB             Aminoacid            yes        https://www.uniprot.org/help/uniprotkb
- UniProtKB/TrEMBL      Aminoacid            yes        https://www.uniprot.org/help/uniprotkb
- UniProtKB/Swiss-Prot  Aminoacid            yes        https://uniprot.org
- NR                    Aminoacid            yes        https://ftp.ncbi.nlm.nih.gov/blast/db/FASTA
- NT                    Nucleotide             -        https://ftp.ncbi.nlm.nih.gov/blast/db/FASTA
- GTDB                  Aminoacid            yes        https://gtdb.ecogenomic.org
- PDB                   Aminoacid              -        https://www.rcsb.org
- PDB70                 Profile                -        https://github.com/soedinglab/hh-suite
- Pfam-A.full           Profile                -        https://pfam.xfam.org
- Pfam-A.seed           Profile                -        https://pfam.xfam.org
- Pfam-B                Profile                -        https://xfam.wordpress.com/2020/06/30/a-new-pfam-b-is-released
- CDD                   Profile                -        https://www.ncbi.nlm.nih.gov/Structure/cdd/cdd.shtml
- eggNOG                Profile                -        http://eggnog5.embl.de
- VOGDB                 Profile                -        https://vogdb.org
- dbCAN2                Profile                -        http://bcb.unl.edu/dbCAN2
- SILVA                 Nucleotide           yes        https://www.arb-silva.de
- Resfinder             Nucleotide             -        https://cge.cbs.dtu.dk/services/ResFinder
- Kalamari              Nucleotide           yes        https://github.com/lskatz/Kalamari

download_sequence_by_accession(
;
    accession,
    outdir,
    database,
    format,
    suffix,
    compressed
)

Download a sequence from NCBI by accession number with flexible format options.

Arguments

• accession::String: NCBI accession number
• outdir::String: Output directory
• database::String: NCBI database ("nuccore", "protein", etc.)
• format::String: Output format ("fasta", "fastacdsna", "fastacdsaa", etc.)
• suffix::String: File suffix to append to accession for filename
• compressed::Bool: Whether to gzip compress the output (default: true)

Returns

• String: Path to the downloaded file

Downloads sequencing reads from NCBI's Sequence Read Archive (SRA).

Downloads reads using fasterq-dump. The function automatically detects whether the data is single-end or paired-end and returns appropriate file paths. Users should apply quality control based on their knowledge of the data type.

Arguments

• srr_identifier: SRA run identifier (e.g., "SRR1234567")
• outdir: Output directory for downloaded files (default: current directory)

Returns

Named tuple with:

• srr_id: The SRA identifier
• outdir: Output directory path
• files: Vector of downloaded file paths (1 file for single-end, 2 for paired-end)
• is_paired: Boolean indicating if data is paired-end

Example

# Download SRA data
result = Mycelia.download_sra_data("SRR1234567", outdir="./data")

# Apply appropriate QC based on data type
if result.is_paired
    # Paired-end data - use paired-end QC
    Mycelia.trim_galore_paired(forward_reads=result.files[1], reverse_reads=result.files[2])
else
    # Single-end data - use single-end QC
    Mycelia.qc_filter_short_reads_fastp(input=result.files[1])
end

download_viroid_reference_data(
    viroid_name::String;
    outdir,
    download_genome,
    download_cds,
    download_protein,
    max_sequences
) -> @NamedTuple{genome_files::Vector{String}, cds_files::Vector{String}, protein_files::Vector{String}, output_directory::String}

Download comprehensive viroid reference data including genome, CDS, and protein sequences.

Arguments

• viroid_name::String: Name or search term for the viroid (e.g., "Potato spindle tuber viroid")
• outdir::String: Output directory for downloaded files (default: current directory)
• download_genome::Bool: Download genomic DNA sequence (default: true)
• download_cds::Bool: Download CDS transcript sequences (default: true)
• download_protein::Bool: Download protein (FAA) sequences (default: true)
• max_sequences::Int: Maximum number of sequences per type (default: 10)

Returns

• NamedTuple: Paths to downloaded files (genomefiles, cdsfiles, protein_files)

Examples

# Download all data for Potato spindle tuber viroid
files = download_viroid_reference_data("Potato spindle tuber viroid", "viroid_data/")

# Download only genome sequences for general viroid search
files = download_viroid_reference_data("viroid", "viroid_genomes/"; 
                                     download_cds=false, download_protein=false)

draw_dendrogram_tree(
    mg::MetaGraphs.MetaDiGraph;
    width,
    height,
    fontsize,
    margins,
    mergenodesize,
    lineweight,
    filename
) -> Luxor.Drawing

Draw a dendrogram visualization of hierarchical clustering results stored in a MetaDiGraph.

Arguments

• mg::MetaGraphs.MetaDiGraph: Graph containing hierarchical clustering results. Must have :hcl in graph properties with clustering data and vertex properties containing :x, :y coordinates.

Keywords

• width::Integer=500: Width of output image in pixels
• height::Integer=500: Height of output image in pixels
• fontsize::Integer=12: Font size for node labels in points
• margins::Float64: Margin size in pixels, defaults to min(width,height)/20
• mergenodesize::Float64=1: Size of circular nodes at merge points
• lineweight::Float64=1: Thickness of dendrogram lines
• filename::String: Output filename, defaults to timestamp with .dendrogram.png extension

Returns

Nothing, but saves dendrogram image to disk and displays preview.

draw_radial_tree(
    mg::MetaGraphs.MetaDiGraph;
    width,
    height,
    fontsize,
    margins,
    mergenodesize,
    lineweight,
    filename
) -> Luxor.Drawing

Draw a radial hierarchical clustering tree visualization and save it as an image file.

Arguments

• mg::MetaGraphs.MetaDiGraph: A meta directed graph containing hierarchical clustering data with required graph properties :hcl containing clustering information.

Keywords

• width::Int=500: Width of the output image in pixels
• height::Int=500: Height of the output image in pixels
• fontsize::Int=12: Font size for node labels
• margins::Float64: Margin size (automatically calculated as min(width,height)/20)
• mergenodesize::Float64=1: Size of the merge point nodes
• lineweight::Float64=1: Thickness of the connecting lines
• filename::String: Output filename (defaults to timestamp with ".radial.png" suffix)

Details

The function creates a radial visualization of hierarchical clustering results where:

• Leaf nodes are arranged in a circle with labels
• Internal nodes represent merge points
• Connections show the hierarchical structure through arcs and lines

The visualization is saved as a PNG file and automatically previewed.

Required Graph Properties

The input graph must have:

• mg.gprops[:hcl].labels: Vector of leaf node labels
• mg.gprops[:hcl].order: Vector of ordered leaf nodes
• mg.gprops[:hcl].merges: Matrix of merge operations
• mg.vprops[v][:x]: X coordinate for each vertex
• mg.vprops[v][:y]: Y coordinate for each vertex

drop_empty_columns!(
    df::DataFrames.AbstractDataFrame
) -> DataFrames.AbstractDataFrame

Identify all columns that have only missing or empty values, and remove those columns from the dataframe in-place.

Returns a modified version of the original dataframe.

See also: dropemptycolumns

drop_empty_columns(df::DataFrames.AbstractDataFrame) -> Any

Identify all columns that have only missing or empty values, and remove those columns from the dataframe.

Returns a modified copy of the dataframe.

See also: dropemptycolumns!

edge_path_to_sequence(kmer_graph, edge_path) -> Any

Converts a path of edges in a kmer graph into a DNA sequence by concatenating overlapping kmers.

Arguments

• kmer_graph: A directed graph where vertices represent kmers and edges represent overlaps
• edge_path: Vector of edges representing a path through the graph

Returns

A BioSequences.LongDNASeq containing the merged sequence represented by the path

Details

The function:

1. Takes the first kmer from the source vertex of first edge
2. For each edge, handles orientation (forward/reverse complement)
3. Verifies overlaps between consecutive kmers
4. Concatenates unique bases to build final sequence

edge_probability(stranded_kmer_graph, edge) -> Any

Calculate the probability of traversing a specific edge in a stranded k-mer graph.

The probability is computed as the ratio of this edge's coverage weight to the sum of all outgoing edge weights from the source vertex.

edge_probability(stranded_kmer_graph, edge) -> Any

Arguments

• stranded_kmer_graph: A directed graph where edges represent k-mer connections
• edge: The edge for which to calculate the probability

Returns

• Float64: Probability in range [0,1] representing likelihood of traversing this edge Returns 0.0 if sum of all outgoing edge weights is zero

Note

Probability is based on the :coverage property of edges, using their length as weights

edgemer_to_vertex_kmers(
    edgemer
) -> Tuple{Kmers.Kmer{BioSequences.DNAAlphabet{2}}, Kmers.Kmer{BioSequences.DNAAlphabet{2}}}

Convert an edgemer to two vertex kmers.

This function takes an edgemer (a sequence of DNA nucleotides) and converts it into two vertex kmers. A kmer is a substring of length k from a DNA sequence. The first kmer is created from the first n-1 elements of the edgemer, and the second kmer is created from the last n-1 elements of the edgemer.

Arguments

• edgemer::AbstractVector{T}: A vector of DNA nucleotides where T is a subtype of BioSequences.DNAAlphabet{2}.

Returns

• Tuple{Kmers.Kmer{BioSequences.DNAAlphabet{2}}, Kmers.Kmer{BioSequences.DNAAlphabet{2}}}: A tuple containing two kmers.

ensure_next_graph(graph) -> Any

Automatically convert legacy graphs to next-generation format if needed.

This is a convenience function that checks the graph type and converts if necessary.

Arguments

• graph: Graph in either legacy or next-generation format

Returns

• MetaGraphsNext.MetaGraph: Graph in next-generation format

equally_spaced_samples(vector, n) -> Any

Sample n equally spaced elements from vector.

Arguments

• vector: Input vector to sample from
• n: Number of samples to return (must be positive)

Returns

A vector containing n equally spaced elements from the input vector.

equivalent_fasta_sequences(fasta_1, fasta_2) -> Bool

Compare two FASTA files to determine if they contain the same set of sequences, regardless of sequence order.

Arguments

• fasta_1::String: Path to first FASTA file
• fasta_2::String: Path to second FASTA file

Returns

• Bool: true if both files contain exactly the same sequences, false otherwise

Details

Performs a set-based comparison of DNA sequences by hashing each sequence. Sequence order differences between files do not affect the result.

error_rate_to_q_value(error_rate) -> Any

Convert a sequencing error probability to a Phred quality score (Q-value).

The calculation uses the standard Phred formula: Q = -10 * log₁₀(error_rate)

Arguments

• error_rate::Float64: Probability of error (between 0 and 1)

Returns

• q_value::Float64: Phred quality score

Analyze k-mer coverage distribution to detect if errors are singletons. Returns true if low-coverage k-mers (likely errors) are well-separated from high-coverage ones.

Estimate copy number for repeat region.

estimate_dense_matrix_memory(nrows::Integer, ncols::Integer)
estimate_dense_matrix_memory(T::DataType, nrows::Integer, ncols::Integer)

Estimate the memory required (in bytes) for a dense matrix.

• If T is provided, estimate memory for a matrix with element type T.
• If T is not provided, defaults to Float64.

estimate_genome_size_from_kmers(
    sequence::Union{AbstractString, BioSequences.LongSequence},
    k::Integer
) -> Dict

Estimate genome size from k-mer analysis using total k-mer count.

This function estimates genome size using the basic relationship: genomesize ≈ totalkmers - k + 1, where total_kmers is the sum of all k-mer counts. This is a simple estimation method; more sophisticated approaches accounting for sequencing depth, repeats, and errors may be more accurate.

Arguments

• sequence::Union{BioSequences.LongSequence, AbstractString}: Input sequence or string
• k::Integer: K-mer size for analysis

Returns

• Dict{String, Any}: Dictionary containing:
  • "unique_kmers": Number of unique k-mers observed
  • "total_kmers": Total k-mer count (sum of all frequencies)
  • "estimatedgenomesize": Estimated genome size
  • "actual_size": Length of input sequence (if provided)

Examples

# Estimate genome size from a sequence
sequence = "ATCGATCGATCGATCG"
result = estimate_genome_size_from_kmers(sequence, 5)

estimate_genome_size_from_kmers(
    records::AbstractArray{T<:Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}, 1},
    k::Integer
)

Estimate genome size from FASTQ/FASTA records using k-mer analysis.

Overload for processing FASTQ or FASTA records directly.

Arguments

• records::AbstractVector{T} where T <: Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}`: Input records
• k::Integer: K-mer size for analysis

Returns

• Dict{String, Any}: Dictionary with k-mer statistics and genome size estimate

Estimate memory usage for a graph with given number of k-mers. Provides rough estimate for memory monitoring.

estimate_sparse_matrix_memory(nrows::Integer, ncols::Integer; nnz=nothing, density=nothing)
estimate_sparse_matrix_memory(T::DataType, nrows::Integer, ncols::Integer; nnz=nothing, density=nothing)

Estimate the memory required (in bytes) for a sparse matrix in CSC format.

• If T is provided, estimate memory for a matrix with element type T.
• If T is not provided, defaults to Float64.
• You must specify either nnz (number of non-zeros) or density (proportion of non-zeros, between 0 and 1).

estimate_transition_probabilities(graph::MetaGraph, sequences::Vector) -> Matrix{Float64}

Estimate transition probabilities from observed sequences in the graph.

evaluate_assembly_agent(policy::DQNPolicy, validation_data::Vector{String}; episodes=10)

Evaluate a trained assembly agent on validation data.

Arguments

• policy::DQNPolicy: Trained policy to evaluate
• validation_data::Vector{String}: Validation FASTQ files
• episodes::Int: Number of evaluation episodes (default: 10)

Returns

• Dict{String, Float64}: Evaluation metrics including mean reward, assembly quality, etc.

Example

# metrics = evaluate_assembly_agent(trained_policy, validation_files)
# println("Mean reward: $(metrics["mean_reward"])")

evaluate_classification(true_labels, pred_labels)

Runs confusionmatrix, precisionrecall_f1, and accuracy. Pretty-prints macro metrics and accuracy. Returns a named tuple with all results and plots.

execute_continue_k_action(env::AssemblyEnvironment, action::AssemblyAction)

Execute a "continue with current k" action by performing error correction iterations.

Arguments

• env::AssemblyEnvironment: Current environment
• action::AssemblyAction: Action specifying correction parameters

Returns

• RewardComponents: Reward breakdown for the action

execute_next_k_action(env::AssemblyEnvironment, action::AssemblyAction)

Execute a "move to next k" action by progressing to the next prime k-mer size.

Arguments

• env::AssemblyEnvironment: Current environment
• action::AssemblyAction: Action specifying progression parameters

Returns

• RewardComponents: Reward breakdown for the action

execute_terminate_action(env::AssemblyEnvironment, action::AssemblyAction)

Execute a "terminate assembly" action and assess final assembly quality.

Arguments

• env::AssemblyEnvironment: Current environment
• action::AssemblyAction: Termination action

Returns

• RewardComponents: Final reward breakdown including assembly quality assessment

export_blast_db(; path_to_db, fasta)

Export sequences from a BLAST database to a gzipped FASTA file.

Arguments

• path_to_db: Path to the BLAST database
• fasta: Output path for the gzipped FASTA file (default: path_to_db * ".fna.gz")

Details

Uses conda's BLAST environment to extract sequences using blastdbcmd. The output is automatically compressed using pigz. If the output file already exists, the function will skip extraction.

export_blast_db_taxonomy_table(; path_to_db, outfile)

Exports a taxonomy mapping table from a BLAST database in seqid2taxid format.

Arguments

• path_to_db::String: Path to the BLAST database
• outfile::String: Output file path (defaults to input path + ".seqid2taxid.txt.gz")

Returns

• String: Path to the created output file

Details

Creates a compressed tab-delimited file mapping sequence IDs to taxonomy IDs. Uses blastdbcmd without GI identifiers for better cross-referencing compatibility. If the output file already exists, returns the path without regenerating.

Dependencies

Requires BLAST+ tools installed via Bioconda.

extract_pacbiosample_information(
    xml
) -> DataFrames.DataFrame

Extract biosample and barcode information from a PacBio XML metadata file.

Arguments

• xml: Path to PacBio XML metadata file

Returns

DataFrame with two columns:

• BioSampleName: Name of the biological sample
• BarcodeName: Associated DNA barcode identifier

extract_typed_sequence(
    record::Union{FASTX.FASTA.Record, FASTX.FASTQ.Record},
    sequence_type::Type{<:BioSequences.BioSequence}
) -> Any

Extract sequence from FASTX record using dynamically determined type.

This function provides type-safe sequence extraction by using the appropriate BioSequence type, avoiding hardcoded sequence types and string conversions.

Arguments

• record::Union{FASTX.FASTA.Record, FASTX.FASTQ.Record}: Input sequence record
• sequence_type::Type{<:BioSequences.BioSequence}: Target BioSequence type

Returns

• BioSequences.BioSequence: Typed sequence from the record

Examples

record = FASTX.FASTQ.Record("read1", "ATCG", "IIII")
seq_type = alphabet_to_biosequence_type(:DNA)
sequence = extract_typed_sequence(record, seq_type)

fasta_and_gff_to_genbank(; fasta, gff, genbank)

Convert FASTA sequence and GFF annotation files to GenBank format using EMBOSS seqret.

Arguments

• fasta::String: Path to input FASTA file containing sequence data
• gff::String: Path to input GFF file containing genomic features
• genbank::String: Path for output GenBank file

Details

Requires EMBOSS toolkit (installed via Bioconda). The function will:

1. Create necessary output directories
2. Run seqret to combine sequence and features
3. Generate a GenBank format file at the specified location

fasta_genome_size(fasta_file) -> Any

Calculate the total size (in bases) of all sequences in a FASTA file.

Arguments

• fasta_file::AbstractString: Path to the FASTA file

Returns

• Int: Sum of lengths of all sequences in the FASTA file

fasta_list_to_dense_kmer_counts(
;
    fasta_list,
    k,
    alphabet,
    temp_dir_parent,
    count_element_type,
    result_file,
    force,
    cleanup_temp
)

Create a dense k-mer counts table for a set of FASTA files, with disk-backed temporary storage, custom element type, robust error handling, and optional output file caching.

fasta_list_to_sparse_kmer_counts(
;
    fasta_list,
    k,
    alphabet,
    temp_dir_parent,
    count_element_type,
    rarefaction_data_filename,
    rarefaction_plot_basename,
    show_rarefaction_plot,
    result_file,
    out_dir,
    force,
    rarefaction_plot_kwargs...
)

Create a sparse kmer counts table (SparseMatrixCSC) from a list of FASTA files using a 3-pass approach. Pass 1 (Parallel): Counts kmers per file and writes to temporary JLD2 files. Pass 2 (Serial): Aggregates unique kmers, max count, nnz per file, and rarefaction data from temp files. Generates and saves a k-mer rarefaction plot. Pass 3 (Parallel): Reads temporary counts again to construct the final sparse matrix.

Optionally, a results filename can be provided to save/load the output. If the file exists and force is false, the result is loaded and returned. If force is true or the file does not exist, results are computed and saved.

Output Directory Behavior

• All auxiliary output files (e.g., rarefaction data, plots) are written to a common output directory.
• By default, this is:
  • The value of out_dir if provided.
  • Otherwise, the directory containing result_file (if provided and has a directory component).
  • Otherwise, the current working directory (pwd()).
• If you provide an absolute path for an output file (e.g. rarefaction_data_filename), that path is used directly.
• If both out_dir and a relative filename are given, the file is written to out_dir.

Arguments

• fasta_list::AbstractVector{<:AbstractString}: A list of paths to FASTA files.
• k::Integer: The length of the kmer.
• alphabet::Symbol: The alphabet type (:AA, :DNA, :RNA).
• temp_dir_parent::AbstractString: Parent directory for creating the temporary working directory. Defaults to Base.tempdir().
• count_element_type::Union{Type{<:Unsigned}, Nothing}: Optional. Specifies the unsigned integer type for the counts. If nothing (default), the smallest UInt type capable of holding the maximum observed count is used.
• rarefaction_data_filename::AbstractString: Filename for the TSV output of rarefaction data. If a relative path, will be written to out_dir.
• rarefaction_plot_basename::AbstractString: Basename for the output rarefaction plots. If a relative path, will be written to out_dir.
• show_rarefaction_plot::Bool: Whether to display the rarefaction plot after generation. Defaults to true.
• result_file::Union{Nothing, AbstractString}: Optional. If provided, path to a file to save/load the full results (kmers, counts, etc) as a JLD2 file.
• out_dir::Union{Nothing, AbstractString}: Optional. Output directory for auxiliary outputs. Defaults as described above.
• force::Bool: If true, recompute and overwrite the output file even if it exists. Defaults to false.
• rarefaction_plot_kwargs...: Keyword arguments to pass to plot_kmer_rarefaction for plot customization.

Returns

• NamedTuple{(:kmers, :counts, :rarefaction_data_path)}:
  • kmers: A sorted Vector of unique kmer objects.
  • counts: A SparseArrays.SparseMatrixCSC{V, Int} storing kmer counts.
  • rarefaction_data_path: Path to the saved TSV file with rarefaction data.

Raises

• ErrorException: If input fasta_list is empty, alphabet is invalid, or required Kmer/counting functions are not found.

fasta_table_to_fasta(fasta_df) -> Any

Convert a DataFrame containing FASTA sequence information into a vector of FASTA records.

Arguments

• fasta_df::DataFrame: DataFrame with columns "identifier", "description", and "sequence"

Returns

• Vector{FASTX.FASTA.Record}: Vector of FASTA records

fasta_to_reference_kmer_counts(; kmer_type, fasta)

Counts k-mer occurrences in a FASTA file, considering both forward and reverse complement sequences.

Arguments

• kmer_type: Type specification for k-mers (e.g., DNAKmer{21})
• fasta: Path to FASTA file containing reference sequences

Returns

• Dict{kmer_type, Int}: Dictionary mapping each k-mer to its total count across all sequences

fasta_to_table(fasta) -> DataFrames.DataFrame

Convert a FASTA file/record iterator to a DataFrame.

Arguments

• fasta: FASTA record iterator from FASTX.jl

Returns

• DataFrame with columns:
  • identifier: Sequence identifiers
  • description: Full sequence descriptions
  • sequence: Biological sequences as strings

fasta_xam_mapping_stats(; fasta, xam)

Calculate mapping statistics by comparing sequence alignments (BAM/SAM) to a reference FASTA.

Arguments

• fasta::String: Path to reference FASTA file
• xam::String: Path to alignment file (BAM or SAM format)

Returns

DataFrame with columns:

• contig: Reference sequence name
• contig_length: Length of reference sequence
• total_aligned_bases: Total number of bases aligned to reference
• mean_depth: Average depth of coverage (totalalignedbases/contig_length)

fastani_list(
;
    query_list,
    reference_list,
    outfile,
    threads,
    force
) -> Union{Nothing, Base.Process}

Run fastani with a query and reference list

Calculate Average Nucleotide Identity (ANI) between genome sequences using FastANI.

Arguments

• query_list::String: Path to file containing list of query genome paths (one per line)
• reference_list::String: Path to file containing list of reference genome paths (one per line)
• outfile::String: Path to output file that will contain ANI results
• threads::Int=Sys.CPU_THREADS: Number of parallel threads to use
• force::Bool=false: If true, rerun analysis even if output file exists

Output

Generates a tab-delimited file with columns:

• Query genome
• Reference genome
• ANI value (%)
• Count of bidirectional fragment mappings
• Total query fragments

Notes

• Requires FastANI to be available via Bioconda
• Automatically sets up required conda environment

fasterq_dump(
;
    outdir,
    srr_identifier
) -> NamedTuple{(:forward_reads, :reverse_reads, :unpaired_reads), <:Tuple{Union{Missing, String}, Union{Missing, String}, Union{Missing, String}}}

Download and compress sequencing reads from the SRA database using fasterq-dump.

Arguments

• outdir::String="": Output directory for the FASTQ files. Defaults to current directory.
• srr_identifier::String="": SRA run accession number (e.g., "SRR12345678")

Returns

Named tuple containing paths to the generated files:

• forward_reads: Path to forward reads file (*_1.fastq.gz) or missing
• reverse_reads: Path to reverse reads file (*_2.fastq.gz) or missing
• unpaired_reads: Path to unpaired reads file (*.fastq.gz) or missing

Outputs

Creates compressed FASTQ files in the output directory:

• {srr_identifier}_1.fastq.gz: Forward reads (for paired-end data)
• {srr_identifier}_2.fastq.gz: Reverse reads (for paired-end data)
• {srr_identifier}.fastq.gz: Unpaired reads (for single-end data)

Dependencies

Requires:

• fasterq-dump from the SRA Toolkit (installed via Conda)
• gzip for compression

Notes

• Skips download if output files already exist
• Uses up to 4 threads or system maximum, whichever is lower
• Allocates 1GB memory for processing
• Skips technical reads
• Handles both paired-end and single-end data automatically

Parallel FASTQ dump for multiple SRA files.

Converts multiple SRA files to FASTQ format in parallel. More efficient than sequential processing for large batches.

Arguments

• srr_identifiers: Vector of SRA run identifiers
• outdir: Output directory for FASTQ files (default: current directory)
• max_parallel: Maximum number of parallel conversions (default: 2)

Returns

Vector of named tuples with conversion results

Example

runs = ["SRR1234567", "SRR1234568"]
results = Mycelia.fasterq_dump_parallel(runs, outdir="./fastq_data")

fastq_record(; identifier, sequence, quality_scores)

Construct a FASTX FASTQ record from its components.

Arguments

• identifier::String: The sequence identifier without the '@' prefix
• sequence::String: The nucleotide sequence
• quality_scores::Vector{Int}: Quality scores (0-93) as raw integers

Returns

• FASTX.FASTQRecord: A parsed FASTQ record

Notes

• Quality scores are automatically capped at 93 to ensure FASTQ compatibility
• Quality scores are converted to ASCII by adding 33 (Phred+33 encoding)
• The record is constructed in standard FASTQ format with four lines:
  1. Header line (@ + identifier)
  2. Sequence
  3. Plus line
  4. Quality scores (ASCII encoded)

fastx2normalized_table(fastx::AbstractString) -> DataFrames.DataFrame

fastx2normalized_table(fastx) -> DataFrames.DataFrame

Read a FASTA or FASTQ file and convert its records into a normalized DataFrames.DataFrame where each row represents a sequence record and columns provide standardized metadata and sequence statistics.

Arguments

• fastx::AbstractString: Path to a FASTA or FASTQ file. The file must exist and be non-empty. The file type is inferred from the filename using Mycelia.FASTA_REGEX and Mycelia.FASTQ_REGEX.

Returns

• DataFrames.DataFrame: A data frame where each row contains information for a record from the input file, and columns include:
  • fastx_path: Basename of the input file.
  • fastx_sha256: Aggregated SHA256 hash of all record SHA256s in the file.
  • record_identifier: Identifier from the record header.
  • record_description: Description from the record header.
  • record_sha256: SHA256 hash of the record sequence.
  • record_quality: Vector of quality scores (Vector{Float64}) for FASTQ, or missing for FASTA.
  • record_alphabet: Sorted, joined string of unique, uppercase characters in the record sequence.
  • record_type: Alphabet type detected by Mycelia.detect_alphabet (e.g., :DNA, :RNA, etc.).
  • mean_record_quality: Mean quality score (for FASTQ), or missing (for FASTA).
  • median_record_quality: Median quality score (for FASTQ), or missing (for FASTA).
  • record_length: Length of the sequence.
  • record_sequence: The sequence string itself.

Notes

• The function asserts that the file exists and is not empty.
• File type is determined by regex matching on the filename.
• For FASTA files, quality-related columns are set to missing.
• For FASTQ files, quality scores are extracted and statistics are computed.
• Record SHA256 hashes are aggregated to compute a file-level SHA256 via Mycelia.metasha256.
• Requires the following namespaces: DataFrames, Statistics, Mycelia, FASTX, and Base.basename.
• The function returns the columns in the order: fastx_path, fastx_sha256, followed by all other record columns.

Example

import DataFrames
import Mycelia
import FASTX

table = fastx2normalized_table("example.fasta")
DataFrames.first(table, 3)

fastx_stats(fastx) -> DataFrames.DataFrame

Calculate basic statistics for FASTQ/FASTA sequence files using seqkit.

Arguments

• fastq::String: Path to input FASTQ/FASTA file

Details

Automatically installs and uses seqkit from Bioconda to compute sequence statistics including number of sequences, total bases, GC content, average length, etc.

Dependencies

• Requires Conda and Bioconda channel
• Installs seqkit package if not present

Returns

Returns a DataFrame of the table

https://bioinf.shenwei.me/seqkit/usage/#stats

fastx_to_contig_lengths(
    fastx
) -> OrderedCollections.OrderedDict

Generate detailed mapping statistics for each reference sequence/contig in a XAM (SAM/BAM/CRAM) file.

Arguments

• xam: Path to XAM file or XAM object

Returns

A DataFrame with per-contig statistics including:

• n_aligned_reads: Number of aligned reads
• total_aligned_bases: Sum of alignment lengths
• total_alignment_score: Sum of alignment scores
• Mapping quality statistics (mean, std, median)
• Alignment length statistics (mean, std, median)
• Alignment score statistics (mean, std, median)
• Percent mismatches statistics (mean, std, median)

Note: Only primary alignments (isprimary=true) and mapped reads (ismapped=true) are considered.

fastx_to_kmer_graph(
    KMER_TYPE,
    fastx::AbstractString
) -> MetaGraphs.MetaGraph

Constructs a k-mer graph from a single FASTX format string.

Arguments

• KMER_TYPE: The k-mer type specification (e.g., DNAKmer{K} where K is k-mer length)
• fastx::AbstractString: Input sequence in FASTX format (FASTA or FASTQ)

Returns

• KmerGraph: A directed graph where vertices are k-mers and edges represent overlaps

fastx_to_kmer_graph(
    KMER_TYPE,
    fastxs::AbstractVector{<:AbstractString}
) -> MetaGraphs.MetaGraph

Create an in-memory kmer-graph that records:

• all kmers
• counts
• all observed edges between kmers
• edge orientations
• edge counts

Construct a kmer-graph from one or more FASTX files (FASTA/FASTQ).

Arguments

• KMER_TYPE: Type for kmer representation (e.g., DNAKmer{K})
• fastxs: Vector of paths to FASTX files

Returns

A MetaGraph where:

• Vertices represent unique kmers with properties:
  • :kmer => The kmer sequence
  • :count => Number of occurrences
• Edges represent observed kmer adjacencies with properties:
  • :orientation => Relative orientation of connected kmers
  • :count => Number of observed transitions

fibonacci_numbers_less_than(
    n::Int64
) -> Union{Vector{Any}, Vector{Int64}}

Generate a sequence of Fibonacci numbers strictly less than the input value.

Arguments

• n::Int: Upper bound (exclusive) for the Fibonacci sequence

Returns

• Vector{Int}: Array containing Fibonacci numbers less than n

filesize_human_readable(f) -> Any

Gets the size of a file and returns it in a human-readable format.

Arguments

• f: The path to the file, either as a String or an AbstractString.

Returns

A string representing the file size in a human-readable format (e.g., "3.40 MB").

Details

This function internally uses filesize(f) to get the file size in bytes, then leverages Base.format_bytes to convert it into a human-readable format with appropriate units (KB, MB, GB, etc.).

Examples

julia> filesize_human_readable("my_image.jpg")
"2.15 MB"

See Also

• filesize: Gets the size of a file in bytes.
• Base.format_bytes: Converts a byte count into a human-readable string.

Finalize assembly by combining information from all k-mer sizes. Phase 5.1b: Enhanced with accuracy metrics and reward tracking.

Finalize cross-validation results with comprehensive summary.

Finalize iterative assembly by combining results from all k-mer sizes and iterations.

Find potential bubble paths from an entry vertex.

Find all connected components in the graph. Returns a vector of vectors of vertex indices.

Specialized method for directed graphs using weakly connected components.

Specialized method for undirected graphs.

find_contigs_next(graph::MetaGraphsNext.MetaGraph, min_contig_length::Int=500) -> Vector{ContigPath}

Extract linear contigs from the assembly graph.

Find Eulerian path starting from a specific vertex using Hierholzer's algorithm.

find_eulerian_paths_next(graph::MetaGraphsNext.MetaGraph) -> Vector{Vector{String}}

Find Eulerian paths in the assembly graph. An Eulerian path visits every edge exactly once.

Find valid starting vertices for Eulerian paths.

find_fasta_files(input_path::String) -> Vector{String}

Find all FASTA files in a directory or return single file if path is a file.

Uses the existing FASTA_REGEX constant to identify FASTA files.

Arguments

• input_path: Path to directory or single FASTA file

Returns

• Vector of FASTA file paths

Example

fasta_files = find_fasta_files("./genomes/")

Find the optimal starting k-mer size using sparsity detection. Only considers prime k-mer sizes for optimal performance.

Find a limited-length path from a starting vertex.

Find a linear path through the graph.

find_matching_prefix(
    filename1::String,
    filename2::String;
    strip_trailing_delimiters
) -> String

Find the longest common prefix between two filenames.

Arguments

• filename1::String: First filename to compare
• filename2::String: Second filename to compare

Keywords

• strip_trailing_delimiters::Bool=true: If true, removes trailing dots, hyphens, and underscores from the result

Returns

• String: The longest common prefix found between the filenames

find_nonempty_columns(df) -> Any

Identify all columns that have only missing or empty values

Returns as a bit array

See also: dropemptycolumns, dropemptycolumns!

Find optimal sequence path through graph using maximum likelihood principles. Returns improved sequence and likelihood improvement score.

Find where two paths converge.

Find all primes in a range (convenience function).

find_quality_weighted_path(
    graph,
    start_vertex;
    max_path_length
) -> Vector

Find a quality-weighted path through a qualmer graph starting from a given vertex. Uses joint probability as the primary weighting factor for path selection.

Arguments

• graph: Qualmer graph (MetaGraphsNext with QualmerVertexData)
• start_vertex: Starting vertex for path traversal
• max_path_length::Int=1000: Maximum path length to prevent infinite loops

Returns

• Vector{Int}: Path as sequence of vertex indices

Details

At each step, selects the unvisited neighbor with the highest joint probability. Terminates when no unvisited neighbors are available or max length is reached.

Find vertices that could be part of repeat regions.

find_resampling_stretches(
;
    record_kmer_solidity,
    solid_branching_kmer_indices
)

Identifies sequence regions that require resampling based on kmer solidity patterns.

Arguments

• record_kmer_solidity::BitVector: Boolean array where true indicates solid kmers
• solid_branching_kmer_indices::Vector{Int}: Indices of solid branching kmers

Returns

• Vector{UnitRange{Int64}}: Array of ranges (start:stop) indicating stretches that need resampling

Details

Finds continuous stretches of non-solid kmers and extends them to the nearest solid branching kmers on either side. These stretches represent regions that need resampling.

If a stretch doesn't have solid branching kmers on both sides, it is excluded from the result. Duplicate ranges are removed from the final output.

find_true_ranges(
    bool_vec::AbstractVector{Bool};
    min_length
) -> Vector

Finds contiguous ranges of true values in a boolean vector.

Arguments

• bool_vec::AbstractVector{Bool}: Input boolean vector to analyze
• min_length=1: Minimum length requirement for a range to be included

Returns

Vector of tuples (start, end) where each tuple represents the indices of a contiguous range of true values meeting the minimum length requirement.

first_of_each_group(
    gdf::DataFrames.GroupedDataFrame{DataFrames.DataFrame}
) -> Any

Return a DataFrame containing the first row from each group in a GroupedDataFrame.

Arguments

• gdf::GroupedDataFrame: A grouped DataFrame created using groupby

Returns

• DataFrame: A new DataFrame containing first row from each group

frequency_matrix_to_bray_curtis_distance_matrix(counts_table)

Pairwise Bray-Curtis distance between columns of counts_table.

frequency_matrix_to_cosine_distance_matrix(probability_matrix)

Pairwise cosine distance between columns of probability_matrix.

frequency_matrix_to_euclidean_distance_matrix(counts_table)

Pairwise Euclidean distance between columns of counts_table.

frequency_matrix_to_jaccard_distance_matrix(matrix)

Thresholds the input matrix at >0 to obtain a binary matrix, then computes pairwise Jaccard distance between columns. Accepts any numeric matrix.

frequency_matrix_to_jensen_shannon_distance_matrix(probability_matrix)

Pairwise Jensen-Shannon divergence between columns of probability_matrix.

Arguments

• probability_matrix: Matrix where each column is a probability distribution (sums to 1.0).

Returns

• Symmetric matrix of Jensen-Shannon divergence values between columns.

gamma_pca_epca(M::AbstractMatrix{<:Real}; k::Int=0)

Perform Gamma EPCA on a matrix M (features × samples).

When to use

Use for positive continuous data, such as rates, times, or strictly positive measurements.

Keyword arguments

• k : desired number of latent dimensions; if k<1 defaults to min(n_samples-1, n_features, 10)

Returns

NamedTuple with fields

• model : the fitted ExpFamilyPCA.GammaEPCA object
• scores : k×n_samples matrix of sample scores
• loadings : k×n_features matrix of feature loadings

gaussian_pca_epca(M::AbstractMatrix{<:Real}; k::Int=0)

Perform Gaussian EPCA on a matrix M (features × samples).

When to use

Use for real-valued continuous data (centered, can be negative or positive), such as normalized or standardized measurements.

Keyword arguments

• k : desired number of latent dimensions; if k<1 defaults to min(n_samples-1, n_features, 10)

Returns

NamedTuple with fields

• model : the fitted ExpFamilyPCA.GaussianEPCA object
• scores : k×n_samples matrix of sample scores
• loadings : k×n_features matrix of feature loadings

genbank_to_codon_frequencies(
    genbank;
    allow_all
) -> Dict{BioSymbols.AminoAcid, Dict{Kmers.Kmer{BioSequences.DNAAlphabet{2}, 3, 1}, Int64}}

Analyze codon usage frequencies from genes in a GenBank file.

Arguments

• genbank: Path to GenBank format file containing genomic sequences and annotations
• allow_all: If true, initializes frequencies for all possible codons with count=1 (default: true)

Returns

Nested dictionary mapping amino acids to their corresponding codon usage counts:

• Outer key: AminoAcid (including stop codon)
• Inner key: DNACodon
• Value: Count of codon occurrences

Details

• Only processes genes marked as ':misc_feature' in the GenBank file
• Analyzes both forward and reverse complement sequences
• Determines coding strand based on presence of stop codons and start codons
• Skips ambiguous sequences that cannot be confidently oriented

genbank_to_fasta(; genbank, fasta, force)

Convert a GenBank format file to FASTA format using EMBOSS seqret.

Arguments

• genbank: Path to input GenBank format file
• fasta: Optional output FASTA file path (defaults to input path with .fna extension)
• force: If true, overwrites existing output file (defaults to false)

Returns

Path to the output FASTA file

Notes

• Requires EMBOSS suite (installed automatically via Conda)
• Will not regenerate output if it already exists unless force=true

generate_all_possible_canonical_kmers(k, alphabet) -> Any

Create distance matrix from a column-major counts matrix (features as rows and entities as columns) where distance is a proportional to total feature count magnitude (size) and cosine similarity (relative frequency)

Generate all possible canonical k-mers of length k from the given alphabet.

For DNA/RNA sequences, returns unique canonical k-mers where each k-mer is represented by the lexicographically smaller of itself and its reverse complement. For amino acid sequences, returns all possible k-mers without canonicalization.

Arguments

• k: Length of k-mers to generate
• alphabet: Vector of BioSymbols (DNA, RNA or AminoAcid)

Returns

• Vector of k-mers, canonicalized for DNA/RNA alphabets

generate_all_possible_kmers(k, alphabet) -> Any

Create distance matrix from a column-major counts matrix (features as rows and entities as columns) where distance is a proportional to total feature count magnitude (size) and cosine similarity (relative frequency)

Generate a sorted list of all possible k-mers for a given alphabet.

Arguments

• k::Integer: Length of k-mers to generate
• alphabet: Collection of symbols (DNA, RNA, or amino acids) from BioSymbols

Returns

• Sorted Vector of Kmers of the appropriate type (DNA, RNA, or amino acid)

Generate alternative qualmer paths through the graph using quality-aware probabilistic sampling.

generate_and_save_kmer_counts(; 
    bioalphabet, 
    fastas, 
    k, 
    output_dir=pwd(),
    filename=nothing
)

Generates and saves k-mer counts for a list of FASTA files for a single k.

Keyword Arguments

• bioalphabet: Alphabet type (e.g., :DNA).
• fastas: List of FASTA file paths.
• k: Value of k (e.g., 9).
• output_dir: (optional) Directory to write output files (default: current directory).
• filename: (optional) Full file name for output (default: "{Mycelia.normalized_current_date()}.{lowercase(string(bioalphabet))}{k}mers.jld2").

Output

Saves a .jld2 file with the specified file name in output_dir if it does not already exist.

Generate a binary (Bernoulli) matrix with given dimensions and probability.

Arguments

• n_features::Int: Number of features (rows)
• n_samples::Int: Number of samples (columns)
• p::Float64: Probability of 1 in the Bernoulli distribution

Returns

• Matrix{Bool}: Binary matrix with dimensions (nfeatures, nsamples)

Generate consensus pangenome from cross-validation results.

Generate sequence for a contig path.

Generate coverage profile for a contig path.

Generate alternative k-mers for improvement attempts using proper k-mer objects.

generate_paired_end_reads(reference_seq, coverage, read_length, insert_size; error_rate=0.01)

Generate realistic paired-end sequencing reads from a reference sequence.

Simulates paired-end Illumina sequencing with realistic insert sizes, read lengths, and optional sequencing errors for assembly benchmarking.

Arguments

• reference_seq: Reference sequence (BioSequences.LongDNA{4})
• coverage: Target sequencing coverage depth
• read_length: Length of each read in base pairs
• insert_size: Insert size between paired reads
• error_rate: Sequencing error rate (default: 0.01)

Returns

• Tuple of (forwardreads, reversereads) as vectors of BioSequences.LongDNA{4}

See Also

• simulate_illumina_paired_reads: For more sophisticated read simulation using ART
• introduce_sequencing_errors: For adding realistic sequencing errors

Generate a Poisson matrix with given dimensions and rate parameter.

Arguments

• n_features::Int: Number of features (rows)
• n_samples::Int: Number of samples (columns)
• λ::Float64: Rate parameter for the Poisson distribution

Returns

• Matrix{Int}: Poisson matrix with dimensions (nfeatures, nsamples)

generate_polished_sequence(states::Vector{ViterbiState}, observations::Vector{String}, 
                          config::ViterbiConfig) -> (String, Vector{Tuple{Int, String, String}})

Generate polished sequence and track corrections made.

Generate sequence of prime k-mer sizes starting from min_k.

Generate reasoning for assembly approach recommendation.

generate_taxa_abundances_plot(
    joint_reads_to_taxon_lineage_table::DataFrames.DataFrame;
    taxa_level::String,
    top_n::Int = 30,
    kwargs...
)

Convenience wrapper function to generate taxa abundance visualization with default parameters and save to a file if requested.

Arguments

• joint_reads_to_taxon_lineage_table: DataFrame with sample_id and taxonomic assignments
• taxa_level: Taxonomic level to analyze
• top_n: Number of top taxa to display individually
• kwargs...: Additional parameters to pass to plottaxaabundances

Returns

• fig: CairoMakie figure object
• ax: CairoMakie axis object
• taxa_colors: Dictionary mapping taxa to their assigned colors

generate_test_fastq_data(n_reads::Int, read_length::Int, filename::String)

Generate test FASTQ data for benchmarking purposes.

Creates a FASTQ file with random DNA sequences and realistic quality scores suitable for performance testing and validation.

Arguments

• n_reads::Int: Number of reads to generate
• read_length::Int: Length of each read in base pairs
• filename::String: Output filename for the FASTQ file

Details

• Generates random DNA sequences using BioSequences.randdnaseq
• Assigns realistic quality scores (Phred+33 encoding, range 20-40)
• Uses existing Mycelia I/O functions for consistency

See Also

• Mycelia.write_fastq: For writing FASTQ records
• Mycelia.fastq_record: For creating FASTQ records
• Mycelia.simulate_illumina_paired_reads: For more sophisticated read simulation

generate_test_genome_with_genes(genome_size, gene_density=0.02)

Generate a test genome with simulated gene positions for annotation benchmarking.

Creates a random DNA sequence with estimated gene positions based on gene density, suitable for testing gene prediction algorithms.

Arguments

• genome_size: Size of the genome in base pairs
• gene_density: Proportion of genome that consists of genes (default: 0.02)

Returns

• Tuple of (genomesequence, genepositions) where gene_positions is a vector of (start, end) tuples

See Also

• random_fasta_record: For generating random FASTA sequences
• save_genome_as_fasta: For saving genomes to FASTA format

generate_test_sequences(genome_size::Int, n_sequences::Int=1)

Generate test DNA sequences for k-mer analysis benchmarking.

Creates random DNA sequences suitable for k-mer counting and analysis performance testing.

Arguments

• genome_size::Int: Size of each generated sequence in base pairs
• n_sequences::Int: Number of sequences to generate (default: 1)

Returns

• Vector of BioSequences.LongDNA{4} sequences

See Also

• random_fasta_record: For generating FASTA records with random sequences
• BioSequences.randdnaseq: For generating random DNA sequences

generate_test_sequences(
    config::Mycelia.BenchmarkConfig
) -> Vector{FASTX.FASTA.Record}

Generate synthetic DNA sequences for benchmarking.

Arguments

• config: BenchmarkConfig with test parameters

Returns

• Vector of FASTA records for testing

generate_training_datasets(; n_datasets=20, genome_sizes=[10000, 50000, 100000], 
                          error_rates=[0.001, 0.01, 0.05], coverage_levels=[20, 30, 50])

Generate diverse training datasets for RL agent training.

This function creates simulated genomic datasets with varying characteristics to provide comprehensive training scenarios for the RL agent.

Arguments

• n_datasets::Int: Total number of datasets to generate (default: 20)
• genome_sizes::Vector{Int}: Range of genome sizes to simulate (default: [10K, 50K, 100K])
• error_rates::Vector{Float64}: Range of sequencing error rates (default: [0.1%, 1%, 5%])
• coverage_levels::Vector{Int}: Range of coverage depths (default: [20x, 30x, 50x])

Returns

• Vector{String}: Paths to generated training FASTQ files

Example

training_files = generate_training_datasets(n_datasets=50, genome_sizes=[50000, 100000, 200000])

generate_transterm_coordinates_from_fasta(fasta) -> Any

Generate minimal coordinate files required for TransTermHP analysis from FASTA sequences.

Creates artificial gene annotations at sequence boundaries to enable TransTermHP to run without real gene annotations. For each sequence in the FASTA file, generates two single-base-pair "genes" at positions 1-2 and (L-1)-L, where L is sequence length.

Arguments

• fasta: Path to input FASTA file containing sequences to analyze

Returns

• Path to generated coordinate file (original path with ".coords" extension)

Format

Generated coordinate file follows TransTermHP format: gene_id start stop chromosome

where chromosome matches FASTA sequence identifiers.

See also: run_transterm

generate_transterm_coordinates_from_gff(gff_file) -> Any

Convert a GFF file to a coordinates file compatible with TransTermHP format.

Arguments

• gff_file::String: Path to input GFF file

Processing

• Converts 1-based to 0-based coordinates
• Extracts gene IDs from the attributes field
• Retains columns: gene_id, start, end, seqid

Returns

• Path to the generated coordinates file (original filename with '.coords' suffix)

Output Format

Space-delimited file with columns: gene_id, start, end, seqid

get_base_extension(filename::String) -> String

Extract the base file extension from a filename, handling compressed files.

For regular files, returns the last extension. For gzipped files, returns the extension before .gz.

get_biosequence_alphabet(s::BioSequences.BioSequence) -> Any

Return the alphabet associated with a BioSequence type.

Arguments

• s::BioSequences.BioSequence: A subtype instance.

Returns

BioSymbols.Alphabet of the sequence type.

get_correct_quality(tech::Symbol, pos::Int, read_length::Int) -> Int

Simulates a Phred quality score (using the Sanger convention) for a correctly observed base. For Illumina, the quality score is modeled to decay linearly from ~40 at the start to ~20 at the end of the read. For other technologies, the score is sampled from a normal distribution with parameters typical for that platform.

Returns an integer quality score.

get_error_quality(tech::Symbol) -> Int

Simulates a Phred quality score (using the Sanger convention) for a base observed with an error. Error bases are assigned lower quality scores than correctly observed bases. For Illumina, scores typically range between 5 and 15; for nanopore and pacbio, slightly lower values are used; and for ultima, a modest quality score is assigned.

Returns an integer quality score.