forked from holtzy/Bioinformatics-Tool-Box
-
Notifications
You must be signed in to change notification settings - Fork 0
/
nucleotidic_distance_between_individuals.py
141 lines (78 loc) · 3.42 KB
/
nucleotidic_distance_between_individuals.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
#!/usr/bin/python
# -*- coding: utf-8 -*-
#-------------------------------------
#
# SCRIPT PYTHON Nucleotidic_distance_between individuals
# A Partir d'une matrice de SNP, calcule la distance entre accessions.
# Les fichiers d'entrées sont :
# - matrice de SNP
#
# Yan Holtz, [email protected] , Alban Besnard
#-------------------------------------
import os
import sys
import re
try:
import argparse
except ImportError:
print"oops, the import /argparse/ didn't work"
My_description="\n\nA partir d\'une matrice de SNP, ce programme va calculer la distance entre individus."+"\n---\n"+"La distance est calculée tel que : distance = (nb de génotype identiques / nb de cas identiques) + nb de cas différents"+"\n---\n"+"Format d'entrée : une ligne d'entete qui comporte le nom des individus de chaque colonne. Et une seule premiere colonne qui donne le nom du SNP ou marqueur."+"\n---\n"+"La matrice de sortie est une demi matrice supérieure\n\n"
parser = argparse.ArgumentParser(description= My_description)
parser.add_argument('-SNP', required=True, help=' fichier de SNP à analyser. Format : une ligne d\'entete avec nom_SNP + tab + nom_individu1 + tab + nom2 etc.... Puis une ligne par SNP, avec le nom du SNP puis le genotype de chaqueindividu en colonne, au format A/T par exemple et des tiret pour les NA')
parser.add_argument('-out', required=True, help=' output')
args = parser.parse_args()
fic_SNP=args.SNP
out=args.out
#------------------------------------------------------------------------------------------------------
#------------------------------------------------------------------------------------------------------
### STEP 0 : Détermination du noombre d'individus
fic = open(fic_SNP,'r')
entete = fic.readline();
nb_indiv=int(len(entete.split("\t")))
fic.close()
print"\n-------\n\nLe fichier de SNP comporte "+str(nb_indiv-1)+" accessions\n-------\n\n"
#------------------------------------------------------------------------------------------------------
#------------------------------------------------------------------------------------------------------
### STEP 1 : Calcul des distances
tmp=open(out , "w")
#Fabrication de m'en tete du fichier de sortie. Je transforme l'entete initiale en une nouvelle "to_add"
entete=entete.strip()
entete=entete.split("\t")
to_add="versus"
for i in entete[1:]:
to_add=to_add+"\t"+i
tmp.write(to_add)
#Je commence les calculs ligne par ligne = accession par accession
for i in range(1,nb_indiv):
print "traitement de l'individu "+str(i)+" sur "+str(nb_indiv-1)+" individus"
indiv=entete[i]
tmp.write("\n" + indiv)
#Pour chaque autre accession
for j in range(1 , nb_indiv) :
#Si on a déja traité l'accession, alors on met juste un tiret
if j<i:
tmp.write("\t-")
#Sinon on fait le calcul de distance
if j>=i :
nb_same=0
nb_different=0
num_ligne=0
for line in open(fic_SNP):
line=line.strip()
num_ligne+=1
if num_ligne>1:
line=line.split("\t")
a=line[i]
b=line[j]
if a!="-" and b!="-":
if a==b:
nb_same+=1
else:
nb_different+=1
# petit conditionnement si les individus ne sont pas assez couverts (nécessaire si il y a trop peu de marqueurs sur un individus)
if float(nb_same+nb_different)!=0 :
distance=float(nb_same)/float(nb_same+nb_different)
tmp.write("\t"+str(distance))
else:
tmp.write("\t-")
tmp.write("\n")