Fixes bgmCompare output handling.

Fixes summary, print, and coef functions.

Author: MaartenMarsman

R/bgmcompare-methods.r

@@ -145,26 +145,54 @@ print.summary.bgmCompare = function(x, digits = 3, ...) {
     print(ind, row.names = FALSE)
     if (nrow(x$indicator) > 6)
       cat("... (use `summary(fit)$indicator` to see full output)\n")
-    cat("\n")
     cat("Note: NA values are suppressed in the print table. They occur when an indicator\n")
     cat("was constant (all 0 or all 1) across all iterations, so sd/mcse/n_eff/Rhat\n")
     cat("are undefined; `summary(fit)$indicator` still contains the NA values.\n\n")
   }
 
-
   if (!is.null(x$main_diff)) {
     cat("Group differences (main effects):\n")
-    print_df(x$main_diff, digits)
+
+    maind <- head(x$main_diff, 6)
+
+    # Only round numeric columns
+    is_num <- vapply(maind, is.numeric, logical(1))
+    maind[is_num] <- lapply(maind[is_num],
+                            function(col) ifelse(is.na(col), "", round(col, digits)))
+
+    print(maind, row.names = FALSE)
+
     if (nrow(x$main_diff) > 6)
       cat("... (use `summary(fit)$main_diff` to see full output)\n")
+
+    if (!is.null(x$indicator)) {
+      cat("Note: NA values are suppressed in the print table. They occur here when an\n")
+      cat("indicator was zero across all iterations, so mcse/n_eff/Rhat are undefined;\n")
+      cat("`summary(fit)$main_diff` still contains the NA values.\n")
+    }
     cat("\n")
   }
 
   if (!is.null(x$pairwise_diff)) {
     cat("Group differences (pairwise effects):\n")
-    print_df(x$pairwise_diff, digits)
+
+    pairwised <- head(x$pairwise_diff, 6)
+
+    # Only round numeric columns
+    is_num <- vapply(pairwised, is.numeric, logical(1))
+    pairwised[is_num] <- lapply(pairwised[is_num],
+                                function(col) ifelse(is.na(col), "", round(col, digits)))
+
+    print(pairwised, row.names = FALSE)
+
     if (nrow(x$pairwise_diff) > 6)
       cat("... (use `summary(fit)$pairwise_diff` to see full output)\n")
+
+    if (!is.null(x$indicator)) {
+      cat("Note: NA values are suppressed in the print table. They occur here when an\n")
+      cat("indicator was zero across all iterations, so mcse/n_eff/Rhat are undefined;\n")
+      cat("`summary(fit)$pairwise_diff` still contains the NA values.\n")
+    }
     cat("\n")
   }
 
