5OWKP3I46JLY7OJC7OQFW7FLYU7Y3KZ3ACD5IS66SE4DVKVSRPWAC
3MMLXI5QRLKAAMSROH3E5L37TQYUOGWHTDOEABUD4DIMFGILKJBAC
RHWQQAAHNHFO3FLCGVB3SIDKNOUFJGZTDNN57IQVBMXXCWX74MKAC
ZV5YCHTZKVBL5ZW7NIQZLI5W6KEXSQ3T4XQNAPUQPJV2OJ6CKNLQC
RNTONMNJKHIFTXWDZ4MCP7UKZYNNAZOC3ELYRVZ47TYMM5Y5NYGQC
BRO5ESII2PGGU4FR7DCWHYC3QK22PL6WARZPKWTRMWBSQBDIYWCAC
PDH2BEBXR6WCCO2GRS3L6HMLQNCC2JU5BOKZLV4DLFPQD2UZFJKQC
FXA3ZBV64FML7W47IPHTAJFJHN3J3XHVHFVNYED47XFSBIGMBKRQC
EBFLN6GQZOT5RSTNMZRBJ2SYHHUOW7QRCJZRLSHAKDWIIIZKXH7AC
23PO3TWBWUMTYBXGMK3QSYX6MYP3LPMJXRKBT6WBQLLO5JJU4NRQC
3BYR6ZW3ZUHTCXNZ3GX6KLHAMNKEKFZVY4VSUYGMMZET47T4AR6QC
OEYACZNB7JKYKXLFHOU6DKLEWARDT53GSEFLBA6TGDLZTNCB4ZXQC
X5CU5SQWDBN2QUXTVD2GXXGFMRGYDAF55SWKJWZ5EX6LU2HXS4UAC
BKY4GUBANZGG3XDKLMG3HPPAMP457E6JMKTVAZM3RQ4NL5ZKKYSAC
CRN5TEHRZUS2JA26WEEHG65MT4ADK6DYAFFY6SCYRRPIF36NRU4AC
F4OH5H3ONZKUVBUI5NKYJ25B66FS22QQ7LRAO53OQX3ZBL2BL6JAC
FUEJT45VG4Q5LJEVCUGSEGMTQ4DUC2JLAEIHVRSFNG6MS53NIWNAC
NCBHHUESPLBDQI2VOWF74N3KMHS6MOSWQPSYBISHOTN2JQLSVLYQC
SYWVRVTSQFTJTEOZOEIIRDGPSJUEIIXT6FLFGXK7MJGL2P2C7MOQC
* <2023-11-30 Thu> Gymnastics
- "Roulade" avant barre symétrique : plutôt en force mais commence à passer
- 540 : travail
- trampoline : travail exercices de base
- champignon: 1 tour !
WES 129.94×
HiSeq2500 SRR1611184 SeqCap EZ Human Exome Lib v3.0 WES 111.90×
Kit acessible ?
**** Résumé
Kit disponible en hg38
| HiSeq 4000 | Agilent SureSelect v7 | SRX11061486 | https://github.com/kevinblighe/agilent |
| NovaSeq 6000 | Agilent SureSelect v7 | SRX11061516 | idem |
Kit disponible en hg19
| HiSeq2000 | SeqCap EZ Human Exome Lib v3.0 | SRR1611178 |http://hgdownload.soe.ucsc.edu/gbdb/hg19/exomeProbesets/
| HiSeq2000 | SeqCap EZ Human Exome Lib v3.0 | SRR1611179 |idem
| HiSeq2500 | SeqCap EZ Human Exome Lib v3.0 | SRR1611183 |idem
| HiSeq2500 | SeqCap EZ Human Exome Lib v3.0 | SRR1611184 |idem
https://emea.support.illumina.com/downloads/truseq-exome-product-files.html
*** Liste de capture
Agilent sureselect v7 hg19 et 38 https://github.com/kevinblighe/agilent
**** UCSCS
- [[http://hgdownload.soe.ucsc.edu/gbdb/hg19/exomeProbesets/][hg19]]
- [[http://hgdownload.soe.ucsc.edu/gbdb/hg38/exomeProbesets/][hg38]]
**** github aztrazeneca
https://github.com/AstraZeneca-NGS/reference_data
- IDT xGen Exome Research Panel v1.0
- Agilent SureSelect Human All Exon V6
- Agilent SureSelect Clinical Research Exome
- Nimblegen SeqCap EZ MedExome
- Nmblegen SeqCap EZ Exome v3
**** Trueseq
https://emea.support.illumina.com/downloads/truseq-exome-product-files.html
*** Exemple de validation avec bcbio:
Télécharge données + bed + liftover avec crossmap
https://github.com/bcbio/bcbio_validation_workflows/blob/master/giab-exome/input/get_data.sh
*** TODO Comment télécharger
**** DONE Tester ligne de commande
CLOSED: [2023-11-29 Wed 23:37] SCHEDULED: <2023-11-28 Tue>
***** KILL Tester aws
CLOSED: [2023-11-28 Tue 23:47] SCHEDULED: <2023-11-28 Tue>
Semble télécharger le .sra vu la taille (manque l'extension)
#+begin_src
aws s3 cp s3://sra-pub-run-odp/sra/SRR1611178/SRR1611178 --no-sign-request .
#+end_src
***** KILL Tester sra faster dump
CLOSED: [2023-11-29 Wed 22:20] SCHEDULED: <2023-11-28 Tue>
Selon la doc https://github.com/ncbi/sra-tools/wiki/08.-prefetch-and-fasterq-dump, il faut faire un "pré" - téléchargement
#+begin_src sh
prefetch SRR1611178
fastqer-dump SRR1611178
#+end_src
Note fasterq-dump créé un répertoire temporaire de la taille de prefetch et le supprime. Les fastq ne sont pas compressés
***** DONE Passer par ENA qui donne un lien vers FTP directement
CLOSED: [2023-11-29 Wed 23:37]
**** TODO Nextflow
***** KILL fromSRA
CLOSED: [2023-11-29 Wed 23:15]
Ne renvoie pas le FTP pour SRR1611178/SRR1611178 même avec clé API
**** TODO DataToolkit.jl
SCHEDULED: <2023-11-28 Tue>
- plusieurs datasets par patient appelé NA12878 par exemple mais avec attributs différents (séquencer, kit, pair1, pair2)
- FTP depuis ENA (FTP)
*** Zone de capture GIAB fourni le .bed pour l'exome . INfo : https://support.illumina.com/sequencing/sequencing_kits/nextera-rapid-capture-exome-kit/downloads.html
*** Valider la méthode
- 1000 genomes + SureSelect human all exon v2 target capture kit : non disponible sur le site d'agilent (V6 ou plus)
https://bmcbioinformatics.biomedcentral.com/articles/10.1186/s12859-019-2928-9
- GIAB + liftover du fichire de capture en hg38
Ce qui est aussi fait par
https://bcbio-nextgen.readthedocs.io/en/stable/contents/germline_variants.html
Mais avec UCSC liftover
** Centogène
https://www.twistbioscience.com/node/23906
Bed non fourni pour exactement cette capture
On prend https://www.twistbioscience.com/resources/data-files/twist-alliance-vcgs-exome-401mb-bed-files
qui content la majeure partie
* Réunion
** <2023-08-10 Thu> Alexis
Ok pour bloquer le développment d'ici mardi prochain
Dév:
- pipeline jusque VEP en T2T + GRCh38
- ok pour valider spip T2T sur quelques variant => à intégrer au pipeline
- annotation :
- ok pour mobidetails hg38
- +OMIM T2T+ non
- +franklin hg38+ non pour le moment
- métriques (fastq a minima) + rapport multiqc
- optionnel
- reformater la sortie
- on abandonne
- XAMScissors ave indel
- parallélisation haplotype caller
- spliceai à la vollée
- pangolin
Test
- GIAB:
- hg38: ok pour refaire les tests NA12878 avec données cento, sinon ok pour "c'est difficile" sur les 3 fichiers de capture
- T2T: ok pour faire des tests rapides mais probablement pas assez de temps !
- patient de synthèse : variant cento confirém par sanger seuls
Résultats
- ok pour scale up bwa mem et haplotyecaller
Manuscrit
- validation de méthode : laisser tomber la version actuelle et faire comme strasbourg (cf ngs diag) dans la présentatino
- a envoyé le powerponit avec les références des différsences articles
- ok pour robo4 si résultat
- architecture cible = VM : 78 coeurs 54Go RAUM et 1To espace disque
Passage en production : ok pour présentation rapide du code
* Nixpkgs :nix:
** DONE GATK
CLOSED: [2023-05-06 Sat 08:51]
*** DONE [[https://github.com/NixOS/nixpkgs/pull/185819][Binaire]]
CLOSED: [2022-09-10 Sat 23:53] SCHEDULED: <2022-08-10 Wed>
/Entered on/ [2022-08-09 Tue 10:57]
PR submitted
*** KILL Corriger code pour utiliser source
CLOSED: [2022-09-11 Sun 22:05]
*** DONE Corriger PATH pour include java et python
CLOSED: [2022-10-11 Tue 11:46]
https://github.com/NixOS/nixpkgs/pull/191548
Review <2022-10-10 Mon> , corrigé dans la journée
*** DONE Update 4.3.0.0
CLOSED: [2023-04-13 Thu 09:01]
** HOLD Nextflow
*** KILL version script seule
CLOSED: [2023-04-01 Sat 18:29]
Fix pour SGE et nextflow
https://github.com/NixOS/nixpkgs/issues/192396
*** KILL Version avec gradle
CLOSED: [2022-10-09 Sun 22:51]
*** HOLD [[https://github.com/NixOS/nixpkgs/issues/192396][Bug report Version 22.10.6]]
**** Notes
Erreur :
ERROR: Cannot download nextflow required file -- make sure you can connect to the internet
Alternatively you can try to download this file:
https://www.nextflow.io/releases/v22.10.6/nextflow-22.10.6-all.jar
and save it as:
.//nix/store/md2b1ah4d7ivj82k8xxap30dmdci00pa-nextflow-22.10.6/bin/.nextflow-wrapped
Dans la mise à jour, il y a la création d'un environnement virtuel qui casse l'exécution de nextflow (besoin de télécharger)
Fix = désactiver
**** KILL Patch NXF_OFFLINE=true
CLOSED: [2023-07-02 Sun 11:02] SCHEDULED: <2023-06-11 Sun>
** WAIT [[https://github.com/NixOS/nixpkgs/pull/249329][Multiqc]]
HG002,sanger-chr20,data/HG002-sanger-inserted-chr20_1.fq.gz,data/HG002-sanger-inserted-chr20_2.fq.gz
** KILL Mutalyzer
CLOSED: [2023-08-16 Wed 19:07] SCHEDULED: <2023-08-13 Sun>
Packaging faisable mais nombreux paquet python
** TODO Variant validator -> hgvs
C'est juste une interface autour d'hgvs mais il faut
- postgresql
- un accès ou télécharger des bases de données
Dépendences
s: wcwidth, pyee, pure-eval, ptyprocess, pickleshare, parsley, parse, fake-useragent, executing, backcall, appdirs, zipp, websockets, w3lib, urllib3, traitlets, tqdm, tabulate, sqlparse, soupsieve, six, pygments, psycopg2, prompt-toolkit, pexpect, parso, lxml, idna, humanfriendly, decorator, cython, cssselect, configparser, charset-normalizer, certifi, attrs, requests, pysam, pyquery, matplotlib-inline, jedi, importlib-metadata, coloredlogs, beautifulsoup4, asttokens, yoyo-migrations, stack-data, pyppeteer, bs4, bioutils, requests-html, ipython, biocommons.seqrepo, hgvs
** TODO SPIP :spip:
*** DONE PR upstream
CLOSED: [2023-08-12 Sat 18:23] SCHEDULED: <2023-08-12 Sat 18:00>
*** DONE Mail R. Lemann :T2T:
CLOSED: [2023-08-12 Sat 18:23] SCHEDULED: <2023-08-12 Sat 18:00>
*** KILL Mise à jour T2T :T2T:
*** WAIT Corriger PR
SCHEDULED: <2023-12-03 Sun>
** TODO VEP :vep:
*** DONE [[https://github.com/NixOS/nixpkgs/pull/185691][BioPerl]]
SCHEDULED: <2022-08-10 Wed>
/Entered on/ [2022-08-09 Tue 10:57]
PR submitted
*** TODO BioDBBBigFile
:PROPERTIES:
:ORDERED: t
:END:
/Entered on/ [2022-08-10 Wed 14:28]
On utilise la dernière version de kent, donc plus de problème.
PRête à être mergé. Rebase faite<2023-07-02 Sun>
**** DONE Venrsion de kent déjà packagée : forcer version 335
CLOSED: [2023-07-02 Sun 11:20]
***** KILL [[https://github.com/NixOS/nixpkgs/pull/206991][Restore building kent 404]]
CLOSED: [2023-05-06 Sat 17:40]
Review faite <2023-03-26 Sun> , atteinte merge]
Relancé <2023-05-06 Sat>
Kent 446 n'a pas ce problème donc PR inutile
***** DONE [[https://github.com/NixOS/nixpkgs/pull/223411][Ajouter les header to package]] (inc folder)
CLOSED: [2023-05-08 Mon 10:18] SCHEDULED: <2023-05-07 Sun>
Review à faire
https://github.com/NixOS/nixpkgs/pull/223411
Corrigé et plus besoin de la PR précédente
***** KILL [[https://github.com/NixOS/nixpkgs/pull/186462][BioDBBBigFile]] avec ces 2 changements
CLOSED: [2023-07-02 Sun 11:20]
**** KILL Version de kent déjà packagée : 404
CLOSED: [2023-03-27 Mon 16:43]
Compile mais les tests de passent pas
**** DONE Modifier selon PR https://github.com/NixOS/nixpkgs/pull/186462
CLOSED: [2023-07-30 Sun 22:01] SCHEDULED: <2023-07-30 Sun 20:00>
:LOGBOOK:
CLOCK: [2023-07-30 Sun 19:13]--[2023-07-30 Sun 20:50] => 1:37
:END:
Modification nécessaire pour kent :
- plus de patch
- suppression d'une boucle dans postPatch
On supprime aussi NIX_BUILD_TOP
**** WAIT Corriger PR biobigfile
SCHEDULED: <2023-12-05 Tue>
/Entered on/ [2023-10-15 Sun 17:21]
*** DONE [[https://github.com/NixOS/nixpkgs/pull/186459][BioDBHTS]]
CLOSED: [2023-05-06 Sat 08:49] SCHEDULED: <2023-04-15 Sat>
/Entered on/ [2022-08-10 Wed 14:28]
Correction pour review faites <2022-10-10 Mon>
*** DONE [[https://github.com/NixOS/nixpkgs/pull/186464][BioExtAlign]]
CLOSED: [2022-10-22 Sat 12:43] SCHEDULED: <2022-08-10 Wed>
/Entered on/ [2022-08-10 Wed 14:28]
Review <2022-10-10 Mon>, correction dans la journée.
Correction 2e pass
e, attente
Impossible de faire marcher les tests Car il ne trouve pas le module Bio::Tools::Align, qui est dans un dossier ailleurs dans le dépôt. Même en compilant tout le dépôt, cela ne fonctionne pas... On skip les tests.
*** TO
DO VEP
** WAIT [[https://github.com/NixOS/nixpkgs/pull/230394][rtg-tools]] :vcfeval:
Soumis
** WAIT Package Spip https://github.com/NixOS/nixpkgs/pull/247476
** TODO Happy :happy:
*** TODO PR python 3 upstream
SCHEDULED: <2023-12-02 Sat>
*** TODO nixpkgs en l'état
SCHEDULED: <2023-12-02 Sat>
** PROJ SpliceAI
** TODO Bamsurgeon
/Entered on/ [2023-05-13 Sat 19:11]
*** TODO Velvet
** TODO PR Picard avec option pour gérer la mémoire
Similaire à
https://github.com/bioconda/bioconda-recipes/blob/master/r
ecipes/picard/picard.sh
* Julia :julia:
** KILL XAM.jl: PR pour modification record :julia:
CLOSED: [2023-05-29 Mon 15:40] SCHEDULED: <2023-05-28 Sun>
/Entered on/ [2023-05-27 Sat 22:39]
** TODO XAMscissors.jl :xamscissors:
Modification de la séquence dans BAM.
*Pas de mise à jour de CIGAR*
On convertit en fastq et on lance le pipeline pour "corriger"
#+begin_src sh
cd /home/alex/code/bisonex/out/63003856/preprocessing/mapped
samtools view 63003856_S135.bam NC_000022.11 -o 63003856_S135_chr22.bam
cd /home/alex/recherche/bisonex/code/BamScissors.jl
cp ~/code/bisonex/out/63003856/preprocessing/mapped/63003856_S135_chr22.bam .
samtools index 63003856_chr22.bam
#+end_src
Le script va modifier le bam, le trier et générer le fastq. !!!
Attention: ne pas oublier l'option -n !!!
#+begin_src sh
time julia --project=.. insertVariant.jl
scp 63003856_S135_chr22_{1,2}.fq.gz meso:/Work/Users/apraga/bisonex/tests/bamscissors/
#+end_src
*** WAIT Implémenter les SNV avec VAF :snv:
Stratégie :
1. calculer la profondeur sur les positions
2. créer un dictionnaire { nom du reads : position dataframe }
3. itérer sur tous les reads et changer ceux marqués
**** DONE VAF = 1
CLOSED: [2023-05-29 Mon 15:34]
**** DONE VAF selon loi normale
CLOSED: [2023-05-29 Mon 15:35]
Tronquée si > 1
**** WAIT Tests unitaires
***** DONE NA12878: 1 gène sur chromosome 22
CLOSED: [2023-05-30 Tue 23:55]
root = "https://ftp-trace.ncbi.nlm.nih.gov/ReferenceSamples/giab/data/NA12878/Garvan_NA12878_HG001_HiSeq_Exome/"
#+begin_src sh
samtools view project.NIST_NIST7035_H7AP8ADXX_NA12878.bwa.markDuplicates.bam chr22 -o project.NIST_NIST7035_H7AP8ADXX_NA12878_chr22.bam
samtools view project.NIST_NIST7035_H7AP8ADXX_NA12878_chr22.bam chr22:19419700-19424000 -o NIST7035_H7AP8ADXX_NA12878_chr22_MRPL40_hg19.bam
#+end_src
***** WAIT Pull request formatspeciment
https://github.com/BioJulia/FormatSpecimens.jl/pull/8
***** DONE Formatspecimens
CLOSED: [2023-05-29 Mon 23:03]
****** DONE 1 read
CLOSED: [2023-05-29 Mon 23:02]
****** DONE VAF sur 1 exon
CLOSED: [2023-05-29 Mon 23:03]
**** DONE [#A] Bug: perte de nombreux reads avec NA12878
CLOSED: [2023-08-19 Sat 20:45] SCHEDULED: <2023-08-18 Fri>
:PROPERTIES:
:ID: 5c1c36f3-f68e-4e6d-a7b6-61dca89abc37
:END:
Ex: chrX:g.124056226 : on passe de 65 reads à 1
Test xamscissors: pas de soucis...
On teste sur cette position +/- 200bp
#+begin_src sh :dir /home/alex/roam/research/bisonex/code/sanger
samtools view /home/alex/code/bisonex/out/2300346867_NA12878-63118093_S260-GRCh38/preprocessing/mapped/2300346867_NA12878-63118093_S260-GRCh38.bam chrX:124056026-124056426 -o chrXsmall.bam
#+end_src
#+RESULTS:
***** DONE Vérifier profondeur avec dernière version :
CLOSED: [2023-08-19 Sat 20:34] SCHEDULED: <2023-08-19 Sat>
****** DONE chr20: profondeur ok
SCHEDULED: <2023-08-19 Sat>
****** DONE toutes les données
CLOSED: [2023-08-19 Sat 20:34] SCHEDULED: <2023-08-19 Sat>
Ok pour 7 variants (IGV) notament chromosome X
*** TODO Implémenter les indel avec VAF :indel:
*** TODO Soumission paquet
* Données
:PROPERTIES:
:CATEGORY: data
:END:
** DONE Remplacer bam par fastq sur mesocentre
CLOSED: [2023-04-16 Sun 16:33]
Commande
*** DONE Supprimer les fastq non "paired"
CLOSED: [2023-04-16 Sun 16:33]
nushell
Liste des fastq avec "paired-end" manquant
#+begin_src nu
ls **/*.fastq.gz | get name | path basename | split column "_" | get column1 | uniq -u | save single.txt
#+end_src
#+RESULTS:
: 62907927
: 62907970
: 62899606
: 62911287
: 62913201
: 62914084
: 62915905
: 62921595
: 62923065
: 62925220
: 62926503
: 62926502
: 62926500
: 62926499
: 62926498
: 62931719
: 62943423
: 62943400
: 62948290
: 62949205
: 62949206
: 62949118
: 62951284
: 62960792
: 62960785
: 62960787
: 62960617
: 62962561
: 62962692
: 62967473
: 62972194
: 62979102
On vérifie
#+begin_src nu
open single.txt | lines | each {|e| ls $"fastq/*_($in)/*" | get 0 }
open single.txt | lines | each {|e| ls $"fastq/*_($in)/*" | get 0.name } | path basename | split column "_" | get column1 | uniq -c
#+end_src
On met tous dans un dossier (pas de suppression )
#+begin_src
open single.txt | lines | each {|e| ls $"fastq/*_($in)/*" | get 0 } | each {|e| ^mv $e.name bad-fastq/}
#+end_src
On vérifie que les dossiier sont videsj
open single.txt | lines | each {|e| ls $"fastq/*_($in)" | get 0.name } | ^ls -l $in
Puis on supprime
open single.txt | lines | each {|e| ls $"fastq/*_($in)" | get 0.name } | ^rm -r $in
*** DONE Supprimer bam qui ont des fastq
CLOSED: [2023-04-16 Sun 16:33]
On liste les identifiants des fastq et bam dans un tableau avec leur type :
#+begin_src
let fastq = (ls fastq/*/*.fastq.gz | get name | parse "{dir}/{full_id}/{id}_{R}_001.fastq.gz" | select dir id | uniq )
let bam = (ls bam/*/*.bam | get name | parse "{dir}/{full_id}/{id}_{S}.bqrt.bam" | select dir id)
#+end_src
On groupe les résultat par identifiant (résultats = liste de records qui doit être convertie en table)
et on trie ceux qui n'ont qu'un fastq ou un bam
#+begin_src
let single = ( $bam | append $fastq | group-by id | transpose id files | get files | where {|x| ($x | length) == 1})
#+end_src
On convertit en table et on récupère seulement les bam
#+begin_src
$single | reduce {|it, acc| $acc | append $it} | where dir == bam | get id | each {|e| ^ls $"bam/*_($e)/*.bam"}
#+end_src
#+RESULTS:
: bam/2100656174_62913201/62913201_S52.bqrt.bam
: bam/2100733271_62925220/62925220_S33.bqrt.bam
: bam/2100738763_62926502/62926502_S108.bqrt.bam
: bam/2100746726_62926498/62926498_S105.bqrt.bam
: bam/2100787936_62931955/62931955_S4.bqrt.bam
: bam/2200066374_62948290/62948290_S130.bqrt.
WES 129.94×
HiSeq2500 SRR1611184 SeqCap EZ Human Exome Lib v3.0 WES 111.90×
Kit acessible ?
**** Résumé
Kit disponible en hg38
| HiSeq 4000 | Agilent SureSelect v7 | SRX11061486 | https://github.com/kevinblighe/agilent |
| NovaSeq 6000 | Agilent SureSelect v7 | SRX11061516 | idem |
Kit disponible en hg19
| HiSeq2000 | SeqCap EZ Human Exome Lib v3.0 | SRR1611178 |http://hgdownload.soe.ucsc.edu/gbdb/hg19/exomeProbesets/
| HiSeq2000 | SeqCap EZ Human Exome Lib v3.0 | SRR1611179 |idem
| HiSeq2500 | SeqCap EZ Human Exome Lib v3.0 | SRR1611183 |idem
| HiSeq2500 | SeqCap EZ Human Exome Lib v3.0 | SRR1611184 |idem
https://emea.support.illumina.com/downloads/truseq-exome-product-files.html
*** Liste de capture
Agilent sureselect v7 hg19 et 38 https://github.com/kevinblighe/agilent
**** UCSCS
- [[http://hgdownload.soe.ucsc.edu/gbdb/hg19/exomeProbesets/][hg19]]
- [[http://hgdownload.soe.ucsc.edu/gbdb/hg38/exomeProbesets/][hg38]]
**** github aztrazeneca
https://github.com/AstraZeneca-NGS/reference_data
- IDT xGen Exome Research Panel v1.0
- Agilent SureSelect Human All Exon V6
- Agilent SureSelect Clinical Research Exome
- Nimblegen SeqCap EZ MedExome
- Nmblegen SeqCap EZ Exome v3
**** Trueseq
https://emea.support.illumina.com/downloads/truseq-exome-product-files.html
*** Exemple de validation avec bcbio:
Télécharge données + bed + liftover avec crossmap
https://github.com/bcbio/bcbio_validation_workflows/blob/master/giab-exome/input/get_data.sh
*** TODO Comment télécharger
**** DONE Tester ligne de commande
CLOSED: [2023-11-29 Wed 23:37] SCHEDULED: <2023-11-28 Tue>
***** KILL Tester aws
CLOSED: [2023-11-28 Tue 23:47] SCHEDULED: <2023-11-28 Tue>
Semble télécharger le .sra vu la taille (manque l'extension)
#+begin_src
aws s3 cp s3://sra-pub-run-odp/sra/SRR1611178/SRR1611178 --no-sign-request .
#+end_src
***** KILL Tester sra faster dump
CLOSED: [2023-11-29 Wed 22:20] SCHEDULED: <2023-11-28 Tue>
Selon la doc https://github.com/ncbi/sra-tools/wiki/08.-prefetch-and-fasterq-dump, il faut faire un "pré" - téléchargement
#+begin_src sh
prefetch SRR1611178
fastqer-dump SRR1611178
#+end_src
Note fasterq-dump créé un répertoire temporaire de la taille de prefetch et le supprime. Les fastq ne sont pas compressés
***** DONE Passer par ENA qui donne un lien vers FTP directement
CLOSED: [2023-11-29 Wed 23:37]
**** TODO Nextflow
***** KILL fromSRA
CLOSED: [2023-11-29 Wed 23:15]
Ne renvoie pas le FTP pour SRR1611178/SRR1611178 même avec clé API
**** TODO DataToolkit.jl
SCHEDULED: <2023-11-28 Tue>
- plusieurs datasets par patient appelé NA12878 par exemple mais avec attributs différents (séquencer, kit, pair1, pair2)
- FTP depuis ENA (FTP)
*** Zone de capture GIAB fourni le .bed pour l'exome . INfo : https://support.illumina.com/sequencing/sequencing_kits/nextera-rapid-capture-exome-kit/downloads.html
*** Valider la méthode
- 1000 genomes + SureSelect human all exon v2 target capture kit : non disponible sur le site d'agilent (V6 ou plus)
https://bmcbioinformatics.biomedcentral.com/articles/10.1186/s12859-019-2928-9
- GIAB + liftover du fichire de capture en hg38
Ce qui est aussi fait par
https://bcbio-nextgen.readthedocs.io/en/stable/contents/germline_variants.html
Mais avec UCSC liftover
** Centogène
https://www.twistbioscience.com/node/23906
Bed non fourni pour exactement cette capture
On prend https://www.twistbioscience.com/resources/data-files/twist-alliance-vcgs-exome-401mb-bed-files
qui content la majeure partie
* Réunion
** <2023-08-10 Thu> Alexis
Ok pour bloquer le développment d'ici mardi prochain
Dév:
- pipeline jusque VEP en T2T + GRCh38
- ok pour valider spip T2T sur quelques variant => à intégrer au pipeline
- annotation :
- ok pour mobidetails hg38
- +OMIM T2T+ non
- +franklin hg38+ non pour le moment
- métriques (fastq a minima) + rapport multiqc
- optionnel
- reformater la sortie
- on abandonne
- XAMScissors ave indel
- parallélisation haplotype caller
- spliceai à la vollée
- pangolin
Test
- GIAB:
- hg38: ok pour refaire les tests NA12878 avec données cento, sinon ok pour "c'est difficile" sur les 3 fichiers de capture
- T2T: ok pour faire des tests rapides mais probablement pas assez de temps !
- patient de synthèse : variant cento confirém par sanger seuls
Résultats
- ok pour scale up bwa mem et haplotyecaller
Manuscrit
- validation de méthode : laisser tomber la version actuelle et faire comme strasbourg (cf ngs diag) dans la présentatino
- a envoyé le powerponit avec les références des différsences articles
- ok pour robo4 si résultat
- architecture cible = VM : 78 coeurs 54Go RAUM et 1To espace disque
Passage en production : ok pour présentation rapide du code
* Nixpkgs :nix:
** DONE GATK
CLOSED: [2023-05-06 Sat 08:51]
*** DONE [[https://github.com/NixOS/nixpkgs/pull/185819][Binaire]]
CLOSED: [2022-09-10 Sat 23:53] SCHEDULED: <2022-08-10 Wed>
/Entered on/ [2022-08-09 Tue 10:57]
PR submitted
*** KILL Corriger code pour utiliser source
CLOSED: [2022-09-11 Sun 22:05]
*** DONE Corriger PATH pour include java et python
CLOSED: [2022-10-11 Tue 11:46]
https://github.com/NixOS/nixpkgs/pull/191548
Review <2022-10-10 Mon> , corrigé dans la journée
*** DONE Update 4.3.0.0
CLOSED: [2023-04-13 Thu 09:01]
** HOLD Nextflow
*** KILL version script seule
CLOSED: [2023-04-01 Sat 18:29]
Fix pour SGE et nextflow
https://github.com/NixOS/nixpkgs/issues/192396
*** KILL Version avec gradle
CLOSED: [2022-10-09 Sun 22:51]
*** HOLD [[https://github.com/NixOS/nixpkgs/issues/192396][Bug report Version 22.10.6]]
**** Notes
Erreur :
ERROR: Cannot download nextflow required file -- make sure you can connect to the internet
Alternatively you can try to download this file:
https://www.nextflow.io/releases/v22.10.6/nextflow-22.10.6-all.jar
and save it as:
.//nix/store/md2b1ah4d7ivj82k8xxap30dmdci00pa-nextflow-22.10.6/bin/.nextflow-wrapped
Dans la mise à jour, il y a la création d'un environnement virtuel qui casse l'exécution de nextflow (besoin de télécharger)
Fix = désactiver
**** KILL Patch NXF_OFFLINE=true
CLOSED: [2023-07-02 Sun 11:02] SCHEDULED: <2023-06-11 Sun>
** WAIT [[https://github.com/NixOS/nixpkgs/pull/249329][Multiqc]]
HG002,sanger-chr20,data/HG002-sanger-inserted-chr20_1.fq.gz,data/HG002-sanger-inserted-chr20_2.fq.gz
** KILL Mutalyzer
CLOSED: [2023-08-16 Wed 19:07] SCHEDULED: <2023-08-13 Sun>
Packaging faisable mais nombreux paquet python
** TODO Variant validator -> hgvs
C'est juste une interface autour d'hgvs mais il faut
- postgresql
- un accès ou télécharger des bases de données
Dépendences
s: wcwidth, pyee, pure-eval, ptyprocess, pickleshare, parsley, parse, fake-useragent, executing, backcall, appdirs, zipp, websockets, w3lib, urllib3, traitlets, tqdm, tabulate, sqlparse, soupsieve, six, pygments, psycopg2, prompt-toolkit, pexpect, parso, lxml, idna, humanfriendly, decorator, cython, cssselect, configparser, charset-normalizer, certifi, attrs, requests, pysam, pyquery, matplotlib-inline, jedi, importlib-metadata, coloredlogs, beautifulsoup4, asttokens, yoyo-migrations, stack-data, pyppeteer, bs4, bioutils, requests-html, ipython, biocommons.seqrepo, hgvs
** TODO SPIP :spip:
*** DONE PR upstream
CLOSED: [2023-08-12 Sat 18:23] SCHEDULED: <2023-08-12 Sat 18:00>
*** DONE Mail R. Lemann :T2T:
CLOSED: [2023-08-12 Sat 18:23] SCHEDULED: <2023-08-12 Sat 18:00>
*** KILL Mise à jour T2T :T2T:
*** WAIT Corriger PR
SCHEDULED: <2023-12-03 Sun>
** TODO VEP :vep:
*** DONE [[https://github.com/NixOS/nixpkgs/pull/185691][BioPerl]]
SCHEDULED: <2022-08-10 Wed>
/Entered on/ [2022-08-09 Tue 10:57]
PR submitted
*** DONE BioDBBBigFile
CLOSED: [2023-11-30 Thu 21:52]
:PROPERTIES:
:ORDERED: t
:END:
/Entered on/ [2022-08-10 Wed 14:28]
On utilise la dernière version de kent, donc plus de problème.
PRête à être mergé. Rebase faite<2023-07-02 Sun>
**** DONE Version de kent déjà packagée : forcer version 335
CLOSED: [2023-07-02 Sun 11:20]
***** KILL [[https://github.com/NixOS/nixpkgs/pull/206991][Restore building kent 404]]
CLOSED: [2023-05-06 Sat 17:40]
Review faite <2023-03-26 Sun> , atteinte merge]
Relancé <2023-05-06 Sat>
Kent 446 n'a pas ce problème donc PR inutile
***** DONE [[https://github.com/NixOS/nixpkgs/pull/223411][Ajouter les header to package]] (inc folder)
CLOSED: [2023-05-08 Mon 10:18] SCHEDULED: <2023-05-07 Sun>
Review à faire
https://github.com/NixOS/nixpkgs/pull/223411
Corrigé et plus besoin de la PR précédente
***** KILL [[https://github.com/NixOS/nixpkgs/pull/186462][BioDBBBigFile]] avec ces 2 changements
CLOSED: [2023-07-02 Sun 11:20]
**** KILL Version de kent déjà packagée : 404
CLOSED: [2023-03-27 Mon 16:43]
Compile mais les tests de passent pas
**** DONE Modifier selon PR https://github.com/NixOS/nixpkgs/pull/186462
CLOSED: [2023-07-30 Sun 22:01] SCHEDULED: <2023-07-30 Sun 20:00>
:LOGBOOK:
CLOCK: [2023-07-30 Sun 19:13]--[2023-07-30 Sun 20:50] => 1:37
:END:
Modification nécessaire pour kent :
- plus de patch
- suppression d'une boucle dans postPatch
On supprime aussi NIX_BUILD_TOP
**** DONE Corriger PR biobigfile
CLOSED: [2023-11-30 Thu 21:52] SCHEDULED: <2023-12-05 Tue>
/Entered on/ [2023-10-15 Sun 17:21]
*** DONE [[https://github.com/NixOS/nixpkgs/pull/186459][BioDBHTS]]
CLOSED: [2023-05-06 Sat 08:49] SCHEDULED: <2023-04-15 Sat>
/Entered on/ [2022-08-10 Wed 14:28]
Correction pour review faites <2022-10-10 Mon>
*** DONE [[https://github.com/NixOS/nixpkgs/pull/186464][BioExtAlign]]
CLOSED: [2022-10-22 Sat 12:43] SCHEDULED: <2022-08-10 Wed>
/Entered on/ [2022-08-10 Wed 14:28]
Review <2022-10-10 Mon>, correction dans la journée.
Correction 2e passe, attente
Impossible de faire marcher les tests Car il ne trouve pas le module Bio::Tools::Align, qui est dans un dossier ailleurs dans le dépôt. Même en compilant tout le dépôt, cela ne fonctionne pas... On skip les tests.
*** TODO VEP
SCHEDULED: <2023-12-07 Thu>
** WAIT [[https://github.com/NixOS/nixpkgs/pull/230394][rtg-tools]] :vcfeval:
Soumis
** WAIT Package Spip https://github.com/NixOS/nixpkgs/pull/247476
** TODO Happy :happy:
*** TODO PR python 3 upstream
SCHEDULED: <2023-12-02 Sat>
*** TODO nixpkgs en l'état
SCHEDULED: <2023-12-02 Sat>
** PROJ SpliceAI
** TODO Bamsurgeon
/Entered on/ [2023-05-13 Sat 19:11]
*** TODO Velvet
** TODO PR Picard avec option pour gérer la mémoire
Similaire à
https://github.com/bioconda/bioconda-recipes/blob/master/recipes/picard/picard.sh
* Julia :julia:
** KILL XAM.jl: PR pour modification record :julia:
CLOSED: [2023-05-29 Mon 15:40] SCHEDULED: <2023-05-28 Sun>
/Entered on/ [2023-05-27 Sat 22:39]
** TODO XAMscissors.jl :xamscissors:
Modification de la séquence dans BAM.
*Pas de mise à jour de CIGAR*
On convertit en fastq et on lance le pipeline pour "corriger"
#+begin_src sh
cd /home/alex/code/bisonex/out/63003856/preprocessing/mapped
samtools view 63003856_S135.bam NC_000022.11 -o 63003856_S135_chr22.bam
cd /home/alex/recherche/bisonex/code/BamScissors.jl
cp ~/code/bisonex/out/63003856/preprocessing/mapped/63003856_S135_chr22.bam .
samtools index 63003856_chr22.bam
#+end_src
Le script va modifier le bam, le trier et générer le fastq. !!!
Attention: ne pas oublier l'option -n !!!
#+begin_src sh
time julia --project=.. insertVariant.jl
scp 63003856_S135_chr22_{1,2}.fq.gz meso:/Work/Users/apraga/bisonex/tests/bamscissors/
#+end_src
*** WAIT Implémenter les SNV avec VAF :snv:
Stratégie :
1. calculer la profondeur sur les positions
2. créer un dictionnaire { nom du reads : position dataframe }
3. itérer sur tous les reads et changer ceux marqués
**** DONE VAF = 1
CLOSED: [2023-05-29 Mon 15:34]
**** DONE VAF selon loi normale
CLOSED: [2023-05-29 Mon 15:35]
Tronquée si > 1
**** WAIT Tests unitaires
***** DONE NA12878: 1 gène sur chromosome 22
CLOSED: [2023-05-30 Tue 23:55]
root = "https://ftp-trace.ncbi.nlm.nih.gov/ReferenceSamples/giab/data/NA12878/Garvan_NA12878_HG001_HiSeq_Exome/"
#+begin_src sh
samtools view project.NIST_NIST7035_H7AP8ADXX_NA12878.bwa.markDuplicates.bam chr22 -o project.NIST_NIST7035_H7AP8ADXX_NA12878_chr22.bam
samtools view project.NIST_NIST7035_H7AP8ADXX_NA12878_chr22.bam chr22:19419700-19424000 -o NIST7035_H7AP8ADXX_NA12878_chr22_MRPL40_hg19.bam
#+end_src
***** WAIT Pull request formatspeciment
https://github.com/BioJulia/FormatSpecimens.jl/pull/8
***** DONE Formatspecimens
CLOSED: [2023-05-29 Mon 23:03]
****** DONE 1 read
CLOSED: [2023-05-29 Mon 23:02]
****** DONE VAF sur 1 exon
CLOSED: [2023-05-29 Mon 23:03]
**** DONE [#A] Bug: perte de nombreux reads avec NA12878
CLOSED: [2023-08-19 Sat 20:45] SCHEDULED: <2023-08-18 Fri>
:PROPERTIES:
:ID: 5c1c36f3-f68e-4e6d-a7b6-61dca89abc37
:END:
Ex: chrX:g.124056226 : on passe de 65 reads à 1
Test xamscissors: pas de soucis...
On teste sur cette position +/- 200bp
#+begin_src sh :dir /home/alex/roam/research/bisonex/code/sanger
samtools view /home/alex/code/bisonex/out/2300346867_NA12878-63118093_S260-GRCh38/preprocessing/mapped/2300346867_NA12878-63118093_S260-GRCh38.bam chrX:124056026-124056426 -o chrXsmall.bam
#+end_src
#+RESULTS:
***** DONE Vérifier profondeur avec dernière version :
CLOSED: [2023-08-19 Sat 20:34] SCHEDULED: <2023-08-19 Sat>
****** DONE chr20: profondeur ok
SCHEDULED: <2023-08-19 Sat>
****** DONE toutes les données
CLOSED: [2023-08-19 Sat 20:34] SCHEDULED: <2023-08-19 Sat>
Ok pour 7 variants (IGV) notament chromosome X
*** TODO Implémenter les indel avec VAF :indel:
*** TODO Soumission paquet
* Données
:PROPERTIES:
:CATEGORY: data
:END:
** DONE Remplacer bam par fastq sur mesocentre
CLOSED: [2023-04-16 Sun 16:33]
Commande
*** DONE Supprimer les fastq non "paired"
CLOSED: [2023-04-16 Sun 16:33]
nushell
Liste des fastq avec "paired-end" manquant
#+begin_src nu
ls **/*.fastq.gz | get name | path basename | split column "_" | get column1 | uniq -u | save single.txt
#+end_src
#+RESULTS:
: 62907927
: 62907970
: 62899606
: 62911287
: 62913201
: 62914084
: 62915905
: 62921595
: 62923065
: 62925220
: 62926503
: 62926502
: 62926500
: 62926499
: 62926498
: 62931719
: 62943423
: 62943400
: 62948290
: 62949205
: 62949206
: 62949118
: 62951284
: 62960792
: 62960785
: 62960787
: 62960617
: 62962561
: 62962692
: 62967473
: 62972194
: 62979102
On vérifie
#+begin_src nu
open single.txt | lines | each {|e| ls $"fastq/*_($in)/*" | get 0 }
open single.txt | lines | each {|e| ls $"fastq/*_($in)/*" | get 0.name } | path basename | split column "_" | get column1 | uniq -c
#+end_src
On met tous dans un dossier (pas de suppression )
#+begin_src
open single.txt | lines | each {|e| ls $"fastq/*_($in)/*" | get 0 } | each {|e| ^mv $e.name bad-fastq/}
#+end_src
On vérifie que les dossiier sont videsj
open single.txt | lines | each {|e| ls $"fastq/*_($in)" | get 0.name } | ^ls -l $in
Puis on supprime
open single.txt | lines | each {|e| ls $"fastq/*_($in)" | get 0.name } | ^rm -r $in
*** DONE Supprimer bam qui ont des fastq
CLOSED: [2023-04-16 Sun 16:33]
On liste les identifiants des fastq et bam dans un tableau avec leur type :
#+begin_src
let fastq = (ls fastq/*/*.fastq.gz | get name | parse "{dir}/{full_id}/{id}_{R}_001.fastq.gz" | select dir id | uniq )
let bam = (ls bam/*/*.bam | get name | parse "{dir}/{full_id}/{id}_{S}.bqrt.bam" | select dir id)
#+end_src
On groupe les résultat par identifiant (résultats = liste de records qui doit être convertie en table)
et on trie ceux qui n'ont qu'un fastq ou un bam
#+begin_src
let single = ( $bam | append $fastq | group-by id | transpose id files | get files | where {|x| ($x | length) == 1})
#+end_src
On convertit en table et on récupère seulement les bam
#+begin_src
$single | reduce {|it, acc| $acc | append $it} | where dir == bam | get id | each {|e| ^ls $"bam/*_($e)/*.bam"}
#+end_src
#+RESULTS:
: bam/2100656174_62913201/62913201_S52.bqrt.bam
: bam/2100733271_62925220/62925220_S33.bqrt.bam
: bam/2100738763_62926502/62926502_S108.bqrt.bam
: bam/2100746726_62926498/62926498_S105.bqrt.bam
: bam/2100787936_62931955/62931955_S4.bqrt.bam
: bam/2200066374_62948290/62948290_S130.bqrt.
| | |
| NC_000020.11 | 62172726 | rs36106901 | G | A | |
| NC_000020.11 | 62172726 | 981031 | G | A | Conflicting_interpretations_of_pathogenicity |
| ------ | | | | | |
| NC_000020.11 | 63349782 | rs1044396 | G | A,C | |
| NC_000020.11 | 63349782 | 93427 | G | A | Benign |
| NC_000020.11 | 63349782 | 857384 | G | C | Conflicting_interpretations_of_pathogenicity |
| ------ | | | | | |
| NC_000020.11 | 63414925 | rs1801545 | G | A,C,T | |
| NC_000020.11 | 63414925 | 194284
| G | A | Conflicting_interpretations_of_pathogenicity |
| NC_000020.11 | 63414925 | 129337 | G | C | Benign |
| NC_000020.11 | 63414925 | 851545 | GG | CA | Uncertain_significance |
| ------ | | | | | |
On a donc plusieurs problèmes :
1. isec devrait fonctionner au moins sur
| NC_000020.11 | 25390747 | rs373200654 | G | C | |
| NC_000020.11 | 25390747 | 338000 | G | C | Conflicting_interpretations_of_pathogenicity |
On teste juste sur cette ligne
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter -i 'POS=25390747' clinvar_chr20.vcf.gz -o clinvar_test.vcf.gz
bcftools filter -i 'POS=25390747' dbSNP_common_chr20.vcf.gz -o dbSNP_test.vcf.gz
#+end_src
On retrouve bien la ligne dans l'intersection...
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter -i 'POS=25390747' clinvar_chr20.vcf.gz -o clinvar_test.vcf.gz
bcftools index dbSNP_test.vcf.gz dbSNP_test.vcf.gz
bcftools index dbSNP_test.vcf.gz clinvar_test.vcf.gz
bcftools isec dbSNP_test.vcf.gz clinvar_test.vcf.gz -p test
#+end_src
#+RESULTS:
2. isec ne semble pas fonctionner sur en cas d'ALT multiples
| NC_000020.11 | 32800145 | rs2424926 | C | G,T | |
| NC_000020.11 | 32800145 | 338173 | C | G | Benign |
| NC_000020.11 | 32800145 | 338174 | C | T | Conflicting_interpretations_of_pathogenicity |
| | | | | | |
3. s'il y a plusieurs variantions à une position, il faut bien vérifier que tous ne sont pas patho.
La version d'Alexis le fait bien
| NC_000020.11 | 3234173 | rs3827075 | T | A,C,G | |
| NC_000020.11 | 3234173 | 262001 | T | G | Conflicting_interpretations_of_pathogenicity |
| NC_000020.11 | 3234173 | 1072511 | T | TGGCGAAGC | Pathogenic |
| NC_000020.11 | 3234173 | 208613 | TGGCGAAGC | G | Pathogenic |
| NC_000020.11 | 3234173 | 1312 | TGGCGAAGC | T | Pathogenic |
****** DONE Voir si isec gère les multiallélique (chr20) : non, impossible de faire marcher
CLOSED: [2022-11-27 Sun 00:37]
******* DONE chr20 en prenant un patho clinvar aussi dans dbSNP
CLOSED: [2022-11-27 Sun 00:37]
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter dbSNP_common_chr20.vcf.gz -i 'POS=10652589' -o test_dbsnp.vcf.gz
bcftools filter clinvar_chr20.vcf.gz -i 'POS=10652589' -o test_clinvar.vcf.gz
bcftools index test_dbsnp.vcf.gz
bcftools index test_clinvar.vcf.gz
#+end_src
#+RESULTS:
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools isec test_dbsnp.vcf.gz test_clinvar.vcf.gz -p tmp
grep '^[^#]' tmp/0002.vcf
grep '^[^#]' tmp/0003.vcf
#+end_src
#+RESULTS:
Même en biallélique, ne fonctionne pas.
Testé en modifiant test_dbsnp !
Fonctionne avec un variant par ligne
****** DONE isec en coupant les sites multialléliques: non
CLOSED: [2022-11-27 Sun 00:37]
******* DONE Exemple simple ok
CLOSED: [2022-11-27 Sun 00:34]
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter -i 'POS=10652589' dbSNP_common_chr20.vcf.gz -o dbsnp_mwi.vcf.gz
bcftools filter -i 'POS=10652589' clinvar_chr20.vcf.gz -o clinvar_mwi.vcf.gz
bcftools index -f dbsnp_mwi.vcf.gz
bcftools index -f clinvar_mwi.vcf.gz
bcftools isec dbsnp_mwi.vcf.gz clinvar_mwi.vcf.gz -n=2
#+end_src
#+RESULTS:
Même en biallélique, ne fonctionne pas.
Chr 20
Avec les fichiers du teste précédent
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools norm -m -any dbsnp_mwi.vcf.gz -o dbsnp_mwi_norm.vcf.gz
bcftools index dbsnp_mwi_norm.vcf.gz
bcftools isec dbsnp_mwi_norm.vcf.gz clinvar_mwi.vcf.gz -n=2
#+end_src
#+RESULTS:
| NC_000020.11 | 10652589 | G | A | 11 |
| NC_000020.11 | 10652589 | G | C | 11 |
******* DONE Sur dbSNP chr20 non
CLOSED: [2023-10-07 Sat 17:57]
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools norm -m -any dbSNP_common_chr20 -o dbSNP_common_chr20_norm.vcf.gz
#+end_src
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools isec -i 'INFO/CLNSIG="Pathogenic"' dbSNP_common_chr20_norm.vcf.gz clinvar_chr20.vcf.gz -p tmp
#+end_src
#+RESULTS:
***** DONE Essai bedtools intersect
#+begin_src sh
bedtools intersect -a dbSNP_common.vcf.gz -b clinvar.vcf.gz
#+end_src
$ wc -l intersect.vcf
220206 intersect.vcf
** TODO Dépendences avec Nix
*** DONE GATK
CLOSED: [2022-10-21 Fri 21:59]
*** WAIT BioDBHTS
Contribuer pull request
*** DONE BioExtAlign
CLOSED: [2022-10-22 Sat 00:38]
*** WAIT BioBigFile
Revoir si on peut utliser kent dernière version
Contribuer pull request
*** HOLD rtg-tools
Convertir clinvar NC
*** DONE simuscop
CLOSED: [2022-12-30 Fri 22:31]
*** DONE Spip
CLOSED: [2022-12-04 Sun 12:49]
Pas de pull request
*** DONE R + packages
CLOSED: [2022-11-19 Sat 21:05]
*** TODO hap.py
https://github.com/Illumina/hap.py
**** DONE Version sans rtgtools avec python 3
CLOSED: [2023-02-02 Thu 22:15]
Procédure pour tester
#+begin_src
nix develop .#hap-py
$ genericBuild
#+end_src
1. Supprimer l’appel à make_dependencies dans cmakelist.txt : on peut tout installer avec nix
2. Patch Roc.cpp pour avoir numeric_limits ( error: 'numeric_limits' is not a member of 'std')
3. ajout de flags de link (essai, error)
set(ZLIB_LIBRARIES -lz -lbz2 -lcurl -lcrypto -llzma)
4. Changer les appels à print en print() dans le code python et suppression de quelques import
[nix-shell:~/source]$ sed -i.orig 's/print \"\(.*\)"/print(\1)/' src/python/*.py
**** DONE Sérialiser json pour écrire données de sorties
CLOSED: [2023-02-17 Fri 19:25]
**** DONE Tester sur example
CLOSED: [2023-02-04 Sat 00:25]
#+begin_src sh
$ cd hap.py
$ ../result/bin/hap.py example/happy/PG_NA12878_chr21.vcf.gz example/happy/NA12878_chr21.vcf.gz -f example/happy/PG_Conf_chr21.bed.gz -o test -r example/chr21.fa
#+end_src
#+RESULTS:
| Type | Filter | TRUTH.TOTAL | TRUTH.TP | TRUTH.FN | QUERY.TOTAL | QUERY.FP | QUERY.UNK | FP.gt | FP.al | METRIC.Recall | METRIC.Precision | METRIC.Frac_NA | METRIC.F1_Score |
| INDEL | ALL | 8937 | 7839 | 1098 | 11812 | 343 | 3520 | 45 | 283 | 0.877140 | 0.958635 | 0.298002 | 0.916079 |
| INDEL | PASS | 8937 | 7550 | 1387 | 9971 | 283 | 1964 | 30 | 242 | 0.844803 | 0.964656 | 0.196971 | 0.900760 |
| SNP | ALL | 52494 | 52125 | 369 | 90092 | 582 | 37348 | 107 | 354 | 0.992971 | 0.988966 | 0.414554 | 0.990964 |
| SNP | PASS | 52494 | 46920 | 5574 | 48078 | 143 | 992
| | |
| NC_000020.11 | 62172726 | rs36106901 | G | A | |
| NC_000020.11 | 62172726 | 981031 | G | A | Conflicting_interpretations_of_pathogenicity |
| ------ | | | | | |
| NC_000020.11 | 63349782 | rs1044396 | G | A,C | |
| NC_000020.11 | 63349782 | 93427 | G | A | Benign |
| NC_000020.11 | 63349782 | 857384 | G | C | Conflicting_interpretations_of_pathogenicity |
| ------ | | | | | |
| NC_000020.11 | 63414925 | rs1801545 | G | A,C,T | |
| NC_000020.11 | 63414925 | 194284 | G | A | Conflicting_interpretations_of_pathogenicity |
| NC_000020.11 | 63414925 | 129337 | G | C | Benign |
| NC_000020.11 | 63414925 | 851545 | GG | CA | Uncertain_significance |
| ------ | | | | | |
On a donc plusieurs problèmes :
1. isec devrait fonctionner au moins sur
| NC_000020.11 | 25390747 | rs373200654 | G | C | |
| NC_000020.11 | 25390747 | 338000 | G | C | Conflicting_interpretations_of_pathogenicity |
On teste juste sur cette ligne
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter -i 'POS=25390747' clinvar_chr20.vcf.gz -o clinvar_test.vcf.gz
bcftools filter -i 'POS=25390747' dbSNP_common_chr20.vcf.gz -o dbSNP_test.vcf.gz
#+end_src
On retrouve bien la ligne dans l'intersection...
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter -i 'POS=25390747' clinvar_chr20.vcf.gz -o clinvar_test.vcf.gz
bcftools index dbSNP_test.vcf.gz dbSNP_test.vcf.gz
bcftools index dbSNP_test.vcf.gz clinvar_test.vcf.gz
bcftools isec dbSNP_test.vcf.gz clinvar_test.vcf.gz -p test
#+end_src
#+RESULTS:
2. isec ne semble pas fonctionner sur en cas d'ALT multiples
| NC_000020.11 | 32800145 | rs2424926 | C | G,T | |
| NC_000020.11 | 32800145 | 338173 | C | G | Benign |
| NC_000020.11 | 32800145 | 338174 | C | T | Conflicting_interpretations_of_pathogenicity |
| | | | | | |
3. s'il y a plusieurs variantions à une position, il faut bien vérifier que tous ne sont pas patho.
La version d'Alexis le fait bien
| NC_000020.11 | 3234173 | rs3827075 | T | A,C,G | |
| NC_000020.11 | 3234173 | 262001 | T | G | Conflicting_interpretations_of_pathogenicity |
| NC_000020.11 | 3234173 | 1072511 | T | TGGCGAAGC | Pathogenic |
| NC_000020.11 | 3234173 | 208613 | TGGCGAAGC | G | Pathogenic |
| NC_000020.11 | 3234173 | 1312 | TGGCGAAGC | T | Pathogenic |
****** DONE Voir si isec gère les multiallélique (chr20) : non, impossible de faire marcher
CLOSED: [2022-11-27 Sun 00:37]
******* DONE chr20 en prenant un patho clinvar aussi dans dbSNP
CLOSED: [2022-11-27 Sun 00:37]
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter dbSNP_common_chr20.vcf.gz -i 'POS=10652589' -o test_dbsnp.vcf.gz
bcftools filter clinvar_chr20.vcf.gz -i 'POS=10652589' -o test_clinvar.vcf.gz
bcftools index test_dbsnp.vcf.gz
bcftools index test_clinvar.vcf.gz
#+end_src
#+RESULTS:
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools isec test_dbsnp.vcf.gz test_clinvar.vcf.gz -p tmp
grep '^[^#]' tmp/0002.vcf
grep '^[^#]' tmp/0003.vcf
#+end_src
#+RESULTS:
Même en biallélique, ne fonctionne pas.
Testé en modifiant test_dbsnp !
Fonctionne avec un variant par ligne
****** DONE isec en coupant les sites multialléliques: non
CLOSED: [2022-11-27 Sun 00:37]
******* DONE Exemple simple ok
CLOSED: [2022-11-27 Sun 00:34]
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools filter -i 'POS=10652589' dbSNP_common_chr20.vcf.gz -o dbsnp_mwi.vcf.gz
bcftools filter -i 'POS=10652589' clinvar_chr20.vcf.gz -o clinvar_mwi.vcf.gz
bcftools index -f dbsnp_mwi.vcf.gz
bcftools index -f clinvar_mwi.vcf.gz
bcftools isec dbsnp_mwi.vcf.gz clinvar_mwi.vcf.gz -n=2
#+end_src
#+RESULTS:
Même en biallélique, ne fonctionne pas.
Chr 20
Avec les fichiers du teste précédent
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools norm -m -any dbsnp_mwi.vcf.gz -o dbsnp_mwi_norm.vcf.gz
bcftools index dbsnp_mwi_norm.vcf.gz
bcftools isec dbsnp_mwi_norm.vcf.gz clinvar_mwi.vcf.gz -n=2
#+end_src
#+RESULTS:
| NC_000020.11 | 10652589 | G | A | 11 |
| NC_000020.11 | 10652589 | G | C | 11 |
******* DONE Sur dbSNP chr20 non
CLOSED: [2023-10-07 Sat 17:57]
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools norm -m -any dbSNP_common_chr20 -o dbSNP_common_chr20_norm.vcf.gz
#+end_src
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools isec -i 'INFO/CLNSIG="Pathogenic"' dbSNP_common_chr20_norm.vcf.gz clinvar_chr20.vcf.gz -p tmp
#+end_src
#+RESULTS:
***** DONE Essai bedtools intersect
#+begin_src sh
bedtools intersect -a dbSNP_common.vcf.gz -b clinvar.vcf.gz
#+end_src
$ wc -l intersect.vcf
220206 intersect.vcf
** TODO Dépendences avec Nix
*** DONE GATK
CLOSED: [2022-10-21 Fri 21:59]
*** WAIT BioDBHTS
Contribuer pull request
*** DONE BioExtAlign
CLOSED: [2022-10-22 Sat 00:38]
*** DONE BioBigFile
CLOSED: [2023-11-30 Thu 21:57]
Revoir si on peut utliser kent dernière version
Contribuer pull request
*** HOLD rtg-tools
Convertir clinvar NC
*** DONE simuscop
CLOSED: [2022-12-30 Fri 22:31]
*** DONE Spip
CLOSED: [2022-12-04 Sun 12:49]
Pas de pull request
*** DONE R + packages
CLOSED: [2022-11-19 Sat 21:05]
*** TODO hap.py
https://github.com/Illumina/hap.py
**** DONE Version sans rtgtools avec python 3
CLOSED: [2023-02-02 Thu 22:15]
Procédure pour tester
#+begin_src
nix develop .#hap-py
$ genericBuild
#+end_src
1. Supprimer l’appel à make_dependencies dans cmakelist.txt : on peut tout installer avec nix
2. Patch Roc.cpp pour avoir numeric_limits ( error: 'numeric_limits' is not a member of 'std')
3. ajout de flags de link (essai, error)
set(ZLIB_LIBRARIES -lz -lbz2 -lcurl -lcrypto -llzma)
4. Changer les appels à print en print() dans le code python et suppression de quelques import
[nix-shell:~/source]$ sed -i.orig 's/print \"\(.*\)"/print(\1)/' src/python/*.py
**** DONE Sérialiser json pour écrire données de sorties
CLOSED: [2023-02-17 Fri 19:25]
**** DONE Tester sur example
CLOSED: [2023-02-04 Sat 00:25]
#+begin_src sh
$ cd hap.py
$ ../result/bin/hap.py example/happy/PG_NA12878_chr21.vcf.gz example/happy/NA12878_chr21.vcf.gz -f example/happy/PG_Conf_chr21.bed.gz -o test -r example/chr21.fa
#+end_src
#+RESULTS:
| Type | Filter | TRUTH.TOTAL | TRUTH.TP | TRUTH.FN | QUERY.TOTAL | QUERY.FP | QUERY.UNK | FP.gt | FP.al | METRIC.Recall | METRIC.Precision | METRIC.Frac_NA | METRIC.F1_Score |
| INDEL | ALL | 8937 | 7839 | 1098 | 11812 | 343 | 3520 | 45 | 283 | 0.877140 | 0.958635 | 0.298002 | 0.916079 |
| INDEL | PASS | 8937 | 7550 | 1387 | 9971 | 283 | 1964 | 30 | 242 | 0.844803 | 0.964656 | 0.196971 | 0.900760 |
| SNP | ALL | 52494 | 52125 | 369 | 90092 | 582 | 37348 | 107 | 354 | 0.992971 | 0.988966 | 0.414554 | 0.990964 |
| SNP | PASS | 52494 | 46920 | 5574 | 48078 | 143 | 992
0.120397 0.844947 3.017198 2.85705 5.560099 2.114633
SNP PASS 15883 15479 404 23597 5277 2841 46 44 0.974564 0.745760 0.120397 0.844947 3.017198 2.85705 5.560099 2.114633
******* DONE Vérifier qu'il ne reste plus de filtre autre que PASS
CLOSED: [2023-07-08 Sat 15:19]
#+begin_src
$ zgrep -c 'PASS' HG001_GRCh38_1_22_v4_lifted_merged.vcf.gz
3730505
$ zgrep -c '^chr' HG001_GRCh38_1_22_v4_lifted_merged.vcf.gz
3730506
#+end_src
****** TODO 1/4 SNP manquant ?
******* DONE Regarder avec Julia si ce sont vraiment des FP: 61/5277 qui ne le sont pas
CLOSED: [2023-07-09 Sun 12:09]
******* DONE Examiner les FP
CLOSED: [2023-07-30 Sun 22:05]
******* DONE Tester un FP
CLOSED: [2023-07-30 Sun 22:05]
2 │ chr1 608765 A G ./.:.:.:.:NOCALL:nocall:. 1/1:FP:.:ti:SNP:homalt:188
liftDown UCSC: rien en GIAB : vrai FP
3 │ chr1 762943 A G ./.:.:.:.:NOCALL:nocall:. 1/1:FP:.:ti:SNP:homalt:287
4 │ chr1 762945 A T ./.:.:.:.:NOCALL:nocall:. 1/1:FP:.:tv:SNP:homalt:287
Remaniements complexes ? Pas dans le gène en HG38
******* DONE La plupart des FP (4705/5566) sont homozygotes: erreur de référence ?
CLOSED: [2023-07-12 Wed 21:10] SCHEDULED: <2023-07-09 Sun>
Sur les 2 premiers variants, ils montrent en fait la différence entre T2T et GRCh38
Erreur à l'alignement ?
******** KILL relancer
l'alignement
CLOSED: [2023-07-09 Sun 17:36]
******** DONE vérifier reads identiques hg38 et T2T: oui
CLOSED: [2023-07-09 Sun 16:36]
T2T CHR1608765
38 chr1:1180168-1180168 (
SRR14724513.24448214
SRR14724513.24448214
******* DONE Vérifier quelques variants sur IGV
CLOSED: [2023-07-09 Sun 17:36]
******* KILL Répartition des FP : cluster ?
CLOSED: [2023-07-09 Sun 17:36]
****** DONE Examiner les FP restant après correction selon séquence de référence
CLOSED: [2023-08-12 Sat 15:57]
****** HOLD Examiner les variants supprimé
****** TODO Enlever les FP qui correspondent à un changement dans le génome
******* Condition:
- pas de variation à la position en GRCh38
- variantion homozygote
- la varation en T2T correspond au changement de pair de base GRC38 -> T2T
pour les SNP:
alt_T2T[i] = DNA_GRC38[j]
avec i la position en T2T et j la position en GRCh38
Note: définir un ID n'est pas correct car les variants peuvent être modifié par happy !
******* Idée
- Pour chaque FP, c'est un "faux" FP si
- REF en hg38 == ALT en T2T
- et REF en hg38 != REF en T2T
- et variant homozygote
Comment obtenir les séquences de réferences ?
1. liftover
2. blat sur la séquence autour du variant
3. identifier quelques reads contenant le variant et regarder leur aligneement en hg38
Après discussion avec Alexis: solution 3
******* Algorithme
1. Extraire les coordonnées en T2T des faux positifs *homozygote*
2. Pour chaque faux positif
1. lister 10 reads contenant le variant
2. pour chacun de ces reads, récupérer la séquence en T2T et GRCh38 via le nom du read dans le bam
3. si la séquence en T2T modifiée par le variant est "identique" à celle en GRCh38, alors on ignore ce faux positif
Note: on ignore les reads qui ont changé de chromosome entre les version
******* DONE Résultat préliminaire
CLOSED: [2023-07-23 Sun 14:30]
cf [[file:~/roam/research/bisonex/code/giab/giab-corrected.csv][script julia]]
3498 faux positifs en moins, soit 0.89 sensibilité
julia> tp=15479
julia> fp=5277
julia> tp/(tp+fp)
0.7457602620928888
julia> tp/(tp+(fp-3498))
0.8969173716537258
On est toujours en dessous des 97%
******* HOLD Corriger proprement VCF ou résultats Happy
******* TODO Adapter pour gérer plusieurs variants par read
****** DONE Méthodologie du pangenome
CLOSED: [2023-10-03 Tue 21:28]
Voir biblio[cite:@liao2023] mais ont aligné sur GRCH38
******* DONE Mail alexis
CLOSED: [2023-10-03 Tue 21:28]
****** DONE Méthodologie T2T
CLOSED: [2023-10-16 Mon 19:42]
Mail alexis
SCHEDULED: <2023-10-04 Wed>
***** TODO Rendre simplement le nombre de vrais positifs
SCHEDULED: <2023-12-02 Sat>
***** KILL Mail Yannis
CLOSED: [2023-07-08 Sat 10:44]
***** DONE Mail GIAB pour version T2T
CLOSED: [2023-07-07 Fri 18:37]
**** KILL HG002 :hg002:T2T:
CLOSED: [2023-11-26 Sun 12:30]
**** KILL HG003 :hg003:T2T:
CLOSED: [2023-11-26 Sun 12:30]
**** KILL HG004 :hg004:T2T:
CLOSED: [2023-11-26 Sun 12:30]
**** DONE Plot : ashkenazim trio :hg38:
CLOSED: [2023-07-30 Sun 16:49] SCHEDULED: <2023-07-30 Sun 15:00>
:LOGBOOK:
CLOCK: [2023-07-30 Sun 16:06]--[2023-07-30 Sun 16:35] => 0:29
CLOCK: [2023-07-30 Sun 15:39]--[2023-07-30 Sun 15:40] => 0:01
:END:
/Entered on/ [2023-04-16 Sun 17:29]
Refaire résultats
**** DONE Mail Paul sur les résultat ashkenazim +/- centogene
CLOSED: [2023-08-06 Sun 20:24] SCHEDULED: <2023-08-06 Sun>
**** DONE Relancer comparaison GIAB avec GATK 4.4.0
CLOSED: [2023-08-12 Sat 15:55]
/Entered on/ [2023-08-03 Thu 12:42]
**** TODO Re-télécharger proprement dans pipeline dédiés
[[*Résumé][Résumé]]
Cf [[*Validation : Quelles données de référence ?][Validation : Quelles données de référence ?]]
https://medium.com/dnanexus/benchmarking-state-of-the-art-secondary-variant-calling-pipelines-5472ca6bace7
Source:
https://trace.ncbi.nlm.nih.gov/Traces/index.html?view=study&acc=SRP047086
https://zenodo.org/records/3597727
Selon https://github.com/genome-in-a-bottle/giab_data_indexes
***** TODO HG001 :hg001:
SCHEDULED: <2023-11-29 Wed>
****** TODO Avec données en hg38
SCHEDULED: <2023-11-29 Wed>
[[*Résumé][Résumé]]
****** TODO Avec données en hg19
SCHEDULED: <2023-12-03 Sun>
Utiliser crossmap ! https://crossmap.readthedocs.io/en/latest/ (inspiré de [[https://github.com/bcbio/bcbio_validation_workflows/blob/master/giab-exome/input/get_data.sh][bcbio]]
pour vérifier
***** TODO HG002 :hg002:
SCHEDULED: <2023-11-30 Thu>
***** TODO HG003 :hg003:
SCHEDULED: <2023-11-30 Thu>
***** TODO HG004 :hg001:
SCHEDULED: <2023-11-30 Thu>
**** TODO Refaire les analyses pour avoir meilleurs résultats
SCHEDULED: <2023-12-03 Sun>
On veut les résultats de https://medium.com/dnanexus/benchmarking-state-of-the-art-secondary-variant-calling-pipelines-5472ca6bace7
***** TODO hap.py avec conda
SCHEDULED: <2023-12-03 Sun>
***** TODO rtgveval
SCHEDULED: <2023-12-03 Sun>
***** TODO Relancer
SCHEDULED: <2023-12-03 Sun>
*** TODO Platinum genome :platinum:
https://emea.illumina.com/platinumgenomes.html
**** TODO Tester sur la zone couverte par l'exome centogène
SCHEDULED: <2023-12-02 Sat>
*** DONE Séquencer NA12878 :cento:hg001:
CLOSED: [2023-10-07 Sat 17:59]
Discussion avec Paul : sous-traitant ne nous donnera pas les données, il faut commander l'ADN
**** DONE ADN commandé
CLOSED: [2023-06-30 Fri 22:29]
**** DONE Sauvegarder les données brutes
CLOSED: [2023-07-30 Sun 14:22] SCHEDULED: <2023-07-19 Wed>
K, scality, S
**** KILL Récupérer le fichier de capture
CLOSED: [2023-07-30 Sun 14:25] SCHEDULED: <2023-07-23 Sun>
Candidats donnés dans publication https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8354858/
#+begin_quote
In short, the Nextera Rapid Capture Exome Kit (Illumina, San Diego, CA), the SureSelect Human All Exon kit (Agilent, Santa Clara, CA) or the Twist Human Core Exome was used for enrichment, and a Nextseq500, HiSeq4000, or Novoseq 6000 (Illumina) instrument was used for the actual sequencing, with the average coverage targeted to at least 100× or at least 98% of the target DNA covered 20×.
#+end_quote
Par défaut, on utilisera https://www.twistbioscience.com/products/ngs/alliance-panels#tab-3
ANnonce récente pour nouveau panel Twist : https://www.centogene.com/news-events/news/newsdetails/twist-bioscience-and-centogene-launch-three-panels-to-advance-rare-disease-and-hereditary-cancer-research-and-support-diagnostics
Masi pas de fichier BED
***** DONE Mail centogène
CLOSED
0.120397 0.844947 3.017198 2.85705 5.560099 2.114633
SNP PASS 15883 15479 404 23597 5277 2841 46 44 0.974564 0.745760 0.120397 0.844947 3.017198 2.85705 5.560099 2.114633
******* DONE Vérifier qu'il ne reste plus de filtre autre que PASS
CLOSED: [2023-07-08 Sat 15:19]
#+begin_src
$ zgrep -c 'PASS' HG001_GRCh38_1_22_v4_lifted_merged.vcf.gz
3730505
$ zgrep -c '^chr' HG001_GRCh38_1_22_v4_lifted_merged.vcf.gz
3730506
#+end_src
****** TODO 1/4 SNP manquant ?
******* DONE Regarder avec Julia si ce sont vraiment des FP: 61/5277 qui ne le sont pas
CLOSED: [2023-07-09 Sun 12:09]
******* DONE Examiner les FP
CLOSED: [2023-07-30 Sun 22:05]
******* DONE Tester un FP
CLOSED: [2023-07-30 Sun 22:05]
2 │ chr1 608765 A G ./.:.:.:.:NOCALL:nocall:. 1/1:FP:.:ti:SNP:homalt:188
liftDown UCSC: rien en GIAB : vrai FP
3 │ chr1 762943 A G ./.:.:.:.:NOCALL:nocall:. 1/1:FP:.:ti:SNP:homalt:287
4 │ chr1 762945 A T ./.:.:.:.:NOCALL:nocall:. 1/1:FP:.:tv:SNP:homalt:287
Remaniements complexes ? Pas dans le gène en HG38
******* DONE La plupart des FP (4705/5566) sont homozygotes: erreur de référence ?
CLOSED: [2023-07-12 Wed 21:10] SCHEDULED: <2023-07-09 Sun>
Sur les 2 premiers variants, ils montrent en fait la différence entre T2T et GRCh38
Erreur à l'alignement ?
******** KILL relancer l'alignement
CLOSED: [2023-07-09 Sun 17:36]
******** DONE vérifier reads identiques hg38 et T2T: oui
CLOSED: [2023-07-09 Sun 16:36]
T2T CHR1608765
38 chr1:1180168-1180168 (
SRR14724513.24448214
SRR14724513.24448214
******* DONE Vérifier quelques variants sur IGV
CLOSED: [2023-07-09 Sun 17:36]
******* KILL Répartition des FP : cluster ?
CLOSED: [2023-07-09 Sun 17:36]
****** DONE Examiner les FP restant après correction selon séquence de référence
CLOSED: [2023-08-12 Sat 15:57]
****** HOLD Examiner les variants supprimé
****** TODO Enlever les FP qui correspondent à un changement dans le génome
******* Condition:
- pas de variation à la position en GRCh38
- variantion homozygote
- la varation en T2T correspond au changement de pair de base GRC38 -> T2T
pour les SNP:
alt_T2T[i] = DNA_GRC38[j]
avec i la position en T2T et j la position en GRCh38
Note: définir un ID n'est pas correct car les variants peuvent être modifié par happy !
******* Idée
- Pour chaque FP, c'est un "faux" FP si
- REF en hg38 == ALT en T2T
- et REF en hg38 != REF en T2T
- et variant homozygote
Comment obtenir les séquences de réferences ?
1. liftover
2. blat sur la séquence autour du variant
3. identifier quelques reads contenant le variant et regarder leur aligneement en hg38
Après discussion avec Alexis: solution 3
******* Algorithme
1. Extraire les coordonnées en T2T des faux positifs *homozygote*
2. Pour chaque faux positif
1. lister 10 reads contenant le variant
2. pour chacun de ces reads, récupérer la séquence en T2T et GRCh38 via le nom du read dans le bam
3. si la séquence en T2T modifiée par le variant est "identique" à celle en GRCh38, alors on ignore ce faux positif
Note: on ignore les reads qui ont changé de chromosome entre les version
******* DONE Résultat préliminaire
CLOSED: [2023-07-23 Sun 14:30]
cf [[file:~/roam/research/bisonex/code/giab/giab-corrected.csv][script julia]]
3498 faux positifs en moins, soit 0.89 sensibilité
julia> tp=15479
julia> fp=5277
julia> tp/(tp+fp)
0.7457602620928888
julia> tp/(tp+(fp-3498))
0.8969173716537258
On est toujours en dessous des 97%
******* HOLD Corriger proprement VCF ou résultats Happy
******* TODO Adapter pour gérer plusieurs variants par read
****** DONE Méthodologie du pangenome
CLOSED: [2023-10-03 Tue 21:28]
Voir biblio[cite:@liao2023] mais ont aligné sur GRCH38
******* DONE Mail alexis
CLOSED: [2023-10-03 Tue 21:28]
****** DONE Méthodologie T2T
CLOSED: [2023-10-16 Mon 19:42]
Mail alexis
SCHEDULED: <2023-10-04 Wed>
***** TODO Rendre simplement le nombre de vrais positifs
SCHEDULED: <2023-12-02 Sat>
***** KILL Mail Yannis
CLOSED: [2023-07-08 Sat 10:44]
***** DONE Mail GIAB pour version T2T
CLOSED: [2023-07-07 Fri 18:37]
**** KILL HG002 :hg002:T2T:
CLOSED: [2023-11-26 Sun 12:30]
**** KILL HG003 :hg003:T2T:
CLOSED: [2023-11-26 Sun 12:30]
**** KILL HG004 :hg004:T2T:
CLOSED: [2023-11-26 Sun 12:30]
**** DONE Plot : ashkenazim trio :hg38:
CLOSED: [2023-07-30 Sun 16:49] SCHEDULED: <2023-07-30 Sun 15:00>
:LOGBOOK:
CLOCK: [2023-07-30 Sun 16:06]--[2023-07-30 Sun 16:35] => 0:29
CLOCK: [2023-07-30 Sun 15:39]--[2023-07-30 Sun 15:40] => 0:01
:END:
/Entered on/ [2023-04-16 Sun 17:29]
Refaire résultats
**** DONE Mail Paul sur les résultat ashkenazim +/- centogene
CLOSED: [2023-08-06 Sun 20:24] SCHEDULED: <2023-08-06 Sun>
**** DONE Relancer comparaison GIAB avec GATK 4.4.0
CLOSED: [2023-08-12 Sat 15:55]
/Entered on/ [2023-08-03 Thu 12:42]
**** TODO Re-télécharger proprement dans pipeline dédiés
[[*Résumé][Résumé]]
Cf [[*Validation : Quelles données de référence ?][Validation : Quelles données de référence ?]]
https://medium.com/dnanexus/benchmarking-state-of-the-art-secondary-variant-calling-pipelines-5472ca6bace7
Source:
https://trace.ncbi.nlm.nih.gov/Traces/index.html?view=study&acc=SRP047086
https://zenodo.org/records/3597727
Selon https://github.com/genome-in-a-bottle/giab_data_indexes
***** TODO HG001 :hg001:
SCHEDULED: <2023-11-29 Wed>
****** TODO Avec données en hg38
SCHEDULED: <2023-11-29 Wed>
[[*Résumé][Résumé]]
****** TODO Avec données en hg19
SCHEDULED: <2023-12-03 Sun>
Utiliser crossmap ! https://crossmap.readthedocs.io/en/latest/ (inspiré de [[https://github.com/bcbio/bcbio_validation_workflows/blob/master/giab-exome/input/get_data.sh][bcbio]]
pour vérifier
***** TODO HG002 :hg002:
SCHEDULED: <2023-12-04 Mon>
***** TODO HG003 :hg003:
SCHEDULED: <2023-12-04 Mon>
***** TODO HG004 :hg001:
SCHEDULED: <2023-12-04 Mon>
**** TODO Refaire les analyses pour avoir meilleurs résultats
SCHEDULED: <2023-12-03 Sun>
On veut les résultats de https://medium.com/dnanexus/benchmarking-state-of-the-art-secondary-variant-calling-pipelines-5472ca6bace7
***** TODO hap.py avec conda
SCHEDULED: <2023-12-03 Sun>
***** TODO rtgveval
SCHEDULED: <2023-12-03 Sun>
***** TODO Relancer
SCHEDULED: <2023-12-03 Sun>
*** TODO Platinum genome :platinum:
https://emea.illumina.com/platinumgenomes.html
**** TODO Tester sur la zone couverte par l'exome centogène
SCHEDULED: <2023-12-02 Sat>
*** DONE Séquencer NA12878 :cento:hg001:
CLOSED: [2023-10-07 Sat 17:59]
Discussion avec Paul : sous-traitant ne nous donnera pas les données, il faut commander l'ADN
**** DONE ADN commandé
CLOSED: [2023-06-30 Fri 22:29]
**** DONE Sauvegarder les données brutes
CLOSED: [2023-07-30 Sun 14:22] SCHEDULED: <2023-07-19 Wed>
K, scality, S
**** KILL Récupérer le fichier de capture
CLOSED: [2023-07-30 Sun 14:25] SCHEDULED: <2023-07-23 Sun>
Candidats donnés dans publication https://www.ncbi.nlm.nih.gov/pmc/articles/PMC8354858/
#+begin_quote
In short, the Nextera Rapid Capture Exome Kit (Illumina, San Diego, CA), the SureSelect Human All Exon kit (Agilent, Santa Clara, CA) or the Twist Human Core Exome was used for enrichment, and a Nextseq500, HiSeq4000, or Novoseq 6000 (Illumina) instrument was used for the actual sequencing, with the average coverage targeted to at least 100× or at least 98% of the target DNA covered 20×.
#+end_quote
Par défaut, on utilisera https://www.twistbioscience.com/products/ngs/alliance-panels#tab-3
ANnonce récente pour nouveau panel Twist : https://www.centogene.com/news-events/news/newsdetails/twist-bioscience-and-centogene-launch-three-panels-to-advance-rare-disease-and-hereditary-cancer-research-and-support-diagnostics
Masi pas de fichier BED
***** DONE Mail centogène
CLOSED
er
open snvs-cento-sanger.csv | select chrom pos | insert pos2 {$in.pos } | to csv --separator="\t" | save snvs-cento-sanger.bed -f
#+end_src
2. Liftover avec UCSC (en ligne)
NB: vérifié sur le premier résultat en cherche le read contenant le variant (samtools view -r puis samtools view | grep en T2T) et avec l'aide d'IGV, on a un variant qui correspond en
chr1:10757746
3. En supposant que l'ordre des variants n'a pas changé, on ajoute simplement REF et ALT avec annotateLifted.jl
Annotation spip *très lente* : 1h13 !
Résultat:
2×3 DataFrame
Row │ variant meanQual depth
│ String Float64 Int64
─────┼──────────────────────────────────────
1 │ chr12:g.13594572 60.0 1
2 │ chr17:g.10204026 60.0 1
144 found over 146
filter depth : another 0 missed variants
filter poly : another 0 missed variants
filter vep : another 0 missed variants
Et on a trop de variants en sortie (7330 !)
**** DONE Mail Paul avec résultats filtre en T2T + nouveau schéma
CLOSED: [2023-09-17 Sun 23:15] SCHEDULED: <2023-09-17 Sun>
** TODO Medically relevant genes
SCHEDULED: <2023-11-30 Thu>
/Entered on/ [2023-10-18 Wed 22:37]
** TODO HG002 en T2T
/Entered on/ [2023-11-25 Sat 17:58]
https://github.com/marbl/HG002
*** TODO Tester les benchmark préliminaires
SCHEDULED: <2023-12-26 Tue>
https://ftp-trace.ncbi.nlm.nih.gov/ReferenceSamples/giab/data/AshkenazimTrio/analysis/NIST_HG002_DraftBenchmark_defrabbV0.011-20230725/
* Ré-interprétation :reanalysis:
** DONE Lancer tests sur données brutes [225/250] <(samples.csv)> <(runs.waiting)>
CLOSED: [2023-10-14 Sat 11:58] SCHEDULED: <2023-10-08 Sun>
- [X] 100222_63015289
- [X] 1600304839_63051311
- [X] 1900007827_62913191
- [X] 1900398899_62999500
- [X] 1900486799_62913197
- [X] 2100422923_62952677
- [X] 2100458888_62933047
- [X] 2100601558_62903840
- [X] 2100609288_62905768
- [X] 2100609501_62905776
- [X] 2100614493_62951074
- [X] 2100622566_62908067
- [X] 2100622601_62908060
- [X] 2100622705_62908063
- [X] 2100640027_62911936
- [X] 2100645285_62913212
- [X] 2100661411_62914081
- [X] 2100661462_62914086
- [X] 2100708257_62921596
- [X] 2100738732_62926501
- [X] 2100738850_62926509
- [X] 2100746751_62926505
- [X] 2100746797_62926506
- [X] 2100782349_62931722
- [X] 2100782416_62931561
- [X] 2100782559_62931718
- [X] 2100799204_62934768
- [X] 2200010202_62940284
- [X] 2200023600_62940631
- [X] 2200024348_62999591
- [X] 2200027505_62942457
- [X] 2200038776_62943412
- [X] 2200041919_62943405
- [X] 2200088014_62951326
- [X] 2200146652_62959388
- [X] 2200151850_62960953
- [X] 2200160014_62959475
- [X] 2200160070_62959478
- [X] 2200201368_62967471
- [X] 2200201400_62967470
- [X] 2200265558_62976332
- [X] 2200265605_62976401
- [X] 2200267046_62975192
- [X] 2200273878_62999530
- [X] 2200279708_62977002
- [X] 2200284408_62979102
- [X] 2200293987_62979116
- [X] 2200294359_62979118
- [X] 2200306299_62982217
- [X] 2200306539_62982193
- [X] 220030671_62982211
- [X] 2200307058_62982231
- [X] 2200307108_62982196
- [X] 2200307136_62982221
- [X] 2200307199_62982239
- [X] 2200307230_62982234
- [X] 2200307262_62982219
- [X] 2200307297_62982227
- [X] 2200324510_62985453
- [X] 2200324549_62985478
- [X] 2200324573_62985445
- [X] 2200324594_62985467
- [X] 2200324606_62985463
- [X] 2200324614_62985459
- [X] 2200338306_62985430
- [X] 2200343880_62989407
- [X] 2200343910_62989460
- [X] 2200343938_62989451
- [X] 2200343966_62989456
- [X] 2200343993_62989440
- [X] 2200344013_62989464
- [X] 2200349749_62989465
- [X] 2200363462_62988848
- [X] 2200377880_62991993
- [X] 2200378032_62991991
- [X] 2200383996_62993828
- [X] 2200384015_62993796
- [X] 2200384046_62993822
- [X] 2200384117_62993808
- [X] 2200384187_62993825
- [X] 2200384231_62992898
- [X] 2200385658_63060260
- [X] 2200394260_62994732
- [X] 2200395817_62994742
- [X] 2200396731_62994737
- [X] 2200424073_62999579
- [X] 2200424207_62999632
- [X] 2200426178_62999630
- [X] 2200426243_62999635
- [X] 2200426466_62999605
- [X] 2200426642_62999627
- [X] 2200427406_62999649
- [X] 2200427512_62999639
- [X] 2200428953_62999572
- [X] 2200428981_62999600
- [X] 2200428999_62999592
- [X] 2200441970_63000868
- [X] 2200441989_63000882
- [X] 2200442135_63000864
- [X] 2200442216_63000886
- [X] 2200442257_63000951
- [X] 2200451801_63003573
- [X] 2200451862_63004218
- [X] 2200451894_63004210
- [X] 2200456165_63051294
- [X] 2200459865_63004933
- [X] 2200459968_63004937
- [X] 2200460073_63004943
- [X] 2200460121_63004684
- [X] 2200467051_63003856
- [X] 2200467225_63004940
- [X] 2200467261_63004930
- [X] 2200467338_63004925
- [X] 2200470099_63004485
- [X] 2200470142_63004480
- [X] 2200471780_63004362
- [X] 2200480910_63006466
- [X] 2200495073_63010427
- [X] 2200495510_63009152
- [X] 2200508677_63060252
- [X] 2200510531_63012582
- [X] 2200510628_63012549
- [X] 2200510657_63012554
- [X] 2200511249_63012533
- [X] 2200511274_63012586
- [X] 2200517952_63060399
- [X] 2200519525_63060439
- [X] 2200524009_63014044
- [X] 2200524609_63014046
- [X] 2200524616_63014048
- [X] 2200533429_63060425
- [X] 2200539735_63060406
- [X] 2200549908_63019339
- [X] 2200549965_63019349
- [X] 2200550414_63019357
- [X] 2200550471_63020031
- [X] 2200550490_63019351
- [X] 2200550505_63019340
- [X] 2200555565_63018614
- [X] 2200559438_63020029
- [X] 2200559682_63020030
- [X] 2200559713_63019623
- [X] 2200559739_63019626
- [X] 2200569969_63019991
- [X] 2200570001_63021580
- [X] 2200570025_63021490
- [X] 2200570035_63021491
- [X] 2200570042_63021493
- [X] 2200570050_63021494
- [X] 2200579897_63024910
- [X] 2200583995_63024866
- [X] 2200584035_63024905
- [X] 2200584069_63024888
- [X] 2200584126_63024810
- [X] 2200589507_63026712
- [X] 2200597365_63027994
- [X] 2200597480_63027988
- [X] 2200597752_63026853
- [X] 2200597778_63027992
- [X] 22005977_63026903
- [X] 2200609031_63026527
- [X] 2200614198_63113928
- [X] 2200620372_63030821
- [X] 2200620442_63030810
- [X] 2200620498_63030816
- [X] 2200620628_63031031
- [X] 2200622310_63030984
- [X] 2200622355_63030956
- [X] 2200625369_63028699
- [X] 2200625410_63028697
- [X] 2200625536_63028694
- [X] 2200630189_63030665
- [X] 2200635149_63033182
- [X] 2200644544_63037731
- [X] 2200644594_63037725
- [X] 2200650089_63038093
- [X] 2200666292_63076568
- [X] 2200669188_63036688
- [X] 2200669320_63040259
- [X] 2200669383_63040254
- [X] 2200669414_63040257
- [X] 2200669446_63040251
- [X] 2200680342_63105271
- [X] 2200694535_63042853
- [X] 2200694789_63042862
- [X] 2200694858_63042702
- [X] 2200694917_63042696
- [X] 2200699290_63043047
- [X] 2200699345_63040238
- [X] 2200699383_63043050
- [X] 2200699412_63040731
- [X] 220071551_63048935
- [X] 2200731515_63048963
- [X] 2200748145_63051198
- [X] 2200748171_63051213
- [X] 2200751046_63051249
- [X] 2200751101_63051234
- [X] 2200766471_63054590
- [X] 2200767731_63054595
- [X] 2200767822_63054464
- [X] 2200775505_63060410
- [X] 2200850441_63019345
- [X] 220597589_63026879
- [X] 2300003253_63060430
- [X] 2300005679_63060370
- [X] 2300009914_63060390
- [X] 2300028784_63060001
- [X] 2300036815_63063357
- [X] 2300055382_63061874
- [X] 2300055421_63061871
- [X] 2300055440_63061880
- [X] 230006894_63064950
- [X] 2300071111_
63070356
- [X] 2300083434_
63071675
- [X] 2300103609_63076239
- [X] 2300104572_63076232
- [X] 2300109602_63076765
- [X] 2300109665_63076770
- [X] 2300119721_63078732
- [X] 2300137773_63078133
- [X] 2300137834_63078123
- [X] 2300167821_63086183
- [X] 2300172698_63113453
- [X] 2300188216_63090609
- [X] 2300188281_63090632
- [ ] 2300188800_63090616
- [ ] 2300193193645_63090623
- [ ] 2300193668_63090611
- [ ] 2300195426_63090608
- [ ] 2300201017_63089636
- [ ] 2300227479_63098330
- [ ] 2300232688_63130821
- [ ] 2300292749_63109239
- [ ] 230029277_63109247
- [ ] 2300294712_63109236
- [ ] 2300308032_63111581
- [ ] 2300323537_63114209
- [ ] 23003
er
open snvs-cento-sanger.csv | select chrom pos | insert pos2 {$in.pos } | to csv --separator="\t" | save snvs-cento-sanger.bed -f
#+end_src
2. Liftover avec UCSC (en ligne)
NB: vérifié sur le premier résultat en cherche le read contenant le variant (samtools view -r puis samtools view | grep en T2T) et avec l'aide d'IGV, on a un variant qui correspond en
chr1:10757746
3. En supposant que l'ordre des variants n'a pas changé, on ajoute simplement REF et ALT avec annotateLifted.jl
Annotation spip *très lente* : 1h13 !
Résultat:
2×3 DataFrame
Row │ variant meanQual depth
│ String Float64 Int64
─────┼──────────────────────────────────────
1 │ chr12:g.13594572 60.0 1
2 │ chr17:g.10204026 60.0 1
144 found over 146
filter depth : another 0 missed variants
filter poly : another 0 missed variants
filter vep : another 0 missed variants
Et on a trop de variants en sortie (7330 !)
**** DONE Mail Paul avec résultats filtre en T2T + nouveau schéma
CLOSED: [2023-09-17 Sun 23:15] SCHEDULED: <2023-09-17 Sun>
** TODO Medically relevant genes
SCHEDULED: <2023-12-07 Thu>
/Entered on/ [2023-10-18 Wed 22:37]
** TODO HG002 en T2T
/Entered on/ [2023-11-25 Sat 17:58]
https://github.com/marbl/HG002
*** TODO Tester les benchmark préliminaires
SCHEDULED: <2023-12-26 Tue>
https://ftp-trace.ncbi.nlm.nih.gov/ReferenceSamples/giab/data/AshkenazimTrio/analysis/NIST_HG002_DraftBenchmark_defrabbV0.011-20230725/
* Ré-interprétation :reanalysis:
** DONE Lancer tests sur données brutes [225/250] <(samples.csv)> <(runs.waiting)>
CLOSED: [2023-10-14 Sat 11:58] SCHEDULED: <2023-10-08 Sun>
- [X] 100222_63015289
- [X] 1600304839_63051311
- [X] 1900007827_62913191
- [X] 1900398899_62999500
- [X] 1900486799_62913197
- [X] 2100422923_62952677
- [X] 2100458888_62933047
- [X] 2100601558_62903840
- [X] 2100609288_62905768
- [X] 2100609501_62905776
- [X] 2100614493_62951074
- [X] 2100622566_62908067
- [X] 2100622601_62908060
- [X] 2100622705_62908063
- [X] 2100640027_62911936
- [X] 2100645285_62913212
- [X] 2100661411_62914081
- [X] 2100661462_62914086
- [X] 2100708257_62921596
- [X] 2100738732_62926501
- [X] 2100738850_62926509
- [X] 2100746751_62926505
- [X] 2100746797_62926506
- [X] 2100782349_62931722
- [X] 2100782416_62931561
- [X] 2100782559_62931718
- [X] 2100799204_62934768
- [X] 2200010202_62940284
- [X] 2200023600_62940631
- [X] 2200024348_62999591
- [X] 2200027505_62942457
- [X] 2200038776_62943412
- [X] 2200041919_62943405
- [X] 2200088014_62951326
- [X] 2200146652_62959388
- [X] 2200151850_62960953
- [X] 2200160014_62959475
- [X] 2200160070_62959478
- [X] 2200201368_62967471
- [X] 2200201400_62967470
- [X] 2200265558_62976332
- [X] 2200265605_62976401
- [X] 2200267046_62975192
- [X] 2200273878_62999530
- [X] 2200279708_62977002
- [X] 2200284408_62979102
- [X] 2200293987_62979116
- [X] 2200294359_62979118
- [X] 2200306299_62982217
- [X] 2200306539_62982193
- [X] 220030671_62982211
- [X] 2200307058_62982231
- [X] 2200307108_62982196
- [X] 2200307136_62982221
- [X] 2200307199_62982239
- [X] 2200307230_62982234
- [X] 2200307262_62982219
- [X] 2200307297_62982227
- [X] 2200324510_62985453
- [X] 2200324549_62985478
- [X] 2200324573_62985445
- [X] 2200324594_62985467
- [X] 2200324606_62985463
- [X] 2200324614_62985459
- [X] 2200338306_62985430
- [X] 2200343880_62989407
- [X] 2200343910_62989460
- [X] 2200343938_62989451
- [X] 2200343966_62989456
- [X] 2200343993_62989440
- [X] 2200344013_62989464
- [X] 2200349749_62989465
- [X] 2200363462_62988848
- [X] 2200377880_62991993
- [X] 2200378032_62991991
- [X] 2200383996_62993828
- [X] 2200384015_62993796
- [X] 2200384046_62993822
- [X] 2200384117_62993808
- [X] 2200384187_62993825
- [X] 2200384231_62992898
- [X] 2200385658_63060260
- [X] 2200394260_62994732
- [X] 2200395817_62994742
- [X] 2200396731_62994737
- [X] 2200424073_62999579
- [X] 2200424207_62999632
- [X] 2200426178_62999630
- [X] 2200426243_62999635
- [X] 2200426466_62999605
- [X] 2200426642_62999627
- [X] 2200427406_62999649
- [X] 2200427512_62999639
- [X] 2200428953_62999572
- [X] 2200428981_62999600
- [X] 2200428999_62999592
- [X] 2200441970_63000868
- [X] 2200441989_63000882
- [X] 2200442135_63000864
- [X] 2200442216_63000886
- [X] 2200442257_63000951
- [X] 2200451801_63003573
- [X] 2200451862_63004218
- [X] 2200451894_63004210
- [X] 2200456165_63051294
- [X] 2200459865_63004933
- [X] 2200459968_63004937
- [X] 2200460073_63004943
- [X] 2200460121_63004684
- [X] 2200467051_63003856
- [X] 2200467225_63004940
- [X] 2200467261_63004930
- [X] 2200467338_63004925
- [X] 2200470099_63004485
- [X] 2200470142_63004480
- [X] 2200471780_63004362
- [X] 2200480910_63006466
- [X] 2200495073_63010427
- [X] 2200495510_63009152
- [X] 2200508677_63060252
- [X] 2200510531_63012582
- [X] 2200510628_63012549
- [X] 2200510657_63012554
- [X] 2200511249_63012533
- [X] 2200511274_63012586
- [X] 2200517952_63060399
- [X] 2200519525_63060439
- [X] 2200524009_63014044
- [X] 2200524609_63014046
- [X] 2200524616_63014048
- [X] 2200533429_63060425
- [X] 2200539735_63060406
- [X] 2200549908_63019339
- [X] 2200549965_63019349
- [X] 2200550414_63019357
- [X] 2200550471_63020031
- [X] 2200550490_63019351
- [X] 2200550505_63019340
- [X] 2200555565_63018614
- [X] 2200559438_63020029
- [X] 2200559682_63020030
- [X] 2200559713_63019623
- [X] 2200559739_63019626
- [X] 2200569969_63019991
- [X] 2200570001_63021580
- [X] 2200570025_63021490
- [X] 2200570035_63021491
- [X] 2200570042_63021493
- [X] 2200570050_63021494
- [X] 2200579897_63024910
- [X] 2200583995_63024866
- [X] 2200584035_63024905
- [X] 2200584069_63024888
- [X] 2200584126_63024810
- [X] 2200589507_63026712
- [X] 2200597365_63027994
- [X] 2200597480_63027988
- [X] 2200597752_63026853
- [X] 2200597778_63027992
- [X] 22005977_63026903
- [X] 2200609031_63026527
- [X] 2200614198_63113928
- [X] 2200620372_63030821
- [X] 2200620442_63030810
- [X] 2200620498_63030816
- [X] 2200620628_63031031
- [X] 2200622310_63030984
- [X] 2200622355_63030956
- [X] 2200625369_63028699
- [X] 2200625410_63028697
- [X] 2200625536_63028694
- [X] 2200630189_63030665
- [X] 2200635149_63033182
- [X] 2200644544_63037731
- [X] 2200644594_63037725
- [X] 2200650089_63038093
- [X] 2200666292_63076568
- [X] 2200669188_63036688
- [X] 2200669320_63040259
- [X] 2200669383_63040254
- [X] 2200669414_63040257
- [X] 2200669446_63040251
- [X] 2200680342_63105271
- [X] 2200694535_63042853
- [X] 2200694789_63042862
- [X] 2200694858_63042702
- [X] 2200694917_63042696
- [X] 2200699290_63043047
- [X] 2200699345_63040238
- [X] 2200699383_63043050
- [X] 2200699412_63040731
- [X] 220071551_63048935
- [X] 2200731515_63048963
- [X] 2200748145_63051198
- [X] 2200748171_63051213
- [X] 2200751046_63051249
- [X] 2200751101_63051234
- [X] 2200766471_63054590
- [X] 2200767731_63054595
- [X] 2200767822_63054464
- [X] 2200775505_63060410
- [X] 2200850441_63019345
- [X] 220597589_63026879
- [X] 2300003253_63060430
- [X] 2300005679_63060370
- [X] 2300009914_63060390
- [X] 2300028784_63060001
- [X] 2300036815_63063357
- [X] 2300055382_63061874
- [X] 2300055421_63061871
- [X] 2300055440_63061880
- [X] 230006894_63064950
- [X] 2300071111_63070356
- [X] 2300083434_63071675
- [X] 2300103609_63076239
- [X] 2300104572_63076232
- [X] 2300109602_63076765
- [X] 2300109665_63076770
- [X] 2300119721_63078732
- [X] 2300137773_63078133
- [X] 2300137834_63078123
- [X] 2300167821_63086183
- [X] 2300172698_63113453
- [X] 2300188216_63090609
- [X] 2300188281_63090632
- [ ] 2300188800_63090616
- [ ] 2300193193645_63090623
- [ ] 2300193668_63090611
- [ ] 2300195426_63090608
- [ ] 2300201017_63089636
- [ ] 2300227479_63098330
- [ ] 2300232688_63130821
- [ ] 2300292749_63109239
- [ ] 230029277_63109247
- [ ] 2300294712_63109236
- [ ] 2300308032_63111581
- [ ] 2300323537_63114209
- [ ] 23003
e "Confirmed in sanger" == "" | length
71
❯ open sangerized.csv | where "Found by bisonex" == "missed" | where "Confirmed in sanger" == "" | length
5
❯ open sangerized.csv | where "Found by bisonex" == "missed" | where "Confirmed in sanger" == "true" | length
0
[[id:cd79a77c-a0b6-4bb1-9e08-fe08dc89e3aa][Résultats finaux]]
*** DONE Regarder 5 variants manqués: 3 explicables, 2 non
CLOSED: [2023-11-09 Thu 00:22] SCHEDULED: <2023-11-05 Sun>
open searched.csv | where "Found by bisonex" == "missed"
62982193 7884996 : haplotypecaller ok... -> filtré car AD=5 <= 10
63012582 102230760 : non présent haplotypcellar mais une délétion en 755 (en 754 CG -> C). Vérifié mobidetails
63019340 50721335 : non présent haplotypecaller (vérifié igv). vérifié mobidetails
63060439 26869324 : filtré car 15 reads
63109239 14358800 : présent haplotypecaller : filtré car DP=29 <= 30
Non présent haplotypecaller avec bcftools mais zgrep ok
zgrep 7884996 call_variant/haplotypecaller/*62982193*/*
zgrep 102230760 call_variant/haplotypecaller/*63012582*/*
zgrep 50721335 call_variant/haplotypecaller/*63019340*/*
zgrep 26869324 call_variant/haplotypecaller/*63060439*/*
zgrep 14358800 call_variant/haplotypecaller/*63109239*/*
*** DONE Flowchart
CLOSED: [2023-11-09 Thu 00:22]
*** DONE Refaire extraction
CLOSED: [2023-11-04 Sat 19:02] SCHEDULED: <2023-11-04 Sat>
*** DONE Refaire annotation avec mobidetails
CLOSED: [2023-11-04 Sat 19:02] SCHEDULED: <2023-11-04 Sat>
*** DONE Refaire annotation avec transcrit non reconnus
CLOSED: [2023-11-04 Sat 20:42] SCHEDULED: <2023-11-04 Sat>
5 transcrits, donnés égalemen tpar
#+begin_src nu
open annotated.csv | where coding != "negatif" | where chrom == ""
#+end_src
| 62676048 | NM_001080420.1 | SHANK3 | référénce non valide |
| 62690893 | NM_001080420.1 | KDM6B | idem |
| 62690893 | NM_001080420.1 | KDM6B | même variant |
| 62795429 | NM_016381.3 | TREX1 | NM_033629.5 |
| 63019340 | NM_001080420.1 | SHANK3 | NM_001372044.2 |
SCHEDULED: <2023-11-01 Wed>
*** DONE Rajouter variant pour 63009152
CLOSED: [2023-11-04 Sat 20:47] SCHEDULED: <2023-11-01 Wed>
*** DONE Regénérer annotation avec NC_
CLOSED: [2023-11-04 Sat 18:59] SCHEDULED: <2023-10-31 Tue>
*** DONE Comparer variants manqué avec sanger: 0 confirmés
CLOSED: [2023-11-06 Mon 23:48] SCHEDULED: <2023-11-04 Sat>
*** DONE Annoter variants avec sanger
CLOSED: [2023-11-08 Wed 23:17] SCHEDULED: <2023-11-07 Tue>
*** DONE Mail paul avec résultats
CLOSED: [2023-11-09 Thu 00:22] SCHEDULED: <2023-11-05 Sun>
*** DONE Vérifier coordonnées des 2 variants manquants
CLOSED: [2023-11-12 Sun 16:53] SCHEDULED: <2023-11-11 Sat>
Les 2 sont des homopolymer
- 1er = même variant mais représenté différement
- SHANK3 ?
**** PITX3: filtrée car AD=8
NB: représentation synonyme
Même séquence
>hg38_dna range=chr10:102230742-102230777 5'pad=2 3'pad=2 strand=+ repeatMasking=none
GGAGCCAGCCCGGGGGGGCCCCCGCCCAGGCCCTG
>hg19_dna range=chr10:103990500-103990534 5'pad=0 3'pad=0 strand=+ repeatMasking=none
GGAGCCAGCCCGGGGGGGCCCCCGCCCAGGCCCTG
Selon IGV:
GGAGCCAGCCC(G)GGGGGGCCCCCGCCCAGGCCCTG
Selon cento
GGAGCCAGCCCGGGGGG(G)CCCCCGCCCAGGCCCTG
#+begin_src sh :dir ~/annex/data/bisonex/
bcftools filter -i 'POS=102230760' call_variant/haplotypecaller/*63012582*/*.vcf.gz
#+end_src
DP ok mais AD trop faible
GT:AD:DP:GQ:PL 0/1:26,8:34:99:146,0,671
**** SHANK3: transcrit supprimé depuis: ok
Retrouvé par ERic: 50721504dup
On vérifie
#+begin_src sh :dir ~/annex/data/bisonex/
bcftools filter -i 'POS=50721504' call_variant/haplotypecaller/*63019340*/*.vcf.gz
#+end_src
#+begin_src sh :dir ~/annex/data/bisonex/
zgrep '50721504' annotate/full/*63019340*.tsv
#+end_src
*** TODO Sanger pour 4 VOUS manqués
SCHEDULED: <2023-12-13 Wed>
/Entered on/ [2023-11-13 Mon 22:40]
** TODO Chercher nouveaux gènes
SCHEDULED: <2023-12-06 Wed>
* Résultats
** TODO Speed-up BWA-mem
SCHEDULED: <2023-12-02 Sat>
** TODO Speed-up Hapotypecaller
SCHEDULED: <2023-12-02 Sat>
** TODO Refaire statistics avec happy+ vcfeval
SCHEDULED: <2023-11-30 Thu>
/Entered on/ [2023-11-18 Sat 20:13]
* Communication
** DONE Mail NGS-diag
CLOSED: [2023-10-06 Fri 08:04] SCHEDULED: <2023-10-06 Fri>
/Entered on/ [2023-10-04 Wed 19:33]
e "Confirmed in sanger" == "" | length
71
❯ open sangerized.csv | where "Found by bisonex" == "missed" | where "Confirmed in sanger" == "" | length
5
❯ open sangerized.csv | where "Found by bisonex" == "missed" | where "Confirmed in sanger" == "true" | length
0
[[id:cd79a77c-a0b6-4bb1-9e08-fe08dc89e3aa][Résultats finaux]]
*** DONE Regarder 5 variants manqués: 3 explicables, 2 non
CLOSED: [2023-11-09 Thu 00:22] SCHEDULED: <2023-11-05 Sun>
open searched.csv | where "Found by bisonex" == "missed"
62982193 7884996 : haplotypecaller ok... -> filtré car AD=5 <= 10
63012582 102230760 : non présent haplotypcellar mais une délétion en 755 (en 754 CG -> C). Vérifié mobidetails
63019340 50721335 : non présent haplotypecaller (vérifié igv). vérifié mobidetails
63060439 26869324 : filtré car 15 reads
63109239 14358800 : présent haplotypecaller : filtré car DP=29 <= 30
Non présent haplotypecaller avec bcftools mais zgrep ok
zgrep 7884996 call_variant/haplotypecaller/*62982193*/*
zgrep 102230760 call_variant/haplotypecaller/*63012582*/*
zgrep 50721335 call_variant/haplotypecaller/*63019340*/*
zgrep 26869324 call_variant/haplotypecaller/*63060439*/*
zgrep 14358800 call_variant/haplotypecaller/*63109239*/*
*** DONE Flowchart
CLOSED: [2023-11-09 Thu 00:22]
*** DONE Refaire extraction
CLOSED: [2023-11-04 Sat 19:02] SCHEDULED: <2023-11-04 Sat>
*** DONE Refaire annotation avec mobidetails
CLOSED: [2023-11-04 Sat 19:02] SCHEDULED: <2023-11-04 Sat>
*** DONE Refaire annotation avec transcrit non reconnus
CLOSED: [2023-11-04 Sat 20:42] SCHEDULED: <2023-11-04 Sat>
5 transcrits, donnés égalemen tpar
#+begin_src nu
open annotated.csv | where coding != "negatif" | where chrom == ""
#+end_src
| 62676048 | NM_001080420.1 | SHANK3 | référénce non valide |
| 62690893 | NM_001080420.1 | KDM6B | idem |
| 62690893 | NM_001080420.1 | KDM6B | même variant |
| 62795429 | NM_016381.3 | TREX1 | NM_033629.5 |
| 63019340 | NM_001080420.1 | SHANK3 | NM_001372044.2 |
SCHEDULED: <2023-11-01 Wed>
*** DONE Rajouter variant pour 63009152
CLOSED: [2023-11-04 Sat 20:47] SCHEDULED: <2023-11-01 Wed>
*** DONE Regénérer annotation avec NC_
CLOSED: [2023-11-04 Sat 18:59] SCHEDULED: <2023-10-31 Tue>
*** DONE Comparer variants manqué avec sanger: 0 confirmés
CLOSED: [2023-11-06 Mon 23:48] SCHEDULED: <2023-11-04 Sat>
*** DONE Annoter variants avec sanger
CLOSED: [2023-11-08 Wed 23:17] SCHEDULED: <2023-11-07 Tue>
*** DONE Mail paul avec résultats
CLOSED: [2023-11-09 Thu 00:22] SCHEDULED: <2023-11-05 Sun>
*** DONE Vérifier coordonnées des 2 variants manquants
CLOSED: [2023-11-12 Sun 16:53] SCHEDULED: <2023-11-11 Sat>
Les 2 sont des homopolymer
- 1er = même variant mais représenté différement
- SHANK3 ?
**** PITX3: filtrée car AD=8
NB: représentation synonyme
Même séquence
>hg38_dna range=chr10:102230742-102230777 5'pad=2 3'pad=2 strand=+ repeatMasking=none
GGAGCCAGCCCGGGGGGGCCCCCGCCCAGGCCCTG
>hg19_dna range=chr10:103990500-103990534 5'pad=0 3'pad=0 strand=+ repeatMasking=none
GGAGCCAGCCCGGGGGGGCCCCCGCCCAGGCCCTG
Selon IGV:
GGAGCCAGCCC(G)GGGGGGCCCCCGCCCAGGCCCTG
Selon cento
GGAGCCAGCCCGGGGGG(G)CCCCCGCCCAGGCCCTG
#+begin_src sh :dir ~/annex/data/bisonex/
bcftools filter -i 'POS=102230760' call_variant/haplotypecaller/*63012582*/*.vcf.gz
#+end_src
DP ok mais AD trop faible
GT:AD:DP:GQ:PL 0/1:26,8:34:99:146,0,671
**** SHANK3: transcrit supprimé depuis: ok
Retrouvé par ERic: 50721504dup
On vérifie
#+begin_src sh :dir ~/annex/data/bisonex/
bcftools filter -i 'POS=50721504' call_variant/haplotypecaller/*63019340*/*.vcf.gz
#+end_src
#+begin_src sh :dir ~/annex/data/bisonex/
zgrep '50721504' annotate/full/*63019340*.tsv
#+end_src
*** TODO Sanger pour 4 VOUS manqués
SCHEDULED: <2023-12-13 Wed>
/Entered on/ [2023-11-13 Mon 22:40]
** TODO Chercher nouveaux gènes
SCHEDULED: <2023-12-06 Wed>
* Résultats
** TODO Speed-up BWA-mem
SCHEDULED: <2023-12-02 Sat>
** TODO Speed-up Hapotypecaller
SCHEDULED: <2023-12-02 Sat>
** TODO Refaire statistics avec happy+ vcfeval
SCHEDULED: <2023-12-04 Mon>
/Entered on/ [2023-11-18 Sat 20:13]
* Communication
** DONE Mail NGS-diag
CLOSED: [2023-10-06 Fri 08:04] SCHEDULED: <2023-10-06 Fri>
/Entered on/ [2023-10-04 Wed 19:33]