# 整个流程可以分为以下几个主要步骤：
   1. PacBio HiFi 测序数据回比对：将长读长（PacBio）测序数据比对回初步组装的基因组上，用于计算测序深度。
   2. 计算测序深度和覆盖度：使用 pbcstat 工具分析比对结果，生成覆盖度统计信息。
   3. 确定深度阈值：使用 calcuts 工具根据深度分布，自动计算用于区分正常序列、重复序列和垃圾序列的深度阈值。
   4. 基因组自身比对：将组装好的基因组序列自己和自己比对，以识别结构上相似的区域（即潜在的重复区）。
   5. 识别并筛选重复序列：结合深度阈值和自身比对结果，使用 purge_dups 主程序识别出哪些是冗余序列。
   6. 提取最终序列：根据上一步的识别结果，从原始组装文件中移除冗余序列，保留主序列，生成最终的纯净版基因组。

1. PacBio HiFi 数据回比对



   1 minimap2 -t 60 -x map-pb final_asm.bp.p_ctg.fa out.fasta.gz | gzip -c - > pb_aln.paf.gz



   * 命令: minimap2
   * 目的: 将 PacBio HiFi 测序数据（out.fasta.gz）比对到初步组装好的基因组（final_asm.bp.p_ctg.fa）上。这是为了后续计算每个contig（拼接序列）的测序深度。
   * 参数:
       * -t 60: 使用 60 个线程进行计算，以加快速度。
       * -x map-pb: 使用 minimap2 预设的参数集，专门用于将 PacBio 测序数据比对到参考基因组。
   * 输入文件:
       * final_asm.bp.p_ctg.fa: 参考序列，即初步组装好的基因组文件。
       * out.fasta.gz: 查询序列，即 PacBio HiFi 测序的 reads 文件（FASTQ 或 FASTA 格式，此处为 .gz 压缩格式）。
   * 输出文件:
       * pb_aln.paf.gz: 比对结果文件，格式为 PAF (Pairwise mApping Format)。通过管道（|）和 gzip 命令进行压缩，节省存储空间。

  2. 计算覆盖度统计



   1 ~/Software/purge_dups/bin/pbcstat pb_aln.paf.gz



   * 命令: pbcstat
   * 目的: 分析上一步生成的 PAF 比对文件，计算基因组中每个位置的测序深度和覆盖度统计信息。
   * 输入文件:
       * pb_aln.paf.gz: minimap2 生成的比对文件。
   * 输出文件:
       * PB.stat: 默认输出文件，包含了每个 contig 的覆盖度统计信息。
       * PB.base.cov: 默认输出文件，包含了基因组每个碱基的覆盖深度。

  3. 计算深度阈值



   1 ~/Software/purge_dups/bin/calcuts PB.stat > cutoffs 2> calcults.log



   * 命令: calcuts
   * 目的: 基于 pbcstat
     生成的统计文件，自动计算和设定测序深度的分界点（cutoffs）。它会根据深度分布图，估算出正常单倍型序列的平均深度，并以此为依据设定一个上限和下限，用于区分正常序列和重复序列（通常深度更高）。
   * 输入文件:
       * PB.stat: pbcstat 生成的统计文件。
   * 输出文件:
       * cutoffs: 包含了计算出的深度阈值。这个文件会被后续的 purge_dups 命令使用。
       * calcults.log: 记录 calcuts 运行过程中的日志信息（通过 2> 重定向标准错误流）。

  4. 基因组自身比对



   1 ~/Software/purge_dups/bin/split_fa final_asm.bp.p_ctg.fa > asm.split
   2 minimap2 -t 60 -xasm5 -DP asm.split asm.split | pigz -c > asm.split.self.paf.gz



   * 第一部分: split_fa
       * 目的: 将原始的基因组 FASTA 文件进行切分或格式化，使其适合后续的自身比对。
       * 输入: final_asm.bp.p_ctg.fa
       * 输出: asm.split (格式化后的 FASTA 文件)
   * 第二部分: minimap2
       * 目的: 将基因组与自身进行比对，以找出所有相似度高的区域。这些区域是潜在的重复序列或单倍型序列。
       * 参数:
           * -t 60: 使用 60 个线程。
           * -x asm5: 使用预设参数集，适用于高质量基因组（错误率 <5%）之间的比对。
           * -DP: -D 标记主链（primary chain），-P 保留次要链（secondary chains），有助于识别所有可能的重叠。
       * 输入文件:
           * asm.split: 两次作为输入，表示将这个文件与它自身进行比对。
       * 输出文件:
           * asm.split.self.paf.gz: 自身比对的结果，以 PAF 格式存储，并用 pigz（一个并行的 gzip 工具）进行压缩。


  5. 识别冗余序列



   1 ~/Software/purge_dups/bin/purge_dups -2 -T cutoffs -c PB.base.cov asm.split.self.paf.gz > dups.bed 2> purge_dups.log



   * 命令: purge_dups
   * 目的: 这是核心步骤。它整合了测序深度信息和自身比对结果，来最终判断哪些序列是冗余的。
   * 参数:
       * -2: 表示使用 purge_dups 的第二版算法（通常是默认或推荐的）。
       * -T cutoffs: 指定包含深度阈值的文件，即 calcuts 命令生成的 cutoffs 文件。
       * -c PB.base.cov: 指定每个碱基的覆盖度文件，即 pbcstat 生成的 PB.base.cov 文件。
   * 输入文件:
       * asm.split.self.paf.gz: 基因组自身比对的结果。
       * cutoffs: 深度阈值文件。
       * PB.base.cov: 碱基覆盖度文件。
   * 输出文件:
       * dups.bed: 一个 BED 格式的文件，列出了所有被识别为冗余（duplication）的序列的坐标。
       * purge_dups.log: 记录 purge_dups 运行过程的日志。

  6. 提取最终的非冗余序列


   * 命令: purge_dups
   * 目的: 这是核心步骤。它整合了测序深度信息和自身比对结果，来最终判断哪些序列是冗余的。
   * 参数:

   1 ~/Software/purge_dups/bin/get_seqs dups.bed final_asm.bp.p_ctg.fa



   * 命令: get_seqs
   * 目的: 根据 purge_dups 生成的 dups.bed 文件，从原始的组装基因组中移除这些被标记为冗余的序列，并输出一个纯净版的基因组。
   * 输入文件:
       * dups.bed: 包含冗余序列坐标的 BED 文件。
       * final_asm.bp.p_ctg.fa: 原始的、未经过滤的组装基因组。
   * 输出文件:
       * 默认情况下，它会生成两个文件：
           * purged.fa: 最终的、去除了冗余序列的基因组文件。这是整个流程最重要的产出文件。
           * hap.fa: 被识别为冗余或单倍型（haplotigs）的序列文件。