This section describes the format of the alignment file and commands
for reading, writing, making, analyzing and using the alignments of
sequences and structures (pairwise and multiple). It also includes a
description of the command for searching a sequence database. For the
underlying dynamic programming methods see Section~\ref{SECTION:comparedesc}.

\vspace{5mm}

\input{align.tex}

\vspace{5mm}

\Command      {READ\_ALIGNMENT}{read sequences and/or their alignment}  

\Options{
\OptLine{FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGNMENT\_FORMAT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{REMOVE\_GAPS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ADD\_SEQUENCE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}, NUMB\_OF\_SEQUENCES, ALIGN\_CODES}

\Description  {This command reads the sequence(s) and/or their alignment from 
a text file. Only sequences with the specified codes are read in;
\K{ALIGN\_CODES} = \Z{all} can be used to read all sequences.

There are several alignment formats:

\begin{enumerate}
\item The \Z{PIR} format resembles that of the PIR sequence database. It
is described in Section~\ref{alignmentformat} and is used for
comparative modeling because it allows for additional data about 
the proteins that are useful for automated access to the atomic coordinates.

\item The \Z{FASTA} format resembles the \Z{PIR} format but has a 
missing second `comment' line and a missing star at the end of each 
sequence.

\item The \Z{PAP} format is nicer to look at but contains less
information and is not used by other programs. When used in 
conjunction with PDB files, the PDB files must contain 
exactly the residues in the sequences in the \Z{PAP} file; \ie, it
is not possible to use only a segment of a PDB file. In addition, the 
\Z{PAP} protein codes must be expandable into proper PDB atom filenames,
as described in Section~\ref{SECTIONFILENAME}. The protein sequence
can now start in any column (this was limited to column 11 before
release 5).

\item The \Z{QUANTA} format can be used to communicate with the 
\QUANTA\ program. You are not supposed to mix \Z{QUANTA} format with any other
format because the \Z{QUANTA} format contains residue numbers which do 
not occur in the other formats and are difficult to guess correctly.  

\item The \Z{INSIGHT} format is very similar to the \Z{PAP} format
and can sometimes be used to communicate with the
\INSIGHT\ program. When used in conjunction with PDB files, the
same rules as for the \Z{PAP} format apply.
\end{enumerate}

If \K{REMOVE\_GAPS} = \V{on}, positions with gaps in all selected
sequences are removed from the alignment.

If \K{ADD\_SEQUENCE} is \V{on}, the new sequences are added to the
current ones, otherwise the old sequences are deleted.
}


\Example{../examples/commands/read_alignment.top}



\Command      {READ\_ALIGNMENT2}{read 2nd alignment}    

\Options{
\OptLine{FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES2}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGNMENT\_FORMAT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{REMOVE\_GAPS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description  {This command reads the sequences and/or their alignment from 
a text file into the second alignment array. The two alignments can be
compared by the \C{COMPARE\_ALIGNMENTS} command.}

\SeeExample{\C{READ\_ALIGNMENT}}




\Command      {CHECK\_ALIGNMENT}{check alignment for modeling}  

\Options {
\OptLine{ATOM\_FILES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Description  {This command evaluates an alignment to be used for
comparative modeling. It uses two criteria:

First, it checks the alignment of the template structures (all but
the last entry in the alignment):
For each pairwise superposition of the templates,
it reports those equivalent pairs of \Calpha\ atoms that are 
more than 6\Ang\ away from each other. Such pairs are almost certainly 
misaligned. The pairwise superpositions are done using the \Calpha\ atoms
and the given alignment.

Second, the command checks the alignment of the target sequence (the last
entry in the alignment) with each of the templates: For all
consecutive pairs of \Calpha\ atoms in the target, it calculates
the distance between the two equivalent \Calpha\ atoms
in each of the templates. If the distance is longer than
8\Ang, it is reported. In such a case, the alignment between 
the template and the target is almost certainly incorrect.
}


\Example{../examples/commands/read_alignment.top}





\Command      {COLOR\_ALN\_MODEL}{color MODEL according to alignment}   

\Description  {This command colors MODEL according to a given alignment between
MODEL and a sequence. MODEL has to be the first protein in the alignment.
The second protein can be any sequence, with or without known structure.
The MODEL can be displayed on the screen, colored by `the fourth parameter' 
and inspected for the structural context of deletions and insertions. This
is useful in optimizing the alignment for comparative modeling. The isotropic 
temperature factors in MODEL are set as follows: 

\begin{itemize}
\item 0, for those regions that have residues in both MODEL and the sequence 
         (blue in \RASMOL; light green in \QUANTA); 
\item 1, for the two residues that span regions occurring in the sequence but 
         not in MODEL (green in \RASMOL; pink in \QUANTA);
\item 2, regions that occur in MODEL but are deleted from the sequence 
         (red in \RASMOL; bright green in \QUANTA).
\end{itemize}
}


\Example{../examples/commands/color_aln_model.top}





\Command      {COMPARE\_ALIGNMENTS}{compare two alignments}     

\Requirements{\C{READ\_ALIGNMENT} \AND\ \C{READ\_ALIGNMENT2}}

\Description  {This command compares two pairwise alignments read by the
\C{READ\_ALIGNMENT} and \C{READ\_\-ALIGNMENT2} commands. The alignment
of the first sequence with the second sequence in ALIGNMENT is
evaluated with respect to ALIGNMENT2. The numbers are not symmetric;
\ie, they will change if the sequences or alignments are swapped.
The output in the \logf\ file is self-explanatory.}


\Example{../examples/commands/compare_alignments.top}




\Command      {SEQUENCE\_TO\_ALI}{copy MODEL sequence and coordinates 
               to alignment}

\Options{
\OptLine{ADD\_SEQUENCE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output       {NUMB\_OF\_SEQUENCES, ALIGN\_CODES, ATOM\_FILES}

\Description  {This command copies the sequence and coordinates of the 
MODEL to the alignment arrays. 

If \K{ADD\_SEQUENCE} is \V{on} the sequence is added to the sequences
that are already in the alignment arrays, otherwise it becomes the
only sequence in those arrays. 

When sequence $i$ is added, the corresponding elements of
\K{ALIGN\_CODES} and \K{ATOM\_FILES} are used to set the protein and
PDB code fields in the alignment file, respectively.
}

\SeeExample{\C{COMPARE\_ALIGNMENTS}}



\Command      {WRITE\_ALIGNMENT}{write sequences and/or their alignment}        

\Options{
\OptLine{FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGNMENT\_FORMAT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_BLOCK}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGNMENT\_FEATURES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{CUT\_OVERHANGS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{HETATM\_IO}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{WATER\_IO}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{HYDROGEN\_IO}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Description  {This command writes the whole alignment to a text file.

The \Z{PAP} format, which corresponds to a relatively nice looking
alignment, has several additional formatting options that can be
selected by the \K{ALIGNMENT\_FEATURES} variable.
This scalar variable can contain any combination of the following keywords:

\begin{itemize}
\item \Z{INDICES}, the alignment position indices;

\item \Z{CONSERVATION}, a star for each absolutely conserved position;

\item \Z{ACCURACY}, the alignment accuracy indices, scaled between 0--9,
as calculated by \C{ALIGN\_\-CONSENSUS};

\item \Z{HELIX}, average content of helical residues for
structures 1 -- \K{ALIGN\_BLOCK} at each position,
0 for 0\% and 9 for 100\%, as calculated by \C{ALIGN2D}.

\item \Z{BETA}, average content of \be-strand residues for
structures 1 -- \K{ALIGN\_BLOCK} at each position,
0 for 0\% and 9 for 100\%, as calculated `by \C{ALIGN2D}.

\item \Z{ACCESSIBILITY}, average relative sidechain buriedness for 
structures 1 -- \K{ALIGN\_BLOCK}, 0 for 0\% (100\% accessibility) and 
9 for 100\% (0\% accessibility), as calculated by \C{ALIGN2D};

\item \Z{STRAIGHTNESS}, average mainchain straightness
structures 1 -- \K{ALIGN\_BLOCK} at each position
0 for 0\% and 9 for 100\%, as calculated by \C{ALIGN2D}.
\end{itemize}

Options \Z{HELIX}, \Z{BETA}, \Z{ACCESSIBILITY}, and \Z{STRAIGHTNESS} 
are valid only after executing command \C{ALIGN2D}, where the 
corresponding quantities are defined. They refer to the 3D profile 
defined for the first \K{ALIGN\_BLOCK}
structures (run \C{ALIGN2D} with \K{FIT} = \V{off} to prepare
these structural data without changing the input alignment).
Similarly, the \Z{ACCURACY} option is valid only after the 
\C{CONSENSUS\_ALIGNMENT} command.

\K{ALIGN\_ALIGNMENT} and \K{ALIGN\_BLOCK} are used to ensure correct 
indication of identical alignment positions, depending on whether
sequences or two blocks of sequences were aligned: For sequences 
(\K{ALIGN\_ALIGNMENT} = \V{off} and \K{ALIGN\_BLOCK} is ignored), a
\Z{*} indicating a conserved position is printed where all sequences
have the same residue type. For blocks (\K{ALIGN\_ALIGNMENT} = \V{on} 
and \K{ALIGN\_BLOCK} indicates the last sequence of the first block),
a \Z{*} is printed only where
the two blocks have the same order of residue types (there has to be
the same number of sequences in both blocks). The blocks option is
useful when comparing two alignments, possibly aligned by the
\C{ALIGN} command.

If \K{CUT\_OVERHANGS} is \V{on}, overhangs longer than \K{OVERHANG}
residues are cut from the alignment. In such a case,
the \K{HETATM\_IO}, \K{WATER\_IO}, \K{HYDROGEN\_IO}, and
\K{ATOM\_\-FILES\_\-DIRECTORY} keywords also apply because the beginning
and ending residue numbers for the {\tt `structure'} entries in the alignment
are renumbered automatically by reading the appropriate atom files.
}

\SeeExample{\C{READ\_ALIGNMENT}}



\Command      {DESCRIBE}{describe proteins}

\Options{
\OptLine{ATOM\_FILES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Requirements {[alignment]}

\Description  {This command outputs basic data about the proteins whose
atom filenames are specified by \K{ATOM\_\-FILES} or \K{ALIGN\_CODES}.
An alternative specification of the proteins to be described can be
provided by the alignment in memory; \ie, \C{READ\_ALIGNMENT}
followed by \C{DESCRIBE} will describe all the proteins in the
alignment. This command is useful for preparation before comparative
modeling because it summarizes disulfides, \cis-prolines, charges,
chain breaks, \etc.  When an alignment is given, results depending
only on the amino acid sequences are still written out even if some 
atom files do not exist.}


\Example{../examples/commands/describe.top}




\Command      {ID\_TABLE}{calculate percentage sequence identities}     

\Options{
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Requirements {alignment}

\Description  {This command calculates percentage residue identities for 
all pairs of sequences in the current alignment.  The percentage
residue identity is defined as the number of identical residues
divided by the length of the shorter sequence.

The \K{ALIGN\_CODES} variable is only used for output, not in
calculations, so it does not have to be set. 

In addition to the output in the \logf\ file, this routine creates file
\K{OUTPUT\_DIRECTORY}/\-\K{MATRIX\_\-FILE} with pairwise sequence
distances that can be used directly as the input to the tree making
programs of the \PHYLIP\ package, such as \KITSCH\ \cite{Fel85}, 
and also for the \C{DENDROGRAM} and \C{PRINCIPAL\_\-COMPONENTS} 
commands. A more general
version of this command, which allows a user specified measure for
residue--residue differences is \C{SEQUENCE\_COMPARISON}.  }


\Example{../examples/commands/id_table.top}





\Command      {SEQUENCE\_COMPARISON}{compare sequences in alignment}

\Options{
\OptLine{RR\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{VARIABILITY\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MAX\_GAPS\_MATCH}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Description  {The pairwise similarity of sequences in the current alignment 
is evaluated using a user specified residue--residue scores file. 

The residue--residue scores, including gap--residue, and gap--gap scores, 
are read
from file \K{RR\_FILE}. The sequence pair score is equal to the average
pairwise residue--residue score for all alignment positions that
have at most \K{MAX\_GAPS\_MATCH} gaps (1 by default). If the 
gap--residue and gap--gap scores are not defined in \K{MATRIX\_FILE},
they are set to the worst and best residue--residue score, respectively.
If \K{MATRIX\_FILE} is a similarity matrix, it is converted into 
a distance matrix ($x' = -x + x_{\mbox{max}}$).

The comparison matrix is written in the \PHYLIP\ format to file
\K{MATRIX\_FILE}.

The family variability as a function of alignment position is
calculated as the \RMS\ deviation of all residue -- residue scores at a given
position, but only for those pairs of residues that have at most 
\K{MAX\_\-GAPS\_\-MATCH} gaps (0, 1, or 2). The variability is written to file 
\K{VARIABILITY\_FILE}, as is the number of pairwise comparisons contributing 
to each positional variability.}

\SeeExample{\C{ID\_TABLE}}


\Command      {DENDROGRAM}{clustering}

\Options{
\OptLine{MATRIX\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Description  {This command calculates a clustering tree 
from the input matrix of pairwise distances. This matrix must be in 
the \PHYLIP\ format and can be produced by the \C{ID\_TABLE},
\C{SEQUENCE\_\-COMPARISON}, or \C{COMPARE} commands. 
The weighted pair-group average clustering method is used.

The tree is written to the \logf\ file.

This command is useful for deciding about which known
3D structures are to be used as templates for comparative modeling.  }

\SeeExample{\C{ID\_TABLE}}



\Command      {PRINCIPAL\_COMPONENTS}{clustering}

\Options{
\OptLine{MATRIX\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE}{TYPE}{VALUES}{DEFAULT}{output file} 
}

\Description  {This command calculates principal components clustering
for the input matrix of pairwise distances. This matrix must be in 
the \PHYLIP\ format and can be produced by the \C{ID\_TABLE},
\C{SEQUENCE\_\-COMPARISON}, or \C{COMPARE} commands.

The projected coordinates $p$ and $q$ are written to file \K{FILE}.
The output file can be used with \ASGL\ to produce a principal components 
plot.

This command is useful for deciding about which known
3D structures are to be used as templates for comparative modeling.  }

\SeeExample{\C{ID\_TABLE}}


\Command      {ALIGN}{align two (blocks of) sequences}

\Options{
\OptLine{RR\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{GAP\_PENALTIES\_1D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_BLOCK}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OFF\_DIAGONAL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_OFFSET}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{LOCAL\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_WHAT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{READ\_WEIGHTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{WRITE\_WEIGHTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{INPUT\_WEIGHTS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT\_WEIGHTS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description{This command aligns two blocks of sequences.

The two blocks of sequences to be aligned are sequences 1 to 
\K{ALIGN\_BLOCK} and \K{ALIGN\_BLOCK}+1 to the last sequence. The
sequences within the two blocks should already be aligned; their
alignment does not change.

{\sloppy The command can do either the global (similar to \cite{NeeWun70};
\K{LOCAL\_\-ALIGNMENT} = \V{off}) or local dynamic programming alignment
(similar to \cite{SmiWat81}; \K{LOCAL\_\-ALIGNMENT} = \V{on}).}

For the global alignment, set overhang length
\K{OVERHANG} to more than 0 so that the corresponding number of
residues at either of the four termini won't be penalized by any gap
penalties (this makes it a pseudo local alignment). 

To speed up the calculation, set \K{OFF\_DIAGONAL} to a number smaller than
the shortest sequence length. The alignments matching residues $i$ and 
$j$ with $|i-j| > \K{OFF\_DIAGONAL}$ are not considered at all in the 
search for the best alignment.

The gap initiation and extension penalties are specified by
\K{GAP\_PENALTIES\_1D}. The default values of -900 -50 for the
\Z{as1.sim.mat} similarity matrix were found to be optimal for
pairwise alignments of sequences that share from 30\% to 45\%
sequence identity (RS and A\v{S}, in preparation).

The residue type -- residue type scores are read from file
\K{RR\_FILE}. The routine automatically determines whether it has to
maximize similarity or minimize distance.

\K{MATRIX\_OFFSET} applies to local alignment only and influences its
length. \K{MATRIX\_OFFSET} should be somewhere between the lowest and
highest residue--residue scores. A smaller value of this parameter
will make the local alignments shorter when distance is minimized, 
and longer when similarity is maximized. This works as follows: The 
recursively constructed dynamic programming comparison matrix is reset to 0 at
position $i,j$ when the current alignment score becomes larger
(distance) or smaller (similarity) than \K{MATRIX\_OFFSET}.  Note that
this is equivalent to the usual shifting of the residue--residue
scoring matrix in the sense that there are two combinations of
\K{GAP\_PENALTIES\_1D} and \K{MATRIX\_OFFSET} values that will give
exactly the same alignments irrespective of whether the matrix is
actually offset (with 0 used to restart local alignments in dynamic
programming) or the matrix is not offset but \K{MATRIX\_OFFSET} is
used as the cutoff for restarting local alignments in dynamic
programming. For the same reason, the matrix offset does not have any
effect on the global alignments if the gap extension penalty is also
shifted for half of the matrix offset.  

The position--position score is an average residue--residue score for
all possible pairwise comparisons between the two blocks ($n \times m$ 
comparisons are done, where $n$ and $m$ are the number of sequences in 
the two blocks, respectively). The first exception to this is when 
\K{ALIGN\_WHAT} is set to \Z{ALIGNMENTS}, in which case
the two alignments defined by \K{ALIGN\_BLOCK} are aligned; \ie, the
score is obtained by comparing only equivalent positions between the
two alignment blocks (only $n$ comparisons are done, where $n$ is the number of
sequences in each of the two blocks). This option is useful in
combination with \C{COMPARE\_ALIGNMENTS} and \C{WRITE\_ALIGNMENT} for
evaluation of various alignment parameters and methods. The second
exception is when \K{ALIGN\_WHAT} is set to \Z{LAST}, in which case
only the last sequences in the two blocks are used to get the scores.
In \Z{block}, \Z{alignment}, and \Z{last} comparisons, penalty for a 
comparison of a gap with a
residue during the calculation of the scoring matrix is obtained from
the score file (gap--gap match should have a score of 0.0). 

Only the 20 standard residue types, plus Asx (changes to Asn) and Glx
(changes to Gln) are recognized. Every other unrecognized residue,
except for a gap and a chain break, changes to Gly for comparison
purposes.

If you receive an error message to increase the \V{MAXRES} constant,
you can try to increase the gap penalties first.  Here and elsewhere
in \MODELLER, \V{MAXRES} is both the maximal number of residues in a 
protein as well as the maximal length of an alignment. If the length 
of the alignment arrays is too small, \K{MODELLER\_STATUS} becomes 1
(Section~\ref{SECTION:breakpoint}).

For the time being, this and the other alignment commands (\C{MALIGN},
\C{ALIGN2D}, \C{ALIGN3D}, and \C{MALIGN3D}) remove chain break information 
from the \V{CALN} array, which means that chain breaks are not retained when 
the alignment is written to a file after executing these commands. }


\Example{../examples/commands/align.top}



\Command      {ALIGN2D}{align sequences with structures}
\label{COMMAND:align2d}

\Options{
\OptLine{RR\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{directory of \K{RR\_FILE}} 
\OptLine{GAP\_PENALTIES\_1D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{GAP\_PENALTIES\_2D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_BLOCK}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OFF\_DIAGONAL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_OFFSET}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{LOCAL\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FIT}{TYPE}{VALUES}{DEFAULT}{whether to align} 
\OptLine{READ\_WEIGHTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{WRITE\_WEIGHTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{INPUT\_WEIGHTS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{OUTPUT\_WEIGHTS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description{This command aligns a block of sequences (second block)
with a block of structures (first block). It is the same as the 
\C{ALIGN} command except that a variable gap {\bf opening} penalty is 
used. This gap penalty depends on the 3D structure of all sequences in
block 1. The variable gap penalty can favor gaps in exposed regions,
avoid gaps within secondary structure elements, favor gaps in
curved parts of the mainchain, and minimize the distance between the two
\Calpha\ positions spanning a gap.  The \C{ALIGN2D} command is
preferred for aligning a sequence with structure(s) in comparative
modeling because it tends to place gaps in a better structural
context. See Section~\ref{SECTION:variablegap} for the dynamic
programming algorithm that implements the variable gap penalty.
\K{GAP\_PENALTIES\_2D} specifies parameters
$\omega_H$, $\omega_S$, $\omega_B$, $\omega_C$, $\omega_D$,
$d_o$, $\gamma$, and $t$ (Section~\ref{SECTION:variablegap}).
The default gap penalties \K{GAP\_PENALTIES\_1D} ($-450, 0$) 
and \K{GAP\_PENALTIES\_2D} (0.35, 1.2, 0.9, 1.2, 0.6, 8.6, 1.2, 0.0) 
as well as the \K{RR\_FILE} substitution
matrix (\Z{as1.sim.mat}) were found to be optimal in pairwise alignments of
structures and sequences sharing from 30\% to 45\% sequence identity
(RS and A\v{S}, in preparation).

--- move to back

The linear gap penalty function for inserting a gap in block 1 of
structures is: 
$g = f_1(H, S, B, C) u + l v$ where $u$ and $v$ are the
usual gap opening and extension penalties, $l$ is gap length, and
$f_1$ is a function that is at least 1, but can be larger to make gap
opening more difficult in the following circumstances: between two consecutive 
(\ie, $i, i+1$) helical positions, two consecutive \be-strand positions, two
consecutive buried positions, or two consecutive positions where 
the mainchain is locally straight. This function is $f_1 = 1 + [\omega_H
H_i H_{i+1} + \omega_S S_i S_{i+1} + \omega_B B_i B_{i+1} +
\omega_C C_i C_{i+1}]$,
$H_i$ is the fraction of helical residues at
position $i$ in block 1, $S_i$ is the fraction of \be-strand residues
at position $i$ in block 1, $B_i$ is the average relative sidechain
buriedness of residues at position $i$ in block 1, and $C_i$ is the
average straightness of residues at position $i$ in block 1. See 
Section~\ref{COMMAND:writedata} for the definition of these features.
The original straightness is modified here by assigning maximal straightness 
of 1 to all residues in a helix or a \be-strand. 

The linear gap penalty function for opening a gap in block 2 of
sequences is: \linebreak $g = f_2(H, S, B, C, D) u + l v$ where $f_2$ is a
function that is at least 1, but can be larger to make the gap
opening in block 2 more difficult in the following circumstances:
when the first gap position is
aligned with a helical residue, a \be-strand residue, a buried residue,
extended mainchain, or when the whole gap in block 2 is spanned by two
residues in block 1 that are far apart in space. This function is $f_2
= 1 + [\omega_H H_i + \omega_S S_i + \omega_B B_i + \omega_C C_i +
\omega_D \sqrt{d-d_o}]$. $d$ is the distance between the two \Calpha\ 
atoms spanning the gap, averaged over all structures in block 1 and
$d_o$ is the distance that is small enough to correspond to no
increase in the opening gap penalty (\eg, 8.6\Ang).

When \K{FIT} is \V{off}, no alignment is done and the routine returns
only the average structural information, which can be written out by the
\C{WRITE\_ALIGNMENT} command.
}


\Example{../examples/commands/align2d.top}



\Command      {MALIGN}{align two or more sequences}

\Options{
\OptLine{RR\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{GAP\_PENALTIES\_1D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OFF\_DIAGONAL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_BLOCK}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_OFFSET}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{LOCAL\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description  {This command performs a multiple sequence alignment. 
The sequences to be aligned are the sequences in the current alignment
arrays.  The command uses the dynamic programming method for the best
sequence alignment, given the gap initiation and extension penalties
specified by \K{GAP\_PENALTIES\_1D}, and residue type weights read from
file \K{RR\_FILE}.  See command \C{ALIGN} for more information.

The algorithm for the multiple alignment is as follows. First,
sequence 2 is aligned with sequence 1 (\ie, block of sequences 
from 1--\K{ALIGN\_BLOCK}). Next, sequence 3 is aligned
with an average of the aligned sequences 1 and 2; \ie, the weight
matrix is an average of the weights 1--3 and 2--3. For this
averaging, the gap--residue and gap--gap weights are obtained from the 
residue--residue weight matrix file, not from gap penalties. If
the corresponding weights are not in the file, they are set to 
the worst and best residue--residue score, respectively.

See instructions for \C{ALIGN} for more details.
}


\Example{../examples/commands/malign.top}





\Command      {ALIGN\_CONSENSUS}{consensus sequence alignment}

\Options{
\OptLine{GAP\_PENALTIES\_1D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_BLOCK}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{READ\_WEIGHTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{WRITE\_WEIGHTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{INPUT\_WEIGHTS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{OUTPUT\_WEIGHTS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
}


\Output {\K{MODELLER\_STATUS} = \Integer{1}} 

\Description {This command is similar to \C{ALIGN} except that a
consensus alignment of two blocks of sequences is produced. A
consensus alignment is obtained from a consensus similarity matrix
using the specified gap penalties and the global dynamic programming
method.  The consensus similarity matrix is obtained by aligning the
two blocks of sequences many times with different parameters and
methods and counting how many times each pair was aligned. This command
is still experimental and no detailed description is given at this time.

This command also produces the alignment accuracy that can be printed
out by the \C{WRITE\_ALIGNMENT} command in the \Z{PAP} format (0
inaccurate, 9 accurate). If the gap initiation penalty is 0, the
gap extension penalty of say 0.4 means that only those positions will
be equivalenced that were aligned in at least 80\% of the individual 
alignments (\ie, 2 times 0.40).}


\Example{../examples/commands/align_consensus.top}





\Command      {SUPERPOSE}{superpose MODEL2 on MODEL given alignment}

\Options{
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FIT}{TYPE}{VALUES}{DEFAULT}{whether to superpose} 
\OptLine{SUPERPOSE\_REFINE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{RMS\_CUTOFFS}{TYPE}{VALUES}{DEFAULT}{only the first element is used for calculating the cutoff RMS and DRMS measures}
\OptLine{REFERENCE\_ATOM}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{REFERENCE\_DISTANCE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{SWAP\_ATOMS\_IN\_RES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
}

\Requirements {MODEL \AND\ MODEL2 [\AND\ alignment]}

\Description  {This command superposes MODEL2 on MODEL, without
changing the alignment.

If an alignment is in memory, it is used to obtain the equivalent
atoms.  MODEL must be the first sequence in the alignment, MODEL2 must
be the second sequence in the alignment.  The equivalent atoms are
those selected atoms (set 1) of the MODEL that have equivalently named atoms
in MODEL2; the atom equivalences are defined in library {\tt
\$ATMEQV\_LIB}.  Use the \C{PICK\_ATOMS} command to select the desired
atoms for superposition.  By default, all atoms are selected.
If there is no alignment, a 1:1 correspondence between the residues is
assumed. 

No fitting is done if \K{FIT} = \V{off}.

The \K{ALIGN\_CODES} variable is used only for output, not in
calculations. 

The \K{RMS\_CUTOFFS}[1] element is the cutoff used in calculating the
cutoff \RMS\ deviations; \ie, those position and distance \RMS\ deviations 
that are defined
on the equivalent atoms which are less than \K{RMS\_\-CUTOFFS}[1] 
angstroms away from each other (as superposed using all aligned positions) 
and those equivalent distances which are less than \K{RMS\_CUTOFFS}[1] 
angstroms different from each other, respectively.  

If \K{SUPERPOSE\_REFINE} is \V{on} the refinement of the superposition
is done by repeating the fitting with only those aligned pairs of
atoms that are within \K{RMS\_CUTOFFS}[1] of each other until there is
no change in the number of equivalent positions. This refinement
can only remove compared positions, not add them like \C{ALIGN3D} can
do. This is useful for comparing equivalent parts of two structures
with a fixed alignment but omitting divergent parts from the
superposition and \RMS\ deviation calculation; \eg, comparing a model with 
the X-ray structure.

If \K{SUPERPOSE\_REFINE} is \V{off} and \K{REFERENCE\_ATOM} is non-blank, only
those pairs of equivalently named selected atoms from aligned residues 
are superposed that come from residues whose \K{REFERENCE\_ATOM} atoms are
closer than \K{REFERENCE\_DISTANCE} \Ang\ to each other.

When MODEL and MODEL2 have exactly the same atoms in the same order, 
one can set \K{SWAP\_ATOMS\_IN\_RES} to any combination of
single character amino acid residue codes in \V{DEFHLNQRVY}.
Certain atoms (see below) in the specified sidechains of MODEL2 are 
then swapped to minimize their \RMS\ deviation relative to MODEL.
The labelling resulting in the lowest \RMS\ deviation is retained. 
The following swaps are attempted: 

\begin{center}
\begin{tabular}{ll} \hline
Residue & Swap(s)  \\ \hline
D       & OD1, OD2 \\
E       & OE1, OE2 \\
F       & CD1, CD2 \\
        & CE1, CE2 \\
H       & ND1, CD2 \\
        & NE2, CE1 \\
N       & OD1, ND2 \\
Q       & OE1, NE2 \\
R       & NH1, NH2 \\
V       & CG1, CG2 \\
Y       & CD1, CD2 \\
        & CE1, CE2 \\
\end{tabular}
\end{center}
}


\Example{../examples/commands/superpose.top}



\Example{../examples/commands/align3d.top}



\Example{../examples/commands/swap_atoms_in_res.top}



\Command      {COMPARE}{compare 3D structures given alignment}  

\Options{
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT}{TYPE}{VALUES}{DEFAULT}{selects output:
\Z{SHORT} \OR\ \Z{LONG}\} \OR\ \Z{RMS}\ \OR\ \Z{DRMS}}
\OptLine{MATRIX\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{COMPARE\_MODE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{RMS\_CUTOFFS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FIT\_ATOMS}{TYPE}{VALUES}{DEFAULT}{whether to superpose before comparison} 
\OptLine{DISTANCE\_ATOMS}{TYPE}{VALUES}{DEFAULT}{atom type used for variability calculations} 
\OptLine{FIT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ASGL\_OUTPUT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Description{This command compares the structures in the given alignment.
It does not make an alignment, but it calculates the \RMS\ and 
\DRMS\ deviations between atomic positions and distances, and class differences
between the mainchain and sidechain dihedral angles. In contrast to the 
\C{SUPERPOSE} command, \C{COMPARE} works with a multiple alignment and 
it writes more information about the pairwise comparisons.

If no alignment is available, it assumes a 1:1 correspondence for the
proteins specified by \K{ATOM\_FILES} or \K{ALIGN\_CODES}. If
\K{ATOM\_FILES} is defined, it is used with the \MODELLER\ file-naming
mechanism (Section~\ref{SECTIONFILENAME}) to find full names for the
atom files. If it is not defined, \K{ALIGN\_CODES}, which is
usually set by the previous \C{READ\_ALIGNMENT} command, is used.
\K{ALIGN\_CODES} does not have to be set if \K{ATOM\_FILES} is set.

\K{OUTPUT} selects short (\Z{SHORT}) or long (\Z{LONG})
form of output to the \logf\ file. If it contains word 
\Z{RMS} or \Z{DRMS} it also outputs the \RMS\ or \DRMS\ deviation matrix
to file \K{MATRIX\_FILE}. This file can be used with the 
\PHYLIP\ program or with the \C{DENDROGRAM} or \C{PRINCIPAL\_COMPONENTS}
commands of \MODELLER\ to calculate a clustering of the structures.
  
\K{COMPARE\_MODE} selects the form of the positional variability
calculated for each position along the sequence:

\begin{itemize}  
\item[1,] for true \RMS\ deviation over all proteins that have a residue at 
the current position. This does not make any sense for periodic quantities 
like dihedral angles.
  
\item[2,] for the average absolute distance over all pairs of 
residues that have a residue at the current position.
  
\item[3,] the same as 2 except that average distance, not its 
absolute value is used (convenient for comparison of 2 structures to 
get the $\pm$ sign of the changes for dihedral angles and distances).
\end{itemize}
  
  
\K{RMS\_CUTOFFS} specifies cutoff values for calculation of the 
position, distance, and dihedral angle \RMS\ deviations for pairwise
overall comparisons. If difference between two equivalent points is
larger than cutoff it is not included in the \RMS\ sum. The order of
cutoffs in this vector is: atomic position, intra-molecular distance,
$\alpha$, $\Phi$, $\Psi$, $\omega$, $\chi_1$, $\chi_2$, $\chi_3$, $\chi_4$, 
and $\chi_5$ (there are 5 dihedrals in a disulfide bridge), where
$\alpha$ is the virtual \Calpha\ dihedral angle between four consecutive
\Calpha\ atoms. These cutoffs do not affect positional variability 
calculations.

\K{FIT\_ATOMS} string specifies all the atom types (including possibly
a generic \Z{ALL}) to be fitted in the least-squares superposition.
These atom types are used in the least-squares superposition, and in
calculation of the position and distance \RMS\ deviations.

\K{DISTANCE\_ATOMS}[1] specifies the atom type that is used for
getting the average structure and \RMS\ deviation at each alignment position in
the \ASGL\ output file \Z{posdif.asgl}. This \ASGL\ file contains the
positional variability of the selected atom type in the family of
compared proteins.  The \ASGL\ output files can then be used with
\ASGL\ scripts \Z{posdif} and \Z{dih} to produce \PostScript\ plots of the 
corresponding variabilities at each alignment position.
\K{ASGL\_OUTPUT} has to be \V{on} to obtain the \ASGL\ output files.
  
If \K{FIT} = \V{on}, a least-squares superposition is done before the
comparisons; otherwise, the orientation of the molecules in the
input atom files is used.}

\SeeExample{\C{MALIGN3D}}




\Command      {ALIGN3D}{align two structures}

\Options{
\OptLine{GAP\_PENALTIES\_3D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FIT\_ATOMS}{TYPE}{VALUES}{DEFAULT}{one atom type used for superposition} 
\OptLine{FIT}{TYPE}{VALUES}{DEFAULT}{whether to align} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT}{TYPE}{VALUES}{DEFAULT}{\Z{SHORT} \OR\ \Z{LONG} \OR\ 
         \Z{VERY\_LONG}}
\OptLine{ALIGN3D\_TRF}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN3D\_REPEAT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OFF\_DIAGONAL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_OFFSET}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{LOCAL\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}


\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description  {This command uses the current alignment as the starting point 
for an iterative least-squares superposition of two 3D structures. This 
results in a new pairwise structural alignment. If no alignment is in memory, 
the initial alignment is the 1:1 alignment. A good initial alignment may
be obtained by sequence alignment (\C{ALIGN}). For superpositions,
only one atom per residue is used, as specified by \K{FIT\_ATOMS}[1].

The alignment algorithm is as follows. First, structure 2 is
least-squares fit on structure 1 using all the equivalent residue
positions in the initial alignment that have the specified atom type.
Next, the residue--residue distance matrix is obtained by calculating
Euclidean distances between all pairs of selected atoms from the two
structures. The alignment of the two structures is then obtained by
the standard dynamic programming optimization based on the
residue--residue distance matrix. 

\K{GAP\_PENALTIES\_3D}[1] is a gap creation penalty (usually 0), and 
\K{GAP\_PENALTIES\_3D}[2] is a gap extension penalty, say 1.75. This 
procedure identifies pairs of residues as equivalent when they have
their selected atoms at most 2 times \K{GAP\_PENALTIES\_3D}[2] angstroms
apart in the current orientation (this is so when the gap initiation
penalty is 0).  The reason is that an equivalence costs the distance
between the two residues while an alternative, the gap--residue and
residue-gap matches, costs twice the gap extension penalty.  

From the dynamic programming run, a new alignment is obtained.  Thus,
structure 2 can be fitted onto structure 1 again, using this new
alignment, and the whole cycle is repeated until there is no change in
the number of equivalent positions and until the difference in the
rotation matrices for the last two superpositions is very small.  At
the end, the framework, that is the alignment positions without gaps,
is written to the \logf\ file. 

If \K{FIT} is \V{off}, no alignment is done. 

If \K{OUTPUT} contains \Z{SHORT}, only the best alignment and its
summary are displayed. If \K{OUTPUT} contains \Z{LONG}, summaries 
are displayed for all initial alignments in each framework cycle.  
If \K{OUTPUT} contains \Z{VERY\_LONG}, all alignments are displayed.

If \K{ALIGN3D\_TRF} is \V{on}, the weights in the weight matrix are
modified distances \cite{SubLauLev93}.

If \K{ALIGN3D\_REPEAT} is \V{on}, three additional initial alignments
are tried and the one resulting in the largest number of equivalent
positions is selected.
}


\Example{../examples/commands/align3d.top}




\Command      {MALIGN3D}{align two or more structures}

\Options{
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ATOM\_FILES\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{GAP\_PENALTIES\_3D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OFF\_DIAGONAL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_OFFSET}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{LOCAL\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FIT\_ATOMS}{TYPE}{VALUES}{DEFAULT}{one atom type for superposition} 
\OptLine{FIT}{TYPE}{VALUES}{DEFAULT}{whether to align} 
\OptLine{OUTPUT}{TYPE}{VALUES}{DEFAULT}{\Z{SHORT} \OR \Z{LONG} \OR
         \Z{VERY\_LONG} \OR\ \Z{NO\_ALIGNMENT}}
\OptLine{WRITE\_FIT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{CURRENT\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{WRITE\_WHOLE\_PDB}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description  {This command uses the current alignment as the starting point 
for an iterative least-squares superposition of two or more 3D structures.
This results in a new multiple structural alignment. If no alignment is in
memory, the initial alignment is the 1:1 alignment.  A good initial
alignment may be obtained by sequence alignment (\C{MALIGN}). For
superpositions, only one atom per residue is used, as specified by
\K{FIT\_ATOMS}. The resulting alignment can be written to a file with
the \C{WRITE\_ALIGNMENT} command. The multiply superposed coordinates
remain in memory and can be used with such commands as \C{TRANSFER\_XYZ}
if \K{ATOM\_FILES} is not changed in the meantime. It is best to use 
the structure that overlaps most with all the other structures as the first 
protein in the alignment. This may prevent an error exit due to too few
equivalent positions during framework construction.

The alignment algorithm is as follows. There are several cycles, each
of which consists of an update of a framework and a calculation of a
new alignment; the new alignment is based on the superposition of the
structures onto the latest framework. The framework in each cycle is
obtained as follows. The initial framework consists of the atoms in
structure 1 that correspond to \K{FIT\_ATOMS}. If there is no specified
atom types in any of the residues at a given position, the coordinates 
for this framework position are approximated by the
neighboring coordinates.  Next, all other structures are fit to this
framework.  The final framework for the current cycle is then obtained
as an average of all the structures, in their fitted orientations, but
only for residue positions that are common to all of them, given the
current alignment. Another result is that all the structures are now
superposed on this framework.  Note that the alignment has not been
changed yet. Next, the multiple alignment itself is re-derived in
$N-1$ dynamic programming runs, where $N$ is the number of structures.
This is done as follows.  First, structure 2 is aligned with structure
1, using the inter-molecular atom--atom distance matrix, for all atoms
of the selected type, as the weight matrix for the dynamic programming
run. Next, structure 3 is aligned with an average of structures 1 and
2 using the same dynamic programming technique. Structure 4 is then
aligned with an average of structures 1--3, and so on.  Averages of
structures $i$--$j$ are calculated for all alignment positions where
there is at least one residue in any of the structures $i$--$j$ (this is
different from a framework which requires that residues from all
structures be present).  Note that in this step, residues out of the
current framework may get aligned and the current framework residues
may get unaligned.
Thus, after the series of $N-1$ dynamic programming runs, a new
multiple alignment is obtained. This is then used in the next cycle to
obtain the next framework and the next alignment. The cycles are
repeated until there is no change in the number of equivalent
positions.  This procedure is best viewed as a way to determine the
framework regions, not the whole alignment. The results from this
command are expected to be similar to the output of program 
\MNYFIT\ \cite{SutHanCar87}.

\K{GAP\_PENALTIES\_3D}[1] is a gap creation penalty (usually 0), and
\K{GAP\_PENALTIES\_3D}[2] is a gap extension penalty, say 1.75.  This
procedure identifies pairs of positions as equivalent when they have
their selected atoms at most 2 times \K{GAP\_PENALTIES\_3D}[2] angstroms
apart in the current superposition (this is so when the gap initiation
penalty is 0), as described for the \C{ALIGN3D} command. 

Argument \K{OUTPUT} can contain the following values:
\begin{itemize}
\item \Z{SHORT}, only the final framework is written to the \logf\ file.  
\item \Z{LONG}, the framework after the alignment stage in
       each cycle is written to the \logf\ file.
\item \Z{VERY\_LONG}, the framework from the framework stage in each cycle
      is also written to the log. 
\end{itemize}

If \K{WRITE\_FIT} is \V{on}, the fitted atom files are written out in
their final fitted orientations. Their filenames are the original
filenames with an extension {\tt .fit}. 

If \K{CURRENT\_DIRECTORY} is \V{on}, the output {\tt .fit} files will
go to the current directory. Otherwise, the output will be in the
directory with the original files.

If \K{WRITE\_WHOLE\_PDB} is \V{on}, the whole PDB files are written out;
otherwise only the parts corresponding to the aligned sequences are output.

If \K{FIT} is \V{off}, the initial alignment is not changed. This is
useful when all the structures have to be superimposed with the
initial alignment (\K{FIT} = \V{off} and \K{WRITE\_FIT} = \V{on}).  }


\Example{../examples/commands/malign3d.top}



\Command      {EXPAND\_ALIGNMENT}{put all models into alignment}

\Options{
\OptLine{ROOT\_NAME}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE\_ID}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{EXPAND\_CONTROL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE\_EXT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {alignment}

\Description  {
\K{ID1}, \K{ID2}, \K{ROOT\_NAME},
\K{FILE\_EXT}, and \K{FILE\_ID} are used to construct atom filenames
for all the models (Section~\ref{SECTIONFILENAME}). Next, all the models 
are added to the alignment, using the last sequence in the input alignment 
as the guide. This allows easy multiple superposition of all the templates
and models after comparative modeling.
}

\Example{../examples/commands/expand_alignment.top}


\Command      {SEQUENCE\_SEARCH}{search for similar sequences}

\Options{
\OptLine{RR\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE}{TYPE}{VALUES}{DEFAULT}{file with the target sequence} 
\OptLine{ALIGN\_CODES}{TYPE}{VALUES}{DEFAULT}{the code of the target sequence} 
\OptLine{DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{GAP\_PENALTIES\_1D}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OFF\_DIAGONAL}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MATRIX\_OFFSET}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OVERHANG}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{LOCAL\_ALIGNMENT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEARCH\_CHAINS\_LIST}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION}
\OptLine{SEARCH\_CHAINS\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEARCH\_GROUP\_LIST}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGNMENT\_FORMAT}{TYPE}{VALUES}{DEFAULT}{sequence file formats; has to be \Z{PIR}} 
\OptLine{ALIGNMENT\_FEATURES}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{REMOVE\_GAPS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEARCH\_TOP\_LIST}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT}{TYPE}{VALUES}{DEFAULT}{\Z{SHORT} \OR\ \Z{LONG}} 
\OptLine{STOP\_ON\_ERROR}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEARCH\_SORT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEARCH\_RANDOMIZATIONS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{RAND\_SEED}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FAST\_SEARCH}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FAST\_SEARCH\_CUTOFF}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{DATA\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SIGNIF\_CUTOFF}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Output {\K{MODELLER\_STATUS} = \Integer{1}}

\Description  {This command searches a sequence database
for proteins that are similar to a given target sequence. 

Target sequence is read from file \K{FILE}.

\K{ALIGN\_CODES} specifies the code of the target sequence in the
\K{FILE} file. If only one sequence is in the file, you can use
\K{ALIGN\_CODES} = \Z{all} to read it without bothering about the actual
sequence code.

\K{SEARCH\_CHAINS\_LIST} specifies a file that contains protein codes 
for the proteins to be compared with the target sequence.

The database sequences specified in \K{SEARCH\_CHAINS\_LIST} file must
occur in the \K{SEARCH\_\-CHAINS\_\-FILE} file.

The command uses the dynamic programming method for the best sequence
alignment, given the gap creation and extension penalties specified by
\K{GAP\_PENALTIES\_1D} and residue type scores read from file
\K{RR\_FILE}. \K{GAP\_PENALTIES\_1D}[1] is a gap creation penalty 
and \K{GAP\_PENALTIES\_1D}[2] is a gap extension penalty.

The \K{SEARCH\_TOP\_LIST} top hits are written to the \logf\ file at the
end. The hits are sorted according to the fractional sequence
identity score obtained by dividing the number of identical residue pairs
by the length of the longer sequence (\K{SEARCH\_SORT} =
\Z{LONGER}) or the shorter sequence (\K{SEARCH\_SORT} = 
\Z{SHORTER}). 

The final list of hits contains three different significance values:

\begin{enumerate}
\item \V{SIGNI}. Z-score from sequence randomizations. This is the most
accurate significance score, but the slowest one to calculate. For each 
pairwise comparison, the two sequences are shuffled a specified number of 
times (\K{SEARCH\_RANDOMIZATIONS}) to obtain the mean and standard deviation 
of ``random'' scores from which the Z-score for an alignment score of a 
given pair of sequences is calculated.

\item \V{SIGNI2}. Z-score for sequence identity from the database scan.
After comparison of the target sequence with all sequences in the 
database is done, the comparisons are sorted by the length of the database
sequence. The pairwise sequence identities of the 20 sequences closest
in length to the target sequence are used to calculate the average and 
standard deviation of the percentage sequence identities for subsequent
calculation of the Z-score for the percentage sequence identity of a 
given pairwise alignment.

\item \V{SIGNI3}. Z-score for alignment score from the database scan.
The procedure is the same as for \V{SIGNI2}, except that the alignment
scores are used instead of the pairwise sequence identities.
\end{enumerate}

The calculation of the Z-scores assumes that the random scores
are distributed according to the Gaussian distribution, instead
of the extreme value distribution \cite{KarAlt90}, which is more
correct.

\K{SEARCH\_RANDOMIZATIONS} specifies how many alignments of
the shuffled sequences are done to calculate the significance score
for the overall sequence similarity. If 0, the significance is not
calculated. If more than 5 randomizations are done, the significance
score, not sequence identity, is used for sorting the hit list.

When \K{FAST\_SEARCH} is \V{on} only those sequences that have a 
database-scan alignment score significance (\V{SIGNI3} in output) 
above \K{FAST\_SEARCH\_CUTOFF} are used for the ``full'' 
randomization-based significance calculation. Since the mean
and the standard deviation of the distribution obtained by 
randomizing the two compared sequences are much more appropriate
than the corresponding quantities for the target/database comparisons,
\K{FAST\_SEARCH} should be \V{on} only when you are in a hurry and
the database is large.

If \K{DATA\_FILE} is \V{on} the final results (list of PDB codes with 
significances, \etc.) are also written to a separate file 
\Z{seqsearch.dat}.

If \K{OUTPUT} is \Z{LONG}, the best alignment for each sequence in
\K{SEARCH\_CHAINS\_FILE} and its various scores are also written to the 
\logf\ file. If \K{OUTPUT} is \Z{VERY\_LONG}, individual scores obtained for
randomized sequences are also written to the \logf\ file (this is almost
never needed).

If the selected significance score is larger than \K{SIGNIF\_CUTOFF}[1]
and not more than \K{SIGNIF\_CUTOFF}[2] units worse than the best
hit, all the members of the same group, as defined in 
\K{SEARCH\_GROUP\_LIST}, are added to the alignment array.
Subsequent \C{MALIGN}, \C{DENDROGRAM} and \C{WRITE\_ALIGNMENT} can then 
be used to write out all related PDB chains aligned to the target sequence.
}


\Example{../examples/commands/sequence_search.top}





\Command      {DELETE\_ALIGNMENT}{delete alignment}     

\Description  {This command deletes an existing alignment from 
the \MODELLER\ memory.  This is useful when a default 1:1
correspondence, such as that between an X-ray structure and its
\MODELLER\ model, is needed. This default alignment is constructed for
the commands that need an alignment only if there is no alignment
already in memory.}

\SeeExample{\C{PATCH}}



\Command      {SEGMENT\_MATCHING}{align segments}  

\Options{
\OptLine{RR\_FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ALIGN\_BLOCK}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEGMENT\_REPORT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEGMENT\_CUTOFF}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEGMENT\_SHIFTS}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEGMENT\_GROWTH\_N}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{SEGMENT\_GROWTH\_C}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{MIN\_LOOP\_LENGTH}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{OUTPUT\_DIRECTORY}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{ROOT\_NAME}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE\_ID}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
\OptLine{FILE\_EXT}{TYPE}{VALUES}{DEFAULT}{DESCRIPTION} 
}

\Requirements{alignment}

\Description  {This command enumerates alignments between 
two blocks of sequences. More precisely, it enumerates
the alignments between the segments in the first block and the
sequences in the second block. The segments can be moved 
to the left and right as well as lengthened and shortened, relative
to the initial alignment. The regions not in segments or
not aligned with segments are left un-aligned, possibly to 
be modeled as insertions. Typically, the first 
block of sequences corresponds to structures, the segments to
secondary structure elements in the first block, and the second block 
to the sequences one of which is to be modeled later on. The command 
is useful for generating many alignments which can then be used 
by another \MODELLER\ script to generate and 
evaluate the corresponding 3D models. 

All the sequences and segments are defined in the alignment 
array. The first block of sequences, the ones with segments,
are the first \K{ALIGN\_BLOCK} sequences. The regions 
corresponding to the segments are defined by the last 
entry in the alignment as contiguous blocks of non-gap 
residues. Any standard single character residue code may be used. The 
segments must be separated by gap residues, `-'. The remaining
sequences from \K{ALIGN\_BLOCK} + 1 to $\mbox{NSEQ}-1$ are the 
second block of sequences. The alignment
of the sequences within the two blocks does not change.
A sample alignment file is

% \verbatimfile{../examples/commands/segment_matching.pap}

The enumeration of alignments explores all possible combinations
of alignments between each segment and the 2nd block of sequences:
The starting position of each segment $i$ is varied relative to
the input alignment in the interval from \K{SEGMENT\_SHIFT}[$2 i - 1$] 
to \K{SEGMENT\_SHIFT}[$2 i$]. There has to be at
least \K{MIN\_LOOP\_LENGTH}[$i$] and \K{MIN\_LOOP\_LENGTH}[$i+1$]
residues that are not in any segment before and after the $i$-th segment, 
respectively. The location of the N-terminus of segment $i$ is varied 
relative to the location in the input alignment in the interval from
\K{SEGMENT\_GROWTH\_N}[$2 i - 1$] to \K{SEGMENT\_GROWTH\_N}[$2 i$].
Similarly, the location of the C-terminus of segment $i$ is varied 
relative to the location in the input alignment in the interval from
\K{SEGMENT\_GROWTH\_C}[$2 i - 1$] to \K{SEGMENT\_GROWTH\_C}[$2 i$].
The shortening and lengthening of the segments may be useful in 
determining the best anchor regions for modeling of a loop. 

Each alignment is scored according to the
similarity scoring matrix specified by filename \K{RR\_FILE}. This matrix may
contain residue---gap scores, the gap being residue type 21; otherwise
the value is set to the smallest value in the matrix. The score
for an alignment is obtained by summing scores only over 
all alignment positions corresponding 
to the segments (no gap penalty is added for loops). When there is more than
one sequence in any of the two blocks, the position score is an average
of all pairwise comparisons between the two blocks of sequences. In the
case where the number of positions in the alignment changes (\ie, the
segments grow or shorten), the scores are not comparable to each other.
It is feasible to enumerate on the order of $10^{10}$ different alignments
in less than one hour of CPU time.

In general, two runs are required. In the first run, the alignments
are scored and a histogram of the scores is written to file
\K{FILE}. Then this file must be inspected to determine the cutoff 
\K{SEGMENT\_CUTOFF}. In the second run, all the alignments with a score 
higher than \K{SEGMENT\_CUTOFF} are written to files in the 
\V{PIR} format, using the standard filenaming convention: 
\K{OUTPUT\_DIRECTORY}/\K{ROOT\_NAME}\K{FILE\_ID}$nnnn0000$\K{FILE\_EXT},
where $nnnn$ is the alignment file counter. In addition, the alignments
are also written out in the \V{PAP} format for easier inspection by
eye. Thus, \K{SEGMENT\_CUTOFF} has to be set to a very large value in
the first run, to avoid writting alignment files. During a run,
a message is written to the \logf\ every \K{SEGMENT\_REPORT} 
aligments; this is useful for knowing what is going on during very
long runs.
}

% \Example{../examples/commands/segment_matching.top}


\vfill
