CSC 15 Lab 8

	      At Least Part I is due next Tuesday 


The target for this lab is to write some procedures for string
matching.  We will apply these procedures to DNA sequence analysis.
Some of the procedures will also involve file io.

                            PART I:

0. Write a procedure that generates a random DNA sequence of length n:

def randomDNA(n):

   Think about the algorithm first: you'll start with an empty string "",
   then run a loop n times.  Inside the loop, you need to generate a 
   randint(0,3) and depending on the value of this random number, you
   will add either an "A", "C", "G" or "T" to the end of your string.
   For example, to add the character "C" to the end of string S, use the
   statement S = S+"C".

   The function should return the string generated.



#######
1. Write a function count(X,S) that returns the number of times that the
character X appears in string S.  For example, count("A","GATCAATC") should
return 3.  Test it on your random DNA sequences.


#######
The following procedures may require procedures to transform strings to
arrays of chars and vice versa:

def stringtoarray(S):
        A = []
        i = 0
        while i<len(S):
            A.append(S[i])
            i = i+1
        return A

def arraytostring(A):
            S = ""
            for c in A:
                S = S+c
            return S

######
2. A DNA sequence mutates when some of its code changes randomly.
   Write a procedure that randomly mutates up to n of the codes in the
   DNA sequence.  For example, mutate("ATTGCGA",2) may return
  "ATAGCGC".  

        Algorithm: first convert the DNA string into an array of chars
        so that it can be easily modified.  Then run a loop n times.
        inside the loop, generate a random array index, and modify the
        array content at the index.  Calling your randomDNA function
        with an argument of 1 should return a random letter: "A", "C",
        "G" or "T".  After the loop is finished, convert the array
        back to a string and return the string.

        Note that it's possible that a mutation could end up with the
        same letter as before.  That's why the function is specified
        to mutate "up to" n of the codes.

    CHALLENGE:  Write another version of this function that mutates EXACTLY
    n codes.


def mutate(dna,n):


######
3. Write a function "splice":

  def splice(A,B,j)

  Where A is a DNA string, B is a DNA string and j is a number.

  The function needs to replace a segment of B with the contents of A.
  j will indicate the starting index in B where this replacement will
  take place.  For example, splice("AAA","ATGTTGC",2) should return the
  string "ATAAAGC": the contents of string B from index 2 (third code) 
  is replaced by string A.  

  Your function also needs to detect error conditions: it's possible that
  the length of A is larger than the length of B.  It's also possible that
  j+len(A) is larger than the length of B.  In either case, your function
  should return B unchanged.
   
  This time, you'll have to come up with the algorithm yourself.


######
4. Using your functions above, generate several dna sequences large and
   small.  mutate and splice the smaller sequences into the larger ones.
   Write the dna sequences to file.  Please do not create DNA sequences 
   with lengths longer than 2**20 (1 megabyte).


You can use the following function to write a string to a file:

# write a string S to a file:
def writefile(S, filename):
    fd = open(filename,"w")
    fd.write(S)
    fd.close()

You need to write a complimentary function that reads the file and return
its contents as a string.  Look at your notes, sample code and tutorial
section on reading files.

def readfile(filename)  # should return a string containing the file contents
    



############################################################################

#   Part II:

After years of painstaking research geneticists have identified the
DNA sequence responsible for causing some computer science students to
doze off during class.

The Dozer Gene:

GTGTGCAGCGTCCCCGATGCTACTCCATTCCTTTCGCGTGTATCTTGAGAGCAATGATGGGATTTACAATACAGATACGTAAACGCAACGTATCAACAGGTAGGGGTTCTGTTTGGAGCCAACTAGTC


Four DNA samples were randomly taken from students at an unnamed university in
Nassau county.  Their files are found on the homepage as sd1.dna, sd2.dna,
sd3.dna, and sd4.dna.

You are to determine the probability that a student carries the dozer gene.
This amounts to finding if one string A is a substring of another string B,
WHILE FORGIVING A NUMBER OF ERRORS.  For example, "abcde" matches "abudg"
with 2 errors.  You need to write a function that takes two strings and
a number indicating the maximum number of errors allowed.  The function
needs to determine the best possible match.  It should compute two values:
the starting position in the longer sequence where the best match occurs,
and the number of errors in the best match.  If there is no match that
contains less than or equal to the maximum number of allowed errors, then
the starting position value should be -1.  For example, calling

 subsequence("defgh","abcefkhdfsfdefgms",3)   // 3 indicates max errors

should return the tuple (11,1): that is, the starting position of the
best match is 11 and the number of errors is 1, indicating the
substring "defgm", which occurs starting at index 11 (idices start
from 0) in the longer string.  Note that the substring "cefkh", which
starts at position 2, is also a match, but has 2 errors.  Your
function needs to report the best possible match.

If the entire dozer gene is found embedded in a dna file, with no
errors, then that person is considered to have a 100% chance of being
a dozer.  For every error, the probability decreases by 5 percent.  So
if there are two errors, the probability that the student is a dozer
becomes 90%. If there are 20 or more errors, then the probability is
0%.

Your program should produce output in the following FORMAT (exact values
will differ)

DNA SAMPLE 1:
A best match with 4 errors was found starting at position 2018.
There is a 80% chance that this student is a dozer


########################################

While conducting their research on the dozer gene, our brilliant
geneticists also took note of another phenomenon: many of the students
who were dozing off during class were also drooling.  The gene that
causes this behavior was also identified:

The Drooler Gene:

ACTAAGGTTGCGTCGATTACAAAAATGGCGCTTGTAGCTTTACAAGCTACGTTTTTGCCGTGCTTGTAGGATAGAGTCAGGCGGAGCGGGCCTAAGAGTGTGAATACATGGGCGAGTCGGGGGATTATCTTTAGGTAAGGCATGATACCAGTTTACTCCTTCTTCTGTCCGATGTGTGACCTAGTCGGTGGAAATGATTGTTTATCGAGACCCATTGTAGACTGAAGCGAAACTTTGGTGTTCTGAGGTCTCAAGT


One geneticist boldly proposed the following, controversial hypothesis:

The "D & D" Hypothesis:  

   "Every carrier of the drooler gene also carries the dozer gene"


You need to modify your experiments to either lend proof to the hypothesis,
or disprove it by finding a counter-example (i.e., a student that drools
while awake). 

Report your findings.

----------------------------------------

EXTRA: DNA sequences may also contain gaps representing unknown code.
Modify your procedures to account for sequences such as "TA-CAG-TGGA".
When finding a match between two sequences, you have a choice on how
to interpret the gap: either as a "wild card", which means that it
matches everything, or report an error (consider it a mismatch).  Your
functions should be parameterized to be able to do BOTH.  That is, add
another parameter, perferably a boolean parameter, that determines how
the gap is to be interpreted.