cobot-io 0.1.3.2 → 0.1.3.3

raw patch · 8 files changed

+293/−196 lines, 8 filesdep ~directorydep ~neat-interpolationdep ~splitPVP: minor bump suggested

API additions: PVP suggests at least a minor version bump

Dependency ranges changed: directory, neat-interpolation, split, vector

API changes (from Hackage documentation)

Files

• ChangeLog.md +7/−0
• cobot-io.cabal +130/−126
• src/Bio/ABI/Clean.hs +11/−6
• src/Bio/ABI/Decode.hs +59/−14
• src/Bio/Sequence/Basecalled/Type.hs +19/−1
• test/ABISpec.hs +16/−4
• test/PDBSpec.hs +45/−39
• test/StructureSpec.hs +6/−6

@@ -2,6 +2,13 @@  ## [Unreleased] +## [0.1.3.3] - 2020-05-07+### Added+- Type and decoder for `ab1` with raw channel data and peak locations.+### Fixed+- Use supplied thresholds instead of default ones in the implementation of `cleanWith`;+- Do not run `IO` in `hspec` `Spec` monad, use `beforeAll`.+ ## [0.1.3.2] - 2020-04-15 ### Changed - Resolver version up.

@@ -1,136 +1,140 @@ cabal-version: 1.12-name: cobot-io-version: 0.1.3.2-license: BSD3-license-file: LICENSE-copyright: 2020, Biocad-maintainer: pavel@yakovlev.me-author: Pavel Yakovlev, Bogdan Neterebskii, Alexander Sadovnikov-homepage: https://github.com/biocad/cobot-io#readme-bug-reports: https://github.com/biocad/cobot-io/issues-synopsis: Biological data file formats and IO-description:-    Please see the README on GitHub at <https://github.com/biocad/cobot-io#readme>-category: Bio-build-type: Simple++-- This file has been generated from package.yaml by hpack version 0.33.0.+--+-- see: https://github.com/sol/hpack+--+-- hash: 476c2b78718477cf61d3cc8b8f752cbb5b4533d9a349a7aa05dc96fa352f2d7e++name:           cobot-io+version:        0.1.3.3+synopsis:       Biological data file formats and IO+description:    Please see the README on GitHub at <https://github.com/biocad/cobot-io#readme>+category:       Bio+homepage:       https://github.com/biocad/cobot-io#readme+bug-reports:    https://github.com/biocad/cobot-io/issues+author:         Pavel Yakovlev, Bogdan Neterebskii, Alexander Sadovnikov+maintainer:     pavel@yakovlev.me+copyright:      2020, Biocad+license:        BSD3+license-file:   LICENSE+build-type:     Simple extra-source-files:     README.md     ChangeLog.md  source-repository head-    type: git-    location: https://github.com/biocad/cobot-io+  type: git+  location: https://github.com/biocad/cobot-io  library-    exposed-modules:-        Bio.ABI-        Bio.ABI.Clean-        Bio.ABI.Decode-        Bio.FASTA-        Bio.FASTA.Parser-        Bio.FASTA.Type-        Bio.FASTA.Writer-        Bio.GB-        Bio.GB.Parser-        Bio.GB.Type-        Bio.GB.Writer-        Bio.MAE-        Bio.MAE.Parser-        Bio.MAE.Type-        Bio.MMTF-        Bio.MMTF.Decode-        Bio.MMTF.Decode.Codec-        Bio.MMTF.Decode.MessagePack-        Bio.MMTF.MessagePack-        Bio.MMTF.Type-        Bio.PDB-        Bio.PDB.BondRestoring-        Bio.PDB.Functions-        Bio.PDB.Parser-        Bio.PDB.Reader-        Bio.PDB.Type-        Bio.Sequence-        Bio.Sequence.Basecalled-        Bio.Sequence.Basecalled.Type-        Bio.Sequence.Class-        Bio.Sequence.Functions.Marking-        Bio.Sequence.Functions.Sequence-        Bio.Sequence.Functions.Weight-        Bio.Sequence.Utilities-        Bio.Structure-        Bio.Structure.Functions-        Bio.Uniprot-        Bio.Uniprot.Parser-        Bio.Uniprot.Type-    hs-source-dirs: src-    other-modules:-        Paths_cobot_io-    default-language: Haskell2010-    default-extensions: DeriveGeneric DeriveFunctor DeriveFoldable-                        DeriveAnyClass FlexibleInstances InstanceSigs MultiParamTypeClasses-                        RecordWildCards ScopedTypeVariables OverloadedStrings-                        TypeApplications TypeFamilies DataKinds ConstraintKinds-                        TypeOperators TemplateHaskell FlexibleContexts-    build-depends:-        array ==0.5.*,-        attoparsec >=0.10 && <0.14,-        base >=4.7 && <5,-        binary >=0.8.3.0 && <1.0,-        bytestring >=0.10.8.1 && <0.11,-        containers >=0.5.7.1 && <0.7,-        data-msgpack >=0.0.9 && <0.1,-        deepseq ==1.4.*,-        http-conduit ==2.3.*,-        hyraxAbif >=0.2.3.15 && <0.2.4.0,-        lens >=4.16 && <5.0,-        linear ==1.20.*,-        mtl >=2.2.1 && <2.3.0,-        split >=0.2.3.4 && <0.3,-        text >=1.2.2.1 && <1.3,-        vector >=0.12.1.2 && <0.13+  exposed-modules:+      Bio.ABI+      Bio.ABI.Clean+      Bio.ABI.Decode+      Bio.FASTA+      Bio.FASTA.Parser+      Bio.FASTA.Type+      Bio.FASTA.Writer+      Bio.GB+      Bio.GB.Parser+      Bio.GB.Type+      Bio.GB.Writer+      Bio.MAE+      Bio.MAE.Parser+      Bio.MAE.Type+      Bio.MMTF+      Bio.MMTF.Decode+      Bio.MMTF.Decode.Codec+      Bio.MMTF.Decode.MessagePack+      Bio.MMTF.MessagePack+      Bio.MMTF.Type+      Bio.PDB+      Bio.PDB.BondRestoring+      Bio.PDB.Functions+      Bio.PDB.Parser+      Bio.PDB.Reader+      Bio.PDB.Type+      Bio.Sequence+      Bio.Sequence.Basecalled+      Bio.Sequence.Basecalled.Type+      Bio.Sequence.Class+      Bio.Sequence.Functions.Marking+      Bio.Sequence.Functions.Sequence+      Bio.Sequence.Functions.Weight+      Bio.Sequence.Utilities+      Bio.Structure+      Bio.Structure.Functions+      Bio.Uniprot+      Bio.Uniprot.Parser+      Bio.Uniprot.Type+  other-modules:+      Paths_cobot_io+  hs-source-dirs:+      src+  default-extensions: DeriveGeneric DeriveFunctor DeriveFoldable DeriveAnyClass FlexibleInstances InstanceSigs MultiParamTypeClasses RecordWildCards ScopedTypeVariables OverloadedStrings TypeApplications TypeFamilies DataKinds ConstraintKinds TypeOperators TemplateHaskell FlexibleContexts+  build-depends:+      array >=0.5 && <0.6+    , attoparsec >=0.10 && <0.14+    , base >=4.7 && <5+    , binary >=0.8.3.0 && <1.0+    , bytestring >=0.10.8.1 && <0.11+    , containers >=0.5.7.1 && <0.7+    , data-msgpack >=0.0.9 && <0.1+    , deepseq >=1.4 && <1.5+    , http-conduit >=2.3 && <2.4+    , hyraxAbif >=0.2.3.15 && <0.2.4.0+    , lens >=4.16 && <5.0+    , linear >=1.20 && <1.21+    , mtl >=2.2.1 && <2.3.0+    , split+    , text >=1.2.2.1 && <1.3+    , vector+  default-language: Haskell2010  test-suite cobot-io-test-    type: exitcode-stdio-1.0-    main-is: Spec.hs-    hs-source-dirs: test-    other-modules:-        ABISpec-        FastaParserSpec-        FASTASpec-        FastaWriterSpec-        GBParserSpec-        GBWriterSpec-        MAEParserSpec-        MAESpec-        MMTFSpec-        PDBParserSpec-        PDBSpec-        SequenceSpec-        StructureSpec-        UniprotSpec-        Paths_cobot_io-    default-language: Haskell2010-    default-extensions: OverloadedStrings TypeFamilies-    ghc-options: -threaded -rtsopts -with-rtsopts=-N-    build-depends:-        QuickCheck >=2.9.2 && <2.14,-        array ==0.5.*,-        attoparsec >=0.10 && <0.14,-        base >=4.7 && <5,-        binary >=0.8.3.0 && <1.0,-        bytestring >=0.10.8.1 && <0.11,-        cobot-io -any,-        containers >=0.5.7.1 && <0.7,-        data-msgpack >=0.0.9 && <0.1,-        deepseq ==1.4.*,-        directory >=1.3.6.0 && <1.4,-        hspec >=2.4.1 && <2.8,-        http-conduit ==2.3.*,-        hyraxAbif >=0.2.3.15 && <0.2.4.0,-        lens >=4.16 && <5.0,-        linear ==1.20.*,-        mtl >=2.2.1 && <2.3.0,-        neat-interpolation ==0.3.*,-        split >=0.2.3.4 && <0.3,-        text >=1.2.2.1 && <1.3,-        vector >=0.12.1.2 && <0.13+  type: exitcode-stdio-1.0+  main-is: Spec.hs+  other-modules:+      ABISpec+      FastaParserSpec+      FASTASpec+      FastaWriterSpec+      GBParserSpec+      GBWriterSpec+      MAEParserSpec+      MAESpec+      MMTFSpec+      PDBParserSpec+      PDBSpec+      SequenceSpec+      StructureSpec+      UniprotSpec+      Paths_cobot_io+  hs-source-dirs:+      test+  default-extensions: OverloadedStrings TypeFamilies+  ghc-options: -threaded -rtsopts -with-rtsopts=-N+  build-depends:+      QuickCheck >=2.9.2 && <2.14+    , array >=0.5 && <0.6+    , attoparsec >=0.10 && <0.14+    , base >=4.7 && <5+    , binary >=0.8.3.0 && <1.0+    , bytestring >=0.10.8.1 && <0.11+    , cobot-io+    , containers >=0.5.7.1 && <0.7+    , data-msgpack >=0.0.9 && <0.1+    , deepseq >=1.4 && <1.5+    , directory+    , hspec >=2.4.1 && <2.8+    , http-conduit >=2.3 && <2.4+    , hyraxAbif >=0.2.3.15 && <0.2.4.0+    , lens >=4.16 && <5.0+    , linear >=1.20 && <1.21+    , mtl >=2.2.1 && <2.3.0+    , neat-interpolation >=0.3+    , split+    , text >=1.2.2.1 && <1.3+    , vector+  default-language: Haskell2010

