使用 "MatrixEQTL "软件包进行实用的 eQTL 分析

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 # 查看结果