TCGA/GTEx泛癌数据任意基因表达量展示

news/2024/11/28 19:25:38/

有了泛癌的数据之后就可以进行各种分析了,当然这些都是在R语言的基础上进行的。如果你不会R语言,也可以通过各种各样的网页工具实现。

我们今天就简单展示下任意基因在泛癌图谱中的表达量情况。

TCGA,GTEx,TCGA+GTEx的泛癌数据都整理好了,大家可以自己通过easyTCGA包实现1行代码整理,也可以直接在公众号后台回复pancancer获取整理好的数据。详情请见:任意基因在泛癌中的表达量展示

GTEx

GTEx的展示比较简单,最常见的就是某个基因在所有组织中的表达量情况。

# 加载数据
load(file="output_pancancer_xena/GTEx_pancancer_mrna_pheno.rdata")
library(tidyverse)
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.4.1     ✔ purrr   1.0.1
## ✔ tibble  3.2.1     ✔ dplyr   1.1.1
## ✔ tidyr   1.3.0     ✔ stringr 1.5.0
## ✔ readr   2.1.4     ✔ forcats 1.0.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()# 简单看下,这几个泛癌数据的详细情况我都给大家有说明,看一下即可
head(colnames(gtex_mrna_pheno))
## [1] "sample_id"    "primary_site" "MT-ATP8"      "MT-ATP6"      "MT-CO2"      
## [6] "MT-CO3"
#table(gtex_mrna_pheno$primary_site)
length(table(gtex_mrna_pheno$primary_site))
## [1] 31

接下来以CXCL1这个基因为例进行展示。

gene <- "CXCL1"# 提取数据就是这么简单
plot_df <- gtex_mrna_pheno %>%select(1:2,all_of(gene))# 画图即可
ggplot2::ggplot(plot_df, aes(fct_reorder(primary_site,CXCL1),CXCL1))+ggplot2::geom_boxplot(aes(fill=primary_site))+ggplot2::labs(x=NULL)+ggplot2::theme_bw()+ggplot2::theme(legend.position = "none",axis.text.x = ggplot2::element_text(angle = 45,hjust = 1))

TCGA

单独使用TCGA泛癌的数据进行展示是花样最多的,你在pubmed中以pan cancer为关键词进行检索,基本上其中的Fig1都是类似的箱线图。

# tcga pancancer,前34列是临床信息
rm(list = ls())
load(file="output_pancancer_xena/TCGA_pancancer_mrna_clin.rdata")head(colnames(tcga_mrna_clin))
## [1] "sample_id"                           "patient_id"                         
## [3] "project"                             "age_at_initial_pathologic_diagnosis"
## [5] "gender"                              "race"

继续以CXCL1这个基因为例进行展示。

gene <- "CXCL1"plot_df <- tcga_mrna_clin %>%select(sample_id,project,all_of(gene)) %>%# 这个分组你可以任意指定,并不一定要tumor、normalmutate(sample_type=ifelse(as.numeric(substr(.$sample_id,14,15))<10,"tumor","normal"))#tcga pancancer中有很多癌种没有normal哦,要注意!
ggplot2::ggplot(plot_df,aes(project,CXCL1))+ggplot2::geom_boxplot(aes(fill=sample_type))+ggplot2::labs(x=NULL,y="expression")+ggplot2::theme_bw()+ggplot2::theme(legend.position = "top",axis.text.x = ggplot2::element_text(angle = 45,hjust = 1))+ggpubr::stat_compare_means(ggplot2::aes(group = sample_type,label = "p.format"),method = "kruskal.test")

或者你也可以展示只在tumor样本中的表达量。

plot_df <- plot_df %>%filter(sample_type=="tumor")ggplot2::ggplot(plot_df,aes(fct_reorder(project,CXCL1),CXCL1))+ggplot2::geom_boxplot(aes(fill=project))+ggplot2::labs(x=NULL,y="expression")+ggplot2::theme_bw()+ggplot2::theme(legend.position = "none",axis.text.x = ggplot2::element_text(angle = 45,hjust = 1))

接下来是大家比较感兴趣的某个基因在泛癌配对样本中的表达。

首先我们把泛癌的表达矩阵(这里应该叫转置后的表达矩阵比较合适,一般我我们说表达矩阵就是指行是基因,列是样本的矩阵)按照project拆分,然后自定义一个可以提取配对样本的函数:

