MutationOrder 0.0.0.1 → 0.0.0.2
raw patch · 4 files changed
+112/−29 lines, 4 filesdep +file-embedPVP ok
version bump matches the API change (PVP)
Dependencies added: file-embed
API changes (from Hackage documentation)
Files
- MutationOrder.cabal +3/−1
- README.md +74/−12
- changelog.md +6/−0
- src/MutationOrder.hs +29/−16
MutationOrder.cabal view
@@ -1,5 +1,5 @@ name: MutationOrder-version: 0.0.0.1+version: 0.0.0.2 author: Maria Beatriz Walter Costa, Christian Hoener zu Siederdissen, 2017 copyright: Christian Hoener zu Siederdissen, 2017 homepage: https://github.com/choener/MutationOrder@@ -111,6 +111,7 @@ build-depends: base , bytestring , cmdargs >= 0.10+ , file-embed >= 0.0.8 , filepath -- , MutationOrder@@ -119,6 +120,7 @@ default-extensions: BangPatterns , DeriveDataTypeable , RecordWildCards+ , TemplateHaskell main-is: MutationOrder.hs default-language:
README.md view
@@ -2,32 +2,94 @@ Determine the most likely order of mutations from one RNA sequence to another. -1. Walter Costa, Maria Beatriz and Hoener zu Siederdissen, Christian and Tulpan, Dan and Stadler, Peter F. and Nowick, Katja + Walter Costa, Maria Beatriz and Hoener zu Siederdissen, Christian and Tulpan, Dan and Stadler, Peter F. and Nowick, Katja *Uncovering the Structural Evolution of the Human Accelerated Region 1* 2017, submitted [preprint](http://www.bioinf.uni-leipzig.de/~choener/pdfs/wal-hoe-2017.pdf) +# General information++Given two RNA sequences, one ancestral, and one extant, we want to determine+the most likely path of evolution under different measures of fitness.++This program produces the (i) maximum-likelihood path, (ii) all end+probabilities, (iii) all start-end probabilities, (iv) all edge probabilities,+and (v) the maximum expected accuracy path for these two RNA sequences.++In detail: +(i) gives the optimal path(s) for the fitness function +(ii) gives for each nucleotide polymorphism, how likely it is, that this mutation was introduced *last* +(iii) looks at all pairs of (first mutation, last mutation) and gives the probability that these two mutations are the begin and end of the chain of mutations +(iv) yields for all pairs of nodes (i -> j) the probability that this path occurs, over the whole ensemble of all possible paths +(v) produces the path of maximal weight using the probabilities produced in (iv) + # Usage instructions -We assume that you have two Fasta files, *from.fa* and *to.fa* but they can be-named however is convenient (say *chimp.fa* and *human.fa*). Each file has to-contain exactly one sequence and both sequences have to be of the same length.+## example usage +We assume that you have two Fasta files, *chimp_118.fa* and *human_118.fa* but+they can be named however is convenient. Each file has to contain exactly one+sequence and both sequences have to be of the same length.++For testing with chimp and human, the provided chimp-human.json.gz database+should be used, otherwise the initial foldings will be recalculated. All+required files are available under 'Binaries' at the bottom of the page.++In case, you don't want or can't use the provided work database, run+./MutationOrder with --verbose+ We then run -```./MutationOrder --workdb from-to.json.gz --scoretype pairdistcen --onlypositive --outputprefix test```+ ./MutationOrder --workdb chimp-human.json.gz --scoretype pairdistcen --onlypositive --outputprefix test chimp_118.fa human_118.fa This will generate ```test.run```, ```test-edge.eps```, and-```test-meaorder.eps```. The ```test.run``` file provides extensive output of-the optimal path, the first-last probabilities, the edge probabilities, and the-mea output. The two ```eps``` files give a graphical representation of the edge+```test-meaorder.eps```.++The ```test.run``` file provides extensive output of the optimal path, the+first-last probabilities, the edge probabilities, and the mea output. This+conforms to (i) -- (v) mentioned above.++The two ```eps``` files give a graphical representation of the edge probabilities, for the ```meaorder``` in order of the path of maximum expected accuracy. -The ```--scoretype``` allows for ```mfe```, ```centroid```, ```pairdistcen```,-and ```pairdistmfe```, which analyse possible evoluationary paths according to-mfe energy, centroid energy, smallest base pair distances for each step in the-```cen```troid or ```mfe``` case.+The work database collects intermediate structures and their folding and is+only created once. The initial run will, however, take some time. I.e. for+'HAR1' this requires 1-4 hours depending on the machine. Further runs complete+*much* faster. In minutes for HAR1.++## Command-line options++ --help provides short help+ --verbose will show folding steps during the initial run++ -w+ --workdb=ITEM the database where to store intermediate foldings+ -t+ --temperature=NUM annealing temperature. Values close to 0 favor optimal paths. The default is 1.0+ --fillweight=FILLWEIGHT provides logarithmic and linear fill styles for probability plots. The full style always fills the box+ --fillstyle=FILLSTYLE normally, boxes are sized, but all in the same color. This changes the opacity of the color as well. Does not work well for eps files+ --cooptcount=INT how many co-optimals to count for (the count in the .run file is produced differently)+ --cooptprint=INT how many co-optimals to actually print+ --outprefix=ITEM how to prefix all output files+ --scoretype=SCORETYPE choose 'mfe', 'centroid', 'pairdistmfe', or 'pairdistcen' for the evaluation of each mutational step+ --positivesquared square positive energies to penalize bad moves+ --onlypositive minimize only over penalties, not energy gains+ --equalstart each possible mutation is selected with equal probability as the initial one+ --posscaled=NUM,NUM in =x,y will exponentiate all numbers >=x by the constant y. For value k>=x, we have k^y+ --lkupfile=ITEM developer option to feed the initial work database with known foldings (usable but very raw and undocumented. needs 5-line rnafold output)+ --showmanual will show this manual++The allowed score types are: ++ mfe+which optimizes based on the minimum free energy of each intermediate sequence+ centroid+which instead looks at the energy of the centroid structure+ pairdistmfe+which minimizes the base pair distance between following mutations using mfe structures+ pairdistcen+which minimizes the base pair distance between following mutations using centroid structures
changelog.md view
@@ -1,3 +1,9 @@+0.0.0.2+-------++- updated readme+- included readme as manual+ 0.0.0.1 -------
src/MutationOrder.hs view
@@ -3,10 +3,12 @@ module Main where -import System.Console.CmdArgs-import System.FilePath+import Control.Monad+import Data.FileEmbed import qualified Data.ByteString.Char8 as BS-import Control.Monad+import System.Console.CmdArgs+import System.Exit (exitSuccess, exitFailure)+import System.FilePath import BioInf.MutationOrder import BioInf.MutationOrder.RNA (createRNAlandscape)@@ -19,7 +21,7 @@ deriving (Show,Data,Typeable) data Options- = Options+ = MutationOrder { infiles :: [FilePath] , workdb :: FilePath , temperature :: Double@@ -34,27 +36,29 @@ , equalStart :: Bool , posscaled :: Maybe (Double,Double) , lkupfile :: Maybe FilePath+ , showmanual :: Bool } | GenSequences { infiles :: [FilePath] } deriving (Show,Data,Typeable) -oOptions = Options+oMutationOrder = MutationOrder { infiles = def &= args , workdb = "work.db" &= help "name of the database to store intermediates in" , temperature = 1.0 &= help "lower temperatures favor the more optimal paths, defaults to 1.0"- , fillweight = FWlog- , fillstyle = FSfull- , cooptcount = 100000- , cooptprint = 2- , outprefix = "tmp"+ , fillweight = FWlog &= help "scale options for probability plots: fwlog, fwlinear, fwfill"+ , fillstyle = FSfull &= help "fill options for probability plots: fsopacitylog, fsopacitylinear, fsfull"+ , cooptcount = 1000 &= help "how many co-optimals to count"+ , cooptprint = 2 &= help "how many co-optimals to actually print out"+ , outprefix = "tmp" &= help "prefix for output files" , scoretype = Centroid &= help "choose 'mfe', 'centroid', 'pairdistmfe', or 'pairdistcen' for the evaluation of each mutational step" , positivesquared = False &= help "square positive energies to penalize worse structures" , onlypositive = False &= help "minimize only over penalties, not energy gains"- , equalStart = False- , posscaled = Nothing- , lkupfile = Nothing+ , equalStart = False &= help "run mea with equal start probabilities"+ , posscaled = Nothing &= help "--posscaled=x,y scale all values >= x by using y as exponent"+ , lkupfile = Nothing &= help "developer option: if an RNAfold file with foldings exists, then use it"+ , showmanual = False &= help "shows the manual" } oGenSequences = GenSequences@@ -63,9 +67,9 @@ main :: IO () main = do- o <- cmdArgs $ modes [oOptions, oGenSequences] &= verbosity+ o <- cmdArgs $ modes [oMutationOrder &= auto, oGenSequences] &= verbosity case o of- Options{} -> mainProgram o+ MutationOrder{} -> mainProgram o GenSequences{} -> genSequences o genSequences o = do@@ -76,8 +80,17 @@ forM_ ls $ \(k,sq) -> BS.putStrLn sq return () +embeddedManual = $(embedFile "README.md")+ mainProgram oOptions = do- let Options{..} = oOptions+ let MutationOrder{..} = oOptions+ when showmanual $ do+ BS.putStrLn embeddedManual+ exitSuccess+ unless (length infiles == 2) $ do+ BS.putStrLn embeddedManual+ putStrLn "\n\n\nThis program expects exactly two equal-length fasta files as input"+ exitFailure isL <- isLoud let fwdScaleFunction = (if positivesquared then squaredPositive else id)