sequence-formats 1.9.0.0 → 1.10.0.0

raw patch · 7 files changed

+90/−93 lines, 7 filesPVP ok

version bump matches the API change (PVP)

API changes (from Hackage documentation)

Files

• Changelog.md +1/−0
• sequence-formats.cabal +1/−1
• src/SequenceFormats/Eigenstrat.hs +36/−42
• src/SequenceFormats/Plink.hs +34/−42
• src/SequenceFormats/Utils.hs +8/−2
• test/SequenceFormats/EigenstratSpec.hs +5/−3
• test/SequenceFormats/PlinkSpec.hs +5/−3

@@ -1,5 +1,6 @@ # Changelog +- V 1.10.0.0: Brought gzip-writing support for Eigenstrat and Plink files back to non-breaking API in `writeEigenstrat` and `writePlink`. Client code can safely update from 1.8.X to 1.10.0.0. - V 1.9.0.0: Added gzip-writing support for Eigenstrat and Plink files. This required a breaking change in `writeEigenstrat` and `writePlink`. - V 1.8.1.0: Added gzip-support (read-only for now) for Plink (bed and bim files) and VCF. - V 1.8.0.1: Allow reading arbitrary letters as reference base in Pileup format. Before only ACTG, N and M were allowed.

@@ -1,5 +1,5 @@ name:                sequence-formats-version:             1.9.0.0+version:             1.10.0.0 synopsis:            A package with basic parsing utilities for several Bioinformatic data formats. description:         Contains utilities to parse and write Eigenstrat, Fasta, FreqSum, VCF, Plink and other file formats used in population genetics analyses. license:             GPL-3

