2a:null 2b:["$","title",null,{"children":"Information Transfer and Decoders | Steve Oney | Observable"}] 30:T706,We can see that while our signal is high, our neuron fires. When our signal is low, it doesn't fire. So how can we try to re-construct our original signal based just on the firing rate? Suppose we have a firing rate ${tex`r(I)`} as a function of input ${tex`I`}. So let's define a weighting term ${tex`\phi{}`} that we can multiply ${tex`r(I)`} by in order to estimate our signal. We will estimate our signal as ${tex`\phi{}r(I)`}. We can calculate the mean squared error (MSE) of our estimate as the actual input (${tex`I`}) minus the estimate: ${tex`E(I) = \frac{1}{2}[I - r(I)\phi{}]^2`} ...but really, we want to minimize our error across multiple values of ${tex`I`} (not just one). So we choose an interval over which we want to minimize error. Let's select ${tex`[-1, 1]`}, so we ultimately can calculate the error as: ${tex` E = \int_{-1}^{1} \frac{1}{2}[x - r(x)\phi{}]^2 \,dx `} To minimize this, we want to find the lowest point of the error curve, which is where ${tex`\frac{dE(x)}{d\phi} = 0`}: ${tex` \frac{dE(x)}{d\phi} = \int_{-1}^{1} -2r(x)\frac{1}{2}[x - r(x)\phi{}] \,dx \newline \frac{dE(x)}{d\phi} = \int_{-1}^{1}[-xr(x) + \phi{}r(x)^2]dx `} If we set this equal to zero, we have: ${tex` \int_{-1}^{1} xr(x)dx = \int_{-1}^{1}\phi{}r(x)^2dx `} Since ${tex`\phi{}`} does not vary with respect to input ${tex`x`}, we can take it outside of the integral and compute its value as: ${tex` \phi{} = \frac{\int_{-1}^{1} xr(x)dx}{\int_{-1}^{1}r(x)^2dx} `} For convenience, we will assign the variables ${tex`\upsilon{} = \int_{-1}^{1}xr(x)dx`} and ${tex`\gamma{} = \int_{-1}^{1}r(x)^2dx`} so we end up with: ${tex` \phi{} = \gamma{}^{-1}\upsilon{} `} We can then estimate the values of ${tex`\upsilon{}`} and ${tex`\gamma{}`} by computing them over several fixed intervals.31:T444,## Neural Populations (Multiple Neurons) So far, we have focused on the dynamics of single neurons but in order to convey information effectively, we need **populations** of neurons---not just one. In other words, we are going to move from using one neuron to using many neurons to better represent (and later on, process) our signals. If we work with a population of neurons, each neuron should have a different tuning curve (if two neurons had the exact same tuning curve, they aren't conveying any more information than one neuron). So each neuron should react to the same input in different ways. This should convey more information than just one neuron, which we can verify by decoding the signals of multiple neurons. So we are going to (1) start using the subscript ${tex`_{i}`} to represent neuron ${tex`i`} and (2) add two additional terms that scale the input of every neuron: a multiplicative factor, ${tex`\alpha_i`} and an additive bias ${tex`b_i`}. This means that given an input signal ${tex`x`}, neuron ${tex`i`} actually reacts as if the input were ${tex`\alpha_ix + b_i`}.32:T7e9,## Computing Scales and Biases (Optional) We could randomly select values of ${tex`\alpha_i`} and ${tex`b_i`} but let's instead compute them from: - A desired maximum firing rate in the context of an interval, which we will call ${tex`r_{max}`} - A desired intercept, where our neuron starts firing, which we will call ${tex`x_{int}`} [We can analytically determine](https://observablehq.com/d/af48ade400490b47?collection=@soney/neuromorphic-computing#cell-216) ${tex`r(I)`} to estimate a firing rate: ${tex` r(I) = \begin{cases} 0, & \text{if } I \leq v_{th}{}\\ \frac{1}{\tau_{ref}-\tau_{rc}\ln{\left(1 - \frac{v_{th}}{I}\right)}}, & \text{otherwise} \end{cases} `} We can turn this equation around to ask what **input current** do we need to achieve a given firing rate. In other words, we can leave ${tex`r(I)`} fixed (at ${tex`r`}) and solve for ${tex`I`}. If we do this on the above equation, we get: ${tex` I = \frac{v_{th}}{1 - e^{\frac{\tau_{ref} - \frac{1}{r}}{\tau_{rc}}}} `} Let's call this ${tex`z(r)`} for short: ${tex` z(r) = \frac{v_{th}}{1 - e^{\frac{\tau_{ref} - \frac{1}{r}}{\tau_{rc}}}} `} Recall that our effective input given input ${tex`x`} is ${tex`I`}: ${tex` I(x) = \alpha x + b `} Let's suppose that we only care about the possible inputs up to ${tex`1`}. This means that our maximum rate, ${tex`r_{max}`} should be achieved at ${tex`x=1`}: ${tex` I_i(1) = 1\alpha_i + b = z(r_{max}) `} or re-arranging: ${tex` b = z(r_{max}) - \alpha `} But we know that at our intercept ${tex`x_{int}`}, we should not fire, meaning that it is **right** at the threshold of our firing voltage, ${tex`v_{th}`}: ${tex`I(x_{int}) = v_{th}`}. So we have: ${tex` I(x_{int}) = v_{th} = \alpha{}x_{int} + b `} or re-arranging: ${tex` b = v_{th} - \alpha{}x_{int} `} Putting both equations for ${tex`b`} together, we get: ${tex` z(r_{max}) - \alpha = v_{th} - \alpha{}x_{int} `} ...re-arranging: ${tex` \alpha = \frac{v_{th} - z(r_{max})}{x_{int} - 1} `} and ${tex`b = v_{th} - \alpha{}x_{int}`}.33:T6ee,## Computing Population Decoders (Optional) When we want to compute the decoders for **multiple** neurons, the procedure is similar as it was for one neuron. However, now we are minimizing the error across several neurons. So now, if we have ${tex`n`} neurons, our decoders ${tex`\vec{\phi}`} will be a vector rather than a single number. We use ${tex`\phi_i`} to refer to the decoder for neuron ${tex`i`}. ${tex` E = \int_{-1}^{1} \frac{1}{2}[x - \sum_{i=1}^{n} r_i(x)\phi_i]^2 \,dx `} We could minimize this procedurally (for example, through gradient descent) or analytically. For large numbers of neurons, it's probably best to do it procedurally but let's instead do it analytically. We want to find the optimal values for each ${tex`\phi_i`} by looking for the point where ${tex`\frac{\partial{}E(j)}{\partial{}\phi_i} = 0`}. ${tex` \frac{\partial{}E(x)}{\partial{}\phi_i} = \int_{-1}^{1} -2r_i(x)\frac{1}{2}[x - \sum_{j=1}^{n} r_j(x)\phi_j] \,dx `} If we set this equal to zero, we have: ${tex` \int_{-1}^{1} xr_i(x)dx = \int_{-1}^{1}[r_i(x)\sum_{j=1}^{n}r_j(x)\phi_j] \,dx `} re-arranging ${tex` \int_{-1}^{1} xr_i(x)dx = \sum_{j=1}^{n}[\int_{-1}^{1}r_i(x)r_j(x) \,dx]\phi_i `} If we define matrices ${tex`\Upsilon{}`} and ${tex`\Gamma{}`} as follows: ${tex` \Gamma_{ij} = \int r_i(x)r_j(x)dx `} and ${tex` \Upsilon_i = \int xr_i(x)dx `} then we get ${tex` \Upsilon = \Gamma \phi `} or ${tex` \phi = \Gamma^{-1} \Upsilon `} Note that [some sources](https://link.springer.com/article/10.1007/BF00962720) recommend adding random noise to be more sure that ${tex`\Gamma`} is invertible. They recommend adding ${tex`\sigma{}`}, the standard deviation of the noise to the diagonal so we will use ${tex`\Gamma + I`} where ${tex`I`} is the identity matrix.34:T505,function plotDecoder2(lifs, J, tau_s, runtime=2) { const color = d3.scaleOrdinal(d3.schemeCategory10); const TIME_STEP = 1e-3; const NUM_STEPS = runtime / TIME_STEP; // Simulate for runtime seconds const phis = getDecoders(lifs); const psps = lifs.neurons.map((lif) => psp(lif, tau_s)); const v_history = []; const historyByNeuron = []; for(let i = 0; i psp.step(j, TIME_STEP)); let out = 0; for(let k = 0; k < lifs.neurons.length; k++) { const output_k = phis[k]*psps[k].output; historyByNeuron[k].push({t, out: output_k}); out += output_k; } v_history.push({t, out, j}); } const marks = historyByNeuron.map((neuronData, i) => { return Plot.line(neuronData, {x: "t", y: "out", stroke: color(i), strokeWidth: 1, opacity: 0.1}); }); const plot = Plot.plot({ y: { grid: true, label: 'j', domain: [-1, 1] }, marks: marks.concat([ Plot.line(v_history, {x: "t", y: "out", stroke: 'black', strokeWidth: .5}), Plot.line(v_history, {x: "t", y: "j", stroke: 'red', strokeWidth: 1}) ]) }); return plot; }35:[] 2e:{"id":"e363ce6cdb928f48","slug":"nmc-information-transfer-and-decoders","title":"Information Transfer and Decoders","thumbnail":"eac39d6a65aa39d29594bbcdecbb61d86b40c606779712c9d4e2906aedb3db11","access_level":"public","create_time":"$D2022-10-24T14:36:45.468Z","update_time":"$D2024-07-18T15:17:10.731Z","view_time":null,"first_public_time":"$D2022-12-07T15:00:10.977Z","fork":false,"stars":0,"starred":false,"owner":{"type":"team","id":"55b87b37da07de37","login":"soney","name":"Steve Oney","avatar_url":"https://avatars-next.observableusercontent.com/565926373d73cf99aa6425cf7c0d7ff79e41635d1a367894a5e3bce0299b18b9","tier":"starter_2024"},"version":1038,"publish_version":null,"role":null,"authored":false,"authors":[{"type":"individual","id":"1d138075b46297f1","login":"soney","name":"Steve Oney","avatar_url":"https://avatars-next.observableusercontent.com/565926373d73cf99aa6425cf7c0d7ff79e41635d1a367894a5e3bce0299b18b9","tier":null}],"members":[{"type":"individual","id":"79581e18d086caf4","login":"ericgallo","name":"","avatar_url":null,"tier":null,"role":"editor"},{"type":"individual","id":"93113dd5fb467076","login":"sharvari-bondre","name":"Sharvari Bondre","avatar_url":"https://avatars-next.observableusercontent.com/754b967a5cf1ce71589f615fed536b94bd68509f3412154d82572c4b6aa44ead","tier":null,"role":"editor"},{"type":"individual","id":"1d138075b46297f1","login":"soney","name":"Steve Oney","avatar_url":"https://avatars-next.observableusercontent.com/565926373d73cf99aa6425cf7c0d7ff79e41635d1a367894a5e3bce0299b18b9","tier":null,"role":"editor"}],"comments":[{"id":"1e494e895830d70d","reply_to_id":null,"content":"Add explanation of why breaking into multiple neurons (representation)","create_time":"$D2023-01-24T17:07:08.929Z","update_time":null,"resolved":true,"user":{"type":"individual","id":"1d138075b46297f1","login":"soney","name":"Steve Oney","avatar_url":"https://avatars-next.observableusercontent.com/565926373d73cf99aa6425cf7c0d7ff79e41635d1a367894a5e3bce0299b18b9","tier":null}},{"id":"f0d04ee600b9881d","reply_to_id":null,"content":"Why do you want to define the maximum firing rate? I thought the maximum firing rate should be defined by \\tau?","create_time":"$D2023-05-07T21:48:37.830Z","update_time":null,"resolved":false,"user":{"type":"individual","id":"a2f9f861fa76becb","login":"a2f9f861fa76becb","name":"Ziqiao Xi","avatar_url":"https://avatars-next.observableusercontent.com/64861e404980b0ca88d5d40b5860ed28e9be4296c77b75d590ccb4e9e6535319","tier":null}},{"id":"646e3035511aed64","reply_to_id":null,"content":"typo: exapt -> except","create_time":"$D2023-06-11T15:07:43.548Z","update_time":null,"resolved":true,"user":{"type":"individual","id":"4132f3acf2c585f1","login":"feiwu","name":"Fei Wu","avatar_url":"https://avatars-next.observableusercontent.com/786fbb8fd38ae9919783bb4368af88fefa862d32ba0a267d664efe2981392be7","tier":null}},{"id":"e4ee32aacbc1424a","reply_to_id":null,"content":"In Priya Panda's talk \"Spiking Neural Networks: Algorithms\", she talks about two coding techniques: \n1. temporal coding: higher pixel values spike first\n2. rate coding: higher pixel values spike more number of times (used more frequently now)\nHere, you are talking about \"rate-coded\". However, in module 4 \"Learning through Spike-Timing Dependent Plasticity (STDP)\", it seems you are using \"temporal-coded\" at paragraph 2. Without related knowledge, readers may think they are the same.","create_time":"$D2023-06-11T16:10:04.230Z","update_time":null,"resolved":false,"user":{"type":"individual","id":"4132f3acf2c585f1","login":"feiwu","name":"Fei Wu","avatar_url":"https://avatars-next.observableusercontent.com/786fbb8fd38ae9919783bb4368af88fefa862d32ba0a267d664efe2981392be7","tier":null}}],"custom_thumbnail":false,"trashed":false,"subscribed":null,"forks":1,"collection":{"id":"ceab76d8d8a0acbc","slug":"neuromorphic-computing","title":"Neuromorphic Tutorials: Spiking Neural Networks (SNNs) from Scratch","owner":{"id":"55b87b37da07de37","login":"soney","name":"Steve Oney"}},"type":"notebook","body":{"title":"Information Transfer and Decoders","cells":[{"id":0,"value":"# Information Transfer and Decoders","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":5,"value":"Neurons' spiking patterns contain **information** about the stimulus that is causing them to fire. However, the form of this information (spike trains---series of spikes over time) can be difficult to interpret. We (as SNN authors) might want to know if these spikes are actually conveying the information that we care about or not. In order to do that, we need to **decode** the series of spikes to re-estimate the original message.\n\n\nThere are many theories of how information is represented in spike trains. Near the periphery, however, it appears that information is \"rate-coded\"---the firing rate of a spike train represents some information about the underlying stimulus. Let's take the example where we have an input in the form of an oscillating sin wave:","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":11,"value":"function signal(t) {\n return 2*Math.sin(t);\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":13,"value":"...and suppose that we have the same LIF neuron model as before.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":15,"value":"class LIF {\n constructor(tau_rc=0.2, v_init=0, tau_ref=0.002, v_th=1) {\n this.tau_rc = tau_rc;\n this.tau_ref = tau_ref;\n this.v = v_init;\n this.v_th = v_th;\n \n this.output = 0;\n this.refractory_time = 0;\n }\n\n step(I, t_step) {\n this.refractory_time -= t_step; // Subtract the amount of time thhat passed from our refractory time\n \n if(this.refractory_time < 0) { // If we aren't in our refractory period, update the voltage\n this.v = this.v * (1 - t_step/this.tau_rc) + (I * t_step / this.tau_rc);\n }\n\n if(this.v >= this.v_th) { // Fired\n this.refractory_time = this.tau_ref; // Set refractory time\n\n this.output = 1 / t_step; // Fire\n this.v = 0; // Reset potential\n } else {\n this.output = 0; // Don't fire\n }\n return this.output;\n }\n \n reset() {\n this.v = this.output = this.refractory_time = 0;\n }\n}\n","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":21,"value":"...and the same synapse model as before","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":23,"value":"function psp(neuron, tau_s=1e-2) {\n return {\n output: 0,\n step: function (t, t_step) {\n neuron.step(t, t_step);\n this.output = this.output * (1 - t_step / tau_s) + neuron.output * t_step / tau_s;\n }\n };\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":47,"value":"So we have our incoming signal, which is a scaled sin wave...","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":44,"value":"plotSignal(signal, Math.PI*4)","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":49,"value":"...and we have our output for that signal for a single neuron.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":28,"value":"plotPotential(signal, 1e-1, Math.PI*4)","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":51,"value":"$30","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":109,"value":"function activity_level(i, lif) {\n const T_STEP = 1e-3;\n const T = 0.2;\n let spikeCount = 0;\n lif.reset();\n for(let t = 0; t<=T; t+=T_STEP) {\n const output = lif.step(i, T_STEP);\n if(output > 0) {\n spikeCount++;\n }\n }\n return spikeCount / T;\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":123,"value":"function computeGamma(lif, range_low, range_high, interval) {\n let gamma = 0;\n const di = interval;\n for(let i = range_low; i<=range_high; i+=interval) {\n gamma += di * Math.pow(activity_level(i, lif), 2);\n }\n\n return gamma;\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":126,"value":"function computeUpsilon(lif, range_low, range_high, interval) {\n let upsilon = 0;\n const di = interval;\n for(let i = range_low; i<=range_high; i+=interval) {\n upsilon += di * i * activity_level(i, lif);\n }\n return upsilon;\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":135,"value":"function computePhi(lif, range_low=-1, range_high=1, interval=1e-2) {\n return computeUpsilon(lif, range_low, range_high, interval) / computeGamma(lif, range_low, range_high, interval);\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":133,"value":"function plotDecoder(I, tau_s, runtime=2) {\n const TIME_STEP = 1e-3;\n const NUM_STEPS = runtime / TIME_STEP; // Simulate for runtime seconds\n\n const neuron = new LIF(/*0.02, 0, 0.002*/);\n const postSynaptic = psp(neuron, tau_s);\n\n const phi = computePhi(neuron, -2, 2);\n const v_history = [ {t: 0, v: neuron.v} ];\n for(let i = 0; i < NUM_STEPS; i++) {\n const inp = I((i+1)*TIME_STEP);\n postSynaptic.step(inp, TIME_STEP);\n v_history.push({t: (i+1)*TIME_STEP, decoded: postSynaptic.output * phi, inp });\n }\n\n const plot = Plot.plot({\n y: {\n grid: true,\n label: 'inp'\n },\n marks: [\n Plot.line(v_history, {x: \"t\", y: \"inp\", stroke: \"red\", strokeWidth: 2}),\n Plot.line(v_history, {x: \"t\", y: \"decoded\", stroke: \"black\", strokeWidth: 1}),\n ]\n });\n return plot\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":148,"value":"plotDecoder(signal, 1e-1, Math.PI*4)","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":179,"value":"We can see that our estimate *kind of* captures some of the highs of our input signal but not much else. Instead of just having **one** neuron, let's see how we can instead represent our signal with multiple neurons to get a better result.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":193,"value":"$31","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":660,"value":"We will define a function `getGainBias` that accepts four arguments---${tex`\\tau_{rc}`}, ${tex`\\tau_{ref}`}, ${tex`v_{th}`}, ${tex`r_{max}`}, and ${tex`x_{int}`}.\n\n${tex`\\tau_{rc}`}, ${tex`\\tau_{ref}`}, and ${tex`v_{th}`} are all neuron parameters that we've seen before. ${tex`r_{max}`} represents the maximum firing rate we want our neuron to have if the input (${tex`x`}) is in the range ${tex`[-1, 1]`} (in Hz). We will assign it randomly in the range ${tex`[200, 400]`}. ${tex`x_{int}`} represents the input (${tex`x`}) where we want our neuron to start firing. We will assign it randomly in the range ${tex`[-1, 1]`}. Our function will return a value for ${tex`\\alpha`} and ${tex`b`} given these five values. [The optional material at the end of the tutorial](https://observablehq.com/d/e363ce6cdb928f48?collection=@soney/neuromorphic-computing#cell-221) explains how `getGainBias` is implemented.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":276,"value":"function getGainBias(tau_rc, tau_ref, thresh, r_max, x_int) {\n const z_rmax = thresh / (1 - Math.exp((tau_ref - 1/r_max)/tau_rc) );\n\n const gain = (thresh - z_rmax) / (x_int - 1);\n const bias = thresh - gain * x_int;\n return [gain, bias];\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":675,"value":"Now, we need to compute the decoders for **multiple** neurons. The procedure is similar to before, except that now we are minimizing the error across several neurons. So now, if we have ${tex`n`} neurons, our decoders ${tex`\\vec{\\phi}`} will be a vector rather than a single number. We use ${tex`\\phi_i`} to refer to the decoder for neuron ${tex`i`}. We will implement code to find ${tex`\\vec{\\phi}`} but [the optional materials at the end of this notebook](https://observablehq.com/d/e363ce6cdb928f48?collection=@soney/neuromorphic-computing#cell-399) detail how it is derived.\n\nNow, we will create a class `NeuronCollection` that contains multiple neurons. It will contain instance variables `this.neurons` to represent an array of its neurons and `this.output` to represent the output of each of those neurons.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":940,"value":"class NeuronCollection {\n constructor(createNeuronFn, n) { // createNeuronFn is a function to create a neuron instance. n is the number of neurons\n this.neurons = []; // Neuron instances\n\n // Create n neurons\n for(let i = 0; i neuron.output); // Assign our output to be a list of neuron outputs\n }\n\n step(inputs, t_step) {\n // Step every individual neuron\n // Assign this.output to the neurons' outputs\n // Return this.output\n return this.output = this.neurons.map((neuron, i) => neuron.step(inputs[i], t_step));\n }\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":944,"value":"Now, we will create a `NeuronContainer` with random values for ${tex`\\alpha`} and ${tex`b`}. We will select intercepts randomly in the range ${tex`[-1, 1]`} and maximum firing rates in the range ${tex`[200, 400]`}.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":299,"value":"function getLIFCollection1(n, tau_rc=0.02, tau_ref=0.002, thresh=1, max_rate_range = [200, 400], intercept_range = [-1, 1]) {\n return new NeuronCollection(() => {\n const maxRate = randomBetween(max_rate_range[0], max_rate_range[1]);\n const intercept = randomBetween(intercept_range[0], intercept_range[1]);\n const [gain, bias] = getGainBias(tau_rc, tau_ref, thresh, maxRate, intercept);\n\n const lif = new LIF(tau_rc, 0, tau_ref, thresh);\n lif.step = (I, t_step) => LIF.prototype.step.call(lif, I * gain + bias, t_step);\n lif.gain = gain;\n lif.bias = bias;\n return lif;\n }, n);\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":314,"value":"Let's look at the tuning curves for these neurons:","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":316,"value":"{\n const color = d3.scaleOrdinal(d3.schemeCategory10);\n const [minJ, maxJ] = [-1, 1];\n const stepSize = .01;\n\n const lifs = getLIFCollection1(10);\n\n const marks = lifs.neurons.map((lif, i) => {\n const data = [];\n for(let j = minJ; j<=maxJ; j+=stepSize) {\n data.push({j, analytical: analytical_activity_level(j, lif), empirical: activity_level(j, lif)});\n }\n return [Plot.line(data, {x: \"j\", y: \"analytical\", stroke: color(i), strokeWidth: 0.5}), Plot.line(data, {x: \"j\", y: \"empirical\", stroke: color(i), strokeWidth: 2})];\n });\n\n const plot = Plot.plot({\n y: {\n grid: true,\n label: 'hz'\n },\n marks\n });\n return plot\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":219,"value":"### Negative Scaling\nWe are going to make one more change to allow ${tex`\\alpha_i`} to be *negative*, meaning that our neuron might fire with negative input values. Although this does not necessarily make sense in the context of a single isolated neuron, (1) biological realism is not our primary goal and (2) it might correspond to a situation where there is some constant background current and an inhibitory background current.\n\nWe will call this our **encoder**. We will call it ${tex`e`} It will randomly be either ${tex`1`} or ${tex`-1`} and we will multiply it by ${tex`\\alpha_i`}:\n\n${tex`\nI(x) = e\\alpha x + b\n`}","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":346,"value":"function getLIFCollection(n, tau_rc=0.02, tau_ref=0.002, thresh=1, max_rate_range = [200, 400], intercept_range = [-1, 1]) {\n return new NeuronCollection(() => {\n const maxRate = randomBetween(max_rate_range[0], max_rate_range[1]);\n const intercept = randomBetween(intercept_range[0], intercept_range[1]);\n\n const [gain, bias] = getGainBias(tau_rc, tau_ref, thresh, maxRate, intercept);\n const encoder = randomItem([-1, 1]);\n\n const lif = new LIF(tau_rc, 0, tau_ref, thresh);\n lif.step = (I, t_step) => LIF.prototype.step.call(lif, I * gain * encoder + bias, t_step);\n lif.gain = gain;\n lif.bias = bias;\n lif.encoder = encoder;\n return lif;\n }, n);\n}","pinned":true,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":685,"value":"Now, if we plot our tuning curves, we should see some that go in the \"opposite\" direction (as in they fire more with negative inputs).","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":694,"value":"lifs = getLIFCollection(10)","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":351,"value":"{\n const color = d3.scaleOrdinal(d3.schemeCategory10);\n const [minJ, maxJ] = [-1, 1];\n const stepSize = .01;\n \n const marks = lifs.neurons.map((lif, i) => {\n const data = [];\n for(let j = minJ; j<=maxJ; j+=stepSize) {\n data.push({j, analytical: analytical_activity_level(j, lif), empirical: activity_level(j, lif)});\n }\n return [ Plot.line(data, {x: \"j\", y: \"analytical\", stroke: color(i), strokeWidth: 0.5}), Plot.line(data, {x: \"j\", y: \"empirical\", stroke: color(i), strokeWidth: 2}) ];\n });\n\n const plot = Plot.plot({\n y: {\n grid: true,\n label: 'hz'\n },\n marks: marks.flat()\n });\n return plot\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":981,"value":"So let's use these neuron collections to start representing different signals. For convenience, we'll also create a class `PSP` to represent the post-synaptic potential of a `NeuronCollection`. The implementation is similar to that of [a single neuron](https://observablehq.com/d/3c638d3aaf9207b4). Our `PSP` class will accept a `NeuronCollection` instance and a time constant (`tau`) as arguments.","pinned":false,"hidden":false,"mode":"md","author":"$undefined","summary":"$undefined"},{"id":974,"value":"class PSP {\n constructor(collection, tau) {\n this.collection = collection; // NeuronCollection instance\n this.tau = tau; // time constant for decay\n this.output = this.collection.output.map(o => o); // Use the collection's output (but make a clone so we don't mutate the original)\n }\n\n step(t_step) {\n // Loop through my output, decay it, and add the output from the incoming collection.\n for(let i = 0; i 0) {\n numSpikes++;\n }\n v_history.push({t: (i+1)*TIME_STEP, v: neuron.v, j, out });\n }\n const plot = Plot.line(v_history, {x: \"t\", y: \"out\"}).plot();\n return plot;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":42,"value":"function plotSignal(J, runtime=2) {\n const TIME_STEP = 1e-3;\n const NUM_STEPS = runtime / TIME_STEP; // Simulate for runtime seconds\n\n const signal_history = [ {t: 0, s: J(0)} ];\n\n let numSpikes = 0;\n for(let i = 0; i < NUM_STEPS; i++) {\n signal_history.push({t: (i+1)*TIME_STEP, s: J((i+1)*TIME_STEP) });\n }\n const plot = Plot.line(signal_history, {x: \"t\", y: \"s\"}).plot();\n return plot;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":355,"value":"function plotPotentials(lifs, J, tau_s, runtime=2) {\n const color = d3.scaleOrdinal(d3.schemeCategory10);\n\n const TIME_STEP = 1e-3;\n const NUM_STEPS = (runtime / TIME_STEP) + 1; // Simulate for runtime seconds\n\n const psp = new PSP(lifs, tau_s);\n\n const v_history = [];\n const historyByNeuron = [];\n for(let i = 0; i {\n historyByNeuron[k].push({t, out});\n });\n }\n\n const marks = historyByNeuron.map((neuronData, i) => {\n return Plot.line(neuronData, {x: \"t\", y: \"out\", stroke: color(i), strokeWidth: 2});\n });\n const plot = Plot.plot({\n y: {\n grid: true,\n label: 'j'\n },\n marks\n });\n return plot;\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":287,"value":"function randomBetween(low, high) {\n return low + sRandom() * (high - low);\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":348,"value":"function randomItem(l) {\n return l[Math.floor(sRandom() * l.length)];\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":289,"value":"function repeat(func, times) {\n const result = [];\n for(let i = 0; i computeTuningCurve(lif));\n const As = tuningCurves.map((tc) => tc[1]);\n const A = _math.matrix(As);\n const Gamma = _math.add(_math.identity(lifs.neurons.length), _math.multiply(A, _math.transpose(A)));\n const GammaInv = _math.inv(Gamma);\n\n const js = tuningCurves[0][0];\n const Upsilon = _math.multiply(_math.matrix(js), _math.transpose(A))\n\n const Phi = _math.multiply(GammaInv, Upsilon);\n\n return Phi.toArray();\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":485,"value":"$34","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":913,"value":"getSeededRandom = function(s){\n let seed = s;\n return function random() {\n let x = Math.sin(seed++) * 10000;\n return x - Math.floor(x);\n }\n}","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":917,"value":"sRandom = getSeededRandom(100)","pinned":false,"hidden":false,"mode":"ojs","author":"$undefined","summary":"$undefined"},{"id":986,"value":"function times(n, funcVal) {\n const result = [];\n const func = (typeof funcVal === 'function') ? funcVal : (() => funcVal);\n\n for(let i = 0; i