前言
關於全基因組關聯分析（GWAS）原理的資料，網上有很多。
這也是我寫了這么多GWAS的軟件教程，卻從來沒有寫過GWAS計算原理的原因。
恰巧之前微博上某位小可愛提問能否寫一下GWAS的計算原理。我一順口就答應了。
后面一直很懶，不願意動筆，但想着既然答應了，不寫說不過去。
我寫這段話的意思是，如果你有任何關於GWAS分析問題或者疑問，希望我能寫一下的，可以跟我說。
如果我認為有價值，寫出來對大家有幫助的話，會寫的。
GWAS所涉及的公式：最小二乘法
首先，我們來一個知識點的回顧：最小二乘法。
看下圖，熟不熟悉！
這可是我們中學時解了很多遍的算術題。
圖片來源：http://kitsprout.blogspot.com/2015/11/algorithm-least-squares.html
公式可以寫為： y = ax + b
y:我們研究的表型
x:基因型數據，這里指每一個SNP
a：SNP的系數
b：殘差，可以是環境變量，或者除了SNP之外的影響表型的因素
來個例子給我們講講唄，公式怎么套進去
如圖所示，假定有一個SNP,叫 rs123: T>C
我們定義C為風險位點，以加性模型為例，一個C=1,T=0
那么CC=2,CT=1,TT=0
根據上面的公式：
SNP對應的值x分別為：2,2,1,2,1,1,0,2
對應的表型y分別為10,7,6,8,5,4,2,6
回顧我們前面提到的公式：y = ax + b
現在我們有：
10= 2a+b
7= 2a+b
6= 1a+b
8= 2a+b
5= 1a+b
4= 1a+b
2= 0+b
6= 2a+b
轉化一下，就是：
2a+b - 10 = 0
2a+b - 7 = 0
1a+b - 6 = 0
2a+b - 8 =0
1a+b - 5 = 0
1a+b - 4 = 0
0+b -2 = 0
2a+b -6 = 0
我們的任務就是，找到合適的a，b使得
（2a+b - 10)^2 + (2a+b - 7)^2 + (1a+b - 6)^2 + (2a+b - 8)^2 + (1a+b - 5)^2 + (1a+b - 4)^2 + (0+b -2)^2 + (2a+b -6)^2 最小。
a，b的求值涉及到最小二乘法的推導，感興趣的看這篇文章：https://zhuanlan.zhihu.com/p/53556591
用公式表示就是：
b = cor(y,x)*Sd(y)/Sd(x)
a = (10+7+6+8+5+4+2+6)/8 - ((2+2+1+2+1+1+0+2)/8)*b
cor(y,x)表示x和y的相關系數
Sd(y),Sd(x)分別表示y和x的標准差
可以自己手算一下，也可以借助R語言：
x=c(2,2,1,2,1,1,0,2)
y=c(10,7,6,8,5,4,2,6)
Ex=mean(x);Ex
Ey=mean(y);Ey
Sx=sd(x);Sx
Sy=sd(y);Sy
corn=cor(y,x) ; corn
b=corn*Sy/Sx ; b
a=Ey-b*Ex ; a
最后擬合的結果是：a的最優化為 2.8387， b的最優化為 2.0968 ，公式 y = 2.8387* x + 2.0968
R語言的lm函數也可以計算a和b，完全不需要我們自己一個個手動推導。lm函數結果的Intercept即為b值，x所在行對應的Estimate值即為a值
回歸到我們的全基因組關聯分析，這里的a即為beta（OR）值
所以搞明白全基因組關聯分析的值是怎么來的了嗎，這個就是它的計算原理
其他變量呢
上面列出來的公式只是簡單的計算基因型與表型之間的相關性。
實際上，我們在計算的時候，會加入其他的變量，比如性別，年齡，品系等。
這些因素也是影響表型的變量。
因此，當考慮其他變量存在時，計算公式會稍微改變一下：y = ax + zβ + b
y:我們研究的表型
x:基因型數據，這里指每一個SNP
a：SNP的系數
z:年齡，性別等因素
β：年齡，性別等因素的系數
b：殘差，除了我們納入的SNP，性別年齡等因素外等其他可能影響表型的因素；
計算原理同上。