本发明涉及微生物检测模型领域，具体的说是一种单颗牙齿微生物群落检测模型的构建方法及其应用。
背景技术：
寄生于人体的大量共生菌群作为后天禀赋的主要承载者与人体的健康息息相关。目前认为共生微生物可作为人体的第二基因组，其遗传信息的总和被称为微生物组(microbiome)，赋予人类不依赖于自身进化而获得的复杂个体特征。因此，全面认识人体共生菌群的基础上，才能深度揭示其对人体健康或疾病状态的影响，从而构建微生物群落存在及变化情况与特定的疾病之间的联系。
口腔系统是连通人体内外的交通枢纽，为人体共生菌群非常重要的栖息位点，其微生物群落结构复杂，这些微生物以细菌为主，也有少量真菌、病毒、古细菌的存在。而由于不同牙齿位点解剖结构、功能、氧化梯度、ph等差异，造成了口腔不同牙齿具有显著不同的菌群。口腔菌群失调与口腔疾病如龋齿、牙龈炎、牙周炎等发生发展有关，甚至与全身性疾病如糖尿病、心血管疾病、胰腺癌等具有显著相关。维持口腔菌群结构和功能的健康平衡状态，对于人体健康具有深刻而不容忽视的重大意义。此外，与血液和尿液作为疾病诊断媒介相比，口腔位点采样具有低侵害性、低成本、样品采集和处理简易、快捷等优势。
现有的检测模型构建的数据通常是通过横向研究方法而获得，即比较健康人群和疾病人群的特征而获得特定特征健康和疾病差异。但慢性疾病发展的不同阶段存在变异性；而疾病的同一发展阶段在不同个体间亦存在一定异质性，如可受到宿主基因、免疫、生活习惯等因素影响。对于疾病诊断采样通常集中于单一位点，如血液样品或尿液样品等。但不同位点提供的有效信息存在差异，目前对于不同采样位点对于同一种疾病诊断的是否具有优劣性，是否不同位点对于疾病诊断存在互补性均较少考察。更为重要的是，目前基于全口菌群建立模型的方法可以判定宿主整体的健康疾病状态，但却无法判断出具体是疾病状态存在在哪颗牙齿，更无法预测单个牙齿疾病发生的概率，即无法判定疾病高风险宿主具体哪颗牙齿即将发生疾病。因此，多数诊断方法和模型止步于整体口腔生物因子的筛选后简单建模，未将其精细化和精准化，限制了其临床转化和应用前景。因此，好的诊断方法应精准考察精细位点状态，使其自身具备高敏感性、高特异性，同时在临床应用上应满足高通量、易操作、无侵害性和低成本等要求。
技术实现要素：
针对现有技术中存在的上述不足之处，本发明要解决的技术问题是提供一种口腔微生物群落检测模型的构建方法及其应用。
本发明为实现上述目的所采用的技术方案是：一种单颗牙位的微生物群落检测模型的构建方法，包括以下步骤：一种单牙位微生物群落检测模型的构建方法，包括以下步骤：
数据收集：在待收集人群中收集每颗牙齿各自的口腔微生物群落数据；
模型的初步构建：以口腔微生物群落数据为样本，通过机器学习算法构建模型；
模型的最终构建：根据初步模型获得每个微生物因子对模型性能的贡献程度，将其进行降低排列，并建立基于不同个数微生物因子组合的所有模型，进而筛选出对于模型性能贡献程度最大的一组微生物因子，作为最终的单牙位微生物群落检测模型。
所述数据收集为在待收集人群中区分出患有牙病人群和未患有牙病人群中每颗独立牙齿在的口腔微生物群落数据。
所述口腔微生物群落数据为高通量测序手段获得微生物群落的16srna或全基因组信息。
所述模型的初步构建，包括以下步骤：
(1)以口腔微生物群落数据作为输入；
(2)样本的状态作为输出，包括微生物标记物的宿主状态，某种状态的可能性；
(3)采用随机森林方法得到初步检测模型；
(4)获取每个微生物因子对于模型性能的贡献程度。
所述模型的最终构建，包括以下步骤：
(1)将初步模型构建后获得的各个微生物因子对模型性能的贡献程度进行降序排列；
(2)根据降序排列结果，取前n个贡献程度值对应的微生物因子进行组合，建立基于不同微生物因子组合的所有模型；n为自然数且小于微生物因子个数；
(3)筛选出对模型贡献程度值低于阈值并且准确率低于设定值的微生物因子作为无关微生物因子；
(4)从微生物组数据中剔除无关生物因子；并删除无关微生物因子对应的模型；
(5)建立以剩余贡献程度的微生物因子为输入的单牙位微生物群落检测模型。
所述样本类别表示微生物标记的宿主的状态，包括宿主此时的正常、此时的异常、未来的正常、未来的异常。
所属微生物因子为微生物物种信息。
所述模型性能是指检测模型区分样本状态的准确率。
所述贡献程度，是指某一种微生物因子对于区分样本状态的能力，体现其对于模型准确率的贡献程度。
所述模型用作检测口腔疾病发生的可能性评价。
本发明具有以下优点及有益效果：
1.模型对象采集处理简易、无侵害性、成本低、检测对象相对简单；
2.建模方法易于操作、数据处理高效、检测能力强；
3.模型应用广泛：其应用对象不仅适用于大规模人群筛选，也可针对个体单牙位实现精准地终身监测；其应用形式不仅可检测口腔微生物宿主此时状态，也可预测宿主未来可能发生状态。
附图说明
图1为本发明方法的实施流程图；
图2为本发明实施提供的实验设计图；
图3为本发明实施提供的口腔微生物群落结构特征图；
图4为微生物因子数量与模型准确率关系图。
图5为本发明实施提供的通过随机森林回归方法展示不同微生物因子数量情况下，si-mic模型的性能示意图。
图6为本发明实施提供的基于随机森林二分类方法最终的的疾病状态检测模型的准确率示意图。
图7为本发明实施提供的基于随机森林二分类方法最终的单牙位菌群检测模型(si-mic)的准确率示意图。
具体实施方式
下面结合附图及实施例对本发明做进一步的详细说明。
数据收集：在待收集人群中收集所有单颗牙齿微生物群落数据，并随着待收集人群的年龄增长，收集多个时间点的口腔微生物群落数据，并在所述待收集人群中区分出在患有牙病人群和未患有牙病人群在多个时间点上的单牙位微生物群落数据；
模型的初步构建：利用机器学习随机森算法得到的口腔微生物群落数据，初步构建口腔微生物群落检测模型；
模型的最终建立：获得对模型准确率贡献程度最大的一组微生物因子，最终建立单牙位龋病微生物指数(single-toothmicrobialindicatorofcaries，si-mic)模型。准确率＝(真阳性样本个数 真阴性样本个数)/样本总个数，真阳性样本指的是实际为阳性且预测为阳性的样本，真阴性样本指的实际为阴性且预测为阴性的样本。如图4所示。
数据筛选：在所述未患有牙病人群在多个时间点上的口腔微生物群落数据中，筛选出于对模型效能无作用、甚至负作用的微生物群落；
数据剔除：在所述患有牙病人群在多个时间点上的口腔微生物群落数据中，剔除数据筛选步骤得到的微生物群落的数据成分；
所述口腔微生物群落数据为高通量测序手段获得微生物群落的16srna或全基因组信息。
数据筛选，具体为：利用随机森林回归分析方法处理健康组健康牙齿和疾病发生组中龋病牙齿多时间点牙菌斑微生物群落数据，分别筛选出不具备区分宿主状态能力或者较低区分能力的微生物标记物。
所述模型的初步构建，包括以下步骤：
(1)去除数据筛选步骤所得的于宿主状态无关的微生物数据；
(2)利用随机森林二分类方法处理已进行剔除步骤的患有牙病人群单牙位多时间点的口腔微生物群落数据，建立基于所有微生物组群落信息的初步的检测模型(si-mic)。
所述模型的最终确定，包括以下步骤：
(1)按照每个微生物因子对于单牙位健康疾病状态的区分能力排序；
(2)逐步加入微生物因子数量，分别建立不同数量微生物因子的检测模型；
(3)在不降低模型性能前提获得最少数量的微生物因子组合，从而确定最终模型。
所述种微生物因子包括所有可注释的种水平微生物。
所述模型用作检测口腔疾病发生的可能性评价。
本发明利用单颗牙齿微生物群落检测模型初步构建和最终确立儿童单颗牙齿状态的检测模型作为实施例(图1)，包括下列内容：
(1)收集儿童口腔健康状态临床信息：
对青岛市爱弥儿幼儿园全日制儿童的口腔健康进行追踪调查，每两个检查一次，持续四次检查，之后再间隔五个月进行检查，根据五次调查记录的儿童dmfs(龋，失，补牙数)指数，根据本研究目的选择具有下述两类口腔健康变化特征的儿童纳入此课题研究：1)健康对照组(h2h组):口腔龋病状况始终保持健康组；2)健康宿主龋病发生组(h2c组):健康宿主口腔状况经历从健康到龋病新发过程；3)龋病宿主龋病发生组(c2c组):：龋病宿主单牙位经历健康宿主口腔状况经历从健康到龋病新发过程。入选标准包括：年龄约5岁，20颗乳牙全部萌出，排除标准包括：有全身系统性疾病和牙周、口臭等口腔疾患，三个月未服用抗生素。选取所有入选儿童的三次口腔检查时所取的龈上牙菌斑和唾液样品共计284个，按照dmfs指数进行如下分类：低龋样品(lc,0<dmfs<6)，重症婴幼儿龋样品(sc,dmfs≥6)，健康对照样品(h,dmfs＝0)。其中，龋病状态(c)是绝对，健康(h)有相对和绝对区分。具体是：龋病样本组包括：h2c-c和c2c-c；而对于观测时健康牙位，根据其样本分组、现在和未来状态，可将其分为以下三种情况：1)绝对健康样本(h2h-hh)：在健康组中，绝对健康样品是指两个时间点均为健康的前一个时间点的样本；对应地，2)相对健康样品：①龋病发生组中保持健康状态样本(h2c-hh和c2c-hh)；②龋病发生组中经历从健康到疾病状态的样本(h2c-hc和c2c-hc)。就整个实验流程各项细节及以后的数据公布等事宜征得志愿者监护人同意，并签署知情同意书。
调查方法:由两名牙体牙髓专科医生以视诊结合探诊的方式进行检查，检查器械高温高压消毒，必要时借助棉签去除软垢。检查前统一认识、方法和标准，标准一致性检验的kappa值均大于0.92。采用世界卫生组织《口腔健康调查基本方法》(1997年)对龋病的诊断标准。冠龋诊断标准：牙齿的窝沟点隙或光滑面有明显龋洞、或明显釉质下破坏、或明确可探及软化洞底或洞壁的病损记为龋齿，包括有充填物或已窝沟封闭同时有龋者。有下列表现而缺乏其他阳性症状时不列入龋齿记录范围：①白色或白垩色斑点；②探诊无软化的着色或粗糙斑点；③釉质点隙或窝沟着色，但无明显釉质下潜行破坏；④中到重度氟斑牙，有光泽、质硬、有小凹陷；⑤根据分布或病史，结合触诊、视诊观察因磨损而造成病损龋齿。
(2)收集单个牙位的龈上菌斑样本：
取样前一小时受试者避免进食及饮水，每次取样均在早上9:100-12:00，取样时儿童保持轻仰头、闭眼、直立座位。使用无菌牙刷采集萌出的上乳前牙以及所有乳后牙龈上菌斑1分钟，将粘附于牙刷上菌斑转移至盛有10ml双蒸水的50ml离心管，取样时避免触碰黏膜等口腔其他位点。对样品分别编号并置于-80℃保存待提取dna。
(3)基因组dna提取和pcr扩增16srrna基因片段
采用dneasybloodtissuekit试剂盒法，分别提纯分离样本中细菌基因组dna。具体步骤如下：1)加入180ulenzymaticlysisbuffer(溶菌酶缓冲液)中重悬，37℃，孵育30min。2)加入25ul蛋白酶k,200ulbufferal56℃,孵育30min。3)加入200ul乙醇(96％～100％)涡旋至均匀。4)移取以上步骤的混合液至离心柱上，离心柱放在2ml收集管上(8000rpm离心1分钟)弃收集管/液。5)离心柱放在新的2ml收集管上，加入500ul缓冲液aw1>6000xg(8000rpm1分钟，弃收集管/液。6)离心柱放在新的2ml收集管上，加入500ul缓冲液aw220000xg(14000rpm)3分钟，弃收集管/液。7)离心柱放在新的1.5-2ml收集管上，200ulbufferae加入到离心柱中心膜上，室温孵育1分钟.8)重复步骤7,以增加dna的产量。采用qubit超微量分光光度仪定量dna浓度，电泳检测dna完整性。提取后的dna保存于-20℃。约15ngdna用于构建16s扩增文库。
为获得相对准确的种系发育信息，选取16srrna片段上v3-v4高变区作为pcr扩增目标片段。确定pcr上游引物(5’-nnnnnnn-gtactcctacgggaggcagca-3’)及下游引物(5’-nnnnnnn-gtggactachvgggtwtctaat-3’)，nnnnnnn即idtag，是为区别不同样品来源而设计的随机组合的七个碱基，分别加入上下游引物的5’端，利用该多样品平行标记技术完成多个样品同时在测序仪上测序。
每个样品进行三次pcr扩增，pcr反应体系(25μl)包含：1)12.5μl的gotaghotstart聚合酶,2)1μl上游引物(浓度5pm)，3)下游引物(浓度5pm)，4)1μl基因组dna(5ngμl-1)，5)9.5μlpcr级别无菌水。整个反应体系需要在在thermocyclerpcrsystem进行反应。反应条件设定为：95℃预变性2min，94℃变性30s，退火56℃25s，72℃延伸25s，共25个循环，最后72℃延伸5min[4]。
(4)illuminami-seq测序
(5)将获得的高通量数据转换成具体的微生物群落数据
序列质量控制：454高质量序列分析流程主要基于qiimeii平台，设定质量控制规范，符合标准的序列片段被视为高质量序列，予以保留。基于16s数据库的种系发育信息分析：采用qiimeii分类方法针对人类口腔核心微生物16s数据库(core)进行从门到种水平细菌种系信息划归，分别统计各个样品在每个分类水平上各物种的序列数，并与该样品总体获得的序列数计算比值，从而获取每个门类各物种的相对丰度。
(6)不同因素对于口腔菌群分布的影响：
以jensen-shannon矩阵为基础的群落结构计算方法：其除了样品间的进化距离外，还可调查样品菌种水平物种丰度或上的区别。样品中细菌种的丰度分布可以看作是物种的概率分布，可以利用样品间这种概率分布的互信息熵(jensen-shannondivergence，jsd)来度量样品间的微生物组的区别。样品间的距离d(a,b)的计算公式如下：
pa和pb分别代表样品a和样品b中的丰度分布。jsd(x,y)定义了两个样品中不同的概率分布x和y间的互信息熵(jensen-shannondivergence)。
kld是x和y间的kullback-leibler离散度，具体的计算方法如下：
xi表示某个i样品的相对丰度，yi表示其他非i样品的相对丰度。
非监督的主坐标分析：将未加权的jensen-shannon矩阵进行主坐标分析(pcoa:principalcoordinatesanalysis)以展示不同样本间差异，pcoa将各个物种信息视为互相独立不关联的变量，以样本×变量相对丰度的矩阵进行分析，以在不考虑环境因子影响的前提下，无偏见、整体的观察样本的内在菌群结果，发现一个或多个潜在的变量(主坐标，principalcoordinate，pc)以最大程度的在较低维度上最好的解释样本内在的变异，每一个主坐标代表在此维度下可解释的整体结构变异程度，从而达到数据降维处理并对样品排序的目的，其中样本的得分(score)是物种得分的线性组合。
置换多元统计分析结果显示决定口腔微生物菌群分布的重要因素，结果提示我们(图1,图2)：
①前后位点是决定菌群分布的最重要因素，即无论在疾病还是健康个体内或个体间，前牙和后牙牙菌斑菌群物种组成显著不同(图1)。
②其他因素的对菌群结构的影响程度，降序顺序为：单牙位龋病状态(包括牙齿的健康状态、牙齿患龋的严重程度)、个体异质性、宿主龋病状态(包括牙齿的健康状态、牙齿患龋的严重程度)(图1)，但性别无影响。
在进一步根据宿主以及宿主单颗牙位现在和未来状态进行分组，发现来自健康组的健康牙齿(h2h-hh)菌群不同于疾病发生组的健康牙齿菌群(h2c-hh和c2c-hh)，且来自疾病发生组中未来会发生疾病的健康牙齿菌群(h2c-hc和c2c-hc)更趋向于该组中的疾病牙齿样本菌群(c2c-c和h2c-c)(图2)：
(7)初步建立单颗牙位的口腔状态检测模型的方法
在机器学习中，随机森林方法是一个包含多个决策树的分类器，并且其输出的类别是由个别树输出的类别的众数而定。随机森林方法不但可以建立分类模型，同时可确定区分特定状态或标签的变量，并可通过其贡献程度(importancescore)以判断其区分能力的大小。在本实例中，随机森林方法利用统计画图软件r(https://www.r-project.org/)的randomforest软件包实现，建立5000棵树，其他均为默认设置。以输入数据的2/3作为训练数据集，以输入数据的1/3作为测试数据集，随机进行100次实验以降低误差。具体的是：
利用机械学习的随机森林方法，以龋病组患龋牙齿(h2c-c和c2c组中)样本为疾病菌群和健康组绝对健康样本(h2h-hh)为健康菌群，将其二者所有微生物组信息作为输入变量，将其与绝对健康样本(h2h-hh)和龋病样本(c2c-c和h2c-c)二分类的输出变量对应，初步建立分类模型，获得可甄别单颗牙齿当前状态的单牙位菌群模型(single-toothmicrobialindictorofcaries,si-mic)。其中，分类模型si-mic结果以对于样本状态的分类结果形式在：1)当指数结果＝50％时，则判定为随机事件，无法进一步判断单个牙位状态；2)当指数结果>50％时，则判定为阳性事件，可以进一步判定单个牙位此时状态为龋病；3)当指数结果<50％时,则判定为阴性事件，可以进一步判定单个牙位此时状态为健康。建立初步模型同时，还可以获得所有变量(微生物因子)对于模型性能/准确性的贡献程度。
(8)筛选低区分能力/重要程度的微生物标记物：
根据变量其对模型重要性排序，进一步变量数量逐步增加，评估不同si-mic模型的性能。其中，si-mic的总数量应为变量的总数(图3)。
(9)基于不同微生物因子的si-mic模型的方法及其性能
根据不同变量数量组合的模型性能，优选出随着变量减少而模型性能未显著改变的14个变量组合作为最终模型变量(图5)，并获得了高到93.4％准确率的模型(图6)。
(10)最终si-mic模型对单牙位疾病风险程度/未来状态的预测
以si-mic为模型，相对健康样本中的未来发生疾病的相对加健康样本(h2c-hc和c2c-hc)作为带入变量，获得评估状态的可能性指数。其中，1)当指数结果＝50％时，则判定为随机事件，无法进一步判断单个牙位状态；2)当指数结果>50％时，则判定为阳性事件，可以认为未来风险预测正确；3)当指数结果<50％时,则判定为阴性事件，可以进一步判定单个牙位未来状态为风险预测错误。最终的模型可以获得高达83.9％的预测准确率(图7)。
本发明所述的随机森林回归分析方法和二分类方法可参见breimanl(2001)randomforests.machlearn45:5–32.)和(knightsd,costelloek,knightr.supervisedclassificationofhumanmicrobiota.femsmicrobiolrev.2011mar；35(2):343-59.doi:10.1111/j.1574-6976.2010.00251.x.epub2010oct7.review.pubmedpmid:21039646。
当然，上述说明并非是对本发明的限制，本发明也并不限于上述举例，本技术领域的普通技术人员，在本发明的实施范围内，做出的变化、改型、添加或替换，都应属于本发明的保护范围。
1.一种单牙位微生物群落检测模型的构建方法，其特征在于，包括以下步骤：
数据收集：在待收集人群中收集每颗牙齿各自的口腔微生物群落数据；
模型的初步构建：以口腔微生物群落数据为样本，通过机器学习算法构建模型；
模型的最终构建：根据初步模型获得每个微生物因子对模型性能的贡献程度，将其进行降低排列，并建立基于不同个数微生物因子组合的所有模型，进而筛选出对于模型性能贡献程度最大的一组微生物因子，作为最终的单牙位微生物群落检测模型。
2.根据权利要求1所述的口腔微生物群落检测模型的构建方法，其特征在于，所述数据收集为在待收集人群中区分出患有牙病人群和未患有牙病人群中每颗独立牙齿在的口腔微生物群落数据。
3.根据权利要求1所述的口腔微生物群落检测模型的构建方法，其特征在于，所述口腔微生物群落数据为高通量测序手段获得微生物群落的16srna或全基因组信息。
4.根据权利要求1所述的口腔微生物群落检测模型的构建方法，其特征在于，所述模型的初步构建，包括以下步骤：
(1)以口腔微生物群落数据作为输入；
(2)样本的状态作为输出，包括微生物标记物的宿主状态，某种状态的可能性；
(3)采用随机森林方法得到初步检测模型；
(4)获取每个微生物因子对于模型性能的贡献程度。
5.根据权利要求1所述的单牙位口腔微生物群落检测模型的构建方法，其所述模型的最终构建，包括以下步骤：
(1)将初步模型构建后获得的各个微生物因子对模型性能的贡献程度进行降序排列；
(2)根据降序排列结果，取前n个贡献程度值对应的微生物因子进行组合，建立基于不同微生物因子组合的所有模型；n为自然数且小于微生物因子个数；
(3)筛选出对模型贡献程度值低于阈值并且准确率低于设定值的微生物因子作为无关微生物因子；
(4)从微生物组数据中剔除无关生物因子；并删除无关微生物因子对应的模型；
(5)建立以剩余贡献程度的微生物因子为输入的单牙位微生物群落检测模型。
6.根据权利要求4所述的单牙位口腔微生物群落检测模型的构建方法，其特征在于，样本类别表示微生物标记的宿主的状态，包括宿主此时的正常、此时的异常、未来的正常、未来的异常。
7.根据权利要求1所述的单牙位口腔微生物群落检测模型的构建方法，其特征在于，所属微生物因子为微生物物种信息。
8.根据权利要求1所述的口腔微生物群落检测模型的构建方法，其所述模型性能是指检测模型区分样本状态的准确率。
9.根据权利要求1所述的口腔微生物群落检测模型的构建方法，其所述贡献程度，是指某一种微生物因子对于区分样本状态的能力，体现其对于模型准确率的贡献程度。
10.根据前述权利要求任一项的口腔微生物群落检测模型的应用，其特征在于，所述模型用作检测口腔疾病发生的可能性评价。
技术总结