Code
--iwdn worm
--version 1.2 30.04.2021
--Script to perform local wiggles/mutates/shakes along the protein.
--What it does:
--Performs local wiggle/mutate/shake along a protein.
--If the user has selected residues before this recipe is run, then the worm will only run along these segments. Otherwise all non-locked segments will be treated.
--The advantage of this worm-implementation is that it will not only act on the sidechains on the protein but also on adjacent sidechains for example on a target-binder. Set this via the "Radius" slider.
--Also, score- and objective-bonus-loss is monitored and it is reverted to a previous state so that no score- and especially objective-bonus-loss is possible.
--This recipe can sometimes gain many points in short time while keeping the structure intact.
--Just play around with it and see what it can do for you.
--Happy folding to all!
--version history
--1.0:
--Initial version. Only for internal use. Imlementation of basic functionality.
--1.1:
--Brushed up for release. Changed a few settings.
--1.2:
--Added options in the GUI to choose any combination of Objectives to monitor for the worm in terms of bonus-/socre-degradation. Like this the user can choose which target to aim for in terms of optimization.
--Added option to choose the qs-slot in the GUI.
--Better handling of ending/errors. The initial solution and the best/last found valid worm-solution are loaded to Undo-stack. Use ctrl+y/z to flip through.
--Added some more text-output (summary of GUI-options at start etc.)
--Added option to wiggle only sidechains. This will leave your secondary structure intact.
--Added option to use more sophisticated distance-check to beta-carbon. This is mandatory for symmetry-puzzles since the symmetric copies are considered as well.
--Improved distance-check so that only those segments are considered which are actually useful for the selected mode (w/m/s).
--Cleaned up the code in general
--Variables and Constants
local version = "1.2"
local dia = {} --Global GUI
local toWorkOn = {}
local toWorkOnIdx = 1
local qs_idx_start = 33 --Index to start for quicksave slots
local qs_idx = qs_idx_start
local guardInt01 = 1
local selectedSegments = {}
local residues = structure.GetCount()
local residues_w = {} --Only these residues need to be considered in distance-checks for WIGGLE.
local residues_m = {} --Only these residues need to be considered in distance-checks for MUTATION.
local residues_s = {} --Only these residues need to be considered in distance-checks for SHAKE.
--Functions
function wait(seconds)
local start = os.time()
repeat until os.time() > start + seconds
end
--For DEBUG only:
--selection.DeselectAll()
--for m=1,#residues_m do
--local testSeg = residues_m[m]
--print("Segment " .. testSeg .. " is mutable: " .. tostring(structure.IsMutable(testSeg)))
--selection.Select(testSeg)
--end
--wait(1000)
function round(inNum)
return inNum-inNum%0.001
end
function findToWorkOn()
--Determine which sidechains shall be processed. The user needs to select them before running this recipe!
for m=1,residues do
if (selection.IsSelected(m) == true) then
table.insert(toWorkOn, m)
end
end
if (#toWorkOn > 0) then
selectedSegments = toWorkOn
end
--print("Selected segments1: " .. #selectedSegments)
--Now THIS is strange!
--If I would just put in: "selectedSegments = toWorkOn" here (without the if), then this would act like a pointer! So if toWorkOn changes afterwards, then the contents of selectedSegments will change to this as well even though it would not be reassigned again!
--Investigate if this is a special Lua-subject. If not then this can lead to all kinds of errors which may explain constant crashes.
--Seems like this is standard behavior of Lua. So I will have to be careful with this in the future.
--If user has not selected anything beforehand then choose all NON-LOCKED residues.
if (#toWorkOn == 0) then
for m=1,residues do
if (structure.IsLocked(m) == false) then
table.insert(toWorkOn, m)
end
end
end
--If still toWorkOn is empty then this is an error.
if (#toWorkOn == 0) then
print("Error. No residues to work on could be found!")
return
end
end
function mysplit(inputstr, sep)
if sep == nil then
sep = "%s"
end
local t={}
for str in string.gmatch(inputstr, "([^"..sep.."]+)") do
table.insert(t, str)
end
return t
end
function checkFilters(cf_toMonitor)
--Check the objectives and store for later reference.
--The input table cf_toMonitor is a list of filter-names to monitor. A special case is "Global Score" which is not a filter but the overall score.
--The output are the actual points/bonus for each filter.
local cf_toMonitor_pts = {}
for m=1,#cf_toMonitor do
if (string.find(cf_toMonitor[m], "Global Score") ~= nil) then
--Treatment of Global Score.
cf_toMonitor_pts[m] = current.GetEnergyScore() --This does of course contain the bonuses of all filters as well.
else
--This means that we shall check an actual filter
cf_toMonitor_pts[m] = filter.GetBonus(cf_toMonitor[m])
end
end
return cf_toMonitor_pts
end
function getDistances(gd_testArr, gd_center)
--Determines reliably the minimal distances between the sidechain of the center-element and ALL others
--Distance refers to the beta-carbon. Only bands can be used in symmetry puzzles for this since structure.GetDistance() cannot be used for symmetric copies.
--For symmetry puzzles, only the minimum distance from center to any of the copies is output. So only one value for each residue in gd_testArr.
local gd_symCount = structure.GetSymCount() --0...
--local tmpBdIdxs = {}
--local tmpBdsIdxs = {}
local tmpDists = {}
--local tmpDist2 = {}
local outArr = {}
--local tmpBdIdx_last = {}
for m=1,#gd_testArr do
if (gd_testArr[m] == gd_center) then
outArr[m] = -1 --Just set an illegal value if center is in gd_testArr so that it will be ignored.
else
local tmpBdIdxs = {}
local tmpDists = {}
for n=1,gd_symCount+1 do
local tmpBdIdx = band.AddBetweenSegments(gd_center, gd_testArr[m], 5, 5, n-1) --Distance-check to beta-carbon. It seems like this will also work well for Glycine even if it has no beta-carbon. It will fall-back to alpha-carbon?! Check this!
table.insert(tmpBdIdxs, tmpBdIdx)
local tmpDist = band.GetLength(tmpBdIdx)
table.insert(tmpDists, tmpDist)
end
--Now all bands from center to the symmetric copies of the current segment are drawn.
--Determine the shortest one and output this as distance.
local tmpMinDist = tmpDists[1]
local tmpMinDistIdx = 1
for n=2,#tmpDists do
if (tmpDists[n] < tmpMinDist) then
tmpMinDist = tmpDists[n]
tmpMinDistIdx = n
end
end
--local tmpString = ""
--for n=1,#tmpDists do
-- tmpString = tmpString .. round(tmpDists[n]) .. " "
--end
--wait(5) --for debug
--Delete all bands. Somehow this needs to be done in reverse since the indices need to be consecutive?
for n=#tmpDists,1,-1 do
band.Delete(tmpBdIdxs[n])
end
--Debug-check:
--local tmpBdIdx = band.AddBetweenSegments(gd_center, gd_testArr[m], 5, 5, tmpMinDistIdx-1) --Distance-check to beta-carbon. Attention! Glycine will be ignored here!
--print("Seg " .. gd_testArr[m] .. " dists: " .. tmpString .. " taken: " .. round(band.GetLength(tmpBdIdx)))
--wait(5) --for debug
--band.Delete(tmpBdIdx)
outArr[m] = tmpDists[tmpMinDistIdx]
end
end
return outArr
end
function worm(w_toWorkOn, w_modus, w_startIdx, w_step, w_radius, w_iters, w_toMonitor)
--Instead of wiggling globally, this function wiggles locally multiple times at different positions.
--Objectives are monitored so that none degrades. Especially core and loops. Uses qs33.
--Changed this so that this function always only does 1 turn. If you want to run multiple turns then just call this function again.
--w_prot=0: unprotected. 1: Objectives only are monitored. 2: Score is monitored. 3: Objectives and score are monitored.
local w_startFilters = checkFilters(w_toMonitor)
--Store solution for later restore
save.Quicksave(qs_idx) --This incorporates the selection of the segments to work on.
--Start the worm-loop
--Go through the residues to process and apply worm-action in the vicinity.
local w_currCenter = nil
local w_currVici = {}
for n=w_startIdx,#w_toWorkOn,w_step do
w_currCenter = w_toWorkOn[n]
w_currVici = {}
--Determine which segments to use for distance-check.
--This depends on the actual mode that is used.
local w_toCheck = {}
if (w_modus == 1 or w_modus == 10) then --wiggle
w_toCheck = residues_w
elseif (w_modus == 2) then --mutate
w_toCheck = residues_m
else --shake
w_toCheck = residues_s
end
--Determine which segments are in the vicinity of the center
if (dia.cbx03c.value == false) then --Standard distance-check. May get incorrect results on symmetry-puzzles since the symmetric copies are not considered because GetDistance() cannot do this.
w_currVici[1] = w_currCenter
for o=1,#w_toCheck do
if (w_toCheck[o] ~= w_currCenter) then
w_currDist = structure.GetDistance(w_currCenter, w_toCheck[o])
if (w_currDist <= w_radius) then
table.insert(w_currVici, w_toCheck[o])
end
end
end --End: Determine vicinity
else
--More sophisticated distance-check with bands. Mandatory for symmetry-puzzles!
local tmpDistances = getDistances(w_toCheck, w_currCenter)
for o=1,#tmpDistances do
if (tmpDistances[o] <= w_radius) then
table.insert(w_currVici, w_toCheck[o])
end
end
end
--Now apply the desired operation onto all segments in the vicinity.
selection.DeselectAll()
for o=1,#w_currVici do
selection.Select(w_currVici[o]) --select residues in the vicinity
end
if (w_modus == 1) then --wiggle
structure.LocalWiggleSelected(w_iters, true, true) --backbone, sidechains
elseif (w_modus == 10) then --wiggle sidechains only
structure.LocalWiggleSelected(w_iters, false, true) --backbone, sidechains
elseif (w_modus == 2) then --mutate
structure.MutateSidechainsSelected(w_iters)
else --shake
structure.ShakeSidechainsSelected(w_iters)
end
--Restore previous solution if ANY of the filters in w_toMonitor decreases.
local w_currFilters = checkFilters(w_toMonitor)
local w_doRestore = false
for o=1,#w_currFilters do
if (w_currFilters[o] < w_startFilters[o]) then
w_doRestore = true
break
end
end
--if w_doRestore == false at this point, then we have a new, better solution.
--This result is updated so that in further steps we cannot loose this anymore.
if (w_doRestore == false) then
w_startFilters = w_currFilters
save.Quicksave(qs_idx)
else
save.Quickload(qs_idx)
end
end --End: Main worm-loop
end --End: worm
function cleanup(err)
print("")
if (string.find(string.lower(err), "cancelled") ~= nil) then
print("User Cancelled. Cleaning up...")
print("Loading initial and recent best solution to Undo-stack...")
undo.SetUndo(true)
for n=1,qs_idx-qs_idx_start+1 do
save.Quickload(qs_idx_start+n-1) --Puts found solutions into restore. Use ctrl+z/y to flip through.
wait(guardInt01) --This seems to be necessary in order for things to settle.
print("Solution in qs" .. qs_idx_start+n-1 .. " (" .. n .. "/" .. qs_idx-qs_idx_start+1 .. ") loaded.")
end
--Restore initial selection. May come in handy for further actions.
selection.DeselectAll()
for o=1,#selectedSegments do
selection.Select(selectedSegments[o])
end
print("iwdn worm version " .. version .. " has finished!")
else
print("ERROR: " .. err)
print("iwdn worm version " .. version .. " has finished with an error!")
end
end
function main()
print("Running iwdn worm version " .. version .. "...")
print("Initial Score is: " .. round(current.GetEnergyScore()) .. " Pts.")
--Identify possible segments to work on. Exclude fully locked segments (bb and sc) etc.
print("Identifying possible segments to work on (please wait)...")
for m=1,residues do
--If the segment is mutable, then it needs to be considered in every case since then it will be shakeable and wiggleable as well. In linker- or interface-puzzles there are typically segments which have locked backbone but are still mutable. These will be incorporated here as well.
--Interesting: On binder-puzzles, IsMutable() seems to work correctly. the residues on the binder are found as false. I remember that I had problems with it on other puzzles. Or was it "CanMutate()"?
--But does this work as well for linker or interface where the sidechains are mutable? Yes, for linker it does!
if (structure.IsMutable(m) == true) then
--segments belonging to this category: general protein-segments, locked or non-locked mutable segments on linker- or interface-targets
table.insert(residues_w, m)
table.insert(residues_m, m)
table.insert(residues_s, m) --For this we could exclude glycine and alanine because they have no alternative shake-position. But this would mean a time-consuming check with GetAminoAcid().
else
--Segments belonging to this category can either be segments in prediction-puzzles (non-locked and non-mutable), completely locked segments (bb and sc) in binder-puzzles, locked segments (only bb) in binder-puzzles which are not mutable.
--Fully locked segments shall be ignored here while all others shall be considered.
--For wiggling and shaking the bb-locked but sc-non-locked shall be considered. For mutate not since they must be non-mutable in this routine!
if (structure.IsLocked(m) == true) then
--In this case we know that the BACKBONE is locked. We don't know from this if the SC is locked as well! See doc for IsLocked()!
--Find out if the sidechain is locked as well. This cannot safely be done in all cases since for example glycine has no sidechain!
--The proposed way is to check if the sidechain has multiple rotamers. Glycine and alanine have only 1?!
local tmpRotamers = rotamer.GetCount(m)
if (tmpRotamers > 1) then
--In this case we know for sure that the sidechain is non-locked since we can shake it to different rotamers.
table.insert(residues_w, m)
table.insert(residues_s, m)
else
--In this case we either have a fully locked segment (bb and sc) or a locked glycine or alanine or a sidechain which has only one rotamer (this can possible happen in certain cases if the environment around the sidechain does prohibit rotamer-configs for this sidechain)
--We still want to wiggle segments that are not glycine but shake will not make sense on them since they only have one rotamer.
if (string.lower(structure.GetAminoAcid(m)) ~= "g") then
table.insert(residues_w, m)
end
end
else
--In this case we have prediction-segments which are non-locked and non-mutable. But we still want to shake and wiggle them but not mutate.
--But this may not make sense for glycine and alanine to shake them??? But it should also do no harm.
table.insert(residues_w, m)
table.insert(residues_s, m)
end
end
end
--Find actual segments to work on
findToWorkOn() --The variable 'toWorkOn' is subsequently required in the GUI
if (#toWorkOn > 0) then
print("Found " .. #toWorkOn .. " segments to work on.")
else
print("iwdn worm version " .. version .. " has finished with an error.")
end
--print("Selected segments2: " .. #selectedSegments)
--Determine which filters are available and have bonus.
local filterNames = filter.GetNames()
local bonusFilterNames = {} --Stores only those filters that actually can provide a bonus.
for m=1,#filterNames do
if (filter.HasBonus(filterNames[m])) then
table.insert(bonusFilterNames, filterNames[m])
end
end
--Now all filters that have bonus are stored with their names in bonusFilterNames.
--Only these filters need to be queried in the future for bonus-change.
--Show GUI
dia = dialog.CreateDialog("Iwdn worm " .. version)
dia.lbl01 = dialog.AddLabel("Worm-mode (just choose one of the next 3):")
dia.cbx01 = dialog.AddCheckbox("Do Wiggle", true)
dia.cbx02 = dialog.AddCheckbox("Do Mutate", false)
dia.cbx03 = dialog.AddCheckbox("Do Shake", false)
dia.lbl02a = dialog.AddLabel("")
dia.lbl02b = dialog.AddLabel("Worm-config:")
dia.sli01 = dialog.AddSlider("Iters", 6, 1, 100, 0)
dia.sli02 = dialog.AddSlider("Radius", 8, 1, 20, 0)
dia.sli05 = dialog.AddSlider("Step", 3, 1, 20, 0)
dia.sli03 = dialog.AddSlider("Total turns", 1, 1, 30, 0)
dia.sli04 = dialog.AddSlider("Qicksave startslot", qs_idx_start, 1, 99, 0)
dia.cbx03b = dialog.AddCheckbox("Wiggle sidechains only", false)
dia.cbx03c = dialog.AddCheckbox("More accurate distance-check", false)
dia.lbl03a = dialog.AddLabel("")
dia.lbl03b = dialog.AddLabel("Settings for Objectives:")
dia.cbx03a = dialog.AddCheckbox("Protect \"Score\"", true)
--For each bonus-objective provide a checkbox
local objCbxs = {} --Stores the names to dia[objCbxs[x]]
for m=1,#bonusFilterNames do
objCbxs[m] = "objCbx" .. m
local filterName = filterNames[m]
local tickIt = false
if (string.find(string.lower(filterName), "core") ~= nil and string.find(string.lower(filterName), "existence") ~= nil) then
tickIt = true
elseif (string.find(string.lower(filterName), "ideal") ~= nil and string.find(string.lower(filterName), "loop") ~= nil) then
tickIt = true
elseif (string.find(string.lower(filterName), "hydrogen") ~= nil and string.find(string.lower(filterName), "network") ~= nil) then
tickIt = true
elseif (string.find(string.lower(filterName), "disulfide") ~= nil and string.find(string.lower(filterName), "count") ~= nil) then
tickIt = true
elseif (string.find(string.lower(filterName), "buried") ~= nil and string.find(string.lower(filterName), "unsats") ~= nil) then
tickIt = true
end
dia[objCbxs[m]] = dialog.AddCheckbox("Protect " .. "\"" .. bonusFilterNames[m] .. "\"", tickIt)
end
dia.but01 = dialog.AddButton("OK", 1)
dia.but02 = dialog.AddButton("Cancel", 0)
tmpBack = dialog.Show(dia)
if (tmpBack == 0) then
print("Iwdn worm " .. version .. " has finished.")
return
end
--Determine the chosen options for Objectives and provide appropriate vector for the worm-function
local toMonitor = {}
if (dia.cbx03a.value == true) then --Global Score
toMonitor[1] = "Global Score" --If all other objectives did not degrade but Score did, then do not accept this solution in the worm.
end
for m=1,#bonusFilterNames do
if(dia[objCbxs[m]].value == true) then
table.insert(toMonitor, bonusFilterNames[m])
end
end
--Now the filter names that shall be monitored - along with flag for score - are stored in "toMonitor".
--We are ready to run the worm with it.
--Summary of GUI-options
print("")
print("Summary of chosen GUI-options:")
local modisActive = 0
if (dia.cbx01.value == true) then
modisActive = modisActive+1
end
if (dia.cbx02.value == true) then
modisActive = modisActive+2
end
if (dia.cbx03.value == true) then
modisActive = modisActive+4
end
if (modisActive == 1) then
--Only wiggle is selected
print("Wiggle-modus active.")
elseif (modisActive == 2) then
print("Mutate-modus active.")
elseif (modisActive == 4) then
print("Shake-modus active.")
else
--Multiple modi were selected. This is an error/warning
print("Warning: only tick 1 of the 3 modi!")
if (modisActive == 3 or modisActive == 5 or modisActive == 7) then
print("Wiggle-modus active.")
else
print("Mutate-modus active.")
end
end
print("Worm iterations: " .. dia.sli01.value)
print("Worm radius: " .. dia.sli02.value)
print("Worm step: " .. dia.sli05.value)
print("Worm total turns: " .. dia.sli03.value)
print("Quicksave startslot: " .. dia.sli04.value)
qs_idx_start = dia.sli04.value
qs_idx = qs_idx_start
print("Wiggle sidechains only (backbone is not wiggled): " .. tostring(dia.cbx03b.value))
print("More accurate distance-check: " .. tostring(dia.cbx03c.value))
for m=1,#toMonitor do
print("Filter to monitor: " .. toMonitor[m])
end
print("Storing initial solution in qs" .. qs_idx_start .. " for later reference.")
save.Quicksave(qs_idx_start)
qs_idx = qs_idx+1
print("Worm uses qs-slot qs" .. qs_idx .. ".")
--Run the worm
print("")
print("Running worm with chosen options (please wait)...")
local w_stIdx = 1
for m=1,dia.sli03.value do --Total turns loop
if (dia.cbx01.value == true) then --Wiggle
print("Wiggle worm protected (" .. dia.sli01.value .. " iters) turn " .. m .. "/" .. dia.sli03.value .. "...")
if (dia.cbx03b.value == false) then
worm(toWorkOn, 1, w_stIdx, dia.sli05.value, dia.sli02.value, dia.sli01.value, toMonitor)
else
worm(toWorkOn, 10, w_stIdx, dia.sli05.value, dia.sli02.value, dia.sli01.value, toMonitor)
end
elseif (dia.cbx02.value == true) then --Mutate
print("Mutate worm protected (" .. dia.sli01.value .. " iters) turn " .. m .. "/" .. dia.sli03.value .. "...")
worm(toWorkOn, 2, w_stIdx, dia.sli05.value, dia.sli02.value, dia.sli01.value, toMonitor)
elseif (dia.cbx03.value == true) then --Shake
print("Shake worm protected (" .. dia.sli01.value .. " iters) turn " .. m .. "/" .. dia.sli03.value .. "...")
worm(toWorkOn, 3, w_stIdx, dia.sli05.value, dia.sli02.value, dia.sli01.value, toMonitor)
end
w_stIdx = w_stIdx+1 --Wrap the start-index
if (w_stIdx >= dia.sli05.value) then
w_stIdx = 1
end
end
print("Worm is done.")
print("Loading initial and recent best solution to Undo-stack...")
undo.SetUndo(true)
for n=1,qs_idx-qs_idx_start+1 do
save.Quickload(qs_idx_start+n-1) --Puts found solutions into restore. Use ctrl+z/y to flip through.
wait(guardInt01) --This seems to be necessary in order for things to settle.
print("Solution in qs" .. qs_idx_start+n-1 .. " (" .. n .. "/" .. qs_idx-qs_idx_start+1 .. ") loaded.")
end
--Restore initial selection. May come in handy for further actions.
selection.DeselectAll()
for o=1,#selectedSegments do
selection.Select(selectedSegments[o])
end
print("Iwdn worm version " .. version .. " has finished.")
end --End: main
xpcall(main, cleanup)
--main()