Contents

1 Introduction

This package provides a feature selection method for single-cell RNA-seq data. It encodes the category number and calculate the spearman correlation coefficient.

The encoding the category number (CAEN) method is developed for selecting feature genes, which is differentially expressed between classes. We have implemented the CAEN method via a set of R functions. We make a R package named CAEN, and give a tutorial for the package. The method consist three steps.

Step 1: Data Pre-processing;

Step 2: Encoding the category number; Calculating the spearman correlation coefficient between the gene and category number;

Step 3: Calculate the classification error using the genes selected.

We employ a simulation dataset to illustrate the usage of the CAEN package. The programs can run under the Linux system and windows 10 system. The R versions should larger than 4.0.

2 Preparations

To install the CAEN package into your R environment, start R and enter:

install.packages("BiocManager")
BiocManager::install("CAEN")

Then, the CAEN package is ready to load.

library(CAEN)
library(SummarizedExperiment)

3 Data format

In order to show the CAEN workflow, the package includes the example data sets, which is a real dataset from Woo, Y., Kim, D., Koo, N., Kim, Y., Lee, S., Ko, J., et al. (2017). Profiling of mirnas and target genes related to cystogenesis in adpkd mouse models. scientific reports 7, 14151.
The next we will show the dataset in the package.

The dataset is in the data subdirectory of CAEN package, to consistent with standard Bioconductor representations, we transform the format of dataset as SummarizedExperiment, please refer R package * SummarizedExperiment* for details. In the dataset, gene expression is compared between two mouse models, Pkd1f/f: HoxB7-cre mice and Pkd2f/f: HoxB7-cre mice. Each group includes 18 samples, and there are total 29996 transcripts in this dataset. The next we will load the data and use CAEN function to find the most differentially expressed genes.

data("realData")
realData
## class: SummarizedExperiment 
## dim: 29996 36 
## metadata(0):
## assays(1): ''
## rownames(29996): NM_001001144 NM_001001152 ... NM_213730 NM_213733
## rowData names(0):
## colnames(36): Pkd1f_P1_KO1 Pkd1f_P1_KO2 ... Pkd2f_P7_WT2 Pkd2f_P7_WT3
## colData names(0):

readData includes 2 groups, each group includes 18 samples, and there are total 29996 transcripts. * the columns are samples of dataset; * the rows are transcripts of dataset.

x <- realData                  
y <- c(rep(1,18),rep(2,18))      
xte <- realData                

prob <- estimatep(x = x, y = y, xte = x, beta = 1, type = "mle", prior = NULL)     
prob0 <- estimatep(x = x, y = y, xte = xte, beta = 1, type = "mle", 
                   prior = NULL)   
myscore <- CAEN(dataTable = x, y = y, K = 2, gene_no_list = 100)

The output of the function CAEN is: A list of computed correlation coefficient and the first some differentially expressed genes , where “r” represents correlation coefficient between gene and category number, and “np” represents the top differential feature label.

4 Calculate classification error rate using genes selected with CAEN method

Getting the important gene, we Calculate classification error rate using genes selected. The step is as follows:

ddd <- myscore$np
datx <- t(assay(x)[ddd,])
datxte <- t(assay(xte)[ddd,])
probb <- prob[ddd,]
probb0 <- prob0[ddd,]

zipldacv.out <- ZIPLDA.cv(x = datx, y = y, prob0 = t(probb))
## 12345
ZIPLDA.out <- ZIPLDA(x = datx, y = y, xte = datxte, transform = FALSE, 
                     prob0 = t(probb0),rho = zipldacv.out$bestrho)
classResult <- ZIPLDA.out$ytehat

In order to reproduce the presented CAEN workflow, the package also includes the simulated data sets, which is generated by function newCountDataSet(). The next we will give an example for how to generate simulation dataset:

dat <- newCountDataSet(n = 100, p = 500, K = 4, param = 10, sdsignal = 2, 
                       drate = 0.2)

