1简介

Milo是一种用于分析由复制的多条件实验生成的复杂单细胞数据集的工具，它可以检测不同条件之间成分的变化。虽然差异丰度(DA)通常在离散的细胞簇中量化，但Milo在KNN图上使用部分重叠的细胞邻域。从忠实地概括细胞群生物学的图表开始，米洛分析包括3个步骤:

1. 有代表性的社区抽样
2. 测试所有社区条件的差异丰度
3. 使用加权FDR程序对多个假设检验进行解释，该程序解释了邻居的重叠

在这篇短文中，我们将详细说明这些步骤是如何在miloR包中。

2加载数据

对于这个演示，我们将使用一个模拟发展轨迹的合成数据集dyntoy．

data("sim_trajectory"， package = "miloR") ##提取singlecel实验对象traj_sce <- sim_trajectory[['SCE']] ##提取样本元数据用于测试traj_meta <- sim_trajectory[["meta"]] ##添加元数据到colData槽colData(traj_sce) <- DataFrame(traj_meta)

3.预处理

对于DA分析，我们需要构造单细胞的无向KNN图。在PCA中，标准单细胞分析管道通常在距离上进行此操作。我们用归一化和计算主成分嘘．我还运行UMAP以实现可视化。

logcounts(traj_sce) <- log(counts(traj_sce) + 1) traj_sce <- runPCA(traj_sce, ncomponents=30) traj_sce <- runUMAP(traj_sce) plotUMAP(traj_sce)

4创建一个Milo对象

对于图邻域的差异丰度分析，我们首先构造一个米洛对象。这扩展了SingleCellExperiment类来存储KNN图上关于邻居的信息。

traj_milo <- Milo(traj_sce) reducedDim(traj_milo， "UMAP") <- reducedDim(traj_sce， "UMAP") traj_milo

##类:Milo ## dim: 500 500 ##元数据(0):## assays(2): counts logcounts ## rownames(500): G1 G2…G499 G500 ## rowData names(0): ## colnames: NULL ## colData names(5): cell_id group_id条件复制样本## reducedDimNames(2): PCA UMAP ## mainExpName: NULL ## altExpNames(0): ## nhood维度(2):1 1 # nhoodCounts维度(2):1 1 # nhood距离维度(1):0 ## graph names(0): # nhoodIndex names(1): 0 ## nhoodExpression维度(2):1 1 # nhoodReducedDim names(0): ## nhoodGraph names(0): ## nhoodAdjacency维度(2):1 1 1

