Package 'WithinHostExperiment'

Description

Provides a RangedSummarizedExperiment-class-based container and workflow functions for analysing intra-host single nucleotide variants (iSNV) from pathogen deep sequencing data.

Details

The package fills an infrastructure gap between upstream variant callers (iVar, LoFreq, Freebayes) and downstream transmission bottleneck workflows. Its core contribution is replicate-aware, uncertainty-preserving quality control – variants are flagged rather than deleted, preserving full provenance for reproducible analysis. Standardised export functions provide interoperability with tools such as ViralBottleneck.

Key function groups

Import

readWithinHost (from VCF), readWithinHostTable (from TSV/CSV)

Quality control

ISNVFilter (configure thresholds), flagISNV (apply standard QC), flagReplicateDiscordance (replicate-aware QC), flagPrimerSites (amplicon primer masking), qcSummary (QC report), passedISNV (extract clean subset)

Diversity

calcDiversity, shannonISNV, piISNV, wattersonISNV, compareDiversity

Transmission analysis

transmissionPairs, exportPairFrequencies, asViralBottleneckInput, calcSharedVariants, quickBottleneck

Neutrality tests

tajimaD, fusFs, calcNeutralityTests, simpsonISNV, chao1ISNV

Annotation

annotateFromGFF (positional gene annotation from GFF3 files), annotateCodonChange (CDS-aware codon translation and amino acid change prediction from GFF3 + reference FASTA)

Selection analysis

dndsWithinHost, estimateSelectionCoefficient

Multi-caller consensus

consensusISNV, calcCallerConcordance, benchmarkCallers

Longitudinal analysis

trackFrequency, detectEmergingVariants, plotFrequencyTrajectory

Exact bottleneck estimation

exactBottleneck

Strand bias

sciStrandBias, flagStrandBiasSCI

Temporal QC

flagTemporalInconsistency (longitudinal consistency filtering – flag transient iSNVs that appear at only one timepoint)

Cross-contamination detection

detectCrossContamination (identify index- hopping signatures between co-sequenced samples)

Adaptive filtering

adaptiveFilter

Diversity modelling

fitDiversityModel, validateDiversity

Site frequency spectrum

WithinHostSFS-class (S4 class), buildSFS and buildSFSList (construct from WHE), correctSFSBias (ascertainment bias correction), neutralityFromSFS (Tajima's D + Fu's Fs + Fay & Wu's H from a single SFS), compareSFS (chi-squared SFS comparison), plotSFS (publication-quality SFS plot)

Visualisation

plotFrequencySpectrum, plotPairScatter, plotQCDashboard

Scope and limitations

Single-site iSNV resolution

This package operates at the level of individual single nucleotide variants (iSNVs). It does not perform haplotype phasing or quasispecies reconstruction, which require read-level (BAM) analysis that is outside the scope of RangedSummarizedExperiment-based containers. For haplotype-aware within-host analysis, consider CliqueSNV (Knyazev et al. 2021 Nucleic Acids Res 49:e102), ShoRAH (Zagordi et al. 2011 BMC Bioinformatics 12:119), or similar tools. Their output can be imported into WithinHostExperiment for downstream diversity and QC analysis.

Co-infection vs diversity

Distinguishing genuine within-host evolution from co-infection (superinfection) requires phylogenetic analysis of full haplotypes. The multi-caller consensus and replicate-concordance modules reduce artefacts but do not resolve this distinction.

Selection coefficient estimation

estimateSelectionCoefficient assumes a simple deterministic model. For stochastic (drift-aware) estimates, consider dedicated population-genetics frameworks.

Design philosophy

1. Flag, don't delete. QC functions mark variants in the qcPass assay rather than removing rows, preserving uncertainty and enabling re-analysis with different thresholds.

2. Replicate-aware. Technical replicate concordance is a first-class QC criterion, not a post-hoc check.

3. Interoperable, don't reinvent. The package provides standardised export functions compatible with existing downstream tools rather than reimplementing their methods.

4. Bioconductor-native. Built on RangedSummarizedExperiment, interoperable via VRanges coercion.

Author(s)

Maintainer: Cuiwei Gao [email protected]

References

Roder AE et al. (2023). Optimized quantification of intra-host viral diversity in SARS-CoV-2 and influenza virus sequence data. mBio 14:e01046-23. doi:10.1128/mbio.01046-23

Bendall EE et al. (2023). Rapid transmission and tight bottlenecks constrain the evolution of highly transmissible SARS-CoV-2 variants. Nat. Commun. 14:272. doi:10.1038/s41467-023-36001-5

Sobel Leonard A et al. (2017). Transmission Bottleneck Size Estimation from Pathogen Deep-Sequencing Data. J Virol 91:e00171-17. doi:10.1128/JVI.00171-17

Mostefai F et al. (2024). Refining SARS-CoV-2 intra-host variation by leveraging large-scale sequencing data. NAR Genomics and Bioinformatics 6:lqae145. doi:10.1093/nargab/lqae145

Zheng B, Johnson PCD, Hughes J (2025). ViralBottleneck: an R package for estimating viral transmission bottlenecks from deep sequencing data using multiple methods. Virus Evolution 11:veaf071. doi:10.1093/ve/veaf071

See Also

Useful links:

• https://github.com/CuiweiG/WithinHostExperiment

• Report bugs at https://github.com/CuiweiG/WithinHostExperiment/issues

Examples

# List exported functions
ls("package:WithinHostExperiment")

Description

Returns an ISNVFilter with default thresholds tuned for the error profile and read characteristics of the specified sequencing technology.

Usage

adaptiveFilter(technology = c("illumina", "ont", "pacbio_hifi", "generic"))

Arguments

Details

Each technology has distinct error modes that motivate different default thresholds:

Illumina (minFreq = 0.03, minDepth = 100)

Illumina short reads have low per-base error rates (~0.1–0.5%) but systematic context-dependent errors. A 3% frequency threshold safely exceeds the noise floor while retaining true low-frequency variants. Grubaugh et al. (2019) measured iSNVs above 3% accurately from at least 1,000 RNA copies sequenced to at least 400x, so 100x is a permissive floor; low-frequency false positives inflate diversity estimates markedly (McCrone & Lauring 2016).

Oxford Nanopore (ONT) (minFreq = 0.10, minDepth = 60, maxStrandBias = Inf)

ONT reads have higher per-base error rates (~1–5% with R10 chemistry), requiring a more conservative frequency threshold of 10% to avoid false positives; Bull et al. (2020) found that ONT sequencing does not accurately detect variants at low read-count frequencies, and reported over 99% sensitivity and precision for single-nucleotide variants above about 60-fold depth, which sets the depth threshold. Strand bias filtering is disabled (Inf) rather than given an untested threshold; supply maxStrandBias explicitly when the data support one.

PacBio HiFi (minFreq = 0.05, minDepth = 50)

HiFi consensus reads reach about 99.8% read accuracy via circular consensus sequencing (Wenger et al. 2019), allowing a frequency threshold below ONT but above Illumina. Depth of 50x reflects the lower throughput but high per-read accuracy.

Generic

Falls back to the standard ISNVFilter defaults (minFreq = 0.03, minDepth = 100), suitable for Illumina or when the technology is unknown.

Value

An ISNVFilter object.

References

McCrone JT, Lauring AS (2016). Measurements of intrahost viral diversity are extremely sensitive to systematic errors in variant calling. J Virol 90:6884-6895. doi:10.1128/JVI.00667-16

Grubaugh ND et al. (2019). An amplicon-based sequencing framework for accurately measuring intrahost virus diversity using PrimalSeq and iVar. Genome Biology 20:8. doi:10.1186/s13059-018-1618-7

Bull RA et al. (2020). Analytical validity of nanopore sequencing for rapid SARS-CoV-2 genome analysis. Nature Communications 11:6272. doi:10.1038/s41467-020-20075-6

Wenger AM et al. (2019). Accurate circular consensus long-read sequencing improves variant detection and assembly of a human genome. Nature Biotechnology 37:1155–1162. doi:10.1038/s41587-019-0217-9

See Also

ISNVFilter for the underlying constructor, flagISNV for applying filters.

Examples

# Illumina defaults
adaptiveFilter("illumina")

# ONT with strand bias disabled
filt <- adaptiveFilter("ont")
slot(filt, "maxStrandBias")  # Inf

# PacBio HiFi
adaptiveFilter("pacbio_hifi")

# Generic (same as ISNVFilter())
adaptiveFilter("generic")

Description

For each SNV in a WithinHostExperiment, determines the affected codon and amino acid change using CDS coordinates from a GFF3 file and the reference genome sequence from a FASTA file. Adds GFF_FEATURE, REF_CODON, ALT_CODON, REF_AA, ALT_AA, and AA_CLASS (Synonymous / Nonsynonymous / Nonsense) columns to mcols(rowRanges).

Usage

annotateCodonChange(whe, gff, refFasta, nameCol = "gene")

Arguments

Details

This is a pure-R implementation using the standard genetic code (NCBI translation table 1). No external annotation tools are required. For complex genomes with overlapping reading frames, ribosomal frameshift sites, or non-standard codon usage, use SnpEff or ANNOVAR upstream.

The algorithm:

1. Parse CDS features from GFF3 (requires column 8 = phase).

2. For each variant overlapping a CDS, compute the codon position within the CDS (accounting for strand and phase).

3. Extract the reference codon from the FASTA, substitute the variant base, translate both codons.

4. Classify: Synonymous (same AA), Nonsynonymous (different AA), or Nonsense (stop codon gained/lost).

Only single-nucleotide substitutions are annotated. Indels and multi-nucleotide variants are marked as NA.

Value

A WithinHostExperiment with additional columns in mcols(rowRanges): GFF_FEATURE, REF_CODON, ALT_CODON, REF_AA, ALT_AA, AA_CLASS.

Examples

## Requires a reference FASTA and GFF3 with CDS entries
ref <- system.file("extdata", "test_ref.fa",
    package = "WithinHostExperiment")
if (file.exists(ref)) {
    gff_file <- tempfile(fileext = ".gff3")
    writeLines(c("##gff-version 3",
        "test_segment\t.\tCDS\t1\t999\t.\t+\t0\tgene=testGene"),
        gff_file)
    vcf <- system.file("extdata", "test_donor.vcf",
        package = "WithinHostExperiment")
    whe <- readWithinHost(vcf,
        colData = S4Vectors::DataFrame(sample_id = "donor"),
        caller = "ivar")
    whe <- annotateCodonChange(whe, gff_file, ref)
    mc <- S4Vectors::mcols(SummarizedExperiment::rowRanges(whe))
    head(mc[, c("GFF_FEATURE", "REF_AA", "ALT_AA", "AA_CLASS")])
}

Description

Reads a GFF3 file and adds gene/feature annotations to mcols(rowRanges) of a WithinHostExperiment. This automates the annotation step required by dndsWithinHost, which expects REF_AA and ALT_AA columns.

Usage

annotateFromGFF(whe, gff, featureType = "gene", nameCol = "Name")

Arguments

Details

The function overlaps variant positions with GFF features of the specified type (default: "gene") and adds a GFF_FEATURE column containing the gene name or ID. For coding-region annotation (amino acid changes), a reference FASTA is additionally required; this function provides the positional overlap only.

The GFF3 file is parsed with base R (no external dependency). Only lines matching featureType in column 3 are used. Attribute parsing follows the GFF3 specification (key=value pairs separated by semicolons).

For full codon-aware annotation (predicting amino acid changes), use dedicated tools such as SnpEff or ANNOVAR upstream, then import the annotated VCF via readWithinHost.

Value

A WithinHostExperiment with a GFF_FEATURE column added to mcols(rowRanges).

Examples

## Create a minimal GFF3 for demonstration
gff_file <- tempfile(fileext = ".gff3")
writeLines(c(
  "##gff-version 3",
  "test_segment\t.\tgene\t1\t500\t.\t+\t.\tName=geneA",
  "test_segment\t.\tgene\t501\t1000\t.\t+\t.\tName=geneB"
), gff_file)
vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
whe <- annotateFromGFF(whe, gff_file)
S4Vectors::mcols(SummarizedExperiment::rowRanges(whe))$GFF_FEATURE

Description

Converts a WithinHostExperiment to the input format expected by the ViralBottleneck package.

Usage

asViralBottleneckInput(
  whe,
  pairId = NULL,
  threshold = 0.03,
  usePassedOnly = TRUE
)

Arguments

Details

For comprehensive bottleneck analysis with multiple methods (presence-absence, KL divergence, binomial, beta-binomial approximate/exact, and Wright-Fisher), use the ViralBottleneck package (Zheng et al. 2025, Virus Evolution 11:veaf071). This function prepares a simplified input suitable for the beta-binomial methods.

Value

A data.frame with columns: pos, donor_freq, recipient_freq – a simplified format for beta-binomial MLE. For the full ViralBottleneck input format (per-base read counts), prepare data externally using the read count assays (altCount, refCount).

Examples

vcf1 <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
vcf2 <- system.file("extdata", "test_recipient.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(c(vcf1, vcf2),
    colData = S4Vectors::DataFrame(
        sample_id = c("donor", "recipient"),
        role = c("donor", "recipient"),
        pair_id = c("p1", "p1")),
    caller = "ivar")
asViralBottleneckInput(whe, pairId = "p1")

Description

If a truth WithinHostExperiment is provided, computes true positives, false positives, false negatives, precision, recall, and F1 for each caller. If no truth is provided, computes pairwise agreement using calcCallerConcordance.

Usage

benchmarkCallers(whe_list, truth = NULL)

Arguments

Value

A DataFrame. When truth is provided:

caller

Caller name.

n_total

Total sites called by this caller.

TP

True positives.

FP

False positives.

FN

False negatives.

precision

TP / (TP + FP).

recall

TP / (TP + FN).

F1

Harmonic mean of precision and recall.

When truth is NULL, the result of calcCallerConcordance is returned.

Examples

gr1 <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200, 300), width = 1))
S4Vectors::mcols(gr1)$ref <- c("A", "C", "G")
S4Vectors::mcols(gr1)$alt <- c("T", "G", "A")
whe1 <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.1, 0.2, 0.3), ncol = 1)),
    rowRanges = gr1,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
