from __future__ import print_function

import sys
import os
import argparse
import math
import errno

import subprocess
from constants import PAM_TEMPLATE_SEQUENCE, DNA_ALPHABET

# Nucleotides we want to mutate are located at positions 5,6,7 in chain D (PAM,NGG)
# and 8,7,6 in chain C (target strand, NCC).
positionIndex = { 0 : "5", 1 : "6", 2 : "7", 3 : "8"}
# Create a dictionary mapping corresponding positions to each other.
positionPairs = { "5" : "8", "6" : "7", "7" : "6", "8" : "5"}

# Create dictionary mapping valid base pairs to each other.
basePairs = {'a' : 't', 'c' : 'g', 'g' : 'c', 't' : 'a'}

# Function to mutate base pairs
def mutateNucleotide(fileNameOriginal,PAM):
    mutationString = createMutationString(PAM)
    print mutationString
    subprocess.call(["mutate_bases", mutationString, fileNameOriginal, "4UN3." + PAM + ".pdb"])

# Function to determine command line output for mutations and positions for
# call call of mutate_base
def createMutationString(PAM):
    # Set up some variables
    mutationString = ""
    # Loop through the PAM sequence and find positions to mutate
    for n in range(len(PAM)):
        pos = positionIndex[n]
        base = PAM[n]
        complementaryBase = basePairs[base]
        complementaryPosition = positionPairs[pos]
        # If the first base is T, don't mutate it
        if base== 't' and n == 0:
            continue
        # If the 2nd or 3rd bases are G, don't mutate them
        if all([base == 'g', any([n == 1, n == 2])]):
            continue
        mutationString = mutationString + "c=d s=" + pos + " m=D" + base + \
                "; c=c s=" + complementaryPosition + " m=D" + complementaryBase + "; "
    mutationString = mutationString[0:-2]
    return mutationString
# Going to assume file"4UN3.clean.pdb" is always used and in the same directory.
# This file has a TGG PAM sequence already present in it
fileNameOriginal = "4UN3.original.pdb"

# Iterate through all possible PAM sites
listOfBases = ['a','c','g','t']
for i in listOfBases:
    for j in listOfBases:
        for k in listOfBases:
            for l in listOfBases:
                PAM = i+j+k+l
                # TGG PAM site already exists, let's skip it
                mutateNucleotide(fileNameOriginal,PAM)
                print "mutate PAM to " + PAM
                # Save and close all files

@@ -2001,16 +2001,1361 @@
 g','t'%5D%0A
+%0Aif __name__ == '__main__':%0A    # create parser and parse arguments%0A    parser = argparse.ArgumentParser(description='Generate PDBs with new PAM sites based on input PDB with Cas9 variant')%0A    parser.add_argument('-n', '--num_pams', metavar='N', type=str, # string type because of how this script must be run by Chimera%0A                        help='how many PAMs in total to run. 64 = all PAMs of length 3, 256 = all PAMs of length 4')%0A    parser.add_argument('-i', '--input_pdb', metavar='F', type=str,%0A                        help='input PDB file containing Cas9 mutant of interest')%0A    parser.add_argument('-o', '--output_dir', metavar='D', type=str,%0A                        help='path to output directory for new PDBs')%0A    args = parser.parse_args()%0A%0A    assert args.num_pams is not None%0A    assert args.input_pdb is not None%0A    assert args.output_dir is not None%0A    args.num_pams = int(args.num_pams) # convert string from command line to int%0A    pam_length = int(math.log(args.num_pams, 4))%0A    assert 64 == args.num_pams or 256 == args.num_pams%0A    assert os.path.isfile(args.input_pdb)%0A%0A    try: # check existence again to handle concurrency problems%0A        os.makedirs(args.output_dir)%0A    except OSError as exc:%0A        if exc.errno == errno.EEXIST and os.path.isdir(args.output_dir):%0A            pass%0A        else: raise%0A%0A%0A    
 for i in
@@ -3364,24 +3364,28 @@
 istOfBases:%0A
+    
     for j in
@@ -3402,24 +3402,28 @@
 es:%0A        
+  
+  
 for k in lis
@@ -3440,24 +3440,28 @@
             
+   
+ 
 for l in lis
@@ -3486,16 +3486,20 @@
         
+    
 PAM = i+
@@ -3520,16 +3520,20 @@
         
+ 
+   
 # TGG PA
@@ -3581,24 +3581,28 @@
             
+  
+  
 mutateNucleo
@@ -3640,24 +3640,28 @@
             
+    
 print %22mutat
@@ -3677,16 +3677,20 @@
 %22 + PAM%0A
+    
         

Input arguments like chimera_generate_pams.py