diff --git a/tests/testthat/test-meanAdjAgree_index.mean.R b/tests/testthat/test-meanAdjAgree_index.mean.R
new file mode 100644
index 0000000000000000000000000000000000000000..21e9f93d7bd4bf8259e686e00b3b6f7edf53d555
--- /dev/null
+++ b/tests/testthat/test-meanAdjAgree_index.mean.R
@@ -0,0 +1,27 @@
+test_that("mean.index always returns a numeric vector of length(index.variables)", {
+  list.res <- rlist::list.flatten(
+    list(
+      list(`1` = c(0.0, 0.5, 0.5, 0.5)),
+      list(`1` = c(0.0, 0.6, 0.7, 0.8)),
+      list(`5` = c(0.4, 0.0, 0.4, 0.4)),
+      list(`100` = c(0.9, 0.3, 0.0, 0.9)),
+      list(`10000` = c(0.7, 0.8, 0.9, 0.0)),
+      list(`5` = c(0.3, 0.0, 0.5, 0.6), `100` = c(0.7, 0.6, 0.0, 0.9))
+    )
+  )
+  index.variables = c(1, 5, 100, 10000)
+
+  # 1 would not go out of bounds, thus should always pass
+  testthat::expect_length(mean.index(1, list.res, index.variables), length(index.variables))
+  testthat::expect_length(mean.index(2, list.res, index.variables), length(index.variables))
+  testthat::expect_length(mean.index(3, list.res, index.variables), length(index.variables))
+  testthat::expect_length(mean.index(4, list.res, index.variables), length(index.variables))
+  # 5 is out of bounds of index.variables, returning only NA of correct length
+  testthat::expect_length(mean.index(5, list.res, index.variables), length(index.variables))
+})
+
+test_that("mean.index returns NA vector when i is out of bounds of index.variables", {
+  list.res <- list()
+  index.variables = 1:4
+  testthat::expect_equal(mean.index(5, list.res, index.variables), rep(NA, length(index.variables)))
+})