benchmarkCallers(list(caller1 = whe1), truth = whe1)

Description

Constructs a WithinHostSFS object from the allele frequencies in a WithinHostExperiment. One SFS is built per sample.

Usage

buildSFS(
  whe,
  sampleIdx = 1L,
  nBins = 20L,
  fold = TRUE,
  usePassedOnly = TRUE,
  threshold = 0.03,
  genomeLength = NULL
)

Arguments

Value

A WithinHostSFS object.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
sfs <- buildSFS(whe, genomeLength = 1000L)
sfs

Description

Convenience wrapper that calls buildSFS for each sample and returns a named list of WithinHostSFS objects.

Usage

buildSFSList(
  whe,
  nBins = 20L,
  fold = TRUE,
  usePassedOnly = TRUE,
  threshold = 0.03,
  genomeLength = NULL
)

Arguments

Value

A named list of WithinHostSFS objects.

Examples

data(example_whe)
sfs_list <- buildSFSList(example_whe, genomeLength = 1000L)
length(sfs_list)

Description

For a list of WithinHostExperiment objects (one per variant caller), computes pairwise concordance statistics including shared sites, caller-unique sites, and Jaccard index.

Usage

calcCallerConcordance(whe_list)

Arguments

Value

A DataFrame with columns:

caller1

First caller name.

caller2

Second caller name.

n_shared

Number of sites detected by both callers.

n_only1

Sites unique to caller1.

n_only2

Sites unique to caller2.

jaccard

Jaccard similarity index (shared / union).

References

Roder AE et al. (2023). Optimized quantification of intra-host viral diversity in SARS-CoV-2 and influenza virus sequence data. mBio 14:e01046-23. doi:10.1128/mbio.01046-23

Examples

gr1 <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200, 300), width = 1))
S4Vectors::mcols(gr1)$ref <- c("A", "C", "G")
S4Vectors::mcols(gr1)$alt <- c("T", "G", "A")
whe1 <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.1, 0.2, 0.3), ncol = 1)),
    rowRanges = gr1,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
gr2 <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200, 400), width = 1))
S4Vectors::mcols(gr2)$ref <- c("A", "C", "T")
S4Vectors::mcols(gr2)$alt <- c("T", "G", "C")
whe2 <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.12, 0.22, 0.15), ncol = 1)),
    rowRanges = gr2,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
calcCallerConcordance(list(ivar = whe1, lofreq = whe2))

Description

Generic for computing diversity metrics from within-host variant data. Methods should return a DataFrame with one row per sample.

Calculates diversity metrics per sample from a WithinHostExperiment. Implements standard population genetics indices adapted for deep sequencing iSNV data, following Lauring (2020).

Usage

calcDiversity(x, ...)

## S4 method for signature 'WithinHostExperiment'
calcDiversity(
  x,
  indices = c("shannon", "pi", "watterson", "richness", "simpson", "chao1"),
  genomeLength = NULL,
  usePassedOnly = TRUE,
  ...
)

Arguments

Value

A DataFrame of diversity metrics.

A DataFrame with one row per sample and columns for each requested index.

References

Lauring AS (2020). Within-Host Viral Diversity: A Window into Viral Evolution. Annu Rev Virol 7:63-81. doi:10.1146/annurev-virology-010320-061642

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
calcDiversity(whe, genomeLength = 1000L)

Description

Generic for computing neutrality tests from within-host variant data. Methods should return a DataFrame with one row per sample.

Calculates neutrality test statistics per sample from a WithinHostExperiment. Implements Tajima's D (Tajima 1989) and Fu's Fs (Fu 1997), adapted for deep sequencing iSNV data.

Usage

calcNeutralityTests(x, ...)

