||
#!/usr/bin/Rscript
#calculate how many down and up -regulated genes are hyper and hypo at promoter
#dealing with down-regulated genes
load("/home/zhanghl/workshop/paper1/methylation/anno/fullannotInd.rda")
# fullannot # colnames rownames
# TSS1500Ind # TSS1500Ind$SID names(TSS1500Ind$SID) TSS1500Ind$SID[[1]] TSS1500Ind$SID[["ZYX"]]
com_down_genes <- scan("/home/zhanghl/workshop/paper1/expression/analysis_common_up_down_genes/survival_gene_list_down",what="")
tumors <- c("BRCA","KIRC","PRAD","THCA","UCEC")
#readin methy data
for (tumor in tumors){
spc_file <- paste("/home/zhanghl/workshop/paper1/expression/analysis_common_up_down_genes/",tumor,".sig.down.genes",sep="")
spc_genes <- scan(spc_file,what="")
com_genes <- intersect(spc_genes,com_down_genes)
com_genes_symbol <- sub("\|[a-zA-Z0-9?]*","",com_genes,fixed=FALSE)
inter_genes_anno <- intersect(names(TSS1500Ind$SID),com_genes_symbol)
methy_file <- paste("/home/zhanghl/workshop/paper1/methylation/methylaiton_data_v3/",tumor,".paired.methylation2",sep="")
methy <- read.table(methy_file,header=TRUE,row.names=1,sep="t",stringsAsFactors=FALSE)
methy_patients <- unique(substr(colnames(methy),1,12))
t_names <- paste(methy_patients,".01",sep="")
n_names <- paste(methy_patients,".11",sep="")
sf <- matrix(ncol=3,nrow=length(inter_genes_anno))
sf <- as.data.frame(sf)
colnames(sf) <- c("normal","Tumor","hypo_hyper")
rownames(sf) <- inter_genes_anno
for (i in inter_genes_anno){
probes <- TSS1500Ind$SID[[i]]
methy_ele_t_media <- apply(methy[probes,t_names],2,median,na.rm=TRUE)
methy_ele_n_median <- apply(methy[probes,n_names],2,median,na.rm=TRUE)
sf[i,1] <- median(methy_ele_n_median,na.rm=TRUE)
sf[i,2] <- median(methy_ele_t_media,na.rm=TRUE)
test <- wilcox.test(methy_ele_t_media,methy_ele_n_median,paired=TRUE,exact = NULL, correct = TRUE)
sig <- test$p.value
sf[i,3] <- ifelse(sig>0.05,"NULL",ifelse(sf[i,1]>sf[i,2],"hypo","hyper"))
}
write.table(sf,paste(tumor,".down.genes",sep=""),row.names=TRUE,col.names=TRUE,sep="t",quote=FALSE)
hyper_num <- sum(sf[,3]=="hyper")
hypo_num <- sum(sf[,3]=="hypo")
null_num <- sum(sf[,3]=="NULL")
print (paste("hyper_num:",hyper_num,"hypo_num:",hypo_num,"null_null:",null_num,sep=" "))
}
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2024-10-6 16:28
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社