Repository: JuliaStats/HypothesisTests.jl
Pull request #145 (2019-01-26 06:52:45)

Add Hotelling's T^2 test
Imported from https://github.com/ararslan/MultivariateTests.jl, updated for the modern era (originally written for 0.6). The test values were validated against SAS when I wrote it.

ararslan: It matches SAS and examples from a textbook so I'd say let's go for it. If something is wrong, people will let us know.

Patch:
diff --git a/Project.toml b/Project.toml
index 1b50c99..9ca7863 100644
--- a/Project.toml
+++ b/Project.toml
@@ -12,7 +12,8 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 
 [extras]
+DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
-test = ["Test"]
+test = ["DelimitedFiles", "Test"]
diff --git a/docs/make.jl b/docs/make.jl
index 8011367..dd0e5cf 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -8,7 +8,8 @@ makedocs(
         "methods.md",
         "parametric.md",
         "nonparametric.md",
-        "time_series.md"
+        "time_series.md",
+        "multivariate.md",
     ]
 )
 
diff --git a/docs/src/index.md b/docs/src/index.md
index 8af6846..bdad7e1 100644
--- a/docs/src/index.md
+++ b/docs/src/index.md
@@ -3,6 +3,6 @@
 This package implements several hypothesis tests in Julia.
 
 ```@contents
-Pages = ["methods.md", "parametric.md", "nonparametric.md", "time_series.md"]
+Pages = ["methods.md", "parametric.md", "nonparametric.md", "time_series.md", "multivariate.md"]
 Depth = 2
 ```
diff --git a/docs/src/multivariate.md b/docs/src/multivariate.md
new file mode 100644
index 0000000..ce7a0a2
--- /dev/null
+++ b/docs/src/multivariate.md
@@ -0,0 +1,9 @@
+# Multivariate tests
+
+## Hotelling's ``T^2`` test
+
+```@docs
+OneSampleHotellingT2
+EqualCovHotellingT2
+UnequalCovHotellingT2
+```
diff --git a/src/common.jl b/src/common.jl
index 9cd5b30..67f73e9 100644
--- a/src/common.jl
+++ b/src/common.jl
@@ -57,3 +57,11 @@ function tiedrank_adj!(ord::AbstractVector, v::AbstractArray)
 end
 
 tiedrank_adj(v::AbstractArray) = tiedrank_adj!(Vector{Float64}(undef, length(v)), v)
+
+# Pool covariance matrices, overwriting the first
+function poolcov!(Sx::AbstractMatrix, nxm1::Int, Sy::AbstractMatrix, nym1::Int)
+    @inbounds for i = eachindex(Sx, Sy)
+        Sx[i] = (Sx[i] * nxm1 + Sy[i] * nym1) / (nxm1 + nym1)
+    end
+    Sx
+end
diff --git a/src/hotelling.jl b/src/hotelling.jl
new file mode 100644
index 0000000..e2ea2f1
--- /dev/null
+++ b/src/hotelling.jl
@@ -0,0 +1,165 @@
+# Hotelling T² test
+
+export OneSampleHotellingT2, EqualCovHotellingT2, UnequalCovHotellingT2
+
+abstract type HotellingT2Test <: HypothesisTest end
+
+default_tail(::HotellingT2Test) = :right
+pvalue(T::HotellingT2Test; tail=:right) = pvalue(FDist(dof(T)...), T.F, tail=tail)
+
+function show_params(io::IO, T::HotellingT2Test, indent="")
+    println(io, indent, "number of observations: ", nobs(T))
+    println(io, indent, "number of variables:    ", T.p)
+    println(io, indent, "T² statistic:           ", T.T²)
+    println(io, indent, "transformed statistic:  ", T.F)
+    println(io, indent, "degrees of freedom:     ", dof(T))
+    println(io, indent, "covariance matrix:")
+    Base.print_matrix(io, T.S, indent^2)
+end
+
+## Utility functions
+
+function checkdims(X::AbstractMatrix, Y::AbstractMatrix)
+    nx, p = size(X)
+    ny, q = size(Y)
+    p == q || throw(DimensionMismatch("Inconsistent number of variables: $p, $q"))
+    (nx > 0 && ny > 0) || throw(ArgumentError("Inputs must be non-empty"))
+    return (p, nx, ny)
+end
+
+# Helper function for computing A'B⁻¹A when B is a covariance matrix
+At_Binv_A(A::AbstractArray, B::AbstractArray) = A'*(B \ A)
+
+## One sample
+
+struct OneSampleHotellingT2 <: HotellingT2Test
+    T²::Real
+    F::Real
+    n::Int
+    p::Int
+    μ₀::Vector
+    x̄::Vector
+    S::Matrix
+end
+
+StatsBase.nobs(T::OneSampleHotellingT2) = T.n
+StatsBase.dof(T::OneSampleHotellingT2) = (T.p, T.n - T.p)
+
+"""
+    OneSampleHotellingT2(X::AbstractMatrix, μ₀=<zero vector>)
+
+Perform a one sample Hotelling's ``T^2`` test of the hypothesis that the vector of
+column means of `X` is equal to `μ₀`.
+"""
+function OneSampleHotellingT2(X::AbstractMatrix{T},
+                              μ₀::AbstractVector=fill(middle(zero(T)), size(X, 2))) where T
+    n, p = size(X)
+    p == length(μ₀) ||
+        throw(DimensionMismatch("Number of variables does not match number of means"))
+    n > 0 || throw(ArgumentError("The input must be non-empty"))
+    x̄ = vec(mean(X, dims=1))
+    S = cov(X)
+    T² = n * At_Binv_A(x̄ .- μ₀, S)
+    F = (n - p) * T² / (p * (n - 1))
+    return OneSampleHotellingT2(T², F, n, p, μ₀, x̄, S)
+end
+
+"""
+    OneSampleHotellingT2(X::AbstractMatrix, Y::AbstractMatrix, μ₀=<zero vector>)
+
+Perform a paired Hotelling's ``T^2`` test of the hypothesis that the vector of mean
+column differences between `X` and `Y` is equal to `μ₀`.
+"""
+function OneSampleHotellingT2(X::AbstractMatrix{T}, Y::AbstractMatrix{S},
+                              μ₀::AbstractVector=fill(middle(zero(T), zero(S)), size(X, 2))) where {T,S}
+    p, nx, ny = checkdims(X, Y)
+    nx == ny || throw(DimensionMismatch("Inconsistent number of observations: $nx, $ny"))
+    p == length(μ₀) ||
+        throw(DimensionMismatch("Number of variables does not match number of means"))
+    return OneSampleHotellingT2(X .- Y, μ₀)
+end
+
+testname(::OneSampleHotellingT2) = "One sample Hotelling's T² test"
+population_param_of_interest(T::OneSampleHotellingT2) = ("Mean vector", T.μ₀, T.x̄)
+
+## Two sample
+
+## Two sample: equal covariance
+
+struct EqualCovHotellingT2 <: HotellingT2Test
+    T²::Real
+    F::Real
+    nx::Int
+    ny::Int
+    p::Int
+    Δ::Vector
+    S::Matrix
+end
+
+"""
+    EqualCovHotellingT2(X::AbstractMatrix, Y::AbstractMatrix)
+
+Perform a two sample Hotelling's ``T^2`` test of the hypothesis that the difference in
+the mean vectors of `X` and `Y` is zero, assuming that `X` and `Y` have equal covariance
+matrices.
+"""
+function EqualCovHotellingT2(X::AbstractMatrix, Y::AbstractMatrix)
+    p, nx, ny = checkdims(X, Y)
+    Δ = vec(mean(X, dims=1) .- mean(Y, dims=1))
+    S = poolcov!(cov(X), nx - 1, cov(Y), ny - 1) .* (inv(nx) .+ inv(ny))
+    T² = At_Binv_A(Δ, S)
+    F = T² * (nx + ny - p - 1) / (p * (nx + ny - 2))
+    return EqualCovHotellingT2(T², F, nx, ny, p, Δ, S)
+end
+
+StatsBase.nobs(T::EqualCovHotellingT2) = (T.nx, T.ny)
+StatsBase.dof(T::EqualCovHotellingT2) = (T.p, T.nx + T.ny - T.p - 1)
+
+testname(::EqualCovHotellingT2) =
+    "Two sample Hotelling's T² test (equal covariance matrices)"
+population_param_of_interest(T::EqualCovHotellingT2) =
+    ("Difference in mean vectors", zeros(eltype(T.Δ), T.p), T.Δ)
+
+## Two sample: unequal covariance
+
+# Store the denominator degrees of freedom in the type, since the computation
+# is expensive and we don't want to redo it every time the user calls dof
+struct UnequalCovHotellingT2 <: HotellingT2Test
+    T²::Real
+    F::Real
+    nx::Int
+    ny::Int
+    p::Int
+    ν::Int
+    Δ::Vector
+    S::Matrix
+end
+
+"""
+    UnequalCovHotellingT2(X::AbstractMatrix, Y::AbstractMatrix)
+
+Perform a two sample Hotelling's ``T^2`` test of the hypothesis that the difference in
+the mean vectors of `X` and `Y` is zero, without assuming that `X` and `Y` have equal
+covariance matrices.
+"""
+function UnequalCovHotellingT2(X::AbstractMatrix, Y::AbstractMatrix)
+    p, nx, ny = checkdims(X, Y)
+    Δ = vec(mean(X, dims=1) .- mean(Y, dims=1))
+    Sx = cov(X) ./ nx
+    Sy = cov(Y) ./ ny
+    ST = Sx .+ Sy
+    T² = At_Binv_A(Δ, ST)
+    F = (nx + ny - p - 1) * T² / (p * (nx + ny - 2))
+    tmp = Symmetric(ST) \ Δ
+    iν = (dot(tmp, Sx * tmp) / T²)^2 / (nx - 1) + (dot(tmp, Sy * tmp) / T²)^2 / (ny - 1)
+    ν = trunc(Int, inv(iν))
+    return UnequalCovHotellingT2(T², F, nx, ny, p, ν, Δ, ST)
+end
+
+StatsBase.nobs(T::UnequalCovHotellingT2) = (T.nx, T.ny)
+StatsBase.dof(T::UnequalCovHotellingT2) = (T.p, T.ν)
+
+testname(::UnequalCovHotellingT2) =
+    "Two sample Hotelling's T² test (unequal covariance matrices)"
+population_param_of_interest(T::UnequalCovHotellingT2) =
+    ("Difference in mean vectors", zeros(eltype(T.Δ), T.p), T.Δ)


