在GSE46224文件夹里，定量以后，获得(15+8)个fc.txt文件，这个是RNA-seq测序数据的表达矩阵，RNA-seq的count数据符合负二项分布，因此用SVA包时，要用ComBat-Seq函数(针对RNA-seq等离散型数据)，不能用combat函数(针对芯片等连续型变量数据)。在RUVSeq包里，用RUVg函数。
分析流程如下：
1.R
读取二十三个样本、合并，获取未进行id转换的原始表达矩阵：
allgeneid
然后用biomaRt包整理基因注释信息，将geneid转换为gene名，获得genematrix表达矩阵
genematrix
并储存。
save(genematrix,file = "genematrix.data.Rdata")
2.R
在genematrix里，行名为阿拉伯数字，第一列为gene名，要想进行后续分析，要把行名设置为基因名，但是获得的表达矩阵里，基因名有很多重复，如果设置为行名会报错，因此就要去重复。
去重复后获得genematrix2表达矩阵，行数从67103到38570，储存。
genematrix2
save(genematrix2,file = "genematrix2.data.Rdata")
3.R
这一步主要是基于RUVseq这个R包去批次以及差异分析
原始数据为genematrix2，这个矩阵是原始矩阵经过去重复及ID转换后的。
过滤数据：#去除表达量为0的，行数从38570到18561，得到genematrix3表达矩阵。
save(genematrix3,file = "genematrix3.data.Rdata")
标准化：使用上四分位数法，标准化前后分别画PCA图
用edgeR的方法找到DEG，挑选排名靠前的DEG，然后将排名靠后的DEG作为negative genes进行后续RUVg的normalization。
这时再画PCA图，HF和normal的样本已经被区分的非常开了。
最后再用edgeR进行GLM建模分析DEG，最后获得差异基因的矩阵DEG_result46224。
DEG_result46224
将行名变为第一列，最后变为下列形式：
DEG_result46224.jpg
4.R
这一步的主要目的：画出差异最大的几个基因；对任意基因画在两组里表达量的箱式图
所需输入文件：
1.分组文件：
group_list46224
2.差异基因
差异基因
3.表达矩阵并进行归一化
表达矩阵
利用这三者去画图
BNP表达量
网友评论