## S4 method for signature 'WithinHostExperiment'
calcNeutralityTests(
  x,
  genomeLength = NULL,
  indices = c("tajimaD", "fusFs"),
  usePassedOnly = TRUE,
  ...
)

Arguments

Value

A DataFrame of neutrality statistics.

A DataFrame with one row per sample and columns sample_id plus one column per requested test.

References

Tajima F (1989). Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics 123:585-595.

Fu YX (1997). Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics 147:915-925.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
calcNeutralityTests(whe, genomeLength = 1000L)

Description

Calculate shared variant statistics

Usage

calcSharedVariants(whe, pairId)

Arguments

Value

A named list: n_donor_only, n_recipient_only, n_shared, n_total, sharing_rate.

Examples

vcf1 <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
vcf2 <- system.file("extdata", "test_recipient.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(c(vcf1, vcf2),
    colData = S4Vectors::DataFrame(
        sample_id = c("donor", "recipient"),
        role = c("donor", "recipient"),
        pair_id = c("p1", "p1")),
    caller = "ivar")
calcSharedVariants(whe, pairId = "p1")

Description

Estimates total species (variant) richness using the Chao1 estimator (Chao 1984):

$\hat{S}_{Chao1} = S_{obs} + \frac{f_1^2}{2 f_2}$

where $S_{obs}$ is the observed number of variants, $f_1$ is the number of singletons (variants with altCount == 1), and $f_2$ is the number of doubletons (variants with altCount == 2). When $f_2 = 0$, the bias-corrected form is used: $S_{obs} + f_1 (f_1 - 1) / 2$.

Usage

chao1ISNV(counts, detected)

Arguments

Value

Numeric scalar. Chao1 richness estimate.

References

Chao A (1984). Nonparametric estimation of the number of classes in a population. Scand J Statist 11:265-270.

Examples

chao1ISNV(c(1, 1, 2, 5, 10), detected = 5)

Description

Coerce WithinHostExperiment to data.frame or VRanges

Value

A data.frame or VRanges object containing the variant data from the WithinHostExperiment object.

Examples

data(example_whe)
df <- as(example_whe, "data.frame")
head(df)

vr <- as(example_whe, "VRanges")
vr

Description

Compare diversity between sample groups

Usage

compareDiversity(whe, groupBy, index = "shannon", genomeLength = NULL)

Arguments

Value

A DataFrame with test results.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
## Single-group, returns diversity table:
compareDiversity(whe, groupBy = "sample_id", index = "shannon")

Description

Computes summary statistics and a chi-squared test for homogeneity between two WithinHostSFS objects (e.g., before vs. after QC, or two timepoints).

Usage

compareSFS(sfs1, sfs2)

Arguments

Details

The test compares raw variant counts, so apply it to uncorrected spectra: the counts of a bias-corrected spectrum are rescaled and no longer multinomial. Bins where both spectra are empty are dropped, and when only two bins remain chisq.test applies Yates' continuity correction. The chi-squared approximation is unreliable when an expected count falls below 5; instead of a warning, that condition is reported in chisq_approx_ok, and the p-value should then be read with caution.

Value

A named list with components:

chisq_stat

Chi-squared statistic.

chisq_p

Chi-squared p-value.

df

Degrees of freedom.

min_expected

Smallest expected count under homogeneity, or NA when fewer than two bins remain and no test is run.

chisq_approx_ok

TRUE when every expected count is at least 5, the condition under which chisq.test does not warn; NA alongside the other NA diagnostics when no test is run.

n1

Total variants in sfs1.

n2

Total variants in sfs2.

proportions1

Proportion in each bin, sfs1.

proportions2

Proportion in each bin, sfs2.

Examples

data(example_whe)
whe1 <- example_whe
sfs1 <- buildSFS(whe1, genomeLength = 1000L)
sfs2 <- buildSFS(whe1, sampleIdx = 2L, genomeLength = 1000L)
compareSFS(sfs1, sfs2)

Description

Takes a named list of WithinHostExperiment objects (one per variant caller) and identifies variant sites detected by multiple callers. Sites are matched by genomic position (seqnames + start) and alternative allele. Only sites detected by at least min_callers with a frequency spread no greater than freq_tolerance are retained.

Usage

consensusISNV(whe_list, min_callers = 2L, freq_tolerance = 0.05)

Arguments

Details

When technical replicates are unavailable, using a multi-caller consensus is recommended for reducing false positive iSNV calls (Roder et al. 2023).

Value

A WithinHostExperiment containing only consensus sites. mcols(rowRanges) includes:

n_callers

Number of callers detecting the site.

caller_names

Comma-separated caller names.

freq_spread

Range of allele frequencies across callers.

The altFreq assay contains the mean frequency across callers for each consensus site. It is the only assay of the result: read counts, depths and qcPass are not carried over, so run flagISNV on the individual callers' objects rather than on the consensus.

References

Roder AE et al. (2023). Optimized quantification of intra-host viral diversity in SARS-CoV-2 and influenza virus sequence data. mBio 14:e01046-23. doi:10.1128/mbio.01046-23

Examples

gr <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200, 300), width = 1))
S4Vectors::mcols(gr)$ref <- c("A", "C", "G")
S4Vectors::mcols(gr)$alt <- c("T", "G", "A")
whe1 <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.10, 0.20, 0.30), ncol = 1)),
    rowRanges = gr, colData = S4Vectors::DataFrame(sample_id = "S1"))
gr2 <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200, 400), width = 1))
S4Vectors::mcols(gr2)$ref <- c("A", "C", "T")
S4Vectors::mcols(gr2)$alt <- c("T", "G", "C")
whe2 <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.12, 0.18, 0.15), ncol = 1)),
    rowRanges = gr2, colData = S4Vectors::DataFrame(sample_id = "S1"))
consensus <- consensusISNV(list(caller1 = whe1, caller2 = whe2))
consensus

Description

Deep sequencing imposes a detection threshold (e.g., 3%): sites with true frequency below this threshold are unobservable. This inflates the apparent fraction of moderate-frequency variants relative to rare variants.

Usage

correctSFSBias(sfs, method = c("truncation", "binomial"))

Arguments

Details

This function applies a correction based on the expected number of unobserved sites under a neutral beta distribution model: given the observed SFS and depth-dependent detectability, it estimates the missing mass in the lowest frequency bins.

The "truncation" method assumes a uniform density of variant frequencies and multiplies every bin by $f_{max} / (f_{max} - \mathrm{threshold})$, the reciprocal of the observable fraction of the frequency range. It therefore rescales the spectrum without changing its shape. The "binomial" method uses the stored depth to compute $P(\mathrm{detect} | \nu, n) = P(\mathrm{Alt} \ge k_{\min} | n, \nu)$ where $k_{\min} = \lceil n \times \mathrm{threshold} \rceil$, and reweights each bin by the inverse of its detectability. Detectability is floored at 0.01, so a bin is inflated at most one hundredfold.

Value

A corrected WithinHostSFS with corrected = TRUE. Counts are no longer raw variant counts, and nSites changes accordingly, which also changes the number of segregating sites used by neutralityFromSFS. Use uncorrected spectra for compareSFS.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
sfs <- buildSFS(whe, genomeLength = 1000L)
sfs_corr <- correctSFSBias(sfs, method = "binomial")
sfs_corr

Description

Identifies iSNV sites where a high-frequency variant in one sample appears at suspiciously low frequency in another sample from the same sequencing run, consistent with index hopping or sample cross-contamination.

Usage

detectCrossContamination(
  whe,
  runCol = NULL,
  sourceMinFreq = 0.2,
  sinkMaxFreq = 0.05,
  sinkMinFreq = 0.005,
  usePassedOnly = TRUE
)

Arguments

Details

Index hopping occurs when adapters are exchanged between library molecules during pooled sequencing, causing reads from one sample to be misassigned to another. The signature: a variant at high frequency (e.g., 80%) in sample A appears at very low frequency (e.g., 0.5-2%) in sample B sequenced on the same flow cell, at a ratio consistent with the expected hopping rate (typically 0.1-2% for Illumina platforms).

This function does not require external metadata beyond what is in colData. For finer control, use the output DataFrame to implement study-specific contamination criteria.

Value

A DataFrame with columns:

source_sample

Sample where the variant is at high frequency.

sink_sample

Sample where the variant appears at suspiciously low frequency.

chrom

Chromosome/segment.

position

Genomic position.

ref

Reference allele.

alt

Alternative allele.

source_freq

Frequency in the source sample.

sink_freq

Frequency in the sink sample.

freq_ratio

sink_freq / source_freq – lower values are more suspicious.

run

Sequencing run (if runCol provided).

References

Bendall EE et al. (2022). SARS-CoV-2 genomic diversity in households highlights the challenges of sequence-based transmission inference. mSphere 7:e00400-22. doi:10.1128/msphere.00400-22

Costello M et al. (2018). Characterization and remediation of sample index swaps by non-redundant dual indexing on massively parallel sequencing platforms. BMC Genomics 19:332. doi:10.1186/s12864-018-4703-0

Examples

gr <- GenomicRanges::GRanges("seg",
    IRanges::IRanges(c(100, 200, 300), width = 1))
