CN113380324A - 一种T细胞受体序列motif组合识别检测方法、存储介质及设备 - Google Patents
一种T细胞受体序列motif组合识别检测方法、存储介质及设备 Download PDFInfo
- Publication number
- CN113380324A CN113380324A CN202110536816.6A CN202110536816A CN113380324A CN 113380324 A CN113380324 A CN 113380324A CN 202110536816 A CN202110536816 A CN 202110536816A CN 113380324 A CN113380324 A CN 113380324A
- Authority
- CN
- China
- Prior art keywords
- motif
- population
- tumor
- matrix
- data
- 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.)
- Granted
Links
Images
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/50—Mutagenesis
-
- 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
-
- 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
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
-
- 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)
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Biophysics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Biotechnology (AREA)
- Evolutionary Biology (AREA)
- General Health & Medical Sciences (AREA)
- Theoretical Computer Science (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Analytical Chemistry (AREA)
- Chemical & Material Sciences (AREA)
- Bioethics (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Data Mining & Analysis (AREA)
- Databases & Information Systems (AREA)
- Epidemiology (AREA)
- Evolutionary Computation (AREA)
- Public Health (AREA)
- Software Systems (AREA)
- Genetics & Genomics (AREA)
- Molecular Biology (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Investigating Or Analysing Biological Materials (AREA)
Abstract
本发明公开了一种T细胞受体序列motif组合识别检测方法、存储介质及设备,构建Tumor‑Health矩阵和Tissue‑Blood矩阵作为输入矩阵;以motif组合的类内类间距离可分性作为优化目标函数;设计双种群遗传算法,对输入矩阵中的每一motif列和label列进行相关性分析,采用点二列相关性作为衡量motif列和label列是否相关的指标;以优化目标函数作为算法优化目标,利相关性指标对输入矩阵的motif进行初步筛选,过滤掉不相关的motif列,通过遗传算法进行迭代计算,得到两组种群的优势解,取两组种群中适应度排前三的染色体作为最优解,并解码成为对应的motif集合,取对应集合的交集作为最终挖掘出的motif,完成识别检测。本发明能够为免疫治疗提供可能的生物标志物研究方向。
Description
技术领域
本发明属于数据分析技术领域,具体涉及一种T细胞受体序列motif组合识别检测方法、存储介质及设备。
背景技术
免疫组库是指某个个体在任何特定时间点其循环系统中所有功能多样性B淋巴细胞和T淋巴细胞的总和。T细胞受体(T cell receptor,TCR)是细胞识别抗原的媒介,可反映肿瘤的发病原因,肿瘤进展和免疫应答情况,其中CDR3(Complementarity determiningregion,CDR)可直接与抗原复合物相互作用且高度可变,经常用于确定T细胞克隆类型。
近年来免疫组库高通量测序技术的发展和不断成熟,为研究人员研究T细胞受体序列多样性开辟了新的途径,使得T细胞受体序列能够用于定量分析不同样本的适应性免疫反应,目前已有的研究和方法大致可分为两类,第一类通过描述性方法寻找共享的T细胞受体序列或过表达的序列簇,CMV-T细胞受体算法挖掘CMV-positive和CVM-negative(巨细胞病毒)的T细胞公共克隆型,使用Fisher检验来识别CMV相关的克隆并在概率分类模型上使用找到的公共克隆型来预测个体的CMV状态;RECOLD算法通过对成对序列的比对,将高维序列的计算结果映射到低维空间,从而比较样本间的免疫系统相似性;第二类方法通过对T细胞受体序列进行motif解构来进行样本分类或挖掘关键motif,LR-MIL选择将T细胞受体序列解构为长度为4的motif片段,并通过结合Atchley因子采用多示例学习和逻辑回归模型识别出区分motif。
然而,现有方法存在以下问题:
1)由于T细胞受体序列多样性的特征,这些方法都存在过拟合、泛化性能较差的缺陷,一些特殊的降维方法也使得特征的可解释性降低,不利于后续个体免疫系统图谱和个性化疫苗的研究。
2)目前的研究使用的都是单一的外周血测序数据或组织测序数据,并未考虑到样本在两组数据中也会存在差异,由于肿瘤的高度异质性和持续进化性,组织测序数据会存在取样偏差,相较之下外周血测序具有无创、快速、全面的特点,能够在一定程度上克服组织测序的劣势,但仍然缺乏足够的试验数据表明其有效性和实用性,且存在组织测序和外周血测序结果不一致的现象。
发明内容
本发明所要解决的技术问题在于针对上述现有技术中的不足,提供一种T细胞受体序列motif组合识别检测方法、存储介质及设备,在进行T细胞受体序列数据分析时采用motif解构的方法缓解维度问题,同时关注两种数据的差异,即不仅以健康样本和肿瘤样本的T细胞受体数据比对分析为结果,同时还要考虑肿瘤样本本身在外周血和组织上的差异;根据已有的研究成果提出合理的推断:存在T细胞受体序列的motif片段在健康样本和肿瘤样本之间存在差异,同时排除个性化motif的可能,这些片段在肿瘤样本的外周血和组织测序中也存在差异,符合这种规律的片段更具有特异性。
本发明采用以下技术方案:
一种T细胞受体序列motif组合识别检测方法,包括以下步骤:
S1、将健康样本和肿瘤样本的基因组DNA比对到胚系序列,鉴定对应的CDR3区段并将CDR3区段翻译成氨基酸序列;进行氨基酸序列剪切并将健康样本和肿瘤样本的CDR3区段对应的氨基酸序列解构为多个相同长度的连续氨基酸片段;利用健康样本的外周血数据、肿瘤样本的组织数据和肿瘤样本外周血数据解构后的信息构建Tumor-Health矩阵,利用肿瘤样本组织数据和肿瘤样本外周血数据构建Tissue-Blood矩阵,Tumor-Health矩阵和Tissue-Blood矩阵的最后一列为label列,代表样本对应的标签,剩余列为motif列,代表样本对应氨基酸片段的克隆数,将Tumor-Health矩阵和Tissue-Blood矩阵共同作为输入矩阵;
S2、根据步骤S1得出的两个输入矩阵分别构建求解对应的优化目标函数,以motif组合的类内类间距离可分性作为优化目标计算方法,将两个优化目标函数加和得到总的优化目标函数;
S3、根据步骤S1得到的输入矩阵和步骤S2得到的总优化目标函数设计双种群遗传算法并进行求解,确定编码方法和初始种群,采用锦标赛选择法加精英保留制度选择算子对种群进行选择,保留最优个体;采取种群内和种群间的混合交叉方式确定交叉算子;采用基本位变异的操作确定变异算子;
S4、对步骤S1构建的输入矩阵中的每一motif列和label列进行相关性分析,采用点二列相关性作为衡量motif列和label列是否相关的指标,根据计算得到的相关性指标ρ判断显著水平,确定相关性指标ρ>0.05为相关性不显著;
S5、以步骤S2构建的优化目标函数作为算法优化目标,利用步骤S4计算的相关性指标对步骤S1构建的输入矩阵进行初步筛选,过滤掉不相关的motif列,通过步骤S3设计的双种群遗传算法进行迭代计算,得到两组种群的优势解,取两组种群中适应度排前三的染色体作为最优解,并解码成为对应的motif集合,最终得到六组motif集合,两两对应,取对应集合的交集作为最终挖掘出的motif,完成识别检测。
具体的,步骤S1中,构建Tumor-Health矩阵和Tissue-Blood矩阵具体为:
肿瘤样本数为M,健康样本数为N,肿瘤样本外周血和组织测序的和为2M,第一个用于输入的矩阵Tumor-Health共有N+2M行,代表肿瘤样本外周血、组织的测序结果和健康样本外周血测序结果,第二个输入的矩阵Tissue-Blood共有2M行,代表肿瘤样本外周血和组织的测序结果,两个矩阵的每一行均包含一个向量C={C1,C2,...Cs,label},其中,s=8000,Ci代表样本中motif出现的计数,如果Ci=0则代表样本无对应motif,label代表样本对应的标签,在Tumor-Health矩阵中0代表正常的样本数据,1代表非正常的样本数据,在Tissue-blood矩阵中0代表非正常外周血数据,1代表非正常组织数据。
具体的,步骤S2中,总优化目标函数min F(x)为:
min F(x)=F1(x)+F2(x)
其中,F1(x)为优化目标1,在Tumor-Health矩阵上进行类内类间距离可分计算,F2(x)为优化目标2,在Tissue-blood矩阵上进行类内类间距离可分计算。
具体的,步骤S3中,确定编码方法和初始种群具体为:
采用二进制编码,使用二值符号集{0,1}构成种群中每个染色体的基因型,每一个二进制位对应一种motif,1表示染色体对应的可行解中包括对应motif,0表示不包含对应motif;对于两组数据采用随机生成的方式分别生成两个种群,然后通过交叉让两组数据中的优势基因进行互换,完成初始化。
具体的,步骤S3中,采用锦标赛选择法加精英保留制度选择算子对种群进行筛选具体为:
S3021、计算种群中所有染色体的适应度值,将最好的染色体直接保留至下一代;
S3022、确定每次选择的个体数量k;
S3023、从种群中随机选择k个个体构成组,根据其适应度值选择最好的个体进入子代种群;
S3024、重复操作直到新一代种群中染色体个数达到要求,选择算子筛选产生一组新的种群解用于混合交叉过程计算。
具体的,步骤S3中,采用种群内和种群间的混合交叉方式确定交叉算子具体为:
生成一个范围在(0,1)之内的随机数,计算染色体适应度值和自适应交叉率,将随机数与适应度最大的染色体对应的自适应交叉率进行比较,如果自适应交叉率大则进行交叉操作生成新个体,自适应交叉率具体计算如下:
其中,fmax代表当前种群中的最大值,f'代表交叉的两个染色体中较大的适应度值,favg代表种群的平均适应度值,k1代表计算常数。
具体的,步骤S3中,确定变异算子具体为:
S3041、比较随机数与自适应变异率判断对应染色体是否需要进行变异,如果变异率大则进行变异,自适应变异率计算如下:
其中,fmax代表当前种群中的最大值,f'代表要变异的染色体适应度值,favg代表种群的平均适应度值,k2代表计算常数,取k2=0.1。
S3042、根据染色体包含motif数量选择变异类型
S3043、对每一条染色体,统计所有变异位点对应的相关系数,根据系数比例计算每个位点对应的变异概率p(xi),同时计算累计概率q(xi),即每个个体之前所有个体的选择概率之和:
S3044、随机生成s,范围[0,1],若q(xi-1)<s<q(xi),选择位点xi作为变异算子进行变异算子操作,经过变异算子操作更新种群中的染色体得到新个体。
具体的,步骤S4中,每一motif列和label列相关性系数ρ具体为:
本发明的另一个技术方案是,一种存储一个或多个程序的计算机可读存储介质,所述一个或多个程序包括指令,所述指令当由计算设备执行时,使得所述计算设备执行所述的方法中的任一方法。
本发明的另一个技术方案是,一种计算设备,包括:
一个或多个处理器、存储器及一个或多个程序,其中一个或多个程序存储在所述存储器中并被配置为所述一个或多个处理器执行,所述一个或多个程序包括用于执行所述的方法中的任一方法的指令。
与现有技术相比,本发明至少具有以下有益效果:
本发明一种T细胞受体序列motif组合识别检测方法,T细胞受体序列motif组合识别检测方法,在进行T细胞受体序列数据分析时采用motif解构的方法缓解维度问题,关注两种数据的差异,同时考虑肿瘤数据本身在外周血和组织上的差异,采用序列剪切和双种群遗传算法的优势在于能够缓解序列异质性较高的问题,同时在查找满足条件的motif片段时能够进行交叉加快收敛速度,提高效率,最终通过优化方法有效挖掘出肿瘤特异性的motif组合。
进一步的,将样本数据处理为算法需要的形式,即Tumor-Health矩阵和Tissue-Blood,便于后续算法在数据上迭代计算的效率,同时为每一条数据赋予标签,是为了后续计算优化目标函数,每个矩阵都包含对应的motif列和label列,label列即为每个样本的标签,motif列即为每个样本在对应氨基酸片段下的克隆数,两个矩阵是后续计算的输入。
进一步的,根据数据的输入和求解目标,以motif组合的类内类间距离可分性作为优化目标函数,根据输入数据的形式,设置了两个优化目标函数,并将其整合为算法整体求解的总目标函数,有助于迭代计算结果的呈现。
进一步的,确定遗传算法中染色体的编码形式为二进制编码,这种编码形式符和我们最终求解的目标,即motif集合,{0,1}能够很好地代表一个motif的状态,即1代表包含在结果集内,而0代表不包含在结果集内,采用随机生成的方式通过编码得到一组初始种群,用于后续的计算和迭代。
进一步的,采用锦标赛选择法加精英保留制度作为选择算子。在迭代计算过程中,选择算子的作用就是从父代中选出部分个体作为遗传信息传递给下一代,通常以适应度作为淘汰指标。如果一个选择算法选择多样性降低,便会导致种群过早的收敛到局部最优点,导致“早熟”,而选择策略过于发散则会导致算法难以收敛到最优点。锦标赛选择法增加了种群多样性,同时精英保留制度能够保存已有的寻优结果,通过选择算子淘汰部分个体,保留优质个体作为新的种群进行后续计算。
进一步的,采取种群内和种群间的混合交叉方式,种群内确定交叉的起始和终止位置,对范围内的基因组进行交叉,种群间采用单点交叉的方式,在不破坏两组数据中找到的优势解的前提下,对两组种群中对应的优势motif进行互换,这种交叉策略的好处在于既增加了种群内部的染色体组的多样性,同时引入其他种群的优势基因从而加速算法在两个种群上找出符合共同条件的解。自适应交叉率能够根据已有结果集更好的进行交叉操作,混合交叉算子操作后改变了部分个体的结构,生成一组新的种群用于探索更多集合的效果。
进一步的,采用基本位变异的操作,变异也采用自适应变异率来确定是否进行变异操作,同时选择突变基因位时采用轮盘置操作,使得突变频率更高的位点更容易发生变异操作,变异算子操作的目的同样是为了改变部分个体的结构,增加解集合的丰富度,生成更优质的种群。
进一步的,采用相关性分析作为迭代计算初始优化方法,因为存在一些motif,他们只出现在少部分样本中,具有极高的样本特异性,这些motif片段的计算不会提高适应度函数值,同时在后续寻优算法的实现中还浪费了计算资源。同时该指标作为突变频率影响变异操作,使得结果集更加符合寻优要求,具体针对motif列数据和label列数据特征采用了点二列相关性计算,将与label列不显著相关的motif删除,确定了初始解对应的总的motif集合。
综上所述,本发明能够在维度问题和数据差异的问题下求解得到肿瘤相关性motif组合,为免疫治疗提供可能的生物标志物研究方向。
下面通过附图和实施例,对本发明的技术方案做进一步的详细描述。
附图说明
图1为MDPGA算法数据流程处理图;
图2为MDPGA算法输入矩阵示例图;
图3为MDPGA算法输出结果流程图;
图4为不同K值下的一类错误率和二类错误了比对结果图,其中,(a)为K=30,(b)为K=40,(c)为K=50,(d)为K=60,(e)为K=70;
图5为与BOW_SVM和1-DBC算法对比结果图;
图6为不同样本中motif的平均计数统计图,其中,(a)为算法寻优结果的motif片段计数在tumor和health样本中的差异,(b)为算法寻优结果的motif片段计数在tissue和blood样本中的差异。
具体实施方式
本发明一种T细胞受体序列motif组合识别检测方法,包括以下步骤:
S1、分别取健康样本和肿瘤样本,使用QIAGEN Multiplex PCR Kit引物将健康样本和肿瘤样本基因组DNA分别进行多重PCR扩增,包括32对正向V基因引物和13对反向J基因引物;扩增得到的产物继续用Illumina通用引物进行第二轮扩增,产生插入片段大小为100bp的文库;在Illumina NovaSeq平台上进行双端测序,对获得的测序数据去除测序接头序列,过滤低质量测序序列,将过滤后的数据使用fastp软件进行Q20、Q30、GC含量、N含量、平均读长长度和Clean_base_ratio的筛选,筛选基于下列设定阈值:Q20>90%、Q30>85%、GC含量>40%并且<60%、N含量<10.00%、平均读长长度>90bp并且≤110bp和Clean_base_ratio>80%;用Pear软件将经质控的高质量pair-end读段合并为contigs;用MiXCR软件将合并后的contigs比对到胚系序列,鉴定出TCR克隆的种类以及CDR3区段,并将CDR3区段翻译成氨基酸序列;选择k=3进行序列剪切,将健康样本和肿瘤样本的CDR3区段对应的氨基酸序列解构为多个相同长度的连续氨基酸片段,对健康样本和肿瘤样本CDR3区段的序列进行motif解构,同时考虑每个样本对应的每种序列出现的计数;利用康样本外周血数据和肿瘤样本所有数据解构后的信息构建Tumor-Health矩阵,利用肿瘤样本组织和外周血数据构建Tissue-Blood矩阵,两个矩阵的最后一列为label列,代表样本对应的标签,除最后一列的其他列为motif列,代表样本对应氨基酸片段的克隆数,两个矩阵共同作为输入矩阵,用于后续步骤的求解;
假设肿瘤样本数为M,健康样本数为N,2M则代表肿瘤样本外周血和组织测序的和,第一个用于算法输入的矩阵Tumor-Health共有N+2M行,代表肿瘤样本外周血、组织的测序结果和健康样本外周血测序结果,第二个矩阵Tissue-Blood共有2M行,代表肿瘤样本外周血和组织的测序结果,两个矩阵的每一行均包含一个向量C={C1,C2,...Cs,label},其中,s=8000,因为20种氨基酸组成的长度为3的motif共有8000种,Ci代表样本中该motif出现的计数,如果Ci=0则代表该样本无对应motif,label代表样本对应的标签,在Tumor-Health矩阵中0代表正常的样本数据,1代表非正常的样本数据,在Tissue-blood矩阵中0代表非正常外周血数据,1代表非正常组织数据。数据处理流程和算法输入形式如图1和图2所示。
S2、根据步骤S1得出的两个输入矩阵分别构建求解对应的优化目标函数,以motif组合的类内类间距离可分性作为优化目标计算方法,将两个优化目标加和得到总的优化目标函数,类内类间距离可分性是用来评估特征用于分类好坏的方法,相较于基于熵和基于概率的可分性判别,这种方法在计算上简便,直观概念清楚,直接从样本距离之间进行计算适合应用于研究的数据。在本发明中,motif在不同类别的样本中具有类内类间可分的特点,具体表现为数据在类内离散度较小,而在类间离散度较大,这样就代表在不同类别间的差异越大,令Sw为类内离散度矩阵,Sb为类间离散度矩阵,其计算方式如下:
其中,i代表类别,j代表第i类的样本,l代表类别总数,x代表样本对应的motif计数,μ代表i类样本在对应motif上的均值,其公式如下:
定义优化的目标函数:
优化目标1:
优化目标2:
总优化目标:
min F(x)=F1(x)+F2(x)
其中,tr代表矩阵迹的计算,优化目标1主要计算第一个矩阵的相关数据,优化目标2主要计算第二个矩阵的相关数据,将问题转化为一个双目标优化问题。
S3、设计双种群遗传算法
S301、编码方法和初始种群
采用二进制编码,使用二值符号集{0,1}构成种群中每个染色体的基因型,其中每一个二进制位对应一种motif,该位为1则表示该染色体对应的可行解中包括对应motif,0则表示不包含对应motif。对于两组数据采用随机生成的方式分别生成两个种群,目的是为了能尽快在两组数据中找到最优的染色体,然后通过交叉让两组数据中的优势基因进行互换,有利于更快的寻找到同时符合两组数据条件的解。
S302、选择算子
在迭代计算过程中,为了符合“适应度大代表染色体更优”的原则,取负优化目标1和负优化目标2作为两组种群的适应度函数;选择算子的作用就是从父代中选出部分个体作为遗传信息传递给下一代,通常以适应度作为淘汰指标。如果一个选择算法选择多样性降低,便会导致种群过早的收敛到局部最优点,导致“早熟”,而选择策略过于发散则会导致算法难以收敛到最优点。
本发明采用锦标赛选择法加精英保留制度选择算子对种群进行筛选,具体包含四个步骤:
S3021、计算种群中所有染色体的适应度值,将最好的染色体直接保留至下一代;
S3022、确定每次选择的个体数量k;
S3023、从种群中随机选择k个个体构成组,根据其适应度值选择最好的个体进入子代种群;
S3024、重复操作直到新一代种群中染色体个数达到要求,保留优质个体。
S303、交叉算子;
采取种群内+种群间的混合交叉方式对染色体进行交叉,种群内采用Partial-Mapped Crossover的算子模式,即确定交叉的起始和终止位置,对范围内的基因组进行交叉,种群间采用单点交叉的方式,在不破坏两组数据中找到的优势解的前提下,对两组种群中对应的优势motif进行互换,这种交叉策略的好处在于既增加了种群内部的染色体组的多样性,同时引入其他种群的优势基因从而加速算法在两个种群上找出符合共同条件的解。采用自适应交叉率来确定是否进行交叉操作,生成新个体,扩展种群丰富度,具体计算如下:
其中,fmax代表当前种群中的最大值,f'代表交叉的两个染色体中较大的适应度值,favg代表种群的平均适应度值,k1代表计算常数,取k1=0.8,由公式可知,自适应交叉率的优势在于可以根据染色体自身的情况判断是否进行交叉,在保留了优势染色体的同时对适应度较低的染色体进行改进,从而达到促进整个种群进化的效果。
S304、变异算子;
采用基本位变异的操作对染色体进行变异,使用相关性计算中得到的motif相关性统计指标作为先验知识影响变异过程,以该指标的绝对值作为每个位点(motif)的突变频率,相关性较高的突变频率更高,且为了限制寻优结果中包含过多的无关motif,通过判断已入选motif的个数来判断突变的类型,如染色体C1包含motif小于阈值p,则从其未选中的motif中选择变异并对相应位点置‘1’,染色体C2包含motif大于阈值q,则从其已选中的motif中选择变异并对相应位点置‘0’,其余情况则随机变异。变异也采用自适应变异率来确定是否进行变异操作,同时选择突变基因位时采用轮盘置操作,使得突变频率更高的位点更容易发生变异操作,生成新个体,扩展种群丰富度,具体步骤如下:
S3041、自适应交叉率确定是否进行变异操作;
其中,fmax代表当前种群中的最大值,f'代表要变异的染色体适应度值,favg代表种群的平均适应度值,k2代表计算常数,取k2=0.1。
S3042、根据染色体包含motif数量选择变异类型;
S3043、对每一条染色体,统计所有变异位点对应的相关系数,根据系数比例计算每个位点对应的变异概率p(xi),同时计算累计概率q(xi),即每个个体之前所有个体的选择概率之和
S3044、随机生成s,范围[0,1],若q(xi-1)<s<q(xi),则选择位点xi进行变异操作。
S4、相关性分析
输入数据的矩阵中存在一些motif,他们只出现在少部分样本中,具有极高的样本特异性,如果直接在两个矩阵上进行同步寻优,这些motif片段的计算不会提高适应度函数值,同时在后续寻优算法的实现中还浪费了计算资源。因此,在进行遗传算法寻优之前首先对矩阵中的每一motif列和label列进行相关性分析,采用点二列相关性作为衡量和label的指标,点二列相关的一列变量是连续变量,另一列是真实二分变量,即样本性质,计算方式具体如下:
在筛选过程中通过假设检验结果ρ值判断显著水平,取ρ=0.05作为划分阈值,即ρ>0.05则认为相关性不显著,将对应motif过滤。同时该相关性指标也将作为变异操作的先验信息影响变异过程。
S5、过滤算法输出
通过遗传算法的迭代计算,得到两组种群的优势解,取两组种群中适应度排前三的染色体作为最优解,并解码成为对应的motif集合,最终得到六组motif集合,两两对应,取对应集合的交集作为算法最终挖掘出的motif,并分析结果,具体流程如图3所示。
为使本发明实施例的目的、技术方案和优点更加清楚,下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例是本发明一部分实施例,而不是全部的实施例。通常在此处附图中的描述和所示的本发明实施例的组件可以通过各种不同的配置来布置和设计。因此,以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围,而是仅仅表示本发明的选定实施例。基于本发明中的实施例,本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
下面通过半仿真数据和真实数据实验及对比结果证明本方法的实用性、准确性和高效性:
1)半仿真数据生成、评估指标设计和实验结果
在原数据的基础上生成新的数据,目的在于扩大数据样本以评估算法的稳定性。半仿真数据的具体生成步骤如下:以样本数N=1000为例,首先随机生成一组取值在{0,1}之间且数量为1000的数组代表样本的标签,n_tumor代表肿瘤样本数量,n_health代表健康样本数量,其次在n种motif中随机选取m个作为预先设置的motif集合,根据真实样本的统计数据范围,分别以随机数生成的方式在blood、tissue和health中填入数据,而对于未选中的motif,则根据health和tumor的综合统计范围分别生成对应的随机数,这样得到了用于寻优的两组数据。同时通过控制n_sample样本总数,Kmotif组合中的motif数量两个参数来生成多组实验数据,其中n_sample取值范围如下:{1000,1500,2000,2500,3000},K值则根据已有研究结果设置为{30,40,50,60,70},共生成25组对应的算法输入数据。
对于半仿真数据选择统计学中的I、II类错误率作为评估指标,其中,I类错误指的是在原假设H0为真时否定原假设,II类错误指的是在原假设为假时未能否定原假设。假设H0:该motif与肿瘤之间无显著性关系,则I类错误具体指的是错误的认为该motif与肿瘤之间有关系(实际无关),即解空间中包含了与肿瘤无关的motif;II类错误具体指的是错误的认为该motif与甲状腺癌之间无关系(实际有关),即解空间中遗漏了部分与肿瘤有关的motif。由此,定义I类错误率的计算公式:
其中,Selected_error_num代表与肿瘤无关的motif被错误的认为相关的数量,Selected_all代表由实验得出的与肿瘤相关的motif的总数量。II类错误率的计算公式如下:
其中,Selected_miss代表与肿瘤相关的motif被错误的认为无关的数量,Ans代表预先设定的相关motif总数量。
经过迭代计算,得到25组数据的结果如下:
由实验结果可知,在所有数据集上,与II类错误率相比,I类错误率相对较低,一半以上的数据结果显示I类错误率达到0。在寻优的过程中,目的是在尽可能降低未知motif包含在染色体解中的情况下尽可能多的识别正确的motif,原因是如果寻优结果中包含无关motif,可能会导致错误的肿瘤免疫研究生物标记物方向,而如果识别结果中包含的motif大多为肿瘤相关motif,即使丢失了一部分,寻优找到的这一部分也具有一定意义。实验结果符合预期目标。
请参阅图4,图中可以更明显的表现I类错误和II类错误之间的差距,在不同的K值下,I类错误率均达到相对较低的水平,其中最高点在K=40,N=2500时取到,I类错误率达到10.85%,所有实验数据的平均I类错误率为2.55%,虽然存在部分数据中最后的寻优结果里仍然有错误的motif,这可能与根据实际数据统计结果随机生成数据有关;II类错误率的水平维持在12.50%~38.67%,最高点在K=50,N=2500时取到,所有实验数据的平均II类错误率为24.53%。
目前在样本间T细胞受体序列重叠研究中通常采用T细胞受体序列或片段作为特征构建模型,由于维度较高的特点,研究者会做特征选择,即通过计算序列或片段的指标作为特征筛选的参考,该过程会选择出重要的序列或片段。为了突出本发明和使用数据的有效性,将BOW_SVM和1-DBC运用于肿瘤样本和健康样本外周血数据上,同时将结果与MDPGA结果进行比对,以Type II error作为比对指标,比较结果如图5所示。
图5中的实验结果表明,BOW_SVM和1-DBC算法的二类错误率比较高,集中在70%~80%之间,而MDPGA的二类错误率集中在17%~30%之间,效果优于对比算法,由此对比结果可以证明本发明在肿瘤特异性motif的挖掘上使用两组对比数据的效果更优,同时在对比数据上的效果比单一数据在其他算法上的效果更好,证明了算法和数据使用的有效性。
2)真实数据检验指标和实验结果
真实数据包含包括85个甲状腺癌样本的组织和血液T细胞受体测序数据以及260个健康样本血液T细胞受体序列数据。按照方法中所述的数据预处理方式,将所有样本数据处理为两个用于寻优的矩阵,然后采用算法的基本参数设置,经过1500次迭代后通过筛选得到11个与肿瘤相关的motif,为了验证这些motif的有效性,采用Mann–Whitney U Test检测该motif下不同性质样本的差异,以p<0.05为判断条件,得到的效果均比较显著。检验效果如下:
请参阅图6,表示不同性质样本中不同motif平均出现次数,可以看出组合中的motif在肿瘤样本中出现的次数更高,同时其出现在组织中的次数也比血液中高,寻优结果符合预期。
综上所述,本发明一种T细胞受体序列motif组合识别检测方法、存储介质及设备,使用的数据不仅包含肿瘤与非肿瘤对照组的数据,同时基于免疫异质性的特点,通过对比肿瘤患者组织和外周血的T细胞受体数据进一步识别肿瘤特异性motif组合。具体通过T细胞受体序列剪切和双种群遗传算法的过程得出组合结果,本发明在仿真数据和真实数据集上均得到验证。
本领域内的技术人员应明白,本申请的实施例可提供为方法、系统、或计算机程序产品。因此,本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且,本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本申请是参照根据本申请实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器,使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中,使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品,该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上,使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理,从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上内容仅为说明本发明的技术思想,不能以此限定本发明的保护范围,凡是按照本发明提出的技术思想,在技术方案基础上所做的任何改动,均落入本发明权利要求书的保护范围之内。
Claims (10)
1.一种T细胞受体序列motif组合识别检测方法,其特征在于,包括以下步骤:
S1、将健康样本和肿瘤样本的基因组DNA比对到胚系序列,鉴定对应的CDR3区段并将CDR3区段翻译成氨基酸序列;进行氨基酸序列剪切并将健康样本和肿瘤样本的CDR3区段对应的氨基酸序列解构为多个相同长度的连续氨基酸片段;利用健康样本的外周血数据、肿瘤样本的组织数据和肿瘤样本外周血数据解构后的信息构建Tumor-Health矩阵,利用肿瘤样本组织数据和肿瘤样本外周血数据构建Tissue-Blood矩阵,Tumor-Health矩阵和Tissue-Blood矩阵的最后一列为label列,代表样本对应的标签,剩余列为motif列,代表样本对应氨基酸片段的克隆数,将Tumor-Health矩阵和Tissue-Blood矩阵共同作为输入矩阵;
S2、根据步骤S1得出的两个输入矩阵分别构建求解对应的优化目标函数,以motif组合的类内类间距离可分性作为优化目标计算方法,将两个优化目标函数加和得到总的优化目标函数;
S3、根据步骤S1得到的输入矩阵和步骤S2得到的总优化目标函数设计双种群遗传算法并进行求解,确定编码方法和初始种群,采用锦标赛选择法加精英保留制度选择算子对种群进行选择,保留最优个体;采取种群内和种群间的混合交叉方式确定交叉算子;采用基本位变异的操作确定变异算子;
S4、对步骤S1构建的输入矩阵中的每一motif列和label列进行相关性分析,采用点二列相关性作为衡量motif列和label列是否相关的指标,根据计算得到的相关性指标ρ判断显著水平,确定相关性指标ρ>0.05为相关性不显著;
S5、以步骤S2构建的优化目标函数作为算法优化目标,利用步骤S4计算的相关性指标对步骤S1构建的输入矩阵进行初步筛选,过滤掉不相关的motif列,通过步骤S3设计的双种群遗传算法进行迭代计算,得到两组种群的优势解,取两组种群中适应度排前三的染色体作为最优解,并解码成为对应的motif集合,最终得到六组motif集合,两两对应,取对应集合的交集作为最终挖掘出的motif,完成识别检测。
2.根据权利要求1所述的方法,其特征在于,步骤S1中,构建Tumor-Health矩阵和Tissue-Blood矩阵具体为:
肿瘤样本数为M,健康样本数为N,肿瘤样本外周血和组织测序的和为2M,第一个用于输入的矩阵Tumor-Health共有N+2M行,代表肿瘤样本外周血、组织的测序结果和健康样本外周血测序结果,第二个输入的矩阵Tissue-Blood共有2M行,代表肿瘤样本外周血和组织的测序结果,两个矩阵的每一行均包含一个向量C={C1,C2,...Cs,label},其中,s=8000,Ci代表样本中motif出现的计数,如果Ci=0则代表样本无对应motif,label代表样本对应的标签,在Tumor-Health矩阵中0代表正常的样本数据,1代表非正常的样本数据,在Tissue-blood矩阵中0代表非正常外周血数据,1代表非正常组织数据。
3.根据权利要求1所述的方法,其特征在于,步骤S2中,总优化目标函数min F(x)为:
min F(x)=F1(x)+F2(x)
其中,F1(x)为优化目标1,在Tumor-Health矩阵上进行类内类间距离可分计算,F2(x)为优化目标2,在Tissue-blood矩阵上进行类内类间距离可分计算。
4.根据权利要求1所述的方法,其特征在于,步骤S3中,确定编码方法和初始种群具体为:
采用二进制编码,使用二值符号集{0,1}构成种群中每个染色体的基因型,每一个二进制位对应一种motif,1表示染色体对应的可行解中包括对应motif,0表示不包含对应motif;对于两组数据采用随机生成的方式分别生成两个种群,然后通过交叉让两组数据中的优势基因进行互换,完成初始化。
5.根据权利要求1所述的方法,其特征在于,步骤S3中,采用锦标赛选择法加精英保留制度选择算子对种群进行筛选具体为:
S3021、计算种群中所有染色体的适应度值,将最好的染色体直接保留至下一代;
S3022、确定每次选择的个体数量k;
S3023、从种群中随机选择k个个体构成组,根据其适应度值选择最好的个体进入子代种群;
S3024、重复操作直到新一代种群中染色体个数达到要求,选择算子筛选产生一组新的种群解用于混合交叉过程计算。
7.根据权利要求1所述的方法,其特征在于,步骤S3中,确定变异算子具体为:
S3041、比较随机数与自适应变异率判断对应染色体是否需要进行变异,如果变异率大则进行变异,自适应变异率计算如下:
其中,fmax代表当前种群中的最大值,f'代表要变异的染色体适应度值,favg代表种群的平均适应度值,k2代表计算常数,取k2=0.1;
S3042、根据染色体包含motif数量选择变异类型
S3043、对每一条染色体,统计所有变异位点对应的相关系数,根据系数比例计算每个位点对应的变异概率p(xi),同时计算累计概率q(xi),即每个个体之前所有个体的选择概率之和:
S3044、随机生成s,范围[0,1],若q(xi-1)<s<q(xi),选择位点xi作为变异算子进行变异算子操作,经过变异算子操作更新种群中的染色体得到新个体。
9.一种存储一个或多个程序的计算机可读存储介质,其特征在于,所述一个或多个程序包括指令,所述指令当由计算设备执行时,使得所述计算设备执行根据权利要求1至8所述的方法中的任一方法。
10.一种计算设备,其特征在于,包括:
一个或多个处理器、存储器及一个或多个程序,其中一个或多个程序存储在所述存储器中并被配置为所述一个或多个处理器执行,所述一个或多个程序包括用于执行根据权利要求1至8所述的方法中的任一方法的指令。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110536816.6A CN113380324B (zh) | 2021-05-17 | 2021-05-17 | 一种T细胞受体序列motif组合识别检测方法、存储介质及设备 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110536816.6A CN113380324B (zh) | 2021-05-17 | 2021-05-17 | 一种T细胞受体序列motif组合识别检测方法、存储介质及设备 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN113380324A true CN113380324A (zh) | 2021-09-10 |
CN113380324B CN113380324B (zh) | 2023-06-27 |
Family
ID=77571199
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110536816.6A Active CN113380324B (zh) | 2021-05-17 | 2021-05-17 | 一种T细胞受体序列motif组合识别检测方法、存储介质及设备 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN113380324B (zh) |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20080167850A1 (en) * | 2007-01-10 | 2008-07-10 | International Business Machines Corporation | Method and apparatus for detecting consensus motifs in data sequences |
US20100016748A1 (en) * | 2008-07-17 | 2010-01-21 | Syed Zeeshan H | Motif Discovery in Physiological Datasets: A Methodology for Inferring Predictive Elements |
CN108319818A (zh) * | 2018-02-07 | 2018-07-24 | 中国科学院生物物理研究所 | 一种预测影响长非编码rna生物学功能的snp位点的方法 |
CN110309822A (zh) * | 2019-06-18 | 2019-10-08 | 哈尔滨工程大学 | 基于量子进化粒子群算法的高光谱图像波段选择方法 |
CN111511911A (zh) * | 2017-12-24 | 2020-08-07 | 诺伊尔免疫生物科技株式会社 | 表达特异性地识别人间皮素的细胞表面分子、il-7和ccl19的免疫活性细胞 |
CN111914470A (zh) * | 2020-06-17 | 2020-11-10 | 西安交通大学 | 一种能源化工生产系统多元监测时间序列回归预测方法 |
CN112307678A (zh) * | 2020-11-05 | 2021-02-02 | 湖南科技大学 | 基于混沌非支配排序遗传算法的机器人多目标搜索方法 |
-
2021
- 2021-05-17 CN CN202110536816.6A patent/CN113380324B/zh active Active
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20080167850A1 (en) * | 2007-01-10 | 2008-07-10 | International Business Machines Corporation | Method and apparatus for detecting consensus motifs in data sequences |
US20100016748A1 (en) * | 2008-07-17 | 2010-01-21 | Syed Zeeshan H | Motif Discovery in Physiological Datasets: A Methodology for Inferring Predictive Elements |
CN111511911A (zh) * | 2017-12-24 | 2020-08-07 | 诺伊尔免疫生物科技株式会社 | 表达特异性地识别人间皮素的细胞表面分子、il-7和ccl19的免疫活性细胞 |
CN108319818A (zh) * | 2018-02-07 | 2018-07-24 | 中国科学院生物物理研究所 | 一种预测影响长非编码rna生物学功能的snp位点的方法 |
CN110309822A (zh) * | 2019-06-18 | 2019-10-08 | 哈尔滨工程大学 | 基于量子进化粒子群算法的高光谱图像波段选择方法 |
CN111914470A (zh) * | 2020-06-17 | 2020-11-10 | 西安交通大学 | 一种能源化工生产系统多元监测时间序列回归预测方法 |
CN112307678A (zh) * | 2020-11-05 | 2021-02-02 | 湖南科技大学 | 基于混沌非支配排序遗传算法的机器人多目标搜索方法 |
Non-Patent Citations (2)
Title |
---|
AL MUTTAKIN ETAL,: "Motif discovery in unaligned DNA sequences using genetic algorithm", 《HTTPS://WWW.RESEARCHGATE.NET/PUBLICATION/322515123》 * |
蒙冰等: "基于非支配排序遗传算法求解启动子识别问题", 《基因组学与应用生物学》 * |
Also Published As
Publication number | Publication date |
---|---|
CN113380324B (zh) | 2023-06-27 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Fischer et al. | High-definition reconstruction of clonal composition in cancer | |
Chen et al. | Random forests for genomic data analysis | |
CN109767810B (zh) | 高通量测序数据分析方法及装置 | |
US20220130488A1 (en) | Methods for detecting copy-number variations in next-generation sequencing | |
CN109994151B (zh) | 基于复杂网络与机器学习方法的肿瘤驱动基因预测系统 | |
CN106033502B (zh) | 鉴定病毒的方法和装置 | |
US20220310199A1 (en) | Methods for identifying chromosomal spatial instability such as homologous repair deficiency in low coverage next- generation sequencing data | |
CN111462823A (zh) | 一种基于dna测序数据的同源重组缺陷判定方法 | |
WO2024027032A1 (zh) | 一种肿瘤形成风险与肿瘤组织来源的评估方法及系统 | |
Wijethilake et al. | Survival prediction and risk estimation of Glioma patients using mRNA expressions | |
Guo et al. | Genome-wide interaction-based association of human diseases-a survey | |
CN108460248B (zh) | 一种基于Bionano平台检测长串联重复序列的方法 | |
Yang et al. | Catfish Taguchi-based binary differential evolution algorithm for analyzing single nucleotide polymorphism interactions in chronic dialysis | |
US20230073973A1 (en) | Deep learning based system and method for prediction of alternative polyadenylation site | |
CN113380324A (zh) | 一种T细胞受体序列motif组合识别检测方法、存储介质及设备 | |
US20190108311A1 (en) | Site-specific noise model for targeted sequencing | |
Riley et al. | Interpreting generative adversarial networks to infer natural selection from genetic data | |
CN111951889B (zh) | 一种rna序列中m5c位点的识别预测方法及系统 | |
CN115035951A (zh) | 一种突变签名的预测方法、装置、终端设备及存储介质 | |
KR102376212B1 (ko) | 신경망 기반의 유전자 선택 알고리즘을 이용한 유전자 발현 마커 선별 방법 | |
US20230407405A1 (en) | Method for diagnosing cancer and predicting type of cancer based on single nucleotide variant in cell-free dna | |
CN114242158B (zh) | ctDNA单核苷酸变异位点检测方法、装置、存储介质及设备 | |
Durvasula et al. | Recovering signals of ghost archaic admixture in the genomes of | |
CN117746982A (zh) | 基于外源dna插入突变导致的优势克隆的检测方法 | |
CN115713107A (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 |