Test for the difference of two cluster means after k-means clustering, for matrix normal model with arbitrary scale matrices.

Usage

test.clusters.km(
  X,
  U = NULL,
  Sigma = NULL,
  Y = NULL,
  UY = NULL,
  precUY = NULL,
  NC,
  clusters,
  itermax = 10,
  tol = 1e-06
)

Arguments

X: A \(n \times p\) matrix drawn from a \(n \times p\) matrix normal distribution \(\mathcal{MN}(\)M, U, Sigma\()\). X must have \(n\) rows and \(p\) columns.
U: A \(n \times n\) positive-definite matrix describing the dependence structure between the rows in X. If NULL, observations are considered independent and U is set to the \(n \times n\) identity matrix.
Sigma: A \(p \times p\) positive-definite matrix describing the dependence structure between the columns in X. If NULL, Sigma is over-estimated (in the sense of the Loewner partial order).
Y: If Sigma is NULL, an i.i.d. copy of X allowing its estimation. Y must have the same number of columns as X.
UY: If Sigma is NULL, a positive-definite matrix describing the dependence structure between the rows in Y. If NULL and its inverse is not provided, set to the identity matrix by default.
precUY: The inverse matrix of UY, that can be provided to increase computational efficiency. If UY is not NULL and precUY is NULL, precUY is obtained by inverting UY.
NC: The number of clusters to choose.
clusters: A vector of two integers from 1 to NC indicating the pair of clusters whose means have to be compared.
itermax: The iter.max parameter of the k-means algorithm in kmeans_estimation function of KmeansInference package.
tol: The tol_eps parameter of the k-means algorithm in kmeans_estimation function of KmeansInference package.

Value

pvalue - The p-value for the difference of cluster means.
stat - The test statistic.
km - The partition of the n observations retrieved by the clustering algorithm.

References

[1] L. L. Gao, J. Bien, and D. Witten. Selective inference for hierarchical clustering. Journal of the American Statistical Association, 0(0):1–11, 2022. [2] Y. T. Chen and D. M. Witten. Selective inference for k-means clustering, 2022. arXiv:2203.15267.

Examples


n <- 50
p <- 20
M <- Matrix::Matrix(0, nrow = n , ncol = p) # Mean matrix
Sigma <- stats::toeplitz(seq(1, 0.1, length = p)) # Sigma: dependence between features
U <- matrixNormal::I(n) # U: dependence between observations
X <- matrixNormal::rmatnorm(s = 1, M, U, Sigma)
Y <- matrixNormal::rmatnorm(s = 1, M, U, Sigma) # i.i.d. copy of X

