### simu_new example workflow Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/simu_new.html Complete example of using simu_new with example_sce data. It constructs data, fits marginal and copula models, extracts parameters, and simulates new counts. The output includes progress messages from fit_copula and simu_new. ```R data(example_sce) my_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "1" ) my_marginal <- fit_marginal( data = my_data, mu_formula = "s(pseudotime, bs = 'cr', k = 10)", sigma_formula = "1", family_use = "nb", n_cores = 1, usebam = FALSE ) my_copula <- fit_copula( sce = example_sce, assay_use = "counts", marginal_list = my_marginal, family_use = c(rep("nb", 5), rep("zip", 5)), copula = "vine", n_cores = 1, input_data = my_data$dat ) #> Convert Residuals to Uniform #> Converting End #> Copula group 1 starts #> Vine Copula Estimation Starts #> Time difference of 0.1432807 secs #> Vine Copula Estimation Ends my_para <- extract_para( sce = example_sce, marginal_list = my_marginal, n_cores = 1, family_use = c(rep("nb", 5), rep("zip", 5)), new_covariate = my_data$new_covariate, data = my_data$dat ) my_newcount <- simu_new( sce = example_sce, mean_mat = my_para$mean_mat, sigma_mat = my_para$sigma_mat, zero_mat = my_para$zero_mat, quantile_mat = NULL, copula_list = my_copula$copula_list, n_cores = 1, family_use = c(rep("nb", 5), rep("zip", 5)), input_data = my_data$dat, new_covariate = my_data$new_covariate, important_feature = my_copula$important_feature, filtered_gene = my_data$filtered_gene ) #> Use Copula to sample a multivariate quantile matrix #> Sample Copula group 1 starts ``` -------------------------------- ### Example for ba() printing "No example" Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/ba.html Shows the only example provided for the ba() function, which simply prints the string "No example". The output is displayed as a comment with the result. ```r print("No example") #> [1] "No example" ``` -------------------------------- ### Construct input data with example_sce Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/construct_data.html Example usage of construct_data with example_sce, specifying assay_use, celltype, pseudotime, and corr_by. The example does not set ncell, so newCovariate will be NULL. ```r data(example_sce) my_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "1" ) ``` -------------------------------- ### Run scDesign3 simulation with scdesign3() Source: https://songdongyuan1994.github.io/scDesign3/docs/index.html Quick start example of running the scDesign3 simulator. The scdesign3() function takes a SingleCellExperiment object with cell covariates stored in colData and returns a simulated dataset. This example uses a negative binomial family, pseudotime as a smooth term, and a Gaussian copula. ```r example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime, k = 10, bs = 'cr')", sigma_formula = "s(pseudotime, k = 5, bs = 'cr')", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", fastmvn = FALSE, DT = TRUE, pseudo_obs = FALSE, family_set = c("gauss", "indep"), important_feature = "all", nonnegative = TRUE, return_model = FALSE, nonzerovar = FALSE, parallelization = "mcmapply", BPPARAM = NULL, trace = FALSE ) ``` -------------------------------- ### Install scDesign3 from GitHub with devtools Source: https://songdongyuan1994.github.io/scDesign3/docs/index.html Installs the development version of scDesign3 from GitHub using devtools. Requires the devtools package, which is installed if missing. ```r if (!require("devtools", quietly = TRUE)) install.packages("devtools") devtools::install_github("SONGDONGYUAN1994/scDesign3") ``` -------------------------------- ### Example: perform_lrt with constructed data and fitted marginal models Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/perform_lrt.html Complete example from the documentation. It constructs two datasets with different correlation structures, fits marginal models with and without a pseudotime smooth term, extracts the fits, and performs the likelihood ratio test. The output includes a progress bar from the second construct_data call. ```r data(example_sce) my_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "cell_type" ) my_data2 <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "pseudotime", ncell = 10000 ) #> fitting ... #> | | | 0% | |====== | 9% | |============= | 18% | |=================== | 27% | |========================= | 36% | |================================ | 45% | |====================================== | 55% | |============================================= | 64% | |=================================================== | 73% | |========================================================= | 82% | |================================================================ | 91% | |======================================================================| 100% my_marginal1 <- fit_marginal( data = my_data, mu_formula = "1", sigma_formula = "1", family_use = "nb", n_cores = 1, usebam = FALSE ) my_marginal2 <- fit_marginal( data = my_data, mu_formula = "s(pseudotime, bs = 'cr', k = 10)", sigma_formula = "1", family_use = "nb", n_cores = 1, usebam = FALSE ) my_fit1 <- lapply(my_marginal1, function(x)x$fit) my_fit2 <- lapply(my_marginal2, function(x)x$fit) my_pvalue <- perform_lrt(my_fit2, my_fit1) ``` -------------------------------- ### Read and inspect example_sce from figshare Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-conditionEffect-vignette.html Downloads the example PBMC dataset from figshare and prints its structure. The dataset is a SingleCellExperiment with 1000 genes and 13999 cells, containing condition information in colData. ```r example_sce <- readRDS((url("https://figshare.com/ndownloader/files/40581977"))) print(example_sce) #> class: SingleCellExperiment #> dim: 1000 13999 #> metadata(0): #> assays(2): counts logcounts #> rownames(1000): FTL HBB ... AC022182.3 SCIN #> rowData names(0): #> colnames(13999): AAACATACATTTCC.1 AAACATACCAGAAA.1 ... TTTGCATGCTAAGC.1 #> TTTGCATGGGACGA.1 #> colData names(8): orig.ident nCount_RNA ... cell_type condition #> reducedDimNames(0): #> mainExpName: RNA #> altExpNames(0): ``` -------------------------------- ### Fit vine copula with example data Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/fit_copula.html Example workflow that constructs data, fits marginal models, and then fits a vine copula. The output shows the expected console messages, including the time difference for the vine copula estimation. ```r data(example_sce) my_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "1" ) my_marginal <- fit_marginal( data = my_data, mu_formula = "s(pseudotime, bs = 'cr', k = 10)", sigma_formula = "1", family_use = "nb", n_cores = 1, usebam = FALSE ) my_copula <- fit_copula( sce = example_sce, assay_use = "counts", marginal_list = my_marginal, family_use = c(rep("nb", 5), rep("zip", 5)), copula = "vine", n_cores = 1, input_data = my_data$dat ) #> Convert Residuals to Uniform #> Converting End #> Copula group 1 starts #> Vine Copula Estimation Starts #> Time difference of 0.1842949 secs #> Vine Copula Estimation Ends ``` -------------------------------- ### extract_para example with construct_data, fit_marginal, and fit_copula Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/extract_para.html Example workflow for extract_para. It constructs data, fits marginal models, fits a vine copula, and then calls extract_para with the fitted marginal list and new covariates. The output shows the copula fitting progress messages. ```r data(example_sce) my_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "1" ) my_marginal <- fit_marginal( data = my_data, mu_formula = "s(pseudotime, bs = 'cr', k = 10)", sigma_formula = "1", family_use = "nb", n_cores = 1, usebam = FALSE ) my_copula <- fit_copula( sce = example_sce, assay_use = "counts", marginal_list = my_marginal, family_use = c(rep("nb", 5), rep("zip", 5)), copula = "vine", n_cores = 1, input_data = my_data$dat ) #> Convert Residuals to Uniform #> Converting End #> Copula group 1 starts #> Vine Copula Estimation Starts #> Time difference of 0.1801219 secs #> Vine Copula Estimation Ends my_para <- extract_para( sce = example_sce, marginal_list = my_marginal, n_cores = 1, family_use = c(rep("nb", 5), rep("zip", 5)), new_covariate = my_data$new_covariate, data = my_data$dat ) ``` -------------------------------- ### Example: fitting marginal models with fit_marginal Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/fit_marginal.html Example usage of fit_marginal with a constructed data object. It uses the 'nb' family, a smooth term for pseudotime, and one core. The example assumes example_sce is available. ```r data(example_sce) my_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, corr_by = "1" ) my_marginal <- fit_marginal( data = my_data, mu_formula = "s(pseudotime, bs = 'cr', k = 10)", sigma_formula = "1", family_use = "nb", n_cores = 1, usebam = FALSE ) ``` -------------------------------- ### Simulate with vine copula using scdesign3 Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-copulaCompare-vignette.html Simulates a new dataset using the vine copula. Uses the same parameters as the Gaussian example but with copula = "vine". The result is stored in Zhengmix4eq_simu_sce_vine. ```R set.seed(123) Zhengmix4eq_simu_sce_vine <- scdesign3(sce = Zhengmix4eq_sce_sub, celltype = 'cell_type', pseudotime = NULL, spatial = NULL, other_covariates = NULL, corr_formula = "cell_type", mu_formula = "cell_type", sigma_formula = "cell_type", n_cores = 2, copula = "vine", assay_use = "normcounts", family_use = "nb", pseudo_obs = TRUE, return_model = TRUE) ``` -------------------------------- ### Inspect colData with batch information Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-batchEffect-vignette.html Displays the first six rows of the colData from the example SingleCellExperiment object, showing the batch column and other metadata such as cell type and experiment. ```r head(colData(example_sce)) #> DataFrame with 6 rows and 13 columns #> phenoid orig.ident nCount_RNA #> #> pbmc2_10X_V2_AAACCTGAGATGGGTC B cell pbmc2 2360 #> pbmc2_10X_V2_AAACCTGAGCGTAATA B cell pbmc2 1888 #> pbmc2_10X_V2_AAACCTGAGCTAGGCA Cytotoxic T cell pbmc2 3456 #> pbmc2_10X_V2_AAACCTGAGGGTCTCC Dendritic cell pbmc2 3802 #> pbmc2_10X_V2_AAACCTGGTCCGAACC CD4+ T cell pbmc2 3826 #> pbmc2_10X_V2_AAACCTGTCGTCCGTT CD4+ T cell pbmc2 2345 #> nFeature_RNA nGene nUMI #> #> pbmc2_10X_V2_AAACCTGAGATGGGTC 1044 1044 2360 #> pbmc2_10X_V2_AAACCTGAGCGTAATA 803 803 1888 #> pbmc2_10X_V2_AAACCTGAGCTAGGCA 1372 1372 3456 #> pbmc2_10X_V2_AAACCTGAGGGTCTCC 1519 1519 3802 #> pbmc2_10X_V2_AAACCTGGTCCGAACC 1451 1451 3826 #> pbmc2_10X_V2_AAACCTGTCGTCCGTT 931 931 2345 #> percent.mito Cluster Experiment #> #> pbmc2_10X_V2_AAACCTGAGATGGGTC 0.0419491525423729 2 pbmc2 #> pbmc2_10X_V2_AAACCTGAGCGTAATA 0.0413135593220339 2 pbmc2 #> pbmc2_10X_V2_AAACCTGAGCTAGGCA 0.0353009259259259 1 pbmc2 #> pbmc2_10X_V2_AAACCTGAGGGTCTCC 0.0420831141504471 6 pbmc2 #> pbmc2_10X_V2_AAACCTGGTCCGAACC 0.0371144798745426 0 pbmc2 #> pbmc2_10X_V2_AAACCTGTCGTCCGTT 0.0652452025586354 0 pbmc2 #> Method ident batch #> #> pbmc2_10X_V2_AAACCTGAGATGGGTC 10x Chromium (v2) pbmc2 10x Chromium (v2) #> pbmc2_10X_V2_AAACCTGAGCGTAATA 10x Chromium (v2) pbmc2 10x Chromium (v2) #> pbmc2_10X_V2_AAACCTGAGCTAGGCA 10x Chromium (v2) pbmc2 10x Chromium (v2) #> pbmc2_10X_V2_AAACCTGAGGGTCTCC 10x Chromium (v2) pbmc2 10x Chromium (v2) #> pbmc2_10X_V2_AAACCTGGTCCGAACC 10x Chromium (v2) pbmc2 10x Chromium (v2) #> pbmc2_10X_V2_AAACCTGTCGTCCGTT 10x Chromium (v2) pbmc2 10x Chromium (v2) #> cell_type #> #> pbmc2_10X_V2_AAACCTGAGATGGGTC B cell #> pbmc2_10X_V2_AAACCTGAGCGTAATA B cell #> pbmc2_10X_V2_AAACCTGAGCTAGGCA Cytotoxic T cell #> pbmc2_10X_V2_AAACCTGAGGGTCTCC Dendritic cell #> pbmc2_10X_V2_AAACCTGGTCCGAACC CD4+ T cell #> pbmc2_10X_V2_AAACCTGTCGTCCGTT CD4+ T cell ``` -------------------------------- ### Load example spatial SingleCellExperiment from figshare Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-DEanalysis-vignette.html Reads a remote RDS file from figshare containing a SingleCellExperiment with 1000 genes and 2696 cells, with counts and logcounts assays and cell_type column. The example output shows the class, dimensions, and assay information. ```R example_sce <- readRDS((url("https://figshare.com/ndownloader/files/40582019"))) print(example_sce) #> class: SingleCellExperiment #> dim: 1000 2696 #> metadata(0): #> assays(2): counts logcounts #> rownames(1000): Calb2 Gng4 ... Fndc5 Gda #> rowData names(0): #> colnames(2696): AAACAAGTATCTCCCA-1 AAACACCAATAACTGC-1 ... #> TTGTTTCATCCAGG-1 TTGTTTCCATACAACT-1 #> colData names(12): orig.ident nCount_Spatial ... spatial2 cell_type #> reducedDimNames(0): #> mainExpName: NULL #> altExpNames(0): ``` -------------------------------- ### Read example SingleCellExperiment data from figshare Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-parallelization-vignette.html Loads the preprocessed pancreatic endocrinogenesis dataset (top 1000 highly variable genes, single trajectory) from a figshare URL using readRDS and prints it. The output shows the SingleCellExperiment structure with 1000 genes and 2087 cells. ```r example_sce <- readRDS((url("https://figshare.com/ndownloader/files/40581992"))) print(example_sce) ``` -------------------------------- ### Print cell type table from example_sce Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-conditionEffect-vignette.html Prints the cell type counts in the example dataset. The output shows the number of cells per cell type, with CD14 Mono being the most abundant and Eryth the least. ```r print(table(colData(example_sce)$cell_type)) #> #> CD14 Mono CD4 Naive T CD4 Memory T CD16 Mono B CD8 T #> 4362 2504 1762 1044 978 814 #> T activated NK DC B Activated Mk pDC #> 633 619 472 388 236 132 #> Eryth #> 55 ``` -------------------------------- ### Construct data with construct_data() Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-DEanalysis-vignette.html Initializes the data structure for simulation using the example_sce object. Sets a random seed for reproducibility. ```R set.seed(123) example_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = NULL, spatial = NULL, other_covariates = NULL, corr_by = "1" ) ``` -------------------------------- ### Load required libraries Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-spatial-deconvolution.html Loads scDesign3, SingleCellExperiment, ggplot2, dplyr, viridis, IOBR, and scatterpie for the tutorial. ```R library(scDesign3) library(SingleCellExperiment) library(ggplot2) library(dplyr) library(viridis) library(IOBR) library(scatterpie) theme_set(theme_bw()) ``` -------------------------------- ### Load required libraries for scDesign3 DE benchmarking Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-DEanalysis-vignette.html Loads all R packages needed for the tutorial, including scDesign3, SingleCellExperiment, Seurat, scran, DESeq2, tradeSeq, and spatial DE tools. Run this before any other code in the vignette. ```R library(scDesign3) library(SingleCellExperiment) library(ggplot2) library(DuoClustering2018) library(Seurat) library(SeuratObject) library(scran) library(parallel) library(monocle3) library(DESeq2) library(BiocParallel) library(NBAMSeq) library(PseudotimeDE) library(dplyr) library(tradeSeq) library(reshape2) library(spatialDE) library(SPARK) library(scales) library(tidyr) theme_set(theme_bw()) ``` -------------------------------- ### Subset to top 100 genes Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-parallelization-vignette.html Reduces the example SingleCellExperiment to the first 100 genes to save computational time before running scdesign3. ```R example_sce <- example_sce[1:100, ] ``` -------------------------------- ### Subset to B cells and regulatory T cells Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-DEanalysis-vignette.html Filters the example SingleCellExperiment to only B cells and regulatory T cells, and adds a 'cell_type' factor column to the colData for downstream analysis. ```R selected_cells <- which(colData(example_sce)$phenoid %in% c("b.cells","regulatory.t")) example_sce <- example_sce[,selected_cells] colData(example_sce)$cell_type <- as.factor(colData(example_sce)$phenoid) example_sce head(colData(example_sce)) ``` -------------------------------- ### Run scdesign3 with default parallelization (mcmapply) Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-parallelization-vignette.html Uses the default parallelization from the 'parallel' package. No need to set the parallelization parameter. ```R set.seed(123) example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime, k = 10, bs = 'cr')", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE ) ``` -------------------------------- ### Run scdesign3 with example_sce and Gaussian copula Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/scdesign3.html Runs scdesign3 on the example_sce dataset with a negative binomial family, a pseudotime smooth term for mu, and a Gaussian copula. The output shows the progress messages and the sequence of fitting steps, from input data construction through marginal and copula fitting to new data generation. ```r data(example_sce) my_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime, bs = 'cr', k = 10)", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, edf_flexible = FALSE, corr_formula = "pseudotime", copula = "gaussian", if_sparse = TRUE, DT = TRUE, pseudo_obs = FALSE, ncell = 1000, return_model = FALSE ) #> Input Data Construction Start #> fitting ... #> | | | 0% | |====== | 9% | |============= | 18% | |=================== | 27% | |========================= | 36% | |================================ | 45% | |====================================== | 55% | |============================================= | 64% | |=================================================== | 73% | |========================================================= | 82% | |================================================================ | 91% | |======================================================================| 100% #> Input Data Construction End #> Start Marginal Fitting #> Marginal Fitting End #> Start Copula Fitting #> Convert Residuals to Multivariate Gaussian #> Converting End #> Copula group 1 starts #> Copula group 2 starts #> Copula group 3 starts #> Copula group 4 starts #> Copula Fitting End #> Start Parameter Extraction #> Parameter #> Extraction End #> Start Generate New Data #> Use Copula to sample a multivariate quantile matrix #> Sample Copula group 1 starts #> Sample Copula group 2 starts #> Sample Copula group 3 starts #> Sample Copula group 4 starts #> New Data Generating End ``` -------------------------------- ### Construct input dataset for step-by-step simulation Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-cellType-vignette.html Uses construct_data to prepare the input dataset for step-by-step simulation, specifying assay_use='counts', celltype='cell_type', and corr_by='cell_type'. ```r example_data <- construct_data( sce = sce, assay_use = "counts", celltype = "cell_type", pseudotime = NULL, spatial = NULL, other_covariates = NULL, corr_by = "cell_type" ) ``` -------------------------------- ### Read and print reference scATAC-seq data Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-scATACseq-vignette.html Downloads the example SingleCellExperiment object from figshare and prints its structure. The object contains 1133 peaks and 7034 cells with assays 'counts' and 'logcounts'. ```R example_sce <- readRDS((url("https://figshare.com/ndownloader/files/40581962"))) print(example_sce) #> class: SingleCellExperiment #> dim: 1133 7034 #> metadata(0): #> assays(2): counts logcounts #> rownames(1133): chr6-44025105-44028184 chr2-113581628-113594911 ... #> chr2-233174300-233175622 chr3-13015237-13015864 #> rowData names(0): #> colnames(7034): AAACGAAAGAGCGAAA-1 AAACGAAAGAGTTTGA-1 ... #> TTTGTGTTCTACTTTG-1 TTTGTGTTCTTGTGCC-1 #> colData names(31): orig.ident nCount_peaks ... ident cell_type #> reducedDimNames(2): LSI UMAP #> mainExpName: peaks #> altExpNames(0) ``` -------------------------------- ### Inspect colData of example_sce Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-conditionEffect-vignette.html Displays the first six rows of the colData, showing the condition column with values CTRL and the cell type annotations for the subsetted dataset. ```r head(colData(example_sce)) #> DataFrame with 6 rows and 8 columns #> orig.ident nCount_RNA nFeature_RNA stim #> #> AAACATACATTTCC.1 IMMUNE_CTRL 3017 248 CTRL #> AAACATACCAGAAA.1 IMMUNE_CTRL 2481 713 CTRL #> AAACATACCTCGCT.1 IMMUNE_CTRL 3420 850 CTRL #> AAACATACGGCATT.1 IMMUNE_CTRL 1581 557 CTRL #> AAACATTGCTTCGC.1 IMMUNE_CTRL 2536 669 CTRL #> AAACGCACTCGCCT.1 IMMUNE_CTRL 3563 908 CTRL #> seurat_annotations ident cell_type condition #> #> AAACATACATTTCC.1 CD14 Mono IMMUNE_CTRL CD14 Mono CTRL #> AAACATACCAGAAA.1 CD14 Mono IMMUNE_CTRL CD14 Mono CTRL #> AAACATACCTCGCT.1 CD14 Mono IMMUNE_CTRL CD14 Mono CTRL #> AAACATACGGCATT.1 CD14 Mono IMMUNE_CTRL CD14 Mono CTRL #> AAACATTGCTTCGC.1 CD14 Mono IMMUNE_CTRL CD14 Mono CTRL #> AAACGCACTCGCCT.1 CD14 Mono IMMUNE_CTRL CD14 Mono CTRL ``` -------------------------------- ### Run scdesign3 with pbmcapply parallelization Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-parallelization-vignette.html Sets parallelization to 'pbmcmapply' to show a progress bar during parallel execution. ```R set.seed(123) example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime, k = 10, bs = 'cr')", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE, parallelization = "pbmcmapply" ) ``` -------------------------------- ### Subset example_sce to top genes and two cell types Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-conditionEffect-vignette.html Subsets the example dataset to the top 100 genes and only the CD14 Mono and B cell types, then removes unused cell type factor levels to save computational time. ```r example_sce <- example_sce[1:100, colData(example_sce)$cell_type %in% c("CD14 Mono", "B")] ## Remove unused cell type levels colData(example_sce)$cell_type <- droplevels(colData(example_sce)$cell_type) ``` -------------------------------- ### Load required R packages for pseudotime GOF evaluation Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-pseudotimeGOF-vignette.html Loads the required R packages for the tutorial: scDesign3, dyngen, SingleCellExperiment, ggplot2, and dplyr. Sets the ggplot2 theme to black and white. ```r library(scDesign3) library(dyngen) library(SingleCellExperiment) library(ggplot2) library(dplyr) theme_set(theme_bw()) ``` -------------------------------- ### Subset Zhengmix4eq to top HVGs and TFs Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-copulaCompare-vignette.html Selects the top 100 highly variable genes, adds human transcription factors, and filters out housekeeping genes and genes starting with RP or TMSB. Subsets the SingleCellExperiment to these genes and prints the resulting object, showing 119 genes and 3555 cells. ```r humantfs <- read_csv("http://humantfs.ccbr.utoronto.ca/download/v_1.01/TF_names_v_1.01.txt", col_names = FALSE) stats <- modelGeneVar(Zhengmix4eq_sce) chosen <- getTopHVGs(stats, n = 100) ### Filter out some HKGs. chosen <- union(chosen, humantfs$X1) chosen <- chosen[!stringr::str_starts(chosen, "RP")] chosen <- chosen[!stringr::str_starts(chosen, "TMSB")] chosen <- chosen[!chosen %in% c("B2M", "MALAT1", "ACTB", "ACTG1", "GAPDH", "FTL", "FTH1")] Zhengmix4eq_sce_sub <- Zhengmix4eq_sce[rownames(Zhengmix4eq_sce) %in% chosen, ] print(Zhengmix4eq_sce_sub) #> class: SingleCellExperiment #> dim: 119 3555 #> metadata(1): log.exprs.offset #> assays(3): counts logcounts normcounts #> rownames(119): TPT1 CD74 ... PHF1 IRF7 #> rowData names(10): id ... total_counts log10_total_counts #> colnames(3555): b.cells6276 b.cells6144 ... regulatory.t1084 #> regulatory.t9696 #> colData names(16): dataset barcode ... sizeFactor cell_type #> reducedDimNames(2): PCA TSNE #> mainExpName: NULL #> altExpNames(0): ``` -------------------------------- ### Read and inspect reference spatial data (example_sce) Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-spatial-vignette.html Downloads and reads a SingleCellExperiment object from figshare containing Visium spatial transcriptomics data with 1000 genes and 2696 cells. Prints the object to show its structure. ```R example_sce <- readRDS((url("https://figshare.com/ndownloader/files/40582019"))) print(example_sce) ``` -------------------------------- ### Initialize q-values matrix for DE methods Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-DEanalysis-vignette.html Creates a matrix to store BH-corrected p-values from four DE methods: Monocle3-DE, tradeSeq, NBAMseq, and PseudotimeDE. ```r qvals <- matrix(0, ncol = 4, nrow = dim(simu_sce)[1]) colnames(qvals) <- c("Monocle3-DE","tradeSeq","NBAMseq","PseudotimeDE") rownames(qvals) <- rownames(simu_sce) ``` -------------------------------- ### Construct synthetic spatial data with construct_data() Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-DEanalysis-vignette.html Builds a data structure for spatial simulation using the cell type and spatial coordinates from example_sce. Sets the seed to 1 for reproducibility and uses "counts" as the assay. ```R set.seed(1) example_data <- construct_data( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = NULL, spatial = c("spatial1", "spatial2"), other_covariates = NULL, corr_by = "1" ) ``` -------------------------------- ### Run scdesign3 with BiocParallel using SnowParam Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-parallelization-vignette.html Uses pbmapply from BiocParallel with SnowParam for parallelization. ```R set.seed(123) example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime, k = 10, bs = 'cr')", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE, parallelization = "pbmapply", BPPARAM = SnowParam() ) ``` -------------------------------- ### Load required R packages for copula comparison Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-copulaCompare-vignette.html Loads the required R packages for the tutorial, including scDesign3, DuoClustering2018, scater, tidygraph, ggrepel, scran, dplyr, igraph, tidyverse, and rvinecopulib. Sets the default ggplot2 theme to black and white. ```r library(scDesign3) library(DuoClustering2018) library(scater) library(tidygraph) library(ggrepel) library(scran) library(dplyr) library(igraph) library(tidyverse) library(rvinecopulib) theme_set(theme_bw()) ``` -------------------------------- ### Generate reference dataset with dyngen and store ground truth pseudotime Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-pseudotimeGOF-vignette.html Generates a reference dataset with ground truth pseudotime using dyngen. Sets a random seed for reproducibility, creates a linear backbone, initializes the model with 500 cells and specified numbers of transcription factors, target genes, and housekeeping genes, then generates a SingleCellExperiment object and stores the true pseudotime in its colData. ```r set.seed(123) backbone <- backbone_linear_simple() config <- initialise( backbone = backbone, num_cells = 500, num_tfs = nrow(backbone$module_info), num_targets = 100, num_hks = 50, verbose = FALSE ) out <- generate_dataset( config, format = "sce", make_plots = FALSE ) example_sce <- out$dataset colData(example_sce)$pseudotime <- out$model$experiment$cell_info$time ``` -------------------------------- ### Load scDesign3 and ggplot2 libraries Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-introduction-vignette.html Loads the required packages and sets a black-and-white theme for plots. ```R library(scDesign3) library(ggplot2) theme_set(theme_bw()) ``` -------------------------------- ### Run scdesign3 with BiocParallel using MulticoreParam Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-parallelization-vignette.html Uses pbmapply from BiocParallel with MulticoreParam for parallelization. ```R set.seed(123) example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = "pseudotime", spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime, k = 10, bs = 'cr')", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE, parallelization = "pbmapply", BPPARAM = MulticoreParam() ) ``` -------------------------------- ### Load example_sce dataset Source: https://songdongyuan1994.github.io/scDesign3/docs/reference/example_sce.html Loads the example_sce dataset, a SingleCellExperiment object containing both cell type and pseudotime information. The dataset has 10 rows (genes) and 1289 columns (cells). ```r data("example_sce") ``` -------------------------------- ### Show session information Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-librarySize-vignette.html Display the R session information, including version, platform, and locale. ```R sessionInfo() #> R version 4.3.1 (2023-06-16) #> Platform: x86_64-pc-linux-gnu (64-bit) #> Running under: Ubuntu 20.04.6 LTS #> #> Matrix products: default #> BLAS: /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 #> LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/liblapack.so.3; LAPACK version 3.9.0 #> #> locale: #> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C #> [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 #> [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_TELEPHONE=C #> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C #> #> time zone: America/Los_Angeles #> tzcode source: system (glibc) #> ``` -------------------------------- ### Run scDesign3 simulation with spatial coordinates Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-spatial-vignette.html Calls scdesign3() to simulate new spatial data. Uses a Gaussian process smooth term for spatial coordinates in mu_formula, negative binomial family, and parallelization with pbmcapply. Sets a seed for reproducibility. ```R set.seed(123) example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = NULL, spatial = c("spatial1", "spatial2"), other_covariates = NULL, mu_formula = "s(spatial1, spatial2, bs = 'gp', k= 400)", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE, parallelization = "pbmcapply" ) ``` -------------------------------- ### Load required packages and set theme Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-librarySize-vignette.html Load required R packages and set the default ggplot2 theme to black and white. ```R library(scDesign3) library(SingleCellExperiment) library(ggplot2) library(DuoClustering2018) library(scran) library(tidyverse) theme_set(theme_bw()) ``` -------------------------------- ### Prepare training data with top 100 HVGs and PCAKmeans results Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-clusterGOF-vignette.html Extracts the PCAKmeans clustering results, selects the top 100 highly variable genes, and splits the data into a training set of half the cells with a fixed seed. ```r kmeans_res <- Filter(function(x){all(x$method == "PCAKmeans")}, res_sub_list) ncell <- ncol(Zhengmix4eq_sce) ngene <- 100 zheng_sce <- modelGeneVar(Zhengmix4eq_sce) chosen <- getTopHVGs(zheng_sce, n = ngene) sce <- Zhengmix4eq_sce[chosen, ] ntrain <- round(ncell/2) set.seed(123) train_index <- sample(seq_len(ncell), ntrain, replace = FALSE) train_sce <- sce[, train_index] ``` -------------------------------- ### Load required libraries for scDesign3 spatial simulation Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-spatial-vignette.html Loads the scDesign3, SingleCellExperiment, ggplot2, dplyr, and viridis packages and sets the ggplot2 theme to theme_bw(). Required for the tutorial. ```R library(scDesign3) library(SingleCellExperiment) library(ggplot2) library(dplyr) library(viridis) theme_set(theme_bw()) ``` -------------------------------- ### Simulate data with scDesign3 using library size offset Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-librarySize-vignette.html Run scDesign3 simulation with library size as an offset in the mu formula. The library size is log-transformed because the link function for mu of the negative binomial distribution in GAMLSS is log. ```R set.seed(123) example_simu <- scdesign3( sce = sce, assay_use = "counts", celltype = "cell_type", pseudotime = NULL, spatial = NULL, other_covariates = "library", mu_formula = "cell_type + offset(log(library))", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE, parallelization = "pbmcmapply", important_feature = "auto" ) ``` -------------------------------- ### Simulate multiple-lineage data with scdesign3 Source: https://songdongyuan1994.github.io/scDesign3/docs/articles/scDesign3-multipleLineages-vignette.html Runs scdesign3 with a mu_formula containing two smooth terms for pseudotime1 and pseudotime2, each weighted by lineage indicators l1 and l2. Uses negative binomial family, Gaussian copula, and DT=TRUE. Set seed 123 for reproducibility. ```R set.seed(123) example_simu <- scdesign3( sce = example_sce, assay_use = "counts", celltype = "cell_type", pseudotime = c("pseudotime1", "pseudotime2", "l1", "l2"), spatial = NULL, other_covariates = NULL, mu_formula = "s(pseudotime1, k = 10, by = l1, bs = 'cr') + s(pseudotime2, k = 10, by = l2, bs = 'cr')", sigma_formula = "1", family_use = "nb", n_cores = 2, usebam = FALSE, corr_formula = "1", copula = "gaussian", DT = TRUE, pseudo_obs = FALSE, return_model = FALSE, nonzerovar = FALSE ) ```