Added util funcs & getVariantsFromRegion

daynefiler · Dec 3, 2021 · a71fdc7 · a71fdc7
1 parent 18b0ae1
commit a71fdc7
Show file tree

Hide file tree

Showing 10 changed files with 313 additions and 1 deletion.
diff --git a/DESCRIPTION b/DESCRIPTION
@@ -16,7 +16,10 @@ LazyData: true
 Roxygen: list(markdown = TRUE)
 RoxygenNote: 7.1.2
 Imports: 
-    data.table
+    data.table,
+    ghql,
+    glue,
+    jsonlite
 Suggests: 
     spelling
 Language: en-US

diff --git a/NAMESPACE b/NAMESPACE
@@ -1,5 +1,15 @@
 # Generated by roxygen2: do not edit by hand
 
+export(DatasetIds)
+export(ReferenceGenomeIds)
+export(apiUrl)
+export(compatibleGenomeDataset)
+export(getDatasets)
+export(getGenomes)
+export(getVariantsFromRegion)
+export(validDatasets)
+export(validGenomes)
+import(ghql)
 importFrom(data.table,":=")
 importFrom(data.table,.BY)
 importFrom(data.table,.EACHI)
@@ -9,3 +19,6 @@ importFrom(data.table,.N)
 importFrom(data.table,.NGRP)
 importFrom(data.table,.SD)
 importFrom(data.table,data.table)
+importFrom(glue,glue)
+importFrom(glue,glue_collapse)
+importFrom(jsonlite,fromJSON)
diff --git a/R/getVariantsFromRegion.R b/R/getVariantsFromRegion.R
@@ -0,0 +1,47 @@
+#' @title Get variants from region
+#' @inheritParams pkgParams
+#' @description
+#' Returns the variant_id, ref, alt, rsids fields for all variants in the given
+#' region(s).
+#'
+#' @return List of data.frames for each region given
+#'
+#' @examples
+#' \dontrun{
+#' getVariantsFromRegion(genomes = 'GRCh37',
+#'                       chroms = c('1'),
+#'                       starts = c(89388944))
+#' }
+#'
+#' @import ghql
+#' @importFrom glue glue glue_collapse
+#' @importFrom jsonlite fromJSON
+#' @export
+
+getVariantsFromRegion <- function(genomes, chroms, starts, stops = starts) {
+  stopifnot(validGenomes(genomes))
+  gmCon <- GraphqlClient$new(url = apiUrl())
+  datasets <- getDatasets(genomes)
+  tmp <-
+  '
+  {genomes}_{chroms}_{starts}_{stops}: region(chrom: "{chroms}",
+           start: {starts},
+           stop: {stops},
+           reference_genome: {genomes}) {{
+      variants(dataset: {datasets}) {{
+        variant_id
+        ref
+        alt
+        rsids
+      }}
+    }}
+  '
+  qryBody <- glue_collapse(glue(tmp), sep = "\n")
+  qry <- Query$new()$query('getRegion',
+                           glue('query getRegion {{ {qryBody} }}'))
+  tryres <- try(jsn <- gmCon$exec(qry$getRegion), silent = TRUE)
+  if (is(tryres, 'try-error')) stop(qfailmessage)
+  res <- lapply(fromJSON(jsn)$data, "[[", 1)
+  res
+}
+
diff --git a/R/gnomadR-package.R b/R/gnomadR-package.R
@@ -13,3 +13,18 @@
 #' @importFrom data.table data.table
 ## usethis namespace: end
 NULL
