Calculating geodesic distances in hemibrain dataset

[Pasha Davoudian]

Hello all,

I've been using the nat-verse predominantly for NBLAST (which has been extremely helpful), but have a few questions about using the anatomical tools on individual skeletons. I am interested in calculating the geodesic distance from specific postsynaptic locations to the soma of neuron in the hemibrain dataset, reading the nat paper it seems this can at least be done in principal (supp. fig 3).

I have been able to download the .swc files for neuron of interest, as well as the x,y,z information for the location of the synapses using the neuprint python API. However, I am having trouble understanding from a neuron-graph perspective how best to map these synapses onto the .swc file to calculate geodesic distance. Is such a task still possible in the hemibrain dataset? Or is it limited by the format of the .swc files for the neurons?

thanks,

Pasha 

[Gregory Jefferis]

Dear Pasha,

Quickly:

1. Use http://natverse.org/neuprintr/reference/neuprint_read_neurons.html to read neurons with connector information included. 

2. The connectors are stored inside each neuron. You can access them via the $connectors field

3. Each synapse is associated with the nearest skeleton node by Euclidean distance. Gory details here:

https://github.com/natverse/neuprintr/blob/7196c8bca15a15e34cb639d5e1fd37a190d38e01/R/neurons.R#L157

4. If the neuro has a known soma in the dataset, the skeleton will be rooted on the soma. 

5. You can use the ngraph representation to find the geodesic distance from the synapse nodes to the root. See the second half 

http://natverse.org/nat/articles/neurons-as-graph.html

Alex Bates May have an example handy. 

Best,

Greg. 

Sent from my iPhone

On 30 Apr 2020, at 20:00, Pasha Davoudian <pashada...@gmail.com> wrote:



--
You received this message because you are subscribed to the Google Groups "nat-user" group.
To unsubscribe from this group and stop receiving emails from it, send an email to nat-user+u...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/nat-user/285920a7-67fe-41ec-b9e5-7073e0bd77f4%40googlegroups.com.

[Pasha Davoudian]

Greg,

Thank you for the information, this was extremely helpful! I was able to calculate geodesic length more easily when using R to grab the hemibrain data. 

Relatedly, it seems I should just transition most of my hemibrain analysis to R since I am using the natvers already, I see that your group is working on a 'hemibrainR' package on GitHub. Is the current 0.10 version a stable release?

best,

Pasha 

To unsubscribe from this group and stop receiving emails from it, send an email to nat-...@googlegroups.com.

[Greg Jefferis]

Dear Pasha,

Thank you for the information, this was extremely helpful! I was able to calculate geodesic length more easily when using R to grab the hemibrain data. 

Great. I'm glad that worked out. Feel free to share a short example if you have one.

Relatedly, it seems I should just transition most of my hemibrain analysis to R since I am using the natvers already, I see that your group is working on a 'hemibrainR' package on GitHub. Is the current 0.10 version a stable release?