@@ -13,27 +13,30 @@ import           SequenceFormats.Utils            (Chrom (..),                                                    SeqFormatException (..),                                                    consumeProducer,+                                                   deflateFinaliser,+                                                   gzipConsumer,                                                    readFileProdCheckCompress,-                                                   word, gzipConsumer, writeFromPopper)+                                                   word)  import           Control.Applicative              ((<|>)) import           Control.Exception                (throw) import           Control.Monad                    (forM_, void) import           Control.Monad.Catch              (MonadThrow) import           Control.Monad.IO.Class           (MonadIO, liftIO)+import           Control.Monad.Trans.Class        (lift) import qualified Data.Attoparsec.ByteString.Char8 as A import qualified Data.ByteString.Char8            as B import           Data.List                        (isSuffixOf) import qualified Data.Streaming.Zlib              as Z import           Data.Vector                      (Vector, fromList, toList) import           Pipes                            (Consumer, Pipe, Producer,-                                                   cat, for, yield, (>->), runEffect)+                                                   cat, for, yield, (>->)) import qualified Pipes.ByteString                 as PB import qualified Pipes.Prelude                    as P-import           Pipes.Safe                       (MonadSafe)+import           Pipes.Safe                       (MonadSafe, register) import qualified Pipes.Safe.Prelude               as PS-import           System.IO                        (Handle, IOMode (..),-                                                   hPutStrLn, withFile)+import           System.IO                        (IOMode (..), hPutStrLn,+                                                   withFile)  -- |A datatype to represent a single genomic SNP. The constructor arguments are: -- Chromosome, Position, Reference Allele, Alternative Allele.@@ -164,31 +167,37 @@         Unknown -> "U"  -- |Function to write an Eigenstrat Snp File. Returns a consumer expecting EigenstratSnpEntries.-writeEigenstratSnp :: (MonadIO m) => Maybe Z.Deflate -- ^If Nothing, then no compression-    -> Handle -- ^The Eigenstrat Snp File Handle.+writeEigenstratSnp :: (MonadSafe m) => FilePath -- ^The Eigenstrat Snp File.     -> Consumer EigenstratSnpEntry m () -- ^A consumer to read EigenstratSnpEntries-writeEigenstratSnp maybeDeflate snpFileH =-    let snpOutTextConsumer = case maybeDeflate of-            Nothing -> PB.toHandle snpFileH-            Just def -> gzipConsumer def snpFileH-        toTextPipe = P.map (\(EigenstratSnpEntry chrom pos gpos gid ref alt) ->+writeEigenstratSnp snpFile = do+    (_, snpFileH) <- lift $ PS.openFile snpFile WriteMode+    snpOutTextConsumer <- if ".gz" `isSuffixOf` snpFile then do+            def <- liftIO $ Z.initDeflate 6 (Z.WindowBits 31)+            _ <- register (deflateFinaliser def snpFileH)+            return $ gzipConsumer def snpFileH+        else+            return $ PB.toHandle snpFileH+    let toTextPipe = P.map (\(EigenstratSnpEntry chrom pos gpos gid ref alt) ->             let snpLine = B.intercalate "\t" [gid, unChrom chrom, B.pack (show gpos),                     B.pack (show pos), B.singleton ref, B.singleton alt]             in  snpLine <> "\n")-    in  toTextPipe >-> snpOutTextConsumer+    toTextPipe >-> snpOutTextConsumer  -- |Function to write an Eigentrat Geno File. Returns a consumer expecting Eigenstrat Genolines.-writeEigenstratGeno :: (MonadIO m) => Maybe Z.Deflate -- ^If Nothing, then no compression-    -> Handle -- ^The Genotype file handle+writeEigenstratGeno :: (MonadSafe m) => FilePath -- ^The Genotype file     -> Consumer GenoLine m () -- ^A consumer to read Genotype entries.-writeEigenstratGeno maybeDeflate genoFileH =-    let genoOutTextConsumer = case maybeDeflate of-            Nothing -> PB.toHandle genoFileH-            Just def -> gzipConsumer def genoFileH-        toTextPipe = P.map (\genoLine ->+writeEigenstratGeno genoFile = do+    (_, genoFileH) <- lift $ PS.openFile genoFile WriteMode+    genoOutTextConsumer <- if ".gz" `isSuffixOf` genoFile then do+            def <- liftIO $ Z.initDeflate 6 (Z.WindowBits 31)+            _ <- register (deflateFinaliser def genoFileH)+            return $ gzipConsumer def genoFileH+        else+            return $ PB.toHandle genoFileH+    let toTextPipe = P.map (\genoLine ->             let genoLineStr = B.concat . map (B.pack . show . toEigenStratNum) . toList $ genoLine             in  genoLineStr <> "\n")-    in  toTextPipe >-> genoOutTextConsumer+    toTextPipe >-> genoOutTextConsumer   where     toEigenStratNum c = case c of         HomRef  -> 2 :: Int@@ -201,25 +210,10 @@                 -> FilePath -- ^The Snp File                 -> FilePath -- ^The Ind file                 -> [EigenstratIndEntry] -- ^The list of individual entries-                -> Producer (EigenstratSnpEntry, GenoLine) m () -- ^A consumer to read joint Snp/Genotype entries.-                -> m ()-writeEigenstrat genoFile snpFile indFile indEntries prod = do-    snpDeflate  <- if ".gz" `isSuffixOf` snpFile  then fmap Just . liftIO $ Z.initDeflate 6 (Z.WindowBits 31) else return Nothing-    genoDeflate <- if ".gz" `isSuffixOf` genoFile then fmap Just . liftIO $ Z.initDeflate 6 (Z.WindowBits 31) else return Nothing+                -> Consumer (EigenstratSnpEntry, GenoLine) m () -- ^A consumer to read joint Snp/Genotype entries.+writeEigenstrat genoFile snpFile indFile indEntries = do     liftIO $ writeEigenstratIndFile indFile indEntries-    (_, snpFileH)  <- PS.openFile snpFile  WriteMode-    (_, genoFileH) <- PS.openFile genoFile WriteMode-    let snpOutConsumer  = writeEigenstratSnp  snpDeflate  snpFileH-        genoOutConsumer = writeEigenstratGeno genoDeflate genoFileH-    runEffect $ prod >-> P.tee (P.map fst >-> snpOutConsumer) >-> P.map snd >-> genoOutConsumer-    case snpDeflate of-        Nothing -> return ()-        Just def -> do-            let finalPop = Z.finishDeflate def-            writeFromPopper finalPop snpFileH-    case genoDeflate of-        Nothing -> return ()-        Just def -> do-            let finalPop = Z.finishDeflate def-            writeFromPopper finalPop genoFileH-    +    let snpOutConsumer  = writeEigenstratSnp  snpFile+        genoOutConsumer = writeEigenstratGeno genoFile+    P.tee (P.map fst >-> snpOutConsumer) >-> P.map snd >-> genoOutConsumer+

@@ -17,13 +17,16 @@                                                    GenoEntry (..), GenoLine,                                                    Sex (..)) import           SequenceFormats.Utils            (Chrom (..), consumeProducer,+                                                   deflateFinaliser,+                                                   gzipConsumer,                                                    readFileProdCheckCompress,-                                                   word, gzipConsumer, writeFromPopper)+                                                   word, writeFromPopper)  import           Control.Applicative              ((<|>)) import           Control.Monad                    (forM_, void) import           Control.Monad.Catch              (MonadThrow, throwM) import           Control.Monad.IO.Class           (MonadIO, liftIO)+import           Control.Monad.Trans.Class        (lift) import           Control.Monad.Trans.State.Strict (runStateT) import qualified Data.Attoparsec.ByteString       as AB import qualified Data.Attoparsec.ByteString.Char8 as A@@ -34,14 +37,14 @@ import qualified Data.Streaming.Zlib              as Z import           Data.Vector                      (fromList, toList) import           Data.Word                        (Word8)-import           Pipes                            (Consumer, Producer, (>->), runEffect)+import           Pipes                            (Consumer, Producer, (>->)) import           Pipes.Attoparsec                 (ParsingError (..), parse) import qualified Pipes.ByteString                 as PB import qualified Pipes.Prelude                    as P-import           Pipes.Safe                       (MonadSafe)+import           Pipes.Safe                       (MonadSafe, register) import qualified Pipes.Safe.Prelude               as PS-import           System.IO                        (Handle, IOMode (..),-                                                   hPutStrLn, withFile)+import           System.IO                        (IOMode (..), hPutStrLn,+                                                   withFile)  -- see https://www.cog-genomics.org/plink/2.0/formats#fam data PlinkFamEntry = PlinkFamEntry {@@ -168,18 +171,21 @@     return (indEntries, P.zip snpProd genoProd)  -- |Function to write a Bim file. Returns a consumer expecting EigenstratSnpEntries.-writeBim :: (MonadIO m) => Maybe Z.Deflate -- ^If Nothing, then no Compression-    -> Handle -- ^The Eigenstrat Snp File handle.+writeBim :: (MonadSafe m) => FilePath -- ^The Plink Bim File.     -> Consumer EigenstratSnpEntry m () -- ^A consumer to read EigenstratSnpEntries-writeBim maybeDeflate snpFileH =-    let snpOutTextConsumer = case maybeDeflate of-            Nothing -> PB.toHandle snpFileH-            Just def -> gzipConsumer def snpFileH-        toTextPipe = P.map (\(EigenstratSnpEntry chrom pos gpos gid ref alt) ->+writeBim bimFile = do+    (_, bimFileH) <- lift $ PS.openFile bimFile WriteMode+    bimOutTextConsumer <- if ".gz" `isSuffixOf` bimFile then do+            def <- liftIO $ Z.initDeflate 6 (Z.WindowBits 31)+            _ <- register (deflateFinaliser def bimFileH)+            return $ gzipConsumer def bimFileH+        else+            return $ PB.toHandle bimFileH+    let toTextPipe = P.map (\(EigenstratSnpEntry chrom pos gpos gid ref alt) ->             let bimLine = B.intercalate "\t" [unChrom chrom, gid, B.pack (show gpos),                     B.pack (show pos), B.singleton ref, B.singleton alt]             in  bimLine <> "\n")-    in  toTextPipe >-> snpOutTextConsumer+    toTextPipe >-> bimOutTextConsumer  -- |Function to write a Plink Fam file. writeFam :: (MonadIO m) => FilePath -> [PlinkFamEntry] -> m ()@@ -193,20 +199,21 @@         Female  -> "2"         Unknown -> "0" --- |Function to write an Eigentrat Geno File. Returns a consumer expecting Eigenstrat Genolines.-writeBed :: (MonadIO m) => Maybe Z.Deflate -- ^If Nothing, then no Compression-  -> Handle -- ^The Bed file handle+-- |Function to write a Plink Bed File. Returns a consumer expecting Eigenstrat Genolines.+writeBed :: (MonadSafe m) => FilePath -- ^The Bed file handle   -> Consumer GenoLine m () -- ^A consumer to read Genotype entries.-writeBed maybeDeflate bedFileH = do+writeBed bedFile = do+    (_, bedFileH) <- lift $ PS.openFile bedFile WriteMode     let stickyBytes = BB.pack [0b01101100, 0b00011011, 0b00000001]-    bedOutConsumer <- case maybeDeflate of-        Nothing  -> do-            liftIO $ BB.hPut bedFileH stickyBytes-            return $ PB.toHandle bedFileH-        Just def -> do+    bedOutConsumer <- if ".gz" `isSuffixOf` bedFile then do+            def <- liftIO $ Z.initDeflate 6 (Z.WindowBits 31)+            _ <- register (deflateFinaliser def bedFileH)             pop <- liftIO (Z.feedDeflate def stickyBytes)             liftIO (writeFromPopper pop bedFileH)             return $ gzipConsumer def bedFileH+        else do+            liftIO $ BB.hPut bedFileH stickyBytes+            return $ PB.toHandle bedFileH     let toPlinkPipe = P.map (BB.pack . genoLineToBytes)     toPlinkPipe >-> bedOutConsumer   where@@ -233,24 +240,9 @@                 -> FilePath -- ^The Bim File                 -> FilePath -- ^The Fam file                 -> [PlinkFamEntry] -- ^The list of individual entries-                -> Producer (EigenstratSnpEntry, GenoLine) m () -- ^A consumer to read joint Snp/Genotype entries.-                -> m ()-writePlink bedFile bimFile famFile indEntries prod = do-    bimDeflate <- if ".gz" `isSuffixOf` bimFile then fmap Just . liftIO $ Z.initDeflate 6 (Z.WindowBits 31) else return Nothing-    bedDeflate <- if ".gz" `isSuffixOf` bedFile then fmap Just . liftIO $ Z.initDeflate 6 (Z.WindowBits 31) else return Nothing+                -> Consumer (EigenstratSnpEntry, GenoLine) m () -- ^A consumer to read joint Snp/Genotype entries.+writePlink bedFile bimFile famFile indEntries = do     liftIO $ writeFam famFile indEntries-    (_, bimFileH) <- PS.openFile bimFile WriteMode-    (_, bedFileH) <- PS.openFile bedFile WriteMode-    let bimOutConsumer = writeBim bimDeflate bimFileH-        bedOutConsumer = writeBed bedDeflate bedFileH-    runEffect $ prod >-> P.tee (P.map fst >-> bimOutConsumer) >-> P.map snd >-> bedOutConsumer-    case bimDeflate of-        Nothing -> return ()-        Just def -> do-            let finalPop = Z.finishDeflate def-            writeFromPopper finalPop bimFileH-    case bedDeflate of-        Nothing -> return ()-        Just def -> do-            let finalPop = Z.finishDeflate def-            writeFromPopper finalPop bedFileH+    let bimOutConsumer = writeBim bimFile+        bedOutConsumer = writeBed bedFile+    P.tee (P.map fst >-> bimOutConsumer) >-> P.map snd >-> bedOutConsumer

