New gpuCorPerm(): permutation p-values for every
correlation between the columns of x and y,
and p-values adjusted for testing them all by Westfall and Young’s maxT.
Pearson or Spearman.
New gpuRowTtest(): Welch’s or Student’s t test of
every row of a matrix between two groups of columns, as
t.test gives them, with permutation and maxT p-values
computed on the GPU.
New gpuKmeans(): Lloyd’s k-means on the GPU, step
for step as kmeans(algorithm = "Lloyd"), from the same
random starts, returning a "kmeans" object. A million
points into a thousand clusters takes seconds.
New gpuAssoc(): an association scan of every column
of x against every column of y (SNPs against
genes, in an eQTL study), linear regression with covariates, reporting
the associations below a p-value or each outcome’s most significant
ones, with slopes, t statistics and false discovery rates. It computes
what Matrix eQTL’s linear model does, on the GPU, and agrees with
lm().
New gpuCorPairs(): the strongly correlated pairs
among many variables (every pair over a threshold, or each variable’s
strongest partners), Pearson or Spearman, without forming the
correlation matrix, which it computes on the GPU a block at a time. For
tens of thousands of variables, such as a gene co-expression
network.
The benchmarks vignette now says at what data size each function
starts to beat base R, and what a loop of small calls costs;
bench/crossover.R repeats the measurement.
Every method was checked against its published or standard
definition; review/algorithms.md in the source repository
has the details. Fixed:
gpuLm, gpuLsfit and gpuGlm
did not detect collinear predictors: where lm() reports a
coefficient as NA, they returned huge values. They now use
the QR decomposition lm() uses (LINPACK’s
dqrdc2) and alias the same columns.gpuLm, gpuLsfit, gpuGlm and
gpuQr treated data in small units (around 1e-8) as rank 0,
returning NA coefficients; their rank tests are now
relative.gpuCor(method = "kendall") was wrong whenever the data
had ties; it now computes tau-b, as cor() does, with exact
pair counts.gpuHclust(method = "flexible group") computed single
linkage; it now computes flexible UPGMA.gpuCor()$ts was 10000 for every correlation of 0.999 or
more, and the t-scores lost precision in double; gpuTtest
gave NaN at |r| = 1.gpuGranger’s result matrices are labelled
cause x effect, and their diagonal is
NA (a series against itself is no test) rather than F = 0,
p = 1.New:
gpuSolve computes on the GPU with NVIDIA’s cuSOLVER: LU
decomposition for a square x, as solve() does,
and QR least squares for an overdetermined one. It used to factor
x with base R on the CPU, which made it slower than
solve(). y may now be a matrix, and results
are named as solve() names them. Building gputools now also
needs cuSOLVER, which comes with the CUDA toolkit.gpuCor computes Spearman’s coefficient
(method = "spearman"), and Kendall’s now honours
use = "pairwise.complete.obs" (it warned and ignored it).
With test = TRUE every method returns p-values
($p.value) as cor.test() gives them – exactly
for Pearson’s, as cor.test(exact = FALSE) for the rank
coefficients – and Kendall’s its z statistics ($z). Results
carry cor()’s dimnames.use = "everything") are one matrix product on the
GPU, and Spearman’s columns are ranked on the GPU too: 25 million
correlations (500 x 5,000) in about 0.5 s, where cor()
takes 7 s.gpuPrcomp and gpuSvd: principal components
and the singular value decomposition as prcomp() and
svd() compute them, on cuSOLVER; gpuPrcomp
returns a "prcomp" object.gpuHclust and gpuDistClust take
beta, the parameter of the "flexible" and
"flexible group" methods, which was fixed at 0.5. The
default stays 0.5; Lance and Williams (1967) recommend -0.25.These gave wrong answers in earlier versions:
gpuDist(method = "canberra") used |a - b| / |a + b|; it
now uses R’s |a - b| / (|a| + |b|), leaving out 0/0 terms as
dist() does.gpuHclust(method = "ward") swapped the two cluster
sizes in its update, so merges of unequal-sized clusters got the wrong
height.gpuGranger(x, y) with a single column in x
tested every y column but the last against the wrong
series.gpuQr chose its pivots from wrongly updated column
norms, so it could report too low a rank.gpuGranger p-values below about 1e-7 were rounded to
0.gpuGranger documentation described
fStatistics the other way round from what the code
computes; it now says what the code does: fStatistics[i, j]
tests whether x[, i] Granger-causes x[, j] (or
y[, j]).And these crashed or failed:
qr.X() and other base-R functions given a
gpuQr() result could corrupt memory and crash R;
gpuTtest, gpuSolve and gpuHclust
results carried the same hazard.gpuTtest could not be called at all, and
gpuGlm stopped with “argument 16 is empty” (#23).gpuMatMult read integer matrices as doubles.gpuDistClust leaked the GPU memory of its distance
matrix.precision = "double", or
options(gputools.precision = "double") for the session.
useSingle = FALSE in gpuLm,
gpuLsfit and gpuGlm, which used to stop with
“not yet implemented”, now works.gpuAvailable() says whether a usable CUDA device is
present.configure finds CUDA through CUDA_HOME,
CUDA_PATH, nvcc on the PATH,
conda environments or /usr/local/cuda.gpuHclust on 3000 points went from 26 s to 0.3 s..Call and registered
routines; about 1500 lines of unused code are gone.added some linear and nonlinear regression functions
added some configure args –with-r-lib and –with-r-include (see INSTALL notes section)
fixed a bug in gpuMi
added a demo of CULA use: fastICA with CULA’s svd plugged in.
Mac os x and nvidia cards with compute capability less than 1.3 are now supported. Systems lacking a gpu altogether also receive some support. See the file INSTALL for details.
gpuHclust, gpuDist, and gpuQr now output the same classes as R’s base hclust, dist, and qr functions.
The c code now uses r’s api. Host memory allocation and error reporting have received a facelift.
gpuGranger now reports the correct p-values and is slightly faster for small data sets.
Some bugs with gpuCor’s pearson method have been fixed.