class cleanup

DESCRIPTION

@@ -140,7 +140,6 @@ Collate:
     'MultiStudyQtlDataset.R'
     'QtlFineMappingResult.R'
     'SldscData.R'
-    'TupleRangesView.R'
     'causalInferencePipeline.R'
     'colocPipeline.R'
     'qtlSumStats.R'

NAMESPACE

@@ -4,10 +4,8 @@ S3method(as.data.frame,GwasSumStats)
 S3method(dplyr::arrange,RangedTupleList)
 S3method(dplyr::filter,RangedTupleList)
 S3method(dplyr::group_by,RangedTupleList)
-S3method(dplyr::group_by,TupleRangesView)
 S3method(dplyr::mutate,RangedTupleList)
 S3method(dplyr::select,RangedTupleList)
-S3method(dplyr::select,TupleRangesView)
 S3method(dplyr::slice,RangedTupleList)
 S3method(dplyr::summarise,RangedTupleList)
 S3method(postprocessFinemappingFit,mvsusie)
@@ -309,9 +307,12 @@ exportClasses(AnnotationMatrix)
 exportClasses(ColocBoostResult)
 exportClasses(ColocResult)
 exportClasses(ColocResultBase)
+exportClasses(CtwasResult)
+exportClasses(CtwasResultEntry)
 exportClasses(FineMappingResultBase)
 exportClasses(FineMappingRow)
 exportClasses(GwasFineMappingResult)
+exportClasses(GwasSumStats)
 exportClasses(H2Estimate)
 exportClasses(LdData)
 exportClasses(LdEigen)
@@ -321,17 +322,18 @@ exportClasses(MashPrior)
 exportClasses(MultiStudyQtlDataset)
 exportClasses(QtlDataset)
 exportClasses(QtlFineMappingResult)
+exportClasses(QtlSumStats)
 exportClasses(RangedTupleList)
 exportClasses(SldscData)
 exportClasses(SumStatsBase)
-exportClasses(TupleRangesView)
 exportClasses(TwasWeights)
 exportClasses(TwasWeightsRow)
 exportMethods("$")
 exportMethods("$<-")
 exportMethods("[")
 exportMethods("[[<-")
 exportMethods(as.data.frame)
+exportMethods(bindROWS)
 exportMethods(colnames)
 exportMethods(colocboostPipeline)
 exportMethods(computeLdScores)
