Back to index

python-biopython  1.60
Functions
Bio.SeqUtils.MeltingTemp Namespace Reference

Functions

def Tm_staluc
def _test

Function Documentation

def Bio.SeqUtils.MeltingTemp._test ( ) [private]
Run the module's doctests (PRIVATE).

Definition at line 176 of file MeltingTemp.py.

00176 
00177 def _test():
00178     """Run the module's doctests (PRIVATE)."""
00179     import doctest
00180     print "Runing doctests..."
00181     doctest.testmod()
00182     print "Done"

def Bio.SeqUtils.MeltingTemp.Tm_staluc (   s,
  dnac = 50,
  saltc = 50,
  rna = 0 
)
Returns DNA/DNA tm using nearest neighbor thermodynamics.

dnac is DNA concentration [nM]
saltc is salt concentration [mM].
rna=0 is for DNA/DNA (default), use 1 for RNA/RNA hybridisation.

For DNA/DNA, see Allawi & SantaLucia (1997), Biochemistry 36: 10581-10594
For RNA/RNA, see Xia et al (1998), Biochemistry 37: 14719-14735

Example:

>>> print "%0.2f" % Tm_staluc('CAGTCAGTACGTACGTGTACTGCCGTA')
59.87
>>> print "%0.2f" % Tm_staluc('CAGTCAGTACGTACGTGTACTGCCGTA', rna=True)
68.14

You can also use a Seq object instead of a string,

>>> from Bio.Seq import Seq
>>> from Bio.Alphabet import generic_nucleotide
>>> s = Seq('CAGTCAGTACGTACGTGTACTGCCGTA', generic_nucleotide)
>>> print "%0.2f" % Tm_staluc(s)
59.87
>>> print "%0.2f" % Tm_staluc(s, rna=True)
68.14

Definition at line 10 of file MeltingTemp.py.

