************************************************************ 差异基因：DESeq2/edgeR/limma ************************************************************

安装数据分析包包

if (!requireNamespace("BiocManager", quietly = TRUE))
  install.packages("BiocManager")
BiocManager::install("DESeq2")
BiocManager::install("edgeR")
BiocManager::install("limma")

设定目录及病例对照个数

setwd("C:\\Users\\acer\\Desktop\\测序数据处理\\2v2\\")
controlSampleNum=2
treatSampleNum=2

用DESeq2进行分析

library("DESeq2")
library("limma")

读取数据及分组

database_all <- read.csv(file = "gene_count_all.csv",header = T, row.names=1,check.names=F)
head(database_all)
data = database_all[,1:(controlSampleNum+treatSampleNum)]
data <- round(as.matrix(data))
head(data)
#ann = database_all[,(controlSampleNum+treatSampleNum+1):ncol(database_all)]
#head(ann)
group <- factor(c(rep("control",controlSampleNum), rep("treat",treatSampleNum)), levels = c("control","treat"))
group
col_data <- data.frame(row.names = colnames(data), group)
col_data

进行差异分析

dds <- DESeqDataSetFromMatrix(data, col_data, design = ~group)
keep = rowSums(counts(dds)) > 2
dds = dds[keep,]
DDS <- DESeq(dds)
result = results(DDS, contrast=c("group", "treat","control"))
result = result[order(result$padj),]
head(result)

结果输出

write.csv(result, file="Diff_all_gene_DESeq2.csv")
#write.csv(ann, file="ann.csv")
diff_DESeq2 = read.csv("Diff_all_gene_DESeq2.csv", header = TRUE)
ann = read.csv("gene_ann.csv", header = TRUE)
result = merge(diff_DESeq2,ann,by = "gene_id")
write.csv(result, file="Diff_all_gene_DESeq2_F.csv")

使用edgeR进行分析

library("edgeR")
library("limma")

读取数据及分组

database_all <- read.table(file = "miRNA_Matrix.txt",header = T, row.names = 1,check.names=F)
head(database_all)
data = database_all[,1:(controlSampleNum+treatSampleNum)]
data <- round(as.matrix(data))
head(data)
ann = database_all[,(controlSampleNum+treatSampleNum+1):ncol(database_all)]
head(ann)
group <- factor(c(rep("LC_1",controlSampleNum), rep("LC_2",treatSampleNum)), levels = c("LC_1","LC_2"))

进行差异分析

dge <- DGEList(counts=data,group=group)
keep <- rowSums(cpm(dge)>1) >= 2###过滤
dge <- dge[keep, , keep.lib.sizes=FALSE]
dge <- calcNormFactors(dge)###归一化
design <- model.matrix(~group)
DGE <- estimateDisp(dge,design)
DGE <- exactTest(DGE)
result = topTags(DGE, n = nrow(DGE$table))$table
head(result)

结果输出

write.table(result, "diff_edgeR.csv", col.names = NA, sep = "," )
write.csv(ann, file="ann.csv")
diff_edgeR = read.csv("diff_edgeR.csv", header = TRUE)
ann = read.csv("ann.csv", header = TRUE)
result = merge(diff_edgeR,ann,by = "X")
write.csv(result, file="result_edgeR.csv")