Subsetting and executing a small subset of the Allen V1 model (layer23) using cpp intermediate objects where possible.
27 Jan
2023
27 Jan
'23
7:17 a.m.
Hi there,
I would like to subset and execute the subset of the Allen V1 model.
I am also hoping to store intermediate simulation configuration files not as Python code,
but if possible as cpp code. I am wondering if there is a way of storing the matrices in a
cpp native way?
Forinstance I have noticed the class
https://github.com/nest/nest-simulator/blob/master/nestkernel/connection_ma… I
think that this would class may have
I acquired the config.json file from DropBox and all of the assoicated hdf5, and csv files
from here:
https://www.dropbox.com/sh/w5u31m3hq6u2x5m/AAB0B6OLRyLJwxCKjuWrG26pa/simula…
I have found the nodes in the model pertaining to only layer 23 by wrangling data frames
in Julia (see the code below). Calling the method `extract_subset` is slower than I would
like. I suspect it would be faster to use the SL or /nest_kernel more directly, but I
don't know what interventions I would need to do to get the PyNest instance to dump
its intermediate objects in a cpp native way. I wonder if any one has any hints on how to
do this?
subset0 = get_layer_sub_net("i23Htr3a")
subset1 = get_layer_sub_net("i23Pvalb")
subset2 = get_layer_sub_net("i23Sst")
subset3 = get_layer_sub_net("e23Cux2")
For instance I would like to be able to use the nestkernel CLI to save
using Pkg
using Conda
using PyCall
ENV["PYTHON"] = $HOME/.julia/conda/3/lib/python3.10/
Pkg.build("PyCall")
py"""
def init_sim_to_call_from_julia():#slice_start,slice_end,nest=None,sim=None)
import os, sys
import math
import numpy as np
from bmtk.simulator import pointnet
from bmtk.simulator.pointnet.pyfunction_cache import synaptic_weight
from bmtk.simulator.pointnet.io_tools import io
import nest
import h5py
try:
nest.Install('glifmodule')
except Exception as e:
pass
@synaptic_weight
def DirectionRule_others(edges, src_nodes, trg_nodes):
src_tuning = src_nodes['tuning_angle'].values
tar_tuning = trg_nodes['tuning_angle'].values
sigma = edges['weight_sigma'].values
nsyn = edges['nsyns'].values
syn_weight = edges['syn_weight'].values
delta_tuning_180 = np.abs(np.abs(np.mod(np.abs(tar_tuning - src_tuning), 360.0) -
180.0) - 180.0)
w_multiplier_180 = np.exp(-(delta_tuning_180 / sigma) ** 2)
return syn_weight * w_multiplier_180 * nsyn
@synaptic_weight
def DirectionRule_EE(edges, src_nodes, trg_nodes):
src_tuning = src_nodes['tuning_angle'].values
tar_tuning = trg_nodes['tuning_angle'].values
x_tar = trg_nodes['x'].values
x_src = src_nodes['x'].values
z_tar = trg_nodes['z'].values
z_src = src_nodes['z'].values
sigma = edges['weight_sigma'].values
nsyn = edges['nsyns'].values
syn_weight = edges['syn_weight'].values
delta_tuning_180 = np.abs(np.abs(np.mod(np.abs(tar_tuning - src_tuning), 360.0) -
180.0) - 180.0)
w_multiplier_180 = np.exp(-(delta_tuning_180 / sigma) ** 2)
delta_x = (x_tar - x_src) * 0.07
delta_z = (z_tar - z_src) * 0.04
theta_pref = tar_tuning * (np.pi / 180.)
xz = delta_x * np.cos(theta_pref) + delta_z * np.sin(theta_pref)
sigma_phase = 1.0
phase_scale_ratio = np.exp(- (xz ** 2 / (2 * sigma_phase ** 2)))
# To account for the 0.07 vs 0.04 dimensions. This ensures the horizontal neurons
are scaled by 5.5/4 (from the
# midpoint of 4 & 7). Also, ensures the vertical is scaled by 5.5/7. This was
a basic linear estimate to get the
# numbers (y = ax + b).
theta_tar_scale = abs(abs(abs(180.0 - np.mod(np.abs(tar_tuning), 360.0)) - 90.0) -
90.0)
phase_scale_ratio = phase_scale_ratio * (5.5 / 4.0 - 11.0 / 1680.0 *
theta_tar_scale)
return syn_weight * w_multiplier_180 * phase_scale_ratio * nsyn
def main(config_file):
configure = pointnet.Config.from_json(config_file)
configure.build_env()
graph = pointnet.PointNetwork.from_config(configure)
sim = pointnet.PointSimulator.from_config(configure, graph)
node_info = list(graph.get_node_populations())
dfv1 = node_info[0].nodes_df()
pop_names = { k:v for k,v in zip(dfv1.index,dfv1["pop_name"].values) }
dfv1e = dfv1[dfv1["ei"] == "e"].index
dfv1i = dfv1[dfv1["ei"] == "i"].index
dfv1e_pop = dfv1[dfv1["ei"] ==
"e"]["pop_name"].values
dfv1i_pop = dfv1[dfv1["ei"] ==
"i"]["pop_name"].values
lgn = node_info[1].nodes_df()
lgne = lgn[lgn["ei"] == "e"].index
return nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names #,lgni
nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names =
main('config.json') # ,lgni
return (nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names)
"""
(nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names) =
py"init_sim_to_call_from_julia"()
function get_layer_sub_net(target)
subset = Int64[]
for (k,v) in pairs(pop_names)
if v==target
println(k,v)
append!(subset,k)
end
end
subset
end
subset0 = get_layer_sub_net("i23Htr3a")
subset1 = get_layer_sub_net("i23Pvalb")
subset2 = get_layer_sub_net("i23Sst")
subset3 = get_layer_sub_net("e23Cux2")
layer_23_ids = []
layer_23_ids = vcat(layer_23_ids,subset0)
layer_23_ids = vcat(layer_23_ids,subset1)
layer_23_ids = vcat(layer_23_ids,subset2)
layer_23_ids = vcat(layer_23_ids,subset3)
py"""def get_conns_from_py_nest(nest,sim,slice_start,nodes=None):
conns = nest.GetConnections(nest.NodeCollection([slice_start]))
nodes = nest.NodeCollection([slice_start])
node_ = nodes.get()
dict_ = conns.get()
return
(dict_["source"],dict_["target"],dict_["weight"])
"""
function extract_subset(nest,sim,subset)
N = 230924+1#+17400#230901
exc_matrix = spzeros(N,N)
@inbounds for srcv in ProgressBars.ProgressBar(subset)
result = py"get_conns_from_py_nest"(nest,sim,srcv)
@inbounds for (x,y,z) in
zip(Vector{Int64}(result[1]),Vector{Int64}(result[2]),Vector{Float64}(result[3]))
exc_matrix[y,x] = z
end
end
exc_matrix
end
exc_matrix = extract_subset(nest,sim,layer_23_ids)
Thanks for your help :)
Russell.
28 Jan
28 Jan
4:50 a.m.
New subject: Update: Ablating and executing a small subset of the Allen V1 model (layer23) using cpp intermediate objects where possible.
Hi Nest forum,
In retrospect, I probably should have just asked for tips about how I could change the
simulation configuration, maybe even before the job is loaded into memory.
The machine I would like to load the model onto only has 64GB of memory (its not a
cluster), I am not interested in any one layer more than another (there was no reason for
me to choose to model layer23 over other layers for example). I would be just as happy to
model all V1 layers, but to diminish the number of cell counts in every cell population by
1,000 fold or more. I think that it's possible that someone has already done a related
smaller scale simulation when they where debugging and developing the full scale
simulation for the publication.
It would be great if anyone knows a code approach or repository that illustrates how to
scale down a simulation specified in a config.json, but without doing 24 hours of
dataframe wrangling.
Thanks again for any help or tips.
Regards,
Russell
________________________________________
From: Russell Jarvis [R.Jarvis(a)westernsydney.edu.au]
Sent: 27 January 2023 17:17
To: users(a)nest-simulator.org
Subject: [External] [NEST Users] Subsetting and executing a small subset of the Allen V1
model (layer23) using cpp intermediate objects where possible.
Hi there,
I would like to subset and execute the subset of the Allen V1 model.
I am also hoping to store intermediate simulation configuration files not as Python code,
but if possible as cpp code. I am wondering if there is a way of storing the matrices in a
cpp native way?
Forinstance I have noticed the class
https://github.com/nest/nest-simulator/blob/master/nestkernel/connection_ma…
I think that this would class may have
I acquired the config.json file from DropBox and all of the assoicated hdf5, and csv files
from here:
https://www.dropbox.com/sh/w5u31m3hq6u2x5m/AAB0B6OLRyLJwxCKjuWrG26pa/simula…
I have found the nodes in the model pertaining to only layer 23 by wrangling data frames
in Julia (see the code below). Calling the method `extract_subset` is slower than I would
like. I suspect it would be faster to use the SL or /nest_kernel more directly, but I
don't know what interventions I would need to do to get the PyNest instance to dump
its intermediate objects in a cpp native way. I wonder if any one has any hints on how to
do this?
subset0 = get_layer_sub_net("i23Htr3a")
subset1 = get_layer_sub_net("i23Pvalb")
subset2 = get_layer_sub_net("i23Sst")
subset3 = get_layer_sub_net("e23Cux2")
For instance I would like to be able to use the nestkernel CLI to save
using Pkg
using Conda
using PyCall
ENV["PYTHON"] = $HOME/.julia/conda/3/lib/python3.10/
Pkg.build("PyCall")
py"""
def init_sim_to_call_from_julia():#slice_start,slice_end,nest=None,sim=None)
import os, sys
import math
import numpy as np
from bmtk.simulator import pointnet
from bmtk.simulator.pointnet.pyfunction_cache import synaptic_weight
from bmtk.simulator.pointnet.io_tools import io
import nest
import h5py
try:
nest.Install('glifmodule')
except Exception as e:
pass
@synaptic_weight
def DirectionRule_others(edges, src_nodes, trg_nodes):
src_tuning = src_nodes['tuning_angle'].values
tar_tuning = trg_nodes['tuning_angle'].values
sigma = edges['weight_sigma'].values
nsyn = edges['nsyns'].values
syn_weight = edges['syn_weight'].values
delta_tuning_180 = np.abs(np.abs(np.mod(np.abs(tar_tuning - src_tuning), 360.0) - 180.0) -
180.0)
w_multiplier_180 = np.exp(-(delta_tuning_180 / sigma) ** 2)
return syn_weight * w_multiplier_180 * nsyn
@synaptic_weight
def DirectionRule_EE(edges, src_nodes, trg_nodes):
src_tuning = src_nodes['tuning_angle'].values
tar_tuning = trg_nodes['tuning_angle'].values
x_tar = trg_nodes['x'].values
x_src = src_nodes['x'].values
z_tar = trg_nodes['z'].values
z_src = src_nodes['z'].values
sigma = edges['weight_sigma'].values
nsyn = edges['nsyns'].values
syn_weight = edges['syn_weight'].values
delta_tuning_180 = np.abs(np.abs(np.mod(np.abs(tar_tuning - src_tuning), 360.0) - 180.0) -
180.0)
w_multiplier_180 = np.exp(-(delta_tuning_180 / sigma) ** 2)
delta_x = (x_tar - x_src) * 0.07
delta_z = (z_tar - z_src) * 0.04
theta_pref = tar_tuning * (np.pi / 180.)
xz = delta_x * np.cos(theta_pref) + delta_z * np.sin(theta_pref)
sigma_phase = 1.0
phase_scale_ratio = np.exp(- (xz ** 2 / (2 * sigma_phase ** 2)))
# To account for the 0.07 vs 0.04 dimensions. This ensures the horizontal neurons are
scaled by 5.5/4 (from the
# midpoint of 4 & 7). Also, ensures the vertical is scaled by 5.5/7. This was a basic
linear estimate to get the
# numbers (y = ax + b).
theta_tar_scale = abs(abs(abs(180.0 - np.mod(np.abs(tar_tuning), 360.0)) - 90.0) - 90.0)
phase_scale_ratio = phase_scale_ratio * (5.5 / 4.0 - 11.0 / 1680.0 * theta_tar_scale)
return syn_weight * w_multiplier_180 * phase_scale_ratio * nsyn
def main(config_file):
configure = pointnet.Config.from_json(config_file)
configure.build_env()
graph = pointnet.PointNetwork.from_config(configure)
sim = pointnet.PointSimulator.from_config(configure, graph)
node_info = list(graph.get_node_populations())
dfv1 = node_info[0].nodes_df()
pop_names = { k:v for k,v in zip(dfv1.index,dfv1["pop_name"].values) }
dfv1e = dfv1[dfv1["ei"] == "e"].index
dfv1i = dfv1[dfv1["ei"] == "i"].index
dfv1e_pop = dfv1[dfv1["ei"] == "e"]["pop_name"].values
dfv1i_pop = dfv1[dfv1["ei"] == "i"]["pop_name"].values
lgn = node_info[1].nodes_df()
lgne = lgn[lgn["ei"] == "e"].index
return nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names #,lgni
nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names =
main('config.json') # ,lgni
return (nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names)
"""
(nest,sim,node_info,dfv1e,dfv1i,lgne,dfv1e_pop,dfv1i_pop,pop_names) =
py"init_sim_to_call_from_julia"()
function get_layer_sub_net(target)
subset = Int64[]
for (k,v) in pairs(pop_names)
if v==target
println(k,v)
append!(subset,k)
end
end
subset
end
subset0 = get_layer_sub_net("i23Htr3a")
subset1 = get_layer_sub_net("i23Pvalb")
subset2 = get_layer_sub_net("i23Sst")
subset3 = get_layer_sub_net("e23Cux2")
layer_23_ids = []
layer_23_ids = vcat(layer_23_ids,subset0)
layer_23_ids = vcat(layer_23_ids,subset1)
layer_23_ids = vcat(layer_23_ids,subset2)
layer_23_ids = vcat(layer_23_ids,subset3)
py"""def get_conns_from_py_nest(nest,sim,slice_start,nodes=None):
conns = nest.GetConnections(nest.NodeCollection([slice_start]))
nodes = nest.NodeCollection([slice_start])
node_ = nodes.get()
dict_ = conns.get()
return (dict_["source"],dict_["target"],dict_["weight"])
"""
function extract_subset(nest,sim,subset)
N = 230924+1#+17400#230901
exc_matrix = spzeros(N,N)
@inbounds for srcv in ProgressBars.ProgressBar(subset)
result = py"get_conns_from_py_nest"(nest,sim,srcv)
@inbounds for (x,y,z) in
zip(Vector{Int64}(result[1]),Vector{Int64}(result[2]),Vector{Float64}(result[3]))
exc_matrix[y,x] = z
end
end
exc_matrix
end
exc_matrix = extract_subset(nest,sim,layer_23_ids)
Thanks for your help :)
Russell.
_______________________________________________
NEST Users mailing list -- users(a)nest-simulator.org
To unsubscribe send an email to users-leave(a)nest-simulator.org
451
days inactive
452
days old
1 comments
1 participants
participants (1)
-
Russell Jarvis