Skip to content

Helixer commands

Jump to bottom
Mariano Olivera Fedi edited this page Apr 25, 2024 · 19 revisions

Helixer commands

Pre-Processing

Go to the directory for this activity.

cd /home/{user}/train/Annotation_workshop/Helixer

Reformat the fasta file with the genome to an h5 format to be used by Helixer.

fasta2h5.py --species A_thaliana --h5-output-path Format_conversion/A_thaliana.h5 --fasta-path Inputs/Genome/Athaliana_447_TAIR10.Chr3.fa --subsequence-length 106920

Analysis of model effect in the output.

Use all the models available to predict protein coding genes.

for f in Inputs/Models/* ; do tag=$(basename $f | rev | cut -f2- -d "." | rev) ; HybridModel.py --load-model-path $f --test-data Format_conversion/A_thaliana.h5  --overlap-offset 53460 --core-length 106920  --cpus 2 --overlap --val-test-batch-size 32 -v -p Predictions/${tag}_predictions.h5 ; done

Transform the predictions into a gff format file annotation.

for f in Predictions/* ; do tag=$(basename $f | rev | cut -f2- -d "." | rev) ; helixer_post_bin Format_conversion/A_thaliana.h5 ${f} 100 0.1 0.8 60 Model_analysis/Annotations/A_thaliana.helixer_${tag}.gff ; done

Use Mikado to generate metrics for each annotation.

for f in Model_analysis/Annotations/* ; do tag=$(basename $f | rev | cut -f2- -d "." | rev) ; mikado util stats $f Model_analysis/Mikado_stats/${tag}.stats ; done

Sumarise the statistics generated by Mikado.

for f in Model_analysis/Annotations/* ; do tag=$(basename $f | rev | cut -f2- -d "." | rev) ; parse_mikado_stats Model_analysis/Mikado_stats/${tag}.stats > Mikado_stats/${tag}.stats.summary ; done

Concat all the summary files into one general file.

Scripts/paste_mikado_summary_stats.sh Mikado_stats/*stats.summary  > Model_analysis/Mikado_stats/all_gff.summary

Visualise the final output.

cat Model_analysis/Mikado_stats/all_gff.summary
	Mikado_stats/A_thaliana.helixer_fungi_v0.3_a_0100.stats.summary	Mikado_stats/A_thaliana.helixer_invertebrate_v0.3_m_0100.stats.summary	Mikado_stats/A_thaliana.helixer_land_plant_v0.3_a_0080.stats.summary	Mikado_stats/A_thaliana.helixer_vertebrate_v0.3_m_0080.stats.summary	
Number of genes                   7740.00   5709.00   5454.00   3143.00
Number of Transcripts             7740.00   5709.00   5454.00   3143.00
Transcripts per gene              1.00      1.00      1.00      1.00
Number of monoexonic genes        1874.00   833.00    1015.00   256.00
Monoexonic transcripts            1874.00   833.00    1015.00   256.00
Transcript mean size cDNA (bp)    1310.96   1597.91   1583.62   1596.03
Transcript median size cDNA (bp)  1146.00   1452.00   1419.00   1366.00
Min cDNA                          62.00     62.00     71.00     62.00
Max cDNA                          15907.00  15832.00  15581.00  15678.00
Total exons                       36231.00  33331.00  29635.00  25066.00
Exons per transcript              4.68      5.84      5.43      7.98
Exon mean size (bp)               280.06    273.69    291.45    200.12
CDS mean size (bp)                238.22    230.80    239.68    166.95
Transcript mean size CDS (bp)     1049.20   1222.30   1236.45   1244.39
Transcript median size CDS (bp)   858.00    1044.00   1065.00   999.00
Min CDS                           60.00     60.00     60.00     60.00
Max CDS                           15495.00  15294.00  15234.00  15234.00
Intron mean size (bp)             113.34    172.90    159.33    229.17
5'UTR mean size (bp)              106.44    132.71    133.69    130.11
3'UTR mean size (bp)              155.32    242.89    213.48    221.53

Use Mikado to compare the annotations to the reference annotation.

for f in Model_analysis/Annotations/* ; do  tag=$(basename $f | rev | cut -f2- -d "." | rev) ; mikado compare -r Inputs/Annotation/Athaliana_447_Araport11.Chr3.gene_exons.gff3  -p $f  -o Model_analysis/Mikado_Compare/${tag} ; done

### Summarise all comparisions.

mikado util collect_compare -fmt tsv -l "all" -o Model_analysis/Mikado_Compare/Athaliana_helixer_models Model_analysis/Mikado_Compare/A_thaliana.helixer_*.stats

Visualise the the F1 metrics of all the comparisons.

less -S Land_plant_analysis/Mikado_Compare/Athaliana_helixer_models.f1.tsv

Analyse the effect of post-bin parameter.

Generate multiple annotation files using a different peak-threshold parameter.

for f in {0.8,0.85,0.9,0.95,0.975} ; do  helixer_post_bin Format_conversion/A_thaliana.h5 Predictions/land_plant_v0.3_a_0080_predictions.h5 100 0.1 ${f} 60 Land_plant_analysis/Annotations/A_thaliana.helixer_${tag}_land_plant_v0.3_a_0080.gff" ; done

Use Mikado to generate metrics for each annotation.

for f in Land_plant_analysis/Annotations/* ; do tag=$(basename $f | rev | cut -f2- -d "." | rev) ; mikado util stats $f Land_plant_analysis/Mikado_stats/${tag}.stats ; done

Sumarise the statistics generated by Mikado.

Scripts/paste_mikado_summary_stats.sh Land_plant_analysis/Mikado_stats/*stats.summary > Land_plant_analysis/Mikado_stats/all_gff.summary

Visualise the final output.

cat Land_plant_analysis/Mikado_stats/all_gff.summary
	Land_plant_analysis/Mikado_stats/A_thaliana.helixer_0.85_land_plant_v0.3_a_0080.stats.summary	Land_plant_analysis/Mikado_stats/A_thaliana.helixer_0.8_land_plant_v0.3_a_0080.stats.summary	Land_plant_analysis/Mikado_stats/A_thaliana.helixer_0.95_land_plant_v0.3_a_0080.stats.summary	Land_plant_analysis/Mikado_stats/A_thaliana.helixer_0.975_land_plant_v0.3_a_0080.stats.summary	Land_plant_analysis/Mikado_stats/A_thaliana.helixer_0.9_land_plant_v0.3_a_0080.stats.summary	
Number of genes                   5395.00   5454.00   5133.00   4927.00   5295.00
Number of Transcripts             5395.00   5454.00   5133.00   4927.00   5295.00
Transcripts per gene              1.00      1.00      1.00      1.00      1.00
Number of monoexonic genes        998.00    1015.00   905.00    804.00    955.00
Monoexonic transcripts            998.00    1015.00   905.00    804.00    955.00
Transcript mean size cDNA (bp)    1592.21   1583.62   1631.69   1661.66   1609.17
Transcript median size cDNA (bp)  1425.00   1419.00   1463.00   1493.00   1440.00
Min cDNA                          71.00     71.00     71.00     71.00     71.00
Max cDNA                          15581.00  15581.00  15581.00  15581.00  15581.00
Total exons                       29444.00  29635.00  28626.00  28195.00  29144.00
Exons per transcript              5.46      5.43      5.58      5.72      5.50
Exon mean size (bp)               291.74    291.45    292.58    290.37    292.36
CDS mean size (bp)                239.68    239.68    239.88    238.37    239.99
Transcript mean size CDS (bp)     1241.99   1236.45   1271.50   1296.89   1254.65
Transcript median size CDS (bp)   1071.00   1065.00   1095.00   1119.00   1083.00
Min CDS                           60.00     60.00     60.00     60.00     60.00
Max CDS                           15234.00  15234.00  15234.00  15234.00  15234.00
Intron mean size (bp)             159.25    159.33    159.15    159.07    159.41
5'UTR mean size (bp)              134.94    133.69    139.33    141.75    136.81
3'UTR mean size (bp)              215.29    213.48    220.86    223.02    217.72

Use Mikado to compare the annotations to the reference annotation.

for f in Land_plant_analysis/Annotations/* ; do  tag=$(basename $f | rev | cut -f2- -d "." | rev) ; mikado compare -r Inputs/Annotation/Athaliana_447_Araport11.Chr3.gene_exons.gff3  -p $f  -o Land_plant_analysis/Mikado_Compare/${tag} ; done

### Summarise all comparisions.

mikado util collect_compare -fmt tsv -l "all" -o Land_plant_analysis/Mikado_Compare/Athaliana_helixer_models Land_plant_analysis/Mikado_Compare/A_thaliana.helixer_*.stats

Visualise the F1 metrics of all the comparisons.

less -S Land_plant_analysis/Mikado_Compare/Athaliana_helixer_models.f1.tsv

Clone this wiki locally