|
|
| (5 intermediate revisions by the same user not shown) |
| Line 1: |
Line 1: |
| − | I use Jmol to study protein folding. Here is a script I wrote that accepts an amino acid sequence (1 letter encoding: "AAC...FYW" for example) and generates an alpha helix using the Model Kit: | + | What I had here before belongs on the user page - so I moved it... |
| − | <pre># RIBOZOME - Jmol script by Ron Mignery
| |
| − | # v1.0 beta 10/19/2013
| |
| − | #
| |
| − | # RIBOZOME takes a string message encoding an amino acid (aa) sequence
| |
| − | # and generates a corresponding alpha helix one aa at a time from the
| |
| − | # N terminus to the C terminus rotating the emerging helix as it goes.
| |
| − | #
| |
| − | # The message is a string entered by the user at a prompt.
| |
| − | # It may be typed in or pasted in and be of any length
| |
| − | #
| |
| − | # The IUPAC/IUBMB 1 letter code is used:
| |
| − | # A=ALAnine B=GLUtamine?* C=CYSteine D=ASPartate E=GLUtamate
| |
| − | # F=PHEnylalanine G=GLYcine H=HIStidine I=IsoLEucine K=LYSine
| |
| − | # L=LEUcine M=METhionine N=ASparagiNe O=PYrroLysine*** P=PROline
| |
| − | # Q=GLutamiNe R=ARGinine S=SERine T=THReonine U=SElenoCysteine
| |
| − | # V=VALine W=TRyPtophan X=UNKnown**** Y=TYRosine Z=ASparagiNe?**
| |
| − | # *Either GLU or GLN: drawn as GLN with chi3 flipped
| |
| − | # **Either ASP or ASN: drawn as ASN with chi3 flipped
| |
| − | # ***Not supported
| |
| − | # ****Unknown aa will be drawn as ALA
| |
| − | #
| |
| − | | |
| − | # The following constant values determine the pitch of the alpha helix
| |
| − | var PHI = -50 # Dihedral angle of N-CA bond (nominally -50)
| |
| − | var PSI = -60 # Dihedral angle of CA-C bond (nominally -60)
| |
| − | var OMEGA = 180 # Dihedral angle of the peptide bond (nominally 180)
| |
| − | var PEPTIDE_ANGLE = 110 # C-N-CA angle (nominally 110)
| |
| − | var PRO_BUMP = -15 # Psi angle change to N-3psi when N is Pro
| |
| − | | |
| − | # Lookup 3 letter code from 1 letter code
| |
| − | function get3from1(c) {
| |
| − | ret = ""
| |
| − | switch (c) {
| |
| − | case 'A':
| |
| − | case 'X':
| |
| − | ret = "ALA";
| |
| − | break;
| |
| − | case 'C':
| |
| − | ret = "CYS";
| |
| − | break;
| |
| − | case 'D':
| |
| − | ret = "ASP";
| |
| − | break;
| |
| − | case 'E':
| |
| − | ret = "GLU";
| |
| − | break;
| |
| − | case 'F':
| |
| − | ret = "PHE";
| |
| − | break;
| |
| − | case 'G':
| |
| − | ret = "GLY";
| |
| − | break;
| |
| − | case 'H':
| |
| − | ret = "HIS";
| |
| − | break;
| |
| − | case 'I':
| |
| − | ret = "ILE";
| |
| − | break;
| |
| − | case 'K':
| |
| − | ret = "LYS";
| |
| − | break;
| |
| − | case 'L':
| |
| − | ret = "LEU";
| |
| − | break;
| |
| − | case 'M':
| |
| − | ret = "MET";
| |
| − | break;
| |
| − | case 'N':
| |
| − | case 'Z':
| |
| − | ret = "ASN";
| |
| − | break;
| |
| − | case 'P':
| |
| − | ret = "PRO";
| |
| − | break;
| |
| − | case 'B':
| |
| − | case 'Q':
| |
| − | ret = "GLN";
| |
| − | break;
| |
| − | case 'R':
| |
| − | ret = "ARG";
| |
| − | break;
| |
| − | case 'S':
| |
| − | ret = "SER";
| |
| − | break;
| |
| − | case 'T':
| |
| − | ret = "THR";
| |
| − | break;
| |
| − | case 'U':
| |
| − | ret = "SEC";
| |
| − | break;
| |
| − | case 'V':
| |
| − | ret = "VAL";
| |
| − | break;
| |
| − | case 'W':
| |
| − | ret = "TRP";
| |
| − | break;
| |
| − | case 'Y':
| |
| − | ret = "TYR";
| |
| − | break;
| |
| − | }
| |
| − | return ret
| |
| − | };
| |
| − | | |
| − | # Generate PDM atom record
| |
| − | function genAtom(n, e, aa, i, xyz) {
| |
| − | a = format("ATOM %5d %3s %3s A", n, e, aa )
| |
| − | a += format("%4d %8.3f", i, xyz[1] )
| |
| − | a += format("%8.3f%8.3f\n", xyz[2], xyz[3] )
| |
| − | return a
| |
| − | };
| |
| − | | |
| − | # Generate an amino acid record set
| |
| − | function genAA(i, aa, x) {
| |
| − | n = x
| |
| − | | |
| − | # From constructed AAs
| |
| − | N0 = [0.0, 0.0, 0.0]
| |
| − | CA = [ 0.200, 1.174, 0.911 ]
| |
| − | C = [ -1.110, 1.668, 1.425 ] #[ -1.129, 1.783, 1.241 ]
| |
| − | O = [ -1.320, 1.693, 2.62 ] #[ -1.241, 1.967, 2.726 ]
| |
| − | CB = [ 1.062, 2.1950, 0.230 ]
| |
| − |
| |
| − | G1 = [ 2.359, 1.607, -0.344] #2.396, 1.588, -0.091 ]
| |
| − | G2 = [ 1.363, 3.336, 1.157 ] #0.680, 3.652, 0.423]
| |
| − | D1 = [ 3.222, 2.656, -1.048 ] #[ 3.225, 2.340, -1.096]
| |
| − | D2 = [ 3.143, 0.904, 0.725 ] #[ 3.189, 1.093, 1.087]
| |
| − | E1 = [ 3.645, 3.749, -0.167 ] #[ 3.652, 3.503, -0.111 ]
| |
| − | E2 = [ 2.491, 3.234, -2.249 ] #[ 4.342, 1.591, -1.456 ]
| |
| − | Z = [ 4.470, 4.717, -0.885 ] #[ 4.115, 3.339, 1.403 ]
| |
| − | H1 = [ 4.450, 6.006, -0.220 ] #[4.087, 4.572, 2.139]
| |
| − | H2 = [5.833, 4.228, -0.984 ] #[5.469, 2.866, 1.296]
| |
| − |
| |
| − | Gp = [ 2.008, 1.24, -0.46 ]
| |
| − | Dp = [1.022, 0.213, -1.031 ]
| |
| − |
| |
| − | Gfy = [ 2.368, 1.471, -0.0152 ]
| |
| − | D1fy = [ 3.346, 1.524, 0.921 ]
| |
| − | D2fy = [ 2.493, 0.516, -1.151 ]
| |
| − | E1fy = [ 4.513, 0.615, 0.8244 ]
| |
| − | E2fy = [ 3.528, -0.336, -1.206 ]
| |
| − | Zfy = [ 4.588, -0.285, -0.168 ]
| |
| − | Hfy = [ 5.738, -1.245, -0.233 ]
| |
| − |
| |
| − | Ghw = [ 2.406, 1.626, -0.134 ]
| |
| − | D1hw = [3.498, 1.936, 0.675]
| |
| − | D2hw = [ 2.713, 0.901, -1.211 ]
| |
| − | E1hw = [ 4.160, 0.518, -1.178 ]
| |
| − | E2hw = [ 4.622, 1.160, 0.0816 ]
| |
| − | E3hw = [ 3.789, 2.523, 1.944 ]
| |
| − | Z2hw = [ 5.973, 1.177, 0.689 ]
| |
| − | Z3hw = [ 5.014, 2.550, 2.503 ]
| |
| − | H2hw = [ 6.153, 1.846, 1.844 ]
| |
| − |
| |
| − | N1 = [ 2.069, -2.122, -0.554] #[ -1.965, 2.307, 0.206 ]
| |
| − |
| |
| − | # Build PDB atom records common to all AAs
| |
| − | a3 = get3from1(seq[i])
| |
| − | if (a3 == "") {
| |
| − | a3 = "UNK"
| |
| − | }
| |
| − | print format("+ %s%d", a3, i)
| |
| − | a = genAtom(n++, "N ", a3, i, N0)
| |
| − | a += genAtom(n++, "CA ", a3, i, CA)
| |
| − | a += genAtom(n++, "C ", a3, i, C)
| |
| − | a += genAtom(n++, "O ", a3, i, O)
| |
| − | if (seq[i] != 'G') {
| |
| − | a += genAtom(n++, "CB ", a3, i, CB)
| |
| − | }
| |
| − | | |
| − | # Now add AA specific atom records
| |
| − | switch (aa) {
| |
| − | case 'A' :
| |
| − | case 'X' :
| |
| − | break;
| |
| − | case 'B' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "CD ", a3, i, D1)
| |
| − | a += genAtom(n++, "OE1", a3, i, E2) # GLN with Es switched
| |
| − | a += genAtom(n++, "NE2", a3, i, E1)
| |
| − | break;
| |
| − | case 'C' :
| |
| − | a += genAtom(n++, "SG ", a3, i, G2)
| |
| − | break;
| |
| − | case 'D' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "OD1", a3, i, D1)
| |
| − | a += genAtom(n++, "OD2", a3, i, D2)
| |
| − | break;
| |
| − | case 'E' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "CD ", a3, i, D1)
| |
| − | a += genAtom(n++, "OE1", a3, i, E1)
| |
| − | a += genAtom(n++, "OE2", a3, i, E2)
| |
| − | break;
| |
| − | case 'F' :
| |
| − | a += genAtom(n++, "CG ", a3, i, Gfy)
| |
| − | a += genAtom(n++, "CD1", a3, i, D1fy)
| |
| − | a += genAtom(n++, "CD2", a3, i, D2fy)
| |
| − | a += genAtom(n++, "CE1", a3, i, E1fy)
| |
| − | a += genAtom(n++, "CE2", a3, i, E2fy)
| |
| − | a += genAtom(n++, "CZ ", a3, i, Zfy)
| |
| − | break;
| |
| − | case 'G' :
| |
| − | break;
| |
| − | case 'H' :
| |
| − | a += genAtom(n++, "CG ", a3, i, Ghw)
| |
| − | a += genAtom(n++, "ND1", a3, i, D1hw)
| |
| − | a += genAtom(n++, "CD2", a3, i, D2hw)
| |
| − | a += genAtom(n++, "CE1", a3, i, E2hw)
| |
| − | a += genAtom(n++, "NE2", a3, i, E1hw)
| |
| − | break;
| |
| − | case 'I' :
| |
| − | a += genAtom(n++, "CG1", a3, i, G1)
| |
| − | a += genAtom(n++, "CG2", a3, i, G2)
| |
| − | a += genAtom(n++, "CD1", a3, i, D1)
| |
| − | break;
| |
| − | case 'K' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "CD ", a3, i, D1)
| |
| − | a += genAtom(n++, "CE ", a3, i, E1)
| |
| − | a += genAtom(n++, "NZ ", a3, i, Z)
| |
| − | break;
| |
| − | case 'L' :
| |
| − | a += genAtom(n++, "CG1", a3, i, G1)
| |
| − | a += genAtom(n++, "CD1", a3, i, D1)
| |
| − | a += genAtom(n++, "CD2", a3, i, D2)
| |
| − | break;
| |
| − | case 'M' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "SD ", a3, i, D1)
| |
| − | a += genAtom(n++, "CE ", a3, i, E1)
| |
| − | break;
| |
| − | case 'N' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "OD1", a3, i, D1)
| |
| − | a += genAtom(n++, "ND2", a3, i, D2)
| |
| − | break;
| |
| − | case 'P' :
| |
| − | a += genAtom(n++, "CG ", a3, i, GP)
| |
| − | a += genAtom(n++, "CD ", a3, i, DP)
| |
| − | break;
| |
| − | case 'Q' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "CD ", a3, i, D1)
| |
| − | a += genAtom(n++, "OE1", a3, i, E1)
| |
| − | a += genAtom(n++, "NE2", a3, i, E2)
| |
| − | break;
| |
| − | case 'R' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "CD ", a3, i, D1)
| |
| − | a += genAtom(n++, "NE ", a3, i, E1)
| |
| − | a += genAtom(n++, "CZ ", a3, i, Z)
| |
| − | a += genAtom(n++, "NH1", a3, i, H1)
| |
| − | a += genAtom(n++, "NH2", a3, i, H2)
| |
| − | break;
| |
| − | case 'S' :
| |
| − | a += genAtom(n++, "OG ", a3, i, G1)
| |
| − | break;
| |
| − | case 'T' :
| |
| − | a += genAtom(n++, "OG1", a3, i, G1)
| |
| − | a += genAtom(n++, "CG2", a3, i, G2)
| |
| − | break;
| |
| − | case 'U' :
| |
| − | a += genAtom(n++, "SeG", a3, i, G1)
| |
| − | break;
| |
| − | case 'V' :
| |
| − | a += genAtom(n++, "CG1", a3, i, G1)
| |
| − | a += genAtom(n++, "CG2", a3, i, G2)
| |
| − | break;
| |
| − | case 'W' :
| |
| − | a += genAtom(n++, "CG ", a3, i, Ghw)
| |
| − | a += genAtom(n++, "CD1", a3, i, D1hw)
| |
| − | a += genAtom(n++, "CD2", a3, i, D2hw)
| |
| − | a += genAtom(n++, "CE2", a3, i, E2hw)
| |
| − | a += genAtom(n++, "NE1", a3, i, E1hw)
| |
| − | a += genAtom(n++, "CE3", a3, i, E3hw)
| |
| − | a += genAtom(n++, "CZ2", a3, i, Z2hw)
| |
| − | a += genAtom(n++, "CZ3", a3, i, Z3hw)
| |
| − | a += genAtom(n++, "CH2", a3, i, H2hw)
| |
| − | break;
| |
| − | case 'Y' :
| |
| − | a += genAtom(n++, "CG ", a3, i, Gfy)
| |
| − | a += genAtom(n++, "CD1", a3, i, D1fy)
| |
| − | a += genAtom(n++, "CD2", a3, i, D2fy)
| |
| − | a += genAtom(n++, "CE1", a3, i, E1fy)
| |
| − | a += genAtom(n++, "CE2", a3, i, E2fy)
| |
| − | a += genAtom(n++, "CZ ", a3, i, Zfy)
| |
| − | a += genAtom(n++, "OH ", a3, i, Hfy)
| |
| − | break;
| |
| − | case 'Z' :
| |
| − | a += genAtom(n++, "CG ", a3, i, G1)
| |
| − | a += genAtom(n++, "OD1", a3, i, D2) # ASN with Ds switched
| |
| − | a += genAtom(n++, "ND2", a3, i, D1)
| |
| − | break;
| |
| − | default :
| |
| − | break;
| |
| − | }
| |
| − | | |
| − | return a
| |
| − | };
| |
| − | | |
| − | function rotateNward (a1, a2, angle) {
| |
| − | select atomno<a1
| |
| − | fix atomno>=a2
| |
| − | rotateselected {atomno=a1} {atomno=a2} @angle
| |
| − | #print format("a1=%d a2=%d angle=%d", a1, a2, angle) #DEBUG
| |
| − | | |
| − | };
| |
| − | | |
| − | # GenAlph
| |
| − | function genAlpha(seq, PHI, PSI, OMEGA, PEPTIDE_ANGLE, PRO_BUMP) {
| |
| − | | |
| − | # For each aa
| |
| − | set appendnew false
| |
| − | n = 1
| |
| − | pn = 0
| |
| − | pc= 0 # previous C
| |
| − | ca1 = 0; ca2 = 0; ca3 = 0
| |
| − | for (var i = 1; i <= seq.count; i++) {
| |
| − | | |
| − | # Step previous N
| |
| − | pn = n
| |
| − | | |
| − | # Move polypeptide C to bond distance from new AA N
| |
| − | select all
| |
| − | fix none
| |
| − | translateselected {2.069, -2.122, -0.554 } #N1
| |
| − | | |
| − | # Gen AA
| |
| − | a = "data \"append aa\"\n"
| |
| − | a += genAA(i, seq[i], n);
| |
| − | a += "end \"append aa\""
| |
| − | script inline @{a}
| |
| − | | |
| − | # If PRO
| |
| − | # Adjust i-3psi as rigid PRO D4 bumps Oi-4 as i-3psi
| |
| − | # is the most torqued psi and psis are the easiest to twist
| |
| − | if (seq[i] == 'P') {
| |
| − | rotateNward (ca3, ca3+1, PRO_BUMP)
| |
| − | }
| |
| − | | |
| − | # If not first AA
| |
| − | if (pc > 0) {
| |
| − | | |
| − | # Gen axis on previous n perpendicular to the plane
| |
| − | # containing atoms pc, pn and pn+1
| |
| − | v1={atomno=pc}.xyz - {atomno=pn}.xyz
| |
| − | v2={atomno=@{pn+1}}.xyz - {atomno=pn}.xyz
| |
| − | axis = cross(v1, v2)
| |
| − |
| |
| − | # Center on atom previous n
| |
| − | axis += {atomno=pn}.xyz
| |
| − | | |
| − | # Rotate the polypeptide N on the new AA C to the
| |
| − | # desired angle (nominally 110)
| |
| − | select atomno<pn
| |
| − | fix atomno>=pn
| |
| − | rotateselected @axis {atomno=pn} @{PEPTIDE_ANGLE - 77}
| |
| − | | |
| − | # Make omega dihedral OMEGA (nominally 180)
| |
| − | rotateselected {atomno=pc} {atomno=pn} @{OMEGA - 123}
| |
| − | | |
| − | # Make the phi PHI (nominally -50)
| |
| − | rotateselected {atomno=pn} {atomno=@{pn+1}} @{PHI - 30}
| |
| − | | |
| − | # Make the psi PSI (nominally -60)
| |
| − | fix atomno>pca
| |
| − | select atomno<pn and (atomno != @{pc+1})
| |
| − | rotateselected {atomno=pc} {atomno=@{pc-1}} @{PSI + 60}
| |
| − | | |
| − | # If aromatic go trans on chi 1
| |
| − | select atomno>@{pn+4} and atomno<n
| |
| − | if ((seq[i] == 'H') || (seq[i] == 'W') || (seq[i] == 'F') || (seq[i] == 'Y')) {
| |
| − | #rotateselected {atomno=@{pn+1}} {atomno=@{pn+4}} -60
| |
| − | }
| |
| − | | |
| − | # Save last three CAs for proline bumps
| |
| − | ca3 = ca2; ca2 = ca1; ca1 = pn + 1
| |
| − | }
| |
| − |
| |
| − | # Step previous C
| |
| − | pc =pn + 2
| |
| − | | |
| − | # Make the peptide bond
| |
| − | connect
| |
| − | }
| |
| − |
| |
| − | # Clean up
| |
| − | select all
| |
| − | fix none
| |
| − | }
| |
| − | | |
| − | echo Generating Alpha Helix
| |
| − | | |
| − | # Get the sequence from the user
| |
| − | seq = prompt("*** Any existing will be cleared ***\nEnter AA sequence (1 char coded)", "")%9999
| |
| − | if (seq.count > 0) {
| |
| − | zap # disable to keep existing structures
| |
| − | print format ("seq=%s phi=%d psi=%d", seq, PHI, PSI)
| |
| − | print format ("phi=%d psi=%d peptide angle=%d pro bump=%d", PEPTIDE_ANGLE, PRO_BUMP)
| |
| − | genAlpha(seq, PHI,PSI, OMEGA, PEPTIDE_ANGLE, PRO_BUMP) # defined at top of scripts
| |
| − | }
| |
| − | </pre>
| |
