Difference between revisions of "User talk:Remig"

From Jmol
Jump to navigation Jump to search
(Alpha helix generater)
 
(Replaced content with "What I had here before belongs on the user page - so I moved it...")
 
(6 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>function get3from1(c) {
 
ret = ""
 
switch (c) {
 
case 'A':
 
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':
 
ret = "ASN";
 
break;
 
case 'P':
 
ret = "PRO";
 
break;
 
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
 
};
 
 
 
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
 
};
 
 
 
function genAA(i, aa, x) {
 
n = x
 
 
 
N0 = [0.0, 0.0, 0.0]
 
CA = [ 0.200, 1.174, 0.911 ]
 
C = [ -1.129, 1.783, 1.241 ]
 
O = [ -1.241, 1.967, 2.726 ]
 
CB = [ 1.062, 2.1950, 0.230 ]
 
G1 = [ 2.396, 1.588, -0.091 ]
 
G2 = [ 0.680, 3.652, 0.423]
 
Gfy = [ 2.368, 1.471, -0.0152 ]
 
Ghw = [ 2.406, 1.626, -0.134 ]
 
D1 = [ 3.225, 2.340, -1.096]
 
D1hw = [3.498, 1.936, 0.675]
 
D1fy = [ 3.346, 1.524, 0.921 ]
 
D2 = [ 3.189, 1.093, 1.087]
 
D2hw = [ 2.713, 0.901, -1.211 ]
 
D2fy = [ 2.493, 0.516, -1.151 ]
 
E1fy = [ 4.513, 0.615, 0.8244 ]
 
E1hw = [ 4.160, 0.518, -1.178 ]
 
E2fy = [ 3.528, -0.336, -1.206 ]
 
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 ]
 
Zfy = [ 4.588, -0.285, -0.168 ]
 
Hfy = [ 5.738, -1.245, -0.233 ]
 
N1 = [ -1.965, 2.307, 0.206 ]
 
Gp = [ 2.008, 1.24, -0.46 ]
 
Dp = [1.022, 0.213, -1.031 ]
 
E1 = [ 3.652, 3.503, -0.111 ]
 
E2 = [ 4.342, 1.591, -1.456 ]
 
Z = [ 4.115, 3.339, 1.403 ]
 
H1 = [4.087, 4.572, 2.139]
 
H2 = [5.469, 2.866, 1.296]
 
 
a3 = get3from1(seq[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)
 
}
 
 
 
switch (aa) {
 
case 'A' :
 
break;
 
case 'C' :
 
a += genAtom(n++, "SG ", a3, i, G1)
 
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;
 
default :
 
break;
 
}
 
 
 
return a
 
}
 
 
 
# GenAlph
 
function genAlpha(seq, pa) {
 
 
 
# For each aa
 
set appendnew false
 
n = 1
 
pn = 0
 
pc= 0 # previous C
 
for (var i = 1; i <= seq.count; i++) {
 
 
 
pn = n
 
m = pn + 1
 
 
 
# Anticipate PRO
 
dpsi = 0
 
dphi = 0
 
if ((seq.count - i) >= 4) {
 
if (seq[i + 4] == 'P') {
 
dpsi = =9
 
}
 
}
 
if ((seq.count - i) >= 3) {
 
if (seq[i + 3] == 'P') {
 
dphi = 10
 
}
 
}
 
if ((seq.count - i) >= 2) {
 
if (seq[i + 2] == 'P') {
 
dpsi = -11
 
#dphi = -0
 
}
 
}
 
if ((seq.count - i) >= 1) {
 
if (seq[i + 1] == 'P') {
 
dpsi = 25
 
#dphi = -81
 
}
 
}
 
 
# Move polypeptide C to bond distance to next AA N
 
select all
 
fix none
 
translateselected { 1.955, -1.993, 0.000 }
 
 
 
# Gen AA
 
a = "data \"append aa\"\n"
 
a += genAA(i, seq[i], n);
 
a += "end \"append aa\""
 
script inline @{a}
 
 
 
# If not first AA
 
if (pc > 0) {
 
v1={atomno=pc}.xyz - {atomno=pn}.xyz
 
v2={atomno=m}.xyz - {atomno=pn}.xyz
 
# Gen normalized axis perpendicular to the plane contaoning atoms pc, pn and m
 
axis = cross(v1, v2) - {atomno=pn}.xyz
 
 
 
# Center on atom pn
 
axis += {atomno=pn}.xyz
 
 
 
# Rotate the polypeptide on the new AA
 
select atomno<pn
 
fix atomno>=pn
 
 
 
# pull peptide  angle trigonal
 
rotateselected @axis {atomno=pn} @{pa-62}
 
 
 
# Make peptide bond dihedral 180
 
rotateselected {atomno=pc} {atomno=pn} 52
 
 
 
# Make the phi -60
 
rotateselected {atomno=pn} {atomno=m} @{dphi - 68}
 
 
 
# Make the psi -45
 
pca = pc - 1
 
fix atomno>pca
 
select atomno<pn and (atomno != @{pc+1})
 
rotateselected {atomno=pc} {atomno=pca} @{dpsi}
 
 
 
# 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}} -120
 
}
 
 
 
}
 
pc =pn + 2
 
 
 
connect
 
}
 
select all
 
fix none
 
}
 
 
 
 
 
echo Generating Alpha Helix
 
zap
 
 
 
# Get the sequence from the user
 
seq = prompt("Enter AA sequence (1 char coded)", "")%9999
 
pa = prompt("Enter the peptide angle", "120")
 
print format ("seq=%s at %d deg", seq, pa)
 
genAlpha(seq , pa)
 
 
 
 
 
 
 
</pre>
 

Latest revision as of 21:32, 31 October 2013

What I had here before belongs on the user page - so I moved it...

Contributors

Remig