pre-Bötzinger complex variability (Fietkiewicz et al. 2016)

 Download zip file   Auto-launch 
Help downloading and running models
Accession:184732
" ... Based on experimental observations, we developed a computational model that can be embedded in more comprehensive models of respiratory and cardiovascular autonomic control. Our simulation results successfully reproduce the variability we observed experimentally. The in silico model suggests that age-dependent variability may be due to a developmental increase in mean synaptic conductance between preBötC neurons. We also used simulations to explore the effects of stochastic spiking in sensory relay neurons. Our results suggest that stochastic spiking may actually stabilize modulation of both respiratory rate and its variability when the rate changes due to physiological demand. "
Reference:
1 . Fietkiewicz C, Shafer GO, Platt EA, Wilson CG (2016) Variability in respiratory rhythm generation: In vitro and in silico models Communications in Nonlinear Science and Numerical Simulation 32:158-168
Citations  Citation Browser
Model Information (Click on a link to find other models with that property)
Model Type: Realistic Network;
Brain Region(s)/Organism:
Cell Type(s): Respiratory column neuron;
Channel(s):
Gap Junctions:
Receptor(s):
Gene(s):
Transmitter(s):
Simulation Environment: NEURON;
Model Concept(s): Noise Sensitivity; Development;
Implementer(s):
/* Use run.hoc to start simulation. Creates a raster plot similar to Figure 2C in "Variability in respiratory rhythm generation: in vitro and in silico models", C Fietkiewicz, GO Shafer, EA Platt, CG Wilson, Communications in Nonlinear Science and Numerical Simulation, Volume 32, March 2016, Pages 158–168 (doi: 10.1016/j.cnsns.2015.08.018).
NOTE: To get more bursts, increase the simulation time.
Author contact: chris.fietkiewicz@case.edu
Author website: http://filer.case.edu/cxf47
*/

/*
Parameter descriptions:

"seed" is a seed for random number generators.

tonicToPMweight: Synaptic weight from tonic to PM cells.
tonicToPMprob: Probability of any tonic cell to be connected to a PM cell.
PMtoPMweight: Synaptic weight between PM cells.
PMtoPMprob: Probability of any 2 PM cells to be connected.
tonicNoise: Controls the "randomness" of the tonic spiking. 0 = none. 1 = 100%
tonicPeriodMean: Tonic cell interspike interval.

Comments: The connection weights and probabilities compliment each other. If you 
increase the weight to a PM cell, the input to the cell will increase and the cell
will brust faster. However the input will decrease if you also decrease the 
probability because there will be fewer connections from other cells.
*/

if (setStop < 0) { // For negative sim time...
	setStop = tonicPeriodMean * 1000 // ...estimate sim time necessary for 100 bursts
}

// Other parameters
recordMode = 0 // 0 = Spike times, 1 = Input currents under v-clamp, 2 = Membrane voltages under i-clamp
recordTonic = 0 // Include(1)/exclude(0) tonic cells in spike raster
printParams = 0
weightSeed = seed
connectSeed = seed
cellSeed = seed

load_file("cellTemplates.hoc")

//***************************Network specification interface****************

objref nclist, netcon
objref pmcells, tonicecells
objref ranlist, cells

cells = new List()
nclist = new List()
tonicecells = new List()
pmcells= new List()

func tonicecell_append() {
	cells.append($o1)
    tonicecells.append($o1) 
	return tonicecells.count - 1
}

func pmcell_append() {
	cells.append($o1)
    pmcells.append($o1) 
	return pmcells.count - 1
}

//************************NEW CELLS*********************************
objref randomDelay
randomDelay = new Random(tonicSeed)
randomDelay.uniform(0, tonicPeriodMean)

proc createTonicE() {
    tonicecells.remove_all()
    for i=0, $1-1 {
		tonicPeriod = tonicPeriodMean
		tonicDelay = randomDelay.repick()
		tonicecell_append(new NoisyTonicECells_NetStim(tonicDelay, tonicNoise, tonicSeed, tonicPeriod)) 
	}	
}
proc createPM() {
	pmcells.remove_all()
    for i=0, $1-1 {
		pmcell_append(new PM())
	}     
}
//***********************************CONNECTIONS***************************
objref randomConnection, randomWeight, leakRecord, napRecord  // Random number generator and lists of random values
randomWeight = new Random(weightSeed)
randomConnection = new Random(connectSeed)
randomConnection.uniform(0, 1)
leakRecord = new Vector(numPM)
napRecord = new Vector(numPM)

proc wireTonic() { 
	for i = 0, tonicecells.count()-1 {
		for j = 0, pmcells.count()-1 {
			randomWeight.uniform(tonicToPMweight*(1-tonicToPMrange), tonicToPMweight*(1+tonicToPMrange))
			r = randomConnection.repick()
			w = randomWeight.repick()
			if (r < tonicToPMprob) {
				netcon = new NetCon (tonicecells.object(i).pp, pmcells.object(j).synlist.object(0), -10, 0, w)
				nclist.append(netcon)
			}
		}
	}
}

// Choose random PMtoPMweight (mean value)
variance = PMtoPMrange * PMtoPMrange
randomWeight.normal(PMtoPMweight, variance) // PMtoPMweight is initially the mean of randomly selected weights
PMtoPMweight = randomWeight.repick()  // PMtoPMweight becomes the actual weight as in previous schemes
print "PMtoPMweight = "
print PMtoPMweight
	
