Skip to content

Commit 448255b

Browse files
committed
Merge branch 'alternative_split' of https://github.com/openproblems-bio/task_cyto_batch_integration into alternative_split
2 parents 016ea73 + 9af2c4e commit 448255b

8 files changed

Lines changed: 513 additions & 5 deletions

File tree

CHANGELOG.md

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -56,7 +56,9 @@
5656

5757
* Added CytoNorm correction to a goal batch (PR #92).
5858
* Added cyCombine correction to a reference batch (PR #90).
59-
* Added `metrics/bras`
59+
* Added `metrics/bras` (PR #91).
60+
61+
* Added Seurat rPCA (PR #95).
6062

6163

6264
## MAJOR CHANGES
Lines changed: 104 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,104 @@
1+
# The API specifies which type of component this is.
2+
# It contains specifications for:
3+
# - The input/output files
4+
# - Common parameters
5+
# - A unit test
6+
__merge__: ../../api/comp_method.yaml
7+
8+
# A unique identifier for your component (required).
9+
# Can contain only lowercase letters or underscores.
10+
name: rpca_to_goal
11+
# A relatively short label, used when rendering visualisations (required)
12+
label: Seurat RPCA (to-goal)
13+
14+
# A one sentence summary of how this method works (required). Used when
15+
# rendering summary tables.
16+
summary: "Batch integrate data to a goal batch using mutual nearest neighbors identified via Seurat reciprocal PCA."
17+
18+
# A multi-line description of how this component works (required). Used
19+
# when rendering reference documentation.
20+
description: |
21+
Seurat RPCA performs batch integration by projecting each query dataset into the PCA space
22+
of a goal batch, and identifying anchors using reciprocal PCA (RPCA).
23+
RPCA identifies mutual nearest neighbors between the goal batch and the remaining batches
24+
in their shared low-dimensional space (PCA), which are used to align
25+
and integrate the batches into the goal batch.
26+
27+
We ran Seurat RPCA implemented in Seurat v4.4.0, available on their github.
28+
This is because subsequent version (>= v5) does not support getting corrected
29+
count matrix, only corrected PC space.
30+
See: https://github.com/satijalab/seurat/issues/8551.
31+
32+
We varied the number of PCs and nearest neighbours considered when running RPCA.
33+
34+
references:
35+
doi:
36+
- 10.1016/j.cell.2021.04.048
37+
bibtex:
38+
- |
39+
@article{hao2021integrated,
40+
title={Integrated analysis of multimodal single-cell data},
41+
author={Hao, Yuhan and Hao, Stephanie and Andersen-Nissen, Erica and Mauck, William M and Zheng, Shiwei and Butler, Andrew and Lee, Madison J and Wilk, Aaron J and Darby, Charlotte and Zager, Michael and others},
42+
journal={Cell},
43+
volume={184},
44+
number={13},
45+
pages={3573--3587},
46+
year={2021},
47+
publisher={Elsevier}
48+
}
49+
50+
links:
51+
# URL to the documentation for this method (required).
52+
documentation: https://satijalab.org/seurat/articles/integration_rpca.html
53+
# URL to the code repository for this method (required).
54+
repository: https://github.com/satijalab/seurat
55+
56+
argument_groups:
57+
- name: Parameters
58+
arguments:
59+
- type: integer
60+
name: --npcs
61+
info:
62+
optimize:
63+
type: linear
64+
lower: 10
65+
upper: 20
66+
default: 10
67+
description: The number of principal component (PC) dimensions to use when computing anchors between batches.
68+
- type: integer
69+
name: --n_neighbours
70+
info:
71+
optimize:
72+
type: linear
73+
lower: 5
74+
upper: 50
75+
default: 5
76+
description: The number of mutual nearest neighbors that are used when identifying anchors between batches
77+
78+
# Resources required to run the component
79+
resources:
80+
# The script of your component (required)
81+
- type: r_script
82+
path: script.R
83+
84+
engines:
85+
# Specifications for the Docker image for this component.
86+
- type: docker
87+
image: openproblems/base_r:1
88+
# Add custom dependencies here (optional). For more information, see
89+
# https://viash.io/reference/config/engines/docker/#setup .
90+
setup:
91+
- type: r
92+
packages: [ anndata ]
93+
url: [
94+
https://cran.r-project.org/src/contrib/Archive/Seurat/Seurat_4.4.0.tar.gz,
95+
https://cran.r-project.org/src/contrib/Archive/SeuratObject/SeuratObject_4.1.4.tar.gz
96+
]
97+
98+
runners:
99+
# This platform allows running the component natively
100+
- type: executable
101+
# Allows turning the component into a Nextflow module / pipeline.
102+
- type: nextflow
103+
directives:
104+
label: [midtime,midmem,midcpu]

src/methods/rpca_to_goal/script.R

Lines changed: 147 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,147 @@
1+
requireNamespace("anndata", quietly = TRUE)
2+
requireNamespace("Seurat", quietly = TRUE)
3+
4+
## VIASH START
5+
par <- list(
6+
input = "resources_test/task_cyto_batch_integration/mouse_spleen_flow_cytometry_subset/unintegrated_censored.h5ad",
7+
output = "resources_test/task_cyto_batch_integration/mouse_spleen_flow_cytometry_subset/output.h5ad",
8+
npcs = 10,
9+
n_neighbours = 50
10+
)
11+
meta <- list(
12+
name = "rpca_to_goal"
13+
)
14+
## VIASH END
15+
16+
cat("Reading input files\n")
17+
input_adata <- anndata::read_h5ad(par[["input"]])
18+
19+
cat("Preparing input Anndata\n")
20+
input_adata$obs$batch <- as.factor(input_adata$obs$batch)
21+
22+
adata_to_correct <- input_adata[, input_adata$var$to_correct]
23+
markers_to_correct <- input_adata$var_names[input_adata$var$to_correct]
24+
25+
cat("Creating Seurat object and preprocess\n")
26+
27+
# create one seurat object per batch
28+
batches <- unique(input_adata$obs$batch)
29+
30+
seurat_objs <- lapply(batches, function(batch) {
31+
32+
cat(paste("Processing batch", batch))
33+
34+
adata_batch <- input_adata[
35+
input_adata$obs$batch == batch,
36+
input_adata$var$to_correct
37+
]
38+
# batch <- batches[1]
39+
mat <- Matrix::as.matrix(adata_batch$layers["preprocessed"])
40+
41+
# have to transpose so cells are columns..
42+
mat <- Matrix::t(mat)
43+
44+
# convert to sparse matrix
45+
mat <- Matrix::Matrix(mat, sparse = TRUE)
46+
47+
seurat_obj <- Seurat::CreateSeuratObject(
48+
counts = mat,
49+
data = mat,
50+
assay = "cyto",
51+
meta.data = adata_batch$obs
52+
)
53+
54+
# save RAM
55+
rm(mat)
56+
57+
# scale all features/markers
58+
seurat_obj <- Seurat::ScaleData(
59+
object = seurat_obj,
60+
features = markers_to_correct,
61+
assay = "cyto",
62+
verbose = FALSE
63+
)
64+
65+
# run pca. mandatory
66+
# if num pcs is more than number of markers, it'll be capped at
67+
# the number of markers
68+
# not using approximate pca as we don't have many markers
69+
seurat_obj <- Seurat::RunPCA(
70+
object = seurat_obj,
71+
features = markers_to_correct,
72+
assay = "cyto",
73+
npcs = par[["npcs"]],
74+
approx = FALSE,
75+
verbose = FALSE
76+
)
77+
78+
return(seurat_obj)
79+
})
80+
81+
names(seurat_objs) <- batches
82+
83+
cat("Finding anchors\n")
84+
85+
# get how many PCs we have calculated
86+
# if the number of PCs is more than how many markers
87+
# seurat set that to how many markers.
88+
# hence we can't just use par[["npcs"]] below.
89+
npcs_computed <- dim(seurat_objs[[1]][["pca"]])[2]
90+
91+
anchors <- Seurat::FindIntegrationAnchors(
92+
object.list = seurat_objs,
93+
anchor.features = markers_to_correct,
94+
dims = seq(npcs_computed),
95+
k.anchor = par[["n_neighbours"]],
96+
reduction = "rpca",
97+
verbose = FALSE,
98+
reference = which(names(seurat_objs) == "1")
99+
)
100+
101+
cat("Batch correct\n")
102+
103+
# Warning will say Layer counts isn't present in the assay object; returning NULL
104+
# Even though the original assay has counts layer.
105+
# Not sure why. But the object has data layer.
106+
batch_corrected_seurat_obj <- Seurat::IntegrateData(
107+
anchorset = anchors,
108+
features = markers_to_correct,
109+
features.to.integrate = markers_to_correct,
110+
dims = seq(npcs_computed),
111+
verbose = FALSE
112+
)
113+
# just to be sure!
114+
Seurat::DefaultAssay(batch_corrected_seurat_obj) <- "integrated"
115+
116+
cat("Creating output AnnData\n")
117+
118+
batch_corrected_mat <- Matrix::t(
119+
Matrix::as.matrix(Seurat::GetAssayData(batch_corrected_seurat_obj))
120+
)
121+
# cbind corrected matrix to matrix containing markers not corrected
122+
batch_corrected_mat <- cbind(
123+
batch_corrected_mat,
124+
input_adata[, !input_adata$var$to_correct]$layers[["preprocessed"]]
125+
)
126+
127+
# make sure the row and column orders are matching
128+
# between input adata and the batch corrected matrix
129+
batch_corrected_mat <- batch_corrected_mat[
130+
input_adata$obs_names, input_adata$var_names
131+
]
132+
133+
cat("Write output AnnData to file\n")
134+
output <- anndata::AnnData(
135+
obs = input_adata$obs[, integer(0)],
136+
var = input_adata$var[colnames(batch_corrected_mat), integer(0)],
137+
layers = list(integrated = batch_corrected_mat),
138+
uns = list(
139+
dataset_id = input_adata$uns$dataset_id,
140+
method_id = meta$name,
141+
parameters = list(
142+
"npcs" = par[["npcs"]],
143+
"n_neighbours" = par[["n_neighbours"]]
144+
)
145+
)
146+
)
147+
output$write_h5ad(par[["output"]], compression = "gzip")
Lines changed: 106 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,106 @@
1+
# The API specifies which type of component this is.
2+
# It contains specifications for:
3+
# - The input/output files
4+
# - Common parameters
5+
# - A unit test
6+
__merge__: ../../api/comp_method.yaml
7+
8+
# A unique identifier for your component (required).
9+
# Can contain only lowercase letters or underscores.
10+
name: rpca_to_mid
11+
# A relatively short label, used when rendering visualisations (required)
12+
label: Seurat RPCA (to-middle)
13+
14+
# A one sentence summary of how this method works (required). Used when
15+
# rendering summary tables.
16+
summary: "Batch integrate data to a midpoint using mutual nearest neighbors identified via Seurat reciprocal PCA."
17+
18+
# A multi-line description of how this component works (required). Used
19+
# when rendering reference documentation.
20+
description: |
21+
Seurat RPCA performs batch integration by identifying mutual nearest neighbors (anchors)
22+
between all batches using reciprocal PCA (RPCA).
23+
Instead of merging batches into a single PCA space, each batch is projected into the
24+
PCA space of the others, and anchors are found where cells are mutual nearest neighbors
25+
across these projections.
26+
These anchors are then used to compute batch correction vectors,
27+
aligning the batches into a shared corrected space.
28+
29+
We ran Seurat RPCA implemented in Seurat v4.4.0, available on their github.
30+
This is because subsequent version (>= v5) does not support getting corrected
31+
count matrix, only corrected PC space.
32+
See: https://github.com/satijalab/seurat/issues/8551
33+
34+
We varied the number of PCs and nearest neighbours considered when running RPCA.
35+
36+
references:
37+
doi:
38+
- 10.1016/j.cell.2021.04.048
39+
bibtex:
40+
- |
41+
@article{hao2021integrated,
42+
title={Integrated analysis of multimodal single-cell data},
43+
author={Hao, Yuhan and Hao, Stephanie and Andersen-Nissen, Erica and Mauck, William M and Zheng, Shiwei and Butler, Andrew and Lee, Madison J and Wilk, Aaron J and Darby, Charlotte and Zager, Michael and others},
44+
journal={Cell},
45+
volume={184},
46+
number={13},
47+
pages={3573--3587},
48+
year={2021},
49+
publisher={Elsevier}
50+
}
51+
52+
links:
53+
# URL to the documentation for this method (required).
54+
documentation: https://satijalab.org/seurat/articles/integration_rpca.html
55+
# URL to the code repository for this method (required).
56+
repository: https://github.com/satijalab/seurat
57+
58+
argument_groups:
59+
- name: Parameters
60+
arguments:
61+
- type: integer
62+
name: --npcs
63+
info:
64+
optimize:
65+
type: linear
66+
lower: 10
67+
upper: 20
68+
default: 10
69+
description: The number of principal component (PC) dimensions to use when computing anchors between batches.
70+
- type: integer
71+
name: --n_neighbours
72+
info:
73+
optimize:
74+
type: linear
75+
lower: 5
76+
upper: 50
77+
default: 5
78+
description: The number of mutual nearest neighbors that are used when identifying anchors between batches
79+
80+
# Resources required to run the component
81+
resources:
82+
# The script of your component (required)
83+
- type: r_script
84+
path: script.R
85+
86+
engines:
87+
# Specifications for the Docker image for this component.
88+
- type: docker
89+
image: openproblems/base_r:1
90+
# Add custom dependencies here (optional). For more information, see
91+
# https://viash.io/reference/config/engines/docker/#setup .
92+
setup:
93+
- type: r
94+
packages: [ anndata ]
95+
url: [
96+
https://cran.r-project.org/src/contrib/Archive/Seurat/Seurat_4.4.0.tar.gz,
97+
https://cran.r-project.org/src/contrib/Archive/SeuratObject/SeuratObject_4.1.4.tar.gz
98+
]
99+
100+
runners:
101+
# This platform allows running the component natively
102+
- type: executable
103+
# Allows turning the component into a Nextflow module / pipeline.
104+
- type: nextflow
105+
directives:
106+
label: [midtime,midmem,midcpu]

0 commit comments

Comments
 (0)