跳转至

生信面试 R 语言统计题精选

一句话概述:生信面试中 R 语言题集中在数据可视化(ggplot2)、差异分析(DESeq2/edgeR)、统计检验(t检验/ANOVA/非参)和生存分析四个方向,掌握 tidyverse 生态就能应对大部分考题。


核心知识点速查表

概念白话解释
tidyverseR语言数据科学工具集(dplyr+ggplot2+tidyr等)
ggplot2R的王牌绑图包,图层叠加语法
dplyr数据操作五大动词:filter/select/mutate/arrange/summarise
DESeq2RNA-seq差异表达分析的标准工具
BioconductorR的生物信息专用包管理平台
正态检验判断数据是否符合正态分布(决定用什么检验)
参数检验假设数据符合某种分布(如t检验假设正态)
非参检验不假设分布(如Wilcoxon检验)

一、R 基础操作(快速回顾)

# ========== R基础操作 ==========

# 数据类型
x <- 5                                # 赋值(推荐用<-)
name <- "TP53"                         # 字符串
flag <- TRUE                           # 逻辑值
vec <- c(1, 2, 3, 4, 5)               # 向量(一维数组)
mat <- matrix(1:12, nrow=3, ncol=4)    # 矩阵(二维数组)
df <- data.frame(                      # 数据框(最常用!类似Excel表)
    gene = c("TP53", "BRCA1", "EGFR"),
    logFC = c(2.5, -1.8, 3.2),
    padj = c(0.001, 0.03, 0.0001)
)
lst <- list(name="John", age=25)       # 列表(可以混合类型)

# 向量操作
length(vec)                            # 长度
vec[1]                                 # 第1个元素(R从1开始!不是0!)
vec[c(1,3)]                            # 第1和第3个元素
vec[vec > 3]                           # 筛选>3的元素

# 数据框操作
df$gene                                # 提取gene列
df[1, ]                                # 第1行
df[, 2]                                # 第2列
df[df$padj < 0.05, ]                   # 筛选padj<0.05的行
nrow(df)                               # 行数
ncol(df)                               # 列数
str(df)                                # 查看结构
summary(df)                            # 统计摘要

# 读写文件
df <- read.csv("data.csv")             # 读CSV
df <- read.delim("data.tsv")           # 读TSV(tab分隔)
write.csv(df, "output.csv", row.names=FALSE)  # 写CSV

二、数据处理(dplyr,高频考点)

# ========== dplyr数据处理五大动词 ==========
library(dplyr)                         # 加载dplyr

# 示例数据
expr <- data.frame(
    gene = paste0("Gene", 1:100),      # 100个基因
    logFC = rnorm(100, mean=0, sd=2),  # 随机log2FC
    padj = runif(100, 0, 1),           # 随机padj
    baseMean = rpois(100, lambda=500)  # 随机表达量
)

# filter — 行筛选(选行)
sig_genes <- expr %>%                  # %>% 管道符(把结果传给下一步)
    filter(abs(logFC) > 1 & padj < 0.05)  # 显著差异基因
nrow(sig_genes)                        # 多少个显著基因

# select — 列选择(选列)
expr %>% select(gene, logFC, padj)     # 只保留3列
expr %>% select(-baseMean)             # 去掉baseMean列

# mutate — 新增/修改列
expr <- expr %>%
    mutate(
        direction = ifelse(logFC > 0, "Up", "Down"),  # 新增方向列
        neg_log_padj = -log10(padj)    # 计算-log10(padj)
    )

# arrange — 排序
expr %>% arrange(padj)                 # 按padj升序
expr %>% arrange(desc(abs(logFC)))     # 按|logFC|降序

# summarise + group_by — 分组统计
expr %>%
    filter(padj < 0.05) %>%            # 先筛选显著的
    group_by(direction) %>%            # 按方向分组
    summarise(
        count = n(),                   # 每组计数
        mean_logFC = mean(logFC),      # 每组平均logFC
        median_padj = median(padj)     # 每组中位padj
    )

# 面试经典:合并两个表格
gene_info <- data.frame(
    gene = paste0("Gene", 1:50),
    pathway = sample(c("apoptosis","cell_cycle","metabolism"), 50, replace=TRUE)
)
merged <- expr %>%
    inner_join(gene_info, by="gene")   # 按gene列合并(只保留交集)
