Known issues with Traub et al 2005.

-----------------------------------

This is a quite complex and detailed model and as discussed in the [original paper](http://www.ncbi.nlm.nih.gov/pubmed/15525801?dopt=Abstract)

> Any model, even of a small bit of cortex, is subject to difficulties and hazards: limited data, large numbers of parameters, criticisms that models with complexity comparable to the modeled system cannot be scientifically useful, the expense and slowness of the necessary computations, and serious uncertainties as to how a complex model can be compared with experiment and shown to be predictive.

> The above difficulties and hazards are too real to be dismissed readily. In our opinion, the only way to proceed is through a state of denial that any of the difficulties need be fatal. The reader must then judge whether the results, preliminary as they must be, help our understanding.

Even the published Fortran version of this model was acknowledged to be incomplete. Each conversion of this model will deviate to a small or large extent from this version.

### Questions about physiological properties of model

**Dependence on Fast Regular Bursting cells for oscillatory behaviour**

**Prevalence of gap junctions**

**High current threshold for deep pyramidal firing**

**Not tested with external synaptic input**

### Limitations of the conversion of the model to NEURON

It is useful to read the [notes on conversion of this model to NEURON from Fortran](http://senselab.med.yale.edu/modeldb/ShowModel.asp?model=82894&file=\nrntraub\README) by Tom Morse and Michael Hines

**Slightly different method of running the simulation** (e.g. in Neuron information about spike is sent immediately, in Fortran every 0.1 ms )

**Sum of transmembrane currents in every single cells sums up to 0 only if you use cvode\_active\*

**Diffrent behaviour of NMDA synapse when thalamus is disconnected\* (some bug in Fortran version?)

In Fortran code:

     z = 0.d0  ! thalamus disconnected

     gAMPA_TCR_to_suppyrRS = z * gAMPA_TCR_to_suppyrRS

     gNMDA_TCR_to_suppyrRS = z * gNMDA_TCR_to_suppyrRS

     gAMPA_TCR_to_suppyrFRB = z * gAMPA_TCR_to_suppyrFRB

     gNMDA_TCR_to_suppyrFRB = z * gNMDA_TCR_to_suppyrFRB

    ...

gNMDA\_TCR\_to\_suppyrFRB becomes 0. Then when you compute NMDA activation

from TCR to suppyrFRB

    ....

    ! NMDA part

            if (delta.le.5.d0) then

           gNMDA_suppyrFRB(k,L) = gNMDA_suppyrFRB(k,L) +

         &  gNMDA_TCR_to_suppyrFRB * delta * 0.2d0

            else

           dexparg = (delta - 5.d0)/tauNMDA_TCR_to_suppyrFRB

              if (dexparg.le.5.d0) then

              z = dexptablesmall (int(dexparg*1000.d0))

             else if (dexparg.le.100.d0) then

              z = dexptablebig (int(dexparg*10.d0))

             else

              z = 0.d0

             endif

           gNMDA_suppyrFRB(k,L) = gNMDA_suppyrFRB(k,L) +

         &  gNMDA_TCR_to_suppyrFRB * z

            endif

    c Test for NMDA saturation

           z = NMDA_saturation_fact * gNMDA_TCR_to_suppyrFRB

           if (gNMDA_suppyrFRB(k,L).gt.z)

         &  gNMDA_suppyrFRB(k,L) = z

    ! end NMDA part

    ....

It seems that this piece of code, more precisely the last three lines:

<pre>

c Test for NMDA saturation

 z = NMDA\_saturation\_fact \* gNMDA\_TCR\_to\_suppyrFRB

 if (gNMDA\_suppyrFRB(k,L).gt.z)

 & gNMDA\_suppyrFRB(k,L) = z

</pre>

kills completely NMDA activation of suppyrFRB cells from all the other populations, not just TCR (except from nontuftRS cells, nontuftRS - suppyrFRB NMDA conductance is calculated after this block). In Neuron version there is no such behaviour.

An **updated version** of this model in NEURON is being worked on [here](https://github.com/hglabska/Thalamocortical/tree/Neuron_version_simplified_groucho_file/Neuron). The version allows to modify easily the network, e.g. to add new population (version commited on 26 June 2013 and later), replace one template by another e.g. tuftIB Traub cell by [Hay cell](http://senselab.med.yale.edu/ModelDb/ShowModel.asp?model=139653) ( version commited on 04 July 2013 or later). The main groucho.hoc file is simpler and much shorter (about 10 times), parameters like AMPA, GABA, NMDA conductances, connections between populations are defined in separated files.

#### Tests for [Neuron](http://senselab.med.yale.edu/ModelDb/ShowModel.asp?model=82894) and [Fortran](https://github.com/hglabska/Thalamocortical/tree/master/Fortran_ifc) version . Trying to reproduce results from the [article](http://www.ncbi.nlm.nih.gov/pubmed/15525801)

Remark: In Fortran version, compilation flag ~~finit-local-zero , seems to be important!

Thanks to kindness of Roger Traub, who sent us parameters which were used to generate figures 2. and 7. in the [article](http://www.ncbi.nlm.nih.gov/pubmed/15525801) , [we](http://www.opensourcebrain.org/groups/71) were able to test how well we can reproduce the results on different version of the model.

h5. Single Cell

Results from Appendix A~~ activity of single cells after applying some current to the soma, were reproduce reasonable well in Neuron version. For more data look [here](http://figshare.com/articles/Neuron_single_cell/861118) .

To compare the single cell result in Neuron with Fortran version you can use [this](https://github.com/hglabska/Thalamocortical/tree/master/Fortran_ifc) code with makefile.single\_cell instead of makefile. This version contains additional 14 programs to simulate single cell from every of 14 populations.

The biggest challenge in Appendix A is to reproduce fig A4C: applying some pulse current in apical dendrite caused somatic burst.

![](A4C.png)

First difficulties is to estimate the amplitude of the current (It is not describe in article).

*I =3\* /10)) \* /20)) nA,*

 looks reasonable well:

![](pulse.png)

but Neuron result doesn’t look similar like the result in the article (colours: green D1, black D2, red soma):

**Neuron**

![](tuftIB_Neuron_voltage.png)

also [Fortran version](https://github.com/hglabska/Thalamocortical/tree/master/Fortran_ifc) (tuftIB.f ) of the model failed to reproduce the somatic burst with the same stimulus.

**Fortran**

![](voltage_tutfIB_fortran.png)

It is possible to obtain this somatic spikes ( in both Neuron and Fortran version) after depolarizing the soma by 1nA current and increasing the apical stimulus 3 times. Decreasing the depolarizing somatic currents two times (0.5 nA) , or using the apical stimulus like at the beginning (*I =3\* /10)) \* /20))* ), caused that the somatic spikes disappear.

