cycle_npreg_outsample: Predict test-sample ordering using training labels (no...

In peco: A Supervised Approach for **P**r**e**dicting **c**ell Cycle Pr**o**gression using scRNA-seq data

Description

Apply the estimates of cycle_npreg_insample to another gene expression dataset to infer an angle or cell cycle phase for each cell.

Usage

cycle_npreg_outsample(Y_test, sigma_est, funs_est,
  method.trend = c("trendfilter", "loess", "bspline"), normed = TRUE,
  polyorder = 2, method.grid = "uniform", ncores = 2, grids = 100,
  get_trend_estimates = FALSE)

Arguments

Y_test	A SingleCellExperiment object.
sigma_est	Input from training data. A vector of gene-specific standard error of the cyclic trends.
funs_est	Input fron training data. A vector of cyclic functions estimating cyclic trends.
method.trend	Varous methods that can be applied to estimate cyclic trend of gene expression levels.
normed	Is the data already normalized? TRUE or FALSE.
polyorder	We estimate cyclic trends of gene expression levels using nonparamtric trend filtering. The default fits second degree polynomials.
method.grid	Method for defining bins along the circle.
ncores	We use doParallel package for parallel computing.
grids	number of bins to be selected along 0 to 2pi.
get_trend_estimates	To re-estimate the cylic trend based on the predicted cell cycle phase or not (T or F). Default FALSE. This step calls trendfilter and is computationally intensive.

Value

A list with the following elements:

Y	The input gene expression marix.
cell_times_est	Inferred angles or cell cycle phases, NOT ordered.
loglik_est	Log-likelihood estimates for each gene.
cell_times_reordered	The inferred angles reordered (in ascending order).
Y_reorded	The input gene expression matrix reordered by cell_times_reordered.
sigma_reordered	Estimated standard error of the cyclic trend for each gene, reordered by cell_times_reordered.
funs_reordered	A list of functions for approximating the cyclic trends of gene express levels for each gene, reordered by cell_times_reordered.
mu_reordered	Estimated cyclic trend of gene expression values for each gene, reordered by cell_times_reordered.
prob_per_cell_by_celltimes	Probabilities of each cell belong to each bin.

See Also

cycle_npreg_insample for obtaining parameteres for cyclic functions from training data, cycle_npreg_loglik for log-likehood at angles between 0 to 2pi, initialize_grids for selecting angles in cycle_npreg_loglik, cycle_npreg_mstep for estimating cyclic functions given inferred phases from cycle_npreg_loglik

Other peco classifier functions: cycle_npreg_insample, cycle_npreg_loglik, cycle_npreg_mstep, initialize_grids

Examples

# import data
library(SingleCellExperiment)
data(sce_top101genes)

# select top 5 cyclic genes
sce_top5 <- sce_top101genes[order(rowData(sce_top101genes)$pve_fucci,
                                  decreasing=TRUE)[1:5],]

# Select samples from NA18511 for our prediction example
coldata <- colData(sce_top5)
which_samples_train <- rownames(coldata)[coldata$chip_id != "NA18511"]
which_samples_predict <- rownames(coldata)[coldata$chip_id == "NA18511"]

# learning cyclic functions of the genes using our training data
sce_top5 <- data_transform_quantile(sce_top5)
expr_quant <- assay(sce_top5, "cpm_quantNormed")
Y_train <- expr_quant[, colnames(expr_quant) %in% which_samples_train]
theta_train <-
    coldata$theta_shifted[rownames(coldata) %in% which_samples_train]
names(theta_train) <-
    rownames(coldata)[rownames(coldata) %in% which_samples_train]

# obtain cyclic function estimates
model_5genes_train <- cycle_npreg_insample(Y = Y_train,
                                           theta = theta_train,
                                           polyorder=2,
                                           ncores=2,
                                           method.trend="trendfilter")

# predict cell cycle
model_5genes_predict <- cycle_npreg_outsample(
  Y_test=sce_top5[,colnames(sce_top5) %in% which_samples_predict],
  sigma_est=model_5genes_train$sigma_est,
  funs_est=model_5genes_train$funs_est,
  method.trend="trendfilter",
  ncores=2,
  get_trend_estimates=FALSE)

# estimate cyclic gene expression levels given cell cycle for each gene
predict_cyclic <-
    fit_cyclical_many(Y=assay(model_5genes_predict$Y,"cpm_quantNormed"),
                      theta=colData(model_5genes_predict$Y)$cellcycle_peco)
all.equal(names(predict_cyclic[[2]]), colnames(predict_cyclic[[1]]))

par(mfrow=c(2,3), mar=c(4,4,3,1))
for (g in seq_along(rownames(model_5genes_predict$Y))) {
  plot(assay(model_5genes_predict$Y,"cpm_quantNormed")[
      rownames(model_5genes_predict$Y)[g],],
       x=colData(model_5genes_predict$Y)$cellcycle_peco, axes=FALSE,
       xlab="FUCCI phase",
       ylab="Predicted phase")
  points(y=predict_cyclic$cellcycle_function[[
           rownames(model_5genes_predict$Y)[g]]](
    seq(0, 2*pi, length.out = 100)),
    x=seq(0, 2*pi, length.out = 100),
    pch=16, col="royalblue")
  axis(2); axis(1,at=c(0,pi/2, pi, 3*pi/2, 2*pi),
                labels=c(0,expression(pi/2), expression(pi),
                         expression(3*pi/2), expression(2*pi)))
  abline(h=0, lty=1, col="black", lwd=.7)
  title(rownames(model_5genes_predict$Y_reordered)[g])
}
title("Predicting cell cycle phase for NA18511", outer=TRUE)

peco documentation built on Nov. 8, 2020, 8:16 p.m.