Skip to content

Commit

Permalink
sandwich.hmm+ works
Browse files Browse the repository at this point in the history
  • Loading branch information
Andrew Gallant (Ocelot) authored and Andrew Gallant (Ocelot) committed Jul 10, 2012
1 parent c0c0a81 commit 2e9b5cc
Show file tree
Hide file tree
Showing 3 changed files with 34 additions and 30 deletions.
40 changes: 22 additions & 18 deletions Beta.hs
Expand Up @@ -36,29 +36,33 @@ import MRFTypes


getBetaStrands :: [StrandPair] -> [BetaStrand]
getBetaStrands ss = addIndexInfo
$ mkBetaStrands
$ mergeStrands
$ addPairings
$ mkBetaResidues ss
getBetaStrands = addIndexInfo
. mkBetaStrands
. mergeStrands
. addPairings
. mkBetaResidues

-- Decorates the *pairs* with index information so that
-- Beta scoring is efficient later on.
addIndexInfo :: [BetaStrand] -> [BetaStrand]
addIndexInfo betas = addIndexInfo' betas
where addIndexInfo' [] = []
addIndexInfo' (b:bs) = b' : addIndexInfo' bs
where b' = b { residues = map addToResidue $ residues b }
addIndexInfo betas = map addIndexInfo' betas
where addIndexInfo' b = b { residues = map addToResidue $ residues b }

addToResidue :: BetaResidue -> BetaResidue
addToResidue r = r { pairs = map (addToPair betas) $ pairs r }
where addToPair [] _ = error "could not find beta position in list of beta strands"
addToPair (b:bs) p =
case find (\r' -> pairPosition p == resPosition r') $ residues b of
Just r' -> p { pairStrandSerial = resStrandSerial r'
, residueInd = maybe (error "whoa there") id $ elemIndex r' $ residues b
}
Nothing -> addToPair bs p

addToPair :: [BetaStrand] -> BetaPair -> BetaPair
addToPair [] _ =
error "could not find beta position in list of beta strands"
addToPair (b:bs) p =
case find (\r' -> pairPosition p == resPosition r') $ residues b of
Just r' ->
p { pairStrandSerial = resStrandSerial r'
, residueInd = maybe (error "unreachable") id
$ elemIndex r'
$ residues b
}
Nothing -> addToPair bs p

-- Creates the BetaStrand structures and adds serial
-- information to each residue.
Expand Down Expand Up @@ -102,7 +106,7 @@ mergeStrands (s:ss) = s' : (mergeStrands $ filter (emptyIntersection s') ss)
-- haven't been made yet. Effectively, a strand is a list
-- of beta nodes.
addPairings :: [[BetaResidue]] -> [[BetaResidue]]
addPairings strands = map (map (\n -> foldl annotate n strands)) strands
addPairings strands = sort $ map (map (\n -> foldl annotate n strands)) strands
where -- Find one or zero equivalent nodes in a beta strand
-- and add any additional pairs found to 'r'
annotate :: BetaResidue -> [BetaResidue] -> BetaResidue
Expand All @@ -118,7 +122,7 @@ mkBetaResidues :: [StrandPair] -> [[BetaResidue]]
mkBetaResidues [] = []
mkBetaResidues (sp:sps) = residues1 : residues2 : mkBetaResidues sps
where (residues1, residues2) =
mkResidues (zip3 (unElist $ exposure sp) [s1..] secondIndexing) [] []
mkResidues (zip3 (exposure sp) [s1..] secondIndexing) [] []

secondIndexing = case parallel sp of Parallel -> [s2..]
Antiparallel -> [s2, s2-1..]
Expand Down
4 changes: 2 additions & 2 deletions HMMPlus.hs
Expand Up @@ -138,7 +138,7 @@ p_int = decimal
p_direction = liftA mkDirection (optSpaces *> dir <* optSpaces <?> "direction")
where dir = oneStringOf [ "1", "-1" ]

p_exposure = EList <$> liftA (map mkExposure)
p_exposure = liftA (map mkExposure)
(optSpaces *> many1 (oneOf "io") <* optSpaces <?> "exposure")


Expand Down Expand Up @@ -216,4 +216,4 @@ logProb = (Score <$> fractional) <|> (negLogZero <$ char '*') <?> "logProb"



eoH = string "HMM" <* manyTill anyChar eol <?> "hmmSeparator"
eoH = string "HMM" <* manyTill anyChar eol <?> "hmmSeparator"
20 changes: 10 additions & 10 deletions MRFTypes.hs
Expand Up @@ -6,8 +6,7 @@ module MRFTypes
, matchEmissions, insertionEmissions
, StrandPair(..)
, Helix(..)
, Exposure(..), mkExposure, unElist
, ExposureList(..)
, Exposure(..), mkExposure
, Direction(..), mkDirection
, BetaStrand(..), BetaPosition, BetaResidue(..), BetaPair(..)
, TProb(..), TProbs, m_m, m_i, m_d, i_m, i_i, d_m, d_d, b_m, m_e
Expand Down Expand Up @@ -124,7 +123,7 @@ instance Arbitrary StrandPair where
where (d, e) = (arbitrary, arbitrary)
(first, second) = (choose (1, max), choose (1, max))
len = choose (1, max)
max = 1000
max = 100

data Helix = Helix { startRes :: Int
, helixLength :: Int
Expand All @@ -144,13 +143,14 @@ data HMMHeader = HMMHeader { betas :: [BetaStrand] -- TODO add alphas (sst
data MRF = MRF { hmmHeader :: HMMHeader, hmm :: HMM }
deriving (Show)

newtype ExposureList = EList [ Exposure ] deriving Show

unElist :: ExposureList -> [Exposure]
unElist (EList exposures) = exposures

instance Arbitrary ExposureList where
arbitrary = fmap EList $ listOf (arbitrary :: Gen Exposure)
type ExposureList = [Exposure]
-- newtype ExposureList = EList [ Exposure ] deriving Show
--
-- unElist :: ExposureList -> [Exposure]
-- unElist (EList exposures) = exposures
--
-- instance Arbitrary ExposureList where
-- arbitrary = fmap EList $ listOf (arbitrary :: Gen Exposure)


data Exposure = Buried -- 'i'
Expand Down

0 comments on commit 2e9b5cc

Please sign in to comment.