GEO数据库高通量测序R语言差异分析

在R语言中对GEO（Gene Expression Omnibus）数据库中的高通量测序数据进行差异表达分析，通常包括以下几个步骤：

1. **数据获取**：利用`geoquery`包从GEO下载实验设计文件（GPL、GSM），然后提取相应的表达矩阵（counts或normalized counts）。

library(geoquery)
gse <- getGEO("GSE accession", destdir=".", GSEMatrix=TRUE)
exp_data <- exprs(gse$elements[[1]])

2. **数据预处理**：处理缺失值和标准化数据，这可能包括log变换、归一化（如TMM或RLE），以及滤除低表达或高度变异性基因。

library(edgeR)
keep <- rowSums(cpm(exp_data) > 1) >= samples(gse)$nRow / 2 # 留下表达水平较高的基因
exp_data <- exp_data[keep,]

3. **设置模型和拟合**：对于配对样品或多组设计，可以建立负二项分布或一般线性模型（GLM）。例如，使用`edgeR`的`DGEList`和`glmQLFit`。

dge <- DGEList(counts=exp_data, group=samples(gse)$condition)
fit <- glmQLFit(dge, design=full.matrix(~ condition))

4. **检测差异表达**：使用`topTags`函数找出显著差异的基因。

tagged genes <- topTags(fit, n=nrow(exp_data), adjust.method="BH")$table

5. **结果可视化和报告**：使用`ggplot2`绘制 volcano plot 和 MA plot，展示差异表达结果。

library(gplots)
volcanoPlot(tagged_genes)