Identify cluster specific network: Difference between revisions

From ZhangLabWiki
Jump to navigation Jump to search
>Shicheng
(Created page with "library(Seurat) setwd("/media/Home_Raid1/zhl002/NAS1/RNA_seq/hiseq_020617/seurat_analysis") load("./ipsnt_33k_latest3.RData") cluster<-as.matrix(pbmc33k.merged@data[,which(pbm...")
 
>Shicheng
No edit summary
 
Line 1: Line 1:
library(Seurat)
library(Seurat)
setwd("/media/Home_Raid1/zhl002/NAS1/RNA_seq/hiseq_020617/seurat_analysis")
setwd("/media/Home_Raid1/zhl002/NAS1/RNA_seq/hiseq_020617/seurat_analysis")
load("./ipsnt_33k_latest3.RData")
load("./ipsnt_33k_latest3.RData")
cluster<-as.matrix(pbmc33k.merged@data[,which(pbmc33k.merged@ident==1)])
cluster<-as.matrix(pbmc33k.merged@data[,which(pbmc33k.merged@ident==1)])
COR<-list()
COR<-list()
for(i in 1:26) {
for(i in 1:26) {
  cluster<-as.matrix(pbmc33k.merged@data[,which(pbmc33k.merged@ident==i)])
cluster<-as.matrix(pbmc33k.merged@data[,which(pbmc33k.merged@ident==i)])
  data<-(na.omit(cluster[match(as.character(pbmc33k.merged@var.genes),rownames(cluster)),]))
data<-(na.omit(cluster[match(as.character(pbmc33k.merged@var.genes),rownames(cluster)),]))
  COR[[i]]<-cor(t(data),method="spearman")
COR[[i]]<-cor(t(data),method="spearman")
  print(i)
print(i)
}
}
 
gene<-c()
gene<-c()
for(i in 1:nrow(COR[[1]])){
for(i in 1:nrow(COR[[1]])){
for(j in 1:nrow(COR[[1]])){
  for(j in 1:nrow(COR[[1]])){
for(z in 1:26){
    for(z in 1:26){
SD=sd(c(COR[[1]][i][j],COR[[2]][i][j],COR[[3]][i][j],
      SD=sd(c(COR[[1]][i][j],COR[[2]][i][j],COR[[3]][i][j],
COR[[4]][i][j],COR[[5]][i][j],COR[[6]][i][j],
            COR[[4]][i][j],COR[[5]][i][j],COR[[6]][i][j],
COR[[7]][i][j],COR[[8]][i][j],COR[[9]][i][j],
            COR[[7]][i][j],COR[[8]][i][j],COR[[9]][i][j],
COR[[10]][i][j],COR[[11]][i][j],COR[[12]][i][j],
            COR[[10]][i][j],COR[[11]][i][j],COR[[12]][i][j],
COR[[13]][i][j],COR[[14]][i][j],COR[[15]][i][j],
            COR[[13]][i][j],COR[[14]][i][j],COR[[15]][i][j],
COR[[16]][i][j],COR[[17]][i][j],COR[[18]][i][j],
            COR[[16]][i][j],COR[[17]][i][j],COR[[18]][i][j],
COR[[19]][i][j],COR[[20]][i][j],COR[[21]][i][j],
            COR[[19]][i][j],COR[[20]][i][j],COR[[21]][i][j],
COR[[22]][i][j],COR[[23]][i][j],COR[[24]][i][j],
            COR[[22]][i][j],COR[[23]][i][j],COR[[24]][i][j],
COR[[25]][i][j],COR[[26]][i][j],na.rm=T),na.rm=T)
            COR[[25]][i][j],COR[[26]][i][j],na.rm=T),na.rm=T)
if(SD>0.6){
      if(SD>0.6){
gene<-c(gene,i,j)
        gene<-c(gene,i,j)
}
      }
}
    }
}
       
}
  }
### example
}
#gene<-c()
 
#for(i in 1:nrow(COR[[1]])){
### example
# print(i)
#gene<-c()
# for(j in 1:nrow(COR[[1]])){
#for(i in 1:nrow(COR[[1]])){
# for(z in 1:4){
# print(i)
# SD=try(sd(c(COR[[1]][i,j],COR[[2]][i,j],COR[[3]][i,j],COR[[4]][i,j]),na.rm=T))
# for(j in 1:nrow(COR[[1]])){
# if(!is.na(SD) && SD>0.6){
#   for(z in 1:4){
# gene<-c(gene,i,j)
#     SD=try(sd(c(COR[[1]][i,j],COR[[2]][i,j],COR[[3]][i,j],COR[[4]][i,j]),na.rm=T))
# }
#     if(!is.na(SD) && SD>0.6){
# }
#       gene<-c(gene,i,j)
#
#     }
# }
#   }
#}
#  
# }
#}

Latest revision as of 02:08, 16 May 2017

library(Seurat)
setwd("/media/Home_Raid1/zhl002/NAS1/RNA_seq/hiseq_020617/seurat_analysis")
load("./ipsnt_33k_latest3.RData")
cluster<-as.matrix(pbmc33k.merged@data[,which(pbmc33k.merged@ident==1)])
COR<-list()
for(i in 1:26) {
cluster<-as.matrix(pbmc33k.merged@data[,which(pbmc33k.merged@ident==i)])
data<-(na.omit(cluster[match(as.character(pbmc33k.merged@var.genes),rownames(cluster)),]))
CORi<-cor(t(data),method="spearman")
print(i)
}
gene<-c()
for(i in 1:nrow(COR1)){
for(j in 1:nrow(COR1)){
for(z in 1:26){
SD=sd(c(COR1[i][j],COR2[i][j],COR3[i][j],
COR4[i][j],COR5[i][j],COR6[i][j],
COR7[i][j],COR8[i][j],COR9[i][j],
COR10[i][j],COR11[i][j],COR12[i][j],
COR13[i][j],COR14[i][j],COR15[i][j],
COR16[i][j],COR17[i][j],COR18[i][j],
COR19[i][j],COR20[i][j],COR21[i][j],
COR22[i][j],COR23[i][j],COR24[i][j],
COR25[i][j],COR26[i][j],na.rm=T),na.rm=T)
if(SD>0.6){
gene<-c(gene,i,j)
}
}
}
}
### example
#gene<-c()
#for(i in 1:nrow(COR1)){
# print(i)
# for(j in 1:nrow(COR1)){
# for(z in 1:4){
# SD=try(sd(c(COR1[i,j],COR2[i,j],COR3[i,j],COR4[i,j]),na.rm=T))
# if(!is.na(SD) && SD>0.6){
# gene<-c(gene,i,j)
# }
# }
#
# }
#}