From 913a99826077bbb87e8eafd7cfe87c3440e7a524 Mon Sep 17 00:00:00 2001 From: Moritz Kern <92092328+Moritz-Alexander-Kern@users.noreply.github.com> Date: Tue, 24 Oct 2023 09:53:42 +0200 Subject: [PATCH] [Enh] add example and handle one-dimensional arrays in CSD.generate lfp (#594) * add example to documentation * add handling for 1-dimensional arrays * add regression test, testing for one- and two-dimensional arrays --------- Co-authored-by: Michael Denker Co-authored-by: Alexander Kleinjohann <33096371+Kleinjohann@users.noreply.github.com> --- elephant/current_source_density.py | 50 +++++++++++++++---- ..._csd.py => test_current_source_density.py} | 25 ++++++++++ 2 files changed, 64 insertions(+), 11 deletions(-) rename elephant/test/{test_csd.py => test_current_source_density.py} (88%) diff --git a/elephant/current_source_density.py b/elephant/current_source_density.py index ad46f11cd..fd73d0733 100644 --- a/elephant/current_source_density.py +++ b/elephant/current_source_density.py @@ -71,7 +71,7 @@ def estimate_csd(lfp, coordinates='coordinates', method=None, Parameters ---------- - lfp : neo.AnalogSignal + lfp : :class:`neo.core.AnalogSignal` Positions of electrodes can be added as an array annotation coordinates : array-like Quantity or string Specifies the corresponding spatial coordinates of the electrodes. @@ -79,8 +79,8 @@ def estimate_csd(lfp, coordinates='coordinates', method=None, with dimension of space, where M is the number of signals in 'lfp', and N is equal to the dimensionality of the method. Alternatively, if coordinates is a string, the function will fetch the - coordinates, supplied in the same format, as annotation of 'lfp' by that - name. + coordinates, supplied in the same format, as annotation of 'lfp' by + that name. Default: 'coordinates' method : string Pick a method corresponding to the setup, in this implementation @@ -182,14 +182,14 @@ def estimate_csd(lfp, coordinates='coordinates', method=None, # All iCSD methods explicitly assume a source # diameter in contrast to the stdCSD that # implicitly assume infinite source radius - raise ValueError("Parameter diam must be specified for iCSD \ - methods: {}".format(", ".join(icsd_methods))) + raise ValueError(f"Parameter diam must be specified for iCSD " + f"methods: {', '.join(icsd_methods)}") if 'f_type' in kwargs: if (kwargs['f_type'] != 'identity') and \ (kwargs['f_order'] is None): - raise ValueError("The order of {} filter must be \ - specified".format(kwargs['f_type'])) + raise ValueError(f"The order of {kwargs['f_type']} filter must" + f" be specified") csd_method = getattr(icsd, method) # fetch class from icsd.py file csd_estimator = csd_method(lfp=lfp.T.magnitude * lfp.units, @@ -228,12 +228,15 @@ def generate_lfp(csd_profile, x_positions, y_positions=None, z_positions=None, 2D : large_source_2D and small_source_2D 3D : gauss_3d_dipole x_positions : np.ndarray - Positions of the x coordinates of the electrodes + A 2D column vector (N x 1 array) containing the positions of the x + coordinates of the electrodes y_positions : np.ndarray, optional - Positions of the y coordinates of the electrodes + A 2D column vector (N x 1 array) containing the positions of the y + coordinates of the electrodes Defaults to None, use in 2D or 3D cases only z_positions : np.ndarray, optional - Positions of the z coordinates of the electrodes + A 2D column vector (N x 1 array) containing the positions of the z + coordinates of the electrodes Defaults to None, use in 3D case only x_limits : list, optional A list of [start, end]. @@ -253,10 +256,31 @@ def generate_lfp(csd_profile, x_positions, y_positions=None, z_positions=None, Returns ------- - LFP : neo.AnalogSignal + LFP : :class:`neo.core.AnalogSignal` The potentials created by the csd profile at the electrode positions. The electrode positions are attached as an annotation named 'coordinates'. + + Examples + -------- + >>> import numpy as np + >>> from elephant.current_source_density import generate_lfp, estimate_csd + >>> from elephant.current_source_density_src.utility_functions import gauss_1d_dipole # noqa + >>> # 1. Define an array xs to x coordinate values with a length of 2304 + >>> xs=np.linspace(0, 10, 2304) + + >>> # 2. Run generate_lfp(gauss_1d_dipole, xs) + >>> lfp = generate_lfp(gauss_1d_dipole, xs) + + >>> # 3. Run estimate_csd(lfp, method="StandardCSD") + >>> csd = estimate_csd(lfp, method="StandardCSD") #doctest: +ELLIPSIS + discrete ... + >>> # 4. Print the results + >>> print(f"LFPs: {lfp}") + LFPs: [[-0.01483716 -0.01483396 -0.01483075 ... 0.01219233 0.0121911 + 0.01218986]] mV + >>> print(f"CSD estimate: {csd}") #doctest: +ELLIPSIS + CSD estimate: [[-1.00025592e-04 -6.06684588e-05 ... """ def integrate_1D(x0, csd_x, csd, h): @@ -297,6 +321,10 @@ def integrate_3D(x, y, z, csd, xlin, ylin, zlin, X, Y, Z): sigma = 1.0 h = 50. if dim == 1: + # Handle one dimensional case, + # see https://github.com/NeuralEnsemble/elephant/issues/546 + if len(x_positions.shape) == 1: + x_positions = np.expand_dims(x_positions, axis=1) chrg_x = x csd = csd_profile(chrg_x) pots = integrate_1D(x_positions, chrg_x, csd, h) diff --git a/elephant/test/test_csd.py b/elephant/test/test_current_source_density.py similarity index 88% rename from elephant/test/test_csd.py rename to elephant/test/test_current_source_density.py index f020553c0..b17ba3f89 100644 --- a/elephant/test/test_csd.py +++ b/elephant/test/test_current_source_density.py @@ -15,6 +15,7 @@ import quantities as pq from elephant import current_source_density as csd import elephant.current_source_density_src.utility_functions as utils +from elephant.current_source_density import generate_lfp available_1d = ['StandardCSD', 'DeltaiCSD', 'StepiCSD', 'SplineiCSD', 'KCSD1D'] available_2d = ['KCSD2D', 'MoIKCSD'] @@ -153,5 +154,29 @@ def test_kcsd2d_init(self): self.assertEqual(len(result.times), 1) +class GenerateLfpTestCase(unittest.TestCase): + @classmethod + def setUpClass(cls) -> None: + cls.one_dimensional = np.linspace(0, 10, 2304) + cls.two_dimensional = np.linspace(0, 10, 2304 + ).reshape(2304, 1) + + def test_generate_lfp_one_dimensional_array(self): + """ + Regression test for Issue #546, + see: https://github.com/NeuralEnsemble/elephant/issues/546 + """ + # this should raise NOT an error + generate_lfp(utils.gauss_1d_dipole, self.one_dimensional) + + def test_generate_lfp_two_dimensional_array(self): + """ + Regression test for Issue #546, + see: https://github.com/NeuralEnsemble/elephant/issues/546 + """ + # this should NOT raise an error + generate_lfp(utils.gauss_1d_dipole, self.two_dimensional) + + if __name__ == '__main__': unittest.main()