## Copyright (C) 2016 John W. Eaton ## ## This file is part of Octave. ## ## Octave is free software; you can redistribute it and/or modify it ## under the terms of the GNU General Public License as published by ## the Free Software Foundation; either version 3 of the License, or (at ## your option) any later version. ## ## Octave is distributed in the hope that it will be useful, but ## WITHOUT ANY WARRANTY; without even the implied warranty of ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU ## General Public License for more details. ## ## You should have received a copy of the GNU General Public License ## along with Octave; see the file COPYING. If not, see ## . ## -*- texinfo -*- ## @deftypefn {} {@var{R} =} corrcoef (@var{X}) ## @deftypefnx {} {@var{R} =} corrcoef (@var{X}, @var{Y}) ## @deftypefnx {} {[@var{R}, @var{P}] =} corrcoef (@dots{}) ## Compute a matrix of correlation coefficients. ## ## @var{X} is an array where each column contains a variable. ## ## If provided, @var{Y} is a vector of the same size than @var{X} and ## corrcoef([X(:), Y(:)]) is returned. ## ## @var{R} is a matrix of Pearson's product moment correlation ## coefficients for each pair of variables. ## ## @var{P} is a matrix of pair-wise p-values testing for the null ## hypothesis of a correlation coefficient of zero. ## ## @seealso{corr,cor_test} ## @end deftypefn ## Author: jwe function [R, P] = corrcoef (X, varargin) if (nargin == 0) print_usage (); endif alpha = 0.05; rows = "all"; if (nargin > 1) if (isnumeric (varargin{1})) X = [X(:) varargin{1}(:)]; else for i=1:2:numel(varargin) if (ischar (varargin{i})) parameter = varargin{i}; else error ("Wrong type argument '%s'.",class (varargin{i})); endif if (nargin < i+2) error ("Parameter '%s' missing value.", parameter); else value = varargin{i+1}; endif switch (lower (parameter)) case "alpha" if (isnumeric(value) && numel(value) == 1 && value >=0 && value <= 1) alpha = value; else error ("'alpha' is a number between 0 and 1."); endif case "rows" if (ischar(value)) switch lower(value) case {"all", "complete", "pairwise"} rows = value; otherwise error ("'rows' is either 'all', 'complete' or 'pairwise'"); end else error ("'rows' is a string."); endif error ("Option 'rows' not implemented."); otherwise error ("Unknown option '%s'.", parameter); endswitch endfor endif endif N = size (X,2); R = eye (N); if (nargout > 1) P = eye (N); endif for i=1:N for j=i+1:N R(i,j) = corr (X(:,i), X(:,j)); R(j,i) = R(i,j); if (nargout > 1) T = cor_test (X(:,i), X(:,j), "!=", "pearson"); P(i,j) = T.pval; P(j,i) = P(i,j); endif endfor endfor endfunction %!test %! X = rand (5); %! R = corrcoef (X); %! assert (size (R) == [5, 5]); %!test %! X = rand (5); %! [R, P] = corrcoef (X); %! assert (size (R) == [5, 5] && size(P) == [5 5]); %!test %! X = rand (5,1); %! Y = rand (5,1); %! R1 = corrcoef (X, Y); %! R2 = corrcoef ([X, Y]); %! assert (R1, R2, sqrt (eps));