Contents
Introduction
Reference Guide
ICM options
Editing
Graph.Controls
Alignment Editor
Constants
Subsets
Molecules
Selections
Regexp
Tree cluster
Arithmetics
Flow control
MolObjects
Energy Terms
Integers
Reals
Logicals
Strings
Preferences
Tables
Other
Chemical
Smiles
Soap
Gui programming
Commands
Functions
Abs
Acc
Acos
Acosh
Align
Angle
Area
Asin
Asinh
Ask
Askg
Atan
Atan2
Atanh
Atom
Augment
Axis
Bfactor
Boltzmann
Box
Bracket
Cad
Ceil
Cell
Charge
Chemical
Cluster
Color
Consensus
Corr
Cos
Cosh
Count
Date
Deletion
Det
Disgeo
Distance
Eigen
Energy
Error
Error soap
Exist
Existenv
Extension
Exp
Field
File
Find
Floor
Formula
Getarg
Getenv
Gradient
Grob
Group
Header
Histogram
Iarray
IcmSequence
Image
Index
Indexx
Insertion
Info
Info image
Integer
Integral
Interrupt
Label
Length
LinearFit
LinearModel
Link
Log
Map
Mass
Moment
Match
Matrix
Max
MaxHKL
Mean
Min
Money
Mod
Mol
Name
Namex
Next
Nof
Norm
Normalize
NotInList
Obj
Occupancy
Path
Parray
Pattern
Pi
Potential
Power
Predict
Probability
Profile
Property
Putarg
Putenv
Radius
Random
Rarray
Real
Remainder
Reference
Replace
Res
Resali
Resolution
Rfactor
Rfree
Rmsd
Rot
Sarray
Score
Select
Sequence
Shuffle
Sign
Sin
Sinh
Site
Slide
Smiles
Smooth
Sql
Sqrt
Sphere
SoapMessage
Sort
Split
Srmsd
String
Sstructure
Sum
Symgroup
Table
Tan
Tanh
Tensor
Temperature
Time
Tointeger
Tolower
Toreal
Torsion
Tostring
Toupper
Tr123
Tr321
Trace
Trans
Transform
Transpose
Trim
Trim sequence
Turn
Type
Unique
Unix
Value
Value soap
Vector
Version
Volume
View
Warning
Xyz
Macros
Files
User's Guide
References
Glossary
Index

[ Abs | Acc | Acos | Acosh | Align | Angle | Area | Asin | Asinh | Ask | Askg | Atan | Atan2 | Atanh | Atom | Augment | Axis | Bfactor | Boltzmann | Box | Bracket | Cad | Ceil | Cell | Charge | Chemical | Cluster | Color | Consensus | Corr | Cos | Cosh | Count | Date | Deletion | Det | Disgeo | Distance | Eigen | Energy | Error | Error soap | Exist | Existenv | Extension | Exp | Field | File | Find | Floor | Formula | Getarg | Getenv | Gradient | Grob | Group | Header | Histogram | Iarray | IcmSequence | Image | Index | Indexx | Insertion | Info | Info image | Integer | Integral | Interrupt | Label | Length | LinearFit | LinearModel | Link | Log | Map | Mass | Moment | Match | Matrix | Max | MaxHKL | Mean | Min | Money | Mod | Mol | Name | Namex | Next | Nof | Norm | Normalize | NotInList | Obj | Occupancy | Path | Parray | Pattern | Pi | Potential | Power | Predict | Probability | Profile | Property | Putarg | Putenv | Radius | Random | Rarray | Real | Remainder | Reference | Replace | Res | Resali | Resolution | Rfactor | Rfree | Rmsd | Rot | Sarray | Score | Select | Sequence | Shuffle | Sign | Sin | Sinh | Site | Slide | Smiles | Smooth | Sql | Sqrt | Sphere | SoapMessage | Sort | Split | Srmsd | String | Sstructure | Sum | Symgroup | Table | Tan | Tanh | Tensor | Temperature | Time | Tointeger | Tolower | Toreal | Torsion | Tostring | Toupper | Tr123 | Tr321 | Trace | Trans | Transform | Transpose | Trim | Trim sequence | Turn | Type | Unique | Unix | Value | Value soap | Vector | Version | Volume | View | Warning | Xyz ]

 
 show Mean(Random(1.,3.,10))  

Abs ( map ) - returns map of absolute values of the source map.

