prototyping.Rmd
The distcomp
package allows one to perform computations on data that is row-partitioned data among several sites. Implemented methods include stratified Cox regression and singular value decomposition (SVD). The distributed computations rely on a REST API for R, such as that provided by Opencpu for example.
In order to prototype new methods for distcomp
, it is convenient to write code that does all computations locally, at least during the development phase. This allows for easy debugging and code changes. The step to using a REST call is merely one step off. For example, in the context of distributed multiple imputation, the distributed SVD is a building block. So one has to be able to be able to run the svd operation not only in a modular fashion, but also quickly to test code.
Previous versions of distcomp
assumed all computations were done using opencpu
remote calls. This meant that the results were time-consuming to check; for reasonably sized problems one would have to wait for hours.
Version 1.1 of distcomp
addresses this by allowing local computations, in essence faking the opencpu
calls. All calculations are done locally at normal speeds (rather than a https POST
request) while ensuring that the sites are only privy to their own data.
## [1] "QueryCount" "StratifiedCoxModel" "RankKSVD"
We will perform an SVD on a 1500 by 20 matrix where the matrix is distributed across three sites each having a 500 by 20 matrix. We will extract the first 10 components using this distributed algorithm.
We first define the computation.
svdDef <- data.frame(compType = "RankKSVD",
rank = 20L,
ncol = 20L,
id = "SVD",
stringsAsFactors = FALSE)
We generate some synthetic data.
set.seed(12345)
## Three sites
nSites <- 3
siteData <- lapply(seq.int(nSites), function(i) matrix(rnorm(10000), ncol=20))
We now create local objects to represent the sites, a list of three items, each with a name and its specific data.
sites <- lapply(seq.int(nSites),
function(x) list(name = paste0("site", x),
worker = makeWorker(defn = svdDef, data = siteData[[x]])
))
We are now ready to create a master object and add sites to it.
master <- makeMaster(svdDef)
for (site in sites) {
master$addSite(name = site$name, worker = site$worker)
}
When any of the sites specify a worker
argument in the call to addSite
, it is an indication that the computation is to be performed locally. Otherwise, a URL will have to be specified using the url
argument to addSite
. (Only one of worker
or url
can be specified.)
Finally, we are ready to run the decomposition, specifying a maximum number of iterations, so that things don’t go haywire. This takes less than 5 seconds on my laptop.
## user system elapsed
## 6.124 3.543 8.325
The result is a named list of two things, the \(d\) diagonal values and the \(v\) matrix.
We can compare the results obtained here with the actual.
Let’s print it side-by-side.
truth | distcomp |
---|---|
41.82263 | 41.82263 |
41.65242 | 41.65242 |
41.32027 | 41.32027 |
40.88629 | 40.88629 |
40.41849 | 40.41849 |
40.27631 | 40.27631 |
40.12698 | 40.12698 |
39.44775 | 39.44775 |
39.24815 | 39.24815 |
39.07074 | 39.07074 |
38.37076 | 38.37076 |
38.21350 | 38.21350 |
37.61047 | 37.61047 |
37.14286 | 37.14286 |
36.98763 | 36.98763 |
36.24778 | 36.24778 |
36.04962 | 36.04962 |
35.66050 | 35.66050 |
35.21343 | 35.21343 |
34.56329 | 34.56329 |
We can also compare the \(v\) matrix, taking into account the signs may be different.
## [1] 1.124008e-05
The approach for the stratified Cox model is similar. We start with a definition.
coxDef <- data.frame(compType = "StratifiedCoxModel",
formula = "Surv(time, censor) ~ age + becktota + ndrugfp1 + ndrugfp2 + ivhx3 + race + treat",
projectName = "STCoxTest",
projectDesc = "STCox Project Desc",
stringsAsFactors = FALSE)
We now create two sites with appropriate data
## Two sites
siteDataFiles <- file.path(system.file("ex", package="distcomp"), c("uis-site1.csv", "uis-site2.csv"))
siteData <- lapply(siteDataFiles, read.csv)
We now create local objects to represent the two sites, each with a name and its specific data.
sites <- lapply(seq_along(siteData),
function(x) list(name = paste0("site", x),
worker = makeWorker(defn = coxDef, data = siteData[[x]])
))
We are now ready to create a master object and add sites to it.
master <- makeMaster(coxDef)
for (site in sites) {
master$addSite(name = site$name, worker = site$worker)
}
Now we can run the model.
## coef exp(coef) se(coef) z p
## 1 -0.028076 0.9723 0.008131 -3.453 5.542e-04
## 2 0.009146 1.0092 0.004991 1.832 6.691e-02
## 3 -0.521973 0.5933 0.124424 -4.195 2.727e-05
## 4 -0.194178 0.8235 0.048252 -4.024 5.717e-05
## 5 0.263634 1.3017 0.108243 2.436 1.487e-02
## 6 -0.240021 0.7866 0.115632 -2.076 3.792e-02
## 7 -0.212616 0.8085 0.093747 -2.268 2.333e-02
The above result will be the same as what we obtain from a Cox regression fit.
coxFit <- survival::coxph(formula=Surv(time, censor) ~ age + becktota + ndrugfp1 + ndrugfp2 +
ivhx3 + race + treat + strata(site),
data = rbind(siteData[[1]], siteData[[2]]))
summary(coxFit)$coefficients
## coef exp(coef) se(coef) z Pr(>|z|)
## age -0.028075893 0.9723146 0.008130685 -3.453078 5.542280e-04
## becktota 0.009145528 1.0091875 0.004991421 1.832250 6.691425e-02
## ndrugfp1 -0.521973047 0.5933487 0.124423881 -4.195119 2.727278e-05
## ndrugfp2 -0.194177573 0.8235117 0.048252289 -4.024215 5.716573e-05
## ivhx3TRUE 0.263634281 1.3016521 0.108243388 2.435569 1.486837e-02
## race -0.240020862 0.7866115 0.115632433 -2.075723 3.791961e-02
## treat -0.212616368 0.8084662 0.093747124 -2.267978 2.333058e-02
## R version 3.6.2 (2019-12-12)
## Platform: x86_64-apple-darwin19.2.0 (64-bit)
## Running under: macOS Catalina 10.15.2
##
## Matrix products: default
## BLAS/LAPACK: /usr/local/Cellar/openblas/0.3.7/lib/libopenblasp-r0.3.7.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] distcomp_1.3 survival_3.1-8
##
## loaded via a namespace (and not attached):
## [1] gmp_0.5-13.5 Rcpp_1.0.3 highr_0.8 compiler_3.6.2
## [5] pillar_1.4.3 later_1.0.0 tools_3.6.2 digest_0.6.23
## [9] jsonlite_1.6 evaluate_0.14 memoise_1.1.0 tibble_2.1.3
## [13] lattice_0.20-38 pkgconfig_2.0.3 rlang_0.4.2 Matrix_1.2-18
## [17] shiny_1.4.0 yaml_2.2.0 pkgdown_1.4.1 xfun_0.11
## [21] fastmap_1.0.1 httr_1.4.1 stringr_1.4.0 dplyr_0.8.3
## [25] knitr_1.26 desc_1.2.0 fs_1.3.1 rprojroot_1.3-2
## [29] grid_3.6.2 tidyselect_0.2.5 glue_1.3.1 R6_2.4.1
## [33] rmarkdown_2.0 homomorpheR_0.2-4 purrr_0.3.3 magrittr_1.5
## [37] backports_1.1.5 promises_1.1.0 htmltools_0.4.0 MASS_7.3-51.5
## [41] splines_3.6.2 assertthat_0.2.1 mime_0.8 xtable_1.8-4
## [45] httpuv_1.5.2 stringi_1.4.3 sodium_1.1 crayon_1.3.4