block diagonal covariance matrix by group of variable

Question

If this is my dataset, 3 subjects, measured at two time points t1,t2 , at each time point each subject is measured twice

     ID         Time     Observation    Y
     1          t1       1              3.1
     1          t1       2              4.5
     1          t2       1              4.2
     1          t2       2              1.7

     2          t1       1              2.3
     2          t1       2              3.1
     2          t2       1              1.5
     2          t2       2              2.2

     3          t1       1              4.1
     3          t1       2              4.9
     3          t2       1              3.5
     3          t2       2              3.2

What is the block diagonal covariance matrix for this dataset ? I am assuming the matrix would have three blocks, each block representing 2x2 variance-covariance matrix of subject i with t1 and t2. Is this correct ?

Answer 1

The code below computes the covariance matrices for each ID, then creates a block matrix out of them.

make_block_matrix <- function(x, y, fill = NA) {
  u <- matrix(fill, nrow = dim(x)[1], ncol = dim(y)[2])
  v <- matrix(fill, nrow = dim(y)[1], ncol = dim(x)[2])
  u <- cbind(x, u)
  v <- cbind(v, y)
  rbind(u, v)
}

# compute covariance matrices for each ID
sp <- split(df1, df1$ID)
cov_list <- lapply(sp, \(x) {
  y <- reshape(x[-1], direction = "wide", idvar = "Time", timevar = "Observation")
  cov(y[-1])
})
# no longer needed, tidy up
rm(sp)

# fill with NA's, simple call
# res <- Reduce(make_block_matrix, cov_list)

# fill with zeros
res <- Reduce(\(x, y) make_block_matrix(x, y, fill = 0), cov_list)
# finally, the results' row and column names
nms <- paste(sapply(cov_list, rownames), names(cov_list), sep = ".")
dimnames(res) <- list(nms, nms)
res
#>        Y.1.1 Y.2.2 Y.1.3 Y.2.1 Y.1.2 Y.2.3
#> Y.1.1  0.605 -1.54  0.00 0.000  0.00 0.000
#> Y.2.2 -1.540  3.92  0.00 0.000  0.00 0.000
#> Y.1.3  0.000  0.00  0.32 0.360  0.00 0.000
#> Y.2.1  0.000  0.00  0.36 0.405  0.00 0.000
#> Y.1.2  0.000  0.00  0.00 0.000  0.18 0.510
#> Y.2.3  0.000  0.00  0.00 0.000  0.51 1.445

^{Created on 2023-02-20 with reprex v2.0.2}

---

Data

df1 <- "ID         Time     Observation    Y
     1          t1       1              3.1
     1          t1       2              4.5
     1          t2       1              4.2
     1          t2       2              1.7

     2          t1       1              2.3
     2          t1       2              3.1
     2          t2       1              1.5
     2          t2       2              2.2

     3          t1       1              4.1
     3          t1       2              4.9
     3          t2       1              3.5
     3          t2       2              3.2"
df1 <- read.table(text = df1, header = TRUE)

^{Created on 2023-02-20 with reprex v2.0.2}