00010 
00011 def Tm_staluc(s,dnac=50,saltc=50,rna=0):
00012     """Returns DNA/DNA tm using nearest neighbor thermodynamics.
00013 
00014     dnac is DNA concentration [nM]
00015     saltc is salt concentration [mM].
00016     rna=0 is for DNA/DNA (default), use 1 for RNA/RNA hybridisation.
00017 
00018     For DNA/DNA, see Allawi & SantaLucia (1997), Biochemistry 36: 10581-10594
00019     For RNA/RNA, see Xia et al (1998), Biochemistry 37: 14719-14735
00020 
00021     Example:
00022 
00023     >>> print "%0.2f" % Tm_staluc('CAGTCAGTACGTACGTGTACTGCCGTA')
00024     59.87
00025     >>> print "%0.2f" % Tm_staluc('CAGTCAGTACGTACGTGTACTGCCGTA', rna=True)
00026     68.14
00027 
00028     You can also use a Seq object instead of a string,
00029 
00030     >>> from Bio.Seq import Seq
00031     >>> from Bio.Alphabet import generic_nucleotide
00032     >>> s = Seq('CAGTCAGTACGTACGTGTACTGCCGTA', generic_nucleotide)
00033     >>> print "%0.2f" % Tm_staluc(s)
00034     59.87
00035     >>> print "%0.2f" % Tm_staluc(s, rna=True)
00036     68.14
00037 
00038     """
00039 
00040     #Credits: 
00041     #Main author: Sebastian Bassi <sbassi@genesdigitales.com>
00042     #Overcount function: Greg Singer <singerg@tcd.ie>
00043     #Based on the work of Nicolas Le Novere <lenov@ebi.ac.uk> Bioinformatics.
00044     #17:1226-1227(2001)
00045 
00046     #This function returns better results than EMBOSS DAN because it uses
00047     #updated thermodynamics values and takes into account inicialization
00048     #parameters from the work of SantaLucia (1998).
00049     
00050     #Things to do:
00051     #+Detect complementary sequences. Change K according to result.
00052     #+Add support for heteroduplex (see Sugimoto et al. 1995).
00053     #+Correction for Mg2+. Now supports only monovalent ions.
00054     #+Put thermodinamics table in a external file for users to change at will
00055     #+Add support for danglings ends (see Le Novele. 2001) and mismatches.
00056     
00057     dh = 0 #DeltaH. Enthalpy
00058     ds = 0 #deltaS Entropy
00059 
00060     def tercorr(stri):
00061         deltah = 0
00062         deltas = 0
00063         if rna==0:
00064             #DNA/DNA
00065             #Allawi and SantaLucia (1997). Biochemistry 36 : 10581-10594
00066             if stri.startswith('G') or stri.startswith('C'):
00067                 deltah -= 0.1
00068                 deltas += 2.8
00069             elif stri.startswith('A') or stri.startswith('T'):
00070                 deltah -= 2.3
00071                 deltas -= 4.1
00072             if stri.endswith('G') or stri.endswith('C'):
00073                 deltah -= 0.1
00074                 deltas += 2.8
00075             elif stri.endswith('A') or stri.endswith('T'):
00076                 deltah -= 2.3
00077                 deltas -= 4.1
00078             dhL = dh + deltah
00079             dsL = ds + deltas
00080             return dsL,dhL
00081         elif rna==1:
00082             #RNA
00083             if stri.startswith('G') or stri.startswith('C'):
00084                 deltah -= 3.61
00085                 deltas -= 1.5
00086             elif stri.startswith('A') or stri.startswith('T') or \
00087                  stri.startswith('U'):
00088                 deltah -= 3.72
00089                 deltas += 10.5
00090             if stri.endswith('G') or stri.endswith('C'):
00091                 deltah -= 3.61
00092                 deltas -= 1.5
00093             elif stri.endswith('A') or stri.endswith('T') or \
00094                  stri.endswith('U'):
00095                 deltah -= 3.72
00096                 deltas += 10.5
00097             dhL = dh + deltah
00098             dsL = ds + deltas
00099             # print "delta h=",dhL
00100             return dsL,dhL
00101         else:
00102             raise ValueError("rna = %r not supported" % rna)
00103 
00104     def overcount(st,p):
00105         """Returns how many p are on st, works even for overlapping"""
00106         ocu = 0
00107         x = 0
00108         while True:
00109             try:
00110                 i = st.index(p,x)
00111             except ValueError:
00112                 break
00113             ocu += 1
00114             x = i + 1
00115         return ocu
00116 
00117     R = 1.987 # universal gas constant in Cal/degrees C*Mol
00118     sup = str(s).upper() #turn any Seq object into a string (need index method)
00119     vsTC, vh = tercorr(sup)
00120     vs = vsTC
00121     
00122     k = (dnac/4.0)*1e-9
00123     #With complementary check on, the 4.0 should be changed to a variable.
00124     
00125     if rna==0:
00126         #DNA/DNA
00127         #Allawi and SantaLucia (1997). Biochemistry 36 : 10581-10594
00128         vh = vh + (overcount(sup,"AA"))*7.9 + (overcount(sup,"TT"))*\
00129         7.9 + (overcount(sup,"AT"))*7.2 + (overcount(sup,"TA"))*7.2 \
00130         + (overcount(sup,"CA"))*8.5 + (overcount(sup,"TG"))*8.5 + \
00131         (overcount(sup,"GT"))*8.4 + (overcount(sup,"AC"))*8.4
00132         vh = vh + (overcount(sup,"CT"))*7.8+(overcount(sup,"AG"))*\
00133         7.8 + (overcount(sup,"GA"))*8.2 + (overcount(sup,"TC"))*8.2
00134         vh = vh + (overcount(sup,"CG"))*10.6+(overcount(sup,"GC"))*\
00135         9.8 + (overcount(sup,"GG"))*8 + (overcount(sup,"CC"))*8
00136         vs = vs + (overcount(sup,"AA"))*22.2+(overcount(sup,"TT"))*\
00137         22.2 + (overcount(sup,"AT"))*20.4 + (overcount(sup,"TA"))*21.3
00138         vs = vs + (overcount(sup,"CA"))*22.7+(overcount(sup,"TG"))*\
00139         22.7 + (overcount(sup,"GT"))*22.4 + (overcount(sup,"AC"))*22.4
00140         vs = vs + (overcount(sup,"CT"))*21.0+(overcount(sup,"AG"))*\
00141         21.0 + (overcount(sup,"GA"))*22.2 + (overcount(sup,"TC"))*22.2
00142         vs = vs + (overcount(sup,"CG"))*27.2+(overcount(sup,"GC"))*\
00143         24.4 + (overcount(sup,"GG"))*19.9 + (overcount(sup,"CC"))*19.9
00144         ds = vs
00145         dh = vh        
00146     elif rna==1:
00147         #RNA/RNA hybridisation of Xia et al (1998)
00148         #Biochemistry 37: 14719-14735         
00149         vh = vh+(overcount(sup,"AA"))*6.82+(overcount(sup,"TT"))*6.6+\
00150         (overcount(sup,"AT"))*9.38 + (overcount(sup,"TA"))*7.69+\
00151         (overcount(sup,"CA"))*10.44 + (overcount(sup,"TG"))*10.5+\
00152         (overcount(sup,"GT"))*11.4 + (overcount(sup,"AC"))*10.2
00153         vh = vh + (overcount(sup,"CT"))*10.48 + (overcount(sup,"AG"))\
00154         *7.6+(overcount(sup,"GA"))*12.44+(overcount(sup,"TC"))*13.3
00155         vh = vh + (overcount(sup,"CG"))*10.64 + (overcount(sup,"GC"))\
00156         *14.88+(overcount(sup,"GG"))*13.39+(overcount(sup,"CC"))*12.2
00157         vs = vs + (overcount(sup,"AA"))*19.0 + (overcount(sup,"TT"))*\
00158         18.4+(overcount(sup,"AT"))*26.7+(overcount(sup,"TA"))*20.5
00159         vs = vs + (overcount(sup,"CA"))*26.9 + (overcount(sup,"TG"))*\
00160         27.8 + (overcount(sup,"GT"))*29.5 + (overcount(sup,"AC"))*26.2
00161         vs = vs + (overcount(sup,"CT"))*27.1 + (overcount(sup,"AG"))*\
00162         19.2 + (overcount(sup,"GA"))*32.5 + (overcount(sup,"TC"))*35.5
00163         vs = vs + (overcount(sup,"CG"))*26.7 + (overcount(sup,"GC"))\
00164         *36.9 + (overcount(sup,"GG"))*32.7 + (overcount(sup,"CC"))*29.7
00165         ds = vs
00166         dh = vh
00167     else:
00168         raise ValueError("rna = %r not supported" %rna)
00169 
00170     ds = ds-0.368*(len(s)-1)*math.log(saltc/1e3)
00171     tm = ((1000* (-dh))/(-ds+(R * (math.log(k)))))-273.15
00172     # print "ds="+str(ds)
00173     # print "dh="+str(dh)
00174     return tm
00175