From 8263f1d70995d5815c30691d2b14ef4a840ff356 Mon Sep 17 00:00:00 2001 From: Weina Ji Date: Mon, 17 Jul 2023 11:28:33 +0000 Subject: [PATCH] [BBPBGLIB-1036] Implement SONATA synapse reader for gap junction ## Context Recently we have dropped the synapse-tool dependency for reading synapse parameters. SynReaderNRN is used to read the old format nrn edges files, and SonataReader is to read sonata edges. This PR is to implement a sonata reader to read gap junction synapse parameters. ## Scope In `gap_junection.py`, a GapJunctionSonata reader derived from the SonataReader. The Sonata parameters for gap junction are `efferent_junction_id` and `afferent_junction_id`. ## Testing Add a unit test to test GapJunctionSynapseReader for both NRN and SONATA edges files. ## Review * [x] PR description is complete * [x] Coding style (imports, function length, New functions, classes or files) are good * [x] Unit/Scientific test added * [ ] Updated Readme, in-code, developer documentation --- neurodamus/gap_junction.py | 48 ++++++------------ neurodamus/io/synapse_reader.py | 36 ++++++------- .../mini_thalamus_sonata/gapjunction/README | 2 + .../mini_thalamus_sonata/gapjunction/edges.h5 | Bin 0 -> 169312 bytes tests/test_synapse_reader.py | 48 ++++++++++++++++++ 5 files changed, 84 insertions(+), 50 deletions(-) create mode 100644 tests/simulations/mini_thalamus_sonata/gapjunction/README create mode 100644 tests/simulations/mini_thalamus_sonata/gapjunction/edges.h5 create mode 100644 tests/test_synapse_reader.py diff --git a/neurodamus/gap_junction.py b/neurodamus/gap_junction.py index 9d5dd67e..314568b1 100644 --- a/neurodamus/gap_junction.py +++ b/neurodamus/gap_junction.py @@ -8,7 +8,8 @@ from .connection_manager import ConnectionManagerBase from .core.configuration import ConfigurationError -from .io.synapse_reader import SynapseReader, SynReaderNRN, SynapseParameters, FormatNotSupported +from .io.synapse_reader import SynapseReader, SynReaderNRN, SonataReader, SynapseParameters, \ + FormatNotSupported from .utils import compat from .utils.logging import log_verbose @@ -17,40 +18,9 @@ class GapJunctionConnParameters(SynapseParameters): # Attribute names of synapse parameters, consistent with the normal synapses _synapse_fields = ("sgid", "isec", "offset", "weight", "D", "F", "ipt", "location") - # Actual fields to read from conectivity files - _gj_v1_fields = [ - "connected_neurons_pre", - "morpho_section_id_post", - "morpho_offset_segment_post", - "conductance", - "junction_id_pre", - "junction_id_post", - "morpho_segment_id_post", - - ] - _gj_v2_fields = [ - "connected_neurons_pre", - "morpho_section_id_post", - "morpho_section_fraction_post", # v2 field - "conductance", - "junction_id_pre", - "junction_id_post" - ] - # SONATA fields, see conversion map in spykfunc/schema.py and - # synapse-tool Naming::get_property_mapping - _gj_sonata_fields = [ - "connected_neurons_pre", - "morpho_section_id_post", - "morpho_section_fraction_post", - "conductance", - "efferent_junction_id", - "afferent_junction_id" - ] - @classmethod def create_array(cls, length): npa = np.recarray(length, cls.dtype) - npa.ipt = -1 npa.location = 0.5 return npa @@ -65,6 +35,10 @@ class GapJunctionSynapseReader(SynapseReader): def create(cls, syn_src, conn_type, population=None, *args, **kw): """Instantiates a synapse reader, giving preference to GapJunctionSynToolReader """ + if fn := cls._get_sonata_circuit(syn_src): + log_verbose("[SynReader] Using SonataReader.") + return GapJunctionSonataReader(fn, conn_type, population, **kw) + # Syn2 support dropped (July 2023) if not ospath.isdir(syn_src) and not syn_src.endswith(".h5"): raise FormatNotSupported("File: {syn_src}") @@ -72,6 +46,16 @@ def create(cls, syn_src, conn_type, population=None, *args, **kw): return SynReaderNRN(syn_src, conn_type, None, *args, **kw) +class GapJunctionSonataReader(SonataReader): + Parameters = GapJunctionConnParameters + parameter_mapping = { + "weight": "conductance", + "D": "efferent_junction_id", + "F": "afferent_junction_id", + } + # "isec", "ipt", "offset" are custom parameters as in base class + + class GapJunctionManager(ConnectionManagerBase): """ The GapJunctionManager is similar to the SynapseRuleManager. It will diff --git a/neurodamus/io/synapse_reader.py b/neurodamus/io/synapse_reader.py index 17a928c9..968b1eb3 100644 --- a/neurodamus/io/synapse_reader.py +++ b/neurodamus/io/synapse_reader.py @@ -12,23 +12,6 @@ from ..utils.logging import log_verbose -def _get_sonata_circuit(path): - """Returns a SONATA edge file in path if present - """ - if os.path.isdir(path): - filename = os.path.join(path, "edges.sonata") - if os.path.exists(filename): - return filename - elif path.endswith(".sonata"): - return path - elif path.endswith(".h5"): - import h5py - f = h5py.File(path, 'r') - if "edges" in f: - return path - return None - - def _constrained_hill(K_half, y): K_half_fourth = K_half**4 y_fourth = y**4 @@ -185,6 +168,23 @@ def has_property(self, field_name): """Checks whether source data has the given additional field. """ + @staticmethod + def _get_sonata_circuit(path): + """Returns a SONATA edge file in path if present + """ + if os.path.isdir(path): + filename = os.path.join(path, "edges.sonata") + if os.path.exists(filename): + return filename + elif path.endswith(".sonata"): + return path + elif path.endswith(".h5"): + import h5py + f = h5py.File(path, 'r') + if "edges" in f: + return path + return None + @classmethod def create(cls, syn_src, conn_type=SYNAPSES, population=None, *args, **kw): """Instantiates a synapse reader, giving preference to SonataReader @@ -194,7 +194,7 @@ def create(cls, syn_src, conn_type=SYNAPSES, population=None, *args, **kw): return cls(syn_src, conn_type, population, *args, **kw) kw["verbose"] = (MPI.rank == 0) - if fn := _get_sonata_circuit(syn_src): + if fn := cls._get_sonata_circuit(syn_src): log_verbose("[SynReader] Using SonataReader.") return SonataReader(fn, conn_type, population, **kw) diff --git a/tests/simulations/mini_thalamus_sonata/gapjunction/README b/tests/simulations/mini_thalamus_sonata/gapjunction/README new file mode 100644 index 00000000..6853bf31 --- /dev/null +++ b/tests/simulations/mini_thalamus_sonata/gapjunction/README @@ -0,0 +1,2 @@ +Extract connections with target gid = 0 from the gap junction edges file in +/gpfs/bbp.cscs.ch/project/proj12/SIT/thalamus_sonata/networks/edges/thalamus_neurons__thalamus_neurons__electrical_synapse/edges.h5 diff --git a/tests/simulations/mini_thalamus_sonata/gapjunction/edges.h5 b/tests/simulations/mini_thalamus_sonata/gapjunction/edges.h5 new file mode 100644 index 0000000000000000000000000000000000000000..b15e58e5e62c8097fbd4ad287bb351141ddb356d GIT binary patch literal 169312 zcmeI*4|G)3oxt%s`5!_g5vjrUlptVYku~i?n?*D8()%cBvA{vjHt2GOG&4zX>150w zH&Ilyh_MSrx7C97R5-i!(4)tLt+@W_QaA@~a9b9RyXvB=o;|h(*Y(dHw~9yCeQ)l) z2{X)0oDe26lkY_4H}}5#-n;kn?!3A0`X;%f$^J(9lnbT^_fb+J3PeTZrFr^PKi1uy z7s)bdvwWR?;n8tG$33dJ@I-!usM7hBdi$~ER{{!xFOiu)3v^}WrM%FKhzTL4vzG9 zbDlNZyK?yo+wH@q-&N~=JN*8L>BnY0VzfY_^T5htRW7XG-=m!e-Yr+Zf&x6>xs?-&1I@hefs^a&I1Jr$Fm(5s`8OR zjq}p_GnF>|V(`yZy2Koy%T#*Mq*apVEEO>4fdhK{`l3FrbzeKX4N~WU;^FZ?o$j}$ zTgCQcOIN3(CuX1>Q%%#4`+Ad?tJCW>cPB5`q05!rsLGl8)#>!gn^k^=E>}M-rb02v zYMpK}E5s~MHR%yW>Ka2=%%iVxJsy>DuHK-zj?H(jd)+1O?$ub#ij*jHQf_`;zS>-g z-cEjAUY{>gLJT=WG&6h?Q$`0e|I15KUEX8 zs(P*Rn=93ZRi331XZV~3kD)eI$Ft_4der!MER`#+i1_OMcB;SIjM!1Vc{k|dr!7}a z7h>-y{oIrzB7gt_2q1s}0tg_000IagfWX-o@Q6Y+|DrGaN1q?6ij=Rcp5@$w;XbOx z`6c4~n{`!_%o8>LP?d4#MbxU1^nhM7GlH6TEEt~1YM!0EUe$W>-Z{z3)$7Fv<|Z$< zQZL?gQSx$`zBG;fs`deekY2x((XWp@;t3B8= zl6LQx@%QwHgMN3Z5Rr6yS7&cnWt^H01%iEM`|Cr2^{vjGHfCPaG7D0SsUJT*kM2xs ztC?-J=OpLH!GBYAikdeT;?5Lnnx5ik;&JN#5lP}v|H+7lnRa$^a>ZlJpQ9b+sx#-< z`EKO>@>%=Gyb*^S>>nw6bS%ngk1@Yrnm!(J$=Gq{c%QXXKR=JzcgK<}@HON54h?;0 zO>*%~dhxE>c0&2b%aWgIl(@fkn&n;shR$idJ3F`tL~em&xnv7dA1 zJUpo1#dQmkznx#xiw7=EUhZAJct>6Ga^^n8Jzq{hZTaPr)v~nHY)9@7F8D=3gJ1e@$cHF_JA#%MW+DUw1k^-<{Km&Ts#oTVFNJb82~&-`UpI9qbXUouNQ)Yfnq4)m__} z4>IkHso!vZO#M#F*Y!IsU)S%H{HXe!mapqKlHc3a*5b@;Y53FfHT-G$8ooQcJ{0tC z2!`5wI^1@uE-1v1o5gTsHe{NQtr{2vOJ)fdJH%mmd?f{Q}CHcXx*X1e}B`>#9 zm#hA2@^T%z+|l{T%Wcx--u;W@<#yWM*Kgj(qk@tOte8X3vZf=o|`s>}p;{*Nc333aC)*yfY0tg_000IagfB*sr zAb{@g8P_D8-i+CSZQ zk^RA&1@=GY&9LjwFSTp#DYJiJ)!6l0X4`KaonrqiSY==G$GP@_8AWzu={ff8-z>1( zo+-C;(S|l5fB*srAb!jBH!`fEswSBmb;ETEpIQ|E_Z%Ng#j#0tg_0KrRU!5{JdB z;uYuTS#exEFMcmhh}*;#ald$4+#y~N$DHLn5I_I{1Q0+Vw*(T+7wnoAznwc~#Gfxv zchkpThg<{@KmY**5Xf$UlDUt0KdLGAc9mPc*?+vp`=yf?`wmk_dVrX+<)Br-pfAgy>A_~ zrh9i<`?lZWTXg6zN3M87(w`Rp}&5mIF z`TUaj#mPkg0R#|00D&A3xOb=fIPA-mm;zgE_uFORVW!FmAbgw6_NakW0jFO|8*6c7Eh0)&xt$>u~5gQI$ojUavfLexJJjTbX=<=0RaRM zKmY**5I_I{1Q0*~0R#|0009ILKmY**5I_I{1Q0*~fr$`Ebiduv;`p`cdo=$2c6srO zlZyZX2q1s}0y!XXLfj^{i2KFU;tuhGI3^B@L(acf#VgLwv*NgTUi>}>T1gWSKmY** z5J2Emfkg8KiT=LrogQA0R#|00D;^Rcy%D?y}x;<_rABM`<68d-=z=snWz^)(zm;SIdL*jMV`=l!Ur*xOZZ`8qH9iFezZ7yAxe8umW?@f6>x{nGcu z=A+)+QlT{nAb$UP*#~S1>R&15K=Uyit zS$w^`yR1?EJU=X#e$X%f^B*pg|9&zo10OcXzpLLWi)UXaD?fL=eBi5Jmo1CKvSxq3 zd}a0JGLzawT@XM30R#|0009ILKmY**5Xd%xROj=n7BZioZO54UB7gt_2q18#0`nHt z*q6LF+kWh_GP`^4OnYtla*_B{lXhPt3Oe=Ylf3wtJ?1 z=%#sg`%jDQKQv6U-`QVmZ+WTQKGS-zJ^~0JfB*srAbhPI{ZCtm>y5ru4$TE{hL;wK<5I_I{1Q0*~0R#|00D-Xyc$~j{?{Qa+ zJh`DNvLTO0Wt^)wYhu2rh!m|T5yJV;d03I>LMP>}k(ZYjS+PWKC*RY2gZe)JcSg#M zJjJF*mQ3WUv+pr9IXs+$_}ThL-^Og}9~0R(**5|RATXW+pJv!&Jli&&+tR}$VeV$g z*BdX;u~5e%9gB5TpMgS@>hu&H%XFNoW4Vs!=r~Qs3LPtT^f;p|w=QawCPbw(9&;ac zH>#Sgyjex}J`Hy?7khXkg`KB?TUEJglR=zR0!>RASGo;u)~&ABF|7NKPKxf%-mcc5 zzo*mR)6&%*>@oQ&<>sprCf~eIs!c{rOw@nh)%&f`{c^Pa``TvuZw&3EN+#;Rb-Mo; zj8JDF=wBajJ8D-;s6E)N7eo?nZZOd97BJ;U>%UhnUtxQkPXTjYAJqNUtQXsltM!KK zx2yhY%B|DsO?Rn$v;WBwJ)$V`Ii{Z6PpaBT(>-eXK$y=-qZ6ARo`-PO`fS&kwOiD- zc101NtJijq%fsJ%#OVj&sO9NO#P+MB-#0CfUxR9$UHil2<>IfM`?=S>~>r-`PXCYpMhXzFRAsi%pi zo+g@lnrP~kY4tJnGJd9hCf`J}+(a|JHqnf;O*G?e6V15WL^J+2(Tu}QG~;m-V~?k) zClLVz5I`V30_h)D4rcVYGST&R?DfZJ*V`pOQX`S-^c+U^l3@aguD5r^T)&R@^>+H- t&mFg`p1v#U2xRJY;)gN*QI&QQ)jDJ4kIj`tY>CO4PJh2D8OY>*_5ZvtCb9qk literal 0 HcmV?d00001 diff --git a/tests/test_synapse_reader.py b/tests/test_synapse_reader.py new file mode 100644 index 00000000..b44d3feb --- /dev/null +++ b/tests/test_synapse_reader.py @@ -0,0 +1,48 @@ +import numpy as np +import numpy.testing as npt +import os +import pytest +from neurodamus.gap_junction import GapJunctionSynapseReader +from pathlib import Path + + +pytestmark = [ + pytest.mark.forked, + pytest.mark.skipif( + not os.environ.get("NEURODAMUS_NEOCORTEX_ROOT"), + reason="Test requires loading a neocortex model to run" + ) +] + + +def test_gapjunction_synreaderNRN(): + nrn_file = "/gpfs/bbp.cscs.ch/project/proj12/jenkins/cellular/circuit-scx-v5-gapjunctions/gap_junctions/nrn_gj.h5" # noqa + nrn_reader = GapJunctionSynapseReader.create(nrn_file, 1) + syn_params_nrn = nrn_reader._load_synapse_parameters(100124) + # check reading of sgid, junction_id_pre and junction_id_post + ref_sgids = np.array([94669., 94723., 95634., 95823., 96581., + 97338., 97455., 98139., 98432., 100725., + 101360., 101506., 101696., 101696., 191567.]) + ref_junction_id_pre = np.array([735., 736., 29., 36., 51., + 77., 744., 134., 148., 286., + 322., 337., 355., 356., 681.]) + ref_junction_id_post = np.array([1251., 1259., 617., 1354., 1002., + 1756., 1027., 924., 709., 624., + 1050., 521., 592., 593., 590.]) + npt.assert_allclose(syn_params_nrn.sgid, ref_sgids) + npt.assert_allclose(syn_params_nrn.D, ref_junction_id_pre) + npt.assert_allclose(syn_params_nrn.F, ref_junction_id_post) + + +def test_gapjunction_sonata_reader(): + SIM_DIR = Path(__file__).parent.absolute() / "simulations" + sonata_file = str(SIM_DIR / "mini_thalamus_sonata/gapjunction/edges.h5") + sonata_reader = GapJunctionSynapseReader.create(sonata_file, 1) + syn_params_sonata = sonata_reader._load_synapse_parameters(1) + ref_junction_id_pre = np.array([10257., 43930., 226003., 298841., 324744., + 1094745., 1167632., 1172523., 1260104.]) + ref_junction_id_post = np.array([14., 52., 71., 76., 78., 84., 89., 90., 93.]) + ref_weight = np.array([0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2, 0.2]) + npt.assert_allclose(syn_params_sonata.D, ref_junction_id_pre) + npt.assert_allclose(syn_params_sonata.F, ref_junction_id_post) + npt.assert_allclose(syn_params_sonata.weight, ref_weight)