Code
-- Globals
-- Version 2.0 Filter: added filters management by BitsPawn
-- Version 3.0 Optional filter and cuts (Bruno Kestemont)
recipename="Helix twister cut 3.0"
original_secondary_structure = {}
helix_starts = {}
helix_ends = {}
use_helix = {}
n_helices = 0
n_residues = 0
loop_spacer = 3
clk_twist = true
kBandLength = 5
band_strengths = {}
target_loss_on_twist = 150
kSlackTolerance = 1.8
mutate_not_shake = false
test_for_credit = false
keep_trying = true
DISABLEFILTERS= false
-- Quicksave slots
kOriginalStructureOrNewBest = 1 -- the starting structure, or any subsequent improvement, will be stored in this quicksave slot
-- function to copy class/table
function CopyTable(orig)
local copy = {}
for orig_key, orig_value in pairs(orig) do
copy[orig_key] = orig_value
end
return copy
end
-- functions for filters
function FiltersOn()
if behavior.GetSlowFiltersDisabled() then
behavior.SetSlowFiltersDisabled(false)
end
end
function FiltersOff()
if behavior.GetSlowFiltersDisabled()==false then
behavior.SetSlowFiltersDisabled(true)
end
end
--Cut functions
allowcut=true
CUTWORST= false -- not used here
function cut(seg) -- 18/2/2016, adapted for HR (n_residues)
if (allowcut or CUTWORST) and not freeze.IsFrozen(seg) then
if seg==n_residues then seg=n_residues-1 end
if seg<1 then seg=1 end
structure.InsertCut(seg)
end -- else nothing !
end
function uncut() -- 18/2/2016, deletes all cuts ! (unless seg is frozen)
if (allowcut or CUTWORST) then
for i=1, n_residues do
structure.DeleteCut(i)
end
end -- else nothing !
end
--end cut functions
-- function to overload a function
function mutFunction(func)
local currentfunc = func
local function mutate(func, newfunc)
local lastfunc = currentfunc
currentfunc = function(...) return newfunc(lastfunc, ...) end
end
local wrapper = function(...) return currentfunc(...) end
return wrapper, mutate
end
-- function to overload a class
-- to do: set the name of function
classes_copied = 0
myclcp = {}
function MutClass(cl, filters)
classes_copied = classes_copied+1
myclcp[classes_copied] = CopyTable(cl)
local mycl =myclcp[classes_copied]
for orig_key, orig_value in pairs(cl) do
myfunc, mutate = mutFunction(mycl[orig_key])
if filters==true then
mutate(myfunc, function(...)
FiltersOn()
if table.getn(arg)>1 then
-- first arg is self (function pointer), we pack from second argument
local arguments = {}
for i=2,table.getn(arg) do
arguments[i-1]=arg[i]
end
return mycl[orig_key](unpack(arguments))
else
--print("No arguments")
return mycl[orig_key]()
end
end)
cl[orig_key] = myfunc
else
mutate(myfunc, function(...)
FiltersOff()
if table.getn(arg)>1 then
local arguments = {}
for i=2, table.getn(arg) do
arguments[i-1]=arg[i]
end
return mycl[orig_key](unpack(arguments))
else
return mycl[orig_key]()
end
end)
cl[orig_key] = myfunc
end
end
end
-- how to use:
--MutClass(structure, false)
--MutClass(band, false)
--MutClass(current, true)
function r3 ( i )
-- printing convenience
return i - i % 0.001
end
function GetScore ()
score = current.GetEnergyScore ()
return score
end
function IsPuzzleMutable ()
for i = 1 , n_residues do
if ( structure.IsMutable ( i ) ) then
return true
end
end
return false
end
function get_helices ()
within_helix = false
for i = 1, n_residues do
if ( structure.GetSecondaryStructure ( i ) == "H" ) then
if ( within_helix == false ) then
-- start of a new helix
within_helix = true
n_helices = n_helices + 1
helix_starts [ n_helices ] = i
end
elseif ( within_helix == true ) then
-- end of a helix
within_helix = false
helix_ends [ n_helices ] = i -1
end
end -- for i
if ( within_helix == true ) then
helix_ends [ n_helices ] = n_residues
end
end
function GetClosestHelix ( i )
-- From the list of selected helices, choose the one closest to helix i
-- Very crude, looks at the distance between helix centres
mid_i = ( helix_starts [ i ] + helix_ends [ i ] ) / 2
closest_j = -1
closest_distance = 99999
for j = 1 , n_helices do
if ( j ~=i ) then
mid_j = ( helix_starts [ j ] + helix_ends [ j ] ) / 2
d = structure.GetDistance ( mid_i , mid_j )
if ( d < closest_distance ) then
closest_distance = d
closest_j = j
end
end
end
return closest_j
end
function IsWithinASelectedHelix ( residue )
for i = 1 , n_helices do
if ( use_helix [ i ] == true ) then
if ( ( residue >= helix_starts [ i ] ) and ( residue >= helix_ends [ i ] ) ) then
return true
end
end
end
return false
end
function ConvertToLoop ( first , last )
for k = first , last do
if ( original_secondary_structure [ k ] ~= "L" ) then
if ( IsWithinASelectedHelix ( k ) == false ) then
selection.SelectRange ( k , k )
structure.SetSecondaryStructureSelected ( "L" )
end
end
end
end
function ConstructSpacers ( start_helix , end_helix )
local i
local ok
start_idx = start_helix - loop_spacer
if ( start_idx < 1 ) then
start_idx = 1
end
i = start_helix - 1
ok = true
while ( i >= start_idx and ok == true ) do
if ( structure.IsLocked ( i ) == true ) then
ok = false
start_idx = i + 1
else
i = i -1
end
end
if ( start_helix > 1 ) then
ConvertToLoop ( start_idx , start_helix - 1 )
end
end_idx = end_helix + loop_spacer
if ( end_idx > n_residues ) then
end_idx = n_residues
end
i = end_helix + 1
ok = true
while ( i <= end_idx and ok == true ) do
if ( structure.IsLocked ( i ) == true ) then
ok = false
end_idx = i - 1
else
i = i +1
end
end
if ( end_helix < n_residues ) then
ConvertToLoop ( end_helix + 1 , end_idx )
end
return start_idx , end_idx
end
function Wiggle ( step , threshold )
kMaxIterations = 10 -- Sanity check
n_it = 0
new_score = GetScore ()
repeat
prev_score = new_score
structure.WiggleSelected ( step )
new_score = GetScore ()
n_it = n_it + 1
until ( ( math.abs ( new_score - prev_score ) < threshold ) or ( n_it > kMaxIterations ) )
end
function AdjustBandStrength ( factor )
n_bands = band.GetCount ()
for i = 1 , n_bands do
strength = band.GetStrength ( i )
band.SetStrength ( i , strength * factor )
strength = band.GetStrength ( i )
if ( i == 1 ) then
undo.SetUndo ( false )
end
end
undo.SetUndo ( true )
end
function WiggleBandedProtein ( helix_index )
-- We are given a protein with a bunch of bands between segments: the intent is to wiggle until
-- until the reduction from the start position is reasonably close to the desired value.
-- Wiggle once to get an approximate idea of whether bands need to be strengthened, weakened,
-- or whether they're just right.
init_score = GetScore ()
structure.WiggleSelected ( 2 )
score = GetScore ()
loss_after_banded_wiggle = init_score - score
--print ( "WBP : init reduction " .. r3 ( loss_after_banded_wiggle ) .. " strength " .. r3 ( band_strengths [ helix_index ] ) )
if ( loss_after_banded_wiggle > target_loss_on_twist ) then
factor = 1.0 - 0.2 * math.random ()
else
factor = 1.0 + 0.2 * math.random ()
end
band_strengths [ helix_index ] = band_strengths [ helix_index ] * factor
band_strengths [ helix_index ] = math.max ( band_strengths [ helix_index ] , 0.2 )
band_strengths [ helix_index ] = math.min ( band_strengths [ helix_index ] , 2 )
end
function ConstructBands ( helix_index , start_helix , end_helix )
for i = math.max ( start_helix , 2 ) , math.min ( end_helix , n_residues - 1 ) do
if ( clk_twist ) then
band.Add ( i , i + 1 , i - 1 , kBandLength , math.pi/2 , math.rad ( 315 ) )
else
band.Add ( i , i - 1 , i + 1 , kBandLength , math.pi/2 , math.rad ( 315 ) )
end
if ( i == 1 ) then
undo.SetUndo ( false )
end
end
if ( math.abs ( band_strengths [ helix_index ] - 1.0 ) > 1e-3 ) then
n_bands = band.GetCount ()
for i = 1 , n_bands do
band.SetStrength ( i , band_strengths [ helix_index ] )
end
end
undo.SetUndo ( true )
end
function TwistHelix ( helix_index , start_helix , end_helix )
if ( clk_twist == true ) then
print ( "H " .. helix_index .. " (c)" )
else
print ( "H " .. helix_index .. " (a)" )
end
improvement_made = false
save.Quickload ( kOriginalStructureOrNewBest )
best_score = GetScore ()
selection.DeselectAll ()
start_idx , end_idx = ConstructSpacers ( start_helix , end_helix )
ConstructBands ( helix_index , start_helix , end_helix )
--cut before and after
cut(start_idx)
cut(end_idx)
selection.SelectRange ( start_idx , end_idx )
WiggleBandedProtein ( helix_index )
uncut()
score_after_wiggle= GetScore ()
if ( math.abs ( score_after_wiggle - best_score ) < 1e-5 ) then
-- If the scores are still equal we probably have a locked segment
print ( "Possible issue with helix " .. helix_index .. " (" .. start_helix .. "-" .. end_helix .. ")" )
-- Trying to unfreeze: if successful rebuild will work next time this function is called
for i = start_idx , end_idx do
freeze.Unfreeze ( i , true , true )
end
return false
end
band.DeleteAll ()
selection.SelectAll ()
if ( mutate_not_shake == true ) then
structure.MutateSidechainsSelected ( 1 )
else
structure.ShakeSidechainsSelected ( 1 )
end
structure.WiggleSelected ( 2 , false , true ) -- sidechains
Wiggle ( 4 , 0.1 )
selection.SelectAll ()
curr_score = GetScore ()
if ( ( curr_score > best_score ) and
( test_for_credit == false or creditbest.AreConditionsMet () == true ) ) then
improvement_made = true
best_score = curr_score
save.LoadSecondaryStructure ()
save.Quicksave ( kOriginalStructureOrNewBest )
print ( "Improvement to " .. r3 ( best_score ) )
end
save.LoadSecondaryStructure ()
return true
end
function get_parameters ()
local dlog = dialog.CreateDialog ( recipename )
dlog.loop_spacer = dialog.AddSlider ( "Loop spacer" , loop_spacer , 0 , 6 , 0 )
dlog.tlot = dialog.AddSlider ( "Target loss" , target_loss_on_twist , 0 , 1000 , 0 )
for i = 1 , n_helices do
dlog [ "helix" .. i ] = dialog.AddCheckbox ( helix_starts [ i ] .. "-" .. helix_ends [ i ] , true )
end
design_puzzle = IsPuzzleMutable ()
if ( design_puzzle == true ) then
dlog.mutate= dialog.AddCheckbox ( "Mutate not shake" , false )
dlog.DISABLEFILTERS=dialog.AddCheckbox ( "Disable filters unless for score" , false )
end
dlog.clk = dialog.AddCheckbox ( "Clockwise " , clk_twist )
dlog.credit_check = dialog.AddCheckbox ( "Check conditions met " , test_for_credit )
dlog.allowcut = dialog.AddCheckbox ( "Cut helix out " , allowcut )
dlog.keep_trying = dialog.AddCheckbox ( "Keep trying " , keep_trying )
dlog.ok = dialog.AddButton ( "OK" , 1 )
dlog.cancel = dialog.AddButton ( "Cancel" , 0 )
if ( dialog.Show ( dlog ) > 0 ) then
loop_spacer = dlog.loop_spacer.value
target_loss_on_twist = dlog.tlot.value
for i = 1 , n_helices do
use_helix [ i ] = dlog [ "helix" .. i ].value
end
clk_twist = dlog.clk.value
if ( design_puzzle == true ) then
mutate_not_shake = dlog.mutate.value
if dlog.DISABLEFILTERS.value then
MutClass(structure, false)
MutClass(band, false)
MutClass(current, true)
end
end
test_for_credit = dlog.credit_check.value
allowcut=dlog.allowcut.value
keep_trying = dlog.keep_trying.value
return true
else
return false
end
end
function main ()
print ( recipename )
band.DeleteAll ()
n_residues = structure.GetCount ()
for i = 1, n_residues do
original_secondary_structure [ i ] = structure.GetSecondaryStructure ( i )
end
save.SaveSecondaryStructure ()
save.Quicksave ( kOriginalStructureOrNewBest )
behavior.SetClashImportance ( 1 )
best_score = GetScore ()
print ( "Start score : " .. r3 ( best_score ) )
get_helices ()
print ( n_helices .. " helices" )
qs10_best_score = 0
if ( get_parameters () == false ) then
print ( "error" )
error ()
end
print ( "Band length : " .. r3 ( kBandLength ) )
-- Delete unused helices
for i = n_helices , 1 , -1 do
if ( use_helix [ i ] == false ) then
table.remove ( helix_starts , i )
table.remove ( helix_ends , i )
end
end
n_helices = #helix_starts
if ( n_helices == 0 ) then
print ( "No helix selected" )
error ()
end
for i = 1 , n_helices do
print ( "Helix " .. i .. " : (" .. helix_starts [ i ] .. "-" .. helix_ends [ i ] .. ")" )
end
for i = 1 , n_helices do
band_strengths [ i ] = 1.0
end
for i = 1 , n_helices do
TwistHelix ( i , helix_starts [ i ] , helix_ends [ i ] )
end
clk_twist = not clk_twist
for i = 1 , n_helices do
TwistHelix ( i , helix_starts [ i ] , helix_ends [ i ] )
end
if ( keep_trying == true ) then
while ( true ) do
if loop_spacer< 6 then loop_spacer=loop_spacer+1 else loop_spacer=1 end
clk_twist = true
for i = 1 , n_helices do
TwistHelix ( i , helix_starts [ i ] , helix_ends [ i ] )
end
clk_twist = false
for i = 1 , n_helices do
TwistHelix ( i , helix_starts [ i ] , helix_ends [ i ] )
end
end
end
band.DeleteAll ()
save.Quickload ( kOriginalStructureOrNewBest )
end
function cleanup ()
print ( "Cleaning up" )
behavior.SetClashImportance ( 1.0 )
save.Quickload ( kOriginalStructureOrNewBest )
band.EnableAll ()
undo.SetUndo ( true )
end
--main ()
xpcall ( main , cleanup )