S4Vectors::mcols(gr)$ref <- c("A", "C", "G")
S4Vectors::mcols(gr)$alt <- c("T", "G", "A")
## Sample A has pos 100 at 80%; sample B has same pos at 1%
freq <- matrix(c(0.80, 0.05, 0.30,   # sample A
                 0.01, 0.40, 0.00),   # sample B
               nrow = 3, ncol = 2)
whe <- WithinHostExperiment(
    assays = list(altFreq = freq),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(
        sample_id = c("A", "B"),
        run = c("run1", "run1")))
detectCrossContamination(whe, runCol = "run")

Description

Identifies variants whose allele frequency crosses a threshold between consecutive time points within each host. A variant is considered emerging if it was below freqThreshold at one time point and increased by at least minIncrease at the next.

Usage

detectEmergingVariants(
  whe,
  hostCol = "host_id",
  timeCol = "timepoint",
  freqThreshold = 0.05,
  minIncrease = 0.05,
  usePassedOnly = TRUE
)

Arguments

Value

A DataFrame with columns:

host_id

Host identifier.

variant_key

Variant identifier (chrom:pos:alt).

time_from

Earlier time point.

time_to

Later time point.

freq_from

Frequency at the earlier time point.

freq_to

Frequency at the later time point.

delta_freq

Change in frequency.

Examples

gr <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200), width = 1))
S4Vectors::mcols(gr)$ref <- c("A", "C")
S4Vectors::mcols(gr)$alt <- c("T", "G")
freq <- matrix(c(0.02, 0.10, 0.12, 0.15), nrow = 2, ncol = 2)
whe <- WithinHostExperiment(
    assays = list(altFreq = freq),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(
        sample_id = c("S1_t1", "S1_t2"),
        host_id = c("H1", "H1"),
        timepoint = c(1, 2)))
detectEmergingVariants(whe)

Description

Estimates the ratio of nonsynonymous to synonymous divergence among the iSNVs of each sample with the counting method of Nei and Gojobori (1986). The numbers of synonymous ($S$) and nonsynonymous ($N$) sites are counted over the complete codons of every CDS in gff, using the reference sequence in refFasta and the standard genetic code. With $S_d$ synonymous and $N_d$ nonsynonymous iSNVs present in a sample,

$p_S = S_d / S, \quad p_N = N_d / N,$

each proportion is corrected for multiple hits with Jukes and Cantor (1969), $d = -\frac{3}{4}\log(1 - \frac{4}{3}p)$, and $d_N/d_S$ is their ratio. Under neutrality the expected ratio is 1; values below 1 suggest purifying selection and values above 1 diversifying selection, but within-host ratios computed from few iSNVs are highly variable and should be read with their counts.

Usage

dndsWithinHost(
  whe,
  gff,
  refFasta,
  gffFeatureCol = "GFF_FEATURE",
  refAACol = "REF_AA",
  altAACol = "ALT_AA",
  usePassedOnly = TRUE,
  nameCol = "gene"
)

Arguments

Details

Each iSNV present in a sample (non-missing, non-zero frequency, and passing QC when usePassedOnly = TRUE) is counted once, irrespective of its frequency. Amino acid classes come from annotateCodonChange, which should be run with the same gff and refFasta. Variants outside a CDS, or without an amino acid annotation, are not counted. Nonsense changes count as nonsynonymous. Codons shared by overlapping CDS features of the same gene are counted once.

Value

A DataFrame with one row per sample and gene carrying at least one counted iSNV (one row with gene = NA when a sample has none):

sample_id

Sample identifier.

nS, nN

Synonymous and nonsynonymous iSNVs in the sample.

S_sites, N_sites

Synonymous and nonsynonymous sites over all CDS features.

pS, pN

Proportions $S_d/S$ and $N_d/N$.

dNdS

Jukes-Cantor corrected $d_N/d_S$ (NA when $d_S$ is zero or a proportion reaches 0.75).

gene, gene_nS, gene_nN, gene_S_sites, gene_N_sites, gene_dNdS

The same quantities for each gene.

References

Nei M, Gojobori T (1986). Simple methods for estimating the numbers of synonymous and nonsynonymous nucleotide substitutions. Mol Biol Evol 3:418-426. doi:10.1093/oxfordjournals.molbev.a040410

Jukes TH, Cantor CR (1969). Evolution of protein molecules. In Munro HN (ed.), Mammalian Protein Metabolism, pp. 21-132. Academic Press, New York.

See Also

annotateCodonChange

Examples

ref <- tempfile(fileext = ".fa")
writeLines(c(">seg1", "ATGCTGAAAGGGTAA"), ref)
gff <- tempfile(fileext = ".gff3")
writeLines(c("##gff-version 3",
    "seg1\t.\tCDS\t1\t15\t.\t+\t0\tgene=geneA"), gff)
gr <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(1, 6), width = 1))
S4Vectors::mcols(gr)$ref <- c("A", "G")
S4Vectors::mcols(gr)$alt <- c("G", "A")
whe <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.1, 0.2), ncol = 1)),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
whe <- annotateCodonChange(whe, gff, ref)
dndsWithinHost(whe, gff = gff, refFasta = ref)

Description

Given a vector of allele frequencies observed over equally-spaced time points, estimates the selection coefficient $s$ using the change in logit-transformed frequency per unit time:

$s = \frac{\Delta \mathrm{logit}(p)}{\Delta t}$

where $\mathrm{logit}(p) = \log(p / (1 - p))$.

Usage

estimateSelectionCoefficient(freqs, times = NULL, generation_time = 1)

Arguments

Details

This assumes a simple deterministic selection model. Each consecutive pair of time points yields one estimate of $s$.

Value

Numeric vector of selection coefficient estimates, one per consecutive time interval. Length is length(freqs) - 1.

References

Feder AF et al. (2016). More effective drugs lead to harder selective sweeps in the evolution of drug resistance in HIV-1. eLife 5:e10670. doi:10.7554/eLife.10670

Examples

freqs <- c(0.05, 0.10, 0.25, 0.50)
estimateSelectionCoefficient(freqs)
estimateSelectionCoefficient(freqs, times = c(0, 3, 7, 14))

Description

Estimates transmission bottleneck size using one of three methods:

exact_bb

Beta-binomial likelihood over the number of transmitted variant virions, which models finite recipient read depth, in the spirit of Sobel Leonard et al. (2017). Requires altCount and totalDepth assays. Two deviations from the published method are deliberate: when all or none of the transmitted virions carry the variant, the recipient read counts are modelled with a fixed sequencing error probability of threshold / 2, and the recipient detection threshold is not integrated out.

wright_fisher

Approximate score in which one generation of binomial drift is combined with a beta kernel of concentration max(Nb, 2) around the drifted frequency. It mixes a density for detected variants with probabilities for lost variants, so it is a heuristic rather than a likelihood; prefer exact_bb or presence_absence for inference.

presence_absence

Presence-absence likelihood using $P(\mathrm{detect} | Nb, \nu) = 1 - (1-\nu)^{Nb}$.

Usage

exactBottleneck(
  whe,
  pairId,
  maxNb = 500L,
  method = c("exact_bb", "wright_fisher", "presence_absence"),
  nboot = 0L,
  threshold = 0.03
)

Arguments

Value

A named list with components:

Nb

Numeric. Maximum likelihood estimate of bottleneck size over the grid 1:maxNb. A warning is issued when the estimate reaches maxNb, because the value is then censored by the grid.

ci

Numeric vector of length 2. 95% confidence interval (likelihood-ratio or bootstrap percentile), also bounded by maxNb.

loglik

Numeric. Log-likelihood at the MLE.

method

Character. Method used.

n_variants

Integer. Number of donor variants used.

pair_id

Character. Pair identifier.

References

Sobel Leonard A et al. (2017). Transmission Bottleneck Size Estimation from Pathogen Deep-Sequencing Data. J Virol 91:e00171-17. doi:10.1128/JVI.00171-17

Examples

vcf1 <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
vcf2 <- system.file("extdata", "test_recipient.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(c(vcf1, vcf2),
    colData = S4Vectors::DataFrame(
        sample_id = c("donor", "recipient"),
        role = c("donor", "recipient"),
        pair_id = c("p1", "p1")),
    caller = "ivar")
exactBottleneck(whe, pairId = "p1", maxNb = 50L,
    method = "presence_absence")

Description

A small synthetic WithinHostExperiment with 3 samples (donor, recipient, independent) containing iSNV data at positions 33 to 945 of a synthetic genome segment. The donor and recipient form a transmission pair (pair_id = "pair_1").

Usage

data(example_whe)

Format

A WithinHostExperiment with 28 variant sites and 3 samples. Assays: altFreq, totalDepth, refCount, altCount, qcPass.

Source

Synthetic. Built by inst/scripts/save_example_data.R from test_donor.vcf, test_recipient.vcf and test_independent.vcf in inst/extdata, read with readWithinHost.

Examples

data(example_whe)
example_whe
transmissionPairs(example_whe)

Description

