如何进行差异表达基因分析-基于DESeq2包
差异表达分析是做RNA-seq分析必经的一步,这里给大家介绍一下。# 1 DEseq2安装
在命令行中输入以下内容:
```
if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")
BiocManager::install("DESeq2")
```
# 2 准备文件
在做差异表达基因分析时需要两个输入文件,一个是表达矩阵文件,一个是分组信息,将两个文件整理为“.txt”格式。
**表达矩阵文件(第一列为基因ID),行为样本信息,列为count值(一般不用FPKM值):**
| gene_id | NIP_1 | NIP_2 | NIP_3 | MUT_1 | MUT_2 | MUT_3 |
| ------------ | ----- | ----- | ----- | ----- | ----- | ----- |
| Os01g0100100 | 643 | 589 | 591 | 644 | 778 | 650 |
| Os01g0100200 | 9 | 7 | 7 | 13 | 5 | 13 |
| Os01g0100400 | 187 | 171 | 130 | 189 | 182 | 184 |
| Os01g0100500 | 1131| 1084| 1033| 1226| 1165| 1160|
| Os01g0100600 | 280 | 332 | 299 | 310 | 415 | 406 |
| Os01g0100700 | 2135| 1883| 1531| 1913| 1593| 2304|
| Os01g0100800 | 400 | 389 | 288 | 305 | 325 | 405 |
| Os01g0100900 | 1384| 1724| 1356| 1629| 1678| 1480|
| Os01g0101200 | 81 | 69 | 71 | 64 | 86 | 83 |
| Os01g0101300 | 57 | 32 | 46 | 29 | 48 | 62 |
| Os01g0101600 | 1445| 1190| 1221| 1437| 1814| 1345|
**分组文件,第一列为样本名,第二列为分组名:**
| sample | group |
| ------ | ----- |
| NIP_1| NIP |
| NIP_2| NIP |
| NIP_3| NIP |
| MUT_1| MUT |
| MUT_2| MUT |
| MUT_3| MUT |
# 3 使用R语言进行分析,话不多说,直接上脚本
```
#---导入R包---
library(DESeq2)
library(dplyr)
#---导入文件---
countdata <- read.table("C:/Users/bf/Desktop/新建文件夹/shootOEvsNIP.txt", header = T,row.names = 1) #表达矩阵文件,注意更换为自己的文件
countdata = countdata
coldata<- read.table("C:/Users/bf/Desktop/新建文件夹/sampleOE_NIP.txt", header = T, row.names = 1) #样本分组文件,注意更换为自己的文件
all(rownames(coldata) %in% colnames(countdata))
all(rownames(coldata) == colnames(countdata))
#---制作差异矩阵---
dds <-DESeqDataSetFromMatrix(countData = countdata,colData = coldata,design = ~ group)
dim(dds)
#---过滤---
dds <- dds
nrow(dds)
#---差异比较---
dep <- DESeq(dds)
res <- results(dep)
diff = res
diff <- na.omit(diff)## 去除NA
dim(diff)
#---导出文件---
write.csv(diff,"C:/Users/bf/Desktop/shootOEvsNIP.xls")
```
# 4 输出文件
输出文件格式如下:
| | baseMean | log2FoldChange | lfcSE | stat | pvalue | padj |
| ------------ | -------- | -------------- | ----- | ------ | ------ | ---- |
| Os01g0100100 | 1164.78| 0.43 | 0.11| 3.90 | 0.00 | 0.00 |
| Os01g0100200 | 26.25 | 1.68 | 0.51| 3.31 | 0.00 | 0.00 |
| Os01g0100400 | 256.76 | -0.09 | 0.24| -0.39| 0.69 | 0.76 |
| Os01g0100500 | 1463.07| 0.02 | 0.12| 0.17 | 0.87 | 0.90 |
| Os01g0100600 | 380.22 | -0.34 | 0.16| -2.17| 0.03 | 0.05 |
| Os01g0100700 | 2196.86| 0.96 | 0.11| 8.81 | 0.00 | 0.00 |
| Os01g0100800 | 473.50 | -0.11 | 0.17| -0.64| 0.52 | 0.60 |
| Os01g0100900 | 1786.58| -0.36 | 0.11| -3.31| 0.00 | 0.00 |
| Os01g0101200 | 686.24 | -10.36 | 0.67| -15.54 | 0.00 | 0.00 |
| Os01g0101300 | 68.55 | 0.02 | 0.32| 0.07 | 0.94 | 0.96 |
| Os01g0101600 | 1886.67| -4.21 | 0.13| -33.60 | 0.00 | 0.00 |
| Os01g0101800 | 89.00 | -0.27 | 0.32| -0.85| 0.40 | 0.48 |
| Os01g0102000 | 219.78 | 0.62 | 0.25| 2.45 | 0.01 | 0.03 |
· gene_id: 基因或转录本的标识符。
· baseMean: 平均表达量。
· log2FoldChange: 对数变换后的 fold change,表示在不同条件之间的表达倍数变化。
· lfcSE: log2 fold change 的标准误差。
· tat: 统计检验值。
· pvalue: 未经校正的 p 值。
· padj: 经过多重假设检验校正后的调整 p 值(通常使用 FDR 校正),用于控制假阳性发现率。
页:
[1]