CN113066530B - 一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法 - Google Patents
一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法 Download PDFInfo
- Publication number
- CN113066530B CN113066530B CN202110346625.3A CN202110346625A CN113066530B CN 113066530 B CN113066530 B CN 113066530B CN 202110346625 A CN202110346625 A CN 202110346625A CN 113066530 B CN113066530 B CN 113066530B
- Authority
- CN
- China
- Prior art keywords
- snp
- gene
- cis
- trans
- txt
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 46
- 238000000034 method Methods 0.000 title claims abstract description 25
- 108090000623 proteins and genes Proteins 0.000 claims abstract description 83
- 210000000349 chromosome Anatomy 0.000 claims description 52
- 240000002791 Brassica napus Species 0.000 claims description 12
- 235000011293 Brassica napus Nutrition 0.000 claims description 12
- 108020004999 messenger RNA Proteins 0.000 claims description 8
- 230000000717 retained effect Effects 0.000 claims description 6
- 108700028369 Alleles Proteins 0.000 claims description 5
- 241000196324 Embryophyta Species 0.000 claims description 3
- 238000011144 upstream manufacturing Methods 0.000 claims description 3
- 238000012217 deletion Methods 0.000 claims 1
- 230000037430 deletion Effects 0.000 claims 1
- ZVQOOHYFBIDMTQ-UHFFFAOYSA-N [methyl(oxido){1-[6-(trifluoromethyl)pyridin-3-yl]ethyl}-lambda(6)-sulfanylidene]cyanamide Chemical compound N#CN=S(C)(=O)C(C)C1=CC=C(C(F)(F)F)N=C1 ZVQOOHYFBIDMTQ-UHFFFAOYSA-N 0.000 abstract 1
- 238000004364 calculation method Methods 0.000 abstract 1
- 230000014509 gene expression Effects 0.000 description 8
- 241000894007 species Species 0.000 description 5
- 238000012098 association analyses Methods 0.000 description 3
- 230000002759 chromosomal effect Effects 0.000 description 3
- 201000010099 disease Diseases 0.000 description 2
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 2
- 230000002068 genetic effect Effects 0.000 description 2
- 239000000463 material Substances 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 101150072531 10 gene Proteins 0.000 description 1
- 241000219194 Arabidopsis Species 0.000 description 1
- 206010071602 Genetic polymorphism Diseases 0.000 description 1
- 240000007594 Oryza sativa Species 0.000 description 1
- 235000007164 Oryza sativa Nutrition 0.000 description 1
- 240000008042 Zea mays Species 0.000 description 1
- 235000016383 Zea mays subsp huehuetenangensis Nutrition 0.000 description 1
- 235000002017 Zea mays subsp mays Nutrition 0.000 description 1
- 238000009412 basement excavation Methods 0.000 description 1
- -1 chromosome Proteins 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000010219 correlation analysis Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 230000007614 genetic variation Effects 0.000 description 1
- 238000003205 genotyping method Methods 0.000 description 1
- 235000009973 maize Nutrition 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 239000003550 marker Substances 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 230000009456 molecular mechanism Effects 0.000 description 1
- 239000002773 nucleotide Substances 0.000 description 1
- 125000003729 nucleotide group Chemical group 0.000 description 1
- 238000011176 pooling Methods 0.000 description 1
- 230000001105 regulatory effect Effects 0.000 description 1
- 235000009566 rice Nutrition 0.000 description 1
- 238000012163 sequencing technique Methods 0.000 description 1
- 238000007619 statistical method Methods 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
- 238000012070 whole genome sequencing analysis Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
- G16B20/20—Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B50/00—ICT programming tools or database systems specially adapted for bioinformatics
Landscapes
- Health & Medical Sciences (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Evolutionary Biology (AREA)
- Medical Informatics (AREA)
- Biophysics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biotechnology (AREA)
- Theoretical Computer Science (AREA)
- General Health & Medical Sciences (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Databases & Information Systems (AREA)
- Bioethics (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Genetics & Genomics (AREA)
- Molecular Biology (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明公开了一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法。本发明所提供的合并eQTL分析结果中存在连锁不平衡SNP的方法基于SNP与靶基因的位置信息将其分为顺式和反式SNP，基于SNP位置信息合并相邻的SNP为一个SNP cluster，基于SNP cluster间的连锁不平衡程度进一步合并结果。本发明的脚本由python3语言写成，速度快，灵活性高，可靠性强，实现了批量化、自动化和流程化计算。本发明将在eQTL分析结果的简化上发挥重要作用。
Description
技术领域
本发明属于生物技术领域，涉及一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法。
背景技术
关联分析(association analysis)是一种以位点间的连锁不平衡(linkagedisequilibrium)为基础，采用统计方法检测遗传多态性与性状之间关联的分析方法。全基因组关联分析(genome-wide association study，GWAS)最初应用于人类疾病相关的研究中，这些研究对人类理解相关疾病的遗传基础和分子机制具有显著的贡献。近年来，随着高密度SNP基因分型芯片和全基因组测序等技术的发展，GWAS广泛应用于作物复杂性状遗传结构的解析。与传统的连锁分析相比，GWAS有多个优势：可以利用自然群体，无需针对特定性状构建作图群体，花费时间少；能同时检测群体同一基因的多个等位基因，有利于优良等位基因的挖掘；当群体连锁不平衡程度低并且标记覆盖度高时，定位精度高，可以达到单基因水平。
eQTL(expression quantitative trait loci)分析是将每个基因的表达量作为表型，将它们与基因组变异位点进行关联分析，从而发掘调控基因表达的遗传变异的分析方法。通过eQTL分析，不但可以获取基因的表达调控位点，还可以构建基因之间的调控网络，对于解析基因的调控机制具有很好的指导作用。目前，eQTL分析已经成功地应用到多个物种中，包括拟南芥、玉米和水稻等。根据eQTL位点与调控的目标基因的位置关系可以将eQTL分成顺式eQTL(cis eQTL)和反式eQTL(trans eQTL)。顺式eQTL与目标基因物理距离较近，表明是该基因本身的差别引起mRNA水平的差别，反式eQTL与目标基因物理距离较近或位于不同的染色体，表明是其他基因的差别控制该基因mRNA水平的差异。
做eQTL分析前首先要对基因的表达谱进行正态化，其次利用正态化的表达谱计算隐性因素，利用基因型数据计算群体结构，最终以群体结构和隐性因素作为协变量，对群体的基因表达谱和基因型数据进行关联分析。由于得到的eQTL结果非常多，常常数以十万、百万计，而相邻的SNP间往往连锁不平衡(linkage disequilibrium，LD)程度高，因此需要合并冗余SNP简化结果。手动合并存在LD的SNP非常耗时耗力，但目前没有软件能提供批量合并的功能，这成为一个亟待解决的问题。
发明内容
本发明提供的是一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法，具体包括如下步骤：
(1)在Windows操作系统下创建工作目录eqtl_analysis及其子文件夹gene_info，将待分析植物基因注释信息文件***.gff3和脚本abstract_gene_info.py放在gene_info文件夹下，运行“python abstract_gene_info.py***.gff3”命令，得到每条染色体各自的基因信息文件，记为G数据集。
所述“***.gff3”为研究物种的基因注释文件，油菜Darmor-bzh参考基因组对应的文件为Brassica_napus.annotation_v5.gff3，油菜中双11参考基因组对应的文件为ZS11.annotation.gff3。
G数据集文件命名方式为染色体名+“_gene_info.txt”，不保留标题行，以油菜为例，染色体A01基因型文件为A01_gene_info.txt。文件包括5列，分别为基因名、染色体、基因起始位置、基因中止位置和正负链信息(图1)。
(2)将待分析植物的eQTL结果文件记为A数据集，格式为eQTL分析常用软件MatrixeQTL的结果文件格式；脚本separate_cis_trans.py参考G数据集提供的基因物理位置，分析SNP与基因的染色体和物理距离，将所有SNP划分为两类，顺式SNP和反式SNP。A数据集和脚本eqtl_cis_trans.py均放在工作目录eqtl_analysis下，运行命令“pythonseparate_cis_trans.py XXX1.txt cis_dis”，得到“XXX1_cis.txt”和“XXX1_trans.txt”两个文件。
“XXX1.txt”代表所述A数据集的文件名，格式与eQTL分析常用软件MatrixeQTL的结果文件格式相同，包含6列“SNP”，“gene”，“beta”，“t-stat”，“p-value”和“FDR”(图2)，脚本利用其中“SNP”，“gene”和“p-value”三列信息。文件需按“gene”和“SNP”两列信息进行排序。SNP的命名方式用染色体名+物理位置，染色体为3位或者10位，如A03和C03_random；物理位置为8位数，不足的位数用0补全。以油菜为例，SNP的染色体为A03，物理位置为41437，SNP命名为A0300041437；SNP的染色体为C03_random，物理位置为44754，SNP命名为A01_random00044754。
“cis_dis”为划分SNP为顺式SNP和反式SNP的距离阈值，默认设定为24,000bp。
所述“XXX1_cis.txt”为包含所有顺式SNP的文件名，记为B1数据集；“XXX1_trans.txt”为包含所有反式SNP的文件名，记为B2数据集。
(3)脚本combine_near_snp.py合并相邻的显著SNP，得到SNP cluster，并用其中最显著、物理位置小的SNP作为代表。将脚本combine_near_snp.py放在工作目录eqtl_analysis下，针对B1、B2数据集分别运行命令“python combine_near_snp.py XXX1_cis.txt part_dis”和“python combine_near_snp.py XXX1_trans.txt part_dis”，得到“XXX1_cis_median.txt”和“XXX1_trans_median.txt”两个文件。
“part_dis”为合并相邻SNP的距离阈值，默认设定为10,000bp。所述“XXX1_cis_median.txt”为合并相邻的顺式SNP后得到的结果文件，记为C1数据集；“XXX1_trans_median.txt”为合并相邻的反式SNP后得到的结果文件，记为C2数据集。
(4)为计算SNP cluster间的LD系数r2，手动创建各染色体的关联群体基因型文件，记为M数据集，放在新创建的eqtl_analysis子目录snp_info下。脚本combine_ld_snp.py参考M数据集，计算SNP cluster间的r2，若r2>0.1，则保留更显著、物理位置更小的SNP cluster。将脚本combine_ld_snp.py放在工作目录eqtl_analysis下，针对C1、C2数据集分别运行命令“python combine_ld_snp.py XXX1_cis_median.txt”和“pythoncombine_ld_snp.py XXX1_trans_median.txt”，得到“XXX1_cis_final.txt”和“XXX1_trans_final.txt”两个文件。
M数据集各染色体基因型文件命名方式为染色体名+“_snp_info.txt”，不保留标题行，以油菜为例，染色体A03基因型文件为A03_snp_info.txt。文件格式为SNP+基因型信息，SNP只能包含两个等位基因，分别用0和2表示，杂合、缺失用NA表示(图3)。
所述“XXX1_cis_final.txt”为合并所有相邻、存在LD的顺式SNP的最终结果文件，记为D1数据集；所述“XXX1_trans_final.txt”为合并所有相邻、存在LD的反式SNP的最终结果文件，记为D2数据集。
在上述方法步骤(1)中，所述脚本abstract_gene_info.py关于提取基因信息的内容，是基于如下原理进行编程的：提取基因注释文件中第三列为“mRNA”的所有行，提取第九列“Parent”后内容作为基因名，第一列内容作为染色体信息，第四、五和七列作为起始、中止和正负链信息。
在上述方法步骤(2)中，所述脚本separate_cis_trans.py关于划分SNP为顺式和反式的内容，是基于如下原理进行编程的：根据G数据集，提取基因染色体和物理位置；若SNP与基因位于不同染色体，划分为反式SNP；SNP与基因位于相同染色体时，若SNP位于基因上下游24kb内，划分为顺式SNP，否则划分为反式SNP。
在上述方法步骤(3)中，所述脚本combine_near_snp.py关于合并相邻SNP为SNPcluster的内容，是基于如下原理进行编程的：判断同一染色体物理位置最大和最小的SNP距离是否超过10kb；距离小于10kb时，合并所有SNP为一个cluster，用cluster中最显著、物理位置最小的SNP作为代表；距离大于10kb时，以10kb为间距，将SNP分为多组，每组分别进行合并，得到多个cluster。
在上述方法步骤(4)中，所述脚本combine_ld_snp.py关于合并SNP cluster的内容，是基于如下原理进行编程的：参考M数据集，计算SNP cluster间的LD系数r2；若r2>0.1，则保留更显著、物理位置更小的SNP cluster，否则两者均保留。
在本发明的一个实施例中，所述脚本abstract_gene_info.py具体为：
在本发明的一个实施例中，所述脚本separate_cis_trans.py具体为：
在本发明的一个实施例中，所述脚本combine_near_snp.py具体为：
/>
/>
在本发明的一个实施例中，所述脚本combine_ld_snp.py具体为：
/>
/>
/>
/>
/>
/>
/>
在所述方法中，所述基因注释文件可以通过下载已公开的物种基因注释文件获得。
具体地，本发明所述物种具体为甘蓝型油菜，所述基因注释文件记录于法国GenoScope网站(https://www.genoscope.cns.fr/brassicanapus/)。
本发明的优点在于：1.可批量合并eQTL分析结果中存在连锁不平衡SNP，策略为先合并相邻位点再合并存在LD的位点，速度快，效率高；2.脚本可以自定义合并相邻SNP的区间大小和划分顺式/反式SNP的距离阈值，灵活性高；3.当输入文件很大时，脚本需要按分段读取文件，会造成部分跨段且存在LD的SNP无法合并，导致结果的冗余。本方法考虑到该问题并在脚本中添加了相应代码予以解决，确保分析结果的代表性。
附图说明
图1为G数据集文件格式，包括5列，分别为基因名、染色体、基因起始位置、基因中止位置和正负链信息。
图2为A数据集文件格式，包含6列“SNP”，“gene”，“beta”，“t-stat”，“p-value”和“FDR”。
图3为M数据集文件格式，包含SNP名和基因型信息。
图4为本发明批量合并eQTL分析结果中存在连锁不平衡SNP的方法的流程图。
具体实施方式
下述实施例中所使用的实验方法如无特殊说明，均为常规方法。
下述实施例中所用的材料，如无特殊说明，均可从商业途径得到。
实施例1、批量合并eQTL分析结果中存在连锁不平衡SNP方法的建立
本发明所提供的批量合并eQTL分析结果中存在连锁不平衡SNP方法的流程图见图4，具体包括如下步骤：
(1)在Windows操作系统下创建工作目录eqtl_analysis及其子文件夹gene_info，将基因注释信息文件***.gff3和脚本abstract_gene_info.py放在gene_info文件夹下，运行“python abstract_gene_info.py***.gff3”命令，得到每条染色体各自的基因信息文件，记为G数据集。
所述“***.gff3”为研究物种的基因注释文件，油菜Darmor-bzh参考基因组对应的文件为Brassica_napus.annotation_v5.gff3，油菜中双11参考基因组对应的文件为ZS11.annotation.gff3。
G数据集文件命名方式为染色体名+“_gene_info.txt”，不保留标题行，以油菜为例，染色体A01基因型文件为A01_gene_info.txt。文件包括5列，分别为基因名、染色体、基因起始位置、基因中止位置和正负链信息。其中，所属脚本abstract_gene_info.py具有如下特点：提取基因注释文件中第三列为“mRNA”的所有行，提取第九列“Parent”后内容作为基因名，第一列内容作为染色体信息，第四、五和七列作为起始、中止和正负链信息。
/>
/>
(2)将待分析的eQTL结果文件记为A数据集，格式为eQTL分析常用软件MatrixeQTL的结果文件格式；脚本separate_cis_trans.py参考G数据集提供的基因物理位置，分析SNP与基因的染色体和物理距离，将所有SNP划分为两类，顺式SNP和反式SNP。A数据集和脚本eqtl_cis_trans.py均放在工作目录eqtl_analysis下，运行命令“python separate_cis_trans.py XXX1.txt cis_dis”，得到“XXX1_cis.txt”和“XXX1_trans.txt”两个文件。所述“XXX1.txt”代表所述A数据集的文件名，格式与eQTL分析常用软件MatrixeQTL的结果文件格式相同，包含6列“SNP”，“gene”，“beta”，“t-stat”，“p-value”和“FDR”，脚本利用其中“SNP”，“gene”和“p-value”三列信息。文件需按“gene”和“SNP”两列信息进行排序。SNP的命名方式用染色体名+物理位置，染色体为3位或者10位，如A03和C03_random；物理位置为8位数，不足的位数用0补全。以油菜为例，SNP的染色体为A03，物理位置为41437，SNP命名为A0300041437；SNP的染色体为C03_random，物理位置为44754，SNP命名为A01_random00044754。所述“cis_dis”为划分SNP为顺式SNP和反式SNP的距离阈值，默认设定为24,000bp。
所述“XXX1_cis.txt”为包含所有顺式SNP的文件名，记为B1数据集；“XXX1_trans.txt”为包含所有反式SNP的文件名，记为B2数据集。其中，所述脚本separate_cis_trans.py具有如下特点：根据G数据集，提取基因染色体和物理位置；若SNP与基因位于不同染色体，划分为反式SNP；SNP与基因位于相同染色体时，若SNP位于基因上下游24kb内，划分为顺式SNP，否则划分为反式SNP。
/>
/>
/>
(3)脚本combine_near_snp.py合并相邻的显著SNP，得到SNP cluster，并用其中最显著、物理位置小的SNP作为代表。将脚本combine_near_snp.py放在工作目录eqtl_analysis下，针对B1、B2数据集分别运行命令“python combine_near_snp.py XXX1_cis.txt part_dis”和“python combine_near_snp.py XXX1_trans.txt part_dis”，得到“XXX1_cis_median.txt”和“XXX1_trans_median.txt”两个文件。
所述“part_dis”为合并相邻SNP的距离阈值，默认设定为10,000bp。
所述“XXX1_cis_median.txt”为合并相邻的顺式SNP后得到的结果文件，记为C1数据集；“XXX1_trans_median.txt”为合并相邻的反式SNP后得到的结果文件，记为C2数据集。
其中，所述脚本combine_near_snp.py具有如下特点：判断同一染色体物理位置最大和最小的SNP距离是否超过10kb；距离小于10kb时，合并所有SNP为一个cluster，用cluster中最显著、物理位置最小的SNP作为代表；距离大于10kb时，以10kb为间距，将SNP分为多组，每组分别进行合并，得到多个cluster。
/>
/>
/>
/>
(4)为计算SNP cluster间的LD系数r2，手动创建各染色体的关联群体基因型文件，记为M数据集，放在新创建的eqtl_analysis子目录snp_info下。脚本combine_ld_snp.py参考M数据集，计算SNP cluster间的r2，若r2>0.1，则保留更显著、物理位置更小的SNP cluster。将脚本combine_ld_snp.py放在工作目录eqtl_analysis下，针对C1、C2数据集分别运行命令“python combine_ld_snp.py XXX1_cis_median.txt”和“pythoncombine_ld_snp.py XXX1_trans_median.txt”，得到“XXX1_cis_final.txt”和“XXX1_trans_final.txt”两个文件。
M数据集各染色体基因型文件命名方式为染色体名+“_snp_info.txt”，不保留标题行，以油菜为例，染色体A03基因型文件为A03_snp_info.txt。文件格式为SNP+基因型信息，SNP只能包含两个等位基因，分别用0和2表示，杂合、缺失用NA表示。
所述“XXX1_cis_final.txt”为合并所有相邻、存在LD的顺式SNP的最终结果文件，记为D1数据集；所述“XXX1_trans_final.txt”为合并所有相邻、存在LD的反式SNP的最终结果文件，记为D2数据集。
其中，所述脚本combine_ld_snp.py具有如下特点：参考M数据集，计算SNPcluster间的LD系数r2；若r2>0.1，则保留更显著、物理位置更小的SNP cluster，否则两者均保留。
/>
/>
/>
/>
/>
/>
/>
实施例2、利用实施例1建立的方法批量合并甘蓝型油菜eQTL分析结果中存在连锁不平衡SNP
进入GenoScope网站下载甘蓝型油菜Darmor-bzh的基因注释文件Brassica_napus.annotation_v5.gff3。我们之前对29个甘蓝型油菜的分枝进行转录组测序，获得基因的表达谱数据，结合已获得的基因型数据，利用MatrixeQTL软件进行eQTL分析，结果文件包含8,489个基因及其5，616，239个显著SNP，命名为“29_all_snp.txt”(A数据集)。利用本发明中的脚本对该文件进行SNP的合并。分。具体包括如下步骤：
1)参照实施例1的步骤(1)进行。
采用abstract_gene_info.py脚本从甘蓝型油菜基因注释文件Brassica_napus.annotation_v5.gff3中提取每条染色体的基因信息，分别生成文件(G数据集)保存在gene_info子文件夹中。
2)参照实施例1的步骤(2)进行。
利用separate_cis_trans.py脚本对A数据集进行分析，参考G数据集提供的基因物理位置，分析SNP与基因的染色体和物理距离，将所有SNP划分为两类，顺式SNP(B1数据集)和反式SNP(B2数据集)。最终检测到203,868顺式SNP和5,412,371反式SNP。
3)参照实施例1的步骤(3)进行。
采用combine_near_snp.py脚本分别对B1和B2数据集进行相邻SNP(<10kb)的合并，得到6,324顺式(C1数据集)和490,971反式SNP cluster(C2数据集)。
4)参照实施例1的步骤(4)进行。
提取C1和C2数据集所有SNP的名称，从29份材料的基因型文件(vcf格式)中提取这些SNP的基因型信息，按说明书部分的描述整理成相应的格式，存放在snp_info子文件夹中(M数据集)。脚本combine_near_snp.py根据M数据集信息，分别计算C1和C2数据集中同一靶基因不同SNP cluster间的LD系数r2，进一步合并存在LD的SNP cluster，最终得到2,674顺式(D1数据集)和42,631反式SNP cluster(D2数据集)。
随机挑选10个基因SNP的eQTL结果进行验证，与脚本分析得到的结果完全相同，这证实了本方法的可靠性。
Claims (6)
1.一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法，其特征在于，步骤如下：
(1)在Windows操作系统下创建工作目录eqtl_analysis及其子文件夹gene_info，将待分析植物的基因注释信息文件***.gff3和脚本abstract_gene_info.py放在gene_info文件夹下，运行“python abstract_gene_info.py***.gff3”命令，得到每条染色体各自的基因信息文件，记为G数据集；
G数据集文件命名方式为染色体名+“_gene_info.txt”，不保留标题行，文件包括5列，分别为基因名、染色体、基因起始位置、基因中止位置和正负链信息；
(2)将待分析植物的eQTL结果文件记为A数据集，格式为eQTL分析常用软件MatrixeQTL的结果文件格式；脚本separate_cis_trans.py参考G数据集提供的基因物理位置，分析SNP与基因的染色体和物理距离，将所有SNP划分为两类，顺式SNP和反式SNP；A数据集和脚本eqtl_cis_trans.py均放在工作目录eqtl_analysis下，运行命令“python separate_cis_trans.py XXX1.txt cis_dis”，得到“XXX1_cis.txt”和“XXX1_trans.txt”两个文件；
“XXX1.txt”代表所述A数据集的文件名，格式与eQTL分析常用软件MatrixeQTL的结果文件格式相同，包含6列“SNP”，“gene”，“beta”，“t-stat”，“p-value”和“FDR”，脚本利用其中“SNP”，“gene”和“p-value”三列信息；文件按“gene”和“SNP”两列信息进行排序；SNP的命名方式用染色体名+物理位置，染色体为3位或者10位，物理位置为8位数，不足的位数用0补全；
“cis_dis”为划分SNP为顺式SNP和反式SNP的距离阈值，默认设定为24,000bp；
所述“XXX1_cis.txt”为包含所有顺式SNP的文件名，记为B1数据集；所述“XXX1_trans.txt”为包含所有反式SNP的文件名，记为B2数据集；
(3)脚本combine_near_snp.py合并相邻的显著SNP，得到SNP cluster，并用其中最显著、物理位置小的SNP作为代表，将脚本combine_near_snp.py放在工作目录eqtl_analysis下，针对B1、B2数据集分别运行命令“python combine_near_snp.py XXX1_cis.txt part_dis”和“python combine_near_snp.py XXX1_trans.txt part_dis”，得到“XXX1_cis_median.txt”和“XXX1_trans_median.txt”两个文件；
“part_dis”为合并相邻SNP的距离阈值，默认设定为10,000bp；所述“XXX1_cis_median.txt”为合并相邻的顺式SNP后得到的结果文件，记为C1数据集；“XXX1_trans_median.txt”为合并相邻的反式SNP后得到的结果文件，记为C2数据集；
(4)为计算SNP cluster间的LD系数r2，手动创建各染色体的关联群体基因型文件，记为M数据集，放在新创建的eqtl_analysis子目录snp_info下；脚本combine_ld_snp.py参考M数据集，计算SNP cluster间的r2，若r2>0.1，则保留更显著、物理位置更小的SNP cluster；将脚本combine_ld_snp.py放在工作目录eqtl_analysis下，针对C1、C2数据集分别运行命令“python combine_ld_snp.py XXX1_cis_median.txt”和“python combine_ld_snp.pyXXX1_trans_median.txt”，得到“XXX1_cis_final.txt”和“XXX1_trans_final.txt”两个文件；
M数据集各染色体基因型文件命名方式为染色体名+“_snp_info.txt”，不保留标题行，文件格式为SNP+基因型信息，SNP包含两个等位基因，分别用0和2表示，杂合、缺失用NA表示；
所述“XXX1_cis_final.txt”为合并所有相邻、存在LD的顺式SNP的最终结果文件，记为D1数据集；所述“XXX1_trans_final.txt”为合并所有相邻、存在LD的反式SNP的最终结果文件，记为D2数据集。
2.根据权利要求1所述的方法，其特征在于：步骤(1)中，所述脚本abstract_gene_info.py关于提取基因信息的内容，是基于如下原理进行编程的：提取基因注释文件中第三列为“mRNA”的所有行，提取第九列“Parent”后内容作为基因名，第一列内容作为染色体信息，第四、五和七列作为起始、中止和正负链信息。
3.根据权利要求1所述的方法，其特征在于：步骤(2)中，所述脚本separate_cis_trans.py关于划分SNP为顺式和反式的内容，是基于如下原理进行编程的：根据G数据集，提取基因染色体和物理位置；若SNP与基因位于不同染色体，划分为反式SNP；SNP与基因位于相同染色体时，若SNP位于基因上下游24kb内，划分为顺式SNP，否则划分为反式SNP。
4.根据权利要求1所述的方法，其特征在于：步骤(3)中，所述脚本combine_near_snp.py关于合并相邻SNP为SNP cluster的内容，是基于如下原理进行编程的：判断同一染色体物理位置最大和最小的SNP距离是否超过10kb；距离小于10kb时，合并所有SNP为一个cluster，用cluster中最显著、物理位置最小的SNP作为代表；距离大于10kb时，以10kb为间距，将SNP分为多组，每组分别进行合并，得到多个cluster。
5.根据权利要求1中所述的方法，其特征在于：步骤(4)中，所述脚本combine_ld_snp.py关于合并SNP cluster的内容，是基于如下原理进行编程的：参考M数据集，计算SNPcluster间的LD系数r2；若r2>0.1，则保留更显著、物理位置更小的SNP cluster，否则两者均保留。
6.根据权利要求1-5中任一所述的方法，其特征在于：所述待分析植物为甘蓝型油菜。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110346625.3A CN113066530B (zh) | 2021-03-31 | 2021-03-31 | 一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110346625.3A CN113066530B (zh) | 2021-03-31 | 2021-03-31 | 一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113066530A CN113066530A (zh) | 2021-07-02 |
| CN113066530B true CN113066530B (zh) | 2024-05-10 |
Family
ID=76565091
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110346625.3A Active CN113066530B (zh) | 2021-03-31 | 2021-03-31 | 一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113066530B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN116564410A (zh) * | 2023-05-23 | 2023-08-08 | 浙江大学 | 一种预测突变位点顺式调控基因的方法、设备和介质 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102369531A (zh) * | 2009-02-06 | 2012-03-07 | 先正达参股股份有限公司 | 用于选择统计上确认的候选基因的方法 |
| CN110349623A (zh) * | 2019-01-17 | 2019-10-18 | 哈尔滨工业大学 | 基于改进孟德尔随机化的老年痴呆病基因及位点筛选方法 |
| CN110534157A (zh) * | 2019-07-26 | 2019-12-03 | 江苏省农业科学院 | 一种批量提取基因组基因信息并翻译比对分析序列的方法 |
Family Cites Families (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20200051661A1 (en) * | 2016-10-18 | 2020-02-13 | Arizona Board Of Regents On Behalf Of The University Of Arizona | Pharmacogenomics of Intergenic Single-Nucleotide Polymorphisms and in Silico Modeling for Precision Therapy |
-
2021
- 2021-03-31 CN CN202110346625.3A patent/CN113066530B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102369531A (zh) * | 2009-02-06 | 2012-03-07 | 先正达参股股份有限公司 | 用于选择统计上确认的候选基因的方法 |
| CN110349623A (zh) * | 2019-01-17 | 2019-10-18 | 哈尔滨工业大学 | 基于改进孟德尔随机化的老年痴呆病基因及位点筛选方法 |
| CN110534157A (zh) * | 2019-07-26 | 2019-12-03 | 江苏省农业科学院 | 一种批量提取基因组基因信息并翻译比对分析序列的方法 |
Non-Patent Citations (4)
| Title |
|---|
| 分析SNP位点：里阿奴搜不平衡-可视化R包LDheatmap;程凉皮儿;《https://www.jianshu.com/p/10c9c26cfd57》;20200110;全文 * |
| 基于eQTL数据的生物网络构建与分析方法;彭启迪;《中国优秀硕士学位论文全文数据库基础科学辑》;20210228;A006-567 * |
| 孙程明 ; 陈锋 ; 陈松 ; 彭琦 ; 张维 ; 易斌 ; 张洁夫 ; 傅廷栋.甘蓝型油菜每角粒数的全基因组关联分析.《作物学报》.2019,第147-153页. * |
| 生物信息学与eQTL作图:问题和展望;李宏;;生物技术通报;20120726(07);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113066530A (zh) | 2021-07-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Jayakodi et al. | The barley pan-genome reveals the hidden legacy of mutation breeding | |
| Ling et al. | Genome sequence of the progenitor of wheat A subgenome Triticum urartu | |
| Daccord et al. | High-quality de novo assembly of the apple genome and methylome dynamics of early fruit development | |
| Allen et al. | Characterization of a Wheat Breeders’ Array suitable for high‐throughput SNP genotyping of global accessions of hexaploid bread wheat (Triticum aestivum) | |
| Grattapaglia et al. | High-throughput SNP genotyping in the highly heterozygous genome of Eucalyptus: assay success, polymorphism and transferability across species | |
| Eckert et al. | Association genetics of coastal Douglas fir (Pseudotsuga menziesii var. menziesii, Pinaceae). I. Cold-hardiness related traits | |
| Mascher et al. | Anchoring and ordering NGS contig assemblies by population sequencing (POPSEQ) | |
| Choi et al. | A soybean transcript map: gene distribution, haplotype and single-nucleotide polymorphism analysis | |
| Hyten et al. | A high density integrated genetic linkage map of soybean and the development of a 1536 universal soy linkage panel for quantitative trait locus mapping | |
| Zhou et al. | Targeted enrichment of the black cottonwood (Populus trichocarpa) gene space using sequence capture | |
| Paschold et al. | Nonsyntenic genes drive highly dynamic complementation of gene expression in maize hybrids | |
| Orr et al. | A phylogenomic approach reveals a low somatic mutation rate in a long-lived plant | |
| Velmurugan et al. | An ultra-high density genetic linkage map of perennial ryegrass (Lolium perenne) using genotyping by sequencing (GBS) based on a reference shotgun genome assembly | |
| Mousavi et al. | De novo SNP discovery and genetic linkage mapping in poplar using restriction site associated DNA and whole-genome sequencing technologies | |
| Liu et al. | Development of genome-wide insertion and deletion markers for maize, based on next-generation sequencing data | |
| Bradeen et al. | Genetics, genomics and breeding of potato | |
| Cottage et al. | Heterozygosity and diversity analysis using mapped single nucelotide polymorphisms in a faba bean inbreeding programme | |
| Shirasawa et al. | Target amplicon sequencing for genotyping genome‐wide single nucleotide polymorphisms identified by whole‐genome resequencing in peanut | |
| CN113066530B (zh) | 一种批量合并eQTL分析结果中存在连锁不平衡SNP的方法 | |
| Li et al. | Fine mapping of the sex locus in Salix triandra confirms a consistent sex determination mechanism in genus Salix | |
| Gallavotti et al. | Positional cloning in maize (Zea mays subsp. mays, Poaceae) | |
| Ashrafi et al. | A long‐read transcriptome assembly of cotton (Gossypium hirsutum L.) and intraspecific single nucleotide polymorphism discovery | |
| Vendramin et al. | Genomic tools for durum wheat breeding: de novo assembly of Svevo transcriptome and SNP discovery in elite germplasm | |
| Kerwin et al. | Rampant misexpression in a Mimulus (monkeyflower) introgression line caused by hybrid sterility, not regulatory divergence | |
| Ali et al. | Gene-specific sex-linked genetic markers in date palm (Phoenix dactylifera L.) |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |