Previous short-read-based genome assembly had a limitation in that short read lengths (tens to hundreds of bp) made alignment and assembly difficult in repetitive sequences or complex genomic regions. As a result, genomes were often assembled in a fragmented manner, or important regions were missed. In contrast, long-read sequencing technologies (such as PacBio and Oxford Nanopore) can generate long reads extending from several kb to tens of kb, allowing for the alignment of entire repetitive sections (Figure 1). This reduces gaps in the assembly process and enables more accurate reconstruction of the genome structure.

A representative case is the project that completely sequenced the human genome from telomere-to-telomere (T2T), which clearly demonstrates the advantages of long-read technology. Furthermore, long-read-based assembly is being used as a core tool in advanced genome research for major crops and livestock.

To obtain genetic variants, we used various bioinformatics software and followed a sequence of steps: quality control and trimming > initial assembly > final assembly (Figure 2). Additionally, using the completed assembly, we constructed a pangenome following the method of minigraph-cactus (github).

1. Long-Read Data Quality Control (Read QC)

   A representative case is the project that completely sequenced the human genome from telomere-to-telomere (T2T), which clearly demonstrates the advantages of long-read technology. Furthermore, long-read-based assembly is being used as a core tool in advanced genome research for major crops and livestock.

   bash Code Copy

   ccs input.bam output.ccs.bam 
   porechop -i KPR_ont_reads.fastq -o KPR_ont_trimmed.fastq

2. Adapter filtering/trim

   For PacBio data, residual adapter sequences were removed using HiFiAdapter. For ONT data, adapter trimming was already applied during the preceding QC step.

   bash Code Copy

   hifiadapterfilt —in KPR_input.ccs.fastq —out KPR_filtered.fastq

3. Initial assembly

   An initial assembly was performed using platform-specific assemblers. For PacBio HiFi reads, assembly was conducted using HiFiasm. For ONT PromethION reads, which have a higher error rate, assemblers such as Flye, Canu, Shasta were applied.

   bash Code Copy

   hifiasm \
       -o KPR.asm \
       --dual-scaf \
       --telo-m CCCTAA \
       --h1 KPR_read1.fq.gz \
       --h2 KPR_read2.fq.gz \
       --ul KPR_ul.fq.gz \
       KPR-HiFi.fq.gz
   
   flye \
       --nano-raw KPR_ont_trimmed.fastq \
       --out-dir KPR_flye_out \
       --genome-size 3g \
       --threads 32
   
   canu \
       -p asm \
       -d KPR_canu_out \
       genomeSize=3g \
       -nanopore-raw KPR_ont_trimmed.fastq \
       useGrid=false \
       maxThreads=32
   
   shasta \
       --input KPR_ont_trimmed.fastq \
       --assemblyDirectory shasta_out \
       --threads 32

4. Long-read polishing

   The initial assembly result underwent a polishing process for error correction. For PacBio, residual errors were corrected using Inspector. For ONT reads, a multi-step polishing was performed by combining tools like Racon, Nanopolish, and Inspector.

   bash Code Copy

   racon KPR_ont_reads.fastq KPR_ont_aln.sam flye_out/KPR_assembly.fasta > racon_KPR_polished.fasta
   nanopolish variants --consensus polished.fa \
       --reads KPR_ont_trimmed.fastq \
       --bam KPR_aln.sorted.bam \
       --genome flye_out/assembly.fasta \
       --threads 32
   inspector -i KPR_output.asm.p_ctg.fa -b KPR_aln.bam -o KPR_inspector_out