+
+#' @name pkgParams
+#' @title Common function parameters
+#' @description Common function parameters
+#' @param dataset Character scalar, gnomAD DatasetId
+#' @param datasets Character, gnomAD DatasetId(s)
+#' @param genome Character scalar, gnomAD ReferenceGenomeId
+#' @param genomes Character, gnomAD ReferenceGenomeId(s)
+#' @param chroms Character, the chromosome(s) ('1', '2', ..., 'X', 'Y')
+#' @param chrom Character scalar, the chromosome ('1', '2', ..., 'X', 'Y')
+#' @param starts Integer, starting position(s)
+#' @param start Integer scalar, starting position
+#' @param stop Integer scalar, ending position
+#'
+NULL
diff --git a/R/utils.R b/R/utils.R
@@ -0,0 +1,96 @@
+#' @title Get gnomAD API URL
+#' @description
+#' Returns the gnomAD API URL as a character. Checks if the option
+#' 'gnomadr.apiurl' has been set to allow for hypothetical parameterization
+#' in the future
+#' @return Character with gnomAD API URL
+#' @export
+
+apiUrl <- function() {
+  getOption("gnomadr.apiurl", "https://gnomad.broadinstitute.org/api")
+}
+
+qfailmessage <-
+  c("Query failed. This is likely due to the size of the query. At this time ",
+    "the gnomAD API has very low query limits. This package is designed to ",
+    "accomodate large queries in the event gnomAD increases query limits in ",
+    "the future.")
+
+#' @name gnomad-ids
+#' @title gnomAD IDs and associated helper functions
+#' @title Validate/associate reference genome build ID and gnomAD data set ID
+#' @description
+#' Description here.
+
+NULL
+
+#' @rdname gnomad-ids
+#' @importFrom glue glue
+#' @export
+
+ReferenceGenomeIds <- c("GRCh38", "GRCh37")
+attr(ReferenceGenomeIds, "defDatasetOpt") <-
+  glue("gnomadr.{ReferenceGenomeIds}.dataset")
+attr(ReferenceGenomeIds, "defDataset") <- c("gnomad_r3", "gnomad_r2_1")
+
+#' @rdname gnomad-ids
+#' @export
+
+DatasetIds <- c("gnomad_r3",
+                "gnomad_r3_controls_and_biobanks",
+                "gnomad_r3_non_cancer",
+                "gnomad_r3_non_neuro",
+                "gnomad_r3_non_topmed",
+                "gnomad_r3_non_v2",
+                "gnomad_r2_1",
+                "gnomad_r2_1_controls",
+                "gnomad_r2_1_non_neuro",
+                "gnomad_r2_1_non_cancer",
+                "gnomad_r2_1_non_topmed",
+                "exac")
+attr(DatasetIds, "associatedGenome") <- rep(ReferenceGenomeIds, each = 6)
+
+#' @rdname gnomad-ids
+#' @inheritParams pkgParams
+#' @export
+
+validGenomes <- function(genomes) all(genomes %in% ReferenceGenomeIds)
+
+#' @rdname gnomad-ids
+#' @inheritParams pkgParams
+#' @export
+
+validDatasets <- function(datasets) all(datasets %in% DatasetIds)
+
+#' @rdname gnomad-ids
+#' @inheritParams pkgParams
+#' @export
+
+compatibleGenomeDataset <- function(datasets, genomes) {
+  stopifnot(validDatasets(datasets))
+  stopifnot(validGenomes(genomes))
+  identical(genomes, getGenomes(datasets))
+}
+
+#' @rdname gnomad-ids
+#' @inheritParams pkgParams
+#' @export
+
+getGenomes <- function(datasets) {
+  stopifnot(validDatasets(datasets))
+  attr(DatasetIds, "associatedGenome")[match(datasets, DatasetIds)]
+}
+
+#' @rdname gnomad-ids
+#' @inheritParams pkgParams
+#' @export
+
+getDatasets <- function(genomes) {
+  stopifnot(validGenomes(genomes))
+  rgi <- match(genomes, ReferenceGenomeIds)
+  dsopt <- attr(ReferenceGenomeIds, "defDatasetOpt")[rgi]
+  dsdef <- attr(ReferenceGenomeIds, "defDataset")[rgi]
+  ds <- mapply(getOption, x = dsopt, default = dsdef, USE.NAMES = FALSE)
+  stopifnot(compatibleGenomeDataset(datasets = ds, genomes = genomes))
+  ds
+}
diff --git a/man/apiUrl.Rd b/man/apiUrl.Rd
diff --git a/man/getVariantsFromRegion.Rd b/man/getVariantsFromRegion.Rd
diff --git a/man/gnomad-ids.Rd b/man/gnomad-ids.Rd
diff --git a/man/gnomadR-package.Rd b/man/gnomadR-package.Rd
diff --git a/man/pkgParams.Rd b/man/pkgParams.Rd