From bb2f7a1738bc158ce2da778aaffaadc734223296 Mon Sep 17 00:00:00 2001 From: ehennestad Date: Mon, 21 Sep 2026 23:18:16 +0200 Subject: [PATCH 1/5] fix(roisignals): leave roi pixels out of neuropil masks estimateNeuropilMasks ignores the pixels in pixelsToIgnore when it grows a neuropil mask. With an image mask it computed pixelsToIgnore as footprintAllRois & ~imageMask, which is only true for roi pixels that lie outside the image mask. extractF and SignalExtractor replace an empty imageMask with true(imageSize), so pixelsToIgnore was false everywhere: excludeRoiFromNeuropil had no effect, and the neuropil signal of a roi included the pixels of neighboring rois. Pixels outside a user-defined image mask also stayed in the neuropil masks. Combine the two conditions with |, so that a pixel is left out if it belongs to a roi or if it lies outside the image mask. Neuropil signals that were extracted through extractF or SignalExtractor before this change include pixels of neighboring rois. Co-Authored-By: Claude Fable 5.1 --- code/+nansen/+processing/+roi/estimateNeuropilMasks.m | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/code/+nansen/+processing/+roi/estimateNeuropilMasks.m b/code/+nansen/+processing/+roi/estimateNeuropilMasks.m index 3dc2d0bb8..e2591e6d5 100644 --- a/code/+nansen/+processing/+roi/estimateNeuropilMasks.m +++ b/code/+nansen/+processing/+roi/estimateNeuropilMasks.m @@ -25,7 +25,7 @@ end if ~isempty(params.imageMask) - pixelsToIgnore = footprintAllRois & ~params.imageMask; + pixelsToIgnore = footprintAllRois | ~params.imageMask; else pixelsToIgnore = footprintAllRois; end From da56a8bcc653fa74f8110d22490947233f6b65ad Mon Sep 17 00:00:00 2001 From: ehennestad Date: Mon, 21 Sep 2026 23:18:26 +0200 Subject: [PATCH 2/5] fix(roisignals): return the options that extractF used extractF parses its options into params and resolves roiInd, imageMask, extractFcn and roiMaskFormat in updateParameters. The second output was P, the struct of default values that the options are parsed against, so callers always received the defaults. The extractSignals session method saves this output as OptionsSignalExtraction, which therefore did not record options that differed from the defaults. Return params as the second output. Co-Authored-By: Claude Fable 5.1 --- code/+nansen/+twophoton/+roisignals/extractF.m | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/code/+nansen/+twophoton/+roisignals/extractF.m b/code/+nansen/+twophoton/+roisignals/extractF.m index 86a327359..11f9a0b3b 100644 --- a/code/+nansen/+twophoton/+roisignals/extractF.m +++ b/code/+nansen/+twophoton/+roisignals/extractF.m @@ -1,4 +1,4 @@ -function [signalArray, P] = extractF(imageData, roiData, varargin) +function [signalArray, params] = extractF(imageData, roiData, varargin) %extractF Extract RoI signals from an image stack in chunkwise manner % % signals = extractF(imageStack, roiArray) extracts signals for rois in From 52d44a688d4c568398e5024503fe4fa41530dc29 Mon Sep 17 00:00:00 2001 From: ehennestad Date: Mon, 21 Sep 2026 23:18:39 +0200 Subject: [PATCH 3/5] refactor(roisignals): take neuropil mask options as name-value arguments estimateNeuropilMasks parsed varargin with utility.parsenvpairs against all signal extraction parameters, although it uses five of them. It now declares roiInd, imageMask, excludeRoiFromNeuropil, neuropilExpansionFactor and cellNeuropilSeparation in an arguments block. An option that is not given takes its default value from nansen.twophoton.roisignals.extract.getDefaultParameters, which stays the only place where the default values are defined. - roiInd="all" is resolved to all rois. It was used as a character index before, so the function only worked when the caller had resolved it. - An imageMask whose size differs from the roi masks throws NANSEN:RoiSignals:ImageMaskSizeMismatch. - createMasks passes the five options by name. - The help text describes how a mask is grown and what each option does. EstimateNeuropilMasksTest covers the options. Two of its tests fail if the | in pixelsToIgnore is changed back to &. Co-Authored-By: Claude Fable 5.1 --- code/+nansen/+processing/+roi/createMasks.m | 7 +- .../+processing/+roi/estimateNeuropilMasks.m | 133 +++++++++++++-- .../+roisignals/EstimateNeuropilMasksTest.m | 157 ++++++++++++++++++ 3 files changed, 279 insertions(+), 18 deletions(-) create mode 100644 tests/+nansen/+unittest/+roisignals/EstimateNeuropilMasksTest.m diff --git a/code/+nansen/+processing/+roi/createMasks.m b/code/+nansen/+processing/+roi/createMasks.m index 8b50e4fa0..7be319e74 100644 --- a/code/+nansen/+processing/+roi/createMasks.m +++ b/code/+nansen/+processing/+roi/createMasks.m @@ -50,7 +50,12 @@ return end - neuropilMaskArray = estimateNeuropilMasks(roiMaskArrayOrig, params); + neuropilMaskArray = estimateNeuropilMasks(roiMaskArrayOrig, ... + roiInd=params.roiInd, ... + imageMask=params.imageMask, ... + excludeRoiFromNeuropil=params.excludeRoiFromNeuropil, ... + neuropilExpansionFactor=params.neuropilExpansionFactor, ... + cellNeuropilSeparation=params.cellNeuropilSeparation); if params.numNeuropilSlices > 1 neuropilMasks = splitNeuropilMasks(neuropilMaskArray, ... diff --git a/code/+nansen/+processing/+roi/estimateNeuropilMasks.m b/code/+nansen/+processing/+roi/estimateNeuropilMasks.m index e2591e6d5..a8858191f 100644 --- a/code/+nansen/+processing/+roi/estimateNeuropilMasks.m +++ b/code/+nansen/+processing/+roi/estimateNeuropilMasks.m @@ -1,36 +1,105 @@ -function neuropilMasks = estimateNeuropilMasks(roiMasksIn, varargin) -%estimateNeuropilMasks Estimate neuropil masks for RoIs. +function neuropilMasks = estimateNeuropilMasks(roiMasks, options) +%estimateNeuropilMasks - Create a mask of the neuropil around each roi +% neuropilMasks = estimateNeuropilMasks(roiMasks) creates one neuropil +% mask for each roi mask. roiMasks is a logical array of size +% imageHeight-by-imageWidth-by-numRois and neuropilMasks is a logical +% array of the same size. A neuropil mask is grown from its roi mask in +% steps of one pixel. The first step adds the four edge neighbors of +% each mask pixel, the second step adds the four corner neighbors, and +% the following steps alternate between the two. Growing stops when +% the mask holds at least neuropilExpansionFactor times as many pixels +% as the roi, or after 101 steps. The pixels of a roi are never part of +% its own neuropil mask. Pixels that imageMask or excludeRoiFromNeuropil +% leave out of the neuropil masks are not counted. % +% neuropilMasks = estimateNeuropilMasks(roiMasks,Name=VALUE) specifies +% options using one or more name-value arguments. An option that is not +% given takes its default value from +% nansen.twophoton.roisignals.extract.getDefaultParameters: +% roiInd - Indices of the rois to create neuropil +% masks for, or "all" (default). The +% third dimension of neuropilMasks +% follows the order of roiInd. +% imageMask - Logical matrix of size +% imageHeight-by-imageWidth. Pixels that +% are false are left out of all neuropil +% masks. The default is empty, which +% leaves out no pixels. +% excludeRoiFromNeuropil - True (default) leaves the pixels of +% all rois in roiMasks out of all +% neuropil masks. +% neuropilExpansionFactor - Number of pixels in a neuropil mask +% relative to the number of pixels in its +% roi. The default is 4. +% cellNeuropilSeparation - Number of one-pixel steps that the +% rois are grown by before their pixels +% are left out of the neuropil masks. It +% has no effect if excludeRoiFromNeuropil +% is false. The default is 0. +% +% Example: Neuropil mask around one of two neighboring rois +% roiMasks = false(64,64,2); +% roiMasks(20:30,20:30,1) = true; +% roiMasks(20:30,34:44,2) = true; +% neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... +% roiMasks,roiInd=1); +% imagesc(neuropilMasks + 2*any(roiMasks,3)) +% +% See also nansen.processing.roi.createMasks, +% nansen.processing.roi.prepareRoiMasks, +% nansen.twophoton.roisignals.extractF, +% nansen.twophoton.roisignals.extract.getDefaultParameters + % TODO - % [ ] Instead of using max count of 100 to stop, use a parameter - % called max neuropil radius.... - % [ ] Add doc - % [ ] test the max radius - % [ ] add warning if mask can't be created (i.e not enough neuropil available) + % [ ] Stop growing at the maxNeuropilRadius parameter instead of + % after a fixed number of steps. + % [ ] Add warning if mask can't be created (i.e not enough neuropil + % available) + + arguments + roiMasks logical + options.roiInd {mustBeRoiIndices} + options.imageMask (:,:) logical + options.excludeRoiFromNeuropil (1,1) logical + options.neuropilExpansionFactor (1,1) {mustBeNumeric, mustBeReal, mustBePositive} + options.cellNeuropilSeparation (1,1) {mustBeNumeric, mustBeReal, mustBeInteger, mustBeNonnegative} + end + + [height, width, numRois] = size(roiMasks); - % Get default parameters and assertion functions. - [P, V] = nansen.twophoton.roisignals.extract.getDefaultParameters(); + params = addDefaultValues(options); - % Parse potential parameters from input arguments - params = utility.parsenvpairs(P, V, varargin{:}); + if isIndexOfAllRois(params.roiInd) + params.roiInd = 1:numRois; + end + + hasImageMask = ~isempty(params.imageMask); + if hasImageMask && ~isequal(size(params.imageMask), [height, width]) + error("NANSEN:RoiSignals:ImageMaskSizeMismatch", ... + "The size of imageMask is %s, but the size of the roi masks is %s. " + ... + "Give an imageMask with one element for each pixel of a roi mask.", ... + mat2str(size(params.imageMask)), mat2str([height, width])) + end % Create mask of pixel to ignore when creating neuropil masks if params.excludeRoiFromNeuropil % Create an array with all rois overlapping. - footprintAllRois = sum(roiMasksIn, 3) >= 1; + footprintAllRois = sum(roiMasks, 3) >= 1; - % Use binary erosion to remove pixels in the vicinity of main rois + % Use binary dilation to add pixels in the vicinity of main rois for i = 1:params.cellNeuropilSeparation footprintAllRois = imdilateoddeven(footprintAllRois, i); end - if ~isempty(params.imageMask) + % Ignore a pixel if it belongs to a roi or if it is excluded by + % the image mask + if hasImageMask pixelsToIgnore = footprintAllRois | ~params.imageMask; else pixelsToIgnore = footprintAllRois; end else - if ~isempty(params.imageMask) + if hasImageMask pixelsToIgnore = ~params.imageMask; else pixelsToIgnore = []; @@ -38,14 +107,13 @@ end % Allocate array for neuropil masks. - [height, width, ~] = size(roiMasksIn); neuropilMasks = false( height, width, numel(params.roiInd) ); % Loop through rois for n = 1:numel(params.roiInd) currentRoiInd = params.roiInd(n); - currentRoiMask = roiMasksIn(:, :, currentRoiInd); + currentRoiMask = roiMasks(:, :, currentRoiInd); thisMask = currentRoiMask; origArea = sum(currentRoiMask(:)); % original area @@ -93,6 +161,37 @@ end end +function params = addDefaultValues(options) +%addDefaultValues - Add the default value of each option that is not given + defaults = nansen.twophoton.roisignals.extract.getDefaultParameters(); + + optionNames = ["roiInd", "imageMask", "excludeRoiFromNeuropil", ... + "neuropilExpansionFactor", "cellNeuropilSeparation"]; + + params = options; + for name = optionNames + if ~isfield(params, name) + params.(name) = defaults.(name); + end + end +end + +function tf = isIndexOfAllRois(roiInd) +%isIndexOfAllRois - Check if roiInd is the text "all" + tf = (ischar(roiInd) || isstring(roiInd)) && isequal(string(roiInd), "all"); +end + +function mustBeRoiIndices(value) +%mustBeRoiIndices - Validate that value is "all" or a vector of indices + isIndexVector = isnumeric(value) && isvector(value) && isreal(value) ... + && all(value == round(value)) && all(value >= 1); + + if ~(isIndexOfAllRois(value) || isIndexVector) + error("NANSEN:RoiSignals:InvalidRoiIndices", ... + "Value must be ""all"" or a vector of positive integers.") + end +end + function mask = imdilateoddeven(mask, count) %imdilateoddeven imdilate using different nhood on odd and even iterations % diff --git a/tests/+nansen/+unittest/+roisignals/EstimateNeuropilMasksTest.m b/tests/+nansen/+unittest/+roisignals/EstimateNeuropilMasksTest.m new file mode 100644 index 000000000..da67f516b --- /dev/null +++ b/tests/+nansen/+unittest/+roisignals/EstimateNeuropilMasksTest.m @@ -0,0 +1,157 @@ +classdef EstimateNeuropilMasksTest < matlab.unittest.TestCase + %EstimateNeuropilMasksTest Unit tests for estimateNeuropilMasks. + % + % Tests nansen.processing.roi.estimateNeuropilMasks on three square + % roi masks in a 64-by-64 image. Roi 1 and roi 2 are 3 pixels apart, + % so that roi 2 lies where the neuropil mask of roi 1 grows. Roi 3 + % has no neighbors. + % + % testMasksLeaveOutRoiPixels with the allTrue parameter is a + % regression test: an image mask that includes all pixels used to + % switch off excludeRoiFromNeuropil, and extractF always passes such + % a mask. + % + % Run tests: + % runtests('nansen.unittest.roisignals.EstimateNeuropilMasksTest') + + properties (Constant, Access = private) + ImageSize = [64, 64] + end + + properties (TestParameter) + % Image masks that leave out no pixels + fullImageMask = struct( ... + 'empty', {logical.empty}, ... + 'allTrue', {true(64, 64)}) + end + + methods (Test) + + function testMasksLeaveOutRoiPixels(testCase, fullImageMask) + roiMasks = createRoiMasks(testCase.ImageSize); + + neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, imageMask=fullImageMask); + + isRoiPixel = any(roiMasks, 3); + testCase.verifyFalse(any(neuropilMasks & isRoiPixel, "all")) + end + + function testNeighborRoiIsKeptWhenRoisAreNotExcluded(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, excludeRoiFromNeuropil=false); + + maskOfRoi1 = neuropilMasks(:, :, 1); + testCase.verifyTrue(any(maskOfRoi1 & roiMasks(:, :, 2), "all"), ... + "The neuropil mask of roi 1 should reach into roi 2.") + testCase.verifyFalse(any(maskOfRoi1 & roiMasks(:, :, 1), "all"), ... + "A neuropil mask should never hold pixels of its own roi.") + end + + function testMasksLeaveOutPixelsOutsideImageMask(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + imageMask = true(testCase.ImageSize); + imageMask(:, 1:19) = false; % The columns to the left of roi 1 + + neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, imageMask=imageMask); + + testCase.verifyFalse(any(neuropilMasks & ~imageMask, "all")) + testCase.verifyTrue(any(neuropilMasks(:, :, 1), "all")) + end + + function testImageMaskAppliesWhenRoisAreNotExcluded(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + imageMask = true(testCase.ImageSize); + imageMask(:, 1:19) = false; + + neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, imageMask=imageMask, excludeRoiFromNeuropil=false); + + testCase.verifyFalse(any(neuropilMasks & ~imageMask, "all")) + end + + function testMaskAreaReachesExpansionFactor(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + expansionFactor = 2; + + neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, neuropilExpansionFactor=expansionFactor); + + roiArea = squeeze(sum(roiMasks, [1, 2])); + neuropilArea = squeeze(sum(neuropilMasks, [1, 2])); + testCase.verifyGreaterThanOrEqual(neuropilArea, expansionFactor * roiArea) + end + + function testLargerExpansionFactorGivesLargerMasks(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + smallMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, neuropilExpansionFactor=2); + largeMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, neuropilExpansionFactor=4); + + smallArea = squeeze(sum(smallMasks, [1, 2])); + largeArea = squeeze(sum(largeMasks, [1, 2])); + testCase.verifyGreaterThan(largeArea, smallArea) + end + + function testSeparationKeepsMasksAwayFromRois(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + neuropilMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, cellNeuropilSeparation=2); + + % Two steps cover the eight neighbors of every roi pixel + isNextToRoi = imdilate(any(roiMasks, 3), ones(3)); + testCase.verifyFalse(any(neuropilMasks & isNextToRoi, "all")) + end + + function testRoiIndSelectsAndOrdersMasks(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + allMasks = nansen.processing.roi.estimateNeuropilMasks(roiMasks); + selectedMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, roiInd=[3, 1]); + + testCase.verifyEqual(selectedMasks, allMasks(:, :, [3, 1])) + end + + function testRoiIndAllEqualsDefault(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + defaultMasks = nansen.processing.roi.estimateNeuropilMasks(roiMasks); + allMasks = nansen.processing.roi.estimateNeuropilMasks( ... + roiMasks, roiInd="all"); + + testCase.verifySize(allMasks, size(roiMasks)) + testCase.verifyEqual(allMasks, defaultMasks) + end + + function testImageMaskOfWrongSizeThrowsError(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + testCase.verifyError( ... + @() nansen.processing.roi.estimateNeuropilMasks(roiMasks, imageMask=true(10, 10)), ... + "NANSEN:RoiSignals:ImageMaskSizeMismatch") + end + + function testInvalidRoiIndThrowsError(testCase) + roiMasks = createRoiMasks(testCase.ImageSize); + + testCase.verifyError( ... + @() nansen.processing.roi.estimateNeuropilMasks(roiMasks, roiInd=0), ... + "NANSEN:RoiSignals:InvalidRoiIndices") + end + end +end + +function roiMasks = createRoiMasks(imageSize) +%createRoiMasks Two neighboring rois and one roi without neighbors + roiMasks = false([imageSize, 3]); + roiMasks(20:30, 20:30, 1) = true; + roiMasks(20:30, 34:44, 2) = true; + roiMasks(45:50, 8:13, 3) = true; +end From e93f050986a3ff65ceb1ed05c727dcab166b9da1 Mon Sep 17 00:00:00 2001 From: ehennestad Date: Mon, 21 Sep 2026 23:18:48 +0200 Subject: [PATCH 4/5] refactor(roisignals): validate extractF options in an arguments block extractF parsed varargin with utility.parsenvpairs, which accepts a struct, name-value pairs with char names, or both in either order, and silently ignores names it does not know. Names given as strings, and therefore all options given with Name=Value syntax, were also ignored. extractF now takes an optional options struct followed by name-value arguments that are declared in an arguments block: - A misspelled option name throws an error. - Options in the struct are passed on as name-value arguments, so they are validated in the same way. Fields that are not options are ignored, because callers pass the struct from getDefaultParameters, which also holds the fields that configure the options editor. - An option given by name replaces the field of the same name in the struct, as before. - No option has a default value in the arguments block. Default values are added afterwards from getDefaultParameters, because a default value in the block could not be told apart from a value given by name and would replace the value in the struct. The struct has to come before the name-value arguments. The extractSignals session method passed it last and is updated. All other callers already pass the struct first. The help text documents the size and layout of signalArray, every option, the second output and the behavior for stacks with more than one current channel. The unused variable method is removed. ExtractFTest compares the signals with the mean pixel values of hand-made masks and covers the calling conventions and the second output. Co-Authored-By: Claude Fable 5.1 --- .../+nansen/+twophoton/+roisignals/extractF.m | 241 +++++++++++++++--- .../+signalExtraction/extractSignals.m | 2 +- .../+unittest/+roisignals/ExtractFTest.m | 227 +++++++++++++++++ 3 files changed, 430 insertions(+), 40 deletions(-) create mode 100644 tests/+nansen/+unittest/+roisignals/ExtractFTest.m diff --git a/code/+nansen/+twophoton/+roisignals/extractF.m b/code/+nansen/+twophoton/+roisignals/extractF.m index 11f9a0b3b..b37a7289a 100644 --- a/code/+nansen/+twophoton/+roisignals/extractF.m +++ b/code/+nansen/+twophoton/+roisignals/extractF.m @@ -1,33 +1,105 @@ -function [signalArray, params] = extractF(imageData, roiData, varargin) -%extractF Extract RoI signals from an image stack in chunkwise manner +function [signalArray, params] = extractF(imageData, roiData, optionsStruct, options) +%extractF - Extract roi signals from an image stack +% signalArray = extractF(imageData,roiData) computes the mean pixel +% value of each roi in roiData for every frame of imageData. imageData +% is a nansen.stack.ImageStack or a numeric array of size +% imageHeight-by-imageWidth-by-numFrames. roiData is an array of RoI +% objects or a roimanager.roiGroup, and the imagesize property of the +% rois must equal [imageHeight, imageWidth]. signalArray is of size +% numFrames-by-numSubregions-by-numRois. signalArray(:,1,:) holds the +% signals of the rois and signalArray(:,2:end,:) holds the signals of +% the neuropil masks that are created around each roi. Frames are read +% and processed in chunks if imageData is a virtual stack or if there +% are 100 rois or more. The chooseChunkLength method of +% nansen.stack.ImageStack sets the number of frames in a chunk. % -% signals = extractF(imageStack, roiArray) extracts signals for rois in -% roiArray from image frames in imageData using default settings. -% roiArray must be a nansen.RoiArray object and imageData can be a -% numeric array of three or more dimensions or a nansen.ImageStack -% object. +% signalArray = extractF(imageData,roiData,optionsStruct) specifies +% options as fields of the scalar struct optionsStruct. The field names +% are the names of the name-value arguments below. Fields with other +% names are ignored, so that optionsStruct can be the struct returned +% by nansen.twophoton.roisignals.extract.getDefaultParameters. % -% signals = extractF(imageStack, roiArray, params) extracts signals -% using options specified in params. Params gan be given as a -% struct or as a cell array of name-value pairs. +% signalArray = extractF(...,Name=VALUE) specifies options using one or +% more name-value arguments. An option given by name replaces the field +% of the same name in optionsStruct. An option that is not given takes +% its default value from +% nansen.twophoton.roisignals.extract.getDefaultParameters: +% roiInd - Indices of the rois to extract signals +% for, or "all" (default). +% imageMask - Logical matrix of size +% imageHeight-by-imageWidth. Pixels that +% are false are left out of the neuropil +% masks. The roi masks are not affected. +% The default includes all pixels. +% pixelComputationMethod - Statistic that is computed over the +% pixels of a mask, "mean" (default) or +% "median". "median" is always computed +% by serialExtract. +% excludeRoiOverlaps - True (default) leaves pixels that +% belong to more than one roi out of the +% roi masks. +% createNeuropilMask - True (default) creates a neuropil mask +% around each roi and extracts its +% signal. +% excludeRoiFromNeuropil - True (default) leaves the pixels of +% all rois out of the neuropil masks. +% neuropilExpansionFactor - Number of pixels in a neuropil mask +% relative to the number of pixels in +% its roi. The default is 4. +% cellNeuropilSeparation - Number of one-pixel steps that the +% rois are grown by before their pixels +% are left out of the neuropil masks. +% The default is 0. +% numNeuropilSlices - Number of slices that each neuropil +% mask is split into. Each slice is one +% subregion in signalArray. The default +% is 1. +% maxNeuropilRadius - Not used. +% extractFcn - Handle of the function that extracts +% the signals from a chunk of frames. +% The default is serialExtract for fewer +% than 200 rois and batchExtract for 200 +% rois or more. Both are in the package +% nansen.twophoton.roisignals.extract. +% roiMaskFormat - Format of the prepared roi masks, +% "struct" for serialExtract and +% "sparse" for batchExtract. It is set +% to match extractFcn. +% signalDataType - Class of signalArray, "single" +% (default) or "double". +% channelNum - Channels to extract signals from if +% imageData has more than one current +% channel. The signals of the channels +% are concatenated along the third +% dimension of signalArray. If roiData +% is an array of roimanager.roiGroup +% objects, roiData(i) is used for +% channel i. The default is 1. +% verbose - True prints a message for each chunk +% of frames and the elapsed time. The +% default is false. +% showTimer - True prints the elapsed time. The +% default is false. % -% Parameters: +% [signalArray,PARAMS] = extractF(...) also returns the options that +% were used, as a scalar struct with one field for each name-value +% argument. PARAMS.roiInd holds the indices of the rois, +% PARAMS.imageMask the image mask that was applied, and +% PARAMS.extractFcn and PARAMS.roiMaskFormat the extraction function +% and the mask format that were selected. If imageData has more than +% one current channel, these four fields hold the values that were +% given, because they are resolved separately for each channel. % -% roiInd : integer vector in range [1, numRois] -% List of roi indices (use for extraction of signals from a subset -% of rois). Default : extract signals from all rois. -% imageMask : logical matrix (imageHeight x imageWidth) -% Use for excluding regions of image. Include pixels that are true -% and exclude pixels that are false. Default : include all pixels -% excludeRoiOverlaps : logical scalar -% Exclude pixels where rois are overlapping. Default = true -% createNeuropilMask : logical scalar -% Create mask (and extract signals) for surrounding neuropil regions -% Default = true -% excludeRoiFromNeuropil : logical scalar -% Exclude rois from neuropil regions. Default = true +% Example: Extract the signal of one roi from a numeric array +% imageData = rand(64,64,100); +% roiArray = RoI("Circle",[32,32,5],[64,64]); +% signalArray = nansen.twophoton.roisignals.extractF( ... +% imageData,roiArray,createNeuropilMask=false); +% plot(signalArray) % -% See also nansen.twophoton.roisignals.extract.getDefaultParameters +% See also nansen.twophoton.roisignals.computeDff, +% nansen.twophoton.roisignals.extract.getDefaultParameters, +% nansen.processing.roi.prepareRoiMasks, nansen.stack.ImageStack, RoI % Todo: % [x] Implement a standard set of options across signal extraction @@ -37,9 +109,9 @@ % versions of signals???? % % [x] Make sure it works to extract from subset of rois.... - % [ ] Test that it works as expected. I.e do we get signals from - % requested rois, and are masks created according to options. Use - % imviewer/roimanager gui for this... + % [x] Test that it works as expected. I.e do we get signals from + % requested rois, and are masks created according to options. See + % nansen.unittest.roisignals.ExtractFTest % % [ ] Extraction function should depend on number of rois. If nRois < % 140 should use looping over rois and array cropping, otherwise @@ -50,16 +122,44 @@ % % [ ] Support for multiple planes. - % Get default parameters and assertion functions. + arguments + imageData + roiData + optionsStruct struct {mustBeScalarOrEmpty} = struct.empty + options.roiInd {mustBeRoiIndices} + options.imageMask (:,:) logical + options.pixelComputationMethod {mustBeTextScalar, ... + mustBeMember(options.pixelComputationMethod, ["mean", "median"])} + options.excludeRoiOverlaps (1,1) logical + options.createNeuropilMask (1,1) logical + options.excludeRoiFromNeuropil (1,1) logical + options.neuropilExpansionFactor (1,1) {mustBeNumeric, mustBeReal, mustBePositive} + options.cellNeuropilSeparation (1,1) {mustBeNumeric, mustBeReal, mustBeInteger, mustBeNonnegative} + options.numNeuropilSlices (1,1) {mustBeNumeric, mustBeReal, mustBeInteger, mustBeNonnegative} + options.maxNeuropilRadius (1,1) {mustBeNumeric, mustBeReal, mustBeInteger, mustBeNonnegative} + options.roiMaskFormat {mustBeTextScalar, ... + mustBeMember(options.roiMaskFormat, ["sparse", "struct"])} + options.extractFcn {mustBeFunctionHandleOrEmpty} + options.showTimer (1,1) logical + options.verbose (1,1) logical + options.signalDataType {mustBeTextScalar, ... + mustBeMember(options.signalDataType, ["single", "double"])} + options.channelNum {mustBeNumeric, mustBeReal, mustBeInteger, mustBePositive, mustBeNonempty} + end - [P, V] = nansen.twophoton.roisignals.extract.getDefaultParameters(); - P.showTimer = false; V.showTimer = @(x) assert(islogical(x), 'Value must be logical'); - P.verbose = false; V.verbose = @(x) assert(islogical(x), 'Value must be logical'); - P.signalDataType = 'single'; V.signalDataType = @(x) assert(any(strcmp(x, {'single', 'double'})), 'Value must be ''single'' or ''double'''); - P.channelNum = 1; V.channelNum = @(x) not(isempty(x)); + % Options in a struct are passed on as name-value arguments, so that + % they are validated in the same way as options given by name. No + % option has a default value in the arguments block, because a default + % value could not be told apart from a value given by name, and would + % replace the value in the struct. + if ~isempty(optionsStruct) + nvPairs = combineOptions(optionsStruct, options); % Local function + [signalArray, params] = nansen.twophoton.roisignals.extractF( ... + imageData, roiData, nvPairs{:}); + return + end - % Parse potential parameters from input arguments - params = utility.parsenvpairs(P, V, varargin{:}); + params = addDefaultValues(options); % Local function % Validate roidata if isa(roiData, 'roimanager.roiGroup') @@ -79,13 +179,14 @@ % If multiple channels are selected, extract from each individually numChannels = numel( imageStack.CurrentChannel ); if numChannels > 1 + nvPairs = namedargs2cell(options); if ~isempty(params.channelNum) currentChannel = imageStack.CurrentChannel; imageStack.CurrentChannel = params.channelNum; - signalArray = extractFromMultipleChannels(imageStack, roiData, varargin{:}); + signalArray = extractFromMultipleChannels(imageStack, roiData, nvPairs{:}); imageStack.CurrentChannel = currentChannel; else - signalArray = extractFromMultipleChannels(imageStack, roiData, varargin{:}); + signalArray = extractFromMultipleChannels(imageStack, roiData, nvPairs{:}); end return end @@ -121,7 +222,6 @@ [elapsedTimeLoad, elapsedTimeExtract] = deal( 0 ); signalExtractionFcn = params.extractFcn; - method = params.pixelComputationMethod; if params.verbose fprintf('Image stack is split in %d parts for signal extraction\n', numParts) @@ -151,6 +251,69 @@ end end +function params = addDefaultValues(options) +%addDefaultValues - Add the default value of each option that is not given + params = nansen.twophoton.roisignals.extract.getDefaultParameters(); + + % Remove the fields that configure the options editor. The name of + % such a field is the name of a parameter followed by an underscore. + names = string(fieldnames(params)); + isConfigField = endsWith(names, "_") & ... + ismember(extractBefore(names, strlength(names)), names); + params = rmfield(params, names(isConfigField)); + + % Options of extractF that are not parameters of the extraction + % functions + params.showTimer = false; + params.verbose = false; + params.signalDataType = "single"; + params.channelNum = 1; + + names = fieldnames(options); + for i = 1:numel(names) + params.(names{i}) = options.(names{i}); + end +end + +function nvPairs = combineOptions(optionsStruct, options) +%combineOptions - Combine an options struct with name-value arguments +% Fields of optionsStruct that are not options of extractF are left out. +% A field of options replaces the field of the same name in +% optionsStruct. + optionNames = fieldnames(addDefaultValues(struct())); + + fieldNames = fieldnames(optionsStruct); + isOption = ismember(fieldNames, optionNames); + combined = rmfield(optionsStruct, fieldNames(~isOption)); + + names = fieldnames(options); + for i = 1:numel(names) + combined.(names{i}) = options.(names{i}); + end + + nvPairs = namedargs2cell(combined); +end + +function mustBeRoiIndices(value) +%mustBeRoiIndices - Validate that value is "all" or a vector of indices + isAll = (ischar(value) || isstring(value)) && isequal(string(value), "all"); + isIndexVector = isnumeric(value) && isvector(value) && isreal(value) ... + && all(value == round(value)) && all(value >= 1); + + if ~(isAll || isIndexVector) + error("NANSEN:RoiSignals:InvalidRoiIndices", ... + "Value must be ""all"" or a vector of positive integers.") + end +end + +function mustBeFunctionHandleOrEmpty(value) +%mustBeFunctionHandleOrEmpty - Validate that value is a function handle or empty + if ~(isempty(value) || isa(value, "function_handle")) + error("NANSEN:RoiSignals:InvalidExtractFcn", ... + "Value must be a function handle or empty.") + end +end + function validateInputDimensions(imageStack, roiArray) %validateInputDimensions Check that ImageStack and roiArray has matching % dimensions. diff --git a/code/modules/+nansen/+module/+ophys/+twophoton/+sessionmethod/+process/+signalExtraction/extractSignals.m b/code/modules/+nansen/+module/+ophys/+twophoton/+sessionmethod/+process/+signalExtraction/extractSignals.m index 6b9b7e6df..df6002255 100644 --- a/code/modules/+nansen/+module/+ophys/+twophoton/+sessionmethod/+process/+signalExtraction/extractSignals.m +++ b/code/modules/+nansen/+module/+ophys/+twophoton/+sessionmethod/+process/+signalExtraction/extractSignals.m @@ -67,7 +67,7 @@ function runMethod(obj) roiArray = sessionData.RoiArray; extractF = @nansen.twophoton.roisignals.extractF; - [signalArray, P] = extractF(imageStack, roiArray, 'verbose', true, obj.Options); + [signalArray, P] = extractF(imageStack, roiArray, obj.Options, verbose=true); % Todo: Save results... obj.saveData('RoiSignals_MeanF', squeeze(signalArray(:, 1, :)) ) diff --git a/tests/+nansen/+unittest/+roisignals/ExtractFTest.m b/tests/+nansen/+unittest/+roisignals/ExtractFTest.m new file mode 100644 index 000000000..ef54c3783 --- /dev/null +++ b/tests/+nansen/+unittest/+roisignals/ExtractFTest.m @@ -0,0 +1,227 @@ +classdef ExtractFTest < matlab.unittest.TestCase + %ExtractFTest Unit tests for nansen.twophoton.roisignals.extractF. + % + % The extracted signals are compared with the mean pixel values of + % hand-made masks. Three square rois lie in a 48-by-64 image: roi 1 + % and roi 2 overlap in three columns and roi 3 has no neighbors. + % + % Two tests are regression tests. In + % testNeuropilSignalLeavesOutNeighboringRois the rois are brighter + % than the background, so the neuropil signal only equals the + % background if no neuropil mask holds a pixel of any roi. The + % testSecondOutput* tests check that the second output holds the + % options that were used and not the default options. + % + % Run tests: + % runtests('nansen.unittest.roisignals.ExtractFTest') + + properties (Constant, Access = private) + ImageSize = [48, 64] + NumFrames = 25 + AbsTol = 1e-9 % Signals are requested as double + end + + methods (TestMethodSetup) + function seedRandomNumberGenerator(testCase) + originalState = rng(); + testCase.addTeardown(@() rng(originalState)) + rng(1, "twister") + end + end + + methods (Test) + + function testRoiSignalIsMeanOfNonOverlappingPixels(testCase) + [rois, roiMasks] = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + isInOneRoi = sum(roiMasks, 3) == 1; + + signalArray = nansen.twophoton.roisignals.extractF( ... + imageData, rois, signalDataType="double"); + + expected = meanInMasks(imageData, roiMasks & isInOneRoi); + testCase.verifySize(signalArray, [testCase.NumFrames, 2, 3]) + testCase.verifyEqual(squeeze(signalArray(:, 1, :)), expected, ... + AbsTol=testCase.AbsTol) + end + + function testRoiSignalUsesOverlappingPixelsOnRequest(testCase) + [rois, roiMasks] = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + + signalArray = nansen.twophoton.roisignals.extractF( ... + imageData, rois, excludeRoiOverlaps=false, signalDataType="double"); + + expected = meanInMasks(imageData, roiMasks); + testCase.verifyEqual(squeeze(signalArray(:, 1, :)), expected, ... + AbsTol=testCase.AbsTol) + end + + function testNeuropilSignalLeavesOutNeighboringRois(testCase) + [rois, roiMasks] = createRois(testCase.ImageSize); + backgroundValue = 1000; + imageData = createTwoLevelStack(any(roiMasks, 3), ... + testCase.NumFrames, backgroundValue, 6000); + + signalArray = nansen.twophoton.roisignals.extractF( ... + imageData, rois, signalDataType="double"); + + expected = backgroundValue * ones(testCase.NumFrames, 3); + testCase.verifyEqual(squeeze(signalArray(:, 2, :)), expected, ... + AbsTol=testCase.AbsTol) + end + + function testImageMaskLeavesPixelsOutOfNeuropilSignal(testCase) + rois = createRois(testCase.ImageSize); + imageMask = true(testCase.ImageSize); + imageMask(:, 1:8) = false; % Columns to the left of roi 1 + backgroundValue = 1000; + imageData = createTwoLevelStack(~imageMask, ... + testCase.NumFrames, backgroundValue, 9000); + + signalArray = nansen.twophoton.roisignals.extractF( ... + imageData, rois, imageMask=imageMask, signalDataType="double"); + + expected = backgroundValue * ones(testCase.NumFrames, 3); + testCase.verifyEqual(squeeze(signalArray(:, 2, :)), expected, ... + AbsTol=testCase.AbsTol) + end + + function testBatchExtractEqualsSerialExtract(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + + serialSignals = nansen.twophoton.roisignals.extractF(imageData, rois, ... + extractFcn=@nansen.twophoton.roisignals.extract.serialExtract, ... + roiMaskFormat="struct", signalDataType="double"); + batchSignals = nansen.twophoton.roisignals.extractF(imageData, rois, ... + extractFcn=@nansen.twophoton.roisignals.extract.batchExtract, ... + roiMaskFormat="sparse", signalDataType="double"); + + testCase.verifyEqual(batchSignals, serialSignals, AbsTol=testCase.AbsTol) + end + + function testSubsetOfRoisEqualsColumnsOfFullExtraction(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + + allSignals = nansen.twophoton.roisignals.extractF(imageData, rois); + selectedSignals = nansen.twophoton.roisignals.extractF( ... + imageData, rois, roiInd=[3, 1]); + + testCase.verifyEqual(selectedSignals, allSignals(:, :, [3, 1])) + end + + function testSecondOutputHoldsOptionsGivenByName(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + + [~, params] = nansen.twophoton.roisignals.extractF( ... + imageData, rois, createNeuropilMask=false, roiInd=[3, 1]); + + testCase.verifyFalse(params.createNeuropilMask) + testCase.verifyEqual(params.roiInd, [3, 1]) + end + + function testSecondOutputHoldsResolvedDefaults(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + + [~, params] = nansen.twophoton.roisignals.extractF(imageData, rois); + + testCase.verifyEqual(params.roiInd, 1:3) + testCase.verifyEqual(params.imageMask, true(testCase.ImageSize)) + testCase.verifyEqual(params.extractFcn, ... + @nansen.twophoton.roisignals.extract.serialExtract) + testCase.verifyTrue(params.createNeuropilMask) + end + + function testOptionsStructEqualsNameValueArguments(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + optionsStruct = struct('createNeuropilMask', false, 'roiInd', 2); + + [signalsFromStruct, paramsFromStruct] = ... + nansen.twophoton.roisignals.extractF(imageData, rois, optionsStruct); + [signalsFromNames, paramsFromNames] = ... + nansen.twophoton.roisignals.extractF(imageData, rois, ... + createNeuropilMask=false, roiInd=2); + + testCase.verifySize(signalsFromStruct, [testCase.NumFrames, 1]) + testCase.verifyEqual(signalsFromStruct, signalsFromNames) + testCase.verifyEqual(paramsFromStruct, paramsFromNames) + end + + function testNameValueArgumentReplacesStructField(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + optionsStruct = struct('createNeuropilMask', true); + + signalArray = nansen.twophoton.roisignals.extractF( ... + imageData, rois, optionsStruct, createNeuropilMask=false); + + testCase.verifySize(signalArray, [testCase.NumFrames, 1, 3]) + end + + function testDefaultParameterStructIsAccepted(testCase) + % The struct holds fields such as roiInd_ that configure the + % options editor and are not options of extractF. + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + optionsStruct = nansen.twophoton.roisignals.extract.getDefaultParameters(); + + signalsFromStruct = nansen.twophoton.roisignals.extractF( ... + imageData, rois, optionsStruct); + signalsFromDefaults = nansen.twophoton.roisignals.extractF(imageData, rois); + + testCase.verifyEqual(signalsFromStruct, signalsFromDefaults) + end + + function testInvalidValueInOptionsStructThrowsError(testCase) + rois = createRois(testCase.ImageSize); + imageData = createRandomStack(testCase.ImageSize, testCase.NumFrames); + optionsStruct = struct('pixelComputationMethod', 'max'); + + testCase.verifyError( ... + @() nansen.twophoton.roisignals.extractF(imageData, rois, optionsStruct), ... + "MATLAB:validators:mustBeMember") + end + end +end + +function [rois, roiMasks] = createRois(imageSize) +%createRois Two overlapping rois and one roi without neighbors + roiMasks = false([imageSize, 3]); + roiMasks(10:18, 10:18, 1) = true; + roiMasks(10:18, 16:24, 2) = true; + roiMasks(30:38, 40:48, 3) = true; + + rois = [RoI('Mask', roiMasks(:, :, 1)), ... + RoI('Mask', roiMasks(:, :, 2)), ... + RoI('Mask', roiMasks(:, :, 3))]; +end + +function imageData = createRandomStack(imageSize, numFrames) +%createRandomStack Stack of uniformly distributed uint16 pixel values + imageData = uint16(randi([100, 4000], [imageSize, numFrames])); +end + +function imageData = createTwoLevelStack(isBright, numFrames, backgroundValue, brightValue) +%createTwoLevelStack Stack with one value in bright pixels and one elsewhere + frame = backgroundValue * ones(size(isBright), 'single'); + frame(isBright) = brightValue; + imageData = repmat(frame, 1, 1, numFrames); +end + +function signals = meanInMasks(imageData, masks) +%meanInMasks Mean pixel value per frame (rows) and mask (columns) + numFrames = size(imageData, 3); + numMasks = size(masks, 3); + pixelValues = reshape(double(imageData), [], numFrames); + + signals = zeros(numFrames, numMasks); + for i = 1:numMasks + mask = masks(:, :, i); + signals(:, i) = mean(pixelValues(mask(:), :), 1)'; + end +end From 23425dc7088a2ec9e8102c68e400eecba76e288b Mon Sep 17 00:00:00 2001 From: "github-actions[bot]" <41898282+github-actions[bot]@users.noreply.github.com> Date: Mon, 21 Sep 2026 21:45:27 +0000 Subject: [PATCH 5/5] Update GitHub badges [skip ci] --- .github/badges/code_issues.svg | 2 +- .github/badges/tests.svg | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/.github/badges/code_issues.svg b/.github/badges/code_issues.svg index edc4000eb..e0d2558a1 100644 --- a/.github/badges/code_issues.svg +++ b/.github/badges/code_issues.svg @@ -1 +1 @@ -code issues: 1290code issues1290 +code issues: 1289code issues1289 diff --git a/.github/badges/tests.svg b/.github/badges/tests.svg index fd00e9b16..783fa2a5e 100644 --- a/.github/badges/tests.svg +++ b/.github/badges/tests.svg @@ -1 +1 @@ -tests: 369 passedtests369 passed +tests: 393 passedtests393 passed