For each transmission pair, creates a table of variant frequencies in the donor and recipient.

Usage

exportPairFrequencies(whe, pairId = NULL, usePassedOnly = TRUE)

Arguments

Value

A DataFrame with columns: chrom, position, ref, alt, donor_freq, recipient_freq, shared (logical).

Examples

vcf1 <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
vcf2 <- system.file("extdata", "test_recipient.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(c(vcf1, vcf2),
    colData = S4Vectors::DataFrame(
        sample_id = c("donor", "recipient"),
        role = c("donor", "recipient"),
        pair_id = c("p1", "p1")),
    caller = "ivar")
exportPairFrequencies(whe)

Description

Fits a linear model (lm) relating a diversity metric to one or more predictors. When group_col is supplied it is included as a fixed effect (a full mixed-effects model would require the nlme dependency).

Usage

fitDiversityModel(div_df, formula, group_col = NULL)

Arguments

Value

A list with components:

model

The fitted lm object.

summary_df

A data.frame of coefficient estimates, standard errors, t-values, and p-values.

residuals

Numeric vector of model residuals.

Examples

df <- data.frame(
    shannon = c(1.2, 1.5, 0.8, 1.1, 1.9, 0.7),
    timepoint = c(1, 2, 3, 1, 2, 3),
    group = c("A", "A", "A", "B", "B", "B")
)
res <- fitDiversityModel(df, shannon ~ timepoint)
res_grouped <- fitDiversityModel(df, shannon ~ timepoint,
    group_col = "group")

Description

Generic for applying quality control filters to within-host variant data. The default method operates on WithinHostExperiment objects. Other packages can define methods for alternative containers.

Applies quality control filters to a WithinHostExperiment. Instead of removing variants, this function sets the corresponding entries in the qcPass assay to FALSE. Each filtering step is recorded in the qcLog for full reproducibility.

Usage

flagISNV(x, filter = ISNVFilter(), ...)

## S4 method for signature 'WithinHostExperiment'
flagISNV(x, filter = ISNVFilter(), ...)

Arguments

Details

Filtering steps are applied in order:

1. Minimum read depth (minDepth)

2. Minimum alternative allele frequency (minFreq)

3. Maximum alternative allele frequency (maxFreq)

4. Minimum alternative allele read count (minAltReads)

5. Strand bias (maxStrandBias) – only if fwdAltCount and revAltCount assays exist

The minBaseQual, minMapQual and replicateConc slots of the filter are recorded for provenance but not applied here; a non-default value raises a warning. Use flagReplicateDiscordance for replicate concordance.

Value

The filtered object (same class as input).

A WithinHostExperiment with updated qcPass assay and qcLog.

See Also

ISNVFilter for filter configuration, passedISNV for extracting passed variants, qcSummary for QC report, flagReplicateDiscordance for replicate QC.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
whe <- flagISNV(whe, ISNVFilter(minDepth = 500L))
qcLog(whe)
passedISNV(whe)

Description

For amplicon-based sequencing data, marks variants that fall within primer binding sites in the qcPass assay.

Usage

flagPrimerSites(whe, primers, action = c("flag", "remove"))

Arguments

Value

A WithinHostExperiment.

Examples

library(GenomicRanges)
vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
primers <- GRanges("test_segment", IRanges::IRanges(
    start = c(40, 880), width = c(25, 30)))
whe <- flagPrimerSites(whe, primers)
qcLog(whe)

Description

Identifies variant sites where technical replicates disagree and flags them in the qcPass assay.

Usage

flagReplicateDiscordance(whe, freqTolerance = 0.05, requireBothDetected = TRUE)

Arguments

Details

For each replicate group, the function compares altFreq between replicates and flags sites where the difference exceeds freqTolerance. Sites are flagged in ALL samples of the replicate group.

Roder et al. (2023) reported that the choice of variant caller and the use of replicate sequencing had the largest effects on single-nucleotide variant discovery. An alternative approach is replicate intersection (requiring detection in both replicates); use requireBothDetected = TRUE with freqTolerance = Inf to achieve this.

Value

A WithinHostExperiment with updated qcPass and qcLog.

References

Roder AE et al. (2023). Optimized quantification of intra-host viral diversity in SARS-CoV-2 and influenza virus sequence data. mBio 14:e01046-23. doi:10.1128/mbio.01046-23

See Also

flagISNV for standard QC filtering.

Examples

rep1 <- system.file("extdata", "test_donor_rep1.vcf",
    package = "WithinHostExperiment")
rep2 <- system.file("extdata", "test_donor_rep2.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(c(rep1, rep2),
    colData = S4Vectors::DataFrame(
        sample_id = c("rep1", "rep2"),
        replicate_group = c("donor_reps", "donor_reps"),
        replicate_id = c("1", "2")),
    caller = "ivar")
whe <- flagReplicateDiscordance(whe, freqTolerance = 0.05)
qcLog(whe)

Description

Applies Strand Coverage Imbalance filtering to a WithinHostExperiment. Sites whose absolute SCI score exceeds maxSCI are flagged in the qcPass assay.

Usage

flagStrandBiasSCI(whe, maxSCI = 1.5)

Arguments

Details

Unlike the simple fold-difference strand bias in flagISNV, this function uses the SCI score, which normalises by the global strand coverage ratio (see sciStrandBias). This is preferred for amplicon-based sequencing where forward and reverse strand coverage is inherently unbalanced.

Value

A WithinHostExperiment with updated qcPass and qcLog.

References

Mostefai F et al. (2024). Refining SARS-CoV-2 intra-host variation by leveraging large-scale sequencing data. NAR Genomics and Bioinformatics 6:lqae145. doi:10.1093/nargab/lqae145

See Also

sciStrandBias for the underlying computation, flagISNV for standard QC.

Examples

library(SummarizedExperiment)
library(GenomicRanges)
gr <- GRanges("seg", IRanges::IRanges(start = 1:3, width = 1))
assays <- list(
    altFreq      = matrix(c(0.1, 0.05, 0.2), ncol = 1),
    totalDepth   = matrix(c(1000L, 1000L, 1000L), ncol = 1),
    qcPass       = matrix(TRUE, nrow = 3, ncol = 1),
    fwdAltCount  = matrix(c(50L, 2L, 30L), ncol = 1),
    revAltCount  = matrix(c(45L, 40L, 28L), ncol = 1),
    fwdTotalCount = matrix(c(500L, 500L, 500L), ncol = 1),
    revTotalCount = matrix(c(480L, 480L, 480L), ncol = 1)
)
whe <- WithinHostExperiment(
    rowRanges = gr, assays = assays,
    colData = S4Vectors::DataFrame(sample_id = "s1"))
whe <- flagStrandBiasSCI(whe, maxSCI = 1.5)
qcLog(whe)

Description

For longitudinal samples from the same host, identifies iSNVs that appear at only one timepoint (transient) versus those detected at multiple consecutive timepoints (persistent). Transient iSNVs are more likely to be sequencing artefacts; persistent ones are more likely genuine within-host variants.

Usage

flagTemporalInconsistency(
  whe,
  hostCol = "host_id",
  timeCol = "timepoint",
  minTimepoints = 2L,
  minConsecutive = NULL,
  usePassedOnly = TRUE
)

Arguments

Details

This check applies only where a host was sampled more than once. The biological rationale: genuine within-host variants arise from viral replication and should persist (or change smoothly) across consecutive timepoints. An iSNV that appears once at 5% and is absent at all other timepoints is more likely a sequencing or library preparation artefact.

The function records the filtering step in qcLog with step name "temporal_consistency".

Value

A WithinHostExperiment with updated qcPass assay and a new column temporal_class in mcols(rowRanges). Sites are classified within each host; a transient site is flagged only in the samples of the hosts where it is transient. The site-level temporal_class is "persistent" if the site is persistent in at least one host, "transient" if it is transient in every host where it was detected, and NA where it was never detected.

References

Farjo M et al. (2024). Within-host evolutionary dynamics and tissue compartmentalization during acute SARS-CoV-2 infection. J. Virol. 98:e01618-23. doi:10.1128/jvi.01618-23

Examples

gr <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200, 300), width = 1))
S4Vectors::mcols(gr)$ref <- c("A", "C", "G")
S4Vectors::mcols(gr)$alt <- c("T", "G", "A")
freq <- matrix(c(0.1, NA,  0.2,   # var1: 2 timepoints
                 0.05, NA,  NA,    # var2: 1 timepoint (transient)
                 0.15, 0.2, 0.1),  # var3: 3 timepoints
               nrow = 3, ncol = 3)
whe <- WithinHostExperiment(
    assays = list(altFreq = freq),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(
        sample_id = c("t1", "t2", "t3"),
        host_id = c("H1", "H1", "H1"),
        timepoint = c(1, 2, 3)))
whe <- flagTemporalInconsistency(whe, minTimepoints = 2L)
S4Vectors::mcols(SummarizedExperiment::rowRanges(whe))$temporal_class

Description