# 拆分
cancer_list <- split(tcga_mrna_clin,tcga_mrna_clin$project)# 自定义函数
get_paired_sample <- function(exprset){# get paired samplessample_group <- ifelse(as.numeric(substr(exprset$sample_id,14,15))<10,"tumor","normal")tmp <- data.frame(sample_group = sample_group, sample_id=exprset$sample_id,project=exprset$project)tmp_nor <- tmp[tmp$sample_group=="normal",]tmp_tum <- tmp[tmp$sample_group=="tumor",]#每一个normal都有配对的tumor吗?并不是keep <- intersect(substr(tmp_tum$sample_id,1,12),substr(tmp_nor$sample_id,1,12))tmp_tum <- tmp_tum[substr(tmp_tum$sample_id,1,12) %in% keep,]tmp_tum <- tmp_tum[!duplicated(substr(tmp_tum$sample_id,1,12)),]tmp_nor <- tmp_nor[substr(tmp_nor$sample_id,1,12) %in% keep,]tmp_nor <- tmp_nor[!duplicated(substr(tmp_nor$sample_id,1,12)),]tmp_pair <- rbind(tmp_tum,tmp_nor)
}

接下来就是把这个函数应用于33种癌症中,然后提取CXCL1这个基因的画图数据即可:

paired_samples <- do.call(rbind,lapply(cancer_list,get_paired_sample))plot_df <- paired_samples %>%left_join(tcga_mrna_clin[,c(gene,"sample_id")]) %>%mutate(sample_id=substr(sample_id,1,12))
## Joining with `by = join_by(sample_id)`

接下来画图就是基本功了,ggplot2搞定一切,下面这个interaction的用法在《R数据可视化手册》中有讲过,我强烈呼吁大家赶紧买本书看看吧!别再天天问图怎么画了!

ggplot(plot_df, aes(interaction(sample_group,project),CXCL1,color=sample_group))+ggplot2::geom_point(size=3,position = position_dodge(0.9))+ggplot2::geom_line(aes(group=interaction(sample_id,project)),color="grey70")+ggplot2::scale_color_manual(values = c("#028EA1","#F2AA9D"))+ggplot2::scale_x_discrete(labels = rep(unique(plot_df$project),each=2))+ggplot2::theme_bw()+ggplot2::theme(legend.position = "top",axis.text.x = element_text(angle = 45,hjust = 1))

当然还有分面的画法:

#下面是分面
ggplot(plot_df, aes(interaction(sample_group,project),CXCL1,color=sample_group))+ggplot2::geom_point(size=3)+ggplot2::geom_line(aes(group=interaction(sample_id,project)),color="grey70")+ggplot2::scale_color_manual(values = c("#028EA1","#F2AA9D"))+ggplot2::scale_x_discrete(name = NULL)+ggplot2::facet_grid(~project,scales="free_x",switch = "x")+ggplot2::theme_bw()+ggplot2::theme(legend.position = "top",axis.text.x = element_blank(),strip.background = element_blank(),axis.ticks.x = element_blank(),panel.border = element_blank())

当然如果你看了书也搞不明白,也可以通过万能的网络解决一切,比如上面这种图,你可以通过关键词搜索:ggplot2 paired line multi groups, 实现方式非常多,任你选择。

TCGA+GTEx

TCGA+GTEx就没有配对展示了,除此之外都和TCGA的泛癌展示方式差不多。

rm(list = ls())
load(file="output_pancancer_xena/TCGA_GTEx_pancancer_mRNA_pheno.rdata")# 前4列是样本信息
head(colnames(tcga_gtex_mrna_pheno))
## [1] "sample_id"    "sample_type"  "project"      "primary_site" "MT-ATP6"     
## [6] "MT-CO2"
table(tcga_gtex_mrna_pheno$sample_type)
## 
## GTEx_normal TCGA_normal  TCGA_tumor 
##        7568         712        9784
table(tcga_gtex_mrna_pheno$project)
## 
##  ACC BLCA BRCA CESC CHOL COAD DLBC ESCA  GBM HNSC KICH KIRC KIRP LAML  LGG LIHC 
##  205  435 1390  319   45  637  491  848 1317  564  119  631  349  243 1674  531 
## LUAD LUSC MESO   OV PAAD PCPG PRAD READ SARC SKCM STAD TGCT THCA THYM UCEC  UCS 
##  862  836   87  515  350  185  648  410  264 1282  624  302  850  565  272  135 
##  UVM 
##   79

