Generate random DNA sequences with variable %identity based on a reference sequence
Hello there!
I have a query sequence (~500 bases long) and I want to generate multiple sequences which have a specific % identity to this query sequence. For example, I would like to have 3 sequences which are 80%, 60%, and 20% similar to this query sequence. Are there any tools that could generate them?
Thank you!
• 1,174 views
•
link
1 answer
This short python script will read a fasta file and randomly mutate the requested fraction of nucleotides. Save it as mutate_dna.py and run:
python mutate_dna.py dna.fa 0.2 > mutated_dna.fa
You will need to run it multiple times with different mutation rates to get what you want.
Starting with this sequence:
>dna
AGGGATCAGCGGATGTTGCTTATAGGACTCCGCTGGCACCTTATGAGAAA
TCAAAGTTTTTGGGTTCCGGgGGGaGtATGGTCGCaAGGcTGAAACTTAA
AGGAATTGACGGAAGGGCACCACCAGGAGTGGAGCCTGCGGCTTAATTTG
ACGCAACACGGGGAAACTTACCAGGTCCAGACATAGTAAGGATTGACAGA
CTGAGAGCTCTTTCTTGATTCTATGGGTGGTGGTGCATGGCCGTTCTTAG
TTGGTGGAGCGATTTGTCTGGTTAATTCCGTTAACGAACGAGACCTCAGC
CTGCTAACTAGCTACGTGGAGGCATCCCTTCACGGCCGGCTTCTTAGAGG
GACTATGGCCGTTTAGGCCAAGGAAGTTTGAGGCAATAACAGGTCTGTGA
TGCCCTTAGATGTTCTGGGCCGCACGCGCGCTACACTGATGTATTCAACG
AGTTCACACCTTGGCCGACAGGCCCGGGTAATCTTTGAAATTTCATCGTG
ATGGGGATAGATCATTGCAATTGTTGGTCTTCAACGAGGAATTCCTAGTA
AGCGCGAGTCATCAGCTCGCGTTGACTACGTCCCTGCCCTTTGTACACAC
CGCCCGTCGCTCCTACCGATTGAATGATCCGGTGAAGTGTTCGGATCGCG
GCGACGTGGGTGGTTCGCCGCCCGCGACGTCGCGAGAAGTCCACTAAACC
TTATCATTTAGAGGAAGGAGAAGTCGTAACAAGGTTTCCGTAGGTGAACC
TGCGG
It makes this sequence:
>mutated_0.2_dna
TGGGATCAGCGGATGTTGCTTCTAGGTCCCCACTGGCACCTTATGAGAATGGAAAGTTTTTGGATTCCGG
gGGGaGtATGGTGGCaTGGcTGATACTTAAATGAATTGTCGGGAGGGCTCCACCCGGAGTGGATCTTGCG
TCTTAATTCGACTTAAGAGGGGGAACCTTACCCGGTCCAGACAGAGTAAGTATTGACAAACTGAGGGCTC
TTTCTAGAAACCATGGGTGGTGCAGCATGCTCGTTCTTTGTTCGTGGAGCGATTTGTCTGGTTAATTTCC
TTGCCGAACGAGACCTTAGCCAGCTTACTAGCGACGTAGAGCCGTCCCTGGACGTCGGGCATCTTAGAAG
GACTAAGGCCGTTTGGACCAAGGAAGTTGAAGGCATTAACAGGTCTGTGATCACCATCGATGTTCGGGGC
CGCACCCGCGCTACTCTGATGTATTTAACGAGTTCCCGTCCTGGCACTCTGGCCCAGGTAATCCTAGATT
TTTCATCATGCTGGGGTTAAAGCATTCCAGTTGTAGGTCTTATGCGAGGGATTTCTAGTAAGCACGAGTC
TTCAGTTCACGTAGACTATGTCCTTGCCCGTTCTACACACCACCCGGCGCTTCTACCGTATGAATGTTCG
ACTGAAGTGTTCGCATCTCGGCGACGTGGGCGGTTCGGCGCCCCCGACGTCGCTAGAAGTCTACGAAACC
TTATCATCTAAAAGAAGGACAAGTCGTATAAAGGCTTCCGTAGGTGAGCCTGCGG
import os
import sys
from random import randint
import textwrap
import numpy as np
from Bio import SeqIO
from Bio.Seq import MutableSeq
def mutate(letter):
pos = randint(0, 2)
if (letter == 'A') or (letter == 'a'):
bases = ['C', 'G', 'T']
if (letter == 'C') or (letter == 'c'):
bases = ['A', 'G', 'T']
if (letter == 'G') or (letter == 'g'):
bases = ['A', 'C', 'T']
if (letter == 'T') or (letter == 't'):
bases = ['A', 'C', 'G']
return str(bases[pos])
if os.access(sys.argv[1], os.R_OK):
FastaFile = open(sys.argv[1], 'rU')
else:
print('\n !!! Input file "%s" does not exist in this directory !!!\n' %
sys.argv[1])
sys.exit(1)
mut_rate = float(sys.argv[2])
if (mut_rate > 1) or (mut_rate < 0):
print('\n Sequence mutation rate must be in 0-1 range !\n')
sys.exit(1)
for rec in SeqIO.parse(FastaFile, 'fasta'):
name = rec.id
seq = MutableSeq(list(rec.seq))
seq_len = len(rec)
mut_len = int(seq_len * mut_rate)
mut_pos = np.random.randint(0, seq_len, mut_len)
for x in range(mut_len):
seq[mut_pos[x]] = mutate(seq[mut_pos[x]])
rec.seq = seq
print('>mutated_%s_%s' % (str(mut_rate), name))
for lines in textwrap.wrap(str(rec.seq)):
print(lines)
FastaFile.close()
• 0 views
•
link
Log in to answer this question.