commons-dev mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From Mikkel Meyer Andersen <m...@mikl.dk>
Subject Re: svn commit: r1053836 - in /commons/proper/math/trunk/src: main/java/org/apache/commons/math/stat/inference/ test/java/org/apache/commons/math/stat/inference/
Date Thu, 30 Dec 2010 10:00:47 GMT
Hi,

This commit consists of new files, so I put
@version $Revision: $ $Date: $
in the top. I thought SVN would then auto update these, but after
checking out again to see if it worked, I see that it didn't. What
should I do to make it work? Fire a svn propXXX, use other variables,
or?

Cheers, Mikkel.

2010/12/30  <mikl@apache.org>:
> Author: mikl
> Date: Thu Dec 30 09:52:00 2010
> New Revision: 1053836
>
> URL: http://svn.apache.org/viewvc?rev=1053836&view=rev
> Log:
> Starting point for fixing MATH-431. Still some validation tests and documentation are
missing.
>
> Added:
>    commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTest.java
>    commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTestImpl.java
>    commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTest.java
>    commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestImpl.java
>    commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/MannWhitneyUTestTest.java
>    commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestTest.java
>
> Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTest.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTest.java?rev=1053836&view=auto
> ==============================================================================
> --- commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTest.java
(added)
> +++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTest.java
Thu Dec 30 09:52:00 2010
> @@ -0,0 +1,96 @@
> +/*
> + * Licensed to the Apache Software Foundation (ASF) under one or more
> + * contributor license agreements.  See the NOTICE file distributed with
> + * this work for additional information regarding copyright ownership.
> + * The ASF licenses this file to You under the Apache License, Version 2.0
> + * (the "License"); you may not use this file except in compliance with
> + * the License.  You may obtain a copy of the License at
> + *
> + *      http://www.apache.org/licenses/LICENSE-2.0
> + *
> + * Unless required by applicable law or agreed to in writing, software
> + * distributed under the License is distributed on an "AS IS" BASIS,
> + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
> + * See the License for the specific language governing permissions and
> + * limitations under the License.
> + */
> +package org.apache.commons.math.stat.inference;
> +
> +import org.apache.commons.math.MathException;
> +
> +/**
> + * An interface for Mann-Whitney U test (also called Wilcoxon rank-sum test).
> + *
> + * @version $Revision: $ $Date: $
> + */
> +public interface MannWhitneyUTest {
> +
> +    /**
> +     * Computes the <a
> +     * href="http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U"> Mann-Whitney
> +     * U statistic</a> comparing mean for two independent samples possibly of
> +     * different length.
> +     * <p>
> +     * This statistic can be used to perform a Mann-Whitney U test evaluating
> +     * the null hypothesis that the two independent samples has equal mean.
> +     * </p>
> +     * <p>
> +     * Let X<sub>i</sub> denote the i'th individual of the first sample
and
> +     * Y<sub>j</sub> the j'th individual in the second sample. Note that
the
> +     * samples would often have different length.
> +     * </p>
> +     * <p>
> +     * <strong>Preconditions</strong>:
> +     * <ul>
> +     * <li>All observations in the two samples are independent.</li>
> +     * <li>The observations are at least ordinal (continuous are also ordinal).</li>
> +     * </ul>
> +     * </p>
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @return mannWhitneyU statistic
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     */
> +    double mannWhitneyU(final double[] x, final double[] y)
> +            throws IllegalArgumentException;
> +
> +    /**
> +     * Returns the asymptotic <i>observed significance level</i>, or <a
href=
> +     * "http://www.cas.lancs.ac.uk/glossary_v1.1/hyptest.html#pvalue">
> +     * p-value</a>, associated with a <a
> +     * href="http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U"> Mann-Whitney
> +     * U statistic</a> comparing mean for two independent samples.
> +     * <p>
> +     * Let X<sub>i</sub> denote the i'th individual of the first sample
and
> +     * Y<sub>j</sub> the j'th individual in the second sample. Note that
the
> +     * samples would often have different length.
> +     * </p>
> +     * <p>
> +     * <strong>Preconditions</strong>:
> +     * <ul>
> +     * <li>All observations in the two samples are independent.</li>
> +     * <li>The observations are at least ordinal (continuous are also ordinal).</li>
> +     * </ul>
> +     * </p>
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @param exactPValue
> +     *            if the exact p-value is wanted (only works for x.length <=
30,
> +     *            if true and x.length > 30, this is ignored because
> +     *            calculations may take too long)
> +     * @return asymptotic p-value
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     * @throws MathException
> +     *             if an error occurs computing the p-value
> +     */
> +    double mannWhitneyUTest(final double[] x, final double[] y)
> +            throws IllegalArgumentException, MathException;
> +}
>
> Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTestImpl.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTestImpl.java?rev=1053836&view=auto
> ==============================================================================
> --- commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTestImpl.java
(added)
> +++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/MannWhitneyUTestImpl.java
Thu Dec 30 09:52:00 2010
> @@ -0,0 +1,199 @@
> +/*
> + * Licensed to the Apache Software Foundation (ASF) under one or more
> + * contributor license agreements.  See the NOTICE file distributed with
> + * this work for additional information regarding copyright ownership.
> + * The ASF licenses this file to You under the Apache License, Version 2.0
> + * (the "License"); you may not use this file except in compliance with
> + * the License.  You may obtain a copy of the License at
> + *
> + *      http://www.apache.org/licenses/LICENSE-2.0
> + *
> + * Unless required by applicable law or agreed to in writing, software
> + * distributed under the License is distributed on an "AS IS" BASIS,
> + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
> + * See the License for the specific language governing permissions and
> + * limitations under the License.
> + */
> +package org.apache.commons.math.stat.inference;
> +
> +import org.apache.commons.math.MathException;
> +import org.apache.commons.math.distribution.NormalDistributionImpl;
> +import org.apache.commons.math.stat.ranking.NaNStrategy;
> +import org.apache.commons.math.stat.ranking.NaturalRanking;
> +import org.apache.commons.math.stat.ranking.TiesStrategy;
> +import org.apache.commons.math.util.FastMath;
> +
> +/**
> + * An implementation of the Mann-Whitney U test (also called Wilcoxon rank-sum
> + * test).
> + *
> + * @version $Revision: $ $Date: $
> + */
> +public class MannWhitneyUTestImpl implements MannWhitneyUTest {
> +    private NaturalRanking naturalRanking;
> +
> +    /**
> +     * Create a test instance using where NaN's are left in place and ties get
> +     * the average of applicable ranks. Use this unless you are very sure of
> +     * what you are doing.
> +     */
> +    public MannWhitneyUTestImpl() {
> +        naturalRanking = new NaturalRanking(NaNStrategy.FIXED,
> +                TiesStrategy.AVERAGE);
> +    }
> +
> +    /**
> +     * Create a test instance using the given strategies for NaN's and ties.
> +     * Only use this if you are sure of what you are doing.
> +     *
> +     * @param nanStrategy
> +     *            specifies the strategy that should be used for Double.NaN's
> +     * @param tiesStrategy
> +     *            specifies the strategy that should be used for ties
> +     */
> +    public MannWhitneyUTestImpl(NaNStrategy nanStrategy,
> +            TiesStrategy tiesStrategy) {
> +        naturalRanking = new NaturalRanking(nanStrategy, tiesStrategy);
> +    }
> +
> +    /**
> +     * Ensures that the provided arrays fulfills the assumptions.
> +     *
> +     * @param x
> +     * @param y
> +     * @throws IllegalArgumentException
> +     *             if assumptions are not met
> +     */
> +    private void ensureDataConformance(final double[] x, final double[] y)
> +            throws IllegalArgumentException {
> +        if (x == null) {
> +            throw new IllegalArgumentException("x must not be null");
> +        }
> +
> +        if (y == null) {
> +            throw new IllegalArgumentException("y must not be null");
> +        }
> +
> +        if (x.length == 0) {
> +            throw new IllegalArgumentException(
> +                    "x must contain at least one element");
> +        }
> +
> +        if (y.length == 0) {
> +            throw new IllegalArgumentException(
> +                    "y must contain at least one element");
> +        }
> +    }
> +
> +    private double[] concatinateSamples(final double[] x, final double[] y) {
> +        final double[] z = new double[x.length + y.length];
> +
> +        System.arraycopy(x, 0, z, 0, x.length);
> +        System.arraycopy(y, 0, z, x.length, y.length);
> +
> +        return z;
> +    }
> +
> +    /**
> +     * {@inheritDoc}
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @return mannWhitneyU statistic U (maximum of U<sup>x</sup> and
U<sup>y</sup>)
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     */
> +    public double mannWhitneyU(final double[] x, final double[] y)
> +            throws IllegalArgumentException {
> +
> +        ensureDataConformance(x, y);
> +
> +        final double[] z = concatinateSamples(x, y);
> +        final double[] ranks = naturalRanking.rank(z);
> +
> +        double sumRankX = 0;
> +
> +        /*
> +         * The ranks for x is in the first x.length entries in ranks because x
> +         * is in the first x.length entries in z
> +         */
> +        for (int i = 0; i < x.length; ++i) {
> +            sumRankX += ranks[i];
> +        }
> +
> +        /*
> +         * U1 = R1 - (n1 * (n1 + 1)) / 2 where R1 is sum of ranks for sample 1,
> +         * e.g. x, n1 is the number of observations in sample 1.
> +         */
> +        final double U1 = sumRankX - (x.length * (x.length + 1)) / 2;
> +
> +        /*
> +         * It can be shown that U1 + U2 = n1 * n2
> +         */
> +        final double U2 = x.length * y.length - U1;
> +
> +        return FastMath.max(U1, U2);
> +    }
> +
> +    /**
> +     * @param Umin
> +     *            smallest Mann-Whitney U value
> +     * @param N
> +     *            number of subjects (corresponding to x.length)
> +     * @return two-sided asymptotic p-value
> +     * @throws MathException
> +     *             if an error occurs computing the p-value
> +     */
> +    private double calculateAsymptoticPValue(final double Umin, final int n1,
> +            final int n2) throws MathException {
> +
> +        final int n1n2prod = n1 * n2;
> +
> +        // http://en.wikipedia.org/wiki/Mann%E2%80%93Whitney_U#Normal_approximation
> +        final double EU = (double) n1n2prod / 2.0;
> +        final double VarU = (double) (n1n2prod * (n1 + n2 + 1)) / 12.0;
> +
> +        final double z = (Umin - EU) / FastMath.sqrt(VarU);
> +
> +        final NormalDistributionImpl standardNormal = new NormalDistributionImpl(
> +                0, 1);
> +
> +        return 2 * standardNormal.cumulativeProbability(z);
> +    }
> +
> +    /**
> +     * Ties give rise to biased variance at the moment. See e.g. <a
> +     * href="http://mlsc.lboro.ac.uk/resources/statistics/Mannwhitney.pdf"
> +     * >http://mlsc.lboro.ac.uk/resources/statistics/Mannwhitney.pdf</a>.
> +     *
> +     * {@inheritDoc}
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @param exactPValue
> +     *            if the exact p-value is wanted (only for x.length <= 50)
> +     * @return asymptotic p-value (biased for samples with ties)
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     * @throws MathException
> +     *             if an error occurs computing the p-value
> +     */
> +    public double mannWhitneyUTest(final double[] x, final double[] y)
> +            throws IllegalArgumentException, MathException {
> +
> +        ensureDataConformance(x, y);
> +
> +        final double Umax = mannWhitneyU(x, y);
> +
> +        /*
> +         * It can be shown that U1 + U2 = n1 * n2
> +         */
> +        final double Umin = x.length * y.length - Umax;
> +
> +        return calculateAsymptoticPValue(Umin, x.length, y.length);
> +    }
> +}
>
> Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTest.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTest.java?rev=1053836&view=auto
> ==============================================================================
> --- commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTest.java
(added)
> +++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTest.java
Thu Dec 30 09:52:00 2010
> @@ -0,0 +1,107 @@
> +/*
> + * Licensed to the Apache Software Foundation (ASF) under one or more
> + * contributor license agreements.  See the NOTICE file distributed with
> + * this work for additional information regarding copyright ownership.
> + * The ASF licenses this file to You under the Apache License, Version 2.0
> + * (the "License"); you may not use this file except in compliance with
> + * the License.  You may obtain a copy of the License at
> + *
> + *      http://www.apache.org/licenses/LICENSE-2.0
> + *
> + * Unless required by applicable law or agreed to in writing, software
> + * distributed under the License is distributed on an "AS IS" BASIS,
> + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
> + * See the License for the specific language governing permissions and
> + * limitations under the License.
> + */
> +package org.apache.commons.math.stat.inference;
> +
> +import org.apache.commons.math.MathException;
> +
> +/**
> + * An interface for Wilcoxon signed-rank test.
> + *
> + * @version $Revision: $ $Date: $
> + */
> +public interface WilcoxonSignedRankTest {
> +
> +    /**
> +     * Computes the <a
> +     * href="http://en.wikipedia.org/wiki/Wilcoxon_signed-rank_test">
> +     * Wilcoxon signed ranked statistic</a> comparing mean for two related
> +     * samples or repeated measurements on a single sample.
> +     * <p>
> +     * This statistic can be used to perform a Wilcoxon signed ranked test
> +     * evaluating the null hypothesis that the two related samples or repeated
> +     * measurements on a single sample has equal mean.
> +     * </p>
> +     * <p>
> +     * Let X<sub>i</sub> denote the i'th individual of the first sample
and
> +     * Y<sub>i</sub> the related i'th individual in the second sample.
Let
> +     * Z<sub>i</sub> = Y<sub>i</sub> - X<sub>i</sub>.
> +     * </p>
> +     * <p>
> +     * <strong>Preconditions</strong>:
> +     * <ul>
> +     * <li>The differences Z<sub>i</sub> must be independent.</li>
> +     * <li>Each Z<sub>i</sub> comes from a continuous population
(they must be
> +     * identical) and is symmetric about a common median.</li>
> +     * <li>The values that X<sub>i</sub> and Y<sub>i</sub>
represent are
> +     * ordered, so the comparisons greater than, less than, and equal to are
> +     * meaningful.</li>
> +     * </ul>
> +     * </p>
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @return wilcoxonSignedRank statistic
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     */
> +    double wilcoxonSignedRank(final double[] x, final double[] y)
> +            throws IllegalArgumentException;
> +
> +    /**
> +     * Returns the <i>observed significance level</i>, or <a href=
> +     * "http://www.cas.lancs.ac.uk/glossary_v1.1/hyptest.html#pvalue">
> +     * p-value</a>, associated with a <a
> +     * href="http://en.wikipedia.org/wiki/Wilcoxon_signed-rank_test">
> +     * Wilcoxon signed ranked statistic</a> comparing mean for two related
> +     * samples or repeated measurements on a single sample.
> +     * <p>
> +     * Let X<sub>i</sub> denote the i'th individual of the first sample
and
> +     * Y<sub>i</sub> the related i'th individual in the second sample.
Let
> +     * Z<sub>i</sub> = Y<sub>i</sub> - X<sub>i</sub>.
> +     * </p>
> +     * <p>
> +     * <strong>Preconditions</strong>:
> +     * <ul>
> +     * <li>The differences Z<sub>i</sub> must be independent.</li>
> +     * <li>Each Z<sub>i</sub> comes from a continuous population
(they must be
> +     * identical) and is symmetric about a common median.</li>
> +     * <li>The values that X<sub>i</sub> and Y<sub>i</sub>
represent are
> +     * ordered, so the comparisons greater than, less than, and equal to are
> +     * meaningful.</li>
> +     * </ul>
> +     * </p>
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @param exactPValue
> +     *            if the exact p-value is wanted (only works for x.length <=
30,
> +     *            if true and x.length > 30, this is ignored because
> +     *            calculations may take too long)
> +     * @return p-value
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     * @throws MathException
> +     *             if an error occurs computing the p-value
> +     */
> +    double wilcoxonSignedRankTest(final double[] x, final double[] y,
> +            boolean exactPValue) throws IllegalArgumentException,
> +            MathException;
> +}
>
> Added: commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestImpl.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestImpl.java?rev=1053836&view=auto
> ==============================================================================
> --- commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestImpl.java
(added)
> +++ commons/proper/math/trunk/src/main/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestImpl.java
Thu Dec 30 09:52:00 2010
> @@ -0,0 +1,271 @@
> +/*
> + * Licensed to the Apache Software Foundation (ASF) under one or more
> + * contributor license agreements.  See the NOTICE file distributed with
> + * this work for additional information regarding copyright ownership.
> + * The ASF licenses this file to You under the Apache License, Version 2.0
> + * (the "License"); you may not use this file except in compliance with
> + * the License.  You may obtain a copy of the License at
> + *
> + *      http://www.apache.org/licenses/LICENSE-2.0
> + *
> + * Unless required by applicable law or agreed to in writing, software
> + * distributed under the License is distributed on an "AS IS" BASIS,
> + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
> + * See the License for the specific language governing permissions and
> + * limitations under the License.
> + */
> +package org.apache.commons.math.stat.inference;
> +
> +import org.apache.commons.math.MathException;
> +import org.apache.commons.math.distribution.NormalDistributionImpl;
> +import org.apache.commons.math.stat.ranking.NaNStrategy;
> +import org.apache.commons.math.stat.ranking.NaturalRanking;
> +import org.apache.commons.math.stat.ranking.TiesStrategy;
> +import org.apache.commons.math.util.FastMath;
> +
> +/**
> + * An implementation of the Wilcoxon signed-rank test.
> + *
> + * @version $Revision: $ $Date: $
> + */
> +public class WilcoxonSignedRankTestImpl implements WilcoxonSignedRankTest {
> +    private NaturalRanking naturalRanking;
> +
> +    /**
> +     * Create a test instance where NaN's are left in place and ties get
> +     * the average of applicable ranks. Use this unless you are very sure
> +     * of what you are doing.
> +     */
> +    public WilcoxonSignedRankTestImpl() {
> +        naturalRanking = new NaturalRanking(NaNStrategy.FIXED,
> +                TiesStrategy.AVERAGE);
> +    }
> +
> +    /**
> +     * Create a test instance using the given strategies for NaN's and ties.
> +     * Only use this if you are sure of what you are doing.
> +     *
> +     * @param nanStrategy
> +     *            specifies the strategy that should be used for Double.NaN's
> +     * @param tiesStrategy
> +     *            specifies the strategy that should be used for ties
> +     */
> +    public WilcoxonSignedRankTestImpl(NaNStrategy nanStrategy,
> +            TiesStrategy tiesStrategy) {
> +        naturalRanking = new NaturalRanking(nanStrategy, tiesStrategy);
> +    }
> +
> +    /**
> +     * Ensures that the provided arrays fulfills the assumptions.
> +     *
> +     * @param x
> +     * @param y
> +     * @throws IllegalArgumentException
> +     *             if assumptions are not met
> +     */
> +    private void ensureDataConformance(final double[] x, final double[] y)
> +            throws IllegalArgumentException {
> +        if (x == null) {
> +            throw new IllegalArgumentException("x must not be null");
> +        }
> +
> +        if (y == null) {
> +            throw new IllegalArgumentException("y must not be null");
> +        }
> +
> +        if (x.length != y.length) {
> +            throw new IllegalArgumentException(
> +                    "x and y must contain the same number of elements");
> +        }
> +
> +        if (x.length == 0) {
> +            throw new IllegalArgumentException(
> +                    "x and y must contain at least one element");
> +        }
> +    }
> +
> +    /**
> +     * Calculates y[i] - x[i] for all i
> +     *
> +     * @param x
> +     * @param y
> +     * @throws IllegalArgumentException
> +     *             if assumptions are not met
> +     */
> +    private double[] calculateDifferences(final double[] x, final double[] y)
> +            throws IllegalArgumentException {
> +
> +        final double[] z = new double[x.length];
> +
> +        for (int i = 0; i < x.length; ++i) {
> +            z[i] = y[i] - x[i];
> +        }
> +
> +        return z;
> +    }
> +
> +    /**
> +     * Calculates |z[i]| for all i
> +     *
> +     * @param z
> +     * @throws IllegalArgumentException
> +     *             if assumptions are not met
> +     */
> +    private double[] calculateAbsoluteDifferences(final double[] z)
> +            throws IllegalArgumentException {
> +        if (z == null) {
> +            throw new IllegalArgumentException("z must not be null");
> +        }
> +
> +        if (z.length == 0) {
> +            throw new IllegalArgumentException(
> +                    "z must contain at least one element");
> +        }
> +
> +        final double[] zAbs = new double[z.length];
> +
> +        for (int i = 0; i < z.length; ++i) {
> +            zAbs[i] = FastMath.abs(z[i]);
> +        }
> +
> +        return zAbs;
> +    }
> +
> +    /**
> +     * {@inheritDoc}
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @return wilcoxonSignedRank statistic (the larger of W+ and W-)
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met
> +     */
> +    public double wilcoxonSignedRank(final double[] x, final double[] y)
> +            throws IllegalArgumentException {
> +
> +        ensureDataConformance(x, y);
> +
> +        // throws IllegalArgumentException if x and y are not correctly
> +        // specified
> +        final double[] z = calculateDifferences(x, y);
> +        final double[] zAbs = calculateAbsoluteDifferences(z);
> +
> +        final double[] ranks = naturalRanking.rank(zAbs);
> +
> +        double Wplus = 0;
> +
> +        for (int i = 0; i < z.length; ++i) {
> +            if (z[i] > 0) {
> +                Wplus += ranks[i];
> +            }
> +        }
> +
> +        final int N = x.length;
> +        final double Wminus = (((double) (N * (N + 1))) / 2.0) - Wplus;
> +
> +        return FastMath.max(Wplus, Wminus);
> +    }
> +
> +    /**
> +     * Algorithm inspired by
> +     * http://www.fon.hum.uva.nl/Service/Statistics/Signed_Rank_Algorihms.html#C
> +     * by Rob van Son, Institute of Phonetic Sciences & IFOTT,
> +     * University of Amsterdam
> +     *
> +     * @param Wmax largest Wilcoxon signed rank value
> +     * @param N number of subjects (corresponding to x.length)
> +     * @return two-sided exact p-value
> +     */
> +    private double calculateExactPValue(final double Wmax, final int N) {
> +
> +        // Total number of outcomes (equal to 2^N but a lot faster)
> +        final int m = 1 << N;
> +
> +        int largerRankSums = 0;
> +
> +        for (int i = 0; i < m; ++i) {
> +            int rankSum = 0;
> +
> +            // Generate all possible rank sums
> +            for (int j = 0; j < N; ++j) {
> +
> +                // (i >> j) & 1 extract i's j-th bit from the right
> +                if (((i >> j) & 1) == 1) {
> +                    rankSum += j + 1;
> +                }
> +            }
> +
> +            if (rankSum >= Wmax) {
> +                ++largerRankSums;
> +            }
> +        }
> +
> +        /*
> +         * largerRankSums / m gives the one-sided p-value, so it's multiplied
> +         * with 2 to get the two-sided p-value
> +         */
> +        return 2 * ((double) largerRankSums) / ((double) m);
> +    }
> +
> +    /**
> +     * @param Wmin smallest Wilcoxon signed rank value
> +     * @param N number of subjects (corresponding to x.length)
> +     * @return two-sided asymptotic p-value
> +     * @throws MathException if an error occurs computing the p-value
> +     */
> +    private double calculateAsymptoticPValue(final double Wmin, final int N) throws
MathException {
> +
> +        final double ES = (double) (N * (N + 1)) / 4.0;
> +
> +        /* Same as (but saves computations):
> +         * final double VarW = ((double) (N * (N + 1) * (2*N + 1))) / 24;
> +         */
> +        final double VarS = ES * ((double) (2 * N + 1) / 6.0);
> +
> +        // - 0.5 is a continuity correction
> +        final double z = (Wmin - ES - 0.5) / FastMath.sqrt(VarS);
> +
> +        final NormalDistributionImpl standardNormal = new NormalDistributionImpl(0,
1);
> +
> +        return 2*standardNormal.cumulativeProbability(z);
> +    }
> +
> +    /**
> +     * {@inheritDoc}
> +     *
> +     * @param x
> +     *            the first sample
> +     * @param y
> +     *            the second sample
> +     * @param exactPValue
> +     *            if the exact p-value is wanted (only for x.length <= 30)
> +     * @return p-value
> +     * @throws IllegalArgumentException
> +     *             if preconditions are not met or exact p-value is wanted for
> +     *             when x.length > 30
> +     * @throws MathException
> +     *             if an error occurs computing the p-value
> +     */
> +    public double wilcoxonSignedRankTest(final double[] x, final double[] y,
> +            boolean exactPValue) throws IllegalArgumentException,
> +            MathException {
> +
> +        ensureDataConformance(x, y);
> +
> +        final int N = x.length;
> +        final double Wmax = wilcoxonSignedRank(x, y);
> +
> +        if (exactPValue && N > 30) {
> +            throw new IllegalArgumentException("Exact test can only be made for
N <= 30.");
> +        }
> +
> +        if (exactPValue) {
> +            return calculateExactPValue(Wmax, N);
> +        } else {
> +            final double Wmin = ( (double)(N*(N+1)) / 2.0 ) - Wmax;
> +            return calculateAsymptoticPValue(Wmin, N);
> +        }
> +    }
> +}
>
> Added: commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/MannWhitneyUTestTest.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/MannWhitneyUTestTest.java?rev=1053836&view=auto
> ==============================================================================
> --- commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/MannWhitneyUTestTest.java
(added)
> +++ commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/MannWhitneyUTestTest.java
Thu Dec 30 09:52:00 2010
> @@ -0,0 +1,101 @@
> +/*
> + * Licensed to the Apache Software Foundation (ASF) under one or more
> + * contributor license agreements.  See the NOTICE file distributed with
> + * this work for additional information regarding copyright ownership.
> + * The ASF licenses this file to You under the Apache License, Version 2.0
> + * (the "License"); you may not use this file except in compliance with
> + * the License.  You may obtain a copy of the License at
> + *
> + *      http://www.apache.org/licenses/LICENSE-2.0
> + *
> + * Unless required by applicable law or agreed to in writing, software
> + * distributed under the License is distributed on an "AS IS" BASIS,
> + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
> + * See the License for the specific language governing permissions and
> + * limitations under the License.
> + */
> +package org.apache.commons.math.stat.inference;
> +
> +import junit.framework.TestCase;
> +
> +/**
> + * Test cases for the ChiSquareTestImpl class.
> + *
> + * @version $Revision: $ $Date: $
> + */
> +
> +public class MannWhitneyUTestTest extends TestCase {
> +
> +    protected MannWhitneyUTest testStatistic = new MannWhitneyUTestImpl();
> +
> +    public MannWhitneyUTestTest(String name) {
> +        super(name);
> +    }
> +
> +    public void testMannWhitneyUSimple() throws Exception {
> +        /* Target values computed using R version 2.11.1
> +         * x <- c(19, 22, 16, 29, 24)
> +         * y <- c(20, 11, 17, 12)
> +         * wilcox.test(x, y, alternative = "two.sided", mu = 0, paired = FALSE,
exact = FALSE, correct = FALSE)
> +         * W = 17, p-value = 0.08641
> +         */
> +        final double x[] = {19, 22, 16, 29, 24};
> +        final double y[] = {20, 11, 17, 12};
> +
> +        assertEquals(17, testStatistic.mannWhitneyU(x, y), 1e-10);
> +        assertEquals(0.08641, testStatistic.mannWhitneyUTest(x, y), 1e-5);
> +    }
> +
> +
> +    public void testMannWhitneyUInputValidation() throws Exception {
> +        /* Samples must be present, i.e. length > 0
> +         */
> +        try {
> +            testStatistic.mannWhitneyUTest(new double[] { }, new double[] { 1.0
});
> +            fail("x does not contain samples (exact), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.mannWhitneyUTest(new double[] { 1.0 }, new double[]
{ });
> +            fail("y does not contain samples (exact), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        /*
> +         * x and y is null
> +         */
> +        try {
> +            testStatistic.mannWhitneyUTest(null, null);
> +            fail("x and y is null (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.mannWhitneyUTest(null, null);
> +            fail("x and y is null (asymptotic), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        /*
> +         * x or y is null
> +         */
> +        try {
> +            testStatistic.mannWhitneyUTest(null, new double[] { 1.0 });
> +            fail("x is null (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.mannWhitneyUTest(new double[] { 1.0 }, null);
> +            fail("y is null (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +    }
> +}
>
> Added: commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestTest.java
> URL: http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestTest.java?rev=1053836&view=auto
> ==============================================================================
> --- commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestTest.java
(added)
> +++ commons/proper/math/trunk/src/test/java/org/apache/commons/math/stat/inference/WilcoxonSignedRankTestTest.java
Thu Dec 30 09:52:00 2010
> @@ -0,0 +1,180 @@
> +/*
> + * Licensed to the Apache Software Foundation (ASF) under one or more
> + * contributor license agreements.  See the NOTICE file distributed with
> + * this work for additional information regarding copyright ownership.
> + * The ASF licenses this file to You under the Apache License, Version 2.0
> + * (the "License"); you may not use this file except in compliance with
> + * the License.  You may obtain a copy of the License at
> + *
> + *      http://www.apache.org/licenses/LICENSE-2.0
> + *
> + * Unless required by applicable law or agreed to in writing, software
> + * distributed under the License is distributed on an "AS IS" BASIS,
> + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
> + * See the License for the specific language governing permissions and
> + * limitations under the License.
> + */
> +package org.apache.commons.math.stat.inference;
> +
> +import junit.framework.TestCase;
> +
> +/**
> + * Test cases for the ChiSquareTestImpl class.
> + *
> + * @version $Revision: $ $Date: $
> + */
> +
> +public class WilcoxonSignedRankTestTest extends TestCase {
> +
> +    protected WilcoxonSignedRankTest testStatistic = new WilcoxonSignedRankTestImpl();
> +
> +    public WilcoxonSignedRankTestTest(String name) {
> +        super(name);
> +    }
> +
> +    public void testWilcoxonSignedRankSimple() throws Exception {
> +        /* Target values computed using R version 2.11.1
> +         * x <- c(1.83, 0.50, 1.62, 2.48, 1.68, 1.88, 1.55, 3.06, 1.30)
> +         * y <- c(0.878, 0.647, 0.598, 2.05, 1.06, 1.29, 1.06, 3.14, 1.29)
> +         */
> +        final double x[] = {1.83, 0.50, 1.62, 2.48, 1.68, 1.88, 1.55, 3.06, 1.30};
> +        final double y[] = {0.878, 0.647, 0.598, 2.05, 1.06, 1.29, 1.06, 3.14, 1.29};
> +
> +        /* EXACT:
> +         * wilcox.test(x, y, alternative = "two.sided", mu = 0, paired = TRUE, exact
= TRUE, correct = FALSE)
> +         * V = 40, p-value = 0.03906
> +         *
> +         * Corresponds to the value obtained in R.
> +         */
> +        assertEquals(40, testStatistic.wilcoxonSignedRank(x, y), 1e-10);
> +        assertEquals(0.03906, testStatistic.wilcoxonSignedRankTest(x, y, true),
1e-5);
> +
> +        /* ASYMPTOTIC:
> +         * wilcox.test(x, y, alternative = "two.sided", mu = 0, paired = TRUE, exact
= FALSE, correct = FALSE)
> +         * V = 40, p-value = 0.03815
> +         *
> +         * This is not entirely the same due to different corrects,
> +         * e.g. http://mlsc.lboro.ac.uk/resources/statistics/wsrt.pdf
> +         * and src/library/stats/R/wilcox.test.R in the R source
> +         */
> +        assertEquals(40, testStatistic.wilcoxonSignedRank(x, y), 1e-10);
> +        assertEquals(0.0329693812, testStatistic.wilcoxonSignedRankTest(x, y, false),
1e-10);
> +    }
> +
> +    public void testWilcoxonSignedRankInputValidation() throws Exception {
> +        /*
> +         * Exact only for sample size <= 30
> +         */
> +        final double[] x1 = new double[30];
> +        final double[] x2 = new double[31];
> +        final double[] y1 = new double[30];
> +        final double[] y2 = new double[31];
> +        for (int i = 0; i < 30; ++i) {
> +            x1[i] = x2[i] = y1[i] = y2[i] = i;
> +        }
> +
> +        // Exactly 30 is okay
> +        testStatistic.wilcoxonSignedRankTest(x1, y1, true);
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(x2, y2, true);
> +            fail("More than 30 samples and exact chosen, IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        /* Samples must be present, i.e. length > 0
> +         */
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { }, new double[]
{ 1.0 }, true);
> +            fail("x does not contain samples (exact), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { }, new double[]
{ 1.0 }, false);
> +            fail("x does not contain samples (asymptotic), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { 1.0 }, new double[]
{ }, true);
> +            fail("y does not contain samples (exact), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { 1.0 }, new double[]
{ }, false);
> +            fail("y does not contain samples (asymptotic), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        /* Samples not same size, i.e. cannot be pairred
> +         */
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { 1.0, 2.0 }, new
double[] { 3.0 }, true);
> +            fail("x and y not same size (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { 1.0, 2.0 }, new
double[] { 3.0 }, false);
> +            fail("x and y not same size (asymptotic), IllegalArgumentException
expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        /*
> +         * x and y is null
> +         */
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(null, null, true);
> +            fail("x and y is null (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(null, null, false);
> +            fail("x and y is null (asymptotic), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        /*
> +         * x or y is null
> +         */
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(null, new double[] { 1.0 }, true);
> +            fail("x is null (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(null, new double[] { 1.0 }, false);
> +            fail("x is null (asymptotic), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { 1.0 }, null, true);
> +            fail("y is null (exact), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +
> +        try {
> +            testStatistic.wilcoxonSignedRankTest(new double[] { 1.0 }, null, false);
> +            fail("y is null (asymptotic), IllegalArgumentException expected");
> +        } catch (IllegalArgumentException ex) {
> +            // expected
> +        }
> +    }
> +}
>
>
>

---------------------------------------------------------------------
To unsubscribe, e-mail: dev-unsubscribe@commons.apache.org
For additional commands, e-mail: dev-help@commons.apache.org


Mime
View raw message