... | ... |
@@ -1,7 +1,7 @@ |
1 | 1 |
Package: lisaClust |
2 | 2 |
Type: Package |
3 | 3 |
Title: lisaClust: Clustering of Local Indicators of Spatial Association |
4 |
-Version: 1.15.1 |
|
4 |
+Version: 1.15.2 |
|
5 | 5 |
Authors@R: c( |
6 | 6 |
person("Ellis", "Patrick", , "[email protected]", role = c("aut", "cre")), |
7 | 7 |
person("Nicolas", "Canete", , "[email protected]", role = "aut"), |
... | ... |
@@ -49,5 +49,5 @@ Suggests: |
49 | 49 |
rmarkdown, |
50 | 50 |
SpatialDatasets, |
51 | 51 |
testthat (>= 3.0.0) |
52 |
-RoxygenNote: 7.3.1 |
|
52 |
+RoxygenNote: 7.3.2 |
|
53 | 53 |
Config/testthat/edition: 3 |
... | ... |
@@ -1,3 +1,4 @@ |
1 |
+ |
|
1 | 2 |
#' Generate local indicators of spatial association |
2 | 3 |
#' |
3 | 4 |
#' @param cells A SingleCellExperiment, SpatialExperiment or data frame that contains at least the |
... | ... |
@@ -71,11 +72,11 @@ lisa <- |
71 | 72 |
} |
72 | 73 |
|
73 | 74 |
cellSummary <- spicyR:::getCellSummary(cells, bind = FALSE) |
74 |
- |
|
75 |
+ |
|
75 | 76 |
if (is.null(Rs)) { |
76 | 77 |
Rs <- c(20, 50, 100) |
77 | 78 |
} |
78 |
- |
|
79 |
+ |
|
79 | 80 |
BPimage <- BPcellType <- BiocParallel::SerialParam() |
80 | 81 |
if (whichParallel == "imageID") { |
81 | 82 |
BPimage <- BPPARAM |
... | ... |
@@ -83,10 +84,10 @@ lisa <- |
83 | 84 |
if (whichParallel == "cellType") { |
84 | 85 |
BPcellType <- BPPARAM |
85 | 86 |
} |
86 |
- |
|
87 |
- |
|
87 |
+ |
|
88 |
+ |
|
88 | 89 |
message("Generating local L-curves.") |
89 |
- |
|
90 |
+ |
|
90 | 91 |
curveList <- |
91 | 92 |
BiocParallel::bplapply( |
92 | 93 |
cellSummary, |
... | ... |
@@ -99,11 +100,11 @@ lisa <- |
99 | 100 |
lisaFunc = lisaFunc, |
100 | 101 |
BPPARAM = BPimage |
101 | 102 |
) |
102 |
- |
|
103 |
+ |
|
103 | 104 |
curvelist <- lapply(curveList, as.data.frame) |
104 | 105 |
curves <- as.matrix(dplyr::bind_rows(curvelist)) |
105 | 106 |
rownames(curves) <- as.character(unlist(lapply(cellSummary, function(x) x$cellID))) |
106 |
- |
|
107 |
+ |
|
107 | 108 |
curves[is.na(curves)] <- 0 |
108 | 109 |
return(curves) |
109 | 110 |
} |
... | ... |
@@ -120,7 +121,7 @@ pppGenerate <- function(cells, window, window.length) { |
120 | 121 |
window = ow, |
121 | 122 |
marks = cells$cellType |
122 | 123 |
) |
123 |
- |
|
124 |
+ |
|
124 | 125 |
pppCell |
125 | 126 |
} |
126 | 127 |
|
... | ... |
@@ -133,7 +134,7 @@ makeWindow <- |
133 | 134 |
data <- data.frame(data) |
134 | 135 |
ow <- |
135 | 136 |
spatstat.geom::owin(xrange = range(data$x), yrange = range(data$y)) |
136 |
- |
|
137 |
+ |
|
137 | 138 |
if (window == "convex") { |
138 | 139 |
p <- spatstat.geom::ppp(data$x, data$y, ow) |
139 | 140 |
ow <- spatstat.geom::convexhull(p) |
... | ... |
@@ -155,8 +156,8 @@ makeWindow <- |
155 | 156 |
})) |
156 | 157 |
ch <- |
157 | 158 |
concaveman::concaveman(bigDat, |
158 |
- length_threshold = window.length, |
|
159 |
- concavity = 1 |
|
159 |
+ length_threshold = window.length, |
|
160 |
+ concavity = 1 |
|
160 | 161 |
) |
161 | 162 |
poly <- as.data.frame(ch[nrow(ch):1, ]) |
162 | 163 |
colnames(poly) <- c("x", "y") |
... | ... |
@@ -186,7 +187,7 @@ borderEdge <- function(X, maxD) { |
186 | 187 |
areas <- unlist(lapply(circs, spatstat.geom::area)) / (pi * maxD^2) |
187 | 188 |
e[inB] <- areas |
188 | 189 |
} |
189 |
- |
|
190 |
+ |
|
190 | 191 |
e |
191 | 192 |
} |
192 | 193 |
|
... | ... |
@@ -213,13 +214,13 @@ generateCurves <- |
213 | 214 |
window = ow, |
214 | 215 |
marks = data$cellType |
215 | 216 |
) |
216 |
- |
|
217 |
+ |
|
217 | 218 |
if (!is.null(sigma)) { |
218 | 219 |
d <- spatstat.explore::density.ppp(p1, sigma = sigma) |
219 | 220 |
d <- d / mean(d) |
220 | 221 |
} |
221 |
- |
|
222 |
- |
|
222 |
+ |
|
223 |
+ |
|
223 | 224 |
locIJ <- |
224 | 225 |
BiocParallel::bplapply(as.list(levels(p1$marks)), function(j) { |
225 | 226 |
locI <- lapply(as.list(levels(p1$marks)), function(i) { |
... | ... |
@@ -227,7 +228,7 @@ generateCurves <- |
227 | 228 |
jID <- data$cellID[p1$marks == j] |
228 | 229 |
locR <- matrix(NA, length(iID), length(Rs)) |
229 | 230 |
rownames(locR) <- iID |
230 |
- |
|
231 |
+ |
|
231 | 232 |
if (length(jID) > 1 & length(iID) > 1) { |
232 | 233 |
if (!is.null(sigma)) { |
233 | 234 |
dFrom <- d * (sum(p1$marks == i) - 1) / spatstat.geom::area(ow) |
... | ... |
@@ -266,7 +267,7 @@ generateCurves <- |
266 | 267 |
} |
267 | 268 |
colnames(locR) <- |
268 | 269 |
paste(j, round(Rs, 2), sep = "_") |
269 |
- |
|
270 |
+ |
|
270 | 271 |
locR |
271 | 272 |
}) |
272 | 273 |
do.call("rbind", locI) |
... | ... |
@@ -279,16 +280,16 @@ sqrtVar <- function(x) { |
279 | 280 |
len <- 1000 |
280 | 281 |
lambda <- (seq(1, 300, length.out = len) / 100)^x |
281 | 282 |
mL <- max(lambda) |
282 |
- |
|
283 |
- |
|
283 |
+ |
|
284 |
+ |
|
284 | 285 |
V <- NULL |
285 | 286 |
for (i in 1:len) { |
286 | 287 |
V[i] <- var(sqrt(rpois(10000, lambda[i]))) |
287 | 288 |
} |
288 |
- |
|
289 |
+ |
|
289 | 290 |
lambda <- lambda^(1 / x) |
290 | 291 |
V <- V |
291 |
- |
|
292 |
+ |
|
292 | 293 |
f <- loess(V ~ lambda, span = 0.1) |
293 | 294 |
} |
294 | 295 |
|
... | ... |
@@ -297,25 +298,25 @@ sqrtVar <- function(x) { |
297 | 298 |
#' @importFrom spatstat.geom nearest.valid.pixel area marks |
298 | 299 |
weightCounts <- function(dt, X, maxD, lam) { |
299 | 300 |
maxD <- as.numeric(as.character(maxD)) |
300 |
- |
|
301 |
+ |
|
301 | 302 |
# edge correction |
302 | 303 |
e <- borderEdge(X, maxD) |
303 |
- |
|
304 |
+ |
|
304 | 305 |
# lambda <- as.vector(e%*%t(maxD^2*lam*pi)) |
305 | 306 |
# pred <- predict(fit,lambda^(1/4)) |
306 | 307 |
# pred[lambda < 0.001] = (lambda - 4*lambda^2)[lambda < 0.001] |
307 | 308 |
# pred[lambda > mL^(1/4)] = 0.25 |
308 | 309 |
# V <- e%*%t(maxD^2*lam*pi) |
309 | 310 |
# V[] <- pred |
310 |
- |
|
311 |
- |
|
311 |
+ |
|
312 |
+ |
|
312 | 313 |
lambda <- as.vector(maxD^2 * lam * pi) |
313 | 314 |
names(lambda) <- names(lam) |
314 | 315 |
LE <- (e) %*% t(lambda) |
315 | 316 |
mat <- apply(dt, 2, function(x) x) |
316 | 317 |
mat <- ((mat) - (LE)) |
317 | 318 |
mat <- mat / sqrt(LE) |
318 |
- |
|
319 |
+ |
|
319 | 320 |
# # plot(apply(mat,2,sd)) |
320 | 321 |
# # plot(apply(mat,2,mean)) |
321 | 322 |
colnames(mat) <- paste(maxD, colnames(mat), sep = "_") |
... | ... |
@@ -380,42 +381,42 @@ inhomLocalK <- |
380 | 381 |
window = ow, |
381 | 382 |
marks = data$cellType |
382 | 383 |
) |
383 |
- |
|
384 |
+ |
|
384 | 385 |
if (is.null(Rs)) { |
385 | 386 |
Rs <- c(20, 50, 100, 200) |
386 | 387 |
} |
387 | 388 |
if (is.null(sigma)) { |
388 | 389 |
sigma <- 100000 |
389 | 390 |
} |
390 |
- |
|
391 |
+ |
|
391 | 392 |
maxR <- min(ow$xrange[2] - ow$xrange[1], ow$yrange[2] - ow$yrange[1]) / 2.01 |
392 | 393 |
Rs <- unique(pmin(c(0, sort(Rs)), maxR)) |
393 |
- |
|
394 |
+ |
|
394 | 395 |
den <- spatstat.explore::density.ppp(X, sigma = sigma) |
395 | 396 |
den <- den / mean(den) |
396 | 397 |
den$v <- pmax(den$v, minLambda) |
397 |
- |
|
398 |
+ |
|
398 | 399 |
p <- spatstat.geom::closepairs(X, max(Rs), what = "ijd") |
399 | 400 |
n <- X$n |
400 | 401 |
p$j <- data$cellID[p$j] |
401 | 402 |
p$i <- data$cellID[p$i] |
402 |
- |
|
403 |
+ |
|
403 | 404 |
cT <- data$cellType |
404 | 405 |
names(cT) <- data$cellID |
405 |
- |
|
406 |
+ |
|
406 | 407 |
p$d <- cut(p$d, Rs, labels = Rs[-1], include.lowest = TRUE) |
407 |
- |
|
408 |
+ |
|
408 | 409 |
# inhom density |
409 | 410 |
np <- spatstat.geom::nearest.valid.pixel(X$x, X$y, den) |
410 | 411 |
w <- den$v[cbind(np$row, np$col)] |
411 | 412 |
names(w) <- data$cellID |
412 | 413 |
p$wt <- 1 / w[p$j] * mean(w) |
413 | 414 |
rm(np) |
414 |
- |
|
415 |
+ |
|
415 | 416 |
lam <- table(data$cellType) / spatstat.geom::area(X) |
416 |
- |
|
417 |
- |
|
418 |
- |
|
417 |
+ |
|
418 |
+ |
|
419 |
+ |
|
419 | 420 |
p$cellTypeJ <- cT[p$j] |
420 | 421 |
p$cellTypeI <- cT[p$i] |
421 | 422 |
p$i <- factor(p$i, levels = data$cellID) |
... | ... |
@@ -424,21 +425,21 @@ inhomLocalK <- |
424 | 425 |
colnames(edge) <- Rs[-1] |
425 | 426 |
edge$i <- data$cellID |
426 | 427 |
edge <- tidyr::pivot_longer(edge, -i, names_to = "d") |
427 |
- |
|
428 |
+ |
|
428 | 429 |
p <- dplyr::left_join(as.data.frame(p), edge, c("i", "d")) |
429 | 430 |
p$d <- factor(p$d, levels = Rs[-1]) |
430 |
- |
|
431 |
- |
|
432 |
- |
|
431 |
+ |
|
432 |
+ |
|
433 |
+ |
|
433 | 434 |
p <- as.data.frame(p) |
434 |
- |
|
435 |
+ |
|
435 | 436 |
if (lisaFunc == "K") { |
436 | 437 |
r <- getK(p, lam) |
437 | 438 |
} |
438 | 439 |
if (lisaFunc == "L") { |
439 | 440 |
r <- getL(p, lam) |
440 | 441 |
} |
441 |
- |
|
442 |
+ |
|
442 | 443 |
as.matrix(r[data$cellID, ]) |
443 | 444 |
} |
444 | 445 |
|
... | ... |
@@ -458,11 +459,11 @@ getK <- |
458 | 459 |
r$wt <- (r$wt - E) / sqrt(E) |
459 | 460 |
r <- r[, value := NULL] |
460 | 461 |
r <- data.table::dcast(r, i ~ d + cellTypeJ, value.var = "wt") |
461 |
- |
|
462 |
+ |
|
462 | 463 |
r <- as.data.frame(r) |
463 | 464 |
rownames(r) <- r$i |
464 | 465 |
r <- r[, -1] |
465 |
- |
|
466 |
+ |
|
466 | 467 |
r |
467 | 468 |
} |
468 | 469 |
|
... | ... |
@@ -482,11 +483,11 @@ getL <- |
482 | 483 |
r$wt <- sqrt(r$wt) - sqrt(E) |
483 | 484 |
r <- r[, value := NULL] |
484 | 485 |
r <- data.table::dcast(r, i ~ d + cellTypeJ, value.var = "wt") |
485 |
- |
|
486 |
+ |
|
486 | 487 |
r <- as.data.frame(r) |
487 | 488 |
rownames(r) <- r$i |
488 | 489 |
r <- r[, -1] |
489 |
- |
|
490 |
+ |
|
490 | 491 |
r |
491 | 492 |
} |
492 | 493 |
|
... | ... |
@@ -535,16 +536,16 @@ regionMap <- function(cells, type = "bubble", cellType = "cellType", region = "r |
535 | 536 |
if (is.data.frame(cells)) { |
536 | 537 |
df <- cells[, c(cellType, region)] |
537 | 538 |
} |
538 |
- |
|
539 |
+ |
|
539 | 540 |
if (is(cells, "SingleCellExperiment") | is(cells, "SpatialExperiment")) { |
540 | 541 |
df <- as.data.frame(SummarizedExperiment::colData(cells))[, c(cellType, region)] |
541 | 542 |
} |
542 |
- |
|
543 |
+ |
|
543 | 544 |
tab <- table(df[, cellType], df[, region]) |
544 | 545 |
tab <- tab / rowSums(tab) %*% t(colSums(tab)) * sum(tab) |
545 |
- |
|
546 |
+ |
|
546 | 547 |
ph <- pheatmap::pheatmap(pmax(pmin(tab, limit[2]), limit[1]), cluster_cols = FALSE, silent = TRUE, ...) |
547 |
- |
|
548 |
+ |
|
548 | 549 |
if (type == "bubble") { |
549 | 550 |
p1 <- tab |> |
550 | 551 |
as.data.frame() |> |
... | ... |
@@ -554,9 +555,9 @@ regionMap <- function(cells, type = "bubble", cellType = "cellType", region = "r |
554 | 555 |
ggplot2::scale_colour_gradient2(low = "#4575B4", mid = "grey90", high = "#D73027", midpoint = 1, guide = "legend") + |
555 | 556 |
ggplot2::theme_minimal() + |
556 | 557 |
ggplot2::labs(x = "Region", y = "Cell-type", colour = "Relative\nFrequency", size = "Relative\nFrequency") |
557 |
- |
|
558 |
+ |
|
558 | 559 |
return(p1) |
559 | 560 |
} |
560 |
- |
|
561 |
+ |
|
561 | 562 |
pheatmap::pheatmap(pmax(pmin(tab, limit[2]), limit[1]), cluster_cols = FALSE, ...) |
562 | 563 |
} |
... | ... |
@@ -14,7 +14,6 @@ lisa( |
14 | 14 |
sigma = NULL, |
15 | 15 |
lisaFunc = "K", |
16 | 16 |
minLambda = 0.05, |
17 |
- fast = TRUE, |
|
18 | 17 |
spatialCoords = c("x", "y"), |
19 | 18 |
cellType = "cellType", |
20 | 19 |
imageID = "imageID" |
... | ... |
@@ -42,9 +41,6 @@ the cellType.} |
42 | 41 |
|
43 | 42 |
\item{minLambda}{Minimum value for density for scaling when fitting inhomogeneous L-curves.} |
44 | 43 |
|
45 |
-\item{fast}{A logical describing whether to use a fast approximation of the |
|
46 |
-inhomogeneous local L-curves.} |
|
47 |
- |
|
48 | 44 |
\item{spatialCoords}{The columns which contain the x and y spatial coordinates.} |
49 | 45 |
|
50 | 46 |
\item{cellType}{The column which contains the cell types.} |
... | ... |
@@ -18,8 +18,7 @@ lisaClust( |
18 | 18 |
whichParallel = "imageID", |
19 | 19 |
sigma = NULL, |
20 | 20 |
lisaFunc = "K", |
21 |
- minLambda = 0.05, |
|
22 |
- fast = TRUE |
|
21 |
+ minLambda = 0.05 |
|
23 | 22 |
) |
24 | 23 |
} |
25 | 24 |
\arguments{ |
... | ... |
@@ -53,9 +52,6 @@ the cellType.} |
53 | 52 |
\item{lisaFunc}{Either "K" or "L" curve.} |
54 | 53 |
|
55 | 54 |
\item{minLambda}{Minimum value for density for scaling when fitting inhomogeneous L-curves.} |
56 |
- |
|
57 |
-\item{fast}{A logical describing whether to use a fast approximation of the |
|
58 |
-inhomogeneous local L-curves.} |
|
59 | 55 |
} |
60 | 56 |
\value{ |
61 | 57 |
A matrix of LISA curves |