Computes Fu's Fs (Fu 1997), a test of neutrality based on the probability of observing at least the observed number of alleles given the estimate of theta from nucleotide diversity. Under the infinite-sites model:

$F_s = \ln\!\left(\frac{S'}{1 - S'}\right)$

where $S'$ is the probability of observing a number of alleles greater than or equal to the observed value, given $\theta = \hat{\pi}$. A large negative Fs indicates an excess of alleles (rare variants), consistent with population expansion or genetic hitchhiking.

Usage

fusFs(S, n, theta_pi)

Arguments

Details

The probability is approximated using Ewens' sampling formula for the expected number of alleles.

Value

Numeric scalar. Fu's Fs, or NA with a warning when $S + 1 > n$: the statistic uses $S + 1$ as the number of haplotypes, which cannot exceed the number of sequences sampled.

References

Fu YX (1997). Statistical tests of neutrality of mutations against population growth, hitchhiking and background selection. Genetics 147:915-925.

Examples

fusFs(S = 10, n = 100, theta_pi = 5)

Description

Constructs a quality control filter configuration for use with flagISNV.

Usage

ISNVFilter(
  minDepth = 100L,
  minFreq = 0.03,
  maxFreq = 0.5,
  minAltReads = 10L,
  maxStrandBias = 10,
  minBaseQual = 20L,
  minMapQual = 20L,
  replicateConc = FALSE
)

Arguments

Value

An ISNVFilter object.

Examples

ISNVFilter()
ISNVFilter(minDepth = 500L, minFreq = 0.05)

Description

An S4 class holding parameters for filtering intra-host single nucleotide variants. Used by flagISNV.

Usage

## S4 method for signature 'ISNVFilter'
show(object)

Arguments

Details

Default values are based on commonly used thresholds in the within-host pathogen diversity literature:

minDepth = 100

Minimum coverage for variant detection. This is a permissive floor: Grubaugh et al. (2019, Genome Biology) measured iSNVs above 3% accurately from at least 1,000 RNA copies sequenced to at least 400x.

minFreq = 0.03

3% minor allele frequency threshold, the lower limit of accurate iSNV measurement in Grubaugh et al. (2019); Popa et al. (2020, Science Translational Medicine) used 2%.

maxFreq = 0.50

Variants above 50% are considered consensus-level changes, not intra-host variants.

minAltReads = 10

Minimum supporting reads to reduce false positives from sequencing errors.

maxStrandBias = 10

Maximum fold-difference between forward and reverse strand alt allele counts. Excessive strand bias suggests systematic error.

Value

An object of class ISNVFilter. Objects are created with the constructor ISNVFilter; the show method prints the filter thresholds and returns NULL invisibly.

Slots

minDepth

Integer scalar. Minimum read depth at iSNV site (default: 100).

minFreq

Numeric scalar. Minimum alternative allele frequency (default: 0.03).

maxFreq

Numeric scalar. Maximum alternative allele frequency; variants above this are excluded as near-fixed (default: 0.50).

minAltReads

Integer scalar. Minimum reads supporting the alternative allele (default: 10).

maxStrandBias

Numeric scalar. Maximum strand bias as fold-difference between strands (default: 10.0).

minBaseQual

Integer scalar. Minimum mean base quality (default: 20). Recorded for provenance only: flagISNV does not apply it, because the object model has no per-site quality assay. Filter on quality in the variant caller instead.

minMapQual

Integer scalar. Minimum mean mapping quality (default: 20). Recorded for provenance only, as for minBaseQual.

replicateConc

Logical scalar. Whether replicate concordance is required (default: FALSE). Recorded for provenance only; apply it with flagReplicateDiscordance.

References

Roder AE et al. (2023). Optimized quantification of intra-host viral diversity in SARS-CoV-2 and influenza virus sequence data. mBio 14:e01046-23. doi:10.1128/mbio.01046-23

McCrone JT, Lauring AS (2016). Measurements of intrahost viral diversity are extremely sensitive to systematic errors in variant calling. J Virol 90:6884-6895. doi:10.1128/JVI.00667-16

Grubaugh ND et al. (2019). An amplicon-based sequencing framework for accurately measuring intrahost virus diversity using PrimalSeq and iVar. Genome Biology 20:8. doi:10.1186/s13059-018-1618-7

Mostefai F et al. (2024). Refining SARS-CoV-2 intra-host variation by leveraging large-scale sequencing data. NAR Genomics and Bioinformatics 6:lqae145. doi:10.1093/nargab/lqae145

See Also

ISNVFilter for the constructor, flagISNV for applying filters.

Examples

# Default filter
ISNVFilter()

# Stringent filter
ISNVFilter(minDepth = 500L, minFreq = 0.05)

Description

Thin wrapper around p.adjust for convenient multiple testing correction of p-value vectors.

Usage

multitestCorrection(pvalues, method = "BH")

Arguments

Value

Numeric vector of adjusted p-values (same length as pvalues).

Examples

raw_p <- c(0.001, 0.04, 0.03, 0.5)
multitestCorrection(raw_p)
multitestCorrection(raw_p, method = "bonferroni")

Description

Computes Tajima's D, Fu's Fs, and Fay and Wu's H from a single WithinHostSFS object. This is the recommended interface for neutrality testing: build the SFS once, then derive all statistics from it.

Usage

neutralityFromSFS(sfs, n = NULL)

Arguments

Details

Fay and Wu's H measures an excess of high-frequency derived variants, which can indicate positive selection or genetic hitchhiking (Fay and Wu 2000). It is defined as:

$H = \hat{\pi} - \hat{\theta}_H$

where $\hat{\theta}_H = \frac{2}{n(n-1)} \sum_i i^2 \xi_i$ and $\xi_i$ is the count of sites with $i$ derived alleles. A negative H indicates an excess of high-frequency variants.

Value

A named list with components:

S

Number of segregating sites.

n

Effective sample size used.

pi

Nucleotide diversity (per site).

theta_W

Watterson's theta.

theta_H

Fay and Wu's theta_H.

tajimaD

Tajima's D.

fusFs

Fu's Fs.

fayWuH

Fay and Wu's H.

sfs_corrected

Whether the input SFS was bias-corrected.

References

Fay JC, Wu CI (2000). Hitchhiking under positive Darwinian selection. Genetics 155:1405-1413.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
sfs <- buildSFS(whe, genomeLength = 1000L, fold = FALSE)
neutralityFromSFS(sfs)

Description

Non-parametric permutation test comparing a summary statistic between two groups. The observed test statistic is compared to a null distribution obtained by randomly permuting group labels.

Usage

permutationTest(values, groups, n_perm = 999L, statistic = "mean_diff")

Arguments

Value

A list with components:

observed

Numeric scalar. The observed test statistic.

p_value

Numeric scalar. Two-sided permutation p-value.

n_perm

Integer. Number of permutations performed.

Examples

set.seed(42)
vals <- c(rnorm(10, 2), rnorm(10, 3))
grps <- rep(c("A", "B"), each = 10)
permutationTest(vals, grps, n_perm = 499)

Description

Estimates nucleotide diversity following Nei & Li (1979), adapted for deep sequencing allele frequencies:

$\pi = \frac{n}{n-1} \cdot \frac{2}{L} \sum_i p_i(1 - p_i)$

where $p_i$ is the alternative allele frequency at site $i$, $L$ is the genome length, and $n$ is the number of reads sampled at a site, for which mean read depth stands in. The $n/(n-1)$ correction removes the finite-sample bias of estimating a frequency from $n$ draws (Nei 1987, Eq. 10.5). This is a different use of "sample size" from the one in tajimaD, where reads must not stand in for independent lineages. When meanDepth is not provided, the uncorrected estimator is returned (equivalent to infinite sample size).

Usage

piISNV(freq, genomeLength, meanDepth = NULL)

Arguments

Value

Numeric scalar. Nucleotide diversity.

References

Nei M, Li WH (1979). Mathematical model for studying genetic variation in terms of restriction endonucleases. Proc Natl Acad Sci USA 76:5269-5273. doi:10.1073/pnas.76.10.5269

Nei M (1987). Molecular Evolutionary Genetics. Columbia University Press, New York. Chapter 10.

Examples

piISNV(c(0.1, 0.3, 0.5), genomeLength = 13588)
piISNV(c(0.1, 0.3, 0.5), genomeLength = 13588, meanDepth = 1000)

Description

Creates a publication-quality histogram or density plot of alternative allele frequencies, optionally faceted by sample or metadata group.

Usage

plotFrequencySpectrum(
  whe,
  type = c("histogram", "density"),
  facetBy = NULL,
  usePassedOnly = TRUE,
  threshold = NULL,
  freqRange = NULL
)

Arguments

Value

A ggplot object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE)) {
    vcf <- system.file("extdata", "test_donor.vcf",
        package = "WithinHostExperiment")
    whe <- readWithinHost(vcf,
        colData = S4Vectors::DataFrame(sample_id = "d"),
        caller = "ivar")
    plotFrequencySpectrum(whe, threshold = 0.03)
}

Description

