2 % Y2RES evaluates basic statistics of a data series
5 % several statistics are estimated from each column of y
8 % R.N number of samples, NaNs are not counted
11 % R.STD standard deviation
15 % ... and many more including:
16 % MEDIAN, Quartiles, Variance, standard error of the mean (SEM),
17 % Coefficient of Variation, Quantization (QUANT), TRIMEAN, SKEWNESS,
18 % KURTOSIS, Root-Mean-Square (RMS), ENTROPY
21 % $Id: y2res.m 5090 2008-06-05 08:12:04Z schloegl $
22 % Copyright (C) 1996-2005,2008 by Alois Schloegl <a.schloegl@ieee.org>
23 % This is part of the TSA-toolbox
24 % http://octave.sourceforge.net/
25 % http://www.dpmi.tugraz.at/~schloegl/matlab/tsa/
27 % This program is free software: you can redistribute it and/or modify
28 % it under the terms of the GNU General Public License as published by
29 % the Free Software Foundation, either version 3 of the License, or
30 % (at your option) any later version.
32 % This program is distributed in the hope that it will be useful,
33 % but WITHOUT ANY WARRANTY; without even the implied warranty of
34 % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
35 % GNU General Public License for more details.
37 % You should have received a copy of the GNU General Public License
38 % along with this program. If not, see <http://www.gnu.org/licenses/>.
41 [R.SUM, R.N, R.SSQ] = sumskipnan(Y,1);
42 %R.S3P = sumskipnan(Y.^3,1);
43 R.S4P = sumskipnan(Y.^4,1);
44 %R.S5P = sumskipnan(Y.^5,1);
49 R.SSQ0 = R.SSQ-R.SUM.*R.MEAN; % sum square of mean removed
51 if 1,%flag_implicit_unbiased_estim,
52 n1 = max(R.N-1,0); % in case of n=0 and n=1, the (biased) variance, STD and STE are INF
57 R.VAR = R.SSQ0./n1; % variance (unbiased)
58 R.STD = sqrt(R.VAR); % standard deviation
59 R.SEM = sqrt(R.SSQ0./(R.N.*n1)); % standard error of the mean
60 R.SEV = sqrt(n1.*(n1.*R.S4P./R.N+(R.N.^2-2*R.N+3).*(R.SSQ./R.N).^2)./(R.N.^3)); % standard error of the variance
61 R.Coefficient_of_variation = R.STD./R.MEAN;
68 %R.NormEntropy = log2(sqrt(2*pi*exp(1)))+log2(R.STD);
70 Q0500=repmat(nan,1,size(Y,2));
76 Q0250(k) = flix(tmp,R.N(k)/4 + 0.75);
77 Q0500(k) = flix(tmp,R.N(k)/2 + 0.50);
78 Q0750(k) = flix(tmp,R.N(k)*3/4 + 0.25);
81 pdf = diff([0; find(tmp>0); R.N(k)])/R.N(k); % empirical probability distribution
82 R.ENTROPY(k) = -sumskipnan(pdf.*log(pdf));
91 qerror = max(abs(tmp-round(tmp)));
96 if any(R.Qerror*1e6>R.QUANT)
97 warning('(Y2RES) Quantization might not be equidistant')
101 R.Quartiles = [Q0250; Q0750];
102 % R.IQR = H_spread = [Q0750 - Q0250];
103 R.TRIMEAN = [Q0250 + 2*Q0500 + Q0750]/4;
105 Y = Y - repmat(R.MEAN,size(Y)./size(R.MEAN));
106 R.CM3 = sumskipnan(Y.^3,1)./n1;
107 R.CM4 = sumskipnan(Y.^4,1)./n1;
108 %R.CM5 = sumskipnan(Y.^5,1)./n1;
110 R.SKEWNESS = R.CM3./(R.STD.^3);
111 R.KURTOSIS = R.CM4./(R.VAR.^2)-3;
113 %R.Skewness.Fisher = (R.CM3)./(R.STD.^3); %%% same as R.SKEWNESS
115 %R.Skewness.Pearson_Mode = (R.MEAN-R.MODE)./R.STD;
116 %R.Skewness.Pearson_coeff1 = (3*R.MEAN-R.MODE)./R.STD;
117 R.Skewness.Pearson_coeff2 = (3*R.MEAN-R.MEDIAN)./R.STD;
118 R.Skewness.Bowley = (Q0750+Q0250 - 2*Q0500)./(Q0750-Q0250); % quartile skewness coefficient
120 R.datatype = 'STAT Level 4';