基因课FTP地址:ftp://http://gsx.genek.tv/2020-3-10%E7%9B%B4%E6%92%AD%E4%B8%80%E4%B8%AA%E5%AE%8C%E6%95%B4%E7%9A%84%E8%BD%AC%E5%BD%95%E7%BB%84%E9%A1%B9%E7%9B%AE/
听张旭东老师的课
加载PCAtools
library(PCAtools)
主成分分析
p <- pca(gene_exp, metadata = sample_info, removeVar = 0.1) # 一般matadata就指样本信息表, removeVar去掉表达量变化最小的10%的基因
- removeVar 该参数可以对数据进行过滤
- 一般要算Z score
- 对基因降维分析
主成分分析的降维对应关系
- pca_loadings <- p$loadings # 提取原来基因与24个主成分的关系
loadings为主成分负载,即线性变换中的加权数 - pca_rotated <- p$rotated # 提取每个基因对样本的描述值
结果解释
- screeplot(p)
该图描述的是主成分对样本差异的解释度 - biplot(p,
x = 'PC1',
y = 'PC2',
colby = 'strain',
colkey = c('KID' = 'red', 'BLO' = 'yellow'),
shape = 'stage', l
egendPosition = 'top')
该图是用pca_rotated画的
colby指定按照颜色分的组(列)
colkey指定颜色
shape指定性状分的组(列)
legendPosition图例指定位置
PCA图加圈要用到ggplot2 - plotloadings(p)
查看各个主成分中哪些个基因对主成分贡献大
注
- 主成分分析中PC的数目由样本数决定
- 如果发现PC1可以区分重要的性状,想要研究PC1中重要的基因,则打开pca_loading表,找出所占比重最高的几个基因进行研究 —— 重点基因
- PC有可能只是将取样批次区分开,即为批次效应
- 批次效应的消除,用到的R软件包是SVA(表达芯片常用,表达测序数据也可以用)
- 在哪一步消除批次效应?
reads_count矩阵 → TPM → TMM → 消除批次效应
TPM消除组内基因长度、不同测序深度影响
TMM消除组间,可部分消除批次效应 - 如何能观察到批次效应?
样本聚类、相关性分析、PCA、箱线图(boxplot)
boxplot(gene_exp) 或 boxplot(log10(gene_exp + 1))# 箱线图可以展示,批次效应表现为:同一批次中位数在一个水平 - 芯片数据要用limma等标准化处理
网友评论