-- | Encoding of Watson-Crick and Wobble Pairs in the Vienna RNA package style. module Biobase.Secondary.Vienna where import Data.Aeson import Data.Binary import Data.Ix import Data.Primitive.Types import Data.Serialize (Serialize(..)) import Data.Tuple (swap) import Data.Vector.Fusion.Stream.Monadic (map,Step(..)) import Data.Vector.Unboxed.Deriving import GHC.Base (remInt,quotInt) import GHC.Generics (Generic) import Prelude hiding (map) import qualified Data.Vector.Generic as VG import qualified Data.Vector.Generic.Mutable as VGM import qualified Data.Vector.Unboxed as VU import qualified Prelude as P import Data.PrimitiveArray hiding (Complement(..),map) import Biobase.Primary.Letter import Biobase.Primary.Nuc import Biobase.Primary.Nuc.RNA -- | Use machine Ints internally newtype ViennaPair = ViennaPair { unViennaPair :: Int } deriving (Eq,Ord,Generic,Ix) instance Binary (ViennaPair) instance Serialize (ViennaPair) instance FromJSON (ViennaPair) instance ToJSON (ViennaPair) instance Index ViennaPair where linearIndex _ _ (ViennaPair p) = p {-# Inline linearIndex #-} smallestLinearIndex _ = error "still needed?" {-# Inline smallestLinearIndex #-} largestLinearIndex (ViennaPair p) = p {-# Inline largestLinearIndex #-} size _ (ViennaPair h) = h+1 {-# Inline size #-} inBounds (ViennaPair l) (ViennaPair h) (ViennaPair p) = l <= p && p <= h {-# Inline inBounds #-} instance IndexStream z => IndexStream (z:.ViennaPair) where streamUp (ls:.ViennaPair l) (hs:.ViennaPair h) = flatten mk step $ streamUp ls hs where mk z = return (z,l) step (z,k) | k > h = return $ Done | otherwise = return $ Yield (z:.ViennaPair k) (z,k+1) {-# Inline [0] mk #-} {-# Inline [0] step #-} {-# Inline streamUp #-} streamDown (ls:.ViennaPair l) (hs:.ViennaPair h) = flatten mk step $ streamDown ls hs where mk z = return (z,h) step (z,k) | k < l = return $ Done | otherwise = return $ Yield (z:.ViennaPair k) (z,k-1) {-# Inline [0] mk #-} {-# Inline [0] step #-} {-# Inline streamDown #-} instance IndexStream ViennaPair where streamUp l h = map (\(Z:.k) -> k) $ streamUp (Z:.l) (Z:.h) {-# Inline streamUp #-} streamDown l h = map (\(Z:.k) -> k) $ streamDown (Z:.l) (Z:.h) {-# Inline streamDown #-} {- instance (Shape sh,Show sh) => Shape (sh :. ViennaPair) where rank (sh:._) = rank sh + 1 zeroDim = zeroDim:.ViennaPair 0 unitDim = unitDim:.ViennaPair 1 -- TODO does this one make sense? intersectDim (sh1:.n1) (sh2:.n2) = intersectDim sh1 sh2 :. min n1 n2 addDim (sh1:.ViennaPair n1) (sh2:.ViennaPair n2) = addDim sh1 sh2 :. ViennaPair (n1+n2) -- TODO will not necessarily yield a valid ViennaPair size (sh1:.ViennaPair n) = size sh1 * n sizeIsValid (sh1:.ViennaPair n) = sizeIsValid (sh1:.n) toIndex (sh1:.ViennaPair sh2) (sh1':.ViennaPair sh2') = toIndex (sh1:.sh2) (sh1':.sh2') fromIndex (ds:.ViennaPair d) n = fromIndex ds (n `quotInt` d) :. ViennaPair r where r | rank ds == 0 = n | otherwise = n `remInt` d inShapeRange (sh1:.n1) (sh2:.n2) (idx:.i) = i>=n1 && i error "empty list in shapeOfList/Primary" x:xs -> shapeOfList xs :. ViennaPair x deepSeq (sh:.n) x = deepSeq sh (n `seq` x) {-# INLINE rank #-} {-# INLINE zeroDim #-} {-# INLINE unitDim #-} {-# INLINE intersectDim #-} {-# INLINE addDim #-} {-# INLINE size #-} {-# INLINE sizeIsValid #-} {-# INLINE toIndex #-} {-# INLINE fromIndex #-} {-# INLINE inShapeRange #-} {-# INLINE listOfShape #-} {-# INLINE shapeOfList #-} {-# INLINE deepSeq #-} instance (Eq sh, Shape sh, Show sh, ExtShape sh) => ExtShape (sh :. ViennaPair) where subDim (sh1:.ViennaPair n1) (sh2:.ViennaPair n2) = subDim sh1 sh2 :. (ViennaPair $ n1-n2) rangeList (sh1:.ViennaPair n1) (sh2:.ViennaPair n2) = [sh:.ViennaPair n | sh <- rangeList sh1 sh2, n <- [n1 .. (n1+n2)]] rangeStream (fs:.ViennaPair f) (ts:.ViennaPair t) = VM.flatten mk step Unknown $ rangeStream fs ts where mk sh = return (sh :. f) step (sh :. k) | k>t = return $ VM.Done | otherwise = return $ VM.Yield (sh :. ViennaPair k) (sh :. k +1) {-# INLINE [1] mk #-} {-# INLINE [1] step #-} {-# INLINE rangeStream #-} -} pattern NP = ViennaPair 0 :: ViennaPair pattern CG = ViennaPair 1 :: ViennaPair pattern GC = ViennaPair 2 :: ViennaPair pattern GU = ViennaPair 3 :: ViennaPair pattern UG = ViennaPair 4 :: ViennaPair pattern AU = ViennaPair 5 :: ViennaPair pattern UA = ViennaPair 6 :: ViennaPair pattern NS = ViennaPair 7 :: ViennaPair pattern Undef = ViennaPair 8 :: ViennaPair class MkViennaPair a where mkViennaPair :: a -> ViennaPair fromViennaPair :: ViennaPair -> a instance MkViennaPair (Letter RNA, Letter RNA) where mkViennaPair = \case (C,G) -> CG (G,C) -> GC (G,U) -> GU (U,G) -> UG (A,U) -> AU (U,A) -> UA _ -> NS {-# INLINE mkViennaPair #-} fromViennaPair = \case CG -> (C,G) GC -> (G,C) GU -> (G,U) UG -> (U,G) AU -> (A,U) UA -> (U,A) _ -> error "non-standard pairs can't be backcasted" {-# INLINE fromViennaPair #-} isViennaPair :: Letter RNA -> Letter RNA -> Bool isViennaPair l r = l==C && r==G || l==G && r==C || l==A && r==U || l==U && r==A || l==G && r==U || l==U && r==G {-# INLINE isViennaPair #-} viennaPairTable :: Unboxed (Z:.Letter RNA:.Letter RNA) ViennaPair viennaPairTable = fromAssocs (Z:.N:.N) (Z:.U:.U) NS [ (Z:.C:.G , CG) , (Z:.G:.C , GC) , (Z:.G:.U , GU) , (Z:.U:.G , UG) , (Z:.A:.U , AU) , (Z:.U:.A , UA) ] {-# NOINLINE viennaPairTable #-} instance Enum ViennaPair where toEnum x | x>=0 && x<=7 = ViennaPair x | otherwise = error $ "can't make to enum" ++ show x fromEnum (ViennaPair x) = x {-# INLINE toEnum #-} {-# INLINE fromEnum #-} instance Bounded ViennaPair where minBound = NP maxBound = NS instance Show ViennaPair where show x | Just s <- x `lookup` pairToString = s | otherwise = "??" instance Read ViennaPair where readsPrec p [] = [] readsPrec p [x] = [] readsPrec p (x:y:xs) | x ==' ' = readsPrec p (y:xs) | Just n <- (x:y:[]) `lookup` s2p = [(n,xs)] | otherwise = [] where s2p = (P.map swap pairToString) -- | reverse a vienna pair revPair :: ViennaPair -> ViennaPair revPair = \case CG -> GC GC -> CG GU -> UG UG -> GU AU -> UA UA -> AU NP -> NP NS -> NS -- * Convenience structures cguaP = [CG .. UA] cgnsP = [CG .. NS] pairToString = [(CG,"CG"),(GC,"GC"),(UA,"UA"),(AU,"AU"),(GU,"GU"),(UG,"UG"),(NS,"NS"),(NP,"NP")] derivingUnbox "ViennaPair" [t| ViennaPair -> Int |] [| unViennaPair |] [| ViennaPair |]