From 95965f8c22a9206b84670290e06d52eeec1a40d3 Mon Sep 17 00:00:00 2001 From: sfmig <33267254+sfmig@users.noreply.github.com> Date: Fri, 6 Sep 2024 09:31:28 +0100 Subject: [PATCH] Simplify and expand kinematic tests for bboxes (#265) * Simplify and expand kinematic tests * Suggestion to rename internal method for clarity * Update docstrings for API reference docs * Add test for values * Refactor test for uniform linear motion * Add notes to dosctrings * Fix kinematics tests * Add fixture with uniform linear motion for poses * Add poses dataset to linear uniform motion test * Add test for dataset with nans * Edits to docstrings * Remove circular fixture * Small edits to fixture comments * Edits to comments in tests * Small edits * Clarify vector vs array in docstrings and make consistent where required * Add missing docstring in test and small edits * Remove TODOs * Fix offset in fixture for uniform linear motion poses * Apply suggestions from code review Co-authored-by: Chang Huan Lo * Differentiation method Co-authored-by: Chang Huan Lo * Docstrings fixes * :py:meth: to :meth: * Combine into one paragraph * Add uniform linear motion to doscstring of fixtures * Simplify valid_poses_array_uniform_linear_motion with suggestions * kinematic_variable --> kinematic_array * Simplify test_kinematics_uniform_linear_motion with suggestion from review * Update tests/test_unit/test_kinematics.py Co-authored-by: Chang Huan Lo * Update tests/test_unit/test_kinematics.py Co-authored-by: Chang Huan Lo * Cosmetic edits to test * Change docstring to time-derivative --------- Co-authored-by: Chang Huan Lo --- movement/analysis/kinematics.py | 107 +++++++---- tests/conftest.py | 129 +++++++++++-- tests/test_unit/test_kinematics.py | 278 ++++++++++++++++++----------- 3 files changed, 365 insertions(+), 149 deletions(-) diff --git a/movement/analysis/kinematics.py b/movement/analysis/kinematics.py index 15375a53..ed2b4b30 100644 --- a/movement/analysis/kinematics.py +++ b/movement/analysis/kinematics.py @@ -6,24 +6,40 @@ def compute_displacement(data: xr.DataArray) -> xr.DataArray: - """Compute displacement between consecutive positions. + """Compute displacement array in cartesian coordinates. - Displacement is the difference between consecutive positions - of each keypoint for each individual across ``time``. - At each time point ``t``, it is defined as a vector in - cartesian ``(x,y)`` coordinates, pointing from the previous - ``(t-1)`` to the current ``(t)`` position. + The displacement array is defined as the difference between the position + array at time point ``t`` and the position array at time point ``t-1``. + + As a result, for a given individual and keypoint, the displacement vector + at time point ``t``, is the vector pointing from the previous + ``(t-1)`` to the current ``(t)`` position, in cartesian coordinates. Parameters ---------- data : xarray.DataArray - The input data containing position information, with - ``time`` as a dimension. + The input data array containing position vectors in cartesian + coordinates, with ``time`` as a dimension. Returns ------- xarray.DataArray - An xarray DataArray containing the computed displacement. + An xarray DataArray containing displacement vectors in cartesian + coordinates. + + Notes + ----- + For the ``position`` array of a ``poses`` dataset, the ``displacement`` + array will hold the displacement vectors for every keypoint and every + individual. + + For the ``position`` array of a ``bboxes`` dataset, the ``displacement`` + array will hold the displacement vectors for the centroid of every + individual bounding box. + + For the ``shape`` array of a ``bboxes`` dataset, the + ``displacement`` array will hold vectors with the change in width and + height per bounding box, between consecutive time points. """ _validate_time_dimension(data) @@ -33,52 +49,73 @@ def compute_displacement(data: xr.DataArray) -> xr.DataArray: def compute_velocity(data: xr.DataArray) -> xr.DataArray: - """Compute the velocity in cartesian ``(x,y)`` coordinates. + """Compute velocity array in cartesian coordinates. - Velocity is the first derivative of position for each keypoint - and individual across ``time``, computed with the second order - accurate central differences. + The velocity array is the first time-derivative of the position + array. It is computed by applying the second-order accurate central + differences method on the position array. Parameters ---------- data : xarray.DataArray - The input data containing position information, with - ``time`` as a dimension. + The input data array containing position vectors in cartesian + coordinates, with ``time`` as a dimension. Returns ------- xarray.DataArray - An xarray DataArray containing the computed velocity. + An xarray DataArray containing velocity vectors in cartesian + coordinates. + + Notes + ----- + For the ``position`` array of a ``poses`` dataset, the ``velocity`` array + will hold the velocity vectors for every keypoint and every individual. + + For the ``position`` array of a ``bboxes`` dataset, the ``velocity`` array + will hold the velocity vectors for the centroid of every individual + bounding box. See Also -------- - :py:meth:`xarray.DataArray.differentiate` : The underlying method used. + :meth:`xarray.DataArray.differentiate` : The underlying method used. """ return _compute_approximate_time_derivative(data, order=1) def compute_acceleration(data: xr.DataArray) -> xr.DataArray: - """Compute acceleration in cartesian ``(x,y)`` coordinates. + """Compute acceleration array in cartesian coordinates. - Acceleration is the second derivative of position for each keypoint - and individual across ``time``, computed with the second order - accurate central differences. + The acceleration array is the second time-derivative of the + position array. It is computed by applying the second-order accurate + central differences method on the velocity array. Parameters ---------- data : xarray.DataArray - The input data containing position information, with - ``time`` as a dimension. + The input data array containing position vectors in cartesian + coordinates, with``time`` as a dimension. Returns ------- xarray.DataArray - An xarray DataArray containing the computed acceleration. + An xarray DataArray containing acceleration vectors in cartesian + coordinates. + + Notes + ----- + For the ``position`` array of a ``poses`` dataset, the ``acceleration`` + array will hold the acceleration vectors for every keypoint and every + individual. + + For the ``position`` array of a ``bboxes`` dataset, the ``acceleration`` + array will hold the acceleration vectors for the centroid of every + individual bounding box. See Also -------- - :py:meth:`xarray.DataArray.differentiate` : The underlying method used. + :meth:`xarray.DataArray.differentiate` : The underlying method used. """ return _compute_approximate_time_derivative(data, order=2) @@ -87,24 +124,26 @@ def compute_acceleration(data: xr.DataArray) -> xr.DataArray: def _compute_approximate_time_derivative( data: xr.DataArray, order: int ) -> xr.DataArray: - """Compute the derivative using numerical differentiation. + """Compute the time-derivative of an array using numerical differentiation. - This function uses :py:meth:`xarray.DataArray.differentiate`, - which differentiates the array with the second order - accurate central differences. + This function uses :meth:`xarray.DataArray.differentiate`, + which differentiates the array with the second-order + accurate central differences method. Parameters ---------- data : xarray.DataArray - The input data containing ``time`` as a dimension. + The input data array containing ``time`` as a dimension. order : int - The order of the derivative. 1 for velocity, 2 for - acceleration. Value must be a positive integer. + The order of the time-derivative. For an input containing position + data, use 1 to compute velocity, and 2 to compute acceleration. Value + must be a positive integer. Returns ------- xarray.DataArray - An xarray DataArray containing the derived variable. + An xarray DataArray containing the time-derivative of the + input data. """ if not isinstance(order, int): @@ -113,7 +152,9 @@ def _compute_approximate_time_derivative( ) if order <= 0: raise log_error(ValueError, "Order must be a positive integer.") + _validate_time_dimension(data) + result = data for _ in range(order): result = result.differentiate("time") diff --git a/tests/conftest.py b/tests/conftest.py index f2c77bed..272e5eaa 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -238,12 +238,18 @@ def valid_bboxes_arrays_all_zeros(): # --------------------- Bboxes dataset fixtures ---------------------------- @pytest.fixture -def valid_bboxes_array(): - """Return a dictionary of valid non-zero arrays for a - ValidBboxesDataset. - - Contains realistic data for 10 frames, 2 individuals, in 2D - with 5 low confidence bounding boxes. +def valid_bboxes_arrays(): + """Return a dictionary of valid arrays for a + ValidBboxesDataset representing a uniform linear motion. + + It represents 2 individuals for 10 frames, in 2D space. + - Individual 0 moves along the x=y line from the origin. + - Individual 1 moves along the x=-y line line from the origin. + + All confidence values are set to 0.9 except the following which are set + to 0.1: + - Individual 0 at frames 2, 3, 4 + - Individual 1 at frames 2, 3 """ # define the shape of the arrays n_frames, n_individuals, n_space = (10, 2, 2) @@ -276,22 +282,21 @@ def valid_bboxes_array(): "position": position, "shape": shape, "confidence": confidence, - "individual_names": ["id_" + str(id) for id in range(n_individuals)], } @pytest.fixture def valid_bboxes_dataset( - valid_bboxes_array, + valid_bboxes_arrays, ): - """Return a valid bboxes dataset with low confidence values and - time in frames. + """Return a valid bboxes dataset for two individuals moving in uniform + linear motion, with 5 frames with low confidence values and time in frames. """ dim_names = MovementDataset.dim_names["bboxes"] - position_array = valid_bboxes_array["position"] - shape_array = valid_bboxes_array["shape"] - confidence_array = valid_bboxes_array["confidence"] + position_array = valid_bboxes_arrays["position"] + shape_array = valid_bboxes_arrays["shape"] + confidence_array = valid_bboxes_arrays["confidence"] n_frames, n_individuals, _ = position_array.shape @@ -409,12 +414,110 @@ def valid_poses_dataset(valid_position_array, request): @pytest.fixture def valid_poses_dataset_with_nan(valid_poses_dataset): """Return a valid pose tracks dataset with NaN values.""" + # Sets position for all keypoints in individual ind1 to NaN + # at timepoints 3, 7, 8 valid_poses_dataset.position.loc[ {"individuals": "ind1", "time": [3, 7, 8]} ] = np.nan return valid_poses_dataset +@pytest.fixture +def valid_poses_array_uniform_linear_motion(): + """Return a dictionary of valid arrays for a + ValidPosesDataset representing a uniform linear motion. + + It represents 2 individuals with 3 keypoints, for 10 frames, in 2D space. + - Individual 0 moves along the x=y line from the origin. + - Individual 1 moves along the x=-y line line from the origin. + + All confidence values for all keypoints are set to 0.9 except + for the keypoints at the following frames which are set to 0.1: + - Individual 0 at frames 2, 3, 4 + - Individual 1 at frames 2, 3 + """ + # define the shape of the arrays + n_frames, n_individuals, n_keypoints, n_space = (10, 2, 3, 2) + + # define centroid (index=0) trajectory in position array + # for each individual, the centroid moves along + # the x=+/-y line, starting from the origin. + # - individual 0 moves along x = y line + # - individual 1 moves along x = -y line + # They move one unit along x and y axes in each frame + frames = np.arange(n_frames) + position = np.empty((n_frames, n_individuals, n_keypoints, n_space)) + position[:, :, 0, 0] = frames[:, None] # reshape to (n_frames, 1) + position[:, 0, 0, 1] = frames + position[:, 1, 0, 1] = -frames + + # define trajectory of left and right keypoints + # for individual 0, at each timepoint: + # - the left keypoint (index=1) is at x_centroid, y_centroid + 1 + # - the right keypoint (index=2) is at x_centroid + 1, y_centroid + # for individual 1, at each timepoint: + # - the left keypoint (index=1) is at x_centroid - 1, y_centroid + # - the right keypoint (index=2) is at x_centroid, y_centroid + 1 + offsets = [ + [(0, 1), (1, 0)], # individual 0: left, right keypoints (x,y) offsets + [(-1, 0), (0, 1)], # individual 1: left, right keypoints (x,y) offsets + ] + for i in range(n_individuals): + for kpt in range(1, n_keypoints): + position[:, i, kpt, 0] = ( + position[:, i, 0, 0] + offsets[i][kpt - 1][0] + ) + position[:, i, kpt, 1] = ( + position[:, i, 0, 1] + offsets[i][kpt - 1][1] + ) + + # build an array of confidence values, all 0.9 + confidence = np.full((n_frames, n_individuals, n_keypoints), 0.9) + # set 5 low-confidence values + # - set 3 confidence values for individual id_0's centroid to 0.1 + # - set 2 confidence values for individual id_1's centroid to 0.1 + idx_start = 2 + confidence[idx_start : idx_start + 3, 0, 0] = 0.1 + confidence[idx_start : idx_start + 2, 1, 0] = 0.1 + + return {"position": position, "confidence": confidence} + + +@pytest.fixture +def valid_poses_dataset_uniform_linear_motion( + valid_poses_array_uniform_linear_motion, +): + """Return a valid poses dataset for two individuals moving in uniform + linear motion, with 5 frames with low confidence values and time in frames. + """ + dim_names = MovementDataset.dim_names["poses"] + + position_array = valid_poses_array_uniform_linear_motion["position"] + confidence_array = valid_poses_array_uniform_linear_motion["confidence"] + + n_frames, n_individuals, _, _ = position_array.shape + + return xr.Dataset( + data_vars={ + "position": xr.DataArray(position_array, dims=dim_names), + "confidence": xr.DataArray(confidence_array, dims=dim_names[:-1]), + }, + coords={ + dim_names[0]: np.arange(n_frames), + dim_names[1]: [f"id_{i}" for i in range(1, n_individuals + 1)], + dim_names[2]: ["centroid", "left", "right"], + dim_names[3]: ["x", "y"], + }, + attrs={ + "fps": None, + "time_unit": "frames", + "source_software": "test", + "source_file": "test_poses.h5", + "ds_type": "poses", + }, + ) + + # -------------------- Invalid datasets fixtures ------------------------------ @pytest.fixture def not_a_dataset(): diff --git a/tests/test_unit/test_kinematics.py b/tests/test_unit/test_kinematics.py index 1f75a824..66822bfd 100644 --- a/tests/test_unit/test_kinematics.py +++ b/tests/test_unit/test_kinematics.py @@ -1,5 +1,3 @@ -from contextlib import nullcontext as does_not_raise - import numpy as np import pytest import xarray as xr @@ -7,106 +5,180 @@ from movement.analysis import kinematics -class TestKinematics: - """Test suite for the kinematics module.""" - - @pytest.fixture - def expected_dataarray(self, valid_poses_dataset): - """Return a function to generate the expected dataarray - for different kinematic properties. - """ - - def _expected_dataarray(property): - """Return an xarray.DataArray with default values and - the expected dimensions and coordinates. - """ - # Expected x,y values for velocity - x_vals = np.array( - [1.0, 2.0, 4.0, 6.0, 8.0, 10.0, 12.0, 14.0, 16.0, 17.0] - ) - y_vals = np.full((10, 2, 2, 1), 4.0) - if property == "acceleration": - x_vals = np.array( - [1.0, 1.5, 2.0, 2.0, 2.0, 2.0, 2.0, 2.0, 1.5, 1.0] - ) - y_vals = np.full((10, 2, 2, 1), 0) - elif property == "displacement": - x_vals = np.array( - [0.0, 1.0, 3.0, 5.0, 7.0, 9.0, 11.0, 13.0, 15.0, 17.0] - ) - y_vals[0] = 0 - - x_vals = x_vals.reshape(-1, 1, 1, 1) - # Repeat the x_vals to match the shape of the position - x_vals = np.tile(x_vals, (1, 2, 2, 1)) - return xr.DataArray( - np.concatenate( - [x_vals, y_vals], - axis=-1, - ), - dims=valid_poses_dataset.dims, - coords=valid_poses_dataset.coords, - ) - - return _expected_dataarray - - kinematic_test_params = [ - ("valid_poses_dataset", does_not_raise()), - ("valid_poses_dataset_with_nan", does_not_raise()), - ("missing_dim_poses_dataset", pytest.raises(ValueError)), - ] +@pytest.mark.parametrize( + "valid_dataset_uniform_linear_motion", + [ + "valid_poses_dataset_uniform_linear_motion", + "valid_bboxes_dataset", + ], +) +@pytest.mark.parametrize( + "kinematic_variable, expected_kinematics", + [ + ( + "displacement", + [ + np.vstack([np.zeros((1, 2)), np.ones((9, 2))]), # Individual 0 + np.multiply( + np.vstack([np.zeros((1, 2)), np.ones((9, 2))]), + np.array([1, -1]), + ), # Individual 1 + ], + ), + ( + "velocity", + [ + np.ones((10, 2)), # Individual 0 + np.multiply( + np.ones((10, 2)), np.array([1, -1]) + ), # Individual 1 + ], + ), + ( + "acceleration", + [ + np.zeros((10, 2)), # Individual 0 + np.zeros((10, 2)), # Individual 1 + ], + ), + ], +) +def test_kinematics_uniform_linear_motion( + valid_dataset_uniform_linear_motion, + kinematic_variable, + expected_kinematics, # 2D: n_frames, n_space_dims + request, +): + """Test computed kinematics for a uniform linear motion case. + + Uniform linear motion means the individuals move along a line + at constant velocity. + + We consider 2 individuals ("id_0" and "id_1"), + tracked for 10 frames, along x and y: + - id_0 moves along x=y line from the origin + - id_1 moves along x=-y line from the origin + - they both move one unit (pixel) along each axis in each frame + + If the dataset is a poses dataset, we consider 3 keypoints per individual + (centroid, left, right), that are always in front of the centroid keypoint + at 45deg from the trajectory. + """ + # Compute kinematic array from input dataset + position = request.getfixturevalue( + valid_dataset_uniform_linear_motion + ).position + kinematic_array = getattr(kinematics, f"compute_{kinematic_variable}")( + position + ) - @pytest.mark.parametrize("ds, expected_exception", kinematic_test_params) - def test_displacement( - self, ds, expected_exception, expected_dataarray, request - ): - """Test displacement computation.""" - ds = request.getfixturevalue(ds) - with expected_exception: - result = kinematics.compute_displacement(ds.position) - expected = expected_dataarray("displacement") - if ds.position.isnull().any(): - expected.loc[ - {"individuals": "ind1", "time": [3, 4, 7, 8, 9]} - ] = np.nan - xr.testing.assert_allclose(result, expected) - - @pytest.mark.parametrize("ds, expected_exception", kinematic_test_params) - def test_velocity( - self, ds, expected_exception, expected_dataarray, request - ): - """Test velocity computation.""" - ds = request.getfixturevalue(ds) - with expected_exception: - result = kinematics.compute_velocity(ds.position) - expected = expected_dataarray("velocity") - if ds.position.isnull().any(): - expected.loc[ - {"individuals": "ind1", "time": [2, 4, 6, 7, 8, 9]} - ] = np.nan - xr.testing.assert_allclose(result, expected) - - @pytest.mark.parametrize("ds, expected_exception", kinematic_test_params) - def test_acceleration( - self, ds, expected_exception, expected_dataarray, request - ): - """Test acceleration computation.""" - ds = request.getfixturevalue(ds) - with expected_exception: - result = kinematics.compute_acceleration(ds.position) - expected = expected_dataarray("acceleration") - if ds.position.isnull().any(): - expected.loc[ - {"individuals": "ind1", "time": [1, 3, 5, 6, 7, 8, 9]} - ] = np.nan - xr.testing.assert_allclose(result, expected) - - @pytest.mark.parametrize("order", [0, -1, 1.0, "1"]) - def test_approximate_derivative_with_invalid_order(self, order): - """Test that an error is raised when the order is non-positive.""" - data = np.arange(10) - expected_exception = ( - ValueError if isinstance(order, int) else TypeError + # Build expected data array from the expected numpy array + expected_array = xr.DataArray( + np.stack(expected_kinematics, axis=1), + # Stack along the "individuals" axis + dims=["time", "individuals", "space"], + ) + if "keypoints" in position.coords: + expected_array = expected_array.expand_dims( + {"keypoints": position.coords["keypoints"].size} ) - with pytest.raises(expected_exception): - kinematics._compute_approximate_time_derivative(data, order=order) + expected_array = expected_array.transpose( + "time", "individuals", "keypoints", "space" + ) + + # Compare the values of the kinematic_array against the expected_array + np.testing.assert_allclose(kinematic_array.values, expected_array.values) + + +@pytest.mark.parametrize( + "valid_dataset_with_nan", + [ + "valid_poses_dataset_with_nan", + "valid_bboxes_dataset_with_nan", + ], +) +@pytest.mark.parametrize( + "kinematic_variable, expected_nans_per_individual", + [ + ("displacement", [5, 0]), # individual 0, individual 1 + ("velocity", [6, 0]), + ("acceleration", [7, 0]), + ], +) +def test_kinematics_with_dataset_with_nans( + valid_dataset_with_nan, + kinematic_variable, + expected_nans_per_individual, + helpers, + request, +): + """Test kinematics computation for a dataset with nans. + + We test that the kinematics can be computed and that the number + of nan values in the kinematic array is as expected. + + """ + # compute kinematic array + valid_dataset = request.getfixturevalue(valid_dataset_with_nan) + position = valid_dataset.position + kinematic_array = getattr(kinematics, f"compute_{kinematic_variable}")( + position + ) + + # compute n nans in kinematic array per individual + n_nans_kinematics_per_indiv = [ + helpers.count_nans(kinematic_array.isel(individuals=i)) + for i in range(valid_dataset.dims["individuals"]) + ] + + # expected nans per individual adjusted for space and keypoints dimensions + expected_nans_adjusted = [ + n + * valid_dataset.dims["space"] + * valid_dataset.dims.get("keypoints", 1) + for n in expected_nans_per_individual + ] + # check number of nans per individual is as expected in kinematic array + np.testing.assert_array_equal( + n_nans_kinematics_per_indiv, expected_nans_adjusted + ) + + +@pytest.mark.parametrize( + "invalid_dataset, expected_exception", + [ + ("not_a_dataset", pytest.raises(AttributeError)), + ("empty_dataset", pytest.raises(AttributeError)), + ("missing_var_poses_dataset", pytest.raises(AttributeError)), + ("missing_var_bboxes_dataset", pytest.raises(AttributeError)), + ("missing_dim_poses_dataset", pytest.raises(ValueError)), + ("missing_dim_bboxes_dataset", pytest.raises(ValueError)), + ], +) +@pytest.mark.parametrize( + "kinematic_variable", + [ + "displacement", + "velocity", + "acceleration", + ], +) +def test_kinematics_with_invalid_dataset( + invalid_dataset, + expected_exception, + kinematic_variable, + request, +): + """Test kinematics computation with an invalid dataset.""" + with expected_exception: + position = request.getfixturevalue(invalid_dataset).position + getattr(kinematics, f"compute_{kinematic_variable}")(position) + + +@pytest.mark.parametrize("order", [0, -1, 1.0, "1"]) +def test_approximate_derivative_with_invalid_order(order): + """Test that an error is raised when the order is non-positive.""" + data = np.arange(10) + expected_exception = ValueError if isinstance(order, int) else TypeError + with pytest.raises(expected_exception): + kinematics._compute_approximate_time_derivative(data, order=order)