    001    /*
    002     * Licensed to the Apache Software Foundation (ASF) under one or more
    003     * contributor license agreements.  See the NOTICE file distributed with
    004     * this work for additional information regarding copyright ownership.
    005     * The ASF licenses this file to You under the Apache License, Version 2.0
    006     * (the "License"); you may not use this file except in compliance with
    007     * the License.  You may obtain a copy of the License at
    008     *
    009     *      http://www.apache.org/licenses/LICENSE-2.0
    010     *
    011     * Unless required by applicable law or agreed to in writing, software
    012     * distributed under the License is distributed on an "AS IS" BASIS,
    013     * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    014     * See the License for the specific language governing permissions and
    015     * limitations under the License.
    016     */
    017    package org.apache.commons.math3.stat.inference;
    018
    019    import org.apache.commons.math3.distribution.NormalDistribution;
    020    import org.apache.commons.math3.exception.ConvergenceException;
    021    import org.apache.commons.math3.exception.MaxCountExceededException;
    022    import org.apache.commons.math3.exception.NoDataException;
    023    import org.apache.commons.math3.exception.NullArgumentException;
    024    import org.apache.commons.math3.stat.ranking.NaNStrategy;
    025    import org.apache.commons.math3.stat.ranking.NaturalRanking;
    026    import org.apache.commons.math3.stat.ranking.TiesStrategy;
    027    import org.apache.commons.math3.util.FastMath;
    028
    029    /**
    030     * An implementation of the Mann-Whitney U test (also called Wilcoxon rank-sum test).
    031     *
    032     * @version $Id: MannWhitneyUTest.java 1416643 2012-12-03 19:37:14Z tn $
    033     */
    034    public class MannWhitneyUTest {
    035
    036        /** Ranking algorithm. */
    037        private NaturalRanking naturalRanking;
    038
    039        /**
    040         * Create a test instance using where NaN's are left in place and ties get
    041         * the average of applicable ranks. Use this unless you are very sure of
    042         * what you are doing.
    043         */
    044        public MannWhitneyUTest() {
    045            naturalRanking = new NaturalRanking(NaNStrategy.FIXED,
    046                    TiesStrategy.AVERAGE);
    047        }
    048
    049        /**
    050         * Create a test instance using the given strategies for NaN's and ties.
    051         * Only use this if you are sure of what you are doing.
    052         *
    053         * @param nanStrategy
    054         *            specifies the strategy that should be used for Double.NaN's
    055         * @param tiesStrategy
    056         *            specifies the strategy that should be used for ties
    057         */
    058        public MannWhitneyUTest(final NaNStrategy nanStrategy,
    059                                final TiesStrategy tiesStrategy) {
    060            naturalRanking = new NaturalRanking(nanStrategy, tiesStrategy);
    061        }
    062
    063        /**
    064         * Ensures that the provided arrays fulfills the assumptions.
    065         *
    066         * @param x first sample
    067         * @param y second sample
    068         * @throws NullArgumentException if {@code x} or {@code y} are {@code null}.
    069         * @throws NoDataException if {@code x} or {@code y} are zero-length.
    070         */
    071        private void ensureDataConformance(final double[] x, final double[] y)
    072            throws NullArgumentException, NoDataException {
    073
    074            if (x == null ||
    075                y == null) {
    076                throw new NullArgumentException();
    077            }
    078            if (x.length == 0 ||
    079                y.length == 0) {
    080                throw new NoDataException();
    081            }
    082        }
    083
    084        /** Concatenate the samples into one array.
    085         * @param x first sample
    086         * @param y second sample
    087         * @return concatenated array
    088         */
    089        private double[] concatenateSamples(final double[] x, final double[] y) {
    090            final double[] z = new double[x.length + y.length];
    091
    092            System.arraycopy(x, 0, z, 0, x.length);
    093            System.arraycopy(y, 0, z, x.length, y.length);
    094
    095            return z;
    096        }
    097
    098        /**
    099         * Computes the <a
    100         * href="http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U"> Mann-Whitney
    101         * U statistic</a> comparing mean for two independent samples possibly of
    102         * different length.
    103         * <p>
    104         * This statistic can be used to perform a Mann-Whitney U test evaluating
    105         * the null hypothesis that the two independent samples has equal mean.
    106         * </p>
    107         * <p>
    108         * Let X<sub>i</sub> denote the i'th individual of the first sample and
    109         * Y<sub>j</sub> the j'th individual in the second sample. Note that the
    110         * samples would often have different length.
    111         * </p>
    112         * <p>
    113         * <strong>Preconditions</strong>:
    114         * <ul>
    115         * <li>All observations in the two samples are independent.</li>
    116         * <li>The observations are at least ordinal (continuous are also ordinal).</li>
    117         * </ul>
    118         * </p>
    119         *
    120         * @param x the first sample
    121         * @param y the second sample
    122         * @return Mann-Whitney U statistic (maximum of U<sup>x</sup> and U<sup>y</sup>)
    123         * @throws NullArgumentException if {@code x} or {@code y} are {@code null}.
    124         * @throws NoDataException if {@code x} or {@code y} are zero-length.
    125         */
    126        public double mannWhitneyU(final double[] x, final double[] y)
    127            throws NullArgumentException, NoDataException {
    128
    129            ensureDataConformance(x, y);
    130
    131            final double[] z = concatenateSamples(x, y);
    132            final double[] ranks = naturalRanking.rank(z);
    133
    134            double sumRankX = 0;
    135
    136            /*
    137             * The ranks for x is in the first x.length entries in ranks because x
    138             * is in the first x.length entries in z
    139             */
    140            for (int i = 0; i < x.length; ++i) {
    141                sumRankX += ranks[i];
    142            }
    143
    144            /*
    145             * U1 = R1 - (n1 * (n1 + 1)) / 2 where R1 is sum of ranks for sample 1,
    146             * e.g. x, n1 is the number of observations in sample 1.
    147             */
    148            final double U1 = sumRankX - (x.length * (x.length + 1)) / 2;
    149
    150            /*
    151             * It can be shown that U1 + U2 = n1 * n2
    152             */
    153            final double U2 = x.length * y.length - U1;
    154
    155            return FastMath.max(U1, U2);
    156        }
    157
    158        /**
    159         * @param Umin smallest Mann-Whitney U value
    160         * @param n1 number of subjects in first sample
    161         * @param n2 number of subjects in second sample
    162         * @return two-sided asymptotic p-value
    163         * @throws ConvergenceException if the p-value can not be computed
    164         * due to a convergence error
    165         * @throws MaxCountExceededException if the maximum number of
    166         * iterations is exceeded
    167         */
    168        private double calculateAsymptoticPValue(final double Umin,
    169                                                 final int n1,
    170                                                 final int n2)
    171            throws ConvergenceException, MaxCountExceededException {
    172
    173            /* long multiplication to avoid overflow (double not used due to efficiency
    174             * and to avoid precision loss)
    175             */
    176            final long n1n2prod = (long) n1 * n2;
    177
    178            // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
    179            final double EU = n1n2prod / 2.0;
    180            final double VarU = n1n2prod * (n1 + n2 + 1) / 12.0;
    181
    182            final double z = (Umin - EU) / FastMath.sqrt(VarU);
    183
    184            // No try-catch or advertised exception because args are valid
    185            final NormalDistribution standardNormal = new NormalDistribution(0, 1);
    186
    187            return 2 * standardNormal.cumulativeProbability(z);
    188        }
    189
    190        /**
    191         * Returns the asymptotic <i>observed significance level</i>, or <a href=
    192         * "http://www.cas.lancs.ac.uk/glossary_v1.1/hyptest.html#pvalue">
    193         * p-value</a>, associated with a <a
    194         * href="http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U"> Mann-Whitney
    195         * U statistic</a> comparing mean for two independent samples.
    196         * <p>
    197         * Let X<sub>i</sub> denote the i'th individual of the first sample and
    198         * Y<sub>j</sub> the j'th individual in the second sample. Note that the
    199         * samples would often have different length.
    200         * </p>
    201         * <p>
    202         * <strong>Preconditions</strong>:
    203         * <ul>
    204         * <li>All observations in the two samples are independent.</li>
    205         * <li>The observations are at least ordinal (continuous are also ordinal).</li>
    206         * </ul>
    207         * </p><p>
    208         * Ties give rise to biased variance at the moment. See e.g. <a
    209         * href="http://mlsc.lboro.ac.uk/resources/statistics/Mannwhitney.pdf"
    210         * >http://mlsc.lboro.ac.uk/resources/statistics/Mannwhitney.pdf</a>.</p>
    211         *
    212         * @param x the first sample
    213         * @param y the second sample
    214         * @return asymptotic p-value
    215         * @throws NullArgumentException if {@code x} or {@code y} are {@code null}.
    216         * @throws NoDataException if {@code x} or {@code y} are zero-length.
    217         * @throws ConvergenceException if the p-value can not be computed due to a
    218         * convergence error
    219         * @throws MaxCountExceededException if the maximum number of iterations
    220         * is exceeded
    221         */
    222        public double mannWhitneyUTest(final double[] x, final double[] y)
    223            throws NullArgumentException, NoDataException,
    224            ConvergenceException, MaxCountExceededException {
    225
    226            ensureDataConformance(x, y);
    227
    228            final double Umax = mannWhitneyU(x, y);
    229
    230            /*
    231             * It can be shown that U1 + U2 = n1 * n2
    232             */
    233            final double Umin = x.length * y.length - Umax;
    234
    235            return calculateAsymptoticPValue(Umin, x.length, y.length);
    236        }
    237
    238    }