estreps 0.1 → 0.3
raw patch · 6 files changed
+41/−26 lines, 6 filesnew-uploader
Files
- README +3/−2
- estreps.cabal +6/−5
- src/RSelect.lhs +21/−11
- src/Repeats.lhs +2/−1
- src/Unigene.hs +5/−4
- src/Util.hs +4/−3
README view
@@ -15,11 +15,12 @@ USAGE - rselect n [m] input.seq+ rselect [-r] n [m] input.seq Selects n sequences from the file input.seq. If the optional m is given, this limits the selection to happen only from the first- m sequences in the file, which may be more efficient.+ m sequences in the file, which may be more efficient. If -r is given,+ the sequences will be reoriented randomly. The selected sequences are written to standard output, so you probably want to redirect them to a file.
estreps.cabal view
@@ -1,5 +1,5 @@ Name: estreps-Version: 0.1+Version: 0.3 License: GPL License-file: LICENSE Author: Ketil Malde@@ -7,12 +7,13 @@ Stability: Beta Category: Bioinformatics Synopsis: Repeats from ESTs-Description: rselect - select a random set of sequences from a FASTA file.- reps - extract exact k-word repeats based that occur in- sequences grouped in different clusters.+Description: * rselect - select a random set of sequences from a FASTA file, optinally+ with random orientation (forward/reverse complement).+ * reps - extract exact k-word repeats based that occur in+ sequences grouped in different clusters. . The Darcs repository is at: <http://malde.org/~ketil/biohaskell/estreps>.-Homepage: http://malde.org/~ketil/+Homepage: http://blog.malde.org/ Tested-With: GHC==6.8.2 Build-Type: Simple
src/RSelect.lhs view
@@ -8,26 +8,28 @@ module Main where -import Data.Set+import Data.Set hiding (partition) import System.Environment (getArgs) import System.Random import System.IO+import Data.List (isPrefixOf,partition) -import Bio.Sequence.Fasta+import Bio.Sequence.Fasta (countSeqs)+import Bio.Sequence import Unigene main :: IO () main = do- as <- getArgs- case as of [m,f] -> select (read m) Nothing f- [m,n,f] -> select (read m) (Just (read n)) f- _ -> error "Usage: rselect n [m] input.seq > output.seq"+ (opts,as) <- partition (isPrefixOf "-") `fmap` getArgs+ case as of [m,f] -> select opts (read m) Nothing f+ [m,n,f] -> select opts (read m) (Just (read n)) f+ _ -> error "Usage: rselect [-r] n [m] input.seq > output.seq" -- | Select m out of maybe n sequences from f -- Not specifying a maximum increases memory consumption-select :: Integer -> Maybe Integer -> FilePath -> IO ()-select m mn f = do+select :: [String] -> Integer -> Maybe Integer -> FilePath -> IO ()+select os m mn f = do seed <- newStdGen hPutStrLn stderr ("Random seed: "++show seed) ss <- ugRead f@@ -35,14 +37,21 @@ n' <- countSeqs f return (fromIntegral n') Just n' -> return n'- let ds = if m>n then error "rselect: m is larger than n, that won't work."++ rds <- randoms `fmap` newStdGen :: IO [Bool]+ let randomize_dirs = if "-r" `elem` os then zipWith maybe_reverse rds else id++ ds = if m>n then error "rselect: m is larger than n, that won't work." else if m==n then take (fromIntegral n) $ repeat 1 else if m<=n `div` 2 then diffs $ pick m $ randomRs (1,n) seed else diffs $ invert n 1 $ pick (n-m) $ randomRs (1,n) seed- out = Prelude.filter (not.Prelude.null) $ choose ds ss+ out = Prelude.filter (not.Prelude.null) $ Prelude.map randomize_dirs $ choose ds ss n `seq` hPutStrLn stderr ("N = "++show n++", M = "++show m)- mapM_ (\o -> putStrLn "# cluster delimiter" >> hWriteFasta stdout o) out+ mapM_ (\o -> putStrLn "# cluster delimiter" >> hWriteFasta stdout o) out +maybe_reverse :: Bool -> Sequence -> Sequence+maybe_reverse b = if b then revcompl else id+ -- | choose elements according to a list of deltas choose :: Show a => [Integer] -> [[a]] -> [[a]] choose = choose' []@@ -62,6 +71,7 @@ diffs' (a:b:cs) = if b<=a then error ("diffs needs a monotonically increasing sequence"++show (take 10 (a:b:cs))) else (b-a) : diffs' (b:cs) diffs' [_] = []+ diffs' [] = [] -- | convert a list to the list of missing values invert :: Integer -> Integer -> [Integer] -> [Integer]
src/Repeats.lhs view
@@ -33,10 +33,11 @@ import Data.Maybe import System.Environment (getArgs) -import Util+-- import Util -- type Word = Int +main :: IO () main = do args <- getArgs when (length args < 2 || length args > 3)
src/Unigene.hs view
@@ -16,10 +16,11 @@ ugRead f = do s <- B.readFile f let ls = dropWhile B.null $ B.lines $ s- when (B.head (head ls) /= '#') (hPutStrLn stderr ("Warning: '"++f++"' does not look like Unigene format"))+ when (B.head (head ls) /= '#') (hPutStrLn stderr ("'"++f++"' does not look like Unigene format")) return (clusters $ ls) clusters :: [B.ByteString] -> [[Sequence]]-clusters (l:ls) = let ls' = if B.head l == '#' then (l:ls) else (B.pack "# dummy UG hdr" :l:ls)- ss = groupBy (const (('#' /=) . B.head)) $ filter (not . B.null) ls'- in map (mkSeqs . tail) ss+clusters (l:ls) = if B.head l == '#' + then map (mkSeqs . tail) $ groupBy (const (('#' /=) . B.head)) $ filter (not . B.null) (l:ls)+ else [mkSeqs (l:ls)]+clusters [] = []
src/Util.hs view
@@ -3,10 +3,11 @@ import System.IO import System.IO.Unsafe +countIO :: Int -> String -> [a] -> IO [a] countIO step msg xs = sequence $ map unsafeInterleaveIO $ ct 0 xs- where ct s xs = let (a,b) = splitAt (step-1) xs- next = s+step- output c = hPutStr stderr ((take (length (show c)+2+length msg) $ repeat ' ')++'\r':msg++show c) >> hFlush stderr+ where ct s xs' = let (a,b) = splitAt (step-1) xs'+ next = s+step+ output c = hPutStr stderr ((take (length (show c)+2+length msg) $ repeat ' ')++'\r':msg++show c) >> hFlush stderr in if null b then map return a ++ [] else map return a ++ [output next >> return (head b)] ++ ct next (tail b)