Proposal for multi-compartmental models

I don't have the time at the moment to go into this in as much detail (or with as much thinking through) as I would like but I wanted to get it out there so it has a chance to make the agenda for the upcoming meeting because it will probably require some workshopping.

Building on my last proposal in issue #4 for sub component containers, what do you think of this format?

<ComponentClass name="MultiCompartmentHodgkinHuxleyContainer" type="dynamics_container">
  <Domain name="soma">
    <SubComponent target="Na">
      <PortConnection receiver="V">
        <Sender component="memb" name="V"/>
      </PortConnection>
    </SubComponent>
    <SubComponent target="K">
      <PortConnection receiver="V">
        <Sender component="memb" name="V"/>
      </PortConnection>
    </SubComponent>
    <SubComponent target="leak">
      <PortConnection receiver="V">
        <Sender component="memb" name="V"/>
      </PortConnection>
    </SubComponent>
    <SubComponent target="compartmentBase">
      <PortConnection receiver="membrane_i">
        <Sender component="na" name="i"/>
        <Sender component="k" name="i"/>
        <Sender component="leak" name="i">
      </PortConnection>
      <PortConnection receiver="neighbouring_v">
        <Sender component="compartmentBase" name="v" domain="soma">
        <Sender component="compartmentBase" name="v" domain="dendrites">
      </PortConnection>
      <PortConnection receive="neighbouring_axonal_g">
        <Sender component="compartmentBase" name="axonal_g" domain="soma">
        <Sender component="compartmentBase" name="axonal_g" domain="dendrites">
      </PortConnection>
    </SubComponent>
  </Domain>
  <Domain name="dendrites">
    <SubComponent target="leak">
      <PortConnection receiver="V">
        <Sender component="memb" name="V"/>
      </PortConnection>
    </SubComponent>
    <SubComponent target="compartmentBase">
      <PortConnection receiver="membrane_i">
        <Sender component="leak" name="i"/>
      </PortConnection>
      <PortConnection receive="neighbouring_v">
        <Sender component="compartmentBase" name="v" domain="soma">
        <Sender component="compartmentBase" name="v" domain="dendrites">
      </PortConnection>
      <PortConnection receive="neighbouring_axonal_g">
        <Sender component="compartmentBase" name="axonal_g" domain="soma">
        <Sender component="compartmentBase" name="axonal_g" domain="dendrites">
      </PortConnection>
    </SubComponent>
  </Domain>
  <Compartments xmlns="http://www.neuroml.org/neuromlv2">
    ... compartment specification goes here, potentially using NeuroML2 to specify a 
        predetermined tree or a component specifying an algorithm to generate the 
        compartment tree ... 
  </Compartments>
  <Compartment2DomainMapping xmlns="http://www.neuroml.org/neuromlv2">
    ... again this could be a list of mappings as in NeuroML2/MorphML (I think Padraig
        mentioned a way to ensure the segment groups are exclusive) or an algorithm to
        map the compartments to exclusive domains ... 
  </Compartment2DomainMapping>
</ComponentClass>

In this format, if the "domain" attribute format is not provided to the "Sender" element then the component is in the same compartment, otherwise if the domain attribute is provided, any connections from that domain are connected to the port. This would be a case that would require the array port connections I suggested in issue #13. I am not completely happy with this, particularly having to pass the axonal conductance as a separate port but maybe it is a start that we can work with.

It might be useful to also introduce the notion of compartments without a spatial extent to support reduced compartmental models, such as e.g. Pinsky-Rinzel CA3 model. I think these would fit well into the domain/subcomponent object model but without requiring the full machinery for cable equation models.

Hi Ivan, yes I have been thinking about that too, whether the morphology tree should only define the compartment connectivity and the compartment dimensions become a property of the domain (which could still be specified per compartment or even from a distribution).

I forgot to mention that this multi-compartmental format is meant to to be a extension namespace that could be flattened first to a sub-component container and then to a "core NineML" representation for simulators that don't internally support multi-compartmental modelling (i.e. not NEURON or GENESIS/MOOSE)