Creates a ggplot2 spaghetti plot of allele frequency over time for each variant within a specified host. Each line represents one variant, coloured with the Okabe-Ito palette (Wong 2011); with more than eight variants the colours repeat.

Usage

plotFrequencyTrajectory(
  whe,
  hostId,
  hostCol = "host_id",
  timeCol = "timepoint",
  usePassedOnly = TRUE
)

Arguments

Value

A ggplot object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE)) {
    gr <- GenomicRanges::GRanges("seg1",
        IRanges::IRanges(c(100, 200), width = 1))
    S4Vectors::mcols(gr)$ref <- c("A", "C")
    S4Vectors::mcols(gr)$alt <- c("T", "G")
    freq <- matrix(c(0.05, 0.10, 0.15, 0.20), nrow = 2, ncol = 2)
    whe <- WithinHostExperiment(
        assays = list(altFreq = freq),
        rowRanges = gr,
        colData = S4Vectors::DataFrame(
            sample_id = c("S1_t1", "S1_t2"),
            host_id = c("H1", "H1"),
            timepoint = c(1, 2)))
    plotFrequencyTrajectory(whe, hostId = "H1")
}

Description

Creates a publication-quality scatter plot comparing variant frequencies between a donor and recipient, highlighting shared and unique variants with distinct shapes and colours.

Usage

plotPairScatter(whe, pairId, threshold = 0.03)

Arguments

Value

A ggplot object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE)) {
    vcf1 <- system.file("extdata", "test_donor.vcf",
        package = "WithinHostExperiment")
    vcf2 <- system.file("extdata", "test_recipient.vcf",
        package = "WithinHostExperiment")
    whe <- readWithinHost(c(vcf1, vcf2),
        colData = S4Vectors::DataFrame(
            sample_id = c("donor", "recipient"),
            role = c("donor", "recipient"),
            pair_id = c("p1", "p1")),
        caller = "ivar")
    plotPairScatter(whe, pairId = "p1")
}

Description

Creates a publication-quality four-panel dashboard showing frequency distribution, read depth, QC filtering waterfall, and per-sample pass/fail summary.

Usage

plotQCDashboard(whe)

Arguments

Value

A patchwork composite ggplot object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE) &&
    requireNamespace("patchwork", quietly = TRUE)) {
    vcf <- system.file("extdata", "test_donor.vcf",
        package = "WithinHostExperiment")
    whe <- readWithinHost(vcf,
        colData = S4Vectors::DataFrame(sample_id = "d"),
        caller = "ivar")
    whe <- flagISNV(whe, ISNVFilter())
    plotQCDashboard(whe)
}

Description

Creates a publication-quality bar plot of a WithinHostSFS object.

Usage

plotSFS(sfs, normalize = FALSE)

Arguments

Value

A ggplot object.

Examples

if (requireNamespace("ggplot2", quietly = TRUE)) {
    vcf <- system.file("extdata", "test_donor.vcf",
        package = "WithinHostExperiment")
    whe <- readWithinHost(vcf,
        colData = S4Vectors::DataFrame(sample_id = "donor"),
        caller = "ivar")
    sfs <- buildSFS(whe, genomeLength = 1000L)
    plotSFS(sfs)
}

Description

Returns the quality control log recording all filtering steps applied to a WithinHostExperiment.

Returns a subset of the WithinHostExperiment containing only variant sites where at least one sample passed QC.

Identifies donor-recipient pairs based on role and pair_id columns in colData.

Usage

qcLog(x)

qcFilters(x)

passedISNV(x)

transmissionPairs(x)

qcLog(x) <- value

qcFilters(x) <- value

## S4 method for signature 'WithinHostExperiment'
qcLog(x)

## S4 method for signature 'WithinHostExperiment'
qcFilters(x)

## S4 replacement method for signature 'WithinHostExperiment'
qcLog(x) <- value

## S4 replacement method for signature 'WithinHostExperiment'
qcFilters(x) <- value

## S4 method for signature 'WithinHostExperiment'
passedISNV(x)

## S4 method for signature 'WithinHostExperiment'
transmissionPairs(x)

Arguments

Value

A DataFrame with columns: step, parameter, value, n_flagged, n_passed, timestamp.

A list of ISNVFilter objects.

A WithinHostExperiment (subset of rows).

A DataFrame with columns: donor (sample_id), recipient (sample_id), pair_id.

The modified WithinHostExperiment.

The modified WithinHostExperiment.

Examples

