apraga/org - Change WHKJYYZX7Y2VYHEX4NAKKV4DHI3RBHY3XJMIGPYBPCD3DNWZUHDQC

Bisonex: isec does not manage properly multiple alternate allelel...

dbSNP identifiants are not always filled in in clinvar

Created by Alexis Praga on November 20, 2022

WHKJYYZX7Y2VYHEX4NAKKV4DHI3RBHY3XJMIGPYBPCD3DNWZUHDQC

Dependencies

In channels

main

Change contents

Insertion in projects.org at line 26 [5.123895]

[3.8274]

[4.2241]

** Quillbot
Possible d'utiliser selenium tout le long
Mais si on copy lese cookies et headers depuis le navigateur:
- on peut bypasser cloudflare (vu qu'on a le résultat)
- mais au lieu d'avoir un json,  on a du binaire...

Replacement in projects.org at line 1611 [5.123895]

B:BD[6.15] → [7.2320:2408]

∅:D[7.2408] → [6.67:116]

B:BD[6.67] → [6.67:116]

******* DONE Test: variants dbSNP étiquettés "clinvar"
CLOSED: [2022-10-23 Sun 18:55]
$ zgrep CLNSIG dbSNP_common.vcf.gz | wc -l
57997

[7.2319]

[6.116]

******* KILL Utiliser directement le numéro dbSNP dans clinvar ? Non
CLOSED: [2022-11-20 Sun 19:51]
Ex: chr20
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools query -f 'rs%INFO/RS \n' -i 'INFO/RS != "." & INFO/CLNSIG="Pathogenic"' clinvar_chr20.vcf.gz | sort > ID_clinvar_patho.txt
bcftools query -f '%ID\n' dbSNP_common_chr20.vcf.gz | sort > ID_of_common_snp.txt
comm -23 ID_of_common_snp.txt ID_clinvar_patho.txt > ID_of_common_snp_not_clinvar_patho.txt
wc -l ID_of_common_snp_not_clinvar_patho.txt
# sort ID
#+end_src

Replacement in projects.org at line 1622 [5.123895]

B:BD[6.117] → [6.117:161]

$ zgrep "RS=" clinvar.vcf.gz | wc -l
834352

[6.117]

[8.267]

#+RESULTS:
: 518846 ID_of_common_snp_not_clinvar_patho.txt

Replacement in projects.org at line 1625 [5.123895]

B:BD[8.268] → [9.237:288]

∅:D[9.288] → [6.216:286]

B:BD[6.216] → [6.216:286]

Donc il faut vraiment faire l'intersection des VCF
soit avec vcftools (naïve) soit avec rtgtools (gères les synonymes)

[8.268]

[10.99]

Version d'alexis
#+begin_src sh :dir ~/code/bisonex/test_isec
snp=dbSNP_common_chr20.vcf.gz
clinvar=clinvar_chr20_notremapped.vcf.gz
python ../script/pythonScript/clinvar_sbSNP.py \
    --clinvar $clinvar \
    --chrm_name_table ../database/RefSeq/refseq_to_number_only_consensual.txt \
    --dbSNP $snp --output prod.txt
wc -l prod.txt
zgrep '^NC' dbSNP_common_chr20.vcf.gz | wc -l
#+end_src
#+RESULTS:
| 518832 | prod.txt |
| 518846 |          |

Replacement in projects.org at line 1653 [5.123895]
B:BD[10.663] → [11.1537:1570]
```
******* DONE Essai bcftools isec
```
[10.663]
[7.2470]
```
******* DONE bcftools isec : non
```

Replacement in projects.org at line 1724 [5.123895]

B:BD[7.4817] → [7.4817:4855]

∅:D[7.4855] → [10.1087:1216]

B:BD[10.1087] → [10.1087:1216]

B:BD[10.1216] → [11.1571:1650]

B:BD[11.1650] → [12.160:237]

******* 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
******* DONE Valider avec Alexis : bcftool isec
CLOSED: [2022-11-07 Mon 21:42]
******* TODO Pourquoi nombre de lignes différentes avec la version d'Alexis

[7.4817]

[12.237]

******** DONE Valider avec Alexis : bcftool isec
CLOSED: [2022-11-07 Mon 21:42   ]
******** DONE Pourquoi nombre de lignes différentes avec la version d'Alexis -> isec ne gère pas plusieurs ALT
CLOSED: [2022-11-20 Sun 17:56]

Replacement in projects.org at line 1736 [5.123895]
∅:D[2.47] → [12.551:604]
B:BD[12.551] → [12.551:604]
```
******** Clinvar classe 4 ? Moins mais toujours trop
```
[2.47]
[12.604]
```
********* Clinvar classe 4 ? Moins mais toujours trop
```

Replacement in projects.org at line 1741 [5.123895]

B:BD[12.669] → [12.669:719]

******** Comparer les ID et regarder ceux en plus

[12.669]

[12.719]

********* Comparer les ID et regarder ceux en plus

Insertion in projects.org at line 1758 [5.123895]

[2.152]

[13.14]

********* Tester sur chromosome 20
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools view --regions NC_000020.11 ../database/dbSNP/dbSNP_common.vcf.gz -o dbSNP_common_chr20.vcf.gz
bcftools view --regions 20 ../database/clinvar/clinvar.vcf.gz -o clinvar_chr20.vcf.gz
tabix dbSNP_common_chr20.vcf.gz
tabix clinvar_chr20.vcf.gz
#+end_src
#+RESULTS:
Attention à bien renommer clinvar !
#+begin_src sh :dir ~/code/bisonex/test_isec
mv clinvar_chr20.vcf.gz clinvar_chr20_notremapped.vcf.gz
bcftools annotate --rename-chrs chromosome_mapping.txt clinvar_chr20_notremapped.vcf.gz -o clinvar_chr20.vcf.gz
#+end_src
#+RESULTS:
On teste l'intersection dbsnp et clinvar patho ainsi que le complémentaire
*ATTENTION*: sans indexer les vcf, les fichiers seront *VIDES*
*ATTENTION*: il faut include tous les variants de SNP avec loption -i-  (par défaut les filtres s'appliquent sur les 2)
#+begin_src sh :dir ~/code/bisonex/test_isec
snp=dbSNP_common_chr20.vcf.gz
clinvar=clinvar_chr20.vcf.gz
bcftools isec -i 'INFO/CLNSIG="Pathogenic"' -p test $snp $clinvar
for i in test/*; do echo $i; grep '^NC' $i | wc -l; done
#+end_src
#+RESULTS:
|   test/0000.vcf |
|               0 |
|   test/0001.vcf |
|            1787 |
|   test/0002.vcf |
|               0 |
|   test/0003.vcf |
|               0 |
| test/README.txt |
|               0 |
|  test/sites.txt |
|            1787 |
A priori, pas de clinvar patho dans dbSNP sur ce chromosome...
Clinvar stocke le numéro dbSNP donc on pourrait vérifier mais il faudrait aller voir si le SNP est bien "common"
Si on utilise la version d'Alexis
#+begin_src sh :dir ~/code/bisonex/test_isec
snp=dbSNP_common_chr20.vcf.gz
clinvar=clinvar_chr20_notremapped.vcf.gz
python ../script/pythonScript/clinvar_sbSNP.py \
    --clinvar $clinvar \
    --chrm_name_table ../database/RefSeq/refseq_to_number_only_consensual.txt \
    --dbSNP $snp --output prod.txt
wc -l prod.txt
zgrep '^NC' dbSNP_common_chr20.vcf.gz | wc -l
#+end_src
#+RESULTS:
| 518832 | prod.txt |
| 518846 |          |
On devrait donc avoir 14 common snp clinvar patho...
Test :
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools query -f '%ID\n' dbSNP_common_chr20.vcf.gz | sort > chr20_id.txt
sort prod.txt > prod_sorted.txt
comm -23 chr20_id.txt prod_sorted.txt
#+end_src
#+RESULTS:
| rs1044396   |
| rs1131695   |
| rs146917730 |
| rs1799990   |
| rs1801138   |
| rs1801545   |
| rs181943893 |
| rs2424926   |
| rs35873579  |
| rs35938843  |
| rs36106901  |
| rs3827075   |
| rs373200654 |
| rs79338570  |
Si on prend le premier
#+begin_src sh :dir ~/code/bisonex/test_isec
bcftools query -f '%CHROM %POS %ID %REF %ALT\n' -i 'ID="rs1044396"' dbSNP_common_chr20.vcf.gz
bcftools query -f '%CHROM %POS %ID %REF %ALT\n' -i 'INFO/RS="1044396"' clinvar_chr20.vcf.gz
#+end_src
#+RESULTS:
| NC_000020.11 | 63349782 | rs1044396 | G | A,C |
| NC_000020.11 | 63349782 |     93427 | G | A   |
| NC_000020.11 | 63349782 |    857384 | G | C   |
isec ne gère donc pas la représentation des allèles alternatifs...
#+begin_src sh :dir ~/code/bisonex/test_isec
snp=dbSNP_common_chr20.vcf.gz
clinvar=clinvar_chr20.vcf.gz
bcftools isec -i- -i 'INFO/CLNSIG="Pathogenic"' -c none -p test $snp $clinvar
for i in test/*; do echo $i; grep '^NC' $i | wc -l; done
#+end_src
#+RESULTS:
|   test/0000.vcf |
|          518846 |
|   test/0001.vcf |
|            1787 |
|   test/0002.vcf |
|               0 |
|   test/0003.vcf |
|               0 |
| test/README.txt |
|               0 |
|  test/sites.txt |
|          520633 |
******* 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

Replacement in projects.org at line 1892 [5.123895]
B:BD[13.185] → [12.1538:1559]
```
***** KILL rtg-tools
```
[13.185]
[12.1559]
```
***** KILL 'rtg-tools
```
Replacement in projects.org at line 1930 [5.123895]
B:BD[14.778] → [12.1788:1831]
```
**** TODO Tester version d'alexis avec Nix
```
[14.778]
[12.1831]
```
**** STRT Tester version d'alexis avec Nix
```

Deletion in projects.org at line 1933 [5.123895]

B:BD[15.94] → [12.1859:1963]

***** TODO Corriger numérotation chromosome pour clinav
***** STRT Ajouter common snp not clinvar path

Replacement in projects.org at line 1956 [5.123895]

B:BD[12.2391] → [12.2391:2577]

***** TODO ID commo snp not clinvar patho
Résoudre [[*Pourquoi nombre de lignes différentes avec la version d'Alexis][Pourquoi nombre de lignes différentes avec la version d'Alexis]]

[12.2391]

[12.2577]

***** STRT ID commo snp not clinvar patho
- isec : [[*Pourquoi nombre de lignes différentes avec la version d'Alexis -> isec ne gère pas plusieurs ALT][Pourquoi nombre de lignes différentes avec la version d'Alexis -> isec ne gère pas plusieurs ALT]]
- utiliser directement l'identifiant dbSNP dans clinvar ?
- rtgtools ?