获取示例H5AD文件。Example_h5ad <- system。文件(“extdata”、“krumsiek11。h5ad", package = "zellkonverter") example_h5ad_sce <- readH5AD(example_h5ad) example_h5ad_milo <- Milo(example_h5ad_sce)

library(Seurat) data("pbmc_small") pbmc_small_sce <- as. singlecel实验(pbmc_small) pbmc_small_milo <- Milo(pbmc_small_sce)

5构造KNN图

我们需要将KNN图添加到Milo对象。存储在图槽,igraph格式。的miloR包包含从存储在reducedDim槽。

traj_milo <- buildGraph(traj_milo, k = 10, d = 30)

构造k:10的kNN图

进展:我们正在完善功能，以添加一个预先计算的KNN图(例如用Seurat或scanpy构造)到图槽使用邻接矩阵。

61.界定有代表性的社区

我们将一个单元的邻域(索引)定义为KNN图中由一条边连接到索引单元的单元组。为了提高效率，我们不测试每个细胞附近的DA，但我们采样作为指标的代表性细胞的子集，使用的KNN采样算法Gut et al. 2015．

为了进行采样，你需要定义几个参数:

• 道具:开始时随机抽样的细胞比例(通常0.1 - 0.2就足够了)
• k:用于KNN优化的k(我们建议使用与KNN图构建相同的k)
• d:用于KNN优化的降维数(我们建议使用与KNN图构建相同的d)
• 精制指示您是想使用采样细化算法，还是只是随机选择单元格。默认和推荐的方法是使用细化。唯一可以考虑使用随机的情况是，如果您已经使用基于图形的校正算法批量校正了数据，例如BBKNN，但DA测试的结果会不理想。

traj_milo <- makeNhoods(traj_milo, prop = 0.1, k = 10, d=30, refined = TRUE)

##检查有效对象

##使用reduced_dim运行精细采样

##警告在.refined_sampling(random_vertices, X_reduced_dims, k): Rownames not ## set on reducedDims -设置为行索引

一旦我们定义了邻域，最好看看邻域有多大(即每个邻域有多少个单元格)。这会影响DA测试的能力。我们可以用plotNhoodSizeHist函数。根据经验，我们发现最好的分布峰值在50到100之间。否则，您可以考虑重新运行makeNhoods增加k和/或道具(这里的分布看起来很可笑，因为它是一个小数据集)。

plotNhoodSizeHist (traj_milo)

7计数小区里的细胞

现在我们必须计算每个样本在每个邻域中有多少个细胞。我们需要使用单元格元数据并指定包含示例信息的列。

traj_milo <- countCells(traj_milo, meta。data = data.frame(colData(traj_milo)), samples="Sample")

##检查meta。数据的有效性

##计数邻居的细胞

这增加了米洛对象一个N \乘以m矩阵，其中n是邻域数和\ \(米)是实验样本的数量。值表示在邻域中计数的每个样本的细胞数。此计数矩阵将用于DA测试。

头(nhoodCounts (traj_milo))

## 6 x 6稀疏矩阵类“dgCMatrix”## B_R1 A_R1 A_R2 B_R2 B_R3 A_R3 ## 1 8 . .14 28 1 ## 2 7 6 5 11 16 3 ## 3 15。1 22 29 1 ## 4 5 3 3 3 3 4 ## 5 2 13 10 3 ## 6 6 3 3 7 5

8差异丰度测试

现在我们都准备好测试不同社区的差异丰度了。我们在广义线性模型(GLM)框架中实现了这个假设检验，特别是在中使用负二项式GLM实现刨边机．

我们首先要考虑实验设计。设计矩阵应使样本与感兴趣的条件相匹配。在这种情况下条件是我们要检验的协变量。

traj_design <- data.frame(colData(traj_milo))[，c("Sample"， "Condition")] traj_design <- distinct(traj_design) rownames(traj_design) <- traj_design$Sample ##重新排列rownames以匹配nhoodCounts(milo)的列traj_design <- traj_design[colnames(nhoodCounts(traj_milo))，， drop=FALSE] traj_design

样本条件## B_R1 B_R1 B ## A_R1 A_R1 A ## A_R2 A_R2 A ## B_R2 B_R2 B_R3 B_R3 B ## A_R3 A_R3 A

Milo使用了空间FDR校正的改编cydar，这说明了社区之间的重叠。具体而言，每个假设检验p值均以第k个最近邻距离的倒数加权。要使用这个统计数据，我们首先需要存储Milo对象中最近邻居之间的距离。

traj_milo <- calcNhoodDistance(traj_milo, d=30)

现在我们可以进行测试，明确地定义我们的实验设计。

rownames(traj_design) <- traj_design$Sample da_results <- testNhoods(traj_milo, design = ~条件，设计。Df = traj_design)

##使用TMM标准化

使用k-距离加权执行空间FDR校正

这将计算每个邻域的折叠变化和校正的p值，这表明条件之间是否存在显著的差异丰度。

da_results %>% arrange(- SpatialFDR) %>% head()

## logFC logCPM F PValue FDR Nhood SpatialFDR ## 19 0.01880599 15.41161 0.001386078 0.9731932 0.9731932 19 0.9731932 ## 9 -0.25556996 15.71117 0.266199971 0.6069875 0.6312670 9 0.6311887 ## 11 -0.32274234 15.44814 0.471935719 0.5503505 0.5962130 11 0.5957746 ## 4 - 0.532274234 15.17063 0.605934505 0.4380901 0.4952323 4 0.4954764 ## 15 0.39376174 15.95905 1.174088971 0.3807867 0.4500207 15 0.4486690 ## 6 -0.54039326 15.51855 1.356294066 0.3688407 6 0.3672821

9可视化显示DA的社区

为了使DA结果可视化，我们构建了一个可以叠加在单细胞嵌入上的邻域抽象图。

traj_milo <- buildNhoodGraph(traj_milo) plotUMAP(traj_milo) + plotNhoodGraphDA(traj_milo, da_results, alpha=0.05) + plot_layout(guides="collect")

sessionInfo ()

## R版本4.2.0 RC (2022-04-19 r82224) ##平台:x86_64-pc-linux-gnu(64位)##运行在Ubuntu 20.04.4 LTS ## ##矩阵产品:默认## BLAS: /home/biocbuild/bbs-3.15-bioc/R/lib/libRblas。/home/biocbuild/bbs-3.15-bioc/R/lib/libRlapack。所以## ## locale: ## [1] LC_CTYPE=en_US。UTF-8 LC_NUMERIC= c# # [3] LC_TIME=en_GB LC_COLLATE= c# # [5] LC_MONETARY=en_US。utf - 8 LC_MESSAGES = en_US。UTF-8 ## [7] LC_PAPER=en_US。UTF-8 LC_NAME= c# # [9] LC_ADDRESS=C lc_phone = c# # [11] LC_MEASUREMENT=en_US。UTF-8 LC_IDENTIFICATION=C ## ##附加的基本包:## [1]stats4 stats graphics grDevices utils datasets methods ##[8]基础## ##其他附加包:# # # # [1] MouseThymusAgeing_1.3.0 patchwork_1.1.1 [3] dplyr_1.0.8 scran_1.24.0 # # [5] scater_1.24.0 ggplot2_3.3.5 # # [7] scuttle_1.6.0 SingleCellExperiment_1.18.0 # # [9] SummarizedExperiment_1.26.0 Biobase_2.56.0 # # [11] GenomicRanges_1.48.0 GenomeInfoDb_1.32.0 # # [13] IRanges_2.30.0 S4Vectors_0.34.0 # # [15] BiocGenerics_0.42.0 MatrixGenerics_1.8.0 # # [17] matrixStats_0.62.0 miloR_1.4.0 # # [19] edgeR_3.38.0 limma_3.52.0 # # [21] BiocStyle_2.24.0 # # # #通过加载一个名称空间(而不是附加):# # # # [1] AnnotationHub_3.4.0 BiocFileCache_2.4.0 [3] igraph_1.3.1 splines_4.2.0 # # [5] BiocParallel_1.30.0 digest_0.6.29 # # [7] htmltools_0.5.2 magick_2.7.3 # # [9] viridis_0.6.2 fansi_1.0.3 # # [11] magrittr_2.0.3 memoise_2.0.1 # # [13] ScaledMatrix_1.4.0 cluster_2.1.3 # # [15] Biostrings_2.64.0 graphlayouts_0.8.0 # # [17] colorspace_2.0-3 blob_1.2.3 # # [19] rappdirs_0.3.3 ggrepel_0.9.1 # # [21] xfun_0.30 crayon_1.5.1 # # [23] rcurl_1.98 - 1.6 jsonlite_1.8.0 # # [25] glue_1.6.2 polyclip_1.10-0 # # [27][39] rsvd_1.0.5 metapod_1.4.0 ## [41] httr_1.4.2 FNN_1.1.3 ## [43] RColorBrewer_1.1-3 ellipsis_0.3.2 ## [45] pkgconfig_2.0.3 farver_2.1.0 ## [47] uwot_0.1.11 sass_0.4.1 ## [51] utf8_1. 1.1 locfit_1. 1.5-9.5 ## [51] utf8_1.2.2 labeling_0.0.2 ## [55] ## [53] tidyselect_1.1.2 rcpp_1 .2.0 ## # [55][63] RSQLite_2.2.12 ExperimentHub_2.4.0 ## [65] evaluate_0.15 string_1 .4.0 ## [67] fastmap_1.1.0 yaml_2.3.5 ## [69] knitr_1.38 bit64_4.0.5 ## [71] tidygraph_1.2.1 purrr_0.3.4 ## [73] KEGGREST_1.36.0 ggraph_2.0.5 ## [75] sparseMatrixStats_1.8.0 mime_0.12 ## [77] compiler_4.2.0 png_0.1-7 ## [81] beeswarm_0.4.0 filelock_1.0.2 ##[83] tibble_1 .1.6 statmod_1.4.36 ## [85] tweenr_1. 0.6 bslib_3.1 ## [87] stringi_1.7.6 highr_0.9 ## [89] RSpectra_0.16-1 lattice_0.20-45 ## [91] bluster_1.6.0 Matrix_1.4-1 ## [93] vctrs_0.4.1 pillar_1.7.0 ## [95] lifecycle_1.0.1 BiocManager_1.30.17 ## [97] jquerylib_0.1.4 BiocNeighbors_1.14.0 ## [99] cowplot_1.1.1 bitops_1.0-7 ## [101] irlba_2.3.5 httpuv_1.6.5 ## [103] R6_2.5.1 bookdown_0.26 ## [105] promises_1.2.0.1 gridExtra_2.3 ## [107] vipor_0.4.5 MASS_7.3-57 ## [109] gtools_3.9.2 assertthat_0.2.1 ## [111] withr_2.5.0 GenomeInfoDbData_1.2.8 ## [113] parallel_4.2.0 grid_4.2.0 ## [115] beachmat_2.12.0 tidyr_1.2.0 ## [117] rmarkdown_2.14 DelayedMatrixStats_1.18.0 ## [119] ggforce_0.3.3 shiny_1.7.1 ## [121] ggbeeswarm_0.6.0