1 GEO数据挖掘
数据类型
(1) 基因表达芯片
(2) 转录组
(3) 单细胞
(4) 突变、甲基化、拷贝数变异……
怎样筛选基因--分析方法
数据下载（DEO、TCGA）-差异分析（芯片与转录组不相同）-WGCNA（加权共表达网络）-富集分析（ORA、GSEA）-PPI网络-预后分析（影响生存的疾病）
1.1
1.1.1 热图
输入数值为数值型矩阵/数据框
以颜色变化代表数值大小
#聚类树：根据基因相似程度进行排序分类，与原表达矩阵基因顺序不同
1.1.2 散点图和箱线图
可以用箱线图代替散点图，显示整体差异
箱线图：
以连续型向量为纵坐标；有重复值的离散型向量为横坐标
箱线图的五条线
max - 75% - median#中位数 - 25% - min
最大值和最小值以外可能存在离群值#离群点
#用于单个基因在几组之间的表达差异
###多基因 --> 差异分析
1.1.3 火山图
两个数值：logFC、P.Value
logFC（横坐标）
Foldchange(FC)：处理组均值/对照组均值
log2Foldchange(logFC)：Foldchange取log2
#实际运算中先取log再相减
#logFC表示处理组和对照组相比的基因表达差异倍数
#存在负值，表示表达降低
#基因的上调/下调，指基因表达量显著上升/下降
--> P.Value
芯片差异分析的起点是一个取过log的表达矩阵（0-20）；
若未进行该操作，数值将非常大，需要先取log
通常设置阈值，例如：
上调基因：logFC > 1, p < 0.01
下调基因：logFC < -1, p < 0.01
#logFC常见阈值：1, 2, 1.2, 1.5, 2.2, 0.585 = log2(1.5)
P.Value
会进行调整将其增大
-->
-log10(P.Value)
P.Value越小，-log10(P.Value)越大，差异越大的置信度越高
1.1.4 主成分分析
PCA样本聚类图降维
点与点之间的相对距离表示相似程度
横、纵坐标：Dimension（Dim1、2）——主成分（综合指标）
几个基因组合到一起成为一个主成分
例如：BMI
#括号内的数字越大越好，没有具体要求
#图中最大的点为聚类的中心点，不是样本点
#至少四个样本点才能在图中形成一簇
#将权重最高的两个主成分作为横、纵坐标，而非全部主成分
#用于简单查看组间是否存在差异
2 GEO背景知识及芯片表达分析思路
2.1 表达数据实验设计
实验目的：通过基因表达量数据的差异分析和富集分析来解释生物学现象
2.2 数据库介绍
Gene Expression Omibus
#Tools-Analyze a Study with GEO2R-代码，可以辅助完成一些操作
Sample：用户提交给GEO的样本数据（GSM）
Series：一个完整的研究，提供了整个研究的描述，包括对数据的描述、总结、分析（GSE）
Platform：用户测定表达量使用的芯片/平台（GPL）
基因表达芯片的原理：探针的表达量代表基因的表达量
#分析思路
找数据，找到GSE编号 -->
下载并读取数据 -->
#表达矩阵 #临床（分组）信息 #GPL编号（探针注释）
数据探索 -->
#分组间是否存在差异，PCA、热图
差异分析并可视化 -->
#P.Value, logFC #火山图、热图
富集分析
#KEGG #GO
为什么不画全部基因的热图
1* 数据太大
2* 并不是所有基因都存在差异
2.3 表达矩阵
行名：探针id #需要转换为gene symbol
列名：GSM，样本编号 #需要分组信息
3 代码分析流程
芯片差异分析所需输入数据
表达矩阵
#数据分布范围0-20
#无异常值，如NA、INF、负值
#无异常样本
分组信息
#同一分组对应同一关键词
#顺序与表达矩阵的列一一对应
#因子，对照组的levels在前
探针注释
#根据GPL编号查找
#探针与基因之间的对应关系
#只能有两列，且均为字符型
#列名必须是probe_id和symbol
批量装包代码
options("repos"="https://mirrors.ustc.edu.cn/CRAN/") if(!require("BiocManager")) install.packages("BiocManager",update = F,ask = F) options(BioC_mirror="https://mirrors.westlake.edu.cn/bioconductor") cran_packages <- c('tidyr', 'tibble', 'dplyr', 'stringr', 'ggplot2', 'ggpubr', 'factoextra', 'FactoMineR', 'devtools', 'cowplot', 'patchwork', 'basetheme', 'paletteer', 'AnnoProbe', 'ggthemes', 'VennDiagram', 'survminer', "tinyarray") Biocductor_packages <- c('GEOquery', 'GO.db', 'hgu133plus2.db', 'ggnewscale', "limma", "impute", "GSEABase", "GSVA", "clusterProfiler", "org.Hs.eg.db", "preprocessCore", "enrichplot")
for (pkg in cran_packages){ if (! require(pkg,character.only=T,quietly = T) ) { install.packages(pkg,ask = F,update = F) require(pkg,character.only=T) } }
for (pkg in Biocductor_packages){ if (! require(pkg,character.only=T,quietly = T) ) { BiocManager::install(pkg,ask = F,update = F) require(pkg,character.only=T) } }
#前面的所有提示和报错都先不要管。主要看这里
for (pkg in c(Biocductor_packages,cran_packages)){ require(pkg,character.only=T) } #没有任何提示就是成功了，如果有warning xx包不存在，用library检查一下。 #library报错，就单独安装。
实战代码
rm(list = ls()) #打破下载时间的限制,改前60秒，改后10w秒 options(timeout = 100000) options(scipen = 20)#不要以科学计数法表示
#传统下载方式
library(GEOquery) eSet = getGEO("GSE7305", destdir = '.', getGPL = F) #网速太慢，下不下来怎么办 #1.从网页上下载/发链接让别人帮忙下，放在工作目录里 #2.试试geoChina(),只能下载2019年前的表达芯片数据 library(AnnoProbe) eSet = geoChina("GSE7305") #选择性代替getGEO() #研究一下这个eSet class(eSet) length(eSet)
eSet = eSet[1] class(eSet) 1 "ExpressionSet" attr(,"package") 1 "Biobase" #“ExpressionSet”为来自R包“Biobase”中的一个对象
#(1)提取表达矩阵exp
exp <- exprs(eSet)
#⭐第一个要检查的地方
