From e4ec55ed3a7bbd4cc6453345fbf653297996f2bc Mon Sep 17 00:00:00 2001
From: Alan O'Callaghan <alan.ocallaghan@outlook.com>
Date: Mon, 11 Mar 2024 16:47:48 +0000
Subject: [PATCH] Add working na.rm option and test

---
 DESCRIPTION                   |  4 ++--
 R/Methods.R                   | 13 ++++++++-----
 tests/testthat/test_methods.R | 11 +++++++++++
 3 files changed, 21 insertions(+), 7 deletions(-)

diff --git a/DESCRIPTION b/DESCRIPTION
index 2135539..0073718 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,8 +1,8 @@
 Package: BASiCS
 Type: Package
 Title: Bayesian Analysis of Single-Cell Sequencing data
-Version: 2.15.2
-Date: 2023-12-22
+Version: 2.15.3
+Date: 2024-03-11
 Authors@R: c(person("Catalina", "Vallejos", role=c("aut", "cre"),
         email="catalina.vallejos@igmm.ed.ac.uk",
         comment=c(ORCID = "0000-0003-3638-1960")), 
diff --git a/R/Methods.R b/R/Methods.R
index 38a6c40..f2ca969 100644
--- a/R/Methods.R
+++ b/R/Methods.R
@@ -157,11 +157,14 @@ setMethod("Summary",
         nrow = ncol(n),
         dimnames = list(colnames(n), c("median", "lower", "upper"))
       )
-      HPD[,1] <- colMedians(n)
-      HPD[!is.na(HPD[,1]),2:3] <- coda::HPDinterval(
-        coda::mcmc(n[,!is.na(HPD[,1])]),
-        prob = prob
-      )
+      HPD[, 1] <- colMedians(n, na.rm = na.rm)
+      ind_not_na <- !is.na(HPD[, 1])
+      HPD[ind_not_na, 2:3] <- apply(n[, ind_not_na], 2, function(col) {
+        if (na.rm) {
+          col <- na.omit(col)
+        }
+        coda::HPDinterval(coda::mcmc(col), prob = prob)
+      })
       HPD
     }
   )
diff --git a/tests/testthat/test_methods.R b/tests/testthat/test_methods.R
index b0191db..70ded1c 100644
--- a/tests/testthat/test_methods.R
+++ b/tests/testthat/test_methods.R
@@ -10,6 +10,17 @@ test_that("show", {
   expect_output(show(v), "An object of class BASiCS_ResultVG.")
 })
 
+test_that("Summary na.rm", {
+  data(ChainSC)
+  ChainSC@parameters$nu[1, 1] <- NA
+  ChainSC@parameters$nu[1, 4] <- NA
+  s1 <- Summary(ChainSC, na.rm = TRUE)
+  s2 <- Summary(ChainSC, na.rm = FALSE)
+  expect_false(any(is.na(s1@parameters$nu)))
+  expect_true(all(is.na(s2@parameters$nu[1, ])))
+  expect_true(all(is.na(s2@parameters$nu[4, ])))
+})
+
 test_that("subset", {
   data(ChainSC)
   data(ChainSCReg)