<commit_before>#!/usr/bin/env python

# python libs
import sys, time, os, random, math

# rasmus libs
from rasmus import env, treelib, alignlib
from rasmus import fasta, util
from rasmus import sindirlib, stats

# scipy libs
from scipy import matrixmultiply as mm
from scipy import transpose, array, dot
from scipy.linalg import expm, logm, eig



options = [
"""\
    Simulates sequence evolution using SINDIR's model. 
""",

    ["s:", "stree=", "stree", "<species tree>",
        {"single": True}],
    ["S:", "smap=", "smap", "<gene2species map>"],
    ["p:", "params=", "params", "<params>",
        {"single": True}],
    ["r:", "baserate=", "baserate", "<baserate>",
        {"single": True,
         "default": None,
         "parser": float}],
        
    "Sequence generation",
    ["n:", "ntrees=", "ntrees", "<number of trees to produce>",
        {"single": True,
         "default": 1,
         "parser": int}],    
    ["l:", "seqlen=", "seqlen", "<length of starting sequence>",
        {"single": True,
         "default": 0,
         "parser": int}],
    ["q:", "seqfile=", "seqfile", "<FASTA file of starting sequences>",
        {"default": []}],
    ["G:", "gc=", "gc", "<GC content>",
        {"single": True,
         "default": .5,
         "parser": float}],
    
    "Duplication/Loss configuration",
    ["D:", "dupprob=", "dupprob", "<duplication probability>",
        {"single": True,
         "default": 0,
         "parser": float}],
    ["P:", "predup=", "predup", "<number of duplication before speciations>",
        {"single": True,
         "default": 0,
         "parser": int}],
    ["L:", "lossprob=", "lossprob", "<loss probability>",
        {"single": True,
         "default": 0,
         "parser": float}],
    
    "Kimura model parameters",
    ["a:", "alpha=", "alpha", "<kimura alpha>",
        {"single": True,
         "default": .25,
         "parser": float}],
    ["b:", "beta=", "beta", "<kimura beta>",
        {"single": True,
         "default": .25,
         "parser": float}],
    
    "Codon simulation",
    ["", "csm=", "csm", "<Codon Substitution Matrix>",
        {"single": True}],
    ["", "csmbranchlen=", "csmbranchlen", "<branch length for species in CSM>",
        {"single": True,
         "parser": float}],
    
    "Output extensions",
    ["T:", "outtree=", "outtree", "<output tree filename extension>",
        {"single": True,
         "default": ".tree"}],
    ["A:", "outalign=", "outalign", "<output align filename extension>",
        {"single": True,
         "default": ".align"}],
#    ["d:", "outdist=", "outdist", "<output distance filename extension>",
#        {"single": True,
#         "default": ".dist"}],
    
]


conf = util.parseOptions(sys.argv, options, quit=True)


BASE2INT = {
    "A": 0,
    "C": 1,
    "G": 2,
    "T": 3
}

INT2BASE = ["A", "C", "G", "T"]

KIMURA_MATRIX = [
    ['r', 's', 'u', 's'],
    ['s', 'r', 's', 'u'],
    ['u', 's', 'r', 's'],
    ['s', 'u', 's', 'r']
]