I am coming to the conclusion that in order for this multi-compartmental proposal to work we will need to introduce the concept of array ports (see #13, although maybe "vector ports" would be a better name) to handle connections between multiple compartments (i.e. forks). This in turn would probably need the introduction of some sort of vector calculus to handle them, which although it could be described in MathML may be a bridge too far for us at this stage...

Some more thoughts I have had on the proposal I put forward are:

The port connections that are received from other compartments will probably need to use a specific "RemoteSender" tag to distinguish them from local connections. Also the end of the compartment they connect to may need to be specified in some cases i.e.

    <SubComponent target="cylindricalCompartment">
      <PortConnection receiver="membrane_i">
        <Sender component="leak" name="i"/>
      </PortConnection>
      <PortConnection receiver="proximalV">
        <RemoteSender component="cylindricalCompartment" name="v" connectedTo="proximal"/>
      </PortConnection>
      <PortConnection receiver="distalV">
        <RemoteSender component="cylindricalCompartment" name="v" connectedTo="distal"/>
      </PortConnection>
      <PortConnection receiver="proximalG">
        <RemoteSender component="cylindricalCompartment" name="axonalG" connectedTo="proximal"/>
      </PortConnection>
      <PortConnection receiver="distalG">
        <RemoteSender component="cylindricalCompartment" name="axonalG" connectedTo="distal"/>
      </PortConnection>
    </SubComponent>

If the port is only connected to compartments of a specific domain(s) (particularly important if some domains do not define the sender port in question) then this could be specified like

<PortConnection receiver="axialDiffusion">
  <RemoteSender component="Ca" name="concentration">
    <FromDomain>soma</FromDomain>
    <FromDomain>proximalDendrites</FromDomain>
  </RemoteSender>
</PortConnection>

Also, in this scheme port connections that receive remote senders must be either a "reduce" port or an "array" port so if there are no connections that meet the specified criteria the received values are still defined.

I just had a brief look at the "spatial" SBML package and found that it uses almost exactly the same names as I was suggesting here (i.e. Domain, Compartment, CompartmentMapping, etc). Don't know whether this is a good thing or a bad thing though (or neither), as they will refer to slightly different beasts.

I will have a closer look at the spatial package and see if there are any other ideas that could be good to include.

SBML has a proposed package called "arrays" (http://sbml.org/Documents/Specifications/SBML_Level_3/Packages/Arrays_and_Sets_%28arrays%29), which can handle vectors, arrays and simple operations (although no inverse), which might be good to look at.

Okay, so I have been doing some more thinking about this issue and have come up with an alternative structure I am quite happy with, inspired by the syntax we came up with for Projections and the SWC format. What do you think of this?

<?xml version="1.0" encoding="UTF-8"?>
<NineML xmlns="http://nineml.net/9ML/1.0">
  <Dimension name="specificCapacitance" m="-1" l="-4" t="4" i="2"/>
  <Dimension name="length" l="1"/>
  <Dimension name="resistance" m="1" l="2" t="-3" i="-2"/>
  <Dimension name="concentration" l="-3" n="1"/>
  <ComponentClass name="ExampleMultiCompartmentalModel">
    <Parameter name="C" dimension="specificCapacitance"/>
    <AnalogSendPort name="somaVoltage" dimension="voltage"/>
    <AnalogReceivePort name="T" dimension="temperature"/>
    <EventSendPort name="outgoingSpike"/>
    <MultiCompartmentDynamics>
      <BranchStructure>
        <!-- Similar to "parent node" column in SWC --> 
        <ExternalArrayValue ...>
        <Annotations>
          <3DPoints>
            <!-- Same as the X, Y, Z coordinates in SWC -->    
            <XCoords><ExternalArrayValue ...></XCoords>
            <YCoords><ExternalArrayValue ...></YCoords>
            <ZCoords><ExternalArrayValue ...></ZCoords>
          </3DPoints>
        </Annotations>
      </BranchStructure>
      <Mapping>
        <!-- Same as "type" column in SWC with ids that correspond to domain map_ids -->
        <ExternalArrayValue ...>
      </Mapping>
      <Domain name="soma" map_id="0">
        <ComponentClass name="somaDynamics">
          <!-- Dynamics component, typically mult-component container -->
          <Parameter name="length" dimension="length"/>  <!-- Could be derived from 3D points -->
          <Parameter name="diameter" dimension="length"/> <!-- Could be read from SWC file -->
          <Parameter name="axialResistance" dimension="resistance"/>
          <Parameter name="C" dimension="specificCapacitance">C</Parameter>
          <AnalogReceivePort name="proximalVoltage" dimension="voltage"/>
          <AnalogArrayReceivePort name="distalVoltage" dimension="voltage"/>
          <EventSendPort name="toDendritePort"/>
          ...
        </ComponentClass>
        <FromProximal>
          <PortConnection port="proximalVoltage">
            <Sender port="voltage"/>
          </PortConnection>
        </FromProximal>
        <FromDistal>
          <PortConnection port="distalVoltage">
            <Sender port="voltage"/>
          </PortConnection>
        </FromDistal>        
      </Domain>
      <Domain name="dendrites" map_id="2">
        <ComponentClass name="dendriteDynamics">
          <!-- Dynamics component, typically multi-component container -->
          <Parameter name="length" dimension="length"/>
          <Parameter name="diameter" dimension="length"/>
          <Parameter name="axialResistance" dimension="resistance"/>
          <Parameter name="C" dimension="specificCapacitance">C</Parameter>
          <AnalogReceivePort name="proximalVoltage" dimension="voltage"/>
          <AnalogArrayReceivePort name="distalVoltage" dimension="voltage"/>
          <EventReceivePort name="fromSomaPort"/>
          <AnalogSendPort name="withinDendritesSendPort" dimension="concentration"/>
          <AnalogArrayReceivePort name="withinDendritesReceivePort" dimension="concentration"/>
          ...
        </ComponentClass>
        <FromProximal>
          <PortConnection port="proximalVoltage">
            <Sender port="voltage"/>
          </PortConnection>
          <PortConnection port="fromSomaPort">
            <Sender name="soma" port="toDendritePort"/>
          </PortConnection>
        </FromProximal>
        <FromDistal>
          <PortConnection port="distalVoltage">
            <Sender port="voltage"/>
          </PortConnection>
          <PortConnection port="withinDendritesReceivePort">
            <Sender name="dendrites" port="withinDendritesSendPort"/>
          </PortConnection>
        </FromDistal>        
      </Domain>
    </MultiCompartmentDynamics>
  </ComponentClass>
</NineML>

The BranchStructure is an array containing the indices of the parent node of the segment. In the future instead of an ArrayValue (or ExternalArrayValue) this could be substituted for a component that generates the branch structure. The format is similar to the parent column in the SWC format with the caveats:
- the root node's ID is 0 and is omitted
- the IDs of the remaining nodes (which correspond to segments) start from 1 and contiguously increase
- (NB: this allows IDs to be implied from the array index so they don't need to be defined separately).
The mapping array is a list of indices corresponding to the map_id's of the domains. Again this is very similar to the SWC "type" column. In fact it would be nearly be possible to use an SWC file as an ExternalArrayValue file, with the exception that it would need column headers and the root node would need to be omitted. This could also be substituted for a component that algorithmically maps the segments to different domains.
Domain component classes would contain all the ion channels and concentration dynamics for each "class" of compartment (or potentially another more abstract type of compartment such as an Izhikevich one).
- Like Projections the parameters of the components could either be SingleValues, ArrayValues or Components. The length of the arrays could either be the same as the number of compartments in the model, in which case only the values that map to the domain would be used, or the same as the number of compartments in the domain.
- PortConnection>Sender tags without a name attribute accept ports from any compartments that are connected to it (assumes ports have the same name in different domains). With a name attribute they are only connected if the compartment is of the corresponding domain (probably would need to be a reduce port).
- FromDistal>PortConnections should be connected to AnalogReducePorts, EventReceivePorts or AnalogArrayReceivePorts (this is where they come in)
3D points could be optionally included as annotations to allow plotting and other analysis but would not effect the length parameter in the dynamics component or the interpretation of the dynamics in general (although could possibly be used for optimisations such as in IBM's Neural Tissue Simulator).

The plan is for this to be an "extended-format", which could be first flattened down to a MultiComponentDynamics model and then further flattened to a vanilla Dynamics component for implementation in NEST, GeNN and other simulators with "unstructured" neuronal models, or translated directly from the MultiCompartmentDynamics format into simulators with "structured" neuronal models such as NEURON, GENESIS, etc.

Feedback welcome. If you guys like it I will write it up into a proposal along with the multi-component dynamics proposed syntax (see #4).

Proposal in dev2.0 spec:

<?xml version="1.0" encoding="UTF-8"?>
<NineML xmlns="http://nineml.net/9ML/2.0">
  <MultiCompartmental name="ExampleMultiCompartmentModel">
    <Definition>ExampleMultiCompartmentClass</Definition>
    <Branches>
      <ExternalArrayValue url="./example_compartments.txt"
        mimetype="vnd.net.nineml.arrayvalue.text" columnName="parentID"/>
    </Branches>
    <Mapping>
      <ExternalArrayValue url="./example_compartments.txt"
        mimetype="vnd.net.nineml.arrayvalue.text" columnName="domain"/>
    </Mapping>
    <Domain name="soma">
      <MultiComponent name="somaDynamics">
        <SubComponent name="PyramidalNa">
          <Component>
            <Prototype url="http://nineml.net/9ML/2.0/catalog/channels"
              >HHPyramidalNa</Prototype>
            <Property name="gbar" units="uS">
              <SingleValue>0.0036</SingleValue>
            </Property>
          </Component>
          <ReceiveConnection port="V">
            <FromSister component="membrane" port="V"/>
          </ReceiveConnection>
        </SubComponent>
        <SubComponent name="PyramidalK">
          <Component>
            <Prototype url="http://nineml.net/9ML/2.0/catalog/channels"
              >HHPyramidalK</Prototype>
            <Property name="gbar" units="uS">
              <SingleValue>0.0025</SingleValue>
            </Property>
          </Component>
          <ReceiveConnection port="V">
            <FromSister component="membrane" port="V"/>
          </ReceiveConnection>
        </SubComponent>
        <SubComponent name="Kv4">
          <Component>
            <Prototype url="http://nineml.net/9ML/2.0/catalog/channels"
              >HHPyramidalK</Prototype>
            <Property name="gbar" units="uS">
              <SingleValue>0.0025</SingleValue>
            </Property>
          </Component>
          <ReceiveConnection port="V">
            <FromSister component="membrane" port="V"/>
          </ReceiveConnection>
        </SubComponent>
        <SubComponent name="leak">
          <Component>
            <Definition url="http://nineml.net/9ML/2.0/catalog/channels"
              >Leak</Definition>
            <Property name="g" units="uS">
              <SingleValue>0.001</SingleValue>
            </Property>
            <Property name="e_rev" units="mV">
              <SingleValue>-21</SingleValue>
            </Property>
          </Component>
          <ReceiveConnection port="V">
            <FromSister component="membrane" port="V"/>
          </ReceiveConnection>
        </SubComponent>
        <SubComponent name="membrane">
          <Component>
            <Definition url="http://nineml.net/9ML/2.0/catalog/channels/basicCompartment.9ml"
              >BasicCompartment</Definition>
            <Property name="length" units="um">
              <SingleValue>15.0</SingleValue>
            </Property>
            <Property name="axialResistance" units="ohm_per_cm">
              <SingleValue>100.0</SingleValue>
            </Property>
            <Property name="diameter" units="um">
              <SingleValue>15.0</SingleValue>
            </Property>
            <Property name="C" units="uF">
              <SingleValue>1.0</SingleValue>
            </Property>
          </Component>
          <ReceiveConnection port="i">
            <FromSister component="Na" name="i"/>
            <FromSister component="K" name="i"/>
            <FromSister component="leak" name="i"/>
          </ReceiveConnection>
        </SubComponent>
        <PortExposure name="proximalVoltage" component="membrane" port="proximalV"/>
        <PortExposure name="distalVoltage" component="membrane" port="distalV"/>
      </MultiComponent>
      <ReceiveConnection port="proximalVoltage">
        <FromProximal port="voltage"/>
      </ReceiveConnection>
      <ReceiveConnection port="distalVoltage">
        <FromDistal port="voltage"/>
      </ReceiveConnection>
    </Domain>
    <Domain name="dendrites">
      <Reference>DendriteDynamics</Reference>
      <ReceiveConnection port="proximalVoltage">
        <FromProximal port="voltage"/>
      </ReceiveConnection>
      <ReceiveConnection port="distalVoltage">
        <FromDistal port="voltage"/>
      </ReceiveConnection>
      <ReceiveConnection port="withinDendritesReducePort">
        <FromProximal domain="dendrites" port="withinDendritesSendPort"/>
        <FromDistal domain="dendrites" port="withinDendritesSendPort"/>
      </ReceiveConnection>
    </Domain>
    <Property name="xCoords">
      <ExternalArrayValue url="./example_compartments.txt"
        mimetype="vnd.net.nineml.arrayvalue.text" columnName="X"/>
    </Property>
    <Property name="yCoords">
      <ExternalArrayValue url="./example_compartments.txt"
         mimetype="vnd.net.nineml.arrayvalue.text" columnName="Y"/>
    </Property>
    <Property name="zCoords">
       <ExternalArrayValue url="./example_compartments.txt"
         mimetype="vnd.net.nineml.arrayvalue.text" columnName="Z"/>
    </Property>
  </MultiCompartmental>
  <Dimension name="time" t="1"/>
  <Dimension name="voltage" m="1" l="2" t="-3" i="-1"/>
  <Dimension name="conductance" m="-1" t="3" l="-2" i="2"/>
  <Dimension name="area" l="2"/>
  <Dimension name="length" l="1"/>
  <Unit symbol="mV" dimension="voltage" power="-3"/>
  <Unit symbol="um" dimension="length" power="-3"/>
  <Unit symbol="ms" dimension="time" power="-3"/>
  <Unit symbol="cm_square" dimension="area" power="-4"/>
  <Unit symbol="mS" dimension="conductance" power="-3"/>
</NineML>

The committee has agreed that both multicomponent and multicompartmental models can be supported in the Core using the Core specification outlined in Issue #46. The image below shows the proposal for multicompartmental models. 20150402_111632

The layout that was agreed upon for the core structure layer with ComponentArrays and ConnectionGroups (see #46) is planned to be used for multi-component and multi-compartmental models. This has the nice benefit that it would make it easier to optimally load balance simulations by spreading large cell models over different processing cores, while the same optimisation algorithm could be used to group small cells connected by gap junctions (since gap junctions and components are both communicated via analog ports without delay and would therefore look the same to the simulator).

Under the hood in the core a multi-compartmental model would therefore look like (see attached image above):

ComponentArrays for each component type (it may be practical to split these along domain groups if some components are present in different domains)
ConnectionGroups forming an explicit connection matrices (flattened into a single vector) to connect each of the trans-compartment components, i.e. those that communicate between compartments (typically only membrane dynamics for standard multi-compartmental models)
ConnectionGroups connecting components within compartments (typically one-to-one), the same as is proposed for multi-component models.

The question now is how best to represent this in a more user-friendly way. One option is to have a connectivity matrix, which could possibly be created by a selection rule as well as loaded from file (will have to see how difficult this is though), to represent to the compartments and their neighbouring relationship to each other (e.g. for an unbranching cylinder it would be a row of ones offset from the diagonal). This is how Herman Cuntz's "trees toolbox" represents morphologies I believe, so that would be kind of nice. Trans-compartment components could then slice up this matrix to get the connectivity rule for their inter-compartment ConnectionGroup during the translation to the core.

When it comes to specifying the distribution of components over the compartments perhaps the selection syntax could be reused. However, this could be difficult and will require generic properties specifying things like distance from the soma, which might be difficult to derive in this format. Also, I am not sure whether the distribution of components should be grouped into domains (groups compartments with consistent component sets) or should be free to map onto any segment (as they are in NeuroML). Free seems to fit better with the selection rule approach but then it is not possible to specify/check that all the receive ports are connected when the components are created in the Structure Layer.

INCF / nineml-spec

Proposal for multi-compartmental models #30