packages feed

bioinformatics-toolkit 0.2.4 → 0.3.0

raw patch · 15 files changed

+124/−348 lines, 15 filesdep +HsHTSLibdep −binarydep −colourdep −optparse-applicativedep ~basedep ~bytestringdep ~textbinary-added

Dependencies added: HsHTSLib

Dependencies removed: binary, colour, optparse-applicative, palette, shelly

Dependency ranges changed: base, bytestring, text

Files

− app/mergeMotifs.hs
@@ -1,229 +0,0 @@-{-# LANGUAGE OverloadedStrings #-}--import           AI.Clustering.Hierarchical        hiding (drawDendrogram)-import           Control.Monad                     (forM, forM_)-import qualified Data.ByteString.Char8             as B-import           Data.Default.Class-import           Data.Double.Conversion.ByteString (toShortest)-import           Data.List-import           Data.List.Split                   (splitOn)-import           Data.Ord-{--import           Diagrams.Backend.Cairo-import           Diagrams.Plots.Dendrogram-import           Diagrams.Prelude           (dims2D, strutX, (|||))--}-import           Options.Applicative-import           System.IO-import           Text.Printf--import           Bio.Data.Fasta-import           Bio.Motif-import           Bio.Motif.Alignment-import           Bio.Motif.Merge-import           Bio.Seq                           (toBS)-import           Bio.Utils.Functions--data CMD = Merge { mergeInput  :: FilePath-                 , mode        :: String-                 , thres       :: Double-                 , alignOption :: AlignOption-                 , mergeOutput :: FilePath-                 }-         | Dist { inputA      :: FilePath-                , inputB      :: Maybe FilePath-                , alignOption :: AlignOption-                }-         | Cluster { clusterInput :: !FilePath-                   , cutoff :: !Double-                   , alignOption :: !AlignOption-                   }--mergeParser :: Parser CMD-mergeParser = Merge-    <$> argument str (metavar "INPUT")-    <*> strOption-        ( long "mode"-       <> short 'm'-       <> value "iter"-       <> metavar "MODE"-       <> help "Merging algorithm, could be iter or tree, default is iter" )-    <*> option auto-        ( long "thres"-       <> short 't'-       <> value 0.2-       <> metavar "THRESHOLD"-       <> help "two motifs that have distance belowing threshold would be merged, default is 0.2" )-    <*> alignParser-    <*> strOption-        ( long "output"-       <> short 'o'-       <> value "merged_output.meme"-       <> metavar "OUTPUT" )--distParser :: Parser CMD-distParser = Dist-    <$> strOption-           ( short 'a'-          <> metavar "INPUT_A" )-    <*> (optional . strOption)-           ( short 'b'-          <> metavar "INPUT_B" )-    <*> alignParser--clusterParser :: Parser CMD-clusterParser = Cluster-    <$> argument str (metavar "INPUT")-    <*> option auto-        ( long "height"-       <> short 'h'-       <> value 0.2-       <> metavar "HEIGHT"-       <> help "Cut hierarchical tree at given height. Default: 0.2" )-    <*> alignParser--data AlignOption = AlignOption-    { gap     :: Double-    , gapMode :: String-    , avgMode :: String-    }--alignParser :: Parser AlignOption-alignParser = AlignOption-     <$> option auto-           ( long "gap"-          <> short 'g'-          <> value 0.05-          <> metavar "GAP_PENALTY"-          <> help "Gap penalty, default: 0.05" )-     <*> strOption-           ( long "gap_mode"-          <> value "exp"-          <> metavar "GAP_MODE"-          <> help "Gap penalty mode, one of linear, quadratic, cubic, and exp. default: exp." )-     <*> strOption-           ( long "avg_mode"-          <> value "l1"-          <> metavar "AVERAGE_MODE"-          <> help "Averaging function, one of l1, l2, l3, max. default: l1." )--treeMerge :: Double -> String -> [Motif] -> AlignOption-          -> ([Motif], Dendrogram Motif)-treeMerge th pre ms alignOpt = (zipWith f [0::Int ..] $ map merge $ tree `cutAt` th, tree)-  where-    f i (suffix, pwm) = Motif ((B.pack $ pre ++ "_" ++ show i ++ "_" ++ show (toIUPAC pwm))-                                 `B.append` "("-                                 `B.append` suffix-                                 `B.append` ")" ) pwm-    merge tr = ( B.intercalate "+" $ map _name $ flatten tr-               , dilute $ mergeTreeWeighted align tr)-    tree = buildTree align ms-    align = mkAlignFn alignOpt-{-# INLINE treeMerge #-}--getSuffix :: String -> String-getSuffix = last . splitOn "."-{-# INLINE getSuffix #-}--mkAlignFn :: AlignOption -> AlignFn-mkAlignFn (AlignOption gap gapMode avgMode) = alignmentBy jsd (pFn gap) avgFn-  where-    pFn = case gapMode of-        "linear" -> linPenal-        "quadratic" -> quadPenal-        "cubic" -> cubPenal-        "exp" -> expPenal-        _ -> error "Unknown gap mode"-    avgFn = case avgMode of-        "l1" -> l1-        "l2" -> l2-        "l3" -> l3-        "max" -> lInf-        _ -> error "Unknown average mode"-{-# INLINE mkAlignFn #-}--readMotif :: FilePath -> IO [Motif]-readMotif fl = case getSuffix fl of-    "fasta" -> readFasta' fl-    _ -> readMEME fl--writeMotif :: FilePath -> [Motif] -> IO ()-writeMotif fl motifs = case getSuffix fl of-    "fasta" -> writeFasta fl motifs-    _ -> writeMEME fl motifs def---defaultMain :: CMD -> IO ()-defaultMain (Dist a b alignOpt) = do-    motifsA <- readMotif a-    pairs <- case b of-        Nothing -> return $ comb motifsA-        Just b' -> do-            motifsB <- readMotif b'-            return [(x,y) | x <- motifsA, y <- motifsB]-    forM_ pairs $ \(x,y) -> do-        let (d,_) = alignFn (_pwm x) $ (_pwm y)-        B.putStrLn $ B.intercalate "\t" [_name x, _name y, toShortest d]-  where-    alignFn = mkAlignFn alignOpt--defaultMain (Merge inFl m th alignOpt outFl)= do-    motifs <- readMotif inFl-    let motifNumber = length motifs--    hPutStrLn stderr $ printf "Merging Mode: %s" m-    hPutStrLn stderr $ printf "Read %d motifs" motifNumber--    motifs' <- case m of-        "tree" -> do-            let (newMotifs, tree) = treeMerge th "" motifs alignOpt-                fn x = B.unpack (_name x) ++ ": " ++ B.unpack (toBS $ toIUPAC $ _pwm x)--            {--            case svg of-                Just fl -> do-                    let w = 80-                        h = 5 * fromIntegral motifNumber-                    renderCairo fl (dims2D (10*w) (10*h)) $ drawDendrogram w h th tree fn ||| strutX 40-                    return newMotifs-                    -}-            return newMotifs-        "iter" -> do-            let rs = iterativeMerge (mkAlignFn alignOpt) th motifs-            forM rs $ \(nm, pwm, ws) -> do-                let pwm' = dilute (pwm, ws)-                return $ Motif (B.intercalate "+" nm) pwm-         -- _ -> error "Unkown mode!"--    hPutStrLn stderr $ printf "Write %d motifs" (length motifs')--    writeMotif outFl motifs'--defaultMain (Cluster inFl c alignOpt) = do-    motifs <- readMotif inFl-    let tree = buildTree align motifs-        align = mkAlignFn alignOpt-    forM_ (tree `cutAt` c) $ \t ->-        B.putStrLn $ B.intercalate "\t" $ map _name $ flatten t-{-# INLINE defaultMain #-}--comb :: [a] -> [(a,a)]-comb (y:ys) = zip (repeat y) ys ++ comb ys-comb _ = []-{-# INLINE comb #-}--main :: IO ()-main = execParser opts >>= defaultMain-  where-    opts = info (helper <*> parser)-            ( fullDesc-           <> header (printf "Compare, align and merge motifs, version %s" v))-    v = "0.2.0" :: String-    parser = subparser $ (-        command "merge" (info (helper <*> mergeParser) $-            fullDesc <> progDesc "Merge motifs")-     <> command "dist" (info (helper <*> distParser) $-            fullDesc <> progDesc "Align and compute pairwise distances")-     <> command "cluster" (info (helper <*> clusterParser) $-            fullDesc <> progDesc "Perform hierarchical clustering on motifs")-     )
− app/mkindex.hs
@@ -1,11 +0,0 @@-module Main where--import Bio.Seq.IO (mkIndex)-import Shelly-import qualified Data.Text as T-import System.Environment--main :: IO ()-main = do [inDir, outF] <- getArgs-          fls <- shelly . lsT . fromText . T.pack $ inDir-          mkIndex (map T.unpack fls) outF
− app/viewSeq.hs
@@ -1,15 +0,0 @@-module Main where--import qualified Data.ByteString.Char8 as B-import Bio.Seq-import Bio.Seq.IO-import System.Environment--main :: IO ()-main = do-    [fl, chr, start, end] <- getArgs-    withGenome fl $ \g -> do-        s <- getSeq g (B.pack chr, read start, read end) :: IO (Either String (DNA IUPAC))-        case s of-            Left err -> error err-            Right x -> B.putStrLn . toBS $ x
bioinformatics-toolkit.cabal view
@@ -1,8 +1,5 @@--- Initial bioinformatics-toolkit.cabal generated by cabal init.  For--- further documentation, see http://haskell.org/cabal/users-guide/- name:                bioinformatics-toolkit-version:             0.2.4+version:             0.3.0 synopsis:            A collection of bioinformatics tools description:         A collection of bioinformatics tools license:             MIT@@ -17,7 +14,7 @@ data-files:   tests/data/*.bed   tests/data/*.sorted.bed---  tests/data/*.bam+  tests/data/*.bam --  tests/data/*.bam.bai   tests/data/*.meme   tests/data/*.fasta@@ -30,8 +27,7 @@     Bio.ChIPSeq     Bio.ChIPSeq.FragLen     Bio.Data.Bed---    Bio.Data.BigWig---    Bio.Data.Bam+    Bio.Data.Bam     Bio.Data.Fasta     Bio.GO     Bio.GO.Parser@@ -53,79 +49,38 @@     Bio.Utils.Overlap     Bio.Utils.Types -  -- other-modules:-  -- other-extensions:   build-depends:       base >=4.8 && <5.0     , aeson     , aeson-pretty---    , bbi-    , binary     , bytestring >=0.10     , bytestring-lexing >=0.5     , case-insensitive     , clustering-    , colour     , conduit-combinators     , containers >=0.5     , data-default-class     , double-conversion+    , HsHTSLib     , http-conduit >=2.1.8     , hexpat+    , IntervalMap >=0.5.0.0     , matrices >=0.4.3     , mtl >=2.1.3.1     , math-functions     , parallel >=3.2     , primitive-    , palette     , split     , statistics >=0.13.2.1     , text >=0.11     , transformers >=0.3.0.0     , unordered-containers >=0.2+    , word8     , vector     , vector-algorithms-    , word8-    , IntervalMap >=0.5.0.0    default-language:    Haskell2010 -executable mkindex-  hs-source-dirs:   app-  main-is: mkindex.hs-  build-depends:-      base >=4.8 && <5.0-    , bioinformatics-toolkit-    , shelly-    , text-  default-language:    Haskell2010--executable viewSeq-  hs-source-dirs:   app-  main-is: viewSeq.hs-  build-depends:-      base >=4.8 && <5.0-    , bioinformatics-toolkit-    , bytestring-  default-language:    Haskell2010--executable mergeMotifs-  hs-source-dirs:   app-  main-is: mergeMotifs.hs-  build-depends:-      base >=4.6 && <5.0-    , bioinformatics-toolkit-    , bytestring-    , clustering-    , data-default-class-    , double-conversion-    , split-    , optparse-applicative---    , haskell-plot >=0.2.0.0---    , diagrams-cairo >=1.3---    , diagrams-lib >=1.3-  default-language:    Haskell2010- benchmark bench   type: exitcode-stdio-1.0   main-is: benchmarks/bench.hs@@ -148,6 +103,7 @@   main-is: test.hs   other-modules:       Tests.Bed+    , Tests.Bam     , Tests.Motif     , Tests.Seq     , Tests.GREAT
+ src/Bio/Data/Bam.hs view
@@ -0,0 +1,45 @@+{-# LANGUAGE FlexibleContexts #-}+module Bio.Data.Bam+    ( Bam+    , HeaderState+    , runBam+    , readBam+    , writeBam+    , bamToBed+    ) where++import           Bio.Data.Bed+import           Bio.HTS+import           Bio.HTS.Types             (Bam, FileHeader (..))+import           Conduit+import           Control.Monad.State (get, lift)++bamToBed :: Conduit Bam HeaderState BED+bamToBed = mapMC f =$= concatC+  where+    f bam = do+        BamHeader hdr <- lift get+        case getChr hdr bam of+            Just chr ->+                let start = fromIntegral $ position bam+                    end = fromIntegral $ endPos bam+                    nm = Just $ qName bam+                    strand = Just $ not $ isRev bam+                in return $ Just $ BED chr start end nm Nothing strand+            _ -> return Nothing+{-# INLINE bamToBed #-}++{-+viewBam :: IdxHandle -> (B.ByteString, Int, Int) -> Source IO Bam1+viewBam handle (chr, s, e) = case lookupTarget (idxHeader handle) chr of+    Nothing -> return ()+    Just chrId -> do+        q <- lift $ query handle chrId (fromIntegral s,fromIntegral e)+        go q+  where+    go q' = do r <- lift $ next q'+               case r of+                   Nothing -> return ()+                   Just bam -> yield bam >> go q'+{-# INLINE viewBam #-}+-}
src/Bio/Motif.hs view
@@ -63,31 +63,32 @@ size :: PWM -> Int size (PWM _ mat) = M.rows mat --- | information content of a poistion in pwm+-- | Information content of a poistion in pwm. (Not implemented) ic :: PWM -> Int -> Double ic = undefined --- | extract sub-PWM given starting position and length, zero indexed+-- | Extract sub-PWM given starting position and length, zero indexed. subPWM :: Int -> Int -> PWM -> PWM subPWM i l (PWM n mat) = PWM n $ M.subMatrix (i,0) (i+l,3) mat {-# INLINE subPWM #-} --- | reverse complementary of PWM+-- | Reverse complementary of PWM. rcPWM :: PWM -> PWM rcPWM (PWM n mat) = PWM n . M.fromVector d . U.reverse . M.flatten $ mat   where     d = M.dim mat {-# INLINE rcPWM #-} --- | GC content of PWM+-- | GC content of PWM. gcContentPWM :: PWM -> Double gcContentPWM (PWM _ mat) = loop 0 0 / fromIntegral m   where     m = M.rows mat-    loop !acc !i | i >= m = acc-                 | otherwise =-                     let acc' = acc + M.unsafeIndex mat (i,1) + M.unsafeIndex mat (i,2)-                     in loop acc' (i+1)+    loop !acc !i+        | i >= m = acc+        | otherwise =+            let acc' = acc + M.unsafeIndex mat (i,1) + M.unsafeIndex mat (i,2)+            in loop acc' (i+1)  data Motif = Motif     { _name :: !B.ByteString@@ -95,13 +96,13 @@     } deriving (Show, Read)  -- | background model which consists of single nucletide frequencies, and di-nucletide--- frequencies+-- frequencies. newtype Bkgd = BG (Double, Double, Double, Double)  instance Default Bkgd where     def = BG (0.25, 0.25, 0.25, 0.25) --- | convert pwm to consensus sequence, see D. R. Cavener (1987).+-- | Convert pwm to consensus sequence, see D. R. Cavener (1987). toIUPAC :: PWM -> DNA IUPAC toIUPAC (PWM _ pwm) = fromBS . B.pack . map f $ M.toRows pwm   where@@ -121,7 +122,7 @@     sort' (x, y) | x > y = (y, x)                  | otherwise = (x, y) --- | get scores of a long sequences at each position+-- | Get scores of a long sequences at each position. scores :: Bkgd -> PWM -> DNA a -> [Double] scores bg p@(PWM _ pwm) dna = go $! toBS dna   where@@ -130,7 +131,7 @@     len = M.rows pwm {-# INLINE scores #-} --- | a streaming version of scores+-- | A streaming version of scores. scores' :: Monad m => Bkgd -> PWM -> DNA a -> Source m Double scores' bg p@(PWM _ pwm) dna = go 0   where@@ -146,7 +147,7 @@ score bg p dna = scoreHelp bg p $! toBS dna {-# INLINE score #-} --- | the best possible score for a pwm+-- | The best possible matching score of a pwm. optimalScore :: Bkgd -> PWM -> Double optimalScore (BG (a, c, g, t)) (PWM _ pwm) = foldl' (+) 0 . map f . M.toRows $ pwm   where f xs = let (i, s) = U.maximumBy (comparing snd) .@@ -190,7 +191,8 @@         pseudoCount = 0.0001 {-# INLINE scoreHelp #-} --- | the probability that a kmer is generated by background (0-orderd model)+-- | The probability that a kmer is generated by background based on a+-- 0-orderd Markov model. pBkgd :: Bkgd -> B.ByteString -> Double pBkgd (BG (a, c, g, t)) = B.foldl' f 1   where@@ -204,14 +206,20 @@               _ -> undefined {-# INLINE pBkgd #-} +-- | calculate the minimum motif mathching score that would produce a kmer with+-- p-Value less than the given number. This score would then be used to search+-- for motif occurrences with significant p-Value+pValueToScore :: Double -> Bkgd -> PWM -> Double+pValueToScore p bg pwm = cdf' (scoreCDF bg pwm) $ 1 - p --- | unlike pValueToScore, this version compute the exact score but much slower--- and is inpractical for long motifs+-- | Unlike pValueToScore, this version compute the exact score but much slower+-- and is inpractical for long motifs. pValueToScoreExact :: Double -- ^ desirable p-Value               -> Bkgd               -> PWM               -> Double-pValueToScoreExact p bg pwm = go 0 0 . sort' . map ((scoreHelp bg pwm &&& pBkgd bg) . B.pack) . replicateM l $ "ACGT"+pValueToScoreExact p bg pwm = go 0 0 $ sort' $+    map ((scoreHelp bg pwm &&& pBkgd bg) . B.pack) $ replicateM l "ACGT"   where     sort' xs = U.create $ do v <- U.unsafeThaw . U.fromList $ xs                              I.sortBy (flip (comparing fst)) v@@ -220,13 +228,7 @@                    | otherwise = go (acc + snd (vec U.! i)) (i+1) vec     l = size pwm --- | calculate the minimum motif mathching score that would produce a kmer with--- p-Value less than the given number. This score would then be used to search--- for motif occurrences with significant p-Value-pValueToScore :: Double -> Bkgd -> PWM -> Double-pValueToScore p bg pwm = cdf' (scoreCDF bg pwm) $ 1 - p---- | The cumulative distribution function in the form of (x, P(X <= x))+-- | The cumulative distribution function in the form of (x, P(X <= x)). newtype CDF = CDF (U.Vector (Double, Double)) deriving (Read, Show)  -- P(X <= x)@@ -245,7 +247,7 @@     n = U.length v {-# INLINE cdf #-} --- the inverse of cdf+-- | The inverse of cdf. cdf' :: CDF -> Double -> Double cdf' (CDF v) p     | p > 1 || p < 0 = error "p must be in [0,1]"@@ -268,7 +270,7 @@ truncateCDF x (CDF v) = CDF $ U.filter ((>=x) . snd) v {-# INLINE truncateCDF #-} --- approximate the cdf of motif matching scores+-- | Approximate the cdf of motif matching scores scoreCDF :: Bkgd -> PWM -> CDF scoreCDF (BG (a,c,g,t)) pwm = toCDF $ loop (U.singleton 1, const 0) 0   where@@ -319,11 +321,11 @@            | otherwise = log x {-# INLINE scoreCDF #-} --- | get pwm from a matrix+-- | Get pwm from a matrix. toPWM :: [B.ByteString] -> PWM toPWM x = PWM Nothing . M.fromLists . map (map readDouble.B.words) $ x --- | pwm to bytestring+-- | Convert pwm to bytestring. fromPWM :: PWM -> B.ByteString fromPWM = B.unlines . map (B.unwords . map toShortest) . M.toLists . _mat 
src/Bio/Motif/Merge.hs view
@@ -9,6 +9,7 @@     , mergeTreeWeighted     , iterativeMerge     , buildTree+    , cutTreeBy     )where  import           AI.Clustering.Hierarchical    hiding (size)
src/Bio/RealWorld/BioGRID.hs view
@@ -78,7 +78,7 @@ -- | retreive first 10,000 records fetchByGeneNames :: [String] -> IO [TAB2] fetchByGeneNames genes = do-    initReq <- parseUrl $ intercalate "&" [url, geneList, tax, accessKey]+    initReq <- parseRequest $ intercalate "&" [url, geneList, tax, accessKey]     let request = initReq { method = "GET"                           , requestHeaders = [("Content-type", "text/plain")]                           }
src/Bio/RealWorld/ENCODE.hs view
@@ -78,7 +78,7 @@ -- search ["chip", "sp1"] ["type=experiment"] search :: KeyWords -> IO (Either String [Value]) search kw = do-    initReq <- parseUrl url+    initReq <- parseRequest url     let request = initReq { method = "GET"                           , requestHeaders = [("accept", "application/json")]                           }@@ -116,7 +116,7 @@  openUrl :: String -> String -> IO B.ByteString openUrl url datatype = do-    initReq <- parseUrl url+    initReq <- parseRequest url     let request = initReq { method = "GET"                           , requestHeaders = [("accept", BS.pack datatype)]                           }
src/Bio/RealWorld/Ensembl.hs view
@@ -27,7 +27,7 @@  lookupHelp :: [EnsemblID] -> IO (Either String Object) lookupHelp xs = do-    initReq <- parseUrl url+    initReq <- parseRequest url     let request = initReq { method = "POST"                           , requestHeaders = [("Content-type", "application/json")]                           , requestBody = body
src/Bio/Seq.hs view
@@ -1,6 +1,6 @@+{-# LANGUAGE BangPatterns          #-}+{-# LANGUAGE FlexibleContexts      #-} {-# LANGUAGE MultiParamTypeClasses #-}-{-# LANGUAGE BangPatterns #-}-{-# LANGUAGE FlexibleContexts #-} module Bio.Seq     (     -- * Alphabet@@ -19,11 +19,11 @@     , nucleotideFreq     ) where -import Prelude hiding (length) import qualified Data.ByteString.Char8 as B-import qualified Data.HashMap.Strict as M-import qualified Data.HashSet as S-import Data.Char8 (toUpper)+import           Data.Char8            (toUpper)+import qualified Data.HashMap.Strict   as M+import qualified Data.HashSet          as S+import           Prelude               hiding (length)  -- | Alphabet defined by http://www.chem.qmul.ac.uk/iupac/ -- | Standard unambiguous alphabet@@ -103,7 +103,7 @@         f x | x `S.member` alphabet (undefined :: RNA Basic) = x             | otherwise = error $ "Bio.Seq.fromBS: unknown character: " ++ [x] --- | O(n) Reverse complementary of DNA sequence+-- | O(n) Reverse complementary of DNA sequence. rc :: DNA alphabet -> DNA alphabet rc (DNA s) = DNA . B.map f . B.reverse $ s   where@@ -114,7 +114,7 @@         'T' -> 'A'         _ -> x --- | O(n) Compute GC content+-- | O(n) Compute GC content. gcContent :: DNA alphabet -> Double gcContent = (\(a,b) -> a / fromIntegral b) . B.foldl' f (0.0,0::Int) . toBS   where@@ -133,7 +133,7 @@                 _ -> x + 0.5     -- "NMKYR"         in (x', n+1) --- | O(n) Compute single nucleotide frequency+-- | O(n) Compute single nucleotide frequency. nucleotideFreq :: BioSeq DNA a => DNA a -> M.HashMap Char Int nucleotideFreq dna = B.foldl' f m0 . toBS $ dna   where
src/Bio/Seq/IO.hs view
@@ -57,8 +57,11 @@ withGenome fl fn = bracket (openGenome fl) closeGenome fn {-# INLINE withGenome #-} -type Query = (B.ByteString, Int, Int) -- (chr, start, end), zero-based index, half-close-half-open+-- | A query is represented by a tuple: (chr, start, end) and is+-- zero-based index, half-close-half-open+type Query = (B.ByteString, Int, Int) +-- | Retrieve sequence. getSeq :: BioSeq s a => Genome -> Query -> IO (Either String (s a)) getSeq (Genome h index headerSize) (chr, start, end) = case M.lookup chr index of     Just (chrStart, chrSize) ->@@ -71,6 +74,7 @@     _ -> return $ Left $ "Bio.Seq.getSeq: Cannot find " ++ show chr {-# INLINE getSeq #-} +-- | Retrieve whole chromosome. getChrom :: Genome -> B.ByteString -> IO (Either String (DNA IUPAC)) getChrom g chr = case lookup chr chrSize of     Just s -> getSeq g (chr, 0, s)@@ -79,6 +83,7 @@     chrSize = getChrSizes g {-# INLINE getChrom #-} +-- | Retrieve chromosome size information. getChrSizes :: Genome -> [(B.ByteString, Int)] getChrSizes (Genome _ table _) = map (\(k, (_, l)) -> (k, l)) . M.toList $ table {-# INLINE getChrSizes #-}
+ tests/Tests/Bam.hs view
@@ -0,0 +1,20 @@+{-# LANGUAGE OverloadedStrings #-}++module Tests.Bam (tests) where++import           Bio.Data.Bam+import           Bio.Data.Bed+import           Conduit+import           Test.Tasty+import           Test.Tasty.HUnit++tests :: TestTree+tests = testGroup "Test: Bio.Data.Bam"+    [ testCase "bamToBed" bamToBedTest+    ]++bamToBedTest :: Assertion+bamToBedTest = do+    bed <- readBed' "tests/data/example.bed"+    bed' <- runBam $ readBam "tests/data/example.bam" =$= bamToBed $$ sinkList+    bed @=? bed'
+ tests/data/example.bam view

binary file changed (absent → 689964 bytes)

tests/test.hs view
@@ -1,4 +1,5 @@ import qualified Tests.Bed as Bed+import qualified Tests.Bam as Bam import qualified Tests.Motif as Motif import qualified Tests.Seq as Seq import qualified Tests.GREAT as GREAT@@ -8,6 +9,7 @@ main :: IO () main = defaultMain $ testGroup "Main"     [ Bed.tests+    , Bam.tests     , Seq.tests     , Motif.tests     , GREAT.tests