生信分析工具包生信基础知识生信

GTEx: 如何找组织特异表达的基因

2019-12-25  本文已影响0人  TOP生物信息

官网--https://gtexportal.org/home/
数据下载链接--https://gtexportal.org/home/datasets

下文方法参考文献:Tissue-specific transcription reprogramming promotes liver metastasis of colorectal cancer

For each gene, we ranked the median expression value for each tissue or cell type in decreasing order. Genes defined as tissue-specific needed to
meet two criteria:

  1. gene expression ranked in the top 5 among all tissue and cell types
  2. also highly expressed (>90th percentile of all genes) in particular tissues and cell types
tmp <- read.table("gene_median_tpm.txt",header = T,sep = "\t",stringsAsFactors = F)
Name_Description <- tmp[,c("Name","Description")]
rownames(tmp) <- tmp$Name
gene_median_tpm <- tmp[,c(-1,-2)]
tissue_gene <- as.data.frame(t(gene_median_tpm))
# #可以用小数据集测试
# small_tissue_gene <- tissue_gene[,1:100]

#判断矩阵
judge_matrix <- matrix(nrow = length(colnames(gene_median_tpm)),ncol = length(rownames(gene_median_tpm)))

#先求分位数,后面直接调用,一定程度上加快速度
vec_quan <- c()
for (tissue_num in seq(1,length(rownames(tissue_gene)),1)) {
  quan=quantile(tissue_gene[tissue_num,],0.9)
  vec_quan <- append(vec_quan,as.numeric(quan))
}

#两层循环遍历,并判断
for (col_i in seq(1,length(colnames(tissue_gene)),1)) {
 is_top5=rank(tissue_gene[,col_i])>=length(rownames(tissue_gene))-(5-1)
  for (row_j in seq(1,length(colnames(gene_median_tpm)),1)) {
    if (is_top5[row_j]) {
      quan=vec_quan[row_j]
      if (tissue_gene[row_j,col_i] > quan) {
        judge_matrix[row_j,col_i]=1
      } else {
        judge_matrix[row_j,col_i]=0
      }
    } else {
      judge_matrix[row_j,col_i]=0
    }
  }
}

#将判断矩阵转换为数据框,并修改行名、列名
judge_matrix <- as.data.frame(judge_matrix)
colnames(judge_matrix) <- colnames(tissue_gene)
rownames(judge_matrix) <- rownames(tissue_gene)
#转置回来
gene_tissue_specific <- as.data.frame(t(judge_matrix))
#下面这个差不多就是一张总表
gene_tissue_specific <- gene_tissue_specific[!(rowSums(gene_tissue_specific)==0),]

比如找肝特异表达的基因

Liver <- as.data.frame(gene_tissue_specific$Liver)
colnames(Liver) <- c("specific")
rownames(Liver) <- rownames(gene_tissue_specific)
Liver$Name <- rownames(Liver)
Liver <- filter(Liver,Liver$specific == 1)
Liver <- inner_join(Liver,Name_Description,by="Name")
Liver <- Liver[,-1]

乙型结肠特异表达的基因

Colon_Sigmoid <- as.data.frame(gene_tissue_specific$Colon...Sigmoid)
colnames(Colon_Sigmoid) <- c("specific")
rownames(Colon_Sigmoid) <- rownames(gene_tissue_specific)
Colon_Sigmoid$Name <- rownames(Colon_Sigmoid)
Colon_Sigmoid <- filter(Colon_Sigmoid,Colon_Sigmoid$specific == 1)
Colon_Sigmoid <- inner_join(Colon_Sigmoid,Name_Description,by="Name")
Colon_Sigmoid <- Colon_Sigmoid[,-1]

韦恩图比较一下

library(VennDiagram)
venn.diagram(
  x = list(Liver$Description,Colon_Sigmoid$Description),
  category.names = c("Liver" , "Colon_Sigmoid"),
  filename = "Liver_vs_Colon_Sigmoid.png",
  output=T
)

和预期一样,交集很少,这些表示两个部位相对于其他部位是特异表达的

上一篇下一篇

猜你喜欢

热点阅读