Recipe: Local Mutate Shake or Rebuild 1.0

created by Chris Klassen

Profile

• Name: Local Mutate Shake or Rebuild 1.0
• ID: 103496
• Shared with: Public
• Parent: None
• Children: None
• Created on: May 14, 2020 at 15:29 PM UTC
• Updated on: May 14, 2020 at 15:29 PM UTC
• Description:

  Based on Local Mutate 3.2.0 allows shaking or rebuild as options.

Code

min_residue = 1
max_residue =  999
n_residues = 0

delta_residue = 9
kMaxDelta = 20

best_score = 0
clash_importance_during_modify = 0.09

kWiggleIterations = 15 

modify_count = 10
discard_limit = 0

original_slow_filters_setting = 0

score_type = 1

loss_from_modify = {}

original_residues = {}
new_residues = {}

--scores = {}
ids = {}

start_time = 0

wiggle_sidechains_first = false

--    Amino acids. Hydrophobics first.
three_letter_code_table = {g = "Gly", a = "Ala", v = "Val",
                           l = "Leu", i = "Ile", m = "Met",
                           f = "Phe", w = "Trp", p = "Pro",
                           s = "Ser", t = "Thr", c = "Cys",
                           y = "Tyr", n = "Asn", q = "Gln",
                           d = "Asp", e = "Glu", k = "Lys",
                           r = "Arg", h = "His"}


function aa1_to_aa3(aa1)
    local aa3 = three_letter_code_table[aa1]
    if not aa3 then
       aa3 = "Unk"
    end
    return aa3
end


function get_aa3 ( i )
    local aa1 = structure.GetAminoAcid (  i )
    return three_letter_code_table[aa1]
end

-- Use of quicksave slots

kOriginalStructureOrNewBest = 10
base_primary_sequence = {}
n_sequence_changes = 0


function r3 ( x )
    -- Round to 3 decimal places
    local t  = 10 ^ 3
    return math.floor ( x*t + 0.5 ) /  t 
end


function GetScore()
    if ( score_type == 1 ) then
        score = current.GetScore () 
    elseif ( score_type == 2 ) then
        behavior.SetFiltersDisabled ( false )
        score = current.GetScore () 
        behavior.SetFiltersDisabled ( true )
    end

    return score
end


function GetBasePrimarySequence()
    for seg_idx = 1, n_residues do
        base_primary_sequence[ seg_idx ] = structure.GetAminoAcid( seg_idx )
    end
end


function GetNumberOfPrimarySequenceChanges()
    local n_changes = 0
    for seg_idx = 1, n_residues do
        if (base_primary_sequence[ seg_idx ] ~= structure.GetAminoAcid( seg_idx )) then
            n_changes = n_changes + 1
        end
    end

    return n_changes
end


function PostModify()
    behavior.SetClashImportance( 1.0 )

    selection.SelectAll ()
    if ( wiggle_sidechains_first == true ) then
        structure.WiggleSelected (1, false, true )    -- sidechains
    end
 
    score_before_wiggle = GetScore ()   
    structure.WiggleSelected ( kWiggleIterations )        
    score_after_wiggle = GetScore ()
end


function SortResiduesByDistanceFrom(from_idx)
    local distances = {}
    for to_idx = 1, structure.GetCount() do
        distance = structure.GetDistance (from_idx, to_idx)
        table.insert(distances, {distance = distance, to_idx = to_idx})
    end
       
    table.sort(distances, function (k1, k2) return k1.distance < k2.distance end)
    
    local indexes_in_distance_order = {}
    
    table.foreachi(distances, function(key, val) table.insert(indexes_in_distance_order, val.to_idx) end)
    return indexes_in_distance_order
end


function PrintResidueChanges()
    local 	seg_idx
    for seg_idx = 1, n_residues do
        if ( base_primary_sequence [ seg_idx ] ~= structure.GetAminoAcid ( seg_idx ) ) then
            print ( "         " .. seg_idx .. "  "  .. aa1_to_aa3 (  base_primary_sequence [ seg_idx ] ) .. " -> " .. aa1_to_aa3 ( structure.GetAminoAcid ( seg_idx ) ) )
        end
    end
end


function MutateResidue( seg_idx )
    local	id_idx

    save.Quickload ( kOriginalStructureOrNewBest )
    
    ids = SortResiduesByDistanceFrom(seg_idx)

--    print("For " .. i)
    selection.DeselectAll ()
    for id_idx = 1, modify_count do