Test patch:
diff --git a/src/HypothesisTests.jl b/src/HypothesisTests.jl
index 44eba09..1e1d455 100644
--- a/src/HypothesisTests.jl
+++ b/src/HypothesisTests.jl
@@ -177,4 +177,5 @@ include("augmented_dickey_fuller.jl")
 include("jarque_bera.jl")
 include("durbin_watson.jl")
 include("permutation.jl")
+include("hotelling.jl")
 end
diff --git a/test/data/nutrient.txt b/test/data/nutrient.txt
new file mode 100644
index 0000000..16de1d2
--- /dev/null
+++ b/test/data/nutrient.txt
@@ -0,0 +1,737 @@
+  1  522.29 10.188  42.561   349.13  54.141
+  2  343.32  4.113  67.793   266.99  24.839
+  3  858.26 13.741  59.933   667.90 155.455
+  4  575.98 13.245  42.215   792.23 224.688
+  5 1927.50 18.919 111.316   740.27  80.961
+  6  607.58  6.800  45.785   165.68  13.050
+  7 1046.19 18.433 116.418  1119.59 158.986
+  8  181.21 12.762  64.156    78.49  26.942
+  9  327.08  8.693  49.161   568.38  49.977
+ 10  383.09 13.667 103.844  1029.20   8.404
+ 11 1227.58 16.810 107.698   623.32  37.063
+ 12  845.44  7.417  56.519   273.16  21.692
+ 13  460.52 12.375  52.126   204.02 139.515
+ 14  349.56 16.074  50.679   597.67  53.576
+ 15   74.46  6.838  77.878     6.48  69.714
+ 16  280.31 14.548  88.798   556.27 100.920
+ 17  422.36 21.638 148.713   418.03 263.039
+ 21  552.90 21.517  91.083  7621.54 172.804
+ 22  588.77  8.950  74.493   567.05   3.589
+ 23  254.65  4.801  28.299   194.07  34.455
+ 24  546.66 20.961 125.557   396.50  65.010
+ 25   86.11  3.348  27.069     0.00  15.133
+ 26   84.39  3.526  16.494    28.58  10.645
+ 27 1131.79 23.181 133.082   544.73 105.937
+ 28  232.62  2.918  28.563   195.81  14.585
+ 29  130.89  4.194  24.770   540.55  61.401
+ 30  795.17  6.474  38.373   395.05  27.459
+ 31 1070.89  9.775  42.505   965.58 113.389
+ 32  357.08  6.023  47.933  1176.82 219.791
+ 33 1262.89 11.751 121.262  1268.65 327.188
+ 34  529.73  6.150  31.428   367.52 173.030
+ 35  491.31 12.072  59.397   524.03  30.565
+ 36  632.10  7.914  43.823   158.26  16.655
+ 37  574.12  9.118  52.544   801.05  40.690
+ 38  414.88 22.005 134.802   418.10 164.954
+ 39  329.71  6.039  57.454   293.76  67.949
+ 47  985.79 24.438 131.299   706.54 234.272
+ 48  326.54  9.460  71.244   375.80 110.426
+ 49  497.68 16.652  41.467  1339.46  48.611
+ 50  403.37 18.091 116.169   201.03 135.351
+ 51  664.97  7.776  71.124   216.34  32.214
+ 52  149.77  6.013  56.819   114.38  23.596
+ 53  831.59 13.038  45.940   578.64  26.859
+ 54  506.87 11.219  88.101   813.67  65.340
+ 55  975.62 13.474  77.293  1016.01  37.206
+ 56  693.27 17.210  88.751   855.60 143.399
+ 57  703.91 17.545 175.164   130.33  50.645
+ 58  419.26 11.112  73.502   179.56 112.504
+ 59  165.79  5.368  33.946   200.28  27.096
+ 60  580.68 14.963 151.601   763.37  20.320
+ 61  266.10  8.233  28.146   150.36  14.341
+ 62  147.71  3.669  34.919   211.89  59.382
+ 63 1089.24 12.961  95.119   760.81  36.846
+ 64  140.66  3.821  29.177    55.32   4.558
+ 65  230.43 23.226  58.614 11261.03 104.315
+ 66  554.24  7.505  56.937  1485.65  88.893
+ 67   69.34  2.033  16.294    12.08   3.605
+ 68  205.56  9.769  98.640   209.61  58.889
+ 69  598.93  8.662  52.446   190.34  13.240
+ 70  253.11  3.176  25.933   104.47   6.836
+ 71  929.92  8.831  63.722   566.31  39.440
+ 72  675.65 15.307  57.685  2385.88  29.351
+ 73  117.18  3.413  19.727    19.38  12.480
+ 74  702.91  8.845  24.706  1008.96 158.215
+ 75 1073.44 20.065 127.713   777.40 144.730
+ 76  618.03  7.329  58.499   272.14  47.659
+ 77  630.08 14.054  84.008  1077.78  78.124
+ 78  401.65  4.853  31.660   308.97   8.355
+ 79  228.63  8.624  30.721   704.79   0.739
+ 80 1263.04 31.079  70.538  2703.62 180.383
+ 81  751.93  5.038  39.375   357.07  23.659
+ 82  877.63 10.553  75.774   753.73 233.362
+ 83 1403.64  6.341  67.204   442.47  10.493
+ 89  155.09  7.850  59.013   528.01  56.255
+ 90  269.72  9.761  45.294   578.92  12.152
+ 91  901.57 15.994  89.959   519.05  26.802
+ 92 1330.62 20.873 114.956   962.20  75.658
+ 98  736.18  8.753  58.770  3500.77  85.589
+ 99  286.23 10.471  34.227    26.34  16.317
+100 1028.47 22.381  73.545  4100.43  76.401
+101  678.47 11.309  63.641   890.96  16.444
+102  732.50 18.467  77.839  1222.24 301.499
+103  252.96 11.295  60.346   189.64  14.782
+104 1354.99 20.229  71.032  1143.18  80.419
+105  442.90  7.136  56.949   211.64  63.467
+106  324.87  8.722  50.468   234.72  72.610
+107  478.60  9.076  73.568   370.25  16.946
+108  889.35  5.447  44.462  1405.41  89.725
+109  487.72  9.598  38.651  1487.75  13.472
+110  627.99 11.916  63.813   450.83 111.607
+111 2072.56  7.214  91.488  1048.08  80.585
+112  239.78  8.507  59.680   301.83  70.687
+113 1043.12 11.758 108.370   631.20  98.783
+114  202.50  9.295  56.303   388.57   4.501
+115  369.43  8.742  24.237  1038.84 239.563
+119  568.65  9.012  30.910    96.12  32.026
+120 1334.58 21.438  94.075  2023.71 139.911
+121  475.10 10.932  83.035   219.59  32.740
+122  604.65 11.452  59.100   713.93  22.055
+123  449.60  9.979  43.877   458.67  16.596
+124  462.60  8.796  25.363   554.59  69.345
+125  559.75  7.185  38.871   533.11 172.841
+126  421.33 12.728  59.518  1381.50 106.216
+127  464.37 11.501  75.722   701.62  46.667
+128  710.70  7.519  44.908   328.25  28.965
+129  877.36 30.032  92.248  1886.60  24.682
+130  313.06  5.705  22.330   166.93  35.614
+131  250.82  8.688  53.577   227.49  29.526
+132 1067.69 11.516  93.199   586.39  95.171
+135  646.59  5.260  40.527  2410.92  47.187
+136 1301.09  4.109  54.479   378.89  30.361
+137  790.08 12.183  82.020  1147.65  45.483
+138  854.41 14.208  93.167   737.77 158.333
+139  318.57  7.003  38.986   246.42 155.471
+140  557.23 11.313  55.922   785.11 140.147
+141  152.08  4.799  34.256   783.90  83.358
+142  171.76  4.231  29.826   205.27  48.691
+143  799.95  5.906  45.523   596.39  58.665
+144  818.47  8.408  63.296   380.53  11.885
+145  548.72  7.864  48.538   537.65  46.064
+146  570.27  8.311  48.982   258.24  29.926
+147  530.57  9.266  55.630   347.69  56.884
+148  556.49 10.479  56.242   438.91  10.030
+149  321.52  6.843  65.080   148.97 218.006
+150  987.90 15.477 107.859   775.94  60.203
+151 1211.30  8.558  85.494   598.98  77.093
+152  579.31 10.196  64.391   300.13  33.274
+153  865.37  8.967  56.490   500.68 218.795
+154  768.31 15.748  86.159  1553.63 332.144
+155  337.82  9.139  47.638   314.55 314.953
+156  548.29 13.262  77.479   322.48  37.937
+157 1155.33 10.289  96.108   773.13 157.918
+158  470.16 13.040 102.323   699.25 120.230
+169  275.17  9.665  60.209   276.54  49.553
+170  165.22  3.800  30.699   232.77 144.951
+171 1101.17 16.411  72.207  1962.16 231.591
+172  261.76 12.533  78.137   515.21  22.300
+173  151.46  6.486  63.829   150.90 109.632
+174  955.82 10.210  51.690   455.03  19.011
+177  920.13  7.255  63.792   584.93 137.243
+178  260.31  2.937  34.589   122.08   8.327
+179 1193.43 24.361 108.140   635.10  65.531
+180  419.22  6.541  38.275   300.00 206.349
+181  807.79  9.011  66.522   464.46  42.797
+182  365.51  8.401  43.626  1749.54  62.071
+183  765.73 15.727  76.945  1008.31 209.337
+184  668.73 22.796 115.905   641.90 163.300
+185  284.05  7.194  49.006   286.38  19.142
+186  523.08 13.789  45.911   375.54  21.236
+187  159.85  5.123  25.912   117.71  43.609
+188  788.02  8.068  36.242   772.25  46.413
+189  678.79  7.476  55.729  1384.16  81.321
+190  958.35 15.437  98.701   752.97  15.030
+191  774.57  6.149  51.593  1350.19  87.282
+192  583.28 20.376  86.681   288.96 108.666
+193 1498.07 13.167  83.049   688.55 215.302
+194 1159.71 10.370 103.114   729.12  16.641
+199 1548.99 13.211 126.912   865.00  61.681
+200  689.02  7.213  88.373   258.87  61.368
+201  265.39  6.651  44.480   304.97  20.624
+202   69.52  2.115   9.591     3.96   0.775
+203 1642.99 14.250 173.639  2084.53 228.976
+204  469.82  7.339  40.354   353.53  26.370
+205  223.67  6.813  55.891   215.23  20.774
+206  543.12  9.456  34.839  1693.80 135.222
+207  875.10 19.039 109.533  1017.56 169.567
+208 1494.06 11.409  79.329   750.70  53.322
+209  467.96  5.950  28.249   188.96   2.446
+210  913.72 13.809  76.899  1079.92  68.302
+211  995.86 10.763  62.814  1318.09 145.113
+212  100.49  5.777  16.134   583.58  20.444
+213  270.82  9.568  72.462   401.20 142.010
+214  429.21  8.546  41.772   650.13 155.567
+215  417.14  6.151  36.837   100.27  62.940
+216 1075.88 11.040  44.366   978.00 229.689
+217  928.45 10.220  66.601  1213.55  70.096
+218  210.26  3.742  46.609  1921.06  51.285
+219  476.67 11.502  72.512   731.03 114.865
+220  822.13 16.062  91.673   912.40 121.049
+221  248.41 11.673  80.901  1155.94 192.167
+222  633.67 18.765  89.889   698.85  45.960
+223  589.53 16.564  96.035   441.83 128.170
+224  279.22  9.619  69.606   764.07  54.272
+225 1022.66  9.614  51.015   484.20 210.494
+226  486.86 11.675  59.398   334.46  80.731
+227  233.19  9.596  91.456   444.90  13.629
+228  574.18  6.514  25.175   284.65 100.295
+229 1117.33 14.316 112.210  1843.32 295.685
+230  780.66 10.606 114.030   540.24  35.461
+231  833.71 16.263 117.516   906.56  74.429
+232  759.17 16.540 150.012   869.56 117.763
+233 1184.02 18.001 130.253  2037.92 187.449
+234  386.37 18.684  82.192   448.12 230.552
+235 1866.24 42.166  81.637  4053.28 139.901
+236  573.34 15.776  86.978  6947.59 118.381
+237  360.64 10.005  45.260   297.23  15.904
+238  550.55  6.895  38.062   505.12 187.999
+239  914.80 12.906  75.585   656.30  61.892
+240 1012.64 18.110  74.178  2001.68 129.127
+241  650.29 19.487  63.444  2262.19 275.854
+242  562.15  8.926  35.479   543.11 103.833
+243  649.77 10.484  59.042  1347.39 145.166
+251  633.09  8.639  89.188   458.03 117.720
+252  416.25  7.973 108.883   297.49 113.779
+253  594.13 11.199  56.152   923.74 168.288
+254  153.59  4.877  53.263   126.42  44.356
+255  338.60  9.253  55.854   153.39  72.161
+256  445.62  2.833  30.111   235.89  22.015
+257  550.36  7.096  44.028   284.51  30.020
+258  594.72  5.917  66.588   394.52  23.135
+259  506.07  9.766  70.207   247.54  19.229
+260  148.56  4.895  24.587    46.30   3.094
+261  456.49  5.807  45.231   428.10  21.876
+262  927.41 11.294  74.209   506.15 137.930
+263 2032.65 20.075 171.933  1446.40  59.105
+264  814.44  4.652  41.975  1652.00  24.238
+265  417.42 10.101  52.613   385.47  15.306
+269  626.62  7.200  60.476   494.30  79.594
+270 1292.16 16.377 104.208   761.59  63.527
+271  613.75 12.322  75.184   846.97 234.990
+272  135.41  9.907  61.477   229.10  50.425
+273  178.24  9.670  52.457   261.10  29.566
+274  239.86 16.802  70.630   328.65  48.834
+275  569.05 16.737  66.308  7134.56  46.127
+276  540.63 10.293  79.052   406.81  65.094
+277  697.13 11.049  87.404   631.74  76.784
+278  214.22  7.225  27.447   120.76  17.191
+279  687.12  7.560  54.985   430.00   5.381
+280  385.31 16.826  88.831   520.79  67.459
+281  843.66 10.419  57.813   475.72 134.750
+282  386.14  8.869  55.029   289.38  68.044
+283  329.10 12.891  86.920   177.10  71.505
+284  290.04  7.460  45.753   314.55  44.304
+285  444.18  9.821  45.168   118.82  11.218
+286  916.75 11.938  61.887  1325.03  53.659
+287  584.16 17.944 109.303  1028.35 160.421
+288  360.69  6.115  21.880   465.85  13.393
+291  169.91  8.114  48.395   465.02  77.967
+292  486.24 12.979 126.142   314.50  24.337
+297  248.88 10.674  64.624   299.82  61.508
+298  577.38 10.541  93.247  1089.25  26.994
+299  711.06  9.030  54.869   336.57  47.492
+300  482.36  9.419  66.579   465.57  20.689
+301  255.51 10.163  21.653   623.29 300.359
+302  479.40 17.586  44.560   961.71  69.112
+303  719.12 14.373  45.229   350.06  14.494
+304  405.75  8.742  69.566   768.66 335.266
+305  603.80  9.687  58.025   774.62  66.770
+306  587.60 10.084  79.705  1605.45  44.104
+307 1263.44 10.162  82.003   374.27  17.132
+308  468.46 11.753  65.704  1235.12  35.792
+309  830.75  9.061  65.750  1322.00 242.297
+311  174.55  6.822  44.461   486.67  29.476
+312 1076.46 10.469  72.306   560.77  92.225
+313  997.59  7.243  55.101   742.48  15.579
+314 1059.54 13.864  85.078   683.86  72.116
+315 1164.95 12.745  82.654   472.99  11.635
+316  670.69 15.619  54.608   435.80  17.258
+317  817.78 27.980  75.203  1253.65  32.479
+318 1006.57 11.238  71.041   387.36  79.886
+319   99.19  4.592  41.378    73.21  46.687
+320  409.12  7.907  51.574   150.80   8.645
+321  458.61  5.925  31.610   169.11  16.553
+322  695.89  6.740  29.873   428.12 263.056
+323  437.36  9.918  26.055   702.05 143.770
+324 1051.12 20.518 122.557  1296.97 112.873
+325  452.46  5.918  53.444   174.15 118.896
+332  252.40  8.054  34.028  1360.71  11.015
+333  847.70 16.791  68.929  1473.92  75.588
+334  835.88  7.138  57.385   385.04  36.804
+335  947.06  5.198  53.146   512.80  23.418
+336  170.19  6.890  24.008   271.51  35.549
+337  595.39 10.205  46.306   883.23 202.941
+338  384.21 13.711  34.934   858.83  42.399
+339  557.09  8.235  47.803   179.83  42.380
+340 1008.06 13.818 120.149   403.02  48.428
+341 1065.27  8.188 121.029   818.59  30.118
+342  219.65  3.917  26.785  1730.16 104.708
+343  341.32 13.743 146.195   654.34  45.254
+344  451.58  4.735  38.307   265.68  91.066
+345  441.56 15.927  72.558   697.94  28.229
+346  653.97  8.750  52.695   157.20  37.630
+347  768.16 20.343  73.579   980.63 149.628
+348  507.82  3.974  38.344   200.55   2.450
+349  365.58  9.901  46.542   268.51 137.272
+350 1006.46  7.555  61.519   506.49  36.318
+351  566.99 12.728  94.574   297.89  63.138
+352  704.14  9.418  49.784   844.98 190.638
+353   88.29  2.723  17.907   323.90  13.466
+354  205.19  5.089  36.171    68.87   0.000
+355  771.19 22.190  86.994   982.75 161.672
+356  109.53  0.645   3.145    20.15   3.010
+357  281.45  5.715  30.939   274.96  29.273
+369  618.49 10.790  53.055   511.45 161.383
+370 1516.42 16.367 124.987  1183.26 161.340
+371  750.26  8.844  79.664   436.10  18.948
+372 1324.00 11.391  72.936   697.17 128.035
+373  787.28  7.086  41.884   987.68  44.833
+374  498.49  8.947  47.448  1042.67  34.022
+375  908.82 12.255  65.811   747.58  21.943
+376  534.40 11.200  66.974   258.70 127.453
+377  728.71  8.213  49.328   458.40  89.510
+378 1414.60 35.077 251.012  2240.57 283.232
+379 1015.38 14.434  59.383  2247.53  82.720
+380  865.08 22.585  89.235  1173.10  25.941
+381  732.80 16.072  27.789   943.23  45.720
+382 1486.75  8.196  70.799  3448.27  15.285
+383  681.09  6.516  39.920   258.27 111.464
+384  324.95  7.551  34.932   301.38  68.223
+385  312.69  5.843  33.385   553.29  19.010
+386 1516.08 19.303  97.230   839.56 178.042
+387  323.13  7.885  39.491   329.32 125.150
+388 1243.94 16.133 122.051   786.23 215.674
+389  764.49  7.183  48.018   340.27  87.327
+390 1642.75 30.548  86.613  2882.30 264.257
+391  431.04 12.331  73.713   341.44 162.025
+392  427.50  9.603  69.799   228.64  28.732
+393  401.10  5.134  25.979   249.16  12.548
+394  325.55 12.459  81.624  7937.76  21.738
+395  278.58 15.320  88.484   351.52  35.792
+396  845.98 14.511  81.977  1006.73 102.249
+397  947.66 19.399 115.936  1084.06 274.953
+398 1149.45 18.676  78.217   477.36 102.365
+399  244.57 11.629  38.805   618.54  38.735
+400  820.52  8.433  54.309   677.47 110.086
+401  473.23 14.134  81.916  1579.73 291.251
+402  508.25 14.222  67.970  1127.83   6.120
+403  162.92  5.461  18.424   165.28  97.219
+404  968.21 15.386  76.665  1023.40  20.837
+405  903.69 10.149  64.593   317.15  53.594
+406  937.14 17.823 109.122   910.57 110.730
+411 1473.85  4.575  48.200   455.20  38.977
+412  708.88 19.365  65.076   134.33  28.242
+413  451.13  5.282  55.110   307.95  29.993
+414  276.22  8.924  36.378   106.27   5.079
+415  274.99  8.488  59.360   191.85  22.749
+416  684.75 12.548  98.943   579.34  94.408
+417  303.78  5.778  40.731   283.80  26.489
+418  436.92  5.776  30.578   576.22  52.726
+419  538.82 19.367  61.070   307.29 112.821
+420  541.22 11.220  84.971   583.61 104.378
+421  443.43  6.855  44.750  1212.51  20.512
+422 2465.78  8.347  85.069   899.03  30.034
+423  956.57 12.073  81.467  1218.89  76.728
+424  544.64 13.400 100.618   526.52  13.124
+425  517.41 12.606  61.450   646.50  89.709
+426  449.19  6.645  44.515   206.72 175.868
+427  406.78  8.136  72.176   281.57  41.472
+428 1032.01  7.153  59.208   549.88  36.201
+429  947.54 24.677  93.018  1814.73  30.287
+430  858.67  7.775  66.604  3226.59 119.937
+431  789.27 17.244  86.418   812.10 269.757
+441  650.74  9.304  58.789   513.18  47.816
+442  897.90 58.668  72.829  2814.08  63.290
+443  545.14  8.865  48.670  1449.52 214.063
+444  893.80  9.263  65.254  5344.78  87.468
+445  938.60  0.367  25.578   433.65   7.350
+446  663.61 17.583  59.775  1337.15  95.386
+447  180.55  7.613  54.831    85.64  27.660
+448 1174.69  8.588  80.349   497.59 112.676
+449  348.43  5.484  89.307   277.34   1.838
+450  772.58 16.368  56.331   175.74  17.002
+451   62.59  1.773  15.901     2.27   4.495
+452  123.80  9.495  64.961   310.07  16.459
+453  449.07 10.178  69.824   725.93  59.100
+454  880.32 17.379 114.932   543.58 183.523
+455  236.66  9.476 112.514   537.21  41.638
+456  461.24 14.698  74.222   881.32 101.691
+457  587.12  8.345  61.072   794.43 183.122
+458  971.88 11.037  52.035  1589.83 256.408
+459  402.27 16.250 112.054   188.00  46.088
+461  386.77  4.814  47.000   277.54  50.926
+462  600.22 12.067  76.146   463.07  50.130
+463  420.74 15.297 138.823   353.06 199.276
+464  683.08 11.976  77.189   489.15  42.179
+465  209.43  9.494  73.298   290.18  63.103
+466  743.88 16.294  74.204   782.27 108.185
+467  332.04 11.568  92.977   585.31 140.124
+468  540.89  9.012  54.589  1467.10 238.763
+469  765.60  5.669  37.220  1489.29 106.154
+472 2866.44  8.457 115.734  1223.09  28.613
+473  374.46 14.337  83.377   496.55  89.432
+474  457.52  8.344  47.120   219.96  49.239
+475 1217.15 13.119  75.173   731.21  51.165
+476  816.75 15.194  95.733   375.93  36.042
+477  741.24 17.971  81.994  1325.20  77.275
+478  515.09 11.029  66.082   424.52  41.574
+479 1922.19 13.608 109.143  1374.98 211.107
+480  878.18 27.791  64.793  2318.95 311.100
+481 1189.61  7.611  55.131  2395.71 202.105
+482  480.45  9.764  45.557  1227.99 103.639
+483  608.03 14.288  90.952   741.97  28.160
+484  479.05 13.823  70.949   220.27  53.709
+485 1184.69 15.688 130.512   671.63  67.333
+486  512.72 41.354 150.625  1597.81  82.446
+487  754.71 10.494  96.384   543.59 142.185
+488  509.40 12.988  63.562   566.94 205.563
+489  611.13 17.636  67.592  2297.44 261.317
+490  533.58 10.974 111.063  1644.52  62.835
+491  228.87  6.751  46.503   140.95  11.630
+492  327.88 12.191  54.865   119.18   7.176
+493  172.72  5.881  57.954   252.64 101.139
+494  471.31 14.453  74.784   291.11  61.071
+495  564.83 18.689 163.442  2594.52  65.405
+496  100.01  5.876  15.257   107.34   8.320
+497  599.88 13.510  56.382   253.72  47.616
+498  667.87 11.619  80.774 10405.32 407.999
+499  228.36  7.238  51.296  1042.57  30.141
+500  278.55  3.850  25.096    20.10  23.552
+501  540.87  8.639  54.918   286.98  99.360
+502 1132.23 11.053  78.618  1314.56  80.919
+503   34.88  1.661  19.985    43.09   4.060
+504  383.59  4.815  37.733   419.43  88.450
+505  719.77 10.522  44.981   192.67  61.392
+519 1196.60 10.447  67.429   948.22  26.563
+520  626.21 15.259  83.776   255.41  24.942
+521  213.11  7.703  56.202   100.00  33.653
+522  787.00  8.692  76.505   563.13 161.141
+523   99.69  3.888  46.776   900.52  24.914
+524  728.82 11.629  73.267   511.28  36.793
+525  123.00  2.493  12.339   289.46  17.399
+526  599.69 12.932  66.696   354.25  64.768
+527  217.19  3.046  28.565   117.29   0.595
+528  235.70  7.553  50.196   236.27  12.874
+529  305.12  9.397  45.924   351.32  49.329
+530  236.46 11.326  73.112   234.33  44.213
+531  257.84  5.521  30.744   190.67  64.561
+533  857.13 12.337  52.620  1684.24  31.303
+534  686.08 14.926 117.282   683.42  64.835
+535  104.70  7.547  52.354    82.12  45.957
+536  400.62 12.513  43.301   498.47 168.978
+537  567.65 13.986  63.661   576.87  60.170
+545  944.50  9.642 103.542   831.69 155.485
+546  685.54 11.606  89.352   480.34  52.869
+547  846.25 19.692 156.897   892.52  19.181
+548  594.05 15.052  72.342   481.90  18.544
+549 1002.46 11.338 120.557   714.79  22.361
+550  619.39  8.505  36.565   246.43 260.671
+551  340.87 14.285  60.183   484.77  63.850
+552  880.27 12.512  62.486  1209.87  45.770
+553 1067.58 18.409  97.031  1277.55 100.782
+554  479.02 15.811  96.324   228.54  59.037
+555 1197.55 11.385  80.052   467.07 177.245
+556  488.93 13.365 100.102   931.42 106.409
+557  127.17  9.324  57.306   656.81  37.433
+558 1863.51 14.900  94.805   917.23  64.341
+559  281.26 11.476  80.394   614.50  54.303
+564  146.63  6.582  45.839   266.79  10.515
+565  636.22  8.420  40.671   174.50  17.984
+566  822.57  8.951  69.750  1070.88 166.254
+567  487.06  8.291  54.301   268.62  31.260
+568  517.68 22.203 140.320  1417.20  41.843
+569  621.91  2.882  28.238   631.58  29.716
+570  464.03 13.443  72.652    73.32   5.547
+571  213.43  7.897  53.109  2189.24  40.268
+572 1263.81 16.683 157.152   924.12  34.638
+573  313.44  8.469  58.588   269.70 414.802
+574 1415.72  9.629  78.613   987.31  14.332
+575  931.02 16.145 114.863  2176.71 232.231
+576  158.25  6.496  33.525   113.44   7.535
+577  508.79  6.498  30.418   612.02 121.792
+578  206.90 10.313  67.492   645.37  10.416
+579  824.29  5.387  71.710   506.37  18.254
+583  729.54 13.773  74.483   864.48  25.574
+584  397.36 10.883  73.267    94.65  42.141
+585 1234.79 18.214  76.503  1199.67  66.745
+586  589.15  9.628  64.466  1161.13  74.315
+587  249.23 12.467  59.446  1241.82 242.861
+588  231.53  9.736  62.768   442.45  53.573
+603  186.53  8.503  37.810   770.50  71.291
+604  595.08  7.542  44.353  1897.13 106.132
+605  336.10 13.836  67.867   576.08  93.716
+606 2111.83 13.637 141.407  1470.88 127.426
+607 1210.33 18.690  89.235  1029.11 197.375
+608  528.20  9.185  45.770   491.74  75.639
+609  550.04  6.336  57.761   889.11  45.934
+610  460.25  5.214  37.891   199.66   5.399
+611  644.49  9.089  54.807  2323.11 151.486
+612  993.95  9.011  74.015  1441.71 107.484
+613  631.12  9.149  76.233   752.33 173.391
+614  555.35  9.400  80.646   813.01  75.863
+615 1060.13 17.798  88.169   636.83  61.182
+616  487.37 15.267  49.771  1017.40 166.354
+617 1150.35  4.906  46.928   515.33  98.252
+618  286.03  8.266  37.834  2029.51  19.396
+619  235.90  7.321  64.098   840.41  68.845
+620  457.28  8.306  43.859  1362.81  20.665
+621 1376.88 15.758  80.707   975.79  47.132
+622  620.32  7.729  69.706   444.75  87.028
+623  949.86  8.208  55.713   938.21  19.321
+624  365.10 31.562  60.583   465.54 122.157
+625  897.43  5.091  43.390   432.41   4.744
+626  686.45  7.808  67.588   902.39 120.945
+627 1512.47 11.721  77.809   550.75  19.195
+628  202.86 10.551  53.690   484.13  38.452
+629  697.41  6.581  56.199   444.66  19.416
+630 1855.69 18.822  79.813  1137.06 199.147
+631  535.29 23.269 103.510   487.86 174.170
+633  498.81  7.293  28.168   244.68  15.584
+634 1290.96 16.930  98.897  1602.27 141.626
+635  345.41 12.606  62.432   219.70 129.331
+636  447.62 29.296  73.915   438.45  62.734
+637  433.74 11.769  51.607  2318.90 138.455
+638  364.21  9.687  57.027   351.17  22.317
+642  671.48 14.527  64.401  1237.32 144.169
+643  261.69 15.783  50.589  1366.99  97.999
+644  259.41  6.901  21.918   276.63 234.122
+645  308.64  7.967  47.096   269.79  16.533
+646  511.85  9.085  60.676  1035.28 141.570
+647  435.20 18.272  83.204  1090.14 224.956
+648  457.51 11.071  37.999   490.74  86.484
+649 1488.64 18.071  85.809  2252.22 100.520
+650  547.62 10.510  44.012   220.82  17.992
+651 1059.48 23.169  71.898   858.45 231.785
+652  817.89  3.529  45.412   602.88  29.553
+653  277.73  3.720  15.845   236.65  97.760
+654  282.64 10.634  53.594   156.02  15.999
+655  356.37 10.744  28.261   884.50   6.727
+656 1102.78  7.438  84.046  1107.42  24.195
+657  564.08  8.109  52.129  2239.96 116.086
+658  351.41 15.108  43.117   547.64 174.539
+659 1522.47 13.251 114.235  1374.74 168.930
+660  524.91 22.220  59.655   211.59  96.489
+664   97.53  4.550  20.384    58.48  24.002
+665  424.69  4.086  23.637    25.84  39.332
+666  687.29  5.559  41.379   430.59  16.356
+667  224.38  5.298  32.088  1792.18  50.884
+668  243.53 10.770  99.552  2233.47  79.535
+669  561.53 11.375  58.778   557.65 106.435
+670  794.07  6.260  44.526   506.90 157.319
+671  799.21  5.732  49.307  1316.13  20.109
+672 1509.41 15.171  99.245   541.44  71.521
+673  309.91  6.972  32.071   595.54   8.209
+674  252.80  9.233  45.632   250.49  49.069
+675  352.45  5.062  31.666   296.57  38.034
+676  325.95  8.926  52.104   949.73 333.161
+677  572.58  4.599  35.471   148.54 137.901
+678 1640.31  1.957  61.813   776.94  38.536
+679  538.69 16.623  78.774   707.83  89.039
+680  291.33  7.248  35.986  1516.67  22.704
+681  220.66  3.887  13.414   146.06  55.251
+682  841.10 10.857  62.530   940.99  28.507
+683  271.54  9.903  49.265   474.93  18.245
+684  194.44  7.430  80.657    78.27  70.525
+685 1370.42 18.881 115.729   675.39 318.132
+686  672.34 12.419  63.703  1915.80  34.413
+689  781.54  6.993  47.404  1316.73 102.383
+690  700.75 10.183  46.479   314.60  36.455
+691  659.02 16.799  84.772 34434.27 123.317
+692  291.52  7.710  37.824   527.50  65.786
+693  253.62  7.704  23.758     5.97  29.741
+694 1323.85 25.263  79.677  1079.41  64.339
+704  340.31 13.370  80.990   296.52  30.969
+705  617.94 26.150  79.519   330.76  23.300
+706 1373.01 14.304  89.944   894.43  28.275
+707  261.19 10.905  49.731   358.29  28.075
+708  267.70  4.859  23.125   189.87  18.952
+709  363.19 18.780 138.225  1043.11  34.634
+710  372.70  9.317  92.041   149.97 148.961
+711  501.58 13.565  79.505   666.59  74.661
+712 1255.55 17.395  95.566  1344.65  67.985
+713  238.48  9.484  37.078   112.40  74.659
+714 1282.07 55.323 229.617  2849.82 150.567
+715  619.64  9.936 119.922   449.03  70.098
+716  391.72 14.117  81.885  2300.43  43.829
+717  601.66 11.162  72.726   522.34  54.902
+718  367.43 10.941  51.935   124.17  18.061
+719  547.22  9.099  65.496   487.93  37.756
+720  732.02 12.391  70.122   424.53  52.668
+721 1972.22 12.796 125.967  1812.19  60.283
+722  800.64 16.147  83.801  4521.66 113.853
+723  183.13  4.399  25.596   297.55  48.197
+724  760.65  7.483  57.328   329.65  74.989
+737  699.31  6.681  46.669   175.31  41.318
+738  367.38  4.297  83.144  2276.14  22.067
+739  255.63 11.409  65.336    40.41 108.508
+740 1081.85 15.683 119.992  1462.74 102.017
+741  956.11 12.976  51.975   321.58  20.664
+742  242.43 10.237  37.238   284.60   0.037
+743  195.82  4.817  60.089   434.37  64.580
+744  226.19  7.903  30.640   239.04  29.072
+745  732.61  8.773  55.200   715.15 167.408
+746  310.69  3.880  12.056   352.36 210.842
+747 1033.37 18.866  73.027  1250.31 106.413
+748 1773.76 18.718 136.846  1609.34  53.722
+749  361.50 16.550  90.565   191.35   6.904
+750  873.97 11.145  96.631   697.53  43.896
+751 1246.59  9.662  71.457   827.41 135.447
+752  547.43 10.753  56.469   304.96  23.490
+753  233.23  8.093  51.051   191.09   2.371
+754  342.00 13.015  74.201   532.65   1.591
+755 1092.71  7.136  52.177  1509.30  32.311
+756  633.92  7.314  61.407   878.85 139.823
+757  411.87  7.111  39.830  2752.42  27.235
+758 1017.64 15.444  82.017   845.90  59.828
+759  310.13 12.089  45.024   922.08  38.449
+760 1611.92 18.191 107.108  1428.36 219.804
+761  675.33 10.322  55.490   417.03  16.090
+762  821.16 12.456  59.553   789.31 189.411
+763  140.90  6.960  23.711   105.61   3.360
+764 1051.02 12.746  68.216  1690.47  72.404
+765  552.62  4.914  46.890   242.47   8.924
+766  244.15  4.974  50.495   465.14  79.004
+767  139.86  5.462  32.528   148.33  17.308
+768  992.74 10.453  74.105   486.61 172.719
+769 1058.71  9.518  56.709   402.99  26.135
+771  101.33 10.725  85.218   150.98  17.147
+772  935.25 19.854 116.062   833.72 110.106
+773  203.64  7.812  42.246   126.80 102.859
+774  436.73 22.902  70.770  1360.94  78.877
+783 1732.46  8.865  69.129   680.41  22.411
+784  411.25  4.148  53.623    61.44   3.248
+785  731.00 13.915  65.335   971.74 232.932
+786  826.51 11.237  59.431   348.56  36.762
+787 1531.85 10.609  91.431   974.99  44.999
+788 1357.81  8.598  92.639   784.67  12.197
+789  609.27  8.015  45.522  1549.32 148.512
+790 1945.91 29.047  91.068  2513.04 198.842
+791  971.12  9.157  72.576   589.13  53.703
+792  316.81  9.054  30.588   230.12  21.729
+793 1108.78  9.300  59.468   575.19  75.768
+794 1379.65 12.434  73.026   543.02 126.612
+795  510.41 16.954  77.844   573.02  57.574
+796  623.93  6.511  52.873   276.27   7.474
+801  631.42 21.396 123.871   815.48  36.433
+802  985.94  9.741  58.867   667.98  42.270
+803  307.93  8.353  52.243   426.99  14.696
+804  750.83 12.645 101.535   578.96 112.144
+805  788.70  9.401  62.752   396.96  22.635
+806  512.21 11.258  45.282   901.78  99.683
+807   13.53  1.318   2.255     0.00   0.000
+808  493.53 10.662  56.605  1375.55  72.377
+814  726.56 10.325  75.775   446.48  36.889
+815  293.64  8.950 115.235   287.85  22.690
+816  570.90 15.999  60.150  1472.86 169.385
+817  617.85 14.392 151.284   330.78  66.900
+818  363.73  8.110  67.695  1809.34 218.513
+819  495.59 12.906 113.826   271.80 134.593
+820  120.61  1.699   2.724   208.86   7.388
+821  451.12  9.490  48.687   414.32 158.923
+822  971.41 13.767  60.609   892.35  41.256
+823  292.63  7.690  73.193   165.90  11.765
+824  608.13 11.616  80.352   246.84  25.775
+825  645.88 12.956  64.858   418.79  23.610
+826  110.92  4.828  24.064   205.59 217.388
+827 1065.18 11.456  80.007   600.52  42.723
+828 1300.44 19.703  96.373  2998.92 150.237
+829  582.17 10.439  48.180   118.96  26.386
+830  224.20  9.070 111.301   112.73   0.069
+831  462.54  9.082  56.896    19.35  46.925
+832  280.34  7.949  52.360   157.46  32.440
+833  646.71  7.891  64.771   485.89 360.579
+834  342.16 14.216 151.425   466.09  27.175
+835  878.67 10.457  92.405  1357.32  25.263
+836  850.46  8.539  85.120   664.67 106.954
+839   91.20  0.785   4.320     0.00   0.000
+840  459.11 14.874  66.964  1288.74  62.139
+841  180.29  3.364  20.517    55.93  15.667
+842  527.23 12.577  70.805   308.76  46.759
+843  671.90 11.098  81.530   922.64  94.855
+844  361.49 12.090  59.525   253.93  58.111
+845  266.63  9.139  70.261   650.66  49.119
+846  171.82  3.963  22.490   320.47  28.622
+847  124.19  3.060  23.367    29.56  23.470
+848   94.52  4.578  56.053   185.69  26.392
+849 1154.61  9.400  76.772   767.03  30.167
+850  929.06 16.052  70.141   299.32  41.399
+851  491.22  6.069  42.136   223.25  81.747
+852  857.78 13.512  74.782   396.93   8.365
+853  362.87 10.155  67.332   990.48  41.077
+854  417.02  8.857  80.715   430.22  22.120
+855  267.88  4.665  45.225   136.26  26.839
+856  293.38 11.072 111.507   293.34  36.759
+859  362.78  7.536  56.140   595.41  41.563
+860  642.42  8.773  43.377   219.42  53.466
+861  609.60 11.049  73.557   221.71  23.918
+862  656.71  7.833  52.561    99.07 131.008
+863  577.56 12.669  56.288   592.04  67.365
+864  368.41  7.962  38.008   742.46  21.467
+865 1412.43 13.747 118.721  1472.62  88.206
+866  209.54  6.353  56.540   251.96  19.270
+867  649.23  8.910  56.640   680.42  19.871
+868  439.59 14.081  85.010   228.44  22.560
+869  412.77 11.128  55.835   283.86  99.254
+870  704.07 13.895  72.930   637.09 148.321
+871  689.53 11.427  92.504   737.40  92.642
+881  224.41 11.393  61.745   164.82  14.560
+882  130.54  2.828  11.529   117.33  66.478
+883  500.56 12.241  67.546 12174.29  30.646
+884    7.44  0.000   0.000     0.00   0.000
+885 1192.92  7.119  65.628   470.33  14.876
+888  344.64 10.870  54.609   259.51  37.523
+889 1036.45 12.747  94.604  3731.88 251.547
+890  367.12  9.759  95.265   420.29  13.160
+891  297.93  6.906  34.484   162.19  78.284
+892  397.05 13.117  58.140   534.25  48.759
+893  446.47  8.107  25.008  2559.24  50.959
+894  291.02 11.612  65.328   244.80 116.371
+895  905.89  8.446  73.206   780.59  46.062
+896  317.87  5.533  29.837   455.88  36.205
+897  180.16  9.910  97.089    72.61   0.683
+898  308.93 10.384  58.482   267.32  35.136
+899  302.96  6.426  49.110   134.30  67.251
+900  697.12 38.454  39.274  3391.02 433.339
+901  895.43  8.599  56.645  1145.89  34.370
+902  380.78  6.378  40.420   128.93  14.636
+903  288.06  4.972  31.619   235.38  30.235
+904  638.43 15.322  78.394   601.27  59.407
+905  176.55 10.033  50.787   252.20 157.730
+906  592.78  7.157  48.977   431.96  99.534
+907  481.68  7.374  45.965   166.91  22.532
+908  546.49 11.096  43.835   634.21  78.979
+909  458.21 10.580  75.056   301.22  36.023
+917  991.76 25.060 117.567  1284.47  80.684
+918  601.67 13.579  63.903   697.36  43.069
+919  262.69 12.960  73.291   291.81  88.297
+920  242.22 12.555  81.500   211.38  26.191
+921  296.91  8.642  36.016  1160.50  35.212
+922  201.87  6.211  32.162   113.35  13.319
+923 1302.47 14.635  70.047  1288.40  92.270
+924  138.65  8.798  45.243    77.66  16.282
+925  209.81 11.042 107.465   466.92  56.784
+926  670.33  4.811  49.567   352.04  24.907
+927  455.03 10.216  55.159   224.12  23.700
+929  350.01  9.985  41.474  1454.53  14.653
+930  331.18 16.598  56.147   663.45  41.931
+931  281.41  4.701  22.349   331.78  14.083
+932  123.06  2.407  22.636    87.75  13.364
+933  602.45 13.920  77.954   732.61  67.533
+934  246.31 11.851  62.797   316.42  15.530
+935 1635.28 23.839 110.801  2021.84 232.025
+936  219.01  4.242  17.434   191.17  99.920
+937  217.96  6.447  37.294   260.16  23.616
+938  416.19  7.156  36.991   552.76  18.908
+940   93.43  6.484  35.467   264.63  48.105
+941  419.04  7.545  16.894    96.79 320.558
+942 1156.49 17.252  69.045  1087.98 152.650
\ No newline at end of file
diff --git a/test/data/spouse.txt b/test/data/spouse.txt
new file mode 100644
index 0000000..43f295c
--- /dev/null
+++ b/test/data/spouse.txt
@@ -0,0 +1,30 @@
+2 3 5 5 4 4 5 5
+5 5 4 4 4 5 5 5
+4 5 5 5 4 4 5 5
+4 3 4 4 4 5 5 5
+3 3 5 5 4 4 5 5
+3 3 4 5 3 3 4 4
+3 4 4 4 4 3 5 4
+4 4 5 5 3 4 5 5
+4 5 5 5 4 4 5 4
+4 4 3 3 3 4 4 4
+4 4 5 5 4 5 5 5
+5 5 4 4 5 5 5 5
+4 4 4 4 4 4 5 5
+4 3 5 5 4 4 4 4
+4 4 5 5 4 4 5 5
+3 3 4 5 3 4 4 4
+4 5 4 4 5 5 5 5
+5 5 5 5 4 5 4 4
+5 5 4 4 3 4 4 4
+4 4 4 4 5 3 4 4
+4 4 4 4 5 3 4 4
+4 4 4 4 4 5 4 4
+3 4 5 5 2 5 5 5
+5 3 5 5 3 4 5 5
+5 5 3 3 4 3 5 5
+3 3 4 4 4 4 4 4
+4 4 4 4 4 4 5 5
+3 3 5 5 3 4 4 4
+4 4 3 3 4 4 5 4
+4 4 5 5 4 4 5 5
\ No newline at end of file
diff --git a/test/data/swiss3.txt b/test/data/swiss3.txt
new file mode 100644
index 0000000..80b07aa
--- /dev/null
+++ b/test/data/swiss3.txt
@@ -0,0 +1,200 @@
+real 214.8 131.0 131.1  9.0  9.7 141.0
+real 214.6 129.7 129.7  8.1  9.5 141.7
+real 214.8 129.7 129.7  8.7  9.6 142.2
+real 214.8 129.7 129.6  7.5 10.4 142.0
+real 215.0 129.6 129.7 10.4  7.7 141.8
+real 215.7 130.8 130.5  9.0 10.1 141.4
+real 215.5 129.5 129.7  7.9  9.6 141.6
+real 214.5 129.6 129.2  7.2 10.7 141.7
+real 214.9 129.4 129.7  8.2 11.0 141.9
+real 215.2 130.4 130.3  9.2 10.0 140.7
+real 215.3 130.4 130.3  7.9 11.7 141.8
+real 215.1 129.5 129.6  7.7 10.5 142.2
+real 215.2 130.8 129.6  7.9 10.8 141.4
+real 214.7 129.7 129.7  7.7 10.9 141.7
+real 215.1 129.9 129.7  7.7 10.8 141.8
+real 214.5 129.8 129.8  9.3  8.5 141.6
+real 214.6 129.9 130.1  8.2  9.8 141.7
+real 215.0 129.9 129.7  9.0  9.0 141.9
+real 215.2 129.6 129.6  7.4 11.5 141.5
+real 214.7 130.2 129.9  8.6 10.0 141.9
+real 215.0 129.9 129.3  8.4 10.0 141.4
+real 215.6 130.5 130.0  8.1 10.3 141.6
+real 215.3 130.6 130.0  8.4 10.8 141.5
+real 215.7 130.2 130.0  8.7 10.0 141.6
+real 215.1 129.7 129.9  7.4 10.8 141.1
+real 215.3 130.4 130.4  8.0 11.0 142.3
+real 215.5 130.2 130.1  8.9  9.8 142.4
+real 215.1 130.3 130.3  9.8  9.5 141.9
+real 215.1 130.0 130.0  7.4 10.5 141.8
+real 214.8 129.7 129.3  8.3  9.0 142.0
+real 215.2 130.1 129.8  7.9 10.7 141.8
+real 214.8 129.7 129.7  8.6  9.1 142.3
+real 215.0 130.0 129.6  7.7 10.5 140.7
+real 215.6 130.4 130.1  8.4 10.3 141.0
+real 215.9 130.4 130.0  8.9 10.6 141.4
+real 214.6 130.2 130.2  9.4  9.7 141.8
+real 215.5 130.3 130.0  8.4  9.7 141.8
+real 215.3 129.9 129.4  7.9 10.0 142.0
+real 215.3 130.3 130.1  8.5  9.3 142.1
+real 213.9 130.3 129.0  8.1  9.7 141.3
+real 214.4 129.8 129.2  8.9  9.4 142.3
+real 214.8 130.1 129.6  8.8  9.9 140.9
+real 214.9 129.6 129.4  9.3  9.0 141.7
+real 214.9 130.4 129.7  9.0  9.8 140.9
+real 214.8 129.4 129.1  8.2 10.2 141.0
+real 214.3 129.5 129.4  8.3 10.2 141.8
+real 214.8 129.9 129.7  8.3 10.2 141.5
+real 214.8 129.9 129.7  7.3 10.9 142.0
+real 214.6 129.7 129.8  7.9 10.3 141.1
+real 214.5 129.0 129.6  7.8  9.8 142.0
+real 214.6 129.8 129.4  7.2 10.0 141.3
+real 215.3 130.6 130.0  9.5  9.7 141.1
+real 214.5 130.1 130.0  7.8 10.9 140.9
+real 215.4 130.2 130.2  7.6 10.9 141.6
+real 214.5 129.4 129.5  7.9 10.0 141.4
+real 215.2 129.7 129.4  9.2  9.4 142.0
+real 215.7 130.0 129.4  9.2 10.4 141.2
+real 215.0 129.6 129.4  8.8  9.0 141.1
+real 215.1 130.1 129.9  7.9 11.0 141.3
+real 215.1 130.0 129.8  8.2 10.3 141.4
+real 215.1 129.6 129.3  8.3  9.9 141.6
+real 215.3 129.7 129.4  7.5 10.5 141.5
+real 215.4 129.8 129.4  8.0 10.6 141.5
+real 214.5 130.0 129.5  8.0 10.8 141.4
+real 215.0 130.0 129.8  8.6 10.6 141.5
+real 215.2 130.6 130.0  8.8 10.6 140.8
+real 214.6 129.5 129.2  7.7 10.3 141.3
+real 214.8 129.7 129.3  9.1  9.5 141.5
+real 215.1 129.6 129.8  8.6  9.8 141.8
+real 214.9 130.2 130.2  8.0 11.2 139.6
+real 213.8 129.8 129.5  8.4 11.1 140.9
+real 215.2 129.9 129.5  8.2 10.3 141.4
+real 215.0 129.6 130.2  8.7 10.0 141.2
+real 214.4 129.9 129.6  7.5 10.5 141.8
+real 215.2 129.9 129.7  7.2 10.6 142.1
+real 214.1 129.6 129.3  7.6 10.7 141.7
+real 214.9 129.9 130.1  8.8 10.0 141.2
+real 214.6 129.8 129.4  7.4 10.6 141.0
+real 215.2 130.5 129.8  7.9 10.9 140.9
+real 214.6 129.9 129.4  7.9 10.0 141.8
+real 215.1 129.7 129.7  8.6 10.3 140.6
+real 214.9 129.8 129.6  7.5 10.3 141.0
+real 215.2 129.7 129.1  9.0  9.7 141.9
+real 215.2 130.1 129.9  7.9 10.8 141.3
+real 215.4 130.7 130.2  9.0 11.1 141.2
+real 215.1 129.9 129.6  8.9 10.2 141.5
+real 215.2 129.9 129.7  8.7  9.5 141.6
+real 215.0 129.6 129.2  8.4 10.2 142.1
+real 214.9 130.3 129.9  7.4 11.2 141.5
+real 215.0 129.9 129.7  8.0 10.5 142.0
+real 214.7 129.7 129.3  8.6  9.6 141.6
+real 215.4 130.0 129.9  8.5  9.7 141.4
+real 214.9 129.4 129.5  8.2  9.9 141.5
+real 214.5 129.5 129.3  7.4 10.7 141.5
+real 214.7 129.6 129.5  8.3 10.0 142.0
+real 215.6 129.9 129.9  9.0  9.5 141.7
+real 215.0 130.4 130.3  9.1 10.2 141.1
+real 214.4 129.7 129.5  8.0 10.3 141.2
+real 215.1 130.0 129.8  9.1 10.2 141.5
+real 214.7 130.0 129.4  7.8 10.0 141.2
+fake 214.4 130.1 130.3  9.7 11.7 139.8
+fake 214.9 130.5 130.2 11.0 11.5 139.5
+fake 214.9 130.3 130.1  8.7 11.7 140.2
+fake 215.0 130.4 130.6  9.9 10.9 140.3
+fake 214.7 130.2 130.3 11.8 10.9 139.7
+fake 215.0 130.2 130.2 10.6 10.7 139.9
+fake 215.3 130.3 130.1  9.3 12.1 140.2
+fake 214.8 130.1 130.4  9.8 11.5 139.9
+fake 215.0 130.2 129.9 10.0 11.9 139.4
+fake 215.2 130.6 130.8 10.4 11.2 140.3
+fake 215.2 130.4 130.3  8.0 11.5 139.2
+fake 215.1 130.5 130.3 10.6 11.5 140.1
+fake 215.4 130.7 131.1  9.7 11.8 140.6
+fake 214.9 130.4 129.9 11.4 11.0 139.9
+fake 215.1 130.3 130.0 10.6 10.8 139.7
+fake 215.5 130.4 130.0  8.2 11.2 139.2
+fake 214.7 130.6 130.1 11.8 10.5 139.8
+fake 214.7 130.4 130.1 12.1 10.4 139.9
+fake 214.8 130.5 130.2 11.0 11.0 140.0
+fake 214.4 130.2 129.9 10.1 12.0 139.2
+fake 214.8 130.3 130.4 10.1 12.1 139.6
+fake 215.1 130.6 130.3 12.3 10.2 139.6
+fake 215.3 130.8 131.1 11.6 10.6 140.2
+fake 215.1 130.7 130.4 10.5 11.2 139.7
+fake 214.7 130.5 130.5  9.9 10.3 140.1
+fake 214.9 130.0 130.3 10.2 11.4 139.6
+fake 215.0 130.4 130.4  9.4 11.6 140.2
+fake 215.5 130.7 130.3 10.2 11.8 140.0
+fake 215.1 130.2 130.2 10.1 11.3 140.3
+fake 214.5 130.2 130.6  9.8 12.1 139.9
+fake 214.3 130.2 130.0 10.7 10.5 139.8
+fake 214.5 130.2 129.8 12.3 11.2 139.2
+fake 214.9 130.5 130.2 10.6 11.5 139.9
+fake 214.6 130.2 130.4 10.5 11.8 139.7
+fake 214.2 130.0 130.2 11.0 11.2 139.5
+fake 214.8 130.1 130.1 11.9 11.1 139.5
+fake 214.6 129.8 130.2 10.7 11.1 139.4
+fake 214.9 130.7 130.3  9.3 11.2 138.3
+fake 214.6 130.4 130.4 11.3 10.8 139.8
+fake 214.5 130.5 130.2 11.8 10.2 139.6
+fake 214.8 130.2 130.3 10.0 11.9 139.3
+fake 214.7 130.0 129.4 10.2 11.0 139.2
+fake 214.6 130.2 130.4 11.2 10.7 139.9
+fake 215.0 130.5 130.4 10.6 11.1 139.9
+fake 214.5 129.8 129.8 11.4 10.0 139.3
+fake 214.9 130.6 130.4 11.9 10.5 139.8
+fake 215.0 130.5 130.4 11.4 10.7 139.9
+fake 215.3 130.6 130.3  9.3 11.3 138.1
+fake 214.7 130.2 130.1 10.7 11.0 139.4
+fake 214.9 129.9 130.0  9.9 12.3 139.4
+fake 214.9 130.3 129.9 11.9 10.6 139.8
+fake 214.6 129.9 129.7 11.9 10.1 139.0
+fake 214.6 129.7 129.3 10.4 11.0 139.3
+fake 214.5 130.1 130.1 12.1 10.3 139.4
+fake 214.5 130.3 130.0 11.0 11.5 139.5
+fake 215.1 130.0 130.3 11.6 10.5 139.7
+fake 214.2 129.7 129.6 10.3 11.4 139.5
+fake 214.4 130.1 130.0 11.3 10.7 139.2
+fake 214.8 130.4 130.6 12.5 10.0 139.3
+fake 214.6 130.6 130.1  8.1 12.1 137.9
+fake 215.6 130.1 129.7  7.4 12.2 138.4
+fake 214.9 130.5 130.1  9.9 10.2 138.1
+fake 214.6 130.1 130.0 11.5 10.6 139.5
+fake 214.7 130.1 130.2 11.6 10.9 139.1
+fake 214.3 130.3 130.0 11.4 10.5 139.8
+fake 215.1 130.3 130.6 10.3 12.0 139.7
+fake 216.3 130.7 130.4 10.0 10.1 138.8
+fake 215.6 130.4 130.1  9.6 11.2 138.6
+fake 214.8 129.9 129.8  9.6 12.0 139.6
+fake 214.9 130.0 129.9 11.4 10.9 139.7
+fake 213.9 130.7 130.5  8.7 11.5 137.8
+fake 214.2 130.6 130.4 12.0 10.2 139.6
+fake 214.8 130.5 130.3 11.8 10.5 139.4
+fake 214.8 129.6 130.0 10.4 11.6 139.2
+fake 214.8 130.1 130.0 11.4 10.5 139.6
+fake 214.9 130.4 130.2 11.9 10.7 139.0
+fake 214.3 130.1 130.1 11.6 10.5 139.7
+fake 214.5 130.4 130.0  9.9 12.0 139.6
+fake 214.8 130.5 130.3 10.2 12.1 139.1
+fake 214.5 130.2 130.4  8.2 11.8 137.8
+fake 215.0 130.4 130.1 11.4 10.7 139.1
+fake 214.8 130.6 130.6  8.0 11.4 138.7
+fake 215.0 130.5 130.1 11.0 11.4 139.3
+fake 214.6 130.5 130.4 10.1 11.4 139.3
+fake 214.7 130.2 130.1 10.7 11.1 139.5
+fake 214.7 130.4 130.0 11.5 10.7 139.4
+fake 214.5 130.4 130.0  8.0 12.2 138.5
+fake 214.8 130.0 129.7 11.4 10.6 139.2
+fake 214.8 129.9 130.2  9.6 11.9 139.4
+fake 214.6 130.3 130.2 12.7  9.1 139.2
+fake 215.1 130.2 129.8 10.2 12.0 139.4
+fake 215.4 130.5 130.6  8.8 11.0 138.6
+fake 214.7 130.3 130.2 10.8 11.1 139.2
+fake 215.0 130.5 130.3  9.6 11.0 138.5
+fake 214.9 130.3 130.5 11.6 10.6 139.8
+fake 215.0 130.4 130.3  9.9 12.1 139.6
+fake 215.1 130.3 129.9 10.3 11.5 139.7
+fake 214.8 130.3 130.4 10.6 11.1 140.0
+fake 214.7 130.7 130.8 11.2 11.2 139.4
+fake 214.3 129.9 129.9 10.2 11.5 139.6
\ No newline at end of file
diff --git a/test/hotelling.jl b/test/hotelling.jl
new file mode 100644
index 0000000..3878d58
--- /dev/null
+++ b/test/hotelling.jl
@@ -0,0 +1,89 @@
+using HypothesisTests
+using Test
+using DelimitedFiles
+
+@testset "Hotelling utility functions" begin
+    HT = HypothesisTests
+
+    @test HT.checkdims(rand(3, 2), rand(4, 2)) == (2, 3, 4)
+    # Mismatched number of variables
+    @test_throws DimensionMismatch HT.checkdims(rand(4, 3), rand(4, 6))
+    # Empty observations
+    @test_throws ArgumentError HT.checkdims(rand(0, 4), rand(2, 4))
+
+    Sx = [ 3.0 -1.5 0.0
+          -1.5  1.0 0.5
+           0.0  0.5 1.0]
+    Sy = fill(4.5, (3, 3))
+    out = [3.5 0.5      1.5
+           0.5 2.166667 1.833333
+           1.5 1.833333 2.166667]
+    @test HT.poolcov!(Sx, 2, Sy, 1) ≈ out atol=1e-6
+    # Input is modified
+    @test Sx ≈ out atol=1e-6
+
+    # Positive semi-definite but not positive definite
+    P = [ 1.0000  0.7426  0.1601 -0.7000 0.5500
+          0.7426  1.0000 -0.2133 -0.5818 0.5000
+          0.1601 -0.2133  1.0000 -0.1121 0.1000
+         -0.7000 -0.5818 -0.1121  1.0000 0.4500
+          0.5500  0.5000  0.1000  0.4500 1.0000]
+    # Positive definite
+    Q = [ 1.0 -0.5  0.0
+         -0.5  1.0 -0.5
+          0.0 -0.5  1.0]
+    @test @inferred(HT.At_Binv_A(ones(5), P)) ≈ -0.8008792 atol=1e-6
+    @test @inferred(HT.At_Binv_A(ones(3), Q)) ≈ 10.0
+end
+
+@testset "One sample Hotelling's T²" begin
+    # Columns are calcium, iron, protein, vitamin A, vitamin C
+    nutrient = readdlm(joinpath(@__DIR__, "data", "nutrient.txt"))[:,2:end]
+
+    t = OneSampleHotellingT2(nutrient, [1000, 15, 60, 800, 75])
+    @test nobs(t) == 737
+    @test dof(t) == (5, 732)
+    @test pvalue(t) ≈ 0.0 atol=eps()
+    @test t.T² ≈ 1758.5413137 atol=1e-6
+    @test t.F ≈ 349.7968048 atol=1e-6
+    let out = sprint(show, t)
+        @test occursin("reject h_0", out) && !occursin("fail to", out)
+    end
+
+    # Columns are survey answers: husbands' to questions 1-4, wives' to questions 1-4
+    spouse = readdlm(joinpath(@__DIR__, "data", "spouse.txt"))
+
+    # Paired
+    p = OneSampleHotellingT2(spouse[:,1:4], spouse[:,5:end])
+    @test nobs(p) == 30
+    @test dof(p) == (4, 26)
+    @test pvalue(p) ≈ 0.039369144 atol=1e-6
+    @test p.T² ≈ 13.127840261 atol=1e-6
+    @test p.F ≈ 2.942446955 atol=1e-6
+    let out = sprint(show, p)
+        @test occursin("reject h_0", out) && !occursin("fail to", out)
+    end
+end
+
+@testset "Two sample Hotelling's T²" begin
+    # Columns are type, length, left, right, bottom, top, diag
+    swiss = readdlm(joinpath(@__DIR__, "data", "swiss3.txt"))
+    genuine = convert(Matrix{Float64}, swiss[view(swiss, :, 1) .== "real", 2:end])
+    counterfeit = convert(Matrix{Float64}, swiss[view(swiss, :, 1) .== "fake", 2:end])
+
+    eq = EqualCovHotellingT2(genuine, counterfeit)
+    @test nobs(eq) == (100, 100)
+    @test dof(eq) == (6, 193)
+    @test pvalue(eq) ≈ 0.0 atol=eps()
+    @test eq.T² ≈ 2412.4506855 atol=1e-6
+    @test eq.F ≈ 391.9217023 atol=1e-6
+    @test occursin("reject h_0", sprint(show, eq))
+
+    un = UnequalCovHotellingT2(genuine, counterfeit)
+    @test nobs(un) == (100, 100)
+    @test dof(un) == (6, 193)
+    @test pvalue(un) ≈ 0.0 atol=eps()
+    @test un.T² ≈ 2412.4506855 atol=1e-6
+    @test un.F ≈ 391.9217023 atol=1e-6
+    @test occursin("reject h_0", sprint(show, un))
+end
diff --git a/test/runtests.jl b/test/runtests.jl
index 4cc4277..2c6ecbd 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -31,3 +31,4 @@ include("show.jl")
 include("t.jl")
 include("wilcoxon.jl")
 include("z.jl")
+include("hotelling.jl")