Corrections to functions that the earlier reviews had not examined. Each has a test that fails on 0.99.1. Several change results; each entry says what the function did and what it does now.
detectEmergingVariants() returns one more column, qc_from: whether
a call at the earlier timepoint of the interval passed QC.combineRows() and combineCols() are an error when the objects do
not hold the same samples, or the same sites (see Corrections).trackFrequency(), detectEmergingVariants() and
plotFrequencyTrajectory() stop on a timepoint column of text that has
no time order, such as "day10", as flagTemporalInconsistency()
already did, and leave out samples without a timepoint, with a warning._pkgdown.yml no longer give
https://cuiweig.github.io/WithinHostExperiment/, which does not
resolve; the repository address remains.qcPass of NA counts as not passed wherever an analysis uses only
passing calls. Ten masks wrote freq[!qc] <- NA, which R skips at an
NA index, so such a call was used as passing by calcDiversity(),
compareDiversity(), buildSFS(), calcNeutralityTests(),
pNpSWithinHost(), exportPairFrequencies(), calcSharedVariants(),
the bottleneck inputs, detectCrossContamination(),
trackFrequency(), detectEmergingVariants(),
plotFrequencyTrajectory() and plotFrequencySpectrum(). No package
function writes NA into qcPass; one set by hand now counts as
failed, as it already did in flagTemporalInconsistency() and the
consensus.benchmarkCallers() and calcCallerConcordance() counted it as a
called site, so a caller with two such rows scored precision 0.6
instead of 1. qcSummary() counted it in total_sites, as flagged
once flagISNV() had run, and in median_freq: four calls among ten
assessed sites gave a pass rate of 0.4 and a median frequency of 0.
Every column of qcSummary() now describes the calls, median_depth
included, and the help page defines each one.
plotFrequencySpectrum() and panels a and d of plotQCDashboard()
drew and counted it as an iSNV; panel a now leaves out a sample with
no call, so a sample of recorded zeros does not stop the plot.plotFrequencySpectrum(),
plotQCDashboard() and plotPairScatter() are rounded, where
as.integer(x * 100) truncated them: a break at 4% could read 3, and
one at 2.5% read 2.trackFrequency(), detectEmergingVariants() and
plotFrequencyTrajectory() read the timepoint column in time order,
as flagTemporalInconsistency() does. Numbers given as text were
sorted as text, so "1", "2", "10" gave an interval from 10 to 2, and
an NA timepoint was sorted last and formed an interval ending in
NA, with no warning.trackFrequency() reports each call's QC state in qc_passed. It read
qcPass only with usePassedOnly = TRUE, which then removed the
failing calls, so the column was TRUE on every row.detectEmergingVariants() combines replicate samples taken at one
timepoint by the mean of their calls before it forms intervals. Two
replicates formed an interval from a timepoint to itself, a gap
between them was reported as emergence, and the answer depended on the
order of the columns.detectEmergingVariants() forms intervals on every call and applies
QC to their end points. It removed the calls that failed QC first, so
after the default flagISNV() a start below minFreq (0.03), a
recorded 0 included, could never be reported, and a timepoint that
failed QC was skipped, joining the timepoints either side of it. The
later point must now have a call that passed QC; the earlier point
uses its passing calls, or all its calls when none passed, since a
frequency below minFreq is the low start that emergence is measured
from, and qc_from says which.detectCrossContamination() keeps a source call that QC failed for
its frequency alone. After the default flagISNV(), whose maxFreq is
0.5, every source above 50% was removed, the 80% source of the help
page's own example included. A failed source call above the strictest
maxFreq of the filters applied is now used when it has the depth,
alternative reads and strand balance that each of them asks for.correctSFSBias() rounds the corrected spectrum as a total and shares
its sites out among the bins by largest remainder. Rounding each bin
left every bin of 7 or fewer iSNVs unchanged in a folded spectrum at
threshold 0.03, and every bin of 16 or fewer in an unfolded one, so the
truncation correction of the help page's example spectrum kept its 15
sites while corrected became TRUE; it now gives 16, and nSites is
the corrected total. The binomial method, which the example runs, still
gives 15: at its median depth of 1,600 reads every bin is detectable.combineRows() matches samples by column name, and combineCols()
matches sites by sequence name, position, ref and alt. Both bound
by position, so an object holding the same samples or sites in another
order had its values written to the wrong ones without a message, and
one holding others was relabelled. The same samples or sites in another
order are now put in the first object's order; others are an error
that names the argument and what differs.exactBottleneck(method = "exact_bb") derives a missing recipient
alternative count at a known depth as round(altFreq * totalDepth), as
readWithinHostTable(format = "generic") does. It took the count as 0,
so an iVar table without its ALT_DP column made every recipient call
an absence: on the shipped pair without the recipient's counts it gave
Nb 5 [3, 7] where the counts give 95 [21, 272]. A recipient with no
frequency there still counts as 0 reads.qcSummary(),
correctSFSBias(), trackFrequency() and detectEmergingVariants();
the truncation correction of its donor spectrum gives 15 variants
rather than 14.Changes made in response to the Bioconductor review. The package has not
been released, so functions are removed and renamed without a deprecation
cycle, and objects saved with 0.99.0 must be rebuilt: an ISNVFilter
saved then has a logical replicateConc, and a WithinHostSFS lacks the
polarization and maxFreq slots. The WithinHostExperiment validity
checks every filter in qcFilters, so validObject(), subsetting and
combining report a filter saved with 0.99.0, with a message to rebuild it
with ISNVFilter(), and flagISNV() given one names its replicateConc
slot.
annotateFromGFF() is removed, together with the package's own GFF3
and FASTA readers. annotateGenes(x, genes) labels each site with the
gene it falls in, adding a GENE column to mcols(rowRanges(x)).
genes is a GRanges, such as a GFF3 or GTF file read with
rtracklayer::import(), or a TxDb, from which
GenomicFeatures::genes() takes the ranges; the overlap is
GenomicRanges::findOverlaps(). nameCol selects the label column, and
an unknown one is an error that lists those available; collapse keeps
the longest overlapping gene or all of them. When the sites and the
genes share no sequence name the function stops and prints both sets,
so a RefSeq/GenBank naming mismatch can no longer return a column of
NA.annotateCodonChange() is an S4 generic, annotateCodonChange(x, cds, ref), taking the coding sequences as a GRanges of CDS records, a
GRangesList or a TxDb, and the reference as a DNAStringSet; a file
path for either is refused with a message naming the Bioconductor reader
to use. The consequence is computed by
VariantAnnotation::predictCoding(), and the output columns take its
names: TXNAME, GENEID, CONSEQUENCE, CDSLOC, PROTEINLOC,
REFCODON, VARCODON, REFAA and VARAA replace GFF_FEATURE,
REF_CODON, ALT_CODON, REF_AA, ALT_AA and AA_CLASS. A site in
several overlapping coding sequences is resolved by collapse, the
longest transcript by default; a nucleotide that one transcript reads
twice, at a -1 ribosomal frameshift, keeps its first reading, and a
message says so. Three classifications change with the delegation:
losing a stop codon is nonsynonymous rather than nonsense, an
alternative start codon (CTG or TTG) at the first base of a coding
sequence reads as M, and a site whose reference codon holds an
ambiguous base is left unannotated, with a warning. The reference
length is recorded in seqlengths(rowRanges(x)). The GRanges and
GRangesList routes build the TxDb with txdbmaker; a TxDb built
with txdbmaker::makeTxDbFromGFF() can lack transcripts the GFF3 holds
(for NC_045512.2 it drops the ORF1ab gene, so sites in ORF1a and ORF1ab
are left unannotated), so the man page advises passing the CDS records
instead.dndsWithinHost() is replaced by pNpSWithinHost(x, cds, ref), which
takes the same inputs as annotateCodonChange() and calls it for the
reference and variant codon of each iSNV. The quantity is a ratio of
nonsynonymous to synonymous polymorphism within a host, not of
substitutions between sequences, so the Jukes-Cantor correction is no
longer applied. Nei-Gojobori sites are counted over the distinct codons
of the genome, so codons shared by same-frame overlapping coding
sequences, such as ORF1a and ORF1ab, are counted once, and a warning
names the overlapping transcripts. Per gene, a shared codon counts in
each gene that reads it, and so does an iSNV in it, each by that gene's
own transcript. On the README's SARS-CoV-2 call set the pooled ratio
moves from 0.2857 to 0.2867 with the correction removed. The man page
says why the sites are counted here rather than by MSA2dist:
MSA2dist::codonmat2pnps() and dnastring2kaks() compare aligned
sequence pairs and cannot take per-site calls; for each of the 61 sense
codons the sites equal those of codonmat2pnps(), which a test checks.consensusSequences(x, ref, minFreq = 0.5) is new. It writes each
sample's substitutions at or above minFreq into the reference and
returns one sequence per sample as a DNAStringSet; insertions and
deletions are not applied, so the sequences stay aligned to the
reference and to one another. region returns one coding sequence,
joined across its parts and read in its own direction. dN/dS between
hosts can then be estimated with MSA2dist::dnastring2kaks(), which is
suggested, rather than with a model of this package's own. When ref
and the sites share no sequence name it stops and names both sets, and
a region on a sequence missing from ref is reported as such.readWithinHost() keeps its arguments but reads files with
VariantAnnotation::readVcf() instead of parsing VCF text, and splits
records with several alternative alleles with expand(). A per-allele
field that the header does not declare Number=A on such a record stops
the import, since readVcf() would give every allele the first allele's
value; a field declared Number=R, whose first value is the reference
allele's, stops it on any file, where that value was read as the
alternative allele's count or frequency. genome takes a
DNAStringSet, not a file path; its lengths, or those of ##contig
header lines, are recorded in seqlengths() for every sequence named,
with a call or without. A genome that names none of the sequences the
calls lie on, or is shorter than a called position, is an error, since
every per-site statistic would divide by its length. A file holding more
than one sample is refused; as(vcf, "WithinHostExperiment") converts a
VCF object with all of its samples. Coercion to and from VRanges
keeps every recorded sequence and its length, so the genome length
survives a round trip, and a sequence name may contain ':', which it
used to be split on.calcDiversity() and compareDiversity() lose the genomeLength
argument. The length that pi and watterson divide by is read from
seqlengths(rowRanges(x)), summed over every sequence recorded for
the reference, so that a segment of a segmented genome without calls
still counts, and a length that is not recorded is an error saying
how to record it. calcDiversity() gains n_sites_assessed, the sites
at which a sample has any call before QC, and detection_limit, the
strictest minFreq of the filters applied, so that a diversity of zero
after QC can be told apart from a sample with no calls. example_whe
records its 1000-base segment length; its values are unchanged.calcNeutralityTests() loses the genomeLength argument it required.
Tajima's D and Fu's Fs use the nucleotide diversity summed over the
genome, which the frequencies give directly, so the length was divided
in and multiplied back out and changed no output beyond rounding
error; passing it is now an error.validateDiversity() is renamed checkDiversityPlausibility(). Its
ranges are the package's own heuristics, not literature values, and the
man page lists every one. Every sample of a calcDiversity() table is
checked, where only the first row was checked before, and the result
names the sample; a missing value is reported as not in range.ISNVFilter() loses the minBaseQual and minMapQual slots and
arguments: the object holds no per-site quality, so they could be
recorded but never applied, and quality filtering belongs in the
variant caller.replicateConc holds a ConcordanceRule, or NULL for none, where
it held a logical that was recorded but never applied.
ConcordanceRule(freqTolerance, requireBothDetected) takes the two
settings of flagReplicateDiscordance(), and flagISNV() applies it
after the per-sample thresholds, flagging what that function flags
and recording the same two qcLog rows, so the rule is stored with the
other thresholds in qcFilters. A rule with no replicate_group of
two or more samples to compare is an error; samples without a
replicate are named in a warning. A logical replicateConc is
refused, and so is a rule that would flag nothing: a freqTolerance of
1 or more with requireBothDetected = FALSE, since frequencies lie in
[0, 1]. ConcordanceRule() refuses a factor, string or logical
tolerance and a requireBothDetected that is not TRUE or FALSE,
rather than convert them, and flagReplicateDiscordance() now checks
its arguments in the same way. flagISNV() checks its filter with
validObject() first, so a threshold or rule changed with @<- into a
value its class refuses is an error rather than a step applied.consensusISNV() keeps the seqinfo() of its inputs, merged, and
stops when a sequence's length, circularity or genome differs between
them. It takes each sample's colData row from the first input that
holds the sample, and keeps a column for every input sample, NA where
no site passed, with a warning naming any sample that fewer than
min_callers inputs hold, since such a sample cannot reach consensus.
Sequence lengths, host, role and pair columns used to be dropped, so
calcDiversity() stopped for want of a genome length and the
transmission functions found no pairs.quickBottleneck() and exactBottleneck() gain
n_excluded, the donor variants at or above threshold left out of
the estimate, by reason (see below).WithinHostSFS gains a polarization slot: "none" for a folded
spectrum, "reference" for an unfolded spectrum of
alternative-allele frequencies, which treats the reference allele as
ancestral, or "outgroup". Validity ties it to folded, and show()
prints it and labels the stored depth a median, which it is.neutralityFromSFS() stops on a folded spectrum, where it printed a
note and returned a value, because Fay and Wu's H needs derived-allele
frequencies; its result records the polarisation it used. It also
stops on a spectrum without a genome length, where it returned values
summed over the genome under the label "per site". buildSFS() now
takes that length from seqlengths(rowRanges(whe)) when it is not
given, as calcDiversity() does.WithinHostSFS gains a maxFreq slot, the highest frequency a call
could keep through QC. buildSFS() sets it to the smallest maxFreq
among the filters applied with flagISNV() when it uses QC-passed
calls, and to 1 otherwise. neutralityFromSFS() stops on an unfolded
spectrum whose cap is below 1: the default ISNVFilter() removes every
call above 0.5, which empties the classes Fay and Wu's H weights most,
so H could not reveal an excess of high-frequency variants whatever
the data. The spectrum must be built from the unfiltered object
flagged with maxFreq = 1. show() prints the cap when it is below 1.neutralityFromSFS() computes pi and theta_H on the same derived-allele
classes, round(mid * n) held within 1 to n - 1, with Fay and Wu's
normalisation. Before, pi came from unrounded midpoints without the
n/(n-1) factor, which biased H by about -theta/n under neutrality even
on a complete spectrum. H is now centred on 0 only when every class
from 1 to n - 1 is observed and each bin holds one class; the detection
threshold and wide bins still shift it, to about -0.03 theta on the
expected neutral spectrum under the defaults (20 bins, threshold 0.03,
n = 100), much as before. pi, Tajima's D, Fu's Fs and H from this
function change.compareSFS() stops unless both spectra share their bin breaks, their
folding, their polarisation and their QC frequency cap, and neither is
bias-corrected; it used to compare any two spectra with the same number
of bins. A spectrum from before QC can therefore be compared with one
from after it only when the filters kept maxFreq = 1. When either
spectrum is empty it runs no test and returns NA diagnostics, where
it returned NaN on a positive df.n_sites, S and theta_pi in
calcNeutralityTests(), to richness, watterson and chao1 in
calcDiversity(), and to the sites buildSFS() bins. A frequency of 0,
or of 1, a fixed difference that iVar reports, used to be counted
without adding to pi, which pulled Tajima's D down.tajimaD() returns NA, with a warning, below three segregating
sites, where it returned 0 at none and a number at one or two, and
when the sample size rounds below four, where the variance of D is
zero and it returned Inf or NaN. fusFs() returns NA, with a
warning, when there is no segregating site or theta_pi is 0, where it
returned 0, a number where none is defined; above n = 500 it computes
the Ewens tail exactly, where a normal approximation made the value
jump at the switch. calcNeutralityTests() reports n_sites, the
number of sites each statistic rests on, documents every rule under
which D or Fs is NA and warns once per cause. A sample with no calls
gives zero diversity and zero segregating sites rather than an
uninformative error, in calcNeutralityTests() and in
calcDiversity().fusFs(), calcNeutralityTests() and
neutralityFromSFS() say that the haplotype count is set to S + 1, the
infinite-sites maximum. The value therefore carries no haplotype
information and is pushed below 0 under neutrality, increasingly so as
theta grows (-7.00 at theta = 5 with n = 100, -0.37 at theta = 1), so a
negative Fs is not evidence of expansion.pNpSWithinHost() counts an iSNV only through a reading whose codon
carries sites. A reading in a stop codon, including the loss of a stop
codon, is not counted, since a stop codon carries no sites; where coding
sequences overlap, such an iSNV is counted through a transcript that
reads it in a sense codon, and one that every transcript reads in a stop
codon enters neither n_S nor n_N. Each call is
classified by translating its reference and variant codons with the
standard genetic code and no initiator rule, so ATG to CTG or TTG
at the first base of a coding sequence is nonsynonymous here, as ATG
carries no synonymous site, although annotateCodonChange() reports it
as synonymous. The rule applies to the totals and to each gene's
counts. The vignette donor's ratio moves from 1.11 (3 synonymous, 11
nonsynonymous) to 1.36 (2 and 9); the README figure values do not
change.calcSharedVariants() gains minDepth, by default the strictest
minDepth applied with flagISNV() (100 when none was), and splits the
donor-only count into n_donor_only_adequate_depth,
n_donor_only_low_depth and n_donor_only_unknown_depth, with
n_donor_only_recipient_failed_qc for the sites where the recipient
made a call that QC removed, which are not absences; only the first is
evidence that a variant was lost in transmission. n_donor_only is
still their sum.exportPairFrequencies() and calcSharedVariants() take a host to
carry an allele only where its frequency is known and above 0. A site
that a table read with readWithinHostTable(format = "generic") lists
with aaf = 0 counted as shared when no QC had been applied; it is now
an absence, split by the recipient's depth. recipient_freq keeps the
recorded 0, and sites that neither host carries are no longer listed;
a pair with no such site gives an empty table.asViralBottleneckInput(), quickBottleneck() and every method of
exactBottleneck() leave out, with a message, a donor variant at a site
where QC removed a recipient call at or above threshold. Such a call,
for example one above maxFreq where the variant rose towards fixation,
was written as a frequency of 0 and scored as a variant lost in
transmission, while exact_bb used its raw read counts. A removed call
below threshold is an absence by the estimators' own definition and
is still written as 0, and a recorded 0 is never taken for a call.exactBottleneck(method = "exact_bb") leaves out, with a message, donor
variants whose recipient depth is unknown, which includes every site
the recipient's VCF file does not list. They added the same constant to
the log-likelihood at every Nb, so the estimate and the
likelihood-ratio interval do not change, but n_variants counted them
as used and a bootstrap interval resampled them; they are now counted
in n_excluded, and the manual page says that exact_bb then sees no
loss in transmission.quickBottleneck() and the
presence_absence and wright_fisher methods of exactBottleneck(),
use donor variants from threshold to 1 - threshold, as the
published methods do. A donor variant at frequency 1 that the recipient
lacked gave a log-likelihood of -Inf at every Nb and an estimate of 1
without a warning. Both functions now stop when the log-likelihood is
not finite at any Nb of the grid.quickBottleneck() checks maxNb, where maxNb = 0 returned an
estimate of length 0 with two base-R warnings. It and
exactBottleneck() refuse a maxNb that is not a single whole number
from 1 to the integer limit; exactBottleneck() truncated 2.5 to 2 and
accepted text and logical values. Both refuse a threshold that is not
a single number above 0 and below 0.5: at 0 every recipient frequency
counts as present, and from 0.5 no donor variant lies between
threshold and 1 - threshold. quickBottleneck() did not check it,
and exactBottleneck() accepted any value from 0 to 1.flagReplicateDiscordance(), and the replicate step of flagISNV(),
let a frequency difference equal to freqTolerance pass, as
documented. A decimal subtraction such as 0.20 - 0.15 rounds to just
above 0.05, so 39 of the 96 two-decimal pairs 0.05 apart were flagged;
the comparison now allows for rounding error, and so does the frequency
spread of consensusISNV(). The vignette's counts are unchanged. The
replicate step also counts the cells it flags when qcPass holds
NA, where the qcLog recorded n_flagged as NA.flagTemporalInconsistency() counts the distinct timepoints at which a
site is detected within a host, not the samples, so technical replicates
taken at one timepoint no longer make a variant persistent;
minConsecutive runs over those timepoints, and a host sampled at a
single timepoint is not classified. Timepoints may be numbers, dates,
date-times, time differences, a factor in time order, or text that reads
as numbers; other text is an error, and a blank is a missing timepoint.
A qcPass cell that is NA counts as not passed, where an NA cell
with a frequency stopped the function with an uninformative error, and
is left NA when the site is flagged. Samples with no timepoint are not
counted, with a warning.flagStrandBiasSCI() computes each sample's global strand ratio over
the sites that sample called. The ratio was summed without removing
missing values, so in a multi-sample object any sample with an absent
call scored NA throughout and nothing in it was flagged.
sciStrandBias() scores NA at a site with a missing count and leaves
it out of the sums, and it refuses negative, infinite and non-numeric
counts, where a negative count returned -Inf.buildSFS() checks sampleIdx, nBins, threshold and fold and
names the one at fault; two samples are no longer pooled, and a filter
with a maxFreq of 0, which leaves no call to count, is an error giving
that value. It folds a sample with no passing
call to an empty spectrum, where it failed and took buildSFSList()
down with it. piISNV(), shannonISNV() and simpsonISNV() stop on
text and on frequencies outside [0, 1], where piISNV() returned a
negative diversity.flagPrimerSites() stops, naming both sets, when the primers and the
sites share no sequence name, where it reported that no site overlapped
a primer and flagged nothing.exactBottleneck(method = "wright_fisher") says in a message that its
result is a heuristic score, not a likelihood; exact_bb and
presence_absence are the estimators.readWithinHostTable() reads the sequence and allele columns as text.
read.delim() and read.csv() turned a column that held only "T" into
the logical TRUE, so an iVar table whose calls were all C>T was
imported with the allele TRUE.consensusISNV() counts each caller once per site and sample: a caller
that listed the same site and allele twice supplied two of the
min_callers votes and was named twice in caller_names. It no longer
splits a sequence name at its colons, which put every site on a contig
such as HLA-A*01:01 at one wrong position.detectCrossContamination() stops, listing the available columns, when
runCol is not a single column of colData; a misspelt name used to
compare every sample with every other without a word. Samples whose
run is NA are left out with a warning giving their number.plotPairScatter() draws a donor variant whose recipient frequency is
a recorded 0 as donor-only, where it drew it as shared, and says there
is nothing to plot for a pair in which neither host carries an allele.consensusISNV() builds one table of calls instead of looping over
every cell: 3,000 sites in 20 samples from two callers took about four
and a half minutes and take half a second, with the same result.detectEmergingVariants() orders the trajectories once and compares
each timepoint with the next in one step, giving the same rows in the
same order: 3,000 sites in 20 samples took about a minute and take a
fraction of a second. Samples with no host are left out with a warning,
where they stopped the function.LongRead is removed from biocViews.qcLog and qcFilters slots that
WithinHostExperiment adds to RangedSummarizedExperiment are
described one by one, and installation is from Bioconductor only.readWithinHost() and readWithinHostTable() list
what the import does not read. The FILTER column of a VCF and iVar's
PASS column are not read, so records the caller rejected are imported
as calls; per-strand counts are not read, so the strand-bias step of
flagISNV() does not run on an imported object; and insertions and
deletions are imported as width-1 rows in the caller's notation,
without being told apart from substitutions. Each page shows how to
drop rejected records, or keep substitutions only, before QC.WithinHostExperiment class page states that assays must be
ordinary R matrices, and a qcPass assay held as a DelayedArray or
another backend is refused with a message saying so.inst/extdata.inst/scripts/make_extdata.R, which ships, records how each file in
inst/extdata was made and with which seed.readme_figure_output/ is removed from the repository root. The values
and checksums the README is checked against are kept in
inst/scripts/readme_values/, and inst/scripts/reproduce_all.R
rebuilds every derived file from the raw data and checks it against
the committed copy.tests/testthat/test-corrections.R holds a regression test for each
item of 'Corrections made before submission' under 0.99.0 that still
applies, named after the bug it guards against, and the README's
statement that QC cannot raise pi as a test. The removed codon
annotator is kept in tests/testthat/helper-reference-codon.R as a
reference implementation, and annotateCodonChange() is tested
against it.flagISNV() is tested at a tolerance other than the default, with
requireBothDetected = FALSE, on calls that fail a threshold, at the
tolerance boundary from both sides and on a group of three replicates,
with hard-coded expected flags.Before submission every exported function was read against its
documentation and run on fixtures built to break it, and each
discrepancy found was fixed, documented as implemented, or both. The
items below record what changed. In 0.99.1 each item that still applies
was pinned by a regression test in tests/testthat/test-corrections.R,
named after the bug it guards against (for example exact_bb matches sites by chrom, pos and allele (was: pos only)), so that none can
return unnoticed.
compareSFS() suppressed every warning from chisq.test(), so a p-value
resting on expected counts below 5 looked as trustworthy as any other. It
now returns the smallest expected count and chisq_approx_ok, which is
FALSE in that case, and silences only that one warning; the statistic and
p-value are unchanged.flagTemporalInconsistency() classified each site with the last host
processed and flagged transient sites in every sample. Sites are now
classified within each host and flagged only in the samples of hosts where
they are transient.exactBottleneck(method = "exact_bb") matched donor variants to rows by
position alone, which is wrong for segmented genomes, and paired donor
frequencies with recipient read counts by index when some sites did not
match. Sites are now matched by chromosome, position and allele. An estimate
that reaches maxNb triggers a warning.exactBottleneck() states how exact_bb departs from
Sobel Leonard et al. (2017) and that wright_fisher is a heuristic score.flagISNV(), calcDiversity() and calcNeutralityTests() stop on
arguments they do not use.adaptiveFilter("ont") uses minDepth = 60, the depth above which Bull et
al. (2020) found over 99% sensitivity and precision for single-nucleotide
variants; it was 30.correctSFSBias(), compareSFS(), the WithinHostSFS bins and depth slot,
consensusISNV(), flagReplicateDiscordance() and validateDiversity() are
documented as implemented: truncation rescales the whole spectrum, bins are
right-closed, the stored depth is a median, the consensus object holds only
allele frequencies, and the diversity ranges are heuristic plausibility
bounds rather than literature values.ISNVFilter() documents minBaseQual, minMapQual and replicateConc as
recorded but not applied, and flagISNV() warns when they are set, instead
of ignoring them silently. WithinHostExperiment() requires an altFreq
assay, as documented.trackFrequency() builds its table in one vectorised step, which is much
faster on genome-wide data and gives the same rows.plotFrequencySpectrum(), plotPairScatter() and plotQCDashboard() clip the
view with coordinate limits instead of removing data outside the axis range,
which dropped jittered points with a warning; plotFrequencyTrajectory()
uses the Okabe-Ito palette, as described, and plotSFS() recycles it for
more than eight samples.WithinHostExperiment S4 class extending
RangedSummarizedExperiment for within-host pathogen iSNV data,
with setValidity enforcing structural invariants.readWithinHost() reads VCFs in the LoFreq,
Freebayes, generic and iVar-style field layouts with auto-detection
(ivar variants writes a TSV, so an iVar-style VCF comes from a
converter). readWithinHostTable() handles the iVar TSV directly,
and generic CSV.flagISNV() applies configurable quality
filters without removing data. flagReplicateDiscordance()
identifies discordant technical replicates with full provenance
logging. flagPrimerSites() masks amplicon primer regions.shannonISNV(), piISNV() (with optional
finite-sample correction per Nei 1987), wattersonISNV(),
simpsonISNV(), chao1ISNV(), and calcDiversity() for
per-sample within-host diversity.tajimaD(), fusFs(), and
calcNeutralityTests() for site-frequency-spectrum-based
neutrality testing with deep-sequencing-aware sample size
capping.WithinHostSFS S4 class stores the within-host SFS with
metadata (sample, genome length, depth, detection threshold).
buildSFS() / buildSFSList() construct from a WHE;
correctSFSBias() applies ascertainment bias correction for
the detection threshold (truncation or binomial method);
neutralityFromSFS() computes Tajima's D, Fu's Fs, and
Fay and Wu's H (new) from a single SFS object;
compareSFS() tests SFS homogeneity between samples;
plotSFS() provides publication-quality SFS visualisation.flagTemporalInconsistency() identifies transient iSNVs that
appear at only one timepoint (likely artefacts) versus persistent
ones detected across multiple consecutive timepoints. Supports
configurable minTimepoints and minConsecutive thresholds.detectCrossContamination() identifies index-hopping signatures
where a high-frequency variant in one sample appears at
suspiciously low frequency in a co-sequenced sample, consistent
with adapter exchange during pooled sequencing. Returns a
structured DataFrame of suspect events with source/sink
frequencies and run metadata.annotateFromGFF() overlaps variant sites with
gene models from a GFF3 file, adding GFF_FEATURE to
mcols(rowRanges) for use with dndsWithinHost().
annotateCodonChange() extends this with CDS-aware codon
translation: given a GFF3 and reference FASTA, it determines
the affected codon, translates both reference and alternative
alleles using the standard genetic code, and classifies each
SNV as Synonymous, Nonsynonymous, or Nonsense. (0.99.1 replaces
both with functions that take Bioconductor objects and delegates the
consequence call to VariantAnnotation::predictCoding().)dndsWithinHost() estimates within-host dN/dS
with the Nei-Gojobori counting method, normalising by the synonymous
and nonsynonymous sites of the CDS features in a GFF3 and reference
FASTA, with Jukes-Cantor correction; estimateSelectionCoefficient()
for time-series selection coefficient inference.consensusISNV() for deriving consensus
variant calls across callers, calcCallerConcordance() for
inter-caller agreement metrics, benchmarkCallers() for caller
performance comparison.trackFrequency() for tracking variant
allele frequencies over time, detectEmergingVariants() for
identifying newly arising variants, plotFrequencyTrajectory()
for publication-quality trajectory plots.transmissionPairs(),
exportPairFrequencies(), calcSharedVariants(), and
quickBottleneck() (approximate beta-binomial MLE, Sobel
Leonard et al. 2017) for bottleneck estimation.exactBottleneck() offers a
beta-binomial likelihood over the transmitted virions that models
finite recipient depth, in the spirit of Sobel Leonard et al. (2017)
with two stated departures, a presence-absence likelihood, and a
Wright-Fisher score, which is a heuristic rather than a likelihood.sciStrandBias() and
flagStrandBiasSCI() for strand bias filtering that accounts
for unequal strand coverage.adaptiveFilter() for technology-specific
QC presets (Illumina, ONT, PacBio HiFi).fitDiversityModel() for regression
modelling of diversity measures, multitestCorrection() for
multiple testing correction, permutationTest() for
non-parametric hypothesis testing.validateDiversity() compares diversity estimates with
wide heuristic plausibility ranges per pathogen, not with literature
values.asViralBottleneckInput() exports
standardised data for the ViralBottleneck package (Zheng
et al. 2025).as(whe, "VRanges"),
as(vr, "WithinHostExperiment"), and
as(whe, "data.frame") for seamless integration.combineRows() and
combineCols() generics merge experiments together with their QC
logs and filters.plotFrequencySpectrum(),
plotPairScatter(), plotQCDashboard(), and
plotFrequencyTrajectory() with a colour-blind-safe palette
(Wong 2011).