Added capabilty to run CodeML on multiple threads.

Author: icwells

AlignmentProcessor.py

@@ -112,10 +112,14 @@ def splitFasta(fasta, outdir):
     if sf.returncode == 0:
         return True
 
-def rmHeader(outdir):
+def rmHeader(outdir, commonNames):
     '''Removes header information and changes genome build to common name.'''
     print("Changing header...")
-    rh = Popen(split("python bin/02_RemoveHeader.py " + outdir))
+    if commonNames == False:
+        rh = Popen(split("python bin/02_RemoveHeader.py " + outdir))
+    elif commonNames == True:
+        rh = Popen(split("python bin/02_RemoveHeader.py " + outdir +
+                         "--changeNames"))
     rh.wait()
     if rh.returncode == 0:
         return True
@@ -193,14 +197,15 @@ def phylipConvert(outdir, starttime, codeml):
             print("Total runtime: ", elapsedtime)
         return True
 
-def runcodeml(outdir, retainStops, starttime):
+def runcodeml(cpu, outdir, retainStops, starttime):
     '''Runs codeml on a directory.'''
     print("Running codeml...")
     if retainStops == False:
-        cm = Popen(split("python bin/07_CodeMLonDir.py " + outdir))
+        cm = Popen(split("python bin/07_CodeMLonDir.py -i " + outdir
+                         + " -n " + cpu))
     elif retainStops == True:
-        cm = Popen(split("python bin/07_CodeMLonDir.py " + outdir +
-                         "--retainStops"))
+        cm = Popen(split("python bin/07_CodeMLonDir.py -i " + outdir +
+                         " --retainStops -n " + cpu))
     cm.wait()
     if cm.returncode == 0:
         elapsedtime = datetime.now() - starttime
@@ -240,12 +245,15 @@ def helplist():
 def main():
     starttime = datetime.now()
     # Set optional parameters to False:
-    axt = False
+    commonNames = False
     retainStops = False
+    axt = False
     phylip = False
     kaks = False
     codeml = False
     conv = False
+    # Popen requires string input, so the following are typecast as strings
+    cpu = "1"
     percent = "0.5"
     ref = "void"
     # Extract values from command line:
@@ -283,6 +291,10 @@ def main():
             codeml = True
         elif i == "--ucsc":
             conv = True
+        elif i == "-n":
+            cpu = str(argv[argv.index(i) + 1])
+        elif i == "--changeNames":
+            commonNames = True
         elif i == "--retainStops":
             retainStops = True
     # Check inout commands prior to running:
@@ -316,7 +328,7 @@ def main():
         else:
             sf = splitFasta(fasta, outdir)
         if sf == True:
-            rh = rmHeader(outdir)
+            rh = rmHeader(outdir, commonNames)
         if rh == True:
             cf = checkFrame(outdir, ref)
         if cf == True:
@@ -340,7 +352,7 @@ def main():
             # Run codeml
             if codeml == True and kaks == False:
                 if pc == True:
-                    runcodeml(outdir, retainStops, starttime)
+                    runcodeml(cpu, outdir, retainStops, starttime)
 
 if __name__ == "__main__":
     main()

AlignmentProcessorReadMe.txt

@@ -162,14 +162,29 @@ in the bin/ directory.
 
 	-i	the path to your input fasta alignment.
 	-o	the path to your working/output directory.
-	-r	the common name of your reference species (more below).
+	-r	the build or common name of your reference species (more below).
 
 # Optional Arguments:
 
 	--ucsc	This will invoke 00_convertHeader.py, which will convert the 
 		headers from UCSC fasta files so they only contain build
 		names and gene IDs. This does not need to be run on Stich Gene
 		Blocks output.
+
+	--retainStops	This will tell AlignmentProcessor to retain sequences
+			that contain internal stop codons. By default,
+			sequences with internal stop codons will be removed 
+			from the analysis as they may bias the results.
+
+	--changeNames	Tells the program to change genome build names to 
+			commom names (more below).
+
+	-%	a decimal value specifying the minimum percentage of reads 
+		that must remain after replacing unknown codons with gaps 
+		(Default = 0.5). You may wish to use a lower threshold
+		for highly diverged species or for low quality genomes.
+
+	# Ka/Ks_Calculator/CodeML related arguments
 
 	--axt/phylip	Specifies which format to convert the files into 
 			(axt format for KaKs_Calculator; phylip for CodeML 
@@ -186,15 +201,10 @@ in the bin/ directory.
 			This file must be titled "codeml.ctl" (the default 
 			name given by PAML). 
 
-	--retainStops	This will tell AlignmentProcessor to retain sequences
-			that contain internal stop codons. By default,
-			sequences with internal stop codons will be removed 
-			from the analysis as they may bias the results.
-
-	-%	a decimal value specifying the minimum percentage of reads 
-		that must remain after replacing unknown codons with gaps 
-		(Default = 0.5). You may wish to use a lower threshold
-		for highly diverged species or for low quality genomes.
+	-n	if "--codeml" is selected, you may specify the number of CPUs 
+		to run CodeML. CodeML itself cannot be parallelized, but 
+		AlignmentProcessor can call multiple instances of CodeML to
+		shorten overall run time. (Default = 1)
 
 # Additional Commands
 
@@ -431,7 +441,7 @@ memory.
 Change directory into the AlignmentProcessor folder. Paste the followig into
 a terminal:
 
-python AlignmentProcessor.py --axt --kaks --ucsc -r Green_anole \
+python AlignmentProcessor.py --axt --kaks --ucsc -r anoCar2 \
 -i kaksTest.fa -o test/
 
 This will return a text file with 11 lines.
@@ -441,7 +451,7 @@ The test directory already contains sample CodeML control and tree files, so
 all you need  to do is change into the AlignmentProcessor direcotry and paste
 the following:
 
-python AlignmentProcessor.py --phylip --codeml --ucsc -r Green_anole \
+python AlignmentProcessor.py --phylip --codeml --ucsc -n 2 -r anoCar2 \
 -i codemlTest.fa -o test/
 
 There should be 8 .mlc files in the 07_codeml directory.

bin/02_RemoveHeader.py

@@ -10,7 +10,7 @@
 from sys import argv
 from glob import glob
 
-def removeHeader(path):
+def removeHeader(path, commonNames):
     # Open all input files in the directory
     inpath = path + "01_splitFastaFiles/" + "*.fa"   
     files = glob(inpath)   
@@ -25,7 +25,11 @@ def removeHeader(path):
                     if line[0] == ">":
 			# Split header to get build name and remove extra info
                         header = line.split(".")
-                        name = changeNames(header[0])
+                        # Replace build name with common name
+                        if commonNames == True:
+                            name = changeNames(header[0])
+                        elif commonNames == False:
+                            name = header[0][1:]
 			# Write common name to file
                         output.write(">" + str(name) + "\n")
                     else:
@@ -56,9 +60,15 @@ def main():
 <path to inut and output directories>")
         quit()
     else:
+        commonNames = False
         # Set directory names and add a trailing "/" if necessary
         path = argv[1]
-        removeHeader(path)
+        try:
+            if argv[2] == "--changeNames":
+                commonNames = True
+        except IndexError:
+            pass
+        removeHeader(path, commonNames)
 
 if __name__ == "__main__":
     main()

bin/02_nameList.txt

@@ -1,4 +1,4 @@
-Aapl1	Anolis_apletophallus
+Aapl1	Slender_anole
 Aaur1	Grass_anole
 Afre1	Bridled_anole
 ailMel1	Panda