Find local maxima of each column in matrix

I have the columns of this 526x26 matrix normalized and plotted here:

n=hkcl;
for ii=[2 4 6 8 10 12 14 16 18 20 22 26]
    n(:,ii)=(hkcl(:,ii)-min(hkcl(:,ii)))/(max(hkcl(:,ii))-min(hkcl(:,ii)));
    plot(hkcl(:,1),n(:,ii))
end

Now I want to find the local maxima of the listed ii columns in matrix n in my case. It would be 2 data points per column. I tried some things with max() and localmax() but no luck.

[Walter Roberson]

islocalmax

[Mathieu NOE]

hello Claire

can you share your data ? what issue did you had ?

[Mathieu NOE]

hello @CLAIRE

problem solved ?

[Constantino Carlos Reyes-Aldasoro]

Do you want to find the 2 highest values per column (which could be next to each other and part of a single peak)? Or find peaks and then take the 2 highest ones?

Most probably findpeaks would find what you need:

https://uk.mathworks.com/help/signal/ref/findpeaks.html

[Ken Garrard]

If you do not have the signal processing toolbox consider using morphological image dialation to find minima and/or maxima in a signal or a surface. Here is my fndpeaks function which is also up to 1000 times faster than fndpeaks.

function [ex,p,v,q,u] = fndpeaks(x,L,M)
% FNDPEAKS  Find local extreme values in a vector or matrix
%   [ex,p,v,q,u] = fndpeaks(x,L,M)
% Input
%    x     vector or matrix of data
%    L,M   size of structuring element for finding local for maxima, minima
% Output
%    ex    matrix with local maxima (col 1) and local minima (col 2)
%    p,v   row    index vectors into x at local peak,valley locations
%    q,u   column index vectors into x at local peak,valley locations
%
% Morphological dilation, is used to find centered peaks and valleys on a
% surface with a rectangular two dimensional structuring element of size L x M.
% If L and M are omitted the global maximum and minimum are returned.
% If the number of peaks and valleys are not equal, the shorter column of ex is
% padded with NaNs.
% If x is a vector, then only x and L need to be given as arguments and only ex,
% p and v specified for return values (q,u will be column vectors of ones).
%
% Execute fndpeaks without input arguments to run these two examples.
% 3D surface
%   n = 201;
%   [xx,yy] = meshgrid(linspace(-10,10,n),linspace(-10,10,n));
%   zz      = pi*cos(xx) + sqrt(2)*sin(yy) - xx/2 - (yy.^2)/8;
%   [ex,p,v,q,u] = fndpeaks(zz,20,20);
%   figure; hold on; grid on;
%   mesh(yy(1:4:end,1:4:end),xx(1:4:end,1:4:end),zz(1:4:end,1:4:end));
%   plot3(xx(1,p),yy(q,1),ex(1:length(p),1),'bo',...
%         'MarkerSize',8,'MarkerFaceColor',[0 0 1 ]);
%   plot3(xx(1,v),yy(u,1),ex(1:length(v),2),'ro',...
%         'MarkerSize',8,'MarkerFaceColor',[1 0 0 ]);
%   set(gca,'Color',[.75 .75 .75]);
%   view([-46,52]);
%
% 2D vector
%   sig = diag(zz);
%   [ex,p,v] = fndpeaks(sig,20);
%   figure; hold on; grid on;
%   plot(sig,'y-','LineWidth',2); set(gca,'Color',[.75 .75 .75]);
%   plot(p,sig(p),'bo','MarkerSize',8,'MarkerFaceColor',[0 0 1]);
%   plot(v,sig(v),'ro','MarkerSize',8,'MarkerFaceColor',[1 0 0]);
if nargin < 1
    n = 201;
    [xx,yy] = meshgrid(linspace(-10,10,n),linspace(-10,10,n));
    zz      = pi*cos(xx) + sqrt(2)*sin(yy) - xx/2 - (yy.^2)/8;
    [ex,p,v,q,u] = fndpeaks(zz,20,20);
    figure; hold on; grid on;
    title(sprintf('%s example',mfilename));
    mesh(yy(1:4:end,1:4:end),xx(1:4:end,1:4:end),zz(1:4:end,1:4:end));
    plot3(xx(1,p),yy(q,1),ex(1:length(p),1),'bo',...
        'MarkerSize',8,'MarkerFaceColor',[0 0 1 ]);
    plot3(xx(1,v),yy(u,1),ex(1:length(v),2),'ro',...
        'MarkerSize',8,'MarkerFaceColor',[1 0 0 ]);
    set(gca,'Color',[.75 .75 .75]);
    view([-46,52]);
    sig = diag(zz);
    [ex,p,v] = fndpeaks(sig,20); %#ok<ASGLU>
    figure; hold on; grid on;
    title(sprintf('%s example',mfilename));
    plot(sig,'y-','LineWidth',2); set(gca,'Color',[.75 .75 .75]);
    plot(p,sig(p),'bo','MarkerSize',8,'MarkerFaceColor',[0 0 1]);
    plot(v,sig(v),'ro','MarkerSize',8,'MarkerFaceColor',[1 0 0]);
    % display help text
    error(['No input arguments given\n'...
        'Please consult the help text and the example plots\n'...
        '--------\n%s'],help(mfilename));
end
% no structure element specified, find global extrema and return
if nargin < 2 || (isempty(L) && isempty(M))
    [ex(1),p] = max(x(:));      % max value and index
    [ex(2),v] = min(x(:));      % min value and index
    [p,q] = ind2sub(size(x),p);
    [v,u] = ind2sub(size(x),v);
    return
end
% empty or missing structure element arguments
if nargin < 2 || isempty(L), L = 1; end      % one unit tall (rows)
if nargin < 3 || isempty(M), M = 1; end      % one unit wide (cols)
% check vector input direction
[r,c] = size(x);               % size of input data
if r==1                        % row vector ?
    x = x.';                    %   yes, transpose to column vector
    [L,M,r,c] = deal(M,L,c,r);  %        swap scalars
end
if c==1 && min(L,M)==1         % column vector ?
    [L,M] = deal(max(L,M),1);   %   yes, L is length of structure element
end
% create 1D structuring element for dilation operation
s = ones(L,M);
% find local maxima
d = imdilate(x,s);             % local maxima of data
p = find(x==d);                % find 1D indices of local maxima
% find local minima
y = mean(x(:))-x;              % flip matrix about its mean
d = imdilate(y,s);             % maxima of flipped data = local minima
v = find(y==d);                % find 1D indices of local minima
ln = max(length(p),length(v)); % length of output matrix
% return local maxima and minima as an (n x 2) matrix
ex = horzcat([x(p); ones(ln-length(p),1)*NaN],...
    [x(v); ones(ln-length(v),1)*NaN]);
% reshape 1D indices as row,column location vectors
[p,q] = ind2sub([r,c],p);    % maxima, ex(:,1) == x(p,q)
[v,u] = ind2sub([r,c],v);    % minima, ex(:,2) == x(v,u)