GroupStats: GroupStats

View source: R/2_buildSOM.R

GroupStatsR Documentation

GroupStats

Description

Calculate statistics between 2 groups based on the GetFeatures output

Usage

GroupStats(features, groups)

Arguments

features

Feature matrix as generated by GetFeatures, e.g. a percentages matrix

groups

Named list with file or patient IDs per group (should match with the rownames of the matrix).

Value

Matrix with the medians per group, the p-values (the raw, Benjamini Hochberg corrected one and the -log10) that resulted from a Wilcox test and the fold and log10 fold changes between the medians of the 2 groups

Examples

# Build FlowSom result
 fileName <- system.file("extdata", "68983.fcs", package = "FlowSOM")
 ff <- flowCore::read.FCS(fileName)
 ff <- flowCore::compensate(ff, flowCore::keyword(ff)[["SPILL"]])
 ff <- flowCore::transform(ff,
         flowCore::transformList(colnames(flowCore::keyword(ff)[["SPILL"]]),
                                flowCore::logicleTransform()))
 flowSOM.res <- FlowSOM(ff, scale = TRUE, colsToUse = c(9, 12, 14:18),
                          nClus = 10)
   
# Create new data
# To illustrate the output, we here generate new FCS files (with more 
# cells in metaclusters 1 and 9).
# In practice you would not generate any new file but use your different
# files from your different groups
 flowCore::write.FCS(ff[sample(1:nrow(ff), 1000), ], file = "ff_tmp1.fcs")
 flowCore::write.FCS(ff[sample(1:nrow(ff), 1000), ], file = "ff_tmp2.fcs")
 flowCore::write.FCS(ff[sample(1:nrow(ff), 1000), ], file = "ff_tmp3.fcs")
 ff_tmp <- ff[c(1:1000,
                which(flowSOM.res$map$mapping[, 1] %in% 
                     which(flowSOM.res$metaclustering == 9)),
                which(flowSOM.res$map$mapping[, 1] %in% 
                     which(flowSOM.res$metaclustering == 1))), ]
 flowCore::write.FCS(ff_tmp[sample(1:nrow(ff_tmp), 1000), ],
                     file = "ff_tmp4.fcs")
 flowCore::write.FCS(ff_tmp[sample(1:nrow(ff_tmp), 1000), ], 
                     file = "ff_tmp5.fcs")
 
# Get the count matrix
 percentages <- GetFeatures(fsom = flowSOM.res, 
                            files = c("ff_tmp1.fcs", 
                                      "ff_tmp2.fcs", 
                                      "ff_tmp3.fcs",
                                      "ff_tmp4.fcs", 
                                      "ff_tmp5.fcs"), 
                            type = "percentages")
                       
  
# Perform the statistics
groups <- list("Group 1" = c("ff_tmp1.fcs", "ff_tmp2.fcs", "ff_tmp3.fcs"), 
               "Group 2" = c("ff_tmp4.fcs", "ff_tmp5.fcs"))
MC_stats <- GroupStats(percentages[["metacluster_percentages"]], groups)
C_stats <- GroupStats(percentages[["cluster_percentages"]], groups)

# Process the fold changes vector
fold_changes <- C_stats["fold changes", ]
fold_changes <- factor(ifelse(fold_changes < -3, 
                              "Underrepresented compared to Group 1",
                              ifelse(fold_changes > 3, 
                                     "Overrepresented compared to Group 1",
                                      "--")), 
                        levels = c("--", 
                                   "Underrepresented compared to Group 1",
                                   "Overrepresented compared to Group 1"))
fold_changes[is.na(fold_changes)] <- "--"

# Show in figure
## Fold change
gr_1 <- PlotStars(flowSOM.res, 
                  title = "Group 1", 
                  nodeSizes = C_stats["medians Group 1", ], 
                  list_insteadof_ggarrange = TRUE)
gr_2 <- PlotStars(flowSOM.res, title = "Group 2", 
            nodeSizes = C_stats["medians Group 2", ], 
            backgroundValues = fold_changes,
            backgroundColors = c("white", "red", "blue"), 
            list_insteadof_ggarrange = TRUE)
p <- ggpubr::ggarrange(plotlist = c(list(gr_1$tree), gr_2),
                       heights = c(3, 1))
ggplot2::ggsave("Groups_foldchanges.pdf", p, width = 10)

## p values
p <- PlotVariable(flowSOM.res, title = "Wilcox test group 1 vs. group 2",
variable = C_stats["p values", ])
ggplot2::ggsave("Groups_pvalues.pdf", p)

## volcano plot
p <- ggplot2::ggplot(data.frame("-log10 p values" = c(C_stats[4, ], 
                                                      MC_stats[4, ]), 
                                "log10 fold changes" = c(C_stats[7, ],
                                                         MC_stats[7, ]), 
check.names = FALSE), ggplot2::aes(x = `log10 fold changes`, 
                                   y = `-log10 p values`)) +
ggplot2::xlim(-3, 3) +
ggplot2::ylim(0, 3) +
ggplot2::geom_point() 


saeyslab/FlowSOM documentation built on July 6, 2024, 10:59 a.m.