#本教程仿自于B站的15天入门生物信息视频,若有疑惑请以视频为准。若侵请联。
#上一个教程我们得到了genes.counts.matrix.OE_vs_WT.DESeq2.DE_results文件,接下来对差异基因绘制火山图,我们这一步在Windows系统上的R进行即可,若需要安装Windows的R请看教程“https://www.cnblogs.com/liangjinghui/p/17855049.html”。
#1,我们把之前算表达量的一个文件genes.TMM.EXPR.matrix给弄到工作文件夹中,并将第一列弄成这样,就是在前面第一列标注为“ID”,genes.counts.matrix.OE_vs_WT.DESeq2.DE_results也是:
#2,制作group.txt,tab分隔:
Sample group
WT1 WT
WT2 WT
WT3 WT
OE1 OE
OE2 OE
OE3 OE
#3,进入R,设置工作文件夹,加载数据,并做数据关联
library(readr) gene_exp <- read_delim("有参转录组实战6-差异基因绘制火山图\\genes.TMM.EXPR.matrix") sample_info <- read_delim("有参转录组实战6-差异基因绘制火山图\\group.txt") de_result <- read_delim("有参转录组实战6-差异基因绘制火山图\\genes.counts.matrix.OE_vs_WT.DESeq2.DE_results") library (dplyr) library(tidyverse) de_result <- dplyr::select(de_result, ID, log2FoldChange, pvalue, padj)%>% mutate(direction=if_else(abs(log2FoldChange)<1|padj>0.05,'ns', if_else(log2FoldChange>=1,'up','down')))%>% left_join(gene_exp,by=c('ID'='ID'))
#4,计算上调下调基因的数量
group_by(de_result, direction) %>% summarise(count=n())
#5,绘制火山图
library(EnhancedVolcano)
EnhancedVolcano(de_result, lab =de_result$ID, selectLab = F, x = 'log2FoldChange', y = 'padj',
title = 'Volcano Plot', subtitle= 'OE vs WT', xlim = c(-10, 10), pCutoff = 0.05, FCcutoff = 1)
#6,火山图没什么用,最多放附件,自己还可以试着用ggplot2画。
#蜂蜜与四叶草
网友评论