The output of the function newCountDataSet() includes: * “sim_train_data” represents training data of qn data matrix.
“sim_test_data” represents test data of qn data matrix.
The colnames of this two matrix are class labels for the n observations. May have q<p because features with 0 total counts are removed. The q features are those with >0 total counts in dataset. So q <= p.
“truesf” denotes size factors for training observations.
* “isDE” represnts the differential gene label.

Calculate the spearman correlation coefficient for data For the category number, we need to consider not only the difference between class but also the Intra-category difference. Therefore, we propose CAEN method, by encoding the category number within class, it get the optimal category number and select the most important genes used for classification.

x <- t(assay(dat$sim_train_data))                  
y <- as.numeric(colnames(dat$sim_train_data))      
xte <- t(assay(dat$sim_test_data))                 

prob <- estimatep(x = x, y = y, xte = x, beta = 1, 
                  type = c("mle","deseq","quantile"),
                  prior = NULL)      
prob0 <- estimatep(x = x, y = y, xte = xte, beta = 1, 
                   type = c("mle","deseq","quantile"),
                   prior = NULL)   
myscore <- CAEN(dataTable = assay(dat$sim_train_data), 
                y = as.numeric(colnames(dat$sim_train_data)), K = 4,
                gene_no_list = 100)

The output of the function CAEN is: A list of computed correlation coefficient and the first some differentially expressed genes , where “r” represents correlation coefficient between gene and category number, and “np” represents the top differential feature label.

Calculate classification error rate using genes selected with CAEN method Getting the important gene, we Calculate classification error rate using genes selected. The step is as follows:

ddd <- myscore$np
datx <- x[,ddd]
datxte <- xte[,ddd]
probb <- prob[ddd,]
probb0 <- prob0[ddd,]

zipldacv.out <- ZIPLDA.cv(x = datx, y = y, prob0 = t(probb))
## 12345
ZIPLDA.out <- ZIPLDA(x = datx, y = y,
                     xte = datxte, transform = FALSE, prob0 = t(probb0),
                     rho = zipldacv.out$bestrho)
classResult <- ZIPLDA.out$ytehat
sessionInfo()
## R Under development (unstable) (2023-10-22 r85388)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 22.04.3 LTS
## 
## Matrix products: default
## BLAS:   /home/biocbuild/bbs-3.19-bioc/R/lib/libRblas.so 
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.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/New_York
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats4    stats     graphics  grDevices utils     datasets  methods  
## [8] base     
## 
## other attached packages:
##  [1] SummarizedExperiment_1.33.0 Biobase_2.63.0             
##  [3] GenomicRanges_1.55.0        GenomeInfoDb_1.39.0        
##  [5] IRanges_2.37.0              S4Vectors_0.41.0           
##  [7] BiocGenerics_0.49.0         MatrixGenerics_1.15.0      
##  [9] matrixStats_1.0.0           CAEN_1.11.0                
## [11] BiocStyle_2.31.0           
## 
## loaded via a namespace (and not attached):
##  [1] Matrix_1.6-1.1          jsonlite_1.8.7          compiler_4.4.0         
##  [4] BiocManager_1.30.22     crayon_1.5.2            bitops_1.0-7           
##  [7] jquerylib_0.1.4         yaml_2.3.7              fastmap_1.1.1          
## [10] lattice_0.22-5          R6_2.5.1                XVector_0.43.0         
## [13] S4Arrays_1.3.0          knitr_1.44              DelayedArray_0.29.0    
## [16] bookdown_0.36           GenomeInfoDbData_1.2.11 bslib_0.5.1            
## [19] rlang_1.1.1             cachem_1.0.8            xfun_0.40              
## [22] sass_0.4.7              PoiClaClu_1.0.2.1       SparseArray_1.3.0      
## [25] cli_3.6.1               zlibbioc_1.49.0         digest_0.6.33          
## [28] grid_4.4.0              evaluate_0.22           abind_1.4-5            
## [31] RCurl_1.98-1.12         rmarkdown_2.25          tools_4.4.0            
## [34] htmltools_0.5.6.1