# left_join 保留左表所有行
# right_join 保留右表所有行
# full_join 保留所有行

三、数据可视化(ggplot2,必考)

火山图(Volcano Plot)

# ========== 火山图(差异表达可视化经典) ==========
library(ggplot2)                       # 加载ggplot2

# 准备数据
expr <- expr %>%
    mutate(
        significance = case_when(      # 根据条件分类
            padj < 0.05 & logFC > 1 ~ "Up",        # 显著上调
            padj < 0.05 & logFC < -1 ~ "Down",     # 显著下调
            TRUE ~ "NS"                             # 不显著
        )
    )

# 画火山图
ggplot(expr, aes(x=logFC, y=-log10(padj),         # X轴logFC,Y轴-log10(padj)
                 color=significance)) +             # 按显著性着色
    geom_point(alpha=0.6, size=1.5) +              # 散点图
    scale_color_manual(values=c(                    # 自定义颜色
        "Up"="red", "Down"="blue", "NS"="grey"
    )) +
    geom_hline(yintercept=-log10(0.05),            # padj=0.05水平线
               linetype="dashed", color="grey50") +
    geom_vline(xintercept=c(-1, 1),                # logFC=±1垂直线
               linetype="dashed", color="grey50") +
    labs(title="Volcano Plot",                      # 标题
         x="log2(Fold Change)",                     # X轴标签
         y="-log10(adjusted p-value)") +            # Y轴标签
    theme_classic()                                 # 经典主题

ggsave("volcano_plot.pdf", width=8, height=6)      # 保存为PDF

热图(Heatmap)

# ========== 热图(表达谱可视化) ==========
library(pheatmap)                      # 加载pheatmap

# 模拟表达矩阵
set.seed(42)                           # 设置随机种子(可重现)
expr_matrix <- matrix(rnorm(200),      # 随机数据
                       nrow=20, ncol=10)
rownames(expr_matrix) <- paste0("Gene", 1:20)   # 行名=基因
colnames(expr_matrix) <- paste0("Sample", 1:10) # 列名=样本

# 注释信息
annotation_col <- data.frame(
    Group = rep(c("Control", "Treatment"), each=5),  # 分组
    row.names = colnames(expr_matrix)
)

# 画热图
pheatmap(expr_matrix,
         scale = "row",                # 按行标准化(Z-score)
         clustering_method = "ward.D2", # 聚类方法
         annotation_col = annotation_col,  # 列注释
         color = colorRampPalette(c("blue","white","red"))(100),  # 颜色
         fontsize_row = 8,             # 行字体大小
         fontsize_col = 10,            # 列字体大小
         main = "Gene Expression Heatmap")  # 标题

箱线图 + 统计检验

# ========== 箱线图 + 统计检验 ==========
library(ggplot2)
library(ggpubr)                        # 加载ggpubr(添加统计检验)

# 模拟数据
data <- data.frame(
    group = rep(c("Control", "Treatment"), each=30),  # 两组各30个
    expression = c(rnorm(30, mean=5, sd=1),           # 对照组
                   rnorm(30, mean=7, sd=1.5))          # 处理组
)

# 箱线图 + t检验p值
ggplot(data, aes(x=group, y=expression, fill=group)) +
    geom_boxplot(alpha=0.7) +                          # 箱线图
    geom_jitter(width=0.2, alpha=0.3) +                # 散点抖动
    stat_compare_means(method="t.test",                # 添加t检验
                       label="p.format") +             # 显示p值
    labs(title="Gene Expression Comparison",
         y="Expression (TPM)") +
    theme_classic() +
    theme(legend.position="none")                      # 隐藏图例

ggsave("boxplot_comparison.pdf", width=5, height=5)

四、统计检验(核心考点)

# ========== 统计检验完整流程 ==========

# === 第 1 步:正态性检验(决定用参数还是非参) ===
# Shapiro-Wilk检验(p>0.05 = 正态分布)
shapiro.test(data$expression[data$group == "Control"])    # 检验对照组
shapiro.test(data$expression[data$group == "Treatment"])  # 检验处理组

# === 第 2 步:方差齐性检验 ===
# Levene检验(p>0.05 = 方差齐)
car::leveneTest(expression ~ group, data=data)

# === 第 3 步:选择检验方法 ===

