octave_packages/nan-2.5.5/ttest2.m

   1 function [h, pval, ci, stats, df] = ttest2 (x, y, alpha, tail, vartype, DIM)
   2 % TTEST2 (unpaired) t-test
   3 %     For two samples x and y from normal distributions with unknown
   4 %     means and unknown equal variances, perform a two-sample t-test of
   5 %     the null hypothesis of equal means.  Under the null, the test
   6 %     statistic T follows a Student distribution with DF degrees of
   7 %     freedom.
   8 %
   9 %     TTEST2 treads NaNs as "Missing values" and ignores these.
  10 %
  11 % H = ttest2(x,y)
  12 % H = ttest2([x;y],C,W)
  13 % H = ttest2(x,y,alpha)
  14 % H = ttest2(x,y,alpha,tail)
  15 % H = ttest2(x,y,alpha,tail,vartype)
  16 % H = ttest2(x,y,alpha,tail,vartype,DIM)
  17 % [H,PVAL] = ttest2(...)
  18 % [h,p,ci,stats] = ttest2(...)
  19 %
  20 %     H=1 indicates a rejection of the Null-hypothesis at a significance
  21 %     level of alpha (default alpha = 0.05).
  22 %
  23 %     With the optional argument string TAIL, the Alternative of interest
  24 %     can be selected.  If TAIL is '!=' or '<>' or 'both', the null is tested
  25 %     against the two-sided Alternative `mean (X) ~= mean (Y)'.  If TAIL
  26 %     is '>' or 'right', the one-sided Alternative `mean (X) > mean (Y)' is used.
  27 %     Similarly for '<' or 'left', the one-sided Alternative `mean (X) < mean
  28 %     (Y)' is used.  The default is the two-sided case.
  29 %
  30 %     vartype support only 'equal' (default value); the value 'unequal' is not supported.
  31 %
  32 %     H returns whether the Null-Hypotheses must be rejected.
  33 %     The p-value of the test is returned in PVAL.
  34 %
  35 %     TTEST2 works on the first non-singleton dimension or on DIM.
  36 %
  37 %     If no output argument is given, the p-value of the test is
  38 %     displayed.
  39 %
  40
  41 %%% not supported yet
  42 % [h,p,ci] = ttest2(...)
  43 % [h,p,ci,stats] = ttest2(...)
  44
  45 %       $Id$
  46 %       Copyright (C) 1995, 1996, 1997, 1998, 2000, 2002, 2005, 2006, 2007
  47 %               Kurt Hornik
  48 %       Copyright (C) 2010 by Alois Schloegl <alois.schloegl@gmail.com>
  49 %       This function is part of the NaN-toolbox
  50 %       http://pub.ist.ac.at/~schloegl/matlab/NaN/
  51
  52 %    This program is free software: you can redistribute it and/or modify
  53 %    it under the terms of the GNU General Public License as published by
  54 %    the Free Software Foundation, either version 3 of the License, or
  55 %    (at your option) any later version.
  56 %
  57 %    This program is distributed in the hope that it will be useful,
  58 %    but WITHOUT ANY WARRANTY; without even the implied warranty of
  59 %    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  60 %    GNU General Public License for more details.
  61 %
  62 %    You should have received a copy of the GNU General Public License
  63 %    along with this program.  If not, see <http://www.gnu.org/licenses/>.
  64
  65   if ((nargin < 2) || (nargin > 6) || nargout > 4)
  66         print_usage ;
  67   end
  68
  69   if (nargin < 3) || isempty(alpha)
  70     alpha = .05;
  71     W = [];
  72   elseif size(x,1)==size(y,1) && size(x,1)==size(alpha,1) && size(y,2)==1 && all(y==1 | y==2)
  73     c = y;
  74     y = x(c==2,:);
  75     x = x(c==1,:);
  76     Wy = alpha(c==2);
  77     Wx = alpha(c==1);
  78     alpha = .05;
  79   end
  80
  81   if (nargin < 4) || isempty(tail)
  82     tail = '~=';
  83   end
  84   if (~ ischar (tail))
  85     error ('ttest2: tail must be a string');
  86   end
  87   if (nargin < 5) || isempty(vartype)
  88     vartype = 'equal';
  89   end
  90   if ~strcmp(vartype,'equal')
  91     error ('test: vartype not supported')
  92   end
  93   if nargin<6,
  94        DIM = find(size(x)>1,1);
  95   end;
  96   if isempty(DIM), DIM=1; end;
  97
  98   szx = size(x);
  99   szy = size(y);
 100   szy(DIM) = 1;
 101   szx(DIM) = 1;
 102
 103   if (any(szx-szy))
 104     error ('ttest2: dimension of X and Y do not fit');
 105   end
 106
 107   [SX, NX] = sumskipnan(x, DIM, Wx);
 108   [SY, NY] = sumskipnan(y, DIM, Wy);
 109   stats.df = NX + NY - 2;
 110   MX = SX ./ NX;
 111   MY = SY ./ NY;
 112
 113   if any(size(x)==0) || any(size(y)==0)
 114         v = NaN;
 115   else
 116         v = sumsq(x-repmat(MX,size(x)./size(MX))) + sumsq(y-repmat(MY,size(y)./size(MY)));
 117   end;
 118   stats.sd    = sqrt(v/stats.df);
 119   stats.tstat = (MX - MY) .* sqrt ((NX .* NY .* stats.df) ./ (v .* (NX + NY)));
 120   cdf         = tcdf (stats.tstat, stats.df);
 121
 122   if (strcmp (tail, '~=') || strcmp (tail, '!=') || strcmp (tail, '<>')) || strcmp(tail,'both'),
 123     pval = 2 * min (cdf, 1 - cdf);
 124   elseif strcmp (tail, '>') || strcmp(tail,'right'),
 125     pval = 1 - cdf;
 126   elseif strcmp (tail, '<') || strcmp(tail,'left'),
 127     pval = cdf;
 128   else
 129     error ('ttest2: option %s not recognized', tail);
 130   end
 131
 132   h = pval < alpha;
 133
 134   if (nargout == 0)
 135     fprintf(1,'  pval: %g\n', pval);
 136   end
 137