# k-means under the global null hypothesis
test.clusters.km(X, U, Sigma, NC = 3, clusters = sample(1:3, 2))
#> $pvalue
#> [1] 0.08818978
#> 
#> $stat
#> [1] 6.895642
#> 
#> $km
#>  [1] 3 3 2 2 1 1 2 2 3 3 3 2 3 1 1 2 1 2 2 1 1 3 3 3 1 3 3 3 1 3 2 1 1 3 2 1 2 2
#> [39] 1 1 3 2 3 1 3 3 2 1 1 3
#> 
#> $Sigma
#> 20 x 20 Matrix of class "dsyMatrix"
#>            [,1]      [,2]      [,3]      [,4]      [,5]      [,6]      [,7]
#>  [1,] 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579 0.7157895
#>  [2,] 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579
#>  [3,] 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263
#>  [4,] 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947
#>  [5,] 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632
#>  [6,] 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316
#>  [7,] 0.7157895 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000
#>  [8,] 0.6684211 0.7157895 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316
#>  [9,] 0.6210526 0.6684211 0.7157895 0.7631579 0.8105263 0.8578947 0.9052632
#> [10,] 0.5736842 0.6210526 0.6684211 0.7157895 0.7631579 0.8105263 0.8578947
#> [11,] 0.5263158 0.5736842 0.6210526 0.6684211 0.7157895 0.7631579 0.8105263
#> [12,] 0.4789474 0.5263158 0.5736842 0.6210526 0.6684211 0.7157895 0.7631579
#> [13,] 0.4315789 0.4789474 0.5263158 0.5736842 0.6210526 0.6684211 0.7157895
#> [14,] 0.3842105 0.4315789 0.4789474 0.5263158 0.5736842 0.6210526 0.6684211
#> [15,] 0.3368421 0.3842105 0.4315789 0.4789474 0.5263158 0.5736842 0.6210526
#> [16,] 0.2894737 0.3368421 0.3842105 0.4315789 0.4789474 0.5263158 0.5736842
#> [17,] 0.2421053 0.2894737 0.3368421 0.3842105 0.4315789 0.4789474 0.5263158
#> [18,] 0.1947368 0.2421053 0.2894737 0.3368421 0.3842105 0.4315789 0.4789474
#> [19,] 0.1473684 0.1947368 0.2421053 0.2894737 0.3368421 0.3842105 0.4315789
#> [20,] 0.1000000 0.1473684 0.1947368 0.2421053 0.2894737 0.3368421 0.3842105
#>            [,8]      [,9]     [,10]     [,11]     [,12]     [,13]     [,14]
#>  [1,] 0.6684211 0.6210526 0.5736842 0.5263158 0.4789474 0.4315789 0.3842105
#>  [2,] 0.7157895 0.6684211 0.6210526 0.5736842 0.5263158 0.4789474 0.4315789
#>  [3,] 0.7631579 0.7157895 0.6684211 0.6210526 0.5736842 0.5263158 0.4789474
#>  [4,] 0.8105263 0.7631579 0.7157895 0.6684211 0.6210526 0.5736842 0.5263158
#>  [5,] 0.8578947 0.8105263 0.7631579 0.7157895 0.6684211 0.6210526 0.5736842
#>  [6,] 0.9052632 0.8578947 0.8105263 0.7631579 0.7157895 0.6684211 0.6210526
#>  [7,] 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579 0.7157895 0.6684211
#>  [8,] 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579 0.7157895
#>  [9,] 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579
#> [10,] 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263
#> [11,] 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947
#> [12,] 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632
#> [13,] 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316
#> [14,] 0.7157895 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000
#> [15,] 0.6684211 0.7157895 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316
#> [16,] 0.6210526 0.6684211 0.7157895 0.7631579 0.8105263 0.8578947 0.9052632
#> [17,] 0.5736842 0.6210526 0.6684211 0.7157895 0.7631579 0.8105263 0.8578947
#> [18,] 0.5263158 0.5736842 0.6210526 0.6684211 0.7157895 0.7631579 0.8105263
#> [19,] 0.4789474 0.5263158 0.5736842 0.6210526 0.6684211 0.7157895 0.7631579
#> [20,] 0.4315789 0.4789474 0.5263158 0.5736842 0.6210526 0.6684211 0.7157895
#>           [,15]     [,16]     [,17]     [,18]     [,19]     [,20]
#>  [1,] 0.3368421 0.2894737 0.2421053 0.1947368 0.1473684 0.1000000
#>  [2,] 0.3842105 0.3368421 0.2894737 0.2421053 0.1947368 0.1473684
#>  [3,] 0.4315789 0.3842105 0.3368421 0.2894737 0.2421053 0.1947368
#>  [4,] 0.4789474 0.4315789 0.3842105 0.3368421 0.2894737 0.2421053
#>  [5,] 0.5263158 0.4789474 0.4315789 0.3842105 0.3368421 0.2894737
#>  [6,] 0.5736842 0.5263158 0.4789474 0.4315789 0.3842105 0.3368421
#>  [7,] 0.6210526 0.5736842 0.5263158 0.4789474 0.4315789 0.3842105
#>  [8,] 0.6684211 0.6210526 0.5736842 0.5263158 0.4789474 0.4315789
#>  [9,] 0.7157895 0.6684211 0.6210526 0.5736842 0.5263158 0.4789474
#> [10,] 0.7631579 0.7157895 0.6684211 0.6210526 0.5736842 0.5263158
#> [11,] 0.8105263 0.7631579 0.7157895 0.6684211 0.6210526 0.5736842
#> [12,] 0.8578947 0.8105263 0.7631579 0.7157895 0.6684211 0.6210526
#> [13,] 0.9052632 0.8578947 0.8105263 0.7631579 0.7157895 0.6684211
#> [14,] 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579 0.7157895
#> [15,] 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263 0.7631579
#> [16,] 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947 0.8105263
#> [17,] 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632 0.8578947
#> [18,] 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316 0.9052632
#> [19,] 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000 0.9526316
#> [20,] 0.7631579 0.8105263 0.8578947 0.9052632 0.9526316 1.0000000
#> 
# k-means under the global null hypothesis and over-estimation of Sigma
test.clusters.km(X, U, Sigma = NULL, Y = Y, NC = 3, clusters = sample(1:3, 2))
#> Sigma not provided: plugging an over-estimate.
#> $pvalue
#> [1] 0.1231035
#> 
#> $stat
#> [1] 7.116161
#> 
#> $km
#>  [1] 2 2 1 1 3 1 1 1 2 1 2 1 2 3 3 1 3 1 1 1 3 1 1 2 1 2 1 1 3 2 1 3 3 2 1 3 1 1
#> [39] 3 3 1 1 2 1 2 2 1 3 3 1
#> 
#> $Sigma
#> 20 x 20 Matrix of class "dgeMatrix"
#>            [,1]      [,2]      [,3]      [,4]      [,5]      [,6]      [,7]
#>  [1,] 0.9603016 0.9305723 0.9159794 0.9233182 0.8375435 0.7510176 0.7252088
#>  [2,] 0.9305723 0.9807935 0.9851779 0.9888548 0.9155179 0.8430032 0.8188754
#>  [3,] 0.9159794 0.9851779 1.0830863 1.0763937 1.0102268 0.9221182 0.8937114
#>  [4,] 0.9233182 0.9888548 1.0763937 1.1619349 1.0919877 0.9960759 0.9779506
#>  [5,] 0.8375435 0.9155179 1.0102268 1.0919877 1.1410679 1.0614338 1.0502664
#>  [6,] 0.7510176 0.8430032 0.9221182 0.9960759 1.0614338 1.1034975 1.0649232
#>  [7,] 0.7252088 0.8188754 0.8937114 0.9779506 1.0502664 1.0649232 1.1072610
#>  [8,] 0.7159366 0.8076853 0.8934622 0.9860960 1.0662553 1.0929081 1.1399145
#>  [9,] 0.6726643 0.7612316 0.8160893 0.9335137 0.9987469 1.0257534 1.0810622
#> [10,] 0.5964272 0.6747457 0.7230974 0.8326702 0.9019117 0.9196906 0.9907687
#> [11,] 0.5624121 0.6222424 0.6803674 0.7859979 0.8725333 0.8645810 0.9307502
#> [12,] 0.5437777 0.5895832 0.6814992 0.7960742 0.8679777 0.8424430 0.9163504
#> [13,] 0.5208889 0.5605890 0.6484986 0.7489808 0.8064976 0.7546865 0.8428468
#> [14,] 0.5339751 0.5695625 0.6587417 0.7534454 0.8152375 0.7427396 0.8149990
#> [15,] 0.4857298 0.5135660 0.6152519 0.6960051 0.7542550 0.6702716 0.7487294
#> [16,] 0.3923730 0.4266917 0.5446956 0.6190627 0.7031704 0.6289224 0.7118074
#> [17,] 0.3401214 0.3616328 0.4850441 0.5666769 0.6482051 0.5780694 0.6732186
#> [18,] 0.2792818 0.3134447 0.4350400 0.5008451 0.5763099 0.5053886 0.6162580
#> [19,] 0.2799745 0.3096147 0.4153897 0.4778372 0.5665693 0.5182635 0.6249115
#> [20,] 0.2467922 0.2929188 0.3906825 0.4412738 0.5666690 0.5362469 0.6388518
#>            [,8]      [,9]     [,10]     [,11]     [,12]     [,13]     [,14]
#>  [1,] 0.7159366 0.6726643 0.5964272 0.5624121 0.5437777 0.5208889 0.5339751
#>  [2,] 0.8076853 0.7612316 0.6747457 0.6222424 0.5895832 0.5605890 0.5695625
#>  [3,] 0.8934622 0.8160893 0.7230974 0.6803674 0.6814992 0.6484986 0.6587417
#>  [4,] 0.9860960 0.9335137 0.8326702 0.7859979 0.7960742 0.7489808 0.7534454
#>  [5,] 1.0662553 0.9987469 0.9019117 0.8725333 0.8679777 0.8064976 0.8152375
#>  [6,] 1.0929081 1.0257534 0.9196906 0.8645810 0.8424430 0.7546865 0.7427396
#>  [7,] 1.1399145 1.0810622 0.9907687 0.9307502 0.9163504 0.8428468 0.8149990
#>  [8,] 1.2631818 1.1898429 1.0910838 1.0345943 1.0439247 0.9630220 0.9467901
#>  [9,] 1.1898429 1.2021147 1.1033170 1.0377979 1.0371751 0.9612791 0.9412960
#> [10,] 1.0910838 1.1033170 1.0928598 1.0208504 1.0188849 0.9285940 0.9069454
#> [11,] 1.0345943 1.0377979 1.0208504 1.0342524 1.0333823 0.9389519 0.9317146
#> [12,] 1.0439247 1.0371751 1.0188849 1.0333823 1.1142905 1.0328720 1.0250947
#> [13,] 0.9630220 0.9612791 0.9285940 0.9389519 1.0328720 1.0565618 1.0347314
#> [14,] 0.9467901 0.9412960 0.9069454 0.9317146 1.0250947 1.0347314 1.1055911
#> [15,] 0.8865048 0.8717669 0.8521642 0.8712581 0.9694873 0.9924925 1.0634564
#> [16,] 0.8575011 0.8325855 0.8229976 0.8408441 0.9591151 0.9671134 1.0227154
#> [17,] 0.7901468 0.8048192 0.7849075 0.8079824 0.9166237 0.9192575 0.9579884
#> [18,] 0.7234925 0.7273194 0.6825019 0.7038466 0.8002181 0.8186977 0.8411513
#> [19,] 0.7252500 0.7322990 0.6961550 0.7190832 0.8078791 0.8052253 0.8252309
#> [20,] 0.7444814 0.7401713 0.7222783 0.7390451 0.8016566 0.7791165 0.8033116
#>           [,15]     [,16]     [,17]     [,18]     [,19]     [,20]
#>  [1,] 0.4857298 0.3923730 0.3401214 0.2792818 0.2799745 0.2467922
#>  [2,] 0.5135660 0.4266917 0.3616328 0.3134447 0.3096147 0.2929188
#>  [3,] 0.6152519 0.5446956 0.4850441 0.4350400 0.4153897 0.3906825
#>  [4,] 0.6960051 0.6190627 0.5666769 0.5008451 0.4778372 0.4412738
#>  [5,] 0.7542550 0.7031704 0.6482051 0.5763099 0.5665693 0.5666690
#>  [6,] 0.6702716 0.6289224 0.5780694 0.5053886 0.5182635 0.5362469
#>  [7,] 0.7487294 0.7118074 0.6732186 0.6162580 0.6249115 0.6388518
#>  [8,] 0.8865048 0.8575011 0.7901468 0.7234925 0.7252500 0.7444814
#>  [9,] 0.8717669 0.8325855 0.8048192 0.7273194 0.7322990 0.7401713
#> [10,] 0.8521642 0.8229976 0.7849075 0.6825019 0.6961550 0.7222783
#> [11,] 0.8712581 0.8408441 0.8079824 0.7038466 0.7190832 0.7390451
#> [12,] 0.9694873 0.9591151 0.9166237 0.8002181 0.8078791 0.8016566
#> [13,] 0.9924925 0.9671134 0.9192575 0.8186977 0.8052253 0.7791165
#> [14,] 1.0634564 1.0227154 0.9579884 0.8411513 0.8252309 0.8033116
#> [15,] 1.0996231 1.0777985 1.0025336 0.8882201 0.8619984 0.8547137
#> [16,] 1.0777985 1.1452911 1.0851325 0.9681028 0.9307063 0.9160037
#> [17,] 1.0025336 1.0851325 1.1405418 1.0181164 1.0069858 0.9566722
#> [18,] 0.8882201 0.9681028 1.0181164 0.9712650 0.9373330 0.8841727
#> [19,] 0.8619984 0.9307063 1.0069858 0.9373330 0.9870705 0.9440803
#> [20,] 0.8547137 0.9160037 0.9566722 0.8841727 0.9440803 1.0095941
#>