def main(conf):
    env.addEnvPaths("DATAPATH")
    
    # parse options
    if conf["dupprob"] >= 1:
        conf["dupnum"] = int(conf["dupprob"])
    else:
        conf["dupnum"] = 0
        
    if conf["lossprob"] >= 1:
        conf["lossnum"] = int(conf["lossprob"])
    else:
        conf["lossnum"] = 0
    
    
    
    # adjustment for baserate (empirically determined)
    #conf["baserate"] = conf["baserate"] / .9
    
    
    stree = treelib.readTree(env.findFile(conf["stree"]))
    params = sindirlib.readParams(conf["params"])
    
    # read species2gene
    for f in conf["smap"]:
        smap = util.readDelim(env.findFile(f))
    species2gene = {}
    
    for prefix, genome in smap:
        species2gene[genome] = prefix.replace("*", "")
    
    # read starting sequences
    rootseqs = []
    for f in conf["seqfile"]:
        rootseqs.extend(fasta.readFasta(f).values())
    rootseqs = [s.replace("-", "") for s in rootseqs]
    
    
    # load CSM
    if "csm" in conf:
        labels, conf["csm"], conf["codon_usage"] = loadCsm(conf["csm"])
        
        # remove AA from labels
        labels = [x[:x.index("/")] for x in labels]
        conf["csmlabels"] = labels
        
        # make csm lookup
        conf["csmlookup"] = {}
        for i, label in enumerate(labels):
            conf["csmlookup"][label] = i
        
        # find rate matrix
        conf["csmrate"] = logm(conf["csm"]) / conf["csmbranchlen"]
        
    
    
    # simulate
    util.tic("simulating %d sequences" % conf["ntrees"])
    for i in range(conf["ntrees"]):
        util.log("simulating", i)
        
        # create initial sequence
        if conf["seqlen"] > 0:
            rootseq = createRandomSeq(conf["seqlen"], 
                                      conf["csmlabels"],<commit_msg>*** empty log message ***<commit_after>#!/usr/bin/env python

# python libs
import sys, time, os, random, math

# rasmus libs
from rasmus import env, treelib, alignlib
from rasmus import fasta, util
from rasmus import sindirlib, stats

# scipy libs
from scipy import matrixmultiply as mm
from scipy import transpose, array, dot
from scipy.linalg import expm, logm, eig



options = [
"""\
    Simulates sequence evolution using SINDIR's model. 
""",

    ["s:", "stree=", "stree", "<species tree>",
        {"single": True}],
    ["S:", "smap=", "smap", "<gene2species map>"],
    ["p:", "params=", "params", "<params>",
        {"single": True}],
    ["r:", "baserate=", "baserate", "<baserate>",
        {"single": True,
         "default": None,
         "parser": float}],
        
    "Sequence generation",
    ["n:", "ntrees=", "ntrees", "<number of trees to produce>",
        {"single": True,
         "default": 1,
         "parser": int}], 
    ["", "start=", "start", "<starting number>",
        {"single": True,
         "default": 0,
         "parser": int}],
    ["l:", "seqlen=", "seqlen", "<length of starting sequence>",
        {"single": True,
         "default": 0,
         "parser": int}],
    ["q:", "seqfile=", "seqfile", "<FASTA file of starting sequences>",
        {"default": []}],
    ["G:", "gc=", "gc", "<GC content>",
        {"single": True,
         "default": .5,
         "parser": float}],
    
    "Duplication/Loss configuration",
    ["D:", "dupprob=", "dupprob", "<duplication probability>",
        {"single": True,
         "default": 0,
         "parser": float}],
    ["P:", "predup=", "predup", "<number of duplication before speciations>",
        {"single": True,
         "default": 0,
         "parser": int}],
    ["L:", "lossprob=", "lossprob", "<loss probability>",
        {"single": True,
         "default": 0,
         "parser": float}],
    
    "Kimura model parameters",
    ["a:", "alpha=", "alpha", "<kimura alpha>",
        {"single": True,
         "default": .25,
         "parser": float}],
    ["b:", "beta=", "beta", "<kimura beta>",
        {"single": True,
         "default": .25,
         "parser": float}],
    
    "Codon simulation",
    ["", "csm=", "csm", "<Codon Substitution Matrix>",
        {"single": True}],
    ["", "csmbranchlen=", "csmbranchlen", "<branch length for species in CSM>",
        {"single": True,
         "parser": float}],
    
    "Output extensions",
    ["T:", "outtree=", "outtree", "<output tree filename extension>",
        {"single": True,
         "default": ".tree"}],
    ["A:", "outalign=", "outalign", "<output align filename extension>",
        {"single": True,
         "default": ".align"}],
#    ["d:", "outdist=", "outdist", "<output distance filename extension>",
#        {"single": True,
#         "default": ".dist"}],
    
]


conf = util.parseOptions(sys.argv, options, quit=True)


BASE2INT = {
    "A": 0,
    "C": 1,
    "G": 2,
    "T": 3
}

INT2BASE = ["A", "C", "G", "T"]

KIMURA_MATRIX = [
    ['r', 's', 'u', 's'],
    ['s', 'r', 's', 'u'],
    ['u', 's', 'r', 's'],
    ['s', 'u', 's', 'r']
]


def main(conf):
    env.addEnvPaths("DATAPATH")
    
    # parse options
    if conf["dupprob"] >= 1:
        conf["dupnum"] = int(conf["dupprob"])
    else:
        conf["dupnum"] = 0
        
    if conf["lossprob"] >= 1:
        conf["lossnum"] = int(conf["lossprob"])
    else:
        conf["lossnum"] = 0
    
    
    
    # adjustment for baserate (empirically determined)
    #conf["baserate"] = conf["baserate"] / .9
    
    
    stree = treelib.readTree(env.findFile(conf["stree"]))
    params = sindirlib.readParams(conf["params"])
    
    # read species2gene
    for f in conf["smap"]:
        smap = util.readDelim(env.findFile(f))
    species2gene = {}
    
    for prefix, genome in smap:
        species2gene[genome] = prefix.replace("*", "")
    
    # read starting sequences
    rootseqs = []
    for f in conf["seqfile"]:
        rootseqs.extend(fasta.readFasta(f).values())
    rootseqs = [s.replace("-", "") for s in rootseqs]
    
    
    # load CSM
    if "csm" in conf:
        labels, conf["csm"], conf["codon_usage"] = loadCsm(conf["csm"])
        
        # remove AA from labels
        labels = [x[:x.index("/")] for x in labels]
        conf["csmlabels"] = labels
        
        # make csm lookup
        conf["csmlookup"] = {}
        for i, label in enumerate(labels):
            conf["csmlookup"][label] = i
        
        # find rate matrix
        conf["csmrate"] = logm(conf["csm"]) / conf["csmbranchlen"]
        
    
    
    # simulate
    util.tic("simulating %d sequences" % conf["ntrees"])
    for i in range(conf["start"], conf["ntrees"]):
        util.log("simulating", i)
        
        # create initial sequence
        if conf["seqlen"] > 0:
            rootseq = createRandomSeq(conf["seqlen"], 
                                      conf["csmlabels"],
{'fraction_non_alphanumeric': 0.19936351503952365, 'fraction_numerical': 0.005440919823426753, 'mean_word_length': 2.324914675767918, 'pattern_counts': {'":': 89, '<': 40, '<?xml version=': 0, '>': 46, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2709466', 'n_tokens_mistral': 3298, 'n_tokens_neox': 3196, 'n_words': 832}