octave_packages/m/statistics/tests/cor_test.m

   1 ## Copyright (C) 1995-2012 Kurt Hornik
   2 ##
   3 ## This file is part of Octave.
   4 ##
   5 ## Octave is free software; you can redistribute it and/or modify it
   6 ## under the terms of the GNU General Public License as published by
   7 ## the Free Software Foundation; either version 3 of the License, or (at
   8 ## your option) any later version.
   9 ##
  10 ## Octave is distributed in the hope that it will be useful, but
  11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 ## General Public License for more details.
  14 ##
  15 ## You should have received a copy of the GNU General Public License
  16 ## along with Octave; see the file COPYING.  If not, see
  17 ## <http://www.gnu.org/licenses/>.
  18
  19 ## -*- texinfo -*-
  20 ## @deftypefn {Function File} {} cor_test (@var{x}, @var{y}, @var{alt}, @var{method})
  21 ## Test whether two samples @var{x} and @var{y} come from uncorrelated
  22 ## populations.
  23 ##
  24 ## The optional argument string @var{alt} describes the alternative
  25 ## hypothesis, and can be @code{"!="} or @code{"<>"} (non-zero),
  26 ## @code{">"} (greater than 0), or @code{"<"} (less than 0).  The
  27 ## default is the two-sided case.
  28 ##
  29 ## The optional argument string @var{method} specifies which
  30 ## correlation coefficient to use for testing.  If @var{method} is
  31 ## @code{"pearson"} (default), the (usual) Pearson's product moment
  32 ## correlation coefficient is used.  In this case, the data should come
  33 ## from a bivariate normal distribution.  Otherwise, the other two
  34 ## methods offer nonparametric alternatives.  If @var{method} is
  35 ## @code{"kendall"}, then Kendall's rank correlation tau is used.  If
  36 ## @var{method} is @code{"spearman"}, then Spearman's rank correlation
  37 ## rho is used.  Only the first character is necessary.
  38 ##
  39 ## The output is a structure with the following elements:
  40 ##
  41 ## @table @var
  42 ## @item pval
  43 ## The p-value of the test.
  44 ##
  45 ## @item stat
  46 ## The value of the test statistic.
  47 ##
  48 ## @item dist
  49 ## The distribution of the test statistic.
  50 ##
  51 ## @item params
  52 ## The parameters of the null distribution of the test statistic.
  53 ##
  54 ## @item alternative
  55 ## The alternative hypothesis.
  56 ##
  57 ## @item method
  58 ## The method used for testing.
  59 ## @end table
  60 ##
  61 ## If no output argument is given, the p-value is displayed.
  62 ## @end deftypefn
  63
  64 ## Author: FL <Friedrich.Leisch@ci.tuwien.ac.at>
  65 ## Adapted-by: KH <Kurt.Hornik@wu-wien.ac.at>
  66 ## Description: Test for zero correlation
  67
  68 function t = cor_test (x, y, alt, method)
  69
  70   if ((nargin < 2) || (nargin > 4))
  71     print_usage ();
  72   endif
  73
  74   if (!isvector (x) || !isvector (y) || length (x) != length (y))
  75     error ("cor_test: X and Y must be vectors of the same length");
  76   endif
  77
  78   if (nargin < 3)
  79     alt = "!=";
  80   elseif (! ischar (alt))
  81     error ("cor_test: ALT must be a string");
  82   endif
  83
  84   if (nargin < 4)
  85     method = "pearson";
  86   elseif (! ischar (method))
  87     error ("cor_test: METHOD must be a string");
  88   endif
  89
  90   n = length (x);
  91   m = method (1);
  92
  93   if (m == "p")
  94     r = corr (x, y);
  95     df = n - 2;
  96     t.method = "Pearson's product moment correlation";
  97     t.params = df;
  98     t.stat = sqrt (df) .* r / sqrt (1 - r.^2);
  99     t.dist = "t";
 100     cdf  = tcdf (t.stat, df);
 101   elseif (m == "k")
 102     tau = kendall (x, y);
 103     t.method = "Kendall's rank correlation tau";
 104     t.params = [];
 105     t.stat = tau / sqrt ((2 * (2*n+5)) / (9*n*(n-1)));
 106     t.dist = "stdnormal";
 107     cdf = stdnormal_cdf (t.stat);
 108   elseif (m == "s")
 109     rho = spearman (x, y);
 110     t.method = "Spearman's rank correlation rho";
 111     t.params = [];
 112     t.stat = sqrt (n-1) * (rho - 6/(n^3-n));
 113     t.dist = "stdnormal";
 114     cdf = stdnormal_cdf (t.stat);
 115   else
 116     error ("cor_test: METHOD `%s' not recognized", method);
 117   endif
 118
 119   if (strcmp (alt, "!=") || strcmp (alt, "<>"))
 120     t.pval = 2 * min (cdf, 1 - cdf);
 121   elseif (strcmp (alt, ">"))
 122     t.pval = 1 - cdf;
 123   elseif (strcmp (alt, "<"))
 124     t.pval = cdf;
 125   else
 126     error ("cor_test: alternative `%s' not recognized", alt);
 127   endif
 128
 129   t.alternative = alt;
 130
 131   if (nargout == 0)
 132     printf ("pval: %g\n", t.pval);
 133   endif
 134
 135 endfunction