5. Contaminated & mis-join contig filtering
   1. Contamination detection

      Assembly results often contain external contaminant sequences or mis-joined contigs. Since these can cause fatal errors in downstream analysis, a filtering process was performed.

      bash Code Copy

      blobtools create --fasta KPR_assembly.fasta --bam KPR_aln.sorted.bam --taxid-mapping taxdump
      kraken2 --db kraken_db --threads 32 --output KPR_kraken.out --report KPR_kraken.report assembly.fasta

   2. Mis-join contig detection

      We verified the quality of the genome assembly using several tools. We used QUAST to quantitatively evaluate the assembly quality with metrics like N50 and mis-assembly rate by comparing it to a reference genome. Additionally, we used Inspector or Tigmint to align long-read or linked-read data and find regions with abnormal read support. Finally, we utilized Hi-C data along with 3D-DNA or Juicer to identify and correct chromosome-scale mis-joins.

      bash Code Copy

      quast.py KPR_assembly.fasta -r reference.fa -o quast_out
      tigmint-make tigmint draft=assembly reads=KPR_ont_trimmed.fastq
      juicer.sh -d juicer_out -D /path/to/juicer_tools -z KPR_assembly.fasta -p chrom.sizes -s MboI
      run-asm-pipeline.sh KPR_assembly.fasta juicer_out/aligned/KPR_merged_nodups.txt

6. Construction of the Pangenome
   1. Preparing the Input Seqfile

      First, prepare the assembled fasta files from the previous steps. Create a Seqfile by entering the sample's structure as follows (using KPR's KOR-101 sample as an example).

      bash Code Copy

      KOR-101.1   paternal.fa
      KOR-101.2   maternal.fa
      CHM13.0   chm13.fa
      GRCh38    grch38.fa

   2. Fixing problem contigs

      Use commands such as the following to ensure all contigs are in the same format.

      bash Code Copy

      wget -q $(grep KOR-102\.1 KPR-mc.seqfile | tail -1 | awk '{print $2}') -O KOR-102.1.fa.gz
      gzip -d KOR-102.1.fa.gz
      samtools faidx KOR-102.1.fa

   3. Adding the Y-chromosome to the CHM13 reference used together

      Since CHM13 lacks chrY, we add chrY. Next, we also update the Seqfile.

      bash Code Copy

      wget -q $(grep CHM13 KPR-mc.seqfile | tail -1 | awk '{print $2}') -O CHM13.Y.v1.1.fa.gz
      wget -q $(grep GRCh38 KPR-mc.seqfile | tail -1 | awk '{print $2}') -O GRCh38.fa.gz
      gzip -d GRCh38.fa.gz
      
      samtools faidx GRCh38.fa chrY chrEBV | bgzip --threads 8 >> CHM13.Y.KPR.fa.gz
      aws s3 cp CHM13.Y.KPR.fa.gz ${MYBUCKET}/fasta/
      
      grep -v CHM13 KPR-mc.seqfile > t
      printf "CHM13_1_1\t${MYBUCKET}/fasta/CHM13.Y.KPR.fa.gz\n" > KPR-mc.seqfile
      cat t >> KPR-mc.seqfile

   4. Standardizing Contig Naming Rules

      Using cactus-preprocess, we corrected the contigs to the SAMPLE.HAPLOTYPE#CONTIG naming convention.

      bash Code Copy

      cactus-preprocess ${MYJOBSTORE} KPR-mc.seqfile KPR-mc.pp.seqfile \
          --configFile ./config_cut_hash.xml --realTimeLogging \
          --batchSystem mesos --provisioner aws --defaultPreemptable \
          --nodeType r5.8xlarge --nodeStorage 500 --maxNodes 2 \
          --skipMasking --logFile KPR-mc-grch38.pp.log

   5. Construction of Pangenome

      A pangenome was constructed by generating a minigraph using cactus-minigraph.

      bash Code Copy

      cactus-minigraph "${MYJOBSTORE}" KPR-mc.pp.seqfile \
          "${MYBUCKET}/KPR-minigraph-grch38.gfa.gz" \
          --reference GRCh38 \
          --realTimeLogging \
          --batchSystem mesos --provisioner aws --defaultPreemptable \
          --nodeType r5.8xlarge:1.5 --nodeStorage 650 --maxNodes 25 \
          --betaInertia 0 --targetTime 1 \
          --logFile KPR-minigraph-grch38.log