[MRG] compute yscale from cell points #327

jasmainak · 2021-05-07T13:45:11Z

closes #326

I tried not to touch too much so others can continue with their PRs

jasmainak · 2021-05-07T13:46:48Z

hnn_core/cell.py

@@ -293,18 +306,28 @@ def move_to_pos(self):
    # 2. a list needs to be created with a Dipole (Point Process) in each
    #    section at position 1
    # In Cell() and not Pyr() for future possibilities
-    def insert_dipole(self, yscale):
+    def insert_dipole(self, p_secs, dipole_sec):


it's becoming clear that p_secs should be an attribute of Cell but I'm holding off on that in this PR.

I don't see why p_secs needs to be an attribute of Cell since it isn't used anywhere outside Cell.build() and here. Since Cell.insert_dipole() is called after build(), why not just access the points via the NEURON Section objects?

I'm imagining that plot_cell_morphology will move from Network to Cell and then self.p_secs will be used for the plot. I'm a little wary of relying on NEURON objects because of : 1) numerical issues we saw in the past 2) unpicklability. Since p_secs is picklable we can make it accessible to the user to inspect/modify etc.

Fair enough!

I agree that to enable Cell.plot we should keep HOC-objects out of ours for as long as possible.

codecov-commenter · 2021-05-07T17:55:44Z

Codecov Report

Merging #327 (b8b75b3) into master (580fee6) will increase coverage by 0.02%.
The diff coverage is 93.93%.

@@            Coverage Diff             @@
##           master     #327      +/-   ##
==========================================
+ Coverage   90.31%   90.33%   +0.01%     
==========================================
  Files          13       13              
  Lines        2416     2421       +5     
==========================================
+ Hits         2182     2187       +5     
  Misses        234      234

Impacted Files	Coverage Δ
hnn_core/viz.py	`87.50% <33.33%> (ø)`
hnn_core/__init__.py	`100.00% <100.00%> (ø)`
hnn_core/cell.py	`97.15% <100.00%> (+0.13%)`	⬆️
hnn_core/cells_default.py	`97.82% <100.00%> (ø)`
hnn_core/params_default.py	`100.00% <100.00%> (ø)`

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 580fee6...b8b75b3. Read the comment docs.

hnn_core/cell.py

hnn_core/params_default.py

doc/whats_new.rst

hnn_core/cell.py

cjayb · 2021-05-08T04:31:38Z

hnn_core/cell.py

 from neuron import h, nrn

 # Units for e: mV
 # Units for gbar: S/cm^2


