case/control的关联分析,本质是寻找在两组间基因型分布有差异的SNP位点,这些位点就是候选的关联信号,常用的分析方法有以下几种
- 卡方检验
- 费舍尔精确检验
- 逻辑回归
卡方检验是一种用途广泛的假设检验,属于非参数的检验一种,适合针对分类变量的分析。从形式上看,数据是由行和列对应的两个分类变量构成的表格,示意如下
对于case/control的关联分析,我们有两个分类变量,第一个就是样本的分组, 有case和control两组;第二个是Allel或者基因型的类别,对于Allele而言有两种,major和minor allele。对于基因型而言, 在上图中有AA, Aa, aa3种,当然在实际分析中,还会考虑遗传模型进一步对基因型的类别进行划分,常用的遗传模型有以下几种
- domanant model, 显性遗传模型,只要有突变位点就会致病,所以杂合突变和纯合突变归位一类,基因型就划分为两类,第一类为AA和Aa, 第二类为aa
- recessive model, 隐性模型, 只有纯合突变会致病,基因型同样划分为两类,第一类为纯合突变AA, 第二类为非纯合突变,Aa和aa
- additive model, 相加模型,突变位点的个数会影响性状的表型值,而且是累加关系,纯合突变的突变位点个数是杂合突变的2倍,对应的性状是不同的,基因型划分为3类, AA,Aa, aa
- multiplicative model, 相乘模型,突变位点的个数会影响性状的表型值,而且是相乘关系,纯合突变的突变位点个数是杂合突变的4倍,对应的性状是不同的,基因型划分为3类, AA,Aa, aa
以上模型根据划分的类别可以分为3大类,第一类是显性遗传模型,第二类是隐性遗传模型,第三类是additive, multiplicative model和常规的基因型分类,这三种模型都是划分为了3种基因型。
对于卡方检验,首先需要根据表格中的频数分布计算卡方统计量,公式如下
A表示实际频数,T表示理论频数,从公式可以看到,卡方统计量代表的是实际值与理论值之间的差异。看一个具体的例子
Genotype | AA | Aa | aa |
---|---|---|---|
Case | 30 | 15 | 55 |
Control | 28 | 12 | 60 |
上图表示的是两组实际观测到的基因型频数分布,对应的频率分布如下
Genotype | AA | Aa | aa |
---|---|---|---|
Case | 30% | 15% | 55% |
Control | 28% | 12% | 60% |
从数值上看,直观的可以看两组间分布有差异,但是这个差异是由抽样导致的误差还是真实存在的差异不知道。先假设两组间没有差异,合并样本,再次统计对应的频率,分别为29%, 13.5%,57.5% ,这3个数值就是理论频率, 根据这个频率来计算理论频数
Genotype | AA | Aa | aa |
---|---|---|---|
Case | 100 x 29% | 100 x 13.5% | 100 x 57.5% |
Control | 100 x 29% | 100 x 13.5 % | 100 x 57.5 % |
然后通过公式来计算卡方值,最终的计算结果为0.61969, 对应的R代码如下
从上图可以看到,对于卡方检验,除了卡方值X-squared之外,还有df和p-value两个值。df表示自由度,取值为(行数 - 1) X (列数 - 1), 上述数据为2X3的表格,自由度为2。为什么要考虑自由度呢?
这就要从卡方分布的定义说起,对于N个符合标准正态分布的变量,其平方和服从卡方分布,自由度指的就是这里的N, 不同自由度卡方分布是不同的,如下图所示
上图所示是不同自由度下卡方值的密度分布,不同自由度之间差别很大,所以我们需要先明确对应的自由度才可以利用卡方值来做出判断。利用自由度和卡方值,我们需要去查询卡方值分布表,获得对应的p值。在R中对应的操作代码如下
代码语言:javascript复制1 - pchisq(0.6196902, df = 2)
[1] 0.7335606
pchisq代表是卡方值的累计分布函数,代表卡方值小于0.6196902的概率。卡方分布表中为大于阈值的概率,示意如下
卡方值越小,对应的概率越大。自由度为2,P=0.05对应的卡方临界值为5.99, 上述示例的卡方值小于该临界值,说明发生的概率大于0.05,拒绝原假设,case/control组间差异不显著。
卡方检验虽然使用范围广泛,但还是有一些限制,样本量必须大于40, 而且最小的频数不能小于5, 这里的频数指的是理论频数
对于2X2的数据,当不满足要求时,推荐使用费舍尔精确检验来进行分析。
·end·