Recipe: Local Mutate 4.1 Neural Mutate Fixed

created by Joanna_H

Profile

• Name: Local Mutate 4.1 Neural Mutate Fixed
• ID: 108181
• Shared with: Public
• Parent: Local Mutate 4.1 Neural Mutate
• Children:
  • analyse aminozuren puzzle S. AUREUS
• Created on: January 19, 2023 at 15:51 PM UTC
• Updated on: January 19, 2023 at 16:02 PM UTC
• Description:

  Local Mutate with Neural Mutate instead of standard Mutate. After mutate in post modify it does an initial shake on the changes to ensure clashes are removed. So far seems to score better than the standard mutate function.

Code

min_residue = 1
max_residue =  999
n_residues = 0
n_mutations = 0

delta_residue = 9
kMaxDelta = 20

best_score = 0

clash_importance_during_modify = 0.09
packing_importance_during_modify = 1.0
hiding_importance_during_modify = 1.0
pairwise_importance_during_modify = 1.0
backboneHbond_importance_during_modify = 1.0
sidechainHbond_importance_during_modify = 1.0

kWiggleIterations = 15 

modify_count = 10

score_type = 1

loss_from_modify = {}

original_residues = {}
new_residues = {}

scores = {}
ids = {}

start_time = 0

wiggle_sidechains_first = false

--    Amino acids. Hydrophobics first.

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

        -- Use of quicksave slots

kOriginalStructureOrNewBest = 10
base_primary_sequence = {}
n_sequence_changes = 0

function r3 ( x )

-- Round to 3 decimal places

    t  = 10 ^ 3
    return math.floor ( x*t + 0.5 ) /  t 
end

function GetScore ()

    if ( score_type == 1 ) then

        score = current.GetScore () 

    end

    return score

end

function aa1_to_aa3 ( k )

  for  j = 1, 20 do
        if ( k == single_letter_codes [  j ] ) then
            return three_letter_codes [  j  ]
        end
    end

    return "Unk"

end

function get_aa3 ( i )

    k = structure.GetAminoAcid (  i )

    for  j = 1, 20 do
        if ( k == single_letter_codes [  j ] ) then
            return three_letter_codes [  j  ]
        end
    end

    return "Unk"

end

function GetBasePrimarySequence ()

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

end

function GetNumberOfPrimarySequenceChanges ()

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

    return n_changes
end


function ShellSort ()

-- Adapted from Numerical Recipes in C

    inc =  1
    
    repeat 
        inc = inc * 3 + 1
    until  inc > n_residues

    repeat

        inc = inc  / 3 
        inc = inc - inc % 1

        for i = inc + 1 ,  n_residues do

            v =  scores [  i  ]
            w = ids [ i ]

            j = i

            flag = false 

            while (  flag == false and scores [  j - inc ] > v ) do
                scores [  j  ] = scores [  j  - inc ]
                ids [ j ] = ids [ j - inc ]

                j = j - inc
                
                if ( j <= inc ) then
                    flag = true
                end

            end

            scores [ j ]  = v
            ids [ j ] = w
        end

    until  inc  <= 1

end


function PostModify ( i )

    behavior.SetClashImportance ( 1.0 )
    behavior.SetPackingImportance ( 1.0 )
    behavior.SetHidingImportance ( 1.0 )
    behavior.SetPairwiseImportance ( 1.0 )
    behavior.SetBackboneHBondImportance ( 1.0 )
    behavior.SetSidechainHBondImportance ( 1.0 )
 
    structure.ShakeSidechainsSelected ( 1 )

    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 MutateResidue ( i  )

    local	j
    local 	k

    save.Quickload ( kOriginalStructureOrNewBest )

     for  j = 1 , n_residues do
        scores [ j ] = structure.GetDistance ( i , j )
        ids [ j ] = j
    end 

    ShellSort ()

    selection.DeselectAll ()
    for j = 1 , modify_count do
        selection.Select ( ids [ j ] )
    end

    behavior.SetClashImportance ( clash_importance_during_modify )
    behavior.SetPackingImportance ( packing_importance_during_modify )
    behavior.SetHidingImportance ( hiding_importance_during_modify )
    behavior.SetPairwiseImportance ( pairwise_importance_during_modify )
    behavior.SetBackboneHBondImportance ( backboneHbond_importance_during_modify )
    behavior.SetSidechainHBondImportance ( sidechainHbond_importance_during_modify )
 
    structure.NeuralNetMutateSidechainsSelected ()  

    score = GetScore ()
    table.insert ( loss_from_modify , best_score - score )

    n_mutations = GetNumberOfPrimarySequenceChanges ()

 --   print ( "Hid " .. r3 ( scores [ modify_count ] ) .. " lim " .. r3 ( best_score - score ) .. "  nc " .. n_mutations )
    print (  "      lom " .. r3 ( best_score - score ) .. "  nc " .. n_mutations )


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


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

        score = GetScore ()
        if ( score > best_score ) then

            save.Quicksave ( kOriginalStructureOrNewBest )
            GetBasePrimarySequence ()
            best_score = score
           print ( "Improvement to " .. r3 ( best_score ) )
       end

    end

end

function MutateAll ()

    n = 1

    for start_idx = 1 , delta_residue do

        i = min_residue - 1 + start_idx
        
        while ( i <= max_residue ) do

            print ( get_aa3 (i) .. " " .. i .. " (" .. n .. "/" .. max_residue - min_residue + 1 .. ")" )

            MutateResidue ( i )
            i = i + delta_residue
            n = n + 1
        end

    end