**Neuron**

![](neuron_burst.png)

**Fortran**

![](fortran_burst.png)

**Remark :** Look at the difference in the somatic membrane potential after the burst, between Fortran and Neuron versions.

For more data look here ( EPSP means apical stimulus with amplitude *I =3\* /10)) \* /20))*):

Fortran

[EPSP](http://figshare.com/articles/tuftIB_cell_strong_dendritic_input_additional_somatick_input_Fortran_version/861127)

[somatic current 1nA + 3\*EPSP](http://figshare.com/articles/tuftIB_cell_strong_dendritic_input_additional_somatick_input_Fortran_version/861127)

[somatic current 1nA + EPSP](http://figshare.com/articles/_tuftIB_cell_dendritic_input_somatic_input_0_5_nA_Fortran_version/861133)

[somatic current 0.5 nA + 3\*EPSP](http://figshare.com/articles/tuftIB_cell_strong_dendritic_input_small_somatic_input_Fortran_version/861130)

Neuron

[EPSP, 3 \* EPSP ](http://figshare.com/articles/tuftIB_cell_dendritic_input_Neuron_version/861138)

[somatic current 1 nA + EPSP, 3\*EPSP](http://figshare.com/articles/_tuftIB_cell_dendritic_input_somatic_current_1_nA_Neuron_version/861145)

[somatic current 0.5 nA + EPSP, 3\*EPSP](http://figshare.com/articles/_tuftIB_cell_dendritic_input_somatic_current_0_5_nA_Neuron_version_/861153)

##### Figure 2

“Simulation of kainate-induced gamma oscillations”

![](test2_labels.png)

The results in both Neuron and Fortran version looks quite similar. Only be aware that activity of suppyrRS differs much between individual cells. One questionable issue is appearance of the burst after about 1500 ms in Fortran and nearly 1200 ms in Neuron version (not shown here), which they didn’t report in the [article](http://www.ncbi.nlm.nih.gov/pubmed/15525801).

You can download the data (+ rasterplot) for Fig 2 from Fortran and Neuron simulation: [Fortran data](http://figshare.com/articles/2_Fortran/858844) and [Neuron data](http://figshare.com/articles/2_Neuron_use_traubexac_0/858878).

For Neuron simulation you can also compare the result with simulation using the “traub\_exact()” algorithm: [Neuron traub\_excat() data](http://figshare.com/articles/2_Neuron_use_traubexac_1/858893). More about “traub\_excat()” algorithm you can read in [notes on conversion of this model to NEURON from Fortran](http://senselab.med.yale.edu/modeldb/ShowModel.asp?model=82894&file=\nrntraub\README) by Tom Morse and Michael Hines.

##### Figure 7

“Effects of disinhibition in model (cortex only, with thalamus disconnected), when there are open gap junctions between the axons of the respective principal cell populations (superficial pyramids, spiny stellates, layer 5 pyramids, layer 6 pyramids), and spiny stellates are strongly interconnected by AMPA receptors .”

Figure 7 from the [article](http://www.ncbi.nlm.nih.gov/pubmed/15525801)

![](7paper.png)

**7A**

In the [article](http://www.ncbi.nlm.nih.gov/pubmed/15525801) they raported about consisting of 17 burst complexes that terminate spontaneously. The last 5 of the bursts are shown. Results from the [Fortran](https://github.com/hglabska/Thalamocortical/tree/master/Fortran_ifc) version are very similar but only 14 bursts appears. In Neuron version the result is much different.

![](7A_small_labels.png)

You can download the data (+ rasterplot) for Fig 7A from Fortran and Neuron simulation: [Fortran data](http://figshare.com/articles/7A/855456) and [Neuron data](http://figshare.com/articles/7A_Neuron_use_traubexac_0/855486).

For Neuron simulation you can also compare the result with simulation using the “traub\_exact()” algoritm: [Neuron traub\_excat() data](http://figshare.com/articles/7A_Neuron_use_traubexac_1/856743). More about “traub\_excat()” algoritm you can read in [notes on conversion of this model to NEURON from Fortran](http://senselab.med.yale.edu/modeldb/ShowModel.asp?model=82894&file=\nrntraub\README) by Tom Morse and Michael Hines.

**7B**

Results from the [Fortran](https://github.com/hglabska/Thalamocortical/tree/master/Fortran_ifc) version looks again very similar, although gives much more complex bursts, at least 6, when prolong the simulation up to 2000 ms (results not shown here - [download](http://figshare.com/articles/7B_Fortran_long/858794) ) .

![](7B_small_labels.png)

You can download the data (+ rasterplot) for Fig 7B from Fortran and Neuron simulation: [Fortran data](http://figshare.com/articles/7B_Fortran/855478) and [Neuron data](http://figshare.com/articles/7B_Neuron_use_traubexac_0/856699) compare with [Neuron traub\_excat() data](http://figshare.com/articles/7B_Neuron_use_traubexac_1/856753) .

**7C**

![](7C_small_labels.png)

You can download the data (+ rasterplot) for Fig 7C from Fortran and Neuron simulation. [Fortran data](http://figshare.com/articles/7C/855462) and [Neuron data](http://figshare.com/articles/7C_Neuron_use_traubexac_0/856722) compare with [Neuron traub\_excat() data](http://figshare.com/articles/7C_Neuron_use_traubexac_1/858769) .

**7D**

![](7D_small_labels.png)

You can download the data (+ rasterplot) for Fig 7D from Fortran and Neuron simulation. [Fortran data](http://figshare.com/articles/7D_Fortran/855470) and [Neuron data](http://figshare.com/articles/7D_Neuron_use_traubexac_0/856732) compare with [Neuron traub\_excat() data](http://figshare.com/articles/7D_Neuron_use_traubexact_1/858779)

#### Response to simple stimulus - comparison between Fortran and Neuron versions.

**Gap junctions are closed,** thalamus is connected with cortex.

Stimulus: current injection to thalamic (TCR) somas . Current delay 300 ms, duration 2 ms, amplitude 3 nA.

**Fortran**

In *normal* case there is small, short response in layers 2/3, 4 and inhibitory neurons in layers 5/6. The answer is much better visible if we decrease GABA conductances, but still there is no response in layer 5 and 6 in pyramidal cells (except ectopic spikes).

**Neuron**

No response in the cortex in *normal* case. Answer in layers 2/3, 4 and inhibitory neurons in layers 5/6 after decreasing GABA conductances, but activity in layers 2/3 is shorter than in Fortran case, single spike in pyramidal cells layer 6 and no response in layer 5 (only ectopic spikes).

![](thalamus_awake0.png)

Applying additional current tu somas in pyramidal cells in layer 5 (1 nA) and 6 (0.75 nA) ( awake = 1 in Neuron version).

The additional stimulus is to big, a lot of spontaneous burst in every case. In Fortran version response in layer 2/3 lasts longer.

![](thalamus_awake1.png)

Additional current to somas; 0.5 nA in pyramidal cells in layer 5 and 0.375 nA in somas of pyramids in layer 6.

**Fortran**

The additional stimulus is still to big in Fortran version ( a lot of spontaneous burst).

**Neuron**

Spontaneous burst still exist but there are very seldom (not shown on the picture). Now can observe response in every layer.

![](thalamus_awake05.png)

Additional current to somas; 0.2 nA in pyramidal cells in layer 5 and 0.15 nA in somas of pyramids in layer 6.

**Fortran**

All layers answer to stimulus, only response in pyramids layer 5 and 6 is quite late.

**Neuron**

Again no response in layers 5 and single spike or no response in layer 6 pyramids.

![](thalamus_awake02.png)

**Conclusions/Remarks:**

* Fortran and Neuron code doesn’t generate the soma output even when gap junctions are closed

* in Fortran version response in layer 2/3 is more complex, 3 bursts versus 1 (why? )

* when gap junctions are closed in *normal* condition when GABA conductance is not decreased, response in layer 2/3 pyramids last extremely short (single spikes)

### Limitations of the conversion of the model to MOOSE

TODO…

### Limitations of the conversion of the model to NeuroML

**Optimal spatial discretisation for each cell needs to be investigated**

Important details of the process of conversion of the cell models to NeuroML, and matching cell behaviour across simulators is present in the [2010 NeuroML paper](http://www.ploscompbiol.org/article/info%3Adoi%2F10.1371%2Fjournal.pcbi.1000815).

The spatial discretisation of the cells influenced precise spike timing. Changing the number of compartments/points used to calculate the membrane potential changed the timing of the cell (e.g. changing the value of nseg in NEURON on all sections). See below for an example of how 3 cells with differing numbers of compartments converged at different rates. A) Nucleus reticularis thalami (nRT) cell; B) Superficial Low Threshold spiking (LTS) cell; C) Layer 6 Non-tufted Regular Spiking pyramidal cell. Traces for NEURON (black) and MOOSE (green) and GENESIS (red).

![](http://www.opensourcebrain.org/attachments/download/114/converge.png)

**NMDA conductance wave form**

The NMDA synapse model used in the network has an unconventional form, with a scaling factor rising lineally between 0 and 5ms, and decaying exponentially. This can probably be approximated by a double exponential synapse (coupled with v & [Mg] dependent blocking mechanism).

**Firing rate vs. injected current of cells**

Many of the cells show unusual F/I curves.

![](/attachments/download/113/ifcurve.png)

**Support in NeuroML**

All model elements from the neuroConstruct generated network can be exported to valid NeuroML v1.8.1.

Model can be exported to [(mostly valid) NeuroML 2](https://github.com/OpenSourceBrain/Thalamocortical/tree/master/neuroConstruct/generatedNeuroML2), but there is not yet an application that can handle such detailed NML2 models (but we’re [working on it](https://github.com/NeuroML/org.neuroml.export/blob/development/src/main/java/org/neuroml/export/neuron/NeuronWriter.java)).