Created
May 8, 2020 17:42
-
-
Save cplaisier/35aa2a89c67dfdf52fdf1d5653be0bc4 to your computer and use it in GitHub Desktop.
Example code to conduct AM and MDM analyses:
This file contains hidden or 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
| # Set working directory | |
| #setwd('C:/Users/cplaisie/Dropbox/Urdahl/Ag85B_ESAT6_ts/network') | |
| #setwd('/home/cplaisier/Dropbox (ASU)/Urdahl/macs/cMonkey') | |
| setwd('C:/Users/plais/Dropbox (ASU)/Urdahl/macs/cMonkey') | |
| #setwd('C:/Users/Christopher/Dropbox/Urdahl/Ag85B_ESAT6_ts/network') | |
| # Get eigengenes | |
| d1 = read.csv('eigengenes/biclusterEigengenes_tfbs_db.csv',header=T, row.names=1) | |
| rownames(d1) = paste('TFBS_DB_',rownames(d1),sep='') | |
| d2 = read.csv('eigengenes/biclusterEigengenes_pita.csv',header=T, row.names=1) | |
| rownames(d2) = paste('PITA_',rownames(d2),sep='') | |
| d3 = read.csv('eigengenes/biclusterEigengenes_targetscan.csv',header=T, row.names=1) | |
| rownames(d3) = paste('TargetScan_',rownames(d3),sep='') | |
| dAll = rbind(d1,d2,d3) | |
| # Get expression data | |
| ge1 = read.table('../gexp_cMonkey_norm_cv8000_preProcessed.tsv', header=T, row.names=1, sep='\t') | |
| # Ordering samples | |
| #os1 = c(c(1:4,9:12,16:19,23:25),c(5:8,13:15,20:22,26:28)) | |
| #ge2 = ge1[,os1] | |
| #dAll2 = dAll[,os1] | |
| ge2 = ge1 | |
| dAll2 = dAll[,colnames(ge1)] | |
| # Plot Bicluster vs. conditions | |
| library(pheatmap) | |
| library(RColorBrewer) | |
| library(gplots) | |
| sets = c(rep(1,4),rep(2,3),rep(3,4),rep(4,3)) | |
| names1 = c(rep('Uninf_AM',4),rep('Inf_AM',3),rep('Uninf_Mono',4),rep('Inf_Mono',4)) | |
| #tfs1 = c('uc007ldr.2', rep('uc007xlk.1',2),rep('uc009orz.1',7),rep('uc007ifo.1',2)) | |
| #tfs2 = c(c('TargetScan_331'), c('TargetScan_82', 'PITA_56'), c('TargetScan_2','PITA_241','PITA_312','TFBS_DB_169','PITA_234','TFBS_DB_274','PITA_47'),c('TargetScan_188','PITA_334')) | |
| tfNames2 = tfs2 | |
| pdf('tfBiclusters_noScale_lines.pdf') | |
| #par(mar=c(7,4,4,2),mfrow=c(1,2)) | |
| for(tf2 in 1:length(tfs2)) { | |
| par(mar=c(7,4,4,2),mfrow=c(1,2)) | |
| tmp2 = t(sapply(1:max(sets), function(x) { mean(as.numeric(ge2[tfs1[tf2],which(sets==x)])) } )) | |
| sdTmp2 = t(sapply(1:max(sets), function(x) { sd(as.numeric(ge2[tfs1[tf2],which(sets==x)])) } )) | |
| tmp1 = t(sapply(1:max(sets), function(x) { mean(as.numeric(dAll2[tfs2[tf2],which(sets==x)])) } )) | |
| sdTmp1 = t(sapply(1:max(sets), function(x) { sd(as.numeric(dAll2[tfs2[tf2],which(sets==x)])) } )) | |
| colnames(tmp1) = names1 | |
| ylim1 = c(min(c(tmp1,tmp2)),max(c(tmp1,tmp2))) | |
| plotCI(y=tmp1[1:4],x=c(0,2,4,8),uiw=sdTmp1[o1],gap=0,ylab='Bicluster Relative Expression',pch=19,col='black',main=paste(tfs1[tf2],' vs. ',tfs2[tf2],'\nESAT6',sep=''),xlab='Time Post Stim.',ylim=ylim1) | |
| lines(c(0,2,4,8),tmp1[1:4],lwd=1,col=rgb(0,0,0,0.5),lty=2) | |
| par(new=T) | |
| plotCI(y=tmp2[1:4],x=c(0,2,4,8),uiw=sdTmp2[o1],gap=0,ylab='',xaxt='n',xlab='',pch=19,col='red',yaxt='n') | |
| lines(c(0,2,4,8),tmp2[1:4],lwd=2,col=rgb(1,0,0,0.5)) | |
| plotCI(y=tmp1[5:8],x=c(0,2,4,8),uiw=sdTmp1[o1],gap=0,ylab='Bicluster Relative Expression',pch=19,col=rgb(0,0,0,1),main=paste(tfs1[tf2],' vs. ',tfs2[tf2],'\nAg85B',sep=''),xlab='Time Post Stim.',ylim=ylim1) | |
| lines(c(0,2,4,8),tmp1[5:8],lwd=1,col=rgb(0,0,0,0.5),lty=2) | |
| par(new=T) | |
| plotCI(y=tmp2[5:8],x=c(0,2,4,8),uiw=sdTmp2[o1],gap=0,ylab='',xaxt='n',xlab='',pch=19,col=rgb(1,0,0,1),yaxt='n') | |
| lines(c(0,2,4,8),tmp2[5:8],lwd=2,col=rgb(1,0,0,0.5)) | |
| } | |
| dev.off() | |
| #### GOOD ONE #### | |
| library(pheatmap) | |
| library(RColorBrewer) | |
| sets1 = c(rep(1,4),rep(2,4),rep(3,4),rep(4,3),rep(5,4),rep(6,3),rep(7,3),rep(8,3)) | |
| names1 = c('ESAT6_0H','ESAT6_2H','ESAT6_4H','ESAT6_8H','Ag85B_0H','Ag85B_2H','Ag85B_4H','Ag85B_8H') | |
| tfs1 = c('uc007ldr.2', rep('uc007xlk.1',2),rep('uc009orz.1',7),rep('uc007ifo.1',2),rep('uc012gkc.1',3)) | |
| tfs2 = c(c('TargetScan_331'), c('TargetScan_82', 'PITA_56'), c('TargetScan_2','PITA_241','PITA_312','TFBS_DB_169','PITA_234','TFBS_DB_274','PITA_47'),c('TargetScan_188','PITA_334'),c('TFBS_DB_71','PITA_50','TFBS_DB_309')) | |
| tfNames2 = tfs2 | |
| pdf('tfBiclusters_noScale.pdf') | |
| #par(mar=c(7,4,4,2)) | |
| for(tf2 in 1:length(tfs2)) { | |
| par(mfrow=c(2,1),mar=c(7,4,4,2)) | |
| tmp2 = t(sapply(1:max(sets1), function(x) { mean(as.numeric(ge2[tfs1[tf2],which(sets1==x)])) } )) | |
| colnames(tmp2) = names1 | |
| c1 = cor.test(as.numeric(ge2[tfs1[tf2],]),as.numeric(dAll2[tfs2[tf2],])) | |
| c2_ESAT6 = cor.test(as.numeric(ge2[tfs1[tf2],])[1:15],as.numeric(dAll2[tfs2[tf2],])[1:15]) | |
| c2_Ag85B = cor.test(as.numeric(ge2[tfs1[tf2],])[16:28],as.numeric(dAll2[tfs2[tf2],])[16:28]) | |
| boxplot(as.numeric(dAll2[tfs2[tf2],]) ~ sets1, names=names1, las=2,main=paste(tfNames2[tf2],':',tfs1[tf2],'\nR.all = ',round(c1$estimate,2),'; R.ESAT6 = ',round(c2_ESAT6$estimate,2),'; R.Ag85B = ',round(c2_Ag85B$estimate,2)), col=brewer.pal(8, "Paired")) | |
| par(new=T) | |
| plot(1:4,tmp2[1:4],col=rgb(1,0,0,0.3),type='l',lty=1,lwd=2,axes=F,xlim=c(0.5,8.5),ylab='',xlab='') | |
| par(new=T) | |
| points(1:4,tmp2[1:4],col=rgb(1,0,0,0.8),pch=19,xlim=c(0.5,8.5),ylab='',xlab='') | |
| par(new=T) | |
| plot(5:8,tmp2[5:8],col=rgb(1,0,0,0.3),type='l',lty=1,lwd=2,axes=F,xlim=c(0.5,8.5),ylab='',xlab='') | |
| par(new=T) | |
| points(5:8,tmp2[5:8],col=rgb(1,0,0,0.8),pch=19,xlim=c(0.5,8.5),ylab='',xlab='') | |
| abline(v=4.5,col='black') | |
| boxplot(t(as.matrix(ge2[biclustMembership[[tfs2[tf2]]],]))~sets1,names=names1,las=2,main=paste('Genes: ',tfNames2[tf2]),col=brewer.pal(8, "Paired")) | |
| tmp1 = t(sapply(1:max(sets1), function(x) { rowMeans(ge2[biclustMembership[[tfs2[tf2]]],which(sets1==x)]) } )) | |
| rownames(tmp1) = names1 | |
| matlines(1:4,tmp1[1:4,],type='l',col=rgb(0.3,0.3,0.3,0.5),lty=1,lwd=1) | |
| matlines(5:8,tmp1[5:8,],type='l',col=rgb(0.3,0.3,0.3,0.5),lty=1,lwd=1) | |
| par(new=T) | |
| plot(1:4,tmp2[1:4],col=rgb(1,0,0,0.3),type='l',lty=1,lwd=2,axes=F,xlim=c(0.5,8.5),ylab='',xlab='') | |
| par(new=T) | |
| points(1:4,tmp2[1:4],col=rgb(1,0,0,0.8),pch=19,xlim=c(0.5,8.5),ylab='',xlab='') | |
| par(new=T) | |
| plot(5:8,tmp2[5:8],col=rgb(1,0,0,0.3),type='l',lty=1,lwd=2,axes=F,xlim=c(0.5,8.5),ylab='',xlab='') | |
| par(new=T) | |
| points(5:8,tmp2[5:8],col=rgb(1,0,0,0.8),pch=19,xlim=c(0.5,8.5),ylab='',xlab='') | |
| abline(v=4.5,col='black') | |
| pheatmap(t(tmp1),cluster_cols=F,scale='row') | |
| } | |
| dev.off() | |
| # Select out clusters of | |
| library(GMD) | |
| tm1 = sapply(unique(sets), function(x) { rowMeans(dAll[,which(sets==x)]) }) | |
| hc1 = hclust(dist(dAll)) | |
| css1 = css.hclust(dist(dAll),hc1) | |
| #hc1 = hclust(dist(cor(t(d1[,names(col.rug)[f1]])))) | |
| #css1 = css.hclust(dist(cor(t(d1[,names(col.rug)[f1]]))),hc1) | |
| elbow1 = elbow.batch(css1) | |
| elbow1 = elbow.batch(css1,inc.thres=0.1,ev.thres=c(0.9,0.8)) | |
| cut1 = cutree(hc1, k=elbow1$k) | |
| #cut1 = cutree(hc1, k=5) | |
| #write.csv(cbind('Cluster Index'=cut1),'clusters.csv') | |
| # | |
| tm1 = sapply(unique(sets), function(x) { rowMeans(dAll[,which(sets==x)]) }) | |
| hc1 = hclust(dist(tm1)) | |
| css1 = css.hclust(dist(tm1),hc1) | |
| #hc1 = hclust(dist(cor(t(d1[,names(col.rug)[f1]])))) | |
| #css1 = css.hclust(dist(cor(t(d1[,names(col.rug)[f1]]))),hc1) | |
| elbow1 = elbow.batch(css1) | |
| elbow1 = elbow.batch(css1,inc.thres=0.1,ev.thres=c(0.9,0.8)) | |
| cut1 = cutree(hc1, k=elbow1$k) | |
| library(pheatmap) | |
| #ann1 = data.frame(Clusters=as.factor(cut2), Subtypes=Subtypes) | |
| #ann1 = data.frame(Subtypes=Subtypes) | |
| #rownames(ann1) = names(cut2) | |
| #Clusters = c('orange','red','black','cadetblue') | |
| #names(Clusters) = c(1,2,3,4) | |
| #SubColors = c('black','orange','cadetblue','red','green','white') | |
| #names(SubColors) = c('classical','mesenchymal','neural','proneural','control','Unknown') | |
| #ann_colors = list(Clusters=Clusters,Subtypes=SubColors) | |
| #ann_colors = list(Subtypes=SubColors) | |
| #pheatmap(d2[c('tfbs_db_230','pita_157','targetscan_530','hsa_miR_34a','hsa_miR_210','CEBPE','ETV6','FOSL1','IRF1','MAFF','NFKB1','TSC22D1'),names(cut2)[order(Subtypes)]],color=colorpanel(256,'blue','black','yellow'),annotation=ann1,show_rownames=T,show_colnames=F, annotation_colors=ann_colors,scale='row',cluster_cols=F) | |
| #pheatmap(dAll,color=colorpanel(256,'blue','black','yellow'),annotation=ann1,show_rownames=T,show_colnames=F, annotation_colors=ann_colors,scale='row',cluster_cols=F) | |
| pheatmap(dAll,color=colorpanel(256,'blue','black','yellow'),show_rownames=F,show_colnames=F,cluster_cols=F) | |
| tm1 = sapply(unique(sets), function(x) { rowMeans(dAll[,which(sets==x)]) }) | |
| pheatmap(tm1,color=colorpanel(256,'blue','black','yellow'),show_rownames=F,show_colnames=F,cluster_cols=F) | |
| m3 = t(sapply(1:length(unique(cut1)), function(x) {apply(dAll[names(which(cut1==x)),],2,mean)})) | |
| heatmap.2(as.matrix(m3),trace='none',col=colorpanel(256,'blue','black','yellow'),Colv=F,Rowv=T,dendrogram='row',density.info='none') #,ColSideColors=as.character(col.rug),scale='row') | |
| #dev.off() | |
| library(miscTools) | |
| tm3 = t(sapply(unique(cut1), function(x) { colMedians(tm1[names(which(cut1==x)),]) } )) | |
| heatmap.2(as.matrix(tm3),trace='none',col=colorpanel(256,'blue','black','yellow'),Colv=F,Rowv=T,dendrogram='row',density.info='none') #,ColSideColors=as.character(col.rug),scale='row') | |
| o1 = order(table(cut1)) | |
| #o1 = order(tm3[,3],decreasing=T) | |
| heatmap.2(as.matrix(m3)[o1,],trace='none',col=colorpanel(256,'blue','black','yellow'),Colv=F,Rowv=F,dendrogram='row',density.info='none') #,ColSideColors=as.character(col.rug),scale='row') | |
| # Resorting clusters | |
| cut2 = cut1 | |
| cut2[which(cut1==o1[1])] = 1 | |
| cut2[which(cut1==o1[2])] = 2 | |
| cut2[which(cut1==o1[3])] = 3 | |
| cut2[which(cut1==o1[4])] = 4 | |
| cut2[which(cut1==o1[5])] = 5 | |
| cut2[which(cut1==o1[6])] = 6 | |
| cut2[which(cut1==o1[7])] = 7 | |
| cut2[which(cut1==o1[8])] = 8 | |
| cut2[which(cut1==o1[9])] = 9 | |
| cut2[which(cut1==o1[10])] = 10 | |
| cut2[which(cut1==o1[11])] = 11 | |
| cut2[which(cut1==o1[12])] = 12 | |
| cut2[which(cut1==o1[13])] = 13 | |
| cut2[which(cut1==o1[14])] = 14 | |
| # Write out bicluster memberhips | |
| write.csv(cut2,'clustersOfBiclusters_VALL.csv') | |
| ######################################## | |
| # Heatmap of profiles for each cluster # | |
| ######################################## | |
| m3 = t(sapply(1:14, function(x) {apply(dAll[names(which(cut2==x)),],2,mean)})) | |
| pdf('tCellNetworkClustersOfBiclusters_trends_VALL.pdf') | |
| heatmap.2(as.matrix(m3),trace='none',col=colorpanel(256,'blue','black','yellow'),Colv=F,Rowv=F,dendrogram='none',RowSideColors=rainbow(length(unique(cut2))),density.info='none') #,scale='row') | |
| dev.off() | |
| ####################### | |
| ## Getting Gene Sets ## | |
| ####################### | |
| g1 = read.csv('genesInBiclusters/cluster.members.genes_tfbs_db.txt',header=F) | |
| biclustMembership = list() | |
| for(j in 1:length(g1[,1])) { | |
| biclustMembership[[paste('TFBS_DB_',j,sep='')]] = strsplit(as.character(g1[j,]),split=' ')[[1]][-1] | |
| } | |
| g2 = read.csv('genesInBiclusters/cluster.members.genes_pita.txt',header=F) | |
| for(j in 1:length(g2[,1])) { | |
| biclustMembership[[paste('PITA_',j,sep='')]] = strsplit(as.character(g2[j,]),split=' ')[[1]][-1] | |
| } | |
| g3 = read.csv('genesInBiclusters/cluster.members.genes_targetScan.txt',header=F) | |
| for(j in 1:length(g3[,1])) { | |
| biclustMembership[[paste('TargetScan_',j,sep='')]] = strsplit(as.character(g3[j,]),split=' ')[[1]][-1] | |
| } | |
| # Take out | |
| clusterMembership = list() | |
| for(i in 1:18) { | |
| tmp = c() | |
| for(j in names(cut2[which(cut2==i)])) { | |
| tmp = c(tmp,biclustMembership[[j]]) | |
| } | |
| clusterMembership[[i]] = unique(tmp) | |
| } | |
| # Write out clusterMembership | |
| dClustMemb = cbind(meta_cluster=1:length(clusterMembership),genes=unlist(lapply(clusterMembership,paste,collapse=';'))) | |
| write.csv(dClustMemb,'metaClusterMembership.csv',row.names=F) | |
| # Number of genes per cluster | |
| # 105 686 602 376 48 297 247 1118 470 929 | |
| sapply(1:length(clusterMembership), function(x) { length(clusterMembership[[x]]) }) | |
| # Functional enrichment | |
| library(topGO) | |
| # Read in GO mappings to affymetrix probe ids | |
| library('org.Mm.eg.db') | |
| # Make a map of Entrez IDs to UCSC transcript IDs | |
| d1 = read.csv('entrez2ucsc.csv',header=F) | |
| ucsc2entrez = list() | |
| for(j in 1:length(d1[,1])) { | |
| tmp = strsplit(as.character(d1[j,2]),split=';')[[1]] | |
| for(k in tmp) { | |
| ucsc2entrez[[as.character(k)]] = as.character(d1[j,1]) | |
| } | |
| } | |
| # Convert IDs to entrez | |
| biclustMembership.eg = list() | |
| for(j in 1:length(clusterMembership)) { | |
| biclustMembership.eg[[j]] = as.character(unlist(ucsc2entrez[clusterMembership[[j]]])) | |
| } | |
| biclustMembership.eg = clusterMembership | |
| # Setup for analysis | |
| xx <- annFUN.org("BP", mapping = "org.Mm.eg.db", ID = "entrez") | |
| #geneNames <- unique(unlist(xx)) | |
| geneNames <- intersect(unique(unlist(biclustMembership.eg)), unique(unlist(xx))) | |
| tmp1 <- biclustMembership.eg[[1]] | |
| geneList <- factor(as.integer(geneNames %in% tmp1)) | |
| names(geneList) <- geneNames | |
| # Make Biological Process GOData object | |
| GOdata.BP <- new("topGOdata", ontology='BP', allGenes = geneList, annot = annFUN.org, mapping = 'org.Mm.eg.db', ID = 'entrez') | |
| m1.BP <- matrix(nrow = length(GOdata.BP@graph@nodes), ncol = length(biclustMembership.eg), dimnames = list(GOdata.BP@graph@nodes, 1:length(biclustMembership.eg))) | |
| m2.BP <- matrix(nrow = length(biclustMembership.eg), ncol=2, dimnames = list(1:length(biclustMembership.eg),c('Top10.Terms.BP','BH.sig.GO.Ids.BP'))) | |
| # Make Molecular Function GOData object | |
| GOdata.MF <- new("topGOdata", ontology='MF', allGenes = geneList, annot = annFUN.org, mapping = 'org.Mm.eg.db') | |
| m1.MF <- matrix(nrow = length(GOdata.MF@graph@nodes), ncol = length(biclustMembership.eg), dimnames = list(GOdata.MF@graph@nodes, 1:length(biclustMembership.eg))) | |
| m2.MF <- matrix(nrow = length(biclustMembership.eg), ncol=2, dimnames = list(1:length(biclustMembership.eg),c('Top10.Terms.MF','BH.sig.GO.Ids.MF'))) | |
| # Make Cellular Component GOData object | |
| GOdata.CC <- new("topGOdata", ontology='CC', allGenes = geneList, annot = annFUN.org, mapping = 'org.Mm.eg.db') | |
| m1.CC <- matrix(nrow = length(GOdata.CC@graph@nodes), ncol = length(biclustMembership.eg), dimnames = list(GOdata.CC@graph@nodes, 1:length(biclustMembership.eg))) | |
| m2.CC <- matrix(nrow = length(biclustMembership.eg), ncol=2, dimnames = list(1:length(biclustMembership.eg),c('Top10.Terms.CC','BH.sig.GO.Ids.CC'))) | |
| for( clust in (1:length(biclustMembership.eg)) ) { | |
| # Expand gene list and change factor in GOdata | |
| clusterGenes <- biclustMembership.eg[[clust]] | |
| GOdata.BP@allScores <- factor(as.integer(geneNames %in% clusterGenes)) | |
| GOdata.MF@allScores <- factor(as.integer(geneNames %in% clusterGenes)) | |
| GOdata.CC@allScores <- factor(as.integer(geneNames %in% clusterGenes)) | |
| # Biological process | |
| r1.BP = runTest(GOdata.BP, algorithm = 'classic', statistic = 'fisher') | |
| m1.BP[,clust] = r1.BP@score | |
| m2.BP[clust,1] = gsub(',','_',paste(GenTable(GOdata.BP, r1.BP)[,2],collapse=';')) | |
| m2.BP[clust,2] = paste(names(which(p.adjust(r1.BP@score,method='BH')<=0.05)),collapse=';') | |
| #m2.BP[clust,2] = paste(names(which(r1.BP@score<=0.05)),collapse=';') | |
| r1.MF = runTest(GOdata.MF, algorithm = 'classic', statistic = 'fisher') | |
| m1.MF[,clust] = r1.MF@score | |
| m2.MF[clust,1] = paste(GenTable(GOdata.MF, r1.MF)[,2],collapse=':') | |
| m2.MF[clust,2] = paste(names(which(p.adjust(r1.MF@score,method='BH')<=0.05)),collapse=';') | |
| #m2.MF[clust,2] = paste(names(which(r1.MF@score<=0.05)),collapse=';') | |
| r1.CC = runTest(GOdata.CC, algorithm = 'classic', statistic = 'fisher') | |
| m1.CC[,clust] = r1.CC@score | |
| m2.CC[clust,1] = paste(GenTable(GOdata.CC, r1.CC)[,2],collapse=';') | |
| m2.CC[clust,2] = paste(names(which(p.adjust(r1.CC@score,method='BH')<=0.05)),collapse=';') | |
| #m2.CC[clust,2] = paste(names(which(r1.CC@score<=0.05)),collapse=';') | |
| } | |
| write.csv(m2.BP,'clusterEnrichment_GOBP_0.05.csv') | |
| write.csv(m2.MF,'clusterEnrichment_GOMF_0.05.csv') | |
| write.csv(m2.CC,'clusterEnrichment_GOCC_0.05.csv') | |
This file contains hidden or 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
| # Read in clusters of biclusters | |
| clusters = {} | |
| inFile = open('clustersOfBiclusters_VALL.csv','r') | |
| inFile.readline() # Get rid of header | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.strip().split(',') | |
| if not splitUp[1] in clusters: | |
| clusters[splitUp[1]] = [] | |
| clusters[splitUp[1]].append(splitUp[0].strip('"').lower()) | |
| inFile.close() | |
| # Read in post processed files | |
| #Bicluster - 0 Genes - 1 Conditions - 2 Norm. Residual - 3 Norm. Residual Perm. P-Value - 4 Var. Exp. First PC - 5 Var. Exp. First PC Perm. P-Value - 6 MEME Motif1 E-Value - 7 Up.MEME Motif1 Consensus - 8 Up.MEME Motif1 Matches - 9 Up.MEME Motif1 Expanded Matches - 10 Up.MEME Motif1 Correlated Matches - 11 Up.MEME Motif1 Original/Expanded - 12 Up.MEME Motif1 Minimum Correlated - 13 Up.MEME Motif2 E-Value - 14 Up.MEME Motif2 Consensus - 15 Up.MEME Motif2 Matches - 16 Up.MEME Motif2 Expanded Matches - 17 Up.MEME Motif2 Correlated Matches - 18 Up.MEME Motif2 Original/Expanded - 19 Up.MEME Motif2 Minimum Correlated - 20 Up.WEEDER Motif1 Score - 21 Up.WEEDER Motif1 Consensus - 22 Up.WEEDER Motif1 Matches - 23 Up.WEEDER Motif1 Expanded Matches - 24 Up.WEEDER Motif1 Correlated Matches - 25 Up.WEEDER Motif1 Original/Expanded - 26 Up.WEEDER Motif1 Minimum Correlated - 27 Up.WEEDER Motif2 Score - 28 Up.WEEDER Motif2 Consensus - 29 Up.WEEDER Motif2 Matches - 30 Up.WEEDR Motif2 Expanded Matches - 31 Up.WEEDER Motif2 Correlated Matches - 32 Up.WEEDER Motif2 Original/Expanded - 33 Up.WEEDER Motif2 Minimum Correlated - 34 TFBS_DB.TFs - 35 TFBS_DB.percTargets - 36 TFBS_DB.pValue - 37 TFBS_DB.Exapnded Matches - 38 TFBS_DB.Correlated Matches - 39 TFBS_DB.Original/Expanded - 40 TFBS_DB.Minimum Correlated - 41 3pUTR.WEEDER Motif1 E-Value - 42 3pUTR.WEEDER Motif1 Consensus - 43 3pUTR.WEEDER Motif1 Matches - 44 3pUTR.WEEDER Motif1 Model - 45 3pUTR.WEEDER Motif2 E-Value - 46 3pUTR.WEEDER Motif2 Consensus - 47 3pUTR.WEEDER Motif2 Matches - 48 3pUTR.WEEDER Motif2 Model - 49 3pUTR_pita.miRNAs - 50 3pUTR_pita.percTargets - 51 3pUTR_pita.pValue - 52 3pUTR_targetScan.miRNAs - 53 3pUTR_targetScan.percTargets - 54 3pUTR_targetScan.pValue - 55 GO_Term_BP - 56 | |
| biclusters = {} | |
| biclusterHeader = '' | |
| for i in ['pita','targetscan','tfbs_db']: | |
| inFile = open('postProcessed/postProcessed_urdahl_alvMac_'+i+'.csv','r') | |
| biclusterHeader = inFile.readline() # Get rid of header | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.split(',') | |
| biclusters[i+'_'+splitUp[0]] = line.strip() | |
| inFile.close() | |
| inFile = open('genesInBiclusters/cluster.members.genes_'+i+'.txt','r') | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.strip().split(' ') | |
| j = splitUp.pop(0) | |
| biclusters[i+'_'+j] = biclusters[i+'_'+j] + ',' + ';'.join(splitUp) | |
| inFile.close() | |
| inFile = open('genesInBiclusters/cluster.members.conditions_'+i+'.txt','r') | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.strip().split(' ') | |
| j = splitUp.pop(0) | |
| biclusters[i+'_'+j] = biclusters[i+'_'+j] + ',' + ';'.join(splitUp) | |
| inFile.close() | |
| # Dump combined postProcessed file with clusters IDs | |
| #Bicluster - 0 Genes - 1 Conditions - 2 Norm. Residual - 3 Norm. Residual Perm. P-Value - 4 Var. Exp. First PC - 5 Var. Exp. First PC Perm. P-Value - 6 MEME Motif1 E-Value - 7 Up.MEME Motif1 Consensus - 8 Up.MEME Motif1 Matches - 9 Up.MEME Motif1 Expanded Matches - 10 Up.MEME Motif1 Correlated Matches - 11 Up.MEME Motif1 Original/Expanded - 12 Up.MEME Motif1 Minimum Correlated - 13 Up.MEME Motif2 E-Value - 14 Up.MEME Motif2 Consensus - 15 Up.MEME Motif2 Matches - 16 Up.MEME Motif2 Expanded Matches - 17 Up.MEME Motif2 Correlated Matches - 18 Up.MEME Motif2 Original/Expanded - 19 Up.MEME Motif2 Minimum Correlated - 20 Up.WEEDER Motif1 Score - 21 Up.WEEDER Motif1 Consensus - 22 Up.WEEDER Motif1 Matches - 23 Up.WEEDER Motif1 Expanded Matches - 24 Up.WEEDER Motif1 Correlated Matches - 25 Up.WEEDER Motif1 Original/Expanded - 26 Up.WEEDER Motif1 Minimum Correlated - 27 Up.WEEDER Motif2 Score - 28 Up.WEEDER Motif2 Consensus - 29 Up.WEEDER Motif2 Matches - 30 Up.WEEDR Motif2 Expanded Matches - 31 Up.WEEDER Motif2 Correlated Matches - 32 Up.WEEDER Motif2 Original/Expanded - 33 Up.WEEDER Motif2 Minimum Correlated - 34 TFBS_DB.TFs - 35 TFBS_DB.percTargets - 36 TFBS_DB.pValue - 37 TFBS_DB.Exapnded Matches - 38 TFBS_DB.Correlated Matches - 39 TFBS_DB.Original/Expanded - 40 TFBS_DB.Minimum Correlated - 41 3pUTR.WEEDER Motif1 E-Value - 42 3pUTR.WEEDER Motif1 Consensus - 43 3pUTR.WEEDER Motif1 Matches - 44 3pUTR.WEEDER Motif1 Model - 45 3pUTR.WEEDER Motif2 E-Value - 46 3pUTR.WEEDER Motif2 Consensus - 47 3pUTR.WEEDER Motif2 Matches - 48 3pUTR.WEEDER Motif2 Model - 49 3pUTR_pita.miRNAs - 50 3pUTR_pita.percTargets - 51 3pUTR_pita.pValue - 52 3pUTR_targetScan.miRNAs - 53 3pUTR_targetScan.percTargets - 54 3pUTR_targetScan.pValue - 55 GO_Term_BP - 56 | |
| tfs = {} | |
| correlatedTfs = {} | |
| bestTfs = {} | |
| miRNAs = {} | |
| goBP = {} | |
| writeMe = ['Meta-Cluster,bicluster,'+biclusterHeader.strip()+',Gene Names,Condition Names'] | |
| for i in clusters: | |
| tfs[i] = {'MEME':[], 'WEEDER':[], 'TFBS_DB':[]} | |
| correlatedTfs[i] = {'MEME':[], 'WEEDER':[], 'TFBS_DB':[]} | |
| bestTfs[i] = {'MEME':[], 'WEEDER':[], 'TFBS_DB':[]} | |
| miRNAs[i] = [] | |
| goBP[i] = [] | |
| for j in clusters[i]: | |
| writeMe.append(i+','+j+','+biclusters[j]) | |
| # Process out tfs and miRNAs | |
| splitUp = dict(zip(biclusterHeader.strip().split(','),biclusters[j].split(','))) | |
| # MEME1 | |
| if (not splitUp['MEME Motif1 E-Value']=='NA') and float(splitUp['MEME Motif1 E-Value'])<=0.05 and (not splitUp['Up.MEME Motif1 Matches']=='NA'): # and len(splitUp[10].split(' '))<10: | |
| tfs[i]['MEME'] += splitUp['Up.MEME Motif1 Matches'].split(' ') | |
| if not splitUp['Up.MEME Motif1 Correlated Matches']=='NA': | |
| correlatedTfs[i]['MEME'] += [k.split(':')[0] for k in splitUp['Up.MEME Motif1 Correlated Matches'].split(' ')] | |
| bestTfs[i]['MEME'] += [splitUp['Up.MEME Motif1 Minimum Correlated'].split(':')[0]] | |
| # MEME2 | |
| if (not splitUp['Up.MEME Motif2 E-Value']=='NA') and float(splitUp['Up.MEME Motif2 E-Value'])<=0.05 and (not splitUp['Up.MEME Motif2 Matches']=='NA'): # and len(splitUp[14].split(' '))<10: | |
| tfs[i]['MEME'] += splitUp['Up.MEME Motif2 Matches'].split(' ') | |
| if not splitUp['Up.MEME Motif2 Correlated Matches']=='NA': | |
| correlatedTfs[i]['MEME'] += [k.split(':')[0] for k in splitUp['Up.MEME Motif2 Correlated Matches'].split(' ')] | |
| bestTfs[i]['MEME'] += [splitUp['Up.MEME Motif2 Minimum Correlated'].split(':')[0]] | |
| # Weeder | |
| if (not splitUp['Up.WEEDER Motif1 Matches']=='NA'): # and len(splitUp[17].split(' '))<10: | |
| tfs[i]['WEEDER'] += splitUp['Up.WEEDER Motif1 Matches'].split(' ') | |
| if not splitUp['Up.WEEDER Motif1 Correlated Matches']=='NA': | |
| correlatedTfs[i]['WEEDER'] += [k.split(':')[0] for k in splitUp['Up.WEEDER Motif1 Correlated Matches'].split(' ')] | |
| bestTfs[i]['WEEDER'] += [splitUp['Up.WEEDER Motif1 Minimum Correlated'].split(':')[0]] | |
| # Weeder | |
| if (not splitUp['Up.WEEDER Motif2 Matches']=='NA'): # and len(splitUp[20].split(' '))<10: | |
| tfs[i]['WEEDER'] += splitUp['Up.WEEDER Motif2 Matches'].split(' ') | |
| if not splitUp['Up.WEEDER Motif2 Correlated Matches']=='NA': | |
| correlatedTfs[i]['WEEDER'] += [k.split(':')[0] for k in splitUp['Up.WEEDER Motif2 Correlated Matches'].split(' ')] | |
| bestTfs[i]['WEEDER'] += [splitUp['Up.WEEDER Motif2 Minimum Correlated'].split(':')[0]] | |
| # TFBS_DB | |
| if (not splitUp['TFBS_DB.percTargets']=='NA') and float((splitUp['TFBS_DB.percTargets'].split(' '))[0])>=0.1 and float((splitUp['TFBS_DB.pValue'].split(' '))[0])<=0.05: | |
| tfs[i]['TFBS_DB'] += splitUp['TFBS_DB.TFs'].split(' ') | |
| if not splitUp['TFBS_DB.Correlated Matches']=='NA': | |
| correlatedTfs[i]['TFBS_DB'] += [k.split(':')[0] for k in splitUp['TFBS_DB.Correlated Matches'].split(' ')] | |
| bestTfs[i]['TFBS_DB'] += [splitUp['TFBS_DB.Minimum Correlated'].split(':')[0]] | |
| # WEEDER:miRvestigator miRNA? | |
| if splitUp['3pUTR.WEEDER Motif1 Model']=='8mer': | |
| miRNAs[i] += [splitUp['3pUTR.WEEDER Motif1 Matches']] | |
| if splitUp['3pUTR.WEEDER Motif2 Model']=='8mer': | |
| miRNAs[i] += [splitUp['3pUTR.WEEDER Motif2 Matches']] | |
| # PITA miRNA? | |
| if (not splitUp['3pUTR_pita.percTargets']=='NA') and float((splitUp['3pUTR_pita.percTargets'].split(' '))[0])>=0.1 and float((splitUp['3pUTR_pita.pValue'].split(' '))[0])<=0.05: # and len(splitUp[36].split(' '))<10: | |
| miRNAs[i] += splitUp['3pUTR_pita.miRNAs'].split(' ') | |
| # TargetScan miRNA? | |
| if (not splitUp['3pUTR_targetScan.percTargets']=='NA') and float((splitUp['3pUTR_targetScan.percTargets'].split(' '))[0])>=0.1 and float((splitUp['3pUTR_targetScan.pValue'].split(' '))[0])<=0.05: # and len(splitUp[39].split(' '))<10: | |
| miRNAs[i] += splitUp['3pUTR_targetScan.miRNAs'].split(' ') | |
| # GO Biological Processes | |
| if not splitUp['GO_Term_BP']=='NA': | |
| goBP[i] += splitUp['GO_Term_BP'].split(';') | |
| outFile = open('postProcessed_clustersOfBiclusters.csv','w') | |
| outFile.write('\n'.join(writeMe)) | |
| outFile.close() | |
| # Dump out TFs to test | |
| outFile = open('tf_regulators.csv','w') | |
| outFile.write('cluster,MEME,WEEDER,TFBS_DB,Number Predicted\n') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(tfs[i]['MEME'])))+','+';'.join(list(set(tfs[i]['WEEDER'])))+','+';'.join(list(set(tfs[i]['TFBS_DB'])))+','+str(len(list(set(list(set(tfs[i]['MEME']))+list(set(tfs[i]['WEEDER']))+list(set(tfs[i]['TFBS_DB'])))))) for i in tfs])) | |
| outFile.close() | |
| # Dump out TFs to test | |
| outFile = open('correlated_tf_regulators.csv','w') | |
| outFile.write('cluster,MEME,WEEDER,TFBS_DB,Number Predicted\n') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(correlatedTfs[i]['MEME'])))+','+';'.join(list(set(correlatedTfs[i]['WEEDER'])))+','+';'.join(list(set(correlatedTfs[i]['TFBS_DB'])))+','+str(len(list(set(list(set(correlatedTfs[i]['MEME']))+list(set(correlatedTfs[i]['WEEDER']))+list(set(correlatedTfs[i]['TFBS_DB'])))))) for i in correlatedTfs])) | |
| outFile.close() | |
| # Dump out TFs to test | |
| outFile = open('best_tf_regulators.csv','w') | |
| outFile.write('cluster,MEME,WEEDER,TFBS_DB,Number Predicted\n') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(bestTfs[i]['MEME'])))+','+';'.join(list(set(bestTfs[i]['WEEDER'])))+','+';'.join(list(set(bestTfs[i]['TFBS_DB'])))+','+str(len(list(set(list(set(bestTfs[i]['MEME']))+list(set(bestTfs[i]['WEEDER']))+list(set(bestTfs[i]['TFBS_DB'])))))) for i in bestTfs])) | |
| outFile.close() | |
| # Dump out miRNAs to test | |
| outFile = open('miRNA_regulators.csv','w') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(miRNAs[i])))+','+str(len(list(set(miRNAs[i])))) for i in miRNAs])) | |
| outFile.close() | |
| # Dump out GO BPs to test | |
| outFile = open('GO_BPs_metaClusters.csv','w') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(goBP[i])))+','+str(len(list(set(goBP[i])))) for i in goBP])) | |
| outFile.close() | |
| # Read in clusters of biclusters | |
| clusters = {} | |
| inFile = open('clustersOfBiclusters_VALL.csv','r') | |
| inFile.readline() # Get rid of header | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.strip().split(',') | |
| if not splitUp[1] in clusters: | |
| clusters[splitUp[1]] = [] | |
| clusters[splitUp[1]].append(splitUp[0].strip('"').lower()) | |
| inFile.close() | |
| # Read in post processed files | |
| #Bicluster - 0 Genes - 1 Conditions - 2 Norm. Residual - 3 Norm. Residual Perm. P-Value - 4 Var. Exp. First PC - 5 Var. Exp. First PC Perm. P-Value - 6 MEME Motif1 E-Value - 7 Up.MEME Motif1 Consensus - 8 Up.MEME Motif1 Matches - 9 Up.MEME Motif1 Expanded Matches - 10 Up.MEME Motif1 Correlated Matches - 11 Up.MEME Motif1 Original/Expanded - 12 Up.MEME Motif1 Minimum Correlated - 13 Up.MEME Motif2 E-Value - 14 Up.MEME Motif2 Consensus - 15 Up.MEME Motif2 Matches - 16 Up.MEME Motif2 Expanded Matches - 17 Up.MEME Motif2 Correlated Matches - 18 Up.MEME Motif2 Original/Expanded - 19 Up.MEME Motif2 Minimum Correlated - 20 Up.WEEDER Motif1 Score - 21 Up.WEEDER Motif1 Consensus - 22 Up.WEEDER Motif1 Matches - 23 Up.WEEDER Motif1 Expanded Matches - 24 Up.WEEDER Motif1 Correlated Matches - 25 Up.WEEDER Motif1 Original/Expanded - 26 Up.WEEDER Motif1 Minimum Correlated - 27 Up.WEEDER Motif2 Score - 28 Up.WEEDER Motif2 Consensus - 29 Up.WEEDER Motif2 Matches - 30 Up.WEEDR Motif2 Expanded Matches - 31 Up.WEEDER Motif2 Correlated Matches - 32 Up.WEEDER Motif2 Original/Expanded - 33 Up.WEEDER Motif2 Minimum Correlated - 34 TFBS_DB.TFs - 35 TFBS_DB.percTargets - 36 TFBS_DB.pValue - 37 TFBS_DB.Exapnded Matches - 38 TFBS_DB.Correlated Matches - 39 TFBS_DB.Original/Expanded - 40 TFBS_DB.Minimum Correlated - 41 3pUTR.WEEDER Motif1 E-Value - 42 3pUTR.WEEDER Motif1 Consensus - 43 3pUTR.WEEDER Motif1 Matches - 44 3pUTR.WEEDER Motif1 Model - 45 3pUTR.WEEDER Motif2 E-Value - 46 3pUTR.WEEDER Motif2 Consensus - 47 3pUTR.WEEDER Motif2 Matches - 48 3pUTR.WEEDER Motif2 Model - 49 3pUTR_pita.miRNAs - 50 3pUTR_pita.percTargets - 51 3pUTR_pita.pValue - 52 3pUTR_targetScan.miRNAs - 53 3pUTR_targetScan.percTargets - 54 3pUTR_targetScan.pValue - 55 GO_Term_BP - 56 | |
| biclusters = {} | |
| biclusterHeader = '' | |
| for i in ['pita','targetscan','tfbs_db']: | |
| inFile = open('postProcessed/postProcessed_urdahl_alvMac_'+i+'.csv','r') | |
| biclusterHeader = inFile.readline() # Get rid of header | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.split(',') | |
| biclusters[i+'_'+splitUp[0]] = line.strip() | |
| inFile.close() | |
| inFile = open('genesInBiclusters/cluster.members.genes_'+i+'.txt','r') | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.strip().split(' ') | |
| j = splitUp.pop(0) | |
| biclusters[i+'_'+j] = biclusters[i+'_'+j] + ',' + ';'.join(splitUp) | |
| inFile.close() | |
| inFile = open('genesInBiclusters/cluster.members.conditions_'+i+'.txt','r') | |
| while 1: | |
| line = inFile.readline() | |
| if not line: | |
| break | |
| splitUp = line.strip().split(' ') | |
| j = splitUp.pop(0) | |
| biclusters[i+'_'+j] = biclusters[i+'_'+j] + ',' + ';'.join(splitUp) | |
| inFile.close() | |
| # Dump combined postProcessed file with clusters IDs | |
| #Bicluster - 0 Genes - 1 Conditions - 2 Norm. Residual - 3 Norm. Residual Perm. P-Value - 4 Var. Exp. First PC - 5 Var. Exp. First PC Perm. P-Value - 6 MEME Motif1 E-Value - 7 Up.MEME Motif1 Consensus - 8 Up.MEME Motif1 Matches - 9 Up.MEME Motif1 Expanded Matches - 10 Up.MEME Motif1 Correlated Matches - 11 Up.MEME Motif1 Original/Expanded - 12 Up.MEME Motif1 Minimum Correlated - 13 Up.MEME Motif2 E-Value - 14 Up.MEME Motif2 Consensus - 15 Up.MEME Motif2 Matches - 16 Up.MEME Motif2 Expanded Matches - 17 Up.MEME Motif2 Correlated Matches - 18 Up.MEME Motif2 Original/Expanded - 19 Up.MEME Motif2 Minimum Correlated - 20 Up.WEEDER Motif1 Score - 21 Up.WEEDER Motif1 Consensus - 22 Up.WEEDER Motif1 Matches - 23 Up.WEEDER Motif1 Expanded Matches - 24 Up.WEEDER Motif1 Correlated Matches - 25 Up.WEEDER Motif1 Original/Expanded - 26 Up.WEEDER Motif1 Minimum Correlated - 27 Up.WEEDER Motif2 Score - 28 Up.WEEDER Motif2 Consensus - 29 Up.WEEDER Motif2 Matches - 30 Up.WEEDR Motif2 Expanded Matches - 31 Up.WEEDER Motif2 Correlated Matches - 32 Up.WEEDER Motif2 Original/Expanded - 33 Up.WEEDER Motif2 Minimum Correlated - 34 TFBS_DB.TFs - 35 TFBS_DB.percTargets - 36 TFBS_DB.pValue - 37 TFBS_DB.Exapnded Matches - 38 TFBS_DB.Correlated Matches - 39 TFBS_DB.Original/Expanded - 40 TFBS_DB.Minimum Correlated - 41 3pUTR.WEEDER Motif1 E-Value - 42 3pUTR.WEEDER Motif1 Consensus - 43 3pUTR.WEEDER Motif1 Matches - 44 3pUTR.WEEDER Motif1 Model - 45 3pUTR.WEEDER Motif2 E-Value - 46 3pUTR.WEEDER Motif2 Consensus - 47 3pUTR.WEEDER Motif2 Matches - 48 3pUTR.WEEDER Motif2 Model - 49 3pUTR_pita.miRNAs - 50 3pUTR_pita.percTargets - 51 3pUTR_pita.pValue - 52 3pUTR_targetScan.miRNAs - 53 3pUTR_targetScan.percTargets - 54 3pUTR_targetScan.pValue - 55 GO_Term_BP - 56 | |
| tfs = {} | |
| correlatedTfs = {} | |
| bestTfs = {} | |
| miRNAs = {} | |
| goBP = {} | |
| writeMe = ['Meta-Cluster,bicluster,'+biclusterHeader.strip()+',Gene Names,Condition Names'] | |
| for i in clusters: | |
| tfs[i] = {'MEME':[], 'WEEDER':[], 'TFBS_DB':[]} | |
| correlatedTfs[i] = {'MEME':[], 'WEEDER':[], 'TFBS_DB':[]} | |
| bestTfs[i] = {'MEME':[], 'WEEDER':[], 'TFBS_DB':[]} | |
| miRNAs[i] = [] | |
| goBP[i] = [] | |
| for j in clusters[i]: | |
| writeMe.append(i+','+j+','+biclusters[j]) | |
| # Process out tfs and miRNAs | |
| splitUp = dict(zip(biclusterHeader.strip().split(','),biclusters[j].split(','))) | |
| # MEME1 | |
| if (not splitUp['MEME Motif1 E-Value']=='NA') and float(splitUp['MEME Motif1 E-Value'])<=0.05 and (not splitUp['Up.MEME Motif1 Matches']=='NA'): # and len(splitUp[10].split(' '))<10: | |
| tfs[i]['MEME'] += splitUp['Up.MEME Motif1 Matches'].split(' ') | |
| if not splitUp['Up.MEME Motif1 Correlated Matches']=='NA': | |
| correlatedTfs[i]['MEME'] += [k.split(':')[0] for k in splitUp['Up.MEME Motif1 Correlated Matches'].split(' ')] | |
| bestTfs[i]['MEME'] += [splitUp['Up.MEME Motif1 Minimum Correlated'].split(':')[0]] | |
| # MEME2 | |
| if (not splitUp['Up.MEME Motif2 E-Value']=='NA') and float(splitUp['Up.MEME Motif2 E-Value'])<=0.05 and (not splitUp['Up.MEME Motif2 Matches']=='NA'): # and len(splitUp[14].split(' '))<10: | |
| tfs[i]['MEME'] += splitUp['Up.MEME Motif2 Matches'].split(' ') | |
| if not splitUp['Up.MEME Motif2 Correlated Matches']=='NA': | |
| correlatedTfs[i]['MEME'] += [k.split(':')[0] for k in splitUp['Up.MEME Motif2 Correlated Matches'].split(' ')] | |
| bestTfs[i]['MEME'] += [splitUp['Up.MEME Motif2 Minimum Correlated'].split(':')[0]] | |
| # Weeder | |
| if (not splitUp['Up.WEEDER Motif1 Matches']=='NA'): # and len(splitUp[17].split(' '))<10: | |
| tfs[i]['WEEDER'] += splitUp['Up.WEEDER Motif1 Matches'].split(' ') | |
| if not splitUp['Up.WEEDER Motif1 Correlated Matches']=='NA': | |
| correlatedTfs[i]['WEEDER'] += [k.split(':')[0] for k in splitUp['Up.WEEDER Motif1 Correlated Matches'].split(' ')] | |
| bestTfs[i]['WEEDER'] += [splitUp['Up.WEEDER Motif1 Minimum Correlated'].split(':')[0]] | |
| # Weeder | |
| if (not splitUp['Up.WEEDER Motif2 Matches']=='NA'): # and len(splitUp[20].split(' '))<10: | |
| tfs[i]['WEEDER'] += splitUp['Up.WEEDER Motif2 Matches'].split(' ') | |
| if not splitUp['Up.WEEDER Motif2 Correlated Matches']=='NA': | |
| correlatedTfs[i]['WEEDER'] += [k.split(':')[0] for k in splitUp['Up.WEEDER Motif2 Correlated Matches'].split(' ')] | |
| bestTfs[i]['WEEDER'] += [splitUp['Up.WEEDER Motif2 Minimum Correlated'].split(':')[0]] | |
| # TFBS_DB | |
| if (not splitUp['TFBS_DB.percTargets']=='NA') and float((splitUp['TFBS_DB.percTargets'].split(' '))[0])>=0.1 and float((splitUp['TFBS_DB.pValue'].split(' '))[0])<=0.05: | |
| tfs[i]['TFBS_DB'] += splitUp['TFBS_DB.TFs'].split(' ') | |
| if not splitUp['TFBS_DB.Correlated Matches']=='NA': | |
| correlatedTfs[i]['TFBS_DB'] += [k.split(':')[0] for k in splitUp['TFBS_DB.Correlated Matches'].split(' ')] | |
| bestTfs[i]['TFBS_DB'] += [splitUp['TFBS_DB.Minimum Correlated'].split(':')[0]] | |
| # WEEDER:miRvestigator miRNA? | |
| if splitUp['3pUTR.WEEDER Motif1 Model']=='8mer': | |
| miRNAs[i] += [splitUp['3pUTR.WEEDER Motif1 Matches']] | |
| if splitUp['3pUTR.WEEDER Motif2 Model']=='8mer': | |
| miRNAs[i] += [splitUp['3pUTR.WEEDER Motif2 Matches']] | |
| # PITA miRNA? | |
| if (not splitUp['3pUTR_pita.percTargets']=='NA') and float((splitUp['3pUTR_pita.percTargets'].split(' '))[0])>=0.1 and float((splitUp['3pUTR_pita.pValue'].split(' '))[0])<=0.05: # and len(splitUp[36].split(' '))<10: | |
| miRNAs[i] += splitUp['3pUTR_pita.miRNAs'].split(' ') | |
| # TargetScan miRNA? | |
| if (not splitUp['3pUTR_targetScan.percTargets']=='NA') and float((splitUp['3pUTR_targetScan.percTargets'].split(' '))[0])>=0.1 and float((splitUp['3pUTR_targetScan.pValue'].split(' '))[0])<=0.05: # and len(splitUp[39].split(' '))<10: | |
| miRNAs[i] += splitUp['3pUTR_targetScan.miRNAs'].split(' ') | |
| # GO Biological Processes | |
| if not splitUp['GO_Term_BP']=='NA': | |
| goBP[i] += splitUp['GO_Term_BP'].split(';') | |
| outFile = open('postProcessed_clustersOfBiclusters.csv','w') | |
| outFile.write('\n'.join(writeMe)) | |
| outFile.close() | |
| # Dump out TFs to test | |
| outFile = open('tf_regulators.csv','w') | |
| outFile.write('cluster,MEME,WEEDER,TFBS_DB,Number Predicted\n') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(tfs[i]['MEME'])))+','+';'.join(list(set(tfs[i]['WEEDER'])))+','+';'.join(list(set(tfs[i]['TFBS_DB'])))+','+str(len(list(set(list(set(tfs[i]['MEME']))+list(set(tfs[i]['WEEDER']))+list(set(tfs[i]['TFBS_DB'])))))) for i in tfs])) | |
| outFile.close() | |
| # Dump out TFs to test | |
| outFile = open('correlated_tf_regulators.csv','w') | |
| outFile.write('cluster,MEME,WEEDER,TFBS_DB,Number Predicted\n') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(correlatedTfs[i]['MEME'])))+','+';'.join(list(set(correlatedTfs[i]['WEEDER'])))+','+';'.join(list(set(correlatedTfs[i]['TFBS_DB'])))+','+str(len(list(set(list(set(correlatedTfs[i]['MEME']))+list(set(correlatedTfs[i]['WEEDER']))+list(set(correlatedTfs[i]['TFBS_DB'])))))) for i in correlatedTfs])) | |
| outFile.close() | |
| # Dump out TFs to test | |
| outFile = open('best_tf_regulators.csv','w') | |
| outFile.write('cluster,MEME,WEEDER,TFBS_DB,Number Predicted\n') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(bestTfs[i]['MEME'])))+','+';'.join(list(set(bestTfs[i]['WEEDER'])))+','+';'.join(list(set(bestTfs[i]['TFBS_DB'])))+','+str(len(list(set(list(set(bestTfs[i]['MEME']))+list(set(bestTfs[i]['WEEDER']))+list(set(bestTfs[i]['TFBS_DB'])))))) for i in bestTfs])) | |
| outFile.close() | |
| # Dump out miRNAs to test | |
| outFile = open('miRNA_regulators.csv','w') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(miRNAs[i])))+','+str(len(list(set(miRNAs[i])))) for i in miRNAs])) | |
| outFile.close() | |
| # Dump out GO BPs to test | |
| outFile = open('GO_BPs_metaClusters.csv','w') | |
| outFile.write('\n'.join([i+','+';'.join(list(set(goBP[i])))+','+str(len(list(set(goBP[i])))) for i in goBP])) | |
| outFile.close() | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment