zoukankan      html  css  js  c++  java
  • 全基因组关联分析(GWAS)的计算原理

    前言

    关于全基因组关联分析(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之外的影响表型的因素

    来个例子给我们讲讲呗,公式怎么套进去

    Kntqd1.png

    如图所示,假定有一个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,性别年龄等因素外等其他可能影响表型的因素;

    计算原理同上。

  • 相关阅读:
    对QR码的初步研究(附:在博客里放上博客文章的QR码)
    EonerCMS——做一个仿桌面系统的CMS(十四)
    终于病了
    【HoorayOS】开源的Web桌面应用框架(第二版 v120311)
    【HoorayOS】开源的Web桌面应用框架——EonerCMS更名为HoorayOS
    一句代码实现 HTML5 语音搜索
    HTML5 拖拽上传图片实例
    【HoorayOS】开源的Web桌面应用框架
    【HoorayOS】开源的Web桌面应用框架(文件夹功能分析)
    从源码分析常见的基于Array的数据结构动态扩容机制
  • 原文地址:https://www.cnblogs.com/chenwenyan/p/11721698.html
Copyright © 2011-2022 走看看