@@ -1,6 +1,8 @@ module Bio.ABI.Clean   ( Cleanable (..)+   , Thresholds (..)+  , defaultThresholds   ) where  import           Bio.Sequence            (mean, meanInRange)@@ -18,22 +20,25 @@   cleanWith :: Thresholds -> a -> Maybe a  -- | Thresholds to clean the data.+-- -- ABI file contains sequence with quality. -- By design of sanger sequencing method start and end of the sequence have bad quality. -- Moreover, internal part of the sequence can also has bad quality. -- To clean the data we make 2 steps. -- -- Step 1. Clean edges:+-- --   * take frame with @frameSize@ and go through the sequence; --   * on each step evaluate mean value;---   * if mean value less than @edgeThreshold@, go further;---   * if mean value more than @edgeThreshold@, stop and cut the sequence from END of this frame;+--   * if mean value less than 'edgeThreshold', go further;+--   * if mean value more than 'edgeThreshold' stop and cut the sequence from END of this frame; --   * repeat this algorithm for the right edge. -- -- Step 2. Evaluate quality:+-- --   * for cropped sequence evaluate mean value;---   * if mean value less then @innerThreshold@, sequence is bad;---   * if mean value more then @innerThreshold@, sequence is acceptable.+--   * if mean value less then 'innerThreshold', sequence is bad;+--   * if mean value more then 'innerThreshold', sequence is acceptable. -- -- Logic of this algorithm and 'defaultThresholds' were obtained by taking experiments with read ABI files. --@@ -54,10 +59,10 @@                           then fromBoth                           else Nothing     where-      fromLeft = cutEdge defaultThresholds input+      fromLeft = cutEdge thr input       fromBoth =  fmap S.reverse                .  join-               $  cutEdge defaultThresholds+               $  cutEdge thr                .  S.reverse               <$> fromLeft 

