-
Notifications
You must be signed in to change notification settings - Fork 3
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
- Loading branch information
Showing
28 changed files
with
806 additions
and
394 deletions.
There are no files selected for viewing
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| Original file line number | Diff line number | Diff line change |
|---|---|---|
|
|
@@ -553,6 +553,7 @@ CreatePRECASTObject <- function(seuList, project = "PRECAST", gene.number=2000 | |
|
|
||
|
|
||
| seulist <- lapply(seuList, function(x) x[genelist, ]) | ||
| ## reset the variable features to the genelist. | ||
| if(verbose) | ||
| message("Filter spots and features from SVGs(HVGs) count data...") | ||
| tstart <- Sys.time() | ||
|
|
@@ -669,200 +670,8 @@ dfList2df <- function(dfList){ | |
| return(df) | ||
| } | ||
|
|
||
| get_correct_exp <- function(XList, RfList, houseKeep, covariateList=NULL, q_unwanted=10){ | ||
|
|
||
|
|
||
| if(!all(sapply(XList, is.matrix))){ | ||
| XList <- lapply(XList, as.matrix) | ||
| } | ||
| XList_sub <- pbapply::pblapply(XList, function(x) x[,houseKeep]) | ||
| M0 <- wpca(matlist2mat(XList_sub), q=q_unwanted, FALSE)$PCs | ||
|
|
||
|
|
||
| Rf <- matlist2mat(RfList) | ||
| colnames(Rf) <- paste0("Rf",1:ncol(Rf)) | ||
| rm(RfList) | ||
| if(!is.null(covariateList)){ | ||
| # covariates <- matlist2mat(covariateList) | ||
| # covariates <- as.matrix(covariates) | ||
| covarites_df <- dfList2df(covariateList) | ||
| covariates <- model.matrix.lm(object = ~.+1, data = covarites_df, na.action = "na.pass") | ||
| rm(covariateList, covarites_df) | ||
| Rf <- cbind(Rf, covariates[,-1]) | ||
| rm(covariates) | ||
| } | ||
| ### XList <- lapply(XList, scale, scale=FALSE) | ||
| X0 <- matlist2mat(XList) | ||
| rm(XList) | ||
| nc_M0 <- ncol(M0) | ||
| lm1 <- lm(X0~ 0+ cbind(M0, Rf)) | ||
| coefmat <- coef(lm1)[c(1:nc_M0),] | ||
| #row.names(coef(lm1)) | ||
| rm(lm1) | ||
| hX <- X0 - M0 %*% coefmat | ||
| return(hX) | ||
| } | ||
|
|
||
|
|
||
| get_correct_mean_exp <- function(XList, hVList, covariateList=NULL){ | ||
|
|
||
| ## XList <- lapply(XList, scale, scale=FALSE) | ||
|
|
||
| if(!all(sapply(XList, is.matrix))){ | ||
| XList <- lapply(XList, as.matrix) | ||
| } | ||
|
|
||
| r_max <- length(XList) | ||
| X0 <- XList[[1]] | ||
| hV0 <- hVList[[1]] | ||
| if(r_max>1){ | ||
| for(r in 2:r_max){ | ||
| X0 <- rbind(X0, XList[[r]]) | ||
| hV0 <- rbind(hV0, hVList[[r]]) | ||
|
|
||
| } | ||
| } | ||
|
|
||
| rm(XList) | ||
| if(!is.null(covariateList)){ | ||
| # covariates <- matlist2mat(covariateList) | ||
| # covariates <- as.matrix(covariates) | ||
| # covariates <- cbind(1, covariates) | ||
| covarites_df <- dfList2df(covariateList) | ||
| covariates <- model.matrix.lm(object = ~.+1, data = covarites_df, na.action = "na.pass") | ||
| rm(covariateList, covarites_df) | ||
| }else{ | ||
| covariates <- matrix(1, nrow=nrow(hV0), ncol=1) | ||
| } | ||
|
|
||
| nc_M0 <- ncol(hV0) | ||
| lm1 <- lm(X0~ 0+ cbind(hV0, covariates)) | ||
| coefmat <- coef(lm1)[c(1:nc_M0),] | ||
| # row.names(coef(lm1)) | ||
| rm(lm1) | ||
| X0 - hV0 %*% coefmat | ||
|
|
||
|
|
||
| } | ||
|
|
||
| IntegrateSpaData <- function(PRECASTObj, species="Human", custom_housekeep=NULL, covariates_use=NULL){ | ||
| # suppressMessages(require(Matrix)) | ||
| # suppressMessages(require(Seurat)) | ||
|
|
||
| verbose <- PRECASTObj@parameterList$verbose | ||
| if(!inherits(PRECASTObj, "PRECASTObj")) | ||
| stop("IntegrateSpaData: Check the argument: PRECASTObj! PRECASTObj must be a PRECASTObj object.") | ||
| if(is.null(PRECASTObj@seulist)) | ||
| stop("IntegrateSpaData: Check the argument: PRECASTObj! The slot seulist in PRECASTObj is NULL!") | ||
|
|
||
|
|
||
| if(!tolower(species) %in% c("human", "mouse", "unknown")) | ||
| stop("IntegrateSpaData: Check the argument: species! it must be one of 'Human', 'Mouse' and 'Unknown'!") | ||
|
|
||
| defAssay_vec <- sapply(PRECASTObj@seulist, DefaultAssay) | ||
| if(any(defAssay_vec!=defAssay_vec[1])) warning("IntegrateSpaData: there are different default assays in PRECASTObj@seulist that will be used to integrating!") | ||
| n_r <- length(defAssay_vec) | ||
|
|
||
| XList <- lapply(1:n_r, function(r) Matrix::t(GetAssayData(PRECASTObj@seulist[[r]], assay = defAssay_vec[r], slot='data'))) | ||
|
|
||
| if(!is.null(covariates_use)){ | ||
| # covariateList <- lapply(PRECASTObj@seulist, function(x) [email protected][, covariates_use]) | ||
| covariateList <- lapply(PRECASTObj@seulist, function(x) x@meta.data[covariates_use]) | ||
| }else{ | ||
| covariateList <- NULL | ||
| } | ||
|
|
||
| if(tolower(species) =='mouse'){ | ||
| for(r in 1:n_r){ | ||
| colnames(XList[[r]]) <- firstup(colnames(XList[[r]])) | ||
| } | ||
| if(!is.null(custom_housekeep)){ | ||
| custom_housekeep <- firstup(custom_housekeep) | ||
| } | ||
| } | ||
| if(tolower(species) =='human'){ | ||
| for(r in 1:n_r){ | ||
| colnames(XList[[r]]) <- toupper(colnames(XList[[r]])) | ||
| } | ||
| if(!is.null(custom_housekeep)){ | ||
| custom_housekeep <- toupper(custom_housekeep) | ||
| } | ||
| } | ||
|
|
||
| barcodes_all <- lapply(XList, row.names) | ||
| if(any(duplicated(unlist(barcodes_all)))){ | ||
|
|
||
| for(r in 1:n_r){ | ||
| row.names(XList[[r]]) <- paste0(row.names(XList[[r]]), r) | ||
| } | ||
| } | ||
| genelist <- colnames(XList[[1]]) | ||
| lower_species <- tolower(species) | ||
| houseKeep <- switch (lower_species, | ||
| human = { | ||
| # data(Human_HK_genes) | ||
| intersect(toupper(genelist), PRECAST::Human_HK_genes$Gene) | ||
| }, | ||
| mouse={ | ||
| #data(Mouse_HK_genes) | ||
| intersect(firstup(genelist), PRECAST::Mouse_HK_genes$Gene) | ||
| }, | ||
| unknown={ | ||
| character() | ||
| } | ||
| ) | ||
| houseKeep <- c(houseKeep, custom_housekeep) | ||
| houseKeep <- intersect(houseKeep, colnames(XList[[1]])) | ||
|
|
||
| tstart <- Sys.time() | ||
| if(length(houseKeep) < 5){ | ||
| if(verbose){ | ||
| message("Using only PRECAST results to obtain the batch corrected gene expressions since species is unknown or the genelist in PRECASTObj has less than 5 overlapp with the housekeeping genes of given species.") | ||
| message("Start integration...") | ||
| } | ||
|
|
||
| hX <- get_correct_mean_exp(XList,PRECASTObj@resList$hV, covariateList=covariateList) | ||
| }else{ | ||
| if(verbose){ | ||
| message("Using bouth housekeeping gene and PRECAST results to obtain the batch corrected gene expressions.") | ||
| message("Start integration...") | ||
| } | ||
|
|
||
| hX <- get_correct_exp(XList, PRECASTObj@resList$Rf, houseKeep=houseKeep, q_unwanted=min(10, length(houseKeep)), covariateList=covariateList) | ||
| } | ||
| .logDiffTime(sprintf(paste0("%s Data integration finished!"), "*****"), t1 = tstart, verbose = verbose) | ||
|
|
||
| if(verbose) | ||
| message("Put the data into a new Seurat object...") | ||
| tstart <- Sys.time() | ||
| meta_data <- data.frame(batch=factor(get_sampleID(XList)), cluster= factor(unlist(PRECASTObj@resList$cluster))) | ||
| row.names(meta_data) <- row.names(hX) | ||
| rm(XList) | ||
| count <- sparseMatrix(i=1,j=1, x=0, dims=dim(t(hX))) | ||
| row.names(count) <- colnames(hX) | ||
| colnames(count) <- row.names(hX) | ||
| seuInt <- CreateSeuratObject(counts = count, assay = 'PRE_CAST', meta.data=meta_data) | ||
| if(inherits(seuInt[['PRE_CAST']], "Assay5") ){ | ||
|
|
||
| seuInt <- SetAssayData(object = seuInt, slot='data', assay = "PRE_CAST", new.data = t(hX)) | ||
| }else{ | ||
| seuInt[['PRE_CAST']]@data <- t(hX) | ||
| } | ||
|
|
||
| rm(hX) | ||
|
|
||
| # seuInt <- CreateSeuratObject(assay, meta.data=meta_data, assay = 'PRECAST') | ||
| seuInt <- Add_embed(matlist2mat(PRECASTObj@resList$hZ), seuInt, embed_name = 'PRECAST', assay='PRE_CAST') | ||
| posList <- lapply(PRECASTObj@seulist, function(x) cbind(x$row, x$col)) | ||
| seuInt <- Add_embed(matlist2mat(posList), seuInt, embed_name = 'position', assay='PRE_CAST') | ||
| Idents(seuInt) <- factor(meta_data$cluster) | ||
|
|
||
| .logDiffTime(sprintf(paste0("%s New Seurat object is generated!"), "*****"), t1 = tstart, verbose = verbose) | ||
|
|
||
|
|
||
| return(seuInt) | ||
| } | ||
|
|
||
|
|
||
| # plot cluster spatial heatmap, or tSNE RGB plot, or UMAP RGB plot | ||
| gg_color_hue <- function(n) { | ||
|
|
||
Oops, something went wrong.