+def _get_yscale(p_secs, dipole_sec):
+    """Get cos theta to compute dipole in z direction."""
+    a = (np.array(p_secs[dipole_sec]['sec_pts'][1]) -


Here we have the hard assumption that our cells live in the xy-plane. Not much to do about it, just a note.

Not really. This just projects along the apical dendrite which may or may not be oriented along the y direction (Neuron coordinates ... z direction in physical coordinates).

There are the following two alternatives which I considered and that does assume the xy-plane:

y_long = (h.y3d(1, sec=sect) - h.y3d(0, sec=sect)) * pos_all

and

y_long = (p_secs[sec_name]['sec_pts'][1][1] - p_secs['apical_trunk']['sec_pts'][0][1]) * pos_all

and get rid of yscale altogether. That didn't work though because of numerical issues I suspect. The code uses a mix of Neuron and Python. At some point we should figure out what kind of rounding/precision Neuron actually uses. We could deliberately typecast in Python then to make things match and then remove it when we're confident it's correct

Not really. This just projects along the apical dendrite which may or may not be oriented along the y direction

But isn't index 1 in p_secs[dipole_sec]['sec_pts'][1] explicitly the y-coordinate? And 0 on the next line x. That's what I mean by "hard assumptions"---we're doing 2D projection, implicitly ignoring the z-orientation, which is assumed to be degenerate (flat).

But isn't index 1 in p_secs[dipole_sec]['sec_pts'][1] explicitly the y-coordinate?

actually no, the default Neuron has only two points and this is the second point. To get the y-coordinate of the second point you would have to do:

p_secs[dipole_sec]['sec_pts'][1][1]

Ah, my bad. All clear now!

cjayb · 2021-05-08T04:33:01Z

hnn_core/cell.py

@@ -293,18 +306,28 @@ def move_to_pos(self):
    # 2. a list needs to be created with a Dipole (Point Process) in each
    #    section at position 1
    # In Cell() and not Pyr() for future possibilities
-    def insert_dipole(self, yscale):
+    def insert_dipole(self, p_secs, dipole_sec):


I agree that to enable Cell.plot we should keep HOC-objects out of ours for as long as possible.

hnn_core/cell.py

cjayb · 2021-05-08T04:43:45Z

hnn_core/cell.py

+            The properties of syn are specified in p_syn.
+            mech is a dict with keys as the mechanism names. The
+            values are dictionaries with properties of the mechanism.
+        dipole_sec : str


Along the lines of my other comments, could this be 'apical_section'?

The name of a section defining the apical orientation, onto which dipole moments are projected.

isn't "apical" specific to the default neuron? I like the use of term "projection" though. I'll try to incorporate it

In my mind, "apical" refers to the single large dendrite sprouting from the top of the pyramid-formed soma (from the "apex" of the pyramid). Some L2/3 pyramidal cells have a bifurcation very close to the soma, but even then the two apical dendrites run in parallel, towards L1. The "textbook answer" is that the apical dendrite(s) is a defining characteristic of all neocortical (indeed, also hippocampal) pyramidal cells. The terms 'proximal' and 'distal' are much more ambiguous in comparison.

humm ... so are you saying that you'd always want to compute the dipole along an apical dendrite? I'll change the name in that case

Yes. Even if we one day included a more morphologically detailed model, the dipole moment would still be calculated along the primary symmetry axis defined by the apical dendrite. You know how in MNE we have the "loose orientation constraint"? The "looseness" is fudge factor to cover for coreg and surface recon errors. The underlying strong assumption is that there's one dominant dipole orientation per source space point. And this is the normal to the cortical (WM) surface because of the apical dendrites.

jasmainak · 2021-05-08T12:33:17Z

hnn_core/cell.py

-            sect_name = sect.name().split('_', 1)[1]
-            y_scale = (yscale[sect_name] * sect.L) * pos_all
-            # y_long = (h.y3d(1, sec=sect) - h.y3d(0, sec=sect)) * pos
+            pos_proj = sect.L * cos_thetas[sect_name] * pos_all


this scaling is super mysterious to me. We're doing:

diff(section length * position vector of segment)
= section length * segment lengths

Yeah I was trying to make sense of that too... I think it's very ugly shorthand. The L and theta are constant, so the diff is effectively on the position, which is in units of section length ([0, 1]). This may have something to do with calculating the dipole moment of current, but it all depends on what Dipole.mod and Dipole_pp.mod are actually doing... I'm digging into the details of HOC (trying to figure out my LFP-woes). I have plans to visit the "Dipole"-calculations, too, hoping this part of our code will make more sense!

oh I get it now. sect.L * np.diff(pos_all) is the fractional length (in units of um) of that segment. Multiply that by the current in the segment and you get the dipole moment of that segment

Indeed, but that's what I mean by having to understand the MOD-files better. Because then they must not be calculating the dipole moment, right? But rather the primary current density...

I think you're right. I'm guessing Qsum is the the sum of all radial current entering the dendrite at a given point in a nrn section?

cjayb

Final nitpicks, good to merge!

cjayb · 2021-05-09T05:47:59Z

hnn_core/cell.py

+    """Get cos(theta) to compute dipole along the apical dendrite."""
+    a = (np.array(p_secs[sec_name_apical]['sec_pts'][1]) -
+         np.array(p_secs[sec_name_apical]['sec_pts'][0]))
+    yscales = dict()


Suggested change

yscales = dict()

cos_theta = dict()

cjayb · 2021-05-09T05:48:14Z

hnn_core/cell.py

+    yscales = dict()
+    for sec_name, p_sec in p_secs.items():
+        b = np.array(p_sec['sec_pts'][1]) - np.array(p_sec['sec_pts'][0])
+        yscales[sec_name] = np.dot(a, b) / (norm(a) * norm(b))


Suggested change

yscales[sec_name] = np.dot(a, b) / (norm(a) * norm(b))

cos_theta[sec_name] = np.dot(a, b) / (norm(a) * norm(b))

cjayb · 2021-05-09T05:48:27Z

hnn_core/cell.py

+    for sec_name, p_sec in p_secs.items():
+        b = np.array(p_sec['sec_pts'][1]) - np.array(p_sec['sec_pts'][0])
+        yscales[sec_name] = np.dot(a, b) / (norm(a) * norm(b))
+    return yscales


Suggested change

return yscales

return cos_theta

cjayb · 2021-05-09T05:49:10Z

hnn_core/cell.py

+            mech is a dict with keys as the mechanism names. The
+            values are dictionaries with properties of the mechanism.
+        sec_name_apical : str
+            The name of the section along which dipole is to


Suggested change

The name of the section along which dipole is to

The name of the section along which dipole moment is calculated.

WDYT?

jasmainak · 2021-05-09T13:41:34Z

Addressed! Thanks for the careful reviews :)

