diff --git a/DESCRIPTION b/DESCRIPTION index 00c09073..0b751a56 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -34,7 +34,8 @@ Imports: nat, readobj, memoise, - dplyr + dplyr, + bit64 Remotes: natverse/drvid Encoding: UTF-8 diff --git a/R/connectivity.R b/R/connectivity.R index f3be16c0..5698ffdb 100644 --- a/R/connectivity.R +++ b/R/connectivity.R @@ -16,9 +16,7 @@ neuprint_get_adjacency_matrix <- function(bodyids, dataset = NULL, all_segments "WHERE n.bodyId IN input AND m.bodyId IN input", "RETURN n.bodyId AS upstream, m.bodyId AS downstream, c.weight AS weight, n.%s AS upName, m.%s AS downName" ), - jsonlite::toJSON(as.numeric(unique(unlist( - bodyids - )))), + id2json(bodyids), all_segments.json, namefield, namefield @@ -72,7 +70,7 @@ neuprint_connection_table <- function(bodyids, prepost = c("PRE","POST"), roi = "UNWIND %s AS k", "RETURN a.bodyId AS %s, b.bodyId AS %s, k AS roi,", "apoc.convert.fromJsonMap(c.roiInfo)[k].post AS weight"), - jsonlite::toJSON(unique(as.numeric(unlist(bodyids)))), + id2json(bodyids), all_segments.json, all_segments.json, ifelse(prepost=="POST","a","b"), @@ -123,7 +121,7 @@ neuprint_common_connectivity <- function(bodyids, statuses = NULL, all_segments = ifelse(all_segments,"true","false") Payload = noquote(sprintf('{"dataset":"%s","neuron_ids":%s,"statuses":%s,"find_inputs":%s,"all_segments":%s}', dataset, - jsonlite::toJSON(unique(as.numeric(unlist(bodyids)))), + id2json(bodyids), ifelse(is.null(statuses),jsonlite::toJSON(list()),jsonlite::toJSON(statuses)), find_inputs, all_segments)) @@ -184,7 +182,7 @@ neuprint_simple_connectivity <- function(bodyids, } Payload = noquote(sprintf('{"dataset":"%s","neuron_ids":%s,"find_inputs":%s}', dataset, - jsonlite::toJSON(as.numeric(unique(unlist(bodyids)))), + id2json(bodyids), find_inputs)) class(Payload) = "json" simp.conn = neuprint_fetch(path = 'api/npexplorer/commonconnectivity', body = Payload, conn = conn, simplifyVector = TRUE, ...) @@ -259,8 +257,8 @@ neuprint_get_paths <- function(body_pre, body_post, n, weightT=5, roi=NULL, "ALL (x in relationships(p) WHERE x.weight >= %s %s)", "RETURN length(p) AS `length(path)`,[n in nodes(p) | [n.bodyId, n.instance, n.type]] AS path,[x in relationships(p) | x.weight] AS weights" ), - jsonlite::toJSON(unique(as.numeric(unlist(body_pre)))), - jsonlite::toJSON(unique(as.numeric(unlist(body_post)))), + id2json(body_pre), + id2json(body_post), all_segments.json, n[1]-1, n[2], @@ -325,8 +323,8 @@ neuprint_get_shortest_paths <- function(body_pre,body_post,weightT=5,roi=NULL,da "ALL (x in relationships(p) WHERE x.weight >= %s %s)", "RETURN length(p) AS `length(path)`,[n in nodes(p) | [n.bodyId, n.instance, n.type]] AS path,[x in relationships(p) | x.weight] AS weights" ), - jsonlite::toJSON(unique(as.numeric(unlist(body_pre)))), - jsonlite::toJSON(unique(as.numeric(unlist(body_post)))), + id2json(body_pre), + id2json(body_post), all_segments.json, all_segments.json, weightT, diff --git a/R/dump.R b/R/dump.R index 67e0d20a..9be4edb1 100644 --- a/R/dump.R +++ b/R/dump.R @@ -32,7 +32,7 @@ neuprint_dump <- function(dir, bodyids = NULL, roi = NULL, preprocess = NULL, co inroi = neuprint_bodies_in_ROI( roi = roi, dataset = dataset, conn=conn, ...) inroi = subset(inroi, voxels>voxel.thresh) - bodyids = as.numeric(unique(c(unlist(inroi$bodyid),bodyids))) + bodyids = unique(c(id2bit64(inroi$bodyid),id2bit64(bodyids))) } # Fetch neuron data message("Reading neurons from ", conn$server, " for dataset: ", dataset) diff --git a/R/ids.R b/R/ids.R new file mode 100644 index 00000000..a0db6302 --- /dev/null +++ b/R/ids.R @@ -0,0 +1,29 @@ +# private function to convert 1 or more 64 bit ids (e.g. body ids) to JSON +id2json <- function(x, uniqueids=TRUE, ...) { + bx=id2bit64(x) + if(isTRUE(uniqueids)) bx=unique(bx) + jsonlite::toJSON(bx, ...) +} + +# ... this bit might also be useful +id2bit64 <- function(x) { + x=unlist(x, use.names = FALSE) + if(is.factor(x)) { + x=as.character(x) + } + if(isTRUE(is.character(x)) || bit64::is.integer64(x)) { + # do nothing + } else if(is.integer(x)){ + BIGGESTINT=2147483647L + if(any(x>BIGGESTINT)) + stop("Some ids cannot be exactly represented! Please use character or bit64!") + } else if(is.double(x)) { + # biggest int that can be represented as double 2^(mantissa bits + 1) + BIGGESTINT=2^53+1 + if(any(x>BIGGESTINT)) + stop("Some ids cannot be exactly represented! Please use character or bit64!") + } else { + stop("Unexpected data type for id. Use character, bit64, or numeric!") + } + bit64::as.integer64(x) +} diff --git a/R/name.R b/R/name.R index bbecb0d8..06a74b84 100644 --- a/R/name.R +++ b/R/name.R @@ -13,7 +13,7 @@ neuprint_get_neuron_names <- function(bodyids, dataset = NULL, all_segments = FA all_segments.json = ifelse(all_segments,"Segment","Neuron") cypher = sprintf("WITH %s AS bodyIds UNWIND bodyIds AS bodyId MATCH (n:`%s`) WHERE n.bodyId=bodyId RETURN n.instance AS name", - jsonlite::toJSON(as.numeric(unlist(bodyids))), + id2json(bodyids), all_segments.json) nc = neuprint_fetch_custom(cypher=cypher, conn = conn, dataset = dataset, ...) @@ -41,7 +41,7 @@ neuprint_get_meta <- function(bodyids, dataset = NULL, all_segments = FALSE, con "MATCH (n:`%s`) WHERE n.bodyId=bodyId", "RETURN n.bodyId AS bodyid, n.%s AS name, n.type AS type, n.status AS status, n.size AS voxels, n.pre AS pre, n.post AS post, n.cropped AS cropped" ), - jsonlite::toJSON(as.numeric(unlist(bodyids))), + id2json(bodyids), all_segments, neuprint_name_field(conn) ) @@ -64,7 +64,7 @@ neuprint_get_roiInfo <- function(bodyids, dataset = NULL, all_segments = FALSE, all_segments = ifelse(all_segments,"Segment","Neuron") cypher = sprintf( "WITH %s AS bodyIds UNWIND bodyIds AS bodyId MATCH (n:`%s`) WHERE n.bodyId=bodyId RETURN n.bodyId AS bodyid, n.roiInfo AS roiInfo", - jsonlite::toJSON(as.numeric(unlist(bodyids))), + id2json(bodyids), all_segments ) nc = neuprint_fetch_custom(cypher=cypher, dataset = dataset, conn = conn, ...) diff --git a/R/neurons.R b/R/neurons.R index 408b8be4..e803c61d 100644 --- a/R/neurons.R +++ b/R/neurons.R @@ -44,8 +44,8 @@ neuprint_read_neurons <- function(bodyids, dataset = NULL, resample = FALSE, conn = NULL, - OmitFailures = TRUE, ...){ - neurons = nat::nlapply(unique(bodyids),function(bodyid) + OmitFailures = TRUE, ...) { + neurons = nat::nlapply(id2bit64(bodyids),function(bodyid) neuprint_read_neuron(bodyid=bodyid, nat=nat, drvid=drvid, @@ -64,7 +64,10 @@ neuprint_read_neurons <- function(bodyids, } names(neurons) = unlist(sapply(neurons,function(n) n$bodyid)) if(meta){ - attr(neurons,"df") = neuprint_get_meta(bodyids = as.numeric(names(neurons)), dataset = dataset, all_segments = all_segments, conn = conn, ...) + attr(neurons,"df") = neuprint_get_meta(bodyids = names(neurons), + dataset = dataset, + all_segments = all_segments, + conn = conn, ...) }else{ attr(neurons,"df") = data.frame(bodyid=names(neurons)) } @@ -92,7 +95,8 @@ neuprint_read_neuron <- function(bodyid, n = drvid::read.neuron.dvid(bodyid) d = n$d }else{ - n = neuprint_read_neuron_simple(as.numeric(bodyid),dataset=dataset,conn = conn, heal = FALSE,...) + n = neuprint_read_neuron_simple(id2bit64(bodyid), dataset=dataset, + conn = conn, heal = FALSE,...) } if(heal){ n = suppressWarnings( heal_skeleton(x = n) ) @@ -201,8 +205,11 @@ heal_skeleton <- function(x, ...){ #' @rdname neuprint_read_neurons neuprint_read_neuron_simple <- function(bodyid, dataset=NULL, conn=NULL, heal=TRUE, ...) { dataset <- check_dataset(dataset) + bodyid=as.character(id2bit64(bodyid)) if(length(bodyid)>1) { - return(nat::nlapply(bodyid, neuprint_read_neuron_simple, dataset=dataset, conn=conn, ...)) + fakenl=structure(bodyid, .Names=bodyid) + nl=nat::nlapply(fakenl, neuprint_read_neuron_simple, dataset=dataset, conn=conn, ...) + return(nl) } path=file.path("api/skeletons/skeleton", dataset, bodyid) res=neuprint_fetch(path, conn=conn, simplifyVector = TRUE, include_headers = FALSE, ...) diff --git a/R/soma.R b/R/soma.R index 25804e04..f6a4f0bd 100644 --- a/R/soma.R +++ b/R/soma.R @@ -14,7 +14,7 @@ neuprint_locate_soma <- function(bodyids, dataset = NULL, all_segments = TRUE, c conn=neuprint_login(conn) all_segments.json = ifelse(all_segments,"Segment","Neuron") cypher = sprintf("WITH %s AS bodyIds UNWIND bodyIds AS bodyId MATCH (n:`%s`) WHERE n.bodyId=bodyId RETURN n.bodyId AS bodyId, n.somaLocation AS soma", - jsonlite::toJSON(as.numeric(unique(bodyids))), + id2json(bodyids), all_segments.json) nc = neuprint_fetch_custom(cypher=cypher, conn = conn, dataset = dataset, ...) if(length(nc$data)==0){ diff --git a/R/synapses.R b/R/synapses.R index 4604b2ab..ca4b0af9 100644 --- a/R/synapses.R +++ b/R/synapses.R @@ -45,7 +45,7 @@ neuprint_get_synapses <- function(bodyids, roi = NULL, progress = FALSE, dataset if(progress){ d = do.call(rbind, pbapply::pblapply(bodyids, function(bi) tryCatch(neuprint_get_synapses( - bodyids = as.numeric(bi), + bodyids = bi, roi = roi, progress = FALSE, dataset = dataset, @@ -61,7 +61,7 @@ neuprint_get_synapses <- function(bodyids, roi = NULL, progress = FALSE, dataset "RETURN DISTINCT id(s) AS connector_id,", "s.type AS prepost, s.location.x AS x ,s.location.y AS y, s.location.z AS z,", "s.confidence AS confidence, a.bodyId AS bodyid, b.bodyId AS partner"), - jsonlite::toJSON(as.numeric(unlist(bodyids))), + id2json(bodyids), prefixed_seg, prefixed_seg, roi) @@ -71,7 +71,7 @@ neuprint_get_synapses <- function(bodyids, roi = NULL, progress = FALSE, dataset "RETURN DISTINCT id(s) AS connector_id,", "s.type AS prepost, s.location.x AS x ,s.location.y AS y, s.location.z AS z,", "s.confidence AS confidence, a.bodyId AS bodyid, b.bodyId AS partner"), - jsonlite::toJSON(as.numeric(unlist(bodyids))), + id2json(bodyids), prefixed_seg, prefixed_seg, roi) diff --git a/tests/testthat/test-ids.R b/tests/testthat/test-ids.R new file mode 100644 index 00000000..7be753c7 --- /dev/null +++ b/tests/testthat/test-ids.R @@ -0,0 +1,28 @@ +test_that("id conversion works", { + bigid="9223372036854775806" + big.json="[9223372036854775806]" + big.json2="[9223372036854775806,9223372036854775806]" + + medid="223372036854775806" + med.json="[223372036854775806]" + + expect_equal(as.character(id2json(bigid)), + big.json) + expect_equal(as.character(id2json(bit64::as.integer64(bigid))), + big.json) + expect_error(id2json(as.numeric(bigid))) + expect_error(id2json(NA)) + expect_equal(as.character(id2json(medid)), + med.json) + expect_equal(as.character(id2json(bit64::as.integer64(medid))), + med.json) + expect_equal(as.character(id2json(as.factor(medid))), + med.json) + + expect_equal(as.character(id2json(c(bigid, bigid))), + big.json) + expect_equal(as.character(id2json(c(bigid, bigid), uniqueids=FALSE)), + big.json2) + expect_equal(as.character(id2json(list(bigid, bigid))), + big.json) +})