## ----install, eval = FALSE------------------------------------------------- # if (!requireNamespace("BiocManager", quietly = TRUE)) # install.packages("BiocManager") # # BiocManager::install("TOAST") ## ----quick_start, eval = FALSE--------------------------------------------- # Design_out <- makeDesign(design, Prop) # fitted_model <- fitModel(Design_out, Y_raw) # fitted_model$all_coefs # list all phenotype names # fitted_model$all_cell_types # list all cell type names # # coef should be one of above listed phenotypes # # cell_type should be one of above listed cell types # res_table <- csTest(fitted_model, coef = "age", # cell_type = "Neuron", contrast_matrix = NULL) # head(res_table) ## ----loadData-------------------------------------------------------------- library(TOAST) data("RA_100samples") Y_raw <- RA_100samples$Y_raw Pheno <- RA_100samples$Pheno Blood_ref <- RA_100samples$Blood_ref ## ----checkData------------------------------------------------------------- dim(Y_raw) Y_raw[1:4,1:4] ## ----checkPheno------------------------------------------------------------ dim(Pheno) head(Pheno, 3) ## ----checkRef-------------------------------------------------------------- dim(Blood_ref) head(Blood_ref, 3) ## ----SelFeature------------------------------------------------------------ refinx <- findRefinx(Y_raw, nmarker = 1000) ## ----Subset---------------------------------------------------------------- Y <- Y_raw[refinx,] Ref <- as.matrix(Blood_ref[refinx,]) ## ----DB2------------------------------------------------------------------- library(EpiDISH) outT <- epidish(beta.m = Y, ref.m = Ref, method = "RPC") estProp_RB <- outT$estF ## ---- DB3------------------------------------------------------------------ refinx = findRefinx(Y_raw, nmarker=1000, sortBy = "var") ## ----DF, message=FALSE, warning=FALSE-------------------------------------- library(RefFreeEWAS) ## ----DF2------------------------------------------------------------------- refinx <- findRefinx(Y_raw, nmarker = 1000) Y <- Y_raw[refinx,] ## ---- DF3, results='hide', message=FALSE, warning=FALSE-------------------- K <- 6 outT <- RefFreeCellMix(Y, mu0=RefFreeCellMixInitialize(Y, K = K)) estProp_RF <- outT$Omega ## ----compareRFRB----------------------------------------------------------- # first we align the cell types from RF # and RB estimations using pearson's correlation estProp_RF <- assignCellType(input=estProp_RF, reference=estProp_RB) mean(diag(cor(estProp_RF, estProp_RB))) ## ----IRB-RFCM1, message=FALSE, warning=FALSE------------------------------- library(TOAST) ## ----IRB-RFCM2, results='hide', message=FALSE, warning=FALSE--------------- K=6 set.seed(1234) outRF1 <- csDeconv(Y_raw, K, TotalIter = 30, bound_negative = TRUE) ## ----IRB-RFCM3, message=FALSE, warning=FALSE------------------------------- ## check the accuracy of deconvolution estProp_RF_improved <- assignCellType(input=outRF1$estProp, reference=estProp_RB) mean(diag(cor(estProp_RF_improved, estProp_RB))) ## ----initFeature, eval = FALSE--------------------------------------------- # refinx <- findRefinx(Y_raw, nmarker = 1000, sortBy = "cv") # InitNames <- rownames(Y_raw)[refinx] # csDeconv(Y_raw, K = 6, nMarker = 1000, # InitMarker = InitNames, TotalIter = 30) ## ---- eval = FALSE--------------------------------------------------------- # mydeconv <- function(Y, K){ # if (is(Y, "SummarizedExperiment")) { # se <- Y # Y <- assays(se)$counts # } else if (!is(Y, "matrix")) { # stop("Y should be a matrix # or a SummarizedExperiment object!") # } # # if (K<0 | K>ncol(Y)) { # stop("K should be between 0 and N (samples)!") # } # outY = RefFreeEWAS::RefFreeCellMix(Y, # mu0=RefFreeEWAS::RefFreeCellMixInitialize(Y, # K = K)) # Prop0 = outY$Omega # return(Prop0) # } # set.seed(1234) # outT <- csDeconv(Y_raw, K, FUN = mydeconv, bound_negative = TRUE) ## ----csDE2----------------------------------------------------------------- head(Pheno, 3) design <- data.frame(disease = as.factor(Pheno$disease)) Prop <- estProp_RF_improved colnames(Prop) <- colnames(Ref) ## columns of proportion matrix should have names ## ----csDE3----------------------------------------------------------------- Design_out <- makeDesign(design, Prop) ## ----csDE4----------------------------------------------------------------- fitted_model <- fitModel(Design_out, Y_raw) # print all the cell type names fitted_model$all_cell_types # print all phenotypes fitted_model$all_coefs ## -------------------------------------------------------------------------- res_table <- csTest(fitted_model, coef = "disease", cell_type = "Gran") head(res_table, 3) Disease_Gran_res <- res_table ## ---- eval = FALSE--------------------------------------------------------- # res_table <- csTest(fitted_model, # coef = "disease", # cell_type = "joint") # head(res_table, 3) ## ----joint, eval = FALSE--------------------------------------------------- # res_table <- csTest(fitted_model, # coef = "disease", # cell_type = NULL) # lapply(res_table, head, 3) # # ## this is exactly the same as # res_table <- csTest(fitted_model, coef = "disease") ## ----general2-------------------------------------------------------------- design <- data.frame(age = Pheno$age, gender = as.factor(Pheno$gender), disease = as.factor(Pheno$disease)) Prop <- estProp_RF_improved colnames(Prop) <- colnames(Ref) ## columns of proportion matrix should have names ## ----general3-------------------------------------------------------------- Design_out <- makeDesign(design, Prop) ## ----general4-------------------------------------------------------------- fitted_model <- fitModel(Design_out, Y_raw) # print all the cell type names fitted_model$all_cell_types # print all phenotypes fitted_model$all_coefs ## ----general5-------------------------------------------------------------- res_table <- csTest(fitted_model, coef = "age", cell_type = "Gran") head(res_table, 3) ## ----general6-------------------------------------------------------------- res_table <- csTest(fitted_model, coef = "disease", cell_type = "Bcell") head(res_table, 3) ## -------------------------------------------------------------------------- res_table <- csTest(fitted_model, coef = c("gender", 2, 1), cell_type = "CD4T") head(res_table, 3) ## ---- eval = FALSE--------------------------------------------------------- # res_table <- csTest(fitted_model, # coef = "age", # cell_type = "joint") # head(res_table, 3) ## ---- eval = FALSE--------------------------------------------------------- # res_table <- csTest(fitted_model, # coef = "age", # cell_type = NULL) # lapply(res_table, head, 3) # # ## this is exactly the same as # res_table <- csTest(fitted_model, # coef = "age") ## ----crossCellType2-------------------------------------------------------- design <- data.frame(age = Pheno$age, gender = as.factor(Pheno$gender), disease = as.factor(Pheno$disease)) Prop <- estProp_RF_improved colnames(Prop) <- colnames(Ref) ## columns of proportion matrix should have names ## ----crossCellType3, eval = FALSE------------------------------------------ # design <- data.frame(disease = as.factor(rep(0,100))) ## ----crossCellType4-------------------------------------------------------- Design_out <- makeDesign(design, Prop) ## ----crossCellType5-------------------------------------------------------- fitted_model <- fitModel(Design_out, Y_raw) # print all the cell type names fitted_model$all_cell_types # print all phenotypes fitted_model$all_coefs ## -------------------------------------------------------------------------- test <- csTest(fitted_model, coef = c("disease", 1), cell_type = c("CD8T", "Bcell"), contrast_matrix = NULL) head(test, 3) ## -------------------------------------------------------------------------- test <- csTest(fitted_model, coef = c("disease", 0), cell_type = c("CD8T", "Bcell"), contrast_matrix = NULL) head(test, 3) ## -------------------------------------------------------------------------- test <- csTest(fitted_model, coef = "joint", cell_type = c("Gran", "CD4T"), contrast_matrix = NULL) head(test, 3) ## -------------------------------------------------------------------------- test <- csTest(fitted_model, coef = NULL, cell_type = c("Gran", "CD4T"), contrast_matrix = NULL) lapply(test, head, 3) ## -------------------------------------------------------------------------- test <- csTest(fitted_model, coef = "disease", cell_type = c("Gran", "CD4T"), contrast_matrix = NULL) head(test, 3) ## -------------------------------------------------------------------------- test <- csTest(fitted_model, coef = "gender", cell_type = c("Gran", "CD4T"), contrast_matrix = NULL) head(test, 3) ## ----sessionInfo, echo=FALSE----------------------------------------------- sessionInfo()