Skip to content

Commit

Permalink
Merge pull request #403 from natverse/prune
Browse files Browse the repository at this point in the history
Prune
  • Loading branch information
jefferis authored Jan 31, 2020
2 parents a4ed90f + 51a6958 commit 7cefd30
Show file tree
Hide file tree
Showing 12 changed files with 535 additions and 7 deletions.
1 change: 1 addition & 0 deletions DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -64,6 +64,7 @@ Collate:
'graph-nodes.R'
'hxsurf.R'
'im3d.R'
'interactive.R'
'nat-data.R'
'nat-package.R'
'ndigest.R'
Expand Down
9 changes: 9 additions & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -125,6 +125,10 @@ S3method(prune,default)
S3method(prune,dotprops)
S3method(prune,neuron)
S3method(prune,neuronlist)
S3method(prune_in_volume,neuron)
S3method(prune_in_volume,neuronlist)
S3method(prune_online,neuron)
S3method(prune_online,neuronlist)
S3method(remotesync,default)
S3method(remotesync,neuronlistfh)
S3method(resample,neuron)
Expand Down Expand Up @@ -213,6 +217,7 @@ export(cmtkparams2affmat)
export(cmtkreg)
export(cmtkreglist)
export(coord2ind)
export(correct_root)
export(dotprops)
export(endpoints)
export(fileformats)
Expand Down Expand Up @@ -241,6 +246,7 @@ export(is.neuronlist)
export(is.nrrd)
export(is.swc)
export(is.vaa3draw)
export(make_model)
export(makeboundingbox)
export(mask)
export(materials)
Expand Down Expand Up @@ -269,6 +275,8 @@ export(progress_natprogress)
export(projection)
export(prune)
export(prune_edges)
export(prune_in_volume)
export(prune_online)
export(prune_strahler)
export(prune_vertices)
export(read.amiramesh)
Expand Down Expand Up @@ -298,6 +306,7 @@ export(seglengths)
export(seglist)
export(seglist2swc)
export(segmentgraph)
export(select_points)
export(setdiff)
export(simplify_reglist)
export(smooth_neuron)
Expand Down
2 changes: 1 addition & 1 deletion R/dotprops.R
Original file line number Diff line number Diff line change
Expand Up @@ -535,7 +535,7 @@ subset.dotprops<-function(x, subset, invert=FALSE, ...){
#' @param ... Additional arguments for methods (eventually passed to
#' \code{prune.default})
#' @seealso \code{\link{prune_strahler}}, \code{\link{spine}},
#' \code{\link{prune_vertices}}
#' \code{\link{prune_vertices}}, \code{\link{subset.neuron}}
#' @examples
#' ## prune single neurons
#' \donttest{
Expand Down
252 changes: 252 additions & 0 deletions R/interactive.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,252 @@
### Functions for interactively working with neurons ###

#' Prune a neuron interactively in an rgl window
#'
#' @description Remove points from a neuron, keeping the root node intact.
#' @details The neuron is plotted initially with all nodes selected (and shown
#' with black points). You can interactively select points to remove (they
#' will now be plotted in red). You can also add points back again (they will
#' return to black). When you are finished, press [e] to exit and then indicate
#' that you have finished (yes).
#' @param ... Additional methods passed to \code{prune_vertices}
#' @inheritParams prune
#' @return A pruned neuron/neuronlist object
#' @examples
#' \dontrun{
#' ## Interactively choose which bit of the neuron you wish to keep
#' pruned.as.you.like.it = prune_online(Cell07PNs[1:2])
#' }
#' @seealso \code{\link{as.neuron.ngraph}}, \code{\link{subset.neuron}},
#' \code{\link{prune.neuron}}
#' @export
#' @rdname prune_online
prune_online <-function(x, ...) UseMethod("prune_online")
#' @export
#' @rdname prune_online
prune_online.neuron <- function(x, ...){
continue = "no"
ids=integer()
while(!continue%in%c("y","yes")){
selected = select_points(xyzmatrix(x), clear_plot_on_exit = TRUE)
v = match(data.frame(t(selected)), data.frame(t(xyzmatrix(x))))
neuron = prune_vertices(x, verticestoprune=v,invert=TRUE)
pop3d(id=unlist(ids))
ids=rgl::plot3d(neuron, col ="black")
continue = readline("Finished with this neuron? yes/no ")
pop3d(id=unlist(ids))
}
neuron
}
#' @export
#' @rdname prune_online
prune_online.neuronlist <- function(x, ...){
nlapply(x,prune_online.neuron)
}


#' Interactively re-root neurons (usually to their soma)
#'
#' @description Cycle through and manually re-root neurons using an rgl window.
#' This will typically be used for manual identification of the soma of a
#' neuron.
#'
#' @param someneuronlist a neuron/neuronlist object
#' @param brain for context, plot some other objects (e.g. a brain from the
#' nat.templatebrains package such as FCWB, or any object that may be plotted
#' using plot3d)
#' @examples
#' \dontrun{
#' ## NB these neurons actually have their somata chopped off
#' correctedsomas = correct_root(Cell07PNs[1:3])
#' plot3d(correctedsomas, soma=TRUE)
#' }
#' @return a matrix of 3D points
#' @export
correct_root <- function(someneuronlist, brain = NULL){
correctedsomas = someneuronlist
nopen3d()
for (n in 1:length(someneuronlist)){
message("neuron ", n, " of ", length(someneuronlist))
w = someneuronlist[[n]]
print(names(someneuronlist[n]))
if(!is.null(brain)){rgl::plot3d(brain)}
rgl::plot3d(w, soma = T)
eps.df=w$d[endpoints(w),]
progress =F
while(!isTRUE(progress)){
cat ("Rotate brain and then hit [enter] to continue")
line <- readline()
message("Select new root from highlighted endpoints")
selection <- select3d()
selected.point <- selection(xyzmatrix(eps.df))
selected.pointno <- eps.df$PointNo[selected.point]
if (length(selected.pointno)!=1){
message("You must select exactly one end point. Try again!")
} else{
corrected = as.neuron(as.ngraph(w), origin = selected.pointno)
rgl::plot3d(corrected, soma = T, col = "blue")
progress = tolower(readline(prompt="Good enough? [y/n] "))=="y"
}
}
rgl::clear3d()
correctedsomas[[n]]=corrected
}
correctedsomas
}


#' Generate a 3D model from connector and/or tree node data
#'
#' @description Generate a mesh3d model based on points contained in a
#' neuronlist or neuron object, or another object that consists of 3D points.
#'
#' @param x a neuronlist or neuron object, or another object that consists of 3D points
#' @param substrate whether to make the model based on the 3D location of
#' connectors, neuron cable or both. Connectors are pre-synapse locations,
#' e.g. the pre-synapses of a \code{catmaidneuron} from the R package \code{catmaid})
#' @param alpha a single value or vector of values for alpha, fed to
#' \code{alphashape3d::ashape3d}. Selection is subsequently interactive.
#' @param auto.selection logical, whether or not to try and remove points based on interactively
#' choosing simple values for clustering.
#' @details Interactive function that allows a users to select points in 3D space from neuronlist/neuron objects,
#' or another object that is coercible in 3D points using \code{\link{xyzmatrix}}. Points can first be automatically chosen, by
#' selecting an integer number of nearest neighbours to find for each point using \code{nabor::\link{knn}}, and then a
#' maximum distance at which nodes can be part of a cluster. Next, \code{\link{select_points}} is used to manually pick desired 3D
#' points. Lastly, \code{alphashape3d::ashape3d} is used to create an alphashape around these points. The user can trial different values for
#' alpha until they get their desired result.
#' @examples
#' \dontrun{
#' # Make a model based off of fly olfactory projection neuron arbours
#' PN_blob = make_model(Cell07PNs)
#' }
#' @seealso \code{\link{prune_online}}
#' @return A mesh3d object
#' @export
make_model <- function(x,
substrate = c("cable", "connectors", "both"),
alpha = 30,
auto.selection = TRUE){
# need alphashape3d
if(!requireNamespace("alphashape3d", quietly = TRUE))
stop("Please install the suggested alphashape3d package in order to use make_model!")

# get starting points
substrate <- match.arg(substrate)
if (is.neuron(x)|is.neuronlist(x)) {
if (substrate=="connectors"){substrate.points <- xyzmatrix(do.call(rbind, lapply(x, function(x) x$connectors)))
} else if(substrate =="cable"){substrate.points <- xyzmatrix(x)
} else if (substrate == "both"){substrate.points <- rbind(xyzmatrix(x), do.call(rbind, lapply(x, function(x) xyzmatrix(x$connectors))))}
}else{
substrate.points <- xyzmatrix(x)
}

# automatic point selection
if (auto.selection == TRUE){
progress = "n"
while (progress == "n"){
groupsize <- as.numeric (readline(prompt="Each node must have this number of near neighbours (e.g. 10): "))
maxdistance <- as.numeric (readline(prompt="And all the neighbours must be this close in Euclidean space (e.g. 1): "))
if(groupsize<1){
groupsize <- 1
}
if(maxdistance<0){
maxdistance <- 0
}
neighbours <- nabor::knn(substrate.points, substrate.points, k = groupsize)
loose <- apply(neighbours$nn.dists, 1, function(x) {(any(as.numeric(x[1:ncol(neighbours$nn.dists)]) > maxdistance))})
keep <- c(neighbours$nn.idx[,1][!loose])
selected.points <- substrate.points[keep,]
message("Selected points in black, deselected in red")
ids <- rgl::points3d(selected.points, cl = 'black')
ids <- union(ids, rgl::points3d(substrate.points, col = 'red') )
progress <- readline(prompt="Continue to manual selection/deselection? y/n ")
pop3d(id=ids)
}
}

# manual point selection
progress = "p"
while (progress != "c"){
if (progress == 'r'){
removed.points <- picked.points(selected.points)
selected.points <- subset(selected.points, !removed.points)
}
if (progress == 'a'){
added.points <- subset(substrate.points, picked.points(substrate.points))
selected.points <- rbind(selected.points, added.points)
}
ids <- rgl::points3d(selected.points, cl = 'black')
ids <- union(ids, rgl::points3d(substrate.points, col = 'red') )
progress <- readline(prompt="Selected points in black. Remove (r) or add (a) points? Or continue to alphashape generation (c)? ")
if(progress != "c"){
picked.points <- rgl::select3d()
pop3d(id=ids)
}
}

# make alphashape
ids <- rgl::points3d(selected.points, cl = 'black')
ids <- union(ids, rgl::points3d(substrate.points, col = 'red') )
progress = "n"
while (progress == "n"){
message("alpha is ", alpha)
alphashape <- alphashape3d::ashape3d(unique(selected.points), alpha = alpha)
mesh3d <- as.mesh3d(alphashape)
ids=rgl::plot3d(mesh3d, alpha = 0.5, col = "purple", add = TRUE)
progress <- readline(prompt="Continue? y/n ")
if(progress == "n"){
alpha <- as.numeric (readline(prompt="Select a numeric value for alpha "))
}
}

# return
mesh3d
}

#' Interactively select 3D points in space
#'
#' @description Plot a set of 3D points in space and select a subset of them
#' interactively, using an rgl window
#'
#' @param points a matrix of 3D points to plot (or an object for which
#' \code{\link{xyzmatrix}} can extract 3D points).
#' @param clear_plot_on_exit Whether to remove points from the rgl scene when
#' selection has been completed.
#' @examples
#' \dontrun{
#' # Select points from 3 olfactory projection neurons
#' selected_points = select_points(Cell07PNs[1:3])
#' }
#' @seealso \code{\link{prune_online}}
#' @return A matrix describing selected 3D points
#' @export
select_points <- function (points, clear_plot_on_exit=FALSE) {
selected.points <- points <- xyzmatrix(points)
ids=rgl::points3d(selected.points)
progress = readline(prompt = "Selected points in black. Add (a) or remove (r) points, or continue (c)? ")
while (progress != "c") {
if (progress == "a") {
keeps = rgl::select3d()
keep.points <- keeps(points)
keep.points = subset(points, keep.points)
selected.points = rbind(selected.points, keep.points)
}
if (progress == "r") {
remove.points <- select3d()
removed.points <- remove.points(selected.points)
selected.points = subset(selected.points, !removed.points)
}
pop3d(id=ids)

ids=integer()
if (length(selected.points) > 0) {
ids=rgl::points3d(selected.points)
}
ids=union(ids, rgl::points3d(points, col = "red"))
progress = readline(prompt = "Add (a) or remove (r) points, or continue (c)? ")
}
if(clear_plot_on_exit)
pop3d(id=ids)
return(selected.points)
}

53 changes: 53 additions & 0 deletions R/ngraph.R
Original file line number Diff line number Diff line change
Expand Up @@ -595,3 +595,56 @@ EdgeListFromSegList<-function(SegList){
starts=unlist(lapply(sl,function(x) x[-length(x)]))
cbind(starts,ends)
}

#' Prune neuron(s) within a volume defined by a 3D mesh
#'
#' @details Prune a neuron to be within, or to exclude arbour within, a 3D
#' object that can be coerced into the mesh3d data structure
#'
#' @param x a \code{neuron} object
#' @param surf An \code{\link{hxsurf}} or \code{\link[rgl]{mesh3d}} object, or
#' any object coercible into \code{\link[rgl]{mesh3d}} by
#' \code{\link{as.mesh3d}}.
#' @param neuropil Character vector specifying a subset of the \code{surf}
#' object. This is only relevant when \code{surf} is of class
#' \code{\link{hxsurf}}. If NULL (default), then the full object given as
#' \code{surf} will be used for the pruning.
#' @param invert Logical when \code{TRUE} indicates that points inside the mesh
#' are kept.
#' @param ... Additional arguments for methods (eventually passed to
#' \code{\link{prune_vertices}}) surface should be pruned.
#' @inherit prune seealso
#' @examples
#' \dontrun{
#' ### Example requires the package nat.flybrains
#' LH_arbour = prune_in_volume(x = Cell07PNs, surf = nat.flybrains::IS2NP.surf,
#' neuropil = "LH_L", OmitFailures = TRUE)
#' }
#' @inheritParams prune
#' @return A pruned neuron/neuronlist object
#' @seealso \code{\link{as.neuron.ngraph}}, \code{\link{subset.neuron}},
#' \code{\link{prune.neuron}}, \code{\link{prune}}
#' @export
#' @rdname prune_in_volume
prune_in_volume <-function(x, surf, neuropil = NULL, invert = TRUE, ...)
UseMethod("prune_in_volume")

#' @export
#' @rdname prune_in_volume
prune_in_volume.neuron <- function(x, surf, neuropil = NULL, invert = TRUE, ...) {
if(is.null(neuropil)){
mesh = as.mesh3d(surf)
} else {
stopifnot(inherits(surf, "hxsurf"))
mesh = as.mesh3d(subset(surf, neuropil))
}
v = which(pointsinside(xyzmatrix(x), surf = mesh))
neuron = prune_vertices(x, verticestoprune=v, invert=invert, ...)
neuron
}

#' @export
#' @rdname prune_in_volume
prune_in_volume.neuronlist <- function(x, surf, neuropil = NULL, invert = TRUE, ...){
nlapply(x, prune_in_volume.neuron, surf = surf, neuropil = neuropil, invert = invert, ...)
}
Loading

0 comments on commit 7cefd30

Please sign in to comment.