@@ -4,7 +4,7 @@  module SequenceFormats.Utils (liftParsingErrors,                               consumeProducer, readFileProd, readFileProdCheckCompress,-                              SeqFormatException(..),+                              SeqFormatException(..), deflateFinaliser,                               Chrom(..), word, gzipConsumer, writeFromPopper, Z.Deflate) where  import           Control.Error                    (readErr)@@ -17,7 +17,8 @@ import           Data.Char                        (isSpace) import           Data.List                        (isSuffixOf) import qualified Data.Streaming.Zlib              as Z-import           Pipes                            (Consumer, Producer, next, await)+import           Pipes                            (Consumer, Producer, await,+                                                   next) import           Pipes.Attoparsec                 (ParsingError (..), parsed) import qualified Pipes.ByteString                 as PB import           Pipes.GZip                       (decompress)@@ -103,3 +104,8 @@       Z.PRNext bs -> do          liftIO $ B.hPut h bs          writeFromPopper pop h++deflateFinaliser :: (MonadIO m) => Z.Deflate -> Handle -> m ()+deflateFinaliser def h = do+    let finalPop = liftIO $ Z.finishDeflate def+    writeFromPopper finalPop h

@@ -4,7 +4,7 @@ import           Control.Foldl              (list, purely) import           Control.Monad.IO.Class     (liftIO) import           Data.Vector                (fromList)-import           Pipes                      (each)+import           Pipes                      (each, runEffect, (>->)) import qualified Pipes.Prelude              as P import           Pipes.Safe                 (runSafeT) import           SequenceFormats.Eigenstrat (EigenstratIndEntry (..),@@ -83,7 +83,8 @@             testDatSnpProd = each mockDatEigenstratSnp             testDatGenoProd = each mockDatEigenstratGeno             testDatJointProd = P.zip testDatSnpProd testDatGenoProd-        runSafeT $ writeEigenstrat tmpGeno tmpSnp tmpInd mockDatEigenstratInd testDatJointProd+            cons = writeEigenstrat tmpGeno tmpSnp tmpInd mockDatEigenstratInd+        runSafeT . runEffect $ testDatJointProd >-> cons         (indEntries, esProd) <- liftIO . runSafeT $ readEigenstrat tmpGeno tmpSnp tmpInd         indEntries `shouldBe` mockDatEigenstratInd         snpGenoEntries <- liftIO . runSafeT $ purely P.fold list esProd@@ -99,7 +100,8 @@             testDatSnpProd = each mockDatEigenstratSnp             testDatGenoProd = each mockDatEigenstratGeno             testDatJointProd = P.zip testDatSnpProd testDatGenoProd-        runSafeT $ writeEigenstrat tmpGeno tmpSnp tmpInd mockDatEigenstratInd testDatJointProd+            cons = writeEigenstrat tmpGeno tmpSnp tmpInd mockDatEigenstratInd+        runSafeT . runEffect $ testDatJointProd >-> cons         (indEntries, esProd) <- liftIO . runSafeT $ readEigenstrat tmpGeno tmpSnp tmpInd         indEntries `shouldBe` mockDatEigenstratInd         snpGenoEntries <- liftIO . runSafeT $ purely P.fold list esProd

@@ -16,7 +16,7 @@ import           Control.Foldl              (list, purely) import           Control.Monad.IO.Class     (liftIO) import           Data.Vector                (fromList)-import           Pipes                      (each)+import           Pipes                      (each, runEffect, (>->)) import qualified Pipes.Prelude              as P import           Pipes.Safe                 (runSafeT) import           Test.Hspec@@ -118,7 +118,8 @@             testDatSnpProd = each mockDatEigenstratSnp             testDatGenoProd = each mockDatPlinkBed             testDatJointProd = P.zip testDatSnpProd testDatGenoProd-        runSafeT $ writePlink tmpGeno tmpSnp tmpInd mockDatPlinkFam testDatJointProd+            cons = writePlink tmpGeno tmpSnp tmpInd mockDatPlinkFam+        runSafeT . runEffect $ testDatJointProd >-> cons         (indEntries, esProd) <- liftIO . runSafeT $ readPlink tmpGeno tmpSnp tmpInd         indEntries `shouldBe` mockDatPlinkFam         snpGenoEntries <- liftIO . runSafeT $ purely P.fold list esProd@@ -134,7 +135,8 @@             testDatSnpProd = each mockDatEigenstratSnp             testDatGenoProd = each mockDatPlinkBed             testDatJointProd = P.zip testDatSnpProd testDatGenoProd-        runSafeT $ writePlink tmpBed tmpBim tmpFam mockDatPlinkFam testDatJointProd+            cons = writePlink tmpBed tmpBim tmpFam mockDatPlinkFam+        runSafeT . runEffect $ testDatJointProd >-> cons         (indEntries, esProd) <- liftIO . runSafeT $ readPlink tmpBed tmpBim tmpFam         indEntries `shouldBe` mockDatPlinkFam         snpGenoEntries <- liftIO . runSafeT $ purely P.fold list esProd