CN117524301B - 一种拷贝数变异的检测方法、装置以及计算机可读介质 - Google Patents
一种拷贝数变异的检测方法、装置以及计算机可读介质 Download PDFInfo
- Publication number
- CN117524301B CN117524301B CN202410010138.3A CN202410010138A CN117524301B CN 117524301 B CN117524301 B CN 117524301B CN 202410010138 A CN202410010138 A CN 202410010138A CN 117524301 B CN117524301 B CN 117524301B
- Authority
- CN
- China
- Prior art keywords
- cnv
- sample
- calculating
- tumor
- copy number
- 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
- 238000001514 detection method Methods 0.000 title claims abstract description 83
- 206010028980 Neoplasm Diseases 0.000 claims abstract description 125
- 238000009826 distribution Methods 0.000 claims abstract description 103
- 238000012163 sequencing technique Methods 0.000 claims abstract description 91
- 230000035945 sensitivity Effects 0.000 claims abstract description 54
- 238000012070 whole genome sequencing analysis Methods 0.000 claims abstract description 8
- 210000000349 chromosome Anatomy 0.000 claims description 60
- 238000000034 method Methods 0.000 claims description 41
- 108090000623 proteins and genes Proteins 0.000 claims description 39
- 238000004364 calculation method Methods 0.000 claims description 36
- 230000003321 amplification Effects 0.000 claims description 27
- 238000003199 nucleic acid amplification method Methods 0.000 claims description 27
- 238000012217 deletion Methods 0.000 claims description 26
- 230000037430 deletion Effects 0.000 claims description 24
- 210000000265 leukocyte Anatomy 0.000 claims description 13
- 238000012360 testing method Methods 0.000 claims description 13
- 239000003550 marker Substances 0.000 claims description 12
- 238000012937 correction Methods 0.000 claims description 9
- 230000006870 function Effects 0.000 claims description 9
- 238000004458 analytical method Methods 0.000 claims description 8
- 238000007781 pre-processing Methods 0.000 claims description 6
- 238000012935 Averaging Methods 0.000 claims description 5
- 238000000605 extraction Methods 0.000 claims description 5
- 238000005070 sampling Methods 0.000 claims description 5
- 208000001490 Dengue Diseases 0.000 claims description 3
- 206010012310 Dengue fever Diseases 0.000 claims description 3
- 230000002759 chromosomal effect Effects 0.000 claims description 3
- 238000010276 construction Methods 0.000 claims description 3
- 208000025729 dengue disease Diseases 0.000 claims description 3
- 230000009466 transformation Effects 0.000 claims 1
- 210000001519 tissue Anatomy 0.000 description 18
- 239000008280 blood Substances 0.000 description 7
- 108700015053 epidermal growth factor receptor activity proteins Proteins 0.000 description 7
- 102000052116 epidermal growth factor receptor activity proteins Human genes 0.000 description 7
- YOHYSYJDKVYCJI-UHFFFAOYSA-N n-[3-[[6-[3-(trifluoromethyl)anilino]pyrimidin-4-yl]amino]phenyl]cyclopropanecarboxamide Chemical compound FC(F)(F)C1=CC=CC(NC=2N=CN=C(NC=3C=C(NC(=O)C4CC4)C=CC=3)C=2)=C1 YOHYSYJDKVYCJI-UHFFFAOYSA-N 0.000 description 7
- 210000004369 blood Anatomy 0.000 description 6
- 206010018338 Glioma Diseases 0.000 description 5
- 239000012634 fragment Substances 0.000 description 5
- 208000010507 Adenocarcinoma of Lung Diseases 0.000 description 4
- 208000032612 Glial tumor Diseases 0.000 description 4
- 230000002159 abnormal effect Effects 0.000 description 4
- 201000007983 brain glioma Diseases 0.000 description 4
- 210000001175 cerebrospinal fluid Anatomy 0.000 description 4
- 201000005249 lung adenocarcinoma Diseases 0.000 description 4
- 208000024891 symptom Diseases 0.000 description 4
- 108010009392 Cyclin-Dependent Kinase Inhibitor p16 Proteins 0.000 description 3
- 101001043809 Homo sapiens Interleukin-7 receptor subunit alpha Proteins 0.000 description 3
- 102100021593 Interleukin-7 receptor subunit alpha Human genes 0.000 description 3
- 208000005016 Intestinal Neoplasms Diseases 0.000 description 3
- 206010027476 Metastases Diseases 0.000 description 3
- 210000004556 brain Anatomy 0.000 description 3
- 210000004027 cell Anatomy 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 201000002313 intestinal cancer Diseases 0.000 description 3
- 230000009401 metastasis Effects 0.000 description 3
- 238000007481 next generation sequencing Methods 0.000 description 3
- 206010009944 Colon cancer Diseases 0.000 description 2
- 208000001333 Colorectal Neoplasms Diseases 0.000 description 2
- 102000009508 Cyclin-Dependent Kinase Inhibitor p16 Human genes 0.000 description 2
- 101150039808 Egfr gene Proteins 0.000 description 2
- 101001030211 Homo sapiens Myc proto-oncogene protein Proteins 0.000 description 2
- 102100038895 Myc proto-oncogene protein Human genes 0.000 description 2
- 238000000692 Student's t-test Methods 0.000 description 2
- 238000012512 characterization method Methods 0.000 description 2
- 238000013461 design Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 108700021358 erbB-1 Genes Proteins 0.000 description 2
- 238000011528 liquid biopsy Methods 0.000 description 2
- 238000012544 monitoring process Methods 0.000 description 2
- 239000002773 nucleotide Substances 0.000 description 2
- 125000003729 nucleotide group Chemical group 0.000 description 2
- 238000012353 t test Methods 0.000 description 2
- 210000004881 tumor cell Anatomy 0.000 description 2
- 208000031639 Chromosome Deletion Diseases 0.000 description 1
- 108020004414 DNA Proteins 0.000 description 1
- 230000004544 DNA amplification Effects 0.000 description 1
- 238000000018 DNA microarray Methods 0.000 description 1
- 206010061818 Disease progression Diseases 0.000 description 1
- 108091028043 Nucleic acid sequence Proteins 0.000 description 1
- 108091081062 Repeated sequence (DNA) Proteins 0.000 description 1
- 241000700605 Viruses Species 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000001574 biopsy Methods 0.000 description 1
- 201000011510 cancer Diseases 0.000 description 1
- 238000005119 centrifugation Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 239000003153 chemical reaction reagent Substances 0.000 description 1
- 201000010099 disease Diseases 0.000 description 1
- 230000005750 disease progression Effects 0.000 description 1
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 description 1
- 230000014509 gene expression Effects 0.000 description 1
- 238000003780 insertion Methods 0.000 description 1
- 230000037431 insertion Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000035772 mutation Effects 0.000 description 1
- 108700024542 myc Genes Proteins 0.000 description 1
- 108020004707 nucleic acids Proteins 0.000 description 1
- 150000007523 nucleic acids Chemical class 0.000 description 1
- 102000039446 nucleic acids Human genes 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 230000010076 replication Effects 0.000 description 1
- 238000010998 test method 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/10—Ploidy or copy number 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
- 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
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
- G16B30/10—Sequence alignment; Homology search
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Biotechnology (AREA)
- Medical Informatics (AREA)
- Biophysics (AREA)
- Theoretical Computer Science (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Bioinformatics & Computational Biology (AREA)
- Chemical & Material Sciences (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Analytical Chemistry (AREA)
- Molecular Biology (AREA)
- Genetics & Genomics (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明提供一种拷贝数变异的检测方法、装置以及计算机可读介质，属于生物信息学技术领域，包括对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20‑30×，评估肿瘤纯度最低检出限及其对应检测性能水平(灵敏性和特异性)，并计算相应理论t值；对待测样本进行低深度全基因组测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比t统计量与理论t值，判定CNV状态。避免由基线样本引入额外的误差，并节省了建立基线样本的成本。
Description
技术领域
本发明属于生物信息学技术领域，具体涉及一种拷贝数变异的检测方法、装置以及计算机可读介质。
背景技术
人类基因组中存在大量的变异，根据碱基数目，其可分为单核苷酸变异(SingleNucleotide Variant，SNV)和结构变异(Structural Variation，SV)。拷贝数变异(CopyNumber Variation，CNV)是结构变异的一种形式，是指与参照基因组相比，大小在50bp到数十Mb DNA片段的缺失、插入、复制和复杂多位点变异。近年来的研究表明，基因组片段的CNV通过改变基因剂量或染色体构象影响基因的表达，在疾病进展、表型多态性和进化研究中占据重要的地位。目前在全基因组范围内寻找CNV主要基于两种技术，分别是基因芯片技术(DNA chip)及新一代测序技术(Next Generation Sequencing，NGS)。
现有技术，如专利公告号为CN108427864B的中国发明专利，公开了一种拷贝数变异的检测方法、装置及计算机可读介质，通过采用T分布拟合获得平均测序深度，进而计算拷贝数基线及分析存在异常的基因片段。该分析方法可有效消除各个检测区域因NGS试验方法导致的测序深度的波动性对检测结果的影响和不同panel间检测结果不一致的问题。但该方法需要大量正常样本构成基线，使得检测效率和难度大大增加。
发明内容
为解决上述问题，本发明提供一种拷贝数变异的检测方法、装置以及计算机可读介质。避免由基线样本引入额外的误差，并节省了建立基线样本的成本，且基于低深度全基因组测序，检测标志物不受到引物设计的影响，更容易形成标准化检测。
为了实现上述目的，本发明所采用的技术方案如下：
一方面，本发明提供一种拷贝数变异的检测方法，包括：对肿瘤样本和正常样本取样，其中正常样本需包括有与所述肿瘤样本严格配对的白细胞，进行全基因组高深度测序，测序深度为20-30×，计算出肿瘤纯度最低检出限，并基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的理论t值；对待测样本进行低深度测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比t统计量与理论t值，判定CNV状态。
进一步的，肿瘤样本的取样例数不少于10例，正常样本的取样例数不小于20例。
进一步的，计算肿瘤纯度最低检出限，包括：通过所述肿瘤样本以及与之配对的正常样本计算出肿瘤纯度；基于所述肿瘤纯度和正常样本按照预设纯度针对每一例肿瘤样本进行reads抽取和混样，并计算t分布曲线的标准差；基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，并计算在一定特异性范围下所对应的特异性t值，并基于所述特异性t值计算灵敏性；将同一预设纯度的混合样本的灵敏性取均值，并与所述特异性范围进行配对组合；通过约登指数将不同预设纯度的特异性与灵敏性进行计算，得到满足预设条件下的肿瘤纯度最低检出限。
进一步的，混样后，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb；采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；针对基因水平的拷贝数变异分析采用bin区来构建t分布曲线；针对染色体水平的拷贝数变异采用segment来构建t分布曲线。
进一步的，t分布曲线的标准差的计算方法为：
logR＝log2A
其中A为通过矫正后的reads数；
M＝median(Var(logR1),Var(logR2),…,Var(logRN))
其中Var(logRN)指第N个segment内logR对应的方差值，M为待测样本在全基因组范围内的噪声；
针对染色体水平，特定长度的segment，t分布的标准差为：
其中N代表特定segment的bin个数；
针对基因水平而言，t分布的标准差为：
进一步的，基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，包括：确定CNV阴性分布拷贝数C1和CNV阳性分布拷贝数C2，对于染色体扩增和缺失而言，C1确定为2，C2分别为3和1；对于基因的高倍扩增而言，C1确定为4，C2被设置为一系列梯度{5，6，7，8，9，10}，对于基因的纯合缺失而言，C1设置为1，C2设置为0；
计算CNV阳性分布均值X1和CNV阴性分布X2：
其中T′为预设纯度。
进一步的，所述灵敏性的计算方法为：
通过特异性t值t1计算灵敏性t值t2，
基于所述灵敏性t值t2，通过R语言的pt函数计算出灵敏性。
进一步的，针对染色体水平，确定约登指数等于98％的肿瘤纯度作为该染色体的肿瘤纯度最低检出限；针对基因水平，采用约登指数增长的饱和点确定肿瘤纯度最低检出限，规定肿瘤纯度上升5％，约登指数上升首次小于2％的点作为约登指数增长的饱和点。
进一步的，计算待测样本CNV阳性状态分布和CNV阴性状态分布的t统计量，
当t3>t5并且t4>t6，判定扩增CNV；
当t3<t5并且t4<t6，判定缺失CNV；
其中，t3和t4为CNV阳性状态分布或CNV阴性状态分布的t统计量，t5和t6为使用R语言的qt函数基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值。
本发明还提供了一种拷贝数变异的检测装置，包括：肿瘤纯度最低检出限计算模块，用于对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20-30×，计算出肿瘤纯度最低检出限，并基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值；CNV判定模块，对待测样本进行低深度测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比t统计量与理论t值，判定CNV状态。
进一步的，所述肿瘤纯度最低检出限计算模块包括：样品高深度测序单元，用于对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20-30×，其中正常样本中有与所述肿瘤样本相同数量的样本进行配对；测序数据预处理单元，用于对测序数据比对、排序、去重、去除低质量序列以及建立索引；肿瘤纯度估计单元，基于所述肿瘤样本和与之配对的相同数量的正常样本进行分析，估计肿瘤样本的肿瘤纯度；模拟样本构建单元，针对每例肿瘤样本和正常样本，按照预设纯度进行reads抽取和混样；覆盖度信息统计单元，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb，对bin区内reads数目进行矫正，并对矫正之后的reads数进行log2转化，记为logR，采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；模型方差计算单元，计算全基因组范围内的噪声，基于所述噪声计算染色体水平及基因水平的t分布的标准差；肿瘤纯度最低检出限计算单元，基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，并计算在一定特异性范围下所对应的特异性t值，并基于所述特异性t值计算灵敏性，将同一预设纯度的混合样本的灵敏性取均值，并与所述特异性范围进行配对组合，通过约登指数将不同预设纯度的特异性与灵敏性进行计算，得到满足预设条件下的肿瘤纯度最低检出限。
进一步的，模型方差计算单元中，
M＝median(Var(logR1),Var(logR2),…,Var(logRN))
其中Var(logRN)指第N个segment内logR对应的方差值，M为待测样本在全基因组范围内的噪声，
针对染色体水平，特定长度的segment，t分布的标准差为：
其中N代表特定segment的bin个数；
针对基因水平而言，t分布的标准差为：
进一步的，所述CNV判定模块包括：待测样本的低深度测序单元，用于对待测样品进行低深度全基因组测序，测序深度至少为1×；测序数据预处理单元，用于对测序数据比对、排序、去重、去除低质量序列以及建立索引；覆盖度信息统计单元，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb，对bin区内reads数目进行矫正，并对矫正之后的reads数进行log2转化，记为logR，采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；构建待检测CNV标志物的假设检验单元，基于拷贝数计算CNV阳性状态和CNV阴性状态t分布曲线的logR平均值，计算全基因组范围内的噪声，基于所述噪声计算染色体水平及基因水平的t分布的标准差，构建检测CNV标志物的假设验证模型；t检测判断待检测CNV状态单元，用于计算CNV阳性状态和CNV阴性状态t分布曲线的t统计量，基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值；通过对比所述t统计量与所述t值判定待检测样本CNV状态；输出样本待检测CNV的状态单元，确定样本所有待检测CNV状态，并输出。
本发明还提供了计算机可读介质，记载有可以运行上述拷贝数变异的检测方法的程序。
本发明实施例提供的技术方案带来的有益效果包括：
1.本发明不需要大量正常样本构成的基线，避免由基线样本引入额外的误差，并节省了建立基线样本的成本；2.本发明不受样本类型的限制，可以应用于组织和液态活检的CNV检测；3.由于本发明基于低深度全基因组测序，检测标志物不受到引物设计的影响，更容易形成标准化检测。
附图说明
为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
图1为本发明实施例所提供的肿瘤纯度最低检出限推断的原理图；
图2为本发明实施例所提供的CNV判定原理示意图；
图3为本发明实施例1中对脑胶质瘤1p/19q共缺失的检测结果，向下箭头表示染色体发生拷贝数缺失，向右箭头表示染色体没有发生拷贝数变异；
图4为本发明实施例1中对脑胶质瘤7号染色体扩增和10号染色体缺失的检测结果，向上箭头表示染色体发生拷贝数扩增，向下箭头表示染色体发生拷贝数缺失，向右箭头表示染色体没有发生拷贝数变异；
图5为本发明实施例1中对脑胶质瘤EGFR基因高倍扩增的检测结果，三角形代表EGFR发生高倍扩增，正方形代表EGFR未发生高倍扩增；
图6为本发明实施例1中对脑胶质瘤CDKN2A/B纯合缺失的检测结果，三角形代表基因CDKN2A/B纯合缺失，正方形代表CDKN2A/B未发生纯合缺失；
图7为本发明实施例2中在肠癌组织和血液cfDNA的CNV检测结果，向上箭头表示染色体发生拷贝数扩增，向下箭头表示染色体发生拷贝数缺失，向右箭头表示染色体没有发生拷贝数变异；
图8为本发明实施例3中症状没有好转的动态监测样本的CNV检测结果，向上箭头表示染色体发生拷贝数扩增，向下箭头表示染色体发生拷贝数缺失，向右箭头表示染色体没有发生拷贝数变异，三角形代表基因发生高倍扩增，正方形代表基因未发生高倍扩增，基因从左到右分别为：TERT、IL7R、EGFR、MYC；
图9为本发明实施例3中症状明显好转的动态监测样本的CNV检测结果，向上箭头表示染色体发生拷贝数扩增，向下箭头表示染色体发生拷贝数缺失，向右箭头表示染色体没有发生拷贝数变异，三角形代表基因发生高倍扩增，正方形代表基因未发生高倍扩增。
具体实施方式
下面通过具体实施方式对本发明作进一步详细说明。但本领域技术人员将会理解，下列实施例仅用于说明本发明，而不应视为限定本发明的范围。实施例中未注明具体技术或条件者，按照本领域内的文献所描述的技术或条件或者按照产品说明书进行。所用试剂或仪器未注明生产厂商者，均为可以通过市购获得的常规产品。
本文使用的词语“包括”、“包含”、“具有”或其任何其他变体意欲涵盖非排它性的包括。例如，包括列出要素的工艺、方法、物品或设备不必受限于那些要素，而是可以包括其他没有明确列出或属于这种工艺、方法、物品或设备固有的要素。除非上下文明确规定，否则单数形式“一个/种”和“所述(该)”包括复数个讨论对象。
本发明提供的检测方法仅仅是用于通过测序结果判定是否存在着拷贝数变异现象，并非是用于疾病的诊断和治疗目的。
本发明中“待测样本”是指需要进行检测，并判定该样本上的一个或者多个区域或基因是否存在有拷贝数变异，可以通过组织或液态活检方式获取。“肿瘤样本”是指通过活检，检测出具有肿瘤细胞的样本，“正常样本”是指血液样本中离心得到的白细胞。
如本文所用的，术语“比对”是指将测序序列与参考基因组序列进行比较并且由此确定该测序序列在参考基因组中具体位置的过程。如果该参考序列含有该测序序列，则可以将其映射至参考序列中的某个特定位置。作为结果，比对可得到测序序列对应的基因组位置，并可进而判别其是否属于某一特定区域。术语“参考基因组”或“参考序列”是指生物体或病毒的已知基因组序列(无论是部分的或完整的)，它可以用于对来自受试者的测序序列进行参比。例如，用于人类受试者以及许多其他生物体的参考基因组可见于美国国家生物技术信息中心(ncbi.nlm.nih.gov)，对于人类样品来说，参照序列可以是基因组hg18或hg19版本。目前hg19的相关数据库相对较多且hg19对应的碱基量比hg18要多，所以优先选择hg19。
术语“测序序列(read)”是指来自核酸样品经测序后所获得的短片段序列。典型地，测序序列代表样品中的相邻碱基对的短序列。测序序列通过碱基对序列(ATCG符号)代表。它可以存储在存储设备中，且经过适当处理后可确定其是否与参考基因组序列匹配。测序序列可以直接地从测序装置中或者间接地从涉及样品的存储序列信息中获得。一般而言，测序序列是有一定长度(例如，至少约30bp)的DNA序列，可以通过比对来确定其在染色体对应的具体位置。
每个目标区域上的测序深度信息是比对结果中位于区域内的测序片段数目；位点的测序深度信息是比对结果中包含该位点的测序片段数目。
本发明提供的拷贝数变异的检测方法，包括：
S1.对肿瘤样本和正常样本取样，其中正常样本中有与所述肿瘤样本相同数量的样本进行配对，进行全基因组高深度测序，测序深度为20-30×，计算出肿瘤纯度最低检出限，并基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值。
1.样本高深度测序。为了得到具有足够的样本，肿瘤样本的取样例数不少于10例，正常样本的取样例数不小于20例。本实施例中，选取10例肿瘤组织样本和20例白细胞样本，进行全基因组测序，得到样本的测序数据。在20例白细胞样本中，其中10例白细胞样本是与肿瘤组织配对的，可用于估计对应样本肿瘤组织的纯度。
2.测序数据预处理。对所述10例肿瘤组织样本和20例白细胞样本的测序数据进行预处理，包括了比对、排序、去重、去除低质量序列、建立索引五个步骤。
比对：使用BWA-MAM(v.0.7.12-r1039)将待测样本的fastq文件比对至hg19参考基因组，比对参数为BWA-MAM默认参数，生成待测样本原始bam文件。
排序：使用samtools(v.1.2)软件sort命令对原始bam文件进行排序。
去重：使用Picard(v.1.124)软件MarkDuplicates命令对排序后的bam文件进行重复序列标记和过滤。
去除低质量序列：使用samtools去除MAPQ<20的低质量序列。
建立索引：使用samtools软件的index命令对去除低质量序列后的bam文件构建索引。
3.肿瘤纯度估计。通过所述肿瘤样本以及与之配对的正常样本计算出肿瘤纯度。使用Sclust(v.1.0)软件默认参数对10组肿瘤组织样本和与之配对的白细胞样本的bam文件进行分析，估计肿瘤组织的肿瘤纯度值T。
4.构建模拟样本。基于所述肿瘤纯度和正常样本按照预设纯度针对每一例肿瘤样本进行reads抽取和混样，并计算t分布曲线的标准差。具体的，使用samtools软件的命令"samtools view–s 12.Y"和"samtools view–s 12.Z"分别对每一例肿瘤组织和20例白细胞样本的bam文件按照预设纯度T'∈{0.025,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6}进行reads抽取。因此，得到每一个预设纯度下具有200个混样。
肿瘤细胞抽取的reads比例Y按照式1进行计算：
T'为上述预设纯度梯度对应的任意值；T为肿瘤纯度值；R1为由预处理后测序数据得到的肿瘤样本的reads数，通过命令“samtools view–c肿瘤样本的bam文件”进行计算；21500000代表1×测序条件下的理论的reads数。
白细胞样本抽取reads比例Z按照式2进行计算：
式2中的R2代表预处理后的测序数据中白细胞样本的reads数，计算方法与R1相同；T'和T与式1中保持一致。
使用samtools merge将抽取后的每一例肿瘤样本分别和抽取后的20例白细胞样本逐一混合。由此在某一特定预设纯度T'条件下10例肿瘤样本可以产生200个模拟样本，共计200*13个模拟样本。
以肿瘤纯度为0.6的样本(样本1)和正常白细胞样本(样本2)混合为预设纯度为0.025的样本(样本3)为例：假设样本1的reads总数为430000000(测序深度为20×)，样本2的reads数目为645000000(测序深度为30×)。首先，使用步骤4的抽样方法，抽取样本1的reads比例为：
即为0.21％，抽取样本2的reads数目的比例为：
即为3.19％。使用samtools软件的"samtools view-s 12.0021样本1.bam"和"samtools view-s 12.0319样本2.bam"，分别对肿瘤样本和白细胞样本进行抽取，使用步骤4的混合方法将抽取的样本1和样本2的reads混合，得到预设纯度为0.025，reads数目为21500000的样本3。
5.覆盖度信息统计。混样后，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb；采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；针对基因水平的拷贝数变异分析采用bin区覆盖度表征值来构建t分布曲线；针对染色体水平的拷贝数变异采用segment覆盖度表征值来构建t分布曲线。
具体的，使用HMMcopy(v0.99.0)软件对样本的bam文件进行计数，间隔窗口为50-500kb均可，本实施例优选中间隔窗口选择190kb，每个窗口称为一个bin区。
使用R语言的loess函数的默认参数分别对全基因组bin区内reads数目进行GC含量和比对率两步矫正，并对矫正后的reads数进行log2转化，记为logR。
logR＝log2A (3)
其中，A为矫正后的reads数。
利用DNAcopy(v1.64.0)软件包的CBS算法的默认参数，按照发生CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分成同一区间，记为segment。
后续计算过程中，所有基因层面的拷贝数变异分析均根据bin区来计算，染色体层面的拷贝数变异均根据segment来计算。
6.模型方差计算：t分布曲线的标准差的计算方法为：用待计算样本所有segment的logR方差中位数代表样本全基因组范围的噪声(每个样本独立计算一个噪音)。
M＝median(Var(logR1),Var(logR2),…,Var(logRN)) (4)
其中Var(logRN)指第N个segment内logR对应的方差值，M为待测样本在全基因组范围内的噪声；
t分布的标准差估计：对于染色体水平而言，特定长度的segment，采用抽样分布原理对其标准差进行估计，计算公式如式5：
其中M代表样本全基因组范围的噪声，N代表特定segment的bin个数。例如：对某个样本一号染色体长臂(1p)的标准差计算，N代表特定长度染色体的bin数目。例如：1p染色体长度为125000000bp，以190kb的bin区划分，1p的bin数目为658个。
对于基因水平而言，bin长度设定为1。因此，其分布的标准差即是该样本全基因组范围的噪声的平方根，即为：
7.肿瘤纯度最低检出限的计算，肿瘤纯度最低检出限被定义为在1×测序深度条件下时仍稳定检出CNV的肿瘤纯度最低值。本发明使用每个肿瘤纯度条件下的200个模拟样本进行肿瘤纯度最低检出限的计算。
(1)构建肿瘤纯度最低检出限的计算模型：如图1所示，本发明是基于两个t分布来进行CNV检出的，分别被定义为CNV阳性状态和CNV阴性状态分布，分别对应特定CNV事件发生的有无。其中长度为N bin的segment其t分布自由度设置为N-1，基因水平的自由度设置为2.1。当模型达到一定的特异性和灵敏性时，CNV阳性分布和CNV阴性分布会处于如图1所示的临界状态。临界状态下的特异性和灵敏性通过后续步骤进行计算。
(2)模型阴性分布和阳性分布拷贝数的确定：阴性状态分布拷贝数是指没有发生指定CNV事件的拷贝数，阳性状态分布拷贝数是指发生指定CNV事件的拷贝数。本发明旨在判定染色体的扩增/缺失，基因的高倍扩增/纯合缺失。因此，需预先确定每个标志物CNV阴性分布拷贝数(C1)和阳性分布拷贝数(C2)。对于染色体扩增而言，C1确定为2，C2确定为3，对于染色体缺失而言，C1确定为2，C2确定为1。对于基因的高倍扩增而言，C1确定为4，C2被设置为一系列梯度{5，6，7，8，9，10}，分别用于判断拷贝数大于或等于5，6，7，8，9，10的情况。对于基因的纯合缺失而言，C1设置为1，C2设置为0，拷贝数为0属于基因纯合缺失。
(3)t分布均值计算：如图1所示，阳性状态分布均值(X1)和阴性分布(X2)均值为logR的理论值。计算方法如下式6和式7所示：
公式中C1和C2分别代表如上述所示的阴性分布和阳性分布的拷贝数值；预设纯度T'∈{0.025,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6}。
(4)计算灵敏性：首先，在特异性∈{0.900,0.905,0.910,...,0.99}条件下，使用R语言的“qt”函数计算特异性对应的特异性t值t1。然后，使用式8计算灵敏性t值t2：
公式中的X1和X2使用步骤7的第(3)步进行计算，标准差(SD)使用步骤6进行计算，最后，使用R语言“pt”函数根据灵敏性t值(t2)计算灵敏性，得到19个特异性(0.900,0.905,0.910,...,0.99)所对应的灵敏性。
通过上述方法计算每个CNV标志物的特定预设纯度和特异性条件下200样本的灵敏性，并将200个样本的灵敏性取均值，此时，每个预设纯度条件下，均存在19对灵敏性和特异性的组合。
(5)计算约登指数：通过计算每个预设纯度条件下19对灵敏性和特异性组合的约登指数，确定每个预设纯度条件下最优特异性和灵敏性组合。每个预设纯度条件下的约登指数如式9所示：
YI＝max (Spec+Sens-1) (9)
YI代表约登指数，Spec代表特异性，Sens代表灵敏性；
染色体水平扩增和缺失的约登指数计算结果如表1所示，基因高倍扩增和纯合缺失的约登指数计算结果如表2所示。
表1.染色体水平的约登指数统计
表2基因水平的约登指数统计
/>
(6)肿瘤纯度最低检出限的确定：对于染色体水平而言，确定约登指数等于98％的预设纯度作为该染色体的肿瘤纯度最低检出限，通过计算所有染色体均可稳定检出的预设纯度的最小值作为染色体CNV的肿瘤纯度最低检出限。对于基因水平而言，采用约登指数增长的饱和点确定肿瘤纯度最低检出限，规定预设纯度上升5％，约登指数上升首次小于2％的点作为约登指数增长的饱和点。所有CNV的肿瘤纯度最低检出限如表3所示。
表3.不同类型CNV的肿瘤纯度最低检出限
S2.对待测样本进行低深度测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比t统计量与理论t值，判定CNV状态。
计算待测样本CNV阳性状态分布和CNV阴性状态分布的t统计量，
当t3>t5并且t4>t6，判定扩增CNV；
当t3<t5并且t4<t6，判定缺失CNV；
其中，t3和t4为CNV阳性状态分布和CNV阴性状态分布的t统计量，t5和t6为使用R语言的qt函数基于所述肿瘤纯度最低检出限计算相应灵敏性(表3)和特异性(表3)所对应的t值。
具体的，包括：
1.待测样本的低深度测序。对待检验的样本进行低深度全基因组测序，测序深度至少为1×。
2.测序数据预处理：用于对测序数据比对、排序、去重、去除低质量序列以及建立索引，与S1中的预处理相同。
3.覆盖度信息统计：与S1中的计算方式相同。
4.构建待检测CNV标志物的假设检验模型：如图2所示：在进行扩增或者缺失判断时，分别构建两组t分布，即CNV阳性状态分布和CNV阴性状态t分布。两个分布平均数和标准差计算方法如下：
(1)计算CNV标志物的CNV阴性状态分布和CNV阳性状态分布的logR平均值(X1和X2)，计算方法与S1中的步骤7保持一致，其中参数T'在此步骤为待检测CNV标志物经评估后确定的肿瘤纯度最低检出限值；
(2)待检测样本的标准差计算的计算方法与S1中的步骤6一致，这里采用待检测样本实际segment划分来进行计算；
(3)构建待检测CNV标志物的假设检验模型：通过计算(1)(2)两步的参数，构建如图2所示的假设检验模型。
5.t检验判断待检测CNV状态：根据待检测CNV标志物构建的假设检验模型，通过两次t检验，判定某个segment/bin是否发生了拷贝数变异。
(1).计算t统计量：分别计算CNV阳性状态分布和CNV阴性状态分布两个分布的t统计量，计算公式如式10和式11所示：
其中R代表待检测segment/bin对应的logR值，X1代表CNV阴性分布的均值，X2代表CNV阳性分布的均值，SD代表t分布的标准差。X1、X2、SD均在S2的步骤4中计算。
(2).判断CNV状态：CNV判断通过两步t检验完成，当判定扩增CNV时，要求：
t3>t5并且t4>t6
当判定缺失CNV时，要求：
t3<t5并且t4<t6
t3和t4被S2的步骤5中计算，t5和t6为使用R语言函数“qt”计算表3中特异性和灵敏性所对应的t值。
6.输出样本的CNV状态：通过上述计算，确定样本所有待检测CNV的状态。
本发明还保护了一种拷贝数变异的检测装置，包括：肿瘤纯度最低检出限计算模块，用于对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20-30×，计算出CNV稳定检出所对应的肿瘤纯度最低检出限，并基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的理论t值；CNV判定模块，对待测样本进行低深度测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比t统计量与理论t值，判定CNV状态。
所述肿瘤纯度最低检出限计算模块包括：
样品高深度测序单元，用于对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20-30×，其中正常样本中有与所述肿瘤样本相同数量的样本进行配对；
测序数据预处理单元，用于对测序数据比对、排序、去重、去除低质量序列以及建立索引；
肿瘤纯度估计单元，基于所述肿瘤样本和与之配对的相同数量的正常样本进行分析，估计肿瘤样本的肿瘤纯度；
模拟样本构建单元，针对每例肿瘤样本和正常样本，按照预设纯度进行reads抽取和混样；
覆盖度信息统计单元，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb，对bin区内reads数目进行矫正，并对矫正之后的reads数进行log2转化，记为logR，采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；
模型方差计算单元，计算全基因组范围内的噪声，基于所述噪声计算染色体水平及基因水平的t分布的标准差；
肿瘤纯度最低检出限计算单元，基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，并计算在一定特异性范围下所对应的特异性t值，并基于所述特异性t值计算灵敏性，将同一预设纯度的混合样本的灵敏性取均值，并与所述特异性范围进行配对组合，通过约登指数将不同预设纯度的特异性与灵敏性进行计算，得到满足预设条件下的肿瘤纯度最低检出限。
模型方差计算单元中，
M＝median(Var(logR1),Var(logR2),…,Var(logRN))
其中Var(logRN)指第N个segment内logR对应的方差值，M为待测样本在全基因组范围内的噪声，
针对染色体水平，特定长度的segment，t分布的标准差为：
其中N代表特定segment的bin个数；
针对基因水平而言，t分布的标准差为：
所述CNV判定模块包括：
待测样本的低深度测序单元，用于对待测样品进行低深度全基因组测序，测序深度至少为1×；
测序数据预处理单元，用于对测序数据比对、排序、去重、去除低质量序列以及建立索引；
覆盖度信息统计单元，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb，对bin区内reads数目进行矫正，并对矫正之后的reads数进行log2转化，记为logR，采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；
构建待检测CNV标志物的假设检验单元，基于拷贝数计算CNV阳性状态和CNV阴性状态t分布曲线的logR平均值，计算全基因组范围内的噪声，基于所述噪声计算染色体水平及基因水平的t分布的标准差，构建检测CNV标志物的假设验证模型；
t检测判断待检测CNV状态单元，用于计算CNV阳性状态和CNV阴性状态t分布曲线的t统计量，基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值；通过对比所述t统计量与所述t值判定待检测样本CNV状态；
输出样本待检测CNV的状态单元，确定样本所有待检测CNV状态，并输出。
本发明还保护了一种计算机可读介质，记载有可以运行上述方法的拷贝数变异的检测方法的程序。
显然，本领域的技术人员应该明白，上述的本发明的各模块或各步骤可以用通用的计算装置来实现，它们可以集中在单个的计算装置上，或者分布在多个计算装置所组成的网络上，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。这样，本发明不限制于任何特定的硬件和软件结合。
实施例1
胶质瘤组织样本的CNV判定。选取1p/19q共缺失、7号染色体扩增和10号染色体缺失(+7/-10)、EGFR基因高倍扩增、CDKN2A/B纯合缺失的阴性样本和阳性样本各1例，并对这些样本进行FISH检测，FISH检测严格按照试剂盒说明书执行，最终确定样本的CNV状态如表4所示。
表4.脑胶质瘤FISH检测的样本信息表。
对8例脑胶质瘤患者组织样本进行低深度测序，并参考实施例1中CNV判定模块的步骤，对脑胶质瘤用于分级分型的CNV标志物进行检测，其中EGFR高倍扩增C2参数被设置为10。最终，本发明检测结果与FISH检测的结果比较如图3-6所示：8例样本的CNV状态与FISH检测的结果完全一致。
实施例2
肠癌组织和cfDNA样本的CNV判定。选取一例结直肠癌组织和血液cfDNA配对样本进行低深度测序。使用本发明所提供的拷贝数变异的检测方法对组织和血液cfDNA中已被报道的结直肠癌高频出现的CNV进行判定，包括扩增：7p、7q、8q、13q、20p、20q；缺失：8p、17p、18p、18q。
如图7所示，上述高频出现的CNV在组织(病例9-组织)和血液cfDNA(病例9-血液)中均发生，且组织和cfDNA判定的CNV状态完全一致，说明本申请方法可以应用于血液cfDNA的CNV判定。
实施例3
肺腺癌脑转移cfDNA样本的CNV判定。取两例肺腺癌脑转移病人动态监测的脑脊液样本进行低深度测序，其中病例10为两次随访病情稳定且无进展的样本，其两次随访样本脑脊液异型细胞占比差距较小，分别为：7％和12％。病例11为第二次随访症状好转的样本，两次随访脑脊液异型细胞占比差距较大，分别为：37％和0％。
选取在肺腺癌中高频出现的CNV区域进行判断，包括扩增：1q、5p、7p、8q；缺失：6q、9p、13q、15q；基因高倍扩增：TERT、IL7R、EGFR、MYC。其中TERT、IL7R、EGFR、MYC基因高倍扩增使用参数C2为5进行CNV判断。
使用本发明所提供的拷贝数变异的检测方法对两次随访病情稳定且无进展样本的低深度数据进行CNV判定，结果如图8所示：两次随访中的所有高频CNV表现一致。
本发明所提供的拷贝数变异的检测方法对两次随访症状有好转样本的低深度数据进行CNV判定，结果如图9所示：随访1表现出多个CNV，包括：1q、5p、7p、8p的扩增，6q、15q的缺失和TERT、EGFR的高倍扩增，而对于随访2而言，本申请方法并未鉴定出CNV事件发生。
综合以上四个实施例，本申请方法可以应用于脑胶质瘤、肠癌、肺腺癌脑转移等多个癌症类型包括组织和cfDNA在内的CNV检测。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。
Claims (13)
1.一种拷贝数变异的检测方法，其特征在于，包括：
对肿瘤样本和正常样本取样，其中正常样本需包括有与所述肿瘤样本严格配对的白细胞，进行全基因组高深度测序，测序深度为20-30×，计算出肿瘤纯度最低检出限，并基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的理论t值；
对待测样本进行低深度测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比t统计量与理论t值，判定CNV状态；
计算肿瘤纯度最低检出限包括：
通过所述肿瘤样本以及与之配对的正常样本计算出肿瘤纯度；
基于所述肿瘤纯度和正常样本按照预设纯度针对每一例肿瘤样本进行reads抽取和混样，并计算t分布曲线的标准差；
基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，并计算在一定特异性范围下所对应的特异性t值，并基于所述特异性t值计算灵敏性；
将同一预设纯度的混合样本的灵敏性取均值，并与所述特异性范围进行配对组合；
通过约登指数将不同预设纯度的特异性与灵敏性进行计算，得到满足预设条件下的肿瘤纯度最低检出限。
2.根据权利要求1所述的拷贝数变异的检测方法，其特征在于，肿瘤样本的取样例数不少于10例，正常样本的取样例数不小于20例。
3.根据权利要求2所述的拷贝数变异的检测方法，其特征在于，混样后，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb；
采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；
针对基因水平的拷贝数变异分析采用bin区来构建t分布曲线；
针对染色体水平的拷贝数变异采用segment来构建t分布曲线。
4.根据权利要求3所述的拷贝数变异的检测方法，其特征在于，t分布曲线的标准差的计算方法为：
logR＝log2A
其中A为通过矫正后的reads数；
M＝median(Var(logR1)，Var(logR2)，...，Var(logRN))
其中Var(logRN)指第N个segment内logR对应的方差值，M为待测样本在全基因组范围内的噪声；
针对染色体水平，特定长度的segment，t分布的标准差为：
其中N代表特定segment的bin个数；
针对基因水平而言，t分布的标准差为：
5.根据权利要求4所述的拷贝数变异的检测方法，其特征在于，基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，包括：
确定CNV阴性分布拷贝数C1和CNV阳性分布拷贝数C2，对于染色体扩增和缺失而言，C1确定为2，C2分别为3和1；对于基因的高倍扩增而言，C1确定为4，C2被设置为一系列梯度{5，6，7，8，9，10}，对于基因的纯合缺失而言，C1设置为1，C2设置为0；
计算CNV阳性分布均值X1和CNV阴性分布X2：
其中T′为预设纯度。
6.根据权利要求5所述的拷贝数变异的检测方法，其特征在于，所述灵敏性的计算方法为：
通过特异性t值t1计算灵敏性t值t2，
基于所述灵敏性t值t2，通过R语言的pt函数计算出灵敏性。
7.根据权利要求2所述的拷贝数变异的检测方法，其特征在于，
针对染色体水平，确定约登指数等于98％的肿瘤纯度作为该染色体的肿瘤纯度最低检出限；
针对基因水平，采用约登指数增长的饱和点确定肿瘤纯度最低检出限，规定肿瘤纯度上升5％，约登指数上升首次小于2％的点作为约登指数增长的饱和点。
8.根据权利要求1所述的拷贝数变异的检测方法，其特征在于，
计算待测样本CNV阳性状态分布和CNV阴性状态分布的t统计量，
当t3>t5并且t4>t6，判定扩增CNV；
当t3<t5并且t4<t6，判定缺失CNV；
其中，t3和t4为CNV阳性状态分布或CNV阴性状态分布的t统计量，t5和t6为使用R语言的qt函数基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值。
9.一种拷贝数变异的检测装置，其特征在于，包括：
肿瘤纯度最低检出限计算模块，用于对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20-30×，计算出肿瘤纯度最低检出限，并基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值；
CNV判定模块，对待测样本进行低深度测序，测序深度至少为1×，构建CNV阳性状态和CNV阴性状态t分布曲线，针对待确定CNV状态的logR值，分别计算其在上述t分布对应的t统计量，通过对比计算t统计量与理论t值，判定CNV状态；
所述肿瘤纯度最低检出限计算模块包括：
样品高深度测序单元，用于对肿瘤样本和正常样本进行全基因组高深度测序，测序深度为20-30×，其中正常样本中有与所述肿瘤样本相同数量的样本进行配对；
肿瘤纯度估计单元，基于所述肿瘤样本和与之配对的相同数量的正常样本进行分析，估计肿瘤样本的肿瘤纯度；
模拟样本构建单元，针对每例肿瘤样本和正常样本，按照预设纯度进行reads抽取和混样；
模型方差计算单元，计算全基因组范围内的噪声，基于所述噪声计算染色体水平及基因水平的t分布的标准差；
肿瘤纯度最低检出限计算单元，基于拷贝数构建CNV阳性状态和CNV阴性状态t分布曲线，并计算在一定特异性范围下所对应的特异性t值，并基于所述特异性t值计算灵敏性，将同一预设纯度的混合样本的灵敏性取均值，并与所述特异性范围进行配对组合，通过约登指数将不同预设纯度的特异性与灵敏性进行计算，得到满足预设条件下的肿瘤纯度最低检出限。
10.根据权利要求9所述的拷贝数变异的检测装置，其特征在于，所述肿瘤纯度最低检出限计算模块还包括：
测序数据预处理单元，用于对测序数据比对、排序、去重、去除低质量序列以及建立索引；
覆盖度信息统计单元，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb，对bin区内reads数目进行矫正，并对矫正之后的reads数进行log2转化，记为logR，采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment。
11.根据权利要求10所述的拷贝数变异的检测装置，其特征在于，模型方差计算单元中，
M＝medianVar(logR1)，Var(logR2)，...，Var(logRN))
其中Var(logRN)指第N个segment内logR对应的方差值，M为待测样本在全基因组范围内的噪声，
针对染色体水平，特定长度的segment，t分布的标准差为：
其中N代表特定segment的bin个数；
针对基因水平而言，t分布的标准差为：
12.根据权利要求9所述的拷贝数变异的检测装置，其特征在于，所述CNV判定模块包括：
待测样本的低深度测序单元，用于对待测样品进行低深度全基因组测序，测序深度至少为1×；
测序数据预处理单元，用于对测序数据比对、排序、去重、去除低质量序列以及建立索引；
覆盖度信息统计单元，计算固定间隔窗口内的reads数，所述固定间隔窗口称为一个bin区，间隔窗口的长度为50-500kb，对bin区内reads数目进行矫正，并对矫正之后的reads数进行log2转化，记为logR，采用CBS算法对CNV变异的断点位置对染色体进行分割，将具有相同拷贝数的连续bin区划分为同一个segment；
构建待检测CNV标志物的假设检验单元，基于拷贝数计算CNV阳性状态和CNV阴性状态t分布曲线的logR平均值，计算全基因组范围内的噪声，基于所述噪声计算染色体水平及基因水平的t分布的标准差，构建检测CNV标志物的假设验证模型；
t检测判断待检测CNV状态单元，用于计算CNV阳性状态和CNV阴性状态t分布曲线的t统计量，基于所述肿瘤纯度最低检出限计算相应灵敏性和特异性所对应的t值；通过对比所述t统计量与所述t值判定待检测样本CNV状态；
输出样本待检测CNV的状态单元，确定样本所有待检测CNV状态，并输出。
13.一种计算机可读介质，其特征在于，记载有可以运行权利要求1～8任一项所述的拷贝数变异的检测方法的程序。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202410010138.3A CN117524301B (zh) | 2024-01-04 | 2024-01-04 | 一种拷贝数变异的检测方法、装置以及计算机可读介质 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202410010138.3A CN117524301B (zh) | 2024-01-04 | 2024-01-04 | 一种拷贝数变异的检测方法、装置以及计算机可读介质 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117524301A CN117524301A (zh) | 2024-02-06 |
| CN117524301B true CN117524301B (zh) | 2024-04-09 |
Family
ID=89742324
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202410010138.3A Active CN117524301B (zh) | 2024-01-04 | 2024-01-04 | 一种拷贝数变异的检测方法、装置以及计算机可读介质 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117524301B (zh) |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106055923A (zh) * | 2016-05-13 | 2016-10-26 | 万康源(天津)基因科技有限公司 | 一种基因拷贝数变异分析方法 |
| CN110029157A (zh) * | 2018-01-11 | 2019-07-19 | 北京大学 | 一种检测肿瘤单细胞基因组单倍体拷贝数变异的方法 |
| CN112669901A (zh) * | 2020-12-31 | 2021-04-16 | 北京优迅医学检验实验室有限公司 | 基于低深度高通量基因组测序的染色体拷贝数变异检测装置 |
| CN112802548A (zh) * | 2021-01-07 | 2021-05-14 | 深圳吉因加医学检验实验室 | 单样本全基因组预测等位基因特异性拷贝数变异的方法 |
| CN114724631A (zh) * | 2022-04-19 | 2022-07-08 | 上海仁东医学检验所有限公司 | 染色体拷贝数变异程度评估模型、方法及应用 |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| GB201318369D0 (en) * | 2013-10-17 | 2013-12-04 | Univ Leuven Kath | Methods using BAF |
| EP4227947A1 (en) * | 2013-10-21 | 2023-08-16 | Verinata Health, Inc. | Method for improving the sensitivity of detection in determining copy number variations |
-
2024
- 2024-01-04 CN CN202410010138.3A patent/CN117524301B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106055923A (zh) * | 2016-05-13 | 2016-10-26 | 万康源(天津)基因科技有限公司 | 一种基因拷贝数变异分析方法 |
| CN110029157A (zh) * | 2018-01-11 | 2019-07-19 | 北京大学 | 一种检测肿瘤单细胞基因组单倍体拷贝数变异的方法 |
| CN112669901A (zh) * | 2020-12-31 | 2021-04-16 | 北京优迅医学检验实验室有限公司 | 基于低深度高通量基因组测序的染色体拷贝数变异检测装置 |
| CN112802548A (zh) * | 2021-01-07 | 2021-05-14 | 深圳吉因加医学检验实验室 | 单样本全基因组预测等位基因特异性拷贝数变异的方法 |
| CN114724631A (zh) * | 2022-04-19 | 2022-07-08 | 上海仁东医学检验所有限公司 | 染色体拷贝数变异程度评估模型、方法及应用 |
Non-Patent Citations (2)
| Title |
|---|
| Sensitive detection of CSF-derived cell-free DNA at the CNV level by shallow whole-genome sequencing for diagnosis and characterization of CNS metastases;Li Yubin.et.;《Journal of Clinical Oncology》;20230531;第41卷(第16期);e21185 * |
| 基于固相增菌法和纳米生物传感技术的沙门氏菌可视化高效检测技术的研究;唐锋;《中国博士学位论文全文数据库医药卫生科技辑》;20170815(第8期);E060-26 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117524301A (zh) | 2024-02-06 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US11371074B2 (en) | Method and system for determining copy number variation | |
| CN107423578B (zh) | 检测体细胞突变的装置 | |
| CN110910957B (zh) | 一种基于单肿瘤样本高通量测序微卫星不稳定性探测位点筛选方法 | |
| CN109767810B (zh) | 高通量测序数据分析方法及装置 | |
| CN104302781B (zh) | 一种检测染色体结构异常的方法及装置 | |
| CN104846089B (zh) | 一种孕妇外周血中胎儿游离dna比例的定量方法 | |
| CN110993029B (zh) | 一种检测染色体异常的方法及系统 | |
| CN107480470B (zh) | 基于贝叶斯与泊松分布检验的已知变异检出方法和装置 | |
| CN115064211B (zh) | 一种基于全基因组甲基化测序的ctDNA预测方法及装置 | |
| CN111081315A (zh) | 一种同源假基因变异检测的方法 | |
| CN110060733B (zh) | 基于单样本的二代测序肿瘤体细胞变异检测装置 | |
| CN110016497B (zh) | 一种检测肿瘤单细胞基因组拷贝数变异的方法 | |
| CN108804876B (zh) | 用于计算癌症样本纯度和染色体倍性的方法和装置 | |
| CN110093417B (zh) | 一种检测肿瘤单细胞体细胞突变的方法 | |
| KR102405245B1 (ko) | 전장유전체 시퀀싱 기반의 염색체 이상 검출 방법 및 그 용도 | |
| CN108304694B (zh) | 基于二代测序数据分析基因突变的方法 | |
| CN111321209A (zh) | 一种用于循环肿瘤dna测序数据双端矫正的方法 | |
| WO2019213811A1 (zh) | 检测染色体非整倍性的方法、装置及系统 | |
| CN109461473B (zh) | 胎儿游离dna浓度获取方法和装置 | |
| CN105483210A (zh) | 一种rna编辑位点的检测方法 | |
| Simonyan et al. | HIVE-heptagon: a sensible variant-calling algorithm with post-alignment quality controls | |
| CN117524301B (zh) | 一种拷贝数变异的检测方法、装置以及计算机可读介质 | |
| CN110373458B (zh) | 一种地中海贫血检测的试剂盒及分析系统 | |
| CN109033752B (zh) | 一种基于长读长测序的多基因融合检测方法 | |
| WO2019213810A1 (zh) | 检测染色体非整倍性的方法、装置及系统 |
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 |