在基因工程中的重组DNA构建的实验中，往往需要用限制性核酸内切酶对DNA分子进行切割，已得到目的基因，但DNA分子上可能含有多个酶切位点，因此需要认为的计算酶切后得到的DNA片段（fragements）大小的分子量，以便后续电泳中可以根据Marker得到目的基因的片段，但对于酶切后得到的片段较多时，人工的计算片段的分子量较为耗时，因此这里特异根据《生物科学家的Python指南》一书，写了这篇文章以计算限制性核酸内切酶切割DNA分子后得到的片段分子量。
1.如何计算DNA片段的分子量？
这里用for循环写了一个函数，然后用一段DNA序列（
gcgatgctaggatccgc
）演示了函数的用法，具体代码如下：
#定义计算DNA片段的分子量大小的函数def oligoMolecularWeight(seqence):dnaMolecularWeight = {'a': 313.2, 'c': 289.2, 't': 304.2, 'g': 329.2}molecularWeight = 0.0for base in seqence:molecularWeight += dnaMolecularWeight[base]return molecularWeight#示例print(oligoMolecularWeight('gcgatgctaggatccgc'))
最后计算得到的结果为5273.399999999999。
2.find函数查找酶切位点在DNA上的位置，下面列了find函数运行后的结果：
mySequence = 'gcgatgctaggatccgcgatcgcgtacgatcgtacgcggtacggacggatccttctc'print(mySequence.find('ggatcc'))#运行结果：9print(mySequence.find('ggatttt'))#运行结果-1
由示例可以看出find函数在搜索字符串（‘ggatttt’） 时，因序列中无ggatttt，因此放回的是-1，表示未找到子字符串。
PS：find函数不能对“没有发现”返回零，因为这实际上是Python字符串中第一个位置的索引，如果返回零，则表示在该位置找到了这个子字符串。
3.while语句遍历整个DNA序列
mySequence = 'gcgatgctaggatccgcgatcgcgtacgatcgtacgcggtacggacggatccttctc'found = 0searchFrom = foundwhile found != -1:found = mySequence.find('ggatcc',searchFrom)if found != -1:print('Substring found at:',found)searchFrom = found +1=====================================================================运行结果：Substring found at: 9Substring found at: 46
4.合并计算分子量和酶切得到DNA片段的函数
#创建限制性核酸内切酶的切割位点restrictionEnzymes = {}restrictionEnzymes['bamH1'] = ['ggatcc',0]restrictionEnzymes['sma1'] = ['cccggg',2]#定义计算DNA片段的分子量大小的函数def oligoMolecularWeight(seqence):dnaMolecularWeight = {'a': 313.2, 'c': 289.2, 't': 304.2, 'g': 329.2}molecularWeight = 0.0for base in seqence:molecularWeight += dnaMolecularWeight[base]return molecularWeight#定义酶切函数def restrictionDigest(sequence,enzyme):#将限制性核酸内切酶的切割序列和切割位点分别赋值给motif和cutPositionmotif = restrictionEnzymes[enzyme][0]cutPosition = restrictionEnzymes[enzyme][1]#提前创建用来存储切割后得到的酶切片段序列及其分子量的列表fragments = []#found表示上次找到的限制性核酸内切酶位点序列的起始位置found = 0#lastCut用来存储上次切割DNA的位置lastCut = found#searchFrom用来告诉我们下次从哪开始找出酶切序列searchFrom = lastCut#while循环遍历DNA序列while found != -1:found = sequence.find(motif,searchFrom)if found != -1:fragment = sequence[lastCut:found+cutPosition]wnt = oligoMolecularWeight(fragment)fragments.append((fragment,wnt))else:fragment = sequence[lastCut:]wnt = oligoMolecularWeight(fragment)fragments.append((fragment,wnt))lastCut = found+cutPositionsearchFrom = lastCut+1return fragments#用法示例digestSquence = 'gcgatgctaggatccgcgatcgcgtacgatcgtacgcggtacggacggatccttctc'print(restrictionDigest(digestSquence,'bamH1'))
总结：大家只需要根据自己的需求添加所需的限制性核酸内切酶到字典，然后根据自己的序列进行运行就可以得到结果。