@@ -1,19 +1,55 @@ {-# OPTIONS_GHC -fno-warn-orphans #-} -module Bio.ABI.Decode () where+module Bio.ABI.Decode+  ( decodeRawSequence+  , decodeRawSequence'+  ) where -import           Bio.Sequence               (SequenceDecodable (..),-                                             weightedSequence)-import           Bio.Sequence.Basecalled    (BasecalledSequence)+import           Control.Applicative        (many)+import           Data.Bifunctor             (bimap)+import           Data.Binary.Get            (getInt16be, runGetOrFail) import           Data.ByteString            as BS (ByteString) import           Data.ByteString.Lazy       as BSL (ByteString, fromStrict) import           Data.ByteString.Lazy.Char8 as BSL8 (unpack) import           Data.Char                  (ord)+import           Data.Functor               ((<&>))+import           Data.Int                   (Int16) import           Data.List                  (find)-import           Data.Text                  (Text)+import           Data.Text                  (Text, pack)+import           Data.Vector                (Vector, fromList) import           Hyrax.Abif                 (Abif (..), Directory (..)) import           Hyrax.Abif.Read            (getAbif) +import           Bio.Sequence               (SequenceDecodable (..),+                                             weightedSequence)+import           Bio.Sequence.Basecalled    (BasecalledSequence,+                                             BasecalledSequenceWithRawData (..))++-- | Decode ABIF file with additional raw data that may be required for later processing.+decodeRawSequence :: BSL.ByteString -> Either Text BasecalledSequenceWithRawData+decodeRawSequence bs = do+    abif <- getAbif bs+    sequence' <- extractSequence abif+    quality' <- extractQuality abif+    bsSequence <- weightedSequence sequence' quality'++    bsRawG <- findDataByDirectory "DATA" 9 abif >>= decodeShortArray+    bsRawA <- findDataByDirectory "DATA" 10 abif >>= decodeShortArray+    bsRawT <- findDataByDirectory "DATA" 11 abif >>= decodeShortArray+    bsRawC <- findDataByDirectory "DATA" 12 abif >>= decodeShortArray++    bsPeakLocations <- findDataByDirectory "PLOC" 2 abif >>= decodeShortArray <&> fmap fromIntegral++    return BasecalledSequenceWithRawData{..}++-- | Same as 'decodeRawSequence', for strict @ByteString@.+decodeRawSequence' :: BS.ByteString -> Either Text BasecalledSequenceWithRawData+decodeRawSequence' = decodeRawSequence . BSL.fromStrict++-- | Discards raw data information.+instance SequenceDecodable BasecalledSequenceWithRawData BasecalledSequence where+  sequenceDecode = pure . bsSequence+ -- | Converts 'Data.ByteString.Lazy.ByteString' (that should be content of ABI file) -- into 'BasecalledSequence'. --@@ -39,13 +75,13 @@ -- | Extracts sequence from ABI file. -- extractSequence :: Abif -> Either Text String-extractSequence abif = findDataByDirectory "PBAS" abif >>= checkACGT+extractSequence abif = findDataByDirectory "PBAS" 1 abif <&> BSL8.unpack . dData >>= checkACGT  -- | Extracts quality from ABI file. -- Number are encoded with letters, thus we have function @fromIntegral . ord@. -- extractQuality :: Abif -> Either Text [Double]-extractQuality abif = map (fromIntegral . ord) <$> findDataByDirectory "PCON" abif+extractQuality abif = map (fromIntegral . ord) . BSL8.unpack . dData <$> findDataByDirectory "PCON" 1 abif  -- | Checks that all chars are from alphabet ACGT --@@ -58,13 +94,22 @@  -- | Looks into ABI file and extract data by 'Directory' name. ---findDataByDirectory :: Text -> Abif -> Either Text String-findDataByDirectory dirName abif =-    let directoryM = find (\Directory{..} -> dTagName == dirName) . aDirs $ abif-    in maybe (Left errorMsg) (Right . getData) directoryM+findDataByDirectory+  :: Text -- ^ Directory name+  -> Int  -- ^ 1-based directory index+  -> Abif+  -> Either Text Directory+findDataByDirectory dirName dirIndex abif =+    let directoryM = find (\Directory{..} -> dTagName == dirName && dTagNum == dirIndex) . aDirs $ abif+    in maybe (Left errorMsg) Right directoryM   where     errorMsg :: Text-    errorMsg = "Bio.ABI.Decode: could not find directory " <> dirName+    errorMsg = "Bio.ABI.Decode: could not find directory " <> dirName <> " with index " <> pack (show dirIndex) -    getData :: Directory -> String-    getData = BSL8.unpack . dData+decodeShortArray :: Directory -> Either Text (Vector Int16)+decodeShortArray =+  bimap+    (\(_, _, msg) -> "Data.ABI.Decode: could not decode short array: " <> pack msg)+    (\(_, _, lst) -> fromList lst)+  . runGetOrFail (many getInt16be)+  . dData

@@ -1,13 +1,32 @@ module Bio.Sequence.Basecalled.Type   ( BasecalledSequence (..)+  , BasecalledSequenceWithRawData (..)   ) where  import           Bio.Sequence (IsSequence (..), WeightedSequence) import           Data.Coerce  (coerce)+import           Data.Int     (Int16)+import           Data.Vector  (Vector)  newtype BasecalledSequence = BasecalledSequence (WeightedSequence Double Char)   deriving (Eq, Show) +data BasecalledSequenceWithRawData+  = BasecalledSequenceWithRawData+      { bsSequence      :: BasecalledSequence+        -- A sequence with quality.+      , bsRawG          :: Vector Int16+      , bsRawA          :: Vector Int16+      , bsRawT          :: Vector Int16+      , bsRawC          :: Vector Int16+      , bsPeakLocations :: Vector Int+        -- ^ Same length as 'bsSequence'.+        --+        -- For every base in 'bsSequence' corresponding element in 'bsPeakLocations' is+        -- an index in one of the raw vectors.+      }+  deriving (Eq, Show)+ instance IsSequence BasecalledSequence where   type Element BasecalledSequence = Char   type Marking BasecalledSequence = ()@@ -15,4 +34,3 @@    toSequence = coerce   fromSequence = coerce-

@@ -1,13 +1,18 @@+{-# LANGUAGE RecordWildCards #-}+ module ABISpec where -import           Bio.ABI                 (Cleanable (..))-import           Bio.Sequence            (SequenceDecodable (..))-import qualified Bio.Sequence            as S (getWeights, length, toList)-import           Bio.Sequence.Basecalled (BasecalledSequence) import           Data.ByteString.Lazy    as BSL (readFile) import           Data.Text               (Text) import           Test.Hspec +import           Bio.ABI                 (Cleanable (..))+import           Bio.ABI.Decode          (decodeRawSequence)+import           Bio.Sequence            (SequenceDecodable (..))+import qualified Bio.Sequence            as S (getWeights, length, toList)+import           Bio.Sequence.Basecalled (BasecalledSequence,+                                          BasecalledSequenceWithRawData (..))+ abiExtractSpec :: Spec abiExtractSpec =   describe "ABI decode" $ do@@ -21,6 +26,13 @@       datM <- readData "test/ABI/not_ab1.txt"       datM `shouldBe` Left "Error reading root: not enough bytes" +    it "decode with raw data" $ do+      Right BasecalledSequenceWithRawData{..} <- decodeRawSequence <$> BSL.readFile "test/ABI/test.ab1"+      S.length bsSequence `shouldBe` length bsPeakLocations+      length bsRawG `shouldSatisfy` (>0)+      length bsRawA `shouldSatisfy` (>0)+      length bsRawT `shouldSatisfy` (>0)+      length bsRawC `shouldSatisfy` (>0)  abiCleanSpec :: Spec abiCleanSpec =

@@ -27,33 +27,34 @@  rawPDBToModelConversionSingleChainSpec :: SpecWith () rawPDBToModelConversionSingleChainSpec = describe "Cobot Model from raw single chain PDB" $ do-  modelFromPDB <- runIO $ firstPDBModel "test/PDB/1PPE_I.pdb"-  modelFromMae <- runIO $ firstMaeModel "test/PDB/1PPE_I.mae"-  let (pdbBondCount, _, pdbChainCount, pdbAtomCount) = getStats modelFromPDB-  let (maeBondCount, _, _, _) = getStats modelFromMae--  it "Should have correct number of atoms" $ pdbAtomCount `shouldBe` 436-  it "Should have correct number of chains" $ pdbChainCount `shouldBe` 1-  it "Should restore bonds correctly" $ pdbBondCount `shouldBe` maeBondCount+    let pdbStats = getStats <$> firstPDBModel "test/PDB/1PPE_I.pdb"+    let maeStats = getStats <$> firstMaeModel "test/PDB/1PPE_I.mae"+    beforeAll ((,) <$> pdbStats <*> maeStats) $ do+        it "Should have correct number of atoms" $ \((_, _, _, pdbAtomCount), _) ->+          pdbAtomCount `shouldBe` 436+        it "Should have correct number of chains" $ \((_, _, pdbChainCount, _), _) ->+          pdbChainCount `shouldBe` 1+        it "Should restore bonds correctly" $ \((pdbBondCount, _, _, _), (maeBondCount, _, _, _)) ->+          pdbBondCount `shouldBe` maeBondCount --- tripeptides are not checked as in checkBiggerMolecule +-- tripeptides are not checked as in checkBiggerMolecule -- because there are inconsistencies in atom numbers between pdb and mae bondsRestoringTripeptideSpec :: SpecWith () bondsRestoringTripeptideSpec = describe "Bonds should be restored correctly in tripeptides" $   sequence_ $ checkTripeptide <$> tripeptides     where       tripeptides :: [String]-      tripeptides = ["ALA_3", "ARG_3", "ASN_3", "ASP_3", "CYS_3", "GLN_3", "GLU_3", "GLY_3", "HID_3", "HIE_3", "HIP_3", +      tripeptides = ["ALA_3", "ARG_3", "ASN_3", "ASP_3", "CYS_3", "GLN_3", "GLU_3", "GLY_3", "HID_3", "HIE_3", "HIP_3",                      "ILE_3", "LEU_3", "LYS_3", "MET_3", "PHE_3", "PRO_3", "SER_3", "THR_3", "TRP_3", "TYR_3", "VAL_3"]       checkTripeptide :: String -> SpecWith (Arg Expectation)-      checkTripeptide tripeptideName = do-        modelFromMae <- runIO . firstMaeModel $ "test/PDB/BondsRestoring/" ++ tripeptideName ++ ".mae"-        modelFromPDB <- runIO . firstPDBModel $ "test/PDB/BondsRestoring/" ++ tripeptideName ++ ".pdb"+      checkTripeptide tripeptideName = it (tripeptideName ++ " equal bond count in Mae and PDB") $ do+          modelFromMae <- firstMaeModel $ "test/PDB/BondsRestoring/" ++ tripeptideName ++ ".mae"+          modelFromPDB <- firstPDBModel $ "test/PDB/BondsRestoring/" ++ tripeptideName ++ ".pdb" -        let pdbBondCount = V.length $ modelBonds modelFromPDB-        let maeBondCount = V.length $ modelBonds modelFromMae+          let pdbBondCount = V.length $ modelBonds modelFromPDB+          let maeBondCount = V.length $ modelBonds modelFromMae -        it (tripeptideName ++ " equal bond count in Mae and PDB") $ pdbBondCount `shouldBe` maeBondCount+          pdbBondCount `shouldBe` maeBondCount  bondsRestoringBiggerMoleculesSpec :: SpecWith () bondsRestoringBiggerMoleculesSpec = describe "Bonds should be restored correctly in bigger molecules" $ do@@ -62,31 +63,36 @@   checkBiggerMolecule "4dn4_ag_b"   where     checkBiggerMolecule moleculeName = do-      modelFromPDB <- runIO . firstPDBModel $ "test/PDB/BondsRestoring/" ++ moleculeName ++ ".pdb"-      modelFromMae <- runIO . firstMaeModel $ "test/PDB/BondsRestoring/" ++ moleculeName ++ ".mae"-      let (pdbGlobalBondCount, pdbLocalBondCount, _, _) = getStats modelFromPDB-      let (maeGlobalBondCount, maeLocalBondCount, _, _) = getStats modelFromMae--      it (moleculeName ++ " equal global bond count in Mae and PDB") $ pdbGlobalBondCount `shouldBe` maeGlobalBondCount-      it (moleculeName ++ " equal local bond count in Mae and PDB") $ pdbLocalBondCount `shouldBe` maeLocalBondCount--      it (moleculeName ++ " no dublicate bonds") $ length (doubleBonds modelFromPDB) `shouldBe` 0--      let _globalBondSetPDB = globalBondSet modelFromPDB-      let _globalBondSetMae = globalBondSet modelFromMae-      let diffMaePDBGlobal = S.difference _globalBondSetMae _globalBondSetPDB-      let diffPDBMaeGlobal = S.difference _globalBondSetPDB _globalBondSetMae+      let pdbIO = firstPDBModel $ "test/PDB/BondsRestoring/" ++ moleculeName ++ ".pdb"+      let maeIO = firstMaeModel $ "test/PDB/BondsRestoring/" ++ moleculeName ++ ".mae"+      beforeAll ((,) <$> pdbIO <*> maeIO) $ do+        it (moleculeName ++ " equal global bond count in Mae and PDB") $ \(modelFromPDB, modelFromMae) -> do+            let (pdbGlobalBondCount, _, _, _) = getStats modelFromPDB+            let (maeGlobalBondCount, _, _, _) = getStats modelFromMae+            pdbGlobalBondCount `shouldBe` maeGlobalBondCount+        it (moleculeName ++ " equal local bond count in Mae and PDB") $ \(modelFromPDB, modelFromMae) -> do+            let (_, pdbLocalBondCount, _, _) = getStats modelFromPDB+            let (_, maeLocalBondCount, _, _) = getStats modelFromMae+            pdbLocalBondCount `shouldBe` maeLocalBondCount -      it (moleculeName ++ " difference in Mae and PDB global bond sets") $ S.size diffMaePDBGlobal `shouldBe` 0-      it (moleculeName ++ " difference in PDB and Mae global bond sets") $ S.size diffPDBMaeGlobal `shouldBe` 0+        it (moleculeName ++ " no dublicate bonds") $ \(modelFromPDB, _) ->+          length (doubleBonds modelFromPDB) `shouldBe` 0 -      let _localBondSetPDB = localBondSet modelFromPDB-      let _localBondSetMae = localBondSet modelFromMae-      let diffMaePDBLocal = S.difference _localBondSetMae _localBondSetPDB-      let diffPDBMaeLocal = S.difference _localBondSetPDB _localBondSetMae+        it (moleculeName ++ " difference in Mae and PDB global bond sets") $ \(modelFromPDB, modelFromMae) -> do+            let _globalBondSetPDB = globalBondSet modelFromPDB+            let _globalBondSetMae = globalBondSet modelFromMae+            let diffPDBMaeGlobal = S.difference _globalBondSetPDB _globalBondSetMae+            let diffMaePDBGlobal = S.difference _globalBondSetMae _globalBondSetPDB+            S.size diffPDBMaeGlobal `shouldBe` 0+            S.size diffMaePDBGlobal `shouldBe` 0 -      it (moleculeName ++ " difference in Mae and PDB local bond sets") $ S.size diffMaePDBLocal `shouldBe` 0-      it (moleculeName ++ " difference in PDB and Mae local bond sets") $ S.size diffPDBMaeLocal `shouldBe` 0+        it (moleculeName ++ " difference in Mae and PDB local bond sets") $ \(modelFromPDB, modelFromMae) -> do+            let _localBondSetPDB = localBondSet modelFromPDB+            let _localBondSetMae = localBondSet modelFromMae+            let diffMaePDBLocal = S.difference _localBondSetMae _localBondSetPDB+            let diffPDBMaeLocal = S.difference _localBondSetPDB _localBondSetMae+            S.size diffMaePDBLocal `shouldBe` 0+            S.size diffPDBMaeLocal `shouldBe` 0      localBondSet :: Model -> Set (Text, Int, Int, Int) -- (ChainID, ResidueNumber, LocalFrom, LocalTo)     localBondSet Model{..} = S.fromList $ do@@ -127,7 +133,7 @@  getStats :: Model -> (Int, Int, Int, Int) getStats model = (globalBondCount, localBondCount, chainCount, atomCount)-  where +  where     globalBondCount = V.length $ modelBonds model     localBondCount = V.length . V.concatMap (V.concatMap resBonds . chainResidues) $ modelChains model     chainCount = V.length $ modelChains model

@@ -23,12 +23,12 @@  structureSpec :: Spec structureSpec = describe "Structure spec." $ do-    m <- runIO $ V.toList . modelsOf <$> fromFile "test/MAE/Capri.mae" >>= \[x] -> pure x--    it "atoms filtering works correctly. only N, CA, C"  $ checkFiltering m $ (`elem` ["N", "CA", "C"]) . atomName-    it "atoms filtering works correctly. only CA"        $ checkFiltering m $ (== "CA") . atomName-    it "atoms filtering works correctly. no atoms"       $ checkFiltering m $ const False-    it "atoms filtering works correctly. all atoms"      $ checkFiltering m $ const True+    let mIO = V.toList . modelsOf <$> fromFile "test/MAE/Capri.mae" >>= \[x] -> pure x+    beforeAll mIO $ do+        it "atoms filtering works correctly. only N, CA, C"  $ \m -> checkFiltering m $ (`elem` ["N", "CA", "C"]) . atomName+        it "atoms filtering works correctly. only CA"        $ \m -> checkFiltering m $ (== "CA") . atomName+        it "atoms filtering works correctly. no atoms"       $ \m -> checkFiltering m $ const False+        it "atoms filtering works correctly. all atoms"      $ \m -> checkFiltering m $ const True   where     checkFiltering :: Model -> (Atom -> Bool) -> Expectation     checkFiltering m p = do