Recipe: Helix twister 2.0 (NC)

created by spvincent

Profile

• Name: Helix twister 2.0 (NC)
• ID: 48501
• Shared with: Public
• Parent: None
• Children:
  • Helix twister cut 3.0
  • Helix twister 2.01 (NC)
  • Helix twister 2.0 Filter
• Created on: March 11, 2014 at 18:06 PM UTC
• Updated on: March 11, 2014 at 18:06 PM UTC
• Description:

  See first comment

Code

-- Globals

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 = false

-- Quicksave slots

kOriginalStructureOrNewBest = 1	  -- the starting structure, or any subsequent improvement, will be stored in this quicksave slot


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 )

    selection.SelectRange ( start_idx , end_idx )

    WiggleBandedProtein ( helix_index )

    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 ( "Helix twister 2.0" )

    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 )  
    end

    dlog.clk = dialog.AddCheckbox ( "Clockwise " , clk_twist )  

    dlog.credit_check = dialog.AddCheckbox ( "Check conditions met " , test_for_credit )  

    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
        end

        test_for_credit = dlog.credit_check.value
        keep_trying = dlog.keep_trying.value

        return true
    else
       return false
     end

end

function main ()

    print (  "Helix twister 2.0" )

    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

            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 )

Comments

[spvincent]

This is an attempt to mimic some of the functionality of the tweak tool using bands. It doesn't take long to run and often produces a few points in the endgame but doesn't seem to be a big point getter.

Changes since 1.0: generally made more NC friendly: suggest using at MWP. Also eliminated some useless options and needless complexity.

Algorithm outline

For each helix
    Bracket it with loops of the specified size
    For each residue of the helix, create a band from the alpha carbon into empty space parallel to the backbone (more or less: not sure if I've understood the spherical coordinate system correctly)
    Wiggle a region making up the helix plus spacers on other side.
    Wiggle
    Delete bands
    Shake or mutate
    Wiggle sidechains 
    Wiggle
    Record any improvements as appropriate
end

Repeat in other direction (clockwise/anticlockwise)

Parameters

Loop spacer

Just as when using the Tweak tool you have to bracket a sheet or helix with loop segments; so you can do something similar here.

Target loss

On banded wiggle, aim to lose this amount of points. Band strength will be adjusted if its wildly off but it won't make extraordinary efforts to get

Helix list

All selected by default</b>

Mutate not shake

For design puzzles only. If set will replace Shake after banded wiggle by Mutate as indicated in the algorithm outline above.

Clockwise

Whether a clockwise twist is done first or an anti-clockwise one. Kind of unnecessary.

Check conditions met

When set, will only accept a score improvement if all conditions are met. For puzzles where some restrictions are imposed.

Keep trying

Repeats everything indefinitely with minor random changes to band length and band strength.