1 function Q=quantile(Y,q,DIM,method)
2 % QUANTILE calculates the quantiles of histograms and sample arrays.
5 % Q = quantile(Y,q,DIM)
6 % returns the q-th quantile along dimension DIM of sample array Y.
7 % size(Q) is equal size(Y) except for dimension DIM which is size(Q,DIM)=length(Q)
10 % returns the q-th quantile from the histogram HIS.
11 % HIS must be a HISTOGRAM struct as defined in HISTO2 or HISTO3.
12 % If q is a vector, the each row of Q returns the q(i)-th quantile
14 % see also: HISTO2, HISTO3, PERCENTILE
17 % $Id: quantile.m 9601 2012-02-09 14:14:36Z schloegl $
18 % Copyright (C) 1996-2003,2005,2006,2007,2009,2011 by Alois Schloegl <alois.schloegl@gmail.com>
19 % This function is part of the NaN-toolbox
20 % http://pub.ist.ac.at/~schloegl/matlab/NaN/
22 % This program is free software; you can redistribute it and/or modify
23 % it under the terms of the GNU General Public License as published by
24 % the Free Software Foundation; either version 3 of the License, or
25 % (at your option) any later version.
27 % This program is distributed in the hope that it will be useful,
28 % but WITHOUT ANY WARRANTY; without even the implied warranty of
29 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 % GNU General Public License for more details.
32 % You should have received a copy of the GNU General Public License
33 % along with this program; If not, see <http://www.gnu.org/licenses/>.
39 DIM = find(size(Y)>1,1);
40 if isempty(DIM), DIM = 1; end;
48 [q, rix] = sort(q(:)'); % sort quantile values
49 [tmp,rix] = sort(rix); % generate reverse index
52 if SW, SW = isfield(Y,'datatype'); end;
53 if SW, SW = strcmp(Y.datatype,'HISTOGRAM'); end;
56 Q = repmat(nan,length(q),yc);
63 h = full(Y.H(tmp,k1));
64 t = Y.X(tmp,min(size(Y.X,2),k1));
67 t2(1:2:2*length(t)) = t;
68 t2(2:2:2*length(t)) = t;
71 x(2:2:2*length(t)) = x2;
72 x(3:2:2*length(t)) = x2(1:end-1);
74 % Q(q < 0 | 1 < q,:) = NaN; % already done at initialization
78 for k2 = find( (0 < q) & (q < 1) )
79 while (q(k2)*N > x(n)),
84 % mean of upper and lower bound
85 Q(k2,k1) = (t2(n)+t2(n+1))/2;
90 Q = Q(rix,:); % order resulting quantiles according to original input q
97 sz = [sz,ones(1,DIM-length(sz))];
100 f = zeros(1,length(q));
101 f( (q < 0) | (1 < q) ) = NaN;
102 D1 = prod(sz(1:DIM-1));
103 D3 = prod(sz(DIM+1:length(sz)));
104 Q = repmat(nan,[sz(1:DIM-1),length(q),sz(DIM+1:length(sz))]);
107 xi = k + l * D1*sz(DIM) + 1 ;
108 xo = k + l * D1*length(q) + 1;
109 t = Y(xi:D1:xi+D1*sz(DIM)-1);
117 t2(1:2:2*length(t)) = t;
118 t2(2:2:2*length(t)) = t;
119 x = floor((1:2*length(t))/2);
120 %f(q < 0 | 1 < q) = NaN; % for efficiency its defined outside loop
125 for k2 = find( (0 < q) & (q < 1) )
126 while (q(k2)*N > x(n)),
131 % mean of upper and lower bound
132 f(k2) = (t2(n) + t2(n+1))/2;
138 Q(xo:D1:xo + D1*length(q) - 1) = f(rix);
143 fprintf(2,'Error QUANTILES: invalid input argument\n');
149 %!assert(quantile(1:10,[.2,.5]),[2.5, 5.5])