异源的人类杂质染色体之间的重组

　　我们分析了HPRC4产生的47 T2T分阶段二倍体DE从头组合（总共94个单倍型）。我们包括T2T-CHM13版本23和GRCH38。 　　我们首先使用全部映射来构建HPRCY1 pangenome中同源关系的无参考模型。这将完整的HPRCY1建模为映射图，其中节点是重叠群，边缘表示它们之间的映射。为了构建HPRCY1映射图，我们基于95％平均核苷酸同一性的50-kb种子的链生成同源映射（我们希望支持同源重组17），以至于（n-1）= 93个替代映射。我们首先将此过程应用于所有其他所有其他38,325个HPRCY1重叠群，以获取38,036个映射，其中覆盖了总组装序列的99.9％。这表明，在38,325（99.2％）的HPRCY1组件重叠群中，有38,036个与至少1个其他重叠群同源。用于构建HPRCY1的组装图中的复杂缠结会产生简短的重叠群，并倾向于导致较高的错误率4。因此，为了简化以后的分析并专注于整理良好的组件区域，我们缩小了重点，以考虑仅考虑涉及至少1 MB长的16118个重叠群的映射，覆盖了总组装序列的98.72％。 　　然后，我们构建了一个图形，节点是重叠群，边缘表示它们之间的映射 - “映射图”。此映射图中的边缘的权重等于估计的序列身份乘以映射的长度。为了推断每个重叠群代表的染色体，我们对T2T-CHM13和GRCH38参考的所有重叠群绘制了所有重叠群，并根据此映射为它们分配了染色体身份。该映射图非常密集，每个重叠群最多93次映射，因此很难直接使用现有方法可视化。为了在此图中开发有关模式的直觉，我们相反查看了一个从每个重叠群段的最佳三个映射构建的简化映射图，并用其参考依赖性分配标记每个重叠群（图1A）。杂技集群（图1b）通常与我们先前对Chr之间更高相似性的期望相匹配。13和chr。21，介于chr之间。14和chr。2218,19。 　　为了量化这些模式的重要性，我们将社区检测算法56应用于完整的映射图。该算法将每个重叠群分配给一个社区，以使总分配最大化模块化，这可以理解为与社区之间的链接相比，社区内部（加权）链接的密度。该过程产生了31个社区（补充文件1）。我们假设每个簇代表一个染色体或染色体臂。大约一半的染色体（n = 11）由一个社区代表。染色体1、2、3、6和18分别在两个社区中代表与它们的短和长臂相对应的社区，这可能是由于频繁的组装在其centromeres上频繁断裂（图1C，D）。染色体X和Y的重叠群在同一社区中落下，尽管X染色体的短臂在两个社区中有代表（图1D）。SAAC构成了社区最不同的染色体和大多数重叠群（1,706个包含3.91％的HPRCY1序列的重叠群），由属于所有五个肢体中心染色体的短臂的重叠群组成。21q and Chr。22Q（图1C，D）。chr。13q，chr。14q and Chr。15q每个都有自己的社区。与染色体13、14和15相比，由P臂重叠群组成的社区中染色体21和22的Q臂可能与它们的短长度有关。当我们提高映射的灵敏度时，我们获得了相似的结果（补充图29）。 　　在HPRCY1的同源映射图中，只有杂技体和性别染色体形成了包含多个染色体的组合群落。性染色体社区反映了X和Y57上的PARS，这是端粒区域，在这些区域中，这些原本非同源的染色体可以重新组合，就好像它们是同源物一样。我们假设以阶级为中心的社区也可能反映了正在进行的伪同源重组 　　我们使用WFMASH序列Aligner58（提交AD8AEBA）对47个T2T分阶段二倍体de Novo组件进行了所有重叠群的成对映射。我们设置以下参数： 　　我们使用了50 kb（）的段种子长度，需要至少约250 kb（）和估计的核苷酸认同，至少为〜95％（）。总共有94张单倍型，我们为每个重叠群（）保留了多达93个映射。此外，当查询和目标在“＃”字符之前具有相同的前缀时，我们跳过了映射，那就是涉及相同的单倍型。为了正确地通过重复区域绘制，只有0.001％的最常见的公共忽略了（）。我们跳过了基础对齐（）。我们还以相同的参数生成了成对映射，但是使用10 kb的段种子长度，并且需要至少约50 kb的同源区域。 　　从由此产生的映射中，我们排除了涉及重叠群比1 MB短的映射，以减少虚假比赛的可能性。然后，我们使用Python脚本（在PGGB存储库中传递）来构建结果的图表表示（映射图），分别代表它们之间的节点和边缘，分别代表它们之间的重叠群和映射。 　　该脚本会产生一个代表边缘的文件，代表边缘权重的文件以及将图形节点映射到序列名称的文件。边缘的重量由长度的乘积和相应映射的核苷酸身份（较高的权重与较高的映射相关）。最后，我们使用Python脚本（在PGGB存储库中交付）应用了IGRaph Tools59中实现的Leiden AlgorithM56，以检测映射图中的基础社区。 　　为了确定每个社区中代表哪些染色体，我们通过将它们映射到具有WFMASH的T2T-CHM13V1.1和GRCH38人类参考基因组中，对所有重叠群进行了映射，这一次需要同源区域至少150 kb的长度和至少90％的核苷标识。 　　我们在映射过程中禁用了重叠群拆分（），以获得覆盖整个重叠群的同源区域。对于未绘制的重叠群，我们用相同的参数重复了映射，但允许重叠群拆分（无需指定）。我们将重叠群的“ P”或“ Q”标记，具体取决于它们覆盖其属于的短臂还是长臂。完全跨越中心粒的重叠群标记为“ PQ”。我们使用此类标签来识别未经参考序列获得的映射图中检测到的社区的染色体组成，并注释映射图中的节点。 　　为了获得对HPRC组件之间同源关系关系的清晰可视化，我们通过使用用于主图的相同参数生成了一个更简单的映射图，但是为每个重叠群保留多达3个映射，并添加T2T-CHM13参考基因组版本2，其中还包括完整的HG002 Chromosome Y（https://www.ncbi.nlm.nih.gov/assembly/gcf_009914755.1）： 　　我们将窗口大小设置为素描等于5000（），以减少较少的kmers来减少运行时。我们使用Python脚本构建映射图，然后使用Gephi60（0.9.4版）可视化它。我们通过运行“随机布局”，然后是“ Yifan Hu”算法来计算映射图布局。 　　对于47个T2T分阶段的二倍体二倍体组件中的每个组件，我们用WFMASH序列对准器（commit ad8aeaeba）绘制了所有重叠群针对T2T-CHM13人参考基因组。对于HG002样本，我们包括了两个组件：HG002-HPRCY1从Novo组装组合（使用HIFI读取）和一个基于HIFI和ONT读取的NOVO DE NOVO组装，并使用Verkko Assember构建。我们设置以下参数： 　　我们使用了50 kb（）的段种子长度，需要至少约150 kb（）和估计的核苷酸认同，至少为〜90％（）。我们仅保留每个重叠群（）的映射（最好的映射）。为了正确地通过重复区域绘制，只有0.001％的最常见的公共忽略了（）。我们跳过了基础对齐（）。对于HG002-HPRCY1重叠群，我们禁用重叠群拆分（）。 　　然后，我们确定了源自杂技染色体的重叠群，并涵盖了它们所属的染色体的短和长臂。我们仅在手臂上至少有1 kb的映射重叠群，距离丝粒至少1 MB。我们称此类重叠群为“ P – Q Acrentric Contigs”。对于HG002-HPRCY1，无论覆盖归属染色体的两个臂，都考虑了仅比重叠群更长或等于300 kb。 　　最后，我们通过应用PGGB22（提交A4A6668）构建了所有P – Q Acrentic重叠群以及T2T-CHM13和GRCH38人参考基因组的Pangenome图。我们设置以下参数： 　　我们使用了50 kb（）的段种子长度，需要至少约250 kb（）和估计的核苷酸认同，至少为〜98％（）。在输入中具有142个P – Q族裔重叠（来自HG002-HPRCY1和HG002-VERKKO的10个）以及10个来自T2T-CHM13和GRCH38参考基因组T2T-CHM13和GRCH38参考基因组加49 HG002-HPRCY1-HPRCY1对其他10个Acrocentric and arpercentric Haplotypes（5个杂志）的10个杂种（5个杂志）（5个杂货）（5个杂志）的10个杂种（5个杂志）（5个杂志）（5个杂货）（5 pere）（5个杂志）（5个杂志）（5个杂志）（5 per）（5 per）（142+10+10）对于每个重叠群（）。为了正确地通过重复区域绘制，只有0.001％的最常见的公共忽略了（）。我们过滤了比对匹配的匹配短于311 bp，以消除由短反复的同源物引起的可能的伪造关系（）。我们为两轮的图归一化（和）设置了大目标序列长度和一个小序列填充。为了可视化中心的pangenome图，我们使用ODGI Layout26（提交E2DE6CB）构建了图形布局，并使用gfaestus61（提交50FE37A）进行了可视化。这将小变体的序列和链条作为线性结构，而由分段重复，反转和其他结构变体引起的重复趋向于形成循环。 　　我们通过应用ODGI Untangle（提交E2DE6CB）来解开pangenome图。实际上，我们将图表投影到一组查询（HPRCY1重叠群）和参考（T2T-CHM13）序列之间的对齐中。我们设置以下参数： 　　我们将图表分为〜50 kb（）的常规大小区域，合并区域短于1 kb（）。我们报告了每个查询段（）的第100个最佳目标映射（），而没有为Jaccard相似性应用任何阈值（）。我们将图中的所有路径用作查询，并将其投影针对T2T-CHM13基因组（）的五个中心染色体。此外，我们发出用于分割图形（）的切割点。 　　对于每个查询段，如果针对不同的目标有多次命中（即具有相同，最高的Jaccard相似性的命中），我们将其作为查询原始染色体的第一个命中（从重叠群的染色体分配获得）。 　　我们重复了五次缠绕的图形，但限制了该算法仅使用T2T-CHM13的杂技染色体之一作为目标（）（）并返回最佳匹配的命中（）。 　　当使用T2T-CHM13的所有杂技染色体作为靶标（）时，我们使用了生成的切口。这样，所有无缠结的运行（总共六个）都使用了相同的切割点，用于段边界。 　　最后，我们将所有中心染色体产生的无污点输出“接地”作为目标：更详细地，将每个无区域性的查询段用于针对特定的中心型染色体（不仅是匹配最佳的匹配）（不仅是匹配最佳的染色体），该片段使用限制于单个目标的无骨输出。我们将结果分割为杂技染色体，仅保留靶标的P臂和Q臂的查询。此外，我们删除了在集会中标记为不可靠（即存在覆盖范围问题）的查询段重叠区域。Flagger是一种基于HIFI读取的管道，它通过识别跨母体和父亲单倍型的读取映射覆盖范围的不一致来检测分阶段二倍体组件中不同类型的错误组件。为了关注更相似的查询目标命中，我们使用Jaccard指标来通过应用参考文献中报道的校正公式来估计序列身份。62，仅保留以大于90％的估计身份的结果。为了分析以SST1数组为中心的分段重复中的HPRCY1重叠群的方向状态，我们生成了一个带有ODGI FLIP的新的Pangenome图（COMPL 0B21B35）。更详细地说，如果相对于pangenome图，我们首先将路径绕过反向补体方向。这导致对HPRCY1重叠群具有均匀的方向，这都是相对于图形的正向方向。然后，我们以与上述相同的方式解开了翻转图。我们显示了配备以下软件包的R开发环境（版本3.6.3）的每个杂技染色体的无缠结结果：Tidyverse（版本1.3.0），Rcolorbrewer（版本1.1.2），GGPLOT2（版本1.1.2），GGPLOT2（版本3.3.3）和GGREPEL（版本0.9.1.1）。 　　对于锚定在T2T-CHM13 Actentric Q-ARM上的每组HPRCY1-ACRO重叠群，我们计算了具有最佳打击的重叠群的数量。特别是，对于每组T2T-CHM13肢体染色体的每个基础位置，我们量化了每个杂技的每一个杂技表中的每一个杂物，在pangenome图中显示出最佳打击。我们仅认为最佳命中率至少为90％。 　　为了量化无障碍结果中的混乱程度，我们计算了以最佳攻击为主的不同杂技院之间的多样性熵。更详细地，我们将每个HPRCY1 Acrentric P – Q重叠群投影，以通过Q-ARM锚定的T2T-CHM13 Accentric，并将每个参考基础位置与相应的中心最佳命中相关联。我们仅认为最佳命中率至少为90％。然后，我们计算了50 kb的窗口中的香农多样性指数（SDI），长距离长方体。我们在重叠群不匹配任何目标的区域中使用-1作为缺少的SDI值。对于每组重叠群，我们通过计算每个参考基础位置的平均值（忽略丢失的SDI值）来汇总SDI结果。我们将此指标称为位置同源性熵，它显示出可以将重叠群描述为不同参考染色体的镶嵌物的区域。但是，它无法区分有不同的参考染色体相似性（可能指示重组交换）的区域，从每个重叠群与T2T-CHM13之间存在区域多样性的地方。如果T2T-CHM13本身在杂技杂志之间包含罕见的重组，或者古代同源性可能会导致毫无障碍的对齐方式“噪声”，则可能会发生后一种情况，因为重叠群从两个同样好的替代映射中挑选。为了避免这些陷阱并建立更严格的图形重组度量，然后我们扩展了毫无障碍的多样性指标，以在多个映射上运行。 　　要考虑到除第一个命令之外的其他命中之外，我们将多样性熵指标概括为在前5个毫无障碍的命中的命令中工作，并共同考虑所有重叠群。对于每个参考段，我们为HPRCY1 Acrentric P – Q重叠群中的每一个收集了相应的最佳5个无缠绕命中。这是可能的，因为参考段在所有重叠群中都是稳定的。我们仅认为最佳命中率至少为90％。为了避免通过以PVG结构建模的分段重复（如Chr。13q，Chr。13Q，Chr。15Q和Chr。22Q;图2;图2）引起的分段重复引起的染色体内相似性，以驱动脱离熵的熵，我们忽略了连续的重复目标命中，我们忽略了连续的重复目标命中 - 其他单词，我们采取了独特的参考目标的顺序。我们采用了独特的参考目标。当多个重叠群段以相同的参考段为基础时，我们认为第一个重叠群段具有最佳的接地，即在针对当前参考段时具有最高的估计身份。然后，我们根据估计的相似性对五个最佳命中进行了排名。最后，对于每个参考段，我们在所有可用的五个最佳命中订单中计算了SDI。我们将-1用作参考区域中缺少SDI值，而没有任何重叠匹配。我们还将输出保存在输出中，以了解每个参考段中有多少hprcy1 accentric p – q重叠群有助于熵计算。这产生了位置同源性熵。 　　为了获得PHR，我们通过考虑具有大于0的位置同源性熵的区域来汇总最终结果，并由至少1个重叠群支撑，与BedTools63合并的区域距离小于30 kb，并删除比30 kb的合并区域。 　　我们显示了每个中心染色体的汇总结果。我们使用T2T-CHM13V2.0 UCSC TrackHub（https://genome.ucsc.edu/cgi-bin/cgi-bin/hgtracks?db = hub_3671779_hs1），使用T2T-CHM13V2.0 UCSC TrackHub（使用T2T-CHM13V2.0）使用基因组注释。我们在https://github.com/pangenome/chromosome_communities/tree/main/main/scripts上使用了脚本的数字。为了绘制图形，我们使用了配备以下软件包的R（版本3.6.3）：Tidyverse（版本1.3.0），Rcolorbrewer（版本1.1.2），GGPLOT2（版本3.3.3）和GGREPEL（版本0.9.1）。最后，我们使用Inkscape（https://inkscape.org/）根据结果来构成主要文本数字，并提供这些方法直接产生的补充图。 　　我们在使用PGGB4与GRCH38和T2T-CHM13构建的完整HPRCY1 Pangenome中显示了双重SNP密度。为此，我们从已发布的VCF文件中与每个染色体提取了双重SNP，均为参考（GET_BISNP.SH）。由于在HPRCY1中使用的T2T-CHM13版本1.1没有Y染色体，因此我们使用了GRCH38（包括蒙版PAR1和PAR2区域）的T2T-CHM13版本。我们使用R（版本4.1.1）和Tidyverse（版本1.3.1）软件包（plot_bisnp_dens.r）显示了100 kb的垃圾箱中的双重SNP密度。 　　我们从参考文献中绘制了BAC克隆。8针对具有WFMASH序列对准器的T2T-CHM13人参考基因组（COMM aD8aeba）。我们只保留覆盖杂技染色体的映射，估计的身份至少为90％。为了绘制图形，我们使用了配备以下软件包的R（版本3.6.3）：Tidyverse（版本1.3.0），Rcolorbrewer（版本1.1.2），GGPLOT2（版本3.3.3）和GGREPEL（版本0.9.1）。根据参考文献，我们为BAC克隆的映射着色。8。 　　我们使用基于最佳拟合替代模型（Kimura 2-参数 +G，参数= 5.5047）的最大似然方法进行了系统发育分析，该模型由Jmodeltest264推断出具有1,000个自举重复的Jmodeltest264。在每个节点的底部指示高于75的引导值。 　　我们从参考文献中获得了人类PRDM9结合基序（总共17个）。41并使用FIMO65扫描其在T2T-CHM13V2.0人类参考基因组中的发生： 　　FIMO相对于每个序列位置计算每个基序的对数 - 基础比率评分，并使用动态编程（假设零阶空模型将这些得分转换为p值，在该模型中，在该模型中以用户指定的每个字体背景频率随机生成序列，然后估算错误发现率65。每个基序都与它代表PRDM9的真正结合目标的可能性相关联。我们保留了下游分析，仅该量度至少为70％（17个中的14个）的基序。对于每个图案，我们使用BedTools63来计算每个T2T-CHM13V2.0染色体中窗口20 kb长的发生次数的数量。 　　为了绘制图形，我们使用了配备以下软件包的R（版本3.6.3）：Tidyverse（版本1.3.0），Rcolorbrewer（版本1.1.2），GGPLOT2（版本3.3.3）和GGREPEL（版本0.9.1）。 　　我们使用VG DeConstruct5： 　　我们称关于T2T-CHM13参考基因组（）的变体，每种HPRCY1杂技p – q重叠群的报告变体（and）。我们仅考虑与图（）中的路径（即重叠群）相对应的遍历（即），还报告了嵌套变体（）。从变体集合中，我们仅考虑单核苷酸变体。我们通过使用PLINK V1.966在单倍体集合时使用PLINK V1.966估算了在70 Kb内的链接不平衡，并保留R2> 0的所有值（plot_ld_1.r）。最后，我们使用R（3.6.3版）的置信区间生成了BINNED连锁不平衡的衰减图，重点是相距不到4 kb的对。 　　我们应用了Verkko（Beta 1，Commit VD3F0B941B5FACF5807C303B0C0171202D83B7C74）的较早版本，以使用HIFI（105X）和ONT（85X）和ONT（85X）和DACEDS reads readS构建HG002单元线的二倍体组装图。所得组装图将每个单倍型的单个重叠群中的近端连接解析为多型碱基碱基，而远端连接仍待解决。我们使用了来自父母Illumina读取的均聚物压缩标记，以将单位物分配给母亲，父亲单倍型或模棱两可的标记，如果不足的标记不支持任何一种单倍型。For estimating the number of times a unitig has to be visited, we aligned HiFi and ONT reads to the assembly graph using GraphAligner with the following parameters: –seeds-mxm-length 30–seeds-mem-count 10000 -b 15–multimap-score-fraction 0.99–precise-clipping 0.85–min-alignment-score 5000–clip-ambiguous-ends100 –超局部不兼容切割0.15 –玛克斯 - 跟踪计数5 – hpc-collapse-reads-discard-cigar68。四个远端连接连接到rDNA阵列，并带有模棱两可的节点连接母体和父亲节点，并支撑它们属于同一染色体。两个远端连接单位，一个母体和父亲，彼此断开，但连接到rDNA阵列，这些阵列分配给了同一染色体。使用标记和ONT对齐，我们确定了图中的路径，并根据最支持的单倍型分配了路径。如果单倍型分配的单元之间仅存在模棱两可的节点，而没有ont读取路径，则将节点随机分配到一个单倍型以构建重叠群。确定所有路径后，我们使用了共识。现在，在最新的Verkko（V1）中，以更具自动化的方式获得了整个生产父母标记，根据其单倍型标记其单倍型的标记。 　　为了提供HIFI重叠组件的交叉验证及其对它们的分析，我们比较了同一样品的两个组件（HG002）的毫无用处。一个是用HPRCY1管道制成的，而另一个使用Verkko二倍体T2T汇编器。Verkko在基于HIFI的组装图中采用了ONT到无障碍区域，该区域最初是在T2T-CHM13生产中开发的。Verkko的组装从ONT汇总了信息，因此提供了使用替代测序和组装方法进行分析的单个集成目标。 　　我们通过比较了两个HG002组件的最佳命中来验证了pangenome毫无障碍的结果，该组件是使用HIFI读取的一个，另一个基于HIFI和ONT读取，并使用Verkko Assembler构建，并基于Verkko Assembler67。对于每个T2T-CHM13杂技染色体的每个基本位置，我们比较了HG002-HPRCY1重叠群的无缠绕最佳打击与由HG002-Verkko重叠支持的最佳打击。我们仅认为最佳命中率至少为90％。当两个HG002组件都支持与最佳命中相同的T2T-CHM13 Accentric时，我们将参考区域定义为一致。我们分别处理了两种单倍型（母亲和父母）。 　　我们观察到两种方法在染色体同源性镶嵌图的水平上的高度一致性。在HG002-Verkko组装中，最佳打击的无缠结显示与HG002-HPRCY1中的模式相似（补充图9-13）。但是，在HG002-HPRCY1中，一些SAAC单倍型似乎组装不力。在Q臂上，我们在HG002-HPRCY1和HG002-VERKKO毫无疑问的结果之间达到了99.93％的一致性，但P臂上只有87.45％的一致性（补充文件8）。由于其多个重复的序列（包括PHR），卫星阵列和rDNA，因此该较低的水平与组装SAAC的难度更大。我们发现这种不和谐是由单个染色体单倍型驱动的：虽然大多数P臂达到了90％的一致性，但HG002-HPRCY1 14P-MARTNAL在组件中表现出很高的不一致性（66.19％的一致性）（补充图19）。尽管此n = 1验证仅着眼于10个单倍型，但它包含了HG002-Verkko中使用的深HIFI（105倍）和ONT（85X）数据提供的许多独立读取。因此，我们比较了在SAACS读取单分子读取中观察到的结构与代表产生我们pangenome的仅HIFI纯组装过程的HG002-HPRCY1组件的结构之间的一致性。 　　但是，该分析应被视为在我们的过程准确性上呈现下限。我们正在考虑所有映射到accentrics的HG002-HPRCY1重叠群，而不仅仅是那些可以满足我们的丝粒划线要求的HG002-HPRCY1本身比我们选择的其他为acro-PVG4选择的组件更为分散。其重叠群的碎片性质（只有一个来自Chr。22满足我们的P – Q映射要求）可能会引入与HG002-Verkko的其他分歧。总体结果表明，在仅HIFI纯组件中观察到的大多数模式可能受到同一样品的近-T2T组件的自动化支持。 　　有关研究设计的更多信息可在与本文有关的自然投资组合报告摘要中获得。