继续以CXCL1这个基因为例进行展示。

用的最多的肯定还是任意基因在不同组别中的表达:

gene <- "CXCL1"plot_df <- tcga_gtex_mrna_pheno %>%select(1:4,all_of(gene)) %>%filter(sample_type %in% c("GTEx_normal","TCGA_tumor"))ggplot(plot_df,aes(project,CXCL1))+geom_boxplot(aes(fill=sample_type))+theme(legend.position = "top")+ggplot2::labs(x=NULL,y="expression")+ggplot2::theme_bw()+ggplot2::theme(legend.position = "top",axis.text.x = ggplot2::element_text(angle = 45,hjust = 1))

扩展

其实任何类似于这个数据的格式都能像这样展示。

比如你可以通过ssGSEA对泛癌进行免疫浸润分析,这样每个样本都可以有一个得分,这样你就可以展示某个细胞在不同组别中的得分情况。

大家一定要多看文献,多积累不同的方法,以及一些好用的网站、图表等,说不定以后就用到了!

后面可能会安排几篇图表复现的推文,敬请期待。


http://www.ppmy.cn/news/966964.html

相关文章

GATK变异检测

变异检测 GATK 变异检测 GATK是Genome Analysis Toolkit的缩写&#xff0c;是用来处理高通量测序数据的一套软件。最初&#xff0c;GATK被设计用来分析人类基因组和外显子&#xff0c;主要用来寻找SNP和indel。后开&#xff0c;GATK的功能越来越丰富&#xff0c;增加了short v…

四、肿瘤全基因组学体细胞点突变特征(The repertoire of mutational signatures in human cancer)

全文链接 一、肿瘤突变特征&#xff1a;碱基置换及插入、缺失突变 单碱基置换&#xff08;49种特征类型&#xff0c;single-base-substitution&#xff0c;SBS&#xff09; 双碱基置换&#xff08;11种特征类型&#xff0c;doublet-base-substitution&#xff0c;DBS&#xf…

TCGA_联合GTEx分析2_查看批次效应

在 TCGA_联合GTEx分析1_得到表达矩阵.tpm_老实人谢耳朵的博客-CSDN博客 中&#xff0c;获取了TCGA和GTEx中样本的表达矩阵数据&#xff0c;数据格式均为tpm。本文对二者进行合并后&#xff0c;通过PCA分析、绘制内参箱线图等方法&#xff0c;查看是否存在批次效应。 关于批次效…

TCGA 亚型突变负荷代码

#1、准备文件/数据并加载相应的包 #1.1下载并加载相应的包&#xff0c;有就直接加载&#xff0c;没有就下载后再加载。 install.packages("pacman") library(pacman) p_load(TCGAbiolinks,DT,tidyverse) BiocManager::install("TCGAbiolinks") library(t…

基因组变异检测SNPcalling(GATK)

基因组变异检测SNPcalling&#xff08;GATK&#xff09; 第一步&#xff0c;将数据软连接到自己的工作路径下第二步&#xff0c;用BWA index功能为ref文件做index&#xff0c;用Samtools软件为ref做index第三步&#xff0c;用BWA软件做序列比对&#xff0c;得到序列比对的sam文…

提取TCGA 中体细胞突变数据的表达矩阵

#因为之前的命令调用GDCquery_Maf 发现用不了 #故找到了一些其他的方法&#xff0c;并且自己试着将其弄成了一个表达矩阵。 #代码如下 #1、下载加载相应的包 install.packages("pacman") library(pacman) p_load(TCGAbiolinks,DT,tidyverse) BiocManager::insta…

maftools|TCGA肿瘤突变数据的汇总,分析和可视化

之前介绍了使用maftools | 从头开始绘制发表级oncoplot&#xff08;瀑布图&#xff09; R-maftools包绘制组学突变结果&#xff08;MAF&#xff09;的oncoplot或者叫“瀑布图”&#xff0c;以及一些细节的更改和注释。 本文继续介绍maftools对于MAF文件的其他应用&#xff0c;为…

新版TCGA的突变SNP数据添加临床信息

文章目录 加载数据和R包读取数据 今天给大家演示下如何用自己的数据完成maftools的分析&#xff0c;主要是snp文件和临床信息的制作&#xff0c;其实很简单&#xff0c;但是网络上的教程都说的不清楚。 这次我们直接用之前TCGA-COAD和TCGA-READ合并后的数据演示&#xff0c;合…