R Under development (unstable) (2023-12-09 r85665 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2023 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library("matrixStats")
> 
> set.seed(1)
> 
> x <- matrix(rnorm(20), nrow = 5L, ncol = 4L)
> print(x)
           [,1]       [,2]       [,3]        [,4]
[1,] -0.6264538 -0.8204684  1.5117812 -0.04493361
[2,]  0.1836433  0.4874291  0.3898432 -0.01619026
[3,] -0.8356286  0.7383247 -0.6212406  0.94383621
[4,]  1.5952808  0.5757814 -2.2146999  0.82122120
[5,]  0.3295078 -0.3053884  1.1249309  0.59390132
> 
> # To check names attribute
> dimnames <- list(letters[1:5], LETTERS[1:4])
> 
> # Weighted row variances (uniform weights - all w = 1)
> # Non-weighted row variances
> x_est0 <- rowVars(x)
> w <- rep(1, times = ncol(x))
> x_est1 <- rowWeightedVars(x, w = w)
> print(x_est1)
[1] 1.11767161 0.05022969 0.83582537 2.76819528 0.35351857
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedVars(t(x), w = w)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedVars(x, w = w, useNames = FALSE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowVars(x, useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedVars(x, w = w, useNames = NA)
+   x_est2 <- colWeightedVars(t(x), w = w, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }
> x_est1 <- rowWeightedVars(x, w = w, useNames = TRUE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> 
> # Weighted row variances (uniform weights - all w = 3)
> x3 <- cbind(x, x, x)
> x_est0 <- rowVars(x3)
> w <- rep(3, times = ncol(x))
> x_est1 <- rowWeightedVars(x, w = w)
> print(x_est1)
[1] 0.91445859 0.04109702 0.68385712 2.26488705 0.28924246
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedVars(t(x), w = w)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedVars(x, w = w, useNames = FALSE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x3 <- cbind(x, x, x)
> x_est0 <- rowVars(x3, useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedVars(x, w = w, useNames = NA)
+   x_est2 <- colWeightedVars(t(x), w = w, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }
> x_est1 <- rowWeightedVars(x, w = w, useNames = TRUE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> 
> # Weighted row variances (excluding some columns)
> w <- c(1, 1, 0, 1)
> x_est0 <- rowVars(x[, (w == 1), drop = FALSE])
> x_est1 <- rowWeightedVars(x, w = w)
> print(x_est1)
[1] 0.16287693 0.06430861 0.94767651 0.28313135 0.21361984
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedVars(t(x), w = w)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedVars(x, w = w, useNames = FALSE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowVars(x[, (w == 1), drop = FALSE], useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedVars(x, w = w, useNames = NA)
+   x_est2 <- colWeightedVars(t(x), w = w, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }
> x_est1 <- rowWeightedVars(x, w = w, useNames = TRUE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> 
> # Weighted row variances (excluding some columns)
> w <- c(0, 1, 0, 0)
> x_est0 <- rowVars(x[, (w == 1), drop = FALSE])
> x_est1 <- rowWeightedVars(x, w = w)
> #stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedVars(t(x), w = w)
> stopifnot(all.equal(x_est2, x_est1))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedVars(x, w = w, useNames = FALSE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowVars(x[, (w == 1), drop = FALSE], useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedVars(x, w = w, useNames = NA)
+   x_est2 <- colWeightedVars(t(x), w = w, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }
> x_est1 <- rowWeightedVars(x, w = w, useNames = TRUE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> 
> # Weighted row variances (all zero weights)
> w <- c(0, 0, 0, 0)
> x_est0 <- rowVars(x[, (w == 1), drop = FALSE])
> x_est1 <- rowWeightedVars(x, w = w)
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedVars(t(x), w = w)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedVars(x, w = w, useNames = FALSE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowVars(x[, (w == 1), drop = FALSE], useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedVars(x, w = w, useNames = NA)
+   x_est2 <- colWeightedVars(t(x), w = w, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }
> x_est1 <- rowWeightedVars(x, w = w, useNames = TRUE)
> x_est2 <- colWeightedVars(t(x), w = w, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> # Weighted variances by rows and columns
> w <- 1:4
> # Test with and without dimnames on x
> for (setDimnames in c(TRUE, FALSE)) {
+   if (setDimnames) dimnames(x) <- dimnames
+   else dimnames(x) <- NULL    
+   # Check names attribute
+   for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
+     x_est1 <- rowWeightedVars(x, w = w, useNames = useNames)
+     print(x_est1)
+     x_est2 <- colWeightedVars(t(x), w = w, useNames = useNames)
+     stopifnot(all.equal(x_est2, x_est1))
+   }
+ }
         a          b          c          d          e 
0.90747094 0.05040463 0.65089440 2.28697479 0.27902734 
[1] 0.90747094 0.05040463 0.65089440 2.28697479 0.27902734
[1] 0.90747094 0.05040463 0.65089440 2.28697479 0.27902734
[1] 0.90747094 0.05040463 0.65089440 2.28697479 0.27902734
> 
> 
> x[sample(length(x), size = 0.3 * length(x))] <- NA
> print(x)
           [,1]      [,2]       [,3]        [,4]
[1,] -0.6264538        NA  1.5117812 -0.04493361
[2,]  0.1836433        NA         NA -0.01619026
[3,] -0.8356286 0.7383247 -0.6212406  0.94383621
[4,]  1.5952808        NA -2.2146999  0.82122120
[5,]  0.3295078        NA         NA  0.59390132
> 
> # Non-weighted row variances with missing values
> x_est0 <- rowVars(x, na.rm = TRUE)
> x_est1 <- rowWeightedVars(x, w = rep(1, times = ncol(x)), na.rm = TRUE)
> print(x_est1)
[1] 1.22226258 0.01996673 0.83582537 4.05532299 0.03495197
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedVars(t(x), w = rep(1, times = ncol(x)), na.rm = TRUE)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedVars(x, w = rep(1, times = ncol(x)), na.rm = TRUE, useNames = FALSE)
> x_est2 <- colWeightedVars(t(x), w = rep(1, times = ncol(x)), na.rm = TRUE, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowVars(x, na.rm = TRUE, useNames = TRUE)
> x_est1 <- rowWeightedVars(x, w = rep(1, times = ncol(x)), na.rm = TRUE, useNames = TRUE)
> x_est2 <- colWeightedVars(t(x), w = rep(1, times = ncol(x)), na.rm = TRUE, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> 
> # Weighted row variances with missing values
> # Test with and without dimnames on x
> for (setDimnames in c(TRUE, FALSE)) {
+   if (setDimnames) dimnames(x) <- dimnames
+   else dimnames(x) <- NULL    
+   # Check names attribute
+   for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
+     x_est1 <- rowWeightedVars(x, w = w, na.rm = TRUE, useNames = useNames)
+     print(x_est1)
+     x_est2 <- colWeightedVars(t(x), w = w, na.rm = TRUE, useNames = useNames)
+     stopifnot(all.equal(x_est2, x_est1))
+   }
+ }
          a           b           c           d           e 
0.788377504 0.007986693 0.650894398 2.795470241 0.013980790 
[1] 0.788377504 0.007986693 0.650894398 2.795470241 0.013980790
[1] 0.788377504 0.007986693 0.650894398 2.795470241 0.013980790
[1] 0.788377504 0.007986693 0.650894398 2.795470241 0.013980790
> 
> 
> # Weighted variances by rows and columns
> w <- 1:4
> # Test with and without dimnames on x
> for (setDimnames in c(TRUE, FALSE)) {
+   if (setDimnames) dimnames(x) <- dimnames
+   else dimnames(x) <- NULL    
+   # Check names attribute
+   for (useNames in c(if (!matrixStats:::isUseNamesNADefunct()) NA, TRUE, FALSE)) {
+     x_est1 <- rowWeightedVars(x, w = w, na.rm = TRUE, useNames = useNames)
+     print(x_est1)
+     x_est2 <- colWeightedVars(t(x), w = w, na.rm = TRUE, useNames = useNames)
+     stopifnot(all.equal(x_est2, x_est1))
+   }
+ }
          a           b           c           d           e 
0.788377504 0.007986693 0.650894398 2.795470241 0.013980790 
[1] 0.788377504 0.007986693 0.650894398 2.795470241 0.013980790
[1] 0.788377504 0.007986693 0.650894398 2.795470241 0.013980790
[1] 0.788377504 0.007986693 0.650894398 2.795470241 0.013980790
> 
> 
> # Weighted row standard deviation (excluding some columns)
> w <- c(1, 1, 0, 1)
> ## FIXME: rowVars()/rowSds() needs na.rm = FALSE (wrong default)
> x_est0 <- rowSds(x[, (w == 1), drop = FALSE], na.rm = FALSE)
> x_est1 <- rowWeightedSds(x, w = w)
> print(x_est1)
[1]        NA        NA 0.9734868        NA        NA
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedSds(t(x), w = w, na.rm = FALSE)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedSds(x, w = w, na.rm = FALSE, useNames = FALSE)
> x_est2 <- colWeightedSds(t(x), w = w, na.rm = FALSE, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowSds(x[, (w == 1), drop = FALSE], na.rm = FALSE, useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedSds(x, w = w, na.rm = FALSE, useNames = NA)
+   x_est2 <- colWeightedSds(t(x), w = w, na.rm = FALSE, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }  
> x_est1 <- rowWeightedSds(x, w = w, na.rm = FALSE, useNames = TRUE)
> x_est2 <- colWeightedSds(t(x), w = w, na.rm = FALSE, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> 
> # Weighted row MADs (excluding some columns)
> w <- c(1, 1, 0, 1)
> x_est0 <- rowMads(x[, (w == 1), drop = FALSE])
> x_est1 <- rowWeightedMads(x, w = w)
> print(x_est1)
[1]        NA        NA 0.3046914        NA        NA
> stopifnot(all.equal(x_est1, x_est0))
> x_est2 <- colWeightedMads(t(x), w = w)
> stopifnot(all.equal(x_est2, x_est0))
> # Check names attribute
> dimnames(x) <- dimnames
> x_est1 <- rowWeightedMads(x, w = w, useNames = FALSE)
> x_est2 <- colWeightedMads(t(x), w = w, useNames = FALSE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> x_est0 <- rowMads(x[, (w == 1), drop = FALSE], useNames = TRUE)
> if (!matrixStats:::isUseNamesNADefunct()) {
+   x_est1 <- rowWeightedMads(x, w = w, useNames = NA)
+   x_est2 <- colWeightedMads(t(x), w = w, useNames = NA)
+   stopifnot(all.equal(x_est1, x_est0))
+   stopifnot(all.equal(x_est2, x_est0))
+ }
> x_est1 <- rowWeightedMads(x, w = w, useNames = TRUE)
> x_est2 <- colWeightedMads(t(x), w = w, useNames = TRUE)
> stopifnot(all.equal(x_est1, x_est0))
> stopifnot(all.equal(x_est2, x_est0))
> dimnames(x) <- NULL
> 
> proc.time()
   user  system elapsed 
   0.18    0.04    0.21