call the program below like:

python mutate.py some.fasta 0.01 > mutated.fasta

And it will change 1% (0.01) of the bases.

from random import random, choice
import sys
from itertools import groupby

def fasta_iter(fasta_name):
    """
    given a fasta file. yield tuples of header, sequence
    """
    fh = open(fasta_name)
    # ditch the boolean (x[0]) and just keep the header or sequence since
    # we know they alternate.
    faiter = (x[1] for x in groupby(fh, lambda line: line[0] == ">"))
    for header in faiter:
        # drop the ">"
        header = header.next()[1:].strip()
        # join all sequence lines to one.
        seq = "".join(s.strip() for s in faiter.next())
        yield header, seq

def main(fasta_name, mutation_freq):

    for header, seq in fasta_iter(fasta_name):
        seq = list(seq)
        for i, s in enumerate(seq):
            val = random()
            if val < mutation_freq:
                # choose a random nucleotide that's different.
                seq[i] = choice([x for x in "ACTG" if x != s.upper()])
        print ">%s\n%s" % (header, "".join(seq))

if __name__ == "__main__":
    main(sys.argv[1], float(sys.argv[2]))

The EMBOSS program msbar will do this. E.g. to introduce 10 point mutations

msbar -count 10 -point 4 seq.fasta

The argument '4' allows only base changes, but msbar can produce other types of mutation. There is also an argument to ensure that the mutation does not create another specific sequence (e.g. the one you started with). You can use infoseq to obtain the sequence length so that you know how many mutations to request.

If this effort is to support research on human genetics, for example, I would add statements that the random variation mimic known frequencies of transitions (T <-> C, or A <-> G) and transversions (A -> C or T, or G -> C or T, C -> G or A, T -> G or A) in the human genome. The most common single nucleotide change in human is C <-> T. Thus, if T is picked as the base to change, then it should change more often to a C than to another base. This would be important if the project wants to align with a biological situation.

The frequency with which a given nucleotide, selected at random by one of the scripts above, changes to any of the other three nucleotides can be taken by scanning the literature for such a table or by parsing dbSNP. I just don't have number on hand.

You can do this with Biopieces www.biopieces.org) like this:

read_fasta -i test.fna | mutate_seq -p 1 | write_fasta -x

There is also indel_seq for insertions and deletions.

http://code.google.com/p/biopieces/wiki/mutate_seq

http://code.google.com/p/biopieces/wiki/indel_seq

Cheers,

Martin