cjayb · 2021-05-09T15:22:24Z

Thanks for the much more readable code @jasmainak, merged after green light from @rythorpe

jasmainak added 2 commits May 7, 2021 09:21

ENH: cos theta from section points

ab4a1c3

ENH: _get_yscale in cell.py

40068a0

jasmainak changed the title ~~[WIP] compute yscale from cell points~~ [MRG] compute yscale from cell points May 7, 2021

jasmainak requested review from ntolley and rythorpe May 7, 2021 13:45

jasmainak commented May 7, 2021

View reviewed changes

DOC: update whats_new.rst

9a1f276

jasmainak force-pushed the auto_yscale branch from b389132 to 9a1f276 Compare May 7, 2021 13:54

Flake8

d1faf0d

rythorpe reviewed May 7, 2021

View reviewed changes

hnn_core/cell.py Outdated Show resolved Hide resolved

rythorpe reviewed May 7, 2021

View reviewed changes

hnn_core/params_default.py Show resolved Hide resolved

ENH: some more simplificatoin

75d1bbf

jasmainak force-pushed the auto_yscale branch from a962177 to 03910e0 Compare May 7, 2021 23:14

DOC: show as code block + fix links"

c4d9eb9

jasmainak force-pushed the auto_yscale branch from 03910e0 to c4d9eb9 Compare May 7, 2021 23:50

cjayb approved these changes May 8, 2021

View reviewed changes

ENH: clarify code

d719935

jasmainak commented May 8, 2021

View reviewed changes

DOC: better variable names

cad01e0

jasmainak force-pushed the auto_yscale branch from ae4f1bf to cad01e0 Compare May 9, 2021 03:35

MAINT: tiny improvement

ad13f8e

cjayb approved these changes May 9, 2021

View reviewed changes

Chris comments

b8b75b3

rythorpe approved these changes May 9, 2021

View reviewed changes

cjayb merged commit eaa7373 into jonescompneurolab:master May 9, 2021

cjayb mentioned this pull request May 11, 2021

section points and lengths don't match #330

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

[MRG] compute yscale from cell points #327

[MRG] compute yscale from cell points #327

jasmainak commented May 7, 2021

jasmainak May 7, 2021

rythorpe May 7, 2021

jasmainak May 7, 2021

rythorpe May 7, 2021

cjayb May 8, 2021

codecov-commenter commented May 7, 2021 •

edited

cjayb May 8, 2021

jasmainak May 8, 2021 •

edited

cjayb May 8, 2021

jasmainak May 8, 2021

cjayb May 9, 2021

cjayb May 8, 2021

cjayb May 8, 2021

jasmainak May 8, 2021

cjayb May 8, 2021

jasmainak May 8, 2021

cjayb May 9, 2021

jasmainak May 8, 2021

cjayb May 8, 2021

jasmainak May 8, 2021

cjayb May 9, 2021

rythorpe May 9, 2021

cjayb left a comment

cjayb May 9, 2021

cjayb May 9, 2021

cjayb May 9, 2021

cjayb May 9, 2021

jasmainak commented May 9, 2021

cjayb commented May 9, 2021

	yscales[sec_name] = np.dot(a, b) / (norm(a) * norm(b))
	cos_theta[sec_name] = np.dot(a, b) / (norm(a) * norm(b))

	The name of the section along which dipole is to
	The name of the section along which dipole moment is calculated.

[MRG] compute yscale from cell points #327

[MRG] compute yscale from cell points #327

Conversation

jasmainak commented May 7, 2021

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

codecov-commenter commented May 7, 2021 • edited

Codecov Report

Choose a reason for hiding this comment

jasmainak May 8, 2021 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

cjayb left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jasmainak commented May 9, 2021

cjayb commented May 9, 2021

codecov-commenter commented May 7, 2021 •

edited

jasmainak May 8, 2021 •

edited