Skip to content

Assembly facilitated by GEP2

Jump to bottom
Diego De Panis edited this page Oct 17, 2026 · 3 revisions

1. processing reads for contiging and profiling

This is a genome I assembled for SIB under the ERGA project… The reads were already on ENA, so we can use GEP2 to get the reads and do a quick pre-processing and reads QC, and also a genome profiling based on the reads alone...

  • The table should only contain the species name and the read types and files (asm_id and asm_files columns empty). The google table (or whatever you are using) should look like this:
sp_name asm_id skip asm_files read_type read_files
Rosalia alpina hifi ERR16037030, ERR16037031
Rosalia alpina hic ERR16029358

  • On the control_panel.yaml file, everything should be off except the related to the # Read processing and # K-mer analysis sections

  • Run GEP2 and check the results.

2. assembly with hifiasm

After checking genoscope2 and the quality of the reads, we can run hifiasm using the trimmed reads ;)

Hifiasm is not part of the pipeline but we use it for correcting ONT reads -optionally-, and therefore is included in the container, so you could run it like:

cd GEP2
source containers_vars
apptainer exec -B {path/to/data}:{path/to/data} $GEP2_CONTAINER hifiasm -h

3. assembly qc

For this assembly, I ran hifiasm, on three different modes (L1, L2 and L3) and with the --primary flag too, so for each mode I have both pri/alt and Hap1/Hap2. The idea now is to check which one is better to decide how to continue with the assembly.

So the table should be updated with each combination of mode and assmebly type (e.g., L1 pri/alt, L1 hap1/hap2, L2 pri/alt...), like this:

sp_name asm_id skip asm_files read_type read_files
Rosalia alpina hifiasm_l1_pa /hifiasm/l1/icRosAlpi1_l1.hic.p_ctg.fasta.gz, /hifiasm/l1/icRosAlpi1_l1.hic.a_ctg.fasta.gz hifi ERR16037030, ERR16037031
Rosalia alpina hifiasm_l1_pa /hifiasm/l1/icRosAlpi1_l1.hic.p_ctg.fasta.gz, /hifiasm/l1/icRosAlpi1_l1.hic.a_ctg.fasta.gz hic ERR16029358
Rosalia alpina hifiasm_l1_h /hifiasm/l1/icRosAlpi1_l1.hic.hap1.p_ctg.fasta.gz, /hifiasm/l1/icRosAlpi1_l1.hic.hap2.p_ctg.fasta.gz hifi ERR16037030, ERR16037031
Rosalia alpina hifiasm_l1_h /hifiasm/l1/icRosAlpi1_l1.hic.hap1.p_ctg.fasta.gz, /hifiasm/l1/icRosAlpi1_l1.hic.hap2.p_ctg.fasta.gz hic ERR16029358
Rosalia alpina hifiasm_l2_pa /hifiasm/l2/icRosAlpi1_l2.hic.p_ctg.fasta.gz, /hifiasm/l2/icRosAlpi1_l2.hic.a_ctg.fasta.gz hifi ERR16037030, ERR16037031
Rosalia alpina hifiasm_l2_pa /hifiasm/l2/icRosAlpi1_l2.hic.p_ctg.fasta.gz, /hifiasm/l2/icRosAlpi1_l2.hic.a_ctg.fasta.gz hic ERR16029358
Rosalia alpina hifiasm_l2_h /hifiasm/l2/icRosAlpi1_l2.hic.hap1.p_ctg.fasta.gz, /hifiasm/l2/icRosAlpi1_l2.hic.hap2.p_ctg.fasta.gz hifi ERR16037030, ERR16037031
Rosalia alpina hifiasm_l2_h /hifiasm/l2/icRosAlpi1_l2.hic.hap1.p_ctg.fasta.gz, /hifiasm/l2/icRosAlpi1_l2.hic.hap2.p_ctg.fasta.gz hic ERR16029358
Rosalia alpina hifiasm_l3_pa /hifiasm/l3/icRosAlpi1_l3.hic.p_ctg.fasta.gz, /hifiasm/l3/icRosAlpi1_l3.hic.a_ctg.fasta.gz hifi ERR16037030, ERR16037031
Rosalia alpina hifiasm_l3_pa /hifiasm/l3/icRosAlpi1_l3.hic.p_ctg.fasta.gz, /hifiasm/l3/icRosAlpi1_l3.hic.a_ctg.fasta.gz hic ERR16029358
Rosalia alpina hifiasm_l3_h /hifiasm/l3/icRosAlpi1_l3.hic.hap1.p_ctg.fasta.gz, /hifiasm/l3/icRosAlpi1_l3.hic.hap2.p_ctg.fasta.gz hifi ERR16037030, ERR16037031
Rosalia alpina hifiasm_l3_h /hifiasm/l3/icRosAlpi1_l3.hic.hap1.p_ctg.fasta.gz, /hifiasm/l3/icRosAlpi1_l3.hic.hap2.p_ctg.fasta.gz hic ERR16029358
  • On the control_panel.yaml file, now switch RUN_COMPL: on too.

