在上一期内容中,米老鼠和大家介绍了eQTL的相关概念和分析原理,今天我就带大家用“MatrixEQTL”进行一下实战演练。
这里我们使用的是该包提供的内置数据集,代码如下:
代码语言:javascript复制install.packages("MatrixEQTL") # 安装R包
library("MatrixEQTL") # 加载R包
base.dir = find.package("MatrixEQTL") # 获取该包的路径信息
useModel = modelLINEAR # 三种模型可选(modelANOVA or modelLINEAR or modelLINEAR_CROSS)
SNP_file_name = paste(base.dir, "/data/SNP.txt", sep="") # 获取SNP文件位置
SNP_file = data.table::fread(SNP_file_name, header=T) # 读取SNP文件,可以在R中查看
expression_file_name = paste(base.dir, "/data/GE.txt", sep="") # 获取基因表达量文件位置
expression_file = data.table::fread(expression_file_name, header=T) # 读取基因表达量文件,可以在R中查看
covariates_file_name = paste(base.dir, "/data/Covariates.txt", sep="") # 读取协变量文件
covariates_file = data.table::fread(covariates_file_name, header=T) # 读取协变量文件,可在R中查看
output_file_name = tempfile() # 将输出文件设置为临时文件
pvOutputThreshold = 1e-2 # 定义gene-SNP associations的显著性P值
errorCovariance = numeric() # 定义误差项的协方差矩阵 (用的很少)
snps = SlicedData$new() # 创建SNP文件为S4对象(S4对象是该包的默认输入方式)
snps$fileDelimiter = "t" # 指定SNP文件的分隔符
snps$fileOmitCharacters = "NA" # 定义缺失值
snps$fileSkipRows = 1 # 跳过第一行(适用于第一行是列名的情况)
snps$fileSkipColumns = 1 # 跳过第一列(适用于第一列是SNP ID的情况
snps$fileSliceSize = 2000 # 每次读取2000条数据
snps$LoadFile( SNP_file_name ) # 载入SNP文件
# 下面的代码同snps的那部分一致
gene = SlicedData$new()
gene$fileDelimiter = "t"
gene$fileOmitCharacters = "NA"
gene$fileSkipRows = 1
gene$fileSkipColumns = 1
gene$fileSliceSize = 2000
gene$LoadFile( expression_file_name )
cvrt = SlicedData$new()
cvrt$fileDelimiter = "t"
cvrt$fileOmitCharacters = "NA"
cvrt$fileSkipRows = 1
cvrt$fileSkipColumns = 1
cvrt$fileSliceSize = 2000
cvrt$LoadFile( covariates_file_name )
# 文件的输入部分结束
me = Matrix_eQTL_engine( # 这是进行eQTL分析的主要函数
snps = snps, # 指定SNP 文件
gene = gene, # 指定基因表达量文件
cvrt = cvrt, # 指定协变量文件
output_file_name = output_file_name, # 指定输出文件
pvOutputThreshold = pvOutputThreshold, # 指定显著性P值
useModel = useModel, # 指定使用的计算模型
errorCovariance = errorCovariance, # 指定误差项的协方差矩阵
verbose = TRUE,
pvalue.hist = TRUE,
min.pv.by.genesnp = FALSE,
noFDRsaveMemory = FALSE)
res <- me$all$eqtls # 把eQTL的显著结果存储到变量res里
res # 查看结果
关于“MatrxiEQTL”包的用法就简单介绍到这里,感兴趣的朋友可以深入学习一下,掌握cis-eQTL和trans-eQTL的分析方法。
生信与临床