// Actual PM wiring
proc wirePM() {
	for i = 0, pmcells.count()-1 {
		for j = 0, pmcells.count()-1 {
			randomWeight.uniform(PMtoPMweight*(1-0.4), PMtoPMweight*(1+0.4))
			r = randomConnection.repick()
			w = randomWeight.repick()
			if (r < PMtoPMprob && i != j) {
				netcon = pmcells.object(i).connect2target(pmcells.object(j).synlist.object(0))
				netcon.threshold = 0
				netcon.delay = 0
				netcon.weight = w   
				nclist.append(netcon)
			}
		}
	}
}


//***********************ADDING HETEROGENIETY*******************************
objref randomCells // Random number generator for cell heterogeneity
randomCells = new Random(cellSeed)

proc randomizeCells () { 
	for i = 0, pmcells.count()-1 {
		// Set gmax_leak
		gLimitFactor = 0.40 // Factor for computing hard limit (+/-)
		gMean = 0.00005
		randomCells.uniform(gMean - gMean * gLimitFactor, gMean + gMean * gLimitFactor) // Set distribution
		r = randomCells.repick() // Pick value
		PM[i].soma gmax_leak = r
		leakRecord.x[i] = r
		// Set gmax_nap
		gLimitFactor = 0.15 // Factor for computing hard limit (+/-)
		gMean = 0.00009
		randomCells.uniform(gMean - gMean * gLimitFactor, gMean + gMean * gLimitFactor) // Set distribution
		r = randomCells.repick() // Pick value
		PM[i].soma gmax_nap = r
		napRecord.x[i] = r
	}
}

//***************** Changing Extracellular K ****************
proc setPotentials() {
	E_NA = 50
	E_K = $1
	E_LEAK = $2
    for i = 0, pmcells.count()-1 {
		PM[i].soma.e_na = E_NA
		PM[i].soma.e_nap = E_NA
		PM[i].soma.e_k = E_K
		PM[i].soma.e_leak = E_LEAK
     }
}

//********************** RECORD *********************
objref netcon, nil
objref spikeTimeVec, spikeIdVec, currentVec[100], vclamp[100], voltageVec[100]
spikeTimeVec = new Vector()
spikeIdVec = new Vector()
nRecord = 0 // Total number of cells. Used for raster plot.
vclampQty = 5 // Quantity of cells for vclamp recording
iclampQty = 5 // Quantity of cells for vclamp recording

proc prepareRecording() {
	if (recordMode == 0) {
		if (recordTonic == 1) {
			// Record tonic E cells
			for i=0, tonicecells.count - 1 {
				netcon = new NetCon(tonicecells.object(i).pp, nil, 0, 1, 0)    
				netcon.record(spikeTimeVec, spikeIdVec, nRecord+i+1)
			}
			nRecord = nRecord + tonicecells.count
		}
		// Record PM cells
		for i=0, pmcells.count()-1 {
			pmcells.object(i).soma netcon = new NetCon(&v(0.5), nil, 0, 1, 0)    
			netcon.record(spikeTimeVec, spikeIdVec, nRecord+i+1)
		}
		nRecord = nRecord + pmcells.count()
	}
	if (recordMode == 1) {
		for i=0, vclampQty-1 {
			// Create v-clamps
			pmcells.object(i).soma vclamp[i] = new SEClamp(0.5)
			vclamp[i].dur1 = setStop
			vclamp[i].rs = 0.01
			vclamp[i].amp1 = -71

			// Set up recordings
			currentVec[i] = new Vector()
			currentVec[i].record(&vclamp[i].i)
		}
	}
	if (recordMode == 2) {
		for i=0, iclampQty-1 {
			// Set up recordings
			voltageVec[i] = new Vector()
			voltageVec[i].record(&PM[i].soma.v(0.5))
		}
	}
}

//****************File Output*************
objref f, state
objref graster

proc writeFiles() {  
	if (recordMode == 0) {
		// Save i-clamp spike times
		f = new File()
		strdef filename
		sprint(filename, "spikeTime_%s.txt", $s1)
		f.wopen(filename)
		spikeTimeVec.printf(f)
		f.close()

		// Save i-clamp spike IDs
		f = new File()
		sprint(filename, "spikeID_%s.txt", $s1)
		f.wopen(filename)
		spikeIdVec.printf(f)
		f.close()
	}
	if (recordMode == 1) {
		strdef fileName
		f = new File()
		for i=0, vclampQty-1 {
			sprint(fileName, "in%d_%s.bin", i + 1, $s1)
			f.wopen(fileName)
			currentVec[i].vwrite(f, 3)
			f.close()
		}
	}
	if (recordMode == 2) {
		strdef fileName
		f = new File()
		for i=0, iclampQty-1 {
			sprint(fileName, "out%d_%s.bin", i + 1, $s1)
			f.wopen(fileName)
			voltageVec[i].vwrite(f, 3)
			f.close()
		}
	}
}

proc writeParams() {
	if (printParams == 1) {
		f = new File()
		sprint(filename, "params_%s.txt", $s1)
		f.wopen(filename)
		leakRecord.printf(f)
		napRecord.printf(f)
		f.close()
	}
}

//*************************INITIALIZATION************************
objref f, state

proc init() {
	celsius = 27
	prepareRecording()
	finitialize(-70)
}

//**************** Setup and Run ****************
createTonicE(numTonic)
createPM(numPM)
randomizeCells()
setPotentials(setEK, setEleak)
wireTonic()
wirePM()
writeParams(fileSuffix)
init()
batch_run(setStop, 0.1)

//**************** Output ****************
if (outputMode == 0) {
	writeFiles(fileSuffix)
		quit()
} else { // Raster plot
	graster = new Graph(0)
	graster.view(0 ,0 , setStop, nRecord, 0, 0, 900, 600)
	spikeIdVec.mark(graster, spikeTimeVec, "|", 6)
}