Case of N GWAS = 2

[danielcoral]

There seems to be a bug if attempting to run a cross-trait GWAS with only 2 traits. It seems to be in step3, in lines 161-175, when doing:

tresm <- foreach(i = 1:nrow(pairma),
                   .combine = "rbind",
                   .inorder = T) %dopar%
    CorE.ICE(pairma[i, 1], pairma[i, 2], resinfm, maSpltRw, maSpltN, cgwasenv)
  tresm <- signif(tresm, 7)

  corm <- cbind(cgwasenv$.TRAIT_NAME[pairma[,1]],
                cgwasenv$.TRAIT_NAME[pairma[,2]],
                signif(resinfm[pairma[,1], 4]-1, 7),
                signif(resinfm[pairma[,2], 4]-1, 7),
                tresm)

  colnames(corm) <- c("GWAS1", "GWAS2",
                      "T1Eff", "T2Eff",
                      "StatCor", "Psi", "EffCov", "allPi", "T1sEff", "T2sEff", "EffsCov", "sigPi")

It fails at naming the columns because the corm object does not end up with the desired number of columns.

With the following changes the function works:

step3 <- function(cgwasenv) {
  logOutput("========== C-GWAS step 3 : GetPsi ==========\n\n", cgwasenv = cgwasenv)
  logOutput("Estimating background correlation for ",
            nrow(pairma <- t(combn(seq_len(cgwasenv$.TRAIT_NUM), 2))),
            " GWAS pairs ..\n", cgwasenv = cgwasenv)

  minsnpn = 1e5
  maSpltRw = ceiling(cgwasenv$.SNP_N / minsnpn)
  maSpltN = floor(cgwasenv$.SNP_N / maSpltRw) * maSpltRw
  resinfm <- as.matrix(
               read.table(file.path(cgwasenv$.CGWAS_DETAIL_PATH, "SummaryGetI.txt"),
                          header = T)[,-1])

  threadNCur <- min(cgwasenv$.PARAL_NUM, nrow(pairma))
  cl <- makeCluster(threadNCur)
  registerDoParallel(cl)
  # globalVariables(c('i'))
  i <- 1 # assign parallel control variants

  tresm <- foreach(i = 1:nrow(pairma),
                   .combine = "rbind",
                   .inorder = T) %dopar%
    CorE.ICE(pairma[i, 1], pairma[i, 2], resinfm, maSpltRw, maSpltN, cgwasenv)
  tresm <- signif(tresm, 7)
  if(is.null(dim(tresm))){
    tresm <- matrix(tresm, 1)
  }
  corm <- cbind(cgwasenv$.TRAIT_NAME[pairma[,1]],
                cgwasenv$.TRAIT_NAME[pairma[,2]],
                signif(resinfm[pairma[,1], 4]-1, 7),
                signif(resinfm[pairma[,2], 4]-1, 7),
                tresm)
  colnames(corm) <- c("GWAS1", "GWAS2",
                      "T1Eff", "T2Eff",
                      "StatCor", "Psi", "EffCov", "allPi", "T1sEff", "T2sEff", "EffsCov", "sigPi")
  write.table(corm,
              file.path(cgwasenv$.CGWAS_DETAIL_PATH, "BCCorrelationStat.txt"),
              row.names = F, quote = F, sep = "\t")
  write.table(corm[, c(1,2,5,6,8,12), drop = FALSE],
              file.path(cgwasenv$.CGWAS_DETAIL_PATH, "SummaryGetPsi.txt"),
              row.names = F, quote = F, sep = "\t")
  logOutput("Summary of GetPsi written to Details/SummaryGetPsi.txt\n", cgwasenv = cgwasenv)

  stopCluster(cl)

  logOutput("\nC-GWAS step 3 completed\n", cgwasenv = cgwasenv)
}