Update utils.R

xuranw · Jul 24, 2018 · 27a51b7 · 27a51b7
1 parent 8dfe5da
commit 27a51b7
Showing 1 changed file with 24 additions and 23 deletions.
diff --git a/R/utils.R b/R/utils.R
@@ -240,54 +240,55 @@ music_prop = function(bulk.eset, sc.eset, markers = NULL, clusters, samples, sel
 #' @return Estimated proportions from MuSiC with cluster information.
 #' @seealso
 #' \code{\link{music_basis}}; \code{\link{music_prop}}
-music_prop.cluster = function(bulk.eset, sc.eset, group.markers, groups, clusters, samples, clusters.type,  
+music_prop.cluster = function(bulk.eset, sc.eset, group.markers, groups, clusters, samples, clusters.type,
                               verbose = TRUE, iter.max = 1000, nu = 0.0001, eps = 0.01, inter.as.bseq = F, normalize = F, ... ){
   bulk.gene = rownames(bulk.eset)[rowMeans(exprs(bulk.eset)) != 0]
   bulk.eset = bulk.eset[bulk.gene, , drop = FALSE]
   select.ct = unlist(clusters.type)
-  
+
   if(length(setdiff(names(group.markers), names(clusters.type))) > 0 || length(setdiff(names(clusters.type), names(group.markers))) > 0){
     stop("Cluster number is not matching!")
   }else{
     group.markers = group.markers[names(clusters.type)]
   }
-  
-  
+
+
   if(verbose){message('Start: cluster estimations...')}
   cluster.sc.basis = music_basis(sc.eset, non.zero = TRUE, markers = NULL, clusters = groups, samples = samples, select.ct = names(clusters.type), verbose = verbose)
   if(verbose){message('Start: cell type estimations...')}
   sc.basis = music_basis(sc.eset, non.zero = TRUE, markers = NULL, clusters = clusters, samples = samples, select.ct = select.ct, verbose = verbose)
   cm.gene = intersect(rownames(sc.basis$Disgn.mtx), bulk.gene)
-  
+
   if(length(cm.gene)< 0.2*min(length(bulk.gene), nrow(sc.eset)) ){
     stop("Too few common genes!")
   }
-  
+
   if(verbose){message(paste('Used', length(cm.gene), 'common genes...'))}
-  
+
   m.sc = match(cm.gene, rownames(sc.basis$Disgn.mtx)); m.bulk = match(cm.gene, bulk.gene)
   group.markers = lapply(group.markers, intersect, cm.gene)
-  
+
   D1 = sc.basis$Disgn.mtx[m.sc,]; M.S = sc.basis$M.S; Sigma = sc.basis$Sigma[m.sc, ]
-  
+
   cluster.select = setdiff(rownames(D1), unique(unlist(group.markers)))
   cluster.diff = unique(unlist(group.markers))
-  
+
   D1.cluster = cluster.sc.basis$Disgn.mtx[cluster.select, ]; M.S.cluster = cluster.sc.basis$M.S;
   Yjg = relative.ab(exprs(bulk.eset)[m.bulk, ]); N.bulk = ncol(bulk.eset);
-  
+
   Sigma.cluster = cluster.sc.basis$Sigma[cluster.select, ];
-  
+
   D1.sub = cluster.sc.basis$Disgn.mtx[cluster.diff, ]; Sigma.sub = cluster.sc.basis$Sigma[cluster.diff, ];
-  
+
   Est.prop.weighted.cluster = NULL
   for(i in 1:N.bulk){
     if(sum(Yjg[, i] == 0) > 0){
-      D1.cluster.temp = D1.cluster[Yjg[, i]!=0, ];
-      D1.sub.temp = D1.sub[Yjg[, i]!=0, ];
+      name.temp = rownames(Yjg)[Yjg[, i] != 0]
+      D1.cluster.temp = D1.cluster[rownames(D1.cluster)%in%name.temp, ];
+      D1.sub.temp = D1.sub[rownames(D1.sub)%in%name.temp, ];
       Yjg.temp = Yjg[Yjg[, i]!=0, i];
-      Sigma.cluster.temp = Sigma.cluster[Yjg[, i]!=0, ];
-      Sigma.sub.temp = Sigma.sub[Yjg[, i]!=0, ];
+      Sigma.cluster.temp = Sigma.cluster[rownames(Sigma.cluster)%in%name.temp, ];
+      Sigma.sub.temp = Sigma.sub[rownames(Sigma.sub)%in%name.temp, ];
       if(verbose) message(paste(colnames(Yjg)[i], 'has common genes', sum(Yjg[, i] != 0), '...') )
     }else{
       D1.cluster.temp = D1.cluster;
@@ -309,22 +310,22 @@ music_prop.cluster = function(bulk.eset, sc.eset, group.markers, groups, cluster
           p.weight = c(p.weight, rep(0, length(clusters.type[[j]])))
         }else{
           c.marker = intersect(group.markers[[j]], names(Yjg.temp))
-          Y.sub = D1.sub.temp[c.marker, j]*p.cluster.weight[j] + 
+          Y.sub = D1.sub.temp[c.marker, j]*p.cluster.weight[j] +
             (Yjg.temp[c.marker] - D1.sub.temp[c.marker, ]%*% p.cluster.weight) * p.cluster.weight[j]
           names(Y.sub) = c.marker
           Y.sub = Y.sub[Y.sub > 0]
-          lm.D1.sub = music.iter(Y.sub, D1[c.marker, clusters.type[[j]]], M.S[clusters.type[[j]]], 
-                                 Sigma[c.markers, clusters.type[[j]]], select.ct = clusters.type[[j]])
-          p.weight = c(p.weight, p.cluster.weight * lm.D1.sub$p.weight)
+          lm.D1.sub = music.iter(Y.sub, D1[c.marker, clusters.type[[j]]], M.S[clusters.type[[j]]],
+                                 Sigma[c.marker, clusters.type[[j]]])
+          p.weight = c(p.weight, p.cluster.weight[j] * lm.D1.sub$p.weight)
         }
       }
     }
     Est.prop.weighted.cluster = rbind(Est.prop.weighted.cluster, p.weight)
   }
-  
+
   colnames(Est.prop.weighted.cluster) = unlist(clusters.type)
   rownames(Est.prop.weighted.cluster) = colnames(Yjg)
-  
+
   return(list(Est.prop.weighted.cluster = Est.prop.weighted.cluster))
 }