Code
-- Random Squeezes
-- 12/17/12 jeff101 adapted Random Squeezes v1.7.0 to 1.7.1 to vary ci less randomly
-- last updated 100pm 12/18/12
-- Generic functions and objects
function repeat_string(s, n)
return n >= .5 and s..repeat_string(s, n-1) or ""
end
function width_of_char(c)
local chars_by_len = {
[0.4] = "ijl",
[0.468] = "|",
[0.5] = " .,;:![]Ift",
[0.599] = "-()",
[0.6] = "{}r",
[0.7] = "*",
[0.843] = "^",
[0.9] = "Jcksvxyz",
[1.05] = "+<>",
[1.099] = "FTZ",
[1.2] = "&ABEKPSVXY",
[1.3] = "CDHNRUw",
[1.4] = "GOQ",
[1.499] = "Mm",
[1.6] = "%",
[1.699] = "W",
[1.828] = "@",
}
for len,chars in pairs(chars_by_len) do
if chars:find(c, nil, true) then
return len
end
end
return 1
end
function width_of_string(s)
local sum = 0
for i=1,s:len() do
sum = sum + width_of_char(s:sub(i,i))
end
return sum
end
alignment_spill = 0
function alignment_reset()
alignment_spill = 0
end
function align(text, width, left)
local fill = (width - width_of_string(text)) / width_of_char(" ") + alignment_spill
local spaces = math.max(0, math.floor(fill+.5))
alignment_spill = fill - spaces
if alignment == "left" then
return text..repeat_string(" ", spaces)
elseif alignment == "right" then
return repeat_string(" ", spaces)..text
else
return repeat_string(" ", math.floor(spaces/2))..text..repeat_string(" ", math.ceil(spaces/2))
end
end
function format(value, style)
-- style fields:
-- integer digits: number of digits before decimal point
-- integer decimals: number of digits after decimal point, default 0
-- boolean left: aligned, default true for text, false for number
-- boolean sign: sign even if positive, default false
local decimals = style.decimals or 0
local width = (style.sign and width_of_char("+") or 0) + (style.digits+decimals) * width_of_char("0") + (decimals > 0 and width_of_char(".") or 0)
if type(value) == "number" then
local fmt = "%"..(style.sign and "+" or "").."."..decimals.."f"
local text = string.format(fmt, value)
return align(text, width, style.left and "left" or "right")
elseif value == "-" then
return align(value, width, "center")
else
return align(value or "", width, "left")
end
end
function is_flexible_segment(si)
if freeze.IsFrozen(si) then
return nil -- might be flexible but frozen by user
else
freeze.Freeze(si, true, false)
if freeze.IsFrozen(si) then
freeze.Unfreeze(si, true, false)
return true -- flexible
else
return false -- intrisically frozen (as in ligand puzzles)
end
end
end
function has_flexible_backbone()
-- If we test every segment until we find one flexible,
-- we rapidly flush the undo buffer in ligand puzzles.
-- So only test the last segment and make assumptions.
local f = is_flexible_segment(structure.GetCount())
if f then
-- At least one is flexible, plain and simple.
return true
elseif f == nil then
-- One segment is frozen by user; assume that there are
-- other segments around not frozen by user and not
-- Intrinsically frozen.
return false
else
-- intrinsically frozen segments tend to be numbered first;
-- so if the last one is Intrinsically frozen, assume all are.
return false
end
end
function has_mutable_sidechain()
for si = 1, structure.GetCount() do
if structure.IsMutable(si) then
return true
end
end
return false
end
band_mgr = { bands_enabled={}, num_bands_enabled=0, num_bands_disabled=0 }
for bi = 1, band.GetCount() do
e = band.IsEnabled(bi)
band_mgr.bands_enabled[bi] = e
band_mgr.num_bands_enabled = band_mgr.num_bands_enabled + (e and 1 or 0)
band_mgr.num_bands_disabled = band_mgr.num_bands_disabled + (e and 0 or 1)
end
function band_mgr.is_any_enabled()
for bi = 1, band.GetCount() do
if band.IsEnabled(bi) then
return true
end
end
return false
end
function band_mgr.set_enabled(enabled)
if band_mgr.is_any_enabled() ~= enabled then
if enabled then
band_mgr.restore_bands()
else
band.DisableAll()
end
end
end
function band_mgr.restore_bands()
if band_mgr.num_bands_enabled > band_mgr.num_bands_disabled then
band.EnableAll()
for bi, e in pairs(band_mgr.bands_enabled) do
if not e then
band.Disable(bi)
end
end
else
for bi, e in pairs(band_mgr.bands_enabled) do
if e then
band.Enable(bi)
end
end
end
end
-- Constants
min_worthy_score = 1e-4
bottom_score = -1e5
slot_mission, slot_dive, slot_wiggle = 1, 2, 3
buckets = 2
time_width = 4
style_iter = { digits=3 }
style_ci = { digits=1, decimals=6 }
style_depth = { digits=1 }
style_em = { digits=1, decimals=5 }
style_score = { digits=5, decimals=3 }
style_gain = { digits=4, sign=true } -- decimals determined by state of desparation
style_relative_score = { digits=5, decimals=4, sign=true }
-- Main
initial_CI = behavior.GetClashImportance()
is_exploration = current.GetExplorationMultiplier() > 0
is_mutable = has_mutable_sidechain()
print('Running Random Squeezes v1.7.1 at '..(os.date())..'\n on puzzle '..(puzzle.GetName())..'\n with Puzzle ID '..(puzzle.GetPuzzleID())..', '..(structure.GetCount())..' residues, and '..(band.GetCount())..' bands:\n ')
dbox = dialog.CreateDialog("Mode of operation")
dbox.minci = dialog.AddSlider("Min CI", 1e-2, 0, 1, 2)
dbox.maxpw = dialog.AddSlider("Max D", 1, 1, 9, 0)
dbox.labl1 = dialog.AddLabel("D = depth = maximum #iterations of perturbing wiggle\n i.e. wiggle with ci < 1.")
dbox.shake = dialog.AddCheckbox((is_mutable and "Mutate" or "Shake").." after substantial gain (otherwise never)", true)
if is_exploration then
dbox.embound = dialog.AddCheckbox("Reject lower exploration multiplier", true)
end
if band_mgr.num_bands_enabled > 0 then
dbox.bands_down = dialog.AddCheckbox("Respect bands while perturbing", true)
dbox.bands_up = dialog.AddCheckbox("Respect bands while recovering", true)
end
dbox.cibins = dialog.AddCheckbox("Vary ci less randomly", false)
dbox.basebin = dialog.AddSlider("Base Bin", 2, 2, 10, 0)
dbox.labl2 = dialog.AddLabel("Base Bin = first # of bins when vary ci less randomly.")
dbox.verbose = dialog.AddCheckbox("Give more detailed Recipe Output", false)
--dbox.model = dialog.AddLabel("Happy with gain:")
dbox.mode1 = dialog.AddButton("Normal", 0)
dbox.mode2 = dialog.AddButton("Desperate", 3)
dbox.cancl = dialog.AddButton("Cancel", -99)
style_gain.decimals = dialog.Show(dbox)
if style_gain.decimals == -99 then return end
if dbox.cibins.value then
print("Below will vary ci less randomly than usual.\n ")
else
print("Below will vary ci the usual way.\n ")
end
enable_bands_down = band_mgr.num_bands_enabled > 0 and dbox.bands_down.value
enable_bands_up = band_mgr.num_bands_enabled > 0 and dbox.bands_up.value
minem = is_exploration and dbox.embound.value and current.GetExplorationMultiplier() or 0
initial_best_score = absolutebest.GetScore()
score_granting_repeat = 10 ^ -style_gain.decimals
score_granting_shake = math.max(10 * score_granting_repeat, 1)
initial_score = current.GetScore()
bottom_ci = dbox.minci.value
initial_time = os.time()
iteration = 0
function header()
local time = align("Time", time_width)
local iter = format("Iter", style_iter)
local ci = format("CI", style_ci)
local depth = dbox.maxpw.value > 1 and " "..format("D", style_depth) or ""
local em = is_exploration and " "..format("Stability", style_em) or ""
local gain = format("Gain", style_gain)
print(time.." "..iter.." "..ci..depth..em.." "..gain.."Score versus track best")
report(0, "-", 0, 0, "")
end
function report(iteration, ci, depth, gain, tip)
local last_reported_time = os.time()
local minutes = os.difftime(last_reported_time, initial_time) / 60
local time = align(string.format("%i", minutes/60)..":"..string.format("%02i", minutes%60), time_width)
local iter = format(iteration, style_iter)
ci = format(ci, style_ci)
depth = dbox.maxpw.value > 1 and " "..format(depth, style_depth) or ""
if gain ~= 0 then
gain = format(gain, style_gain)
else
gain = format("-", { digits=style_gain.digits, decimals=style_gain.decimals, center=true})
end
local score = current.GetScore()
local relative_score = score - initial_best_score
local em = is_exploration and " "..format(current.GetExplorationMultiplier(), style_em) or ""
alignment_reset()
print(time.." "..iter.." "..ci..depth..em.." "..format(gain, style_gain).." "..format(trunc3(score), style_score)..format(relative_score, style_relative_score).." "..tip)
end
dive = {}
function dive.start(ci, depth)
save.Quicksave(slot_dive)
dive.ref_score = current.GetScore()
dive.ci = ci
dive.depth = depth
end
function dive.abort()
recentbest.Restore()
end
recentbest.Save()
function dive.wiggle_up(bands_enabled)
band_mgr.set_enabled(bands_enabled)
local score = current.GetScore()
--save.Quicksave(slot_wiggle)
local min_win = min_worthy_score
if score < dive.ref_score then
min_win = dive.ref_score - score
end
local its_per_depth = { 1, 3, 9, 27 }
local its_divisor = 32 -- should be less than sum of its_per_depth
local depth = 0
local prev_score
repeat
depth = depth + 1
prev_score = score
structure.WiggleAll(its_per_depth[depth], true, true)
local score1 = current.GetScore()
recentbest.Restore()
local score2 = current.GetScore()
if score2 > prev_score then
score = score2
else
score = score1
--save.Quickload(slot_wiggle)
band_mgr.set_enabled(bands_enabled)
end
until depth == #its_per_depth or score < prev_score + min_win * its_per_depth[depth] / its_divisor
end
function dive.finish()
if dive.ref_score ~= nil then
local em = current.GetExplorationMultiplier()
local gain = current.GetScore() - dive.ref_score
local tip
if em < minem or gain < 0 then
tip = ":("
save.Quickload(slot_dive)
elseif gain == 0 then
tip = ":|"
else
tip = ":)"
save.Quicksave(slot_dive)
minem = em
end
report(iteration, dive.ci, dive.depth, gain, tip)
dive.ref_score = nil
return gain
end
end
-- below will show the sequence and secondary structure
function showss()
local j,aa,ss,j1s,j10s
local numstr=" "
local aastr="AAs="
local ssstr="SSs="
local tot=structure.GetCount()
print('Puzzle has '..tot..' residues with AAs and SSs as below:')
for j=1,tot do
aa=structure.GetAminoAcid(j)
ss=structure.GetSecondaryStructure(j)
j10s=math.floor(j/10)
j1s=j-10*j10s
numstr=(numstr..j1s)
aastr=(aastr..aa)
ssstr=(ssstr..ss)
if j1s==0 and j<tot then
numstr=(numstr..' ')
aastr=(aastr..' ') -- add blank space after every 10th one
ssstr=(ssstr..' ')
end -- if j
end -- for j
print(numstr)
print(aastr)
print(ssstr)
end -- showss()
-- below rounds val down to nearest 0.001
function trunc3(val)
return val-(val % 0.001)
end -- trunc3()
-- below rounds to 2 decimal places
function round2(numi)
local numo=round(numi*100)/100
return numo
end -- round2()
-- below rounds val to the nearest integer to get res
function round(val)
local res=math.floor(val)
if val-res>=0.5 then
res=res+1
end -- if val
return res
end -- round()
-- below returns an array named res of length dim with all elements equal to 0
function zeros(dim)
local res = {}
local num
for num=1, dim do
res[num] = 0
end -- for num
return res
end -- function zeros
-- below determines what bin (binval) a certain ci value (valci) lies in
function getbin(valci, min_ci, ncibins, pflag)
local binval = math.ceil(ncibins * (valci - min_ci) / (1 - min_ci))
-- ideally binval is always from 1 to ncibins, but what if valci==min_ci or valci==1 ?
if binval<1 then
binval=1
end
if binval>ncibins then
binval=ncibins
end
if pflag==1 and dbox.verbose.value then
print('ci='..format(valci, style_ci)..' is in bin '..binval..' of '..ncibins..' for ci='..min_ci..' to 1.')
end -- if pflag
return binval
end -- function getbin
-- below lists how well present conformation obeys user-supplied bands
function listdiff()
local nbands, lens, glens, strs, diffs, ison, tots, avgdiff, rmsdiff, btype, mindiff, maxdiff
local maxstr, minstr, num, str, diff, diff2, indx, btype, typenams
nbands=band.GetCount()
lens={}
glens={}
strs={}
diffs={}
ison={}
tots={}
avgdiff={}
rmsdiff={}
for btype=0,2 do
avgdiff[btype]=0
rmsdiff[btype]=0
tots[btype]=0
end -- for btype
mindiff={}
maxdiff={}
maxstr={}
minstr={}
for num=1,nbands do
ison[num]=band.IsEnabled(num)
lens[num]=band.GetLength(num)
glens[num]=band.GetGoalLength(num)
str=band.GetStrength(num)
strs[num]=str
diff=math.abs(lens[num]-glens[num])
diffs[num]=diff
diff2=diff*diff
if ison[num]==true then
btypes={1,2} -- for enabled bands
else
btypes={0,2} -- for disabled bands
end -- if ison
for indx=1,2 do
btype=btypes[indx]
avgdiff[btype]=avgdiff[btype]+diff
rmsdiff[btype]=rmsdiff[btype]+diff2
tots[btype]=tots[btype]+1
if tots[btype]==1 then
mindiff[btype]=diff
maxdiff[btype]=diff
maxstr[btype]=str
minstr[btype]=str
else
if diff>maxdiff[btype] then
maxdiff[btype]=diff
elseif diff<mindiff[btype] then
mindiff[btype]=diff
end
if str>maxstr[btype] then
maxstr[btype]=str
elseif str<minstr[btype] then
minstr[btype]=str
end
end -- if tots
end -- for indx
end -- for num
for btype=0,2 do
avgdiff[btype]=avgdiff[btype]/tots[btype]
rmsdiff[btype]=math.sqrt(rmsdiff[btype]/tots[btype])
end
print('For '..nbands..' user-supplied bands and score '..trunc3(current.GetScore())..' we have:')
typenams={}
typenams[0]='disabled'
typenams[1]=' enabled'
typenams[2]='any-type'
for btype=0,2 do
if tots[btype]>0 then
print(' '..tots[btype]..' '..typenams[btype]..' bands with diffs='..round2(mindiff[btype])..'-'..round2(maxdiff[btype])..', avgdiff='..round2(avgdiff[btype])..', rmsdiff='..round2(rmsdiff[btype])..', and strengths='..round2(minstr[btype])..'-'..round2(maxstr[btype])..'.')
end -- if tots
end -- for btype
end -- listdiff()
math.randomseed(os.time())
save.Quicksave(slot_mission)
save.Quicksave(slot_dive)
-- save.SaveSecondaryStructure()
mission = {}
function mission.finalize(why)
dive.abort() -- needed on error or cancel
band_mgr.restore_bands()
-- save.LoadSecondaryStructure()
behavior.SetClashImportance(initial_CI)
local win = current.GetScore() - initial_score
if win > 0 then
-- push single step undo solutions
save.Quicksave(slot_dive)
save.Quickload(slot_mission)
save.Quickload(slot_dive)
print(why..string.format(". Gained %.6f over initial %.3f", win, trunc3(initial_score)))
else
print(why..". Sorry, couldn't help...")
save.Quickload(slot_mission) -- different if bands were disabled
end
end
function mission.Main()
if dbox.verbose.value then
showss()
print(' ')
if band.GetCount()>0 then
listdiff()
print(' ')
end
end
header()
if disable_bands then
band.DisableAll()
end
local last_shaken_score = initial_score
local min_ci = bottom_ci
local fruitless_iterations = 0
local bucket = 1
local ncibins, cibindone, oldci, noldci, isnewciok, ci, cibin, oldcinum
repeat -- start of main loop
-- below picks next ci value to use
if dbox.cibins.value then -- vary ci less randomly
if fruitless_iterations == 0 then
ncibins = dbox.basebin.value -- the number of bins to break ci's range into
-- above makes ncibins be dbox.basebin.value raised to some power
cibindone = zeros(ncibins) -- cibindone is 0 if bin hasn't been done, 1 if it has
oldci = {} -- this will list the ci's used since the last gain in score
noldci = 0 -- the number of ci's used since last gain in score
end -- if fruitless_iterations
isnewciok = 0
repeat -- below tries ci at random until it finds one in a range not done yet
ci = min_ci + (math.random()) * (1 - min_ci)
cibin = getbin(ci, min_ci, ncibins, 0) -- use 0 so won't print
if cibindone[cibin] == 0 then -- it found a ci from a bin that hasn't been done yet
isnewciok = 1 -- this is 1 if the new ci is ok, 0 if the new ci is not ok
cibindone[cibin] = 1
noldci = noldci + 1
oldci[noldci] = ci
cibin = getbin(ci, min_ci, ncibins, 1) -- use 1 so it prints
if noldci == ncibins then
ncibins = ncibins * dbox.basebin.value -- make ncibins be a power of dbox.basebin.value
cibindone = zeros(ncibins)
for oldcinum = 1, noldci do
cibin = getbin(oldci[oldcinum], min_ci, ncibins, 1) -- use 1 so it prints
cibindone[cibin] = 1
end -- for oldcinum
end -- if noldci
end -- if cibindone
until isnewciok == 1
else -- use original way
bucket = 1 - bucket -- have bucket alternate between 0 and 1 when buckets equals 2
ci = min_ci + (bucket + math.random()) * (1 - min_ci) / buckets
end -- if dbox.cibins.value
-- above picks next ci value to use
local depth = math.random(dbox.maxpw.value)
local best_score = current.GetScore()
local gain = 0
repeat -- start of inner loop
iteration = iteration + 1
dive.start(ci, depth)
behavior.SetClashImportance(ci)
band_mgr.set_enabled(enable_bands_down)
structure.WiggleAll(depth, true, false)
behavior.SetClashImportance(1)
score = current.GetScore()
if score < bottom_score then
dive.abort()
report(iteration, dive.ci, dive.depth, 0, "ouch!")
min_ci = math.max(ci + bottom_ci, ci * 1.01)
else
dive.wiggle_up(enable_bands_up)
gain = dive.finish()
local score = current.GetScore()
if score >= last_shaken_score + score_granting_shake then
if dbox.shake.value then
band.DisableAll()
if is_mutable then
structure.MutateSidechainsAll(1)
else
structure.ShakeSidechainsAll(1)
end
band_mgr.restore_bands()
end
last_shaken_score = current.GetScore()
if dbox.shake.value then
report("", is_mutable and "mutate" or "shake", 0, last_shaken_score-score, "")
end
min_ci = math.max(bottom_ci, min_ci * 0.9921875)
end
-- depth = depth + 1
end
until gain < score_granting_repeat -- end of inner loop
if current.GetScore() > best_score then
fruitless_iterations = 0
if dbox.verbose.value and band.GetCount()>0 then
print(' ')
listdiff()
print(' ')
end
else
fruitless_iterations = fruitless_iterations + 1
end
if min_ci == 1 then break end
until min_ci == 1 or fruitless_iterations > 999 -- end of main loop
mission.finalize("Petered out")
end
function mission.Catch(what) -- on cancel or error
if what:match("Cancelled$") then
dive.finish()
mission.finalize("Cancelled")
what = nil -- value doesn't matter, somehow suppresses further error message
else
print(what)
mission.finalize("Error")
end
end
-- selection.SelectAll()
-- structure.SetSecondaryStructureSelected("l")
xpcall(mission.Main, mission.Catch)
