Code
-- Base on Local Quake 1.0 by spvinvent with helep by thom
-- Bands radiate from a structure (helix, loop or sheet) to spatially close neighbours
-- one round with the strcuture frozen one round not
-- Foregoes Fuses and generally tries to optimize for speed
-- "lua v1 in v2" library by rav3n_pl
--just add it in front of your v1 code and use v2 and v1 code in v2 scripts :)
-- print(arg1[,...,argN]) no change :)
function are_conditions_met()
return current.AreConditionsMet()
end
function band_add_segment_segment(sgi1, sgi2)
band.AddBetweenSegments(sgi1, sgi2)
end
function band_delete(bndIdx)
if bndIdx~= nil then
band.Delete(bndIdx)
else
band.DeleteAll()
end
end
function band_disable(bndIdx)
if bndIdx~=nil then
band.Disable(bndIdx)
else
band.DisableAll()
end
end
function band_enable(bndIdx)
if bndIdx~=nil then
band.Enable(bndIdx)
else
band.EnableAll()
end
end
function band_set_length(bndIdx, len)
band.SetGoalLength(bndIdx, len)
end
function band_set_strength(bndIdx, str)
band.SetStrength(bndIdx, str)
end
function get_band_count()
return band.GetCount()
end
function deselect_all()
selection.DeselectAll()
end
function deselect_index(sgn)
selection.Deselect(sgn)
end
function select_all()
selection.SelectAll()
end
function select_index(sgn)
selection.Select(sgn)
end
function select_index_range(sg1,sg2)
selection.SelectRange(sg1,sg2)
end
function do_freeze(bbone,schain)
freeze.FreezeSelected(bbone,schain)
end
function do_unfreeze_all()
freeze.UnfreezeAll()
end
function do_global_wiggle_all(iters)
structure.WiggleAll(iters,true,true)
end
function do_global_wiggle_backbone(iters)
structure.WiggleAll(iters, true,false)
end
function do_global_wiggle_sidechains(iters)
structure.WiggleAll(iters,false,true)
end
function do_local_rebuild(iters)
structure.RebuildSelected(iters)
end
function do_local_wiggle(iters)
structure.LocalWiggleSelected(iters,true,true)
end
function do_mutate(iters)
structure.MutateSidechainsSelected(iters)
end
function do_shake(iters)
structure.ShakeSidechainsSelected(iters)
end
function do_sidechain_snap(sgn, snap)
rotamer.SetRotamer(sgn, snap)
end
function get_sidechain_snap_count(sgn)
return rotamer.GetCount(sgn)
end
function load_structure()
save.LoadSecondaryStructure()
end
function save_structure()
save.SaveSecondaryStructure()
end
function quickload(slot)
save.Quickload(slot)
end
function quicksave(slot)
save.Quicksave(slot)
end
function get_exploration_score()
return current.GetExplorationMultiplier()
end
function get_ranked_score()
return current.GetScore()
end
function get_score()
return current.GetEnergyScore()
end
function get_segment_distance(sg1,sg2)
return structure.GetDistance(sg1,sg2)
end
function get_aa(sn)
return structure.GetAminoAcid(sn)
end
function get_segment_count()
return structure.GetCount()
end
function get_ss(sn)
return structure.GetSecondaryStructure(sn)
end
function is_hydrophobic(sn)
return structure.IsHydrophobic(sn)
end
function replace_aa(aa)
for i=1,structure.GetCount() do
if selection.IsSelected(i) then
structure.SetAminoAcid(i, aa)
end
end
end
function replace_ss(ss)
for i=1,structure.GetCount() do
if selection.IsSelected(i) then
structure.SetSecondaryStructure(i,ss)
end
end
end
function get_segment_score(sg)
return current.GetSegmentEnergyScore(sg)
end
function get_segment_score_part(score_part,sg)
return current.GetSegmentEnergySubscore(sg,score_part)
end
function reset_puzzle()
puzzle.StartOver()
end
function restore_abs_best()
absolutebest.Restore()
end
function restore_credit_best()
creditbest.Restore()
end
function reset_recent_best()
recentbest.Save()
end
function restore_recent_best()
recentbest.Restore()
end
function set_behavior_clash_importance(ci)
behavior.SetClashImportance(ci)
end
-- end of library
-- 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
-- function to overload a funtion
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)
use_ranked_score = 1 -- if 1, use the score being used for ranking (best normally). If 0, use the
-- old style scoring system
max_bands_per_residue = 6 -- Increase this for more bands at each segment
max_segment_length = 999 -- Don't create bands greater than this length. Set to 999 or something if you're
-- not worried about the potential for long bands creeping in,
min_exclusion = 999 -- These can be set to establish an "exclusion zone", within which no bands will
-- be created. This might be used to prevent banding to a wayward loop.
-- If min_exclusion > max_exclusion no exclusion will be set.
max_exclusion = 1
target_wiggle_reduction = 150 -- On wiggle, try and reduce the score by this amount after creating the bands.
-- If it's wildly out, reduce/increase band strength
-- for the next try (don't try to adjust for the current scramble).
max_segment_length = 10 -- Don't create bands greater than this length. Set to 999 or something if you're
-- not worried about the potential for long bands creeping in,
n_quakes = 50000 -- Go crazy
band_to_minima = true -- Band to minima or maxima in the distance. See bands_from_segment for meaning.
do_second_shake = false -- After creating the bands and pulling with them do a Shake/Wiggle/Shake
-- if do_second_shake is false (best I think) don't do the last shake.
-- Speeds things up quite a bit
scale_band_strength = true -- If true, band strength will be scaled inversely with length
monit = false
---------- A few more globals
candidate_ids_start = {}
candidate_ids_end = {}
candidate_distances = {}
n_candidates = 0
best_score = 0
function GetScore ()
if ( use_ranked_score == 1 ) then
return get_ranked_score ( true )
else
return get_score ()
end
end
function Coprime ( n )
if ( n < 31 ) then
return 1
elseif ( n >= 31*37 ) then
return 1
elseif ( n%31 ~= 0 ) then
return 31
else
return 37
end
end
function pseudorandom ( score )
-- Returns a pseudorandom number >= 0 and < 1.
-- Based on the fractional part of the score
if ( score >= 0 ) then
return score % 1
else
return (-score ) % 1
end
end
function delete_all_bands ()
band_delete ()
--band_count = get_band_count ()
-- for i = band_count , 1 , -1 do
-- band_delete ( i )
--end
end
function ShellSort ( ids , distances , n )
-- Adapted from Numerical Recipes in C
local inc = 1
repeat
inc = inc * 3 + 1
until inc > n
repeat
inc = inc / 3
inc = inc - inc % 1
for i = inc + 1 , n do
v = distances [ i ]
w = ids [ i ]
j = i
flag = false
while ( flag == false and distances [ j - inc ] > v ) do
distances [ j ] = distances [ j - inc ]
ids [ j ] = ids [ j - inc ]
j = j - inc
if ( j <= inc ) then
flag = true
end
end
distances [ j ] = v
ids [ j ] = w
end
until inc <= 1
end
function test_for_minimum ( i , seg , seg_distances )
if ( ( i ~= seg ) and
( seg_distances [ i ] <= seg_distances [ i - 1 ] ) and
( seg_distances [ i ] <= seg_distances [ i - 2 ] ) and
( seg_distances [ i ] <= seg_distances [ i + 1 ] ) and
( seg_distances [ i ] <= seg_distances [ i + 2 ] ) ) then
return true
else
return false
end
end
function test_for_maximum ( i , seg , seg_distances )
if ( ( seg_distances [ i ] >= seg_distances [ i - 1 ] ) and
( seg_distances [ i ] >= seg_distances [ i - 2 ] ) and
( seg_distances [ i ] >= seg_distances [ i + 1 ] ) and
( seg_distances [ i ] >= seg_distances [ i + 2 ] ) ) then
return true
else
return false
end
end
function get_band_segment_data ( seg )
seg_distances = {}
-- Get the function that looks like this:
-- x = residue number
-- y = distance of that segment from seg
-- ignore the degenerate minimum
n_segments = 0
segment_ids = {}
segment_distances = {}
for i = 1, n_residues do
segment_distances [ i ] = get_segment_distance ( seg , i )
end
for i = 3 , n_residues - 3 do
k = segment_distances [ i ]
if ( ( band_to_minima == true ) and
( test_for_minimum ( i , seg , segment_distances ) == true ) ) then
n_segments = n_segments + 1
segment_ids [ n_segments ] = i
segment_distances [ n_segments ] = segment_distances [ i ]
elseif ( ( band_to_minima == false ) and
( test_for_maximum ( i , seg , segment_distances ) == true ) ) then
n_segments = n_segments + 1
segment_ids [ n_segments ] = i
segment_distances [ n_segments ] = segment_distances [ i ]
end
end
ShellSort ( segment_ids , segment_distances , n_segments )
if ( n_segments >= max_bands_per_residue ) then
n_c = max_bands_per_residue
else
n_c = n_segments
end
for i = 1 , n_c do
if ( ( segment_distances [ i ] <= max_segment_length ) and
( seg < min_exclusion or seg > max_exclusion ) and
( segment_ids [ i ] < min_exclusion or segment_ids [ i ] > max_exclusion ) ) then
n_candidates = n_candidates + 1
candidate_ids_start [ n_candidates ] = seg
candidate_ids_end [ n_candidates ] = segment_ids [ i ]
candidate_distances [ n_candidates ] = segment_distances [ i ]
end
end
end
function create_bands ()
if ( scale_band_strength == true ) then
total_band_length = 0
for i = 1 , n_candidates do
total_band_length = total_band_length + candidate_distances [ i ]
end
average_band_length = total_band_length / n_candidates
end
for i = 1 , n_candidates do
band_add_segment_segment ( candidate_ids_start [ i ] , candidate_ids_end [ i ] )
if ( scale_band_strength == true ) then
band_set_strength ( i , band_strength * average_band_length / candidate_distances [ i ] )
else
band_set_strength ( i , band_strength )
end
end
end
function Quake ( seg )
select_all ()
init_score = GetScore ()
do_global_wiggle_backbone ( 1 )
score_after_wiggle = GetScore ()
score_reduction_after_wiggle = init_score - score_after_wiggle
if ( monit ) then
print ( "seg ", seg , " score loss " , score_reduction_after_wiggle )
end
-- Adaptively change band strength based on how much the previous wiggle overshot/undershot the target. Could
-- be made much more elaborate.
if ( score_reduction_after_wiggle > target_wiggle_reduction * 1.2 ) then
band_strength = band_strength * ( 0.85 + 0.1 *pseudorandom ( score_after_wiggle ) )
elseif ( score_reduction_after_wiggle < target_wiggle_reduction * 0.8 ) then
band_strength = band_strength + ( 1 - band_strength ) * ( 0.05 + 0.1 *pseudorandom ( score_after_wiggle ) )
end
if ( monit ) then
print ( " band_strength ", band_strength )
end
delete_all_bands ()
do_shake ( 1 )
gas = GetScore () - score_after_wiggle
if ( monit ) then
print ( " gas ", gas )
end
do_global_wiggle_all ( 25 )
if ( do_second_shake == true ) then
score_a2w = GetScore ()
do_shake ( 1 )
score_a2s = GetScore ()
if ( score_a2s - score_a2w > 1 ) then
do_global_wiggle_all ( 4 )
end
end
end
n_residues = get_segment_count ()
band_strength = 0.5
set_behavior_clash_importance ( 1 )
reset_recent_best ()
best_score = GetScore ()
print ( " init score ", best_score )
idx = 1
gel=1
for i = 1, n_quakes do
restore_recent_best ()
-- Need to record any improvements from the previous iteration
score = GetScore()
if ( score > best_score ) then
best_score = score
print ( "Improvement to " , score )
end
do_unfreeze_all()
delete_all_bands ()
n_candidates = 0
debut=0
if (idx==1) then debut=1
else if ( get_ss(idx)~=get_ss(idx-1) ) then debut =idx end end
if ( debut>0) then
while (idx < n_residues+1 and get_ss(idx)==get_ss(debut)) do
get_band_segment_data ( idx )
if gel ==1 then
deselect_all()
select_index(idx)
do_freeze(true,false)
end -- of gel
idx=idx+1
end
end
select_all()
if ( n_candidates > 0 ) then
create_bands ()
print ( "i: " , i , " seg ", debut," " ,idx-1 )
Quake ( idx )
end
if idx>n_residues then idx=1 if gel==1 then gel=0 else gel=1 end end
end
-- to see every line