Merge branch 'master' of https://github.com/KatyBrown/CIAlign

Author: KatyBrown

CIAlign/CIAlign.py

@@ -326,7 +326,7 @@ def main():
         print("Error! Your input alignmnent path could not be found.")
         exit()
 
-    arr, nams = utilityFunctions.FastaToArray(args.infile, args.outfile_stem)
+    arr, nams = utilityFunctions.FastaToArray(args.infile, log, args.outfile_stem)
 
     # check if at least names are unique
     if len(nams) > len(set(nams)):

CIAlign/utilityFunctions.py

@@ -43,7 +43,7 @@ def unAlign(arr):
     return (arr)
 
 
-def FastaToArray(infile, outfile_stem=None):
+def FastaToArray(infile, log, outfile_stem=None):
     '''
     Convert an alignment into a numpy array.
 
@@ -63,19 +63,26 @@ def FastaToArray(infile, outfile_stem=None):
         List of sequence names in the same order as in the input file
     '''
 
+    formatErrorMessage = "The MSA file needs to be in FASTA format."
     nams = []
     seqs = []
     nam = ""
     seq = ""
     with open(infile) as input:
         for line in input:
             line = line.strip()
+            if len(line) == 0:
+                continue # todo: test!
             if line[0] == ">":
                 seqs.append([s.upper() for s in seq])
                 nams.append(nam)
                 seq = []
                 nam = line.replace(">", "")
             else:
+                if len(nams) == 0:
+                    log.error(formatErrorMessage)
+                    print(formatErrorMessage)
+                    exit()
                 seq += list(line)
     seqs.append(np.array([s.upper() for s in seq]))
     nams.append(nam)
@@ -198,6 +205,7 @@ def seqType(arr):
     '''
     Detects if an alignment is of nucleotides or amino acids using pre-built
     dictionarys of amino acid and nucleotide codes.
+    Checks if arr contains characters that are not in the dictionary (not IUPAC)
 
     Parameters
     ----------
@@ -209,27 +217,36 @@ def seqType(arr):
     str
     'aa' for amino acid and 'nt for nucleotide
     '''
-    seq1 = arr[0]
-    nucs = set(list(getNtColours().keys()))
-    aas = set(list(getAAColours().keys()))
-    n = 0
-    a = 0
-    x = 0
-    for s in seq1:
-        s = s.upper()
-        if s in nucs:
-            n += 1
-        if s in aas:
-            a += 1
-        if s not in aas and s not in nucs:
-            x += 1
-    counts = n, a, x
-    if n == max(counts):
+    nt_count = 0
+    aa_count = 0
+    for seq in arr:
+        nucs = set(list(getNtColours().keys()))
+        aas = set(list(getAAColours().keys()))
+        n = 0
+        a = 0
+        x = 0
+        for s in seq:
+            s = s.upper()
+            if s in nucs:
+                n += 1
+            if s in aas:
+                a += 1
+            if s not in aas and s not in nucs:
+                x += 1
+        counts = n, a, x
+        if n == len(seq):
+            nt_count += 1
+        elif a == len(seq):
+            aa_count += 1
+        else:
+            print("Unknown nucleotides or amino acids detected. Please fix your MSA.")
+            exit()
+    if nt_count == len(arr):
         return "nt"
-    elif a == max(counts):
+    elif aa_count == len(arr):
         return "aa"
     else:
-        print("Majority of positions are not known nucleotides or amino acids")
+        print("MSA type couldn't be established. Please fix your MSA.")
         exit()