# 情况1:正态+两组 → t检验
t.test(expression ~ group, data=data)              # 独立样本t检验
t.test(expression ~ group, data=data, var.equal=TRUE)  # 方差齐时用等方差t检验
t.test(before, after, paired=TRUE)                  # 配对t检验

# 情况2:不正态+两组 → Wilcoxon检验(非参数)
wilcox.test(expression ~ group, data=data)          # Mann-Whitney U检验

# 情况3:正态+三组及以上 → ANOVA
# 模拟三组数据
data3 <- data.frame(
    group = rep(c("A", "B", "C"), each=20),
    value = c(rnorm(20,5,1), rnorm(20,7,1), rnorm(20,6,1))
)
aov_result <- aov(value ~ group, data=data3)        # 单因素方差分析
summary(aov_result)                                  # 查看结果
TukeyHSD(aov_result)                                 # 事后两两比较

# 情况4:不正态+三组及以上 → Kruskal-Wallis检验
kruskal.test(value ~ group, data=data3)              # Kruskal-Wallis检验
pairwise.wilcox.test(data3$value, data3$group,       # 事后两两比较
                      p.adjust.method="BH")          # BH多重检验校正

# === 多重检验校正 ===
pvalues <- c(0.01, 0.04, 0.03, 0.001, 0.06)         # 原始p值
p.adjust(pvalues, method="BH")                        # BH校正(FDR)
p.adjust(pvalues, method="bonferroni")                 # Bonferroni校正(更保守)

# 面试要点:
# BH(Benjamini-Hochberg)= 控制FDR(假发现率),生信最常用
# Bonferroni = 控制FWER(家族错误率),最保守
# 为什么要多重检验校正?因为同时做很多检验,偶然发现显著的概率增加

五、差异表达分析(DESeq2)

# ========== DESeq2差异表达分析流程 ==========
# BiocManager::install("DESeq2")       # 安装
library(DESeq2)                        # 加载

# 准备输入:原始counts矩阵 + 样本信息
# counts矩阵:行=基因,列=样本,值=raw read counts
count_matrix <- as.matrix(read.csv("counts.csv", row.names=1))

# 样本信息
coldata <- data.frame(
    condition = factor(c(rep("Control",3),     # 3个对照
                          rep("Treatment",3)),  # 3个处理
                       levels=c("Control","Treatment")),  # Control为参考
    row.names = colnames(count_matrix)          # 行名=样本名
)

# 创建DESeq2对象
dds <- DESeqDataSetFromMatrix(
    countData = count_matrix,                   # counts矩阵
    colData = coldata,                          # 样本信息
    design = ~ condition                        # 设计公式
)

# 预过滤低表达基因(加速运算)
keep <- rowSums(counts(dds)) >= 10              # 至少总计10个reads
dds <- dds[keep, ]                              # 过滤

# 运行差异分析
dds <- DESeq(dds)                               # 核心函数(做了3步)
# 1. estimateSizeFactors — 标准化
# 2. estimateDispersions — 估计离散度
# 3. nbinomWaldTest — 负二项分布Wald检验

# 提取结果
res <- results(dds,
               contrast=c("condition","Treatment","Control"),  # 比较组
               alpha=0.05)                       # 显著性阈值
res <- res[order(res$padj), ]                    # 按padj排序
summary(res)                                     # 结果摘要

# 保存结果
write.csv(as.data.frame(res), "DESeq2_results.csv")

# 可视化:MA图
plotMA(res, ylim=c(-5,5),                        # MA图
       main="DESeq2 MA Plot")

# 可视化:PCA图
vsd <- vst(dds, blind=FALSE)                     # 方差稳定化转换
plotPCA(vsd, intgroup="condition") +               # PCA图
    theme_classic()

六、生存分析

# ========== Kaplan-Meier生存分析 ==========
library(survival)                      # 生存分析核心包
library(survminer)                     # 漂亮的生存曲线

# 模拟临床数据
clinical <- data.frame(
    patient = paste0("P", 1:100),
    time = rexp(100, rate=0.02),       # 生存时间(指数分布)
    status = sample(0:1, 100,          # 1=事件(死亡), 0=删失
                    replace=TRUE, prob=c(0.3,0.7)),
    gene_group = sample(c("High","Low"), 100,  # 基因高/低表达组
                        replace=TRUE)
)