gr <- GenomicRanges::GRanges("seg", IRanges::IRanges(1:3, width = 1))
whe <- WithinHostExperiment(
    assays = list(altFreq = matrix(runif(3), ncol = 1)),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
qcLog(whe)
gr <- GenomicRanges::GRanges("seg", IRanges::IRanges(1:3, width = 1))
whe <- WithinHostExperiment(
    assays = list(altFreq = matrix(runif(3), ncol = 1)),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
passedISNV(whe)

Description

Produces per-sample summary statistics from a WithinHostExperiment.

Usage

qcSummary(whe)

Arguments

Value

A DataFrame with columns: sample_id, total_sites, passed_sites, flagged_sites, pass_rate, median_depth, median_freq.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
whe <- flagISNV(whe, ISNVFilter(minDepth = 500L))
qcSummary(whe)

Description

Estimates transmission bottleneck size using the approximate beta-binomial method of Sobel Leonard et al. (2017), which does not model sequencing depth explicitly. For the exact beta-binomial method (which accounts for finite read depth) or other approaches (presence-absence, KL divergence, Wright-Fisher), use the ViralBottleneck package (Zheng et al. 2025, Virus Evolution 11:veaf071) via asViralBottleneckInput.

Usage

quickBottleneck(whe, pairId, maxNb = 200L, threshold = 0.03)

Arguments

Value

A named list: Nb (MLE), ci (95% CI), loglik, n_variants, pair_id.

References

Sobel Leonard A et al. (2017). Transmission Bottleneck Size Estimation from Pathogen Deep-Sequencing Data. J Virol 91:e00171-17. doi:10.1128/JVI.00171-17

Examples

vcf1 <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
vcf2 <- system.file("extdata", "test_recipient.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(c(vcf1, vcf2),
    colData = S4Vectors::DataFrame(
        sample_id = c("donor", "recipient"),
        role = c("donor", "recipient"),
        pair_id = c("p1", "p1")),
    caller = "ivar")
quickBottleneck(whe, pairId = "p1", maxNb = 50L)

Description

Reads one or more VCF files produced by variant callers (LoFreq, Freebayes, or generic) and constructs a WithinHostExperiment.

Usage

readWithinHost(vcfFiles, colData, caller = "auto", genome = NULL)

Arguments

Details

Each VCF file corresponds to one sample. The colData argument provides sample metadata and must contain a sample_id column whose length matches the number of files.

caller = "ivar" refers to a VCF carrying iVar's per-sample fields (REF_DP, ALT_DP, ALT_FREQ, ALT_QUAL). ivar variants itself writes a TSV, not a VCF, so such files come from a converter, for instance ivar_variants_to_vcf.py in the nf-core/viralrecon pipeline. An iVar TSV is read directly by readWithinHostTable(format = "ivar").

Records with several alternative alleles are split into one site per allele. Per-allele values (INFO AF; FORMAT AO, AF, ALT_DP, ALT_FREQ and ALT_QUAL) are taken for that allele, and per-site values such as total depth and reference counts are shared.

When technical replicates are not available, Roder et al. (2023, mBio 14:e01046-23) recommend using a combination of multiple variant callers with stringent cutoffs to reduce false positives. To implement multi-caller consensus, import VCFs from each caller separately and intersect the resulting objects using standard GRanges overlap operations.

Value

A WithinHostExperiment.

Examples

vcf <- system.file("extdata", "test_donor.vcf",
    package = "WithinHostExperiment")
whe <- readWithinHost(vcf,
    colData = S4Vectors::DataFrame(sample_id = "donor"),
    caller = "ivar")
whe

Description

Reads pre-computed iSNV data from TSV or CSV files and constructs a WithinHostExperiment.

Usage

readWithinHostTable(files, format = c("ivar", "generic"), colData = NULL)

Arguments

Value

A WithinHostExperiment.

Examples

tsv <- system.file("extdata", "test_donor_ivar.tsv",
    package = "WithinHostExperiment")
whe <- readWithinHostTable(tsv, format = "ivar")
whe

Description

Computes a strand bias score that normalises the per-site forward/reverse alternative-allele ratio by the global strand coverage ratio, making it suitable for amplicon data where strand coverage is inherently unequal. Mostefai et al. (2024) showed that strand bias filters which ignore this imbalance can discard genuine variants; their own metric is a per-strand binomial likelihood.

Usage

sciStrandBias(fwd_alt, rev_alt, fwd_total, rev_total)

Arguments

Details

The SCI is defined as:

$\mathrm{SCI} = \log_2 \!\left( \frac{f_{\mathrm{alt}} / r_{\mathrm{alt}}} {F_{\mathrm{total}} / R_{\mathrm{total}}} \right)$

where $F_{\mathrm{total}}$ and $R_{\mathrm{total}}$ are the sums of forward and reverse total counts across all sites. A pseudocount of 1 is added to strand counts to avoid division by zero.

Value

A numeric vector of SCI scores (log2 scale). Values near zero indicate no strand bias; large absolute values indicate strand-specific enrichment or depletion.

References

Mostefai F et al. (2024). Refining SARS-CoV-2 intra-host variation by leveraging large-scale sequencing data. NAR Genomics and Bioinformatics 6:lqae145. doi:10.1093/nargab/lqae145

Examples

sci <- sciStrandBias(
    fwd_alt   = c(50, 2, 30),
    rev_alt   = c(45, 40, 28),
    fwd_total = c(500, 500, 500),
    rev_total = c(480, 480, 480)
)
sci

Description

Calculates per-site biallelic Shannon entropy, summed across all iSNV sites:

$H = -\sum_i [p_i \log(p_i) + (1-p_i) \log(1-p_i)]$

where $p_i$ is the alternative allele frequency at site $i$. This formulation treats each iSNV as an independent biallelic locus, following McCrone & Lauring (2016) and Popa et al. (2020). Higher values indicate greater within-host diversity.

Usage

shannonISNV(freq)

Arguments

Value

Numeric scalar. Shannon entropy.

References

Lauring AS (2020). Within-Host Viral Diversity: A Window into Viral Evolution. Annu Rev Virol 7:63-81. doi:10.1146/annurev-virology-010320-061642

Examples

shannonISNV(c(0.05, 0.15, 0.30, 0.45))
shannonISNV(0.5) # log(2)

Description

Calculates Simpson diversity summed over biallelic iSNV sites:

$D = \sum_i \left[1 - p_i^2 - (1 - p_i)^2\right] = \sum_i 2 p_i (1 - p_i)$

where $p_i$ is the alternative allele frequency at site $i$. Each site contributes the probability that two alleles drawn at random at that site differ (at most 0.5), so the index is the expected number of differing sites. Higher values indicate greater diversity.

Usage

simpsonISNV(freq)

Arguments

Value

Numeric scalar. Simpson diversity index.

Examples

simpsonISNV(c(0.05, 0.15, 0.30, 0.45))
simpsonISNV(0.5) # maximum per-site diversity = 0.5

Description

Computes Tajima's D (Tajima 1989) from summary statistics of a within-host pathogen population. The statistic measures the difference between nucleotide diversity ($\hat{\pi}$) and Watterson's estimator ($\hat{\theta}_W = S / a_1$), normalised by its expected standard deviation under neutrality:

$D = \frac{\hat{\pi} - S/a_1}{\sqrt{e_1 S + e_2 S(S-1)}}$

A significantly negative D indicates an excess of rare variants (e.g., population expansion or purifying selection), while a positive D suggests balancing selection or population contraction.

Usage

tajimaD(S, n, pi_hat, genomeLength)

Arguments

Details

Deep sequencing caveat: The classical Tajima (1989) formula was designed for independent haplotype samples (e.g., Sanger-sequenced clones or individual genomes). In deep sequencing, reads are not independent lineages; they are technical replicates of a viral population whose effective size is much smaller than the read depth. Using raw read depth (e.g., n = 2000) inflates the harmonic number $a_1 = \sum 1/i$, making $\theta_W$ artificially small and biasing D toward positive values.

We recommend capping n at a biologically meaningful effective sample size. Lauring (2020 Annu Rev Virol) suggests n = 100 as a conservative proxy. The high-level calcNeutralityTests method applies this cap automatically (n = min(depth, 100)).

Interpretation warning: Deep-sequencing D values are not directly comparable to classical population-genetic D computed from independently sampled haplotypes. We recommend reporting D as a descriptive statistic of the site frequency spectrum (e.g., "the majority of samples show negative D, consistent with an excess of rare variants"), rather than as a formal neutrality test with p-value interpretation. For rigorous neutrality testing on deep sequencing data, consider simulation-based approaches that account for the read-sampling process.

Value

Numeric scalar. Tajima's D.

References

Tajima F (1989). Statistical method for testing the neutral mutation hypothesis by DNA polymorphism. Genetics 123:585-595.

Lauring AS (2020). Within-Host Viral Diversity: A Window into Viral Evolution. Annu Rev Virol 7:63-81. doi:10.1146/annurev-virology-010320-061642

Examples

tajimaD(S = 10, n = 1000, pi_hat = 0.001, genomeLength = 10000)

Description

Extracts per-variant allele frequency trajectories across longitudinal time points within each host. Returns a long-format DataFrame suitable for downstream analysis or plotting.

Usage

trackFrequency(
  whe,
  hostCol = "host_id",
  timeCol = "timepoint",
  usePassedOnly = TRUE
)

Arguments

Value

A DataFrame with columns:

host_id

Host identifier.

variant_key

Unique variant identifier (chrom:pos:alt).

timepoint

Time point value.

frequency

Alternative allele frequency at that time point.

qc_passed

Logical. Whether the variant passed QC at that time point (TRUE if no QC assay present).

Examples

gr <- GenomicRanges::GRanges("seg1",
    IRanges::IRanges(c(100, 200), width = 1))
S4Vectors::mcols(gr)$ref <- c("A", "C")
S4Vectors::mcols(gr)$alt <- c("T", "G")
freq <- matrix(c(0.05, 0.10, 0.15, 0.20), nrow = 2, ncol = 2)
whe <- WithinHostExperiment(
    assays = list(altFreq = freq),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(
        sample_id = c("S1_t1", "S1_t2"),
        host_id = c("H1", "H1"),
        timepoint = c(1, 2)))
trackFrequency(whe)

Description

Compares observed diversity metrics against wide plausibility ranges for common pathogens, as a sanity check on units and orders of magnitude rather than a reference interval. The bounds are the package's own heuristics, not estimates from a specific study, and are deliberately generous: $\pi$ is compared on a per-site scale and Shannon entropy on the scale returned by calcDiversity. A value inside the range is not evidence that an estimate is correct.

Usage

validateDiversity(
  div_df,
  pathogen = c("sars_cov_2", "influenza", "hiv", "tb", "generic")
)

Arguments

Value

A DataFrame with columns:

metric

Diversity metric name.

observed

Observed value.

expected_min

Lower bound of the plausibility range.

expected_max

Upper bound of the plausibility range.

in_range

Logical; whether the observed value falls within the range.

Examples

div <- S4Vectors::DataFrame(pi = 1e-4, shannon = 5.0)
validateDiversity(div, pathogen = "sars_cov_2")

Description

Estimates Watterson's theta (Watterson 1975):

$\hat{\theta}_W = S / (a_n \cdot L)$

where $S$ is the number of segregating sites, $a_n = \sum_{i=1}^{n-1} 1/i$ is the harmonic number, and $L$ is the genome length, with $n$ the number of sampled viral genomes.

Usage

wattersonISNV(nSites, genomeLength, meanDepth)

Arguments

Details

As tajimaD sets out, raw read depth should not be passed as $n$ for deep sequencing data: reads are technical replicates of a much smaller viral population, and a depth-sized $n$ inflates $a_n$ and makes $\theta_W$ correspondingly small. Pass an effective sample size, or a depth capped at one, as calcNeutralityTests does.

Value

Numeric scalar. Watterson's theta estimate.

References

Watterson GA (1975). On the number of segregating sites in genetical models without recombination. Theor Popul Biol 7:256-276. doi:10.1016/0040-5809(75)90020-9

Examples

## n is an effective sample size, not the raw depth
wattersonISNV(15, genomeLength = 13588, meanDepth = 100)

Description

Methods for the combineRows and combineCols generics. combineRows combines objects vertically (adding variant sites); all objects must have the same samples (columns). combineCols combines objects horizontally (adding samples); all objects must have the same variant sites (rows). QC logs and filters from all objects are merged.

Usage

## S4 method for signature 'WithinHostExperiment'
combineRows(x, ...)

## S4 method for signature 'WithinHostExperiment'
combineCols(x, ..., use.names = TRUE)

Arguments

Value

A WithinHostExperiment.

Examples

data(example_whe)
whe1 <- example_whe[1:10, ]
whe2 <- example_whe[11:nrow(example_whe), ]
merged <- combineRows(whe1, whe2)
nrow(merged) == nrow(example_whe)

whe_d <- example_whe[, 1]
whe_r <- example_whe[, 2]
merged <- combineCols(whe_d, whe_r)
ncol(merged)

Description

Constructs a WithinHostExperiment from assay matrices, genomic ranges, and sample metadata.

Usage

WithinHostExperiment(assays, rowRanges, colData, metadata = list())

Arguments

Value

A WithinHostExperiment object.

Examples

library(GenomicRanges)
gr <- GRanges("seg1", IRanges::IRanges(c(100, 200, 300), width = 1))
whe <- WithinHostExperiment(
    assays = list(altFreq = matrix(c(0.1, 0.2, 0.3), ncol = 1)),
    rowRanges = gr,
    colData = S4Vectors::DataFrame(sample_id = "S1"))
whe