001/* 002 * Copyright (C) 2012 The Guava Authors 003 * 004 * Licensed under the Apache License, Version 2.0 (the "License"); you may not use this file except 005 * in compliance with the License. You may obtain a copy of the License at 006 * 007 * http://www.apache.org/licenses/LICENSE-2.0 008 * 009 * Unless required by applicable law or agreed to in writing, software distributed under the License 010 * is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express 011 * or implied. See the License for the specific language governing permissions and limitations under 012 * the License. 013 */ 014 015package com.google.common.math; 016 017import static com.google.common.base.Preconditions.checkArgument; 018import static com.google.common.base.Preconditions.checkNotNull; 019import static com.google.common.base.Preconditions.checkState; 020import static java.lang.Double.NaN; 021import static java.lang.Double.doubleToLongBits; 022import static java.lang.Double.isNaN; 023 024import com.google.common.annotations.Beta; 025import com.google.common.annotations.GwtIncompatible; 026import com.google.common.base.MoreObjects; 027import com.google.common.base.Objects; 028import java.io.Serializable; 029import java.nio.ByteBuffer; 030import java.nio.ByteOrder; 031import org.checkerframework.checker.nullness.qual.Nullable; 032 033/** 034 * An immutable value object capturing some basic statistics about a collection of paired double 035 * values (e.g. points on a plane). Build instances with {@link PairedStatsAccumulator#snapshot}. 036 * 037 * @author Pete Gillin 038 * @since 20.0 039 */ 040@Beta 041@GwtIncompatible 042public final class PairedStats implements Serializable { 043 044 private final Stats xStats; 045 private final Stats yStats; 046 private final double sumOfProductsOfDeltas; 047 048 /** 049 * Internal constructor. Users should use {@link PairedStatsAccumulator#snapshot}. 050 * 051 * <p>To ensure that the created instance obeys its contract, the parameters should satisfy the 052 * following constraints. This is the callers responsibility and is not enforced here. 053 * 054 * <ul> 055 * <li>Both {@code xStats} and {@code yStats} must have the same {@code count}. 056 * <li>If that {@code count} is 1, {@code sumOfProductsOfDeltas} must be exactly 0.0. 057 * <li>If that {@code count} is more than 1, {@code sumOfProductsOfDeltas} must be finite. 058 * </ul> 059 */ 060 PairedStats(Stats xStats, Stats yStats, double sumOfProductsOfDeltas) { 061 this.xStats = xStats; 062 this.yStats = yStats; 063 this.sumOfProductsOfDeltas = sumOfProductsOfDeltas; 064 } 065 066 /** Returns the number of pairs in the dataset. */ 067 public long count() { 068 return xStats.count(); 069 } 070 071 /** Returns the statistics on the {@code x} values alone. */ 072 public Stats xStats() { 073 return xStats; 074 } 075 076 /** Returns the statistics on the {@code y} values alone. */ 077 public Stats yStats() { 078 return yStats; 079 } 080 081 /** 082 * Returns the population covariance of the values. The count must be non-zero. 083 * 084 * <p>This is guaranteed to return zero if the dataset contains a single pair of finite values. It 085 * is not guaranteed to return zero when the dataset consists of the same pair of values multiple 086 * times, due to numerical errors. 087 * 088 * <h3>Non-finite values</h3> 089 * 090 * <p>If the dataset contains any non-finite values ({@link Double#POSITIVE_INFINITY}, {@link 091 * Double#NEGATIVE_INFINITY}, or {@link Double#NaN}) then the result is {@link Double#NaN}. 092 * 093 * @throws IllegalStateException if the dataset is empty 094 */ 095 public double populationCovariance() { 096 checkState(count() != 0); 097 return sumOfProductsOfDeltas / count(); 098 } 099 100 /** 101 * Returns the sample covariance of the values. The count must be greater than one. 102 * 103 * <p>This is not guaranteed to return zero when the dataset consists of the same pair of values 104 * multiple times, due to numerical errors. 105 * 106 * <h3>Non-finite values</h3> 107 * 108 * <p>If the dataset contains any non-finite values ({@link Double#POSITIVE_INFINITY}, {@link 109 * Double#NEGATIVE_INFINITY}, or {@link Double#NaN}) then the result is {@link Double#NaN}. 110 * 111 * @throws IllegalStateException if the dataset is empty or contains a single pair of values 112 */ 113 public double sampleCovariance() { 114 checkState(count() > 1); 115 return sumOfProductsOfDeltas / (count() - 1); 116 } 117 118 /** 119 * Returns the <a href="http://mathworld.wolfram.com/CorrelationCoefficient.html">Pearson's or 120 * product-moment correlation coefficient</a> of the values. The count must greater than one, and 121 * the {@code x} and {@code y} values must both have non-zero population variance (i.e. {@code 122 * xStats().populationVariance() > 0.0 && yStats().populationVariance() > 0.0}). The result is not 123 * guaranteed to be exactly +/-1 even when the data are perfectly (anti-)correlated, due to 124 * numerical errors. However, it is guaranteed to be in the inclusive range [-1, +1]. 125 * 126 * <h3>Non-finite values</h3> 127 * 128 * <p>If the dataset contains any non-finite values ({@link Double#POSITIVE_INFINITY}, {@link 129 * Double#NEGATIVE_INFINITY}, or {@link Double#NaN}) then the result is {@link Double#NaN}. 130 * 131 * @throws IllegalStateException if the dataset is empty or contains a single pair of values, or 132 * either the {@code x} and {@code y} dataset has zero population variance 133 */ 134 public double pearsonsCorrelationCoefficient() { 135 checkState(count() > 1); 136 if (isNaN(sumOfProductsOfDeltas)) { 137 return NaN; 138 } 139 double xSumOfSquaresOfDeltas = xStats().sumOfSquaresOfDeltas(); 140 double ySumOfSquaresOfDeltas = yStats().sumOfSquaresOfDeltas(); 141 checkState(xSumOfSquaresOfDeltas > 0.0); 142 checkState(ySumOfSquaresOfDeltas > 0.0); 143 // The product of two positive numbers can be zero if the multiplication underflowed. We 144 // force a positive value by effectively rounding up to MIN_VALUE. 145 double productOfSumsOfSquaresOfDeltas = 146 ensurePositive(xSumOfSquaresOfDeltas * ySumOfSquaresOfDeltas); 147 return ensureInUnitRange(sumOfProductsOfDeltas / Math.sqrt(productOfSumsOfSquaresOfDeltas)); 148 } 149 150 /** 151 * Returns a linear transformation giving the best fit to the data according to <a 152 * href="http://mathworld.wolfram.com/LeastSquaresFitting.html">Ordinary Least Squares linear 153 * regression</a> of {@code y} as a function of {@code x}. The count must be greater than one, and 154 * either the {@code x} or {@code y} data must have a non-zero population variance (i.e. {@code 155 * xStats().populationVariance() > 0.0 || yStats().populationVariance() > 0.0}). The result is 156 * guaranteed to be horizontal if there is variance in the {@code x} data but not the {@code y} 157 * data, and vertical if there is variance in the {@code y} data but not the {@code x} data. 158 * 159 * <p>This fit minimizes the root-mean-square error in {@code y} as a function of {@code x}. This 160 * error is defined as the square root of the mean of the squares of the differences between the 161 * actual {@code y} values of the data and the values predicted by the fit for the {@code x} 162 * values (i.e. it is the square root of the mean of the squares of the vertical distances between 163 * the data points and the best fit line). For this fit, this error is a fraction {@code sqrt(1 - 164 * R*R)} of the population standard deviation of {@code y}, where {@code R} is the Pearson's 165 * correlation coefficient (as given by {@link #pearsonsCorrelationCoefficient()}). 166 * 167 * <p>The corresponding root-mean-square error in {@code x} as a function of {@code y} is a 168 * fraction {@code sqrt(1/(R*R) - 1)} of the population standard deviation of {@code x}. This fit 169 * does not normally minimize that error: to do that, you should swap the roles of {@code x} and 170 * {@code y}. 171 * 172 * <h3>Non-finite values</h3> 173 * 174 * <p>If the dataset contains any non-finite values ({@link Double#POSITIVE_INFINITY}, {@link 175 * Double#NEGATIVE_INFINITY}, or {@link Double#NaN}) then the result is {@link 176 * LinearTransformation#forNaN()}. 177 * 178 * @throws IllegalStateException if the dataset is empty or contains a single pair of values, or 179 * both the {@code x} and {@code y} dataset must have zero population variance 180 */ 181 public LinearTransformation leastSquaresFit() { 182 checkState(count() > 1); 183 if (isNaN(sumOfProductsOfDeltas)) { 184 return LinearTransformation.forNaN(); 185 } 186 double xSumOfSquaresOfDeltas = xStats.sumOfSquaresOfDeltas(); 187 if (xSumOfSquaresOfDeltas > 0.0) { 188 if (yStats.sumOfSquaresOfDeltas() > 0.0) { 189 return LinearTransformation.mapping(xStats.mean(), yStats.mean()) 190 .withSlope(sumOfProductsOfDeltas / xSumOfSquaresOfDeltas); 191 } else { 192 return LinearTransformation.horizontal(yStats.mean()); 193 } 194 } else { 195 checkState(yStats.sumOfSquaresOfDeltas() > 0.0); 196 return LinearTransformation.vertical(xStats.mean()); 197 } 198 } 199 200 /** 201 * {@inheritDoc} 202 * 203 * <p><b>Note:</b> This tests exact equality of the calculated statistics, including the floating 204 * point values. Two instances are guaranteed to be considered equal if one is copied from the 205 * other using {@code second = new PairedStatsAccumulator().addAll(first).snapshot()}, if both 206 * were obtained by calling {@code snapshot()} on the same {@link PairedStatsAccumulator} without 207 * adding any values in between the two calls, or if one is obtained from the other after 208 * round-tripping through java serialization. However, floating point rounding errors mean that it 209 * may be false for some instances where the statistics are mathematically equal, including 210 * instances constructed from the same values in a different order... or (in the general case) 211 * even in the same order. (It is guaranteed to return true for instances constructed from the 212 * same values in the same order if {@code strictfp} is in effect, or if the system architecture 213 * guarantees {@code strictfp}-like semantics.) 214 */ 215 @Override 216 public boolean equals(@Nullable Object obj) { 217 if (obj == null) { 218 return false; 219 } 220 if (getClass() != obj.getClass()) { 221 return false; 222 } 223 PairedStats other = (PairedStats) obj; 224 return xStats.equals(other.xStats) 225 && yStats.equals(other.yStats) 226 && doubleToLongBits(sumOfProductsOfDeltas) == doubleToLongBits(other.sumOfProductsOfDeltas); 227 } 228 229 /** 230 * {@inheritDoc} 231 * 232 * <p><b>Note:</b> This hash code is consistent with exact equality of the calculated statistics, 233 * including the floating point values. See the note on {@link #equals} for details. 234 */ 235 @Override 236 public int hashCode() { 237 return Objects.hashCode(xStats, yStats, sumOfProductsOfDeltas); 238 } 239 240 @Override 241 public String toString() { 242 if (count() > 0) { 243 return MoreObjects.toStringHelper(this) 244 .add("xStats", xStats) 245 .add("yStats", yStats) 246 .add("populationCovariance", populationCovariance()) 247 .toString(); 248 } else { 249 return MoreObjects.toStringHelper(this) 250 .add("xStats", xStats) 251 .add("yStats", yStats) 252 .toString(); 253 } 254 } 255 256 double sumOfProductsOfDeltas() { 257 return sumOfProductsOfDeltas; 258 } 259 260 private static double ensurePositive(double value) { 261 if (value > 0.0) { 262 return value; 263 } else { 264 return Double.MIN_VALUE; 265 } 266 } 267 268 private static double ensureInUnitRange(double value) { 269 if (value >= 1.0) { 270 return 1.0; 271 } 272 if (value <= -1.0) { 273 return -1.0; 274 } 275 return value; 276 } 277 278 // Serialization helpers 279 280 /** The size of byte array representation in bytes. */ 281 private static final int BYTES = Stats.BYTES * 2 + Double.SIZE / Byte.SIZE; 282 283 /** 284 * Gets a byte array representation of this instance. 285 * 286 * <p><b>Note:</b> No guarantees are made regarding stability of the representation between 287 * versions. 288 */ 289 public byte[] toByteArray() { 290 ByteBuffer buffer = ByteBuffer.allocate(BYTES).order(ByteOrder.LITTLE_ENDIAN); 291 xStats.writeTo(buffer); 292 yStats.writeTo(buffer); 293 buffer.putDouble(sumOfProductsOfDeltas); 294 return buffer.array(); 295 } 296 297 /** 298 * Creates a {@link PairedStats} instance from the given byte representation which was obtained by 299 * {@link #toByteArray}. 300 * 301 * <p><b>Note:</b> No guarantees are made regarding stability of the representation between 302 * versions. 303 */ 304 public static PairedStats fromByteArray(byte[] byteArray) { 305 checkNotNull(byteArray); 306 checkArgument( 307 byteArray.length == BYTES, 308 "Expected PairedStats.BYTES = %s, got %s", 309 BYTES, 310 byteArray.length); 311 ByteBuffer buffer = ByteBuffer.wrap(byteArray).order(ByteOrder.LITTLE_ENDIAN); 312 Stats xStats = Stats.readFrom(buffer); 313 Stats yStats = Stats.readFrom(buffer); 314 double sumOfProductsOfDeltas = buffer.getDouble(); 315 return new PairedStats(xStats, yStats, sumOfProductsOfDeltas); 316 } 317 318 private static final long serialVersionUID = 0; 319}