bio-0.3.5: A bioinformatics librarySource codeContentsIndex
Bio.Sequence
Contents
Data structures etc (Bio.Sequence.SeqData)
Accessor functions
Converting to and from String.
Nucleotide functionality.
Protein sequence functionality
File formats
The Fasta file format (Bio.Sequence.Fasta)
Quality data
The FastQ format (Bio.Sequence.FastQ)
The phd file format (Bio.Sequence.Phd)
TwoBit file format support (Bio.Seqeunce.TwoBit)
Hashing functionality (Bio.Sequence.HashWord)
Entropy calculations
Description

This is a meta-module importing and re-exporting sequence-related stuff.

It encompasses the Bio.Sequence.SeqData, Bio.Sequence.Fasta, and Bio.Sequence.TwoBit modules.

Synopsis
data Sequence = Seq !SeqData !SeqData !(Maybe QualData)
type Offset = Int64
type SeqData = ByteString
type Qual = Word8
type QualData = ByteString
seqlength :: Sequence -> Offset
seqlabel :: Sequence -> SeqData
seqheader :: Sequence -> SeqData
seqdata :: Sequence -> SeqData
seqqual :: Sequence -> QualData
(!) :: Sequence -> Offset -> Char
appendHeader :: Sequence -> String -> Sequence
setHeader :: Sequence -> String -> Sequence
fromStr :: String -> SeqData
toStr :: SeqData -> String
compl :: Char -> Char
revcompl :: Sequence -> Sequence
data Amino
= Ala
| Arg
| Asn
| Asp
| Cys
| Gln
| Glu
| Gly
| His
| Ile
| Leu
| Lys
| Met
| Phe
| Pro
| Ser
| Thr
| Tyr
| Trp
| Val
| STP
| Asx
| Glx
| Xle
| Xaa
translate :: Sequence -> Offset -> [Amino]
fromIUPAC :: SeqData -> [Amino]
toIUPAC :: [Amino] -> SeqData
readFasta :: FilePath -> IO [Sequence]
hReadFasta :: Handle -> IO [Sequence]
writeFasta :: FilePath -> [Sequence] -> IO ()
hWriteFasta :: Handle -> [Sequence] -> IO ()
readQual :: FilePath -> IO [Sequence]
writeQual :: FilePath -> [Sequence] -> IO ()
hWriteQual :: Handle -> [Sequence] -> IO ()
readFastaQual :: FilePath -> FilePath -> IO [Sequence]
writeFastaQual :: FilePath -> FilePath -> [Sequence] -> IO ()
hWriteFastaQual :: Handle -> Handle -> [Sequence] -> IO ()
readFastQ :: FilePath -> IO [Sequence]
writeFastQ :: FilePath -> [Sequence] -> IO ()
hReadFastQ :: Handle -> IO [Sequence]
hWriteFastQ :: Handle -> [Sequence] -> IO ()
readPhd :: FilePath -> IO Sequence
hReadPhd :: Handle -> IO Sequence
decode2Bit :: ByteString -> [Sequence]
read2Bit :: FilePath -> IO [Sequence]
hRead2Bit :: Handle -> IO [Sequence]
data HashF k = HF {
hash :: SeqData -> Offset -> Maybe k
hashes :: SeqData -> [(k, Offset)]
ksort :: [k] -> [k]
}
contigous :: Integral k => Int -> HashF k
rcontig :: Integral k => Int -> HashF k
rcpacked :: Integral k => Int -> HashF k
class KWords s where
kwords :: Int -> s -> [s]
entropy :: (Ord str, KWords str) => Int -> str -> Double
Data structures etc (Bio.Sequence.SeqData)
data Sequence Source
A sequence consists of a header, the sequence data itself, and optional quality data.
Constructors
Seq !SeqData !SeqData !(Maybe QualData)header and actual sequence
show/hide Instances
type Offset = Int64Source
An offset, index, or length of a SeqData
type SeqData = ByteStringSource
The basic data type used in Sequences
type Qual = Word8Source
Basic type for quality data. Range 0..255. Typical Phred output is in the range 6..50, with 20 as the line in the sand separating good from bad.
type QualData = ByteStringSource
Quality data is a Qual vector, currently implemented as a ByteString.
Accessor functions
seqlength :: Sequence -> OffsetSource
Return sequence length.
seqlabel :: Sequence -> SeqDataSource
Return sequence label (first word of header)
seqheader :: Sequence -> SeqDataSource
Return full header.
seqdata :: Sequence -> SeqDataSource
Return the sequence data.
seqqual :: Sequence -> QualDataSource
Return the quality data, or error if none exist. Use hasqual if in doubt.
(!) :: Sequence -> Offset -> CharSource
Read the character at the specified position in the sequence.
appendHeader :: Sequence -> String -> SequenceSource
setHeader :: Sequence -> String -> SequenceSource
Modify the header by appending text, or by replacing all but the sequence label (i.e. first word).
Converting to and from String.
fromStr :: String -> SeqDataSource
Convert a String to SeqData
toStr :: SeqData -> StringSource
Convert a SeqData to a String
Nucleotide functionality.
compl :: Char -> CharSource
Complement a single character. I.e. identify the nucleotide it can hybridize with. Note that for multiple nucleotides, you usually want the reverse complement (see revcompl for that).
revcompl :: Sequence -> SequenceSource
Calculate the reverse complement. This is only relevant for the nucleotide alphabet, and it leaves other characters unmodified.
Protein sequence functionality
data Amino Source
Constructors
Ala
Arg
Asn
Asp
Cys
Gln
Glu
Gly
His
Ile
Leu
Lys
Met
Phe
Pro
Ser
Thr
Tyr
Trp
Val
STP
Asx
Glx
Xle
Xaa
show/hide Instances
translate :: Sequence -> Offset -> [Amino]Source
Translate a nucleotide sequence into the corresponding protein sequence. This works rather blindly, with no attempt to identify ORFs or otherwise QA the result.
fromIUPAC :: SeqData -> [Amino]Source
Convert a sequence in IUPAC format to a list of amino acids.
toIUPAC :: [Amino] -> SeqDataSource
Convert a list of amino acids to a sequence in IUPAC format.
File formats
The Fasta file format (Bio.Sequence.Fasta)
readFasta :: FilePath -> IO [Sequence]Source
Lazily read sequences from a FASTA-formatted file
hReadFasta :: Handle -> IO [Sequence]Source
Lazily read sequence from handle
writeFasta :: FilePath -> [Sequence] -> IO ()Source
Write sequences to a FASTA-formatted file. Line length is 60.
hWriteFasta :: Handle -> [Sequence] -> IO ()Source
Write sequences in FASTA format to a handle.
Quality data
Not part of the Fasta format, and treated separately.
readQual :: FilePath -> IO [Sequence]Source
Read quality data for sequences to a file.
writeQual :: FilePath -> [Sequence] -> IO ()Source
Write quality data for sequences to a file.
hWriteQual :: Handle -> [Sequence] -> IO ()Source
readFastaQual :: FilePath -> FilePath -> IO [Sequence]Source
Read sequence and associated quality. Will error if the sequences and qualites do not match one-to-one in sequence.
writeFastaQual :: FilePath -> FilePath -> [Sequence] -> IO ()Source
Write sequence and quality data simulatnously This may be more laziness-friendly.
hWriteFastaQual :: Handle -> Handle -> [Sequence] -> IO ()Source
The FastQ format (Bio.Sequence.FastQ)
readFastQ :: FilePath -> IO [Sequence]Source
writeFastQ :: FilePath -> [Sequence] -> IO ()Source
hReadFastQ :: Handle -> IO [Sequence]Source
hWriteFastQ :: Handle -> [Sequence] -> IO ()Source
The phd file format (Bio.Sequence.Phd)
These contain base (nucleotide) calling information, and are generated by phred.
readPhd :: FilePath -> IO SequenceSource
Parse a .phd file, extracting the contents as a Sequence
hReadPhd :: Handle -> IO SequenceSource
Parse .phd contents from a handle
TwoBit file format support (Bio.Seqeunce.TwoBit)
Used by BLAT and related tools.
decode2Bit :: ByteString -> [Sequence]Source
Parse a (lazy) ByteString as sequences in the 2bit format.
read2Bit :: FilePath -> IO [Sequence]Source
Extract sequences from a file in 2bit format.
hRead2Bit :: Handle -> IO [Sequence]Source
Extract sequences in the 2bit format from a handle.
Hashing functionality (Bio.Sequence.HashWord)
Packing words from sequences into integral data types
data HashF k Source
This is a struct for containing a set of hashing functions
Constructors
HF
hash :: SeqData -> Offset -> Maybe kcalculates the hash at a given offset in the sequence
hashes :: SeqData -> [(k, Offset)]calculate all hashes from a sequence, and their indices
ksort :: [k] -> [k]for sorting hashes
contigous :: Integral k => Int -> HashF kSource
Contigous constructs an int/eger from a contigous k-word.
rcontig :: Integral k => Int -> HashF kSource
Like contigous, but returns the same hash for a word and its reverse complement.
rcpacked :: Integral k => Int -> HashF kSource
Like rcontig, but ignoring monomers (i.e. arbitrarily long runs of a single nucelotide are treated the same a single nucleotide.
Entropy calculations
class KWords s whereSource
Methods
kwords :: Int -> s -> [s]Source
show/hide Instances
KWords ([] a)
entropy :: (Ord str, KWords str) => Int -> str -> DoubleSource
Produced by Haddock version 2.4.2