Code
rebuild_length = 6 -- Length of sections to rebuild.
number_of_rebuilds = 8 -- It'll look for the best score from this number of rebuilds.
start_idx = 0 -- if this is a residue id, use it as the start point. If not
-- use a semi-random number computed from the score
loops_only = true -- If true, will still allow a single non-loop residue in the rebuild segment
-- If false, non-loop residues will be temporarily converted to loops
min_residue = 1 -- Change these if you want to rebuild part of a protein only. For example, to
-- rebuild 53-64, set min_residue to 53 and max_residue to 64. To restore, set to
-- 1 and 9999 (or any large number) respectively
max_residue = 999
kWorstAllowableAngle = 75 -- ignore rebuilds with backbone angles less than this (degrees)
kMaxSphereDistance = 14 -- Residues further than this away won't be shaken
-- Use of quicksave slots
kOriginalStructureOrNewBest = 1
kBestPostRebuild = 2
kTemp = 3
-- Globals
original_secondary_structure = {}
fixed_residues = {}
--wiggle_accuracy = 0
n_residues = 0
best_score = 0
cos_worst_allowable_angle = 0
-- NudgeWiggle : adaptive thresholds. Change the thresholds based on current protein.
kInitNudgeWiggleThreshold = 25
nudge_wiggle_gains = {}
kLowestNudgeWiggleThreshold = 8
-- Performance monitoring. Result intended for import into Excel as comma-delimited text
monit = false
monit_init_worst_angle = 0
monit_worst_angle = 0
monit_score_after_rebuild = 0
monit_score_after_post_rebuild_shake = 0
monit_gain_from_possible_initial_nudge_wiggle = 0
monit_gain_from_wiggle = 0
monit_gain_from_shake = 0
monit_gain_from_nudge_wiggle = 0
monit_total_gain = 0
-- Begin
function r3 ( x )
-- Round to 3 decimal places
t = 10 ^ 3
return math.floor ( x*t + 0.5 ) / t
end
function GetScore ()
score = current.GetEnergyScore ()
return score
end
function Coprime ( n )
-- returns a number that has no common factors with n
-- an ugly hack
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 IsSegmentALoop ( first , last )
-- Returns number of non-loop residues (helices or sheets)
n_non_loops = 0
for i = first , last do
if ( structure.GetSecondaryStructure ( i ) ~= "L" ) then
n_non_loops = n_non_loops + 1
end
end
return n_non_loops
end
function ConvertToLoop ( first , last )
for k = first , last do
if ( original_secondary_structure [ k ] ~= "L" ) then
selection.DeselectAll ( )
selection.Select ( k )
structure.SetSecondaryStructureSelected ( "L" )
end
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 GetListOfFixedResidues ()
for i = 1 , n_residues do
fixed_residues [ i ] = false
end
for i = 1 , n_residues do
seg_backbone_frozen, seg_sidechain_frozen = freeze.IsFrozen ( i )
if ( seg_backbone_frozen == true ) then
fixed_residues [ i ] = true
end
s = structure.IsLocked ( i )
if ( s == true ) then
fixed_residues [ i ] = true
end
end
end
function DoesSectionContainFixedResidues ( start_idx , end_idx )
for i = start_idx , end_idx do
if ( fixed_residues [ i ] == true ) then
return true
end
end
return false
end
function GetNudgeWiggleThreshold ()
n = #nudge_wiggle_gains
k = math.max ( n - 5 , 1 )
local x = 0
for i = k , n do
x = x + nudge_wiggle_gains [ i ]
end
av = x / ( n - k + 1 )
return math.max ( kLowestNudgeWiggleThreshold , 1.4*av )
end
function MiniFuse ()
selection.SelectAll ()
score_in_kTemp = GetScore()
save.Quicksave ( kTemp)
selection.SelectAll ()
behavior.SetClashImportance ( 0.07 )
structure.ShakeSidechainsSelected ( 1 )
behavior.SetClashImportance ( 1 )
structure.WiggleSelected ( 8 )
score = GetScore()
if ( score < score_in_kTemp ) then
save.Quickload ( kTemp )
end
end
function NudgeWiggle ( score )
selection.SelectAll ()
score_in_kTemp = GetScore()
save.Quicksave ( kTemp)
-- behavior.SetWiggleAccuracy ( 2 )
behavior.SetClashImportance ( 0.4 + 0.2 * pseudorandom ( score ) )
structure.WiggleSelected ( 1 )
-- behavior.SetWiggleAccuracy ( wiggle_accuracy )
behavior.SetClashImportance ( 1 )
structure.WiggleSelected ( 6 )
score = GetScore()
gain_from_nw = score - score_in_kTemp
if ( monit ) then
monit_gain_from_nudge_wiggle = gain_from_nw
end
if ( score < score_in_kTemp ) then
save.Quickload ( kTemp )
end
end
function WiggleShakeNudgeWiggle ()
-- returns the score after doing a WSNW: short-circuiting if things don't look promising
score_after_rebuild_shake = GetScore ()
structure.WiggleSelected ( 10 )
score_after_wiggle = GetScore ()
if ( monit ) then
monit_gain_from_wiggle = score_after_wiggle - score_after_rebuild_shake
end
if ( score_after_wiggle - score_after_rebuild_shake < 10 ) then
-- The wiggle got stuck and didn't achieve anything
NudgeWiggle ( score_after_wiggle )
score_after_nudge_wiggle = GetScore ()
if ( monit ) then
monit_gain_from_possible_initial_nudge_wiggle = score_after_nudge_wiggle - score_after_rebuild_shake
end
score_after_wiggle = score_after_nudge_wiggle
end
nwt = GetNudgeWiggleThreshold ()
if ( best_score - score_after_wiggle > 1.5*nwt ) then
return score_after_wiggle
end
structure.ShakeSidechainsSelected ( 1 )
score_after_second_shake = GetScore ()
if ( monit ) then
monit_gain_from_shake = score_after_second_shake - score_after_wiggle
end
if ( best_score - score_after_second_shake < nwt ) then
-- nwt = GetNudgeWiggleThreshold ()
start_score = GetScore()
NudgeWiggle ( score_after_second_shake )
MiniFuse ()
-- print ( "nwt " .. r3 ( nwt ) )
score = GetScore ()
if ( score - start_score > 0 ) then
table.insert ( nudge_wiggle_gains , score - start_score )
end
else
score = score_after_second_shake
end
return score
end
function post_rebuild ()
save.Quickload ( kBestPostRebuild )
if ( monit ) then
monit_score_after_post_rebuild_shake = GetScore ()
end
curr_score = WiggleShakeNudgeWiggle ()
if ( monit ) then
monit_total_gain = curr_score - best_score
end
if ( curr_score > best_score ) then
best_score = curr_score
save.LoadSecondaryStructure ()
save.Quicksave ( kOriginalStructureOrNewBest )
print ( "Improvement to ".. r3 ( best_score ) )
end
end
function get_cos_backbone_angle ( i )
if ( i == 1 or i == n_residues ) then
return -1
else
-- Cosine rule
a = structure.GetDistance ( i - 1 , i )
b = structure.GetDistance ( i , i + 1 )
c = structure.GetDistance ( i - 1 , i + 1 )
cos_angle = ( a*a + b*b - c*c ) / ( 2*a*b )
return cos_angle
end
end
function get_worst_angle ( start_idx , end_idx )
-- For monitoring only
cos_worst_angle = -1
for i = start_idx , end_idx do
cos_worst_angle = math.max ( cos_worst_angle , get_cos_backbone_angle ( i ) )
end
return math.deg ( math.acos ( cos_worst_angle ) )
end
function is_backbone_reasonably_angled ( start_idx , end_idx )
for i = start_idx , end_idx do
if ( get_cos_backbone_angle ( i ) > cos_worst_allowable_angle ) then
return false
end
end
return true
end
function ConstructListOfResiduesWithinRange ( SphereSegList , start_idx , end_idx )
mid_idx = math.floor ( ( end_idx + start_idx ) / 2 )
-- Look at start_idx, end_idx, and mid_idx only to save time at the expense of a small amount of accuracy
for i = 1 , n_residues do
if ( ( i >= start_idx and i <= end_idx ) or
( structure.GetDistance ( i , mid_idx ) < kMaxSphereDistance ) or
( structure.GetDistance ( i , end_idx ) < kMaxSphereDistance ) or
( structure.GetDistance ( i , start_idx ) < kMaxSphereDistance ) ) then
table.insert ( SphereSegList , i )
end
end
-- print ( "n side chains " , #SphereSegList )
end
function SelectSphere ( SphereSegList )
selection.DeselectAll ()
for i = 1 , #SphereSegList do
selection.Select ( SphereSegList [ i ] )
end
end
function rebuild_seg ( start_idx , end_idx )
-- return codes
--
-- 0 rebuild had no effect: seg probably locked
-- 1 normal
best_score_after_rebuild_and_shake = -99999
good_structure_found = false
SphereSegList = {} -- Speed things up when shaking
ConstructListOfResiduesWithinRange ( SphereSegList , start_idx , end_idx )
for i = 1 , number_of_rebuilds do
save.Quickload ( kOriginalStructureOrNewBest )
ConvertToLoop ( start_idx , end_idx )
selection.DeselectAll ()
selection.SelectRange ( start_idx , end_idx )
if ( monit ) then
monit_init_worst_angle = get_worst_angle ( start_idx , end_idx )
monit_worst_angle = 0
monit_score_after_rebuild = 0
monit_score_after_post_rebuild_shake = 0
end
-- for some reason, structure.RebuildSelected ( 1) often leaves the protein unchanged
k = 0
repeat
k = k + 1
structure.RebuildSelected ( k )
score_after_rebuild = GetScore ()
until math.abs ( score_after_rebuild - best_score ) >= 1e-5 or k > 2
if ( monit ) then
monit_score_after_rebuild = score_after_rebuild
monit_worst_angle = get_worst_angle ( start_idx , end_idx )
end
if ( ( math.abs ( score_after_rebuild - best_score ) >= 1e-5 ) and
( is_backbone_reasonably_angled ( start_idx , end_idx ) ) ) then
good_structure_found = true
SelectSphere ( SphereSegList )
structure.ShakeSidechainsSelected ( 1 )
selection.SelectAll ()
score_after_shake = GetScore()
if ( monit ) then
monit_score_after_post_rebuild_shake = score_after_shake
end
if ( score_after_shake > best_score_after_rebuild_and_shake ) then
best_score_after_rebuild_and_shake = score_after_shake
save.Quicksave ( kBestPostRebuild )
end
end
if ( monit ) then
print ( r3 ( best_score ) .. " , " ..
start_idx .. " , " ..
end_idx .. " , " ..
r3 ( monit_score_after_rebuild ) .. " , " ..
r3 ( monit_score_after_post_rebuild_shake ) .. " , " ..
" , " ..
" , " ..
" , " ..
" , " ..
" , " ..
r3 ( monit_init_worst_angle ) .. " , " ..
r3 ( monit_worst_angle ) .. " , " ..
r3 ( get_worst_angle ( start_idx , end_idx ) ) )
end
end -- for i
if ( good_structure_found == true ) then
return 1
else
return 0
end
end
function rebuild_all_at_specified_length ( rebuild_length )
n_possible_segs = n_residues - rebuild_length + 1
inc = Coprime ( n_possible_segs )
-- The point of this Coprime, start_idx business is so we go through the protein in a
-- nonuniform way. Rather than rebuilding 1-4, 2-5 , 3-6 for example, we might do 21 4, 53-57, 3-6
-- or something like that
for i = 1 , n_possible_segs do
end_idx = start_idx + rebuild_length - 1
if ( ( start_idx >= min_residue and end_idx <= max_residue ) and
( DoesSectionContainFixedResidues ( start_idx , end_idx ) == false ) and
( loops_only == false or
( loops_only == true and IsSegmentALoop ( start_idx , end_idx ) < 2 ) ) ) then
if ( monit ) then
monit_init_worst_angle = get_worst_angle ( start_idx , end_idx )
monit_worst_angle = 0
monit_score_after_rebuild = 0
monit_score_after_post_rebuild_shake = 0
monit_gain_from_possible_initial_nudge_wiggle = 0
monit_gain_from_wiggle = 0
monit_gain_from_shake = 0
monit_gain_from_nudge_wiggle = 0
monit_total_gain = 0
end
rb = rb + 1
if ( monit == false ) then
print ( "rb ".. rb .. " " .. start_idx .."-" .. end_idx .. " (" ..i .. "/" .. n_possible_segs .. ")" )
end
if ( rebuild_seg ( start_idx , end_idx ) == 1 ) then
post_rebuild ()
end
if ( monit ) then
print ( r3 ( best_score ) .. " , " ..
start_idx .. " , " ..
end_idx .. " , " ..
" , " ..
r3 ( monit_score_after_post_rebuild_shake ) .. " , " ..
r3 ( monit_gain_from_wiggle ) .. " , " ..
r3 ( monit_gain_from_possible_initial_nudge_wiggle ) .. " , " ..
r3 ( monit_gain_from_shake ) .. " , " ..
r3 ( monit_gain_from_nudge_wiggle ) .. " , " ..
r3 ( monit_total_gain ) .. " , " ..
r3 ( monit_init_worst_angle ) .. " , " ..
r3 ( monit_worst_angle ) .. " , " ..
r3 ( get_worst_angle ( start_idx , end_idx ) ) .. " , " ..
r3 ( GetNudgeWiggleThreshold () ) )
end
save.LoadSecondaryStructure ()
end
start_idx = start_idx + inc
if ( start_idx > n_possible_segs ) then
start_idx = start_idx - n_possible_segs
end
end -- for i
end
function get_parameters ()
local dlog = dialog.CreateDialog ( "Loop rebuild 3" )
dlog.rb_length = dialog.AddSlider ( "Rebuild length" , rebuild_length , 3 , 10 , 0 )
dlog.n_rebuilds = dialog.AddSlider ( "Number of rebuilds" , number_of_rebuilds , 1 , 50 , 0 )
dlog.min_residue = dialog.AddSlider ( "Min residue" , 1 , 1 , n_residues , 0 )
dlog.max_residue = dialog.AddSlider ( "Max residue" , n_residues , 1 , n_residues , 0 )
dlog.loops = dialog.AddCheckbox ( "Loops only" , loops_only )
dlog.ok = dialog.AddButton ( "OK" , 1 )
dlog.cancel = dialog.AddButton ( "Cancel" , 0 )
if ( dialog.Show ( dlog ) > 0 ) then
rebuild_length = dlog.rb_length.value
number_of_rebuilds = dlog.n_rebuilds.value
min_residue = dlog.min_residue.value
max_residue = dlog.max_residue.value
loops_only = dlog.loops.value
return true
else
return false
end
end
function main()
band.DisableAll ()
n_residues = structure.GetCount ()
for i = 1, n_residues do
original_secondary_structure [ i ] = structure.GetSecondaryStructure ( i )
end
save.SaveSecondaryStructure ()
rb = 0
save.Quicksave ( kOriginalStructureOrNewBest )
behavior.SetClashImportance ( 1 )
--wiggle_accuracy = behavior.GetWiggleAccuracy ()
--if ( wiggle_accuracy < 3 ) then
--wiggle_accuracy = 3
--behavior.SetWiggleAccuracy ( 3 )
--print ( "Wiggle accuracy set to 3" )
--end
if ( get_parameters () == false ) then
return -- graceful exit
end
best_score = GetScore ()
print ( "Start score : " .. r3 ( best_score ) )
if ( loops_only == true ) then
print ( "Loops only" )
else
print ( "Loops, sheets, and helices" )
end
max_residue = math.min ( n_residues , max_residue )
min_residue = math.max ( 1 , min_residue )
if ( max_residue < min_residue + rebuild_length - 1 ) then
print ( "Rebuild range error." )
print ( "Max residue (" .. max_residue .. ") must be greater than or equal to Min residue (" .. min_residue .. ")" )
print ( "plus Rebuild Length (" .. rebuild_length .. ")" )
return
end
print ( "Rebuild range " .. min_residue .. " to " .. max_residue )
GetListOfFixedResidues ()
if ( start_idx < 1 or start_idx > n_residues ) then
r = pseudorandom ( best_score )
start_idx = 1 + n_residues * r
start_idx = start_idx - start_idx % 1
end
cos_worst_allowable_angle = math.cos ( math.rad ( kWorstAllowableAngle ) )
table.insert ( nudge_wiggle_gains , kInitNudgeWiggleThreshold )
if ( monit ) then
print ( "Start score" .. " , " ..
"Start idx" .. " , " ..
"End idx" .. " , " ..
"Score after rebuild" .. " , " ..
"Score after shake" .. " , " ..
"Gain from wiggle" .. " , " ..
"Gain from possible nudge wiggle" .. " , " ..
"Gain from shake" .. " , " ..
"Gain from nudge wiggle" .. " , " ..
"Total gain" .. " , " ..
"Initial worst angle" .. " , " ..
"Post rebuild worst angle" .. " , " ..
"Final worst angle" .. " , " ..
"Nw Threshold" )
end
while true do
rebuild_all_at_specified_length ( rebuild_length )
end
save.Quickload ( kOriginalStructureOrNewBest )
band.EnableAll ()
end
function cleanup ()
print ( "Cleaning up" )
behavior.SetClashImportance ( 1.0 )
save.Quickload ( kOriginalStructureOrNewBest )
save.LoadSecondaryStructure ()
band.EnableAll ()
end
xpcall ( main , cleanup )
rebuild_length = 6 -- Length of sections to rebuild.
number_of_rebuilds = 8 -- It'll look for the best score from this number of rebuilds.
start_idx = 0 -- if this is a residue id, use it as the start point. If not
-- use a semi-random number computed from the score
loops_only = true -- If true, will still allow a single non-loop residue in the rebuild segment
-- If false, non-loop residues will be temporarily converted to loops
min_residue = 1 -- Change these if you want to rebuild part of a protein only. For example, to
-- rebuild 53-64, set min_residue to 53 and max_residue to 64. To restore, set to
-- 1 and 9999 (or any large number) respectively
max_residue = 999
kWorstAllowableAngle = 75 -- ignore rebuilds with backbone angles less than this (degrees)
kMaxSphereDistance = 14 -- Residues further than this away won't be shaken
-- Use of quicksave slots
kOriginalStructureOrNewBest = 1
kBestPostRebuild = 2
kTemp = 3
-- Globals
original_secondary_structure = {}
fixed_residues = {}
--wiggle_accuracy = 0
n_residues = 0
best_score = 0
cos_worst_allowable_angle = 0
-- NudgeWiggle : adaptive thresholds. Change the thresholds based on current protein.
kInitNudgeWiggleThreshold = 25
nudge_wiggle_gains = {}
kLowestNudgeWiggleThreshold = 8
-- Performance monitoring. Result intended for import into Excel as comma-delimited text
monit = false
monit_init_worst_angle = 0
monit_worst_angle = 0
monit_score_after_rebuild = 0
monit_score_after_post_rebuild_shake = 0
monit_gain_from_possible_initial_nudge_wiggle = 0
monit_gain_from_wiggle = 0
monit_gain_from_shake = 0
monit_gain_from_nudge_wiggle = 0
monit_total_gain = 0
-- Begin
function r3 ( x )
-- Round to 3 decimal places
t = 10 ^ 3
return math.floor ( x*t + 0.5 ) / t
end
function GetScore ()
score = current.GetEnergyScore ()
return score
end
function Coprime ( n )
-- returns a number that has no common factors with n
-- an ugly hack
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 IsSegmentALoop ( first , last )
-- Returns number of non-loop residues (helices or sheets)
n_non_loops = 0
for i = first , last do
if ( structure.GetSecondaryStructure ( i ) ~= "L" ) then
n_non_loops = n_non_loops + 1
end
end
return n_non_loops
end
function ConvertToLoop ( first , last )
for k = first , last do
if ( original_secondary_structure [ k ] ~= "L" ) then
selection.DeselectAll ( )
selection.Select ( k )
structure.SetSecondaryStructureSelected ( "L" )
end
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 GetListOfFixedResidues ()
for i = 1 , n_residues do
fixed_residues [ i ] = false
end
for i = 1 , n_residues do
seg_backbone_frozen, seg_sidechain_frozen = freeze.IsFrozen ( i )
if ( seg_backbone_frozen == true ) then
fixed_residues [ i ] = true
end
s = structure.IsLocked ( i )
if ( s == true ) then
fixed_residues [ i ] = true
end
end
end
function DoesSectionContainFixedResidues ( start_idx , end_idx )
for i = start_idx , end_idx do
if ( fixed_residues [ i ] == true ) then
return true
end
end
return false
end
function GetNudgeWiggleThreshold ()
n = #nudge_wiggle_gains
k = math.max ( n - 5 , 1 )
local x = 0
for i = k , n do
x = x + nudge_wiggle_gains [ i ]
end
av = x / ( n - k + 1 )
return math.max ( kLowestNudgeWiggleThreshold , 1.4*av )
end
function MiniFuse ()
selection.SelectAll ()
score_in_kTemp = GetScore()
save.Quicksave ( kTemp)
selection.SelectAll ()
behavior.SetClashImportance ( 0.07 )
structure.ShakeSidechainsSelected ( 1 )
behavior.SetClashImportance ( 1 )
structure.WiggleSelected ( 8 )
score = GetScore()
if ( score < score_in_kTemp ) then
save.Quickload ( kTemp )
end
end
function NudgeWiggle ( score )
selection.SelectAll ()
score_in_kTemp = GetScore()
save.Quicksave ( kTemp)
-- behavior.SetWiggleAccuracy ( 2 )
behavior.SetClashImportance ( 0.4 + 0.2 * pseudorandom ( score ) )
structure.WiggleSelected ( 1 )
-- behavior.SetWiggleAccuracy ( wiggle_accuracy )
behavior.SetClashImportance ( 1 )
structure.WiggleSelected ( 6 )
score = GetScore()
gain_from_nw = score - score_in_kTemp
if ( monit ) then
monit_gain_from_nudge_wiggle = gain_from_nw
end
if ( score < score_in_kTemp ) then
save.Quickload ( kTemp )
end
end
function WiggleShakeNudgeWiggle ()
-- returns the score after doing a WSNW: short-circuiting if things don't look promising
score_after_rebuild_shake = GetScore ()
structure.WiggleSelected ( 10 )
score_after_wiggle = GetScore ()
if ( monit ) then
monit_gain_from_wiggle = score_after_wiggle - score_after_rebuild_shake
end
if ( score_after_wiggle - score_after_rebuild_shake < 10 ) then
-- The wiggle got stuck and didn't achieve anything
NudgeWiggle ( score_after_wiggle )
score_after_nudge_wiggle = GetScore ()
if ( monit ) then
monit_gain_from_possible_initial_nudge_wiggle = score_after_nudge_wiggle - score_after_rebuild_shake
end
score_after_wiggle = score_after_nudge_wiggle
end
nwt = GetNudgeWiggleThreshold ()
if ( best_score - score_after_wiggle > 1.5*nwt ) then
return score_after_wiggle
end
structure.ShakeSidechainsSelected ( 1 )
score_after_second_shake = GetScore ()
if ( monit ) then
monit_gain_from_shake = score_after_second_shake - score_after_wiggle
end
if ( best_score - score_after_second_shake < nwt ) then
-- nwt = GetNudgeWiggleThreshold ()
start_score = GetScore()
NudgeWiggle ( score_after_second_shake )
MiniFuse ()
-- print ( "nwt " .. r3 ( nwt ) )
score = GetScore ()
if ( score - start_score > 0 ) then
table.insert ( nudge_wiggle_gains , score - start_score )
end
else
score = score_after_second_shake
end
return score
end
function post_rebuild ()
save.Quickload ( kBestPostRebuild )
if ( monit ) then
monit_score_after_post_rebuild_shake = GetScore ()
end
curr_score = WiggleShakeNudgeWiggle ()
if ( monit ) then
monit_total_gain = curr_score - best_score
end
if ( curr_score > best_score ) then
best_score = curr_score
save.LoadSecondaryStructure ()
save.Quicksave ( kOriginalStructureOrNewBest )
print ( "Improvement to ".. r3 ( best_score ) )
end
end
function get_cos_backbone_angle ( i )
if ( i == 1 or i == n_residues ) then
return -1
else
-- Cosine rule
a = structure.GetDistance ( i - 1 , i )
b = structure.GetDistance ( i , i + 1 )
c = structure.GetDistance ( i - 1 , i + 1 )
cos_angle = ( a*a + b*b - c*c ) / ( 2*a*b )
return cos_angle
end
end
function get_worst_angle ( start_idx , end_idx )
-- For monitoring only
cos_worst_angle = -1
for i = start_idx , end_idx do
cos_worst_angle = math.max ( cos_worst_angle , get_cos_backbone_angle ( i ) )
end
return math.deg ( math.acos ( cos_worst_angle ) )
end
function is_backbone_reasonably_angled ( start_idx , end_idx )
for i = start_idx , end_idx do
if ( get_cos_backbone_angle ( i ) > cos_worst_allowable_angle ) then
return false
end
end
return true
end
function ConstructListOfResiduesWithinRange ( SphereSegList , start_idx , end_idx )
mid_idx = math.floor ( ( end_idx + start_idx ) / 2 )
-- Look at start_idx, end_idx, and mid_idx only to save time at the expense of a small amount of accuracy
for i = 1 , n_residues do
if ( ( i >= start_idx and i <= end_idx ) or
( structure.GetDistance ( i , mid_idx ) < kMaxSphereDistance ) or
( structure.GetDistance ( i , end_idx ) < kMaxSphereDistance ) or
( structure.GetDistance ( i , start_idx ) < kMaxSphereDistance ) ) then
table.insert ( SphereSegList , i )
end
end
-- print ( "n side chains " , #SphereSegList )
end
function SelectSphere ( SphereSegList )
selection.DeselectAll ()
for i = 1 , #SphereSegList do
selection.Select ( SphereSegList [ i ] )
end
end
function rebuild_seg ( start_idx , end_idx )
-- return codes
--
-- 0 rebuild had no effect: seg probably locked
-- 1 normal
best_score_after_rebuild_and_shake = -99999
good_structure_found = false
SphereSegList = {} -- Speed things up when shaking
ConstructListOfResiduesWithinRange ( SphereSegList , start_idx , end_idx )
for i = 1 , number_of_rebuilds do
save.Quickload ( kOriginalStructureOrNewBest )
ConvertToLoop ( start_idx , end_idx )
selection.DeselectAll ()
selection.SelectRange ( start_idx , end_idx )
if ( monit ) then
monit_init_worst_angle = get_worst_angle ( start_idx , end_idx )
monit_worst_angle = 0
monit_score_after_rebuild = 0
monit_score_after_post_rebuild_shake = 0
end
-- for some reason, structure.RebuildSelected ( 1) often leaves the protein unchanged
k = 0
repeat
k = k + 1
structure.RebuildSelected ( k )
score_after_rebuild = GetScore ()
until math.abs ( score_after_rebuild - best_score ) >= 1e-5 or k > 2
if ( monit ) then
monit_score_after_rebuild = score_after_rebuild
monit_worst_angle = get_worst_angle ( start_idx , end_idx )
end
if ( ( math.abs ( score_after_rebuild - best_score ) >= 1e-5 ) and
( is_backbone_reasonably_angled ( start_idx , end_idx ) ) ) then
good_structure_found = true
SelectSphere ( SphereSegList )
structure.ShakeSidechainsSelected ( 1 )
selection.SelectAll ()
score_after_shake = GetScore()
if ( monit ) then
monit_score_after_post_rebuild_shake = score_after_shake
end
if ( score_after_shake > best_score_after_rebuild_and_shake ) then
best_score_after_rebuild_and_shake = score_after_shake
save.Quicksave ( kBestPostRebuild )
end
end
if ( monit ) then
print ( r3 ( best_score ) .. " , " ..
start_idx .. " , " ..
end_idx .. " , " ..
r3 ( monit_score_after_rebuild ) .. " , " ..
r3 ( monit_score_after_post_rebuild_shake ) .. " , " ..
" , " ..
" , " ..
" , " ..
" , " ..
" , " ..
r3 ( monit_init_worst_angle ) .. " , " ..
r3 ( monit_worst_angle ) .. " , " ..
r3 ( get_worst_angle ( start_idx , end_idx ) ) )
end
end -- for i
if ( good_structure_found == true ) then
return 1
else
return 0
end
end
function rebuild_all_at_specified_length ( rebuild_length )
n_possible_segs = n_residues - rebuild_length + 1
inc = Coprime ( n_possible_segs )
-- The point of this Coprime, start_idx business is so we go through the protein in a
-- nonuniform way. Rather than rebuilding 1-4, 2-5 , 3-6 for example, we might do 21 4, 53-57, 3-6
-- or something like that
for i = 1 , n_possible_segs do
end_idx = start_idx + rebuild_length - 1
if ( ( start_idx >= min_residue and end_idx <= max_residue ) and
( DoesSectionContainFixedResidues ( start_idx , end_idx ) == false ) and
( loops_only == false or
( loops_only == true and IsSegmentALoop ( start_idx , end_idx ) < 2 ) ) ) then
if ( monit ) then
monit_init_worst_angle = get_worst_angle ( start_idx , end_idx )
monit_worst_angle = 0
monit_score_after_rebuild = 0
monit_score_after_post_rebuild_shake = 0
monit_gain_from_possible_initial_nudge_wiggle = 0
monit_gain_from_wiggle = 0
monit_gain_from_shake = 0
monit_gain_from_nudge_wiggle = 0
monit_total_gain = 0
end
rb = rb + 1
if ( monit == false ) then
print ( "rb ".. rb .. " " .. start_idx .."-" .. end_idx .. " (" ..i .. "/" .. n_possible_segs .. ")" )
end
if ( rebuild_seg ( start_idx , end_idx ) == 1 ) then
post_rebuild ()
end
if ( monit ) then
print ( r3 ( best_score ) .. " , " ..
start_idx .. " , " ..
end_idx .. " , " ..
" , " ..
r3 ( monit_score_after_post_rebuild_shake ) .. " , " ..
r3 ( monit_gain_from_wiggle ) .. " , " ..
r3 ( monit_gain_from_possible_initial_nudge_wiggle ) .. " , " ..
r3 ( monit_gain_from_shake ) .. " , " ..
r3 ( monit_gain_from_nudge_wiggle ) .. " , " ..
r3 ( monit_total_gain ) .. " , " ..
r3 ( monit_init_worst_angle ) .. " , " ..
r3 ( monit_worst_angle ) .. " , " ..
r3 ( get_worst_angle ( start_idx , end_idx ) ) .. " , " ..
r3 ( GetNudgeWiggleThreshold () ) )
end
save.LoadSecondaryStructure ()
end
start_idx = start_idx + inc
if ( start_idx > n_possible_segs ) then
start_idx = start_idx - n_possible_segs
end
end -- for i
end
function get_parameters ()
local dlog = dialog.CreateDialog ( "Loop rebuild 3" )
dlog.rb_length = dialog.AddSlider ( "Rebuild length" , rebuild_length , 3 , 10 , 0 )
dlog.n_rebuilds = dialog.AddSlider ( "Number of rebuilds" , number_of_rebuilds , 1 , 50 , 0 )
dlog.min_residue = dialog.AddSlider ( "Min residue" , 1 , 1 , n_residues , 0 )
dlog.max_residue = dialog.AddSlider ( "Max residue" , n_residues , 1 , n_residues , 0 )
dlog.loops = dialog.AddCheckbox ( "Loops only" , loops_only )
dlog.ok = dialog.AddButton ( "OK" , 1 )
dlog.cancel = dialog.AddButton ( "Cancel" , 0 )
if ( dialog.Show ( dlog ) > 0 ) then
rebuild_length = dlog.rb_length.value
number_of_rebuilds = dlog.n_rebuilds.value
min_residue = dlog.min_residue.value
max_residue = dlog.max_residue.value
loops_only = dlog.loops.value
return true
else
return false
end
end
function main()
band.DisableAll ()
n_residues = structure.GetCount ()
for i = 1, n_residues do
original_secondary_structure [ i ] = structure.GetSecondaryStructure ( i )
end
save.SaveSecondaryStructure ()
rb = 0
save.Quicksave ( kOriginalStructureOrNewBest )
behavior.SetClashImportance ( 1 )
--wiggle_accuracy = behavior.GetWiggleAccuracy ()
--if ( wiggle_accuracy < 3 ) then
--wiggle_accuracy = 3
--behavior.SetWiggleAccuracy ( 3 )
--print ( "Wiggle accuracy set to 3" )
--end
if ( get_parameters () == false ) then
return -- graceful exit
end
best_score = GetScore ()
print ( "Start score : " .. r3 ( best_score ) )
if ( loops_only == true ) then
print ( "Loops only" )
else
print ( "Loops, sheets, and helices" )
end
max_residue = math.min ( n_residues , max_residue )
min_residue = math.max ( 1 , min_residue )
if ( max_residue < min_residue + rebuild_length - 1 ) then
print ( "Rebuild range error." )
print ( "Max residue (" .. max_residue .. ") must be greater than or equal to Min residue (" .. min_residue .. ")" )
print ( "plus Rebuild Length (" .. rebuild_length .. ")" )
return
end
print ( "Rebuild range " .. min_residue .. " to " .. max_residue )
GetListOfFixedResidues ()
if ( start_idx < 1 or start_idx > n_residues ) then
r = pseudorandom ( best_score )
start_idx = 1 + n_residues * r
start_idx = start_idx - start_idx % 1
end
cos_worst_allowable_angle = math.cos ( math.rad ( kWorstAllowableAngle ) )
table.insert ( nudge_wiggle_gains , kInitNudgeWiggleThreshold )
if ( monit ) then
print ( "Start score" .. " , " ..
"Start idx" .. " , " ..
"End idx" .. " , " ..
"Score after rebuild" .. " , " ..
"Score after shake" .. " , " ..
"Gain from wiggle" .. " , " ..
"Gain from possible nudge wiggle" .. " , " ..
"Gain from shake" .. " , " ..
"Gain from nudge wiggle" .. " , " ..
"Total gain" .. " , " ..
"Initial worst angle" .. " , " ..
"Post rebuild worst angle" .. " , " ..
"Final worst angle" .. " , " ..
"Nw Threshold" )
end
while true do
rebuild_all_at_specified_length ( rebuild_length )
end
save.Quickload ( kOriginalStructureOrNewBest )
band.EnableAll ()
end
function cleanup ()
print ( "Cleaning up" )
behavior.SetClashImportance ( 1.0 )
save.Quickload ( kOriginalStructureOrNewBest )
save.LoadSecondaryStructure ()
band.EnableAll ()
end
xpcall ( main , cleanup )