001: /*
002: * Copyright 2003-2004 The Apache Software Foundation.
003: *
004: * Licensed under the Apache License, Version 2.0 (the "License");
005: * you may not use this file except in compliance with the License.
006: * You may obtain a copy of the License at
007: *
008: * http://www.apache.org/licenses/LICENSE-2.0
009: *
010: * Unless required by applicable law or agreed to in writing, software
011: * distributed under the License is distributed on an "AS IS" BASIS,
012: * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
013: * See the License for the specific language governing permissions and
014: * limitations under the License.
015: */
016: package org.apache.commons.math.stat.descriptive.moment;
017:
018: import org.apache.commons.math.stat.descriptive.AbstractStorelessUnivariateStatistic;
019:
020: /**
021: * Computes the Kurtosis of the available values.
022: * <p>
023: * We use the following (unbiased) formula to define kurtosis:
024: * <p>
025: * kurtosis = { [n(n+1) / (n -1)(n - 2)(n-3)] sum[(x_i - mean)^4] / std^4 } - [3(n-1)^2 / (n-2)(n-3)]
026: * <p>
027: * where n is the number of values, mean is the {@link Mean} and std is the
028: * {@link StandardDeviation}
029: * <p>
030: * Note that this statistic is undefined for n < 4. <code>Double.Nan</code>
031: * is returned when there is not sufficient data to compute the statistic.
032: * <p>
033: * <strong>Note that this implementation is not synchronized.</strong> If
034: * multiple threads access an instance of this class concurrently, and at least
035: * one of the threads invokes the <code>increment()</code> or
036: * <code>clear()</code> method, it must be synchronized externally.
037: *
038: * @version $Revision: 348519 $ $Date: 2005-11-23 12:12:18 -0700 (Wed, 23 Nov 2005) $
039: */
040: public class Kurtosis extends AbstractStorelessUnivariateStatistic {
041:
042: /** Serializable version identifier */
043: private static final long serialVersionUID = 2784465764798260919L;
044:
045: /**Fourth Moment on which this statistic is based */
046: protected FourthMoment moment;
047:
048: /**
049: * Determines whether or not this statistic can be incremented or cleared.
050: * <p>
051: * Statistics based on (constructed from) external moments cannot
052: * be incremented or cleared.
053: */
054: protected boolean incMoment;
055:
056: /**
057: * Construct a Kurtosis
058: */
059: public Kurtosis() {
060: incMoment = true;
061: moment = new FourthMoment();
062: }
063:
064: /**
065: * Construct a Kurtosis from an external moment
066: *
067: * @param m4 external Moment
068: */
069: public Kurtosis(final FourthMoment m4) {
070: incMoment = false;
071: this .moment = m4;
072: }
073:
074: /**
075: * @see org.apache.commons.math.stat.descriptive.StorelessUnivariateStatistic#increment(double)
076: */
077: public void increment(final double d) {
078: if (incMoment) {
079: moment.increment(d);
080: } else {
081: throw new IllegalStateException(
082: "Statistics constructed from external moments cannot be incremented");
083: }
084: }
085:
086: /**
087: * @see org.apache.commons.math.stat.descriptive.StorelessUnivariateStatistic#getResult()
088: */
089: public double getResult() {
090: double kurtosis = Double.NaN;
091: if (moment.getN() > 3) {
092: double variance = moment.m2 / (double) (moment.n - 1);
093: if (moment.n <= 3 || variance < 10E-20) {
094: kurtosis = 0.0;
095: } else {
096: double n = (double) moment.n;
097: kurtosis = (n * (n + 1) * moment.m4 - 3 * moment.m2
098: * moment.m2 * (n - 1))
099: / ((n - 1) * (n - 2) * (n - 3) * variance * variance);
100: }
101: }
102: return kurtosis;
103: }
104:
105: /**
106: * @see org.apache.commons.math.stat.descriptive.StorelessUnivariateStatistic#clear()
107: */
108: public void clear() {
109: if (incMoment) {
110: moment.clear();
111: } else {
112: throw new IllegalStateException(
113: "Statistics constructed from external moments cannot be cleared");
114: }
115: }
116:
117: /**
118: * @see org.apache.commons.math.stat.descriptive.StorelessUnivariateStatistic#getN()
119: */
120: public long getN() {
121: return moment.getN();
122: }
123:
124: /* UnvariateStatistic Approach */
125:
126: /**
127: * Returns the kurtosis of the entries in the specified portion of the
128: * input array.
129: * <p>
130: * See {@link Kurtosis} for details on the computing algorithm.
131: * <p>
132: * Throws <code>IllegalArgumentException</code> if the array is null.
133: *
134: * @param values the input array
135: * @param begin index of the first array element to include
136: * @param length the number of elements to include
137: * @return the kurtosis of the values or Double.NaN if length is less than
138: * 4
139: * @throws IllegalArgumentException if the input array is null or the array
140: * index parameters are not valid
141: */
142: public double evaluate(final double[] values, final int begin,
143: final int length) {
144: // Initialize the kurtosis
145: double kurt = Double.NaN;
146:
147: if (test(values, begin, length) && length > 3) {
148:
149: // Compute the mean and standard deviation
150: Variance variance = new Variance();
151: variance.incrementAll(values, begin, length);
152: double mean = variance.moment.m1;
153: double stdDev = Math.sqrt(variance.getResult());
154:
155: // Sum the ^4 of the distance from the mean divided by the
156: // standard deviation
157: double accum3 = 0.0;
158: for (int i = begin; i < begin + length; i++) {
159: accum3 += Math.pow((values[i] - mean), 4.0);
160: }
161: accum3 /= Math.pow(stdDev, 4.0d);
162:
163: // Get N
164: double n0 = length;
165:
166: double coefficientOne = (n0 * (n0 + 1))
167: / ((n0 - 1) * (n0 - 2) * (n0 - 3));
168: double termTwo = ((3 * Math.pow(n0 - 1, 2.0)) / ((n0 - 2) * (n0 - 3)));
169:
170: // Calculate kurtosis
171: kurt = (coefficientOne * accum3) - termTwo;
172: }
173: return kurt;
174: }
175:
176: }
|