# 拟合生存模型
surv_obj <- Surv(time=clinical$time,   # 生存时间
                  event=clinical$status) # 事件状态

fit <- survfit(surv_obj ~ gene_group,  # 按基因分组
               data=clinical)
print(fit)                              # 打印中位生存时间

# 画Kaplan-Meier生存曲线
ggsurvplot(fit,
           data=clinical,
           pval=TRUE,                  # 显示log-rank p值
           risk.table=TRUE,            # 风险表
           conf.int=TRUE,              # 置信区间
           palette=c("red","blue"),    # 颜色
           xlab="Time (months)",       # X轴
           ylab="Survival Probability",# Y轴
           title="Kaplan-Meier Survival Curve",
           ggtheme=theme_classic())    # 主题

# Cox回归(多变量)
cox_model <- coxph(surv_obj ~ gene_group + age + stage,
                    data=clinical)      # 多因素Cox回归
summary(cox_model)                      # 查看HR、p值、CI

七、常见报错与解决

问题原因解决方法
object not found变量名拼写错误或未赋值检查变量名,ls() 列出所有变量
package not found包未安装install.packages("pkg") 或 BiocManager
non-numeric argument对非数值数据做数学运算as.numeric() 转换,检查NA
Error in plot数据中有NA或Infna.omit()is.finite() 过滤
factor levels错误分组变量顺序不对factor(x, levels=c("Control","Treatment"))
中文乱码编码问题Sys.setlocale("LC_ALL","UTF-8")

八、面试高频问题

Q1:R 中 data.frame 和 matrix 的区别?

data.frame 可以混合不同类型(数值+字符),matrix 只能存同一类型。data.frame 有列名可以用 $ 访问,适合存储样本信息等异质数据。matrix 运算速度更快,适合存储表达矩阵等同质数据。

Q2:什么时候用 t 检验,什么时候用 Wilcoxon?

数据正态分布 + 方差齐 → t检验(参数检验,统计效力更强)。数据不正态或样本量很小 → Wilcoxon(非参数检验,不依赖分布假设)。先做 Shapiro-Wilk 正态性检验再决定。

Q3:DESeq2 和 edgeR 的区别?

都用负二项分布建模counts数据。DESeq2用中位数比值法(median of ratios)标准化,Wald检验;edgeR用TMM标准化,精确检验或QLF检验。样本量<3时edgeR更稳定,大样本两者结果接近。一般推荐DESeq2(更保守、假阳性更低)。

Q4:p值和FDR(q值)的区别?

p值是单次检验犯错的概率。当同时检验上万个基因时,即使p<0.05也会有大量假阳性(5%×20000=1000个假阳性)。FDR(False Discovery Rate)控制的是"你说显著的里面有多少比例是假的",padj<0.05意味着预期假阳性不超过5%。生信中必须用FDR校正后的padj。


九、速查表

# === 生信 R 语言速查 ===

# 数据操作 (dplyr)
df %>% filter(padj < 0.05)            # 筛选行
df %>% select(gene, logFC)            # 选择列
df %>% mutate(new_col = logFC * 2)    # 新增列
df %>% arrange(padj)                  # 排序
df %>% group_by(group) %>% summarise(n=n())  # 分组统计

# 可视化 (ggplot2)
ggplot(df, aes(x, y, color=group)) +  # 基本语法
    geom_point() +                     # 散点
    geom_boxplot() +                   # 箱线图
    geom_bar(stat="identity") +        # 柱状图
    theme_classic()                    # 主题

# 统计检验选择
# 2组正态 → t.test()
# 2组非正态 → wilcox.test()
# 3+组正态 → aov() + TukeyHSD()
# 3+组非正态 → kruskal.test()
# 多重校正 → p.adjust(method="BH")

# DESeq2
dds <- DESeqDataSetFromMatrix(counts, coldata, ~condition)
dds <- DESeq(dds)
res <- results(dds, alpha=0.05)

# 生存分析
fit <- survfit(Surv(time, status) ~ group, data)
ggsurvplot(fit, pval=TRUE, risk.table=TRUE)

参考资料:R for Data Science (Wickham 2023)、DESeq2 vignette、ggplot2: Elegant Graphics (Wickham 2024)、Bioconductor