⚠️ IMPORTANT: if using google tables, every time you update the table you need to remove the data_{}.csv file in the GEP2_results folder so the pipeline will check the google sheet again and see your changes.

  • Run GEP2, the results will be in GEP2_results/Rosalia_alpina/

4. scaffolding files generation

After checking the nice md reports, the best assembly ended being hifiasm_l2_h. The idea now is to use GEP2 to generate the files necessary for doing the scaffolding with YaHS.

  • We want to run things only on that assembly, so we need to play with the skip flags in the table and the config. On the table, add in the on in the skip column for all the other assemblies, only hifiasm_l2_h should be off (or empty).

  • On the control_panel.yaml file, turn RUN_HIC: on and SKIP_HIC: on (now, that analysis will be skipped for everything that has skip activated on the table!). Set HIC_MAPQ: 10 (we are going to use the files generated for scaffolding in the next step, so MAPQ 0 would not be good for that use!). turn all the tracks off.

  • Run GEP2. The pipeline will generate .pairs.gz files for each one of the two assemblies using chromap.

  • The next step is use those files to generate the input needed for the scaffolder yahs using the script GEP2/workflow/scripts/get_files_for_yahs.sh. You can run it like (check the help!):

cd GEP2
source containers_vars
apptainer exec -B {path/to/data}:{path/to/data} $HIC_CONTAINER workflow/scripts/get_files_for_yahs.sh -h

The actual command is something like (this is only for hap1, you need to run the same for hap2):

apptainer exec -B /scratch/ddepanis/:/scratch/ddepanis/ $HIC_CONTAINER bash /home/ddepanis/GEP2/workflow/scripts/get_files_for_yahs.sh /scratch/ddepanis/ASSEMBLIES/Rosalia_alpina/hifiasm/l2/icRosAlpi1_l2.hic.hap1.p_ctg.fasta.gz 151 /scratch/ddepanis/ASSEMBLIES/Rosalia_alpina/GEP2_results/Rosalia_alpina/hifiasm_l2_h/hic/icRosAlpi1_l2.hic.hap1.p_ctg/icRosAlpi1_l2.hic.hap1.p_ctg.pairs.gz /scratch/ddepanis/ASSEMBLIES/Rosalia_alpina/yahs 2

When the script finishes, it points how you should run yahs (which is also included in the container, together with purge_dups that you also are welcome to use if needed).

⚠️ IMPORTANT: A neat idea if you have hap-phased assemblies, is to combine the scaffolded assemblies into a single fasta to enable a dual-hap curation later... Do something like sed '/^>/s/$/_Hap1/‘ inFile > outFile for renaming before combining!

5. qc and curation files generation

Now that we have the scaffolded hap assemblies and the combined scaffolded haps, we need to update the table with those two

  • Remember setting skip on to everything except to the new rows (those should be off).

  • On the control_panel.yaml file, set HIC_MAPQ: 0 (this is better for curation, and esencial if dual-hap). Turn on all the tracks, turn on inspector and the contamination stuff, also turn skip on to all those.

  • Run GEP2 and enjoy the qc and curation files

GEP2 wiki

Clone this wiki locally