Source code that demonstrates the solution

Here is the source code for two files that, taken together, demonstrate the solution. I suggest that you copy this to a couple of text files that you save on your own computer, then print out a hard copy of each to examine while reviewing comments on the implementation of the solution.

File 1: ranstream.hoc

This file defines the RandomStream class. It comes straight out of the NEURON/common directory in the source code for entry 83319 from ModelDB http://modeldb.yale.edu/

There are two ways to produce independent random streams with MCellRan4:

1. Give each Random instance a highindex that differs from all the other Random instances' highindexes by an amount that is at least as large as the number of events that its NetStim will generate.
2. Give each Random instance a lowindex that differs from all other Random instances' lowindices.

The code in ranstream.hoc uses the first approach, but it would be easy to change so that it uses the second approach.

random_stream_offset_ = 1000

begintemplate RandomStream
public r, repick, start, stream
external random_stream_offset_
objref r
proc init() {
	stream = $1
	r = new Random()
	start()
}
func start() {
	return r.MCellRan4(stream*random_stream_offset_ + 1)
}
func repick() {
	return r.repick()
}
endtemplate RandomStream

File 2: init.hoc

This file shows how to take advantage of the RandomStream class in your own programs, so that your NetStims can generate independent random streams produced by the MCellRan4 generator. For the most part, this file is quite similar to initn.hoc, which demonstrated one of the problems that may occur when NetStims share a common random number generator. The next section in this tutorial--"Discussion of solution source code"--will focus on just those aspects of init.hoc that deserve special comment.

NSNUM = 5 // how many NetStims to create
MAXSTREAM = 1000 // max # events in a NetStim's stream
  // before it begins to repeat values already generated
  // by a different stream.
  // set to 0 and all NetStims will produce identical streams

ISI = 10 // default NetStim parameters
NUM = 50
START = 1
NOISE = 0

TSTOP = 20 // length of simulation

load_file("nrngui.hoc")
load_file("ranstream.hoc")

///// model setup

objref nslist, rslist
nslist = new List()
rslist = new List()

proc makenetstims() { local i  localobj ns, rs
  nslist = new List()
  rslist = new List()
  random_stream_offset_ = MAXSTREAM
  for i = 0, $1-1 {
    ns = new NetStim()
    nslist.append(ns)
    rs = new RandomStream(i)
    rslist.append(rs)
    ns.noiseFromRandom(rs.r)
    rs.r.negexp(1) // must specify negexp distribution with mean = 1
    rs.start()
  }
}

makenetstims(NSNUM)

proc setparams() { local i
  for i = 0,nslist.count()-1 {
    nslist.o(i).interval = ISI
    nslist.o(i).number = NUM
    nslist.o(i).start = START
    nslist.o(i).noise = NOISE
  }
}

setparams() // default is deterministic spike times (NOISE = 0)

///// instrumentation

objref tvec, idvec, nil
tvec = new Vector()
idvec = new Vector()

proc instrumentation() { local i  localobj nc
  tvec = new Vector()
  idvec = new Vector()
  for i = 0,nslist.count()-1 {
    nc = new NetCon(nslist.o(i), nil)
    nc.record(tvec, idvec, i)
  }
}

instrumentation()

///// simulation control

// spike times will be at 0.5*ISI + sample from negexp that has mean 0.5*ISI
NOISE = 0.5
setparams()

tstop = TSTOP

proc restart() { local i
  for i = 0, rslist.count()-1 rslist.o(i).start()
}

objref g

proc myrun() {
  run()
  // make a raster plot of spike times
  idvec.add(1) // so cell i spikes are marked at y = i+1
  objref g
  g = new Graph(0)
  g.size(0,tstop,0,NSNUM)
  g.view(0, 0, tstop, NSNUM, 47, 109, 300.48, 200.32)
  idvec.mark(g, tvec, "|", 10, 2, 1)
  // print out spike times
  print " time         cell"
  for i=0, tvec.size-1 printf("%7.3f \t%d\n", tvec.x(i), idvec.x(i))
}

/////

print "============================================================"
print "Control:  all have interval 10, number 50, start 10, noise 1"
print "============================================================"

myrun()

print " "
print "Take a good look at the graph.  The marks at y = 1"
print "indicate the spikes generated by NetStim[0]."
print "When ready for next test, type"
print "test()"
print "at the oc> prompt, then press Return"
print " "

proc test() {
  print "============================================================"
  print "Test:  NetStim 0 (i.e. cell 1) interval 1, number 3"
  print "============================================================"

  nslist.o(0).interval = 1
  nslist.o(0).number = 3
  restart() // restore each generator to the start of its stream

  myrun()

  print " "
  print "Note that only the marks at y = 1 changed location in time."
  print "Compare the NetStim spike times in this run"
  print "against the spike times from the control run."
  print "Note that only the spikes from \"cell 1\" (NetStim[0]) are different."
}