diff --git a/README.md b/README.md index 92a7464..bb5eb48 100644 --- a/README.md +++ b/README.md @@ -55,14 +55,14 @@ Evaluate the model before and after training: Optional: Maybe evaluate in depth - for ((i=0;i<20;i+=2)); do + for ((i=1;i<2;i+=2)); do echo "Evaluating at $i" py neurosim/main.py eval $WDIR --resume_tidx=$i done Run all evaluation: - WDIR=results/20210713 + WDIR=results/20210805 py neurosim/tools/evaluate.py frequency $WDIR --timestep 10000 py neurosim/tools/evaluate.py medians $WDIR py neurosim/tools/evaluate.py rewards $WDIR diff --git a/config.json b/config.json index 8005743..cc0edbb 100644 --- a/config.json +++ b/config.json @@ -4,9 +4,9 @@ "savemp4": 0, "render": 0, "observation_map": [ - {"type": "rf_normal", "bins": 5, "mean": -0.00645789, "std": 0.06894682}, - {"type": "rf_normal", "bins": 5, "mean": -0.03634136, "std": 0.53989509}, - {"type": "rf_normal", "bins": 10, "mean": 0.00529333, "std": 0.09129902}, + {"type": "rf_normal", "bins": 20, "mean": -0.00645789, "std": 1.0}, + {"type": "rf_normal", "bins": 20, "mean": -0.03634136, "std": 0.53989509}, + {"type": "rf_normal", "bins": 20, "mean": 0.00529333, "std": 0.09129902}, {"type": "rf_normal", "bins": 20, "mean": 0.04752025, "std": 0.82070579} ], "episodes": 21, @@ -24,9 +24,13 @@ "rewardsToKeep": 100, "critic": { "angv_bias": 0.2, - "total_gain": 5.0, - "max_reward": 2.0, - "posRewardBias": 2.0 + "total_gain": 6.0, + "max_reward": 3.0, + "posRewardBias": 2.0, + "modulation": { + "type": "linear", + "steps": 25 + } }, "verbose": 0, "stayStepLim": 0, @@ -35,12 +39,13 @@ "ResumeSimFromFile": "...pkl" }, "sim": { - "duration": 1000, + "duration": 2000, "dt": 0.6, "verbose": 0, "recordStep": 1.0, "recordStim": 0, "recordWeightStepSize": 1000, + "normalizeStepSize": 2333, "RLFakeUpRule": 0, "RLFakeDownRule": 0, "RLFakeStayRule": 0, @@ -114,15 +119,43 @@ "I": 0.0 }, "allpops": { - "ES": 40, + "ES": 80, "EA": 40, "IA": 10, "IAL": 10, - "EM": 40, - "IM": 10, - "IML": 10 + "EM": 60, + "IM": 15, + "IML": 15 + }, + "out_pop_sort_tag": "x", + "topographic": { + "stimModES": { + "xnormRange": [0.4, 0.6], + "ynormRange": [0.0, 0.1] + }, + "ES": { + "ynormRange": [0.1, 0.3] + }, + "EA": { + "ynormRange": [0.35, 0.45] + }, + "IA": { + "ynormRange": [0.50, 0.55] + }, + "IAL": { + "ynormRange": [0.55, 0.60] + }, + "EM": { + "ynormRange": [0.70, 0.90] + }, + "IM": { + "ynormRange": [0.90, 0.95] + }, + "IML": { + "ynormRange": [0.95, 1.0] + } }, "inputPop": "ES", "STDPconns": { @@ -161,18 +194,18 @@ "IA":{"GA2": 2.25,"p": 0.2}, "IAL":{"GA2": 5.5,"p": 0.1}}, "EM":{ - "EM":{"AM2": 0.01,"NM2": 0.001, "conv": 1}, - "IM":{"AM2": 1.95,"NM2":0.0195, "conv": 8}, - "IML":{"AM2": 0.98,"NM2":0.098, "conv": 8}, - "EA":{"AM2": 0.5,"NM2": 0.01,"conv": 6}}, + "EM":{"AM2": 0.01,"NM2": 0.001, "p": "0.2*exp(-dist_xnorm*3)"}, + "IM":{"AM2": 1.95,"NM2":0.0195, "p": "0.4*(1.0-exp(-dist_xnorm*3))"}, + "IML":{"AM2": 0.98,"NM2":0.098, "p": "0.4*(1.0-exp(-dist_xnorm*3))"}, + "EA":{"AM2": 0.5,"NM2": 0.01,"conv": 4}}, "IM":{ - "EM":{"GA": 9.0,"p": 0.4}, - "IM":{"GA": 4.5,"p": 0.1}, - "IML":{"GA": 4.5,"p": 0.2}}, + "EM":{"GA": 9.0,"p": "0.5*exp(-dist_xnorm*3)"}, + "IM":{"GA": 4.5,"p": "0.2*(1.0-exp(-dist_xnorm*3))"}, + "IML":{"GA": 4.5,"p": "0.05*(1.0-exp(-dist_xnorm*3))"}}, "IML":{ - "EM":{"GA2": 2.5,"p": 0.4}, - "IM":{"GA2": 2.25,"p": 0.2}, - "IML":{"GA2": 5.5,"p": 0.1}} + "EM":{"GA2": 2.5,"p": "0.5*exp(-dist_xnorm*3)"}, + "IM":{"GA2": 2.25,"p": "0.2*(1.0-exp(-dist_xnorm*3))"}, + "IML":{"GA2": 5.5,"p": "0.05*(1.0-exp(-dist_xnorm*3))"}} }, "weightVar": 0.5, "delays": { @@ -241,38 +274,5 @@ "softthresh": 0, "verbose": 0 } - }, - "noise": { - "EMLEFT": { - "AM": {"w":0, "rate":0}, - "NM": {"w":0, "rate":0}, - "GA": {"w":1.875, "rate": 100}, - "AM2": {"w": 3.938, "rate": 200}, - "NM2": {"w": 0, "rate": 0}, - "GA2": {"w": 1.875, "rate": 100}}, - "EMRIGHT": {"AM": {"w":0, "rate":0}, - "NM": {"w":0, "rate":0}, - "GA": {"w":1.875, "rate": 100}, - "AM2": {"w": 3.938, "rate": 200}, - "NM2": {"w": 0, "rate": 0}, - "GA2": {"w": 1.875, "rate": 100}}, - "EA": {"AM": {"w":0, "rate":0}, - "NM": {"w":0, "rate":0}, - "GA": {"w":1.875, "rate": 100}, - "AM2": {"w": 3.938, "rate": 200}, - "NM2": {"w": 0, "rate": 0}, - "GA2": {"w": 1.875, "rate": 100}}, - "IA": {"AM": {"w":1.875, "rate": 100}, - "NM": {"w":0, "rate":0}, - "GA": {"w":0, "rate":0}, - "AM2": {"w": 4.125, "rate": 200}, - "NM2": {"w": 0, "rate": 0}, - "GA2": {"w": 1.875, "rate": 100}}, - "IAL": {"AM": {"w":1.875, "rate": 100}, - "NM": {"w":0, "rate":0}, - "GA": {"w":0, "rate":0}, - "AM2": {"w": 3.0, "rate": 200}, - "NM2": {"w": 0, "rate": 0}, - "GA2": {"w": 1.875, "rate": 100}} } } \ No newline at end of file diff --git a/neurosim/critic.py b/neurosim/critic.py index f20b5b4..7003310 100644 --- a/neurosim/critic.py +++ b/neurosim/critic.py @@ -11,9 +11,10 @@ def _modulate_linear(mod_steps, reward, min_steps=1): k_pos = len([st for st in mod_steps if st > EPS]) k_neg = len([st for st in mod_steps if st < -EPS]) if reward >= 0: - return reward * (M - k_pos) / M + new_reward = reward * (M - k_pos + 1) / (M+1) else: - return reward * (M - k_neg) / M + new_reward = reward * (M - k_neg + 1) / (M+1) + return new_reward class Critic: diff --git a/neurosim/sim.py b/neurosim/sim.py index 35a6ba0..d6b043d 100644 --- a/neurosim/sim.py +++ b/neurosim/sim.py @@ -15,7 +15,7 @@ from cells import intf7 from game_interface import GameInterface from critic import Critic -from utils.conns import getconv, getSec, getInitDelay, getDelay, setdminID, setrecspikes +from utils.conns import getconn, getSec, getInitDelay, getDelay, setdminID, setrecspikes from utils.plots import plotRaster, plotWeights, saveActionsPerEpisode from utils.sync import syncdata_alltoall from utils.weights import saveSynWeights, readWeights, getWeightIndex, getInitSTDPWeight @@ -62,6 +62,8 @@ def outpath(fname): return os.path.join(dconf['sim']['outdir'], fname) sim.saveMotionFields = dconf['sim']['saveMotionFields'] sim.saveObjPos = 1 # save ball and paddle position to file self.recordWeightStepSize = dconf['sim']['recordWeightStepSize'] + self.normalizeStepSize = dconf['sim']['normalizeStepSize'] + self.normalizationMeans = None # time step per action (in ms) self.tstepPerAction = dconf['sim']['tstepPerAction'] @@ -81,6 +83,9 @@ def outpath(fname): return os.path.join(dconf['sim']['outdir'], fname) netParams = specs.NetParams() # spike threshold, 10 mV is NetCon default, lower it for all cells netParams.defaultThreshold = 0.0 + netParams.sizeX = 100 + netParams.sizeY = 100 + netParams.sizeZ = 1 self.netParams = netParams # object of class SimConfig to store simulation configuration @@ -131,13 +136,8 @@ def outpath(fname): return os.path.join(dconf['sim']['outdir'], fname) self.ECellModel = dconf['net']['ECellModel'] self.ICellModel = dconf['net']['ICellModel'] - - # Population parameters - for ty in self.allpops: - self.netParams.popParams[ty] = { - 'cellType': ty, - 'numCells': self.dnumc[ty], - 'cellModel': self.ECellModel if isExc(ty) else self.ICellModel} + self.topographic = dconf['net']['topographic'] if 'topographic' in dconf['net'] else None + self.dCellGids = None self.makeECellModel() self.makeICellModel() @@ -200,6 +200,12 @@ def makeECellModel(self): self.netParams.popParams[ty][k] = v self.PlastWeightIndex = intf7.dsyn['AM2'] + # Add topographic: + if self.topographic: + for ty in self.allpops: + if isExc(ty) and ty in self.topographic: + self.netParams.popParams[ty].update(self.topographic[ty]) + def makeICellModel(self): # create rules for inhibitory neuron models if self.ICellModel == 'IntFire4': @@ -228,6 +234,12 @@ def makeICellModel(self): for k, v in icell.dparam.items(): self.netParams.popParams[ty][k] = v + # Add topographic: + if self.topographic: + for ty in self.allpops: + if isInh(ty) and ty in self.topographic: + self.netParams.popParams[ty].update(self.topographic[ty]) + def readSTDPParams(self): lsy = ['AMPA', 'NMDA', 'AMPAI'] gains = [self.cfg.Gain[gainType] for gainType in ['EE', 'EE', 'EI']] @@ -262,6 +274,8 @@ def setupStimMod(self): 'rate': 'variable', 'noise': 0, 'start': 0} + if self.topographic and stimty in self.topographic: + self.netParams.popParams[stimty].update(self.topographic[stimty]) blist = [[i, i] for i in range(self.dnumc[poty])] self.netParams.connParams[stimty+'->'+poty] = { 'preConds': {'pop': stimty}, @@ -273,11 +287,12 @@ def setupStimMod(self): return lstimty def setupNoiseStim(self): + if 'noise' not in self.dconf: + return {} lnoisety = [] - dnoise = self.dconf['noise'] # setup noisy NetStim sources (send random spikes) if self.ECellModel == 'IntFire4' or self.ECellModel == 'INTF7': - for poty, dpoty in dnoise.items(): + for poty, dpoty in self.dconf['noise'].items(): for sy, dsy in dpoty.items(): damp = self.dconf['net']['noiseDamping']['E' if isExc(poty) else 'I'] Weight, Rate = dsy['w'] * damp, dsy['rate'] @@ -326,7 +341,6 @@ def setupSTDPWeights(self): self.netParams.connParams[k] = { 'preConds': {'pop': prety}, 'postConds': {'pop': poty}, - 'convergence': getconv(self.cmat, prety, poty, self.dnumc[prety]), 'weight': weight, 'delay': getDelay(self.dconf, prety, poty, sy, sec), 'synMech': synToMech[sy], @@ -335,6 +349,8 @@ def setupSTDPWeights(self): 'weightIndex': getWeightIndex( sy, self.ICellModel if isInh(poty) else self.ECellModel) } + self.netParams.connParams[k].update( + getconn(self.cmat, prety, poty, self.dnumc[prety])) # Setup STDP plasticity rules if ct in stdpConns and stdpConns[ct] and self.dSTDPparams[synToMech[sy]]['RLon']: print('Setting RL-STDP on {} ({})'.format(k, weight)) @@ -378,6 +394,39 @@ def recordWeights(self, sim, t): conn['hObj'].weight[self.PlastWeightIndex])] sim.allSTDPWeights.append(weightItem) + def weightsMean(self, sim): + pop_sizes = {} + for cell in sim.net.cells: + for conn in cell.conns: + if 'hSTDP' in conn: + pops = [pop for pop,dpop in sim.net.pops.items() if cell.gid in dpop.cellGids] + pop = pops[0] + if pop not in pop_sizes: + pop_sizes[pop] = [] + weight = conn['hObj'].weight[self.PlastWeightIndex] + if weight > 0: + # once the weight reaches 0, we stop counting that weight towards the mean + # then it the mean of the _firing_ neurons will remain the same. + pop_sizes[pop].append(weight) + norm_means = dict([(pop, np.mean(weights)) for pop, weights in pop_sizes.items()]) + return norm_means + + + def normalizeWeights(self, sim): + pop_means = self.weightsMean(sim) + norm_means = self.normalizationMeans + print('\n!Normalizing from means: {} to {}!'.format(pop_means, norm_means)) + for cell in sim.net.cells: + for conn in cell.conns: + if 'hSTDP' in conn: + pops = [pop for pop,dpop in sim.net.pops.items() if cell.gid in dpop.cellGids] + pop = pops[0] + new_weight = conn['hObj'].weight[self.PlastWeightIndex] - pop_means[pop] + norm_means[pop] + if new_weight < 0: + new_weight = 0.0 + conn['hObj'].weight[self.PlastWeightIndex] = new_weight + + def getSpikesWithInterval(self, trange=None, neuronal_pop=None): if len(neuronal_pop) < 1: return 0.0 @@ -458,6 +507,30 @@ def getAllSTDPObjects(self, sim): dSTDPmech[pop].append(STDPmech) return dSTDPmech + def cellGidsMap(self, sim): + if not self.dCellGids: + # Check if the self.dCellGids is cached. if not, compute it only once + pop_to_moves = self.dconf['pop_to_moves'] + cgids_map = {} + for p, pop_moves in pop_to_moves.items(): + if type(pop_moves) == str: + pop_moves = [pop_moves] + cgids = sim.net.pops[p].cellGids + if 'out_pop_sort_tag' in self.dconf['net'] and self.dconf['net']['out_pop_sort_tag']: + # sort the cgids by `sort_tag` dimension (easiest: 'x') + sort_tag = self.dconf['net']['out_pop_sort_tag'] + cgids, _ = list(zip(*sorted( + [(gid, sim.net.cells[gid].tags[sort_tag]) for gid in cgids], + key=lambda x:x[1]))) + # Create map of cell gids to env movement. + # For example: {150: 'LEFT', 151: 'RIGHT', ...} + cells_per_move = math.floor(len(cgids) / len(pop_moves)) + for idx, cgid in enumerate(cgids): + cgids_map[cgid] = pop_moves[math.floor(idx / cells_per_move)] + self.dCellGids = cgids_map + + return self.dCellGids + def getActions(self, sim, t, moves, pop_to_moves): # Get move frequencies move_freq = {} @@ -466,16 +539,7 @@ def getActions(self, sim, t, moves, pop_to_moves): for ts in range(int(self.dconf['actionsPerPlay'])): ts_beg = t-self.tstepPerAction*(self.dconf['actionsPerPlay']-ts-1) ts_end = t-self.tstepPerAction*(self.dconf['actionsPerPlay']-ts) - cgids_map = {} - for p, pop_moves in pop_to_moves.items(): - if type(pop_moves) == str: - pop_moves = [pop_moves] - cells_per_move = math.floor( - len(sim.net.pops[p].cellGids) / len(pop_moves)) - for idx, cgid in enumerate(sim.net.pops[p].cellGids): - cgids_map[cgid] = pop_moves[math.floor(idx / cells_per_move)] - - freq.append(self.getSpikesWithInterval([ts_end, ts_beg], cgids_map)) + freq.append(self.getSpikesWithInterval([ts_end, ts_beg], self.cellGidsMap(sim))) for move in moves: freq_move = [q[move] for q in freq] sim.pc.allreduce(vec.from_python(freq_move), 1) # sum @@ -619,11 +683,16 @@ def trainAgent(self, t): print('Weights Recording Time:', t, 'NBsteps:', self.NBsteps, 'recordWeightStepSize:', self.recordWeightStepSize) self.recordWeights(sim, t) + if self.normalizeStepSize and not self.normalizationMeans: + # first step count averages per population + self.normalizationMeans = self.weightsMean(sim) + if self.normalizeStepSize and self.NBsteps % self.normalizeStepSize == 0: + self.normalizeWeights(sim) t5 = datetime.now() - t5 if random.random() < 0.001: - print(t, [round(tk.microseconds / 1000, 0) - for tk in [t1, t2, t3, t4, t5]]) + print('\n', t, [round(tk.microseconds / 1000, 0) + for tk in [t1, t2, t3, t4, t5]]) if 'sleeptrial' in dconf['sim'] and dconf['sim']['sleeptrial']: time.sleep(dconf['sim']['sleeptrial']) @@ -690,7 +759,7 @@ def run(self): self.InitializeNoiseRates(sim) # Plot 2d net - # sim.analysis.plot2Dnet(saveFig='data/net.png', showFig=False) + sim.analysis.plot2Dnet(saveFig=self.outpath('net.png'), showFig=False) sim.analysis.plotConn( saveFig=self.outpath('connsCells.png'), showFig=False, groupBy='cell', feature='weight') diff --git a/neurosim/utils/conns.py b/neurosim/utils/conns.py index e2b1054..9bd2ab7 100644 --- a/neurosim/utils/conns.py +++ b/neurosim/utils/conns.py @@ -8,14 +8,16 @@ def prob2conv(prob, npre): return int(0.5 + prob * npre) -def getconv(cmat, prety, poty, npre): +def getconn(cmat, prety, poty, npre): # get convergence value from cmat dictionary # (uses convergence if specified directly, otherwise uses p to calculate) if 'conv' in cmat[prety][poty]: - return cmat[prety][poty]['conv'] + conv = cmat[prety][poty]['conv'] elif 'p' in cmat[prety][poty]: - return prob2conv(cmat[prety][poty]['p'], npre) - return 0 + if type(cmat[prety][poty]['p']) == str: + return {'probability': cmat[prety][poty]['p']} + conv = prob2conv(cmat[prety][poty]['p'], npre) + return {"convergence": conv} def getInitDelay(dconf, sec):