end

function GetOptionalParameters ()

    local dlog = dialog.CreateDialog ( "Options" )

    dlog.pidm = dialog.AddSlider ( "Packing importance" , packing_importance_during_modify , 0 , 3.0 , 2 )  
    dlog.hidm = dialog.AddSlider ( "Hiding importance" , hiding_importance_during_modify , 0 , 3.0 , 2 )  
    dlog.widm = dialog.AddSlider ( "Pairwise importance" , pairwise_importance_during_modify , 0 , 3.0 , 2 )  
    dlog.bidm = dialog.AddSlider ( "Backbone hbond importance" , backboneHbond_importance_during_modify , 0 , 3.0 , 2 )  
    dlog.sidm = dialog.AddSlider ( "Sidechain hbond importance" , sidechainHbond_importance_during_modify , 0 , 3.0 , 2 )  

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

    if ( dialog.Show ( dlog ) > 0 ) then
    
        packing_importance_during_modify = dlog.pidm.value
        hiding_importance_during_modify = dlog.hidm.value
        pairwise_importance_during_modify = dlog.widm.value
        backboneHbond_importance_during_modify = dlog.bidm.value
        sidechainHbond_importance_during_modify = dlog.sidm.value
 
        return true
      else
        return false
    end

end

function GetParameters ()

    local dlog = dialog.CreateDialog ( "Local Mutate 4.1 Neural Mutate" )

    dlog.min_residue = dialog.AddSlider ( "Min residue" , min_residue , 1 , n_residues , 0 )  
    dlog.max_residue = dialog.AddSlider ( "Max residue" , max_residue , 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.ok = dialog.AddButton ( "OK" , 1 )
    dlog.cancel = dialog.AddButton ( "Cancel" , 0 )
    dlog.other_options = dialog.AddButton ( "Other options" , 2 )

    return_code = dialog.Show ( dlog )

    if ( return_code > 0 ) then
    
        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
 
    end
    
    return return_code

end

function main ()

    band.DisableAll ()
    n_residues = structure.GetCount ()
 
    -- Trim off locked terminal sequences as a UI convenience
    
    min_residue = 1
    while ( ( structure.IsLocked ( min_residue ) == true ) and ( min_residue <= n_residues ) ) do
      min_residue = min_residue + 1
    end
    
    max_residue = n_residues
    while ( ( structure.IsLocked ( max_residue ) == true ) and ( max_residue >= 1 ) ) do
      max_residue = max_residue - 1
    end  

    save.Quicksave ( kOriginalStructureOrNewBest )
    GetBasePrimarySequence ()

    repeat
       dialog_code = GetParameters ()
       if ( dialog_code == 2 ) then
           GetOptionalParameters ()
       end
    until ( dialog_code < 2 )

    if ( dialog_code == 0 ) then
        return
    end

    print ( "Local Mutate 4.1 Neural Mutate" )
 
    max_residue = math.min ( n_residues , max_residue )
    min_residue = math.max ( 1 , min_residue )

    print  ( "Mutate range " .. min_residue .. " to " .. max_residue )
    if ( max_residue < min_residue ) then
        print ( "Mutate range error." )
        return
    end
 
    print  ( "Delta " .. delta_residue )
    
    print ( "Clashing importance " .. r3 ( clash_importance_during_modify ) )
    print ( "Packing importance " .. r3 ( packing_importance_during_modify ) )
    print ( "Hiding importance " .. r3 ( hiding_importance_during_modify ) )
    print ( "Pairwise importance " .. r3 ( pairwise_importance_during_modify ) )
    print ( "Backbone Hbond importance " .. r3 ( backboneHbond_importance_during_modify ) )
    print ( "Sidechain Hbond importance " .. r3 ( sidechainHbond_importance_during_modify ) )
    
    print  ( "CI on modify " .. r3 ( clash_importance_during_modify ) )
    print  ( "Modify count " .. modify_count )
    print  ( "" )

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

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

    if ( ( math.abs ( packing_importance_during_modify ) - 1 > 1e-3 ) or 
         ( math.abs ( hiding_importance_during_modify ) - 1 > 1e-3 ) or 
         ( math.abs ( pairwise_importance_during_modify ) - 1 > 1e-3 ) or 
         ( math.abs ( backboneHbond_importance_during_modify ) - 1 > 1e-3 ) or 
         ( math.abs ( sidechainHbond_importance_during_modify ) - 1 > 1e-3 ) ) then 
       print ( "" )
       print ( "MAKE SURE Recipe Score Modding IS CHECKED" )
       print ( "" )
    end
    
    start_time = os.time ()

    MutateAll ()

    cleanup ()

end

function cleanup ()
   
    print ( "Cleaning up" )
   
    behavior.SetClashImportance ( 1.0 )
    behavior.SetPackingImportance ( 1.0 )
    behavior.SetHidingImportance ( 1.0 )
    behavior.SetPairwiseImportance ( 1.0 )
    behavior.SetBackboneHBondImportance ( 1.0 )
    behavior.SetSidechainHBondImportance ( 1.0 )

    save.Quickload ( kOriginalStructureOrNewBest )

    selection.SelectAll ()
    band.EnableAll ()

    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

    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

--main ()
xpcall ( main , cleanup )