Only the filename is correc. Conversion to org has not been done properly. Org-roam id are not generated.
M6GJ5MQ7PKVD5SKOL3VMTZAEGONBDDRQIXW3LGTJAVCQX77ADRCQC
***** Choix du génome de référence
- inclusion de zones hypervariables (HMC = histocmopatibility locas), PAR du chrY -> perte de l'aliggenent donc diminution sensibilité
- contigé non placés et non localisité : permet d'éviter le mauvais alignement de ces régions donc diminue les faux positif
- reco : référence "primaire" (chr, mitochondriae, contig non placé et non localisé)
- Exome = la plus populaire et la plus rentable (prix proportionnel aux nombre de variants)
- Liste appels de variants
***** Alignement
Privilégier sensbilité à la spécificité (ex: GATK)
***** difficulté : duplicats PCR
= reads du même fragment ADN pendant amplification PCR. Solution : pas d'amplification si suffisament de matériel, utilisation de "barcode" (UMI) si peu de materiel. Sinon on marque les reads qui partagent les même coordonnées comment potentientls doublons (mais risque de trop corriger dans les régions dupliquée) ou manquant les vrais doublons (régions répétée)
***** Algorithms
#+caption: Variant caller. Abbreviations: AS, assembly; CNV, copy-number variants; GATK, Genome Analysis ToolKit; MEI, mobile element insertions; RD, read depth; RP, read pairing; SR, split-read; SV, structural variants.
| Category | Name | Algorithm | Use | Paper |
|---------------------+----------------------+-------------------------------------+----------------------+-------------------------------------|
| Small variants | GATK Haplotypecaller | Local reassembly of haplotypes | Germline, MNPs | (Poplin, Ruano-Rubio, et al., 2018) |
| | BCFtools | Positional, pileups | Germline | (Danecek et al., 2021) |
| | FreeBayes | Haplotype-based, Bayesian model | Germline, MNPs | (Garrison & Marth, 2012) |
| | GATK Mutect2 | Local reassembly | Somatic | (Cibulskis et al., 2013) |
| | Strelka2 | Tiered haplotype model | Germline, somatic | (Kim et al., 2018) |
|---------------------+----------------------+-------------------------------------+----------------------+-------------------------------------|
| Structural variants | Delly2 | RP, SR, RD | Germline SVs | (Rausch et al., 2012) |
| | Pindel | SR, RP | Germline SVs | (Ye et al., 2018) |
| | Manta | SR, RP, AS | Germline, somatic | (Chen et al., 2016) |
| | GRIDSS2 | AS, SV Breakpoint | Somatic | (Cameron et al., 2021) |
| | Varscan2 | RD, Circular Binary Segmentation | Exome, somatic, CNVs | (Koboldt et al., 2012) |
| | EXCAVATOR2 | RD, In-,Off-target | Exome, CNVs | (D'Aurizio et al., 2016) |
| | ExomeDepth | RD, beta-binomial | Exome, CNVs | (Plagnol et al., 2012) |
|---------------------+----------------------+-------------------------------------+----------------------+-------------------------------------|
| Other | Mobster | RP, clipped reads | MEIs | (Thung et al., 2014) |
| | Expansion-Hunter | Reads spanning, flanking, in-repeat | Repeat expansion | (Dolzhenko et al., 2017) |
| | sideRETRO | SR, RP at insert | GRIP | (Miller et al., 2021) |
| | Harpak et al. (2017) | HMM model | NAGC | (Harpak et al., 2017) |
| | Li et al. (2019) | k-mer count, MDS | NUMT | (W. Li et al., 2019) |
**** Méthode:
- Génome : NA12878 (TruSeq), synthetic diploid (HiSeq) et simulé (NEAT)
- NEAT: 1 version avec mutations aléatoires (moyenne 0.002), 1 version défini à partir du (depuis le "golden truth")
- GRCh37
- 4 pipeline: bwamem+GATK, bwamem+deepvariant, Dragen+dragen, dragen+Deepvariant
*mais en GRch37*
- 64 coeurs
- comparaison avec happy + vcfeval
**** Résultat
- bonne concordance des pipeline :
- SNV: 91.7-9.6% retrouvés par tous les pipeline et 95.3-99.95% par au moins 2
- indel : 83.5-99.4%
- meilleur performance dragen vs gatk mais pas forcément mieux que deepvariant
- temps d'exécuton:
- dragen >> autres en total (facteur 5 et 30min à 1h30 d'exécution sur HPC)
- pour bwa-mem/gatk ou deepvariant, alignement domine appel de variant. Deepvariant ~ gatk pour appel de variant
Limite:
- résultat différent du fda challenge
- GRch37
Vérification de la syntaxe des variants selon HGVS
* Méthode
- 126 variants vérifés manuellement
-- 50 clinvar, dbsnp, cosmic, mycancer genome, LOVD, Emory genetics lboratory
-- 76 variants générés (mutalizer + variant viewer)
-- difficiles à annoter (selon expérience du laboratoire)
- Données clinvar (106 110) et COSMIC (2 215 000)
- left-justified (vt-normalize)
- snpeff, vep, variant reporter (et snpeff, vep sur données COSMIC et clinvar)
- grch37
- snpeff avec NCBI homo sapiens 105, vep avec refseq/ensembl
* Vérification de la syntaxe
Même transcrit et même syntax pour codant/protéine (ou "équivalent" si synonyme, ex c.482del et c482delT)
* Résultats
** Base de référence
- Transcrit différents: 4 pour SNPeff et 5 pour VEP
- Variation sur l'annotation proétique surtout (cf figure) mais si on considère les "synonymes", équivalent
- discordance sur indel mais insertion ou deletion ~ SNV (sic)
- sur 1 variant à un site d'épissage (splice site junction) : le variant est aligné sur la droite (c.1200 ou c:1201 au lieu de de c.1200-1)
- 4 variants avec discordance entre transcrit Refseq et séquence génomique (leur argument : refseq est plus "curé") -> VEP et SNPEFF les rendent comme des variants (fauxsense ou deletion)
- annotation ok pour substitution à 2 nucléotides
- effet prédit (si plusieurs effets, ok si au moins un correspond) : bons résultats (cf figure)
** Clinvar
- très bonne concordance snpeff et vep (mais principalement snv)
- plus mauvaise performance = insertion suprennament
- souci semble être l'échec du right-shift
- très grande variabilité pour syntaxe proétique (16-78% différence)
-- majorité = erreur systématique et correctible (ex: substitution au lieu d'indel p.AspAla2625GluPro versus p.Asp2625_Ala2626delins) et une petite partie annotent mal les délétion (ex: pIle55fs au lieu de p.Ile55Ter)
GluPro for rs267606668
- outils globalement bons
- s'ils sont d'accords, a permi de mettre en évidence ererur dans clinvar
* Conclusion
- outils globalement bons pour SNV mais vérifier les insertion ou délétion manuellement
- variabilité mais vep semble meilleur que snpeff
:PROPERTIES:
:ID: 5dbb8b67-ab7b-455a-8ed3-d657ee152941
:ROAM_REFS: @Yang2015
:END:
* 3 types d'annotations
- gène:
-- relation par rapport au gène : variant exoniqe, intronique
-- si exonique, rôle fonctionnel (synonmique, frameshift...) + changement acide aminé
- région: utile pour savoir si un variant est dans une région d'intérêt
-- segdup
-- cytogénétique
-- site d'attache pour transcription
-- microARN, snoARN
-- GFF3
* filtre : 1000G, dbsnp, prédicion, ESP6500, exac, clivan, GWAS
* Limites
1. filtre par défaut ("disease model") peut ne pas être appropripadapt
2. Parfois trop d'information (ex: 10 scores de pathogénétique) -> score MétaSVM. Recommandation: SIFT + metascore (metaSNVM) + CADD pour variant condant. Pour noncodant (phyloP, CADD)
3. pas de nomenclature uniforme
*
*
* Annotation
https://annovar.openbioinformatics.org/en/latest/
- Sur le gène Refseq par défaut (ENSEMBLE, UCSC, GENCODE ou GFF3 "maison")
- Basé sur la région:
--- régions conservée sur plusieurs espécèse
--- site prédiction de liaison pour facteur de transcription
--- bandes cytogénétique
--- Variants impactant microARN et snoARN (small nucleolar ARN)
--- Variants impactant site d'accroches prédits pour microARN
--- duplications segmentaire
--- variants structurels
--- variants publis dans études GWAS (genome-wide association studies)
--- zones ENCODE
--- variants non codant touchant enhancer, répsesseur, promoteur (prédit par chromHMM par exemple)
- peut annoter variants mitochondrial
- génome : hg19,hg38 (humain), souris, vers, mouche, levure + autres
- intronique, exonique, intergénique, 5'/3' UTR, épissage, upstream/downstram à 1kb
- si intronique, retourne les 2 gènes les plus proches
- Filtres :
- dbSNP, férquence allélique (1000 genomes, exac, gnomad)
--- score :
- SIFT, polyphen, lrt, muttaionteaster, mutation assesor, fathmm, metasnv, metalr
- CADD, GERP++
* Méthode
- processus en plusieurs passes
- génes "à filtrer" si contiennent variants causant un codon stop mais fréquent (MAF combinée > 1% dans 1000 génomes (par exemple: un variant avec MAF à 0.5 et un autre à 0.8%)
* Test
- 4 exomes avec 4 variants connus sur /MYH3/
Exemple d'utilisation de GATK best practive.
**** Hard-filter
:PROPERTIES:
:ID: 17aa8b97-b6d1-465e-bcdd-f9dafa9e8f2c
:END:
quand VQSR n'est pas possible. Exemple donné par
“QD < 2.0 || FS > 60.0 || MQ < 40.0 ||
HaplotypeScore > 13.0 || MappingQualityRankSum < -12.5 ||
ReadPosRankSum < -8.
- motivation: SNP = variant la plus commune -> besoin d'un catalogue pour étude d'assocation, fonctionnel, pharmacogénomique...
- single nucleatide substitution (majorité), polymorphism petit indel, région invariant (important aussi !), STR
- peut inclure des variants pathègène
- seulemen homosapiens avec 3,341,554,567 soumissions https://www.ncbi.nlm.nih.gov/projects/SNP/snp_summary.cgi, accédé 17 décembre 2023
- majrité = homo sapiens
Selon https://www.ncbi.nlm.nih.gov/books/NBK21088/,
- partie humaine = soumission du consortium SNP, extraction du génome (human gnome project), laboratoires (sans précisions)
- contexte fonctionnel (codant, intro, près gène...) calculé à partir de la séquence adjacente
- position dans la protéine si changement d'acide aminé
- calcul hétérozygosité moyenne (= probabilité que les 2 allèles soit dans un individu diploîde) selon allèles soumies
Rappel historique génome:
Différence par rapport aux autres assemblage: longue contig et scaffold, précision importante par base pair, représentation robuste région répétée et degmentally duplicately
- lié aux clones "BAC" (> 150kbp) et avec technique complémentaires (Radiation Hybrid Mapping = chromosome fragmenté par radiation =>_ si 2 marqueurs sont "loins", ils seront probablement sur des fragments différents; étude héritabilité (Genetic Linkage Maps); comparaison de "fingerprint" sur les clones [si similaire, viennent probablement de régions qui se recouvrents] Fingerprint Maps)
- cherche à représenter un "pan-génome" avec plusieurs donneurs
GRch37 :
- utilisé par le projet 1000 Genomes project pour représenter la diversité
- ajout de scaffold pour locus aletrnatif: représentation alternativ pour régions fortement variable MHC + haplotype divergenge aux locus /MAPT/ et /UGT2D/ tout en conservant la linéairité des chromosomes
- aprè 13 patches, mise à jour des coordénes -> GRCh38
GRCh38 = mosaic haplotype :
- réslution de 1000 problèmes (mettre fig 1)
- 75mb nouvelle séquence (2.3% total) avec 5Mb supprimé
- > 95% séquence totale, 98% séquence non centrmérique
- odnneur principale (probablement africcan-européan ) légèment diminué
- augmentation des "trous" mais lié au remplacement des centromère (=1 trou) en un scaffold
- augmentation nombre de gènes 41 722
- amélioration de l'alignement des transcrit
- test : > 1000 gèes perte de fonction pour autisme et un kit de > 4600 gène d'exome (plus utilisé) https://www.genomeweb.com/sequencing/emory-chop-harvard-develop-medical-exome-kit-complete-coverage-5k-disease-associ ->
- ajout de séquence alternative (pour augmen), qui améliorent l'alignement
- tests sur clinvar:
- test 2 : sur 113 000 variants clinvar, ~200 unique à GRCh37. Sont sur 9 régions curées. Parmi ces variants, exemple d'un variant patho de KCNE1 (rs1805128) qui en 38 s'aligne sur plusieurs endroits avec l'ajout d'une région paralogue....
En 2023, il est effectivement rendu "conflicting" https://www.ncbi.nlm.nih.gov/clinvar/variation/13479/?oq=((28518[AlleleID]))&m=NM_000219.6(KCNE1):c.253G%3EA%20(p.Asp85Asn) avec 23 soumissions
- test 1: alignement (bwa mem) + appel de variant pour NA24143
- 10 emplacement où la profondeur a été améliorée
- 135 emplacements oà la couverture a été "perdue" ! (région paralogue/alléliques)
(NB: paralogue = duplication d'un gène, avec séquence et fonction pouvant différer)
- impact sur aligement: 99.2% alignéen 38 mais 64% reads non mappés le sont. Ceux restant: 23% sont mappés sur les régions alternatives
Approches classique:
1. aligner reads sur référence et examine différences -> sensibe, peu coûteux mais souvent ne se concentre que un 1 type de variant (errours sur indel/large variant), peut échoeur sur régions qui divergent fortement, plutôt pour SNP (cuteux indel)
2. Assemblage sans réference : graph (dit "de Bruijn") et on cherche dans ce graphes des polymorphosimes. Avantage: gère les régions fortement divergentes, fonctionne pour haplotype au niveau local (plutôt qu'un variant), évite alignement. Mais coûteux, moins sensible, limité sur séquence répété (quand reads sont "coupés" in k-remrs)
3.Utilisation d'échantillions apparentés : diminue faux négatif
Apport de Platypus (mais vieux...) = combine ces 3 approches
- assemblage local et génération d'haplotyage
- variants candidats a partir des reads, de l'alignement et de base de données (polymorphismes connus)
- graphe de Bruijn à partir des reads puis extraction des chemin qui commencent et finisshes sur la séquence -> ce sont les allèle alternative candidatet (normalement non affectés par régions répéties, bloucles)
- répartir dans différents intervalles puis génération hapoltype (si trop d'haplotype, on ne garde que 256 les plus "prometteurs")
- estimation de la fréquence des haplotypes : matrice avc les probabilités pour chaque reads et chaque haplotype
- alignement selon un model qui utilise des probablités d'indel et SNP à partir des scores de qualité
- les haplotypes sont "coupés" en leurs variants et on calcule génotype ed probabilité en utilisant les variants de la région (somme au niveau proba, cf https://www.quora.com/What-is-marginalization-in-probability )
- filtre sur allèle, biais de brin, qualité alignement, qualité de la paire de base, contexet, probabilité après théorème de Bayess
#+caption: Platypus algorithm [cite:@Rimmer_2014]
[file:./img/platypus.png]
- plupart des "trous" insolubles à cause de de polymorphismes de structures aux extrémitées et nombreux éléments répétés/polymorphes non finis ou mal assemblé
- 151Mbp de séquencce inconnues:
- péri/sous-télomorique, duplication segmentaire, ADN ribosomal
- bras court des 4 chromosomes acrocentrique (13,14,15,21,22), HSat (human satellite repeat array)
- certaines régions sont artificielles : dans le centromère, certaines zones du centromère dite alpha satellite (Fig. \ref{fig:alpha-satellite} [fn:35] sont représenté sous forme d'un modèle informatique[fn:34]
- régions mal assemblées : bras court du chr21
Méthode:
- utilisation d'une mode hydatiform (Complete HM soit CHM) qui ont un caryotype 46XX (perte du complémente maternel avec duplication paternel)
- plus part est d'origine européen avec régions introgression nanderthale
- pas plus de variant perte de fonction que échantillions du projet 1000 génome, ni d'allèles isolée
- ultra-long read pour régions répétées (nanopore) et pacbio hifi pour read 20kb (compromis avec taux d'erreur de 0.1%), + autre (illumina pcr-free, illumina hi-c pour interaction de chromatine, bionano map optinque single cell Strand-seq)
Résultat:
- 8% du génome: ajout de 238Mbp (4.5% dont 182 exclusif T2T
- pas de trous, pas de contig
- +5% gène
- sur 1000 génomes : diminue faux positifs et faux négatifs
- manque 5 régions des ADN ribsomaux les plus répété (9.9Mbp)
Motivation
#+begin_quote
the widely used GATK uses logistic regression to model base errors, hidden Markov mod- els to compute read likelihoods, and naive Bayes classification to identify variants, which are then filtered to remove likely false positives using a Gaussian mixture model with hand-crafted features capturing common error modes5."
#+end_quote
GATK donc orienté Illumina mais algorithmes difficiles à porter sur autres technologies
Algorithme = deep learning
1. SNP et indel candidats avec techniques classiques (sensible mais peu spécifique) [comparaison à la réference ?]
2. Probabilité à chaque locus en utilisant les pileup de la référénce et des reads autour de cette position (cf figure \ref{fig:deepvariant})
3. Avec le modèle entrainé, donne les variants (entrainement du model se fait sur sur génotype connu -> version "fixée"
GATK supose que les erreurs de reads sont indépendantes. Le réseau de neurones (convolutional) prend en compte les dépendences complexes donc plus performant
Contient algorithme haplotypecaller
Figure: https://gatk.broadinstitute.org/hc/en-us/articles/360035531412-HaplotypeCaller-in-a-nutshell
Principe : réalignement local
1. Régions d'intérêt
À partir de l'alignement, compare à l'allèle de réference (en tenant compte des mismatch, indel, softclib) -> Calcule une "probabilité" pour chaque position (NB: en lissant avec un kernel gaussien).
Une région est donc un ensemble ed position telles que la probabilité > seuil (0.002 par défaut). NB: utilise des reads qui sont sont à +/-100pb autour de cette région mais qui ne contribuent pas au calcul de probabilité. Taille mi = 50bp, max 300bp (si supérieur, coupée en 2)
2. Construction d'un graphe pour tous les haplotype possibles (NB: si trop de k-mers (65 ), région exclues). Haplotypes avec seulement 2 reads exclus. Si le chemin dans le graphe ne "revient" pas à l'haplotype de référence, essai de fusion avec smith-waterman. Si echec, suppression du chemin.
*Chaque haplotype du graphe est réaligné avec smith-waterman*
3. Pour chaque read, calcul d'une probabilité pour chaque haplotype (paried Hidden Markov Model). Utilise BQSR (si on a suivi les best pratice) donc permet de corriger les erreur PCR
4. probabilité par read transformée en variant (modèle Bayeson).
- proba d'un génotype * probabilité "brute" d'un génotype / somme de probabilité de tous les génotype pour les allles d'un variant
- proba brutes d'un génotype = calcul sur tous les reads en faisant le produit des probabilité des différentes allèles pour ce génote pour chaque read
NB: plus de détails sur le calcul de probabilité avec un exemple simple https://eriqande.github.io/eca-bioinf-handbook/variant-calling.html
Apport : DeepVariant et Sention non testé
**** Méthodes
- Alignement : BWA-MEM
- Données :
- germline = génome(?)
- HG001...7 = Hiseq2500 avec 30x (sauf 40x pour lymphocytbe B HG001)
- HG001, 2 et 5 sur pacbio à 30x
- somatique hiseq2500 à 100X
- Génome de référence : hs37d5 (utilisé par 1000 génomes)
- 11 appel de variants
- Sentieon: Germline (DNAseq, DNAscope) + somatique (TNscop, TNseq)
- GATK : Germnline (haplotypecaller) + somatique (mutect2)
- Deepvariant : somatique
- neusomatic, varscan2 et strelka2 pour somatique
NB: third generation sequening pour DNAseq (sention) et haplotypecaller (GATK)
- référence : GIAB pour germline, données de synthèse pour somatique (mélange germline)
- rtgtool
**** Résultats (germline)
***** Next-generation sequencing
- Bonne concordance pour germline: F1 > 0.99 pour SNV et > 0.98 pour indel (cohérent avec biblio)
- Downsampling (2,5,10,15,30x)
- 30x: score F1 identique -> 30x recommandé
- < 15x précision semblable mais recall bas
- Hard filter (voir [[id:17aa8b97-b6d1-465e-bcdd-f9dafa9e8f2c][Hard-filter]]): diminution légère F1 avec diminution recall 0.001% et augmentation precision 0.001
- Modification de 2 paramtre (emit_conf, call_conf): F1 identique avec diminution precision et augmentation recall
- Impact BQSR
***** Third generation sequencing
- SNV : bonne concordance (F1 > 0.99) avec très faible différence, en accord avec [cite:@wenger2019]
- indel : variabilité importante deepvariant > dnaseq >> haplotypecaller
- 0.9902, 0.9927, 0.9924 for TGS001, TGS002 and TGS005,
- followed by DNAseq mode of Sentieon (0.9433, 0.9390, 0.9393),
- HaplotypeCaller from GATK (V4.0.7) )
- hard filter : nette amélioration F1 score pour indel
***** Comparaison next- et third generation sequencing
- 3Rd-gen : +2.13% de SNP, 3.89% indel
- confirmé pour régions GIAB: high-confidenc, filtered et en dehors
- Idem dans régions riches en GC e pauvre en GC
- Idem dans régions hautement répétées
- Confirmation manuelle sur IGV
***** Coût
Test 4 CPUS: facteur 5:
- dnaseq/score < 50 heures CPU , haplotypecaller ~ 150 et deepvariant ~200 (sic)
Méthodes
- Alignement sur génome de champignon
- BWA, Bowtie2, MUMmer4, STAR, HiSAt2
- Dual Xeon E5-2643 (six cores and 12 threads for each processor) with 512 GB of RAM
- timed as CPU time (user + system).
- données RNA-seq
RésultatS
- Qualité:
- taux d'alignement : BWA et bowtie2 = meilleur. MUMmer, STAR intermédiaire
- runtime par read : HISAT2 = le plus rapide, intermédiaire = BWA, bowtie2, star
- parallisésation: bowtie2 (0.999) et MUMmer4 (0.981= meilleur, bwa correcte avec (0.833), hisat 0.876
Utile pour scaleup
Les zones alpha satellites sont représentées à partir d'un modèle crée à partir de séquencage du génome. Cf Fig 1
* Contexte
- Interprétation souvent basée sur impact sur transcrit ou protéine
- Donc dépend de l'annotation du transcrit et de la localisationo du variant par rapport aux régions codantes pour des protéines ou 2
- 2 sources principales:
1. GENCODE : dépend du génome de référence, vise à représenter tous les isoformes pour tous les tissus et étapes du développement -> en moyenne 4 transcrit par gène codant. Parfois plusieur dizaines de transcrit pour un variant
2. Refseq (indépendant)
- Attention aux haplotype alternatifs dans le génome de référence : 1 variant peut être associé à plusieurs ALT mais n'avoir un codon que sur l'un d'entre eux...
Exemple: https://www.ncbi.nlm.nih.gov/snp/rs150580082#hgvs_tab C>G,T: la référence est C => codon stop en G seulement. Or il y a G et T sur des transcrits alternatifs
HGVS: basé sur trascnit donc confusion possible + nombreuses annotation possibles
* VEP
- annotation pour SNVi, indel, substiution multiple bp, microsatellite, tandem repeat + structural variant > 50 nucléotide (not. CNV)
- GRCh38 + 37 + T2T (unstable)
- effet sur trnascrit, protéine et région régulatrice
- si variant connu, allele fréquence + information (phénotype)
- utilisable dès qu'il y a un génome + ensemble de gène annoté
- chaque version est liée à une release d'Ensembl
- open source + libre d'utilisation
Annotation
** Annotation de transcrit
- Annotation selon
- Ensembl = GENCODE (fusion des prédictions de transcrits d'Ensemble basé sur des preuves + annotation manuel, pour L'homme)
- Refseq
- Une ligne par allèle et par caractéristique génome (transcrit, zone régulatrice...)
- Problème : pas de consensus sur filtre -> à faire par l'utilisateur
Sortie de VEP pour gène et transcrit (Table 2 article)
- identification du gène affecté : identifiant Ensembl, nom "commun" du gène (ex: HGNC)
- identifiant du transcript: Ensemble, NCBI Refseq
- identifiant pour CCDS (Consensus coding sequence)
- biotype selon GENCODE (codant pour protéine, pseudogène... voir https://www.gencodegenes.org/pages/biotypes.html )
- coordonnées du variant : cDNA, processed coding sequence (CDS)
- distance au transcrit s'il est en dehors des bornes
- conséquence sur transcrit
- nombre d'exons et introns touchés
- Transcript Support Level (TLS) qui indique la fiabilité des modèles de transcrit
- Annotation principle splice isoformos (APPRIS) pour l'annotation de transcrit sur épissage alternatif (modèles informnatique)
** Annotation protéique
Liste:
- identifiant proétiqiue : Ensembl, Refseq, UniProt (généré automatiquement, nettoyé la main ou combiné)
- coordonées protéique
- codon de référence et alternatif
- acides aminés de référence et alternatif
- score SIFT et PolyPhen2 (prédictif de pathogénicitié)
- domaines protéiques
- notation HGVS
Score pré-calculé pour toutes les combinaisons, mise à jour quand nécessaifre
Autres scores (FATHMM, mutationtaster... disponibles via plugins)
** Annotation variants non-codants
Impact si éléments régulateurs de la transcription/traduction -> VEP donne les ARNs non codant, zone régulatrices génomique ou motifs d'attache pour facteur de transcription (transcription factor binding motif) + changement sur score consensus.
Source = ensembl regulator build (ENCODE + BLUEPRINT + NIH epigenomics roadmap).
Autres scores disponibles via plugin: CADD
** Fréquence, phénotype
dbSNP + autres source not. COSMIC (somatique, HGMD + variants structures et CNV de "Database of Genomic Variants"
fréquence 10000 génome, exac
omim, orphanet, gwas
# Input, output
- input: VCF, identifiants (dbSNP, HGVS, refseq, ensembl)
- output: HTML, text (vcf, tsv, json)
# Performance
4 coeurs: 4millions de variants : 62min
Plus lent que la compétition:
- snpeff est écrit en java
- annovar : moins d'annotation, Perl
| Type | Temps d'exécution |
| ----- | ------------------|
| Annovar | 21min50 s (3415 v/s)
| SnpEff | 46min39 s (1598 v/s)
| SnpEff (threaded)* | 10min28 s (7121 v/s)
| VEP | 62min9 s (1200 v/s)
#+begin_quote
In brief, Minigraph builds a pangenome by starting from a reference assembly,
here GRCh38, and iteratively and progressively adds in additional assemblies,
recording only SVs ≥ 50 bases.
#+end_quote
**** Validation avec GIAB
1. Assemblage (Trio-Hifiasm)
2. Alignement sur GRch38
3. Appel de variant
4. Exclusion des régions "structurally variant": alignement reads pacbio en GRCh38, appel de variant puis fusion des VCFs (il faut au moins 2 calleur) puis dipcall (pas clair, mais voir méthode rDefining confident regions for variant benchmarking)
5. Comparaison avec GIAB
#+begin_quote
Evaluation using the GIAB benchmark
The small variant calls were evaluated using HG001–HG007 with the GIAB benchmark
(v.4.2.1)124, on HG002 in challenging medically relevant autosomal genes47, and
on HG002 using a preliminary draft assembly-based benchmark. For the draft
assembly-based benchmark, we used Dipcall38 to align a scaffolded, high-coverage
Trio-Hifiasm assembly21,38 to GRCh38 and call variants, and then we excluded
structurally variant regions from the dip.bed file as described above for the
benchmarking of small variants from the pangenome graph.
#+end_quote
Différente GRCh38, GRCH37: 1.5% SNV, 2% indel discondants
- fusion d'interprétations soumises par utilisateur pour signification clinique/fonctionnelle de variant
- germline + somatic, toutel se régions, toute les types , toutes les tailles, SNVs, CNV, cyto
- contexte recherche ou clinique
- OMIM + genereviews + laboratoire clinique, + recherche + uniprot
- archive: permet de voir l'évolution d'un variant au cours du temps
- équipe clinvar va vérifier HGS, relation gene-maladie mais pas li'tenprértation -> création de panels d'experts qui vont "nettoyer" les interprétations
- peut soumetter : laboratoire clinique, chercheurs, base de données, groupes d'expertes / ceux étabilssant les guideline
- agrégation des interprétations
- gratuit
#+caption: Statistique clinvar au 17 décembre 2023
#+name: tab:clinvar
| Records with assertion criteria | 3310698 |
| Total genes represented | 92082 |
| Unique variation records | 2425225 |
| Unique variation records with interpretations | 2414095 |
| Unique variation records with assertion criteria | 2309591 |
| Unique variation records with practice guidelines (4 stars) | 663 |
| Unique variation records from expert panels (3 stars) | 15598 |
| Unique variation records with assertion criteria, multiple submitters, and no conflicts (2 stars) | 345874 |
| Unique variation records with assertion criteria (1 star) | 1841230 |
| Unique variation records with assertion criteria and a conflict (1 star) | 106226 |
| Unique variation records with conflicting interpretations | 106544 |
| Genes with variants specific to one protein-coding gene | 17210 |
| Genes included in a variant spanning more than one gene | 92411 |
| Variants affecting overlapping genes | 35523 |
| | |
Total submitters 2704
**** Méthode :
20 pipelines
- 4 appels de variants: deepvariant, samtools, freebayes, gatk
- 5 aligneurs: bowtie, bwa, novoalign, SOAP,mosaik
- GIAB: NA12878, NA24385, NA24631
- données simulées avec ART
- GIAB 37 *et 38* pour NA12878
**** Biblliography
En désaccord avec article plus récents apparements:
#+begin_quote
SAMtools is best for Ion Proton data [6], and GATK is best for Illumina data [7].
#+end_quote
#+begin_quote
They have also shown low concordance when examining the same set of sequencing data. Thus the accuracy of the variant callers is still not adequate
#+end_quote
#+begin_quote
Applying multiple tools can result in more misleading output
#+end_quote
#+begin_quote
It has also been reported that read aligners influence the accuracy of variant de- tection
#+end_quote
**** Résultats
NA12878 (réel + simulé, 37 + 38):
- 4 meilleurs sont BWA_DeepVariant, Novoa- lign_DeepVariant, BWA_SAMtools and Novoalign [F-score en 0.97 et 0.99]
- pipelines basés sur deepvariant > GATK
- indel: bwa+deep variant et novoalign+deep variant sont les meilleurs, puis gatk+bwa et gatk+novoalign
- profondeur:
- profil similair SNV et indel.
- plupart des SNVs détectés à 150x
- Genotype quality (GQ): meilleure preformance avec augmentation GQ mais deepvariant + {bwa, novoalign} ont de bonnes performances même bas GQ
- heterozygote/homozygote supérieur pour SNV que pour indel (1.6-1.5 vs 1.2-1.3)
- transitio/transversion: 3.4-3.3[fn:4]
- détection indel : deepvariant et gatk avec bwa/novoalign mais attention variable selon la longeur (échec de détection certaines longeur)
- amélioration des performances en fusionnant des 4 meilleurs pipeline : 99-98% précision pour SNs, 96-98% pour indel
- faux négatifs : 0.5-1.5% (SNV) ou 0.5-4% pour indel. Ils sont dans des zones avec profondeur < 30x et < 10GQ
- NA24385 and NA24631 = F-score identique, deevariant reste le meilleur avec bwa/novoalign
**** discussion : qualité : deepvariant > samtools > gatk
- indel moins performant: données d’exome (loupe grand indel ?)
- gold standard na12878: deepvariant > autre (contradiction avec bibliograhie)
- Influence aligneur :
- bwa = temps d’exécution, mémoire et précision sont équilibré
- novo align = lent et consommateur ed mémoire mais mappe mieux
- appel de variant a plus d’impact
- meilleure performance en gRCh38 que 37 (plus de vrai positifs ?). Faux négatifs : -8%) et -20% pour SNV et indel respectivement
**** Problème
ART ne simule que du génome. Ils ont du restreinte le fasta (mais non précisé dans le github ...)
**** Données
- SRR pour données d'exomes
- Scripts pour les différents pipeline en bash : https://github.com/bharani-lab/WES-Benchmarking-Pipeline_Manoj
Fait l'alignement + appel de variant.
Peu de détails sur l'appel de variant. Le meilleur alignement de chaque read est examiné pour changement. Les variants sur plusieurs reads sont fusionnés.
Il y a des seuils sur la couverture, qualité, fréquence, nombre de reads définit par défaut.
Combine 2 appproches
1. algorithme de Freebayes [cite:@garrison2012]
2. réassemblage local (même idée que GATK apparement)
Valeur ajoutées
- filtre "maison" pour enlever le "bruit" du séquenceur ("novel heuristic filter for sequencer phasing noise")
- évite le coût d'une "pair Hidden Markov model" avec
1. évite de faire tous les alignements locaux possibles
2. utilise une probabilité maximal locale (comment ?) pour la probabilité "marginale" (qui est le produit de toutes les combinaisons)
- score "maison" pour le variant basé sur du machine learning (pour filtre ou priorisitation)
- contient proba du génotype, qualité de l'alignement (root-mean-square), biais de brin, fruction de reads correspondant aux modèle et complexité de l'alignement (compressiblitié, longeur homopolymer
- pré-entrainé sur NA12878 (génome)
#+begin_quote
Strelka2’s haplotype model uses an efficient tiered scheme for haplotype discovery, combining the advantages simple model based on input read alignments3 with those of a more complex model using local assembly2,5,6 , with the appropriate method selected according to the properties of each variant locus.
#+end_quote
Limite de GRCh38:
- toujours incomplet : 150 000 "trous" -> 995
- ne représente pas la diversité du génome (70% = 1 individu et le reste = 19 autres)
- difficulté : duplication > taille du reads, SNV et variant structurelle mal capturés par short-read : enzymes utilisées ont des difficultés avec structures complexe (régions riches en GC pour les promoteur)
- problème de la réprésentation haploid (= gène représenté une seule fois):
- variants de différents haplotype sont regroupés dans le même (pseud-haplotype))
- fausse duplication: soit par un mauvais alignement de régions hétérozygotes qui apparemment comme des paralogue (gènes provenant d'une duplication initiale) (= "false heterotype duplication"), soit quand des régions qui devraient être fusionnées ne le sont pas et sont identifiées comme des duplications (="homotype false duplications), par exemple erreur de séquencage sur long reads
- une représentation diploide (= gènes des 2 parents) pour chromosome X et Y et gènes soumis à empreinte
Ici T2T semi-automatisé appliqué à HG002
Suite de Hwang 2015
Aligneur + appel de variant mais utile surtout pour appel de variant
Méthodes
- 70 combinaisons: 7 aligneurs avec 10 appel de variant
- 1 génome européen et 1 african
- variants de référénces :
- GIAB et platinum genomes pour NA12878
- 1000 genompes pour les 2
- distance de Jaccard [fn:2] entre 2 pipeline sur 4 figures (SNP et indel, NA12878 et NA19240 respectivement)
- avec les référence donc 70 +1 ou +3
**** Résultats
- haplotypecaller est assez groupé (probablement car réaligne localement)
- plus de différence entre le type de variant (SNP vs indel) qu'entre les individus
- *majorité des variants sont appelé par plupart des pipeline*
- discordance expliquée par la profondeur ?
- pas pour SNP + indel de NA12878 (plus important pour les concordances !). Idem pour NA19240 sauf pour indel homozygote
- discordance expliquée par la balance allélique ?
- oui pour les indels. Non SNP (supérieur pour SNP Cconcordant sur NA19240)
- plus de discordance entre les pipeline pour patient NA19240: analyse statistque :
- minor allele frequnce++, impact fonctionnel prédit, éléments éléments, contenu GC, profondeur, mapping quality
- plus de discordance pour faible fréquence (MAF[fn:1] 0.5-5%), variants rares (MAF < 0.5)
- plus de discordance pour variants avce un impact fonctionnel important prédit par VEP
- zone répétes effet déléter pour SINE, simple repeat, low complexité. Notamment indel dans zone répétes par rapport aux indel.
- GC influence la concordance (dans quel sens ?)
- profondeur et MAPQ des short reads influence significativement :
- régions bien couverte ont une concordance x 1.08-1.28 [logicque]
- MAPQ : x1.76
- Selon la référence :
- 1000 genomes:
- glftool : bonne sensibilité et VPP pour SNP
- haplotypecaller : bonne performance (?) pour indel
- GIAB et platinum genome :
- haplotypecaller : bonne sensibilité et VPP (sauf pour SNV platinum)
- différence due à la région génomique ? GIAB = 90% du génome et 97% pour platinum. 1000 genome : tout -> plus de variabilité dans les performances
- combinaison d’appel de variant = performance supérieure dans la bibilé. Test ici :
- méthode: variants filtré
- soit par "call condordance" (i.e appelé par plusieurs pipeline ?)
- soit par 6 facteurs qui l’influent : MAF, impact fonctionnel prédit, éléments répété, contenu GC, profondeur, MAPQ (sélection sur régression logistique)
- résultat:
- pas de vrai gain de performance pour GIAB et platinum (sauf SNV de platinum).
- Pour GIAB et indel platinum = haplotypecaller a une performance supérieure ou égale que les ensemble.
- 1000 genome : ensemble de pipeline = meilleur.
- Le meilleur pipeline solo = gsnap + glfsingle (2 samples)
- pour indel : bwa mem + haplotypcaller (2 samples)
- -> certains pipeline exploitent les régions hors "zones de confiance" de GIAB
**** Conclusion
- *influence appel de variant semble > aligneur* (mais peut être lié aux différentes options)
- discordance pour variants rare et nouveau et pour SNP rare
- probable diminution des performances pour population non-européennes (plus de SNP et indel rare)
- *BWA + haplotypecaller = VPP et bonne sensibilité pour NA12878*
- limitations
- comparaison de vCF = left-normalized seulement (différentes représentation d’un variant complexe non considérée)
- pas el même séquencage pour les 2 patients (profondeur 49 et 72, short-read length 101 vs 250)
- pas d’optimisation des paramètre
*BWA-mem + haplotype caller n’est pas inférieur à une combinaison de pipeline pour la plupart des région*
Conclusion: méthode solide
Apport:
Article précédant hwang2019
Différents pipeline sur données d’exome NA12878 (3 aligneur et 4 appels de variants)
Context : résultat différent pour diffents outils d’appels de variant
**** Méthodes
- NA1287 sur Hiseq200 (7 données), Hiseq2500 (4 données) et ion proton (1 seul):
| HiSeq2000 | SRR1611178 | SeqCap EZ Human Exome Lib v3.0 | WES | 79.93x |
| HiSeq2000 | SRR1611179 | SeqCap EZ Human Exome Lib v3.0 | WES | 79.84x |
| HiSeq2000 | SRR292250 | SeqCap EZ Exome SeqCap v2 | WES | 116.06x |
| HiSeq2000 | SRR515199 | SureSelect v4 | WES | 298.45x |
| HiSeq2000 | SRR098401 | SureSelect v2 | WES | 116.84x |
| HiSeq2500 | SRR1611183 | SeqCap EZ Human Exome Lib v3.0 | WES | 129.94x |
| HiSeq2500 | SRR1611184 | SeqCap EZ Human Exome Lib v3.0 | WES | 111.90x |
| HiSeq2000 | ERR194147 | UCSC Known gene | WGS | 45.68x |
| HiSeq2000 | SRX485062 | UCSC Known gene | WGS | 56.60x |
| HiSeq2500 | SRX515284 | UCSC Known gene | WGS | 56.87x |
| HiSeq2500 | SRX516752 | UCSC Known gene | WGS | 43.61x |
| IonProton | NA12878_combine | UCSC Known gene | WGS | 9.87 |
| | | | | |
- alignement sur GRCh37
- vcflib pour normalisation
- restriction aux zones de capture pour données d’exomes
- comparaison avec GIAB: precision et recall calculé "à la main" (article avant hap.py)
- aire sous la courbe precision-recall (!= aire ROC)
- 13 pipeline : combinaison
- aligneur : Bwa-mem, bowtie2, novoalign
- appel de variant : haplotypcaller, samtools, freebays, ioproton variant caller
**** Résultats
CLOSED: [2023-10-23 lun. 13:33]
- Appel de variant : variabilité suivant les échantillons...
- Impact de l’aligneur et appel de variant : appel de variant >> alignement pour SNP et indel (surtout indel): average standard deviation: 3.46e-3 (aligneur) and 4.02e-3 (appel de variant) for SNPs and 0.72e-2 (aligneur) and 7.2e-2 (appel de variant)
- illumina:
- SNP : samtools + bwa-mem
- indel : haplotypecaller + n’importe quel aligneur
- ion proton : samtools
- concordance: metter figure 3
https://www.nature.com/articles/srep17875/figures/3
- illumina 92% (intersection des 3 appels): explication possible: différence de version, pipeline
- variabilité sur données 82-97% opevrlap !
- biais:
- ignorer allèle de référence (hétérozygote considéré comme homozygote): 7290/19851 -> freebayes sur illumina
- ignorer allèle alternative (homozygote considéré comme hétérozygote ): 9917/19851 -> haplotypecaller
- autre: 2644
- donc attention aux hétérozygote avec haplotypecaller
- Filtre sur la qualité (qualité < 20) peut impacter la performance. HaplotypeCaller ne rende que les > 30
- Profondeur n’affecte pas les résultats (significativement)
Revue plus courte
Short-read: Illumina (MiSeq, NextSeq, NextSeq...) et Thermofirschser (Ion Torrent).
- amplification clonale = soit avec des "billes" (beads) (Io torrent) our sur des flow cells. L'ampliication se fait par "emulsion PCR" (ion torrent) ou "bridging" (Illumina) ()"solid phase"
- ion torrent: limites = homopolymère mais taux d'erreur par bp < 0.1%
- illumina: taux d'erreur 0.1% mais temps d'exécution long sur miseq, fortement dimunié sur génération suivantes (https://www.illumina.com/systems/sequencing-platforms.html)
# TODO image de [cite:@goodwin2016]
Long-read: > 10kb mais avec des taux d'erreurs initialeent trè élevé.
- Pacbio (Pacific Biosiciences)): + long read (5% > 135kb) (utile pour HLA), données en temps réals, pas d'amplification PCR. Taux d'erreur initialemen très élevé (14%) mais combiné avec du short-read, diminution. Coût.
- Oxford nanopore (ONT)): reads > 1mb, même avantages que PacBio mais moins cher, plus petits et portables.
Vieil article mais plus exhaustif
Méthodes
- Bowtie, Bowtie2, BWA, SOAP2, MAQ, RMAP, GSNAP, Novoalign, and mrsFAST (mrFAST)
- données génomiques de synthèse: wgsim (erreur de distirbution uniforme) et ART (biais illumina)
- données RNA-seq
- génome de souris ("similair behaviour" que sur le génome humain
Résultat:
- performance = taux de reads alignés
- throughput (bp mappé par seconde)
- SNP: GSNAP > bowtie2 = bwa
- comparaison autre génomes animaux : performances semblables entre les outils
- sensibilité : taux de reads aligné correctement : bowtie (93)>BWA > bowtie2 (85%)
- scalabilité :
- 8 threads : bowtie < bowtie2 > bwa
- 8 coeurs : bowtie = bwa = SOA = FANGS
Conclusion : pas de meilleur aligneur pour toutes les métriques. Bowtie = meilleur throughput. BWA = meilleur pour reads plus long
Intro (pour thèse)
- séquencage initial du génome = long et coûteux
- séquencage haut-début milieu des années 2000 = révolution sur le temps et côut
- coût génome < 1000$ [cite:@dnacost]
- limitations : taux d'erreur 0.1-1.5% et reads de petite taille (35-700bp) par rapport au Sanger
Technologies
Short-read
Domination d'Ill%umina avec précision > 99.5% mais
- sous représendtation région riche AT-GC et tendance aux erreurs de substitution
**** Principe
Théorème de Bayes pour détecter des haplotypes à partir des reads.
Principe: sur un locus donné pour R_1... R_n, on calcule la probabilité d'un génotype G à partir d'hypothèse sur la distribution d'alèlle a priori P(G) et sur la qualité de séquencage P(R_1... R):
P(G|R_1...) = P(G)P(R_...|G)/P(R_1..)
P(R_i|G_i) est la probabilité d'un génotype à partir d'un ensebmel de reads. Calculé à partir de la probabilité d'avoir R_i échantillons à partir de G multilié par les erreurs de séquencage (calculé à partir des scores de qualités)
P(G) est la probabilité d'échantilloner le génotype en étant donné la fréquence de l'allèle, multiplié par la probabilité d'observer cet allèle dans la population (Théorême de Bayes également)
Le second terme est donné par la loi d'Ewens (correspond à la distribution de probabilité des décompositions d'un entier, voir https://en.wikipedia.org/wiki/Ewens%27s_sampling_formula)
**** Implémentation
1. Si les variations sont "assez proches" (moins de n bp selon un seuil défini), combiné dans un seul haplotype
2. Définition d'un intervalle.
On filtre les allèles selon le nombre de reads portant ALT et sur la somme de qualité des pb dans un échantillon. Puis processure itératif : à partir de l'allèle la plus grande, pour tous les alignements qui sont complètement contenus dans l'intervalle on va chercher les haptolypes. On prend le plus grand haplotype la plus grande et l'intervalle est étiré vers la droite pour inclure cet allèle. On réitère jusqu'à ce que la borne à droite ne contienent aucun haplotype passant les critères de qualité
3. On applique le principe précedent mais pour calculer la probabilité a postério à uene position, on fait une recherche de gradient sur tous les génotypes (donc plusieurs échantillons). Marginal probability.
4. Estimation de la qualité du génotypage définie approxivmativement comme la probabilité d'un polymorphismes
5. estimalation de la qualité d'_un_ génotype en somme mais seulement dans un espace restreint
BWA MEM est peusdo-random et si les reads sont mélangés, un read ambigu peut changer d'aligenment !
A noter aussi qu'HaplotypeCaller est le plus sensible à ce mélange (1%)..
17 aligneurs. Comparaison récente d'aligneurs pour Illumina et Ion torrent.
**** Notes
Données
- simulées (génomique):
- humain (GRCh38.p13) et souris
- read 50bp et 150bp
- génération de reads avec erreur aléatoire selon modèle statistique basé sur couverture, erreur de séquencage, distribution mutation, GC
- ART,DWGSIM, WGSIM, MASON, CURIM
- réelles : 4 échantillons
- 1 RNA-seq, 1 genome (rétinine pigmentaire) sur Io torrent
- 1 RNA seq et 1 exome (ostéomyélite) sur illumina hiset 2500 (paired end)
Méthode:
- calcul sur données simulée
- sensibilité = reads mappé correctement/reads mappé incorrect (seuil fixé)
- scores
- MAPQ = 10log_10 proba que la psoition soit fausse
- Alignement Scores = similarité entre query et référence : proportionel au nombre de matches et inv. prop au nombre de mismaps et gaps
- grand AS et petit MAPQ = alignement parfait à plusieurs position.
- petit AS et grand MAPQ = alignement avec discordance mais la positdion est la plus probable
- mapping efficinency = nombre de reads aligné.
- dépend de la longueur des reads, de leur qualité, absence contamination, logiciel, génome de référence +/- librairie
- note: clipping dears reads = les extrémitiés qui ne s'alignent pas sont ignées
- distribution selon les régions GENCODE
- distribution indel
- GC biais. Théoriquement chez les humaine, gaussienne autour de 41%
- duplicats (préparation échantillion ou 1 cluster amplification) -> les supprimer n'implique pas forcément une augmentation de la précision
- 6 coeurs i7 32Gob RAM
Résultat:
- données simulée en paired-end: meilleur = BBMAP, BWA-MEM, Novoalign, DNASTAR, YARA,
Segemehl and TopHat2
- Échantillons :
- accuracy (si chaque read a bien été aligné: comment est-il calculé sur données non simulées TODO ???): novoalign, DNASTAR.
- Meilleur résultat pour exon mapping (nombre d'exons marqué comme tels par kb) = Segemehl
- *Pire résultat pour exon alignement : BWA*. Meilleur = variable selon écdhantillons mais RUM, Segemehl pour WES, novolagin pour WGS
- CLC, BWA-MEM, GEM and Magic-BLAST ont le plus grand nombre de reads alignés
- clipped reads : bowtie2, esegemehld et tophat2 = peu de clipping sur exome/génome
- délétion : distribution augmentée pour génome sur toute la longueur du read pour tous mapper
- GC: génome : aucun aligner n'a une distribution romale !! (biais de séquencage). Exome= double pic pour bowtie2, bwa, gema, minimap2, star, subrdea
- duplicats : faible pour dna
Performances : dépend de la longeur du reads et du nombre de reads
- mémoire : milimu = tophat2, bwa, bwa-mem, hisat2, bowtie2, GEM < 10Gb
-> utiliser plutôt [cite:@alser2021] pour données brutes (figure trop petite, mauvaise méthodo)
Conclusion
- single-end = moins de redas mappé mais diminué nombre d'alignement multiples
- paired-end = plus de reads mappé, meilleur efficacaté, meilleurs reads uniquué mappé (exome) => compromose efficacicé et qualité
- magic-blast = aligne le plus de redas mais peu précise. novoalign et dnastar = meilleur qualité. Segemehl et novoalign pour exon.
- le plus rapide = subread puis minimap2. .segemelhl, tophat2 et novoalign = les lpus lents
* Description
- SNP, indel, "multiple nucleotide polymorphisms"
- annotation sur position (intronic, UTR, upstream, downsstream, site d'épissage, intergénique)
- pour gènes non codant, annotation + biotype si information disponible
- prédiction par rapport aux gènes codant pour les protéines (synonyme ou non, gain/pert d'un codon start ou stop, frameshift)
- différence avec annovar : plus de version du génome, un peu plus rapide
@table
| Warnings | Any warnings or errors. |
| Gene_ID | Gene ID (usually ENSEMBL) |
| Gene_name | Gene name |
| Bio_type | BioType, as reported by ENSEMBL. |
| Trancript_ID | Transcript ID (usually ENSEMBL) |
| Exon_ID | Exon ID (usually ENSEMBL) |
| Exon_Rank | Exon number on a transcript |
| Effect | Effect of this variant. See details below. |
| old_AA/new_AA | Amino acid change |
| old_codon/new_codon | Codon change |
| Codon_Num(CDS) | Codon number in CDS |
| Codon_degenaracy | Codon degenaracy |
| CDS_size | CDS size in bases |
| Custom_interval_ID | If any custom interval was used, add the IDs here (may be more than one). |
@end
- Testé sur 356 000 SNPs entre 2 strans de drosophile
Apport: plusieurs séquencage de NA12878 BGI+ illumina
**** Méthodes
- NA12878 sur plusieurs séquenceur (exome + genome)
| Sequencing Samples | Type | Bases (Gbp) | Read (x10^6) | Clean rare | >Q20 | >Q30 | GC content | Mean coverage |
| BGISEQ500 | WES | 29.41 | 294.30 | 0.41% | 96.72% | 89.14% | 49.75% | 328.49X |
| MGISEQ2000 | WES | 16.34 | 163.55 | 0.25% | 98.18% | 92.08% | 49.71% | 129.40X |
| HiSeq4000 | WES | 41.93 | 283.70 | 4.46% | 97.36% | 93.01% | 50.63% | 395.17X |
| NovaSeq | WES | 25.88 | 178.87 | 2.25% | 95.33% | 92.67% | 49.73% | 241.52X |
| BGISEQ500 | WGS | 126.86 | 1270.02 | 1.76% | 93.73% | 83.33% | 41.76% | 41.03X |
| MGISEQ2000 | WGS | 137.36 | 1374.87 | 0.21% | 96.17% | 88.19% | 41.76% | 45.13X |
| HiSeq4000 | WGS | 191.00 | 1276.10 | 8.25% | 95.90% | 90.11% | 41.69% | 58.00X |
| NovaSeq | WGS | 98.30 | 657.45 | 1.28% | 95.89% | 93.86% | 41.61% | 28.96X |
| HiSeq Xten | WGS | 134.00 | 894.58 | 7.29% | 94.50% | 87.63% | 40.71% | 38.93X |
>Q20 = moins de 1% d'erreur. >Q30 = moins de 1‰
- 1 aligneur (bwa [mem non spécifiqu]) et 3 appel de variant : haplotypcaller, strelka, samtools varscan2
- filtering et trimming
**** Résultat
- qualité :
- exome 95% avec moins de 1% d'erreur et 89% avec moins de 1‰
- genome 92% et 93% respectivement
- profondeur :
| | | exome | genome |
| BGI | bgiseq500 | 328X | 41x |
| BGI | mgiseq2000 | 129 | 45x |
| illumina | hiseq | 395 | 58 |
| | novaseq | 241 | 28.96 (*) |
(*) 92.1% > 20x
- Appel de variant
- Exome:
- comparaison entre pipeliness
- meilleur précision quelque soit le seuil de recall pour exome. Mais en regardant la figure, GATK4 identique ou légèrement meilleur sur indel sur ihseq400 et novaseq (NB: l'article parle d'un "weaker recall" pour SK par rapport à GATK4 ??)
- F-score : divergences importance (0.75-0.91)
- SNP: SK > GAK > sv quelque soit la plateforme.
- indel : SK > autre pour une plateforme donnée. Pour BIG, Sk est meilleur mais GATK ~ SK
- très bonne concordance :
- SNP : > 97.13% des variants détectés par 10 combinaisons (90.92% retrouvé par les 12 combinaisons)
- meilleur score : novaseq + samtools varscan. Pire = hiseq400-SK2
- indel : > 93.2% indel détecté par au moins 10 combinaisons
- idem (meilleur score : novaseq + samtools varscan. Pire = hiseq400-SK2)
- temps d'exécution :
- méthode: downsampling -> 100x puis test avec donwnsampling 20x, 40x, 60X et 80x
- 88GB mémoire, 24 cpus
- résultat: strelka2 plus rapide que GATK (facteur 6-8 selon séquenceur) et que SV (facteur 67)
- Genome
- comparaison entre pipelines
- SK2 = meilleur précision.
- performance similaire exome
- Pour SNP, bonne performance pour toutes les combinaison. SK > GATK pour gbi et illumina (F-score).
- indel : gatk et sk2 ont meilleur performance. Illumina: SK2 > gatk, SV. Idem BGI (f-score)
- Conclusion : SK2 meilleur pour bg/illumina sur exomegenome
- très bonne concordance:
- 94.22% SNV détecté par au moins 10 combinaison
- 90.63% indel détecté par au moins 13 (plus grande divergence, comme pour les exomes)
- temps d'exécution (même méthode et même config): sk2 42-45x plus rapide que GATK et 62-69x que SV (illumina-bgi)
Conclusion: recommande strelka2
Croissance anormale de cellules trophoblastique. Pas de tissu embryonnaire ou foetal.
Pas de contribution maternelle, seulement paternelle => homozyoget.
Ovum vide sans noyau fertilisé par du sperm 23X, qui va se réplicer avec caryotype 46XX (isodisomy paternelle uniparentelle)o
NB: si fertilisation avec 2 seperme, il peut y avoir 46XX ou 46XY. Très rarement tétraploidie.
- GATK = le plus utilisé
- pas le plus performant : deepvariant = precisionFDA en 2016 (préision SNP) et en 2020 Illumina dragen
- quid de la vitesse et précision ?
- apport : pas de données simulée, test dragen et teste l'aliggnement
Méthode:
- alignement : dragen vs GATK avec BWA mem2
- appel de variant : deep variant et dragen contre GATK
- HG002 séqunecé 70x et trio x3
- 64 coeurs par pipeline
Résultat
- dragen plus rapide (x10 alignement)
- chr20-22: meilleur F1 avec dragen en alignemnt (varianlt calling = résultat moins nets)
- chr20-22: idem régions compliqué
- deletion: GATK pas moins inntéressér selon la proportion
Limite:
- chr20-22 seulement !!
- ne reproduit pas les résultat de Hwang et al.15 "who demonstrated an excellent performance of BWA-MEM for mapping/alignment when combined with the GATK Haplotypecaller. "
- dragen = illumina, non open source
La plus récente
- 7 patients de référence (exome + genome) GIIAB: trio NA1287, 91 et 92, trio chinois, trio ashkenazi
- genome + exome pour chaque
- exome : agilent sreselect v5, 100-200x
- Mais l'alignement se fait sur GRCh37...
- différents filtre : 3 pour haplotypecaller, 2 pour octopus (standard, random forest)
Aligneurs testé :
- BWA MEM v.0.7.17
- Bowtie2 v.2.3.5.1 [6],
- Novoalign v. 4.02.01
- Isaac
Variant caller
FreeBayes v. 1.3.1 [35], GATK HaplotypeCaller (HC) v. 4.2.3 [8, 9], Strelka2 v. 2.9.10 [13], DeepVariant v. 1.2.0 [10], Clair3 v. 0.1-r8 [25], and Octopus v. 0.7.4 [26].
**** Résultats
Comparatifs
- régions CDS: appel de variant plus important que aligner (confirmé par test statstique)
- bwamem + deepvariant = meilleur score F1
- deepvariant = meilleur quelque soit aligneur SNP et indel
- GATK bonne performance seulement avec filtre CNNScorevVariant CNN reference-based (autre chois = reads bases et hard filtering mais plus mauvais )
- mettre figure 1.a
- éviter bowtie2 comme aligneur, freebayes comme appel de variant
- précision semblable quelque soit le pipeline mais variabilité recall => non lié aux faux positifs mais aux faux négatifs (par définition...)
- effet filtre sur variant dépend du type de données et aligneur
Facteurs
- exomes ~ genome pour SNP mais plus mauvais performance pour indel
- test sur régions au-delà des CDs (+/- 50bp): diminution faible pour SNP et forte diminution pour indel en précision (meme à partir de 25pb)
- faible couverture = diminution performance mais aussi pour certains appels de variants en cas de couverture très importante ! (GTAK + filtre réseau de neurone)
- GC : pas d'effet très important mais si fortement enrichi/pauvre en GC, sera mal couver par exome. Ici: confirmé pour fortement enrichi en GC
- régions difficile à mapper : algorithmes basé sur haplotype (gatk, freebayes) y sont moins sensible, contrairement à des modèles basés sur le machine learning (deepvariant)
Patients hors GIAB
- métrique: concordance des appels de variants
- ont vérifié que c'était bien le cas : les meilleurs sont GIAB sont bien concordants, les plus mauvais sont les moins concordants. Ont vérifié aussi que FN = manqué par un seul outils mais FP = seulement rapporté par un autils
- 3 individus nigerien (Yoruba) (genome), 3 russes (exome)
- résultats identique
Résultat de GIAB : genome et exome (différent kit de capture)
- gatk4, strelka2, octopus, deepvariant
- comparé au gold truth GAB: HG0001, 5,6,7 sur grch37 (car lifté) et 2,3,4 en 38
- pour exome : en 37, ont fait l'union des capture et des gènes refseq. EN 38, exons refesq +100bp de chahaque côté
- novaseq et hiseq4000 pour exome
***** Résultat
- genome :
- deepvariant meilleur (priciosn + recall).
- indel : choix de la méthode plus important que couverture
- note: les version avec PCR sont moins précise sur les indel (probablement du au homopolymérie)
- exome :
- précision diminué, cohérent avec études sanger[fn:29]
- exons refseq : recall diminué (problème de couverture), précision aussi (régions difficiles à capturée)
- deepvariant plus précis (soit à couvertures plus faible, soit car entraîné sur les erreur des séquencage)
- choix de la capture: aligent: meilleur recall, IDT = meilleur précision[fn:30]
- novaseq plus précs hiseq4000
- précision diminuée pour indel
Tutorial + présentation des outils. Référence recommandées par GATK
*** DONE Best practice
CLOSED: [2023-11-22 Wed 23:03]
**** MarkDuplicates
Duplicats =
- soit lié à l'amplification PCR lors de la préparation (si elle est utilisées)
- version récente : amplification utiles donc les duplicats vont venir de problme optiques (1 cluster lu comme 2 sur illumina donc traité comme 2 reads)
- stratégie
- cf figure 6.2 : indel ou reads plus long = ne correspondent pas mais softclip oui
- groupe reads et pour chaque groupe marque le read qui a la meilleure qualité de séquencage. Les autres reads sont marqué comme "à ignorer"
NB: étape coûteuse. version parallélisée avec MarkDuplicatesSpark (même sur architecture non spark) mais nécessite de nombreuses écriture/lectures. Voir
https://gatk.broadinstitute.org/hc/en-us/articles/360037224932-MarkDuplicatesSpark
**** Recalibration base qualité score
- Détection et correction 'erreurs systématiques dans les score de qualité de paires de base (données par séquences).
- Biais possibles : durant préparation de libraire, séquencage, défaut de la flowcell ou du séquenceur
- 1. Récupération de toutes les statistiques "covariées" des bp
- 2. Détermination d'un modèle statistique: voir aussi https://gatk.broadinstitute.org/hc/en-us/articles/360035890531
- ! dépend de la librairie et du séquenceur donc à déterminer à chaque fois...A
- nécessite des positions de SNPs connus pour les exclure (typiquement dbSNP)
- calculé à partir du group de reads, quality score, cycle machine et base + base précédente
- 3. Recalibration
*** Haplotypecaller
- Continuité de A. MkKusick Mendelian INheritance in Man
- publié en ligne 1987
- source = veille scientifique journaux biomédicaux peer-reviewed (45) + pubmed
- objectifs: gene-phenotype avec 2 types d'entrée : gene et phénotype
- chaque entrée est formattée de la même manière + liexns experten (genome(ensemble), clinique (gene reviews), variant (clinvar)...)
- identifiants OMIM sont largement utilisé dans la communautié + littérature (6 chiffres: le epremier est 1,2,6 si autosomoqiue, 3 si lié à l'X, 4 si lié à l'Y et 5 si mtochondrial. Astérisque devant el chiffre si gène)
- gène : éléments régulateurs, gène codnant et autres éléments fonctionnels
- phénotype(préfixe par # si base moléculaire connue, % si mendelien mais sans base connuée ou rien): maladies monogénique, phénotypes, sensibilité à réaction médicaments, sensibiilté à infection, sensibilité constit à un cancer, syndrome del et duplication
- gene map et morbid map : seulement les données OMIM avec emplacement cytogénétique publié
- critère pour
- inclusion mutation gene-phénotype : 1ere mutation découverte, fréquence dans la population, phénotype distinctif, signification historique, mécanisme de mutation/pathogénétique inhabituation ou hérité de manière atypique. Plupart causent des maladies mais il y a des SNP
- relation gène phénotype
- variant chez plusieurs individus non apparenté
- ségrège avec le phénotype dans plusieurs familles
- de novo dans un nombre sigificatif de famille
- "qualifié" (qualified) si 1 famille avec plusieurs variants dans 1 gène et les variants ségrègent avec le phénotype + peuve fonctionnelle (NB: 1 seul patient peut suffir si l'argumentaire est solide)
#+caption: OMIM Stats (19 décembre 2023)
| MIM Number Prefix | Autosomal | X Linked | Y Linked | Mitochondrial | Totals |
| Gene description * | 16,309 | 769 | 51 | 37 | 17,166 |
| Gene and phenotype, combined + | 21 | 0 | 0 | 0 | 21 |
| Phenotype description, molecular basis known # | 6,344 | 383 | 5 | 34 | 6,766 |
| Phenotype description or locus, molecular basis unknown % | 1,390 | 111 | 4 | 0 | 1,505 |
| Other, mainly phenotypes with suspected mendelian basis | 1,640 | 100 | 3 | 0 | 1,743 |
| Totals | 25,704 | 1,363 | 63 | 71 | 27,201 |
Revue des algorithmes
**** Notes
Introduction
- Brute-force = non utilisable sur machine modernes (centains millions de read)
- Méthode actuelle
- Alignement = étape cruciale
- 2 possibilité pour un génome : assemblage /de novo/ ou sur un génomne de rééfrence. La première solution n'est pas encore utilisée (génome complexe et reads petite taille, problème non résolu par long read) car trop lente et plus coûteuse
Problèmes
- génome de réference incomplet (reads non aligné/incorrectement)
- régions répétées (1 reads peut aligner sur différente région) Ex: BWA mem qui aligne "au hasard" ([cite:@firtina2016]: si les reads sont mélangés, un read ambigu est aligné à un endroit différent)
- Doit tolérer SNV ou variants structuraux
- Doit gérer les brins forward et reverse pour éviter d’avoir des génotypes différents (strand bias)
***** Liste complète d'aligner jusque 2020: Table 1
Aligneur gérant l'ADN après 2012
(W) = wrapper (utilise d'autres outils)
****** Initiale
| Aligner | Paper | publication | Application | Indexing | GP fix | GP Spaced speed | GP seed chaining | Pairwise alignment | Max. read length |
| | | | | | | | | | tested (bp) |
|----------------+-------+-------------+----------------+--------------------+--------+-----------------+------------------+-------------------------+------------------|
| [[https://github.com/lh3/bwa][BWA]] | [54] | 2009 | DNA | BWT-FM | N | N | N | Semi-Global | 125 |
| [[http://bowtie-bio.sourceforge.net/manual.shtml][Bowtie]] | [55] | 2009 | DNA | BWT-FM | Y | N | N | HD | 76 |
| [[https://sourceforge.net/projects/cloudburst-bio/][CloudBurst]] | [56] | 2009 | DNA | Hashing | Y | N | N | Landau-Vishkin | 36 |
| [[https://github.com/byucsl/gnumap][GNUMAP]] | [57] | 2009 | DNA | Hashing | Y | N | Y | NW | 36 |
| [[http://1001genomes.org/software/genomemapper_singleref.html][GenomeMapper]] | [58] | 2009 | DNA | Hashing | Y | N | Y | NW | 200 |
| [[https://github.com/hugheaves/MOM][MOM]] | [59] | 2009 | DNA | Hashing | Y | N | N | HD | 40 |
| [[http://pass.cribi.unipd.it/cgi-bin/pass.pl][PASS]] | [60] | 2009 | DNA | Hashing | Y | N | Y | NW | 32 |
| [[https://code.google.com/archive/p/perm/downloads][PerM]] | [61] | 2009 | DNA | Hashing | Y | Y | N | HD | 47 |
| [[https://github.com/seqan/seqan/tree/master/apps/razers][RazerS]] | [62] | 2009 | DNA | Hashing | Y | Y | Y | Myers Bit Vector | 76 |
| [[http://compbio.cs.toronto.edu/shrimp/][SHRiMP]] | [63] | 2009 | DNA | Hashing | N | N | N | SW | 35 |
| [[https://github.com/ShujiaHuang/SOAPaligner][SOAP2]] | [64] | 2009 | DNA | BWT-FM | Y | N | N | SW | 44 |
| [[http://www.bcgsc.ca/platform/bioinfo/software/slider][Slider]] | [65] | 2009 | DNA | Hashing | Y | N | N | HD | 36 |
| [[https://www.bioinf.uni-leipzig.de/Software/segemehl/][segemehl]] | [66] | 2009 | DNA | Suffix array | N | N | Y | SW | 35 |
| [[https://github.com/lh3/bwa][BWA-SW]] | [54] | 2010 | DNA | BWT-FM | N | N | N | SW | 10000 |
| [[http://www.irisa.fr/symbiose/projects/gassst/][GASSST]] | [35] | 2010 | DNA | Hashing | Y | Y | Y | Semi-Global | 500 |
| [[https://github.com/juliangehring/GMAP-GSNAP][GSNAP]] | [37] | 2010 | DNA | Hashing | Y | N | Y | Non-DP Heuristic | 100 |
| [[https://github.com/rcallahan/smalt][SMALT]] | [69] | 2010 | DNA | Hashing | Y | N | Y | SW | 150 |
| [[http://www.bcgsc.ca/platform/bioinfo/software/SliderII][SliderII]] (W) | [70] | 2010 | DNA | Hashing | Y | N | N | HD | 42 |
| [[http://www.vmatch.de/][VMATCH]] (W) | [71] | 2010 | DNA | Suffix array | Y | N | Y | SW | N/A |
| [[https://github.com/sfu-compbio/mrsfast][mrsFAST]] | [72] | 2010 | DNA | Hashing | Y | N | N | HD | 100 |
| [[http://last.cbrc.jp/][LAST]] | [78] | 2011 | DNA/BS-Seq/RNA | Suffix array | N | Y | N | SW & NW | 105 |
| [[https://dl.acm.org/citation.cfm?id=2147845&dl=ACM&coll=DL][DynMap]] | [79] | 2011 | DNA | Hashing | Y | N | N | NW | 52 |
| [[http://compbio.cs.toronto.edu/shrimp/][SHRiMP2]] | [80] | 2011 | DNA | Hashing | Y | Y | Y | SW | 75 |
| [[http://snap.cs.berkeley.edu/][SNAP]] | [81] | 2011 | DNA | Hashing | Y | N | N | NW | 10000 |
| [[https://www.well.ox.ac.uk/project-stampy][Stampy]] | [82] | 2011 | DNA | Hashing | Y | N | N | NW | 4500 |
| [[https://github.com/iontorrent/TS/tree/master/Analysis/TMAP][TMAP]] | ? | 2011 | DNA | BWT-FM | N | N | Y | SW | N/A |
| [[http://grimmond.imb.uq.edu.au/X-MATE/][X-Mate]] | [83] | 2011 | DNA | Hashing | N | N | N | Non-DP Heuristic | 50 |
| [[https://github.com/mchaisso/blasr/][BLASR]] | [85] | 2012 | DNA | Suffix array | Y | N | Y | NW | 8000 |
| [[https://code.google.com/archive/p/batmis/][Batmis]] | [86] | 2012 | DNA | BWT-ST | Y | N | N | HD | 100 |
| [[http://bowtie-bio.sourceforge.net/bowtie2][Bowtie2]] | [87] | 2012 | DNA | BWT-FM | Y | N | Y | SW & NW | 400 |
| [[https://github.com/smarco/gem3-mapper][GEM]] | [88] | 2012 | DNA | BWT-FM | N | N | Y | SW & NW | 150 |
| [[https://github.com/seqan/seqan/tree/master/apps/razers3][RazerS3]] | [89] | 2012 | DNA | Hashing | Y | Y | Y | Banded Myers Bit Vector | 800 |
| [[https://web.stanford.edu/group/wonglab/seqalto/][SeqAlto]] | [90] | 2012 | DNA | Hashing | Y | N | N | NW | 200 |
| [[https://github.com/seqan/seqan/blob/master/apps/splazers/README][SplazerS]] | [91] | 2012 | DNA | Hashing | Y | N | Y | Banded Myers Bit Vector | 150 |
| [[http://pages.cs.wisc.edu/~jignesh/wham/][WHAM]] | [92] | 2012 | DNA | Hashing | Y | N | N | NW | 74 |
| [[https://github.com/GregoryFaust/yaha][YAHA]] | [93] | 2012 | DNA | Hashing | Y | N | Y | SW | 10000 |
| [[http://subread.sourceforge.net/][Subread]] | [97] | 2013 | DNA/RNA-Seq | Hashing | Y | Y | Y | SW | 202 |
| [[https://github.com/lh3/bwa][BWA-MEM]] | [98] | 2013 | DNA | BWT-FM | N | N | Y | SW & NW | 650 |
| [[http://www.seqan.de/projects/masai][Masai]] | [99] | 2013 | DNA | Suffix tree | N | N | Y | Banded Myers Bit Vector | 150 |
| [[http://cibiv.github.io/NextGenMap/][NextGenMap]] | [100] | 2013 | DNA | Hashing | Y | N | N | SW & NW | 250 |
| [[http://www.umsl.edu/~wongch/software.html][SRmapper]] | [101] | 2013 | DNA | Hashing | Y | N | N | HD | 100 |
| [[https://github.com/BilkentCompGen/mrfast][mrFAST]] | [102] | 2013 | DNA | Hashing | Y | N | N | Semi-Global | 180 |
| [[http://bwa-pssm.binf.ku.dk/][BWA-PSSM]] (W) | [107] | 2014 | DNA | BWT-FM | Y | N | N | SW | 100 |
| [[http://cushaw3.sourceforge.net/homepage.htm#latest][CUSHAW3]] | [108] | 2014 | DNA | BWT-FM | Y | N | Y | SW & Semi-Global | 100 |
| [[https://hobbes.ics.uci.edu/download.shtml][Hobbes2]] | [109] | 2014 | DNA | Hashing | Y | N | Y | Banded Myers Bit Vector | |
| [[https://github.com/wanpinglee/MOSAIK][MOSAIK]] | [110] | 2014 | DNA | Hashing | Y | N | N | SW | 100 |
| [[https://github.com/opencb/hpg-aligner][hpg-Aligner]] | [111] | 2014 | DNA | Suffix array | N | N | Y | SW | 5000 |
| [[https://github.com/sfu-compbio/mrsfast][mrsFAST-Ultra]] | [112] | 2014 | DNA | Hashing | Y | N | N | HD | 100 |
| [[http://erne.sourceforge.net/][ERNE2]] | [116] | 2016 | DNA/BS-Seq | BWT-FM +hashing | Y | N | N | HD | 100 |
| [[https://github.com/isovic/graphmap][GraphMap]] | [117] | 2016 | DNA | Hashing | Y | Y | Y | Semi-global | 9000 |
| [[https://github.com/ruhulsbu/NanoBLASTer][NanoBLASTer]] | [118] | 2016 | DNA | Hashing | Y | N | Y | NW | 7040 |
| [[https://github.com/lh3/minimap][minimap]] | [119] | 2016 | DNA | Hashing | Y | N | N | N/A | 13000 |
| [[https://github.com/dfguan/rHAT][rHAT]] | [120] | 2016 | DNA | Hashing | Y | N | Y | SW | 8000 |
| [[https://github.com/hsinnan75/KART][KART]] | [121] | 2017 | DNA | BWT-FM | N | N | Y | NW | 7118 |
| [[https://github.com/hitbc/LAMSA][LAMSA]] (W) | [122] | 2017 | DNA | BWT-FM + hashing | Y | N | Y | Sparse DP | 1000 |
| [[https://github.com/lh3/minimap2][minimap2]] | [124] | 2018 | DNA/RNA-Seq | Hashing | Y | N | Y | NW | 11628 |
| [[https://gitlab.com/pirovc/dream_yara/][DREAM-Yara]] (W) | [125] | 2018 | DNA | BWT-FM | Y | N | N | Banded Myers Bit Vector | |
| [[https://github.com/mummer4/mummer][MUMmer4]] (W) | [126] | 2018 | DNA | Suffix array | Y | N | Y | SW | 7821 |
| [[https://github.com/philres/ngmlr][NGMLR]] | [127] | 2018 | DNA | Hashing | Y | N | Y | SW | 50000 |
| [[https://github.com/vpc-ccg/lordfast][lordFAST]] | [128] | 2018 | DNA | BWT-FM + hashing | N | N | Y | SW & NW | 35489 |
| [[https://github.com/lbcb-sci/graphmap2][GraphMap2]] | [130] | 2019 | DNA/RNA-Seq | Hashing | Y | Y | Y | Semi-global | 9000 |
| [[https://github.com/ncbi/magicblast][Magic-BLAST]] | [131] | 2019 | DNA/RNA-Seq | Hashing | Y | N | N | Non-DP Heuristic | 90000 |
| [[https://github.com/bwa-mem2/bwa-mem2][BWA-MEM2]] | [132] | 2019 | DNA | BWT-FM | N | N | Y | SW | 650 |
| [[https://ccb.jhu.edu/software/hisat2/index.shtml][HISAT2]] | [133] | 2019 | DNA | BWT-FM | Y | N | N | Non-DP Heuristic | 100 |
| [[https://www.dropbox.com/s/3jcu4i240kyu2tc/source%20code%20conLSH_bio.tar.gz?dl=0][conLSH]] | [135] | 2020 | DNA | Hashing | Y | N | Y | Sparse DP | 8000 |
Étapes
1. Index du génome
2. Liste endroits possible avec cet index
3. Pour chaque position possible, calcul de similarité avec le read
****** Avec citations
#+begin_src julia
using DataFramesMeta, CSV
info = CSV.read("biblio-aligner-init.tsv", DataFrame, delim="\t")
cite = CSV.read("citation-aligner-init.tsv", DataFrame, delim="\t")
d = innerjoin(info, cite, on = :ID)
CSV.write("biblio-aligner.tsv", d, delim="\t")
#+end_src
| Link | Aligner | ID | publication | Application | Indexing | GP fix | GP spaced speed | GP seed chaining | Pairwise alignment | Max. read length | Doi | Pubmed |
| http://www.irisa.fr/symbiose/projects/gassst/ | GASSST | 35 | 2010 | DNA | Hashing | Y | Y | Y | Semi-Global | 500 | 10.1093/bioinformatics/btq485 | 20739310 |
| https://github.com/juliangehring/GMAP-GSNAP | GSNAP | 37 | 2010 | DNA | Hashing | Y | N | Y | Non-DP Heuristic | 100 | 10.1093/bioinformatics/btq057 | 20147302 |
| https://github.com/lh3/bwa | BWA | 54 | 2009 | DNA | BWT-FM | N | N | N | Semi-Global | 125 | 10.1093/bioinformatics/btp698 | 20080505 |
| https://github.com/lh3/bwa | BWA-SW | 54 | 2010 | DNA | BWT-FM | N | N | N | SW | 10000 | 10.1093/bioinformatics/btp698 | 20080505 |
| http://bowtie-bio.sourceforge.net/manual.shtml | Bowtie | 55 | 2009 | DNA | BWT-FM | Y | N | N | HD | 76 | 10.1186/gb-2009-10-3-r25 | 19261174 |
| https://sourceforge.net/projects/cloudburst-bio/ | CloudBurst | 56 | 2009 | DNA | Hashing | Y | N | N | Landau-Vishkin | 36 | 10.1093/bioinformatics/btp236 | 19357099 |
| https://github.com/byucsl/gnumap | GNUMAP | 57 | 2009 | DNA | Hashing | Y | N | Y | NW | 36 | 10.1093/bioinformatics/btp614 | 19861355 |
| http://1001genomes.org/software/genomemapper_singleref.html | GenomeMapper | 58 | 2009 | DNA | Hashing | Y | N | Y | NW | 200 | 10.1186/gb-2009-10-9-r98 | 19761611 |
| https://github.com/hugheaves/MOM | MOM | 59 | 2009 | DNA | Hashing | Y | N | N | HD | 40 | 10.1093/bioinformatics/btp092 | 19228804 |
| http://pass.cribi.unipd.it/cgi-bin/pass.pl | PASS | 60 | 2009 | DNA | Hashing | Y | N | Y | NW | 32 | 10.1093/bioinformatics/btp087 | 19218350 |
| https://code.google.com/archive/p/perm/downloads | PerM | 61 | 2009 | DNA | Hashing | Y | Y | N | HD | 47 | 10.1093/bioinformatics/btp486 | 19675096 |
| https://github.com/seqan/seqan/tree/master/apps/razers | RazerS | 62 | 2009 | DNA | Hashing | Y | Y | Y | Myers Bit Vector | 76 | 10.1101/gr.088823.108 | 19592482 |
| http://compbio.cs.toronto.edu/shrimp/ | SHRiMP | 63 | 2009 | DNA | Hashing | N | N | N | SW | 35 | 10.1371/journal.pcbi.1000386 | 19461883 |
| https://github.com/ShujiaHuang/SOAPaligner | SOAP2 | 64 | 2009 | DNA | BWT-FM | Y | N | N | SW | 44 | 10.1093/bioinformatics/btp336 | 19497933 |
| http://www.bcgsc.ca/platform/bioinfo/software/slider | Slider | 65 | 2009 | DNA | Hashing | Y | N | N | HD | 36 | 10.1093/bioinformatics/btn565 | 18974170 |
| https://www.bioinf.uni-leipzig.de/Software/segemehl/ | segemehl | 66 | 2009 | DNA | Suffix array | N | N | Y | SW | 35 | 10.1371/journal.pcbi.1000502 | 19750212 |
| https://github.com/rcallahan/smalt | SMALT | 69 | 2010 | DNA | Hashing | Y | N | Y | SW | 150 | | |
| http://www.bcgsc.ca/platform/bioinfo/software/SliderII | SliderII (W) | 70 | 2010 | DNA | Hashing | Y | N | N | HD | 42 | 10.1093/bioinformatics/btq092 | 20190250 |
| https://github.com/sfu-compbio/mrsfast | mrsFAST | 72 | 2010 | DNA | Hashing | Y | N | N | HD | 100 | 10.1038/nmeth0810-576 | 20676076 |
| http://last.cbrc.jp/ | LAST | 78 | 2011 | DNA/BS-Seq/RNA | Suffix array | N | Y | N | SW & NW | 105 | 10.1101/gr.113985.110 | 21209072 |
| https://dl.acm.org/citation.cfm?id=2147845&dl=ACM&coll=DL | DynMap | 79 | 2011 | DNA | Hashing | Y | N | N | NW | 52 | 10.1145/2147805.2147845 | |
| http://compbio.cs.toronto.edu/shrimp/ | SHRiMP2 | 80 | 2011 | DNA | Hashing | Y | Y | Y | SW | 75 | 10.1093/bioinformatics/btr046 | 21278192 |
| http://snap.cs.berkeley.edu/ | SNAP | 81 | 2011 | DNA | Hashing | Y | N | N | NW | 10000 | | |
| https://www.well.ox.ac.uk/project-stampy | Stampy | 82 | 2011 | DNA | Hashing | Y | N | N | NW | 4500 | 10.1101/gr.111120.110 | 20980556 |
| http://grimmond.imb.uq.edu.au/X-MATE/ | X-Mate | 83 | 2011 | DNA | Hashing | N | N | N | Non-DP Heuristic | 50 | 10.1093/bioinformatics/btq698 | 21216778 |
| https://github.com/mchaisso/blasr/ | BLASR | 85 | 2012 | DNA | Suffix array | Y | N | Y | NW | 8000 | 10.1186/1471-2105-13-238 | |
| https://code.google.com/archive/p/batmis/ | Batmis | 86 | 2012 | DNA | BWT-ST | Y | N | N | HD | 100 | 10.1093/bioinformatics/bts339 | 22689389 |
| http://bowtie-bio.sourceforge.net/bowtie2 | Bowtie2 | 87 | 2012 | DNA | BWT-FM | Y | N | Y | SW & NW | 400 | 10.1038/nmeth.1923 | 22388286 |
| https://github.com/smarco/gem3-mapper | GEM | 88 | 2012 | DNA | BWT-FM | N | N | Y | SW & NW | 150 | 10.1038/nmeth.2221 | 23103880 |
| https://github.com/seqan/seqan/tree/master/apps/razers3 | RazerS3 | 89 | 2012 | DNA | Hashing | Y | Y | Y | Banded Myers Bit Vector | 800 | 10.1093/bioinformatics/bts505 | 22923295 |
| https://web.stanford.edu/group/wonglab/seqalto/ | SeqAlto | 90 | 2012 | DNA | Hashing | Y | N | N | NW | 200 | 10.1093/bioinformatics/bts450 | 22811546 |
| https://github.com/seqan/seqan/blob/master/apps/splazers/README | SplazerS | 91 | 2012 | DNA | Hashing | Y | N | Y | Banded Myers Bit Vector | 150 | 10.1093/bioinformatics/bts019 | 22238266 |
| http://pages.cs.wisc.edu/~jignesh/wham/ | WHAM | 92 | 2012 | DNA | Hashing | Y | N | N | NW | 74 | 10.1145/1989323.1989370 | |
| https://github.com/GregoryFaust/yaha | YAHA | 93 | 2012 | DNA | Hashing | Y | N | Y | SW | 10000 | 10.1093/bioinformatics/bts456 | 22829624 |
| http://subread.sourceforge.net/ | Subread | 97 | 2013 | DNA/RNA-Seq | Hashing | Y | Y | Y | SW | 202 | 10.1093/nar/gkt214 | 23558742 |
| https://github.com/lh3/bwa | BWA-MEM | 98 | 2013 | DNA | BWT-FM | N | N | Y | SW & NW | 650 | | |
| http://www.seqan.de/projects/masai | Masai | 99 | 2013 | DNA | Suffix tree | N | N | Y | Banded Myers Bit Vector | 150 | 10.1093/nar/gkt005 | 23358824 |
| http://cibiv.github.io/NextGenMap/ | NextGenMap | 100 | 2013 | DNA | Hashing | Y | N | N | SW & NW | 250 | 10.1093/bioinformatics/btt468 | 23975764 |
| http://www.umsl.edu/~wongch/software.html | SRmapper | 101 | 2013 | DNA | Hashing | Y | N | N | HD | 100 | 10.1093/bioinformatics/bts712 | 23267171 |
| https://github.com/BilkentCompGen/mrfast | mrFAST | 102 | 2013 | DNA | Hashing | Y | N | N | Semi-Global | 180 | 10.1038/ng.437 | 19718026 |
| http://bwa-pssm.binf.ku.dk/ | BWA-PSSM (W) | 107 | 2014 | DNA | BWT-FM | Y | N | N | SW | 100 | 10.1186/1471-2105-15-100 | |
| http://cushaw3.sourceforge.net/homepage.htm#latest | CUSHAW3 | 108 | 2014 | DNA | BWT-FM | Y | N | Y | SW & Semi-Global | 100 | 10.1371/journal.pone.0086869 | 24466273 |
| https://hobbes.ics.uci.edu/download.shtml | Hobbes2 | 109 | 2014 | DNA | Hashing | Y | N | Y | Banded Myers Bit Vector | | 10.1186/1471-2105-15-42 | |
| https://github.com/wanpinglee/MOSAIK | MOSAIK | 110 | 2014 | DNA | Hashing | Y | N | N | SW | 100 | 10.1371/journal.pone.0090581 | 24599324 |
| https://github.com/opencb/hpg-aligner | hpg-Aligner | 111 | 2014 | DNA | Suffix array | N | N | Y | SW | 5000 | 10.1093/bioinformatics/btu553 | 25143289 |
| https://github.com/sfu-compbio/mrsfast | mrsFAST-Ultra | 112 | 2014 | DNA | Hashing | Y | N | N | HD | 100 | 10.1093/nar/gku370 | 24810850 |
| http://erne.sourceforge.net/ | ERNE2 | 116 | 2016 | DNA/BS-Seq | BWT-FM +hashing | Y | N | N | HD | 100 | 10.1186/s12859-016-0910-3 | |
| https://github.com/isovic/graphmap | GraphMap | 117 | 2016 | DNA | Hashing | Y | Y | Y | Semi-global | 9000 | 10.1038/ncomms11307 | 27079541 |
| https://github.com/ruhulsbu/NanoBLASTer | NanoBLASTer | 118 | 2016 | DNA | Hashing | Y | N | Y | NW | 7040 | 10.1109/iccabs.2016.7802776 | |
| https://github.com/lh3/minimap | minimap | 119 | 2016 | DNA | Hashing | Y | N | N | N/A | 13000 | 10.1093/bioinformatics/btw152 | 27153593 |
| https://github.com/dfguan/rHAT | rHAT | 120 | 2016 | DNA | Hashing | Y | N | Y | SW | 8000 | 10.1093/bioinformatics/btv662 | 26568628 |
| https://github.com/hsinnan75/KART | KART | 121 | 2017 | DNA | BWT-FM | N | N | Y | NW | 7118 | 10.1093/bioinformatics/btx189 | 28379292 |
| https://github.com/hitbc/LAMSA | LAMSA (W) | 122 | 2017 | DNA | BWT-FM + hashing | Y | N | Y | Sparse DP | 1000 | 10.1093/bioinformatics/btw594 | 27667793 |
| https://github.com/lh3/minimap2 | minimap2 | 124 | 2018 | DNA/RNA-Seq | Hashing | Y | N | Y | NW | 11628 | 10.1093/bioinformatics/bty191 | 29750242 |
| https://gitlab.com/pirovc/dream_yara/ | DREAM-Yara (W) | 125 | 2018 | DNA | BWT-FM | Y | N | N | Banded Myers Bit Vector | | 10.1093/bioinformatics/bty567 | 30423080 |
| https://github.com/mummer4/mummer | MUMmer4 (W) | 126 | 2018 | DNA | Suffix array | Y | N | Y | SW | 7821 | 10.1371/journal.pcbi.1005944 | 29373581 |
| https://github.com/philres/ngmlr | NGMLR | 127 | 2018 | DNA | Hashing | Y | N | Y | SW | 50000 | 10.1038/s41592-018-0001-7 | 29713083 |
| https://github.com/vpc-ccg/lordfast | lordFAST | 128 | 2018 | DNA | BWT-FM + hashing | N | N | Y | SW & NW | 35489 | 10.1093/bioinformatics/bty544 | 30561550 |
| https://github.com/bwa-mem2/bwa-mem2 | BWA-MEM2 | 132 | 2019 | DNA | BWT-FM | N | N | Y | SW | 650 | 10.1109/ipdps.2019.00041 | |
| https://ccb.jhu.edu/software/hisat2/index.shtml | HISAT2 | 133 | 2019 | DNA | BWT-FM | Y | N | N | Non-DP Heuristic | 100 | 10.1038/s41587-019-0201-4 | 31375807 |
| https://www.dropbox.com/s/3jcu4i240kyu2tc/source%20code%20conLSH_bio.tar.gz?dl=0 | conLSH | 135 | 2020 | DNA | Hashing | Y | N | Y | Sparse DP | 8000 | 10.1016/j.compbiolchem.2020.107206 | 32000034 |
| | | | | | | | | | | | | |
***** Index
****** Par hash
Pour des séquences courtes (seed), stocke la liste des positions (k-mer). Pour un read, on extrait certains seed et donc la liste de positions.
****** suffix tree
Permet de faire des correspondances partielles en groupant les sous-séquences communes. Cf figure 1.b. BWT utilise une approche similaire mais en diminuant le stockage.
Diminution des performances avec l s erreurs de sequencage ou avec la dissimilarité avec la référence.
****** Comparaison
Hash = index de grande taille mais recherche efficace (O(1) !) et coût d'indexage faible.
Suffix = recherche partielle possible mais recherche lente et coût d'indexage élevé (cf Table 2)
****** Performance
- Run sequential genome
- Pas de différence stastiquement significative entre index et hash pour CPU mais significative pour mémoire (loique)
- BWA-FMT = 3.8 moins de ressources
- BWA, bowtie et bowtie2 proches pour runtime + coût mémoire (cf figure)
Popularité : mesure par le nombre de citation de l'article initial : BLAST > BWA,bowtie,bowtie2
****** Endroits possibles
On prend quelques seeds et on génère une liste de possibilités (voisinage de chaque seed). Pour des seeds courts, liste trop grande donc heuristiques pour réduire le nombre de candidats. Augmenter la taille des seeds diminue ma sensibilité
Majorité utilisent des seeds de taille fixe (autre utilise des suffix car le hash nécessiterait d’être recalculé)
****** Vérification
soit dynamique (local ou global) soit non dynamique (distance de Hamming… probablement comparaison indépendamment du contexte). Pour substitution/indel, dynamic. Alignement local quand on ne veut qu’une correspondance partielle (structural variant) -> 38% outil.
Comparaison quadratique en temps et en espace malgré années recherche …
****** Note :
Long-read
indexing : les 2 méthodes existent
Peuvent couper les long read en short read
Difficulté : plus d’erreur et plus de seed. Solution : calculer moins de seed mais plus représentatives du read
Comparaison avec short read
plus d’erreur
Moins de reads (throughput)
Moins d’ambiguïté car read plus longs
Comparaison : plus facile avec short read car moins d’erreur
SNP plus facile à détecter avec short read (moins d’erreur) mais variants structurels avec long read
***** Autres
RNAseq : problème = aligner des reads sur des zones non contigùes (à cause de l'épissage)
***** Performances
Meilleure que [cite:@donato2021]
10 génomes
- Exome !! haéplotyypecaller vs dragen
- Comparaison en hg19, novaseq6000
- n'utilise pas happy mais un outil illumina
- bwa-mem2 + dragen gatk comparé à bwa mem + haplotypecaller
- Performance variable selon échantillon, mais nouvelle version domine pour l'un des échantillions
NB: 229Mb sont "non-syntenique" = ne s'alignent pas de manière linéaire en GRCh38 sur un intervalle de 1Mbp (ex: inversion)
**** Résultats
- Héritage
- GRCh38: 1 indivu représente 72:6% du génome (56% africain, 28.1% Européan), le second plus gros contributeur est 5.5% majoritairement asie east
- T2T : principalement européenne
- les 2 ont des "introgression" de Néanderthal
- problème de GRCh38 : aux bornes des clones utilisés, il y a des haplotypes de structures anormales -> représente des haplotypes rare pour un individu mais à une fréquence non représentative de la population:
- discordance GRCh38 et 1000 genomes et beaucoup plus rares en T2T
- et ce sont bien aux bornes des clones BAC
- correction false duplications
- méthode : examen de "cluster" de variants hétérozygotes (CHM est homozygote...) en alignant sur GRCh38
- ces zones sont associées à des sgemental duplications, centromère ou problèmes dans gRCh38
- elles contiennent également des variants marqué par gnomad comme "beaucoup trop hétérozygote" (23%)
- après liftover T2T: 48 gènes codant (dont 14 complement contenu) dont DUSP22 (régulation immune), KMT2C (syndrome d Kleefstrah)
- identifié zones marqué comme duplications mais non T2T (22 gènes codants)
- liftover clinvar : 99.8% soit > 800 000 variants (dont 99.6% des patho/probablement patho). Ceux qui ont échoué sont du à des indel différents entre GRCh38 et T2T ou bien à un liftover difficile
- test sur 1000 genomes
- aligner (bwa mem) 0.97% en plus avec taux d'erreur diminué et couverture plus uniforme et charactéristiques. NB: africains = taux de mismatch le plus élevé
- appel de variant (haplytype caller):
- moins de variants par échantillions, attribué à une diminution des allèles rares, erreurs de consense et structurels (surtout pour non-africains peut-ête parce que 70% de gRCh38 vient d'un individu d'ascendace afro-européeenne et les africains ont plus de variants rares)
- ne proposent pas de faire l'appel de variant en T2T et lifter en 38 (si l'allèle de référence n'est pas dans l'échantilon par exemple)
- analyse statitisque : moins de variants de mauvaise qualité, moins de variant discordants au niveau mendelien (chez les enfants mais pas chez les parents, homozygote chez parent mais pas chez les enfants)
- variation de la distribution des AF : surtout rare < 0.05, intermediate et fixé ou quasi fixé. Ces derniers sont du à un variant "privé" d'un donneur ou d'une erreur de GRCh38 donc normal qu'il y en ait moins en T2T. Ceux ou tous les donneur sont hétérozygote : surtout du à des correction dans les segmental duplication fusionné (enrichi en paralag spécific variant hétérozygote). Ces 2 catégories sont la majeure partie
- long read: améliore mapping, balance indel apparent et aide de novo SV ou SV dans séquence non résolue
- régions non résolues précédemment
- région no synthenize : 73-78% des SNVs avec illumina sont concordants avec long-read
- tests de ces régions sur le trio de GIAB et le trio du personal genome project en comparant illumian avec pacbio hifi
- : recall semblable et bas (21-28, 21-25) mais restreint aux duplications fusonnées de 38: 98-99 vs 64-67).
- correction des faux duplicats : 1% -> 57-68% en T2T (recall) et 76-95 -> 98-99 pour précision
- impact clinique
- il existe des variants "délétère" dans le génome de référence qui peut donc altérer l'interprétation d'un variant avec des efforts pour les suppriméer (NB: les citations n'ont pas l'air de le dire clairement...)
- test : alignement avec dipcall:: 210 variants potentiellement perte de fonction[fn:6] identifiés sur 31 gènes d'intérêt clinique[fn:7]. 158 sur 1 individu et la plupart AF=0.47 donc probablement bien toléré. Autres = indel plus grand ou allèle rares. Les 10 varintas sur gènes d'intérêt clinique diréves de paralogue dupliqué
- test sur 4964 gènes d'intérêt clinique ([cite:@Wagner_2022])
- 28 gènes sur régions non résolue /non synthénique
- 756 gènes touché par allèles rares ou structurelles fauesse (306 T2T) dont plupart semble corrigé en T2T
- exemple TNNT3 (arthrogrypose): remaniement complexe supposé (GRC a détemriné qu'il y a un problème) -> T2T: région upstream transposé "inversement" proximal du gène, ce qui le place près de TNNI2
- 17 sur régions probableemnt fusionnée ou dupliquée par erreur. Ex KCNE1 (fausse duplication) qui a une couverture moindre (erreur d'alignement), KCNJ18 (duplication fusonée): mielleur couverture. Discordance entre les 2 référence sur ces gènes
- benchmark maison en GRCh38 et T2T 269 gènes d'intérêt clinique mais difficiles: dminution FP et FN avec T2T