--        print("    For " .. seg_idx .. select " .. ids[id_idx]  .. " distance " .. structure.GetDistance (seg_idx , ids[id_idx] ))
        selection.Select ( ids[id_idx] )
    end

    behavior.SetClashImportance ( clash_importance_during_modify )

    -- Take action
    if (action_type == 1) then
        structure.MutateSidechainsSelected ( 1 )
    elseif (action_type == 2) then
        structure.ShakeSidechainsSelected ( 1 )
    elseif (action_type == 3) then
        structure.RebuildSelected( 1 )
    end
 
    local score = GetScore ()
    table.insert ( loss_from_modify , best_score - score )

    local n_mutations = GetNumberOfPrimarySequenceChanges ()

    print (  "      lom " .. r3 ( best_score - score ) .. "  nc " .. n_mutations )

    PrintResidueChanges()

    if ( math.abs ( score - best_score ) > 0.5 ) then
        PostModify()

        score = GetScore()
        if ( score > best_score + discard_limit ) then
            save.Quicksave ( kOriginalStructureOrNewBest )
            GetBasePrimarySequence ()
            best_score = score
           print ( "Improvement to " .. r3 ( best_score ) )
       end
    end
end


function MutateAll ()
    local nth = 1
    local seg_idx

    for start_idx = 1, delta_residue do
        seg_idx = min_residue - 1 + start_idx
        
        while ( seg_idx <= max_residue ) do

            print ( get_aa3 (seg_idx) .. " " .. seg_idx .. " (" .. nth .. "/" .. max_residue - min_residue + 1 .. ")" )

            MutateResidue ( seg_idx )
            seg_idx = seg_idx + delta_residue
            nth = nth + 1
        end
    end
end


function GetParameters ()
    local dlog = dialog.CreateDialog ( "Local Mutate Shake Rebuild 3.2" )

    dlog.action_type = dialog.AddSlider ( "Action type" , 1 , 1 , 3 , 0 )  
    dlog.tp1 = dialog.AddLabel ( "1 = Mutate, 2 = Shake, 3 = Rebuild"  )

    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.delta_residue = dialog.AddSlider ( "Delta residue" , delta_residue  , 1 , kMaxDelta , 0 )  
    dlog.cidm = dialog.AddSlider ( "Clash importance" , clash_importance_during_modify , 0 , 1.0 , 2 )  
    dlog.modify_count = dialog.AddSlider ( "Modify count" , modify_count , 1 , 40 , 0 )
    dlog.discard_limit = dialog.AddSlider ( "Discard limit" , discard_limit , 0 , 25 , 0 )

    dlog.score_type = dialog.AddSlider ( "Score type" , score_type , 1 , 2 , 0 )  
    dlog.tp2 = dialog.AddLabel ( "1 = Normal : 2 = Normal for Filters"  )

    dlog.ok = dialog.AddButton ( "OK" , 1 )
    dlog.cancel = dialog.AddButton ( "Cancel" , 0 )

    if ( dialog.Show ( dlog ) > 0 ) then
    
        action_type = dlog.action_type.value
        min_residue = dlog.min_residue.value
        max_residue = dlog.max_residue.value
        delta_residue = dlog.delta_residue.value
        clash_importance_during_modify = dlog.cidm.value
        modify_count = dlog.modify_count.value
        discard_limit = dlog.discard_limit.value
        score_type = dlog.score_type.value

        return true
      else
        return false
    end
end


function PrintParameters()
    print ( "Local Mutate 3.2" )

    print ( "Action type " .. action_type )
 

    print ( "Mutate range " .. min_residue .. " to " .. max_residue )
 
    print  ( "Delta " .. delta_residue )
    print  ( "CI on modify " .. r3 ( clash_importance_during_modify ) )
    print  ( "Modify count " .. modify_count )
    print  ( "Discard limit " .. discard_limit )

    print  ( "Score type " .. score_type  )

end


function Main()

    band.DisableAll ()
    n_residues = structure.GetCount ()
 
    save.Quicksave ( kOriginalStructureOrNewBest )
    GetBasePrimarySequence ()

    original_slow_filters_setting = behavior.GetSlowFiltersDisabled ()

    if ( GetParameters () == false ) then
         return                -- graceful exit
    end
    
    max_residue = math.min (n_residues , max_residue)
    min_residue = math.max (1 , min_residue)
    if ( max_residue < min_residue ) then
        print ( "Mutate range error." )
        return
    end

    PrintParameters()
    
    print  ( "" )

    best_score = GetScore ()
    print (  "Start  score : " .. r3 ( best_score ) )

    for i = 1 , n_residues do
        original_residues [ i ] = structure.GetAminoAcid ( i )
    end

    start_time = os.time ()

    MutateAll()

    Cleanup()
end


function Cleanup ()
    print ( "Cleaning up" )
    behavior.SetClashImportance ( 1.0 )
    save.Quickload ( kOriginalStructureOrNewBest )

    selection.SelectAll ()
    band.EnableAll ()

    behavior.SetSlowFiltersDisabled ( original_slow_filters_setting )

    if ( #loss_from_modify > 1 ) then
        table.sort ( loss_from_modify )
        mid_pt = math.ceil ( #loss_from_modify / 2 )
        print ( "Median loss from modify " .. r3 ( loss_from_modify [ mid_pt ] ) )
     end

    local new_residues = {}
    for i = 1, n_residues do
        new_residues[ i ] = structure.GetAminoAcid ( i )
    end

    for i = 1 , n_residues do
        if ( new_residues[ i ] ~= original_residues [ i ] ) then
            print ( i .. " : " .. aa1_to_aa3 ( original_residues  [ i ]  ) .. " -> " .. aa1_to_aa3 ( new_residues  [ i ] ) )
        end
    end

    end_time = os.time ()
    print ( "Elapsed time : " .. os.difftime ( end_time , start_time ) .. " secs" )

end


function ErrorHandler(err)
    print("Error: " .. err)
    Cleanup()
end

--main ()
xpcall ( Main , ErrorHandler )