Browse code

some more diagnostics

Greg Finak authored on 18/09/2020 22:06:41
Showing 4 changed files

... ...
@@ -1,7 +1,7 @@
1 1
 Package: COMPASS
2 2
 Type: Package
3 3
 Title: Combinatorial Polyfunctionality Analysis of Single Cells
4
-Version: 1.27.1
4
+Version: 1.27.2
5 5
 Date: 2014-07-11
6 6
 Author: Lynn Lin, Kevin Ushey, Greg Finak, Ravio Kolde (pheatmap)
7 7
 Authors@R: c( person("Lynn", "Lin", role="aut", email="[email protected]"),
... ...
@@ -62,4 +62,4 @@ LazyLoad: yes
62 62
 LazyData: yes
63 63
 biocViews: ImmunoOncology, FlowCytometry
64 64
 Encoding: UTF-8
65
-RoxygenNote: 7.1.0
65
+RoxygenNote: 7.1.1
... ...
@@ -29,6 +29,7 @@ export(COMPASS)
29 29
 export(COMPASSContainer)
30 30
 export(COMPASSContainerFromGatingSet)
31 31
 export(COMPASSDescription)
32
+export(COMPASSMCMCDiagnosis)
32 33
 export(COMPASSfitToCountsTable)
33 34
 export(CellCounts)
34 35
 export(Combinations)
... ...
@@ -154,3 +154,25 @@ checkCOMPASSConvergence<-function(mlist,ncores=1){
154 154
   }
155 155
 }
156 156
 
157
+#'@title Diagnostic of a set of COMPASS Models.
158
+#' @param x a list of compass model fits of the same data with the same number of iterations, different seeds.
159
+#' Run some mcmc diagnostics on a series of COMPASS model fits.
160
+#' Assuming the input is a list of model fits for the same data with the same number of iterations and different seeds.
161
+#' Run Gelman's Rhat diagnostics on the alpha_s and alpha_u hyperparameter chains, treating each model as an independent chain.
162
+#' Rhat should be near 1 but rarely are in practice. Very large values may be a concern.
163
+#' The method returns an average model, by averaging the mean_gamma matrices (equally weighted since each input has the same number of iterations).
164
+#' This mean model should be better then any of the individual models.
165
+#' It can be plotted via "plot(result$mean_model)".
166
+#' @export
167
+COMPASSMCMCDiagnosis<-function(x){
168
+    require(coda)
169
+    require(abind)
170
+    diag<-list()
171
+    diag$alpha_s<-coda::gelman.diag(Map(function(x)coda::as.mcmc(x$fit$alpha_s),x))
172
+    diag$alpha_u<-coda::gelman.diag(Map(function(x)coda::as.mcmc(x$fit$alpha_u),x))
173
+  mean_gamma <- apply(Map(function(x) abind(x, along = 3), Map(function(x) x$fit$mean_gamma, x))[[1]], 1:2, mean)
174
+  mean_model <- x[[1]]
175
+  mean_model$fit$mean_gamma <- mean_gamma
176
+    return(list(diag=diag,mean_model=mean_model))
177
+}
178
+
157 179
new file mode 100644
... ...
@@ -0,0 +1,21 @@
1
+% Generated by roxygen2: do not edit by hand
2
+% Please edit documentation in R/Rhat.R
3
+\name{COMPASSMCMCDiagnosis}
4
+\alias{COMPASSMCMCDiagnosis}
5
+\title{Diagnostic of a set of COMPASS Models.}
6
+\usage{
7
+COMPASSMCMCDiagnosis(x)
8
+}
9
+\arguments{
10
+\item{x}{a list of compass model fits of the same data with the same number of iterations, different seeds.
11
+Run some mcmc diagnostics on a series of COMPASS model fits.
12
+Assuming the input is a list of model fits for the same data with the same number of iterations and different seeds.
13
+Run Gelman's Rhat diagnostics on the alpha_s and alpha_u hyperparameter chains, treating each model as an independent chain.
14
+Rhat should be near 1 but rarely are in practice. Very large values may be a concern.
15
+The method returns an average model, by averaging the mean_gamma matrices (equally weighted since each input has the same number of iterations).
16
+This mean model should be better then any of the individual models.
17
+It can be plotted via "plot(result$mean_model)".}
18
+}
19
+\description{
20
+Diagnostic of a set of COMPASS Models.
21
+}