@@ -204,29 +232,36 @@ coef.bgmCompare <- function(object, ...) {
   is_ordinal     <- as.logical(args$is_ordinal_variable)
   num_groups     <- as.integer(args$num_groups)
   num_variables  <- as.integer(args$num_variables)
-  projection     <- args$projection   # matrix [num_groups x (num_groups-1)]
+  projection     <- args$projection   # [num_groups x (num_groups-1)]
 
-  # helper: combine chains into array3d [iter, chain, param]
+  # ---- helper: combine chains into [iter, chain, param], robust to vectors/1-col
   to_array3d <- function(xlist) {
-    nchains <- length(xlist)
-    niter   <- nrow(xlist[[1]])
-    nparam  <- ncol(xlist[[1]])
-    arr <- array(NA_real_, dim = c(niter, nchains, nparam))
-    for (c in seq_len(nchains)) arr[, c, ] <- xlist[[c]]
+    if (is.null(xlist)) return(NULL)
+    stopifnot(length(xlist) >= 1)
+    mats <- lapply(xlist, function(x) {
+      m <- as.matrix(x)
+      if (is.null(dim(m))) m <- matrix(m, ncol = 1L)
+      m
+    })
+    niter  <- nrow(mats[[1]])
+    nparam <- ncol(mats[[1]])
+    arr <- array(NA_real_, dim = c(niter, length(mats), nparam))
+    for (c in seq_along(mats)) arr[, c, ] <- mats[[c]]
     arr
   }
 
   # ============================================================
   # ---- main effects ----
   array3d_main <- to_array3d(object$raw_samples$main)
-  mean_main    <- apply(array3d_main, 3, mean)
+  stopifnot(!is.null(array3d_main))
+  mean_main <- apply(array3d_main, 3, mean)
 
-  num_main <- length(mean_main) / num_groups
-  main_mat <- matrix(mean_main,
-                     nrow = num_main, ncol = num_groups,
-                     byrow = FALSE)
+  stopifnot(length(mean_main) %% num_groups == 0L)
+  num_main <- as.integer(length(mean_main) / num_groups)
 
-  # row names
+  main_mat <- matrix(mean_main, nrow = num_main, ncol = num_groups, byrow = FALSE)
+
+  # row names in sampler row order
   rownames(main_mat) <- unlist(lapply(seq_len(num_variables), function(v) {
     if (is_ordinal[v]) {
       paste0(var_names[v], "(c", seq_len(num_categories[v]), ")")
@@ -235,15 +270,13 @@ coef.bgmCompare <- function(object, ...) {
         paste0(var_names[v], "(quadratic)"))
     }
   }))
+  colnames(main_mat) <- c("baseline", paste0("diff", seq_len(num_groups - 1L)))
 
-  # column names: baseline + diffs
-  colnames(main_mat) <- c("baseline", paste0("diff", seq_len(num_groups - 1)))
-
-  # compute group effects
+  # group-specific main effects: baseline + P %*% diffs
   main_effects_groups <- matrix(NA_real_, nrow = num_main, ncol = num_groups)
   for (r in seq_len(num_main)) {
     baseline <- main_mat[r, 1]
-    diffs    <- main_mat[r, -1]
+    diffs    <- main_mat[r, -1, drop = TRUE]
     main_effects_groups[r, ] <- baseline + as.vector(projection %*% diffs)
   }
   rownames(main_effects_groups) <- rownames(main_mat)
@@ -252,51 +285,63 @@ coef.bgmCompare <- function(object, ...) {
   # ============================================================
   # ---- pairwise effects ----
   array3d_pair <- to_array3d(object$raw_samples$pairwise)
-  mean_pair    <- apply(array3d_pair, 3, mean)
+  stopifnot(!is.null(array3d_pair))
+  mean_pair <- apply(array3d_pair, 3, mean)
+
+  stopifnot(length(mean_pair) %% num_groups == 0L)
+  num_pair <- as.integer(length(mean_pair) / num_groups)
 
-  num_pair <- length(mean_pair) / num_groups
-  pairwise_mat <- matrix(mean_pair,
-                         nrow = num_pair, ncol = num_groups,
-                         byrow = FALSE)
+  pairwise_mat <- matrix(mean_pair, nrow = num_pair, ncol = num_groups, byrow = FALSE)
 
-  # row names
+  # row names in sampler row order (upper-tri i<j)
   pair_names <- character()
-  for (i in 1:(num_variables - 1)) {
-    for (j in (i + 1):num_variables) {
-      pair_names <- c(pair_names, paste0(var_names[i], "-", var_names[j]))
+  if (num_variables >= 2L) {
+    for (i in 1L:(num_variables - 1L)) {
+      for (j in (i + 1L):num_variables) {
+        pair_names <- c(pair_names, paste0(var_names[i], "-", var_names[j]))
+      }
     }
   }
   rownames(pairwise_mat) <- pair_names
-  colnames(pairwise_mat) <- c("baseline", paste0("diff", seq_len(num_groups - 1)))
+  colnames(pairwise_mat) <- c("baseline", paste0("diff", seq_len(num_groups - 1L)))
 
-  # compute group effects
+  # group-specific pairwise effects
   pairwise_effects_groups <- matrix(NA_real_, nrow = num_pair, ncol = num_groups)
   for (r in seq_len(num_pair)) {
     baseline <- pairwise_mat[r, 1]
-    diffs    <- pairwise_mat[r, -1]
+    diffs    <- pairwise_mat[r, -1, drop = TRUE]
     pairwise_effects_groups[r, ] <- baseline + as.vector(projection %*% diffs)
   }
   rownames(pairwise_effects_groups) <- rownames(pairwise_mat)
   colnames(pairwise_effects_groups) <- paste0("group", seq_len(num_groups))
 
   # ============================================================
-  # ---- indicators ----
+  # ---- indicators (present only if selection was used) ----
+  indicators <- NULL
   array3d_ind <- to_array3d(object$raw_samples$indicator)
-  mean_ind    <- apply(array3d_ind, 3, mean)
-
-  indicators <- matrix(0, nrow = num_variables, ncol = num_variables,
-                       dimnames = list(var_names, var_names))
-
-  diag(indicators) <- mean_ind[seq_len(num_variables)]
-
-  counter <- num_variables + 1
-  for (i in 1:(num_variables - 1)) {
-    for (j in (i + 1):num_variables) {
-      val <- mean_ind[counter]
-      indicators[i, j] <- val
-      indicators[j, i] <- val
-      counter <- counter + 1
+  if (!is.null(array3d_ind)) {
+    mean_ind <- apply(array3d_ind, 3, mean)
+
+    # reconstruct VxV matrix using the sampler’s interleaved order:
+    # (1,1),(1,2),...,(1,V),(2,2),...,(2,V),...,(V,V)
+    V <- num_variables
+    stopifnot(length(mean_ind) == V * (V + 1L) / 2L)
+
+    ind_mat <- matrix(0, nrow = V, ncol = V,
+                      dimnames = list(var_names, var_names))
+    pos <- 1L
+    for (i in seq_len(V)) {
+      # diagonal (main indicator)
+      ind_mat[i, i] <- mean_ind[pos]; pos <- pos + 1L
+      if (i < V) {
+        for (j in (i + 1L):V) {
+          val <- mean_ind[pos]; pos <- pos + 1L
+          ind_mat[i, j] <- val
+          ind_mat[j, i] <- val
+        }
+      }
     }
+    indicators <- ind_mat
   }
 
   # ============================================================
@@ -308,4 +353,4 @@ coef.bgmCompare <- function(object, ...) {
     pairwise_effects_groups = pairwise_effects_groups,
     indicators              = indicators
   )
-}
+}

R/bgms-methods.R

@@ -100,8 +100,16 @@ print.summary.bgms <- function(x, digits = 3, ...) {
 
   if (!is.null(x$pairwise)) {
     cat("Pairwise interactions:\n")
-    print(round(head(x$pairwise, 6), digits = digits))
+    pair <- head(x$pairwise, 6)
+    pair[] <- lapply(pair, function(col) ifelse(is.na(col), "", round(col, digits)))
+    print(pair)
+    #print(round(head(x$pairwise, 6), digits = digits))
     if (nrow(x$pairwise) > 6) cat("... (use `summary(fit)$pairwise` to see full output)\n")
+    if (!is.null(x$indicator)) {
+      cat("Note: NA values are suppressed in the print table. They occur here when an \n")
+      cat("indicator was zero across all iterations, so mcse/n_eff/Rhat are undefined;\n")
+      cat("`summary(fit)$pairwise` still contains the NA values.\n")
+    }
     cat("\n")
   }
 
@@ -112,7 +120,6 @@ print.summary.bgms <- function(x, digits = 3, ...) {
     print(ind, row.names = FALSE)
     if (nrow(x$indicator) > 6)
       cat("... (use `summary(fit)$indicator` to see full output)\n")
-    cat("\n")
     cat("Note: NA values are suppressed in the print table. They occur when an indicator\n")
     cat("was constant (all 0 or all 1) across all iterations, so sd/mcse/n_eff/Rhat\n")
     cat("are undefined; `summary(fit)$indicator` still contains the NA values.\n\n")