|
| 1 | +### ========================================================================= |
| 2 | +### Download data from OGRDB thru the OGRDB API |
| 3 | +### ------------------------------------------------------------------------- |
| 4 | +### |
| 5 | +### Nothing in this file is exported. |
| 6 | +### |
| 7 | + |
| 8 | + |
| 9 | +### See short introduction to OGRDB REST API at |
| 10 | +### https://wordpress.vdjbase.org/index.php/ogrdb_news/downloading-germline-sets-from-the-command-line-or-api/ |
| 11 | +.OGRDB_API_URL <- paste0(OGRDB_URL, "api") |
| 12 | + |
| 13 | +.encode_OGRDB_API_query <- function(query) |
| 14 | +{ |
| 15 | + stopifnot(is.character(query)) |
| 16 | + if (length(query) != 0L) |
| 17 | + stopifnot(!is.null(names(query))) |
| 18 | + vapply(query, encode_OGRDB_URL_component, character(1)) |
| 19 | +} |
| 20 | + |
| 21 | + |
| 22 | +### - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - |
| 23 | +### fetch_germline_set_from_OGRDB() |
| 24 | +### |
| 25 | +### No longer used. Superseded by download_OGRDB_germline_set() defined |
| 26 | +### in R/OGRDB-utils.R. |
| 27 | +### |
| 28 | +### Retrieves: |
| 29 | +### /germline/set/ |
| 30 | +### {species}/{set_name}/{release_version}/{format} |
| 31 | +### or |
| 32 | +### {species}/{species_subgroup}/{set_name}/{release_version}/{format} |
| 33 | +### as documented at https://ogrdb.airr-community.org/api/ |
| 34 | +### |
| 35 | +### Typical usage: |
| 36 | +### |
| 37 | +### fetch_germline_set_from_OGRDB("Human", set_name="IGH_VDJ") |
| 38 | +### fetch_germline_set_from_OGRDB("Mouse", species_subgroup="C57BL/6", |
| 39 | +### set_name="C57BL/6 IGH") |
| 40 | +### |
| 41 | +### Note that passing "Homo_sapiens" or "Mus musculus" instead of "Human" |
| 42 | +### or "Mouse" above also works and produces the same results. |
| 43 | +### Returns a character vector containing the nucleotide sequences in FASTA |
| 44 | +### format, except when 'format' is set to "airr" or "airr_ex". |
| 45 | +### |
| 46 | +### About the "airr" and "airr_ex" formats |
| 47 | +### -------------------------------------- |
| 48 | +### |
| 49 | +### When 'format' is set to "airr" or "airr_ex", the function returns a |
| 50 | +### named list containing a bunch of information about the germline set |
| 51 | +### in adition to the nucleotide sequences. |
| 52 | +### If 'res' is the results obtained with 'format="airr"' then all the |
| 53 | +### sequences in the germline set are described in 'res$allele_descriptions' |
| 54 | +### which is a data.frame with 1 row per sequence and dozens of columns. |
| 55 | +### However, for some mysterious reason, this data.frame has more rows than |
| 56 | +### the number of sequences returned when 'format' is not "airr" or "airr_ex". |
| 57 | +### For example, 'res$allele_descriptions' has 245 rows for Human/IGH_VDJ |
| 58 | +### but using 'format="ungapped"' returns only 236 sequences! |
| 59 | +### |
| 60 | +### Assuming 'alleles' is the data.frame obtained by removing the extra |
| 61 | +### sequences from 'res$allele_descriptions', then the seq ids, species, |
| 62 | +### species subgroups, loci, nucleotide sequences (ungapped and gapped), |
| 63 | +### and region types, can be extracted with: |
| 64 | +### |
| 65 | +### - seqids <- alleles$label |
| 66 | +### |
| 67 | +### - species <- alleles$species$label |
| 68 | +### |
| 69 | +### - species_subgroups <- alleles$species_subgroup |
| 70 | +### |
| 71 | +### - loci <- alleles$locus |
| 72 | +### |
| 73 | +### - region_types <- alleles$sequence_type |
| 74 | +### |
| 75 | +### - ungapped_seqs <- alleles$coding_sequence |
| 76 | +### Note that for sequences of type V, this should match: |
| 77 | +### sapply(alleles$v_gene_delineations, |
| 78 | +### function(x) if (is.null(x)) NA_character_ else x$unaligned_sequence) |
| 79 | +### |
| 80 | +### - gapped_seqs <- sapply(alleles$v_gene_delineations, |
| 81 | +### function(x) if (is.null(x)) NA_character_ else x$aligned_sequence) |
| 82 | +### Note that only sequences of type V can have gaps. |
| 83 | + |
| 84 | +fetch_germline_set_from_OGRDB <- |
| 85 | + function(species, species_subgroup=NULL, set_name, |
| 86 | + release_version="published", |
| 87 | + format=c("gapped", "ungapped", "airr", |
| 88 | + "gapped_ex", "ungapped_ex", "airr_ex")) |
| 89 | +{ |
| 90 | + stopifnot(isSingleNonWhiteString(species), |
| 91 | + is.null(species_subgroup) || |
| 92 | + isSingleNonWhiteString(species_subgroup), |
| 93 | + isSingleNonWhiteString(set_name), |
| 94 | + isSingleNonWhiteString(release_version)) |
| 95 | + format <- match.arg(format) |
| 96 | + |
| 97 | + query <- c(species=species, species_subgroup=species_subgroup, |
| 98 | + set_name=set_name, release_version=release_version, |
| 99 | + format=format) |
| 100 | + query <- .encode_OGRDB_API_query(query) |
| 101 | + url <- paste(c(.OGRDB_API_URL, "germline/set", query), collapse="/") |
| 102 | + content <- getUrlContent(url, type="text", encoding="UTF-8") |
| 103 | + ### Is it possible that 'content' will be NA? See R/REST_API.R in |
| 104 | + ### the UCSC.utils package for more info. |
| 105 | + stopifnot(isSingleString(content)) |
| 106 | + |
| 107 | + if (!(format %in% c("airr", "airr_ex"))) { |
| 108 | + fasta_lines <- strsplit(content, split="\n", fixed=TRUE)[[1L]] |
| 109 | + return(fasta_lines) # nucleotide sequences in FASTA format |
| 110 | + } |
| 111 | + |
| 112 | + parsed_json <- fromJSON(content, simplifyDataFrame=FALSE) |
| 113 | + ## Sanity checks. |
| 114 | + stopifnot(is.list(parsed_json), |
| 115 | + length(parsed_json) == 1L, |
| 116 | + identical(names(parsed_json), "GermlineSet"), |
| 117 | + is.list(parsed_json[[1L]]), |
| 118 | + length(parsed_json[[1L]]) == 1L) |
| 119 | + ans <- parsed_json[[1L]][[1L]] |
| 120 | + stopifnot(is.list(ans), !is.null(names(ans)), |
| 121 | + "allele_descriptions" %in% names(ans)) |
| 122 | + ans$allele_descriptions <- |
| 123 | + jsonlite:::simplifyDataFrame(ans$allele_descriptions, |
| 124 | + flatten=FALSE, simplifyMatrix=TRUE) |
| 125 | + ans # named list |
| 126 | +} |
| 127 | + |
0 commit comments