原文下载
数据及代码下载
Korte A, Vilhjálmsson B J, Segura V, et al. A mixed-model approach for genome-wide association studies of correlated traits in structured populations[J]. Nature Genetics, 2012, 44(9):1066.
摘要
GWAS是研究标记和性状关系的通用方法。其中一个比较重要的一点是GWAS需要考虑数据之间的关系，包括位点间和个体间的。
混合线性模型因其自身特点，能够灵活的定义GWAS研究中的群体结构和数据关系。
这里我们考虑到GWAS中性状间的关系，扩展了混合线性模型的应用范围，提出了MTMM（多性状混合线性模型），他可以考虑性状间和性状内的变异，因此非常适合多性状的分析。
我们用人类的数据作为比较分析，结果表明，MTMM模型可以显著提高检验的效率和准确性，而且可以分析位点与环境之间的关系。
背景及前言模型
大部分GWAS的研究都是单性状和SNP的分析，但是很多性状都是多基因控制的，这导致如果不同群体结构，一些位点缺乏独立性（LD），从而影响结果。
当个体有多个测量性状时，性状之间可能是有关联的，就像处于同一环境或者LD状态，这些因素通常被忽视。但是这些因素是很重要的，如果我们在分析时考虑在内，会降低误差，从而提高检测的功效。
比如田间试验中，对不同环境自交系表型值的测量，这些表型值在不同环境中不是独立的，它们与环境存在着互作。
在数量遗传学中，多性状分析的研究有很长的历史，但是很少在GWAS的研究应用。在这篇文章中，我们通过混合线性模型研究动物模型中多性状的关系。群体结构和亲缘关系kinship在GWAS中会被考虑到，这里我们又增加了多性状之间的关系，它可以考虑性状间和性状内的方差组分，用来进行GWAS的分析。我们的结果表明，它可以提高对位点的检测效率，同时保证较低的假阳性率。
模型
用三种模型进行比较：
- MTMM
- marginal
- single-trait analysis
测试模型用三种方法：
- full test，主要是比较：marker效应和交互效应 VS 无效应
- interaction effect test：交互效应 VS 无效应
- common effect test：主要是比较有一个基因型marker效应 VS无效应
程序运行
MTMM - A mixed-model approach for genome-wide association studies of correlated traits in structured populations
Introduction
The MTMM function as published in Nature Genetics currently don't support estimates on missing data and replicates.
This is work in progress and will be accordingly updated here.
For questions and comments feel free to contact me: arthur.korte@gmi.oeaw.ac.at
How to use
# Load libraries and source needed functions
# The AsREML package needs a valid license that can be obtained at http://www.vsni.co.uk/software/asreml
library(lattice)
library(asreml)
# msm and nadiv librarys are used to estimate SE of the correlation estimates, only used if run=FALSE
#library(msm)
#library(nadiv)
source('mtmm_function.r')
source('emma.r')
# load your data (Phenotype(Y),Genotype(X) and Kinship(K))
# note you can calculate K using the emma package K<-emma.kinship(t(X)), make sure to set colnames(K)=rownames(K)=rownames(X)
# alternativley load the sample data
load('data/MTMM_SAMPLE_DATA.Rdata')
# different options include method(default or errorcorrelation, include.single.analysis, calculate.effect.size (if TRUE, ###analysis is more time consuming) default for X is binary coding of 0 and 1, if your data are code 0,1 and 2 use ###gen.data='heterozygot', run=FALSE will not perform the GWAS, but only output the correlation estimates (fast)
mtmm(Y,X,K,method='default',include.single.analysis=T,calculate.effect.size=T,gen.data='binary',exclude=T,run=T)
# the function outputs a list called results ($phenotype ,$pvals, $statistics, $kinship)
output<-results$pvals
# manhattan plots
# default plots for include.single.analysis=T
par(mfrow=c(5,1),mar=c(3, 4, 1, 4))
plot_gwas(output,h=8)
plot_gwas(output,h=9)
plot_gwas(output,h=10)
plot_gwas(output,h=11)
plot_gwas(output,h=12)
#qq plots
par(mfrow=c(1,1),mar=c(3, 4, 1, 4))
qq_plot_all(output)
