diff --git a/RNAlien.cabal b/RNAlien.cabal --- a/RNAlien.cabal +++ b/RNAlien.cabal @@ -2,7 +2,7 @@ -- see http://haskell.org/cabal/users-guide/ name: RNAlien -version: 1.1.3 +version: 1.2.4 synopsis: Unsupervized construction of RNA family models description: RNAlien is a tool for automatic construction of RNAfamily models from a single sequence. . @@ -53,8 +53,8 @@ source-repository this type: git - location: https://github.com/eggzilla/RNAlien/tree/1.1.3 - tag: 1.1.3 + location: https://github.com/eggzilla/RNAlien/tree/1.2.4 + tag: 1.2.4 executable RNAlien Hs-Source-Dirs: ./src/Bio/ @@ -74,10 +74,22 @@ ghc-options: -Wall -O2 build-depends: base >=4.5 && <5, cmdargs, either-unwrap, RNAlien +executable cmsearchToBed + Hs-Source-Dirs: ./src/Bio/ + main-is: cmsearchToBED.hs + ghc-options: -Wall -O2 + build-depends: base >=4.5 && <5, cmdargs, either-unwrap, RNAlien, bytestring, text + +executable RNAcentralHTTPRequest + Hs-Source-Dirs: ./src/Bio/ + main-is: RNAcentralHTTPRequest.hs + ghc-options: -Wall -O2 + build-depends: base >=4.5 && <5, cmdargs, either-unwrap, RNAlien + Library Hs-Source-Dirs: ./src/ ghc-options: -Wall -O2 -fno-warn-unused-do-bind - build-depends: base >=4.5 && <5, cmdargs, ViennaRNAParser>=1.2.9, process, directory, blastxml>=0.3.2, biofasta, parsec, random, BlastHTTP, biocore, bytestring, Taxonomy >= 1.0.1, either-unwrap, containers, ClustalParser>=1.1.0, EntrezHTTP>=1.0.1, vector, edit-distance, cassava, matrix, hierarchical-clustering, filepath, HTTP, http-conduit, hxt, network, aeson, text, transformers >= 0.4 && <0.5, pureMD5 + build-depends: base >=4.5 && <5, cmdargs, ViennaRNAParser>=1.3.1, process, directory, blastxml>=0.3.2, biofasta, parsec, random, BlastHTTP, biocore, bytestring, Taxonomy >= 1.0.2, either-unwrap, containers, ClustalParser>=1.1.0, EntrezHTTP>=1.0.2, vector, edit-distance, cassava, matrix, hierarchical-clustering, filepath, HTTP, http-conduit, hxt, network, aeson, text, transformers, pureMD5 Exposed-Modules: Bio.RNAlienData Bio.RNAlienLibrary Bio.RNAcentralHTTP diff --git a/src/Bio/RNAcentralHTTP.hs b/src/Bio/RNAcentralHTTP.hs --- a/src/Bio/RNAcentralHTTP.hs +++ b/src/Bio/RNAcentralHTTP.hs @@ -6,6 +6,7 @@ -- module Bio.RNAcentralHTTP (rnaCentralHTTP, buildSequenceViaMD5Query, + buildStringViaMD5Query, getRNACentralEntries, showRNAcentralAlienEvaluation, RNAcentralEntryResponse, @@ -29,9 +30,9 @@ -- | Data structure for RNAcentral entry response data RNAcentralEntryResponse = RNAcentralEntryResponse { - _count :: Int, - _next :: Maybe Text, - _previous :: Maybe Text, + count :: Int, + next :: Maybe Text, + previous :: Maybe Text, results :: [RNAcentralEntry] } deriving (Show, Eq, Generic) @@ -44,13 +45,13 @@ data RNAcentralEntry = RNAcentralEntry { - _url :: Text, + url :: Text, rnacentral_id :: Text, md5 :: Text, - _sequence :: Text, + sequence :: Text, length :: Int, - _xrefs :: Text, - _publications :: Text + xrefs :: Text, + publications :: Text } deriving (Show, Eq, Generic) @@ -65,13 +66,18 @@ startSession query' = do requestXml <- withSocketsDo $ sendQuery query' + --putStr (L8.unpack requestXml) let eitherErrorResponse = eitherDecode requestXml :: Either String RNAcentralEntryResponse return eitherErrorResponse -- | Send query and return response XML sendQuery :: String -> IO L8.ByteString -sendQuery query' = simpleHttp ("http://rnacentral.org/api/v1/rna/" ++ query') - +sendQuery query' = do + let address = "http://rnacentral.org/api/v1/rna/" + let request = address ++ query' + --putStrLn request + simpleHttp request + -- | Function for querying the RNAcentral REST interface. rnaCentralHTTP :: String -> IO (Either String RNAcentralEntryResponse) rnaCentralHTTP query' = do @@ -88,19 +94,29 @@ responses <- mapM delayedRNACentralHTTP queries return responses +--Build a query from a input sequence buildSequenceViaMD5Query :: Sequence -> String buildSequenceViaMD5Query s = qString where querySequence = unSD (seqdata s) querySequenceUreplacedwithT = L8.map bsreplaceUT querySequence - md5Sequence = M.md5 querySequenceUreplacedwithT - qString = "?md5=" ++ (show md5Sequence) + querySequenceU2Twolb = L8.filter (\a -> a /= '\n') querySequenceUreplacedwithT + md5Sequence = M.md5 querySequenceU2Twolb + qString = "?md5=" ++ show md5Sequence +--Build a query from a input string +buildStringViaMD5Query :: String -> String +buildStringViaMD5Query s = qString + where querySequenceUreplacedwithT = L8.map bsreplaceUT (L8.pack s) + querySequenceU2Twolb = L8.filter (\a -> a /= '\n') querySequenceUreplacedwithT + md5Sequence = M.md5 querySequenceU2Twolb + qString = "?md5=" ++ show md5Sequence + showRNAcentralAlienEvaluation :: [(Either String RNAcentralEntryResponse)] -> String showRNAcentralAlienEvaluation responses = output where resultEntries = Prelude.concatMap results (rights responses) resulthead = "rnacentral_id\tmd5\tlength\n" resultentries = Prelude.concatMap showRNAcentralAlienEvaluationLine resultEntries - output = if resultentries == [] then resulthead ++ "No matching sequences found in RNAcentral\n" else resulthead ++ resultentries + output = if resultentries == [] then "No matching sequences found in RNAcentral\n" else resulthead ++ resultentries showRNAcentralAlienEvaluationLine :: RNAcentralEntry -> String showRNAcentralAlienEvaluationLine entry = unpack (rnacentral_id entry) ++ "\t" ++ unpack (md5 entry) ++ "\t" ++ show (Bio.RNAcentralHTTP.length entry) ++"\n" @@ -109,3 +125,4 @@ bsreplaceUT a | a == 'U' = 'T' | otherwise = a + diff --git a/src/Bio/RNAcentralHTTPRequest.hs b/src/Bio/RNAcentralHTTPRequest.hs new file mode 100644 --- /dev/null +++ b/src/Bio/RNAcentralHTTPRequest.hs @@ -0,0 +1,28 @@ +{-# LANGUAGE RecordWildCards #-} +{-# LANGUAGE DeriveDataTypeable #-} + +-- | RNAcentralHTTPRequest +-- Testcommand: dist/build/RNAcentralHTTPRequest/RNAcentralHTTPRequest -i ATACTTACCTGGCACAGGGGATACCACGATCACCAAGGTGGTTCCCCCAAGACGAGGCTCACCATTGCACTCCGGTGGCGCTGACCCTTGCAATGACCCCAAATGTGGGTTACTCGGGTGTGTAATTTCTGTTAGCTGGGGACTGCGTTCGCGCTTTCCCCTT +module Main where + +import System.Console.CmdArgs +import Bio.RNAlienLibrary +import Data.Either.Unwrap +import Bio.RNAcentralHTTP + +data Options = Options + { inputSequence :: String + } deriving (Show,Data,Typeable) + +options :: Options +options = Options + { inputSequence = def &= name "i" &= help "input sequence" + } &= summary "RNAcentralHTTPRequest" &= help "Florian Eggenhofer 2016" &= verbosity + +main :: IO () +main = do + Options{..} <- cmdArgs options + let query = buildStringViaMD5Query inputSequence + rnacentralentries <- getRNACentralEntries [query] + print rnacentralentries + diff --git a/src/Bio/RNAlien.hs b/src/Bio/RNAlien.hs --- a/src/Bio/RNAlien.hs +++ b/src/Bio/RNAlien.hs @@ -15,6 +15,8 @@ import Data.Either.Unwrap import Data.Time import qualified System.FilePath as FP +import Paths_RNAlien (version) +import Data.Version (showVersion) data Options = Options { inputFastaFilePath :: String, @@ -51,7 +53,7 @@ threads = 1 &= name "c" &= help "Number of available cpu slots/cores. Default: 1", taxonomyRestriction = Nothing &= name "r" &= help "Restrict search space to taxonomic kingdom (bacteria,archea,eukaryia). Default: not set", sessionIdentificator = Nothing &= name "d" &= help "Optional session id that is used instead of automatically generated one." - } &= summary "RNAlien version 1.1.2" &= help "Florian Eggenhofer, Ivo L. Hofacker, Christian Höner zu Siederdissen - 2013 - 2016" &= verbosity + } &= summary ("RNAlien " ++ alienVersion) &= help "Florian Eggenhofer, Ivo L. Hofacker, Christian Höner zu Siederdissen - 2013 - 2016" &= verbosity main :: IO () main = do @@ -64,7 +66,7 @@ createDirectoryIfMissing False temporaryDirectoryPath createDirectoryIfMissing False (temporaryDirectoryPath ++ "log") -- Create Log files - writeFile (temporaryDirectoryPath ++ "Log") ("RNAlien 1.1.2" ++ "\n") + writeFile (temporaryDirectoryPath ++ "Log") ("RNAlien " ++ alienVersion ++ "\n") writeFile (temporaryDirectoryPath ++ "log/warnings") ("") logMessage ("Timestamp: " ++ (show timestamp) ++ "\n") temporaryDirectoryPath logMessage ("Temporary Directory: " ++ temporaryDirectoryPath ++ "\n") temporaryDirectoryPath @@ -103,3 +105,6 @@ writeFile (temporaryDirectoryPath ++ "result.csv") resultTaxonomyRecordsCSVTable resultSummary modelConstructionResults staticOptions writeFile (temporaryDirectoryPath ++ "done") "" + +alienVersion :: String +alienVersion = showVersion version diff --git a/src/Bio/RNAlienLibrary.hs b/src/Bio/RNAlienLibrary.hs --- a/src/Bio/RNAlienLibrary.hs +++ b/src/Bio/RNAlienLibrary.hs @@ -1507,7 +1507,7 @@ let clustalFilepath = evaluationDirectoryFilepath ++ "result.clustal" let reformatedClustalPath = evaluationDirectoryFilepath ++ "result.clustal.reformated" let cmFilepath = tempDirPath staticOptions ++ "result.cm" - let resultSequences = map nucleotideSequence (concatMap sequenceRecords (taxRecords mCResult)) + let resultSequences = (inputFasta mCResult):map nucleotideSequence (concatMap sequenceRecords (taxRecords mCResult)) let resultNumber = length resultSequences + 1 let rnaCentralQueries = map buildSequenceViaMD5Query resultSequences rnaCentralEntries <- getRNACentralEntries rnaCentralQueries @@ -1563,7 +1563,7 @@ rnaCodeEntries = concatMap showRNACodeHits (RC.rnacodeHits rnaCode) showRNACodeHits :: RC.RNAcodeHit -> String -showRNACodeHits rnacodeHit = show (RC.hss rnacodeHit) ++ "\t" ++ show (RC.frame rnacodeHit) ++ "\t" ++ show (RC.length rnacodeHit) ++ "\t"++ show (RC.from rnacodeHit) ++ "\t" ++ show (RC.to rnacodeHit) ++ "\t" ++ (RC.name rnacodeHit) ++ "\t" ++ show (RC.start rnacodeHit) ++ "\t" ++ show (RC.end rnacodeHit) ++ "\t" ++ show (RC.score rnacodeHit) ++ show (RC.pvalue rnacodeHit) ++ "\n" +showRNACodeHits rnacodeHit = show (RC.hss rnacodeHit) ++ "\t" ++ show (RC.frame rnacodeHit) ++ "\t" ++ show (RC.hitLength rnacodeHit) ++ "\t"++ show (RC.from rnacodeHit) ++ "\t" ++ show (RC.to rnacodeHit) ++ "\t" ++ (RC.name rnacodeHit) ++ "\t" ++ show (RC.start rnacodeHit) ++ "\t" ++ show (RC.end rnacodeHit) ++ "\t" ++ show (RC.score rnacodeHit) ++ show (RC.pvalue rnacodeHit) ++ "\n" -- | Call for external preprocessClustalForRNAz preprocessClustalForRNAzExternal :: String -> String -> IO (Either String String) diff --git a/src/Bio/cmsearchToBED.hs b/src/Bio/cmsearchToBED.hs new file mode 100644 --- /dev/null +++ b/src/Bio/cmsearchToBED.hs @@ -0,0 +1,178 @@ +{-# LANGUAGE RecordWildCards #-} +{-# LANGUAGE DeriveDataTypeable #-} + +-- | Convert cmsearch output to Browser Extensible Data (BED) format +-- Testcommand: cmsearchToBED -i /path/to/test.clustal +module Main where +import Prelude +import System.Console.CmdArgs +import Bio.RNAlienLibrary +import Data.Either.Unwrap +import qualified Data.ByteString.Lazy.Char8 as L +import qualified Data.Text as T +import Data.List + +data Bed = Bed + { browserPostition :: T.Text, + browserSettings :: T.Text, + bedName :: T.Text, + bedDescription :: T.Text, + bedVisibility :: Int, + bedItemRgb :: Bool, + bedEntries :: [BedEntry] + } deriving (Eq, Read) + +instance Show Bed where + show (Bed _browserPostition _browserSettings _bedName _bedDescription _bedVisibility _bedItemRgb _bedEntries) = a ++ b ++ c ++ d ++ e ++ f ++ g + where a = "browser position " ++ (T.unpack _browserPostition) ++ "\n" + b = (T.unpack _browserSettings) ++ "\n" + c = "track name=\"" ++ T.unpack _bedName ++ "\" " + d = "description=\"" ++ T.unpack _bedDescription ++ "\" " + e = "visibility=" ++ show _bedVisibility ++ " " + f = "itemRgb=\"" ++ itemRbg ++ "\"\n" + itemRbg = if _bedItemRgb then "On" else "Off" + g = concatMap show _bedEntries + + +data BedEntry = BedEntry + { chrom :: T.Text, + chromStart :: Int, + chromEnd :: Int, + chromName :: Maybe T.Text, + score :: Maybe Int, + strand :: Maybe Char, + thickStart :: Maybe Int, + thickEnd :: Maybe Int, + color :: Maybe T.Text, + blockCount :: Maybe Int, + blockSizes :: Maybe [Int], + blockStarts :: Maybe [Int] + } deriving (Eq, Read) + +instance Show BedEntry where + show (BedEntry _chrom _chromStart _chromEnd _chromName _score _strand _thickStart _thickEnd _color _blockCount _blockSizes _blockStarts) = a ++ b ++ c ++ d ++ e ++ f ++ g ++ h ++ i ++ j ++ k ++ l + where a = T.unpack _chrom ++ "\t" + b = show _chromStart ++ "\t" + c = show _chromEnd ++ "\t" + d = maybe "" T.unpack _chromName ++ "\t" + e = maybe "" show _score ++ "\t" + f = maybe "" (\s -> [s]) _strand ++ "\t" + g = maybe "" show _thickStart ++ "\t" + h = maybe "" show _thickEnd ++ "\t" + i = maybe "" T.unpack _color ++ "\t" + j = maybe "" show _blockCount ++ "\t" + k = maybe "" (\sizes -> intercalate "," (map show sizes)) _blockSizes ++ "\t" + l = maybe "" (\starts -> intercalate "," (map show starts)) _blockStarts ++ "\n" + +data Options = Options + { cmsearchPath :: String, + inputBrowserSettings :: String, + inputBedVisibility :: Int, + inputTrackName :: String, + inputTrackDescription :: String, + inputItemRgb :: Bool, + inputTrackColor :: String, + sortBed :: Bool, + withHeader :: Bool + } deriving (Show,Data,Typeable) + +options :: Options +options = Options + { cmsearchPath = def &= name "i" &= help "Path to input cmsearch file", + inputBrowserSettings = "browser hide all" &= name "b" &= help "Browser settings. Default: browser hide all", + inputBedVisibility = (2 :: Int) &= name "y" &= help "Visibility setting of track. Default: 2", + inputTrackName = "PredictedRNA" &= name "n" &= help "Name of the track Default: PredictedRNA", + inputTrackDescription = "RNA loci predicted by cmsearch" &= name "d" &= help "Description of the track. Default: RNA loci predicted by cmsearch", + inputItemRgb = True &= name "r" &= help "RGB Color of the track. Default: True", + inputTrackColor = "255,0,0" &= name "c" &= help "RGB Color of the track. Default: 255,0,0", + sortBed = True &= name "s" &= help "Sort entries of Bed file by start end end cooridinates. Default: True", + withHeader = True &= name "w" &= help "Output contains bed header. Default: True" + } &= summary "cmsearchToBED - Converts cmsearch file hits to BED file entries" &= help "Florian Eggenhofer 2016" &= verbosity + +main :: IO () +main = do + Options{..} <- cmdArgs options + parsedCmsearch <- readCMSearch cmsearchPath + if (isRight parsedCmsearch) + then do + let outputBED = convertcmSearchToBED (fromRight parsedCmsearch) inputBrowserSettings inputTrackName inputTrackDescription inputTrackColor inputBedVisibility inputItemRgb sortBed + if (isRight outputBED) + then + if withHeader + then print (fromRight outputBED) + else do + let output = concatMap show (bedEntries (fromRight outputBED)) + putStr output + else putStr (fromLeft outputBED) + else (putStr ("A problem occured converting from cmsearch to BED format:\n " ++ show (fromLeft parsedCmsearch))) + +--convertcmSearchToBED :: CMsearch -> String -> String -> Either String String +--convertcmSearchToBED inputcmsearch trackName trackColor +-- | null cmHits = Left "cmsearch file contains no hits" +-- | otherwise = Right (bedHeader ++ bedEntries) +-- where cmHits = cmsearchHits inputcmsearch +-- bedHeader = "browser position " ++ browserPosition ++ "\nbrowser hide all\ntrack name=\"cmsearch hits\" description=\"cmsearch hits\" visibility=2 itemRgb=\"On\"\n" +-- bedEntries = concatMap (cmsearchHitToBEDentry trackName trackColor) cmHits +-- browserPosition = L.unpack (hitSequenceHeader firstHit) ++ ":" ++ entryStart firstHit ++ "-" ++ entryEnd firstHit +-- firstHit = (head cmHits) + +convertcmSearchToBED :: CMsearch -> String -> String -> String -> String -> Int -> Bool -> Bool -> Either String Bed +convertcmSearchToBED inputcmsearch inputBrowserSettings trackName trackDescription trackColor inputBedVisibility inputItemRgb sortBed + | null cmHits = Left "cmsearch file contains no hits" + | otherwise = Right bed + where cmHits = cmsearchHits inputcmsearch + --bedHeader = "browser position " ++ browserPosition ++ "\nbrowser hide all\ntrack name=\"cmsearch hits\" description=\"cmsearch hits\" visibility=2 itemRgb=\"On\"\n" + bedEntries = map (cmsearchHitToBEDentry trackName trackColor) cmHits + sortedBedEntries = if sortBed then sortBy orderBedEntry bedEntries else bedEntries + currentBrowserPosition = T.unpack (chrom firstEntry) ++ ":" ++ show (chromStart firstEntry) ++ "-" ++ show (chromEnd firstEntry) + firstEntry = (head sortedBedEntries) + bed = Bed (T.pack currentBrowserPosition) (T.pack inputBrowserSettings) (T.pack trackName) (T.pack trackDescription) inputBedVisibility inputItemRgb sortedBedEntries + +cmsearchHitToBEDentry :: String -> String -> CMsearchHit -> BedEntry +cmsearchHitToBEDentry hitName hitColor cmHit = entry + where entry = BedEntry chromosome entrystart entryend (Just (T.pack hitName)) entryscore entrystrand thickstart thickend entrycolor blocks blockSize blockStart + chromosome = T.pack (L.unpack (hitSequenceHeader cmHit)) + --entryline = L.unpack (hitSequenceHeader cmHit) ++ "\t" ++ entryStart cmHit ++ "\t" ++ entryEnd cmHit++ "\t" ++ (hitName) ++ "\t" ++ "0" ++ "\t" ++ [(hitStrand cmHit)] ++ "\t" ++ show (hitStart cmHit) ++ "\t" ++ show (hitEnd cmHit) ++ "\t" ++ hitColor ++ "\n" + entrystart = if hitStrand cmHit == '+' then hitStart cmHit else hitEnd cmHit + entryend = if hitStrand cmHit == '+' then hitEnd cmHit else hitStart cmHit + entryscore = Just (0 :: Int) + entrystrand = Just (hitStrand cmHit) + thickstart = Just entrystart + thickend = Just entryend + entrycolor = Just (T.pack hitColor) + blocks = Just (1 :: Int) + blockSize = Just [(entryend - entrystart)] + blockStart = Just [(0 :: Int)] + + +--cmsearchHitToBEDentry :: String -> String -> CMsearchHit -> String +--cmsearchHitToBEDentry hitName hitColor cmHit = entryline +-- where entryline = L.unpack (hitSequenceHeader cmHit) ++ "\t" ++ entryStart cmHit ++ "\t" ++ entryEnd cmHit++ "\t" ++ (hitName) ++ "\t" ++ "0" ++ "\t" ++ [(hitStrand cmHit)] ++ "\t" ++ show (hitStart cmHit) ++ "\t" ++ show (hitEnd cmHit) ++ "\t" ++ hitColor ++ "\n" + --entrystart = if (hitStrand cmHit) == '+' then show (hitStart cmHit) else show (hitEnd cmHit) + --entryend = if (hitStrand cmHit) == '+' then show (hitEnd cmHit) else show (hitStart cmHit) + +entryStart :: CMsearchHit -> String +entryStart cmHit + | (hitStrand cmHit) == '+' = show (hitStart cmHit) + | otherwise = show (hitEnd cmHit) + +entryEnd :: CMsearchHit -> String +entryEnd cmHit + | (hitStrand cmHit) == '+' = show (hitEnd cmHit) + | otherwise = show (hitStart cmHit) + +orderBedEntry :: BedEntry -> BedEntry -> Ordering +orderBedEntry firstHit secondHit + | start1 > start2 = GT + | start1 < start2 = LT + | otherwise = orderBedEntryEnd firstHit secondHit + where start1 = chromStart firstHit + start2 = chromStart secondHit + +orderBedEntryEnd :: BedEntry -> BedEntry -> Ordering +orderBedEntryEnd firstHit secondHit + | end1 > end2 = GT + | end1 < end2 = LT + | otherwise = EQ + where end1 = chromEnd firstHit + end2 = chromEnd secondHit