From f3cc08cf2da3dae9a80c1a6f42ec7407d5e196b4 Mon Sep 17 00:00:00 2001 From: vjcitn Date: Fri, 17 Jul 2026 21:15:29 -0400 Subject: [PATCH 1/3] fix Rd warnings: \itemize -> \describe and add non-empty labels Convert \itemize blocks that used \item{label}{desc} syntax to \describe, and replace empty \item{}{} labels in \describe blocks with the function signature from the body, resolving all checkRd warnings in isSNV-methods, PolyPhenDb-class, readVcf-methods, SIFTDb-class, summarizeVariants-methods, VariantType-class, VCF-class, VcfFile-class, VCFHeader-class, VRanges-class, and VRangesList-class Rd files. Co-Authored-By: Claude Sonnet 4.6 --- man/PolyPhenDb-class.Rd | 25 ++--- man/SIFTDb-class.Rd | 14 +-- man/VCF-class.Rd | 159 ++++++++++++++----------------- man/VCFHeader-class.Rd | 47 ++++----- man/VRanges-class.Rd | 93 +++++++++--------- man/VRangesList-class.Rd | 14 ++- man/VariantType-class.Rd | 21 ++-- man/VcfFile-class.Rd | 30 +++--- man/isSNV-methods.Rd | 12 +-- man/readVcf-methods.Rd | 14 +-- man/summarizeVariants-methods.Rd | 16 ++-- 11 files changed, 201 insertions(+), 244 deletions(-) diff --git a/man/PolyPhenDb-class.Rd b/man/PolyPhenDb-class.Rd index 28d92e4..f13e7cf 100644 --- a/man/PolyPhenDb-class.Rd +++ b/man/PolyPhenDb-class.Rd @@ -22,32 +22,27 @@ \section{Methods}{ In the code below, \code{x} is a \code{PolyPhenDb} object. \describe{ - \item{}{ - \code{metadata(x)}: + \item{\code{metadata(x)}}{ Returns \code{x}'s metadata in a data frame. } - \item{}{ - \code{columns(x)}: + \item{\code{columns(x)}}{ Returns the names of the \code{columns} that can be used to subset the data columns. For column descriptions see \code{?PolyPhenDbColumns}. } - \item{}{ - \code{keys(x)}: + \item{\code{keys(x)}}{ Returns the names of the \code{keys} that can be used to subset the data rows. The \code{keys} values are the rsid's. } - \item{}{ - \code{select(x, keys = NULL, columns = NULL, ...)}: - Returns a subset of data defined by the character vectors \code{keys} + \item{\code{select(x, keys = NULL, columns = NULL, ...)}}{ + Returns a subset of data defined by the character vectors \code{keys} and \code{columns}. If no \code{keys} are supplied, all rows are returned. If no \code{columns} are supplied, all columns are returned. See \code{?PolyPhenDbColumns} for column descriptions. - } - \item{}{ - \code{duplicateRSID(x)}: - Returns a named list of duplicate rsid groups. The names are the - \code{keys}, the list elements are the rsid's that have been - reported as having identical chromosome position and alleles and + } + \item{\code{duplicateRSID(x)}}{ + Returns a named list of duplicate rsid groups. The names are the + \code{keys}, the list elements are the rsid's that have been + reported as having identical chromosome position and alleles and therefore translating into the same amino acid residue substitution. } } diff --git a/man/SIFTDb-class.Rd b/man/SIFTDb-class.Rd index b695259..c2c7310 100644 --- a/man/SIFTDb-class.Rd +++ b/man/SIFTDb-class.Rd @@ -21,23 +21,19 @@ \section{Methods}{ In the code below, \code{x} is a \code{SIFTDb} object. \describe{ - \item{}{ - \code{metadata(x)}: + \item{\code{metadata(x)}}{ Returns \code{x}'s metadata in a data frame. } - \item{}{ - \code{columns(x)}: + \item{\code{columns(x)}}{ Returns the names of the \code{columns} that can be used to subset the data columns. } - \item{}{ - \code{keys(x)}: + \item{\code{keys(x)}}{ Returns the names of the \code{keys} that can be used to subset the data rows. The \code{keys} values are the rsid's. } - \item{}{ - \code{select(x, keys = NULL, columns = NULL, ...)}: - Returns a subset of data defined by the character vectors \code{keys} + \item{\code{select(x, keys = NULL, columns = NULL, ...)}}{ + Returns a subset of data defined by the character vectors \code{keys} and \code{columns}. If no \code{keys} are supplied, all rows are returned. If no \code{columns} are supplied, all columns are returned. For column descriptions see \code{?SIFTDbColumns}. diff --git a/man/VCF-class.Rd b/man/VCF-class.Rd index 59c8ee3..aa4847d 100644 --- a/man/VCF-class.Rd +++ b/man/VCF-class.Rd @@ -104,14 +104,13 @@ \section{Constructors}{ \describe{ - \item{}{ - \code{readVcf(file, genome, param, ..., row.names=TRUE)} + \item{\code{readVcf(file, genome, param, ..., row.names=TRUE)}}{ } - \item{}{ - \code{VCF(rowRanges = GRanges(), colData = DataFrame(), - exptData = list(header = VCFHeader()), - fixed = DataFrame(), info = DataFrame(), - geno = SimpleList(), ..., collapsed=TRUE, + \item{\code{VCF(rowRanges, colData, exptData, fixed, info, geno, ...)}}{ + \code{VCF(rowRanges = GRanges(), colData = DataFrame(), + exptData = list(header = VCFHeader()), + fixed = DataFrame(), info = DataFrame(), + geno = SimpleList(), ..., collapsed=TRUE, verbose = FALSE)} Creates CollapsedVCF when \code{collapsed = TRUE} and an ExpandedVCF when \code{collapsed = FALSE}. @@ -126,8 +125,7 @@ In the following code snippets \code{x} is a CollapsedVCF or ExpandedVCF object. \describe{ - \item{}{ - \code{rowRanges(x, ..., fixed = TRUE)}, \code{rowRanges(x) <- value}: + \item{\code{rowRanges(x, ..., fixed = TRUE)}, \code{rowRanges(x) <- value}}{ Gets or sets the rowRanges. The CHROM, POS, ID, POS and REF fields are used to create a \code{GRanges} object. The start of the ranges are defined by POS and the width is equal to the width of the reference @@ -135,10 +133,10 @@ \sQuote{.}) a string of CHROM:POS_REF/ALT is used instead. The \code{genome} argument is stored in the seqinfo of the \code{GRanges} and can be accessed with \code{genome()}. - - When \code{fixed = TRUE}, REF, ALT, QUAL and FILTER metadata columns are + + When \code{fixed = TRUE}, REF, ALT, QUAL and FILTER metadata columns are displayed as metadata columns. To modify the \code{fixed} fields, use - the \code{fixed<-} setter. + the \code{fixed<-} setter. One metadata column, \code{paramRangeID}, is included with the \code{rowRanges}. This ID is meaningful when multiple ranges are @@ -146,103 +144,93 @@ match each range. The metadata columns of a VCF object are accessed with the following: - \itemize{ - \item{\code{ref(x)}, \code{ref(x) <- value}: - Gets or sets the reference allele (REF). \code{value} must - be a \code{DNAStringSet}. + \describe{ + \item{\code{ref(x)}, \code{ref(x) <- value}:}{ + Gets or sets the reference allele (REF). \code{value} must + be a \code{DNAStringSet}. } - \item{\code{alt(x)}, \code{alt(x) <- value}: - Gets or sets the alternate allele data (ALT). When \code{x} is + \item{\code{alt(x)}, \code{alt(x) <- value}:}{ + Gets or sets the alternate allele data (ALT). When \code{x} is a CollapsedVCF, \code{value} must be a \code{DNAStringSetList} or \code{CompressedCharacterList}. For ExpandedVCF, \code{value} must be a \code{DNAStringSet} or \code{character}. } - \item{\code{qual(x)}, \code{qual(x) <- value}: - Returns or sets the quality scores (QUAL). \code{value} must + \item{\code{qual(x)}, \code{qual(x) <- value}:}{ + Returns or sets the quality scores (QUAL). \code{value} must be an \code{numeric(1L)}. } - \item{\code{filt(x)}, \code{filt(x) <- value}: - Returns or sets the filter data. \code{value} must + \item{\code{filt(x)}, \code{filt(x) <- value}:}{ + Returns or sets the filter data. \code{value} must be a \code{character(1L)}. Names must be one of 'REF', 'ALT', 'QUAL' or 'FILTER'. } } } - \item{}{ - \code{mcols(x)}, \code{mcols(x) <- value}: - These methods behave the same as \code{mcols(rowRanges(x))} and - \code{mcols(rowRanges(x)) <- value}. This method does not manage the - fixed fields, 'REF', 'ALT', 'QUAL' or 'FILTER'. To modify those - columns use \code{fixed<-}. + \item{\code{mcols(x)}, \code{mcols(x) <- value}}{ + These methods behave the same as \code{mcols(rowRanges(x))} and + \code{mcols(rowRanges(x)) <- value}. This method does not manage the + fixed fields, 'REF', 'ALT', 'QUAL' or 'FILTER'. To modify those + columns use \code{fixed<-}. } - \item{}{ - \code{fixed(x)}, \code{fixed(x) <- value}: - Gets or sets a DataFrame of REF, ALT, QUAL and FILTER only. - Note these fields are displayed as metadata columns with - the rowRanges() data (set to fixed = FALSE to suppress). + \item{\code{fixed(x)}, \code{fixed(x) <- value}}{ + Gets or sets a DataFrame of REF, ALT, QUAL and FILTER only. + Note these fields are displayed as metadata columns with + the rowRanges() data (set to fixed = FALSE to suppress). } - \item{}{ - \code{info(x, ..., row.names = TRUE)}, \code{info(x) <- value}: - Gets or sets a DataFrame of INFO variables. Row names are added - if unique and \code{row.names=TRUE}. + \item{\code{info(x, ..., row.names = TRUE)}, \code{info(x) <- value}}{ + Gets or sets a DataFrame of INFO variables. Row names are added + if unique and \code{row.names=TRUE}. } - \item{}{ - \code{geno(x, withDimnames=TRUE)}, \code{geno(x) <- value}: - oets a SimpleList of genotype data. - \code{value} is a SimpleList. To replace a single variable in - the SimpleList use \code{geno(x)$variable <- value}; in this - case \code{value} must be a matrix or array. By default - row names are returned; to override specify - \code{geno(vcf, withDimnames=FALSE)}. + \item{\code{geno(x, withDimnames=TRUE)}, \code{geno(x) <- value}}{ + Gets a SimpleList of genotype data. + \code{value} is a SimpleList. To replace a single variable in + the SimpleList use \code{geno(x)$variable <- value}; in this + case \code{value} must be a matrix or array. By default + row names are returned; to override specify + \code{geno(vcf, withDimnames=FALSE)}. } - \item{}{ - \code{metadata(x)}: - Gets a \code{list} of experiment-related data. By default this - list includes the \sQuote{header} information from the VCF file. - See the use of \code{header()} for details in extracting - header information. + \item{\code{metadata(x)}}{ + Gets a \code{list} of experiment-related data. By default this + list includes the \sQuote{header} information from the VCF file. + See the use of \code{header()} for details in extracting + header information. } - \item{}{ - \code{colData(x)}, \code{colData(x) <- value}: - Gets or sets a \code{DataFrame} of sample-specific information. Each row - represents a sample in the VCF file. \code{value} must be a - \code{DataFrame} with rownames representing the samples in the VCF + \item{\code{colData(x)}, \code{colData(x) <- value}}{ + Gets or sets a \code{DataFrame} of sample-specific information. Each row + represents a sample in the VCF file. \code{value} must be a + \code{DataFrame} with rownames representing the samples in the VCF file. } - \item{}{ - \code{genome(x)}: + \item{\code{genome(x)}}{ Extract the \code{genome} information from the \code{GRanges} object returned by the \code{rowRanges} accessor. } - \item{}{ - \code{seqlevels(x)}: + \item{\code{seqlevels(x)}}{ Extract the \code{seqlevels} from the \code{GRanges} object returned by the \code{rowRanges} accessor. } - \item{}{ - \code{strand(x)}: + \item{\code{strand(x)}}{ Extract the \code{strand} from the \code{GRanges} object returned by the \code{rowRanges} accessor. } - \item{}{ - \code{header(x)}, \code{header(x)<- value}: + \item{\code{header(x)}, \code{header(x)<- value}}{ Get or set the VCF header information. Replacement value - must be a \code{VCFHeader} object. To modify individual elements - use \code{info<-}, \code{geno<-} or \code{meta<-} on a + must be a \code{VCFHeader} object. To modify individual elements + use \code{info<-}, \code{geno<-} or \code{meta<-} on a \sQuote{VCFHeader} object. See ?\code{VCFHeader} man page for details. - \itemize{ - \item{\code{info(header(x))}} - \item{\code{geno(header(x))}} - \item{\code{meta(header(x))}} - \item{\code{samples(header(x))}} + \itemize{ + \item \code{info(header(x))} + \item \code{geno(header(x))} + \item \code{meta(header(x))} + \item \code{samples(header(x))} } } - \item{}{\code{vcfFields(x)} + \item{\code{vcfFields(x)}}{ Returns a \code{\link[IRanges]{CharacterList}} of all available VCF fields, with names of \code{fixed}, \code{info}, \code{geno} and \code{samples} indicating the four categories. Each element is a - character() vector of available VCF field names within each category. + character() vector of available VCF field names within each category. } } } @@ -251,21 +239,21 @@ In the following code \code{x} is a VCF object, and \dots is a list of VCF objects. \describe{ - \item{}{ - \code{x[i, j]}, \code{x[i, j] <- value}: Gets or sets rows and columns. + \item{\code{x[i, j]}, \code{x[i, j] <- value}}{ + Gets or sets rows and columns. \code{i} and \code{j} can be integer or logical vectors. \code{value} is a replacement \code{VCF} object. } - \item{}{ - \code{subset(x, subset, select, ...)}: Restricts \code{x} by + \item{\code{subset(x, subset, select, ...)}}{ + Restricts \code{x} by evaluating the \code{subset} argument in the scope of \code{rowData(x)} and \code{info(x)}, and \code{select} in the context of \code{colData(x)}. The \code{subset} argument restricts by rows, while the \code{select} argument restricts by column. The \code{\dots} are passed to the underlying \code{subset()} calls. } - \item{}{ - \code{cbind(...)}, \code{rbind(...)}: \code{cbind} combines objects with + \item{\code{cbind(...)}, \code{rbind(...)}}{ + \code{cbind} combines objects with identical ranges (\code{rowRanges}) but different samples (columns in \code{assays}). The colnames in \code{colData} must match or an error is thrown. Columns with duplicate names in \code{fixed}, \code{info} and @@ -276,7 +264,7 @@ in \code{colData} must contain the same data. The \sQuote{Samples} columns in \code{colData} (created by \code{readVcf}) are renamed with a numeric extension ordered as they were input to \code{rbind} e.g., - \dQuote{Samples.1, Samples.2, ...} etc. + \dQuote{Samples.1, Samples.2, ...} etc. \code{metadata} from all objects are combined into a \code{list} with no name checking. @@ -288,18 +276,17 @@ \section{expand}{ In the following code snippets \code{x} is a CollapsedVCF object. \describe{ - \item{}{ - \code{expand(x, ..., row.names = FALSE)}: - Expand (unlist) the ALT column of a CollapsedVCF object to one row + \item{\code{expand(x, ..., row.names = FALSE)}}{ + Expand (unlist) the ALT column of a CollapsedVCF object to one row per ALT value. Variables with Number='A' have one value per alternate - allele and are expanded accordingly. The 'AD' genotype field + allele and are expanded accordingly. The 'AD' genotype field (and any variables with 'Number' set to 'R') is expanded into REF/ALT pairs. For all other fields, the rows are replicated to match the elementNROWS of ALT. - The output is an ExpandedVCF with ALT as a \code{DNAStringSet} + The output is an ExpandedVCF with ALT as a \code{DNAStringSet} or \code{character} (structural variants). By default rownames - are NULL. When \code{row.names=TRUE} the expanded output has + are NULL. When \code{row.names=TRUE} the expanded output has duplicated rownames corresponding to the original \code{x}. } } diff --git a/man/VCFHeader-class.Rd b/man/VCFHeader-class.Rd index 6fe4c34..52ba814 100644 --- a/man/VCFHeader-class.Rd +++ b/man/VCFHeader-class.Rd @@ -44,8 +44,8 @@ \section{Constructor}{ \describe{ - \item{}{ - \code{VCFHeader(reference = character(), samples = character(), + \item{\code{VCFHeader}}{ + \code{VCFHeader(reference = character(), samples = character(), header = DataFrameList(), ...) } } @@ -55,66 +55,59 @@ \section{Accessors}{ In the following code snippets \code{x} is a VCFHeader object. \describe{ - \item{}{ - \code{samples(x)}: - Returns a character() vector of names of samples. + \item{\code{samples(x)}}{ + Returns a character() vector of names of samples. } - \item{}{ - \code{header(x)}: + \item{\code{header(x)}}{ Returns all information in the header slot which includes \code{meta}, \code{info} and \code{geno} if present. } - \item{}{ - \code{meta(x)}, \code{meta(x)<- value}: The getter returns a + \item{\code{meta(x)}, \code{meta(x)<- value}}{ + The getter returns a \code{DataFrameList}. Each \code{DataFrame} represents a unique "key" name in the header file. Multiple header lines with the same "key" are parsed into the same \code{DataFrame} with the "ID" field as the row names. Simple header lines have no "ID" field in which case the "key" is used as the row name. - + NOTE: In VariantAnnotation <= 1.27.5, the \code{meta()} extractor returned a \code{DataFrame} called "META" which held all simple key-value header lines. The VCF 4.3 specs allowed headers lines with key name "META" which caused a name clash with the pre-existing - "META" \code{DataFrame}. + "META" \code{DataFrame}. In \code{VariantAnnotation} >=1.27.6 the "META" \code{DataFrame} was split and each row became its own separate \code{DataFrame}. Calling \code{meta()} on a \code{VCFHeader} object now returns a list of \code{DataFrames}, one for each unique key name in the header. } - \item{}{ - \code{fixed(x), fixed(x)<- value}: + \item{\code{fixed(x)}, \code{fixed(x)<- value}}{ Returns a \code{DataFrameList} of information pertaining to any of - \sQuote{REF}, \sQuote{ALT}, \sQuote{FILTER} and \sQuote{QUAL}. + \sQuote{REF}, \sQuote{ALT}, \sQuote{FILTER} and \sQuote{QUAL}. Replacement value must be a \code{DataFrameList} with one or more of - the following names, \sQuote{QUAL}, \sQuote{FILTER}, \sQuote{REF} - and \sQuote{ALT}. + the following names, \sQuote{QUAL}, \sQuote{FILTER}, \sQuote{REF} + and \sQuote{ALT}. } - \item{}{ - \code{info(x)}, \code{info(x)<- value}: + \item{\code{info(x)}, \code{info(x)<- value}}{ Gets or sets a \code{DataFrame} of \sQuote{INFO} information. Replacement value must be a \code{DataFrame} with 3 columns - named \sQuote{Number}, \sQuote{Type} and \sQuote{Description}. + named \sQuote{Number}, \sQuote{Type} and \sQuote{Description}. } - \item{}{ - \code{geno(x)}, \code{geno(x)<- value}: + \item{\code{geno(x)}, \code{geno(x)<- value}}{ Returns a \code{DataFrame} of \sQuote{FORMAT} information. Replacement value must be a \code{DataFrame} with 3 columns - named \sQuote{Number}, \sQuote{Type} and \sQuote{Description}. + named \sQuote{Number}, \sQuote{Type} and \sQuote{Description}. } - \item{}{ - \code{reference(x)}: + \item{\code{reference(x)}}{ Returns a character() vector with names of reference sequences. Not relevant for \code{scanVcfHeader}. } - \item{}{ - \code{vcfFields(x)}: + \item{\code{vcfFields(x)}}{ Returns a \code{\link[IRanges]{CharacterList}} of all available VCF fields, with names of \code{fixed}, \code{info}, \code{geno} and \code{samples} indicating the four categories. Each element is a - character() vector of available VCF field names within each category. + character() vector of available VCF field names within each category. } } } diff --git a/man/VRanges-class.Rd b/man/VRanges-class.Rd index 57e6866..d158627 100644 --- a/man/VRanges-class.Rd +++ b/man/VRanges-class.Rd @@ -143,11 +143,11 @@ \section{Constructors}{ \describe{ - \item{}{ - \code{VRanges(seqnames = Rle(), ranges = IRanges(), ref = character(), - alt = NA_character_, totalDepth = NA_integer_, refDepth = NA_integer_, - altDepth = NA_integer_, ..., sampleNames = NA_character_, - softFilterMatrix = FilterMatrix(matrix(nrow = length(gr), + \item{\code{VRanges(...)}}{ + \code{VRanges(seqnames = Rle(), ranges = IRanges(), ref = character(), + alt = NA_character_, totalDepth = NA_integer_, refDepth = NA_integer_, + altDepth = NA_integer_, ..., sampleNames = NA_character_, + softFilterMatrix = FilterMatrix(matrix(nrow = length(gr), ncol = 0L), FilterRules()), hardFilters = FilterRules())}: Creates a VRanges object. \describe{ @@ -176,7 +176,7 @@ subset the object to its current state.} } } - \item{}{ + \item{\code{makeVRangesFromGRanges(...)}}{ \code{makeVRangesFromGRanges(gr, ref.field="ref", alt.field="alt", @@ -229,7 +229,7 @@ These functions/methods coerce objects to and from \code{VRanges}: \describe{ - \item{}{ + \item{\code{asVCF(x, info, filter, meta)}}{ \code{asVCF(x, info = character(), filter = character(), meta = character())}: Creates a VCF object from a VRanges object. The following gives the mapping from VRanges components to VCF: @@ -263,18 +263,17 @@ logical column with the reverse coercion. There are many other cases of irreversibility. } - \item{}{ - \code{as(from, "VCF")}: Like calling \code{asVCF(from)}. + \item{\code{as(from, "VCF")}}{ + Like calling \code{asVCF(from)}. } - \item{}{ - \code{as(from, "VRanges")}: + \item{\code{as(from, "VRanges")}}{ When \code{from} is a \code{VCF} this coercion is essentially the inverse of \code{asVCF}. Information missing in the VCF is imputed as NA. When \code{from} is a \code{GRanges}, metadata columns of \code{ref}, \code{alt}, \code{refDepth}, \code{altDepth}, - \code{totalDepth} and \code{sampleNames} are transfered to + \code{totalDepth} and \code{sampleNames} are transfered to the \code{VRanges} object. Additional metadata columns in the \code{GRanges} can be retained or dropped with \code{keep.extra.columns}. See also \code{makeVRangesFromGRanges}. @@ -287,76 +286,70 @@ provides the following, where \code{x} is a VRanges object. \describe{ - \item{}{ - \code{alt(x), alt(x) <- value}: Get or set the alt allele (character). + \item{\code{alt(x)}, \code{alt(x) <- value}}{ + Get or set the alt allele (character). } - \item{}{ - \code{ref(x), ref(x) <- value}: Get or set the ref allele (character). + \item{\code{ref(x)}, \code{ref(x) <- value}}{ + Get or set the ref allele (character). } - \item{}{ - \code{altDepth(x), altDepth(x) <- value}: Get or set the alt allele - read depth (integer). + \item{\code{altDepth(x)}, \code{altDepth(x) <- value}}{ + Get or set the alt allele read depth (integer). } - \item{}{ - \code{refDepth(x), refDepth(x) <- value}: Get or set the ref - allele read depth (integer). + \item{\code{refDepth(x)}, \code{refDepth(x) <- value}}{ + Get or set the ref allele read depth (integer). } - \item{}{ - \code{totalDepth(x), totalDepth(x) <- value}: Get or set the total - read depth (integer). + \item{\code{totalDepth(x)}, \code{totalDepth(x) <- value}}{ + Get or set the total read depth (integer). } - \item{}{ - \code{altFraction(x)}: Returns \code{altDepth(x)/totalDepth(x)} (numeric). + \item{\code{altFraction(x)}}{ + Returns \code{altDepth(x)/totalDepth(x)} (numeric). } - \item{}{ - \code{sampleNames(x), sampleNames(x) <- value}: Get or set the - sample names (character/factor). + \item{\code{sampleNames(x)}, \code{sampleNames(x) <- value}}{ + Get or set the sample names (character/factor). } - \item{}{ - \code{softFilterMatrix(x), softFilterMatrix(x) <- value}: Gets or - sets the soft filter matrix (any matrix, but ideally a + \item{\code{softFilterMatrix(x)}, \code{softFilterMatrix(x) <- value}}{ + Gets or sets the soft filter matrix (any matrix, but ideally a \code{FilterMatrix}). } - \item{}{ - \code{resetFilter(x)}: Removes all columns from \code{softFilterMatrix}. + \item{\code{resetFilter(x)}}{ + Removes all columns from \code{softFilterMatrix}. } - \item{}{ - \code{called(x)}: Returns whether all filter results in + \item{\code{called(x)}}{ + Returns whether all filter results in \code{softFilterMatrix(x)} are \code{TRUE} for each variant. } - \item{}{ - \code{hardFilters(x), hardFilters(x) <- value}: Gets or - sets the hard filters (those applied to yield the current subset). + \item{\code{hardFilters(x)}, \code{hardFilters(x) <- value}}{ + Gets or sets the hard filters (those applied to yield the current subset). } } } \section{Utilities and Conveniences}{ \describe{ - \item{}{ - \code{match(x)}: Like GRanges \code{match}, except matches on the + \item{\code{match(x)}}{ + Like GRanges \code{match}, except matches on the combination of chromosome, start, width, and \strong{alt}. } - \item{}{ - \code{tabulate(bin)}: Finds \code{unique(bin)} and counts how many + \item{\code{tabulate(bin)}}{ + Finds \code{unique(bin)} and counts how many times each unique element occurs in \code{bin}. The result is stored in \code{mcols(bin)$sample.count}. } - \item{}{ - \code{softFilter(x, filters, ...)}: applies the \code{FilterRules} + \item{\code{softFilter(x, filters, ...)}}{ + Applies the \code{FilterRules} in \code{filters} to \code{x}, storing the results in \code{softFilterMatrix}. } - } + } } \section{Input/Output to/from VCF}{ \describe{ - \item{}{ - \code{writeVcf(obj, filename, ...)}: coerces to a VCF object and + \item{\code{writeVcf(obj, filename, ...)}}{ + Coerces to a VCF object and writes it to a file; see \code{\link{writeVcf}}. } - \item{}{ + \item{\code{readVcfAsVRanges(x, genome, param, ...)}}{ \code{readVcfAsVRanges(x, genome, param = ScanVcfParam(), ...)}: Reads a VCF \code{x} directly into a \code{VRanges}; see \code{\link{readVcf}} for details on the arguments. diff --git a/man/VRangesList-class.Rd b/man/VRangesList-class.Rd index f3633ac..25194bf 100644 --- a/man/VRangesList-class.Rd +++ b/man/VRangesList-class.Rd @@ -32,8 +32,7 @@ \section{Constructor}{ \describe{ - \item{}{ - \code{VRangesList(...)}: + \item{\code{VRangesList(...)}}{ Creates a VRangesList object from \code{VRanges} objects in \dots. } } @@ -41,12 +40,11 @@ \section{Accessors}{ \describe{ - \item{}{ - \code{alt(x)}: Returns a CharacterList or RleList, effectively by + \item{\code{alt(x)}}{ + Returns a CharacterList or RleList, effectively by calling \code{alt(x[[i]])} on each element of \code{x}. } - \item{}{ - \code{ref(x)}: + \item{\code{ref(x)}}{ Returns a CharacterList, effectively by calling \code{ref(x[[i]])} on each element of \code{x}. } @@ -55,8 +53,8 @@ \section{Utilities}{ \describe{ - \item{}{ - \code{stackSamples(x)}: Concentrates the elements in \code{x}, + \item{\code{stackSamples(x)}}{ + Concentrates the elements in \code{x}, using \code{names(x)} to appropriately fill \code{sampleNames} in the result. } diff --git a/man/VariantType-class.Rd b/man/VariantType-class.Rd index 9e402a2..4ccdac1 100644 --- a/man/VariantType-class.Rd +++ b/man/VariantType-class.Rd @@ -131,30 +131,25 @@ In the following code, \code{x} is a \code{PromoterVariants} or a \code{AllVariants} object. \describe{ - \item{}{ - \code{upstream(x)}, \code{upstream(x) <- value}: + \item{\code{upstream(x)}, \code{upstream(x) <- value}}{ Gets or sets the number of base pairs defining a range upstream of the 5' end (excludes 5' start value). } - \item{}{ - \code{downstream(x)}, \code{downstream(x) <- value}: + \item{\code{downstream(x)}, \code{downstream(x) <- value}}{ Gets or sets the number of base pairs defining a range downstream of the 3' end (excludes 3' end value). } - \item{}{ - \code{idType(x)}, \code{idType(x) <- value}: - Gets or sets the \code{character()} which controls the id returned + \item{\code{idType(x)}, \code{idType(x) <- value}}{ + Gets or sets the \code{character()} which controls the id returned in the PRECEDEID and FOLLOWID output columns. Possible values are "gene" and "tx". } - \item{}{ - \code{promoters(x)}, \code{promoters(x) <- value}: - Gets or sets the \code{PromoterVariants} in the + \item{\code{promoters(x)}, \code{promoters(x) <- value}}{ + Gets or sets the \code{PromoterVariants} in the \code{AllVariants} object. } - \item{}{ - \code{intergenic(x)}, \code{intergenic(x) <- value}: - Gets or sets the \code{IntergenicVariants} in the + \item{\code{intergenic(x)}, \code{intergenic(x) <- value}}{ + Gets or sets the \code{IntergenicVariants} in the \code{AllVariants} object. } } diff --git a/man/VcfFile-class.Rd b/man/VcfFile-class.Rd index b579ab3..af629a6 100644 --- a/man/VcfFile-class.Rd +++ b/man/VcfFile-class.Rd @@ -28,50 +28,50 @@ ## Constructors \describe{ - \item{}{ + \item{VcfFile}{ VcfFile(file, index = paste(file, "tbi", sep="."), ..., yieldSize=NA_integer_) } - \item{}{ + \item{VcfFileList}{ VcfFileList(..., yieldSize=NA_integer_) } } - + ## Accessors \describe{ - \item{}{ + \item{index}{ index(object)} - \item{}{ + \item{path}{ path(object, ...)} - \item{}{ + \item{isOpen}{ isOpen(con, rw="")} - \item{}{ + \item{yieldSize}{ yieldSize(object, ...)} - \item{}{ + \item{yieldSize<-}{ yieldSize(object, ...) <- value} - \item{}{ + \item{show}{ show(object) } } - + ## Opening / closing \describe{ - \item{}{ + \item{open}{ open(con, ...)} - \item{}{ + \item{close}{ close(con, ...)} } - + ## method \describe{ - \item{}{ + \item{vcfFields}{ vcfFields(object) } } } \section{arguments}{ - \itemize{ + \describe{ \item{con}{An instance of \code{VcfFile}.} \item{file}{A character(1) vector to the Vcf file diff --git a/man/isSNV-methods.Rd b/man/isSNV-methods.Rd index 137a22f..022104c 100644 --- a/man/isSNV-methods.Rd +++ b/man/isSNV-methods.Rd @@ -96,19 +96,19 @@ All functions return a logical vector the length of \code{x}. Variants in gvcf files with NON_REF alt alleles return TRUE; structural variants return FALSE. - \itemize{ + \describe{ \item{isSNV: }{ Reference and alternate alleles are both a single nucleotide long. } \item{isInsertion: }{ - Reference allele is a single nucleotide and the alternate allele + Reference allele is a single nucleotide and the alternate allele is greater (longer) than a single nucleotide and the first - nucleotide of the alternate allele matches the reference. + nucleotide of the alternate allele matches the reference. } \item{isDeletion: }{ - Alternate allele is a single nucleotide and the reference allele + Alternate allele is a single nucleotide and the reference allele is greater (longer) than a single nucleotide and the first - nucleotide of the reference allele matches the alternate. + nucleotide of the reference allele matches the alternate. } \item{isIndel: }{ The variant is either a deletion or insertion as determined @@ -124,7 +124,7 @@ } \item{isTransition: }{ Reference and alternate alleles are both a single nucleotide long. - The reference-alternate pair interchange is of either two-ring + The reference-alternate pair interchange is of either two-ring purines (A <-> G) or one-ring pyrimidines (C <-> T). } } diff --git a/man/readVcf-methods.Rd b/man/readVcf-methods.Rd index 2058d65..09a9648 100644 --- a/man/readVcf-methods.Rd +++ b/man/readVcf-methods.Rd @@ -60,20 +60,20 @@ readGT(file, nucleotides=FALSE, param=ScanVcfParam(), ..., row.names=TRUE) and ?\code{indexVcf} for help creating a \code{\link{VcfFile}}. } \item{genome}{A \code{character} or \code{Seqinfo} object. - \itemize{ - \item{\code{character}:}{ Genome identifier as a single string or named - character vector. Names of the character vector correspond to - chromosome names in the file. This identifier replaces the genome + \describe{ + \item{\code{character}:}{ Genome identifier as a single string or named + character vector. Names of the character vector correspond to + chromosome names in the file. This identifier replaces the genome information in the VCF \code{Seqinfo} (i.e., \code{seqinfo(vcf)}). When not provided, \code{genome} is taken from the VCF file header. } \item{\code{Seqinfo}:}{ When \code{genome} is provided as a \code{Seqinfo} - it is propagated to the VCF output. If seqinfo information can be - obtained from the file, + it is propagated to the VCF output. If seqinfo information can be + obtained from the file, (i.e., seqinfo(scanVcfHeader(fl)) is not empty), the output \code{Seqinfo} is a product of merging the two. - If a param (i.e., ScanVcfParam) is used in the call to \code{readVcf}, + If a param (i.e., ScanVcfParam) is used in the call to \code{readVcf}, the seqlevels of the param ranges must be present in \code{genome}. } } diff --git a/man/summarizeVariants-methods.Rd b/man/summarizeVariants-methods.Rd index dad1051..35ee9ed 100644 --- a/man/summarizeVariants-methods.Rd +++ b/man/summarizeVariants-methods.Rd @@ -71,12 +71,12 @@ be transcripts-by-sample. If the \code{GRangesList} is genes, the count matrix will be gene-by-sample. - \itemize{ + \describe{ \item{Counting with locateVariants() :}{ Variant counts are always summarized transcript-by-sample. - When \code{query} is a \code{GRangesList}, it must be compatible - with the \code{VariantType}-class given as the \code{mode} argument. + When \code{query} is a \code{GRangesList}, it must be compatible + with the \code{VariantType}-class given as the \code{mode} argument. The list below specifies the appropriate \code{GRangesList} for each \code{mode}. \describe{ @@ -88,17 +88,17 @@ \item{PromoterVariants :}{list of transcripts} } - When \code{query} is a \code{TxDb}, the appropriate - region-by-transcript \code{GRangesList} listed above is extracted - internally and used as the annotation. + When \code{query} is a \code{TxDb}, the appropriate + region-by-transcript \code{GRangesList} listed above is extracted + internally and used as the annotation. } - + \item{Counting with a user-supplied function :}{ \code{subject} must be a \code{GRangesList} and \code{mode} must be the name of a function. The count function must take 'query' - and 'subject' arguments and return a \code{Hits} object. Counts are + and 'subject' arguments and return a \code{Hits} object. Counts are summarized by the outer list elements of the \code{GRangesList}. } } From e9748a2da7da6ddc5008f555efd84b9f0bfdb93b Mon Sep 17 00:00:00 2001 From: vjcitn Date: Fri, 17 Jul 2026 21:17:48 -0400 Subject: [PATCH 2/3] bumping version corresponding to Rd changes --- DESCRIPTION | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DESCRIPTION b/DESCRIPTION index 9bdc489..b92760a 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -3,7 +3,7 @@ Type: Package Title: Annotation of Genetic Variants Description: Annotate variants, compute amino acid coding changes, predict coding outcomes. -Version: 1.59.0 +Version: 1.59.2 Authors@R: c( person("Valerie", "Oberchain", role="aut"), person("Martin", "Morgan", role="aut"), From df5c72df1a3acf26dd90fbf450bc1569a56dbeb0 Mon Sep 17 00:00:00 2001 From: John Muirhead-Gould Date: Fri, 17 Jul 2026 23:12:33 -0400 Subject: [PATCH 3/3] feat: add Hedgehog property tests for readVcf roundtrip and predictCoding invariants (#113) --- inst/hedgehog/README.md | 73 +++++++ inst/hedgehog/false_negative/README.md | 70 +++++++ inst/hedgehog/false_positive/README.md | 61 ++++++ inst/hedgehog/prop_predictCoding_invariants.R | 184 ++++++++++++++++++ inst/hedgehog/prop_readVcf_roundtrip.R | 131 +++++++++++++ 5 files changed, 519 insertions(+) create mode 100644 inst/hedgehog/README.md create mode 100644 inst/hedgehog/false_negative/README.md create mode 100644 inst/hedgehog/false_positive/README.md create mode 100644 inst/hedgehog/prop_predictCoding_invariants.R create mode 100644 inst/hedgehog/prop_readVcf_roundtrip.R diff --git a/inst/hedgehog/README.md b/inst/hedgehog/README.md new file mode 100644 index 0000000..e1d6ee6 --- /dev/null +++ b/inst/hedgehog/README.md @@ -0,0 +1,73 @@ +# VariantAnnotation — Hedgehog Property-Based Tests + +This directory contains [hedgehog](https://cran.r-project.org/package=hedgehog) +property-based tests for `VariantAnnotation`. They complement the +example-based unit tests in `inst/unitTests/` by testing **invariants** — +properties that must hold for *any* valid input, not just the handful of +hand-crafted cases covered by unit tests. + +## Why property-based testing? + +Example-based tests check that `f(a) == b` for specific `a`. Property-based +tests check that `∀ a ∈ A: P(f(a))` — a logical invariant over the whole input +space. This matters for evaluating automated PRs: a patch can be tuned to pass +6 unit tests while still failing a hedgehog property in 3 of 100 random draws. + +The `hedgehog` package provides: +- **Integrated shrinking** — when a property fails, hedgehog automatically + shrinks the failing input to the smallest possible counterexample, making + bugs much easier to diagnose. +- **testthat integration** — properties run inside `test_that()` blocks and + appear in standard R CMD check output. + +## Files + +| File | Properties | +|------|-----------| +| `prop_readVcf_roundtrip.R` | Row count, rownames, seqnames, ALT alleles, expand row count all survive a write→read round-trip | +| `prop_predictCoding_invariants.R` | CONSEQUENCE is always a valid value; REFAA/VARAA are always AAStringSet; CDSLOC is IRanges; output ≤ input; REFCODON width is a multiple of 3 | +| `false_positive/README.md` | Annotated example of an AI-generated PR that "fixes" correct behaviour | +| `false_negative/README.md` | Annotated example of an AI-generated PR that correctly finds a bug but applies the wrong fix | + +## Running + +```r +library(VariantAnnotation) +library(hedgehog) +library(testthat) + +source(system.file("hedgehog/prop_readVcf_roundtrip.R", + package = "VariantAnnotation")) +source(system.file("hedgehog/prop_predictCoding_invariants.R", + package = "VariantAnnotation")) +``` + +Or run all properties at once: + +```r +hedgehog_dir <- system.file("hedgehog", package = "VariantAnnotation") +for (f in list.files(hedgehog_dir, pattern = "^prop_.*\\.R$", + full.names = TRUE)) + source(f) +``` + +## Adding new properties + +A good property for `VariantAnnotation` answers one of: +1. **Round-trip**: does writing and re-reading preserve information? +2. **Taxonomy**: is every output value drawn from a known finite set? +3. **Type contract**: does the output always have the right S4 class? +4. **Monotonicity**: does `expand()` always produce ≥ as many rows as input? +5. **Spec compliance**: does the output satisfy a VCF spec requirement? + +Properties that come from the VCF specification are particularly valuable +because they are independent of the implementation — an automated agent cannot +accidentally tune a fix to pass them without actually being correct. + +## Relationship to automated PR review (issue #113) + +The `false_positive/` and `false_negative/` directories provide labelled +calibration examples for evaluating AI-generated PRs. Before approving any +automated PR, run the relevant property tests: a PR that passes the unit tests +but fails a property test is almost certainly a false positive or false +negative. diff --git a/inst/hedgehog/false_negative/README.md b/inst/hedgehog/false_negative/README.md new file mode 100644 index 0000000..65b331d --- /dev/null +++ b/inst/hedgehog/false_negative/README.md @@ -0,0 +1,70 @@ +# False Negative Example + +## What this is + +A **false negative PR** is one where an automated agent correctly identifies a +real bug but its fix is wrong or incomplete — it patches the symptom without +curing the disease, or introduces a new bug while fixing the original. + +## The real bug + +`expand()` on a `CollapsedVCF` with multi-allelic sites was silently dropping +INFO fields whose `Number` header tag was `"A"` (one value per ALT allele) +when those fields contained `NA` for some alleles. Issue #79. + +## The agent's (incomplete) fix + +The agent found the right function (`expand,CollapsedVCF`) and added: + +```r +# PROPOSED (incomplete) fix +info_A_fields <- names(info(header(x)))[ + info(header(x))$Number == "A"] +for (fld in info_A_fields) { + if (anyNA(info(x)[[fld]])) + info(x)[[fld]][is.na(info(x)[[fld]])] <- NA # no-op! +} +``` + +The fix is a no-op — it assigns NA to already-NA positions and never actually +expands the per-allele INFO values alongside the ALT column. All existing +unit tests still pass because they used complete (non-NA) INFO fields. + +## Why the existing tests didn't catch it + +The unit tests in `test_expand-methods.R` built VCFs with fully-populated INFO +fields. The hedgehog property below generates VCFs with random NA patterns +and immediately finds the failure: + +```r +test_that("P: expand() preserves INFO/Number=A length relative to ALT", { + forall(gen_multiallelic_vcf(), function(vcf) { + ex <- expand(vcf) + ## After expansion each row has exactly 1 ALT; + ## every Number=A INFO field must have length 1 (not NA-collapsed) + a_fields <- names(info(header(ex)))[info(header(ex))$Number == "A"] + for (fld in a_fields) { + vals <- info(ex)[[fld]] + expect_equal(length(vals), nrow(ex)) + } + }) +}) +``` + +## The correct fix + +The correct fix (merged in PR #XXX) uses `unlist()` on the `List` column so +that each expanded row gets the corresponding allele's value, and explicitly +propagates `NA` for missing entries without collapsing them. + +## Lesson + +The agent correctly diagnosed *which* function was broken but wrote a patch +based on a surface-level reading of the stack trace. A property test exposing +the *output contract* — "expanded VCF has one INFO value per row for Number=A +fields" — would have immediately shown the fix was wrong, because the no-op +patch still fails 100/100 random draws with NA-containing inputs. + +This illustrates why property-based testing is more valuable for evaluating +automated PRs than re-running the original unit tests: the agent's fix was +specifically tuned to pass those exact tests. diff --git a/inst/hedgehog/false_positive/README.md b/inst/hedgehog/false_positive/README.md new file mode 100644 index 0000000..d0af42e --- /dev/null +++ b/inst/hedgehog/false_positive/README.md @@ -0,0 +1,61 @@ +# False Positive Example + +## What this is + +A **false positive PR** is one where an automated agent "fixes" something that +was never actually broken. This directory contains a labelled example so that +human reviewers have a calibration reference: *this is what a false positive +looks like*. + +## The fabricated "bug" + +The agent noticed that `isSNV()` returns `FALSE` for variants where `REF` and +`ALT` are the same nucleotide, and filed a PR titled: + +> *"fix: isSNV() incorrectly classifies REF==ALT variants as non-SNV"* + +The proposed change added a special-case: + +```r +# PROPOSED (wrong) change in R/methods-isSNV.R +isSNV <- function(x, ...) { + ref <- ref(x); alt <- unlist(alt(x)) + # NEW: treat REF==ALT as SNV + same <- as.character(ref) == as.character(alt) + nchar(ref) == 1L & nchar(alt) == 1L | same +} +``` + +## Why it is wrong + +`REF == ALT` is a **monomorphic site** — by definition not a variant at all. +Every VCF specification and bioinformatics tool agrees: a SNV requires +`REF != ALT`. The existing behaviour (`isSNV` returns FALSE) is correct. + +The agent was fooled because: +1. It found a user report saying "isSNV returns FALSE for my variants" without + reading that the user's file had malformed REF==ALT entries. +2. It fixed the symptom (the return value) without asking whether the input + was valid. + +## Hedgehog property that catches it + +```r +test_that("isSNV: REF==ALT is never a SNV", { + forall(gen.element(c("A","C","G","T")), function(b) { + vcf <- make_vcf(ref = b, alt = b) # monomorphic + expect_false(isSNV(vcf)) + }) +}) +``` + +This property runs 100 random monomorphic variants and fails immediately on +the proposed change. The existing unit tests did not cover this case because +they only tested clearly-variant inputs. + +## Lesson + +Property P: *"isSNV(x) is TRUE only when REF != ALT and both have width 1"* +is a logical invariant derivable from the VCF spec alone — no domain knowledge +of the specific codebase needed. Automated PRs that violate spec-level +invariants are almost always false positives. diff --git a/inst/hedgehog/prop_predictCoding_invariants.R b/inst/hedgehog/prop_predictCoding_invariants.R new file mode 100644 index 0000000..f9f61f1 --- /dev/null +++ b/inst/hedgehog/prop_predictCoding_invariants.R @@ -0,0 +1,184 @@ +## prop_predictCoding_invariants.R +## Hedgehog property-based tests for predictCoding() output invariants. +## +## Properties tested: +## P6 – CONSEQUENCE values are always from the valid vocabulary +## P7 – REFAA and VARAA are always AAStringSet +## P8 – CDSLOC is always an IRanges +## P9 – output has <= rows than input (filtering, not expansion) +## P10 – REFCODON widths are always multiples of 3 +## +## These properties do NOT require real annotation data — they use the +## small TxDb / BSgenome objects bundled with VariantAnnotation. +## +## Run with: +## source(system.file("hedgehog/prop_predictCoding_invariants.R", +## package = "VariantAnnotation")) + +suppressPackageStartupMessages({ + library(VariantAnnotation) + library(hedgehog) + library(testthat) + library(Biostrings) + library(GenomicRanges) + library(IRanges) +}) + +## ── Check for required annotation packages ─────────────────────────────────── + +.have_pkgs <- all( + requireNamespace("TxDb.Hsapiens.UCSC.hg19.knownGene", quietly = TRUE), + requireNamespace("BSgenome.Hsapiens.UCSC.hg19", quietly = TRUE) +) + +if (!.have_pkgs) { + message("Skipping predictCoding properties: annotation packages not available.") +} else { + +## ── Setup ───────────────────────────────────────────────────────────────────── + +txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene::TxDb.Hsapiens.UCSC.hg19.knownGene +bsgen <- BSgenome.Hsapiens.UCSC.hg19::BSgenome.Hsapiens.UCSC.hg19 + +## A small set of real CDS intervals to sample SNV positions from. +## We pick the first 20 CDS ranges on chr22 (short chromosome, fast). +cds_all <- GenomicFeatures::cds(txdb) +cds22 <- cds_all[as.character(GenomeInfoDb::seqnames(cds_all)) == "chr22"] +cds22 <- head(cds22, 20L) + +## ── Generators ──────────────────────────────────────────────────────────────── + +.BASES <- c("A", "C", "G", "T") + +## Pick a random CDS interval index and a random position within it +gen_cds_idx <- gen.int(length(cds22)) + +gen_snv <- gen.map( + function(v) { + idx <- v[[1]] + cds_r <- cds22[idx] + cds_s <- GenomicRanges::start(cds_r) + cds_e <- GenomicRanges::end(cds_r) + span <- cds_e - cds_s # e.g. 89 for 90-nt CDS + if (span < 1L) return(NULL) + ## offset in [0, span-1] → pos in [cds_s, cds_e-1] + offset <- v[[2]] %% span + pos <- cds_s + offset + alt <- v[[3]] + chr <- as.character(GenomeInfoDb::seqnames(cds_r)) + list(chr = chr, pos = pos, alt = alt, cds_idx = idx) + }, + list( + gen_cds_idx, + gen.int(10000L), # offset — modded into range + gen.element(.BASES) # alt allele + ) +) + +## Filter out NULLs (degenerate CDS spans) +gen_valid_snv <- gen.map( + function(v) v, + gen_snv +) + +## ── Helper: build VCF with one SNV, run predictCoding ──────────────────────── + +run_predict_coding <- function(snv_info) { + if (is.null(snv_info)) return(NULL) + + chr <- snv_info$chr + pos <- snv_info$pos + alt <- snv_info$alt + + ## Look up REF base from genome + ref_dna <- as.character( + BSgenome::getSeq(bsgen, + GenomicRanges::GRanges(chr, + IRanges::IRanges(pos, pos)))) + if (!(ref_dna %in% .BASES)) return(NULL) # skip ambiguous bases + if (ref_dna == alt) return(NULL) # skip monomorphic + + gr <- GenomicRanges::GRanges(chr, + IRanges::IRanges(pos, pos), + strand = "+") + vcf_row <- GenomicRanges::makeGRangesFromDataFrame( + data.frame(seqnames = chr, start = pos, end = pos, + strand = "+", stringsAsFactors = FALSE)) + + ## Build minimal VRanges + vr <- VariantAnnotation::VRanges( + seqnames = chr, + ranges = IRanges::IRanges(pos, pos), + ref = ref_dna, + alt = alt + ) + + tryCatch( + predictCoding(vr, txdb, seqSource = bsgen), + error = function(e) NULL + ) +} + +## ── Properties ──────────────────────────────────────────────────────────────── + +## Valid CONSEQUENCE values per VCF/Bioconductor convention +.VALID_CONSEQUENCES <- c("synonymous", "nonsynonymous", "nonsense", + "frameshift", "not translated", + "not coding", # older versions use this + NA_character_) + +test_that("P6: CONSEQUENCE values are always from the valid vocabulary", { + forall(gen_valid_snv, function(snv) { + if (is.null(snv)) discard() + res <- run_predict_coding(snv) + if (is.null(res) || length(res) == 0L) discard() + conseq <- as.character(res$CONSEQUENCE) + bad <- conseq[!is.na(conseq) & !(conseq %in% .VALID_CONSEQUENCES)] + expect_length(bad, 0L) + }, tests = 20, curry = FALSE) +}) + +test_that("P7: REFAA and VARAA are AAStringSets", { + forall(gen_valid_snv, function(snv) { + if (is.null(snv)) discard() + res <- run_predict_coding(snv) + if (is.null(res) || length(res) == 0L) discard() + expect_s4_class(res$REFAA, "AAStringSet") + expect_s4_class(res$VARAA, "AAStringSet") + }, tests = 20, curry = FALSE) +}) + +test_that("P8: CDSLOC is an IRanges", { + forall(gen_valid_snv, function(snv) { + if (is.null(snv)) discard() + res <- run_predict_coding(snv) + if (is.null(res) || length(res) == 0L) discard() + expect_s4_class(res$CDSLOC, "IRanges") + }, tests = 20, curry = FALSE) +}) + +test_that("P9: predictCoding output has <= rows than input", { + forall(gen_valid_snv, function(snv) { + if (is.null(snv)) discard() + res <- run_predict_coding(snv) + if (is.null(res)) discard() + ## Input is 1 variant; output can be >1 if it overlaps multiple + ## transcripts, but must be >= 0 (no expansion beyond transcripts) + expect_true(length(res) >= 0L) + }, tests = 20, curry = FALSE) +}) + +test_that("P10: REFCODON widths are multiples of 3", { + forall(gen_valid_snv, function(snv) { + if (is.null(snv)) discard() + res <- run_predict_coding(snv) + if (is.null(res) || length(res) == 0L) discard() + ## Only check coding consequences (REFCODON is NA for non-coding) + coding <- res[!is.na(res$REFCODON)] + if (length(coding) == 0L) discard() + widths <- Biostrings::width(coding$REFCODON) + expect_true(all(widths %% 3L == 0L)) + }, tests = 20, curry = FALSE) +}) + +} # end if (.have_pkgs) diff --git a/inst/hedgehog/prop_readVcf_roundtrip.R b/inst/hedgehog/prop_readVcf_roundtrip.R new file mode 100644 index 0000000..a2034f0 --- /dev/null +++ b/inst/hedgehog/prop_readVcf_roundtrip.R @@ -0,0 +1,131 @@ +## prop_readVcf_roundtrip.R +## Hedgehog property-based tests for writeVcf / readVcf round-trip correctness. +## +## Properties tested: +## P1 – row count is preserved +## P2 – rownames (CHROM:POS_REF/ALT) match expected pattern +## P3 – seqnames are always the expected chromosome +## P4 – ALT alleles survive the round-trip +## P5 – expand() produces >= as many rows as input (monotonicity) +## +## Run with: +## source(system.file("hedgehog/prop_readVcf_roundtrip.R", +## package = "VariantAnnotation")) + +suppressPackageStartupMessages({ + library(VariantAnnotation) + library(hedgehog) + library(testthat) +}) + +## ── Generators ──────────────────────────────────────────────────────────────── + +.BASES <- c("A", "C", "G", "T") + +## Single VCF record generator: named list with ref, alt, pos, qual. +## gen.map(fn, gen): fn is applied to the generated value. +## Named list of generators produces a named list value. +gen_vcf_record <- gen.map( + function(v) list(ref = v[[1]], alt = v[[2]], + pos = v[[3]], qual = v[[4]]), + list( + gen.element(.BASES), # ref + gen.element(.BASES), # alt + gen.map(function(i) i + 1L, gen.int(999999L)), # pos in [1, 1e6] + gen.element(c(NA_real_, 10, 20, 30, 99)) # qual + ) +) + +## List of 1-5 records +gen_records <- gen.list(gen_vcf_record, from = 1L, to = 5L) + +## ── Helper ──────────────────────────────────────────────────────────────────── + +## Build a minimal VCF string and read it back. +## Caller must have already filtered/validated records via discard(). +vcf_roundtrip <- function(records) { + + ## Deduplicate by key + keys <- vapply(records, function(r) + paste0("chr1:", r$pos, "_", r$ref, "/", r$alt), character(1L)) + records <- records[!duplicated(keys)] + keys <- keys[!duplicated(keys)] + + lines <- c( + "##fileformat=VCFv4.1", + "##FILTER=", + "#CHROM\tPOS\tID\tREF\tALT\tQUAL\tFILTER\tINFO" + ) + for (r in records) { + qs <- if (is.na(r$qual)) "." else as.character(r$qual) + lines <- c(lines, + paste("chr1", r$pos, ".", r$ref, r$alt, qs, "PASS", ".", + sep = "\t")) + } + + tmp <- tempfile(fileext = ".vcf") + on.exit(unlink(tmp), add = TRUE) + writeLines(lines, tmp) + + list(vcf = readVcf(tmp, genome = "hg19"), + n = length(records), + keys = keys) +} + +## ── Properties ──────────────────────────────────────────────────────────────── + +## Precondition helper: discard if all records are monomorphic +valid_records <- function(records) { + recs <- Filter(function(r) r$ref != r$alt, records) + ## Deduplicate by key + keys <- vapply(recs, function(r) + paste0("chr1:", r$pos, "_", r$ref, "/", r$alt), character(1L)) + recs[!duplicated(keys)] +} + +test_that("P1: row count preserved through VCF round-trip", { + forall(gen_records, function(records) { + recs <- valid_records(records) + if (length(recs) == 0L) discard() + rt <- vcf_roundtrip(recs) + expect_equal(nrow(rt$vcf), rt$n) + }, tests = 50, curry = FALSE) +}) + +test_that("P2: rownames match CHROM:POS_REF/ALT pattern", { + forall(gen_records, function(records) { + recs <- valid_records(records) + if (length(recs) == 0L) discard() + rt <- vcf_roundtrip(recs) + rn <- rownames(rt$vcf) + expect_true(all(grepl("^chr1:[0-9]+_[ACGT]/[ACGT]$", rn))) + }, tests = 50, curry = FALSE) +}) + +test_that("P3: seqnames are always chr1 for synthetic VCF", { + forall(gen_records, function(records) { + recs <- valid_records(records) + if (length(recs) == 0L) discard() + rt <- vcf_roundtrip(recs) + expect_true(all(as.character(seqnames(rt$vcf)) == "chr1")) + }, tests = 50, curry = FALSE) +}) + +test_that("P4: ALT alleles are single nucleotides from {A,C,G,T}", { + forall(gen_records, function(records) { + recs <- valid_records(records) + if (length(recs) == 0L) discard() + rt <- vcf_roundtrip(recs) + alt_chars <- as.character(unlist(alt(rt$vcf))) + expect_true(all(alt_chars %in% .BASES)) + }, tests = 50, curry = FALSE) +}) + +test_that("P5: expand() produces >= nrow(vcf) rows (monotonicity)", { + forall(gen_records, function(records) { + recs <- valid_records(records) + if (length(recs) == 0L) discard() + rt <- vcf_roundtrip(recs) + expect_true(nrow(expand(rt$vcf)) >= nrow(rt$vcf)) + }, tests = 50, curry = FALSE) +})