研究生数学建模竞赛优秀论文基于假设检验与关联分析的多性状致病位点与致病基因方法研究24.docx
《研究生数学建模竞赛优秀论文基于假设检验与关联分析的多性状致病位点与致病基因方法研究24.docx》由会员分享，可在线阅读，更多相关《研究生数学建模竞赛优秀论文基于假设检验与关联分析的多性状致病位点与致病基因方法研究24.docx（26页珍藏版）》请在冰豆网上搜索。
研究生数学建模竞赛优秀论文基于假设检验与关联分析的多性状致病位点与致病基因方法研究24
1.1问题背景
1问题重述
人体的每条染色体携带一个DNA分子，人的遗传密码由人体中的DNA携带。
DNA是由分别带有A，T，C，G四种碱基的脱氧核苷酸链接组成的双螺旋长链分子。
在这条双螺旋的长链中，共有约30亿个碱基对，而基因则是DNA长链中有遗传效应的一些片段。
在组成DNA的数量浩瀚的碱基对（或对应的脱氧核苷酸）中，有一些特定位置的单个核苷酸经常发生变异引起DNA的多态性，我们称之为位点。
在DNA长链中，位点个数约为碱基对个数的1/1000。
由于位点在DNA长链中出现频繁，多态性丰富，近年来成为人们研究DNA遗传信息的重要载体，被称为人类研究遗传学的第三类遗传标记。
大量研究表明，人体的许多表型性状差异以及对药物和疾病的易感性等都可能与某些位点相关联，或和包含有多个位点的基因相关联。
因此，定位与性状或疾病相关联的位点在染色体或基因中的位置，能帮助研究人员了解性状和一些疾病的遗传机理，也能使人们对致病位点加以干预，防止一些遗传病的发生。
近年来，研究人员大都采用全基因组的方法来确定致病位点或致病基因，具体做法是：
招募大量志愿者（样本），包括具有某种遗传病的人和健康的人，通常用1表示病人，0表示健康者。
对每个样本，采用碱基（A，T，C，G）的编码方式来获取每个位点的信息（因为染色体具有双螺旋结构，所以用两个碱基的组合表示一个位点的信息）；例如本题中在位点rs100015位置，不同样本的编码都是T和C的组合，有三种不同编码方式TT，TC和CC。
类似地其他的位点虽然碱基的组合不同，但也只有三种不同编码。
研究人员可以通过对样本的健康状况和位点编码的对比分析来确定致病位点，从而发现遗传病或性状的遗传机理。
1.2问题分析
为了通过全基因组的方法来确定致病位点或致病基因，本题提出了下面四个问题，下面我们对问题进行重述，并简要叙述解决方案。
问题一、请用适当的方法，把genotype.dat中每个位点的碱基（A，T，C，G）编码方式转化成数值编码方式，便于进行数据分析。
解决方案：
通过对每个位点所在列进行观察，我们发现每一列只有三种编码方式，因此我们可以采用0/1/2的方式对每一个位点进行编码，编码时我们考虑到了次脱氧核苷酸（编码为0）和主脱氧核苷酸（编码为2）的问题，为后面的分析提供了方便。
问题二、根据附录中1000个样本在某条有可能致病的染色体片段上的9445个位点的编码信息（见genotype.dat）和样本患有遗传疾病A的信息（见phenotype.txt文件）。
设计或采用一个方法，找出某种疾病最有可能的一个或几个致病位点，并给出相关的理论依据。
解决方案：
这里我们通过构造假设检验的方式来找出可能的致病位点，特别的为了保证所得结果的可靠性，我们采用费舍尔精确检验、Cochran-Armitagetrend检验、logistics回归中β是否为0，这三种检验得到的p值进行分析，同时我们也考虑到了多重假设检验所带来的问题，对所得的p值进行了修正。
问题三、同上题中的样本患有遗传疾病A的信息（phenotype.txt文件）。
现有300个基因，每个基因所包含的位点名称见文件夹geneinfo中的300个dat文件，每个dat文件列出了对应基因所包含的位点（位点信息见文件genotype.dat）。
由于可以把基因理解为若干个位点组成的集合，遗传疾病与基因的关联性可以由基因中包含的位点的全集或其子集合表现出来请找出与疾病最有可能相关的一个或几个基因，并说明理由。
解决方案：
因为每个基因含有多个位点，为了充分利用每一个位点的信息，我们通过费舍尔方法，分别对第二问中三种假设检验得到的p值进行融合，通过融合后得到的检验统计量对每一个基因与疾病的关联性进行假设检验，找出可能的致病基因。
问题四、在问题二中，已知9445个位点，其编码信息见genotype.dat文件。
在实际的研究中，科研人员往往把相关的性状或疾病看成一个整体，然后来探寻与它们相关的位点或基因。
试根据multiphenos.txt文件给出的1000个样本的10个相关联性状的信息及其9445个位点的编码信息（见genotype.dat），找出与multiphenos.txt中10个性状有关联的位点。
解决方案：
我们使用了惩罚典型相关分析，并且通过最优缩放将位点的离散数值编码转换为连续变量，便于进行关联计算。
结果显示“rs351617”位点的致病可能性最大，“rs7538876”位点和“rs4949516”位点致病可能性次之，“rs780983”位点、“rs12746773”位点、“rs406985”位点和“rs716325”位点也有致病的可能性。
2.1符号说明
2
模型假设与符号说明
表1符号说明
I
样本数量
J
位点数量
D1
主脱氧核苷酸
D2
次脱氧核苷酸
sij
第i个样本第j个位点的数值编码
P-value
假定值
pj
位点j的假定值
X2
k
自由度为k服从X2的假设检验
gjk
第k个基因是否包含第j个位点
tj
第j个位点的自然底数值的二倍
χ2
m
自由度为m的卡方分布
σ2
第k个基因对应的位点协方差之和
Ω
t的协方差矩阵
ωT
i
第i个特征的权重向量
ρu,v
u和v的Pearson相关系数
Σij
第i个向量与第j个向量的协方差
λ
拉格朗日系数
θ
拉格朗日系数
aˆk
第k次迭代所估计的X的权重向量
ˆbk
第k次迭代所估计的Y的权重向量
x∗
转换后得到的连续型变量
J
转换函数
Gj
指示矩阵
cj
变量j的分类量化结果
3问题一模型建立与求解
DNA即脱氧核糖核酸，是一种生物大分子，可组成遗传指令，引导生物发育与生命机能运作。
主要功能是信息储存，可比喻为“蓝图”或“食谱”。
其中包含的指令，是建构细胞内其他的化合物，如蛋白质与核糖核酸所需。
带有蛋白质编码的DNA片段称为基因。
DNA是一种长链聚合物，组成单位为四种脱氧核苷酸，即：
腺嘌呤脱氧核苷酸（A）、胸腺嘧啶脱氧核苷酸（T）、胞嘧啶脱氧核苷酸（C）、鸟嘌呤脱氧核苷酸（G）。
四种碱基两两配对，即（A-T，C-G相互作用）形成DNA单体以及编码遗传信息的化学结构。
在组成DNA的数量浩瀚的碱基对（或对应的脱氧核苷酸）中，有一些特定位置的单个核苷酸经常发生变异引起DNA的多态性，我们称之为位点。
在DNA长链中，位点个数约为碱基对个数的1/1000。
由于位点在DNA长链中出现频繁，多态性丰富，近年来成为人们研究DNA遗传信息的重要载体，被称为人类研究遗传学的第三类遗传标记。
通过对1000个样本的位点信息的数据分析，我们可知，每个位点都包含两个碱基对，即可用两个脱氧核苷酸表示，并且其编码方式有三种，即两个碱基对的排列方式为2*2-1。
如位点“rs3094315”可用C与T的组合方式表示，其三种编码为CC、TT与TC。
通过对每一个位点1000个样本的两种脱氧核苷酸的数量进行统计，我们可得到该位点的主脱氧核苷酸（数量多的脱氧核苷酸）与次脱氧核苷酸（数量少的脱氧核苷酸）。
再如“rs3094315”位点，1000个样本中，TT有664个，TC有293个，CC有43个，因此T的数量为1621个，C的数量为379个，所以在该位点中T为主脱氧核苷酸，C为次脱氧核苷酸。
文献[1]中对基因型的编码方式，本文对位点碱基编码方式进行如下的数值编码方式的转换。
≫
假定数据为从I个个体中采集的J个位点的信息，并且J远远大于I，即JI，本文中，I=1000，J=9445。
假设第j个位点包含的两种碱基分别为D1，D2，其中D1为主脱氧核苷酸，D2为次脱氧核苷酸。
则第i个个体的第j个位点的数值sij如公式
（1）所示
sij=
0,ifD2D2
1,ifD1D2orD2D1
2,ifD1D1.
（1）
也即如果该位点包含两个主脱氧核苷酸则该位点编码为2，如果该位点为杂合的（包含一个主脱氧核苷酸和一个次脱氧核苷酸），该位点编码为1，如果该位点包含两个次脱氧核苷酸，该位点编码为0。
因此所有位点均为0，1，2三种数值码。
一段数值编码的示意图如图1所示。
图1数值编码示意图
图2多重检验各类情况图列联表
4问题二模型建立与求解
我们收集到了1000个样本中患A病的情况，其中500人为A病患组，另500人为对照组。
本题要根据1000个样本的患病情况以及与其对应的9445个位点信息确定一个或几个可能致病的位点。
为解决此问题，本文采用了几种基于假设检验p值的方法，这里首先介绍一些假设检验的基本知识。
4.1假设检验概述
−−
−
设有M个假设检验，H1,H2,...,HM，相应的p值为p1,p2,...,pM，通常情况下，如果p值小于0.05，我们称这个假设检验是显著的。
然而此处我们的问题是多重假设检验问题以0.05作为p值的比较阈值是不合适的。
具体原因如下：
以0.05为阈值时，对于M个假设检验，所有检验都正确的概率是（10.05）M，即我们至少犯一个错误的概率为1（10.05）M，当M变大时该概率会很大，在我们的问题中M为9445，该值近似为1。
另一方面，由大数定律可知，在M=9445阈值为0.05条件下，平均意义下应该有9445∗0.05≈472个检验被拒绝。
{
因此我们应对p值进行修正，一种方法采用比较常用的方法的Bonferroni修正，即将阈值改为0.05/M=0.05/9445≈5.3×10−6.另一种是控制错误发现率FDR,
E（V/R）,ifR>00,ifR=0.
（2）
其中V,R分别为M个假设检验中被错误拒绝的检验的个数和拒绝的总个数，如图2所示：
我们可以通过Benjamini-Hochberg过程来达到控制FDR的目的：
•固定错误发现率，如α=0.05,对p值排序p
（1）≤p
（2）≤...≤p（M）
•
jL=max{j:
p（j）<αM}
•拒绝所有满足p（j）≤p（L）的Hj.
通过上面两种修正，我们便可以解决多重假设检验的问题。
下面我们将分别采用费舍尔精确检验、Cochran-Armitagetrend检验、logistics回归中β是否为0，这三种检验得到的p值进行分析。
4.2费舍尔精确检验
费舍尔精确检验[2]是一种通过分析列联表来判断统计显著性的检验方法。
尽管经常在样本数量较小的情况下应用，但它在任意样本数量下都是有效的。
它是以它的发明者RonaldFisher命名的，是精确检验方法的一种。
之所以被称为精确检验是因为原
表2某一统计量的示例列联表
Class1
Class2
RowTotal
Case
a
b
a+b
Control
c
d
c+d
ColumnTotal
a+c
b+d
a+b+c+d=n
表3
升级会员 