@@ -492,6 +494,7 @@ importFrom(S4Vectors,
   "mcols<-",
   DataFrame,
   SimpleList,
+  bindROWS,
   endoapply,
   mcols,
   queryHits,

R/AllClasses.R

@@ -79,7 +79,8 @@ setClassUnion("LdMixtureWeights", c("numeric", "NULL"))
 # -----------------------------------------------------------------------------
 # Shared parent of the QTL and GWAS summary statistics collections.
 # Concrete subclasses (QtlSumStats, GwasSumStats) inherit from
-# RangedTupleList and share the ldSketch / genome / qcInfo slots. Each element
+# RangedTupleList and share the ldSketch / qcInfo slots (the genome build
+# lives in seqinfo, not a slot). Each element
 # is one tuple's per-variant GRanges: x[[i]], formerly x$entry[[i]].
 #
 # getZ / getN / getMaf / nSnps are
@@ -92,7 +93,8 @@ setClassUnion("LdMixtureWeights", c("numeric", "NULL"))
 #' @description Virtual base class for QTL and GWAS summary statistics
 #'   collections. Concrete subclasses (\code{QtlSumStats}, \code{GwasSumStats})
 #'   inherit from \code{\linkS4class{RangedTupleList}} and share the
-#'   \code{ldSketch} / \code{genome} / \code{qcInfo} slots.
+#'   \code{ldSketch} / \code{qcInfo} slots, and the genome build in
+#'   \code{seqinfo()}.
 #'
 #'   Each element is the per-variant \code{GRanges} of one tuple, so
 #'   \code{x[[i]]} is that tuple's summary statistics and the identity columns
@@ -102,7 +104,6 @@ setClassUnion("LdMixtureWeights", c("numeric", "NULL"))
 #'   or \code{NULL}. Optional: LD-free workflows (e.g. mash, which operates
 #'   across conditions per variant) carry \code{NULL}; pipelines that need LD
 #'   validate its presence when they consume the collection.
-#' @slot genome Character, genome build label.
 #' @slot qcInfo A \code{list} recording which QC steps ran. Empty \code{list()}
 #'   on construction; populated by \code{summaryStatsQc()} with a per-step audit
 #'   record (filter names, drop counts, liftover target, RAISS settings, etc.).
@@ -115,7 +116,6 @@ setClass(
     contains = c("VIRTUAL", "RangedTupleList"),
     representation(
         ldSketch = "LdSketchOrNULL",
-        genome = "character",
         qcInfo = "list"
     )
 )
@@ -176,12 +176,51 @@ setClass(
     gr[onWindowChrom & IRanges::overlapsAny(gr, win)]
 }
 
+# The build recorded in seqinfo must be exactly one non-NA value. A missing
+# build is an error rather than a default: every downstream liftover / LD
+# join keys on it, and silently guessing hg38 is how a mismatched panel gets
+# through. Mixed builds mean the parts were never comparable.
+# @noRd
+.sumStatsCheckGenome <- function(object) {
+    # seqinfo records the build per SEQLEVEL, so a collection spanning none
+    # has nowhere to keep one -- whether it has no elements at all or only
+    # empty ones (a PIP-screened region emptied by summaryStatsQc). That is
+    # the one case where a missing build is not a defect: there is nothing
+    # for it to describe. A subset that empties an existing collection keeps
+    # its seqinfo (see .rtlRebuild), so this exempts only what was built with
+    # no ranges in the first place.
+    if (length(GenomeInfoDb::seqlevels(object)) == 0L) {
+        return(NULL)
+    }
+    build <- unique(GenomeInfoDb::genome(object))
+    build <- build[!is.na(build)]
+    if (length(build) == 1L && str_length(build) > 0L) {
+        return(NULL)
+    }
+    if (length(build) == 0L) {
+        return("no genome build in seqinfo(); set one with genome(x) <- ...")
+    }
+    str_c(
+        "seqinfo() names more than one genome build (",
+        str_flatten(build, ", "),
+        ")"
+    )
+}
+
 #' @rdname getGenome
 #' @examples
 #' data(qtlSumStatsExample)
 #' getGenome(qtlSumStatsExample)
 #' @export
-setMethod("getGenome", "SumStatsBase", function(x, ...) x@genome)
+setMethod("getGenome", "SumStatsBase", function(x, ...) {
+    # The build lives in seqinfo, exactly as it does on LdStatistic: a
+    # GRangesList already has somewhere to keep it, and a parallel `genome`
+    # slot went stale against it (getGenome() said hg19 while genome(x) said
+    # NA, so every Bioconductor path that reads genome(x) saw nothing).
+    build <- unique(GenomeInfoDb::genome(x))
+    build <- build[!is.na(build)]
+    if (length(build) == 0L) NA_character_ else build[[1L]]
+})
 
 #' @rdname getQcInfo
 #' @examples

R/CtwasResult.R

@@ -14,6 +14,20 @@
 #' @include AllGenerics.R tupleSelectors.R
 NULL
 
+#' @title cTWAS Result Collection
+#' @description S4 collection of cTWAS runs keyed by the identity tuple
+#'   \code{(gwasStudy, study, context, method)}. Each row holds a
+#'   \code{\linkS4class{CtwasResultEntry}} payload -- fine-mapping posteriors,
+#'   the jointly-estimated group priors, and region metadata -- for one run.
+#' @details Unlike the QTL family, \code{trait} is not part of the key: a cTWAS
+#'   run is multi-gene, so genes live inside the payload. The optional
+#'   \code{jointStudies} / \code{jointContexts} columns tag rows born from a
+#'   multi-study or multi-context run and participate in the uniqueness key,
+#'   exactly as in the \code{\linkS4class{TwasWeights}} and fine-mapping
+#'   families.
+#' @seealso \code{\link{CtwasResult}} for the constructor and
+#'   \code{\link{ctwasPipeline}} for the pipeline that builds one.
+#' @export
 setClass("CtwasResult", contains = "DFrame", validity = function(object) {
     .validateCtwasResult(object)
 })

R/CtwasResultEntry.R

@@ -10,6 +10,23 @@
 #' @include AllGenerics.R
 NULL
 
+#' @title cTWAS Per-Run Payload
+#' @description Per-run cTWAS payload: the fine-mapping posterior table, the
+#'   full per-effect susie alpha table, the jointly-estimated group prior(s),
+#'   and per-region metadata. One entry sits in every row of a
+#'   \code{\linkS4class{CtwasResult}} collection.
+#' @slot finemap The per-gene (and, when SNPs are retained, per-SNP) posterior
+#'   summary table (\code{ctwas::finemap_regions} \code{finemap_res} shape),
+#'   or \code{NULL}.
+#' @slot susieAlpha The per-effect susie alpha table
+#'   (\code{ctwas::finemap_regions} \code{susie_alpha_res} shape) -- the
+#'   fuller cTWAS output retained so the raw run is reconstructable, or
+#'   \code{NULL}.
+#' @slot param The estimated \code{group_prior} / \code{group_prior_var} for
+#'   this run, or \code{NULL}.
+#' @slot regionInfo Per-region metadata, or \code{NULL}.
+#' @seealso \code{\link{CtwasResultEntry}} for the constructor.
+#' @export
 setClass(
     "CtwasResultEntry",
     representation(

R/QtlDataset.R

@@ -46,7 +46,6 @@ NULL
 #' @slot study Character (length 1). Study identifier; used in collection
 #'   classes to tag downstream \code{FineMappingResult} / \code{TwasWeights}
 #'   entries.
-#' @slot genotypes The genotype source for lazy access to dosages.
 #'   The \code{genotype} experiment's assay reads through this handle; the
 #'   extraction accessors read it directly, so that QC can be applied per
 #'   block.
@@ -78,7 +77,6 @@ setClass(
     contains = "MultiAssayExperiment",
     representation(
         study = "character",
-        genotypes = "GenotypeHandle",
         scaleResiduals = "logical",
         mafCutoff = "numeric",
         macCutoff = "numeric",
@@ -403,7 +401,6 @@ QtlDataset <- function(
         "QtlDataset",
         .qtlRestrictSamples(mae, keepSamples),
         study = as.character(study),
-        genotypes = handle,
         scaleResiduals = isTRUE(scaleResiduals),
         mafCutoff = as.numeric(mafCutoff),
         macCutoff = as.numeric(macCutoff),
@@ -575,19 +572,15 @@ setMethod("longForm", "QtlDataset", function(object, ..., genotype = FALSE) {
     x[, which(is_in(ids, as.character(keepSamples))), ]
 }
 
-# Replace the genotype handle, rebuilding the genotype experiment along with
-# it. The handle is deliberately held twice -- once as a slot, once inside
-# the assay's seed -- because that is what lets the assay read lazily. Moving
-# one without the other would leave the dosages describing a different panel
-# from the one the extraction accessors read, so nothing may set the slot
-# directly.
+# Replace the genotype handle by rebuilding the genotype experiment around
+# it. The handle lives in exactly one place -- the assay's seed -- so there
+# is no second copy to keep in step; getGenotypeHandle() reads it back.
 # @noRd
 .qtlWithGenotypeHandle <- function(x, handle) {
     exps <- MultiAssayExperiment::experiments(x)
     gCov <- .qtlColDataMatrix(exps[[.QTL_GENO_EXPERIMENT]])
     exps[[.QTL_GENO_EXPERIMENT]] <- .genotypeExperiment(handle, gCov)
     MultiAssayExperiment::experiments(x) <- exps
-    x@genotypes <- handle
     validObject(x)
     x
 }
@@ -673,7 +666,15 @@ setMethod("getScaleResiduals", "QtlDataset", function(x) x@scaleResiduals)
 
 #' @rdname getGenotypeHandle
 #' @keywords internal
-setMethod("getGenotypeHandle", "QtlDataset", function(x) x@genotypes)
+setMethod("getGenotypeHandle", "QtlDataset", function(x) {
+    # Derived, not stored. The handle already lives inside the genotype
+    # assay's seed -- that is what lets the dosages read lazily -- so a
+    # parallel slot was a second copy that had to be kept in step by hand.
+    # Reading it back removes the invariant instead of policing it.
+    .ldSketchHandle(
+        MultiAssayExperiment::experiments(x)[[.QTL_GENO_EXPERIMENT]]
+    )
+})
 
 #' @rdname qtlDatasetFilters
 #' @export
@@ -2069,8 +2070,9 @@ setMethod("show", "QtlDataset", function(object) {
         .trim = FALSE
     ))
     cat(glue("  {totalTraits} unique traits across contexts\n", .trim = FALSE))
+    gh <- getGenotypeHandle(object)
     cat(glue(
-        "  Genotypes: {object@genotypes@format} @ {object@genotypes@path}\n",
+        "  Genotypes: {getFormat(gh)} @ {getPath(gh)}\n",
         .trim = FALSE
     ))
     cat(glue(