We are very actively developing the hemibrainr package (https://github.com/flyconnectome/hemibrainr) and it forms part of our own analysis of hemibrain data. It is designed to complement neuprintr by offering some functions that are specific to the hemibrain dataset rather than any neuprint database. It is currently incubating in the flyconnectome rather than the natverse GitHub organisation, but I would say at this point it should soon go from experimental to maturing. Rapid change still, but most of the functions that now exist are stable and you can certainly use if for serious analysis.

With best wishes,

Greg.

[Pasha Davoudian]

Greg,

Awesome, thanks again for your help (and for all the other tools as well). 

Below is my short example. Hopefully someone else will find it useful, even if it may not be the most well written code:

pre_synaptic_neuron = (hemibrain_id)

target_neuron_list = neuprint_read_neurons(hemibrain_id)

target_neuron = target_neuron_list[[1]]

  

 #get node IDs of each synapse and soma node ID

pre_to_post = filter(target_neuron$connectors, partner==pre_synaptic_neuron)

pre_to_post_IDs = (pre_to_post$treenode_id)

soma_ID = target_neuron$StartPoint

  

#convert neuron to ngraph object for calculation

g = as.ngraph(target_neuron, weights = TRUE)

  

#go through and get cable length from each synapse to the soma

all_cable_lengths <-c()

  for (i in pre_to_post_IDs) {

    synapse_to_soma <- igraph::distances(g, v=soma_ID, to=i, mode='all')

    cable_length = synapse_to_soma[1]

    all_cable_lengths <-c(all_cable_lengths, cable_length)

    i<-i+1

  }

  

#get average value and convert from nm to 1um

average_cable_length = (mean(all_cable_lengths)/(10^3))

[Greg Jefferis]

Thanks a lot, Pasha. Here is a sightlier fancier version starting from your code. The core difference is that the distances are computed in one go, which is rather fast. However it is necessary to check that there are no duplicates in the the requested start nodes because the graph distances doesn't like that.

I will probably add something like this to nat soon.

All the best,

Greg.

#' Find geodesic distance along neuronal tree from the soma to every synapse

#'

#' @description Find the geodesic distance (i.e. path length along) the neuron

#'   from the the root (usually soma) to all synapses or other selected nodes.

#' @param x A \code{\link{neuron}}

#' @param root The index of the node to use as the root for the calculation.

#'   Defaults to the soma.

#' @param nodes Optional argument specifying to which nodes you want to measure

#'   distances. By default this will be all connectors in the neuron, but you

#'   can also specify all nodes, or pass in a any data frame whose first columm

#'   corresponds to (a subset) of the node identifiers of the neuron.

#' @param rval What kind of object to return

#' @return A vector of distances, a data.frame or a \code{\link{neuron}} object

#'   with an extra column in the connectors data.frame entitled

#'   \code{rootdistance}.

#' @export

#'

#' @examples

#' # nb pick first neuron out of neuron list returned by neuprint_read_neurons

#' dl4=neuprint_read_neurons('DL4')[[1]]

#' hist(geodesic_distances(dl4), br=100)

#' df=geodesic_distances(dl4, rval='dataframe')

#' # the correction factor of 125 scales to microns

#' plot(dl4/125, col='black', WithNodes=F)

#' points(xyzmatrix(df)[,c(1,2)]/125, col=topo.colors(10)[cut(df$rootdistance, 10)], pch=19)

#' # plot the soma position

#' points(xyzmatrix(dl4)[rootpoints(dl4),,drop=FALSE]/125, col='red', pch=19)

geodesic_distances <- function(x, root=rootpoints(x),

                               nodes=c("synapses", "nodes"),

                               rval=c("distance", "dataframe", 'neuron')) {

  rval=match.arg(rval)

  if(is.character(nodes)) {

    nodes=match.arg(nodes)

    nodes <- if(nodes=="synapses") x$connectors else x$d

  }

  if(is.null(nodes))

    return(NULL)

  if(!is.data.frame(nodes))

    stop("I cannot understand the nodes argument. See docs for help!")

  g=as.ngraph(x, weights=TRUE)

  # Some nodes may make multiple connections, igraph::distances won't cope

  havedups=any(duplicated(nodes[[1]]))

  tns <- if(havedups) unique(nodes[[1]]) else nodes[[1]]

  # nb converting the to nodes to character handles both the cases where they are

  # 64 bit ids (typical of catmaid neurons) and the case where they are indices

  # into the point array (typical of neuprint neurons)

  dd=igraph::distances(g, v=root, to=as.character(tns), mode='all')

  if(havedups) {

    # re-expand to match the incoming data.frame of connectors

    dd=dd[match(nodes[[1]], tns)]

  }

  if(rval=='distance')

    return(dd)

  nodes$rootdistance=dd

  if(rval=='dataframe') nodes else {x$connectors=nodes; x}

[Pasha Davoudian]

Greg,

Thanks, the speed improvement is notable especially when scaling this up to measure several neurons. I think it would be great addition to nat.

One last question: I noticed you used a correction factor of 125 to convert to um, while I used 10^3. Is there a reason why I should be using 125 instead as well? Maybe I am missing something with how ngraph() converts measurements when weights = TRUE is passed?

best,

Pasha

[Gregory Jefferis]

The factor of / 125 is the product of 8 / 1000. This is because neurons in the hemibrain are in a raw EM coordinate space. The EM voxels are 8nm, so multiply by 8. The / 1000 is then  to convert from nm to microns. All the best, Greg.

Sent from my iPhone

--
You received this message because you are subscribed to the Google Groups "nat-user" group.

To unsubscribe from this group and stop receiving emails from it, send an email to nat-user+u...@googlegroups.com.
To view this discussion on the web, visit https://groups.google.com/d/msgid/nat-user/82ea0cfe-322c-4eca-90e6-69f79aea468f%40googlegroups.com.