Export a subset of loaded morphology

When Python is the interpreter, what is a good
design for the interface to the basic NEURON
concepts.

Moderator: hines

Post Reply
ahwillia
Posts: 17
Joined: Wed Apr 23, 2014 2:17 pm

Export a subset of loaded morphology

Post by ahwillia »

I have a use case where I import a (large, complex) morphology from an swc file.

Then I have a procedure for pruning branches to make a simpler morphology. At the end of the day, I have a python list with of NEURON section objects (they are all connected and the list contains the root node).

Now I want to export a new morphology file with just this subset of sections. The exported file could also be swc, though I'm not attached to this (anything that is easily reimported back to NEURON would work fine for me). Is there an easy function/procedure to do this, or do I need to write my own swc parser? Thanks!
ramcdougal
Posts: 267
Joined: Fri Nov 28, 2008 3:38 pm
Location: Yale School of Public Health

Re: Export a subset of loaded morphology

Post by ramcdougal »

To the best of my knowledge, there is no built-in 1-line solution, but:

The key functions and classes are h.parent_connection, h.section_orientation, and h.SectionRef. It's a little tricky to do this right because of the distinction between a "true parent" and a mere "parent".

There are two functions in neuron.rxd.morphology that allow you to ignore this issue: parent(sec) and parent_loc(sec, trueparent).

SWC has the advantage of being portable, but it loses some information (in particular, the section names)... and some parsers may expect the parents to already be defined before their children.

Here's some code that generates JSON containing the morphology and connectivity, along with a demo:

Code: Select all

from neuron import h
from neuron.rxd.morphology import parent, parent_loc
import json

def morphology_to_dict(sections):
    section_map = {sec: i for i, sec in enumerate(sections)}
    result = []
    h.define_shape()
    for sec in sections:
        my_parent = parent(sec)
        my_parent_loc = -1 if my_parent is None else parent_loc(sec, my_parent)
        my_parent = -1 if my_parent is None else section_map[my_parent]
        n3d = int(h.n3d(sec=sec))
        result.append({
            'section_orientation': h.section_orientation(sec=sec),
            'parent': my_parent,
            'parent_loc': my_parent_loc,
            'x': [h.x3d(i, sec=sec) for i in xrange(n3d)],
            'y': [h.y3d(i, sec=sec) for i in xrange(n3d)],
            'z': [h.z3d(i, sec=sec) for i in xrange(n3d)],
            'diam': [h.diam3d(i, sec=sec) for i in xrange(n3d)],
            'name': sec.hname()           
        })
    return result

if __name__ == '__main__':

    s = [h.Section(name='s[%d]' % i) for i in xrange(13)]

    """
        Create the tree
        
              s0
        s1    s2         s3
        s4           s5      s6
        s7         s8 s9       s10 
    """
    for p, c in [[0, 1], [0, 2], [0, 3], [1, 4], [4, 7], [3, 5], [3, 6], [5, 8], [5, 9], [6, 10]]:
        s[c].connect(s[p])
    
    import json
    print json.dumps(morphology_to_dict([s[3], s[5], s[8], s[0], s[1], s[4], s[7]]), indent=2)
ahwillia
Posts: 17
Joined: Wed Apr 23, 2014 2:17 pm

Re: Export a subset of loaded morphology

Post by ahwillia »

Ok nice - and then you would need a `load_json`function, something like:

Code: Select all

def load_json(morphfile):

	with open(morphfile, 'r') as f:
		secdata = json.load(morphfile)

	seclist = []
	for sd in secdata:
		# make section
		sec = h.Section(name=sd['name'])
		seclist.append(sec)

		# make 3d morphology
		for x,y,z,d in zip(sd['x'], sd['y'], sd['z'], sd('diam')): 
			h.pt3dadd(x, y, z, d, sec=sec)

	# connect children to parent compartments
	for sec,sd in zip(seclist,secdata):
		if sd['parent_loc'] >= 0:
			parent_sec = sec[sd['parent_loc']]
			sec.connect(parent_sec)

	return seclist
ramcdougal
Posts: 267
Joined: Fri Nov 28, 2008 3:38 pm
Location: Yale School of Public Health

Re: Export a subset of loaded morphology

Post by ramcdougal »

Almost.

In my example, parent is the index of the parent section; parent_loc is the normalized position within that section. Thus:
  • parent_sec = sec_list[sd['parent']] # not parent_loc, grab from sec_list not sec
  • sec.connect(parent_sec(sd['parent_loc']), sd['section_orientation'])
This last line is important because otherwise you would be implicitly assuming section_orientation is 0 and parent_loc is 1. (In particular, parent_loc is often not 1 at the soma.)
Post Reply