DESeq2分析预处理

一、在R中去除低表达基因的一个常用方法是设定一个阈值，比如，一个基因在一定数量的样本中至少需要达到某个计数值。以下是一个示例R代码，展示如何去除在少于一半的样本中计数少于10的基因：

# Load necessary libraries
library(DESeq2)

# Read the count data from the CSV file without setting row names
counts_data <- read.csv("merged_count.csv", header = TRUE)

# Assume the first column contains the gene identifiers, so we skip it when applying the threshold
# Here we use a threshold of 10 counts, and the gene must be expressed at this level in at least half the samples
threshold <- 10
minimum_samples <- ceiling(ncol(counts_data) / 2)
filtered_counts <- counts_data[rowSums(counts_data[, -1] >= threshold) >= minimum_samples, ]

# Save the filtered data to a new CSV file
write.csv(filtered_counts, "filtered_merged_count.csv", row.names = FALSE)

在上面的代码中，rowSums(counts_data[, -1] >= threshold) 计算每一行（基因）有多少个样本满足计数值大于或等于阈值，然后检查这个数量是否至少为总样本数的一半。根据您的具体情况，您可能需要调整阈值和样本数的比例。

请注意，这个方法是基于一个简单的阈值，而实际中可能需要一个更复杂的方法来确定低表达基因，这可能涉及到基因的生物学知识和数据的具体情况。如果您的数据处理需要更专业的生物信息学方法，请咨询相关专家。

二、将csv文件中第一列的ENSG ID中的小数点及小数点后数字去掉

# 读取数据
counts_data <- read.csv("filtered_merged_count.csv")

# 使用sub函数去除第一列中的小数点和随后的数字
# 这里假设第一列的名称是"Geneid"
counts_data$Geneid <- sub("//..*", "", counts_data$Geneid)

# 保存修改后的数据到新的CSV文件
write.csv(counts_data, "modified_filtered_merged_count.csv", row.names = FALSE)

三、将这个文件第一列的ENSG ID转化为GeneSymbol和GeneID，分别将每行转化后的结果放在本行的第二列和第三列

#######记得在packages中把biocManager和biomaRt加载

# 首先，安装并加载必要的包。如果您还没有安装biomaRt，请取消注释并运行下面的两行代码。
# if (!requireNamespace("BiocManager", quietly = TRUE))
#     install.packages("BiocManager")
# BiocManager::install("biomaRt")

# 加载biomaRt包
library(biomaRt)

# 设置ensembl mart
ensembl <- useEnsembl(biomart = "genes", dataset = "hsapiens_gene_ensembl")

# 读取CSV文件
counts_data <- read.csv("modified_filtered_merged_count.csv", stringsAsFactors = FALSE)

# 提取ENSG IDs
ensg_ids <- counts_data[, 1]

# 查询Ensembl数据库获取Gene Symbols 和 Entrez Gene IDs
annotations <- getBM(attributes = c('ensembl_gene_id', 'hgnc_symbol', 'entrezgene_id'),
                     filters = 'ensembl_gene_id', values = ensg_ids, mart = ensembl)

# 因为有时候annotations的顺序可能和原始的ensg_ids不一致，所以我们需要以ensg_ids为基础进行左连接
counts_data_annotated <- merge(counts_data, annotations, by.x = 'Geneid', by.y = 'ensembl_gene_id', all.x = TRUE, sort = FALSE)

# 将Gene Symbol和Gene ID插入到第一列后面
counts_data_annotated <- counts_data_annotated[, c(1, ncol(counts_data_annotated)-1, ncol(counts_data_annotated), 2:(ncol(counts_data_annotated)-2))]

# 为新插入的列重命名
names(counts_data_annotated)[2] <- "GeneSymbol"
names(counts_data_annotated)[3] <- "GeneID"

# 将新的数据框保存为CSV文件
write.csv(counts_data_annotated, "updated_modified_filtered_merged_count.csv", row.names = FALSE)

四、检测这个文件第一列各行有没有重复内容，并将重复行的内容只保留一行，把剩余重复行删除

# 读取CSV文件
data <- read.csv("/mnt/data/updated.csv", header = TRUE)

# 查找第一列中的重复项，并排除第一次出现的重复项
data_unique <- data[!duplicated(data[[1]]), ]

# 将处理后的数据写回到一个新的CSV文件
write.csv(data_unique, "/mnt/data/updated_unique.csv", row.names = FALSE)

五、在R中使用 DESeq2 包进行差异表达分析通常涉及以下步骤：

加载必要的库（如 DESeq2 和其他数据处理工具）。
读取数据和样本信息。
创建 DESeq2 数据集对象。
运行差异表达分析。
提取结果并将其写入文件。
下面是一个基本的R脚本，包含执行上述步骤的代码：

# 加载所需的库
library(DESeq2)
library(tibble)

# 假设有两组样品，每组3个生物学重复
# 读取数据文件
data <- read.csv("/mnt/data/updated_unique.csv", header = TRUE)

# 样品信息，根据你的样品布局修改colData的内容
# 这里假设文件的列格式是：第一列为基因ID，后面六列分别为三个重复的两组样品
# 样本名称为 Sample1_rep1, Sample1_rep2, Sample1_rep3, Sample2_rep1, Sample2_rep2, Sample2_rep3
colData <- data.frame(
  sampleName = colnames(data)[-1],
  condition = rep(c("Control", "Treatment"), each = 3)
)

# 为DESeq2创建矩阵，行为基因，列为样品
dds <- DESeqDataSetFromMatrix(
  countData = data[, -1],
  colData = colData,
  design = ~ condition
)

# 运行DESeq2差异表达分析
dds <- DESeq(dds)

# 获取结果
res <- results(dds)

# 将结果保存为CSV文件
write.csv(as.data.frame(res), file = "/mnt/data/deseq2_results.csv")

在运行此代码之前，请确保 updated_unique.csv 文件的格式符合要求。第一列应该是基因或特征的标识符，后面的列应该是计数数据。

你需要在你的R环境中执行这段代码，例如在RStudio中。如果你希望我在这里运行这段代码，请确保我理解了样本的布局，并提供必要的信息以正确设置样本信息（colData）。