MsaMultipleAnlignmentClasses {msa} | R Documentation |
MsaAAMultipleAlignment
,
MsaDNAMultipleAlignment
, and MsaRNAMultipleAlignment
S4 classes for storing multiple alignments of amino acid, DNA, and RNA sequences along with algorithm metadata
Objects of these classes are returned by the multiple sequence
alignment algorithms msaClustalW
,
msaClustalOmega
, msaMuscle
, and the
wrapper function msa
, all of which are
provided by the msa package.
The class MsaAAMultipleAlignment
extends the
AAMultipleAlignment
class, the class
MsaDNAMultipleAlignment
extends the
DNAMultipleAlignment
class, and the class
MsaRNAMultipleAlignment
extends the
RNAMultipleAlignment
class. All three classes
extend their parent classes by the slots contained in the
MsaMetaData
, i.e. all three classes are class
unions of the aforementioned parent classes and the class
MsaMetaData
.
print(x, show=c("alignment", "version", "call"),
showNames=TRUE, showConsensus=TRUE, halfNrow=9, nameWidth=20, ...)
:prints information about the object x
; the show
argument allows for determining what should be printed.
The show
must be a character vector and may contain any
combination of the following strings:
if show
contains "alignment"
, the multiple
sequence alignment is printed in a way similar to the
corresponding method from the Biostrings package
(except for the consensus sequence, see below).
If show
contains "complete"
, the entire width of
the alignment is printed by splitting it over multiple blocks of
lines if necessary. This overrules "alignment"
if both
are contained in the show
argument.
If show
contains "version"
,
the version
slot is shown. If show
contains
"call"
, the call
slot is shown.
If show
contains "standardParams"
, the
settings of the parameters that are common to all three
multiple sequence alignment algorithms are shown. If show
contains "algParams"
, the
algorithm-specific parameters are shown.
The order in which the strings are placed in the show
argument does not have an effect on the order in which
data are printed. The default is
show=c("alignment", "version", "call")
, i.e. by default,
the multiple sequence alignment is shown along with version and
call information. If show
contains "all"
, the
complete alignment is shown along with version information,
call, and the complete set of parameters.
As said above, by default, printing alignments is similar to
the standard print
method provided by the Biostrings
package, whereas including "complete"
in the argument
show
prints the entire width of the alignment.
Unlike the method from the Biostrings
package, the appearance can be customized: by default,
the consensus sequence is appended below the alignment. To switch
this off, use showConsensus=FALSE
. Whether or not sequence
names should be printed can be controlled via the
showNames
argument. The width reserved for the sequence
names can be adjusted using the nameWidth
argument;
the default is 20 like in the Biostrings method.
If the number of sequences in the alignment is large, output
can become quite lengthy. That is why only the first
halfNrow
and the last halfNrow
sequences are
shown. To show all sequences, set halfNrow
to NA
or -1. Note that print
can also handle masked objects,
where the masked sequences/positions are shown as hash marks.
However, the consensus sequences are computed from the
complete, unmasked alignment and displayed as such.
Additional arguments are passed on to
msaConsensusSequence
for customizing how the
consensus sequence is computed.
show(object)
:displays the alignment along with
metadata; synonymous to calling print
with default
arguments.
version(object)
:displays the algorithm with which
the multiple alignment has been computed along with its
version number (see also MsaMetaData
).
params(x)
:accessor to the params
slot (see
also MsaMetaData
)
Enrico Bonatesta, Christoph Horejs-Kainrath, and Ulrich Bodenhofer <msa@bioinf.jku.at>
http://www.bioinf.jku.at/software/msa
U. Bodenhofer, E. Bonatesta, C. Horejs-Kainrath, and S. Hochreiter (2015). msa: an R package for multiple sequence alignment. Bioinformatics 31(24):3997-3999. DOI: 10.1093/bioinformatics/btv494.
msa
, msaClustalW
,
msaClustalOmega
, msaMuscle
,
MsaMetaData
## read sequences filepath <- system.file("examples", "exampleAA.fasta", package="msa") mySeqs <- readAAStringSet(filepath) ## simple call with default values myAlignment <- msaClustalOmega(mySeqs) ## show the algorithm version with which the results were created version(myAlignment) ## show the results show(myAlignment) ## print the results print(myAlignment, show="alignment") print(myAlignment, show="alignment", showConsensus=FALSE) print(myAlignment, show=c("alignment", "version")) print(myAlignment, show="standardParams") print(myAlignment, show="algParams") print(myAlignment, show=c("call", "version")) ## print results with custom consensus sequence print(myAlignment, show="complete", type="upperlower", thresh=c(50, 20)) ## show the params params(myAlignment)