Examples:

 
 a=Abs(-5.)                 # a=5.  
 print Abs({-2.,0.1,-3.})   # prints rarray {2., 0.1, 3.}  
 if (Abs({-3, 1}=={3 1}) print "ok" 

 
                                   # let us select interface residues  
 read object s_icmhome+"complex" 
                                   # display all surface residues 
 show surface area 
 display Acc( a_/* ) 
                                   # now let us show the interface residues 
 display a_1,2 
 color a_1 yellow 
 color a_2 blue 
 show surface area a_1 a_1         # calculate surface of  
                                   # the first molecule only  
 
                                   # select interface residues  
                                   # of the first molecule  
 color red Sphere(a_2/* a_1/* 4.) & Acc(a_1/*)      
 
 read object "crn" 
 show energy "sf" 
 display 
 display cpk Acc(a_//* 0.1)   # display accessible atoms  
 
 show surface area            # prior to invoking Acc function  
                              # surface area should be calculated  
 color Acc(a_/*) red          # color residues with relative  
                              # accessibility > 25% red  

 
 print Acos(1.)               # equal to 0.  
 print Acos(1)                # the same  
 
 print Acos({-1., 0., 1.})    # returns {180. 90. 0.}  

 
 print Acosh(1.)               # returns 0  
 print Acosh(1)                # the same  
 
 print Acosh({1., 10., 100.})  # returns {0., 2.993223, 5.298292}  

[ sequence | structural alignment | sub_alignments ]

• sequence alignment (with the Needleman and Wunsch algorithm with zero gap end penalties ( ZEGA ),
• structural alignment, or
• subalignment extraction

Pairwise sequence alignment or sequence-structure alignment

Align ( [ sequence1, sequence2 ] [ area ] [ M_scores ] )

- returns ZEGA- alignment. If no arguments are given, the function aligns the first two sequences in the sequence list. For sequence alignments, the ZEGA-statistics of structural significance ( Abagyan, Batalov, 1997) is given and can be additionally evaluated with the Probability function. The reported pP value is -Log(Probability,10).

Returned variables:

• i_out - the number of identical residues in the alignment
• r_out - contains Log( Probability_of_structural_dissimilarity ) only for pairwise alignments
• r_2out - percent identity of the alignment.

 
   read sequences s_icmhome+"sh3.seq" # read 3 sequences 
   print Align(Fyn,Spec)              # align two of them  
   Align( )                           # the first two 
   a=Align( sequence[1] sequence[3] ) # 1st and 3rd 
   if(r_out > 5.) print "Sequences are struct. related"  

 
 read pdb "1lbd" 
 show surface area  
 make sequence 
  Info> sequence  1lbd_m  extracted  
 1lbd_m                       # see the relative areas 
 read pdb sequence "1fm6.a/"  # does not have areas 
 Info> 1 sequence 1fm6_a  read from /data/pdb/fm/pdb1fm6.ent.Z    
 ali3d = Align( 1lbd_m 1fm6_a area )  

 
 set area seq1 R_weights  # must be > 0. and less than 2.37 
 Align( seq1 seq2 area ) 

 
 seq1 = Sequence("WEARSLTTGETGYIPSA") 
 seq2 = Sequence("WKVEVNDRQGFVPAAY") 
 Align() 
 # Consensus   W.#.  .~~.~G%#P^ 
           WEARSLTTGETGYIPS-- 
           WKVE--VNDRQGFVPAAY 
 m = Matrix(17,16,0.) 
 m[10,4] = 3. # reward alignment of E in seq1[10]  and E in seq2[4] 
 Align(seq1 seq2 m ) 
 # Consensus   W.#      E    ~G%#P^ 
           WEARSLTTGE----TGYIPS-- 
           WKV------EVNDRQGFVPAAY 

See also: Probability( ali .. ) for local alignment reliability.

Local pairwise structural alignment

[ alignment-from-tethers ]

Two types of structural alignments or mixed sequence/structural alignments can be performed with the Align function.
Align( seq_1 seq_2 distance [ i_window ] [ r_seq_weight ] ) - performs local structural alignment, using distance RMSD as structural fitness criterion. The RMSD is calculated in a window i_window and the dynamic programming algorithm combines structural scores with sequence alignment scores if r_seq_weight > 0.,
Align( seq_1 seq_2 superimpose [ i_window ] [ r_seq_weight ] ) - performs local structural alignment, using superposition followed by coordinate RMSD calculation as structural fitness criterion. The RMSD is calculated in a window i_window and the dynamic programming algorithm combines structural scores with sequence alignment scores if r_seq_weight > 0.,
In both cases the function uses the dynamic algorithm to find the alignment of the locally structurally similar backbone conformations.
The alignment based on optimal structural superposition of two 3D structures may be different from purely sequence alignment
Preconditions:

• sequences must be linked to 3D molecules to access the coordinate information;
• two 3D structures must have a superimposable subsets

See also: align ms1 ms2 function

Deriving an alignment from tethers between two 3D objects
Align ( ms_ ) - returns alignment between sequences of the specified molecule and the template molecule to which it is tethered. The alignment is deduced from the tethers imposed.

Example:

 
 build string "se ala his leu gly trp ala" "a"  # obj. a  
 build string "se his val gly trp gly ala" "b"  # obj. b  
 set tether a_2./1:3 a_1./2:4 align                 # impose tethers  
 show Align(a_2.1)         # derive alignment from tethers  
 write Align(a_2.1) "aa"   # save it to a file  

 
             # 14 sequences  
 read alignment msf s_icmhome + "azurins"  
             # extract a pairwise alignment by names  
 aa = Align(azurins,Azu2_Metj,Azur_Alcde)  
             # reordered sub-alignment extracted by numbers  
 bb = Align(azurins,{2 5 3 4 10 11 12})    

 
 macro reorderAlignmentSeq( ali_ ) 
  nn=Name(ali_)  # names in the alignment order 
  ii=Iarray(Nof(nn)) 
  j=0 
  for i=1,Nof(sequence)  # the original order 
   ipos = Index( nn, Name(sequence[i] ) ) 
   if ipos >0 then 
     j=j+1 
     ii[j] = ipos 
   endif 
  endfor 
  ali_new = Align( ali_ ii ) 
  keep ali_new 
 endmacro 

 
 print Angle(  # and then click the atom of interest. 

 
 d=Angle( a_/4/c )                         # d equals N-Ca-C angle  
 print Angle( a_/4/ca a_/5/ca a_/6/ca )    # virtual Ca-Ca-Ca planar angle  

The rotation angle corresponding to a transformation vector is returned as r_out by the Axis( R_12 ) function.

[ Area contact matrix ]

Area( grob ) ⇒ r

Area( as_ | rs_ ) ⇒ R_atomAreas|R_resAreas # needs surface calculation beforehand

Area( rs_ type ) ⇒ R_maxAreas_in_GLY_X_GLY

Area( as_ R_typeEyPerArea energy ) ⇒ R_atomEnergies

Area( seq ) ⇒ R_relAreasPerResidue

Area( s_icmResType ) ⇒ r

Area( rs_ rs_2 ) ⇒ M_contactAreas

Note that if an atom selection is provided as an argument the surface area needs to be computed beforehand with the show area or show energy "sf" command. The detailed description can be found below:

Area ( grob ) - returns real surface area of a solid graphics object.
See also: the Volume( grob ) function, the split command and How to display and characterize protein cavities section.
Area ( as_ [ [ R_userSolvationDensities ] [ energy ] ] ) - returns rarray of pre-calculated solvent accessible areas or energies for selected atoms `as_ . This areas are set by the show area surface|skin of show energy "sf" commands. Make sure to clean up the areas with the set area a_//* 0. command before computing the areas with show energy command since the command ignores hydrogens.

With option energy returns the product of the individual atomic accessibilities by the atomic surface energy density. The values of the density depend on the surfaceMethod preference and are stored in the icm.hdt file. The "contant tension" value of the preference is a trivial case in which all areas are multiplied by the surfaceTension parameter. For the "atomic solvation" and "apolar" styles, the densities depend on atom types. Normally the atomic solvation densities are taken from the icm.hdt file where the density values are listed for each hydration atom type for "atomic solvation" and "apolar" styles. However, you can provide your own array of n values R_userSolvationDensities with the number of elements less or equal to the number of types to overwrite the first n types.

Examples:

set area a_//* 0.
surfaceMethod = "apolar"
show energy "sf" # only heavy atoms
Area( a_/15:30/* ) # areas of this atoms
#
# Now let us redefine the first three solvation parameters
# of icm.hdt and calculated E*A contributions of selected atoms
#
Area( a_/15:30/* {10., 20. 30.} energy) 

 
read object "1crn" 
show surface area 
a=Area(a_/* )      # absolute conformation dependent residue accessilities  
b=Area(a_/* type ) # maximal residue accessilities in the extended conformation 
c = a/b            # relative (normalized) accessibilities 

Area ( resCode ) ⇒ r_standard_area

- returns the real value of solvent accessible area for the specified residue type in the standard "exposed" conformation surrounded by the Gly residues, e.g. Area("ala"). It is the same value as the Area( .. type ) function.
Area( seq ) ⇒ R_relAreasPerResidue

- returns an array of relative areas per residue stored with the sequence by the make sequence command from molecules in which the areas had been computed beforehand. Note that the sequence keeps only a very limited accuracy areas. Example:

read pdb "1crn"
show area surface 
make sequence  # 1crn_a now has relative areas
group table Sarray( a_/* residue) Area(1crn_a)  Area(a_/*)/Area(a_/* type) name="t"
show t

Important : "pre-calculated" above means that before invoking this function, you should calculate the surface by show area surface , show area skin or show energy "sf" commands.
Examples:

 
 build                   # build a molecule according to the sequence 
                         # from file def.se (default) 
 show area surface       # calculate surface area  
 a = Area(a_//o*)        # individual accessibilities of oxygens  
 
 stdarea = Area("lys")   # standard accessibility of lysine  
 
# More curious example 
 read object "crn"           
 show energy "sf"        # calculate the surface energy contribution 
                         # (hence, the accessibilities are 
                         # also calculated) 
 
 assign sstructure a_/* "_"     
                         # remove current secondary structure assignment 
                         # for tube representation 
 display ribbon 
                         # calculate smoothed relative accessibilities  
                         # and color tube representation accordingly  
 color ribbon a_/* Smooth(Area(a_/*)/Area(a_/* type) 5) 
                         # plot residue accessibility profile 
 plot Count(1 Nof(a_/*)) Smooth(Area(a_/*)/Area(a_/* type) 5) display 

Area contact matrix

Area ( rs_1 rs_2 ) - returns rarray of areas of contact between selected residues. You can do it for intramolecular residue contacts, in which case both selections should be the same, i.e. Area(a_1/* a_1/*) ; or, alternatively, you can analyze intermolecular residue contacts, for example, Area(a_1/A a_2/A). See also the Cad function, and example in plot area in which a contact matrix is calculated via interatomic Ca-Ca distances. The table of the pairwise contact area differences is written to the s_out string which can later be read into a proper table via: read column group name="aa" input=s_out and sorted by the area (see below).
Example:

 
 read object "crn"  # good old crambin  
 s=String(Sequence(a_/A)) 
 PLOT.rainbowStyle="blue/rainbow/red" 
 plot area Area(a_/A, a_/A) comment=s//s color={-50.,50.} \ 
    link transparent={0., 2.} ds 
 
 read object "complex" 
 plot area Area(a_1/A, a_2/A) grid color={-50.,50.} \ 
    link transparent={0., 2.} ds 

 
 print Asin(1.)               # equal to 90 degrees  
 print Asin(1)                # the same  
 
 print Asin({-1., 0., 1.})    # returns {-90., 0., 90.}  

 
 print Asinh(1.)               # returns 0.881374  
 print Asinh(1)                # the same 
 
 print Asinh({-1., 0., 1.})    # returns {-0.881374, 0., 0.881374} 

 
 windowSize=Ask("Enter window size",windowSize) 
 s_mask=Ask("Enter alignment mask","xxx----xxx") 
 
 grobName=Ask("Enter grob name","xxx") 
 display $grobName 
 
 show Ask("Enter string, it will be interpreted by ICM:", "")  
          #e.g. Consensus( myAlignm ) 
 
 show Ask("Enter string:", "As Is",simple)     
          #your input taken directly as a string 

See also: Askg

interactive input function that generates a GUI dialog. Return entered text Askg( s_prompt, i_default ) ⇒ s_returnsTheInputString

E.g.

Askg( "Enter your name", ""  ) # empty default 
Askg( "Enter your name", "Michael"  )

Return the pressed button. Askg( s_Question, "Reply1/Reply2/.." simple ) ⇒ s_theReply

Makes a GUI dialog with the question and several alternatives separated by a slash. This dialog returns one of the string selected ,e.g. "Yes", "No" , or "Cancel" for the "Yes/No/Cancel" argument. Example:

s = Askg("Do you like bananas?","Yes/No/Fried only",simple)
if s=="Fried only" print "Impressive" 

Creating a special chemical dialog for library enumeration.This one is very specialized and is used in combi-chem generator.

Askg( chem_scaffold , enumerate ) ⇒ s_makeLib_React_Args

Askg( chem_reaction , enumerate ) ⇒ s_makeLib_React_Args

prompts for arguments for the enumerate library or make reaction commands to create a combinatorial library. To use this function you need to have the chemical array objects with Markush-scaffolds or reactions, plus the building blocks loaded into ICM. The function returns a string with the agruments for the enumerate library or make reaction commands. E.g.

args = Askg( scaff1 enumerate )
enumerate library scaff1 $args

Askg( s_dialogDeclaration ) ⇒ "yes"/"no"

Generates a dialog from GUI dialog description text. Values from each input field can be accessed either by :

$field_num

or

Getarg( i_field_num gui )

buf  = "#dialog{\"Select InSilco Models\"}\n"
buf += "#1 l_Passive_GUT_Absorption (yes)\n"
buf += "#2 l_ToxCheck (no)\n"
buf += "#3 l_hERG_QSAR (yes)\n"
buf += "#4 s_Comment_Here ()\n"

Askg(buf)

print $1, $2, Getarg( 3 gui ), $4

Using Askg in shell, html-docs and table tool panels. These variants of the Askg function can also be used as a part of an ICM script in dialogs generated from built-in html documents, or in actions associated with tables.

See also : gui programming

Atan

arctangent trigonometric function Returned values are in degrees.
Atan ( real | integer )
- returns the real arctangent of its real or integer argument.
Atan ( rarray )
- returns the rarray of arctangents of rarray elements.
Examples:

 
 print Atan(1.)                    # equal to 45.  
 print Atan(1)                     # the same.  
 
 print Atan({-1., 0., 1.})         # returns {-45., 0., 45.}  

 
 print Atan2(1.,-1.)                        # equal to 135.  
 print Atan2({-1., 0., 1.},{-0.3, 1., 0.3}) # returns phases {-106.7 0. 73.3}  

 
 print Atanh(0.)                    # returns 0.  
 print Atanh(1.)                    # returns error  
 
 print Atanh({-0.9999, 0., .9999})  # returns { -4.951719, 0., -4.951719 }  

 
 asel=Acc(a_2/his)             # select accessible His residues of  
                               # the second molecule  
 show Atom(asel)               # show atoms of these residues 
 show Atom( v_//phi )          # carbonyl Cs 

 
 a1  a2  a3  | a4 
 a5  a6  a7  | a8 
 a9  a10 a11 | a12 
 ------------+---- 
 0.  0.  0.  | 1. 

 
 M_inv = Power(Augment(R_12direct),-1)) 
 R_12inv = Vector(M_inv) 

 
 display a__crn.        # load and display crambin: P21 group  
 obl = Augment(Cell( ))  # extract oblique matrix  
 A = obl[1:3,1]         # vectors A, B, C 
 B = obl[1:3,2] 
 C = obl[1:3,3] 
 g1=Grob("cell",Cell( )) # first cell 
 g2=g1+ (-A)            # second cell 
 display g1 g2          

• r_out rotation angle (in degrees);
• r_2out helix rise;
• R_out 3-rarray with a point on the axis.

 
 color ribbon a_/ Trim(Bfactor( a_/ simple ),-0.5,3.)//-0.5//3. # or 
 color a_// Trim(Bfactor( a_// simple ),-0.5,3.)//-0.5//3.  # for atoms 

 
 avB=Min(Bfactor(a_//ca))     # minimal B-factor of Ca-atoms  
 show Bfactor(a_//!h*)        # array of B-factors of heavy atoms  
 color a_//* Bfactor(a_//*)   # color previously displayed atoms  
                              # according to their B-factor  
 color ribbon a_/A Bfactor(a_/A) # color the whole residue by mean B-fac. 

 
 deltaE = Boltzmann*temperature  # energy  

 
 build string "se ala his"  # a peptide  
 display box Box(a_/2 1.2)  # surround the a_/2 by a box with 1.2A margin 
 color a_//* & Box( )  

 
 make map potential "gh,gc,gb,ge,gs" a_1 Box() 
 m_ge = Bracket(m_ge, Box( a_1/15:18,33:47 )) # redefine m_ge 

[ Cad1 | Cadalign ]

• between two conformations of the same protein based on full atom residue-residue contact area calculation, Cad(..)
• between two conformations of the same protein based on Cbeta-Cbeta distance evaluation (`Cad1{Cad}(.. distance ) .ICM uses an empirically derived ContactStrength( Cb-distance ) function.
• between two homologous structures based preservation of the residue contacts through the alignment ( Cad (.. alignment )) . The contact strength in this case is also derived from the interresidue distances.

• by a single selection rs_A1 : all pairs between selected residues (is equivalent to rs_A1 rs_A1 )
• by two residue selections rs_A1 rs_A2: cross pairs between two sets of selected residues (e.g. the contacts between two subunits)

• distance option allows to compare approximations of the inter-residue contact areas by the Ca and Cbeta positions. This allows to calculated deformations between two homologous proteins which is not possible in the default mode in which two chemically identical molecules are compared. The residue pairs in two homologs are equivalenced according to the alignments linked to the molecules. Residues deleted in a homologue are considered to have zero contact.

 
# Ab initio structure prediction, Overall models by homology  
   read pdb "cnf1"   # one conformation of a protein  
   read pdb "cnf2"   # another conformation of the same protein 
   show 1.8*Cad(a_1. a_2.)   # CAD=0. - identical; =100. different 
   show 1.8*Cad(a_1.1 a_2.1) # CAD between the 1st molecules (domains) 
   show 1.8*Cad(a_1.1/2:10 a_2.1/2:10) # CAD in a window 
   PLOT.rainbowStyle = 2 
   plot area grid M_out comment=String(Sequence(a_1,2.1)) link display 
 
# Loop prediction: 0% - identical; ~100% totally different 
#    CAD for loop 10:20 and its interactions with the environment 
   show 1.8*Cad(a_1.1/10:20 a_1.1/* a_2.1/10:20 a_2.1/*)   
#    CAD for loop 10:20 itself 
   show 1.8*Cad(a_1.1/10:20 a_1.1/10:20 a_2.1/10:20 a_2.1/10:20)   
 
# Evaluation of docking solutions: 0% - identical; 100% totally different 
   read pdb "expr"   # one conformation of a complex  
   read pdb "pred"   # another conformation of the same complex 
   show Cad(a_1.1 a_1.2 a_2.1 a_2.2)  # CAD between two docking solutions  
# 
# ANOTHER EXAMPLE: the most changed contacts  
   read object "crn" 
   copy a_ "crn2" 
   randomize v_ 5. 
   Cad(a_1. a_2.) 
   show s_out 
   read column group input= s_out name="cont" 
   sort cont.1 
   show cont 
# the table looks like this (the diffs can be both + and -): 
#>T cont 
#>-1-----------2-----------3---------- 
   -39.        a_crn.m/38  a_crn.m/1 
   -36.4       a_crn.m/46  a_crn.m/4 
   -32.1       a_crn.m/46  a_crn.m/5 
   -29.8       a_crn.m/30  a_crn.m/9 
   -25.2       a_crn.m/37  a_crn.m/1 
... 
   42.5        a_crn.m/43  a_crn.m/5 
   45.1        a_crn.m/44  a_crn.m/6 
   45.2        a_crn.m/43  a_crn.m/6 
   55.3        a_crn.m/46  a_crn.m/7 
   56.         a_crn.m/45  a_crn.m/7 

• formal : return formal charges
• mmff : return formal charges calculated according to mmff atom types and rules. Note: do not confuse this option with a function to return the mmff charges.

 
 build string "ala his glu lys arg asp" 
 show Charge(a_1)      # charge per molecule 
 show Charge( a_1/* )  # charge per residue 
 show Charge( a_1//* ) # charge per atom  
 
 avC=Charge(a_/15)        # total electric charge of 15th residue  
 avC=Sum(Charge(a_/15/*)) # another way to calculate it 
 show Charge(a_//o*)      # array of oxygen charges  
 
# to return mmff charges: 
 set type mmff 
 set charge mmff 
 Charge( a_//* ) 
 
# to return total charges per molecular object: 
 read mol s_icmhome+"ex_mol.mol" 
 set type mmff 
 set charge mmff 
 Charge( a_*. ) 

Converting 3D objects to chemical arrays.

Chemical( ms|os [exact] [hydrogen] [unique] [pharmacophore] )

returns an array of chemicals from a molecular selection of 3D molecular objects, e.g. a_H for hetero-molecules By default the selected molecules will be converted to 2D graphs. However with the exact option the original 3D coordinates will be retained in the elements of the chemical array. If you want to preserve explicitly drawn hydrogens hydrogen option should be used. Note that the number of chemicals in the array will be determined by the selection level. At the object level multiple molecules of the same object will be merged into one array element. With unique option duplicates will be excluded from the result.

Example:

read pdb "1ch8"
group table t_2D Chemical(a_H)   # convert to 2D chemical table 
group table t_3D Chemical(a_H exact)  # convert to 3D chemical table 
group table t_3D_hyd Chemical(a_H exact hydrogen) # convert to 3D chemical table and preserve hydrogens

With pharmacophore option the function generates pharmacophore points for the input selection.

Example:

read object s_icmhome + "biotin.ob" name="biotin"
read mol input = String(  Chemical(a_ pharmacophore )) name="biotin_ph4"
display xstick 
display wire a_biotin.

To display supported pharmacophore types and use show pharmacophore type command

Converting smiles to chemical arrays:

Chemical( S_smiles|s_smiles )

returns an array of chemicals from a string arrays of smiles.

Example: <> add column t Chemical({"N[C@@](F)(C)C(=O)O", "C[C@H]1CCCCO1"}) <>

See also: chemical functions

Cluster

Cluster( I_NxM_NearestNeighb i_M_totalNofNearNeighbors i_minNofCommonNeighbors ) ⇒ I_N_clusterNumbers
function returns iarray of cluster numbers for each or N points.
The input to the first function is an array of M nearest neighbors (defined by the second argument i_M_totalNofNearNeighbors) for each of N points. For example for an array for 5 points, and i_M_totalNofNearNeighbors = 3 it can be an array like this: {3,4,5, 1,3,4 1,2,5 2,3,5 1,2,3} . The points will be grouped into the same cluster if the number of neighbors they share is larger or equal than i_minNofCommonNeighbors . This clustering algorithm is adaptive to the cluster density and does not depend on absolute distance threshold. In other words it will identify both very sparse clusters and very dense ones. The nearest neighbor array can be calculated by the with the Link ( I_bitkeys , nBits, nNearestNeighbors ) function.
Cluster( M_NxNdist r_maxDist ) ⇒ I_N_clusterNumbers
This function identifies the i_totalNofNeighbors nearest neighbors from the full distance matrix M_NxNdist for each point and assembles points sharing the specified number of common neighbors in clusters.
All singlets (a single item not in any cluster) are placed in a special cluster number 0 . Other items are assigned to a cluster starting from 1.
Example with a distance matrix:

 
# let us make a distance matrix D 
# we will cook it from 5 vectors {0. 0. 0.} 
 m=Matrix(5,3)  # initialize 5 vectors 
 m[2,1:3]={1. 0. 0.}   # v2 
 m[3,1:3]={1. 1. 0.}   # v3 
 m[4,1:3]={1. 1. 1.}   # v4 
 m[5,1:3]={1. 0.1 0.1} # v5 close to v2 
 
 D = Distance( m )  # 5x5 distance matrix created 
 
 Cluster( D , 0.2 ) # v2 and v5 are assigned to cluster 1 
 Cluster( D , 0.1 ) # radius too small. All items are singlets 
 Cluster( D , 2.  ) # radius too large. All items are in cluster 1 

[ Color from gradient | Color image ]

 
 build string "se his" 
 display xstick 
 make grob image name="g_" 
 display g_ only smooth 
 M_clr = Color( g_ ) 
 for i=1,20    # shineStyle = "color" makes it disappear completely 
   color g_ (1.-i/20.)*M_clr 
 endfor 
 color g_ M_clr  

- returns sarray of color names approximating the rgb values in the matrix. The color names and definitions are taken from the icm.clr file. Example:

m = Matrix(3)
Color(m) # returns {"red","blue1","green"}

Color( background )
- returns rarray of three RGB components of the background color.

Interpolating colors by gradient

Color( r_value s_gradient [ r_from r_to ] ) - returns 3-element rarray with RGB components describing the color.

Color( R_N_values s_gradient [ r_from r_to ] ) - returns matrix with N rows and 3 columns where each row is the RGB representation of the interpolated color for the respective value in the R_N_values array.

Examples:

s = "red/lime/blue"
Color( 0.  s 0. 1. ) 
Color( 0.5 s 0. 1. )
Color( 1.0 s 0. 1. )
Color( 0.1 s 0. 1. )
Color( 0.8 s 0. 1. )
Color( {0.1 0.8} s 0. 1. )
Color( {1. 8.} s 0. 10. )

Color( 0.1 "red/lime/blue,0:1" )
Color( {0.1 0.8} "red/lime/blue,0:1" )
Color( {1. 8.} "red/lime/blue,0:10" )

Image color functions

Color( imageArray_ background )

returns sarray with background colors of the images in imageArray_.The color of the top left pixel of the image is returned as the background color currently.

See also: Image, image parray

Consensus

Consensus ( ali_ ) ⇒ s_consensus
- returns the string consensus of alignment ali_. The consensus characters are these: # hydrophobic; + RK; - DE; ^ ASGS; % FYW; ~ polar. In the selections by consensus a letter code (h,o,n,s,p,a) is used.
Consensus ( ali_ { i_seq | seq_ } )
- returns the string consensus of alignment ali_ as projected to the sequence.
Sequence can be specified by its order number in the alignment or by name.
Example displaying conserved residues:

 
 read alignment "sx"  # load alignment  
 read pdb "x"         # structure  
 display ribbon 
        # multiply rs_ by a mask like "  A C   N  .." 
 cnrv = a_/A & Replace(Consensus(sx cd59),"[.^~#]"," ")  
 display cnrv red  
 display residue label cnrv 

Note: this function returns R , not R² . Taking it to the 2nd power can be a humbling experience.
Examples:

 
 r=Corr(a,b)                 # two vectors a and b  
 if (Abs(r_out) < 0.3) print "it is actually as good as no correlation" 

 
 show Cos(60.)                           # returns 0.5  
 show Cos(60)                            # the same  
 
 rho={3.2 1.4 2.3}                       # structure factors  
 phi={60. 30. 180.}                      # phases  
 show rho phi rho*Cos(phi) rho*Sin(phi)  # show in columns rho, phi,  
                                         # Re, Im  

 
 show Cosh(1.)                           # 1.543081  
 show Cosh(1)                            # the same  
 
 show Cosh({-1., 0., 1.})                # returns {1.543081, 1., 1.543081}  

 
 show Count(-2,1)         # returns {-2,-1,0,1}  
 show Count(4)            # returns {1,2,3,4}  

Date ( os_ all ) ⇒ I_dates

returns the date of the pdb file creation in an integer (or iarray ) format. The date read from the HEADER record of a pdb file and is stored with the object. The format of the integer date (e.g. 19920115 for Jan 15, 1992) enables convinient numerical comparison and sorting. Example:

read pdb "1crn"
if Date(a_) > 19820000 print "released after 1992"

Date ( [ s_format ] )

returns string with the current week day, month, year and time in the given format. The default format is "%m %DD %Y" which returns something like "Oct 08 2004" .
The allowed format specifications are the following:

format	description	example
%D	day	1 or 12 or 30
%DD	day with compulsary two digits	01 03 12 etc.
%M	month	1
%MM	see above at %DD	01
%m	brief month (3 letter)	dec oct
%mm	full month	december
%Y	four letter year	2004
%YY	two letter year	04
%y	two letter year (damn with Y3K)	03
%W	week day	Saturday
%w	three letter day	sat
%z	time zone
%HH	always aligned by two chars	03
%H	one or two chars	3
%hh	the same with in pm/am style,	03pm
%h	the same in pm/am style	3pm
%pm	the string "pm" or "am", empty otherwise.	pm
%UU	minutes	01
%U	minutes	1, 12
%SS	seconds	04
%S	seconds	1, 4,12

Example:

 
Date()  # uses the default format 
Date("%m %DD,%Y %h%pm") # today returns 
 Oct 08,2002 5pm 

• "all" (default: no string option) select deletions of all types
• "nter" select only N-terminal fragments
• "cter" select only C-terminal fragments
• "loop" select only the internal zones of deleted loops

 
 a=Rot({0. 0. 1.}, 30.)     # Z-rotation matrix by 30 degrees  
 print Det(a)               # naturally, it is equal to 1.  

 
 group sequence se1 se2 se2 se4 mySeqs 
 align mySeqs 
 distMatr=Distances(mySeqs) 

 
 read sequences s_icmhome+"zincFing"   # read sequences from the file, 
 list sequences              # see them, then ...  
 group sequence alZnFing     # group them, then ...  
 align alZnFing              # align them, then ...  
 a=Distance(alZnFing)        # a matrix of pairwise distances  
 n=Nof(a)                    # number of points  
 b=Disgeo(a)                 # calculate principal components  
 corMat=b[1:n,1:n-1]         # coordinate matrix [n,n-1] of n points  
 eigenV=b[1:n,n]             # vector with n sorted eigenvalues  
 xplot= corMat[1:n,1] 
 yplot= corMat[1:n,2] 
 plot xplot yplot CIRCLE display # call plot a 2D distribution  

[ Distance iarray | Distance rarray | Distance as_ | Distance as_ rarray | Distance matrix | Distance Tanimoto | Distance 2 matrixes | Distance tether | Distance Dayhoff | Distance in alignment | Distance 2 alignments | Distance tree | Distance chemical ]

• Display your molecule
• type Dist , press TAB
• Ctrl-RightMB click on the atom you want (or double click for a residue) and press RETURN

Distance( I_keys1 I_keys2 nBits | R_nBitWeights [simple] ) ⇒ M_distances
- returns the matrix of Tanimoto distances between two arrays of bit-strings. Each array of N-strings is represented by an iarray I_keys of N*( nBits/32 ) elements (e.g. if nBits is 32 , each integer represents 1 bit-string, if nBits i 64, I_keys1 has two integers for each bit string, etc.). The returned matrix dimensions are N1 x N2 . The distance is defined as 1. - similarity , where The Tanimoto similarity between bitstrings is defined as follows: The number of the on-bits in-common between two strings divided by the number of the on-bits in either bit-string.
You can provide a relative weight for each bit in a bit-string as a rarray R_weights. In this case the weighted Tanimoto distance is calculated as follows:

 
  distWeighted = 1. - Sum( Wi_of_common_On_Bits ) / Sum( Wi_of_On_Bits ) 

 
Distance({3} {1} 32 simple ) # returns 0. 
Distance({1} {3} 32 simple ) # returns 0.5 

 
Distance({1 2 3},{1 2 3},32) 
 #>M 
 0. 1. 0.5 
 1. 0. 0.5 
 0.5 0.5 0. 

See also:

• Iarray-bits-to-integers{ Iarray({1 0 0 1 1 0 ..} key ) } to generate compressed integer bit vectors

read pdb "1crn" convert tether # keeps tethers to the pdb original deviations = Distance( a_//!h*,vt* , 9.9) perResDevs = Group( deviations, a_//!h*,vt* ,"max") # find max.devs per residue display ribbon color ribbon a_/* perResDevs # Another example Distance( a_//T ) # selects only tethered atoms #>R 1.677 1.493 1.386 1.435 1.645 1.570 2.165 1.399

Distance( seqArr[n]> <seq ) ⇒ R_n

- returns the real measure of similarity between two aligned sequences. Zero distance means 100% identity. The distance is calculated by the following two steps:

1. d1 = 1.0 - (nResidueIdentities/Min(Length(Seq1), Length(Seq2)) (d1 belongs to [0.,1.] range)
2. if there is no identity option the distance is corrected: Distance(Seq1,Seq2) = DayhoffTransformation( d1 )

Transformation practically does not change small distances d1, whereas large distances, especially above 0.9 (10% sequence identity) are increased to take occasional reversals into account. Distances d1 within [0.9,1.0] are transformed to [5.17, 10.] range.

Distance between sequences or alignment sequences

Distance ( alignment )

Distance( seqArr[n]> ) ⇒ <M_nn

Distance( seqArr[n]> <seqArr[m]> ) ⇒ <M_nm

- returns matrix of pairwise sequence-sequence distances in the alignment. These distances are calculated with the fast option as follows

 
 1.-(nResidueIdentities-gapPenalty)/Min(Length(Seq1), Length(Seq2)) 

 
 read alignment msf "azurins"            # read azurins.msf  
 NormCoord = Disgeo(Distance(azurins))   # 2D sequence diversity in  

- returns the real distance between two alignments formed by the same sequences.
The distance is defined as a number of non-gap columns identical between two alignments.
Two different normalizations are available:
The default normalization is to the shorter alignment. ( Distance ( ali_1 ali_2 ) ). In this case the number of equivalent pairs is calculated and is divided by the total number of aligned pairs in the shorter alignment. This method detects alignment shifts but does not penalize un-alignment of previously aligned residue pairs. D = (La_min - N_commonPairs)/La_min In the following alignment the residue pairs which are aligned in both alignments are the same, therefore the distance is 0.

 
 show a1   # La1 = 3 
 ABC---XYZ 
 ABCDEF--- 
 show a2   # La2 = 6 
 ABCXYZ 
 ABCDEF 
 Distance(a1,a2)   # a1 is a subalignment of a2, distance is 0. 
 0. 

 
 Distance(a1,a2,exact)  # returns 0.5 for the above a1 and a2 

 
 read sequence msf "azurins.msf"  
 gapOpen =2.2 
 a=Align(Azu2_Metj  Azup_Alcfa)   # the first alignment 
 gapOpen =1.9                     # smaller gap penalty and .. 
 b=Align(Azu2_Metj  Azup_Alcfa)   # the alignment changes 
 show 100*Distance(a b )          # 20% difference 
 show 100*Distance(a b exact )    # 21.7% difference 
 show a b   

The distance of the cluster splitting level

Distance( treeArr i_at separator )

- return the current value of the cluster splitting level set by split command.

Chemical similarity distance

Distance( chemarray )

- return square matrix of chemical distacnes.

Example:

Distance( Chemical( { "CCC", "CCO"} ) )

Distance( chemarray1 chemarray2 )

- return a MxN matrix where M is number of elements in chemarray1 and N is number of elements in chemarray2

Example:

Distance( Chemical({ "CCC", "CCO"}) Chemical("CC" ))

See also: find table find molcart other chemical functions

Eigen

eigenvalues/eigenvectors function.
Eigen ( M_ )
- returns the square matrix of eigenvector columns of the input symmetric square matrix M_ . Eigenvalues sorted by their values are stored in the R_out rarray.
Example:

 
 A = Matrix(3, 3, 0.)      # create a zero square matrix... 
 A[1:3,1] = {1.,-2.,-1.}   # and set its elements 
 A[2,2] =  4. 
 for i = 1, 3-1            # the matrix must be symmetric 
   for j = i+1,