WGCNA，cpm值和分组信息分别保存在2个txt文件中，bulk RNA-seq分析R语言代码？

以下是WGCNA和cpm值和分组信息分别保存在2个txt文件中的bulk RNA-seq分析的R语言代码：

首先，我们需要安装WGCNA包：

install.packages("WGCNA")
library(WGCNA)

然后，我们需要读取cpm值和分组信息并进行预处理：

# 读取cpm值文件
cpm_data <- read.table("cpm_data.txt", header=T, row.names=1)

# 读取分组信息文件
group_data <- read.table("group_data.txt", header=T, row.names=1)

# 转置矩阵
cpm_data <- t(cpm_data)

# 删除无用基因
keep <- rowSums(cpm_data)>1
cpm_data <- cpm_data[keep,]

# 对表达矩阵进行log2转换和标准化
cpm_data <- log2(cpm_data+1)
cpm_data <- scale(cpm_data)

# 将分组信息转换为因子变量
group_data$group <- as.factor(group_data$group)

接下来，我们使用WGCNA来构建共表达网络和进行模块识别：

# 构建共表达网络
powers <- c(1:10)
sft <- pickSoftThreshold(cpm_data, powerVector=powers, verbose=5)
power <- sft$powerEstimate

# 构建共表达网络
net <- blockwiseModules(cpm_data, power, TOMType="signed", 
                        minModuleSize=30, reassignThreshold=0, 
                        mergeCutHeight=0.25, numericLabels=TRUE, 
                        pamRespectsDendro=FALSE, 
                        saveTOMs=TRUE, saveTOMFileBase="TOM")

# 绘制模块-属性关联图
MEs <- net$MEs
plotMEs(MEs, main="Module Eigengene Associations")

最后，我们根据模块的基因表达模式进行功能注释和富集分析：

# 将模块的基因表达模式与外部信息进行关联
geneInfo <- data.frame(Module = net$colors, row.names = colnames(cpm_data))
group <- group_data$group

# 功能注释和富集分析
for (module in unique(geneInfo$Module)) {
  moduleGenes <- names(which(geneInfo$Module == module))
  moduleData <- cpm_data[moduleGenes,]
  moduleAnnotation <- geneInfo[geneInfo$Module == module,]
  
  # 进行富集分析
  moduleEnrichment <- enrichGO(moduleData, 
                               universe = rownames(cpm_data),
                               OrgDb = "org.Hs.eg.db",
                               ont = "BP",
                               pvalueCutoff = 0.05,
                               qvalueCutoff = 0.1,
                               readable = TRUE)
  
  # 输出结果
  cat(paste("Module", module, "enrichment analysis:\n"))
  print(moduleEnrichment)
  
  # 进行差异分析
  module_data <- data.frame(t(moduleData))
  module_data$group <- group
  module_data$group <- relevel(module_data$group, ref = "control")
  contrasts <- c("group", "treated", "control")
  design <- model.matrix(~group, data=module_data)
  colnames(design) <- contrasts
  fit <- glmFit(module_data, design)
  lrt <- glmLRT(fit, contrast = makeContrasts(treated - control, levels=contrasts))
  table <- topTags(lrt, n=nrow(moduleData))
  table <- table$table
  table$padj <- p.adjust(table$PValue, method="fdr")
  table <- table[table$padj < 0.05,]
  
  # 输出结果
  cat(paste("Module", module, "differential expression analysis:\n"))
  print(table)
}

这些代码可以帮助您开始使用WGCNA进行cpm值和分组信息分别保存在2个txt文件中的bulk RNA-seq数据分析。